(0.95.3) Reset model clock and skip time_step!ing if next actuation…

… time is tiny (#3606)

* reset model clock and skip time_step if next_actuation_time is really tiny

* move next_actuation_time definition to after IterationInterval is defined

* bugfix

* add minimum_relative_step to Simulation

* apply logic

* bump minor version

* Update src/Simulations/simulation.jl

Co-authored-by: Gregory L. Wagner <[email protected]>

* added minimum_time_step to docstring

* better logic

* add test

* fix test

* update docstring call signature

* don't use next_actuation_time for setting next_time

* Update src/Simulations/run.jl

Co-authored-by: Simone Silvestri <[email protected]>

---------

Co-authored-by: Gregory L. Wagner <[email protected]>
Co-authored-by: Simone Silvestri <[email protected]>

src/Simulations/run.jl

@@ -102,11 +102,24 @@ end
 
 const ModelCallsite = Union{TendencyCallsite, UpdateStateCallsite}
 
+
+function time_step_or_skip!(sim)
+    model_callbacks = Tuple(cb for cb in values(sim.callbacks) if cb.callsite isa ModelCallsite)
+    Δt = aligned_time_step(sim, sim.Δt)
+    if Δt < sim.minimum_relative_step * sim.Δt
+        next_time = sim.model.clock.time + Δt
+        @warn "Resetting clock to $next_time and skipping time step of size Δt = $Δt"
+        sim.model.clock.time = next_time
+    else
+        time_step!(sim.model, Δt, callbacks=model_callbacks)
+    end
+end
+
+
 """ Step `sim`ulation forward by one time step. """
 function time_step!(sim::Simulation)
 
     start_time_step = time_ns()
-    model_callbacks = Tuple(cb for cb in values(sim.callbacks) if cb.callsite isa ModelCallsite)
 
     if !(sim.initialized) # execute initialization step
         initialize!(sim)
@@ -118,8 +131,7 @@ function time_step!(sim::Simulation)
                 start_time = time_ns()
             end
 
-            Δt = aligned_time_step(sim, sim.Δt)
-            time_step!(sim.model, Δt, callbacks=model_callbacks)
+            time_step_or_skip!(sim)
 
             if sim.verbose
                 elapsed_initial_step_time = prettytime(1e-9 * (time_ns() - start_time))
@@ -130,8 +142,7 @@ function time_step!(sim::Simulation)
         end
 
     else # business as usual...
-        Δt = aligned_time_step(sim, sim.Δt)
-        time_step!(sim.model, Δt, callbacks=model_callbacks)
+        time_step_or_skip!(sim)
     end
 
     # Callbacks and callback-like things

src/Simulations/simulation.jl

@@ -9,26 +9,28 @@ import Oceananigans.TimeSteppers: reset!
 default_progress(simulation) = nothing
 
 mutable struct Simulation{ML, DT, ST, DI, OW, CB}
-              model :: ML
-                 Δt :: DT
-     stop_iteration :: Float64
-          stop_time :: ST
-    wall_time_limit :: Float64
-        diagnostics :: DI
-     output_writers :: OW
-          callbacks :: CB
-      run_wall_time :: Float64
-            running :: Bool
-        initialized :: Bool
-            verbose :: Bool
+                    model :: ML
+                       Δt :: DT
+           stop_iteration :: Float64
+                stop_time :: ST
+          wall_time_limit :: Float64
+              diagnostics :: DI
+           output_writers :: OW
+                callbacks :: CB
+            run_wall_time :: Float64
+                  running :: Bool
+              initialized :: Bool
+                  verbose :: Bool
+    minimum_relative_step :: Float64
 end
 
 """
     Simulation(model; Δt,
                verbose = true,
                stop_iteration = Inf,
                stop_time = Inf,
-               wall_time_limit = Inf)
+               wall_time_limit = Inf,
+               minimum_relative_step = 0)
 
 Construct a `Simulation` for a `model` with time step `Δt`.
 
@@ -44,12 +46,16 @@ Keyword arguments
 
 - `wall_time_limit`: Stop the simulation if it's been running for longer than this many
                      seconds of wall clock time.
+- `minimum_relative_step`: time steps smaller than `Δt * minimum_relative_step` will be skipped.
+                           This avoids extremely high values when writing the pressure to disk.
+                           Default value is 0. See github.com/CliMA/Oceananigans.jl/issues/3593 for details.
 """
 function Simulation(model; Δt,
                     verbose = true,
                     stop_iteration = Inf,
                     stop_time = Inf,
-                    wall_time_limit = Inf)
+                    wall_time_limit = Inf,
+                    minimum_relative_step = 0)
 
    if stop_iteration == Inf && stop_time == Inf && wall_time_limit == Inf
        @warn "This simulation will run forever as stop iteration = stop time " *
@@ -88,7 +94,8 @@ function Simulation(model; Δt,
                      0.0,
                      false,
                      false,
-                     verbose)
+                     verbose,
+                     Float64(minimum_relative_step))
 end
 
 function Base.show(io::IO, s::Simulation)
@@ -100,6 +107,7 @@ function Base.show(io::IO, s::Simulation)
                      "├── Stop time: $(prettytime(s.stop_time))", "\n",
                      "├── Stop iteration: $(s.stop_iteration)", "\n",
                      "├── Wall time limit: $(s.wall_time_limit)", "\n",
+                     "├── Minimum relative step: ", prettysummary(s.minimum_relative_step), "\n",
                      "├── Callbacks: $(ordered_dict_show(s.callbacks, "│"))", "\n",
                      "├── Output writers: $(ordered_dict_show(s.output_writers, "│"))", "\n",
                      "└── Diagnostics: $(ordered_dict_show(s.diagnostics, "│"))")

src/Utils/schedules.jl

@@ -108,6 +108,8 @@ IterationInterval(interval; offset=0) = IterationInterval(interval, offset)
 
 (schedule::IterationInterval)(model) = (model.clock.iteration - schedule.offset) % schedule.interval == 0
 
+next_actuation_time(schedule::IterationInterval) = Inf
+
 #####
 ##### WallTimeInterval
 #####

test/test_simulations.jl

@@ -156,6 +156,17 @@ function run_basic_simulation_tests(arch)
     @show called_at
     @test all(called_at .≈ 0.0:schedule.interval:simulation.stop_time)
 
+
+    # Test that minimum_relative_step is running correctly
+    final_time = 1.0 + 1e-11
+    model.clock.time = 0
+    model.clock.iteration = 0
+    simulation = Simulation(model, Δt=1, stop_time=final_time, minimum_relative_step=1e-10)
+    run!(simulation)
+
+    @test simulation.model.clock.time == final_time
+    @test simulation.model.clock.iteration == 1
+
     return nothing
 end