PartMC  2.2.1
aero_dist.F90
Go to the documentation of this file.
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_dist module.
00007 
00008 !> The aero_dist_t structure and associated subroutines.
00009 !!
00010 !! The initial size distributions are computed as number densities, so
00011 !! they can be used for both sectional and particle-resolved
00012 !! simulations. The routine dist_to_n() converts a number concentration
00013 !! distribution to an actual number of particles ready for a
00014 !! particle-resolved simulation.
00015 !!
00016 !! Initial distributions should be normalized so that <tt>sum(n_den) =
00017 !! 1/log_width</tt>.
00018 module pmc_aero_dist
00019 
00020   use pmc_bin_grid
00021   use pmc_util
00022   use pmc_constants
00023   use pmc_spec_file
00024   use pmc_aero_data
00025   use pmc_aero_mode
00026   use pmc_mpi
00027   use pmc_rand
00028 #ifdef PMC_USE_MPI
00029   use mpi
00030 #endif
00031 
00032   !> A complete aerosol distribution, consisting of several modes.
00033   type aero_dist_t
00034      !> Number of modes.
00035      integer :: n_mode
00036      !> Internally mixed modes [length \c n_mode].
00037      type(aero_mode_t), pointer :: mode(:)
00038   end type aero_dist_t
00039 
00040 contains
00041 
00042 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00043 
00044   !> Allocates an aero_dist.
00045   subroutine aero_dist_allocate(aero_dist)
00046 
00047     !> Aerosol distribution.
00048     type(aero_dist_t), intent(out) :: aero_dist
00049 
00050     integer :: i
00051 
00052     aero_dist%n_mode = 0
00053     allocate(aero_dist%mode(0))
00054 
00055   end subroutine aero_dist_allocate
00056 
00057 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00058 
00059   !> Allocates an aero_dist of the given size.
00060   subroutine aero_dist_allocate_size(aero_dist, n_mode, n_spec)
00061 
00062     !> Aerosol distribution.
00063     type(aero_dist_t), intent(out) :: aero_dist
00064     !> Number of modes.
00065     integer, intent(in) :: n_mode
00066     !> Number of species.
00067     integer, intent(in) :: n_spec
00068 
00069     integer :: i
00070 
00071     aero_dist%n_mode = n_mode
00072     allocate(aero_dist%mode(n_mode))
00073     do i = 1,n_mode
00074        call aero_mode_allocate_size(aero_dist%mode(i), n_spec)
00075     end do
00076 
00077   end subroutine aero_dist_allocate_size
00078 
00079 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00080 
00081   !> Free all storage.
00082   subroutine aero_dist_deallocate(aero_dist)
00083 
00084     !> Aerosol distribution.
00085     type(aero_dist_t), intent(inout) :: aero_dist
00086 
00087     integer :: i
00088 
00089     do i = 1,aero_dist%n_mode
00090        call aero_mode_deallocate(aero_dist%mode(i))
00091     end do
00092     deallocate(aero_dist%mode)
00093 
00094   end subroutine aero_dist_deallocate
00095 
00096 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00097 
00098   !> Copy an aero_dist.
00099   subroutine aero_dist_copy(aero_dist_from, aero_dist_to)
00100 
00101     !> Aero_dist original.
00102     type(aero_dist_t), intent(in) :: aero_dist_from
00103     !> Aero_dist copy.
00104     type(aero_dist_t), intent(inout) :: aero_dist_to
00105 
00106     integer :: n_spec, i
00107 
00108     if (aero_dist_from%n_mode > 0) then
00109        n_spec = size(aero_dist_from%mode(1)%vol_frac)
00110     else
00111        n_spec = 0
00112     end if
00113     call aero_dist_deallocate(aero_dist_to)
00114     call aero_dist_allocate_size(aero_dist_to, aero_dist_from%n_mode, n_spec)
00115     do i = 1,aero_dist_from%n_mode
00116        call aero_mode_copy(aero_dist_from%mode(i), aero_dist_to%mode(i))
00117     end do
00118 
00119   end subroutine aero_dist_copy
00120 
00121 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00122 
00123   !> Returns the total number concentration of a distribution. (#/m^3)
00124   real(kind=dp) function aero_dist_total_num_conc(aero_dist)
00125 
00126     !> Aerosol distribution.
00127     type(aero_dist_t), intent(in) :: aero_dist
00128 
00129     integer :: i_mode
00130 
00131     aero_dist_total_num_conc = 0d0
00132     do i_mode = 1,aero_dist%n_mode
00133        aero_dist_total_num_conc = aero_dist_total_num_conc &
00134             + aero_mode_total_num_conc(aero_dist%mode(i_mode))
00135     end do
00136 
00137   end function aero_dist_total_num_conc
00138 
00139 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00140 
00141   !> Returns the total number of particles of a distribution.
00142   real(kind=dp) function aero_dist_number(aero_dist, aero_weight)
00143 
00144     !> Aerosol distribution.
00145     type(aero_dist_t), intent(in) :: aero_dist
00146     !> Aerosol weight.
00147     type(aero_weight_t), intent(in) :: aero_weight
00148 
00149     integer :: i_mode
00150 
00151     aero_dist_number = 0d0
00152     do i_mode = 1,aero_dist%n_mode
00153        aero_dist_number = aero_dist_number &
00154             + aero_mode_number(aero_dist%mode(i_mode), aero_weight)
00155     end do
00156 
00157   end function aero_dist_number
00158 
00159 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00160 
00161   !> Return the binned number concentration for an aero_dist.
00162   subroutine aero_dist_num_conc(aero_dist, bin_grid, aero_data, &
00163        num_conc)
00164 
00165     !> Aero dist for which to compute number concentration.
00166     type(aero_dist_t), intent(in) :: aero_dist
00167     !> Bin grid.
00168     type(bin_grid_t), intent(in) :: bin_grid
00169     !> Aerosol data.
00170     type(aero_data_t), intent(in) :: aero_data
00171     !> Number concentration (#(ln(r))d(ln(r))).
00172     real(kind=dp), intent(out) :: num_conc(bin_grid%n_bin)
00173 
00174     integer :: i_mode
00175     real(kind=dp) :: mode_num_conc(size(num_conc, 1))
00176 
00177     num_conc = 0d0
00178     do i_mode = 1,aero_dist%n_mode
00179        call aero_mode_num_conc(aero_dist%mode(i_mode), bin_grid, &
00180             aero_data, mode_num_conc)
00181        num_conc = num_conc + mode_num_conc
00182     end do
00183 
00184   end subroutine aero_dist_num_conc
00185 
00186 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00187 
00188   !> Return the binned per-species volume concentration for an
00189   !> aero_dist.
00190   subroutine aero_dist_vol_conc(aero_dist, bin_grid, aero_data, &
00191        vol_conc)
00192 
00193     !> Aero dist for which to compute volume concentration.
00194     type(aero_dist_t), intent(in) :: aero_dist
00195     !> Bin grid.
00196     type(bin_grid_t), intent(in) :: bin_grid
00197     !> Aerosol data.
00198     type(aero_data_t), intent(in) :: aero_data
00199     !> Volume concentration (V(ln(r))d(ln(r))).
00200     real(kind=dp), intent(out) :: vol_conc(bin_grid%n_bin, aero_data%n_spec)
00201 
00202     integer :: i_mode
00203     real(kind=dp) :: mode_vol_conc(size(vol_conc, 1), size(vol_conc, 2))
00204 
00205     vol_conc = 0d0
00206     do i_mode = 1,aero_dist%n_mode
00207        call aero_mode_vol_conc(aero_dist%mode(i_mode), bin_grid, &
00208             aero_data, mode_vol_conc)
00209        vol_conc = vol_conc + mode_vol_conc
00210     end do
00211 
00212   end subroutine aero_dist_vol_conc
00213 
00214 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00215 
00216   !> Whether any of the modes are of the given type.
00217   elemental logical function aero_dist_contains_aero_mode_type(aero_dist, &
00218        aero_mode_type)
00219 
00220     !> Aerosol distribution.
00221     type(aero_dist_t), intent(in) :: aero_dist
00222     !> Aerosol mode type to test for.
00223     integer, intent(in) :: aero_mode_type
00224 
00225     aero_dist_contains_aero_mode_type &
00226          = any(aero_dist%mode%type == aero_mode_type)
00227 
00228   end function aero_dist_contains_aero_mode_type
00229 
00230 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00231 
00232   !> Determine the current aero_dist and rate by interpolating at the
00233   !> current time with the lists of aero_dists and rates.
00234   subroutine aero_dist_interp_1d(aero_dist_list, time_list, &
00235          rate_list, time, aero_dist, rate)
00236 
00237     !> Gas states.
00238     type(aero_dist_t), intent(in) :: aero_dist_list(:)
00239     !> Times (s).
00240     real(kind=dp), intent(in) :: time_list(size(aero_dist_list))
00241     !> Rates (s^{-1}).
00242     real(kind=dp), intent(in) :: rate_list(size(aero_dist_list))
00243     !> Current time (s).
00244     real(kind=dp), intent(in) :: time
00245     !> Current gas state.
00246     type(aero_dist_t), intent(inout) :: aero_dist
00247     !> Current rate (s^{-1}).
00248     real(kind=dp), intent(out) :: rate
00249 
00250     integer :: n, p, n_bin, n_spec, i, i_new
00251     real(kind=dp) :: y, alpha
00252 
00253     n = size(aero_dist_list)
00254     p = find_1d(n, time_list, time)
00255     if (p == 0) then
00256        ! before the start, just use the first state and rate
00257        call aero_dist_copy(aero_dist_list(1), aero_dist)
00258        rate = rate_list(1)
00259     elseif (p == n) then
00260        ! after the end, just use the last state and rate
00261        call aero_dist_copy(aero_dist_list(n), aero_dist)
00262        rate = rate_list(n)
00263     else
00264        ! in the middle, use the previous dist
00265        call aero_dist_copy(aero_dist_list(p), aero_dist)
00266        rate = rate_list(p)
00267     end if
00268 
00269   end subroutine aero_dist_interp_1d
00270 
00271 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00272 
00273   !> Read continuous aerosol distribution composed of several modes.
00274   subroutine spec_file_read_aero_dist(file, aero_data, aero_dist)
00275 
00276     !> Spec file to read data from.
00277     type(spec_file_t), intent(inout) :: file
00278     !> Aero_data data.
00279     type(aero_data_t), intent(inout) :: aero_data
00280     !> Aerosol dist.
00281     type(aero_dist_t), intent(inout) :: aero_dist
00282 
00283     type(aero_mode_t), pointer :: new_aero_mode_list(:)
00284     type(aero_mode_t) :: aero_mode
00285     integer :: i, j
00286     logical :: eof
00287 
00288     ! note that the <p> is needed below to force correct paragraph
00289     ! breaking by doxygen
00290     
00291     !> \page input_format_aero_dist Input File Format: Aerosol Distribution
00292     !!
00293     !! <p>An aerosol distribution file consists of zero or more modes,
00294     !! each in the format described by \subpage input_format_aero_mode
00295     !!
00296     !! See also:
00297     !!   - \ref spec_file_format --- the input file text format
00298     !!   - \ref input_format_aero_mode --- the format for each mode
00299     !!     of an aerosol distribution
00300 
00301     call aero_dist_deallocate(aero_dist)
00302     call aero_dist_allocate(aero_dist)
00303     call aero_mode_allocate(aero_mode)
00304     call spec_file_read_aero_mode(file, aero_data, aero_mode, eof)
00305     do while (.not. eof)
00306        aero_dist%n_mode = aero_dist%n_mode + 1
00307        allocate(new_aero_mode_list(aero_dist%n_mode))
00308        do i = 1,aero_dist%n_mode
00309           call aero_mode_allocate(new_aero_mode_list(i))
00310        end do
00311        call aero_mode_copy(aero_mode, &
00312             new_aero_mode_list(aero_dist%n_mode))
00313        do i = 1,(aero_dist%n_mode - 1)
00314           call aero_mode_copy(aero_dist%mode(i), &
00315                new_aero_mode_list(i))
00316           call aero_mode_deallocate(aero_dist%mode(i))
00317        end do
00318        deallocate(aero_dist%mode)
00319        aero_dist%mode => new_aero_mode_list
00320        nullify(new_aero_mode_list)
00321        call spec_file_read_aero_mode(file, aero_data, aero_mode, eof)
00322     end do
00323     call aero_mode_deallocate(aero_mode)
00324 
00325   end subroutine spec_file_read_aero_dist
00326 
00327 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00328 
00329   !> Read an array of aero_dists with associated times and rates from
00330   !> the given file.
00331   subroutine spec_file_read_aero_dists_times_rates(file, aero_data, &
00332        times, rates, aero_dists)
00333 
00334     !> Spec file to read data from.
00335     type(spec_file_t), intent(inout) :: file
00336     !> Aero data.
00337     type(aero_data_t), intent(inout) :: aero_data
00338     !> Times (s).
00339     real(kind=dp), pointer :: times(:)
00340     !> Rates (s^{-1}).
00341     real(kind=dp), pointer :: rates(:)
00342     !> Aero dists.
00343     type(aero_dist_t), pointer :: aero_dists(:)
00344 
00345     type(spec_line_t) :: aero_dist_line
00346     type(spec_file_t) :: aero_dist_file
00347     integer :: n_time, i_time
00348     character(len=SPEC_LINE_MAX_VAR_LEN), pointer :: names(:)
00349     real(kind=dp), pointer :: data(:,:)
00350 
00351     !> \page input_format_aero_dist_profile Input File Format: Aerosol Distribution Profile
00352     !!
00353     !! An aerosol distribution profile input file must consist of
00354     !! three lines:
00355     !! - the first line must begin with \c time and should be followed
00356     !!   by \f$N\f$ space-separated real scalars, giving the times (in
00357     !!   s after the start of the simulation) of the aerosol
00358     !!   distrbution set points --- the times must be in increasing
00359     !!   order
00360     !! - the second line must begin with \c rate and should be
00361     !!   followed by \f$N\f$ space-separated real scalars, giving the
00362     !!   scalings at the corresponding times
00363     !! - the third line must begin with \c dist and should be followed
00364     !!   by \f$N\f$ space-separated filenames, each specifying an
00365     !!   aerosol distribution in the format \ref input_format_aero_dist
00366     !!   at the corresponding time
00367     !!
00368     !! The units of the \c rate line depend on the type of aerosol
00369     !! distribution profile:
00370     !! - Emissions aerosol profiles have rates with units m/s ---
00371     !!   the aerosol distribution number concentrations are multiplied
00372     !!   by the rate to give an emission rate with unit #/(m^2 s)
00373     !!   which is then divided by the current mixing layer height
00374     !!   to give a per-volume emission rate.
00375     !! - Background aerosol profiles have rates with units \f$\rm
00376     !!   s^{-1}\f$, which is the dilution rate between the background
00377     !!   and the simulated air parcel. That is, if the simulated
00378     !!   number concentration is \f$N\f$ and the background number
00379     !!   concentration is \f$N_{\rm back}\f$, then dilution is modeled
00380     !!   as \f$\dot{N} = r N_{\rm back} - r N\f$, where \f$r\f$ is the
00381     !!   rate.
00382     !!
00383     !! Between the specified times the aerosol profile is interpolated
00384     !! step-wise and kept constant at its last value. That is, if the
00385     !! times are \f$t_i\f$, the rates are \f$r_i\f$, and the aerosol
00386     !! distributions are \f$a_i\f$ (all with \f$i = 1,\ldots,n\f$),
00387     !! then between times \f$t_i\f$ and \f$t_{i+1}\f$ the aerosol
00388     !! state is constant at \f$r_i a_i\f$. Before time \f$t_1\f$ the
00389     !! aerosol state is \f$r_1 a_1\f$, while after time \f$t_n\f$ it
00390     !! is \f$r_n a_n\f$.
00391     !!
00392     !! Example: an emissions aerosol profile could be:
00393     !! <pre>
00394     !! time  0          600        1800       # time (in s) after simulation start
00395     !! rate  1          0.5        1          # scaling factor in m/s
00396     !! dist  dist1.dat  dist2.dat  dist3.dat  # aerosol distribution files
00397     !! </pre>
00398     !! Here the emissions between 0&nbsp;min and 10&nbsp;min are given
00399     !! by <tt>dist1.dat</tt> (with the number concentration
00400     !! interpreted as having units 1/(m^2 s)), the emissions between
00401     !! 10&nbsp;min and 30&nbsp;min are given by <tt>dist2.dat</tt>
00402     !! (scaled by 0.5), while the emissions after 30&nbsp;min are
00403     !! given by <tt>dist3.dat</tt>.
00404     !!
00405     !! See also:
00406     !!   - \ref spec_file_format --- the input file text format
00407     !!   - \ref input_format_aero_data --- the aerosol species list
00408     !!     and material data
00409     !!   - \ref input_format_aero_dist --- the format of the
00410     !!     instantaneous aerosol distribution files
00411 
00412     ! read the data from the file
00413     allocate(names(0))
00414     allocate(data(0,0))
00415     call spec_file_read_real_named_array(file, 2, names, data)
00416     call spec_line_allocate(aero_dist_line)
00417     call spec_file_read_line_no_eof(file, aero_dist_line)
00418     call spec_file_check_line_name(file, aero_dist_line, "dist")
00419     call spec_file_check_line_length(file, aero_dist_line, size(data, 2))
00420 
00421     ! check the data size
00422     if (trim(names(1)) /= 'time') then
00423        call die_msg(570205795, 'row 1 in ' // trim(file%name) &
00424             // ' must start with: time not: ' // trim(names(1)))
00425     end if
00426     if (trim(names(2)) /= 'rate') then
00427        call die_msg(221270915, 'row 2 in ' // trim(file%name) &
00428             // ' must start with: rate not: ' // trim(names(1)))
00429     end if
00430     n_time = size(data, 2)
00431     if (n_time < 1) then
00432        call die_msg(457229710, 'each line in ' // trim(file%name) &
00433             // ' must contain at least one data value')
00434     end if
00435 
00436     ! copy over the data
00437     do i_time = 1,size(aero_dists)
00438        call aero_dist_deallocate(aero_dists(i_time))
00439     end do
00440     deallocate(aero_dists)
00441     deallocate(times)
00442     deallocate(rates)
00443     allocate(aero_dists(n_time))
00444     allocate(times(n_time))
00445     allocate(rates(n_time))
00446     do i_time = 1,n_time
00447        call spec_file_open(aero_dist_line%data(i_time), aero_dist_file)
00448        call aero_dist_allocate(aero_dists(i_time))
00449        call spec_file_read_aero_dist(aero_dist_file, &
00450             aero_data, aero_dists(i_time))
00451        call spec_file_close(aero_dist_file)
00452        times(i_time) = data(1,i_time)
00453        rates(i_time) = data(2,i_time)
00454     end do
00455     deallocate(names)
00456     deallocate(data)
00457     call spec_line_deallocate(aero_dist_line)
00458 
00459   end subroutine spec_file_read_aero_dists_times_rates
00460 
00461 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00462 
00463   !> Determines the number of bytes required to pack the given value.
00464   integer function pmc_mpi_pack_size_aero_dist(val)
00465 
00466     !> Value to pack.
00467     type(aero_dist_t), intent(in) :: val
00468 
00469     integer :: i, total_size
00470 
00471     total_size = pmc_mpi_pack_size_integer(val%n_mode)
00472     do i = 1,size(val%mode)
00473        total_size = total_size + pmc_mpi_pack_size_aero_mode(val%mode(i))
00474     end do
00475     pmc_mpi_pack_size_aero_dist = total_size
00476 
00477   end function pmc_mpi_pack_size_aero_dist
00478 
00479 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00480 
00481   !> Packs the given value into the buffer, advancing position.
00482   subroutine pmc_mpi_pack_aero_dist(buffer, position, val)
00483 
00484     !> Memory buffer.
00485     character, intent(inout) :: buffer(:)
00486     !> Current buffer position.
00487     integer, intent(inout) :: position
00488     !> Value to pack.
00489     type(aero_dist_t), intent(in) :: val
00490 
00491 #ifdef PMC_USE_MPI
00492     integer :: prev_position, i
00493 
00494     prev_position = position
00495     call pmc_mpi_pack_integer(buffer, position, val%n_mode)
00496     do i = 1,size(val%mode)
00497        call pmc_mpi_pack_aero_mode(buffer, position, val%mode(i))
00498     end do
00499     call assert(440557910, &
00500          position - prev_position <= pmc_mpi_pack_size_aero_dist(val))
00501 #endif
00502 
00503   end subroutine pmc_mpi_pack_aero_dist
00504 
00505 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00506 
00507   !> Unpacks the given value from the buffer, advancing position.
00508   subroutine pmc_mpi_unpack_aero_dist(buffer, position, val)
00509 
00510     !> Memory buffer.
00511     character, intent(inout) :: buffer(:)
00512     !> Current buffer position.
00513     integer, intent(inout) :: position
00514     !> Value to pack.
00515     type(aero_dist_t), intent(inout) :: val
00516 
00517 #ifdef PMC_USE_MPI
00518     integer :: prev_position, i
00519 
00520     call aero_dist_deallocate(val)
00521     prev_position = position
00522     call pmc_mpi_unpack_integer(buffer, position, val%n_mode)
00523     allocate(val%mode(val%n_mode))
00524     do i = 1,size(val%mode)
00525        call aero_mode_allocate(val%mode(i))
00526        call pmc_mpi_unpack_aero_mode(buffer, position, val%mode(i))
00527     end do
00528     call assert(742535268, &
00529          position - prev_position <= pmc_mpi_pack_size_aero_dist(val))
00530 #endif
00531 
00532   end subroutine pmc_mpi_unpack_aero_dist
00533 
00534 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00535 
00536 end module pmc_aero_dist