有限单元法程序代码.doc

弹性力学平面问题有限元程序#define NODE_NUM 300 #define ELE_NUM 200 #define NODE2 NODE_NUM*2 #include<> #include<> #include<> int nj,ne,nm,nz,npj,jm[ELE_NUM][4],nzc[60],nj2,bw,lxm; float ea,e,um,oul,te,xy[NODE_NUM][4],EA[5],p[NODE2],pj,F[ELE_NUM],etm[5][5]; float ek[7][7]; double **K,*f; int BAND2(int n,int m, double **K,double *f); double **TwoArrayAlloc(int,int); void TwoArrayFree(double **); FILE *outf; /*-------------------------------------------*/ double **TwoArrayAlloc(int r,int c) { double *x,**y; int n; x=(double*)calloc(r*c,sizeof(double)); y=(double**)calloc(r,sizeof(double*)); for(n=0;n<=r-1;++n) y[n]=&x[c*n]; return(y); } /*-------------------------------------------*/ void TwoArrayFree(double **x) { free(x[0]);free(x); } /*------------------------------------------*/ int BAND2(int n,int m, double **K,double *f) { int i,j,t,ij,ji,it,jt,tm,m1; double s,w; m1=m+1; for(i=1;i<=n;i++) { if(K[i-1][m1-1]<=0) return(0); w=; if(i>m1) tm=i-m;else tm=1; for(j=tm;j<=i;j++) { s=; ij=j-i+m1; for(t=tm;t<=j-1;t++) { it=t-i+m1; jt=t-j+m1; s=s+K[i-1][it-1]*K[j-1][jt-1]/K[t-1][m1-1]; } K[i-1][ij-1]=K[i-1][ij-1]-s; if(j==i) f[i-1]=f[i-1]-w; else w=w+K[i-1][ij-1]*f[j-1]/K[j-1][m1-1]; }} for(i=n;i>=1;i--) { s=; if(i>n-m1) tm=n; else tm=i+m; for(j=i+1;j<=tm;j++) { ji=i-j+m1; s=s+K[j-1][ji-1]*f[j-1]; } f[i-1]=(f[i-1]-s)/K[i-1][m1-1]; } return 1;} /*-----------------------------------------*/ void input() { int jj,j,i,nh,nl; float dx,dy; FILE *infile; char name1[30],name2[30]; R:printf("please enter data-filename\n"); scanf("%s",name1); if((infile=fopen(name1,"r"))==NULL) { printf("the data-file not exit!"); goto R;} printf("Please enter result-filename\n"); scanf("%s",name2); outf=fopen(name2,"w"); fscanf(infile,"%d%d%d%d%d%d5d",&nj,&ne,&nm,&nz,&npj,&lxm); for(i=1;i<=nj;i++) { for(j=1;j<3;j++) fscanf(infile,"%f",&xy[i][j]); } for(i=1;i<=ne;i++) { for(j=0;j<4;j++) fscanf(infile,"%d",&jm[i][j]); } for(i=1;i<=nm;i++) { for(j=1;j<=4;j++) fscanf(infile,"%f",&etm[i][