Accurate projection methods for the incompressible Navier—Stokes equations

TL;DR: In this article, the authors consider the accuracy of projection method approximations to the initial-boundary-value problem for the incompressible Navier-Stokes equations and present an improved projection algorithm which is fully second-order accurate.

About: This article is published in Journal of Computational Physics.The article was published on 2001-04-10 and is currently open access. It has received 841 citations till now. The article focuses on the topics: Projection method & Pressure-correction method.

Summary

1. INTRODUCTION

• This paper considers the accuracy of projection method approximations to the initial– boundary-value problem for the incompressible Navier–Stokes equations.
• In recent years, projection methods which exactly enforce a discrete divergence constraint, or “exact” projection methods, have often been replaced with “approximate” projection methods (e.g., [3, 4, 26, 28]), which are similar to pressure-Poisson methods in that the velocity satisfies a discrete divergence constraint only to within the truncation error of the method.
• Additionally, as with all fractional step methods, a crucial issue is how boundary conditions are determined for some or all of the intermediate variables.
• In the process the authors show that several existing methods fall short of second-order accuracy up to the boundary precisely because this coupling was not considered.
• In particular, the analysis shows that second-order accuracy in both the velocity and the pressure are obtainable with the correct choice of boundary conditions and pressure-update equations.

2. COMMENTS ON SOME EXISTING METHODS

• In this section the authors make brief comments about some of the methods mentioned earlier viewed in the context established in the Introduction.
• The authors also comment on reported results that have contributed to the debate on the topic.
• The authors demonstrate later in this paper that for this method, u∗ differs at most byO(1t2) from the correct velocity un+1, justifying the use of the velocity boundary condition for u∗.
• Later, the authors show using normal mode analysis that this is also a necessary condition for second-order accuracy for this method.
• If both the pressure and φ are evaluated at the same time level (i.e., if the right-hand side of Eq. (24) is set equal to pn+1), the resulting pressure is only first-order accurate, as reported by Strikwerda and Lee [36].

3. BOUNDARY CONDITIONS

• The numerical methods presented in the last section require the solution of implicit equations for which boundary conditions must be imposed.
• Besides the implicit momentum Eqs. (8) and (19), the implementation of a projection also requires a boundary condition.
• ∇χn+1 can be approximated by extrapolating the values from previous time steps as proposed by E and Liu [15].
• This boundary condition must be chosen so that when u∗ is projected to yield un+1, the tangential boundary condition on un+1 is satisfied.
• These choices will be analyzed in detail in the following section.

4. NORMAL MODE ANALYSIS

• The original Dirichlet problem as stated in Eqs. (1)–(3) requires only a condition on the velocity u on the boundary.
• The authors consider the second-order-in-time formulation and include in the analysis the extrapolation of the boundary values of χ .
• In order to obtain an accurate solution to the incompressible Navier–Stokes equations using the projection methods described by Eqs. (8)–(12), one either must devise a procedure for accurately approximating the boundary conditions u∗ −1t∇φ = (α, β)T or reformulate the problem in such a way that u∗ is a sufficiently accurate approximation to u.
• By examining the size of the remaining terms, the following results are evident: PmI uses q = pn−1/2 and L = I .
• As explained before, the pressure is given by p̂ = κ 3/2 κ − 1 φ̂, which includes the coefficient of the spurious mode κ.

5. THE NUMERICAL METHODS

• This section describes the numerical methods that will be applied to the full Navier– Stokes equations.
• 3. A Projection Method without Pressure Gradient.
• In the numerical methods presented above, extrapolation in time is used to compute the time-centered advective derivatives as well as the tangential boundary conditions for the implicit treatment of the momentum equation.
• Since these terms cannot be extrapolated at the first time step, an iterative procedure is employed.
• A concern relating to the fact that the tangential component of the velocity boundary condition is not satisfied exactly remains to be addressed.

6. NUMERICAL RESULTS

• In this section numerical examples are presented which confirm the validity of the normal mode analysis presented in Section 4 for the gauge and projection methods.
• Since the first-order temporal error terms for the pressure in the normal mode analysis are scaled by the viscosity, it is important that the viscosity be large enough compared to the grid size so that 1x2 ¿ ν1t .
• To illustrate this, the forced flow problem was recomputed using a different tangential boundary condition for u∗ (or mn+1) for each method.
• For this test, the lagged value τ̂ · ∇χn is used instead, which results in a loss of accuracy in both the velocities and pressure.
• The order of the method is not changed.

