Print preview

1D Navier-Stokes equations - Stationary shock wave

Problem description

The exact solution of a stationary shock wave for the one dimensional Navier-stokes equations has been used as an AMROC verification problem. In a shock stationary frame of reference, the 1D Navier-Stokes equations for constant viscosity and thermal conductivity are

$\frac{\partial( \rho u)}{\partial x} = 0$

$\frac{\partial( \rho u)}{\partial x} = -\frac{\partial p}{\partial x}} + \frac{4}{3} \frac{\partial}{\partial x} \left( \mu \frac{\partial u}{\partial x} \right)$

$\frac{\partial( \rho e_t u)}{\partial x} = -\frac{\partial (p u)}{\partial x}} + \frac{\partial}{\partial x}\left( u \frac{4}{3} \mu \frac{\partial u}{\partial x} + k \frac{\partial T}{\partial x}\right).$

which is supplemented by the ideal gas law,

$T =\frac{p}{R \rho}.$

This solution has a Prandtl number equal to 3/4, which changes the energy equation to an algebraic equation:

$h + \frac{u^2}{2} = h_0 + \frac{u_0^2}{2},$

using the constraints,

$Pr = 3/4 = \frac{c_p \mu}{k}, \quad c_p = \frac{\gamma R_u}{(\gamma-1) W} , \quad k = \frac{4 \mu c_p}{3},$

with the following non dimensionalization,

$\bar{u} = \frac{u}{u_0}, \quad \bar{\rho} = \frac{\rho}{\rho_0}, \quad \bar{p} = \frac{p}{p_0}, \quad \xi = \frac{x}{\lambda_0}.$

The analytical solution is formulated in non-dimensional form, with the upstream density, pressure, and velocity in addition to an equivalent perfect gas mean free path are used to non-dimensionalize the solution. The mean free path is quite abritrary, and for the solution used is

$\lambda_0 = \frac{8 \mu}{ 5} \sqrt{\frac{2}{\pi \rho_0 p_0}}$

All the above simplifies the 3 PDE's of the 1D Navier-Stokes equations to one ODE, which is solved analytically. See Kramer et al[1] for the implicit solution, which is expressed as a function of the Mach number and specific heat ratio and relates the non-dimensional velocity and position. The density, pressure, and hence the total energy are found with the relations,

$\rho u = \rho_0 u_0$

$p/p_0 = \frac{-M^2 ((\gamma-1)\bar{u}^2-(\gamma-1)))-2)}{2 \bar{u}}.$

This derivation is carried out for for the monoatomic gas case, $\gamma = 5/3$ , with different arbitrary perfect gas equivalent mean free path in the textbook by Zel'dovich and Raizer[3].

Initial / Boundary Conditions

This specific exact solution is available in the weno and clawpack: applications/euler/1D/ViscousShock directories. For the following tests, we use the computational domain, x = [-30,30], inflow/outflow boundary conditions and the parameters

$Pr=3/4,\quad \mu=1, \quad R=1, \quad R_{u}=1, \quad \gamma = 1.4, \quad M=2, \quad \rho_0=1, \quad p_0 =1, \quad t_{final}=0.5$

where the 0 subsript denotes the upstream flow with Mach number equal to two. The domain is large enough and the final time short enough, to minimize the boundary errors.

CFL number definition:

By using this exact solution as a verification problem an equivalent CFL number was derived. "CFL number" in this context, is the Navier-Stokes CFL number (including both convective and diffusive operators and analogous to the Euler equations case) that is used in a CFL condition. This CFL condition garrantees a stable solution for which the CFL number is less than or equal to one.

Results:

A stable solution and the order of accuracy for both Clawpack's 2nd order slope-limited Roe solvers and WENO was verified, with viscous and heat conduction terms for a perfect gas. In the process, the differences in the calculations of the CFL number were found in the new Navier-Stokes modification for Clawpack with that of the current WENO-TCD implementation. In actuality, the CFL number for stability of an Explicit method is very complicated for the full one and three dimensional Navier-Stokes equations. However, using a condition which is based on a one-dimensional advection-diffusion equation [2] (which virtually adds the convective, diffusive, and thermal conditions together rather than taking a maximum over all) is robust and stable up to CFL = 1, . The current WENO implementation is stable and non-ossillatory for this verification problem only up to CFL=~0.6, as it used a less robust CFL number.

The newly implemented CFL number evaluation for Clawpack in 1D is

$CFL = max \left( \frac{2 \Delta t}{\Delta x^2}*max\left(\frac{4 \mu}{3 \rho},\frac{k}{\rho c_v} \right) + max\left(u + c\right)*\frac{\Delta t}{\Delta x} \right),$

where the maximum is take over all cell centers. This is derived by using Von Neumann analysis and ignoring the $u \frac{\partial^2 u}{\partial x^2}$ term.

[1] R.M.J. Kramer, C. Patano, D.I. Pullin. A class of energy stable, high-order finite-difference interface schemes suitable for adaptive mesh refinement of hyperbolic problems. J. Comput. Phys., 226:1458-1484, 2007.

[2] A.C. Hindmarsh, P.M. Gresho, and D.F. Griffiths. The stability of explicit Euler time-integration for certain finite difference approximations of the multi-dimensional advection-diffusion equation. Inter. J. Num. Methods in Fluids, 4:853-897, 1984.

[3] Y.B. Zel'dovich and Y.P. Raizer. Physics of shock waves and high-temperature hydrodynamic phenomena. Dover. 2002.

-- JackZiegler^? - 29 Jul 2008

End of topic
Skip to action links | Back to top

Attachment	Size	Date	Who	Comment
velocity.ps	40.5 K	29 Jul 2008 - 12:54	JackZiegler^?	velocity
velocity.png	4.2 K	29 Jul 2008 - 13:12	JackZiegler^?	velocity
entropy.gif	4.1 K	29 Jul 2008 - 13:19	JackZiegler^?	entropy

You are here: Amroc > ViscousShockWaveConvergenceStudy > ViscousShock

to top