CFD Direct

The Architects of OpenFOAM

Modular Solvers in OpenFOAM

In August 2022, CFD Direct introduced modular solvers to the development version of OpenFOAM (OpenFOAM-dev). Modular solvers are written as classes, in contrast to the traditional application solvers, integral to OpenFOAM since icoFoam in 1993.  They are simpler to use and maintain than application solvers. Their source code is easier to navigate, promoting better understanding. They are more flexible; in particular, modules for different fluids and solids can be coupled across multiple regions, notably for conjugate heat transfer (CHT) with any type of flow, e.g. multiphase.  Modular solvers are deployed using the foamRun or foamMultiRun applications, which contain a generic solution algorithm for single and multiple regions, respectively. Additional modules and applications replace existing tools for data processing and case configuration.

Application Solvers

OpenFOAM was originally conceived as a general development platform for computational fluid dynamics (CFD). Its creator, Henry Weller, exploited object orientation and generic programming in C++ to manage the complex science/engineering, mathematics and computing encountered in CFD.  These programming paradigms minimise code duplication for more efficient maintenance. Data and functions are packaged in classes, with the functions providing a convenient interface to manipulate the data.

Most of the code in OpenFOAM is object-oriented and generic, i.e. contained within C++ classes and templates. It is divided across numerous libraries (~120 in OpenFOAM-v10), from the large core OpenFOAM and finiteVolume libraries to small libraries of physical models. By contrast, the application solvers contained a series of instructions that solved the governing equations within an iterative solution algorithm, and made function calls to various physical models, e.g. turbulence. The code was procedural, not object-oriented. 

The Rise and Fall of Application Solvers

Over time, individual application solvers were customised for different categories of flow, e.g. incompressible, compressible, multiphase, etc.  There were baseline solvers within each category, often with separate solvers for steady and transient flow.  There were also variants of the baseline solvers which included combinations of: 1) additional modelling, e.g. porous media, dynamic meshes, Boussinesq’s approximation for buoyancy, single and multiple rotating frames (SRF and MRF), etc.; 2) numerics, e.g. local time stepping, different implementations of equation, the “consistent” option in pressure-velocity coupling, etc.

Between the release of OpenFOAM v1.0 in 2004 and  OpenFOAM v2.1 in 2011, the number of application solvers doubled from 41 to 82. While this established OpenFOAM across a broad range of CFD, the proliferation of solvers created its own problems.  Users often struggled to choose between the different solvers. They then were forced to switch between solvers as they introduced new features into their simulations.  They  faced the daunting task of writing their own solver when none existed that included the particular combination of options they required. Code duplication was widespread, creating inconsistencies when a change to one application was not propagated through duplicate code in others.

To overcome these problems, CFD Direct generalised the optional modelling and numerics to move it out of the solver code. The number of solvers then decreased — along with the number of outstanding bugs — with users able to select the combination of options they wanted from their case configuration files. Some large components of OpenFOAM were significantly redesigned in the process, e.g. thermophysical transport, general models and constraints, and mesh motion, topology change and redistribution.

General Algorithm

When OpenFOAM v10 was released in 2022, the number of solvers has fallen to 48, close to the level in 2004.  However, v10 solvers are very different to those from v1.0, since they contain modelling and numerics that did not exist in v1.0. They are much more complex, making them more difficult to maintain, and still including a lot of code duplication. It was natural, therefore, to turn to object orientation as the tried and tested approach of managing complexity and reducing code duplication.

The concept was to have one application, foamRun, to describe a general solution algorithm with functions that could be defined across a collection of modular solvers. Development of the idea was motivated by a need for a general method to couple solutions in multiple regions, containing different types of fluid and/or solid. The foamRun application for a single mesh region could be extended to a foamMultiRun application which executed the same functions across multiple mesh regions. Both applications contain the same algorithm with a time loop whose key functions are shown below (from foamRun.C)

