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

Last change on this file since 2176 was 2005, checked in by sglmd, 6 years ago

in a old comit, the particle size for NH3 was fixed to 1 micrometer ; now read in callphys.def
File size: 12.5 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                pres_bottom_tropo,pres_top_tropo,size_tropo,        &
12                pres_bottom_strato,size_strato
13     
14      real, save ::  rad_h2o
15      real, save ::  rad_h2o_ice
16      real, save ::  Nmix_h2o
17      real, save ::  Nmix_h2o_ice
18!$OMP THREADPRIVATE(rad_h2o,rad_h2o_ice,Nmix_h2o,Nmix_h2o_ice)
19      real, parameter ::  coef_chaud=0.13
20      real, parameter ::  coef_froid=0.09
21
22
23contains
24
25
26!==================================================================
27   subroutine su_aer_radii(ngrid,nlayer,reffrad,nueffrad)
28!==================================================================
29!     Purpose
30!     -------
31!     Compute the effective radii of liquid and icy water particles
32!
33!     Authors
34!     -------
35!     Jeremy Leconte (2012)
36!
37!==================================================================
38      use ioipsl_getin_p_mod, only: getin_p
39      use radinc_h, only: naerkind
40      use aerosol_mod, only: iaero_back2lay, iaero_co2, iaero_dust, &
41                             iaero_h2o, iaero_h2so4,iaero_nh3,iaero_aurora
42      use callkeys_mod, only: size_nh3_cloud
43
44      Implicit none
45
46      integer,intent(in) :: ngrid
47      integer,intent(in) :: nlayer
48
49      real, intent(out) :: reffrad(ngrid,nlayer,naerkind)      !aerosols radii (K)
50      real, intent(out) :: nueffrad(ngrid,nlayer,naerkind)     !variance     
51
52      logical, save :: firstcall=.true.
53!$OMP THREADPRIVATE(firstcall)
54      integer :: iaer   
55     
56      print*,'enter su_aer_radii'
57          do iaer=1,naerkind
58!     these values will change once the microphysics gets to work
59!     UNLESS tracer=.false., in which case we should be working with
60!     a fixed aerosol layer, and be able to define reffrad in a
61!     .def file. To be improved!
62
63            if(iaer.eq.iaero_co2)then ! CO2 ice
64               reffrad(1:ngrid,1:nlayer,iaer) = 1.e-4
65               nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
66            endif
67
68            if(iaer.eq.iaero_h2o)then ! H2O ice
69               reffrad(1:ngrid,1:nlayer,iaer) = 1.e-5
70               nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
71            endif
72
73            if(iaer.eq.iaero_dust)then ! dust
74               reffrad(1:ngrid,1:nlayer,iaer) = 1.e-5
75               nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
76            endif
77 
78            if(iaer.eq.iaero_h2so4)then ! H2O ice
79               reffrad(1:ngrid,1:nlayer,iaer) = 1.e-6
80               nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
81            endif
82           
83            if(iaer.eq.iaero_back2lay)then ! Two-layer aerosols
84               reffrad(1:ngrid,1:nlayer,iaer) = 2.e-6
85               nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
86            endif
87
88
89            if(iaer.eq.iaero_nh3)then ! Nh3 cloud
90               reffrad(1:ngrid,1:nlayer,iaer) = size_nh3_cloud
91               nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
92            endif
93
94            if(iaer.eq.iaero_aurora)then ! Auroral aerosols
95               reffrad(1:ngrid,1:nlayer,iaer) = 3.e-7
96               nueffrad(1:ngrid,1:nlayer,iaer) = 0.1
97            endif
98
99
100            if(iaer.gt.5)then
101               print*,'Error in callcorrk, naerkind is too high (>5).'
102               print*,'The code still needs generalisation to arbitrary'
103               print*,'aerosol kinds and number.'
104               call abort
105            endif
106
107         enddo
108
109
110         if (radfixed) then
111
112            write(*,*)"radius of H2O water particles:"
113            rad_h2o=13. ! default value
114            call getin_p("rad_h2o",rad_h2o)
115            write(*,*)" rad_h2o = ",rad_h2o
116
117            write(*,*)"radius of H2O ice particles:"
118            rad_h2o_ice=35. ! default value
119            call getin_p("rad_h2o_ice",rad_h2o_ice)
120            write(*,*)" rad_h2o_ice = ",rad_h2o_ice
121
122         else
123
124            write(*,*)"Number mixing ratio of H2O water particles:"
125            Nmix_h2o=1.e6 ! default value
126            call getin_p("Nmix_h2o",Nmix_h2o)
127            write(*,*)" Nmix_h2o = ",Nmix_h2o
128
129            write(*,*)"Number mixing ratio of H2O ice particles:"
130            Nmix_h2o_ice=Nmix_h2o ! default value
131            call getin_p("Nmix_h2o_ice",Nmix_h2o_ice)
132            write(*,*)" Nmix_h2o_ice = ",Nmix_h2o_ice
133         endif
134
135      print*,'exit su_aer_radii'
136
137   end subroutine su_aer_radii
138!==================================================================
139
140
141!==================================================================
142   subroutine h2o_reffrad(ngrid,nlayer,pq,pt,reffrad,nueffrad)
143!==================================================================
144!     Purpose
145!     -------
146!     Compute the effective radii of liquid and icy water particles
147!
148!     Authors
149!     -------
150!     Jeremy Leconte (2012)
151!
152!==================================================================
153      use watercommon_h, Only: T_h2O_ice_liq,T_h2O_ice_clouds,rhowater,rhowaterice
154      use comcstfi_mod, only: pi
155      Implicit none
156
157      integer,intent(in) :: ngrid
158      integer,intent(in) :: nlayer
159
160      real, intent(in) :: pq(ngrid,nlayer) !water ice mixing ratios (kg/kg)
161      real, intent(in) :: pt(ngrid,nlayer) !temperature (K)
162      real, intent(out) :: reffrad(ngrid,nlayer)      !aerosol radii
163      real, intent(out) :: nueffrad(ngrid,nlayer) ! dispersion     
164
165      integer :: ig,l
166      real zfice ,zrad,zrad_liq,zrad_ice
167      real,external :: CBRT           
168     
169
170      if (radfixed) then
171         do l=1,nlayer
172            do ig=1,ngrid
173               zfice = 1.0 - (pt(ig,l)-T_h2O_ice_clouds) / (T_h2O_ice_liq-T_h2O_ice_clouds)
174               zfice = MIN(MAX(zfice,0.0),1.0)
175               reffrad(ig,l)= rad_h2o * (1.-zfice) + rad_h2o_ice * zfice
176               nueffrad(ig,l) = coef_chaud * (1.-zfice) + coef_froid * zfice
177            enddo
178         enddo
179      else
180         do l=1,nlayer
181            do ig=1,ngrid
182               zfice = 1.0 - (pt(ig,l)-T_h2O_ice_clouds) / (T_h2O_ice_liq-T_h2O_ice_clouds)
183               zfice = MIN(MAX(zfice,0.0),1.0)
184               zrad_liq  = CBRT( 3*pq(ig,l)/(4*Nmix_h2o*pi*rhowater) )
185               zrad_ice  = CBRT( 3*pq(ig,l)/(4*Nmix_h2o_ice*pi*rhowaterice) )
186               nueffrad(ig,l) = coef_chaud * (1.-zfice) + coef_froid * zfice
187               zrad = zrad_liq * (1.-zfice) + zrad_ice * zfice
188
189               reffrad(ig,l) = min(max(zrad,1.e-6),1000.e-6)
190               enddo
191            enddo     
192      end if
193
194   end subroutine h2o_reffrad
195!==================================================================
196
197
198!==================================================================
199   subroutine h2o_cloudrad(ngrid,nlayer,pql,reffliq,reffice)
200!==================================================================
201!     Purpose
202!     -------
203!     Compute the effective radii of liquid and icy water particles
204!
205!     Authors
206!     -------
207!     Jeremy Leconte (2012)
208!
209!==================================================================
210      use watercommon_h, Only: rhowater,rhowaterice
211      use comcstfi_mod, only: pi
212      Implicit none
213
214      integer,intent(in) :: ngrid
215      integer,intent(in) :: nlayer
216
217      real, intent(in) :: pql(ngrid,nlayer) !condensed water mixing ratios (kg/kg)
218      real, intent(out) :: reffliq(ngrid,nlayer),reffice(ngrid,nlayer)     !liquid and ice water particle radii (m)
219
220      real,external :: CBRT           
221      integer :: i,k
222
223      if (radfixed) then
224         reffliq(1:ngrid,1:nlayer)= rad_h2o
225         reffice(1:ngrid,1:nlayer)= rad_h2o_ice
226      else
227         do k=1,nlayer
228           do i=1,ngrid
229             reffliq(i,k) = CBRT(3*pql(i,k)/(4*Nmix_h2o*pi*rhowater))
230             reffliq(i,k) = min(max(reffliq(i,k),1.e-6),1000.e-6)
231           
232             reffice(i,k) = CBRT(3*pql(i,k)/(4*Nmix_h2o_ice*pi*rhowaterice))
233             reffice(i,k) = min(max(reffice(i,k),1.e-6),1000.e-6)
234           enddo
235         enddo
236      endif
237
238   end subroutine h2o_cloudrad
239!==================================================================
240
241
242
243!==================================================================
244   subroutine co2_reffrad(ngrid,nlayer,nq,pq,reffrad)
245!==================================================================
246!     Purpose
247!     -------
248!     Compute the effective radii of co2 ice particles
249!
250!     Authors
251!     -------
252!     Jeremy Leconte (2012)
253!
254!==================================================================
255      USE tracer_h, only:igcm_co2_ice,rho_co2
256      use comcstfi_mod, only: pi
257      Implicit none
258
259      integer,intent(in) :: ngrid,nlayer,nq
260
261      real, intent(in) :: pq(ngrid,nlayer,nq) !tracer mixing ratios (kg/kg)
262      real, intent(out) :: reffrad(ngrid,nlayer)      !co2 ice particles radii (m)
263
264      integer :: ig,l
265      real :: zrad   
266      real,external :: CBRT           
267           
268     
269
270      if (radfixed) then
271         reffrad(1:ngrid,1:nlayer) = 5.e-5 ! CO2 ice
272      else
273         do l=1,nlayer
274            do ig=1,ngrid
275               zrad = CBRT( 3*pq(ig,l,igcm_co2_ice)/(4*Nmix_co2*pi*rho_co2) )
276               reffrad(ig,l) = min(max(zrad,1.e-6),100.e-6)
277            enddo
278         enddo     
279      end if
280
281   end subroutine co2_reffrad
282!==================================================================
283
284
285
286!==================================================================
287   subroutine dust_reffrad(ngrid,nlayer,reffrad)
288!==================================================================
289!     Purpose
290!     -------
291!     Compute the effective radii of dust particles
292!
293!     Authors
294!     -------
295!     Jeremy Leconte (2012)
296!
297!==================================================================
298      Implicit none
299
300      integer,intent(in) :: ngrid
301      integer,intent(in) :: nlayer
302
303      real, intent(out) :: reffrad(ngrid,nlayer)      !dust particles radii (m)
304           
305      reffrad(1:ngrid,1:nlayer) = 2.e-6 ! dust
306
307   end subroutine dust_reffrad
308!==================================================================
309
310
311!==================================================================
312   subroutine h2so4_reffrad(ngrid,nlayer,reffrad)
313!==================================================================
314!     Purpose
315!     -------
316!     Compute the effective radii of h2so4 particles
317!
318!     Authors
319!     -------
320!     Jeremy Leconte (2012)
321!
322!==================================================================
323      Implicit none
324
325      integer,intent(in) :: ngrid
326      integer,intent(in) :: nlayer
327
328      real, intent(out) :: reffrad(ngrid,nlayer)      !h2so4 particle radii (m)
329               
330      reffrad(1:ngrid,1:nlayer) = 1.e-6 ! h2so4
331
332   end subroutine h2so4_reffrad
333!==================================================================
334
335!==================================================================
336   subroutine back2lay_reffrad(ngrid,reffrad,nlayer,pplev)
337!==================================================================
338!     Purpose
339!     -------
340!     Compute the effective radii of particles in a 2-layer model
341!
342!     Authors
343!     -------
344!     Sandrine Guerlet (2013)
345!
346!==================================================================
347 
348      use aerosol_mod   !! Particle sizes and boundaries of aerosol layers defined there
349     Implicit none
350
351      integer,intent(in) :: ngrid
352
353      real, intent(out) :: reffrad(ngrid,nlayer)      ! particle radii (m)
354      REAL,INTENT(IN) :: pplev(ngrid,nlayer+1) ! inter-layer pressure (Pa)
355      INTEGER,INTENT(IN) :: nlayer ! number of atmospheric layers
356      REAL :: expfactor
357      INTEGER l,ig
358           
359      reffrad(:,:)=1e-6  !!initialization, not important
360          DO ig=1,ngrid
361            DO l=1,nlayer-1
362              IF (pplev(ig,l) .le. pres_bottom_tropo .and. pplev(ig,l) .ge. pres_top_tropo) THEN
363                reffrad(ig,l) = size_tropo
364              ELSEIF (pplev(ig,l) .lt. pres_top_tropo .and. pplev(ig,l) .gt. pres_bottom_strato) THEN
365                expfactor=log(size_strato/size_tropo) / log(pres_bottom_strato/pres_top_tropo)
366                reffrad(ig,l)= size_tropo*((pplev(ig,l)/pres_top_tropo)**expfactor)
367              ELSEIF (pplev(ig,l) .le. pres_bottom_strato) then
368                reffrad(ig,l) = size_strato
369              ENDIF
370            ENDDO
371          ENDDO
372
373   end subroutine back2lay_reffrad
374!==================================================================
375
376
377end module radii_mod
378!==================================================================
Note: See TracBrowser for help on using the repository browser.