Soil Models

AbstractSoilModel{FT} <: ClimaLand.AbstractImExModel{FT}

The abstract type for all soil models.

Currently, we only have plans to support a RichardsModel, simulating the flow of liquid water through soil via the Richardson-Richards equation, and a fully integrated soil heat and water model, with phase change.

RichardsModel

A model for simulating the flow of water in a porous medium by solving the Richardson-Richards Equation.

A variety of boundary condition types are supported, including FluxBC, RichardsAtmosDrivenFluxBC, MoistureStateBC, and FreeDrainage (only for the bottom of the domain).

If you wish to simulate soil hydrology under the context of a prescribed precipitation volume flux (m/s) as a function of time, the RichardsAtmosDrivenFluxBC type should be chosen. Please see the documentation for more details.

• parameters: the parameter set

• domain: the soil domain, using ClimaCore.Domains

• boundary_conditions: the boundary conditions, of type AbstractSoilBoundaryConditions

• sources: A tuple of sources, each of type AbstractSoilSource

• lateral_flow: A boolean flag which, when false, turns off the horizontal flow of water

RichardsModel{FT}(;
    parameters::RichardsParameters,
    domain::D,
    boundary_conditions::NamedTuple,
    sources::Tuple,
    lateral_flow::Bool = true
) where {FT, D}

A constructor for a RichardsModel, which sets the default value of lateral_flow to be true.

RichardsModel{FT}(domain, forcing;
                     runoff =  ClimaLand.Soil.Runoff.TOPMODELRunoff{FT}(f_over = FT(3.28),
                                                                        R_sb = FT(1.484e-4 / 1000),
                                                                        f_max = topmodel_fmax(domain.space.surface,FT),
                                                                        ),
                     retention_parameters = soil_vangenuchten_parameters(domain.space.subsurface, FT),
                     S_s = ClimaCore.Fields.zeros(domain.space.subsurface) .+ 1e-3,
                     )

Creates a RichardsModel model with the given float type FT, domain and forcing. Here, forcing should be a NamedTuple containing a field atmos with the atmospheric forcing.

Default spatially varying parameters (for retention curve parameters and specific storativity) are provided but can be changed with keyword arguments.

The runoff parameterization is also provided and can be changed via keyword argument.

EnergyHydrology <: AbstractSoilModel

A model for simulating the flow of water and heat in a porous medium by solving the Richardson-Richards equation and the heat equation, including terms for phase change.

A variety of boundary condition types are supported, including FluxBC, MoistureStateBC/TemperatureStateBC, FreeDrainage (only for the bottom of the domain), and an AtmosDrivenFluxBC (under which radiative fluxes and turbulent surface fluxes are computed and used as boundary conditions). Please see the documentation for this boundary condition type for more details.

• parameters: The parameter sets

• domain: the soil domain, using ClimaCore.Domains

• boundary_conditions: the boundary conditions for RRE and heat, of type AbstractSoilBoundaryConditions

• sources: A tuple of sources, each of type AbstractSoilSource

• lateral_flow: A boolean flag which, when false, turns off the horizontal flow of water and heat

EnergyHydrology{FT}(;
    parameters::PS
    domain::D,
    boundary_conditions::NamedTuple,
    sources::Tuple,
    lateral_flow::Bool = false,
) where {FT, D, PS}

A constructor for a EnergyHydrology model, which sets the default value of the lateral_flow flag to false.

EnergyHydrology{FT}(domain, forcing, toml_dict::CP.ParamDict;
                     prognostic_land_components = (:soil),
                     albedo = CLMTwoBandSoilAlbedo{FT}(; clm_soil_albedo_parameters(domain.space.surface)...),
                     runoff::Runoff.AbstractRunoffModel = Runoff.TOPMODELRunoff(toml_dict,
                                                            f_max = topmodel_fmax(domain.space.surface, FT),
                                                        ),
                     retention_parameters = soil_vangenuchten_parameters(domain.space.subsurface, FT),
                     composition_parameters = soil_composition_parameters(domain.space.subsurface, FT),
                     S_s = ClimaCore.Fields.zeros(domain.space.subsurface) .+ 1e-3,
                     z_0m = toml_dict["soil_momentum_roughness_length"],
                     z_0b = toml_dict["soil_scalar_roughness_length"],
                     emissivity = toml_dict["emissivity_bare_soil"],
                     additional_sources = (),
                     ) where {FT <: AbstractFloat}

Creates a EnergyHydrology model with the given float type FT, domain, toml_dict, forcing, and prognostic land components.

