H1 projection from fine mesh to coarse mesh

General
1

Hello folks, hope you all are doing well. I would appreciate if you can offer any help regarding this issue. I was trying to do a H1 projection from fine mesh to coarse mesh but after computing fine mesh sol in quadrature space, fenicsx doesn’t allow to calculate gradient of such function, I would appreciate if you can provide any insight egarding how to achieve this. Here is the minimal code example and error.

import ufl
import dolfinx
from mpi4py import MPI
import basix.ufl
import dolfinx.fem.petsc
from ufl import inner, grad, dx
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,
                                                                 V_quadrature.element,
                                                                 V.mesh, 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 = inner(u, v)*dx + inner(grad(u),grad(v))*dx 
L_coarse = inner(q_func, v)*dx + inner(grad(q_func),grad(v))*dx
problem = dolfinx.fem.petsc.LinearProblem(a_coarse, L_coarse)
u_coarse = problem.solve()

error

ValueError                                Traceback (most recent call last)
<ipython-input-10-f496e4efcf70> in <cell line: 38>()
     36 a_coarse = inner(u, v)*dx + inner(grad(u),grad(v))*dx
     37 L_coarse = inner(q_func, v)*dx + inner(grad(q_func),grad(v))*dx
---> 38 problem = dolfinx.fem.petsc.LinearProblem(a_coarse, L_coarse)
     39 u_coarse = problem.solve()

16 frames
/usr/local/lib/python3.10/dist-packages/basix/ufl.py in tabulate(self, nderivs, points)
   1532         """
   1533         if nderivs > 0:
-> 1534             raise ValueError("Cannot take derivatives of Quadrature element.")
   1535 
   1536         if points.shape != self._points.shape:

ValueError: Cannot take derivatives of Quadrature element
2

Can anyone please help with this, how to resolve the issue. I have been trying different approaches but nothing worked so far.

3

Please note that this is an open-source forum, where people volunteer. You posted a message on a Sunday, and people have other priorities in life as well.

For your particular query, what I would do would be to:

  1. Interpolate grad(u) into a suitable function space in its respective grid (ie on mesh). For your use case that would be a DG-0 space.
  2. Interpolate the function containing grad(u) into a coarser space on the quadrature mesh.
2 Likes
4

Dr @dokken, Thanks for taking your time to respond to my query. I actually tried to do something similar on those lines and I am getting the following shapes do not match error. I once again appreciate your kind reponse.

from types import CellType
import ufl
import dolfinx
from mpi4py import MPI
import basix.ufl
from basix.ufl import element
from dolfinx import default_scalar_type
import dolfinx.fem.petsc
from ufl import inner, grad, dx
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])
# compute gradient and interpolate into DG-0 space
gdim = mesh.geometry.dim
dg = dolfinx.fem.functionspace(mesh,('DG',0,(gdim,)))
grad_int = dolfinx.fem.Function(dg)
grad_uh = dolfinx.fem.Expression(grad(uh),dg.element.interpolation_points())
grad_int.interpolate(grad_uh)

# 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)
nmmid = dolfinx.fem.create_nonmatching_meshes_interpolation_data(V_quadrature.mesh,
                                                                 V_quadrature.element,
                                                                 V.mesh, padding=1e-5)

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

# Compute gradient 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)
nmmid1 = dolfinx.fem.create_nonmatching_meshes_interpolation_data(V_quadrature.mesh,
                                                                 V_quadrature.element,
                                                                 dg.mesh, padding=1e-5)
g_func = dolfinx.fem.Function(V_quadrature)
g_func.interpolate(grad_int, nmm_interpolation_data=nmmid1)

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

error

ValueError                                Traceback (most recent call last)
<ipython-input-49-90b75379c062> in <cell line: 55>()
     53 v = ufl.TestFunction(V_coarse)
     54 a_coarse = inner(u, v)*dx + inner(grad(u),grad(v))*dx
---> 55 L_coarse = inner(q_func, v)*dx + inner(g_func,grad(v))*dx
     56 problem = dolfinx.fem.petsc.LinearProblem(a_coarse, L_coarse)
     57 u_coarse = problem.solve()

1 frames
/usr/local/lib/python3.10/dist-packages/ufl/tensoralgebra.py in __new__(cls, a, b)
    164         ash, bsh = a.ufl_shape, b.ufl_shape
    165         if ash != bsh:
--> 166             raise ValueError(f"Shapes do not match: {ufl_err_str(a)} and {ufl_err_str(b)}")
    167 
    168         # Simplification

ValueError: Shapes do not match: <Coefficient id=138824788492416> and <Grad id=138824788885376>
5

You cannot use the same space for u and grad_u in quadrature , as you need to set the Qe to be a vector element;

I.e. add shape to Qe here, as done in

1 Like
6

Thank you so much Dr. @dokken, got it work. Again, I appreciate your response and being patient with me.