MODULE calldrag_noro_mod IMPLICIT NONE CONTAINS SUBROUTINE calldrag_noro(ngrid,nlayer,ptimestep,pplay,pplev,pt,pu,pv, & !output pdtgw,pdugw,pdvgw) !===================================================================================================== ! MODULE designed to call SUBROUTINE drag_noro_mod Interface for sub-grid scale orographic scheme ! The purpose of this subroutine is: ! 1) Make some initial calculation at first call ! 2) Split the calculation in several sub-grid ("sub-domain") to save memory and be able run on ! a workstation at high resolution.The sub-grid size is defined in dimradmars_mod. ! 3) Transfer the increment output of drag_noro_mod into tendencies ! Christophe Hourdin + Francois Forget ! ! VERSION Update: This version uses the variable's splitting, which can be usefull when performing ! very high resolution simulation like LES. ! ! J.-B. Madeleine 10W12 ! ! UPDATE Rewirten into .F90 J.Liu 3/3/2022 Translate the code into .F90. The name of the ! variables are uniformed. Comments are made. ! Unused variables are deleted. !======================================================================================================= use surfdat_h, only: zstd, zsig, zgam, zthe !sub-grid scale standard devitation,slope,anisotropy, and priciple axes angle ! which will be coded as pvar, psig,pgam,pthe in the called subroutines. use dimradmars_mod, only: ndomainsz use drag_noro_mod, only: drag_noro USE ioipsl_getin_p_mod, ONLY : getin_p IMPLICIT NONE ! 0.1 Declarations : ! 0.1.0 Input INTEGER, intent(in):: ngrid ! Number of atmospheric columns INTEGER, intent(in):: nlayer ! Number of atmospheric layers REAL, intent(in):: ptimestep ! Time step of the Physics(s) REAL, intent(in):: pplev(ngrid,nlayer+1) ! Pressure at 1/2 levels(Pa) REAL, intent(in):: pplay(ngrid,nlayer) ! Pressure at full levels(Pa) REAL, intent(in):: pt(ngrid,nlayer) ! Temp at full levels(K) REAL, intent(in):: pu(ngrid,nlayer) ! Zonal wind at full levels(m/s) REAL, intent(in):: pv(ngrid,nlayer) ! Meridional wind at full levels(m/s) ! 0.1.1 Output REAL, intent(out):: pdtgw(ngrid,nlayer) ! Tendency on temperature (K/s) due to orographic gravity waves REAL, intent(out):: pdugw(ngrid,nlayer) ! Tendency on zonal wind (m/s/s) due to orographic gravity waves REAL, intent(out):: pdvgw(ngrid,nlayer) ! Tendency on meridional wind (m/s/s) due to orographic gravity waves !0.1.2 Local variables : REAL sigtest(nlayer+1) ! sigma coordinate at 1/2 level INTEGER kgwd ! Number of points at which the scheme is called INTEGER kgwdim ! kgwdIM=MAX(1,kgwd) INTEGER ktest(ngrid) ! Map of calling points INTEGER kdx(ndomainsz) ! Points at which to call the scheme INTEGER ndomain ! Parameter (ndomain = (ngrid-1) / ndomainsz + 1) INTEGER l,ig INTEGER jd,ig0,nd REAL,save:: zstdthread ! Detecting points concerned by the scheme by height !$OMP THREADPRIVATE(zstdthread) REAL zulow(ngrid) ! Low level zonal wind REAL zvlow(ngrid) ! Low level Meridional wind REAL zustr(ngrid) ! Low level zonal stress REAL zvstr(ngrid) ! Low level meridional stress REAL zplev(ndomainsz,nlayer+1) ! pplev in sub domain REAL zplay(ndomainsz,nlayer) ! pplay in sub domain REAL zt(ndomainsz,nlayer) ! pt in sub domain REAL zu(ndomainsz,nlayer) ! pu in sub domain REAL zv(ndomainsz,nlayer) ! pv in sub domain REAL d_t(ndomainsz,nlayer) ! Temperature increment(K) from DRAG_NORO REAL d_u(ndomainsz,nlayer) ! Zonal wind increment(K=m/s) from DRAG_NORO REAL d_v(ndomainsz,nlayer) ! Meridional wind increment(K=m/s) from DRAG_NORO logical,save:: ll=.false. !$OMP THREADPRIVATE(ll) LOGICAL,save:: firstcall=.true. !$OMP THREADPRIVATE(firstcall) !----------------------------------------------------------------------------------------------------------------------- ! 1. INITIALISATIONS !----------------------------------------------------------------------------------------------------------------------- IF (firstcall) THEN write(*,*) "calldrag_noro: orographic GW scheme is active!" zstdthread = 50.0 ! Mars' value !! call getin_p("oro_gwd_zstdthread", zstdthread) write(*,*) "oro_gwd_zstdthread: zstdthread=", zstdthread do l=1,nlayer+1 sigtest(l)=pplev(1,l)/pplev(1,1) enddo call sugwd(nlayer,sigtest) if (ngrid .EQ. 1) then IF (ndomainsz .NE. 1) THEN print* print*,'ATTENTION !!!' print*,'pour tourner en 1D, meme pour drag_noro ' print*,'fixer ndomainsz=1 dans phymars/dimradmars_mod' print* call exit(1) ENDIF endif firstcall=.false. END IF !firstcall ! Initialization the tendency of this scheme into zero ! pdtgw(ngrid,nlayer)=0.0 ! Tendency on temperature (K/s) due to orographic gravity waves ! pdugw(ngrid,nlayer)=0.0 ! Tendency on zonal wind (m/s/s) due to orographic gravity waves ! pdvgw(ngrid,nlayer)=0.0 ! Tendency on meridional wind (m/s/s) due to orographic gravity waves ! AS: moved out of firstcall to allow nesting+evoluting horiz domain ndomain = (ngrid-1) / ndomainsz + 1 !----------------------------------------------------------------------------------------------------------------------- ! 2. Starting loop on sub-domain !----------------------------------------------------------------------------------------------------------------------- ! start the loop from the first sub-domain to the ndomain, in which the orographic gravity waves induced ! wind/temperature tendencies are calculated DO jd=1,ndomain ig0=(jd-1)*ndomainsz if (jd.eq.ndomain) then nd=ngrid-ig0 else nd=ndomainsz endif ! 2.1 Detecting points concerned by the scheme kgwd=0 DO ig=ig0+1,ig0+nd ktest(ig)=0 ! only the points have a zstd greater than the thread value(default is 50) are concerned. ll=zstd(ig).gt.zstdthread IF(ll) then ! The map of calling is set 1 (which means the orographic map in this points is called?) ktest(ig)=1 ! The numbers of the calling points added kgwd=kgwd+1 ! The position of the calling points are picked up kdx(kgwd)=ig - ig0 ENDIF ENDDO kgwdIM=MAX(1,kgwd) ! 2.2 Spliting input variable in sub-domain input variables in each level DO l=1,nlayer+1 do ig = 1,nd zplev(ig,l) = pplev(ig0+ig,l) enddo end DO DO l=1,nlayer do ig = 1,nd zplay(ig,l) = pplay(ig0+ig,l) zt(ig,l) = pt(ig0+ig,l) zu(ig,l) = pu(ig0+ig,l) zv(ig,l) = pv(ig0+ig,l) enddo end DO ! 2.3 Calling gravity wave and subgrid scale topo parameterization(all physics in this subroutine) to get ! zonal wind, meridional wind, and temperature increament call drag_noro (nd,nlayer,ptimestep,zplay,zplev,zstd(ig0+1),zsig(ig0+1),zgam(ig0+1),zthe(ig0+1),& kgwd,kgwdim,kdx,ktest(ig0+1),zt, zu, zv, & ! output zulow(ig0+1),zvlow(ig0+1),zustr(ig0+1),zvstr(ig0+1),d_t,d_u,d_v) ! 2.4 Un-spliting output variable from sub-domain input variables (and devide by ptimestep -> true tendencies) DO l=1,nlayer do ig = 1,nd pdtgw(ig0+ig,l) = d_t(ig,l)/ptimestep pdugw(ig0+ig,l) = d_u(ig,l)/ptimestep pdvgw(ig0+ig,l) = d_v(ig,l)/ptimestep enddo end DO end DO ! (jd=1, ndomain) ! 2. ending loop on sub-domain !****************************************************************************** END SUBROUTINE calldrag_noro ! end of the subroutine END MODULE calldrag_noro_mod ! end of the module