Runge-Kutta 4/5 scheme requirement for extra Variable ‘next_dt’`

[selipot]

I would like to test the Runge-Kutta 4 advection scheme against the Runge-Kutta 4/5 scheme. The docstring for the kernel AdvectionRK45 states that:

Note that this kernel requires a Particle Class that has an extra Variable ‘next_dt’

Is there more information in the documentation about what that extra variable corresponds to? I am unable to find it. Any help or pointer would be greatly appreciated.

[erikvansebille]

HI @selipot, you can see an example of how that's done (using parcels v3) below:

RK45Particle = JITParticle.add_variable('next_dt', dtype=np.float32, initial=dt) 

And then use RK45Particle instead of JITParticle in your pset creation

For older Parcels versions, you don't need to set this next_dt explicitly

Note that @daanreijnders found an issue with the RK45 implementation at #1516, that's on my list to explore/fix in the coming days. So keep an eye on that Issue

[selipot]

Hi @erikvansebille thank you for your reply. I am currently using version 2.4.1.dev9 so if I understand properly what you wrote I should simply define my kernel as kernel = pset.Kernel(AdvectionRK45) or is it something else? The same code works with kernel = pset.Kernel(AdvectionRK4) but fails with AdvectionRK45. I am getting a zarr.errors.BoundsCheckError: index out of bounds for dimension with length 2

(I'd like to take note of the issue that @daanreijnders posted at #1516 but the link to the notebook there seems outdated.)

[daanreijnders]

Hmm, that is interesting. Did you play around with the tolerance? The culprit was that Parcels didn't automatically convert the tolerance from degrees to meters, so in spherical setups, it allowed for very large errors. Erik's fixed this in #1539.

[selipot]

I did not change the tolerance at all. Maybe the fact that my velocity field was at approximately 1/25 degree resolution helped?

[daanreijnders]

I think the resolution wouldn't have influenced this. The bug was in the tolerance being read in degrees rather than meters, which was determined by whether the data was provided in spherical or flat coordinates (the latter being the case for the unit tests, which is why the issue wasn't spotted before).

I'm a bit puzzled as to why you didn't experience this issue at all. @erikvansebille you verified the results of my run too, right?
It might be worthwhile double checking your run with RK45 in the latest version of parcels, for just 10 days or so...

[erikvansebille]

@erikvansebille you verified the results of my run too, right?

No, I didn't rerun your simulations from #1516. I wrote some unit tests to verify if the new code worked as expected, but haven't tested it in the wild.
The link to the zip you prepared has expired, can you try to rerun it, @daanreijnders?

[daanreijnders]

I reran my test with the new changes and it was fine, though much slower than RK4.

What I meant to ask was whether you reproduced the error I had before the fix was implemented. I just find it strange that @selipot didn't run into it, so I'm now doubting whether something else was wrong with my setup.