An evaluation of cell type finite difference methods for solving viscous flow problems

Computers & Fluids, Vol. 1, pp. 3 - 1 . Pergamon Press, 1973. Printed in Great Britain.

AN EVALUATION OF CELL TYPE FINITE DIFFERENCE METHODS FOR SOLVING VISCOUS FLOW PROBLEMS T. D. TAYLOR Fluid Mechanics Department, The Aerospace Corporation, E1 Segundo, Los Angeles, Calif. 90045, U.S.A.

(Received 6 October 1972)

Abstract--Cell (integral) type equations are examined and compared with grid point techniques for I-D compressible flow problems. An analysis of Burger's equation for weak wave propagation is made in order to establish the nature of the cell scheme for both constant cell properties and for linear variations of cell properties. The results of the analysis yield insight into the behavior of the Godunov cell method for inviscid flows as well as the Lax-Wendroff technique. This understanding is employed to extend Godunov's method to viscous flow calculations. Tests of the method for the complete 1-D compressible viscous flow equations are conducted.

INTRODUCTION

The rapid progress in the development of high speed computers in the past 10 yr has made numerical methods for solving complicated transport problems practical. Two principal developments which are of long range importance have occurred. The first is application of numerical calculations for unsteady transport problems to obtain steady state results. This approach has been applied primarily in the solution of fluid flow problems. Examples of this technique appear in refs. [1-18]. The second and possibly the most important advancement in the numerical methods field is the introduction of the method of splitting by Bagrinovskii and Godunov[19] and Yanenko[20]. This scheme permits reduction of an "N" dimensional unsteady transport problem to "N"--substeps composed of one-dimensional unsteady transport problems. As a result, methods which have been developed for solving one-dimensional unsteady problems can be applied to complicated multi-dimensional problems. The combination of the unsteady approach and the method of splitting has made solution of many complicated problems feasible. The principal limitation, however, lies in the differencing method for solving the 1-D unsteady transport problem. The differencing methods for the partial differential equations describing transport processes take two forms. The first is a grid point approach in which the difference equations are derived from truncated Taylor's series expansions of the differential equations. The second method is the cell scheme in which the difference equations are derived from the integral conservation laws. Schemes based on the integral laws tend to express average values of the unknowns over a control volume or cell rather than a specific value at a point. This averaging property of cell schemes often permits their use in situations where a grid scheme may encounter difficulty; for example, where corners may introduce multi-valued functions. Cell schemes, however, have not been investigated

4

T . D . TAYLOR

to the point where their relative merits can be assessed. In this paper, an attempt is made to improve the present understanding of cell schemes and their relation to the grid point approach. In order to quantitatively evaluate cell techniques, it is necessary to examine problems in which analytic solutions to the flow equations can be constructed. For compressible viscous flow, the available analytic solutions appear to be limited to one-dimensional wave propagation problems. As a result, the analysis in this paper is limited to this class of problems. In this study, the weak wave propagation as described by Burger's equation and the strong wave propagation described by the more general 1-D viscous flow equations are considered. The one-dimensional limitation of the study prohibits one from drawing sweeping conclusions regarding multi-dimensional problems. If, however, one can apply the method of splitting to reduce a multi-dimensional problem to one-dimensional unsteady operators then the results can provide direct insight into the behavior of a cell differencing method.

FUNDAMENTAL

APPROACH

The principal idea of cell finite difference approximations is to consider the flow region to be divided into a number of finite elements or cells. The conservation equations are then applied to each of these finite cells. The result is a set of equations which relate the mean flow quantities in the cell to the fluxes of mass, momentum or energy at the boundaries. In order to demonstrate this idea, consider the one-dimensional propagation of a weak wave as described by Burger's equation which has the form ~u gu 1 gZu 8t + u d~ = Re ~x 2

(i)

where u is the velocity and Re denotes the parameter simulating the Reynolds number per unit length, defined as Uo/6, in which u0 is the characteristic velocity and ~ is the diffusivity of sound. This equation has the same type of nonlinearity as the more general viscous flow equations and in limited situations described the excess velocity of a weak wave propagating into a constant state region.J21] Burger's equation along with its exact solution will now be used to attempt to gain information regarding the cell difference method. If Burger's equation is integrated over a space interval Ax = x 2 - xl and a time interval At = t2 - tl the result is ut~';> - ut<~> + u"22 - u~2 - q~2 - q"<'>

