Float precision results in invalid time stamps

[finsberg]

Consider the following example where we want to save some function after a time loop and then later we would like to read in the same function (this could also be in a different postprocessing script)

from pathlib import Path
from mpi4py import MPI
import adios4dolfinx
import dolfinx


N = 1
comm = MPI.COMM_WORLD
domain = dolfinx.mesh.create_unit_square(comm, N, N)
V = dolfinx.fem.functionspace(domain, ("CG", 1))
u = dolfinx.fem.Function(V, name="u")
t = 0.0
dt = 0.01  # Time step size
T_final = 1.0  # Final time

checkpoint_fname = Path("checkpoint.bp")
adios4dolfinx.write_mesh(checkpoint_fname, domain)
while t < T_final:
    t += dt
print("Save at time:", t)
print(abs(t - T_final))

adios4dolfinx.write_function(checkpoint_fname, u, time=t, name="u")

new_mesh = adios4dolfinx.read_mesh(filename=checkpoint_fname, comm=MPI.COMM_WORLD)
new_V = dolfinx.fem.functionspace(new_mesh, ("CG", 1))
new_u = dolfinx.fem.Function(new_V, name="u")
adios4dolfinx.read_function(checkpoint_fname, new_u, name="u", time=1.0)
print("u at time 1.0:", new_u.x.array)

The output of the program is as follows

$ python mwe.py
Save at time: 1.0000000000000007
6.661338147750939e-16
Traceback (most recent call last):
  File "/home/shared/examples-fenicsx/mms/mwe.py", line 28, in <module>
    adios4dolfinx.read_function(checkpoint_fname, new_u, name="u", time=1.0)
  File "/dolfinx-env/lib/python3.12/site-packages/adios4dolfinx/checkpointing.py", line 446, in read_function
    input_array, starting_pos = read_array(
                                ^^^^^^^^^^^
  File "/dolfinx-env/lib/python3.12/site-packages/adios4dolfinx/adios2_helpers.py", line 308, in read_array
    raise KeyError(
KeyError: 'No data associated with u_time=1.0 found in checkpoint.bp'

Now, one way to mitigate this is to use time_stamps stored in the file i.e

import numpy as np

time_stamps = adios4dolfinx.read_timestamps(
    checkpoint_fname, comm=comm, function_name="u"
)
time = 1.0
index = np.argmin(np.abs(time_stamps - time))
assert np.isclose(time_stamps[index], time)
adios4dolfinx.read_function(checkpoint_fname, new_u, name="u", time=time_stamps[index])
print("u at time 1.0:", new_u.x.array)

which works. However, I wonder whether we should do this within adios4dolfinx instead, since this would be a very common usecase which might trip off many users.

There are two other options to solve this which I think are also reasonable, so one of the following fixes should be implemented I think (ordered by preference)

1. Handle float precision within adios4dolfinx and only raise KeyError if a time stamp is not within a given tolerance of the available time stamps (this tolerance should also be possible to specify). In this case we need to update the API for read_function to also take an optional tolerance which could have a default value close to the float precision.
2. Print a more informative error message saying that the KeyError might be due to a float precision error and link to documentation of how to solve this (i.e the fix above). Then this needs to also be documented
3. Only allow integer time stamps to avoid floating number precision errors.