while (pimple.run(runTime))
{
    solver.preSolve();
    adjustDeltaT(runTime, solver);
    runTime++;

    while (pimple.loop())
    {
        solver.moveMesh();
        solver.fvModels().correct();
        solver.prePredictor();
        solver.momentumPredictor();
        solver.thermophysicalPredictor();
        solver.pressureCorrector();
        solver.postCorrector();
    }

    solver.postSolve();
    runTime.write();
}

The functions in the algorithm are principally responsible for the following actions.

  • preSolve(): changes the mesh topology and maps between meshes.
  • fvModels().correct(): updates the fvModels.
  • moveMesh(): moves the mesh.
  • prePredictor(): updates the fluid properties.
  • momentumPredictor(): solves the momentum equation.
  • thermophysicalPredictor(): solves an energy equation and updates thermodynamics.
  • pressureCorrector(): solves a pressure equation.
  • postCorrector(): solves the momentum transport, e.g. turbulence modelling, and thermophysical transport.
  • postSolve(): cleans up, e.g. clearing temporary data.

Modular Solvers

The solver described in the code above forms the base class for a hierarchy of modular solvers shown below. Solver classes are derived from the base class, e.g. solid, fluidSolver. As the hierarchy is traversed from left to right, functions within the general algorithm (which begin as virtual functions in the solver class) are defined for different descriptions of fluids, e.g. incompressible, isothermal, multicomponent and multiphase. Sometimes the function is empty, e.g. thermophysicalPredictor() in the incompressibleFluid solver, since it does not include temperature. A lot of functionality is inherited by a derived class from its base class; in some case, base classes are created to define a general category of fluid, e.g. VoFSolver, whose functionality is inherited by four derived solvers based on the volume of fluid (VoF) method.

At the time of writing [updated: 2 May 2023], there are 13 fluid solvers, 2 film solvers, 2 solid solvers and 2 other solvers. These modular solvers replace the application solvers as follows.

A modular solver can specified with foamRun through a solver keyword entry in the case controlDict file, e.g.

solver incompressibleFluid;

Alternatively, the solver can be specified using the -solver option when executing foamRun, i.e.

foamRun -solver incompressibleFluid

In addition to the modular solvers, two new applications have been created to replace functionality relating to application solvers.

  • foamPostProcess utility : replicates the postProcess utility and the solver -postProcess option with modular solvers.
  • foamToC utility : replaces the -list... options to locate models list their options.

Working with Modular Solvers

Through modular solvers, the benefits of object orientation are now extended to the solver component of OpenFOAM.  Whereas application solvers presented equations, algorithm and model interface in a procedural manner, a modular solver is built around a generalised algorithm contained within the foamRun/foamMultiRun applications. Equations and interface are introduced at different levels of the class hierarchy. The solver for a particular fluid type includes little more than the code relating to that particular type. For example, the fluid solver is characterised by solving an energy equation in the thermophysicalPredictor() function, but otherwise inherits its functionality from the isothermalFluid solver.

It becomes clearer when using object orientation that a CFD solution involves time, space (mesh) and fields/properties that define the fluid/solid state.  The foamRun/foamMultiRun applications initialise time and space and evolve time with the time loop. The moveMesh() function evolves space. Each modular solver represents a type of fluid/solid and is responsible for initialising and evolving the fluid/solid state within the pre-defined functions in foamRun/foamMultiRun. Developments focus on evolving the state properties specific to the particular fluid type, e.g. phase fraction in multiphase flow.

Modular solvers provide a distinct advantage for coupled solutions across multiple fluid and/or solid regions. Rather than writing customised application solvers for every combination of fluid and solid types, foamMultiRun allows any modular solver to be applied on any number of mesh regions. It is rapidly being adopted to calculate conjugate heat transfer with a variety of multiphase flows. Users no longer have an overwhelming choice of solvers, as their number decreases. Switching between steady and transient solutions no longer requires a change of solver. Physical models can be optionally selected through case configuration files, including buoyancy which is automatically activated by the presence of the p_rgh pressure field in the start-time directory.

Finally, the solvers in OpenFOAM capture the essence of its software design, 30 years after the first solver application, icoFoam, was written. As Andy Warhol once said: “They always say time changes things, but you actually have to change them yourself.”

Modular Solvers in OpenFOAM