The argument forcing should be a NamedTuple containing two fields: atmos and radiation.

When running the soil model in standalone mode, prognostic_land_components = (:soil,), while for running integrated land models, this should be a list of the component models. This value of this argument must be the same across all components in the land model.

Default spatially varying parameters (for retention curve parameters, composition, and specific storativity) are provided but can be changed with keyword arguments. Note that these parameters must all be of the same type: either FT or ClimaCore Fields. By default they are Fields read in from data, so in practice this means if some values are provided as Floats, all of these parameter defaults must be overwritten as Floats.

retention_parameters should be a NamedTuple with the following fields:

• ν: soil porosity
• hydrology_cm: the hydrology closure model being used
• K_sat: saturated hydraulic conductivity
• θ_r: residual water content

composition_parameters should be a NamedTuple with the following fields:

• ν_ss_om: soil organic matter volume fraction
• ν_ss_quartz: soil quartz volume fraction
• ν_ss_gravel: soil gravel volume fraction

The runoff and albedo parameterizations are also provided and can be changed via keyword argument; additional sources may be required in your model if the soil model will be composed with other component models.

Roughness lengths and soil emissivity are currently treated as constants; these can be passed in as Floats by kwarg; otherwise the default values are used.

RichardsParameters{F <: Union{<: AbstractFloat, ClimaCore.Fields.Field}, C <: AbstractSoilHydrologyClosure}

A struct for storing parameters of the RichardsModel.

• ν: The porosity of the soil (m^3/m^3)

• hydrology_cm: The hydrology closure model: vanGenuchten or BrooksCorey

• K_sat: The saturated hydraulic conductivity (m/s)

• S_s: The specific storativity (1/m)

• θ_r: The residual water fraction (m^3/m^3)

RichardsParameters(;
    hydrology_cm::C,
    ν::F,
    K_sat::F,
    S_s::F,
    θ_r::F,
) where {F <: Union{<:AbstractFloat, ClimaCore.Fields.Field}, C}

Construct a RichardsParameters.

EnergyHydrologyParameters{
        FT <: AbstractFloat,
        F <: Union{<:AbstractFloat, ClimaCore.Fields.Field},
        AP,
        C,
        PSE,
    }

A parameter structure for the integrated soil water and energy equation system.

• soil composition and retention model parameters defined in the interior
• an albedo parameterization
• Scalar parameters: emissivity, α, β, γ, γT_ref, Ω,

roughness lengths z0, dds ) (FT)

• κ_dry: The dry soil thermal conductivity, W/m/K

• κ_sat_frozen: The saturated thermal conductivity of frozen soil, W/m/K

• κ_sat_unfrozen: The saturated thermal conductivity of unfrozen soil, W/m/K

• ρc_ds: The volumetric heat capacity of dry soil, J/m^3/K (per volume dry soil, not per volume soil solids)

• ν: The porosity of the soil (m^3/m^3)

• ν_ss_om: The volumetric fraction of the soil solids in organic matter (m^3/m^3)

• ν_ss_quartz: The volumetric fraction of the soil solids in quartz (m^3/m^3)

• ν_ss_gravel: The volumetric fraction of the soil solids in gravel (m^3/m^3)

• α: The parameter α used in computing Kersten number, unitless

• β: The parameter β used in computing Kersten number, unitless

• hydrology_cm: The soil hydrology closure model: van Genuchten or Brooks and Corey

• K_sat: The saturated hydraulic conductivity (m/s)

• S_s: The specific storativity (1/m)

• θ_r: The residual water fraction (m^3/m^3)

• Ω: Ice impedance factor for the hydraulic conductivity

• γ: Coefficient of viscosity factor for the hydraulic conductivity

• γT_ref: Reference temperature for the viscosity factor

• albedo: Soil Albedo Parameterization

• emissivity: Soil Emissivity

• z_0m: Roughness length for momentum

• z_0b: Roughness length for scalars

• d_ds: Maximum dry soil layer thickness under evaporation (m)

• earth_param_set: Physical constants and clima-wide parameters

EnergyHydrologyParameters(
    toml_dict;
    ν,
    ν_ss_om,
    ν_ss_quartz,
    ν_ss_gravel,
    hydrology_cm,
    K_sat,
    S_s,
    θ_r,
    albedo = Soil.ConstantTwoBandSoilAlbedo{FT}(),
    kwargs...,)

