rzb.txt

人工合成地震波的一个源程序

 !////////////////////////////////////////////////////////////////////////////////
 !////////人工合成地震波程序//////////////////////////////////////////////////////
 !////////////////////////////////////////////////////////////////////////////////
 ! G为重力加速度，S_G(400)为目标反应谱，W(400)为圆频率						//////
 ! DW为圆频率步长,S(400)为加速度功率谱，A(400)为幅值谱，AA(50)为平稳随机过程//////
 !T为时间划分数组，A_G为人工合成的加速度									//////
 !F1为地震影响系数曲线函数													//////
 !F2为强度包络线函数 														//////
 !JS_C为人工合成地震波的加速度反应谱    									//////
 !SIMPSON为复化辛普森积分子过程
 !设计反应谱采用教材235页中的影响系数曲线，其中地震烈度选用6度，设计阶段为第一阶段
 !包络线采用教材45页的分段函数型，其中的持时10秒
 !主程序： 
 IMPLICIT NONE
 INTEGER  N,M,I,J,K	!N为频率域等分数，M为时间域等分数
 REAL  G,S_G(400),W(400),DW,TT,PI,S(400),A1,A2,A(400),AA(50),DT,T(50),A_G(50),F1,F2,Q,P,QQ
 REAL  TS,AS,BS,JS_C(400),SMAX
 G=9.8
 N=400
 M=50
 W(1)=1.0
 W(400) =50.0
 DW=(W(400)-W(1))/N
 PI=3.1415926
 DT=10/50.0
 AS=0
 OPEN(1,FILE='EARTH1.TXT')
 OPEN(2,FILE='EARTH2.TXT')
 DO I=1,M
    T(I)=I*DT
 ENDDO

 DO I=2,N-1	  
   W(I)=W(1)+(I-0.5)*DW
 ENDDO

 DO I=1,N							 !求解目标反应谱
	TT=2*PI/W(I)
	S_G(I)=G*F1(TT)
 endDO
 WRITE(1,'(4X,10A)')'圆频率为：    目标反应谱为：'
 WRITE(1,10)(W(I),S_G(I),I=1,N)
 10 FORMAT(1X,F10.4,4X,F15.6)
 CLOSE(1)
 DO I=1,N							 !求解加速度功率谱
    A1=0.05/(PI*W(I))
	A2=LOG( (-PI*LOG(1-0.85)/(W(I)*10) )**(-1) )
	S(I)=A1*S_G(I)*S_G(I)/A2
 ENDDO
									 !求解幅值谱
DO I=1,N
   A(I)=SQRT(4*S(I)*DW)
ENDDO
DO I=1,M
   DO J=1,N
      K=I+J
	  CALL RANDOM_SEED(K)
      CALL RANDOM_NUMBER(Q)
	  Q=Q*2*PI
	  AA(I)=AA(I)+A(J)*COS(W(J)*T(I)-Q)	  !求解平稳随机过程
   ENDDO
ENDDO
DO I=1,M
   A_G(I)=F2(T(I))*AA(I)
ENDDO
WRITE(2,'(4X,15A)') '时刻点：    人工地震波：'		  !人工合成的加速度
WRITE(2,20)	(T(I),A_G(I),I=1,M)
20 FORMAT(1X,F10.4,F15.6)
CLOSE(2)

DO I=1,N											  !人工地震波的加速度反应谱
   P=W(I)
   SMAX=0.0
   DO J=1,M
	  QQ=A_G(J)
	  BS=T(M)
	  CALL 	SIMPSON(TS,AS,BS,P,QQ)
	  IF(ABS(TS)>SMAX) SMAX=ABS(TS)
   ENDDO
   JS_C(I)=SMAX*(P**2)
ENDDO
OPEN(3,FILE='ESRTH3.TXT')
WRITE(3,'(4X,10A)')'人工地震波的加速度反应谱为：'
WRITE(3,30)(JS_C(I),I=1,N)							  
30 FORMAT(1X,F15.8)
CLOSE(3) 
END

!//////////////////////////////////////
!//////////////////////////////////////
FUNCTION F1(T)                        !地震影响系数曲线
   IMPLICIT NONE
   REAL	 F1,T
   IF(T>=0 .AND. T<0.1)	 THEN
      F1=0.32*T+0.018
   ELSEIF(T>=0.1 .AND. T<0.2)  THEN
      F1=0.04
   ELSEIF(T>=0.2 .AND. T<2.0)	THEN
      F1=((0.2/T)**0.9)*0.04
   ELSE
      F1=0.004
  ENDIF
END FUNCTION F1 

 
FUNCTION F2(T)						 !强度包络线函数
  IMPLICIT NONE
  REAL	F2,T
  IF(T>=0 .AND. T<1.0 )	  THEN
     F2=T*T
  ELSEIF(T>=1 .AND. T<7.0)	THEN
     F2=1.0
  ELSEIF(T>=7.0 .AND. T<=10.0)	  THEN
  F2=EXP(-1.15*(T-7))
  ENDIF
END FUNCTION F2

SUBROUTINE SIMPSON(TS,AS,BS,P,QQ)						!为复化辛普森积分子过程
  IMPLICIT NONE
  REAL TS,AS,BS,P,QQ,EPS,T_N,H_N,T1,T2,H
  INTEGER M,N,I
  M=100
  EPS=0.00001
  N=M
  H=(BS-AS)/N
  T_N=(F(AS)+F(BS))/2
  DO I=1,N-1
     T_N=T_N+F(AS+H*I)
  ENDDO   
  T_N=T_N*H
  T2=T_N
  T1=T2+100
  DO WHILE(ABS(T1-T2)>EPS)
     T1=T2
     H_N=0
    DO I=0,N-1
      H_N=H_N+F(AS+H*I+H/2)
    ENDDO
    H_N=H_N*H
    T2=(T1+H_N)/2
    H=H/2
    N=N*2
  ENDDO
  TS=T2
CONTAINS
  FUNCTION  F(X)
  REAL X,F
  F=QQ*EXP(-0.05*P*(BS-X))*SIN(P*(BS-X)) /P
  END FUNCTION
END SUBROUTINE