NASA Contnetor Report 178042

lCASE REPORT NO. 86-3 NASA-CR-I78042

19860011681

leASE

A RIEMANN SOLVER BASED ON A GLOBAL EXISTENCE PROOF

FOR THE RIEMANN PROBLEM

Pravir Dutt .-- -.- . '- :.> ,~-:-, ..,. .;.... "':" .. .,'1(

I, • ,

.-,,,'" .'

Contract Nos. NASl-17070 and NASl-18107

January 1986

'. -\ • _ _.- •• 1.

INSTITUTE FOR CO~WUTER APPLICATIONS IN SCIENCE AND ENGINEERING

NASA Langley Research Center, Hampton, Virginia 23665

Operated by the Universities Space Research Association

NI\SI\

National Aeronautics and

Space Administration

Langley Research Center

Hampton, Virginia 23665

A RIEHANN SOLVER BASED ON A GLOBAL EXISTENCE PROOF

FOR THE RIEMANN PROBLEM

Pravir Dutt

Institute for Computer Applications in Science and Engineering

ABSTRACT

Godunov's method and several other methods for computing solutions to the

equations of gas dynamics use Riemann solvers to resolve discontinuities at

the interface between cells. A new method is proposed here for solving the

Riemann problem based on a global existence proof for the solution to the

Riemann problem. The method is found to be very reliable and computationally

efficient.

Research was supported by the National Aeronautics and Space Administration

under NASA Contract Nos. NASI-17070 and NAS1-18107 while the author was in

residence at ICASE, NASA Langley Research Center, Hampton, VA 23665-5225.

i

I. INTRODUCTION

The introduction of a finite-difference scheme based on Riemann solvers by

Godunov permitted the development of robust numerical methods for problems

involving strong shocks. He proposed to consider all quantities in a

computational cell as given by cell averages and to resolve the resulting

discontinuities at the cell edges by the solution to a Riemann problem. The

method has the important advantage that it automatically handles strong shocks

and interactions and is able to predict cavitation should it occur.

The use of Godunov's method has attracted particular interest in the field

of aerodynamics. Here one expects to encounter only ideal gas (polytropic)

equations of state. A Riemann solver for the equations of gas dynamics was

first proposed by Godunov. It was subsequently modified and improved by

Chorin [2].

In this paper we take advantage of a _lobal existence proof for the

Riemann problem to propose an altogether different method for solving the

Riemann problem. This method has the advantage that it is computationally

very efficient and is guaranteed to converge to the solution of the Riemann

problem.

2. THE RIEMANN PROBLEM IN GAS DYNAMICS

The gas dynamics equations in Eulerian coordinates are

Pt + (PU)x= 0

(PU)t + (pu2 + P)x = 0
-2-

e t + ((e + p)U)x = 0

where p is the density of the gas, u is the velocity, pu is the momentum,

e is the energy per unit volume, and p is the pressure.

We assume that the gas is ideal and polytropic so that the equation of

state is given by

p = pRT

where p is the pressure, p the density, T the temperature, and R the

universal gas constant. Further, we assume that y > i, where y is the

ratio of specific heats.

We seek a solution to the gas dynamics equations with initial data

H[ I = for x < 0LPLJ

and

= for x > O.

Lp j

The solution to the Riemann problem is illustrated in Figure 1 on a

x - t graph.
-3-

2-contactdiscontinuity

/

# u /

/

p_ / .

/ PR'u,

/

/ P _ckor

/ rarefaction /

PL'UL' /

/

PL PR'UR'

PR

:_:='-XFigure |: Solution to the Riemann problem.

The initial discontinuity is resolved by a system of waves consisting of a

l-shock or rarefaction, a 2-contact discontinuity, and a 3-shock or rarefaction. If cavitation is excluded the pressure and velocity on either side

of the contact discontinuity are equal and are given by p and u .

3. THE GLOBAL EXISTENCETHEOREM

A parametric representation of the solution to the Riemann problem and

global existence proof are discussed in Smoller [3]. We outline the

discussion to relate it to our numerical scheme. Let

B =Y + 1

y - 1

and

"f- 1

T - 27 •
-4-

We denote the sound speed by c, where c = (YP)I/2 . Thus cL is the

sound speed for the state (PL,UL,PL). The solution to the Riemann problem is

obtained by the following three one-parameter family of curves.

1 - family, for x I E ]P_

