The WENO reconstruction in the Godunov method for modeling hydrodynamic flows with shock waves

In the paper, we describe one simple WENO reconstruction for the Godunov method that allows one to get a version of the low-dissipation method. We describe in detail the procedure for reconstructing the ”right” and ”left” values of physical variables, which are used as arguments for an exact solution of the Riemann problem in the Godunov method. In the discontinuity decay test, we verify the quality of the developed numerical method and study its accuracy order. As a model task, the problem of multiple explosion of white dwarfs as a result of their high-speed collision in a three-dimensional problem statement will be considered.


Introduction
The Godunov method has been known for over 60 years and it is successfully used to solve hydrodynamic problems with discontinuous solutions and shock waves. After the original method had been developed many modifications of the Godunov method, aimed to reduce its numerical dissipation, have been proposed. The main modifications are MUSCL-like schemes [1,2] and Kolgan solver [3], piecewise-parabolic method [4] and its compact implementation [5,6] with extension to operator splitting approach [7], equations of magnetohydrodynamics [8] and relativistic hydrodynamics [9].
The main idea of all modifications is the use of piecewise polynomial representations of the solution. It brought with the development in the form of the WENO schemes [10,11,12,13] with piecewise cubic representation of the solution [14,15] up to the 17th order of precision [16]. Note that schemes like "Harten -Lax -van Leer" [17] or Lax-Wendroff [18] had been as the base solver. In this paper, we propose a new version of the modern implementation of the Godunov scheme [19] using the WENO reconstruction to obtain the low-dissipation property of the numerical method.
In the second section, we describe in detail the construction of the modern version of the Godunov method and the reconstruction of the solution to obtain the low-dissipation property. The third section deals with the verification of the numerical method on the Sod problem and the study of the convergence order of the developed method. The fourth section demonstrates the application of the developed numerical method to simulate multiple explosions of white dwarfs during their high-speed collision in the form of the Iax supernova explosion. The fifth section is a conclusion.

Numerical Method
Let us consider the equations of hydrodynamics in a one-dimensional case: where ρ is a density, u is a normal velocity, p is a pressure, γ is an adiabatic exponent , c = γp ρ is a sound speed.
In the computational domain, we introduce a uniform grid with a spatial step h. The time step τ is computed from the Courant condition: where CF L is Courant-Friedrichs-Lewy number. Also, to describe the numerical scheme, we introduce the value of the total mechanical energy: The values of conservative (ρ, ρu, ε) and physical (ρ, u, p) variables are defined at the cell centers with a half-integer index. The values of flows through the boundaries of cells (ρu, ρu 2 +p, [ε+p]u) are defined at nodes with an integer index. The scheme of Godunov is written as: where R n i , U n i , P n i is a solution of the discontinuity decay problem or the Riemann problem. Dirichlet-type conditions are used as boundary conditions. Let us consider the solution of the linearized decay of the discontinuity for two neighboring cells (the index L denotes the left cell, the index R stands for the right cell) in accordance with the work [19]. In case of a supersonic flow on the left, when the condition u L > c L is satisfied, the solution of the Riemann problem is: In the case of a supersonic flow on the right, that is, when the condition u R < −c R is satisfied, the solution to the Riemann problem is: and the following condition is satisfied on it: The speed of the right wave is computed by the formula: and there is the next condition on it: As a result, the values of the velocity U and pressure P at the boundary of the cells are computed by the equations: In the paper [19] there had been proposed the original approach based on the constant value of p − ρc 2 at the characteristics dx/dt = u ± c: As a result, the sign of the velocity U depends on the solution for the density R at the boundary of the cells, calculated by the equation: Using the above formulas, fluxes across the boundary of the corresponding physical variables are computed in the Godunov method. For WENO, the reconstructions of the physical variables ρ, u and p (denoted by the function f ) on the interface i are computed by the following equations: To compute the coefficients ω L,i and ω R,i , the following equations are used: where σ L,i and σ R,i are calculated as: where = 10 −36 has been used. Values of β L,i and β R,i are computed by the equations: After the WENO reconstruction, we solve the Riemann problem for the values: Note that the choice of some coefficients was based on computational experiments. To estimate the convergence of this method, we study the behavior of the L 1 norm errors:

Verification
where u(x i ) is the exact solution at x i , u i is the numerical result, and h is spacing of an uniform grid. The behavior of the L 1 norm for the Sod problem can be seen from the table (1). From the table (1) it can be seen that the convergence for the density function drops almost to its half value, then increases and becomes of the convergence order of ∼ 0.6. For pressure and velocity functions, the behavior of the convergence order is similar. Such behavior of the convergence order for a discontinuous solution took place for the classical Godunov scheme, also.

Astrophysics Simulation
As a model astrophysical problem, let us consider the multiple expolsion of white dwarfs of the solar mass as a result of their high-speed collision in the three-dimensional formulation [20]. The mathematical model of the evolution of white dwarfs is based on solving the overdetermined system of gravitational hydrodynamics [21]. To close the gravitational hydrodynamics, the adaptation of the stellar equation of state [22] is used. It consists of the contribution of the pressure of a nondegenerate hot gas, pressure due to radiation and a degenerate gas. In the case of a degenerate gas, the relativistic and nonrelativistic regimes are considered. As a net of nuclear systems, we consider the α-network [23]. In the computational experiment, the temperature of the dwarfs reached the value of T = 10 8 K. At the distance of 200 km from the point of explosion, many satellite bubbles had appeared. The figure (2) shows the density isosurface at the time of t = 5 seconds. It can be seen from the figure (2) that the combustion fronts are correctly  reproduced due to the subsonic turbulent combustion of carbon. These results confirm the conclusions that ignition and transition to the detonation combustion are not required to obtain sufficiently powerful explosions. In the mathematical model, we use the ultimate adiabatic model of the state equation for a degenerate gas, which limits our possibilities for a more realistic account of the physics of the explosion in terms of chemical composition. However, the state equation applied in this work and allowing to describe enough the hydrodynamics of the evolution of white dwarfs and the supernova explosions is also widely used.

Conclusion
A simple low-dissipation WENO reconstruction of the Godunov scheme has been proposed. In the paper there has been described in detail the procedure for reconstructing the "right" and "left" values of physical variables, which are used as arguments for the exact solution of the Riemann problem in the Godunov method. The numerical method has been verified on the Sod problem both in terms of the solution quality and accuracy order of the scheme. The supernova explosion of Iax type based on a high-speed collision of white dwarfs has been used as an astrophysical application to verify the numerical method.