An unfactored implicit method for the unsteady, turbulent Navier-Stokes equations is implemented in parallel. The algorithm uses the alternating direction implicit (ADI) method as a preconditioner for a conjugate gradient solution of the linear system at each time step. The sequential nature of the ADI process is tackled by using a transposition method to solve complete ADI sweeps on single processors. Results are presented for a standard pitching aerofoil case solved on the iPSC

Badcock, K.J.

English

Enlighten: Publications

DEPARTMENT OF
AEROSPACE
ENGINEERING
imm
.F-CGS-P:An Unfactored Method in
Parallel
)r Turbulent Pitching Aerofoil Flows
Engineuring
periodicals
K.J. Badcock 
Aero Report 9323
AF-CGS-P:An Unfactored Method in
Parallel
for Turbulent Pitching Aerofoil Flows
K.J. Badcock 
Aero Report 9323
Engh-’ering
PERIODICALS
September 14, 1993
AF-CGS-P:An Unfactored Method in 
Parallel for Turbulent Pitching Aerofoil
Flows
K.J. Badcock *
Aerospace Engineering Department 
University of Glasgow,
Glasgow, G12 8QQ, U.K.
EMAIL: giiaaSG &uk. ac. giasgow. udcf
September 14, 1993
Abstract
An unfactored implicit method for the unsteady, turbulent Navier- 
Stokes equations is implemented in parallel. The algorithm uses the alter­
nating direction implicit (ADI) method as a preconditioner for a conjugate 
gradient solution of the linear system at each time step. The sequential 
nature of the ADI process is tackled by using a transposition method to 
solve complete ADI sweeps on single processors. Results are presented for 
a standard pitching aerofoil case solved on the iPSC.
‘This work wa.s carried out under SERC contract.
Contents
1 Introduction
2 AF-CGS Method
3 The Solution of the linear system
4 Parallel Implementation
5 Test Cases
6 Conclusions
3
4
5
6 
9
10
1 Introduction
The solution of the Navier-Stokes equations is becoming increasingly com­
mon for aeroelastic problems (see [1] for a review). The major computa­
tional work in such an analysis is the solution of the unsteady aerodynamic 
model. CFD research at Glasgow has centred on using upwind schemes 
and fully implicit methods for hypersonic and transonic viscous flows. This 
work is being extended to develop an unsteady 3-D Navier-Stokes solver 
for incorporation into an aeroelastic simulation code.
Two main algorithms have been under development for aerofoil flows. 
Parallel implementation has been a major consideration due to the ex­
tensive memory requirements of 3-D implicit schemes in general and the 
intensive CPU needs of a 3-D unsteady Navier-Stokes solution.
The first algorithm is partially implicit. Normal-direction terms are 
treated imphcitly to help to remove the stability limits arising from the 
fine mesh spacing needed to resolve boundary layers. The terms in the 
streamwise direction are differenced exphcitly. This yields an algorithm 
which has very efficient parallel properties in that the implicit part can 
proceed indepently on separate nodes without the need for any communi­
cation during this phase of the calculation. Parallel efficiencies of around 
ninety percent have been achieved on 16 nodes of the iPSC. However, the 
stability restrictions arising from the streamwise direction can prove re­
strictive and for aerofoil problems the algorithm is too slow unless a large 
number of nodes are available to take advantage of its good parallel prop­
erties. The method has been tested for aerofoil and viscous shocktube 
flows and fuU details can be found in [2].
To achieve a more efficient method it is necessary to examine ways of 
increasing the allowable time step. The conventional way of achieving this 
is to include more implicitness in the time-stepping method. By treating 
all of the spatial terms implicitly the stability limits can be removed en­
tirely and the time step can then be chosen from accuracy considerations.
However, the large sparse linear system that needs to be solved at 
each step represents a significant computational obstacle. Traditionally for 
Navier-Stokes codes this system is solved by an approximate factorisation 
method. There is an inherently sequential nature to this process which 
limits the utility of parallel implementations. This was tackled in [3] l)y 
using a special communication of the data known as transposition.
The factorisation error introduces a practical stability limit and also 
degrades the accuracy of the solution. An unfactored method based on 
an approximate factorisation as a preconditioner for a conjugate gradient 
method was discussed in [4] and [1]. The method proved to be successful in
serial and in this paper a parallel version of the AF-CGS code is developed 
using the transposition idea.
The AF-CGS method is briefly restated in the following two sections. 
The reader is referred to previous publications for a more detailed descrip­
tion of the method and a more thorough demonstration of its capabilities. 
Then the parallel version of the algorithm is given and results for standard 
pitching aerofoil flows are given.
2 AF-CGS Method
The thin-layer Navier-Stokes equations in generalised co-ordinates are 
given by
dt 0^ 01]
ds
where
p
II
1 pu
pv
e
j = .r1
s = •/
pU
PUU +^j!:P 
PVV + iyp
U{e + p)-^tP 
0
,g = •/ -1
(1)
pV
puV -f i]xp
pvV + 1]yP
V(e + p)- i]tp
Here.
_! pmiu^ -f {p/i)m-2^x
pniiv^ + (p/3)in2^y 
pmims -f (p/'i)rn2{^,xu + ^yV)
U = G -|- u 
V = i]t -t- ih-u + Vx V 
mi = C + ^2y
m2 — Vx + 1ly
m3 = ( «2 + t’2 )/2 + Pr-1 (7 ~ 1 )_1 (‘■2 k
and J is the determinant of the Jacobian of the transformation ,r = 
.r(^, 7/,t) and y = yii. pj). Here p, u, v, e, p. Pr. Re. c. 7 denote density, 
the two components of velocity, energy, pressure, the Prandtl numl)er, the 
Reynolds number, the speed of sound and the constant ratio of the speciflc 
heats respectively. The viscosity is composed of a part due to the natural 
viscosity of the fluid and a term to account for turbulence. Sutherland's 
law is used to describe the variation of the fluid viscosity with teiiiperatuie. 
The Baldwin-Lomax model is used to provide a value for the turbulent vis­
cosity. Since none of the flows examined herein involve massive separation 
no modification of the turbulence model is used.
To solve this system of partial differential equations a finite volume 
scheme is used which has various features. For the spatial terms Osher's 
method is used. For general geometries the details of Osher s method are 
described in [5]. A MUSCL interpolation is used to provide second or third 
order accuracy and the Von Albada limiter prevents spurious oscillations 
from occuring around shock waves. Central differencing is employed for 
the viscous terms. Far-held boundary conditions are imposed by Riemann- 
invariants and no vortex correction is applied due to the unsteadyness of 
the how.
The temporal discretisation is based on the backward Euler method. 
An efficient mixed analytic and hnite difference procedure is used to gen­
erate the required .lacobian of the spatial discretisation. The linear system 
obtained is solved by the conjugate gradient squared (CGS) method with 
the alternating direction implicit approximate factorisation providing a 
preconditioner. The reader is refered to [4] [1] for full details. A summary 
is given in the next section.
The problems we consider herein relate to pitching aerofoils. Since 
no deformation of the aerofoils is considered the mesh can be rotated 
rigidly with the aerofoil. Dehning the starting mesh points by {.r = x0}ij, 
{y = yo}i,j the transformed mesh at time tj is given by
{.r = rocos(a{ti) - a0) + 1{1 - co.sfao))},^
{y = rosin(a(ti) - a0) -|- lsin{a0)}i,j
where the aerofoil pitches about the point (1,0), the initial incidence is a0 
and r'n = -vl -f tJn ■ The mesh velocities can be easilv evaluated given 
an analytical form for a. No-slip boundary conditions on the aerofod can 
then be imposed.
3 The Solution of the linear system
One time step of the method can be denoted by
(I + At^ + At^)6 w = -At{ -f i?,,:
r)tn r))r
OR
aw ow
where R^ and Rv are terms arising from the spatial discretisation in the 
^ and 1] directions respectively. The derivatives required on the left-hand 
side of 2 are calculated by an efficient mixed analytic/finite difference 
procedure [6]. The alternating direction implicit (ADI) factorisation of 2 
is
dRr,
dw dw
)Sw
{I+At^)(I+ At^ )6w = -At{R^ + R,,). 
,1...dw
(3)
Equation 3 is much easier to solve than 2 because each of the factors is 
typically block tridiagonal or block pentadiagonal. However, the solution 
to the ADI factored system is not an exact one for equation 2 and in 
practice the factorisation error leads to a stability limit on the time step 
and introduces another source of error to the calculation.
The approach adopted in the present work involves solving equation 2 
to a required tolerance by using a conjugate gradient method. Denoting 
the unfactored linear system 2 to be solved at each time step by
Ax = b (4)
we seek an approximation to A-1 ss C'_1 which yields a system
C~1Ax = C'-^h (5)
more amenable to conjugate gradient methods. The ADI method gives a 
fast way of calculating an approximate solution to 4 or, restating this , of 
forming the matrix vector product
c-1b = X. (6)
Hence, if we use the inverse of the ADI factorisation as the preconditioner 
then multiplying a vector by the preconditioner can be achieved simply by 
solving a linear system with the right-hand side given by the multiplicand 
and the left hand side given the approximate factorisation. The factors 
in C can be diagonalised once at each time step with the row operations 
being stored for use at each midtiplication by the preconditioner.
The exact form of the algorithm for one step of the Navier-Stokes 
solution is
• calculate matrices and diagonalise ADI factors
• calculate updated solution by ADI
• use this solution as starting solution for CGS
• perform CGS iterations until 4 has been solved to required tolerance
4 Parallel Implementation
The major obstacle to an efficient parallel implementation of the AF-CGS 
method is the inherently sequential nature of the ADI procedure. This was 
overcome in [3] by using a transposition of the data to allow complete ADI 
sweeps to proceed independently on each processor. We use this approach 
here although extra communication is required for the present method 
because of the matrix-vector products required in the CGS algorithm.
6
The computational space in generalised coordinates is mapped onto the 
nodes by grouping several mesh lines in both the ^ ^^d the i] directions 
onto one node. Care has to be taken to make sure that ^ hues separated 
only by the wake cut are mapped to the same processor. The computation 
then falls into three phases. First, the matrix is generated and the factors 
are diagonalised. The next phase is the multiplication of a vector by the 
matrix and finally we have multiplication of a vector by the preconditioner 
which reduces to back substititutiou on the diagoiialised factors of the ADI 
factorisation. For each phase data is held on a node for complete lines 
in one direction and the entire computation relating to that direction is 
completed. The data is then communicated so that information for lines 
in the other direction is held on each node and the computation for that 
direction proceeds.
For the calculation of the matrix, the solution from the previous time 
step is initially held on nodes by lines in the ?/ direction. The fluxes and 
derivatives are calculated for each line and the contributions to the residual 
are formed. At this stage the factor
I + d w (7)
is diagonalised. The solution and the residual are then transposed so that 
each node holds information on complete lines in the ^ direction. No 
communication of the components of the matrix dRv/dtu is made. They 
are stored on the nodes on which they are calculated. The calculation 
proceeds so that we have the left hand side of 2 stored in complete ^ and 
1] lines on a node and the solution and residual stored on complete ^ lines 
only.
Next, one step of ADI is performed to obtain a starting solution for 
CGS. This is achieved by performing ADI sweeps in the ^ direction first, 
then transposing the solution obatined and performing sweeps in the // 
direction to yield the updated ADI solution. A multiplication by the 
preconditioner is performed in the same way. To ensure that the CGS 
starting solution and the right hand side of the hnear system are stored in 
the same direction the residual is transposed.
The CGS solution proceeds with inner products and matrix multipli­
cations being recpiired. Contributions to an inner product are calculated 
locally on each node and are then transmitted to a chosen node for sum­
ming. The results is then transmitted back to all the nodes for use. A 
matrix multipUcation is achieved by finding the product in the current di­
rection (?/), transposing the multiplicand and the partial product and then 
finding the product in the other direction. This has the effect of changing 
the direction of the storage. The matrix multipfication is always followed
Xi fluxes and matrix — (eta) 
i
Eta fluxes and matrix 
i
ADI sweep eta direction
l . . (tp)
(xi)
(xi)
(eta) (eta)vector ops
(eta)
(et»- xi)
C(AX) (xi* eta)
vector ops (eta)
(eta) test solution (eta)
CGS solution
Calculate lift/ drag/turbulence (eta)
1
restart
Figure 1: Schf inatic diagram of the imrallel AF-CGS method. Here the operation 
is denoted on the left hand side, the direction of storage is indicated in parenthesis 
and TP denotes when a transposition of the data takes place. The box inset shows 
the flow of the main parts of the CGS algorithm computation.
by a multiplication by the preconditioiier as described above and this re­
verses the direction back again. As a result of this all of the C GS vectors 
lie in the direction which is current on entry to the CGS subroutine. The 
complete sequence of transpositions is shown in figure 1.
Finally, the turbulent viscosity and the lift and drag coeffients are 
calculated in parallel before the next step proceeds.
The parallel code was tested on the Intel Hypercube at the SERC 
Daresbury Laboratory. The test case considered was to calculate thiee _ 
non-dimensional time units of the flow for the pitching aerofoil case 2 de­
fined below. The CPU times are shown in table 1 and indicate that the 
transpositions degrade the parallel efficiency. However, this is the case for 
most parallel preconditioned conjugate gradient solvers and the efficien­
cies are good compared to the incomplete ILU-conjugate gradient solvei 
of [7]. It is anticipated that the efficiency will increase on finer meshes as 
the amount of computation increases on each node. A significant speed 
up is still achieved in the present case. In addition, the large storage re­
quirements of the algorithm have been divided amongst the processors.
The breakdown of the time taken by each main phase of the calculation is
Machine CPU time efficiency
SPARC 10 338 -
4 nodes 187 1.00
8 nodes 106 0.88
16 nodes 68 0.69
Table 1: CPU times in seconds on various machines. The parallel efficiency is 
based on the CPU time for 4 nodes since the problem is too large in terms of 
memory to fit onto one node.
Machine total Jacobian calculation 
and ADI diagonalisation
ADI step CGS solution
SPARC 10 338 8.00 0.31 1.82
iPSC-8 nodes 106^ 2.27 0.21 0.58
speed up 3.2 3.5 1.48 3.14
Table 2; Breakdown of calculation CPU times in .seconds on a parallel and a 
serial machine.
shown in table 2. The net speed up achieved by using the parallel machine 
is 3.2. The main computational work is involved with calculating the jaco- 
bian and diagonalising the ADI factors It is clear that the parallelisation is 
most effective for this phase of the calculation which takes up about eighty 
percent of the total work. The second most intensive part of the work is 
the CGS solution and this too is effectively parallelised even though the 
ADI part of the CGS solution slightly degrades the efficiency. Finally, 
an ADI step is relatively inefficient due to the large communication time 
relative to the computational time involved.
5 Test Cases
The AGARD sub-committee on aeroelasticity defined test cases to act 
as standard flows for computer code evaluation and verification. In this 
paper we present results for two of these test cases for pitching aerofoils.
The motion is defined by the angle of attack as a function of time and the 
centre of rotation .vc which is given herein as a distance along the chord 
as a percentage of the chord length. The angle of attack is defined as
o(r) = am + c\Q.sin{ kr) («)
where r = is the non-dimensional time.
number aerofoil 3/inf Be X lO6 ^ m c>o k ■>‘c
1 NACA0012 0.60 4.8 4.86 2.44 0.1620 0.25
2 NACA64A010 0.796 12.56 0. 1.0 0.204 0.248
Table 3: Test cases examined in this paper.
The cases considered, which are Usted in table 3. were selected because 
detailed pressure distributions are available at a number of points during 
the cycle. All of the results were obtained on all by 33 C-mesh which was 
generated by the Eagle grid generation package. The far field was set at 10 
chords. The AFCGS method was used throughout with the CGS tolerance 
set at 10_1. The number of time steps per cycle for each case is 150. The 
unsteady computations are started impulsively from the converged steady 
solution at the mean angle of attack and the results are shown after two 
cycles of the motion. Detailed comparisons for all the cases are shown for 
the pressure distril)utions at eight separate times during one cycle. Data 
for case 2 was only available for the upper surface. For a more detailed 
discussion of the results the reader is referred to [1]. The detailed pressure 
distribution comparisons are shown in figures 2-3. Excellent agreement 
with the experimental data of [8] is noted.
6 Conclusions
It was demonstrated in [4] and [1] that the AF-CGS method can produce 
efficient and accurate results to attached aerofoil flows. A parallel version 
of the algorithm was presented in this paper and an efficiency comparable 
with other preconditioned conjugate gradient algorithms was achieved. 
Accurate results for standard pitching aerofoil flow's can be achieved in 
less than 15 minutes per period using eight nodes of the iPSC. A drawback 
of implicit methods is the substantial of memory required. The parallel 
version tackles this problem by splitting the memory needs between nodes 
without duplication.
Tackling three-dimensional unsteady aerodynamic problems on the iPSC 
with the AF-CGS method now seems feasible and will for the next major 
task of this work. It is also intended to examine some two-dimensional 
aeroelastic apphcations in collaboration with British Aerospace.
10
a = 5.94°
O' = 6.59°
o = 3.51°
o = 2.67°
a = 6.96°
a = 5.12°
o = 2.44°
Q = 4.27{
Figure 2: —CP vs x/c for case 1. x-experiment. line - computed
11
o = 1.00° Q = 0.74°
a = -1.06° a = -0.73°
a = -1.01° a = -0.59°
a = 0.21° a = 0-870
Figure 3: -Cp vs x/c for case 2. x-experiment, line - computed
12
References
[1] K.J.Badcock. Computation of turbulent pitching aerofoil flows. Tech- 
irical report. G.U. Aero report 9-322, 1993.
[2] K.J.Badcock. A parallelisable partially implicit method for unsteady 
viscous aerofoil flows. Technical report, G.U. Aero report 9312, 1993.
[3] R. Pelz T. Chyczewski, F. Marconi and E. Churchitser. Solution of 
the Euler and Navier-Stokes eciuations on a parallel processor using 
a transposed/Thomas ADIalgorithm. In 11th AIAA Computational 
Fluid Dynamics Conference. AIAA, 1993.
[4] K.J.Badcock. An efficient unfactored imphcit method for unsteady 
aerofoil flows. Technical report, G.U. Aero report 9313, 1993.
[-5] S.Osher and S.R.Chakravarthy. Upwind schemes and boundary con­
ditions with apphcations to Euler equations in general coordinates. ./. 
Cejinp. Phys.. .50:447-481. 1983.
[6] K.J.Badcock. Newton’s method for laminar aerofoil flows. Technical 
report, G.U. Aero report 9310, 1993.
[7] P.J. Wesson. Parallel Algorithms for Systems of Equations. PhD thesis. 
Oxford University Computing Laboratory, 1992.
[8] J.J. Olsen. Compendium of unsteady aerodynamic measurements. 
Technical Report 702, AGARD, 1982.
13


AF-CGS-P: an Unfactored Method in Parallel for Turbulent Pitching Aerofoil Flows. G.U. Aero Report 9323

http://eprints.gla.ac.uk/183010/1/183010.pdf

AF-CGS-P: an Unfactored Method in Parallel for Turbulent Pitching Aerofoil Flows. G.U. Aero Report 9323

Abstract

Similar works

Full text

Available Versions

Enlighten: Publications