Soil heat with SFCC
In this example, we use the preset SoilHeatTile
to construct a Tile
consisting of a soil column with heat conduction forced using air temperatures from Samoylov Island. We use the SFCC formulation of Painter and Karra (2014). For the purpose of demonstration, we use the apparent heat capacity form of the heat equation in this example (i.e. Heat.Diffusion1D(:T)
).
using CryoGrid
First we set up the model:
forcings = loadforcings(CryoGrid.Presets.Forcings.Samoylov_ERA5_fitted_daily_1979_2020);
grid = CryoGrid.Presets.DefaultGrid_5cm
soilprofile, tempprofile = CryoGrid.Presets.SamoylovDefault
initT = initializer(:T, tempprofile)
sfcc = PainterKarra(swrc=VanGenuchten(α=0.1, n=2.0))
upperbc = TemperatureBC(forcings.Tair, NFactor(nf=0.6))
lowerbc = GeothermalHeatFlux(0.053u"W/m^2")
tile = CryoGrid.Presets.SoilHeatTile(upperbc, lowerbc, soilprofile, initT; grid=grid, freezecurve=sfcc)
tspan = (DateTime(2010,10,30),DateTime(2011,10,30))
u0, du0 = initialcondition!(tile, tspan)
prob = CryoGridProblem(tile, u0, tspan, saveat=3*3600.0, savevars=(:T,));
... then solve it with the built-in forward Euler integrator.
sol = @time solve(prob);
out = CryoGridOutput(sol)
CryoGridOutput with 2920 time steps (2010-10-30T00:00:00 to 2011-10-30T00:00:00) and 2 variables:
H => DimArray of Quantity{Float64, 𝐌 𝐋^-1 𝐓^-2, Unitful.FreeUnits{(J, m^-3), 𝐌 𝐋^-1 𝐓^-2, nothing}} with dimensions (218, 2920)
T => DimArray of Quantity{Float64, 𝚯, Unitful.FreeUnits{(K,), 𝚯, Unitful.Affine{-5463//20}}} with dimensions (218, 2920)
Finally, plot the resulting temperatures.
import Plots
zs = [5,10,15,20,25,30,40,50,100,500]u"cm"
Diagnostics.plot_at_depths(:T, out, zs, ylabel="Temperature (°C)")
This page was generated using Literate.jl.