📄 struct2d.h
字号:
}
for(i=0; i<F2ConElem_num; i++)
{
xyytemp1=(*(FTnode0+*(F2CElemNode_Info+i*4))+*(FTnode0+*(F2CElemNode_Info+i*4+1))+*(FTnode0+*(F2CElemNode_Info+i*4+2))+*(FTnode0+*(F2CElemNode_Info+i*4+3)))/4.0;
*(CElemtemp0+i)=xyytemp1;
}
R2Weiyi<<"恒载升温的恒载段的分析,荷载步:"<<Pdiv<<endl;
R2Sectwei<<"恒载升温的恒载段的分析,荷载步:"<<Pdiv<<endl;
R2Weiyi<<0<<","<<UAxial0<<","<<CUy0<<","<<CUz0<<endl;
R2Sectwei<<0<<curtime0<<","<<e0<<","<<wany<<","<<wanz<<endl;
WResult<<"恒温加载段分析,荷载步数为:"<<Pdiv<<endl;
cout<<"恒温加载段分析,荷载步数为:"<<Pdiv<<endl;
//****************************************************
//****************************************************
//**********************恒温加载段的分析******************************
//****************************************************
for(i=0; i<Pdiv; i++)
{
cout<<" 恒温加载段的第"<<i+1<<" 子步\n";
WResult<<" 恒温加载段的第"<<i+1<<" 子步\n";
//形成混凝土的应力切线模量矩阵
curtime0=0.0;
for(j=0; j<F2ConElem_num; j++)
{
//先求单元的平均温度
xyytemp1=(*(FTnode0+*(F2CElemNode_Info+j*4))+*(FTnode0+*(F2CElemNode_Info+j*4+1))+*(FTnode0+*(F2CElemNode_Info+j*4+2))+*(FTnode0+*(F2CElemNode_Info+j*4+3)))/4.0;
*(ConEstr+j*F2ConElem_num+j)=Ecst(Csresta,*(CElemstress+j),*(CElemstrain+j),xyytemp1,curtime0);
}
//形成钢筋的应力切线模量矩阵
for(j=0;j<F2Longbar_num;j++)
{
Wendubu=*(SJZ_Wendubu+j); //钢筋温度补偿时间(小时)
WendubuDT=*(SJZ_WendubuDT+j);
//先求单元的平均温度
xyytemp1=*(STnode0+j);
*(SteEstr+j*F2Longbar_num+j)=Esst(Ssresta,*(SElemstress+j),*(SElemstrain+j),xyytemp1,curtime0,*(F2LbarInfo+5*j+3),F2LbarInfo[5*j+3]/F2LbarInfo[5*j+4]);
}
//下面形成总刚度矩阵
for(j=0;j<F2ConElem_num;j++)//计算混凝土对总刚度矩阵的贡献
{
*(CEsmAC+j*F2ConElem_num+j)=*(ConEstr+j*F2ConElem_num+j)*(*(ConAc+j*F2ConElem_num+j));
//printf("%.5f\n",*(CEsmAC+j*F2ConElem_num+j));
}
//下面为两矩阵相乘;
for(j=0; j<3; j++)//先将混凝土单元形函数转置
for(k=0; k<F2ConElem_num; k++)
*(JuZhenC1+j*F2ConElem_num+k)=*(NCSect+k*3+j);
matmpl(JuZhenC1,3,F2ConElem_num,CEsmAC,F2ConElem_num,NCmES);
matmpl(NCmES,3,F2ConElem_num,NCSect,3,F2KZ);
//计算钢筋对总刚度矩阵的贡献
for(j=0; j<F2Longbar_num; j++)
SEsmAs[j*F2Longbar_num+j]=*(SteEstr+j*F2Longbar_num+j)*(*(SteAs+j*F2Longbar_num+j));
//下面为两矩阵相乘
for(j=0; j<3; j++)//先将钢筋单元形函数转置
for(k=0; k<F2Longbar_num; k++)
*(JuZhenS1+j*F2Longbar_num+k)=*(NSSect+k*3+j);
matmpl(JuZhenS1,3,F2Longbar_num,SEsmAs,F2Longbar_num,NSmES);
matmpl(NSmES,3,F2Longbar_num,NSSect,3,F2KZtemp);
for(j=0; j<3*3; j++)//总刚度矩阵形成完毕
{
*(F2KZ+j)=*(F2KZ+j)+*(F2KZtemp+j);
//printf("%.5f\n",*(F2KZ+j));
*(F2KZtemp+j)=*(F2KZ+j);
}
if(F2method==1) //按Lie的方法
{
CUy0=-(ComLength*ComLength*wanz/12.0+u0y); //上一迭代步结束时y方向的侧移
CUz0=-(ComLength*ComLength*wany/12.0+u0z); //上一迭代步结束时z方向的侧移
}
if(F2method==2) //按韩林海的方法
{
CUy0=-(ComLength*ComLength*wanz/pi/pi+u0y); //上一迭代步结束时y方向的侧移
CUz0=-(ComLength*ComLength*wany/pi/pi+u0z); //上一迭代步结束时z方向的侧移
}
if(PHxiaozheng==1) //采用平衡校正的迭代算法
{
//计算迭代步开始时截面总内力
//先求混凝土对内力的贡献
matmpl(ConAc,F2ConElem_num,F2ConElem_num,CElemstress,1,JuZhenC2);
//先将混凝土单元形函数转置;
for(j=0;j<3;j++)
for(k=0;k<F2ConElem_num;k++)
*(JuZhenC1+j*F2ConElem_num+k)=*(NCSect+k*3+j);
matmpl(JuZhenC1,3,F2ConElem_num,JuZhenC2,1,FcSect);
//printf("截面实际内力：%.3f %.3f %.3f\n",*(FcSect+0),*(FcSect+1),*(FcSect+2));
//计算钢筋对内力的贡献;
matmpl(SteAs,F2Longbar_num,F2Longbar_num,SElemstress,1,JuZhenS2);
//先将钢筋单元形函数转置
for(j=0;j<3;j++)
for(k=0;k<F2Longbar_num;k++)
*(JuZhenS1+j*F2Longbar_num+k)=*(NSSect+k*3+j);
matmpl(JuZhenS1,3,F2Longbar_num,JuZhenS2,1,JuZhens31);
//printf("钢筋截面实际内力：%.3f %.3f %.3f\n",*(JuZhens31+0),*(JuZhens31+1),*(JuZhens31+2));
for(j=0;j<3;j++)
*(FcSect+j)=*(FcSect+j)+*(JuZhens31+j);
//printf("截面实际内力：%d %.3f %.3f %.3f\n",i,*(FcSect+0),*(FcSect+1),*(FcSect+2));
//计算不平衡力
*(RFnbanlan+0)=i*DelP-*(FcSect+0);
*(RFnbanlan+1)=-i*DelP*CUy0-*(FcSect+1);
*(RFnbanlan+2)=-i*DelP*CUz0-*(FcSect+2);
}
//下面形成荷载增量向量
if(F2method==1) //按Lie的方法
{
CUy0=-(ComLength*ComLength*wanz/12.0+u0y);
CUz0=-(ComLength*ComLength*wany/12.0+u0z);
*(Zload+0)=DelP+*(RFnbanlan+0); //轴压荷载增量
*(Zload+1)=-DelP*(CUy0)+*(RFnbanlan+1);
*(Zload+2)=-DelP*(CUz0)+*(RFnbanlan+2);
//下面修正总刚度系数;
*(F2KZtemp+3*1+1)=*(F2KZtemp+3*1+1)-(i+1)*DelP*ComLength*ComLength/12.0;
*(F2KZtemp+3*2+2)=*(F2KZtemp+3*2+2)-(i+1)*DelP*ComLength*ComLength/12.0;
for(j=0;j<3;j++)
*(Zxewywz+j)=*(Zload+j);
//printf("未知量：%.8f %.8f %.8f\n",*(Zxewywz+0),*(Zxewywz+1),*(Zxewywz+2));
//下面求解位移增量;
if(agaus(F2KZtemp,Zxewywz,3)==0)
{
cout<<"按Lie简化方法求解时总刚度矩阵出现病态!\n";
WResult<<" ******按Lie简化方法求解时总刚度矩阵出现病态!******\n";
cout<<" ******按Lie简化方法求解结束!******\n";
WResult<<" ******按Lie简化方法求解结束!******\n";
cout<<" ******恒载升温的恒载段的第 "<<(i+1)<<" 分析子步求解结束******\n";
WResult<<" ******恒载升温的恒载段的第 "<<(i+1)<<" 分析子步求解结束******\n";
goto s102;
}
}
if(F2method==2) //按韩林海的方法
{
CUy0=-(ComLength*ComLength*wanz/pi/pi+u0y);
CUz0=-(ComLength*ComLength*wany/pi/pi+u0z);
*(Zload+0)=DelP+*(RFnbanlan+0); //轴压荷载增量
*(Zload+1)=-DelP*(CUy0)+*(RFnbanlan+1);
*(Zload+2)=-DelP*(CUz0)+*(RFnbanlan+2);
//下面修正总刚度系数;
*(F2KZtemp+3*1+1)=*(F2KZtemp+3*1+1)-(i+1)*DelP*ComLength*ComLength/pi/pi;
*(F2KZtemp+3*2+2)=*(F2KZtemp+3*2+2)-(i+1)*DelP*ComLength*ComLength/pi/pi;
for(j=0;j<3;j++)
*(Zxewywz+j)=*(Zload+j);
//下面求解位移增量;
if(agaus(F2KZtemp,Zxewywz,3)==0)
{
cout<<"按韩林海简化方法求解总刚度矩阵出现病态!\n";
WResult<<" ******按韩林海简化方法求解总刚度矩阵出现病态!******\n";
cout<<" ******按韩林海简化方法求解结束!******\n";
WResult<<" ******按韩林海简化方法求解结束!******\n";
cout<<" ******恒载升温的恒载段的第 "<<(i+1)<<" 分析子步求解结束******\n";
WResult<<" ******恒载升温的恒载段的第 "<<(i+1)<<" 分析子步求解结束******\n";
goto s102;
}
}
//截面未知位移更新
e0=e0+*(Zxewywz+0);
wanz=wanz+*(Zxewywz+1);
wany=wany+*(Zxewywz+2);
//printf("未知量：%.8f %.8f %.8f\n",e0,wany,wanz);
//构件位移更新
if(F2method==1) //按Lie方法
{
xewywz[0][0]=e0;
xewywz[1][0]=wanz;
xewywz[2][0]=wany;
CUy1=-(wanz*ComLength*ComLength/12.0+u0y);
CUz1=-(wany*ComLength*ComLength/12.0+u0z);
UAxial1=e0*CFirLength+(1/2.0+e0/2.0)*((CUy1+u0y)*(CUy1+u0y)+(CUz1+u0z)*(CUz1+u0z))*24.0/5.0/ComLength;
}
if(F2method==2) //按韩林海方法
{
xewywz[0][0]=e0;
xewywz[1][0]=wanz;
xewywz[2][0]=wany;
CUy1=-(wanz*ComLength*ComLength/pi/pi+u0y);
CUz1=-(wany*ComLength*ComLength/pi/pi+u0z);
UAxial1=e0*CFirLength+(1/2.0+e0/2.0)*((CUy1+u0y)*(CUy1+u0y)+(CUz1+u0z)*(CUz1+u0z))*pi*pi/2.0/ComLength;
}
for(j=0; j<3; j++)
{
*(Zxewywz0+j)=*(Zxewywz+j); //截面基本未知量增量备份
//printf("截面基本未知量：%.14f\n",*(Zxewywz0+j));
}
//下面形成总应变增量向量
matmpl(NCSect,F2ConElem_num,3,Zxewywz,1,ZCZstrain); //混凝土的总应变增量
matmpl(NSSect,F2Longbar_num,3,Zxewywz,1,ZSZstrain); //钢筋的总应变增量
//形成应力应变增量向量,这时候没有温度和时间应变的影响
for(j=0; j<F2ConElem_num; j++)
{
*(ZCStstrain+j)=*(ZCZstrain+j);
}
for(j=0; j<F2Longbar_num; j++)
{
*(ZSStstrain+j)=*(ZSZstrain+j);
}
//迭代步中的应变
for(j=0;j<F2ConElem_num;j++)
{
*(CElemstrain+j)=*(CElemstrain+j)+*(ZCStstrain+j);
}
for(j=0;j<F2Longbar_num;j++)
{
*(SElemstrain+j)=*(SElemstrain+j)+*(ZSStstrain+j);
}
//求混凝土应力向量
for(j=0; j<F2ConElem_num; j++)
{
//先求单元的平均温度
xyytemp1=(*(FTnode0+*(F2CElemNode_Info+j*4))+*(FTnode0+*(F2CElemNode_Info+j*4+1))+*(FTnode0+*(F2CElemNode_Info+j*4+2))+*(FTnode0+*(F2CElemNode_Info+j*4+3)))/4.0;
*(CElemstress+j)=CStress(Csresta,*(CElemstrain+j),xyytemp1);
}
//求钢筋应力向量
for(j=0;j<F2Longbar_num;j++)
{
xyytemp1=*(STnode0+j); //钢筋温度
*(SElemstress+j)=SStress(Ssresta,*(SElemstrain+j),xyytemp1,*(F2LbarInfo+5*j+3),*(F2LbarInfo+5*j+3)/(*(F2LbarInfo+5*j+4)));
}
matmpl(ConAc,F2ConElem_num,F2ConElem_num,CElemstress,1,JuZhenC2);
for(j=0;j<3;j++)//先将混凝土单元形函数转置
for(k=0;k<F2ConElem_num;k++)
*(JuZhenC1+j*F2ConElem_num+k)=*(NCSect+k*3+j);
matmpl(JuZhenC1,3,F2ConElem_num,JuZhenC2,1,FcSect);
//printf("截面实际内力：%.3f %.3f %.3f\n",*(FcSect+0),*(FcSect+1),*(FcSect+2));
//计算钢筋对内力的贡献
matmpl(SteAs,F2Longbar_num,F2Longbar_num,SElemstress,1,JuZhenS2);
//先将钢筋单元形函数转置
for(j=0;j<3;j++)
for(k=0;k<F2Longbar_num;k++)
*(JuZhenS1+j*F2Longbar_num+k)=*(NSSect+k*3+j);
matmpl(JuZhenS1,3,F2Longbar_num,JuZhenS2,1,JuZhens31);
for(j=0;j<3;j++)
*(FcSect+j)=*(FcSect+j)+*(JuZhens31+j);
//计算不平衡力
*(RFnbanlan+0)=(i+1)*DelP-*(FcSect+0);
*(RFnbanlan+1)=-(i+1)*DelP*CUy1-*(FcSect+1);
*(RFnbanlan+2)=-(i+1)*DelP*CUz1-*(FcSect+2);
liang2=0;
do
{
liang2=liang2+1;
for(j=0;j<3*3;j++)
*(F2KZtemp+j)=*(F2KZ+j);
//下面形成荷载增量向量
if(F2method==1) //按Lie的方法
{
*(Zload+0)=*(RFnbanlan+0); //轴压荷载增量
*(Zload+1)=*(RFnbanlan+1);
*(Zload+2)=*(RFnbanlan+2);
//下面修正总刚度系数
*(F2KZtemp+3*1+1)=*(F2KZtemp+3*1+1)-(i+1)*DelP*ComLength*ComLength/12.0;
*(F2KZtemp+3*2+2)=*(F2KZtemp+3*2+2)-(i+1)*DelP*ComLength*ComLength/12.0;
for(j=0;j<3;j++)
*(Zxewywz+j)=*(Zload+j);
//下面求解位移增量
if (agaus(F2KZtemp,Zxewywz,3)==0)
{
cout<<"按Lie简化方法求解时总刚度矩阵出现病态!"<<endl;
WResult<<" ******按Lie简化方法求解时总刚度矩阵出现病态!******"<<endl;
cout<<" ******按Lie简化方法求解结束!******\n";
WResult<<" ******按Lie简化方法求解结束!******\n";
cout<<" ******恒载升温的恒载段的第"<<(i+1)<<" 分析子步求解结束******\n";
WResult<<" ******恒载升温的恒载段的第"<<(i+1)<<" 分析子步求解结束******\n";
goto s102;
}
}
if(F2method==2) //按韩林海的方法
{
*(Zload+0)=*(RFnbanlan+0); //轴压荷载增量
*(Zload+1)=*(RFnbanlan+1);
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -