Context Navigation

watercloud_mod.F @ 2800

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

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

Download in other formats:

Original Format