use timestepping

include/aspect/time_stepping/interface.h

@@ -234,6 +234,36 @@ namespace aspect
         void
         parse_parameters (ParameterHandler &prm);
 
+        /**
+         * Go through the list of all mesh refinement strategies that have been selected
+         * in the input file (and are consequently currently active) and return
+         * true if one of them has the desired type specified by the template
+         * argument.
+         *
+         * This function can only be called if the given template type (the first template
+         * argument) is a class derived from the Interface class in this namespace.
+         */
+        template <typename PluginType,
+                  typename = typename std::enable_if_t<std::is_base_of<Interface<dim>,PluginType>::value>>
+        bool
+        has_matching_plugin () const;
+
+        /**
+         * Go through the list of all mesh refinement strategies that have been selected
+         * in the input file (and are consequently currently active) and see
+         * if one of them has the type specified by the template
+         * argument or can be casted to that type. If so, return a reference
+         * to it. If no mesh refinement strategy is active that matches the
+         * given type, throw an exception.
+         *
+         * This function can only be called if the given template type (the first template
+         * argument) is a class derived from the Interface class in this namespace.
+         */
+        template <typename PluginType,
+                  typename = typename std::enable_if_t<std::is_base_of<Interface<dim>,PluginType>::value>>
+        const PluginType &
+        get_matching_plugin () const;
+
         /**
          * For the current plugin subsystem, write a connection graph of all of the
          * plugins we know about, in the format that the
@@ -247,7 +277,6 @@ namespace aspect
         void
         write_plugin_graph (std::ostream &output_stream);
 
-
         /**
          * A function that is used to register time stepping model objects in such
          * a way that the Manager can deal with all of them without having to
@@ -308,6 +337,50 @@ namespace aspect
         std::list<std::unique_ptr<Interface<dim>>> active_plugins;
     };
 
+
+
+    template <int dim>
+    template <typename PluginType, typename>
+    inline
+    bool
+    Manager<dim>::has_matching_plugin () const
+    {
+      for (typename std::list<std::unique_ptr<Interface<dim>>>::const_iterator
+           p = active_plugins.begin();
+           p != active_plugins.end(); ++p)
+        if (Plugins::plugin_type_matches<PluginType>(*(*p)))
+          return true;
+
+      return false;
+    }
+
+
+
+    template <int dim>
+    template <typename PluginType, typename>
+    inline
+    const PluginType &
+    Manager<dim>::get_matching_plugin () const
+    {
+      AssertThrow(has_matching_plugin<PluginType> (),
+                  ExcMessage("You asked TimeStepping::Manager::get_matching_plugin() for a "
+                             "plugin of type <" + boost::core::demangle(typeid(PluginType).name()) + "> "
+                             "that could not be found in the current model. Activate this "
+                             "plugin in the input file."));
+
+      for (typename std::list<std::unique_ptr<Interface<dim>>>::const_iterator
+           p = active_plugins.begin();
+           p != active_plugins.end(); ++p)
+        if (Plugins::plugin_type_matches<PluginType>(*(*p)))
+          return Plugins::get_plugin_as_type<PluginType>(*(*p));
+
+      // We will never get here, because we had the Assert above. Just to avoid warnings.
+      typename std::list<std::unique_ptr<Interface<dim>>>::const_iterator plugin;
+      return Plugins::get_plugin_as_type<PluginType>(*(*plugin));
+    }
+
+
+
     /**
      * Given a class name, a name, and a description for the parameter file, register it with the
      * aspect::TimeStepping::Manager class.

source/simulator/core.cc

@@ -37,6 +37,7 @@
 #include <aspect/simulator/assemblers/interface.h>
 #include <aspect/geometry_model/initial_topography_model/zero_topography.h>
 #include <aspect/material_model/rheology/elasticity.h>
+#include <aspect/time_stepping/repeat_on_nonlinear_fail.h>
 
 #include <deal.II/base/index_set.h>
 #include <deal.II/base/conditional_ostream.h>
@@ -2049,14 +2050,7 @@ namespace aspect
                         throw;
                       }
 
-                    // reduce timestep size and update time to "go back":
-                    const double new_time_step = 0.5*time_step;
-                    const double multiplier = (parameters.convert_to_years ? 1./year_in_seconds : 1.0);
-                    pcout << "\tReducing timestep by " << new_time_step *multiplier << '\n' << std::endl;
-
-                    time -= time_step - new_time_step;
-                    time_step = new_time_step;
-                    continue; // skip the rest of the main loop that would advance in time
+                    time_stepping_manager.template get_matching_plugin<TimeStepping::RepeatOnNonlinearFail<dim>>().nonlinear_solver_has_failed();
                   }
                 else if (parameters.nonlinear_solver_failure_strategy
                          == Parameters<dim>::NonlinearSolverFailureStrategy::abort_program)

source/time_stepping/repeat_on_nonlinear_fail.cc

@@ -35,9 +35,18 @@ namespace aspect
 
 
 
+    template <int dim>
+    void
+    RepeatOnNonlinearFail<dim>::nonlinear_solver_has_failed() const
+    {
+      nonlinear_solver_just_failed = true;
+    }
+
+
+
     template <int dim>
     std::pair<Reaction, double>
-    RepeatOnNonlinearFail<dim>::determine_reaction (const TimeStepInfo &info)
+    RepeatOnNonlinearFail<dim>::determine_reaction (const TimeStepInfo &/*info*/)
     {
       if (nonlinear_solver_just_failed)
         {