A coarse-grid correction algorithm has been implemented into an implicit upwind Euler solver and tested for transonic airfoil problems. The Euler solver uses split-flux formulation and penta-diagonal scalar equations, respectively, for the explicit and implicit operators. The multigrid sequence starts at the fine grid level, then steps down to each coarse grid level to smooth error components using implicit operators. Estimate of residuals can be obtained by two approaches, which differ in the level where the residuals are collected. Both approaches will lead to a work reduction factor of 12 for a Mach 0.75 flow at 2 degrees incidence on a 65x26 grid. The work reduction factor is found to increase proportional to the number of grid levels

Li, Chien-Peng

English

NASA Technical Reports Server

.<  .' 
NASA Technical Memorandum 58276 
Euler Solutions Using an 
Implicit Multigrid Technique 
Chien-peng Li 
November 1986 
(NASA-TM-58276) EULER 3GLUZICNS U S I N G  A N  N87-16792 
I B P L l C c I l  H U L T 1 G E I . D  TECt lNICUE ( N A S A )  1 9  p 
CSCL 01A 
U n c l a s  
G3/02 43750 
NASA 
National Aeronautics and 
Space Ad ministration 
Lyndon B. Johnson Space Center 
Houston, Texas 
https://ntrs.nasa.gov/search.jsp?R=19870007359 2020-03-20T12:53:49+00:00Z
NASA Technical Memorandum 58276 
Euler Solutions Using an Implicit 
Multigrid Technique 
Chien-peng Li 
Lyndon B. Johnson Space Center 
Houston, Texas 
N A S A  
National Aeronautics and 
Space Administration 
Scientific and Technical 
Information Branch 
1986 
CONTENTS 
Section Page 
ABSTRACT . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  1 
INTRODUCTION . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  1 
SPLIT-FLUX EULER FORMULATION . . . . . . . . . . . . . . . . . . . . .  1 
UPWIND FACTORIZATION TECHNIQUE . . . . . . . . . . . . . . . . . . .  3 
MULTIGRID ALGORITHM . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  3 
DISCUSSION . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  4 
CONCLUSIONS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  5 
REFERENCES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  5 
... 
111 
FIGURES 
Figure 
1 
2 
3 
4 
5 
6 
Page 
Schematic of a complete MG cycle. . . . . . . . . . . . . . . . . . . . . . . . . . . 
Comparisons of single- and multigrid results. . . . . . . . . . . . . . . . . . . . 
8 
9 
Grid and Mach contours for M = 0.8, a = 1.25" calculation . . . . . . . . . . . 10 
Comparison of convergence history between single- and multigrid 
calculations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 
Comparison of convergence rate between two estimates 
of residuals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 
Comparisons of MG results using different types of 
residuals collections . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 
iv 
ABSTRACT 
A coarse-grid correction algorithm has been implemented into an  implicit upwind Euler 
solver and tested for transonic airfoil problems. The Euler solver uses split-flux formulation and 
pentadiagonal scalar equations, respectively, for the explicit and implicit operators. The multigrid 
sequence starts at the fine grid level, then steps down to each coarse grid level to smooth error 
components using implicit operators. An estimate of residuals can be obtained by two approaches, 
which differ in the level a t  which the residuals are collected. Both approaches will lead to a work 
reduction factor of 12 for a Mach 0.75 flow at 2" incidence angle on a 65 X 26 grid. The work 
reduction factor is found to increase in proportion to the number of grid levels. 
INTRODUCTION 
The multigrid (MG) method was first advocated in 1980 by Brandt (ref. 1) as a potentially 
efficient technique to solve fluid dynamics equations. Since then, Ni has developed an  innovative 
explicit scheme to solve the Euler equations, and Jameson and his coworkers (ref. 3) have attempted 
to combine the MG strategy with explicit as  well as with implicit schemes. Jameson's latest work 
has shown an  impressive work reduction factor (WRF) as great as 15 for a nonlifting transonic 
airfoil flow on a 128 X 32 grid. These applications share two prominent features; namely, the finite- 
volume formulation and the centered difference approximation to the conventional flux terms. 
Thereby, the error components can be removed and expelled out the computation domain very 
effectively on a succession of grid levels and smoothed by the numerical damping added to the 
difference operators. I t  is of interest to determine whether the MG strategy will work equally well if 
the flux terms are split into subflux terms according to their eigenvalues and approximated by 
upwind schemes. Because the directional preference of signal propagation might be lost by 
collecting residuals on coarse grids and the numerical damping is generally not required to stabilize 
the results near the shock, the basic mechanisms attributable to the removal of error components 
seem to be missing. To resolve the unsettling issues, the author has incorporated the multigrid 
sequencing into an  AD1 factorization Euler solver and reported the findings in a 1984 paper (ref 4). 
The MG algorithm is different from other versions in using the corrections at each level rather than 
the residuals at the fine grid level. Although the algorithm is stable, as well as easy to implement, 
and requires no additional storage for intermediate variables, the improvement in WRF is 
relatively small compared to others. It has been found in this study that a much higher WRF can be 
obtained by introducing local time increment and implicit damping to the algorithm. 
SPLIT-FLUX EULER FORMULATION 
The conservation-law form of the inviscid flow equations is given in generalized coordinates €, 
and q. 
u , + F t + G q = O  (1) 
Each part of the flux term is composed of three subflux terms which are associated with local 
characteristics AS, Aq. They are of the following form. 
1 
+ lel Ye 1 
1/2 - 
ii = t,,u + t Y v ,  ij = qxu + qyv ,  t,, = C;<, c= (t,: + 5) ,etc. 
- -  
and t,,, t, ,Ti,, Ti are treated as local invariants. The three parts corresponding to the characteristic 
values can be cast alternately in the following form. Y Y 
Standard notation is used for flow variables; viz, the density p, the pressurep, and the velocity 
components u and v in Cartesian coordinates, the total internal energy e and (q  = u2 + ~2112, and 
the internal energy e = CUT, which relates t op  and p by the equation of state. The sonic speed is 
c = ( y p / p P ;  y is the ratio of specific heats. The integer is used to simplify the equation form, in 
which each characteristics component can be identified by e = -1 ,O,  1. The conventional matrix of 
characteristics has four components: 
In accordance with the sign of the characteristics, second-order one-sided difference equations are 
used upstream or downstream a t  each grid point. The order of accuracy is the same on the boundary 
points as those in the field. Details of boundary treatment are given in reference 4 and a 
forthcoming publication. 
2 
UPWIND FACTORIZATION TECHNIQUE 
Let Av and A u  denote the unknown correction vectors abbreviated for (Ap, Au, Av, Ae)T and 
(Ap, Apu, Apv, ApdT, respectively. The solution procedure for implicit calculation can be described 
in the following four steps: 
( I +  At6(At - d i 6 ~ t ) A w i , j = ( T - 1 A v ~ .  b J  . (3) 
( I +  At6?A9 - d . 6  I 99 ( S - l T A W ) i , j  
Av. . = (S A w ' ) ~ , ~  
l .J  
where r refers to the residual vector obtained from the explicit operation in equations (1) and (2) 
after replacing the derivatives of flux terms by their difference formulas. The subscripts i and j 
denote the spatial location of a grid point network ranging from i = 1 to imax a n d j  = 2 tojnax.  
The operators 6t and 6 are for upwind difference and 6t5' 6,, are for centered difference associated 
with a damping coefficient, d,. The difference formulas stop switching sides if the magnitude of the 
characteristics is smaller than a prescribed tolerance in order to reduce the truncation errors near 
the stagnation and sonic points. The notation used and the expanded form of the second and third 
equations in equation (3) are given in reference 5. The step increment is generally At = CFL - min 
(AEjlAJe, Aq/IA,le) where CFL is the Courant number. It must be equal to or less than 0.25 if 
implicit operations are skipped. A global At implies that one value of time increment is used for the 
entire domain; a local At is obtained at each point and used at that point. 
9 
MULTIGRID ALGORITHM 
The implicit procedure for single-level grid calculations may be summarized as follows: 
k + 1  k 
l * J  1.J 
L A v .  . = r .  . 
vk?' = vk . + Au. . k+ 1 
b J  1.J 1.J 
where L is the operator representing equation (3) and ri , j ,  Av,, are called the residual and 
correction vectors, respectively. On a coarse grid with spacing twice as large as in all coordinates, 
equation (4) becomes 
L Av2, = W, 2h rh or 
(4) 
2h = W, Av, (5) 
3 
On a still coarser grid, 
etc., where W is an operator representing the weighted average of the residuals and corrections on 
the grid points surrounding the coarser grid point, and T is a n  operator transferring directly the 
values from fine to coarse grid. The correction obtained a t  each level is interpolated back to the fine 
grid by 
h vh = vh + Izh A u ~ ~  (6) 
The implicit operator Z is a linear function on the boundaries and bilinear inside the domain. The 
weighted operator W serves a very important role to eliminate high frequency modes before passing 
the residuals to coarser grids. Equation ( 5 )  is in a general form so that the residuals may be 
obtained from either r or Av, or from both. For example, letting 
4h - 4h etc. W, 2h = T Y  and W Z h  - TZh 
The present MG strategy uses the explicit residuals only once on the fine grid level, then uses the 
implicit corrections on other grid levels. A complete cycle of the MG algorithm is illustrated in 
figure 1. 
DISCUSSION 
To help evaluate the performance of the implicit MG algorithm, three parameters (Ap/p)rnax, 
number of supersonic points NSUP, and stagnation pressure Pt  are used for monitoring the 
convergence rate, whereas Pt and Mach contours are used for checking solution accuracy. Several 
grids (113 X 34,97 X 34,65 X 26,65 X 22) have beeen considered, but most of the results shown 
here are based on a 65 X 26 grid and for the NACA 0012 airfoil. The baseline calculation uses 
global At with CFL = 10, as it encounters stability difficulties with local At even with CFL = 2. 
Figure 2 shows the convergence history and results of a non-lifting case. The magnitude of (Aplp) in 
the MG calculation is typified by a plateau two orders higher than that in the single-grid (SG) 
calculation. Nevertheless, the MG calculation is at  least four times as  efficient as the single-grid 
calculation. The WRF is close to 7 ifoptimal values of CFL = 4 and di = 1 are used instead. Figure 
3 shows the 65 X 26 grid and the Mach contours for a Mach 0.8 flow at 1.25" incidence angle. Figure 
4 is a comparison of the convergence rate between the calculations made from the MG algorithm 
with CFL = 2 local A t  and from the single-grid algorithm with CFL = 10. The WRF is estimated to 
be 7 for this case and could have been 13 if CFL = 4 had been used. Note the sharpness of the shock 
above the airfoil and the compression wave underneath. The accuracy of the current results seems 
to be on a par with the conventional SG solution on a 128 X 32 grid. 
4 
If the MG algorithm uses the explicit residuals on coarse grids; viz, 
Pmaz 
Pmin 
NSUP 
Lhv = W, 2h rh 
2h 
39.46 40.11 38.96 38.72 
15.01 14.59 14.78 
86 88 88 
the convergence rate is nearly equivalent to that of the current algorithm. Shown in figure 5 is the 
comparison between the two. A comparison of different types of residual collection is shown in 
figure 6. It is seen that the use of local time step increment substantially accelerates the 
convergence rate. 
In the following table the converged values of pressures on the surface and of NSUP are 
listed. The slight differences seen are partly caused by the peculiar behavior of the solution 
technique (eq. (3)), whereby the factorization error affects the converged solutions. 
Work unit Wu 1000 200 116 
M = 0.75, a = 2"; 65 X 26,33 X 14,17 X 8 grids 
CONCLUSIONS 
The current MG algorithm is capable of enhancing the stability and the efficiency of the AD1 
factorization technique, and of retaining the excellent shock-capturing capability rendered by the 
split-flux formulation. Satisfactory results are obtained from a 65 X 26 grid for Mach 0.75 flow at 2" 
incidence angle with a work reduction factor equal to 12. The apparent shortcoming of the 
algorithm is in the difficulty of lowering the maximum incremental vector further after a 
reasonable number of iterations. 
REFERENCES 
1. Brandt, A.: Multi-Level Adaptive Computations in Fluid Dynamics. AIAA J., vol. 18, 
no. 10, Oct. 1980, pp. 1165-1172. 
2. Ni, R. H.: A Multiple Grid Scheme for Solving the Euler Equations. AIAA J., vol. 20, 
no. 11, Nov. 1982, pp. 1565-1571. 
3. Jameson, A.; and Yoon, S.: Multigrid Solution of the Euler Equations Using Implicit Schemes. 
AIAA Paper 85-0293, AIAA 23rd Aerospace Sciences Meeting, Jan. 14-17,1985. 
4. Li, C. P.: A Multigrid Factorization Technique for the Flux-Split Euler Equations. 
Proceedings of the Ninth Inter. Conf. on Numerical Methods for Fluid Dynamics, Lectures 
in Physics, vol. 218, Springer-Verlag, Jan. 1985, pp. 375-379. 
5 
5. Li, C. P.: A Finite-Difference Method for Solving Transient Viscous Flow Problems. 
AIAA J., vol. 23, no. 5,1985, pp. 659-668 . 
6 
wu 0.5 0.5 0.125 
E = EXPLICIT OPERATOR 
I IMPLICIT OPERATOR 
Figure 1.- Schematic of a complete MG cycle. 
7 
10-1 
10-2 
- 
8 
P 
2 10-2 
3 
10-4 
lo-' 
?A = 0.75 a = 0' 
65 x 26 GRID 
3 GRID LEVELS 
- LOCAL CN = 2 
\ P - - - - 1 GRID LEVEL 
AP 
GLOBAL CN = 10 \ 
\ 
\ 
\ '-,- 4 #-.# -- \ 
I 
I \ 
100 
n 
z 
10 2 
1 
0 200 400 600 800 1000 
wu 
XICHORD 
-l.50,0 0;2 0.,4 0;6 0;s 1;O 
-1.0 3 LEVELS 
XICHORD 
0;2 0: 0;6 0;s 1;O 
-1.0 1 
1 LEVEL 
-0.5 
8 0.0 
0.5 
1 .o 
1.5 1 
n 
l.ol 1.5 
Figure 2.- Comparisons of single-grid and multigrid results. 
8 
Figure 3.- Grid and Mach contours for M = 0.8, Q = 1.25" calculation. 
9 
XICHORD 
- l , 5 i  0.0 0;2 0;4 0; 0.8 , 1 .o ,
-1.0 
-0.5 
n 
0 
0.0 
0.5 
1 .o 
-1.0 
-0.5 
n u 
0.0 
0.5 
1 .o 
wu = 100 
XICHORD 
0.0 0.2 0.4 0.6 0.8 1 .o 
-1.0 
-0.5 
a. 
0 
0.0 
0.5 
1 .o 
-1.5 
-1.0 
-0.5 
n 
0 
0.0 
0.5 
1 .o 
XICHORD 
1 0.2 0.4 0.6 0.8 1 .o 
1 1 1 J 
Wu = 500 
Figure 4.- Comparison of convergence history between single-grid and multigrid calculations. 
10 
a 
a 
a 
a 
XICHORD 
-1.5 0.0 1 4 ' 1  0.4 0.6 0.8 1 .o 
-1 .o 
-0.5 
a 
0 
0.0 
0.5 
1 .o 
XICHORD 
0.0 "2 0;4 0;6 Oi8 1;o 
-1.5 
1.5 ' 
XICHORD 
0.0 
-1.0 
-0.5 
n 
0 
0.0 
0.5 
1 .o 
r Wu 900 
Wu = 167 
LOCAL CFL = 2 
MULTIGRID 
1.5 J 
Figure 4.- Concluded 
wu = 1200 
GLOBAL CFL = 10 
SINGLE GRID 
11 
XlCHORD 
0: 0.; 0;s 0;s 1; 
X/CHORD 
-1.5 0.0 0;2 0: 0.,6 0;8 1.; 
-1 .o 
-0.5 
a 
0 
0.0 
0.5 
1.0 y wu = 100 
XICHORD 
0.0 0.2 0.4 0.6 0.8 1 .o 
-1.5 i 1 i 1 1 1 
-1.0 
-0.5 
n u 
0.0 
1.0 1 Wu = 50 
X/CHORD 
0.0 0.6 0.8 1 .o 
-1.0 
-0.5 
a 
0 
0.0 
0.5 
1.0 f wu = 100 
Figure 5.- Comparison of convergence rate between two estimates of residuals. 
12 
XICHORD 
-1.5 0.0 0; 0; 0;s 0;s 1.; 
XICHORD 
0.6 0.8 1 .o 
1 J 
-1.0 
-0.5 
a 
0 
0.0 
0.5 
1 .o 
-1.5 
-1 .o 
-0.5 
8 0.0 
0.5 
1 .o 
1.5 
-1.0 
-0.5 
a 
0 
0.0 
0.5 
1 .o 
XICHORD 
I 0.2 0.4 0.6 0.8 1 .o 
1 I 1 1 b 
XICHORD 
-1.5 Of 0;2 0.4 , 0: 0;8 1 , .o 
wu = 200 
AV 
-1.0 
-0.5 
8 0.0 
0.5 
1 .o 
wu = 200 
r 
1.5 
Figure 5.- Concluded 
13 
X/CHORD 
0.0 0.2 0.4 0.6 0.8 1 .o 
-1.5 I I I 1 1 
Wu = 133 
1.5 ’ 
X/CHORD 
0.0 0.2 0.4 0.6 0.8 1 .o 
-1.5 
-1.0 
-0.5 
8 0.0 
0.5 
1 .o 
1.5 
Wu = 466 
(b) r2,, = T F  huh; CFL = 10, global. (a) r2h = W r  Avh; CFL = 4, local. 
Figure 6.- Comparisons of MG results using different types of residual collections. 
14 
1. Report No. 
NASA TM-58276 
Lyndon B. Johnson Space Center 
Houston, Texas 77058 
1. Government Actcuion No. 3. R h h t ' 8  CrUbg No. 
11. Contract or Grant No. I 
4. Titlo md SuMitlc 
EULER SOLUTIONS USING AN IMPLICIT MULTIGRID TECHNIQUE 
7. Authorlr) 
Chien-peng Li 
9. hrforming Organization N a m  and Address 
13. Type of Repon md Period Covered 
Tec hn i cal Memorandum 12. Sponsoring Agency Name m d  Address 
National Aeronautics and Space Administration 
Washington, D.C. 20546 14. Sponsoring Apncy codc 
5. Report Date 
6. M w m i n g  Orwniratia, co6 
November 1986 
S-556 
8. Performing Orpniratim R ~ p a t  No. 
10. Work Unit No. 
I 
15. Supplementary Notes 
D. Scfuriry Classif. (of this report) 
Unclassified 
Presented at the Second Copper Mountain Nultigrid Conference 
March 31-April 3,  1985 
20. Security Classif. (of this page) 21. NO. of Pages 22. Rice' 
Unclassified 18 
IS. Abstract 
A coarse-grid correction algorithm has been implemented into an implicit upwind Euler solver 
and tested for transonic airfoil problems. The Euler solver uses split-flux formulation and 
penta-diagonal scalar equations, respectively, for the explicit and implicit operators. The 
multigrid sequence starts at the fine grid level, then steps down to each coarse grid level 
to smooth error components using implicit operators. Estimate of residuals can be obtained 
by two approaches, which differ in the level where the residuals are collected. Both 
approaches will lead to a work reduction factor of 12 for a Mach 0.75 flow at 2 degrees 
incidence on a 65x26 grid. The work reduction factor is found to increase proportional to 
the number of grid levels. 
7. Key Words (Suggested by Authorls)) I 18. Distriburion Statement 
Finite-difference 
Mu1 tigrid 
Imp1 icit method 
Transonic airfoil flow 
*For sale by the National Technical Informatioh Service, Springfield. Virginia 22161 


Euler solutions using an implicit multigrid technique

Abstract

Similar works

Full text

Available Versions

NASA Technical Reports Server