Skip to content

Removing repeatdt #1387

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Closed
erikvansebille wants to merge 6 commits into from
Closed

Removing repeatdt #1387

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
f1f2ef9
Removing support for 'repeatdt' argument in particleset
erikvansebille Jun 23, 2023
3173042
Updating tutorials to reflect removal of repeatdt argument
erikvansebille Jun 23, 2023
60687e2
Adding deprecation errors to particleset.__init__ that repeatdt is no…
erikvansebille Jun 23, 2023
9b9a487
Merge branch 'master' into removing_repeatdt
erikvansebille Jun 23, 2023
5057425
Fixing pre-commit error (D202)
erikvansebille Jun 23, 2023
93f8115
Fixing bug in test_pset_create_fromparticlefile without repeatdt
erikvansebille Jun 23, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
13 changes: 7 additions & 6 deletions docs/examples/documentation_MPI.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

8 changes: 3 additions & 5 deletions docs/examples/example_stommel.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -87,7 +87,7 @@ def AgeP(particle, fieldset, time):


def stommel_example(npart=1, mode='jit', verbose=False, method=AdvectionRK4, grid_type='A',
outfile="StommelParticle.zarr", repeatdt=None, maxage=None, write_fields=True, pset_mode='soa'):
outfile="StommelParticle.zarr", maxage=None, write_fields=True, pset_mode='soa'):
timer.fieldset = timer.Timer('FieldSet', parent=timer.stommel)
fieldset = stommel_fieldset(grid_type=grid_type)
if write_fields:
Expand All @@ -105,7 +105,7 @@ class MyParticle(ParticleClass):
p_start = Variable('p_start', dtype=np.float32, initial=fieldset.P)
age = Variable('age', dtype=np.float32, initial=0.)

pset = pset_type[pset_mode]['pset'].from_line(fieldset, size=npart, pclass=MyParticle, repeatdt=repeatdt,
pset = pset_type[pset_mode]['pset'].from_line(fieldset, size=npart, pclass=MyParticle,
start=(10e3, 5000e3), finish=(100e3, 5000e3), time=0)

if verbose:
Expand Down Expand Up @@ -166,8 +166,6 @@ def main(args=None):
help='Numerical method used for advection')
p.add_argument('-o', '--outfile', default='StommelParticle.zarr',
help='Name of output file')
p.add_argument('-r', '--repeatdt', default=None, type=int,
help='repeatdt of the ParticleSet')
p.add_argument('-a', '--maxage', default=None, type=int,
help='max age of the particles (after which particles are deleted)')
p.add_argument('-psm', '--pset_mode', choices=('soa', 'aos'), default='soa',
Expand All @@ -179,7 +177,7 @@ def main(args=None):
timer.args.stop()
timer.stommel = timer.Timer('Stommel', parent=timer.root)
stommel_example(args.particles, mode=args.mode, verbose=args.verbose, method=method[args.method],
outfile=args.outfile, repeatdt=args.repeatdt, maxage=args.maxage, write_fields=args.write_fields)
outfile=args.outfile, maxage=args.maxage, write_fields=args.write_fields)
timer.stommel.stop()
timer.root.stop()
timer.root.print_tree()
Expand Down
2,096 changes: 530 additions & 1,566 deletions docs/examples/tutorial_delaystart.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

127 changes: 57 additions & 70 deletions docs/examples/tutorial_diffusion.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

283 changes: 16 additions & 267 deletions docs/examples/tutorial_sampling.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

59 changes: 10 additions & 49 deletions parcels/particleset/baseparticleset.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -31,22 +31,10 @@ class BaseParticleSet(NDCluster):
interaction_kernel = None
fieldset = None
time_origin = None
repeat_starttime = None
repeatlon = None
repeatlat = None
repeatdepth = None
repeatpclass = None
repeatkwargs = None

