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

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
[38]	29	IMPLICIT NONE
	30
[633]	31
[38]	32	c=======================================================================
[358]	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]	38	c
[633]	39	c There is a time loop specific to cloud formation
	40	c due to timescales smaller than the GCM integration timestep.
	41	c
[358]	42	c Authors: Franck Montmessin, Francois Forget, Ehouarn Millour,
[522]	43	c J.-B. Madeleine, Thomas Navarro
[38]	44	c
[633]	45	c 2004 - 2012
[38]	46	c=======================================================================
	47
	48	c-----------------------------------------------------------------------
	49	c declarations:
	50	c -------------
	51
[1964]	52	include "callkeys.h"
[38]	53
[1976]	54	c Inputs/outputs:
[38]	55	c ------
	56
[1976]	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.
[38]	66
[1976]	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)
[38]	69
[1976]	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)
[38]	73
[1976]	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
[633]	76	! a la chaleur latente
[1976]	77	REAL, INTENT(INOUT) :: rice(ngrid,nlay) ! Ice mass mean radius (m)
[38]	78	! (r_c in montmessin_2004)
[1976]	79	REAL, INTENT(OUT) :: nuice(ngrid,nlay) ! Estimated effective variance
[38]	80	! of the size distribution
[1976]	81	REAL, INTENT(OUT) :: rsedcloud(ngrid,nlay) ! Cloud sedimentation radius
	82	REAL, INTENT(OUT) :: rhocloud(ngrid,nlay) ! Cloud density (kg.m-3)
[38]	83
[1711]	84	REAL, INTENT(INOUT):: totcloudfrac(ngrid) ! Cloud fraction (A. Pottier 2013)
[1976]	85
	86	c Locals:
[38]	87	c ------
[1976]	88
[633]	89	! for ice radius computation
	90	REAL Mo,No
	91	REAl ccntyp
	92
	93	! for time loop
	94	INTEGER microstep ! time subsampling step variable
[1969]	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
[633]	98
[2616]	99	!$OMP THREADPRIVATE(imicro,microtimestep)
	100	!$OMP THREADPRIVATE(microtimestep_prev)
[633]	101	! tendency given by clouds (inside the micro loop)
	102	REAL subpdqcloud(ngrid,nlay,nq) ! cf. pdqcloud
	103	REAL subpdtcloud(ngrid,nlay) ! cf. pdtcloud
[38]	104
[633]	105	! global tendency (clouds+physics)
[1909]	106	REAL sum_subpdq(ngrid,nlay,nq) ! cf. pdqcloud
	107	REAL sum_subpdt(ngrid,nlay) ! cf. pdtcloud
[633]	108
[1467]	109	! no supersaturation when option supersat is false
	110	REAL zt(ngrid,nlay) ! local value of temperature
	111	REAL zqsat(ngrid,nlay) ! saturation
[633]	112
	113	INTEGER iq,ig,l
[38]	114	LOGICAL,SAVE :: firstcall=.true.
[2616]	115
	116	!$OMP THREADPRIVATE(firstcall)
[38]	117
[1711]	118	! Representation of sub-grid water ice clouds A. Pottier 2013
[1880]	119	REAL :: ztclf(ngrid, nlay)
	120	REAL :: zqclf(ngrid, nlay,nq)
[1711]	121	REAL :: zdelt
	122	REAL :: norm
	123	REAL :: ponder
	124	REAL :: tcond(ngrid,nlay)
[1880]	125	REAL :: zqvap(ngrid,nlay)
[1711]	126	REAL :: zqice(ngrid,nlay)
	127	REAL :: spant ! delta T for the temperature distribution
	128	! REAL :: zqsat(ngrid,nlay) ! saturation
[1909]	129	REAL :: pteff(ngrid, nlay)! effective temperature in the cloud,neb
[1711]	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.
[2616]	134
	135	!$OMP THREADPRIVATE(flagcloud)
[2312]	136
[38]	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
[2311]	145	call abort_physic("watercloud","nq.gt.nqmx",1)
[38]	146	endif
	147
[358]	148	write(,) "watercloud: igcm_h2o_vap=",igcm_h2o_vap
	149	write(,) " igcm_h2o_ice=",igcm_h2o_ice
[633]	150
	151	write(,) "time subsampling for microphysic ?"
	152	#ifdef MESOSCALE
	153	imicro = 2
	154	#else
[951]	155	imicro = 30
[633]	156	#endif
[2311]	157	call getin_p("imicro",imicro)
[1969]	158	write(,)"watercloud: imicro = ",imicro
[633]	159
[38]	160	firstcall=.false.
	161	ENDIF ! of IF (firstcall)
