PartMC  2.3.0
aero_weight_array.F90
Go to the documentation of this file.
1 ! Copyright (C) 2012 Jeff Curtis
2 ! Copyright (C) 2012-2015 Matthew West
3 ! Licensed under the GNU General Public License version 2 or (at your
4 ! option) any later version. See the file COPYING for details.
5 
6 !> \file
7 !> The pmc_aero_weight_array module.
8 
9 !> The aero_weight_array_t structure and associated subroutines.
11 
12  use pmc_util
13  use pmc_constants
14  use pmc_rand
15  use pmc_spec_file
17  use pmc_netcdf
18  use pmc_mpi
19  use pmc_aero_weight
20 #ifdef PMC_USE_MPI
21  use mpi
22 #endif
23 
24  !> An array of aerosol size distribution weighting functions.
25  !!
26  !! Each particle has a weight group and a weight class, which
27  !! defines the weighting function associated with it. To determine
28  !! the actual weight for a particle, the harmonic mean of the weight
29  !! from each group is computed (for the given class). All particles
30  !! of a given size within a weight class will thus have the same
31  !! weight, irrespective of which group they are in. The group is
32  !! only important when doubling or halving occurs, which takes place
33  !! for a specific group and class.
35  !> Aero weight array, <tt>n_group x n_class</tt>.
36  type(aero_weight_t), pointer :: weight(:, :)
37  end type aero_weight_array_t
38 
39 contains
40 
41 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
42 
43  !> Allocates an aero_weight_array.
44  subroutine aero_weight_array_allocate(aero_weight_array)
45 
46  !> Aerosol weight array.
47  type(aero_weight_array_t), intent(out) :: aero_weight_array
48 
49  allocate(aero_weight_array%weight(0, 0))
50 
51  end subroutine aero_weight_array_allocate
52 
53 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
54 
55  !> Allocates an \c aero_weight_array to the given size.
56  subroutine aero_weight_array_allocate_size(aero_weight_array, n_group, &
57  n_class)
58 
59  !> Aerosol weight array.
60  type(aero_weight_array_t), intent(out) :: aero_weight_array
61  !> Number of weight groups.
62  integer,intent(in) :: n_group
63  !> Number of weight classes.
64  integer,intent(in) :: n_class
65 
66  allocate(aero_weight_array%weight(n_group, n_class))
67  call aero_weight_allocate(aero_weight_array%weight)
68 
69  end subroutine aero_weight_array_allocate_size
70 
71 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
72 
73  !> Allocates an \c aero_weight_array as flat weightings.
74  subroutine aero_weight_array_allocate_flat(aero_weight_array, n_class)
75 
76  !> Aerosol weight array.
77  type(aero_weight_array_t), intent(out) :: aero_weight_array
78  !> Number of weight classes.
79  integer,intent(in) :: n_class
80 
81  call aero_weight_array_allocate_size(aero_weight_array, 1, n_class)
82  aero_weight_array%weight%type = aero_weight_type_none
83  aero_weight_array%weight%magnitude = 1d0
84  aero_weight_array%weight%exponent = 0d0
85 
86  end subroutine aero_weight_array_allocate_flat
87 
88 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
89 
90  !> Allocates an \c aero_weight_array as power weightings.
91  subroutine aero_weight_array_allocate_power(aero_weight_array, n_class, &
92  exponent)
93 
94  !> Aerosol weight array.
95  type(aero_weight_array_t), intent(out) :: aero_weight_array
96  !> Number of weight classes.
97  integer,intent(in) :: n_class
98  !> Exponent for power-law.
99  real(kind=dp), intent(in) :: exponent
100 
101  call aero_weight_array_allocate_size(aero_weight_array, 1, n_class)
102  aero_weight_array%weight%type = aero_weight_type_power
103  aero_weight_array%weight%magnitude = 1d0
104  aero_weight_array%weight%exponent = exponent
105 
106  end subroutine aero_weight_array_allocate_power
107 
108 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
109 
110  !> Allocates an \c aero_weight_array as joint flat/power-3 weightings..
111  subroutine aero_weight_array_allocate_nummass(aero_weight_array, n_class)
112 
113  !> Aerosol weight array.
114  type(aero_weight_array_t), intent(out) :: aero_weight_array
115  !> Number of weight classes.
116  integer,intent(in) :: n_class
117 
118  call aero_weight_array_allocate_size(aero_weight_array, 2, n_class)
119  aero_weight_array%weight(1, :)%type = aero_weight_type_none
120  aero_weight_array%weight(1, :)%magnitude = 1d0
121  aero_weight_array%weight(1, :)%exponent = 0d0
122  aero_weight_array%weight(2, :)%type = aero_weight_type_power
123  aero_weight_array%weight(2, :)%magnitude = 1d0
124  aero_weight_array%weight(2, :)%exponent = -3d0
125 
127 
128 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
129 
130  !> Free all storage.
131  subroutine aero_weight_array_deallocate(aero_weight_array)
132 
133  !> Structure to deallocate.
134  type(aero_weight_array_t), intent(inout) :: aero_weight_array
135 
136  call aero_weight_deallocate(aero_weight_array%weight)
137  deallocate(aero_weight_array%weight)
138 
139  end subroutine aero_weight_array_deallocate
140 
141 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
142 
143  !> Zeros the contents of the \c aero_weight_array.
144  subroutine aero_weight_array_zero(aero_weight_array)
145 
146  !> Aerosol weight array.
147  type(aero_weight_array_t), intent(inout) :: aero_weight_array
148 
149  call aero_weight_zero(aero_weight_array%weight)
150 
151  end subroutine aero_weight_array_zero
152 
153 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
154 
155  !> Normalizes the \c aero_weight_array to a non-zero value.
156  elemental subroutine aero_weight_array_normalize(aero_weight_array)
157 
158  !> Aerosol weight array.
159  type(aero_weight_array_t), intent(inout) :: aero_weight_array
160 
161  call aero_weight_normalize(aero_weight_array%weight)
162 
163  end subroutine aero_weight_array_normalize
164 
165 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
166 
167  !> Copy an \c aero_weight_array.
168  subroutine aero_weight_array_copy(aero_weight_array_from, &
169  aero_weight_array_to)
170 
171  !> Origin aerosol weight array.
172  type(aero_weight_array_t), intent(in) :: aero_weight_array_from
173  !> Destination aerosol weight array.
174  type(aero_weight_array_t), intent(inout) :: aero_weight_array_to
175 
176  if (any(shape(aero_weight_array_to%weight) &
177  /= shape(aero_weight_array_from%weight))) then
178  deallocate(aero_weight_array_to%weight)
179  call aero_weight_array_allocate_size(aero_weight_array_to, &
180  size(aero_weight_array_from%weight, 1), &
181  size(aero_weight_array_from%weight, 2))
182  end if
183  call aero_weight_copy(aero_weight_array_from%weight, &
184  aero_weight_array_to%weight)
185 
186  end subroutine aero_weight_array_copy
187 
188 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
189 
190  !> Return the number of weight groups.
191  integer function aero_weight_array_n_group(aero_weight_array)
192 
193  !> Aerosol weight array.
194  type(aero_weight_array_t), intent(in) :: aero_weight_array
195 
196  aero_weight_array_n_group = size(aero_weight_array%weight, 1)
197 
198  end function aero_weight_array_n_group
199 
200 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
201 
202  !> Return the number of weight classes.
203  integer function aero_weight_array_n_class(aero_weight_array)
204 
205  !> Aerosol weight array.
206  type(aero_weight_array_t), intent(in) :: aero_weight_array
207 
208  aero_weight_array_n_class = size(aero_weight_array%weight, 2)
209 
210  end function aero_weight_array_n_class
211 
212 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
213 
214  !> Scale the weights by the given factor, so
215  !> <tt>new_weight = old_weight * factor</tt>.
216  subroutine aero_weight_array_scale(aero_weight_array, factor)
217 
218  !> Aerosol weight array to scale.
219  type(aero_weight_array_t), intent(inout) :: aero_weight_array
220  !> Factor to scale by.
221  real(kind=dp), intent(in) :: factor
222 
223  call aero_weight_scale(aero_weight_array%weight, factor)
224 
225  end subroutine aero_weight_array_scale
226 
227 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
228 
229  !> Combine \c aero_weight_array_delta into \c aero_weight_array
230  !> with a harmonic mean.
231  subroutine aero_weight_array_combine(aero_weight_array, &
232  aero_weight_array_delta)
233 
234  !> Aerosol weight array to combine into.
235  type(aero_weight_array_t), intent(inout) :: aero_weight_array
236  !> Aerosol weight array to combine from.
237  type(aero_weight_array_t), intent(in) :: aero_weight_array_delta
238 
239  call aero_weight_combine(aero_weight_array%weight, &
240  aero_weight_array_delta%weight)
241 
242  end subroutine aero_weight_array_combine
243 
244 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
245 
246  !> Adjust source and destination weights to reflect moving \c
247  !> sample_prop proportion of particles from
248  !> \c aero_weight_array_from to \c aero_weight_array_to.
249  subroutine aero_weight_array_shift(aero_weight_array_from, &
250  aero_weight_array_to, sample_prop, overwrite_to)
251 
252  !> Aerosol weight array to shift from.
253  type(aero_weight_array_t), intent(inout) :: aero_weight_array_from
254  !> Aerosol weight array to shift to.
255  type(aero_weight_array_t), intent(inout) :: aero_weight_array_to
256  !> Proportion of particles being transfered.
257  real(kind=dp), intent(in) :: sample_prop
258  !> Whether to overwrite the destination weight (default: no).
259  logical, intent(in), optional :: overwrite_to
260 
261  call aero_weight_shift(aero_weight_array_from%weight, &
262  aero_weight_array_to%weight, sample_prop, overwrite_to)
263 
264  end subroutine aero_weight_array_shift
265 
266 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
267 
268  !> Compute the number concentration for a particle (m^{-3}).
269  real(kind=dp) function aero_weight_array_single_num_conc(aero_weight_array, &
270  aero_particle)
271 
272  !> Aerosol weight array.
273  type(aero_weight_array_t), intent(in) :: aero_weight_array
274  !> Aerosol particle to compute number concentration for.
275  type(aero_particle_t), intent(in) :: aero_particle
276 
278  aero_weight_array%weight(aero_particle%weight_group, &
279  aero_particle%weight_class), aero_particle)
280 
282 
283 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
284 
285  !> Compute the total number concentration at a given radius (m^3).
286  real(kind=dp) function aero_weight_array_num_conc_at_radius( &
287  aero_weight_array, i_class, radius)
288 
289  !> Aerosol weight array.
290  type(aero_weight_array_t), intent(in) :: aero_weight_array
291  !> Weight class number.
292  integer, intent(in) :: i_class
293  !> Radius to compute number concentration at (m).
294  real(kind=dp), intent(in) :: radius
295 
296  integer :: i_group
297  real(kind=dp) :: num_conc(size(aero_weight_array%weight, 1))
298 
299  do i_group = 1,size(aero_weight_array%weight, 1)
300  num_conc(i_group) = aero_weight_num_conc_at_radius( &
301  aero_weight_array%weight(i_group, i_class), radius)
302  end do
303  ! harmonic mean (same as summing the computational volumes)
304  aero_weight_array_num_conc_at_radius = 1d0 / sum(1d0 / num_conc)
305 
307 
308 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
309 
310  !> Compute the number concentration for a particle (m^{-3}).
311  real(kind=dp) function aero_weight_array_num_conc(aero_weight_array, &
312  aero_particle)
313 
314  !> Aerosol weight array.
315  type(aero_weight_array_t), intent(in) :: aero_weight_array
316  !> Aerosol particle to compute number concentration for.
317  type(aero_particle_t), intent(in) :: aero_particle
318 
320  aero_weight_array, aero_particle%weight_class, &
321  aero_particle_radius(aero_particle))
322 
323  end function aero_weight_array_num_conc
324 
325 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
326 
327  !> Check whether a given aero_weight array is flat in total.
328  logical function aero_weight_array_check_flat(aero_weight_array)
329 
330  !> Aerosol weight array.
331  type(aero_weight_array_t), intent(in) :: aero_weight_array
332 
333  integer :: i_group, i_class
334 
335  ! check that all weights have correctly set exponents
336  do i_group = 1,size(aero_weight_array%weight, 1)
337  do i_class = 1,size(aero_weight_array%weight, 2)
339  aero_weight_array%weight(i_group, i_class))
340  end do
341  end do
342 
343  if (all(abs(sum(aero_weight_array%weight%exponent, 1)) &
344  < 1d-20 * sum(abs(aero_weight_array%weight%exponent), 1))) then
346  else
348  end if
349 
350  end function aero_weight_array_check_flat
351 
352 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
353 
354  !> Determine whether all weight functions in an array are monotone
355  !> increasing, monotone decreasing, or neither.
356  subroutine aero_weight_array_check_monotonicity(aero_weight_array, &
357  monotone_increasing, monotone_decreasing)
358 
359  !> Aerosol weight array.
360  type(aero_weight_array_t), intent(in) :: aero_weight_array
361  !> Whether all weights are monotone increasing.
362  logical, intent(out) :: monotone_increasing
363  !> Whether all weights are monotone decreasing.
364  logical, intent(out) :: monotone_decreasing
365 
366  integer :: i_group, i_class
367  logical :: mono_increasing_array(size(aero_weight_array%weight, 1), &
368  size(aero_weight_array%weight, 2))
369  logical :: mono_decreasing_array(size(aero_weight_array%weight, 1), &
370  size(aero_weight_array%weight, 2))
371 
372  do i_group = 1,size(aero_weight_array%weight, 1)
373  do i_class = 1,size(aero_weight_array%weight, 2)
375  aero_weight_array%weight(i_group, i_class), &
376  mono_increasing_array(i_group, i_class), &
377  mono_decreasing_array(i_group, i_class))
378  end do
379  end do
380 
381  monotone_increasing = all(mono_increasing_array)
382  monotone_decreasing = all(mono_decreasing_array)
383 
385 
386 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
387 
388  !> Compute the maximum and minimum number concentrations between the
389  !> given radii.
390  subroutine aero_weight_array_minmax_num_conc(aero_weight_array, i_class, &
391  radius_1, radius_2, num_conc_min, num_conc_max)
392 
393  !> Aerosol weight array.
394  type(aero_weight_array_t), intent(in) :: aero_weight_array
395  !> Weight class.
396  integer, intent(in) :: i_class
397  !> First radius.
398  real(kind=dp), intent(in) :: radius_1
399  !> Second radius.
400  real(kind=dp), intent(in) :: radius_2
401  !> Minimum number concentration.
402  real(kind=dp), intent(out) :: num_conc_min
403  !> Maximum number concentration.
404  real(kind=dp), intent(out) :: num_conc_max
405 
406  real(kind=dp) :: num_conc_1, num_conc_2
407  logical :: monotone_increasing, monotone_decreasing
408 
409  call aero_weight_array_check_monotonicity(aero_weight_array, &
410  monotone_increasing, monotone_decreasing)
411  call assert(857727714, monotone_increasing .or. monotone_decreasing)
412 
413  num_conc_1 = aero_weight_array_num_conc_at_radius(aero_weight_array, &
414  i_class, radius_1)
415  num_conc_2 = aero_weight_array_num_conc_at_radius(aero_weight_array, &
416  i_class, radius_2)
417  num_conc_min = min(num_conc_1, num_conc_2)
418  num_conc_max = max(num_conc_1, num_conc_2)
419 
420  end subroutine aero_weight_array_minmax_num_conc
421 
422 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
423 
424  !> Choose a random group at the given radius, with probability
425  !> inversely proportional to group weight at that radius.
426  integer function aero_weight_array_rand_group(aero_weight_array, i_class, &
427  radius)
428 
429  !> Aerosol weight array.
430  type(aero_weight_array_t), intent(in) :: aero_weight_array
431  !> Weight class to select from.
432  integer, intent(in) :: i_class
433  !> Radius to sample group at (m).
434  real(kind=dp), intent(in) :: radius
435 
436  real(kind=dp) :: comp_vols(size(aero_weight_array%weight, 1))
437  integer :: i_group
438 
439  do i_group = 1,size(aero_weight_array%weight, 1)
440  comp_vols(i_group) = 1d0 / aero_weight_num_conc_at_radius( &
441  aero_weight_array%weight(i_group, i_class), radius)
442  end do
444 
445  end function aero_weight_array_rand_group
446 
447 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
448 
449  !> Read an aero_weight_array from a spec file.
450  subroutine spec_file_read_aero_weight_array(file, aero_weight_array)
451 
452  !> Spec file.
453  type(spec_file_t), intent(inout) :: file
454  !> Aerosol weight array.
455  type(aero_weight_array_t), intent(inout) :: aero_weight_array
456 
457  integer :: i_group, i_class
458 
459  do i_group = 1,size(aero_weight_array%weight, 1)
460  do i_class = 1,size(aero_weight_array%weight, 2)
461  call spec_file_read_aero_weight(file, &
462  aero_weight_array%weight(i_group, i_class))
463  end do
464  end do
465 
466  end subroutine spec_file_read_aero_weight_array
467 
468 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
469 
470  !> Determines the number of bytes required to pack the given value.
472 
473  !> Value to pack.
474  type(aero_weight_array_t), intent(in) :: val
475 
476  integer :: i_group, i_class
477 
479  pmc_mpi_pack_size_integer(size(val%weight, 1)) &
480  + pmc_mpi_pack_size_integer(size(val%weight, 2))
481 
482  do i_group = 1,size(val%weight, 1)
483  do i_class = 1,size(val%weight, 2)
486  + pmc_mpi_pack_size_aero_weight(val%weight(i_group, i_class))
487  end do
488  end do
489 
491 
492 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
493 
494  !> Packs the given value into the buffer, advancing position.
495  subroutine pmc_mpi_pack_aero_weight_array(buffer, position, val)
496 
497  !> Memory buffer.
498  character, intent(inout) :: buffer(:)
499  !> Current buffer position.
500  integer, intent(inout) :: position
501  !> Value to pack.
502  type(aero_weight_array_t), intent(in) :: val
503 
504 #ifdef PMC_USE_MPI
505  integer :: prev_position, i_group, i_class
506 
507  prev_position = position
508  call pmc_mpi_pack_integer(buffer, position, size(val%weight, 1))
509  call pmc_mpi_pack_integer(buffer, position, size(val%weight, 2))
510  do i_group = 1,size(val%weight, 1)
511  do i_class = 1,size(val%weight, 2)
512  call pmc_mpi_pack_aero_weight(buffer, position, &
513  val%weight(i_group, i_class))
514  end do
515  end do
516  call assert(84068036, &
517  position - prev_position <= pmc_mpi_pack_size_aero_weight_array(val))
518 #endif
519 
520  end subroutine pmc_mpi_pack_aero_weight_array
521 
522 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
523 
524  !> Unpacks the given value from the buffer, advancing position.
525  subroutine pmc_mpi_unpack_aero_weight_array(buffer, position, val)
526 
527  !> Memory buffer.
528  character, intent(inout) :: buffer(:)
529  !> Current buffer position.
530  integer, intent(inout) :: position
531  !> Value to pack.
532  type(aero_weight_array_t), intent(inout) :: val
533 
534 #ifdef PMC_USE_MPI
535  integer :: prev_position, n_group, n_class, i_group, i_class
536 
538  prev_position = position
539  call pmc_mpi_unpack_integer(buffer, position, n_group)
540  call pmc_mpi_unpack_integer(buffer, position, n_class)
541  call aero_weight_array_allocate_size(val, n_group, n_class)
542  do i_group = 1,size(val%weight, 1)
543  do i_class = 1,size(val%weight, 2)
544  call pmc_mpi_unpack_aero_weight(buffer, position, &
545  val%weight(i_group, i_class))
546  end do
547  end do
548  call assert(321022868, &
549  position - prev_position <= pmc_mpi_pack_size_aero_weight_array(val))
550 #endif
551 
552  end subroutine pmc_mpi_unpack_aero_weight_array
553 
554 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
555 
556  !> Write the \c aero_weight_group dimension to the given NetCDF file
557  !> if it is not already present and in any case return the
558  !> associated dimid.
559  subroutine aero_weight_netcdf_dim_aero_weight_group(aero_weight_array, &
560  ncid, dimid_aero_weight_group)
561 
562  !> Aero_weight structure array.
563  type(aero_weight_array_t), intent(in) :: aero_weight_array
564  !> NetCDF file ID, in data mode.
565  integer, intent(in) :: ncid
566  !> Dimid of the species dimension.
567  integer, intent(out) :: dimid_aero_weight_group
568 
569  integer :: status, i_group, n_group
570  integer :: varid_aero_weight_group
571  integer :: aero_weight_group_centers(size(aero_weight_array%weight, 1))
572 
573  ! try to get the dimension ID
574  status = nf90_inq_dimid(ncid, "aero_weight_group", dimid_aero_weight_group)
575  if (status == nf90_noerr) return
576  if (status /= nf90_ebaddim) call pmc_nc_check(status)
577 
578  ! dimension not defined, so define it
579  call pmc_nc_check(nf90_redef(ncid))
580 
581  n_group = size(aero_weight_array%weight, 1)
582 
583  call pmc_nc_check(nf90_def_dim(ncid, "aero_weight_group", n_group, &
584  dimid_aero_weight_group))
585  call pmc_nc_check(nf90_def_var(ncid, "aero_weight_group", nf90_int, &
586  dimid_aero_weight_group, varid_aero_weight_group))
587  call pmc_nc_check(nf90_put_att(ncid, varid_aero_weight_group, &
588  "description", "dummy dimension variable (no useful value)"))
589 
590  call pmc_nc_check(nf90_enddef(ncid))
591 
592  do i_group = 1,n_group
593  aero_weight_group_centers(i_group) = i_group
594  end do
595  call pmc_nc_check(nf90_put_var(ncid, varid_aero_weight_group, &
596  aero_weight_group_centers))
597 
599 
600 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
601 
602  !> Write the \c aero_weight_class dimension to the given NetCDF file
603  !> if it is not already present and in any case return the
604  !> associated dimid.
605  subroutine aero_weight_netcdf_dim_aero_weight_class(aero_weight_array, &
606  ncid, dimid_aero_weight_class)
607 
608  !> Aero_weight structure array.
609  type(aero_weight_array_t), intent(in) :: aero_weight_array
610  !> NetCDF file ID, in data mode.
611  integer, intent(in) :: ncid
612  !> Dimid of the species dimension.
613  integer, intent(out) :: dimid_aero_weight_class
614 
615  integer :: status, i_class, n_class
616  integer :: varid_aero_weight_class
617  integer :: aero_weight_class_centers(size(aero_weight_array%weight, 2))
618 
619  ! try to get the dimension ID
620  status = nf90_inq_dimid(ncid, "aero_weight_class", dimid_aero_weight_class)
621  if (status == nf90_noerr) return
622  if (status /= nf90_ebaddim) call pmc_nc_check(status)
623 
624  ! dimension not defined, so define it
625  call pmc_nc_check(nf90_redef(ncid))
626 
627  n_class = size(aero_weight_array%weight, 2)
628 
629  call pmc_nc_check(nf90_def_dim(ncid, "aero_weight_class", n_class, &
630  dimid_aero_weight_class))
631  call pmc_nc_check(nf90_def_var(ncid, "aero_weight_class", nf90_int, &
632  dimid_aero_weight_class, varid_aero_weight_class))
633  call pmc_nc_check(nf90_put_att(ncid, varid_aero_weight_class, &
634  "description", "dummy dimension variable (no useful value)"))
635 
636  call pmc_nc_check(nf90_enddef(ncid))
637 
638  do i_class = 1,n_class
639  aero_weight_class_centers(i_class) = i_class
640  end do
641  call pmc_nc_check(nf90_put_var(ncid, varid_aero_weight_class, &
642  aero_weight_class_centers))
643 
645 
646 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
647 
648  !> Write full aero_weight_array.
649  subroutine aero_weight_array_output_netcdf(aero_weight_array, ncid)
650 
651  !> Aero weight array to write.
652  type(aero_weight_array_t), intent(in) :: aero_weight_array
653  !> NetCDF file ID, in data mode.
654  integer, intent(in) :: ncid
655 
656  integer :: dimid_aero_weight_group, dimid_aero_weight_class
657 
658  !> \page output_format_aero_weight_array Output File Format: Aerosol Weighting Functions
659  !!
660  !! The aerosol weighting function NetCDF dimensions are:
661  !! - \b aero_weight_group: number of aerosol weighting groups
662  !! - \b aero_weight_class: number of aerosol weighting classes
663  !!
664  !! The aerosol weighting function NetCDF variables are:
665  !! - \b weight_type (no unit, dim \c aero_weight): the type of
666  !! each weighting function, with 0 = invalid weight, 1 = no
667  !! weight (\f$w(D) = 1\f$), 2 = power weight (\f$w(D) =
668  !! (D/D_0)^\alpha\f$), 3 = MFA weight (\f$w(D) =
669  !! (D/D_0)^{-3}\f$)
670  !! - \b weight_magnitude (unit m^{-3}, dim \c aero_weight): the
671  !! number concentration magnitude associated with each
672  !! weighting function
673  !! - \b weight_exponent (no unit, dim \c aero_weight): for each
674  !! weighting function, specifies the exponent \f$\alpha\f$ for
675  !! the power \c weight_type, the value -3 for the MFA \c
676  !! weight_type, and zero for any other \c weight_type
677 
678  call aero_weight_netcdf_dim_aero_weight_group(aero_weight_array, ncid, &
679  dimid_aero_weight_group)
680  call aero_weight_netcdf_dim_aero_weight_class(aero_weight_array, ncid, &
681  dimid_aero_weight_class)
682 
683  call pmc_nc_write_integer_2d(ncid, aero_weight_array%weight%type, &
684  "weight_type", &
685  (/ dimid_aero_weight_group, dimid_aero_weight_class /), &
686  description="type of each aerosol weighting function: 0 = invalid, " &
687  // "1 = none (w(D) = 1), 2 = power (w(D) = (D/D_0)^alpha), " &
688  // "3 = MFA (mass flow) (w(D) = (D/D_0)^(-3))")
689  call pmc_nc_write_real_2d(ncid, aero_weight_array%weight%magnitude, &
690  "weight_magnitude", &
691  (/ dimid_aero_weight_group, dimid_aero_weight_class /), &
692  unit="m^{-3}", &
693  description="magnitude for each weighting function")
694  call pmc_nc_write_real_2d(ncid, aero_weight_array%weight%exponent, &
695  "weight_exponent", &
696  (/ dimid_aero_weight_group, dimid_aero_weight_class /), unit="1", &
697  description="exponent alpha for the power weight_type, " &
698  // "set to -3 for MFA, and zero otherwise")
699 
700  end subroutine aero_weight_array_output_netcdf
701 
702 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
703 
704  !> Read full aero_weight_array.
705  subroutine aero_weight_array_input_netcdf(aero_weight_array, ncid)
706 
707  !> Aero weight array to read.
708  type(aero_weight_array_t), intent(inout) :: aero_weight_array
709  !> NetCDF file ID, in data mode.
710  integer, intent(in) :: ncid
711 
712  integer :: dimid_aero_weight_group, dimid_aero_weight_class, n_group
713  integer :: n_class
714  character(len=1000) :: name
715  integer, allocatable :: type(:, :)
716  real(kind=dp), allocatable :: magnitude(:, :), exponent(:, :)
717 
718  call pmc_nc_check(nf90_inq_dimid(ncid, "aero_weight_group", &
719  dimid_aero_weight_group))
720  call pmc_nc_check(nf90_inq_dimid(ncid, "aero_weight_class", &
721  dimid_aero_weight_class))
722  call pmc_nc_check(nf90_inquire_dimension(ncid, &
723  dimid_aero_weight_group, name, n_group))
724  call pmc_nc_check(nf90_inquire_dimension(ncid, &
725  dimid_aero_weight_class, name, n_class))
726  call assert(719221386, n_group < 1000)
727  call assert(520105999, n_class < 1000)
728 
729  allocate(type(n_group, n_class))
730  allocate(magnitude(n_group, n_class))
731  allocate(exponent(n_group, n_class))
732 
733  call pmc_nc_read_integer_2d(ncid, type, "weight_type")
734  call pmc_nc_read_real_2d(ncid, magnitude, "weight_magnitude")
735  call pmc_nc_read_real_2d(ncid, exponent, "weight_exponent")
736 
737  call assert(309191498, size(magnitude) == size(type))
738  call assert(588649520, size(magnitude) == size(exponent))
739 
740  call aero_weight_array_deallocate(aero_weight_array)
741  call aero_weight_array_allocate_size(aero_weight_array, n_group, n_class)
742 
743  aero_weight_array%weight%type = type
744  aero_weight_array%weight%magnitude = magnitude
745  aero_weight_array%weight%exponent = exponent
746 
747  deallocate(type)
748  deallocate(magnitude)
749  deallocate(exponent)
750 
751  end subroutine aero_weight_array_input_netcdf
752 
753 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
754 
755 end module pmc_aero_weight_array
real(kind=dp) function aero_weight_array_single_num_conc(aero_weight_array, aero_particle)
Compute the number concentration for a particle (m^{-3}).
real(kind=dp) function aero_weight_num_conc_at_radius(aero_weight, radius)
Compute the number concentration at a given radius (m^{-3}).
subroutine aero_weight_array_minmax_num_conc(aero_weight_array, i_class, radius_1, radius_2, num_conc_min, num_conc_max)
Compute the maximum and minimum number concentrations between the given radii.
subroutine pmc_mpi_unpack_aero_weight(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
An input file with extra data for printing messages.
Definition: spec_file.F90:59
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:347
elemental subroutine aero_weight_normalize(aero_weight)
Sets the aero_weight to a non-zero normalized value.
Definition: aero_weight.F90:94
subroutine aero_weight_array_check_monotonicity(aero_weight_array, monotone_increasing, monotone_decreasing)
Determine whether all weight functions in an array are monotone increasing, monotone decreasing...
subroutine aero_weight_netcdf_dim_aero_weight_class(aero_weight_array, ncid, dimid_aero_weight_class)
Write the aero_weight_class dimension to the given NetCDF file if it is not already present and in an...
subroutine aero_weight_array_allocate_size(aero_weight_array, n_group, n_class)
Allocates an aero_weight_array to the given size.
subroutine aero_weight_array_deallocate(aero_weight_array)
Free all storage.
elemental subroutine aero_weight_copy(aero_weight_from, aero_weight_to)
Copy an aero_weight.
The aero_weight_t structure and associated subroutines.
Definition: aero_weight.F90:9
Physical constants.
Definition: constants.F90:9
The aero_particle_t structure and associated subroutines.
elemental subroutine aero_weight_allocate(aero_weight)
Allocates an aero_weight.
Definition: aero_weight.F90:46
elemental subroutine aero_weight_shift(aero_weight_from, aero_weight_to, sample_prop, overwrite_to)
Adjust source and destination weights to reflect moving sample_prop proportion of particles from aero...
subroutine aero_weight_array_scale(aero_weight_array, factor)
Scale the weights by the given factor, so new_weight = old_weight * factor.
integer function pmc_mpi_pack_size_aero_weight_array(val)
Determines the number of bytes required to pack the given value.
subroutine pmc_nc_check(status)
Check the status of a NetCDF function call.
Definition: netcdf.F90:22
subroutine aero_weight_array_allocate_flat(aero_weight_array, n_class)
Allocates an aero_weight_array as flat weightings.
subroutine spec_file_read_aero_weight_array(file, aero_weight_array)
Read an aero_weight_array from a spec file.
integer function sample_cts_pdf(pdf)
Sample the given continuous probability density function.
Definition: rand.F90:489
subroutine aero_weight_check_valid_exponent(aero_weight)
Ensures that a weight function exponent is valid.
elemental subroutine aero_weight_scale(aero_weight, factor)
Scale the weight by the given fraction, so new_weight = old_weight * factor.
elemental real(kind=dp) function aero_particle_radius(aero_particle)
Total radius of the particle (m).
integer function aero_weight_array_rand_group(aero_weight_array, i_class, radius)
Choose a random group at the given radius, with probability inversely proportional to group weight at...
An aerosol size distribution weighting function.
Definition: aero_weight.F90:32
Random number generators.
Definition: rand.F90:9
real(kind=dp) function aero_weight_array_num_conc(aero_weight_array, aero_particle)
Compute the number concentration for a particle (m^{-3}).
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:534
integer function pmc_mpi_pack_size_aero_weight(val)
Determines the number of bytes required to pack the given value.
elemental subroutine aero_weight_combine(aero_weight, aero_weight_delta)
Combine aero_weight_delta into aero_weight with a harmonic mean.
subroutine pmc_mpi_pack_aero_weight(buffer, position, val)
Packs the given value into the buffer, advancing position.
Common utility subroutines.
Definition: util.F90:9
subroutine aero_weight_array_input_netcdf(aero_weight_array, ncid)
Read full aero_weight_array.
subroutine aero_weight_array_allocate(aero_weight_array)
Allocates an aero_weight_array.
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:771
subroutine aero_weight_netcdf_dim_aero_weight_group(aero_weight_array, ncid, dimid_aero_weight_group)
Write the aero_weight_group dimension to the given NetCDF file if it is not already present and in an...
An array of aerosol size distribution weighting functions.
subroutine pmc_nc_read_real_2d(ncid, var, name, must_be_present)
Read a simple real 2D array from a NetCDF file.
Definition: netcdf.F90:288
Wrapper functions for MPI.
Definition: mpi.F90:13
subroutine aero_weight_array_allocate_power(aero_weight_array, n_class, exponent)
Allocates an aero_weight_array as power weightings.
real(kind=dp) function aero_weight_array_num_conc_at_radius(aero_weight_array, i_class, radius)
Compute the total number concentration at a given radius (m^3).
The aero_weight_array_t structure and associated subroutines.
subroutine pmc_mpi_unpack_aero_weight_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine pmc_mpi_pack_aero_weight_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
integer function aero_weight_array_n_class(aero_weight_array)
Return the number of weight classes.
Single aerosol particle data structure.
subroutine pmc_nc_write_real_2d(ncid, var, name, dimids, dim_name_1, dim_name_2, unit, long_name, standard_name, description)
Write a simple real 2D array to a NetCDF file.
Definition: netcdf.F90:602
Reading formatted text input.
Definition: spec_file.F90:43
integer function aero_weight_array_n_group(aero_weight_array)
Return the number of weight groups.
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_write_integer_2d(ncid, var, name, dimids, dim_name_1, dim_name_2, unit, long_name, standard_name, description)
Write a simple integer 2D array to a NetCDF file.
Definition: netcdf.F90:656
subroutine aero_weight_array_allocate_nummass(aero_weight_array, n_class)
Allocates an aero_weight_array as joint flat/power-3 weightings..
subroutine aero_weight_check_monotonicity(aero_weight, monotone_increasing, monotone_decreasing)
Determine whether a weight function is monotone increasing, monotone decreasing, or neither...
logical function aero_weight_array_check_flat(aero_weight_array)
Check whether a given aero_weight array is flat in total.
real(kind=dp) function aero_weight_num_conc(aero_weight, aero_particle)
Compute the number concentration for a particle (m^{-3}).
subroutine aero_weight_array_zero(aero_weight_array)
Zeros the contents of the aero_weight_array.
subroutine pmc_nc_read_integer_2d(ncid, var, name, must_be_present)
Read a simple integer 2D array from a NetCDF file.
Definition: netcdf.F90:323
elemental subroutine aero_weight_zero(aero_weight)
Zeros the contents of the aero_weight.
Definition: aero_weight.F90:80
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition: util.F90:102
subroutine aero_weight_array_shift(aero_weight_array_from, aero_weight_array_to, sample_prop, overwrite_to)
Adjust source and destination weights to reflect moving sample_prop proportion of particles from aero...
subroutine aero_weight_array_combine(aero_weight_array, aero_weight_array_delta)
Combine aero_weight_array_delta into aero_weight_array with a harmonic mean.
subroutine aero_weight_array_copy(aero_weight_array_from, aero_weight_array_to)
Copy an aero_weight_array.
elemental subroutine aero_weight_deallocate(aero_weight)
Free all storage.
Definition: aero_weight.F90:70
elemental subroutine aero_weight_array_normalize(aero_weight_array)
Normalizes the aero_weight_array to a non-zero value.