The Visual Room

1. Poisson Equation for Pressure

«  4. von Neumann Stability Analysis   ::   Contents   ::   2. The Marker and Cell Method  »

1. Poisson Equation for Pressure

For compressible flow, pressure and velocity can be coupled with the Equation of State. But for incompressible flow, there is no obvious way to couple pressure and velocity.

So, take the divergence of the momentum equation and use the continuity equation to get a Poisson equation for pressure.

1.1. Poisson’s Equation in Continuous Domain

1.1.2. The Divergence of the Momentum Equation in Partial Notation

Now take the divergence of the momentum equations, if \(\mathbf{M}\) is the “vector” of momentum equations, the divergence is:

\[\nabla \cdot \mathbf{M} = {\partial \over {\partial x}} M_x + {\partial \over {\partial y}} M_y\]
\[{\partial \over {\partial x}} M_x = {\partial \over {\partial x}} {{\partial u} \over {\partial t}} + {{\partial u} \over {\partial x}} {{\partial u} \over {\partial x}} + u {{\partial^2 u} \over {\partial x^2}} + {{\partial v} \over {\partial x}} {{\partial u} \over {\partial y}} + v {{\partial^2 u} \over {\partial x \partial y}} = -{1 \over \rho} {{\partial^2 p} \over {\partial x^2}} + \nu \left( {{\partial^3 u} \over {\partial x^3}} + {{\partial^3 u} \over {\partial x \partial y^2}} \right)\]
\[{\partial \over {\partial y}} M_y = {\partial \over {\partial y}} {{\partial v} \over {\partial t}} + {{\partial u} \over {\partial y}} {{\partial v} \over {\partial x}} + u {{\partial^2 v} \over {\partial x \partial y}} + {{\partial v} \over {\partial y}} {{\partial v} \over {\partial y}} + v {{\partial^2 y} \over {\partial y^2}} = -{1 \over \rho} {{\partial^2 p} \over {\partial y^2}} + \nu \left( {{\partial^3 v} \over {\partial x^2 \partial y}} + {{\partial^3 v} \over {\partial y^3}} \right)\]

1.1.2.1. Add the LHS

\[{\partial \over {\partial x}} M_x + {\partial \over {\partial y}} M_y = {\partial \over {\partial t}} \left( {{\partial u} \over {\partial x}} + {{\partial v} \over {\partial y}} \right) + \left( {{\partial u} \over {\partial x}} \right)^2 + {{\partial u} \over {\partial y}} {{\partial v} \over {\partial x}}+ u {{\partial^2 u} \over {\partial x^2}} + u {{\partial^2 v} \over {\partial x \partial y}} + {{\partial v} \over {\partial x}} {{\partial u} \over {\partial y}} + \left( {{\partial v} \over {\partial y}} \right)^2 + v {{\partial^2 u} \over {\partial x \partial y}} + v {{\partial^2 v} \over {\partial y^2}} = RHS\]

Re-arrange:

\[{\partial \over {\partial x}} M_x + {\partial \over {\partial y}} M_y = {\partial \over {\partial t}} \left( {{\partial u} \over {\partial x}} + {{\partial v} \over {\partial y}} \right) + \left( {{\partial u} \over {\partial x}} \right)^2 + 2 {{\partial u} \over {\partial y}} {{\partial v} \over {\partial x}}+ u {\partial \over {\partial x}} \left( {{\partial u} \over {\partial x}} + {{\partial v} \over {\partial y}} \right) + \left( {{\partial v} \over {\partial y}} \right)^2 + v {\partial \over {\partial y}} \left( {{\partial u} \over {\partial x}} + {{\partial v} \over {\partial y}} \right) = RHS\]

Apply Continuity, so \(\nabla \cdot \mathbf{u} = 0\). Hence:

\[{\partial \over {\partial x}} M_x + {\partial \over {\partial y}} M_y = \left( {{\partial u} \over {\partial x}} \right)^2 + 2 {{\partial u} \over {\partial y}} {{\partial v} \over {\partial x}}+ \left( {{\partial v} \over {\partial y}} \right)^2 = RHS\]

1.1.2.2. Add the RHS

\[-{1 \over \rho} \left( {{\partial^2 p} \over {\partial x^2}} + {{\partial^2 p} \over {\partial y^2}} \right)+ \nu \left( {{\partial^3 u} \over {\partial x^3}} + {{\partial^3 u} \over {\partial x \partial y^2}} + {{\partial^3 v} \over {\partial x^2 \partial y}} + {{\partial^3 v} \over {\partial y^3}} \right) = LHS\]

Re-arrange:

