**4. FOHS for Navier-Stokes equations**

The FOHS method is straightforward for model equations by introducing the gradients of the unknowns as auxiliary variables. However, for compressible NS equations, the construction of the FOHS becomes trickier. A HNS14 method is first introduced by including viscous stresses and heat fluxed as additional variables [19], but it turns out to obtain reduced order of accuracy in velocity gradients. Later, a HNS 17 method is developed to improve the scheme by including the velocity gradients scaled by viscosity and heat fluxes [20], resulting in the designed order of accuracy in velocity gradients. However, the loss of designed accuracy in density gradients is observed. In order to overcome this issue, an artificial hyperbolic density diffusion term is added to HNS17 to construct HNS20 [12]. The additional term is designed to be small enough such that the scheme can provide expected accurate gradients for density while the continuity equation is not affected. At this point, accurate gradients can be obtained for all primary variables by HNS20 in the pseudo-steady state. However, this efficient construction is not straightforward for other conservative or primary variables. To make the recycling procedure trivial, a new formulation HNS20G [21] was developed recently. Unlike the above-mentioned methods, HNS20G uses the gradients of the primary variables, that is, density, velocity, and temperature, as auxiliary variables to the first-order hyperbolic system. Moreover, with the utilization of the reconstruction method, HNS20G is able to deliver a more accurate solution and *Hyperbolic Navier-Stokes with Reconstructed Discontinuous Galerkin Method DOI: http://dx.doi.org/10.5772/intechopen.109605*

gradient while remaining the same degrees of freedom as the conventional DG (P1) method.

#### **4.1 Navier-Stokes equations**

The non-dimensionalized Navier–Stokes equations governing unsteady compressible viscous flows can be expressed as

$$\begin{cases} \frac{\partial \rho}{\partial t} + \frac{\partial \rho v\_k}{\partial \mathbf{x}\_k} = \mathbf{0}, \\\\ \frac{\rho v\_i}{\partial t} + \frac{\partial (\rho v\_i v\_k + p \delta\_{ik} - \tau\_{ik})}{\partial \mathbf{x}\_k} = \mathbf{0}, \\\\ \frac{\rho \epsilon}{\partial t} + \frac{\partial \left[ v\_j \rho \epsilon + v\_i (p \delta\_{ik} - \tau\_{ik}) + q\_k \right]}{\partial \mathbf{x}\_k} = \mathbf{0}, \end{cases} \tag{23}$$

where the summation convention (*k* = 1, 2, 3) has been used as *x*<sup>1</sup> ¼ *x*, *x*<sup>2</sup> ¼ *y*, and *x*<sup>3</sup> ¼ *z*. The symbols *ρ*, *p,* and *e* denote the density, pressure, and specific total energy of the fluid, respectively, and *vx* ¼ *u*, *vy* ¼ *v*, and *vz* ¼ *w* are the velocity components of the flow in the coordinate direction *x*, *y,* and *z, respectively*. The symbol *δik* is the Kronecker delta function. The pressure can be computed from the equation of state for a perfect gas,

$$p = (\boldsymbol{\chi} - \mathbf{1})\rho \left(\boldsymbol{e} - \frac{\mathbf{1}}{2} (\boldsymbol{u}^2 + \boldsymbol{v}^2 + \boldsymbol{w}^2)\right),\tag{24}$$

where the ratio of the specific heats, is assumed to be constant, that is, *γ* ¼ 1*:*4. Moreover, the specific total enthalpy *h* is defined as

$$h = e + \frac{p}{\rho}.\tag{25}$$

The viscous stress tensor *τ* can be found through

$$
\pi = \begin{pmatrix}
\pi\_{\text{xx}} & \pi\_{\text{xy}} & \pi\_{\text{xx}} \\
\pi\_{\text{yx}} & \pi\_{\text{yy}} & \pi\_{\text{yx}} \\
\pi\_{\text{xx}} & \pi\_{\text{xy}} & \pi\_{\text{xx}}
\end{pmatrix}.
\tag{26}
$$

The Newtonian fluid with the Stokes hypothesis is valid under the current framework. Thus, *τ* is a symmetric tensor, which is a linear function of the velocity gradients

