Finite Volume Method for the Relativistic Burgers Model on a (1+1)-Dimensional de Sitter Spacetime

Several generalizations of the relativistic models of Burgers equations have recently been established and developed on different spacetime geometries. In this work, we take into account the de Sitter spacetime geometry, introduce our relativistic model by a technique based on the vanishing pressure Euler equations of relativistic compressible fluids on a (1+1)-dimensional background and construct a second order Godunov type finite volume scheme to examine numerical experiments within an analysis of the cosmological constant. Numerical results demonstrate the efficiency of the method for solutions containing shock and rarefaction waves.


Introduction
A fundamental relativistic Burgers equation has recently been derived and generalized to different spacetime geometries by LeFloch and his collaborators [1,2,3,4]. The first model has been obtained by establishing a hyperbolic balance law satisfying the Lorentz invariance property for flat geometry and then its relativistic generalizations were extended to the Schwarzschild, Friedmann-Lemaitre-Robertson-Walker (FLRW), de Sitter and Schwarzschild-de Sitter spacetimes [5,1,2,3,4].
In this article, we consider the "de Sitter" (dS) spacetime which is a member of the family of the FLRW geometry. The line elements of the FLRW and the dS geometries share certain common properties. In particular, these metrics are solutions to the Einstein field equations. Furthermore, the Minkowski metric can be derived from particular cases of both the FLRW and the dS metrics. The inspiration of derivation of the relativistic Burgers equation on the dS spacetimes is based on the common properties of these geometries. We follow the instruction via [1,4,6] and take into account the relativistic Euler equations on a given curved background M by the following general formulation of balance laws: T αβ (ρ, u) = ρc 2 u α u β + p(ρ) u α u β + g αβ , (1.1) where ∇ is the covariant derivative operator, T αβ is the energy-momentum tensor of perfect fluids, ρ ≥ 0 denotes the mass-energy density of the fluid, c is the light speed and the unit timelike vector field u α represents the velocity of the fluid so that g αβ u α u β = −1.

Motivation of the Paper
One of the basic nonlinear hyperbolic conservation laws: is the (inviscid) Burgers equation which has various applications in many physical fields. Indeed, this equation is a simplified version of the Euler equations of compressible fluids ∂ρ ∂t + ∂ ∂x (ρv) = 0, ∂ ∂t (ρv) + ∂ ∂x (ρv 2 + p(ρ)) = 0, (1.3) where v denotes the velocity, ρ the density, and p the pressure of the fluid. We recover the inviscid Burgers equation from this system by first imposing that p ≡ 0 in Equation (1.3). Then, we take partial derivatives for both equations, combine them and finally deduce the Equation (1.2). This analysis was generalized for the relativistic models in the articles [1,2,4] in order to get the relativistic Burgers equations on the related geometries.
In the present article, we apply this analysis by considering the relativistic Euler equations on a smooth, time-oriented, curved (1+1)-dimensional dS spacetime in order to derive the relativistic Burgers model. We pursue and cite the article by Ceylan and Okutmustur [2] for the details of derivation of the model on a dS background. Once we obtain our model, we then construct a second order Godunov type finite volume scheme in order to inspect the stability and efficiency of the method. Following the papers by LeFloch et al. [1,4], the relativistic Burgers equations are derived from Equation (1.1) on a (1+1)-dimensional background by imposing a vanishing pressure to Equation (1.1). The derived equations satisfy the Lorentz invariance property and in the case where c → ∞, the classical (non-relativistic) Burgers equation is recovered directly from the given system. In paper [4], the curved spacetime is considered to be flat and then generalized to a Schwarzschild background; whereas in [1], an FLRW background is taken into account. In both [1,4], the derived relativistic models are used for the numerical experiments via a second-order accurate finite volume method based on the geometric formulation Equation (1.1). A similar work is attained by Ceylan and Okutmustur for the Schwarzschild-de Sitter geometry [3]. For the convergence and geometric formulation of finite volume schemes on Lorentzian (curved) manifolds, we follow the articles by LeFloch et al. [7,16]. The theory of conservation laws on the flat space is introduced by Kruzkov [8]. We refer the reader to the following work by Groah, Smoller and Temple [9] for the shock wave interactions and further details in general relativity. On the other hand, we cite a recent work on balance laws posed on (1+1)-dimensional spacetimes by Gosse [10] where a locally inertial Godunov scheme is presented numerically. Moreover, further details of second order Godunov schemes are presented in the book of Guinot [11] and the article by Van Leer [12]. The reader can find a more general instruction for the theory of general relativity and related concepts about Einstein's theory in the book of Wald [13].
An outline of this paper is as follows. We introduce the basic features of a dS geometry and its metric in the first part. We then derive the Euler equations by means of the Christoffel symbol where we refer the reader to the article [2] for further details. The next issue is deriving the particular cases of the relativistic Burgers equations depending on the cosmological constant parameter Λ. Then, we present a geometric formulation of finite volume approximation on a curved background which follows by a construction of finite volume scheme on local coordinates. We then establish a Godunov type second order scheme for the given model. We also compare the non-relativistic and relativistic Burgers equations depending on different values of Λ. The effects of the Λ on the numerical scheme are studied in the numerical experiments part. We finish the article by investigating the efficiency of the scheme with different Λ values. The results demonstrate that the scheme is consistent with the conservative form of our model which yields correct computations of weak solutions containing shock and rarefaction waves.

