Context Navigation

concvl.F @ 1085

Last change on this file since 1085 was 987, checked in by Laurent Fairhead, 16 years ago
Du nettoyage sur le parallelisme, inclusion de nouvelles interfaces pour OPA9 et ORCHIDEE YM LF
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 14.8 KB

Line
1
2	!
3	! $Header$
4	!
5	SUBROUTINE concvl (iflag_con,iflag_clos,
6	. dtime,paprs,pplay,
7	. t,q,t_wake,q_wake,u,v,tra,ntra,
8	. ALE,ALP,work1,work2,
9	. d_t,d_q,d_u,d_v,d_tra,
10	. rain, snow, kbas, ktop, sigd,
11	. upwd,dnwd,dnwdbis,Ma,mip,Vprecip,
12	. cape,cin,tvp,Tconv,iflag,
13	. pbase,bbase,dtvpdt1,dtvpdq1,dplcldt,dplcldr,
14	. qcondc,wd,pmflxr,pmflxs,
15	. da,phi,mp,dd_t,dd_q,lalim_conv,wght_th)
16	***************************************************************
17	* *
18	* CONCVL *
19	* *
20	* *
21	* written by : Sandrine Bony-Lena, 17/05/2003, 11.16.04 *
22	* modified by : *
23	***************************************************************
24	*
25	c
26	USE dimphy
27	IMPLICIT none
28	c======================================================================
29	c Auteur(s): S. Bony-Lena (LMD/CNRS) date: ???
30	c Objet: schema de convection de Emanuel (1991) interface
31	c======================================================================
32	c Arguments:
33	c dtime--input-R-pas d'integration (s)
34	c s-------input-R-la valeur "s" pour chaque couche
35	c sigs----input-R-la valeur "sigma" de chaque couche
36	c sig-----input-R-la valeur de "sigma" pour chaque niveau
37	c psolpa--input-R-la pression au sol (en Pa)
38	C pskapa--input-R-exponentiel kappa de psolpa
39	c h-------input-R-enthalpie potentielle (CpT/P*kappa)
40	c q-------input-R-vapeur d'eau (en kg/kg)
41	c
42	c work*: input et output: deux variables de travail,
43	c on peut les mettre a 0 au debut
44	c ALE-----input-R-energie disponible pour soulevement
45	c ALP-----input-R-puissance disponible pour soulevement
46	c
47	C d_h-----output-R-increment de l'enthalpie potentielle (h)
48	c d_q-----output-R-increment de la vapeur d'eau
49	c rain----output-R-la pluie (mm/s)
50	c snow----output-R-la neige (mm/s)
51	c upwd----output-R-saturated updraft mass flux (kg/m**2/s)
52	c dnwd----output-R-saturated downdraft mass flux (kg/m**2/s)
53	c dnwd0---output-R-unsaturated downdraft mass flux (kg/m**2/s)
54	c Ma------output-R-adiabatic ascent mass flux (kg/m2/s)
55	c mip-----output-R-mass flux shed by adiabatic ascent (kg/m2/s)
56	c Vprecip-output-R-vertical profile of precipitations (kg/m2/s)
57	c Tconv---output-R-environment temperature seen by convective scheme (K)
58	c Cape----output-R-CAPE (J/kg)
59	c Cin ----output-R-CIN (J/kg)
60	c Tvp-----output-R-Temperature virtuelle d'une parcelle soulevee
61	c adiabatiquement a partir du niveau 1 (K)
62	c deltapb-output-R-distance entre LCL et base de la colonne (<0 ; Pa)
63	c Ice_flag-input-L-TRUE->prise en compte de la thermodynamique de la glace
64	c dd_t-----output-R-increment de la temperature du aux descentes precipitantes
65	c dd_q-----output-R-increment de la vapeur d'eau du aux desc precip
66	c======================================================================
67	c
68	#include "dimensions.h"
69	cccccc#include "dimphy.h"
70	c
71	integer NTRAC
72	PARAMETER (NTRAC=nqmx-2)
73	c
74	INTEGER iflag_con,iflag_clos
75	c
76	REAL dtime, paprs(klon,klev+1),pplay(klon,klev)
77	REAL t(klon,klev),q(klon,klev),u(klon,klev),v(klon,klev)
78	REAL t_wake(klon,klev),q_wake(klon,klev)
79	REAL tra(klon,klev,ntrac)
80	INTEGER ntra
81	REAL work1(klon,klev),work2(klon,klev),ptop2(klon)
82	REAL pmflxr(klon,klev+1),pmflxs(klon,klev+1)
83	REAL ALE(klon),ALP(klon)
84	c
85	REAL d_t(klon,klev),d_q(klon,klev),d_u(klon,klev),d_v(klon,klev)
86	REAL dd_t(klon,klev),dd_q(klon,klev)
87	REAL d_tra(klon,klev,ntrac)
88	REAL rain(klon),snow(klon)
89	c
90	INTEGER kbas(klon),ktop(klon)
91	REAL em_ph(klon,klev+1),em_p(klon,klev)
92	REAL upwd(klon,klev),dnwd(klon,klev),dnwdbis(klon,klev)
93	REAL Ma(klon,klev), mip(klon,klev),Vprecip(klon,klev)
94	real da(klon,klev),phi(klon,klev,klev),mp(klon,klev)
95	REAL cape(klon),cin(klon),tvp(klon,klev)
96	REAL Tconv(klon,klev)
97	c
98	cCR:test: on passe lentr et alim_star des thermiques
99	INTEGER lalim_conv(klon)
100	REAL wght_th(klon,klev)
101	REAL em_sig1feed ! sigma at lower bound of feeding layer
102	REAL em_sig2feed ! sigma at upper bound of feeding layer
103	REAL em_wght(klev) ! weight density determining the feeding mixture
104	con enleve le save
105	c SAVE em_sig1feed,em_sig2feed,em_wght
106	c
107	INTEGER iflag(klon)
108	REAL rflag(klon)
109	REAL pbase(klon),bbase(klon)
110	REAL dtvpdt1(klon,klev),dtvpdq1(klon,klev)
111	REAL dplcldt(klon),dplcldr(klon)
112	REAL qcondc(klon,klev)
113	REAL wd(klon)
114	REAL Plim1(klon),Plim2(klon),asupmax(klon,klev)
115	REAL supmax0(klon),asupmaxmin(klon)
116	c
117	REAL sigd(klon)
118	REAL zx_t,zdelta,zx_qs,zcor
119	c
120	! INTEGER iflag_mix
121	! SAVE iflag_mix
122	INTEGER noff, minorig
123	INTEGER i,k,itra
124	REAL qs(klon,klev),qs_wake(klon,klev)
125	cLF REAL cbmf(klon)
126	cLF SAVE cbmf
127	REAL, SAVE, ALLOCATABLE :: cbmf(:)
128	c$OMP THREADPRIVATE(cbmf)!
129	REAL cbmflast(klon)
130	INTEGER ifrst
131	SAVE ifrst
132	DATA ifrst /0/
133	c$OMP THREADPRIVATE(ifrst)
134
135	c
136	C Variables supplementaires liees au bilan d'energie
137	c Real paire(klon)
138	cLF Real ql(klon,klev)
139	c Save paire
140	cLF Save ql
141	cLF Real t1(klon,klev),q1(klon,klev)
142	cLF Save t1,q1
143	c Data paire /1./
144	REAL, SAVE, ALLOCATABLE :: ql(:,:), q1(:,:), t1(:,:)
145	c$OMP THREADPRIVATE(ql, q1, t1)
146	c
147	C Variables liees au bilan d'energie et d'enthalpi
148	REAL ztsol(klon)
149	REAL h_vcol_tot, h_dair_tot, h_qw_tot, h_ql_tot
150	$ , h_qs_tot, qw_tot, ql_tot, qs_tot , ec_tot
151	SAVE h_vcol_tot, h_dair_tot, h_qw_tot, h_ql_tot
152	$ , h_qs_tot, qw_tot, ql_tot, qs_tot , ec_tot
153	c$OMP THREADPRIVATE(h_vcol_tot, h_dair_tot, h_qw_tot, h_ql_tot)
154	c$OMP THREADPRIVATE(h_qs_tot, qw_tot, ql_tot, qs_tot , ec_tot)
155	REAL d_h_vcol, d_h_dair, d_qt, d_qw, d_ql, d_qs, d_ec
156	REAL d_h_vcol_phy
157	REAL fs_bound, fq_bound
158	SAVE d_h_vcol_phy
159	c$OMP THREADPRIVATE(d_h_vcol_phy)
160	REAL zero_v(klon)
161	CHARACTER*15 ztit
162	INTEGER ip_ebil ! PRINT level for energy conserv. diag.
163	SAVE ip_ebil
164	DATA ip_ebil/2/
165	c$OMP THREADPRIVATE(ip_ebil)
166	INTEGER if_ebil ! level for energy conserv. dignostics
167	SAVE if_ebil
168	DATA if_ebil/2/
169	c$OMP THREADPRIVATE(if_ebil)
170	c+jld ec_conser
171	REAL d_t_ec(klon,klev) ! tendance du a la conersion Ec -> E thermique
172	REAL ZRCPD
173	c-jld ec_conser
174	cLF
175	INTEGER nloc
176	logical, save :: first=.true.
177	c$OMP THREADPRIVATE(first)
178	INTEGER, SAVE :: itap, igout
179	c$OMP THREADPRIVATE(itap, igout)
180	c
181	#include "YOMCST.h"
182	#include "YOMCST2.h"
183	#include "YOETHF.h"
184	#include "FCTTRE.h"
185	#include "iniprint.h"
186	c
187	if (first) then
188	c Allocate some variables LF 04/2008
189	c
190	allocate(cbmf(klon))
191	allocate(ql(klon,klev))
192	allocate(t1(klon,klev))
193	allocate(q1(klon,klev))
194	itap=0
195	igout=klon/2+1/klon
196	endif
197	c Incrementer le compteur de la physique
198	itap = itap + 1
199
200	c Copy T into Tconv
201	DO k = 1,klev
202	DO i = 1,klon
203	Tconv(i,k) = T(i,k)
204	ENDDO
205	ENDDO
206	c
207	IF (if_ebil.ge.1) THEN
208	DO i=1,klon
209	ztsol(i) = t(i,1)
210	zero_v(i)=0.
211	Do k = 1,klev
212	ql(i,k) = 0.
213	ENDDO
214	END DO
215	END IF
216	c
217	cym
218	snow(:)=0
219
220	c IF (ifrst .EQ. 0) THEN
221	c ifrst = 1
222	if (first) then
223	first=.false.
224	c
225	C===========================================================================
226	C READ IN PARAMETERS FOR THE CLOSURE AND THE MIXING DISTRIBUTION
227	C===========================================================================
228	C
229	if (iflag_con.eq.3) then
230	c CALL cv3_inicp()
231	CALL cv3_inip()
232	endif
233	c
234	C===========================================================================
235	C READ IN PARAMETERS FOR CONVECTIVE INHIBITION BY TROPOS. DRYNESS
236	C===========================================================================
237	C
238	cc$$$ open (56,file='supcrit.data')
239	cc$$$ read (56,*) Supcrit1, Supcrit2
240	cc$$$ close (56)
241	c
242	print*, 'supcrit1, supcrit2' ,supcrit1, supcrit2
243	C
244	C===========================================================================
245	C Initialisation pour les bilans d'eau et d'energie
246	C===========================================================================
247	IF (if_ebil.ge.1) d_h_vcol_phy=0.
248	c
249	DO i = 1, klon
250	cbmf(i) = 0.
251	sigd(i) = 0.
252	ENDDO
253	ENDIF !(ifrst .EQ. 0)
254
255	DO k = 1, klev+1
256	DO i=1,klon
257	em_ph(i,k) = paprs(i,k) / 100.0
258	pmflxs(i,k)=0.
259	ENDDO
260	ENDDO
261	c
262	DO k = 1, klev
263	DO i=1,klon
264	em_p(i,k) = pplay(i,k) / 100.0
265	ENDDO
266	ENDDO
267	c
268	!
269	! Feeding layer
270	!
271	em_sig1feed = 1.
272	em_sig2feed = 0.97
273	c em_sig2feed = 0.8
274	! Relative Weight densities
275	do k=1,klev
276	em_wght(k)=1.
277	end do
278	cCRtest: couche alim des tehrmiques ponderee par a*
279	c DO i = 1, klon
280	c do k=1,lalim_conv(i)
281	c em_wght(k)=wght_th(i,k)
282	c print*,'em_wght=',em_wght(k),wght_th(i,k)
283	c end do
284	c END DO
285
286	if (iflag_con .eq. 4) then
287	DO k = 1, klev
288	DO i = 1, klon
289	zx_t = t(i,k)
290	zdelta=MAX(0.,SIGN(1.,rtt-zx_t))
291	zx_qs= MIN(0.5 , r2es * FOEEW(zx_t,zdelta)/em_p(i,k)/100.0)
292	zcor=1./(1.-retv*zx_qs)
293	qs(i,k)=zx_qs*zcor
294	ENDDO
295	DO i = 1, klon
296	zx_t = t_wake(i,k)
297	zdelta=MAX(0.,SIGN(1.,rtt-zx_t))
298	zx_qs= MIN(0.5 , r2es * FOEEW(zx_t,zdelta)/em_p(i,k)/100.0)
299	zcor=1./(1.-retv*zx_qs)
300	qs_wake(i,k)=zx_qs*zcor
301	ENDDO
302	ENDDO
303	else ! iflag_con=3 (modif de puristes qui fait la diffce pour la convergence numerique)
304	DO k = 1, klev
305	DO i = 1, klon
306	zx_t = t(i,k)
307	zdelta=MAX(0.,SIGN(1.,rtt-zx_t))
308	zx_qs= r2es * FOEEW(zx_t,zdelta)/em_p(i,k)/100.0
309	zx_qs= MIN(0.5,zx_qs)
310	zcor=1./(1.-retv*zx_qs)
311	zx_qs=zx_qs*zcor
312	qs(i,k)=zx_qs
313	ENDDO
314	DO i = 1, klon
315	zx_t = t_wake(i,k)
316	zdelta=MAX(0.,SIGN(1.,rtt-zx_t))
317	zx_qs= r2es * FOEEW(zx_t,zdelta)/em_p(i,k)/100.0
318	zx_qs= MIN(0.5,zx_qs)
319	zcor=1./(1.-retv*zx_qs)
320	zx_qs=zx_qs*zcor
321	qs_wake(i,k)=zx_qs
322	ENDDO
323	ENDDO
324	endif ! iflag_con
325	c
326	C------------------------------------------------------------------
327
328	C Main driver for convection:
329	C iflag_con=3 -> nvlle version de KE (JYG)
330	C iflag_con = 30 -> equivalent to convect3
331	C iflag_con = 4 -> equivalent to convect1/2
332	c
333	c
334	if (iflag_con.eq.30) then
335
336	CALL cv_driver(klon,klev,klev+1,ntra,iflag_con,
337	: t,q,qs,u,v,tra,
338	$ em_p,em_ph,iflag,
339	$ d_t,d_q,d_u,d_v,d_tra,rain,
340	$ pmflxr,cbmf,work1,work2,
341	$ kbas,ktop,
342	$ dtime,Ma,upwd,dnwd,dnwdbis,qcondc,wd,cape,
343	$ da,phi,mp)
344
345	else
346
347	cLF necessary for gathered fields
348	nloc=klon
349	CALL cva_driver(klon,klev,klev+1,ntra,nloc,
350	$ iflag_con,iflag_mix,iflag_clos,dtime,
351	: t,q,qs,t_wake,q_wake,qs_wake,u,v,tra,
352	$ em_p,em_ph,
353	. ALE,ALP,
354	. em_sig1feed,em_sig2feed,em_wght,
355	. iflag,d_t,d_q,d_u,d_v,d_tra,rain,kbas,ktop,
356	$ cbmf,work1,work2,ptop2,sigd,
357	$ Ma,mip,Vprecip,upwd,dnwd,dnwdbis,qcondc,wd,
358	$ cape,cin,tvp,
359	$ dd_t,dd_q,Plim1,Plim2,asupmax,supmax0,
360	$ asupmaxmin,lalim_conv)
361	endif
362	C------------------------------------------------------------------
363
364	DO i = 1,klon
365	rain(i) = rain(i)/86400.
366	rflag(i)=iflag(i)
367	ENDDO
368
369	DO k = 1, klev
370	DO i = 1, klon
371	d_t(i,k) = dtime*d_t(i,k)
372	d_q(i,k) = dtime*d_q(i,k)
373	d_u(i,k) = dtime*d_u(i,k)
374	d_v(i,k) = dtime*d_v(i,k)
375	ENDDO
376	ENDDO
377	c
378	if (iflag_con.eq.30) then
379	DO itra = 1,ntra
380	DO k = 1, klev
381	DO i = 1, klon
382	d_tra(i,k,itra) =dtime*d_tra(i,k,itra)
383	ENDDO
384	ENDDO
385	ENDDO
386	endif
387
388	DO k = 1, klev
389	DO i = 1, klon
390	t1(i,k) = t(i,k)+ d_t(i,k)
391	q1(i,k) = q(i,k)+ d_q(i,k)
392	ENDDO
393	ENDDO
394	c
395	cc IF (if_ebil.ge.2) THEN
396	cc ztit='after convect'
397	cc CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
398	cc e , t1,q1,ql,qs,u,v,paprs,pplay
399	cc s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
400	cc call diagphy(paire,ztit,ip_ebil
401	cc e , zero_v, zero_v, zero_v, zero_v, zero_v
402	cc e , zero_v, rain, zero_v, ztsol
403	cc e , d_h_vcol, d_qt, d_ec
404	cc s , fs_bound, fq_bound )
405	cc END IF
406	C
407	c
408	c les traceurs ne sont pas mis dans cette version de convect4:
409	if (iflag_con.eq.4) then
410	DO itra = 1,ntra
411	DO k = 1, klev
412	DO i = 1, klon
413	d_tra(i,k,itra) = 0.
414	ENDDO
415	ENDDO
416	ENDDO
417	endif
418	c print*, 'concvl->: dd_t,dd_q ',dd_t(1,1),dd_q(1,1)
419
420	DO k = 1, klev
421	DO i = 1, klon
422	dtvpdt1(i,k) = 0.
423	dtvpdq1(i,k) = 0.
424	ENDDO
425	ENDDO
426	DO i = 1, klon
427	dplcldt(i) = 0.
428	dplcldr(i) = 0.
429	ENDDO
430	c
431	if(prt_level.GE.20) THEN
432	DO k=1,klev
433	! print*,'physiq apres_add_con i k it d_u d_v d_t d_q qdl0',igout
434	! .,k,itap,d_u_con(igout,k) ,d_v_con(igout,k), d_t_con(igout,k),
435	! .d_q_con(igout,k),dql0(igout,k)
436	! print*,'phys apres_add_con itap Ma cin ALE ALP wak t q undi t q'
437	! .,itap,Ma(igout,k),cin(igout),ALE(igout), ALP(igout),
438	! . t_wake(igout,k),q_wake(igout,k),t_undi(igout,k),q_undi(igout,k)
439	! print*,'phy apres_add_con itap CON rain snow EMA wk1 wk2 Vpp mip'
440	! .,itap,rain_con(igout),snow_con(igout),ema_work1(igout,k),
441	! .ema_work2(igout,k),Vprecip(igout,k), mip(igout,k)
442	! print*,'phy apres_add_con itap upwd dnwd dnwd0 cape tvp Tconv '
443	! .,itap,upwd(igout,k),dnwd(igout,k),dnwd0(igout,k),cape(igout),
444	! .tvp(igout,k),Tconv(igout,k)
445	! print*,'phy apres_add_con itap dtvpdt dtvdq dplcl dplcldr qcondc'
446	! .,itap,dtvpdt1(igout,k),dtvpdq1(igout,k),dplcldt(igout),
447	! .dplcldr(igout),qcondc(igout,k)
448	! print*,'phy apres_add_con itap wd pmflxr Kpmflxr Kp1 Kpmflxs Kp1'
449	! .,itap,wd(igout),pmflxr(igout,k),pmflxr(igout,k+1),pmflxs(igout,k)
450	! .,pmflxs(igout,k+1)
451	! print*,'phy apres_add_con itap da phi mp ftd fqd lalim wgth',
452	! .itap,da(igout,k),phi(igout,k,k),mp(igout,k),ftd(igout,k),
453	! . fqd(igout,k),lalim_conv(igout),wght_th(igout,k)
454	ENDDO
455	endif !(prt_level.EQ.20) THEN
456	c
457	RETURN
458	END
459

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

Download in other formats:

Original Format