remove 3 vars from canopy cache (#1179)

Author: kmdeck

lib/ClimaLandSimulations/src/utilities/climaland_output_dataframe.jl

@@ -45,10 +45,9 @@ function make_output_df(
     output_list = vcat(
         (1, :SW_u),
         (1, :LW_u),
-        (1, :canopy, :conductance, :gs),
+        (1, :canopy, :conductance, :r_stomata_canopy),
         (1, :canopy, :autotrophic_respiration, :Ra),
         (1, :canopy, :photosynthesis, :GPP),
-        (1, :canopy, :hydraulics, :β),
         (1, :canopy, :hydraulics, :area_index, :leaf),
         (1, :canopy, :turbulent_fluxes, :transpiration),
         (1, :canopy, :turbulent_fluxes, :lhf),

lib/ClimaLandSimulations/src/utilities/makie_plots.jl

@@ -197,19 +197,11 @@ function timeseries_H2O_fig(
         color = :blue,
     )
 
-    # Moisture stress
-    p_MS = lines!(
-        ax_MS,
-        datetime2unix.(climaland.DateTime),
-        climaland.β,
-        color = :green,
-    ) # not sure about units
-
     # Stomatal conductance
     p_SC = lines!(
         ax_SC,
         datetime2unix.(climaland.DateTime),
-        climaland.gs,
+        climaland.r_stomata_canopy,
         color = :green,
     )
 
@@ -244,7 +236,7 @@ function timeseries_H2O_fig(
                 datetime2unix(climaland.DateTime[1]),
                 datetime2unix(climaland.DateTime[end]),
             ),
-        ) for axes in [ax_H2O, ax_H2O_rain, ax_MS, ax_SC]
+        ) for axes in [ax_H2O, ax_H2O_rain, ax_SC]
     ]
 
     fig

src/ClimaLand.jl

@@ -2,7 +2,7 @@ module ClimaLand
 using DocStringExtensions
 
 using ClimaCore
-using LazyBroadcast: @lazy
+using LazyBroadcast: lazy
 import ClimaCore: Fields, Spaces
 
 include("shared_utilities/Parameters.jl")

src/diagnostics/Diagnostics.jl

@@ -14,6 +14,7 @@ import ..LandModel
 
 import ..Soil: EnergyHydrology
 
+import ..Canopy: medlyn_conductance, medlyn_term, moisture_stress
 import ..Domains:
     top_center_to_surface, AbstractDomain, SphericalShell, HybridBox
 

src/diagnostics/land_compute_methods.jl

@@ -82,7 +82,50 @@ end
 @diagnostic_compute "autotrophic_respiration" Union{SoilCanopyModel, LandModel} p.canopy.autotrophic_respiration.Ra
 
 # Canopy - Conductance
-@diagnostic_compute "stomatal_conductance" Union{SoilCanopyModel, LandModel} p.canopy.conductance.gs
+function compute_stomatal_conductance!(
+    out,
+    Y,
+    p,
+    t,
+    land_model::Union{SoilCanopyModel, LandModel},
+)
+    canopy = land_model.canopy
+    (; g1, g0, Drel) = canopy.conductance.parameters
+    earth_param_set = canopy.parameters.earth_param_set
+    thermo_params = LP.thermodynamic_parameters(earth_param_set)
+
+    if isnothing(out)
+        return medlyn_conductance.(
+            g0,
+            Drel,
+            medlyn_term.(
+                g1,
+                p.drivers.T,
+                p.drivers.P,
+                p.drivers.q,
+                thermo_params,
+            ),
+            p.canopy.photosynthesis.An,
+            p.drivers.c_co2,
+        )
+
+
+    else
+        @. out = medlyn_conductance(
+            g0,
+            Drel,
+            medlyn_term(
+                g1,
+                p.drivers.T,
+                p.drivers.P,
+                p.drivers.q,
+                thermo_params,
+            ),
+            p.canopy.photosynthesis.An,
+            p.drivers.c_co2,
+        )
+    end
+end
 
 function compute_canopy_transpiration!(
     out,
@@ -132,7 +175,32 @@ end
     LandModel,
 } p.canopy.hydraulics.fa_roots
 @diagnostic_compute "leaf_area_index" Union{SoilCanopyModel, LandModel} p.canopy.hydraulics.area_index.leaf
