diff --git a/docs/Manifest-v1.11.toml b/docs/Manifest-v1.11.toml
index b5600e673c..3c3dcc84a5 100644
--- a/docs/Manifest-v1.11.toml
+++ b/docs/Manifest-v1.11.toml
@@ -2,7 +2,7 @@
 
 julia_version = "1.11.1"
 manifest_format = "2.0"
-project_hash = "90782a2a9c550ec4050c1b9541950a8a615c2564"
+project_hash = "375d00e829a5d8c49c86dce3a9f254044ec00898"
 
 [[deps.ADTypes]]
 git-tree-sha1 = "eea5d80188827b35333801ef97a40c2ed653b081"
@@ -25,16 +25,6 @@ git-tree-sha1 = "574baf8110975760d391c710b6341da1afa48d8c"
 uuid = "a4c015fc-c6ff-483c-b24f-f7ea428134e9"
 version = "0.0.1"
 
-[[deps.About]]
-deps = ["InteractiveUtils", "JuliaSyntaxHighlighting", "PrecompileTools", "StyledStrings"]
-git-tree-sha1 = "efbf5b623b7ee2a41ce5aed6299aa62ab7f2d5b9"
-uuid = "69d22d85-9f48-4c46-bbbe-7ad8341ff72a"
-version = "1.0.1"
-weakdeps = ["Pkg"]
-
-    [deps.About.extensions]
-    PkgExt = "Pkg"
-
 [[deps.AbstractFFTs]]
 deps = ["LinearAlgebra"]
 git-tree-sha1 = "d92ad398961a3ed262d8bf04a1a2b8340f915fef"
@@ -1464,17 +1454,6 @@ git-tree-sha1 = "af433a10f3942e882d3c671aacb203e006a5808f"
 uuid = "9c1d0b0a-7046-5b2e-a33f-ea22f176ac7e"
 version = "0.2.1+0"
 
-[[deps.JuliaSyntax]]
-git-tree-sha1 = "937da4713526b96ac9a178e2035019d3b78ead4a"
-uuid = "70703baa-626e-46a2-a12c-08ffd08c73b4"
-version = "0.4.10"
-
-[[deps.JuliaSyntaxHighlighting]]
-deps = ["JuliaSyntax", "StyledStrings"]
-git-tree-sha1 = "19ecee1ea81c60156486a92b062e443b6bba60b7"
-uuid = "ac6e5ff7-fb65-4e79-a425-ec3bc9c03011"
-version = "0.1.0"
-
 [[deps.JuliaVariables]]
 deps = ["MLStyle", "NameResolution"]
 git-tree-sha1 = "49fb3cb53362ddadb4415e9b73926d6b40709e70"
diff --git a/docs/Project.toml b/docs/Project.toml
index 1702430d56..83032395f1 100644
--- a/docs/Project.toml
+++ b/docs/Project.toml
@@ -1,5 +1,4 @@
 [deps]
-About = "69d22d85-9f48-4c46-bbbe-7ad8341ff72a"
 AbstractTrees = "1520ce14-60c1-5f80-bbc7-55ef81b5835c"
 CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
 CairoMakie = "13f3f980-e62b-5c42-98c6-ff1f3baf88f0"
diff --git a/docs/tutorials/standalone/Soil/phase_change_analytic.jl b/docs/tutorials/standalone/Soil/phase_change_analytic.jl
index 38e1471e9e..933f7b1ce0 100644
--- a/docs/tutorials/standalone/Soil/phase_change_analytic.jl
+++ b/docs/tutorials/standalone/Soil/phase_change_analytic.jl
@@ -9,15 +9,17 @@
 # volumetric ice fraction `θ_i` for all time, everywhere in the domain.
 # In this example, we add in a source term to the right hand side for both `θ_i`
 # and `ϑ_l` which models freezing and thawing and conserves water mass during the process.
-# The equations are
+
+# To facilitate comparison to an analytically tractable problem (the Stefan problem), we set
+# K_sat of the soil to zero. Then the simplified equations we solve are
 
 
 # ``
-# \frac{∂ ρe_{int}}{∂ t} =  ∇ ⋅ κ(θ_l, θ_i; ν, ...) ∇T + ∇ ⋅ ρe_{int_{liq}} K (T,θ_l, θ_i; ν, ...) \nabla h( ϑ_l, z; ν, ...)
+# \frac{∂ ρe_{int}}{∂ t} =  ∇ ⋅ κ(θ_l, θ_i; ν, ...) ∇T
 # ``
 
 # ``
-# \frac{ ∂ ϑ_l}{∂ t} = ∇ ⋅ K (T,θ_l, θ_i; ν, ...) ∇h( ϑ_l, z; ν, ...) -\frac{F_T}{ρ_l}
+# \frac{ ∂ ϑ_l}{∂ t} =  -\frac{F_T}{ρ_l}
 # ``
 
 # ``
@@ -36,12 +38,6 @@
 
 # ``κ`` is the thermal conductivity (W/m/K),
 
-# ``ρe_{int_{liq}}`` is the volumetric internal energy of liquid water (J/m^3),
-
-# ``K`` is the hydraulic conductivity (m/s),
-
-# ``h`` is the hydraulic head (m),
-
 # ``ϑ_l`` is the augmented volumetric liquid water fraction,
 
 # ``θ_i`` is the volumetric ice fraction,
@@ -113,10 +109,6 @@ boundary_fluxes = (;
     ),
 );
 
-# Create the source term instance. Our phase change model requires
-# knowledge of the vertical spacing, so we pass
-# that information in via an attribute of the
-# [`PhaseChange`](@ref ClimaLand.Soil.PhaseChange) structure.
 # Sources are added as elements of a list of sources. Here we just add freezing
 # and thawing.
 
@@ -141,7 +133,8 @@ soil = Soil.EnergyHydrology{FT}(;
 Y, p, coords = initialize(soil);
 
 # After which, we can specify the initial condition
-# function, and initialze the variables:
+# function, and initialze the variables. We chose these to match
+# the initial conditions of the Stefan problem:
 
 function init_soil!(Ysoil, z, params)
     ν = params.ν
@@ -207,7 +200,7 @@ sv = (;
     t = Array{Float64}(undef, length(saveat)),
     saveval = Array{NamedTuple}(undef, length(saveat)),
 )
-saving_cb = ClimaLand.NonInterpSavingCallback(sv, saveat)
+saving_cb = ClimaLand.NonInterpSavingCallback(sv, saveat);
 # Now we can solve the problem.
 sol = SciMLBase.solve(
     prob,