Context Navigation

watercloud_mod.F @ 2578

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

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

Download in other formats:

Original Format