NASA/TM--2002-211583

Towards Arbitrary Accuracy Inviscid

Surface Boundary Conditions

Rodger W. Dyson

Glenn Research Center, Cleveland, Ohio

Ray Hixon

Institute for Computational Mechanics in Propulsion, Cleveland, Ohio

and

University of Toledo, Toledo, Ohio

ICOMP-2002-04

AIAA-2002-2438

July 2002
The NASA STI Program Office... in Profile

Since its founding, NASA has been dedicated to

the advancement of aeronautics and space

science. The NASA Scientific and Technical

Information (STI) Program Office plays a key part

in helping NASA maintain this important role.

The NASA STI Program Office is operated by

Langley Research Center, the Lead Center for

NASA's scientific and technical information. The

NASA STI Program Office provides access to the

NASA STI Database, the largest collection of

aeronautical and space science STI in the world.

The Program Office is also NASA's institutional

mechanism for disseminating the results of its

research and development activities. These results

are published by NASA in the NASA STI Report

Series, which includes the following report types:

TECHNICAL PUBLICATION. Reports of

completed research or a major significant

phase of research that present the results of

NASA programs and include extensive data

or theoretical analysis. Includes compilations

of significant scientific and technical data and

information deemed to be of continuing

reference value. NASA's counterpart of peerreviewed formal professional papers but

has less stringent limitations on manuscript

length and extent of graphic presentations.

TECHNICAL MEMORANDUM. Scientific

and technical findings that are preliminary or

of specialized interest, e.g., quick release

reports, working papers, and bibliographies

that contain minimal annotation. Does not

contain extensive analysis.

CONTRACTOR REPORT. Scientific and

technical findings by NASA-sponsored

contractors and grantees.

CONFERENCE PUBLICATION. Collected

papers from scientific and technical

conferences, symposia, seminars, or other

meetings sponsored or cosponsored by

NASA.

SPECIAL PUBLICATION. Scientific,

technical, or historical information from

NASA programs, projects, and missions,

often concerned with subjects having

substantial public interest.

TECHNICAL TRANSLATION. Englishlanguage translations of foreign scientific

and technical material pertinent to NASA's

mission.

Specialized services that complement the STI

Program Office's diverse offerings include

creating custom thesauri, building customized

data bases, organizing and publishing research

results.., even providing videos.

For more information about the NASA STI

Program Office, see the following:

• Access the NASA STI Program Home Page

at http:llwww.sti.nasa.gov

• E-mail your question via the Intemet to

help@sti.nasa.gov

• Fax your question to the NASA Access

Help Desk at 301-621-0134

• Telephone the NASA Access Help Desk at

301-621-0390

Write to:

NASA Access Help Desk

NASA Center for AeroSpace Information

7121 Standard Drive

Hanover, MD 21076
NASA/TM--2002-211583

Towards Arbitrary Accuracy Inviscid

Surface Boundary Conditions

Rodger W. Dyson

Glenn Research Center, Cleveland, Ohio

Ray Hixon

Institute for Computational Mechanics in Propulsion, Cleveland, Ohio

and

University of Toledo, Toledo, Ohio

ICOMP-2002-04

AIAA-2002-2438

Prepared for the

Eighth Aeroacoustics Conference

cosponsored by the American Institute of Aeronautics and Astronautics

and the Confederation of European Aerospace Societies

Breckenridge, Colorado, June 17-19, 2002

National Aeronautics and

Space Administration

Glenn Research Center

July 2002
This report is a formal draft or working

paper, intended to solicit comments and

ideas from a technical peer group.

This report contains preliminary

findings, subject to revision as

analysis proceeds.

This report is a preprint of a paper intended for presentation at a conference. Because

of changes that may be made before formal publication, this preprint is made

available with the understanding that it will not be cited or reproduced without the

permission of the author.

NASA Center for Aerospace Information

7121 Standard Drive

Hanover, MD 21076

Available from

National Technical Information Service

5285 Port Royal Road

Springfield, VA 22100

Available electronically at http: / / gltrs.grc._aasa.gov/GLTRS
Towards Arbitrary Accuracy Inviscid Surface

Boundary Conditions

Rodger W. Dyson*

NASA Glenn Research Center at Lewis Field, Cleveland, OH _135

Ray Hixon t

Institute for Computational Mechanics in Propulsion, University of Toledo, Toledo, OH 33606

Invlscid nonlinear surface boundary conditions are currently limited to 3 rd order accuracy in time for non-moving surfaces and actually reduce to 1st order in time when the

surfaces move. For steady-state calculations it may be possible to achieve higher accuracy in space, but high accuracy in time is required for efficient simulation of multiscale

unsteady phenomena. A surprisingly simple technique is shown here that can be used

to correct the normal pressure derivatives of the flow at a surface on a Cartesian grid

so that arbitrarily high order time accuracy is achieved in idealized cases. This work

demonstrates that nonlinear high order time accuracy at a solid surface is possible and

deMrable_ but it also shows that the current practice of only correcting the pressure is

inadequate.

I. Introduction

'IGH order accuracy is required for efficient .aeroacoustic simulations. 1-3 Recently, it has

been demonstrated 4' 5 that there is an optimal range

of design accuracies in which the best efficiency and

solution accuracy is obtained. The optimal value depends upon the length of the simulation and the level

of accuracy required. For example, nonlinear inviscid aeroacoutics simulations with 6 th order accuracy

in space and time is optimal for propagating a wave

five wavelengths with absolute error less than 10 -6 .

Unfortunately, this is not typically achieved in practice because of the solid boundary conditions 6 which

may have high order spatial accuracy, but rarely exceed 2 nd order time accuracy. This severely limits the

efficiency of time marching methods.

This paper demonstrates why most surface boundary approaches are typically only 1st order accurate in

time and it presents a surprisingly simple procedure for

correcting the normal pressure derivatives to achieve

any order of time accuracy in special cases. The basic idea follows the work of Tam 7's in which walls

are modeled by pressure corrections; and follows the

work of Hixon 9 in which the pressure corrections are

treated separately from flow variables to greatly simplify their nonlinear implementation; and finally, those

contributions are combined with Goodrich's boundary

treatment t° to produce arbitrary time accuracy at a

*Aerospace Engineer, Acoustics Branch, Structures and

Acoustics Division, Research and Technology Directorate.

?Senior Research Associate, Member AIAA

Copyright @ 2002 by the American Institute of Aeronautics and

Astronautics, Inc. No copyright is asserted in the United States

under Title 17, U.S. Code. The U.S. Government has a royaltyfree license to exercise all rights under the copyright claimed herein

for Governmental Purposes. All other rights are reserved by the

copyright owner.

surface.

The usefulness of these ideas are currently limited to

cases in which the pressure corrections are much larger

than the density and velocity corrections. Nonetheless,

the procedures shown here demonstrate the feasibility

and efficiency of increasing the time accuracy at surfaces.

II. Problem Description

The following form of the two-dimensional nonlinear

Euler equations are used:

OP + u_x + V_ + p Ou -_x + Or) = Qt (1)

apO--_+ U-_x ap + V-_v ap + Tpy (ou Ov) -_x + -_y =Q2 (2)

O(puv) + O(p 2+ p) + - Q3 (3)

Ot Ox Oy

O(pv) + O(puv) + O(pv 2 + p) _ Q4 (4)

Ot Ox Oy

At the solid boundary we also have the following

infinite set of boundary conditions: 1°

a'_(? '_) - a"Q_w : (_ _>o) (5)

Ot_ Ot_

Ordinarily, the source terms (Qt, Q2, Q3, Q4, Q5)

are zero. But, for testing purposes, these terms will

be modified so that the following analytical solution is

defined: 5

p(x,y, t) -- at cos(kx_-x) cos(ky_y) cos(kt_rt) + ct (6)

u(x,y,t) = a2 cos(k_x) cos(ky_y) cos(kraft) + c2 (7)

v(x,y, t) = a3 cos(kzrx) cos(kyTry) cos(kt_ct) + c3 (8)

p(x, y, t) = a4 cos(kxrx) cos(k_ry) cos(ktrt) + c4 (9)

NASA/T_2002-211583 1
Duetothelargenumberofderivativesusedinwhat

follows,thefollowingshorthandnotationwillbeused:

x s 0°+b+kl(x, y, t) (10)

Ca,b, k : Oxaybt k

where X, when used, will indicate: updated (U), new

(N), old (O), or correction (C).

III. Maintaining High Accuracy at

Boundary

This work uses high accuracy Hermitian Modified

Expansion Solution Approximation (MESA) methods

because high time accuracy at a surface requires correcting the spatial cross derivatives of the flow variables, many of which the MESA scheme explicitly includes on the grid. Any method could be used as long

as the effect of cross derivative corrections are consistently applied. In Fig. 1, the absolute error of the

numerical wavenumber for the derivative of a harmonic

function, f(x) = exp ikx, is plotted by points per wavelength (PPW) as done in Ref. m The extraordinary

resolution of the two-point MESA schemes allows one

to use fewer grid points and larger time steps, resulting

in significant computational efficiency for problems requiring high physical accuracy as shown in Fig. 1.5' 15

However, maintaining high accuracy at surfaces, 12 particularly when waves are being scattered, 13 is required

to maintain this high efficiency.

Unfortunately, most inviscid surface boundary conditions are at most 2 nd order accurate in time since

only the following conditions:

Ot a

- 0Vc_: (a = 0, I) (11)

are used when in fact Eq. (5) is true for any a.

As an example, the numerical time advancement of

the pressure on a surface expressed as a Taylor series

is:

p(x, y, t + At) =

C po,0,0 +C_,o,lAt-_f:P-- "_0,0,2T_t _ --I-.. •

(12)

And from the governing equations Eqs. (1, 2, 3, 4)

we can convert the series' time derivatives into space

derivatives. For example:

_C p _-

0,0,1

P v P P u v

C&o,oC_,o,o+ C;,o,oC_,_,o+ _Co,o,o(C_,o,o+ C&l,O)

