Configuring models

Overview

Terrarium models are constructed by passing a grid and optional keyword arguments that customize component choices and parameters. This page shows how to build SoilModel, LandModel, and VegetationModel with various configurations, and how to inspect their state variables.

Setting up a SoilModel

All models in Terrarium are constructed first and foremost from a grid representing the spatial discretization. As with many of our other examples, we will use here a simple ColumnGrid representing a single vertical column:

num_columns = 1
grid = ColumnGrid(ExponentialSpacing(N = 10), num_columns)

ColumnGrid{Float64} on CPU() with
1×1×10 RectilinearGrid{Float64, Periodic, Flat, Bounded} on CPU with 1×0×3 halo
├── Periodic x ∈ [0.0, 1.0)      regularly spaced with Δx=1.0
├── Flat y                       
└── Bounded  z ∈ [-175.387, 0.0] variably spaced with min(Δz)=0.05, max(Δz)=100.0

We will start here by constructing a simple SoilModel using its defaults:

model = SoilModel(grid)

SoilModel{Float64} on CPU()
├── grid:  ColumnGrid{Float64, CPU} with dimensions (1, 1, 10)
├── soil:  SoilEnergyWaterCarbon{Float64, HomogeneousStratigraphy{Float64, ConstantSoilPorosity{Float64}}, SoilEnergyBalance{Float64, Terrarium.ExplicitTwoPhaseHeatConduction, SoilEnergyTemperatureClosure, SoilThermalProperties{Float64, FreeWater, InverseQuadratic}}, SoilHydrology{Float64, NoFlow, SoilSaturationPressureClosure, SoilHydraulicsSURFEX{Float64, BrooksCorey{FreezeCurves.SoilWaterVolume{Float64, Float64, Float64}, Float64, Float64}, UnsatKLinear{Float64}}, Nothing}, ConstantSoilCarbonDensity{Float64}}
├── constants:  PhysicalConstants{Float64}
├── initializer:  DefaultInitializer{Float64}

Calling variables(model) will give us a more detailed look at the state variables defined by this model:

variables(model)

Variables
├─ Prognostic: 
├── internal_energy [J m^-3] on XYZ{Center, Center, Center}
├─ Auxiliary: 
├── temperature [°C] on XYZ{Center, Center, Center}
├── liquid_water_fraction [-] on XYZ{Center, Center, Center}
├── ground_temperature [°C] on XY{Center, Center}
├── saturation_water_ice [-] on XYZ{Center, Center, Center}
├── water_table [m] on XY{Center, Center}
├── hydraulic_conductivity [m s^-1] on XYZ{Center, Center, Face}
├─ Inputs: 
├─ Namespaces:

By default SoilModel uses a NoFlow soil hydrology scheme that treats soil water/ice as constant in time. Let's try changing that. The (hopefully) easiest way to learn how to do this would be to go look at the documentation page for soil hydrology. But let's be lazy and try to figure it ourselves.

We can see above that model has a property soil (see also the page for SoilModel). Let's inspect that:

model.soil

SoilEnergyWaterCarbon{Float64}
├── Processes:
├──── energy:  SoilEnergyBalance{Float64, Terrarium.ExplicitTwoPhaseHeatConduction, SoilEnergyTemperatureClosure, SoilThermalProperties{Float64, FreeWater, InverseQuadratic}}
├──── hydrology:  SoilHydrology{Float64, NoFlow, SoilSaturationPressureClosure, SoilHydraulicsSURFEX{Float64, BrooksCorey{FreezeCurves.SoilWaterVolume{Float64, Float64, Float64}, Float64, Float64}, UnsatKLinear{Float64}}, Nothing}
├──── biogeochem:  ConstantSoilCarbonDensity{Float64}
├── Properties:
├──── strat:  HomogeneousStratigraphy{Float64, ConstantSoilPorosity{Float64}}

We can see that soil is a coupled process type SoilEnergyWaterCarbon with processes energy, hydrology, and biogeochem. Looking at model.soil.hydrology:

model.soil.hydrology

SoilHydrology{Float64}
├── vertical_flow:  NoFlow
├── closure:  SoilSaturationPressureClosure
├── hydraulic_properties:  SoilHydraulicsSURFEX{Float64, BrooksCorey{FreezeCurves.SoilWaterVolume{Float64, Float64, Float64}, Float64, Float64}, UnsatKLinear{Float64}}
├── vwc_forcing:  Nothing

we can see that it has a property vertical_flow. As a general rule, most process and parameterization types in Terrarium declare abstract base types which can help us figure out what our choices are. We can use some built-in Julia tools for this.

# First get the type of vertical_flow
vftype = typeof(model.soil.hydrology.vertical_flow)
# Then get the "super" type (i.e. the abstract type it inherits from)
suptype = supertype(vftype)
# Then list all of the subtypes of the abstract type
subtypes(suptype)