-@diagnostic_compute "moisture_stress_factor" Union{SoilCanopyModel, LandModel} p.canopy.hydraulics.β
+# Canopy - Hydraulics
+function compute_moisture_stress_factor!(
+    out,
+    Y,
+    p,
+    t,
+    land_model::Union{SoilCanopyModel, LandModel},
+)
+    canopy = land_model.canopy
+
+    hydraulics = canopy.hydraulics
+    n_stem = hydraulics.n_stem
+    n_leaf = hydraulics.n_leaf
+    n = n_stem + n_leaf
+
+    earth_param_set = canopy.parameters.earth_param_set
+    grav = LP.grav(earth_param_set)
+    ρ_l = LP.ρ_cloud_liq(earth_param_set)
+    (; sc, pc) = canopy.photosynthesis.parameters
+    ψ = p.canopy.hydraulics.ψ
+    if isnothing(out)
+        return @. moisture_stress(ψ.:($$n) * ρ_l * grav, sc, pc)
+    else
+        @. out = moisture_stress(ψ.:($$n) * ρ_l * grav, sc, pc)
+    end
+end
 @diagnostic_compute "root_area_index" Union{SoilCanopyModel, LandModel} p.canopy.hydraulics.area_index.root
 @diagnostic_compute "stem_area_index" Union{SoilCanopyModel, LandModel} p.canopy.hydraulics.area_index.stem
 

src/standalone/Snow/Snow.jl

@@ -3,7 +3,7 @@ module Snow
 using DocStringExtensions
 import ...Parameters as LP
 using ClimaCore
-using LazyBroadcast: @lazy
+using LazyBroadcast: lazy
 using Thermodynamics
 using ClimaLand
 using ClimaLand:

src/standalone/Snow/snow_parameterizations.jl

@@ -498,7 +498,7 @@ function volumetric_energy_flux_falling_snow(atmos, p, parameters)
     _LH_f0 = LP.LH_f0(parameters.earth_param_set)
     _ρ_liq = LP.ρ_cloud_liq(parameters.earth_param_set)
     ρe_snow = -_LH_f0 * _ρ_liq
-    return @lazy @. ρe_snow * p.drivers.P_snow # per unit vol of liquid water
+    return @. lazy(ρe_snow * p.drivers.P_snow) # per unit vol of liquid water
 end
 
 """
@@ -513,8 +513,10 @@ CoupledAtmosphere, and the energy flux of the falling rain is passed in the
 cache `p`.  In that case, this should specify `atmos::PrescribedAtmosphere`.
 """
 function volumetric_energy_flux_falling_rain(atmos, p, parameters)
-    return @lazy @. volumetric_internal_energy_liq(p.drivers.T, parameters) *
-                    p.drivers.P_liq
+    return @. lazy(
+        volumetric_internal_energy_liq(p.drivers.T, parameters) *
+        p.drivers.P_liq,
+    )
 end
 
 

src/standalone/Vegetation/Canopy.jl

@@ -2,6 +2,7 @@ module Canopy
 using DocStringExtensions
 using Thermodynamics
 using ClimaLand
+using LazyBroadcast: lazy
 using ClimaCore
 using ClimaCore.MatrixFields
 import ClimaCore.MatrixFields: @name, ⋅
