Context Navigation

watercloud.F @ 1528

Last change on this file since 1528 was 1467, checked in by tnavarro, 10 years ago
new option supersat to allow for supersaturation of water
File size: 16.3 KB

Rev	Line
[633]	1	SUBROUTINE watercloud(ngrid,nlay,ptimestep,
	2	& pplev,pplay,pdpsrf,pzlay,pt,pdt,
[626]	3	& pq,pdq,pdqcloud,pdtcloud,
[358]	4	& nq,tau,tauscaling,rdust,rice,nuice,
	5	& rsedcloud,rhocloud)
[633]	6	! to use 'getin'
	7	USE ioipsl_getincom
[740]	8	USE updaterad
[1226]	9	USE comcstfi_h
[1036]	10	use tracer_mod, only: nqmx, igcm_h2o_vap, igcm_h2o_ice,
	11	& igcm_dust_mass, igcm_dust_number,
	12	& igcm_ccn_mass, igcm_ccn_number,
	13	& rho_dust, nuice_sed, nuice_ref
[1246]	14	use dimradmars_mod, only: naerkind
[38]	15	IMPLICIT NONE
	16
[633]	17
[38]	18	c=======================================================================
[358]	19	c Water-ice cloud formation
	20	c
	21	c Includes two different schemes:
	22	c - A simplified scheme (see simpleclouds.F)
	23	c - An improved microphysical scheme (see improvedclouds.F)
[38]	24	c
[633]	25	c There is a time loop specific to cloud formation
	26	c due to timescales smaller than the GCM integration timestep.
	27	c
[358]	28	c Authors: Franck Montmessin, Francois Forget, Ehouarn Millour,
[522]	29	c J.-B. Madeleine, Thomas Navarro
[38]	30	c
[633]	31	c 2004 - 2012
[38]	32	c=======================================================================
	33
	34	c-----------------------------------------------------------------------
	35	c declarations:
	36	c -------------
	37
	38	#include "callkeys.h"
	39
	40	c Inputs:
	41	c ------
	42
	43	INTEGER ngrid,nlay
[633]	44	INTEGER nq ! nombre de traceurs
[38]	45	REAL ptimestep ! pas de temps physique (s)
	46	REAL pplev(ngrid,nlay+1) ! pression aux inter-couches (Pa)
	47	REAL pplay(ngrid,nlay) ! pression au milieu des couches (Pa)
[633]	48	REAL pdpsrf(ngrid) ! tendence surf pressure
[38]	49	REAL pzlay(ngrid,nlay) ! altitude at the middle of the layers
	50	REAL pt(ngrid,nlay) ! temperature at the middle of the layers (K)
[633]	51	REAL pdt(ngrid,nlay) ! tendence temperature des autres param.
[38]	52
	53	real pq(ngrid,nlay,nq) ! traceur (kg/kg)
[633]	54	real pdq(ngrid,nlay,nq) ! tendence avant condensation (kg/kg.s-1)
[38]	55
[1047]	56	REAL tau(ngrid,naerkind) ! Column dust optical depth at each point
	57	REAL tauscaling(ngrid) ! Convertion factor for dust amount
	58	real rdust(ngrid,nlay) ! Dust geometric mean radius (m)
[38]	59
	60	c Outputs:
	61	c -------
	62
[633]	63	real pdqcloud(ngrid,nlay,nq) ! tendence de la condensation H2O(kg/kg.s-1)
	64	REAL pdtcloud(ngrid,nlay) ! tendence temperature due
	65	! a la chaleur latente
[38]	66
	67	REAL rice(ngrid,nlay) ! Ice mass mean radius (m)
	68	! (r_c in montmessin_2004)
	69	REAL nuice(ngrid,nlay) ! Estimated effective variance
	70	! of the size distribution
[1047]	71	real rsedcloud(ngrid,nlay) ! Cloud sedimentation radius
	72	real rhocloud(ngrid,nlay) ! Cloud density (kg.m-3)
