source: trunk/LMDZ.GENERIC/libf/phystd/radii_mod.F90 @ 2937

Last change on this file since 2937 was 2831, checked in by emillour, 3 years ago

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