$$
\pi\_{\vec{\eta}} = \mu \left( \frac{\partial u\_i}{\partial \mathbf{x}\_j} + \frac{\partial u\_j}{\partial \mathbf{x}\_i} \right) - \frac{2}{3} \mu \frac{\partial u\_k}{\partial \mathbf{x}\_k} \mathfrak{G}\_{\vec{\eta}}, \tag{27}
$$

where *μ* is the molecular viscosity coefficient or dynamic viscosity as in many literatures, which can be determined through Sutherland's law

$$\frac{\mu}{\mu\_0} = \left(\frac{T}{T\_0}\right)^{\frac{3}{2}} \frac{T\_0 + \mathcal{S}}{T + \mathcal{S}},\tag{28}$$

where *μ*<sup>0</sup> represents the viscosity coefficient at the reference temperature *T*<sup>0</sup> and *S* ¼ 110*K*. The temperature of the fluid *T* is determined by

$$T = \frac{P}{\rho \mathcal{R}},\tag{29}$$

where *R* denotes the universal gas constant for a perfect gas. According to Fourier's law, the heat flux vector *qj* is given by

$$q\_j = -\lambda \frac{\partial T}{\partial \mathbf{x}\_j},\tag{30}$$

where *λ* is the thermal conductivity coefficient and expressed as

$$
\lambda = \frac{\mu c\_p}{\text{Pr}},
\tag{31}
$$

where *cp* is the specific heat capacity at constant pressure and Pr is the nondimensional laminar Prandtl number, which is taken as 0.7 for air.

#### **4.2 First-order hyperbolic Navier–Stokes system: HNS20G**

The gradients of density, velocity, and temperature are first introduced as auxiliary variables

$$\begin{aligned} \mathbf{r} &= \begin{pmatrix} r\_x \\ r\_\mathcal{I} \\ r\_\mathcal{I} \\ r\_x \end{pmatrix} = \nabla \rho, \quad \mathbf{g}\_u = \begin{pmatrix} \mathbf{g}\_{ux} \\ \mathbf{g}\_{uy} \\ \mathbf{g}\_{ux} \end{pmatrix} = \nabla u, \quad \mathbf{g}\_v = \begin{pmatrix} \mathbf{g}\_{vx} \\ \mathbf{g}\_{v\mathcal{I}} \\ \mathbf{g}\_{vx} \end{pmatrix} = \nabla v, \\\ \mathbf{g}\_{uv} &= \begin{pmatrix} \mathbf{g}\_{ux} \\ \mathbf{g}\_{uv} \\ \mathbf{g}\_{uv} \end{pmatrix} = \nabla w, \quad \mathbf{h} = \begin{pmatrix} h\_x \\ h\_\mathcal{I} \\ h\_\mathcal{I} \\ h\_z \end{pmatrix} = \nabla T. \end{aligned} \tag{32}$$

Then, the governing equations of compressible viscous flows are reformulated as

$$\begin{cases} \frac{\partial \rho}{\partial t} + \frac{\partial \rho v\_k}{\partial x\_k} - \nu\_r r\_k = 0, \\ \frac{\rho v\_i}{\partial t} + \frac{\partial (\rho v\_i v\_k + p \delta\_{ik} - \tau\_{ik})}{\partial x\_k} = 0, \\ \frac{\rho e}{dt} + \frac{\partial \left[v\_j \rho e + v\_i (p \delta\_{ik} - \tau\_{ik}) + q\_k\right]}{\partial x\_k} = 0, \\ -\frac{\partial v\_i}{\partial x\_k} = -g\_{v\_i k}, \\ -\frac{\partial T}{\partial x\_k} = -h\_k, \\ -\frac{\partial \rho}{\partial x\_k} = -r\_k, \end{cases} \tag{33}$$

*Hyperbolic Navier-Stokes with Reconstructed Discontinuous Galerkin Method DOI: http://dx.doi.org/10.5772/intechopen.109605*

where *ν<sup>r</sup>* is introduced as an additional term to the continuity equation as artificial viscosity associated with the artificial hyperbolic mass diffusion. In the presented work, we take *<sup>ν</sup><sup>r</sup>* <sup>¼</sup> min 10�12, *<sup>h</sup>*<sup>3</sup> � �, where *<sup>h</sup>* is the scale of cell size. A pseudo time *<sup>τ</sup>* is then introduced into Eq. (33) to make the system hyperbolic, and the resulting hyperbolic system is written as