EnergyHydrologyParameters has two constructors: float-type and toml dict based. Additional parameters must be added manually: ν, ν_ss_om, ν_ss_quartz, ν_ss_gravel, hydrology_cm,Ksat,Ss, andθ_r. All parameters can be manually overridden via keyword arguments. Note, however, that certain parameters must have the same type (e.g, if a field is supplied for porosity, it must be supplied for all other parameters defined in the interior of the domain). Some parameters are defined only on the surface of the domain (e.g albedo), while other are defined everywhere (e.g. porosity). These are indicated with typesFandSF`. If both dry/wet albedos and general albedos are given as keywords, the dry/wet albedos will override the general albedos.

Please see the EnergyHydrologyParameters documentation for a complete list.

volumetric_liquid_fraction(ϑ_l::FT, ν_eff::FT, θ_r::FT) where {FT}

A pointwise function returning the volumetric liquid fraction given the augmented liquid fraction and the effective porosity. The output is guaranteed to be in (θ_r, ν_eff].

For Richards model, ν_eff = ν; and the clipping below is not required, (and will do nothing; ν_eff_safe = ν_eff), but we leave it in for a simpler interface.

pressure_head(
    cm::vanGenuchten{FT},
    θ_r::FT,
    ϑ_l::FT,
    ν_eff::FT,
    S_s::FT,
) where {FT}

A point-wise function returning the pressure head in variably saturated soil, using the van Genuchten matric potential if the soil is not saturated, and an approximation of the positive pressure in the soil if the soil is saturated.

pressure_head(
    cm::BrooksCorey{FT},
    θ_r::FT,
    ϑ_l::FT,
    ν_eff::FT,
    S_s::FT,
) where {FT}

A point-wise function returning the pressure head in variably saturated soil, using the Brooks and Corey matric potential if the soil is not saturated, and an approximation of the positive pressure in the soil if the soil is saturated.

hydraulic_conductivity(cm::vanGenuchten{FT}, K_sat::FT, S::FT) where {FT}

A point-wise function returning the hydraulic conductivity, using the van Genuchten formulation.

 hydraulic_conductivity(cm::BrooksCorey{FT}, K_sat::FT, S::FT) where {FT}

A point-wise function returning the hydraulic conductivity, using the Brooks and Corey formulation.

impedance_factor(
    f_i::FT,
    Ω::FT
) where {FT}

Returns the multiplicative factor reducing conductivity when a fraction of ice f_i is present.

Only for use with the EnergyHydrology model.

viscosity_factor(
    T::FT,
    γ::FT,
    γT_ref::FT,
) where {FT}

Returns the multiplicative factor which accounts for the temperature dependence of the conductivity.

Only for use with the EnergyHydrology model.

effective_saturation(ν_eff::FT, ϑ_l::FT, θr::FT) where {FT}

A point-wise function computing the effective saturation given the effective porosity, augmented liquid fraction, and residual water fraction as input.

For Richards model, or any other parameterization where ice is not relevant, νeff = ν; and the clipping below is not required, (and will do nothing; νeffsafe = νeff), but we leave it in for a simpler interface.

 matric_potential(cm::vanGenuchten{FT}, S::FT) where {FT}

A point-wise function returning the matric potential, using the van Genuchten formulation.

 matric_potential(cm::BrooksCorey{FT}, S::FT) where {FT}

A point-wise function returning the matric potential, using the Brooks and Corey formulation.

dψdϑ(cm::vanGenuchten{FT}, ϑ, ν_eff, θ_r, S_s)

Computes and returns the derivative of the pressure head with respect to ϑ for the van Genuchten formulation.

dψdϑ(cm::BrooksCorey{FT}, ϑ, νeff, θr, S_s)

Computes and returns the derivative of the pressure head with respect to ϑ for the Brooks and Corey formulation.

inverse_matric_potential(cm::vanGenuchten{FT}, ψ::FT) where {FT}

A point-wise function returning the effective saturation, given the matric potential, using the van Genuchten formulation.

inverse_matric_potential(cm::BrooksCorey{FT}, ψ::FT) where {FT}

A point-wise function returning the effective saturation, given the matric potential, using the Brooks and Corey formulation.

AbstractSoilHydrologyClosure{FT <: AbstractFloat}

The abstract type of soil hydrology closure, of which vanGenuchten{FT} and BrooksCorey{FT} are the two supported concrete types.

To add a new parameterization, methods are required for:

• matric_potential,
• inversematricpotential,
• pressure_head,
• dψdϑ,
• hydraulic_conductivity.

vanGenuchten{FT} <: AbstractSoilHydrologyClosure{FT}

The van Genuchten soil hydrology closure, chosen when the hydraulic conductivity and matric potential are modeled using the van Genuchten parameterization (van Genuchten 1980; see also Table 8.2 of G. Bonan 2019).

• α: The inverse of the air entry potential (1/m)

• n: The van Genuchten pore-size distribution index (unitless)

• m: The van Genuchten parameter m = 1 - 1/n (unitless)

• S_c: A derived parameter: the critical saturation at which capillary flow no longer replenishes the surface

BrooksCorey{FT} <: AbstractSoilHydrologyClosure{FT}

The Brooks and Corey soil hydrology closure, chosen when the hydraulic conductivity and matric potential are modeled using the Brooks and Corey parameterization (Brooks and Corey, 1964, 1966; see also Table 8.2 of G. Bonan 2019).

• c: The pore-size distribution index (unitless)

• ψb: The air entry matric potential, when S=1 (m)

• S_c: A derived parameter: the critical saturation at which capillary flow no longer replenishes the surface

volumetric_heat_capacity(
    θ_l::FT,
    θ_i::FT,
    ρc_ds::FT,
    earth_param_set::EP,
) where {FT,EP}

Compute the expression for volumetric heat capacity.

κ_solid(ν_ss_om::FT,
        ν_ss_quartz::FT,
        κ_om::FT,
        κ_quartz::FT,
        κ_minerals::FT) where {FT}

Computes the thermal conductivity of the solid material in soil. The _ss_ subscript denotes that the volumetric fractions of the soil components are referred to the soil solid components, not including the pore space.

function κ_sat_frozen(
    κ_solid::FT,
    ν::FT,
    κ_ice::FT
) where {FT}

Computes the thermal conductivity for saturated frozen soil.

function κ_sat_unfrozen(
    κ_solid::FT,
    ν::FT,
    κ_l::FT
) where {FT}

Computes the thermal conductivity for saturated unfrozen soil.

κ_sat(
    θ_l::FT,
    θ_i::FT,
    κ_sat_unfrozen::FT,
    κ_sat_frozen::FT
) where {FT}

Compute the expression for saturated thermal conductivity of soil matrix.

function κ_dry(ρp::FT,
               ν::FT,
               κ_solid::FT,
               κ_air::FT;
               a = FT(0.053)) where {FT}

Computes the thermal conductivity of dry soil according to the model of Balland and Arp.

kersten_number(
    θ_i::FT,
    S_r::FT,
    α::FT,
    β::FT,
    ν_ss_om::FT,
    ν_ss_quartz::FT,
    ν_ss_gravel::FT,
    ) where {FT}

Compute the expression for the Kersten number, using the Balland and Arp model.

relative_saturation(
        θ_l::FT,
        θ_i::FT,
        ν::FT
) where {FT}

Compute the expression for relative saturation. This is referred to as θ_sat in Balland and Arp's paper.

volumetric_internal_energy(θ_i::FT, ρc_s::FT, T::FT,
                             earth_param_set::EP) where {FT, EP}

A pointwise function for computing the volumetric internal energy of the soil, given the volumetric ice content, volumetric heat capacity, and temperature.

volumetric_internal_energy_liq(T::FT, earth_param_set::EP) where {FT, EP}

A pointwise function for computing the volumetric internal energy of the liquid water in the soil, given the temperature T.

temperature_from_ρe_int(ρe_int::FT, θ_i::FT, ρc_s::FT
                        earth_param_set::EP) where {FT, EP}

A pointwise function for computing the temperature from the volumetric internal energy, volumetric ice content, and volumetric heat capacity of the soil.

thermal_conductivity(
    κ_dry::FT,
    K_e::FT,
    κ_sat::FT
) where {FT}

Compute the expression for thermal conductivity of soil matrix.

phase_change_source(
    θ_l::FT,
    θ_i::FT,
    T::FT,
    τ::FT,
    ν::FT,
    θ_r::FT,
    hydrology_cm::CM,
    _ρ_i::FT,
    _ρ_l::FT,
    _LH_f0::FT,
    _T_freeze::FT,
    _grav::FT
    ) where {FT, CM}

Returns the source term (1/s) used for converting liquid water and ice into each other during phase changes. Note that there are unitless prefactors multiplying this term in the equations.

Note that these equations match what is in Dall'Amico (for θstar, ψ(T), ψw0). We should double check them in the case where we have ϑl > θl, but they should be very close to the form we want regardless.

thermal_time(ρc::FT, Δz::FT, κ::FT) where {FT}

Returns the thermal timescale for temperature differences across a typical thickness Δz to equilibrate.