def __init__(self, fieldset=None, pclass=None, lon=None, lat=None,
depth=None, time=None, repeatdt=None, lonlatdepth_dtype=None, pid_orig=None, **kwargs):
depth=None, time=None, lonlatdepth_dtype=None, pid_orig=None, **kwargs):
self._collection = None
self.repeat_starttime = None
self.repeatlon = None
self.repeatlat = None
self.repeatdepth = None
self.repeatpclass = None
self.repeatkwargs = None
self.kernel = None
self.interaction_kernel = None
self.fieldset = None
Expand Down Expand Up @@ -110,7 +98,7 @@ def __create_progressbar(self, starttime, endtime):
return pbar

@classmethod
def from_list(cls, fieldset, pclass, lon, lat, depth=None, time=None, repeatdt=None, lonlatdepth_dtype=None, **kwargs):
def from_list(cls, fieldset, pclass, lon, lat, depth=None, time=None, lonlatdepth_dtype=None, **kwargs):
"""Initialise the ParticleSet from lists of lon and lat.

Parameters
Expand All @@ -127,8 +115,6 @@ def from_list(cls, fieldset, pclass, lon, lat, depth=None, time=None, repeatdt=N
Optional list of initial depth values for particles. Default is 0m
time :
Optional list of start time values for particles. Default is fieldset.U.time[0]
repeatdt :
Optional interval (in seconds) on which to repeat the release of the ParticleSet (Default value = None)
lonlatdepth_dtype :
Floating precision for lon, lat, depth particle coordinates.
It is either np.float32 or np.float64. Default is np.float32 if fieldset.U.interp_method is 'linear'
Expand All @@ -137,10 +123,10 @@ def from_list(cls, fieldset, pclass, lon, lat, depth=None, time=None, repeatdt=N
**kwargs :
Keyword arguments passed to the particleset constructor.
"""
return cls(fieldset=fieldset, pclass=pclass, lon=lon, lat=lat, depth=depth, time=time, repeatdt=repeatdt, lonlatdepth_dtype=lonlatdepth_dtype, **kwargs)
return cls(fieldset=fieldset, pclass=pclass, lon=lon, lat=lat, depth=depth, time=time, lonlatdepth_dtype=lonlatdepth_dtype, **kwargs)

@classmethod
def from_line(cls, fieldset, pclass, start, finish, size, depth=None, time=None, repeatdt=None, lonlatdepth_dtype=None):
def from_line(cls, fieldset, pclass, start, finish, size, depth=None, time=None, lonlatdepth_dtype=None):
"""Create a particleset in the shape of a line (according to a cartesian grid).

Initialise the ParticleSet from start/finish coordinates with equidistant spacing
Expand All @@ -163,8 +149,6 @@ def from_line(cls, fieldset, pclass, start, finish, size, depth=None, time=None,
Optional list of initial depth values for particles. Default is 0m
time :
Optional start time value for particles. Default is fieldset.U.time[0]
repeatdt :
Optional interval (in seconds) on which to repeat the release of the ParticleSet (Default value = None)
lonlatdepth_dtype :
Floating precision for lon, lat, depth particle coordinates.
It is either np.float32 or np.float64. Default is np.float32 if fieldset.U.interp_method is 'linear'
Expand All @@ -174,7 +158,7 @@ def from_line(cls, fieldset, pclass, start, finish, size, depth=None, time=None,
lat = np.linspace(start[1], finish[1], size)
if type(depth) in [int, float]:
depth = [depth] * size
return cls(fieldset=fieldset, pclass=pclass, lon=lon, lat=lat, depth=depth, time=time, repeatdt=repeatdt, lonlatdepth_dtype=lonlatdepth_dtype)
return cls(fieldset=fieldset, pclass=pclass, lon=lon, lat=lat, depth=depth, time=time, lonlatdepth_dtype=lonlatdepth_dtype)

@classmethod
@abstractmethod
Expand All @@ -200,7 +184,7 @@ def monte_carlo_sample(cls, start_field, size, mode='monte_carlo'):
pass

@classmethod
def from_field(cls, fieldset, pclass, start_field, size, mode='monte_carlo', depth=None, time=None, repeatdt=None, lonlatdepth_dtype=None):
def from_field(cls, fieldset, pclass, start_field, size, mode='monte_carlo', depth=None, time=None, lonlatdepth_dtype=None):
"""Initialise the ParticleSet randomly drawn according to distribution from a field.

