Merge branch 'master' into v5

Author: pancetta

pySDC/core/Controller.py

@@ -100,6 +100,19 @@ def __setup_custom_logger(level=None, log_to_file=None, fname=None):
         else:
             pass
 
+    def welcome_message(self):
+        out = 'Welcome to the one and only, really very astonishing and 87.3% bug free' + \
+              '\n                                 _____ _____   _____ ' + \
+              '\n                                / ____|  __ \ / ____|' + \
+              '\n                    _ __  _   _| (___ | |  | | |     ' + \
+              '\n                   | \'_ \| | | |\___ \| |  | | |     ' + \
+              '\n                   | |_) | |_| |____) | |__| | |____ ' + \
+              '\n                   | .__/ \__, |_____/|_____/ \_____|' + \
+              '\n                   | |     __/ |                     ' + \
+              '\n                   |_|    |___/                      ' + \
+              '\n                                                     '
+        self.logger.info(out)
+
     def dump_setup(self, step, controller_params, description):
         """
         Helper function to dump the setup used for this controller
@@ -110,6 +123,7 @@ def dump_setup(self, step, controller_params, description):
             description (dict): description of the problem
         """
 
+        self.welcome_message()
         out = 'Setup overview (--> user-defined) -- BEGIN'
         self.logger.info(out)
         out = '----------------------------------------------------------------------------------------------------\n\n'

pySDC/implementations/convergence_controller_classes/adaptivity.py

@@ -67,7 +67,7 @@ def get_new_step_size(self, controller, S):
         """
         raise NotImplementedError('Please implement a rule for updating the step size!')
 
-    def compute_optimatal_step_size(self, beta, dt, e_tol, e_est, order):
+    def compute_optimal_step_size(self, beta, dt, e_tol, e_est, order):
         """
         Compute the optimal step size for the current step based on the order of the scheme.
         This function can be called from `get_new_step_size` for various implementations of adaptivity, but notably not
@@ -199,9 +199,8 @@ def get_new_step_size(self, controller, S):
             order = S.status.iter  # embedded error estimate is same order as time marching
 
             e_est = self.get_local_error_estimate(controller, S)
-            L.status.dt_new = self.compute_optimatal_step_size(
-                self.params.beta, L.params.dt, self.params.e_tol, e_est, order
-            )
+            L.status.dt_new = self.compute_optimal_step_size(self.params.beta, L.params.dt, self.params.e_tol, e_est,
+                                                             order)
             self.log(f'Adjusting step size from {L.params.dt:.2e} to {L.status.dt_new:.2e}', S)
 
         return None
@@ -220,6 +219,55 @@ def get_local_error_estimate(self, controller, S, **kwargs):
         return S.levels[0].status.error_embedded_estimate
 
 
+class AdaptivityRK(Adaptivity):
+    '''
+    Adaptivity for Runge-Kutta methods. Basically, we need to change the order in the step size update
+    '''
+    def check_parameters(self, controller, params, description):
+        '''
+        Check whether parameters are compatible with whatever assumptions went into the step size functions etc.
+        For adaptivity, we need to know the order of the scheme.
+
+        Args:
+            controller (pySDC.Controller): The controller
+            params (dict): The params passed for this specific convergence controller
+            description (dict): The description object used to instantiate the controller
+
+        Returns:
+            bool: Whether the parameters are compatible
+            str: The error message
+        '''
+        if 'update_order' not in params.keys():
+            return False, 'Adaptivity needs an order for the update rule! Please set some up in \
+description[\'convergence_control_params\'][\'update_order\']!'
+
+        return super(AdaptivityRK, self).check_parameters(controller, params, description)
+
+    def get_new_step_size(self, controller, S):
+        '''
+        Determine a step size for the next step from an embedded estimate of the local error of the current step.
+
+        Args:
+            controller (pySDC.Controller): The controller
+            S (pySDC.Step): The current step
+
+        Returns:
+            None
+        '''
+        # check if we performed the desired amount of sweeps
+        if S.status.iter == S.params.maxiter:
+            L = S.levels[0]
+
+            # compute next step size
+            order = self.params.update_order
+            e_est = self.get_local_error_estimate(controller, S)
+            L.status.dt_new = self.compute_optimal_step_size(self.params.beta, L.params.dt, self.params.e_tol, e_est,
+                                                             order)
+            self.log(f'Adjusting step size from {L.params.dt:.2e} to {L.status.dt_new:.2e}', S)
+
+        return None
+
+
 class AdaptivityResidual(AdaptivityBase):
     """
     Do adaptivity based on residual.

pySDC/implementations/convergence_controller_classes/estimate_embedded_error.py

@@ -3,6 +3,8 @@
 from pySDC.core.ConvergenceController import ConvergenceController, Pars
 from pySDC.implementations.convergence_controller_classes.store_uold import StoreUOld
 
+from pySDC.implementations.sweeper_classes.Runge_Kutta import RungeKutta
+
 
 class EstimateEmbeddedError(ConvergenceController):
     """