\[-{1 \over \rho} \left( {{\partial^2 p} \over {\partial x^2}} + {{\partial^2 p} \over {\partial y^2}} \right)+ \nu \left( {{\partial^2} \over {\partial x^2}} \left( {{\partial u} \over {\partial x}} + {{\partial v} \over {\partial y}} \right) + {{\partial^2} \over {\partial y^2}} \left( {{\partial u} \over {\partial x}} + {{\partial v} \over {\partial y}} \right) \right) = LHS\]

Apply Continuity, so \(\nabla \cdot \mathbf{u} = 0\). Hence:

\[-{1 \over \rho} \left( {{\partial^2 p} \over {\partial x^2}} + {{\partial^2 p} \over {\partial y^2}} \right) = LHS\]

1.1.3. The Poisson Equation in Vector Notation

1.1.3.1. Equate LHS and RHS

\[-{1 \over \rho} \left( {{\partial^2 p} \over {\partial x^2}} + {{\partial^2 p} \over {\partial y^2}} \right) = \left( {{\partial u} \over {\partial x}} \right)^2 + 2 {{\partial u} \over {\partial y}} {{\partial v} \over {\partial x}}+ \left( {{\partial v} \over {\partial y}} \right)^2\]

1.1.3.2. In Vector Form

\[\nabla^2 p = -f\]
  • A Poisson Equation for pressure, which ensures that continuity is satisfied.
  • Now pressure and velocity are coupled in the continuous domain.

1.2. Poisson’s Equation in Numerical Domain

We have shown that the Poisson’s Equation is valid in the continuous domain, but in the numerical domain, we use discretisation

1.2.1. Momentum Equations in Vector Form

1.2.1.1. The Momentum Equation for an incompressible fluid:

\[{{\partial \mathbf{u}} \over {\partial t}} + \mathbf{u} \cdot \nabla \mathbf{u} = {-{1 \over \rho} \nabla p} + {\nu \nabla^2 \mathbf{u}}\]

1.2.2. Discretised Momentum Equations in Vector Form

Discretise in time: FD in time, with pressure at time \(n+1\) (the pressure that corresponds with the velocity at \(n+1\))

\[\mathbf{u}^{n+1} = \mathbf{u}^n + \Delta t \left(-\mathbf{u}^n \cdot \nabla \mathbf{u}^n - {1 \over \rho} \nabla p^{n+1} + \nu \nabla^2 \mathbf{u}^n \right)\]

1.2.3. The Divergence of the Momentum Equations in Vector Form

Now take the divergence of the momentum equations, if \(\mathbf{M}\) is the “vector” of momentum equations, the divergence is:

\[\nabla \cdot \mathbf{u}^{n+1} = \nabla \cdot \mathbf{u}^n + \Delta t \left(-\nabla \cdot (\mathbf{u}^n \cdot \nabla \mathbf{u}^n) - {1 \over \rho} \nabla^2 p^{n+1} + \nu \nabla^2 (\nabla \cdot \mathbf{u}^n) \right)\]
  • In the numerical scheme, we want a divergence free velocity at the next step, i.e. \(\nabla \cdot \mathbf{u}^{n+1} = 0\) to satisfy continuity
  • But at the current step we may have \(\nabla \cdot \mathbf{u}^{n} \ne 0\) due to numerical errors
  • So we can’t cancel out the divergence of velocity in the numerical domain (although we could in the continuous domain)

This is a Fractional Step approach:

  • Solve the momentum equation for velocity, numerically (the divergence of the velocity might not be zero) \(\mathbf{u}^{n+1/2}\)
  • Solve the Poisson equation for pressure, forcing the divergence to be zero
  • Correct the velocity to satisfy continuity

1.2.4. The Poisson Equation in Vector Notation

Poisson Equation for \(p\) at time \(n+1\) and forcing \(\nabla \cdot \mathbf{u}^{n+1} = 0\)

\[\nabla^2 p^{n+1} = \rho {{\nabla \cdot \mathbf{u}^n} \over {\Delta t}}- \rho \nabla \cdot (\mathbf{u}^n \cdot \nabla \mathbf{u}^n)+ \mu \nabla^2 (\nabla \cdot \mathbf{u}^n)\]
  • In the numerical domain the velocity field we are producing with the Navier Stokes equations is not completely divergence free
  • Think of the velocity obtained from the Navier Stokes as being at an intermediate step \(\mathbf{u}^{n+1/2}\)
  • And \(\nabla \cdot \mathbf{u}^{n+1/2} \ne 0\)
  • We need \(p^{n+1}\) so that continuity is satisfied

«  4. von Neumann Stability Analysis   ::   Contents   ::   2. The Marker and Cell Method  »