Repository: GazzolaLab/PyElastica

Issue:
Current-step passed to callback is determined by floating point division.
<a href="https://github.com/skim0119"><img src="https://avatars.githubusercontent.com/u/3798023?v=4" align="left" width="96" height="96" hspace="10"></img></a> **Issue by [skim0119](https://github.com/skim0119)**
_Friday Dec 31, 2021 at 21:31 GMT_
_Originally opened as https://github.com/GazzolaLab/elastica-python/issues/116_

----

https://github.com/GazzolaLab/PyElastica/blob/d5c18c3750fa64681150c73124fd349f533f2c48/elastica/timestepper/symplectic_steppers.py#L108

During the operation `int(time/dt)`, we loose up to 1 step error. If user write the callback such that the data is recorded every 100 steps, some of the data will be recorded either 99 steps.


Discussion:

Patch:
diff --git a/elastica/timestepper/symplectic_steppers.py b/elastica/timestepper/symplectic_steppers.py
index c268abe..b87483b 100644
--- a/elastica/timestepper/symplectic_steppers.py
+++ b/elastica/timestepper/symplectic_steppers.py
@@ -1,6 +1,7 @@
 __doc__ = """Symplectic time steppers and concepts for integrating the kinematic and dynamic equations of rod-like objects.  """
 
 import numpy as np
+import math
 
 # from elastica._elastica_numba._timestepper._symplectic_steppers import (
 #     SymplecticStepperTag,
@@ -105,7 +106,7 @@ class _SystemCollectionStepper:
         SystemCollection.constrain_values(time)
 
         # Call back function, will call the user defined call back functions and store data
-        SystemCollection.apply_callbacks(time, int(time / dt))
+        SystemCollection.apply_callbacks(time, round(time / dt))
 
         # Zero out the external forces and torques
         for system in SystemCollection._memory_blocks:
diff --git a/elastica/wrappers/base_system.py b/elastica/wrappers/base_system.py
index 392086f..c4396d4 100644
--- a/elastica/wrappers/base_system.py
+++ b/elastica/wrappers/base_system.py
@@ -131,6 +131,10 @@ class BaseSystemCollection(MutableSequence):
 
         """
 
+        # FIXME: This is probably the most bizarre way to collect the functions. It is definitely
+        # impressive and even surprising that it is working, but it is far from readable and maintainable.
+        # The code is difficult to debug, because the behavior is un-interpretable except during run-time.
+        # We need a lot more documentation on this part: clear explanation and reasoning.
         def get_methods_from_feature_classes(method_name: str):
             methods = [
                 [v for (k, v) in cls.__dict__.items() if k.endswith(method_name)]
@@ -146,7 +150,9 @@ class BaseSystemCollection(MutableSequence):
         self._features_that_constrain_rates = get_methods_from_feature_classes(
             "_constrain_rates"
         )
-        self._callback_features = get_methods_from_feature_classes("_callBack")
+        self._callback_features = get_methods_from_feature_classes(
+            "_callback_execution"
+        )
         finalize_methods = get_methods_from_feature_classes("_finalize")
 
         # construct memory block
diff --git a/elastica/wrappers/callbacks.py b/elastica/wrappers/callbacks.py
index d5d51d8..c476b18 100644
--- a/elastica/wrappers/callbacks.py
+++ b/elastica/wrappers/callbacks.py
@@ -16,12 +16,12 @@ class CallBacks:
 
         Attributes
         ----------
-        _callbacks: list
+        _callback_list: list
             List of call back classes defined for rod-like objects.
     """
 
     def __init__(self):
-        self._callbacks = []
+        self._callback_list = []
         super(CallBacks, self).__init__()
 
     def collect_diagnostics(self, system):
@@ -43,7 +43,7 @@ class CallBacks:
 
         # Create _Constraint object, cache it and return to user
         _callbacks = _CallBack(sys_idx)
-        self._callbacks.append(_callbacks)
+        self._callback_list.append(_callbacks)
 
         return _callbacks
 
@@ -54,9 +54,9 @@ class CallBacks:
         # dev : the first index stores the rod index to collect data.
         # Technically we can use another array but it its one more book-keeping
         # step. Being lazy, I put them both in the same array
-        self._callbacks[:] = [
+        self._callback_list[:] = [
             (callback.id(), callback(self._systems[callback.id()]))
-            for callback in self._callbacks
+            for callback in self._callback_list
         ]
 
         # Sort from lowest id to highest id for potentially better memory access
@@ -65,13 +65,12 @@ class CallBacks:
         # [(0, MyCallBack), (1, MyVelocityCallBack), ... ]
         # Thus using lambda we iterate over the list of tuples and use rod number (x[0])
         # to sort callbacks.
-        self._callbacks.sort(key=lambda x: x[0])
+        self._callback_list.sort(key=lambda x: x[0])
 
-        self._callBack(time=0.0, current_step=0)
+        self._callback_execution(time=0.0, current_step=0)
 
-    # TODO: same as above naming of _callBack function
-    def _callBack(self, time, current_step: int, *args, **kwargs):
-        for sys_id, callback in self._callbacks:
+    def _callback_execution(self, time, current_step: int, *args, **kwargs):
+        for sys_id, callback in self._callback_list:
             callback.make_callback(
                 self._systems[sys_id], time, current_step, *args, **kwargs
             )
@@ -97,6 +96,7 @@ class _CallBack:
         Parameters
         ----------
         sys_idx: int
+            rod object index
         """
         self._sys_idx = sys_idx
         self._callback_cls = None
@@ -112,10 +112,6 @@ class _CallBack:
         ----------
         callback_cls: object
             User defined callback class.
-        *args
-            Variable length argument list
-        **kwargs
-            Arbitrary keyword arguments.
 
         Returns
         -------
@@ -134,18 +130,8 @@ class _CallBack:
     def id(self):
         return self._sys_idx
 
-    def __call__(self, *args, **kwargs):
-        """Constructs a callback functions after checks
-
-        Parameters
-        ----------
-        args
-        kwargs
-
-        Returns
-        -------
-
-        """
+    def __call__(self, *args, **kwargs) -> CallBackBaseClass:
+        """Constructs a callback functions after checks"""
         if not self._callback_cls:
             raise RuntimeError(
                 "No callback provided to act on rod id {0}"


Test patch:
diff --git a/tests/test_wrappers/test_base_system.py b/tests/test_wrappers/test_base_system.py
index cbd9288..4003801 100644
--- a/tests/test_wrappers/test_base_system.py
+++ b/tests/test_wrappers/test_base_system.py
@@ -203,7 +203,7 @@ class TestBaseSystemWithFeaturesUsingCosseratRod:
         simulator_class.collect_diagnostics(rod).using(legal_callback)
         simulator_class.finalize()
         # After finalize check if the created callback object is instance of the class we have given.
-        assert isinstance(simulator_class._callbacks[-1][-1], legal_callback)
+        assert isinstance(simulator_class._callback_list[-1][-1], legal_callback)
 
         # TODO: this is a dummy test for apply_callbacks find a better way to test them
         simulator_class.apply_callbacks(time=0, current_step=0)
diff --git a/tests/test_wrappers/test_callbacks.py b/tests/test_wrappers/test_callbacks.py
index 5e94ec9..02238bf 100644
--- a/tests/test_wrappers/test_callbacks.py
+++ b/tests/test_wrappers/test_callbacks.py
@@ -131,7 +131,7 @@ class TestCallBacksMixin:
         scwc.append(mock_rod)
 
         _mock_callback = scwc.collect_diagnostics(mock_rod)
-        assert _mock_callback in scwc._callbacks
+        assert _mock_callback in scwc._callback_list
         assert _mock_callback.__class__ == _CallBack
 
     from elastica.callback_functions import CallBackBaseClass
@@ -165,7 +165,7 @@ class TestCallBacksMixin:
 
         scwc._finalize()
 
-        for (x, y) in scwc._callbacks:
+        for (x, y) in scwc._callback_list:
             assert type(x) is int
             assert type(y) is callback_cls
 
@@ -177,6 +177,6 @@ class TestCallBacksMixin:
 
         # this is allowed to fail (not critical)
         num = -np.inf
-        for (x, _) in scwc._callbacks:
+        for (x, _) in scwc._callback_list:
             assert num < x
             num = x