small fixes

- remove deprecated PModelDrivers constructor in test file
- allow users to change some PModelConstants by placing them in the kwargs
- add some more documentation to pmodel.jl

Author: ychnli

src/standalone/Vegetation/pmodel.jl

@@ -99,31 +99,56 @@ Base.broadcastable(m::PModelConstants) = tuple(m)
     PModelConstants(FT)
 
 Creates a `PModelConstants` object with default values for the P-model constants.
+See Stocker et al. (2020) Table A2 and references within for more information. 
 """
-function PModelConstants(FT)
-    return PModelConstants(
-        R = LP.gas_constant(LP.LandParameters(FT)),
-        Kc25 = FT(39.97),
-        Ko25 = FT(27480),
-        To = FT(298.15),
-        ΔHkc = FT(79430),
-        ΔHko = FT(36380),
-        Drel = FT(1.6),
-        ΔHΓstar = FT(37830),
-        Γstar25 = FT(4.332),
-        Ha_Vcmax = FT(71513),
-        Hd_Vcmax = FT(200000),
-        aS_Vcmax = FT(668.39),
-        bS_Vcmax = FT(1.07),
-        Ha_Jmax = FT(49884),
-        Hd_Jmax = FT(200000),
-        aS_Jmax = FT(659.70),
-        bS_Jmax = FT(0.75),
-        Mc = FT(0.0120107),
-        oi = FT(0.2095),
-        aRd = FT(0.1012),
-        bRd = FT(-0.0005),
-        fC3 = FT(0.015),
+function PModelConstants{FT}(;
+    R = LP.gas_constant(LP.LandParameters(FT)),
+    Kc25 = FT(39.97), # Pa (see note on line 38 above)
+    Ko25 = FT(27480), # Pa (see note on line 38 above)
+    To = FT(298.15),
+    ΔHkc = FT(79430),
+    ΔHko = FT(36380),
+    Drel = FT(1.6),
+    ΔHΓstar = FT(37830),
+    Γstar25 = FT(4.332), # Pa (see note on line 38 above)
+    Ha_Vcmax = FT(71513),
+    Hd_Vcmax = FT(200000),
+    aS_Vcmax = FT(668.39),
+    bS_Vcmax = FT(1.07),
+    Ha_Jmax = FT(49884),
+    Hd_Jmax = FT(200000),
+    aS_Jmax = FT(659.70),
+    bS_Jmax = FT(0.75),
+    Mc = FT(0.0120107),
+    oi = FT(0.2095),
+    aRd = FT(0.1012),
+    bRd = FT(-0.0005),
+    fC3 = FT(0.015),
+) where {FT <: AbstractFloat}
+
+    return PModelConstants{FT}(
+        R,
+        Kc25,
+        Ko25,
+        To,
+        ΔHkc,
+        ΔHko,
+        Drel,
+        ΔHΓstar,
+        Γstar25,
+        Ha_Vcmax,
+        Hd_Vcmax,
+        aS_Vcmax,
+        bS_Vcmax,
+        Ha_Jmax,
+        Hd_Jmax,
+        aS_Jmax,
+        bS_Jmax,
+        Mc,
+        oi,
+        aRd,
+        bRd,
+        fC3,
     )
 end
 
@@ -132,6 +157,13 @@ end
                 OPFT <: PModelParameters{FT},
                 OPCT <: PModelConstants{FT}
                 } <: AbstractPhotosynthesisModel{FT}
+
+    An implementation of the optimality photosynthesis model "P-model v1.0" of Stocker et al. (2020). 
+
+    Stocker, B. D., Wang, H., Smith, N. G., Harrison, S. P., Keenan, T. F., Sandoval, D., Davis, T., 
+        and Prentice, I. C.: P-model v1.0: an optimality-based light use efficiency model for simulating 
+        ecosystem gross primary production, Geosci. Model Dev., 13, 1545–1581, 
+        https://doi.org/10.5194/gmd-13-1545-2020, 2020.
 """
 struct PModel{FT, OPFT <: PModelParameters{FT}, OPCT <: PModelConstants{FT}} <:
        AbstractPhotosynthesisModel{FT}
@@ -231,6 +263,9 @@ function compute_full_pmodel_outputs(
     ) = constants
 
     # Compute intermediate values
+
+    # The quantum yield ϕ0 can either be taken as a constant or computed from the temperature
+    # dependent function. If ϕ0 is not NaN, it is used directly. 
     ϕ0 = isnan(ϕ0) ? intrinsic_quantum_yield(T_canopy, ϕc, ϕa0, ϕa1, ϕa2) : ϕ0
 
     Γstar = co2_compensation_p(T_canopy, To, P_air, R, ΔHΓstar, Γstar25)

test/standalone/Vegetation/test_pmodel.jl

@@ -53,36 +53,6 @@ function percent_difference(a, b)
     return abs(a - b) / abs(b) * 100
 end
 
-"""
-    create_pmodel_drivers(inputs::Dict{String, Any}, FT)
-
-Create P-Model driver variables from test input data.
-Converts input data from the test CSV files into the appropriate driver structure
-for the P-Model, handling unit conversions and optional soil moisture stress.
-"""
-function PModelDrivers(inputs::Dict{String, Any}, FT)
-    T_canopy = FT(inputs["tc"] + 273.15)  # Convert from Celsius to Kelvin
-    VPD = FT(inputs["vpd"])
-    ca = FT(inputs["co2"]) * FT(1e-6) * FT(101325.0)  # Convert ppm to Pa
-    P_air = FT(get(inputs, "patm", 101325.0))
-
-    # Calculate I_abs directly from fapar and ppfd
-    I_abs = FT(inputs["fapar"]) * FT(inputs["ppfd"])
-
-    βm =
-        Bool(inputs["do_soilmstress"]) ?
-        quadratic_soil_moisture_stress(FT(inputs["soilm"])) : FT(1.0)
-
-    return ClimaLand.Canopy.PModelDrivers(
-        T_canopy = T_canopy,
-        I_abs = I_abs,
-        ca = ca,
-        P_air = P_air,
-        VPD = VPD,
-        βm = βm,
-    )
-end
-
 """
     create_pmodel_parameters(inputs::Dict{String, Any}, FT)
 
@@ -191,7 +161,7 @@ end
                     println("Running test case: $testcase_name with FT = $FT")
 
                 # prepare constants, parameters, and drivers for the current FT
-                constants = PModelConstants(FT)
+                constants = PModelConstants{FT}()
                 parameters = PModelParameters(inputs, FT)
 
                 T_canopy = FT(inputs["tc"] + 273.15)  # Convert from Celsius to Kelvin