source: LMDZ6/branches/Amaury_dev/libf/phylmd/lmdz_ratqs_main.F90 @ 5411

Last change on this file since 5411 was 5226, checked in by abarral, 3 months ago

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