SIAM News Blog

The FEniCS Project

Towards Automated Scientific Computing

By Robert C. Kirby

The numerical approximation of differential equations to solve problems in science and engineering was a driving force in the development of early computers. Classical and modern techniques alike work to produce solutions by the repetition of vast amounts of arithmetic, but performing this arithmetic by human effort is a tedious and error-prone task. The world now takes for granted the revolution begun by Alan Turing, John von Neumann, and other pioneers of computer science; machines can carry out simple computations with speed and accuracy far surpassing even a team of humans, the original parallel computers [2].

However, programming these early computers in machine code proved to be tedious and error-prone as well. This led to automated programming by assemblers, and later by FORTRAN compilers. In the 1950s, a machine-independent language with high-level data structures, such as arrays, offered scientists new freedom. Ongoing improvements in programming languages and compilers, machine architecture, and algorithms led scientists to tackle bigger challenges that further stretched boundaries in all these areas.

Moving through several generations of progress to the present day, one may rightly say that implementing algorithms for inverse problems dealing with finite element discretizations of coupled systems of nonlinear partial differential equations (PDEs) in unstructured three-dimensional geometry is a tedious and error-prone task. Can this be helped by computers?

The resolution of our (recursive) conundrum begs the development of even higher-level descriptions of our tasks, along with tools for reducing these descriptions into lower-level code. Much like how FORTRAN compilers turn mathematical formulas and array syntax into assembly code, higher-level, domain-specific languages can describe higher-level mathematical structures. MATLAB is a quintessential domain-specific language, concealing implementation details of linear algebra behind a custom, operator-overloaded syntax. Similarly, one could “embed” a domain-specific language in existing high-level languages such as C++ or Python. For example, looking at expression template libraries such as Blitz++ for use of C++ facilities can help implement linear algebra syntax.

Several modern software projects seek to bridge the even greater distance between machines and PDEs via domain-specific languages. A leading example of this approach is the FEniCS project. For example, the “Hello World” program in FEniCS offers a complete numerical solution of the weak form of Poisson’s equation (see Figure 1).

Figure 1. Complete working Python script for solving Poisson’s equation with FEniCS.
Many things happen when this deceptively-simple Python program is run. The expression inner(grad(u), grad(v))*dx is a symbolic representation of the weak form of the Laplace operator. When the solve function is called, a special-purpose just-in-time compiler generates low-level code based on this expression, links it into the environment, builds the sparse stiffness matrix, and invokes a linear solver. The matrix format and solver library support many options, including PETSc and Trilinos. With a bit more code, it is possible for users to configure exactly how these libraries are utilized.

FEniCS was created in 2003 to achieve automated computational mathematical modeling [1]. The two original FEniCS components were DOLFIN, a C++ library originally written by Anders Logg and Johan Hoffman for describing meshes, solvers, and bilinear forms, and FIAT, a novel if somewhat esoteric Python package for automatically constructing high-order bases for Raviart-Thomas and other complicated finite element spaces. Ridgway Scott brought prior experience to the project with a Scheme-based system called Analysa [2].

FEniCS progressed rapidly over the next several years. The FEniCS Form Compiler (FFC) [3] enabled the generation of low-level code for element matrix and vector construction linkable against DOLFIN. Later, its input language was lifted out and extended to become the Unified Form Language (UFL). Garth Wells joined the project early on, expanding FEniCS’s boundaries by tackling new problems and building new capabilities as needed. When Simula Research Laboratory received Norwegian Center of Excellence status and built the Center for Biomedical Computing, FEniCS advanced to the forefront of the Center’s modeling efforts in problems related to blood flow, aneurisms, and cardiac mechanics and electrophysiology.

While FEniCS is now a global, visible project, we cannot overlook the work of similar projects. Even older than FEniCS is Sundance [4], a C++ library developed by Kevin Long that provides a high-level description of variational forms and their symbolic derivatives. Rather than generating customized code, it relies on a highly-optimized interpreter for bilinear forms. In contrast to these embedded languages, the FreeFEM++ project of comparable vintage provides an entire domain-specific language and user environment. This offers a fairly complete user environment (more like MATLAB), but requires the recreation of many language features, making interaction with external codes more challenging. Deal.II, a more traditional C++ library from the 1990s that won the 2007 Wilkinson Prize, delivers a rich set of finite element classes with special attention to adaptive methods but lacks an automation layer present in these other projects. The Deal.II community has produced an incredible set of documentation and tutorials for learning finite element methods through the software’s many features.

Also worth noting is a new generation of projects inspired by the vision and success of those before them. Employing FEniCS as a point of departure, Firedrake attempts to provide compatible and enhanced functionality, engineered accordingly to different design principles. A notable result of the Firedrake/FEniCS collaboration is the dolfin-adjoint project, featuring lead developers of both software. Using UFL, the team, which won the 2015 Wilkinson prize, is able to automatically derive and discretize adjoints and other operations required for inverse problems.

It is an exciting time for numerical software. Much as early computers and FORTRAN lifted us above the slide rule and assembly code respectively, projects like FEniCS allow us to think of variational forms rather than arrays. Enduring challenges such as uncertainty quantification, accelerator-based architectures, and even more challenging applications will continue to spur future developments. The FEniCS project welcomes new users and developers and invites users, developers, colleagues, competitors, and other interested parties to its annual meeting.

[1] Dupont, T., Hoffman, J., Johnson, C., Kirby, R.C., Larson, M., Logg, A., & Scott, L.R. (2003). The FEniCS Project. Preprint 2003-21. Chalmers Finite Element Center: Chalmers University of Technology.

[2] Grier, D.A. (2013). When Computers Were Human. Princeton, NJ: Princeton University Press.

[3] Kirby, R.C., & Logg, A. (2006). A compiler for variational forms. ACM Transactions on Mathematical Software, 32(3), 417-444.

[4] Long, K.R. (2003). Sundance Rapid Prototyping Tool for Parallel PDE Optimization. In L.T. Biegler (Ed.), Large-Scale PDE-Constrained Optimization (pp. 331-341). Berlin: Springer-Verlag.

Robert C. Kirby earned his Ph.D. in computational and applied mathematics from the University of Texas at Austin. He is currently on the faculty of the Department of Mathematics at Baylor University.

Randall J. LeVeque ([email protected]) of the University of Washington, Seattle, is the editor of the Software and Programming column.

blog comments powered by Disqus