source: trunk/LMDZ.MARS/libf/phymars/calldrag_noro_mod.F90 @ 2997

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