SIAM News Blog

FreeFem++: A High Level MultiPhysics Finite Element Software

By Olivier Pironneau

The finite element method (FEM) was invented shortly after computers as a natural framework for solid mechanics; the success of NASTRAN—a finite element analysis program—in the 1960s is well known. Variational methods were then popular among mathematicians in the analysis of partial differential equations (PDEs), which rendered FEM easily adaptable to other domains of physics, such as thermodynamics, fluid mechanics, and electromagnetism. However, the preparation of data at the time was a nightmare, and still is for many 3D applications. Thus, teaching FEM for PDEs was a challenge: one could easily spend half the course on triangulation and graphics methods.

Niklaus Wirth and the Finite Element Method

Dominique Bernardi and I wrote FreeFem+ in the 1990s at the University of Paris VI to ease the teaching and prototyping of PDE algorithms. Upon Frédéric Hecht’s introduction of the advanced automatic mesh generator, I realized the power of that module and asked him to partake in the FreeFem++ project. We prioritized the user interface’s ability to define the PDE and hide the numerical methods. We aimed for the highest possible level, nearest to the mathematical statements. For instance, a Stokes system in

\[-\Delta\mathbf{u} + \nabla p=0, \quad \nabla\cdot\mathbf{u}=0, \quad \int_\Omega p=0, \quad \mathbf{u} \enspace \textrm {given on the boundary},\]

discretized with the Hood-Taylor element and regularization of the pressure constraint, would be defined and solved by the code in Figure 1. The only non-intuitive statement is “on(3,u=x*(1-x),v=0),” referring as it does to the numbering convention of the four sides of the square, 3 being the top side. The last instruction generates Figure 2. 

Figure 1. FreeFem++ code for the cavity-driven Stokes problem in 2D.

For nonlinear, time-dependent multiphysics systems, it is the user’s responsibility to specify his/her scheme. For instance, Navier-Stokes equations in a 3D cube—up to time 5, discretized with an implicit Euler scheme, and with semi-linearization of the convection terms—would need a problem block—which only defines the PDE—to replace the solve block, which defines and solves the PDE, for reusability in a time loop (see Figure 3a). This yields 3b. Notice the use of macro for clarity and the approach’s generality. For instance, adding a temperature equation to simulate natural convection would cost only a few extra lines of code.

Figure 2. Freefem++ code for the Navier-Stokes equations in 3D.

Parsing regular expression, such as x*(1-x), requires techniques that are at the core of compilation. Algorithms and Data Structures by Niklaus Wirth explains it well [2].

The syntax of FreeFem++ scripts obeys the rule of LL(1) computer languages. It can be parsed in one pass, generating and interpreting a stack of instructions at execution time, much like Pascal. It does not generate bytecode, meaning it is not possible to embed this Navier-Stokes solver in another application written in MATLAB, for instance, without including the whole FreeFem++ environment. Incorporating MATLAB into FreeFem++ may be easier.

Figure 3. Results of Figures 1 and 2. 3a. Stokes flow: velocity vectors and pressure lines. 3b. Navier-Stokes flow: iso-pressure surfaces. Image credit: Olivier Pironneau.

The Advantages of C++ to FreeFem++

As the years passed, more and more instructions were added to FreeFem++, which now has a considerable subset of C++ instructions in its syntax (useful to open and close files, manipulate matrices, etc.). Actually, the name FreeFem++ comes from the fact that the software was entirely rewritten in C++ in 2002 and again in 2010. It makes extensive use of operator overloading and templates, something known as generic programming. Consequently, including complex numbers and adding 3D problems—vector problems like Stokes’—was easy, since the syntax is the same as for the heat equation.

FreeFem++ is a collaborative open-source effort like its cousin FEniCS,1 (Ridgeway Scott, the latter’s co-creator, is a friend), but Hecht does 95 percent of the development, and I am still wondering how one man can do so much. I’d be tempted to say “try both and tell us which is best,” but that is difficult; these are easy to try but time-consuming to master. Both software packages are huge, with FreeFem++ incorporating 130K lines of code. Most well-known libraries of numerical analysis are integrated and can be called from within FreeFem++, including blas, mpi, Linpack, Arpack, GSL, PETSc, MUMPS, UFMPACK, several optimization modules like CMAES and IPOPT, and the 3D mesh generator TetGen. The documentation contains 418 pages that nobody wants to read, and you truly don’t have to for simple cases.

Strengths and Weaknesses

By comparison to its commercial competitor COMSOL, FreeFem++ is primitive for graphics but interfaced with paraview and others. Nonetheless, it does things that COMSOL cannot do. Both run on Mac OS, Windows, and Linux in integrated environments; FreeFem++ runs on smartphones thanks to Antoine Le Hyaric. However, we haven’t found much use for this extension yet.

Figure 4. Complex 3D examples. 4a. Iso surfaces of pressure in aortic blood flow. Image credit: Olivier Pironneau. More details in [1]. 4b. Navier-Stokes flow: iso-pressure surfaces. Blood flow reconstruction in the brain by simulation of Maxwell equations. Top: Target permitivity of the brain generated by high-precision imaging. Bottom: Computer-reconstructed image from noisy microwave measurements. The arrows indicate the region of stroke. Image credit: Frédéric Nataf.
FreeFem++ has the best open-source module for mesh adaptivity. One can also easily interpolate functions on different meshes, simplifying domain decomposition methods, for example. It has a robust implementation of the characteristic Galerkin method for convection terms, and a large library of triangular finite elements—P0. . .P4, RT, HCT, Edge, mixed—adapted to discontinuous Galerkin methods (DG) and a posteriori mesh refinements. The one thing FreeFem++ does not handle well is conservation laws with shocks. We have not yet been able to write a fast Riemann solver for a general hyperbolic system of any size; one may solve such PDEs with DG but doing so requires considerable expertise, which is not in the spirit of FreeFem++.

Many C++/Fortran90 toolboxes are potentially just as powerful, such as feel++, deal.II, and Elmer, among others. But in general, these require that users have a high level of programming ability. Field-specific languages are coming; the meta language scala makes such a promise.

FreeFem++ has at least 3,000 active users. Ph.D. students worldwide seem to love it; in France, almost all courses on numerical analysis for PDEs use FreeFem++. But it is also increasingly popular as a research tool because execution time is comparable with handwritten C++ codes, and development time is so much faster.

For instance, Figure 4a shows the iso-pressure surfaces of blood flow in an aorta using a mesh constructed from an MRI. A simulation of Maxwell equations in 3D (see Figure 4b), entirely written with FreeFem++, won the French supercomputing Atos-Fourier prize because it scales perfectly on a parallel machine up to 10,000 cores. The next challenge is hiding parallel computing instructions from the user and automatically employing the computer’s resource. Frédéric Nataf and Pierre Jolivet are working hard at it, with Hecht’s help, of course!

1 Read a previous SIAM News article by Robert C. Kirby about FEniCS.

[1] Rebollo, T.C., Girault, V., Murat, F., & Pironneau, O. (2016). Analysis of a Coupled Fluid-Structure Model with Applications to Hemodynamics. SIAM J. Numer. Anal., 54(2), 994-1019.
[2] Wirth, N. (1985). Algorithms and Data Structures. Upper Saddle River, NJ: Prentice Hall.

Olivier Pironneau is emeritus professor at the Université Pierre et Marie Curie in Paris, France. He is also a member of the French Académie des Sciences and co-author of eight books, including Computational Methods for Option Pricing (SIAM, 2005). Visit his website for more information.

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