L2 projection of fine mesh solution to a coarse mesh

dolfinx
1

I am working on a project, where I came across doing L2 projection of a fine mesh solution to a coarse mesh in fenicsx. I have tried it multiple times but it is throwing this error "not supported between instances of ‘Mesh’ and ‘Mesh’ ". I am compeletely new to fenicsx, I was wondering if there is any alternate project function in fenicsx like in fenics or any example of doing L2 projection from fine mesh to a coarse mesh. Also, please let me know if it is supported in fenicsx at all or any prospective work is going on in this direction.

2

L_2 projection between nonmatching meshes is entirely possible with FEniCS. There is no built in implementation because the operation is dependent on the users’ needs. There is no “catch all” procedure.

3

I understand there doesnt need to be built in implementation but I tried to solve variational problem for projection like this:
domain = create_unit_square(MPI.COMM_WORLD, 50, 50)
V_coarse = FunctionSpace(domain, (“CG”, 1))
u_coarse = TrialFunction(V_coarse)
v_coarse = TestFunction(V_coarse)

u_fine = fem.Function(V_fine)
u_fine.interpolate(uh)

a = u_coarsev_coarsedx
L = u_finev_coarsedx

error is like this

TypeError Traceback (most recent call last)
Cell In[28], line 2
1 a = u_coarsev_coarsedx
----> 2 L = u_finev_coarsedx

File /usr/lib/python3/dist-packages/ufl/measure.py:409, in Measure.rmul(self, integrand)
407 domain = self.ufl_domain()
408 if domain is None:
→ 409 domains = extract_domains(integrand)
410 if len(domains) == 1:
411 domain, = domains

File /usr/lib/python3/dist-packages/ufl/domain.py:340, in extract_domains(expr)
338 for t in traverse_unique_terminals(expr):
339 domainlist.extend(t.ufl_domains())
→ 340 return sorted(join_domains(domainlist))

TypeError: ‘<’ not supported between instances of ‘Mesh’ and ‘Mesh’

4

I would suggest interpolating the function from one mesh to the other by using a Quadrature space. The then interpolated function can in turn be used in a L2 projection of the given quadrature order

5

Can you please elaborate bit more on using quadrature space. Here ‘‘uh’’ is a solution function corresponding to a pde solved on fine mesh(say 100 by 100). It needs to live in same space(fine space) if I interpolate it first into coarse space( on a mesh say 50 by 50) then do projection that doesnt make sense.

6

At some point you have to evaulate the function at the quadrature points of the coarse grid when doing the projection. The following code illlustrates how to compute the L2 projection from a fine grid to a coarse grid by:

  1. Get values of the fine function at the quadrature points of the coarse grid.
  2. Then do the projection
# Illustration of projection from one mesh to another
# Author: Jørgen S. Dokken
# SPDX License: MIT
import ufl
import dolfinx
from mpi4py import MPI
import basix.ufl
import dolfinx.fem.petsc
mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 10, 10)

# Compute function on fine mesh
V = dolfinx.fem.functionspace(mesh, (("Lagrange", 1)))
uh = dolfinx.fem.Function(V)
uh.interpolate(lambda x: x[0]+2*x[1])

# Create coarse mesh
coarser_mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 5, 5)
V_coarse = dolfinx.fem.functionspace(coarser_mesh, ("Lagrange", 1))


# Compute function for fine grid at quadrature points of the coarse grid
degree = 4
Qe = basix.ufl.quadrature_element(
    coarser_mesh.topology.cell_name(), degree=degree)
V_quadrature = dolfinx.fem.functionspace(coarser_mesh, Qe)
q_e = dolfinx.fem.Function(V_quadrature)
nmmid = dolfinx.fem.create_nonmatching_meshes_interpolation_data(V_quadrature.mesh._cpp_object,
                                                                 V_quadrature.element,
                                                                 V.mesh._cpp_object, padding=1e-5)

q_func = dolfinx.fem.Function(V_quadrature)
q_func.interpolate(uh, nmm_interpolation_data=nmmid)

# Project fine function at quadrature points to coarse grid
u = ufl.TrialFunction(V_coarse)
v = ufl.TestFunction(V_coarse)
a_coarse = ufl.inner(u, v) * ufl.dx
L_coarse = ufl.inner(q_func, v)*ufl.dx
problem = dolfinx.fem.petsc.LinearProblem(a_coarse, L_coarse)
u_coarse = problem.solve()

with dolfinx.io.VTXWriter(coarser_mesh.comm, "u_coarse.bp", [u_coarse], engine="BP4") as bp:
    bp.write(0.0)
1 Like
7

Thank you dr @dokken for your prompt response, one last question I would like to ask is I don’t have exact expression for the fine mesh solution, can I directly interpolate my solution function(uh) then compute it at quadrature points of the coarse grid and then project or I need to evaluate my fine mesh solution (by using eval_method) then compute it at quadrature points of the coarse grid and then project. I am truly grateful for your kind response.

8

Not really sure what you are trying to say here.
I used an analytical expression as an illustration:

This does of course not have to be how uh is defined.
Any function coming from the definition

would work (if it is obtained by solving a PDE does not matter).

9

Thank you dr @dokken, for your help. I appreciate it.