Update to latest version of dolfinx

.cspell_dict.txt

@@ -3,6 +3,8 @@ allreduce
 argsort
 astype
 atol
+basix
+cellname
 cofac
 dirichletbc
 dofs
@@ -11,6 +13,7 @@ fenics
 fenicsx
 fenicsx_pulse
 finsberg
+functionspace
 gradu
 Holzapfel
 holzapfelogden

demo/unit_cube.ipynb

@@ -192,6 +192,7 @@
     ") -> list[dolfinx.fem.bcs.DirichletBC]:\n",
     "    V, _ = state_space.sub(0).collapse()\n",
     "    facets = geo.facet_tags.find(1)  # Specify the marker used on the boundary\n",
+    "    mesh.topology.create_connectivity(mesh.topology.dim - 1, mesh.topology.dim)\n",
     "    dofs = dolfinx.fem.locate_dofs_topological((state_space.sub(0), V), 2, facets)\n",
     "    u_fixed = dolfinx.fem.Function(V)\n",
     "    u_fixed.x.array[:] = 0.0\n",

src/fenicsx_pulse/mechanicsproblem.py

@@ -5,6 +5,7 @@
 import dolfinx.fem.petsc
 import dolfinx.nls.petsc
 import ufl
+import basix
 
 from . import kinematics
 from .boundary_conditions import BoundaryConditions
@@ -13,7 +14,7 @@
 
 
 @dataclass(slots=True)
-class MechanicsProblem:
+class BaseMechanicsProblem:
     model: CardiacModel
     geometry: Geometry
     bcs: BoundaryConditions = field(default_factory=BoundaryConditions)
@@ -22,49 +23,33 @@ class MechanicsProblem:
     state_space: dolfinx.fem.FunctionSpace = field(init=False, repr=False)
     state: dolfinx.fem.Function = field(init=False, repr=False)
     test_state: dolfinx.fem.Function = field(init=False, repr=False)
-    _virtual_work: ufl.form.Form = field(init=False, repr=False)
-    _dirichlet_bc: typing.Sequence[dolfinx.fem.bcs.DirichletBC] = field(
+    virtual_work: ufl.form.Form = field(init=False, repr=False)
+    _dirichlet_bc: typing.Sequence[dolfinx.fem.bcs.DirichletBC] | None = field(
+        default=None,
         init=False,
         repr=False,
     )
 
+    @property
+    def dirichlet_bc(self) -> typing.Sequence[dolfinx.fem.bcs.DirichletBC]:
+        return self._dirichlet_bc or []
+
+    @dirichlet_bc.setter
+    def dirichlet_bc(self, value: typing.Sequence[dolfinx.fem.bcs.DirichletBC]) -> None:
+        self._dirichlet_bc = value
+        self._set_dirichlet_bc()
+        self._init_solver()
+
     def __post_init__(self):
         self._init_space()
         self._init_form()
         self._init_solver()
 
-    def _init_space(self) -> None:
-        P2 = ufl.VectorElement("Lagrange", self.geometry.mesh.ufl_cell(), 2)
-        P1 = ufl.FiniteElement("Lagrange", self.geometry.mesh.ufl_cell(), 1)
-
-        self.state_space = dolfinx.fem.FunctionSpace(self.geometry.mesh, P2 * P1)
-        self.state = dolfinx.fem.Function(self.state_space)
-        self.test_state = ufl.TestFunction(self.state_space)
-
-    def _init_form(self) -> None:
-        u, p = ufl.split(self.state)
-        v, _ = ufl.split(self.test_state)
-
-        self.model.compressibility.register(p)
-
-        F = kinematics.DeformationGradient(u)
-        psi = self.model.strain_energy(F, p)
-        self._virtual_work = ufl.derivative(
-            psi * self.geometry.dx,
-            coefficient=self.state,
-            argument=self.test_state,
-        )
-        external_work = self._external_work(u, v)
-        if external_work is not None:
-            self._virtual_work += external_work
-
-        self._set_dirichlet_bc()
-
     def _init_solver(self) -> None:
         self._problem = dolfinx.fem.petsc.NonlinearProblem(
-            self._virtual_work,
+            self.virtual_work,
             self.state,
-            self._dirichlet_bc,
+            self.dirichlet_bc,
         )
         self._solver = dolfinx.nls.petsc.NewtonSolver(
             self.geometry.mesh.comm,
@@ -114,3 +99,43 @@ def _set_dirichlet_bc(self) -> None:
 
     def solve(self):
         return self._solver.solve(self.state)
+
+
+@dataclass(slots=True)
+class MechanicsProblem(BaseMechanicsProblem):
+    def _init_space(self) -> None:
+        P2 = basix.ufl.element(
+            family="Lagrange",
+            cell=self.geometry.mesh.ufl_cell().cellname(),
+            degree=2,
+            shape=(self.geometry.mesh.ufl_cell().topological_dimension(),),
+        )
+        P1 = basix.ufl.element(
+            family="Lagrange",
+            cell=self.geometry.mesh.ufl_cell().cellname(),
+            degree=1,
+        )
+        element = basix.ufl.mixed_element([P2, P1])
+
+        self.state_space = dolfinx.fem.functionspace(self.geometry.mesh, element)
+        self.state = dolfinx.fem.Function(self.state_space)
+        self.test_state = ufl.TestFunction(self.state_space)
+
+    def _init_form(self) -> None:
+        u, p = ufl.split(self.state)
+        v, _ = ufl.split(self.test_state)
+
+        self.model.compressibility.register(p)
+
+        F = kinematics.DeformationGradient(u)
+        psi = self.model.strain_energy(F, p)
+        self.virtual_work = ufl.derivative(
+            psi * self.geometry.dx,
+            coefficient=self.state,
+            argument=self.test_state,
+        )
+        external_work = self._external_work(u, v)
+        if external_work is not None:
+            self.virtual_work += external_work
+
+        self._set_dirichlet_bc()

tests/conftest.py

@@ -10,12 +10,12 @@ def mesh():
 
 @pytest.fixture(scope="session")
 def P1(mesh):
-    return dolfinx.fem.FunctionSpace(mesh, ("Lagrange", 1))
+    return dolfinx.fem.functionspace(mesh, ("Lagrange", 1))
 
 
 @pytest.fixture(scope="session")
 def P2(mesh):
-    return dolfinx.fem.VectorFunctionSpace(mesh, ("Lagrange", 2))
+    return dolfinx.fem.functionspace(mesh, ("Lagrange", 2, (mesh.geometry.dim,)))
 
 
 @pytest.fixture

tests/test_mechanicsproblem.py

@@ -37,6 +37,7 @@ def dirichlet_bc(
     ) -> list[dolfinx.fem.bcs.DirichletBC]:
         V, _ = state_space.sub(0).collapse()
         facets = geo.facet_tags.find(1)
+        mesh.topology.create_connectivity(mesh.topology.dim - 1, mesh.topology.dim)
         dofs = dolfinx.fem.locate_dofs_topological((state_space.sub(0), V), 2, facets)
         u_fixed = dolfinx.fem.Function(V)
         u_fixed.x.array[:] = 0.0