Issue with time values written to XDMF file and reloading H5 as TimeSeries: "Sample points are not strictly monotone in series"


125
views
0
4 months ago by
I'm attempting to do something similar as in the reaction-advection-diffusion example, where I need to load a h5 file as a TimeSeries() using

timeseries_w = TimeSeries('timeseries')​


This gives me the following error

*** -------------------------------------------------------------------------
*** Error:   Unable to read time series from file.
*** Reason:  Sample points for mesh data are not strictly monotone in series "./data/pressure_ref/p1".
*** Where:   This error was encountered inside TimeSeries.cpp.
*** Process: 0
*** 
*** DOLFIN version: 2017.2.0
*** Git changeset:  4c59bbdb45b95db2f07f4e3fd8985c098615527f
*** -------------------------------------------------------------------------


In the corresponding XDMF file time points were written as 

<Time Value="0" />
...
<Time Value="0.01" />
...
<Time Value="0.02" />
...
<Time Value="0.029999999999999999" />
...
<Time Value="0.040000000000000001" />
...


and so forth even though my time steps are always 0.1

t = 0.0
dt = 0.1

for n in range(nt):
    t += dt
    ...
    file.write(u, t)


Why does this occur and how can I avoid it? I tried using

t = round(t, 2)


but that didn't do the trick.

MWE:

from dolfin import *

# Create mesh and define function space
mesh = UnitSquareMesh(32, 32)
V = FunctionSpace(mesh, "Lagrange", 1)

# Define Dirichlet boundary (x = 0 or x = 1)
def boundary(x):
    return x[0] < DOLFIN_EPS or x[0] > 1.0 - DOLFIN_EPS

t = 0.0
dt = 0.1
nt = 20

# Define boundary condition
u0 = Constant(0.0)
bc = DirichletBC(V, u0, boundary)

# Define variational problem
u = TrialFunction(V)
v = TestFunction(V)
f = Expression("10*exp(-(pow(x[0] - 0.5, 2) + pow(x[1] - 0.5, 2)) / 0.02)", degree=2)
g = Expression("sin(5*x[0])", degree=2)
a = inner(grad(u), grad(v))*dx
L = f*v*dx + g*v*ds

file = XDMFFile(mesh.mpi_comm(), 'u.xdmf')
file.write(mesh)

for n in range(nt):
	# Compute solution
	u = Function(V)
	solve(a == L, u, bc)

	file.write(u, t)
	t += dt

file.close()
Community: FEniCS Project

3 Answers


3
4 months ago by
If using TimeSeries is not fundamental, you could read from file in this way:
u = Function(V)
fileh5 = HDF5File(MPI.comm_world, "u.h5", "r")
for i in range(10):
    attribute_name = "VisualisationVector/%d" % i
    fileh5.read(u, attribute_name)
​
1
Thanks Eleonora, your suggestion works!
written 4 months ago by alexdiem  
0
4 months ago by
I don't understand why this would cause an error regarding monotonicity, but I'll answer your question represented in the MWE.

When the time step size is constant, as in your case, you should instead set

    t = dt*n

so that you do not accumulate floating point errors. Note that you should set t in this case before writing in the loop, and remove your previous initialization of t = 0.
That unfortunately didn't make a difference. If it did I believe my attempt with

t = round(t, 2)​


should have been successful. But if that shouldn't normally cause a monotonicity issue, do you have any ideas what else could be the problem? I get the following messages from HDF5-DIAG just before the Fenics error message (sorry for not adding them earlier, I just thought I'd identified the source of the problem)

EDIT: Message removed, too long...

written 4 months ago by alexdiem  
Would you please provide a MWE that reproduces the error?
written 4 months ago by Alexander G. Zimmerman  
I have provided one above and after running that I just need to run 

from fenics import *
tw1 = TimeSeries('u')​


and I get the error

written 4 months ago by alexdiem  
0
4 months ago by
So why there is this long floating number? It is the way how decimal number is converted into string in XDMFFile.cpp, see https://bitbucket.org/fenics-project/dolfin/src/5c18533c286f81a177e6f2a3a27ca4ca57e07f8e/dolfin/io/XDMFFile.cpp?at=master&fileviewer=file-view-default#XDMFFile.cpp-450
Boost::lexical_cast just produces this value, if you wonder, play with the following:
import dolfin as d

base_code = '''
#include <pybind11/pybind11.h>
#include <boost/lexical_cast.hpp>
namespace py = pybind11;

#include <iostream>

void to_double(double number)
{
    std::cout << boost::lexical_cast<std::string>(number) << std::endl;
};

PYBIND11_MODULE(SIGNATURE, m)
{
    m.def("to_double", &to_double);
}
'''

code = d.compile_cpp_code(base_code)
code.to_double(0.03)​

The effect is combination of both IEEE-754 Floating Point representation and boost casting.

Why TimeSeries is not working? I do not know, but I'd recommend not to use it. It is an old piece of code, no longer maintained and will be dropped. Use HDF5File or XDMFFile's checkpointing functionality, see e.g. https://www.allanswered.com/post/ovxg/
Thanks! I was using TimeSeries because it was used in the reaction-advection-diffusion chapter in the tutorial. So it appeared to me that that was the recommended way of loading data from a previous solution.
written 4 months ago by alexdiem  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »