explicit and implicit sources

Author: kmdeck

docs/tutorials/standalone/Usage/model_tutorial.jl

@@ -131,7 +131,7 @@ import ClimaLand:
     make_exp_tendency,
     make_compute_exp_tendency,
     make_update_aux,
-    make_update_boundary_fluxes,
+    make_update_explicit_boundary_fluxes,
     prognostic_vars,
     prognostic_types,
     prognostic_domain_names,
@@ -257,7 +257,9 @@ function ClimaLand.make_update_aux(model::RichardsTutorialModel)
     return update_aux!
 end;
 
-function ClimaLand.make_update_boundary_fluxes(model::RichardsTutorialModel)
+function ClimaLand.make_update_explicit_boundary_fluxes(
+    model::RichardsTutorialModel,
+)
     function update_boundary_fluxes!(p, Y, t)
         FT = ClimaLand.FTfromY(Y)
         p.soil.top_flux .= model.F_sfc
@@ -271,7 +273,7 @@ end;
 # ```julia
 # function make_exp_tendency(model::AbstractModel)
 #     update_aux! = make_update_aux(model)
-#     update_boundary_fluxes! = make_update_boundary_fluxes(model)
+#     update_boundary_fluxes! = make_update_explicit_boundary_fluxes(model)
 #     compute_exp_tendency! = make_compute_exp_tendency(model)
 #     function exp_tendency!(dY,Y,p,t)
 #         update_aux!(p,Y,t)

src/ClimaLand.jl

@@ -142,7 +142,7 @@ function make_imp_tendency(land::AbstractLandModel)
     compute_imp_tendency_list =
         map(x -> make_compute_imp_tendency(getproperty(land, x)), components)
     update_aux! = make_update_aux(land)
-    update_boundary_fluxes! = make_update_boundary_fluxes(land)
+    update_boundary_fluxes! = make_update_implicit_boundary_fluxes(land)
     function imp_tendency!(dY, Y, p, t)
         update_aux!(p, Y, t)
         update_boundary_fluxes!(p, Y, t)
@@ -158,7 +158,7 @@ function make_exp_tendency(land::AbstractLandModel)
     compute_exp_tendency_list =
         map(x -> make_compute_exp_tendency(getproperty(land, x)), components)
     update_aux! = make_update_aux(land)
-    update_boundary_fluxes! = make_update_boundary_fluxes(land)
+    update_boundary_fluxes! = make_update_explicit_boundary_fluxes(land)
     function exp_tendency!(dY, Y, p, t)
         update_aux!(p, Y, t)
         update_boundary_fluxes!(p, Y, t)
@@ -181,16 +181,32 @@ function make_update_aux(land::AbstractLandModel)
     return update_aux!
 end
 
-function make_update_boundary_fluxes(land::AbstractLandModel)
+function make_update_explicit_boundary_fluxes(land::AbstractLandModel)
     components = land_components(land)
-    update_fluxes_function_list =
-        map(x -> make_update_boundary_fluxes(getproperty(land, x)), components)
-    function update_boundary_fluxes!(p, Y, t)
+    update_fluxes_function_list = map(
+        x -> make_update_explicit_boundary_fluxes(getproperty(land, x)),
+        components,
+    )
+    function update_explicit_boundary_fluxes!(p, Y, t)
         for f! in update_fluxes_function_list
             f!(p, Y, t)
         end
     end
-    return update_boundary_fluxes!
+    return update_explicit_boundary_fluxes!
+end
+
+function make_update_implicit_boundary_fluxes(land::AbstractLandModel)
+    components = land_components(land)
+    update_fluxes_function_list = map(
+        x -> make_update_implicit_boundary_fluxes(getproperty(land, x)),
+        components,
+    )
+    function update_implicit_boundary_fluxes!(p, Y, t)
+        for f! in update_fluxes_function_list
+            f!(p, Y, t)
+        end
+    end
+    return update_implicit_boundary_fluxes!
 end
 
 function make_compute_jacobian(land::AbstractLandModel)

src/integrated/land.jl

@@ -277,6 +277,7 @@ lsm_aux_domain_names(m::LandModel) = (
     :surface,
 )
 
