Context Navigation

← Previous Revision
Next Revision →
Normal
Revision Log

watercloud_mod.F

Last change on this file was 2984, checked in by jnaar, 18 months ago
Adaptative timestep working and computed directly in improvedclouds. Simpleclouds working again. JN
File size: 27.4 KB

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

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

Download in other formats:

Original Format