有限元编程算例(Fortran)本程序通过Fortran语言编写，程序在Intel Parallel Studio XE 2013 with VS2013中成功运行，程序为《计算力学》（龙述尧等编）一书中的源程序，仅作研究学习使用，省去了敲写的麻烦。

源程序为：!Page149COMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AECOMMON/X2/JM(100,3),NZC(50),CJZ(100,2),PJ(100,2),B(3,6),D(3,3),S(3,6),TKZ(200,20),EKE(6,6),P(200) OPEN(5,FILE='DATAIN')!OPEN(6,FILE='DATAOUT',STATUS='NEW')CALL DATAIF(IND.EQ.0)GOTO 10EO=EO/(1.0-UN*UN)UN=UN/(1.0-UN)10 CALL TOTSTICALL LOADCALL SUPPORCALL SOLVEQCALL STRESSPAUSE!STOPENDSUBROUTINE DATACOMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AECOMMON/X2/JM(100,3),NZC(50),CJZ(100,2),PJ(100,2),B(3,6),D(3,3),S(3,6),TKZ(200,20),EKE(6,6),P(200)READ(5,*)NJ,NE,NZ,NDD,NPJ,INDNJ2=NJ*2NPJ1=NPJ+1READ(5,*)EO,UN,GAMA,TEREAD(5,*)((JM(I,J),J=1,3),I=1,NE)READ(5,*)((CJZ(I,J),J = 1,2),I=1,NJ)!Page150READ(5,*)(NZC(I),I=1,NZ)READ(5,*)((PJ(I,J),J=1,2),I=1,NPJ1)WRITE(6,10)(I,(CJZ(I,J),J=1,2),I=1,NJ)10 FORMA T(4X,2HNO,6X,1HX,6X,1HY/(I6,2X,F7.2,F7.2))RETURNENDSUBROUTINE ELEST(MEO,IASK)COMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AECOMMON/X2/JM(100,3),NZC(50),CJZ(100,2),PJ(100,2),B(3,6),D(3,3),S(3,6),TKZ(200,20),EKE(6,6),P(200) IE=JM(MEO,1)JE=JM(MEO,2)ME=JM(MEO,3)CM=CJZ(JE,1)-CJZ(IE,1)BM=CJZ(IE,2)-CJZ(JE,2)CJ=CJZ(IE,1)-CJZ(ME,1)BJ=CJZ(ME,2)-CJZ(IE,2)AE=(BJ*CM-BM*CJ)/2.0IF(IASK.LE.1) GOTO 50DO 10 I=1,3DO 10 J=1,6B(I,J)=0.010 CONTINUEB(1,1)=-BJ-BMB(1,3)=BJB(1,5)=BMB(2,2)=-CJ-CMB(2,4)=CJB(2,6)=CMB(3,1)=B(2,2)B(3,2)=B(1,1)B(3,3)=B(2,4)B(3,4)=B(1,3)B(3,5)=B(2,6)!Page151B(3,6)=B(1,5)DO 20 I=1,3DO 20 J=1,6B(I,J)=B(I,J)/(2.0*AE) 20 CONTINUED(1,1)=EO/(1.0-UN*UN)D(1,2)=EO*UN/(1.0-UN*UN) D(2,1)=D(1,2)D(2,2)=D(1,1)D(1,3)=0.0D(2,3)=0.0D(3,1)=0.0D(3,2)=0.0D(3,3)=EO/(2.0*(1.0+UN))DO 30 I=1,3DO 30 J=1,6S(I,J)=0.0DO 30 K=1,3S(I,J)=S(I,J)+D(I,K)*B(K,J)30 CONTINUEIF(IASK.LE.2) GOTO 50DO 40 I=1,6DO 40 J=1,6EKE(I,J)=0.0DO 40 K=1,3!**********************************Exchange B And S***********************************************EKE(I,J)=EKE(I,J)+B(K,I)*S(K,J)*AE*TE40 CONTINUE50 CONTINUERETURNENDSUBROUTINE TOTSTICOMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AECOMMON/X2/JM(100,3),NZC(50),CJZ(100,2),PJ(100,2),B(3,6),D(3,3),S(3,6),TKZ(200,20),EKE(6,6),P(200)!Page152DO 20 I=1,NJ2DO 20 J=1,NDDTKZ(I,J)=0.020 CONTINUE!*************Not Understanded*****************************DO 30 MEO=1,NECALL ELEST(MEO,3)DO 30 I=1,3DO 30 II=1,2LH=2*(I-1)+IILDH=2*(JM(MEO,I)-1)+IIDO 30 J=1,3DO 30 JJ=1,2L=2*(J-1)+JJLZ=2*(JM(MEO,J)-1)+JJLD=LZ-LDH+1IF(LD.LE.0) GOTO 30TKZ(LDH,LD)=TKZ(LDH,LD)+EKE(LH,L)30 CONTINUERETURNENDSUBROUTINE LOADCOMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AECOMMON/X2/JM(100,3),NZC(50),CJZ(100,2),PJ(100,2),B(3,6),D(3,3),S(3,6),TKZ(200,20),EKE(6,6),P(200)DO 10 I=1,NJ2P(I)=0.010 CONTINUEIF(NPJ.EQ.0) GOTO 30DO 20 I=1,NPJI1=I+1J=IFIX(PJ(I1,2))P(J)=PJ(I1,1)20 CONTINUE30 IF(GAMA.LE.0.0) GOTO 50!Page153DO 40 MEO=1,NECALL ELEST(MEO,1)PE=-GAMA*AE*TE/3.0IE=JM(MEO,1)JE=JM(MEO,2)ME=JM(MEO,3)P(2*IE)=P(2*IE)+PEP(2*JE)=P(2*JE)+PEP(2*ME)=P(2*ME)+PE40 CONTINUE50 CONTINUERETURNENDSUBROUTINE SUPPORCOMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AECOMMON/X2/JM(100,3),NZC(50),CJZ(100,2),PJ(100,2),B(3,6),D(3,3),S(3,6),TKZ(200,20),EKE(6,6),P(200)DO 60 I=1,NZMZ=NZC(I)TKZ(MZ,1)=1.0DO 10 J=2,NDDTKZ(MZ,J)=0.010 CONTINUEIF(MZ-NDD)20,20,3020 JO=MZGOTO 4030 JO=NDD40 DO 50 J = 2,JOJ1=MZ-JTKZ(J1+1,J)=0.050 CONTINUEP(MZ)=0.060 CONTINUERETURNEND!Page154SUBROUTINE SOLVEQCOMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AECOMMON/X2/JM(100,3),NZC(50),CJZ(100,2),PJ(100,2),B(3,6),D(3,3),S(3,6),TKZ(200,20),EKE(6,6),P(200)NJ1=NJ2-1DO 50 K=1,NJ1IF(NJ2-K-NDD+1)10,10,2010 IM=NJ2GOTO 3020 IM=K+NDD-130 K1=K+1DO 50 I=K1,IML=I-K+1C=TKZ(K,L)/TKZ(K,1)LD1=NDD-L+1DO 40 J=1,LD1M=J+I-KTKZ(I,J)=TKZ(I,J)-C*TKZ(K,M)40 CONTINUEP(I)=P(I)-C*P(K)50 CONTINUEP(NJ2)=P(NJ2)/TKZ(NJ2,1)DO 100 I1 = 1,NJ1I=NJ2-I1!************************************************************************下面一行可能出错IF(NDD-NJ2+I-1)60,60,7060 JO=NDDGOTO 8070 JO=NJ2-I+180 DO 90 J=2,JOLH=J+I-1P(I)=P(I)-TKZ(I,J)*P(LH)90 CONTINUEP(I)=P(I)/TKZ(I,1)100 CONTINUE!Page155WRITE(6,110)(I,P(2*I-1),P(2*I),I=1,NJ)!************************************************************************************ 110 FORMA T(2X,3HJD=,3X,2HU=,12X,2HV=/(I4,3X,F16.7,3X,F16.7))RETURNENDSUBROUTINE STRESSCOMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AECOMMON/X2/JM(100,3),NZC(50),CJZ(100,2),PJ(100,2),B(3,6),D(3,3),S(3,6),TKZ(200,20),EKE(6,6),P(200) DIMENSION WY(6),YL(3)DO 60 MEO=1,NECALL ELEST(MEO,2)DO 10 I=1,3DO 10 J=1,2LH=2*(I-1)+JLDH=2*(JM(MEO,I)-1)+JWY(LH)=P(LDH)10 CONTINUEDO 20 I=1,3YL(I)=0.0DO 20 J=1,6YL(I)=YL(I)+S(I,J)*WY(J)20 CONTINUESIGX=YL(1)SIGY=YL(2)TOXY=YL(3)PYL=(SIGX+SIGY)/2.0SIG=(SIGX-SIGY)**2/4.0+TOXY*TOXYRYL=SQRT(SIG)SIG1=PYL+RYLSIG2=PYL-RYLIF(SIGY.EQ.SIG2) GOTO 30CETA1=TOXY/(SIGY-SIG2)CETA=90.0-57.29578*ATAN(CETA1)GOTO 40!Page15630 CETA=0.040 WRITE(6,50)MEO,SIGX,SIGY,TOXY,SIG1,SIG2,CETA50FORMA T(4X,2HE=,I3/2X,3HSX=,F11.3,3X,3HSY=,F11.3,3X,4HTAU=,F11.3/2X,3HS1=,F11.3,3X,3HS2=,F11. 3,3X,4HCET=,F11.3)!50FORMA T(4X,2HE=,I3/2X,3HSX=,Fll.3,3X,3HSY=,F11.3,3X,4HTAU=,F11.3/2X,3HSl=,Fll.3,3X,3HS2=,F11.3,3 X,4HCET=,F11.3)60 CONTINUERETURNEND输入文件为datain28,36,9,10,4,01,0.17,0,11,5,22,5,62,6,33,6,73,7,44,7,85,9,66,9,106,10,77,10,117,11,88,11,129,13,1010,13,1410,14,11 11,14,15 11,15,12 12,15,16 13,17,14 14,17,18 14,18,15 15,18,19 15,19,16 16,19,20 17,21,18 18,21,22 18,22,19 19,22,23 19,23,20 20,23,24 21,25,22 22,25,26 22,26,23 23,26,27 23,27,24 24,27,28 0,61,62,63,60,51,52,53,50,41,42,43,40,31,32,33,30,21,22,23,20,11,12,13,10,01,02,03,07,15,23,31,39,47,49,50,550,0-5E4,2-10E4,4-10E4,6-5E4,8输出结果为：DATAOUTNO X Y1 0.00 6.002 1.00 6.003 2.00 6.004 3.00 6.005 0.00 5.006 1.00 5.007 2.00 5.008 3.00 5.009 0.00 4.0010 1.00 4.0011 2.00 4.0012 3.00 4.0013 0.00 3.0014 1.00 3.0015 2.00 3.0016 3.00 3.0017 0.00 2.0018 1.00 2.0019 2.00 2.0020 3.00 2.0021 0.00 1.0022 1.00 1.0023 2.00 1.0024 3.00 1.0025 0.00 0.0026 1.00 0.0027 2.00 0.0028 3.00 0.00JD= U= V=1 -29766.873 -1173917.7502 -14003.185 -1174018.8753 -3753.270 -1179518.1254 0.000 -1181719.7505 -26382.471 -1072681.5006 -10746.993 -1073615.0007 -2064.593 -1082360.7508 0.000 -1085873.2509 -13536.995 -964010.12510 3372.794 -970055.12511 7268.415 -989269.12512 0.000 -998401.81213 7816.581 -835383.43814 27176.234 -861713.93815 22063.230 -905726.12516 0.000 -927165.18817 29514.479 -665602.87518 53419.637 -747340.43819 34876.832 -839806.81220 0.000 -881219.12521 29580.273 -416288.71922 52944.918 -632601.12523 17504.195 -803765.68824 0.000 -859481.93825 0.000 0.00026 -120102.820 -583505.37527 -76202.375 -787347.18828 0.000 -829170.812E= 1SX= -1489.530 SY=-101489.383 TAU= -1489.531 S1= -1467.348 S2=-101511.562 CET= 179.147 E= 2SX= -1475.844 SY=-100654.875 TAU= -1790.500 S1= -1443.531 S2=-100687.188 CET= 178.966 E= 3SX= -7021.670 SY=-101597.672 TAU= -3741.688 S1= -6873.875 S2=-101745.469 CET= 177.738 E= 4SX= -8067.500 SY= -98528.750 TAU= -4459.156 S1= -7848.227 S2= -98748.023 CET= 177.185 E= 5SX= -13143.328 SY= -99391.750 TAU= -1662.500 S1= -13111.293 S2= -99423.781 CET= 178.896 E= 6SX= -14652.781 SY= -98337.500 TAU= -1501.062S1= -14625.867 S2= -98364.414 CET= 178.973 E= 7SX= -2923.122 SY=-109168.297 TAU= -5888.469 S1= -2597.762 S2=-109493.656 CET= 176.837 E= 8SX= -716.078 SY=-103681.562 TAU= -8617.406 S1= 0.148 S2=-104397.789 CET= 175.249 E= 9SX= -9188.316 SY=-105121.867 TAU= -9771.594 S1= -8203.125 S2=-106107.062 CET= 174.243 E= 10SX= -12285.000 SY= -95180.250 TAU= -12199.594 S1= -10526.887 S2= -96938.359 CET= 171.799 E= 11SX= -14170.516 SY= -95500.750 TAU= -5489.531 S1= -13801.664 S2= -95869.602 CET= 176.156 E= 12SX= -22797.406 SY= -91347.000 TAU= -3902.844 S1= -22575.914 S2= -91568.492 CET= 176.752 E= 13SX= -5104.269 SY=-129494.438 TAU= -11708.750 S1= -4011.727 S2=-130586.977 CET= 174.669 E= 14SX= 969.672 SY=-108176.375 TAU= -21424.750 S1= 5024.582 S2=-112231.281 CET= 169.283 E= 15SX= -14954.572 SY=-110883.469 TAU= -18383.531 S1= -11552.273 S2=-114285.766 CET= 169.515 E= 16SX= -19890.141 SY= -86924.312 TAU= -25131.188 S1= -11514.844 S2= -95299.609 CET= 161.569 E= 17SX= -22109.688 SY= -87301.625 TAU= -10225.406 S1= -20543.453 S2= -88867.859 CET= 171.292 E= 18SX= -35190.453 SY= -77219.000 TAU= -9162.000 S1= -33280.023 S2= -79129.430 CET= 168.222 E= 19SX= -9785.850 SY=-171444.172 TAU= -20524.969 S1= -7220.594 S2=-174009.422 CET= 172.876 E= 20SX= 4594.438 SY=-113592.375 TAU= -46145.688 S1= 20477.398 S2=-129475.336 CET= 161.007 E= 21SX= -25287.307 SY=-118672.312 TAU= -30023.750 S1= -16467.512 S2=-127492.109 CET= 163.629 E= 22SX= -30634.422 SY= -71127.188 TAU= -44991.469 S1= -1543.715 S2=-100217.891 CET= 147.114 E= 23SX= -34259.609 SY= -71743.438 TAU= -14637.906 S1= -29220.699 S2= -76782.344 CET= 161.005 E= 24SX= -43958.047 SY= -53418.938 TAU= -17697.562 S1= -30369.627 S2= -67007.359 CET= 142.482 E= 25SX= -19028.160 SY=-252549.000 TAU= -34958.688 S1= -13907.055 S2=-257670.094 CET= 171.666 E= 26SX= 3973.812 SY=-114063.750 TAU= -92238.344 S1= 54459.047 S2=-164548.984 CET= 151.307 E= 27SX= -39180.809 SY=-121400.055 TAU= -39312.688 S1= -23409.074 S2=-137171.781 CET= 158.140 E= 28SX= -42804.766 SY= -43317.938 TAU= -65723.062 S1= 22662.211 S2=-108784.914 CET= 135.112 E= 29SX= -42224.094 SY= -43219.188 TAU= -10273.375 S1= -32436.225 S2= -53007.055 CET= 136.386 E= 30SX= -21830.422 SY= -25448.312 TAU= -23810.344 S1= 239.594 S2= -47518.328 CET= 137.172 E= 31SX= -48815.199 SY=-424587.344 TAU= -79800.078 S1= -32570.844 S2=-440831.688 CET= 168.494 E= 32SX=-132271.750 SY= -71582.000 TAU=-175409.250 S1= 76087.781 S2=-279941.531 CET= 130.093 E= 33SX= -45090.102 SY= -56761.105 TAU= 804.781 S1= -45034.867 S2= -56816.336 CET= 3.926 E= 34SX= 42332.711 SY= -9221.938 TAU= -47066.344 S1= 70218.328 S2= -37107.555 CET= 149.354 E= 35SX= -20899.344 SY= -19971.375 TAU= 16235.219 S1= -4193.512 S2= -36677.207 CET= 45.819E= 36SX= 73163.914 SY= -17873.250 TAU= -17873.344 S1= 76547.250 S2= -21256.586 CET= 169.281。