Changing mesh coordinates

Users may want to change the coordinates of an existing mesh object for certain reasons. The coordinates can be accessed as a Function through mesh.coordinates where mesh is a mesh object. For example,

mesh.coordinates.dat.data[:, 1] *= 2.0

streches the mesh in the y-direction. Another possibility is to use assign():

Vc = mesh.coordinates.function_space()
x, y = SpatialCoordinate(mesh)
f = Function(Vc).interpolate(as_vector([x, y*2.0]))
mesh.coordinates.assign(f)

This can also be used if f is a solution to a PDE.

Warning

Features which rely on the coordinates field of a mesh’s PETSc DM (usually a DMPlex) such as VertexOnlyMesh() and MeshHierarchy() will not work as expected if the mesh.coordinates field has been modified: at present, the this does not correspondingly update the coordinates field of the DM. This will be fixed in a future Firedrake update.

Changing the coordinate function space

For more complicated situations, one might wish to replace the mesh coordinates with a field which lives on a different FunctionSpace (e.g. higher-order meshes).

Note

Re-assigning the coordinates property of a mesh used to be an undocumented feature. However, this no longer works:

mesh.coordinates = f  # Raises an exception

Instead of re-assigning the coordinates of a mesh, one can create new mesh object from a field f:

new_mesh = Mesh(f)

new_mesh has the same mesh topology as the original mesh, but its coordinate values and coordinate function space are from f. The coordinate function space must be a rank-1 FunctionSpace, constructed either with VectorFunctionSpace(), or by providing a VectorElement to FunctionSpace(). For efficiency, the new mesh object shares data with f. That is, changing the values of f will change the coordinate values of the mesh, and vice versa. If this behaviour is undesired, one should explicitly copy:

g = Function(f)  # creates a copy of f
new_mesh = Mesh(g)

Or simply:

new_mesh = Mesh(Function(f))

Immersing a mesh in higher dimensional space

A useful special case of creating a mesh on modified coordinates is to immerse a lower dimensional mesh in a higher dimension, for example to create a mesh of a two-dimensional manifold immersed in 3D.

This is accomplished by setting the value dimension of the new VectorFunctionSpace() to that of the space in which it should be immersed. For example, a mesh of square bent into a sine wave using linear (flat) elements can be created with:

mesh = UnitSquareMesh(10, 10)
x, y = SpatialCoordinate(mesh)
coord_fs = VectorFunctionSpace(mesh, "CG", 1, dim=3)
new_coord = assemble(interpolate(as_vector([x, y, sin(2*pi*x)]), coord_fs))
new_mesh = Mesh(new_coord)

A sine-wave shaped triangle mesh immersed in three-dimensional space.

Replacing the mesh geometry of an existing function

Creating a new mesh geometry object, as described above, leaves any existing Functions untouched – they continue to live on their original mesh geometries. One may wish to move these functions over to the new mesh. To move f over to mesh, use:

g = Function(functionspaceimpl.WithGeometry.create(f.function_space(), mesh),
             val=f.topological)

This creates a Function g which shares data with f, but its mesh geometry is mesh.

Warning

The example above uses Firedrake internal APIs, which might change in the future.