Context Navigation

calltherm.F90 @ 979

Last change on this file since 979 was 973, checked in by lmdzadmin, 16 years ago
Initialisations : concvl, cv3_routines, cva_driver, physiq Correction bug i0 + ajout tests : cv3p1_closure Ajout sorties : ale, alp, cin, wape Ajout variables wake : phyetat0, phyredem IM
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 12.2 KB

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

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

Download in other formats:

Original Format