Re-enable ITime in ocean and sea-ice models

config/ci_configs/cmip_oceananigans_climaseaice.yml

@@ -37,7 +37,7 @@ surface_setup: "PrescribedSurface"
 t_end: "3hours"
 topo_smoothing: true
 topography: "Earth"
-use_itime: false
+use_itime: true
 viscous_sponge: true
 z_elem: 39
 z_max: 60000.0

config/ci_configs/cmip_oceananigans_climaseaice_bucket.yml

@@ -35,7 +35,7 @@ surface_setup: "PrescribedSurface"
 t_end: "3hours"
 topo_smoothing: true
 topography: "Earth"
-use_itime: false
+use_itime: true
 viscous_sponge: true
 z_elem: 39
 z_max: 60000.0

config/longrun_configs/cmip_edonly_bucket.yml

@@ -30,7 +30,7 @@ surface_setup: "PrescribedSurface"
 t_end: "366days"
 topo_smoothing: true
 topography: "Earth"
-use_itime: false
+use_itime: true
 extra_atmos_diagnostics:
   - short_name: [hur, hus, ta]
     period: 30days

config/longrun_configs/cmip_edonly_land.yml

@@ -37,4 +37,4 @@ surface_setup: "PrescribedSurface"
 t_end: "366days"
 topo_smoothing: true
 topography: "Earth"
-use_itime: false
+use_itime: true

config/longrun_configs/cmip_progedmf_land.yml

@@ -39,4 +39,4 @@ surface_setup: "PrescribedSurface"
 t_end: "186days"
 topo_smoothing: true
 topography: "Earth"
-use_itime: false
+use_itime: true

experiments/test/compare_cmip.jl

@@ -7,6 +7,7 @@ import ClimaCore as CC
 import Oceananigans as OC
 import ClimaSeaIce as CSI
 using ClimaSeaIce.SeaIceThermodynamics.HeatBoundaryConditions: IceWaterThermalEquilibrium
+using Oceananigans.Architectures: on_architecture, CPU
 import NCDatasets
 
 """
@@ -20,8 +21,10 @@ function _error(arr1::AbstractArray, arr2::AbstractArray; ABS_TOL = 100eps(eltyp
     # There are some parameters, e.g. Obukhov length, for which Inf
     # is a reasonable value (implying a stability parameter in the neutral boundary layer
     # regime, for instance). We account for such instances with the `isfinite` function.
-    arr1 = Array(arr1) .* isfinite.(Array(arr1))
-    arr2 = Array(arr2) .* isfinite.(Array(arr2))
+    arr1 = Array(arr1)
+    arr2 = Array(arr2)
+    arr1 = arr1 .* isfinite.(arr1)
+    arr2 = arr2 .* isfinite.(arr2)
     diff = abs.(arr1 .- arr2)
     denominator = abs.(arr1)
     error = ifelse.(denominator .> ABS_TOL, diff ./ denominator, diff)
@@ -71,7 +74,7 @@ function compare(
     },
 }
     pass = true
-    return _compare(pass, v1, v2; name, ignore)
+    _compare(pass, v1, v2; name, ignore)
 end
 
 function _compare(pass, v1::T, v2::T; name, ignore) where {T}
@@ -152,7 +155,10 @@ function _compare(
     name,
     ignore,
 ) where {FT <: AbstractFloat}
-    return pass && _compare(collect(v1), collect(v2); name, ignore)
+    # move to CPU and flatten to 1D so that we don't run into issues with halo offsets (i.e. `OffsetArrays`)
+    h1 = vec(on_architecture(CPU(), v1))
+    h2 = vec(on_architecture(CPU(), v2))
+    return pass && _compare(h1, h2; name, ignore)
 end
 
 function _compare(

experiments/test/restart_cmip.yml

@@ -25,7 +25,7 @@ rayleigh_sponge: true
 simple_ocean: true
 save_cache: true
 start_date: "20100101"
-use_itime: false
+use_itime: true
 surface_setup: "PrescribedSurface"
 t_end: "480secs"
 topo_smoothing: true

ext/ClimaCouplerCMIPExt/clima_seaice.jl

@@ -85,7 +85,15 @@ function ClimaSeaIceSimulation(
     arch = OC.Architectures.architecture(grid)
 
     advection = ocean.ocean.model.advection.T
-    ice = CO.SeaIces.sea_ice_simulation(grid, ocean.ocean; Δt = float(dt), advection)
+
+    ice = CO.SeaIces.sea_ice_simulation(
+        grid,
+        ocean.ocean;
+        clock = deepcopy(ocean.ocean.model.clock),
+        stop_time = ocean.ocean.stop_time,
+        Δt = float(dt),
+        advection,
+    )
 
     ocean_ice_flux_formulation = CO.InterfaceComputations.ThreeEquationHeatFlux(ice)
     interface_temperature = OC.Field{OC.Center, OC.Center, Nothing}(grid)

ext/ClimaCouplerCMIPExt/oceananigans.jl

@@ -41,17 +41,22 @@ function Interfacer.OceanSimulation(::Type{FT}, ::Val{:oceananigans}; kwargs...)
 end
 
 """
-    tripolar_ocean_simulation(arch; active_cells_map, kwargs...)
+    tripolar_ocean_simulation(arch; clock, stop_time, active_cells_map, kwargs...)
 
 One-degree tripolar ocean setup matching
 `ClimaOcean.OceanConfigurations.one_degree_tripolar_ocean`, but exposing
 `active_cells_map`. Tripolar immersed grids with `active_cells_map = true`
 exceed the 4 KiB CUDA kernel parameter limit on sm_60 (P100) even with a
 minimal closure; disabling the map keeps the tripolar grid while staying
 within the limit.
+
+`clock` and `stop_time` are required and forwarded to `ClimaOcean.ocean_simulation`;
+the coupler builds them from `tspan` so the ocean steps in calendar time when coupled.
 """
 function tripolar_ocean_simulation(
     arch;
+    clock,
+    stop_time,
     active_cells_map = true,
     zstar = true,
     Nz = 32,
@@ -79,6 +84,8 @@ function tripolar_ocean_simulation(
 
     return CO.ocean_simulation(
         grid;
+        clock,
+        stop_time,
         momentum_advection,
         tracer_advection,
         free_surface,
@@ -209,24 +216,17 @@ function OceananigansSimulation(
 
     if tspan[1] isa ITime
         # create a model clock that uses DateTime, for compatibility with ITime.
-        model_clock = OC.TimeSteppers.Clock(time = start_date)
+        clock = OC.TimeSteppers.Clock(time = start_date)
         stop_time = Dates.DateTime(tspan[2])
     elseif tspan[1] isa Float64
-        model_clock = OC.TimeSteppers.Clock{Float64}(time = tspan[1])
+        clock = OC.TimeSteppers.Clock{Float64}(time = tspan[1])
         stop_time = tspan[2]
     else
         error("Unsupported time type: $(typeof(tspan[1]))")
     end
 
-    ocean = ocean_simulation(
-        arch,
-        ocean_grid;
-        simple_ocean,
-        closure,
-        clock = model_clock,
-        depth,
-    )
-    ocean.stop_time = stop_time
+    ocean =
+        ocean_simulation(arch, ocean_grid; clock, stop_time, simple_ocean, closure, depth)
     ocean.Δt = float(dt)
 
     # Set initial condition to EN4 state estimate at start_date
@@ -262,7 +262,7 @@ function OceananigansSimulation(
         ocean_properties,
         remapping,
         ice_concentration,
-        ocean.Δt,
+        dt,
     )
 
     add_ocean_diagnostics!(