@@ -14,7 +16,7 @@ class EstimateEmbeddedError(ConvergenceController):
 
     def setup(self, controller, params, description):
         """
-        Add a default value for control order to the parameters
+        Add a default value for control order to the parameters and check if we are using a Runge-Kutta sweeper
 
         Args:
             controller (pySDC.Controller): The controller
@@ -24,11 +26,12 @@ def setup(self, controller, params, description):
         Returns:
             dict: Updated parameters
         """
-        return {'control_order': -80, **params}
+        sweeper_type = 'RK' if RungeKutta in description['sweeper_class'].__bases__ else 'SDC'
+        return {'control_order': -80, 'sweeper_type': sweeper_type, **params}
 
     def dependencies(self, controller, description):
         """
-        Load the convergence controller that stores the solution of the last sweep
+        Load the convergence controller that stores the solution of the last sweep unless we are doing Runge-Kutta
 
         Args:
             controller (pySDC.Controller): The controller
@@ -37,9 +40,32 @@ def dependencies(self, controller, description):
         Returns:
             None
         """
-        controller.add_convergence_controller(StoreUOld, description=description)
+        if RungeKutta not in description['sweeper_class'].__bases__:
+            controller.add_convergence_controller(StoreUOld, description=description)
         return None
 
+    def estimate_embedded_error_serial(self, L):
+        """
+        Estimate the serial embedded error, which may need to be modified for a parallel estimate.
+
+        Depending on the type of sweeper, the lower order solution is stored in a different place.
+
+        Args:
+            L (pySDC.level): The level
+
+        Returns:
+            dtype_u: The embedded error estimate
+        """
+        if self.params.sweeper_type == 'RK':
+            # lower order solution is stored in the second to last entry of L.u
+            return abs(L.u[-2] - L.u[-1])
+        elif self.params.sweeper_type == 'SDC':
+            # order rises by one between sweeps, making this so ridiculously easy
+            return abs(L.uold[-1] - L.u[-1])
+        else:
+            raise NotImplementedError(f'Don\'t know how to estimate embedded error for sweeper type \
+{self.params.sweeper_type}')
+
 
 class EstimateEmbeddedErrorNonMPI(EstimateEmbeddedError):
     def __init__(self, controller, params, description):
@@ -85,10 +111,9 @@ def post_iteration_processing(self, controller, S):
 level'
             )
 
-        if S.status.iter > 1:
+        if S.status.iter > 1 or self.params.sweeper_type == 'RK':
             for L in S.levels:
-                # order rises by one between sweeps, making this so ridiculously easy
-                temp = abs(L.uold[-1] - L.u[-1])
+                temp = self.estimate_embedded_error_serial(L)
                 L.status.error_embedded_estimate = max([abs(temp - self.buffers.e_em_last), np.finfo(float).eps])
 
             self.buffers.e_em_last = temp * 1.0