Merge pull request #682 from KulaginVladimir/main

Fix bug with calculations on small timescales

festim/exports/derived_quantities/derived_quantities.py

@@ -123,7 +123,7 @@ def is_export(self, t, final_time, nb_iterations):
             nb_its_between_exports = self.nb_iterations_between_exports
             if nb_its_between_exports is None:
                 # export at the end
-                return np.isclose(t, final_time)
+                return np.isclose(t, final_time, atol=0)
             else:
                 # export every N iterations
                 return nb_iterations % nb_its_between_exports == 0

festim/exports/txt_export.py

@@ -47,7 +47,7 @@ def is_it_time_to_export(self, current_time):
         if self.times is None:
             return True
         for time in self.times:
-            if np.isclose(time, current_time):
+            if np.isclose(time, current_time, atol=0):
                 return True
         return False
 

festim/generic_simulation.py

@@ -320,7 +320,7 @@ def run_transient(self):
         #  Time-stepping
         print("Time stepping...")
         while self.t < self.settings.final_time and not np.isclose(
-            self.t, self.settings.final_time
+            self.t, self.settings.final_time, atol=0
         ):
             self.iterate()
 
@@ -370,7 +370,10 @@ def display_time(self):
         msg = "{:.1f} %        ".format(simulation_percentage)
         msg += "{:.1e} s".format(self.t)
         msg += "    Ellapsed time so far: {:.1f} s".format(elapsed_time)
-        if not np.isclose(self.t, self.settings.final_time) and self.log_level == 40:
+        if (
+            not np.isclose(self.t, self.settings.final_time, atol=0)
+            and self.log_level == 40
+        ):
             print(msg, end="\r")
         else:
             print(msg)

festim/stepsize.py

@@ -95,7 +95,7 @@ def adapt(self, t, nb_it, converged):
         next_milestone = self.next_milestone(t)
         if next_milestone is not None:
             if t + float(self.value) > next_milestone and not np.isclose(
-                t, next_milestone
+                t, next_milestone, atol=0
             ):
                 self.value.assign((next_milestone - t))
 

test/system/test_misc.py

@@ -121,20 +121,25 @@ def test_wrong_value_for_bc_field(field):
         sim.initialise()
 
 
-def test_txt_export_desired_times(tmp_path):
+@pytest.mark.parametrize(
+    "final_time,stepsize,export_times",
+    [(1, 0.1, [0.2, 0.5]), (1e-7, 1e-9, [1e-8, 1.5e-8, 2e-8])],
+)
+def test_txt_export_desired_times(tmp_path, final_time, stepsize, export_times):
     """
     Tests that TXTExport can be exported at desired times
+    Also catches the bug #682
     """
     my_model = F.Simulation()
 
     my_model.mesh = F.MeshFromVertices(np.linspace(0, 1))
     my_model.materials = F.Material(1, 1, 0)
-    my_model.settings = F.Settings(1e-10, 1e-10, final_time=1)
+    my_model.settings = F.Settings(1e-10, 1e-10, final_time=final_time)
     my_model.T = F.Temperature(500)
-    my_model.dt = F.Stepsize(0.1)
+    my_model.dt = F.Stepsize(stepsize)
 
     my_export = F.TXTExport(
-        "solute", times=[0.2, 0.5], filename="{}/mobile_conc.txt".format(tmp_path)
+        "solute", times=export_times, filename="{}/mobile_conc.txt".format(tmp_path)
     )
     my_model.exports = [my_export]
 
@@ -252,3 +257,20 @@ def test_derived_quantities_exported_last_timestep_with_small_stepsize(tmp_path)
     my_model.run()
 
     assert os.path.exists(f"{tmp_path}/out.csv")
+
+
+def test_small_timesteps_final_time_bug():
+    """
+    Test to catch the bug #682
+    Runs a sim on small timescales and checks that the final time is reached
+    """
+    my_model = F.Simulation(log_level=40)
+    my_model.mesh = F.MeshFromVertices(np.linspace(0, 1, 10))
+    my_model.materials = F.Material(1, 1, 0.1)
+    my_model.T = F.Temperature(1000)
+    my_model.settings = F.Settings(1e-10, 1e-10, final_time=1e-7)
+    my_model.dt = F.Stepsize(1e-9)
+    my_model.initialise()
+    my_model.run()
+
+    assert np.isclose(my_model.t, my_model.settings.final_time, atol=0.0)