Simulation gets stuck on first timestep since I installed Parcels v3.0

[martoconnh]

Hi !

So I have just installed the latest parcels version, v3.0. Before that, I had been getting this warning with Parcels v2.4.2 :

<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast

I know this has been fixed in the v3.0 (already asked in a discussion) , but now I am having a new issue. Even though I got that warning in v2.4.2, my simulation worked, but with this new version, my simulation gets stuck on the first time step. No errors nor warnings, the only thing that I am getting on screen is the following:

INFO: Output files are stored in output_arousa1.zarr.
  0%|          | 900.0/432000.0 [00:30<2:35:45, 46.13it/s]

I don't know if this is a common issue, maybe I should change something in my code which is different from v2.4.2 to v3.0 ?

[michaeldenes]

Hi @martoconnh,

It's quite difficult to diagnose a problem without all the information. Would you be able to share a minimal example of code that reproduces this issue and we can try and help from there?

Cheers,
Michael

[martoconnh]

Sure!

Here it is my code. It simulates only 2 particles through 5 days. I still get the same problem even with this one.

# -*- coding: utf-8 -*-

""" Brief example """

import os
from datetime import timedelta as delta
from parcels import (AdvectionRK4,FieldSet,JITParticle,ParticleSet,StatusCode)

# Parameters
total_days = 5  # Max of 23 days (since 7 to 30 september)

#%% FieldSet creation

# List with our .nc4 files to create the FieldSet
dataset_folder = r'C:\Users\Marto\Documents\UNIVERSIDADE\TFG\FICHEIROS_AROUSA\hidrodinamicos_desde_20170907'
arquivos_nc4 = [os.path.join(dataset_folder, archivo) for archivo in os.listdir(dataset_folder) if archivo.endswith('.nc4')]

# Variables and dimensions
variables = {'U': 'u', 'V': 'v'}
dimensions = {'time': 'time', 'lon': 'lon', 'lat': 'lat'}

# FieldSet Creation
fieldset = FieldSet.from_netcdf(arquivos_nc4, variables, dimensions)

#%%  Coordinates and ParticleSet definition

class Particula(JITParticle):
    pass

# Only 2 particles from 2 separate points
lats = [42.635,42.5]
lons = [-8.78,-8.83]

pset = ParticleSet(fieldset=fieldset, pclass=Particula, lon=lons, lat=lats)

#%%  Execution

total_runtime = delta(days=total_days)
dt = delta(minutes=15)

# Delete particle if out of bounds
def DeleteParticle(particle, fieldset, time):
    if particle.state == StatusCode.ErrorOutOfBounds:
        particle.delete()

kernels = [AdvectionRK4,DeleteParticle]

pset.execute(kernels,
             runtime=total_runtime,
             dt=dt,
             output_file=pset.ParticleFile(name='output_example.zarr'))

Thanks!

Martinho

[erikvansebille]

Hmm, it's strange that this code hangs. And what if you make these two changes to further simplify the code:

1. Use pset = ParticleSet(fieldset=fieldset, pclass=JITParticle, lon=lons, lat=lats) (so use a simple JITParticle)
2. Use kernels = [AdvectionRK4] (so don't delete the particles)

Does either of these make a difference?

[martoconnh]

Hi !

I have just tried, but it does not make any difference. Maybe there is something else that changed from v2.4.2 to v3.0.0 that I'm not using propperly??

[erikvansebille]

OK, good to know it's not either of these two options at #1461 (comment)

Can you then try to add, immediately after you created the fieldset

fieldset.computeTimeChunk(0,1)

If your script then hangs on that command, it clearly is related to the loading in of the fieldset data from the NectCDF files. But that would be strange, because changes to the code for the fieldset creation were minimal going from 2.4.2 to 3.0.0. Still, good to check it before we dig deeper.

[martoconnh]

Hi!

Just tried that, but it does not hang on that command. The output is 3600 (don't know if that's useful somehow)

[erikvansebille]

OK, final test: what if you remove the output_file=pset.ParticleFile(name='output_example.zarr') line in the pset.execute()? Does that make a difference?

[martoconnh]

Oh, that did actually work. So maybe I should store the simulation data in some other way??

[martoconnh]

Oh okay, just realized I should add this to the output_file line:
output_file=pset.ParticleFile(name='output_example.zarr' , outputdt=timedelta(hours=1))

I have tried this and now everything is working just fine. Maybe this has changed from v2.4.2 to v3.0? I looked it up in https://docs.oceanparcels.org/en/latest/examples/parcels_tutorial.html but i don't recall this extra argument in the output_file was necessary in the former version (maybe I'm wrong).

Anyway, thank you so much for your help (and patience) !!

Cheers,

Martinho Rial

[erikvansebille]

Ah, good to know that adding outputdt fixed your problem. This was already necessary in v2.4.2 too, but perhaps the code would back then be bit more lenient if you had forgotten it. Glad it's solved!