PartMC  2.3.0
aero_binned.F90
Go to the documentation of this file.
1 ! Copyright (C) 2005-2015 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 the.
4 
5 !> \file
6 !> The pmc_aero_binned module.
7 
8 !> The aero_binned_t structure and associated subroutines.
10 
11  use pmc_bin_grid
13  use pmc_spec_file
14  use pmc_util
15  use pmc_bin_grid
16  use pmc_aero_dist
17  use pmc_mpi
18  use pmc_aero_data
19 #ifdef PMC_USE_MPI
20  use mpi
21 #endif
22 
23  !> Aerosol number and volume distributions stored per bin.
24  !!
25  !! These quantities are densities both in volume (per m^3) and in
26  !! radius (per log_width). The total concentration per volume is computed as
27  !! sum(aero_binned%%num_conc * bin_grid%%log_width).
28  !!
29  !! An aero_binned_t is similar to an aero_dist_t in that they both
30  !! store binned aerosol distributions. The difference is that an
31  !! aero_dist_t has the same composition in every bin, whereas an
32  !! aero_binned_t can have aerosol composition that varies per bin.
34  !> Number concentration per bin (#/m^3/log_width).
35  !! Array length is typically \c bin_grid%%n_bin.
36  real(kind=dp), pointer :: num_conc(:)
37  !> Volume concentration per bin and per species (m^3/m^3/log_width).
38  !! Array size is typically \c bin_grid%%n_bin x \c aero_data%%n_spec.
39  real(kind=dp), pointer :: vol_conc(:,:)
40  end type aero_binned_t
41 
42 contains
43 
44 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
45 
46  !> Allocate an aero_binned_t.
47  subroutine aero_binned_allocate(aero_binned)
48 
49  !> Structure to be allocated.
50  type(aero_binned_t), intent(out) :: aero_binned
51 
52  allocate(aero_binned%num_conc(0))
53  allocate(aero_binned%vol_conc(0, 0))
54 
55  end subroutine aero_binned_allocate
56 
57 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
58 
59  !> Allocate an aero_binned_t of the given size.
60  subroutine aero_binned_allocate_size(aero_binned, n_bin, n_spec)
61 
62  !> Structure to be allocated.
63  type(aero_binned_t), intent(out) :: aero_binned
64  !> Number of aerosol bins to allocate (typically \c bin_grid%%n_bin).
65  integer, intent(in) :: n_bin
66  !> Number of aerosol species to allocate (typically
67  !> \c aero_data%%n_spec).
68  integer, intent(in) :: n_spec
69 
70  allocate(aero_binned%num_conc(n_bin))
71  allocate(aero_binned%vol_conc(n_bin, n_spec))
72  call aero_binned_zero(aero_binned)
73 
74  end subroutine aero_binned_allocate_size
75 
76 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
77 
78  !> Free internal memory in an aero_binned_t structure.
79  subroutine aero_binned_deallocate(aero_binned)
80 
81  !> Structure to free.
82  type(aero_binned_t), intent(inout) :: aero_binned
83 
84  deallocate(aero_binned%num_conc)
85  deallocate(aero_binned%vol_conc)
86 
87  end subroutine aero_binned_deallocate
88 
89 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
90 
91  !> Set all internal data in an aero_binned_t structure to zero.
92  subroutine aero_binned_zero(aero_binned)
93 
94  !> Structure to zero.
95  type(aero_binned_t), intent(inout) :: aero_binned
96 
97  aero_binned%num_conc = 0d0
98  aero_binned%vol_conc = 0d0
99 
100  end subroutine aero_binned_zero
101 
102 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
103 
104  !> Add two aero_binned_t structures together.
105  !!
106  !! Symbolically does aero_binned = aero_binned + aero_binned_delta.
107  subroutine aero_binned_add(aero_binned, aero_binned_delta)
108 
109  !> Base aero_binned_t structure that will be added to.
110  type(aero_binned_t), intent(inout) :: aero_binned
111  !> Structure to add to aero_binned.
112  type(aero_binned_t), intent(in) :: aero_binned_delta
113 
114  aero_binned%num_conc = aero_binned%num_conc + aero_binned_delta%num_conc
115  aero_binned%vol_conc = aero_binned%vol_conc + aero_binned_delta%vol_conc
116 
117  end subroutine aero_binned_add
118 
119 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
120 
121  !> Add a scaled \c aero_binned_t structure to an existing one.
122  !!
123  !! Symbolically does aero_binned = aero_binned + alpha * aero_binned_delta.
124  subroutine aero_binned_add_scaled(aero_binned, aero_binned_delta, alpha)
125 
126  !> Base aero_binned_t structure that will be added to.
127  type(aero_binned_t), intent(inout) :: aero_binned
128  !> Structure to add to aero_binned.
129  type(aero_binned_t), intent(in) :: aero_binned_delta
130  !> Scale factor.
131  real(kind=dp), intent(in) :: alpha
132 
133  aero_binned%num_conc = aero_binned%num_conc &
134  + alpha * aero_binned_delta%num_conc
135  aero_binned%vol_conc = aero_binned%vol_conc &
136  + alpha * aero_binned_delta%vol_conc
137 
138  end subroutine aero_binned_add_scaled
139 
140 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
141 
142  !> Subtract one aero_binned_t structure from another.
143  !!
144  !! Symbolically does aero_binned = aero_binned - aero_binned_delta.
145  subroutine aero_binned_sub(aero_binned, aero_binned_delta)
146 
147  !> Base aero_binned_t structure that will be subtracted from.
148  type(aero_binned_t), intent(inout) :: aero_binned
149  !> Structure to subtract from aero_binned.
150  type(aero_binned_t), intent(in) :: aero_binned_delta
151 
152  aero_binned%num_conc = aero_binned%num_conc - aero_binned_delta%num_conc
153  aero_binned%vol_conc = aero_binned%vol_conc - aero_binned_delta%vol_conc
154 
155  end subroutine aero_binned_sub
156 
157 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
158 
159  !> Scale an aero_binned_t by a real number.
160  !!
161  !! Symbolically does aero_binned = aero_binned * alpha.
162  subroutine aero_binned_scale(aero_binned, alpha)
163 
164  !> Base aero_binned to scale.
165  type(aero_binned_t), intent(inout) :: aero_binned
166  !> Scale factor.
167  real(kind=dp), intent(in) :: alpha
168 
169  aero_binned%num_conc = aero_binned%num_conc * alpha
170  aero_binned%vol_conc = aero_binned%vol_conc * alpha
171 
172  end subroutine aero_binned_scale
173 
174 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
175 
176  !> Scales an aero_binned_t element-wise by an array of reals.
177  subroutine aero_binned_scale_by_array(aero_binned, alpha_array)
178 
179  !> Base aero_binned_t structure that will be scaled.
180  type(aero_binned_t), intent(inout) :: aero_binned
181  !> Structure to scale aero_binned.
182  real(kind=dp), allocatable, intent(in) :: alpha_array(:)
183 
184  integer :: i
185 
186  do i = 1, size(aero_binned%num_conc)
187  aero_binned%num_conc(i) = alpha_array(i)*aero_binned%num_conc(i)
188  aero_binned%vol_conc(i,:) = alpha_array(i)*aero_binned%vol_conc(i,:)
189  end do
190 
191  end subroutine aero_binned_scale_by_array
192 
193 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
194 
195  !> Copy one aero_binned_t structure to another.
196  !!
197  !! Symbolically does aero_binned_to = aero_binned_from.
198  subroutine aero_binned_copy(aero_binned_from, aero_binned_to)
199 
200  !> Base aero_binned_t structure to copy from.
201  type(aero_binned_t), intent(in) :: aero_binned_from
202  !> Structure to copy to.
203  type(aero_binned_t), intent(inout) :: aero_binned_to
204 
205  integer :: n_bin, n_spec
206 
207  n_bin = size(aero_binned_from%vol_conc, 1)
208  n_spec = size(aero_binned_from%vol_conc, 2)
209  call aero_binned_deallocate(aero_binned_to)
210  call aero_binned_allocate_size(aero_binned_to, n_bin, n_spec)
211  aero_binned_to%num_conc = aero_binned_from%num_conc
212  aero_binned_to%vol_conc = aero_binned_from%vol_conc
213 
214  end subroutine aero_binned_copy
215 
216 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
217 
218  !> Add an aero_dist_t to an aero_binned_t.
219  !!
220  !! Symbolically does aero_binned = aero_binned + aero_dist.
221  subroutine aero_binned_add_aero_dist(aero_binned, bin_grid, &
222  aero_data, aero_dist)
223 
224  !> Base aero_binned_t structure to add to.
225  type(aero_binned_t), intent(inout) :: aero_binned
226  !> Bin grid.
227  type(bin_grid_t), intent(in) :: bin_grid
228  !> Aerosol data.
229  type(aero_data_t), intent(in) :: aero_data
230  !> The aero_dist_t structure to add.
231  type(aero_dist_t), intent(in) :: aero_dist
232 
233  real(kind=dp) :: dist_num_conc(size(aero_binned%num_conc, 1))
234  real(kind=dp) :: dist_vol_conc(size(aero_binned%vol_conc, 1), &
235  size(aero_binned%vol_conc, 2))
236 
237  call aero_dist_num_conc(aero_dist, bin_grid, aero_data, &
238  dist_num_conc)
239  call aero_dist_vol_conc(aero_dist, bin_grid, aero_data, &
240  dist_vol_conc)
241  aero_binned%num_conc = aero_binned%num_conc + dist_num_conc
242  aero_binned%vol_conc = aero_binned%vol_conc + dist_vol_conc
243 
244  end subroutine aero_binned_add_aero_dist
245 
246 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
247 
248  !> Determine the number of bytes required to pack the structure.
249  !!
250  !! See pmc_mpi for usage details.
251  integer function pmc_mpi_pack_size_aero_binned(val)
252 
253  !> Structure to pack.
254  type(aero_binned_t), intent(in) :: val
255 
257  pmc_mpi_pack_size_real_array(val%num_conc) &
258  + pmc_mpi_pack_size_real_array_2d(val%vol_conc)
259 
260  end function pmc_mpi_pack_size_aero_binned
261 
262 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
263 
264  !> Pack the structure into the buffer and advance position.
265  !!
266  !! See pmc_mpi for usage details.
267  subroutine pmc_mpi_pack_aero_binned(buffer, position, val)
268 
269  !> Memory buffer.
270  character, intent(inout) :: buffer(:)
271  !> Current buffer position.
272  integer, intent(inout) :: position
273  !> Structure to pack.
274  type(aero_binned_t), intent(in) :: val
275 
276 #ifdef PMC_USE_MPI
277  integer :: prev_position
278 
279  prev_position = position
280  call pmc_mpi_pack_real_array(buffer, position, val%num_conc)
281  call pmc_mpi_pack_real_array_2d(buffer, position, val%vol_conc)
282  call assert(348207873, &
283  position - prev_position <= pmc_mpi_pack_size_aero_binned(val))
284 #endif
285 
286  end subroutine pmc_mpi_pack_aero_binned
287 
288 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
289 
290  !> Unpack the structure from the buffer and advance position.
291  !!
292  !! See pmc_mpi for usage details.
293  subroutine pmc_mpi_unpack_aero_binned(buffer, position, val)
294 
295  !> Memory buffer.
296  character, intent(inout) :: buffer(:)
297  !> Current buffer position.
298  integer, intent(inout) :: position
299  !> Structure to unpack into (must not be allocated).
300  type(aero_binned_t), intent(inout) :: val
301 
302 #ifdef PMC_USE_MPI
303  integer :: prev_position
304 
305  prev_position = position
306  call pmc_mpi_unpack_real_array(buffer, position, val%num_conc)
307  call pmc_mpi_unpack_real_array_2d(buffer, position, val%vol_conc)
308  call assert(878267066, &
309  position - prev_position <= pmc_mpi_pack_size_aero_binned(val))
310 #endif
311 
312  end subroutine pmc_mpi_unpack_aero_binned
313 
314 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
315 
316  !> Computes the average of the structure across all processes,
317  !> storing the result on the root process.
318  subroutine pmc_mpi_reduce_avg_aero_binned(val, val_avg)
319 
320  !> Per-process value to average.
321  type(aero_binned_t), intent(in) :: val
322  !> Averaged result (only valid on root process).
323  type(aero_binned_t), intent(inout) :: val_avg
324 
325  call pmc_mpi_reduce_avg_real_array(val%num_conc, val_avg%num_conc)
326  call pmc_mpi_reduce_avg_real_array_2d(val%vol_conc, val_avg%vol_conc)
327 
328  end subroutine pmc_mpi_reduce_avg_aero_binned
329 
330 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
331 
332  !> Write full state.
333  subroutine aero_binned_output_netcdf(aero_binned, ncid, bin_grid, &
334  aero_data)
335 
336  !> Aero_binned to write.
337  type(aero_binned_t), intent(in) :: aero_binned
338  !> NetCDF file ID, in data mode.
339  integer, intent(in) :: ncid
340  !> bin_grid structure.
341  type(bin_grid_t), intent(in) :: bin_grid
342  !> aero_data structure.
343  type(aero_data_t), intent(in) :: aero_data
344 
345  integer :: dimid_aero_diam, dimid_aero_species
346  real(kind=dp) :: mass_conc(bin_grid%n_bin, aero_data%n_spec)
347  integer :: i_bin
348 
349  !> \page output_format_aero_binned Output File Format: Aerosol Binned Sectional State
350  !!
351  !! The aerosol size distributions (number and mass) are stored on
352  !! a logarmithmic grid (see the \ref output_format_diam_bin_grid
353  !! section). To compute the total number or mass concentration,
354  !! compute the sum over \c i of <tt>aero_number_concentration(i) *
355  !! aero_diam_widths(i)</tt>, for example.
356  !!
357  !! The aerosol binned sectional state uses the \c aero_species
358  !! NetCDF dimension as specified in the \ref
359  !! output_format_aero_data section, as well as the \c aero_diam
360  !! NetCDF dimension specified in the \ref
361  !! output_format_diam_bin_grid section.
362  !!
363  !! The aerosol binned sectional state NetCDF variables are:
364  !! - \b aero_number_concentration (unit 1/m^3, dim \c aero_diam): the
365  !! number size distribution for the aerosol population,
366  !! \f$ dN(r)/d\ln r \f$, per bin
367  !! - \b aero_mass_concentration (unit kg/m^3, dim
368  !! <tt>dimid_aero_diam x dimid_aero_species</tt>): the mass size
369  !! distribution for the aerosol population,
370  !! \f$ dM(r,s)/d\ln r \f$, per bin and per species
371 
372  do i_bin = 1,bin_grid%n_bin
373  mass_conc(i_bin,:) = aero_binned%vol_conc(i_bin,:) * aero_data%density
374  end do
375 
376  call bin_grid_netcdf_dim(bin_grid, ncid, "aero_diam", "m", &
377  dimid_aero_diam, "aerosol diameter", scale=2d0)
378  call aero_data_netcdf_dim_aero_species(aero_data, ncid, dimid_aero_species)
379 
380  call pmc_nc_write_real_1d(ncid, aero_binned%num_conc, &
381  "aero_number_concentration", (/ dimid_aero_diam /), &
382  unit="1/m^3", &
383  long_name="aerosol number size concentration distribution", &
384  description="logarithmic number size concentration, " &
385  // "d N(r)/d ln D --- multiply by aero_diam_widths(i) " &
386  // "and sum over i to compute the total number concentration")
387  call pmc_nc_write_real_2d(ncid, mass_conc, &
388  "aero_mass_concentration", &
389  (/ dimid_aero_diam, dimid_aero_species /), unit="kg/m^3", &
390  long_name="aerosol mass size concentration distribution", &
391  description="logarithmic mass size concentration per species, " &
392  // "d M(r,s)/d ln D --- multiply by aero_diam_widths(i) " &
393  // "and sum over i to compute the total mass concentration of " &
394  // "species s")
395 
396  end subroutine aero_binned_output_netcdf
397 
398  ! output_format_diam_bin_grid is here, as this is the only place it's used
399 
400  ! this belongs in the subroutine above, but is outside because
401  ! Doxygen 1.8.7 doesn't resolve references when multiple \page
402  ! blocks are in one subroutine
403 
404  !> \page output_format_diam_bin_grid Output File Format: Diameter Bin Grid Data
405  !!
406  !! The aerosol diameter bin grid data NetCDF dimensions are:
407  !! - \b aero_diam: number of bins (grid cells) on the diameter axis
408  !! - \b aero_diam_edges: number of bin edges (grid cell edges) on
409  !! the diameter axis --- always equal to <tt>aero_diam + 1</tt>
410  !!
411  !! The aerosol diameter bin grid data NetCDF variables are:
412  !! - \b aero_diam (unit m, dim \c aero_diam): aerosol diameter axis
413  !! bin centers --- centered on a logarithmic scale from the edges, so
414  !! that <tt>aero_diam(i) / aero_diam_edges(i) =
415  !! sqrt(aero_diam_edges(i+1) / aero_diam_edges(i))</tt>
416  !! - \b aero_diam_edges (unit m, dim \c aero_diam_edges): aersol
417  !! diameter axis bin edges (there is one more edge than center)
418  !! - \b aero_diam_widths (dimensionless, dim \c aero_diam):
419  !! the base-e logarithmic bin widths --- <tt>aero_diam_widths(i)
420  !! = ln(aero_diam_edges(i+1) / aero_diam_edges(i))</tt>, so
421  !! all bins have the same width
422  !!
423  !! See also:
424  !! - \ref input_format_diam_bin_grid --- the corresponding input format
425 
426 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
427 
428  !> Read full state.
429  subroutine aero_binned_input_netcdf(aero_binned, ncid, bin_grid, &
430  aero_data)
431 
432  !> Aero_binned to write.
433  type(aero_binned_t), intent(inout) :: aero_binned
434  !> NetCDF file ID, in data mode.
435  integer, intent(in) :: ncid
436  !> bin_grid structure.
437  type(bin_grid_t), intent(in) :: bin_grid
438  !> aero_data structure.
439  type(aero_data_t), intent(in) :: aero_data
440 
441  real(kind=dp) :: mass_conc(bin_grid%n_bin, aero_data%n_spec)
442  integer :: i_bin
443 
444  call aero_binned_deallocate(aero_binned)
445  call aero_binned_allocate_size(aero_binned, bin_grid%n_bin, &
446  aero_data%n_spec)
447 
448  call pmc_nc_read_real_1d(ncid, aero_binned%num_conc, &
449  "aero_number_concentration")
450  call pmc_nc_read_real_2d(ncid, mass_conc, "aero_mass_concentration")
451 
452  do i_bin = 1,bin_grid%n_bin
453  aero_binned%vol_conc(i_bin,:) = mass_conc(i_bin,:) / aero_data%density
454  end do
455 
456  end subroutine aero_binned_input_netcdf
457 
458 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
459 
460 end module pmc_aero_binned
subroutine aero_binned_allocate_size(aero_binned, n_bin, n_spec)
Allocate an aero_binned_t of the given size.
Definition: aero_binned.F90:60
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
The aero_data_t structure and associated subroutines.
Definition: aero_data.F90:9
subroutine pmc_mpi_pack_real_array_2d(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:742
subroutine pmc_mpi_unpack_aero_binned(buffer, position, val)
Unpack the structure from the buffer and advance position.
The aero_dist_t structure and associated subroutines.
Definition: aero_dist.F90:18
subroutine aero_binned_sub(aero_binned, aero_binned_delta)
Subtract one aero_binned_t structure from another.
The aero_particle_t structure and associated subroutines.
subroutine aero_binned_deallocate(aero_binned)
Free internal memory in an aero_binned_t structure.
Definition: aero_binned.F90:79
subroutine aero_binned_allocate(aero_binned)
Allocate an aero_binned_t.
Definition: aero_binned.F90:47
subroutine bin_grid_netcdf_dim(bin_grid, ncid, dim_name, unit, dimid, long_name, scale)
Write a bin grid dimension to the given NetCDF file if it is not already present and in any case retu...
Definition: bin_grid.F90:438
subroutine pmc_mpi_reduce_avg_real_array(val, val_avg)
Computes the average of val across all processes, storing the result in val_avg on the root process...
Definition: mpi.F90:1163
A complete aerosol distribution, consisting of several modes.
Definition: aero_dist.F90:33
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 aero_binned_add(aero_binned, aero_binned_delta)
Add two aero_binned_t structures together.
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_nc_read_real_2d(ncid, var, name, must_be_present)
Read a simple real 2D array from a NetCDF file.
Definition: netcdf.F90:288
subroutine aero_binned_scale(aero_binned, alpha)
Scale an aero_binned_t by a real number.
Wrapper functions for MPI.
Definition: mpi.F90:13
subroutine pmc_mpi_pack_aero_binned(buffer, position, val)
Pack the structure into the buffer and advance position.
integer function pmc_mpi_pack_size_real_array_2d(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:511
1D grid, either logarithmic or linear.
Definition: bin_grid.F90:33
The bin_grid_t structure and associated subroutines.
Definition: bin_grid.F90:9
subroutine pmc_mpi_unpack_real_array_2d(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:983
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
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
The aero_binned_t structure and associated subroutines.
Definition: aero_binned.F90:9
subroutine pmc_mpi_reduce_avg_aero_binned(val, val_avg)
Computes the average of the structure across all processes, storing the result on the root process...
subroutine aero_binned_scale_by_array(aero_binned, alpha_array)
Scales an aero_binned_t element-wise by an array of reals.
subroutine aero_binned_input_netcdf(aero_binned, ncid, bin_grid, aero_data)
Read full state.
subroutine aero_binned_zero(aero_binned)
Set all internal data in an aero_binned_t structure to zero.
Definition: aero_binned.F90:92
subroutine pmc_mpi_pack_real_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:688
integer function pmc_mpi_pack_size_aero_binned(val)
Determine the number of bytes required to pack the structure.
subroutine aero_binned_add_aero_dist(aero_binned, bin_grid, aero_data, aero_dist)
Add an aero_dist_t to an aero_binned_t.
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_binned_copy(aero_binned_from, aero_binned_to)
Copy one aero_binned_t structure to another.
subroutine pmc_mpi_reduce_avg_real_array_2d(val, val_avg)
Computes the average of val across all processes, storing the result in val_avg on the root process...
Definition: mpi.F90:1190
subroutine aero_binned_add_scaled(aero_binned, aero_binned_delta, alpha)
Add a scaled aero_binned_t structure to an existing one.
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
Aerosol number and volume distributions stored per bin.
Definition: aero_binned.F90:33