PartMC
2.2.1
|
00001 ! Copyright (C) 2005-2011 Nicole Riemer and Matthew West 00002 ! Licensed under the GNU General Public License version 2 or (at your 00003 ! option) any later version. See the file COPYING for details. 00004 00005 !> \file 00006 !> The pmc_aero_particle module. 00007 00008 !> The aero_particle_t structure and associated subroutines. 00009 module pmc_aero_particle 00010 00011 use pmc_util 00012 use pmc_aero_data 00013 use pmc_spec_file 00014 use pmc_mpi 00015 #ifdef PMC_USE_MPI 00016 use mpi 00017 #endif 00018 00019 !> Single aerosol particle data structure. 00020 !! 00021 !! The \c vol array stores the total volumes of the different 00022 !! species that make up the particle. This array must have length 00023 !! equal to aero_data%%n_spec, so that \c vol(i) is the volume (in 00024 !! m^3) of the i'th aerosol species. 00025 type aero_particle_t 00026 !> Constituent species volumes [length aero_data%%n_spec] (m^3). 00027 real(kind=dp), pointer :: vol(:) 00028 !> Number of original particles from each source that coagulated 00029 !> to form this one [length aero_data%%n_source]. 00030 integer, pointer :: n_orig_part(:) 00031 !> Weighting function group number. 00032 integer :: weight_group 00033 !> Absorption cross-section (m^2). 00034 real(kind=dp) :: absorb_cross_sect 00035 !> Scattering cross-section (m^2). 00036 real(kind=dp) :: scatter_cross_sect 00037 !> Asymmetry parameter (1). 00038 real(kind=dp) :: asymmetry 00039 !> Refractive index of the shell (1). 00040 complex(kind=dc) :: refract_shell 00041 !> Refractive index of the core (1). 00042 complex(kind=dc) :: refract_core 00043 !> Volume of the core (m^3). 00044 real(kind=dp) :: core_vol 00045 !> Water hysteresis curve section (0 = lower, 1 = upper) 00046 integer :: water_hyst_leg 00047 !> Unique ID number. 00048 integer :: id 00049 !> First time a constituent was created (s). 00050 real(kind=dp) :: least_create_time 00051 !> Last time a constituent was created (s). 00052 real(kind=dp) :: greatest_create_time 00053 end type aero_particle_t 00054 00055 !> Next unique ID number to use for a particle. 00056 integer, save :: next_id = 1 00057 00058 contains 00059 00060 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00061 00062 !> Allocates memory in an aero_particle_t. 00063 subroutine aero_particle_allocate(aero_particle) 00064 00065 !> Particle to init. 00066 type(aero_particle_t), intent(out) :: aero_particle 00067 00068 allocate(aero_particle%vol(0)) 00069 allocate(aero_particle%n_orig_part(0)) 00070 call aero_particle_zero(aero_particle) 00071 00072 end subroutine aero_particle_allocate 00073 00074 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00075 00076 !> Allocates an aero_particle_t of the given size. 00077 subroutine aero_particle_allocate_size(aero_particle, n_spec, n_source) 00078 00079 !> Particle to init. 00080 type(aero_particle_t), intent(out) :: aero_particle 00081 !> Number of species. 00082 integer, intent(in) :: n_spec 00083 !> Number of sources. 00084 integer, intent(in) :: n_source 00085 00086 allocate(aero_particle%vol(n_spec)) 00087 allocate(aero_particle%n_orig_part(n_source)) 00088 call aero_particle_zero(aero_particle) 00089 00090 end subroutine aero_particle_allocate_size 00091 00092 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00093 00094 !> Deallocates memory associated with an aero_particle_t. 00095 subroutine aero_particle_deallocate(aero_particle) 00096 00097 !> Particle to free. 00098 type(aero_particle_t), intent(inout) :: aero_particle 00099 00100 deallocate(aero_particle%vol) 00101 deallocate(aero_particle%n_orig_part) 00102 00103 end subroutine aero_particle_deallocate 00104 00105 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00106 00107 !> Copies a particle. 00108 subroutine aero_particle_copy(aero_particle_from, aero_particle_to) 00109 00110 !> Reference particle. 00111 type(aero_particle_t), intent(in) :: aero_particle_from 00112 !> Destination particle (already alloced on entry). 00113 type(aero_particle_t), intent(inout) :: aero_particle_to 00114 00115 integer :: n_spec, n_source 00116 00117 n_spec = size(aero_particle_from%vol) 00118 n_source = size(aero_particle_from%n_orig_part) 00119 if ((n_spec /= size(aero_particle_to%vol)) & 00120 .or. (n_source /= size(aero_particle_to%n_orig_part))) then 00121 call aero_particle_deallocate(aero_particle_to) 00122 call aero_particle_allocate_size(aero_particle_to, n_spec, n_source) 00123 end if 00124 call assert(651178226, size(aero_particle_from%vol) & 00125 == size(aero_particle_to%vol)) 00126 call assert(105980551, size(aero_particle_from%n_orig_part) & 00127 == size(aero_particle_to%n_orig_part)) 00128 aero_particle_to%vol = aero_particle_from%vol 00129 aero_particle_to%n_orig_part = aero_particle_from%n_orig_part 00130 aero_particle_to%weight_group = aero_particle_from%weight_group 00131 aero_particle_to%absorb_cross_sect = aero_particle_from%absorb_cross_sect 00132 aero_particle_to%scatter_cross_sect = & 00133 aero_particle_from%scatter_cross_sect 00134 aero_particle_to%asymmetry = aero_particle_from%asymmetry 00135 aero_particle_to%refract_shell = aero_particle_from%refract_shell 00136 aero_particle_to%refract_core = aero_particle_from%refract_core 00137 aero_particle_to%core_vol = aero_particle_from%core_vol 00138 aero_particle_to%water_hyst_leg = aero_particle_from%water_hyst_leg 00139 aero_particle_to%id = aero_particle_from%id 00140 aero_particle_to%least_create_time = aero_particle_from%least_create_time 00141 aero_particle_to%greatest_create_time = & 00142 aero_particle_from%greatest_create_time 00143 00144 end subroutine aero_particle_copy 00145 00146 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00147 00148 !> Shift data from one aero_particle_t to another and free the first 00149 !> one. 00150 !! 00151 !! This is roughly equivalent to aero_particle_copy(from, to) 00152 !! followed by aero_particle_deallocate(from), but faster and with 00153 !! different memory allocation requirements. 00154 subroutine aero_particle_shift(aero_particle_from, aero_particle_to) 00155 00156 !> Reference particle (will be deallocated on return). 00157 type(aero_particle_t), intent(inout) :: aero_particle_from 00158 !> Destination particle (not allocated on entry). 00159 type(aero_particle_t), intent(inout) :: aero_particle_to 00160 00161 aero_particle_to%vol => aero_particle_from%vol 00162 nullify(aero_particle_from%vol) 00163 aero_particle_to%n_orig_part => aero_particle_from%n_orig_part 00164 nullify(aero_particle_from%n_orig_part) 00165 aero_particle_to%weight_group = aero_particle_from%weight_group 00166 aero_particle_to%absorb_cross_sect = aero_particle_from%absorb_cross_sect 00167 aero_particle_to%scatter_cross_sect = & 00168 aero_particle_from%scatter_cross_sect 00169 aero_particle_to%asymmetry = aero_particle_from%asymmetry 00170 aero_particle_to%refract_shell = aero_particle_from%refract_shell 00171 aero_particle_to%refract_core = aero_particle_from%refract_core 00172 aero_particle_to%core_vol = aero_particle_from%core_vol 00173 aero_particle_to%water_hyst_leg = aero_particle_from%water_hyst_leg 00174 aero_particle_to%id = aero_particle_from%id 00175 aero_particle_to%least_create_time = aero_particle_from%least_create_time 00176 aero_particle_to%greatest_create_time = & 00177 aero_particle_from%greatest_create_time 00178 00179 end subroutine aero_particle_shift 00180 00181 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00182 00183 !> Resets an aero_particle to be zero. 00184 subroutine aero_particle_zero(aero_particle) 00185 00186 !> Particle to zero. 00187 type(aero_particle_t), intent(inout) :: aero_particle 00188 00189 aero_particle%vol = 0d0 00190 aero_particle%n_orig_part = 0 00191 aero_particle%weight_group = 0 00192 aero_particle%absorb_cross_sect = 0d0 00193 aero_particle%scatter_cross_sect = 0d0 00194 aero_particle%asymmetry = 0d0 00195 aero_particle%refract_shell = (0d0, 0d0) 00196 aero_particle%refract_core = (0d0, 0d0) 00197 aero_particle%core_vol = 0d0 00198 aero_particle%water_hyst_leg = 0 00199 aero_particle%id = 0 00200 aero_particle%least_create_time = 0d0 00201 aero_particle%greatest_create_time = 0d0 00202 00203 end subroutine aero_particle_zero 00204 00205 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00206 00207 !> Assigns a globally-unique new ID number to the particle. 00208 subroutine aero_particle_new_id(aero_particle) 00209 00210 !> Particle to set ID for. 00211 type(aero_particle_t), intent(inout) :: aero_particle 00212 00213 aero_particle%id = (next_id - 1) * pmc_mpi_size() + pmc_mpi_rank() + 1 00214 next_id = next_id + 1 00215 00216 end subroutine aero_particle_new_id 00217 00218 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00219 00220 !> Sets the creation times for the particle. 00221 subroutine aero_particle_set_create_time(aero_particle, create_time) 00222 00223 !> Particle to set time for. 00224 type(aero_particle_t), intent(inout) :: aero_particle 00225 !> Creation time. 00226 real(kind=dp), intent(in) :: create_time 00227 00228 aero_particle%least_create_time = create_time 00229 aero_particle%greatest_create_time = create_time 00230 00231 end subroutine aero_particle_set_create_time 00232 00233 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00234 00235 !> Sets the aerosol particle volumes. 00236 subroutine aero_particle_set_vols(aero_particle, vols) 00237 00238 !> Particle. 00239 type(aero_particle_t), intent(inout) :: aero_particle 00240 !> New volumes. 00241 real(kind=dp), intent(in) :: vols(size(aero_particle%vol)) 00242 00243 aero_particle%vol = vols 00244 00245 end subroutine aero_particle_set_vols 00246 00247 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00248 00249 !> Sets the aerosol particle source. 00250 subroutine aero_particle_set_source(aero_particle, i_source) 00251 00252 !> Particle. 00253 type(aero_particle_t), intent(inout) :: aero_particle 00254 !> Source number for the particle. 00255 integer, intent(in) :: i_source 00256 00257 aero_particle%n_orig_part(i_source) = 1 00258 00259 end subroutine aero_particle_set_source 00260 00261 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00262 00263 !> Sets the aerosol particle weight group. 00264 subroutine aero_particle_set_group(aero_particle, i_group) 00265 00266 !> Particle. 00267 type(aero_particle_t), intent(inout) :: aero_particle 00268 !> Weight group number for the particle. 00269 integer, intent(in) :: i_group 00270 00271 aero_particle%weight_group = i_group 00272 00273 end subroutine aero_particle_set_group 00274 00275 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00276 00277 !> Total mass of the particle (kg). 00278 elemental real(kind=dp) function aero_particle_mass(aero_particle, aero_data) 00279 00280 !> Particle. 00281 type(aero_particle_t), intent(in) :: aero_particle 00282 !> Aerosol data. 00283 type(aero_data_t), intent(in) :: aero_data 00284 00285 aero_particle_mass = sum(aero_particle%vol * aero_data%density) 00286 00287 end function aero_particle_mass 00288 00289 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00290 00291 !> Mass of a single species in the particle (kg). 00292 elemental real(kind=dp) function aero_particle_species_mass(aero_particle, & 00293 i_spec, aero_data) 00294 00295 !> Particle. 00296 type(aero_particle_t), intent(in) :: aero_particle 00297 !> Species number to find mass of. 00298 integer, intent(in) :: i_spec 00299 !> Aerosol data. 00300 type(aero_data_t), intent(in) :: aero_data 00301 00302 aero_particle_species_mass = aero_particle%vol(i_spec) & 00303 * aero_data%density(i_spec) 00304 00305 end function aero_particle_species_mass 00306 00307 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00308 00309 !> Total moles in the particle (1). 00310 elemental real(kind=dp) function aero_particle_moles(aero_particle, aero_data) 00311 00312 !> Particle. 00313 type(aero_particle_t), intent(in) :: aero_particle 00314 !> Aerosol data. 00315 type(aero_data_t), intent(in) :: aero_data 00316 00317 aero_particle_moles = sum(aero_particle%vol * aero_data%density & 00318 / aero_data%molec_weight) 00319 00320 end function aero_particle_moles 00321 00322 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00323 00324 !> Total volume of the particle (m^3). 00325 elemental real(kind=dp) function aero_particle_volume(aero_particle) 00326 00327 !> Particle. 00328 type(aero_particle_t), intent(in) :: aero_particle 00329 00330 aero_particle_volume = sum(aero_particle%vol) 00331 00332 end function aero_particle_volume 00333 00334 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00335 00336 !> Total dry volume of the particle (m^3). 00337 elemental real(kind=dp) function aero_particle_dry_volume(aero_particle, & 00338 aero_data) 00339 00340 !> Particle. 00341 type(aero_particle_t), intent(in) :: aero_particle 00342 !> Aerosol data. 00343 type(aero_data_t), intent(in) :: aero_data 00344 00345 integer :: i_spec 00346 00347 aero_particle_dry_volume = 0d0 00348 do i_spec = 1,aero_data%n_spec 00349 if (i_spec /= aero_data%i_water) then 00350 aero_particle_dry_volume = aero_particle_dry_volume & 00351 + aero_particle%vol(i_spec) 00352 end if 00353 end do 00354 00355 end function aero_particle_dry_volume 00356 00357 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00358 00359 !> Total volume (dry or wet) of the particle (m^3). 00360 elemental real(kind=dp) function aero_particle_volume_maybe_dry( & 00361 aero_particle, aero_data, dry_volume) 00362 00363 !> Particle. 00364 type(aero_particle_t), intent(in) :: aero_particle 00365 !> Aerosol data. 00366 type(aero_data_t), intent(in) :: aero_data 00367 !> Whether the desired volume is dry (otherwise wet). 00368 logical, intent(in) :: dry_volume 00369 00370 if (dry_volume) then 00371 aero_particle_volume_maybe_dry & 00372 = aero_particle_dry_volume(aero_particle, aero_data) 00373 else 00374 aero_particle_volume_maybe_dry = aero_particle_volume(aero_particle) 00375 end if 00376 00377 end function aero_particle_volume_maybe_dry 00378 00379 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00380 00381 !> Total radius of the particle (m). 00382 elemental real(kind=dp) function aero_particle_radius(aero_particle) 00383 00384 !> Particle. 00385 type(aero_particle_t), intent(in) :: aero_particle 00386 00387 aero_particle_radius = vol2rad(aero_particle_volume(aero_particle)) 00388 00389 end function aero_particle_radius 00390 00391 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00392 00393 !> Total dry radius of the particle (m). 00394 elemental real(kind=dp) function aero_particle_dry_radius(aero_particle, & 00395 aero_data) 00396 00397 !> Particle. 00398 type(aero_particle_t), intent(in) :: aero_particle 00399 !> Aerosol data. 00400 type(aero_data_t), intent(in) :: aero_data 00401 00402 aero_particle_dry_radius = vol2rad( & 00403 aero_particle_dry_volume(aero_particle, aero_data)) 00404 00405 end function aero_particle_dry_radius 00406 00407 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00408 00409 !> Total diameter of the particle (m). 00410 elemental real(kind=dp) function aero_particle_diameter(aero_particle) 00411 00412 !> Particle. 00413 type(aero_particle_t), intent(in) :: aero_particle 00414 00415 aero_particle_diameter = 2d0 * aero_particle_radius(aero_particle) 00416 00417 end function aero_particle_diameter 00418 00419 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00420 00421 !> Total dry diameter of the particle (m). 00422 elemental real(kind=dp) function aero_particle_dry_diameter(aero_particle, & 00423 aero_data) 00424 00425 !> Particle. 00426 type(aero_particle_t), intent(in) :: aero_particle 00427 !> Aerosol data. 00428 type(aero_data_t), intent(in) :: aero_data 00429 00430 aero_particle_dry_diameter & 00431 = 2d0 * aero_particle_dry_radius(aero_particle, aero_data) 00432 00433 end function aero_particle_dry_diameter 00434 00435 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00436 00437 !> Average density of the particle (kg/m^3). 00438 real(kind=dp) function aero_particle_density(aero_particle, aero_data) 00439 00440 !> Particle. 00441 type(aero_particle_t), intent(in) :: aero_particle 00442 !> Aerosol data. 00443 type(aero_data_t), intent(in) :: aero_data 00444 00445 aero_particle_density = aero_particle_mass(aero_particle, aero_data) & 00446 / aero_particle_volume(aero_particle) 00447 00448 end function aero_particle_density 00449 00450 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00451 00452 !> Returns the volume-average of the non-water elements of quantity. 00453 real(kind=dp) function aero_particle_average_solute_quantity( & 00454 aero_particle, aero_data, quantity) 00455 00456 !> Aerosol particle. 00457 type(aero_particle_t), intent(in) :: aero_particle 00458 !> Aerosol data. 00459 type(aero_data_t), intent(in) :: aero_data 00460 !> Quantity to average. 00461 real(kind=dp), intent(in) :: quantity(:) 00462 00463 real(kind=dp) :: ones(aero_data%n_spec) 00464 00465 ones = 1d0 00466 aero_particle_average_solute_quantity = & 00467 aero_particle_total_solute_quantity(aero_particle, & 00468 aero_data, quantity) & 00469 / aero_particle_total_solute_quantity(aero_particle, aero_data, ones) 00470 00471 end function aero_particle_average_solute_quantity 00472 00473 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00474 00475 !> Returns the volume-total of the non-water elements of quantity. 00476 real(kind=dp) function aero_particle_total_solute_quantity(aero_particle, & 00477 aero_data, quantity) 00478 00479 !> Aerosol particle. 00480 type(aero_particle_t), intent(in) :: aero_particle 00481 !> Aerosol data. 00482 type(aero_data_t), intent(in) :: aero_data 00483 !> Quantity to total. 00484 real(kind=dp), intent(in) :: quantity(:) 00485 00486 real(kind=dp) total 00487 integer i 00488 00489 total = 0d0 00490 do i = 1,aero_data%n_spec 00491 if (i /= aero_data%i_water) then 00492 total = total + aero_particle%vol(i) * quantity(i) 00493 end if 00494 end do 00495 aero_particle_total_solute_quantity = total 00496 00497 end function aero_particle_total_solute_quantity 00498 00499 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00500 00501 !> Returns the water element of quantity. 00502 real(kind=dp) function aero_particle_average_water_quantity(aero_particle, & 00503 aero_data, quantity) 00504 00505 !> Aerosol particle. 00506 type(aero_particle_t), intent(in) :: aero_particle 00507 !> Aerosol data. 00508 type(aero_data_t), intent(in) :: aero_data 00509 !> Quantity to average. 00510 real(kind=dp), intent(in) :: quantity(:) 00511 00512 call assert(420016623, aero_data%i_water > 0) 00513 aero_particle_average_water_quantity = quantity(aero_data%i_water) 00514 00515 end function aero_particle_average_water_quantity 00516 00517 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00518 00519 !> Returns the volume-total of the water element of quantity. 00520 real(kind=dp) function aero_particle_total_water_quantity(aero_particle, & 00521 aero_data, quantity) 00522 00523 !> Aerosol particle. 00524 type(aero_particle_t), intent(in) :: aero_particle 00525 !> Aerosol data. 00526 type(aero_data_t), intent(in) :: aero_data 00527 !> Quantity to total. 00528 real(kind=dp), intent(in) :: quantity(:) 00529 00530 call assert(223343210, aero_data%i_water > 0) 00531 aero_particle_total_water_quantity & 00532 = aero_particle%vol(aero_data%i_water) & 00533 * quantity(aero_data%i_water) 00534 00535 end function aero_particle_total_water_quantity 00536 00537 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00538 00539 !> Returns the water molecular weight. 00540 !> (kg/mole) 00541 real(kind=dp) function aero_particle_water_molec_weight(aero_data) 00542 00543 !> Aerosol data. 00544 type(aero_data_t), intent(in) :: aero_data 00545 00546 call assert(772012490, aero_data%i_water > 0) 00547 aero_particle_water_molec_weight & 00548 = aero_data%molec_weight(aero_data%i_water) 00549 00550 end function aero_particle_water_molec_weight 00551 00552 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00553 00554 !> Returns the average of the solute molecular weight (kg/mole). 00555 real(kind=dp) function aero_particle_solute_molec_weight(aero_particle, & 00556 aero_data) 00557 00558 !> Aerosol particle. 00559 type(aero_particle_t), intent(in) :: aero_particle 00560 !> Aerosol data. 00561 type(aero_data_t), intent(in) :: aero_data 00562 00563 aero_particle_solute_molec_weight & 00564 = aero_particle_average_solute_quantity(aero_particle, & 00565 aero_data, aero_data%molec_weight) 00566 00567 end function aero_particle_solute_molec_weight 00568 00569 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00570 00571 !> Returns the average of the solute ion number (1). 00572 real(kind=dp) function aero_particle_solute_num_ions(aero_particle, & 00573 aero_data) 00574 00575 !> Aerosol data. 00576 type(aero_data_t), intent(in) :: aero_data 00577 !> Aerosol particle. 00578 type(aero_particle_t), intent(in) :: aero_particle 00579 00580 aero_particle_solute_num_ions & 00581 = aero_particle_average_solute_quantity(aero_particle, & 00582 aero_data, real(aero_data%num_ions, kind=dp)) 00583 00584 end function aero_particle_solute_num_ions 00585 00586 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00587 00588 !> Returns the water density (kg/m^3). 00589 real(kind=dp) function aero_particle_water_density(aero_data) 00590 00591 !> Aerosol data. 00592 type(aero_data_t), intent(in) :: aero_data 00593 00594 call assert(235482108, aero_data%i_water > 0) 00595 aero_particle_water_density = aero_data%density(aero_data%i_water) 00596 00597 end function aero_particle_water_density 00598 00599 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00600 00601 !> Returns the average of the solute densities (kg/m^3). 00602 real(kind=dp) function aero_particle_solute_density(aero_particle, & 00603 aero_data) 00604 00605 !> Aerosol data. 00606 type(aero_data_t), intent(in) :: aero_data 00607 !> Aerosol particle. 00608 type(aero_particle_t), intent(in) :: aero_particle 00609 00610 aero_particle_solute_density & 00611 = aero_particle_average_solute_quantity(aero_particle, & 00612 aero_data, aero_data%density) 00613 00614 end function aero_particle_solute_density 00615 00616 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00617 00618 !> Returns the water mass (kg). 00619 real(kind=dp) function aero_particle_water_mass(aero_particle, aero_data) 00620 00621 !> Aerosol data. 00622 type(aero_data_t), intent(in) :: aero_data 00623 !> Aerosol particle. 00624 type(aero_particle_t), intent(in) :: aero_particle 00625 00626 call assert(888636139, aero_data%i_water > 0) 00627 aero_particle_water_mass = aero_particle%vol(aero_data%i_water) & 00628 * aero_data%density(aero_data%i_water) 00629 00630 end function aero_particle_water_mass 00631 00632 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00633 00634 !> Returns the total solute mass (kg). 00635 real(kind=dp) function aero_particle_solute_mass(aero_particle, aero_data) 00636 00637 !> Aerosol data. 00638 type(aero_data_t), intent(in) :: aero_data 00639 !> Aerosol particle. 00640 type(aero_particle_t), intent(in) :: aero_particle 00641 00642 aero_particle_solute_mass & 00643 = aero_particle_total_solute_quantity(aero_particle, & 00644 aero_data, aero_data%density) 00645 00646 end function aero_particle_solute_mass 00647 00648 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00649 00650 !> Returns the total solute volume (m^3). 00651 real(kind=dp) function aero_particle_solute_volume(aero_particle, & 00652 aero_data) 00653 00654 !> Aerosol data. 00655 type(aero_data_t), intent(in) :: aero_data 00656 !> Aerosol particle. 00657 type(aero_particle_t), intent(in) :: aero_particle 00658 00659 real(kind=dp) :: ones(aero_data%n_spec) 00660 00661 ones = 1d0 00662 aero_particle_solute_volume & 00663 = aero_particle_total_solute_quantity(aero_particle, & 00664 aero_data, ones) 00665 00666 end function aero_particle_solute_volume 00667 00668 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00669 00670 !> Returns the total solute radius (m). 00671 real(kind=dp) function aero_particle_solute_radius(aero_particle, & 00672 aero_data) 00673 00674 !> Aerosol data. 00675 type(aero_data_t), intent(in) :: aero_data 00676 !> Aerosol particle. 00677 type(aero_particle_t), intent(in) :: aero_particle 00678 00679 aero_particle_solute_radius & 00680 = vol2rad(aero_particle_solute_volume(aero_particle, aero_data)) 00681 00682 end function aero_particle_solute_radius 00683 00684 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00685 00686 !> Returns the average of the solute kappas (1). 00687 real(kind=dp) function aero_particle_solute_kappa(aero_particle, aero_data) 00688 00689 !> Aerosol data. 00690 type(aero_data_t), intent(in) :: aero_data 00691 !> Aerosol particle. 00692 type(aero_particle_t), intent(in) :: aero_particle 00693 00694 real(kind=dp) :: kappa(aero_data%n_spec), M_w, rho_w, M_a, rho_a 00695 integer :: i_spec 00696 00697 M_w = aero_particle_water_molec_weight(aero_data) 00698 rho_w = aero_particle_water_density(aero_data) 00699 do i_spec = 1,aero_data%n_spec 00700 if (i_spec == aero_data%i_water) then 00701 kappa(i_spec) = 0d0 00702 elseif (aero_data%num_ions(i_spec) > 0) then 00703 call assert_msg(123681459, aero_data%kappa(i_spec) == 0d0, & 00704 "species has nonzero num_ions and kappa: " & 00705 // trim(aero_data%name(i_spec))) 00706 M_a = aero_data%molec_weight(i_spec) 00707 rho_a = aero_data%density(i_spec) 00708 kappa(i_spec) = M_w * rho_a / (M_a * rho_w) & 00709 * real(aero_data%num_ions(i_spec), kind=dp) 00710 else 00711 kappa(i_spec) = aero_data%kappa(i_spec) 00712 end if 00713 end do 00714 aero_particle_solute_kappa & 00715 = aero_particle_average_solute_quantity(aero_particle, & 00716 aero_data, kappa) 00717 00718 end function aero_particle_solute_kappa 00719 00720 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00721 00722 !> Coagulate two particles together to make a new one. The new 00723 !> particle will not have its ID set. 00724 subroutine aero_particle_coagulate(aero_particle_1, & 00725 aero_particle_2, aero_particle_new) 00726 00727 !> First particle. 00728 type(aero_particle_t), intent(in) :: aero_particle_1 00729 !> Second particle. 00730 type(aero_particle_t), intent(in) :: aero_particle_2 00731 !> Result particle. 00732 type(aero_particle_t), intent(inout) :: aero_particle_new 00733 00734 call assert(203741686, size(aero_particle_1%vol) & 00735 == size(aero_particle_new%vol)) 00736 call assert(586181003, size(aero_particle_2%vol) & 00737 == size(aero_particle_new%vol)) 00738 call assert(483452167, size(aero_particle_1%n_orig_part) & 00739 == size(aero_particle_new%n_orig_part)) 00740 call assert(126911035, size(aero_particle_2%n_orig_part) & 00741 == size(aero_particle_new%n_orig_part)) 00742 aero_particle_new%vol = aero_particle_1%vol + aero_particle_2%vol 00743 aero_particle_new%n_orig_part = aero_particle_1%n_orig_part & 00744 + aero_particle_2%n_orig_part 00745 aero_particle_new%weight_group = 0 00746 aero_particle_new%absorb_cross_sect = 0d0 00747 aero_particle_new%scatter_cross_sect = 0d0 00748 aero_particle_new%asymmetry = 0d0 00749 aero_particle_new%refract_shell = (0d0, 0d0) 00750 aero_particle_new%refract_core = (0d0, 0d0) 00751 aero_particle_new%core_vol = 0d0 00752 if ((aero_particle_1%water_hyst_leg == 1) & 00753 .and. (aero_particle_2%water_hyst_leg == 1)) then 00754 aero_particle_new%water_hyst_leg = 1 00755 else 00756 aero_particle_new%water_hyst_leg = 0 00757 end if 00758 aero_particle_new%id = 0 00759 aero_particle_new%least_create_time = & 00760 min(aero_particle_1%least_create_time, & 00761 aero_particle_2%least_create_time) 00762 aero_particle_new%greatest_create_time = & 00763 max(aero_particle_1%greatest_create_time, & 00764 aero_particle_2%greatest_create_time) 00765 00766 end subroutine aero_particle_coagulate 00767 00768 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00769 00770 !> Determines the number of bytes required to pack the given value. 00771 integer function pmc_mpi_pack_size_aero_particle(val) 00772 00773 !> Value to pack. 00774 type(aero_particle_t), intent(in) :: val 00775 00776 pmc_mpi_pack_size_aero_particle = & 00777 pmc_mpi_pack_size_real_array(val%vol) & 00778 + pmc_mpi_pack_size_integer_array(val%n_orig_part) & 00779 + pmc_mpi_pack_size_integer(val%weight_group) & 00780 + pmc_mpi_pack_size_real(val%absorb_cross_sect) & 00781 + pmc_mpi_pack_size_real(val%scatter_cross_sect) & 00782 + pmc_mpi_pack_size_real(val%asymmetry) & 00783 + pmc_mpi_pack_size_complex(val%refract_shell) & 00784 + pmc_mpi_pack_size_complex(val%refract_core) & 00785 + pmc_mpi_pack_size_real(val%core_vol) & 00786 + pmc_mpi_pack_size_integer(val%water_hyst_leg) & 00787 + pmc_mpi_pack_size_integer(val%id) & 00788 + pmc_mpi_pack_size_real(val%least_create_time) & 00789 + pmc_mpi_pack_size_real(val%greatest_create_time) 00790 00791 end function pmc_mpi_pack_size_aero_particle 00792 00793 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00794 00795 !> Packs the given value into the buffer, advancing position. 00796 subroutine pmc_mpi_pack_aero_particle(buffer, position, val) 00797 00798 !> Memory buffer. 00799 character, intent(inout) :: buffer(:) 00800 !> Current buffer position. 00801 integer, intent(inout) :: position 00802 !> Value to pack. 00803 type(aero_particle_t), intent(in) :: val 00804 00805 #ifdef PMC_USE_MPI 00806 integer :: prev_position 00807 00808 prev_position = position 00809 call pmc_mpi_pack_real_array(buffer, position, val%vol) 00810 call pmc_mpi_pack_integer_array(buffer, position, val%n_orig_part) 00811 call pmc_mpi_pack_integer(buffer, position, val%weight_group) 00812 call pmc_mpi_pack_real(buffer, position, val%absorb_cross_sect) 00813 call pmc_mpi_pack_real(buffer, position, val%scatter_cross_sect) 00814 call pmc_mpi_pack_real(buffer, position, val%asymmetry) 00815 call pmc_mpi_pack_complex(buffer, position, val%refract_shell) 00816 call pmc_mpi_pack_complex(buffer, position, val%refract_core) 00817 call pmc_mpi_pack_real(buffer, position, val%core_vol) 00818 call pmc_mpi_pack_integer(buffer, position, val%water_hyst_leg) 00819 call pmc_mpi_pack_integer(buffer, position, val%id) 00820 call pmc_mpi_pack_real(buffer, position, val%least_create_time) 00821 call pmc_mpi_pack_real(buffer, position, val%greatest_create_time) 00822 call assert(810223998, position - prev_position & 00823 <= pmc_mpi_pack_size_aero_particle(val)) 00824 #endif 00825 00826 end subroutine pmc_mpi_pack_aero_particle 00827 00828 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00829 00830 !> Unpacks the given value from the buffer, advancing position. 00831 subroutine pmc_mpi_unpack_aero_particle(buffer, position, val) 00832 00833 !> Memory buffer. 00834 character, intent(inout) :: buffer(:) 00835 !> Current buffer position. 00836 integer, intent(inout) :: position 00837 !> Value to pack. 00838 type(aero_particle_t), intent(inout) :: val 00839 00840 #ifdef PMC_USE_MPI 00841 integer :: prev_position 00842 00843 prev_position = position 00844 call pmc_mpi_unpack_real_array(buffer, position, val%vol) 00845 call pmc_mpi_unpack_integer_array(buffer, position, val%n_orig_part) 00846 call pmc_mpi_unpack_integer(buffer, position, val%weight_group) 00847 call pmc_mpi_unpack_real(buffer, position, val%absorb_cross_sect) 00848 call pmc_mpi_unpack_real(buffer, position, val%scatter_cross_sect) 00849 call pmc_mpi_unpack_real(buffer, position, val%asymmetry) 00850 call pmc_mpi_unpack_complex(buffer, position, val%refract_shell) 00851 call pmc_mpi_unpack_complex(buffer, position, val%refract_core) 00852 call pmc_mpi_unpack_real(buffer, position, val%core_vol) 00853 call pmc_mpi_unpack_integer(buffer, position, val%water_hyst_leg) 00854 call pmc_mpi_unpack_integer(buffer, position, val%id) 00855 call pmc_mpi_unpack_real(buffer, position, val%least_create_time) 00856 call pmc_mpi_unpack_real(buffer, position, val%greatest_create_time) 00857 call assert(287447241, position - prev_position & 00858 <= pmc_mpi_pack_size_aero_particle(val)) 00859 #endif 00860 00861 end subroutine pmc_mpi_unpack_aero_particle 00862 00863 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00864 00865 end module pmc_aero_particle