Merge pull request festim-dev#698 from KulaginVladimir/Stepsize

`max_stepsize` attribute for `festim.Stepsize`

docs/source/userguide/stepsize.rst

@@ -28,6 +28,14 @@ To cap the stepsize after some time, the parameters ``t_stop`` and ``stepsize_st
 
     my_stepsize = F.Stepsize(initial_value=1.2, stepsize_change_ratio=1.5, dt_min=1e-6, t_stop=10, stepsize_stop_max=1.5)
 
+.. warning::
+
+    Please note that the parameters ``t_stop`` and ``stepsize_stop_max`` will be deprecated in a future release. To set the maximal value of the stepsize, consider using the ``max_stepsize`` parameter:
+
+    .. code-block:: python
+
+        my_stepsize = F.Stepsize(initial_value=1.2, stepsize_change_ratio=1.5, dt_min=1e-6, max_stepsize=lambda t: 1 if t < 1 else 2)
+
 The ``milestones`` argument can be used to make sure the simulation passes through specific times.
 This will modify the stepsize as needed.
 

festim/stepsize.py

@@ -1,5 +1,6 @@
 import fenics as f
 import numpy as np
+import warnings
 
 
 class Stepsize:
@@ -14,6 +15,8 @@ class Stepsize:
             stops. Defaults to None.
         stepsize_stop_max (float, optional): Maximum stepsize after
             t_stop. Defaults to None.
+        max_stepsize (float or callable, optional): Maximum stepsize.
+            Can be a function of festim.t. Defaults to None.
         dt_min (float, optional): Minimum stepsize below which an error is
             raised. Defaults to None.
         milestones (list, optional): list of times by which the simulation must
@@ -24,6 +27,15 @@ class Stepsize:
         value (fenics.Constant): value of dt
         milestones (list): list of times by which the simulation must
             pass.
+
+    Example::
+
+        my_stepsize = Stepsize(
+            initial_value=0.5,
+            stepsize_change_ratio=1.1,
+            max_stepsize=lambda t: None if t < 1 else 2,
+            dt_min=1e-05
+        )
     """
 
     def __init__(
@@ -32,15 +44,21 @@ def __init__(
         stepsize_change_ratio=None,
         t_stop=None,
         stepsize_stop_max=None,
+        max_stepsize=None,
         dt_min=None,
         milestones=None,
     ) -> None:
         self.adaptive_stepsize = None
         if stepsize_change_ratio is not None:
+            if t_stop or stepsize_stop_max:
+                warnings.warn(
+                    "stepsize_stop_max and t_stop attributes will be deprecated in a future release, please use max_stepsize instead",
+                    DeprecationWarning,
+                )
+                max_stepsize = lambda t: stepsize_stop_max if t >= t_stop else None
             self.adaptive_stepsize = {
                 "stepsize_change_ratio": stepsize_change_ratio,
-                "t_stop": t_stop,
-                "stepsize_stop_max": stepsize_stop_max,
+                "max_stepsize": max_stepsize,
                 "dt_min": dt_min,
             }
         self.initial_value = initial_value
@@ -75,8 +93,8 @@ def adapt(self, t, nb_it, converged):
         if self.adaptive_stepsize:
             change_ratio = self.adaptive_stepsize["stepsize_change_ratio"]
             dt_min = self.adaptive_stepsize["dt_min"]
-            stepsize_stop_max = self.adaptive_stepsize["stepsize_stop_max"]
-            t_stop = self.adaptive_stepsize["t_stop"]
+            max_stepsize = self.adaptive_stepsize["max_stepsize"]
+
             if not converged:
                 self.value.assign(float(self.value) / change_ratio)
                 if float(self.value) < dt_min:
@@ -86,10 +104,11 @@ def adapt(self, t, nb_it, converged):
             else:
                 self.value.assign(float(self.value) / change_ratio)
 
-            if t_stop is not None:
-                if t >= t_stop:
-                    if float(self.value) > stepsize_stop_max:
-                        self.value.assign(stepsize_stop_max)
+            if callable(max_stepsize):
+                max_stepsize = max_stepsize(t)
+            if max_stepsize is not None:
+                if float(self.value) > max_stepsize:
+                    self.value.assign(max_stepsize)
 
         # adapt for next milestone
         next_milestone = self.next_milestone(t)

test/unit/test_stepsize.py

@@ -30,13 +30,27 @@ def test_value_is_reduced(self, my_stepsize):
             == old_value / my_stepsize.adaptive_stepsize["stepsize_change_ratio"]
         )
 
-    def test_hit_stepsize_max(self, my_stepsize):
+    def test_hit_stepsize_max_with_float(self, my_stepsize):
+        """
+        Assigns an initial value to the stepsize, then calls adapt at t > 0
+        and checks that the new value is equal to max_stepsize
+        """
         my_stepsize.value.assign(10)
-        my_stepsize.adaptive_stepsize["stepsize_stop_max"] = 1
-        my_stepsize.adaptive_stepsize["t_stop"] = 0
+        my_stepsize.adaptive_stepsize["max_stepsize"] = 1
         my_stepsize.adapt(t=6, converged=True, nb_it=2)
         new_value = float(my_stepsize.value)
-        assert new_value == my_stepsize.adaptive_stepsize["stepsize_stop_max"]
+        assert new_value == my_stepsize.adaptive_stepsize["max_stepsize"]
+
+    def test_hit_stepsize_max_with_callable(self, my_stepsize):
+        """
+        Assigns an initial value to the stepsize and a callable for max_stepsize
+        and checks that the new value is equal to max_stepsize
+        """
+        my_stepsize.value.assign(10)
+        my_stepsize.adaptive_stepsize["max_stepsize"] = lambda t: t
+        my_stepsize.adapt(t=6, converged=True, nb_it=2)
+        new_value = float(my_stepsize.value)
+        assert new_value == my_stepsize.adaptive_stepsize["max_stepsize"](6)
 
 
 def test_milestones_are_hit():
@@ -86,3 +100,33 @@ def test_next_milestone():
             if expected_milestone is not None
             else next_milestone is None
         )
+
+
+def test_DeprecationWarning_t_stop():
+    """A temporary test to check DeprecationWarning in festim.Stepsize"""
+
+    with pytest.deprecated_call():
+        festim.Stepsize(
+            initial_value=1e-8,
+            stepsize_change_ratio=2,
+            dt_min=1,
+            t_stop=0,
+            stepsize_stop_max=1,
+        )
+
+
+@pytest.mark.parametrize("time", (0, 2))
+def test_hit_stepsize_max_with_t_stop(time):
+    """
+    A temporary test to check that when old attributes t_stop and stepsize_stop_max
+    are used their work is re-created with max_stepsize
+    """
+    my_stepsize = festim.Stepsize(
+        initial_value=10,
+        stepsize_change_ratio=2,
+        dt_min=0.1,
+        t_stop=1,
+        stepsize_stop_max=1,
+    )
+    max_stepsize = lambda t: 1 if t >= 1 else None
+    assert my_stepsize.adaptive_stepsize["max_stepsize"](time) == max_stepsize(time)