@@ -438,10 +439,6 @@ function ClimaLand.make_update_aux(
 
         # Other auxiliary variables being updated:
         Ra = p.canopy.autotrophic_respiration.Ra
-        β = p.canopy.hydraulics.β
-        medlyn_factor = p.canopy.conductance.medlyn_term
-        gs = p.canopy.conductance.gs # leaf level
-        rs_canopy = p.canopy.conductance.r_stomata_canopy
         An = p.canopy.photosynthesis.An
         GPP = p.canopy.photosynthesis.GPP
         Rd = p.canopy.photosynthesis.Rd
@@ -461,13 +458,13 @@ function ClimaLand.make_update_aux(
 
         # unpack parameters
         earth_param_set = canopy.parameters.earth_param_set
-        c = FT(LP.light_speed(earth_param_set))
-        planck_h = FT(LP.planck_constant(earth_param_set))
-        N_a = FT(LP.avogadro_constant(earth_param_set))
-        grav = FT(LP.grav(earth_param_set))
-        ρ_l = FT(LP.ρ_cloud_liq(earth_param_set))
-        R = FT(LP.gas_constant(earth_param_set))
-        T_freeze = FT(LP.T_freeze(earth_param_set))
+        c = LP.light_speed(earth_param_set)
+        planck_h = LP.planck_constant(earth_param_set)
+        N_a = LP.avogadro_constant(earth_param_set)
+        grav = LP.grav(earth_param_set)
+        ρ_l = LP.ρ_cloud_liq(earth_param_set)
+        R = LP.gas_constant(earth_param_set)
+        T_freeze = LP.T_freeze(earth_param_set)
         thermo_params = earth_param_set.thermo_params
         (; G_Function, Ω, λ_γ_PAR) = canopy.radiative_transfer.parameters
         energy_per_mole_photon_par = planck_h * c / λ_γ_PAR * N_a
@@ -560,14 +557,14 @@ function ClimaLand.make_update_aux(
         # We update the fa[n_stem+n_leaf] element once we have computed transpiration, below
         # update photosynthesis and conductance terms
         # This should still use T_air, P_air, q_air
-        medlyn_factor .=
-            medlyn_term.(g1, T_air, P_air, q_air, Ref(thermo_params))
+        medlyn_factor =
+            @. lazy(medlyn_term(g1, T_air, P_air, q_air, thermo_params))
         # Anywhere we use an Arrhenius factor, use T_canopy instead T_air
         T_canopy = canopy_temperature(canopy.energy, canopy, Y, p, t)
 
         # update moisture stress
         i_end = n_stem + n_leaf
-        @. β = moisture_stress(ψ.:($$i_end) * ρ_l * grav, sc, pc)
+        β = @. lazy(moisture_stress(ψ.:($$i_end) * ρ_l * grav, sc, pc))
 
         # Update Rd, An, Vcmax25 (if applicable to model) in place
         Vcmax25 = p.canopy.photosynthesis.Vcmax25
@@ -600,8 +597,14 @@ function ClimaLand.make_update_aux(
             par_d,
         )
         @. GPP = compute_GPP(An, extinction_coeff(G_Function, cosθs), LAI, Ω)
-        @. gs = medlyn_conductance(g0, Drel, medlyn_factor, An, c_co2_air)
-        @. rs_canopy = 1 / upscale_leaf_conductance(gs, LAI, T_air, R, P_air)
+        @. p.canopy.conductance.r_stomata_canopy =
+            1 / upscale_leaf_conductance(
+                medlyn_conductance(g0, Drel, medlyn_factor, An, c_co2_air), #conductance, leaf level
+                LAI,
+                T_air,
+                R,
+                P_air,
+            )
         # update autotrophic respiration
         h_canopy = hydraulics.compartment_surfaces[end]
         @. Ra = compute_autrophic_respiration(

src/standalone/Vegetation/PlantHydraulics.jl

@@ -287,13 +287,11 @@ prognostic_vars(model::PlantHydraulicsModel) = (:ϑ_l,)
 
 A function which returns the names of the auxiliary
 variables of the `PlantHydraulicsModel`,
-the transpiration stress factor `β` (unitless),
 the water potential `ψ` (m), the volume flux*cross section `fa` (1/s),
 and the volume flux*root cross section in the roots `fa_roots` (1/s),
 where the cross section can be represented by an area index.
 """
-auxiliary_vars(model::PlantHydraulicsModel) =
-    (:β, :ψ, :fa, :fa_roots, :area_index)
+auxiliary_vars(model::PlantHydraulicsModel) = (:ψ, :fa, :fa_roots, :area_index)
 
 """
     ClimaLand.prognostic_types(model::PlantHydraulicsModel{FT}) where {FT}
@@ -310,14 +308,13 @@ ClimaLand.prognostic_domain_names(::PlantHydraulicsModel) = (:surface,)
 Defines the auxiliary types for the PlantHydraulicsModel.
 """
 ClimaLand.auxiliary_types(model::PlantHydraulicsModel{FT}) where {FT} = (
-    FT,
     NTuple{model.n_stem + model.n_leaf, FT},
     NTuple{model.n_stem + model.n_leaf, FT},
     FT,
     NamedTuple{(:root, :stem, :leaf), Tuple{FT, FT, FT}},
 )
 ClimaLand.auxiliary_domain_names(::PlantHydraulicsModel) =
-    (:surface, :surface, :surface, :surface, :surface)
+    (:surface, :surface, :surface, :surface)
 
 function clip(x::FT, threshold::FT) where {FT}
     x > threshold ? x : FT(0)

src/standalone/Vegetation/stomatalconductance.jl

@@ -38,9 +38,6 @@ end
 
 ClimaLand.name(model::AbstractStomatalConductanceModel) = :conductance
 
-ClimaLand.auxiliary_vars(model::MedlynConductanceModel) =
-    (:medlyn_term, :gs, :r_stomata_canopy)
-ClimaLand.auxiliary_types(model::MedlynConductanceModel{FT}) where {FT} =
-    (FT, FT, FT)
-ClimaLand.auxiliary_domain_names(::MedlynConductanceModel) =
-    (:surface, :surface, :surface)
+ClimaLand.auxiliary_vars(model::MedlynConductanceModel) = (:r_stomata_canopy,)
+ClimaLand.auxiliary_types(model::MedlynConductanceModel{FT}) where {FT} = (FT,)
+ClimaLand.auxiliary_domain_names(::MedlynConductanceModel) = (:surface,)