[38]	73
	74	c local:
	75	c ------
[633]	76
	77	! for ice radius computation
	78	REAL Mo,No
	79	REAl ccntyp
	80
	81	! for time loop
	82	INTEGER microstep ! time subsampling step variable
	83	INTEGER imicro ! time subsampling for coupled water microphysics & sedimentation
	84	SAVE imicro
	85	REAL microtimestep ! integration timestep for coupled water microphysics & sedimentation
	86	SAVE microtimestep
	87
	88	! tendency given by clouds (inside the micro loop)
	89	REAL subpdqcloud(ngrid,nlay,nq) ! cf. pdqcloud
	90	REAL subpdtcloud(ngrid,nlay) ! cf. pdtcloud
[38]	91
[633]	92	! global tendency (clouds+physics)
	93	REAL subpdq(ngrid,nlay,nq) ! cf. pdqcloud
	94	REAL subpdt(ngrid,nlay) ! cf. pdtcloud
	95
[1467]	96	! no supersaturation when option supersat is false
	97	REAL zt(ngrid,nlay) ! local value of temperature
	98	REAL zqsat(ngrid,nlay) ! saturation
[633]	99
	100	INTEGER iq,ig,l
[38]	101	LOGICAL,SAVE :: firstcall=.true.
	102
	103	c ** un petit test de coherence
	104	c --------------------------
	105
	106	IF (firstcall) THEN
	107
	108	if (nq.gt.nqmx) then
	109	write(,) 'stop in watercloud (nq.gt.nqmx)!'
	110	write(,) 'nq=',nq,' nqmx=',nqmx
	111	stop
	112	endif
	113
[358]	114	write(,) "watercloud: igcm_h2o_vap=",igcm_h2o_vap
	115	write(,) " igcm_h2o_ice=",igcm_h2o_ice
[633]	116
	117	write(,) "time subsampling for microphysic ?"
	118	#ifdef MESOSCALE
	119	imicro = 2
	120	#else
[951]	121	imicro = 30
[633]	122	#endif
	123	call getin("imicro",imicro)
	124	write(,)"imicro = ",imicro
	125
	126	microtimestep = ptimestep/real(imicro)
	127	write(,)"Physical timestep is",ptimestep
	128	write(,)"Microphysics timestep is",microtimestep
[38]	129
	130	firstcall=.false.
	131	ENDIF ! of IF (firstcall)
[522]	132
[633]	133	c-----Initialization
	134	subpdq(1:ngrid,1:nlay,1:nq) = 0
	135	subpdt(1:ngrid,1:nlay) = 0
	136
	137	! default value if no ice
	138	rhocloud(1:ngrid,1:nlay) = rho_dust
[38]	139
[633]	140
	141
	142	c------------------------------------------------------------------
	143	c Time subsampling for microphysics
	144	c------------------------------------------------------------------
	145	DO microstep=1,imicro
[522]	146
[633]	147	c-------------------------------------------------------------------
	148	c 1. Tendencies:
	149	c------------------