(13)

NASA/T_2002-211583

10°1

._10 _ ]

¢u

_= 10-6 IE i

==_o-,=1_

_t

P"t-

]0""1-

7-Po nt DRP

_. .... Optimized Compact

t "_\- I ....... Explicit 6th Order

"_-_ [ ....................Compact 6th Order

I _ ..........-_. -- MESAlgthOrder

nll:l,fl ibl

5 10 15

PPW

a) MESA scheme compares well to popular techniques

lo°1-!

_o I-i',_\'.."--.

•'_ 4| _\\ ...... "_ _ _

_I0- Fiil\", \ \ .............. "--___ ..... 1st 3rd Order

=1o-'1..-/_,_\_,,., '.... ........................ "... Z.:.-.:.-.:Zsmo,_

I1,110 r _l\/\\ \.. x. .............. - 11thOrdet

,,,o-,oL =. v. \ \.- _. ---.. ' --.... 15th ,.,,,ord., Order

\', \ \.. "--...... ....... E | "\\,, ', _ ........................ 17thOrde,

.,_ -,_ • \ ,, .... _ ...... 19thOrder

_10 F \(-\'\"\ _ ""'-. ------ 2,st Order

-v* _ \. \ k "'_-- --

, _,_,.iT, , ,_

0 5 10 15

PPW

b) Resolution of Two-Point MESA schemes

l / /' /.."

,.,_ 10' / /./ ......."

¢_ / /" ....< /" .i-"

/_ 1st Order

3rd Order

_10 2 5_ Order

,_ _/_ --.._._-_ 7th Order

_ 9th Order

10o _ 11th Order

13th Order

15th Order

17th Order

(_, 0" _/ . _ ,1.Order 19t, Order

V,,,, i, _ _ i _ _ _ 11111 i I , , , , i

0 5 10 15 20 25

-Log,0(Error Level Required)

c) Efficiency of Nonlinear Euler Solver

Fig. 1
and

-C p =

0,0,2

o + 1,o,o C_,o,oCl,o,o + C_,o,oC_,l,o + Co,o,oC_,o,o

u P u P v P

C_o,o,o(Cl,o,oCl,o,o + C_,o,oC_,o,o+ Cl,o,oC_,1,o+

CLo,oCL,o+ _C_,o,o(c1_,o,o+ cL_,o)+

_C_o,0(c2_0,0+ c_,1,o))+

P u v v v 1/p p

C_,l,o (Co,o,oCl,o,o -+- C_,o,oC_,l,O -}- Co,o,oC_,l,o) -}-

v u P u P v P C_,o,o(cL_,oCl,o,o+ Co,o,oG,l,o+ cLl,oCLl,O+

v P P u v

CLo,oCo,2,o + _cLI,o(Ci,o,o + cL,,o)+

_CP,0,0 (ClU, l,O -t- 6_,2,o) ) +

u v u p v p

-Y(C_,o,o+ cL_,o)(CLo,oC_,o,o+ CLo,oCL_,o+

p v P

_Co,o,o(CLo,oCL_,o)) +

_CLo,o(- (Ci_o,oC1_o,o + Co_o,oC_o,o + Cf,o,oC_,o+

v u lip p lip p

Co,o,oC1,1,o -_ Cl,o,o61,o,o _ Co,o,oC2,o,o_-) ) -_

P _ V _ V _) V

_C0,o, 0 (- (C_,l,oCl,o,o -4- C_,o,oC1,1, 0 -4- C_,l,oC_,l,O-/-

..v ,-_ _ ,.1/ p ,.p ± r.1/ o r:p _

t._O,O,O_0,2, 0 "I"-_0,i,0_0,i,0 T ,JO,O,O_JO,2,0] ]

(14)

These equations clearly show that for 2 nd order time

accuracy we must have the correct first order spatial

derivatives of the pressure and velocities at the surface. In particular, Eq. (5) with a = 1 provides a

constraint on these 1st order spatial derivatives by converting the time derivatives to space derivatives again

using Eqs. (31, 31):

