# Why are upwind schemes stable in convection flow calculation?

It is well known that upwind schemes are stable when calculating convection flows with $|\text{Pe}|>2$, $\text{Pe}$ is the Peclet number. Why is that, and why is central difference unstable?

Is there any intrinsic reason there?

Any explanation, references, links will be helpful.

• This is not a direct answer, so I'll just post a comment, but one article in this field which I found remarkably insightful is this: H.-G. Roos, "Ten ways to generate the Il'in and related schemes" – cfh May 23 '15 at 22:32
• The truthfulness of the assertion in the question depends on how you discretize in time. For instance, (centered differences + backward Euler) is stable. – David Ketcheson May 27 '15 at 18:28

## 4 Answers

The reasoning for the stability of upwind schemes based on an understanding of the characteristics of the hyperbolic equation(s). Characteristics are essentially the finite speeds at which information in a hyperbolic system travel, and are found via decomposing a hyperbolic system into independent hyperbolic PDEs.

Now characteristics are essentially just pushing along the initial conditions of a given hyperbolic equation (though nonlinear equations can distort them). The fact this speed is finite results in a need to be careful with your numerical stencil.

The typical example to illustrate the need to sample carefully is to imagine an initial condition of:

$$u_0(x) = \begin{cases} 1 & x < 2\\ 0 & 2 \leq x \end{cases}$$

Since the characteristic just pushes this initial condition along the space-time domain, its general shape and the center discontinuity remains. Now imagine some time later, you are aiming to evaluate the derivative at some x* location. There's a chance that the x* location is at the center of the discontinuity, where the derivative is technically undefined. To counter this, we sample on one side. Now how do we decide which side to sample on? It should be the side that would change value first.. Which is the side the characteristic would touch first. This means that if the characteristic is moving in the positive x direction, we should sample to the left to ensure we capture any possible changes in the solution.

This sampling of data that is going in the opposite direction of the characteristic is known as upwinded sampling. This helps ensure numerical stability by only sampling data we know we have information for.

An intuitive reason is that upwinding can be viewed as providing additional numerical diffusion, which is typically associated with stability in standard Finite Difference/Finite Element schemes. In finite differences, this means that the first order upwind scheme can be rewritten

$$\frac{u_i-u_{i-1}}{h} = \frac{u_{i+1}-u_{i-1}}{2h} - \frac{h}{2} \frac{u_{i+1}-2u_i + u_{i-1}}{h^2}.$$

In other words, the first order upwind difference can be interpreted as adding additional artificial diffusion relative to the 2nd order central difference scheme.

The reason why central differences is unstable is a little more involved. IMO, it's easier to analyze stability in the Finite Element framework. For a variational problem $b(u,v) = l(v)$, the constant $\gamma_h$ in the discrete coercivity condition (where $\|u\|_h$ is some norm on $u$)

$$b(u,u)\geq \gamma_h \|u\|_h^2$$

gives a measure of the stability of the problem (sort of like the smallest singular value tells you the stability of a matrix equation). For standard Galerkin methods (which are related to 2nd order central finite differencing) and Discontinuous Galerkin with central fluxes, $\gamma_h$ is inversely proportional to the Peclet number, and grows more unstable as Peclet number increases. For example, take the 1D Galerkin scheme for constant advection diffusion $u' - \epsilon u = f$:

$$b(u,v) = \int_0^1 u'v + \epsilon \int_0^1 u'v' = \int f v$$

Assuming $0$ Dirichlet boundary conditions, $b(u,u)$ gives

$$\int_0^1 u'u + \epsilon \int_0^1 (u')^2 = \int_0^1 (u^2/2)' + \epsilon \int_0^1 (u')^2 = [u^2/2]_0^1 + \epsilon\int_0^1 (u')^2 = \epsilon\|u\|_{H^1}^2$$

which implies issues with stability as $\epsilon \rightarrow 0$. Note that this implies this issue with stability is usually worse for boundary layer problems (Dirichlet boundary conditions on both sides).

Upwinded Galerkin (i.e. SUPG, streamline diffusion, etc) and upwind DG restore discrete coercivity with respect to a specific norm $\|u\|_h$, such that $\gamma_h$ does not approach $0$ as Peclet number increases. This is also a part of why adding numerical diffusion tends to stabilizes advective flows - additional numerical diffusion adds a coercive term to the variational form.

• Is there any book or paper talk about this in a similar way like yours? It is very appreciated. – Kozuki May 11 '15 at 5:56
• I took most of these from some class notes that Tom Hughes wrote up previously - I believe he was planning on making it a book, but I'm not sure if it's out yet. More mathematical discussions may be found in "Robust Methods for Singular Perturbation Problems" by Roos, Stynes and Tobiska springer.com/us/book/9783540344667. – Jesse Chan May 11 '15 at 15:01

Also, here is a pretty extensive discussion with results of various numerical schemes:

http://www.holzmann-cfd.de/index.php/en/numerical-schemes

I hope to answer your question in a general matter without any references to mathematical equations.

To be stable, a numerical scheme needs to capture the physics of the equation you are discretizing. Without going into the details, each term in your equation has eigenvalues which define the direction in which the information of your flow travels.

For example, when discretizing the diffusion term in the Navier-Stokes Equations (NSE), you should use a central scheme since diffusion, by definition, happens in all directions. If you were to use a one-sided numerical scheme, you could make your scheme stabilized by satisfying the CFL condition, but your results may become non-physical, since it can't exactly capture expansion waves. Where we should use Godunov's family schemes. Please note that Godunov's family schemes use upwinding with some constraints. Now, if you are solving the NSE for a supersonic flow, where the flow characteristics travel in ONLY one direction, then you MUST use an upwind scheme for stability. Good reference for upwind schemes in high speed flow is a paper by van Leer

• anyone care to comment on the downvote? My answer is pretty much the same as the accepted answer without the charasteristics background – solalito May 28 '15 at 9:14
• Yes, your answer is the simplest one to understand. And it is correct within my knowledge. – Kozuki Jun 6 '15 at 9:31