[38]	150
[633]	151
	152	c------ Temperature tendency subpdt
	153	! Each microtimestep we give the cloud scheme a stepped entry subpdt instead of pdt
	154	! If imicro=1 subpdt is the same as pdt
	155	DO l=1,nlay
	156	DO ig=1,ngrid
	157	subpdt(ig,l) = subpdt(ig,l)
	158	& + pdt(ig,l) ! At each micro timestep we add pdt in order to have a stepped entry
	159	ENDDO
	160	ENDDO
	161	c------ Tracers tendencies subpdq
	162	c------ At each micro timestep we add pdq in order to have a stepped entry
	163	IF (microphys) THEN
	164	DO l=1,nlay
	165	DO ig=1,ngrid
	166	subpdq(ig,l,igcm_dust_mass) =
	167	& subpdq(ig,l,igcm_dust_mass)
	168	& + pdq(ig,l,igcm_dust_mass)
	169	subpdq(ig,l,igcm_dust_number) =
	170	& subpdq(ig,l,igcm_dust_number)
	171	& + pdq(ig,l,igcm_dust_number)
	172	subpdq(ig,l,igcm_ccn_mass) =
	173	& subpdq(ig,l,igcm_ccn_mass)
	174	& + pdq(ig,l,igcm_ccn_mass)
	175	subpdq(ig,l,igcm_ccn_number) =
	176	& subpdq(ig,l,igcm_ccn_number)
	177	& + pdq(ig,l,igcm_ccn_number)
	178	ENDDO
	179	ENDDO
	180	ENDIF
	181	DO l=1,nlay
	182	DO ig=1,ngrid
	183	subpdq(ig,l,igcm_h2o_ice) =
	184	& subpdq(ig,l,igcm_h2o_ice)
	185	& + pdq(ig,l,igcm_h2o_ice)
	186	subpdq(ig,l,igcm_h2o_vap) =
	187	& subpdq(ig,l,igcm_h2o_vap)
	188	& + pdq(ig,l,igcm_h2o_vap)
	189	ENDDO
	190	ENDDO
	191
	192
	193	c-------------------------------------------------------------------
	194	c 2. Main call to the different cloud schemes:
	195	c------------------------------------------------
	196	IF (microphys) THEN
	197	CALL improvedclouds(ngrid,nlay,microtimestep,
	198	& pplay,pt,subpdt,
	199	& pq,subpdq,subpdqcloud,subpdtcloud,
	200	& nq,tauscaling)
	201
	202	ELSE
	203	CALL simpleclouds(ngrid,nlay,microtimestep,
	204	& pplay,pzlay,pt,subpdt,
	205	& pq,subpdq,subpdqcloud,subpdtcloud,
[645]	206	& nq,tau,rice)
[633]	207	ENDIF
	208
	209
	210	c-------------------------------------------------------------------
	211	c 3. Updating tendencies after cloud scheme:
	212	c-----------------------------------------------
	213
	214	IF (microphys) THEN
	215	DO l=1,nlay
	216	DO ig=1,ngrid
	217	subpdq(ig,l,igcm_dust_mass) =
	218	& subpdq(ig,l,igcm_dust_mass)
	219	& + subpdqcloud(ig,l,igcm_dust_mass)
	220	subpdq(ig,l,igcm_dust_number) =
	221	& subpdq(ig,l,igcm_dust_number)
	222	& + subpdqcloud(ig,l,igcm_dust_number)
	223	subpdq(ig,l,igcm_ccn_mass) =
	224	& subpdq(ig,l,igcm_ccn_mass)
	225	& + subpdqcloud(ig,l,igcm_ccn_mass)
	226	subpdq(ig,l,igcm_ccn_number) =
	227	& subpdq(ig,l,igcm_ccn_number)
	228	& + subpdqcloud(ig,l,igcm_ccn_number)
	229	ENDDO
	230	ENDDO
	231	ENDIF
	232	DO l=1,nlay
	233	DO ig=1,ngrid
	234	subpdq(ig,l,igcm_h2o_ice) =
	235	& subpdq(ig,l,igcm_h2o_ice)
	236	& + subpdqcloud(ig,l,igcm_h2o_ice)
	237	subpdq(ig,l,igcm_h2o_vap) =
	238	& subpdq(ig,l,igcm_h2o_vap)
	239	& + subpdqcloud(ig,l,igcm_h2o_vap)
	240	ENDDO
	241	ENDDO
[882]	242
	243	IF (activice) THEN
	244	DO l=1,nlay
	245	DO ig=1,ngrid
	246	subpdt(ig,l) =
	247	& subpdt(ig,l) + subpdtcloud(ig,l)
	248	ENDDO
	249	ENDDO
	250	ENDIF