Errors the Pressure for the Forced Flow Test Problem

• The same lagged boundary condition as above can also be used for projection method II.
• Again this choice decreases the size of the errors somewhat, but the order of the method is not changed.
• For projection method III, the normal mode analysis indicates that using the lagged value τ̂ · ∇n is necessary for second-order accuracy.
• For this test, the less accurate boundary condition τ̂ · u∗ = 0 was used (as is done normally done for PmII) which results in a loss of accuracy in both the velocities and the pressure.
• The same periodic channel geometry is used with zero boundary conditions at the bottom wall,.

Different Boundary Extrapolations Are Used

• For the gauge method, the initial condition m = u is used and the boundary condition m · n̂ = 0 is specified at both top and bottom boundaries throughout the computation.
• Since no exact solution is known, a reference solution was computed on a 1152× 1152 grid, and errors are estimated by the difference from this solution.
• To assure that the reference solution being used is valid, both the impulse method and PmII were used to compute the solution; it was observed that the maximum difference between the two solutions was 1.31× 10−6 in the velocity, 2.23× 10−6 in the pressure, and 8.73× 10−5 in py .

Flow Test Problem

• Significantly smaller than the estimated errors used to compute the convergence rates, using the reference solution is justified.
• It should be noted that the standard Richardson extrapolation techniques commonly employed to estimate convergence rates can be misleading in this context.
• In particular, the pressure gradient computed with projection method I will appear to converge quite nicely at the boundary if only a Richardson procedure is used.
• Since the flow is not forced except by the motion of the top wall, the magnitude of the v-component of the velocity decays rapidly while that of the u-component increases throughout the run at the top wall.
• The errors are estimated at time 0.25 in the u-component of the velocity and the pressure when the maximum value of u is about 0.86, while the maximum of v has dropped to about 0.39.

Errors in the Pressure for the Unforced Flow Test Problem

• Projection method I displays second-order accuracy in the velocity but only first-order accuracy in the pressure.
• Since an approximate projection is being used, the discrete divergence of un will not be zero for any of the methods.
• The size of the boundary layer has decreased an order of magnitude to the size of that in the interior.
• Be made is that although the normal boundary condition for u∗ is mathematically arbitrary, the choice can affect the accuracy of the numerical solution.
• 5. Smoothness of the Pressure Error Despite the fact that projection methods II and III display optimal convergence rates in the pressure, the pressure error is not a completely smooth function near the solid wall boundaries.

Errors in the Unforced Flow Test Problem for Projection Method III Using the Free Boundary Condition

• The slightly irregular shape of the error, the overall size is still converging to zero at a second-order rate.
• The slight irregularities in the pressure error create noticeable irregularities in the error of py .
• Several comments can be made based on the data.
• (Note that PmIII was computed using the modified boundary condition for n̂ · u∗.).
• Both projection methods II and III show a decrease in the observed convergence rate for the pressure gradient measured in the L∞ norm.

Errors in py for the Unforced Flow Test Problem

• The cause of the slightly lower convergence rates for the py can again be traced to the lack of smoothness of the Laplacian term in the pressure-update Eq. (74).
• The fact that the pressure itself is converging at the optimal rate indicates that the drop in convergence rates for the gradient is caused by spatial rather than temporal error.
• Depending on the implementation, the error in the pressure gradient due to a lack of smoothness in the pressure correction terms could potentially be exacerbated by the presence of complex geometries.

7. CONCLUSIONS

• The class of incremental pressure projection methods discussed in this paper is characterized by the choice of three ingredients: the approximation to the pressure gradient term in the momentum equation, the formula used for the global pressure update during the time step, and the boundary conditions.
• If the conditions for u∗ are separated into normal and tangential components, there is apparently some freedom in choosing the normal component since the required boundary condition for the potential φ in the projection step can be adjusted to ensure that n̂ · un+1|∂Ä = n̂ · un+1b .
• As demonstrated by the numerical experiments with PmIII, the choice of the normal boundary condition can affect the smoothness of u∗ near the boundary and therefore can also play a role in the accuracy with which the pressure is recovered.
• The gauge method variable m (equivalent to u∗ during the first time step) is not discarded but used throughout the computation.
• In applications where an accurate representation of the pressure near solid wall boundaries is required, the results in this paper provide an important improvement in accuracy for a popular class of projection methods.

