Context Navigation

lmdz_thermcell_alp.F90 @ 5117

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

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

Download in other formats:

Original Format