source: LMDZ6/trunk/libf/phylmd/lmdz_ratqs_main.f90 @ 5826

Last change on this file since 5826 was 5826, checked in by rkazeroni, 2 months ago

For GPU porting of ratqs_main routine:

  • Simplify refactoring of rev 5820
  • Add forgotten THREADPRIVATE
  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
File size: 10.8 KB
Line 
1!$gpum horizontal klon ngrid
2MODULE lmdz_ratqs_main
3  PRIVATE
4
5  LOGICAL, SAVE :: first=.TRUE.  ! fisrt call to ratqs_main
6  !$OMP THREADPRIVATE(first)
7
8  REAL, SAVE :: resolmax
9  !$OMP THREADPRIVATE(resolmax)
10
11  PUBLIC ratqs_main_first, ratqs_main
12
13CONTAINS
14
15SUBROUTINE ratqs_main_first(klon, cell_area)
16
17  IMPLICIT NONE
18  INTEGER, INTENT(in) :: klon
19  REAL, DIMENSION(klon), INTENT(in) :: cell_area
20
21  IF (first) THEN
22     resolmax=sqrt(maxval(cell_area))
23     first = .FALSE.
24  END IF
25
26END SUBROUTINE ratqs_main_first
27
28SUBROUTINE ratqs_main(klon,klev,nbsrf,prt_level,lunout,       &
29           iflag_ratqs,iflag_con,iflag_cld_th,pdtphys, &
30           ratqsbas,ratqshaut,ratqsp0,ratqsdp, &
31           pctsrf,s_pblh,zstd, &
32           tau_ratqs,fact_cldcon,wake_s, wake_deltaq,   &
33           ptconv,ptconvth,clwcon0th, rnebcon0th,       &
34           paprs,pplay,t_seri,q_seri,                   &
35           qtc_cv, sigt_cv,detrain_cv,fm_cv,fqd,fqcomp,sigd,zqsat,             &
36           omega,tke,tke_dissip,lmix,wprime, &
37           t2m,q2m,fm_therm,entr_therm,detr_therm,cell_area,&
38           ratqs,ratqsc,ratqs_inter_,sigma_qtherm)
39
40
41USE clouds_gno_mod,     ONLY: clouds_gno
42USE lmdz_ratqs_multi,   ONLY: ratqs_inter, ratqs_oro, ratqs_hetero, ratqs_tke
43
44implicit none
45
46!========================================================================
47! Computation of ratqs, the width of the subrid scale water distribution
48! (normalized by the mean value)
49! Various options controled by flags iflag_con and iflag_ratqs
50! F Hourdin 2012/12/06
51!========================================================================
52
53! Declarations
54
55! Input
56integer,intent(in) :: klon,klev,nbsrf,prt_level,lunout
57integer,intent(in) :: iflag_con,iflag_cld_th,iflag_ratqs
58real,intent(in) :: pdtphys,ratqsbas,ratqshaut,fact_cldcon,tau_ratqs
59real,intent(in) :: ratqsp0, ratqsdp
60real, dimension(klon,klev),intent(in) :: omega
61real, dimension(klon,klev+1),intent(in) :: paprs,tke,tke_dissip,lmix,wprime
62real, dimension(klon,klev),intent(in) :: pplay,t_seri,q_seri,zqsat
63real, dimension(klon,klev),intent(in) :: entr_therm,detr_therm,qtc_cv, sigt_cv
64real, dimension(klon,klev) :: detrain_cv,fm_cv,fqd,fqcomp
65real, dimension(klon) :: sigd
66
67real, dimension(klon,klev+1),intent(in) :: fm_therm
68logical, dimension(klon,klev),intent(in) :: ptconv
69real, dimension(klon,klev),intent(in) :: rnebcon0th,clwcon0th
70real, dimension(klon,klev),intent(in) :: wake_deltaq
71real, dimension(klon),intent(in) :: wake_s
72real, dimension(klon,nbsrf),intent(in) :: t2m,q2m
73real, dimension(klon), intent(in) :: cell_area
74real, dimension(klon,nbsrf),intent(in) :: pctsrf
75real, dimension(klon),intent(in) :: s_pblh
76real, dimension(klon),intent(in) :: zstd
77
78! Output
79real, dimension(klon,klev),intent(inout) :: ratqs,ratqsc,ratqs_inter_,sigma_qtherm
80
81logical, dimension(klon,klev),intent(inout) :: ptconvth
82
83! local
84integer i,k
85real, dimension(klon,klev) :: ratqss
86real facteur,zfratqs1,zfratqs2
87real, dimension(klon,klev) :: ratqs_hetero_,ratqs_oro_,ratqs_tke_
88real :: resol, fact
89
90!-------------------------------------------------------------------------
91!  Caclul des ratqs
92!-------------------------------------------------------------------------
93
94!      print*,'calcul des ratqs'
95!   ratqs convectifs a l'ancienne en fonction de q(z=0)-q / q
96!   ----------------
97!   on ecrase le tableau ratqsc calcule par clouds_gno
98      if (iflag_cld_th.eq.1) then
99         do k=1,klev
100         do i=1,klon
101            if(ptconv(i,k)) then
102              ratqsc(i,k)=ratqsbas &
103              +fact_cldcon*(q_seri(i,1)-q_seri(i,k))/q_seri(i,k)
104            else
105               ratqsc(i,k)=0.
106            endif
107         enddo
108         enddo
109
110!-----------------------------------------------------------------------
111!  par nversion de la fonction log normale
112!-----------------------------------------------------------------------
113      else if (iflag_cld_th.eq.4) then
114         ptconvth(:,:)=.false.
115         ratqsc(:,:)=0.
116         if(prt_level.ge.9) print*,'avant clouds_gno thermique'
117         call clouds_gno &
118         (klon,klev,q_seri,zqsat,clwcon0th,ptconvth,ratqsc,rnebcon0th)
119         if(prt_level.ge.9) print*,' CLOUDS_GNO OK'
120       
121       endif
122
123!   ratqs stables
124!   -------------
125
126      if (iflag_ratqs.eq.0) then
127
128! Le cas iflag_ratqs=0 correspond a la version IPCC 2005 du modele.
129         do k=1,klev
130            do i=1, klon
131               ratqss(i,k)=ratqsbas+(ratqshaut-ratqsbas)* &
132               min((paprs(i,1)-pplay(i,k))/(paprs(i,1)-30000.),1.)
133            enddo
134         enddo
135
136! Pour iflag_ratqs=1 ou 2, le ratqs est constant au dessus de
137! 300 hPa (ratqshaut), varie lineariement en fonction de la pression
138! entre 600 et 300 hPa et est soit constant (ratqsbas) pour iflag_ratqs=1
139! soit lineaire (entre 0 a la surface et ratqsbas) pour iflag_ratqs=2
140! Il s'agit de differents tests dans la phase de reglage du modele
141! avec thermiques.
142
143      else if (iflag_ratqs.eq.1) then
144
145         do k=1,klev
146            do i=1, klon
147               if (pplay(i,k).ge.60000.) then
148                  ratqss(i,k)=ratqsbas
149               else if ((pplay(i,k).ge.30000.).and.(pplay(i,k).lt.60000.)) then
150                  ratqss(i,k)=ratqsbas+(ratqshaut-ratqsbas)*(60000.-pplay(i,k))/(60000.-30000.)
151               else
152                  ratqss(i,k)=ratqshaut
153               endif
154            enddo
155         enddo
156
157      else if (iflag_ratqs.eq.2) then
158
159         do k=1,klev
160            do i=1, klon
161               if (pplay(i,k).ge.60000.) then
162                  ratqss(i,k)=ratqsbas*(paprs(i,1)-pplay(i,k))/(paprs(i,1)-60000.)
163               else if ((pplay(i,k).ge.30000.).and.(pplay(i,k).lt.60000.)) then
164                    ratqss(i,k)=ratqsbas+(ratqshaut-ratqsbas)*(60000.-pplay(i,k))/(60000.-30000.)
165               else
166                    ratqss(i,k)=ratqshaut
167               endif
168            enddo
169         enddo
170
171      else if (iflag_ratqs==3) then
172         do k=1,klev
173           ratqss(:,k)=ratqsbas+(ratqshaut-ratqsbas) &
174           *min( ((paprs(:,1)-pplay(:,k))/70000.)**2 , 1. )
175         enddo
176
177      else if (iflag_ratqs==4) then
178         do k=1,klev
179           ratqss(:,k)=ratqsbas+0.5*(ratqshaut-ratqsbas) &
180!          *( tanh( (50000.-pplay(:,k))/20000.) + 1.)
181           *( tanh( (ratqsp0-pplay(:,k))/ratqsdp) + 1.)
182         enddo
183
184
185      else if (iflag_ratqs==5) then
186! Dependency of ratqs on model resolution
187! Audran, Meryl, Lea, Gwendal and Etienne
188! April 2023
189         do k=1,klev
190            do i=1,klon
191              resol=sqrt(cell_area(i))
192              fact = sqrt(resol/resolmax)
193              ratqss(i,k)=ratqsbas*fact+0.5*(ratqshaut-ratqsbas)*fact &
194              *( tanh( (ratqsp0-pplay(i,k))/ratqsdp) + 1.)
195           enddo
196         enddo
197
198
199       else if (iflag_ratqs .GT. 9) then
200 
201       ! interactive ratqs calculations that depend on cold pools, orography, surface heterogeneity and small-scale turbulence
202       ! This should help getting a more realistic ratqs in the low and mid troposphere
203       ! We however need a "background" ratqs to account for subgrid distribution of qt (or qt/qs)
204       ! in the high troposphere
205       
206       ! background ratqs and initialisations
207          do k=1,klev
208             do i=1,klon
209              ratqss(i,k)=ratqsbas+0.5*(ratqshaut-ratqsbas) &
210              *( tanh( (ratqsp0-pplay(i,k))/ratqsdp) + 1.)
211              ratqss(i,k)=max(ratqss(i,k),0.0)
212              ratqs_hetero_(i,k)=0.
213              ratqs_oro_(i,k)=0.
214              ratqs_tke_(i,k)=0.
215              ratqs_inter_(i,k)=0
216             enddo
217          enddo
218     
219          if (iflag_ratqs .EQ. 10) then
220             print*,'avant ratqs_inter'
221            ! interactive ratqs with several sources
222             call ratqs_inter(klon,klev,iflag_ratqs,pdtphys,paprs, &
223                       ratqsbas,wake_deltaq,wake_s,q_seri,qtc_cv, sigt_cv, &
224                       fm_therm,entr_therm,detr_therm,detrain_cv,fm_cv,fqd,fqcomp,sigd, &
225                       ratqs_inter_,sigma_qtherm)
226             ratqss=ratqss+ratqs_inter_
227          else if (iflag_ratqs .EQ. 11) then
228            print*,'avant ratqs_inter'
229            ! interactive ratqs with several sources
230             call ratqs_inter(klon,klev,iflag_ratqs,pdtphys,paprs, &
231                       ratqsbas,wake_deltaq,wake_s,q_seri,qtc_cv, sigt_cv, &
232                       fm_therm,entr_therm,detr_therm,detrain_cv,fm_cv,fqd,fqcomp,sigd, &
233                       ratqs_inter_,sigma_qtherm)
234             ratqss=ratqss+ratqs_inter_
235          else if (iflag_ratqs .EQ. 12) then
236             ! contribution of surface heterogeneities to ratqs
237             call ratqs_hetero(klon,klev,pctsrf,s_pblh,t2m,q2m,t_seri,q_seri,pplay,paprs,ratqs_hetero_)
238             ratqss=ratqss+ratqs_hetero_
239          else if (iflag_ratqs .EQ. 13) then
240             ! contribution of ubgrid orography to ratqs
241             call ratqs_oro(klon,klev,pctsrf,zstd,zqsat,t_seri,pplay,paprs,ratqs_oro_)
242             ratqss=ratqss+ratqs_oro_
243          else if (iflag_ratqs .EQ. 14) then
244             ! effect of subgrid-scale TKE on ratqs (in development)
245             call ratqs_tke(klon,klev,pdtphys,t_seri,q_seri,zqsat,pplay,paprs,omega,tke,tke_dissip,lmix,wprime,ratqs_tke_)     
246             ratqss=ratqss+ratqs_tke_
247          endif
248         
249     
250      endif
251
252
253!  ratqs final
254!  -----------
255
256      if (iflag_cld_th.eq.1 .or.iflag_cld_th.eq.2.or.iflag_cld_th.eq.4) then
257
258! On ajoute une constante au ratqsc*2 pour tenir compte de
259! fluctuations turbulentes de petite echelle
260
261         do k=1,klev
262            do i=1,klon
263               if ((fm_therm(i,k)>1.e-10)) then
264                  ratqsc(i,k)=sqrt(ratqsc(i,k)**2+0.05**2)
265               endif
266            enddo
267         enddo
268
269!   les ratqs sont une combinaison de ratqss et ratqsc
270       if(prt_level.ge.9) write(lunout,*)'PHYLMD NOUVEAU TAU_RATQS ',tau_ratqs
271
272         if (tau_ratqs>1.e-10) then
273            facteur=exp(-pdtphys/tau_ratqs)
274         else
275            facteur=0.
276         endif
277         ratqs(:,:)=ratqsc(:,:)*(1.-facteur)+ratqs(:,:)*facteur
278!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
279! FH 22/09/2009
280! La ligne ci-dessous faisait osciller le modele et donnait une solution
281! assymptotique bidon et d??pendant fortement du pas de temps.
282!        ratqs(:,:)=sqrt(ratqs(:,:)**2+ratqss(:,:)**2)
283!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
284         ratqs(:,:)=max(ratqs(:,:),ratqss(:,:))
285      else if (iflag_cld_th<=6) then
286!   on ne prend que le ratqs stable pour fisrtilp
287         ratqs(:,:)=ratqss(:,:)
288      else
289          zfratqs1=exp(-pdtphys/10800.)
290          zfratqs2=exp(-pdtphys/10800.)
291          do k=1,klev
292             do i=1,klon
293                if (ratqsc(i,k).gt.1.e-10) then
294                   ratqs(i,k)=ratqs(i,k)*zfratqs2+(iflag_cld_th/100.)*ratqsc(i,k)*(1.-zfratqs2)
295                endif
296                ratqs(i,k)=min(ratqs(i,k)*zfratqs1+ratqss(i,k)*(1.-zfratqs1),0.5)
297             enddo
298          enddo
299      endif
300
301
302return
303END SUBROUTINE ratqs_main
304
305END MODULE lmdz_ratqs_main
Note: See TracBrowser for help on using the repository browser.