Many classical and modern finite element spaces are derived by dividing each computational cell into finer pieces. Such macroelements frequently enable the enforcement of mathematically desirable properties such as divergence-free conditions or \(C^{1}\) continuity in a simpler or more efficient manner than elements without the subdivision. Although a few modern software projects provide one-off support for particular macroelements, a general approach facilitating broad-based support has, until now, been lacking. In this work, we describe a major addition to the FIAT project to support a wide range of different macroelements. These enhancements have been integrated into the Firedrake code stack. We provide numerical evaluation of the new macroelement facility.

FIAT (the FInite element Automatic Tabulator) provides a powerful Python library for the generation and evaluation of finite element basis functions on a reference element. This release paper describes recent improvements to FIAT aimed at improving its run time and the accuracy and efficiency of code generated using FIAT-provided information. In the first category, we have greatly streamlined the implementation of orthogonal polynomials out of which finite element bases are built. The second category comprises several more advances. For one, we have built an interface to the recursivenodes package to enable more accurate Lagrange bases at high order. We have also implemented integral-type degrees of freedom for \(H(\operatorname{div})\) and \(H(\operatorname{curl})\) elements, which match the mathematical definitions of the elements more closely and also avoid loss of accuracy in interpolation. More fundamentally, we have included families of simplicial quadrature rules that require many fewer quadrature points than the Stroud rules previously used in FIAT. Finally, FIAT now provides support for fast diagonalization methods, which enable fast solution algorithms at very high order. In each case, we describe the new features in FIAT and illustrate some of the gains obtained through simple numerical tests.

Irksome is a library based on the Unified Form Language (UFL) that enables automated generation of Runge–Kutta methods for time-stepping finite element spatial discretizations of Partial Differential Equations (PDEs). Allowing users to express semidiscrete forms of PDEs, it generates UFL representations for the stage-coupled variational problems to be solved at each timestep. The Firedrake package then generates efficient code for evaluating these variational problems and allows users a wide range of options to deploy efficient algebraic solvers in PETSc. In this article, we describe several recent advances in Irksome . These include alternate formulations of the Runge–Kutta time-stepping methods and optimized support for diagonally implicit (DIRK) methods. Additionally, we present new and improved tools for building preconditioners for the resulting linear and linearized systems, demonstrating that these can lead to efficient approaches for solving fully implicit Runge–Kutta discretizations. The new features are demonstrated through a sequence of computational examples demonstrating the high-level interface and obtained solver performance.

We present an implementation of the trimmed serendipity finite element family, using the open-source finite element package Firedrake. The new elements can be used seamlessly within the software suite for problems requiring H 1 , H (curl), or H (div)-conforming elements on meshes of squares or cubes. To test how well trimmed serendipity elements perform in comparison to traditional tensor product elements, we perform a sequence of numerical experiments including the primal Poisson, mixed Poisson, and Maxwell cavity eigenvalue problems. Overall, we find that the trimmed serendipity elements converge, as expected, at the same rate as the respective tensor product elements, while being able to offer significant savings in the time or memory required to solve certain problems.

While implicit Runge–Kutta (RK) methods possess high order accuracy and important stability properties, implementation difficulties and the high expense of solving the coupled algebraic system at each time step are frequently cited as impediments. We present Irksome , a high-level library for manipulating UFL (Unified Form Language) expressions of semidiscrete variational forms to obtain UFL expressions for the coupled Runge–Kutta stage equations at each time step. Irksome works with the Firedrake package to enable the efficient solution of the resulting coupled algebraic systems. Numerical examples confirm the efficacy of the software and our solver techniques for various problems.

Many classical finite elements such as the Argyris and Bell elements have long been absent from high-level PDE software. Building on recent theoretical work, we describe how to implement very general finite-element transformations in FInAT and hence into the Firedrake finite-element system. Numerical results evaluate the new elements, comparing them to existing methods for classical problems. For a second-order model problem, we find that new elements give smooth solutions at a mild increase in cost over standard Lagrange elements. For fourth-order problems, however, the newly enabled methods significantly outperform interior penalty formulations. We also give some advanced use cases, solving the nonlinear Cahn-Hilliard equation and some biharmonic eigenvalue problems (including Chladni plates) using C 1 discretizations.

We revisit the question of optimizing the construction and application of finite element matrices. By using commuting properties of the reference mappings and duality, we reorganize stiffness matrix construction and matrix-free application so that the bulk of the work can be done by optimized matrix multiplication libraries. We provide examples, including numerical experiments, with the Laplace and curl-curl operators as well as develop a general framework. Our techniques are applicable in general geometry and are not restricted to constant coefficient operators.

The L 2 -orthogonal polynomials used in finite and spectral element methods on nonrectangular elements may be defined in terms of collapsed coordinates, wherein the shapes are mapped to a square or cube by means of a singular change of variables. The orthogonal basis is a product of specific Jacobi polynomials in these new coordinates. Implementations of these polynomials require special handling of the coordinate singularities. We derive new recurrence relations for these polynomials on triangles and tetrahedra that work directly in the original coordinates. These relations, also applicable to pyramids and prisms, do not require any special treatment of singular points. These recurrences are seen to speed up both symbolic and numerical computation of the orthogonal polynomials.

We examine the effect of using complexity-reducing relations [Kirby et al. 2006] to generate optimized code for the evaluation of finite-element variational forms. The optimizations are implemented in a prototype code named FErari, which has been integrated as an optimizing backend to the FEniCS form compiler, FFC [Kirby and Logg 2006; 2007]. In some cases, FErari provides very little speedup, while in other cases we obtain reduced local operation counts by a factor of as much as 7.9 and speedups for the assembly of the global sparse matrix by as much as a factor of 2.8 (see Figure 9).

We investigate the compilation of general multilinear variational forms over affines simplices and prove a representation theorem for the representation of the element tensor (element stiffness matrix) as the contraction of a constant reference tensor and a geometry tensor that accounts for geometry and variable coefficients. Based on this representation theorem, we design an algorithm for efficient pretabulation of the reference tensor. The new algorithm has been implemented in the FEniCS Form Compiler (FFC) and improves on a previous loop-based implementation by several orders of magnitude, thus shortening compile-times and development cycles for users of FFC.

Our previous work on FIAT (Finite Element Automatic Tabulator) developed a “computational representation theory ” that allowed us to construct arbitrary order instances of a wide range of finite elements, many of which are infrequently used owing to their associated code complexity. In our present work, we further hone this theory by rephrasing most of the internal operations as linear transformations over finite-dimensional Banach spaces. This additional insight has led to increased code granularity and allowed the use of level 3 BLAS operations. This is both a conceptual and a practical development; as the run-time performance of FIAT has been improved multiple orders of magnitude.

As a key step towards a complete automation of the finite element method, we present a new algorithm for automatic and efficient evaluation of multilinear variational forms. The algorithm has been implemented in the form of a compiler, the FEniCS Form Compiler (FFC). We present benchmark results for a series of standard variational forms, including the incompressible Navier--Stokes equations and linear elasticity. The speedup compared to the standard quadrature-based approach is impressive; in some cases the speedup is as large as a factor of 1000.

Much of finite element computation is constrained by the difficulty of evaluating high-order nodal basis functions. While most codes rely on explicit formulae for these basis functions, we present a new approach that allows us to construct a general class of finite element basis functions from orthonormal polynomials and evaluate and differentiate them at any points. This approach relies on fundamental ideas from linear algebra and is implemented in Python using several object-oriented and functional programming techniques.