Out of boundaries error #1086
Replies: 3 comments 14 replies
-
|
Dear @mainarabgouveia. You will have to include a recovery kernel to your execution: pset.execute(..., recovery={ErrorCode.ErrorOutOfBounds: DeleteParticle})And then create a Recovery Kernel to delete the particle: def DeleteParticle(particle, fieldset, time):
particle.delete()See also https://github.com/OceanParcels/parcels/blob/master/parcels/examples/tutorial_Agulhasparticles.ipynb for a full example |
Beta Was this translation helpful? Give feedback.
-
|
Hi Erik, I have included the DeleteErrorParticle kernel you mentioned above and included it in my set execute However, the model is still terminating Do you have any insight on why this may be happening? Thanks in advance |
Beta Was this translation helpful? Give feedback.
-
|
It's really important that the kernels = [AdvectionRK4, DiffusionUniformKh, SampleAge, Sample_land, Distance, DeleteErrorParticle]See also the last paragraph of https://docs.oceanparcels.org/en/latest/examples/tutorial_kernelloop.html#Working-with-Status-Codes |
Beta Was this translation helpful? Give feedback.
-
|
Hi Erik, I am getting an out-of-bound error while using the AdvectionAnalytical advection kernel. Even though I did use the above kernel that deletes the particles at the boundaries, the code still gives me the out-of-bound error for U velocity, while the particle state is still not over 40. This is the code: This is the error: The (x, y) at (-79.0750963318386, 31.960647583007816) are the upper bounds of my domain. Do you have any ideas as to why I am getting this error? This is the first time I am using analytical advection to advect particles and I am not sure if I am missing something. Thank you in advance! |
Beta Was this translation helpful? Give feedback.
-
|
I think in this case, the error is in the AnalyticalAdvection scheme. This was coded up as a quick showcase that it would in principle be possible to implement such a scheme within parcels, but wasn't well-tested. I.e. no proper error-handling. You've run into the limits of its capabilities, I fear... However, my next big project the coming weeks/months is to 'properly' implement AnalyticalAdvection into Parcels, also in JIT-mode. I will then make sure that these out-of-bounds errors are properly treated, too. So I suggest you either wait a few more weeks and/or use RungeKutta for the timing being; or write your own removal kernel to remove particles before they even get close to the boundary (if that is even an option for your use-case) |
Beta Was this translation helpful? Give feedback.
-
|
Thank you for letting me know! I might be able to implement a removal kernel after the parcels exit the Straits of Florida -my region is the Gulf of Mexico, so something like that could work. In your future project, is the fully time-dependent scheme going to be implemented or is it going to be a further development of the stepwise scheme? Thank you! |
Beta Was this translation helpful? Give feedback.
-
|
Tks for your help and time. |
Beta Was this translation helpful? Give feedback.
-
|
Hi, I have a similar problem. I am getting def CheckOutOfBounds(particle, fieldset, time):
if particle.state == StatusCode.ErrorOutOfBounds:
particle.delete()
kernels = [parcels.AdvectionRK4,CheckOutOfBounds]
pset.execute(kernels,
runtime=timedelta(hours=72),
dt=timedelta(minutes=5),
output_file=output_file,
)I am still getting another error message: If anyone can help? Thank you. This is my code: import math
from datetime import timedelta as delta
from datetime import datetime,timedelta
from operator import attrgetter
from glob import glob
import matplotlib.pyplot as plt
import numpy as np
import trajan as ta
import xarray as xr
from IPython.display import HTML
from matplotlib.animation import FuncAnimation
import parcels
from parcels import StatusCode
data_folder = '.'
hydrofiles = sorted(glob(f"{data_folder}/Strait_2004*.nc"))
def DeleteParticle(particle, fieldset, time):
particle.delete()
def CheckOutOfBounds(particle, fieldset, time):
if particle.state == StatusCode.ErrorOutOfBounds:
particle.delete()
filenames = {
"U": hydrofiles,
"V": hydrofiles,
"W": hydrofiles,
}
variables = {"U": "UVEL", "V": "VVEL", "W": "WVEL"}
dimensions = {
"U": {"lon": "XG", "lat": "YG", "depth": "Z", "time": "time"},
"V": {"lon": "XG", "lat": "YG", "depth": "Z", "time": "time"},
"W": {"lon": "XG", "lat": "YG", "depth": "Z", "time": "time"},}
fieldset = parcels.FieldSet.from_mitgcm(
filenames, variables, dimensions, mesh="spherical",
allow_time_extrapolation=True,)
advection_duration = 48
ReleaseTime=0
day_start = datetime(2004,1,1,1,0)
day_end = day_start+delta(hours=advection_duration)
nptc=10
npart=8
lonp = [115.00,115.25,115.5,115.75,116,116.25,116.5,116.75]
latp = -8.0*np.ones(npart)
lonp=np.array(lonp)
latp=np.array(latp)
depthp = -1.0*np.ones(npart)
releaseTimes = np.arange(0, advection_duration) * delta(hours=1)
releaseTimes += day_start
particles_t = np.array([])
particles_lat = np.array([])
particles_lon = np.array([])
particles_depth = np.array([])
for i1,time in enumerate(releaseTimes):
for i2 in range(nptc):
particles_t = np.append(particles_t,np.array([(time-day_start).total_seconds() for i in range(len(lonp))]))
particles_lat = np.append(particles_lat,latp)
particles_lon = np.append(particles_lon,lonp)
particles_depth=np.append(particles_depth,depthp)
time_origin = datetime.utcfromtimestamp(np.datetime64(str(fieldset.time_origin)).astype(int)*1e-9)
time_start = day_start - time_origin
particles_t += time_start.total_seconds()
pset = parcels.ParticleSet(fieldset = fieldset,
pclass = parcels.JITParticle,
lon = particles_lon,
lat = particles_lat,
depth=particles_depth,
time = particles_t,)
output_file = pset.ParticleFile(name="test.zarr", outputdt=timedelta(minutes=30))
kernels = [parcels.AdvectionRK4,CheckOutOfBounds]
pset.execute(kernels,
runtime=timedelta(hours=72),
dt=timedelta(minutes=5),
output_file=output_file,
) |
Beta Was this translation helpful? Give feedback.
-
|
@kartadikaria what sort of environment are you using? I'm interested in the following:
I notice that you're running on Python 3.8, which means that you likely have an outdated version of Parcels. Newer versions of Parcels (and a lot of the scientific Python ecosystem) have dropped support for Python 3.8 - I recommend that you update to Python 3.10 at the minimum. |
Beta Was this translation helpful? Give feedback.
-
|
Hi @VeckoTheGecko, Thanks in advance. I experienced the above problem on 2 PCs.
Previously, I used python 3.10 on both machines. However, since I found a similar problem on both of them, then I downgraded the python on one machine. Yet, the problem is still appearing. Here the other libraries' version:
This is the fieldset: fields: name : 'U' grid : RectilinearZGrid(lon=array([ 114.10, 114.11, 114.11, ..., 117.08, 117.09, 117.09], dtype=float32), lat=array([-9.50, -9.49, -9.49, ..., -7.54, -7.54, -7.53], dtype=float32), time=array([ 0.00, 3600.00, 7200.00, ..., 248400.00, 252000.00, 255600.00]), time_origin=2004-01-01T00:30:00.000000000, mesh='spherical') extrapolate time: True time_periodic : False gridindexingtype: 'mitgcm' to_write : False name : 'V' grid : RectilinearZGrid(lon=array([ 114.10, 114.11, 114.11, ..., 117.08, 117.09, 117.09], dtype=float32), lat=array([-9.50, -9.49, -9.49, ..., -7.54, -7.54, -7.53], dtype=float32), time=array([ 0.00, 3600.00, 7200.00, ..., 248400.00, 252000.00, 255600.00]), time_origin=2004-01-01T00:30:00.000000000, mesh='spherical') extrapolate time: True time_periodic : False gridindexingtype: 'mitgcm' to_write : False name : 'W' grid : RectilinearZGrid(lon=array([ 114.10, 114.11, 114.11, ..., 117.08, 117.09, 117.09], dtype=float32), lat=array([-9.50, -9.49, -9.49, ..., -7.54, -7.54, -7.53], dtype=float32), time=array([ 0.00, 3600.00, 7200.00, ..., 248400.00, 252000.00, 255600.00]), time_origin=2004-01-01T00:30:00.000000000, mesh='spherical') extrapolate time: True time_periodic : False gridindexingtype: 'mitgcm' to_write : False name: 'UV' U: V: W: None name: 'UVW' U: V: W: I hope, I can get the solution. :) Thanks. |
Beta Was this translation helpful? Give feedback.
-
|
Hi @kartadikaria ,
Can you elaborate more on what this means? You need to run the script as a file, and not line by line like this. This is a requirement for the Can you try run the Python script? (i.e., Failing that: it would be good to get this to a minimal example so that we can debug on our end. Can you do the following, seeing that you still get the same errors?:
|
Beta Was this translation helpful? Give feedback.
-
|
Hi @VeckoTheGecko, I tried to run block-by-block to understand the process. For instance after pset = parcels.ParticleSet() --> output_file ---> kernels ---> pset.execute. The problem occured during pset. execute. Anyway, after I ran it using "python [my_parcels.py]", as you mentioned.....It had no problem and could finish running. I am still confused, but I accept. Thanks and now I can continue to create other fieldsets in parcels V3. |
Beta Was this translation helpful? Give feedback.
-
|
I'm not entirely sure which editor you're using, or how you're running block by block. We recommend that you develop in Jupyter Notebook/Jupyter Lab or by running the script normally via |
Beta Was this translation helpful? Give feedback.
Uh oh!
There was an error while loading. Please reload this page.
Uh oh!
There was an error while loading. Please reload this page.
-
Hi,
How do I make my released particles that leave the domain or arrive on the land don't cause my simulations to stop?
I'm getting this error:
INFO: Compiled ArrayJITParticleAdvectionRK45AdvectionDiffusionM1 ==> /tmp/parcels-1000/lib77962480406540cb0477e74d24ef8364_0.so
98% (14938200.0 of 15120000.0) |####### | Elapsed Time: 0:15:43 ETA: 0:00:11
OutOfBoundsError Traceback (most recent call last)
/tmp/ipykernel_32755/3033041448.py in
1 # Now run the remaining ParticleSet (using less accurate Euler Forward for speedup)
2
----> 3 pset.execute(pset.Kernel(AdvectionRK45) + AdvectionDiffusionM1, # the kernel (which defines how particles move)
4 runtime=timedelta(days=175), # the total length of the run
5 dt=-timedelta(minutes=5), # the timestep of the kernel (negative timestep for backwards run)
~/anaconda3/envs/py3_parcels/lib/python3.9/site-packages/parcels/particleset/baseparticleset.py in execute(self, pyfunc, pyfunc_inter, endtime, runtime, dt, moviedt, recovery, output_file, movie_background_field, verbose_progress, postIterationCallbacks, callbackdt)
456 # If we don't perform interaction, only execute the normal kernel efficiently.
457 if self.interaction_kernel is None:
--> 458 self.kernel.execute(self, endtime=next_time, dt=dt, recovery=recovery, output_file=output_file,
459 execute_once=execute_once)
460 # Interaction: interleave the interaction and non-interaction kernel for each time step.
~/anaconda3/envs/py3_parcels/lib/python3.9/site-packages/parcels/kernel/kernelsoa.py in execute(self, pset, endtime, dt, recovery, output_file, execute_once)
225 recovery_kernel = recovery_map[p.state]
226 p.set_state(StateCode.Success)
--> 227 recovery_kernel(p, self.fieldset, p.time)
228 if p.isComputed():
229 p.reset_state()
~/anaconda3/envs/py3_parcels/lib/python3.9/site-packages/parcels/tools/statuscodes.py in recovery_kernel_out_of_bounds(particle, fieldset, time)
194 # TODO: JIT does not yet provide the context that created
195 # the exception. We need to pass that info back from C.
--> 196 raise OutOfBoundsError(particle, fieldset)
197 else:
198 error = particle.exception
OutOfBoundsError: 0
Particle P[4](lon=17.006130, lat=-29.650196, depth=0.000000, time=36208800.000000)
Time: 2021-02-23T02:30:00.000000000, timestep dt: -3600.000000
Out-of-bounds sampling by particle at (17.006130, -29.650196, 0.000000)
Beta Was this translation helpful? Give feedback.
All reactions