Figures

TABLE II Errors the Pressure for the Forced Flow Test Problem

TABLE I Errors in the u-Component of Velocity for the Forced Flow Test Problem

FIG. 1. First-order boundary layer error for projection method I. The three graphs correspond to profiles atx locations 3/16, 6/16, and 9/16. Each graph shows the error from the 256× (∗) and 384× 384 (o) runs.

TABLE VII Errors in the Divergence of un for the Unforced Flow Test Problem

TABLE IX Errors in py for the Unforced Flow Test Problem

FIG. 5. Error in thepy for the projection method II on the unforced problem. The three graphs correspond to profiles atx locations 3/16, 6/16, and 9/16. Each graph shows the error from the 96× 96 (o) and 192× 192 (∗) runs.

TABLE V Errors in the u-Component of Velocity for the Unforced Flow Test Problem

TABLE VI Errors in the Pressure for the Unforced Flow Test Problem

FIG. 4. Error in the pressure for the projection method II on the unforced problem. The three graphs correspond to profiles atx locations 3/16, 6/16, and 9/16. Each graph shows the error from the 96× 96 (o) and 192× 192 (∗) runs.

TABLE IV Errors in the Pressure for the Forced Flow Test Problem When Different Boundary Extrapolations Are Used

TABLE III Errors in the u-Component of Velocity for the Forced Flow Test Problem When Different Boundary Extrapolations Are Used

TABLE VIII Errors in the Unforced Flow Test Problem for Projection Method III Using the Free Boundary Condition

Frequently asked questions

Q1. What have the authors contributed in "Accurate projection methods for the incompressible navier–stokes equations" ?

In this paper, the authors consider the second-order convergence of a projection method for the Navier-Stokes equations. 

Q2. What is the reason for the drop in convergence rates for the gradient?

The fact that the pressure itself is converging at the optimal rate indicates that the drop in convergence rates for the gradient is caused by spatial rather than temporal error. 

Q3. How is the inversion of the Laplacian in the projection made efficient?

Since the test problems studied in the next section are all set in a periodic channel, the inversion of the Laplacian in the projection is made efficient by first taking the discrete Fourier transform of the equation in the x-direction. 

Q4. How many errors are estimated in the velocity and pressure?

The errors are estimated at time 0.25 in the u-component of the velocity and the pressure when the maximum value of u is about 0.86, while the maximum of v has dropped to about 0.39. 

Q5. What is the class of incremental pressure projection methods discussed in this paper?

The class of incremental pressure projection methods discussed in this paper is characterized by the choice of three ingredients: the approximation to the pressure gradient term in the momentum equation, the formula used for the global pressure update during the time step, and the boundary conditions. 

Q6. How can the authors approximate the solution of the coupled system?

a fractional step procedure can be used to approximate the solution of the coupled system by first solving an analog to Eq. (5) (without regard to the divergence constraint) for an intermediate quantity u∗, and then projecting this quantity onto the space of divergence-free fields to yield un+1. 

Q7. What is the definition of pressure-Poisson methods?

Methods are often categorized as “pressure-Poisson” or “projection” methods based on which form of the elliptic constraint equation is being used. 

Q8. What is the second-order method proposed by Perot?

The second-order method proposed by Perot uses q = 0 and replaces the pressure-update formula (12) with(I + ν1t 2 ∇2 ) pn+1/2 = φn+1. (25)This method still only obtains first-order convergence in the pressure since n̂ · ∇ p = 0 is the boundary condition used for the elliptic pressure equation. 

Q9. What is the continuity of u in time?

The continuity of ∇φ in time is implied by the fact that u∗ satisfies an elliptic equation with continuous forcing and 1t∇φ is simply (I− P)u∗. 

Q10. What is the difference between projection method II and the other methods?

projection method II has substantially less error in the divergence of un than the other methods, and this error appears to be converging to zero at a higher rate than the other methods. 

Q11. What is the need for the extrapolation of the velocity and pressure?

In the next section it is shown that this extrapolation is necessary for the resulting velocity and pressure to be second-order accurate in the maximum norm.