2-element Vector{Any}:
 NoFlow
 RichardsEq

Aha! We see here a second implementation RichardsEq. This happens to correspond to the configuration option for SoilHydrology that enables vertical water flow governed by the Richardson-Richards equation. We can enable this by changing vertical_flow when constructing the process and then building back up the SoilModel from there.

hydrology = SoilHydrology(eltype(grid), vertical_flow = RichardsEq())
soil = SoilEnergyWaterCarbon(eltype(grid); hydrology)
model = SoilModel(grid; soil)
variables(model)

Variables
├─ Prognostic: 
├── internal_energy [J m^-3] on XYZ{Center, Center, Center}
├── saturation_water_ice [-] on XYZ{Center, Center, Center}
├── surface_excess_water [m] on XY{Center, Center}
├─ Auxiliary: 
├── temperature [°C] on XYZ{Center, Center, Center}
├── liquid_water_fraction [-] on XYZ{Center, Center, Center}
├── pressure_head [m] on XYZ{Center, Center, Center}
├── ground_temperature [°C] on XY{Center, Center}
├── hydraulic_conductivity [m s^-1] on XYZ{Center, Center, Face}
├── water_table [m] on XY{Center, Center}
├─ Inputs: 
├─ Namespaces:

We could also do something similar for the soil thermal properties:

thermal_properties = SoilThermalProperties(eltype(grid))
energy = SoilEnergyBalance(eltype(grid); thermal_properties)
hydrology = SoilHydrology(eltype(grid), RichardsEq()) # also a valid constructor, same as above
biogeochem = ConstantSoilCarbonDensity(eltype(grid))
soil = SoilEnergyWaterCarbon(eltype(grid); energy, hydrology, biogeochem)
model = SoilModel(grid; soil)

SoilModel{Float64} on CPU()
├── grid:  ColumnGrid{Float64, CPU} with dimensions (1, 1, 10)
├── soil:  SoilEnergyWaterCarbon{Float64, HomogeneousStratigraphy{Float64, ConstantSoilPorosity{Float64}}, SoilEnergyBalance{Float64, Terrarium.ExplicitTwoPhaseHeatConduction, SoilEnergyTemperatureClosure, SoilThermalProperties{Float64, FreeWater, InverseQuadratic}}, SoilHydrology{Float64, RichardsEq, SoilSaturationPressureClosure, SoilHydraulicsSURFEX{Float64, BrooksCorey{FreezeCurves.SoilWaterVolume{Float64, Float64, Float64}, Float64, Float64}, UnsatKLinear{Float64}}, Nothing}, ConstantSoilCarbonDensity{Float64}}
├── constants:  PhysicalConstants{Float64}
├── initializer:  DefaultInitializer{Float64}

Setting up a LandModel

LandModel couples soil, surface energy and water fluxes, and optionally vegetation. A key option here is whether or not the land surface has vegetation.

No vegetation (bare soil)

For bare soil simulations, we can simply pass vegetation=nothing:

model = LandModel(grid; vegetation = nothing)
variables(model)

Variables
├─ Prognostic: 
├── internal_energy [J m^-3] on XYZ{Center, Center, Center}
├── skin_temperature [°C] on XY{Center, Center}
├─ Auxiliary: 
├── temperature [°C] on XYZ{Center, Center, Center}
├── liquid_water_fraction [-] on XYZ{Center, Center, Center}
├── ground_temperature [°C] on XY{Center, Center}
├── saturation_water_ice [-] on XYZ{Center, Center, Center}
├── water_table [m] on XY{Center, Center}
├── hydraulic_conductivity [m s^-1] on XYZ{Center, Center, Face}
├── ground_heat_flux [W m^-2] on XY{Center, Center}
├── surface_shortwave_up [W m^-2] on XY{Center, Center}
├── surface_longwave_up [W m^-2] on XY{Center, Center}
├── surface_net_radiation [W m^-2] on XY{Center, Center}
├── sensible_heat_flux [W m^-2] on XY{Center, Center}
├── latent_heat_flux [W m^-2] on XY{Center, Center}
├── rainfall_ground [m s^-1] on XY{Center, Center}
├── evaporation_ground [m s^-1] on XY{Center, Center}
├── surface_runoff [m s^-1] on XY{Center, Center}
├── infiltration [m s^-1] on XY{Center, Center}
├─ Inputs: 
├── air_temperature [°C] on XY{Center, Center}
├── air_pressure [Pa] on XY{Center, Center}
├── windspeed [m s^-1] on XY{Center, Center}
├── specific_humidity [-] on XY{Center, Center}
├── rainfall [m s^-1] on XY{Center, Center}
├── snowfall [m s^-1] on XY{Center, Center}
├── surface_shortwave_down [W m^-2] on XY{Center, Center}
├── surface_longwave_down [W m^-2] on XY{Center, Center}
├── daytime_length [hr] on XY{Center, Center}
├── CO2 [ppm] on XY{Center, Center}
├─ Namespaces:

The bare soil configuration uses simpler hydrology (no canopy interception or canopy evapotranspiration) and focuses on direct soil-atmosphere exchanges. You can still customize soil components the same way as with SoilModel:

hydrology = SoilHydrology(eltype(grid), vertical_flow = RichardsEq())
soil = SoilEnergyWaterCarbon(eltype(grid); hydrology)
model = LandModel(grid; vegetation=nothing, soil)

LandModel{Float64} on CPU()
├── grid:  ColumnGrid{Float64, CPU} with dimensions (1, 1, 10)
├── vegetation:  Nothing
├── soil:  SoilEnergyWaterCarbon{Float64, HomogeneousStratigraphy{Float64, ConstantSoilPorosity{Float64}}, SoilEnergyBalance{Float64, Terrarium.ExplicitTwoPhaseHeatConduction, SoilEnergyTemperatureClosure, SoilThermalProperties{Float64, FreeWater, InverseQuadratic}}, SoilHydrology{Float64, RichardsEq, SoilSaturationPressureClosure, SoilHydraulicsSURFEX{Float64, BrooksCorey{FreezeCurves.SoilWaterVolume{Float64, Float64, Float64}, Float64, Float64}, UnsatKLinear{Float64}}, Nothing}, ConstantSoilCarbonDensity{Float64}}
├── surface_energy_balance:  SurfaceEnergyBalance{Float64, ImplicitSkinTemperature{Float64}, DiagnosedTurbulentFluxes{Float64}, DiagnosedRadiativeFluxes{Float64}, ConstantAlbedo{Float64}}
├── surface_hydrology:  SurfaceHydrology{Float64, NoCanopyInterception{Float64}, BareGroundEvaporation{Float64, SoilMoistureResistanceFactor{Float64}}, DirectSurfaceRunoff{Float64}}
├── atmosphere:  PrescribedAtmosphere{Float64, (:CO2,), RainSnow, LongShortWaveRadiation, Terrarium.SpecificHumidity, Terrarium.ConstantAerodynamics{Float64}, Tuple{TracerGas{Float64, :CO2}}}
├── constants:  PhysicalConstants{Float64}
├── initializer:  DefaultInitializer{Float64}

Customized vegetation processes

You can customize individual vegetation processes by passing them to VegetationCarbon:

photosynthesis = LUEPhotosynthesis(eltype(grid))
stomatal_conductance = MedlynStomatalConductance(eltype(grid))
carbon_dynamics = PALADYNCarbonDynamics(eltype(grid))
vegetation = VegetationCarbon(eltype(grid);
    photosynthesis,
    stomatal_conductance,
    carbon_dynamics
)
model = LandModel(grid; vegetation)

LandModel{Float64} on CPU()
├── grid:  ColumnGrid{Float64, CPU} with dimensions (1, 1, 10)
├── vegetation:  VegetationCarbon{Float64, LUEPhotosynthesis{Float64}, MedlynStomatalConductance{Float64}, PALADYNAutotrophicRespiration{Float64}, PALADYNPhenology{Float64}, PALADYNCarbonDynamics{Float64}, PALADYNVegetationDynamics{Float64}, StaticExponentialRootDistribution{Float64}, FieldCapacityLimitedPAW{Float64}}
├── soil:  SoilEnergyWaterCarbon{Float64, HomogeneousStratigraphy{Float64, ConstantSoilPorosity{Float64}}, SoilEnergyBalance{Float64, Terrarium.ExplicitTwoPhaseHeatConduction, SoilEnergyTemperatureClosure, SoilThermalProperties{Float64, FreeWater, InverseQuadratic}}, SoilHydrology{Float64, RichardsEq, SoilSaturationPressureClosure, SoilHydraulicsSURFEX{Float64, BrooksCorey{FreezeCurves.SoilWaterVolume{Float64, Float64, Float64}, Float64, Float64}, UnsatKLinear{Float64}}, Nothing}, ConstantSoilCarbonDensity{Float64}}
├── surface_energy_balance:  SurfaceEnergyBalance{Float64, ImplicitSkinTemperature{Float64}, DiagnosedTurbulentFluxes{Float64}, DiagnosedRadiativeFluxes{Float64}, ConstantAlbedo{Float64}}
├── surface_hydrology:  SurfaceHydrology{Float64, PALADYNCanopyInterception{Float64}, PALADYNCanopyEvapotranspiration{Float64, ConstantEvaporationResistanceFactor{Float64}}, DirectSurfaceRunoff{Float64}}
├── atmosphere:  PrescribedAtmosphere{Float64, (:CO2,), RainSnow, LongShortWaveRadiation, Terrarium.SpecificHumidity, Terrarium.ConstantAerodynamics{Float64}, Tuple{TracerGas{Float64, :CO2}}}
├── constants:  PhysicalConstants{Float64}
├── initializer:  DefaultInitializer{Float64}

