Context Navigation

watercloud_mod.F @ 2311

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

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

Download in other formats:

Original Format