At

2Ax

(2)

Re. Ax

In this equation, ft2,,2 dt SXfu d x uZ _ .,~-

u --

Ax

At

and

fr2,7 dt q _ jt~-~

At

where q = Ou/Ox. As the definitions indicate superscripts have been used to designate the time and subscripts the position.

An evaluation of cell type finite difference methods for solving viscous flow problems

5

Equation (2) is a difference equation for the mean velocity, u
n÷l

I1

%x> n-I I I I

I

IL

× Fig. 1. First order cell velocity distribution.

The next order of approximation to the intial conditions is a linear fit between cell centers. This approach permits discontinuous gradients but u itself remains continuous. Both the jump conditions between cells and the linear approximation have been employed to develop difference equivalents to Burger's equation. It is important to note that these difference approximations are exact in the special cases where the initial distribution u(x, O) is a jump or linear function. As a r+esult, they should reflect information on the proper approach to differencing a nonlinear equation with viscous and convective effects.

6

T.D. TAYLOR

FIRST ORDER RESULTS The results obtained for an initial step function fit to u(x, O) will now be discussed. The boundary fluxes for the jump approximations to the initial condition are derived by employing Lighthill's solution to Burger's equation[21] subject to the conditions

u=u. forx < 0 a t t = 0 U=Un+X f o r x > 0 a t t = 0 . The solution has the form U--

u,A + u,+ IB A+B

where

I

orfcF-U'q

A=exp ~

