Merge remote-tracking branch 'upstream/master'

pySDC/implementations/convergence_controller_classes/adaptivity.py

@@ -52,7 +52,7 @@ def dependencies(self, controller, description, **kwargs):
         Returns:
             None
         """
-        step_limiter_keys = ['dt_min', 'dt_max', 'dt_slope_min', 'dt_slope_max']
+        step_limiter_keys = ['dt_min', 'dt_max', 'dt_slope_min', 'dt_slope_max', 'dt_rel_min_slope']
         available_keys = [me for me in step_limiter_keys if me in self.params.__dict__.keys()]
 
         if len(available_keys) > 0:

pySDC/implementations/convergence_controller_classes/step_size_limiter.py

@@ -39,7 +39,7 @@ def dependencies(self, controller, description, **kwargs):
         Returns:
             None
         """
-        slope_limiter_keys = ['dt_slope_min', 'dt_slope_max']
+        slope_limiter_keys = ['dt_slope_min', 'dt_slope_max', 'dt_rel_min_slope']
         available_keys = [me for me in slope_limiter_keys if me in self.params.__dict__.keys()]
 
         if len(available_keys) > 0:
@@ -90,7 +90,9 @@ class StepSizeSlopeLimiter(ConvergenceController):
     """
     Class to set limits to adaptive step size computation during run time
 
-    Please supply dt_min or dt_max in the params to limit in either direction
+    Please supply `dt_slope_min` or `dt_slope_max` in the params to limit in either direction.
+    You can also supply `dt_rel_min_slope` in order to keep the old step size in case the relative change is smaller
+    than this minimum.
     """
 
     def setup(self, controller, params, description, **kwargs):
@@ -109,6 +111,7 @@ def setup(self, controller, params, description, **kwargs):
             "control_order": 91,
             "dt_slope_min": 0,
             "dt_slope_max": np.inf,
+            "dt_rel_min_slope": 0,
         }
         return {**defaults, **super().setup(controller, params, description, **kwargs)}
 
@@ -143,5 +146,11 @@ def get_new_step_size(self, controller, S, **kwargs):
                         S,
                     )
                     L.status.dt_new = dt_new
+                elif abs(L.status.dt_new / L.params.dt - 1) < self.params.dt_rel_min_slope:
+                    L.status.dt_new = L.params.dt
+                    self.log(
+                        f"Step size did not change sufficiently to warrant step size change, keeping {L.status.dt_new:.2e}",
+                        S,
+                    )
 
         return None

pySDC/implementations/datatype_classes/mesh.py

