Mystery surrounding order of accuracy in CFD...

Just came back from a trip to Beijing, where I gave two lectures in a Workshop, one on high-order CFD methods, and the other on the quadrature method of moment for multiphase flows. On the last day, I got to see the beautiful Great Wall at Mu Tian Yu. See the top picture for a magnificent view! Personally I believe this is a better tourist attraction than the more famous Great Wall at Ba Da Lin because it is less crowded. I am now paying a price for the trip: I am still recovering from the jet lag.

In the Workshop at Beijing as well as in many other technical conferences, I often heard people talking about and making mistakes on order of accuracy (OOA) in CFD. As often as it is mentioned, there appears to be as much confusion about it too. There are two major mistakes regarding OOA:

1. If a flow solver produces a more accurate solution than another, it must posses  higher order accuracy. 
2. The OOA is the same as the order of the polynomial approximating the (unknown) solution.

I still remember an incident in an AIAA conference (maybe at Reno). A guy believed like a religion that the OOA is the same as the order of the reconstructed solution polynomial in a finite volume scheme. In fact, he accused the finite volume community of lying about the accuracy of their favorite schemes. You cannot help feeling bad for the guy with such a strong (but wrong) conviction. The truth is that the OOA of a finite volume scheme with a degree p solution polynomial is (p+1).

Let's define OOA. A scheme is pth order if the solution error is proportional to (dx)^p, where dx is the mesh size, i.e.,

Error = C(dx)^p,

where C is a constant depending on the particular scheme. Or equivalently,

log(Error) = p*log(dx) + log(C).

In other words, p is the slope of the dx-Error plot in a log-log scale. Strictly speaking, OOA has nothing to do with the magnitude of the solution error, but only related to the rate of error reduction with mesh refinement, i.e., the slope.

How do you actually measure the OOA? If you know the exact solution, you need only perform two simulations: one on a coarse mesh, and the other on a fine mesh with half the mesh size of the coarse mesh. With a simple derivation, you can obtain the order p from

Error_(dx)/Error_(dx/2) = 2^p, or p = log_2[Error_(dx)/Error_(dx/2)]

where Error_(dx) is the error on the coarse mesh, and  Error_(dx/2) the error on the fine mesh.

More later on high-order methods...