From my experience, I can say that convection terms are the most dangereous terms for numerical stability compared to dissipation terms. For example when we approximate the dissipation term on RHS we have obtained a stability condition and as long as we choose the free parameters according to this restriction, the scheme was stable. Also, when we approximate the dissipation term on the LHS, the scheme was stable. Therefore, for dissipation terms we were able to have numerical stability even in explicit case. On the other hand, for convection terms there was no numerical stability in explicit case and always stability in implicit case. This is like choosing between solution and no solution or convergence and divergence. We didnot see it in our model equation (49) before, however, most of the times when convection terms are approximated on the explicit side, the finite difference equations have a CFL number restriction for numerical stability as we will see in the next model equation. At high Reynolds numbers, convection terms dominate over other terms in the equations in magnitude and become even more dangerous for numerical stability. Therefore in developing our numerical method as a rule of thumb we will ALWAYS try to keep the convective terms on the implicit (LHS) side.

Statements 1 and 2 are basicaly saying that implicit equations have better numerical stability characteristics, which is actually an abvious statement.

Let us have a general equation such that

(84)

If we approximate the convection and dissipation terms on the explicit side we get

(85)

If we apply a stability analysis to this equation we can find that for numerical stability the following conditions have to be satisfied

(86)

----------------------------------------------------------

Von Neuman Stability Analysis:

We have the following equation

(87)

The finite difference equation will be

(88)

Then, following the same procedure explained before, the amplification factor will have the following form

(89)

Equation (89) describes an ellipse that is centered on the positive real axis at  and the major and minor axes of this ellipse are  and For absolute stability the following restrictions have to satisfied

(90)

Reference:

Tannehill, Anderson and Pletcher, Computational Fluid Mechanics and Heat Transfer, Taylor & Francis, page 112 (last paragraph)

----------------------------------------------------------

How ever, if we go back to equation (84) and approximate the convection and dissipation terms on the implicit side we get

(91)

or in operator notation

(92)

This equation is unconditionally stable.

----------------------------------------------------------

Von Neuman Stability Analysis:

We have the following equation

(93)

The finite difference equation will be

(94)

Then, the amplification factor will have the following form

(95)

Note that  is complex, multiplying with complex conjugate we get

(96)

The magnitude of the amplification factor, then, will be

(97)

The magnitude of the amplification factor is always less than unity such that , and therefore this scheme is unconditionally stable.

----------------------------------------------------------

From all of these statements (1 and 2) and stability analysis of model equations we have seen up to now, we deduce the following rules that we are going to follow in our numerical method.

Rule 1:

We are going to approximate the dissipation terms on the implicit side. For some reason, if we have to, eventhough we donot want to, approximate dissipation terms on the explicit side this is still OK.

Rule 2:

We are going to approximate the convection terms on the implicit side ONLY. We DO NOT WANT any convection term on the explicit side, at least not any first order () convection terms.

At this point, let us apply these two rules to the governing equations we want to numerically solve, the 2-D Steady Incompressible Navier Stokes equations in streamfunction and vorticity formulation