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