CAMP 1.0.0
Chemistry Across Multiple Phases
aero_rep_modal_binned_mass.F90
Go to the documentation of this file.
1! Copyright (C) 2021 Barcelona Supercomputing Center and University of
2! Illinois at Urbana-Champaign
3! SPDX-License-Identifier: MIT
4
5!> \file
6!> The camp_aero_rep_modal_binned_mass module.
7
8!> \page camp_aero_rep_modal_binned_mass CAMP: Modal/Binned Mass Aerosol Representation
9!!
10!! The modal/binned mass aerosol representation includes a set of sections/bins
11!! that are made up of one or more \ref camp_aero_phase "aerosol phases." The
12!! \c json object for this \ref camp_aero_rep "aerosol representation" has the
13!! following format:
14!! \code{.json}
15!! { "camp-data" : [
16!! {
17!! "name" : "my modal/binned aero rep",
18!! "type" : "AERO_REP_MODAL_BINNED_MASS",
19!! "modes/bins" :
20!! {
21!! "dust" :
22!! {
23!! "type" : "BINNED",
24!! "phases" : [ "insoluble", "organic", "aqueous" ],
25!! "bins" : 8,
26!! "minimum diameter [m]" : 0.8e-9,
27!! "maximum diameter [m]" : 1.0e-6,
28!! "scale" : "LOG"
29!! },
30!! "depeche" :
31!! {
32!! "type" : "MODAL",
33!! "phases" : [ "moody", "listless" ],
34!! "shape" : "LOG_NORMAL",
35!! }
36!! }
37!! },
38!! ...
39!! ]}
40!! \endcode
41!! The key-value pair \b type is required and must be
42!! \b AERO_REP_MODAL_BINNED_MASS. The key-value pair \b modes/bins is also
43!! required and must contain a set of at least one uniquely named mode or
44!! bin-set key-value pair whose value(s) specify a \b type that must be either
45!! \b MODAL or \b BINNED and an array of \b phases that correspond to existing
46!! \ref camp_aero_phase "aerosol phase" objects. Each phase will be present
47!! once within a mode or once within each bin in a bin-set.
48!!
49!! Modes must also specify a distribution \b shape which must be \b LOG_NORMAL
50!! (the available shapes may be expanded in the future). Log-normal sections
51!! have \b geometric \b mean \b diameter (m) and a \b geometric
52!! \b standard \b deviation (unitless) properties that must be set using an
53!! \c aero_rep_update_data_modal_binned_mass_GMD_t or and
54!! \c aero_rep_update_data_modal_binned_mass_GSD_t object and will be used along
55!! with the mass
56!! concentration of species in each phase and their densities to calculate a
57!! lognormal distribution for each mode at runtime.
58!!
59!! Bin sets must specify the number of \b bins, a \b minimum \b diameter (m),
60!! a \b maximum \b diameter (m) and a \b scale, which must be \b LOG or
61!! \b LINEAR. The number concentration will be calculated at run-time based on
62!! the total mass of each bin, the species densities and the diameter of
63!! particles in that bin.
64!!
65!! The GMD and GSD for each mode must be set from an external model using
66!! \c camp_aero_rep_modal_binned_mass::aero_rep_update_data_modal_binned_mass_GMD_t
67!! and
68!! \c camp_aero_rep_modal_binned_mass::aero_rep_update_data_modal_binned_mass_GSD_t
69!! objects.
70
71!> The abstract aero_rep_modal_binned_mass_t structure and associated subroutines.
73
78 use camp_mpi
80 use camp_util, only: dp, i_kind, &
83
84 use iso_c_binding
85
86 implicit none
87 private
88
89#define BINNED 1
90#define MODAL 2
91
92#define NUM_SECTION_ this%condensed_data_int(1)
93#define INT_DATA_SIZE_ this%condensed_data_int(2)
94#define REAL_DATA_SIZE_ this%condensed_data_int(3)
95#define AERO_REP_ID_ this%condensed_data_int(4)
96#define NUM_INT_PROP_ 4
97#define NUM_REAL_PROP_ 0
98#define NUM_ENV_PARAM_ 0
99#define MODE_INT_PROP_LOC_(x) this%condensed_data_int(NUM_INT_PROP_+x)
100#define MODE_REAL_PROP_LOC_(x) this%condensed_data_int(NUM_INT_PROP_+NUM_SECTION_+x)
101#define SECTION_TYPE_(x) this%condensed_data_int(MODE_INT_PROP_LOC_(x))
102
103! For modes, NUM_BINS_ = 1
104#define NUM_BINS_(x) this%condensed_data_int(MODE_INT_PROP_LOC_(x)+1)
105
106! Number of aerosol phases in this mode/bin set
107#define NUM_PHASE_(x) this%condensed_data_int(MODE_INT_PROP_LOC_(x)+2)
108
109! Phase state and model data ids
110#define PHASE_STATE_ID_(x,y,b) this%condensed_data_int(MODE_INT_PROP_LOC_(x)+2+(b-1)*NUM_PHASE_(x)+y)
111#define PHASE_MODEL_DATA_ID_(x,y,b) this%condensed_data_int(MODE_INT_PROP_LOC_(x)+2+NUM_BINS_(x)*NUM_PHASE_(x)+(b-1)*NUM_PHASE_(x)+y)
112
113! Number of Jacobian elements in a phase
114#define PHASE_NUM_JAC_ELEM_(x,y,b) this%condensed_data_int(MODE_INT_PROP_LOC_(x)+2+2*NUM_BINS_(x)*NUM_PHASE_(x)+(b-1)*NUM_PHASE_(x)+y)
115
116! Bin diameter (bins only)
117#define BIN_DP_(x,b) this%condensed_data_real(MODE_REAL_PROP_LOC_(x)+(b-1)*3)
118
119! Real-time number concetration - used for modes and bins - for modes, b=1
120#define NUMBER_CONC_(x,b) this%condensed_data_real(MODE_REAL_PROP_LOC_(x)+(b-1)*3+1)
121
122! Real-time effective radius - for modes, b=1
123#define EFFECTIVE_RADIUS_(x,b) this%condensed_data_real(MODE_REAL_PROP_LOC_(x)+(b-1)*3+2)
124
125! Real-time aerosol phase mass - used for modes and bins - for modes, b=1
126#define PHASE_MASS_(x,y,b) this%condensed_data_real(MODE_REAL_PROP_LOC_(x)+3*NUM_BINS_(x)+(b-1)*NUM_PHASE_(x)+y-1)
127
128! Real-time aerosol phase average MW - used for modes and bins - for modes, b=0
129#define PHASE_AVG_MW_(x,y,b) this%condensed_data_real(MODE_REAL_PROP_LOC_(x)+(3+NUM_PHASE_(x))*NUM_BINS_(x)+(b-1)*NUM_PHASE_(x)+y-1)
130
131 ! Update types (These must match values in aero_rep_modal_binned_mass.c)
132 integer(kind=i_kind), parameter, public :: update_gmd = 0
133 integer(kind=i_kind), parameter, public :: update_gsd = 1
134
138
139 !> Modal mass aerosol representation
140 !!
141 !! Time-invariant data related to a modal/binned mass aerosol representation.
143 !> Mode names (only used during initialization)
144 type(string_t), allocatable :: section_name(:)
145 !> Phase state id (only used during initialization)
146 integer(kind=i_kind), allocatable :: phase_state_id(:)
147 contains
148 !> Initialize the aerosol representation data, validating component data and
149 !! loading any required information from the \c
150 !! aero_rep_data_t::property_set. This routine should be called once for
151 !! each aerosol representation at the beginning of a model run after all
152 !! the input files have been read in. It ensures all data required during
153 !! the model run are included in the condensed data arrays.
154 procedure :: initialize
155 !> Initialize an update data GSD object
156 procedure :: update_data_initialize_gsd => update_data_init_gsd
157 !> Initialize an update data GMD object
158 procedure :: update_data_initialize_gmd => update_data_init_gmd
159 !> Get an id for a mode or bin in the aerosol representation by name for
160 !! use with updates from external modules
161 procedure :: get_section_id
162 !> Get the size of the section of the
163 !! \c camp_camp_state::camp_state_t::state_var array required for this
164 !! aerosol representation.
165 !!
166 !! For a modal/binned mass representation, the size will correspond to the
167 !! the sum of the sizes of a single instance of each aerosol phase
168 !! provided to \c aero_rep_modal_binned_mass::initialize()
169 procedure :: size => get_size
170 !> Get a list of unique names for each element on the
171 !! \c camp_camp_state::camp_state_t::state_var array for this aerosol
172 !! representation. The list may be restricted to a particular phase and/or
173 !! aerosol species by including the phase_name and spec_name arguments.
174 !!
175 !! For a modal/binned mass representation, the unique names for bins are:
176 !! - "bin name.bin #.phase name.species name"
177 !!
178 !! ... and for modes are:
179 !! - "mode name.phase name.species name"
180 procedure :: unique_names
181 !> Get a species id on the \c camp_camp_state::camp_state_t::state_var
182 !! array by its unique name. These are unique ids for each element on the
183 !! state array for this \ref camp_aero_rep "aerosol representation" and
184 !! are numbered:
185 !!
186 !! \f[x_u \in x_f ... (x_f+n-1)\f]
187 !!
188 !! where \f$x_u\f$ is the id of the element corresponding to the species
189 !! with unique name \f$u\f$ on the \c
190 !! camp_camp_state::camp_state_t::state_var array, \f$x_f\f$ is the index
191 !! of the first element for this aerosol representation on the state array
192 !! and \f$n\f$ is the total number of variables on the state array from
193 !! this aerosol representation.
194 procedure :: spec_state_id
195 !> Get the non-unique name of a species by its unique name
196 procedure :: spec_name
197 !> Get the number of instances of an aerosol phase in this representation
199 !> Get the number of Jacobian elements used in calculations of aerosol mass,
200 !! volume, number, etc. for a particular phase
201 procedure :: num_jac_elem
202 !> Returns index_pair_t type with phase_ids of adjacent phases
203 !! for modal/binned representation there are no adjacent phases
204 procedure :: adjacent_phases
205 !> Finalize the aerosol representation
207
209
210 ! Constructor for aero_rep_modal_binned_mass_t
212 procedure :: constructor
213 end interface aero_rep_modal_binned_mass_t
214
215 !> Update GMD object
216 type, extends(aero_rep_update_data_t) :: &
218 private
219 !> Flag indicating whether the update data has been allocated
220 logical :: is_malloced = .false.
221 !> Unique id for finding aerosol representations during initialization
222 integer(kind=i_kind) :: aero_rep_unique_id = 0
223 contains
224 !> Update the GMD
225 procedure :: set_gmd => update_data_set_gmd
226 !> Determine the pack size of the local update data
228 !> Pack the local update data to a binary
230 !> Unpack the local update data from a binary
232 !> Finalize the GMD update data
235
236 !> Update GSD object
237 type, extends(aero_rep_update_data_t) :: &
239 private
240 !> Flag indicating whether the update data has been allocated
241 logical :: is_malloced = .false.
242 !> Unique id for finding aerosol representations during initialization
243 integer(kind=i_kind) :: aero_rep_unique_id = 0
244 contains
245 !> Update the GSD
246 procedure :: set_gsd => update_data_set_gsd
247 !> Determine the pack size of the local update data
249 !> Pack the local update data to a binary
251 !> Unpack the local update data from a binary
253 !> Finalize the GSD update data
256
257 !> Interface to c aerosol representation functions
258 interface
259
260 !> Allocate space for a GMD update object
262 result(update_data) bind (c)
263 use iso_c_binding
264 !> Allocated update data object
265 type(c_ptr) :: update_data
267
268 !> Set a new mode GMD
270 aero_rep_unique_id, section_id, gmd) bind (c)
271 use iso_c_binding
272 !> Update data
273 type(c_ptr), value :: update_data
274 !> Aerosol representation unique id
275 integer(kind=c_int), value :: aero_rep_unique_id
276 !> Section id from
277 !! camp_aero_rep_modal_binned_mass::aero_rep_modal_binned_mass_t::get_section_id
278 integer(kind=c_int), value :: section_id
279 !> New GMD (m)
280 real(kind=c_double), value :: gmd
282
283 !> Allocate space for a GSD update object
285 result(update_data) bind (c)
286 use iso_c_binding
287 !> Allocated update data object
288 type(c_ptr) :: update_data
290
291 !> Set a new mode GSD
293 aero_rep_unique_id, section_id, gsd) bind (c)
294 use iso_c_binding
295 !> Update data
296 type(c_ptr), value :: update_data
297 !> Aerosol representation unique id
298 integer(kind=c_int), value :: aero_rep_unique_id
299 !> Section id from
300 !! camp_aero_rep_modal_binned_mass::aero_rep_modal_binned_mass_t::get_section_id
301 integer(kind=c_int), value :: section_id
302 !> New GSD (m)
303 real(kind=c_double), value :: gsd
305
306 !> Free an update data object
307 pure subroutine aero_rep_free_update_data(update_data) bind (c)
308 use iso_c_binding
309 !> Update data
310 type(c_ptr), value, intent(in) :: update_data
311 end subroutine aero_rep_free_update_data
312
313 end interface
314
315contains
316
317!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
318
319 !> Constructor for aero_rep_modal_binned_mass_t
320 function constructor() result (new_obj)
321
322 !> New aerosol representation
323 type(aero_rep_modal_binned_mass_t), pointer :: new_obj
324
325 allocate(new_obj)
326
327 end function constructor
328
329!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
330
331 !> Initialize the aerosol representation data, validating component data and
332 !! loading any required information from the \c
333 !! aero_rep_data_t::property_set. This routine should be called once for
334 !! each aerosol representation at the beginning of a model run after all
335 !! the input files have been read in. It ensures all data required during
336 !! the model run are included in the condensed data arrays.
337 subroutine initialize(this, aero_phase_set, spec_state_id)
338
339 !> Aerosol representation data
340 class(aero_rep_modal_binned_mass_t), intent(inout) :: this
341 !> The set of aerosol phases
342 type(aero_phase_data_ptr), pointer, intent(in) :: aero_phase_set(:)
343 !> Beginning state id for this aerosol representationin the model species
344 !! state array
345 integer(kind=i_kind), intent(in) :: spec_state_id
346
347 type(property_t), pointer :: sections, section, phases
348 integer(kind=i_kind) :: i_section, i_phase, j_phase, k_phase, &
349 i_bin
350 integer(kind=i_kind) :: curr_spec_state_id
351 integer(kind=i_kind) :: num_phase, num_bin
352 integer(kind=i_kind) :: n_int_param, n_float_param
353 character(len=:), allocatable :: key_name, phase_name, sect_type, str_val
354 real(kind=dp) :: min_dp, max_dp, d_log_dp
355
356 ! Determine the size of the condensed data arrays
357 n_int_param = num_int_prop_
358 n_float_param = num_real_prop_
359
360 ! Get the set of sections/bin-sets
361 key_name = "modes/bins"
362 call assert_msg(877855909, &
363 this%property_set%get_property_t(key_name, sections), &
364 "Missing sections/bins for modal/binned mass aerosol "// &
365 "representation '"//this%rep_name//"'")
366 call assert_msg(894962494, sections%size().gt.0, "No sections or bins "// &
367 "specified for modal/binned mass aerosol representation '"// &
368 this%rep_name//"'")
369
370 ! Allocate space for the mode/bin names
371 allocate(this%section_name(sections%size()))
372
373 ! Loop through the sections, adding names and counting the spaces needed
374 ! on the condensed data arrays, and counting the total phases instances
375 num_phase = 0
376 call sections%iter_reset()
377 do i_section = 1, sections%size()
378
379 ! Get the mode/bin name
380 call assert(867378489, sections%get_key(key_name))
381 call assert_msg(234513113, len(key_name).gt.0, "Missing mode/bin "// &
382 "name in modal/binned mass aerosol representation '"// &
383 this%rep_name//"'")
384 this%section_name(i_section)%string = key_name
385
386 ! Get the mode/bin properties
387 call assert_msg(517138327, sections%get_property_t(val=section), &
388 "Invalid structure for mode/bin '"// &
389 this%section_name(i_section)%string// &
390 "' in modal/binned mass aerosol representation '"// &
391 this%rep_name//"'")
392
393 ! Get the section type
394 key_name = "type"
395 call assert_msg(742404898, section%get_string(key_name, sect_type), &
396 "Missing mode/bin type in mode/bin '"// &
397 this%section_name(i_section)%string// &
398 "' in modal/binned mass aerosol representation '"// &
399 this%rep_name//"'")
400 call assert_msg(618556995, &
401 sect_type.eq."MODAL".or.sect_type.eq."BINNED", &
402 "Invalid mode/bin type '"//sect_type//"' in mode/bin '"// &
403 this%section_name(i_section)%string// &
404 "' in modal/binned mass aerosol representation '"// &
405 this%rep_name//"'")
406
407 ! Get the number of bins (or set to 1 for a mode)
408 num_bin = 1
409 if (sect_type.eq."BINNED") then
410
411 key_name = "bins"
412 call assert_msg(824494286, section%get_int(key_name, num_bin), &
413 "Missing number of bins in bin '"// &
414 this%section_name(i_section)%string// &
415 "' in modal/binned mass aerosol representation '"// &
416 this%rep_name//"'")
417 end if
418
419 ! Add space for the mode/bin type, number of bins, and phase count
420 ! and parameter locations
421 n_int_param = n_int_param + 5
422
423 ! Add space for the bin diameter, number concentration, and
424 ! effective radius
425 n_float_param = n_float_param + 3*num_bin
426
427 ! Get the set of phases
428 key_name = "phases"
429 call assert_msg(815518058, section%get_property_t(key_name, phases), &
430 "Missing phases for mode '"// &
431 this%section_name(i_section)%string// &
432 "' in modal/binned mass aerosol representation '"// &
433 this%rep_name//"'")
434
435 ! Add the phases to the counter
436 call assert_msg(772593427, phases%size().gt.0, &
437 "No phases specified for mode '"// &
438 this%section_name(i_section)%string// &
439 "' in modal/binned mass aerosol representation '"// &
440 this%rep_name//"'")
441 num_phase = num_phase + phases%size() * num_bin
442
443 ! Loop through the phases and make sure they exist
444 call phases%iter_reset()
445 do i_phase = 1, phases%size()
446
447 ! Get the phase name
448 call assert_msg(393427582, phases%get_string(val=phase_name), &
449 "Non-string phase name for mode '"// &
450 this%section_name(i_section)%string// &
451 "' in modal/binned mass aerosol representation '"// &
452 this%rep_name//"'")
453
454 ! Find the aerosol phase and add space for its variables
455 do j_phase = 1, size(aero_phase_set)
456 if (phase_name.eq.aero_phase_set(j_phase)%val%name()) then
457
458 ! Add space for the phase state, model data ids, and number
459 ! of Jacobian elements
460 n_int_param = n_int_param + 3 * num_bin
461
462 ! Add space for total aerosol phase mass and average MW,
463 n_float_param = n_float_param + 2 * num_bin
464
465 exit
466 else if (j_phase.eq.size(aero_phase_set)) then
467 call die_msg(652391420, "Non-existant aerosol phase '"// &
468 phase_name//"' specified for mode '"// &
469 this%section_name(i_section)%string// &
470 "' in modal/binned mass aerosol representation '"// &
471 this%rep_name//"'")
472 end if
473 end do
474
475 call phases%iter_next()
476 end do
477
478 call sections%iter_next()
479 end do
480
481 ! Allocate space for the aerosol phases and species state ids
482 allocate(this%aero_phase(num_phase))
483 allocate(this%aero_phase_is_at_surface(num_phase))
484 allocate(this%phase_state_id(size(this%aero_phase)))
485
486 ! Allocate condensed data arrays
487 allocate(this%condensed_data_int(n_int_param))
488 allocate(this%condensed_data_real(n_float_param))
489 this%condensed_data_int(:) = int(0, kind=i_kind)
490 this%condensed_data_real(:) = real(0.0, kind=dp)
491 int_data_size_ = n_int_param
492 real_data_size_ = n_float_param
493
494 ! Set the number of sections
495 num_section_ = sections%size()
496
497 ! Save space for the environment-dependent parameters (GMD and GSD)
498 this%num_env_params = 2 * num_section_
499
500 ! Loop through the sections, adding names and distribution parameters and
501 ! counting the phases in each section
502 i_phase = 1
503 curr_spec_state_id = spec_state_id
504 n_int_param = num_int_prop_+2*num_section_+1
505 n_float_param = num_real_prop_+1
506 call sections%iter_reset()
507 do i_section = 1, num_section_
508
509 ! Set the data locations for this mode
510 mode_int_prop_loc_(i_section) = n_int_param
511 mode_real_prop_loc_(i_section) = n_float_param
512
513 ! Get the mode/bin properties
514 call assert(394743663, sections%get_property_t(val=section))
515
516 ! Get the mode/bin type
517 key_name = "type"
518 call assert(667058653, section%get_string(key_name, sect_type))
519 if (sect_type.eq."MODAL") then
520 section_type_(i_section) = modal
521 else if (sect_type.eq."BINNED") then
522 section_type_(i_section) = binned
523 else
524 call die_msg(256924433, "Internal error")
525 end if
526
527 ! Get the number of bins (or set to 1 for a mode)
528 num_bins_(i_section) = 1
529 if (section_type_(i_section).eq.binned) then
530 key_name = "bins"
531 call assert(315215287, section%get_int(key_name, num_bins_(i_section)))
532 end if
533
534 ! Get mode parameters
535 if (section_type_(i_section).eq.modal) then
536
537 ! Currently no mode parameters
538
539 effective_radius_(i_section,1) = -9999.9
540 number_conc_(i_section,1) = -9999.9
541
542 ! Get bin parameters
543 else if (section_type_(i_section).eq.binned) then
544
545 ! Get the minimum diameter (m)
546 key_name = "minimum diameter [m]"
547 call assert_msg(548762180, section%get_real(key_name, min_dp), &
548 "Missing minimum diameter for bin '"// &
549 this%section_name(i_section)%string// &
550 "' in modal/binned mass aerosol representation '"// &
551 this%rep_name//"'")
552
553 ! Get the maximum diameter (m)
554 key_name = "maximum diameter [m]"
555 call assert_msg(288632226, section%get_real(key_name, max_dp), &
556 "Missing maximum diameter for bin '"// &
557 this%section_name(i_section)%string// &
558 "' in modal/binned mass aerosol representation '"// &
559 this%rep_name//"'")
560
561 ! Get the scale
562 key_name = "scale"
563 call assert_msg(404761639, section%get_string(key_name, str_val), &
564 "Missing bin scale for bin '"// &
565 this%section_name(i_section)%string// &
566 "' in modal/binned mass aerosol representation '"// &
567 this%rep_name//"'")
568
569 ! Assign the bin diameters
570 if (str_val.eq."LINEAR") then
571 do i_bin = 1, num_bins_(i_section)
572 bin_dp_(i_section,i_bin) = min_dp + &
573 (i_bin-1) * (max_dp-min_dp)/(num_bins_(i_section)-1)
574 end do
575 else if (str_val.eq."LOG") then
576 d_log_dp = (log10(max_dp)-log10(min_dp))/(num_bins_(i_section)-1)
577 do i_bin = 1, num_bins_(i_section)
578 bin_dp_(i_section,i_bin) = 10.0d0**( log10(min_dp) + &
579 (i_bin-1) * d_log_dp )
580 end do
581 else
582 call die_msg(236797392, "Invalid scale specified for bin '"// &
583 this%section_name(i_section)%string// &
584 "' in modal/binned mass aerosol representation '"// &
585 this%rep_name//"'")
586 end if
587 do i_bin = 1, num_bins_(i_section)
588 ! Set the effective radius
589 effective_radius_(i_section,i_bin) = bin_dp_(i_section,i_bin) / 2.0
590 number_conc_(i_section,i_bin) = -9999.9
591 end do
592 end if
593
594 ! Get the set of phases
595 key_name = "phases"
596 call assert(712411046, section%get_property_t(key_name, phases))
597
598 ! Save the number of phases
599 num_phase_(i_section) = phases%size()
600
601 ! Add space for the mode/bin type, number of bins, and number of phases
602 n_int_param = n_int_param + 3
603
604 ! Add space for bin diameter, number concentration and effective radius
605 n_float_param = n_float_param + 3 * num_bins_(i_section)
606
607 ! Loop through the phase names, look them up, and add them to the list
608 call phases%iter_reset()
609 do j_phase = 1, phases%size()
610
611 ! Get the phase name
612 call assert(775801035, phases%get_string(val=phase_name))
613
614 ! Find the aerosol phase and add it to the list
615 do k_phase = 1, size(aero_phase_set)
616 if (phase_name.eq.aero_phase_set(k_phase)%val%name()) then
617
618 ! Loop through the bins
619 do i_bin = 1, num_bins_(i_section)
620
621 ! Add the aerosol phase to the list
622 this%aero_phase(i_phase) = aero_phase_set(k_phase)
623
624 ! No species exist at surface
625 this%aero_phase_is_at_surface(i_phase) = .true.
626
627 ! Save the starting id for this phase on the state array
628 this%phase_state_id(i_phase) = curr_spec_state_id
629 phase_state_id_(i_section, j_phase, i_bin) = curr_spec_state_id
630
631 ! Increment the state id by the size of the phase
632 curr_spec_state_id = curr_spec_state_id + &
633 aero_phase_set(k_phase)%val%size()
634
635 ! Save the phase model data id
636 phase_model_data_id_(i_section, j_phase, i_bin) = k_phase
637
638 i_phase = i_phase + 1
639 end do
640
641 ! Add space for aerosol phase state, model data ids, and
642 ! number of Jacobian elements
643 n_int_param = n_int_param + 3*num_bins_(i_section)
644
645 ! Add space for aerosol phase mass and average MW
646 n_float_param = n_float_param + 2*num_bins_(i_section)
647
648 exit
649 else if (k_phase.eq.size(aero_phase_set)) then
650 call die_msg(652391420, "Internal error.")
651 end if
652 end do
653
654 call phases%iter_next()
655 end do
656
657 call sections%iter_next()
658 end do
659
660 ! Initialize the aerosol representation id
661 aero_rep_id_ = -1
662
663 ! Check the data sizes
664 call assert(951534966, i_phase-1.eq.num_phase)
665 call assert(951534966, n_int_param.eq.int_data_size_+1)
666 call assert(325387136, n_float_param.eq.real_data_size_+1)
667
668 end subroutine initialize
669
670!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
671
672 !> Get an id for a mode or bin by name for use with updates from external
673 !! modules
674 function get_section_id(this, section_name, section_id) result (found)
675
676 !> Flag indicating whether the mode/bin was found
677 logical :: found
678 !> Aerosol representation
679 class(aero_rep_modal_binned_mass_t), intent(in) :: this
680 !> Section name
681 character(len=*), intent(in) :: section_name
682 !> Section id
683 integer(kind=i_kind), intent(out) :: section_id
684
685 integer(kind=i_kind) :: i_section
686
687 call assert_msg(194186171, len(trim(section_name)).gt.0, &
688 "Trying to get section id of unnamed aerosol "// &
689 "representation.")
690
691 found = .false.
692 do i_section = 1, size(this%section_name)
693 if (this%section_name(i_section)%string.eq.trim(section_name)) then
694 found = .true.
695 section_id = i_section
696 return
697 end if
698 end do
699
700 end function get_section_id
701
702!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
703
704 !> Get the size of the section of the
705 !! \c camp_camp_state::camp_state_t::state_var array required for this
706 !! aerosol representation.
707 !!
708 !! For a modal/binned mass representation, the size will correspond to the
709 !! the sum of the sizes of a single instance of each aerosol phase
710 !! provided to \c aero_rep_modal_binned_mass::initialize()
711 function get_size(this) result (state_size)
712
713 !> Size on the state array
714 integer(kind=i_kind) :: state_size
715 !> Aerosol representation data
716 class(aero_rep_modal_binned_mass_t), intent(in) :: this
717
718 integer(kind=i_kind) :: i_phase
719
720 ! Get the total number of species across all phases
721 state_size = 0
722 do i_phase = 1, size(this%aero_phase)
723 state_size = state_size + this%aero_phase(i_phase)%val%size()
724 end do
725
726 end function get_size
727
728!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
729
730 !> Get a list of unique names for each element on the
731 !! \c camp_camp_state::camp_state_t::state_var array for this aerosol
732 !! representation. The list may be restricted to a particular phase and/or
733 !! aerosol species by including the phase_name and spec_name arguments.
734 !!
735 !! For a modal/binned mass representation, the unique names for bins are:
736 !! - "bin name.bin #.phase name.species name"
737 !!
738 !! ... and for modes are:
739 !! - "mode name.phase name.species name"
740 function unique_names(this, phase_name, tracer_type, spec_name, &
741 phase_is_at_surface)
742
743 !> List of unique names
744 type(string_t), allocatable :: unique_names(:)
745 !> Aerosol representation data
746 class(aero_rep_modal_binned_mass_t), intent(in) :: this
747 !> Aerosol phase name
748 character(len=*), optional, intent(in) :: phase_name
749 !> Aerosol-phase species tracer type
750 integer(kind=i_kind), optional, intent(in) :: tracer_type
751 !> Aerosol-phase species name
752 character(len=*), optional, intent(in) :: spec_name
753 !> Aerosol-phase species is at surface
754 logical, optional, intent(in) :: phase_is_at_surface
755
756 integer(kind=i_kind) :: num_spec, i_spec, j_spec, i_phase, j_phase, &
757 i_section, i_bin
758 integer(kind=i_kind) :: curr_tracer_type
759 character(len=:), allocatable :: curr_section_name, curr_phase_name, &
760 curr_bin_str
761 type(string_t), allocatable :: spec_names(:)
762
763 ! Count the number of unique names
764 num_spec = 0
765 do i_phase = 1, size(this%aero_phase)
766
767 ! Filter by phase name
768 if (present(phase_name)) then
769 curr_phase_name = this%aero_phase(i_phase)%val%name()
770 if (phase_name.ne.curr_phase_name) cycle
771 end if
772
773 ! Filter by phase is at surface
774 if (present(phase_is_at_surface)) then
775 if (phase_is_at_surface .neqv. &
776 this%aero_phase_is_at_surface(i_phase)) cycle
777 end if
778
779 ! Filter by spec name and/or tracer type
780 if (present(spec_name).or.present(tracer_type)) then
781 spec_names = this%aero_phase(i_phase)%val%get_species_names()
782 do j_spec = 1, size(spec_names)
783 curr_tracer_type = &
784 this%aero_phase(i_phase)%val%get_species_type( &
785 spec_names(j_spec)%string)
786 if (present(spec_name)) then
787 if (spec_name.ne.spec_names(j_spec)%string) cycle
788 end if
789 if (present(tracer_type)) then
790 if (tracer_type.ne.curr_tracer_type) cycle
791 end if
792 num_spec = num_spec + 1
793 end do
794 else
795 num_spec = num_spec + this%aero_phase(i_phase)%val%size()
796 end if
797
798 end do
799
800 ! Allocate space for the unique names
801 allocate(unique_names(num_spec))
802
803 ! Loop through the modes/bin sets
804 i_phase = 1
805 i_spec = 1
806 do i_section = 1, num_section_
807
808 ! Get the current section name
809 curr_section_name = this%section_name(i_section)%string
810
811 ! Loop through the phases for this mode/bin set
812 do j_phase = 1, num_phase_(i_section)
813
814 ! Set the current phase name
815 curr_phase_name = this%aero_phase(i_phase)%val%name()
816
817 ! Filter by phase name
818 if (present(phase_name)) then
819 if (phase_name.ne.curr_phase_name) then
820 i_phase = i_phase + num_bins_(i_section)
821 cycle
822 end if
823 end if
824
825 ! Filter by phase is at surface
826 if (present(phase_is_at_surface)) then
827 if (phase_is_at_surface .neqv. &
828 this%aero_phase_is_at_surface(i_phase)) then
829 i_phase = i_phase + num_bins_(i_section)
830 cycle
831 end if
832 end if
833
834 ! Get the species names in this phase
835 spec_names = this%aero_phase(i_phase)%val%get_species_names()
836
837 ! Loop through the bins (one iteration for modes)
838 do i_bin = 1, num_bins_(i_section)
839
840 ! Set the current bin label (except for single bins or mode)
841 if (num_bins_(i_section).gt.1) then
842 curr_bin_str = trim(to_string(i_bin))//"."
843 else
844 curr_bin_str = ""
845 end if
846
847 ! Add species from this phase/bin
848 num_spec = this%aero_phase(i_phase)%val%size()
849 do j_spec = 1, num_spec
850
851 ! Filter by species name
852 if (present(spec_name)) then
853 if (spec_name.ne.spec_names(j_spec)%string) cycle
854 end if
855
856 ! Filter by species tracer type
857 if (present(tracer_type)) then
858 curr_tracer_type = &
859 this%aero_phase(i_phase)%val%get_species_type( &
860 spec_names(j_spec)%string)
861 if (tracer_type.ne.curr_tracer_type) cycle
862 end if
863
864 ! Add the unique name for this species
865 unique_names(i_spec)%string = curr_section_name//"."// &
866 curr_bin_str//curr_phase_name//'.'// &
867 spec_names(j_spec)%string
868
869 i_spec = i_spec + 1
870 end do
871
872 ! Move to the next phase instance
873 i_phase = i_phase + 1
874
875 end do
876
877 deallocate(spec_names)
878
879 end do
880 end do
881
882 end function unique_names
883
884!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
885
886 !> Get a species id on the \c camp_camp_state::camp_state_t::state_var
887 !! array by its unique name. These are unique ids for each element on the
888 !! state array for this \ref camp_aero_rep "aerosol representation" and
889 !! are numbered:
890 !!
891 !! \f[x_u \in x_f ... (x_f+n-1)\f]
892 !!
893 !! where \f$x_u\f$ is the id of the element corresponding to the species
894 !! with unique name \f$u\f$ on the \c
895 !! camp_camp_state::camp_state_t::state_var array, \f$x_f\f$ is the index
896 !! of the first element for this aerosol representation on the state array
897 !! and \f$n\f$ is the total number of variables on the state array from
898 !! this aerosol representation.
899 function spec_state_id(this, unique_name) result (spec_id)
900
901 !> Species state id
902 integer(kind=i_kind) :: spec_id
903 !> Aerosol representation data
904 class(aero_rep_modal_binned_mass_t), intent(in) :: this
905 !> Unique name
906 character(len=*), intent(in) :: unique_name
907
908 type(string_t), allocatable :: unique_names(:)
909 integer(kind=i_kind) :: i_spec
910
911 spec_id = 0
912 unique_names = this%unique_names()
913 do i_spec = 1, size(unique_names)
914 if (unique_names(i_spec)%string .eq. unique_name) then
915 spec_id = this%phase_state_id(1) + i_spec - 1
916 return
917 end if
918 end do
919 call die_msg( 105414960, "Cannot find species '"//unique_name//"'" )
920
921 end function spec_state_id
922
923!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
924
925 !> Get the non-unique name of a species by its unique name
926 function spec_name(this, unique_name)
927
928 !> Chemical species name
929 character(len=:), allocatable :: spec_name
930 !> Aerosol representation data
931 class(aero_rep_modal_binned_mass_t), intent(in) :: this
932 !> Unique name of the species in this aerosol representation
933 character(len=*), intent(in) :: unique_name
934
935 ! Indices for iterators
936 integer(kind=i_kind) :: i_spec, j_spec, i_phase
937
938 ! species names in the aerosol phase
939 type(string_t), allocatable :: spec_names(:)
940
941 ! unique name list
942 type(string_t), allocatable :: unique_names(:)
943
944 unique_names = this%unique_names()
945
946 i_spec = 1
947 do i_phase = 1, size(this%aero_phase)
948 spec_names = this%aero_phase(i_phase)%val%get_species_names()
949 do j_spec = 1, this%aero_phase(i_phase)%val%size()
950 if (unique_name.eq.unique_names(i_spec)%string) then
951 spec_name = spec_names(j_spec)%string
952 end if
953 i_spec = i_spec + 1
954 end do
955 deallocate(spec_names)
956 end do
957
958 deallocate(unique_names)
959
960 end function spec_name
961
962!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
963
964 !> Get the number of instances of a specified aerosol phase.
965 function num_phase_instances(this, phase_name, is_at_surface)
966
967 !> Number of instances of the aerosol phase
968 integer(kind=i_kind) :: num_phase_instances
969 !> Aerosol representation data
970 class(aero_rep_modal_binned_mass_t), intent(in) :: this
971 !> Aerosol phase name
972 character(len=*), intent(in) :: phase_name
973 !> Indicates if aerosol phase is at the surface of particle
974 logical, intent(in), optional :: is_at_surface
975
976 integer(kind=i_kind) :: i_phase
977
979 if (present(is_at_surface)) then
980 if (is_at_surface) then
981 do i_phase = 1, size(this%aero_phase)
982 if (this%aero_phase(i_phase)%val%name().eq.phase_name .and. &
983 this%aero_phase_is_at_surface(i_phase)) then
985 end if
986 end do
987 else
988 do i_phase = 1, size(this%aero_phase)
989 if (this%aero_phase(i_phase)%val%name().eq.phase_name .and. &
990 .not. this%aero_phase_is_at_surface(i_phase)) then
991 call die_msg(507144607, "Species must exist at surface "// &
992 "in modal/binned mass aerosol representation '"// &
993 this%rep_name//"'")
994 end if
995 end do
996 end if
997 else
998 do i_phase = 1, size(this%aero_phase)
999 if (this%aero_phase(i_phase)%val%name().eq.phase_name) then
1001 end if
1002 end do
1003 end if
1004
1005 end function num_phase_instances
1006
1007!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1008
1009 !> Get the number of Jacobian elements used in calculations of aerosol mass,
1010 !! volume, number, etc. for a particular phase
1011 function num_jac_elem(this, phase_id)
1012
1013 !> Number of Jacobian elements used
1014 integer(kind=i_kind) :: num_jac_elem
1015 !> Aerosol respresentation data
1016 class(aero_rep_modal_binned_mass_t), intent(in) :: this
1017 !> Aerosol phase id
1018 integer(kind=i_kind), intent(in) :: phase_id
1019
1020 integer(kind=i_kind) :: i_section, i_phase, i_bin, j_phase, &
1021 phase_id_to_find, section_size, &
1022 section_start
1023
1024 phase_id_to_find = phase_id
1025 section_start = 1
1026 do i_section = 1, num_section_
1027 section_size = num_phase_(i_section) * num_bins_(i_section)
1028 if( phase_id_to_find .le. section_size ) then
1029 i_phase = ( phase_id_to_find - 1 ) / num_bins_(i_section) + 1
1030 i_bin = mod( phase_id_to_find - 1, num_bins_(i_section) ) + 1
1031 num_jac_elem = 0
1032 do j_phase = section_start + i_bin - 1, &
1033 section_start + i_bin - 1 + &
1034 ( num_phase_(i_section) - 1 ) * num_bins_(i_section), &
1035 num_bins_(i_section)
1037 this%aero_phase( j_phase )%val%num_jac_elem( )
1038 end do
1039 return
1040 end if
1041 phase_id_to_find = phase_id_to_find - section_size
1042 section_start = section_start + section_size
1043 end do
1044
1045 end function num_jac_elem
1046
1047!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1048
1049 !> Determine is specified phase(s) exist in adjacent layers. Returns array
1050 !! of phase_ids for adjacent phases first and second.
1051
1052 function adjacent_phases(this, phase_name_first, &
1053 phase_name_second) result(index_pairs)
1054 !> Aerosol representation data
1055 class(aero_rep_modal_binned_mass_t), intent(in) :: this
1056 character(len=*), intent(in) :: phase_name_first
1057 character(len=*), intent(in) :: phase_name_second
1058 type(index_pair_t), allocatable :: index_pairs(:)
1059
1060 allocate(index_pairs(0))
1061
1062 end function adjacent_phases
1063
1064!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1065
1066 !> Finalize the aerosol representation
1067 subroutine finalize(this)
1068
1069 !> Aerosol representation data
1070 type(aero_rep_modal_binned_mass_t), intent(inout) :: this
1071
1072 if (allocated(this%rep_name)) deallocate(this%rep_name)
1073 if (allocated(this%aero_phase)) then
1074 ! The core will deallocate the aerosol phases
1075 call this%aero_phase(:)%dereference()
1076 deallocate(this%aero_phase)
1077 end if
1078 if (associated(this%property_set)) &
1079 deallocate(this%property_set)
1080 if (allocated(this%condensed_data_real)) &
1081 deallocate(this%condensed_data_real)
1082 if (allocated(this%condensed_data_int)) &
1083 deallocate(this%condensed_data_int)
1084
1085 end subroutine finalize
1086
1087!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1088
1089 !> Finalize the aerosol representation array
1090 subroutine finalize_array(aero_reps)
1091
1092 !> Aerosol representation array
1093 type(aero_rep_modal_binned_mass_t), intent(inout) :: aero_reps(:)
1094
1095 integer(kind=i_kind) :: i_rep
1096
1097 do i_rep = 1, size(aero_reps)
1098 call finalize(aero_reps(i_rep))
1099 end do
1100
1101 end subroutine finalize_array
1102
1103!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1104
1105 !> Initialize a GMD update object
1106 subroutine update_data_init_gmd(this, update_data, aero_rep_type)
1107
1108 use camp_rand, only : generate_int_id
1109
1110 !> Aerosol representation to update
1111 class(aero_rep_modal_binned_mass_t), intent(inout) :: this
1112 !> Update data object
1113 class(aero_rep_update_data_modal_binned_mass_gmd_t), intent(out) :: &
1114 update_data
1115 !> Aerosol representation id
1116 integer(kind=i_kind), intent(in) :: aero_rep_type
1117
1118 ! If an aerosol representation id has not been generated, do it now
1119 if (aero_rep_id_.eq.-1) then
1120 aero_rep_id_ = generate_int_id()
1121 end if
1122
1123 update_data%aero_rep_unique_id = aero_rep_id_
1124 update_data%aero_rep_type = int(aero_rep_type, kind=c_int)
1125 update_data%update_data = aero_rep_modal_binned_mass_create_gmd_update_data()
1126 update_data%is_malloced = .true.
1127
1128 end subroutine update_data_init_gmd
1129
1130!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1131
1132 !> Set packed update data for mode GMD
1133 subroutine update_data_set_gmd(this, section_id, GMD)
1134
1135 !> Update data
1136 class(aero_rep_update_data_modal_binned_mass_gmd_t), intent(inout) :: this
1137 !> Aerosol section id from
1138 !! camp_aero_rep_modal_binned_mass::aero_rep_modal_binned_mass_t::get_section_id
1139 integer(kind=i_kind), intent(in) :: section_id
1140 !> Updated GMD (m)
1141 real(kind=dp), intent(in) :: gmd
1142
1144 this%aero_rep_unique_id, section_id-1, gmd)
1145
1146 end subroutine update_data_set_gmd
1147
1148!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1149
1150 !> Determine the size of a binary required to pack the reaction data
1151 integer(kind=i_kind) function internal_pack_size_gmd(this, comm) &
1152 result(pack_size)
1153
1154 !> Aerosol representation update data
1155 class(aero_rep_update_data_modal_binned_mass_gmd_t), intent(in) :: this
1156 !> MPI communicator
1157 integer, intent(in) :: comm
1158
1159 pack_size = &
1160 camp_mpi_pack_size_logical(this%is_malloced, comm) + &
1161 camp_mpi_pack_size_integer(this%aero_rep_unique_id, comm)
1162
1163 end function internal_pack_size_gmd
1164
1165!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1166
1167 !> Pack the given value to the buffer, advancing position
1168 subroutine internal_bin_pack_gmd(this, buffer, pos, comm)
1169
1170 !> Aerosol representation update data
1171 class(aero_rep_update_data_modal_binned_mass_gmd_t), intent(in) :: this
1172 !> Memory buffer
1173 character, intent(inout) :: buffer(:)
1174 !> Current buffer position
1175 integer, intent(inout) :: pos
1176 !> MPI communicator
1177 integer, intent(in) :: comm
1178
1179#ifdef CAMP_USE_MPI
1180 integer :: prev_position
1181
1182 prev_position = pos
1183 call camp_mpi_pack_logical(buffer, pos, this%is_malloced, comm)
1184 call camp_mpi_pack_integer(buffer, pos, this%aero_rep_unique_id, comm)
1185 call assert(685522546, &
1186 pos - prev_position <= this%pack_size(comm))
1187#endif
1188
1189 end subroutine internal_bin_pack_gmd
1190
1191!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1192
1193 !> Unpack the given value from the buffer, advancing position
1194 subroutine internal_bin_unpack_gmd(this, buffer, pos, comm)
1195
1196 !> Aerosol representation update data
1197 class(aero_rep_update_data_modal_binned_mass_gmd_t), intent(inout) :: this
1198 !> Memory buffer
1199 character, intent(inout) :: buffer(:)
1200 !> Current buffer position
1201 integer, intent(inout) :: pos
1202 !> MPI communicator
1203 integer, intent(in) :: comm
1204
1205#ifdef CAMP_USE_MPI
1206 integer :: prev_position
1207
1208 prev_position = pos
1209 call camp_mpi_unpack_logical(buffer, pos, this%is_malloced, comm)
1210 call camp_mpi_unpack_integer(buffer, pos, this%aero_rep_unique_id, comm)
1211 call assert(855679450, &
1212 pos - prev_position <= this%pack_size(comm))
1214#endif
1215
1216 end subroutine internal_bin_unpack_gmd
1217
1218!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1219
1220 !> Finalize a GMD update data object
1222
1223 !> Update data object to free
1224 type(aero_rep_update_data_modal_binned_mass_gmd_t), intent(inout) :: this
1225
1226 if (this%is_malloced) call aero_rep_free_update_data(this%update_data)
1227
1228 end subroutine update_data_gmd_finalize
1229
1230!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1231
1232 !> Finalize a GMD update data array
1233 subroutine update_data_gmd_finalize_array(update_data)
1234
1235 !> Update data array to free
1236 type(aero_rep_update_data_modal_binned_mass_gmd_t), intent(inout) :: &
1237 update_data(:)
1238
1239 integer(kind=i_kind) :: i_data
1240
1241 do i_data = 1, size(update_data)
1242 call update_data_gmd_finalize(update_data(i_data))
1243 end do
1244
1245 end subroutine update_data_gmd_finalize_array
1246
1247!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1248
1249 !> Initialize a GSD update data object
1250 subroutine update_data_init_gsd(this, update_data, aero_rep_type)
1251
1252 use camp_rand, only : generate_int_id
1253
1254 !> Aerosol representation to update
1255 class(aero_rep_modal_binned_mass_t), intent(inout) :: this
1256 !> Update data object
1257 class(aero_rep_update_data_modal_binned_mass_gsd_t), intent(out) :: &
1258 update_data
1259 !> Aerosol representation id
1260 integer(kind=i_kind), intent(in) :: aero_rep_type
1261
1262 ! If an aerosol representation id has not been generated, do it now
1263 if (aero_rep_id_.eq.-1) then
1264 aero_rep_id_ = generate_int_id()
1265 end if
1266
1267 update_data%aero_rep_unique_id = aero_rep_id_
1268 update_data%aero_rep_type = int(aero_rep_type, kind=c_int)
1269 update_data%update_data = aero_rep_modal_binned_mass_create_gsd_update_data()
1270 update_data%is_malloced = .true.
1271
1272 end subroutine update_data_init_gsd
1273
1274!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1275
1276 !> Set packed update data for mode GSD
1277 subroutine update_data_set_gsd(this, section_id, GSD)
1278
1279 !> Update data
1280 class(aero_rep_update_data_modal_binned_mass_gsd_t), intent(inout) :: this
1281 !> Aerosol section id from
1282 !! camp_aero_rep_modal_binned_mass::aero_rep_modal_binned_mass_t::get_section_id
1283 integer(kind=i_kind), intent(in) :: section_id
1284 !> Updated GSD (m)
1285 real(kind=dp), intent(in) :: gsd
1286
1288 this%aero_rep_unique_id, section_id-1, gsd)
1289
1290 end subroutine update_data_set_gsd
1291
1292!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1293
1294 !> Determine the size of a binary required to pack the reaction data
1295 integer(kind=i_kind) function internal_pack_size_gsd(this, comm) &
1296 result(pack_size)
1297
1298 !> Aerosol representation update data
1299 class(aero_rep_update_data_modal_binned_mass_gsd_t), intent(in) :: this
1300 !> MPI communicator
1301 integer, intent(in) :: comm
1302
1303 pack_size = &
1304 camp_mpi_pack_size_logical(this%is_malloced, comm) + &
1305 camp_mpi_pack_size_integer(this%aero_rep_unique_id, comm)
1306
1307 end function internal_pack_size_gsd
1308
1309!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1310
1311 !> Pack the given value to the buffer, advancing position
1312 subroutine internal_bin_pack_gsd(this, buffer, pos, comm)
1313
1314 !> Aerosol representation update data
1315 class(aero_rep_update_data_modal_binned_mass_gsd_t), intent(in) :: this
1316 !> Memory buffer
1317 character, intent(inout) :: buffer(:)
1318 !> Current buffer position
1319 integer, intent(inout) :: pos
1320 !> MPI communicator
1321 integer, intent(in) :: comm
1322
1323#ifdef CAMP_USE_MPI
1324 integer :: prev_position
1325
1326 prev_position = pos
1327 call camp_mpi_pack_logical(buffer, pos, this%is_malloced, comm)
1328 call camp_mpi_pack_integer(buffer, pos, this%aero_rep_unique_id, comm)
1329 call assert(295993259, &
1330 pos - prev_position <= this%pack_size(comm))
1331#endif
1332
1333 end subroutine internal_bin_pack_gsd
1334
1335!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1336
1337 !> Unpack the given value from the buffer, advancing position
1338 subroutine internal_bin_unpack_gsd(this, buffer, pos, comm)
1339
1340 !> Aerosol representation update data
1341 class(aero_rep_update_data_modal_binned_mass_gsd_t), intent(inout) :: this
1342 !> Memory buffer
1343 character, intent(inout) :: buffer(:)
1344 !> Current buffer position
1345 integer, intent(inout) :: pos
1346 !> MPI communicator
1347 integer, intent(in) :: comm
1348
1349#ifdef CAMP_USE_MPI
1350 integer :: prev_position
1351
1352 prev_position = pos
1353 call camp_mpi_unpack_logical(buffer, pos, this%is_malloced, comm)
1354 call camp_mpi_unpack_integer(buffer, pos, this%aero_rep_unique_id, comm)
1355 call assert(518724415, &
1356 pos - prev_position <= this%pack_size(comm))
1358#endif
1359
1360 end subroutine internal_bin_unpack_gsd
1361
1362!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1363
1364 !> Finalize a GSD update data object
1366
1367 !> Update data object to free
1368 type(aero_rep_update_data_modal_binned_mass_gsd_t), intent(inout) :: this
1369
1370 if (this%is_malloced) call aero_rep_free_update_data(this%update_data)
1371
1372 end subroutine update_data_gsd_finalize
1373
1374!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1375
1376 !> Finalize a GSD update data array
1377 subroutine update_data_gsd_finalize_array(update_data)
1378
1379 !> Update data array to free
1380 type(aero_rep_update_data_modal_binned_mass_gsd_t), intent(inout) :: &
1381 update_data(:)
1382
1383 integer(kind=i_kind) :: i_data
1384
1385 do i_data = 1, size(update_data)
1386 call update_data_gsd_finalize(update_data(i_data))
1387 end do
1388
1389 end subroutine update_data_gsd_finalize_array
1390
1391!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1392
Get the size of the section of the camp_camp_state::camp_state_t::state_var array required for this a...
Initialize the aerosol representation data, validating component data and loading any required inform...
Extending-type binary pack function (Internal use only)
Extending-type binary unpack function (Internal use only)
Extending-type binary pack size (internal use only)
Get the number of Jacobian elements used in calculations of aerosol mass, volume, number,...
Get the number of instances of a specified aerosol phase.
Get the non-unique name of a chemical species by its unique name.
Get a species id on the camp_camp_state::camp_state_t::state_var array by unique name....
Get a list of unique names for each element on the camp_camp_state::camp_state_t::state_var array for...
Interface for to_string functions.
Definition util.F90:32
The abstract aero_phase_data_t structure and associated subroutines.
subroutine finalize_array(this)
Finalize the aerosol phase data.
type(aero_phase_data_t) function, pointer constructor(phase_name, init_size)
Constructor for aero_phase_data_t.
subroutine finalize(this)
Finalize the aerosol phase data.
integer(kind=i_kind) function pack_size(this, comm)
Determine the size of a binary required to pack the aerosol representation data.
The abstract aero_rep_data_t structure and associated subroutines.
The abstract aero_rep_modal_binned_mass_t structure and associated subroutines.
subroutine update_data_init_gsd(this, update_data, aero_rep_type)
Initialize a GSD update data object.
subroutine internal_bin_unpack_gsd(this, buffer, pos, comm)
Unpack the given value from the buffer, advancing position.
logical function get_section_id(this, section_name, section_id)
Get an id for a mode or bin by name for use with updates from external modules.
subroutine internal_bin_pack_gmd(this, buffer, pos, comm)
Pack the given value to the buffer, advancing position.
subroutine update_data_gmd_finalize_array(update_data)
Finalize a GMD update data array.
subroutine update_data_set_gsd(this, section_id, gsd)
Set packed update data for mode GSD.
type(index_pair_t) function, dimension(:), allocatable adjacent_phases(this, phase_name_first, phase_name_second)
Determine is specified phase(s) exist in adjacent layers. Returns array of phase_ids for adjacent pha...
subroutine internal_bin_unpack_gmd(this, buffer, pos, comm)
Unpack the given value from the buffer, advancing position.
subroutine internal_bin_pack_gsd(this, buffer, pos, comm)
Pack the given value to the buffer, advancing position.
integer(kind=i_kind) function internal_pack_size_gmd(this, comm)
Determine the size of a binary required to pack the reaction data.
subroutine update_data_set_gmd(this, section_id, gmd)
Set packed update data for mode GMD.
integer(kind=i_kind) function internal_pack_size_gsd(this, comm)
Determine the size of a binary required to pack the reaction data.
integer(kind=i_kind), parameter, public update_gsd
subroutine update_data_gsd_finalize(this)
Finalize a GSD update data object.
subroutine update_data_gsd_finalize_array(update_data)
Finalize a GSD update data array.
subroutine update_data_gmd_finalize(this)
Finalize a GMD update data object.
integer(kind=i_kind), parameter, public update_gmd
subroutine update_data_init_gmd(this, update_data, aero_rep_type)
Initialize a GMD update object.
The camp_state_t structure and associated subroutines.
Definition camp_state.F90:9
The chem_spec_data_t structure and associated subroutines.
Wrapper functions for MPI.
Definition mpi.F90:13
subroutine camp_mpi_pack_logical(buffer, position, val, comm)
Packs the given value into the buffer, advancing position.
Definition mpi.F90:792
subroutine camp_mpi_unpack_integer(buffer, position, val, comm)
Unpacks the given value from the buffer, advancing position.
Definition mpi.F90:1023
subroutine camp_mpi_unpack_logical(buffer, position, val, comm)
Unpacks the given value from the buffer, advancing position.
Definition mpi.F90:1131
integer function camp_mpi_pack_size_logical(val, comm)
Determines the number of bytes required to pack the given value.
Definition mpi.F90:484
subroutine camp_mpi_pack_integer(buffer, position, val, comm)
Packs the given value into the buffer, advancing position.
Definition mpi.F90:691
integer function camp_mpi_pack_size_integer(val, comm)
Determines the number of bytes required to pack the given value.
Definition mpi.F90:398
The property_t structure and associated subroutines.
Definition property.F90:9
Random number generators.
Definition rand.F90:9
integer(kind=i_kind) function generate_int_id()
Generate an integer id Ids will be sequential, and can only be generated by the primary process.
Definition rand.F90:435
Common utility subroutines.
Definition util.F90:9
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition util.F90:165
subroutine die_msg(code, error_msg)
Error immediately.
Definition util.F90:196
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
Definition util.F90:130
Pointer type for building arrays.
Abstract aerosol representation data type.
Define index_pair array for adjacent_phases functions.
String type for building arrays of string of various size.
Definition util.F90:53