ball_berry_gs.cpp

Go to the documentation of this file.

#include <cmath>                          // for std::min
#include <limits>                         // for std::numeric_limits
#include <stdexcept>                      // for std::range_error
#include "../framework/constants.h"       // for dr_boundary, eps_zero
#include "../framework/quadratic_root.h"  // for quadratic_root_plus
#include "water_and_air_properties.h"     // for saturation_vapor_pressure
#include "ball_berry_gs.h"
 
using calculation_constants::eps_zero;
using physical_constants::dr_boundary;
 
stomata_outputs ball_berry_gs(
    double assimilation,            // mol / m^2 / s
    double ambient_c,               // mol / mol
    double ambient_rh,              // Pa / Pa
    double bb_offset,               // mol / m^2 / s
    double bb_slope,                // dimensionless from [mol / m^2 / s] / [mol / m^2 / s]
    double gbw,                     // mol / m^2 / s
    double leaf_temperature,        // degrees C
    double ambient_air_temperature  // degrees C
)
{
    // If An < 0, set b1 = 0 to ensure that gsw = b0 in Equation (1) as defined
    // above
    if (assimilation < 0) {
        bb_slope = 0.0;  // mol / m^2 / s
    }
 
    // Determine Cs using Equation (2) as defined above
    const double Cs = ambient_c -
                      (dr_boundary / gbw) * assimilation;  // mol / mol.
 
    // Check for error conditions (Cs = 0 or Cs < 0)
    if (Cs < -eps_zero) {
        throw std::range_error("Thrown in ball_berry_gs: Cs is negative.");
    } else if (Cs <= eps_zero) {
        // Stomatal conductance becomes infinite as Cs approaches zero from the
        // right. In this case, there is no water vapor drawndown across the
        // stomata, so hs becomes 1.
        double const inf = std::numeric_limits<double>::infinity();
        return stomata_outputs{
            /* .cs = */ 0,    // micromol / mol
            /* .hs = */ 1,    // dimensionless
            /* .gsw = */ inf  // mol / m^2 / s
        };
    }
 
    // Calculate some variables that will be used in later equations
    const double acs = assimilation / Cs;  // mol / m^2 / s
 
    const double swvp_ratio =
        saturation_vapor_pressure(ambient_air_temperature) /
        saturation_vapor_pressure(leaf_temperature);  // dimensionless
 
    // Calculate hs using Equation (3) as defined above
    const double a = bb_slope * acs;                              // mol / m^2 / s
    const double b = bb_offset + gbw - a;                         // mol / m^2 / s
    const double c = -ambient_rh * gbw * swvp_ratio - bb_offset;  // mol / m^2 / s
 
    // If hs is calculated to be larger than 1, this indicates dew formation. We
    // do not handle this in BioCro at the moment, so just limit hs to 1.
    const double hs = std::min(1.0, quadratic_root_plus(a, b, c));  // dimensionless
 
    if (hs < 0) {
        throw std::range_error("Thrown in ball_berry_gs: hs is less than 0.");
    }
 
    // Calculate stomatal conductance using Equation (1) above
    double const gswmol = a * hs + bb_offset;  // mol / m^2 / s
 
    return stomata_outputs{
        /* .cs = */ Cs * 1e6,  // micromol / mol
        /* .hs = */ hs,        // dimensionless
        /* .gsw = */ gswmol    // mol / m^2 / s
    };
}