[1774]	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)
[1969]	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
[522]	172
[633]	173	c-----Initialization
[1909]	174	sum_subpdq(1:ngrid,1:nlay,1:nq) = 0
	175	sum_subpdt(1:ngrid,1:nlay) = 0
[633]	176
	177	! default value if no ice
	178	rhocloud(1:ngrid,1:nlay) = rho_dust
[38]	179
[1711]	180	c-------------------------------------------------------------------
	181	c 0. Representation of sub-grid water ice clouds
	182	c------------------
[1880]	183	c-----Initialization
[1909]	184	pteff(1:ngrid,1:nlay) = 0
[1880]	185	pqeff(1:ngrid,1:nlay,1:nq) = 0
	186	DO l=1,nlay
	187	DO ig=1,ngrid
[1909]	188	pteff(ig,l)=pt(ig,l)
[1880]	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
[1711]	198	c-----Tendencies
	199	DO l=1,nlay
[1880]	200	DO ig=1,ngrid
	201	ztclf(ig,l)=pt(ig,l)+ pdt(ig,l)*ptimestep
	202	ENDDO
[1711]	203	ENDDO
	204	DO l=1,nlay
	205	DO ig=1,ngrid
	206	DO iq=1,nq
[1880]	207	zqclf(ig,l,iq)=pq(ig,l,iq)+pdq(ig,l,iq)*ptimestep
[1711]	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
[1880]	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)
[1711]	223	CALL tcondwater(ngrid*nlay,pplay,zqvap+zqice,tcond)
	224	DO l=1,nlay
	225	DO ig=1,ngrid
[1880]	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
[1909]	228	pteff(ig,l)=ztclf(ig,l)
[1711]	229	cloudfrac(ig,l)=1.
[1880]	230	ELSE IF (tcond(ig,l) .le. (ztclf(ig,l)-zdelt)) THEN
[1909]	231	pteff(ig,l)=ztclf(ig,l)-zdelt
[1711]	232	cloudfrac(ig,l)=mincloud
	233	ELSE
[1880]	234	cloudfrac(ig,l)=(tcond(ig,l)-ztclf(ig,l)+zdelt)/
[1711]	235	& (2.0*zdelt)
[1909]	236	pteff(ig,l)=(tcond(ig,l)+ztclf(ig,l)-zdelt)/2.
[1711]	237	END IF
[1909]	238	pteff(ig,l)=pteff(ig,l)-pdt(ig,l)*ptimestep
[1880]	239	IF (cloudfrac(ig,l).le.mincloud) THEN !MV17: replaced .le.0 by .le.mincloud
[1711]	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
[1880]	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(:,:)
[1909]	267	!! CLFvarying outputs
[1880]	268	CALL WRITEDIAGFI(ngrid,'pqeffice','pqeffice',
	269	& 'kg/kg',3,pqeff(:,:,igcm_h2o_ice))
[1909]	270	CALL WRITEDIAGFI(ngrid,'pteff','pteff',
	271	& 'K',3,pteff(:,:))
[1880]	272	CALL WRITEDIAGFI(ngrid,'tcond','tcond',
	273	& 'K',3,tcond(:,:))
	274	CALL WRITEDIAGFI(ngrid,'cloudfrac','cloudfrac',
	275	& 'K',3,cloudfrac(:,:))
[1909]	276	END IF ! end if (CLFvarying)
[633]	277	c------------------------------------------------------------------
[1880]	278	c Time subsampling for microphysics
	279	c------------------------------------------------------------------
[1711]	280	rhocloud(1:ngrid,1:nlay) = rho_dust
[633]	281	DO microstep=1,imicro
[522]	282
[633]	283	c-------------------------------------------------------------------
	284	c 1. Tendencies:
	285	c------------------
[38]	286
[633]	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
[1909]	293	sum_subpdt(ig,l) = sum_subpdt(ig,l)
[633]	294	& + pdt(ig,l) ! At each micro timestep we add pdt in order to have a stepped entry
	295	ENDDO
	296	ENDDO
[1909]	297	c------ Tracers tendencies subpdq are additionned
[633]	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
[1909]	302	sum_subpdq(ig,l,igcm_dust_mass) =
	303	& sum_subpdq(ig,l,igcm_dust_mass)
[633]	304	& + pdq(ig,l,igcm_dust_mass)
[1909]	305	sum_subpdq(ig,l,igcm_dust_number) =
	306	& sum_subpdq(ig,l,igcm_dust_number)
