From 043b50a113932a03faad49443ba8c274f509b0c6 Mon Sep 17 00:00:00 2001 From: Stephan Kramer Date: Thu, 8 Aug 2024 18:49:29 +0100 Subject: [PATCH] Fix NonHydrostaticTimeIntegrator2D for projection of solutions with adaptivity The advance method of a timeintegrator should always advance from the solved-for-fields, which are the same fields that will contain the solution at the end of the timestep _after_ calling advance(). This is because this is the only field(s) that the "user" of a time integrator has control over, and the user may want to adjust/overwrite/project (for adaptivity) at the end of a timestep and expect the timestepper to advance from the adjusted values in the next timestep. Thus, if a time-integrator stores some copy of a "solution_old" it should be set at the _beginning_ of advance() to the current solution, rather than at the end. --- thetis/coupled_timeintegrator_2d.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/thetis/coupled_timeintegrator_2d.py b/thetis/coupled_timeintegrator_2d.py index 518bde275..a5f10e0a0 100644 --- a/thetis/coupled_timeintegrator_2d.py +++ b/thetis/coupled_timeintegrator_2d.py @@ -199,11 +199,12 @@ def initialize(self, solution2d): self.timesteppers.swe2d.initialize(self.fields.solution_2d) if self.nh_options.update_free_surface: self.timesteppers.fs2d.initialize(self.fields.elev_2d) - self.elev_old.assign(self.fields.elev_2d) @PETSc.Log.EventDecorator("thetis.NonHydrostaticTimeIntegrator2D.advance") def advance(self, t, update_forcings=None): """Advances equations for one time step.""" + if self.nh_options.update_free_surface: + self.elev_old.assign(self.fields.elev_2d) if self.serial_advancing: # --- advance in serial --- self.timesteppers.swe2d.advance(t, update_forcings=update_forcings) @@ -213,7 +214,6 @@ def advance(self, t, update_forcings=None): if self.nh_options.update_free_surface: self.fields.elev_2d.assign(self.elev_old) self.timesteppers.fs2d.advance(t, update_forcings=update_forcings) - self.elev_old.assign(self.fields.elev_2d) # update old solution if self.options.swe_timestepper_type == 'SSPIMEX': self.timesteppers.swe2d.erk.solution_old.assign(self.fields.solution_2d) @@ -234,4 +234,3 @@ def advance(self, t, update_forcings=None): elif last_stage: self.fields.elev_2d.assign(self.elev_old) self.timesteppers.fs2d.advance(t, update_forcings=update_forcings) - self.elev_old.assign(self.fields.elev_2d)