Sep 1, 1981 - (NISI-CB-173254) TUDEE-DIMENSIONAL. EULER EQUATION SOLUTIONS USING FL~X. SPLITTING1Mississippi State Unlv.) He A04/MF AO ...
NASA-CR-173254 19840008789
A Reproduced Copy OF
Reproduced for NASA
by the
I'IASA
Scientific and Technical Information facility 111111111111111111111111111111111111111111111
NF01452
I'
SEP 1 1981 LANGLEY 1\·~cSfAf?CH GEN I ER 1.18RArlY, NASA
FFNo 672 AU~J 65
tlAJYIf':ION,
VlltGII':II"
, ,I
.-'
~
.
'. t
THREE-DIHENSIONAL UNSTEADY EULER EQUATION SOLUTIONS USING FLUX VECTOR SPLITTING
(NISI-CB-173254)
TUDEE-DIMENSIONAL
EULER EQUATION SOLUTIONS USING
FL~X
SPLITTING1Mississippi State Unlv.) He A04/MF AO 1
Ui~STIHWY
N84- J 605"1
VECTOR
58 p cseL t2A
Unclns
00555
.. ,..
,
.
.;
.,/ uect::mber 1983
.
'
,
~:.
....
'
-. ,.,,"
.
., ' .. "
.
r,.' ,
.
, iii'"
FOREWARD Thls research was sponsored by the NASA Langley Resenn:h Center under G~ant ~A~A NAG-1-226 with Dr. Frank Thames as TecJmifal Honitor,
and the United States Air Force, Eglin Air Force Base. Florida under Grant F08635-82-K-0409 wi th Dr. Larry Lij m.rski as Tt'chnical:. Officer. The support of NASA and the Air Force, as "Jell as tedmical. interactions with Drs. Thames and Lije\"s1d. in gratefully
"
acknotvledged~
. \
".,
," .
.... '
,./
" " I
,.\
"
..
t '.: . ~' ,
,
"
, ....
ABSTRACT A method for numerically,.soJving the three-d:f.mensional unsteady Euler equations using flu;cvt!"etuj: split.ti.ng is developed.
The equations
are cast in curvilinear coordinates and a finite volume discretization is used.
ru.
explici.t upwlnd second-order predictor-corrector scheme
:L8 used to solve the discreti:r.cd equations.
GFL
num~er
The scheme is st[tble· for a
of 2 and local time stepping is used to accelerate convergence
..
,,",,;
for steady-state problems.
Characteristic variabl~'boundary conditions
lIre developed and used in the far field and at surfaces. dissipation terms are included in the scheme.
No additional
Numerical results are
compared with results from an existing three-dimensional1!:uler code
and experimental data,
.....;
.,
.\
I.
INTRODUCTION
In support of NASA's interest in the use of prop-fans as a propulsion device, acomput'll1,t'l:aliitl mert,Jlod was recently developed for numerically solving the flowfield about a swept. tapered, supercriti.cal '''ing with a propeller producing thrust and
s\~irl
(Ref. 1).
The equations solved
were essentially Euler equations with body force terms included to simulate the propeller.
The Euler solver used was an extension of the
method of Jameson. et. aL (Refs. 2 and 3) xcferred to as FL057.
Although
good agreement was obtained with exped.mental data (Ref. l). some difficulty in convergence, particularly the energy equation, was encountered using the FL057 central-difference scheme.
Recently. convergence
difficulties Here also reported by Swafford (Ref. 4) in solving a sct of hyperbolic equations \'Jith source terms (ltlh:tch can be consj.dered similar to the force terms included in the Euler equations in Ref. 1) using a similar central-difference scheme. the convergence problems in Rcf. 4.
An upwind,. scheme eU.m:l.nated
Horeover, any additional smoothing
added to the upwind scheme was found to be detrimentaJ. with regard to convergence •. (Additional smoot.hing was found to ahmrys be necessary in usi,ng the central-difference scheme in Ref. 4,)
Tlt~e):efore.
it
seemed appr9priate to consi.der an upwind scheme for Sl:!lv:!.ng the thrcedimensional Euler equations. The flux-vector-split form of the equations used is developed in the following section.
The motivation and bacl,ground' of using splitting
is given in the literature {see, for example, Steger and \varming (Ref. 5» and is not repeated here. given in the literature
Horeover. many of the ma.tr:kes needed are (see Refs. 6 and 7); hot-lever. all equations
~
.
"
and matrices needed are developed and included in Section II fo·r clarity and completeness.
Formulation of
t~e
equations for numerical solution
a.nd the algorithm used to solve the discretized equations are discussed in Section In,' 'i1!f~~;; 'Ciu.:ract"~'· ~ siic variable boundary conditions used :l.n the far field and at surfaces are developed in Section IV.
NUmerical
results, including compa:dsons with FL057 solutions and experimental data, are presented and discussed in Section V.
"r'
.,
.,
II.
SPLITTING
The conservation law vector form of the
I:U))~$.!·Jr
equations in
Cartesian coordina:t.#asx .• y. and z are
(2.1) l~here
-
q '" [p, pu, pv, P\'1, el
T
f :: Cpu, pu 2 + p, puv, puw, u(e -7i' p),l g '" [pv, puv, pv
2
h '" [pt., • puw, pvw • P'd
P
=
.
+ p,
(y - 1) [e - JliP(u
pvw, v(e +. p}}
2
T T
+ p. wee +. p.~, JT
222'
+
+ w )1
v
.
'..
Using curvilinear coordinates defined as E:; ::
..
E:;(x,y,z)
.1
n '" n(x,y.z) I; '" l;(x,y,z) T '" t
,. ""
i t is straightforward to transform Eq. (2.1) t.'o
(2.2) "
,,,here
.
. Q" J [p, pu, P\', pw,
el
T
"
,
F" J[pU, puU
+
G = J[pV. puV
+ nxp.
H" J[pW, puW
+ C;)~p, pvW + C;yP, pwW ...
l;xP , pvU+E:yP, pwU pvv
+ nl'
pwv
+
~zP. U(e
+ nzp.
+
p)]
T T'
V(e
'+
p)]
l;zP, W(e
+
p)}T
, J = xE;(ynzr; - znYI;) -YE;(~nzz; - znxr.)
+ z~(xnYr.
- Ynx~)
z y ) n r;
nx 1;
E;
x
Y
.. J .. J
'" J
ny '" J /;
Y
E;,
z
n
, Z
'" J
-1
(zf;Yl; - y :r. )
t:
-1
(z x
n /;
-1
-1
/;z .. J
(Zt;X
n
-1 (Y
-1
- x z ) n r.
- x Z ) t; n
- Yn,cr;)
t;'\ -
(xE;,Y
n
x,y 1;)
- Y x )
t;, n
U .. t; u x
+ t: v + t: zv.' Y
V'" nx u
+ nv + n
W ..
ORiGiNl~t. rt'~c!S r~
OF POOR QUALITY
(Xt;Z I;; - z t; x c; )
= J -1 (x nYr.• = J
zE;Y ) n
(Y E;zn
-1
I;;
'
~
Y
r; x u + r; v + y
.'., ~
\oJ
'(,zw
The strong conservation law of the Euler equati.ons i.n curvilinear coordinates (Eq. (2.2»
can be wr.:ltten in the quasilinear form
.. ,
.
"
j.
, -I' . " ..
(2.3)
"
. ~ ~ .' I
,
"',.
i \ .. "
.
where the matrices At B, and C are given by
'.'.3k '" k xu +
It v Y
+ k zW
'"
Ok
(2.13)
A5 ., 0 I, It. !~
-
C
I17k I
/
ORlGfN.4L P,l~Cf! ffJ
OF POOI~ QUALiTY 1
2
345
,,,here A , A , A , A , and Ak are the eigenvalues of a for k = E;, b for k k k k k = 11. and c for k ::: r;.
Correspondingly,
O]T [k;x' 0, k z • -ky ,
Xl
O)T [It • ·-k z' O. kx ' y
x2
::
x3
= [It , k ,
z
1
x .. 4 x5
the eigenvec tors are
y
[l~
e
Ii
-kx • O. ('IT
(2.14)
-
T Kx ' ky • kz • pc]
1 = -12
P pc] 'I." ['(! • -kx' -ky • -k z'
_
k
where kx =
x
TViZr =
k
x
k
ky
=~T
(2.15)
Sufficient equations have now been d,~veloped to carry out the splitt:l.ng.
The integral conservation
18\11
form of the Euler equations
will be formulated and discretized for numerical solution :i.n the following section.
This formulation requires evaluation of the
vectors F. G, and H at faces of the f:Lni.te volumes.
By splitting
the vectors F, G. and Ii into thE! sum of separate vectors. each having an eigenvalue as a coefficient. the evaluation of each separate vector by extrapolating from the .!lppropdate direction indicated by the sign of its eigenvalue can be carried out.
OH!GiN;\~. i:;:.,~l~ ~r~
OF POOl1 QUt>.U'f'l The ve..:tors F, G, and H are homogeneous func.t:h:vns
(I.e
degree one
in Q; therefore, by E41er's Theorem (Ref. 9) ;2.H)
K = KQ
where K == F and K == A for k K = Hand K '" C for k
=
.=
E;. K
=G
and K
'C
B for It == TJ. and
The matrices K can be· t;rrlt.ten
1;.
(2.17)
where Ak is the diagonal w.atrb: whose diagonal elem.€nt·"I are the eigenvalues of
K
(cr, also, K) given by Eq. (2.13)
The matri.ces K given
by Eq. (2.12) can be Hrittcn
K ==
where .t·he columns of P
k
respective eigenvalues.
P II p-1
are the
(2.18)
k k k
(;.~gel1\
eCLot!:; of
K
corresponding to the
From Eqs. (2.8) -
K '" MKM
-1
(2.19)
Using Eq. (2.18) in Eq. (2.19) (2.20)
Th,en from Eqs. (2.17) and (2.20)
(2.21a) ··1 -1
Pk M
The matrices M and M- 1 are given by Eqs. (2.9) and (2.10). matrix P
'\
k
is
(2.7.1b)
The
. ,
.;
.' , \\
.....:-
, ,
-~ ~,
\
:~ ORIGINAL PI\[~2 ;r:J POOR QUAtnv
OF
fI i \;.
k )C
1,
y
-k Z
k
k
0
-k x
k
0
Y
Ii
'-..... ~'" i
1' • ...,
\'L
P
I
,
,,,/~' /
k
..
z
k _'L
If
"I.
,.>
k
-k
Y
x
Z
If ac
0
0
0
~
k -2S.
0
! !
a
z
2
.
__k!S
Ii __It.l:'
(2.22)
.,If k
--1'2 &
ac
2
where
... ]
Rather than inVel"t P cUrcctly. the mntrix P ' can bu dt\t(~rl\linerl casily k k from the left: eigenvectors of
of P
-1 k
left (~igmlV(lctors of
K
are the
rO\,6
• Le.
k
k
o
.
-kx
0
Y
k _l£.
-..Y-_
k
k
If
fi
Ii
z
k
-
-ky
It
k
z
z
0
-k
Y
o
k
0
x
..
k
It'
'."lli;!
Ko
x
----
/2
x
(2.23)
...' "
k
k
12
1£
.. _'1-
. ~
z
z
I"~
.-:. 4
:,
"., .
.!/. \
ll~ ~
,
'.
OmG~NAl ~AGI! V~ Of' pom~ Q"AU1'Y
\-Ihere P ""
1 .-=---
Ii
~;i~
The mat.rices Tk and T;l can now be determined using Eqs. (2.21). These matrices are k
TIt ""
k It
k Y
uk x
uk
vI-k
0
0,
fJ
0
0
;.,3
0
0
0
0
0
>-,.t.
0
0
0
or
,1:.-
. 2
0
II
It
..
It
I~·
(2.26)
Q
~~•. i
. •.
;
,can be writt:en II
k
c:
,,1 I k
'.l
;It ,r ).4 I., 1,2,3 + k 1,+ + 'fu. ,)..
(2.27)
.
1iJhere 11 2 3 is a matrix \vhich has unity lAS t:he. £:h:st three diagonal
,
elCluents and Clll other elc:mcnt$ zero, 14 has !)A1ity as the fourth diagonal .~lement
and all other ele;nents zero, and IS i\13:'", unity as the fifth
diagonal element and all other elements zerO'.
lirO!t.e that only three
terms are needed on. the right hand side of Eq: •. (2 •.27) instead of five because the first three eigenvalues are I'he same '/;;' and ASr,; • , (. ':I eigenvalues are interrogated. and Similarly for Qi . ~n 1 Qi j f k except A •
, . k' -!2.J.
~n
i.
.. OR/(~/N I..\./
/.,,:.... .. " .. rr)f'.t.:d .... :.-3
. OF POOR QUALITY
Corrector:
(3.10) " "
'::n+l where the Qi+~,j,k' etc, variables are determined: in the same way ~n
as the Qi+'ii,j ,k' etc, varia.bles as described belo\Yi Eq. (3.9
n Qi.l,k' etc, are used in place of Qi,j,k' etc.
)
except
"n
The Qi+lti.j.k' etc, are
determined by
n Ql-1,j,k
,1 4 \5 if the corr:iEsponding A~ • A~, or At; eigenvalue evaluated at i+~,j.k
is :;. 0 or "n n Qi +~, I j , k" 2Qi+l
1
4
if the conr.e:sponding
A~,
eigenvalue evaluated
~t i+~,j,k
5
'.)
I.E;' or I.E;
is < 0
"n
and sim:l.1arly for Q,
I.
l.-'~.
j
,
k'
...n
Again, the same is oY:r.ne for Qi j+l_ k ,
-"'2,
n and Qi,j ,k±.lli by interrogaUng the appropriate nor. t eigenvalues. Although the algorithm given by Eqs. (3.9) and< (3.10) is an extensi.on to thr.ee dimcnsiona of that used in Ref. 13 for two dimensions, the interrogation of ei.genvalues differs significD.l'1.t1y •
Three methods
were tried in Ref. 13 to handle the computation of flow variables at cell faces Hhen e:lgenvalues were of different sign on either- side of
..
~
,/
ORIG1rJP.L FJ\G;;." ~'.; OF POOR QUALn:y the cell face.
The approach taken in the present
sdH~me
is to compute
e:Lgenvalues at cell faces rather than at cell centers as illustrated in F:lg. 2 on a z; ..
constlmt':~lan"',
changing sign is thus eliminated.
The difficulty ass.-r"""'_·""".~~""""'~""'''_''' __~-''''·''_'~'' _ _ _''''''-'''---_''--_'·-
,
.•~-----,.,---.,.... ...
G, I.
Ci~IGli~f"tl P!:.C!-~ tS OF POOR QUALITY
.
a(p-l q) h.,o
aT
a(p-l -I- It
q)
k,o
ak
k
+ s
k,m
.,
0 .
(4.4)
Defining the characteristic vector, W , as k ,I'
.
(4.5)
Eq, (4.4) becomes
(4.6)
The elements of the characteristic vecto.r, H , are -::.alled characteristi!! k variables, and are denoted hy wk,:i.'
., are determLed using Eq. (.f. 5). k ,~
The charac-teristic variables, w The vector q is
q = J [ p, -1 and the matrix PI
( ,0
U,
T v, w, p]
is given by Eqs, (2.23) ,.here the variables p and
c (speed of sound) in Eq. (2.23) are denoted rand c o
a reference
~onditil)n.
Using q and P
-1
k,o
c
to indicate
in Eq. (4.5), the elements
of the characteristic vector (4.7) "
.' are
.
~,
wk,l .,
J
IVkl
[kx(p -
-1-) c
+ kz v
k \~'l y
°
J
Il1k/
(I•. 8b)
c0
··,. ·..... ., ,,,
"
',.'
wk , 3
=_J_. [k (p
IVkl
z
- _1.'_) + k u 2 y c
°
.
"
wk 2 == - - - [ky (p - .J?i) - k z u + l~x 'v]
,
\
(4, Ra )
- kx v)
(4.8e).
'.
,,'
'.
"
,
.
1 . -J - -
(4. 3d ),
.. ---IVkl 12
(4.8e)
IVkl
Ii 1
J
The chat'llcteristic variables correspond in o'l"der to the cis111 in
3a with the flow in '.i'! dir0ction of increaslng
(~omputntl.onal
sign
l~ig.
cO(I:!.'dinHte
k, the first four eigenvalues ure posi.tivo nnd the fifth 113 ncgnth,"Iil.
,.
.,
'"
~
For the subsonic inflow case 5hO\m in Fig. 3b ,·lith the flow in
t.ht~
directioll of decreasing computational coordinate k, the first thl'.'tle nnd fifth eigenvalues are negative. and the fourth eigenvalue is positivu.
For a totally general three··dimensional code eithfell' situation in rlg. 3 Each possiblity in ! _.
1'ho subscript a refers to approaching the
the situatIon in Fig. 3b.
lind
r·· ',..
ilK€< defined by Eqs. (2.15).
Equat:ions (4. 1 2a),
(4. 12b) , (4.12
P.e
could be the ambient static pressure,
, ..
.,..''
·t· . ~ '1
~.: '.
,-' 'I.",'
\
ror a boundary across which there is no
flOlli
\
the first three
eigenvillues given by Eq,AA':'ff) are zero. the fouz:th is positive. and 'the f,ifth is negative.
One' condition must, therefore,. be specified.
'rhe condition specified is that there is no flow across the boundary. The following relations among the characteristic variables are used to determine the boundary conditions
-1)
- k w] + kv z y b
[k (p - ..Q_)
[~(p -
c0
y
[k
._----.
c
(p -
g'
2
=
[k (p •. .....E.) +k.v - k w) z x 2 y r
(4.158)
k u + kxw]b z
[k (p •. _2...) - k. u y 2 z
+ k xw] r
(4.15b)
k u - k v]
[k (p - L) +. 11 u - k v] z Y x r c2 o
(4.15c)
c0
0
.l?.) 2
+
Y
x
b
Co
/
.....
c0
'.
[kxu + ky v + kzw]b = 0
,'
kw) ) z
b
..
[.p IVk I. -+ Pc
o
0
~Jf.",~U '"
... +
k v Y
+ It w») z
r
(4.15c)
.,"f
.. .(
The subscript r refers to a reference vahlc. which is selected as the center of the first cell from the boundary.
The'minus and plus signs
-., ~.
:In Eq. (4. 15e)
corr~spond
to the locntion of the point r.
If the point
, "
lC
is in the positive k direction from the boundary
then the minus
:;1gn is used in Eq. (l1.15e). and i f it is in the In:bt.us direction then
"
.'
.1
the plus sign is used.
""
.
, \
Finite volume codes only require the pressure at an impermeable boundary and consequently Eqs. (4.15a). (4.l5b), and (4.l5c) are not needed.
However, to
fad:;i.~;',;t:;lW.,"','Lhe
handling of
point~'
near boundaries I
and aid in code vectori.zllItion, phantom points are used in the present' version of the code.
The use of phantom points r('!quires information of
variables other than pressure to ensure, for eJtal1lple, an impermeable boundary.
.
l/!Cl'O
flow across
Such information can be obtained from Eqs.
(4.15) • '~
Equations (4.lSd) and (4.1S e ) can be solved for P • b
Equations
(4.15a) - (4.1Sd) can then be solved for the remaining four variables. The solution of Eqs. (4.15), is ~
P b '" P r :j:.' p oc 0 (kx u r
P b '" P r
ub
.. ur
;.
Pb - P
c
r
(4.J.68)
(4.l6b)
2 0
- k (k u
x
+ kv + k Zwr ) y r
x r
+kv + k ".. ) Z r y r
(4 .16c)
.
4,
''..
"
v = v
,- k (k u
+ k Z wr ) y x r +kv y 1'-
(4.16d)
+ k Zt"r ) wb '" wr .- k Z (kxu r +kv y r
(4.16e)
b
r
-. .
. '
"
where the point r 1s the center of the first cell fr,?1n the boundary and the minus sign in Eq. :I '
(4~16~
is used if r is in the positive k
direction from the boundary, and the plus sign is used if r is in the negative k direction from the boundary.
,......
"
. i , I '.
/."
I
I
......
4.3
Phantom Points Phantom points are denoted by the subscript p.
obtained from the relations
/
/ ,/
The points are
",.
:,.'"
(4.17a)
p
U
p
p
:: 2p
::
2u
b
b
- p
- u
(4.17b)
in
(4.17c)
in
(4.l7d)
(4.17e)
where th(,: subscript in refers to the center of the first cell inside the computatj.onal domain and can be any of the points a, in this section.
./ .
For
.t, 'or r used
the phantom points for an impermeable
~xample.
. surface are
+ 2p oc 0 (Itxu r . + Iey v r ~
p
p
Pp
\
\
u •
.
Pr
2(P
. ,.\ .
=:
p
.. p
·r
., u
r
+
k2\~r)
- P )
b +--r2
c
(4.18a) (4.1Ub)
"
.,
"
0
.
x r + kv y r + k:e.\-lr>
(4.18c)
- 2k (k u +kv + k z wr ) y x r y r
(4.1Bd)
2k (k u
x
..
I
v
/
p
.. v
r
./
w .. w p r
/
"
:>'~-~' .I
z
x r
+kv
y r
+kw) z r
.. ~
(4.18e)
The velocity vector components are the same as those used by Jacocks and Kneile
.I
2k
(k u
/:
(R~!f.
14).
,',
..1~
' ,
v (t~
' ....... r.-:.
E t" -< '.:.~
..
ltil~~I~~ ~ I·
I
'It 0 •• :"
'~
-cp
~ ".~
. \
1~ ........
--,.... "'-.~......
,"
,".
---."~,
.....
/_
-
.....
.
'-,
..
'~
\~:
\, .\....... ,
\-
'.\," ..
\ \;
'"
.\
• F ... • •
..
'-,
\,
...~,,:-
~-=-
-.
..........
,.-
.:.!
~~'--~.; --. ...
'>-.
---.~.
.......
!
.\
•• _k
J"'-.. .', • . ,'\. '~.:'
~
".
~
....
'
~
",
-
'-.
"0
"
_
~.~~--~------
.
.....
!!. 1.4
r
ONER4 ~16 WING' M", = 0.84
.1~2 ~ A I.D
c.~
~~~~ I
!
HI
\
I
=
"'1
11
=
0.45
~
---:.
» \~-'e.-~/.
-cp
0
v...,
\
c.~ i 11 13.4
nh
~
-...:.
"\
H-Il
ti
.
:.:'
.
.
:
::::
"""" x LO
1.O
r-i
8
tJ
I='
r-I
/'
,/
;
/'
. j'
'-.J"
'-'
I'..
-' .:::r;
< 1.O
O'l
/
LO r-f
r-i
>< 1.O
;.
II
II ,~
'..0
-
.-...
II
r,
1.O
1.O
X
~
rl
s::
roof
1.O
rl
.-f X
X
,
~O
l.O
cn
cn
'-'
"-"
, I'·
'-. N
• ;.,1
-
t.:I
•
.
100
,
",
\
'"
,
),
'"
0.1
o
o
CD(TOTAL.)
I
\/ \'
I;j
\ \
•
ONEPJ\ N6 WING ~l >: 0.84 01 ~.:.: 3.06 ---- fl057 (48xS)