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 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 !> Add two aero_binned_t structures together. 00105 !! 00106 !! Symbolically does aero_binned = aero_binned + aero_binned_delta. 00107 subroutine aero_binned_add(aero_binned, aero_binned_delta) 00108 00109 !> Base aero_binned_t structure that will be added to. 00110 type(aero_binned_t), intent(inout) :: aero_binned 00111 !> Structure to add to aero_binned. 00112 type(aero_binned_t), intent(in) :: aero_binned_delta 00113 00114 aero_binned%num_conc = aero_binned%num_conc + aero_binned_delta%num_conc 00115 aero_binned%vol_conc = aero_binned%vol_conc + aero_binned_delta%vol_conc 00116 00117 end subroutine aero_binned_add 00118 00119 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00120 00121 !> Subtract one aero_binned_t structure from another. 00122 !! 00123 !! Symbolically does aero_binned = aero_binned - aero_binned_delta. 00124 subroutine aero_binned_sub(aero_binned, aero_binned_delta) 00125 00126 !> Base aero_binned_t structure that will be subtracted from. 00127 type(aero_binned_t), intent(inout) :: aero_binned 00128 !> Structure to subtract from aero_binned. 00129 type(aero_binned_t), intent(in) :: aero_binned_delta 00130 00131 aero_binned%num_conc = aero_binned%num_conc - aero_binned_delta%num_conc 00132 aero_binned%vol_conc = aero_binned%vol_conc - aero_binned_delta%vol_conc 00133 00134 end subroutine aero_binned_sub 00135 00136 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00137 00138 !> Scale an aero_binned_t by a real number. 00139 !! 00140 !! Symbolically does aero_binned = aero_binned * alpha. 00141 subroutine aero_binned_scale(aero_binned, alpha) 00142 00143 !> Base aero_binned to scale. 00144 type(aero_binned_t), intent(inout) :: aero_binned 00145 !> Scale factor. 00146 real(kind=dp), intent(in) :: alpha 00147 00148 aero_binned%num_conc = aero_binned%num_conc * alpha 00149 aero_binned%vol_conc = aero_binned%vol_conc * alpha 00150 00151 end subroutine aero_binned_scale 00152 00153 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00154 00155 !> Copy one aero_binned_t structure to another. 00156 !! 00157 !! Symbolically does aero_binned_to = aero_binned_from. 00158 subroutine aero_binned_copy(aero_binned_from, aero_binned_to) 00159 00160 !> Base aero_binned_t structure to copy from. 00161 type(aero_binned_t), intent(in) :: aero_binned_from 00162 !> Structure to copy to. 00163 type(aero_binned_t), intent(inout) :: aero_binned_to 00164 00165 integer :: n_bin, n_spec 00166 00167 n_bin = size(aero_binned_from%vol_conc, 1) 00168 n_spec = size(aero_binned_from%vol_conc, 2) 00169 call aero_binned_deallocate(aero_binned_to) 00170 call aero_binned_allocate_size(aero_binned_to, n_bin, n_spec) 00171 aero_binned_to%num_conc = aero_binned_from%num_conc 00172 aero_binned_to%vol_conc = aero_binned_from%vol_conc 00173 00174 end subroutine aero_binned_copy 00175 00176 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00177 00178 !> Add an aero_dist_t to an aero_binned_t. 00179 !! 00180 !! Symbolically does aero_binned = aero_binned + aero_dist. 00181 subroutine aero_binned_add_aero_dist(aero_binned, bin_grid, & 00182 aero_data, aero_dist) 00183 00184 !> Base aero_binned_t structure to add to. 00185 type(aero_binned_t), intent(inout) :: aero_binned 00186 !> Bin grid. 00187 type(bin_grid_t), intent(in) :: bin_grid 00188 !> Aerosol data. 00189 type(aero_data_t), intent(in) :: aero_data 00190 !> The aero_dist_t structure to add. 00191 type(aero_dist_t), intent(in) :: aero_dist 00192 00193 real(kind=dp) :: dist_num_conc(size(aero_binned%num_conc, 1)) 00194 real(kind=dp) :: dist_vol_conc(size(aero_binned%vol_conc, 1), 00195 size(aero_binned%vol_conc, 2)) 00196 00197 call aero_dist_num_conc(aero_dist, bin_grid, aero_data, & 00198 dist_num_conc) 00199 call aero_dist_vol_conc(aero_dist, bin_grid, aero_data, & 00200 dist_vol_conc) 00201 aero_binned%num_conc = aero_binned%num_conc + dist_num_conc 00202 aero_binned%vol_conc = aero_binned%vol_conc + dist_vol_conc 00203 00204 end subroutine aero_binned_add_aero_dist 00205 00206 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00207 00208 !> Determine the number of bytes required to pack the structure. 00209 !! 00210 !! See pmc_mpi for usage details. 00211 integer function pmc_mpi_pack_size_aero_binned(val) 00212 00213 !> Structure to pack. 00214 type(aero_binned_t), intent(in) :: val 00215 00216 pmc_mpi_pack_size_aero_binned = & 00217 pmc_mpi_pack_size_real_array(val%num_conc) & 00218 + pmc_mpi_pack_size_real_array_2d(val%vol_conc) 00219 00220 end function pmc_mpi_pack_size_aero_binned 00221 00222 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00223 00224 !> Pack the structure into the buffer and advance position. 00225 !! 00226 !! See pmc_mpi for usage details. 00227 subroutine pmc_mpi_pack_aero_binned(buffer, position, val) 00228 00229 !> Memory buffer. 00230 character, intent(inout) :: buffer(:) 00231 !> Current buffer position. 00232 integer, intent(inout) :: position 00233 !> Structure to pack. 00234 type(aero_binned_t), intent(in) :: val 00235 00236 #ifdef PMC_USE_MPI 00237 integer :: prev_position 00238 00239 prev_position = position 00240 call pmc_mpi_pack_real_array(buffer, position, val%num_conc) 00241 call pmc_mpi_pack_real_array_2d(buffer, position, val%vol_conc) 00242 call assert(348207873, & 00243 position - prev_position <= pmc_mpi_pack_size_aero_binned(val)) 00244 #endif 00245 00246 end subroutine pmc_mpi_pack_aero_binned 00247 00248 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00249 00250 !> Unpack the structure from the buffer and advance position. 00251 !! 00252 !! See pmc_mpi for usage details. 00253 subroutine pmc_mpi_unpack_aero_binned(buffer, position, val) 00254 00255 !> Memory buffer. 00256 character, intent(inout) :: buffer(:) 00257 !> Current buffer position. 00258 integer, intent(inout) :: position 00259 !> Structure to unpack into (must not be allocated). 00260 type(aero_binned_t), intent(inout) :: val 00261 00262 #ifdef PMC_USE_MPI 00263 integer :: prev_position 00264 00265 prev_position = position 00266 call pmc_mpi_unpack_real_array(buffer, position, val%num_conc) 00267 call pmc_mpi_unpack_real_array_2d(buffer, position, val%vol_conc) 00268 call assert(878267066, & 00269 position - prev_position <= pmc_mpi_pack_size_aero_binned(val)) 00270 #endif 00271 00272 end subroutine pmc_mpi_unpack_aero_binned 00273 00274 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00275 00276 !> Computes the average of the structure across all processes, 00277 !> storing the result on the root process. 00278 subroutine pmc_mpi_reduce_avg_aero_binned(val, val_avg) 00279 00280 !> Per-process value to average. 00281 type(aero_binned_t), intent(in) :: val 00282 !> Averaged result (only valid on root process). 00283 type(aero_binned_t), intent(inout) :: val_avg 00284 00285 call pmc_mpi_reduce_avg_real_array(val%num_conc, val_avg%num_conc) 00286 call pmc_mpi_reduce_avg_real_array_2d(val%vol_conc, val_avg%vol_conc) 00287 00288 end subroutine pmc_mpi_reduce_avg_aero_binned 00289 00290 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00291 00292 !> Write full state. 00293 subroutine aero_binned_output_netcdf(aero_binned, ncid, bin_grid, & 00294 aero_data) 00295 00296 !> Aero_binned to write. 00297 type(aero_binned_t), intent(in) :: aero_binned 00298 !> NetCDF file ID, in data mode. 00299 integer, intent(in) :: ncid 00300 !> bin_grid structure. 00301 type(bin_grid_t), intent(in) :: bin_grid 00302 !> aero_data structure. 00303 type(aero_data_t), intent(in) :: aero_data 00304 00305 integer :: dimid_aero_diam, dimid_aero_species 00306 real(kind=dp) :: mass_conc(bin_grid%n_bin, aero_data%n_spec) 00307 integer :: i_bin 00308 00309 !> \page output_format_aero_binned Output File Format: Aerosol Binned Sectional State 00310 !! 00311 !! The aerosol size distributions (number and mass) are stored on 00312 !! a logarmithmic grid (see the \ref output_format_bin_grid 00313 !! section). To compute the total number or mass concentration, 00314 !! compute the sum over \c i of <tt>aero_number_concentration(i) * 00315 !! aero_diam_widths(i)</tt>, for example. 00316 !! 00317 !! The aerosol binned sectional state uses the \c aero_species 00318 !! NetCDF dimension as specified in the \ref 00319 !! output_format_aero_data section, as well as the \c aero_diam 00320 !! NetCDF dimension specified in the \ref output_format_bin_grid 00321 !! section. 00322 !! 00323 !! The aerosol binned sectional state NetCDF variables are: 00324 !! - \b aero_number_concentration (unit 1/m^3, dim \c aero_diam): the 00325 !! number size distribution for the aerosol population, 00326 !! \f$ dN(r)/d\ln r \f$, per bin 00327 !! - \b aero_mass_concentration (unit kg/m^3, dim 00328 !! <tt>dimid_aero_diam x dimid_aero_species</tt>): the mass size 00329 !! distribution for the aerosol population, 00330 !! \f$ dM(r,s)/d\ln r \f$, per bin and per species 00331 00332 do i_bin = 1,bin_grid%n_bin 00333 mass_conc(i_bin,:) = aero_binned%vol_conc(i_bin,:) * aero_data%density 00334 end do 00335 00336 call bin_grid_netcdf_dim_aero_diam(bin_grid, ncid, dimid_aero_diam) 00337 call aero_data_netcdf_dim_aero_species(aero_data, ncid, dimid_aero_species) 00338 00339 call pmc_nc_write_real_1d(ncid, aero_binned%num_conc, & 00340 "aero_number_concentration", (/ dimid_aero_diam /), & 00341 unit="1/m^3", & 00342 long_name="aerosol number size concentration distribution", & 00343 description="logarithmic number size concentration, " & 00344 // "d N(r)/d ln D --- multiply by aero_diam_widths(i) " & 00345 // "and sum over i to compute the total number concentration") 00346 call pmc_nc_write_real_2d(ncid, mass_conc, & 00347 "aero_mass_concentration", & 00348 (/ dimid_aero_diam, dimid_aero_species /), unit="kg/m^3", & 00349 long_name="aerosol mass size concentration distribution", & 00350 description="logarithmic mass size concentration per species, " & 00351 // "d M(r,s)/d ln D --- multiply by aero_diam_widths(i) " & 00352 // "and sum over i to compute the total mass concentration of " & 00353 // "species s") 00354 00355 end subroutine aero_binned_output_netcdf 00356 00357 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00358 00359 !> Read full state. 00360 subroutine aero_binned_input_netcdf(aero_binned, ncid, bin_grid, & 00361 aero_data) 00362 00363 !> Aero_binned to write. 00364 type(aero_binned_t), intent(inout) :: aero_binned 00365 !> NetCDF file ID, in data mode. 00366 integer, intent(in) :: ncid 00367 !> bin_grid structure. 00368 type(bin_grid_t), intent(in) :: bin_grid 00369 !> aero_data structure. 00370 type(aero_data_t), intent(in) :: aero_data 00371 00372 real(kind=dp) :: mass_conc(bin_grid%n_bin, aero_data%n_spec) 00373 integer :: i_bin 00374 00375 call aero_binned_deallocate(aero_binned) 00376 call aero_binned_allocate_size(aero_binned, bin_grid%n_bin, & 00377 aero_data%n_spec) 00378 00379 call pmc_nc_read_real_1d(ncid, aero_binned%num_conc, & 00380 "aero_number_concentration") 00381 call pmc_nc_read_real_2d(ncid, mass_conc, "aero_mass_concentration") 00382 00383 do i_bin = 1,bin_grid%n_bin 00384 aero_binned%vol_conc(i_bin,:) = mass_conc(i_bin,:) / aero_data%density 00385 end do 00386 00387 end subroutine aero_binned_input_netcdf 00388 00389 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00390 00391 end module pmc_aero_binned