The dS Metric and the Euler System
The dS spacetime is a particular case of the Lorentzian manifold and its metric provides a cosmological solution to the Einstein field equations similar to the FLRW metric. The line element (metric) of the dS background involves a constant so called "cosmological constant" and it is denoted by Λ. Einstein firstly introduced this cosmological constant in his field equations. Depending on the sign of the cosmological constant Λ, we distinguish the de Sitter and the Anti-de Sitter spacetimes; namely, if Λ > 0, the geometry is called the "de Sitter spacetime" and, if Λ < 0, it is called the "Anti-de Sitter spacetime". Particularly, if Λ = 0, the metric turns to be a Minkowski metric and thus we have a flat geometry.
The corresponding line element for a (3+1)-dimensional case in terms of time t, the radial r and angular coordinates θ and ϕ is given by: Thus, Equation (2.1) is called the "de Sitter metric" if Λ > 0 and the "Antide Sitter metric" if Λ < 0. Notice that Equation (2.1) can be rewritten in matrix form as with the non-zero "covariant" elements g 00 = −(1 − Λr 2 ), g 11 = 1 1−Λr 2 g 22 = r 2 , g 33 = r 2 sin 2 θ and the inverse of (g ij ) is with the corresponding "contravariant" elements Here, is satisfied where δ i j is the Kronecker's delta function.

Energy-Momentum Tensor
We take into account our spacetime to be (1+1)-dimensional so that the solutions to the Euler system depend only on the time variable t and radial variable r and the angular components θ, ϕ vanish. It follows that and where (u 0 ) 2 and (u 1 ) 2 are unitary vectors. Next, we substitute the covariant terms g 00 and g 11 in the above relation which yields Next, we define a velocity component v depending on u 0 and u 1 as follows: Examining the relations (2.4) and (2.5) together, it follows that It remains to substitute the Equation (2.6) terms into the energy momentum tensor for perfect fluids equation, which is defined by Thus, by using Equations (2.6) and (2.7), the tensors are evaluated as follows:

The Euler System on a (1 + 1)-Dimensional dS Background with Zero Pressure
We are about to write the Euler equations and deduce the presureless Euler system from these equations. We recall that the Euler system (1.1) is given by which is equivalent to where T αβ terms are calculated in the previous subsection. In order to derive the desired system, we start by substituting β = 0 in Equation (2.8) to get Then, we put α, γ ∈ {0, 1, 2, 3} into the above equation which yields and putting α, γ ∈ {0, 1, 2, 3} into the above equation yields Replacing the Christoffel symbols and the tensor terms with their calculated values and imposing zero pressure to this system yields which is the "Euler system with zero pressure" on a (1+1)-dimensional dS background.

Derivation of the Model
We observed in the introduction part that the classical Burgers equation is deduced from the Euler equations by imposing zero pressure to the given system. In this section, we follow the same technique in order to derive the relativistic Burgers equation on the dS geometry. The investigation of static and spatially homogeneous solutions of the given model is also an issue of this section. We start by rewriting the first and second equations of (2.11) with the notation ∂ 0 = ∂ t , ∂ 1 = ∂ r and combining these two equations (by taking partial derivatives of each term of (2.11) and summing these equations after elimination of the vanishing terms) in order to derive the following single equation: which is the desired so called "relativistic Burgers equation on the dS background". Notice that the left-hand side of the model can also be written in the conservative form as follows:

Limiting Case of the Derived Model
The Equation (3.1) may differ depending on different values of Λ. We recall that the particular case of the dS metrics for the value Λ = 0 gives the Minkowski metric. If we substitute Λ = 0 in the relativistic Burgers equation (3.1), we recover the classical (inviscid) Burgers equation: In other words, the limiting case of the relativistic Burgers Equation (3.1) on a dS spacetime yields the classical (non-relativistic) Burgers equation. This common property is also shared by the other relativistic models on the curved spacetimes (such as flat, Schwarzschild and FLRW).

Static and Spatially Homogeneous Solutions for the Model
In this subsection, we search for the static and homogeneous solutions of our model if they exist. To examine this, we need to write the conservative form of the derived model (3.1), which is We start by seeking the static solutions, that is the t-independent solutions, to Equation (3.2). Due to the t-independency, the term ∂ t v vanishes, and we obtain Applying the change of variable where N ∈ (0, c) is a constant parameter. It follows that Hence, we find the "static solutions" described by On the other hand, we follow an analogue process for the homogeneous solutions, that is the r-dependent solutions to the Equation (3.2). It is obvious to conclude that there is no spatially homogeneous (r-independent) solution since the equation depends already on r.
As a result, we conclude that our model (3.1) admits only "static solutions".

Finite Volume Method on a Curved Spacetime
The finite volume method is an important discretization technique for partial differential equations. In the present work, we use a finite volume approximation for general balance laws of hyperbolic partial differential equations on an (n+1)-dimensional manifold following the papers [7,14,4,15,16]. We consider the geometry in (n+1)-dimensional case, where n refers to the space and 1 refers to the time dimension and introduce a geometric formulation of the finite volume method for the general balance laws. We then consider the finite volume approximation on coordinates to work on the particular (1+1)-dimensional spacetime for analyzing numerically the relativistic Burgers Equation (3.1).
In the following, the introduction of a geometric formulation of the finite volume method on a curved spacetime and its particular case for the local coordinates are given.

Geometric Spacetime Finite Volume Method
Following the article [7], we search for a formulation of the finite volume scheme for a hyperbolic balance law given by posed on an (n + 1)-dimensional curved spacetime M , where the unknown function v is a scalar field, div(·) is the divergence operator, T (v) is the flux vector field and S(v) is the scalar field.
We are now ready to introduce the geometric formulation of the finite method for discretizing the balance laws (4.1) on M . First of all, we establish a general triangulation of the spacetimes by for the manifold M such that M is composed of these spacetime elements K ⊂ M satisfying the following conditions: • The boundary of an element K, denoted by ∂K, is given by which is piecewise smooth and includes exactly two "spacelike" faces (having an induced Riemannian type metric) denoted by e + K and e − K , and "timelike" elements (having an induced Lorentzian type metric) denoted by e 0 ∈ ∂ 0 K := ∂K \ e + K , e − K .
• The intersection K ∩ K of two distinct elements K and K in T h is either a common face of K, K or is a smooth submanifold with dimension at most (n − 1).
• Along the timelike faces e + K , e − K , we introduce the outgoing unit normal vector field denoted by n e + K , n e − K , respectively.
Furthermore, we assume that M permits a foliation {H t } t∈(0,∞) , by oriented spacelike hypersurfaces such that the parameter t : M ∪∂M → [0, +∞) provides us with a global time function. This allows us to separate between future-oriented (t increasing) and past-oriented (t decreasing) timelike directions on M . Thus, M is assumed to be foliated by hypersurfaces H t (which are spacelike elements) such that with the time function t if it is determined from a sequence of discrete times t 0 < t 1 < t 2 < · · · < t N such that all spacelike faces are submanifolds of H t = H tn , (n = 0, 1, · · · , N ) and H 0 is an initial slice. It follows that the class of nonlinear hyperbolic Equations (4.1) and (4.2) gives a scalar model on which we can analyze numerical methods. We are ready to define the finite volume approximations by averaging the balanced law (4.1) over each element K ∈ T h of the triangulation. By integrating in space and time, we can write Applying the Stokes formula to the above equation, it follows that Using the given averages v − K and v − K e 0 computed along e − K and e 0 ∈ ∂ 0 K, we need to compute the average v + K along e + K . To this aim, we introduce and where, for each K and e 0 ∈ ∂ 0 K, we selected a numerical flux Q K,e 0 : R 2 → R satisfying natural properties of consistency, conservation, and monotonicity. Therefore, the finite volume method of interest takes the form Notice that a standard CourantFriedrichsLewy (CFL) condition is assumed for the sake of stability of this scheme. For further details, the reader may consult [7,14,15,16].

