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

Last change on this file since 3807 was 3757, checked in by emillour, 2 months ago

Mars PCM:
More code tidying: puting routines in modules and modernizing some old
constructs. Tested to not change results with respect to previous version.
EM

File size: 8.5 KB
Line 
1MODULE calldrag_noro_mod
2     
3IMPLICIT NONE
4     
5CONTAINS
6     
7      SUBROUTINE calldrag_noro(ngrid,nlayer,ptimestep,pplay,pplev,pt,pu,pv, &
8                  !output
9                  pdtgw,pdugw,pdvgw)
10                 
11!=====================================================================================================
12! MODULE designed to call SUBROUTINE drag_noro_mod Interface for sub-grid scale orographic scheme
13! The purpose of this subroutine is:
14! 1) Make some initial calculation at first call
15! 2) Split the calculation in several sub-grid ("sub-domain") to save memory and be able run on
16!   a workstation at high resolution.The sub-grid size is defined in dimradmars_mod.
17! 3) Transfer the increment output of drag_noro_mod into tendencies
18! Christophe Hourdin + Francois Forget
19!
20! VERSION Update: This version uses the variable's splitting, which can be usefull when performing
21! very high resolution simulation like LES. !
22! J.-B. Madeleine 10W12
23!
24! UPDATE   Rewirten into .F90     J.Liu   3/3/2022      Translate the code into .F90. The name of the
25!                                                       variables are uniformed. Comments are made.
26!                                                       Unused variables are deleted.
27!=======================================================================================================
28
29      use surfdat_h, only: zstd, zsig, zgam, zthe   !sub-grid scale standard devitation,slope,anisotropy, and priciple axes angle
30                                                    ! which will be coded as pvar, psig,pgam,pthe in the called subroutines.
31      use dimradmars_mod, only: ndomainsz
32      use drag_noro_mod, only: drag_noro
33      use sugwd_mod, only: sugwd
34      USE ioipsl_getin_p_mod, ONLY : getin_p
35     
36      IMPLICIT NONE
37
38      ! 0.1  Declarations :
39
40      ! 0.1.0 Input
41      INTEGER, intent(in):: ngrid                 ! Number of atmospheric columns
42      INTEGER, intent(in):: nlayer                ! Number of atmospheric layers
43      REAL, intent(in):: ptimestep                ! Time step of the Physics(s)
44      REAL, intent(in):: pplev(ngrid,nlayer+1)    ! Pressure at 1/2 levels(Pa)
45      REAL, intent(in):: pplay(ngrid,nlayer)      ! Pressure at full levels(Pa)
46      REAL, intent(in):: pt(ngrid,nlayer)         ! Temp at full levels(K)
47      REAL, intent(in):: pu(ngrid,nlayer)         ! Zonal wind at full levels(m/s)
48      REAL, intent(in):: pv(ngrid,nlayer)         ! Meridional wind at full levels(m/s)
49     
50      ! 0.1.1 Output   
51      REAL, intent(out):: pdtgw(ngrid,nlayer)     ! Tendency on temperature (K/s) due to orographic gravity waves
52      REAL, intent(out):: pdugw(ngrid,nlayer)     ! Tendency on zonal wind (m/s/s) due to orographic gravity waves
53      REAL, intent(out):: pdvgw(ngrid,nlayer)     ! Tendency on meridional wind (m/s/s) due to orographic gravity waves
54
55      !0.1.2 Local variables :
56      REAL sigtest(nlayer+1)      ! sigma coordinate at 1/2 level
57      INTEGER kgwd                ! Number of points at which the scheme is called
58      INTEGER kgwdim              ! kgwdIM=MAX(1,kgwd)
59      INTEGER ktest(ngrid)        ! Map of calling points
60      INTEGER kdx(ndomainsz)      ! Points at which to call the scheme
61      INTEGER ndomain             ! Parameter (ndomain = (ngrid-1) / ndomainsz + 1)
62      INTEGER l,ig
63      INTEGER jd,ig0,nd
64      REAL,save:: zstdthread      ! Detecting points concerned by the scheme by height
65!$OMP THREADPRIVATE(zstdthread)     
66      REAL zulow(ngrid)              ! Low level zonal wind
67      REAL zvlow(ngrid)              ! Low level Meridional wind
68      REAL zustr(ngrid)              ! Low level zonal stress
69      REAL zvstr(ngrid)              ! Low level meridional stress
70      REAL zplev(ndomainsz,nlayer+1) ! pplev in sub domain
71      REAL zplay(ndomainsz,nlayer)   ! pplay in sub domain
72      REAL zt(ndomainsz,nlayer)      ! pt in sub domain
73      REAL zu(ndomainsz,nlayer)      ! pu in sub domain
74      REAL zv(ndomainsz,nlayer)      ! pv in sub domain
75      REAL d_t(ndomainsz,nlayer)     ! Temperature increment(K) from DRAG_NORO
76      REAL d_u(ndomainsz,nlayer)     ! Zonal wind increment(K=m/s) from DRAG_NORO
77      REAL d_v(ndomainsz,nlayer)     ! Meridional wind increment(K=m/s) from DRAG_NORO
78      logical,save:: ll=.false.
79!$OMP THREADPRIVATE(ll)
80      LOGICAL,save:: firstcall=.true.
81!$OMP THREADPRIVATE(firstcall)
82
83!-----------------------------------------------------------------------------------------------------------------------
84!  1. INITIALISATIONS
85!-----------------------------------------------------------------------------------------------------------------------
86      IF (firstcall) THEN     
87        write(*,*) "calldrag_noro: orographic GW scheme is active!"
88        zstdthread = 50.0  ! Mars' value !!
89        call getin_p("oro_gwd_zstdthread", zstdthread)
90        write(*,*) "oro_gwd_zstdthread: zstdthread=", zstdthread
91       
92        do l=1,nlayer+1
93           sigtest(l)=pplev(1,l)/pplev(1,1)
94        enddo
95         
96        call sugwd(nlayer,sigtest)
97
98        if (ngrid .EQ. 1) then
99           IF (ndomainsz .NE. 1) THEN
100             print*
101             print*,'ATTENTION !!!'
102             print*,'pour tourner en 1D, meme pour drag_noro '
103             print*,'fixer ndomainsz=1 dans phymars/dimradmars_mod'
104             print*
105             call exit(1)
106           ENDIF
107        endif               
108        firstcall=.false.
109      END IF !firstcall
110     
111      ! Initialization the tendency of this scheme into zero
112!      pdtgw(ngrid,nlayer)=0.0     ! Tendency on temperature (K/s) due to orographic gravity waves
113!      pdugw(ngrid,nlayer)=0.0     ! Tendency on zonal wind (m/s/s) due to orographic gravity waves
114!      pdvgw(ngrid,nlayer)=0.0     ! Tendency on meridional wind (m/s/s) due to orographic gravity waves
115     
116      ! AS: moved out of firstcall to allow nesting+evoluting horiz domain
117      ndomain = (ngrid-1) / ndomainsz + 1
118!-----------------------------------------------------------------------------------------------------------------------     
119! 2. Starting loop on sub-domain
120!-----------------------------------------------------------------------------------------------------------------------
121      ! start the loop from the first sub-domain to the ndomain, in which the orographic gravity waves induced
122      ! wind/temperature tendencies are calculated
123      DO jd=1,ndomain
124         ig0=(jd-1)*ndomainsz
125         if (jd.eq.ndomain) then
126            nd=ngrid-ig0
127         else
128            nd=ndomainsz
129         endif
130                     
131         ! 2.1 Detecting points concerned by the scheme
132         kgwd=0
133         DO ig=ig0+1,ig0+nd
134            ktest(ig)=0
135            ! only the points have a zstd greater than the thread value(default is 50) are concerned.
136            ll=zstd(ig).gt.zstdthread
137            IF(ll) then
138               ! The map of calling is set 1 (which means the orographic map in this points is called?)
139               ktest(ig)=1
140               ! The numbers of the calling points added
141               kgwd=kgwd+1
142               ! The position of the calling points are picked up
143               kdx(kgwd)=ig - ig0
144            ENDIF
145         ENDDO
146         kgwdIM=MAX(1,kgwd)
147
148         ! 2.2 Spliting input variable in sub-domain input variables in each level
149         DO l=1,nlayer+1
150            do ig = 1,nd
151               zplev(ig,l) = pplev(ig0+ig,l)
152            enddo
153         end DO
154         DO l=1,nlayer
155            do ig = 1,nd
156               zplay(ig,l) = pplay(ig0+ig,l)
157               zt(ig,l) = pt(ig0+ig,l)
158               zu(ig,l) = pu(ig0+ig,l)
159               zv(ig,l) = pv(ig0+ig,l)
160            enddo
161         end DO
162
163         ! 2.3 Calling gravity wave and subgrid scale topo parameterization(all physics in this subroutine) to get
164         ! zonal wind, meridional wind, and temperature increament
165         call drag_noro (nd,nlayer,ptimestep,zplay,zplev,zstd(ig0+1),zsig(ig0+1),zgam(ig0+1),zthe(ig0+1),&
166              kgwd,kgwdim,kdx,ktest(ig0+1),zt, zu, zv, &
167              ! output
168              zulow(ig0+1),zvlow(ig0+1),zustr(ig0+1),zvstr(ig0+1),d_t,d_u,d_v)
169
170         ! 2.4  Un-spliting output variable from sub-domain input variables (and devide by ptimestep -> true tendencies)
171         DO l=1,nlayer
172            do ig = 1,nd
173                   pdtgw(ig0+ig,l) = d_t(ig,l)/ptimestep
174                   pdugw(ig0+ig,l) = d_u(ig,l)/ptimestep
175                   pdvgw(ig0+ig,l) = d_v(ig,l)/ptimestep
176            enddo
177         end DO
178
179      end DO ! (jd=1, ndomain)     
180!     2. ending loop on sub-domain
181
182!******************************************************************************
183
184      END SUBROUTINE calldrag_noro  ! end of the subroutine   
185     
186END MODULE calldrag_noro_mod  ! end of the module
187
Note: See TracBrowser for help on using the repository browser.