Context Navigation

watercloud_mod.F @ 2538

Last change on this file since 2538 was 2516, checked in by mvals, 4 years ago

Mars GCM:
Correction of the threshold value used to prevent from negative values in watercloud_mod.F: the threshold is set to 1e-16 (=qperemin, an already
existing parameter for the water cycle), instead of 1e-8. This former value was at the origin of a strong and unrealistic drying out of the upper
atmosphere, when the photochemistry was active, during the first part of the martian year (months 3 and 4). (bug reported by Francisco
Gonzalez-Galindo).
MV

File size: 26.7 KB

Line
1	MODULE watercloud_mod
2
3	IMPLICIT NONE
4
5	REAL,SAVE,ALLOCATABLE :: zdqcloud(:,:,:) ! tendencies on pq due to condensation of H2O(kg/kg.s-1)
6	REAL,SAVE,ALLOCATABLE :: zdqscloud(:,:) ! tendencies on qsurf (calculated only by calchim but declared here)
7
8	CONTAINS
9
10	SUBROUTINE watercloud(ngrid,nlay,ptimestep,
11	& pplev,pplay,pdpsrf,pzlay,pt,pdt,
12	& pq,pdq,pdqcloud,pdtcloud,
13	& nq,tau,tauscaling,rdust,rice,nuice,
14	& rsedcloud,rhocloud,totcloudfrac)
15	USE ioipsl_getin_p_mod, ONLY : getin_p
16	USE updaterad, ONLY: updaterdust, updaterice_micro,
17	& updaterice_typ
18	USE improvedclouds_mod, ONLY: improvedclouds
19	USE watersat_mod, ONLY: watersat
20	use tracer_mod, only: nqmx, igcm_h2o_vap, igcm_h2o_ice,
21	& igcm_hdo_vap, igcm_hdo_ice,
22	& igcm_dust_mass, igcm_dust_number,
23	& igcm_ccn_mass, igcm_ccn_number,
24	& rho_dust, nuice_sed, nuice_ref,
25	& qperemin
26	use dimradmars_mod, only: naerkind
27	IMPLICIT NONE
28
29
30	c=======================================================================
31	c Water-ice cloud formation
32	c
33	c Includes two different schemes:
34	c - A simplified scheme (see simpleclouds.F)
35	c - An improved microphysical scheme (see improvedclouds.F)
36	c
37	c There is a time loop specific to cloud formation
38	c due to timescales smaller than the GCM integration timestep.
39	c
40	c Authors: Franck Montmessin, Francois Forget, Ehouarn Millour,
41	c J.-B. Madeleine, Thomas Navarro
42	c
43	c 2004 - 2012
44	c=======================================================================
45
46	c-----------------------------------------------------------------------
47	c declarations:
48	c -------------
49
50	include "callkeys.h"
51
52	c Inputs/outputs:
53	c ------
54
55	INTEGER, INTENT(IN) :: ngrid,nlay
56	INTEGER, INTENT(IN) :: nq ! nombre de traceurs
57	REAL, INTENT(IN) :: ptimestep ! pas de temps physique (s)
58	REAL, INTENT(IN) :: pplev(ngrid,nlay+1) ! pression aux inter-couches (Pa)
59	REAL, INTENT(IN) :: pplay(ngrid,nlay) ! pression au milieu des couches (Pa)
60	REAL, INTENT(IN) :: pdpsrf(ngrid) ! tendence surf pressure
61	REAL, INTENT(IN) :: pzlay(ngrid,nlay) ! altitude at the middle of the layers
62	REAL, INTENT(IN) :: pt(ngrid,nlay) ! temperature at the middle of the layers (K)
63	REAL, INTENT(IN) :: pdt(ngrid,nlay) ! tendence temperature des autres param.
64
65	REAL, INTENT(IN) :: pq(ngrid,nlay,nq) ! traceur (kg/kg)
66	rEAL, INTENT(IN) :: pdq(ngrid,nlay,nq) ! tendence avant condensation (kg/kg.s-1)
67
68	REAL, INTENT(IN) :: tau(ngrid,naerkind) ! Column dust optical depth at each point
69	REAL, INTENT(IN) :: tauscaling(ngrid) ! Convertion factor for dust amount
70	REAL, INTENT(INOUT) :: rdust(ngrid,nlay) ! Dust geometric mean radius (m)
71
72	REAL, INTENT(OUT) :: pdqcloud(ngrid,nlay,nq) ! tendence de la condensation H2O(kg/kg.s-1)
73	REAL, INTENT(OUT) :: pdtcloud(ngrid,nlay) ! tendence temperature due
74	! a la chaleur latente
75	REAL, INTENT(INOUT) :: rice(ngrid,nlay) ! Ice mass mean radius (m)
76	! (r_c in montmessin_2004)
77	REAL, INTENT(OUT) :: nuice(ngrid,nlay) ! Estimated effective variance
78	! of the size distribution
79	REAL, INTENT(OUT) :: rsedcloud(ngrid,nlay) ! Cloud sedimentation radius
80	REAL, INTENT(OUT) :: rhocloud(ngrid,nlay) ! Cloud density (kg.m-3)
81
82	REAL, INTENT(INOUT):: totcloudfrac(ngrid) ! Cloud fraction (A. Pottier 2013)
83
84	c Locals:
85	c ------
86
87	! for ice radius computation
88	REAL Mo,No
89	REAl ccntyp
90
91	! for time loop
92	INTEGER microstep ! time subsampling step variable
93	INTEGER,SAVE :: imicro ! time subsampling for coupled water microphysics & sedimentation
94	REAL,SAVE :: microtimestep ! integration timestep for coupled water microphysics & sedimentation
95	REAL,SAVE :: microtimestep_prev=-999
96
97	! tendency given by clouds (inside the micro loop)
98	REAL subpdqcloud(ngrid,nlay,nq) ! cf. pdqcloud
99	REAL subpdtcloud(ngrid,nlay) ! cf. pdtcloud
100
101	! global tendency (clouds+physics)
102	REAL sum_subpdq(ngrid,nlay,nq) ! cf. pdqcloud
103	REAL sum_subpdt(ngrid,nlay) ! cf. pdtcloud
104
105	! no supersaturation when option supersat is false
106	REAL zt(ngrid,nlay) ! local value of temperature
107	REAL zqsat(ngrid,nlay) ! saturation
108
109	INTEGER iq,ig,l
110	LOGICAL,SAVE :: firstcall=.true.
111
112	! Representation of sub-grid water ice clouds A. Pottier 2013
113	REAL :: ztclf(ngrid, nlay)
114	REAL :: zqclf(ngrid, nlay,nq)
115	REAL :: zdelt
116	REAL :: norm
117	REAL :: ponder
118	REAL :: tcond(ngrid,nlay)
119	REAL :: zqvap(ngrid,nlay)
120	REAL :: zqice(ngrid,nlay)
121	REAL :: spant ! delta T for the temperature distribution
122	! REAL :: zqsat(ngrid,nlay) ! saturation
123	REAL :: pteff(ngrid, nlay)! effective temperature in the cloud,neb
124	REAL :: pqeff(ngrid, nlay, nq)! effective tracers quantities in the cloud
125	REAL :: cloudfrac(ngrid,nlay) ! cloud fraction
126	REAL :: mincloud ! min cloud frac
127	LOGICAL, save :: flagcloud=.true.
128
129	c ** un petit test de coherence
130	c --------------------------
131
132	IF (firstcall) THEN
133
134	if (nq.gt.nqmx) then
135	write(,) 'stop in watercloud (nq.gt.nqmx)!'
136	write(,) 'nq=',nq,' nqmx=',nqmx
137	call abort_physic("watercloud","nq.gt.nqmx",1)
138	endif
139
140	write(,) "watercloud: igcm_h2o_vap=",igcm_h2o_vap
141	write(,) " igcm_h2o_ice=",igcm_h2o_ice
142
143	write(,) "time subsampling for microphysic ?"
144	#ifdef MESOSCALE
145	imicro = 2
146	#else
147	imicro = 30
148	#endif
149	call getin_p("imicro",imicro)
150	write(,)"watercloud: imicro = ",imicro
151
152	firstcall=.false.
153	ENDIF ! of IF (firstcall)
154
155	!! AS: moved out of firstcall to allow nesting+evoluting timestep
156	!! TBD: consider possible diff imicro with domains?
157	microtimestep = ptimestep/real(imicro)
158	if (microtimestep/=microtimestep_prev) then
159	! only tell the world if microtimestep has changed
160	write(,)"watercloud: Physical timestep is ",ptimestep
161	write(,)"watercloud: Microphysics timestep is ",microtimestep
162	microtimestep_prev=microtimestep
163	endif
164
165	c-----Initialization
166	sum_subpdq(1:ngrid,1:nlay,1:nq) = 0
167	sum_subpdt(1:ngrid,1:nlay) = 0
168
169	! default value if no ice
170	rhocloud(1:ngrid,1:nlay) = rho_dust
171
172	c-------------------------------------------------------------------
173	c 0. Representation of sub-grid water ice clouds
174	c------------------
175	c-----Initialization
176	pteff(1:ngrid,1:nlay) = 0
177	pqeff(1:ngrid,1:nlay,1:nq) = 0
178	DO l=1,nlay
179	DO ig=1,ngrid
180	pteff(ig,l)=pt(ig,l)
181	END DO
182	END DO
183	DO l=1,nlay
184	DO ig=1,ngrid
185	DO iq=1,nq
186	pqeff(ig,l,iq)=pq(ig,l,iq)
187	ENDDO
188	ENDDO
189	ENDDO
190	c-----Tendencies
191	DO l=1,nlay
192	DO ig=1,ngrid
193	ztclf(ig,l)=pt(ig,l)+ pdt(ig,l)*ptimestep
194	ENDDO
195	ENDDO
196	DO l=1,nlay
197	DO ig=1,ngrid
198	DO iq=1,nq
199	zqclf(ig,l,iq)=pq(ig,l,iq)+pdq(ig,l,iq)*ptimestep
200	ENDDO
201	ENDDO
202	ENDDO
203	c-----Effective temperature calculation
204	IF (CLFvarying) THEN
205	spant=3.0 ! delta T for the temprature distribution
206	mincloud=0.1 ! min cloudfrac when there is ice
207	IF (flagcloud) THEN
208	write(,) "Delta T", spant
209	write(,) "mincloud", mincloud
210	flagcloud=.false.
211	END IF
212	!CALL watersat(ngrid*nlay,ztclf,pplay,zqsat) !MV17: we dont need zqsat in the CLFvarying scheme
213	zqvap=zqclf(:,:,igcm_h2o_vap)
214	zqice=zqclf(:,:,igcm_h2o_ice)
215	CALL tcondwater(ngrid*nlay,pplay,zqvap+zqice,tcond)
216	DO l=1,nlay
217	DO ig=1,ngrid
218	zdelt=spant !MAX(spant*ztclf(ig,l),1.e-12), now totally in K. Fixed width
219	IF (tcond(ig,l) .ge. (ztclf(ig,l)+zdelt)) THEN
220	pteff(ig,l)=ztclf(ig,l)
221	cloudfrac(ig,l)=1.
222	ELSE IF (tcond(ig,l) .le. (ztclf(ig,l)-zdelt)) THEN
223	pteff(ig,l)=ztclf(ig,l)-zdelt
224	cloudfrac(ig,l)=mincloud
225	ELSE
226	cloudfrac(ig,l)=(tcond(ig,l)-ztclf(ig,l)+zdelt)/
227	& (2.0*zdelt)
228	pteff(ig,l)=(tcond(ig,l)+ztclf(ig,l)-zdelt)/2.
229	END IF
230	pteff(ig,l)=pteff(ig,l)-pdt(ig,l)*ptimestep
231	IF (cloudfrac(ig,l).le.mincloud) THEN !MV17: replaced .le.0 by .le.mincloud
232	cloudfrac(ig,l)=mincloud
233	ELSE IF (cloudfrac(ig,l).gt.1) THEN
234	cloudfrac(ig,l)=1.
235	END IF
236	ENDDO
237	ENDDO
238	c-----Calculation of the total cloud coverage of the column
239	DO ig=1,ngrid
240	totcloudfrac(ig) = 0.
241	norm=0.
242	DO l=1,nlay
243	ponder=zqice(ig,l)*(pplev(ig,l) - pplev(ig,l+1))
244	totcloudfrac(ig) = totcloudfrac(ig)
245	& + cloudfrac(ig,l)*ponder
246	norm=norm+ponder
247	ENDDO
248	totcloudfrac(ig)=MAX(totcloudfrac(ig)/norm,1.e-12) ! min value if NaNs
249	ENDDO
250	c-----Effective tracers quantities in the cloud fraction
251	IF (microphys) THEN
252	pqeff(:,:,igcm_ccn_mass)=pq(:,:,igcm_ccn_mass)/
253	& cloudfrac(:,:)
254	pqeff(:,:,igcm_ccn_number)=pq(:,:,igcm_ccn_number)/
255	& cloudfrac(:,:)
256	END IF ! end if (microphys)
257	pqeff(:,:,igcm_h2o_ice)=pq(:,:,igcm_h2o_ice)/
258	& cloudfrac(:,:)
259	!! CLFvarying outputs
260	CALL WRITEDIAGFI(ngrid,'pqeffice','pqeffice',
261	& 'kg/kg',3,pqeff(:,:,igcm_h2o_ice))
262	CALL WRITEDIAGFI(ngrid,'pteff','pteff',
263	& 'K',3,pteff(:,:))
264	CALL WRITEDIAGFI(ngrid,'tcond','tcond',
265	& 'K',3,tcond(:,:))
266	CALL WRITEDIAGFI(ngrid,'cloudfrac','cloudfrac',
267	& 'K',3,cloudfrac(:,:))
268	END IF ! end if (CLFvarying)
269	c------------------------------------------------------------------
270	c Time subsampling for microphysics
271	c------------------------------------------------------------------
272	rhocloud(1:ngrid,1:nlay) = rho_dust
273	DO microstep=1,imicro
274
275	c-------------------------------------------------------------------
276	c 1. Tendencies:
277	c------------------
278
279
280	c------ Temperature tendency subpdt
281	! Each microtimestep we give the cloud scheme a stepped entry subpdt instead of pdt
282	! If imicro=1 subpdt is the same as pdt
283	DO l=1,nlay
284	DO ig=1,ngrid
285	sum_subpdt(ig,l) = sum_subpdt(ig,l)
286	& + pdt(ig,l) ! At each micro timestep we add pdt in order to have a stepped entry
287	ENDDO
288	ENDDO
289	c------ Tracers tendencies subpdq are additionned
290	c------ At each micro timestep we add pdq in order to have a stepped entry
291	IF (microphys) THEN
292	DO l=1,nlay
293	DO ig=1,ngrid
294	sum_subpdq(ig,l,igcm_dust_mass) =
295	& sum_subpdq(ig,l,igcm_dust_mass)
296	& + pdq(ig,l,igcm_dust_mass)
297	sum_subpdq(ig,l,igcm_dust_number) =
298	& sum_subpdq(ig,l,igcm_dust_number)
299	& + pdq(ig,l,igcm_dust_number)
300	sum_subpdq(ig,l,igcm_ccn_mass) =
301	& sum_subpdq(ig,l,igcm_ccn_mass)
302	& + pdq(ig,l,igcm_ccn_mass)
303	sum_subpdq(ig,l,igcm_ccn_number) =
304	& sum_subpdq(ig,l,igcm_ccn_number)
305	& + pdq(ig,l,igcm_ccn_number)
306	ENDDO
307	ENDDO
308	ENDIF
309	DO l=1,nlay
310	DO ig=1,ngrid
311	sum_subpdq(ig,l,igcm_h2o_ice) =
312	& sum_subpdq(ig,l,igcm_h2o_ice)
313	& + pdq(ig,l,igcm_h2o_ice)
314	sum_subpdq(ig,l,igcm_h2o_vap) =
315	& sum_subpdq(ig,l,igcm_h2o_vap)
316	& + pdq(ig,l,igcm_h2o_vap)
317	IF (hdo) THEN
318	sum_subpdq(ig,l,igcm_hdo_ice) =
319	& sum_subpdq(ig,l,igcm_hdo_ice)
320	& + pdq(ig,l,igcm_hdo_ice)
321	sum_subpdq(ig,l,igcm_hdo_vap) =
322	& sum_subpdq(ig,l,igcm_hdo_vap)
323	& + pdq(ig,l,igcm_hdo_vap)
324	ENDIF !hdo
325	ENDDO
326	ENDDO
327
328	c-------------------------------------------------------------------
329	c 2. Main call to the different cloud schemes:
330	c------------------------------------------------
331	IF (microphys) THEN
332	CALL improvedclouds(ngrid,nlay,microtimestep,
333	& pplay,pteff,sum_subpdt,
334	& pqeff,sum_subpdq,subpdqcloud,subpdtcloud,
335	& nq,tauscaling)
336
337	ELSE
338	CALL simpleclouds(ngrid,nlay,microtimestep,
339	& pplay,pzlay,pteff,sum_subpdt,
340	& pqeff,sum_subpdq,subpdqcloud,subpdtcloud,
341	& nq,tau,rice)
342	ENDIF
343
344	c-------------------------------------------------------------------
345	c 3. Updating tendencies after cloud scheme:
346	c-----------------------------------------------
347
348	IF (microphys) THEN
349	DO l=1,nlay
350	DO ig=1,ngrid
351	sum_subpdq(ig,l,igcm_dust_mass) =
352	& sum_subpdq(ig,l,igcm_dust_mass)
353	& + subpdqcloud(ig,l,igcm_dust_mass)
354	sum_subpdq(ig,l,igcm_dust_number) =
355	& sum_subpdq(ig,l,igcm_dust_number)
356	& + subpdqcloud(ig,l,igcm_dust_number)
357	sum_subpdq(ig,l,igcm_ccn_mass) =
358	& sum_subpdq(ig,l,igcm_ccn_mass)
359	& + subpdqcloud(ig,l,igcm_ccn_mass)
360	sum_subpdq(ig,l,igcm_ccn_number) =
361	& sum_subpdq(ig,l,igcm_ccn_number)
362	& + subpdqcloud(ig,l,igcm_ccn_number)
363	ENDDO
364	ENDDO
365	ENDIF
366	DO l=1,nlay
367	DO ig=1,ngrid
368	sum_subpdq(ig,l,igcm_h2o_ice) =
369	& sum_subpdq(ig,l,igcm_h2o_ice)
370	& + subpdqcloud(ig,l,igcm_h2o_ice)
371	sum_subpdq(ig,l,igcm_h2o_vap) =
372	& sum_subpdq(ig,l,igcm_h2o_vap)
373	& + subpdqcloud(ig,l,igcm_h2o_vap)
374
375	IF (hdo) THEN
376	sum_subpdq(ig,l,igcm_hdo_ice) =
377	& sum_subpdq(ig,l,igcm_hdo_ice)
378	& + subpdqcloud(ig,l,igcm_hdo_ice)
379	sum_subpdq(ig,l,igcm_hdo_vap) =
380	& sum_subpdq(ig,l,igcm_hdo_vap)
381	& + subpdqcloud(ig,l,igcm_hdo_vap)
382	ENDIF ! hdo
383
384	ENDDO
385	ENDDO
386
387	IF (activice) THEN
388	DO l=1,nlay
389	DO ig=1,ngrid
390	sum_subpdt(ig,l) =
391	& sum_subpdt(ig,l) + subpdtcloud(ig,l)
392	ENDDO
393	ENDDO
394	ENDIF
395
396	! !! Example of how to use writediagmicrofi useful to
397	! !! get outputs at each microphysical sub-timestep (better to be used in 1D)
398	! CALL WRITEDIAGMICROFI(ngrid,imicro,microstep,
399	! & microtimestep,'subpdtcloud',
400	! & 'subpdtcloud','K/s',1,subpdtcloud(:,:))
401
402	ENDDO ! of DO microstep=1,imicro
403
404	c-------------------------------------------------------------------
405	c 6. Compute final tendencies after time loop:
406	c------------------------------------------------
407
408	c------ Temperature tendency pdtcloud
409	DO l=1,nlay
410	DO ig=1,ngrid
411	pdtcloud(ig,l) =
412	& sum_subpdt(ig,l)/real(imicro)-pdt(ig,l)
413	ENDDO
414	ENDDO
415
416	c------ Tracers tendencies pdqcloud
417	DO l=1,nlay
418	DO ig=1,ngrid
419	pdqcloud(ig,l,igcm_h2o_ice) =
420	& sum_subpdq(ig,l,igcm_h2o_ice)/real(imicro)
421	& - pdq(ig,l,igcm_h2o_ice)
422	pdqcloud(ig,l,igcm_h2o_vap) =
423	& sum_subpdq(ig,l,igcm_h2o_vap)/real(imicro)
424	& - pdq(ig,l,igcm_h2o_vap)
425	IF (hdo) THEN
426	pdqcloud(ig,l,igcm_hdo_ice) =
427	& sum_subpdq(ig,l,igcm_hdo_ice)/real(imicro)
428	& - pdq(ig,l,igcm_hdo_ice)
429	pdqcloud(ig,l,igcm_hdo_vap) =
430	& sum_subpdq(ig,l,igcm_hdo_vap)/real(imicro)
431	& - pdq(ig,l,igcm_hdo_vap)
432	ENDIF !hdo
433	ENDDO
434	ENDDO
435
436	IF(microphys) THEN
437	DO l=1,nlay
438	DO ig=1,ngrid
439	pdqcloud(ig,l,igcm_ccn_mass) =
440	& sum_subpdq(ig,l,igcm_ccn_mass)/real(imicro)
441	& - pdq(ig,l,igcm_ccn_mass)
442	pdqcloud(ig,l,igcm_ccn_number) =
443	& sum_subpdq(ig,l,igcm_ccn_number)/real(imicro)
444	& - pdq(ig,l,igcm_ccn_number)
445	ENDDO
446	ENDDO
447	ENDIF
448
449	IF(scavenging) THEN
450	DO l=1,nlay
451	DO ig=1,ngrid
452	pdqcloud(ig,l,igcm_dust_mass) =
453	& sum_subpdq(ig,l,igcm_dust_mass)/real(imicro)
454	& - pdq(ig,l,igcm_dust_mass)
455	pdqcloud(ig,l,igcm_dust_number) =
456	& sum_subpdq(ig,l,igcm_dust_number)/real(imicro)
457	& - pdq(ig,l,igcm_dust_number)
458	ENDDO
459	ENDDO
460	ENDIF
461
462	c------- Due to stepped entry, other processes tendencies can add up to negative values
463	c------- Therefore, enforce positive values and conserve mass
464	IF(microphys) THEN
465	DO l=1,nlay
466	DO ig=1,ngrid
467	IF ((pq(ig,l,igcm_ccn_number) +
468	& ptimestep* (pdq(ig,l,igcm_ccn_number) +
469	& pdqcloud(ig,l,igcm_ccn_number)) .le. 1.)
470	& .or. (pq(ig,l,igcm_ccn_mass) +
471	& ptimestep* (pdq(ig,l,igcm_ccn_mass) +
472	& pdqcloud(ig,l,igcm_ccn_mass)) .le. 1.e-20)) THEN
473	pdqcloud(ig,l,igcm_ccn_number) =
474	& - pq(ig,l,igcm_ccn_number)/ptimestep
475	& - pdq(ig,l,igcm_ccn_number) + 1.
476	pdqcloud(ig,l,igcm_dust_number) =
477	& -pdqcloud(ig,l,igcm_ccn_number)
478	pdqcloud(ig,l,igcm_ccn_mass) =
479	& - pq(ig,l,igcm_ccn_mass)/ptimestep
480	& - pdq(ig,l,igcm_ccn_mass) + 1.e-20
481	pdqcloud(ig,l,igcm_dust_mass) =
482	& -pdqcloud(ig,l,igcm_ccn_mass)
483	ENDIF
484	ENDDO
485	ENDDO
486	ENDIF
487
488	IF(scavenging) THEN
489	DO l=1,nlay
490	DO ig=1,ngrid
491	IF ((pq(ig,l,igcm_dust_number) +
492	& ptimestep* (pdq(ig,l,igcm_dust_number) +
493	& pdqcloud(ig,l,igcm_dust_number)) .le. 1.)
494	& .or. (pq(ig,l,igcm_dust_mass) +
495	& ptimestep* (pdq(ig,l,igcm_dust_mass) +
496	& pdqcloud(ig,l,igcm_dust_mass)) .le. 1.e-20)) THEN
497	pdqcloud(ig,l,igcm_dust_number) =
498	& - pq(ig,l,igcm_dust_number)/ptimestep
499	& - pdq(ig,l,igcm_dust_number) + 1.
500	pdqcloud(ig,l,igcm_ccn_number) =
501	& -pdqcloud(ig,l,igcm_dust_number)
502	pdqcloud(ig,l,igcm_dust_mass) =
503	& - pq(ig,l,igcm_dust_mass)/ptimestep
504	& - pdq(ig,l,igcm_dust_mass) + 1.e-20
505	pdqcloud(ig,l,igcm_ccn_mass) =
506	& -pdqcloud(ig,l,igcm_dust_mass)
507	ENDIF
508	ENDDO
509	ENDDO
510	ENDIF
511
512	DO l=1,nlay
513	DO ig=1,ngrid
514
515	IF (pq(ig,l,igcm_h2o_ice) + ptimestep*
516	& (pdq(ig,l,igcm_h2o_ice) + pdqcloud(ig,l,igcm_h2o_ice))
517	& .le. qperemin) THEN
518	pdqcloud(ig,l,igcm_h2o_ice) =
519	& - pq(ig,l,igcm_h2o_ice)/ptimestep - pdq(ig,l,igcm_h2o_ice)
520	pdqcloud(ig,l,igcm_h2o_vap) = -pdqcloud(ig,l,igcm_h2o_ice)
521	if (hdo) then
522	pdqcloud(ig,l,igcm_hdo_ice) =
523	& - pq(ig,l,igcm_hdo_ice)/ptimestep - pdq(ig,l,igcm_hdo_ice)
524	pdqcloud(ig,l,igcm_hdo_vap) = -pdqcloud(ig,l,igcm_hdo_ice)
525	endif
526	ENDIF
527	IF (pq(ig,l,igcm_h2o_vap) + ptimestep*
528	& (pdq(ig,l,igcm_h2o_vap) + pdqcloud(ig,l,igcm_h2o_vap))
529	& .le. qperemin) THEN
530	pdqcloud(ig,l,igcm_h2o_vap) =
531	& - pq(ig,l,igcm_h2o_vap)/ptimestep - pdq(ig,l,igcm_h2o_vap)
532	pdqcloud(ig,l,igcm_h2o_ice) = -pdqcloud(ig,l,igcm_h2o_vap)
533	if (hdo) then
534	pdqcloud(ig,l,igcm_hdo_vap) =
535	& - pq(ig,l,igcm_hdo_vap)/ptimestep - pdq(ig,l,igcm_hdo_vap)
536	pdqcloud(ig,l,igcm_hdo_ice) = -pdqcloud(ig,l,igcm_hdo_vap)
537	endif
538	ENDIF
539
540	ENDDO
541	ENDDO
542
543	c------Update the ice and dust particle size "rice" for output or photochemistry
544	c------Only rsedcloud is used for the water cycle
545
546	IF(scavenging) THEN
547	DO l=1, nlay
548	DO ig=1,ngrid
549
550	call updaterdust(
551	& pq(ig,l,igcm_dust_mass) + ! dust mass
552	& (pdq(ig,l,igcm_dust_mass) + ! dust mass
553	& pdqcloud(ig,l,igcm_dust_mass))*ptimestep, ! dust mass
554	& pq(ig,l,igcm_dust_number) + ! dust number
555	& (pdq(ig,l,igcm_dust_number) + ! dust number
556	& pdqcloud(ig,l,igcm_dust_number))*ptimestep, ! dust number
557	& rdust(ig,l))
558
559	ENDDO
560	ENDDO
561	ENDIF
562
563	IF(microphys) THEN
564
565	! In case one does not want to allow supersatured water when using microphysics.
566	! Not done by default.
567	IF(.not.supersat) THEN
568	zt = pt + (pdt+pdtcloud)*ptimestep
569	call watersat(ngrid*nlay,zt,pplay,zqsat)
570	DO l=1, nlay
571	DO ig=1,ngrid
572	IF (pq(ig,l,igcm_h2o_vap)
573	& + (pdq(ig,l,igcm_h2o_vap) + pdqcloud(ig,l,igcm_h2o_vap))
574	& * ptimestep .ge. zqsat(ig,l)) THEN
575	pdqcloud(ig,l,igcm_h2o_vap) =
576	& (zqsat(ig,l) - pq(ig,l,igcm_h2o_vap))/ptimestep
577	& - pdq(ig,l,igcm_h2o_vap)
578	pdqcloud(ig,l,igcm_h2o_ice) =
579	& -pdqcloud(ig,l,igcm_h2o_vap)
580	! no need to correct ccn_number, updaterad can handle this properly.
581	ENDIF
582	ENDDO
583	ENDDO
584	ENDIF
585
586	DO l=1, nlay
587	DO ig=1,ngrid
588
589	call updaterice_micro(
590	& pq(ig,l,igcm_h2o_ice) + ! ice mass
591	& (pdq(ig,l,igcm_h2o_ice) + ! ice mass
592	& pdqcloud(ig,l,igcm_h2o_ice))*ptimestep, ! ice mass
593	& pq(ig,l,igcm_ccn_mass) + ! ccn mass
594	& (pdq(ig,l,igcm_ccn_mass) + ! ccn mass
595	& pdqcloud(ig,l,igcm_ccn_mass))*ptimestep, ! ccn mass
596	& pq(ig,l,igcm_ccn_number) + ! ccn number
597	& (pdq(ig,l,igcm_ccn_number) + ! ccn number
598	& pdqcloud(ig,l,igcm_ccn_number))*ptimestep, ! ccn number
599	& tauscaling(ig),rice(ig,l),rhocloud(ig,l))
600
601	ENDDO
602	ENDDO
603
604	ELSE ! no microphys
605
606	DO l=1,nlay
607	DO ig=1,ngrid
608
609	call updaterice_typ(
610	& pq(ig,l,igcm_h2o_ice) + ! ice mass
611	& (pdq(ig,l,igcm_h2o_ice) + ! ice mass
612	& pdqcloud(ig,l,igcm_h2o_ice))*ptimestep, ! ice mass
613	& tau(ig,1),pzlay(ig,l),rice(ig,l))
614
615	ENDDO
616	ENDDO
617
618	ENDIF ! of IF(microphys)
619
620
621
622	c A correction if a lot of subliming CO2 fills the 1st layer FF04/2005
623	c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
624	c Then that should not affect the ice particle radius
625	do ig=1,ngrid
626	if(pdpsrf(ig)ptimestep.gt.0.9(pplev(ig,1)-pplev(ig,2)))then
627	if(pdpsrf(ig)ptimestep.gt.0.9(pplev(ig,1)-pplev(ig,3)))
628	& rice(ig,2)=rice(ig,3)
629	rice(ig,1)=rice(ig,2)
630	end if
631	end do
632
633
634	DO l=1,nlay
635	DO ig=1,ngrid
636	rsedcloud(ig,l)=max(rice(ig,l)*
637	& (1.+nuice_sed)(1.+nuice_sed)(1.+nuice_sed),
638	& rdust(ig,l))
639	! rsedcloud(ig,l)=min(rsedcloud(ig,l),1.e-4)
640	ENDDO
641	ENDDO
642
643	! used for rad. transfer calculations
644	! nuice is constant because a lognormal distribution is prescribed
645	nuice(1:ngrid,1:nlay)=nuice_ref
646
647	c------Update tendencies for sub-grid water ice clouds
648	IF (CLFvarying) THEN
649	DO ig=1,ngrid
650	DO l=1,nlay
651	pdqcloud(ig,l,igcm_dust_mass)=pdqcloud(ig,l,igcm_dust_mass)
652	& *cloudfrac(ig,l)
653	pdqcloud(ig,l,igcm_ccn_mass)=pdqcloud(ig,l,igcm_ccn_mass)
654	& *cloudfrac(ig,l)
655	pdqcloud(ig,l,igcm_dust_number)=pdqcloud(ig,l,
656	& igcm_dust_number) *cloudfrac(ig,l)
657	pdqcloud(ig,l,igcm_ccn_number)=pdqcloud(ig,l,
658	& igcm_ccn_number) *cloudfrac(ig,l)
659	pdqcloud(ig,l,igcm_h2o_vap)=pdqcloud(ig,l,
660	& igcm_h2o_vap) *cloudfrac(ig,l)
661	pdqcloud(ig,l,igcm_h2o_ice)=pdqcloud(ig,l,
662	& igcm_h2o_ice) *cloudfrac(ig,l)
663	ENDDO
664	ENDDO
665	pdtcloud(:,:)=pdtcloud(:,:)*cloudfrac(:,:)
666	ENDIF
667	#ifndef MESOSCALE
668	c=======================================================================
669	call WRITEDIAGFI(ngrid,"pdqice2","pdqcloudice apres microphysique"
670	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,igcm_h2o_ice))
671	call WRITEDIAGFI(ngrid,"pdqvap2","pdqcloudvap apres microphysique"
672	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
673	& igcm_h2o_vap))
674	if (hdo) then
675	call WRITEDIAGFI(ngrid,"pdqiceD","pdqiceD apres microphysique"
676	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,igcm_hdo_ice))
677	call WRITEDIAGFI(ngrid,"pdqvapD","pdqvapD apres microphysique"
678	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
679	& igcm_hdo_vap))
680	endif
681	call WRITEDIAGFI(ngrid,"pdqccn2","pdqcloudccn apres microphysique"
682	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
683	& igcm_ccn_mass))
684	call WRITEDIAGFI(ngrid,"pdqccnN2","pdqcloudccnN apres "//
685	& "microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
686	& igcm_ccn_number))
687	call WRITEDIAGFI(ngrid,"pdqdust2", "pdqclouddust apres "//
688	& "microphysique","kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
689	& igcm_dust_mass))
690	call WRITEDIAGFI(ngrid,"pdqdustN2", "pdqclouddustN apres "//
691	& "microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
692	& igcm_dust_number))
693	c=======================================================================
694	#endif
695
696	END SUBROUTINE watercloud
697
698	subroutine ini_watercloud_mod(ngrid,nlayer,nq)
699	implicit none
700
701	integer,intent(in) :: ngrid ! number of atmospheric columns
702	integer,intent(in) :: nlayer ! number of atmospheric layers
703	integer,intent(in) :: nq ! number of tracers
704
705	allocate(zdqcloud(ngrid,nlayer,nq))
706	zdqcloud(:,:,:)=0
707	allocate(zdqscloud(ngrid,nq))
708	zdqscloud(:,:)=0
709
710	end subroutine ini_watercloud_mod
711
712
713	subroutine end_watercloud_mod
714	implicit none
715
716	if (allocated(zdqcloud)) deallocate(zdqcloud)
717	if (allocated(zdqscloud)) deallocate(zdqscloud)
718
719	end subroutine end_watercloud_mod
720
721	END MODULE watercloud_mod

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

Download in other formats:

Original Format