-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathmodal_aero_data.F90
434 lines (389 loc) · 20 KB
/
modal_aero_data.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
module modal_aero_data
!--------------------------------------------------------------
! ... Basic aerosol mode parameters and arrays
!--------------------------------------------------------------
use shr_kind_mod, only: r8 => shr_kind_r8
use constituents, only: pcnst
use radconstants, only: nswbands, nlwbands
implicit none
save
integer, parameter :: maxd_aspectype = 14
! aerosol mode definitions
!
#if ( defined MODAL_AERO_7MODE )
integer, parameter :: ntot_amode = 7
#elif ( defined MODAL_AERO_4MODE )
integer, parameter :: ntot_amode = 4
#elif ( defined MODAL_AERO_3MODE )
integer, parameter :: ntot_amode = 3
#elif ( defined MODAL_AERO_9MODE )
integer, parameter :: ntot_amode = 9 !kzm
#endif
!
! definitions for aerosol chemical components
!
integer, parameter :: ntot_aspectype = 8
character(len=*),parameter :: specname_amode(ntot_aspectype) = (/ 'sulfate ', 'ammonium ', 'nitrate ', &
'p-organic ', 's-organic ', 'black-c ', &
'seasalt ', 'dust ' /)
! set specdens_amode from physprop files via rad_cnst_get_aer_props
!specdens_amode(:ntot_aspectype) = (/1770.0,1770.0,1770.0, 1000.0, 1000.0, 1700.0,1900.0,2600.0 /)
! rce - 06-aug-2007 - changed specmw for almost everything to match mozart
#if ( defined MODAL_AERO_7MODE )
real(r8), parameter :: specmw_amode(ntot_aspectype) = (/ 96.0_r8, 18.0_r8, 62.0_r8, &
12.0_r8, 12.0_r8, 12.0_r8, 58.5_r8, 135.0_r8 /)
#elif ( defined MODAL_AERO_9MODE )
real(r8), parameter :: specmw_amode(ntot_aspectype) = (/ 96.0_r8, 18.0_r8, 62.0_r8, &
12.0_r8, 12.0_r8, 12.0_r8, 58.5_r8, 135.0_r8 /) !kzm
#elif ( defined MODAL_AERO_4MODE )
real(r8), parameter :: specmw_amode(ntot_aspectype) = (/ 115.0_r8, 115.0_r8, 62.0_r8, &
12.0_r8, 12.0_r8, 12.0_r8, 58.5_r8, 135.0_r8 /)
#elif ( defined MODAL_AERO_3MODE )
real(r8), parameter :: specmw_amode(ntot_aspectype) = (/ 115.0_r8, 115.0_r8, 62.0_r8, &
12.0_r8, 12.0_r8, 12.0_r8, 58.5_r8, 135.0_r8 /)
#endif
! input modename_amode, nspec_amode
#if ( defined MODAL_AERO_7MODE )
character(len=*), parameter :: modename_amode(ntot_amode) = (/ &
'accum ', &
'aitken ', &
'primary_carbon ', &
'fine_seasalt ', &
'fine_dust ', &
'coarse_seasalt ', &
'coarse_dust '/)
#elif ( defined MODAL_AERO_9MODE )
character(len=*), parameter :: modename_amode(ntot_amode) = (/ &
'accum ', &
'aitken ', &
'primary_carbon ', &
'fine_seasalt ', &
'fine_dust ', &
'coarse_seasalt ', &
'coarse_dust1 ', &
'coarse_dust2 ', &
'coarse_dust3 '/) !kzm
#elif ( defined MODAL_AERO_4MODE )
character(len=*), parameter :: modename_amode(ntot_amode) = (/ &
'accum ', &
'aitken ', &
'coarse ', &
'primary_carbon '/)
#elif ( defined MODAL_AERO_3MODE )
character(len=*), parameter :: modename_amode(ntot_amode) = (/ &
'accum ', &
'aitken ', &
'coarse '/)
#endif
#if ( defined MODAL_AERO_7MODE )
integer, parameter :: nspec_amode(ntot_amode) = (/ 6, 4, 2, 3, 3, 3, 3 /) ! SS
#elif ( defined MODAL_AERO_9MODE )
integer, parameter :: nspec_amode(ntot_amode) = (/ 6, 4, 2, 3, 3, 3, 3, 3, 3 /) !kzm
#elif ( defined MODAL_AERO_4MODE )
integer, parameter :: nspec_amode(ntot_amode) = (/ 6, 4, 3, 2 /) ! Aitken dust
#elif ( defined MODAL_AERO_3MODE )
integer, parameter :: nspec_amode(ntot_amode) = (/ 6, 3, 3 /)
#endif
integer, parameter :: nspec_amode_max = 6
! input mprognum_amode, mdiagnum_amode, mprogsfc_amode, mcalcwater_amode
#if ( defined MODAL_AERO_7MODE )
integer, parameter :: mprognum_amode(ntot_amode) = (/ 1, 1, 1, 1, 1, 1, 1/)
integer, parameter :: mdiagnum_amode(ntot_amode) = (/ 0, 0, 0, 0, 0, 0, 0/)
integer, parameter :: mprogsfc_amode(ntot_amode) = (/ 0, 0, 0, 0, 0, 0, 0/)
integer, parameter :: mcalcwater_amode(ntot_amode) = (/ 1, 1, 1, 1, 1, 1, 1/)
#elif ( defined MODAL_AERO_9MODE )
integer, parameter :: mprognum_amode(ntot_amode) = (/ 1, 1, 1, 1, 1, 1, 1, 1, 1/)!kzm
integer, parameter :: mdiagnum_amode(ntot_amode) = (/ 0, 0, 0, 0, 0, 0, 0, 0, 0/)!kzm
integer, parameter :: mprogsfc_amode(ntot_amode) = (/ 0, 0, 0, 0, 0, 0, 0, 0, 0/)!kzm
integer, parameter :: mcalcwater_amode(ntot_amode) = (/ 1, 1, 1, 1, 1, 1, 1, 1, 1/)!kzm
#elif ( defined MODAL_AERO_4MODE )
integer, parameter :: mprognum_amode(ntot_amode) = (/ 1, 1, 1, 1/)
integer, parameter :: mdiagnum_amode(ntot_amode) = (/ 0, 0, 0, 0/)
integer, parameter :: mprogsfc_amode(ntot_amode) = (/ 0, 0, 0, 0/)
integer, parameter :: mcalcwater_amode(ntot_amode) = (/ 0, 0, 0, 0/)
#elif ( defined MODAL_AERO_3MODE )
integer, parameter :: mprognum_amode(ntot_amode) = (/ 1, 1, 1/)
integer, parameter :: mdiagnum_amode(ntot_amode) = (/ 0, 0, 0/)
integer, parameter :: mprogsfc_amode(ntot_amode) = (/ 0, 0, 0/)
integer, parameter :: mcalcwater_amode(ntot_amode) = (/ 0, 0, 0/)
#endif
! input dgnum_amode, dgnumlo_amode, dgnumhi_amode (units = m)
real(r8) :: dgnum_amode(ntot_amode)
real(r8) :: dgnumlo_amode(ntot_amode)
real(r8) :: dgnumhi_amode(ntot_amode)
! input sigmag_amode
real(r8) :: sigmag_amode(ntot_amode)
! input crystalization and deliquescence points
real(r8) :: rhcrystal_amode(ntot_amode)
real(r8) :: rhdeliques_amode(ntot_amode)
integer :: msectional = -1
integer & !
lspectype_amode( maxd_aspectype, ntot_amode ), & !
lmassptr_amode( maxd_aspectype, ntot_amode ), & !
lmassptrcw_amode( maxd_aspectype, ntot_amode ), & !
numptr_amode( ntot_amode ), & !
numptrcw_amode( ntot_amode )
real(r8) :: & !
alnsg_amode( ntot_amode ), & !
voltonumb_amode( ntot_amode ), & !
voltonumblo_amode( ntot_amode ), & !
voltonumbhi_amode( ntot_amode ), & !
alnv2n_amode( ntot_amode ), & !
alnv2nlo_amode( ntot_amode ), & !
alnv2nhi_amode( ntot_amode ), & !
specdens_amode( maxd_aspectype ), & !
spechygro( maxd_aspectype )
complex(r8) & !
specrefndxsw( nswbands, maxd_aspectype ), & !
specrefndxlw( nlwbands, maxd_aspectype )
character(len=16) :: cnst_name_cw( pcnst )
character(len=8) :: aodvisname(ntot_amode ), &
ssavisname(ntot_amode )
character(len=48) :: aodvislongname(ntot_amode ), &
ssavislongname(ntot_amode )
character(len=8) :: fnactname(ntot_amode ), &
fmactname(ntot_amode ), &
nactname(ntot_amode )
character(len=48) :: fnactlongname(ntot_amode ), &
fmactlongname(ntot_amode ), &
nactlongname(ntot_amode )
integer & !
lptr_so4_a_amode(ntot_amode), lptr_so4_cw_amode(ntot_amode), & !
lptr_msa_a_amode(ntot_amode), lptr_msa_cw_amode(ntot_amode), & !
lptr_nh4_a_amode(ntot_amode), lptr_nh4_cw_amode(ntot_amode), & !
lptr_no3_a_amode(ntot_amode), lptr_no3_cw_amode(ntot_amode), & !
lptr_pom_a_amode(ntot_amode), lptr_pom_cw_amode(ntot_amode), & !
lptr_soa_a_amode(ntot_amode), lptr_soa_cw_amode(ntot_amode), & !
lptr_bc_a_amode(ntot_amode), lptr_bc_cw_amode(ntot_amode), & !
lptr_nacl_a_amode(ntot_amode), lptr_nacl_cw_amode(ntot_amode),& !
lptr_dust_a_amode(ntot_amode), lptr_dust_cw_amode(ntot_amode),& !
modeptr_accum, modeptr_aitken, & !
modeptr_ufine, modeptr_coarse, & !
modeptr_pcarbon, & !
modeptr_finedust, modeptr_fineseas, & !
modeptr_coardust, modeptr_coarseas, &
!---added by kzm
modeptr_coardust1, modeptr_coardust2, &
modeptr_coardust3
!----------------------
real(r8) :: &
specmw_so4_amode, specdens_so4_amode, &
specmw_nh4_amode, specdens_nh4_amode, &
specmw_no3_amode, specdens_no3_amode, &
specmw_pom_amode, specdens_pom_amode, &
specmw_soa_amode, specdens_soa_amode, &
specmw_bc_amode, specdens_bc_amode, &
specmw_dust_amode, specdens_dust_amode, &
specmw_seasalt_amode, specdens_seasalt_amode
integer species_class(pcnst) ! indicates species class (
! cldphysics, aerosol, gas )
integer spec_class_undefined
parameter ( spec_class_undefined = 0 )
integer spec_class_cldphysics
parameter ( spec_class_cldphysics = 1 )
integer spec_class_aerosol
parameter ( spec_class_aerosol = 2 )
integer spec_class_gas
parameter ( spec_class_gas = 3 )
integer spec_class_other
parameter ( spec_class_other = 4 )
! threshold for reporting negatives from subr qneg3
real(r8) :: qneg3_worst_thresh_amode(pcnst)
integer, private :: qqcw(pcnst)=-1 ! Remaps modal_aero indices into pbuf
contains
subroutine qqcw_set_ptr(index, iptr)
use cam_abortutils, only : endrun
use time_manager, only : is_first_step
integer, intent(in) :: index, iptr
if(index>0 .and. index <= pcnst ) then
qqcw(index)=iptr
else
call endrun('attempting to set qqcw pointer already defined')
end if
end subroutine qqcw_set_ptr
function qqcw_get_field(pbuf, index, lchnk, errorhandle)
use cam_abortutils, only : endrun
use physics_buffer, only : physics_buffer_desc, pbuf_get_field
integer, intent(in) :: index, lchnk
real(r8), pointer :: qqcw_get_field(:,:)
logical, optional :: errorhandle
type(physics_buffer_desc), pointer :: pbuf(:)
logical :: error
nullify(qqcw_get_field)
error = .false.
if (index>0 .and. index <= pcnst) then
if (qqcw(index)>0) then
call pbuf_get_field(pbuf, qqcw(index), qqcw_get_field)
else
error = .true.
endif
else
error = .true.
end if
if (error .and. .not. present(errorhandle)) then
call endrun('attempt to access undefined qqcw')
end if
end function qqcw_get_field
end module modal_aero_data
!----------------------------------------------------------------
!
! maxd_aspectype = maximum allowable number of chemical species
! in each aerosol mode
!
! ntot_amode = number of aerosol modes
! ( ntot_amode_gchm = number of aerosol modes in gchm
! ntot_amode_ccm2 = number of aerosol modes to be made known to ccm2
! These are temporary until multi-mode activation scavenging is going.
! Until then, ntot_amode is set to either ntot_amode_gchm or
! ntot_amode_ccm2 depending on which code is active )
!
! msectional - if positive, moving-center sectional code is utilized,
! and each mode is actually a section.
! msectional_concinit - if positive, special code is used to initialize
! the mixing ratios of all the sections.
!
! nspec_amode(m) = number of chemical species in aerosol mode m
! nspec_amode_ccm2(m) = . . . while in ccm2 code
! nspec_amode_gchm(m) = . . . while in gchm code
! nspec_amode_nontracer(m) = number of "non-tracer" chemical
! species while in gchm code
! lspectype_amode(l,m) = species type/i.d. for chemical species l
! in aerosol mode m. (1=sulfate, others to be defined)
! lmassptr_amode(l,m) = gchm r-array index for the mixing ratio
! (moles-x/mole-air) for chemical species l in aerosol mode m
! that is in clear air or interstitial air (but not in cloud water)
! lmassptrcw_amode(l,m) = gchm r-array index for the mixing ratio
! (moles-x/mole-air) for chemical species l in aerosol mode m
! that is currently bound/dissolved in cloud water
! lwaterptr_amode(m) = gchm r-array index for the mixing ratio
! (moles-water/mole-air) for water associated with aerosol mode m
! that is in clear air or interstitial air
! lkohlercptr_amode(m) = gchm r-array index for the kohler "c" parameter
! for aerosol mode m. This is defined on a per-dry-particle-mass basis:
! c = r(i,j,k,lkohlercptr_amode) * [rhodry * (4*pi/3) * rdry^3]
! numptr_amode(m) = gchm r-array index for the number mixing ratio
! (particles/mole-air) for aerosol mode m that is in clear air or
! interstitial are (but not in cloud water). If zero or negative,
! then number is not being simulated.
! ( numptr_amode_gchm(m) = same thing but for within gchm
! numptr_amode_ccm2(m) = same thing but for within ccm2
! These are temporary, to allow testing number in gchm before ccm2 )
! numptrcw_amode(m) = gchm r-array index for the number mixing ratio
! (particles/mole-air) for aerosol mode m
! that is currently bound/dissolved in cloud water
! lsfcptr_amode(m) = gchm r-array index for the surface area mixing ratio
! (cm^2/mole-air) for aerosol mode m that is in clear air or
! interstitial are (but not in cloud water). If zero or negative,
! then surface area is not being simulated.
! lsfcptrcw_amode(m) = gchm r-array index for the surface area mixing ratio
! (cm^2/mole-air) for aerosol mode m that is currently
! bound/dissolved in cloud water.
! lsigptr_amode(m) = gchm r-array index for sigmag for aerosol mode m
! that is in clear air or interstitial are (but not in cloud water).
! If zero or negative, then the constant sigmag_amode(m) is used.
! lsigptrcw_amode(m) = gchm r-array index for sigmag for aerosol mode m
! that is currently bound/dissolved in cloud water.
! If zero or negative, then the constant sigmag_amode(m) is used.
! lsigptrac_amode(m) = gchm r-array index for sigmag for aerosol mode m
! for combined clear-air/interstial plus bound/dissolved in cloud water.
! If zero or negative, then the constant sigmag_amode(m) is used.
!
! dgnum_amode(m) = geometric dry mean diameter (m) of the number
! distribution for aerosol mode m.
! (Only used when numptr_amode(m) is zero or negative.)
! dgnumlo_amode(m), dgnumhi_amode(m) = lower and upper limits on the
! geometric dry mean diameter (m) of the number distribution
! (Used when mprognum_amode>0, to limit dgnum to reasonable values)
! sigmag_amode(m) = geometric standard deviation for aerosol mode m
! sigmaglo_amode(m), sigmaghi_amode(m) = lower and upper limits on the
! geometric standard deviation of the number distribution
! (Used when mprogsfc_amode>0, to limit sigmag to reasonable values)
! alnsg_amode(m) = alog( sigmag_amode(m) )
! alnsglo_amode(m), alnsghi_amode(m) = alog( sigmaglo/hi_amode(m) )
! voltonumb_amode(m) = ratio of number to volume for mode m
! voltonumblo_amode(m), voltonumbhi_amode(m) = ratio of number to volume
! when dgnum = dgnumlo_amode or dgnumhi_amode, respectively
! voltosfc_amode(m), voltosfclo_amode(m), voltosfchi_amode(m) - ratio of
! surface to volume for mode m (like the voltonumb_amode's)
! alnv2n_amode(m), alnv2nlo_amode(m), alnv2nhi_amode(m) -
! alnv2n_amode(m) = alog( voltonumblo_amode(m) ), ...
! alnv2s_amode(m), alnv2slo_amode(m), alnv2shi_amode(m) -
! alnv2s_amode(m) = alog( voltosfclo_amode(m) ), ...
! rhcrystal_amode(m) = crystalization r.h. for mode m
! rhdeliques_amode(m) = deliquescence r.h. for mode m
! (*** these r.h. values are 0-1 fractions, not 0-100 percentages)
!
! mcalcwater_amode(m) - if positive, water content for mode m will be
! calculated and stored in rclm(k,lwaterptr_amode(m)). Otherwise, no.
! mprognum_amode(m) - if positive, number mixing-ratio for mode m will
! be prognosed. Otherwise, no.
! mdiagnum_amode(m) - if positive, number mixing-ratio for mode m will
! be diagnosed and put into rclm(k,numptr_amode(m)). Otherwise, no.
! mprogsfc_amode(m) - if positive, surface area mixing-ratio for mode m will
! be prognosed, and sigmag will vary temporally and spatially.
! Otherwise, sigmag is constant.
! *** currently surface area is not prognosed when msectional>0 ***
!
! ntot_aspectype = overall number of aerosol chemical species defined (over all modes)
! specdens_amode(l) = dry density (kg/m^3) of aerosol chemical species type l
! specmw_amode(l) = molecular weight (kg/kmol) of aerosol chemical species type l
! specname_amode(l) = name of aerosol chemical species type l
! specrefndxsw(l) = complex refractive index (visible wavelengths)
! of aerosol chemical species type l
! specrefndxlw(l) = complex refractive index (infrared wavelengths)
! of aerosol chemical species type l
! spechygro(l) = hygroscopicity of aerosol chemical species type l
!
! lptr_so4_a_amode(m), lptr_so4_cw_amode(m) = gchm r-array index for the
! mixing ratio for sulfate associated with aerosol mode m
! ("a" and "cw" phases)
! (similar for msa, oc, bc, nacl, dust)
!
! modename_amode(m) = character-variable name for mode m,
! read from mirage2.inp
! modeptr_accum - mode index for the main accumulation mode
! if modeptr_accum = 1, then mode 1 is the main accumulation mode,
! and modename_amode(1) = "accum"
! modeptr_aitken - mode index for the main aitken mode
! if modeptr_aitken = 2, then mode 2 is the main aitken mode,
! and modename_amode(2) = "aitken"
! modeptr_ufine - mode index for the ultrafine mode
! if modeptr_ufine = 3, then mode 3 is the ultrafine mode,
! and modename_amode(3) = "ufine"
! modeptr_coarseas - mode index for the coarse sea-salt mode
! if modeptr_coarseas = 4, then mode 4 is the coarse sea-salt mode,
! and modename_amode(4) = "coarse_seasalt"
! modeptr_coardust - mode index for the coarse dust mode
! if modeptr_coardust = 5, then mode 5 is the coarse dust mode,
! and modename_amode(5) = "coarse_dust"
!
! specdens_XX_amode = dry density (kg/m^3) of aerosol chemical species type XX
! where XX is so4, om, bc, dust, seasalt
! contains same values as the specdens_amode array
! allows values to be referenced differently
! specmw_XX_amode = molecular weight (kg/kmol) of aerosol chemical species type XX
! contains same values as the specmw_amode array
!
!-----------------------------------------------------------------------
!--------------------------------------------------------------
!
! ... aerosol size information for the current chunk
!
!--------------------------------------------------------------
!
! dgncur = current geometric mean diameters (cm) for number distributions
! dgncur_a - for unactivated particles, dry
! (in physics buffer as DGNUM)
! dgncur_awet - for unactivated particles, wet at grid-cell ambient RH
! (in physics buffer as DGNUMWET)
!
! the dgncur are computed from current mass and number
! mixing ratios in the grid cell, BUT are then adjusted to be within
! the bounds defined by dgnumlo/hi_amode
!
! v2ncur = current (number/volume) ratio based on dgncur and sgcur
! (volume in cm^3/whatever, number in particles/whatever)
! == 1.0 / ( pi/6 * dgncur**3 * exp(4.5*((log(sgcur))**2)) )
! v2ncur_a - for unactivated particles
! (currently just defined locally)
!