Pressure

Question: what is the difference between “incompressible” and “constant density?”

Compressible flows are high speed relative to the sound speed.
When changes in fluid velocity are high enough to significantly impact fluid properties, we say a flow is compressible.
From the Bernoulli equation we have
$$\Delta P = \frac{\rho \Delta v^2}{2}$$
Using $c^2 = \gamma P/\rho$, and letting $\Delta v = v$, we have
$$\frac{\Delta P}{P} = \frac{\gamma}{2}Ma^2,$$
where $\gamma = c_p/c_v$.
For v=100 m/s, $\Delta P/P_{atm}$ is 6%.
- 100 m/s is 233 mph and gives $Ma=0.288$ for air at 1 atm.
For v=347 m/s, $\Delta P/P_{atm}$ is 70%.
- v=347 m/s is 776 mph and gives $Ma=1$ for air at 1 atm.

We can have a variable density flow that is considered incompressible.

Density variations may occur due to, e.g. temperature or species composition changes.
Incompressible for $Ma \le 0.3$ (approximately).

Low Ma flows have pressure fluctuations due to the flow that do not significantly affect the flow.

Flow solvers: stable stepsize

Unfortunately, stable explicit numerical solvers have step sizes set by the smallest timescales, even if the corresponding physical process doesn’t significantly affect the flow.

We would like our stepsize to be set by an advective CFL proportional to $\Delta x/u$.
For stability, a compressible solver requires an acoustic CFL proportional to $\Delta x/(u+c)$.
If Ma = u/c is small, then we need many fewer timesteps for stability than we need for accuratly representing the flow. This defines a stiff problem.
The number of stepsizes needed for stability relative to the number needed for accuracy is $(\Delta x/u)/(\Delta x/(u+c)) = 1 + 1/Ma$.

Pressure gradient scaling

For the Euler equations, the wave speeds are $u$, $u-c$, and $u+c$.
In following the derivation of these, $c$ can be replaced with $\sqrt{\beta}c$, where $\beta$ is a constant, if we multiply the pressure gradient in the momentum equation by $\beta$.
- This is equivalent to adjusting the gas constant $R$.
- The gasdyn code illustrates this.
  - The stable timestep size varies approximately as $\sqrt{\beta}$. In the code, $\Delta t\approx 0.088\,\mu s$ for $\beta=1$, and $\Delta t\approx 1.2\,\mu s$ for $\beta=0.005$. This timestep size ratio is 0.0733, which is practically the same as $\sqrt{0.005}=0.0707$.
In this case the relative number of timesteps required compared to that for advection only is $1 + \sqrt{\beta}/Ma$.
- The relative change in pressure is $\gamma Ma^2/(2\beta)$.
If Ma = 0.1, $\beta = 1.0$ the relative steps are 11 and $\Delta P/P = 0.7\%$.
If Ma = 0.1, $\beta = 0.1$ the relative steps are 4 and $\Delta P/P = 7\%$.
Hence, scaling $\beta$ has allowed fewer timesteps, but has artificially increased pressure fluctuations. The pressure fluctuations increase faster than the timesteps decrease.

Pressure Projections

Low Mach flows can also be solved using pressure projection approaches.
References:
- “Computational Methods for Fluid Dynamics,” J. H. Ferziger, M. Peric, 3rd edition, Springer, 2002
- “Numerical Heat Transfer and Fluid Flow”, S. V. Patankar, Hemisphere Pub., 1980.
- “Numerical Heat Transfer and Fluid Flow”, S. V. Patankar, Hemisphere Pub., 1980.
- Saad et al. (2018)

Constant density, unsteady

Following Ferziger

We have the continuity and momentum equations, which are two equations in pressure and velocity:

$$\nabla\cdot v = 0,$$

$$\frac{\partial v}{\partial t} = -\frac{\partial v v}{\partial x} - \nabla\cdot\tau - \frac{1}{\rho}\nabla P.$$

The momentum equation is an equation for velocity.
There is no explicit pressure equation.
The continuity equation provides a constraint on velocity.
The pressure is a scalar field that will be set so that the continuity equation is satisfied. (Note that the pressure field is only specified up to an additive constant, since only $\nabla P$ appears.)

Apply a simple Explicit Euler timestep to momentum:

$$v^{n+1} = \underbrace{v^n + \Delta t\left(-\frac{\partial v v}{\partial x} - \nabla\cdot \tau \right)}_{H^n} - \frac{\Delta t}{\rho}\nabla P.$$

Take the divergence of this (which is zero), and rearrange to get the pressure equation:

$$\nabla^2P = \frac{\rho}{\Delta t}\nabla\cdot H^n.$$

To advance a given step, we solve this equation for P. And then use this in the momentum equation for $v$, which will satisfy continuity by construction.

$$v^{n+1} = H^n - \frac{\Delta t}{\rho}\nabla P.$$

Variable density, unsteady

Saad et al. (2018) extend the treatment above to the case of variable density flows, where the key difference is the treatment of the nonzero velocity divergence.

Steady treatments

For steay flows (as in RANS), variations of the SIMPLE algorithm are common. See Patankar and Ferziger.