C Twave.for C ----------------------------------------------------------------------!------- ! 27 Dec 2010 - Subroutine is working ! 15 Dec 2010 - First try ! 3 Dec 2010 - Copied from Ttanh6.for and Tplume.for C ----------------------------------------------------------------------!------- SUBROUTINE Twave C ----------------------------------------------------------------------!------- C TEMPERATURE PERTURBATION MODEL C Twave is a temperature perturbation model that has height profiles C composed of three segments with hyperbolic tangent transitions C that are the same as those in ttanh5 and ttanh6. The parameters for C the profiles vary harmonically with latitude to emulate a gravity wave. C The lower segment is meant to approximately cancel the gradient of the C background profile so that the total temperature in the lower segment C is constant in both height and latitude. The upper segment is meant to C give approximately zero temperature perturation to the background C temperature. The middle segment gives a sharp transition from the C lower segment to the upper segment. The height of the transition C varies sinusoidally with latitude. C ----------------------------------------------------------------------!------- implicit double precision (a-h,o-z) TYPE :: values_and_derivs_1_2 DOUBLE PRECISION :: value DOUBLE PRECISION :: ptm,pr,pth,pph DOUBLE PRECISION :: ptmtm,prtm,pthtm,pphtm DOUBLE PRECISION :: ptmr,prr,pthr,pphr DOUBLE PRECISION :: ptmth,prth,pthth,pphth DOUBLE PRECISION :: ptmph,prph,pthph,pphph END TYPE TYPE (values_and_derivs_1_2) :: temperature TYPE :: model CHARACTER :: model_name DOUBLE PRECISION :: data_id CHARACTER :: pert_name DOUBLE PRECISION :: pert_id END TYPE model TYPE (model) :: modtg DOUBLE PRECISION :: h, z1, delta1, z2, delta2, c1, c2 TYPE :: value_and_derivs_1_dim DOUBLE PRECISION :: value DOUBLE PRECISION :: p,pp END TYPE TYPE (value_and_derivs_1_dim) :: psiA1, psiB1sg, psiC1 TYPE (value_and_derivs_1_dim) :: psiA2, psiB2sg, psiC2 TYPE :: value_and_derivs_2_dim DOUBLE PRECISION :: value DOUBLE PRECISION :: h,g DOUBLE PRECISION :: hh,hg,gg C g means d/d cos gamma END TYPE TYPE (value_and_derivs_2_dim) :: psi1sg, psi2sg, psisg COMMON/PCONST/CREF,RGAS,GAMMA 1 ,MO2,MN2,MO,MN,MH2O,MO3 DOUBLE PRECISION :: MO2,MN2,MO,MN,MH2O,MO3 COMMON/MCONST/PI,PIT2,PID2,DEGS,RAD,ALN10 C COMMON DECK "RKAM" INSERTED HERE DOUBLE PRECISION :: KR,KTH,KPH COMMON//R,TH,PH,KR,KTH,KPH,RKVARS(14),TPULSE,CSTEP,DRDT(20) C C COMMON DECK "TT" INSERTED HERE DOUBLE PRECISION :: MODT COMMON/TT/MODT(4),T,PTT,PTR,PTTH,PTPH C COMMON DECK "WW" INSERTED HERE COMMON/WW/ dum(10),MAXW,W(1000) EQUIVALENCE (EARTHR,W(1)) C DELTA TEMPERATURE 225-249 EQUIVALENCE (W(225),DTMODEL),(W(226),DTFORM),(W(227),DTID) EQUIVALENCE (h0,W(228)),(deltaH,W(229)) & ,(deltaT,W(230)),(delta1,W(231)) & ,(delta2,W(232)),(slope,W(233)) & ,(A,W(234)),(kgrav0,W(235)),(phase,W(236)) DOUBLE PRECISION :: kgrav0,kgravTH C COMMON DECK "B6" INSERTED HERE INTEGER DTMX,DTNTBL,DTITBL,DTFRMTB,IDSDT(10) COMMON/B6/DTMX,DTNTBL(10),DTITBL(10),DTFRMTB(10),DTGP(611) EQUIVALENCE (DTGP,IDSDT),(ANDT,DTGP(11)) DOUBLE PRECISION :: 1 c(100),dTm(99),z(99),dl(99),DLC(100),ALC(100),DC(100),TC(100) 2 ,pdelTz(0:101),pcz(0:101) INTEGER dTPX ,dTQTBL(10),dTLTBL(10),dTIRMTB(10) DATA dTPX/2/ DATA dTQTBL/1,11,612,7*0/ DATA dTLTBL/1,100,8*0/ DATA dTIRMTB/1,2,8*0/ DATA RECOGdT,N/10.d0,0/ ! ? DATA AQdT/0.d0/ C ----------------------------------------------------------------------!------- ENTRY SETPTMP C ----------------------------------------------------------------------!------- ANdT=AQdT dTMX=dTPX CALL IMOVE(dTNTBL,dTQTBL,10) CALL IMOVE(dTITBL,dTLTBL,10) CALL IMOVE(dTFRMTB,dTIRMTB,10) RETURN C ----------------------------------------------------------------------!------- ENTRY IPTEMP C ----------------------------------------------------------------------!------- IF(W(225).NE.6.) THEN C IF(IDINT(W(225)).NE.6) THEN PRINT*,' Model check number mismatch: ' STOP ' W(225) should be 6 if you are using model Twave.' ENDIF C MODT(3)='Twave' MODT(4)=W(227) z0=0.d0 z(3)=1000.0d0 dTm(1)=-deltaT dTm(2)=0.d0 dTm(3)=0.d0 dT0=dTm(1)+(z0-h0)*slope dl(1)=delta1 dl(2)=delta2 dl(3)=0.d0 kgravTH=earthr*kgrav0 n=2 RETURN C C ----------------------------------------------------------------------!------- ENTRY IPTEMPg(modtg) C ----------------------------------------------------------------------!------- C This subroutine is called once for each new W array. C IF(W(225).NE.6.) THEN C IF(IDINT(W(225)).NE.6) THEN PRINT*,' Model check number mismatch: ' STOP ' W(225) should be 6 if you are using model Twave.' ENDIF C modtg%pert_name = 'Twave' modtg%pert_id = W(227) RETURN C ----------------------------------------------------------------------!------- ENTRY PTEMP C ----------------------------------------------------------------------!------- H=R-EARTHR z(1)=h0-deltaH/2.d0+A*dcos(kgravTH*(pi/2.d0-TH)+phase) z(2)=h0+deltaH/2.d0+A*dcos(kgravTH*(pi/2.d0-TH)+phase) CALL ITANH3(N,Z0,Z,dT0,dTM,DL,C,CAV,DLC,ALC,DC & ,TC,PCZ) CALL ftanh3(h,delT,delpTr,pdelTz,n,z0,z,dT0,dl,c,CAV,DLC,ALC,DC & ,TC,PCZ) pzipTH = A*kgravTH*dsin(kgravTH*(pi/2.d0-TH)+phase) c PRINT *,'delT = ',delT delpTt = 0.d0 delpTth = pzipTH*(pdelTz(1)+pdelTz(2)) delpTph = 0.d0 T=T+delT PTR =PTR + delPTR PTTH=PTTH + delpTTH ! PTPH=PTPH + delPTPH RETURN C ----------------------------------------------------------------------!------- ENTRY PTEMPg(temperature) C ----------------------------------------------------------------------!------- h = r - earthr c NEED TO ADD SOME CALLS TO ITANH3 HERE ! (Actually, for gravity waves, we do need the second derivative of T with ! respect to height and latitude, so we would need to generalize ftanh2 ! and ftanh3. Neither would be enough.) CALL ftanh3 & (h,delT,delpTr,pdelTz,n,z0,z,delT0,dl,c,CAV,DLC,ALC,DC & ,PCZ) c PRINT *,' at tempg, delT0 = ', delT0 temperature%value = temperature%value + delT temperature%pr = temperature%pr + delpTr temperature%prr = temperature%prr + delpTrr c PRINT *,'delT = ',delT c PRINT *,'delpTr = ',delpTr c PRINT *,'delpTrr = ',temperature%prr C Still need to add above formulas plus higher derivatives for gravity version RETURN C ----------------------------------------------------------------------!------- END SUBROUTINE Twave C ----------------------------------------------------------------------!-------