Log in
Log in Simulate Now
Log in Simulate Now

Solve Verification: The Method of Manufactured Solutions

By Jeremy Theler 24 September 2021

An application to the steady-state thermal problem in OnScale Solve

The fact that you are reading this blog probably means that you are using OnScale for simulation or are considering it. And since OnScale may be a new tool in your analysis workflow you want to know that it is accurate. This post describes how we use the Method of Manufactured Solutions (MMS) to verify the correctness of our simulation solvers in a rigorous and comprehensive way.

MMS covers the verification part of a “Verification & Validation” plan commonly used to test simulation correctness. These two V&V steps are roughly defined as

The process of determining that a model implementation accurately represents the developer’s conceptual description of the model and the solution to the model.
Are we solving right the equations?
The process of determining the degree to which a model is an accurate representation of the real world from the perspective of the intended uses of the model.
Are we solving the right equations?

In this post, we are going to focus on the first of these two subjects and use an example to show how we verify the computer code for solving a thermal analysis problem in OnScale. It should be noted that V&V is a wide and complex subject. This is an introductory overview.

Exact and manufactured solutions

The traditional way to check if a simulation program “solves the equations right” is to compare its results with known analytical results. This is called the Method of Exact Solutions (MES). This comparison can only be made by solving problems that do have analytical solutions, which usually involve simple geometries, properties and boundary conditions. That is to say, problems with space-dependent material properties or with non-linear boundary conditions which usually lead to either difficult or impossible analytical solutions. When applying MES, one has to find an analytical expression for the solution of a partial differential equation as a function of space and eventually time. This procedure involves an integration procedure which—when possible—is hard and time consuming. Note that we require exact analytical solutions for quantitative comparison, whereas comparing against experimental results would fall under “validation” of the simulation.

On the other hand, the basic idea of the Method of Manufactured Solutions is to go the other way round: start by proposing a certain solution and find out which source terms are needed to obtain that solution. The source terms can be obtained by simple differentiation. Let us focus on the the steady-state heat conduction problem to illustrate how this is done.

warning: technical content ahead

Say we propose (i.e. manufacture) a certain temperature distribution T_\text{mms}(\mathbf{x}) as an algebraic expression of known functions such as sine and cosine. If we also manufacture an expression for the thermal conductivity k(\mathbf{x}), we can easily find which source term q'''(\mathbf{x}) is needed so as to obtain the proposed temperature distribution as the solution of the equation

\text{div}\Big[ k(\mathbf{x}) \cdot \text{grad} \left[ T_\text{mms}(\mathbf{x}) \right] \Big] = q'''(\mathbf{x})\qquad(1)

Indeed, while the MES involves the integration of the differential equation—which is hard and maybe impossible—the MMS approach involves computing derivatives of known algebraic functions—which is easy and always possible since the choice for T_\text{mms}(\mathbf{x}) is arbitrary but very flexible. This procedure can be done “by hand” or be automated with the aid of symbolic algebra software. Then the derived material properties, volumetric heat source term and boundary conditions are set through algebraic expressions. The temperature distribution T(\mathbf{x}) computed by the simulation can then be compared to the manufactured temperature distribution T_\text{mms}(\mathbf{x}) to see how well they match.

By repeating this procedure with different mesh densities, we can find the rate at which the error decreases as we use smaller elements to solve the problem. That is to say, we can compute the order of convergence p of the computational method. If this order p matches the theoretical order of the method, we can then conclude that the solver solves the equations right.

Application of MMS to the thermal problem in OnScale

Let us solve a thermal conduction problem in three dimensions over a unit cube [0,1]\times[0,1]\times[0,1] with Onscale, using algebraic expressions to set the source term. Let our non-dimensional “manufactured temperature” be

T_\text{mms}(x,y,z) = 1 + \sin ^2\left(2\,x\right)\,\cos ^2\left(3\,y\right)\,\sin \left(4\,z \right)

To keep the mathematical details as simple as possible, we choose a uniform conductivity identically equal to a non-dimensional value of one over the whole cube.

Then, inserting eq. 2 into eq. 1 and computing the gradient and the divergence (with the aid of a symbolic mathematics software such as Maxima or Mathematica) we can compute an analytical expression for the source term q'''(x,y,z).1

By performing a parametric study with successively more refined meshes (meaning smaller elements), we can look at how the difference between the manufactured solution T_\text{mms}(x,y,z) and the simulation solution T(x,y,z) decreases. We plot this error versus the mesh element size h on a log-log scale so that the slope of the line gives us the order of convergence of the method. fig. 1 shows the results of applying the MMS procedure to OnScale’s steady-state thermal analysis.


Figure 1: L_2-norm of the numerical error vs. element size h for different meshes

The results show that as the element size h is reduced, the error also decreases. A monotonic decrease in error with third-order convergence O(h^3) as the mesh is refined as expected for the second-order structured tetrahedral elements used in this analysis. Obtaining the expected convergence rate verifies that the equations (including material behavior, boundary and source terms) are being solved correctly.

The “wiggles” seen in the curves for unstructured grids are also expected, since the tetrahedra might need to accommodate unevenly and hence introduce inhomogeneities in the mesh. The unstructured cases were deliberately created with a very low quality so as to test the solver’s robustness. A comparison with a cutaway view of the solution T(x,y,z) for both structured and unstructured tetrahedra is shown in fig. 2.




Figure 2: Illustration of structured and unstructured tetrahedral meshes in a cube.. a — Structured grid, b — Unstructured grid


This blog post describes part of the rigorous verification processes that OnScale uses to ensure that it’s simulation solvers are accurate. Validation of OnScale solvers (ensuring that they correctly represent real-world behavior) is also very important.

Now that you have made it this far, put us to the test! Try OnScale here, and start simulating easier than ever before, with the accuracy that professional simulation tools deserve.

  1. The full analytical expression of the source term is
    q'''(x,y,z) = -18\,\sin ^2\left(2\,x\right)\,\sin ^2\left(3\,y\right)\,\sin \left(4\,z\right)\\ + 42\,\sin ^2\left(2\,x\right)\,\cos ^2\left(3\,y \right)\,\sin \left(4\,z\right) \\ -8\,\cos ^2\left(2\,x\right)\,\cos ^2 \left(3\,y\right)\,\sin \left(4\,z\right)
Jeremy Theler
Jeremy Theler

Jeremy is a Solver Developer at OnScale, creating the next-generation finite-element solvers' Reflex in C++ using MPI, PETSc, VTK and state-of-the-art cloud and web technologies, i.e., "the way FEM solvers are supposed to be."

Discover how customers like you found success by
leaving traditional engineering simulation behind

Try OnScale following
our simulation guides

Simulate Now

Discuss your engineering
applications with us

Request a Demo