fast_arbitrary_rate was too slow. Also made it non-constant

CMakeLists.txt

@@ -1,5 +1,5 @@
 cmake_minimum_required(VERSION "3.14")
-project(Photochem LANGUAGES Fortran C VERSION "0.3.14")
+project(Photochem LANGUAGES Fortran C VERSION "0.3.15")
 
 include(FortranCInterface)
 FortranCInterface_VERIFY()

photochem/cython/PhotochemVars.pyx

@@ -144,5 +144,13 @@ cdef class PhotochemVars:
     def __set__(self, int val):
       var_pxd.photochemvars_verbose_set(&self._ptr, &val)
 
+  property fast_arbitrary_rate:
+    def __get__(self):
+      cdef double val
+      var_pxd.photochemvars_fast_arbitrary_rate_get(&self._ptr, &val)
+      return val
+    def __set__(self, double val):
+      var_pxd.photochemvars_fast_arbitrary_rate_set(&self._ptr, &val)
+
 
 

photochem/cython/PhotochemVars_pxd.pxd

@@ -50,4 +50,7 @@ cdef extern void photochemvars_equilibrium_time_get(void *ptr, double *val)
 cdef extern void photochemvars_equilibrium_time_set(void *ptr, double *val)
 
 cdef extern void photochemvars_verbose_get(void *ptr, int *val)
-cdef extern void photochemvars_verbose_set(void *ptr, int *val)
+cdef extern void photochemvars_verbose_set(void *ptr, int *val)
+
+cdef extern void photochemvars_fast_arbitrary_rate_get(void *ptr, double *val)
+cdef extern void photochemvars_fast_arbitrary_rate_set(void *ptr, double *val)

photochem/fortran/PhotochemVars_wrapper.f90

@@ -267,4 +267,20 @@ subroutine photochemvars_verbose_set(ptr, val) bind(c)
     call c_f_pointer(ptr, var)
     var%verbose = val
   end subroutine
+
+  subroutine photochemvars_fast_arbitrary_rate_get(ptr, val) bind(c)
+    type(c_ptr), intent(in) :: ptr
+    real(c_double), intent(out) :: val
+    type(PhotochemVars), pointer :: var
+    call c_f_pointer(ptr, var)
+    val = var%fast_arbitrary_rate
+  end subroutine
+
+  subroutine photochemvars_fast_arbitrary_rate_set(ptr, val) bind(c)
+    type(c_ptr), intent(in) :: ptr
+    real(c_double), intent(in) :: val
+    type(PhotochemVars), pointer :: var
+    call c_f_pointer(ptr, var)
+    var%fast_arbitrary_rate = val
+  end subroutine
 

src/atmosphere/photochem_atmosphere_rhs.f90

@@ -18,7 +18,7 @@ subroutine dochem(self, neqs, nsp, np, nsl, nq, nz, trop_ind, nrT, usol, density
     use photochem_enum, only: CondensingParticle
     use photochem_common, only: chempl, chempl_sl
     use photochem_eqns, only: damp_condensation_rate
-    use photochem_const, only: N_avo, pi, small_real, T_crit_H2O, fast_arbitrary_rate
+    use photochem_const, only: N_avo, pi, small_real, T_crit_H2O
 
     class(Atmosphere), target, intent(in) :: self
     integer, intent(in) :: neqs, nsp, np, nsl, nq, nz, trop_ind, nrT
@@ -84,7 +84,7 @@ subroutine dochem(self, neqs, nsp, np, nsl, nq, nz, trop_ind, nrT, usol, density
     if (dat%fix_water_in_trop) then
       do j = 1,trop_ind
         k = dat%LH2O + (j - 1) * dat%nq
-        rhs(k) = rhs(k) + fast_arbitrary_rate*(H2O_sat_mix(j)*H2O_rh(j) - usol(dat%LH2O,j))
+        rhs(k) = rhs(k) + var%fast_arbitrary_rate*(H2O_sat_mix(j)*H2O_rh(j) - usol(dat%LH2O,j))
       enddo
     endif
     ! H2O condensation

src/atmosphere/photochem_atmosphere_utils.f90

@@ -100,7 +100,6 @@ module function atom_conservation(self, atom, err) result(con)
     use photochem_enum, only: VelocityDistributedFluxBC
     use photochem_eqns, only: damp_condensation_rate
     use photochem_types, only: AtomConservation
-    use photochem_const, only: fast_arbitrary_rate
     class(Atmosphere), target, intent(inout) :: self
     character(len=*), intent(in) :: atom
     character(:), allocatable, intent(out) :: err
@@ -184,7 +183,7 @@ module function atom_conservation(self, atom, err) result(con)
 
     if (dat%fix_water_in_trop) then
       do j = 1,var%trop_ind
-        con_evap_rate = fast_arbitrary_rate*(wrk%H2O_sat_mix(j)*wrk%H2O_rh(j) - wrk%usol(dat%LH2O,j)) &
+        con_evap_rate = var%fast_arbitrary_rate*(wrk%H2O_sat_mix(j)*wrk%H2O_rh(j) - wrk%usol(dat%LH2O,j)) &
                         *wrk%density(j)*var%dz(j)*dat%species_composition(kk,dat%LH2O)
         if (con_evap_rate > 0.0_dp) then
           con%in_other = con%in_other + con_evap_rate

src/evoatmosphere/photochem_evoatmosphere_rhs.f90

@@ -23,7 +23,7 @@ subroutine dochem(self, usol, rx_rates, &
     use photochem_enum, only: CondensingParticle
     use photochem_common, only: chempl, chempl_sl
     use photochem_eqns, only: damp_condensation_rate
-    use photochem_const, only: N_avo, pi, small_real, T_crit_H2O, fast_arbitrary_rate
+    use photochem_const, only: N_avo, pi, small_real, T_crit_H2O
 
     class(EvoAtmosphere), target, intent(in) :: self
     real(dp), intent(in) :: usol(:,:)
@@ -89,7 +89,7 @@ subroutine dochem(self, usol, rx_rates, &
     if (dat%fix_water_in_trop) then
       do j = 1,var%trop_ind
         k = dat%LH2O + (j - 1) * dat%nq
-        rhs(k) = rhs(k) + fast_arbitrary_rate*(density(j)*H2O_sat_mix(j)*H2O_rh(j) - usol(dat%LH2O,j))
+        rhs(k) = rhs(k) + var%fast_arbitrary_rate*(density(j)*H2O_sat_mix(j)*H2O_rh(j) - usol(dat%LH2O,j))
       enddo
     endif
     if (dat%water_cond) then

src/photochem_const.f90

@@ -22,5 +22,4 @@ module photochem_const
   real(dp), parameter :: small_real = tiny(1.0_dp)**0.25_dp
   real(dp), parameter :: ln_small_real = log(small_real)
 
-  real(dp), parameter :: fast_arbitrary_rate = 1.0e-5_dp ! arbitrary rate that is fast (1/s)
 end module

src/photochem_types.f90

@@ -441,6 +441,8 @@ subroutine time_dependent_rate_fcn(tn, nz, rate)
     integer(c_int) :: max_order = 5
     real(dp) :: epsj = 1.d-9 ! perturbation for jacobian calculation
     integer :: verbose = 1 ! 0 == no printing. 1 == some printing. 2 == bunch of printing.
+    !> arbitrary rate that is fast (1/s). Used for keeping H2O at saturation in troposphere
+    real(dp) :: fast_arbitrary_rate = 1.0e-2_dp 
   end type
 
   type :: SundialsData