[633]	251
	252
	253	ENDDO ! of DO microstep=1,imicro
	254
	255	c-------------------------------------------------------------------
	256	c 6. Compute final tendencies after time loop:
	257	c------------------------------------------------
	258
	259	c------ Temperature tendency pdtcloud
	260	DO l=1,nlay
	261	DO ig=1,ngrid
	262	pdtcloud(ig,l) =
	263	& subpdt(ig,l)/real(imicro)-pdt(ig,l)
	264	ENDDO
	265	ENDDO
[740]	266
[633]	267	c------ Tracers tendencies pdqcloud
[703]	268	DO l=1,nlay
[633]	269	DO ig=1,ngrid
[703]	270	pdqcloud(ig,l,igcm_h2o_ice) =
	271	& subpdq(ig,l,igcm_h2o_ice)/real(imicro)
	272	& - pdq(ig,l,igcm_h2o_ice)
	273	pdqcloud(ig,l,igcm_h2o_vap) =
	274	& subpdq(ig,l,igcm_h2o_vap)/real(imicro)
	275	& - pdq(ig,l,igcm_h2o_vap)
[740]	276	ENDDO
	277	ENDDO
	278
	279	IF(microphys) THEN
	280	DO l=1,nlay
	281	DO ig=1,ngrid
[703]	282	pdqcloud(ig,l,igcm_ccn_mass) =
	283	& subpdq(ig,l,igcm_ccn_mass)/real(imicro)
	284	& - pdq(ig,l,igcm_ccn_mass)
	285	pdqcloud(ig,l,igcm_ccn_number) =
	286	& subpdq(ig,l,igcm_ccn_number)/real(imicro)
	287	& - pdq(ig,l,igcm_ccn_number)
[633]	288	ENDDO
[740]	289	ENDDO
	290	ENDIF
	291
	292	IF(scavenging) THEN
	293	DO l=1,nlay
	294	DO ig=1,ngrid
	295	pdqcloud(ig,l,igcm_dust_mass) =
	296	& subpdq(ig,l,igcm_dust_mass)/real(imicro)
	297	& - pdq(ig,l,igcm_dust_mass)
	298	pdqcloud(ig,l,igcm_dust_number) =
	299	& subpdq(ig,l,igcm_dust_number)/real(imicro)
	300	& - pdq(ig,l,igcm_dust_number)
	301	ENDDO
	302	ENDDO
	303	ENDIF
[633]	304
	305	c------- Due to stepped entry, other processes tendencies can add up to negative values
	306	c------- Therefore, enforce positive values and conserve mass
	307
[740]	308
[633]	309	IF(microphys) THEN
	310	DO l=1,nlay
	311	DO ig=1,ngrid
[654]	312	IF ((pq(ig,l,igcm_ccn_number) +
[633]	313	& ptimestep* (pdq(ig,l,igcm_ccn_number) +
[654]	314	& pdqcloud(ig,l,igcm_ccn_number)) .le. 1.)
	315	& .or. (pq(ig,l,igcm_ccn_mass) +
	316	& ptimestep* (pdq(ig,l,igcm_ccn_mass) +
	317	& pdqcloud(ig,l,igcm_ccn_mass)) .le. 1.e-20)) THEN
[633]	318	pdqcloud(ig,l,igcm_ccn_number) =
	319	& - pq(ig,l,igcm_ccn_number)/ptimestep
[654]	320	& - pdq(ig,l,igcm_ccn_number) + 1.
[633]	321	pdqcloud(ig,l,igcm_dust_number) =
	322	& -pdqcloud(ig,l,igcm_ccn_number)
	323	pdqcloud(ig,l,igcm_ccn_mass) =
	324	& - pq(ig,l,igcm_ccn_mass)/ptimestep
[654]	325	& - pdq(ig,l,igcm_ccn_mass) + 1.e-20
[633]	326	pdqcloud(ig,l,igcm_dust_mass) =
	327	& -pdqcloud(ig,l,igcm_ccn_mass)
	328	ENDIF
	329	ENDDO
	330	ENDDO
	331	ENDIF
	332
[740]	333	IF(scavenging) THEN
[633]	334	DO l=1,nlay
	335	DO ig=1,ngrid
