有限元分析程序与算例
RF(I)=0.0 DO 30 J=1,8 SKE(I,J)=0.0
30 CONTINUE 40 CONTINUE WRITE(*,*)'NE=',NEE DO 50 I=1,4 JB=NOD(I) DO 50 M=1,2 JJ=2*(I-1)+M NN(JJ)=JR(M,JB) 50 CONTINUE L=N
DO 80 I=1,8
IF(NN(I).EQ.0) GOTO 80 IF(NN(I).LT.L) L=NN(I) 80 CONTINUE
DO 85 M=1,8 JP=NN(M)
IF(JP.EQ.0) GOTO 85
IF(JP-L+1.GT.MA(JP)) MA(JP)=JP-L+1 85 CONTINUE
CALL STIF(AE)
IF(NK.EQ.0) GOTO 100 DO 95 IW=1,NK
READ(15,*) NNN,NFACE WRITE(16,650) NNN,NFACE Y0=WG(1,NNN) GAMA=WG(2,NNN) DO 90 JW=1,4 ND=NFACE(JW)
IF(ND.EQ.0) GOTO 90
CALL SURFOR (JW,AE,PR,Y0,GAMA,1) 90 CONTINUE 95 CONTINUE 100 CONTINUE
IF(NSF.EQ.0) GOTO 200 DO 110 IW=1,NSF READ(15,*) ND,PR WRITE(16,750) ND,PR
CALL SURFOR(ND,AE,PR,Y0,GAMA,0)
- 6 -
有限元分析程序与算例
110 CONTINUE 200 CONTINUE
IF(NST.EQ.0) GOTO 300 DO 210 IW=1,NST READ(15,*) ND,PR WRITE(16,750) ND,PR
CALL SURNST (ND,AE,PR) 210 CONTINUE 300 CONTINUE
!C WRITE (16,800) (RF(I),I=1,8) WRITE (12) NN,SKE
WRITE (7) NEE,NME,NOD,NN,XY
600 FORMAT(3X,' NEE= NME= NET= NK= NSF= NST= NOD='/3X,6I5,4I5) 650 FORMAT(3X,' SURFACE NEWS NNN=',I3,2X,4I1)
750 FORMAT(3X,' SURFACE NEWS ND=',I3,5X, 'PR=',2F12.2) 800 FORMAT(3X,'RF=',/(1X,8F10.4)) RETURN END
!C****************************************************************** SUBROUTINE STIF(AE)
!C IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION AE(4,*),RJX(2,2),BV(8),D(4) COMMON /C1/ NP,NE,NM,NR,NW,NF,NI,NDP COMMON /C2/ N,MX,NH
COMMON /C3/ SKE(8,8),NN(8),RF(8),B(3,8),XY(2,4) COMMON /C4/ NEE,NME,NET,NK,NSF,NST COMMON /C5/ FUN(4),P(2,4),XJR(2,2) COMMON /GAUSS/ RSTG(2)
E=AE(1,NME) U=AE(2,NME)
GAMMA=AE(3,NME) T=AE(4,NME)
D1=E*(1.0-U)/((1.0+U)*(1.0-2.0*U)) D2=E*U/((1.0+U)*(1.0-2.0*U)) D3=E*0.5/(1.0+U)
DO 99 IS=1,2 S=RSTG(IS) DO 98 IR=1,2
- 7 -
有限元分析程序与算例
R=RSTG(IR)
CALL RMSD (DET,R,S,RJX) IF(NET.EQ.0) GOTO 40 DO 30 I=1,4 J=2*I
RF(J)=RF(J)-FUN(I)*DET*GAMMA*T 30 CONTINUE 40 K2=0 DO 50 I=1,4 K2=K2+2 K1=K2-1
BV(K1)=B(1,K1) BV(K2)=B(2,K2) 50 CONTINUE DO 70 I=1,8 DO 60 J=I,8
SKE(I,J)=SKE(I,J)+BV(I)*BV(J)*DET*T 60 CONTINUE 70 CONTINUE 98 CONTINUE 99 CONTINUE
DO 104 I=2,8 M=I-1
DO 102 J=1,M SKE(I,J)=SKE(J,I) 102 CONTINUE 104 CONTINUE KK=-2
DO 120 I=1,4 KK=KK+2 K1=KK+1 K2=K1+1 DO 115 J=1,4 LL=(J-1)*2 L1=LL+1 L2=L1+1 IC=0
DO 110 K=1,2 M=K+KK DO 105 L=1,2 NL=L+LL IC=IC+1
- 8 -
有限元分析程序与算例
D(IC)=SKE(M,NL) 105 CONTINUE 110 CONTINUE
SKE(K1,L1)=D(1)*D1+D(4)*D3 SKE(K2,L2)=D(4)*D1+D(1)*D3 SKE(K1,L2)=D(2)*D2+D(3)*D3 SKE(K2,L1)=D(3)*D2+D(2)*D3 115 CONTINUE 120 CONTINUE
!C WRITE (16,500) ((SKE(I,J),I=1,4),J=1,4) 500 FORMAT (10X,'SKE='/(1X,8F10.2)) RETURN END
!C****************************************************************** SUBROUTINE RMSD(DET,R,S,RJX)
!C IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION RJX(2,2)
COMMON /C4/ NEE,NME,NET,NK,NSF,NST COMMON /C5/ FUN(4),P(2,4),XJR(2,2)
COMMON /C3/ SKE(8,8),NN(8),RF(8),B(3,8),XY(2,4)
CALL FPJD(R,S,DET)
REC=1.0/DET
RJX(1,1)=REC*XJR(2,2) RJX(2,1)=-REC*XJR(2,1) RJX(1,2)=-REC*XJR(1,2) RJX(2,2)=REC*XJR(1,1)
K2=0
DO 30 K=1,4 K2=K2+2 K1=K2-1 DO 20 L=1,2 B(L,K1)=0.0 B(L,K2)=0.0
20 CONTINUE DO 25 I=1,2
B(1,K1)=B(1,K1)+RJX(1,I)*P(I,K) B(2,K2)=B(2,K2)+RJX(2,I)*P(I,K) 25 CONTINUE B(3,K1)=B(2,K2) B(3,K2)=B(1,K1)
- 9 -
有限元分析程序与算例
30 CONTINUE
RETURN END
!C****************************************************************** SUBROUTINE FPJD(R,S,DET)
!C IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION XI(4),ETA(4)
COMMON /C4/ NEE,NME,NET,NK,NSF,NST COMMON /C5/ FUN(4),P(2,4),XJR(2,2)
COMMON /C3/ SKE(8,8),NN(8),RF(8),B(3,8),XY(2,4) DATA XI/-1.0,1.0,1.0,-1./ DATA ETA/-1.0,-1.0,1.0,1.0/ DO 10 I=1,4
FUN(I)=0.25*(1.0+XI(I)*R)*(1.0+ETA(I)*S) 10 CONTINUE
DO 50 I=1,4
P(1,I)=0.25*XI(I)*(1.0+ETA(I)*S) P(2,I)=0.25*ETA(I)*(1.0+XI(I)*R) 50 CONTINUE
DO 80 I=1,2 DO 75 J=1,2 DET=0.
DO 70 K=1,4
DET=DET+P(I,K)*XY(J,K) 70 CONTINUE XJR(I,J)=DET
75 CONTINUE 80 CONTINUE
DET=XJR(1,1)*XJR(2,2)-XJR(1,2)*XJR(2,1) !C WRITE (0,500) DET IF (DET.LT.1.0D-5) GOTO 100 RETURN
500 FORMAT (10X,'DET=',F15.3) 100 WRITE(16,600) NEE,R,S,DET
600 FORMAT(1X,'ERROR OF NEGTIVE OR ZERO', & 'JACOBIAN DETERMINANT CONPUTED FOR',& 'ELEMENT (',I5,')'/5X,'R=S=DET=',3F10.5) RETURN END
!C******************************************************************
- 10 -