Context Navigation

calltherm.F90 @ 4134

Last change on this file since 4134 was 1399, checked in by idelkadi, 14 years ago

Modifications relatives a l'inclusion des nouveaux nuages d'Arnaud Jam

calltherm.F90 : remontee de variables suplementaires.

Passage de ztv, zpspsk, zthla, zthl en argument

fisrtilp.F : le programme de nuages avec appel a cloudth

Ajout des arguments : zqta, fraca, ztv, zpspsk, ztla, ztlh,
iflag_cldcon

physiq.F : appel modifie a calltherm et fisrtilp
thermcell_main.F90 : changement dans les flag_thermals_ed
thermcell_plume.F90 : version modifiee entrainement et detrainement

(on retrouve la precendente pour iflag_thermals_ed=10)

cloudth.F90 : PDF d'eau bi-gaussiennes

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 13.0 KB

Line
1	!
2	! $Id: calltherm.F90 1399 2010-06-08 08:29:11Z fairhead $
3	!
4	subroutine calltherm(dtime &
5	& ,pplay,paprs,pphi,weak_inversion &
6	& ,u_seri,v_seri,t_seri,q_seri,zqsat,debut &
7	& ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs &
8	& ,fm_therm,entr_therm,detr_therm,zqasc,clwcon0,lmax,ratqscth, &
9	& ratqsdiff,zqsatth,Ale_bl,Alp_bl,lalim_conv,wght_th, &
10	& zmax0,f0,zw2,fraca,ztv,zpspsk,ztla,zthl)
11
12	USE dimphy
13	implicit none
14	#include "dimensions.h"
15	!#include "dimphy.h"
16	#include "thermcell.h"
17	#include "iniprint.h"
18
19	!IM 140508
20	INTEGER itap
21	REAL dtime
22	LOGICAL debut
23	LOGICAL logexpr0, logexpr2(klon,klev), logexpr1(klon)
24	REAL fact(klon)
25	INTEGER nbptspb
26
27	REAL u_seri(klon,klev),v_seri(klon,klev)
28	REAL t_seri(klon,klev),q_seri(klon,klev),qmemoire(klon,klev)
29	REAL weak_inversion(klon)
30	REAL paprs(klon,klev+1)
31	REAL pplay(klon,klev)
32	REAL pphi(klon,klev)
33	real zlev(klon,klev+1)
34	!test: on sort lentr et a* pour alimenter KE
35	REAL wght_th(klon,klev)
36	INTEGER lalim_conv(klon)
37	REAL zw2(klon,klev+1),fraca(klon,klev+1)
38
39	!FH Update Thermiques
40	REAL d_t_ajs(klon,klev), d_q_ajs(klon,klev)
41	REAL d_u_ajs(klon,klev),d_v_ajs(klon,klev)
42	real fm_therm(klon,klev+1)
43	real entr_therm(klon,klev),detr_therm(klon,klev)
44
45	!********************************************************
46	! declarations
47	LOGICAL flag_bidouille_stratocu
48	real fmc_therm(klon,klev+1),zqasc(klon,klev)
49	real zqla(klon,klev)
50	real zqta(klon,klev)
51	real ztv(klon,klev)
52	real zpspsk(klon,klev)
53	real ztla(klon,klev)
54	real zthl(klon,klev)
55	real wmax_sec(klon)
56	real zmax_sec(klon)
57	real f_sec(klon)
58	real detrc_therm(klon,klev)
59	! FH WARNING : il semble que ces save ne servent a rien
60	! save fmc_therm, detrc_therm
61	real clwcon0(klon,klev)
62	real zqsat(klon,klev)
63	real zw_sec(klon,klev+1)
64	integer lmix_sec(klon)
65	integer lmax(klon)
66	real ratqscth(klon,klev)
67	real ratqsdiff(klon,klev)
68	real zqsatth(klon,klev)
69	!nouvelles variables pour la convection
70	real Ale_bl(klon)
71	real Alp_bl(klon)
72	real Ale(klon)
73	real Alp(klon)
74	!RC
75	!on garde le zmax du pas de temps precedent
76	real zmax0(klon), f0(klon)
77	!********************************************************
78
79
80	! variables locales
81	REAL d_t_the(klon,klev), d_q_the(klon,klev)
82	REAL d_u_the(klon,klev),d_v_the(klon,klev)
83	!
84	real zfm_therm(klon,klev+1),zdt
85	real zentr_therm(klon,klev),zdetr_therm(klon,klev)
86	! FH A VERIFIER : SAVE INUTILES
87	! save zentr_therm,zfm_therm
88
89	character (len=20) :: modname='calltherm'
90	character (len=80) :: abort_message
91
92	integer i,k
93	logical, save :: first=.true.
94	!$OMP THREADPRIVATE(first)
95	!********************************************************
96	if (first) then
97	itap=0
98	first=.false.
99	endif
100
101	! Incrementer le compteur de la physique
102	itap = itap + 1
103
104	! Modele du thermique
105	! ===================
106	! print*,'thermiques: WARNING on passe t au lieu de t_seri'
107
108
109	! On prend comme valeur initiale des thermiques la valeur du pas
110	! de temps precedent
111	zfm_therm(:,:)=fm_therm(:,:)
112	zdetr_therm(:,:)=detr_therm(:,:)
113	zentr_therm(:,:)=entr_therm(:,:)
114
115	! On reinitialise les flux de masse a zero pour le cumul en
116	! cas de splitting
117	fm_therm(:,:)=0.
118	entr_therm(:,:)=0.
119	detr_therm(:,:)=0.
120
121	Ale_bl(:)=0.
122	Alp_bl(:)=0.
123	if (prt_level.ge.10) then
124	print*,'thermV4 nsplit: ',nsplit_thermals,' weak_inversion'
125	endif
126
127	! tests sur les valeurs negatives de l'eau
128	logexpr0=prt_level.ge.10
129	nbptspb=0
130	do k=1,klev
131	do i=1,klon
132	! Attention teste abderr 19-03-09
133	! logexpr2(i,k)=.not.q_seri(i,k).ge.0.
134	logexpr2(i,k)=.not.q_seri(i,k).ge.1.e-15
135	if (logexpr2(i,k)) then
136	q_seri(i,k)=1.e-15
137	nbptspb=nbptspb+1
138	endif
139	! if (logexpr0) &
140	! & print*,'WARN eau<0 avant therm i=',i,' k=',k &
141	! & ,' dq,q',d_q_the(i,k),q_seri(i,k)
142	enddo
143	enddo
144	if(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
145
146	zdt=dtime/REAL(nsplit_thermals)
147	do isplit=1,nsplit_thermals
148
149	if (iflag_thermals.eq.1) then
150	CALL thermcell_2002(klon,klev,zdt &
151	& ,pplay,paprs,pphi &
152	& ,u_seri,v_seri,t_seri,q_seri &
153	& ,d_u_the,d_v_the,d_t_the,d_q_the &
154	& ,zfm_therm,zentr_therm &
155	& ,r_aspect_thermals,30.,w2di_thermals &
156	& ,tau_thermals,3)
157	else if (iflag_thermals.eq.2) then
158	CALL thermcell_sec(klon,klev,zdt &
159	& ,pplay,paprs,pphi,zlev &
160	& ,u_seri,v_seri,t_seri,q_seri &
161	& ,d_u_the,d_v_the,d_t_the,d_q_the &
162	& ,zfm_therm,zentr_therm &
163	& ,r_aspect_thermals,30.,w2di_thermals &
164	& ,tau_thermals,3)
165	else if (iflag_thermals.eq.3) then
166	CALL thermcell(klon,klev,zdt &
167	& ,pplay,paprs,pphi &
168	& ,u_seri,v_seri,t_seri,q_seri &
169	& ,d_u_the,d_v_the,d_t_the,d_q_the &
170	& ,zfm_therm,zentr_therm &
171	& ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
172	& ,tau_thermals,3)
173	else if (iflag_thermals.eq.10) then
174	CALL thermcell_eau(klon,klev,zdt &
175	& ,pplay,paprs,pphi &
176	& ,u_seri,v_seri,t_seri,q_seri &
177	& ,d_u_the,d_v_the,d_t_the,d_q_the &
178	& ,zfm_therm,zentr_therm &
179	& ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
180	& ,tau_thermals,3)
181	else if (iflag_thermals.eq.11) then
182	abort_message = 'cas non prevu dans calltherm'
183	CALL abort_gcm (modname,abort_message,1)
184
185	! CALL thermcell_pluie(klon,klev,zdt &
186	! & ,pplay,paprs,pphi,zlev &
187	! & ,u_seri,v_seri,t_seri,q_seri &
188	! & ,d_u_the,d_v_the,d_t_the,d_q_the &
189	! & ,zfm_therm,zentr_therm,zqla &
190	! & ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
191	! & ,tau_thermals,3)
192	else if (iflag_thermals.eq.12) then
193	CALL calcul_sec(klon,klev,zdt &
194	& ,pplay,paprs,pphi,zlev &
195	& ,u_seri,v_seri,t_seri,q_seri &
196	& ,zmax_sec,wmax_sec,zw_sec,lmix_sec &
197	& ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
198	& ,tau_thermals)
199	else if (iflag_thermals==13.or.iflag_thermals==14) then
200	CALL thermcellV0_main(itap,klon,klev,zdt &
201	& ,pplay,paprs,pphi,debut &
202	& ,u_seri,v_seri,t_seri,q_seri &
203	& ,d_u_the,d_v_the,d_t_the,d_q_the &
204	& ,zfm_therm,zentr_therm,zdetr_therm,zqasc,zqla,lmax &
205	& ,ratqscth,ratqsdiff,zqsatth &
206	& ,r_aspect_thermals,l_mix_thermals &
207	& ,tau_thermals,Ale,Alp,lalim_conv,wght_th &
208	& ,zmax0,f0,zw2,fraca)
209	else if (iflag_thermals==15.or.iflag_thermals==16) then
210
211	! print*,'THERM iflag_thermas_ed=',iflag_thermals_ed
212	CALL thermcell_main(itap,klon,klev,zdt &
213	& ,pplay,paprs,pphi,debut &
214	& ,u_seri,v_seri,t_seri,q_seri &
215	& ,d_u_the,d_v_the,d_t_the,d_q_the &
216	& ,zfm_therm,zentr_therm,zdetr_therm,zqasc,zqla,lmax &
217	& ,ratqscth,ratqsdiff,zqsatth &
218	& ,r_aspect_thermals,l_mix_thermals &
219	& ,tau_thermals,iflag_thermals_ed,Ale,Alp,lalim_conv,wght_th &
220	& ,zmax0,f0,zw2,fraca,ztv,zpspsk &
221	& ,ztla,zthl)
222	if (prt_level.gt.10) write(lunout,*)'Apres thermcell_main OK'
223	else
224	abort_message = 'Cas des thermiques non prevu'
225	CALL abort_gcm (modname,abort_message,1)
226	endif
227
228	flag_bidouille_stratocu=iflag_thermals.eq.14.or.iflag_thermals.eq.16
229
230	fact(:)=0.
231	DO i=1,klon
232	logexpr1(i)=flag_bidouille_stratocu.or.weak_inversion(i).gt.0.5
233	IF(logexpr1(i)) fact(i)=1./REAL(nsplit_thermals)
234	ENDDO
235
236	DO k=1,klev
237	! transformation de la derivee en tendance
238	d_t_the(:,k)=d_t_the(:,k)dtimefact(:)
239	d_u_the(:,k)=d_u_the(:,k)dtimefact(:)
240	d_v_the(:,k)=d_v_the(:,k)dtimefact(:)
241	d_q_the(:,k)=d_q_the(:,k)dtimefact(:)
242	fm_therm(:,k)=fm_therm(:,k) &
243	& +zfm_therm(:,k)*fact(:)
244	entr_therm(:,k)=entr_therm(:,k) &
245	& +zentr_therm(:,k)*fact(:)
246	detr_therm(:,k)=detr_therm(:,k) &
247	& +zdetr_therm(:,k)*fact(:)
248	ENDDO
249	fm_therm(:,klev+1)=0.
250
251
252
253	! accumulation de la tendance
254	d_t_ajs(:,:)=d_t_ajs(:,:)+d_t_the(:,:)
255	d_u_ajs(:,:)=d_u_ajs(:,:)+d_u_the(:,:)
256	d_v_ajs(:,:)=d_v_ajs(:,:)+d_v_the(:,:)
257	d_q_ajs(:,:)=d_q_ajs(:,:)+d_q_the(:,:)
258
259	! incrementation des variables meteo
260	t_seri(:,:) = t_seri(:,:) + d_t_the(:,:)
261	u_seri(:,:) = u_seri(:,:) + d_u_the(:,:)
262	v_seri(:,:) = v_seri(:,:) + d_v_the(:,:)
263	qmemoire(:,:)=q_seri(:,:)
264	q_seri(:,:) = q_seri(:,:) + d_q_the(:,:)
265	if (prt_level.gt.10) write(lunout,*)'Apres apres thermcell_main OK'
266
267	DO i=1,klon
268	if(prt_level.GE.10) print*,'calltherm i Alp_bl Alp Ale_bl Ale',i,Alp_bl(i),Alp(i),Ale_bl(i),Ale(i)
269	fm_therm(i,klev+1)=0.
270	Ale_bl(i)=Ale_bl(i)+Ale(i)/REAL(nsplit_thermals)
271	! write(22,*)'ALE CALLTHERM',Ale_bl(i),Ale(i)
272	Alp_bl(i)=Alp_bl(i)+Alp(i)/REAL(nsplit_thermals)
273	! write(23,*)'ALP CALLTHERM',Alp_bl(i),Alp(i)
274	ENDDO
275
276	!IM 060508 marche pas comme cela !!! enddo ! isplit
277
278	! tests sur les valeurs negatives de l'eau
279	nbptspb=0
280	DO k = 1, klev
281	DO i = 1, klon
282	logexpr2(i,k)=.not.q_seri(i,k).ge.0.
283	if (logexpr2(i,k)) then
284	q_seri(i,k)=1.e-15
285	nbptspb=nbptspb+1
286	! if (prt_level.ge.10) then
287	! print*,'WARN eau<0 apres therm i=',i,' k=',k &
288	! & ,' dq,q',d_q_the(i,k),q_seri(i,k), &
289	! & 'fm=',zfm_therm(i,k),'entr=',entr_therm(i,k)
290	endif
291	ENDDO
292	ENDDO
293	IF(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
294	! tests sur les valeurs de la temperature
295	nbptspb=0
296	DO k = 1, klev
297	DO i = 1, klon
298	logexpr2(i,k)=t_seri(i,k).lt.50..or.t_seri(i,k).gt.370.
299	if (logexpr2(i,k)) nbptspb=nbptspb+1
300	! if ((t_seri(i,k).lt.50.) .or. &
301	! & (t_seri(i,k).gt.370.)) then
302	! print*,'WARN temp apres therm i=',i,' k=',k &
303	! & ,' t_seri',t_seri(i,k)
304	! CALL abort
305	! endif
306	ENDDO
307	ENDDO
308	IF(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
309	enddo ! isplit
310
311	!
312	!***************************************************************
313	! calcul du flux ascencant conservatif
314	! print*,'<<<<calcul flux ascendant conservatif'
315
316	fmc_therm=0.
317	do k=1,klev
318	do i=1,klon
319	if (entr_therm(i,k).gt.0.) then
320	fmc_therm(i,k+1)=fmc_therm(i,k)+entr_therm(i,k)
321	else
322	fmc_therm(i,k+1)=fmc_therm(i,k)
323	endif
324	detrc_therm(i,k)=(fmc_therm(i,k+1)-fm_therm(i,k+1)) &
325	& -(fmc_therm(i,k)-fm_therm(i,k))
326	enddo
327	enddo
328
329
330	!****************************************************************
331	! calcul de l'humidite dans l'ascendance
332	! print*,'<<<<calcul de lhumidite dans thermique'
333	!CR:on ne le calcule que pour le cas sec
334	if (iflag_thermals.le.11) then
335	do i=1,klon
336	zqasc(i,1)=q_seri(i,1)
337	do k=2,klev
338	if (fmc_therm(i,k+1).gt.1.e-6) then
339	zqasc(i,k)=(fmc_therm(i,k)*zqasc(i,k-1) &
340	& +entr_therm(i,k)*q_seri(i,k))/fmc_therm(i,k+1)
341	!CR:test on asseche le thermique
342	! zqasc(i,k)=zqasc(i,k)/2.
343	! else
344	! zqasc(i,k)=q_seri(i,k)
345	endif
346	enddo
347	enddo
348
349
350	! calcul de l'eau condensee dans l'ascendance
351	! print*,'<<<<calcul de leau condensee dans thermique'
352	do i=1,klon
353	do k=1,klev
354	clwcon0(i,k)=zqasc(i,k)-zqsat(i,k)
355	if (clwcon0(i,k).lt.0. .or. &
356	& (fm_therm(i,k+1)+detrc_therm(i,k)).lt.1.e-6) then
357	clwcon0(i,k)=0.
358	endif
359	enddo
360	enddo
361	else
362	do i=1,klon
363	do k=1,klev
364	clwcon0(i,k)=zqla(i,k)
365	if (clwcon0(i,k).lt.0. .or. &
366	& (fm_therm(i,k+1)+detrc_therm(i,k)).lt.1.e-6) then
367	clwcon0(i,k)=0.
368	endif
369	enddo
370	enddo
371	endif
372	!*******************************************************************
373
374
375	return
376
377	end

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

Download in other formats:

Original Format