[633]	307	& + pdq(ig,l,igcm_dust_number)
[1909]	308	sum_subpdq(ig,l,igcm_ccn_mass) =
	309	& sum_subpdq(ig,l,igcm_ccn_mass)
[633]	310	& + pdq(ig,l,igcm_ccn_mass)
[1909]	311	sum_subpdq(ig,l,igcm_ccn_number) =
	312	& sum_subpdq(ig,l,igcm_ccn_number)
[633]	313	& + pdq(ig,l,igcm_ccn_number)
	314	ENDDO
	315	ENDDO
	316	ENDIF
	317	DO l=1,nlay
	318	DO ig=1,ngrid
[1909]	319	sum_subpdq(ig,l,igcm_h2o_ice) =
	320	& sum_subpdq(ig,l,igcm_h2o_ice)
[633]	321	& + pdq(ig,l,igcm_h2o_ice)
[1909]	322	sum_subpdq(ig,l,igcm_h2o_vap) =
	323	& sum_subpdq(ig,l,igcm_h2o_vap)
[633]	324	& + pdq(ig,l,igcm_h2o_vap)
[2312]	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
[633]	333	ENDDO
[1880]	334	ENDDO
[633]	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,
[1909]	341	& pplay,pteff,sum_subpdt,
	342	& pqeff,sum_subpdq,subpdqcloud,subpdtcloud,
[633]	343	& nq,tauscaling)
	344
	345	ELSE
	346	CALL simpleclouds(ngrid,nlay,microtimestep,
[1909]	347	& pplay,pzlay,pteff,sum_subpdt,
	348	& pqeff,sum_subpdq,subpdqcloud,subpdtcloud,
[645]	349	& nq,tau,rice)
[633]	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
[1909]	359	sum_subpdq(ig,l,igcm_dust_mass) =
	360	& sum_subpdq(ig,l,igcm_dust_mass)
[633]	361	& + subpdqcloud(ig,l,igcm_dust_mass)
[1909]	362	sum_subpdq(ig,l,igcm_dust_number) =
	363	& sum_subpdq(ig,l,igcm_dust_number)
[633]	364	& + subpdqcloud(ig,l,igcm_dust_number)
[1909]	365	sum_subpdq(ig,l,igcm_ccn_mass) =
	366	& sum_subpdq(ig,l,igcm_ccn_mass)
[633]	367	& + subpdqcloud(ig,l,igcm_ccn_mass)
[1909]	368	sum_subpdq(ig,l,igcm_ccn_number) =
	369	& sum_subpdq(ig,l,igcm_ccn_number)
[633]	370	& + subpdqcloud(ig,l,igcm_ccn_number)
	371	ENDDO
	372	ENDDO
	373	ENDIF
	374	DO l=1,nlay
	375	DO ig=1,ngrid
[1909]	376	sum_subpdq(ig,l,igcm_h2o_ice) =
	377	& sum_subpdq(ig,l,igcm_h2o_ice)
[633]	378	& + subpdqcloud(ig,l,igcm_h2o_ice)
[1909]	379	sum_subpdq(ig,l,igcm_h2o_vap) =
	380	& sum_subpdq(ig,l,igcm_h2o_vap)
[633]	381	& + subpdqcloud(ig,l,igcm_h2o_vap)
[2312]	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
[633]	392	ENDDO
	393	ENDDO
[882]	394
	395	IF (activice) THEN
	396	DO l=1,nlay
	397	DO ig=1,ngrid
[1909]	398	sum_subpdt(ig,l) =
	399	& sum_subpdt(ig,l) + subpdtcloud(ig,l)
[882]	400	ENDDO
	401	ENDDO
	402	ENDIF
[2437]	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(:,:))
[633]	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) =
[1909]	420	& sum_subpdt(ig,l)/real(imicro)-pdt(ig,l)
[633]	421	ENDDO
	422	ENDDO
[740]	423
[633]	424	c------ Tracers tendencies pdqcloud
[703]	425	DO l=1,nlay
[633]	426	DO ig=1,ngrid
[703]	427	pdqcloud(ig,l,igcm_h2o_ice) =
[1909]	428	& sum_subpdq(ig,l,igcm_h2o_ice)/real(imicro)
[703]	429	& - pdq(ig,l,igcm_h2o_ice)
	430	pdqcloud(ig,l,igcm_h2o_vap) =
[1909]	431	& sum_subpdq(ig,l,igcm_h2o_vap)/real(imicro)
[703]	432	& - pdq(ig,l,igcm_h2o_vap)
[2312]	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
[740]	441	ENDDO
	442	ENDDO
	443
	444	IF(microphys) THEN
	445	DO l=1,nlay
	446	DO ig=1,ngrid
