有限元分析程序与算例
I=ILJ(M)*3 B(I-2)=PX(M) B(I-1)=PY(M) B(I)=PM(M) 100 CONTINUE RETURN END
! Initialize terminal forces of members SUBROUTINE F0 (NLM,NM,F)
IMPLICIT DOUBLE PRECISION (A-H,O-Z) DIMENSION F(6,NM) WRITE (*,1) NLM
1 FORMAT(/10X,'The Loadings at Members'//17X,'NLM=',I4) DO 200 J=1,NM DO 100 I=1,6 F(I,J)=0.0 100 CONTINUE 200 CONTINUE RETURN END
! Read data of loading at member & calculate fixed-end forces SUBROUTINE IOLMFB (NM,NJ,N,NLM,ILM,ITL,PV,DST,IST,IEN,X,Y,LV,F,B) IMPLICIT DOUBLE PRECISION (A-H,O-Z) DIMENSION
ILM(NLM),ITL(NLM),PV(NLM),DST(NLM),IST(NM),IEN(NM),X(NJ),& Y(NJ),LV(6,NM),F(6,NM),B(N) READ (3,*)(ILM(M),ITL(M),PV(M),DST(M),M=1,NLM) WRITE (*,1)
1 FORMAT(/2X,'ILM',2X,'ITL',11X,'PV',12X,'DST') WRITE (*,2)(ILM(M),ITL(M),PV(M),DST(M),M=1,NLM) 2 FORMAT(1X,I4,I5,F16.4,F16.6) DO 100 M=1,NLM L=ILM(M) CALL RLCS (L,NM,NJ,IST,IEN,X,Y,RL,C,S) D1=DST(M) D2=RL-D1 IF (ITL(M).EQ.1.OR.ITL(M).EQ.3) THEN P1=PV(M)*C P2=-PV(M)*S END IF IF (ITL(M).EQ.2.OR.ITL(M).EQ.4) THEN P1=PV(M)*S P2=PV(M)*C END IF
- 8 -
有限元分析程序与算例
IF (ITL(M).EQ.1.OR.ITL(M).EQ.2) THEN F1=-P1*D2/RL F4=-P1-F1 F2=-P2*D2*D2*(RL+2.0*D1)/(RL*RL*RL) F5=-P2-F2 F3=-P2*D1*D2*D2/(RL*RL) F6=P2*D1*D1*D2/(RL*RL) END IF IF (ITL(M).EQ.3.OR.ITL(M).EQ.4) THEN G=P2*D1*D1/(12.0*RL*RL) F3=-G*((6.0*RL-8.0*D1)*RL+3.0*D1*D1) F6=G*D1*(4.0*RL-3.0*D1) F5=-6.0*G*D1*(2.0-D1/RL) F2=-P2*D1-F5 F4=-P1*D1*D1/(2.0*RL) F1=-P1*D1-F4 END IF F(1,L)=F(1,L)+F1 F(2,L)=F(2,L)+F2 F(3,L)=F(3,L)+F3 F(4,L)=F(4,L)+F4 F(5,L)=F(5,L)+F5 F(6,L)=F(6,L)+F6 B(LV(1,L))=B(LV(1,L))-F1*C+F2*S B(LV(2,L))=B(LV(2,L))-F1*S-F2*C B(LV(3,L))=B(LV(3,L))-F3 B(LV(4,L))=B(LV(4,L))-F4*C+F5*S B(LV(5,L))=B(LV(5,L))-F4*S-F5*C B(LV(6,L))=B(LV(6,L))-F6 100 CONTINUE RETURN END
! Introduce support conditions into load vector B SUBROUTINE BS (NS,N,IS,VS,B)
IMPLICIT DOUBLE PRECISION (A-H,O-Z) DIMENSION IS(NS),VS(NS),B(N) DO 100 M=1,NS I=3*(IS(M)/10)-3+MOD(IS(M),10) B(I)=VS(M)*1D22 100 CONTINUE RETURN END
! print joint displacements SUBROUTINE OJD (NJ,N,B)
- 9 -
有限元分析程序与算例
IMPLICIT DOUBLE PRECISION (A-H,O-Z) DOUBLE PRECISION B(N) WRITE (*,1)
1 FORMAT(/13X,'The Results of Calculation' //10X,'The Joint Displacements'& //1X,'joint',8X,'u',18X,'v',14X,'phi') WRITE (*,2)(M,B(3*M-2),B(3*M-1),B(3*M),M=1,NJ) 2 FORMAT(1X,I4,1P3E18.8) RETURN END
! Calculate & print terminal forces of members SUBROUTINE COTF (E,NM,NJ,N,IST,IEN,AR,RI,X,Y,LV,B,F) IMPLICIT DOUBLE PRECISION (A-H,O-Z) DIMENSION
IST(NM),IEN(NM),AR(NM),RI(NM),X(NJ),Y(NJ),LV(6,NM),F(6,NM),B(N) DOUBLE PRECISION U1,U2,U3,U4,U5,U6 WRITE (*,1)
1 FORMAT(/10X,'The Terminal Forces'//1X,'member',4X,'N(st)',6X,'Q(st)',7X,& 'M(st)', 6X,'N(en)',6X,'Q(en)',7X,'M(en)') DO 100 M=1,NM CALL KEBAR (M,E,NM,NJ,IST,IEN,AR,RI,X,Y,C,S,E1,E2,E3,E4) U1=B(LV(1,M))*C+B(LV(2,M))*S U2=-B(LV(1,M))*S+B(LV(2,M))*C U3=B(LV(3,M)) U4=B(LV(4,M))*C+B(LV(5,M))*S U5=-B(LV(4,M))*S+B(LV(5,M))*C U6=B(LV(6,M)) F(1,M)=F(1,M)+E1*(U1-U4) F(2,M)=F(2,M)+E2*(U2-U5)+E3*(U3+U6) F(3,M)=F(3,M)+E3*(U2-U5)+E4*(U3+0.5*U6) F(4,M)=F(4,M)+E1*(U4-U1) F(5,M)=F(5,M)+E2*(U5-U2)-E3*(U3+U6) F(6,M)=F(6,M)+E3*(U2-U5)+E4*(0.5*U3+U6) WRITE (*,2)M,F(1,M),F(2,M),F(3,M),F(4,M),F(5,M),F(6,M) 2 FORMAT(1X,I4,2(2F12.3,F12.3)) 100 CONTINUE RETURN END
二、数据输入：
********************************** * * * PLANE FRAME EX.1 2007.6.14 * * *
- 10 -
有限元分析程序与算例
********************************** 2.1E8 2 3 3 1
1 2 102E-4 32240E-8 2 3 102E-4 32240E-8 0 0 4 3 8 3 11 0 12 -0.1 32 0 1
2 0 -20 0 2
2 2 -40 2 2 4 -10 4
三、结果输出：
********************************** * * * PLANE FRAME EX.1 2007.6.14 * * * **********************************
The Input Data
The General Information
E NM NJ NS NLC 2.100E+08 2 3 3 1
The Information of Members
member start end A I
1 1 2 1.020000E-02 3.224000E-04 2 2 3 1.020000E-02 3.224000E-04
The Joint Coordinates
joint X Y
1 .000000 .000000 2 4.000000 3.000000 3 8.000000 3.000000
- 11 -
有限元分析程序与算例
The Information of Supports
IS VS
11 .000000 12 -.100000 32 .000000
( NA= 36 ) ( NW= 179 ) Loading Case 1
The Loadings at Joints
NLJ= 1
ILJ PX PY PM 2 .0000 -20.0000 .00000
The Loadings at Members
NLM= 2
ILM ITL PV DST 2 2 -40.0000 2.000000 2 4 -10.0000 4.000000
The Results of Calculation
The Joint Displacements
joint u v phi
1 1.14352090E-34 -1.00000000E-01 7.87516325E-03 2 -2.80901559E-02 -6.26164867E-02 1.23062152E-02 3 -2.80901559E-02 -3.00000175E-21 1.78204134E-02
The Terminal Forces
member N(st) Q(st) M(st) N(en) Q(en) 1 18.000 24.000 .000 -18.000 -24.000 2 .000 10.000 -120.000 .000 70.000
( NA= 36 ) ( NW= 187 )
- 12 -
M(en) 120.000 .000