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 the. 00004 00005 !> \file 00006 !> The pmc_aero_binned module. 00007 00008 !> The aero_binned_t structure and associated subroutines. 00009 module pmc_aero_binned 00010 00011 use pmc_bin_grid 00012 use pmc_aero_particle 00013 use pmc_spec_file 00014 use pmc_util 00015 use pmc_bin_grid 00016 use pmc_aero_dist 00017 use pmc_mpi 00018 use pmc_aero_data 00019 #ifdef PMC_USE_MPI 00020 use mpi 00021 #endif 00022 00023 !> Aerosol number and volume distributions stored per bin. 00024 !! 00025 !! These quantities are densities both in volume (per m^3) and in 00026 !! radius (per log_width). The total concentration per volume is computed as 00027 !! sum(aero_binned\%num_conc * bin_grid\%log_width). 00028 !! 00029 !! An aero_binned_t is similar to an aero_dist_t in that they both 00030 !! store binned aerosol distributions. The difference is that an 00031 !! aero_dist_t has the same composition in every bin, whereas an 00032 !! aero_binned_t can have aerosol composition that varies per bin. 00033 type aero_binned_t 00034 !> Number concentration per bin (#/m^3/log_width). 00035 !! Array length is typically \c bin_grid\%n_bin. 00036 real(kind=dp), pointer :: num_conc(:) 00037 !> Volume concentration per bin and per species (m^3/m^3/log_width). 00038 !! Array size is typically \c bin_grid\%n_bin x \c aero_data\%n_spec. 00039 real(kind=dp), pointer :: vol_conc(:,:) 00040 end type aero_binned_t 00041 00042 contains 00043 00044 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00045 00046 !> Allocate an aero_binned_t. 00047 subroutine aero_binned_allocate(aero_binned) 00048 00049 !> Structure to be allocated. 00050 type(aero_binned_t), intent(out) :: aero_binned 00051 00052 allocate(aero_binned%num_conc(0)) 00053 allocate(aero_binned%vol_conc(0, 0)) 00054 00055 end subroutine aero_binned_allocate 00056 00057 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00058 00059 !> Allocate an aero_binned_t of the given size. 00060 subroutine aero_binned_allocate_size(aero_binned, n_bin, n_spec) 00061 00062 !> Structure to be allocated. 00063 type(aero_binned_t), intent(out) :: aero_binned 00064 !> Number of aerosol bins to allocate (typically \c bin_grid%%n_bin). 00065 integer, intent(in) :: n_bin 00066 !> Number of aerosol species to allocate (typically 00067 !> \c aero_data%%n_spec). 00068 integer, intent(in) :: n_spec 00069 00070 allocate(aero_binned%num_conc(n_bin)) 00071 allocate(aero_binned%vol_conc(n_bin, n_spec)) 00072 call aero_binned_zero(aero_binned) 00073 00074 end subroutine aero_binned_allocate_size 00075 00076 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00077 00078 !> Free internal memory in an aero_binned_t structure. 00079 subroutine aero_binned_deallocate(aero_binned) 00080 00081 !> Structure to free. 00082 type(aero_binned_t), intent(inout) :: aero_binned 00083 00084 deallocate(aero_binned%num_conc) 00085 deallocate(aero_binned%vol_conc) 00086 00087 end subroutine aero_binned_deallocate 00088 00089 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00090 00091 !> Set all internal data in an aero_binned_t structure to zero. 00092 subroutine aero_binned_zero(aero_binned) 00093 00094 !> Structure to zero. 00095 type(aero_binned_t), intent(inout) :: aero_binned 00096 00097 aero_binned%num_conc = 0d0 00098 aero_binned%vol_conc = 0d0 00099 00100 end subroutine aero_binned_zero 00101 00102 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00103 00104 !> Update aero_binned_t structure for the addition of the given 00105 !> particle whose bin is also given. 00106 !! 00107 !! If the bin of the particle is not known the more expensive 00108 !! aero_binned_add_particle() can be used. 00109 subroutine aero_binned_add_particle_in_bin(aero_binned, bin_grid, & 00110 bin, comp_vol, aero_particle) 00111 00112 !> Structure to update with the new particle. 00113 type(aero_binned_t), intent(inout) :: aero_binned 00114 !> Bin grid. 00115 type(bin_grid_t), intent(in) :: bin_grid 00116 !> Bin number that new particle is in (must be correct). 00117 integer, intent(in) :: bin 00118 !> Computational volume (m^3). 00119 real(kind=dp), intent(in) :: comp_vol 00120 !> Particle to add. 00121 type(aero_particle_t), intent(in) :: aero_particle 00122 00123 aero_binned%num_conc(bin) = aero_binned%num_conc(bin) & 00124 + 1d0 / comp_vol / bin_grid%log_width 00125 aero_binned%vol_conc(bin,:) = aero_binned%vol_conc(bin,:) & 00126 + aero_particle%vol / comp_vol / bin_grid%log_width 00127 00128 end subroutine aero_binned_add_particle_in_bin 00129 00130 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00131 00132 !> Update aero_binned_t structure for the addition of the given 00133 !> particle. 00134 !! 00135 !! If the correct bin for the particle is already known then it is 00136 !! cheaper to call aero_binned_add_particle_in_bin(). 00137 subroutine aero_binned_add_particle(aero_binned, bin_grid, & 00138 comp_vol, aero_particle) 00139 00140 !> Structure to update with the new particle. 00141 type(aero_binned_t), intent(inout) :: aero_binned 00142 !> Bin grid. 00143 type(bin_grid_t), intent(in) :: bin_grid 00144 !> Computational volume (m^3). 00145 real(kind=dp), intent(in) :: comp_vol 00146 !> Particle to add. 00147 type(aero_particle_t), intent(in) :: aero_particle 00148 00149 call aero_binned_add_particle_in_bin(aero_binned, bin_grid, & 00150 aero_particle_in_bin(aero_particle, bin_grid), comp_vol, & 00151 aero_particle) 00152 00153 end subroutine aero_binned_add_particle 00154 00155 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00156 00157 !> Update aero_binned_t structure for the removal of the given 00158 !> particle whose bin is also given. 00159 !! 00160 !! If the bin of the particle is not known the more expensive 00161 !! aero_binned_remove_particle() can be used. 00162 subroutine aero_binned_remove_particle_in_bin(aero_binned, bin_grid, & 00163 bin, comp_vol, aero_particle) 00164 00165 !> Structure to remove the particle from. 00166 type(aero_binned_t), intent(inout) :: aero_binned 00167 !> Bin grid. 00168 type(bin_grid_t), intent(in) :: bin_grid 00169 !> Bin number of the aero_particle. 00170 integer, intent(in) :: bin 00171 !> Computational volume (m^3). 00172 real(kind=dp), intent(in) :: comp_vol 00173 !> Particle to remove. 00174 type(aero_particle_t), intent(in) :: aero_particle 00175 00176 aero_binned%num_conc(bin) = aero_binned%num_conc(bin) & 00177 - 1d0 / comp_vol / bin_grid%log_width 00178 aero_binned%vol_conc(bin,:) = aero_binned%vol_conc(bin,:) & 00179 - aero_particle%vol / comp_vol / bin_grid%log_width 00180 00181 end subroutine aero_binned_remove_particle_in_bin 00182 00183 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00184 00185 !> Update the aero_binned_t structure for the removal of the given 00186 !> particle. 00187 !! 00188 !! If the correct bin for the particle is already known then it is 00189 !! cheaper to call aero_binned_remove_particle_in_bin(). 00190 subroutine aero_binned_remove_particle(aero_binned, bin_grid, & 00191 comp_vol, aero_particle) 00192 00193 !> Structure to remove the particle from. 00194 type(aero_binned_t), intent(inout) :: aero_binned 00195 !> Bin grid. 00196 type(bin_grid_t), intent(in) :: bin_grid 00197 !> Computational volume (m^3). 00198 real(kind=dp), intent(in) :: comp_vol 00199 !> Particle to remove. 00200 type(aero_particle_t), intent(in) :: aero_particle 00201 00202 call aero_binned_remove_particle_in_bin(aero_binned, bin_grid, & 00203 aero_particle_in_bin(aero_particle, bin_grid), comp_vol, & 00204 aero_particle) 00205 00206 end subroutine aero_binned_remove_particle 00207 00208 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00209 00210 !> Add two aero_binned_t structures together. 00211 !! 00212 !! Symbolically does aero_binned = aero_binned + aero_binned_delta. 00213 subroutine aero_binned_add(aero_binned, aero_binned_delta) 00214 00215 !> Base aero_binned_t structure that will be added to. 00216 type(aero_binned_t), intent(inout) :: aero_binned 00217 !> Structure to add to aero_binned. 00218 type(aero_binned_t), intent(in) :: aero_binned_delta 00219 00220 aero_binned%num_conc = aero_binned%num_conc + aero_binned_delta%num_conc 00221 aero_binned%vol_conc = aero_binned%vol_conc + aero_binned_delta%vol_conc 00222 00223 end subroutine aero_binned_add 00224 00225 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00226 00227 !> Subtract one aero_binned_t structure from another. 00228 !! 00229 !! Symbolically does aero_binned = aero_binned - aero_binned_delta. 00230 subroutine aero_binned_sub(aero_binned, aero_binned_delta) 00231 00232 !> Base aero_binned_t structure that will be subtracted from. 00233 type(aero_binned_t), intent(inout) :: aero_binned 00234 !> Structure to subtract from aero_binned. 00235 type(aero_binned_t), intent(in) :: aero_binned_delta 00236 00237 aero_binned%num_conc = aero_binned%num_conc - aero_binned_delta%num_conc 00238 aero_binned%vol_conc = aero_binned%vol_conc - aero_binned_delta%vol_conc 00239 00240 end subroutine aero_binned_sub 00241 00242 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00243 00244 !> Scale an aero_binned_t by a real number. 00245 !! 00246 !! Symbolically does aero_binned = aero_binned * alpha. 00247 subroutine aero_binned_scale(aero_binned, alpha) 00248 00249 !> Base aero_binned to scale. 00250 type(aero_binned_t), intent(inout) :: aero_binned 00251 !> Scale factor. 00252 real(kind=dp), intent(in) :: alpha 00253 00254 aero_binned%num_conc = aero_binned%num_conc * alpha 00255 aero_binned%vol_conc = aero_binned%vol_conc * alpha 00256 00257 end subroutine aero_binned_scale 00258 00259 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00260 00261 !> Copy one aero_binned_t structure to another. 00262 !! 00263 !! Symbolically does aero_binned_to = aero_binned_from. 00264 subroutine aero_binned_copy(aero_binned_from, aero_binned_to) 00265 00266 !> Base aero_binned_t structure to copy from. 00267 type(aero_binned_t), intent(in) :: aero_binned_from 00268 !> Structure to copy to. 00269 type(aero_binned_t), intent(inout) :: aero_binned_to 00270 00271 integer :: n_bin, n_spec 00272 00273 n_bin = size(aero_binned_from%vol_conc, 1) 00274 n_spec = size(aero_binned_from%vol_conc, 2) 00275 call aero_binned_deallocate(aero_binned_to) 00276 call aero_binned_allocate_size(aero_binned_to, n_bin, n_spec) 00277 aero_binned_to%num_conc = aero_binned_from%num_conc 00278 aero_binned_to%vol_conc = aero_binned_from%vol_conc 00279 00280 end subroutine aero_binned_copy 00281 00282 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00283 00284 !> Add an aero_dist_t to an aero_binned_t. 00285 !! 00286 !! Symbolically does aero_binned = aero_binned + aero_dist. 00287 subroutine aero_binned_add_aero_dist(aero_binned, bin_grid, & 00288 aero_data, aero_dist) 00289 00290 !> Base aero_binned_t structure to add to. 00291 type(aero_binned_t), intent(inout) :: aero_binned 00292 !> Bin grid. 00293 type(bin_grid_t), intent(in) :: bin_grid 00294 !> Aerosol data. 00295 type(aero_data_t), intent(in) :: aero_data 00296 !> The aero_dist_t structure to add. 00297 type(aero_dist_t), intent(in) :: aero_dist 00298 00299 real(kind=dp) :: dist_num_conc(size(aero_binned%num_conc, 1)) 00300 real(kind=dp) :: dist_vol_conc(size(aero_binned%vol_conc, 1), 00301 size(aero_binned%vol_conc, 2)) 00302 00303 call aero_dist_num_conc(aero_dist, bin_grid, aero_data, & 00304 dist_num_conc) 00305 call aero_dist_vol_conc(aero_dist, bin_grid, aero_data, & 00306 dist_vol_conc) 00307 aero_binned%num_conc = aero_binned%num_conc + dist_num_conc 00308 aero_binned%vol_conc = aero_binned%vol_conc + dist_vol_conc 00309 00310 end subroutine aero_binned_add_aero_dist 00311 00312 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00313 00314 !> Determine the number of bytes required to pack the structure. 00315 !! 00316 !! See pmc_mpi for usage details. 00317 integer function pmc_mpi_pack_size_aero_binned(val) 00318 00319 !> Structure to pack. 00320 type(aero_binned_t), intent(in) :: val 00321 00322 pmc_mpi_pack_size_aero_binned = & 00323 pmc_mpi_pack_size_real_array(val%num_conc) & 00324 + pmc_mpi_pack_size_real_array_2d(val%vol_conc) 00325 00326 end function pmc_mpi_pack_size_aero_binned 00327 00328 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00329 00330 !> Pack the structure into the buffer and advance position. 00331 !! 00332 !! See pmc_mpi for usage details. 00333 subroutine pmc_mpi_pack_aero_binned(buffer, position, val) 00334 00335 !> Memory buffer. 00336 character, intent(inout) :: buffer(:) 00337 !> Current buffer position. 00338 integer, intent(inout) :: position 00339 !> Structure to pack. 00340 type(aero_binned_t), intent(in) :: val 00341 00342 #ifdef PMC_USE_MPI 00343 integer :: prev_position 00344 00345 prev_position = position 00346 call pmc_mpi_pack_real_array(buffer, position, val%num_conc) 00347 call pmc_mpi_pack_real_array_2d(buffer, position, val%vol_conc) 00348 call assert(348207873, & 00349 position - prev_position <= pmc_mpi_pack_size_aero_binned(val)) 00350 #endif 00351 00352 end subroutine pmc_mpi_pack_aero_binned 00353 00354 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00355 00356 !> Unpack the structure from the buffer and advance position. 00357 !! 00358 !! See pmc_mpi for usage details. 00359 subroutine pmc_mpi_unpack_aero_binned(buffer, position, val) 00360 00361 !> Memory buffer. 00362 character, intent(inout) :: buffer(:) 00363 !> Current buffer position. 00364 integer, intent(inout) :: position 00365 !> Structure to unpack into (must not be allocated). 00366 type(aero_binned_t), intent(inout) :: val 00367 00368 #ifdef PMC_USE_MPI 00369 integer :: prev_position 00370 00371 prev_position = position 00372 call pmc_mpi_unpack_real_array(buffer, position, val%num_conc) 00373 call pmc_mpi_unpack_real_array_2d(buffer, position, val%vol_conc) 00374 call assert(878267066, & 00375 position - prev_position <= pmc_mpi_pack_size_aero_binned(val)) 00376 #endif 00377 00378 end subroutine pmc_mpi_unpack_aero_binned 00379 00380 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00381 00382 !> Computes the average of the structure across all processes, 00383 !> storing the result on the root process. 00384 subroutine pmc_mpi_reduce_avg_aero_binned(val, val_avg) 00385 00386 !> Per-process value to average. 00387 type(aero_binned_t), intent(in) :: val 00388 !> Averaged result (only valid on root process). 00389 type(aero_binned_t), intent(inout) :: val_avg 00390 00391 call pmc_mpi_reduce_avg_real_array(val%num_conc, val_avg%num_conc) 00392 call pmc_mpi_reduce_avg_real_array_2d(val%vol_conc, val_avg%vol_conc) 00393 00394 end subroutine pmc_mpi_reduce_avg_aero_binned 00395 00396 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00397 00398 !> Write full state. 00399 subroutine aero_binned_output_netcdf(aero_binned, ncid, bin_grid, & 00400 aero_data) 00401 00402 !> Aero_binned to write. 00403 type(aero_binned_t), intent(in) :: aero_binned 00404 !> NetCDF file ID, in data mode. 00405 integer, intent(in) :: ncid 00406 !> bin_grid structure. 00407 type(bin_grid_t), intent(in) :: bin_grid 00408 !> aero_data structure. 00409 type(aero_data_t), intent(in) :: aero_data 00410 00411 integer :: dimid_aero_diam, dimid_aero_species 00412 real(kind=dp) :: mass_den(bin_grid%n_bin, aero_data%n_spec) 00413 integer :: i_bin 00414 00415 !> \page output_format_aero_binned Output File Format: Aerosol Binned Sectional State 00416 !! 00417 !! The aerosol size distributions (number and mass) are stored on 00418 !! a logarmithmic grid (see the \ref output_format_bin_grid 00419 !! section). To compute the total number or mass concentration, 00420 !! compute the sum over \c i of <tt>aero_number_concentration(i) * 00421 !! aero_diam_widths(i)</tt>, for example. 00422 !! 00423 !! The aerosol binned sectional state uses the \c aero_species 00424 !! NetCDF dimension as specified in the \ref 00425 !! output_format_aero_data section, as well as the \c aero_diam 00426 !! NetCDF dimension specified in the \ref output_format_bin_grid 00427 !! section. 00428 !! 00429 !! The aerosol binned sectional state NetCDF variables are: 00430 !! - \b aero_number_concentration (unit 1/m^3, dim \c aero_diam): the 00431 !! number size distribution for the aerosol population, 00432 !! \f$ dN(r)/d\ln r \f$, per bin 00433 !! - \b aero_mass_concentration (unit kg/m^3, dim 00434 !! <tt>dimid_aero_diam x dimid_aero_species</tt>): the mass size 00435 !! distribution for the aerosol population, 00436 !! \f$ dM(r,s)/d\ln r \f$, per bin and per species 00437 00438 do i_bin = 1,bin_grid%n_bin 00439 mass_den(i_bin,:) = aero_binned%vol_conc(i_bin,:) & 00440 * aero_data%density 00441 end do 00442 00443 call bin_grid_netcdf_dim_aero_diam(bin_grid, ncid, & 00444 dimid_aero_diam) 00445 call aero_data_netcdf_dim_aero_species(aero_data, ncid, & 00446 dimid_aero_species) 00447 00448 call pmc_nc_write_real_1d(ncid, aero_binned%num_conc, & 00449 "aero_number_concentration", (/ dimid_aero_diam /), & 00450 unit="1/m^3", & 00451 long_name="aerosol number size concentration distribution", & 00452 description="logarithmic number size concentration, " & 00453 // "d N(r)/d ln r --- multiply by aero_diam_widths(i) " & 00454 // "and sum over i to compute the total number concentration") 00455 call pmc_nc_write_real_2d(ncid, mass_den, & 00456 "aero_mass_concentration", & 00457 (/ dimid_aero_diam, dimid_aero_species /), unit="kg/m^3", & 00458 long_name="aerosol number size concentration distribution", & 00459 description="logarithmic mass size concentration per species, " & 00460 // "d M(r,s)/d ln r --- multiply by aero_diam_widths(i) " & 00461 // "and sum over i to compute the total mass concentration of " & 00462 // "species s") 00463 00464 end subroutine aero_binned_output_netcdf 00465 00466 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00467 00468 !> Read full state. 00469 subroutine aero_binned_input_netcdf(aero_binned, ncid, bin_grid, & 00470 aero_data) 00471 00472 !> Aero_binned to write. 00473 type(aero_binned_t), intent(inout) :: aero_binned 00474 !> NetCDF file ID, in data mode. 00475 integer, intent(in) :: ncid 00476 !> bin_grid structure. 00477 type(bin_grid_t), intent(in) :: bin_grid 00478 !> aero_data structure. 00479 type(aero_data_t), intent(in) :: aero_data 00480 00481 real(kind=dp) :: mass_den(bin_grid%n_bin, aero_data%n_spec) 00482 integer :: i_bin 00483 00484 call pmc_nc_read_real_1d(ncid, aero_binned%num_conc, & 00485 "aero_number_concentration") 00486 call pmc_nc_read_real_2d(ncid, mass_den, & 00487 "aero_mass_concentration") 00488 00489 do i_bin = 1,bin_grid%n_bin 00490 aero_binned%vol_conc(i_bin,:) = mass_den(i_bin,:) & 00491 / aero_data%density 00492 end do 00493 00494 end subroutine aero_binned_input_netcdf 00495 00496 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00497 00498 end module pmc_aero_binned