5-1.c

简单有限元计算程序。 NJF---结点自由度数 E---弹性模量 NE---单元数 T---板厚度 VM---泊松比 z---实型二维数组

#include <stdio.h>
#include <conio.h>
main()
{ int i,j,je,k,i1,i2,i3,j3,j8,j4,i4,i8,NE,NJ,NZ,NG,KK,NJF,ie[98][3];
float E,VM,T,cc;
float b[3][3],c[3][3],aa[98],a[98],g[2][2],zk[20][20],ZJ[90][2],z[99][2],p[20];

printf("E=\n");
scanf("%f",&E);

printf("VM=\n");
scanf("%f",&VM);

printf("T=\n");
scanf("%f",&T);



printf("NE=\n");
scanf("%d",&NE); 


printf("NJ=\n");
scanf("%d",&NJ); 


printf("NZ=\n");
scanf("%d",&NZ); 


printf("NG=\n");
scanf("%d",&NG); 


printf("KK=\n");
scanf("%d",&KK); 



printf("Input ZJ number is:");
   for(i=0;i<NZ;i++)
   for(j=0;j<2;j++)
   scanf("%f",&ZJ[i][j]);
  
printf("Input ie number is:");
   for(i=0;i<NE;i++)
   for(j=0;j<3;j++)
   scanf("%d",&ie[i][j]);
   
printf("Input Z number is:");
   for (i=0;i<NJ;i++)
   for (j=0;j<2;j++)
   scanf("%f",&z[i][j]);
   
printf("Input p number is:");
   for (i=0;i<NG;i++)
   scanf("%f",&p[i]);
  

if (KK==0)
   E=E;
   VM=VM;
if (KK==1)
   E=E/(1-VM*VM);
   VM=VM/(1-VM);

for(je=0;je<NE;je++)
	 {i1=ie[je][0]-1;
	  i2=ie[je][1]-1;
	  i3=ie[je][2]-1;
	  b[je][0]=z[i2][1]-z[i3][1];
	  b[je][1]=z[i3][1]-z[i1][1];
	  b[je][2]=z[i1][1]-z[i2][1];
	  c[je][0]=z[i3][0]-z[i2][0];
	  c[je][1]=z[i1][0]-z[i3][0];
	  c[je][2]=z[i2][0]-z[i1][0];
	  aa[je]=a[je]=0.5*(b[je][0]*c[je][1]-b[je][1]*c[je][0]);
	  a[je]=E*T/(4.0*(1.0-VM*VM)*a[je]);



 for(i=0;i<3;i++)
	  {for(j=i;j<3;j++)
	  {g[0][0]=b[je][i]*b[je][j]+0.5*(1.0-VM)*c[je][i]*c[je][j];
	   g[1][1]=c[je][i]*c[je][j]+0.5*(1.0-VM)*b[je][i]*b[je][j];
	   g[1][0]=VM*b[je][j]*c[je][i]+0.5*(1.0-VM)*b[je][i]*c[je][j];
	   g[0][1]=VM*c[je][j]*b[je][i]+0.5*(1.0-VM)*c[je][i]*b[je][j];
	   for(i3=0;i3<2;i3++)
	   {if(i>=j) j8=i3; else j8=0;
	   for(j3=j8;j3<2;j3++)
	   {i4=2*(ie[je][i]-1)+i3;
	    j4=2*(ie[je][j]-1)+j3;
		if(i4>j4) {i8=i4;i4=j4;j4=i8;}
		zk[i4][j4]+=g[i3][j3]*a[je];
	 }}}}}
for(i=0;i<8;i++)
{for(j=i;j<8;j++)
{zk[j][i]=zk[i][j];}
}


for(i=0;i<NZ;i++)
    j=ZJ[i][0];
    zk[j][j]=zk[j][j]*1.0e+20;
    p[j]=zk[j][j]*ZJ[i][2];

NJF=2*NJ;



	for(k=0;k<NJF-1;k++)
		{for(i=k+1;i<NJF;i++)
		{cc=zk[k][i]/z[k][k];
		 for(j=0;j<NJF;j++)
			 zk[i][j]-=cc*zk[k][j];
		 p[i]=p[i]-cc*p[k];}
		}
		p[NJF-1]/=zk[NJF-1][NJF-1];
		for(i=NJF-2;i>=0;i--)
		{for(j=i+1;j<NJF;j++)
		p[i]-=zk[i][j]*p[j];
		p[i]/=zk[i][i];
		}

		for(i=0;i<8;i++)
		{
			for(j=0;j<8;j++)
			{
				printf("%d\n",zk[i][j]);

			}
		
		
		}

}