Since the release of v0.8.0 in May there has been 190 pull requests merged into DOLFINx from 12 contributors. This release marks a milestone in DOLFINx development, as several new features have been added.
New features
- Mixed dimensional assembly
- Simpler constructions of block forms with
ufl.MixedFunctionSpace
andufl.extract_blocks
- Data-independent form compilation with
dolfinx.fem.compile_form
- Native Windows compatibility
API changes
- Improving
dolfinx.mesh.refine
- New function for nonmatching interpolation
- Remove
dolfinx.fem.set_bc
- Rename
scale
inapply_lifting
toalpha
- Move PETSc Vector wrapper from
dolfinx.fem.Function.vector
todolfinx.fem.Function.x.petsc_vec
Python API additions
- Expose
DirichletBC.dof_indices
to Python API, DOLFINx #3389 - Expose
NewtonSolver.set_convergence_check
, DOLFINx #3386
Other updates
- All users of
PETSc.SNES
should upgrade from PETSc<3.21.6 - Performance update for form compilation with subdomains
- Improved documentation of
dolfinx.mesh
- Initial support for mixed topology meshes
- Output of tensors of arbitrary size
- Built-in interpolation matrix
A detailed summary of all pull requests can be found at Release: DOLFINx v0.9.0
New features
Mixed dimensional assembly
Main contributors: Joseph Dean and Jørgen S. Dokken
In DOLFINx, we now support assembly of 0th, 1st and 2nd order tensors consisting of meshes from different sub-meshes.
One can now create function space over sub-meshes, and assemble with coefficients, test-functions or trial-functions from the parent mesh (and vice versa).
The feature currently works for mixed assembly of co-dimension 0 (subset of cells) and co-dimension 1 (subset of facets).
An example of usage is shown below, where a TrialFunction
from a parent mesh is combined in a variational form with a TestFunction
from a sub-mesh of all
exterior facets. We use the integration ds
over the parent mesh, to ensure that we can use quantities such as ufl.FacetNormal
.
domain = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 10, 10)
V = dolfinx.fem.functionspace(domain, ("Lagrange", 1))
# Create submesh and function space on submesh
tdim = domain.topology.dim
domain.topology.create_connectivity(tdim - 1, tdim)
facet_domain, sub_to_parent, _, _ = dolfinx.mesh.create_submesh(
domain, tdim - 1, dolfinx.mesh.exterior_facet_indices(domain.topology)
)
Q = dolfinx.fem.functionspace(domain, ("DG", 2))
u = ufl.TrialFunction(V)
q = ufl.TestFunction(Q)
a_mixed = ufl.inner(u, q) * ufl.ds(domain=domain)
# Invert entity map to go from integration domain to sub mesh
facet_imap = msh.topology.index_map(tdim - 1)
num_facets = facet_imap.size_local + facet_imap.num_ghosts
parent_to_sub = np.full(num_facets, -1)
parent_to_sub[sub_to_parent] = np.arange(len(sub_to_parent))
entity_maps = {facet_domain:parent_to_sub}
# Create dolfinx.fem.Form
a_compiled = dolfinx.fem.form(a_mixed)
Relevant pull requests since last release:
- Codim-1 mixed assembly: DOLFINx #3224
- Interior facet assembly with submesh of codimensional 1: DOLFINx #3452
- Submesh of vertices: DOLFINx #3455
- C++ demo: DOLFINx #3436, DOLFINx #3262
- Packing of coefficients for more than 2 sub-meshes: DOLFINx #3260, DOLFINx #3361
- Fix for block vector assembly: DOLFINx #3346
- Add create submesh documentation in Python: DOLFINx #3112
- Consistent orientation of sub-entities, DOLFINx #3209
- Interpolation to and from co-dim 0 meshes, DOLFINx #3114
Simpler construction of block forms
Main contributor: Jørgen S. Dokken
Related to the introduction of mixed dimensional assembly, having an easier way of creating blocked variational forms is crucial for the usability of the feature. This feature can also be used for “standard” block forms, used with dolfinx.fem.petsc.assemble_*_block/nest
.
With ufl.MixedFunctionSpace
and ufl.extract_blocks
the users can now create simple variational formulations
V = dolfinx.fem.functionspace(domain, ("Lagrange", 2))
Q = dolfinx.fem.functionspace(submesh, ("Lagrange", 1))
# Create mixed problem residual F
W = ufl.MixedFunctionSpace(V, Q)
u = dolfinx.fem.Function(V)
p = dolfinx.fem.Function(Q)
v, q = ufl.TestFunctions(W)
F = ufl.inner(u * ufl.grad(u), ufl.grad(v)) * ufl.dx
F += ufl.inner(p, v) * ufl.ds
F += ufl.inner(1, q) * ufl.ds
F_blocked = ufl.extract_blocks(F)
F_compiled = dolfinx.fem.form(F_blocked)
# Compute Jacobian on block form
du, dp = ufl.TrialFunctions(W)
J_blocked = ufl.extract_blocks(ufl.derivative(F, u, du) + ufl.derivative(F, p, dp))
J_compiled = dolfinx.fem.form(J_blocked)
# Assemble block matrix and block vector
import dolfinx.fem.petsc
A = dolfinx.fem.petsc.assemble_matrix_block(J_compiled)
A.assemble()
b = dolfinx.fem.petsc.assemble_vector_block(F_compiled, a=J_compiled)
# Assemble in petsc nest vector to use PC-FieldSplit
A_nest = dolfinx.fem.petsc.assemble_matrix_nest(J_compiled)
A_nest.assemble()
b_nest = dolfinx.fem.petsc.assemble_vector_nest(F_compiled)
Relevant pull requests since last release:
- Simplify DOLFINx demos with
ufl.extract_block
: DOLFINx #3450 - Extract block fixes: UFL #308, UFL #310, UFL #285
Data-independent form compilation
Main contributor: Jørgen S. Dokken
To unify the C++ and Python interface with UFL and Basix, the user can now create data-independent variational forms for Python (similar to what is done for C++ problems), and then in turn attached data to the compiled form.
This means that the same code can be used for both C++ and Python programs.
This approach could also simplify code when working with mesh-refinement of a sequence of meshes through re-meshing.
real_type = np.float64
dtype = np.complex 128
c_el = basix.ufl.element("Lagrange", "triangle", 1, shape=(2,), dtype=real_type)
domain = ufl.Mesh(c_el)
el = basix.ufl.element("Lagrange", "triangle", 2, dtype=real_type)
V = ufl.FunctionSpace(domain, el)
u = ufl.Coefficient(V)
w = ufl.Coefficient(V)
c = ufl.Constant(domain)
e = ufl.Constant(domain)
J = c * e * u * w * ufl.dx(domain=domain)
# Compile form using dolfinx.jit.ffcx_jit
compiled_form = dolfinx.fem.compile_form(
MPI.COMM_WORLD, J, form_compiler_options={"scalar_type": dtype}
)
for N in [2, 4, 8, 16]:
mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, N, N, dtype=real_type)
Vh = dolfinx.fem.functionspace(mesh, u.ufl_element())
uh = dolfinx.fem.Function(Vh, dtype=dtype)
uh.interpolate(lambda x: x[0])
wh = dolfinx.fem.Function(Vh, dtype=dtype)
wh.interpolate(lambda x: x[1])
eh = dolfinx.fem.Constant(mesh, dtype(3.0))
ch = dolfinx.fem.Constant(mesh, dtype(2.0))
# Attach data to create a `dolfinx.fem.Form` that can be used in assembly
form = dolfinx.fem.create_form(compiled_form, [], mesh, {u: uh, w: wh}, {c: ch, e: eh})
Relevant pull requests since last release:
- Form independent compilation:DOLFINx #3263
- Form independent compilation with submeshes: DOLFINx #3262.
Native Windows compatibility
Main contributors: Jack Hale, Chris Richardson and Min Ragan-Kelley
There is now support in DOLFINx for native Windows builds with MPI (without PETSc). You can install DOLFINx with conda on Windows, see: https://github.com/FEniCS/dolfinx/?tab=readme-ov-file#conda for details.
Example use-case for solving PDEs without PETSc can be found in the PYAMG demo by Chris Richardson and solving PDEs with different scalar types by Garth N. Wells
Relevant pull requests since last release:
- DOLFINx Windows support DOLFINx #3198
- FFCx Windows support FFCx #689
- Basix Windows support Basix #819
- Conda feedstock: Conda Feedstock #78
API updates
Improve and unify the mesh refinement interface
Main contributor: Paul Kühner
Work by Paul during his Google Summer of Code internship on multigrid.
As part of this work he improved the dolfinx.mesh.refine
interface to return cell and facet relation, to easily transfer cell and facet tags to the refined mesh.
There is now also the possibility of passing your own partitioner to the refine function if you want to redistribute the mesh.
There is now also support of refinement of 1D meshes.
Example of how to transfer meshtags for cells and facets with the new interface is shown below
# Create and uniformly refine mesh
domain = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 10, 10)
tdim = domain.topology.dim
domain.topology.create_entities(1)
refined_mesh, parent_cell, parent_facet = dolfinx.mesh.refine(
domain, option=dolfinx.mesh.RefinementOption.parent_cell_and_facet
)
# Create cell tag
domain.topology.create_entities(tdim)
cell_indices = np.arange(domain.topology.index_map(tdim).size_local, dtype=np.int32)
values = cell_indices.copy()
meshtag = dolfinx.mesh.meshtags(
domain,
tdim,
cell_indices,
values,
)
# Transfer cell tag
refined_mesh.topology.create_connectivity(tdim, tdim - 1)
refined_meshtag = dolfinx.mesh.transfer_meshtag(meshtag, refined_mesh, parent_cell)
# Create facet tag
domain.topology.create_entities(tdim - 1)
facet_indices = np.arange(
domain.topology.index_map(tdim - 1).size_local, dtype=np.int32
)
values = facet_indices.copy()
facet_meshtag = dolfinx.mesh.meshtags(
domain,
tdim - 1,
facet_indices,
values,
)
# Transfer facet tag to refined mesh
refined_facettag = dolfinx.mesh.transfer_meshtag(
facet_meshtag, refined_mesh, parent_cell, parent_facet
)
Relevant pull requests since last release:
- 1D mesh-refinement: DOLFINx #3314
- Support graph partitioner in refine: DOLFINx #3444
- Unify refine interface DOLFINx #3322
- Various improvements: DOLFINx: 3447, DOLFINx #3360, DOLFINx 3151
Standalone function for non-matching interpolation
Main contributor: Jørgen S. Dokken
With the introduction of sub-meshes, it became clear that having a single entry-point to interpolate
for both interpolation from a sub-mesh with a known map, and a non-matching mesh was not possible.
To fix this, we have decided to create a stand-alone interpolation for non-matching meshes
dolfinx.fem.Function.interpolate_nonmatching
.
Relevant pull requests since last release:
- Rewrite interpolation API for nonmatching meshes: DOLFINx #3177 for details.
Remove dolfinx.fem.set_bc
Main contributor: Garth N. Wells
The underlying code for setting Dirichlet boundary conditions have been simplified, and what was previously done with
dolfinx.fem.set_bc(b.array, bcs)
is now done with
[bc.set(b.array) for bc in bcs]
Note that dolfinx.fem.petsc.set_bc
has not changed.
Relevant pull requests since last release:
Rename scale
in apply_lifting
to alpha
Main contributor: Garth N. Wells
The name scale was misleading, as it only relates to the lifting part of the operation, i.e.
b -= alpha * A (g - x)
Relevant pull requests since last release:
PETSc object memory management
Main Contributors: Umberto Villa, Jack Hale and Garth N. Wells
Removal of dolfinx.fem.Function.vector
To avoid extra MPI communicator duplication when re-using the underlying dolfinx.la.Vector
the PETSc.Vec
wrapper in DOLFINx now lives in dolfinx.fem.Function.x.petsc_vec
.
Note on memory management
A note on the memory management of all PETSc objects created with DOLFINx has been added to relevant functions.
Relevant pull requests
- Associate PETsc vector with
dolfinx.la.Vector
: DOLFINx #3092 - Note on memory management: DOLFINx #3329
- Remove deprecated method: DOLFINx #3435
Other updates
Upgrade from PETSc<3.21.6
A bug was introduced in PETSc 3.21 that means that newtontr
in PETSc.SNES
did not work in serial in DOLFINx: PETSc #1645.
This has been rectified in PETSc 3.21.6.
The new stable docker images uses PETsc 3.22: DOLFINx: 3445
Performance improvement for compilation of forms with multiple subdomains
A fundamental re-design in UFL has been added to reduce the amount of generated integration kernels when using sub-domains, see: UFL: #305 for details
Documentation improvements for dolfinx.mesh
The dolfinx.mesh.Topology
and dolfinx.mesh.Geometry
classes are now fully wrapped in Python and properly documented. This means that the user can inspect the signatures and input/output types of these functions. See: DOLFINx #3406, DOLFINx #3403, DOLFINx #3400
ARM VTK support
The ghcr.io/fenics/dolfinx/lab
images now includes VTK binaries to use with pyvista, thanks to Henrik Finsberg, DOLFINx #3259
Create mixed topology mesh
Initial support by Chris Richardson for mixed topology meshes, see: DOLFINx #3271 for details.
Relevant PRs: DOLFINx #3237, DOLFINx #3240, DOLFINx #3234, DOLFINx #3223
Support arbitrary length tensor in outputting
Previously only up to 3x3 tensors and length 3 vectors would be outputted. DOLFINx #3227 resolves this for VTKFile
, VTXWriter
and XDMFFile
.
Interpolation matrix with DOLFINx native matrices
Previously only supported with PETSc matrices. Now extended with DOLFINx #3226.