Method of Manufactured Solutions

Nowadays, we often use the tools of computational analysis in engineering design. For the results of such analyses to be believable, the tools that are used have to pass rigorous tests. There are two categories of tests involved:

  1. Verification tests where one checks whether the mathematical model has been implemented according to design.
  2. Validation tests where one checks how closely the correctly implemented mathematical model mimics reality.

A method of verification of a numerical code that has received some attention is is the method of manufactured solutions. A systematic application of the method for complex nonlinear constitutive models in solid mechanics is yet to be found in the literature. The fluid mechanics community, however, uses the method on a regular basis.

Last year, at the World Conference on Computational Mechanics in Los Angeles, Prof. R. C. Batra told me that he had been using the method to verify his codes for a while. However, he had never felt the need to publish the approach because of its seemingly trivial nature. In this post I will give you an example of the method of fictitious body forces used by Prof. Batra [1,2].

The method of fictitious body forces involves assuming a displacement field over the body; finding the body force, initial conditions, and boundary conditions that fit that solution; and then using that body force and the computed boundary conditions and initial conditions in the numerical algorithm to arrive at a solution. If the computed solution matches the assumed displacement, then we're good to go.

The example below was meant to provide a starting point for a student. Please let me know if you find it useful or whether the equations are confusing rather that clarifying.

A 1-D example:

Consider a body in which the material point $ \ensuremath{\mathbf{X}}$ moves to $ \ensuremath{\mathbf{x}}$ during a motion. The motion is governed by the equations

$\displaystyle \begin{aligned} &\ensuremath{\ensuremath{\ensuremath{\boldsymbol{... ...ensuremath{\boldsymbol{\nabla}}_0 \ensuremath{\mathbf{u}}})^T \\ \end{aligned}$

where

$\displaystyle \ensuremath{\mathbf{x}}= \ensuremath{\boldsymbol{\varphi}}(\ensur... ...{\sigma}}= J^{-1} \ensuremath{\boldsymbol{P}} \ensuremath{\boldsymbol{F}}^T . $

If we now consider a one-dimensional bar of the material, then

$\displaystyle x_1 = \varphi(X_1, t) ;   F_{11} = 1 + \ensuremath{\frac{\partia... ...{\partial X_1}} ;   J = 1 + \ensuremath{\frac{\partial u_1}{\partial X_1}} . $

The left Cauchy-Green deformation tensor is

$\displaystyle B_{11} = 1 + 2\ensuremath{\frac{\partial u_1}{\partial X_1}} + \left(\ensuremath{\frac{\partial u_1}{\partial X_1}}\right)^2 . $

The Cauchy stress is given by (using a simplified constitutive model):

$\displaystyle \sigma_{11} = \ensuremath{\frac{1}{2}}  C (B_{11} - 1) = \ensure... ... X_1}} + \left(\ensuremath{\frac{\partial u_1}{\partial X_1}}\right)^2\right] $

where $ C$ is an elastic constant. Also, the Cauchy stress and the 1st Piola-Kirchhoff stress are related by

$\displaystyle \sigma_{11} = \frac{1}{F_{11}} P_{11} F_{11} = P_{11} . $

The one-dimensional momentum equation is

$\displaystyle \ensuremath{\frac{\partial P_{11}}{\partial X_1}} + \rho_0 b_1 = \rho_0 \ddot{u_1} . $

 

The method of manufactured solutions:

Let us try to apply the method of fictitious body forces to the bar problem.

The first step is to assume that the displacement field is given by (you can choose any other appropriate functional form)

$\displaystyle u_1(X_1, t) = X_1^2 \sin(\omega t) . $

Therefore,

