From 4d9a0df70421d4177efc709cc51c436d05a85a41 Mon Sep 17 00:00:00 2001 From: Edward Speer Date: Wed, 15 May 2024 11:49:09 -0700 Subject: [PATCH 01/10] Add BigLeaf Hydraulics Model --- docs/src/APIs/canopy/PlantHydraulics.md | 1 + src/integrated/soil_canopy_model.jl | 89 +++++- src/standalone/Vegetation/Canopy.jl | 4 +- src/standalone/Vegetation/PlantHydraulics.jl | 283 ++++++++++++++++++- 4 files changed, 364 insertions(+), 13 deletions(-) diff --git a/docs/src/APIs/canopy/PlantHydraulics.md b/docs/src/APIs/canopy/PlantHydraulics.md index 4e7ab468d3..8be7095dfb 100644 --- a/docs/src/APIs/canopy/PlantHydraulics.md +++ b/docs/src/APIs/canopy/PlantHydraulics.md @@ -7,6 +7,7 @@ CurrentModule = ClimaLand.PlantHydraulics ```@docs ClimaLand.PlantHydraulics.PlantHydraulicsModel +ClimaLand.PlantHydraulics.BigLeafHydraulicsModel ``` ## Plant Hydraulics Diagnostic Variables diff --git a/src/integrated/soil_canopy_model.jl b/src/integrated/soil_canopy_model.jl index 754f9184b3..7f1871cffc 100644 --- a/src/integrated/soil_canopy_model.jl +++ b/src/integrated/soil_canopy_model.jl @@ -224,13 +224,14 @@ lsm_aux_domain_names(m::SoilCanopyModel) = FT, MM <: Soil.Biogeochemistry.SoilCO2Model{FT}, SM <: Soil.RichardsModel{FT}, - RM <: Canopy.CanopyModel{FT} + RM <: Canopy.CanopyModel{FT, PlantHydraulics.PlantHydraulicsModel{FT}} } A method which makes a function; the returned function updates the additional auxiliary variables for the integrated model, as well as updates the boundary auxiliary variables for all component -models. +models. This version is for an integrated canopy model which suports a +multi-layered hydraulics model. This function is called each ode function evaluation, prior to the tendency function evaluation. @@ -241,7 +242,7 @@ function make_update_boundary_fluxes( FT, MM <: Soil.Biogeochemistry.SoilCO2Model{FT}, SM <: Soil.EnergyHydrology{FT}, - RM <: Canopy.CanopyModel{FT}, + RM <: Canopy.CanopyModel{FT, PlantHydraulics.PlantHydraulicsModel{FT}}, } update_soil_bf! = make_update_boundary_fluxes(land.soil) update_soilco2_bf! = make_update_boundary_fluxes(land.soilco2) @@ -300,6 +301,86 @@ function make_update_boundary_fluxes( return update_boundary_fluxes! end +""" + make_update_boundary_fluxes( + land::SoilCanopyModel{FT, MM, SM, RM}, + ) where { + FT, + MM <: Soil.Biogeochemistry.SoilCO2Model{FT}, + SM <: Soil.RichardsModel{FT}, + RM <: Canopy.CanopyModel{FT, PlantHydraulics.BigLeafHydraulicsModel{FT}} + } + +A method which makes a function; the returned function updates the additional +auxiliary variables for the integrated model, as well as updates the boundary +auxiliary variables for all component models. This version is for an integrated +canopy model which supports only a big leaf hydraulics model. + +This function is called each ode function evaluation, prior to the tendency +function evaluation. +""" +function make_update_boundary_fluxes( + land::SoilCanopyModel{FT, MM, SM, RM}, +) where { + FT, + MM <: Soil.Biogeochemistry.SoilCO2Model{FT}, + SM <: Soil.EnergyHydrology{FT}, + RM <: Canopy.CanopyModel{FT, PlantHydraulics.BigLeafHydraulicsModel{FT}}, +} + update_soil_bf! = make_update_boundary_fluxes(land.soil) + update_soilco2_bf! = make_update_boundary_fluxes(land.soilco2) + update_canopy_bf! = make_update_boundary_fluxes(land.canopy) + function update_boundary_fluxes!(p, Y, t) + # update root extraction + z = + ClimaCore.Fields.coordinate_field( + land.soil.domain.space.subsurface, + ).z + (; h_stem, h_leaf, conductivity_model) = land.canopy.hydraulics.parameters + (label, midpoint) = h_stem > 0 ? (:stem, h_stem / 2) : (:leaf, h_leaf / 2) + + area_index = p.canopy.hydraulics.area_index + + above_ground_area_index = getproperty(area_index, label) + + @. p.root_extraction = + PlantHydraulics.flux( + z, + midpoint, + p.soil.ψ, + p.canopy.hydraulics.ψ.:1, + PlantHydraulics.hydraulic_conductivity( + conductivity_model, + p.soil.ψ, + ), + PlantHydraulics.hydraulic_conductivity( + conductivity_model, + p.canopy.hydraulics.ψ.:1, + ), + ) * + (land.canopy.hydraulics.parameters.root_distribution(z)) + @. p.root_energy_extraction = + p.root_extraction * ClimaLand.Soil.volumetric_internal_energy_liq( + p.soil.T, + land.soil.parameters, + ) + + # Radiation + lsm_radiant_energy_fluxes!( + p, + land.canopy.radiative_transfer, + land.canopy, + land.soil, + Y, + t, + ) + update_soil_bf!(p, Y, t) + update_canopy_bf!(p, Y, t) + update_soilco2_bf!(p, Y, t) + end + return update_boundary_fluxes! +end + """ lsm_radiant_energy_fluxes!(p, canopy_radiation::Canopy.AbstractRadiationModel{FT}, @@ -475,7 +556,7 @@ roots and soil at each soil layer. function PlantHydraulics.root_water_flux_per_ground_area!( fa::ClimaCore.Fields.Field, s::PrognosticSoil, - model::Canopy.PlantHydraulics.PlantHydraulicsModel{FT}, + model::Union{Canopy.PlantHydraulics.PlantHydraulicsModel{FT}, Canopy.PlantHydraulics.BigLeafHydraulicsModel{FT}}, Y::ClimaCore.Fields.FieldVector, p::NamedTuple, t, diff --git a/src/standalone/Vegetation/Canopy.jl b/src/standalone/Vegetation/Canopy.jl index 1e4baf9194..18b3fa3a28 100644 --- a/src/standalone/Vegetation/Canopy.jl +++ b/src/standalone/Vegetation/Canopy.jl @@ -488,9 +488,7 @@ function ClimaLand.make_update_aux( # update plant hydraulics aux hydraulics = canopy.hydraulics - n_stem = hydraulics.n_stem - n_leaf = hydraulics.n_leaf - PlantHydraulics.lai_consistency_check.(n_stem, n_leaf, area_index) + PlantHydraulics.lai_consistency_check.(hydraulics, area_index) (; retention_model, conductivity_model, S_s, ν) = hydraulics.parameters # We can index into a field of Tuple{FT} to extract a field of FT # using the following notation: field.:index diff --git a/src/standalone/Vegetation/PlantHydraulics.jl b/src/standalone/Vegetation/PlantHydraulics.jl index 07de4f6909..ff8053e364 100644 --- a/src/standalone/Vegetation/PlantHydraulics.jl +++ b/src/standalone/Vegetation/PlantHydraulics.jl @@ -22,6 +22,7 @@ import ClimaLand: prognostic_domain_names, name export PlantHydraulicsModel, + BigLeafHydraulicsModel, AbstractPlantHydraulicsModel, flux, effective_saturation, @@ -87,9 +88,7 @@ end """ lai_consistency_check( - n_stem::Int64, - n_leaf::Int64, - area_index::NamedTuple{(:root, :stem, :leaf), Tuple{FT, FT, FT}}, + model::PlantHydraulicsModel{FT}, ) where {FT} Carries out consistency checks using the area indices supplied and the number of @@ -106,10 +105,11 @@ Specifically, this checks that: be zero. """ function lai_consistency_check( - n_stem::Int64, - n_leaf::Int64, + model::PlantHydraulicsModel{FT}, area_index::NamedTuple{(:root, :stem, :leaf), Tuple{FT, FT, FT}}, ) where {FT} + n_stem = model.n_stem + n_leaf = model.n_leaf @assert n_leaf > 0 if area_index.leaf > eps(FT) || area_index.stem > eps(FT) @assert area_index.root > eps(FT) @@ -316,7 +316,7 @@ Sets the canopy prescribed fields pertaining to the PlantHydraulics component (the area indices) with their values at time t. """ function ClimaLand.Canopy.set_canopy_prescribed_field!( - component::PlantHydraulicsModel{FT}, + component::Union{PlantHydraulicsModel{FT}, BigLeafHydraulicsModel{FT}}, p, t, ) where {FT} @@ -327,6 +327,25 @@ function ClimaLand.Canopy.set_canopy_prescribed_field!( end +""" + update_canopy_prescribed_field!(component::Union{PlantHydraulics{FT}, BigLeafHydraulicsModel{FT}}, + p, + t, + ) where {FT} + +Updates the canopy prescribed fields pertaining to the PlantHydraulics +component (the LAI only in this case) with their values at time t. +""" +function ClimaLand.Canopy.update_canopy_prescribed_field!( + component::Union{PlantHydraulicsModel{FT}, BigLeafHydraulicsModel{FT}}, + p, + t, +) where {FT} + (; LAIfunction) = component.parameters.ai_parameterization + evaluate!(p.canopy.hydraulics.area_index.leaf, LAIfunction, t) +end + + """ flux( z1, @@ -703,4 +722,256 @@ function transpiration_per_ground_area( @inbounds return p.canopy.conductance.transpiration end +# BIGLEAF HYDRAULICS MODEL + +""" + BigLeafHydraulicsModel{FT, PS, T} <: AbstractPlantHydraulicsModel{FT} + +Defines, and constructs instances of, the BigLeafHydraulicsModel type which is +used for simulation flux of water to/from the soil for the case of a big leaf +model, meaning that the canopy is treated as a single stem and single leaf +compartment. Water flows along roots of different depths, along the stem, to +the leaf, and ultimately being lost from the system through transpiration. Note +that the canopy height is specified as part of the BigLeafHydraulicsModel, along +with the area indices of the leaves, roots, and stems. + +This model can also be combined with the soil model using ClimaLand, in which +case the prognostic soil water content is used to determine root extraction, and +the transpiration is also computed diagnostically. In a global run with patches +of bare soil, you can "turn off" the canopy model (to get zero root extraction, +zero absorption, and emission, zero transpiration, and sensible heat flux from +the canopy) by setting the canopy heights to zero. + +Finally, the model can be used in Canopy standalone mode by prescribing +the soil matric potential at the root tips or flux in the roots. There is also the +option (intendend only for debugging) to use a prescribed transpiration rate. + +$(DocStringExtensions.FIELDS) +""" +struct BigLeafHydraulicsModel{FT, PS, T} <: AbstractPlantHydraulicsModel{FT} + "The height of the stem compartment, in meters" + h_stem::FT + "The height of the leaf compartment, in meters" + h_leaf::FT + "Parameters required by the BigLeafHydraulics model" + parameters::PS + "The transpiration mode, of type `AbstractTranspiration`" + transpiration::T +end + +function BigLeafHydraulicsModel{FT}(; + h_stem::FT, + h_leaf::FT, + parameters::PlantHydraulicsParameters{FT}, + transpiration::AbstractTranspiration{FT} = DiagnosticTranspiration{FT}(), +) where {FT} + args = (parameters, transpiration) + @assert h_stem >= 0 + @assert h_leaf >= 0 + return BigLeafHydraulicsModel{ + FT, + typeof.(args)..., + }( + h_stem, + h_leaf, + args..., + ) +end + +""" + prognostic_vars(model::BigLeafHydraulicsModel) + +A function which returns the names of the prognostic variables of the +`BigLeafHydraulicsModel`. +""" +prognostic_vars(model::BigLeafHydraulicsModel) = (:ϑ_l,) + +""" + auxiliary_vars(model::BigLeafHydraulicsModel) + +A function which returns the names of the auxiliary variables of the +`BigLeafHydraulicsModel`, 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::BigLeafHydraulicsModel) = + (:β, :ψ, :fa, :fa_roots, :area_index) + +""" + ClimaLand.prognostic_types(model::BigLeafHydraulicsModel{FT}) where {FT} + +Defines the prognostic types for the BigLeafHydraulicsModel. +""" +ClimaLand.prognostic_types(model::BigLeafHydraulicsModel{FT}) where {FT} = + (NTuple{2, FT},) +ClimaLand.prognostic_domain_names(::BigLeafHydraulicsModel) = (:surface,) + +""" + ClimaLand.auxiliary_types(model::BigLeafHydraulicsModel{FT}) where {FT} + +Defines the auxiliary types for the BigLeafHydraulicsModel. +""" +ClimaLand.auxiliarty_types(model::BigLeafHydraulicsModel{FT}) where {FT} = ( + FT, + NTuple{2, FT}, + NTuple{2, FT}, + FT, + NamedTuple{(:root, :stem, :leaf), Tuple{FT, FT, FT}}, +) +ClimaLand.auxiliary_domain_names(::BigLeafHydraulicsModel) = + (:surface, :surface, :surface, :surface, :surface) + +""" + make_compute_exp_tendency(model::BigLeafHydraulicsModel, _) + +A function which creates the compute_exp_tendency! function for the BigLeafHydraulicsModel. +The compute_exp_tendency! function must comply with a rhs function of SciMLBase.jl. + +Below, `fa` denotes a flux multiplied by the relevant cross section +(per unit area ground, or area index, AI). The tendency for the +ith compartment can be written then as: +∂ϑ[i]/∂t = 1/(AI*dz)[fa[i]-fa[i+1]). + +Note that if the area_index is zero because no plant is present, +AIdz is zero, and the fluxes `fa` appearing in the numerator are +zero because they are scaled by AI. + +To prevent dividing by zero, we change AI/(AI x dz)" to +"AI/max(AI x dz, eps(FT))" +""" +function make_compute_exp_tendency( + model::BigLeafHydraulicsmodel{FT}, + canopy, +) where {FT} + function compute_exp_tendency!(dY, Y, p, t) + area_index = p.canopy.hydraulics.area_index + fa = p.canopy.hydraulics.fa + fa_roots = p.canopy.hydraulics.fa_roots + labels = model.parameters.h_stem > 0 ? [:stem, :leaf] : [:leaf] + + @inbounds for i in eachindex(labels) + im1 = i - 1 + ip1 = i + 1 + compartment_type = labels[i] + dz = compartment_type == :stem ? model.h_stem : model.h_leaf - model.h_stem + # To prevent dividing by zero, change AI/(AI x dz)" to + # "AI/max(AI x dz, eps(FT))" + AIdz = + max.(getproperty(area_index, labels[i]) * (dz), eps(FT)) + if i == 1 + @inbounds @. dY.canopy.hydraulics.ϑ_l.:($$i) = + 1 / AIdz * (fa_roots - fa.:($$i)) + else + @inbounds @. dY.canopy.hydraulics.ϑ_l.:($$i) = + 1 / AIdz * (fa.:($$im1) - fa.:($$i)) + end + end + end + return compute_exp_tendency! +end + +""" + root_water_flux_per_ground_area!( + fa::ClimaCore.Fields.Field, + s::PrescribedSoil, + model::BigLeafHydraulicsModel{FT}, + Y::ClimaCore.Fields.FieldVector, + p::NamedTuple, + t, + ) where {FT} + +A method which computes the water flux between the soil and the stem, via the +roots, and multiplied by the RAI, in the case of a big leaf model running +without an integrated soil model. + +The returned flux is per unit ground area. This assumes that the stem compartment +is the first element of `Y.canopy.hydraulics.ϑ_l`. +""" +function root_water_flux_per_ground_area!( + fa::ClimaCore.Fields.Field, + s::PrescribedSoil, + model::BigLeafHydraulicsModel{FT}, + Y::ClimaCore.Fields.FieldVector, + p::NamedTuple, + t, +) where {FT} + (; conductivity_model, root_distribution) = model.parameters + area_index = p.canopy.hydraulics,area_index + # We can index into a field of Tuple{FT} to extract a field of FT + # using the following notation: field.:index + ψ_base = p.canopy.hydraulics.ψ.:1 + root_depths = s.root_depths + n_root_layers = length(root_depths) + ψ_soil::FT = s.ψ(t) + fa .= FT(0.0) + labels = model.parameters.h_stem > 0 ? [:stem, :leaf] : [:leaf] + @inbounds for i in 1:n_root_layers + above_ground_area_index = getproperty(area_index, labels[1]) + compartment = labels[i] + midpoint = compartment == :stem ? model.h_stem / 2 : model.h_stem + (model.h_leaf / 2) + if i != n_root_layers + @. fa += + flux( + root_depths[i], + midpoint, + ψ_soil, + ψ_base, + hydraulic_conductivity(conductivity_model, ψ_soil), + hydraulic_conductivity(conductivity_model, ψ_base), + ) * + root_distribution(root_depths[i]) * + (root_depths[i + 1] - root_depths[i]) * + (area_index.root + above_ground_area_index) / 2 + else + @. fa += + flux( + root_depths[i], + midpoint, + ψ_soil, + ψ_base, + hydraulic_conductivity(conductivity_model, ψ_soil), + hydraulic_conductivity(conductivity_model, ψ_base), + ) * + root_distribution(root_depths[i]) * + (model.h_leaf - root_depths[n_root_layers]) * + (area_index.root + above_ground_area_index) / 2 + end + end + + +""" + lai_consistency_check( + model::BigLeafHydraulicsModel{FT}, + ) where {FT} + +Carries out consistency checks using the area indices supplied and the heights +of the stem and leaf compartments. + +Note that it is possible to have a plant with 0 height stem compartment +but with a leaf compartment, and that a plant must have leaf compartments +(even if LAI = 0). + +Specifically, this checks that: +1. n_leaf > 0 +2. if LAI is nonzero or SAI is nonzero, RAI must be nonzero. +3. if SAI > 0, n_stem must be > 0 for consistency. If SAI == 0, n_stem must +be zero. +""" +function lai_consistency_check( + model::PlantHydraulicsModel{FT}, + area_index::NamedTuple{(:root, :stem, :leaf), Tuple{FT, FT, FT}}, +) where {FT} + if area_index.leaf > esp(FT) || area_index.stem > eso(FT) + @assert area_index.root > eps(FT) + end + # If the SAI is 0, then the stem compartment should be height 0 + if area_index.stem < eps(FT) + @assert model.h_stem == FT(0) + else + # if SAI is > 0, n_stem should be > 0 for consistency + @assert model.h_stem > 0 + end +end + end From bd368507a892ebaad0875750ded623f8014b042e Mon Sep 17 00:00:00 2001 From: Edward Speer Date: Wed, 15 May 2024 12:00:34 -0700 Subject: [PATCH 02/10] switch integrated tutorial to BigLeafHydraulics --- docs/tutorials/integrated/soil_canopy_tutorial.jl | 8 +++----- src/standalone/Vegetation/PlantHydraulics.jl | 1 + 2 files changed, 4 insertions(+), 5 deletions(-) diff --git a/docs/tutorials/integrated/soil_canopy_tutorial.jl b/docs/tutorials/integrated/soil_canopy_tutorial.jl index 3309fe6beb..f366f3e34c 100644 --- a/docs/tutorials/integrated/soil_canopy_tutorial.jl +++ b/docs/tutorials/integrated/soil_canopy_tutorial.jl @@ -215,7 +215,7 @@ canopy_component_types = (; radiative_transfer = Canopy.TwoStreamModel{FT}, photosynthesis = Canopy.FarquharModel{FT}, conductance = Canopy.MedlynConductanceModel{FT}, - hydraulics = Canopy.PlantHydraulicsModel{FT}, + hydraulics = Canopy.BigLeafHydraulicsModel{FT}, ); # Then provide arguments to the canopy radiative transfer, stomatal conductance, @@ -272,10 +272,8 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; plant_hydraulics_args = ( parameters = plant_hydraulics_ps, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_midpoints = compartment_midpoints, - compartment_surfaces = compartment_surfaces, + h_stem = h_stem, + h_leaf = h_leaf, ); # We may now collect all of the canopy component argument tuples into one diff --git a/src/standalone/Vegetation/PlantHydraulics.jl b/src/standalone/Vegetation/PlantHydraulics.jl index ff8053e364..5bcf1c6518 100644 --- a/src/standalone/Vegetation/PlantHydraulics.jl +++ b/src/standalone/Vegetation/PlantHydraulics.jl @@ -938,6 +938,7 @@ function root_water_flux_per_ground_area!( (area_index.root + above_ground_area_index) / 2 end end +end """ From 875c2eb408d73dfa5827fe7770b6873d2bdaf9a6 Mon Sep 17 00:00:00 2001 From: Edward Speer Date: Wed, 15 May 2024 12:31:28 -0700 Subject: [PATCH 03/10] autoformat --- .../integrated/soil_canopy_tutorial.jl | 7 +- src/integrated/soil_canopy_model.jl | 18 +- src/standalone/Vegetation/Canopy.jl | 2 +- src/standalone/Vegetation/PlantHydraulics.jl | 154 ++++++++---------- 4 files changed, 80 insertions(+), 101 deletions(-) diff --git a/docs/tutorials/integrated/soil_canopy_tutorial.jl b/docs/tutorials/integrated/soil_canopy_tutorial.jl index f366f3e34c..b88eb2eb18 100644 --- a/docs/tutorials/integrated/soil_canopy_tutorial.jl +++ b/docs/tutorials/integrated/soil_canopy_tutorial.jl @@ -270,11 +270,8 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; retention_model = retention_model, ) -plant_hydraulics_args = ( - parameters = plant_hydraulics_ps, - h_stem = h_stem, - h_leaf = h_leaf, -); +plant_hydraulics_args = + (parameters = plant_hydraulics_ps, h_stem = h_stem, h_leaf = h_leaf); # We may now collect all of the canopy component argument tuples into one # arguments tuple for the canopy component models. diff --git a/src/integrated/soil_canopy_model.jl b/src/integrated/soil_canopy_model.jl index 7f1871cffc..47a4dcd68f 100644 --- a/src/integrated/soil_canopy_model.jl +++ b/src/integrated/soil_canopy_model.jl @@ -336,14 +336,16 @@ function make_update_boundary_fluxes( ClimaCore.Fields.coordinate_field( land.soil.domain.space.subsurface, ).z - (; h_stem, h_leaf, conductivity_model) = land.canopy.hydraulics.parameters - (label, midpoint) = h_stem > 0 ? (:stem, h_stem / 2) : (:leaf, h_leaf / 2) + (; h_stem, h_leaf, conductivity_model) = + land.canopy.hydraulics.parameters + (label, midpoint) = + h_stem > 0 ? (:stem, h_stem / 2) : (:leaf, h_leaf / 2) area_index = p.canopy.hydraulics.area_index above_ground_area_index = getproperty(area_index, label) - @. p.root_extraction = + @. p.root_extraction = PlantHydraulics.flux( z, midpoint, @@ -357,14 +359,13 @@ function make_update_boundary_fluxes( conductivity_model, p.canopy.hydraulics.ψ.:1, ), - ) * - (land.canopy.hydraulics.parameters.root_distribution(z)) + ) * (land.canopy.hydraulics.parameters.root_distribution(z)) @. p.root_energy_extraction = p.root_extraction * ClimaLand.Soil.volumetric_internal_energy_liq( p.soil.T, land.soil.parameters, ) - + # Radiation lsm_radiant_energy_fluxes!( p, @@ -556,7 +557,10 @@ roots and soil at each soil layer. function PlantHydraulics.root_water_flux_per_ground_area!( fa::ClimaCore.Fields.Field, s::PrognosticSoil, - model::Union{Canopy.PlantHydraulics.PlantHydraulicsModel{FT}, Canopy.PlantHydraulics.BigLeafHydraulicsModel{FT}}, + model::Union{ + Canopy.PlantHydraulics.PlantHydraulicsModel{FT}, + Canopy.PlantHydraulics.BigLeafHydraulicsModel{FT}, + }, Y::ClimaCore.Fields.FieldVector, p::NamedTuple, t, diff --git a/src/standalone/Vegetation/Canopy.jl b/src/standalone/Vegetation/Canopy.jl index 18b3fa3a28..a15afabb19 100644 --- a/src/standalone/Vegetation/Canopy.jl +++ b/src/standalone/Vegetation/Canopy.jl @@ -488,7 +488,7 @@ function ClimaLand.make_update_aux( # update plant hydraulics aux hydraulics = canopy.hydraulics - PlantHydraulics.lai_consistency_check.(hydraulics, area_index) + PlantHydraulics.lai_consistency_check(hydraulics, area_index) (; retention_model, conductivity_model, S_s, ν) = hydraulics.parameters # We can index into a field of Tuple{FT} to extract a field of FT # using the following notation: field.:index diff --git a/src/standalone/Vegetation/PlantHydraulics.jl b/src/standalone/Vegetation/PlantHydraulics.jl index 5bcf1c6518..09e562956a 100644 --- a/src/standalone/Vegetation/PlantHydraulics.jl +++ b/src/standalone/Vegetation/PlantHydraulics.jl @@ -86,44 +86,6 @@ function PrescribedSiteAreaIndex{FT}( PrescribedSiteAreaIndex{FT, typeof(LAIfunction)}(LAIfunction, SAI, RAI) end -""" - lai_consistency_check( - model::PlantHydraulicsModel{FT}, - ) where {FT} - -Carries out consistency checks using the area indices supplied and the number of -stem elements being modeled. - -Note that it is possible to have a plant with no stem compartments -but with leaf compartments, and that a plant must have leaf compartments -(even if LAI = 0). - -Specifically, this checks that: -1. n_leaf > 0 -2. if LAI is nonzero or SAI is nonzero, RAI must be nonzero. -3. if SAI > 0, n_stem must be > 0 for consistency. If SAI == 0, n_stem must -be zero. -""" -function lai_consistency_check( - model::PlantHydraulicsModel{FT}, - area_index::NamedTuple{(:root, :stem, :leaf), Tuple{FT, FT, FT}}, -) where {FT} - n_stem = model.n_stem - n_leaf = model.n_leaf - @assert n_leaf > 0 - if area_index.leaf > eps(FT) || area_index.stem > eps(FT) - @assert area_index.root > eps(FT) - end - # If there SAI is zero, there should be no stem compartment - if area_index.stem < eps(FT) - @assert n_stem == FT(0) - else - # if SAI is > 0, n_stem should be > 0 for consistency - @assert n_stem > 0 - end - -end - """ PlantHydraulicsParameters @@ -306,45 +268,42 @@ ClimaLand.auxiliary_domain_names(::PlantHydraulicsModel) = """ - set_canopy_prescribed_field!(component::PlantHydraulics{FT}, - p, - t, - ) where {FT} - - -Sets the canopy prescribed fields pertaining to the PlantHydraulics -component (the area indices) with their values at time t. -""" -function ClimaLand.Canopy.set_canopy_prescribed_field!( - component::Union{PlantHydraulicsModel{FT}, BigLeafHydraulicsModel{FT}}, - p, - t, -) where {FT} - (; LAIfunction, SAI, RAI) = component.parameters.ai_parameterization - evaluate!(p.canopy.hydraulics.area_index.leaf, LAIfunction, t) - @. p.canopy.hydraulics.area_index.stem = SAI - @. p.canopy.hydraulics.area_index.root = RAI -end + lai_consistency_check( + model::PlantHydraulicsModel{FT}, + ) where {FT} +Carries out consistency checks using the area indices supplied and the number of +stem elements being modeled. -""" - update_canopy_prescribed_field!(component::Union{PlantHydraulics{FT}, BigLeafHydraulicsModel{FT}}, - p, - t, - ) where {FT} +Note that it is possible to have a plant with no stem compartments +but with leaf compartments, and that a plant must have leaf compartments +(even if LAI = 0). -Updates the canopy prescribed fields pertaining to the PlantHydraulics -component (the LAI only in this case) with their values at time t. +Specifically, this checks that: +1. n_leaf > 0 +2. if LAI is nonzero or SAI is nonzero, RAI must be nonzero. +3. if SAI > 0, n_stem must be > 0 for consistency. If SAI == 0, n_stem must +be zero. """ -function ClimaLand.Canopy.update_canopy_prescribed_field!( - component::Union{PlantHydraulicsModel{FT}, BigLeafHydraulicsModel{FT}}, - p, - t, +function lai_consistency_check( + model::PlantHydraulicsModel{FT}, + area_index::NamedTuple{(:root, :stem, :leaf), Tuple{FT, FT, FT}}, ) where {FT} - (; LAIfunction) = component.parameters.ai_parameterization - evaluate!(p.canopy.hydraulics.area_index.leaf, LAIfunction, t) -end + n_stem = model.n_stem + n_leaf = model.n_leaf + @assert n_leaf > 0 + if area_index.leaf > eps(FT) || area_index.stem > eps(FT) + @assert area_index.root > eps(FT) + end + # If there SAI is zero, there should be no stem compartment + if area_index.stem < eps(FT) + @assert n_stem == FT(0) + else + # if SAI is > 0, n_stem should be > 0 for consistency + @assert n_stem > 0 + end +end """ flux( @@ -768,14 +727,7 @@ function BigLeafHydraulicsModel{FT}(; args = (parameters, transpiration) @assert h_stem >= 0 @assert h_leaf >= 0 - return BigLeafHydraulicsModel{ - FT, - typeof.(args)..., - }( - h_stem, - h_leaf, - args..., - ) + return BigLeafHydraulicsModel{FT, typeof.(args)...}(h_stem, h_leaf, args...) end """ @@ -803,7 +755,7 @@ auxiliary_vars(model::BigLeafHydraulicsModel) = Defines the prognostic types for the BigLeafHydraulicsModel. """ -ClimaLand.prognostic_types(model::BigLeafHydraulicsModel{FT}) where {FT} = +ClimaLand.prognostic_types(model::BigLeafHydraulicsModel{FT}) where {FT} = (NTuple{2, FT},) ClimaLand.prognostic_domain_names(::BigLeafHydraulicsModel) = (:surface,) @@ -812,7 +764,7 @@ ClimaLand.prognostic_domain_names(::BigLeafHydraulicsModel) = (:surface,) Defines the auxiliary types for the BigLeafHydraulicsModel. """ -ClimaLand.auxiliarty_types(model::BigLeafHydraulicsModel{FT}) where {FT} = ( +ClimaLand.auxiliary_types(model::BigLeafHydraulicsModel{FT}) where {FT} = ( FT, NTuple{2, FT}, NTuple{2, FT}, @@ -841,7 +793,7 @@ To prevent dividing by zero, we change AI/(AI x dz)" to "AI/max(AI x dz, eps(FT))" """ function make_compute_exp_tendency( - model::BigLeafHydraulicsmodel{FT}, + model::BigLeafHydraulicsModel{FT}, canopy, ) where {FT} function compute_exp_tendency!(dY, Y, p, t) @@ -854,11 +806,12 @@ function make_compute_exp_tendency( im1 = i - 1 ip1 = i + 1 compartment_type = labels[i] - dz = compartment_type == :stem ? model.h_stem : model.h_leaf - model.h_stem + dz = + compartment_type == :stem ? model.h_stem : + model.h_leaf - model.h_stem # To prevent dividing by zero, change AI/(AI x dz)" to # "AI/max(AI x dz, eps(FT))" - AIdz = - max.(getproperty(area_index, labels[i]) * (dz), eps(FT)) + AIdz = max.(getproperty(area_index, labels[i]) * (dz), eps(FT)) if i == 1 @inbounds @. dY.canopy.hydraulics.ϑ_l.:($$i) = 1 / AIdz * (fa_roots - fa.:($$i)) @@ -897,7 +850,7 @@ function root_water_flux_per_ground_area!( t, ) where {FT} (; conductivity_model, root_distribution) = model.parameters - area_index = p.canopy.hydraulics,area_index + area_index = p.canopy.hydraulics, area_index # We can index into a field of Tuple{FT} to extract a field of FT # using the following notation: field.:index ψ_base = p.canopy.hydraulics.ψ.:1 @@ -909,7 +862,9 @@ function root_water_flux_per_ground_area!( @inbounds for i in 1:n_root_layers above_ground_area_index = getproperty(area_index, labels[1]) compartment = labels[i] - midpoint = compartment == :stem ? model.h_stem / 2 : model.h_stem + (model.h_leaf / 2) + midpoint = + compartment == :stem ? model.h_stem / 2 : + model.h_stem + (model.h_leaf / 2) if i != n_root_layers @. fa += flux( @@ -960,7 +915,7 @@ Specifically, this checks that: be zero. """ function lai_consistency_check( - model::PlantHydraulicsModel{FT}, + model::BigLeafHydraulicsModel{FT}, area_index::NamedTuple{(:root, :stem, :leaf), Tuple{FT, FT, FT}}, ) where {FT} if area_index.leaf > esp(FT) || area_index.stem > eso(FT) @@ -975,4 +930,27 @@ function lai_consistency_check( end end +# SHARED METHODS FOR BigLeafHydraulics AND PlantHydraulics + +""" + set_canopy_prescribed_field!(component::PlantHydraulics{FT}, + p, + t, + ) where {FT} + + +Sets the canopy prescribed fields pertaining to the PlantHydraulics +component (the area indices) with their values at time t. +""" +function ClimaLand.Canopy.set_canopy_prescribed_field!( + component::Union{PlantHydraulicsModel{FT}, BigLeafHydraulicsModel{FT}}, + p, + t, +) where {FT} + (; LAIfunction, SAI, RAI) = component.parameters.ai_parameterization + evaluate!(p.canopy.hydraulics.area_index.leaf, LAIfunction, t) + @. p.canopy.hydraulics.area_index.stem = SAI + @. p.canopy.hydraulics.area_index.root = RAI +end + end From ae4150733a598a85e2c574abc6498c1c6fca85f8 Mon Sep 17 00:00:00 2001 From: Edward Speer Date: Wed, 15 May 2024 14:21:39 -0700 Subject: [PATCH 04/10] broadcastable PLantHydraulics --- src/standalone/Vegetation/Canopy.jl | 2 +- src/standalone/Vegetation/PlantHydraulics.jl | 3 +++ 2 files changed, 4 insertions(+), 1 deletion(-) diff --git a/src/standalone/Vegetation/Canopy.jl b/src/standalone/Vegetation/Canopy.jl index a15afabb19..18b3fa3a28 100644 --- a/src/standalone/Vegetation/Canopy.jl +++ b/src/standalone/Vegetation/Canopy.jl @@ -488,7 +488,7 @@ function ClimaLand.make_update_aux( # update plant hydraulics aux hydraulics = canopy.hydraulics - PlantHydraulics.lai_consistency_check(hydraulics, area_index) + PlantHydraulics.lai_consistency_check.(hydraulics, area_index) (; retention_model, conductivity_model, S_s, ν) = hydraulics.parameters # We can index into a field of Tuple{FT} to extract a field of FT # using the following notation: field.:index diff --git a/src/standalone/Vegetation/PlantHydraulics.jl b/src/standalone/Vegetation/PlantHydraulics.jl index 09e562956a..33f275612d 100644 --- a/src/standalone/Vegetation/PlantHydraulics.jl +++ b/src/standalone/Vegetation/PlantHydraulics.jl @@ -49,6 +49,9 @@ abstract type AbstractPlantHydraulicsModel{FT} <: AbstractCanopyComponent{FT} en ClimaLand.name(::AbstractPlantHydraulicsModel) = :hydraulics +# Hydraulics models should be broadcastable +Base.broadcastable(x::AbstractPlantHydraulicsModel) = tuple(x) + """ AbstractTranspiration{FT <: AbstractFloat} From 2c6f8c7077d7412b0f78f54e8d260fe4813191f9 Mon Sep 17 00:00:00 2001 From: Edward Speer Date: Sun, 2 Jun 2024 00:29:42 -0700 Subject: [PATCH 05/10] remove PlantHydraulicsModel --- docs/src/APIs/canopy/PlantHydraulics.md | 1 - .../integrated/soil_canopy_tutorial.jl | 4 - .../standalone/Canopy/canopy_tutorial.jl | 12 +- .../fluxnet/US-Ha1/US-Ha1_simulation.jl | 2 - .../fluxnet/US-MOz/US-MOz_simulation.jl | 2 - .../fluxnet/US-NR1/US-NR1_simulation.jl | 2 - .../fluxnet/US-Var/US-Var_simulation.jl | 2 - .../integrated/fluxnet/fluxnet_domain.jl | 3 - experiments/integrated/fluxnet/ozark_pft.jl | 12 +- experiments/integrated/fluxnet/run_fluxnet.jl | 12 +- .../conservation/ozark_conservation_setup.jl | 10 +- .../performance/profile_allocations.jl | 12 +- .../fluxnet_sites/US-Ha1/US-Ha1_simulation.jl | 8 - .../fluxnet_sites/US-MOz/US-MOz_simulation.jl | 8 - .../fluxnet_sites/US-NR1/US-NR1_simulation.jl | 8 - .../fluxnet_sites/US-Var/US-Var_simulation.jl | 8 - .../Fluxnet/fluxnet_utilities/make_setup.jl | 6 - .../src/Fluxnet/run_fluxnet.jl | 12 +- .../src/utilities/make_domain.jl | 11 - src/integrated/soil_canopy_model.jl | 23 +- src/standalone/Vegetation/Canopy.jl | 28 +- src/standalone/Vegetation/PlantHydraulics.jl | 309 +--------------- .../Vegetation/canopy_boundary_fluxes.jl | 6 +- test/standalone/Vegetation/canopy_model.jl | 103 +----- .../Vegetation/plant_hydraulics_test.jl | 343 +----------------- 25 files changed, 97 insertions(+), 850 deletions(-) diff --git a/docs/src/APIs/canopy/PlantHydraulics.md b/docs/src/APIs/canopy/PlantHydraulics.md index 8be7095dfb..aaeecf693d 100644 --- a/docs/src/APIs/canopy/PlantHydraulics.md +++ b/docs/src/APIs/canopy/PlantHydraulics.md @@ -6,7 +6,6 @@ CurrentModule = ClimaLand.PlantHydraulics ## Models ```@docs -ClimaLand.PlantHydraulics.PlantHydraulicsModel ClimaLand.PlantHydraulics.BigLeafHydraulicsModel ``` diff --git a/docs/tutorials/integrated/soil_canopy_tutorial.jl b/docs/tutorials/integrated/soil_canopy_tutorial.jl index b88eb2eb18..120c63368f 100644 --- a/docs/tutorials/integrated/soil_canopy_tutorial.jl +++ b/docs/tutorials/integrated/soil_canopy_tutorial.jl @@ -68,12 +68,8 @@ f_root_to_shoot = FT(3.5) SAI = FT(0.00242) maxLAI = FT(4.2) plant_ν = FT(2.46e-4) -n_stem = Int64(1) -n_leaf = Int64(1) h_stem = FT(9) h_leaf = FT(9.5) -compartment_midpoints = [h_stem / 2, h_stem + h_leaf / 2] -compartment_surfaces = [zmax, h_stem, h_stem + h_leaf] land_domain = Column(; zlim = (zmin, zmax), nelements = nelements); # - We will be using prescribed atmospheric and radiative drivers from the diff --git a/docs/tutorials/standalone/Canopy/canopy_tutorial.jl b/docs/tutorials/standalone/Canopy/canopy_tutorial.jl index e84449b1d7..58b65af0a7 100644 --- a/docs/tutorials/standalone/Canopy/canopy_tutorial.jl +++ b/docs/tutorials/standalone/Canopy/canopy_tutorial.jl @@ -84,12 +84,8 @@ f_root_to_shoot = FT(3.5) SAI = FT(0.00242) maxLAI = FT(4.2) plant_ν = FT(2.46e-4) # kg/m^2 -n_stem = Int64(1) -n_leaf = Int64(1) h_stem = FT(9) h_leaf = FT(9.5) -compartment_midpoints = [h_stem / 2, h_stem + h_leaf / 2] -compartment_surfaces = [zmax, h_stem, h_stem + h_leaf] land_domain = Point(; z_sfc = FT(0.0)) # - We will be using prescribed atmospheric and radiative drivers from the @@ -233,12 +229,10 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; # Define the remaining variables required for the plant hydraulics model. -plant_hydraulics = PlantHydraulics.PlantHydraulicsModel{FT}(; +plant_hydraulics = PlantHydraulics.BigLeafHydraulicsModel{FT}(; parameters = plant_hydraulics_ps, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_surfaces = compartment_surfaces, - compartment_midpoints = compartment_midpoints, + h_stem = h_stem, + h_leaf = h_leaf, ); # Now, instantiate the canopy model, using the atmospheric and radiative diff --git a/experiments/integrated/fluxnet/US-Ha1/US-Ha1_simulation.jl b/experiments/integrated/fluxnet/US-Ha1/US-Ha1_simulation.jl index 6afcd4670d..1eb17ad89d 100644 --- a/experiments/integrated/fluxnet/US-Ha1/US-Ha1_simulation.jl +++ b/experiments/integrated/fluxnet/US-Ha1/US-Ha1_simulation.jl @@ -9,8 +9,6 @@ dz_bottom = FT(1.5) dz_top = FT(0.025) # Stem and leaf compartments and their heights: -n_stem = Int64(1) -n_leaf = Int64(1) h_leaf = FT(12) # m h_stem = FT(14) # m diff --git a/experiments/integrated/fluxnet/US-MOz/US-MOz_simulation.jl b/experiments/integrated/fluxnet/US-MOz/US-MOz_simulation.jl index 993e855acc..115144cce8 100644 --- a/experiments/integrated/fluxnet/US-MOz/US-MOz_simulation.jl +++ b/experiments/integrated/fluxnet/US-MOz/US-MOz_simulation.jl @@ -9,8 +9,6 @@ dz_bottom = FT(1.5) dz_top = FT(0.025) # Stem and leaf compartments and their heights: -n_stem = Int64(1) -n_leaf = Int64(1) h_stem = FT(9) # m h_leaf = FT(9.5) # m diff --git a/experiments/integrated/fluxnet/US-NR1/US-NR1_simulation.jl b/experiments/integrated/fluxnet/US-NR1/US-NR1_simulation.jl index e27c0b4827..3ad13a9079 100644 --- a/experiments/integrated/fluxnet/US-NR1/US-NR1_simulation.jl +++ b/experiments/integrated/fluxnet/US-NR1/US-NR1_simulation.jl @@ -9,8 +9,6 @@ dz_bottom = FT(1.25) dz_top = FT(0.05) # Stem and leaf compartments and their heights: -n_stem = Int64(1) -n_leaf = Int64(1) h_leaf = FT(6.5) # m h_stem = FT(7.5) # m diff --git a/experiments/integrated/fluxnet/US-Var/US-Var_simulation.jl b/experiments/integrated/fluxnet/US-Var/US-Var_simulation.jl index c19746d33c..093b73a16a 100644 --- a/experiments/integrated/fluxnet/US-Var/US-Var_simulation.jl +++ b/experiments/integrated/fluxnet/US-Var/US-Var_simulation.jl @@ -9,8 +9,6 @@ dz_bottom = FT(1.0) dz_top = FT(0.05) # Stem and leaf compartments and their heights: -n_stem = Int64(0) -n_leaf = Int64(1) h_leaf = FT(0.7) # m h_stem = FT(0) # m diff --git a/experiments/integrated/fluxnet/fluxnet_domain.jl b/experiments/integrated/fluxnet/fluxnet_domain.jl index b5d1be8012..0f969ad59c 100644 --- a/experiments/integrated/fluxnet/fluxnet_domain.jl +++ b/experiments/integrated/fluxnet/fluxnet_domain.jl @@ -8,9 +8,6 @@ nelements = 20 zmin = FT(-10) zmax = FT(0) h_canopy = h_stem + h_leaf -compartment_midpoints = - n_stem > 0 ? [h_stem / 2, h_stem + h_leaf / 2] : [h_leaf / 2] -compartment_surfaces = n_stem > 0 ? [zmax, h_stem, h_canopy] : [zmax, h_leaf] land_domain = Column(; zlim = (zmin, zmax), diff --git a/experiments/integrated/fluxnet/ozark_pft.jl b/experiments/integrated/fluxnet/ozark_pft.jl index 1537c3c77b..66accb4e81 100644 --- a/experiments/integrated/fluxnet/ozark_pft.jl +++ b/experiments/integrated/fluxnet/ozark_pft.jl @@ -168,7 +168,7 @@ canopy_component_types = (; radiative_transfer = Canopy.TwoStreamModel{FT}, photosynthesis = Canopy.FarquharModel{FT}, conductance = Canopy.MedlynConductanceModel{FT}, - hydraulics = Canopy.PlantHydraulicsModel{FT}, + hydraulics = Canopy.BigLeafHydraulicsModel{FT}, energy = Canopy.BigLeafEnergyModel{FT}, ) # Individual Component arguments @@ -211,10 +211,8 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; ) plant_hydraulics_args = ( parameters = plant_hydraulics_ps, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_midpoints = compartment_midpoints, - compartment_surfaces = compartment_surfaces, + h_stem = h_stem, + h_leaf = h_leaf, ) energy_args = (parameters = Canopy.BigLeafEnergyParameters{FT}(ac_canopy),) @@ -276,12 +274,12 @@ Y.soil.ρe_int = Y.soilco2.C .= FT(0.000412) # set to atmospheric co2, mol co2 per mol air ψ_stem_0 = FT(-1e5 / 9800) # pressure in the leaf divided by rho_liquid*gravitational acceleration [m] ψ_leaf_0 = FT(-2e5 / 9800) -ψ_comps = n_stem > 0 ? [ψ_stem_0, ψ_leaf_0] : ψ_leaf_0 +ψ_comps = h_stem > 0 ? [ψ_stem_0, ψ_leaf_0] : ψ_leaf_0 S_l_ini = inverse_water_retention_curve.(retention_model, ψ_comps, plant_ν, plant_S_s) -for i in 1:(n_stem + n_leaf) +for i in 1:(h_stem < 0 + 1) Y.canopy.hydraulics.ϑ_l.:($i) .= augmented_liquid_fraction.(plant_ν, S_l_ini[i]) end diff --git a/experiments/integrated/fluxnet/run_fluxnet.jl b/experiments/integrated/fluxnet/run_fluxnet.jl index 321a2f0d46..82c6de378b 100644 --- a/experiments/integrated/fluxnet/run_fluxnet.jl +++ b/experiments/integrated/fluxnet/run_fluxnet.jl @@ -130,7 +130,7 @@ canopy_component_types = (; radiative_transfer = Canopy.TwoStreamModel{FT}, photosynthesis = Canopy.FarquharModel{FT}, conductance = Canopy.MedlynConductanceModel{FT}, - hydraulics = Canopy.PlantHydraulicsModel{FT}, + hydraulics = Canopy.BigLeafHydraulicsModel{FT}, energy = Canopy.BigLeafEnergyModel{FT}, ) # Individual Component arguments @@ -171,10 +171,8 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; ) plant_hydraulics_args = ( parameters = plant_hydraulics_ps, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_midpoints = compartment_midpoints, - compartment_surfaces = compartment_surfaces, + h_stem = h_stem, + h_leaf = h_leaf, ) energy_args = (parameters = Canopy.BigLeafEnergyParameters{FT}(ac_canopy),) @@ -236,12 +234,12 @@ Y.soil.ρe_int = Y.soilco2.C .= FT(0.000412) # set to atmospheric co2, mol co2 per mol air ψ_stem_0 = FT(-1e5 / 9800) # pressure in the leaf divided by rho_liquid*gravitational acceleration [m] ψ_leaf_0 = FT(-2e5 / 9800) -ψ_comps = n_stem > 0 ? [ψ_stem_0, ψ_leaf_0] : ψ_leaf_0 +ψ_comps = h_stem > 0 ? [ψ_stem_0, ψ_leaf_0] : ψ_leaf_0 S_l_ini = inverse_water_retention_curve.(retention_model, ψ_comps, plant_ν, plant_S_s) -for i in 1:(n_stem + n_leaf) +for i in 1:(h_stem > 0 + 1) Y.canopy.hydraulics.ϑ_l.:($i) .= augmented_liquid_fraction.(plant_ν, S_l_ini[i]) end diff --git a/experiments/integrated/performance/conservation/ozark_conservation_setup.jl b/experiments/integrated/performance/conservation/ozark_conservation_setup.jl index ebf0e98ae2..3c308abcd1 100644 --- a/experiments/integrated/performance/conservation/ozark_conservation_setup.jl +++ b/experiments/integrated/performance/conservation/ozark_conservation_setup.jl @@ -6,8 +6,6 @@ dz_bottom = FT(1.5) dz_top = FT(0.025) # Stem and leaf compartments and their heights: -n_stem = Int64(1) -n_leaf = Int64(1) h_stem = FT(9) # m h_leaf = FT(9.5) # m @@ -112,7 +110,7 @@ canopy_component_types = (; radiative_transfer = Canopy.TwoStreamModel{FT}, photosynthesis = Canopy.FarquharModel{FT}, conductance = Canopy.MedlynConductanceModel{FT}, - hydraulics = Canopy.PlantHydraulicsModel{FT}, + hydraulics = Canopy.BigLeafHydraulicsModel{FT}, energy = Canopy.BigLeafEnergyModel{FT}, ) # Individual Component arguments @@ -155,10 +153,8 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; ) plant_hydraulics_args = ( parameters = plant_hydraulics_ps, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_midpoints = compartment_midpoints, - compartment_surfaces = compartment_surfaces, + h_stem = h_stem, + h_leaf = h_leaf, ) # Canopy component args diff --git a/experiments/integrated/performance/profile_allocations.jl b/experiments/integrated/performance/profile_allocations.jl index bcc5e675fd..0506c05e6b 100644 --- a/experiments/integrated/performance/profile_allocations.jl +++ b/experiments/integrated/performance/profile_allocations.jl @@ -86,14 +86,10 @@ dz_bottom = FT(2.0) dz_top = FT(0.2) soil_depth = FT(5) z_sfc = FT(0) -n_stem = Int64(1) -n_leaf = Int64(1) h_stem = FT(9) # m h_leaf = FT(9.5) # m h_canopy = h_stem + h_leaf -compartment_midpoints = [h_stem / 2, h_stem + h_leaf / 2] -compartment_surfaces = [z_sfc, h_stem, h_canopy] land_domain = ClimaLand.Domains.SphericalShell(; radius = FT(6.3781e6), @@ -258,7 +254,7 @@ canopy_component_types = (; radiative_transfer = Canopy.TwoStreamModel{FT}, photosynthesis = Canopy.FarquharModel{FT}, conductance = Canopy.MedlynConductanceModel{FT}, - hydraulics = Canopy.PlantHydraulicsModel{FT}, + hydraulics = Canopy.BigLeafHydraulicsModel{FT}, energy = Canopy.BigLeafEnergyModel{FT}, ) # Individual Component arguments @@ -289,10 +285,8 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; ) plant_hydraulics_args = ( parameters = plant_hydraulics_ps, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_midpoints = compartment_midpoints, - compartment_surfaces = compartment_surfaces, + h_stem = h_stem, + h_leaf = h_leaf, ) # Canopy component args diff --git a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-Ha1/US-Ha1_simulation.jl b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-Ha1/US-Ha1_simulation.jl index bccdbc01ef..455114293f 100644 --- a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-Ha1/US-Ha1_simulation.jl +++ b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-Ha1/US-Ha1_simulation.jl @@ -4,8 +4,6 @@ export harvard_default_args function harvard_default_args(; dz_bottom = FT(1.5), dz_top = FT(0.025), - n_stem = Int64(1), - n_leaf = Int64(1), h_leaf = FT(12), # m h_stem = FT(14), # m t0 = Float64(120 * 3600 * 24), @@ -15,8 +13,6 @@ export harvard_default_args return ( dz_bottom = dz_bottom, dz_top = dz_top, - n_stem = n_stem, - n_leaf = n_leaf, h_stem = h_stem, h_leaf = h_leaf, t0 = t0, @@ -30,8 +26,6 @@ Named Tuple containing default simulation parameter for the Harvard site. function harvard_default_args(; dz_bottom = FT(1.5), dz_top = FT(0.025), - n_stem = Int64(1), - n_leaf = Int64(1), h_leaf = FT(12), # m h_stem = FT(14), # m t0 = Float64(120 * 3600 * 24), @@ -41,8 +35,6 @@ function harvard_default_args(; return ( dz_bottom = dz_bottom, dz_top = dz_top, - n_stem = n_stem, - n_leaf = n_leaf, h_stem = h_stem, h_leaf = h_leaf, t0 = t0, diff --git a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-MOz/US-MOz_simulation.jl b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-MOz/US-MOz_simulation.jl index dcbff11279..23d551438b 100644 --- a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-MOz/US-MOz_simulation.jl +++ b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-MOz/US-MOz_simulation.jl @@ -4,8 +4,6 @@ export ozark_default_args function ozark_default_args(; dz_bottom = FT(1.5), dz_top = FT(0.025), - n_stem = Int64(1), - n_leaf = Int64(1), h_stem = FT(9), h_leaf = FT(9.5), t0 = Float64(120 * 3600 * 24), @@ -15,8 +13,6 @@ export ozark_default_args return ( dz_bottom = dz_bottom, dz_top = dz_top, - n_stem = n_stem, - n_leaf = n_leaf, h_stem = h_stem, h_leaf = h_leaf, t0 = t0, @@ -30,8 +26,6 @@ Named Tuple containing default simulation parameter for the Ozark site. function ozark_default_args(; dz_bottom = FT(1.5), dz_top = FT(0.025), - n_stem = Int64(1), - n_leaf = Int64(1), h_stem = FT(9), h_leaf = FT(9.5), t0 = Float64(120 * 3600 * 24), @@ -41,8 +35,6 @@ function ozark_default_args(; return ( dz_bottom = dz_bottom, dz_top = dz_top, - n_stem = n_stem, - n_leaf = n_leaf, h_stem = h_stem, h_leaf = h_leaf, t0 = t0, diff --git a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-NR1/US-NR1_simulation.jl b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-NR1/US-NR1_simulation.jl index 0d3acf86be..4227a98414 100644 --- a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-NR1/US-NR1_simulation.jl +++ b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-NR1/US-NR1_simulation.jl @@ -4,8 +4,6 @@ export niwotridge_default_args function niwotridge_default_args(; dz_bottom = FT(1.25), dz_top = FT(0.05), - n_stem = Int64(1), - n_leaf = Int64(1), h_leaf = FT(6.5), # m, h_stem = FT(7.5), # m, t0 = Float64(120 * 3600 * 24), # start mid year to avoid snow, @@ -15,8 +13,6 @@ export niwotridge_default_args return ( dz_bottom = dz_bottom, dz_top = dz_top, - n_stem = n_stem, - n_leaf = n_leaf, h_stem = h_stem, h_leaf = h_leaf, t0 = t0, @@ -30,8 +26,6 @@ Named Tuple containing default simulation parameter for the Niwot Ridge site. function niwotridge_default_args(; dz_bottom = FT(1.25), dz_top = FT(0.05), - n_stem = Int64(1), - n_leaf = Int64(1), h_leaf = FT(6.5), # m, h_stem = FT(7.5), # m, t0 = Float64(120 * 3600 * 24), # start mid year to avoid snow, @@ -41,8 +35,6 @@ function niwotridge_default_args(; return ( dz_bottom = dz_bottom, dz_top = dz_top, - n_stem = n_stem, - n_leaf = n_leaf, h_stem = h_stem, h_leaf = h_leaf, t0 = t0, diff --git a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-Var/US-Var_simulation.jl b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-Var/US-Var_simulation.jl index 66dc5735c3..360ec63b2a 100644 --- a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-Var/US-Var_simulation.jl +++ b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-Var/US-Var_simulation.jl @@ -4,8 +4,6 @@ export vairaranch_default_args function vairaranch_default_args(; dz_bottom = FT(1.0), dz_top = FT(0.05), - n_stem = Int64(0), - n_leaf = Int64(1), h_leaf = FT(0.7), # m, h_stem = FT(0), # m, t0 = Float64(21 * 3600 * 24), # start day 21 of the year, @@ -15,8 +13,6 @@ export vairaranch_default_args return ( dz_bottom = dz_bottom, dz_top = dz_top, - n_stem = n_stem, - n_leaf = n_leaf, h_stem = h_stem, h_leaf = h_leaf, t0 = t0, @@ -30,8 +26,6 @@ Named Tuple containing default simulation parameter for the Vaira Ranch site. function vairaranch_default_args(; dz_bottom = FT(1.0), dz_top = FT(0.05), - n_stem = Int64(0), - n_leaf = Int64(1), h_leaf = FT(0.7), # m, h_stem = FT(0), # m, t0 = Float64(21 * 3600 * 24), # start day 21 of the year, @@ -41,8 +35,6 @@ function vairaranch_default_args(; return ( dz_bottom = dz_bottom, dz_top = dz_top, - n_stem = n_stem, - n_leaf = n_leaf, h_stem = h_stem, h_leaf = h_leaf, t0 = t0, diff --git a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_utilities/make_setup.jl b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_utilities/make_setup.jl index 5f80beee14..9e9434004f 100644 --- a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_utilities/make_setup.jl +++ b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_utilities/make_setup.jl @@ -6,8 +6,6 @@ export make_setup setup site specific domain and time stepping: - dz_bottom, the bottom of the soil domain (m) - dz_top, the top of the soil domain (m) -- n_stem, the number of stems -- n_leaf, the number of leaves - h_stem, the height of each stem (m) - h_leaf, the height of each leaf (m) - t0, the start time of the simulation (s) @@ -17,8 +15,6 @@ setup site specific domain and time stepping: function make_setup(; dz_bottom, dz_top, - n_stem, - n_leaf, h_stem, h_leaf, t0, @@ -28,8 +24,6 @@ function make_setup(; return ( dz_bottom = dz_bottom, dz_top = dz_top, - n_stem = n_stem, - n_leaf = n_leaf, h_stem = h_stem, h_leaf = h_leaf, t0 = t0, diff --git a/lib/ClimaLandSimulations/src/Fluxnet/run_fluxnet.jl b/lib/ClimaLandSimulations/src/Fluxnet/run_fluxnet.jl index 0c4a36fa11..d27c702fc1 100644 --- a/lib/ClimaLandSimulations/src/Fluxnet/run_fluxnet.jl +++ b/lib/ClimaLandSimulations/src/Fluxnet/run_fluxnet.jl @@ -93,7 +93,7 @@ function run_fluxnet( radiative_transfer = Canopy.TwoStreamModel{FT}, photosynthesis = Canopy.FarquharModel{FT}, conductance = Canopy.MedlynConductanceModel{FT}, - hydraulics = Canopy.PlantHydraulicsModel{FT}, + hydraulics = Canopy.BigLeafHydraulicsModel{FT}, energy = Canopy.BigLeafEnergyModel{FT}, ) # Individual Component arguments @@ -182,10 +182,8 @@ function run_fluxnet( ) plant_hydraulics_args = ( parameters = plant_hydraulics_ps, - n_stem = setup.n_stem, - n_leaf = setup.n_leaf, - compartment_midpoints = domain.compartment_midpoints, - compartment_surfaces = domain.compartment_surfaces, + h_stem, + h_leaf, ) energy_args = ( @@ -261,7 +259,7 @@ function run_fluxnet( Y.soilco2.C .= FT(0.000412) # set to atmospheric co2, mol co2 per mol air ψ_stem_0 = FT(-1e5 / 9800) # pressure in the leaf divided by rho_liquid*gravitational acceleration [m] ψ_leaf_0 = FT(-2e5 / 9800) - ψ_comps = setup.n_stem > 0 ? [ψ_stem_0, ψ_leaf_0] : ψ_leaf_0 + ψ_comps = setup.h_stem > 0 ? [ψ_stem_0, ψ_leaf_0] : ψ_leaf_0 S_l_ini = inverse_water_retention_curve.( @@ -271,7 +269,7 @@ function run_fluxnet( params.plant_hydraulics.S_s, ) - for i in 1:(setup.n_stem + setup.n_leaf) + for i in 1:(setup.h_stem > 0 + 1) Y.canopy.hydraulics.ϑ_l.:($i) .= augmented_liquid_fraction.(drivers.plant_ν, S_l_ini[i]) end diff --git a/lib/ClimaLandSimulations/src/utilities/make_domain.jl b/lib/ClimaLandSimulations/src/utilities/make_domain.jl index 42755bb4f5..45b6a86cf9 100644 --- a/lib/ClimaLandSimulations/src/utilities/make_domain.jl +++ b/lib/ClimaLandSimulations/src/utilities/make_domain.jl @@ -9,8 +9,6 @@ export make_domain setup the soil and canopy domain of the simulation. Returns: - h_canopy, -- compartment_midpoints, -- compartment_surfaces, - land_domain, - canopy_domain, """ @@ -18,13 +16,6 @@ function make_domain(setup, FT; nelements = 20, zmin = FT(-10), zmax = FT(0)) h_canopy = setup.h_stem + setup.h_leaf - compartment_midpoints = - setup.n_stem > 0 ? [setup.h_stem / 2, setup.h_stem + setup.h_leaf / 2] : - [setup.h_leaf / 2] - - compartment_surfaces = - setup.n_stem > 0 ? [zmax, setup.h_stem, h_canopy] : [zmax, setup.h_leaf] - land_domain = Column(; zlim = (zmin, zmax), nelements = nelements, @@ -38,8 +29,6 @@ function make_domain(setup, FT; nelements = 20, zmin = FT(-10), zmax = FT(0)) zmin = zmin, zmax = zmax, h_canopy = h_canopy, - compartment_midpoints = compartment_midpoints, - compartment_surfaces = compartment_surfaces, land_domain = land_domain, canopy_domain = canopy_domain, ) diff --git a/src/integrated/soil_canopy_model.jl b/src/integrated/soil_canopy_model.jl index 47a4dcd68f..c86b43d39a 100644 --- a/src/integrated/soil_canopy_model.jl +++ b/src/integrated/soil_canopy_model.jl @@ -224,7 +224,7 @@ lsm_aux_domain_names(m::SoilCanopyModel) = FT, MM <: Soil.Biogeochemistry.SoilCO2Model{FT}, SM <: Soil.RichardsModel{FT}, - RM <: Canopy.CanopyModel{FT, PlantHydraulics.PlantHydraulicsModel{FT}} + RM <: Canopy.CanopyModel{FT, PlantHydraulics.BigLeafHydraulicsModel{FT}} } A method which makes a function; the returned function @@ -242,7 +242,7 @@ function make_update_boundary_fluxes( FT, MM <: Soil.Biogeochemistry.SoilCO2Model{FT}, SM <: Soil.EnergyHydrology{FT}, - RM <: Canopy.CanopyModel{FT, PlantHydraulics.PlantHydraulicsModel{FT}}, + RM <: Canopy.CanopyModel{FT, PlantHydraulics.BigLeadHydraulicsModel{FT}}, } update_soil_bf! = make_update_boundary_fluxes(land.soil) update_soilco2_bf! = make_update_boundary_fluxes(land.soilco2) @@ -257,16 +257,24 @@ function make_update_boundary_fluxes( area_index = p.canopy.hydraulics.area_index + labels = land.canopy.hydraulics.parameters.h_stem > 0 ? [:stem, :leaf] : [:leaf] + above_ground_area_index = getproperty( area_index, - land.canopy.hydraulics.compartment_labels[1], + labels[1], ) + h_stem = land.canopy.hydraulics.parameters.h_stem + h_leaf = land.canopy.hydraulics.parameters.h_leaf + + midpoint = + labels[1] == :stem ? h_stem / 2 : model.h_leaf / 2 + @. p.root_extraction = (area_index.root + above_ground_area_index) / 2 * PlantHydraulics.flux( z, - land.canopy.hydraulics.compartment_midpoints[1], + midpoint, p.soil.ψ, p.canopy.hydraulics.ψ.:1, PlantHydraulics.hydraulic_conductivity( @@ -538,7 +546,7 @@ end PlantHydraulics.root_water_flux_per_ground_area!( fa::ClimaCore.Fields.Field, s::PrognosticSoil, - model::Canopy.PlantHydraulics.PlantHydraulicsModel{FT}, + model::Canopy.PlantHydraulics.BigLeafHydraulicsModel{FT}, Y::ClimaCore.Fields.FieldVector, p::NamedTuple, t, @@ -557,10 +565,7 @@ roots and soil at each soil layer. function PlantHydraulics.root_water_flux_per_ground_area!( fa::ClimaCore.Fields.Field, s::PrognosticSoil, - model::Union{ - Canopy.PlantHydraulics.PlantHydraulicsModel{FT}, - Canopy.PlantHydraulics.BigLeafHydraulicsModel{FT}, - }, + model::Canopy.PlantHydraulics.BigLeafHydraulicsModel{FT}, Y::ClimaCore.Fields.FieldVector, p::NamedTuple, t, diff --git a/src/standalone/Vegetation/Canopy.jl b/src/standalone/Vegetation/Canopy.jl index 18b3fa3a28..910772de7f 100644 --- a/src/standalone/Vegetation/Canopy.jl +++ b/src/standalone/Vegetation/Canopy.jl @@ -87,7 +87,7 @@ running the canopy model in standalone mode), or prognostic (of type PrognosticSoil, for running integrated soil+canopy models) Note that the canopy height is specified as part of the -PlantHydraulicsModel, along with the area indices of the leaves, roots, and +BigLeafHydraulicsModel, along with the area indices of the leaves, roots, and stems. Eventually, when plant biomass becomes a prognostic variable (by integrating with a carbon model), some parameters specified here will be treated differently. @@ -143,7 +143,7 @@ end An outer constructor for the `CanopyModel`. The primary constraints this applies are (1) ensuring that the domain is 1d or 2d (a ``surface" domain of a column, box, or sphere) and (2) ensuring -consistency between the PlantHydraulics model and the general canopy model, +consistency between the BigLeafHydraulics model and the general canopy model, since these are built separately. """ function CanopyModel{FT}(; @@ -362,13 +362,13 @@ end <:Union{BeerLambertModel, TwoStreamModel}, <:FarquharModel, <:MedlynConductanceModel, - <:PlantHydraulicsModel,}, + <:BigLeafHydraulicsModel,}, ) where {FT} Creates the `update_aux!` function for the `CanopyModel`; a specific method for `update_aux!` for the case where the canopy model components are of the type in the parametric type signature: `AutotrophicRespirationModel`, `AbstractRadiationModel`, -`FarquharModel`, `MedlynConductanceModel`, and `PlantHydraulicsModel`. +`FarquharModel`, `MedlynConductanceModel`, and `BigLeafHydraulicsModel`. Please note that the plant hydraulics model has auxiliary variables that are updated in its prognostic `compute_exp_tendency!` function. @@ -386,7 +386,7 @@ function ClimaLand.make_update_aux( <:Union{BeerLambertModel, TwoStreamModel}, <:Union{FarquharModel, OptimalityFarquharModel}, <:MedlynConductanceModel, - <:PlantHydraulicsModel, + <:BigLeafHydraulicsModel, <:AbstractCanopyEnergyModel, }, ) where {FT} @@ -504,7 +504,8 @@ function ClimaLand.make_update_aux( # for broadcasted expressions using the macro @. # field.:($index) .= value # works # @ field.:($$index) = value # works - @inbounds for i in 1:(n_stem + n_leaf - 1) + labels = hydraulics.parameters.h_stem > 0 ? [:stem, :leaf] : [:leaf] + @inbounds for i in eachindex(labels) ip1 = i + 1 @. ψ.:($$ip1) = PlantHydraulics.water_retention_curve( retention_model, @@ -513,16 +514,19 @@ function ClimaLand.make_update_aux( S_s, ) - areai = getproperty(area_index, hydraulics.compartment_labels[i]) + areai = getproperty(area_index, labels[i]) areaip1 = - getproperty(area_index, hydraulics.compartment_labels[ip1]) + getproperty(area_index, labels[ip1]) + + midpoint = hydraulics.h_stem / 2 + midpointip1 = hydraulics.h_stem + (model.h_leaf / 2) # Compute the flux*area between the current compartment `i` # and the compartment above. @. fa.:($$i) = PlantHydraulics.flux( - hydraulics.compartment_midpoints[i], - hydraulics.compartment_midpoints[ip1], + midpoint, + midpointip1, ψ.:($$i), ψ.:($$ip1), PlantHydraulics.hydraulic_conductivity( @@ -544,7 +548,7 @@ function ClimaLand.make_update_aux( T_canopy = canopy_temperature(canopy.energy, canopy, Y, p, t) # update moisture stress - i_end = n_stem + n_leaf + i_end = hydraulics.h_stem > 0 + 1 @. β = moisture_stress(ψ.:($$i_end) * ρ_l * grav, sc, pc) # Update Rd, An, Vcmax25 (if applicable to model) in place @@ -591,7 +595,7 @@ function make_compute_exp_tendency( <:Union{BeerLambertModel, TwoStreamModel}, <:Union{FarquharModel, OptimalityFarquharModel}, <:MedlynConductanceModel, - <:PlantHydraulicsModel, + <:BigLeafHydraulicsModel, <:Union{PrescribedCanopyTempModel, BigLeafEnergyModel}, }, ) where {FT} diff --git a/src/standalone/Vegetation/PlantHydraulics.jl b/src/standalone/Vegetation/PlantHydraulics.jl index 33f275612d..410776be42 100644 --- a/src/standalone/Vegetation/PlantHydraulics.jl +++ b/src/standalone/Vegetation/PlantHydraulics.jl @@ -21,8 +21,7 @@ import ClimaLand: auxiliary_domain_names, prognostic_domain_names, name -export PlantHydraulicsModel, - BigLeafHydraulicsModel, +export BigLeafHydraulicsModel, AbstractPlantHydraulicsModel, flux, effective_saturation, @@ -141,173 +140,6 @@ function PlantHydraulicsParameters(; ) end -""" - PlantHydraulicsModel{FT, PS, T, AA} <: AbstractPlantHydraulicsModel{FT} - -Defines, and constructs instances of, the PlantHydraulicsModel type, which is used -for simulation flux of water to/from soil, along roots of different depths, -along a stem, to a leaf, and ultimately being lost from the system by -transpiration. Note that the canopy height is specified as part of the -PlantHydraulicsModel, along with the area indices of the leaves, roots, and -stems. - -This model can also be combined with the soil model using ClimaLand, in which -case the prognostic soil water content is used to determine root extraction, and -the transpiration is also computed diagnostically. In global run with patches -of bare soil, you can "turn off" the canopy model (to get zero root extraction, zero absorption and -emission, zero transpiration and sensible heat flux from the canopy), by setting: -- n_leaf = 1 -- n_stem = 0 -- LAI = SAI = RAI = 0. - -A plant model can have leaves but no stem, but not vice versa. If n_stem = 0, SAI must be zero. - -Finally, the model can be used in Canopy standalone mode by prescribing -the soil matric potential at the root tips or flux in the roots. There is also the -option (intendend only for debugging) to use a prescribed transpiration rate. - -$(DocStringExtensions.FIELDS) -""" -struct PlantHydraulicsModel{FT, PS, T, AA <: AbstractArray{FT}} <: - AbstractPlantHydraulicsModel{FT} - "The number of stem compartments for the plant; can be zero" - n_stem::Int64 - "The number of leaf compartments for the plant; must be >=1" - n_leaf::Int64 - "The height of the center of each leaf compartment/stem compartment, in meters" - compartment_midpoints::AA - "The height of the compartments' top faces, in meters. The canopy height is the last element of the vector." - compartment_surfaces::AA - "The label (:stem or :leaf) of each compartment" - compartment_labels::Vector{Symbol} - "Parameters required by the Plant Hydraulics model" - parameters::PS - "The transpiration model, of type `AbstractTranspiration`" - transpiration::T -end - -function PlantHydraulicsModel{FT}(; - n_stem::Int64, - n_leaf::Int64, - compartment_midpoints::Vector{FT}, - compartment_surfaces::Vector{FT}, - parameters::PlantHydraulicsParameters{FT}, - transpiration::AbstractTranspiration{FT} = DiagnosticTranspiration{FT}(), -) where {FT} - args = (parameters, transpiration) - @assert (n_leaf + n_stem) == length(compartment_midpoints) - @assert (n_leaf + n_stem) + 1 == length(compartment_surfaces) - for i in 1:length(compartment_midpoints) - @assert compartment_midpoints[i] == - ((compartment_surfaces[i + 1] - compartment_surfaces[i]) / 2) + - compartment_surfaces[i] - end - compartment_labels = Vector{Symbol}(undef, n_stem + n_leaf) - for i in 1:(n_stem + n_leaf) - if i <= n_stem - compartment_labels[i] = :stem - else - compartment_labels[i] = :leaf - end - end - return PlantHydraulicsModel{ - FT, - typeof.(args)..., - typeof(compartment_midpoints), - }( - n_stem, - n_leaf, - compartment_midpoints, - compartment_surfaces, - compartment_labels, - args..., - ) -end - -""" - prognostic_vars(model::PlantHydraulicsModel) - -A function which returns the names of the prognostic -variables of the `PlantHydraulicsModel`. -""" -prognostic_vars(model::PlantHydraulicsModel) = (:ϑ_l,) - -""" - auxiliary_vars(model::PlantHydraulicsModel) - -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) - -""" - ClimaLand.prognostic_types(model::PlantHydraulicsModel{FT}) where {FT} - -Defines the prognostic types for the PlantHydraulicsModel. -""" -ClimaLand.prognostic_types(model::PlantHydraulicsModel{FT}) where {FT} = - (NTuple{model.n_stem + model.n_leaf, FT},) -ClimaLand.prognostic_domain_names(::PlantHydraulicsModel) = (:surface,) - -""" - ClimaLand.auxiliary_types(model::PlantHydraulicsModel{FT}) where {FT} - -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) - - -""" - lai_consistency_check( - model::PlantHydraulicsModel{FT}, - ) where {FT} - -Carries out consistency checks using the area indices supplied and the number of -stem elements being modeled. - -Note that it is possible to have a plant with no stem compartments -but with leaf compartments, and that a plant must have leaf compartments -(even if LAI = 0). - -Specifically, this checks that: -1. n_leaf > 0 -2. if LAI is nonzero or SAI is nonzero, RAI must be nonzero. -3. if SAI > 0, n_stem must be > 0 for consistency. If SAI == 0, n_stem must -be zero. -""" -function lai_consistency_check( - model::PlantHydraulicsModel{FT}, - area_index::NamedTuple{(:root, :stem, :leaf), Tuple{FT, FT, FT}}, -) where {FT} - n_stem = model.n_stem - n_leaf = model.n_leaf - @assert n_leaf > 0 - if area_index.leaf > eps(FT) || area_index.stem > eps(FT) - @assert area_index.root > eps(FT) - end - # If there SAI is zero, there should be no stem compartment - if area_index.stem < eps(FT) - @assert n_stem == FT(0) - else - # if SAI is > 0, n_stem should be > 0 for consistency - @assert n_stem > 0 - end - -end - """ flux( z1, @@ -490,135 +322,6 @@ function effective_saturation(ν::FT, ϑ_l::FT) where {FT} return safe_S_l # (m3 m-3) end -""" - make_compute_exp_tendency(model::PlantHydraulicsModel, _) - -A function which creates the compute_exp_tendency! function for the PlantHydraulicsModel. -The compute_exp_tendency! function must comply with a rhs function of SciMLBase.jl. - -Below, `fa` denotes a flux multiplied by the relevant cross section -(per unit area ground, or area index, AI). The tendency for the -ith compartment can be written then as: -∂ϑ[i]/∂t = 1/(AI*dz)[fa[i]-fa[i+1]). - -Note that if the area_index is zero because no plant is present, -AIdz is zero, and the fluxes `fa` appearing in the numerator are -zero because they are scaled by AI. - -To prevent dividing by zero, we change AI/(AI x dz)" to -"AI/max(AI x dz, eps(FT))" -""" -function make_compute_exp_tendency( - model::PlantHydraulicsModel{FT}, - canopy, -) where {FT} - function compute_exp_tendency!(dY, Y, p, t) - area_index = p.canopy.hydraulics.area_index - n_stem = model.n_stem - n_leaf = model.n_leaf - fa = p.canopy.hydraulics.fa - fa_roots = p.canopy.hydraulics.fa_roots - - # Inside of a loop, we need to use a single dollar sign - # for indexing into Fields of Tuples in non broadcasted - # expressions, and two dollar signs for - # for broadcasted expressions using the macro @. - # field.:($index) .= value # works - # @ field.:($$index) = value # works - @inbounds for i in 1:(n_stem + n_leaf) - im1 = i - 1 - ip1 = i + 1 - # To prevent dividing by zero, change AI/(AI x dz)" to - # "AI/max(AI x dz, eps(FT))" - AIdz = - max.( - getproperty(area_index, model.compartment_labels[i]) .* ( - model.compartment_surfaces[ip1] - - model.compartment_surfaces[i] - ), - eps(FT), - ) - if i == 1 - @inbounds @. dY.canopy.hydraulics.ϑ_l.:($$i) = - 1 / AIdz * (fa_roots - fa.:($$i)) - else - @inbounds @. dY.canopy.hydraulics.ϑ_l.:($$i) = - 1 / AIdz * (fa.:($$im1) - fa.:($$i)) - end - end - end - return compute_exp_tendency! -end - -""" - root_water_flux_per_ground_area!( - fa::ClimaCore.Fields.Field, - s::PrescribedSoil, - model::PlantHydraulicsModel{FT}, - Y::ClimaCore.Fields.FieldVector, - p::NamedTuple, - t, - ) where {FT} - -A method which computes the water flux between the soil and the stem, via the roots, -and multiplied by the RAI, in the case of a model running without an integrated -soil model. - -The returned flux is per unit ground area. This assumes that the stem compartment -is the first element of `Y.canopy.hydraulics.ϑ_l`. -""" -function root_water_flux_per_ground_area!( - fa::ClimaCore.Fields.Field, - s::PrescribedSoil, - model::PlantHydraulicsModel{FT}, - Y::ClimaCore.Fields.FieldVector, - p::NamedTuple, - t, -) where {FT} - - (; conductivity_model, root_distribution) = model.parameters - area_index = p.canopy.hydraulics.area_index - # We can index into a field of Tuple{FT} to extract a field of FT - # using the following notation: field.:index - ψ_base = p.canopy.hydraulics.ψ.:1 - root_depths = s.root_depths - n_root_layers = length(root_depths) - ψ_soil::FT = s.ψ(t) - fa .= FT(0.0) - @inbounds for i in 1:n_root_layers - above_ground_area_index = - getproperty(area_index, model.compartment_labels[1]) - - if i != n_root_layers - @. fa += - flux( - root_depths[i], - model.compartment_midpoints[1], - ψ_soil, - ψ_base, - hydraulic_conductivity(conductivity_model, ψ_soil), - hydraulic_conductivity(conductivity_model, ψ_base), - ) * - root_distribution(root_depths[i]) * - (root_depths[i + 1] - root_depths[i]) * - (area_index.root + above_ground_area_index) / 2 - else - @. fa += - flux( - root_depths[i], - model.compartment_midpoints[1], - ψ_soil, - ψ_base, - hydraulic_conductivity(conductivity_model, ψ_soil), - hydraulic_conductivity(conductivity_model, ψ_base), - ) * - root_distribution(root_depths[i]) * - (model.compartment_surfaces[1] - root_depths[n_root_layers]) * - (area_index.root + above_ground_area_index) / 2 - end - end -end - """ PrescribedTranspiration{FT, F <: Function} <: AbstractTranspiration{FT} @@ -807,7 +510,6 @@ function make_compute_exp_tendency( @inbounds for i in eachindex(labels) im1 = i - 1 - ip1 = i + 1 compartment_type = labels[i] dz = compartment_type == :stem ? model.h_stem : @@ -912,9 +614,8 @@ but with a leaf compartment, and that a plant must have leaf compartments (even if LAI = 0). Specifically, this checks that: -1. n_leaf > 0 -2. if LAI is nonzero or SAI is nonzero, RAI must be nonzero. -3. if SAI > 0, n_stem must be > 0 for consistency. If SAI == 0, n_stem must +1. if LAI is nonzero or SAI is nonzero, RAI must be nonzero. +2. if SAI > 0, h_stem must be > 0 for consistency. If SAI == 0, h_stem must be zero. """ function lai_consistency_check( @@ -928,7 +629,7 @@ function lai_consistency_check( if area_index.stem < eps(FT) @assert model.h_stem == FT(0) else - # if SAI is > 0, n_stem should be > 0 for consistency + # if SAI is > 0, h_stem should be > 0 for consistency @assert model.h_stem > 0 end end @@ -946,7 +647,7 @@ Sets the canopy prescribed fields pertaining to the PlantHydraulics component (the area indices) with their values at time t. """ function ClimaLand.Canopy.set_canopy_prescribed_field!( - component::Union{PlantHydraulicsModel{FT}, BigLeafHydraulicsModel{FT}}, + component::BigLeafHydraulicsModel{FT}, p, t, ) where {FT} diff --git a/src/standalone/Vegetation/canopy_boundary_fluxes.jl b/src/standalone/Vegetation/canopy_boundary_fluxes.jl index 5199989ddf..673273e0fd 100644 --- a/src/standalone/Vegetation/canopy_boundary_fluxes.jl +++ b/src/standalone/Vegetation/canopy_boundary_fluxes.jl @@ -110,7 +110,7 @@ end <:Union{BeerLambertModel, TwoStreamModel}, <:Union{FarquharModel,OptimalityFarquharModel}, <:MedlynConductanceModel, - <:PlantHydraulicsModel, + <:BigLeafHydraulicsModel, <:Union{PrescribedCanopyTempModel,BigLeafEnergyModel} }, radiation::PrescribedRadiativeFluxes, @@ -141,7 +141,7 @@ function canopy_boundary_fluxes!( <:Union{BeerLambertModel, TwoStreamModel}, <:Union{FarquharModel, OptimalityFarquharModel}, <:MedlynConductanceModel, - <:PlantHydraulicsModel, + <:BigLeafHydraulicsModel, <:Union{PrescribedCanopyTempModel, BigLeafEnergyModel}, }, radiation::PrescribedRadiativeFluxes, @@ -157,7 +157,7 @@ function canopy_boundary_fluxes!( shf = p.canopy.energy.shf lhf = p.canopy.energy.lhf r_ae = p.canopy.energy.r_ae - i_end = canopy.hydraulics.n_stem + canopy.hydraulics.n_leaf + i_end = canopy.hydraulics.h_stem > 0 + 1 # Compute transpiration, SHF, LHF canopy_turbulent_fluxes = turbulent_fluxes(atmos, canopy, Y, p, t) diff --git a/test/standalone/Vegetation/canopy_model.jl b/test/standalone/Vegetation/canopy_model.jl index bc2e8958d9..006e4c601d 100644 --- a/test/standalone/Vegetation/canopy_model.jl +++ b/test/standalone/Vegetation/canopy_model.jl @@ -143,38 +143,15 @@ for FT in (Float32, Float64) conductivity_model = conductivity_model, retention_model = retention_model, ) - Δz = FT(1.0) # height of compartments - n_stem = Int64(0) # number of stem elements - n_leaf = Int64(1) # number of leaf elements - compartment_centers = - FT.( - Vector( - range( - start = Δz / 2, - step = Δz, - stop = Δz * (n_stem + n_leaf) - (Δz / 2), - ), - ), - ) - compartment_faces = - FT.( - Vector( - range( - start = 0.0, - step = Δz, - stop = Δz * (n_stem + n_leaf), - ), - ) - ) + h_stem = FT(0.0) + h_leaf = FT(1.0) ψ_soil0 = FT(0.0) soil_driver = PrescribedSoil(FT) - plant_hydraulics = PlantHydraulics.PlantHydraulicsModel{FT}(; + plant_hydraulics = PlantHydraulics.BigLeafHydraulicsModel{FT}(; parameters = param_set, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_surfaces = compartment_faces, - compartment_midpoints = compartment_centers, + h_stem = h_stem, + h_leaf = h_leaf, ) canopy = ClimaLand.Canopy.CanopyModel{FT}(; parameters = shared_params, @@ -332,7 +309,7 @@ for FT in (Float32, Float64) ) < 0.5 @test ClimaLand.surface_evaporative_scaling(canopy, Y, p) == FT(1.0) - @test ClimaLand.surface_height(canopy, Y, p) == compartment_faces[1] + @test ClimaLand.surface_height(canopy, Y, p) == h_stem T_sfc = FT.(T_atmos(t0)) @test Array(parent(ClimaLand.surface_temperature(canopy, Y, p, t0))) == [T_sfc] @@ -403,36 +380,13 @@ for FT in (Float32, Float64) conductivity_model = conductivity_model, retention_model = retention_model, ) - Δz = FT(1.0) # height of compartments - n_stem = Int64(0) # number of stem elements - n_leaf = Int64(1) # number of leaf elements - compartment_centers = - FT.( - Vector( - range( - start = Δz / 2, - step = Δz, - stop = Δz * (n_stem + n_leaf) - (Δz / 2), - ), - ), - ) - compartment_faces = - FT.( - Vector( - range( - start = 0.0, - step = Δz, - stop = Δz * (n_stem + n_leaf), - ), - ) - ) - - plant_hydraulics = PlantHydraulics.PlantHydraulicsModel{FT}(; + h_stem = FT(0.0) + h_leaf = FT(1.0) + + plant_hydraulics = PlantHydraulics.BigLeafHydraulicsModel{FT}(; parameters = param_set, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_surfaces = compartment_faces, - compartment_midpoints = compartment_centers, + h_stem = h_stem, + h_leaf = h_leaf, ) t0 = FT(100) @@ -621,29 +575,8 @@ for FT in (Float32, Float64) conductivity_model = conductivity_model, retention_model = retention_model, ) - Δz = FT(1.0) # height of compartments - n_stem = Int64(0) # number of stem elements - n_leaf = Int64(1) # number of leaf elements - compartment_centers = - FT.( - Vector( - range( - start = Δz / 2, - step = Δz, - stop = Δz * (n_stem + n_leaf) - (Δz / 2), - ), - ), - ) - compartment_faces = - FT.( - Vector( - range( - start = 0.0, - step = Δz, - stop = Δz * (n_stem + n_leaf), - ), - ) - ) + h_stem = FT(0.0) + h_leaf = FT(1.0) ψ_soil0 = FT(0.0) T_soil0 = FT(290) @@ -656,12 +589,10 @@ for FT in (Float32, Float64) FT(0.98), ) - plant_hydraulics = PlantHydraulics.PlantHydraulicsModel{FT}(; + plant_hydraulics = PlantHydraulics.BigLeafHydraulicsModel{FT}(; parameters = param_set, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_surfaces = compartment_faces, - compartment_midpoints = compartment_centers, + h_stem = h_stem, + h_leaf = h_leaf, ) autotrophic_parameters = AutotrophicRespirationParameters(FT) autotrophic_respiration_model = diff --git a/test/standalone/Vegetation/plant_hydraulics_test.jl b/test/standalone/Vegetation/plant_hydraulics_test.jl index 341aff58d4..267669dce8 100644 --- a/test/standalone/Vegetation/plant_hydraulics_test.jl +++ b/test/standalone/Vegetation/plant_hydraulics_test.jl @@ -20,42 +20,42 @@ for FT in (Float32, Float64) LAI = FT.([0, 1]) SAI = FT.([0, 1]) RAI = FT.([0, 1]) - n_stem = [0, 1] - n_leaf = [1] + h_stem = FT.([0, 1]) + h_leaf = FT.([1]) for L in LAI for S in SAI for R in RAI - for n_s in n_stem - for n_l in n_leaf + for h_s in h_stem + for h_l in h_leaf area_index = (root = R, stem = S, leaf = L) if n_l == 0 @test_throws AssertionError ClimaLand.Canopy.PlantHydraulics.lai_consistency_check( - n_s, - n_l, + h_s, + h_l, area_index, ) end if (L > eps(FT) || S > eps(FT)) && R < eps(FT) @test_throws AssertionError ClimaLand.Canopy.PlantHydraulics.lai_consistency_check( - n_s, - n_l, + h_s, + h_l, area_index, ) end if S > FT(0) && n_s == 0 @test_throws AssertionError ClimaLand.Canopy.PlantHydraulics.lai_consistency_check( - n_s, - n_l, + h_s, + h_l, area_index, ) end if S < eps(FT) && n_s > 0 @test_throws AssertionError ClimaLand.Canopy.PlantHydraulics.lai_consistency_check( - n_s, - n_l, + h_s, + h_l, area_index, ) end @@ -101,309 +101,6 @@ for FT in (Float32, Float64) ) end - @testset "Plant hydraulics model integration tests, FT = $FT" begin - domains = [ - Point(; z_sfc = FT(0.0)), - Plane(; - xlim = FT.((0, 1)), - ylim = FT.((0, 1)), - nelements = (2, 2), - periodic = (true, true), - npolynomial = 1, - ), - ] - - AR_params = AutotrophicRespirationParameters(FT) - RTparams = BeerLambertParameters(FT) - photosynthesis_params = FarquharParameters(FT, C3()) - stomatal_g_params = MedlynConductanceParameters(FT) - - AR_model = AutotrophicRespirationModel{FT}(AR_params) - stomatal_model = MedlynConductanceModel{FT}(stomatal_g_params) - photosynthesis_model = FarquharModel{FT}(photosynthesis_params) - rt_model = BeerLambertModel{FT}(RTparams) - - earth_param_set = LP.LandParameters(FT) - LAI = (t) -> 1.0 # m2 [leaf] m-2 [ground] - z_0m = FT(2.0) # m, Roughness length for momentum - z_0b = FT(0.1) # m, Roughness length for scalars - h_canopy = FT(20.0) # m, canopy height - h_sfc = FT(20.0) # m, canopy height - h_int = FT(30.0) # m, "where measurements would be taken at a typical flux tower of a 20m canopy" - shared_params = SharedCanopyParameters{FT, typeof(earth_param_set)}( - z_0m, - z_0b, - earth_param_set, - ) - lat = FT(0.0) # degree - long = FT(-180) # degree - - function zenith_angle( - t, - ref_time; - latitude = lat, - longitude = long, - insol_params = earth_param_set.insol_params, - ) - current_datetime = ref_time + Dates.Second(round(t)) - d, δ, η_UTC = - FT.( - Insolation.helper_instantaneous_zenith_angle( - current_datetime, - ref_time, - insol_params, - ) - ) - return Insolation.instantaneous_zenith_angle( - d, - δ, - η_UTC, - longitude, - latitude, - )[1] - end - - function shortwave_radiation( - t; - latitude = lat, - longitude = long, - insol_params = earth_param_set.insol_params, - ) - return 1000 # W/m^2 - end - - function longwave_radiation(t) - return 200 # W/m^2 - end - - u_atmos = t -> 10 #m.s-1 - - liquid_precip = (t) -> 0 # m - snow_precip = (t) -> 0 # m - T_atmos = t -> 290 # Kelvin - q_atmos = t -> 0.001 # kg/kg - P_atmos = t -> 1e5 # Pa - h_atmos = h_int # m - c_atmos = (t) -> 4.11e-4 # mol/mol - ref_time = DateTime(2005) - atmos = PrescribedAtmosphere( - TimeVaryingInput(liquid_precip), - TimeVaryingInput(snow_precip), - TimeVaryingInput(T_atmos), - TimeVaryingInput(u_atmos), - TimeVaryingInput(q_atmos), - TimeVaryingInput(P_atmos), - ref_time, - h_atmos, - earth_param_set; - c_co2 = TimeVaryingInput(c_atmos), - ) - radiation = PrescribedRadiativeFluxes( - FT, - TimeVaryingInput(shortwave_radiation), - TimeVaryingInput(longwave_radiation), - ref_time; - θs = zenith_angle, - ) - Δz = FT(1.0) # height of compartments - n_stem = Int64(5) # number of stem elements - n_leaf = Int64(6) # number of leaf elements - SAI = FT(1) # m2/m2 - RAI = FT(1) # m2/m2 - ai_parameterization = PlantHydraulics.PrescribedSiteAreaIndex{FT}( - TimeVaryingInput(LAI), - SAI, - RAI, - ) - K_sat_plant = 1.8e-8 # m/s. - ψ63 = FT(-4 / 0.0098) # / MPa to m - Weibull_param = FT(4) # unitless - a = FT(0.05 * 0.0098) # 1/m - plant_ν = FT(0.7) # m3/m3 - plant_S_s = FT(1e-2 * 0.0098) # m3/m3/MPa to m3/m3/m - conductivity_model = - PlantHydraulics.Weibull{FT}(K_sat_plant, ψ63, Weibull_param) - retention_model = PlantHydraulics.LinearRetentionCurve{FT}(a) - root_depths = FT.(-Array{FT}(10:-1:1.0) ./ 10.0 * 2.0 .+ 0.2 / 2.0) # 1st element is the deepest root depth - function root_distribution(z::T) where {T} - return T(1.0 / 0.5) * exp(z / T(0.5)) # (1/m) - end - compartment_midpoints = Vector{FT}( - range( - start = Δz / 2, - step = Δz, - stop = Δz * (n_stem + n_leaf) - (Δz / 2), - ), - ) - - compartment_surfaces = Vector{FT}( - range(start = 0.0, step = Δz, stop = Δz * (n_stem + n_leaf)), - ) - - param_set = PlantHydraulics.PlantHydraulicsParameters(; - ai_parameterization = ai_parameterization, - ν = plant_ν, - S_s = plant_S_s, - root_distribution = root_distribution, - conductivity_model = conductivity_model, - retention_model = retention_model, - ) - - function leaf_transpiration(t) - T = FT(1e-8) # m/s - end - - ψ_soil0 = FT(0.0) - transpiration = PrescribedTranspiration{FT}(leaf_transpiration) - - soil_driver = PrescribedSoil(FT) - - plant_hydraulics = PlantHydraulics.PlantHydraulicsModel{FT}(; - parameters = param_set, - transpiration = transpiration, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_surfaces = compartment_surfaces, - compartment_midpoints = compartment_midpoints, - ) - autotrophic_parameters = AutotrophicRespirationParameters(FT) - autotrophic_respiration_model = - AutotrophicRespirationModel{FT}(autotrophic_parameters) - for domain in domains - model = ClimaLand.Canopy.CanopyModel{FT}(; - parameters = shared_params, - domain = domain, - autotrophic_respiration = autotrophic_respiration_model, - radiative_transfer = rt_model, - photosynthesis = photosynthesis_model, - conductance = stomatal_model, - hydraulics = plant_hydraulics, - soil_driver = soil_driver, - atmos = atmos, - radiation = radiation, - ) - # Set system to hydrostatic equilibrium state by setting fluxes to zero, and setting LHS of both ODEs to 0 - function initial_compute_exp_tendency!(F, Y) - AI = (; leaf = LAI(0.0), root = RAI, stem = SAI) - T0A = FT(1e-8) * AI[:leaf] - for i in 1:(n_leaf + n_stem) - if i == 1 - fa = - sum( - flux.( - root_depths, - plant_hydraulics.compartment_midpoints[i], - ψ_soil0, - Y[i], - PlantHydraulics.hydraulic_conductivity( - conductivity_model, - ψ_soil0, - ), - PlantHydraulics.hydraulic_conductivity( - conductivity_model, - Y[i], - ), - ) .* root_distribution.(root_depths) .* ( - vcat(root_depths, [0.0])[2:end] - - vcat(root_depths, [0.0])[1:(end - 1)] - ), - ) * AI[:root] - else - fa = - flux( - plant_hydraulics.compartment_midpoints[i - 1], - plant_hydraulics.compartment_midpoints[i], - Y[i - 1], - Y[i], - PlantHydraulics.hydraulic_conductivity( - conductivity_model, - ψ_soil0, - ), - PlantHydraulics.hydraulic_conductivity( - conductivity_model, - Y[i], - ), - ) * ( - AI[plant_hydraulics.compartment_labels[1]] + - AI[plant_hydraulics.compartment_labels[2]] - ) / 2 - end - F[i] = fa - T0A - end - end - - soln = nlsolve( - initial_compute_exp_tendency!, - Vector{FT}(-0.03:0.01:0.07); - ftol = 1e-11, - ) - - S_l = - inverse_water_retention_curve.( - retention_model, - soln.zero, - plant_ν, - plant_S_s, - ) - - ϑ_l_0 = augmented_liquid_fraction.(plant_ν, S_l) - - Y, p, coords = initialize(model) - if typeof(domain) <: ClimaLand.Domains.Point - @test propertynames(p) == (:canopy, :drivers) - elseif typeof(domain) <: ClimaLand.Domains.Plane - @test propertynames(p) == (:canopy, :dss_buffer_2d, :drivers) - @test typeof(p.dss_buffer_2d) == typeof( - ClimaCore.Spaces.create_dss_buffer( - ClimaCore.Fields.zeros(domain.space.surface), - ), - ) - end - - dY = similar(Y) - for i in 1:(n_stem + n_leaf) - Y.canopy.hydraulics.ϑ_l.:($i) .= ϑ_l_0[i] - p.canopy.hydraulics.ψ.:($i) .= NaN - p.canopy.hydraulics.fa.:($i) .= NaN - dY.canopy.hydraulics.ϑ_l.:($i) .= NaN - end - set_initial_cache! = make_set_initial_cache(model) - set_initial_cache!(p, Y, 0.0) - canopy_exp_tendency! = make_exp_tendency(model) - canopy_exp_tendency!(dY, Y, p, 0.0) - - m = similar(dY.canopy.hydraulics.ϑ_l.:1) - m .= FT(0) - for i in 1:(n_stem + n_leaf) - @. m += sqrt(dY.canopy.hydraulics.ϑ_l.:($$i)^2.0) - end - @test maximum(parent(m)) < 200 * eps(FT) # starts in equilibrium - - - # repeat using the plant hydraulics model directly - # make sure it agrees with what we get when use the canopy model ODE - Y, p, coords = initialize(model) - standalone_dY = similar(Y) - for i in 1:(n_stem + n_leaf) - Y.canopy.hydraulics.ϑ_l.:($i) .= ϑ_l_0[i] - p.canopy.hydraulics.ψ.:($i) .= NaN - p.canopy.hydraulics.fa.:($i) .= NaN - standalone_dY.canopy.hydraulics.ϑ_l.:($i) .= NaN - end - set_initial_cache!(p, Y, 0.0) - standalone_exp_tendency! = - make_compute_exp_tendency(model.hydraulics, model) - standalone_exp_tendency!(standalone_dY, Y, p, 0.0) - - m = similar(dY.canopy.hydraulics.ϑ_l.:1) - m .= FT(0) - for i in 1:(n_stem + n_leaf) - @. m += sqrt(dY.canopy.hydraulics.ϑ_l.:($$i)^2.0) - end - @test maximum(parent(m)) < 10^2.1 * eps(FT) # starts in equilibrium - end - end - @testset "No plant, FT = $FT" begin domain = Point(; z_sfc = FT(0.0)) @@ -500,8 +197,8 @@ for FT in (Float32, Float64) θs = zenith_angle, ) - n_stem = Int64(0) # number of stem elements - n_leaf = Int64(1) # number of leaf elements + h_stem = FT(0.0) # m + h_leaf = FT(0.0) # m SAI = FT(0) # m2/m2 RAI = FT(0) # m2/m2 lai_fun = t -> 0 @@ -523,8 +220,6 @@ for FT in (Float32, Float64) function root_distribution(z::T) where {T} return T(0) # (1/m) end - compartment_midpoints = [h_canopy] - compartment_surfaces = [FT(0.0), h_canopy] param_set = PlantHydraulics.PlantHydraulicsParameters(; ai_parameterization = ai_parameterization, @@ -537,13 +232,11 @@ for FT in (Float32, Float64) transpiration = DiagnosticTranspiration{FT}() soil_driver = PrescribedSoil(FT) - plant_hydraulics = PlantHydraulics.PlantHydraulicsModel{FT}(; + plant_hydraulics = PlantHydraulics.BigLeafHydraulicsModel{FT}(; parameters = param_set, transpiration = transpiration, - n_stem = n_stem, - n_leaf = n_leaf, - compartment_surfaces = compartment_surfaces, - compartment_midpoints = compartment_midpoints, + h_stem = h_stem, + h_leaf = h_leaf, ) autotrophic_parameters = AutotrophicRespirationParameters(FT) @@ -565,7 +258,7 @@ for FT in (Float32, Float64) Y, p, coords = initialize(model) dY = similar(Y) - for i in 1:(n_stem + n_leaf) + for i in 1:(h_stem != 0 + 1) Y.canopy.hydraulics.ϑ_l.:($i) .= FT(0.1) p.canopy.hydraulics.ψ.:($i) .= NaN p.canopy.hydraulics.fa.:($i) .= NaN From ecdd138ab5784aa9ae0dec2d6799bc94170eb59d Mon Sep 17 00:00:00 2001 From: Edward Speer Date: Sun, 2 Jun 2024 00:35:29 -0700 Subject: [PATCH 06/10] reformat --- experiments/integrated/fluxnet/ozark_pft.jl | 7 ++----- experiments/integrated/fluxnet/run_fluxnet.jl | 7 ++----- .../conservation/ozark_conservation_setup.jl | 7 ++----- .../integrated/performance/profile_allocations.jl | 7 ++----- .../src/Fluxnet/fluxnet_utilities/make_setup.jl | 10 +--------- .../src/Fluxnet/run_fluxnet.jl | 6 +----- src/integrated/soil_canopy_model.jl | 14 ++++++-------- src/standalone/Vegetation/Canopy.jl | 3 +-- src/standalone/Vegetation/PlantHydraulics.jl | 2 +- test/standalone/Vegetation/canopy_model.jl | 2 +- 10 files changed, 19 insertions(+), 46 deletions(-) diff --git a/experiments/integrated/fluxnet/ozark_pft.jl b/experiments/integrated/fluxnet/ozark_pft.jl index 66accb4e81..baae661d20 100644 --- a/experiments/integrated/fluxnet/ozark_pft.jl +++ b/experiments/integrated/fluxnet/ozark_pft.jl @@ -209,11 +209,8 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; conductivity_model = conductivity_model, retention_model = retention_model, ) -plant_hydraulics_args = ( - parameters = plant_hydraulics_ps, - h_stem = h_stem, - h_leaf = h_leaf, -) +plant_hydraulics_args = + (parameters = plant_hydraulics_ps, h_stem = h_stem, h_leaf = h_leaf) energy_args = (parameters = Canopy.BigLeafEnergyParameters{FT}(ac_canopy),) diff --git a/experiments/integrated/fluxnet/run_fluxnet.jl b/experiments/integrated/fluxnet/run_fluxnet.jl index 82c6de378b..69cc9a60dd 100644 --- a/experiments/integrated/fluxnet/run_fluxnet.jl +++ b/experiments/integrated/fluxnet/run_fluxnet.jl @@ -169,11 +169,8 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; conductivity_model = conductivity_model, retention_model = retention_model, ) -plant_hydraulics_args = ( - parameters = plant_hydraulics_ps, - h_stem = h_stem, - h_leaf = h_leaf, -) +plant_hydraulics_args = + (parameters = plant_hydraulics_ps, h_stem = h_stem, h_leaf = h_leaf) energy_args = (parameters = Canopy.BigLeafEnergyParameters{FT}(ac_canopy),) diff --git a/experiments/integrated/performance/conservation/ozark_conservation_setup.jl b/experiments/integrated/performance/conservation/ozark_conservation_setup.jl index 3c308abcd1..df788f1b18 100644 --- a/experiments/integrated/performance/conservation/ozark_conservation_setup.jl +++ b/experiments/integrated/performance/conservation/ozark_conservation_setup.jl @@ -151,11 +151,8 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; conductivity_model = conductivity_model, retention_model = retention_model, ) -plant_hydraulics_args = ( - parameters = plant_hydraulics_ps, - h_stem = h_stem, - h_leaf = h_leaf, -) +plant_hydraulics_args = + (parameters = plant_hydraulics_ps, h_stem = h_stem, h_leaf = h_leaf) # Canopy component args canopy_component_args = (; diff --git a/experiments/integrated/performance/profile_allocations.jl b/experiments/integrated/performance/profile_allocations.jl index 0506c05e6b..214cd40227 100644 --- a/experiments/integrated/performance/profile_allocations.jl +++ b/experiments/integrated/performance/profile_allocations.jl @@ -283,11 +283,8 @@ plant_hydraulics_ps = PlantHydraulics.PlantHydraulicsParameters(; conductivity_model = conductivity_model, retention_model = retention_model, ) -plant_hydraulics_args = ( - parameters = plant_hydraulics_ps, - h_stem = h_stem, - h_leaf = h_leaf, -) +plant_hydraulics_args = + (parameters = plant_hydraulics_ps, h_stem = h_stem, h_leaf = h_leaf) # Canopy component args canopy_component_args = (; diff --git a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_utilities/make_setup.jl b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_utilities/make_setup.jl index 9e9434004f..9fb3738ecd 100644 --- a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_utilities/make_setup.jl +++ b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_utilities/make_setup.jl @@ -12,15 +12,7 @@ setup site specific domain and time stepping: - dt, the time step (s) - n, the number of time step between saving outputs """ -function make_setup(; - dz_bottom, - dz_top, - h_stem, - h_leaf, - t0, - dt, - n, -) +function make_setup(; dz_bottom, dz_top, h_stem, h_leaf, t0, dt, n) return ( dz_bottom = dz_bottom, dz_top = dz_top, diff --git a/lib/ClimaLandSimulations/src/Fluxnet/run_fluxnet.jl b/lib/ClimaLandSimulations/src/Fluxnet/run_fluxnet.jl index d27c702fc1..b49356f4a3 100644 --- a/lib/ClimaLandSimulations/src/Fluxnet/run_fluxnet.jl +++ b/lib/ClimaLandSimulations/src/Fluxnet/run_fluxnet.jl @@ -180,11 +180,7 @@ function run_fluxnet( conductivity_model = params.plant_hydraulics.conductivity_model, retention_model = params.plant_hydraulics.retention_model, ) - plant_hydraulics_args = ( - parameters = plant_hydraulics_ps, - h_stem, - h_leaf, - ) + plant_hydraulics_args = (parameters = plant_hydraulics_ps, h_stem, h_leaf) energy_args = ( parameters = Canopy.BigLeafEnergyParameters{FT}( diff --git a/src/integrated/soil_canopy_model.jl b/src/integrated/soil_canopy_model.jl index c86b43d39a..20877e8527 100644 --- a/src/integrated/soil_canopy_model.jl +++ b/src/integrated/soil_canopy_model.jl @@ -242,7 +242,7 @@ function make_update_boundary_fluxes( FT, MM <: Soil.Biogeochemistry.SoilCO2Model{FT}, SM <: Soil.EnergyHydrology{FT}, - RM <: Canopy.CanopyModel{FT, PlantHydraulics.BigLeadHydraulicsModel{FT}}, + RM <: Canopy.CanopyModel{FT, PlantHydraulics.BigLeafHydraulicsModel{FT}}, } update_soil_bf! = make_update_boundary_fluxes(land.soil) update_soilco2_bf! = make_update_boundary_fluxes(land.soilco2) @@ -257,18 +257,16 @@ function make_update_boundary_fluxes( area_index = p.canopy.hydraulics.area_index - labels = land.canopy.hydraulics.parameters.h_stem > 0 ? [:stem, :leaf] : [:leaf] + labels = + land.canopy.hydraulics.parameters.h_stem > 0 ? [:stem, :leaf] : + [:leaf] - above_ground_area_index = getproperty( - area_index, - labels[1], - ) + above_ground_area_index = getproperty(area_index, labels[1]) h_stem = land.canopy.hydraulics.parameters.h_stem h_leaf = land.canopy.hydraulics.parameters.h_leaf - midpoint = - labels[1] == :stem ? h_stem / 2 : model.h_leaf / 2 + midpoint = labels[1] == :stem ? h_stem / 2 : model.h_leaf / 2 @. p.root_extraction = (area_index.root + above_ground_area_index) / 2 * diff --git a/src/standalone/Vegetation/Canopy.jl b/src/standalone/Vegetation/Canopy.jl index 910772de7f..c42dd619d8 100644 --- a/src/standalone/Vegetation/Canopy.jl +++ b/src/standalone/Vegetation/Canopy.jl @@ -515,8 +515,7 @@ function ClimaLand.make_update_aux( ) areai = getproperty(area_index, labels[i]) - areaip1 = - getproperty(area_index, labels[ip1]) + areaip1 = getproperty(area_index, labels[ip1]) midpoint = hydraulics.h_stem / 2 midpointip1 = hydraulics.h_stem + (model.h_leaf / 2) diff --git a/src/standalone/Vegetation/PlantHydraulics.jl b/src/standalone/Vegetation/PlantHydraulics.jl index 410776be42..acc0a80ec5 100644 --- a/src/standalone/Vegetation/PlantHydraulics.jl +++ b/src/standalone/Vegetation/PlantHydraulics.jl @@ -622,7 +622,7 @@ function lai_consistency_check( model::BigLeafHydraulicsModel{FT}, area_index::NamedTuple{(:root, :stem, :leaf), Tuple{FT, FT, FT}}, ) where {FT} - if area_index.leaf > esp(FT) || area_index.stem > eso(FT) + if area_index.leaf > eps(FT) || area_index.stem > eps(FT) @assert area_index.root > eps(FT) end # If the SAI is 0, then the stem compartment should be height 0 diff --git a/test/standalone/Vegetation/canopy_model.jl b/test/standalone/Vegetation/canopy_model.jl index 006e4c601d..b0d9c2d065 100644 --- a/test/standalone/Vegetation/canopy_model.jl +++ b/test/standalone/Vegetation/canopy_model.jl @@ -382,7 +382,7 @@ for FT in (Float32, Float64) ) h_stem = FT(0.0) h_leaf = FT(1.0) - + plant_hydraulics = PlantHydraulics.BigLeafHydraulicsModel{FT}(; parameters = param_set, h_stem = h_stem, From 899fe4269c5a7357a2c6f03350f013c319c16bb3 Mon Sep 17 00:00:00 2001 From: Edward Speer Date: Sun, 2 Jun 2024 01:23:00 -0700 Subject: [PATCH 07/10] h_stem not in parameters --- src/integrated/soil_canopy_model.jl | 6 +++--- src/standalone/Vegetation/Canopy.jl | 2 +- src/standalone/Vegetation/PlantHydraulics.jl | 4 ++-- 3 files changed, 6 insertions(+), 6 deletions(-) diff --git a/src/integrated/soil_canopy_model.jl b/src/integrated/soil_canopy_model.jl index 20877e8527..20b1cd506d 100644 --- a/src/integrated/soil_canopy_model.jl +++ b/src/integrated/soil_canopy_model.jl @@ -258,13 +258,13 @@ function make_update_boundary_fluxes( area_index = p.canopy.hydraulics.area_index labels = - land.canopy.hydraulics.parameters.h_stem > 0 ? [:stem, :leaf] : + land.canopy.hydraulics.h_stem > 0 ? [:stem, :leaf] : [:leaf] above_ground_area_index = getproperty(area_index, labels[1]) - h_stem = land.canopy.hydraulics.parameters.h_stem - h_leaf = land.canopy.hydraulics.parameters.h_leaf + h_stem = land.canopy.hydraulics.h_stem + h_leaf = land.canopy.hydraulics.h_leaf midpoint = labels[1] == :stem ? h_stem / 2 : model.h_leaf / 2 diff --git a/src/standalone/Vegetation/Canopy.jl b/src/standalone/Vegetation/Canopy.jl index c42dd619d8..f057ab69ba 100644 --- a/src/standalone/Vegetation/Canopy.jl +++ b/src/standalone/Vegetation/Canopy.jl @@ -504,7 +504,7 @@ function ClimaLand.make_update_aux( # for broadcasted expressions using the macro @. # field.:($index) .= value # works # @ field.:($$index) = value # works - labels = hydraulics.parameters.h_stem > 0 ? [:stem, :leaf] : [:leaf] + labels = hydraulics.h_stem > 0 ? [:stem, :leaf] : [:leaf] @inbounds for i in eachindex(labels) ip1 = i + 1 @. ψ.:($$ip1) = PlantHydraulics.water_retention_curve( diff --git a/src/standalone/Vegetation/PlantHydraulics.jl b/src/standalone/Vegetation/PlantHydraulics.jl index acc0a80ec5..1c06a82a4d 100644 --- a/src/standalone/Vegetation/PlantHydraulics.jl +++ b/src/standalone/Vegetation/PlantHydraulics.jl @@ -506,7 +506,7 @@ function make_compute_exp_tendency( area_index = p.canopy.hydraulics.area_index fa = p.canopy.hydraulics.fa fa_roots = p.canopy.hydraulics.fa_roots - labels = model.parameters.h_stem > 0 ? [:stem, :leaf] : [:leaf] + labels = model.h_stem > 0 ? [:stem, :leaf] : [:leaf] @inbounds for i in eachindex(labels) im1 = i - 1 @@ -563,7 +563,7 @@ function root_water_flux_per_ground_area!( n_root_layers = length(root_depths) ψ_soil::FT = s.ψ(t) fa .= FT(0.0) - labels = model.parameters.h_stem > 0 ? [:stem, :leaf] : [:leaf] + labels = model.h_stem > 0 ? [:stem, :leaf] : [:leaf] @inbounds for i in 1:n_root_layers above_ground_area_index = getproperty(area_index, labels[1]) compartment = labels[i] From c9ab29bef534c27f086ffac362efe9634817dbeb Mon Sep 17 00:00:00 2001 From: Edward Speer Date: Sun, 2 Jun 2024 01:44:36 -0700 Subject: [PATCH 08/10] fix loop index bug --- src/standalone/Vegetation/Canopy.jl | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/standalone/Vegetation/Canopy.jl b/src/standalone/Vegetation/Canopy.jl index f057ab69ba..11ce5c8950 100644 --- a/src/standalone/Vegetation/Canopy.jl +++ b/src/standalone/Vegetation/Canopy.jl @@ -505,8 +505,9 @@ function ClimaLand.make_update_aux( # field.:($index) .= value # works # @ field.:($$index) = value # works labels = hydraulics.h_stem > 0 ? [:stem, :leaf] : [:leaf] - @inbounds for i in eachindex(labels) - ip1 = i + 1 + if hydraulics.h_stem > 0 + i = 1 + ip1 = 2 @. ψ.:($$ip1) = PlantHydraulics.water_retention_curve( retention_model, PlantHydraulics.effective_saturation(ν, ϑ_l.:($$ip1)), From 328e3b32c200848ee39669165aa4d1f3e7e0f57d Mon Sep 17 00:00:00 2001 From: Edward Speer Date: Sun, 2 Jun 2024 01:47:17 -0700 Subject: [PATCH 09/10] reformat --- src/integrated/soil_canopy_model.jl | 4 +--- src/standalone/Vegetation/Canopy.jl | 2 +- src/standalone/Vegetation/canopy_boundary_fluxes.jl | 2 +- 3 files changed, 3 insertions(+), 5 deletions(-) diff --git a/src/integrated/soil_canopy_model.jl b/src/integrated/soil_canopy_model.jl index 20b1cd506d..92b13171a3 100644 --- a/src/integrated/soil_canopy_model.jl +++ b/src/integrated/soil_canopy_model.jl @@ -257,9 +257,7 @@ function make_update_boundary_fluxes( area_index = p.canopy.hydraulics.area_index - labels = - land.canopy.hydraulics.h_stem > 0 ? [:stem, :leaf] : - [:leaf] + labels = land.canopy.hydraulics.h_stem > 0 ? [:stem, :leaf] : [:leaf] above_ground_area_index = getproperty(area_index, labels[1]) diff --git a/src/standalone/Vegetation/Canopy.jl b/src/standalone/Vegetation/Canopy.jl index 11ce5c8950..f52c9cbd96 100644 --- a/src/standalone/Vegetation/Canopy.jl +++ b/src/standalone/Vegetation/Canopy.jl @@ -548,7 +548,7 @@ function ClimaLand.make_update_aux( T_canopy = canopy_temperature(canopy.energy, canopy, Y, p, t) # update moisture stress - i_end = hydraulics.h_stem > 0 + 1 + i_end = (hydraulics.h_stem > 0 ? 1 : 0) + 1 @. β = moisture_stress(ψ.:($$i_end) * ρ_l * grav, sc, pc) # Update Rd, An, Vcmax25 (if applicable to model) in place diff --git a/src/standalone/Vegetation/canopy_boundary_fluxes.jl b/src/standalone/Vegetation/canopy_boundary_fluxes.jl index 673273e0fd..24c5163e3d 100644 --- a/src/standalone/Vegetation/canopy_boundary_fluxes.jl +++ b/src/standalone/Vegetation/canopy_boundary_fluxes.jl @@ -157,7 +157,7 @@ function canopy_boundary_fluxes!( shf = p.canopy.energy.shf lhf = p.canopy.energy.lhf r_ae = p.canopy.energy.r_ae - i_end = canopy.hydraulics.h_stem > 0 + 1 + i_end = (canopy.hydraulics.h_stem > 0 ? 1 : 0) + 1 # Compute transpiration, SHF, LHF canopy_turbulent_fluxes = turbulent_fluxes(atmos, canopy, Y, p, t) From e4cfba419282c3df3fc22aa53fc16a55a05baaf5 Mon Sep 17 00:00:00 2001 From: Edward Speer Date: Sun, 2 Jun 2024 10:38:54 -0700 Subject: [PATCH 10/10] remove compartment_surfaces --- src/standalone/Vegetation/Canopy.jl | 4 ++-- src/standalone/Vegetation/PlantHydraulics.jl | 2 +- src/standalone/Vegetation/canopy_boundary_fluxes.jl | 6 +++--- 3 files changed, 6 insertions(+), 6 deletions(-) diff --git a/src/standalone/Vegetation/Canopy.jl b/src/standalone/Vegetation/Canopy.jl index f52c9cbd96..9b22bc748e 100644 --- a/src/standalone/Vegetation/Canopy.jl +++ b/src/standalone/Vegetation/Canopy.jl @@ -519,7 +519,7 @@ function ClimaLand.make_update_aux( areaip1 = getproperty(area_index, labels[ip1]) midpoint = hydraulics.h_stem / 2 - midpointip1 = hydraulics.h_stem + (model.h_leaf / 2) + midpointip1 = hydraulics.h_stem + (hydraulics.h_leaf / 2) # Compute the flux*area between the current compartment `i` # and the compartment above. @@ -568,7 +568,7 @@ function ClimaLand.make_update_aux( @. GPP = compute_GPP(An, K, LAI, Ω) @. gs = medlyn_conductance(g0, Drel, medlyn_factor, An, c_co2_air) # update autotrophic respiration - h_canopy = hydraulics.compartment_surfaces[end] + h_canopy = hydraulics.h_stem + hydraulics.h_leaf @. Ra = compute_autrophic_respiration( canopy.autotrophic_respiration, Vcmax25, diff --git a/src/standalone/Vegetation/PlantHydraulics.jl b/src/standalone/Vegetation/PlantHydraulics.jl index 1c06a82a4d..c3b001935c 100644 --- a/src/standalone/Vegetation/PlantHydraulics.jl +++ b/src/standalone/Vegetation/PlantHydraulics.jl @@ -555,7 +555,7 @@ function root_water_flux_per_ground_area!( t, ) where {FT} (; conductivity_model, root_distribution) = model.parameters - area_index = p.canopy.hydraulics, area_index + area_index = p.canopy.hydraulics.area_index # We can index into a field of Tuple{FT} to extract a field of FT # using the following notation: field.:index ψ_base = p.canopy.hydraulics.ψ.:1 diff --git a/src/standalone/Vegetation/canopy_boundary_fluxes.jl b/src/standalone/Vegetation/canopy_boundary_fluxes.jl index 24c5163e3d..690e22aa9b 100644 --- a/src/standalone/Vegetation/canopy_boundary_fluxes.jl +++ b/src/standalone/Vegetation/canopy_boundary_fluxes.jl @@ -16,7 +16,7 @@ model. See Cowan 1968; Brutsaert 1982, pp. 113–116; Campbell and Norman 1998, p. 71; Shuttleworth 2012, p. 343; Monteith and Unsworth 2013, p. 304. """ function ClimaLand.displacement_height(model::CanopyModel{FT}, Y, p) where {FT} - return FT(0.67) * model.hydraulics.compartment_surfaces[end] + return FT(0.67) * (model.hydraulics.h_stem + model.hydraulics.h_leaf) end """ @@ -68,8 +68,8 @@ end A helper function which returns the surface height for the canopy model, which is stored in the parameter struct. """ -function ClimaLand.surface_height(model::CanopyModel, _...) - return model.hydraulics.compartment_surfaces[1] +function ClimaLand.surface_height(model::CanopyModel{FT}, _...) where {FT} + return FT(0.0) end """