有限元分析中《结构力学》矩阵位移法C语言程序(附例题).doc

程序：#include stdafx.h#include stdio.h#include math.h#include stdlib.hvoid main()int loc32=0,ifix6=0; float area3=0.0,fint3=0.0,cx4=0.0,cy4=0.0,f12=0.0,fr12=0.0,fe36=0.0; int nn,ne,nd,nfix; float ea; int i,j,k; FILE *shuru,*shuchu; shuru=fopen(shuru.dat,r); shuchu=fopen(shuchu.dat,w); fscanf(shuru,%d%d%d%d%f,&nn,&ne,&nd,&nfix,&ea); fprintf(shuchu,nn ne nd nfix en%d %d %d %d %fn,nn,ne,nd,nfix,ea); i=0; while(i=ne-1) fscanf(shuru,%d%d%f%f,&loci0,&loci1,&areai,&finti); i+; fprintf(shuchu,element node1 node2 area fintn); i=0; while(i=ne-1) fprintf(shuchu,%d %d %d %f %fn,i+1,loci0,loci1,areai,finti); i+; j=0; while(j=nn-1) fscanf(shuru,%f%f,&cxj,&cyj); j+; fprintf(shuchu,node x-coord y-coordn); j=0; while(j=nn-1) fprintf(shuchu,%d %f %fn,j+1,cxj,cyj); j+; k=0; while(k=nfix-1) fscanf(shuru,%d,&ifixk); k+; fprintf(shuchu,ifix=); k=0; while(k=nfix-1) fprintf(shuchu,%d ,ifixk); k+; fprintf(shuchu,n); void cst(int (*loc)2,int *ifix,float *area,float *fint,float *cx,float *cy,float *f,float *fr,float (*fe)6,FILE *shuru,FILE *shuchu,float ea); cst(loc,ifix,area,fint,cx,cy,f,fr,fe,shuru,shuchu,ea); fprintf(shuchu,node x-disp y-disp thitan); i=0; while(i=3) fprintf(shuchu,%d %f %f %fn,i+1,f3*i,f3*i+1,f3*i+2); i+; fprintf(shuchu,reaction nodal forces from the equationsn); fprintf(shuchu,node x-load y-load momentn); i=0; while(i=3) fprintf(shuchu,%d %f %f %fn,i+1,fr3*i,fr3*i+1,fr3*i+2); i+; fprintf(shuchu,element axi-f shear-q moment-mn); i=0; while(i=ne-1) fprintf(shuchu,%d %f %f %f %f %f %fn,i+1,fei0,fei1,fei2,fei3,fei4,fei5); i+; fclose(shuru); fclose(shuchu);void cst(int (*loc)2,int *ifix,float *area,float *fint,float *cx,float *cy,float *f,float *fr,float (*fe)6,FILE *shuru,FILE *shuchu,float ea)int np,nvd; float p136=0.0,p236=0.0,gk1212=0.0,gk11212=0.0,al3=0.0,tt366=0.0,bkl366=0.0; float t66=0.0,css3=0.0,snn3=0.0,ek66=0.0,ekl66=0.0,ekk366=0.0,xx6=0.0,ba66=0.0; int nn=4,ne=3,nd=12,nfix=6; int ii,jj,i,j,k,l,inode,nodei,idofn,nrows,nrowe,jnode,nodej,jdofn,ncols,ncole; int i1,i2,ie,ix; float x12,y12,q,eal,eil1,eil2,eil3; i=0; for(;i=ne-1;i+) i1=loci0; i2=loci1; x12=cxi2-1-cxi1-1; y12=cyi2-1-cyi1-1; ali=sqrt(pow(x12,2)+pow(y12,2); cssi=x12/ali; snni=y12/ali; fscanf(shuru,%d%d,&np,&nvd); if(np!=0) i=0; for(;i=np-1;i+) fscanf(shuru,%d%f%f%f,&i,&f3*i,&f3*i+1,&f3*i+2); if(nvd!=0) i=0; for(;i=nvd-1;i+) fscanf(shuru,%d%f,&ie,&q); i1=locie-10; i2=locie-11; p1ie-11=q*alie-1/2; p1ie-12=q*alie-1*alie-1/12; p1ie-14=q*alie-1/2; p1ie-15=-q*alie-1*alie-1/12; i=0; for(;i=nd-1;i+) j=0; for(;j=nd-1;j+) gkij=0.0; for(i=0;i=ne-1;i+) j=0; for(;j=5;j+) k=0; for(;k=5;k+) ekljk=0.0; ekjk=0.0; tjk=0.0; eal=ea*areai/ali; eil1=ea*finti/ali; eil2=ea*finti/(ali*ali); eil3=ea*finti/(ali*ali*ali); ekl00=eal; ekl11=12.0*eil3; ekl22=4.0*eil1; ekl33=eal; ekl44=12.0*eil3; ekl55=4.0*eil1; ekl21=6.0*eil2; ekl30=-eal; ekl41=-12.0*eil3; ekl42=-6.0*eil2; ekl51=6.0*eil2; ekl52=2.0*eil1; ekl54=-6.0*eil2; ii=0; for(;ii=4;ii+) jj=ii+1; for(;jj=5;jj+) ekliijj=ekljjii; k=0; for(;k=5;k+) l=0; for(;l=5;l+) ekkikl=eklkl; fprintf(shuchu,%d %d %d %f %fn,i+1,k+1,l+1,eklkl,ekkikl); t00=cssi; t01=-snni; t10=snni; t11=cssi; t22=1.0; t33=cssi; t34=-snni; t43=snni; t44=cssi; t55=1.0; j=0; for(;j=5;j+) k=0; for(;k=5;k+) ttijk=tjk; p2ij=p2ij+tjk*p1ik; ii=0; for(;ii=5;ii+) j=0; for(;j=5;j+) baiij=0.0; k=0; for(;k=5;k+) baiij=baiij+ttiiik*eklkj; ekiij=0.0; ii=0; for(;ii=5;ii+) j=0; for(;j=5;j+) k=0; for(;k=5;k+) ekiij=ekiij+baiik*ttijk; j=0; for(;j=5;j+) ii=0; for(;ii=5;ii+) fprintf(shuchu,i,ii,j,ek,tt=%d %d %d %f %fn,i+1,ii+1,j+1,ekiij,ttiiij); inode=0; while(inode=1) nodei=lociinode; idofn=0; while(idofn=2) nrows=(nodei-1)*3+idofn; nrowe=inode*3+idofn; fnrows=fnrows+p2inrowe; jnode=0; while(jnode=1) nodej=locijnode; jdofn=0; while(jdofn=2) ncols=(nodej-1)*3+jdofn; ncole=jnode*3+jdofn; gknrowsncols=gknrowsncols+eknrowencole; jdofn+; jnode+; idofn+; inode+; i=0; for(;i=nd-1;i+) j=0; for(;j=nd-1;j+) gk1ij=gkij; fprintf(shuchu,nodal forces from applied loadsnode x-load y-load momentn); i=0; for(;i=nn-1;i+) fprintf(shuchu,%d %f %f %fn,i+1,f3*i,f3*i+1,f3*i+2); i=0; for(;i=nd-1;i+) j=0; for(;j=nd-1;j+) fprintf(shuchu,i,j,gk1 %d %d %fn,i+1,j+1,gkij); i=0; for(;i=nfix-1;i+) ix=ifixi; gkix-1ix-1=gkix-1ix-1*1.0e20; void gauss(float (*a)12,float *b,int n); gauss(gk,f,nd); i=0; for(;i=nd-1;i+) fri=0.0; j=0; for(;j=nd-1;j+) fri=fri+gk1ij*fj; fri=fri-fi; for(i=0;i=ne-1;i+) i1=loci0; i2=loci1; xx0=f3*i1-3; xx1=f3*i1-2; xx2=f3*i1-1; xx3=f3*i2-3; xx4=f3*i2-2; xx5=f3*i2-1; j=0; for(;j=5;j+) k=0; for(;k=5;k+) bklijk=0.0; l=0; for(;l=5;l+) bklijk=bklijk+ekkijl*ttikl; j=0; for(;j=5;j+) feij=0.0; k=0; for(;k=5;k+) feij=feij+bklijk*xxk; j=0; for(;j=5;j+) feij=feij-p1ij; void gauss(float (*a)12,float *b,int n)int i,i1,j,m;i=0; for(;i=n-1;i+) i1=i+1; for(j=i1;j=n-1;j+) aij=aij/aii; bi=bi/aii; aii=1.0; for(j=i1;j=n-1;j+) for(m=i1;m=0;i-) j=i+1; for(;j=n-1;j+) bi=bi-aij*bj; 一、如图所示平面刚架的内力，各杆面积A=76.3cm2,惯性矩I=15760cm4，弹性模量E=2105MPa程序：#include stdafx.h#include stdio.h#include math.h#include stdlib.hvoid main()int loc32=0,ifix6=0; float area3=0.0,fint3=0.0,cx4=0.0,cy4=0.0,f12=0.0,fr12=0.0,fe36=0.0; int nn,ne,nd,nfix; float ea; int i,j,k; FILE *shuru,*shuchu; shuru=fopen(shuru.dat,r); shuchu=fopen(shuchu.dat,w); fscanf(shuru,%d%d%d%d%f,&nn,&ne,&nd,&nfix,&ea); fprintf(shuchu,nn ne nd nfix en%d %d %d %d %fn,nn,ne,nd,nfix,ea); i=0; while(i=ne-1) fscanf(shuru,%d%d%f%f,&loci0,&loci1,&areai,&finti); i+; fprintf(shuchu,element node1 node2 area fintn); i=0; while(i=ne-1) fprintf(shuchu,%d %d %d %f %fn,i+1,loci0,loci1,areai,finti); i+; j=0; while(j=nn-1) fscanf(shuru,%f%f,&cxj,&cyj); j+; fprintf(shuchu,node x-coord y-coordn); j=0; while(j=nn-1) fprintf(shuchu,%d %f %fn,j+1,cxj,cyj); j+; k=0; while(k=nfix-1) fscanf(shuru,%d,&ifixk); k+; fprintf(shuchu,ifix=); k=0; while(k=nfix-1) fprintf(shuchu,%d ,ifixk); k+; fprintf(shuchu,n); void cst(int (*loc)2,int *ifix,float *area,float *fint,float *cx,float *cy,float *f,float *fr,float (*fe)6,FILE *shuru,FILE *shuchu,float ea); cst(loc,ifix,area,fint,cx,cy,f,fr,fe,shuru,shuchu,ea); fprintf(shuchu,node x-disp y-disp thitan); i=0; while(i=3) fprintf(shuchu,%d %f %f %fn,i+1,f3*i,f3*i+1,f3*i+2); i+; fprintf(shuchu,reaction nodal forces from the equationsn); fprintf(shuchu,node x-load y-load momentn); i=0; while(i=3) fprintf(shuchu,%d %f %f %fn,i+1,fr3*i,fr3*i+1,fr3*i+2); i+; fprintf(shuchu,element axi-f shear-q moment-mn); i=0; while(i=ne-1) fprintf(shuchu,%d %f %f %f %f %f %fn,i+1,fei0,fei1,fei2,fei3,fei4,fei5); i+; fclose(shuru); fclose(shuchu);void cst(int (*loc)2,int *ifix,float *area,float *fint,float *cx,float *cy,float *f,float *fr,float (*fe)6,FILE *shuru,FILE *shuchu,float ea)int np,nvd; float p136=0.0,p236=0.0,gk1212=0.0,gk11212=0.0,al3=0.0,tt366=0.0,bkl366=0.0; float t66=0.0,css3=0.0,snn3=0.0,ek66=0.0,ekl66=0.0,ekk366=0.0,xx6=0.0,ba66=0.0; int nn=4,ne=3,nd=12,nfix=6; int ii,jj,i,j,k,l,inode,nodei,idofn,nrows,nrowe,jnode,nodej,jdofn,ncols,ncole; int i1,i2,ie,ix; float x12,y12,q,eal,eil1,eil2,eil3; i=0; for(;i=ne-1;i+) i1=loci0; i2=loci1; x12=cxi2-1-cxi1-1; y12=cyi2-1-cyi1-1; ali=sqrt(pow(x12,2)+pow(y12,2); cssi=x12/ali; snni=y12/ali; fscanf(shuru,%d%d,&np,&nvd); if(np!=0) i=0; for(;i=np-1;i+) fscanf(shuru,%d%f%f%f,&i,&f3*i,&f3*i+1,&f3*i+2); if(nvd!=0) i=0; for(;i=nvd-1;i+) fscanf(shuru,%d%f,&ie,&q); i1=locie-10; i2=locie-11; p1ie-11=q*alie-1/2; p1ie-12=q*alie-1*alie-1/12; p1ie-14=q*alie-1/2; p1ie-15=-q*alie-1*alie-1/12; i=0; for(;i=nd-1;i+) j=0; for(;j=nd-1;j+) gkij=0.0; for(i=0;i=ne-1;i+) j=0; for(;j=5;j+) k=0; for(;k=5;k+) ekljk=0.0; ekjk=0.0; tjk=0.0; eal=ea*areai/ali; eil1=ea*finti/ali; eil2=ea*finti/(ali*ali); eil3=ea*finti/(ali*ali*ali); ekl00=eal; ekl11=12.0*eil3; ekl22=4.0*eil1; ekl33=eal; ekl44=12.0*eil3; ekl55=4.0*eil1; ekl21=6.0*eil2; ekl30=-eal; ekl41=-12.0*eil3; ekl42=-6.0*eil2; ekl51=6.0*eil2; ekl52=2.0*eil1; ekl54=-6.0*eil2; ii=0; for(;ii=4;ii+) jj=ii+1; for(;jj=5;jj+) ekliijj=ekljjii; k=0; for(;k=5;k+) l=0; for(;l=5;l+) ekkikl=eklkl; fprintf(shuchu,%d %d %d %f %fn,i+1,k+1,l+1,eklkl,ekkikl); t00=cssi; t01=-snni; t10=snni; t11=cssi; t22=1.0; t33=cssi; t34=-snni; t43=snni; t44=cssi; t55=1.0; j=0; for(;j=5;j+) k=0; for(;k=5;k+) ttijk=tjk; p2ij=p2ij+tjk*p1ik; ii=0; for(;ii=5;ii+) j=0; for(;j=5;j+) baiij=0.0; k=0; for(;k=5;k+) baiij=baiij+ttiiik*eklkj; ekiij=0.0; ii=0; for(;ii=5;ii+) j=0; for(;j=5;j+) k=0; for(;k=5;k+) ekiij=ekiij+baiik*ttijk; j=0; for(;j=5;j+) ii=0; for(;ii=5;ii+) fprintf(shuchu,i,ii,j,ek,tt=%d %d %d %f %fn,i+1,ii+1,j+1,ekiij,ttiiij); inode=0; while(inode=1) nodei=lociinode; idofn=0; while(idofn=2) nrows=(nodei-1)*3+idofn; nrowe=inode*3+idofn; fnrows=fnrows+p2inrowe; jnode=0; while(jnode=1) nodej=locijnode; jdofn=0; while(jdofn=2) ncols=(nodej-1)*3+jdofn; ncole=jnode*3+jdofn; gknrowsncols=gknrowsncols+eknrowencole; jdofn+; jnode+; idofn+; inode+; i=0; for(;i=nd-1;i+) j=0; for(;j=nd-1;j+) gk1ij=gkij; fprintf(shuchu,nodal forces from applied loadsnode x-load y-load momentn); i=0; for(;i=nn-1;i+) fprintf(shuchu,%d %f %f %fn,i+1,f3*i,f3*i+1,f3*i+2); i=0; for(;i=nd-1;i+) j=0; for(;j=nd-1;j+) fprintf(shuchu,i,j,gk1 %d %d %fn,i+1,j+1,gkij); i=0; for(;i=nfix-1;i+) ix=ifixi; gkix-1ix-1=gkix-1ix-1*1.0e20; void gauss(float (*a)12,float *b,int n); gauss(gk,f,nd); i=0; for(;i=nd-1;i+) fri=0.0; j=0; for(;j=nd-1;j+) fri=fri+gk1ij*fj; fri=fri-fi; for(i=0;i=ne-1;i+) i1=loci0; i2=loci1; xx0=f3*i1-3; xx1=f3*i1-2; xx2=f3*i1-1; xx3=f3*i2-3; xx4=f3*i2-2; xx5=f3*i2-1; j=0; for(;j=5;j+) k=0; for(;k=5;k+) bklijk=0.0; l=0; for(;l=5;l+) bklijk=bklijk+ekkijl*ttikl; j=0; for(;j=5;j+) feij=0.0; k=0; for(;k=5;k+) feij=feij+bklijk*xxk; j=0; for(;j=5;j+) feij=feij-p1ij; void gauss(float (*a)12,float *b,int n)int i,i1,j,m;i=0; for(;i=n-1;i+) i1=i+1; for(j=i1;j=n-1;j+) aij=aij/aii; bi=bi/aii; aii=1.0; for(j=i1;j=n-1;j+) for(m=i1;m=0;i-) j=i+1; for(;j=n-1;j+) bi=bi-aij*bj; 程序输出：nn ne nd nfix e4 3 12 6 200000000.000000element node1 node2 area fint1 1 2 0.007630 0.0001582 3 1 0.007630 0.0001583 2 4 0.007630 0.000158node x-coord y-coord1 0.000000 5.0000002 6.400000 5.0000003 0.000000 0.0000004 9.600000 0.000000ifix=7 8 9 10 11 12 1 1 1 238437.500000 238437.5000001 1 2 0.000000 0.0000001 1 3 0.000000 0.0000001 1 4 -238437.500000 -238437.5000001 1 5 0.000000 0.0000001 1 6 0.000000 0.0000001 2 1 0.000000 0.0000001 2 2 1442.871094 1442.8710941 2 3 4617.187500 4617.1875001 2 4 0.000000 0.0000001 2 5 -1442.871094 -1442.8710941 2 6 4617.187500 4617.1875001 3 1 0.000000 0.0000001 3 2 4617.187500 4617.1875001 3 3 19700.000000 19700.0000001 3 4 0.000000 0.0000001 3 5 -4617.187500 -4617.1875001 3 6 9850.000000 9850.0000001 4 1 -238437.500000 -238437.5000001 4 2 0.000000 0.0000001 4 3 0.000000 0.0000001 4 4 238437.500000 238437.5000001 4 5 0.000000 0.0000001 4 6 0.000000 0.0000001 5 1 0.000000 0.0000001 5 2 -1442.871094 -1442.8710941 5 3 -4617.187500 -4617.1875001 5 4 0.000000 0.0000001 5 5 1442.871094 1442.8710941 5 6 -4617.187500 -4617.1875001 6 1 0.000000 0.0000001 6 2 4617.187500 4617.1875001 6 3 9850.000000 9850.0000001 6 4 0.000000 0.0000001 6 5 -4617.187500 -4617.1875001 6 6 19700.000000 19700.000000i,ii,j,ek,tt=1 1 1 238437.500000 1.000000i,ii,j,ek,tt=1 2 1 0.000000 0.000000i,ii,j,ek,tt=1 3 1 0.000000 0.000000i,ii,j,ek,tt=1 4 1 -238437.500000 0.000000i,ii,j,ek,tt=1 5 1 0.000000 0.000000i,ii,j,ek,tt=1 6 1 0.000000 0.000000i,ii,j,ek,tt=1 1 2 0.000000 0.000000i,ii,j,ek,tt=1 2 2 1442.871094 1.000000i,ii,j,ek,tt=1 3 2 4617.187500 0.000000i,ii,j,ek,tt=1 4 2 0.000000 0.000000i,ii,j,ek,tt=1 5 2 -1442.871094 0.000000i,ii,j,ek,tt=1 6 2 4617.187500 0.000000i,ii,j,ek,tt=1 1 3 0.000000 0.000000i,ii,j,ek,tt=1 2 3 4617.187500 0.000000i,ii,j,ek,tt=1 3 3 19700.000000 1.000000i,ii,j,ek,tt=1 4 3 0.000000 0.000000i,ii,j,ek,tt=1 5 3 -4617.187500 0.000000i,ii,j,ek,tt=1 6 3 9850.000000 0.000000i,ii,j,ek,tt=1 1 4 -238437.500000 0.000000i,ii,j,ek,tt=1 2 4 0.000000 0.000000i,ii,j,ek,tt=1 3 4 0.000000 0.000000i,ii,j,ek,tt=1 4 4 238437.500000 1.000000i,ii,j,ek,tt=1 5 4 0.000000 0.000000i,ii,j,ek,tt=1 6 4 0.000000 0.000000i,ii,j,ek,tt=1 1 5 0.000000 0.000000i,ii,j,ek,tt=1 2 5 -1442.871094 0.000000i,ii,j,ek,tt=1 3 5 -4617.187500 0.000000i,ii,j,ek,tt=1 4 5 0.000000 0.000000i,ii,j,ek,tt=1 5 5 1442.871094 1.000000i,ii,j,ek,tt=1 6 5 -4617.187500 0.000000i,ii,j,ek,tt=1 1 6 0.000000 0.000000i,ii,j,ek,tt=1 2 6 4617.187500 0.000000i,ii,j,ek,tt=1 3 6 9850.000000 0.000000i,ii,j,ek,tt=1 4 6 0.000000 0.000000i,ii,j,ek,tt=1 5 6 -4617.187500 0.000000i,ii,j,ek,tt=1 6 6 19700.000000 1.0000002 1 1 305200.000000 305200.0000002 1 2 0.000000 0.0000002 1 3 0.000000 0.0000002 1 4 -305200.000000 -305200.0000002 1 5 0.000000 0.0000002 1 6 0.000000 0.0000002 2 1 0.000000 0.0000002 2 2 3025.919922 3025.9199222 2 3 7564.800293 7564.8002932 2 4 0.000000 0.0000002 2 5 -3025.919922 -3025.9199222 2 6 7564.800293 7564.8002932 3 1 0.000000 0.0000002 3 2 7564.800293 7564.8002932 3 3 25216.001953 25216.0019532 3 4 0.000000 0.0000002 3 5 -7564.800293 -7564.8002932 3 6 12608.000977 12608.0009772 4 1 -305200.000000 -305200.0000002 4 2 0.000000 0.0000002 4 3 0.000000 0.0000002 4 4 305200.000000 305200.0000002 4 5 0.000000 0.0000002 4 6 0.000000 0.0000002 5 1 0.000000 0.0000002 5 2 -3025.919922 -3025.9199222 5 3 -7564.800293 -7564.8002932 5 4 0.000000 0.0000002 5 5 3025.919922 3025.9199222 5 6 -7564.800293 -7564.8002932 6 1 0.000000 0.0000002 6 2 7564.800293 7564.8002932 6 3 12608.000977 12608.0009772 6 4 0.000000 0.0000002 6 5 -7564.800293 -7564.8002932 6 6 25216.001953 25216.001953i,ii,j,ek,tt=2 1 1 3025.919922 0.000000i,ii,j,ek,tt=2 2 1 0.000000 1.000000i,ii,j,ek,tt=2 3 1 -7564.800293 0.000000i,ii,j,ek,tt=2 4 1 -3025.919922 0.000000i,ii,j,ek,tt=2 5 1 0.000000 0.000000i,ii,j,ek,tt=2 6 1 -7564.800293 0.000000i,ii,j,ek,tt=2 1 2 0.000000 -1.000000i,ii,j,ek,tt=2 2 2 305200.000000 0.000000i,ii,j,ek,tt=2 3 2 0.000000 0.000000i,ii,j,ek,tt=2 4 2 0.000000 0.000000i,ii,j,ek,tt=2 5 2 -305200.000000 0.000000i,ii,j,ek,tt=2 6 2 0.000000 0.000000i,ii,j,ek,tt=2 1 3 -7564.800293 0.000000i,ii,j,ek,tt=2 2 3 0.