Merge pull request #1558 from OceanParcels/error_when_dt_too_high_precision

Throwing an error when dt has too high precision

Author: erikvansebille

docs/examples/example_globcurrent.py

@@ -173,14 +173,6 @@ def test_globcurrent_time_extrapolation_error(mode, use_xarray):
     pset.execute(AdvectionRK4, runtime=delta(days=1), dt=delta(minutes=5))
 
 
-@pytest.mark.parametrize('mode', ['scipy', 'jit'])
-@pytest.mark.parametrize('use_xarray', [True, False])
-def test_globcurrent_dt0(mode, use_xarray):
-    fieldset = set_globcurrent_fieldset(use_xarray=use_xarray)
-    pset = ParticleSet(fieldset, pclass=ptype[mode], lon=[25], lat=[-35])
-    pset.execute(AdvectionRK4, dt=0.)
-
-
 @pytest.mark.parametrize('mode', ['scipy', 'jit'])
 @pytest.mark.parametrize('dt', [-300, 300])
 @pytest.mark.parametrize('with_starttime', [True, False])
@@ -260,3 +252,22 @@ def test_globcurrent_pset_fromfile(mode, dt, pid_offset, tmpdir):
 
     for var in ['lon', 'lat', 'depth', 'time', 'id']:
         assert np.allclose([getattr(p, var) for p in pset], [getattr(p, var) for p in pset_new])
+
+
+@pytest.mark.parametrize('mode', ['scipy', 'jit'])
+def test_error_outputdt_not_multiple_dt(mode, tmpdir):
+    # Test that outputdt is a multiple of dt
+    fieldset = set_globcurrent_fieldset()
+
+    filepath = tmpdir.join("pfile_error_outputdt_not_multiple_dt.zarr")
+
+    dt = 81.2584344538292  # number for which output writing fails
+
+    pset = ParticleSet(fieldset, pclass=ptype[mode], lon=[0], lat=[0])
+    ofile = pset.ParticleFile(name=filepath, outputdt=delta(days=1))
+
+    def DoNothing(particle, fieldset, time):
+        pass
+
+    with pytest.raises(ValueError):
+        pset.execute(DoNothing, runtime=delta(days=10), dt=dt, output_file=ofile)

parcels/particleset.py

@@ -878,8 +878,10 @@ def execute(self, pyfunc=AdvectionRK4, pyfunc_inter=None, endtime=None, runtime=
             runtime = runtime.total_seconds()
         if isinstance(dt, delta):
             dt = dt.total_seconds()
-        if dt > 0 and dt <= 1e-6:
+        if abs(dt) <= 1e-6:
             raise ValueError('Time step dt is too small')
+        if (dt * 1e6) % 1 != 0:
+            raise ValueError('Output interval should not have finer precision than 1e-6 s')
         outputdt = output_file.outputdt if output_file else np.infty
         if isinstance(outputdt, delta):
             outputdt = outputdt.total_seconds()

tests/test_interaction.py

@@ -169,7 +169,7 @@ def f(p):
 
 @pytest.mark.parametrize('runtime, dt',
                          [(1, 1e-2),
-                          (1, -2.1234e-3),
+                          (1, -2.123e-3),
                           (1, -3.12452-3)])
 def test_pseudo_interaction(runtime, dt):
     # A linear field where advected particles are moving at