MODULE drag_noro_mod IMPLICIT NONE CONTAINS SUBROUTINE drag_noro (ngrid,nlayer,ptimestep,pplay,pplev,pvar, psig, pgam, pthe, & kgwd,kgwdim,kdx,ktest,t, u, v, & ! output pulow, pvlow, pustr, pvstr, d_t, d_u, d_v) !-------------------------------------------------------------------------------- ! MODULE contains DRAG_NORO SUBROUNTINE for sub-grid scale orographic sheme. ! 1) Initialization the variables ! 2) Overturn the levels and computes geopotential height ! 3) Call the and ORODRAG subroutine to updates the tendencies ! 4) Overturn back to the normal level of the tendencies and transfer it into ! increments and sum the oro-gw stress. ! the scheme has been called. ! Z.X. Li + F.Lott (LMD/CNRS/ENS) 1995-02-01 ! ! UPDATE: J.Liu 03/03/2022 Translate the code into .F90. The name of the ! variables are uniformed. Comments are made. ! Unused variables are deleted. !------------------------------------------------------------------------------- use dimradmars_mod, only: ndomainsz USE orodrag_mod, ONLY: orodrag USE comcstfi_h, ONLY: g, r IMPLICIT none ! 0. DECLARATIONS: ! 0.1 Input: INTEGER, intent(in):: ngrid ! Number of atmospheric columns [only nd] INTEGER, intent(in):: nlayer ! Number of atmospheric layers REAL, intent(in):: ptimestep ! Time step of the Physics(s) real, intent(in):: pplay(ndomainsz,nlayer) ! Pressure at full levels(Pa) real, intent(in):: pplev(ndomainsz,nlayer+1) ! Pressure at 1/2 levels(Pa) REAL, intent(in):: pvar(ndomainsz) ! sub-grid scale standard deviation REAL, intent(in):: psig(ndomainsz) ! sub-grid scale standard slope REAL, intent(in):: pgam(ndomainsz) ! sub-grid scale standard anisotropy REAL, intent(in):: pthe(ndomainsz) ! sub-grid scale principal axes angle REAL, intent(in):: u(ndomainsz,nlayer) ! Zonal wind at full levels(m/s) REAL, intent(in):: v(ndomainsz,nlayer) ! Meridional wind at full levels(m/s) REAL, intent(in):: t(ndomainsz,nlayer) ! Temperature at full levels(m/s) INTEGER, intent(in):: kgwd ! Number of points at which the scheme is called INTEGER, intent(in):: kgwdim ! kgwdim=MAX(1,kgwd) INTEGER, intent(in):: kdx(ndomainsz) ! Points at which to call the scheme INTEGER, intent(in):: ktest(ndomainsz) ! Map of calling points ! 0.2 Output: REAL, intent(out):: pulow(ndomainsz) ! Low level zonal wind REAL, intent(out):: pvlow(ndomainsz) ! Low level meridional wind REAL, intent(out):: pustr(ndomainsz) ! Low level zonal stress REAL, intent(out):: pvstr(ndomainsz) ! Low level meridional stress REAL, intent(out):: d_t(ndomainsz,nlayer) ! Temperature increment(K) due to orographic gravity waves REAL, intent(out):: d_u(ndomainsz,nlayer) ! Zonal increment(m/s) due to orographic gravity waves REAL, intent(out):: d_v(ndomainsz,nlayer) ! Meridional increment(m/s) due to orographic gravity waves ! 0.3 Variables locales: INTEGER i, k REAL inv_pplev(ndomainsz,nlayer+1) ! Inversed (by inverse the pplev pressure) pressure at 1/2 levels REAL inv_pplay(ndomainsz,nlayer) ! Inversed (by inverse the pplay pressure) pressure at full levels REAL zgeom(ndomainsz,nlayer) ! Geopotetial height (Inversed ??) REAL inv_pt(ndomainsz,nlayer) ! Inversed (by inverse the t) temperature (K) at full levels REAL inv_pu(ndomainsz,nlayer) ! Inversed (by inverse the u) zonal wind (m/s) at full levels REAL inv_pv(ndomainsz,nlayer) ! Inversed (by inverse the v) meridional wind (m/s) at full levels REAL pdtdt(ndomainsz,nlayer) ! Temperature tendency outputs from main routine ORODRAG REAL pdudt(ndomainsz,nlayer) ! Zonal winds tendency outputs from main routine ORODRAG REAL pdvdt(ndomainsz,nlayer) ! Meridional winds tendency outputs from main routine ORODRAG !----------------------------------------------------------------------------------------------------------------------- ! 1. INITIALISATIONS !----------------------------------------------------------------------------------------------------------------------- ! 1.1 Initialize the input/output variables (for security purpose) DO i = 1,ngrid pulow(i) = 0.0 pvlow(i) = 0.0 pustr(i) = 0.0 pvstr(i) = 0.0 ENDDO DO k = 1, nlayer DO i = 1, ngrid d_t(i,k) = 0.0 d_u(i,k) = 0.0 d_v(i,k) = 0.0 pdudt(i,k)=0.0 pdvdt(i,k)=0.0 pdtdt(i,k)=0.0 ENDDO ENDDO ! 1.2 Prepare the ininv_put varibales (Attention: The order of vertical levels increases from top to bottom ) ! Here the levels are overturned, that is, the surface becomes to the top and the top becomes to the surface ! and so forth DO k = 1, nlayer DO i = 1, ngrid inv_pt(i,k) = t(i,nlayer-k+1) inv_pu(i,k) = u(i,nlayer-k+1) inv_pv(i,k) = v(i,nlayer-k+1) inv_pplay(i,k) = pplay(i,nlayer-k+1) inv_pplev(i,k) = pplev(i,nlayer+1-k+1) ENDDO endDO DO i = 1, ngrid inv_pplev(i,nlayer+1) = pplev(i,1) ENDDO !calculate g*dz by R*T*[LN(P(z)/P(z+dz))] DO i = 1, ngrid zgeom(i,nlayer) = r * inv_pt(i,nlayer)* LOG(inv_pplev(i,nlayer+1)/inv_pplay(i,nlayer)) ENDDO ! sum g*dz from surface to top to get geopotential height DO k = nlayer-1, 1, -1 DO i = 1, ngrid zgeom(i,k) = zgeom(i,k+1) + r * (inv_pt(i,k)+inv_pt(i,k+1))/2.0 * LOG(inv_pplay(i,k+1)/inv_pplay(i,k)) ENDDO endDO !----------------------------------------------------------------------------------------------------------------------- ! 2. CALL the Main Rountine OORDRAG !----------------------------------------------------------------------------------------------------------------------- ! 2.1 call the main rountine to get the tendencies CALL ORODRAG(ngrid,nlayer,kgwd,kgwdim,kdx,ktest,ptimestep,inv_pplev,inv_pplay,zgeom,& inv_pt, inv_pu, inv_pv, pvar, psig, pgam, pthe, & ! output: pulow,pvlow,pdudt,pdvdt,pdtdt) ! 2.2 Transfer tendency into increment by multiply the physical time steps ! maybe in the future here cancel multiply the ptimestep thus it will not divide ptimestep again in the main rountine ! calldrag_noro_mod.F90 DO k = 1, nlayer DO i = 1, ngrid d_u(i,nlayer+1-k) = ptimestep*pdudt(i,k) d_v(i,nlayer+1-k) = ptimestep*pdvdt(i,k) d_t(i,nlayer+1-k) = ptimestep*pdtdt(i,k) pustr(i) = pustr(i)+g*pdudt(i,k)*(inv_pplev(i,k+1)-inv_pplev(i,k)) pvstr(i) = pvstr(i)+g*pdvdt(i,k)*(inv_pplev(i,k+1)-inv_pplev(i,k)) ENDDO ENDDO END SUBROUTINE drag_noro END MODULE drag_noro_mod