function that sets T and edd to input P, T, and edd profile

photochem/cython/Atmosphere.pyx

@@ -450,6 +450,38 @@ cdef class Atmosphere:
     if len(err.strip()) > 0:
       raise PhotoException(err.decode("utf-8").strip())
 
+  def set_press_temp_edd(self, ndarray[double, ndim=1] P, ndarray[double, ndim=1] T, ndarray[double, ndim=1] edd, trop_p = None):
+    """Given an input P, T, and edd, the code will find the temperature and eddy diffusion profile
+    on the current altitude-grid that matches the inputs.
+    
+    Parameters
+    ----------
+    P : ndarray[double,ndim=1]
+        Pressure (dynes/cm^2)
+    T : ndarray[double,ndim=1]
+        Temperature (K)
+    edd : ndarray[double,ndim=1]
+        Eddy diffusion (cm^2/s)
+    trop_p : float, optional
+        Tropopause pressure (dynes/cm^2). Only necessary if rainout == True, 
+        or fix_water_in_trop == True.
+    """
+    cdef char err[ERR_LEN+1]
+    cdef int P_dim1 = P.size
+    cdef int T_dim1 = T.size
+    cdef int edd_dim1 = edd.size
+
+    cdef double trop_p_ = 0.0
+    cdef bool trop_p_present = False
+    if trop_p != None:
+      trop_p_present = True
+      trop_p_ = trop_p
+
+    a_pxd.atmosphere_set_press_temp_edd_wrapper(&self._ptr, &P_dim1, <double *>P.data, &T_dim1, <double *>T.data, &edd_dim1, <double *>edd.data,
+                                               &trop_p_, &trop_p_present, err)
+    if len(err.strip()) > 0:
+      raise PhotoException(err.decode("utf-8").strip())
+
   def set_rate_fcn(self, str species, object fcn):
     """Sets a function describing a custom rate for a species.
     This could be useful for modeling external processes not in the
@@ -487,7 +519,15 @@ cdef class Atmosphere:
        raise PhotoException(err.decode("utf-8").strip())
 
   def update_vertical_grid(self, double TOA_pressure):
-    """docstring
+    """Re-does the vertical grid so that the pressure at the top of the
+    atmosphere is at `TOA_pressure`. If the TOA needs to be raised above the current
+    TOA, then the function constantly extrapolates mixing ratios, temperature,
+    eddy diffusion, and particle radii.
+
+    Parameters
+    ----------
+    TOA_pressure : float
+        New top of the atmosphere pressure in dynes/cm^2
     """
     cdef char err[ERR_LEN+1]
 

photochem/cython/Atmosphere_pxd.pxd

@@ -41,6 +41,9 @@ cdef extern void atmosphere_evolve_wrapper(void *ptr, char *filename,
                 double *tstart, int *nq, int *nz, double *usol, 
                 int *nt, double *t_eval, bool *overwrite, bool *success, char *err)
 
+cdef extern void atmosphere_set_press_temp_edd_wrapper(void *ptr, int *P_dim1, double *P, int *T_dim1, double *T, int *edd_dim1, double *edd,
+                                                      double *trop_p, bool *trop_p_present, char *err)
+
 cdef extern void atmosphere_set_temperature_wrapper(void *ptr, int *nz, double *temperature, 
                                                     double *trop_alt, bool *trop_alt_present, char *err)
 

photochem/cython/PhotochemVars.pyx

@@ -53,6 +53,20 @@ cdef class PhotochemVars:
       var_pxd.photochemvars_usol_init_get(&self._ptr, &dim1, &dim2, <double *>arr.data)
       return arr
 
+  property trop_alt:
+    "double. Tropopause altitude."
+    def __get__(self):
+      cdef double val
+      var_pxd.photochemvars_trop_alt_get(&self._ptr, &val)
+      return val
+
+  property trop_ind:
+    "int. Tropopause index."
+    def __get__(self):
+      cdef int val
+      var_pxd.photochemvars_trop_ind_get(&self._ptr, &val)
+      return val
+
   property relative_humidity:
     "double. Relative humidity of H2O in the troposphere."
     def __get__(self):

photochem/cython/PhotochemVars_pxd.pxd