@@ -70,12 +70,11 @@ def __abs__(self):
             float: absolute maximum of all mesh values
         """
         # take absolute values of the mesh values
-        local_absval = float(np.amax(np.ndarray.__abs__(self)))
+        local_absval = float(np.max(np.ndarray.__abs__(self)))
 
         if self.comm is not None:
-            if self.comm.Get_size() > 1:
-                global_absval = 0.0
-                global_absval = max(self.comm.allreduce(sendobj=local_absval, op=MPI.MAX), global_absval)
+            if self.comm.size > 1:
+                global_absval = self.comm.allreduce(sendobj=local_absval, op=MPI.MAX)
             else:
                 global_absval = local_absval
         else:

pySDC/implementations/hooks/log_solution.py

@@ -81,13 +81,13 @@ class LogToFile(Hooks):
 
     Keep in mind that the hook will overwrite files without warning!
     You can give a custom file name by setting the ``file_name`` class attribute and give a custom way of rendering the
-    index associated with individual files by giving a different lambda function ``format_index`` class attribute. This
-    lambda should accept one index and return one string.
+    index associated with individual files by giving a different function ``format_index`` class attribute. This should
+    accept one index and return one string.
 
     You can also give a custom ``logging_condition`` function, accepting the current level if you want to log selectively.
 
     Importantly, you may need to change ``process_solution``. By default, this will return a numpy view of the solution.
-    Of course, if you are not using numpy, you need to change this. Again, this is a lambda accepting the level.
+    Of course, if you are not using numpy, you need to change this. Again, this is a function accepting the level.
 
     After the fact, you can use the classmethod `get_path` to get the path to a certain data or the `load` function to
     directly load the solution at a given index. Just configure the hook like you did when you recorded the data
@@ -99,6 +99,7 @@ class LogToFile(Hooks):
 
     path = None
     file_name = 'solution'
+    counter = 0
 
     def logging_condition(L):
         return True
@@ -111,7 +112,6 @@ def format_index(index):
 
     def __init__(self):
         super().__init__()
-        self.counter = 0
 
         if self.path is None:
             raise ValueError('Please set a path for logging as the class attribute `LogToFile.path`!')
@@ -124,20 +124,41 @@ def __init__(self):
         if not os.path.isdir(self.path):
             os.mkdir(self.path)
 
-    def post_step(self, step, level_number):
+    def log_to_file(self, step, level_number, condition, process_solution=None):
         if level_number > 0:
             return None
 
         L = step.levels[level_number]
 
-        if type(self).logging_condition(L):
+        if condition:
             path = self.get_path(self.counter)
-            data = type(self).process_solution(L)
+
+            if process_solution:
+                data = process_solution(L)
+            else:
+                data = type(self).process_solution(L)
 
             with open(path, 'wb') as file:
                 pickle.dump(data, file)
+            self.logger.info(f'Stored file {path!r}')
+
+            type(self).counter += 1
+
+    def post_step(self, step, level_number):
+        L = step.levels[level_number]
+        self.log_to_file(step, level_number, type(self).logging_condition(L))
+
+    def pre_run(self, step, level_number):
+        L = step.levels[level_number]
+        L.uend = L.u[0]
+
+        def process_solution(L):
+            return {
+                **type(self).process_solution(L),
+                't': L.time,
+            }
 
-            self.counter += 1
+        self.log_to_file(step, level_number, True, process_solution=process_solution)
 
     @classmethod
     def get_path(cls, index):
@@ -148,3 +169,22 @@ def load(cls, index):
         path = cls.get_path(index)
         with open(path, 'rb') as file:
             return pickle.load(file)
+
+
+class LogToFileAfterXs(LogToFile):
+    r'''
+    Log to file after certain amount of time has passed instead of after every step
+    '''
+
+    time_increment = 0
+    t_next_log = 0
+
+    def post_step(self, step, level_number):
+        L = step.levels[level_number]
+
+        if self.t_next_log == 0:
+            self.t_next_log = self.time_increment
+
+        if L.time + L.dt >= self.t_next_log and not step.status.restart:
+            super().post_step(step, level_number)
+            self.t_next_log = max([L.time + L.dt, self.t_next_log]) + self.time_increment

pySDC/implementations/sweeper_classes/Runge_Kutta.py

@@ -477,7 +477,7 @@ class BackwardEuler(RungeKutta):
     nodes, weights, matrix = generator.genCoeffs()
 
 
-class CrankNicholson(RungeKutta):
+class CrankNicolson(RungeKutta):
     """
     Implicit Runge-Kutta method of second order, A-stable.
     """

pySDC/playgrounds/Allen_Cahn/AllenCahn_contracting_circle_standard_integrators.py

@@ -112,7 +112,7 @@ def run_imex_Euler(t0, dt, Tend):
     return err, radius, exact_radius
 
 
-def run_CrankNicholson(t0, dt, Tend):
+def run_CrankNicolson(t0, dt, Tend):
     """
     Routine to run particular SDC variant
 
@@ -202,8 +202,8 @@ def main_radius(cwd=''):
     radii['implicit-Euler'] = radius
     _, radius, exact_radius = run_imex_Euler(t0=t0, dt=dt, Tend=Tend)
     radii['imex-Euler'] = radius
-    _, radius, exact_radius = run_CrankNicholson(t0=t0, dt=dt, Tend=Tend)
-    radii['CrankNicholson'] = radius
+    _, radius, exact_radius = run_CrankNicolson(t0=t0, dt=dt, Tend=Tend)
+    radii['CrankNicolson'] = radius
 
     xcoords = [t0 + i * dt for i in range(int((Tend - t0) / dt))]
     plot_radius(xcoords, exact_radius, radii)
@@ -219,8 +219,8 @@ def main_error(cwd=''):
     # errors['implicit-Euler'] = err
     # err, _, _ = run_imex_Euler(t0=t0, dt=0.001/512, Tend=Tend)
     # errors['imex-Euler'] = err
-    err, _, _ = run_CrankNicholson(t0=t0, dt=0.001 / 64, Tend=Tend)
-    errors['CrankNicholson'] = err
+    err, _, _ = run_CrankNicolson(t0=t0, dt=0.001 / 64, Tend=Tend)
+    errors['CrankNicolson'] = err
 
 
 if __name__ == "__main__":

pySDC/projects/DAE/sweepers/rungeKuttaDAE.py

@@ -2,7 +2,7 @@
 from pySDC.implementations.sweeper_classes.Runge_Kutta import (
     RungeKutta,
     BackwardEuler,
-    CrankNicholson,
+    CrankNicolson,
     EDIRK4,
     DIRK43_2,
 )
@@ -171,7 +171,7 @@ class BackwardEulerDAE(RungeKuttaDAE, BackwardEuler):
     pass
 
 
-class TrapezoidalRuleDAE(RungeKuttaDAE, CrankNicholson):
+class TrapezoidalRuleDAE(RungeKuttaDAE, CrankNicolson):
     pass
 
 

pySDC/tests/test_convergence_controllers/test_adaptivity.py

