zero out tendencies

Author: juliasloan25

experiments/ClimaEarth/components/ocean/prescr_seaice.jl

@@ -270,10 +270,13 @@ function ice_rhs!(dY, Y, p, t)
 
     # Calculate the conductive flux, and set it to zero if the area fraction is zero
     F_conductive = @. params.k_ice / (params.h) * (params.T_base - Y.T_sfc) # fluxes are defined to be positive when upward
-    @. F_conductive = ifelse(p.area_fraction ≈ 0, zero(F_conductive), F_conductive)
 
     rhs = @. (-p.F_turb_energy - p.F_radiative + F_conductive) /
        (params.h * params.ρ * params.c)
+
+    # Zero out tendencies where there is no ice, so that ice temperature remains constant there
+    @. rhs = ifelse(p.area_fraction ≈ 0, zero(rhs), rhs)
+
     # If tendencies lead to temperature above freezing, set temperature to freezing
     @. dY.T_sfc = min(rhs, (params.T_freeze - Y.T_sfc) / float(p.dt))
 end

experiments/ClimaEarth/components/ocean/slab_ocean.jl

@@ -199,6 +199,9 @@ function slab_ocean_rhs!(dY, Y, cache, t)
     p, F_turb_energy, F_radiative = cache
     rhs = @. -(F_turb_energy + F_radiative) / (p.h * p.ρ * p.c)
 
+    # Zero out tendencies where there is no ocean, so that temperature remains constant there
+    @. rhs = ifelse(cache.area_fraction ≈ 0, zero(rhs), rhs)
+
     # Note that the area fraction has already been applied to the fluxes,
     #  so we don't need to multiply by it here.
     @. dY.T_sfc = rhs * p.evolving_switch

experiments/ClimaEarth/test/component_model_tests/prescr_seaice_tests.jl

@@ -71,12 +71,12 @@ for FT in (Float32, Float64)
         # check that nothing changes with no fluxes (zero conduction when T_base == initial T_sfc)
         T_base_eq = IceSlabParameters{FT}().T_freeze - FT(5.0)
         dY, Y, p = test_sea_ice_rhs(T_base = T_base_eq)
-        @test all([i for i in extrema(dY)] .≈ [FT(0.0), FT(0.0)])
+        @test all(T -> T == FT(0), Array(parent(dY.T_sfc)))
 
         # check expected dT due to radiative flux only (again set T_base == initial T_sfc)
         dY, Y, p = test_sea_ice_rhs(F_radiative = 1.0, T_base = T_base_eq)
         dT_expected = -1.0 / (p.params.h * p.params.ρ * p.params.c)
-        @test all(extrema(dY) .≈ FT(dT_expected))
+        @test all(T -> T == FT(0) || T ≈ dT_expected, Array(parent(dY.T_sfc)))
 
         # check that tendency will not result in above freezing T
         dY, Y, p = test_sea_ice_rhs(F_radiative = 0.0, T_base = 330.0) # Float32 requires a large number here!