弹性链杆法半衬砌程序.doc

DIMENSION X(31),Y(31),DYCS(30,4),KB(93,94),DRTA(93) REAL KB INTEGER JDCS(31,4),TXFLAG(31),CSFLAG OPEN(1,FILE=INPUT.DAT,STATUS=OLD) OPEN(2,FILE=OUTPUT.DAT,STATUS=UNKNOWN) CALL INPUT(NJ,X,Y,DYCS,JDCS,RZ,E,TKH,TKV,QH,QV,HA,TXFLAG) NTRY=0 11 CSFLAG=0 NTRY=NTRY+1 CALL SUBKB(NJ,DYCS,JDCS,TXFLAG,X,Y,QH,QV,TKH,TKV,E,RZ,HA,KB) N=3*NJ M=N+1 CALL FCQJ(KB,N,M,DRTA) WRITE(2,100) NTRY WRITE(*,100) NTRY 100 FORMAT(1X,I3,STEP) WRITE(2,200) 200 FORMAT(1X,NO,5X,X-DISP,5X,Y-DISP,1X,X-SPRING, * 1X,R-DISP) DO 21 I=1,NJ DISP=SQRT(DRTA(3*I-2)**2+DRTA(3*I-1)**2) WRITE(2,300) I,DRTA(3*I-2),DRTA(3*I-1),TXFLAG(I),DISP 21 CONTINUE 300 FORMAT(1X,I2,2E11.3,I5,E11.3) DO 31 I=1,NJ FLAG=DRTA(3*I-2) IF(FLAG.LE.0.0.AND.TXFLAG(I).EQ.0) THEN TXFLAG(I)=1 CSFLAG=1 END IF IF(FLAG.GT.0.0.AND.TXFLAG(I).EQ.1) THEN TXFLAG(I)=0 CSFLAG=1 END IF 31 CONTINUE IF(CSFLAG.NE.0) GO TO 11 CALL SUBFE(NJ,DYCS,E,DRTA) END SUBROUTINE INPUT(NJ,X,Y,DYCS,JDCS,RZ,E,TKH,TKV, QH,QV,HA,TXFLAG) DIMENSION X(31),Y(31),DYCS(30,4),JDHD(31) INTEGER JDCS(31,4),TXFLAG(31) REAL JDHD READ(1,*) NJ,RZ,E,QH,QV,TKH,TKV,HA DO 10 I=1,NJ READ(1,*) X(I),Y(I),(JDCS(I,J),J=1,4),JDHD(I) TXFLAG(I)=JDCS(I,1) 10 CONTINUE DO 20 I=1,NJ-1 DX=X(I+1)-X(I) DY=Y(I+1)-Y(I) DYCS(I,1)=SQRT(DX*DX+DY*DY) DYCS(I,2)=(JDHD(I+1)+JDHD(I))/2.0 DYCS(I,3)=DY/DYCS(I,1) DYCS(I,4)=DX/DYCS(I,1) 20 CONTINUE END SUBROUTINE SUBKT(NJ,Y,TKH,KT) REAL Y(31),KT(31) DO 10 I=1,NJ KT(I)=0 10 CONTINUE DO 20 I=2,NJ-1 KT(I)=TKH*ABS(Y(I+1)-Y(I-1))/2 20 CONTINUE KT(1)=TKH*ABS(Y(2)-Y(1))/2 KT(NJ)=TKH*ABS(Y(NJ)-Y(NJ-1))/2 END SUBROUTINE SUBKC(HA,TKV,KC) REAL KC(3,3) DO 10 I=1,3 DO 10 J=1,3 KC(I,J)=0.0 10 CONTINUE KC(1,1)=0 KC(2,2)=TKV*HA KC(3,3)=(TKV*HA**3)/12.0 END SUBROUTINE JDHZ(NJ,X,Y,RZ,QH,QV,DYCS,HZ) DIMENSION HZ(31,3),X(31),Y(31),DYCS(30,4) DO 10 I=1,NJ DO 10 J=1,3 HZ(I,J)=0 10 CONTINUE DO 20 I=1,NJ-1 PH=(Y(I+1)-Y(I))*QH/2 HZ(I,1)=HZ(I,1)+PH HZ(I+1,1)=HZ(I+1,1)+PH 20 CONTINUE DO