Finite Volume Schemes in Coordinates
Following the geometric formulation of the finite volume scheme on (n + 1)dimensional spacetime in the previous part, we now take into account the particular case (for n = 1) and assume that the spacetime is described in coordinates (t, r). We divide into equally spaced cells I j = [r j−1/2 , r j+1/2 ] of size ∆r, centered at r j , that is, r j+1/2 = r j−1/2 + ∆r, satisfying r j−1/2 = j∆r, r j = (j + 1/2)∆r.
Furthermore, we denote by ∆t the constant time length and we set t n = n∆t.
In order to introduce the finite volume method in local coordinates, we start by rewriting the hyperbolic balance law (4.1) in (1+1) dimension, that is, ∂ t T 0 (t, r) + ∂ r T 1 (t, r) = S(t, r), (4.5) where T 0 , T 1 are flux fields, and S is the source term. We integrate Equation (4.5) over each grid cell which is equivalent to The numerical flux functions are approximated by S(t, r) dt dr.
Using these approximations in the above integrations, the scheme becomes Since T is convex (see [4]), by taking the inverse of the above relation, the scheme (4.6) is rewritten by v n+1 For further details of triangulations and discretization of finite volume approximations on a curved spacetime, we refer the reader to the papers by LeFloch et al. [7,4,6,16].

Second Order Godunov-Type Scheme
In this section, we construct a second order Godunov type finite volume scheme for the derived model which is the main objective of this paper. Fundamentally, Godunov type algorithms could be summarized by the following steps: • Discretize the cells in finite volumes; • Obtain the profile reconstruction on the cells; • Do specification of the Riemann problems at cell interfaces; • Find solution for the Riemann problems; • Perform computation for fluxes on cell interfaces; • Determine the variables for the next time step.
The second order scheme differs from the original Godunov type scheme in reconstruction on the profile and specification of the Riemann problems. More generally, for the equation a first order scheme can be converted to a second order method by proceeding the cell-boundary values which are used in the numerical flux functions to decide the intermediate time level t n+1/2 = (t n + t n+1 )/2. In other words, the second order Godunov scheme is obtained from the edge-boundary values of the reconstructed profile proceeded by a half time step.
The construction of the second order type of scheme is formulated by the following steps: • Firstly, the variable within the computational cells is reconstructed which gives couples of values (v n i,L , v n i,R ) in each computational cell where the letters L and R refer to left and right sides, respectively. Note that v n i,L lies between v n i−1 and v n i , whereas v n i,R lies between v n i and v n i+1 .
• Secondly, one should advance the solution by half a step in time. The intermediate values v at the cell edges at the time t n+1/2 = (t n + t n+1 )/2 are denoted by (v n i,L , v n i,R ). These values can be calculated by the following: • As a third step, we solve the Riemann problem formed by the intermediate values (v n i,L , v n i,R ). The solution v n+1/2 i+1/2 is used in order to compute the flux f • In the final step, we advance the solution by the time step ∆t from t n using the usual formula We refer to the article by Van Leer [12] and the book by Guinot [11] for construction of the algorithm to the second order Godunov type finite volume method.

