Context Navigation

lmdz_thermcell_alp.F90 @ 5136

Last change on this file since 5136 was 5134, checked in by abarral, 8 weeks ago
Replace academic.h, alpale.h, comdissip.h, comdissipn.h, comdissnew.h by modules Remove unused clesph0.h
File size: 15.0 KB

Line
1	MODULE lmdz_thermcell_alp
2	! $Id: thermcell_main.F90 2351 2015-08-25 15:14:59Z emillour $
3
4	CONTAINS
5
6	SUBROUTINE thermcell_alp(ngrid,nlay,ptimestep & ! in
7	,pplay,pplev & ! in
8	,fm0,entr0,lmax & ! in
9	,pbl_tke,pctsrf,omega,airephy & ! in
10	,zw2,fraca & ! in
11	,pcon,rhobarz,wth3,wmax_sec,lalim,fm,alim_star,zmax & ! in
12
13	,zcong,ale_bl,alp_bl,lalim_conv,wght_th & ! out
14	,zlcl,fraca0,w0,w_conv,therm_tke_max0,env_tke_max0 & ! out
15	,n2,s2,strig,ale_bl_stat & ! out
16	,therm_tke_max,env_tke_max & ! out
17	,alp_bl_det,alp_bl_fluct_m,alp_bl_fluct_tke & ! out
18	,alp_bl_conv,alp_bl_stat & ! out
19	)
20
21	USE indice_sol_mod
22	USE lmdz_thermcell_main, ONLY: thermcell_tke_transport
23	USE lmdz_alpale
24
25	IMPLICIT NONE
26
27	!=======================================================================
28
29	! Auteurs: Catherine Rio
30	! Modifications :
31	! Nicolas Rochetin et Jean-Yves Grandpeix
32	! pour la fermeture stochastique. 2012
33	! Frédéric Hourdin :
34	! netoyage informatique. 2022
35
36	!=======================================================================
37	!-----------------------------------------------------------------------
38	! declarations:
39	! -------------
40
41	INCLUDE "YOMCST.h"
42	INCLUDE "YOETHF.h"
43	INCLUDE "FCTTRE.h"
44
45	! arguments:
46	! ----------
47
48	!------Entrees
49	INTEGER, INTENT(IN) :: ngrid,nlay
50	REAL, INTENT(IN) :: ptimestep
51	REAL, INTENT(IN) :: pplay(ngrid,nlay),pplev(ngrid,nlay+1)
52	INTEGER, INTENT(IN), DIMENSION(ngrid) ::lmax,lalim
53	REAL, INTENT(IN), DIMENSION(ngrid) :: zmax, zcong
54	REAL, INTENT(IN), DIMENSION(ngrid,nlay+1) :: zw2
55	REAL, INTENT(IN), DIMENSION(ngrid,nlay+1) :: fraca
56	REAL, INTENT(IN), DIMENSION(ngrid,nlay) :: wth3
57	REAL, INTENT(IN), DIMENSION(ngrid,nlay) :: rhobarz
58	REAL, INTENT(IN), DIMENSION(ngrid) :: wmax_sec
59	REAL, INTENT(IN), DIMENSION(ngrid,nlay) :: entr0
60	REAL, INTENT(IN), DIMENSION(ngrid,nlay+1) :: fm0,fm
61	REAL, INTENT(IN), DIMENSION(ngrid) :: pcon
62	REAL, INTENT(IN), DIMENSION(ngrid,nlay) :: alim_star
63	REAL, INTENT(IN), DIMENSION(ngrid,nlay+1,nbsrf) :: pbl_tke
64	REAL, INTENT(IN), DIMENSION(ngrid,nbsrf) :: pctsrf
65	REAL, INTENT(IN), DIMENSION(ngrid,nlay) :: omega
66	REAL, INTENT(IN), DIMENSION(ngrid) :: airephy
67	!------Sorties
68	REAL, INTENT(OUT), DIMENSION(ngrid) :: ale_bl,alp_bl
69	REAL, INTENT(OUT), DIMENSION(ngrid,nlay) :: wght_th
70	INTEGER, INTENT(OUT), DIMENSION(ngrid) :: lalim_conv
71	REAL, INTENT(OUT), DIMENSION(ngrid) :: zlcl,fraca0,w0,w_conv
72	REAL, INTENT(OUT), DIMENSION(ngrid) :: therm_tke_max0,env_tke_max0,n2,s2,ale_bl_stat,strig
73	REAL, INTENT(OUT), DIMENSION(ngrid,nlay) :: therm_tke_max,env_tke_max
74	REAL, INTENT(OUT), DIMENSION(ngrid) :: alp_bl_det,alp_bl_fluct_m,alp_bl_fluct_tke
75	REAL, INTENT(OUT), DIMENSION(ngrid) :: alp_bl_conv,alp_bl_stat
76
77	!=============================================================================================
78	!------Local
79	!=============================================================================================
80
81	REAL susqr2pi, reuler
82	INTEGER ig,k,l
83	INTEGER nsrf
84	REAL rhobarz0(ngrid) ! Densité au LCL
85	LOGICAL ok_lcl(ngrid) ! Existence du LCL des thermiques
86	INTEGER klcl(ngrid) ! Niveau du LCL
87	REAL interp(ngrid) ! Coef d'interpolation pour le LCL
88	!--Triggering
89	REAL, parameter :: su_cst=4e4 ! Surface unite: celle d'un updraft élémentaire
90	REAL, parameter :: hcoef=1 ! Coefficient directeur pour le calcul de s2
91	REAL, parameter :: hmincoef=0.3 ! Coefficient directeur pour l'ordonnée à l'origine pour le calcul de s2
92	REAL, parameter :: eps1=0.3 ! Fraction de surface occupée par la population 1 : eps1=n1s1/(fraca0Sd)
93	REAL, DIMENSION(ngrid) :: hmin ! Ordonnée à l'origine pour le calcul de s2
94	REAL, DIMENSION(ngrid) :: zmax_moy ! Hauteur moyenne des thermiques : zmax_moy = zlcl + 0.33 (zmax-zlcl)
95	REAL, parameter :: zmax_moy_coef=0.33
96	REAL, DIMENSION(ngrid) :: depth ! Epaisseur moyenne du cumulus
97	REAL, DIMENSION(ngrid) :: zcong_moy
98	REAL, DIMENSION(ngrid) :: w_max ! Vitesse max statistique
99	REAL, DIMENSION(ngrid) :: s_max(ngrid)
100	!--Closure
101	REAL, DIMENSION(ngrid,nlay) :: pbl_tke_max ! Profil de TKE moyenne
102	REAL, DIMENSION(ngrid) :: pbl_tke_max0 ! TKE moyenne au LCL
103	REAL, DIMENSION(ngrid,nlay) :: w_ls ! Vitesse verticale grande échelle (m/s)
104	REAL, parameter :: coef_m=1. ! On considère un rendement pour alp_bl_fluct_m
105	REAL, parameter :: coef_tke=1. ! On considère un rendement pour alp_bl_fluct_tke
106	REAL :: zdp
107	REAL, DIMENSION(ngrid) :: alp_int,dp_int
108	REAL, DIMENSION(ngrid) :: fm_tot
109
110	!------------------------------------------------------------
111	! Initialize output arrays related to stochastic triggering
112	!------------------------------------------------------------
113	DO ig = 1,ngrid
114	zlcl(ig) = 0.
115	fraca0(ig) = 0.
116	w0(ig) = 0.
117	w_conv(ig) = 0.
118	therm_tke_max0(ig) = 0.
119	env_tke_max0(ig) = 0.
120	n2(ig) = 0.
121	s2(ig) = 0.
122	ale_bl_stat(ig) = 0.
123	strig(ig) = 0.
124	alp_bl_det(ig) = 0.
125	alp_bl_fluct_m(ig) = 0.
126	alp_bl_fluct_tke(ig) = 0.
127	alp_bl_conv(ig) = 0.
128	alp_bl_stat(ig) = 0.
129	ENDDO
130	DO l = 1,nlay
131	DO ig = 1,ngrid
132	therm_tke_max(ig,l) = 0.
133	env_tke_max(ig,l) = 0.
134	ENDDO
135	ENDDO
136
137	!------------Test sur le LCL des thermiques
138	do ig=1,ngrid
139	ok_lcl(ig)=.FALSE.
140	IF ( (pcon(ig) > pplay(ig,nlay-1)) .AND. (pcon(ig) < pplay(ig,1)) ) ok_lcl(ig)=.TRUE.
141	enddo
142
143	!------------Localisation des niveaux entourant le LCL et du coef d'interpolation
144	do l=1,nlay-1
145	do ig=1,ngrid
146	IF (ok_lcl(ig)) THEN
147	!ATTENTION,zw2 calcule en pplev
148	! if ((pplay(ig,l) .ge. pcon(ig)) .AND. (pplay(ig,l+1) .le. pcon(ig))) THEN
149	! klcl(ig)=l
150	! interp(ig)=(pcon(ig)-pplay(ig,klcl(ig)))/(pplay(ig,klcl(ig)+1)-pplay(ig,klcl(ig)))
151	! endif
152	IF ((pplev(ig,l) >= pcon(ig)) .AND. (pplev(ig,l+1) <= pcon(ig))) THEN
153	klcl(ig)=l
154	interp(ig)=(pcon(ig)-pplev(ig,klcl(ig)))/(pplev(ig,klcl(ig)+1)-pplev(ig,klcl(ig)))
155	endif
156	endif
157	enddo
158	enddo
159
160	do ig =1,ngrid
161	!CR:REHABILITATION ZMAX CONTINU
162	IF (ok_lcl(ig)) THEN
163	rhobarz0(ig)=rhobarz(ig,klcl(ig))+(rhobarz(ig,klcl(ig)+1) &
164	-rhobarz(ig,klcl(ig)))*interp(ig)
165	zlcl(ig)=(pplev(ig,1)-pcon(ig))/(rhobarz0(ig)*RG)
166	zlcl(ig)=min(zlcl(ig),zmax(ig)) ! Si zlcl > zmax alors on pose zlcl = zmax
167	else
168	rhobarz0(ig)=0.
169	zlcl(ig)=zmax(ig)
170	endif
171	enddo
172	!!jyg fin
173
174	!------------Calcul des propriétés du thermique au LCL
175	IF ( (iflag_trig_bl>=1) .OR. (iflag_clos_bl>=1) ) THEN
176
177	!-----Initialisation de la TKE moyenne
178	do l=1,nlay
179	do ig=1,ngrid
180	pbl_tke_max(ig,l)=0.
181	enddo
182	enddo
183
184	!-----Calcul de la TKE moyenne
185	do nsrf=1,nbsrf
186	do l=1,nlay
187	do ig=1,ngrid
188	pbl_tke_max(ig,l)=pctsrf(ig,nsrf)*pbl_tke(ig,l,nsrf)+pbl_tke_max(ig,l)
189	enddo
190	enddo
191	enddo
192
193	!-----Initialisations des TKE dans et hors des thermiques
194	do l=1,nlay
195	do ig=1,ngrid
196	therm_tke_max(ig,l)=pbl_tke_max(ig,l)
197	env_tke_max(ig,l)=pbl_tke_max(ig,l)
198	enddo
199	enddo
200
201	!-----Calcul de la TKE transportée par les thermiques : therm_tke_max
202	CALL thermcell_tke_transport(ngrid,nlay,ptimestep,fm0,entr0, & ! in
203	rg,pplev,therm_tke_max) ! out
204	! print *,' thermcell_tke_transport -> ' !!jyg
205
206	!-----Calcul des profils verticaux de TKE hors thermiques : env_tke_max, et de la vitesse verticale grande échelle : W_ls
207	do l=1,nlay
208	do ig=1,ngrid
209	pbl_tke_max(ig,l)=fraca(ig,l)therm_tke_max(ig,l)+(1.-fraca(ig,l))env_tke_max(ig,l) ! Recalcul de TKE moyenne après transport de TKE_TH
210	env_tke_max(ig,l)=(pbl_tke_max(ig,l)-fraca(ig,l)*therm_tke_max(ig,l))/(1.-fraca(ig,l)) ! Recalcul de TKE dans l'environnement après transport de TKE_TH
211	w_ls(ig,l)=-1.omega(ig,l)/(RGrhobarz(ig,l)) ! Vitesse verticale de grande échelle
212	enddo
213	enddo
214	! print *,' apres w_ls = ' !!jyg
215
216	do ig=1,ngrid
217	IF (ok_lcl(ig)) THEN
218	fraca0(ig)=fraca(ig,klcl(ig))+(fraca(ig,klcl(ig)+1) &
219	-fraca(ig,klcl(ig)))*interp(ig)
220	w0(ig)=zw2(ig,klcl(ig))+(zw2(ig,klcl(ig)+1) &
221	-zw2(ig,klcl(ig)))*interp(ig)
222	w_conv(ig)=w_ls(ig,klcl(ig))+(w_ls(ig,klcl(ig)+1) &
223	-w_ls(ig,klcl(ig)))*interp(ig)
224	therm_tke_max0(ig)=therm_tke_max(ig,klcl(ig)) &
225	+(therm_tke_max(ig,klcl(ig)+1)-therm_tke_max(ig,klcl(ig)))*interp(ig)
226	env_tke_max0(ig)=env_tke_max(ig,klcl(ig))+(env_tke_max(ig,klcl(ig)+1) &
227	-env_tke_max(ig,klcl(ig)))*interp(ig)
228	pbl_tke_max0(ig)=pbl_tke_max(ig,klcl(ig))+(pbl_tke_max(ig,klcl(ig)+1) &
229	-pbl_tke_max(ig,klcl(ig)))*interp(ig)
230	IF (therm_tke_max0(ig)>=20.) therm_tke_max0(ig)=20.
231	IF (env_tke_max0(ig)>=20.) env_tke_max0(ig)=20.
232	IF (pbl_tke_max0(ig)>=20.) pbl_tke_max0(ig)=20.
233	else
234	fraca0(ig)=0.
235	w0(ig)=0.
236	!!jyg le 27/04/2012
237	!! zlcl(ig)=0.
238	!!
239	endif
240	enddo
241
242	ENDIF ! IF ( (iflag_trig_bl.ge.1) .OR. (iflag_clos_bl.ge.1) )
243	! print *,'ENDIF ( (iflag_trig_bl.ge.1) .OR. (iflag_clos_bl.ge.1) ) ' !!jyg
244
245	!------------Triggering------------------
246	IF (iflag_trig_bl>=1) THEN
247
248	!-----Initialisations
249	depth(:)=0.
250	n2(:)=0.
251	s2(:)=100. ! some low value, arbitrary
252	s_max(:)=0.
253
254
255
256	!-----Epaisseur du nuage (depth) et détermination de la queue du spectre de panaches (n2,s2) et du panache le plus gros (s_max)
257	do ig=1,ngrid
258	zmax_moy(ig)=zlcl(ig)+zmax_moy_coef*(zmax(ig)-zlcl(ig))
259	depth(ig)=zmax_moy(ig)-zlcl(ig)
260	hmin(ig)=hmincoef*zlcl(ig)
261	IF (depth(ig)>=10.) THEN
262	s2(ig)=(hcoefdepth(ig)+hmin(ig))*2
263	n2(ig)=(1.-eps1)fraca0(ig)airephy(ig)/s2(ig)
264	!!
265	!!jyg le 27/04/2012
266	!! s_max(ig)=s2(ig)*log(n2(ig))
267	!! if (n2(ig) .lt. 1) s_max(ig)=0.
268	s_max(ig)=s2(ig)*log(max(n2(ig),1.))
269	!!fin jyg
270	else
271	n2(ig)=0.
272	s_max(ig)=0.
273	endif
274	enddo
275	! print *,'avant Calcul de Wmax ' !!jyg
276
277	!CR: calcul de strig
278	IF (iflag_strig==0) THEN
279	strig(:)=s_trig
280	ELSE IF (iflag_strig==1) THEN
281	do ig=1,ngrid
282	! zcong_moy(ig)=zlcl(ig)+zmax_moy_coef*(zcong(ig)-zlcl(ig))
283	! strig(ig)=(hcoef(zcong_moy(ig)-zlcl(ig))+hmin(ig))*2
284	strig(ig)=(zcong(ig)-zlcl(ig))**2
285	enddo
286	ELSE IF (iflag_strig==2) THEN
287	do ig=1,ngrid
288	IF (h_trig>zlcl(ig)) THEN
289	strig(ig)=(h_trig-zlcl(ig))**2
290	else
291	strig(ig)=s_trig
292	endif
293	enddo
294	endif
295
296	susqr2pi=su_cstsqrt(2.Rpi)
297	reuler=exp(1.)
298	do ig=1,ngrid
299	IF ( (depth(ig)>=10.) .AND. (s_max(ig)>susqr2pi*reuler) ) THEN
300	w_max(ig)=w0(ig)(1.+sqrt(2.log(s_max(ig)/susqr2pi)-log(2.*log(s_max(ig)/susqr2pi))))
301	ale_bl_stat(ig)=0.5w_max(ig)*2
302	else
303	w_max(ig)=0.
304	ale_bl_stat(ig)=0.
305	endif
306	enddo
307
308	ENDIF ! iflag_trig_bl
309	! print *,'ENDIF iflag_trig_bl' !!jyg
310
311	!------------Closure------------------
312
313	IF (iflag_clos_bl>=2) THEN
314
315	!-----Calcul de ALP_BL_STAT
316	do ig=1,ngrid
317	alp_bl_det(ig)=0.5coef_mrhobarz0(ig)(w0(ig)3)fraca0(ig)(1.-2.fraca0(ig))/((1.-fraca0(ig))**2)
318	alp_bl_fluct_m(ig)=1.5rhobarz0(ig)fraca0(ig)(w_conv(ig)+coef_mw0(ig))* &
319	(w0(ig)**2)
320	alp_bl_fluct_tke(ig)=3.coef_mrhobarz0(ig)w0(ig)fraca0(ig)*(therm_tke_max0(ig)-env_tke_max0(ig)) &
321	+3.rhobarz0(ig)w_conv(ig)*pbl_tke_max0(ig)
322	IF (iflag_clos_bl>=2) THEN
323	alp_bl_conv(ig)=1.5coef_mrhobarz0(ig)fraca0(ig)(fraca0(ig)/(1.-fraca0(ig)))w_conv(ig) &
324	(w0(ig)**2)
325	else
326	alp_bl_conv(ig)=0.
327	endif
328	alp_bl_stat(ig)=alp_bl_det(ig)+alp_bl_fluct_m(ig)+alp_bl_fluct_tke(ig)+alp_bl_conv(ig)
329	enddo
330
331	!-----Sécurité ALP infinie
332	do ig=1,ngrid
333	IF (fraca0(ig)>0.98) alp_bl_stat(ig)=2.
334	enddo
335
336	ENDIF ! (iflag_clos_bl.ge.2)
337
338	!!! fin nrlmd le 10/04/2012
339
340	! PRINT*,'avant calcul ale et alp'
341	!calcul de ALE et ALP pour la convection
342	alp_bl(:)=0.
343	ale_bl(:)=0.
344	! PRINT*,'ALE,ALP ,l,zw2(ig,l),ale_bl(ig),alp_bl(ig)'
345	do l=1,nlay
346	do ig=1,ngrid
347	alp_bl(ig)=max(alp_bl(ig),0.5rhobarz(ig,l)wth3(ig,l) )
348	ale_bl(ig)=max(ale_bl(ig),0.5zw2(ig,l)*2)
349	! PRINT*,'ALE,ALP',l,zw2(ig,l),ale_bl(ig),alp_bl(ig)
350	enddo
351	enddo
352
353	! ale sec (max de wmax/2 sous la zone d'inhibition) dans
354	! le cas iflag_trig_bl=3
355	IF (iflag_trig_bl==3) ale_bl(:)=0.5wmax_sec(:)*2
356
357	!test:calcul de la ponderation des couches pour KE
358	!initialisations
359
360	fm_tot(:)=0.
361	wght_th(:,:)=1.
362	lalim_conv(:)=lalim(:)
363
364	do k=1,nlay
365	do ig=1,ngrid
366	IF (k<=lalim_conv(ig)) fm_tot(ig)=fm_tot(ig)+fm(ig,k)
367	enddo
368	enddo
369
370	! assez bizarre car, si on est dans la couche d'alim et que alim_star et
371	! plus petit que 1.e-10, on prend wght_th=1.
372	do k=1,nlay
373	do ig=1,ngrid
374	IF (k<=lalim_conv(ig).AND.alim_star(ig,k)>1.e-10) THEN
375	wght_th(ig,k)=alim_star(ig,k)
376	endif
377	enddo
378	enddo
379
380	! PRINT*,'apres wght_th'
381	!test pour prolonger la convection
382	do ig=1,ngrid
383	!v1d if ((alim_star(ig,1).lt.1.e-10).AND.(therm)) THEN
384	IF ((alim_star(ig,1)<1.e-10)) THEN
385	lalim_conv(ig)=1
386	wght_th(ig,1)=1.
387	! PRINT*,'lalim_conv ok',lalim_conv(ig),wght_th(ig,1)
388	endif
389	enddo
390
391	!------------------------------------------------------------------------
392	! Modif CR/FH 20110310 : alp integree sur la verticale.
393	! Integrale verticale de ALP.
394	! wth3 etant aux niveaux inter-couches, on utilise d play comme masse des
395	! couches
396	!------------------------------------------------------------------------
397
398	alp_int(:)=0.
399	dp_int(:)=0.
400	do l=2,nlay
401	do ig=1,ngrid
402	IF(l<=lmax(ig)) THEN
403	zdp=pplay(ig,l-1)-pplay(ig,l)
404	alp_int(ig)=alp_int(ig)+0.5rhobarz(ig,l)wth3(ig,l)*zdp
405	dp_int(ig)=dp_int(ig)+zdp
406	endif
407	enddo
408	enddo
409
410	IF (iflag_coupl>=3 .AND. iflag_coupl<=5) THEN
411	do ig=1,ngrid
412	!valeur integree de alp_bl * 0.5:
413	IF (dp_int(ig)>0.) THEN
414	alp_bl(ig)=alp_int(ig)/dp_int(ig)
415	endif
416	enddo!
417	endif
418
419
420	! Facteur multiplicatif sur alp_bl
421	alp_bl(:)=alp_bl_k*alp_bl(:)
422
423	!------------------------------------------------------------------------
424
425
426
427	RETURN
428	END
429	END MODULE lmdz_thermcell_alp

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format