Move to float64 throughout Parcels?

[erikvansebille]

Early on in the design of Parcels (pre-v2), we decided to do most computations in float32. With the development of C-grid interpolation in v2, we required float64 accuracy to properly deal with boundary conditions. Still, we left a lot of float32 support under the hood in Parcels

With the recent release of NumPy v2, some of these float32-implementations seem to be breaking some of the unit tests. This probably has to do with this change to data type promotion in NumPy. See also https://github.com/OceanParcels/parcels/pull/1603#issuecomment-2238974528.

This has been quick-fixed in https://github.com/OceanParcels/parcels/pull/1603/commits/20206c23e7f993d6f846c47ca1bdbf9b5f325d15, but ideally, we may want to reconsider the use of float32 in Parcels altogether.

This Issue is a place to discuss/explore what would happen if we move fully to float64 in Parcels

[VeckoTheGecko]

I think it would be good to go ahead with this migration.

[VeckoTheGecko]

I think there's a distinction here between Field input and output. I think:

• don't cast dtypes for the input field data
• write trajectory output in float64

Looking at fieldfilebuffer.py (ie.., as of fd1dda2e16962522a0fa6ff817bcb4843387fb30) we cast the field data to float32.

[JamiePringle]

Dear all -- for those of us who do very large particle tracking runs, it is very, very helpful to be able to coerce the output to float32. I routinely end up with 10s of terabytes of output, and do not look forward to dealing with 20s of terabytes of output. While internally it can be helpful to avoid rounding errors, for analysis of pathways it usually makes little difference.

To estimate the error in meters for a coordinate system that would cover the earth, take the circumference of the earth in meters, or 4e7 meters, by the precision of float32, finfo(float32).eps, and you find the error is about 5 meters. This number is the smallest number you can add to 4e7 meters and get a different location, and is good for most analysis.

Can we include some way to coerce the output to float32?

Jamie

[VeckoTheGecko]

Thanks for the input @JamiePringle - we'll keep the flexibility of float32 trajectory outputs in v4