[740]	336	IF ((pq(ig,l,igcm_dust_number) +
	337	& ptimestep* (pdq(ig,l,igcm_dust_number) +
	338	& pdqcloud(ig,l,igcm_dust_number)) .le. 1.)
	339	& .or. (pq(ig,l,igcm_dust_mass) +
	340	& ptimestep* (pdq(ig,l,igcm_dust_mass) +
	341	& pdqcloud(ig,l,igcm_dust_mass)) .le. 1.e-20)) THEN
	342	pdqcloud(ig,l,igcm_dust_number) =
	343	& - pq(ig,l,igcm_dust_number)/ptimestep
	344	& - pdq(ig,l,igcm_dust_number) + 1.
	345	pdqcloud(ig,l,igcm_ccn_number) =
	346	& -pdqcloud(ig,l,igcm_dust_number)
	347	pdqcloud(ig,l,igcm_dust_mass) =
	348	& - pq(ig,l,igcm_dust_mass)/ptimestep
	349	& - pdq(ig,l,igcm_dust_mass) + 1.e-20
	350	pdqcloud(ig,l,igcm_ccn_mass) =
	351	& -pdqcloud(ig,l,igcm_dust_mass)
	352	ENDIF
	353	ENDDO
	354	ENDDO
	355	ENDIF
	356
	357	DO l=1,nlay
	358	DO ig=1,ngrid
[633]	359	IF (pq(ig,l,igcm_h2o_ice) + ptimestep*
	360	& (pdq(ig,l,igcm_h2o_ice) + pdqcloud(ig,l,igcm_h2o_ice))
	361	& .le. 1.e-8) THEN
	362	pdqcloud(ig,l,igcm_h2o_ice) =
	363	& - pq(ig,l,igcm_h2o_ice)/ptimestep - pdq(ig,l,igcm_h2o_ice)
	364	pdqcloud(ig,l,igcm_h2o_vap) = -pdqcloud(ig,l,igcm_h2o_ice)
	365	ENDIF
	366	IF (pq(ig,l,igcm_h2o_vap) + ptimestep*
	367	& (pdq(ig,l,igcm_h2o_vap) + pdqcloud(ig,l,igcm_h2o_vap))
	368	& .le. 1.e-8) THEN
	369	pdqcloud(ig,l,igcm_h2o_vap) =
	370	& - pq(ig,l,igcm_h2o_vap)/ptimestep - pdq(ig,l,igcm_h2o_vap)
	371	pdqcloud(ig,l,igcm_h2o_ice) = -pdqcloud(ig,l,igcm_h2o_vap)
	372	ENDIF
	373	ENDDO
	374	ENDDO
	375
	376
	377	c------Update the ice and dust particle size "rice" for output or photochemistry
	378	c------Only rsedcloud is used for the water cycle
	379
	380	IF(scavenging) THEN
	381	DO l=1, nlay
	382	DO ig=1,ngrid
	383
[740]	384	call updaterdust(
	385	& pq(ig,l,igcm_dust_mass) + ! dust mass
	386	& (pdq(ig,l,igcm_dust_mass) + ! dust mass
	387	& pdqcloud(ig,l,igcm_dust_mass))*ptimestep, ! dust mass
	388	& pq(ig,l,igcm_dust_number) + ! dust number
	389	& (pdq(ig,l,igcm_dust_number) + ! dust number
	390	& pdqcloud(ig,l,igcm_dust_number))*ptimestep, ! dust number
	391	& rdust(ig,l))
[633]	392
	393	ENDDO
	394	ENDDO