-_ I C_.o,oCl,o,o u u + C_,o,oC_,l,o v u + (wl/p _o,o,o-_1,0,o rlp ' +

-_y t C_,o,oCl,o,o u v + C_,o,oC_,l,o v v + g."ll/p vo,o,o-_0,1,o g"_p = 0

(15)

which if not true, implies the time accuracy is only 1st

order.

Similarly, the boundary condition Eq. (5) with a =

2 will provide a constraint on the first and second order spatial derivatives (including density), which if not

true, will reduce the numerical method to below 3 rd

order time accuracy. And in general, one must enforce

Eq. (5) for all a up to the desired time accuracy. Note,

however, in the special case of a nonmoving inviscid

wall that the even alpha conditions will automatically

be satisfied if the initial flow conditions and their spatial derivatives are correct at the surface.

A. Time Advancing Boundary

If the spatial derivatives are correctly specified to

maintain the desired time accuracy, then after time

advancing the flow variables we find that the boundary

condition, Eq. (5), with a = 0 will still be satisfied as

shown below (via its Taylor series expansion in time):

_.p(_, y,t + At)_(_, y,t + At)+

nyp(x, y, t + At)v(x, y, t + At) =

= nzp(x, y, t)u(x, y, t) + nyp(X, y, t)v(x, y, t)+

o 1 ( O_(pu) Ok(pv)_

E# \nz Ot------v-+ny_)(At)k

k=l

= p (,_(_, y, t) + _(_, y, t)) +

• -_o "

0 1 (_-_.(k ok-ip [ O'u Oiv]_(At)k

=0

=0

(16)

At the new time step, the boundary condition,

Eq. (5), for a > 0 will in general not be true. For

example, we know that after time advancing the primitive variables that for c_ = 1:

n_(p(x,y,t + At)u(x,y,t + At))t+

ny(p(x,y,t + At)v(x,y,t + At))t =

= -nz (p_u 2 + 2puuz + Pz + pyuv + pUyV + pUVy)

--ny (pyV _ + 2pvv u + Pu + pzuv + puzv + puvz)

# o

(17)

Clearly it m necessary to correct the spatial derivatives at each time step to maintain time accuracy at

the boundary. In this work however, only the spatial

derivatives of the pressure are corrected because this

approach was successful for the I st order corrections

used by Tam 7' 14 and Hixon2 This approach reduces

the seemingly intractible complexity of very high order nonlinear boundary conditions to a few lines of

code. The downside to this approach is entropy and

vorticity waves may be excited as no mechanism for

properly correcting the density and velocity derivatives

is provided here. For example, the one-dimensional

characteristic form (entropy and two acoustic waves)

of the Euler equations is:

o____. + A. o____=0 for i=1,2,3 Ot z Ox

e@_)= e(p) - @,

d(v2) = d(u) + _,

d@_)= d(u)-

Al=u, A2=u+c, )_3=u-c

vi = s = const, for dx = udt

2____¢_c = const, for dx = (u + a)dt

V2 : U -t- 7_cl

va u-_ const, for dx (u-a)dt

(18)

where c is the speed of sound, s is the entropy, d(vi)

is the i th characteristic variable differential which is

integrable for isentropic flows, and Ai is the speed of

the characteristic wave.

The time derivative of the Euler equations provides

a new set of characteristics with the same velocities

NASA/T_2002-211583 3
butincludehigherorderspacederivatives:

_o(_) + (_)_ + _£ (8) =0,i = 1,2,3

-_ [-- i Ox ) -w i-_ [-- i Ot J -- k i)t Ox

d(vl)= d (-u(p=+ 9))

 (v2): I-(u + d(v_)=d -(_ - c)(_=- _))

_l=u, A2=u+c, A3=u-c

(19)

where c is the speed of sound. In this analysis, the

lower order terms are source terms and are assumed

known (using the procedure provided later).

Differentiating once more in time provides a new set

of characteristics with the same velocity as before:

d(vl) = d (u(ux(px - 9) - u(p_ - _P-_)))

/

(20)

In Eq. (18) the entropy characteristic, d(vl) does not

depend on velocity, but the new differentiated characteristic does depend upon velocity in Eq. (19). And,

in Eq. (20) the new entropy characteristic depends on

all flow variable first derivatives. By only correcting

the pressure derivatives, the Euler system is underspecified because the velocity and density derivatives

are incorrect at the surface. Ideally, all the flow variable derivatives would be corrected so that d(vl) and

the outgoing acoustic characteristic are not changed.

This issue is beyond the scope of this paper since it is

currently standard practice to only correct the pressure terms and since physically correct solutions have

been achieved with it (probably because the required

corrections in the other terms were small compared to

the pressure). This will become more important as

higher order time accuracy is required and the procedure shown next could be extended for such cases.

IV. Surface Pressure Correction

Algorithm

The basic idea 7 is to modify the normal pressure

derivatives on a boundary to satisfy the governing

equations and boundary conditions. One-sided derivatives at a surface will generally be incorrect since they

do not account for the walls presence, thus they need

to be corrected. And, contrary to the assertion of only

one physical boundary condition being available for

high order finite-difference schemes, we make use of

the infinite additional conditions 1° in Eq. (5).

For a wall with normal vector, ff = (_x, r/u), tangential vector, F = (Tx, Ty), and a numerical method with

order O accuracy, the pressure derivatives:

oa+bp

_-=----az-_V(a, b): (1 _< a + b _< O)

a_xa_y

(21)

will be modified such that the boundary condition,

Eq. (5), is satisfied while simultaneously insuring the

tangential derivatives of the pressure:

0_P 0_P 0_P V(a): (a < O) (22)

0_-, 0_-IT0_/' • ""' 0_-0_-1_/

do not change on an inviscid wall.

We can represent the new pressure, ()N, as the sum

of the old, ()o, and correction, ()c, pressures:

Oa+bp _ (_ _ + ( °_+_p_ (23)

Op Cp

= C_,b,o + C_,b,o

where

Cp

Ca_i,i, 0

-n_Co,o,= u - nyC_,o,, c + CQo,0,_j _ (-1) iMi (24)

_j-_ Aj(-1)Y Mj

j=0

Mi = (-1)_+i b i% ,_y

b=0

Aj = (26)

-n_ (T_,l,i .... C3,0,o+T_,1o'-1,_,. c 0,o,o+ _ (__ • cxo_ _-

\ 0,0,0 ] ]

C_ C _ ( T_ o,_- _....

-ny (Ta,l,j-l,u,v o,o,o + Ta,l,j-2,u,v o,o,o + 'c 0,0,0 ° ))

Ta,fLi,Lg =

_+i _ r :._c_ )(_+_)

k 0,0,0

(-1)a-i E E IHm'i'a'z'¢ (C" "V'+"+_)

m=0 i=0 L k o,o,o/

(C<oo) k o,o,o1 j

Hm,i,,x,2,i =

(a - 1)!(m + i)!

(27)

(a - 13- 2i - _ - I)!(_- 2m)!(2(i + m) + fl)!i!m!

u/

The updated term, Co,o,_, in Eq. (24), refers to the

value of C _ after the lower order pressure deriva0,0,a

tives have been corrected. For example, we first correct

the 1_t (a = 1) order pressure derivatives, (p,, p_), using Eq. (24). Next, the higher order time derivatives,

U U

C u and C_,o,2 are calculated, but using the new val- 0,0,2

ues of p= and py. These "updated" time derivatives are

then used to calculate new 2_d (a = 2) order pressure

derivatives, pz_, p_y, pyy, with Eq. (24).

This process is repeated until the desired accuracy

is achieved. A complete derivation of Eq. (24) may be

found in the appendix.

NASA/T_2002-211583 4
V. Results

Theadequacyofcorrectingonlythepressurederivativesto achievehightime accuracywastestedby

scatteringtheinitialGaussianpulse:

p(x,y,O)=O.OOlexp( ln2[(x-1)2+Y2])_04 (28)

with four grid points per unit interval from a flat plate

at progressively higher order accuracy. The initial uniform flow conditions were p = 1.4, p = 1.0, u = 0, and

v--0

In Fig. 2, pressure contours of the reflected pulse

are shown for the 3r4 order conditions. The solution is

stable, but the reflected pressure is slightly amplified.

The pressure is over-corrected (amplified) because the

effects of the other uncorrected variables must be accounted for to maintain the no-flow boundary conditions, Eq. (5). As the accuracy increases, the number

and effect of uncorrected variables increases resulting

in an over-correction of the pressure that leads to numerical instability, unless some damping mechanism is

added.

In Fig. 3, the reflected pulse using the 3 rd order

conditions are compared with the current practice of

correcting only the first derivatives, px and py using

3_4 and 1st order design time accuracy at the surface

(though both achieve only 1st order time accuracy as

mentioned). The 3 rd order time accuracy cases in Figs.

3a and 3b, showed little difference here suggesting the

higher order corrections are small, but the 1_t order

time advance in Fig. 3c generated spurious waves and

was unstable.

The 5th order example in Fig. 4 shows more clearly

the importance of correcting the higher derivatives.

Despite only correcting the pressure, the reflected

pulse is qualitatively superior to correcting only the

first order pressure derivatives. The simulation becomes unstable however on the surface after the pulse

has reflected from it at time t = 0.8. Higher order conditions result in larger over-corrections of the pressure,

as expected. Note that the MESA scheme explicitly

uses the 2 nd order spatial derivatives and the effects of

high order boundary corrections are directly observable. Other methods may perform differently for the

cases shown in Fig. 4b and 4c., though their results

remain nonphysical.

As a test of the procedures on curved surfaces with a

Cartesian grid, the same Gaussian pulse was scattered

from a sphere of radius, r = .25, with a grid density

of 128 points per unit interval. In Fig. 5, a 3rd order solution is shown. The curvature results in greater

changes in the velocity and the pressure corrections

must compensate for that, resulting in a dissipated,

but clearly reflected pulse and was numerically stable.

The higher order conditons (> 4) were unstable without adding artificial dissipation for the circle case.

As a more realistic example, observe the entropy

production shown for the 2D stator blade passage in

Fig. 6 in which the height of the grid indicates magnitude of the entropy. Only the 1_t order pressure

derivatives are corrected and for this inviscid problem

the solution should be flat, indicating no change in

entropy. Here, the curvlinear form of the nonlinear

Euler equations are solved using a 6th order compact

scheme 1S for spatial derivatives and a 4th order 5-6

Low Dispersion and Dissipation Runge-Kutta scheme

for time marching. 19 Notice the entropy production at

the suction surface and convected downstream in the

wake.

Overall, the .advantages of correcting the high order spatial derivatives is clear, and as shown in Fig. 1,

the efficiency of high accuracy is evident for demanding long-term low error tolerance simulations. But the

current practice of correcting only the pressure is not

physically correct despite the flow variables satisfying

the governing equations and boundary conditions at

the surface. The practice is still nonphysical (e.g.,

produces entropy and vorticity) because the characteristics are under-specified and therefore the other flow

variables must be corrected as well.

VI. Conclusion

Most surface boundary methods in use today

achieve only 1st order accuracy in time since only the

no flow condition is imposed. Higher time accuracy

requires enforcing the highier order time derivatives of

the no flow condition. In the special case of stationary

surfaces, the even ordered time derivatives of the no

flow condition will be automatically satisfied if they

were correctly specified at the start of the simulation.

The primary advantage of high time accuracy is

better resolution and larger time steps at surfaces.

Achieving high time accuracy requires correcting

cross-derivatives of the flow at a surface. In principle,

any numerical method could be used with the boundary procedure presented here. But a Hermitian MESA

scheme was used because of its extraordinary resolution and potential for efficiency at a surface with high

time accuracy. The exact performance advantages of

high time accuracy will depend upon the physical accuracy required.

The design time accuracy was increased here by

correcting only the pressure as is the currently accepted 1st order practice. Despite the apparently high

complexity in deriving very high order pressure corrections, a simple process for doing this was found that

requires only a few lines of code to implement. However, a characteristic analysis and numerical testing

has revealed that it is important to correct the velocity and density spatial derivatives as well. Therefore,

an important next step is to correct all flow derivatives

to avoid underspecifying the hyperbolic Euler system

and to extend this analysis to curvilinear coordinates

NASA/T_2002-211583 5
a) t_-O b) t----.3 c) t_.6

Fig. 2 Third Order Scattered Pulse From Plate

d) t_l

a) Fully 3 rd Order in Time b) Partially 3 rd Order in Time c) Partially I st Order in Time

