Context Navigation

watercloud.F @ 937

Last change on this file since 937 was 882, checked in by tnavarro, 13 years ago
small correction in makegcm, latent heat in clouds, changed threshold in dust radius
File size: 15.4 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
[38]	9	IMPLICIT NONE
	10
[633]	11
[38]	12	c=======================================================================
[358]	13	c Water-ice cloud formation
	14	c
	15	c Includes two different schemes:
	16	c - A simplified scheme (see simpleclouds.F)
	17	c - An improved microphysical scheme (see improvedclouds.F)
[38]	18	c
[633]	19	c There is a time loop specific to cloud formation
	20	c due to timescales smaller than the GCM integration timestep.
	21	c
[358]	22	c Authors: Franck Montmessin, Francois Forget, Ehouarn Millour,
[522]	23	c J.-B. Madeleine, Thomas Navarro
[38]	24	c
[633]	25	c 2004 - 2012
[38]	26	c=======================================================================
	27
	28	c-----------------------------------------------------------------------
	29	c declarations:
	30	c -------------
	31
	32	#include "dimensions.h"
	33	#include "dimphys.h"
	34	#include "comcstfi.h"
	35	#include "callkeys.h"
	36	#include "tracer.h"
	37	#include "comgeomfi.h"
[358]	38	#include "dimradmars.h"
[38]	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
[633]	56	REAL tau(ngridmx,naerkind) ! Column dust optical depth at each point
	57	REAL tauscaling(ngridmx) ! Convertion factor for dust amount
	58	real rdust(ngridmx,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
[633]	71	real rsedcloud(ngridmx,nlay) ! Cloud sedimentation radius
	72	real rhocloud(ngridmx,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
	96
	97	INTEGER iq,ig,l
[38]	98	LOGICAL,SAVE :: firstcall=.true.
	99
	100	c ** un petit test de coherence
	101	c --------------------------
	102
	103	IF (firstcall) THEN
	104	IF(ngrid.NE.ngridmx) THEN
	105	PRINT*,'STOP dans watercloud'
	106	PRINT*,'probleme de dimensions :'
	107	PRINT*,'ngrid =',ngrid
	108	PRINT*,'ngridmx =',ngridmx
	109	STOP
	110	ENDIF
	111
	112	if (nq.gt.nqmx) then
	113	write(,) 'stop in watercloud (nq.gt.nqmx)!'
	114	write(,) 'nq=',nq,' nqmx=',nqmx
	115	stop
	116	endif
	117
[358]	118	write(,) "watercloud: igcm_h2o_vap=",igcm_h2o_vap
	119	write(,) " igcm_h2o_ice=",igcm_h2o_ice
[633]	120
	121	write(,) "time subsampling for microphysic ?"
	122	#ifdef MESOSCALE
	123	imicro = 2
	124	#else
	125	imicro = 15
	126	#endif
	127	call getin("imicro",imicro)
	128	write(,)"imicro = ",imicro
	129
	130	microtimestep = ptimestep/real(imicro)
	131	write(,)"Physical timestep is",ptimestep
	132	write(,)"Microphysics timestep is",microtimestep
[38]	133
	134	firstcall=.false.
	135	ENDIF ! of IF (firstcall)
[522]	136
[633]	137	c-----Initialization
	138	subpdq(1:ngrid,1:nlay,1:nq) = 0
	139	subpdt(1:ngrid,1:nlay) = 0
	140
	141	! default value if no ice
	142	rhocloud(1:ngrid,1:nlay) = rho_dust
[38]	143
[633]	144
	145
	146	c------------------------------------------------------------------
	147	c Time subsampling for microphysics
	148	c------------------------------------------------------------------
	149	DO microstep=1,imicro
[522]	150
[633]	151	c-------------------------------------------------------------------
	152	c 1. Tendencies:
	153	c------------------
[38]	154
[633]	155
	156	c------ Temperature tendency subpdt
	157	! Each microtimestep we give the cloud scheme a stepped entry subpdt instead of pdt
	158	! If imicro=1 subpdt is the same as pdt
	159	DO l=1,nlay
	160	DO ig=1,ngrid
	161	subpdt(ig,l) = subpdt(ig,l)
	162	& + pdt(ig,l) ! At each micro timestep we add pdt in order to have a stepped entry
	163	ENDDO
	164	ENDDO
	165	c------ Tracers tendencies subpdq
	166	c------ At each micro timestep we add pdq in order to have a stepped entry
	167	IF (microphys) THEN
	168	DO l=1,nlay
	169	DO ig=1,ngrid
	170	subpdq(ig,l,igcm_dust_mass) =
	171	& subpdq(ig,l,igcm_dust_mass)
	172	& + pdq(ig,l,igcm_dust_mass)
	173	subpdq(ig,l,igcm_dust_number) =
	174	& subpdq(ig,l,igcm_dust_number)
	175	& + pdq(ig,l,igcm_dust_number)
	176	subpdq(ig,l,igcm_ccn_mass) =
	177	& subpdq(ig,l,igcm_ccn_mass)
	178	& + pdq(ig,l,igcm_ccn_mass)
	179	subpdq(ig,l,igcm_ccn_number) =
	180	& subpdq(ig,l,igcm_ccn_number)
	181	& + pdq(ig,l,igcm_ccn_number)
	182	ENDDO
	183	ENDDO
	184	ENDIF
	185	DO l=1,nlay
	186	DO ig=1,ngrid
	187	subpdq(ig,l,igcm_h2o_ice) =
	188	& subpdq(ig,l,igcm_h2o_ice)
	189	& + pdq(ig,l,igcm_h2o_ice)
	190	subpdq(ig,l,igcm_h2o_vap) =
	191	& subpdq(ig,l,igcm_h2o_vap)
	192	& + pdq(ig,l,igcm_h2o_vap)
	193	ENDDO
	194	ENDDO
	195
	196
	197	c-------------------------------------------------------------------
	198	c 2. Main call to the different cloud schemes:
	199	c------------------------------------------------
	200	IF (microphys) THEN
	201	CALL improvedclouds(ngrid,nlay,microtimestep,
	202	& pplay,pt,subpdt,
	203	& pq,subpdq,subpdqcloud,subpdtcloud,
	204	& nq,tauscaling)
	205
	206	ELSE
	207	CALL simpleclouds(ngrid,nlay,microtimestep,
	208	& pplay,pzlay,pt,subpdt,
	209	& pq,subpdq,subpdqcloud,subpdtcloud,
[645]	210	& nq,tau,rice)
[633]	211	ENDIF
	212
	213
	214	c-------------------------------------------------------------------
	215	c 3. Updating tendencies after cloud scheme:
	216	c-----------------------------------------------
	217
	218	IF (microphys) THEN
	219	DO l=1,nlay
	220	DO ig=1,ngrid
	221	subpdq(ig,l,igcm_dust_mass) =
	222	& subpdq(ig,l,igcm_dust_mass)
	223	& + subpdqcloud(ig,l,igcm_dust_mass)
	224	subpdq(ig,l,igcm_dust_number) =
	225	& subpdq(ig,l,igcm_dust_number)
	226	& + subpdqcloud(ig,l,igcm_dust_number)
	227	subpdq(ig,l,igcm_ccn_mass) =
	228	& subpdq(ig,l,igcm_ccn_mass)
	229	& + subpdqcloud(ig,l,igcm_ccn_mass)
	230	subpdq(ig,l,igcm_ccn_number) =
	231	& subpdq(ig,l,igcm_ccn_number)
	232	& + subpdqcloud(ig,l,igcm_ccn_number)
	233	ENDDO
	234	ENDDO
	235	ENDIF
	236	DO l=1,nlay
	237	DO ig=1,ngrid
	238	subpdq(ig,l,igcm_h2o_ice) =
	239	& subpdq(ig,l,igcm_h2o_ice)
	240	& + subpdqcloud(ig,l,igcm_h2o_ice)
	241	subpdq(ig,l,igcm_h2o_vap) =
	242	& subpdq(ig,l,igcm_h2o_vap)
	243	& + subpdqcloud(ig,l,igcm_h2o_vap)
	244	ENDDO
	245	ENDDO
[882]	246
	247	IF (activice) THEN
	248	DO l=1,nlay
	249	DO ig=1,ngrid
	250	subpdt(ig,l) =
	251	& subpdt(ig,l) + subpdtcloud(ig,l)
	252	ENDDO
	253	ENDDO
	254	ENDIF
[633]	255
	256
	257	ENDDO ! of DO microstep=1,imicro
	258
	259	c-------------------------------------------------------------------
	260	c 6. Compute final tendencies after time loop:
	261	c------------------------------------------------
	262
	263	c------ Temperature tendency pdtcloud
	264	DO l=1,nlay
	265	DO ig=1,ngrid
	266	pdtcloud(ig,l) =
	267	& subpdt(ig,l)/real(imicro)-pdt(ig,l)
	268	ENDDO
	269	ENDDO
[740]	270
[633]	271	c------ Tracers tendencies pdqcloud
[703]	272	DO l=1,nlay
[633]	273	DO ig=1,ngrid
[703]	274	pdqcloud(ig,l,igcm_h2o_ice) =
	275	& subpdq(ig,l,igcm_h2o_ice)/real(imicro)
	276	& - pdq(ig,l,igcm_h2o_ice)
	277	pdqcloud(ig,l,igcm_h2o_vap) =
	278	& subpdq(ig,l,igcm_h2o_vap)/real(imicro)
	279	& - pdq(ig,l,igcm_h2o_vap)
[740]	280	ENDDO
	281	ENDDO
	282
	283	IF(microphys) THEN
	284	DO l=1,nlay
	285	DO ig=1,ngrid
[703]	286	pdqcloud(ig,l,igcm_ccn_mass) =
	287	& subpdq(ig,l,igcm_ccn_mass)/real(imicro)
	288	& - pdq(ig,l,igcm_ccn_mass)
	289	pdqcloud(ig,l,igcm_ccn_number) =
	290	& subpdq(ig,l,igcm_ccn_number)/real(imicro)
	291	& - pdq(ig,l,igcm_ccn_number)
[633]	292	ENDDO
[740]	293	ENDDO
	294	ENDIF
	295
	296	IF(scavenging) THEN
	297	DO l=1,nlay
	298	DO ig=1,ngrid
	299	pdqcloud(ig,l,igcm_dust_mass) =
	300	& subpdq(ig,l,igcm_dust_mass)/real(imicro)
	301	& - pdq(ig,l,igcm_dust_mass)
	302	pdqcloud(ig,l,igcm_dust_number) =
	303	& subpdq(ig,l,igcm_dust_number)/real(imicro)
	304	& - pdq(ig,l,igcm_dust_number)
	305	ENDDO
	306	ENDDO
	307	ENDIF
[633]	308
	309	c------- Due to stepped entry, other processes tendencies can add up to negative values
	310	c------- Therefore, enforce positive values and conserve mass
	311
[740]	312
[633]	313	IF(microphys) THEN
	314	DO l=1,nlay
	315	DO ig=1,ngrid
[654]	316	IF ((pq(ig,l,igcm_ccn_number) +
[633]	317	& ptimestep* (pdq(ig,l,igcm_ccn_number) +
[654]	318	& pdqcloud(ig,l,igcm_ccn_number)) .le. 1.)
	319	& .or. (pq(ig,l,igcm_ccn_mass) +
	320	& ptimestep* (pdq(ig,l,igcm_ccn_mass) +
	321	& pdqcloud(ig,l,igcm_ccn_mass)) .le. 1.e-20)) THEN
[633]	322	pdqcloud(ig,l,igcm_ccn_number) =
	323	& - pq(ig,l,igcm_ccn_number)/ptimestep
[654]	324	& - pdq(ig,l,igcm_ccn_number) + 1.
[633]	325	pdqcloud(ig,l,igcm_dust_number) =
	326	& -pdqcloud(ig,l,igcm_ccn_number)
	327	pdqcloud(ig,l,igcm_ccn_mass) =
	328	& - pq(ig,l,igcm_ccn_mass)/ptimestep
[654]	329	& - pdq(ig,l,igcm_ccn_mass) + 1.e-20
[633]	330	pdqcloud(ig,l,igcm_dust_mass) =
	331	& -pdqcloud(ig,l,igcm_ccn_mass)
	332	ENDIF
	333	ENDDO
	334	ENDDO
	335	ENDIF
	336
[740]	337	IF(scavenging) THEN
[633]	338	DO l=1,nlay
	339	DO ig=1,ngrid
[740]	340	IF ((pq(ig,l,igcm_dust_number) +
	341	& ptimestep* (pdq(ig,l,igcm_dust_number) +
	342	& pdqcloud(ig,l,igcm_dust_number)) .le. 1.)
	343	& .or. (pq(ig,l,igcm_dust_mass) +
	344	& ptimestep* (pdq(ig,l,igcm_dust_mass) +
	345	& pdqcloud(ig,l,igcm_dust_mass)) .le. 1.e-20)) THEN
	346	pdqcloud(ig,l,igcm_dust_number) =
	347	& - pq(ig,l,igcm_dust_number)/ptimestep
	348	& - pdq(ig,l,igcm_dust_number) + 1.
	349	pdqcloud(ig,l,igcm_ccn_number) =
	350	& -pdqcloud(ig,l,igcm_dust_number)
	351	pdqcloud(ig,l,igcm_dust_mass) =
	352	& - pq(ig,l,igcm_dust_mass)/ptimestep
	353	& - pdq(ig,l,igcm_dust_mass) + 1.e-20
	354	pdqcloud(ig,l,igcm_ccn_mass) =
	355	& -pdqcloud(ig,l,igcm_dust_mass)
	356	ENDIF
	357	ENDDO
	358	ENDDO
	359	ENDIF
	360
	361	DO l=1,nlay
	362	DO ig=1,ngrid
[633]	363	IF (pq(ig,l,igcm_h2o_ice) + ptimestep*
	364	& (pdq(ig,l,igcm_h2o_ice) + pdqcloud(ig,l,igcm_h2o_ice))
	365	& .le. 1.e-8) THEN
	366	pdqcloud(ig,l,igcm_h2o_ice) =
	367	& - pq(ig,l,igcm_h2o_ice)/ptimestep - pdq(ig,l,igcm_h2o_ice)
	368	pdqcloud(ig,l,igcm_h2o_vap) = -pdqcloud(ig,l,igcm_h2o_ice)
	369	ENDIF
	370	IF (pq(ig,l,igcm_h2o_vap) + ptimestep*
	371	& (pdq(ig,l,igcm_h2o_vap) + pdqcloud(ig,l,igcm_h2o_vap))
	372	& .le. 1.e-8) THEN
	373	pdqcloud(ig,l,igcm_h2o_vap) =
	374	& - pq(ig,l,igcm_h2o_vap)/ptimestep - pdq(ig,l,igcm_h2o_vap)
	375	pdqcloud(ig,l,igcm_h2o_ice) = -pdqcloud(ig,l,igcm_h2o_vap)
	376	ENDIF
	377	ENDDO
	378	ENDDO
	379
	380
	381	c------Update the ice and dust particle size "rice" for output or photochemistry
	382	c------Only rsedcloud is used for the water cycle
	383
	384	IF(scavenging) THEN
	385	DO l=1, nlay
	386	DO ig=1,ngrid
	387
[740]	388	call updaterdust(
	389	& pq(ig,l,igcm_dust_mass) + ! dust mass
	390	& (pdq(ig,l,igcm_dust_mass) + ! dust mass
	391	& pdqcloud(ig,l,igcm_dust_mass))*ptimestep, ! dust mass
	392	& pq(ig,l,igcm_dust_number) + ! dust number
	393	& (pdq(ig,l,igcm_dust_number) + ! dust number
	394	& pdqcloud(ig,l,igcm_dust_number))*ptimestep, ! dust number
	395	& rdust(ig,l))
[633]	396
	397	ENDDO
	398	ENDDO
[740]	399	ENDIF
[645]	400
[740]	401
	402	IF(microphys) THEN
	403
	404	DO l=1, nlay
	405	DO ig=1,ngrid
	406
	407	call updaterice_micro(
	408	& pq(ig,l,igcm_h2o_ice) + ! ice mass
	409	& (pdq(ig,l,igcm_h2o_ice) + ! ice mass
	410	& pdqcloud(ig,l,igcm_h2o_ice))*ptimestep, ! ice mass
	411	& pq(ig,l,igcm_ccn_mass) + ! ccn mass
	412	& (pdq(ig,l,igcm_ccn_mass) + ! ccn mass
	413	& pdqcloud(ig,l,igcm_ccn_mass))*ptimestep, ! ccn mass
	414	& pq(ig,l,igcm_ccn_number) + ! ccn number
	415	& (pdq(ig,l,igcm_ccn_number) + ! ccn number
	416	& pdqcloud(ig,l,igcm_ccn_number))*ptimestep, ! ccn number
	417	& tauscaling(ig),rice(ig,l),rhocloud(ig,l))
	418
[645]	419	ENDDO
[740]	420	ENDDO
[645]	421
[740]	422	ELSE ! no microphys
	423
[645]	424	DO l=1,nlay
	425	DO ig=1,ngrid
[740]	426
	427	call updaterice_typ(
	428	& pq(ig,l,igcm_h2o_ice) + ! ice mass
	429	& (pdq(ig,l,igcm_h2o_ice) + ! ice mass
	430	& pdqcloud(ig,l,igcm_h2o_ice))*ptimestep, ! ice mass
[746]	431	& tau(ig,1),pzlay(ig,l),rice(ig,l))
[740]	432
[633]	433	ENDDO
[740]	434	ENDDO
[633]	435
[740]	436	ENDIF ! of IF(microphys)
[633]	437
[740]	438
	439
[358]	440	c A correction if a lot of subliming CO2 fills the 1st layer FF04/2005
	441	c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	442	c Then that should not affect the ice particle radius
	443	do ig=1,ngridmx
	444	if(pdpsrf(ig)ptimestep.gt.0.9(pplev(ig,1)-pplev(ig,2)))then
	445	if(pdpsrf(ig)ptimestep.gt.0.9(pplev(ig,1)-pplev(ig,3)))
	446	& rice(ig,2)=rice(ig,3)
	447	rice(ig,1)=rice(ig,2)
	448	end if
	449	end do
[740]	450
	451
	452	DO l=1,nlay
	453	DO ig=1,ngrid
	454	rsedcloud(ig,l)=max(rice(ig,l)*
	455	& (1.+nuice_sed)(1.+nuice_sed)(1.+nuice_sed),
	456	& rdust(ig,l))
	457	! rsedcloud(ig,l)=min(rsedcloud(ig,l),1.e-4)
	458	ENDDO
	459	ENDDO
	460
	461	! used for rad. transfer calculations
	462	! nuice is constant because a lognormal distribution is prescribed
	463	nuice(1:ngrid,1:nlay)=nuice_ref
[38]	464
[740]	465
	466
[633]	467	c=======================================================================
	468
[38]	469	END
	470

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

Download in other formats:

Original Format