PartMC  2.4.0
aero_data.F90
Go to the documentation of this file.
1 ! Copyright (C) 2005-2012, 2016 Nicole Riemer and Matthew West
2 ! Licensed under the GNU General Public License version 2 or (at your
3 ! option) any later version. See the file COPYING for details.
4 
5 !> \file
6 !> The pmc_aero_data module.
7 
8 !> The aero_data_t structure and associated subroutines.
10 
11  use pmc_spec_file
12  use pmc_mpi
13  use pmc_util
14  use pmc_fractal
15  use pmc_netcdf
16 #ifdef PMC_USE_MPI
17  use mpi
18 #endif
19 
20  integer, parameter :: aero_name_len = 15
21  integer, parameter :: aero_source_name_len = 100
22 
23  !> Aerosol material properties and associated data.
24  !!
25  !! The data in this structure is constant, as it represents physical
26  !! quantities that cannot change over time.
27  !!
28  !! Each aerosol species is identified by an index <tt>i =
29  !! 1,...,aero_data_n_spec(aero_data)</tt>. Then \c name(i) is the name of
30  !! that species, \c density(i) is its density, etc. The ordering of the
31  !! species is arbitrary and should not be relied upon (currently it is the
32  !! order in the species data file). The only exception is that it is
33  !! possible to find out which species is water from the \c i_water
34  !! variable.
35  !!
36  !! The names of the aerosol species are not important to PartMC, as
37  !! only the material properties are used. The names are used for
38  !! input and output, and also for communication with MOSAIC. For the
39  !! MOSAIC interface to work correctly the species must be named the
40  !! same, but without the \c _a suffix.
42  !> Water species number (0 if water is not a species).
43  integer :: i_water
44  !> Length \c aero_data_n_spec(aero_data), species names.
45  character(len=AERO_NAME_LEN), allocatable :: name(:)
46  !> Length \c aero_data_n_spec(aero_data), mosaic_index(i) a positive
47  !> integer giving the mosaic index of species i, or 0 if there is no match.
48  integer, allocatable :: mosaic_index(:)
49  !> Length \c aero_data_n_spec(aero_data), densities (kg m^{-3}).
50  real(kind=dp), allocatable :: density(:)
51  !> Length \c aero_data_n_spec(aero_data), number of ions in solute.
52  integer, allocatable :: num_ions(:)
53  !> Length \c aero_data_n_spec(aero_data), molecular weights (kg mol^{-1}).
54  real(kind=dp), allocatable :: molec_weight(:)
55  !> Length \c aero_data_n_spec(aero_data), kappas (1).
56  real(kind=dp), allocatable :: kappa(:)
57  !> Length \c aero_data_n_source(aero_data), source names.
58  character(len=AERO_SOURCE_NAME_LEN), allocatable :: source_name(:)
59  !> Fractal particle parameters.
60  type(fractal_t) :: fractal
61  end type aero_data_t
62 
63 contains
64 
65 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
66 
67  !> Convert mass-equivalent volume \f$V\f$ (m^3) to geometric radius
68  !> \f$R_{\rm geo}\f$ (m).
69  real(kind=dp) elemental function aero_data_vol2rad(aero_data, v)
70 
71  !> Aero data structure.
72  type(aero_data_t), intent(in) :: aero_data
73  !> Particle mass-equivalent volume (m^3).
74  real(kind=dp), intent(in) :: v
75 
76  aero_data_vol2rad = fractal_vol2rad(aero_data%fractal, v)
77 
78  end function aero_data_vol2rad
79 
80 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
81 
82  !> Convert mass-equivalent volume \f$V\f$ (m^3) to geometric diameter
83  !> \f$D_{\rm geo}\f$ (m).
84  real(kind=dp) elemental function aero_data_vol2diam(aero_data, v)
85 
86  !> Aero data structure.
87  type(aero_data_t), intent(in) :: aero_data
88  !> Particle mass-equivalent volume (m^3).
89  real(kind=dp), intent(in) :: v
90 
91  aero_data_vol2diam = fractal_vol2diam(aero_data%fractal, v)
92 
93  end function aero_data_vol2diam
94 
95 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
96 
97  !> Convert geometric radius \f$R_{\rm geo}\f$ (m) to mass-equivalent volume
98  !> \f$V\f$ (m^3).
99  real(kind=dp) elemental function aero_data_rad2vol(aero_data, r)
101  !> Aero data structure.
102  type(aero_data_t), intent(in) :: aero_data
103  !> Geometric radius (m).
104  real(kind=dp), intent(in) :: r
105 
106  aero_data_rad2vol = fractal_rad2vol(aero_data%fractal, r)
107 
108  end function aero_data_rad2vol
109 
110 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
111 
112  !> Convert geometric diameter \f$D_{\rm geo}\f$ (m) to
113  !> mass-equivalent volume \f$V\f$ (m^3).
114  real(kind=dp) elemental function aero_data_diam2vol(aero_data, d)
116  !> Aero data structure.
117  type(aero_data_t), intent(in) :: aero_data
118  !> Geometric diameter (m).
119  real(kind=dp), intent(in) :: d
120 
121  aero_data_diam2vol = fractal_rad2vol(aero_data%fractal, diam2rad(d))
122 
123  end function aero_data_diam2vol
124 
125 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
126 
127  !> Convert mass-equivalent volume \f$V\f$ (m^3) to number of
128  !> monomers \f$N\f$ in a fractal particle cluster.
129  !!
130  !! Based on Eq. 5 in Naumann [2003].
131  real(kind=dp) elemental function aero_data_vol_to_num_of_monomers( &
132  aero_data, v)
134  !> Aero data structure.
135  type(aero_data_t), intent(in) :: aero_data
136  !> Particle mass-equivalent volume (m^3).
137  real(kind=dp), intent(in) :: v
138 
140  aero_data%fractal, v)
141 
143 
144 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
145 
146  !> Convert mass-equivalent volume \f$V\f$ (m^3) to mobility equivalent
147  !> radius \f$R_{\rm me}\f$ (m).
148  !!
149  !! Based on Eq. 5, 21 and 30 in Naumann [2003].
150  real(kind=dp) function aero_data_vol_to_mobility_rad(aero_data, v, temp, &
151  pressure)
153  !> Aero data structure.
154  type(aero_data_t), intent(in) :: aero_data
155  !> Particle mass-equivalent volume (m^3).
156  real(kind=dp), intent(in) :: v
157  !> Temperature (K).
158  real(kind=dp), intent(in) :: temp
159  !> Pressure (Pa).
160  real(kind=dp), intent(in) :: pressure
161 
163  aero_data%fractal, v, temp, pressure)
164 
165  end function aero_data_vol_to_mobility_rad
166 
167 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
168 
169  !> Convert mobility equivalent radius \f$R_{\rm me}\f$ (m) to
170  !> mass-equivalent volume \f$V\f$ (m^3).
171  real(kind=dp) function aero_data_mobility_rad_to_vol(aero_data, &
172  mobility_rad, temp, pressure)
174  !> Aero data structure.
175  type(aero_data_t), intent(in) :: aero_data
176  !> Mobility equivalent radius (m).
177  real(kind=dp), intent(in) :: mobility_rad
178  !> Temperature (K).
179  real(kind=dp), intent(in) :: temp
180  !> Pressure (Pa).
181  real(kind=dp), intent(in) :: pressure
182 
184  aero_data%fractal, mobility_rad, temp, pressure)
185 
186  end function aero_data_mobility_rad_to_vol
187 
188 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
189 
190  !> Convert mobility equivalent radius \f$R_{\rm me}\f$ (m) to
191  !> geometric radius \f$R_{\rm geo}\f$ (m^3).
192  real(kind=dp) function aero_data_mobility_rad_to_geometric_rad(aero_data, &
193  mobility_rad, temp, pressure)
195  !> Aero data structure.
196  type(aero_data_t), intent(in) :: aero_data
197  !> Mobility equivalent radius (m).
198  real(kind=dp), intent(in) :: mobility_rad
199  !> Temperature (K).
200  real(kind=dp), intent(in) :: temp
201  !> Pressure (Pa).
202  real(kind=dp), intent(in) :: pressure
203 
205  = fractal_mobility_rad_to_geometric_rad(aero_data%fractal, &
206  mobility_rad, temp, pressure)
207 
209 
210 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
211 
212  !> Return the number of aerosol species, or -1 if uninitialized.
213  elemental integer function aero_data_n_spec(aero_data)
215  !> Aero data structure.
216  type(aero_data_t), intent(in) :: aero_data
217 
218  if (allocated(aero_data%name)) then
219  aero_data_n_spec = size(aero_data%name)
220  else
221  aero_data_n_spec = -1
222  end if
223 
224  end function aero_data_n_spec
225 
226 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
227 
228  !> Return the number of aerosol sources, or -1 if uninitialized.
229  elemental integer function aero_data_n_source(aero_data)
231  !> Aero data structure.
232  type(aero_data_t), intent(in) :: aero_data
233 
234  if (allocated(aero_data%source_name)) then
235  aero_data_n_source = size(aero_data%source_name)
236  else
237  aero_data_n_source = -1
238  end if
239 
240  end function aero_data_n_source
241 
242 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
243 
244  !> Returns the number of the species in aero_data with the given name, or
245  !> returns 0 if there is no such species.
246  integer function aero_data_spec_by_name(aero_data, name)
248  !> Aero_data data.
249  type(aero_data_t), intent(in) :: aero_data
250  !> Name of species to find.
251  character(len=*), intent(in) :: name
252 
253  integer i
254  logical found
255 
256  found = .false.
257  do i = 1,aero_data_n_spec(aero_data)
258  if (name == aero_data%name(i)) then
259  found = .true.
260  exit
261  end if
262  end do
263  if (found) then
265  else
267  end if
268 
269  end function aero_data_spec_by_name
270 
271 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
272 
273  !> Returns the number of the source in aero_data with the given name, or
274  !> adds the source if it doesn't exist yet.
275  integer function aero_data_source_by_name(aero_data, name)
277  !> Aero_data data.
278  type(aero_data_t), intent(inout) :: aero_data
279  !> Name of source to find.
280  character(len=*), intent(in) :: name
281 
282  if (.not. allocated(aero_data%source_name)) then
283  aero_data%source_name = [name(1:aero_source_name_len)]
285  return
286  end if
287  aero_data_source_by_name = string_array_find(aero_data%source_name, name)
288  if (aero_data_source_by_name > 0) return
289  aero_data%source_name = [aero_data%source_name, &
290  name(1:aero_source_name_len)]
291  aero_data_source_by_name = size(aero_data%source_name)
292 
293  end function aero_data_source_by_name
294 
295 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
296 
297  !> Fills in aero_data%%i_water.
298  subroutine aero_data_set_water_index(aero_data)
300  !> Aero_data data.
301  type(aero_data_t), intent(inout) :: aero_data
302 
303  integer :: i
304 
305  aero_data%i_water = string_array_find(aero_data%name, "H2O")
306 
307  end subroutine aero_data_set_water_index
308 
309 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
310 
311  !> Fills in aero_data%%mosaic_index.
312  subroutine aero_data_set_mosaic_map(aero_data)
314  !> Aero_data data.
315  type(aero_data_t), intent(inout) :: aero_data
316 
317  integer, parameter :: n_mosaic_spec = 19
318  character(AERO_NAME_LEN), parameter, dimension(n_mosaic_spec) :: &
319  mosaic_spec_name = [ &
320  "SO4 ", "NO3 ", "Cl ", "NH4 ", "MSA ", "ARO1 ", &
321  "ARO2 ", "ALK1 ", "OLE1 ", "API1 ", "API2 ", "LIM1 ", &
322  "LIM2 ", "CO3 ", "Na ", "Ca ", "OIN ", "OC ", &
323  "BC "]
324 
325  integer :: i_spec, i_mosaic_spec, i
326 
327  aero_data%mosaic_index = 0
328  do i_spec = 1,aero_data_n_spec(aero_data)
329  aero_data%mosaic_index(i_spec) = string_array_find(mosaic_spec_name, &
330  aero_data%name(i_spec))
331  end do
332 
333  end subroutine aero_data_set_mosaic_map
334 
335 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
336 
337  !> Read aero_data specification from a spec file.
338  subroutine spec_file_read_aero_data(file, aero_data)
339 
340  !> Spec file to read data from.
341  type(spec_file_t), intent(inout) :: file
342  !> Aero_data data.
343  type(aero_data_t), intent(inout) :: aero_data
344 
345  integer :: n_species, species, i
346  character(len=SPEC_LINE_MAX_VAR_LEN), allocatable :: species_name(:)
347  real(kind=dp), allocatable :: species_data(:,:)
348 
349  !> \page input_format_aero_data Input File Format: Aerosol Material Data
350  !!
351  !! A aerosol material data file must consist of one line per
352  !! aerosol species, with each line having:
353  !! - species name (string)
354  !! - density (real, unit kg/m^3)
355  !! - ions per fully dissociated molecule (integer) - used to
356  !! compute kappa value if the corresponding kappa value is
357  !! zero
358  !! - molecular weight (real, unit kg/mol)
359  !! - kappa hygroscopicity parameter (real, dimensionless) - if
360  !! zero, then inferred from the ions value
361  !!
362  !! This specifies both which species are to be recognized as
363  !! aerosol consituents, as well as their physical properties. For
364  !! example, an aerosol material data file could contain:
365  !! <pre>
366  !! # species dens (kg/m^3) ions (1) molec wght (kg/mole) kappa (1)
367  !! SO4 1800 0 96e-3 0.65
368  !! NO3 1800 0 62e-3 0.65
369  !! Cl 2200 1 35.5e-3 0
370  !! NH4 1800 0 18e-3 0.65
371  !! </pre>
372  !!
373  !! Note that it is an error to specify a non-zero number of ions
374  !! and a non-zero kappa value for a species. If both values are
375  !! zero then that species has zero hygroscopicity parameter. If
376  !! exactly one of kappa or ions is non-zero then the non-zero
377  !! value is used and the zero value is ignored.
378  !!
379  !! See also:
380  !! - \ref spec_file_format --- the input file text format
381  !! - \ref output_format_aero_data --- the corresponding output format
382 
383  call spec_file_read_real_named_array(file, 0, species_name, species_data)
384 
385  ! check the data size
386  n_species = size(species_data, 1)
387  if (.not. ((size(species_data, 2) == 4) .or. (n_species == 0))) then
388  call die_msg(428926381, 'each line in ' // trim(file%name) &
389  // ' should contain exactly 5 values')
390  end if
391 
392  ! allocate and copy over the data
393  call ensure_string_array_size(aero_data%name, n_species)
394  call ensure_integer_array_size(aero_data%mosaic_index, n_species)
395  call ensure_real_array_size(aero_data%density, n_species)
396  call ensure_integer_array_size(aero_data%num_ions, n_species)
397  call ensure_real_array_size(aero_data%molec_weight, n_species)
398  call ensure_real_array_size(aero_data%kappa, n_species)
399  do i = 1,n_species
400  aero_data%name(i) = species_name(i)(1:aero_name_len)
401  aero_data%density(i) = species_data(i,1)
402  aero_data%num_ions(i) = nint(species_data(i,2))
403  aero_data%molec_weight(i) = species_data(i,3)
404  aero_data%kappa(i) = species_data(i,4)
405  call assert_msg(232362742, &
406  (aero_data%num_ions(i) == 0) .or. (aero_data%kappa(i) == 0d0), &
407  "ions and kappa both non-zero for species " &
408  // trim(aero_data%name(i)) // " in " // trim(file%name))
409  if (species_name(i) == "H2O") then
410  aero_data%i_water = i
411  call warn_assert_msg(945800387, &
412  aero_data%density(i) == const%water_density, &
413  "input H2O density not equal to const%water_density (" &
414  // trim(real_to_string(aero_data%density(i))) // " /= " &
415  // trim(real_to_string(const%water_density)) // ")")
416  call warn_assert_msg(233517437, &
417  aero_data%molec_weight(i) == const%water_molec_weight, &
418  "input H2O molec_weight not equal " &
419  // "to const%water_molec_weight (" &
420  // trim(real_to_string(aero_data%molec_weight(i))) // " /= " &
421  // trim(real_to_string(const%water_molec_weight)) // ")")
422  end if
423  end do
424  call aero_data_set_water_index(aero_data)
425  call aero_data_set_mosaic_map(aero_data)
426 
427  end subroutine spec_file_read_aero_data
428 
429 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
430 
431  !> Read a list of species from the given file with the given name.
432  subroutine spec_file_read_species_list(file, name, aero_data, species_list)
434  !> Spec file.
435  type(spec_file_t), intent(inout) :: file
436  !> Name of line.
437  character(len=*), intent(in) :: name
438  !> Aero_data data.
439  type(aero_data_t), intent(in) :: aero_data
440  !> List of species numbers.
441  integer, allocatable :: species_list(:)
442 
443  type(spec_line_t) :: line
444  integer :: i, spec
445 
446  call spec_file_read_line_no_eof(file, line)
447  call spec_file_check_line_name(file, line, name)
448  call ensure_integer_array_size(species_list, size(line%data))
449  do i = 1,size(line%data)
450  spec = aero_data_spec_by_name(aero_data, line%data(i))
451  if (spec == 0) then
452  call spec_file_die_msg(964771833, file, &
453  'unknown species: ' // trim(line%data(i)))
454  end if
455  species_list(i) = spec
456  end do
457 
458  end subroutine spec_file_read_species_list
459 
460 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
461 
462  !> Determines the number of bytes required to pack the given value.
463  integer function pmc_mpi_pack_size_aero_data(val)
465  !> Value to pack.
466  type(aero_data_t), intent(in) :: val
467 
469  pmc_mpi_pack_size_integer(val%i_water) &
470  + pmc_mpi_pack_size_string_array(val%name) &
471  + pmc_mpi_pack_size_integer_array(val%mosaic_index) &
472  + pmc_mpi_pack_size_real_array(val%density) &
473  + pmc_mpi_pack_size_integer_array(val%num_ions) &
474  + pmc_mpi_pack_size_real_array(val%molec_weight) &
475  + pmc_mpi_pack_size_real_array(val%kappa) &
476  + pmc_mpi_pack_size_string_array(val%source_name) &
477  + pmc_mpi_pack_size_fractal(val%fractal)
478 
479  end function pmc_mpi_pack_size_aero_data
480 
481 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
482 
483  !> Packs the given value into the buffer, advancing position.
484  subroutine pmc_mpi_pack_aero_data(buffer, position, val)
486  !> Memory buffer.
487  character, intent(inout) :: buffer(:)
488  !> Current buffer position.
489  integer, intent(inout) :: position
490  !> Value to pack.
491  type(aero_data_t), intent(in) :: val
492 
493 #ifdef PMC_USE_MPI
494  integer :: prev_position
495 
496  prev_position = position
497  call pmc_mpi_pack_integer(buffer, position, val%i_water)
498  call pmc_mpi_pack_string_array(buffer, position, val%name)
499  call pmc_mpi_pack_integer_array(buffer, position, val%mosaic_index)
500  call pmc_mpi_pack_real_array(buffer, position, val%density)
501  call pmc_mpi_pack_integer_array(buffer, position, val%num_ions)
502  call pmc_mpi_pack_real_array(buffer, position, val%molec_weight)
503  call pmc_mpi_pack_real_array(buffer, position, val%kappa)
504  call pmc_mpi_pack_string_array(buffer, position, val%source_name)
505  call pmc_mpi_pack_fractal(buffer, position, val%fractal)
506  call assert(183834856, &
507  position - prev_position <= pmc_mpi_pack_size_aero_data(val))
508 #endif
509 
510  end subroutine pmc_mpi_pack_aero_data
511 
512 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
513 
514  !> Unpacks the given value from the buffer, advancing position.
515  subroutine pmc_mpi_unpack_aero_data(buffer, position, val)
517  !> Memory buffer.
518  character, intent(inout) :: buffer(:)
519  !> Current buffer position.
520  integer, intent(inout) :: position
521  !> Value to pack.
522  type(aero_data_t), intent(inout) :: val
523 
524 #ifdef PMC_USE_MPI
525  integer :: prev_position
526 
527  prev_position = position
528  call pmc_mpi_unpack_integer(buffer, position, val%i_water)
529  call pmc_mpi_unpack_string_array(buffer, position, val%name)
530  call pmc_mpi_unpack_integer_array(buffer, position, val%mosaic_index)
531  call pmc_mpi_unpack_real_array(buffer, position, val%density)
532  call pmc_mpi_unpack_integer_array(buffer, position, val%num_ions)
533  call pmc_mpi_unpack_real_array(buffer, position, val%molec_weight)
534  call pmc_mpi_unpack_real_array(buffer, position, val%kappa)
535  call pmc_mpi_unpack_string_array(buffer, position, val%source_name)
536  call pmc_mpi_unpack_fractal(buffer, position, val%fractal)
537  call assert(188522823, &
538  position - prev_position <= pmc_mpi_pack_size_aero_data(val))
539 #endif
540 
541  end subroutine pmc_mpi_unpack_aero_data
542 
543 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
544 
545  !> Write the aero species dimension to the given NetCDF file if it
546  !> is not already present and in any case return the associated
547  !> dimid.
548  subroutine aero_data_netcdf_dim_aero_species(aero_data, ncid, &
549  dimid_aero_species)
551  !> Aero_data structure.
552  type(aero_data_t), intent(in) :: aero_data
553  !> NetCDF file ID, in data mode.
554  integer, intent(in) :: ncid
555  !> Dimid of the species dimension.
556  integer, intent(out) :: dimid_aero_species
557 
558  integer :: status, i_spec
559  integer :: varid_aero_species
560  integer :: aero_species_centers(aero_data_n_spec(aero_data))
561  character(len=(AERO_NAME_LEN * aero_data_n_spec(aero_data))) :: &
562  aero_species_names
563 
564  ! try to get the dimension ID
565  status = nf90_inq_dimid(ncid, "aero_species", dimid_aero_species)
566  if (status == nf90_noerr) return
567  if (status /= nf90_ebaddim) call pmc_nc_check(status)
568 
569  ! dimension not defined, so define now define it
570  call pmc_nc_check(nf90_redef(ncid))
571 
572  call pmc_nc_check(nf90_def_dim(ncid, "aero_species", &
573  aero_data_n_spec(aero_data), dimid_aero_species))
574  aero_species_names = ""
575  do i_spec = 1,aero_data_n_spec(aero_data)
576  aero_species_names((len_trim(aero_species_names) + 1):) &
577  = trim(aero_data%name(i_spec))
578  if (i_spec < aero_data_n_spec(aero_data)) then
579  aero_species_names((len_trim(aero_species_names) + 1):) = ","
580  end if
581  end do
582  call pmc_nc_check(nf90_def_var(ncid, "aero_species", nf90_int, &
583  dimid_aero_species, varid_aero_species))
584  call pmc_nc_check(nf90_put_att(ncid, varid_aero_species, "names", &
585  aero_species_names))
586  call pmc_nc_check(nf90_put_att(ncid, varid_aero_species, "description", &
587  "dummy dimension variable (no useful value) - read species names " &
588  // "as comma-separated values from the 'names' attribute"))
589 
590  call pmc_nc_check(nf90_enddef(ncid))
591 
592  do i_spec = 1,aero_data_n_spec(aero_data)
593  aero_species_centers(i_spec) = i_spec
594  end do
595  call pmc_nc_check(nf90_put_var(ncid, varid_aero_species, &
596  aero_species_centers))
597 
598  end subroutine aero_data_netcdf_dim_aero_species
599 
600 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
601 
602  !> Write the aero source dimension to the given NetCDF file if it
603  !> is not already present and in any case return the associated
604  !> dimid.
605  subroutine aero_data_netcdf_dim_aero_source(aero_data, ncid, &
606  dimid_aero_source)
608  !> Aero_data structure.
609  type(aero_data_t), intent(in) :: aero_data
610  !> NetCDF file ID, in data mode.
611  integer, intent(in) :: ncid
612  !> Dimid of the source dimension.
613  integer, intent(out) :: dimid_aero_source
614 
615  integer :: status, i_source
616  integer :: varid_aero_source
617  integer :: aero_source_centers(aero_data_n_source(aero_data))
618  character(len=(AERO_SOURCE_NAME_LEN * aero_data_n_source(aero_data))) &
619  :: aero_source_names
620 
621  ! try to get the dimension ID
622  status = nf90_inq_dimid(ncid, "aero_source", dimid_aero_source)
623  if (status == nf90_noerr) return
624  if (status /= nf90_ebaddim) call pmc_nc_check(status)
625 
626  ! dimension not defined, so define now define it
627  call pmc_nc_check(nf90_redef(ncid))
628 
629  call pmc_nc_check(nf90_def_dim(ncid, "aero_source", &
630  aero_data_n_source(aero_data), dimid_aero_source))
631  aero_source_names = ""
632  do i_source = 1,aero_data_n_source(aero_data)
633  aero_source_names((len_trim(aero_source_names) + 1):) &
634  = trim(aero_data%source_name(i_source))
635  if (i_source < aero_data_n_source(aero_data)) then
636  aero_source_names((len_trim(aero_source_names) + 1):) = ","
637  end if
638  end do
639  call pmc_nc_check(nf90_def_var(ncid, "aero_source", nf90_int, &
640  dimid_aero_source, varid_aero_source))
641  call pmc_nc_check(nf90_put_att(ncid, varid_aero_source, "names", &
642  aero_source_names))
643  call pmc_nc_check(nf90_put_att(ncid, varid_aero_source, "description", &
644  "dummy dimension variable (no useful value) - read source names " &
645  // "as comma-separated values from the 'names' attribute"))
646 
647  call pmc_nc_check(nf90_enddef(ncid))
648 
649  do i_source = 1,aero_data_n_source(aero_data)
650  aero_source_centers(i_source) = i_source
651  end do
652  call pmc_nc_check(nf90_put_var(ncid, varid_aero_source, &
653  aero_source_centers))
654 
655  end subroutine aero_data_netcdf_dim_aero_source
656 
657 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
658 
659  !> Write full state.
660  subroutine aero_data_output_netcdf(aero_data, ncid)
661 
662  !> Aero_data to write.
663  type(aero_data_t), intent(in) :: aero_data
664  !> NetCDF file ID, in data mode.
665  integer, intent(in) :: ncid
666 
667  integer :: dimid_aero_species, dimid_aero_source
668 
669  !> \page output_format_aero_data Output File Format: Aerosol Material Data
670  !!
671  !! The aerosol material data NetCDF dimensions are:
672  !! - \b aero_species: number of aerosol species
673  !! - \b aero_source: number of aerosol sources
674  !!
675  !! The aerosol material data NetCDF variables are:
676  !! - \b aero_species (dim \c aero_species): dummy dimension variable
677  !! (no useful value) - read species names as comma-separated values
678  !! from the 'names' attribute
679  !! - \b aero_source (dim \c aero_source): dummy dimension variable
680  !! (no useful value) - read source names as comma-separated values
681  !! from the 'names' attribute
682  !! - \b aero_mosaic_index (dim \c aero_species): indices of species
683  !! in MOSAIC
684  !! - \b aero_density (unit kg/m^3, dim \c aero_species): densities
685  !! of aerosol species
686  !! - \b aero_num_ions (dim \c aero_species): number of ions produced
687  !! when one molecule of each species fully dissociates in water
688  !! - \b aero_molec_weight (unit kg/mol, dim \c aero_species): molecular
689  !! weights of aerosol species
690  !! - \b aero_kappa (unit kg/mol, dim \c aero_species): hygroscopicity
691  !! parameters of aerosol species
692  !! - \b fractal parameters, see \ref output_format_fractal
693  !!
694  !! See also:
695  !! - \ref input_format_aero_data --- the corresponding input format
696 
697  call aero_data_netcdf_dim_aero_species(aero_data, ncid, &
698  dimid_aero_species)
699  call aero_data_netcdf_dim_aero_source(aero_data, ncid, &
700  dimid_aero_source)
701 
702  call pmc_nc_write_integer_1d(ncid, aero_data%mosaic_index, &
703  "aero_mosaic_index", (/ dimid_aero_species /), &
704  long_name="MOSAIC indices of aerosol species")
705  call pmc_nc_write_real_1d(ncid, aero_data%density, &
706  "aero_density", (/ dimid_aero_species /), unit="kg/m^3", &
707  long_name="densities of aerosol species")
708  call pmc_nc_write_integer_1d(ncid, aero_data%num_ions, &
709  "aero_num_ions", (/ dimid_aero_species /), &
710  long_name="number of ions after dissociation of aerosol species")
711  call pmc_nc_write_real_1d(ncid, aero_data%molec_weight, &
712  "aero_molec_weight", (/ dimid_aero_species /), unit="kg/mol", &
713  long_name="molecular weights of aerosol species")
714  call pmc_nc_write_real_1d(ncid, aero_data%kappa, &
715  "aero_kappa", (/ dimid_aero_species /), unit="1", &
716  long_name="hygroscopicity parameters (kappas) of aerosol species")
717  call fractal_output_netcdf(aero_data%fractal, ncid)
718 
719  end subroutine aero_data_output_netcdf
720 
721 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
722 
723  !> Read full state.
724  subroutine aero_data_input_netcdf(aero_data, ncid)
726  !> Aero_data to read.
727  type(aero_data_t), intent(inout) :: aero_data
728  !> NetCDF file ID, in data mode.
729  integer, intent(in) :: ncid
730 
731  integer, parameter :: MAX_SPECIES = 1000
732  integer, parameter :: MAX_SOURCES = 1000
733 
734  character(len=1000) :: name
735  integer :: dimid_aero_species, n_spec, varid_aero_species, i_spec, i
736  integer :: dimid_aero_source, n_source, varid_aero_source, i_source
737  character(len=((AERO_NAME_LEN + 2) * MAX_SPECIES)) :: aero_species_names
738  character(len=((AERO_SOURCE_NAME_LEN + 2) * MAX_SPECIES)) &
739  :: aero_source_names
740 
741  call pmc_nc_check(nf90_inq_dimid(ncid, "aero_species", &
742  dimid_aero_species))
743  call pmc_nc_check(nf90_inquire_dimension(ncid, &
744  dimid_aero_species, name, n_spec))
745  call assert(141013948, n_spec < max_species)
746 
747  call pmc_nc_check(nf90_inq_dimid(ncid, "aero_source", &
748  dimid_aero_source))
749  call pmc_nc_check(nf90_inquire_dimension(ncid, &
750  dimid_aero_source, name, n_source))
751  call assert(739238793, n_source < max_sources)
752 
753  call pmc_nc_read_integer_1d(ncid, aero_data%mosaic_index, &
754  "aero_mosaic_index")
755  call pmc_nc_read_real_1d(ncid, aero_data%density, "aero_density")
756  call pmc_nc_read_integer_1d(ncid, aero_data%num_ions, "aero_num_ions")
757  call pmc_nc_read_real_1d(ncid, aero_data%molec_weight, "aero_molec_weight")
758  call pmc_nc_read_real_1d(ncid, aero_data%kappa, "aero_kappa")
759 
760  call pmc_nc_check(nf90_inq_varid(ncid, "aero_species", &
761  varid_aero_species))
762  call pmc_nc_check(nf90_get_att(ncid, varid_aero_species, "names", &
763  aero_species_names))
764  ! aero_species_names are comma-separated, so unpack them
765  call ensure_string_array_size(aero_data%name, n_spec)
766  do i_spec = 1,aero_data_n_spec(aero_data)
767  i = 1
768  do while ((aero_species_names(i:i) /= " ") &
769  .and. (aero_species_names(i:i) /= ","))
770  i = i + 1
771  end do
772  call assert(852937292, i > 1)
773  aero_data%name(i_spec) = aero_species_names(1:(i-1))
774  aero_species_names = aero_species_names((i+1):)
775  end do
776  call assert(729138192, aero_species_names == "")
777 
778  call pmc_nc_check(nf90_inq_varid(ncid, "aero_source", &
779  varid_aero_source))
780  call pmc_nc_check(nf90_get_att(ncid, varid_aero_source, "names", &
781  aero_source_names))
782  ! aero_source_names are comma-separated, so unpack them
783  call ensure_string_array_size(aero_data%source_name, n_source)
784  do i_source = 1,aero_data_n_source(aero_data)
785  i = 1
786  do while ((aero_source_names(i:i) /= " ") &
787  .and. (aero_source_names(i:i) /= ","))
788  i = i + 1
789  end do
790  call assert(840982478, i > 1)
791  aero_data%source_name(i_source) = aero_source_names(1:(i-1))
792  aero_source_names = aero_source_names((i+1):)
793  end do
794  call assert(377166446, aero_source_names == "")
795 
796  call aero_data_set_water_index(aero_data)
797 
798  call fractal_input_netcdf(aero_data%fractal, ncid)
799 
800  end subroutine aero_data_input_netcdf
801 
802 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
803 
804 end module pmc_aero_data
subroutine pmc_mpi_pack_integer_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:688
real(kind=dp) elemental function fractal_vol_to_num_of_monomers(fractal, v)
Convert mass-equivalent volume (m^3) to number of monomers in a fractal particle cluster...
Definition: fractal.F90:91
type(fractal_t), intent(in) fractal
Definition: fractal.F90:614
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
Definition: util.F90:77
An input file with extra data for printing messages.
Definition: spec_file.F90:59
integer function pmc_mpi_pack_size_integer_array(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:447
subroutine spec_file_read_species_list(file, name, aero_data, species_list)
Read a list of species from the given file with the given name.
Definition: aero_data.F90:433
Wrapper functions for NetCDF. These all take a NetCDF ncid in data mode and return with it again in d...
Definition: netcdf.F90:11
subroutine pmc_nc_read_real_1d(ncid, var, name, must_be_present)
Read a simple real array from a NetCDF file.
Definition: netcdf.F90:219
real(kind=dp) function aero_data_vol_to_mobility_rad(aero_data, v, temp, pressure)
Convert mass-equivalent volume (m^3) to mobility equivalent radius (m).
Definition: aero_data.F90:152
subroutine ensure_real_array_size(x, n, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
Definition: util.F90:1041
real(kind=dp) elemental function aero_data_diam2vol(aero_data, d)
Convert geometric diameter (m) to mass-equivalent volume (m^3).
Definition: aero_data.F90:115
subroutine pmc_mpi_unpack_integer_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:946
subroutine ensure_integer_array_size(x, n, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
Definition: util.F90:1115
real(kind=dp) elemental function aero_data_rad2vol(aero_data, r)
Convert geometric radius (m) to mass-equivalent volume (m^3).
Definition: aero_data.F90:100
subroutine pmc_mpi_pack_string_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:752
subroutine aero_data_netcdf_dim_aero_species(aero_data, ncid, dimid_aero_species)
Write the aero species dimension to the given NetCDF file if it is not already present and in any cas...
Definition: aero_data.F90:550
integer function pmc_mpi_pack_size_real_array(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:477
The fractal_t structure and associated subroutines.
Definition: fractal.F90:35
Fractal data.
Definition: fractal.F90:59
subroutine aero_data_set_mosaic_map(aero_data)
Fills in aero_data%mosaic_index.
Definition: aero_data.F90:313
real(kind=dp) function aero_data_mobility_rad_to_vol(aero_data, mobility_rad, temp, pressure)
Convert mobility equivalent radius (m) to mass-equivalent volume (m^3).
Definition: aero_data.F90:173
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition: util.F90:103
real(kind=dp) elemental function fractal_rad2vol(fractal, r)
Convert geometric radius (m) to mass-equivalent volume (m^3).
Definition: fractal.F90:143
integer function aero_data_source_by_name(aero_data, name)
Returns the number of the source in aero_data with the given name, or adds the source if it doesn&#39;t e...
Definition: aero_data.F90:276
subroutine pmc_mpi_pack_fractal(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: fractal.F90:505
subroutine pmc_nc_check(status)
Check the status of a NetCDF function call.
Definition: netcdf.F90:23
subroutine warn_assert_msg(code, condition_ok, warning_msg)
Prints a warning message if condition_ok is false.
Definition: util.F90:59
subroutine die_msg(code, error_msg)
Error immediately.
Definition: util.F90:134
subroutine pmc_mpi_unpack_aero_data(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: aero_data.F90:516
real(kind=dp) elemental function aero_data_vol_to_num_of_monomers(aero_data, v)
Convert mass-equivalent volume (m^3) to number of monomers in a fractal particle cluster...
Definition: aero_data.F90:133
subroutine ensure_string_array_size(x, n)
Allocate or reallocate the given array to ensure it is of the given size, without preserving data...
Definition: util.F90:1189
subroutine aero_data_input_netcdf(aero_data, ncid)
Read full state.
Definition: aero_data.F90:725
subroutine spec_file_die_msg(code, file, msg)
Exit with an error message containing filename and line number.
Definition: spec_file.F90:74
subroutine pmc_mpi_pack_aero_data(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: aero_data.F90:485
character(len=pmc_util_convert_string_len) function real_to_string(val)
Convert a real to a string format.
Definition: util.F90:782
integer function pmc_mpi_pack_size_fractal(val)
Determines the number of bytes required to pack the given value.
Definition: fractal.F90:490
The aero_data_t structure and associated subroutines.
Definition: aero_data.F90:9
real(kind=dp) elemental function fractal_vol2rad(fractal, v)
Convert mass-equivalent volume (m^3) to geometric radius (m).
Definition: fractal.F90:107
real(kind=dp) elemental function diam2rad(d)
Convert diameter (m) to radius (m).
Definition: util.F90:277
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:561
integer, parameter aero_name_len
Definition: aero_data.F90:20
subroutine pmc_mpi_unpack_fractal(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: fractal.F90:528
Reading formatted text input.
Definition: spec_file.F90:43
subroutine pmc_mpi_unpack_string_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:1010
integer function string_array_find(array, val)
Return the index of the first occurance of the given value in the array, or 0 if it is not present...
Definition: util.F90:1023
subroutine pmc_mpi_pack_real_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:720
subroutine pmc_mpi_unpack_real_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:978
real(kind=dp) function fractal_mobility_rad_to_vol(fractal, mobility_rad, temp, pressure)
Convert mobility equivalent radius (m) to mass-equivalent volume (m^3).
Definition: fractal.F90:468
subroutine aero_data_netcdf_dim_aero_source(aero_data, ncid, dimid_aero_source)
Write the aero source dimension to the given NetCDF file if it is not already present and in any case...
Definition: aero_data.F90:607
subroutine spec_file_read_real_named_array(file, max_lines, names, vals)
Read an array of named lines with real data. All lines must have the same number of data elements...
Definition: spec_file.F90:649
real(kind=dp) elemental function aero_data_vol2rad(aero_data, v)
Convert mass-equivalent volume (m^3) to geometric radius (m).
Definition: aero_data.F90:70
subroutine pmc_nc_write_real_1d(ncid, var, name, dimids, dim_name, unit, long_name, standard_name, description)
Write a simple real array to a NetCDF file.
Definition: netcdf.F90:536
subroutine spec_file_check_line_name(file, line, name)
Check that the name of the line data is as given.
Definition: spec_file.F90:389
subroutine aero_data_set_water_index(aero_data)
Fills in aero_data%i_water.
Definition: aero_data.F90:299
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:345
real(kind=dp) elemental function aero_data_vol2diam(aero_data, v)
Convert mass-equivalent volume (m^3) to geometric diameter (m).
Definition: aero_data.F90:85
subroutine spec_file_read_line_no_eof(file, line)
Read a spec_line from the spec_file. This will always succeed or error out, so should only be called ...
Definition: spec_file.F90:298
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:818
real(kind=dp) elemental function fractal_vol2diam(fractal, v)
Convert mass-equivalent volume (m^3) to geometric diameter (m).
Definition: fractal.F90:128
integer function pmc_mpi_pack_size_aero_data(val)
Determines the number of bytes required to pack the given value.
Definition: aero_data.F90:464
elemental integer function aero_data_n_source(aero_data)
Return the number of aerosol sources, or -1 if uninitialized.
Definition: aero_data.F90:230
elemental integer function aero_data_n_spec(aero_data)
Return the number of aerosol species, or -1 if uninitialized.
Definition: aero_data.F90:214
subroutine pmc_nc_read_integer_1d(ncid, var, name, must_be_present)
Read a simple integer array from a NetCDF file.
Definition: netcdf.F90:260
Aerosol material properties and associated data.
Definition: aero_data.F90:41
integer, parameter aero_source_name_len
Definition: aero_data.F90:21
subroutine fractal_input_netcdf(fractal, ncid)
Read full state.
Definition: fractal.F90:633
real(kind=dp) function fractal_mobility_rad_to_geometric_rad(fractal, mobility_rad, temp, pressure)
Convert mobility equivalent radius (m) to geometric radius (m^3).
Definition: fractal.F90:436
Common utility subroutines.
Definition: util.F90:9
Wrapper functions for MPI.
Definition: mpi.F90:13
subroutine pmc_nc_write_integer_1d(ncid, var, name, dimids, dim_name, unit, long_name, standard_name, description)
Write a simple integer array to a NetCDF file.
Definition: netcdf.F90:585
real(kind=dp) function aero_data_mobility_rad_to_geometric_rad(aero_data, mobility_rad, temp, pressure)
Convert mobility equivalent radius (m) to geometric radius (m^3).
Definition: aero_data.F90:194
integer function aero_data_spec_by_name(aero_data, name)
Returns the number of the species in aero_data with the given name, or returns 0 if there is no such ...
Definition: aero_data.F90:247
integer function pmc_mpi_pack_size_string_array(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:507
real(kind=dp) function fractal_vol_to_mobility_rad(fractal, v, temp, pressure)
Convert mass-equivalent volume (m^3) to mobility equivalent radius (m).
Definition: fractal.F90:317