$\displaystyle \ensuremath{\frac{\partial u_1}{\partial X_1}} = 2 X_1 \sin(\omeg... ...uremath{\frac{\partial^2 u_1}{\partial t^2}} = -X_1^2 \omega^2 \cos(\omega t) $

and the deformation gradient is

$\displaystyle F_{11} = 1 + 2 X_1 \sin(\omega t) . $

The left Cauchy-Green deformation is

$\displaystyle B_{11} = 1 + 4 X_1 \sin(\omega t) + 4 X_1^2 \sin^2(\omega t) . $

The first Piola-Kirchhoff stress is

$\displaystyle P_{11} = \ensuremath{\frac{1}{2}} C [4 X_1 \sin(\omega t) + 4 X_1^2 \sin^2(\omega t)] . $

Therefore,

$\displaystyle \ensuremath{\frac{\partial P_{11}}{\partial X_1}} = 2 C \sin(\omega t)[1 + 2 X_1 \sin(\omega t)] . $

Plugging $ \ensuremath{\frac{\partial P_{11}}{\partial X_1}}$ and $ \ensuremath{\frac{\partial^2 u_1}{\partial t^2}}$ into the momentum equation, we get

$\displaystyle 2 C \sin(\omega t)[1 + 2 X_1 \sin(\omega t)] + \rho_o b_1 = -\rho_0 X_1^2 \omega^2 \cos(\omega t) . $

Therefore, the body force is

$\displaystyle \boxed{ b_1 = - X_1^2 \omega^2 \cos(\omega t) - \cfrac{2 C}{\rho_0}  \sin(\omega t)[1 + 2 X_1 \sin(\omega t)] . } $

If the bar is of length $ L$ , the boundary conditions are

$\displaystyle \boxed{ u_1(0, t) = 0 \qquad\text{and}\qquad u_1(L, t) = L^2 \sin(\omega t) . } $

 

The initial conditions are

$\displaystyle \boxed{ u_1(X_1, 0) = X_1^2 \sin(0) = 0 . } $

When you apply the body force $ b_1(X_1, t)$ and the initial and boundary conditions, you should get a solution that matches the chosen function.

Note that all these have been done in a Lagrangian configuration. You could alternatively do the same assuming a function $ \ensuremath{\mathbf{u}}(\ensuremath{\mathbf{x}}, t)$ and transforming the equations accordingly.

 

References

 

1
R. C. Batra and X. Q. Liang.
Finite dynamic deformations of smart structures.
Computational Mechanics, 20:427-438, 1997.

 

2
R. C. Batra and B. M. Love.
Multiscale analysis of adiabatic shear bands in tungsten heavy alloy particulate composites.
Int. J. Multiscale Computational Engineering, 4(1):95-114, 2006.

Very nice, thanks.

Very nice, thanks. I tested my hyperelasticity codes using simple linear displacement fields or the Poynting effect. The test and the procedure that you have documented here provide very nice addition tests. The more tests the better, because "there is always one more bug" :-).


Manufactured solutions

Nachiket,

Thanks for your comment.   As more and more complex, multiphysics, multiscale computations are attempted questions about their accuracy have become important.  You will often see beautiful movies of simulations of various physical processes.  However, the more one works in computational science the more one realizes that most computational results are to be taken with a large grain of salt.

Unless consistent standards for verification and validation are established (such as the ASTM standards for testing materials) one can never be sure that computational results are correct.   I feel that a library of manufactured solutions will be a step in that direction.   There are numerous complicated nonlinear constitutive models that are used in solid mechanics simulations.  We don't have good ways of deriving manufactured solutions (in 3-D) for a large subset of these models.  In fact, I don't know of any papers that specifically deal with these issues.

If we want our computational results to become the third approach in science (after theory and experiment) then it's time to take verification and validation seriously.

Biswajit 


Mogadalai Gururajan's picture

Open sourcing the codes

Dear Biswajit,

A major boost for computational results to become the third approach in science (in addition to verification and validation), would be to make the code used to generate the results open (In fact, I think like the video category in iMechanica, we should have a repository for codes). It helps remove the bugs fast, and gives enormous confidence in the results presented since anybody who does not believe the answers can always go check it for herself. In the present system, while experiments can be redone, and theoretical calculations can be checked, if there is a problem in reproducing the reported simulation results, it is not clear whether the problem is a bug in the original code or that of the researcher trying to reproduce the results; this leads to lots of frustration, duplication of work, and general lack of confidence in simulation results.


