Context Navigation

calltherm_interface.F90 @ 355

Last change on this file since 355 was 342, checked in by acolaitis, 14 years ago

M 341 libf/phymars/calltherm_interface.F90
D 341 libf/phymars/calltherm_mars.F90
^{---------------- Merged calltherm_mars with calltherm_interface for simplicity.}

Cleaned up variables and code

M 341 libf/phymars/thermcell_dqup.F90
^{---------------- dqup now output derivatives and not tendancies}

File size: 14.2 KB

Line
1	!
2	! AC 2011-01-05
3	!
4	SUBROUTINE calltherm_interface (firstcall, &
5	& long,lati,zzlev,zzlay, &
6	& ptimestep,pu,pv,pt,pq,pdu,pdv,pdt,pdq,q2, &
7	& pplay,pplev,pphi,zpopsk, &
8	& pdu_th,pdv_th,pdt_th,pdq_th,lmax,zmaxth,pbl_dtke,hfmax,wmax)
9
10	USE ioipsl_getincom
11
12	implicit none
13	#include "callkeys.h"
14	#include "dimensions.h"
15	#include "dimphys.h"
16	#include "comcstfi.h"
17
18	!--------------------------------------------------------
19	! Input Variables
20	!--------------------------------------------------------
21
22	REAL, INTENT(IN) :: ptimestep
23	REAL, INTENT(IN) :: pplev(ngridmx,nlayermx+1),pplay(ngridmx,nlayermx)
24	REAL, INTENT(IN) :: pphi(ngridmx,nlayermx)
25	REAL, INTENT(IN) :: pu(ngridmx,nlayermx),pv(ngridmx,nlayermx)
26	REAL, INTENT(IN) :: pt(ngridmx,nlayermx),pq(ngridmx,nlayermx,nqmx)
27	REAL, INTENT(IN) :: zzlay(ngridmx,nlayermx)
28	REAL, INTENT(IN) :: zzlev(ngridmx,nlayermx+1)
29	LOGICAL, INTENT(IN) :: firstcall
30	REAL, INTENT(IN) :: pdu(ngridmx,nlayermx),pdv(ngridmx,nlayermx)
31	REAL, INTENT(IN) :: pdq(ngridmx,nlayermx,nqmx),pdt(ngridmx,nlayermx)
32	REAL, INTENT(IN) :: q2(ngridmx,nlayermx+1)
33	REAL, INTENT(IN) :: long(ngridmx),lati(ngridmx)
34	REAL, INTENT(IN) :: zpopsk(ngridmx,nlayermx)
35
36	!--------------------------------------------------------
37	! Output Variables
38	!--------------------------------------------------------
39
40	REAL, INTENT(OUT) :: pdu_th(ngridmx,nlayermx)
41	REAL, INTENT(OUT) :: pdv_th(ngridmx,nlayermx)
42	REAL, INTENT(OUT) :: pdt_th(ngridmx,nlayermx)
43	REAL, INTENT(OUT) :: pdq_th(ngridmx,nlayermx,nqmx)
44	INTEGER, INTENT(OUT) :: lmax(ngridmx)
45	REAL, INTENT(OUT) :: zmaxth(ngridmx)
46	REAL, INTENT(OUT) :: pbl_dtke(ngridmx,nlayermx+1)
47
48	!--------------------------------------------------------
49	! Thermals local variables
50	!--------------------------------------------------------
51	REAL zu(ngridmx,nlayermx), zv(ngridmx,nlayermx)
52	REAL zt(ngridmx,nlayermx)
53	REAL d_t_ajs(ngridmx,nlayermx)
54	REAL d_u_ajs(ngridmx,nlayermx), d_q_ajs(ngridmx,nlayermx,nqmx)
55	REAL d_v_ajs(ngridmx,nlayermx)
56	REAL fm_therm(ngridmx,nlayermx+1), entr_therm(ngridmx,nlayermx)
57	REAL detr_therm(ngridmx,nlayermx)
58	REAL zw2(ngridmx,nlayermx+1)
59	REAL fraca(ngridmx,nlayermx+1)
60	REAL ztla(ngridmx,nlayermx)
61	REAL q_therm(ngridmx,nlayermx), pq_therm(ngridmx,nlayermx,nqmx)
62	REAL dq_therm(ngridmx,nlayermx), dq_thermdown(ngridmx,nlayermx)
63	REAL q2_therm(ngridmx,nlayermx), dq2_therm(ngridmx,nlayermx)
64	REAL lmax_real(ngridmx)
65	REAL masse(ngridmx,nlayermx)
66	REAL zdz(ngridmx,nlayermx)
67	LOGICAL qtransport_thermals,dtke_thermals
68	INTEGER l,ig,iq
69	CHARACTER (LEN=20) :: modname
70
71	!--------------------------------------------------------
72	! Local variables for sub-timestep
73	!--------------------------------------------------------
74
75	REAL d_t_the(ngridmx,nlayermx), d_q_the(ngridmx,nlayermx,nqmx)
76	REAL d_u_the(ngridmx,nlayermx),d_v_the(ngridmx,nlayermx)
77	REAL dq2_the(ngridmx,nlayermx)
78	INTEGER isplit,nsplit_thermals
79	REAL r_aspect_thermals
80	REAL fact
81	REAL zfm_therm(ngridmx,nlayermx+1),zdt
82	REAL zentr_therm(ngridmx,nlayermx),zdetr_therm(ngridmx,nlayermx)
83	REAL zheatFlux(ngridmx,nlayermx)
84	REAL zheatFlux_down(ngridmx,nlayermx)
85	REAL zbuoyancyOut(ngridmx,nlayermx)
86	REAL zbuoyancyEst(ngridmx,nlayermx)
87	REAL zzw2(ngridmx,nlayermx+1)
88	REAL zmax(ngridmx)
89
90	!--------------------------------------------------------
91	! Diagnostics
92	!--------------------------------------------------------
93
94	REAL heatFlux(ngridmx,nlayermx)
95	REAL heatFlux_down(ngridmx,nlayermx)
96	REAL buoyancyOut(ngridmx,nlayermx)
97	REAL buoyancyEst(ngridmx,nlayermx)
98	REAL hfmax(ngridmx),wmax(ngridmx)
99
100	!---------------------------------------------------------
101
102
103	! **********************************************************************
104	! Initialization
105	! **********************************************************************
106
107	lmax(:)=0.
108	pdu_th(:,:)=0.
109	pdv_th(:,:)=0.
110	pdt_th(:,:)=0.
111	entr_therm(:,:)=0.
112	detr_therm(:,:)=0.
113	q2_therm(:,:)=0.
114	dq2_therm(:,:)=0.
115	dq_therm(:,:)=0.
116	dq_thermdown(:,:)=0.
117	ztla(:,:)=0.
118	pbl_dtke(:,:)=0.
119	fm_therm(:,:)=0.
120	zw2(:,:)=0.
121	fraca(:,:)=0.
122	if (tracer) then
123	pdq_th(:,:,:)=0.
124	end if
125	d_t_ajs(:,:)=0.
126	d_u_ajs(:,:)=0.
127	d_v_ajs(:,:)=0.
128	d_q_ajs(:,:,:)=0.
129	heatFlux(:,:)=0.
130	heatFlux_down(:,:)=0.
131	buoyancyOut(:,:)=0.
132	buoyancyEst(:,:)=0.
133	zmaxth(:)=0.
134	lmax_real(:)=0.
135
136
137	! **********************************************************************
138	! Preparing inputs for the thermals
139	! **********************************************************************
140
141	zu(:,:)=pu(:,:)+pdu(:,:)*ptimestep
142	zv(:,:)=pv(:,:)+pdv(:,:)*ptimestep
143	zt(:,:)=pt(:,:)+pdt(:,:)*ptimestep
144
145	pq_therm(:,:,:)=0.
146	qtransport_thermals=.true. !! default setting
147	!call getin("qtransport_thermals",qtransport_thermals)
148
149	if(qtransport_thermals) then
150	if(tracer) then
151	pq_therm(:,:,:)=pq(:,:,:)+pdq(:,:,:)*ptimestep
152	endif
153	endif
154
155	! dtke_thermals=.false. !! default setting
156	! !call getin("dtke_thermals",dtke_thermals)
157	!
158	! IF(dtke_thermals) THEN
159	! DO l=1,nlayermx
160	! q2_therm(:,l)=0.5*(q2(:,l)+q2(:,l+1))
161	! ENDDO
162	! ENDIF
163
164	! **********************************************************************
165	! **********************************************************************
166	! **********************************************************************
167	! CALLTHERM
168	! **********************************************************************
169	! **********************************************************************
170	! **********************************************************************
171
172	! r_aspect_thermals ! ultimately conrols the amount of mass going through the thermals
173	! decreasing it increases the thermals effect. Tests at gcm resolution
174	! shows that too low values destabilize the model
175	! when changing this value, one should check that the surface layer model
176	! outputs the correct Cdu and Chu through changing the gustiness coefficient beta
177
178	#ifdef MESOSCALE
179	!! valid for timesteps < 200s
180	nsplit_thermals=2
181	r_aspect_thermals=0.7
182	#else
183	nsplit_thermals=35
184	r_aspect_thermals=1.5
185	#endif
186	call getin("nsplit_thermals",nsplit_thermals)
187	call getin("r_aspect_thermals",r_aspect_thermals)
188
189	! **********************************************************************
190	! SUB-TIMESTEP LOOP
191	! **********************************************************************
192
193	zdt=ptimestep/REAL(nsplit_thermals)
194
195	DO isplit=1,nsplit_thermals
196
197	! Initialization of intermediary variables
198
199	zfm_therm(:,:)=0.
200	zentr_therm(:,:)=0.
201	zdetr_therm(:,:)=0.
202	zheatFlux(:,:)=0.
203	zheatFlux_down(:,:)=0.
204	! zbuoyancyOut(:,:)=0.
205	! zbuoyancyEst(:,:)=0.
206	zzw2(:,:)=0.
207	zmax(:)=0.
208	lmax(:)=0.
209	d_t_the(:,:)=0.
210	d_u_the(:,:)=0.
211	d_v_the(:,:)=0.
212	dq2_the(:,:)=0.
213	if (nqmx .ne. 0) then
214	d_q_the(:,:,:)=0.
215	endif
216
217	CALL thermcell_main_mars(zdt &
218	& ,pplay,pplev,pphi,zzlev,zzlay &
219	& ,zu,zv,zt,pq_therm,q2_therm &
220	& ,d_u_the,d_v_the,d_t_the,d_q_the,dq2_the &
221	& ,zfm_therm,zentr_therm,zdetr_therm,lmax,zmax &
222	& ,r_aspect_thermals &
223	& ,zzw2,fraca,zpopsk &
224	& ,ztla,zheatFlux,zheatFlux_down &
225	& ,zbuoyancyOut,zbuoyancyEst)
226
227	fact=1./REAL(nsplit_thermals)
228
229	d_t_the(:,:)=d_t_the(:,:)ptimestepfact
230	! d_u_the(:,:)=d_u_the(:,:)*fact
231	! d_v_the(:,:)=d_v_the(:,:)*fact
232	! dq2_the(:,:)=dq2_the(:,:)*fact
233	! if (nqmx .ne. 0) then
234	! d_q_the(:,:,:)=d_q_the(:,:,:)*fact
235	! endif
236
237	zmaxth(:)=zmaxth(:)+zmax(:)*fact
238	lmax_real(:)=lmax_real(:)+float(lmax(:))*fact
239	fm_therm(:,:)=fm_therm(:,:) &
240	& +zfm_therm(:,:)*fact
241	entr_therm(:,:)=entr_therm(:,:) &
242	& +zentr_therm(:,:)*fact
243	detr_therm(:,:)=detr_therm(:,:) &
244	& +zdetr_therm(:,:)*fact
245
246	heatFlux(:,:)=heatFlux(:,:) &
247	& +zheatFlux(:,:)*fact
248	heatFlux_down(:,:)=heatFlux_down(:,:) &
249	& +zheatFlux_down(:,:)*fact
250	! buoyancyOut(:,:)=buoyancyOut(:,:) &
251	! & +zbuoyancyOut(:,:)*fact
252	! buoyancyEst(:,:)=buoyancyEst(:,:) &
253	! & +zbuoyancyEst(:,:)*fact
254
255	zw2(:,:)=zw2(:,:) + zzw2(:,:)*fact
256
257	! accumulation de la tendance
258
259	d_t_ajs(:,:)=d_t_ajs(:,:)+d_t_the(:,:)
260	! d_u_ajs(:,:)=d_u_ajs(:,:)+d_u_the(:,:)
261	! d_v_ajs(:,:)=d_v_ajs(:,:)+d_v_the(:,:)
262	! d_q_ajs(:,:,:)=d_q_ajs(:,:,:)+d_q_the(:,:,:)
263	! dq2_therm(:,:)=dq2_therm(:,:)+dq2_the(:,:)
264	! incrementation des variables meteo
265
266	zt(:,:) = zt(:,:) + d_t_the(:,:)
267	! zu(:,:) = zu(:,:) + d_u_the(:,:)
268	! zv(:,:) = zv(:,:) + d_v_the(:,:)
269	! pq_therm(:,:,:) = pq_therm(:,:,:) + d_q_the(:,:,:)
270	! q2_therm(:,:) = q2_therm(:,:) + dq2_therm(:,:)
271
272
273	ENDDO ! isplit
274	!****************************************************************
275
276	! Now that we have computed total entrainment and detrainment, we can
277	! advect u, v, and q in thermals. (theta already advected). We can do
278	! that separatly because u,v,and q are not used in thermcell_main for
279	! any thermals-related computation : they are purely passive.
280
281	! mass of cells
282	do l=1,nlayermx
283	masse(:,l)=(pplev(:,l)-pplev(:,l+1))/g
284	enddo
285
286	! thickness of layers
287	do l=1,nlayermx
288	zdz(:,l)=zzlev(:,l+1)-zzlev(:,l)
289	enddo
290
291	modname='momentum'
292	call thermcell_dqup(ngridmx,nlayermx,ptimestep &
293	& ,fm_therm,entr_therm,detr_therm, &
294	& masse,zu,d_u_ajs,modname,zdz)
295
296	call thermcell_dqup(ngridmx,nlayermx,ptimestep &
297	& ,fm_therm,entr_therm,detr_therm, &
298	& masse,zv,d_v_ajs,modname,zdz)
299
300	if (nqmx .ne. 0.) then
301	modname='tracer'
302	DO iq=1,nqmx
303	call thermcell_dqup(ngridmx,nlayermx,ptimestep &
304	& ,fm_therm,entr_therm,detr_therm, &
305	& masse,pq_therm(:,:,iq),d_q_ajs(:,:,iq),modname,zdz)
306
307	ENDDO
308	endif
309
310	DO ig=1,ngridmx
311	hfmax(ig)=MAXVAL(heatFlux(ig,:)+heatFlux_down(ig,:))
312	wmax(ig)=MAXVAL(zw2(ig,:))
313	ENDDO
314
315	lmax(:)=nint(lmax_real(:))
316
317	! **********************************************************************
318	! **********************************************************************
319	! **********************************************************************
320	! CALLTHERM END
321	! **********************************************************************
322	! **********************************************************************
323	! **********************************************************************
324
325
326	! **********************************************************************
327	! Preparing outputs
328	! **********************************************************************
329
330	! Winds and tracers PDU, PDV, and PDQ are DERIVATIVES
331
332	pdu_th(:,:)=d_u_ajs(:,:)
333	pdv_th(:,:)=d_v_ajs(:,:)
334
335	if(qtransport_thermals) then
336	if(tracer) then
337	pdq_th(:,:,:)=d_q_ajs(:,:,:)
338	endif
339	endif
340
341	! IF(dtke_thermals) THEN
342	! DO l=2,nlayermx
343	! pbl_dtke(:,l)=0.5*(dq2_therm(:,l-1)+dq2_therm(:,l))
344	! ENDDO
345	!
346	! pbl_dtke(:,1)=0.5*dq2_therm(:,1)
347	! pbl_dtke(:,nlayermx+1)=0.
348	! ENDIF
349
350
351	! Temperature PDT is a TENDANCY
352	pdt_th(:,:)=d_t_ajs(:,:)/ptimestep
353
354	! **********************************************************************
355	! Diagnostics
356	! **********************************************************************
357
358	if(outptherm) then
359	if (ngridmx .eq. 1) then
360	call WRITEDIAGFI(ngridmx,'entr_therm','entrainement thermique',&
361	& 'kg/m-2',1,entr_therm)
362	call WRITEDIAGFI(ngridmx,'detr_therm','detrainement thermique',&
363	& 'kg/m-2',1,detr_therm)
364	call WRITEDIAGFI(ngridmx,'fm_therm','flux masse thermique',&
365	& 'kg/m-2',1,fm_therm)
366	call WRITEDIAGFI(ngridmx,'zw2','vitesse verticale thermique',&
367	& 'm/s',1,zw2)
368	call WRITEDIAGFI(ngridmx,'heatFlux_up','heatFlux_updraft',&
369	& 'SI',1,heatFlux)
370	call WRITEDIAGFI(ngridmx,'heatFlux_down','heatFlux_downdraft',&
371	& 'SI',1,heatFlux_down)
372	call WRITEDIAGFI(ngridmx,'fraca','fraction coverage',&
373	& 'percent',1,fraca)
374	call WRITEDIAGFI(ngridmx,'buoyancyOut','buoyancyOut',&
375	& 'm.s-2',1,buoyancyOut)
376	call WRITEDIAGFI(ngridmx,'buoyancyEst','buoyancyEst',&
377	& 'm.s-2',1,buoyancyEst)
378	call WRITEDIAGFI(ngridmx,'d_t_th', &
379	& 'tendance temp TH','K',1,d_t_ajs)
380	call WRITEDIAGFI(ngridmx,'zmax', &
381	& 'pbl height','m',0,zmaxth)
382	else
383
384	call WRITEDIAGFI(ngridmx,'entr_therm','entrainement thermique',&
385	& 'kg/m-2',3,entr_therm)
386	call WRITEDIAGFI(ngridmx,'detr_therm','detrainement thermique',&
387	& 'kg/m-2',3,detr_therm)
388	call WRITEDIAGFI(ngridmx,'fm_therm','flux masse thermique',&
389	& 'kg/m-2',3,fm_therm)
390	call WRITEDIAGFI(ngridmx,'zw2','vitesse verticale thermique',&
391	& 'm/s',3,zw2)
392	call WRITEDIAGFI(ngridmx,'heatFlux','heatFlux',&
393	& 'SI',3,heatFlux)
394	call WRITEDIAGFI(ngridmx,'buoyancyOut','buoyancyOut',&
395	& 'SI',3,buoyancyOut)
396	call WRITEDIAGFI(ngridmx,'d_t_th', &
397	& 'tendance temp TH','K',3,d_t_ajs)
398
399	endif
400	endif
401
402	END

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

Download in other formats:

Original Format