Context Navigation

watercloud_mod.F @ 2740

Last change on this file since 2740 was 2616, checked in by romain.vande, 3 years ago
LMDZ_MARS RV : Open_MP; Put all the "save" variables as "!$OMP THREADPRIVATE" in phymars. The code can now be tested, see README for more info
File size: 26.9 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
118	! Representation of sub-grid water ice clouds A. Pottier 2013
119	REAL :: ztclf(ngrid, nlay)
120	REAL :: zqclf(ngrid, nlay,nq)
121	REAL :: zdelt
122	REAL :: norm
123	REAL :: ponder
124	REAL :: tcond(ngrid,nlay)
125	REAL :: zqvap(ngrid,nlay)
126	REAL :: zqice(ngrid,nlay)
127	REAL :: spant ! delta T for the temperature distribution
128	! REAL :: zqsat(ngrid,nlay) ! saturation
129	REAL :: pteff(ngrid, nlay)! effective temperature in the cloud,neb
130	REAL :: pqeff(ngrid, nlay, nq)! effective tracers quantities in the cloud
131	REAL :: cloudfrac(ngrid,nlay) ! cloud fraction
132	REAL :: mincloud ! min cloud frac
133	LOGICAL, save :: flagcloud=.true.
134
135	!$OMP THREADPRIVATE(flagcloud)
136
137	c ** un petit test de coherence
138	c --------------------------
139
140	IF (firstcall) THEN
141
142	if (nq.gt.nqmx) then
143	write(,) 'stop in watercloud (nq.gt.nqmx)!'
144	write(,) 'nq=',nq,' nqmx=',nqmx
145	call abort_physic("watercloud","nq.gt.nqmx",1)
146	endif
147
148	write(,) "watercloud: igcm_h2o_vap=",igcm_h2o_vap
149	write(,) " igcm_h2o_ice=",igcm_h2o_ice
150
151	write(,) "time subsampling for microphysic ?"
152	#ifdef MESOSCALE
153	imicro = 2
154	#else
155	imicro = 30
156	#endif
157	call getin_p("imicro",imicro)
158	write(,)"watercloud: imicro = ",imicro
159
160	firstcall=.false.
161	ENDIF ! of IF (firstcall)
162
163	!! AS: moved out of firstcall to allow nesting+evoluting timestep
164	!! TBD: consider possible diff imicro with domains?
165	microtimestep = ptimestep/real(imicro)
166	if (microtimestep/=microtimestep_prev) then
167	! only tell the world if microtimestep has changed
168	write(,)"watercloud: Physical timestep is ",ptimestep
169	write(,)"watercloud: Microphysics timestep is ",microtimestep
170	microtimestep_prev=microtimestep
171	endif
172
173	c-----Initialization
174	sum_subpdq(1:ngrid,1:nlay,1:nq) = 0
175	sum_subpdt(1:ngrid,1:nlay) = 0
176
177	! default value if no ice
178	rhocloud(1:ngrid,1:nlay) = rho_dust
179
180	c-------------------------------------------------------------------
181	c 0. Representation of sub-grid water ice clouds
182	c------------------
183	c-----Initialization
184	pteff(1:ngrid,1:nlay) = 0
185	pqeff(1:ngrid,1:nlay,1:nq) = 0
186	DO l=1,nlay
187	DO ig=1,ngrid
188	pteff(ig,l)=pt(ig,l)
189	END DO
190	END DO
191	DO l=1,nlay
192	DO ig=1,ngrid
193	DO iq=1,nq
194	pqeff(ig,l,iq)=pq(ig,l,iq)
195	ENDDO
196	ENDDO
197	ENDDO
198	c-----Tendencies
199	DO l=1,nlay
200	DO ig=1,ngrid
201	ztclf(ig,l)=pt(ig,l)+ pdt(ig,l)*ptimestep
202	ENDDO
203	ENDDO
204	DO l=1,nlay
205	DO ig=1,ngrid
206	DO iq=1,nq
207	zqclf(ig,l,iq)=pq(ig,l,iq)+pdq(ig,l,iq)*ptimestep
208	ENDDO
209	ENDDO
210	ENDDO
211	c-----Effective temperature calculation
212	IF (CLFvarying) THEN
213	spant=3.0 ! delta T for the temprature distribution
214	mincloud=0.1 ! min cloudfrac when there is ice
215	IF (flagcloud) THEN
216	write(,) "Delta T", spant
217	write(,) "mincloud", mincloud
218	flagcloud=.false.
219	END IF
220	!CALL watersat(ngrid*nlay,ztclf,pplay,zqsat) !MV17: we dont need zqsat in the CLFvarying scheme
221	zqvap=zqclf(:,:,igcm_h2o_vap)
222	zqice=zqclf(:,:,igcm_h2o_ice)
223	CALL tcondwater(ngrid*nlay,pplay,zqvap+zqice,tcond)
224	DO l=1,nlay
225	DO ig=1,ngrid
226	zdelt=spant !MAX(spant*ztclf(ig,l),1.e-12), now totally in K. Fixed width
227	IF (tcond(ig,l) .ge. (ztclf(ig,l)+zdelt)) THEN
228	pteff(ig,l)=ztclf(ig,l)
229	cloudfrac(ig,l)=1.
230	ELSE IF (tcond(ig,l) .le. (ztclf(ig,l)-zdelt)) THEN
231	pteff(ig,l)=ztclf(ig,l)-zdelt
232	cloudfrac(ig,l)=mincloud
233	ELSE
234	cloudfrac(ig,l)=(tcond(ig,l)-ztclf(ig,l)+zdelt)/
235	& (2.0*zdelt)
236	pteff(ig,l)=(tcond(ig,l)+ztclf(ig,l)-zdelt)/2.
237	END IF
238	pteff(ig,l)=pteff(ig,l)-pdt(ig,l)*ptimestep
239	IF (cloudfrac(ig,l).le.mincloud) THEN !MV17: replaced .le.0 by .le.mincloud
240	cloudfrac(ig,l)=mincloud
241	ELSE IF (cloudfrac(ig,l).gt.1) THEN
242	cloudfrac(ig,l)=1.
243	END IF
244	ENDDO
245	ENDDO
246	c-----Calculation of the total cloud coverage of the column
247	DO ig=1,ngrid
248	totcloudfrac(ig) = 0.
249	norm=0.
250	DO l=1,nlay
251	ponder=zqice(ig,l)*(pplev(ig,l) - pplev(ig,l+1))
252	totcloudfrac(ig) = totcloudfrac(ig)
253	& + cloudfrac(ig,l)*ponder
254	norm=norm+ponder
255	ENDDO
256	totcloudfrac(ig)=MAX(totcloudfrac(ig)/norm,1.e-12) ! min value if NaNs
257	ENDDO
258	c-----Effective tracers quantities in the cloud fraction
259	IF (microphys) THEN
260	pqeff(:,:,igcm_ccn_mass)=pq(:,:,igcm_ccn_mass)/
261	& cloudfrac(:,:)
262	pqeff(:,:,igcm_ccn_number)=pq(:,:,igcm_ccn_number)/
263	& cloudfrac(:,:)
264	END IF ! end if (microphys)
265	pqeff(:,:,igcm_h2o_ice)=pq(:,:,igcm_h2o_ice)/
266	& cloudfrac(:,:)
267	!! CLFvarying outputs
268	CALL WRITEDIAGFI(ngrid,'pqeffice','pqeffice',
269	& 'kg/kg',3,pqeff(:,:,igcm_h2o_ice))
270	CALL WRITEDIAGFI(ngrid,'pteff','pteff',
271	& 'K',3,pteff(:,:))
272	CALL WRITEDIAGFI(ngrid,'tcond','tcond',
273	& 'K',3,tcond(:,:))
274	CALL WRITEDIAGFI(ngrid,'cloudfrac','cloudfrac',
275	& 'K',3,cloudfrac(:,:))
276	END IF ! end if (CLFvarying)
277	c------------------------------------------------------------------
278	c Time subsampling for microphysics
279	c------------------------------------------------------------------
280	rhocloud(1:ngrid,1:nlay) = rho_dust
281	DO microstep=1,imicro
282
283	c-------------------------------------------------------------------
284	c 1. Tendencies:
285	c------------------
286
287
288	c------ Temperature tendency subpdt
289	! Each microtimestep we give the cloud scheme a stepped entry subpdt instead of pdt
290	! If imicro=1 subpdt is the same as pdt
291	DO l=1,nlay
292	DO ig=1,ngrid
293	sum_subpdt(ig,l) = sum_subpdt(ig,l)
294	& + pdt(ig,l) ! At each micro timestep we add pdt in order to have a stepped entry
295	ENDDO
296	ENDDO
297	c------ Tracers tendencies subpdq are additionned
298	c------ At each micro timestep we add pdq in order to have a stepped entry
299	IF (microphys) THEN
300	DO l=1,nlay
301	DO ig=1,ngrid
302	sum_subpdq(ig,l,igcm_dust_mass) =
303	& sum_subpdq(ig,l,igcm_dust_mass)
304	& + pdq(ig,l,igcm_dust_mass)
305	sum_subpdq(ig,l,igcm_dust_number) =
306	& sum_subpdq(ig,l,igcm_dust_number)
307	& + pdq(ig,l,igcm_dust_number)
308	sum_subpdq(ig,l,igcm_ccn_mass) =
309	& sum_subpdq(ig,l,igcm_ccn_mass)
310	& + pdq(ig,l,igcm_ccn_mass)
311	sum_subpdq(ig,l,igcm_ccn_number) =
312	& sum_subpdq(ig,l,igcm_ccn_number)
313	& + pdq(ig,l,igcm_ccn_number)
314	ENDDO
315	ENDDO
316	ENDIF
317	DO l=1,nlay
318	DO ig=1,ngrid
319	sum_subpdq(ig,l,igcm_h2o_ice) =
320	& sum_subpdq(ig,l,igcm_h2o_ice)
321	& + pdq(ig,l,igcm_h2o_ice)
322	sum_subpdq(ig,l,igcm_h2o_vap) =
323	& sum_subpdq(ig,l,igcm_h2o_vap)
324	& + pdq(ig,l,igcm_h2o_vap)
325	IF (hdo) THEN
326	sum_subpdq(ig,l,igcm_hdo_ice) =
327	& sum_subpdq(ig,l,igcm_hdo_ice)
328	& + pdq(ig,l,igcm_hdo_ice)
329	sum_subpdq(ig,l,igcm_hdo_vap) =
330	& sum_subpdq(ig,l,igcm_hdo_vap)
331	& + pdq(ig,l,igcm_hdo_vap)
332	ENDIF !hdo
333	ENDDO
334	ENDDO
335
336	c-------------------------------------------------------------------
337	c 2. Main call to the different cloud schemes:
338	c------------------------------------------------
339	IF (microphys) THEN
340	CALL improvedclouds(ngrid,nlay,microtimestep,
341	& pplay,pteff,sum_subpdt,
342	& pqeff,sum_subpdq,subpdqcloud,subpdtcloud,
343	& nq,tauscaling)
344
345	ELSE
346	CALL simpleclouds(ngrid,nlay,microtimestep,
347	& pplay,pzlay,pteff,sum_subpdt,
348	& pqeff,sum_subpdq,subpdqcloud,subpdtcloud,
349	& nq,tau,rice)
350	ENDIF
351
352	c-------------------------------------------------------------------
353	c 3. Updating tendencies after cloud scheme:
354	c-----------------------------------------------
355
356	IF (microphys) THEN
357	DO l=1,nlay
358	DO ig=1,ngrid
359	sum_subpdq(ig,l,igcm_dust_mass) =
360	& sum_subpdq(ig,l,igcm_dust_mass)
361	& + subpdqcloud(ig,l,igcm_dust_mass)
362	sum_subpdq(ig,l,igcm_dust_number) =
363	& sum_subpdq(ig,l,igcm_dust_number)
364	& + subpdqcloud(ig,l,igcm_dust_number)
365	sum_subpdq(ig,l,igcm_ccn_mass) =
366	& sum_subpdq(ig,l,igcm_ccn_mass)
367	& + subpdqcloud(ig,l,igcm_ccn_mass)
368	sum_subpdq(ig,l,igcm_ccn_number) =
369	& sum_subpdq(ig,l,igcm_ccn_number)
370	& + subpdqcloud(ig,l,igcm_ccn_number)
371	ENDDO
372	ENDDO
373	ENDIF
374	DO l=1,nlay
375	DO ig=1,ngrid
376	sum_subpdq(ig,l,igcm_h2o_ice) =
377	& sum_subpdq(ig,l,igcm_h2o_ice)
378	& + subpdqcloud(ig,l,igcm_h2o_ice)
379	sum_subpdq(ig,l,igcm_h2o_vap) =
380	& sum_subpdq(ig,l,igcm_h2o_vap)
381	& + subpdqcloud(ig,l,igcm_h2o_vap)
382
383	IF (hdo) THEN
384	sum_subpdq(ig,l,igcm_hdo_ice) =
385	& sum_subpdq(ig,l,igcm_hdo_ice)
386	& + subpdqcloud(ig,l,igcm_hdo_ice)
387	sum_subpdq(ig,l,igcm_hdo_vap) =
388	& sum_subpdq(ig,l,igcm_hdo_vap)
389	& + subpdqcloud(ig,l,igcm_hdo_vap)
390	ENDIF ! hdo
391
392	ENDDO
393	ENDDO
394
395	IF (activice) THEN
396	DO l=1,nlay
397	DO ig=1,ngrid
398	sum_subpdt(ig,l) =
399	& sum_subpdt(ig,l) + subpdtcloud(ig,l)
400	ENDDO
401	ENDDO
402	ENDIF
403
404	! !! Example of how to use writediagmicrofi useful to
405	! !! get outputs at each microphysical sub-timestep (better to be used in 1D)
406	! CALL WRITEDIAGMICROFI(ngrid,imicro,microstep,
407	! & microtimestep,'subpdtcloud',
408	! & 'subpdtcloud','K/s',1,subpdtcloud(:,:))
409
410	ENDDO ! of DO microstep=1,imicro
411
412	c-------------------------------------------------------------------
413	c 6. Compute final tendencies after time loop:
414	c------------------------------------------------
415
416	c------ Temperature tendency pdtcloud
417	DO l=1,nlay
418	DO ig=1,ngrid
419	pdtcloud(ig,l) =
420	& sum_subpdt(ig,l)/real(imicro)-pdt(ig,l)
421	ENDDO
422	ENDDO
423
424	c------ Tracers tendencies pdqcloud
425	DO l=1,nlay
426	DO ig=1,ngrid
427	pdqcloud(ig,l,igcm_h2o_ice) =
428	& sum_subpdq(ig,l,igcm_h2o_ice)/real(imicro)
429	& - pdq(ig,l,igcm_h2o_ice)
430	pdqcloud(ig,l,igcm_h2o_vap) =
431	& sum_subpdq(ig,l,igcm_h2o_vap)/real(imicro)
432	& - pdq(ig,l,igcm_h2o_vap)
433	IF (hdo) THEN
434	pdqcloud(ig,l,igcm_hdo_ice) =
435	& sum_subpdq(ig,l,igcm_hdo_ice)/real(imicro)
436	& - pdq(ig,l,igcm_hdo_ice)
437	pdqcloud(ig,l,igcm_hdo_vap) =
438	& sum_subpdq(ig,l,igcm_hdo_vap)/real(imicro)
439	& - pdq(ig,l,igcm_hdo_vap)
440	ENDIF !hdo
441	ENDDO
442	ENDDO
443
444	IF(microphys) THEN
445	DO l=1,nlay
446	DO ig=1,ngrid
447	pdqcloud(ig,l,igcm_ccn_mass) =
448	& sum_subpdq(ig,l,igcm_ccn_mass)/real(imicro)
449	& - pdq(ig,l,igcm_ccn_mass)
450	pdqcloud(ig,l,igcm_ccn_number) =
451	& sum_subpdq(ig,l,igcm_ccn_number)/real(imicro)
452	& - pdq(ig,l,igcm_ccn_number)
453	ENDDO
454	ENDDO
455	ENDIF
456
457	IF(scavenging) THEN
458	DO l=1,nlay
459	DO ig=1,ngrid
460	pdqcloud(ig,l,igcm_dust_mass) =
461	& sum_subpdq(ig,l,igcm_dust_mass)/real(imicro)
462	& - pdq(ig,l,igcm_dust_mass)
463	pdqcloud(ig,l,igcm_dust_number) =
464	& sum_subpdq(ig,l,igcm_dust_number)/real(imicro)
465	& - pdq(ig,l,igcm_dust_number)
466	ENDDO
467	ENDDO
468	ENDIF
469
470	c------- Due to stepped entry, other processes tendencies can add up to negative values
471	c------- Therefore, enforce positive values and conserve mass
472	IF(microphys) THEN
473	DO l=1,nlay
474	DO ig=1,ngrid
475	IF ((pq(ig,l,igcm_ccn_number) +
476	& ptimestep* (pdq(ig,l,igcm_ccn_number) +
477	& pdqcloud(ig,l,igcm_ccn_number)) .le. 1.)
478	& .or. (pq(ig,l,igcm_ccn_mass) +
479	& ptimestep* (pdq(ig,l,igcm_ccn_mass) +
480	& pdqcloud(ig,l,igcm_ccn_mass)) .le. 1.e-20)) THEN
481	pdqcloud(ig,l,igcm_ccn_number) =
482	& - pq(ig,l,igcm_ccn_number)/ptimestep
483	& - pdq(ig,l,igcm_ccn_number) + 1.
484	pdqcloud(ig,l,igcm_dust_number) =
485	& -pdqcloud(ig,l,igcm_ccn_number)
486	pdqcloud(ig,l,igcm_ccn_mass) =
487	& - pq(ig,l,igcm_ccn_mass)/ptimestep
488	& - pdq(ig,l,igcm_ccn_mass) + 1.e-20
489	pdqcloud(ig,l,igcm_dust_mass) =
490	& -pdqcloud(ig,l,igcm_ccn_mass)
491	ENDIF
492	ENDDO
493	ENDDO
494	ENDIF
495
496	IF(scavenging) THEN
497	DO l=1,nlay
498	DO ig=1,ngrid
499	IF ((pq(ig,l,igcm_dust_number) +
500	& ptimestep* (pdq(ig,l,igcm_dust_number) +
501	& pdqcloud(ig,l,igcm_dust_number)) .le. 1.)
502	& .or. (pq(ig,l,igcm_dust_mass) +
503	& ptimestep* (pdq(ig,l,igcm_dust_mass) +
504	& pdqcloud(ig,l,igcm_dust_mass)) .le. 1.e-20)) THEN
505	pdqcloud(ig,l,igcm_dust_number) =
506	& - pq(ig,l,igcm_dust_number)/ptimestep
507	& - pdq(ig,l,igcm_dust_number) + 1.
508	pdqcloud(ig,l,igcm_ccn_number) =
509	& -pdqcloud(ig,l,igcm_dust_number)
510	pdqcloud(ig,l,igcm_dust_mass) =
511	& - pq(ig,l,igcm_dust_mass)/ptimestep
512	& - pdq(ig,l,igcm_dust_mass) + 1.e-20
513	pdqcloud(ig,l,igcm_ccn_mass) =
514	& -pdqcloud(ig,l,igcm_dust_mass)
515	ENDIF
516	ENDDO
517	ENDDO
518	ENDIF
519
520	DO l=1,nlay
521	DO ig=1,ngrid
522
523	IF (pq(ig,l,igcm_h2o_ice) + ptimestep*
524	& (pdq(ig,l,igcm_h2o_ice) + pdqcloud(ig,l,igcm_h2o_ice))
525	& .le. qperemin) THEN
526	pdqcloud(ig,l,igcm_h2o_ice) =
527	& - pq(ig,l,igcm_h2o_ice)/ptimestep - pdq(ig,l,igcm_h2o_ice)
528	pdqcloud(ig,l,igcm_h2o_vap) = -pdqcloud(ig,l,igcm_h2o_ice)
529	if (hdo) then
530	pdqcloud(ig,l,igcm_hdo_ice) =
531	& - pq(ig,l,igcm_hdo_ice)/ptimestep - pdq(ig,l,igcm_hdo_ice)
532	pdqcloud(ig,l,igcm_hdo_vap) = -pdqcloud(ig,l,igcm_hdo_ice)
533	endif
534	ENDIF
535	IF (pq(ig,l,igcm_h2o_vap) + ptimestep*
536	& (pdq(ig,l,igcm_h2o_vap) + pdqcloud(ig,l,igcm_h2o_vap))
537	& .le. qperemin) THEN
538	pdqcloud(ig,l,igcm_h2o_vap) =
539	& - pq(ig,l,igcm_h2o_vap)/ptimestep - pdq(ig,l,igcm_h2o_vap)
540	pdqcloud(ig,l,igcm_h2o_ice) = -pdqcloud(ig,l,igcm_h2o_vap)
541	if (hdo) then
542	pdqcloud(ig,l,igcm_hdo_vap) =
543	& - pq(ig,l,igcm_hdo_vap)/ptimestep - pdq(ig,l,igcm_hdo_vap)
544	pdqcloud(ig,l,igcm_hdo_ice) = -pdqcloud(ig,l,igcm_hdo_vap)
545	endif
546	ENDIF
547
548	ENDDO
549	ENDDO
550
551	c------Update the ice and dust particle size "rice" for output or photochemistry
552	c------Only rsedcloud is used for the water cycle
553
554	IF(scavenging) THEN
555	DO l=1, nlay
556	DO ig=1,ngrid
557
558	call updaterdust(
559	& pq(ig,l,igcm_dust_mass) + ! dust mass
560	& (pdq(ig,l,igcm_dust_mass) + ! dust mass
561	& pdqcloud(ig,l,igcm_dust_mass))*ptimestep, ! dust mass
562	& pq(ig,l,igcm_dust_number) + ! dust number
563	& (pdq(ig,l,igcm_dust_number) + ! dust number
564	& pdqcloud(ig,l,igcm_dust_number))*ptimestep, ! dust number
565	& rdust(ig,l))
566
567	ENDDO
568	ENDDO
569	ENDIF
570
571	IF(microphys) THEN
572
573	! In case one does not want to allow supersatured water when using microphysics.
574	! Not done by default.
575	IF(.not.supersat) THEN
576	zt = pt + (pdt+pdtcloud)*ptimestep
577	call watersat(ngrid*nlay,zt,pplay,zqsat)
578	DO l=1, nlay
579	DO ig=1,ngrid
580	IF (pq(ig,l,igcm_h2o_vap)
581	& + (pdq(ig,l,igcm_h2o_vap) + pdqcloud(ig,l,igcm_h2o_vap))
582	& * ptimestep .ge. zqsat(ig,l)) THEN
583	pdqcloud(ig,l,igcm_h2o_vap) =
584	& (zqsat(ig,l) - pq(ig,l,igcm_h2o_vap))/ptimestep
585	& - pdq(ig,l,igcm_h2o_vap)
586	pdqcloud(ig,l,igcm_h2o_ice) =
587	& -pdqcloud(ig,l,igcm_h2o_vap)
588	! no need to correct ccn_number, updaterad can handle this properly.
589	ENDIF
590	ENDDO
591	ENDDO
592	ENDIF
593
594	DO l=1, nlay
595	DO ig=1,ngrid
596
597	call updaterice_micro(
598	& pq(ig,l,igcm_h2o_ice) + ! ice mass
599	& (pdq(ig,l,igcm_h2o_ice) + ! ice mass
600	& pdqcloud(ig,l,igcm_h2o_ice))*ptimestep, ! ice mass
601	& pq(ig,l,igcm_ccn_mass) + ! ccn mass
602	& (pdq(ig,l,igcm_ccn_mass) + ! ccn mass
603	& pdqcloud(ig,l,igcm_ccn_mass))*ptimestep, ! ccn mass
604	& pq(ig,l,igcm_ccn_number) + ! ccn number
605	& (pdq(ig,l,igcm_ccn_number) + ! ccn number
606	& pdqcloud(ig,l,igcm_ccn_number))*ptimestep, ! ccn number
607	& tauscaling(ig),rice(ig,l),rhocloud(ig,l))
608
609	ENDDO
610	ENDDO
611
612	ELSE ! no microphys
613
614	DO l=1,nlay
615	DO ig=1,ngrid
616
617	call updaterice_typ(
618	& pq(ig,l,igcm_h2o_ice) + ! ice mass
619	& (pdq(ig,l,igcm_h2o_ice) + ! ice mass
620	& pdqcloud(ig,l,igcm_h2o_ice))*ptimestep, ! ice mass
621	& tau(ig,1),pzlay(ig,l),rice(ig,l))
622
623	ENDDO
624	ENDDO
625
626	ENDIF ! of IF(microphys)
627
628
629
630	c A correction if a lot of subliming CO2 fills the 1st layer FF04/2005
631	c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
632	c Then that should not affect the ice particle radius
633	do ig=1,ngrid
634	if(pdpsrf(ig)ptimestep.gt.0.9(pplev(ig,1)-pplev(ig,2)))then
635	if(pdpsrf(ig)ptimestep.gt.0.9(pplev(ig,1)-pplev(ig,3)))
636	& rice(ig,2)=rice(ig,3)
637	rice(ig,1)=rice(ig,2)
638	end if
639	end do
640
641
642	DO l=1,nlay
643	DO ig=1,ngrid
644	rsedcloud(ig,l)=max(rice(ig,l)*
645	& (1.+nuice_sed)(1.+nuice_sed)(1.+nuice_sed),
646	& rdust(ig,l))
647	! rsedcloud(ig,l)=min(rsedcloud(ig,l),1.e-4)
648	ENDDO
649	ENDDO
650
651	! used for rad. transfer calculations
652	! nuice is constant because a lognormal distribution is prescribed
653	nuice(1:ngrid,1:nlay)=nuice_ref
654
655	c------Update tendencies for sub-grid water ice clouds
656	IF (CLFvarying) THEN
657	DO ig=1,ngrid
658	DO l=1,nlay
659	pdqcloud(ig,l,igcm_dust_mass)=pdqcloud(ig,l,igcm_dust_mass)
660	& *cloudfrac(ig,l)
661	pdqcloud(ig,l,igcm_ccn_mass)=pdqcloud(ig,l,igcm_ccn_mass)
662	& *cloudfrac(ig,l)
663	pdqcloud(ig,l,igcm_dust_number)=pdqcloud(ig,l,
664	& igcm_dust_number) *cloudfrac(ig,l)
665	pdqcloud(ig,l,igcm_ccn_number)=pdqcloud(ig,l,
666	& igcm_ccn_number) *cloudfrac(ig,l)
667	pdqcloud(ig,l,igcm_h2o_vap)=pdqcloud(ig,l,
668	& igcm_h2o_vap) *cloudfrac(ig,l)
669	pdqcloud(ig,l,igcm_h2o_ice)=pdqcloud(ig,l,
670	& igcm_h2o_ice) *cloudfrac(ig,l)
671	ENDDO
672	ENDDO
673	pdtcloud(:,:)=pdtcloud(:,:)*cloudfrac(:,:)
674	ENDIF
675	#ifndef MESOSCALE
676	c=======================================================================
677	call WRITEDIAGFI(ngrid,"pdqice2","pdqcloudice apres microphysique"
678	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,igcm_h2o_ice))
679	call WRITEDIAGFI(ngrid,"pdqvap2","pdqcloudvap apres microphysique"
680	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
681	& igcm_h2o_vap))
682	if (hdo) then
683	call WRITEDIAGFI(ngrid,"pdqiceD","pdqiceD apres microphysique"
684	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,igcm_hdo_ice))
685	call WRITEDIAGFI(ngrid,"pdqvapD","pdqvapD apres microphysique"
686	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
687	& igcm_hdo_vap))
688	endif
689	call WRITEDIAGFI(ngrid,"pdqccn2","pdqcloudccn apres microphysique"
690	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
691	& igcm_ccn_mass))
692	call WRITEDIAGFI(ngrid,"pdqccnN2","pdqcloudccnN apres "//
693	& "microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
694	& igcm_ccn_number))
695	call WRITEDIAGFI(ngrid,"pdqdust2", "pdqclouddust apres "//
696	& "microphysique","kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
697	& igcm_dust_mass))
698	call WRITEDIAGFI(ngrid,"pdqdustN2", "pdqclouddustN apres "//
699	& "microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
700	& igcm_dust_number))
701	c=======================================================================
702	#endif
703
704	END SUBROUTINE watercloud
705
706	subroutine ini_watercloud_mod(ngrid,nlayer,nq)
707	implicit none
708
709	integer,intent(in) :: ngrid ! number of atmospheric columns
710	integer,intent(in) :: nlayer ! number of atmospheric layers
711	integer,intent(in) :: nq ! number of tracers
712
713	allocate(zdqcloud(ngrid,nlayer,nq))
714	zdqcloud(:,:,:)=0
715	allocate(zdqscloud(ngrid,nq))
716	zdqscloud(:,:)=0
717
718	end subroutine ini_watercloud_mod
719
720
721	subroutine end_watercloud_mod
722	implicit none
723
724	if (allocated(zdqcloud)) deallocate(zdqcloud)
725	if (allocated(zdqscloud)) deallocate(zdqscloud)
726
727	end subroutine end_watercloud_mod
728
729	END MODULE watercloud_mod

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

Download in other formats:

Original Format