Hi there,
is it possible in fenics to save/load a function space?
Hi there,
is it possible in fenics to save/load a function space?
What do you imply when you say saving/loading a function space? What data structures would you expect to write/load from file?
Hello everyone,
I am using the HDF5 file format to save all relevant data (functions/mesh/meshfunctions) + metadata in a self-contained manner. For example,
from fenics import *
mesh = UnitSquareMesh(4,4)
def WriteH5Sample():
V = FunctionSpace(mesh,"CG",1)
u = project(Expression(("x[0]*x[1]"), degree=1), V)
with HDF5File(MPI.comm_world, 'test.h5', file_mode='w') as hdf5file:
hdf5file.write(u, 'f')
WriteH5Sample() # Write out a simple H5 file
I can then read back in the function using
def ReadH5Func(V):
with HDF5File(MPI.comm_world, 'test.h5', file_mode='r') as hdf5file:
# V = ???
u = Function(V)
hdf5file.read(u, 'f')
return u
# Define FunctionSpace.
# Can this be automatically interpreted from the data?
newV = FunctionSpace(mesh, 'CG', 1)
u = ReadH5Func(newV)
The problem is that we must decide on a FunctionSpace newV
prior to reading the data. ReadH5(V)
takes V
as an argument, so it is not really self-contained.
How can I go about “saving the function space” so that the appropriate V
is generated automatically? Perhaps I can just save information from str(V)
and use it to create newV
?
Thanks for any advice!
This is my attempt at solving my own problem. It is admittedly clunky, but I wanted it to be robust enough to recreate function spaces (scalar/vector/tensor) as well as mixed function spaces from a small piece of metadata.
If anyone would like to revise my code, point out any flaws, or suggest any improvements that would be great.
from fenics import *
mesh = UnitSquareMesh(4,4)
Q = FunctionSpace(mesh, "DG", 2)
V = VectorFunctionSpace(mesh, "CG", 1)
# We will save info about W and use it later to recreate the equivalent functionspace _W.
W = MixedFunctionSpace(Q,V)
u = Function(W)
def FunctionSpaceInfo(V):
E = V.ufl_element()
info = {'family':E.family(),'cell':E.cell(),'degree':E.degree(),
'num_sub_elements':E.num_sub_elements(),'symmetry':E.symmetry(),
'value_shape':E.value_shape(),'value_size':E.value_size()}
return info
def SaveFunctionSpace(V):
saveinfo = []
if V.num_sub_spaces()>0:
for V_i in (V.sub_spaces()):
saveinfo.append(FunctionSpaceInfo(V_i))
else:
saveinfo.append(FunctionSpaceInfo(V))
return saveinfo
def RecreateFunctionSpace(info):
FS_list = []
for FS in info:
if FS['value_size']==1: # scalar
FS_list.append( FunctionSpace(mesh,FS['family'],FS['degree']))
if FS['value_size']==2: # vector
FS_list.append(VectorFunctionSpace(mesh,FS['family'],FS['degree']))
if FS['value_size']==4: # tensor
FS_list.append(TensorFunctionSpace(mesh,FS['family'],FS['degree']))
return FS_list
# Extract info from FunctionSpace(s)
saveinfo = SaveFunctionSpace(W)
# Recreate FunctionSpace
FS_list = RecreateFunctionSpace(saveinfo)
_W = FS_list[0] if len(saveinfo)==1 else MixedFunctionSpace(*FS_list)
_u = Function(_W)
# This is just a reminder to split a Function defined in a MixedFunctionSpace
if _u.num_sub_spaces()>1: print(_u.sub(0).name())
else: print(_u.name())