Parameters
Expand All @@ -219,20 +203,18 @@ def from_field(cls, fieldset, pclass, start_field, size, mode='monte_carlo', dep
Optional list of initial depth values for particles. Default is 0m
time :
Optional start time value for particles. Default is fieldset.U.time[0]
repeatdt :
Optional interval (in seconds) on which to repeat the release of the ParticleSet (Default value = None)
lonlatdepth_dtype :
Floating precision for lon, lat, depth particle coordinates.
It is either np.float32 or np.float64. Default is np.float32 if fieldset.U.interp_method is 'linear'
and np.float64 if the interpolation method is 'cgrid_velocity'
"""
lon, lat = cls.monte_carlo_sample(start_field, size, mode)

return cls(fieldset=fieldset, pclass=pclass, lon=lon, lat=lat, depth=depth, time=time, lonlatdepth_dtype=lonlatdepth_dtype, repeatdt=repeatdt)
return cls(fieldset=fieldset, pclass=pclass, lon=lon, lat=lat, depth=depth, time=time, lonlatdepth_dtype=lonlatdepth_dtype)

@classmethod
@abstractmethod
def from_particlefile(cls, fieldset, pclass, filename, restart=True, restarttime=None, repeatdt=None, lonlatdepth_dtype=None, **kwargs):
def from_particlefile(cls, fieldset, pclass, filename, restart=True, restarttime=None, lonlatdepth_dtype=None, **kwargs):
"""Initialise the ParticleSet from a netcdf ParticleFile.
This creates a new ParticleSet based on locations of all particles written
in a netcdf ParticleFile at a certain time. Particle IDs are preserved if restart=True
Expand All @@ -252,8 +234,6 @@ def from_particlefile(cls, fieldset, pclass, filename, restart=True, restarttime
time at which the Particles will be restarted. Default is the last time written.
Alternatively, restarttime could be a time value (including np.datetime64) or
a callable function such as np.nanmin. The last is useful when running with dt < 0.
repeatdt :
Optional interval (in seconds) on which to repeat the release of the ParticleSet (Default value = None)
lonlatdepth_dtype :
Floating precision for lon, lat, depth particle coordinates.
It is either np.float32 or np.float64. Default is np.float32 if fieldset.U.interp_method is 'linear'
Expand Down Expand Up @@ -482,8 +462,6 @@ def execute(self, pyfunc=AdvectionRK4, pyfunc_inter=None, endtime=None, runtime=

# Derive _starttime and endtime from arguments or fieldset defaults
_starttime = min_rt if dt >= 0 else max_rt
if self.repeatdt is not None and self.repeat_starttime is None:
self.repeat_starttime = _starttime
if runtime is not None:
endtime = _starttime + runtime * np.sign(dt)
elif endtime is None:
Expand Down Expand Up @@ -511,14 +489,8 @@ def execute(self, pyfunc=AdvectionRK4, pyfunc_inter=None, endtime=None, runtime=
moviedt = np.infty
if callbackdt is None:
interupt_dts = [np.infty, moviedt, outputdt]
if self.repeatdt is not None:
interupt_dts.append(self.repeatdt)
callbackdt = np.min(np.array(interupt_dts))
time = _starttime
if self.repeatdt:
next_prelease = self.repeat_starttime + (abs(time - self.repeat_starttime) // self.repeatdt + 1) * self.repeatdt * np.sign(dt)
else:
next_prelease = np.infty if dt > 0 else - np.infty
next_output = time + outputdt if dt > 0 else time - outputdt
next_movie = time + moviedt if dt > 0 else time - moviedt
next_callback = time + callbackdt if dt > 0 else time - callbackdt
Expand All @@ -539,9 +511,9 @@ def execute(self, pyfunc=AdvectionRK4, pyfunc_inter=None, endtime=None, runtime=
verbose_progress = True

if dt > 0:
next_time = min(next_prelease, next_input, next_output, next_movie, next_callback, endtime)
next_time = min(next_input, next_output, next_movie, next_callback, endtime)
else:
next_time = max(next_prelease, next_input, next_output, next_movie, next_callback, endtime)
next_time = max(next_input, next_output, next_movie, next_callback, endtime)

# If we don't perform interaction, only execute the normal kernel efficiently.
if self.interaction_kernel is None:
Expand All @@ -567,17 +539,6 @@ def execute(self, pyfunc=AdvectionRK4, pyfunc_inter=None, endtime=None, runtime=
break
# End of interaction specific code
time = next_time
if abs(time-next_prelease) < tol:
pset_new = self.__class__(
fieldset=self.fieldset, time=time, lon=self.repeatlon,
lat=self.repeatlat, depth=self.repeatdepth,
pclass=self.repeatpclass,
lonlatdepth_dtype=self.collection.lonlatdepth_dtype,
partitions=False, pid_orig=self.repeatpid, **self.repeatkwargs)
for p in pset_new:
p.dt = dt
self.add(pset_new)
next_prelease += self.repeatdt * np.sign(dt)
if abs(time - next_output) < tol or dt == 0:
for fld in self.fieldset.get_fields():
if hasattr(fld, 'to_write') and fld.to_write:
Expand Down
50 changes: 7 additions & 43 deletions parcels/particleset/particlesetaos.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,8 +1,6 @@
import sys
from copy import copy
from ctypes import c_void_p
from datetime import date, datetime
from datetime import timedelta as delta

import numpy as np
import xarray as xr
Expand Down Expand Up @@ -57,8 +55,6 @@ class ParticleSetAOS(BaseParticleSet):
Optional list of initial depth values for particles. Default is 0m
time :
Optional list of initial time values for particles. Default is fieldset.U.grid.time[0]
repeatdt : datetime.timedelta or float, optional
Optional interval on which to repeat the release of the ParticleSet. Either timedelta object, or float in seconds.
lonlatdepth_dtype :
Floating precision for lon, lat, depth particle coordinates.
It is either np.float32 or np.float64. Default is np.float32 if fieldset.U.interp_method is 'linear'
Expand All @@ -72,7 +68,7 @@ class ParticleSetAOS(BaseParticleSet):
Other Variables can be initialised using further arguments (e.g. v=... for a Variable named 'v')
"""

def __init__(self, fieldset=None, pclass=JITParticle, lon=None, lat=None, depth=None, time=None, repeatdt=None, lonlatdepth_dtype=None, pid_orig=None, **kwargs):
def __init__(self, fieldset=None, pclass=JITParticle, lon=None, lat=None, depth=None, time=None, lonlatdepth_dtype=None, pid_orig=None, **kwargs):
super().__init__()

# ==== first: create a new subclass of the pclass that includes the required variables ==== #
Expand Down Expand Up @@ -199,6 +195,10 @@ def ObjectClass_sizeof_forward(self):
# ==== dynamic re-classing completed ==== #
_pclass = object_class

if 'repeatdt' in kwargs.keys():
raise NotImplementedError('The repeatdt option has been deprecated after v2.4.2. '
'See https://docs.oceanparcels.org/en/latest/examples/tutorial_delaystart.html for how to release a particleset over multiple times')

self.fieldset = fieldset
if self.fieldset is None:
logger.warning_once("No FieldSet provided in ParticleSet generation. "
Expand Down Expand Up @@ -246,38 +246,9 @@ def ObjectClass_sizeof_forward(self):
assert lon.size == kwargs[kwvar].size, (
f"{kwvar} and positions (lon, lat, depth) don't have the same lengths.")

self.repeatdt = repeatdt.total_seconds() if isinstance(repeatdt, delta) else repeatdt
if self.repeatdt:
if self.repeatdt <= 0:
raise 'Repeatdt should be > 0'
if time[0] and not np.allclose(time, time[0]):
raise 'All Particle.time should be the same when repeatdt is not None'
self.repeatpclass = pclass
self.repeatkwargs = kwargs

ngrids = fieldset.gridset.size if fieldset is not None else 0
self._collection = ParticleCollectionAOS(_pclass, lon=lon, lat=lat, depth=depth, time=time, lonlatdepth_dtype=lonlatdepth_dtype, pid_orig=pid_orig, partitions=partitions, ngrid=ngrids, **kwargs)

if self.repeatdt:
if len(time) > 0 and time[0] is None:
self.repeat_starttime = time[0]
else:
collect_time = self._collection.time
if collect_time[0] and not np.allclose(collect_time, collect_time[0]):
raise ValueError('All Particle.time should be the same when repeatdt is not None')
self.repeat_starttime = copy(collect_time[0])
self.repeatlon = copy(self._collection.lon)
self.repeatlat = copy(self._collection.lat)
self.repeatdepth = copy(self._collection.depth)
for kwvar in kwargs:
self.repeatkwargs[kwvar] = copy(getattr(self._collection, kwvar))

if self.repeatdt:
if MPI and self._collection.pu_indicators is not None:
mpi_comm = MPI.COMM_WORLD
mpi_rank = mpi_comm.Get_rank()
self.repeatpid = pid_orig[self._collection.pu_indicators == mpi_rank]

self.kernel = None

def __del__(self):
Expand Down Expand Up @@ -546,7 +517,7 @@ def monte_carlo_sample(cls, start_field, size, mode='monte_carlo'):
raise NotImplementedError(f'Mode {mode} not implemented. Please use "monte carlo" algorithm instead.')

@classmethod
def from_particlefile(cls, fieldset, pclass, filename, restart=True, restarttime=None, repeatdt=None, lonlatdepth_dtype=None, **kwargs):
def from_particlefile(cls, fieldset, pclass, filename, restart=True, restarttime=None, lonlatdepth_dtype=None, **kwargs):
"""Initialise the ParticleSet from a zarr ParticleFile.
This creates a new ParticleSet based on locations of all particles written
in a zarr ParticleFile at a certain time. Particle IDs are preserved if restart=True
Expand All @@ -566,20 +537,13 @@ def from_particlefile(cls, fieldset, pclass, filename, restart=True, restarttime
time at which the Particles will be restarted. Default is the last time written.
Alternatively, restarttime could be a time value (including np.datetime64) or
a callable function such as np.nanmin. The last is useful when running with dt < 0.
repeatdt : datetime.timedelta or float, optional
Optional interval on which to repeat the release of the ParticleSet. Either timedelta object, or float in seconds.
lonlatdepth_dtype :
Floating precision for lon, lat, depth particle coordinates.
It is either np.float32 or np.float64. Default is np.float32 if fieldset.U.interp_method is 'linear'
and np.float64 if the interpolation method is 'cgrid_velocity'
**kwargs :
Keyword arguments passed to the particleset constructor.
"""
if repeatdt is not None:
logger.warning(f'Note that the `repeatdt` argument is not retained from {filename}, and that '
'setting a new repeatdt will start particles from the _new_ particle '
'locations.')

pfile = xr.open_zarr(str(filename))
pfile_vars = [v for v in pfile.data_vars]

Expand Down Expand Up @@ -621,7 +585,7 @@ def from_particlefile(cls, fieldset, pclass, filename, restart=True, restarttime

return cls(fieldset=fieldset, pclass=pclass, lon=vars['lon'], lat=vars['lat'],
depth=vars['depth'], time=vars['time'], pid_orig=vars['id'],
lonlatdepth_dtype=lonlatdepth_dtype, repeatdt=repeatdt, **kwargs)
lonlatdepth_dtype=lonlatdepth_dtype, **kwargs)

def __iadd__(self, particles):
"""Add particles to the ParticleSet.
Expand Down
Loading