有限元编程算例.docx
- 文档编号:7821267
- 上传时间:2023-01-26
- 格式:DOCX
- 页数:20
- 大小:146.38KB
有限元编程算例.docx
《有限元编程算例.docx》由会员分享,可在线阅读,更多相关《有限元编程算例.docx(20页珍藏版)》请在冰豆网上搜索。
有限元编程算例
有限元编程算例(Fortran)
本程序通过Fortran语言编写,程序在IntelParallelStudioXE2013with
VS2013中成功运行,程序为《计算力学》(龙述尧等编)一书中的源程序,仅作
研究学习使用,省去了敲写的麻烦。
3.7.4算例
例工夕设谦聲祇受沟布戴荷、如图3-3叙叮所示「假進£=」,泊松比世二①17•不计容重+厚度才=1为平面应力问題,因对称取半边结构计真•結枸支承,跟兀划分,节点绸号如图3,16(b)JPfziip试矗出及y=Sm截啲的竖向位移图応=总mflg面的靳应力分布图,
⑴(山
图玉開燮均布栽荷的简支漂梁
源程序为:
!
Page149
C0MM0N/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,E0,UN,GAMA,TE,AE
COMMON/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)
0PEN(5,FILE='DATAIN')
!
OPEN(6,FILE='DATAOUT',STATUS='NEW')
CALLDATA
IF(IND.EQ.0)GOTO10
EO=EO/(1.0-UN*UN)
UN=UN/(1.0-UN)
10CALLTOTSTI
CALLLOAD
CALLSUPPOR
CALLSOLVEQ
CALLSTRESS
PAUSE
!
STOP
END
SUBROUTINEDATA
COMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AE
COMMON/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,IND
NJ2=NJ*2
NPJ1=NPJ+1
READ(5,*)EO,UN,GAMA,TE
READ(5,*)((JM(I,J),J=1,3),I=1,NE)
READ(5,*)((CJZ(I,J),J=1,2),I=1,NJ)
!
Page150
READ(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)
10FORMAT(4X,2HNO,6X,1HX,6X,1HY/(I6,2X,F7.2,F7.2))
RETURN
END
SUBROUTINEELEST(MEO,IASK)
COMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AE
COMMON/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.0
IF(IASK.LE.1)GOTO50
DO10I=1,3
DO10J=1,6
B(I,J)=0.0
10CONTINUE
B(1,1)=-BJ-BM
B(1,3)=BJ
B(1,5)=BM
B(2,2)=-CJ-CM
B(2,4)=CJ
B(2,6)=CM
B(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)
!
Page151
B(3,6)=B(1,5)
DO20I=1,3
DO20J=1,6
B(I,J)=B(I,J)/(2.0*AE)
20CONTINUE
D(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.0
D(2,3)=0.0
D(3,1)=0.0
D(3,2)=0.0
D(3,3)=EO/(2.0*(1.0+UN))
DO30I=1,3
DO30J=1,6
S(I,J)=0.0
DO30K=1,3
S(I,J)=S(I,J)+D(I,K)*B(K,J)
30CONTINUE
IF(IASK.LE.2)GOTO50
DO40I=1,6
DO40J=1,6
EKE(I,J)=0.0
DO40K=1,3
And
!
**********************************ExchangeB
S***********************************************
EKE(I,J)=EKE(I,J)+B(K,I)*S(K,J)*AE*TE
40CONTINUE
50CONTINUE
RETURN
END
SUBROUTINETOTSTI
COMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AE
COMMON/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)!
Page152
DO20I=1,NJ2
DO20J=1,NDD
TKZ(I,J)=0.0
20CONTINUE
!
*************NotUnderstanded*****************************
DO30MEO=1,NE
CALLELEST(MEO,3)
DO30I=1,3
DO30II=1,2
LH=2*(I-1)+II
LDH=2*(JM(MEO,I)-1)+II
DO30J=1,3
DO30JJ=1,2
L=2*(J-1)+JJ
LZ=2*(JM(MEO,J)-1)+JJ
LD=LZ-LDH+1
IF(LD.LE.0)GOTO30
TKZ(LDH,LD)=TKZ(LDH,LD)+EKE(LH,L)
30CONTINUE
RETURN
END
SUBROUTINELOAD
COMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AE
COMMON/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)
DO10I=1,NJ2
P(I)=0.0
10CONTINUE
IF(NPJ.EQ.0)GOTO30
DO20I=1,NPJ
I1=I+1
J=IFIX(PJ(I1,2))
P(J)=PJ(I1,1)
20CONTINUE
30IF(GAMA.LE.0.0)GOTO50
!
Page153
DO40MEO=1,NE
CALLELEST(MEO,1)
PE=-GAMA*AE*TE/3.0
IE=JM(MEO,1)
JE=JM(MEO,2)
ME=JM(MEO,3)
P(2*IE)=P(2*IE)+PE
P(2*JE)=P(2*JE)+PE
P(2*ME)=P(2*ME)+PE
40CONTINUE50CONTINUE
RETURN
END
SUBROUTINESUPPOR
COMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AE
COMMON/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)
DO60I=1,NZ
MZ=NZC(I)
TKZ(MZ,1)=1.0
DO10J=2,NDD
TKZ(MZ,J)=0.0
10CONTINUE
IF(MZ-NDD)20,20,30
20JO=MZ
GOTO40
30JO=NDD
40DO50J=2,JO
J1=MZ-J
TKZ(J1+1,J)=0.0
50CONTINUE
P(MZ)=0.0
60CONTINUE
RETURN
END
!
Page154
SUBROUTINESOLVEQ
COMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AE
COMMON/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-1
DO50K=1,NJ1
IF(NJ2-K-NDD+1)10,10,20
10IM=NJ2
GOTO30
20IM=K+NDD-1
30K1=K+1
DO50I=K1,IM
L=I-K+1
C=TKZ(K,L)/TKZ(K,1)
LD1=NDD-L+1
DO40J=1,LD1M=J+I-KTKZ(I,J)=TKZ(I,J)-C*TKZ(K,M)
40CONTINUE
P(I)=P(I)-C*P(K)
50CONTINUE
P(NJ2)=P(NJ2)/TKZ(NJ2,1)
DO100I1=1,NJ1
下面一行可能出错
I=NJ2-I1!
************************************************************************
IF(NDD-NJ2+I-1)60,60,7060JO=NDD
GOTO80
70JO=NJ2-I+1
80DO90J=2,JO
LH=J+I-1
P(I)=P(I)-TKZ(I,J)*P(LH)
90CONTINUE
P(I)=P(I)/TKZ(I,1)
100CONTINUE
!
Page155
WRITE(6,110)(I,P(2*I-1),P(2*I),I=1,NJ)
!
************************************************************************************
110FORMAT(2X,3HJD=,3X,2HU=,12X,2HV=/(I4,3X,F16.7,3X,F16.7))
RETURN
END
SUBROUTINESTRESS
COMMON/X1/NJ,NE,NZ,NDD,NPJ,IND,NJ2,EO,UN,GAMA,TE,AE
DIMENSIONWY(6),YL(3)
DO60MEO=1,NE
CALLELEST(MEO,2)
DO10I=1,3
DO10J=1,2
LH=2*(I-1)+J
LDH=2*(JM(MEO,I)-1)+J
WY(LH)=P(LDH)
10CONTINUE
DO20I=1,3
YL(I)=0.0
DO20J=1,6
YL(I)=YL(I)+S(I,J)*WY(J)
20CONTINUE
SIGX=YL
(1)
SIGY=YL
(2)
TOXY=YL(3)
PYL=(SIGX+SIGY)/2.0
SIG=(SIGX-SIGY)**2/4.0+TOXY*TOXY
RYL=SQRT(SIG)
SIG1=PYL+RYL
SIG2=PYL-RYL
IF(SIGY.EQ.SIG2)GOTO30
CETA1=TOXY/(SIGY-SIG2)
CETA=90.0-57.29578*ATAN(CETA1)
GOTO40
!
Page156
30CETA=0.0
40WRITE(6,50)MEO,SIGX,SIGY,TOXY,SIG1,SIG2,CETA
50
FORMAT(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)
!
50
FORMAT(4X,2HE=,I3/2X,3HSX=,Fll.3,3X,3HSY=,F11.3,3X,4HTAU=,F11.3/2X,3HSl=,Fll.3,3X,3HS2=,F11.3,3X,4HCET=,F11.3)
60CONTINUE
RETURN
END
输入文件为datain
28,36,9,10,4,0
1,0.17,0,1
1,5,2
2,5,6
2,6,3
3,6,7
3,7,4
4,7,8
5,9,6
6,9,10
6,10,7
7,10,11
7,11,8
8,11,12
9,13,10
10,13,14
10,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,6
1,6
2,6
3,6
0,5
1,5
2,5
3,5
0,4
1,4
2,4
3,4
0,3
1,3
2,3
3,3
0,2
1,2
2,2
3,2
0,1
1,1
2,1
3,1
0,0
1,0
2,0
3,0
7,15,23,31,39,47,49,50,55
0,0
-5E4,2
-10E4,4
-10E4,6
-5E4,8
输出结果为:
DATAOUT
NO
X
Y
1
0.00
6.00
2
1.00
6.00
3
2.00
6.00
4
3.00
6.00
5
0.00
5.00
6
1.00
5.00
7
2.00
5.00
8
3.00
5.00
9
0.00
4.00
10
1.00
4.00
11
2.00
4.00
12
3.00
4.00
13
0.00
3.00
14
1.00
3.00
15
2.00
3.00
16
3.00
3.00
17
0.00
2.00
18
1.00
2.00
19
2.00
2.00
20
3.00
2.00
21
0.00
1.00
22
1.00
1.00
23
2.00
1.00
24
3.00
1.00
25
0.00
0.00
26
1.00
0.00
27
2.00
0.00
28
3.00
0.00
JD=U=V=
1-29766.873-1173917.750
2
-14003.185
-1174018.875
3
-3753.270
-1179518.125
4
0.000
-1181719.750
5
-26382.471
-1072681.500
6
-10746.993
-1073615.000
7
-2064.593
-1082360.750
8
0.000
-1085873.250
9
-13536.995
-964010.125
10
3372.794
-970055.125
11
7268.415
-989269.125
12
0.000
-998401.812
13
7816.581
-835383.438
14
27176.234
-861713.938
15
22063.230
-905726.125
16
0.000
-927165.188
17
29514.479
-665602.875
18
53419.637
-747340.438
19
34876.832
-839806.812
20
0.000
-881219.125
21
29580.273
-416288.719
22
52944.918
-632601.125
23
17504.195
-803765.688
24
0.000
-859481.938
25
0.000
0.000
26
-120102.820
-583505.375
27
-76202.375
-787347.188
28
0.000
-829170.812
E=1
SX=
-1489.530
SY=-101489.383
TAU=
-1489.531
S1=
-1467.348
S2=-101511.562
CET=
179.147
E=
2
SX=
-1475.844
SY=-100654.875
TAU=
-1790.500
S1=
-1443.531
S2=-100687.188
CET=
178.966
E=
3
SX=
-7021.670
SY=-101597.672
TAU=
-3741.688
S1=
-6873.875
S2=-101745.469
CET=
177.738
E=
4
SX=
-8067.500
SY=-98528.750
TAU=
-4459.156
S1=
-7848.227
S2=-98748.023
CET=
177.185
E=
5
SX=-
13143.328
SY=-99391.750
TAU=
-1662.500
S1=-13111.293
S2=-99423.781
CET=
178.896
E=
6
SX=-14652.781
SY=-98337.500
TAU=
-1501.062
S1=-14625.867
S2=-98364.414
CET=
178.973
E=7
SX=-2923.122
SY=-109168.297
TAU=
-5888.469
S1=-2597.762
S2=-109493.656
CET=
176.837
E=8
SX=-716.078
SY=-103681.562
TAU=
-8617.406
S1=0.148
S2=-104397.789
CET=
175.249
E=9
SX=-9188.316
SY=-105121.867
TAU=
-9771.594
S1=-8203.125
S2=-106107.062
CET=
174.243
E=10
SX=-12285.000
SY=-95180.250
TAU=-12199.594
S1=-10526.887
S2=-96938.359
CET=
171.799
E=11
SX=-14170.516
SY=-95500.750
TAU=
-5489.531
S1=-13801.664
S2=-95869.602
CET=
176.156
E=12
SX=-22797.406
SY=-91347.000
TAU=
-3902.844
S1=-22575.914
S2=-91568.492
CET=
176.752
E=13
SX=-5104.269
SY=-129494.438
TAU=
-11708.750
S1=-4011.727
S2=-130586.977
CET=
174.669
E=14
SX=969.672
SY=-108176.375
TAU=-21424.750
S1=5024.582
S2=-112231.281
CET=
169.283
E=15
SX=-14954.572
SY=-110883.469
TAU=-
18383.531
S1=-11552.273
S2=-114285.766
CET=
169.515
E=16
SX=-19890.141
SY=-86924.312
TAU=-25131.188
S1=-11514.844
S2=-95299.609
CET=
161.569
E=17
SX=-22109.688
SY=-87301.625
TAU=-
10225.406
S1=-20543.453
S2=-88867.859
CET=
171.292
E=18
SX=-35190.453
SY=-77219.000
TAU=
-9162.000
S1=-33280.023
S2=-79129.430
CET=
168.222
E=19
SX=-9785.850
SY=-171444.172
TAU=
-20524.969
S1=-7220.594
S2=-174009.422
CET=
172.876
E=20
SX=4594.438
SY=-113592.375
TAU=-46145.688
S1=20477.398
S2=-129475.336
CET=
161.007
E=21
SX=-25287.307
SY=-118672.312
TAU=
-30023.750
S1=-16467.512
S2=-127492.109
CET=
163.629
E=22
SX=-30634.422
SY=-71127.188
TAU=-44991.469
S1=-1543.715
S2=-100217.891
CET=
147.114
E=23
SX=-34259.609
SY=-71743.438
TAU=-14637.906
S1=-29220.699
S2=-76782.344
CET=
161.005
E=24
SX=-43958
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 有限元 编程