Context Navigation

thermcell_dqupdown.F90 @ 190

Last change on this file since 190 was 185, checked in by acolaitis, 13 years ago

17/06/2011 == AC

Added new settings for the Martian thermals from new LES observations
Revamped thermcell's module variables to allow it's removal
Minor changes in physiq and meso_physiq for the call to thermals
Switched from dynamic to static memory allocation for all thermals variable

to gain computation speed

Thermals now output maximum speed, maximum heat flux, and maximum height in thermal plumes

Property svn:executable set to *

File size: 9.8 KB

Line
1	subroutine thermcell_dqupdown(ngrid,nlayer,ptimestep,fm0,entr0, &
2	& detr0,masse0,q_therm,dq_therm,ztvd,fm_down,ztv,charvar,lmax)
3	implicit none
4
5	! #include "iniprint.h"
6	!=======================================================================
7	!
8	! Calcul du transport verticale dans la couche limite en presence
9	! de "thermiques" explicitement representes
10	! calcul du dq/dt une fois qu'on connait les ascendances
11	! Version modifiee pour prendre les downdrafts a la place de la
12	! subsidence compensatoire
13	!=======================================================================
14
15	#include "dimensions.h"
16	#include "dimphys.h"
17
18	! ============================ INPUTS ============================
19
20	INTEGER, INTENT(IN) :: ngrid,nlayer
21	REAL, INTENT(IN) :: ptimestep
22	REAL, INTENT(IN) :: fm0(ngridmx,nlayermx+1)
23	REAL, INTENT(IN) :: entr0(ngridmx,nlayermx),detr0(ngridmx,nlayermx)
24	REAL, INTENT(IN) :: q_therm(ngridmx,nlayermx)
25	REAL, INTENT(IN) :: fm_down(ngridmx,nlayermx+1)
26	REAL, INTENT(IN) :: ztvd(ngridmx,nlayermx)
27	REAL, INTENT(IN) :: ztv(ngridmx,nlayermx)
28	CHARACTER (LEN=20), INTENT(IN) :: charvar
29	REAL, INTENT(IN) :: masse0(ngridmx,nlayermx)
30	INTEGER, INTENT(IN) :: lmax(ngridmx)
31
32	! ============================ OUTPUTS ===========================
33
34	REAL, INTENT(OUT) :: dq_therm(ngridmx,nlayermx)
35
36	! ============================ LOCAL =============================
37
38	! REAL detr0(ngridmx,nlayermx)
39	REAL detrd(ngridmx,nlayermx)
40	REAL entrd(ngridmx,nlayermx)
41	REAL fmd(ngridmx,nlayermx+1)
42	REAL q(ngridmx,nlayermx)
43	REAL qa(ngridmx,nlayermx)
44	REAL qd(ngridmx,nlayermx)
45	INTEGER ig,k
46	LOGICAL active(ngridmx,nlayermx)
47	INTEGER lmax_down(ngridmx),lmin_down(ngridmx)
48	INTEGER ncorec
49
50	! =========== Init ==============================================
51
52	entrd(:,:)=0.
53	detrd(:,:)=0.
54	qa(:,:)=q_therm(:,:)
55	q(:,:)=q_therm(:,:)
56	qd(:,:)=q_therm(:,:)
57	active(:,:)=.false.
58
59	! previous calculation of zdthl_down uses the divergence of fmd
60	! so it can be negative without problem. Here we include the sign
61	! of fmd in the equations, so it has to be positive
62	!
63	! fmd(:,:)=-fm_down(:,:)
64	!
65	!! ========== Entrainment, Detrainement and Mass =================
66	!
67	!! ========== DOWNDRAFT TRANSPORT DISABLED FOR NOW ===============
68	!
69	! do ig=1,ngridmx
70	! if (ztv(ig,nlayermx)-ztvd(ig,nlayermx) .gt. 0.5) then
71	! print*,"downdraft non nul derniere couche !!! (dqupdown)"
72	! endif
73	! detrd(ig,nlayermx)=0.
74	! entrd(ig,nlayermx)=0.
75	! enddo
76	!
77	! do k=nlayermx-1,1,-1
78	! do ig=1,ngridmx
79	!
80	! if (ztv(ig,k)-ztvd(ig,k) .gt. 0.0001) then
81	! detrd(ig,k)=MAX(0.,(fmd(ig,k+1)*(ztv(ig,k)-ztvd(ig,k+1))) &
82	! & /(ztv(ig,k)-ztvd(ig,k)) - fmd(ig,k))
83	! entrd(ig,k)=MAX(0.,(fmd(ig,k+1)*(ztvd(ig,k)-ztvd(ig,k+1))) &
84	! & /(ztv(ig,k)-ztvd(ig,k)))
85	!
86	! endif
87	! enddo
88	! enddo
89	!
90	!! ======= We have computed entrainment and detrainment from a prescribed
91	!! mass flux and potential temp profile. Due to the way downdraft are parametrized,
92	!! this can yield negative entr and detr. We force it to be positive, but in
93	!! order to conserve tracers, we need to recompute an adequate mass flux
94	!! and modify interface rates, to preserve consistency.
95	! lmax_down(:)=1
96	! lmin_down(:)=1
97	!
98	! do k=1,nlayermx
99	! do ig=1,ngridmx
100	! if ((entrd(ig,k).gt.0.) .or. (detrd(ig,k).gt.0.)) then
101	!! if (entrd(ig,k).gt.detrd(ig,k)) then
102	! lmax_down(ig)=min(k,lmax(ig))
103	! endif
104	! enddo
105	! enddo
106	! do k=nlayermx,1,-1
107	! do ig=1,ngridmx
108	! if ((entrd(ig,k).gt.0.) .or. (detrd(ig,k).gt.0.)) then
109	!! if (detrd(ig,k).gt.entrd(ig,k)) then
110	! lmin_down(ig)=k
111	! endif
112	! enddo
113	! enddo
114	!
115	!
116	! fmd(:,:)=0.
117	!
118	! do ig=1,ngridmx
119	! if ((lmax_down(ig).gt.1) .and. ((lmax_down(ig)-lmin_down(ig)).gt.1)) then
120	!! fmd(ig,lmax_down(ig))=0.
121	!! entrd(ig,lmax_down(ig))=detrd(ig,lmax_down(ig))
122	!! detrd(ig,lmax_down(ig))=0.
123	!! print*,lmin_down(ig),lmax_down(ig),lmax(ig)
124	!
125	! fmd(ig,lmax_down(ig)+1)=0.
126	!
127	! do k=lmax_down(ig),lmin_down(ig)+1,-1
128	! fmd(ig,k)=fmd(ig,k+1)+entrd(ig,k)-detrd(ig,k)
129	! enddo
130	!
131	! fmd(ig,lmin_down(ig))=0.
132	! detrd(ig,lmin_down(ig))=fmd(ig,lmin_down(ig)+1)+entrd(ig,lmin_down(ig))
133	!
134	! else
135	! entrd(ig,:)=0.
136	! detrd(ig,:)=0.
137	! active(ig,:)=.false.
138	! endif
139	!
140	! enddo
141	! ncorec=0
142	! do k=nlayermx,2,-1
143	! do ig=1,ngridmx
144	! if (fmd(ig,k).lt.0.) then
145	!! detrd(ig,k)=max(0.,detrd(ig,k)+fmd(ig,k-1))
146	!! fmd(ig,k-1)=0.
147	!! entrd(ig,k-1)=0.
148	!! detrd(ig,k-1)=0.
149	!! lmin_down(ig)=k-1
150	! fmd(ig,k)=fmd(ig,k+1)
151	! detrd(ig,k)=entrd(ig,k)
152	! ncorec=ncorec+1
153	!! fmd(ig,k)=0.
154	!! detrd(ig,k)=entrd(ig,k)+fmd(ig,k+1)
155	!
156	! endif
157	! enddo
158	! enddo
159	!
160	! if (ncorec .ne. 0) then
161	! print*, 'corrections for negative downward mass flux :',ncorec
162	! endif
163	! print*, lmin_down(:),lmax_down(:)
164	!
165	! do k=2,nlayermx
166	! do ig=1,ngridmx
167	! active(ig,k)=(k.ge.lmin_down(ig)).and.(k.le.lmax_down(ig)) &
168	! & .and.(((fmd(ig,k)+detrd(ig,k))ptimestep).gt.1.e-6masse0(ig,k))
169	! enddo
170	! enddo
171	!
172	!
173	! do ig=1,ngridmx
174	! do k=lmin_down(ig),lmax_down(ig)
175	! if(.not.active(ig,k)) then
176	! active(ig,:)=.false.
177	! endif
178	! enddo
179	! enddo
180	!
181	! if(charvar .eq. 'tke') then
182	! active(:,:)=.false.
183	! endif
184	!
185	!! do ig=1,ngridmx
186	!! active(ig,lmax_down(ig))=(((fmd(ig,lmax_down(ig))+detrd(ig,lmax_down(ig)))* &
187	!! & ptimestep).gt.1.e-8*masse0(ig,lmax_down(ig)))
188	!! enddo
189	!!
190	!! do ig=1,ngridmx
191	!! if (lmax_down(ig).gt.1) then
192	!! do k=lmax_down(ig)-1,lmin_down(ig),-1
193	!! active(ig,k)=(((fmd(ig,k)+detrd(ig,k))* &
194	!! & ptimestep).gt.1.e-8*masse0(ig,k)) &
195	!! & .and. active(ig,k+1)
196	!! enddo
197	!! else
198	!! active(ig,:)=.false.
199	!! endif
200	!! enddo
201	!!
202	!! ========== qa : q in updraft ==================================
203	!
204	do k=2,nlayermx
205	do ig=1,ngridmx
206	if ((fm0(ig,k+1)+detr0(ig,k))*ptimestep.gt. &
207	& 1.e-5*masse0(ig,k)) then
208	qa(ig,k)=(fm0(ig,k)qa(ig,k-1)+entr0(ig,k)q(ig,k)) &
209	& /(fm0(ig,k+1)+detr0(ig,k))
210
211	if ((qa(ig,k).lt.0.) .and. (charvar .ne. 'momentum')) then
212	print*,'qa<0!!!',charvar,ig,k,fm0(ig,k),qa(ig,k-1),entr0(ig,k),q(ig,k),fm0(ig,k+1),detr0(ig,k)
213	print*,'---------> Cancelling qa'
214	qa(ig,k)=q(ig,k)
215	endif
216	else
217	qa(ig,k)=q(ig,k)
218	endif
219	enddo
220	enddo
221
222	! ========== qd : q in downdraft =================================
223	!
224	! do k=nlayermx-1,1,-1
225	! do ig=1,ngridmx
226	! if (active(ig,k)) then
227	! qd(ig,k)=(fmd(ig,k+1)qd(ig,k+1)+entrd(ig,k)q(ig,k)) &
228	! & /(fmd(ig,k)+detrd(ig,k))
229	!
230	! if ((qd(ig,k).lt.0.) .and. (charvar .ne. 'momentum')) then
231	! print*,'qd<0!!!',charvar,ig,k,fmd(ig,k),qd(ig,k),entrd(ig,k),q(ig,k),fmd(ig,k+1),detrd(ig,k),lmin_down(ig),lmax_down(ig)
232	! print*, '---------> cancelling qd, no downdraft for this gridpoint'
233	! qd(ig,k)=q(ig,k)
234	! active(ig,:)=.false.
235	! endif
236	! else
237	! qd(ig,k)=q(ig,k)
238	! endif
239	!! print*,'active,k,entr,detr,q,qd (down) :',active(ig,k),k,entrd(ig,k),detrd(ig,k),q(ig,k),qd(ig,k)
240	! enddo
241	! enddo
242	!
243	!
244	! ====== dq ======================================================
245
246	do ig=1,ngridmx
247	if(active(ig,1)) then
248
249	dq_therm(ig,1)=(detr0(ig,1)qa(ig,1)+detrd(ig,1)qd(ig,1) &
250	& +fm0(ig,2)*q(ig,2) &
251	& -entr0(ig,1)q(ig,1)-entrd(ig,1)q(ig,1) &
252	& -fmd(ig,2)*q(ig,1)) &
253	& *ptimestep/masse0(ig,1)
254
255	else
256	dq_therm(ig,1)=(detr0(ig,1)qa(ig,1)+fm0(ig,2)q(ig,2) &
257	& -entr0(ig,1)*q(ig,1)) &
258	& *ptimestep/masse0(ig,1)
259
260	endif
261	enddo
262
263	do k=2,nlayermx-1
264	do ig=1, ngridmx
265
266	if(active(ig,k)) then
267
268	dq_therm(ig,k)=(detr0(ig,k)qa(ig,k)+detrd(ig,k)qd(ig,k) &
269	& +fm0(ig,k+1)q(ig,k+1)+fmd(ig,k)q(ig,k-1) &
270	& -entr0(ig,k)q(ig,k)-entrd(ig,k)q(ig,k) &
271	& -fm0(ig,k)q(ig,k)-fmd(ig,k+1)q(ig,k)) &
272	& *ptimestep/masse0(ig,k)
273
274
275	else
276	dq_therm(ig,k)=(detr0(ig,k)qa(ig,k)+fm0(ig,k+1)q(ig,k+1) &
277	& -entr0(ig,k)q(ig,k)-fm0(ig,k)q(ig,k)) &
278	& *ptimestep/masse0(ig,k)
279
280
281	endif
282
283	enddo
284	enddo
285
286	do ig=1, ngridmx
287
288	if(active(ig,nlayermx)) then
289
290	dq_therm(ig,nlayermx)=(detr0(ig,nlayermx)qa(ig,nlayermx)+detrd(ig,nlayermx)qd(ig,nlayermx) &
291	& +fmd(ig,nlayermx)*q(ig,nlayermx-1) &
292	& -entr0(ig,nlayermx)q(ig,nlayermx)-entrd(ig,nlayermx)q(ig,nlayermx) &
293	& -fm0(ig,nlayermx)*q(ig,nlayermx)) &
294	& *ptimestep/masse0(ig,nlayermx)
295
296	else
297	dq_therm(ig,nlayermx)=(detr0(ig,nlayermx)*qa(ig,nlayermx) &
298	& -entr0(ig,nlayermx)q(ig,nlayermx)-fm0(ig,nlayermx)q(ig,nlayermx)) &
299	& *ptimestep/masse0(ig,nlayermx)
300	endif
301
302	enddo
303	return
304	end

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

Download in other formats:

Original Format