The fluid flow inside the space shuttle main engine (SSME) traverses through a complex geometrical configuration. The flow is compressible, viscous, and turbulent with pockets of separated regions. Several computer codes are being developed to solve three dimensional Navier-Stokes equations with different turbulence models for analyzing the SSME internal flow. The locally implicit scheme is a computationally efficient scheme which converges rapidly in multi-grid modes for elliptic problems. It has the promise of providing a rapidly converging algorithm for steady-state viscous flow problems

Reddy, K. C.

English

NASA Technical Reports Server

I 
- !''I r , 
A Locally Implicit Method 
for Fluid Flow Problems 
K . C. Reddy 
Reddy and hssociates 
Tullahoma, TN 37388 
R t . 2 ,  B O X  253 
Final R e p o r t  o f  Contract NO. N A S 8 - 3 6 6 4 8  
f o r  Period: July 1985-July 1986 
Submitted to 
George C Marshall Space Flight Center 
Marshall Space Flight Center, A L  35812 
BY 
Reddy and Associates 
Tu!?ahoma, TN 37388 
R t . 2 ,  B O X  253 
November 1986 
(hASA-CR- 178969)  
1585 - Jul. ISEf (Reddy and Asscciates) 
16  P CSCL 20D Dnclas  
B LGCALLY l E E L I C I l  METRCD N87- 1 f 4 4  3 E C R  E U J 1 : D  FLOW & & C B L E E S  F i n a l  h f p o r t .  Jul. 1 
a G3/34 40295 
https://ntrs.nasa.gov/search.jsp?R=19870006010 2020-03-20T12:08:47+00:00Z
I 
A Locally Implicit Method 
For Fluid Flow Problems 
K. C. Reddy 
Yovember 1986 
Table of Contents 
1 . Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  1 
2 . Locally Implicit Scheme . . . . . . . . . . . . . . . . . . . . . . . . .  2 
3 . Multi-Grid Technique . . . . . . . . . . . . . . . . . . . . . . . . .  4 
4 . Three-Dimensional Diffusion Equation . . . . . . . . . . . . . . . . . .  5 
5 . Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  7 
6 . References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  8 
I 
A Locally Implicit M e t h o d  
For Fluid Flow Problems 
1. Introduct ion 
The fluid flow inside the space shuttle main engine (SSME) traverses through a com- 
plex geometrical configuration. The flow is compressible, viscous and turbulent with pock- 
ets of separated regions. Computation of such a flow field requires the application of the 
state-of-the-art CFD technology. Several computer are being developed to solve 
three dimensional Navier-Stokes equations with different turbulence models for analyzing 
the SSME internal flow. Most of these computer codes require extensive computational 
effort and are rather slow to converge if accurate solutions are to be obtained in three 
dimensional domains. This report outlines an algorithm which has the potential of solving 
fluid flow problems with rapid convergence to steady-state solutions. 
The computational methods used in the Yavier-Stokes codes fall into two major cat- 
egories: finite difference and finite element methods. Each have their own strengths and 
weaknesses. Some of the algorithms are designed to solve the unsteady Compressible 
Navier-Stokes equations, either by explicit or by implicit factorization methods, for sev- 
eral hundred or thousands of time steps to reach a steady-state solution asymptotically. 
Other algorithms attempt to  reach the steady-state solution by relaxation methods. A11 of 
them require body-fitting curvilinear grids with sufficient resolution. Grid requirements, 
however, differ greatly with the algorithm used. Implicit factorization methods, typically 
use numerical transformations whereby the transformed grid in the computational space is 
uniform and rectilinear. This requires the grid to have indices which are separable in three 
directions for three dimensional problems, and also be reasonably smooth. However, such 
requirements may introduce grid singularities when complicated domains are discretized. 
Flow solver algorithm will have to deal with such grid singularities. Explicit schemes and 
b i t e  element algorithms have less stringent requirements on the grid structure. However, 
their efficiency and accuracy also depend heavily on the grid structure. 
@ 
Explicit algorithms for solving time-dependent Yavier-S tokes equations are usually 
bound by stability restrictions and take thousands of time steps to converge to a steady- 
state for three-dimensional flow problems. Most of the current implicit algorithms use 
approximate factorization techniques and require index separable grid structure. The 
approximate factorization techniques are reasonably efficient for two-dimensional problems, 
but can be quite slow in convergence for three dimensional problems. We propose a locally 
*implicit dgoiithm, which is based on a relaxation procedure, for solviiig time dependent 
partial differential equations to obtain steady-state solutions. The basic algorithm which 
is locally implicit but globally explicit is computationally efficient and is outlined in the 
next section. It is suitable for the application of convergence acceleration procedures such 
1 
0 
as the multi-grid method which is also outlined. 
0 Test calculations with the current method for solving an  elliptic partial differential 
equation, the diffusion equation, which simulates many features of the subsonic flow prob- 
lems are presented in the third section. Both single grid and multi-grid results are pre- 
sented. Solutions for inviscid and viscous flow problems from subsonic to supersonic flow 
problems are being computed successfully by this method and those results will be pub- 
lished in due course. 
2. Locally Imdicit Scheme 
A relaxation scheme known as the locally implicit scheme is described below. The 
compressible Navier-Stokes equations can be written in the conservation form. 
aG * 
- + V . f = O  
at 
where G is the conservation variable vector and F' = (Fz, Fy, F,) represents the flux vectors, 
including inviscid and viscous terms. Addition of the equation of state makes this system 
complete. .Integrating over space and using the Gauss Theorem, 
V d 
where v is any closed volume element and s is its surface. In the finite volume discretization, 
the computational domain of the problem is covered by a body conforming grid, typically 
of quadrilateral bricks. For the current scheme the grid need not be index separable. We 
only need a grid which provides reasonable resolution of the flow field. Finite element type 
of grids are acceptable as are also the well structured curvilinear grids used by many of 
the finite difference codes. Euler implicit discretization of the equation (2.2) is written as 
- where 
~ j , k , t ( ~ )  = (F.9 j+ 1/2 ,k , i  - ( f . ~  )I - 1/2,k.t  .- ( ~ - 9  )I,k+ 112,e 
- (J ' .g ) j ,k- i /2 , t  + ( F . g ) I , k , e + 1 / 2  - (F .g  ) l , k , e - i j 2  t 
4 - 
z j + l / 2 , k , &  = (6ds)1,1/2,k,L t etc- 
and n' is the unit normal vector pointing in the increasing j direction. In the above, 
Dj,k,t(G) are the artificial dissipation terms similar to those used by Jameson, et al(5). 
They are smaller than the trucation error terms of the scheme and are designed to suppress 
@ 
2 
the nonlinear instabilities which arise due to the central difference approximation of the 
convective flux terms, without affecting the accuracy of the solution. The implicit difference 
equations (2.3) for all the cells of the domain represent a set of coupled non-linear algebraic 
equations. The coupling of the equations for each cell is primarily with the equations of the 
neighboring nodes. Thus one can devise Gauss-Seidel type of iterations where equations 
at a group of neighboring nodes can be solved implicitly with the latest available iterates 
from the neighboring zones. The simplest scheme would be to solve the linearized form of 
the equations (2.3) one at  a time with an inner iteration loop. 
0 
Locallv implicit approximation: The equation (2.3) is linearized and written in the form 
where 
Resy,k,( = -&,k,t(Gn) 7 P ] , k , t ( G n )  (2.i) 
The coefficients C J M .  CJP, etc. are matrices which include the Jacobian matrices on the 
cell faces of the ( J ,  k, e) cell and the dissipation coefficients. DJ.Ci. D J P  etc., are only 
diagonal matrices with artificial dissipation coefficients at  the cell faces of (J ,  I C ,  e ) ,  and CC 
is a matrix which can be approximated by a diagonal matrix which includes the time step 
term and the dissipation coefficients, both real and artificial. The linearized equations (2.4) 
are solved by an inner iteration for each time step by a modified symmetric Gauss-Seidel 
iteration. 
c - ~G:.:I, = Resj”lk,e - L(AG( ) ) j , k , t  (2.8) 
(2.9) A d r n )  = AG(m-*) + dG(m), m = 1,. . . , 8  
The superscript ( ) for A C  inside the terms of the L operator of (2.8) is either (m) or 
(rn - 1) whichever is the latest available value depending on the direction of the iteration 
sweep. The coefficient C is a modified form of the diagonal matrix CC a t  j ,  k, t designed 
3 
to  enhance the stability and convergence characteristics of the scheme for large Courant 
numbers. Local time stepping enhances the convergence of the scheme for obtaining the 
steady state solutions. The eight inner iterations start from each of the eight corners of 
the computational domain to make the scheme symmetric. The initial AG(O) can either be 
zero or can be the AG computed in the previous time step. The symmetric inner iteration 
is stable and removes any possible instabilities which may arise due to sweep direction 
being opposite to the flow direction. 
The equation (2.8) is now a set of scalar equations and we have explicit expressions 
for computing the corrections to the solution a t  each point. Thus the scheme is globally 
explicit while the inner iteration gives it a locally implicit character. Stability analysis of 
the model equations indicates that the scheme is unconditionally stable. In practice inner 
and outer relaxation parameters are introduced for equations (2.5) and (2.9) as follows: 
(2.10) 
(2.11) 
Under relaxation for w,,t in the range of 0.67 to 0.95 and over relaxation for win from 1 
to 1.2 are indicated for good convergence of the scheme at  Courant numbers of order 10 
and local time stepping. e 3. Multi-Grid Techniaue 
One of the most effective methods to speed up the convergence of iterative procedures 
for the solution of partial differential equations is the multi-grid technique. From a chosen 
fine grid, a sequence of coarser grids are formed by combining groups of adjacent cells 
in a systematic manner. For example, by omitting every other mesh line in an even cell 
numbered grid one can produce a coarse grid and repeating the process the next coarser 
level grid can be obtained. A set of iteration equations are solved on each of the grids 
where information from the fine grid is passed to the coarser grids and vice versa in each 
cycle. The highest frequency errors corresponding to each grid are dissipated rapidly in the 
iteration process at  that grid level and thus errors of all frequencies are damped much more 
rapidly in the multi-grid mode than in a single grid iteration. For fluid flow problems, Xi(6) 
and Jameson(') have applied the multi-grid technique for inviscid transonic flow problems. 
A summary of the multi-grid technique is outlined below. 
At a particular grid level denoted by a mesh spacing parameter h, one or more time 
step calciilations are carried out and the solution Is updated. The solution Gh :Is injected 
to a coarser grid level indicated by the subscript 2h by the volume weighted average. 
4 
and the residual is injected by the summation over the fine grid cells, a 
A forcing function F F  is defined on the coarse grid as the difference between the injected 
residual and the Res function (equation (2.7)) calculated on that grid using the injected 
solution. 
FFzh = RIzh - Reszh 
The residual for advancing the solution by one or more time steps on the coarse grid is 
defined by 
R2h = ReSzh -t- FFzh 
(3.3) 
(3.4) 
The linearized equation corresponding to (2.4) on the coarse grid is 
which is solved for one or more time steps by the locally implicit method outlined in the 
previous section. This process is repeated till the coarsest grid solution is completed. The 
coarse grid information is passed back to the finer grids as follows. The difference in the 
injected solution and the computed solution on the coarse grid is interpolated at  the grid 
points of the next finer mesh and is added to the solution which was previously computed 
at  that level. The solution is recomputed at that level for one or more time steps by 
the equations (3.4) and (3.3), using the forcing function FF which already exists at  that 
level. This process is continued till the finest grid level is reached where the R and Res 
functions are the same. This completes one multi-grid cycle. In practice there are many 
variations to this method. One of the effective methods for fluid flows involves the solution 
advancement at  each grid level for only one time step. Multi-grid scheme propagates the 
information from one boundary to another in each cycle and allows it to adjust to the 
boundary conditions quickly to reach a steady-state solution. 
4. Three- Dimensional Diffusion Eauat ion 
As a test case for the performance of the algorithm, a time dependent diffusion equa- 
tion is solved with Dirichlet boundary conditions. 
aU - = vV2u.  i n D :  0 5 z , y , t  5 1 
at 
u = O  on aD 
At t = 0, u = Random numbers in (-1, 1) 
at all grid points 
5 
0 The solution u - 0 as t -+ 00. This problem has components a t  all frequencies supported by the grid and provides the most severe test for the process of convergence of the solution 
to steady state. The discretization chosen is the Euler implicit scheme, given by 
where A is forward difference and 62 is central difference operator. The scheme has no 
limitation on grid sizes, which can be non-uniform in all directions. Constant grid size h 
is indicated here for the presentation of results for the test cases. The scheme (4.2) can be 
rewritten in the A-form as follows: 
uAt 
h2 
R =  -
Locallv imdicit aDDroximation: The implicit equations (4.3) are solved for Au's by one 
or more Gauss-Seidel iterations. This test problem being relatively simple, one Gauss- 
Seidel iteration per time step without any further modification is adequate. It is stable 
for any R and is insensitive to the iteration sweep direction, though 8 symmetric sweeps 
(starting from each of the eight corners of the computational domain) per time step are 
more desirable for difficult problems, particularly with convection. The iteration for (4.3) 
can be written as 
@ 
where the superscript ( ) for Au on the right side is either (m) or (n - 1) depending 
on the iteration sweep direction. A simple iteration would be with rn = I. Au(") = 0 
and the iteration sweep carried out in the directions of increasing indices. The multi-grid 
technique outlined in the last equation is applicable for this problem and is simpler since 
the equation (4.1) is scaiar and linear. In each multi-grid cycle, the number of iterations 
carried out while going from fine to coarse grids are 1,2,3*. etc.. and the same number 
in the reverse process. Figures 1-3 show the convergence of the solution (I from random 
numbers to its asymptotic state u = 0. DELU = M a z  ' u  - uasyl = M a s  jui. Work unit 
is the work required to compute the solution on the finest grid for one time step. For a 
203 mesh, convergence to more than 4 orders of magnitude is achieved in 60 work units 
for a single grid method where as it takes only 18 work units in a multi-grid mode with 
3 grids. With .LO3 mesh it takes 500 iterations for the same level of convergence with a 
single grid while it takes only 21 work units in the multi-grid mode with 4 grids. Even @ 
6 
with 40 x 40 x 80 mesh it takes only 21 work n i t s  with a 4 level multi-grid scheme while 
it takes more than 600 work units by a single grid scheme. Fig. 4 shows the convergence 
of Max /Residue1 against work units for a 403 mesh, which again shows a convergence of 
4 orders of magnitude in 21 work units. 
0 
I 
The locally implicit scheme can also be applied for a group of nodes at  a time instead 
of one node at a time(’). For example, the equations (4.3) a t  a group of eight adjacent 
nodes can be solved with a locally implicit approximation. Fig. 5 shows the convergence 
of such a zonal scheme in both single grid and multi-grid modes. While the convergence 
is slightly better in terms of number of work units than the single point scheme, actual 
computational work is many times larger, since the locally implicit equations for each zone 
have to be solved simultaneously which is time consuming. Thus the single point scheme 
seem to be the cornputationally efficient scheme. 
5. Conclusions 
The locally implicit scheme is a computationally efficient scheme which converges 
rapidly in multi-grid modes for elliptic problems. It has the promise of providing a rapidly 
converging algorithm for steady-state viscous flow problems. 
7 
6. Reference2 
1. Chang, L. C., Kwak, D., Dao, S. C., and Rosen, R., “A Three-Dimensional Incom- 
pressible Flow Simulation Method and Its Application to the Space Shuttle Main 
Engine, Part 1 - Laminar Flow”, AIAA-85-0175, AIAA 23rd Aerospace Sciences 
Meeting, Reno, January 1985. 
2. Spradley, L. W., “Applicatin of Computational iMethods to  Internal and External 
Fluid Flows”, Proceedings of a Workshop on Computational Fluid Dynamics in Aerospace 
Design, UTSI, Tullahoma, T N ,  June 1985. 
3. Mukhenjee, T., Przekwas, A. J., Tam, L. T., Holland, R. L. and Costes, N. C., “A 
Multidomain Global-Model Analysis of SSME, Workshop on Computational Fluid 
Dynamics, YASA Marshall Space Flight Center, Huntsville, AL, April 1986. 
4. Rhie, J., “A General Purpose Navier-Stokes Solver for Internal Propulsion System 
Flows”, AIAA 86-0207, Aerospace Sciences Meeting, January 1986. 
5.  Reddy, K. C., “A Locally Implicit Scheme for Elliptic Problems”, A Workshop on 
Computational Fluid Dynamics, NASA Marshall Space Flight Center, Huntsville, 
AL, April 1986. 
8 
a 
rl 
k 
M 
I 
rl u 
i w 0 
i 
i 
Ll 
M 
I 
d 
U 
3’ 
-I 
w 
11 
I 
? 
I 
'FI 
d 
L 
bb 
I 
d 
U 
h 
M 
P-. d 
ki w 
x 
i 
I 
c 
a 
m 
M 
I I 


A locally implicit method for fluid flow problems

Abstract

Similar works

Full text

Available Versions

NASA Technical Reports Server