We could also customize surface energy balance processes, e.g. prescribed surface energy balance fluxes:

radiative_fluxes = PrescribedRadiativeFluxes(eltype(grid))
turbulent_fluxes = PrescribedTurbulentFluxes(eltype(grid))
seb = SurfaceEnergyBalance(eltype(grid);
    radiative_fluxes,
    turbulent_fluxes
)
vegetation = VegetationCarbon(eltype(grid))
model = LandModel(grid; vegetation, surface_energy_balance = seb)

LandModel{Float64} on CPU()
├── grid:  ColumnGrid{Float64, CPU} with dimensions (1, 1, 10)
├── vegetation:  VegetationCarbon{Float64, LUEPhotosynthesis{Float64}, MedlynStomatalConductance{Float64}, PALADYNAutotrophicRespiration{Float64}, PALADYNPhenology{Float64}, PALADYNCarbonDynamics{Float64}, PALADYNVegetationDynamics{Float64}, StaticExponentialRootDistribution{Float64}, FieldCapacityLimitedPAW{Float64}}
├── soil:  SoilEnergyWaterCarbon{Float64, HomogeneousStratigraphy{Float64, ConstantSoilPorosity{Float64}}, SoilEnergyBalance{Float64, Terrarium.ExplicitTwoPhaseHeatConduction, SoilEnergyTemperatureClosure, SoilThermalProperties{Float64, FreeWater, InverseQuadratic}}, SoilHydrology{Float64, RichardsEq, SoilSaturationPressureClosure, SoilHydraulicsSURFEX{Float64, BrooksCorey{FreezeCurves.SoilWaterVolume{Float64, Float64, Float64}, Float64, Float64}, UnsatKLinear{Float64}}, Nothing}, ConstantSoilCarbonDensity{Float64}}
├── surface_energy_balance:  SurfaceEnergyBalance{Float64, ImplicitSkinTemperature{Float64}, PrescribedTurbulentFluxes{Float64}, PrescribedRadiativeFluxes{Float64}, ConstantAlbedo{Float64}}
├── surface_hydrology:  SurfaceHydrology{Float64, PALADYNCanopyInterception{Float64}, PALADYNCanopyEvapotranspiration{Float64, ConstantEvaporationResistanceFactor{Float64}}, DirectSurfaceRunoff{Float64}}
├── atmosphere:  PrescribedAtmosphere{Float64, (:CO2,), RainSnow, LongShortWaveRadiation, Terrarium.SpecificHumidity, Terrarium.ConstantAerodynamics{Float64}, Tuple{TracerGas{Float64, :CO2}}}
├── constants:  PhysicalConstants{Float64}
├── initializer:  DefaultInitializer{Float64}

Changing parameters

Model parameters are stored as fields in their respective process implementations or the parameterization types therein. One simple way to change parameter values is to simply modify the values when constructing the relevant process types.

As an example, consider again the above example with SoilModel. Suppose we want to change the thermal conductivity of mineral and organic soil material. We could do as follows:

conductivities = SoilThermalConductivities(eltype(grid); mineral = 3.0, organic = 0.8)
thermal_properties = SoilThermalProperties(eltype(grid); conductivities)
energy = SoilEnergyBalance(eltype(grid); thermal_properties)
soil = SoilEnergyWaterCarbon(eltype(grid); energy)
model = SoilModel(grid; soil)

SoilModel{Float64} on CPU()
├── grid:  ColumnGrid{Float64, CPU} with dimensions (1, 1, 10)
├── soil:  SoilEnergyWaterCarbon{Float64, HomogeneousStratigraphy{Float64, ConstantSoilPorosity{Float64}}, SoilEnergyBalance{Float64, Terrarium.ExplicitTwoPhaseHeatConduction, SoilEnergyTemperatureClosure, SoilThermalProperties{Float64, FreeWater, InverseQuadratic}}, SoilHydrology{Float64, NoFlow, SoilSaturationPressureClosure, SoilHydraulicsSURFEX{Float64, BrooksCorey{FreezeCurves.SoilWaterVolume{Float64, Float64, Float64}, Float64, Float64}, UnsatKLinear{Float64}}, Nothing}, ConstantSoilCarbonDensity{Float64}}
├── constants:  PhysicalConstants{Float64}
├── initializer:  DefaultInitializer{Float64}