PartMC  2.3.0
aero_dist.F90
Go to the documentation of this file.
1 ! Copyright (C) 2005-2012 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_dist module.
7 
8 !> The aero_dist_t structure and associated subroutines.
9 !!
10 !! The initial size distributions are computed as number densities, so
11 !! they can be used for both sectional and particle-resolved
12 !! simulations. The routine dist_to_n() converts a number concentration
13 !! distribution to an actual number of particles ready for a
14 !! particle-resolved simulation.
15 !!
16 !! Initial distributions should be normalized so that <tt>sum(n_den) =
17 !! 1/log_width</tt>.
19 
20  use pmc_bin_grid
21  use pmc_util
22  use pmc_constants
23  use pmc_spec_file
24  use pmc_aero_data
25  use pmc_aero_mode
26  use pmc_mpi
27  use pmc_rand
28 #ifdef PMC_USE_MPI
29  use mpi
30 #endif
31 
32  !> A complete aerosol distribution, consisting of several modes.
34  !> Number of modes.
35  integer :: n_mode
36  !> Internally mixed modes [length \c n_mode].
37  type(aero_mode_t), pointer :: mode(:)
38  end type aero_dist_t
39 
40 contains
41 
42 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
43 
44  !> Allocates an aero_dist.
45  subroutine aero_dist_allocate(aero_dist)
46 
47  !> Aerosol distribution.
48  type(aero_dist_t), intent(out) :: aero_dist
49 
50  integer :: i
51 
52  aero_dist%n_mode = 0
53  allocate(aero_dist%mode(0))
54 
55  end subroutine aero_dist_allocate
56 
57 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
58 
59  !> Allocates an aero_dist of the given size.
60  subroutine aero_dist_allocate_size(aero_dist, n_mode, n_spec)
61 
62  !> Aerosol distribution.
63  type(aero_dist_t), intent(out) :: aero_dist
64  !> Number of modes.
65  integer, intent(in) :: n_mode
66  !> Number of species.
67  integer, intent(in) :: n_spec
68 
69  integer :: i
70 
71  aero_dist%n_mode = n_mode
72  allocate(aero_dist%mode(n_mode))
73  do i = 1,n_mode
74  call aero_mode_allocate_size(aero_dist%mode(i), n_spec)
75  end do
76 
77  end subroutine aero_dist_allocate_size
78 
79 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
80 
81  !> Free all storage.
82  subroutine aero_dist_deallocate(aero_dist)
83 
84  !> Aerosol distribution.
85  type(aero_dist_t), intent(inout) :: aero_dist
86 
87  integer :: i
88 
89  do i = 1,aero_dist%n_mode
90  call aero_mode_deallocate(aero_dist%mode(i))
91  end do
92  deallocate(aero_dist%mode)
93 
94  end subroutine aero_dist_deallocate
95 
96 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
97 
98  !> Copy an aero_dist.
99  subroutine aero_dist_copy(aero_dist_from, aero_dist_to)
100 
101  !> Aero_dist original.
102  type(aero_dist_t), intent(in) :: aero_dist_from
103  !> Aero_dist copy.
104  type(aero_dist_t), intent(inout) :: aero_dist_to
105 
106  integer :: n_spec, i
107 
108  if (aero_dist_from%n_mode > 0) then
109  n_spec = size(aero_dist_from%mode(1)%vol_frac)
110  else
111  n_spec = 0
112  end if
113  call aero_dist_deallocate(aero_dist_to)
114  call aero_dist_allocate_size(aero_dist_to, aero_dist_from%n_mode, n_spec)
115  do i = 1,aero_dist_from%n_mode
116  call aero_mode_copy(aero_dist_from%mode(i), aero_dist_to%mode(i))
117  end do
118 
119  end subroutine aero_dist_copy
120 
121 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
122 
123  !> Returns the total number concentration of a distribution. (#/m^3)
124  real(kind=dp) function aero_dist_total_num_conc(aero_dist)
125 
126  !> Aerosol distribution.
127  type(aero_dist_t), intent(in) :: aero_dist
128 
129  integer :: i_mode
130 
132  do i_mode = 1,aero_dist%n_mode
134  + aero_mode_total_num_conc(aero_dist%mode(i_mode))
135  end do
136 
137  end function aero_dist_total_num_conc
138 
139 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
140 
141  !> Returns the total number of particles of a distribution.
142  real(kind=dp) function aero_dist_number(aero_dist, aero_weight)
143 
144  !> Aerosol distribution.
145  type(aero_dist_t), intent(in) :: aero_dist
146  !> Aerosol weight.
147  type(aero_weight_t), intent(in) :: aero_weight
148 
149  integer :: i_mode
150 
151  aero_dist_number = 0d0
152  do i_mode = 1,aero_dist%n_mode
154  + aero_mode_number(aero_dist%mode(i_mode), aero_weight)
155  end do
156 
157  end function aero_dist_number
158 
159 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
160 
161  !> Return the binned number concentration for an aero_dist.
162  subroutine aero_dist_num_conc(aero_dist, bin_grid, aero_data, &
163  num_conc)
164 
165  !> Aero dist for which to compute number concentration.
166  type(aero_dist_t), intent(in) :: aero_dist
167  !> Bin grid.
168  type(bin_grid_t), intent(in) :: bin_grid
169  !> Aerosol data.
170  type(aero_data_t), intent(in) :: aero_data
171  !> Number concentration (#(ln(r))d(ln(r))).
172  real(kind=dp), intent(out) :: num_conc(bin_grid%n_bin)
173 
174  integer :: i_mode
175  real(kind=dp) :: mode_num_conc(size(num_conc, 1))
176 
177  num_conc = 0d0
178  do i_mode = 1,aero_dist%n_mode
179  call aero_mode_num_conc(aero_dist%mode(i_mode), bin_grid, &
180  aero_data, mode_num_conc)
181  num_conc = num_conc + mode_num_conc
182  end do
183 
184  end subroutine aero_dist_num_conc
185 
186 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
187 
188  !> Return the binned per-species volume concentration for an
189  !> aero_dist.
190  subroutine aero_dist_vol_conc(aero_dist, bin_grid, aero_data, &
191  vol_conc)
192 
193  !> Aero dist for which to compute volume concentration.
194  type(aero_dist_t), intent(in) :: aero_dist
195  !> Bin grid.
196  type(bin_grid_t), intent(in) :: bin_grid
197  !> Aerosol data.
198  type(aero_data_t), intent(in) :: aero_data
199  !> Volume concentration (V(ln(r))d(ln(r))).
200  real(kind=dp), intent(out) :: vol_conc(bin_grid%n_bin, aero_data%n_spec)
201 
202  integer :: i_mode
203  real(kind=dp) :: mode_vol_conc(size(vol_conc, 1), size(vol_conc, 2))
204 
205  vol_conc = 0d0
206  do i_mode = 1,aero_dist%n_mode
207  call aero_mode_vol_conc(aero_dist%mode(i_mode), bin_grid, &
208  aero_data, mode_vol_conc)
209  vol_conc = vol_conc + mode_vol_conc
210  end do
211 
212  end subroutine aero_dist_vol_conc
213 
214 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
215 
216  !> Whether any of the modes are of the given type.
217  elemental logical function aero_dist_contains_aero_mode_type(aero_dist, &
218  aero_mode_type)
219 
220  !> Aerosol distribution.
221  type(aero_dist_t), intent(in) :: aero_dist
222  !> Aerosol mode type to test for.
223  integer, intent(in) :: aero_mode_type
224 
226  = any(aero_dist%mode%type == aero_mode_type)
227 
229 
230 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
231 
232  !> Determine the current aero_dist and rate by interpolating at the
233  !> current time with the lists of aero_dists and rates.
234  subroutine aero_dist_interp_1d(aero_dist_list, time_list, &
235  rate_list, time, aero_dist, rate)
236 
237  !> Gas states.
238  type(aero_dist_t), intent(in) :: aero_dist_list(:)
239  !> Times (s).
240  real(kind=dp), intent(in) :: time_list(size(aero_dist_list))
241  !> Rates (s^{-1}).
242  real(kind=dp), intent(in) :: rate_list(size(aero_dist_list))
243  !> Current time (s).
244  real(kind=dp), intent(in) :: time
245  !> Current gas state.
246  type(aero_dist_t), intent(inout) :: aero_dist
247  !> Current rate (s^{-1}).
248  real(kind=dp), intent(out) :: rate
249 
250  integer :: n, p, n_bin, n_spec, i, i_new
251  real(kind=dp) :: y, alpha
252 
253  n = size(aero_dist_list)
254  p = find_1d(n, time_list, time)
255  if (p == 0) then
256  ! before the start, just use the first state and rate
257  call aero_dist_copy(aero_dist_list(1), aero_dist)
258  rate = rate_list(1)
259  elseif (p == n) then
260  ! after the end, just use the last state and rate
261  call aero_dist_copy(aero_dist_list(n), aero_dist)
262  rate = rate_list(n)
263  else
264  ! in the middle, use the previous dist
265  call aero_dist_copy(aero_dist_list(p), aero_dist)
266  rate = rate_list(p)
267  end if
268 
269  end subroutine aero_dist_interp_1d
270 
271 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
272 
273  !> Read continuous aerosol distribution composed of several modes.
274  subroutine spec_file_read_aero_dist(file, aero_data, aero_dist)
275 
276  !> Spec file to read data from.
277  type(spec_file_t), intent(inout) :: file
278  !> Aero_data data.
279  type(aero_data_t), intent(inout) :: aero_data
280  !> Aerosol dist.
281  type(aero_dist_t), intent(inout) :: aero_dist
282 
283  type(aero_mode_t), pointer :: new_aero_mode_list(:)
284  type(aero_mode_t) :: aero_mode
285  integer :: i, j
286  logical :: eof
287 
288  ! note that the <p> is needed below to force correct paragraph
289  ! breaking by doxygen
290 
291  !> \page input_format_aero_dist Input File Format: Aerosol Distribution
292  !!
293  !! <p>An aerosol distribution file consists of zero or more modes,
294  !! each in the format described by \subpage input_format_aero_mode
295  !!
296  !! See also:
297  !! - \ref spec_file_format --- the input file text format
298  !! - \ref input_format_aero_mode --- the format for each mode
299  !! of an aerosol distribution
300 
301  call aero_dist_deallocate(aero_dist)
302  call aero_dist_allocate(aero_dist)
303  call aero_mode_allocate(aero_mode)
304  call spec_file_read_aero_mode(file, aero_data, aero_mode, eof)
305  do while (.not. eof)
306  aero_dist%n_mode = aero_dist%n_mode + 1
307  allocate(new_aero_mode_list(aero_dist%n_mode))
308  do i = 1,aero_dist%n_mode
309  call aero_mode_allocate(new_aero_mode_list(i))
310  end do
311  call aero_mode_copy(aero_mode, &
312  new_aero_mode_list(aero_dist%n_mode))
313  do i = 1,(aero_dist%n_mode - 1)
314  call aero_mode_copy(aero_dist%mode(i), &
315  new_aero_mode_list(i))
316  call aero_mode_deallocate(aero_dist%mode(i))
317  end do
318  deallocate(aero_dist%mode)
319  aero_dist%mode => new_aero_mode_list
320  nullify(new_aero_mode_list)
321  call spec_file_read_aero_mode(file, aero_data, aero_mode, eof)
322  end do
323  call aero_mode_deallocate(aero_mode)
324 
325  end subroutine spec_file_read_aero_dist
326 
327 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
328 
329  !> Read an array of aero_dists with associated times and rates from
330  !> the given file.
331  subroutine spec_file_read_aero_dists_times_rates(file, aero_data, &
332  times, rates, aero_dists)
333 
334  !> Spec file to read data from.
335  type(spec_file_t), intent(inout) :: file
336  !> Aero data.
337  type(aero_data_t), intent(inout) :: aero_data
338  !> Times (s).
339  real(kind=dp), pointer :: times(:)
340  !> Rates (s^{-1}).
341  real(kind=dp), pointer :: rates(:)
342  !> Aero dists.
343  type(aero_dist_t), pointer :: aero_dists(:)
344 
345  type(spec_line_t) :: aero_dist_line
346  type(spec_file_t) :: aero_dist_file
347  integer :: n_time, i_time
348  character(len=SPEC_LINE_MAX_VAR_LEN), pointer :: names(:)
349  real(kind=dp), pointer :: data(:,:)
350 
351  !> \page input_format_aero_dist_profile Input File Format: Aerosol Distribution Profile
352  !!
353  !! An aerosol distribution profile input file must consist of
354  !! three lines:
355  !! - the first line must begin with \c time and should be followed
356  !! by \f$N\f$ space-separated real scalars, giving the times (in
357  !! s after the start of the simulation) of the aerosol
358  !! distrbution set points --- the times must be in increasing
359  !! order
360  !! - the second line must begin with \c rate and should be
361  !! followed by \f$N\f$ space-separated real scalars, giving the
362  !! values at the corresponding times
363  !! - the third line must begin with \c dist and should be followed
364  !! by \f$N\f$ space-separated filenames, each specifying an
365  !! aerosol distribution in the format \ref input_format_aero_dist
366  !! at the corresponding time
367  !!
368  !! The units of the \c rate line depend on the type of aerosol
369  !! distribution profile:
370  !! - Emissions aerosol profiles have rates with units m/s ---
371  !! the aerosol distribution number concentrations are multiplied
372  !! by the rate to give an emission rate with unit #/(m^2 s)
373  !! which is then divided by the current mixing layer height
374  !! to give a per-volume emission rate.
375  !! - Background aerosol profiles have rates with units \f$\rm
376  !! s^{-1}\f$, which is the dilution rate between the background
377  !! and the simulated air parcel. That is, if the simulated
378  !! number concentration is \f$N\f$ and the background number
379  !! concentration is \f$N_{\rm back}\f$, then dilution is modeled
380  !! as \f$\dot{N} = r N_{\rm back} - r N\f$, where \f$r\f$ is the
381  !! rate.
382  !!
383  !! Between the specified times the aerosol profile is interpolated
384  !! step-wise and kept constant at its last value. That is, if the
385  !! times are \f$t_i\f$, the rates are \f$r_i\f$, and the aerosol
386  !! distributions are \f$a_i\f$ (all with \f$i = 1,\ldots,n\f$),
387  !! then between times \f$t_i\f$ and \f$t_{i+1}\f$ the aerosol
388  !! state is constant at \f$r_i a_i\f$. Before time \f$t_1\f$ the
389  !! aerosol state is \f$r_1 a_1\f$, while after time \f$t_n\f$ it
390  !! is \f$r_n a_n\f$.
391  !!
392  !! Example: an emissions aerosol profile could be:
393  !! <pre>
394  !! time 0 600 1800 # time (in s) after sim start
395  !! rate 1 0.5 1 # scaling factor in m/s
396  !! dist dist1.dat dist2.dat dist3.dat # aerosol distribution files
397  !! </pre>
398  !! Here the emissions between 0&nbsp;min and 10&nbsp;min are given
399  !! by <tt>dist1.dat</tt> (with the number concentration
400  !! interpreted as having units 1/(m^2 s)), the emissions between
401  !! 10&nbsp;min and 30&nbsp;min are given by <tt>dist2.dat</tt>
402  !! (scaled by 0.5), while the emissions after 30&nbsp;min are
403  !! given by <tt>dist3.dat</tt>.
404  !!
405  !! See also:
406  !! - \ref spec_file_format --- the input file text format
407  !! - \ref input_format_aero_data --- the aerosol species list
408  !! and material data
409  !! - \ref input_format_aero_dist --- the format of the
410  !! instantaneous aerosol distribution files
411 
412  ! read the data from the file
413  allocate(names(0))
414  allocate(data(0,0))
415  call spec_file_read_real_named_array(file, 2, names, data)
416  call spec_line_allocate(aero_dist_line)
417  call spec_file_read_line_no_eof(file, aero_dist_line)
418  call spec_file_check_line_name(file, aero_dist_line, "dist")
419  call spec_file_check_line_length(file, aero_dist_line, size(data, 2))
420 
421  ! check the data size
422  if (trim(names(1)) /= 'time') then
423  call die_msg(570205795, 'row 1 in ' // trim(file%name) &
424  // ' must start with: time not: ' // trim(names(1)))
425  end if
426  if (trim(names(2)) /= 'rate') then
427  call die_msg(221270915, 'row 2 in ' // trim(file%name) &
428  // ' must start with: rate not: ' // trim(names(1)))
429  end if
430  n_time = size(data, 2)
431  if (n_time < 1) then
432  call die_msg(457229710, 'each line in ' // trim(file%name) &
433  // ' must contain at least one data value')
434  end if
435 
436  ! copy over the data
437  do i_time = 1,size(aero_dists)
438  call aero_dist_deallocate(aero_dists(i_time))
439  end do
440  deallocate(aero_dists)
441  deallocate(times)
442  deallocate(rates)
443  allocate(aero_dists(n_time))
444  allocate(times(n_time))
445  allocate(rates(n_time))
446  do i_time = 1,n_time
447  call spec_file_open(aero_dist_line%data(i_time), aero_dist_file)
448  call aero_dist_allocate(aero_dists(i_time))
449  call spec_file_read_aero_dist(aero_dist_file, &
450  aero_data, aero_dists(i_time))
451  call spec_file_close(aero_dist_file)
452  times(i_time) = data(1,i_time)
453  rates(i_time) = data(2,i_time)
454  end do
455  deallocate(names)
456  deallocate(data)
457  call spec_line_deallocate(aero_dist_line)
458 
459  end subroutine spec_file_read_aero_dists_times_rates
460 
461 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
462 
463  !> Determines the number of bytes required to pack the given value.
464  integer function pmc_mpi_pack_size_aero_dist(val)
465 
466  !> Value to pack.
467  type(aero_dist_t), intent(in) :: val
468 
469  integer :: i, total_size
470 
471  total_size = pmc_mpi_pack_size_integer(val%n_mode)
472  do i = 1,size(val%mode)
473  total_size = total_size + pmc_mpi_pack_size_aero_mode(val%mode(i))
474  end do
475  pmc_mpi_pack_size_aero_dist = total_size
476 
477  end function pmc_mpi_pack_size_aero_dist
478 
479 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
480 
481  !> Packs the given value into the buffer, advancing position.
482  subroutine pmc_mpi_pack_aero_dist(buffer, position, val)
483 
484  !> Memory buffer.
485  character, intent(inout) :: buffer(:)
486  !> Current buffer position.
487  integer, intent(inout) :: position
488  !> Value to pack.
489  type(aero_dist_t), intent(in) :: val
490 
491 #ifdef PMC_USE_MPI
492  integer :: prev_position, i
493 
494  prev_position = position
495  call pmc_mpi_pack_integer(buffer, position, val%n_mode)
496  do i = 1,size(val%mode)
497  call pmc_mpi_pack_aero_mode(buffer, position, val%mode(i))
498  end do
499  call assert(440557910, &
500  position - prev_position <= pmc_mpi_pack_size_aero_dist(val))
501 #endif
502 
503  end subroutine pmc_mpi_pack_aero_dist
504 
505 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
506 
507  !> Unpacks the given value from the buffer, advancing position.
508  subroutine pmc_mpi_unpack_aero_dist(buffer, position, val)
509 
510  !> Memory buffer.
511  character, intent(inout) :: buffer(:)
512  !> Current buffer position.
513  integer, intent(inout) :: position
514  !> Value to pack.
515  type(aero_dist_t), intent(inout) :: val
516 
517 #ifdef PMC_USE_MPI
518  integer :: prev_position, i
519 
520  call aero_dist_deallocate(val)
521  prev_position = position
522  call pmc_mpi_unpack_integer(buffer, position, val%n_mode)
523  allocate(val%mode(val%n_mode))
524  do i = 1,size(val%mode)
525  call aero_mode_allocate(val%mode(i))
526  call pmc_mpi_unpack_aero_mode(buffer, position, val%mode(i))
527  end do
528  call assert(742535268, &
529  position - prev_position <= pmc_mpi_pack_size_aero_dist(val))
530 #endif
531 
532  end subroutine pmc_mpi_unpack_aero_dist
533 
534 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
535 
536 end module pmc_aero_dist
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_dist_interp_1d(aero_dist_list, time_list, rate_list, time, aero_dist, rate)
Determine the current aero_dist and rate by interpolating at the current time with the lists of aero_...
Definition: aero_dist.F90:234
subroutine pmc_mpi_unpack_aero_mode(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: aero_mode.F90:1100
integer function pmc_mpi_pack_size_aero_dist(val)
Determines the number of bytes required to pack the given value.
Definition: aero_dist.F90:464
The aero_dist_t structure and associated subroutines.
Definition: aero_dist.F90:18
subroutine spec_file_open(filename, file)
Open a spec file for reading.
Definition: spec_file.F90:111
Physical constants.
Definition: constants.F90:9
integer function find_1d(n, x_vals, x)
Find the position of a real number in an arbitrary 1D array.
Definition: util.F90:577
subroutine aero_dist_allocate_size(aero_dist, n_mode, n_spec)
Allocates an aero_dist of the given size.
Definition: aero_dist.F90:60
subroutine aero_mode_vol_conc(aero_mode, bin_grid, aero_data, vol_conc)
Return the binned per-species volume concentration for an aero_mode.
Definition: aero_mode.F90:466
subroutine aero_dist_deallocate(aero_dist)
Free all storage.
Definition: aero_dist.F90:82
Random number generators.
Definition: rand.F90:9
integer function pmc_mpi_pack_size_aero_mode(val)
Determines the number of bytes required to pack the given value.
Definition: aero_mode.F90:1046
subroutine spec_file_check_line_length(file, line, length)
Check that the length of the line data is as given.
Definition: spec_file.F90:446
real(kind=dp) function aero_mode_number(aero_mode, aero_weight)
Return the total number of computational particles for an aero_mode.
Definition: aero_mode.F90:547
subroutine aero_mode_num_conc(aero_mode, bin_grid, aero_data, num_conc)
Return the binned number concentration for an aero_mode.
Definition: aero_mode.F90:431
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:534
subroutine aero_mode_allocate(aero_mode)
Allocates an aero_mode.
Definition: aero_mode.F90:99
real(kind=dp) function aero_mode_total_num_conc(aero_mode)
Returns the total number concentration of a mode. (#/m^3)
Definition: aero_mode.F90:176
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
A complete aerosol distribution, consisting of several modes.
Definition: aero_dist.F90:33
Common utility subroutines.
Definition: util.F90:9
An aerosol size distribution mode.
Definition: aero_mode.F90:47
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:771
real(kind=dp) function aero_dist_total_num_conc(aero_dist)
Returns the total number concentration of a distribution. (#/m^3)
Definition: aero_dist.F90:124
subroutine aero_mode_copy(aero_mode_from, aero_mode_to)
Copy an aero_mode.
Definition: aero_mode.F90:146
Wrapper functions for MPI.
Definition: mpi.F90:13
real(kind=dp) function aero_dist_number(aero_dist, aero_weight)
Returns the total number of particles of a distribution.
Definition: aero_dist.F90:142
The aero_mode_t structure and associated subroutines.
Definition: aero_mode.F90:9
1D grid, either logarithmic or linear.
Definition: bin_grid.F90:33
The bin_grid_t structure and associated subroutines.
Definition: bin_grid.F90:9
Reading formatted text input.
Definition: spec_file.F90:43
subroutine spec_file_close(file)
Close a spec file.
Definition: spec_file.F90:134
subroutine aero_dist_num_conc(aero_dist, bin_grid, aero_data, num_conc)
Return the binned number concentration for an aero_dist.
Definition: aero_dist.F90:162
subroutine aero_dist_copy(aero_dist_from, aero_dist_to)
Copy an aero_dist.
Definition: aero_dist.F90:99
elemental logical function aero_dist_contains_aero_mode_type(aero_dist, aero_mode_type)
Whether any of the modes are of the given type.
Definition: aero_dist.F90:217
subroutine spec_line_allocate(spec_line)
Allocates memory for a spec_line.
Definition: spec_line.F90:39
subroutine pmc_mpi_unpack_aero_dist(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: aero_dist.F90:508
subroutine spec_line_deallocate(spec_line)
Frees all storage.
Definition: spec_line.F90:65
subroutine pmc_mpi_pack_aero_dist(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: aero_dist.F90:482
subroutine pmc_mpi_pack_aero_mode(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: aero_mode.F90:1068
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 aero_mode_allocate_size(aero_mode, n_spec)
Allocates an aero_mode of the given size.
Definition: aero_mode.F90:114
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_mode_deallocate(aero_mode)
Free all storage.
Definition: aero_mode.F90:131
subroutine aero_dist_allocate(aero_dist)
Allocates an aero_dist.
Definition: aero_dist.F90:45
subroutine aero_dist_vol_conc(aero_dist, bin_grid, aero_data, vol_conc)
Return the binned per-species volume concentration for an aero_dist.
Definition: aero_dist.F90:190