Source code for fluidsim_ocean.time_stepping

r"""Shallow water time stepper (:mod:`fluidsim_ocean.time_stepping`)
====================================================================

.. autoclass:: TimeSteppingPseudoSpectralSWnL
   :members:
   :private-members:

"""

from math import floor, log10, sqrt

from fluidsim.base.time_stepping.pseudo_spect import TimeSteppingPseudoSpectral

from fluiddyn.util import mpi


[docs]class TimeSteppingPseudoSpectralSWnL(TimeSteppingPseudoSpectral):
[docs] @staticmethod def _complete_params_with_default(params): """This static method is used to complete the *params* container. """ TimeSteppingPseudoSpectral._complete_params_with_default(params) params.time_stepping.deltat0 = 1e3 params.time_stepping.deltat_max = 1e3
def _init_compute_time_step(self): params_ts = self.params.time_stepping if params_ts.USE_CFL: if params_ts.cfl_coef is not None: self.CFL = params_ts.cfl_coef elif params_ts.type_time_scheme == "RK2": self.CFL = 0.4 elif params_ts.type_time_scheme == "RK4": self.CFL = 1.0 else: raise ValueError("Problem name time_scheme") else: self.deltat = params_ts.deltat0 self.deltat = params_ts.deltat0 self.deltat_max = params_ts.deltat_max nb_layers = self.sim.state.nb_layers if nb_layers == 2: self.compute_time_increment_CLF = ( self._compute_time_increment_CLF_2layers ) else: raise NotImplementedError
[docs] def _compute_time_increment_CLF_2layers(self): """Compute the time increment deltat with a CLF condition: .. math:: dt = CFL \times \max( \min(dx,dy)/\sqrt{g (h_1+h_2)} ) """ # barotropic CFL h0 = self.sim.state.get_var("h0") h1 = self.sim.state.get_var("h1") h = h0 + h1 g = self.sim.params.g tmp = sqrt(g * h.max()) / min(self.sim.oper.deltax, self.sim.oper.deltay) if mpi.nb_proc > 1: tmp = mpi.comm.allreduce(tmp, op=mpi.MPI.MAX) if tmp > 0: deltat_CFL = self.CFL / tmp else: deltat_CFL = self.deltat_max maybe_new_dt = min(deltat_CFL, self.deltat_max) # 4 significant numbers tmp = floor(log10(maybe_new_dt)) - 3 maybe_new_dt = floor(maybe_new_dt * 10 ** -tmp) * 10 ** tmp normalize_diff = abs(self.deltat - maybe_new_dt) / maybe_new_dt if normalize_diff > 0.02: self.deltat = maybe_new_dt