I need to store a solution of expensive FEM calculation to use it in further analysis. Browsing through tutorials I have so far discovered, that I may store my results like this:
from fenics import *
mesh = Mesh('mesh/UnitSquare8x8.xml')
V = FunctionSpace(mesh, 'P', 1)
u = TrialFunction(V)
v = TestFunction(V)
f = Constant(-6.0)
a = dot(grad(u), grad(v))*dx
L = f*v*dx
u_D = Expression('1 + x[0]*x[0] + 2*x[1]*x[1]', degree=2)
def boundary(x, on_boundary):
return on_boundary
bc = DirichletBC(V, u_D, boundary)
A = assemble(a)
b = assemble(L)
bc.apply(A, b)
u = Function(V)
solver = KrylovSolver("cg", "ilu")
solver.solve(A, u.vector(), b)
File('solution.xml') << u.vector()
and later load them like this:
from fenics import *
mesh = Mesh('mesh/UnitSquare8x8.xml')
V = FunctionSpace(mesh, 'P', 1)
u = Function(V)
File('solution.xml') >> u.vector()
Unfortunately I hardly know what exactly I am doing here. Is that a proper way of storing and loading calculated results? Is order of elements in u.vector() (for the same mesh file) fixed within/between different FEniCS versions, or it is just an implementation detail which may change any time? If it is unsafe, then what is the proper way of doing so?
I have found another (possibly even more dangerous) solution. I may use VALUES = u.vector().get_local() and u.vector().set_local(VALUES) methods, as VALUES is a numpy array which I may easily store and load.
No, according to answer to Is the order of Vector elements preserved between runs? the order of Vector elements is not guarranted to be preserved.
It is recommended to use XDMFFile.write_checkpoint() and XDMFFile.read_checkpoint() methods instead.