AUSM⁺-up Flux Scheme

The wellbore simulator uses the AUSM⁺-up (Advection Upstream Splitting Method) flux scheme for computing numerical fluxes at cell interfaces. Based on Liou (2006) with area-weighted extensions from Sacconi & Mahgerefteh (2020).

Interface Mach Number

The numerical flux at each cell interface is determined by a split Mach number. Left and right states each contribute through the $\mathcal{M}_4^{\pm}$ polynomials. A pressure-difference dissipation term scaled by $K_\mathrm{p}$ improves stability at low Mach numbers.

$$M_{i+\frac{1}{2}} = \mathcal{M}_4^{+} (M_\mathrm{L}) \;+\; \mathcal{M}_4^{-} (M_\mathrm{R}) \;-\; K_\mathrm{p} \, \max\bigl(1 - \bar{M}^2,\; 0\bigr) \, \frac{p_\mathrm{R}^{\,\mathrm{f}} - p_\mathrm{L}^{\,\mathrm{f}}}{\bar{\rho}\,\bar{a}^2} \tag{1}$$

Well-Balanced Face Pressures

In a static wellbore column with gravity, neighbouring cells have different pressures due to the hydrostatic gradient. Without correction, the $K_\mathrm{p}$ dissipation term sees this difference and generates spurious mass flux.

The well-balanced modification extrapolates cell-centre pressures to the shared face using the hydrostatic increment $\Delta P = \tfrac{1}{2}(\rho_\mathrm{L} + \rho_\mathrm{R})\,g\,\Delta z\,\sin\theta$:

$$p_\mathrm{L}^{\,\mathrm{f}} = p_\mathrm{L} + \frac{\rho_\mathrm{L}}{\rho_\mathrm{L} + \rho_\mathrm{R}}\,\Delta P_{\text{hydro}}, \qquad p_\mathrm{R}^{\,\mathrm{f}} = p_\mathrm{R} - \frac{\rho_\mathrm{R}}{\rho_\mathrm{L} + \rho_\mathrm{R}}\,\Delta P_{\text{hydro}} \tag{2}$$

At hydrostatic equilibrium $p_\mathrm{R} - p_\mathrm{L} = \Delta P$ exactly, so both face values converge to the same pressure and the correction vanishes.

Mach Splitting Polynomials

The fourth-order Mach polynomials smoothly split the interface Mach number into left- and right-running contributions. For supersonic flow ($|M| \geq 1$) they reduce to simple upwinding. The subsonic branch uses $\beta$ to control blending:

$$\mathcal{M}_4^{\pm}(M) = \begin{cases} \dfrac{1}{2}(M \pm |M|) & |M| \geq 1 \\[8pt] \pm\,\dfrac{1}{4}(M \pm 1)^2 \Bigl[1 + 16\,\beta\;\dfrac{1}{4}(M \mp 1)^2\Bigr] & |M| < 1 \end{cases} \tag{3}$$

Interface Pressure

The interface pressure is reconstructed using $\mathcal{P}_5^{\pm}$ polynomials with a velocity-difference dissipation term scaled by $K_\mathrm{u}$:

$$p_{i+\frac{1}{2}} = \mathcal{P}_5^{+}\,p_\mathrm{L}^{\,\mathrm{f}} \;+\; \mathcal{P}_5^{-}\,p_\mathrm{R}^{\,\mathrm{f}} \;-\; K_\mathrm{u}\;\mathcal{P}_5^{+}\,\mathcal{P}_5^{-}\,(\rho_\mathrm{L} + \rho_\mathrm{R})\,\bar{a}\;\frac{A_\mathrm{R}\,u_\mathrm{R} - A_\mathrm{L}\,u_\mathrm{L}}{\max(A_\mathrm{L},\,A_\mathrm{R})} \tag{4}$$

Pressure Splitting Polynomials

$$\mathcal{P}_5^{\pm}(M) = \begin{cases} \dfrac{\mathcal{M}_1^{\pm}}{M} & |M| \geq 1 \\[8pt] \pm\,\mathcal{M}_2^{\pm} \Bigl(2 \mp M - 16\,\alpha\;M\,\mathcal{M}_2^{\mp}\Bigr) & |M| < 1 \end{cases} \tag{5}$$

Mass Flux

The final convective mass flux uses upwind selection based on the interface Mach number:

$$\dot{m}_{i+\frac{1}{2}} = \bar{a}\left[\frac{(\rho A)_\mathrm{L}}{2}\bigl(M_{i+\frac{1}{2}} + |M_{i+\frac{1}{2}}|\bigr) \;+\; \frac{(\rho A)_\mathrm{R}}{2}\bigl(M_{i+\frac{1}{2}} - |M_{i+\frac{1}{2}}|\bigr)\right] \tag{6}$$

References

• Liou, M.-S. (2006). A sequel to AUSM, Part II: AUSM⁺-up for all speeds. Journal of Computational Physics, 214:137–170.
• Sacconi, D. & Mahgerefteh, H. (2020). Wellbore flow modelling for CO₂ injection.