Parallelized Interpolation from mesh to submesh

dolfinx
, ,
1

I am following the interpolation process from mesh to submesh as listed at dolfinx/python/test/unit/fem/test_interpolation.py at e4439ccca81b976d11c6f606d9c612afcf010a31 · FEniCS/dolfinx · GitHub. However, the speed of execution is slow. Is there a way to speed up the interpolation?

2

Are you using the main branch? If so, have a look at:

it fixes a regression that I introduced when making the interpolation more stable.

If you interpolation is from a mesh to a sub mesh and between the same space, there are much faster ways of doing the interpolation than using this operator.

If you provide an example with a built in mesh I can sketch out what I would do

1 Like
3

No, am using conda version 0.7.1. If by builtin mesh you mean mesh created using dolfinx then yeah.

import dolfinx
import numpy as np

from dolfinx.fem import (Expression, Function, FunctionSpace, assemble_scalar,
                         create_nonmatching_meshes_interpolation_data, form)
from dolfinx.geometry import bb_tree, compute_collisions_points
from dolfinx.mesh import (CellType, create_mesh, create_rectangle,
                          create_unit_cube, create_unit_square,
                          locate_entities, locate_entities_boundary, meshtags)
from mpi4py import MPI

xtype = np.float64

mesh1 = dolfinx.mesh.create_rectangle(comm=MPI.COMM_WORLD,
                            points=((0.0, 0.0), (1.0, 1.0 )), n=(2, 2),
                            cell_type=dolfinx.mesh.CellType.quadrilateral,)

def submesh_marker(x):
    return x[0] <= 0.5


submesh_entities = dolfinx.mesh.locate_entities(mesh1, mesh1.topology.dim, marker = lambda x: submesh_marker(x))
mesh2 = dolfinx.mesh.create_submesh(mesh1, mesh1.topology.dim, submesh_entities)[0]

# Test interpolation from mesh1 to mesh2

u1 = Function(FunctionSpace(mesh1, ("CG", 1)), name="u1", dtype=xtype)
u2 = Function(FunctionSpace(mesh2, ("CG", 1)), name="u2", dtype=xtype)

def f_test1(x):
    return 1.0 - x[0] * x[1]

u1.interpolate(f_test1)
u1.x.scatter_forward()

u1_2_u2_nmm_data = \
    create_nonmatching_meshes_interpolation_data(
        u2.function_space.mesh._cpp_object,
        u2.function_space.element,
        u1.function_space.mesh._cpp_object)

u2.interpolate(u1, nmm_interpolation_data=u1_2_u2_nmm_data)
u2.x.scatter_forward()
4

The speed of the function has been improved in the aforementioned PR, and is released in Release v0.7.2 · FEniCS/dolfinx · GitHub

To interpolate first order Lagrange from a parent mesh to a submesh, you can use:

from IPython import embed
# Author: Jørgen S. Dokken
# SPDX-License-Identifier: MIT

from mpi4py import MPI
import numpy as np
import dolfinx


mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 10, 10)


V = dolfinx.fem.functionspace(mesh, ("Lagrange", 1, (2, )))
u = dolfinx.fem.Function(V)

u.interpolate(lambda x: (x[0], np.sin(x[1])))
u.x.scatter_forward()


def submesh_marker(x):
    return x[0] <= 0.5


submesh_entities = dolfinx.mesh.locate_entities(
    mesh, mesh.topology.dim, marker=lambda x: submesh_marker(x))
submesh, cell_map, _, _ = dolfinx.mesh.create_submesh(
    mesh, mesh.topology.dim, submesh_entities)

V_sub = dolfinx.fem.functionspace(submesh, V.ufl_element())
u_sub = dolfinx.fem.Function(V_sub)

num_sub_cells = submesh.topology.index_map(submesh.topology.dim).size_local
for cell in range(num_sub_cells):
    sub_dofs = V_sub.dofmap.cell_dofs(cell)
    parent_dofs = V.dofmap.cell_dofs(cell_map[cell])
    assert V_sub.dofmap.bs == V.dofmap.bs
    for parent, child in zip(parent_dofs, sub_dofs):
        for b in range(V_sub.dofmap.bs):
            u_sub.x.array[child*V_sub.dofmap.bs +
                          b] = u.x.array[parent*V.dofmap.bs+b]

u_sub.x.scatter_forward()
with dolfinx.io.XDMFFile(submesh.comm, "submesh.xdmf", "w") as xdmf:
    xdmf.write_mesh(submesh)
    xdmf.write_function(u_sub)


with dolfinx.io.XDMFFile(submesh.comm, "parentmesh.xdmf", "w") as xdmf:
    xdmf.write_mesh(mesh)
    xdmf.write_function(u)

Note: Some submesh bugs have been fixed on the main branch by @jpdean:

5

Thank you! I noticed the performance improvements!