SURFEX  V8_0
Surfex V8_0 release
 All Classes Files Functions Variables
mode_aer_surf.F90
Go to the documentation of this file.
1 !SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
2 !SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
3 !SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
4 !SFX_LIC for details. version 1.
5 !! ########################
7 !! ########################
8 !!
9 !! MODULE DUST PSD (Particle Size Distribution)
10 !! Purpose: Contains subroutines to convert from transported variables (ppp)
11 !! to understandable aerosol variables, e.g. #/m3, kg/m3, sigma, R_{n}
12 !!
13 !! Modifications
14 !! M.Leriche 2015 : masse molaire Black carbon à 12 g/mol
15 !-------------------------------------------------------------------------------
16 !! MODIFICATIONS
17 !! -------------
18 !!
19 !! J.Escobar 06/2013 for REAL4/8 add EPSILON management
20 !!
21 !-------------------------------------------------------------------------------
23  USE modd_dst_surf, ONLY : xdensity_dst
24 !
25  USE yomhook ,ONLY : lhook, dr_hook
26  USE parkind1 ,ONLY : jprb
27  USE modd_surf_par , ONLY : xsurf_tiny
28 !
29  IMPLICIT NONE
30 !!
31  CONTAINS
32 !!
33 SUBROUTINE init_var(PSV,PFAC,PCTOTA)
34 !
35 IMPLICIT NONE
36 !
37 REAL,DIMENSION(:,:), INTENT(IN) :: psv ! [aerosol concentration]
38 REAL,DIMENSION(:), INTENT(OUT) :: pfac ! M3 / mass conversion factor
39 REAL,DIMENSION(:,:,:), INTENT(OUT) :: pctota
40 !
41 REAL,DIMENSION(NSP+NCARB+NSOA) :: zmi ! [kg/mol] molar weight of aerosol
42 REAL,DIMENSION(NSP+NCARB+NSOA) :: zrhoi ! aerosol density
43 REAL :: zpi
44 INTEGER :: jj
45 REAL, PARAMETER :: zmol = 6.0221367e+11
46 REAL(KIND=JPRB) :: zhook_handle
47 !
48 IF (lhook) CALL dr_hook('MODE_AER_SURF:INIT_VAR',0,zhook_handle)
49 !
50 zrhoi(:) = 1.8e3
51 zrhoi(jp_aer_h2o) = 1.0e3 ! water
52 zrhoi(jp_aer_dst) = xdensity_dst
53 !
54 ! Moments index
55 nm0(1) = 1
56 nm3(1) = 2
57 nm6(1) = 3
58 nm0(2) = 4
59 nm3(2) = 5
60 nm6(2) = 6
61 !
62 !Set density of aerosol, here (kg/m3)
63 ! Aerosol Density
64 ! Cf Ackermann (all to black carbon except water)
65 !
66 !Set molecular weightn g/mol
67 zmi(:) = 250.
68 zmi(jp_aer_so4) = 98.
69 zmi(jp_aer_no3) = 63.
70 zmi(jp_aer_nh3) = 17.
71 zmi(jp_aer_h2o) = 18.
72 zmi(jp_aer_bc) = 12.
73 zmi(jp_aer_dst) = 100.
74 IF (nsoa .EQ. 10) THEN
75  zmi(jp_aer_soa1) = 88.
76  zmi(jp_aer_soa2) = 180.
77  zmi(jp_aer_soa3) = 1.5374857e3
78  zmi(jp_aer_soa4) = 1.9586780e3
79  zmi(jp_aer_soa5) = 195.
80  zmi(jp_aer_soa6) = 195.
81  zmi(jp_aer_soa7) = 165.
82  zmi(jp_aer_soa8) = 195.
83  zmi(jp_aer_soa9) = 270.
84  zmi(jp_aer_soa10) = 210.
85 END IF
86 !
87 zpi = 2.*asin(1.)
88 DO jj=1,nsp+ncarb+nsoa
89  pfac(jj)=(4./3.)*zpi*zrhoi(jj)*1.e-9
90 ENDDO
91 !
92 !* 2 transfer aerosol mass from gas to aerosol variables
93 ! (and conversion of mol.cm-3 --> microgram/m3)
94 !
95 pctota(:,:,:) = 0.
96 ! aerosol phase
97 pctota(:,jp_aer_so4,1) = psv(:,jp_ch_so4i)*zmi(jp_aer_so4)/zmol
98 pctota(:,jp_aer_so4,2) = psv(:,jp_ch_so4j)*zmi(jp_aer_so4)/zmol
99 
100 pctota(:,jp_aer_no3,1) = psv(:,jp_ch_no3i)*zmi(jp_aer_no3)/zmol
101 pctota(:,jp_aer_no3,2) = psv(:,jp_ch_no3j)*zmi(jp_aer_no3)/zmol
102 
103 pctota(:,jp_aer_nh3,1) = psv(:,jp_ch_nh3i)*zmi(jp_aer_nh3)/zmol
104 pctota(:,jp_aer_nh3,2) = psv(:,jp_ch_nh3j)*zmi(jp_aer_nh3)/zmol
105 !
106 ! water
107 pctota(:,jp_aer_h2o,1) = psv(:,jp_ch_h2oi)*zmi(jp_aer_h2o)/zmol
108 pctota(:,jp_aer_h2o,2) = psv(:,jp_ch_h2oj)*zmi(jp_aer_h2o)/zmol
109 !
110 ! primary organic carbon
111 pctota(:,jp_aer_oc,1) = psv(:,jp_ch_oci)*zmi(jp_aer_oc)/zmol
112 pctota(:,jp_aer_oc,2) = psv(:,jp_ch_ocj)*zmi(jp_aer_oc)/zmol
113 !
114 ! primary black carbon
115 pctota(:,jp_aer_bc,1) = psv(:,jp_ch_bci)*zmi(jp_aer_bc)/zmol
116 pctota(:,jp_aer_bc,2) = psv(:,jp_ch_bcj)*zmi(jp_aer_bc)/zmol
117 !
118 !dust
119 pctota(:,jp_aer_dst,1) = psv(:,jp_ch_dsti)*zmi(jp_aer_dst)/6.0221367e+11
120 pctota(:,jp_aer_dst,2) = psv(:,jp_ch_dstj)*zmi(jp_aer_dst)/6.0221367e+11
121 !
122 IF (nsoa .EQ. 10) THEN
123  pctota(:,jp_aer_soa1,1) = psv(:,jp_ch_soa1i)*zmi(jp_aer_soa1)/zmol
124  pctota(:,jp_aer_soa1,2) = psv(:,jp_ch_soa1j)*zmi(jp_aer_soa1)/zmol
125  pctota(:,jp_aer_soa2,1) = psv(:,jp_ch_soa2i)*zmi(jp_aer_soa2)/zmol
126  pctota(:,jp_aer_soa2,2) = psv(:,jp_ch_soa2j)*zmi(jp_aer_soa2)/zmol
127  pctota(:,jp_aer_soa3,1) = psv(:,jp_ch_soa3i)*zmi(jp_aer_soa3)/zmol
128  pctota(:,jp_aer_soa3,2) = psv(:,jp_ch_soa3j)*zmi(jp_aer_soa3)/zmol
129  pctota(:,jp_aer_soa4,1) = psv(:,jp_ch_soa4i)*zmi(jp_aer_soa4)/zmol
130  pctota(:,jp_aer_soa4,2) = psv(:,jp_ch_soa4j)*zmi(jp_aer_soa4)/zmol
131  pctota(:,jp_aer_soa5,1) = psv(:,jp_ch_soa5i)*zmi(jp_aer_soa5)/zmol
132  pctota(:,jp_aer_soa5,2) = psv(:,jp_ch_soa5j)*zmi(jp_aer_soa5)/zmol
133 
134  pctota(:,jp_aer_soa6,1) = psv(:,jp_ch_soa6i)*zmi(jp_aer_soa6)/zmol
135  pctota(:,jp_aer_soa6,2) = psv(:,jp_ch_soa6j)*zmi(jp_aer_soa6)/zmol
136  pctota(:,jp_aer_soa7,1) = psv(:,jp_ch_soa7i)*zmi(jp_aer_soa7)/zmol
137  pctota(:,jp_aer_soa7,2) = psv(:,jp_ch_soa7j)*zmi(jp_aer_soa7)/zmol
138  pctota(:,jp_aer_soa8,1) = psv(:,jp_ch_soa8i)*zmi(jp_aer_soa8)/zmol
139  pctota(:,jp_aer_soa8,2) = psv(:,jp_ch_soa8j)*zmi(jp_aer_soa8)/zmol
140  pctota(:,jp_aer_soa9,1) = psv(:,jp_ch_soa9i)*zmi(jp_aer_soa9)/zmol
141  pctota(:,jp_aer_soa9,2) = psv(:,jp_ch_soa9j)*zmi(jp_aer_soa9)/zmol
142  pctota(:,jp_aer_soa10,1) = psv(:,jp_ch_soa10i)*zmi(jp_aer_soa10)/zmol
143  pctota(:,jp_aer_soa10,2) = psv(:,jp_ch_soa10j)*zmi(jp_aer_soa10)/zmol
144 END IF
145 !
146 IF (lhook) CALL dr_hook('MODE_AER_SURF:INIT_VAR',1,zhook_handle)
147 !
148 END SUBROUTINE init_var
149 !
150 !! ############################################################
151 SUBROUTINE ppp2aero_surf( &
152  psvt, &!I [ppp] input scalar variables (moment of distribution)
153  prhodref, &!I [kg/m3] density of air
154  psig1d, &!O [-] standard deviation of aerosol distribution
155  prg1d, &!O [um] number median diameter of aerosol distribution
156  pn1d, &!O [#/m3] number concentration of aerosols
157  pctota, &!O [ug/m3] mass of each aerosol compounds
158  pm1d &!moments 0, 3 and 6
159  )
160 !! ############################################################
161 !
162 !!
163 !! PURPOSE
164 !! -------
165 !! Translate the three moments M0, M3 and M6 given in ppp into
166 !! Values which can be understood more easily (R, sigma, N, M)
167 !!
168 !! CALLING STRUCTURE NOTE: OPTIONAL VARIABLES
169 !! -------
170 !! CALL PPP2AERO_SURFS(PSVT, PRHODREF, PSIG1D=SIGVAR, &
171 !! PRG1D=RVAR, PN1D=NVAR, PM1D=ZM)
172 !!
173 !! REFERENCE
174 !! ---------
175 !! none
176 !!
177 !! AUTHOR
178 !! ------
179 !! Pierre TULET (LA)
180 !!
181 !! MODIFICATIONS
182 !! -------------
183 !! Alf Grini (CNRM)
184 !!
185 !! EXTERNAL
186 !! --------
187 !!
188  IMPLICIT NONE
189 !!
190 !-------------------------------------------------------------------------------
191 !
192 !* 0. DECLARATIONS
193 ! ------------
194 !
195 !* 0.1 declarations of arguments
196 !
197 !INPUT
198 REAL, DIMENSION(:,:), INTENT(IN) :: psvt !I [#/molec_{air}] first moment
199  !I [molec_{aer}/molec_{air} 3rd moment
200  !I [um6/molec_{air}*(cm3/m3)] 6th moment
201 REAL, DIMENSION(:), INTENT(IN) :: prhodref !I [kg/m3] density of air
202 
203 !OUTPUT
204 REAL, DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: psig1d !O [-] standard deviation
205 REAL, DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: prg1d !O [um] number median diameter
206 REAL, DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: pn1d !O [#/m3] number concentration
207 REAL, DIMENSION(:,:,:),OPTIONAL, INTENT(OUT) :: pctota !O [ug/m3] mass of each component
208 REAL, DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: pm1d !O moments 0,3 and 6
209 !
210 !* 0.2 declarations local variables
211 !
212 REAL,DIMENSION(SIZE(PSVT,1), SIZE(PSVT,2)) :: zsv ! [aerosol concentration]
213 REAL,DIMENSION(SIZE(PSVT,1)) :: zsigma ! [-] standard deviation
214 REAL,DIMENSION(SIZE(PSVT,1),NSP+NCARB+NSOA,JPMODE) :: zctota
215 REAL,DIMENSION(SIZE(PSVT,1),JPMODE*3) :: zm
216 REAL,DIMENSION(JPMODE*3) :: zmmin
217 !
218 REAL,DIMENSION(NSP+NCARB+NSOA) :: zfac ! M3 / mass conversion factor
219 REAL, PARAMETER :: zden2mol = 1e-6 * 6.0221367e+23 / 28.9644e-3
220 INTEGER :: jj, jn ! [idx] loop counters
221 REAL(KIND=JPRB) :: zhook_handle
222 !
223 !-------------------------------------------------------------------------------
224 IF (lhook) CALL dr_hook('MODE_AER_SURF:PPP2AERO_SURF',0,zhook_handle)
225 !
226 ! 1. initialisation
227 !
228 DO jj=1, SIZE(psvt,2)
229  zsv(:,jj) = psvt(:,jj) * zden2mol * prhodref(:)
230  zsv(:,jj) = max(zsv(:,jj),1e-40 * zden2mol * prhodref(:))
231 ENDDO
232 !
233  CALL init_var(zsv,zfac,zctota)
234 !
235 !* 2 calculate moment 3 from total aerosol mass
236 !
237 zm(:,2) = 0.
238 zm(:,5) = 0.
239 DO jj = 1,nsp+ncarb+nsoa
240  zm(:,2) = zm(:,2)+zctota(:,jj,1)/zfac(jj) !==>um3_{aer}/m3_{air} (volume ==> 3rd moment)
241  zm(:,5) = zm(:,5)+zctota(:,jj,2)/zfac(jj) !==>um3_{aer}/m3_{air} (volume ==> 3rd moment)
242 ENDDO
243 !
244 !
245 !* 3 set moment 0
246 !
247 zm(:,1)= max(zsv(:,jp_ch_m0i) * 1e+6, xsurf_tiny) ! molec_{aer}/m3_{air}
248 zm(:,4)= max(zsv(:,jp_ch_m0j) * 1e+6, xsurf_tiny) ! molec_{aer}/m3_{air}
249 WHERE ((zm(:,1) .LT. zmmin(1))) !.OR.(ZM(:,2) .LT. ZMMIN(2)))
250  zm(:,1)= zmmin(1)
251 ! ZM(:,2)= ZMMIN(2)
252 !
253 ! ZCTOTA(:,JP_AER_H2O,1) = 0.
254 ! ZCTOTA(:,JP_AER_NH3,1) = 0.
255 ! ZCTOTA(:,JP_AER_SO4,1) = 0.
256 ! ZCTOTA(:,JP_AER_NO3,1) = 0.
257 ! ZCTOTA(:,JP_AER_BC,1) = 0.5 * ZM(:,2) * ZFAC(JP_AER_BC)
258 ! ZCTOTA(:,JP_AER_OC,1) = 0.5 * ZM(:,2) * ZFAC(JP_AER_OC)
259 END WHERE
260 !!
261 WHERE ((zm(:,4) .LT. zmmin(4))) !.OR.(ZM(:,5) .LT. ZMMIN(5)))
262  zm(:,4)= zmmin(4)
263 ! ZM(:,5)= ZMMIN(5)
264 !
265 ! ZCTOTA(:,JP_AER_H2O,2) = 0.
266 ! ZCTOTA(:,JP_AER_NH3,2) = 0.
267 ! ZCTOTA(:,JP_AER_SO4,2) = 0.
268 ! ZCTOTA(:,JP_AER_NO3,2) = 0.
269 ! ZCTOTA(:,JP_AER_BC,2) = 0.5 * ZM(:,5) * ZFAC(JP_AER_BC)
270 ! ZCTOTA(:,JP_AER_OC,2) = 0.5 * ZM(:,5) * ZFAC(JP_AER_OC)
271 END WHERE
272 !
273 !* 4 set moment 6 ==> um6_{aer}/m3_{air}
274 !
275 IF (lvarsigi) THEN ! set M6 variable standard deviation
276  zm(:,3) = max(zsv(:,jp_ch_m6i), xsurf_tiny)
277 
278  zsigma(:)=zm(:,2)**2/(zm(:,1)*zm(:,3))
279  zsigma(:)=min(1-1e-10,zsigma(:))
280  zsigma(:)=max(1e-10,zsigma(:))
281  zsigma(:)= log(zsigma(:))
282  zsigma(:)= exp(1./3.*sqrt(-zsigma(:)))
283  zm(:,3) = zm(:,1) &
284  * ( (zm(:,2)/zm(:,1))**(1./3.) &
285  * exp(-(3./2.)*log(zsigma(:))**2))**6 &
286  * exp(18.*log(zsigma(:))**2)
287 
288  IF(present(psig1d)) psig1d(:,1) = zsigma(:)
289 
290 ELSE ! fixed standard deviation
291  zm(:,3) = zm(:,1) &
292  * ( (zm(:,2)/zm(:,1))**(1./3.) &
293  * exp(-(3./2.)*log(xemissigi)**2))**6 &
294  * exp(18.*log(xemissigi)**2)
295 
296  IF(present(psig1d)) psig1d(:,1) = xemissigi
297 END IF
298 
299 IF (lvarsigj) THEN ! set M6 variable standard deviation
300  zm(:,6) = max(zsv(:,jp_ch_m6j), xsurf_tiny)
301 
302  zsigma(:)=zm(:,5)**2/(zm(:,4)*zm(:,6))
303  zsigma(:)=min(1-1e-10,zsigma(:))
304  zsigma(:)=max(1e-10,zsigma(:))
305  zsigma(:)= log(zsigma(:))
306  zsigma(:)= exp(1./3.*sqrt(-zsigma(:)))
307 
308  zm(:,6) = zm(:,4) &
309  * ( (zm(:,5)/zm(:,4))**(1./3.) &
310  * exp(-(3./2.)*log(zsigma(:))**2))**6 &
311  * exp(18.*log(zsigma(:))**2)
312 
313  IF(present(psig1d)) psig1d(:,2) = zsigma(:)
314 
315 ELSE ! fixed standard deviation
316  zm(:,6) = zm(:,4) &
317  * ( (zm(:,5)/zm(:,4))**(1./3.) &
318  * exp(-(3./2.)*log(xemissigj)**2))**6 &
319  * exp(18.*log(xemissigj)**2)
320 
321  IF(present(psig1d)) psig1d(:,2) = xemissigj
322 END IF
323 
324 
325 !* 6 calculate modal parameters from moments
326 DO jn=1,jpmode
327  IF(present(pn1d)) pn1d(:,jn) = zm(:,nm0(jn))
328 
329  IF(present(prg1d)) prg1d(:,jn)=(zm(:,nm3(jn))**4. &
330  / (zm(:,nm6(jn))*zm(:,nm0(jn))**3.))**(1./6.)
331 ENDDO
332 
333 IF(present(pctota)) pctota(:,:,:) = zctota(:,:,:)
334 IF(present(pm1d)) pm1d(:,:) = zm(:,:)
335 IF (lhook) CALL dr_hook('MODE_AER_SURF:PPP2AERO_SURF',1,zhook_handle)
336 !
337 END SUBROUTINE ppp2aero_surf
338 !! ############################################################
339 SUBROUTINE aero2ppp_surf( &
340  psvt, &!IO [ppp] input scalar variables (moment of distribution)
341  prhodref, &!I [kg/m3] density of air
342  psig1d, &!I [-] standard deviation of aerosol distribution
343  prg1d &!I [um] number median diameter of aerosol distribution
344  )
345 !! ############################################################
346 !
347 !!
348 !! PURPOSE
349 !! -------
350 !! Translate the aerosol Mass, RG and SIGMA in the three moments M0, M3 and M6 given in ppp
351 !!
352 !! REFERENCE
353 !! ---------
354 !! none
355 !!
356 !! AUTHOR
357 !! ------
358 !! Pierre TULET (LA)
359 !!
360 !! MODIFICATIONS
361 !! -------------
362 !! Alf Grini (CNRM)
363 !!
364 !! EXTERNAL
365 !! --------
366 !!
367 IMPLICIT NONE
368 !!
369 !-------------------------------------------------------------------------------
370 !
371 !* 0. DECLARATIONS
372 ! ------------
373 !
374 !* 0.1 declarations of arguments
375 !
376 !INPUT
377 REAL, DIMENSION(:,:), INTENT(INOUT) :: psvt !I [#/molec_{air}] first moment
378  !I [molec_{aer}/molec_{air} 3rd moment
379  !I [um6/molec_{air}*(cm3/m3)] 6th moment
380 REAL, DIMENSION(:), INTENT(IN) :: prhodref !I [kg/m3] density of air
381 
382 !OUTPUT
383 REAL, DIMENSION(:,:), INTENT(IN) :: psig1d !O [-] standard deviation
384 REAL, DIMENSION(:,:), INTENT(IN) :: prg1d !O [um] number median diameter
385 !
386 !* 0.2 declarations local variables
387 !
388 REAL,DIMENSION(SIZE(PSVT,1),NSP+NCARB+NSOA,JPMODE):: zctota
389 REAL,DIMENSION(SIZE(PSVT,1), JPMODE*3) :: zm ! [aerosol units] local array which goes to output later
390 !
391 REAL,DIMENSION(NSP+NCARB+NSOA) :: zfac ! M3 / mass conversion factor
392 REAL, PARAMETER :: zden2mol = 1e-6 * 6.0221367e+23 / 28.9644e-3
393 INTEGER :: jj ! [idx] loop counters
394 REAL(KIND=JPRB) :: zhook_handle
395 !
396 !-------------------------------------------------------------------------------
397 IF (lhook) CALL dr_hook('MODE_AER_SURF:AERO2PPP_SURF',0,zhook_handle)
398 !
399 ! 1. initialisation
400 !
401 DO jj=1, SIZE(psvt, 2)
402  psvt(:,jj) = psvt(:,jj) * zden2mol * prhodref(:)
403 ENDDO
404 !
405  CALL init_var(psvt,zfac,zctota)
406 !
407 !* 3 calculate moment 3 from total aerosol mass
408 !
409 zm(:,2) = 0.
410 zm(:,5) = 0.
411 DO jj = 1,nsp+ncarb+nsoa
412  zm(:,2) = zm(:,2)+zctota(:,jj,1)/zfac(jj)
413  zm(:,5) = zm(:,5)+zctota(:,jj,2)/zfac(jj)
414 ENDDO
415 !
416 !
417 !* 4 calculate moment 0 from dispersion and mean radius
418 !
419 zm(:,1)= zm(:,2)/ ( (prg1d(:,1)**3)*exp(4.5 * log(psig1d(:,1))**2) )
420 zm(:,4)= zm(:,5)/ ( (prg1d(:,2)**3)*exp(4.5 * log(psig1d(:,2))**2) )
421 !
422 !* 5 calculate moment 6 from dispersion and mean radius
423 !
424 zm(:,3) = zm(:,1)*(prg1d(:,1)**6) * exp(18 *(log(psig1d(:,1)))**2)
425 zm(:,6) = zm(:,4)*(prg1d(:,2)**6) * exp(18 *(log(psig1d(:,2)))**2)
426 !
427 !* 6 return to ppp
428 !
429 psvt(:,jp_ch_m0i) = zm(:,1) * 1e-6
430 psvt(:,jp_ch_m0j) = zm(:,4) * 1e-6
431 !
432 IF (lvarsigi) psvt(:,jp_ch_m6i) = zm(:,3)
433 IF (lvarsigj) psvt(:,jp_ch_m6j) = zm(:,6)
434 !
435 DO jj=1,SIZE(psvt,2)
436  psvt(:,jj) = psvt(:,jj) / (zden2mol * prhodref(:))
437 ENDDO
438 !
439 IF (lhook) CALL dr_hook('MODE_AER_SURF:AERO2PPP_SURF',1,zhook_handle)
440 !
441 END SUBROUTINE aero2ppp_surf
442 
443 END MODULE mode_aer_surf
subroutine ppp2aero_surf(PSVT, PRHODREF, PSIG1D, PRG1D, PN1D, PCTOTA, PM1D)
subroutine aero2ppp_surf(PSVT, PRHODREF, PSIG1D, PRG1D)
subroutine init_var(PSV, PFAC, PCTOTA)