[ 2x/~j

[ 2x/6t _ff

_

U. +

Un + x

2

This equation can be integrated and differentiated accordingly to obtain u and q at cell boundaries are then determined analytically from the solution by setting x = 0 as the cell boundary. The details of the analysis will not be given due to the complexity, but the resulting difference equations which describe the essential features are (1) For small Reynolds numbers based on cell quantities, i.e. (u,. Ax)/6 ~ 1, the difference approximation of equation (1) in conservation form becomes

T

7 n)-

2

I

In'+' l

Ax~/-~ght

[

Ax ~ - t -

(3)

or in convective form

At

+

4

[

2Ax

]

(4)

where the superscript t 2 denotes the new time level and no superscript denotes the old time level. Note that in this equation it has been assumed that u._ 1 and u.+ 1 are the same order of magnitude as u.. Their exact magnitude is not important as long as the Reynolds number formed by u._ xAx/6 or u.+ lax/6 remains less than one. (2) For the limit of large Reynolds numbers (u.. Ax)/6 >>1, the difference equation has the conservation form

An evaluation of cell type finite differencemethods for solving viscousflowproblems

7

or in convective form

Un'--Un [U,+U,IlI,--u, ll

A ~

+

2

Ax

= Ax [u,+~. At

u,.At J'

(6)

These results show that the cell difference scheme depends on the magnitude of the transport coefficient. For small Reynolds numbers, the convective terms appear as a product of a central difference and a characteristic velocity formed by a three point average. In the high Reynolds number case, these terms "switch" to forward differences (donor cell) with a two point average characteristic velocity. The diffusion terms also vary in both cases. For the small Reynolds number case, the term is characteristic of flows in which diffusion plays a significant role since the characteristic length dividing the second difference is the diffusive length x/(6. At). In the high Reynolds number case, note that the characteristic length becomes a convective length (ui. At) which is dependent on the solution at each cell and time step. It is important, however, to note that in order for the first order cell difference equation to remain consistent with the differential equation some serious constraints must hold. For small Reynolds numbers, note that the limit of Axx/(g6 At) must approach one and for large Reynolds numbers limit Ax/u At must approach one. These constraints arise due to the presence of the viscous terms. As a consequence, one is lead to the conclusion that for calculating the viscous terms the jump approximation is not the most practical approach. If, however, one were only to approximate the convective term u(gu/gx) the approximation appears to be useful and in fact it provides information on how the differencing should be adjusted to stabilize the calculation. In addition, it suggests that possibly the first order cell approach could be used for convective terms with a grid approach for diffusive terms, This possibility will be pursued in a following section, but first it is necessary to determine more information about how the convective term is switched from forward to centered differences. This can be understood by examining the linearized Burger's equation.

gu

du

1 ~2U

0t + ~x

Regx 2"

(7)

The solution of this equation subject to the initial conditions u=un+l u = u,

forx>0 for x < 0

is u=(U"+l-u")[l+erf

x-t

When equation (7) is integrated in the same manner as (1) and the values of u <° and q
_

at

+ ux2<'> - ux~<'> Ax

1 erf Re. Ax. At

]

[u,+

~-

2un + u,-1]

(8)

8

T . D . TAYLOR

where ux2

u,+ l + u.

. . . . . . . 2-

uxl <'> - u, + u._ 1 2

-- (u, +l

- u.)F

(u. - u._ 1)F

2

}¢R,. a,)/4 erf(x) xt2 dx. F - Re ~ At -o This result provides an interesting view of the cell approach. Again, it is clear that the first approximation is unsuitable for computing the viscous terms. From the results for the convective terms, one can, however, infer how difference schemes for inviscid compressible flow can be adapted to viscous flows. In order to determine how the convective terms affect the difference equation, one must first examine the forms of ux <'>. This examination reveals that the first term is an average of propertiesbetween cells and the second term is a gradient. The gradient term is an effective artificial diffusive term with a viscosity given by (Ax. F). Note, however, that for small values of the Reynolds number, F --. 0 and the effective artificial viscosity is zero. As the Reynolds number becomes large, the artificial viscosity approaches a magnitude of + A x / 2 and consequently switches the convective term from a central difference to a forward or backward difference. With this understanding of how to "turn off" artificial viscosity, it is possible to adapt essentially any inviscid method for approximating convective terms to viscous calculations. An example of this approach for Godunov's cell method is presented in the last section of the paper. SECOND

ORDER

APPROXIMATION

RESULTS

If one employs a linear extrapolation between cell midpoints to approximate the initial distribution of u(x, 0), then a cell scheme is developed which can be compared with grid point methods. For this approximation, a solution to Burger's equation is developed subject to the conditions Ax u=u. for x < - ~ - , t = 0 u=u.+(u.+l-u.)(X+AxAX/2)

u=u.+l

for

for

-2--AX
Ax x>--. 2

The result is obtained by substituting these conditions into the general solution of Lighthill(21). Unfortunately the details become complicated and consequently they will not be reproduced here. The results are u~ = u2 - u2 t (u.+l - u.) + O[exp(_~2) ] hx UX 1 -~" -Ul - - -Ul t

(u. - A Xu._,) + O[exp(_a2) ]

An evaluation of cell type finite difference methods for solving viscous flow problems

+0[exp(_~2)]

qx<, > _ un - un-1

Ax

qx2(t)

-

Un+ 1 -- Un

Ax

9

+ O[exp(-

where U2 "~- Un+l + Un

Ul -- Un -~ Un-1 2

'

2

'

~ = A x / 2 - ~t 2 ( t / R e ) 1/2

(~ is the max of ~x or u2) subject to the conditions ~t _ Ax/2 and ~ ~ 1-0. The latter conditions one will recognize as essentially: (1) the Courant-Friedrichs-Lewy condition and (2) the viscous stability criteria for explicit viscous flow problems. The expressions for u at the cell boundaries can now be used to compute Uxx 2 and Ux22. Neglecting the error terms of order ( A t / A x ) 2 and higher, the results are

At (u: - u._l) /

UXl2(D ~- Ul 2 1 - Axx

Ux22(,) __. ~22{1 __ ~1XA(Un+ t _ u,)}.

When these results along with the expressions for qx
+

U22 -- Ul 2 2Ax

At

I[_Un+I--2Un+Un_I1

2Ax2 [U22(Un+1 -- un) -- ~12(Un -- Un-1)] = Re [

A~-~

A"

(9) This result is most interesting since one finds upon examination that the convective terms are differenced in essentially the same manner as the first order scheme. The exception being that a different artificial viscosity exists. Note also that the viscous terms are approximated in an appropriate form and the limits imposed by the first order approximation are no longer present. As a consequence, the difference equation appears to provide a practical equation for computing a solution to Burger's equation. This is further substantiated when the resulting equation is compared to other finite difference procedures. In particular, one finds that the convective terms, without neglecting terms of order (At/Ax)2, are identical to those obtained by employing the Lax-Wendroff[25] method to difference the convective terms of Burger's equation. In addition, the cell difference approximation of the viscous terms correspond to the standard second order Taylor's series approximation results. As a result of these facts, one concludes that employing a second order cell method for viscous flow problems should yield results equivalent to those obtained by employing the Lax-Wendroffapproach along with second order approximations for the viscous terms. This equivalence seems not to have been recognized in the past. S U M M A R Y OF B U R G E R ' S E Q U A T I O N R E S U L T S

The study of first and second order cell approximations for Burger's equations indicate a number of important facts. First the approximation of discontinuous initial conditions between cells appears to yield difference formula which provides a satisfactory

10

T.D. TAYLOR

approximation to the inertia terms, but fails to approximate the viscous terms in an optimum fashion. Secondly, the approximation that initial conditions are linear fits between cell centers yields a cell difference method which is equivalent to the second order Lax-Wendroff method. Another interesting result is obtained if one expands the difference form of ½0U2/OX.The result is (u.+l + 2u. + u . - 1 ) ( u . + l 4 2Ax From this one concludes that the three point average of u in the convective terms is appropriate for small cell Reynolds number calculations. One of the most significant results, however, is the general form which appears for approximating the convective (boundary flux) terms. If one denotes the center points of three adjacent cells as n, n - 1 and n + 1, then the difference approximations (either first or second order) yield the result that the flux term OF/dx can be approximated in the form

OF = t?x

Fx2(t)

Fxl

-

Ax

(lo)

where Fx
Fx <,>_F.+I +F.

fl(F.+ 1-F.)

2 Fxl<0 _ F. +2F._l

Ax 13 (F. -AxF._x)

In these expressions, F, is the initial value of F at point n and fl denotes an artificial viscosity function. This difference form seems to occur in a number of inviscid difference schemes for convective terms. For example, the schemes of Godunov[12], Rusanov[6] and Lax-Wendroff[25] all can be written in this form. The form suggests, however, that a more general approach can be adopted when applying the inviscid flow methods to viscous flows. In particular, the artificial viscosity function fl could be adjusted to become small in regions where the flow is viscous dominated and be permitted to increase in regions of inviscid flow. This would permit stability in inviscid flows and eliminate the artificial viscosity effect in the viscous regions. In the next section, the results of an attempt to pursue this approach for a 1-D viscous compressible flow is presented. APPLICATION OF THE CELL APPROACH TO COMPRESSIBLE VISCOUS FLOWS In this section, the observations made earlier regarding the nature of artificial viscosity in difference methods are employed to adopt an inviscid flow technique to viscous flows. The example problem considered is 1-D wave propagation in a viscous compressible fluid. The development of a difference method for integrating the equations of fluid flow can best be demonstrated for the equations of 1-D compressible flow which can be written in the form

1

w, + F~= ~ G~

(11)

An evaluation of cell type finite difference methods for solving flow problems

11

where w=

E=[

,

F=lp+pu 2 , u(E + p)l

and

P ]+(pu2/2)

r-l!

G= 0 + u2/2

O=

~P p. P r . (7 - 1)"

In these expressions, p denotes the pressure, p the density, y the adiabatic gas constant, Pr the Prandtl number and Re the Reynolds number per unit length. The results obtained from studying Burger's equations indicate that an appropriate difference approximation to (11) should be wnt2 - wn Fx2 Gn+ I - 2G~ + G~-I At lAx = Re. Ax 2

(12)

where quantities with no superscript are assumed to be known at an initial time tl. In this equation, the time derivatives and second order space derivatives have been approximated in a form which is consistent with the standard Taylor's series expansion formulas and the linear cell approximation results. The spatial derivatives of the fluxes, F, are not, however, expressed by direct expansions about an initial time, but are chosen to have the general form F~+ 1 + F~ F*2<'> 2 + fl(hn+l - hn)" This expression is generalized beyond the previous results for Burger's equation since (11) is a vector rather than a scalar equation. The second term acts principally as an artificial viscosity term which stabilizes the calculation for large cell Reynolds numbers. As a consequence, h need not be identical to F. The question remains of how to choose the appropriate form of fl(h~ +1 - h,,). At first glance, one may wish to select the Lax-Wendroff method for the artificial viscosity term, but a study of inviscid flow methods[17] indicates that the method of Godunov[12] appears to be an optimum approach for computing inviscid compressible forms. Based on this result, the method of Godunov was selected over the method of Lax-Wendroff as the approach to employ in this study. If one examines the Godunov approach in detail, it becomes apparent that flux quantities Fn<0 a r e calculated by an equation of the form F(OI/2 n+

_

Fn+l + F, 2

(hn+l - h,)

(13)

where h is a function of the flow variables. When this form is used to approximate the flux terms of equation (12), the result is Fn+<°l/2 - F<~/2 = F.+I - F._I Ax 2Ax

(h.+~ - 2hn + h.-l) Ax

(14)

The first term of this result gives a second order approximation for a viscous calculation while the second term acts as an artificial viscous term. Thus, it becomes apparent that for accurate flow calculations the artificial viscosity term should be damped. If we note that the nature of the flow is characterized primarily by the local cell Reynolds

12

T.D. TAYLOR

number, Re = ~(~Ax//0, since the Prandtl number has only small variations, an appropriate form for F,+ <'>1/2 would be

Fn+l/2(t)

_

Fn+12+ F.

(h.+ 1 - h,). fl(Re)

(15)

where fl is a function of the local cell Reynolds number which tends to zero as Re ~ 0, minus one as Re ~ - 0o and plus one as Re ~ ~ . A function having these characteristics is the error function, eft(x). Consequently, a form of -~- - n + 1/2 which yields second order accurate results at small cell Reynolds numbers and the optimum Godunov results at large cell Reynolds numbers is

Fn+l/2(o -- F,,+12+ F,

(h,,+x - h.)erf(Re).

(16)

Using this suggested form, the Godunov formulas for computing fluxes have been modified and the resulting expressions are (1) for supersonic flow at the cell boundary where

(P.+I-P.) > ~ -

lu*l= ux2 = u2

u,+l - u. erf(Re) 2

P~ = P2

P"+' 2- P" erf(Re)

Pxz = P2

P.+I - P. erf(-Re) 2

(17)

(2) for subsonic flow at the cell boundary where lu*l < ux2 = u2

(P. +1 - P.) erf(Re) 2x//~2 P2

Px2 = P2 - [(pu).+l - (pu).] erf(Re)

P~2 = P2

(U.+l 2- u.).

(18)

erf( e).

In these expressions the subscript, 2, denotes the average between states, n, and n + 1 defined by ~ = (f~ + f . + 1)/2 and R--e denotes the cell Reynolds number (P2 u2. Ax/~2). Note that in the Reynolds number definition it is implicit that the values of P2, Us,/~2 are scaled by characteristic values and x is scaled by the characteristic Reynolds number per unit length. If one attempts to retrieve the original Godunov method from the modified expressions it will become apparent that the equation for the density in subsonic flow is not the same. One can show that the results are consistent, however, by expanding the original Godunov density expression for the weak wave limit. As a test to make sure the

An evaluation of cell type finite difference methods for solving viscous flow problems

13

linearized density result was sufficient, the original Godunov density expressions were employed; i.e.

P~+1/2 = (~ + 1)Pn+l/2 + (T - 1)P~-1

(y _ 1 ) p n + l / 2 + (y + 1 ) p n _ l P n - 1

foru* > 0

(19)

(7 -F 1)pn+l/2 -'F (T -- 1)P.+I (~ _ 1)p,+,/2 + (T + 1)p,+, P"+'

foru* < 0

(20)

or

Pn+l/2 :

along with the linearized results to compute viscous flow through a one dimensional wave for Prandtl number of 0.75 and a cell Reynolds number with the limits 0-2 < Re < 4.0. These limits were selected since for Re < 0.2 the erf(Re) is small and the difference scheme approaches second order accuracy and for Re > 4.0, the erf(Re) is one and the flow is then basically inviscid. Two sets of initial conditions were used for the calculation. They are for and for

{i:,48° /i:2.885} /i°.8852

x < 0

x<0

0-59049~ and 0"47008)

1"0 1"0

and

1.0 1.0

for x > 0

3.259 ~ f o r x > 0 . 6.770 )

The first set of conditions represent a weak wave and served as the initial test of the method. In this test the linearized density equations were employed. Using this set of conditions the flow was computed for Re = 0-5, 1.0 and 2.0. The calculations were attempted, using the linear stability criteria At/Ax = 0"5 Re, but to insure stability the condition A t / A x = 0.125 Re had to be introduced for Re < 1. The velocity results of the calculations for the three different cell Reynolds numbers along with the values of (At/Ax) are shown in Figs. 2(a)-(c). Also shown is the exact solution for each case which was computed from the solution of Von Mises[24]. The root mean square error between the computed and exact solution is shown in Table 1. From these results it is seen that for cell Reynolds numbers less than one the error is less than 1 per cent. As the cell Reynolds number exceeds two, the error begins to increase from about 2 per cent and as the Reynolds number tends to the inviscid (Godunov) limit, the error increases to about 15 per cent for the inviscid shock as shown in ref. [17]. Table 1. Root mean square error for velocity solution of 1-D flow (Case (a)) Re - puAx #

RMS

At/Ax

0"5 1"0 2"0

0-0027 0"0075 0-0174

0"0625 0-125 0'125

14

T.D. TAYLOR 1.6

--

uI = 1'480

1.4

--

puAX =2"0

+

Computed solution Exact solution

-~

At

~x

1.2

=0425

1.0

u2= 0 - 8 7 4

0'8

I

I

-O

I

-6

-4

I -2

I

I

I

I

0

2

4

6

Number o f special increments origin a'~ inviscid shock position

Fig. 2(a).

puZ~x 1=l.O

1'6 - -

-I- C o m p u t e d s o l u t i o n E x a c t solution

U =1.480

--

I

1.4 - -

Ax-l.O

_~ 12

I0

--

u2=0.874 08

I

]

-12

-tO

I

I

-8

-6

I -4

I

I

I

I

I

I

-2

0

2

4

6

8

Number of special increments origin at inviscid shock position

Fig. 2(b).

An evaluation of cell type finite difference methods for solving viscous flow problems 1.6 -

puL~x ul = 1.4.80

~

1.4

p.

+

,

--

~+

=0.5

Computed solution Exact solution

~txt = 0.0625 ~'~,+,

~1.2

15

Ax :l.o

--

u = 0.874 0.8--

2

I

I

I

I

I

I

I

I

I

I

I

I

1

I

-16 -14 -12 - I 0 - 8 - 6

-4

-2

0

2

4

6

8

I0

Number o f spociol increments origin a t inviscid shock position

Fig.

2(c).

Fig. 2(a)-(c) Comparison of exact and numerical solution to 1-D compressible viscous flow equations (Case (a)).

The success of the method in computing weak waves suggested a test of the method on the stronger wave as given by the second set of conditions. For this test, both the linearized and strong wave density relations were employed and dual computations were conducted for cell Reynolds numbers of 0.2, 0.5, 1.0, 2.0 and 4.0. The calculations were advanced in time until the wave shape stabilized and the results showed that the nonlinear density relation did not offer any improvement in accuracy over the linearized density results. Figures 3(a)-(c) show results obtained for Re = 0"5, 1.0 and 4.0 along with the exact Van Mises solution. The results at Re = 0-2 and 2.0 did not show any significant deviation for the indicated trends. Table 2 shows the root mean square error between the computed and exact solution for the different values of the cell Reynolds number. These results show that the error from 0.3 per cent for cell Reynolds numbers of 0.2 to 12 per cent for cell Reynolds numbers of 4-0. It was found that these results could be improved slightly for cell Reynolds numbers less than 1"0 by reducing the artificial viscosity beyond that accomplished by erf(Re). A quantitative and stable procedure for accomplishing this was not developed, however. The error in the computation for the large cell Reynolds numbers could lead one to discard the method. It is important to note, however, that the large Reynolds number limit is essentially the inviscid flow limit and the test computation reduces to calculation Table 2. Root mean square error for velocity solution of 1-D flow (Case (b)) Re - puAx 0"2 0"5 1'0 2"0 4"0

RMS

At/Ax

0"003 0"024 04)44 0"099 0"12

0.01 0-025 0-05 0"10 0"10

16

T.D.

TAYLOR

3.0 I~e: 4 . 0 At

2-8 2"6--

ZX---~-o,~

2.4--

+

Computed solution

--

Exact solution

+

2"2-~ 2-0-o 1~ I - 8 - 1,6-1.4 - 1.2-1.0--

I

08

I

I

I

I

I

I~I

I

I

I

-I0 -8 -6 -4 -2 0 2 4 6 8 0 Number o f spotiol increments origin of inviscid shock position

~ig. 3(a).

3"0

--

Re

-I-O

-

Exact solution

2"8 26 2"4

+

~

2"2 'G 2" 0 O 18 I 6

1"4 I 2 I 0

08

I

I

I

I

I

I

I

I

I

I

I

I

I

-10 -8 -6 -4 -2 0 2 4 6 Number of spotiol increments origin ot inviscid shock position

rig. 3(b).

An evaluation of cell type finite difference methods for solving viscous flow problems

17

:5'0 2'8

t~e = 0 . 5

2'6

A* Ax

= 0"025

"~ 2 2

+

Colculated solution

~

--

E x o c t solution

2.4

20 1"8

+\

1.6 1.4 1.2 1.0 0.8

+

I

I

I

I

I

I

I

I

I

-18 -16 -14 -12 - I 0 - 0

-6

-4

-2

0

2

4

6

8

I

I

I

I

I

Number of spotiol increments origin of inviscid shock position

Fig. 3(c). Fig. 3(a)-(c). Comparison of exact and numerical solution to 1-D compressible viscous flow equations (Case (b)).

through an inviscid shock wave. When viewed in this context the error is less than would be obtained from other 1st, 2nd and 3rd order methods as was demonstrated in ref. [17]. It is also important to note that the error indicated when computing through a shock wave is an upper bound and should not be viewed as the expected error when gradients are not extreme. As an example, if the shock wave is replaced by a rarefaction wave in the inviscid case the error drops to about 2 per cent. Consequently, the method which has been outlined should perform reasonably for mixed viscous inviscid calculations provided the gradients do not become severe in regions of large cell Reynolds numbers. In the regions of severe gradients the method will remain stable but it will exhibit an artificial viscosity effect if the cell Reynolds number becomes large.

CONCLUSIONS

The study of Burger's equation indicates that the cell approach and the grid point method yield a similar form of difference approximation for convective terms. The difference between the methods occurs principally in the nature of the artificial viscosity behavior. For the viscous terms the grid approach and the second order (linear fit between cells) cell approach appear to yield the same results for difference representations. The first order 0ump in conditions between cells) technique appears to be unsatisfactory for differencing viscous terms. It does, however, provide useful information regarding convective term differencing. The combination of the convective cell differencing with standard grid point differencing of viscous terms appears to be a practical method since it provides a procedure for controlling artificial viscosity. A test of this approach on 1-D viscous wave propagation using Godunov's method showed encouraging results when compared with the exact solution. The true test of extending the results to multi-dimensions remains as a research task.

18

T.D. TSYLOR REFERENCES

1. I. O. Bohachevsky and R. E. Mates, A direct method for calculation of flow about an axisymmetric blunt body at angle of attack, AIAA Jl, 4, 776 (1966). 2. G. Moretti and M. Abbett, A time-dependent computation method for blunt body flows, AIAA Jl, 4, 2136 (1966). 3. S. M. Scala and P. Gordon, Solution of the time-dependent Navier-Stokes equations for the flow around a circular cylinder, AIAA Jl, 6, 815 (1968). 4. S. Z. Burstein, Numerical methods in multi-dimensional shocked flows, AIAA JI, 2, 2111 (1964). 5. F. N. Frenkiel and K. Stewartson (Eds), High-speed computing in fluid dynamics, (Proceedings of a Symposium), Phys. Fluids Suppl. (1969). 6. V. V. Rusanov, Supersonic flow around a blunt body, J. Comput. Math. math. Phys., (in Russian) 8, 616 (1968). 7. B. S. Masson, T. D. Taylor and R. M. Foster, Application of Godunov's method to blunt body calculations, AIAA Jl, 7, 694 (1969). 8. T. D. Taylor and B. S. Masson, Application of the unsteady method of Godunov to computation of supersonic flows past bell-shaped bodies, J. Comp. Phys., 5, 443 (1970). 9. R. A. Gentry, R. E. Martin and B. J. Daly, An eulerian differencing method for unsteady compressible flow problems, J. Comp. Phys., 1, 87 (1966). 10. F. H. Harlow, The particle-in-cell computing method for fluid dynamics in Methods in Computational Physics, Vol. 3 (Edited by B. Alder, S. Fernbach and M. Rotenberg). Academic Press, New York (1964). ll. W. F. Nob, CEL: a time-dependent two-space-dimensional coupled Eulerian-Lagrange code in Methods in Computational Physics, Vol. 3 (Edited by B. E. Alder, S. Fernbach and M. Rotenberg). Academic Press, New York (1964). 12. S. K. Godunov, Finite difference method for numerical computation of discontinuous solutions of the equations of fluid dynamics, Mat. Sb, 47 (89), 271 (1959); translated by I. Bohacheysky. 13. S. K. Godunov, A. V. Zabrodin and G. P. Prokopov, The difference schemes for two-dimensional unsteady problems in gas dynamics and the calculation of flows with a detached shock wave, English transl., J. Comput. Math. math. Phys., Acad. Sci., USSR, 1, No. 6 (1961). 14. J. S. Allen and S. I. Cheng, Numerical solutions of the compressible Navier-Stokes equations for the laminar near wake, Phys. Fluids, 13, 37 (1970). 15. E. L. Rubin and S. Z. Burstein, Difference methods for the inviscid and viscous equations of a compressible gas, J. Comp. Phys., 2, 178 (1967). 16. E. L. Rubin and P. K. Khosla, A time dependent method for the solution of one-dimensional radiating flow, Z. angew. Math. Phys., 21, 962 (1970). 17. T. D. Taylor, B. S. Masson and E. Ndefo, A study of difference methods for solving viscous and inviscid flow problems, d. Comp. Phys., 9, 99 (1972). 18. J. G. Trulio, Theory and structure of the AFTON codes, Technical Rep. AWW-TR-66-19 (1966). 19. K. A. Bagrinovskii and S. K. Godunov, Difference methods for multi-dimensional problems, Dokl. Akad. SSSR, 115, 431 (1951). 20. N. N. Yanenko, Method of Fractional Steps for Solution of Multi-Dimensional Problems of Mathematical Physics. Academy of Sciences SSR, Siberian Division, Novosibirsk (1967). 21. M. J. Lighthill, Viscosity effects in sound waves of finite amplitude in Surveys in Mechanics (Edited by G. K. Batchelor and R. M. Davies), p. 250. Cambridge University Press (1956). 22. J. D. Cole, On a quasilinear equation occurring in aerodynamics, Q. appl. Math., 9, 225 (1951). 23. E. Hopf, The partial differential equation u, + uux = t~uxx, Communs. pure appl. Math., 3, 201 (1950). 24. R. VonMises, Mathematical Theory of Compressible Fluid Flow, p. 140. Academic Press, New York (1958). 25. P. D. Lax and B. Wendroff, Difference schemes with high order of accuracy for solving hyperbolic equations, NYU Rep. 9759, Courant Inst. of Math. Sci., New York University (July 1962).