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_bin_grid module. 00007 00008 !> The bin_grid_t structure and associated subroutines. 00009 module pmc_bin_grid 00010 00011 use pmc_constants 00012 use pmc_util 00013 use pmc_spec_file 00014 use pmc_mpi 00015 use pmc_netcdf 00016 #ifdef PMC_USE_MPI 00017 use mpi 00018 #endif 00019 00020 !> 1D grid of size bins. 00021 !! 00022 !! The grid of bins is logarithmically spaced in volume, an 00023 !! assumption that is quite heavily incorporated into the code. At 00024 !! some point in the future it would be nice to relax this 00025 !! assumption. 00026 type bin_grid_t 00027 !> Number of bins. 00028 integer :: n_bin 00029 !> Len n_bin, bin center radii (m^3). 00030 real(kind=dp), pointer :: center_radius(:) 00031 !> Len (n_bin + 1), bin edge radii (m^3). 00032 real(kind=dp), pointer :: edge_radius(:) 00033 !> Bin logarithmic width, equal to <tt>log(edge_radius(i+1)) - 00034 !> log(edge_radius(i))</tt> for any \c i (dimensionless). 00035 real(kind=dp) :: log_width 00036 end type bin_grid_t 00037 00038 contains 00039 00040 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00041 00042 !> Allocates a bin_grid. 00043 subroutine bin_grid_allocate(bin_grid) 00044 00045 !> Bin grid. 00046 type(bin_grid_t), intent(out) :: bin_grid 00047 00048 bin_grid%n_bin = 0 00049 allocate(bin_grid%center_radius(0)) 00050 allocate(bin_grid%edge_radius(0)) 00051 00052 end subroutine bin_grid_allocate 00053 00054 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00055 00056 !> Allocates a bin_grid of the given size. 00057 subroutine bin_grid_allocate_size(bin_grid, n_bin) 00058 00059 !> Bin grid. 00060 type(bin_grid_t), intent(out) :: bin_grid 00061 !> Number of bins. 00062 integer, intent(in) :: n_bin 00063 00064 bin_grid%n_bin = n_bin 00065 allocate(bin_grid%center_radius(n_bin)) 00066 allocate(bin_grid%edge_radius(n_bin + 1)) 00067 00068 end subroutine bin_grid_allocate_size 00069 00070 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00071 00072 !> Frees all memory. 00073 subroutine bin_grid_deallocate(bin_grid) 00074 00075 !> Bin_grid to free. 00076 type(bin_grid_t), intent(inout) :: bin_grid 00077 00078 deallocate(bin_grid%center_radius) 00079 deallocate(bin_grid%edge_radius) 00080 00081 end subroutine bin_grid_deallocate 00082 00083 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00084 00085 !> Copies a bin grid. 00086 subroutine bin_grid_copy(bin_grid_from, bin_grid_to) 00087 00088 !> Bin_grid to copy from. 00089 type(bin_grid_t), intent(in) :: bin_grid_from 00090 !> Bin_grid to copy to. 00091 type(bin_grid_t), intent(inout) :: bin_grid_to 00092 00093 if (bin_grid_from%n_bin /= bin_grid_to%n_bin) then 00094 call bin_grid_deallocate(bin_grid_to) 00095 call bin_grid_allocate_size(bin_grid_to, bin_grid_from%n_bin) 00096 end if 00097 bin_grid_to%center_radius = bin_grid_from%center_radius 00098 bin_grid_to%edge_radius = bin_grid_from%edge_radius 00099 bin_grid_to%log_width = bin_grid_from%log_width 00100 00101 end subroutine bin_grid_copy 00102 00103 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00104 00105 !> Returns the size of a bin in a grid. 00106 real(kind=dp) function bin_grid_bin_size(bin_grid, i_bin) 00107 00108 !> Bin grid. 00109 type(bin_grid_t), intent(in) :: bin_grid 00110 !> Bin number to return size of. 00111 integer, intent(in) :: i_bin 00112 00113 bin_grid_bin_size = bin_grid%log_width 00114 00115 end function bin_grid_bin_size 00116 00117 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00118 00119 !> Convert a concentration f(vol)d(vol) to f(ln(r))d(ln(r)) 00120 !> where vol = 4/3 pi r^3. 00121 subroutine vol_to_lnr(r, f_vol, f_lnr) 00122 00123 !> Radius (m). 00124 real(kind=dp), intent(in) :: r 00125 !> Concentration as a function of volume. 00126 real(kind=dp), intent(in) :: f_vol 00127 !> Concentration as a function of ln(r). 00128 real(kind=dp), intent(out) :: f_lnr 00129 00130 f_lnr = f_vol * 4d0 * const%pi * r**3 00131 00132 end subroutine vol_to_lnr 00133 00134 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00135 00136 !> Generates the bin grid given the range and number of bins. 00137 subroutine bin_grid_make(bin_grid, n_bin, r_min, r_max) 00138 00139 !> New bin grid. 00140 type(bin_grid_t), intent(inout) :: bin_grid 00141 !> Number of bins. 00142 integer, intent(in) :: n_bin 00143 !> Minimum radius (m^3). 00144 real(kind=dp), intent(in) :: r_min 00145 !> Minimum radius (m^3). 00146 real(kind=dp), intent(in) :: r_max 00147 00148 integer :: i_bin 00149 00150 call assert_msg(538534122, n_bin > 0, & 00151 "bin_grid requires a positive n_bin, not: " & 00152 // trim(integer_to_string(n_bin))) 00153 call assert_msg(966541762, r_min > 0d0, & 00154 "bin_grid requires a positive r_min, not: " & 00155 // trim(real_to_string(r_min))) 00156 call assert_msg(711537859, r_min < r_max, & 00157 "bin_grid requires r_min < r_max, not: " & 00158 // trim(real_to_string(r_min)) // " and " & 00159 // trim(real_to_string(r_max))) 00160 call bin_grid_deallocate(bin_grid) 00161 call bin_grid_allocate_size(bin_grid, n_bin) 00162 call logspace(r_min, r_max, bin_grid%edge_radius) 00163 do i_bin = 1,n_bin 00164 bin_grid%center_radius(i_bin) & 00165 = exp(0.5d0 * log(bin_grid%edge_radius(i_bin)) & 00166 + 0.5d0 * log(bin_grid%edge_radius(i_bin + 1))) 00167 end do 00168 bin_grid%log_width = (log(r_max) - log(r_min)) / real(n_bin, kind=dp) 00169 00170 end subroutine bin_grid_make 00171 00172 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00173 00174 !> Find the bin number that contains a given particle. 00175 !! 00176 !! This assumes logarithmically spaced bins. If a particle is below 00177 !! the smallest bin then its bin number is 0. If a particle is above 00178 !! the largest bin then its bin number is <tt>n_bin + 1</tt>. 00179 integer function bin_grid_particle_in_bin(bin_grid, radius) 00180 00181 !> Bin_grid. 00182 type(bin_grid_t), intent(in) :: bin_grid 00183 !> Radius of particle. 00184 real(kind=dp), intent(in) :: radius 00185 00186 call assert(448215689, bin_grid%n_bin >= 1) 00187 bin_grid_particle_in_bin = logspace_find(bin_grid%edge_radius(1), & 00188 bin_grid%edge_radius(bin_grid%n_bin + 1), bin_grid%n_bin + 1, & 00189 radius) 00190 00191 end function bin_grid_particle_in_bin 00192 00193 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00194 00195 !> Make a histogram with of the given weighted data, scaled by the 00196 !> bin sizes. 00197 subroutine bin_grid_histogram_1d(x_bin_grid, x_data, weight_data, hist) 00198 00199 !> x-axis bin grid. 00200 type(bin_grid_t), intent(in) :: x_bin_grid 00201 !> Data values on the x-axis. 00202 real(kind=dp), intent(in) :: x_data(:) 00203 !> Data value weights. 00204 real(kind=dp), intent(in) :: weight_data(size(x_data)) 00205 !> Histogram to compute. 00206 real(kind=dp), intent(out) :: hist(x_bin_grid%n_bin) 00207 00208 integer :: i_data, i_bin 00209 00210 hist = 0d0 00211 do i_data = 1,size(x_data) 00212 i_bin = bin_grid_particle_in_bin(x_bin_grid, x_data(i_data)) 00213 if ((i_bin >= 1) .and. (i_bin <= x_bin_grid%n_bin)) then 00214 hist(i_bin) = hist(i_bin) & 00215 + weight_data(i_data) / bin_grid_bin_size(x_bin_grid, i_bin) 00216 end if 00217 end do 00218 00219 end subroutine bin_grid_histogram_1d 00220 00221 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00222 00223 !> Read the specification for a bin_grid from a spec file and 00224 !> generate it. 00225 subroutine spec_file_read_bin_grid(file, bin_grid) 00226 00227 !> Spec file. 00228 type(spec_file_t), intent(inout) :: file 00229 !> Bin grid. 00230 type(bin_grid_t), intent(inout) :: bin_grid 00231 00232 integer :: n_bin 00233 real(kind=dp) :: d_min, d_max 00234 00235 !> \page input_format_bin_grid Input File Format: Diameter Axis Bin Grid 00236 !! 00237 !! The bin grid is logarithmic in diameter, consisting of 00238 !! \f$n_{\rm bin}\f$ bins with centers \f$c_i\f$ (\f$i = 00239 !! 1,\ldots,n_{\rm bin}\f$) and edges \f$e_i\f$ (\f$i = 00240 !! 1,\ldots,(n_{\rm bin} + 1)\f$) such that \f$e_{i+1}/e_i\f$ is a 00241 !! constant and \f$c_i/e_i = \sqrt{e_{i+1}/e_i}\f$. That is, 00242 !! \f$\ln(e_i)\f$ are uniformly spaced and \f$\ln(c_i)\f$ are the 00243 !! arithmetic centers. 00244 !! 00245 !! The diameter axis bin grid is specified by the parameters: 00246 !! - \b n_bin (integer): The number of bins \f$n_{\rm bin}\f$. 00247 !! - \b d_min (real, unit m): The left edge of the left-most bin, 00248 !! \f$e_1\f$. 00249 !! - \b d_max (real, unit m): The right edge of the right-most bin, 00250 !! \f$e_{n_{\rm bin} + 1}\f$. 00251 !! 00252 !! See also: 00253 !! - \ref spec_file_format --- the input file text format 00254 !! - \ref output_format_bin_grid --- the corresponding output format 00255 00256 call spec_file_read_integer(file, 'n_bin', n_bin) 00257 call spec_file_read_real(file, 'd_min', d_min) 00258 call spec_file_read_real(file, 'd_max', d_max) 00259 call bin_grid_make(bin_grid, n_bin, diam2rad(d_min), diam2rad(d_max)) 00260 00261 end subroutine spec_file_read_bin_grid 00262 00263 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00264 00265 !> Determines the number of bytes required to pack the given value. 00266 integer function pmc_mpi_pack_size_bin_grid(val) 00267 00268 !> Value to pack. 00269 type(bin_grid_t), intent(in) :: val 00270 00271 pmc_mpi_pack_size_bin_grid = & 00272 pmc_mpi_pack_size_integer(val%n_bin) & 00273 + pmc_mpi_pack_size_real_array(val%center_radius) & 00274 + pmc_mpi_pack_size_real_array(val%edge_radius) & 00275 + pmc_mpi_pack_size_real(val%log_width) 00276 00277 end function pmc_mpi_pack_size_bin_grid 00278 00279 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00280 00281 !> Packs the given value into the buffer, advancing position. 00282 subroutine pmc_mpi_pack_bin_grid(buffer, position, val) 00283 00284 !> Memory buffer. 00285 character, intent(inout) :: buffer(:) 00286 !> Current buffer position. 00287 integer, intent(inout) :: position 00288 !> Value to pack. 00289 type(bin_grid_t), intent(in) :: val 00290 00291 #ifdef PMC_USE_MPI 00292 integer :: prev_position 00293 00294 prev_position = position 00295 call pmc_mpi_pack_integer(buffer, position, val%n_bin) 00296 call pmc_mpi_pack_real_array(buffer, position, val%center_radius) 00297 call pmc_mpi_pack_real_array(buffer, position, val%edge_radius) 00298 call pmc_mpi_pack_real(buffer, position, val%log_width) 00299 call assert(385455586, & 00300 position - prev_position <= pmc_mpi_pack_size_bin_grid(val)) 00301 #endif 00302 00303 end subroutine pmc_mpi_pack_bin_grid 00304 00305 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00306 00307 !> Unpacks the given value from the buffer, advancing position. 00308 subroutine pmc_mpi_unpack_bin_grid(buffer, position, val) 00309 00310 !> Memory buffer. 00311 character, intent(inout) :: buffer(:) 00312 !> Current buffer position. 00313 integer, intent(inout) :: position 00314 !> Value to pack. 00315 type(bin_grid_t), intent(inout) :: val 00316 00317 #ifdef PMC_USE_MPI 00318 integer :: prev_position 00319 00320 prev_position = position 00321 call pmc_mpi_unpack_integer(buffer, position, val%n_bin) 00322 call pmc_mpi_unpack_real_array(buffer, position, val%center_radius) 00323 call pmc_mpi_unpack_real_array(buffer, position, val%edge_radius) 00324 call pmc_mpi_unpack_real(buffer, position, val%log_width) 00325 call assert(741838730, & 00326 position - prev_position <= pmc_mpi_pack_size_bin_grid(val)) 00327 #endif 00328 00329 end subroutine pmc_mpi_unpack_bin_grid 00330 00331 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00332 00333 !> Check whether all processors have the same value. 00334 logical function pmc_mpi_allequal_bin_grid(val) 00335 00336 !> Value to compare. 00337 type(bin_grid_t), intent(inout) :: val 00338 00339 #ifdef PMC_USE_MPI 00340 if (.not. pmc_mpi_allequal_integer(val%n_bin)) then 00341 pmc_mpi_allequal_bin_grid = .false. 00342 return 00343 end if 00344 00345 if (val%n_bin == 0) then 00346 pmc_mpi_allequal_bin_grid = .true. 00347 return 00348 end if 00349 00350 if (pmc_mpi_allequal_real(val%edge_radius(1)) & 00351 .and. pmc_mpi_allequal_real(val%edge_radius(val%n_bin))) then 00352 pmc_mpi_allequal_bin_grid = .true. 00353 else 00354 pmc_mpi_allequal_bin_grid = .false. 00355 end if 00356 #else 00357 pmc_mpi_allequal_bin_grid = .true. 00358 #endif 00359 00360 end function pmc_mpi_allequal_bin_grid 00361 00362 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00363 00364 !> Write the aero_diam dimension to the given NetCDF file if it is 00365 !> not already present and in any case return the associated dimid. 00366 subroutine bin_grid_netcdf_dim_aero_diam(bin_grid, ncid, & 00367 dimid_aero_diam) 00368 00369 !> Bin_grid structure. 00370 type(bin_grid_t), intent(in) :: bin_grid 00371 !> NetCDF file ID, in data mode. 00372 integer, intent(in) :: ncid 00373 !> Dimid of the aero_diam dimension. 00374 integer, intent(out) :: dimid_aero_diam 00375 00376 integer :: status, varid_aero_diam 00377 integer :: dimid_aero_diam_edges, varid_aero_diam_edges, 00378 varid_aero_diam_widths 00379 real(kind=dp) :: aero_diam_centers(bin_grid%n_bin), 00380 aero_diam_edges(bin_grid%n_bin + 1) 00381 real(kind=dp) :: aero_diam_widths(bin_grid%n_bin) 00382 00383 status = nf90_inq_dimid(ncid, "aero_diam", dimid_aero_diam) 00384 if (status == NF90_NOERR) return 00385 if (status /= NF90_EBADDIM) call pmc_nc_check(status) 00386 00387 ! dimension not defined, so define now define it 00388 call pmc_nc_check(nf90_redef(ncid)) 00389 00390 call pmc_nc_check(nf90_def_dim(ncid, "aero_diam", & 00391 bin_grid%n_bin, dimid_aero_diam)) 00392 call pmc_nc_check(nf90_def_dim(ncid, "aero_diam_edges", & 00393 bin_grid%n_bin + 1, dimid_aero_diam_edges)) 00394 00395 call pmc_nc_check(nf90_enddef(ncid)) 00396 00397 aero_diam_centers = rad2diam(bin_grid%center_radius) 00398 aero_diam_widths = bin_grid%log_width 00399 aero_diam_edges = rad2diam(bin_grid%edge_radius) 00400 00401 call pmc_nc_write_real_1d(ncid, aero_diam_centers, & 00402 "aero_diam", (/ dimid_aero_diam /), unit="m", & 00403 long_name="aerosol diameter axis bin centers", & 00404 description="logarithmically spaced centers of diameter axis grid, " & 00405 // "so that aero_diam(i) / aero_diam_edges(i) = " & 00406 // "0.5 * aero_diam_edges(i+1) / aero_diam_edges(i)") 00407 call pmc_nc_write_real_1d(ncid, aero_diam_edges, & 00408 "aero_diam_edges", (/ dimid_aero_diam_edges /), unit="m", & 00409 long_name="aerosol diameter axis bin edges", & 00410 description="logarithmically spaced edges of diameter axis grid, " & 00411 // "with one more edge than center") 00412 call pmc_nc_write_real_1d(ncid, aero_diam_widths, & 00413 "aero_diam_widths", (/ dimid_aero_diam /), unit="m", & 00414 long_name="aerosol diameter axis bin widths", & 00415 description="base-e logarithmic widths of diameter axis grid, " & 00416 // "so that aero_diam_widths(i)" & 00417 // "= ln(aero_diam_edges(i+1) / aero_diam_edges(i)) and " & 00418 // "all bins have the same width") 00419 00420 end subroutine bin_grid_netcdf_dim_aero_diam 00421 00422 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00423 00424 !> Write full state. 00425 subroutine bin_grid_output_netcdf(bin_grid, ncid) 00426 00427 !> bin_grid to write. 00428 type(bin_grid_t), intent(in) :: bin_grid 00429 !> NetCDF file ID, in data mode. 00430 integer, intent(in) :: ncid 00431 00432 integer :: dimid_aero_diam 00433 00434 !> \page output_format_bin_grid Output File Format: Bin Grid Data 00435 !! 00436 !! The bin grid data NetCDF dimensions are: 00437 !! - \b aero_diam: number of bins (grid cells) on the diameter axis 00438 !! - \b aero_diam_edges: number of bin edges (grid cell edges) on 00439 !! the diameter axis --- always equal to <tt>aero_diam + 1</tt> 00440 !! 00441 !! The bin grid data NetCDF variables are: 00442 !! - \b aero_diam (unit m, dim \c aero_diam): aerosol diameter axis 00443 !! bin centers --- centered on a logarithmic scale from the edges, so 00444 !! that <tt>aero_diam(i) / aero_diam_edges(i) = 00445 !! sqrt(aero_diam_edges(i+1) / aero_diam_edges(i))</tt> 00446 !! - \b aero_diam_edges (unit m, dim \c aero_diam_edges): aersol 00447 !! diameter axis bin edges (there is one more edge than center) 00448 !! - \b aero_diam_widths (dimensionless, dim \c aero_diam): 00449 !! the base-e logarithmic bin widths --- <tt>aero_diam_widths(i) 00450 !! = ln(aero_diam_edges(i+1) / aero_diam_edges(i))</tt>, so 00451 !! all bins have the same width 00452 !! 00453 !! See also: 00454 !! - \ref input_format_bin_grid --- the corresponding input format 00455 00456 call bin_grid_netcdf_dim_aero_diam(bin_grid, ncid, & 00457 dimid_aero_diam) 00458 00459 ! no need to write any more data as it's all contained in the 00460 ! dimension and associated variables 00461 00462 end subroutine bin_grid_output_netcdf 00463 00464 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00465 00466 !> Read full state. 00467 subroutine bin_grid_input_netcdf(bin_grid, ncid) 00468 00469 !> bin_grid to read. 00470 type(bin_grid_t), intent(inout) :: bin_grid 00471 !> NetCDF file ID, in data mode. 00472 integer, intent(in) :: ncid 00473 00474 integer :: dimid_aero_diam, n_bin 00475 character(len=1000) :: name 00476 real(kind=dp), allocatable :: aero_diam_centers(:) 00477 real(kind=dp), allocatable :: aero_diam_edges(:) 00478 real(kind=dp), allocatable :: aero_diam_widths(:) 00479 00480 call pmc_nc_check(nf90_inq_dimid(ncid, "aero_diam", dimid_aero_diam)) 00481 call pmc_nc_check(nf90_Inquire_Dimension(ncid, dimid_aero_diam, name, & 00482 n_bin)) 00483 00484 call bin_grid_deallocate(bin_grid) 00485 call bin_grid_allocate_size(bin_grid, n_bin) 00486 00487 allocate(aero_diam_centers(n_bin)) 00488 allocate(aero_diam_edges(n_bin + 1)) 00489 allocate(aero_diam_widths(n_bin)) 00490 00491 call pmc_nc_read_real_1d(ncid, aero_diam_centers, & 00492 "aero_diam") 00493 call pmc_nc_read_real_1d(ncid, aero_diam_edges, & 00494 "aero_diam_edges") 00495 call pmc_nc_read_real_1d(ncid, aero_diam_widths, & 00496 "aero_diam_widths") 00497 00498 bin_grid%center_radius = diam2rad(aero_diam_centers) 00499 bin_grid%edge_radius = diam2rad(aero_diam_edges) 00500 bin_grid%log_width = aero_diam_widths(1) 00501 00502 deallocate(aero_diam_centers) 00503 deallocate(aero_diam_edges) 00504 deallocate(aero_diam_widths) 00505 00506 end subroutine bin_grid_input_netcdf 00507 00508 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00509 00510 end module pmc_bin_grid