Release notes#

v0.11.0#

Since the 0.10.0 release, there has been 177 merged pull requests from 27 contributors. Below follows a summary of the biggest changes to the Python API. A full diff can be found here.

SuperLU_DIST interface#

Authors: Jack Hale and Chris Richardson

After the support for Windows was added in v0.9.0, there has been a lack of parallel supported solvers to solve the resulting linear problems. Furthermore, PETSc sometimes can feel complicated for simple problems presented during teaching. In this release, we add support for using SuperLU_DIST with native sparse matrices dolfinx.la.MatrixCSR.

The following constructors and classes have been added:

Built-in matrix support#

Authors: Chris Richardson

The ‘real’ element#

Authors: Jørgen S. Dokken and Matthew Scroggs

The real-element has for a long time not been present in DOLFINx. This has mainly been due to the fact that the implementation in legacy FEniCS introduced a ton of special casing within core functionality, which we thought was better to avoid in DOLFINx. However, a prototype implementation of the real-element has been around for a few releases, and has now been implemented in the core libraries. Users can now call

import basix.ufl
import dolfinx.fem

el = basix.ufl.real_element(mesh.basix_cell(), dtype=dtype, shape=(N, ))
R = dolfinx.fem.functionspace(mesh, el)

to create a function space consisting of N values (of data type dtype, which can be a complex type).

Furthermore, to use this space alongside other spaces, for instance for Lagrange multipliers, users are recommended to use ufl.MixedFunctionSpace(V, R, ...) and ufl.extract_blocks() to create blocked systems that can be used in dolfinx.fem.petsc.LinearProblem or dolfinx.fem.petsc.NonlinearProblem.

Threading#

Authors: Chris Richardson, Jørgen S. Dokken and Garth N. Wells

