How are 2D manifolds immersed/embedded in 3D space handled in dolfinx?

dokken · August 8, 2025, 6:42am

If one want’s to work with H-div conforming spaces on orientable manifolds, a slight workaround is needed, see:

[BUG]: Piola mapping on orientable manifolds are broken.

opened 07:37PM - 07 Aug 25 UTC

### Summarize the issue Consider the MWE from legacy DOLFIN, where div conformi…ng spaces works on orientable manifolds: ```python import numpy as np from dolfin import * orient = True meshid = "sphere_ico3" meshname = "meshes/%s.xml.gz" % meshid print("Computing on %s" % meshname) mesh = Mesh(meshname) if orient: mesh.init_cell_orientations(Expression(("x[0]", "x[1]", "x[2]"), degree=1)) el_v = FiniteElement("RT", mesh.ufl_cell(), 1) el_q = FiniteElement("DG", mesh.ufl_cell(), 0) W = FunctionSpace(mesh, MixedElement(el_v, el_q)) (sigma, u) = TrialFunctions(W) (tau, v) = TestFunctions(W) a = (div(sigma) * v) * dx A = as_backend_type(assemble(a)).mat() print(f"{orient=}", np.sum(A[:, :]), np.linalg.norm(A[:, :])) ``` which yields ```python orient=True 0.0 30.983866769659336 orient=False 320.0 30.983866769659336 ``` as one observes, we get identical norms, but signs are changed, which corresponds with what is stated in @meg-simula's paper: https://gmd.copernicus.org/articles/6/2099/2013/ The orientable version is what we should get in DOLFINx as well. Equivalent code: ```python from mpi4py import MPI import meshio import numpy as np import ufl import basix.ufl import dolfinx tmp_data = meshio.dolfin.read("./meshes/sphere_ico3.xml") mesh = dolfinx.mesh.create_mesh( comm=MPI.COMM_WORLD, cells=tmp_data.cells_dict["triangle"], x=tmp_data.points, e=ufl.Mesh(basix.ufl.element("Lagrange", "triangle", 1, shape=(3,))), ) el_v = basix.ufl.element("RT", mesh.basix_cell(), 1) el_q = basix.ufl.element("DG", mesh.basix_cell(), 0) mixed_el = basix.ufl.mixed_element([el_v, el_q]) W = dolfinx.fem.functionspace(mesh, mixed_el) (sigma, u) = ufl.TrialFunctions(W) (tau, v) = ufl.TestFunctions(W) a = (ufl.div(sigma) * v) * ufl.dx A = dolfinx.fem.assemble_matrix(dolfinx.fem.form(a)) A.scatter_reverse() A_d = A.to_dense() print(np.sum(A_d[:, :]), np.linalg.norm(A_d[:, :])) ``` yielding ```bash -32.0 30.983866769659336 ``` Mesh can be gotten from: https://gmd.copernicus.org/articles/6/2099/2013/gmd-6-2099-2013-supplement.zip ### Version main branch

Topic		Replies	Views
Does Dolfinx's solver have an issue or require some preconditioning when handling certain meshes like the Icosahedral mesh? General dolfinx , legacy-dolfin	10	157	September 6, 2025
The problem in a simple example I made and how to 'project' the value assigned to 'FunctionSpace' into a numpy array in dolfinx environment? dolfinx	11	1218	September 8, 2022
[SOLVED] Solve PDEs on a Manifold dolfinx	4	566	February 7, 2025
3D problem (sphere) Navier-Stokes dolfinx	50	5680	February 8, 2022
My feedback thread for website patches General	2	330	October 20, 2023

How are 2D manifolds immersed/embedded in 3D space handled in dolfinx?

Related topics