CTITLESAZM16 -- SETS UP COEFF. IN TRI-DIAGONAL MATRIX FOR ZMPS 00000102 C*********************************************************************** C COPYRIGHT ATLANTIC RICHFIELD COMPANY 1991 * C*********************************************************************** CA 00000200 CA AUTHOR F. G. SHERRILL 00000300 CA DESIGNER F. G. SHERRILL/R. D. KNIGHT 00000400 CA LANGUAGE FORTRAN H 00000500 CA SYSTEM S/370 00000600 C WRITTEN JAN 1983 00000700 C REVISED 00000800 CA 00000900 CA 00001000 CA CALL SAZM16 ( BETA,CVIS,DXQ,IPR,ISRTIO,ISTO,IVCHK,IZ,NX, 00001100 CA * NZ,NZE,NZF,NZS,OM,OMDZ,Q1,Q3,RR1,RR3,RVIS,S,VEL ) 00001200 CA 00001300 CA 00001400 CA IN/OUT ARG TYPE DESCRIPTION 00001500 CA 00001600 CA IN BETA R4 CONSTANT USED IN 2ND DERIVATIVE APPROXIMATION00001700 CA IN CVIS R4 COMPLEX PART OF VISCOSITY CONSTANT 00001800 CA IN DXQ R4 CONSTANT EQUAL TO DZ/4 00001900 CA IN IPR I4 UNIT NUMBER OF SYSOUT FILE 00002000 CA IN ISRTIO I4 SHOT INTERVAL TO RECEIVER INTERVAL RATIO 00002100 CA IN ISTO I4 CONSTANT LATERAL VELOCITY CHECKING ARRAY 00002200 CA IN IVCHK I4 FLAG ARRAY FOR VELOCITIES WHICH ARE CONSTANT 00002300 CA LATERALLY AND ARE EQUAL IN ADJACENT LAYERS 00002400 CA IN IZ I4 INDEX OF THE CURRENT LAYER 00002500 CA IN NX I4 TOTAL TRACES TO MIGRATE 00002600 CA IN NZ I4 TOTAL DEPTH LAYERS 00002700 CA IN NZE I4 INDEX OF THE LAST LAYER 00002800 CA IN NZF I4 TOTAL DOWNWARD CONTINUATION STEPS 00002900 CA IN NZS I4 INDEX OF THE STARTING LAYER 00003000 CA IN OM R4 RADIAL FREQUENCY 00003100 CA IN OMDZ R4 OM*DZ 00003200 CA IN Q1 R4 CONSTANT IN THE 45 DEGREE APPROXIMATION 00003300 CA (NOT SAME AS Q1 IN SAZM13) 00003400 CA IN Q3 R4 CONSTANT IN THE 45 DEGREE APPROXIMATION 00003500 CA IN RR1 R4 FREQUENCY/DEPTH INDEPENDENT VARIABLE 00003600 CA IN RR3 R4 FREQUENCY/DEPTH INDEPENDENT VARIABLE 00003700 CA IN RVIS R4 REAL PART OF NUMERICAL VISCOSITY CONSTANT 00003800 CA IN/OUT S R4 ARRAY IMPLICITY CONTAINING WORK ARRAYS 00003900 CA IN VEL R4 ARRAY CONTAINING VELOCITIES FOR CURRENT 00004000 CA AND PREVIOUS LAYERS 00004100 CA 00004200 CA THIS SUBROUTINE COMPUTES MAIN DIAGONAL(B AND BB) AND OFF 00004300 CA DIAGONAL (A AND AA) ARRAYS IN TRI-DIAGONAL MATRIX. THESE 00004400 CA ARE NEEDED TO COMPUTE THE WAVEFIELD AT THE NEXT Z-LEVEL. 00004500 CA THE PHASE SHIFT ARRAY IS ALSO CALCULATED IN THIS SUBROUTINE 00004600 CA AS ARE SOME OF THE RECURSIVE ARRAYS IF THERE IS A CONSTANT 00004700 CA VELOCITY STRIP. 00004800 CA 00004900 CA SUBROUTINES CALLED: NONE 00005000 CA 00005100 C 00005200 C 00005300 SUBROUTINE SAZM16 (BETA,CVIS,DXQ,IPR,ISRTIO,ISTO,IVCHK,IZ,NX, 00005400 * NZ,NZE,NZF,NZS,OM,OMDZ,Q1,Q3,RR1,RR3,RVIS,S,VEL ) 00005500 C 00005600 C REAL ARRAYS--EXTERNAL 00005700 C 00005800 REAL S ( 01) 00005900 REAL VEL(NX, 2) 00006000 C 00006100 C INTEGER ARRAYS--EXTERNAL 00006200 C 00006300 DIMENSION ISTO( NX) 00006400 DIMENSION IVCHK( NZ) 00006500 C 00006600 C INTEGER VARIABLES--LOCAL 00006700 C 00006800 INTEGER AAI 00006900 INTEGER AAR 00007000 INTEGER AI 00007100 INTEGER AR 00007200 INTEGER BBI 00007300 INTEGER BBR 00007400 INTEGER BI 00007500 INTEGER BR 00007600 INTEGER DENDEN 00007700 INTEGER DENE 00007800 INTEGER DENI 00007900 INTEGER DENR 00008000 INTEGER DI 00008100 INTEGER DR 00008200 INTEGER EDENR 00008300 INTEGER EDENI 00008400 INTEGER EEI 00008500 INTEGER EER 00008600 INTEGER EI 00008700 INTEGER ER 00008800 INTEGER FI 00008900 INTEGER FR 00009000 INTEGER SHIFTI 00009100 INTEGER SHIFTR 00009200 INTEGER TI 00009300 INTEGER TR 00009400 C 00009500 C GENERATE RELATIVE ADDRESSES FOR WORK AREAS 00009600 C 00009700 AAI = 0 00009800 AAR = AAI + NX 00009900 AI = AAR + NX 00010000 AR = AI + NX 00010100 BBI = AR + NX 00010200 BBR = BBI + NX 00010300 BI = BBR + NX 00010400 BR = BI + NX 00010500 DENDEN = BR + NX 00010600 DENE = DENDEN+ NX 00010700 DENI = DENE + NX 00010800 DENR = DENI + NX 00010900 DI = DENR + NX 00011000 DR = DI + NX 00011100 EDENR = DR + NX 00011200 EDENI = EDENR + NX 00011300 EEI = EDENI + NX 00011400 EER = EEI + NX 00011500 EI = EER + NX 00011600 ER = EI + NX 00011700 FI = ER + NX 00011800 FR = FI + NX 00011900 IWAVR = FR + NX 00012000 IWAVI = IWAVR + NX 00012100 SHIFTI = IWAVI + NX 00012200 SHIFTR = SHIFTI+ NX 00012300 TI = SHIFTR+ NX 00012400 TR = TI + NX 00012500 C 00012600 C ==================================================================== 00012700 C FIND THE TRI-D ELEMENTS FOR THIS LAYER 00012800 C ==================================================================== 00012900 C 00013000 IBGN = 1 00013100 C 00013200 DO 20 KK = 1 ,NX 00013300 KEND =ISTO(KK) 00013400 VCT = VEL(KEND, 2) 00013500 RRSET =RR1*VCT 00013600 CC3 = RR3*OM/VCT 00013700 COMP =(RVIS * RRSET) + (CC3 *BETA) 00013800 RREAL = -(RRSET * CVIS) 00013900 ZAAR = Q1-RREAL 00014000 ZAAI = - COMP 00014100 ZAR = ZAAR-Q3 00014200 ZBBR = 2.*RREAL-1. 00014300 ZBBI = 2. * COMP -CC3 00014400 ZBR = ZBBR + 2. 00014500 ZSFTR= COS(OMDZ/VCT) 00014600 ZSFTI= SIN(OMDZ/VCT) 00014700 C 00014800 DO 10 IK = IBGN,KEND 00014900 S(AAR+IK) = ZAAR 00015000 S(AAI+IK) = ZAAI 00015100 S(AR+IK) = ZAR 00015200 S(AI+IK) = ZAAI 00015300 S(BBR+IK) = ZBBR 00015400 S(BBI+IK) = ZBBI 00015500 S(BR+IK) = ZBR 00015600 S(BI+IK) = ZBBI 00015700 S(SHIFTR+IK)= ZSFTR 00015800 10 S(SHIFTI+IK)= ZSFTI 00015900 C 00016000 IF(KEND.EQ.NX) GO TO 30 00016100 20 IBGN =KEND+1 00016200 C 00016300 30 IF(IVCHK(IZ).NE.1) GO TO 70 00016400 IF((IZ.EQ.NZS).OR.(ISRTIO.NE.1).OR.(IZ.EQ.NZE)) GO TO 40 00016500 IF((IVCHK(IZ).EQ.1).AND.(IVCHK(IZ-1).EQ.1).AND. 00016600 * (VEL(1,2).EQ.VEL(1,1)) ) GO TO 70 00016700 C 00016800 40 BAB1= OM/VEL(1, 2) * DXQ 00016900 RA1D = 1. + BAB1*BAB1 00017000 RA1R = (1 - BAB1*BAB1)/RA1D 00017100 RA1I = (BAB1+BAB1)/RA1D 00017200 BCR = S(AR+1)*RA1R - S(AI+1)*RA1I 00017300 BCI = S(AR+1)*RA1I + S(AI+1)*RA1R 00017400 BBCR = S(AAR+1)*RA1R - S(AAI+1)*RA1I 00017500 BBCI = S(AAR+1)*RA1I + S(AAI+1)*RA1R 00017600 BR1 = S(BR+1) + BCR 00017700 BI1 = S(BI+1) + BCI 00017800 BBR1= S(BBR+1) + BBCR 00017900 BBI1= S(BBI+1) + BBCI 00018000 NN = NX-1 00018100 C 00018200 50 DB = BR1 * BR1 + BI1 * BI1 00018300 S(EER+1) =-(BR1 * S(AR+2) + BI1* S(AI+2))/DB 00018400 S(EEI+1) = (BI1 * S(AR+2) - BR1* S(AI+2))/DB 00018500 C 00018600 DO 60 I =2,NN 00018700 S(EDENR+I)=S(BR+I)+S(AR+I-1)*S(EER+I-1)-S(AI+I-1)*S(EEI+ 00018800 * I-1) 00018900 S(EDENI+I)=S(BI+I)+S(AI+I-1)*S(EER+I-1)+S(AR+I-1)*S(EEI+ 00019000 * I-1) 00019100 S(DENE+I)=S(EDENR+I)*S(EDENR+I)+S(EDENI+I)*S(EDENI+I) 00019200 S(EER+I)=-(S(AR+I+1)*S(EDENR+I)+S(AI+I+1)*S(EDENI+I))/S( 00019300 * DENE+I) 00019400 60 S(EEI+I)= (S(AR+I+1)*S(EDENI+I)-S(AI+I+1)*S(EDENR+I))/S( 00019500 * DENE+I) 00019600 C 00019700 70 RETURN 00019800 C 00019900 END 00020000