How to use prism mesh cells in FEniCSx?

dolfinx
1

Hi, I’ve been learning FEniCSx lately. I was trying to use prism cells for an FEM problem, but I got stuck. Following is a minimal code based on the Poisson demo:

from pathlib import Path
from mpi4py import MPI
import ufl
from dolfinx import fem, io, mesh, plot
from dolfinx.fem.petsc import LinearProblem

msh = mesh.create_box(
    comm=MPI.COMM_WORLD,
    points=((0.0, 0.0, 0.0), (2.0, 1.0, 1.0)),
    n=(32, 16, 16),
    cell_type=mesh.CellType.prism,
)

V = fem.functionspace(msh, ("Lagrange", 1))
u = ufl.TrialFunction(V)
v = ufl.TestFunction(V)
x = ufl.SpatialCoordinate(msh)
f = 10 * ufl.exp(-((x[0] - 0.5) ** 2 + (x[1] - 0.5) ** 2) / 0.02)
g = ufl.sin(5 * x[0])
a = ufl.inner(ufl.grad(u), ufl.grad(v)) * ufl.dx
L = ufl.inner(f, v) * ufl.dx + ufl.inner(g, v) * ufl.ds

problem = LinearProblem(
    a,
    L,
    bcs=[],
    petsc_options_prefix="demo_poisson_",
    petsc_options={"ksp_type": "preonly", "pc_type": "lu", "ksp_error_if_not_converged": True},
)
uh = problem.solve()
assert isinstance(uh, fem.Function)

out_folder = Path("out_poisson")
out_folder.mkdir(parents=True, exist_ok=True)
with io.XDMFFile(msh.comm, out_folder / "poisson.xdmf", "w") as file:
    file.write_mesh(msh)
    file.write_function(uh)

This code runs well when I set the cell_type as hexahedron. But when I change the cell_type to prism, it raise an error:
AttributeError: 'LinearProblem' object has no attribute '_solver'
Is there some extra setup I need to get prism cells working properly?

2

The problem with prism cells is that the assembly over the exterior facets are not trivial (i.e. the ds term in L).

If you remove the term

and change the petsc options to:

    petsc_options={"ksp_type": "preonly", "pc_type": "lu", "ksp_error_if_not_converged": True,
                   "pc_factor_mat_solver_type": "mumps"},

the code will successfully run.

I would have to ask @Chris_Richardson for guidance/status on exterior facet integrals on prisms.

3

Thanks for your help.