* -x1

P---= e

PL

-Xl/Y

, e , Xl>_0

PL

PL fl (Xl) - x 1

B+e

xI ' x1<0

1+Be

I _zx 1

u* - y 2 1 (i - e ), Xl_> 0 uL

CL = hi (Xl) - /o -Xl

2(_) i0= 1 - e < 0

--X x 1

_ _ [ (i + Be 1) I/2

2 - family, for x2 €

OR x2

-- = e

PL
-5-

3 - family, for x3 E

PR x3

p, - e

OR 1

, - f3(x3) : f1(x3)

PR

2 Tx3

* y - 1 (e - I), x3 _>0

uR - u

, = h3(x3) -

cR

x3

2(T) I/2 (e - I)

x ' x3 <0.

y - 1 I1 + 8e 3)1/2

To obtain a solution of the Riemann problem, we need to solve the above

system of equations for the parameters Xl, x2, and x3. It should be noted

that x I < 0 corresponds to a 1-shock and x I > 0 to a 1-rarefaction wave.

A similar statement holds for the parameter x3. Using the above representation the following theorem can be proved.

Theorem: Let (OL,UL,PL) an___d (PR,UR,PR) be any two states (not

necessarily close). Then there is a unique solution to the Riemann problem

with these initial states if and only if

2

uR - uL <-y- 1 (cL + cR).

If this condition is violated, then a vacuum is present in the solution.
-6-

4. THE NUMERICAL METHOD

We first check whether

2

uR - uL < _y - I (cL + CR)"

If the condition is violated, then cavitation occurs.

To obtain the values Xl, x2, and x3 explicitly, we define

0R

APL

PR

B =

PL

uR - uL

C -

cL

We first solve for xl, where x I satisfies the equation:

g(x1) m hl(Xl) + (_)1_hl(Xl + log B) - C = 0.

It should be noted that the function g is a monotone function. Hence, the

equation can be solved efficiently by a numerical routine. We first locate an

interval [a,b] such that g(a) < 0 and g(b) > 0. With a and b as our

initial guesses, a few iterations of the regula-falsi method are needed to

obtain an accurate value for xI.

We then have

x3 = x I + log B.
-7-

Finally, we obtain x 2 from the equation:

x2

fl(xl)e f3(x3) = A.

Having obtained x I, x2, and x3 we substitute them into the parametric

representaion for the solution to the Riemann problem and obtain the values

u ' P ' PL' and PR"

5. NUMERICAL RESULTS

We found that over a range of y and on a variety of shock tube problems

our Riemann solver converged to the exact solution within about five

iterations of the regula-falsi method.

6. CONCLUSIONS

A numerical method for solving the Riemann problem in gas dynamics has

been proposed. The method is based on an existence proof for the solution to

the Riemann problem. Its advantages are that it can predict the occurrence of

cavitation and it is computationally very efficient.
-8-

APPENDIX

Computer Program for the Riemann Solver

subroutine riemann

common/input/ul,ur,rhol,rhor,pl,pr

common/output/pstar,ustar

common/parm/gamma,beta,tau,a,b,c

c input:ul,ur--velocltles to the left and right of interface

c rhol,rhor--densities

c pl,pr--pressures

c output:pstar--interface pressure

c ustar--interface velocity

beta=(gamma+l.)/(gamma-l.)

tau=(gamma-l.)/(2.*gamma)

