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