For a long time, DOLFINx has been exclusively using MPI for the distribution of computational load. However, with the computational landscape evolving to more and more heterogeneous systems, the need for additional parallelisation methods are required. In this release, we introduce initial threading support using std::jthread` in the following methods:

New and improved demos#

Authors: Jørgen S. Dokken, Paul T. Kühner

A new demo showcasing how to use PETSc and matrix-free solvers can be found in demo_matrix-free-petsc

The PML demo now shows how to use one-sided interior facet integrals with manual specification of integration entities.

The Biharmonic demo has gone through a major revision, using:

  • A more suitable choice of finite elements (as P2 yields sub-optimal convergence)

  • Better choice of penalty parameter

  • Change of boundary conditions from simply supported to clamped and explaining the effect of different BCs.

  • Verify solution with the method of manufactured solutions and add relevant references

In general, demos now consistently use tdim=mesh.topology.dim, fdim = tdim -1 and gdim = mesh.geometry.dim to avoid confusion for new users.

Further improvements in submesh support#

Authors: Jørgen S. Dokken

A feature that for a long time has existed outside of the FEniCS core is the dolfinx.mesh.transfer_meshtags_to_submesh(), which makes it possible to transfer a meshtag from a parent mesh to a submesh. This function is now part of the core library.

Another feature introduced in this release is the usage of submesh quantities in dolfinx.fem.Expression. You can now pass Expressions containing coefficients and constants from a submesh, combined with geometric quantities, coefficients, and constants of the parent mesh. For example

V = fem.functionspace(mesh, el)
u = fem.Function(V, dtype=dtype)
# Populate `u` ...

mesh.topology.create_connectivity(mesh.topology.dim - 1, mesh.topology.dim)
exterior_facets = exterior_facet_indices(mesh.topology)
submesh, entity_map, _, _ = create_submesh(mesh, mesh.topology.dim - 1, exterior_facets)
u_sub = fem.Function(fem.functionspace(submesh, sub_el), dtype=dtype)
# Populate `u_sub` ...

quadrature_points, _ = basix.make_quadrature(basix.CellType.interval, qdegree)
quadrature_points = quadrature_points.astype(xtype)

n_h = ufl.FacetNormal(mesh)
f = u * n + u_sub * n

mesh.topology.create_connectivity(mesh.topology.dim - 1, mesh.topology.dim)
compiled_expr = fem.Expression(expr, quadrature_points, dtype=dtype, entity_maps=[entity_map])

Extending GMSH and VTKHDF IO#

Authors: Jørgen S. Dokken

With the changes to mesh in v0.10.0 max_facet_to_cell_links was introduced to make it possible to create meshes with joints, branches, etc. This is now exposed in dolfinx.io.gmsh.model_to_mesh().

Furthermore, new cell types are supported for the vtkhdf backend, including all linear and quadratic VTK cell types.

The GMSH interface can now read in meshes on mixed-topology grids using dolfinx.io.gmsh.model_to_mesh() or dolfinx.io.gmsh.read_from_msh(). There is currently no support for reading entity tags on mixed-topology grids through this interface yet.

Exposing tolerances for non-affine pull-backs#

Authors: Jørgen S. Dokken

A set of users have had issues with non-affine geometries, in particular higher order grids, and the tolerance and maximum number of iterations in the dolfinx.fem.CoordinateElement.pull_back(), dolfinx.fem.Function.interpolate_nonmatching() and dolfinx.fem.Function.eval() yielding errors such as:

RuntimeError: Newton method failed to converge for non-affine geometry

If you see this error message, try to increase the maxiter or tol.

Mixed topology meshes, prisms and pyramid cells#

Authors: Jørgen S. Dokken and Chris Richardson

For mixed topology meshes, there are many notions that do not align with the original design of DOLFINx. Examples are the members num_entity_dofs and num_entity_closure_dofs of dolfinx.cpp.fem.ElementDofLayout, as prisms and pyramids do not have the same number of dofs per sub-entity. They have been removed, and users should instead call len(ElementDofLayout.entity_dofs(dim, entity_index))

Furthermore, dolfinx.fem.apply_lifting() and dolfinx.fem.assemble_scalar() now work for mixed-topology meshes.

Meshes#

Authors: Jørgen S. Dokken and Jack Hale

  • New function dolfinx.mesh.create_point_mesh() to create a point cloud mesh with no points shared between the different processes. Useful for reading in point measures or outputting data.

  • Built in mesh-generators such as dolfinx.mesh.create_rectangle() now takes an optional argument gdim that embeds in a larger space. This simplifies the testing process for problems on manifolds.

  • New function dolfinx.fem.interpolate_geometry() that allows users to create a new mesh either raising or lowering the polynomial order of the mesh geometry. The topology is shared between the old and new grid. It is also possible to switch the Lagrange-variant of the underlying coordinate element.

Interpolation#

Authors: Jørgen S. Dokken and Garth N. Wells

A crucial bug interpolating Piola-mapped elements from parent to a codim-0 submesh has been fixed for this release.

Modernizing UFL#

Authors: Paul T. Kühner

UFL is a Python project that has been in development for almost 20 years, and Python has gone through a massive modernization during this time. One of the visually pleasing improvements is the use of @property-decorators. ufl.AbstractCell now uses properties for topological_dimension and cellname, etc. while ufl.Mesh now has geometric_dimension. See UFL PR #385 for more details.

v0.10.0#

Since the 0.9.0 release, there have been 311 merged pull requests from 25 contributors. Below follows a summary of the biggest changes to the Python-API from these pull requests. In addition to the changes below, the ever-lasting quest of improving performance and squashing bugs continues.

PETSc API#

Authors: Jørgen S. Dokken, Francesco Ballarin, Jack Hale and Garth N. Wells

Mapping data between PETSc.Vec and dolfinx.fem.Functions is now trivial for blocked problems by using dolfinx.fem.petsc.assign().

Both solvers and assembly routines interfacing with PETSc has received a drastic make-over to improve usability and maintenance, both for developers and end-users.

Improved non-linear (Newton) solver#

The FEniCS project has for the last 15 years had its own implementation of a Newton solver. We no longer see the need of providing this solver, as the PETSc SNES solver, and equivalent solver for C++ provides more features than our own implementation.

The previously shipped dolfinx.nls.petsc.NewtonSolver is deprecated, in favor of dolfinx.fem.petsc.NonlinearProblem, which now integrates directly with petsc4py.PETSc.SNES.

The non-linear problem object that was sent into dolfinx.nls.petsc.NewtonSolver has been renamed to NewtonSolverNonlinearProblem and is also deprecated.

The new NonlinearProblem has additional support for blocked systems, such as NEST by supplying kind="nest" to its initializer. See the documentation for further information.

Improved dolfinx.fem.petsc.LinearProblem#

  • The dolfinx.fem.petsc.LinearProblem now supports blocked problems, either specified manually or by using ufl.extract_blocks(). By changing the input-argument kind, the user can now decide if they want to use the DOLFINx blocked PETSc implementation (kind="mpi") or the kind="nest".

  • The default behavior for non-blocked systems remains the same as before.

  • The users can now also specify a (blocked) form for preconditioning through the P keyword argument in the constructor.

Assembly routines#

In earlier versions of DOLFINx, there were three assembly routines for PETSc.Vec and PETSc.Mat:

  • assemble_*

  • assemble_*_block

  • assemble_*_nest

This caused a lot of duplicate logic in codes. Therefore, we have unified all these assembly routines under dolfinx.fem.petsc.assemble_vector() and dolfinx.fem.petsc.assemble_matrix(). The input keyword argument kind selects the relevant assembler routine. See for instance the Stokes demo for a detailed introduction. Similar changes has been done to dolfinx.fem.petsc.apply_lifting().

Linear algebra submodule#

There is now a sub-module (dolfinx.la.petsc) containing PETSc LA operations.

Interpolation#

The dolfinx.fem.discrete_curl() operator has been added to DOLFINx, to cater to Hypre Auxiliary-space Divergence Solver

Simplified demos#

Usage of ufl.MixedFunctionSpace and ufl.extract_blocks()#

Authors: Jørgen S. Dokken and Joe Dean

Initially introduced as part of the v0.9.0-release, usage of these two UFL-abstractions has been propagated into the demos, to make it even easier for users to see examples of how to work with blocked problems.

TODO: Add profiling of blocked/mixed-element vs mixedfunction-space.

Usage of ufl.ZeroBaseForm#

Author: Garth N. Wells For a long time, it has not been possible to specify the right hand side of a linear PDE as empty. This means that users often have had to resolve to adding dolfinx.fem.Constant(mesh, 0.0)*v*ufl.dx to ensure that one can use the dolfinx form compilation functions. With the introduction of ufl.ZeroBaseForm this is no longer required. The aforementioned workaround can now be reduced to ufl.ZeroBaseForm((v, )), which avoid extra assembly calls within DOLFINx.

Form compiler and integral types#

Author: Susanne Claus, Paul T. Kühner, and Jørgen S. Dokken

Mesh#

Authors: Paul T. Kühner, Joe Dean, Garth N. Wells, Jørgen S. Dokken and Chris Richardson

  • Uniform mesh refinement of all CellTypes is available through dolfinx.mesh.uniform_refine().

  • Branching meshes (a mesh where a single facet is connected to more than two cells), such as T-joints (3 cells connected to a single facet) are now supported as input meshes to DOLFINx. To ensure proper partitioning in parallel, one should change the default option max_facet_to_cell_links to how many cells a facet can be attached to in dolfinx.io.XDMFFile.read_mesh(), dolfinx.io.vtkhd.read_mesh() and dolfinx.mesh.create_mesh().

  • One can no longer use set_connectivity or set_index_map to modify dolfinx.mesh.Topology objects. Any connectivity that is not (tdim, 0), (tdim, tdim) or (0, 0) should be created with dolfinx.mesh.Topology.create_connectivity(). The aforementioned connections should be attached to the topology when calling dolfinx.cpp.mesh.create_topology().

  • Mixed-dimensional support has been vastly improved by creating dolfinx.mesh.EntityMap, which replaces the numpy arrays used as entity_maps in dolfinx.fem.form() in the previous release. This is a two-way map, meaning that the user no longer has to take care of creating the correct mapping. The two-way map from a sub-mesh to a parent mesh is returned as part of dolfinx.mesh.create_submesh().

Linear Algebra#

Authors: Chris Richardson

The native matrix-format now has a sparse matrix-vector multiplication dolfinx.la.MatrixCSR.mult(). Note that the dolfinx.la.Vector that you multiply with should use the dolfinx.la.MatrixCSR.index_map(1) rather than the one stemming from the dolfinx.fem.FunctionSpace.dofmap.index_map.

Collision detection#

Author: Chris Richardson The collision detection algorithm dolfinx.geometry.compute_distance_gjk() now used multiprecision to ensure proper collision detection. The algorithm has also been improve to work on co-planar convex hulls.

Documentation#

Authors: Paul T. Kühner, Garth N. Wells, Mehdi Slimani and Jørgen S. Dokken

Revised timer logic#

Authors: Garth N. Wells and Paul T. Kühner

Instead of using the Boost timer library, we have opted for the standard timing library std::chrono. The switch is mainly due to some observed inaccuracies in timings with Boost. This removes the notion of wall, system and user time. See or Timer for examples of usage.

IO#

GMSH#

Authors: Paul T. Kühner, Jørgen S. Dokken, Henrik N.T. Finsberg and Pierric Mora

The GMSH interface to DOLFINx has received a major upgrade.

  • An API-breaking change is that the module dolfinx.io.gmshio has been renamed to dolfinx.io.gmsh.

  • Another API-breaking change is the return type of dolfinx.io.gmshio.model_to_mesh() and dolfinx.io.read_from_msh(). Instead of returning the dolfinx.mesh.Mesh, cell and facet dolfinx.mesh.MeshTags, it now returns a dolfinx.io.gmsh.MeshData data-class, that can contain dolfinx.mesh.MeshTags of an sub-entity:

    • Cell (codim 0)

    • Facet (codim 1)

    • Ridge (codim 2)

    • Peak (codim 3)

  • Additional checks and error handing for Physical tags from GMSH has been added to improve the user experience.

VTKHDF5#

Authors: Chris Richardson and Jørgen S. Dokken

As Kitware has stated that VTKHDF is the future format they want to support, we have started the transition to this format. Currently, the following features have been implemented:

VTXWriter#

Authors: Mehdi Slimani and Jørgen S. Dokken

The writer does now support time-dependent DG-0 data, which can be written in the same file as a set of functions from another (unique) function space.

XDMF#

Author: Massimiliano Leoni and Paul T. Kühner

  • When using dolfinx.io.XDMFFIle.read_meshtags() one can now specify the attribute name, if the grid has multiple tags assigned to it.

  • Flushing data to file is now possible with dolfinx.io.XDMFFile.flush(). This is useful when wanting to visualize long-running jobs in Paraview.

Remove Fides backend#

As we unfortunately haven’t seen an expanding set of features for the Fides Reader in Paraview, we have decided to remove it from DOLFINx.

Pyvista#

Pyvista no longer requires pyvista.start_xvfb() if one has installed vtk with OSMesa support.