Context Navigation

radii_mod.F90 @ 2987

Last change on this file since 2987 was 2957, checked in by emillour, 19 months ago

Generic PCM:
Fix a buggy behavior concerning H2O aerosol variance; aeroptproperties is not
designed to handle aerosol variance which is not constant, whereas h2o_reffrad
returns a variance which varies (between 0.09 and 0.13) with location and time.
Revert to a simpler setup where the H2O aerosol variance is uniform and set by
the user (nueff_iaero_h2o flag in callphys.def; default value 0.1)
Also added some "intent()" in optci arguments and increased length of string
to store varaible name in writediagfi.
EM

File size: 15.1 KB

Line
1	!==================================================================
2	module radii_mod
3	!==================================================================
4	! module to centralize the radii calculations for aerosols
5	!==================================================================
6
7	! CO2 cloud properties (initialized in inifis)
8	real,save :: Nmix_co2 ! Number mixing ratio of CO2 ice particles
9	!$OMP THREADPRIVATE(Nmix_co2)
10
11	! flag to specify if we assume a constant fixed radius for particles
12	logical,save :: radfixed ! initialized in inifis
13	!$OMP THREADPRIVATE(radfixed)
14
15	! water cloud optical properties (initialized in su_aer_radii below)
16	real, save :: rad_h2o
17	real, save :: rad_h2o_ice
18	real, save :: Nmix_h2o
19	real, save :: Nmix_h2o_ice
20	!$OMP THREADPRIVATE(rad_h2o,rad_h2o_ice,Nmix_h2o,Nmix_h2o_ice)
21
22	real,save :: nueff_iaero_h2o ! effective variance of H2O aerosol
23	! (initialized in su_aer_radii below)
24	!$OMP THREADPRIVATE(nueff_iaero_h2o)
25	! coefficients for a variable nueff() for h2o aerosol; disabled for now
26	real, parameter :: coef_hot=0.13
27	real, parameter :: coef_cold=0.09
28
29
30	contains
31
32
33	!==================================================================
34	subroutine su_aer_radii(ngrid,nlayer,reffrad,nueffrad)
35	!==================================================================
36	! Purpose
37	! -------
38	! Compute the effective radii of liquid and icy water particles
39	! Jeremy Leconte (2012)
40	! Extended to dust, CO2, NH3, 2-lay,Nlay,auroral aerosols by ??
41	! Added Radiative Generic Condensable Species effective radii
42	! calculations (Lucas Teinturier 2022)
43	!
44	! Authors
45	! -------
46	! Jeremy Leconte (2012)
47	!
48	!==================================================================
49	use mod_phys_lmdz_para, only : is_master
50	use ioipsl_getin_p_mod, only: getin_p
51	use radinc_h, only: naerkind
52	use aerosol_mod, only: iaero_back2lay, iaero_co2, iaero_dust, &
53	iaero_h2o, iaero_h2so4, iaero_nh3, iaero_nlay, &
54	iaero_aurora, iaero_generic, i_rgcs_ice, &
55	iaero_venus1, iaero_venus2, iaero_venus2p, &
56	iaero_venus3, iaero_venusUV
57	use callkeys_mod, only: size_nh3_cloud, nlayaero, aeronlay_size, &
58	aeronlay_nueff,aerogeneric
59	use tracer_h, only: radius, nqtot, is_rgcs
60	Implicit none
61
62	integer,intent(in) :: ngrid
63	integer,intent(in) :: nlayer
64
65	real, intent(out) :: reffrad(ngrid,nlayer,naerkind) !aerosols radii (K)
66	real, intent(out) :: nueffrad(ngrid,nlayer,naerkind) !variance
67
68	logical, save :: firstcall=.true.
69	!$OMP THREADPRIVATE(firstcall)
70	integer :: iaer, ia , iq, i_rad
71
72	do iaer=1,naerkind
73	! these values will change once the microphysics gets to work
74	! UNLESS tracer=.false., in which case we should be working with
75	! a fixed aerosol layer, and be able to define reffrad in a
76	! .def file. To be improved!
77	! \|-> Done in th n-layer aerosol case (JVO 20)
78
79	if(iaer.eq.iaero_co2)then ! CO2 ice
80	reffrad(1:ngrid,1:nlayer,iaer) = 1.e-4
81	nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
82	endif
83
84	if(iaer.eq.iaero_h2o)then ! H2O ice
85	nueff_iaero_h2o=0.1 ! default value for variance of h2o aerosols
86	call getin_p("nueff_iaero_h2o",nueff_iaero_h2o)
87	if (is_master) write(,)" nueff_iaero_h2o = ",nueff_iaero_h2o
88	reffrad(1:ngrid,1:nlayer,iaer) = 1.e-5
89	nueffrad(1:ngrid,1:nlayer,iaer) = nueff_iaero_h2o
90	endif
91
92	if(iaer.eq.iaero_dust)then ! dust
93	reffrad(1:ngrid,1:nlayer,iaer) = 1.e-5
94	nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
95	endif
96
97	if(iaer.eq.iaero_h2so4)then ! H2SO4 ice
98	reffrad(1:ngrid,1:nlayer,iaer) = 1.e-6
99	nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
100	endif
101
102	if(iaer.eq.iaero_back2lay)then ! Two-layer aerosols
103	reffrad(1:ngrid,1:nlayer,iaer) = 2.e-6
104	nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
105	endif
106
107
108	if(iaer.eq.iaero_nh3)then ! Nh3 cloud
109	reffrad(1:ngrid,1:nlayer,iaer) = size_nh3_cloud
110	nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
111	endif
112
113	do ia=1,nlayaero
114	if(iaer.eq.iaero_nlay(ia))then ! N-layer aerosols
115	reffrad(1:ngrid,1:nlayer,iaer) = aeronlay_size(ia)
116	nueffrad(1:ngrid,1:nlayer,iaer) = aeronlay_nueff(ia)
117	endif
118	enddo
119
120	if(iaer.eq.iaero_aurora)then ! Auroral aerosols
121	reffrad(1:ngrid,1:nlayer,iaer) = 3.e-7
122	nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
123	endif
124
125	if(iaer.eq.iaero_venus1)then ! Venus cloud, mode 1, Haus13 model
126	reffrad(1:ngrid,1:nlayer,iaer) = 0.49e-6
127	nueffrad(1:ngrid,1:nlayer,iaer) = 0.21
128	endif
129
130	if(iaer.eq.iaero_venus2)then ! Venus cloud, mode 2, Haus13 model
131	reffrad(1:ngrid,1:nlayer,iaer) = 1.23e-6
132	nueffrad(1:ngrid,1:nlayer,iaer) = 0.067
133	endif
134
135	if(iaer.eq.iaero_venus2p)then ! Venus cloud, mode 2p, Haus13 model
136	reffrad(1:ngrid,1:nlayer,iaer) = 1.56e-6
137	nueffrad(1:ngrid,1:nlayer,iaer) = 0.044
138	endif
139
140	if(iaer.eq.iaero_venus3)then ! Venus cloud, mode 3, Haus13 model
141	reffrad(1:ngrid,1:nlayer,iaer) = 4.25e-6
142	nueffrad(1:ngrid,1:nlayer,iaer) = 0.062
143	endif
144
145	if(iaer.eq.iaero_venusUV)then ! Venus cloud, UV abs, 1 val as in table
146	reffrad(1:ngrid,1:nlayer,iaer) = 0.5e-6
147	nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
148	endif
149
150	do ia=1,aerogeneric ! Radiative Generic Condensable Species
151	if (iaer .eq. iaero_generic(ia)) then
152	i_rad = i_rgcs_ice(ia)
153	reffrad(1:ngrid,1:nlayer,iaer)=radius(i_rad)
154	nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
155	endif
156	enddo ! generic radiative condensable aerosols
157
158	enddo ! iaer=1,naerkind
159
160
161	if (radfixed) then
162
163	if (is_master) write(,)"radius of H2O water particles:"
164	rad_h2o=13. ! default value
165	call getin_p("rad_h2o",rad_h2o)
166	if (is_master) write(,)" rad_h2o = ",rad_h2o
167
168	if (is_master) write(,)"radius of H2O ice particles:"
169	rad_h2o_ice=35. ! default value
170	call getin_p("rad_h2o_ice",rad_h2o_ice)
171	if (is_master) write(,)" rad_h2o_ice = ",rad_h2o_ice
172
173	else
174
175	if (is_master) write(,)"Number mixing ratio of H2O water particles:"
176	Nmix_h2o=1.e6 ! default value
177	call getin_p("Nmix_h2o",Nmix_h2o)
178	if (is_master) write(,)" Nmix_h2o = ",Nmix_h2o
179
180	if (is_master) write(,)"Number mixing ratio of H2O ice particles:"
181	Nmix_h2o_ice=Nmix_h2o ! default value
182	call getin_p("Nmix_h2o_ice",Nmix_h2o_ice)
183	if (is_master) write(,)" Nmix_h2o_ice = ",Nmix_h2o_ice
184	endif
185
186
187	end subroutine su_aer_radii
188	!==================================================================
189
190
191	!==================================================================
192	subroutine h2o_reffrad(ngrid,nlayer,pq,pt,reffrad,nueffrad)
193	!==================================================================
194	! Purpose
195	! -------
196	! Compute the effective radii of liquid and icy water particles
197	!
198	! Authors
199	! -------
200	! Jeremy Leconte (2012)
201	!
202	!==================================================================
203	use watercommon_h, Only: T_h2O_ice_liq,T_h2O_ice_clouds,rhowater,rhowaterice
204	use comcstfi_mod, only: pi
205	Implicit none
206
207	integer,intent(in) :: ngrid
208	integer,intent(in) :: nlayer
209
210	real, intent(in) :: pq(ngrid,nlayer) !water ice mixing ratios (kg/kg)
211	real, intent(in) :: pt(ngrid,nlayer) !temperature (K)
212	real, intent(out) :: reffrad(ngrid,nlayer) !aerosol radii
213	real, intent(out) :: nueffrad(ngrid,nlayer) ! dispersion
214
215	integer :: ig,l
216	real zfice ,zrad,zrad_liq,zrad_ice
217	real,external :: CBRT
218
219
220	if (radfixed) then
221	do l=1,nlayer
222	do ig=1,ngrid
223	zfice = 1.0 - (pt(ig,l)-T_h2O_ice_clouds) / (T_h2O_ice_liq-T_h2O_ice_clouds)
224	zfice = MIN(MAX(zfice,0.0),1.0)
225	reffrad(ig,l)= rad_h2o * (1.-zfice) + rad_h2o_ice * zfice
226	! nueffrad(ig,l) = coef_hot * (1.-zfice) + coef_cold * zfice
227	enddo
228	enddo
229	else
230	do l=1,nlayer
231	do ig=1,ngrid
232	zfice = 1.0 - (pt(ig,l)-T_h2O_ice_clouds) / (T_h2O_ice_liq-T_h2O_ice_clouds)
233	zfice = MIN(MAX(zfice,0.0),1.0)
234	zrad_liq = CBRT( 3pq(ig,l)/(4Nmix_h2opirhowater) )
235	zrad_ice = CBRT( 3pq(ig,l)/(4Nmix_h2o_icepirhowaterice) )
236	! nueffrad(ig,l) = coef_hot * (1.-zfice) + coef_cold * zfice
237	zrad = zrad_liq * (1.-zfice) + zrad_ice * zfice
238
239	reffrad(ig,l) = min(max(zrad,1.e-6),1000.e-6)
240	enddo
241	enddo
242	end if
243
244	! For now only constant nueff is enabled (otherwise some specific handling
245	! of variable nueff is required in aeroptproperties)
246	nueffrad(1:ngrid,1:nlayer)=nueff_iaero_h2o
247
248	end subroutine h2o_reffrad
249	!==================================================================
250
251
252	!==================================================================
253	subroutine h2o_cloudrad(ngrid,nlayer,pql,reffliq,reffice)
254	!==================================================================
255	! Purpose
256	! -------
257	! Compute the effective radii of liquid and icy water particles
258	!
259	! Authors
260	! -------
261	! Jeremy Leconte (2012)
262	!
263	!==================================================================
264	use watercommon_h, Only: rhowater,rhowaterice
265	use comcstfi_mod, only: pi
266	Implicit none
267
268	integer,intent(in) :: ngrid
269	integer,intent(in) :: nlayer
270
271	real, intent(in) :: pql(ngrid,nlayer) !condensed water mixing ratios (kg/kg)
272	real, intent(out) :: reffliq(ngrid,nlayer),reffice(ngrid,nlayer) !liquid and ice water particle radii (m)
273
274	real,external :: CBRT
275	integer :: i,k
276
277	if (radfixed) then
278	reffliq(1:ngrid,1:nlayer)= rad_h2o
279	reffice(1:ngrid,1:nlayer)= rad_h2o_ice
280	else
281	do k=1,nlayer
282	do i=1,ngrid
283	reffliq(i,k) = CBRT(3pql(i,k)/(4Nmix_h2opirhowater))
284	reffliq(i,k) = min(max(reffliq(i,k),1.e-6),1000.e-6)
285
286	reffice(i,k) = CBRT(3pql(i,k)/(4Nmix_h2o_icepirhowaterice))
287	reffice(i,k) = min(max(reffice(i,k),1.e-6),1000.e-6)
288	enddo
289	enddo
290	endif
291
292	end subroutine h2o_cloudrad
293	!==================================================================
294
295
296
297	!==================================================================
298	subroutine co2_reffrad(ngrid,nlayer,nq,pq,reffrad)
299	!==================================================================
300	! Purpose
301	! -------
302	! Compute the effective radii of co2 ice particles
303	!
304	! Authors
305	! -------
306	! Jeremy Leconte (2012)
307	!
308	!==================================================================
309	USE tracer_h, only:igcm_co2_ice,rho_co2
310	use comcstfi_mod, only: pi
311	Implicit none
312
313	integer,intent(in) :: ngrid,nlayer,nq
314
315	real, intent(in) :: pq(ngrid,nlayer,nq) !tracer mixing ratios (kg/kg)
316	real, intent(out) :: reffrad(ngrid,nlayer) !co2 ice particles radii (m)
317
318	integer :: ig,l
319	real :: zrad
320	real,external :: CBRT
321
322
323
324	if (radfixed) then
325	reffrad(1:ngrid,1:nlayer) = 5.e-5 ! CO2 ice
326	else
327	do l=1,nlayer
328	do ig=1,ngrid
329	zrad = CBRT( 3pq(ig,l,igcm_co2_ice)/(4Nmix_co2pirho_co2) )
330	reffrad(ig,l) = min(max(zrad,1.e-6),100.e-6)
331	enddo
332	enddo
333	end if
334
335	end subroutine co2_reffrad
336	!==================================================================
337
338
339
340	!==================================================================
341	subroutine dust_reffrad(ngrid,nlayer,reffrad)
342	!==================================================================
343	! Purpose
344	! -------
345	! Compute the effective radii of dust particles
346	!
347	! Authors
348	! -------
349	! Jeremy Leconte (2012)
350	!
351	!==================================================================
352	Implicit none
353
354	integer,intent(in) :: ngrid
355	integer,intent(in) :: nlayer
356
357	real, intent(out) :: reffrad(ngrid,nlayer) !dust particles radii (m)
358
359	reffrad(1:ngrid,1:nlayer) = 2.e-6 ! dust
360
361	end subroutine dust_reffrad
362	!==================================================================
363
364
365	!==================================================================
366	subroutine h2so4_reffrad(ngrid,nlayer,reffrad)
367	!==================================================================
368	! Purpose
369	! -------
370	! Compute the effective radii of h2so4 particles
371	!
372	! Authors
373	! -------
374	! Jeremy Leconte (2012)
375	!
376	!==================================================================
377	Implicit none
378
379	integer,intent(in) :: ngrid
380	integer,intent(in) :: nlayer
381
382	real, intent(out) :: reffrad(ngrid,nlayer) !h2so4 particle radii (m)
383
384	reffrad(1:ngrid,1:nlayer) = 1.e-6 ! h2so4
385
386	end subroutine h2so4_reffrad
387	!==================================================================
388
389	!==================================================================
390	subroutine back2lay_reffrad(ngrid,reffrad,nlayer,pplev)
391	!==================================================================
392	! Purpose
393	! -------
394	! Compute the effective radii of particles in a 2-layer model
395	!
396	! Authors
397	! -------
398	! Sandrine Guerlet (2013)
399	!
400	!==================================================================
401	use callkeys_mod, only: pres_bottom_tropo,pres_top_tropo,size_tropo, &
402	pres_bottom_strato,size_strato
403
404	Implicit none
405
406	integer,intent(in) :: ngrid
407
408	real, intent(out) :: reffrad(ngrid,nlayer) ! particle radii (m)
409	REAL,INTENT(IN) :: pplev(ngrid,nlayer+1) ! inter-layer pressure (Pa)
410	INTEGER,INTENT(IN) :: nlayer ! number of atmospheric layers
411	REAL :: expfactor
412	INTEGER l,ig
413
414	reffrad(:,:)=1e-6 !!initialization, not important
415	DO ig=1,ngrid
416	DO l=1,nlayer-1
417	IF (pplev(ig,l) .le. pres_bottom_tropo .and. pplev(ig,l) .ge. pres_top_tropo) THEN
418	reffrad(ig,l) = size_tropo
419	ELSEIF (pplev(ig,l) .lt. pres_top_tropo .and. pplev(ig,l) .gt. pres_bottom_strato) THEN
420	expfactor=log(size_strato/size_tropo) / log(pres_bottom_strato/pres_top_tropo)
421	reffrad(ig,l)= size_tropo((pplev(ig,l)/pres_top_tropo)*expfactor)
422	ELSEIF (pplev(ig,l) .le. pres_bottom_strato) then
423	reffrad(ig,l) = size_strato
424	ENDIF
425	ENDDO
426	ENDDO
427
428	end subroutine back2lay_reffrad
429	!==================================================================
430
431	end module radii_mod
432	!==================================================================

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

Download in other formats:

Original Format