Fig. 3 Fully VS. Partially Corrected 3 rd Order

a) Fully 5 th Order in Time b) Partially 5 th Order in Time c) Partially 18t Order in Time

Fig. 4 Fully VS. Partially Corrected 5 th Order

NASA/T_2002-211583 6
a) n=190

Fig. 5

to efficiently resolve curvature effects.

b) n=690 c) n----l190

Third Order Pulse Scattered From Circle

find the following relationship:

VII. Appendix: Derivation of Pressure

Correction

The pressure correction procedure shown in section IV was derived by assuming only the pressure's

spatial derivatives required adjustment to account for

the presence of the wall. After the corrections are applied, the flow variables and their spatial derivatives

must satisfy the nonlinear Euler equations and the

boundary condition at the surface. In addition, the

flow's spatial derivatives must satisfy the time derivatives of these equations to achieve high accuracy in

time.

Starting with the following boundary condition,

from Eq. (5), for any a:

cOa(V'_t) cOau cOav = C Q_ (29)

Ot_ - n_ Ot---- W + ny _ 0,o,_

we want to correct the pressure derivatives so that the

following is satisfied:

(30)

C C

The terms, C u0,0,a and C 0,0,_, v are the corrections to

the time derivatives of the velocities resulting from the

corrections to the spatial derivatives of the pressure.

We need to determine how changes in the pressure

derivatives will effect the velocity time derivatives so

that Eq. (30) is satisfied up to some specified order a.

This relationship can be found using the governing

equations, Eqs. (1, 2, 3, 4) and the following alternate

form of the momentum equations: is

Ou __ (U Ou -l- V Ou -.l- 10p_

or- o / (31)

Or_ _ O_ .v Ov + l op _

From this alternate momentum equation form we

d) n----1410

C

C_,0,a = C

°

U U C v ['_U

- Co,o,oCl,o,_-, + 0,0,0_0,1,_-1 + _1,o,_-1]

C

C_,0,a = C

- \ o,o,o ,o,o-1+ +

(32)

which shows how changes in pressure and velocity

mixed space and time derivatives will affect the boundary condition equation.

The critical and most difficult next step in this

derivation is expressing the right hand side of Eq. (31)

in terms of the corrected pressure's spatial derivatives.

This was done by recursively differentiating all six governing equation forms with computer algebra until a

pattern could be discerned. The following pattern was

observed and inductively verified:

c X2-'cv T •

Cl,o,a-i = L (-/a-J,J,O a,l,2,u,v

j=0

Cu a Cp T

C0,1,a-1 = E Cct-J,J, 0 a,l,j-l,u,v

j--o

C v a C

C1 ,o ,c_-i -_- E C p -3,3, ' -oTcr .... 1 j--luv

j=o

c x--" % T

C_,l,a-1 = LUa-J,J, 0 a,l,/--2,u,v

j=0

% T •

C1,0,a_ 1 = ,2 _ _vt--j,j,O a,O,3,u,v

j=0

% T

C0,1,a-1 = _ (-_c_--j,j,O a,O,c_--j,v,u

j----0

(33)

where T_,Z,;,Lg is defined in Eq. (27).

For each value of a, we can now correct the following

pressure derivatives:

_P__j,j,oVj = O, 1,..., a (34)

NASA/T_2002-211583 7
Fig. 6 Entropy Produced on Surface of Stator Blade On Curvilinear Grid

NASA/T_2002-211583 8
by subsituting Eq. (33) in to Eq. (32) and solving

boundary condition Eq. (30):

cp - n _ (35) Z Ca_j,j,oA j =- cQo,o,a5 zCo,o,a u - nyCo,o,a U v

j=o

with Aj defined in Eq. (26).

In addition, since we do not want corrections of the

pressure derivative to effect the tangential derivatives,

we impose, for each a the following a conditions :

Oap =OVj:j=O,l,...,a-1 O'r_-5-b_ j co_ea_o_

(36)

which can be written as a condition of only the Cartesian x and y pressure derivatives as:

Or"-JOrlJ ] correction.

(37)

c_

x-'% B

/_ t_.--a,a,O a,3

a=0

= 0

b=0

(--1)(b-a) (nx)(j-2b+a) (ny)(_-J+ 2b-_)

since 7_ = yy and ry = -r/_.

For each value of a, we can write the a + 1 conditions

in Eqs. (35,37) in the following matrix form:

A0

Bo,0

Bo,1

Bo,2

Bo,a-1

A1 A,

B1,0 B_,o

B1,1 B_,I

B1,2 B_,2

:

BI,,-1 B,,a-1

Cp

• . . C.,o,o

Cp

•.. C._1,1,0

''' Cp

•. • 6a_2,2, 0

Cp

• " " v, Co,.,0

5

cQ 5 -- r_ uu uv

O,O,a xCo,o,a - TtyCo,o,a

0

0

0

0

(38)

Since the right hand side of this matrix systemis

zero, except the first term, we can efficiently solve this

system with Cramer's rule for the i th element of vector

C:

Ci- det(F) (39)

get(E)

where matrix F equals matrix E with column number

i replace by vector D.

The determinant of a matrix may be solved using

cofactor expansion along a column. 17 So we can take

advantage of all the zeros in F by expanding along

column i to efficiently evaluate the determinant of the

matrix:

D1

_. co.factor

U U "_

det(F) = (Co,o,. Q5 -n _Co,o,. u - nyCo,o,a) v (-1)iMi

(40)

where Mi is the minor of the element in the first row

and ith column of matrix F.

Similarly we can evaluate the determinant of matrix

E by cofactor expansion along the top row:

det(E) = _ Aj(-1)JMj

j=O

2_ _(_-1) ,-i i

= (-1)_(n_ +nv) 2 n_ nyA_

\i=0 /

(41)

where Aj is defined in Eq. (26). The expression for

Mj, from Eq. (25), was derived by observation with

the help of computer algebra and verified inductively•

With Mj known, we can quickly apply Cramer's rule

to solve Eq. (39) for the unknown pressure correction

vector C in Eq. (38). This simple procedure is repeated

for a = 1, 2,..., O until all higher order pressure spatial derivatives are corrected.

References

1Tam, C. K. W., "Computational Aeroacoustics: Issues and Methods", AIAA Journal, Vol. 33, No. 10, 1995, pp.

1788-1796.

2Lele, S. K., "Computational Aeroacoustics: A Review", AIAA Paper 97-0018, Jan. 1997

aMankbadi R. R. "Review of computational aeroacoustics in propulsion systems", Journal of Propulsion and

Power, Vol. 15, No. 4, 1999, pp. 504-512

4Dyson, R.W., Goodrich, J.W. "Automated Approach to Very High-Order Aeroacoustic Computations", AIAA

Journal, Vol• 39, No., 3, pp. 396-406

5Dyson, R.W. "Technique for Very High Order Nonlinear Simulation and Validation", Journal of Computational

Acoustics, Voh 10, No. 2, June 2002.

6Dadone, A., Grossman, B. "Surface Boundary Conditions for the Numerical Solution of the Euler Equations",

AIAA J., Vol. 32, No. 2, February, 1994.

7Tam, C.K.W., and Dong, Z. "Wall Boundary Conditions for High-Order Finite Difference Schemes in Computational Aeroazoustics", AIAA Paper 94-0487, Reno, NV, Jan.

1994.

STare, C.K.W., "Advances in Numerical Boundary

Conditions for Computational Aeroacoustics", J. of Computational Acoustics, Vol. 6, No. 4, 1998, pp. 377-402

9Hixon, R. "Curvilinear Wall Boundary Conditions

for Computational Aeroacoustics", AIAA Paper 99-2395, 35th

AIAA/ASME/SAE/ASEE Joint Propulsion Conference and Exhibit, Los Angeles, CA, June 1999.

1°Goodrich, J.W. "High Order Implementations of Accurate Boundary Conditions", AIAA 99-1942, 1999.

NASA/TM 2002-211583 9
11Hixon, R. "Nonlinear Comparison of High-order and

Optimized Finite-difference Schemes", Int. J. Comp. Fluid Dynamics, Vol. 13, 2000, pp. 259-277

12Gustaffson, B. "The Convergence Rate for Difference

Approximations to Mixed Initial Boundary Value Problems",

Mathematics of Computations, Vol. 29, 1975, pp. 396-406

13Casper, J., and Carpenter, M.H. "Computational

Considerations for the Simulation of Shock-Induced Sound",

Siam Journal on Scientific Computing, Vol. 19 No. 3, May 1998,

pp. 813-828

14Kurbatskii, K.A., Tam, C.K.W. "Cartesian Boundary Treatment of Curved Walls for High-Order Computational

Aeroacoustics Schemes", AIAA J., Vol. 35, No. 1, January 1997,

pp. 133-140

15 Goodrich, J.W. "A Comparison of Numerical Methods for Computational Aeroacoustics', 5th AIAA/CEAS Aeroacoustics Conference, AIAA 99-1943, May 1999.

16Anderson, J.D., Fundamentals of Aerodynamics,

McGraw-Hill Book Company, 1984, pp.472.

17Anton, H. Elementary Linear Algebra, 5th ed., John

Wiley & Sons, 1987, p. 85

lSHixon, R._ "A New Class of Compact Schemes",

AIAA Paper 98-0367, Reno, NV, January 1998

19Stanescu, D., and Habashi, W.G., "2N-Storage Low

Dissipation and Dispersion Runge-Kutta Schemes for Computational Aeroacoustics", Journal of Computational Physics, Vol.

143, pp. 674-681, 1998

NASA/T_2002-211583 10
Form Approved REPORT DOCUMENTATION PAGE

OMB No. 0704-0188

Public reporting burden for this collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources,

gathering and maintainingthe data needed, and completing and reviewingthe collection of information. Send comments regarding this burden estimate or any other aspect of this

collection of information, including suggestions for reducing this burden,to Washington Headquarters Services, Directorate for Information Operations and Reports, 1215 Jefferson

Davis Highway, Suite 1204,Arlington, VA 22202-4302, and to the Office of Management and Budget, Paperwork Reduction Project (0704-0188),Washington, DC 20503.

1. AGENCY USE ONLY (Leave blank) 2. REPORT DATE 3. REPORT TYPE AND DATES COVERED

July 2002 Technical Memorandum

4. TITLE AND SUBTITLE 5. FUNDING NUMBERS

Towards Arbitrary Accuracy Inviscid Surface Boundary Conditions

6. AUTHOR(S)

Rodger W. Dyson and Ray Hixon

7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES)

National Aeronautics and Space Administration

John H. Glenn Research Center at Lewis Field

Cleveland, Ohio 44135-3191

9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES)

National Aeronautics and Space Administration

Washington, DC 20546-0001

WU-780-30-11-00

8. PERFORMING ORGANIZATION

REPORT NUMBER

E-13368

10. SPONSORIN G/MONITORING

AGENCY REPORT NUMBER

NASA TM--2002-211583

ICOMP-2002-04

AIAA-2002-2438

11. SUPPLEMENTARY NOTES

Prepared for the Eighth Aeroacoustics Conference cosponsored by the American Institute of Aeronautics and Astronautics

and the Confederation of European Aerospace Societies, Breckenridge, Colorado, June 17-19, 2002. Rodger W. Dyson,

NASA Glenn Research Center; and Ray Hixon, Senior Research Associate, Institute for Computational Mechanics in

Propulsion, Cleveland, Ohio and University of Toledo, Toledo, Ohio 43606. Responsible person, Rodger W. Dyson,

organization code 5940, 216-433-9083.

12a. DISTRIBUTION/AVAILABILITY STATEMENT

Unclassified - Unlimited

Subject Category: 71 Distribution: Nonstandard

Available electronically at t_ltp:ligltrs.grc.n;,sa.gov/GKIRS

This publication is available from the NASA Center for AeroSpace Information, 301-621-0390.

12b. DISTRIBUTION CODE

13. ABSTRACT (Maximum 200 words)

Inviscid nonlinear surface boundary conditions are currently limited to third order accuracy in time for non-moving

surfaces and actually reduce to first order in time when the surfaces move. For steady-state calculations it may be possible

to achieve higher accuracy in space, but high accuracy in time is required for efficient simulation of multiscale unsteady

phenomena. A surprisingly simple technique is shown here that can be used to correct the normal pressure derivatives of

the flow at a surface on a Cartesian grid so that arbitrarily high order time accuracy is achieved in idealized cases. This

work demonstrates that nonlinear high order time accuracy at a solid surface is possible and desirable, but it also shows

that the current practice of only correcting the pressure is inadequate.

14. SUBJECT TERMS

Nonlinear Euler; MESA; High resolution; Efficiency; Computational aeroacoustics;

Characteristics

17. SECURITY CLASSIFICATION

OF REPORT

Unclassified

NSN 7540-01-280-5500

15. NUMBER OF PAGES

16

16. PRICE CODE

18. SECURITY CLASSIFICATION 19. SECURITY CLASSIFICATION 20. LIMITATION OF ABSTRACT

OF THIS PAGE OF ABSTRACT

Unclassified Unclassified

Standard Form 298 (Rev. 2-89)

Prescribed by ANSI Std. Z39-18

298-102