1、一 题目如图所示固支梁，高 3米，长 15m，承受均布载荷q=10kN/，E=20GPa，u=0.167，厚度 t=1m，忽略自重，按平面应力问题分析。用有限元方法计算梁的变形及应力分布，要求用矩形单元。要求：1.单元数目不得少于 20个。2.采用矩形单元计算求解。3.计算结果并给出变形图、应力分布图、单元划分图。二、力学分析1.题目可以看做是平面应力问题故有 LXM=02.单元划分图q=10kN/三程序框图及程序程序框图：开始插入基本参数插入其他参数形成整体刚度矩阵形成荷载列阵引入支承条件解方程输出位移求应力输出应力结束子程序单元刚度矩阵（ASK=3）计算单元刚度矩阵（ASK=1）计算单元面
2、积（ASK=2）计算 S 矩阵四源程序#include#include#define NE 30 /单元数#define NJ 44 /节点数#define NZ 16 /支承数#define NPJ 11 /节点载荷作用数#define DD 26 /半带宽#define NJ2 88 /节点位移数int LXM=0;double E0=2e7;double MU=0.167;double LOU=0.0;double TE=1 ;double AJZNJ+13=0,0,0,0,0,3,0,1.5,3,0,3,3,0,4.5,3,0,6,3,0,7.5,3,0,9,3,0,10.5,3,0,
3、12,3,0,13.5,3,0,15,3,0,0,2,0,1.5,2,0,3,2,0,4.5,2,0,6,2,0,7.5,2,0,9,2,0,10.5,2,0,12,2,0,13.5,2,0,15,2,0,0,1,0,1.5,1,0,3,1,0,4.5,1,0,6,1,0,7.5,1,0,9,1,0,10.5,1,0,12,1,0,13.5,1,0,15,1,0,0,0,0,1.5,0,0,3,0,0,4.5,0,0,6,0,0,7.5,0,0,9,0,0,10.5,0,0,12,0,0,13.5,0,0,15,0; /节点坐标int JMNE+15=0,0,0,0,0,0,10,11,22,
4、21,0,9,10,21,20,0,8,9,20,19,0,7,8,19,18,0,6,7,18,17,0,5,6,17,16,0,4,5,16,15,0,3,4,15,14,0,2,3,14,13,0,1,2,13,12,0,21,22,33,32,0,20,21,32,31,0,19,20,31,30,0,18,19,30,29,0,17,18,29,28,0,16,17,28,27,0,15,16,27,26,0,14,15,26,25,0,13,14,25,24,0,12,13,24,23,0,32,33,44,43,0,31,32,43,42,0,30,31,42,41,0,29,3
5、0,41,40,0,28,29,40,39,0,27,28,39,38,0,26,27,38,37,0,25,26,37,36,0,24,25,36,35,0,23,24,35,34;int NZCNZ+1=0,1,2,21,22,23,24,43,44,45,46,65,66,67,68,87,88 ; /支撑数组double PJNPJ+12+1=0,0,0,0,-7500,68,0,-15000,70,0,-15000,72,0,-15000,74,0,-15000,76,0,-15000,78,0,-15000,80,0,-15000,82,0,-15000,84,0,-15000,8
6、6,0,-7500,88; /节点载荷数组 /节点载荷数组double AE,KZNJ2+1DD+1,PNJ2+1,S3+18+1,KE8+18+1,SZ3+15*8+1;double JDYLNJ6; /节点应力矩阵double DYYL4NE4; /单元应力矩阵int IE,JE,ME,PE;void DUGD(int,int);FILE *fp1,*fp2,*ab;void main()int NJ1,k,IN,IM,jn,m,i,j,z,JO,ii,jj,h,dh,E,l,zl,dl,r,n,o,f;double F,c,SIG1,SIG2,SIG3,PYL,RYL,MAYL,MIYL
7、,CETA;double WY8+1,YL3+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;i0)KZdhdl=KZdhdl+KEhl;for(i=1;i0)for(i=1;i0)for(E=1;EDD)JO=DD;elseJO=z;for(j=2;jk+DD-1)IM=k+DD-1;elseIM=NJ2;IN=k+1;for(i=IN;i=1;i-)if(DD=NJ2-i+1)JO=NJ2-i+1;elseJO=DD;for(j=2;j1)for(w=1;w=5;w+)for(i=0;i=3;i+)for(j=0;j=8;j+)Bij=0.0;x=xyw1*a;y=xyw2*b;B11=(y-b)/AE;B13=(b-y)/AE;B15=(b+y)/AE;B17=-(b+y)/AE;B22=(x-a)/AE;B24=-(a+x)/AE;B26=(a+x)/AE;B28=(a-x)/AE;B31=(x-a)/AE;B32=(y-b)/AE;B33=-(a+x)/AE;