[740]	395	ENDIF
[1467]	396
[740]	397	IF(microphys) THEN
[1467]	398
	399	! In case one does not want to allow supersatured water when using microphysics.
	400	! Not done by default.
	401	IF(.not.supersat) THEN
	402	zt = pt + (pdt+pdtcloud)*ptimestep
	403	call watersat(ngrid*nlay,zt,pplay,zqsat)
	404	DO l=1, nlay
	405	DO ig=1,ngrid
	406	IF (pq(ig,l,igcm_h2o_vap)
	407	& + (pdq(ig,l,igcm_h2o_vap) + pdqcloud(ig,l,igcm_h2o_vap))
	408	& * ptimestep .ge. zqsat(ig,l)) THEN
	409	pdqcloud(ig,l,igcm_h2o_vap) =
	410	& (zqsat(ig,l) - pq(ig,l,igcm_h2o_vap))/ptimestep
	411	& - pdq(ig,l,igcm_h2o_vap)
	412	pdqcloud(ig,l,igcm_h2o_ice) =
	413	& -pdqcloud(ig,l,igcm_h2o_vap)
	414	! no need to correct ccn_number, updaterad can handle this properly.
	415	ENDIF
	416	ENDDO
	417	ENDDO
	418	ENDIF
[740]	419
	420	DO l=1, nlay
	421	DO ig=1,ngrid
	422
	423	call updaterice_micro(
	424	& pq(ig,l,igcm_h2o_ice) + ! ice mass
	425	& (pdq(ig,l,igcm_h2o_ice) + ! ice mass
	426	& pdqcloud(ig,l,igcm_h2o_ice))*ptimestep, ! ice mass
	427	& pq(ig,l,igcm_ccn_mass) + ! ccn mass
	428	& (pdq(ig,l,igcm_ccn_mass) + ! ccn mass
	429	& pdqcloud(ig,l,igcm_ccn_mass))*ptimestep, ! ccn mass
	430	& pq(ig,l,igcm_ccn_number) + ! ccn number
	431	& (pdq(ig,l,igcm_ccn_number) + ! ccn number
	432	& pdqcloud(ig,l,igcm_ccn_number))*ptimestep, ! ccn number
	433	& tauscaling(ig),rice(ig,l),rhocloud(ig,l))
	434
[645]	435	ENDDO
[740]	436	ENDDO
[645]	437
[740]	438	ELSE ! no microphys
	439
[645]	440	DO l=1,nlay
	441	DO ig=1,ngrid
[740]	442
	443	call updaterice_typ(
	444	& pq(ig,l,igcm_h2o_ice) + ! ice mass
	445	& (pdq(ig,l,igcm_h2o_ice) + ! ice mass
	446	& pdqcloud(ig,l,igcm_h2o_ice))*ptimestep, ! ice mass
[746]	447	& tau(ig,1),pzlay(ig,l),rice(ig,l))
[740]	448
[633]	449	ENDDO
[740]	450	ENDDO
[633]	451
[740]	452	ENDIF ! of IF(microphys)
[633]	453
[740]	454
	455
[358]	456	c A correction if a lot of subliming CO2 fills the 1st layer FF04/2005
	457	c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	458	c Then that should not affect the ice particle radius
[1047]	459	do ig=1,ngrid
[358]	460	if(pdpsrf(ig)ptimestep.gt.0.9(pplev(ig,1)-pplev(ig,2)))then
	461	if(pdpsrf(ig)ptimestep.gt.0.9(pplev(ig,1)-pplev(ig,3)))
	462	& rice(ig,2)=rice(ig,3)
	463	rice(ig,1)=rice(ig,2)
	464	end if
	465	end do
[740]	466
	467
	468	DO l=1,nlay
	469	DO ig=1,ngrid
	470	rsedcloud(ig,l)=max(rice(ig,l)*
	471	& (1.+nuice_sed)(1.+nuice_sed)(1.+nuice_sed),
	472	& rdust(ig,l))
	473	! rsedcloud(ig,l)=min(rsedcloud(ig,l),1.e-4)
	474	ENDDO
	475	ENDDO
	476
	477	! used for rad. transfer calculations
	478	! nuice is constant because a lognormal distribution is prescribed
	479	nuice(1:ngrid,1:nlay)=nuice_ref
[38]	480
[740]	481
	482
[633]	483	c=======================================================================
	484
[38]	485	END
	486

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

Download in other formats:

Original Format