Open source codes

Guru,

Many source codes are available on the web.  But I haven't see any serious extensions of those codes other than by the original developers.  I feel that the main reason is a lack of documentation. 

Even though our code (Uintah) is available upon request for free along with a host of visualization tools, we have not seen much interest in the code outside of a few people that we know.  We do document parts of the code when we have the time and we tried a formal release last year.  Ours is a research code and it changes all the time.  A systematic set of documents that a new developer can read and use in a reasonably short period of time does not exist (and is never a priority).

Also, the mechanics community is relatively small and the number of serious propgrammers in that community is even smaller.  It is hard to get together a team of developers that can add new features and verify/validate them without any compensation.  Writing multiphysics codes needs input from a variety of different disciplines.  Getting people from different fields to agree on the algorithmic framework and on how the different physics should communicate is a nontrivial problem.

But someone has to take the lead and do it.  Why not you?  I've thought about it and calculated that it's a multiyear project just to get a reasonable range of solid mechanics functionalities into an open code if there are five developers.  If you can get more people interested, then the time needed will be shorter.

Regarding experiments:

I  have spent much of the last couple of years poring through papers that present experimental data on stress-strain behavior at high strain rates and high temperatures.   The amount of variability in the data has to be seen to be believed.  As an example, consider the initial yield stress of 6061-T6 aluminum.   If you fit a curve through that data you will find that around 10% of the data fall well outside the 95% confidence interval.  It you look at the rate of hardening and the hardening modulus data you could well be looking ant entirely different materials!   

The question that arises is how do you know which experimental data to trust?  Why are the data so variable?  There are many reasons, for example, the natural variability of the samples, the heat treatment, the rolling direction etc. but, most importantly,  were the experiments done correctly?  We can never know for sure unless the experiments are repeated by other groups and reasonably close results are obtained.  Otherwise, the experiments are just qualitative and cannot be used for computations - particularly where small differences in behavior can have major consequences.

However, if the experiments follow a clearly described standard then you can assume that the results are accurate and that what you're seeing is just material variablity independent of the actual apparatus used.  I would like to see similar standards for computations too so that you can quantify the uncertainty in the results and say things like "it's 80% certain that what this simulation shows is accurate to within 90% of the error measure".

And, finally, it would also help greatly if the raw data generated by researchers were available for public consumption.  Then people could go through the data and come up with interpretations that the original researchers had missed.

Biswajit 

 


Analytical solutions

Your idea of a library of manufactured solutions is indeed a good one. It provides a great way of verifying codes. A library of analytical solutions combined with benchmark computational examples will be a great verification tool. A complete analytical solution, though extremely useful,  is not always necessary -- even something like the maximum strain in this particular problem will not exceed such and such value -- will be useful. I don't know of any such "benchmark suite" relating to mechanics. In contrast, there is a "Die Hard" set of tests for random number generators http://en.wikipedia.org/wiki/Diehard_tests or http://plato.asu.edu/sub/benchm.html for optimization or http://math.nist.gov/MatrixMarket/.  for linear solvers. The only mechanics related collection I am aware of is http://www.columbia.edu/~ma2325/femarket/ -- and it seems to be currently unavailable. 


Saurabh Puri's picture

thanks

HI Biswajit, Its a really good example


MMS for Navier-Stokes Equations

The Method of Manufactured Solutions has been used for fluid mechanics & CFD too.

Wendl, M. C.
(1999)
``A Method for Obtaining Benchmark Navier-Stokes Solutions'',

International Journal for Numerical Methods in Fluids,

31 (3), 657-660.

reviews the history of MMS in CFD and gives a general method to construct Navier-Stokes solutions...