cl=sqrt(gamma*pl/rhol)

cr=sqrt(gamma*pr/rhor)

c check for cavltalton

if((ur-ul).gt.(2.*(cl+cr)/(gamma-1.)))go to 3

a=rhor/rhol

b=pr/pl

c=(ur-ul)/cl

x=0.

c locate interval on which g(x) has a zero

9 if((g(x).gt.0.).or.(g(x+l.).it.0.))then

if(g(x).gt.0.)then

x_x-l.

else

x=x+l.

end if

end if

if((g(x).gt.0.).or.(g(x+l.).it.0.))go to 9

c do regula-falsi iterations

xmax=x+l.

xmln=x

n=l

8 x=xmax-g(xmax)*(xmax-xmin)/(g(xmax)-g(xmln))

if(g(x).it.0.) then

xmin=x

else

xmax=x

end if

n=n+]

if((abs(xmax-xmin).gt.0.0001).and.(n.le.10)) go to 8

xl=xmax

c calculate ustar,pstar

pstar=exp(-l.*xl)*pl

ustar=ul+cl*h1(xl)

return
-9-

3 print *,'cavitation occurs"

return

end

function hl(x)

common/parm/gamma,beta,tau,a,b,c

if(x.gt.O.)then

hl=(2./(gamma-1.))*(1.-exp(-1.*tau*x))

else

h1=2.*sqrt(tau)*(l.-exp(-1.*x))/(gamma-l.)

h1=hl/sqrt(1.+beta*exp(-1.*x))

end if

return

end

function g(x)

external hl

common/parm/gamma,beta,tau,a,b,c

y=x+alog(b)

g=h1(x)+sqrt(b/a)*hl(y)-c

return

end
-10-

_FE_N_S

1. J. K. DUKOWICZ, J. Comput. Phys. 61 (1985), 119.

2. A. J. CHORIN, J. Comput. Phys. 22 (1976), 517.

3. J. SMOLLER, Shock Waves and Reaction-Diffusion Equations, Springer-Verlag,

New York, 1983.



1. Report No. NASA CR-178042 12. Government Accession No. 3. Recipient's Catalog No.

ICASE Report No. 86-3

4. Title and Subtitle 5. Report Date

A RIEMANN SOLVER BASED ON A GLOBAL EXISTENCE January 1986

PROOF FOR THE RIEMANN PROBLEM 6. Performing Organization Code

7. Author(s) 8'ffff!.fjrming Organization Report No. Pravir Dutt

9. Performing Organization Name and Address 10. Work Unit No.

Institute for Computer Applications in Science

and Engineering 11. Contract or Grant No.

Mail Stop 132C, NASA Langley Research Center NASI-17070; NASI-18107 U",~~"nn \7" ,)'l("(.." ."')'),, 13. Type of Report and Period Covered 12. SpOlisoring Agency Name and Address Contractor Report

National Aeronautics and Space Administration 14. Sponsoring Agency Code Washington, D.C. 20546

505-31-83-01

15. Supplementary Notes

Langley Technical Monitor: Submitted to J. Comput. Phys.

J. C. South

Final Report

16. Abstract

Godunov's method and several other methods for computing solutions to the

equations of gas dynamics use Riemann solvers to resolve discontinuities at

the interface between cells. A new method is proposed here for solving the

Riemann problem based on a global existence proof for the solution to the

Riemann problem. The method is found to be very reliable and computationally

efficient.

17. !{ey Words (Suggested by Authors(s)) 18. Distribution Statement

Riemann solver 59 - Mathematics & Computer Science

(General)

6!. - Numerical Analysis

Unclassified - unlimited

10. Security CIa.9sif.(of this report) 120. Security CIa.9sif.(of this page) 21. No. ~f2pages122. P~i~; Unclassified Unclassified

For sale by the National Technical Information Service, Springfield, Virginia 22161

NASA·Langley, HJ8fi