source: codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/radii_mod.F90 @ 222

Last change on this file since 222 was 222, checked in by ymipsl, 10 years ago

Creating temporary dynamico/lmdz/saturn branche

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