Interpolate Function into a higher-order function space

dolfinx
1

Hello,

I’m currently trying to interpolate the solution of my (nonlinear) problem U into a higher order function space U_target. I therefore created a class which contains the solution U and a function which evaluates U at positions x. This function was indented as input for creating U_target but - even if the evaluated nodal values of U are what they should be - the interpolation fails. I attached a minimum working example where a expression is interpolated into a function space of order 2 while this function should be subsequently interpolated into a higher order functionspace.

import numpy as np
from mpi4py import MPI
from petsc4py import PETSc

import dolfinx
import dolfinx.geometry
import dolfinx.io
import ufl
from dolfinx.cpp.mesh import CellType

class RefSol:
    def __init__(self, t):
        self.t = t

    def interpol_displacement(self, x):
        u = np.zeros((2, x.shape[1]), dtype=PETSc.ScalarType)
        u[0] = np.exp(-self.t) * (- np.cos(np.pi * x[0]) * np.sin(np.pi * x[1]))
        u[1] = np.exp(-self.t) * (- np.cos(np.pi * x[0]) * np.sin(np.pi * x[1]))
        return u

class Interpol_PointValue:
    def __init__(self, mesh, function):
        self.function = function
        self.mesh     = mesh

    def eval_pointvalue(self, x):
        "Code from https://jorgensd.github.io/dolfinx-tutorial/chapter1/membrane_code.html"

        # Create searchable data structure of current mesh
        bb_tree = dolfinx.geometry.BoundingBoxTree(self.mesh, self.mesh.topology.dim)

        # Initialize subset arrays
        cells = []
        points_on_proc = []

        # Select points on current processor
        for point in x.T:
            # Find cells that are close to the point
            cell_candidates = dolfinx.geometry.compute_collisions_point(bb_tree, point)
            # Choose one of the cells that contains the point
            cell = dolfinx.geometry.select_colliding_cells(self.mesh, cell_candidates, point, 1)
            # Only use evaluate for points on current processor
            if len(cell) == 1:
                points_on_proc.append(point)
                cells.append(cell[0])

        # Interpolate values at points on current processor
        points_on_proc = np.array(points_on_proc, dtype=np.float64)
        values = self.function.eval(points_on_proc, cells)

        # DEBUG: Compare output of evaluation of Expression with values from function
        #IntFunc    = RefSol(0.5)
        #values_ref = IntFunc.interpol_displacement(x)
        #error      = values_ref - values.transpose()

        return values.transpose()

# Set mesh
mesh = dolfinx.RectangleMesh(MPI.COMM_WORLD, [[0, 0, 0], [1, 1, 0]], [10, 10], cell_type=CellType.triangle)

# Function spaces
Pu       = ufl.VectorElement('CG', mesh.ufl_cell(), 2)
Pu_t     = ufl.VectorElement('CG', mesh.ufl_cell(), 4)
V        = dolfinx.FunctionSpace(mesh, Pu)
V_target = dolfinx.FunctionSpace(mesh, Pu_t)

uh       = dolfinx.Function(V_target)
u_target = dolfinx.Function(V_target)

# Interpolate functions
IntFunc = RefSol(0.5)
uh.interpolate(IntFunc.interpol_displacement)
uh.x.scatter_forward()

PointVal = Interpol_PointValue(mesh, uh)
u_target.interpolate(PointVal.eval_pointvalue)
u_target.x.scatter_forward()

with dolfinx.io.XDMFFile(MPI.COMM_WORLD, "test_uint-uh.xdmf", "w") as xdmf:
    xdmf.write_mesh(mesh)
    xdmf.write_function(uh)

with dolfinx.io.XDMFFile(MPI.COMM_WORLD, "test_uint-utarget.xdmf", "w") as xdmf:
    xdmf.write_mesh(mesh)
    xdmf.write_function(u_target)

Thanks!

Best wishes,
Maximilian

2

Hi Maximilian ,

I’m currently working on a functionality to interpolate between spaces of the same family (Eg: CG to CG) but of a different order (2 to 4, or 4 to 2).
I expect it to land soon in dolfinx.

The proposed interface is something like:

u_target.interpolate(uh)

As soon as I create the PR Ill write a msg here.

3 Likes
3

Hi Igor,

thanks for your reply. Such a functionality would be highly desirable.

I did some additional testing yesterday, and what currently fixed my problem is to directly modify the nodal values of the (higher order) function. Therefore, I used the coordinates of the degrees of freedom, interpolated the values of uh and manipulate the structure of the resulting array before copying them into the higher order function. Maybe this can help someone, until the interpolation function is extended.

Best wishes,
Maximilian

4

Hello @IgorBaratta and @mbr,

Do we have a functionality in Dolfinx now to interpolate from DG0 to CG1 function space?
If yes, what is the corresponding procedure to do so?
Thanks!

5

Yes, simply use

DOLFINx also supports interpolation between nonmatching meshes

2 Likes