《工程地质学程序》由会员分享，可在线阅读，更多相关《工程地质学程序（14页珍藏版）》请在金锄头文库上搜索。
1、READ(5,*) NS !土条总数 READ(5,*)(X2(I),I=1,NS) READ(5,*)(Y2(I),I=1,NS) N=NS+1 DO 302 I=1,N READ(5,*)W(I),RU(I),C(I),F(I),BDA(I),WL(I),FQUH(I),FDIS(I),YSF(I) 302 CONTINUE !读入两端点和土条底中点的数据,其中涉及力的物理量均以单宽土条计 !X2,Y2=坐标,W=土条重量,RU=孔压系数,C=凝聚力=,F=磨擦系数,BDA=土容重, !WL=浸润线 Y 坐标,FQUH,FDIS,YSF=水平地震力的大小,Y 坐标,垂直水平地震力的大小 I
2、F(OPTION.EQ.3)THEN READ(5,*)ASP WRITE(6,*)陆军工程师团法，平均坝坡=,ASP ASP=ASP*3.14159/180. ENDIF WRITE(6,709)RW 709 FORMAT(1X,水容重=,F10.3/) WRITE(6,702) 702 FORMAT(T25,DATA FOR ASSUMED SIDE READ(5,*) NS !土条总数 READ(5,*)(X2(I),I=1,NS) READ(5,*)(Y2(I),I=1,NS) N=NS+1 DO 302 I=1,N READ(5,*)W(I),RU(I),C(I),F(I),BDA(
3、I),WL(I),FQUH(I),FDIS(I),YSF(I) 302 CONTINUE !读入两端点和土条底中点的数据,其中涉及力的物理量均以单宽土条计 !X2,Y2=坐标,W=土条重量,RU=孔压系数,C=凝聚力=,F=磨擦系数,BDA=土容重, !WL=浸润线 Y 坐标,FQUH,FDIS,YSF=水平地震力的大小,Y 坐标,垂直水平地震力的大小 IF(OPTION.EQ.3)THEN READ(5,*)ASP WRITE(6,*)陆军工程师团法，平均坝坡=,ASP ASP=ASP*3.14159/180. ENDIF WRITE(6,709)RW 709 FORMAT(1X,水容重=,
4、F10.3/) WRITE(6,702)DO 208 I=2,NS X(I)=(X2(I)+X2(I-1)/2. Y(I)=(Y2(I)+Y2(I-1)/2. ALF(I)=ATAN(Y2(I)-Y2(I-1)/(X2(I)-X2(I-1) 208 CONTINUE ALF(1)=ALF(2) ALF(N)=0 X(N)=X2(NS)Y(N)=Y2(NS) X(1)=X2(1) Y(1)=Y2(1) WRITE(6,716) 716 FORMAT(T5,NO. ,T14,X,T24,Y,T34,W,T44,RU,T54,C, #T64,F,T74,BDA,T84,WL,T94,QH,T104,
5、YE) DO 901 I=1,N WRITE(6,718) I,X(I),Y(I),W(I),RU(I),C(I),F(I),BDA(I),WL(I) #,FQUH(I),FDIS(I) 901 CONTINUE 718 FORMAT(1X,I5,11F10.3) WRITE(6,719) 719 FORMAT(1X,/NOTE:) WRITE(6,720) 720 FORMAT(5X,W=土条的单宽重量) WRITE(6,721) 721 FORMAT(5X,RU=孔隙水压力系数) WRITE(6,722) 722 FORMAT(5X,C=凝聚力) WRITE(6,723) 723 FORM
6、AT(5X,F=摩擦系数) WRITE(6,724) 724 FORMAT(5X,BDA=土条的平均容重) WRITE(6,725) FORMAT(5X,WL=浸润线的 Y 坐标) SELECT CASE (OPTION) CASE(0) WRITE(6,701) 701 FORMAT(/10X,SPENCER METHOD) CALL SAFE2(X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,4) IWALL=0 DO 301 I=1,80 FXO(I)=0 FX(I)=1 301 CONTINUE !tan(BETA(I)=FXO(I)+ALAM*FX(I) KXYX
7、=0 CALL MP(N,KXYX,KS,AB1,DQ,DM,BN,FX,FXO, &C,F,X2,Y2,X,Y,W,RU,FQUH,WL,FDIS,RW,DB,DF) CASE(1) WRITE(6,704) 704 FORMAT(/20X,BISHOP 法/) CALL SAFE2(X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,OPTION) CASE(2) WRITE(6,705)705 FORMAT(/20X,瑞典法/) CALL SAFE2(X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,OPTION) CASE(3) WRITE(6,731
8、) 731 FORMAT(/10X,工程师团法) CALL SAFE2(X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,4) IWALL=0 DO 311 I=1,80 FXO(I)=0 FX(I)=1 311 CONTINUE !tan(BETA(I)=FXO(I)+ALAM*FX(I) KXYX=1 CALL MP(N,KXYX,KS,AB1,DQ,DM,BN,FX,FXO, &C,F,X2,Y2,X,Y,W,RU,FQUH,WL,FDIS,RW,DB,DF) CASE(4) WRITE(6,789) 789 FORMAT(/20X,简化法法/) CALL SAFE2(
9、X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,OPTION) CASE DEFAULT END SELECT END SUBROUTINE SAFE2(X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,ID) REAL X2(80),F(80),FQUH(80),W(80),RU(80),ALF(80),X(80),C(80) REAL AB1,FDIS(80) CINTEGER N REAL RO1,DM,RM,BX1,CDN,DG(80),G3,DX1,G,GX REAL O,S1,FRIC,S2,S,DRM,P1,P2,DDM,DBX1,DCDN
10、CALL SAFE2(X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,OPTION) CASE(2) WRITE(6,705) 705 FORMAT(/20X,瑞典法/) CALL SAFE2(X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,OPTION) CASE(3) WRITE(6,731) 731 FORMAT(/10X,工程师团法) CALL SAFE2(X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,4) IWALL=0 DO 311 I=1,80 FXO(I)=0 FX(I)=1 311 CONTINUE !tan(BE
11、TA(I)=FXO(I)+ALAM*FX(I) KXYX=1 CALL MP(N,KXYX,KS,AB1,DQ,DM,BN,FX,FXO, &C,F,X2,Y2,X,Y,W,RU,FQUH,WL,FDIS,RW,DB,DF)CASE(4) WRITE(6,789) 789 FORMAT(/20X,简化法法/) CALL SAFE2(X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,OPTION) CASE DEFAULT END SELECT END SUBROUTINE SAFE2(X2,F,FQUH,FDIS,W,RU,N,X,ALF,C,AB1,ID) REAL X2(
12、80),F(80),FQUH(80),W(80),RU(80),ALF(80),X(80),C(80) REAL AB1,FDIS(80) CINTEGER N REAL RO1,DM,RM,BX1,CDN,DG(80),G3,DX1,G,GX REAL O,S1,FRIC,S2,S,DRM,P1,P2,DDM,DBX1,DCDN C 79 CONTINUE DG(I+1)=DG(I)+(F(I)-F(I+1)*ALF(I) 93 G3=G3+(F(I)*ALF(I)+DG(I)*DX1 G3=G3/(X(N)-X(1) DO 202 I=2,NS IF(ID_C.EQ.0)THEN RO1=
13、1 ELSE RO1=(FDIS(I)-CY)/(DS-CY) ENDIF G=ALF(I) GX=G*F(I)-G3+DG(I) DX1=X2(I)-X2(I-1) O=FQUH(I)*SIN(G) S1=W(I)*(COS(G)-RU(I)/COS(G) FRIC=F(I) S2=(S1-O)*FRIC S=C(I)/COS(G) DRM=S+S2 DRM=DRM*DX1 C 79 CONTINUE DG(I+1)=DG(I)+(F(I)-F(I+1)*ALF(I) 93 G3=G3+(F(I)*ALF(I)+DG(I)*DX1 G3=G3/(X(N)-X(1) DO 202 I=2,NS
14、 IF(ID_C.EQ.0)THEN RO1=1 ELSE RO1=(FDIS(I)-CY)/(DS-CY) ENDIF G=ALF(I)GX=G*F(I)-G3+DG(I) DX1=X2(I)-X2(I-1) O=FQUH(I)*SIN(G) S1=W(I)*(COS(G)-RU(I)/COS(G) FRIC=F(I) S2=(S1-O)*FRIC S=C(I)/COS(G) DRM=S+S2 DRM=DRM*DX1 P1=W(I)*SIN(G) P2=FQUH(I)*RO1 DDM=(P1+P2)*DX1 DM=DM+DDM RM=RM+DRM DBX1=-DRM-SIN(G)*GX*W(
15、I)*DX1 BX1=BX1+DBX1 DCDN=DRM*GX CDN=CDN+DCDN 202 CONTINUE AB1=-BX1/DM+CDN/BX1 !通过简化法求得 AB1 的迭代初值 IF(ID.EQ.4)WRITE(6,*)简化法安全系数= ,AB1 IF(ID.EQ.2)AB1=RM/DM IF(ID.EQ.2)WRITE(6,*)瑞典法安全系数= ,AB1 IF(ID.NE.1)RETURN WRITE(6,*)迭代过程 LM=1 41 DT=0. DM=0. RM=0. DO 21 I=2,NS G=ALF(I) DX1=X2(I)-X2(I-1) S1=W(I)*(1-RU
16、(I)*F(I) S2=S1+C(I) S=COS(G)*(1+F(I)*TAN(G)/AB1) DRM=DX1*S2/S RM=RM+DRM P1=W(I)*SIN(G) P2=FQUH(I)*RO1 DDM=(P1+P2)*DX1 DM=DM+DDM DT=DT+DRM*TAN(G)*F(I)/(AB1*2+AB1*F(I)*TAN(G) 21 CONTINUE AB2=RM/DM !BISHOP 方法IF(AB2.LT.0.OR.AB2.GE.15) THEN WRITE(6,*)F.S NOT REASONABLE = ,AB2 AB2=9999.9 RETURN ENDIF WRITE(6,740) LM,AB2 740 FORMAT(1X,I5,F10.5) IF(ABS(AB2-AB1).LT.0.0001)GO TO 49 IF(LM.LT.10)GO TO 9191 WRITE(6,919
