Merge pull request #690 from CliMA/kd/spatially_varying_biogeochemistry

spatially varying and consistent biogeochemistry
CliMA · Jul 11, 2024 · 5928728 · 5928728
2 parents 18eb764 + 42c41a9
commit 5928728
Show file tree

Hide file tree

Showing 26 changed files with 290 additions and 290 deletions.
diff --git a/NEWS.md b/NEWS.md
@@ -3,6 +3,7 @@ ClimaLand.jl Release Notes
 
 main
 --------
+- PR[#690] Use the soil parameters in creating the biogeochemistry SoilMet driver for consistency.
 
 v0.13.0
 --------

diff --git a/docs/tutorials/integrated/soil_canopy_tutorial.jl b/docs/tutorials/integrated/soil_canopy_tutorial.jl
@@ -174,7 +174,7 @@ soil_model_type = Soil.EnergyHydrology{FT}
 # The domain is defined similarly to the soil domain described above.
 soilco2_type = Soil.Biogeochemistry.SoilCO2Model{FT}
 
-soilco2_ps = SoilCO2ModelParameters(FT; ν = soil_ν);
+soilco2_ps = SoilCO2ModelParameters(FT);
 
 # soil microbes args
 Csom = (z, t) -> eltype(z)(5); # kg C m⁻³, this is a guess, not measured at the site
@@ -185,18 +185,11 @@ soilco2_sources = (MicrobeProduction{FT}(),);
 
 soilco2_boundary_conditions = (; top = soilco2_top_bc, bottom = soilco2_bot_bc);
 
-soilco2_drivers = Soil.Biogeochemistry.SoilDrivers(
-    Soil.Biogeochemistry.PrognosticMet{FT}(),
-    Soil.Biogeochemistry.PrescribedSOC{FT}(Csom),
-    atmos,
-);
-
 soilco2_args = (;
     boundary_conditions = soilco2_boundary_conditions,
     sources = soilco2_sources,
     domain = soil_domain,
     parameters = soilco2_ps,
-    drivers = soilco2_drivers,
 );
 
 # Next we need to set up the [`CanopyModel`](https://clima.github.io/ClimaLand.jl/dev/APIs/canopy/Canopy/#Canopy-Model-Structs).
@@ -306,7 +299,7 @@ canopy_model_args = (; parameters = shared_params, domain = canopy_domain);
 # atmospheric and radiative flux conditions from the observations at the Ozark
 # site as was done in the previous tutorial.
 
-land_input = (atmos = atmos, radiation = radiation)
+land_input = (atmos = atmos, radiation = radiation, soil_organic_carbon = Csom)
 
 land = SoilCanopyModel{FT}(;
     soilco2_type = soilco2_type,

diff --git a/experiments/benchmarks/land.jl b/experiments/benchmarks/land.jl
@@ -400,10 +400,7 @@ function setup_prob(t0, tf, Δt; nelements = (101, 15))
     z0_b = FT(0.1) * z0_m
 
 
-    soilco2_ps = Soil.Biogeochemistry.SoilCO2ModelParameters(
-        FT;
-        ν = 1.0,# INCORRECT!
-    )
+    soilco2_ps = Soil.Biogeochemistry.SoilCO2ModelParameters(FT)
 
     soil_args = (domain = domain, parameters = soil_params)
     soil_model_type = Soil.EnergyHydrology{FT}
@@ -422,18 +419,11 @@ function setup_prob(t0, tf, Δt; nelements = (101, 15))
     soilco2_boundary_conditions =
         (; top = soilco2_top_bc, bottom = soilco2_bot_bc)
 
-    soilco2_drivers = Soil.Biogeochemistry.SoilDrivers(
-        Soil.Biogeochemistry.PrognosticMet{FT}(),
-        Soil.Biogeochemistry.PrescribedSOC{FT}(Csom),
-        atmos,
-    )
-
     soilco2_args = (;
         boundary_conditions = soilco2_boundary_conditions,
         sources = soilco2_sources,
         domain = domain,
         parameters = soilco2_ps,
-        drivers = soilco2_drivers,
     )
 
     # Now we set up the canopy model, which we set up by component:
@@ -538,7 +528,12 @@ function setup_prob(t0, tf, Δt; nelements = (101, 15))
     )
 
     # Integrated plant hydraulics and soil model
-    land_input = (atmos = atmos, radiation = radiation, runoff = runoff_model)
+    land_input = (
+        atmos = atmos,
+        radiation = radiation,
+        runoff = runoff_model,
+        soil_organic_carbon = Csom,
+    )
     land = SoilCanopyModel{FT}(;
         soilco2_type = soilco2_type,
         soilco2_args = soilco2_args,

diff --git a/experiments/integrated/fluxnet/ozark_pft.jl b/experiments/integrated/fluxnet/ozark_pft.jl
@@ -135,10 +135,7 @@ soil_model_type = Soil.EnergyHydrology{FT}
 # Soil microbes model
 soilco2_type = Soil.Biogeochemistry.SoilCO2Model{FT}
 
-soilco2_ps = SoilCO2ModelParameters(
-    FT;
-    ν = soil_ν, # same as soil
-)
+soilco2_ps = SoilCO2ModelParameters(FT)
 
 # soil microbes args
 Csom = (z, t) -> eltype(z)(5.0)
@@ -150,18 +147,11 @@ soilco2_sources = (MicrobeProduction{FT}(),)
 
 soilco2_boundary_conditions = (; top = soilco2_top_bc, bottom = soilco2_bot_bc)
 
-soilco2_drivers = Soil.Biogeochemistry.SoilDrivers(
-    Soil.Biogeochemistry.PrognosticMet{FT}(),
-    Soil.Biogeochemistry.PrescribedSOC{FT}(Csom),
-    atmos,
-)
-
 soilco2_args = (;
     boundary_conditions = soilco2_boundary_conditions,
     sources = soilco2_sources,
     domain = soil_domain,
     parameters = soilco2_ps,
-    drivers = soilco2_drivers,
 )
 
 # Now we set up the canopy model, which we set up by component:
@@ -242,7 +232,7 @@ shared_params = SharedCanopyParameters{FT, typeof(earth_param_set)}(
 canopy_model_args = (; parameters = shared_params, domain = canopy_domain)
 
 # Integrated plant hydraulics and soil model
-land_input = (atmos = atmos, radiation = radiation)
+land_input = (atmos = atmos, radiation = radiation, soil_organic_carbon = Csom)
 land = SoilCanopyModel{FT}(;
     soilco2_type = soilco2_type,
     soilco2_args = soilco2_args,

diff --git a/experiments/integrated/fluxnet/run_fluxnet.jl b/experiments/integrated/fluxnet/run_fluxnet.jl
@@ -94,10 +94,7 @@ soil_model_type = Soil.EnergyHydrology{FT}
 # Soil microbes model
 soilco2_type = Soil.Biogeochemistry.SoilCO2Model{FT}
 
-soilco2_ps = SoilCO2ModelParameters(
-    FT;
-    ν = soil_ν, # same as soil
-)
+soilco2_ps = SoilCO2ModelParameters(FT)
 
 # soil microbes args
 Csom = (z, t) -> eltype(z)(5.0)
@@ -109,18 +106,11 @@ soilco2_sources = (MicrobeProduction{FT}(),)
 
 soilco2_boundary_conditions = (; top = soilco2_top_bc, bottom = soilco2_bot_bc)
 
-soilco2_drivers = Soil.Biogeochemistry.SoilDrivers(
-    Soil.Biogeochemistry.PrognosticMet{FT}(),
-    Soil.Biogeochemistry.PrescribedSOC{FT}(Csom),
-    atmos,
-)
-
 soilco2_args = (;
     boundary_conditions = soilco2_boundary_conditions,
     sources = soilco2_sources,
     domain = soil_domain,
     parameters = soilco2_ps,
-    drivers = soilco2_drivers,
 )
 
 # Now we set up the canopy model, which we set up by component:
@@ -199,7 +189,7 @@ shared_params = SharedCanopyParameters{FT, typeof(earth_param_set)}(
 canopy_model_args = (; parameters = shared_params, domain = canopy_domain)
 
 # Integrated plant hydraulics and soil model
-land_input = (atmos = atmos, radiation = radiation)
+land_input = (atmos = atmos, radiation = radiation, soil_organic_carbon = Csom)
 land = SoilCanopyModel{FT}(;
     soilco2_type = soilco2_type,
     soilco2_args = soilco2_args,

diff --git a/experiments/integrated/global/global_parameters.jl b/experiments/integrated/global/global_parameters.jl
@@ -188,8 +188,4 @@ z0_m = FT(0.13) * h_canopy
 z0_b = FT(0.1) * z0_m
 
 
-soilco2_ps = Soil.Biogeochemistry.SoilCO2ModelParameters(
-    FT;
-    ν = 1.0,# INCORRECT! This should be the same as the soil porosity, but
-    # currently, SoilCO2 does not support spatially varying parameters
-)
+soilco2_ps = Soil.Biogeochemistry.SoilCO2ModelParameters(FT)
diff --git a/experiments/integrated/global/global_soil_canopy.jl b/experiments/integrated/global/global_soil_canopy.jl
@@ -190,18 +190,11 @@ soilco2_sources = (Soil.Biogeochemistry.MicrobeProduction{FT}(),)
 
 soilco2_boundary_conditions = (; top = soilco2_top_bc, bottom = soilco2_bot_bc)
 
-soilco2_drivers = Soil.Biogeochemistry.SoilDrivers(
-    Soil.Biogeochemistry.PrognosticMet{FT}(),
-    Soil.Biogeochemistry.PrescribedSOC{FT}(Csom),
-    atmos,
-)
-
 soilco2_args = (;
     boundary_conditions = soilco2_boundary_conditions,
     sources = soilco2_sources,
     domain = domain,
     parameters = soilco2_ps,
-    drivers = soilco2_drivers,
 )
 
 # Now we set up the canopy model, which we set up by component:
@@ -295,7 +288,12 @@ canopy_model_args = (;
 )
 
 # Integrated plant hydraulics and soil model
-land_input = (atmos = atmos, radiation = radiation, runoff = runoff_model)
+land_input = (
+    atmos = atmos,
+    radiation = radiation,
+    runoff = runoff_model,
+    soil_organic_carbon = Csom,
+)
 land = SoilCanopyModel{FT}(;
     soilco2_type = soilco2_type,
     soilco2_args = soilco2_args,

diff --git a/experiments/integrated/performance/conservation/ozark_conservation_setup.jl b/experiments/integrated/performance/conservation/ozark_conservation_setup.jl
@@ -79,7 +79,7 @@ soil_model_type = Soil.EnergyHydrology{FT}
 # Soil microbes model
 soilco2_type = Soil.Biogeochemistry.SoilCO2Model{FT}
 
-soilco2_ps = SoilCO2ModelParameters(FT; ν = soil_ν)
+soilco2_ps = SoilCO2ModelParameters(FT)
 
 # soil microbes args
 Csom = (z, t) -> eltype(z)(5.0)
@@ -91,18 +91,11 @@ soilco2_sources = (MicrobeProduction{FT}(),)
 
 soilco2_boundary_conditions = (; top = soilco2_top_bc, bottom = soilco2_bot_bc)
 
-soilco2_drivers = Soil.Biogeochemistry.SoilDrivers(
-    Soil.Biogeochemistry.PrognosticMet{FT}(),
-    Soil.Biogeochemistry.PrescribedSOC{FT}(Csom),
-    atmos,
-)
-
 soilco2_args = (;
     boundary_conditions = soilco2_boundary_conditions,
     sources = soilco2_sources,
     domain = soil_domain,
     parameters = soilco2_ps,
-    drivers = soilco2_drivers,
 )
 
 # Now we set up the canopy model, which we set up by component:
@@ -180,7 +173,7 @@ shared_params = SharedCanopyParameters{FT, typeof(earth_param_set)}(
 canopy_model_args = (; parameters = shared_params, domain = canopy_domain)
 
 # Integrated plant hydraulics and soil model
-land_input = (atmos = atmos, radiation = radiation)
+land_input = (atmos = atmos, radiation = radiation, soil_organic_carbon = Csom)
 land = SoilCanopyModel{FT}(;
     soilco2_type = soilco2_type,
     soilco2_args = soilco2_args,

diff --git a/experiments/integrated/performance/profile_allocations.jl b/experiments/integrated/performance/profile_allocations.jl
@@ -226,7 +226,7 @@ soil_model_type = Soil.EnergyHydrology{FT}
 # Soil microbes model
 soilco2_type = Soil.Biogeochemistry.SoilCO2Model{FT}
 
-soilco2_ps = SoilCO2ModelParameters(FT; ν = soil_ν)
+soilco2_ps = SoilCO2ModelParameters(FT)
 
 # soil microbes args
 Csom = (z, t) -> eltype(z)(5.0)
@@ -238,18 +238,11 @@ soilco2_sources = (MicrobeProduction{FT}(),)
 
 soilco2_boundary_conditions = (; top = soilco2_top_bc, bottom = soilco2_bot_bc)
 
-soilco2_drivers = Soil.Biogeochemistry.SoilDrivers(
-    Soil.Biogeochemistry.PrognosticMet{FT}(),
-    Soil.Biogeochemistry.PrescribedSOC{FT}(Csom),
-    atmos,
-)
-
 soilco2_args = (;
     boundary_conditions = soilco2_boundary_conditions,
     sources = soilco2_sources,
     domain = soil_domain,
     parameters = soilco2_ps,
-    drivers = soilco2_drivers,
 )
 
 # Now we set up the canopy model, which we set up by component:
@@ -315,7 +308,7 @@ shared_params = SharedCanopyParameters{FT, typeof(earth_param_set)}(
 canopy_model_args = (; parameters = shared_params, domain = canopy_domain)
 
 # Integrated plant hydraulics and soil model
-land_input = (atmos = atmos, radiation = radiation)
+land_input = (atmos = atmos, radiation = radiation, soil_organic_carbon = Csom)
 land = SoilCanopyModel{FT}(;
     soilco2_type = soilco2_type,
     soilco2_args = soilco2_args,

diff --git a/experiments/standalone/Biogeochemistry/experiment.jl b/experiments/standalone/Biogeochemistry/experiment.jl
@@ -71,7 +71,7 @@ for (FT, tf) in ((Float32, 2 * dt), (Float64, tf))
     # Make biogeochemistry model args
     Csom = (z, t) -> eltype(z)(5.0)
 
-    co2_parameters = Soil.Biogeochemistry.SoilCO2ModelParameters(FT; ν = 0.556)
+    co2_parameters = Soil.Biogeochemistry.SoilCO2ModelParameters(FT)
     C = FT(100)
 
     co2_top_bc = Soil.Biogeochemistry.SoilCO2StateBC((p, t) -> 0.0)
@@ -102,23 +102,17 @@ for (FT, tf) in ((Float32, 2 * dt), (Float64, tf))
         earth_param_set;
         c_co2 = TimeVaryingInput(atmos_co2),
     )
-
-    soil_drivers = Soil.Biogeochemistry.SoilDrivers(
-        Soil.Biogeochemistry.PrognosticMet{FT}(),
-        Soil.Biogeochemistry.PrescribedSOC{FT}(Csom),
-        atmos,
-    )
-
     soilco2_args = (;
         boundary_conditions = co2_boundary_conditions,
         sources = co2_sources,
         domain = lsm_domain,
         parameters = co2_parameters,
-        drivers = soil_drivers,
     )
 
     # Create integrated model instance
+    land_args = (atmos = atmos, soil_organic_carbon = Csom)
     model = LandSoilBiogeochemistry{FT}(;
+        land_args = land_args,
         soil_args = soil_args,
         soilco2_args = soilco2_args,
     )

diff --git a/ext/CreateParametersExt.jl b/ext/CreateParametersExt.jl
@@ -484,10 +484,6 @@ SoilCO2ModelParameters(::Type{FT}; kwargs...) where {FT <: AbstractFloat} =
     SoilCO2ModelParameters(CP.create_toml_dict(FT); kwargs...)
 
 function SoilCO2ModelParameters(toml_dict::CP.AbstractTOMLDict; kwargs...)
-    # These parameters have defaults that should not go in ClimaParams
-    θ_a100 = 0.1816
-    b = 4.547
-
     name_map = (;
         :CO2_diffusion_coefficient => :D_ref,
         :soil_C_substrate_diffusivity => :D_liq,
@@ -504,8 +500,6 @@ function SoilCO2ModelParameters(toml_dict::CP.AbstractTOMLDict; kwargs...)
     earth_param_set = LP.LandParameters(toml_dict)
     return SoilCO2ModelParameters{FT, typeof(earth_param_set)}(;
         earth_param_set,
-        θ_a100,
-        b,
         parameters...,
         kwargs...,
     )

diff --git a/lib/ClimaLandSimulations/experiments/ozark.jl b/lib/ClimaLandSimulations/experiments/ozark.jl
@@ -4,7 +4,7 @@ using ClimaLandSimulations.Fluxnet
 # default parameters, except for custom parameter hetero_resp.b
 sv_test, sol_test, Y_test, p_test = run_fluxnet(
     "US-MOz";
-    params = ozark_default_params(; hetero_resp = hetero_resp_ozark(; b = 2)),
+    params = ozark_default_params(; hetero_resp = hetero_resp_ozark()),
 )
 
 # defaults, except start time

diff --git a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-Ha1/US-Ha1_parameters.jl b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-Ha1/US-Ha1_parameters.jl
@@ -31,9 +31,7 @@ function harvard_default_params(;
 end
 
 function hetero_resp_harvard(;
-    θ_a100 = FT(0.1816),
     D_ref = FT(1.39e-5),
-    b = FT(4.547),
     D_liq = FT(3.17),
     α_sx = FT(194e3),
     Ea_sx = FT(61e3),
@@ -44,9 +42,7 @@ function hetero_resp_harvard(;
     p_sx = FT(0.024),
 )
     return HeteroRespP(
-        θ_a100,
         D_ref,
-        b,
         D_liq,
         α_sx,
         Ea_sx,

diff --git a/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-MOz/US-MOz_parameters.jl b/lib/ClimaLandSimulations/src/Fluxnet/fluxnet_sites/US-MOz/US-MOz_parameters.jl
@@ -31,9 +31,7 @@ function ozark_default_params(;
 end
 
 function hetero_resp_ozark(;
-    θ_a100 = FT(0.1816),
     D_ref = FT(1.39e-5),
-    b = FT(4.547),
     D_liq = FT(3.17),
     α_sx = FT(194e3),
     Ea_sx = FT(61e3),
@@ -44,9 +42,7 @@ function hetero_resp_ozark(;
     p_sx = FT(0.024),
 )
     return HeteroRespP(
-        θ_a100,
         D_ref,
-        b,
         D_liq,
         α_sx,
         Ea_sx,