THE NUMERICAL CALCULATION OF THE SUPERSONIC CIRCULATIONOF A CURRENT OF VISCOUS PERFECT GAS ROUND A SPHERE* V.K.MOLODTSOV Moscow (Received
19 December
1968)
1. IN view of the wide use of blunt aerodynamicshapes the construction of a picture of the flow round such bodies taking viscosity into account is of great practical value. In flights at great heights, where the Reynolds numberis small, the influence of viscosity becomes considerable over the whole range of flow and cannot be described by the equations of boundarylayer theory. It therefore becomes necessary to investigate the flow by means of the complete NavierStokes equations. In [ll a statementof the problem is given and a method of calculating the flow in the region of blunting based on simplified Navier-Stokes equations by scheme I of the method of integral relations is described for solving the problem using the complete non-stationary Navier-Stokes equations. To approximatethe Navier-Stokes equations an explicit finite-difference scheme is used. A solution is found by the methodof establishment. A numerical method is given below for calculating the supersonic circulation of a viscous perfect gas round a sphere, based on the complete Navier-Stokes equations by scheme II of the method of integral relations. Witha slight change of functions along the surface of the body this scheme enables us to choose a comparatively small numberof approximationnodes, which is importantfrom the point of view of saving computingtime.
*Zh.
vTchis1. Mat. mat. Fiz.,
9, 5, 1211-1217,
320
1969.
The numerical
calculation
of the supersonic
321
circulation
Notation r is the radius of the sphere; s is the length of arc measured along the contour of the body; n is the normal to the body; u is the component of the velocity along s relative to Urn; u is the component of the velocity along n relative to U,; p is the pressure relative to p,Um2; h is the enthalpy relative to UWz; p is the density relative to p,; p is the coefficient of viscosity relative to p* for h = 1; Re is Reynolds number @wU=J~!-G Pr is Prandtl’s number cpp+ f k; k is the thermal conductivity; cp is the specific enthalpy at constant pressure; y is the adiabatic index; M, is the Mach number of the undist~~d flow; Urn is the velocity of the undisturbed flow. All the measurements of length are relative to r, the radius of the sphere and parameters of the flow at infinity are denoted by the subscript 00and those on the body by the subscript CO. We shall consider the supersonic viscous flow of a perfect gas round a sphere, on the assumption that the Navier-Stokes equations are valid. The complete set of equations describing the motion of the gas are written in the followirjg divergent form: a
(
puz
p -
TJ6 -
3(GlL -
)
___ ._.. _. ..___ +;(~aii-PLIV)+--0s . i@ r a / Z,n-puv
c.-
____
iis (
a
-_
as
.--
2@ pIi
(
-
---
>
UTSS -
-___-
-
-
utsn
- 2pl9 + pz f(z*n
qs
-
1+&H VTnn
-
-
UTsn
-
-
puu)ctg s = 0, (‘4
nnn - qn)+
2%
asn
= OT(l)
qR)+(pUH - mss - man - q*)ctgs
+.
a
s
&%
8 B(pvH
pqctg
2%
- Jj - pvZ)+ 3Tnn_____..
+;@nn V-ban
puv) + (z, R + 27ra -
( > p”
+;MJ,+
2pv + p” ctg s =o, ~
(3)
= 0, (4)
V. K. Moiodtsov
322
T,* = -
2
12 3 Re
2 au --G% ds
=
-&-Fv
(
,
(7)
)
au
1
li Tan
au
au
u
+
-
Re
(8)
g%-
0%
I.8
4n=--,
where 1=
i+ II,R = h + (1~+
PrRe
~2)
/2
ah 6%
is the total enthalpy.
Here (1) and (2) are the equations of the projections of the momentumon the s and n-axes, (3) is the energy equation, (4) is the equation of continuity and (5) is the equation of state. The dependence of the viscosity on the temperatureis given by (6). The equations for the components of the stress tensor and the heat-flow vector are respectively (7)-U). We now write the boundary conditions for the problem. On the body n = 0 8h a=u=O,
h=h,(s)
or
an
at infinity n = M
and on the axes of symmetrys = 0, s = R
au
Ice-=as
ah as
= 0.
I ,=
0,
(12)
The numerical
calculation
of the supersonic
323
circulation
Henceforth we shall assume that the boundary condition (13) may be given at a finite distance R from the body. Thus the problem, reduces to finding a solution of the set of equations (l)-(6)
with given boundary conditions (12)-(14)
a&o <
n <
in the bounded region 0 < s <
R.
3. The construction of the approximate system by scheme II of the method of integral relations is well known [31. In the given problem the region of integration is divided by lines s = const, which are approximate directions, into R + 1 strips and the functions are represented as follows. For even functions (and complexes) an even trigonometric interpolation polynomial is constructed from the values of the function at !Q+ 1 points SirSiE [O,n), i =
0, 1,2,...,
93,
f (s, n) =
5
ak (n) cos ks;
(15)
k=o
and for odd functions (and complexes) an odd trigonometric interpolation polynomial is constructed from the values of the functions at !R points i =
1,2 ,...,
si,
sf E (0, n),
93, 92
cp(s, n)t
2
(16)
pk(n)sinks.
k=l
There is a linear relation between the coefficients the approximations and the values of the functions
ak (n)
and l3k(n)
in
f(si, n) and q(sr, n) at the
interpolation nodes. Equations (l), (8) and (10) together with the approximations (16) are integrated with respect to s from s = 0 to each si (i = 1, 2,. . . , g), and equations (2)-(4), (7), (9) and (11) with approximations (15) are integrated with respecttosfroms=Otoeach
sj (j=1,
2,...,
93+
I),s~+~=JL
As a result we obtain a closed set of ordinary differential equations (in the required functions ui(n), VI(~), pi(n), hi(n), i = 0, 1, 2,..., %),which is approximated by a set of finitedifference equations of the second order of accuracy in n, except for the momentum equation on the normal; the latter is approximated with the first order of accuracy:
324
(h
- pUv)t,j+$
‘(Fi,j+s +Fi,f-l/,)m Ov
2 - ‘Fnnt,j-ll,- (p + pV’)ij + (P + pVz)&j-i + :(O:r:cv, T*nO,i+$$
I + i! 3Trn +(z,,+2zs,)ctg
3Znn+Zln
ctg s
(17)
-((Zsn - P”V)i,L’Iz + 2
8 - (3PUV+ P”2 ctg 4+
&a
+ 0?w/J = 0, (13)
-P
- PfJ2)
I
(pz9 - 2PV”- Pup ctg s1-t ;h
a
1ikt%
- PUV) i fjk’lr,
The numerical
Similar expressions
calculation
of the supersonic
325
circulation
are valid at the nodes (i, j - M).
The values of the required functions at the semi-integral nodes are calculated by the formula fr, j*l/, = 'h(ft3 + fi, jai). Here the accepted notation is
!fj
=
i = 0, 1, 2, . .., 33,
ii=R/W,
f(sf9 nj)9
For the odd functions (and complexes)
zi -
= AB-‘G,
i = 0, i, 2, . . . , St.
C#J we have
B =IlLll,
A = Ila~mll,
3s
where k=0,1,2 akm =
1 .
bl, 5 -sin m
cos msk,
and for even functions (and complexes)
,...,
93,
v=I,Z
1...,
99,
1, m = 1,2, . . . , 9;
ms,,
f
ar -
= CD-fAT
c
6%
=
Ilwnll,
D = Il&nll,
where
xf cim =
sin mq,
k = $2, . . . , 9.
f(o, n)lI,
=Ilf(SAn)-
drm= L(1 - cos msl),
1, m = 1,2, . . . , 33.
m
(ah /an) Iu,= o
The boundary condition the expression h.ii -
hi0 = 0,
is approximated, for instance, by
i = 0, 1, 2, . . . , 9.
w
5. The solution of the set of non-linear algebraic equations (l2), (l3), (1% (21) is constructed as follows. We consider some non-negative form 0 (for example
@ =
x la/ v-l
or
326
V. K. Molodtsov E
a = &‘, where Av is the discrepancy of the z+th equation of the system and Y-i K the total number of equations), which vanishes for the solution of the required algebraic system. The zero minimum of the function Cpis determined by the method of local variations [41. For this some initial approximation is given which, as a result of a variation of each variable in turn, is refined until the quantity Q reaches a given value t. By the scheme given in this paper we carried out calculations of the super sonic flow of a perfect viscous gas round a thermally isolated and cooled sphere for various parameters of the incident flow. We give some results of the calculations. Figures 1 and 2 show the profiles of the components of the velocity, pressure and enthalpy on the rays s = 0 and s = 0.8 for the case of the flow round a cooled sphere with parameters of the incident flow M, = lo, Re = 70, pr = 0.7, v = ii/g; the value of the enthalpy on the surface of the body was taken as 0.029 and the distance R as 0.3. To throw more light on the evaluation of the accuracy and convergence of the scheme, calculations were carried out on various computational 9RXW networks (9~ is the number of nodes along the n-axis and H the number of rays s = const in the region o c s .c n). In Fig. 1 the continuous line gives the results of calculations obtained by the method of integral relations and the broken lines the results obtained in bl. It is obvious from the graph that the curves almost coincide in the neighbourhood of the body and differ somewhat in the region of the shock wave. Figures 2-4 illustrate the convergence of the method as the number of the approximation ‘8 increases, and also the influence of the position of the extreme ray on the parameter of flow in the field and along the surface of the sphere. In Fig. 2 the continuous line gives the results of calculations for 93= 4, 4~ = 60, 8%~ 1.6, and the broken line for R = 3, 9R= 30, w = 1.2. The maximum absolute difference in the calculations of the components of the velocity here consists of a quantity of the order lO_‘, and the relative error in the calculation of the pressure and enthalpy is about 7%. The results are practically
the same on the axis of symmetry;
thus the
The numerical
calculation
of the supersonic
circulation
32’7
s-o.6
FIG.
FIG. 2.
1.
h
1
n
0.8 FIG. 3.
!2
I.6 FIG. 4.
position of the extreme ray and the number of interpolation nodes in s slightly influence the solution in the neighbourhood of the axis of symmetry. The distribution of the pressure ratio P = P(S) / P (01, the coefficient of ct = (k / ~e)au / an and the ratio of the heat flows 4 = q(s) / q(0) resistance along the surface of the sphere is shown in Fig. 3. In the region where s varies from 0 to 0.8 the greatest difference in the results for 93= ;j and 91=.4 is about 5%.
328
V. K. Molodtsov
as
0 FIG. 5.
I
0.5
FIG. 6.
I
I
I
I
44
43
Q2
a1
I
0
FIG. 7.
The calculation of the flow round a thermally isolated sphere with flow parameters M, = 2.6, Re = -XI, or = 0.7, y = 7/S was carried out on the networks z x w: XI x 3 II 30 x 4, K = 1.2, and the coordinates of the extreme rays were chosen as s3 = 1.2, sp zsI 1.6 respectively. Figures 4 and 5 show the curves for the pressure, components on the rays s = 0 and s = 0.8.
enthalpy and velocity
The distributions of pressure p = P(S) / P (O), enthalpy E = h(s) / h (0) and the coefficient of resistance cf along the surface of the sphere are shown in Fig. 6.
16
The numerical
calculation
of the supersonic
329
circulation
In Fig. 7 a comparison is made with the experimental distribution of the density along the axis of symmetry. The continuous line shows the results obtained by the method of integral relations, and the points represent experimental results from r61. From this comparison it is obvious that the calculated data agree well with the experimental results. An analysis of the calculations shows that in the neighbourhood of the axis of symmetry satisfactory results are obtained even with a small number of strips. The distributions of pressure, enthalpy and the coefficient of resistance along the surface of the sphere are quite well defined even for 9p= 3. To construct the solution in the whole region 0 < s < ;n and obtain the best possible evaluations of the accuracy and convergence of the scheme additional investigations are now being carried out. In conclusion the author wishes to express his thanks to 0. M. Belotserkovskii and A. I. Tolstykh for their useful advice and interest, and also to V. H. Nadorov for his help in carrying out the computer calculations. TrunsZated
by H. F. Cleaves
REFERENCES
1.
TOLSTYKH, A. I. The calculation of the supersonic flow of a flow of viscous round blunt bodies. Zh. uychisl. Mat. mat. Fiz., 6, 1, 113-120, 1966,
2.
PAVLOV, B. M. The calculation of the supersonic flow round blunt bodies using complete Navier-Stokes equations. IZV. Akad. Nauk SSSR, Mekh. zhidkosti i gaza, 3, 128-133, 1968.
3.
The flow of a supersonic tel sverkhzvukovym Nauk SSSR, 1967.
gaza).
(Obtekanie zatuplennykh Ed. 0. M. Belotserkovsii, M., VTs Akad.
CHERNOUS’KO, F. L. The method of local variations for the numerical solution variational problems. Zh. uy’chisl. Mat. mat. Fiz., 5, 4, 749-754, 1965.
5.
KAO, H. C.
Hypersonic
6.
the
current of gas round blunt bodies. potokom
4.
Amer. Inst. Aeronaut
gas
viscous
flow near the stagnation J., 2, 11, 1892-1906,
and Astron.
streamline 1964.
of
of a blunt body.
IVANOV, A. V. The structure of a shock wave in air with Mach numbers from 2.6 to 6. Itu. Akad. Nauk SSSR, Mekh. zhidkosti i gaza, 2, 152-155, 1967.