This may not strictly be a fenics question, but suppose I run a time dependent problem, of the sort:
output = XDMFFile("out.xdmf")
for n in range(nmax):
# computational steps
output.write(u,t)
output.close()
If I then read this in with paraview, it’s all as expected. But I was wondering how to read the data back in with python. I saw the post here, but I don’t quite see how to adapt it to what I’m doing.
You should use write_checkpoint and read_checkpoint as described in the previous post to read things back into FEniCS.
Two follow up questions:
write?nx = ny = 8
mesh = UnitSquareMesh(nx, ny)
V = FunctionSpace(mesh, 'P', 1)
alpha = 3; beta = 1.2
f = Expression('1 + x[0]*x[0] + alpha*x[1]*x[1] + beta*t', degree=2, alpha=alpha, beta=beta, t=0)
f_out = XDMFFile("test.xdmf")
for j in range(5):
t = float(j/5)
f.t = t
f_out.write_checkpoint(project(f,V), "f",t)
f_out.close()
f1 = Function(V)
f_in = XDMFFile("test.xdmf")
f_in.read_checkpoint(f1,"f",0) # this works fine
f_in.read_checkpoint(f1,"f",1) # this throws an error
The error that I get is:
---------------------------------------------------------------------------
RuntimeError Traceback (most recent call last)
<ipython-input-48-ab851acfbb36> in <module>
----> 1 f_in.read_checkpoint(f1,"f",1)
RuntimeError:
*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at
***
*** fenics-support@googlegroups.com
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error: Unable to reading data from XDMF file.
*** Reason: Storage format "" is unknown.
*** Where: This error was encountered inside XDMFFile.cpp.
*** Process: 0
***
*** DOLFIN version: 2019.1.0
*** Git changeset:
*** -------------------------------------------------------------------------
write saves the data at the vertices, efficiently projecting the solution to a CG-1 function space. write_checkpoint saves the data with the appropriate function space, and can therefore be read into dolfin again, to the appropriate function space.
You are not saving multiple time steps to your file as you have not changed the append kwarg in write_checkpoint (see: https://fenicsproject.org/docs/dolfin/latest/cpp/de/d71/classdolfin_1_1XDMFFile.html#a9ef9031aec7bfec8f289d5210a1c811c). The following works:
from dolfin import *
nx = ny = 8
mesh = UnitSquareMesh(nx, ny)
V = FunctionSpace(mesh, 'P', 1)
alpha = 3; beta = 1.2
f = Expression('1 + x[0]*x[0] + alpha*x[1]*x[1] + beta*t', degree=2, alpha=alpha, beta=beta, t=0)
f_out = XDMFFile("test.xdmf")
f_out.write_checkpoint(project(f, V), "f", 0, XDMFFile.Encoding.HDF5, False) # Not appending to file
for j in range(1,5):
t = float(j/5)
f.t = t
f_out.write_checkpoint(project(f,V), "f",t, XDMFFile.Encoding.HDF5, True) #appending to file
f_out.close()
f1 = Function(V)
f_in = XDMFFile("test.xdmf")
f_in.read_checkpoint(f1,"f",0)
f_in.read_checkpoint(f1,"f",1)
Yes, that now works. One question: output in generated with write_checkpoint now doesn’t seem to be readable by Paraview. Is this to be expected, or is there some way to make it compatible with both?
then your paraview is out of data. You Need to use Paraview >=5.5.0
is there a way with the method read_checkpoint to read the actual time of the i_th time step ? let’s say for example you want to recreate the exact same xdmf file ?
from fenics import *
file_in = XDMFFile("in.xdmf")
file_out = XDMFFile("out.xdmf")
mesh = Mesh()
file_in.read(mesh)
V = FunctionSpace(mesh, "P", 1)
u = Function(V)
for i in range(0, 10):
file_in.read_checkpoint(u, "u", i)
t = i # ?????
file_out.write_checkpoint(u, "u", i)
export_file.write_checkpoint(u, "u", t)
As you can see from the documentation: https://fenicsproject.org/docs/dolfin/latest/cpp/de/d71/classdolfin_1_1XDMFFile.html#aa86a07b92188a2da32d1a0e128106c91
That is not possible, as read_checpoint only takes in a counter.
I’m using 5.7.0, and I’m seeing that some of the fields are listed as partial. The consequence seems to be that I don’t have the field at every time at which I performed a write_checkpoint. See attached screen capture.
Are you using the append functionality. I see none of the issues with the following code:
from dolfin import *
nx = ny = 8
mesh = UnitSquareMesh(nx, ny)
V = FunctionSpace(mesh, 'P', 1)
alpha = 3; beta = 1.2
f = Expression('1 + x[0]*x[0] + alpha*x[1]*x[1] + beta*t', degree=2, alpha=alpha, beta=beta, t=0)
f_out = XDMFFile("test.xdmf")
f_out.write_checkpoint(project(f,V), "f", 0.0, XDMFFile.Encoding.HDF5, False)
for j in range(1,5):
t = float(j/5)
f.t = t
f_out.write_checkpoint(project(f,V), "f", t, XDMFFile.Encoding.HDF5, True)
f_out.close()
f1 = Function(V)
f_in = XDMFFile("test.xdmf")
f_in.read_checkpoint(f1,"f",0)
f_in.read_checkpoint(f1,"f",1)
f_in.close()
Yeah, my mistake is I was appending to an existing file. Once I put False in the first call, it cleared up.
Hello Jorgen, is there still any way to read the data back in when they were already written using write instead of write_checkpoint?
No, you cannot read in data written with XDMFFile.write.
Hi there,
I have a problem relating to this topic. My aim is to get a time-averaged velocity field for the vortex shedding behind a cylinder. I have run a 2D Navier-Stokes simulation in FEniCS until a quasi-steady shedding regime has been established. Once established, I saved the data into an XDMF file using write_checkpoint for 30 shedding cycles. I’ve then used Matlab to get an average of the ‘vector’ dataset in the hdf5 files over the 30 shedding cycles. I saved this new average ‘vector’ dataset into a separate xdmf file using the same domain data (dofs, mesh geom, mesh topology etc.) that was in the original xdmf file created by initial simulation. All of the groups, datasets and function have the same names as those prescribed in the ‘write_checkpoint’ function, however, I am having trouble loading this time-averaged velocity field back into FEniCS as a function using read_checkpoint. Any help with this topic is much appreciated, is it even possible??
Thanks,
Ben
You would need to make a minimal example, reproducing your error for anyone to be of any guidance. What error do you get once you use read_checkpoint?
Hi Jorgen,
Thanks for the swift response. I’ve put together a minimal example, although it is probably easiest if I attach my mesh file and corresponding mean velocity file, since I have already done the Matlab manipulation mentioned above (is this possible). Here is the what I run:
from dolfin import *
import numpy as np
mesh = Mesh()
with XDMFFile("hex_mesh_bigdom.xdmf") as infile:
infile.read(mesh)
V = VectorFunctionSpace(mesh, 'P', 2)
Q = FunctionSpace(mesh, 'P', 1)
f1 = Function(V)
f_in = XDMFFile('hex_velocity_MEAN_u.xdmf')
umean = f_in.read_checkpoint(f1, "u_mean", 0)
However, I get the following:
>>> type(umean)
<class 'NoneType'>
Even though everything in the xdmf and h5 files containing the mean velocity field has the same format as the files produced when running the original ‘write_checkpoint’ command.
Happy to send the input files over to see if you can reproduce this error.
Thanks,
Ben
I would suggest you put them in a public location (for instance github gist) or google drive and share the link here.
Please supply both the original velocity file, and the manipulated file.
The mesh, original velocity and time average velocity files in the link below:
https://drive.google.com/drive/folders/1s2ItHUF2Y4SaQDQs5KK-Kxt9J3XB8Lix?usp=sharing
It it worth noting the data was only saved every 22 time steps in the original data, this just helped save on memory/simulation time.
You should note that you read data into f1 here, thus you can check the output of f1 pre and post reading the data:
print(f1.vector().get_local())
f_in.read_checkpoint(f1, "u_mean", 0)
print(f1.vector().get_local())
As you can see, the data is successfully read in
Perfect, thank you for your help. I thought it would be reading the function into the variable umean. All works well know.