eld.h

这是等参单元的有限元程序

//计算D矩阵
#ifndef __ELD_
#define __ELD_
void eld(double E,double u,double (*D)[4],int IK);

void eld(double E,double u,double (*D)[4],int IK)
/*参数说明
     E—是杨氏弹性模量；
     u—是泊松比；
  IK—＝0是平面应力，＝1是平面应变，＝2是轴对称
*/
{     int      i,j;
      double t,a,b;
        for(i=0;i<4;i++)
			for(j=0;j<4;j++)
				D[i][j]=0.0;
        if(IK)
		{ t=E*(1-u)/(1+u)/(1-2*u); 
		  a=t*u/(1-u); b=0.5*t*(1-2*u)/(1-u);
          D[0][0]=t;  D[1][1]=t;  D[2][2]=b;  D[3][3]=t; D[0][1]=a; 
          D[1][0]=a;  D[0][3]=a;  D[3][0]=a;  D[1][3]=a; D[3][1]=a;
		}
          else 
		  {   t=E/(1-u*u); D[0][0]=t;  D[1][1]=t;  
              D[2][2]=0.5*t*(1.0-u);  D[0][1]=t*u;  D[1][0]=t*u;
		  }
}

#endif