固支梁有限元毕业设计.doc

 PAGE \* MERGEFORMAT 68 题目 如图所示固支梁，高3米，长15m，承受均布载荷q=10kN/，E=20GPa，u=0.167，厚度t=1m，忽略自重，按平面应力问题分析。用有限元方法计算梁的变形及应力分布，要求用矩形单元。 q=10kN/ 要求： 1.单元数目不得少于20个。 2.采用矩形单元计算求解。 3.计算结果并给出变形图、应力分布图、单元划分图。 二、力学分析 1.题目可以看做是平面应力问题故有LXM=0 2.单元划分图 三．程序框图及程序 程序框图： 开始 插入基本参数 插入其他参数 形成整体刚度矩阵 形成荷载列阵 引入支承条件 解方程输出位移 求应力输出应力 结束 子程序单元刚度矩阵 （ASK=3）计算单元刚度矩阵 （ASK=1）计算单元面积 （ASK=2）计算S矩阵 四．源程序 #includestdio.() { int NJ1,k,IN,IM,jn,m,i,j,z,JO,ii,jj,,o,f; double F,c,SIG1,SIG2,SIG3,PYL,RYL,MAYL,MIYL,CETA; double WY[8+1],YL[3+1]; ab=fopen(节点应力.txt,w); fp1=fopen(节点位移.txt,w); fp2=fopen(单元应力.txt,w); if(LXM!=0) { E0=E0/(1.0-MU*MU); MU=MU/(1.0-MU); } for(i=0;i=NJ2;i++) { for(j=0;j=DD;j++) KZ[i][j]=0.0; } for(E=1;E=NE;E++) { DUGD(E,3); for(i=1;i=4;i++) { for(ii=1;ii=2;ii++) { =DD-l+1; for(j=1;j=jn;j++) { m=j+i-k; KZ[i][j]=KZ[i][j]-c*KZ[k][m]; } P[i]=P[i]-c*P[k]; } } P[NJ2]=P[NJ2]/KZ[NJ2][1]; for(i=NJ1;i=1;i--) { if(DD=NJ2-i+1) JO=NJ2-i+1; else JO=DD; for(j=2;j=JO;j++) { ); printf(JD U V\n); fprintf(fp1,JD U V\n); for(i=1;i=NJ;i++) { printf(%d %-9.6f %-9.6f\n,i,P[2*i-1],P[2*i]); fprintf(fp1,%d %-9.6f %-9.6f\n,i,P[2*i-1],P[2*i]); } for(i=1;i=NJ;i++) { fprintf(fp1,a%d=%-9.6f;b%d=%-9.6f;\n,i,AJZ[i][1]+P[2*i-1],i,AJZ[i][2]+P[2*i]); } for(E=1;E=NE;E++) { DUGD(E,2); for(i=1;i=4;i++) { for(j=1;j=2;j++) { =1;n=5;n++) { for(i=1;i=3;i++) { YL[i]=0; for(j=1;j=8;j++) YL[i]=YL[i]+SZ[i][8*(n-1)+j]*WY[j]; } SIG1=YL[1]; SIG2=YL[2]; SIG3=YL[3]; PYL=(SIG1+SIG2)/2; RYL=sqrt(pow((SIG1-SIG2)/2.0,2)+pow(SIG3,2)); MAYL=PYL+RYL; MIYL=PYL-RYL; if(SIG2==MIYL) CETA=0; else CETA=90-57.29578*atan2(SIG3,(SIG2-MIYL)); printf(\n); printf(E=%d(%d)\n,E,n); printf(sx=%-9.6f sy=%-9.6f tou=%-9.6f\n,SIG1,SIG2,SIG3); printf(s1=%-9.6f