remove two cache variables (#1177)

Author: kmdeck

experiments/integrated/fluxnet/snow_soil/simulation.jl

@@ -253,11 +253,41 @@ _ρ_i = FT(LP.ρ_cloud_ice(earth_param_set))
 _ρ_l = FT(LP.ρ_cloud_liq(earth_param_set))
 fig = Figure(size = (1600, 1200), fontsize = 26)
 ax1 = Axis(fig[2, 1], ylabel = "ΔEnergy (J/A)", xlabel = "Days")
+function compute_surface_energy_fluxes(p)
+    ρe_flux_falling_snow =
+        Snow.volumetric_energy_flux_falling_snow(atmos, p, land.snow.parameters)
+    ρe_flux_falling_rain =
+        Snow.volumetric_energy_flux_falling_rain(atmos, p, land.snow.parameters)
+
+    return @. (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 +
+                  ρe_flux_falling_rain
+              ) +
+              ρe_flux_falling_snow
+end
+
+function compute_surface_water_vol_fluxes(p)
+    return @. 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
+end
+
 ΔE_expected =
     cumsum(
         -1 .* [
             parent(
-                sv.saveval[k].atmos_energy_flux .-
+                compute_surface_energy_fluxes(sv.saveval[k]) .-
                 sv.saveval[k].soil.bottom_bc.heat,
             )[end] for k in 1:1:(length(sv.t) - 1)
         ],
@@ -270,7 +300,7 @@ E_measured = [
     cumsum(
         -1 .* [
             parent(
-                sv.saveval[k].atmos_water_flux .-
+                compute_surface_water_vol_fluxes(sv.saveval[k]) .-
                 sv.saveval[k].soil.bottom_bc.water .+ sv.saveval[k].soil.R_s,
             )[end] for k in 1:1:(length(sv.t) - 1)
         ],

src/integrated/land.jl

@@ -213,8 +213,6 @@ lsm_aux_vars(m::LandModel) = (
     :scratch3,
     :excess_water_flux,
     :excess_heat_flux,
-    :atmos_energy_flux,
-    :atmos_water_flux,
     :ground_heat_flux,
     :effective_soil_sfc_T,
     :sfc_scratch,
@@ -250,8 +248,6 @@ lsm_aux_types(m::LandModel{FT}) where {FT} = (
     FT,
     FT,
     FT,
-    FT,
-    FT,
     NamedTuple{(:PAR, :NIR), Tuple{FT, FT}},
 )
 
@@ -274,8 +270,6 @@ lsm_aux_domain_names(m::LandModel) = (
     :surface,
     :surface,
     :surface,
-    :surface,
-    :surface,
     :subsurface,
     :surface,
     :surface,
@@ -362,34 +356,6 @@ function make_update_boundary_fluxes(
         update_canopy_bf!(p, Y, t)
         # Update soil CO2
         update_soilco2_bf!(p, 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
-
     end
     return update_boundary_fluxes!
 end

src/integrated/soil_snow_model.jl

@@ -107,8 +107,6 @@ included in the integrated Soil-Snow model.
 lsm_aux_vars(m::SoilSnowModel) = (
     :excess_water_flux,
     :excess_heat_flux,
-    :atmos_energy_flux,
-    :atmos_water_flux,
     :ground_heat_flux,
     :effective_soil_sfc_T,
     :sfc_scratch,
@@ -121,26 +119,16 @@ lsm_aux_vars(m::SoilSnowModel) = (
 The types of the additional auxiliary variables that are
 included in the integrated Soil-Snow model.
 """
-lsm_aux_types(m::SoilSnowModel{FT}) where {FT} =
-    (FT, FT, FT, FT, FT, FT, FT, FT, FT)
+lsm_aux_types(m::SoilSnowModel{FT}) where {FT} = (FT, FT, FT, FT, FT, FT, FT)
 
 """
     lsm_aux_domain_names(m::SoilSnowModel)
 
 The domain names of the additional auxiliary variables that are
 included in the integrated Soil-Snow model.
 """
-lsm_aux_domain_names(m::SoilSnowModel) = (
-    :surface,
-    :surface,
-    :surface,
-    :surface,
-    :surface,
-    :surface,
-    :surface,
-    :subsurface,
-    :surface,
-)
+lsm_aux_domain_names(m::SoilSnowModel) =
+    (:surface, :surface, :surface, :surface, :surface, :subsurface, :surface)
 
 """
     make_update_boundary_fluxes(
@@ -191,31 +179,6 @@ function make_update_boundary_fluxes(
         # Now we can update the soil BC, and use the excess fluxes there in order
         # to conserve energy and water
         update_soil_bf!(p, 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.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
         return nothing
     end
     return update_boundary_fluxes!

test/integrated/soil_snow.jl

@@ -95,8 +95,6 @@ p_soil_alone = deepcopy(p)
 for lsm_aux_var in (
     :excess_water_flux,
     :excess_heat_flux,
-    :atmos_energy_flux,
-    :atmos_water_flux,
     :ground_heat_flux,
     :effective_soil_sfc_T,
     :sfc_scratch,
@@ -184,29 +182,6 @@ dY_soil_snow = deepcopy(Y) .* 0
 ClimaLand.source!(dY_soil_snow, src, Y, p, land_model.soil)
 @test all(parent(dY_soil_snow.soil.θ_i) .≈ 0)
 
-# Check conservation
-@test 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
-
-_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
-@test 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
 # Make sure soil boundary flux method also worked
 G = deepcopy(p.ground_heat_flux)
 p_snow_alone = deepcopy(p)