+
 """
     make_update_boundary_fluxes(
         land::LandModel{FT, MM, SM, RM, SnM},
@@ -296,7 +297,7 @@ models.
 This function is called each ode function evaluation, prior to the tendency function
 evaluation.
 """
-function make_update_boundary_fluxes(
+function make_update_explicit_boundary_fluxes(
     land::LandModel{FT, MM, SM, RM, SnM},
 ) where {
     FT,
@@ -305,16 +306,26 @@ function make_update_boundary_fluxes(
     RM <: Canopy.CanopyModel{FT},
     SnM <: Snow.SnowModel{FT},
 }
-    update_soil_bf! = make_update_boundary_fluxes(land.soil)
-    update_soilco2_bf! = make_update_boundary_fluxes(land.soilco2)
-    update_canopy_bf! = make_update_boundary_fluxes(land.canopy)
-    update_snow_bf! = make_update_boundary_fluxes(land.snow)
+    update_soil_bf! = make_update_explicit_boundary_fluxes(land.soil)
+    update_soilco2_bf! = make_update_explicit_boundary_fluxes(land.soilco2)
+    update_canopy_bf! = make_update_explicit_boundary_fluxes(land.canopy)
+    update_snow_bf! = make_update_explicit_boundary_fluxes(land.snow)
 
     function update_boundary_fluxes!(p, Y, t)
         earth_param_set = land.soil.parameters.earth_param_set
         # update root extraction
         update_root_extraction!(p, Y, t, land) # defined in src/integrated/soil_canopy_root_interactions.jl
 
+        # Compute the ground heat flux in place:
+        update_soil_snow_ground_heat_flux!(
+            p,
+            Y,
+            land.soil.parameters,
+            land.snow.parameters,
+            land.soil.domain,
+            FT,
+        )
+
         # Radiation - updates Rn for soil and snow also
         lsm_radiant_energy_fluxes!(
             p,
@@ -330,15 +341,6 @@ function make_update_boundary_fluxes(
             land.soil.parameters.earth_param_set,
         )
 
-        # Compute the ground heat flux in place:
-        update_soil_snow_ground_heat_flux!(
-            p,
-            Y,
-            land.soil.parameters,
-            land.snow.parameters,
-            land.soil.domain,
-            FT,
-        )
         #Now update snow boundary conditions, which rely on the ground heat flux
         update_snow_bf!(p, Y, t)
 
@@ -355,34 +357,74 @@ function make_update_boundary_fluxes(
         update_canopy_bf!(p, Y, t)
         # Update soil CO2
         update_soilco2_bf!(p, Y, t)
+    end
+    return update_boundary_fluxes!
+end
+
+
+"""
+    make_update_implicit_boundary_fluxes(
+        land::LandModel{FT, MM, SM, RM, SnM},
+    ) where {
+        FT,
+        MM <: Soil.Biogeochemistry.SoilCO2Model{FT},
+        SM <: Soil.RichardsModel{FT},
+        RM <: Canopy.CanopyModel{FT}
+        SnM <: Snow.SnowModel{FT}
+        }
+
+A method which makes a function; the returned function
+updates the additional auxiliary variables for the integrated model,
+as well as updates the boundary auxiliary variables for all component
+models.
+
+This function is called each ode function evaluation, prior to the tendency function
+evaluation.
+"""
+function make_update_implicit_boundary_fluxes(
+    land::LandModel{FT, MM, SM, RM, SnM},
+) where {
+    FT,
+    MM <: Soil.Biogeochemistry.SoilCO2Model{FT},
+    SM <: Soil.EnergyHydrology{FT},
+    RM <: Canopy.CanopyModel{FT},
+    SnM <: Snow.SnowModel{FT},
+}
+    update_soil_bf! = make_update_implicit_boundary_fluxes(land.soil)
+    update_soilco2_bf! = make_update_implicit_boundary_fluxes(land.soilco2)
+    update_canopy_bf! = make_update_implicit_boundary_fluxes(land.canopy)
+    update_snow_bf! = make_update_implicit_boundary_fluxes(land.snow)
+
+    function update_boundary_fluxes!(p, Y, t)
+        # Compute the ground heat flux in place:
+        update_soil_snow_ground_heat_flux!(
+            p,
+            Y,
+            land.soil.parameters,
+            land.snow.parameters,
+            land.soil.domain,
+            FT,
+        )
+
+        # Radiation - updates Rn for soil and snow also
+        lsm_radiant_energy_fluxes!(
+            p,
+            land,
+            land.canopy.radiative_transfer,
+            Y,
+            t,
+        )
 
-        # compute net flux with atmosphere, this is useful for monitoring conservation
-        _LH_f0 = FT(LP.LH_f0(earth_param_set))
-        _ρ_liq = FT(LP.ρ_cloud_liq(earth_param_set))
-        ρe_falling_snow = -_LH_f0 * _ρ_liq # per unit vol of liquid water
-        @. p.atmos_energy_flux =
-            (1 - p.snow.snow_cover_fraction) * (
-                p.soil.turbulent_fluxes.lhf + p.soil.turbulent_fluxes.shf -
-                p.soil.R_n
-            ) +
-            p.snow.snow_cover_fraction * (
-                p.snow.turbulent_fluxes.lhf + p.snow.turbulent_fluxes.shf -
-                p.snow.R_n
-            ) +
-            p.drivers.P_snow * ρe_falling_snow +
-            p.canopy.turbulent_fluxes.shf +
-            p.canopy.turbulent_fluxes.lhf - p.canopy.radiative_transfer.SW_n -
-            p.canopy.radiative_transfer.LW_n
-        @. p.atmos_water_flux =
-            p.drivers.P_snow +
-            p.drivers.P_liq +
-            (1 - p.snow.snow_cover_fraction) * (
-                p.soil.turbulent_fluxes.vapor_flux_liq +
-                p.soil.turbulent_fluxes.vapor_flux_ice
-            ) +
-            p.snow.snow_cover_fraction * p.snow.turbulent_fluxes.vapor_flux +
-            p.canopy.turbulent_fluxes.transpiration
+        # Effective (radiative) land properties
+        set_eff_land_radiation_properties!(
+            p,
+            land.soil.parameters.earth_param_set,
+        )
 
+        update_snow_bf!(p, Y, t)
+        update_soil_bf!(p, Y, t)
+        update_canopy_bf!(p, Y, t)
+        update_soilco2_bf!(p, Y, t)
     end
     return update_boundary_fluxes!
 end

src/integrated/pond_soil_model.jl

@@ -156,7 +156,7 @@ function infiltration_capacity(Y::ClimaCore.Fields.FieldVector, p::NamedTuple)
 end
 
 """
-    make_update_boundary_fluxes(
+    make_update_explicit_boundary_fluxes(
         land::LandHydrology{FT, SM, SW},
     ) where {FT, SM <: Soil.RichardsModel{FT}, SW <: Pond.PondModel{FT}}
 
@@ -167,11 +167,36 @@ term (runoff) for the surface water model.
 
 This function is called each ode function evaluation.
 """
-function make_update_boundary_fluxes(
+function make_update_explicit_boundary_fluxes(
     land::LandHydrology{FT, SM, SW},
 ) where {FT, SM <: Soil.RichardsModel{FT}, SW <: Pond.PondModel{FT}}
-    update_soil_bf! = make_update_boundary_fluxes(land.soil)
-    update_pond_bf! = make_update_boundary_fluxes(land.surface_water)
+    update_soil_bf! = make_update_explicit_boundary_fluxes(land.soil)
+    update_pond_bf! = make_update_explicit_boundary_fluxes(land.surface_water)
+    function update_boundary_fluxes!(p, Y, t)
+        update_soil_bf!(p, Y, t)
+        update_pond_bf!(p, Y, t)
+    end
+    return update_boundary_fluxes!
+end
+
+
+"""
+    make_update_implicit_boundary_fluxes(
+        land::LandHydrology{FT, SM, SW},
+    ) where {FT, SM <: Soil.RichardsModel{FT}, SW <: Pond.PondModel{FT}}
+
+A method which makes a function; the returned function
+updates the auxiliary variable `p.soil_infiltration`, which
+is needed for both the boundary condition for the soil model and the source
+term (runoff) for the surface water model.
+
+This function is called each ode function evaluation.
+"""
+function make_update_implicit_boundary_fluxes(
+    land::LandHydrology{FT, SM, SW},
+) where {FT, SM <: Soil.RichardsModel{FT}, SW <: Pond.PondModel{FT}}
+    update_soil_bf! = make_update_implicit_boundary_fluxes(land.soil)
+    update_pond_bf! = make_update_implicit_boundary_fluxes(land.surface_water)
     function update_boundary_fluxes!(p, Y, t)
         i_c = infiltration_capacity(Y, p)
         @. p.soil_infiltration =

src/integrated/soil_canopy_model.jl

@@ -209,7 +209,7 @@ lsm_aux_domain_names(m::SoilCanopyModel) = (
 )
 
 """
-    make_update_boundary_fluxes(
+    make_update_explicit_boundary_fluxes(
         land::SoilCanopyModel{FT, MM, SM, RM},
     ) where {
         FT,
@@ -226,20 +226,21 @@ models.
 This function is called each ode function evaluation, prior to the tendency function
 evaluation.
 """
-function make_update_boundary_fluxes(
+function make_update_explicit_boundary_fluxes(
     land::SoilCanopyModel{FT, MM, SM, RM},
 ) where {
     FT,
     MM <: Soil.Biogeochemistry.SoilCO2Model{FT},
     SM <: Soil.EnergyHydrology{FT},
     RM <: Canopy.CanopyModel{FT},
 }
-    update_soil_bf! = make_update_boundary_fluxes(land.soil)
-    update_soilco2_bf! = make_update_boundary_fluxes(land.soilco2)
-    update_canopy_bf! = make_update_boundary_fluxes(land.canopy)
+    update_soil_bf! = make_update_explicit_boundary_fluxes(land.soil)
+    update_soilco2_bf! = make_update_explicit_boundary_fluxes(land.soilco2)
+    update_canopy_bf! = make_update_explicit_boundary_fluxes(land.canopy)
     function update_boundary_fluxes!(p, Y, t)
         # update root extraction
         update_root_extraction!(p, Y, t, land)
+
         # Radiation
         lsm_radiant_energy_fluxes!(
             p,
@@ -262,6 +263,45 @@ function make_update_boundary_fluxes(
 end
 
 
+
+"""
+    make_update_implicit_boundary_fluxes(
+        land::SoilCanopyModel{FT, MM, SM, RM},
+    ) where {
+        FT,
+        MM <: Soil.Biogeochemistry.SoilCO2Model{FT},
+        SM <: Soil.RichardsModel{FT},
+        RM <: Canopy.CanopyModel{FT}
+        }
+
+A method which makes a function; the returned function
+updates the additional auxiliary variables for the integrated model,
+as well as updates the implicit_boundary auxiliary variables for all component
+models.
+
+This function is called each ode function evaluation, prior to the tendency function
+evaluation.
+"""
+function make_update_implicit_boundary_fluxes(
+    land::SoilCanopyModel{FT, MM, SM, RM},
+) where {
+    FT,
+    MM <: Soil.Biogeochemistry.SoilCO2Model{FT},
+    SM <: Soil.EnergyHydrology{FT},
+    RM <: Canopy.CanopyModel{FT},
+}
+    update_soil_bf! = make_update_implicit_boundary_fluxes(land.soil)
+    update_soilco2_bf! = make_update_implicit_boundary_fluxes(land.soilco2)
+    update_canopy_bf! = make_update_implicit_boundary_fluxes(land.canopy)
+    function update_boundary_fluxes!(p, Y, t)
+        update_soil_bf!(p, Y, t)
+        update_canopy_bf!(p, Y, t)
+        update_soilco2_bf!(p, Y, t)
+    end
+    return update_boundary_fluxes!
+end
+
+
 """
     lsm_radiant_energy_fluxes!(p, land::SoilCanopyModel{FT},
                                 canopy_radiation::Canopy.AbstractRadiationModel{FT},