$$\begin{cases} \frac{\partial \rho}{\partial \tau} + \frac{\partial \rho}{\partial t} + \frac{\partial \rho v\_k}{\partial \mathbf{x}\_k} = \mathbf{0}, \\\\ \frac{\partial \rho v\_i}{\partial \tau} + \frac{\partial \rho v\_i}{\partial t} + \frac{\partial (\rho v\_i v\_k + p \delta\_{ik} - \tau\_{ik})}{\partial \mathbf{x}\_k} = \mathbf{0}, \\\\ \frac{\partial \rho \epsilon}{\partial \tau} + \frac{\partial \rho \epsilon}{\partial t} + \frac{\partial \left[ v\_i \rho \epsilon + v\_i (p \delta\_{ik} - \tau\_{ik}) + q\_k \right]}{\partial \mathbf{x}\_k} = \mathbf{0}, \\\\ T\_v \frac{\partial \mathbf{z}\_{vk}}{\partial \tau} - \frac{\partial v\_i}{\partial \mathbf{x}\_k} = -\mathbf{g}\_{v,k}, \\\\ T\_h \frac{\partial h\_k}{\partial \tau} - \frac{\partial T}{\partial \mathbf{x}\_k} = -h\_k, \\ T\_r \frac{\partial r\_k}{\partial \tau} - \frac{\partial \rho}{\partial \mathbf{x}\_k} = -r\_k. \end{cases} (34)$$

The relaxation times *Tr*, *Tv*, and *Th* are defined as

$$T\_r = \frac{L\_r^2}{\nu\_r}, \quad T\_v = \frac{L\_v^2}{\nu\_v}, \quad T\_h = \frac{L\_h^2}{\nu\_h}, \tag{35}$$

where the length scale is taken as *Lr* ¼ *Lv* ¼ *Lh* ¼ 1*=*2*π*. The normalized viscosity coefficients are

$$
\nu\_v = \frac{\mu\_v}{\rho}, \quad \nu\_h = \frac{\mu\_h}{\rho}, \tag{36}
$$

Note that for high Reynolds number flows, the length scale should be modified according to the Reynolds number as

$$\begin{aligned} L\_r &= \frac{L\_d}{\sqrt{Re \frac{r^{\infty}}{L\_d}}}, \quad Re^{r\alpha}\_{L\_d} = \frac{U\_{\infty} L\_d}{\nu\_r}, \\ L\_v &= \frac{L\_d}{\sqrt{Re \frac{r^{\infty}}{L\_d}}}, \quad Re^{r\alpha}\_{L\_d} = \frac{U\_{\infty} L\_d}{\nu\_{\nu \alpha}}, \\ L\_h &= \frac{L\_d}{\sqrt{Re \frac{h^{\infty}}{L\_d}}}, \quad Re^{h\alpha}\_{L\_d} = \frac{U\_{\infty} L\_d}{\nu\_{\nu \alpha}}, \\ L\_d &= \frac{1}{2\pi}, \end{aligned} \tag{37}$$

where *U*<sup>∞</sup> is the free-stream flow velocity, and the superscript and subscript ∞ refers to the free-stream flow.

After writing (Eq. (34)) in the same formulation as (Eq. (12)), either FV formulation (Eq. (2)), DG formulation (Eq. (5)), or rDG formulation (Eq. (7)) can be used to integrate the first-order spatial operators.

#### **4.3 Laminar flow past a sphere**

The laminar flow past a sphere is given here to compare the present method with experimental data. The free-stream Reynolds number is taken as *Re* <sup>∞</sup> ¼ 100, and the free-stream Mach number is taken as *M*<sup>∞</sup> ¼ 0*:*5. The triangular meshes for the symmetric plane and the spherical surface are shown in **Figures 2** and **3**. Characteristic condition is given at the left, right, upper, and lower boundaries in **Figure 2**. Adiabatic

**Figure 2.** *Triangular meshes of the symmetric plane for laminar past a sphere.*

**Figure 3.** *Triangular meshes of the sphere surface.*

*Hyperbolic Navier-Stokes with Reconstructed Discontinuous Galerkin Method DOI: http://dx.doi.org/10.5772/intechopen.109605*

wall boundary is given at the spherical surface. The vorticity magnitude contours are given in **Figure 4**. The pressure coefficient distribution on the spherical surface is compared to Ref. [22], as shown in **Figure 5**.