[703]	447	pdqcloud(ig,l,igcm_ccn_mass) =
[1909]	448	& sum_subpdq(ig,l,igcm_ccn_mass)/real(imicro)
[703]	449	& - pdq(ig,l,igcm_ccn_mass)
	450	pdqcloud(ig,l,igcm_ccn_number) =
[1909]	451	& sum_subpdq(ig,l,igcm_ccn_number)/real(imicro)
[703]	452	& - pdq(ig,l,igcm_ccn_number)
[633]	453	ENDDO
[740]	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) =
[1909]	461	& sum_subpdq(ig,l,igcm_dust_mass)/real(imicro)
[740]	462	& - pdq(ig,l,igcm_dust_mass)
	463	pdqcloud(ig,l,igcm_dust_number) =
[1909]	464	& sum_subpdq(ig,l,igcm_dust_number)/real(imicro)
[740]	465	& - pdq(ig,l,igcm_dust_number)
	466	ENDDO
	467	ENDDO
	468	ENDIF
[633]	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
[654]	475	IF ((pq(ig,l,igcm_ccn_number) +
[633]	476	& ptimestep* (pdq(ig,l,igcm_ccn_number) +
[654]	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
[633]	481	pdqcloud(ig,l,igcm_ccn_number) =
	482	& - pq(ig,l,igcm_ccn_number)/ptimestep
[654]	483	& - pdq(ig,l,igcm_ccn_number) + 1.
[633]	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
[654]	488	& - pdq(ig,l,igcm_ccn_mass) + 1.e-20
[633]	489	pdqcloud(ig,l,igcm_dust_mass) =
	490	& -pdqcloud(ig,l,igcm_ccn_mass)
	491	ENDIF
	492	ENDDO
	493	ENDDO
	494	ENDIF
	495
[740]	496	IF(scavenging) THEN
[633]	497	DO l=1,nlay
	498	DO ig=1,ngrid
[740]	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
[2312]	522
[633]	523	IF (pq(ig,l,igcm_h2o_ice) + ptimestep*
	524	& (pdq(ig,l,igcm_h2o_ice) + pdqcloud(ig,l,igcm_h2o_ice))
[2516]	525	& .le. qperemin) THEN
[633]	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)
[2312]	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
[633]	534	ENDIF
	535	IF (pq(ig,l,igcm_h2o_vap) + ptimestep*
	536	& (pdq(ig,l,igcm_h2o_vap) + pdqcloud(ig,l,igcm_h2o_vap))
[2516]	537	& .le. qperemin) THEN
[633]	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)
[2312]	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
[633]	546	ENDIF
[2312]	547
[633]	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
[740]	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))
[633]	566
	567	ENDDO
	568	ENDDO
[740]	569	ENDIF
[1467]	570
[740]	571	IF(microphys) THEN
[1467]	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
[740]	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
[645]	609	ENDDO
[740]	610	ENDDO
[645]	611
[740]	612	ELSE ! no microphys
	613
[645]	614	DO l=1,nlay
	615	DO ig=1,ngrid
[740]	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
[746]	621	& tau(ig,1),pzlay(ig,l),rice(ig,l))
[740]	622
[633]	623	ENDDO
[740]	624	ENDDO
[633]	625
[740]	626	ENDIF ! of IF(microphys)
[633]	627
[740]	628
	629
[358]	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
[1047]	633	do ig=1,ngrid
[358]	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
[740]	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
[38]	654
[1711]	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
[1758]	675	#ifndef MESOSCALE
[1711]	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))
[2312]	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
[1711]	689	call WRITEDIAGFI(ngrid,"pdqccn2","pdqcloudccn apres microphysique"
	690	& ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
	691	& igcm_ccn_mass))
[2494]	692	call WRITEDIAGFI(ngrid,"pdqccnN2","pdqcloudccnN apres "//
	693	& "microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
[1711]	694	& igcm_ccn_number))
[2494]	695	call WRITEDIAGFI(ngrid,"pdqdust2", "pdqclouddust apres "//
	696	& "microphysique","kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
[1711]	697	& igcm_dust_mass))
[2494]	698	call WRITEDIAGFI(ngrid,"pdqdustN2", "pdqclouddustN apres "//
	699	& "microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
[1711]	700	& igcm_dust_number))
[633]	701	c=======================================================================
[1758]	702	#endif
[633]	703
[1711]	704	END SUBROUTINE watercloud
	705
[2162]	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
[1711]	729	END MODULE watercloud_mod

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

Download in other formats:

Original Format