PGF90 (Version 12.8) 08/30/2020 00:27:47 page 1 Switches: -noasm -nodclchk -nodebug -nodlines -noline -list -idir ../include -inform warn -opt 1 -nosave -object -noonetrip -depchk on -nostandard -nosymbol -noupcase Filename: sig2p.f90 ( 1) subroutine sig2p ( 2) ( 3) USE mod_vrbls ( 4) USE mod_extra ( 5) implicit none ( 6) include 'param_o.h' ( 7) !! include 'extra.comm' ../include/param_o.h ( 1)* integer,parameter :: im0=401 ( 2)* integer, parameter :: nm=6 ( 3)* integer, parameter :: lm = 50 ( 4)* integer, parameter :: nsub = 10 ( 5)* ( 6)* integer, parameter :: im=im0, jm=im ( 7)* integer, parameter :: im1=im-1, jm1=jm-1 ( 8)* ( 9)* integer, parameter :: lm1 = lm-1, lp1 =lm +1 ( 10)* ( 11)* integer, parameter :: ixm = nsub, jym = ixm, nxy = ixm*jym ( 12)* integer, parameter :: ildom = (im - 1)/ixm, jldom = (jm - 1)/jym ( 13)* integer, parameter :: ilm = (im - 1)/ixm +1, jlm = (jm - 1)/jym +1 ( 14)* ( 15)* logical,parameter::flat=.false. ( 16)* logical,parameter::hstst=.false. ( 17)* ( 18)*! integer,parameter::igm=360, jgm = 181 ( 19)* integer,parameter::igm=(im0-1)*4, jgm = (igm/2)+1 ( 20)* real,parameter::alfa=0.0, beta=0.0, gamm=0.0 ( 21)*! real,parameter::alfa=0.*3.1415926/180.,beta=66.*3.1415926/180. & ( 22)*! ,gamm=175.*3.1415926/180. ( 23)* ( 24)*!GSM integer,parameter::lsm=20 ( 25)* integer,parameter::lsm=31 ( 26)* ( 27)*!RESTART ( 28)* ( 29)* character(len=10):: restartdate='2020080412' !restart date to create appropriate folder name ( 8) include 'mapot.comm' ../include/mapot.comm ( 1)* real,dimension(lsm)::spl,alsl ( 2)* ( 3)* common/mapot/ spl,alsl ( 9) include 'const.h' ../include/const.h ( 1)* logical :: run, first, restrt, subpost ( 2)* integer :: nfcst,nbc,list,ntsd,nddamp,nprec, & ( 3)* nboco,nshde,ncp,ntddmp ( 4)* logical,parameter::lcornerm=.FALSE. ( 5)* ( 6)* real,parameter :: dt=40 ! length of each time step ( 7)* logical,parameter:: sigma=.false. PGF90 (Version 12.8) 08/30/2020 00:27:47 page 2 ( 8)*!GSM integer,parameter::outnum=40 ! number of outputs ( 9)*!GSM integer,parameter::nday=10 ! number of days ( 10)* integer,parameter:: ntstm=3600*24*10/dt ! total time steps ( 11)*! integer,parameter:: idtad=1 ( 12)*! integer,parameter:: ncnvc=45 ( 13)*!test integer,parameter::nphs=45 ( 14)* ( 15)* integer,parameter:: idtad=2 ( 16)* integer,parameter:: ncnvc=6 ( 17)* integer,parameter::nphs=6 ( 18)* ( 19)* ( 20)* integer,parameter::nradsh=1 ( 21)*!GSM integer,parameter::nradlh=2 ( 22)* integer,parameter::nradlh=1 ( 23)* real,parameter::tsph=3600./dt ( 24)* integer,parameter::nrads=tsph*nradsh ( 25)* integer,parameter::nradl=tsph*nradlh ( 26)* real,parameter::weig=0.25 ( 27)* ( 28)*! ( 29)*! parameters for diagnostics ( 30)*! ( 31)* ( 32)* integer,parameter::idgns=10 ( 33)* integer,parameter::dgnstr=24*3600*200/dt ( 34)* integer,parameter::dgnnum=(ntstm-dgnstr)/idgns+1 ( 35)* integer,parameter::ickmm=1 ( 36)* ( 37)*! ( 38)*! parameter for initial values ( 39)*! ( 40)* real(KIND=4),parameter::pt = 2500.0 ( 41)*!!!test Dragan 1dec2015 ( 42)*!! real,parameter::pt=1000 ( 43)* ( 44)*! ( 45)*! date ( 46)*! ( 47)* integer::idat(3) ( 48)*!GSM data idat/07,21,2020/ ( 49)*!GSM integer::ihrst=0 ( 50)* ( 51)*! variable sst ( 52)* logical::lvsst ( 53)* ( 54)*! data assimilation constants ( 55)* ( 56)* integer,parameter::ndassim=3*3600/dt ( 57)*! integer,parameter::ndassimm=8*ndassim ( 58)* integer,parameter::ndassimm=0 ( 59)* integer,parameter::outstr=ndassimm ! step at which starts output ( 60)* integer,parameter::outend= ntstm+ndassimm ! step at which ends output ( 61)*!GSM integer,parameter::iout=(outend-outstr)/outnum ! steps between two output ( 62)* ( 63)*! data choice ( 64)* ( 65)* integer,parameter::sstc=1 ! 1 is NCEP sst data PGF90 (Version 12.8) 08/30/2020 00:27:47 page 3 ( 66)* ! 2 is TMI sst data ( 67)* ( 68)* ( 10) include 'dynam.comm' ( 11) ../include/dynam.comm ( 1)* real::fadv,fadv2, fadt, rd, f4d, ef4t, fkin ( 2)* real, dimension (lm) :: deta, rdeta, aeta, daeta,f4q2 ( 3)* real, dimension (lp1) :: eta, dfl ( 4)* ( 5)* common/dynam/ fadv,fadv2,fadt,rd,f4d,f4q2,ef4t,fkin,deta,rdeta, & ( 6)* eta,dfl,aeta,daeta ( 7)* ( 8)* ( 9)*!!!!!!!!!!!copy from dynam_comm.h ( 10)*! real :: fadv, fadv2,fadt, rd, f4d, ef4t, fkin, fcp ( 11)*! real, dimension (lm) :: deta, rdeta, aeta, daeta,f4q2 ( 12)*! real, dimension (lp1) :: eta, dfl ( 13)*! real, dimension (0:im+1,0:jm+1,nm) :: wpdar, f11, f12, f21, f22, & ( 14)*! p11, p12, p21, p22, fdiv, fddmp,fvdiff & ( 15)*! ,hbmsk,hsinp,hcosp ( 16)*! common /dynam/ fadv,fadv2, fadt, rd, f4d,f4q2, ef4t, fkin, & ( 17)*! deta, rdeta, aeta, eta, dfl, daeta, & ( 18)*! wpdar, f11, f12, f21, f22, & ( 19)*! p11, p12, p21, p22, fcp, fdiv, fddmp,hsinp,hcosp,& ( 20)*! fvdiff,hbmsk ( 12) integer,parameter::im_jm_nm=(im+2)*(jm+2)*nm ( 13) real,parameter::g=9.80616,RD1=287.04,GAMMA=6.5E-3, & ( 14) RGAMOG=RD1*GAMMA/G ( 15) integer::nhold,lp,i,j,n,lsl,nn,l,lmvij,lv ( 16) real,dimension(0:im+1,0:jm+1,nm)::tsl,qsl,fsl ( 17) integer,dimension(0:im+1,0:jm+1,nm)::nl1x ( 18) real,dimension(0:im+1,0:jm+1,nm)::alpetux,alpet2x,usl,vsl ( 19) integer,dimension(im_jm_nm)::ihold,jhold,nnhold ( 20) real::pnl1,B,fac,pu,tu,tabv,pl,tl,tabo,ahf,tblo, & ( 21) petal,petau,alpetl,alpetu,alpet2,fact,alpet1,trf ( 22) ( 23) open(1,file='../data_out/tmpfld.dat',form='unformatted') ( 24) ( 25) lsl=lsm ( 26) DO 310 LP=1,LSL ( 27) NHOLD=0 ( 28) ( 29) do 125 n=1,nm ( 30) DO 125 J=1,jm ( 31) DO 125 I=1,IM ( 32) ! ( 33) TSL(I,J,n)=-1.E6 ( 34) QSL(I,J,n)=-1.E6 ( 35) FSL(I,J,n)=-1.E6 ( 36) TRF=2.*ALSL(LP) ( 37) ! ( 38) !*** LOCATE VERTICAL INDEX OF MODEL INTERFACE JUST BELOW ( 39) !*** THE PRESSURE LEVEL TO WHICH WE ARE INTERPOLATING. ( 40) ! ( 41) DO 115 L=2,LP1 ( 42) if(i.eq.10.and.j.eq.12.and.n.eq.1)then ( 43) ! print *,l,alpint(i,j,n,l),alsl(lp) PGF90 (Version 12.8) 08/30/2020 00:27:47 page 4 ( 44) endif ( 45) IF(ALPINT(I,J,n,L).GE.ALSL(LP))THEN ( 46) NL1X(I,J,n)=L ( 47) NHOLD=NHOLD+1 ( 48) IHOLD(NHOLD)=I ( 49) JHOLD(NHOLD)=J ( 50) nnhold(nhold)=n ( 51) GO TO 125 ( 52) ELSEIF(ALPINT(I,J,n,LP1).LT.ALSL(LP))THEN ( 53) NL1X(I,J,n)=LP1 ( 54) NHOLD=NHOLD+1 ( 55) IHOLD(NHOLD)=I ( 56) JHOLD(NHOLD)=J ( 57) nnhold(nhold)=n ( 58) GO TO 125 ( 59) ENDIF ( 60) 115 CONTINUE ( 61) ! ( 62) 125 CONTINUE ( 63) ( 64) DO 220 NN=1,NHOLD ( 65) I=IHOLD(NN) ( 66) J=JHOLD(NN) ( 67) n=nnHOLD(NN) ( 68) PNL1=PINT(I,J,n,NL1X(I,J,n)) ( 69) ( 70) IF(NL1X(I,J,n).EQ.1)THEN ( 71) !--------------------------------------------------------------------- ( 72) !*** EXTRAPOLATE ABOVE THE TOPMOST MIDLAYER OF THE MODEL ( 73) !--------------------------------------------------------------------- ( 74) ! ( 75) PU=PINT(I,J,n,2) ( 76) TU=0.5*(T(I,J,n,1)+T(I,J,n,2)) ( 77) ( 78) TABV=TU*(SPL(LP)/PU)**RGAMOG ( 79) ( 80) B =TABV ( 81) FAC =0. ( 82) AHF =0. ( 83) ELSEIF(NL1X(I,J,n).EQ.LP1)THEN ( 84) !--------------------------------------------------------------------- ( 85) !*** EXTRAPOLATE BELOW LOWEST MODEL MIDLAYER (BUT STILL ABOVE GROUND) ( 86) !--------------------------------------------------------------------- ( 87) ! ( 88) PL=PINT(I,J,n,LM-1) ( 89) TL=0.5*(T(I,J,n,LM-2)+T(I,J,n,LM-1)) ( 90) TBLO=TL*(SPL(LP)/PL)**RGAMOG ( 91) ( 92) B =TBLO ( 93) FAC =0. ( 94) AHF =0. ( 95) ELSE ( 96) !--------------------------------------------------------------------- ( 97) !*** INTERPOLATION BETWEEN NORMAL LOWER AND UPPER BOUNDS ( 98) !--------------------------------------------------------------------- ( 99) ! ( 100) B =T(I,J,n,NL1X(I,J,n)) ( 101) FAC =2.*ALOG(PT+PDSL(I,J,n)*AETA(NL1X(I,J,n))) PGF90 (Version 12.8) 08/30/2020 00:27:47 page 5 ( 102) AHF =(B-T(I,J,n,NL1X(I,J,n)-1))/ & ( 103) (ALPINT(I,J,n,NL1X(I,J,n)+1)-ALPINT(I,J,n,NL1X(I,J,n)-1)) ( 104) ENDIF ( 105) ( 106) TSL(I,J,n)=B+AHF*(TRF-FAC) ( 107) FSL(I,J,n)=(PNL1-SPL(LP))/(SPL(LP)+PNL1) & ( 108) *((ALSL(LP)+ALPINT(I,J,n,NL1X(I,J,n))-FAC)*AHF+B)*Rd*2. & ( 109) +ZINT(I,J,n,NL1X(I,J,n))*G ( 110) if(i.eq.30.and.j.eq.101.and.n.eq.3)then ( 111) print *,i,j,n,tsl(i,j,n),t(i,j,n,37),t(i,j,n,38),spl(lp),pl ( 112) endif ( 113) 220 CONTINUE ( 114) ( 115) ( 116) do 281 n=1,nm ( 117) DO 281 J=1,jm-1 ( 118) DO 281 I=1,IM-1 ( 119) !NOTE ( 120) !NOTE 29 JANUARY 1993, RUSS TREADON. ( 121) !NOTE - AS FOR THE OTHER FIELDS WE INTERPOLATE ONLY ( 122) !NOTE BETWEEN THE FAL AND THE MODEL TOP. BELOW ( 123) !NOTE SURFACE VALUES ARE FAL VALUES. ( 124) ! ( 125) LMVIJ=LM ( 126) ! ( 127) DO 280 LV=2,LMVIJ+1 ( 128) PETAL=PT+PDVP1(I,J,n)*ETA(LV) ( 129) PETAU=PT+PDVP1(I,J,n)*ETA(LV-1) ( 130) ALPETL=ALOG(PETAL) ( 131) ALPETU=ALOG(PETAU) ( 132) ALPET2=SQRT(0.5*(ALPETL*ALPETL+ALPETU*ALPETU)) ( 133) ! ( 134) !*** SEARCH FOR HIGHEST MID-LAYER MODEL SURFACE ( 135) !*** THAT IS BELOW THE GIVEN STANDARD PRESSURE LEVEL. ( 136) ! ( 137) IF(ALSL(LP).LT.ALPET2)THEN ( 138) NL1X(I,J,n)=LV-1 ( 139) ALPETUX(I,J,n)=ALPETU ( 140) ALPET2X(I,J,n)=ALPET2 ( 141) GO TO 281 ( 142) ENDIF ( 143) 280 CONTINUE ( 144) ! ( 145) NL1X(I,J,n)=LMVIJ+1 ( 146) ALPETUX(I,J,n)=ALPETU ( 147) ALPET2X(I,J,n)=ALPET2 ( 148) 281 CONTINUE ( 149) ! ( 150) !--------------------------------------------------------------------- ( 151) !*** BELOW GROUND USE WIND IN LOWEST LAYER ( 152) !--------------------------------------------------------------------- ( 153) ! ( 154) do 290 n=1,nm ( 155) DO 290 J=1,jm-1 ( 156) DO 290 I=1,IM-1 ( 157) LMVIJ=LM ( 158) IF(NL1X(I,J,n).GT.LMVIJ)THEN ( 159) USL(I,J,n)=U(I,J,n,LMVIJ) PGF90 (Version 12.8) 08/30/2020 00:27:47 page 6 ( 160) VSL(I,J,n)=V(I,J,n,LMVIJ) ( 161) ! ( 162) !--------------------------------------------------------------------- ( 163) !*** IF REQUESTED PRESSURE LEVEL IS NOT BELOW THE LOCAL GROUND ( 164) !*** THEN WE HAVE TWO POSSIBILITIES. IF THE REQUESTED PRESSURE ( 165) !*** LEVEL IS BETWEEN THE LOCAL SURFACE PRESSURE AND TOP OF ( 166) !*** MODEL PRESSURE, VERTICALLY INTERPOLATE BETWEEN NEAREST ( 167) !*** BOUNDING ETA LEVELS TO GET THE WIND COMPONENTS. IF THE ( 168) !*** REQUESTED PRESSURE LEVEL IS ABOVE THE MODEL TOP, USE ( 169) !*** CONSTANT EXTRAPOLATION OF TOP ETA LAYER (L=1) WINDS. ( 170) !--------------------------------------------------------------------- ( 171) ! ( 172) ELSE ( 173) IF(NL1X(I,J,n).GT.1)THEN ( 174) ALPETL=ALPETUX(I,J,n) ( 175) PETAU=PT+PDVP1(I,J,n)*ETA(NL1X(I,J,n)-1) ( 176) ALPETU=ALOG(PETAU) ( 177) ALPET1=SQRT(0.5*(ALPETL*ALPETL+ALPETU*ALPETU)) ( 178) FACT=(ALPET2X(I,J,n)-ALSL(LP))/(ALPET2X(I,J,n)-ALPET1) ( 179) USL(I,J,n)=U(I,J,n,NL1X(I,J,n)) & ( 180) +(U(I,J,n,NL1X(I,J,n)-1)-U(I,J,n,NL1X(I,J,n)))*FACT ( 181) VSL(I,J,n)=V(I,J,n,NL1X(I,J,n))+(V(I,J,n,NL1X(I,J,n)-1) & ( 182) -V(I,J,n,NL1X(I,J,n)))*FACT ( 183) ELSE ( 184) USL(I,J,n)=U(I,J,n,NL1X(I,J,n)) ( 185) VSL(I,J,n)=V(I,J,n,NL1X(I,J,n)) ( 186) ENDIF ( 187) ! ( 188) !--------------------------------------------------------------------- ( 189) !*** ALPET2 IS MID-LAYER ETA SURFACE JUST BELOW STANDARD PRESSURE ( 190) !*** LEVEL AND ALPET1 IS MID-LAYER ETA SURFACE JUST ABOVE. ( 191) !*** NOTE THAT IF THE STANDARD PRESSURE SURFACE IS SUBMERGED, THEN ( 192) !*** ALPET2 AND ALPET1 ARE THE LOWEST AND 2ND LOWEST MID-LAYER ( 193) !*** ETA SURFACES ABOVE THE TOPOGRAPHY (WITH OLDRD=.TRUE., ZJ). ( 194) !--------------------------------------------------------------------- ( 195) ! ( 196) ENDIF ( 197) 290 CONTINUE ( 198) ( 199) call bocoh(fsl,im,jm,nm) ( 200) call bocov(usl,vsl,im,jm,nm) ( 201) write(1)fsl ( 202) write(1)usl ( 203) write(1)vsl ( 204) write(1)tsl ( 205) ( 206) 310 continue ( 207) ( 208) do n=1,nm ( 209) do j=1,jm ( 210) do i=1,im ( 211) ! print *,i,j,n,fsl(i,j,n,1) ( 212) enddo ( 213) enddo ( 214) enddo ( 215) ( 216) close(1) ( 217) PGF90 (Version 12.8) 08/30/2020 00:27:47 page 7 ( 218) end subroutine sig2p ( 219)