Implementation of the Godunov Scheme for the Derived Model
In this part, numerical tests are displayed for the model derived on a dS spacetime based on a second order Godunov scheme. We prepare the second order finite volume approximation for our model based on the construction in the previous section. The method for the intermediate values and the proceeding solutions are formulated as follows. We analyze the given model with a single shock and rarefaction for an initial function considering the Godunov scheme with a local Riemann problem for each grid cell depending on two particular cases of cosmological constant Λ. Note that, for the numerical experiments, we fix r ∈ [0, 1]. We take into consideration the fastest wave at each grid cell since both shock and rarefaction waves are produced in the Riemann problem. We compel transmissive boundary conditions on the scheme. By normalization (taking c = 1) in Equation (3.1), we attain the following model: and the finite volume schemes for this model can be given as Notice that the speed term is and, in order for the stability condition in the scheme to be satisfied, we choose ∆t and ∆r so that the CFL condition satisfied, i.e., ∆t ∆r max The implementation of the second order Godunov scheme is based on the construction given above. Observe that Λ > 0 and Λ = 0 refer to the dS and Minkowski geometries, respectively. We compare the particular cases Λ = 1 > 0 and Λ = 0 for the schemes of the model with shock and rarefaction waves; that is, we analyze the model via the second order Godunov scheme for two different geometries.
The results are sketched in Figures 6.1 and 6.1. The rarefactions are drawn in Figure 6.1 and shocks are given in Figure 6.1 for two particular cases of Λ. The iterations are performed for different numbers of iterations; for instance, for the first figure of Figure 1-2, the number of iterations n = 100, for the second figure n = 400, for the third one n = 600, and for the last one n = 800. Moreover, the CFL condition (6.6) is strictly less than 1 in each test for stability of the scheme. We observe that the numerical solution for the case Λ = 0 (represented by the red curve) moves faster than the case Λ = 1 (represented by the blue curve). In other words, we have a relatively higher speed when the spacetime is flat comparing with the dS geometry. This fact can also be verified numerically by substituting Λ = 0 and Λ = 1 into the speed term given by (6.5), i.e., by the term 1 − Λ(r n+1/2 j±1/2 ) 2 , which feeds the theoretical background of the model with the numerical results.
We also observe that the negative values of Λ yield the velocity to become greater than 1 by Equation (6.5), which is not acceptable since they exceed the speed of light (as we assumed that c = 1 for the numerical calculations). In other words, the sign of the cosmological constant Λ is significant for the numerical results by which the geometry under consideration is identified.

Conclusions
A new nonlinear hyperbolic model, which describes the propagation and interactions of shock waves on the dS and flat spacetimes, is derived in this study. We take into account the relativistic Euler system on a curved background and impose a zero pressure in the statement of the energymomentum tensor for perfect fluids in order to derive our model (3.1). This model involves a cosmological constant Λ, which can be normalized by taking the positive, negative and zero values. According to the positive, zero and negative values of Λ, the geometry of interest becomes the de Sitter, Minkowski and Anti-de Sitter spacetimes, respectively. We have the following remarkable observations: • The standard Burgers equation is recovered when Λ is chosen to vanish; • Various mathematical properties concerning the hyperbolicity, genuine nonlinearity, shock waves, and rarefaction waves are established; • Investigation for the class of static and spatially homogenous solutions is carried out; • Shock wave solutions to the model for the different possible values of the cosmological constant Λ are investigated.
The analysis for two possible values of the cosmological constant Λ relies on a numerical scheme, which applies discontinuous solutions and is based on the finite volume technique. We establish a second-order Godunov scheme for our model which gives the main originality to the current paper.
• The scheme is consistent with the conservative form of our model, which yields correct computations of weak solutions containing shock waves.
• The numerical solution curves for the Anti-de Sitter spacetimes yields the velocity to become greater than one, which is not acceptable as it exceeds the speed of light.
• One of the most obvious findings emerging from this study is that it allows us to make a numerical comparison of the relativistic models on the de Sitter and Anti-de Sitter spacetimes by means of finite volume approximations. Indeed, the method is not efficient for the relativistic model on the Anti-de Sitter spacetimes, whereas it is efficient for the de Sitter spacetime.
• The efficiency and convergency of the numerical method depends on the value of the cosmological constant Λ. The sign of Λ is significant. If Λ is smaller than zero, the method is not efficient and the solution curves do not converge. On the other hand, if Λ = 0 and Λ > 0, the method gives efficient and convergent numerical results.
• The numerical experiments illustrate the convergence, efficiency and robustness of the proposed scheme on the de Sitter background.
To sum up, we highlight that the existence of the model of the concerning equation that permits one to improve and analyze numerical methods for shock capturing schemes and to attain definite conclusions concerning the convergence, robustness and efficiency of the schemes. According to the background geometry, different techniques may be needed to guarantee that certain classes of time-dependent or space-dependent solutions for the model be preserved by the scheme. A strategic perspective can be given as to investigate relativistic Burgers equations on different backgrounds with more complicated points.