PartMC  2.2.0
coag_kernel_constant.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_coag_kernel_constant module.
00007 
00008 !> Constant coagulation kernel.
00009 module pmc_coag_kernel_constant
00010 
00011   use pmc_env_state
00012   use pmc_bin_grid
00013   use pmc_util
00014   use pmc_constants
00015   use pmc_aero_binned
00016   use pmc_aero_data
00017   use pmc_aero_dist
00018   use pmc_aero_data
00019   use pmc_aero_particle
00020 
00021   !> Coefficient for constant kernel.
00022   real(kind=dp), parameter :: beta_0 = 0.25d0 / (60d0 * 2d8)
00023 
00024 contains
00025   
00026 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00027 
00028   !> Constant coagulation kernel.
00029   subroutine kernel_constant(aero_particle_1, aero_particle_2, &
00030        aero_data, env_state, k)
00031 
00032     !> First particle.
00033     type(aero_particle_t), intent(in) :: aero_particle_1
00034     !> Second particle.
00035     type(aero_particle_t), intent(in) :: aero_particle_2
00036     !> Aerosol data.
00037     type(aero_data_t), intent(in) :: aero_data
00038     !> Environment state.
00039     type(env_state_t), intent(in) :: env_state
00040     !> Coagulation kernel.
00041     real(kind=dp), intent(out) :: k
00042 
00043     k = beta_0
00044     
00045   end subroutine kernel_constant
00046   
00047 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00048 
00049   !> Minimum and maximum values of the constant coagulation kernel.
00050   subroutine kernel_constant_minmax(v1, v2, aero_data, env_state, k_min, k_max)
00051 
00052     !> Volume of first particle.
00053     real(kind=dp), intent(in) :: v1
00054     !> Volume of second particle.
00055     real(kind=dp), intent(in) :: v2
00056     !> Aerosol data.
00057     type(aero_data_t), intent(in) :: aero_data
00058     !> Environment state.
00059     type(env_state_t), intent(in) :: env_state
00060     !> Coagulation kernel minimum value.
00061     real(kind=dp), intent(out) :: k_min
00062     !> Coagulation kernel maximum value.
00063     real(kind=dp), intent(out) :: k_max
00064 
00065     k_min = beta_0
00066     k_max = beta_0
00067     
00068   end subroutine kernel_constant_minmax
00069   
00070 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00071 
00072   !> Exact solution with a constant coagulation kernel and an
00073   !> exponential initial condition.
00074   !!
00075   !! Given input paramaters \f$R\f$ and \f$N_0\f$ we let the mean
00076   !! volume be \f$v_\mu = \frac{4\pi}{3} R^3\f$ and define the
00077   !! rescaled time \f$\tau = N_0 \beta_0 t\f$, where \f$\beta_0\f$ is
00078   !! the fixed constant kernel value. We also set the parameter
00079   !! \f$\lambda = 1\f$. Then the solution is
00080   !! \f[
00081   !!     n(D,t) \ {\rm d}\ln D
00082   !!     = \frac{\pi}{2} D^3 \frac{N_0}{v_\mu} \frac{4}{(\tau + 2)^2}
00083   !!     \exp\left(-\frac{v}{v_\mu} \frac{2}{\tau + 2}
00084   !!     \exp(-\lambda \tau) - \lambda \tau\right) {\rm d}\ln D
00085   !! \f]
00086   !! This thus has initial condition
00087   !! \f[
00088   !!     n(D,t) \ {\rm d}\ln D
00089   !!     = \frac{\pi}{2} D^3 \frac{N_0}{v_\mu}
00090   !!     \exp\left(-\frac{v}{v_\mu}\right) {\rm d}\ln D
00091   !! \f]
00092   subroutine soln_constant_exp(bin_grid, aero_data, time, num_conc, &
00093        radius_at_mean_vol, env_state, aero_binned)
00094 
00095     !> Bin grid.
00096     type(bin_grid_t), intent(in) :: bin_grid
00097     !> Aerosol data.
00098     type(aero_data_t), intent(in) :: aero_data
00099     !> Current time.
00100     real(kind=dp), intent(in) :: time
00101     !> Particle number concentration (#/m^3).
00102     real(kind=dp), intent(in) :: num_conc
00103     !> Mean init radius (m).
00104     real(kind=dp), intent(in) :: radius_at_mean_vol
00105     !> Environment state.
00106     type(env_state_t), intent(in) :: env_state
00107     !> Output state.
00108     type(aero_binned_t), intent(inout) :: aero_binned
00109     
00110     real(kind=dp) :: tau, T, rat_v, nn, b, sigma, mean_vol
00111     integer :: k
00112     
00113     real(kind=dp), parameter :: lambda = 1d0
00114 
00115     mean_vol = rad2vol(radius_at_mean_vol)
00116     if (time .eq. 0d0) then
00117        do k = 1,bin_grid%n_bin
00118           aero_binned%num_conc(k) = const%pi / 2d0 &
00119                * (2d0 * bin_grid%center_radius(k))**3 * num_conc / mean_vol &
00120                * exp(-(rad2vol(bin_grid%center_radius(k)) / mean_vol))
00121        end do
00122     else
00123        tau = num_conc * beta_0 * time
00124        do k = 1,bin_grid%n_bin
00125           rat_v = rad2vol(bin_grid%center_radius(k)) / mean_vol
00126           nn = 4d0 * num_conc / (mean_vol * ( tau + 2d0 ) ** 2d0) &
00127                * exp(-2d0*rat_v/(tau+2d0)*exp(-lambda*tau)-lambda*tau)
00128           aero_binned%num_conc(k) = const%pi / 2d0 &
00129                * (2d0 * bin_grid%center_radius(k))**3d0 * nn
00130        end do
00131     end if
00132     
00133     aero_binned%vol_conc = 0d0
00134     do k = 1,bin_grid%n_bin
00135        aero_binned%vol_conc(k,1) = rad2vol(bin_grid%center_radius(k)) &
00136             * aero_binned%num_conc(k)
00137     end do
00138     
00139   end subroutine soln_constant_exp
00140   
00141 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00142   
00143 end module pmc_coag_kernel_constant