@@ -20,6 +20,10 @@ cdef extern void photochemvars_at_photo_equilibrium_get(void *ptr, bool *at_phot
 cdef extern void photochemvars_usol_init_get_size(void *ptr, int *dim1, int *dim2)
 cdef extern void photochemvars_usol_init_get(void *ptr, int *dim1, int *dim2, double *usol_init)
 
+cdef extern void photochemvars_trop_alt_get(void *ptr, double *val)
+
+cdef extern void photochemvars_trop_ind_get(void *ptr, int *val)
+
 cdef extern void photochemvars_relative_humidity_get(void *ptr, double *val)
 cdef extern void photochemvars_relative_humidity_set(void *ptr, double *val)
 

photochem/fortran/Atmosphere_wrapper.f90

@@ -441,6 +441,36 @@ subroutine atmosphere_set_temperature_wrapper(ptr, nz, temperature, &
 
   end subroutine
 
+  subroutine atmosphere_set_press_temp_edd_wrapper(ptr, P_dim1, P, T_dim1, T, edd_dim1, edd, &
+                                                  trop_p, trop_p_present, err) bind(c)
+    type(c_ptr), intent(in) :: ptr
+    integer(c_int), intent(in) :: P_dim1
+    real(c_double), intent(in) :: P(P_dim1)
+    integer(c_int), intent(in) :: T_dim1
+    real(c_double), intent(in) :: T(T_dim1)
+    integer(c_int), intent(in) :: edd_dim1
+    real(c_double), intent(in) :: edd(edd_dim1)
+    real(c_double), intent(in) :: trop_p
+    logical(c_bool), intent(in) :: trop_p_present
+    character(kind=c_char), intent(out) :: err(err_len+1)
+
+    character(:), allocatable :: err_f
+    type(Atmosphere), pointer :: pc
+
+    call c_f_pointer(ptr, pc)
+
+    if (trop_p_present) then
+      call pc%set_press_temp_edd(P, T, edd, trop_p, err_f)
+    else
+      call pc%set_press_temp_edd(P, T, edd, err=err_f)
+    endif
+    err(1) = c_null_char
+    if (allocated(err_f)) then
+      call copy_string_ftoc(err_f, err)
+    endif
+
+  end subroutine
+
   subroutine atmosphere_set_rate_fcn_wrapper(ptr, species_c, fcn_c, err) bind(c)
     use photochem_types, only: time_dependent_rate_fcn
     type(c_ptr), intent(in) :: ptr

photochem/fortran/PhotochemVars_wrapper.f90

@@ -71,6 +71,22 @@ subroutine photochemvars_usol_init_get(ptr, dim1, dim2, usol_init) bind(c)
     usol_init = var%usol_init
   end subroutine
 
+  subroutine photochemvars_trop_alt_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%trop_alt
+  end subroutine
+
+  subroutine photochemvars_trop_ind_get(ptr, val) bind(c)
+    type(c_ptr), intent(in) :: ptr
+    integer(c_int), intent(out) :: val
+    type(PhotochemVars), pointer :: var
+    call c_f_pointer(ptr, var)
+    val = var%trop_ind
+  end subroutine
+
   subroutine photochemvars_relative_humidity_get(ptr, val) bind(c)
     type(c_ptr), intent(in) :: ptr
     real(c_double), intent(out) :: val

src/atmosphere/photochem_atmosphere.f90

@@ -40,6 +40,7 @@ module photochem_atmosphere
     procedure :: set_lower_bc
     procedure :: set_upper_bc
     procedure :: set_temperature
+    procedure :: set_press_temp_edd
     procedure :: set_rate_fcn
     procedure :: update_vertical_grid
   end type
@@ -271,6 +272,17 @@ module subroutine set_temperature(self, temperature, trop_alt, err)
       character(:), allocatable, intent(out) :: err
     end subroutine
 
+    !> Given an input P, T, and edd, the code will find the temperature and eddy diffusion profile
+    !> on the current altitude-grid that matches the inputs.
+    module subroutine set_press_temp_edd(self, P, T, edd, trop_p, err)
+      class(Atmosphere), target, intent(inout) :: self
+      real(dp), intent(in) :: P(:) !! Pressure (dynes/cm^2)
+      real(dp), intent(in) :: T(:) !! Temperature (K)
+      real(dp), intent(in) :: edd(:) !! Eddy diffusion (cm^2/s)
+      real(dp), optional, intent(in) :: trop_p !! Tropopause pressure (dynes/cm^2)
+      character(:), allocatable, intent(out) :: err
+    end subroutine
+
     !> Sets a function describing a custom rate for a species.
     !> This could be useful for modeling external processes not in the
     !> model.
@@ -282,9 +294,13 @@ module subroutine set_rate_fcn(self, species, fcn, err)
       character(:), allocatable, intent(inout) :: err
     end subroutine
 
+    !> Re-does the vertical grid so that the pressure at the top of the
+    !> atmosphere is at `TOA_pressure`. If the TOA needs to be raised above the current
+    !> TOA, then the function constantly extrapolates mixing ratios, temperature,
+    !> eddy diffusion, and particle radii.
     module subroutine update_vertical_grid(self, TOA_pressure, err)
       class(Atmosphere), target, intent(inout) :: self
-      real(dp), intent(in) :: TOA_pressure !! dynes/cm^2
+      real(dp), intent(in) :: TOA_pressure !! New top of atmosphere pressure (dynes/cm^2)
       character(:), allocatable, intent(out) :: err
     end subroutine