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