@@ -13,7 +13,6 @@ def get_controller(dt, num_nodes, useMPI, adaptivity, adaptivity_params, **kwarg
         adaptivity_params (dict): Parameters for convergence controller
 
     Returns:
-       (dict): Stats object generated during the run
        (pySDC.Controller.controller): Controller used in the run
     """
     from pySDC.implementations.problem_classes.polynomial_test_problem import polynomial_testequation

pySDC/tests/test_convergence_controllers/test_step_size_limiter.py

@@ -0,0 +1,112 @@
+import pytest
+
+
+def get_controller(step_size_limier_params):
+    """
+    Runs a single advection problem with certain parameters
+
+    Args:
+        step_size_limier_params (dict): Parameters for convergence controller
+
+    Returns:
+       (pySDC.Controller.controller): Controller used in the run
+    """
+    from pySDC.implementations.problem_classes.polynomial_test_problem import polynomial_testequation
+    from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI
+    from pySDC.implementations.sweeper_classes.generic_implicit import generic_implicit
+    from pySDC.implementations.convergence_controller_classes.step_size_limiter import StepSizeLimiter
+
+    level_params = {}
+    level_params['dt'] = 1.0
+    level_params['restol'] = 1.0
+
+    sweeper_params = {}
+    sweeper_params['quad_type'] = 'GAUSS'
+    sweeper_params['num_nodes'] = 1
+    sweeper_params['do_coll_update'] = True
+
+    problem_params = {'degree': 10}
+
+    step_params = {}
+    step_params['maxiter'] = 0
+
+    controller_params = {}
+    controller_params['logger_level'] = 30
+
+    description = {}
+    description['problem_class'] = polynomial_testequation
+    description['problem_params'] = problem_params
+    description['sweeper_class'] = generic_implicit
+    description['sweeper_params'] = sweeper_params
+    description['level_params'] = level_params
+    description['step_params'] = step_params
+    description['convergence_controllers'] = {StepSizeLimiter: step_size_limier_params}
+
+    controller = controller_nonMPI(num_procs=1, controller_params=controller_params, description=description)
+
+    controller.add_convergence_controller(StepSizeLimiter, description, step_size_limier_params)
+
+    return controller
+
+
+@pytest.mark.base
+def test_step_size_slope_limiter():
+    from pySDC.implementations.convergence_controller_classes.step_size_limiter import StepSizeSlopeLimiter
+
+    params = {'dt_slope_max': 2, 'dt_slope_min': 1e-3, 'dt_rel_min_slope': 1e-1}
+    controller = get_controller(params)
+
+    limiter = controller.convergence_controllers[
+        [type(me) for me in controller.convergence_controllers].index(StepSizeSlopeLimiter)
+    ]
+
+    S = controller.MS[0]
+    S.status.slot = 0
+    L = S.levels[0]
+    L.status.time = 0
+
+    L.params.dt = 1
+    L.status.dt_new = 3
+    limiter.get_new_step_size(controller, S)
+    assert L.status.dt_new == 2
+
+    L.params.dt = 1
+    L.status.dt_new = 0
+    limiter.get_new_step_size(controller, S)
+    assert L.status.dt_new == 1e-3
+
+    L.params.dt = 1
+    L.status.dt_new = 1 + 1e-3
+    limiter.get_new_step_size(controller, S)
+    assert L.status.dt_new == 1
+
+
+@pytest.mark.base
+def test_step_size_limiter():
+    from pySDC.implementations.convergence_controller_classes.step_size_limiter import StepSizeLimiter
+
+    params = {'dt_max': 2, 'dt_min': 0.5}
+    controller = get_controller(params)
+
+    limiter = controller.convergence_controllers[
+        [type(me) for me in controller.convergence_controllers].index(StepSizeLimiter)
+    ]
+
+    S = controller.MS[0]
+    S.status.slot = 0
+    L = S.levels[0]
+    L.status.time = 0
+
+    L.params.dt = 1
+    L.status.dt_new = 3
+    limiter.get_new_step_size(controller, S)
+    assert L.status.dt_new == 2
+
+    L.params.dt = 1
+    L.status.dt_new = 0
+    limiter.get_new_step_size(controller, S)
+    assert L.status.dt_new == 0.5
+
+
+if __name__ == '__main__':
+    test_step_size_slope_limiter()

pySDC/tests/test_hooks/test_log_to_file.py

@@ -32,7 +32,7 @@ def run(hook, Tend=0):
         u0 = prob.u_exact(0)
 
         _, stats = controller.run(u0, 0, Tend)
-        return stats
+        return u0, stats
 
 
 @pytest.mark.base
@@ -68,8 +68,8 @@ def test_logging():
     LogToFile.path = path
     Tend = 2
 
-    stats = run([LogToFile, LogSolution], Tend=Tend)
-    u = get_sorted(stats, type='u')
+    u0, stats = run([LogToFile, LogSolution], Tend=Tend)
+    u = [(0.0, u0)] + get_sorted(stats, type='u')
 
     u_file = []
     for i in range(len(u)):