Context Navigation

calltherm.F90 @ 1550

Last change on this file since 1550 was 1517, checked in by jghattas, 15 years ago

Errors in argument list not corresponding between CALL and SUBROUTINE in some routines. These errors were found with compile option -qextchk at vargas.

Removed option -qextchk from debug compile options at vargas. This option is not usable with mpi codes.

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

File size: 13.2 KB

Line
1	!
2	! $Id: calltherm.F90 1517 2011-05-10 13:07:35Z lguez $
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)
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)
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)
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)
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,iflag_coupl &
220	& ,Ale,Alp,lalim_conv,wght_th &
221	& ,zmax0,f0,zw2,fraca,ztv,zpspsk &
222	& ,ztla,zthl)
223	if (prt_level.gt.10) write(lunout,*)'Apres thermcell_main OK'
224	else
225	abort_message = 'Cas des thermiques non prevu'
226	CALL abort_gcm (modname,abort_message,1)
227	endif
228
229	flag_bidouille_stratocu=iflag_thermals.eq.14.or.iflag_thermals.eq.16
230
231	fact(:)=0.
232	DO i=1,klon
233	logexpr1(i)=flag_bidouille_stratocu.or.weak_inversion(i).gt.0.5
234	IF(logexpr1(i)) fact(i)=1./REAL(nsplit_thermals)
235	ENDDO
236
237	DO k=1,klev
238	! transformation de la derivee en tendance
239	d_t_the(:,k)=d_t_the(:,k)dtimefact(:)
240	d_u_the(:,k)=d_u_the(:,k)dtimefact(:)
241	d_v_the(:,k)=d_v_the(:,k)dtimefact(:)
242	d_q_the(:,k)=d_q_the(:,k)dtimefact(:)
243	fm_therm(:,k)=fm_therm(:,k) &
244	& +zfm_therm(:,k)*fact(:)
245	entr_therm(:,k)=entr_therm(:,k) &
246	& +zentr_therm(:,k)*fact(:)
247	detr_therm(:,k)=detr_therm(:,k) &
248	& +zdetr_therm(:,k)*fact(:)
249	ENDDO
250	fm_therm(:,klev+1)=0.
251
252
253
254	! accumulation de la tendance
255	d_t_ajs(:,:)=d_t_ajs(:,:)+d_t_the(:,:)
256	d_u_ajs(:,:)=d_u_ajs(:,:)+d_u_the(:,:)
257	d_v_ajs(:,:)=d_v_ajs(:,:)+d_v_the(:,:)
258	d_q_ajs(:,:)=d_q_ajs(:,:)+d_q_the(:,:)
259
260	! incrementation des variables meteo
261	t_seri(:,:) = t_seri(:,:) + d_t_the(:,:)
262	u_seri(:,:) = u_seri(:,:) + d_u_the(:,:)
263	v_seri(:,:) = v_seri(:,:) + d_v_the(:,:)
264	qmemoire(:,:)=q_seri(:,:)
265	q_seri(:,:) = q_seri(:,:) + d_q_the(:,:)
266	if (prt_level.gt.10) write(lunout,*)'Apres apres thermcell_main OK'
267
268	DO i=1,klon
269	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)
270	fm_therm(i,klev+1)=0.
271	Ale_bl(i)=Ale_bl(i)+Ale(i)/REAL(nsplit_thermals)
272	! write(22,*)'ALE CALLTHERM',Ale_bl(i),Ale(i)
273	Alp_bl(i)=Alp_bl(i)+Alp(i)/REAL(nsplit_thermals)
274	! write(23,*)'ALP CALLTHERM',Alp_bl(i),Alp(i)
275	ENDDO
276
277	!IM 060508 marche pas comme cela !!! enddo ! isplit
278
279	! tests sur les valeurs negatives de l'eau
280	nbptspb=0
281	DO k = 1, klev
282	DO i = 1, klon
283	logexpr2(i,k)=.not.q_seri(i,k).ge.0.
284	if (logexpr2(i,k)) then
285	q_seri(i,k)=1.e-15
286	nbptspb=nbptspb+1
287	! if (prt_level.ge.10) then
288	! print*,'WARN eau<0 apres therm i=',i,' k=',k &
289	! & ,' dq,q',d_q_the(i,k),q_seri(i,k), &
290	! & 'fm=',zfm_therm(i,k),'entr=',entr_therm(i,k)
291	endif
292	ENDDO
293	ENDDO
294	IF(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
295	! tests sur les valeurs de la temperature
296	nbptspb=0
297	DO k = 1, klev
298	DO i = 1, klon
299	logexpr2(i,k)=t_seri(i,k).lt.50..or.t_seri(i,k).gt.370.
300	if (logexpr2(i,k)) nbptspb=nbptspb+1
301	! if ((t_seri(i,k).lt.50.) .or. &
302	! & (t_seri(i,k).gt.370.)) then
303	! print*,'WARN temp apres therm i=',i,' k=',k &
304	! & ,' t_seri',t_seri(i,k)
305	! CALL abort
306	! endif
307	ENDDO
308	ENDDO
309	IF(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
310	enddo ! isplit
311
312	!
313	!***************************************************************
314	! calcul du flux ascencant conservatif
315	! print*,'<<<<calcul flux ascendant conservatif'
316
317	fmc_therm=0.
318	do k=1,klev
319	do i=1,klon
320	if (entr_therm(i,k).gt.0.) then
321	fmc_therm(i,k+1)=fmc_therm(i,k)+entr_therm(i,k)
322	else
323	fmc_therm(i,k+1)=fmc_therm(i,k)
324	endif
325	detrc_therm(i,k)=(fmc_therm(i,k+1)-fm_therm(i,k+1)) &
326	& -(fmc_therm(i,k)-fm_therm(i,k))
327	enddo
328	enddo
329
330
331	!****************************************************************
332	! calcul de l'humidite dans l'ascendance
333	! print*,'<<<<calcul de lhumidite dans thermique'
334	!CR:on ne le calcule que pour le cas sec
335	if (iflag_thermals.le.11) then
336	do i=1,klon
337	zqasc(i,1)=q_seri(i,1)
338	do k=2,klev
339	if (fmc_therm(i,k+1).gt.1.e-6) then
340	zqasc(i,k)=(fmc_therm(i,k)*zqasc(i,k-1) &
341	& +entr_therm(i,k)*q_seri(i,k))/fmc_therm(i,k+1)
342	!CR:test on asseche le thermique
343	! zqasc(i,k)=zqasc(i,k)/2.
344	! else
345	! zqasc(i,k)=q_seri(i,k)
346	endif
347	enddo
348	enddo
349
350
351	! calcul de l'eau condensee dans l'ascendance
352	! print*,'<<<<calcul de leau condensee dans thermique'
353	do i=1,klon
354	do k=1,klev
355	clwcon0(i,k)=zqasc(i,k)-zqsat(i,k)
356	if (clwcon0(i,k).lt.0. .or. &
357	& (fm_therm(i,k+1)+detrc_therm(i,k)).lt.1.e-6) then
358	clwcon0(i,k)=0.
359	endif
360	enddo
361	enddo
362	else
363	do i=1,klon
364	do k=1,klev
365	clwcon0(i,k)=zqla(i,k)
366	if (clwcon0(i,k).lt.0. .or. &
367	& (fm_therm(i,k+1)+detrc_therm(i,k)).lt.1.e-6) then
368	clwcon0(i,k)=0.
369	endif
370	enddo
371	enddo
372	endif
373	!*******************************************************************
374
375
376	!jyg Protection contre les temperatures nulles
377	do i=1,klon
378	do k=1,klev
379	if (ztla(i,k) .lt. 1.e-10) fraca(i,k) =0.
380	enddo
381	enddo
382
383
384	return
385
386	end

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

Download in other formats:

Original Format