T 3.. .

t, . . A .
:5 Titltf
.u. x .

.2.
t. 4..

..
it .

1 .111.
. 96.1.

fit. :33" o

.

w
m
a.

.3 Jr
. . .1...

.

 

“F

K.)

C)

9‘)

This is to certify that the

dissertation entitled

GRID-QUALITY MEASURES FOR ERROR ESTIMATION
AND SOLUTION-ADAPTIVE MESH REFINEMENT IN CFD

presented by

Xubin Gu

has been accepted towards fulfillment
of the requirements for

Ph.D. Mechanical Engineering

 

 

degree in

Wjor professor '

Date %7 241 2&02

MS U is an Affirmative Action/Equal Opportunity Institution 0-12771

 

LIBRARY
Michigan State
University

 

 

 

PLACE IN RETURN Box to remove this checkout from your record.
TO AVOID FINES return on or before date due.
MAY BE RECALLED with earlier due date if requested.

 

DATE DUE DATE DUE DATE DUE

 

@g’z‘fUfi’flDUd

 

 

 

 

 

 

 

 

 

 

 

 

 

 

6/01 cJCIRC/DateDue.p65~p. 15

GRID-QUALITY MEASURES FOR ERROR ESTIMATION
AND SOLUTION-ADAPTIVE MESH REFINEMENT IN CFD
By

Xubin Gu

A DISSERTATION
Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of
DOCTOR OF PHILOSOPHY
Department of Mechanical Engineering

2002

ABSTRACT
Grid-Quality Measures for Error Estimation
and Solution-Adaptive Mesh Refinement in CPU
by

Xubin Gu

As computational fluid dynamics (CFD) gains increasing acceptance as a design
and analysis tool, so has the demand on its accuracy and reliability. One major
source of error is from the grid or mesh. This research seeks to develop
methods to estimate error obtained on poor-quality meshes and to examine how

best to refine a mesh to obtain grid-independent solutions.

On error estimation, it is hypothesized that grid-induced error at any location is a
function of a set of grid-quality measures. Based on this hypothesis, six grid-
quality measures were developed and evaluated. These measures, which can
be applied to both structured and unstructured meshes, link information about the
geometry of the cells in a mesh to the scalar, vector, and tensor nature of the
solution. Four measures were developed to gauge mesh resolution, and two
measures were developed to gauge mesh smoothness. Of the measures that
gauge mesh resolution, one is based on the directional change in the velocity
vector. The other three are based on the local gradient of flow variable. Of
these, one compares a length scale representing the curvature of the velocity

field and an effective cell dimension in the direction of change. The second

compares the spectral radius of the deformation tensor scaled by the local
velocity and a cell dimension along its eigenvector. The third compares the
relative changes in effective velocity perpendicular to the streamline. For the two
measures that gauge the mesh smoothness, one considers change in effective
grid spacing among all neighbors. The other considers change in effective grid
spacing only along flow direction. The six measures developed were evaluated
by establishing and interrogating a database, made up of three test problems
with complex flow patterns. For each test problem, a series of CPD solutions
were generated by using a variety of meshes that range from “poor" to “good”
quality. This evaluation showed the measures developed to have correlation with

the solution error defined with respect to the grid-independent solution.

On mesh refinement, a series of numerical experiments were performed to study
the effectiveness of popular error indicators, including grid-quality measures
developed in this study. Results show that even when the exact relative and
absolute errors are used as error indicators, h-type mesh refinements could still
be highly inefficient. Results also show that grid refinement is needed in two
regions - a region referred to as error-source (ES) region, where error initiates
but is not yet apparent and a region referred to as error-apparent (EA) region,
where the actual error is high. In this study, a new adaptation strategy was
developed that refines the mesh in both ES and EA regions. This strategy was
shown to be more effective in reducing relative error than when the exact

absolute and relative errors are used as error indicators.

Copyright by
Xubin Gu

2002

ACKNOWLEDGEMENTS

I would like to express my gratitude to my advisor, Prof. Tom Shih, for his
support, patience, and encouragement throughout my graduate studies. His
technical and editorial advice was essential to the completion of this dissertation
and has taught me innumerable lessons and insights on the workings of

academic research in general.

My thanks also go to my co-advisor, Prof. Harold Schock, and the members of
my committee, Prof. Andre Benard, Prof. C. Y. Wang, and Prof. Mei Zhuang for
reading previous drafts of this dissertation and providing many valuable
comments that improved the presentation and contents of this dissertation.
Thanks are also owed to personals at the CFD laboratory for their whole-hearted

cooperation and friendship.

Last, but not least, I would like to thank my parents for their dedication and the
many years of support during my studies that provided the foundation and
motivation for this work. My wife, Wenyue, receives my deepest gratitude and
love for her understanding and love during the past few years. Her support and

encouragement was in the end what made this dissertation possible.

This work is supported in part by Department of Energy. DaimlerChrysler, EPA,

and For Motor Inc.

TABLE OF CONTENTS

LIST OF TABLES ........................................................................... ix
LIST OF FIGURES ........................................................................... xii
NOMENCLATURE ........................................................................ xvi
CHAPTER 1
INTRODUCTION ........................................................................... 1
1.1 Background ........................................................................... 1
1.2 Literature Survey .................................................................. 5
1.2.1 Grid-Quality Measures ......................................................... 5
1.2.2 Error Estimation .................................................................. 9
1.2.3 Grid Adaptation .................................................................. 15
1.3 Objectives and Outline ......................................................... 20
CHAPTER 2
DEVELOPMENT OF GRID-QUALITY MEASURES .............................. 23
2.1 Conventional Grid-Quality Measures ....................................... 23
2.2 Characteristics of Ideal Grid-Quality Measures ........................ . ..... 27
2.2.1 Basic Classification of Measures ....................................... 27
2.2.2 Attributes of Ideal Measures ................................................ 28
2.3 Newly Developed Grid-Quality Measures: Schematic and Formula 29
2.3.1 Directional Resolution ......................................................... 29
2.3.2 Gradient Resolution ......................................................... 31
2.3.3 Defamation Tensor Resolution ....................................... 34
2.3.4 Directional Smoothness ................................................ 36
2.3.5 Streamline Smoothness ................................................ 37
2.3.6 Relative Velocity Difference ................................................ 41
CHAPTER 3
APPROACH TO EVALUATE GRID-QUALITY MEASURES ..................... 42
CHAPTER 4
TEST PROBLEMS AND SOLUTIONS GENERATED .............................. 46
4.1 Basis of Problem Selection ......................................................... 46
4.2 Test Problem: 2-D Flow in a Square Chamber with Two Inlet Ducts
and One Exit Duct .................................................................. 47
4.2.1 Problem Description and Formulate ....................................... 47
4.2.2 Grid-Independent Solution Generation .............................. 49
4.2.3 Imperfect Grids .................................................................. 53

4.3 SD Test Problem I: 3-D Flow in a Cubic Chamber with Four Curved

vi

Inlet Ducts and One Exit Duct ................................................ 55

4.3.1 Problem Description and Formulate ....................................... 55
4.3.2 Grid-Independent Solution Generation .............................. 57
4.3.3 Imperfect Grids .................................................................. 61
4.4 3D Test Problem ll: 3-D Flow in a Cubic Chamber with Two
Opposing Inlet Ducts and Two Exit Ducts .............................. 62
4.4.1 Problem Description and Formulate ....................................... 62
4.4.2 Grid-Independent Solution Generation .............................. 64
4.4.3 Imperfect Grids .................................................................. 67
CHAPTER 5
EVALUATION OF GRID-QULITY MEASURES ....................................... 69
5.1 Database: Construction and Organization .............................. 69
5.2 Evaluation of Error with One Measure at a Time .............................. 75
5.2.1 Evaluation of Directional Resolution ....................................... 75
5.2.1.1 Correlation with Different Errors ....................................... 75
5.2.1.2 Evaluation of Generality for all ZDISD Simulations ............ 79
5.2.1.3 Evaluation of Different Data Descriptor - Mean, Median,
Mode .................................................................. 81
5.2.1.4 Analysis of Statistics in Data Partition - Variance,
Skewness, Kurtosis ................................................ 84
5.2.1.5 Detailed Data Distribution in Each Partition .. .................... 87
5.2.2 Evaluation of Remaining Measures ....................................... 89
5.2.2.1 Resolution-Type Measures ....................................... 90
5.2.2.2 Smoothness-Type Measures ....................................... 92
5.2.2.3 Cell Shape-Type Measures ....................................... 94
5.2.3 Summary of Evaluation One Measure at a Time ..................... 94
5.3 Evaluation of Error with Multi-Measures ....................................... 96
5.3.1 Evaluation of Errors with Two Measures at a Time ............ 96
5.3.1.1 Contour Pattern ......................................................... 96
5.3.1.2 Mean Value in Partition ................................................ 98
5.3.1.3 Distribution Feature in Each Partition .............................. 99
5.3.2 Analysis of Correlation Coefficient ....................................... 101
5.3.3 Evaluation of Errors with Three Measures at a Time ............ 104
5.4 Summary ........................................................................... 106
CHAPTER 6
APPLICATION OF GRID-QUALITY MEASURES AS ERROR INDOCATORS
............................................................................................. 110
6.1 Introduction ........................................................................... 110
6.2 The Error Pattern of Previous 2D Test Problem .............................. 111
6.3 Application of Grid-Quality Measures as Error Indicator ............ 115
6.3.1 20 Moving Wall Cavity Flow ................................................ 116
6.3.1.1 Problem Description and Formulate .............................. 116

vii

6.3.1.2 Grid-Independent Solution Generation .............................. 117

6.3.1.3 Error Generation Check ................................................ 119
6.3.1.4 Example #1: Uniform Coarse Grid .............................. 122
6.3.1.5 Example #2: Non-Uniform Coarse Grid .............................. 127
6.3.2 20 Chamber Flow ......................................................... 131
6.3.2.1 Problem Description and Formulate .............................. 131
6.3.2.2 Grid-Independent Solution Generation .............................. 132
6.3.2.3 Error Generation Check ................................................ 135
6.3.2.4 Example #1: Uniform Coarse Grid .............................. 138
6.3.2.5 Example #2: Non-Uniform Coarse Grid .............................. 142
6.4 Summary ........................................................................... 145
CHAPTER 7
DIFFERENTIATING BETWEEN SOURCE AND LOCATION OF ERROR
FOR SOLUTION-ADAPTIVE MESH REFINEMENT .............................. 146
7.1 Introduction ........................................................................... 146
7.2 Fix-Percentage Adaptation by (Pen and Exact Solution Error ............ 148
7.2.1 20 Moving Wall Problem ............ . ................................... 148
7.2.2 20 Chamber Problem ......................................................... 161
7.3 Error Source Region and Error Apparent Region Model ............ 166
7.4 Error Source Indicator (ESI) ................................................ 170
7.4.1 Two new developed ESI ................................................ 170
7.4.2 Test by 2D Moving Wall Problem ....................................... 173
7.5 Fixed-Criterion Adaptation ........................................... . ...... . ...... 181
7.6 New Adaptation Strategy Using ESR/EAR Concept ..................... 187
7.6.1 20 moving wall case ......................................................... 187
7.6.2 20 Chamber Problem ......................................................... 192
CHAPTER 8
SUMMARY AND CONCLUSION ......................................................... 199
APPENDEX A
FLUENT Code Validation .................................................................. 203
APPENDEX B
PUBLICATIONS ASSOCIATED WITH THIS WORK .............................. 209
BIBLIOGRAPHY ........................................................................... 210

viii

LIST OF TABLES

Table 1. Uniform grids used to get grid-independent solution for 20 test

problem .................................................................. 52
Table 2. Average errors between successive levels for ZD test problem . 52
Table 3. Meshes with different aspect ratio ..................................... 54
Table 4. Meshes with different smoothness ..................................... 54
Table 5. Meshes with different cell Shape ..................................... 54
Table 6. Uniform grids used to get grid-independent solution for 30 test

problem I ......................................................................... 58
Table 7. Average error between successive levels for 30 test problem I . 59
Table 8. Meshes with different aspect ratio ..................................... 61
Table 9. Meshes with different smoothness ..................................... 61
Table 10. Meshes with different cell shape ..................................... 61
Table 11. Uniform grids used to get grid-independent solution for SD test

problem II ......................................................................... 65
Table 12. Average errors between successive levels for 3D problem II . 65
Table 13. Meshes with different aspect ratio ....................................... 67
Table 14. Meshes with different smoothness ....................................... 68
Table 15. Meshes with different cell shape ....................................... 69
Table 16. Correlation coefficient between the relative error and new
Measures ........................................................................... 103

Table 17. Average error level in uniform meshes .............................. 112
Table 18. Uniform grids used to get grid-independent solution ............ 117

ix

Table 19

Table 20

Table 21.

Table 22.

Table 23.

Table 24.

Table 25.

Table 26.

Table 27.

Table 28.

Table 29.

Table 30.

Table 31.

Table 32.

Table 33.

Table 34.

Table 35.

Table 36.

Table 37.

Table 38.

Table 39.

Convergence process to get grid-independent solution ............ 118
Convergence process to get grid-independent solution ............ 120
Uniform grids used to get grid-independent solution ............ 133
Average relative errors between successive levels ..................... 133
Average relative errors on each coarse mesh ..................... 135
Fixed-percentage adaptation process by using (boa ..................... 148
Fixed-percentage adaptation process by using Er.a ..................... 153
Comparison of error reduction ................................................ 155
Fixed-percentage adaptation process by using Eab ..................... 159
Fixed-percentage adaptation process by using ¢oR ..................... 161
Fixed-percentage adaptation process by using IE... ..................... 164
Fixed-percentage adaptation process by using Eab ..................... 165
Fixed-percentage adaptation process by using (I). ..................... 176
Fixed-percentage adaptation process by using ll)" ..................... 179
Fixed-criteria adaptation process by using lion. 1% = 5 ............ 182
Fixed-criteria adaptation process by using Em, Illm = 5.0% ............ 182
Fixed-criteria adaptation process by using Eab, ¢m = 8.0x10‘5 ...183
Fixed-criteria adaptation process by using (1)., W = 5.x1043 ............ 184
Fixed-criteria adaptation process by using on, (I’m = 2.x10‘4 ...184
Alternate adaptation A by using I»: and 4m, fixed-percentage=25°/o
.................................................................................... 1 88

Alternate adaptation B by using 4).. and (bun. set ¢..=1. x104, ¢oR=5

.................................................................................... 190
Table 40. Fixed-percentage adaptation process by using in ..................... 193
Table 41. Fixed-percentage adaptation process by using ll)“ ..................... 195

Table 42. Alternate adaptation A by using Ill" and (ion. fixed-percentage=25%
.................................................................................... 196

Table 43. Alternate adaptation B by using (Pu and can, set ¢..=1. x10“, ¢DR=5
.................................................................................... 1 97

xi

LIST OF FIGURES

Figure 1. Accuracy of CPD results ..............................................
Figure 2. Schematic for directional resolution om ............................
Figure 3. Schematic for gradient resolution than .....................................
Figure 4. Schematic for defamation tensor resolution com ...................
Figure 5. Schematic for directional smoothness ops ............................
Figure 6. Schematic for streamline smoothness <l>ss ............................
Figure 7. Comparison of ¢ss and ops ..............................................
Figure 8. Schematic for error definition ..............................................

Figure 9. Geometry and flow conditions of the 20 test problem ..........

Figure 10. Velocity profile comparison between level k and level k+1

Figure 11. Average error level between successive meshes ...................
Figure 12. Flow pattern of grid-independent solution in 2D case ..........
Figure 13. Geometry description of 30 test problem I ............................
Figure 14. Flow conditions of 3D test problem I .....................................
Figure 15. Average error level between successive meshes ...................
Figure 16. Flow pattern of grid-independent solution for 30 test problem I
Figure 17. Geometry description of 3D test problem ll ............................
Figure 18. Flow conditions of 30 test problem ll .....................................

Figure 19. Average error level between successive meshes ...................

Figure 20. Flow pattern of grid-independent solution for 3D problem ll

xii

.. 31

..34

.. 35

.. 37

.. 39

.. 41

.. 43

.. 48

.51

52

53

56

57

59

60

63

64

66

.67

Figure 21.

Figure 22.

Figure 23.
Figure 24.
Figure 25.
Figure 26.
Figure 27.
Figure 28.
Figure 29.
Figure 30.
Figure 31 .

Figure 32.

Figure 33.
Figure 34.
Figure 35.
Figure 36.
Figure 37.
Figure 38.

Figure 39.

Figure 40

Summary of grid-quality study ....................................... 74
Correlation between (bun and errors, Partition by Error
(a) Eab = F(tlii ). (b) E", = F(tl>i ). (0) ME...) = F01)i ) ..................... 77
Correlation between then and different errors, Partition by ¢DR . 78
Generality in all simulations ................................................ 80
Comparison of different data descriptors .............................. 83
Comparison of different statistical parameters ..................... 86
Detailed data distribution in selected partitions ..................... 88
Evaluation of resolution-type measures .............................. 92
Evaluation of smoothness-type measures .............................. 93
Evaluation of skewness ........ . ....................................... 94
Er = “one, ops) ............. . .................................................... 97
Ere = f(¢on. ops) .................................................................. 99
Detailed data distribution in selected partitions ..................... 100
Detailed data distribution in selected partitions ..................... 105
Average error level in uniform meshes .............................. 112
Error pattern on uniform meshes of ZD case ..................... 114
Schematic of cavity-driven problem studied .............................. 116
Convergence process to get grid-independent solution ............ 118
Flow pattern of grid-independent solution, obtained
by using a mesh with 32,768 aspect ratio unity cells. ............ 119
Convergence process of the 20 moving wall cavity flow ............ 120

xiii

Figure 41.
Figure 42.
Figure 43.
Figure 44.
Figure 45.
Figure 46.
Figure 47.
Figure 48.
Figure 49.
Figure 50.
Figure 51.
Figure 52.
Figure 53.
Figure 54.
Figure 55.
Figure 56.
Figure 57.
Figure 58.
Figure 59.
Figure 60.
Figure 61.

Figure 62.

Error pattern on uniform meshes of 2D cavity case ............ 121
Uniform example mesh ................................................ 122
Measure patterns for example #1 ....................................... 124
Exact relative error pattern for example#l .............................. 126
Non-uniform examples grid ................................................ 127
Measure patterns for example #2 ....................................... 128
Cell aspect ratio pattern for example #2 .............................. 129
Exact relative error pattern for example #2 .............................. 130
Schematic of 20 chamber flow problem .............. . ............... 132
Average error level in uniform meshes .............................. 134
Flow pattern of grid-independent solution .............................. 134
Average error level on coarse meshes .............................. 135
Error pattern on uniform meshes of 2D chamber case ............ 137
Uniform mesh for example #1 ....................................... 138
Measure patterns for example #1 ....................................... 140
Exact error pattern of example #1 mesh .............................. 141
Non-uniform example #2 mesh ....................................... 142
Measure patterns for example #2 ....................................... 143
Exact error patterns for example#2 ....................................... 144
Schematic for isotropic subdivision ....................................... 147
Adaptation by (Don using fixed-percentage method ..................... 149
Adaptation by Era using fixed-percentage method ..................... 154

xiv

Figure 63.
Figure 64.
Figure 65.
Figure 66.
Figure 67.
Figure 68.
Figure 69.
Figure 70.
Figure 71.
Figure 72.
Figure 73.

Figure 74.

Figure 75.

Figure 76.

Figure 77.

Figure 78.

Figure 79.

Figure 80.

Figure 81.
Figure 82.
Figure 83.

Figure 84.

Exact absolute error pattern for solution on starting mesh ............ 158
Adaptation by Eat, using fixed-percentage method ..................... 160
Adaptation process in 20 chamber test case ..................... 162
Exact relative and absolute error ...................................... 163
Definition of terms in on (Eqs. (6.2)-(6.7)). .............................. 172
Evaluation of o. and Ill" ......................................................... 173
4). pattern ........................................................................... 174
Adaptation by (I). using fixed-percentage method ..................... 175
on pattern .................................................................. 176
Adaptation by 4).. using fixed-percentage method ................ 178
Summary of fixed-percentage adaptation .............................. 180
Summary of fixed-criterion approach ....................................... 185
Adaptation by using fixed-criterion method .............................. 186
Adaptation by alternate approach A ....................................... 189
Adaptation by alternate approach B ....................................... 191
Efficiency comparison ......................................................... 192
4). pattern ........................................................................... 193
Adaptation by ¢i . ¢ii ......................................................... 194
(bit pattern .................................................................. 195
Efficiency comparison ......................................................... 196
Adaptation by alternate adaptation A, B .............................. 197
Efficiency comparison ......................................................... 198

XV

NOMENCLATURE

>>

O. Q. Q;

e,E

Cell area (m2)

Deformation tensor

Vector between centers of two cells
Diagonal length of cell (m)

Diameter of chamber inlet or exit duct (m)
Solution error

Edge length of cell

Error Apparent

Error Source

Height of cavity or chamber (m)

Edge length of cell (m)

Length of cavity or chamber (m)
Number of all neighboring cells
lnscribe radii of circles or spheres
Circumscribe radii of circles or spheres
Stress tensor

X-component velocity (m/s)
Y-component velocity (m/s)

Z-component velocity (m/s)

Velocity vector (m/s)

Velocity magnitude (m/s)

xvi

Vector between centers of two cells

Measure or error indicator

Dynamic viscosity (Kg/m-s)

Density (kg/m3)

Spectral radius

Directional change in velocity vectors
Shear stress

Gradient

Second-derivatives

Subscripts and Superscripts

ab
AR
AVD
DR
DS
DTR
EAS
Gl
GR

Absolute error

Aspect ratio

Absolute velocity difference
Directional resolution
Directional smoothness
Deformation tensor resolution
Equi-Angle skew
Grid-independent solution
Gradient resolution

Cell index

xvii

re

RVD

SS
VR

Projection of velocity vector
Relative error

Relative velocity difference
Direction of streamline
Streamline smoothness
Volume Ratio

Direction of the eigenvector

xviii

Chapter 1

Introduction

1.1 Background

Computational Fluid Dynamics, popularly referred to as CFD, is concerned with
using computers to obtain numerical solutions to mathematical equations that
describe fluid flow, heat transfer, and combustion problems. CFD started
developing in the late 1960’s and early 1970’s based on the evolution of
aerodynamics theory, applied mathematics, and computer techniques. With the
development of high-speed computers and the improvements in numerical
algorithms, there has been an extraordinary reduction in numerical computation
cost. Thus, more and more research and industrial CFD services were created
in the 1980’s and got significantly expanded in the 1990’s. The wide variety of
CFD applications includes aerospace, automotive, materials processing,

electronics appliances, bioengineering, and biochemistry.

Despite the wide use of CFD nowadays, one of the major concern is what can be
trusted in the predicted CFD results. The accuracy of a CFD solution is
determined by the difference between the flow problem being studied and
numerical solutions. Basically, the error in CFD originates from three main

sources (see Fig. 1).

The first main source of error arises from limitations in the models used to
describe unsolved physics. This error can be defined as the difference between
the real physics and the exact solution of the physical model. For any thermal-
fluids phenomenon, the model used should describe all main physics relevant to
the problem of interest. Otherwise, there are mistakes in the model. Roache
[1998] summarized the reasons for errors in model as being due to a lack of
understanding of the physical phenomenon and uncertainties of parameters used
in the model. For single-phase flow, the physics with the greatest concerns in
modeling is turbulence. Thus, in order to use CFD to simulate a turbulent flow
problem, a good understanding is needed on the limitations and capabilities of

the available models (Wilcox [1993]).

The second main source of error in CFD solutions can be attributed to the
numerical schemes. When numerical schemes are used to replace the
governing partial differential equations by a set of algebraic expressions in a
discrete domain of space and time, many errors are introduced. Advanced high-
resolution schemes endeavor to represent physics such as convection and wave
propagation as accurately as possible with minimal numerical diffusion,

dispersion, and other spurious modes.

The third source of error in CFD solutions is introduced by grid. The grid, also
referred as to the mesh, must represent the geometry and enable the physical

model and the numerical scheme to resolve the relevant flow physics to the

required accuracy. This requirement on the mesh can be difficult to achieve for
three-dimensional (30) problems with complicated geometry and flow features.
One reason for this difficulty is that CFD is constantly being challenged to
address the most complicated problem possible based on the latest computing
capability. As a result, the number of grid points or cells that can be used in a
mesh is often restricted by either the available computer memory or a need to
have a practical tum-around time in computing a solution. With a limit on the
number of grid points or cells that can be used, accuracy considerations demand
grid points or cells to be placed in regions where they are most needed to resolve
the geometry and flow physics and sparsely distributing them elsewhere.
Unfortunately, this non-uniform distribution can create what are referred to as
poor quality cells, which can induce considerable errors in the computed
solutions. Sometimes, poor quality grids can also induce numerical instability or

divergence.

 

Numerical
Schemes

 
 
 

 
 
 

Physical
Models

   
       

Figure 1. Accuracy of CFD results

As summarized in Fig. 1, the physical model, the numerical schemes, and the
grid/mesh are the three critical components that affect the accuracy of CFD
results. For most practitioners or users of CFD, it is the grid that is the most
problematic. There are several reasons for this. First, for many engineering
problems, it is impractical to generate grid-independent solutions. Second, even
when it is essential to generate grid-independent solutions, it is difficult to know
where the grid quality needs to be improved and where higher resolutions are
needed. This is especially true for unstructured meshes, where it is impossible to
visually inspect the mesh. Thus, methods are needed to evaluate grid quality.
Also, methods are needed to evaluate the confidence levels that can be attached
to conclusions derived from CFD results obtained on imperfect meshes. This

dissertation addresses these two issues.

1.2 therature Survey

1.2.1 Grid Quality Measures

In an effort to obtain better understanding of the relationship between mesh and
solution accuracy, many indicators, also referred to as grid-quality measures,
have been developed to evaluate meshes. These measures can be classified
into two broad groups - those that only involve information about the mesh
geometry and those that involve information about the mesh and the solution

generated on that mesh.

Most grid-quality measures developed so far only account for the geometry of the
cells in a mesh. In early work, many measures were developed to assess the
shapes of two-dimensional element. The aspect ratio and angular measure were
widely used to check how the quality influences the mathematical error bounds
(Mitchell et al. [1971], Wilson et al. [1990], Bhatia and Lawrence [1990]). The
measures of triangles were then applied to quadrilaterals (Robinson [1987,
1988], Yuan et al. [1994]). Robinson showed that the quadrilaterals present
much more complex geometric features than triangles and more complicated
definitions of measures should be used. He summarized a set of shape or
distortion parameters for the four-node quadrilateral. They are aspect ratio,
warpage, skew angle, and tapers. Among these measures, aspect ratio is a
shape parameter and can be described as the ratio of maximum side to minimum

side of a planar quadrilateral mesh. The warpage of a quadrilateral is measured

by its deviation from a flat plane. Skew angle and tapers are also distortion

parameters and measured by its deviation from square.

For the extension to three-dimensional measures, various grid-quality measures
have been extensively studied on tetrahedral elements (Liu and Joe [1994], Lo
[1997]). The series of geometrical measures developed for tetrahedral elements
included shape ratio and regularity of tetrahedron (Babuska and Aziz [1976],
Dannelongue and Tanguy [1991]). For 3D hexahedral elements, Robinson
[1994] developed to a list of definitions for three dimensional meshes from his
two dimensional measures, including aspect ratios, skew ratios, tapers and
warpages. Recent studies were also reported on geometric measures for 3D
hexahedral, pyramid, and wedge shape elements (Owen [1999], Wa and Chen
[2000]). Warsi and Thomson [1990] collected a set of weight measures used for
the variational principle approaches such as smoothness, orthogonality, and

clustering in some manner.

Field [2000] reviewed these geometric measures and gave a definition for a good
measure. He presented that a good measure should possess following four
attributes:
0 An ability to detect different types of elements.
0 Non-dimensionality.
This can make triangle/tetrahedral and quadrilateral/hexahedral have

similar dimensionless values.

o Boundedness.
This can limit the measure with arbitrarily large values.
0 A normalization.
This will let the value of measure take positive value between zero and

one and make it convenient and general for comparison.

Compared to the research on geometric measures, relative fewer investigators
studied measures that can account for both geometry and solution. The early
stage of this study was the work of Zlamal [1968] and Babuska and Aziz [1976],
who studied the minimum/maximum angle in different flow conditions. In the
analysis, they presented the requirements for angle in the triangular element
based on the approximate solution, which was in the form of piecewise linear
functions, piecewise quadratics, etc. in finite element methods. The study of
angle condition was further extended to three-dimensional tetrahedral element
(Krizek [1992]). In the work of highly directional flow problem, researchers
realized that it is very important to consider both the local element shape and the
local solution behavior to evaluate the quality of the unstructured meshes
(Peraire et al. [1976], Simpson [1994]). Though a good triangulation should have
as few anisotropic triangles as possible, Rippa [1992] set up a connection
between the shape of the triangle and the behavior of the approximated solution.
He showed that the highly distorted element aligned with the flow direction might

be very effective for the strongly anisotropic solution situations. Brackbill [1993]

and Knupp [1995] extended the formulation to include adaptation to align the grid

with a vector field.

Sorokin [1994] used a grid-quality measure, which was based on cell shape and
size combined with gradient, to evaluate the quality of both structured and
unstructured grids in 2D. Lepape and Jacquotte [1994] used the cell deformation
relative to rectangular parallelepipeds combined with the local gradient to
evaluate the grid quality in block interfaces and boundary surfaces in structured
grids. This approach was also extended to unstructured grids in their study.
Berzins [1998, 1999, and 2000] developed a grid-quality measure based on both
geometry and solution information for triangular and tetrahedral meshes. His
work showed that an indicator, which accounted for both geometry and solution,
might perform better than the geometric only measures when evaluating the

mesh quality.

All efforts reviewed above provide encouraging evidence of the usefulness of
grid-quality measures for evaluate 2D or 3D structured and unstructured meshes,
especially for measures that can account for both geometry and solution
information. However, a serious deficiency of existing grid-quality measures is
that they do not directly correlate solutions errors with grid quality. Also, these
grid-quality measures are not proved to be problem-independent so that certain
magnitudes of measure value can imply certain errors or confidence levels, at

least for classes of flows.

1.2.2 Error estimation

Research on error indicators that account for the solution has been reported by a
number of investigators. For the general elliptic systems of partial differential
equations, once the approximate solution is obtained on the initial mesh, this
result is the only available data and can be constructed as a posteriori error
estimates to evaluate the accuracy of the initial solution. For earlier work on a
posteriori error estimators on isotropic meshes we refer to Babuska and
Rheinboldt [1978]. Babuska and Rheinboldt presented an approach that is
similar to the well-known residual method. The principle difference is that they
used the norm of negative Sobolev spaces and the error estimates are given in
an asymptotic form for h—>0 where h is the size of the elements. For further
work, of special interest are the indicators of Bank and Weiser [1985] and the

residual based indicators of Verturth [1994].

Within the framework of finite element methods, numerous theories and
approaches were developed for a posteriori error estimation. These strategies
can be classified into following different techniques (Oden et al. [1989], Verfurth
[1994]).
0 Residual based methods

This error estimate can give a strict upper bound on the solution error in the

energy norm. This approach is very popular in recent error estimation in the

finite element method and produces impressive results for a large variety of
problems (Ainsworth and Oden [1992], Strouboulis and Hague [1992]).
Duality methods

This approach was developed based on the Synge’s ‘hypercircle’ in
mathematics physics and Ekeland's study on convex analysis and variational
problems (Synge [1957], Ekeland et al. [1976]). In this method, the duality
theory of convex optimization, which is valid for self-adjoint elliptic problems,
is used to compute the error bounds.

Subdomain-residual methods

The local error in an element can be derived from a patch of elements
surrounding the given element. This approach was presented by Babuska
and Rheinboldt [1978] for bilinear finite element approximations.

Interpolation methods

In this method, the interpolation theory of finite element is used to generate
the estimation of local error for each element based on the leading term of
truncation error.

Using solution of local problem

By using this method, the original problem was transformed to a similar, but
simpler locally discrete problem. Then, the approximate solution of original
problem is computed on the locally simplified problem and used to make error
estimation (Bank and Weiser [1985], Oden et al. [1989], Bank and Welfert
[1990]).

Sharp a priori error estimates method

10

The higher-order difference schemes are used based on the computed
numerical results to improve the accuracy in predicting the derivatives in the
priori error estimators. The new estimators can be called “sharp a priori error

estimator” (Eriksson [1985], Eriksson and Johnson [1988]).

Above results and methods relate mainly to the linear elliptic problems. In the
area of nonlinear problems or problems of other types such as parabolic and
hyperbolic, the algorithms of error estimation are still not well understood. Some
works on the error indicator of finite element method for parabolic equation have
done by Davis and Flaherty [1982], Eriksson and Johnson [1991], and some
recent development was described by Adjerid et al. [1999]. New constructed a
posteriori error analysis of hyperbolic problems and applications of fluid

mechanics can be found in Johnson et al. [1995], Suli [1998], Barth et al. [1999].

The error estimation for the finite difference method and finite volume method is
still in its early development stage compared with the numerous developed error
estimators for finite element method. The study of error estimation for the finite
volume method began with the study of numerical diffusion induced by the
turbulence modeling. Large numerical diffusion can be generated when the
upwind scheme is applied for turbulence model and interfered with the turbulence
diffusion. This numerical diffusion was estimated based on upwind and central

difference scheme (McGuirk and Rodi [1978]).

ll

Taylor series expansions are often used to estimate errors in finite difference
/finite volume (FD/FV) solutions of boundary-value and initial-value problems.
The leading error term in Taylor series can be used as an error indicator. The
other method is based on polynomial approximation (Weierstrass theorem),
which states that one can always find a polynomial to fit any data to any accuracy
desired. Its polynomial coefficients are computed by using Taylor series
expansion. When estimating the bound of high order derivatives or the
polynomial coefficients without using the approximate solution, the error estimate
can be referred as a priori error estimate. As the grid spacing is refined, these a
priori error estimates can show the asymptotic trend of the FD/FV solution
approaching the exact solution. In other words, the a priori error estimates can

only provide qualitative information on error behavior. (Carey [1993])

In order to quantitatively estimate the error for a given problem, a posteriori error
is generally needed. When high order derivatives in Taylor series method used
in a priori estimate can be bounded by a computable quantity based on the
FD/FV solution, the a priori estimate then becomes a posteriori error estimators.
Richardson extrapolation is so far the most popular error estimator used for the
FD/FV methods (Richardson, [1927]). It is a method for obtaining a higher-order
estimation from a series of lower-order discrete values. First, two sets of grids
with two levels of refinement should be generated for a given problem. Then, a
higher-order estimate can be computed from the two level composite solutions by

using extrapolation method. Finally, the error level in the coarse mesh can be

12

estimated by the difference of the lower-order solution on coarse mesh and the
extrapolated higher-order solution. With more levels of refinement, this process
can continue until it reaches the desired accuracy level. If the grids are fine
enough, Richardson extrapolation can be very reliable and is the only method
that can handle the non-linearities in problem (Muzafeijia and Gosman [1997]). It
has been successfully used in extensive problems, including supersonic flows

and incompressible flows.

Haworth et al. [1993] proposed a local error estimation approach based on the
cell-to-cell imbalances in kinetic energy and momentum and used it on a
transient flow problem in an axis-symmetric IC engine flow. Jasak and Gosman
[1998, 1999] presented a residual error estimate for the finite volume method.
This approach borrowed the idea from residual method for finite element method,
in which the error was estimated from the inconsistency between control volume
integration and the face interpolation value. It can also be extended to include

the temporal error.

Though extensive studies have been performed in the error estimation for the
finite volume method, there are many limitations and weaknesses of those
developed approaches. For methods of error estimation that are based on the
leading term of the truncation error (TE) from a Taylor series expansion, they can
reveal information only when the grid spacing or cell size is sufficiently small so

that TE is bounded. The polynomial approximation (Weierstrass theorem) states

13

that one can always find a polynomial to fit any data to any accuracy desired.
But, Weierstrass theorem breaks down when there are round-off errors. For
example, the slightest round-off error in the coefficients of a high-order
polynomial would produce spurious oscillations. Richardson extrapolation, which
is considered as the most successful approach for error estimation in finite
volume method, requires not only the high resolution in mesh but also the
generation of successive meshes, which is not feasible for many industrial
applications. For other approaches based on cell imbalance, residual, etc. are

still in its early age of development.

Also, the main disadvantage of these existing error indicators for use is that they
do not directly link solutions errors to grid quality. They only show where error at
one location may be higher than another region, which is maybe adequate for
mesh refinement but not good enough for estimating grid-induced errors in

practical problems.

14

1.2.3 Grid Adaptation

The accuracy and efficiency of finite difference and finite volume methods used
in CFD are strongly dependent upon the mesh. On accuracy, the mesh must
well represent the geometry and have sufficient good quality in complex flow
regions so that all the relevant flow physics can be resolved and minimum grid-
related errors will be induced. On efficiency, the number of grid points or cells
used should be kept to the minimum needed to achieve a prescribed accuracy
because of the limitation on computational resource and time. The conflicting
requirements of accuracy and efficiency on the mesh imply a need for

optimization.

The goal of solution-adaptive mesh refinement is to produce an optimum mesh
for a given problem with a prescribed accuracy requirement. For any mesh
adaptation process, it is made up of two key elements. The first one is an error
estimator, which is used to provide information about the regions where the grid
should be refined/coarsened. The development of different error estimators has
been given in previous section. The second one is an adaptation procedure,
which is a mechanism to alter the mesh in an appropriate manner. This
procedure can be accomplished in various ways. The h-refinement, p-

refinement, and r-refinement are the three basic categories.

15

The technique h-refinement is referred to as mesh subdivision when used with
unstructured grids (Coelho et al. [1991], Vilsmeier and Hanel [1993]). It is easy
to change the spacing of the computational points by adding points in regions
where refinement is required or removing points where the region should be
coarsened. Extensive studies have been made in this area. New refinement
nodes may be inserted into the baseline mesh by using the Delaunay point
insertion algorithm (Bowyer [1981], Mavriplis [1992]), or utilizing an element sub-
division with local reconnection procedure (Marcum [1995]), or the hierarchical-
based element sub-division approaches (Lohner and Baum [1992], Connell and
Holmes [1994]). A local edge-face swapping techniques sometimes was
followed after the primary refinement in order to optimize the grid. Good
examples also includes that where tetrahedra were divided into eight cells with
special care to remove hanging points (Kallinderis and Vijayan [1993]), where
nodes were inserted at the mid-points of cell edges (Webster et al. [1994]), and

which introduced three types of subdivisions (Biswas and Strawn [1993]).

There are also a number of studies on h-refinement applied to structured meshes
(Schonfeld [1994]). If only hexahedral elements were used throughout the flow
domain, h-refinement produces hanging nodes and they need special
consideration. For cell-centered flow solver, this difficulty can be easily handled
by using a data structure based on the faces of cells (Spragle and Smith [1997]).
For flow solver using vertex-based schemes, these hanging nodes must be

handled differently from the regular nodes (Aftosmis [1994]).

16

Some researchers create a hybrid grid by using triangular elements or using
pyramidal and prismatic cells to treat the hanging nodes, which is called hybrid
structured-unstructured grids or mixed element unstructured grids (Szmelter et
al. [1991], Aftosmis [1994]). Such adaptive grid can reduce complexity and
maybe improve the solution accuracy compared with using fully tetrahedral
meshes. On the other hand, it can also reduce the complexity of handling
hanging nodes and modifying codes compared with using fully hexahedral

meshes.

Other work using h-refinement grid adaptation includes a grid utilizing a set of
overlapping patches with increasing resolution (Caruso [1985]). This approach
was used for transonic and supersonic flows with discontinuities. Each patch
used a structured and orthogonal mesh and was solved by flow solver
individually. The data transfer between patch and the background mesh were
done through the patch boundary conditions. Berger and Oliger [1984]
presented an adaptive technique for 1D and 20 hyperbolic equations. The grid
refinement was based upon Richardson-type estimation of the truncation error
and done locally in both time and space. In 2D grid, the subgrids were
automatically handled as overlapping rectangles in arbitrary directions. Based on
this starting point, Berger and Colella [1989] further improved the automatic
adaptive mesh refinement strategy for solving 2D hyperbolic conservation-laws.

This strategy was used for solving unsteady flows with shocks and the refined

l7

subgrids can be conveniently rotated according to the baseline grid. Thus, it
would allow the mesh geometry to be approximately aligned with the flow
patterns. Chen et al. [1997] presented a multi-level flow adaptive mesh-
refinement strategy used with structured, non-orthogonal finite-volume scheme.
This approach can be applied to laminar and turbulent flows. Based on a coarse
baseline grid, the regions needing finer grids were identified by some parameters
derived from the preliminary solution. This work got much improvement over
earlier work, but it is not suitable for situations requiring large number of
refinement level. Otherwise the number of block is too large to be handled

properly by code.

Mesh movement, or r-refinement, is also an alternative of mesh adaptation
method (Hawken et al. [1991], Dwyer [1984], Dandekar et al. [1993]). In this
method, the mesh topology and number of degrees of freedom is fixed while the
mesh is distorted to reduce the local errors. The advantage of this approach is
that the number of mesh cells does not change and hence the computational
cost does not increase when a new adapted mesh is applied. However, for r-
refinement, there is no guarantee of generating a sufficiently accurate solution
and the procedure of redistributing the grid nodes can easily create highly

distorted cells.

The third adaptation method, the p-refinement, does not involve the mesh.

Indeed, it alters the order of accuracy of the discretization scheme (Zienkiewicz

l8

[1989]). It is particularly suitable for finite element computation and gets
considerable success in structural mechanics and many fluid flow problems.
However, the higher-order formulae are not likely to be much more accurate than
lower-order ones for situations in which flow discontinuities can arise. So p-

refinement should be applied to solutions with “sufficiently smooth”.

Besides above three broad categories of adaptive schemes, a composition of
these methods combines the advantages of h-refinement and p-refinement, and
leads to a hybrid hp scheme (Babuska and Dorr [1981], Oden et al. [1989]). This
approach can generate significantly accurate results at the cost of increasing
program complexity. Other combinations are also in an early stage of
development, like a p-refinement associated with r-refinement strategy or
redistribution perhaps following h-refinement in unstructured grids as a

smoothing operation (Richter [1994]).

Thus, considerable research has been conducted in grid adaptation. Despite the
progress made, there are two concerns. The first is uncertainty about the
accuracy of the error indicators. The second is uncertainty about the most

efficient way to refine the mesh.

19

1.3 Objectives and Outline

In the previous review, a literature survey was given on grid-quality measures,
error estimation, and grid adaptation. In that survey, a number of research needs
were identified. The objective of this dissertation is to address these needs.

More specifically, the objectives are as follows:

1. Develop grid-quality measures that extract and link information about the
mesh and the scalar, vector, and tensor nature of computed CFD solutions.

2. Evaluate and examine how well these measures correlate with grid-induced
errors.

3. Examine how well these measures can be used as error indicators and for

grid refinement.

This dissertation is composed of two parts. Part I, consisting of Chapter 2 to 5,
presents the development and evaluation of grid-quality measures. Chapter 2
summarizes the characteristics of ideal measures and describes a set of newly
developed grid-quality measures. In Chapter 3, the approach used to evaluate
these measures is introduced. Chapter 4 describes the test problems used to
evaluate grid-quality measures and all solutions generated. Chapter 5 addresses
the details of how measures are evaluated and also presents the results on how

well these newly developed measures correlate with errors.

20

Part II, consisting of Chapter 6 & 7, describes the application of these measures
as error indictors and for grid refinement. In Chapter 6, these measures are
applied to two examples and examined how well they can indicate the real error.
Chapter 7 presents the application of measures developed for grid refinement.

The summary and conclusion of this dissertation are presented in Chapter 8.

21

PART I

Development and Evaluation of Grid-Quality Measures

22

Chapter 2

Development of Grid-quality Measures

In this chapter, we present the grid-quality measures developed in this study.
These measures account for the geometry of the mesh and the solution
information. They were developed after extensive numerical experiments in
which numerous measures were developed and evaluated. The approach used

to evaluate these measures is given in chapter 3.

The organization of this chapter is as follows. First, the conventional grid-quality
measures are introduced. Then, we summarize the basic classification of
measures and give the attributes of ideal measures. Finally, we present the

definitions of six newly developed grid-quality measures.

2.1 Conventional Grid-Quality Measures

The conventional approach to check grid quality typically involves three steps.
The first is called “fundamental check”. In this step, those cells with negative
volumes, overlaps, left-handed faces, and topological inconsistencies are
detected and fixed. In the second step, the check focuses on the local grid
quality. The main qualities considered include cell size, smoothness, and cell
shape. The third and last step is to view the grid in details and further confirm its

quality through visual inspection. By setting the threshold of displaying

23

measures, the cells with bad measure values will be indicated. These
highlighted cells can then be fixed or just accepted depending on the local flow

features and user’s experience.

As noted, conventional grid-quality measures mainly focus on cell size,
smoothness, and cell shape. On cell size, since the continuous domain is
discretely defined, the grid spacing in complex flows (e.g., flows with shear
layers, separated regions, and shock waves) must be small enough to resolve all
of the relevant flow physics. In many flow situations, inadequate resolution in
critical regions can significantly alter the local as well as the global flow
characteristics. Smoothness normally means the changes in cell volume
between neighboring cells. Smoothness with bad value can cause large
truncation errors because cancellation of errors will not take place. The cell
shape (including cell aspect ratio and cell skewness) also has effect on solution
accuracy. Aspect ratio is a measure describing the stretching feature of the cell.
For highly sheared flows, high aspect ratios may achieve the desired accuracy
with fewer cells if grid lines can be aligned with the flow. For recirculating flows,
aspect ratio should be near unity. As for skewness, it can be expressed as the
difference between the actual cell shape and the shape of a perfect cell, namely
an equilateral cell, in terms of angle of intersection. For example, optimal
quadrilateral and triangular meshes should have all vertex angles close to 90 and
60 degrees, respectively. Highly skewed cells can decrease accuracy and

destabilize the solution.

24

Some of the well-known and widely used conventional grid-quality measures are

summarized below.

. Cell resolution
2-D: cell area
3-D: cell volume
0 Smoothness
It is used to describe the changes in grid spacing and defined as the
maximum area or volume ratio of the current cell and all its neighbors.

2-D: QSM = max(Ai lAj),j = 1,2,...,N (2.1a)
3-D: QVR = max(Vi /Vj),j = 1,2,...,N (2.1b)

where, A., Vi is the area and volume of current cell, respectively. A], V,- is the
area and volume of its neighbor, respectively. N is the number of all its
neighbors.

o Aspect ratio

It is used to describe the stretching feature of the cell and defined as:
Q AR = )6 (5), triangularel ements
I‘

= %[§), tetrahedral ements
r

_ max[e1,e2] (2'2)

— , , quadrilateral elements
mm[e,,e2]

max e ,e ,e
= [ ‘ 2 3], hexahedral elements

 

min[e1,e2,e3]

25

where, r, R represent the inscribe and circumscribe radii of circles or spheres.
e; is the edge length of cell. OAR = 1, when describing an equilateral element.

Cell equiAngle skew

It is a dimensionless parameter calculated using the normalized angle

deviation method and is defined as

 

QEAS = max{2’-'—m — 9O 90 — 9m—}, for structured mesh

 

180 — 9o ’ 9o

9 9 (2.3)
Q as = max{—"—“" _ 6O , 6O _ "“4— , for unstructured mesh

180 — 6O 60

where, 8mm, and 8min represent the maximum and minimum angles between
the edges of the element. QEAS = 0, when describing an equilateral element.
QEAs = 1, when describing a degenerate element.

Cell warpage

It is the square root of the ratio of the distance between the cell centroid and

cell circumcenter and the circumcenter radius:

 

|'r‘(centroid) - 'r’(circumcenter)

 

(2.4)

 

w a e=
arp g R(circumcenter)

26

2.2 Characteristlcs of Ideal Grld-Ouality Measures

2.2.1 Basic Classification of Measures

Before developing new grid-quality measures, it is important to obtain better
understanding of measures that can relate to grid quality. Basically, the grid-
quality measures can be classified into two broad groups -- those that only
involve information about the mesh geometry and those that involve information
about both the mesh geometry and the solution obtained on that mesh. Most
grid-quality measures developed so far only focus on the geometry of the cells in
a mesh (9.9., the conventional measures stimmarized in the previous section).
Relatively few investigators studied measures that account for both geometry
and solution. These two groups can further be classified as local (involving only
one cell), macro (involving the cell and its neighboring cells), and physical
(involving cell that related according to physics such as cells along the 'same
streamline). Obviously, physical measures require information on both geometry

and solution. This classification of measures is summarized as below:

I Mesh resolution
I Local i
Geometry i [Cell shape
Only I Macro: Change of measure

I Local: Resolution only related to local solution behavior

Geometry 1 Macro: Combination of local and neighboring solution
&Solution I
I Physical: Combination of local and physics-related solution

I
Grid I
I
I

27

Thus, the conventional measures -- cell volume, aspect ratio, and skewness --
can be considered as local and geometry only measures, while the smoothness

is a macro and geometry only measure.

2.2.2 Attributes of Ideal Measures

Before developing and evaluating grid-quality measures, we summarize some

features that we think ideal measures should have. This is not an exhaustive list,

but rather the minimum requirements.

0 Account for the geometry of the cells in a mesh and the solution obtained on
that mesh;

. The definition should be universally applied for all 2D/3D, structured
lunstructured meshes;

o The formulas should be non-dimensional;

0 Have certain correlation with the solution error;

The goal in following section is to present a set of grid-quality measures that

possesses these characteristics.

28

2.3 Newly Developed Grid-Quality Measures: Schematic and Formula

In this section, we present six grid-quality measures that were developed under
this study. These measures were developed through extensive numerical
experiments. As will be shown in Chapter 5, all of these measures correlate with

solution error.

2.3.1 Directional Resolution

The first measure that gauges resolution is based on the directional change in
the velocity vector. This measure, referred to as directional resolution and
denoted as (lion, can be defined as the maximum change in flow direction at a cell

relative to all its neighboring cells. Its mathematical expression (see Fig. 2) is
on, = max{9i }, i= 1,2,...,M (2.1a)

where e, is the change in flow direction at a cell relative to those at neighboring

cells and M is the number of all the neighboring cells. For 2-D problems, M = 3

for triangular mesh and M = 4 for quadrilateral mesh with 8. given by

 

[ uou.+v0vh 2]
6i=COS
f§+v§ *fJu +v,2 (2.1b)

VO = not: + voj (2.10)
V, =“t'7+"ij (2.1d)

For 3-D problems, M = 4 for tetra mesh and M = 6 for hex mesh with 8; given by

29

 

6. = cos" uou, + vov, + wow,
' Ju2+v2+w2*Ju2+v2+w2
0 0 0 i I i

(219)
I70 = uoi’ + v0] + wait (2”)
lZ=uf+vj+wj€ (219)

This measure accounts for both the geometry feature and the local solution
behaviors. So, it is a macro and solution-based grid-quality measure. For a
recirculating flow region of fixed size, it is clear that by decreasing grid spacing,
the maximum change in flow direction from one cell to the next will be reduced.
Thus, one measures the resolution of the directional changes in a vector field and
has clear physical meaning. Additionally, its definition can be applied to
structured and unstructured cells in 20 and 30. Also, this measure is

dimensionless. So, one meets all properties of an ideal measure except for its

correlation with solution error, which we will describe in Chapter 5.

 

 

 

(a) 2-D Triangle mesh (b) 2-D Quadrilateral mesh

30

 

(c) 3-D Tetrahedral mesh

Figure 2. Schematic for directional resolution one

2.3.2 Gradient Resolution

The second measure that gauges resolution is based on a length scale that
represents the curvature of the velocity field and an effective cell dimension in
the direction of change. This measure, referred to as gradient resolution and

denoted as one, is mathematically defined as follows (see Fig. 3):

Lo

 

¢GR = L ii (2.28)
W
where,
LGR=max max{|l. e ,|},i= ..,,Mj=1,2... ,N (2.2s)

 

     

lie length of edge i , d, a length of diagonal j;

31

For 20 problems, M = 3, N = 0 for triangle mesh and M = 4, N = 2 for

quadrilateral mesh.

 

 

Vu+Vv

LW_ V2u+V2v ’
au— au~

Vu=—i+—' 2.2c
8x 3y ( )
av~ av-

Vv=——i+—°
8x By

For 3-D problems, M = 6, N = 0 for tetra mesh and M = 12, N = 16 for hex mesh.

|Vu+Vv+VwI
‘72IV2u+V2v+VwI’
Vu= 3“ z+a—-uj+a—u§
3" By 32 (2.2d)
av~ av. av-

v =—- +—k
V ax'+ay az

Bw- 6W1 8w~
V — — k
w: axi+ayj+37

 

In the above equations, the effective grid spacing Lee is defined as the largest
distance spanned by a cell in a direction perpendicular to the velocity since most
velocity gradients are typically largest in that direction (e.g., boundary-layer
flows). M is the total number of edges that connect the nodes of a cell (denoted
as l.), and N is the total number of diagonals (denoted as di); see Fig. 3.

Diagonals are relevant only to elements with five or more faces.

32

For this measure, the length scale of the curvature in velocity components is
represented by the ratio of first- and second- derivatives with the largest
derivatives dominating. Since resolving a solution depends on the magnitude of

the spatial variation in its gradients, the effective grid spacing is normalized by
LW, given by Eqs. (2.2c) and (2.2d). With this normalization, smallest grid
spacing would typically be needed in regions where the second-derivatives are
the highest or where the first derivatives are changing the most. Note that if the

second derivatives are all zero, then LW = co, and lies: 0, indicating that grid

spacing can be very large. oGR can be considered as a local and solution-based

grid-quality measure because only the information on the given cell is included.
As noted, its expression can readily be generalized to other scalar variables such
as pressure and temperature by defining the effective grid spacing to be in the

direction of the gradient of that scalar field.

 

 

 

 

 

   

(a) 2-D Triangle mesh (b) 2-D Quadrilateral mesh

33

 

 

 

 

 

 

(0) 3-D Tetra mesh (d) 3-D Hex mesh

Figure 3. Schematic for gradient resolution one

2.3.3 Deformation Tensor Resolution

The third measure that gauges resolution, referred to as deformation tensor
resolution and denoted as onTe, is based on the deformation tensor,
A = T/u (2.3a)

where T is the stress tensor for a Newtonian fluid given by

T“ 1,, In
T = I” T” In (2.3b)
In In T
1 = 1 = a_u. + 2V—
xy ’1 8y ax

(2.30)

_ _ 3v 3w
Tn—sz- 5:4"?

With Eq. (2.3), this measure can be defined as follows (see Fig. 4):

34

 

= LG“ (2.4a)

In this equation, Lees is the effective grid spacing in the direction of the
eigenvector, 5, corresponding to the spectral radius of the deformation tensor
given by Eq. (2.3). Lx, the length scale used to normalize Lees, is magnitude of

the local velocity divided by the spectral radius; i.e.,
1., = J..2 + v2 + W2 / e (A) (2.4b)
where u (A) is the spectral radius. For 20 meshes, the w-component in

equations needs to be ignored.

The physical meaning of this measure is that the spectral radius represents the
principle strain or the largest deformation, and the eigenvector is the direction of

that deformation. This measure is also a local and solution-based grid-quality

measure.

 

 

 

(a) 2-D Triangle mesh (b) 2-D Quadrilateral mesh
Figure 4. Schematic for deformation tensor resolution onTe (only 2D meshes is

shown, SD meshes is similar)

35

2.3.4 Directional Smoothness

The fourth measure gauges smoothness of a mesh along the flow direction; i.e.,
how effective grid spacing changes along flow. This measure, referred to

directional smoothness and denoted as ons, is defined as follows (see Fig. 5):

. L . . 1
(905 = max{r,. },z=1,2,...,M I; =—5—‘-, lf 1;. < 1, then r; =— (2.5)
L... r.-
where M is the number of all the neighboring cells. For 2-D problems, M = 3 for
triangular mesh and M = 4 for quadrilateral mesh. For 3-D problems, M = 4 for

tetra mesh and M = 6 for hex mesh.

As noted, this definition of smoothness considers changes in effective grid
spacing along the direction of the velocity. It differs from conventional definitions
of smoothness, which are based on ratios of grid spacing along certain
coordinate directions for structured meshes and ratios of contiguous cell volumes
for unstructured meshes. This measure can be classified as a macro and

solution-based grid-quality measure.

36

  

 

:5V4

 

(a) 2-D Triangle mesh (b) 2-D Quadrilateral mesh
Figure 5. Schematic for directional smoothness ons (only 20 meshes is shown,

30 meshes is similar)

2.3.5 Streamline Smoothness

The fifth measure gauges the relative changes in effective velocity along the
streamline. This measure, referred to as streamline smoothness and denoted as

¢ss. is mathematically expressed as follows (See Fig. 6):

(055 = max{r, },i= 1,2,...,M, r. = ”Vi'Pl-IVOII (2.6)

lVol
where 17”. is called effective velocity, which is the projection of velocity vector 17,

on I70 direction. The M is the total number of neighboring cells along the

streamline. The streamline cell is determined by the procedure below.

For 2D problems,

37

V0 = 1407+ voj (2.7a)

 

 

 

 

6" :COS-l(u0(xi _x(:)+30(yf —y0)) (2.7b)
lees-"Va ~
at; = (xi _ 10);..1'0’.’ _ yo); (2.7C)
For 3D problems,
I70 = uof. + v0] + wol-c. (2.7d)
6'. = cos'1(u°(x" _ x0) + v09“ :yo) + W0(Zi — 10)) (2.76)
lecx IVOI
e.- = (x.- - x0)? + (y. — yo)? + (z. — z, )1? (2.7:)

where, Cn (xn, yn) is the cell center of the given cell. C1 (x1, y1), CIn (xn, ye) are

the cell center of the neighboring cells (see Fig. 6). EU is a vector pointing from
the center of a neighboring cell i (i=1,,..,2 n) to the center of the given cell. 170 is

the velocity vector on cell 00. 8i E (0,180) is the angle between ZCJ. and I70.

In terms of the value of 8., all the neighboring cells can be divided into two
groups. One group of cells is called “streamline cell”, which locates on the same
streamline as the given cell. The other group is called “perpendicular streamline
cell”, which locates perpendicular to the streamline. This can be described as
below (See Fig. 6).

058.<45°, indicate streamline cell, see Ca in (a) and (b);

45°58i<135°, indicate perpendicular streamline cell, see c1 in (a) and 01, on In (D);

135°58i<180°, indicate streamline cell, see 03 in (a) and c4 in (b);

38

 

(a) 2-D Triangle mesh

 

 

 

 

(b) 2-D Quadrilateral mesh

Figure 6. Schematic for streamline smoothness one (only 20 meshes is shown,

SD meshes is similar)

¢ss describes the change of effective velocity magnitude between cells on the
same streamline. So, it is a physical and solution-based grid-quality measure
and can be considered as an improvement of ons. Fig. 7 presents the
advantages of files over one in certain flow situations. For the flow situation of Fig.
7(a), the velocity distribution is uniform but the grid spacing is changing along
flow direction. L51, L52, L53 are the effective grid spacing on each cells. Since
the effective grid spacing changes significantly, ons will have very bad value

according to its definition. While measure ¢ss will indicate perfect value here

39

since velocity keeps unchanged for all cells. According to research experience,
the requirement for grid quality is not so critical for uniform flows. Probably no

high error will be induced even though high smoothness mesh is used. So the

prediction by oss seems more reasonable than that by ons.

In contrast, Fig. 7(b) shows another situation in which the generated mesh is
uniform whereas the velocity is changing. In this kind of flow condition, the local
grid-quality is very important and high error can be easily activated. According to
the definition of one, the effective grid spacing is same for all cells, so ons gets an
ideal value here. As for (has. it can indicate that the local grid quality is not very
good based on the local solution behavior. So oss may be more useful than ons

for this kind of flow.

C6111 C6112 C6113

 

 

 

(a) Uniform velocity on non-uniform mesh

40

Cell] Cell; C6113

 

 

 

(b) Non-uniform velocity on uniform mesh

Figure 7. Comparison of oss and ons

2.3.6 Relative Velocity Difference

The sixth measure gauges the relative changes in effective velocity
perpendicular to the streamline. This measure, referred as velocity difference
and denoted as ¢RVD. is expressed as Eq. (2.8) (See Fig. 6), which is very similar
to the definition of oss. The only difference is that <I>ss describes the relative
changes in effective velocity along the streamline direction, while ¢vo considers
the relative changes in effective velocity perpendicular to the streamline direction.

¢RVD = max{ri 1:1 = 1,2,...,N, r. :11le

. (2.8)
MI

where N is the total number of perpendicular streamline cells. oevn presents bad
values in the regions with high velocity gradients, such as boundary-layer flows

and is also a physical and solution-based grid-quality measure.

41

Chapter 3

Approach to Evaluate Grid-Quality Measures

In our approach to develop and evaluate grid-quality measures, we did not use
Richardson’s extrapolation, Taylor series, or higher-order polynomial
approximations to estimate solution errors. Our approach to develop and

evaluate grid-quality measures involves three major components:

1. Define errors. This task is to define all errors that are introduced by grid
quality in a computed solution. What we did first is to compute a grid-
independent solution. This solution is taken to be accurate baseline or reference
solution. Next, a series of approximate solutions are computed based on
imperfect meshes by using the same flow model and same numerical scheme.
Then, the difference between each approximate solution and the grid-
independent solution can be defined as the error, which induced by the imperfect
grid. The errors examined at an arbitrary cell i are absolute error (Em), relative

error (E,.,,;), and spatial change in relative error (AEm), which defined as follows

(see Fig. 8):
Eabj = Ive” _ Vil (3-131
Em : Ivan — VII/IVGIJI (3-1b)

ABM. = max IBM - E (3.1c)

re.j|

42

In the above equations, I7. is the solution for the velocity vector at the ith cell

obtained on an “imperfect” mesh, and so has errors. ‘75:..- is the grid-independent

solution at the ith cell for the same problem. In Eq. (3.1c), j includes all cells that
share a face with the i'h cell (i.e., all the neighboring cells).
Vi

/

C611] V 61

var Vi

 

 

 

 

 

Figure 8. Schematic for error definition

Though the absolute error is included, it is recognized that this error type is not
meaningful since there is no reference to gauge what is good or bad. The
relative error and the change in relative error, however, are meaningful. The
change in relative error was defined to examine situations in which the grid
quality is good, but the error is high. Such situations occur because the errors
were created elsewhere, where grid qualities are poor. If this is the case, the
change in relative error may reveal information about how a mesh quality is

locally increasing relative error.

2. Specify test problems and generate solutions. Three steady and
incompressible laminar flow problems that have complicated flow features are
selected as the test problems. One problem involves 2-D flow in a square

chamber with two inlet ducts and one exit duct. The second problem describes a

43

3-D flow in a cubic chamber with four curved inlet ducts and one exit duct. The
third one is a 3-D flow in a cubic chamber with two opposing inlet ducts and two
exit ducts. The reason for using rectangular or cubic chamber is that aspect ratio
of unity grids can be generated. Aspect ratio of unity grids are needed because
they have the least grid-induced errors. Based on the excellent grid quality (all
cells are uniform and all have aspect ratios of unity), the grid-independent
solution for each problems can be obtained. Once the grid-independent
solutions for each test problem have been generated, a series of approximate

solutions are generated based on these imperfect meshes.

3. Evaluate grid-quality measures. Develop and interrogate a database to
determine the correlation or functional relationship between grid-quality

measures and solution errors; i.e.,
E = F(¢1’ ¢2v "-t ¢N) (32)

where E is the solution error; oi is the ith grid-quality measure; and F is the

relationship between E and the oi ’s.

The basis of the approach described above for developing and evaluating grid-
quality measures hinges on two premises concerning the correlation in Eq. (3.2).
The first is that the correlation should be monotonic. The reason is that only
monotonic functions have inverses. Having an inverse will enable us to
determine the o values needed to achieve a desired error bound. The second is

that the function F is fairly universal, at least within a class of problems. With

these premises, the goal is to determine the complete set of grid-quality

measures, o,, and its correlation with grid-induced error.

45

Chapter 4

Test Problems and Solutions Generated

As noted in Chapter 3, the grid-quality measures developed in this study are
evaluated by interrogating a database. This database is made up of many
solutions, generated on a variety of meshes with grid quality varying over a wide
range. This database also has the grid-independent solution so that the error at
every cell can be assessed. In this chapter, we described the problems, referred
to as test problems, which are selected to develop the database. We also

describe the solutions generated for each problem.

4.1 Basis of Problem Selection

In this study, the test problem were determined according to the following

considerations:

0 Geometry: Rectangular shape should be used. This is so that it is possible
to generate aspect ratio of unity grids, which have minimal grid-induced
errors.

. Flow features: The flow features should be complex. For example, it should
include shear layers, recirculating flows, and flow separation, etc.

o Formulation: The formulation should be well established. This is so that the

physical meaning of the solutions generated can be interpreted rigorously.

46

Based on these considerations, we formulated three problems, and they are

. 20 flow in a square chamber with two inlet ducts and one exit duct

. SD flow in a cubic chamber with four curved inlet ducts and one exit duct

. 3D flow in a cubic chamber with two opposing inlet ducts and two exit ducts

In the following subsections, we will describe these three test problems and the

solutions generated.

4.2 20 Test Problem: 20 Flow In a Square Chamber with Two Inlet Ducts

and One Exlt Duct

In this section, the description and formulation for this 20 test problem will be first
presented. Also described are the details of how to obtain grid-independent

solution and generate imperfect solutions.
4.2.1 Problem Description and Formulate
A schematic diagram of this test problem is shown in Fig. 9. It is a 20 square
chamber (300mm x 300mm) with two inlet ducts and one exit duct. The inlets

are two straight ducts (200mm x 50mm). The exit duct has the same dimension

as the inlet ducts. This geometry is symmetric about the X-axis.

47

The fluid in the chamber is water with density p = 998.2 Kg/m3 and dynamic
viscosity p. = 1.003E-3 kg/m-s. The inlet velocity is uniform. The flow speed (u =
0.0048m/s) is determined so that the flow remains laminar. The static pressure
at the outlet is fixed at the standard atmosphere pressure Po. At all walls, the no-
slip condition is employed.

No-slip
Wall

   
  
 
  

50

1 50
3L2..__
f Symmetric Plane

:; / Pressure

Velocity Outlet
Inlet
{-— 200 4— 300 4— 200 ——§

Figure 9. Geometry and flow conditions of the 20 test problem

As noted, this flow problem was selected because it possesses complicated flow
features such as flow separation, jets, and recirculation zones. Because of
symmetry in this problem, we consider only one half of the geometry. Also, we
are only interested in the steady-state solution of this incompressible laminar
flow. The equations govern the flow are the continuity and momentum equations
(full Navior—Stokes equation). Since p. st f(T), the energy equation is decoupled
from the continuity and momentum equations and was not used in this study.

The details of governing equations, boundary conditions, numerical schemes,

48

and control parameters can be found in Fluent manual. Solution to the governing

equations was obtained by using the Fluent code.
4.2.2 Grid-Independent Solution Generation

The detailed process used in this study to compute the grid-independent solution

is as follows:

1. Generate an initial grid with unit aspect ratio and then obtain its solution
(referred to as level k grid/solution).

2. Generate a finer uniform grid with unit aspect ratio and obtain its solution
(referred to as level k+1 grid/solution).

3. Project the solution on the level k+1 grid onto the level k grid by using bilinear
interpolation.

4. Compute the average errors between level k and level k+1 solutions by

i[(IVk,i —vk+l.i I) _. ]
— IVk+l,i I

 

 

x __ i=1
E... .. N (4.1)
N -—o «no
_ 2(I ij " Vk+l.i I)
13:, = i=1 N (4.2)

where V“. is the velocity vector at the ith cell obtained on level k mesh. Vb,”

is the level k+1 solution projected on the level k mesh at the ith cell. N is the

total number of nodes in the level k mesh.

49

5. Check the convergence of solution to grid independent. If the E; is less than

the lower bound 8, then, the level k+1 solution is taken to be the grid-

independent solution. Otherwise, repeat steps 2 to 4.

Table 1 summarizes the series of uniform meshes used to generate the grid-

independent solution of this problem. Table 2 lists the average errors E; and

E; between successive levels. From Fig. 11, it can be seen that the average

relative error decreases very fast as the grid spacing diminishes.

In this study, the solution on level 8 grid is taken to be the grid-independent

solution. This level was selected because

1) e is taken to be 0.5%. E; in Table 2 is well below a;

2) The velocity profiles at different planes between successive levels were also
compared. Figure 10(a), (c), and (9) presents the U profiles at X = 0.15m.
The maximum relative errors of the solutions between level k and level k+1
for k = 1, 5, and 7 are equal to 12%, 0.4%, and 0.1%, respectively. Figure
10(b), (d), and (f) presents the V profiles at Y = 0.08m. The maximum relative
errors of the solutions between level k and level k+1 for k = 1, 5, and 7 is

approximately 19%, 1.1%, and 0.4%, respectively.

50

Absolute Error

 

O
O
O
N
I

 

 

 

 

 

hwuk
------ hwuke1
------ AnumneEnor

 

 

 

A IALLAJ

J A l l

.I

0.001

0.0006

 

 

-0.0025 0 0.0025 0.005
Tfifi 1 T_r v v . . r . . ' H—r—r—j
It —— Ievelk
0.1 r ------ levellt+1
- ------ RelatlveError
50.00 F‘
a I} ‘.
J r
0.04r’
\
I\
1 l 1 44:? fl. I 1 l 1 J +4]
-0.005 “0 0.005 0.01
0 (III)

a) U profile at X = 0.15m, k =1

Absolute Error

-0.0025 0 0.0025 0.005

0 1 f
' I

50.08
>4
0.04

-0.005 0
0 (ale)
0) Uprofile atX=0.15m, k=5

 

L L l

LlAL
0.005 0.01

Absolute Error
1 . ‘—v-+T . Um

7 Y V l

0.1x:-

50.08
>1
0.04

U (n/s)
e) U profile atX = 0.15m, k = 7

7 r

7 Vi' rrrTT—T

 

l L I

11111
0.005 0.01

 

 

 

4A1

 

l l l I ‘

ALIAAA

 

10.001

4

1

‘0 0008
O

4

<

4

40.0006

4

0.0004

lll‘l‘

4

€0.0002

 

 

. L_‘L_I s
0.2

(n)

‘0.“

.
.
l o

d) V profile at Y = 0.08m, k = 5

 

 

d
d
‘0

0.001

0.0008

LALIALJAJ

0.0006

A l A

0.0004

0.0002

JILAJLJIL‘A

 

 

f) V profile at Y = 0.08m, k = 7

Figure 10. Velocity profile comparison between level k and level k+1

51

0.0006

0.0002

Absolute lrror Absolute lrror

Absolute lrror

Table 1. Uniform grids used to get grid-independent solution for 2D test problem

 

 

 

 

 

 

 

 

 

 

 

 

 

Level Grid Spacing (m) Total
Cell no.
1 8.333x10':(1/120) 864
2 6.250x10' (1/160) 1,536
3 4.167x117'fi1/240y 3,456
4 3.125x10’3 (1/320) 6,144
5 2.083x10'3 (1/480) 13,824
6 1 .563x10'3 Q/640) 24,576
7 1 .250x10'f(1/800L 38,400
8 1 .042x10'3 (1/960) 55,296

 

Table 2. Average errors between successive levels for 20 test problem

 

 

 

 

 

 

 

 

Levels E... (m/s) E. (%)
1-2 1.144x10“l 7.29
2-3 1 .097x10" 6.70
3-4 3.676x10's 2.33
4-5 3.029x10'S 1 .63
5-6 . 1.221x10'5 0.605
6-7 6.666x10'6 0.330
7-8 4.175x10‘ 0.224

 

 

 

 

 

0.1 -

 

 

. 1 . . 1 . . . . 1.1-11.1. .fi—fij
0 10000 20000 30000 40000
Cell number

 

Figure 11. Average error level between successive meshes

52

Figure 12 shows the grid-independent solution represented by streamlines. The
flow pattern is symmetric about the X-axis because of its geometric and
boundary condition characteristics. As seen from the upper half of figure, a
strong jet-like flow goes through the whole flow field from the inlet duct to the exit.
One big anticlockwise recirculating flow exists at the top right corner and one
small clockwise recirculating flow is located between the jet-like flow and the
symmetric plane. Other complicated flow features include smaller recirculating
flows in corners, and flow separation at chamber inlet and comers. Thus, this

test problem does posses a number of complicate flow characteristics.

 

   

 
  

éc<((¢¢((c
W
m
I“‘(“¢“

3m

 
    

v
'4-

  

 

Figure 12. Flow pattern of grid-independent solution in 20 case

4.2.3 Imperfect Grids

Once the grid-independent solution is available, the next step is to generate a

series of imperfect grids with different cell aspect ratio, cell smoothness, and cell

shape. A total of 13 ‘imperfect' grids were generated. Tables 3 to 5 give the

53

characteristics of these meshes, respectively. In these tables, Xm & Yllr denote
aspect ratio along the X-and Y- axis, respectively. Ase/A. presents the volume
change between cell and its neighbors. X01 & Ya represents the skewed angle of

generated mesh according to X- and Y- axis, respectively.

Table 3. Meshes with different aspect ratio

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Level Xar Y,lr Am/A; No. of Cells / Nodes
l 10. Baseline 1.0 7074 / 7375
II 20. Baseline 1.0 5898 / 6185
III Baseline 3. 1.0 7938 / 8305
IV Baseline 5. 1.0 5670 / 6025
Table 4. Meshes with different smoothness
Coarse Smoothness No. of Cells / Nodes
I AXm/AX. =1 .05 12450 / 12815
ll AY.+1/AY.= 1.2 9742 / 9715
III AYi.1/AY.= 1.3 7362 / 7725
Table 5. Meshes with different cell Shape
Coarse X01 Ya No. of Cells / Nodes
I 45. 0. 9762 / 10095
II 45. 0. 6906 / 7205
III 0. 45. 10746/ 11125
IV 0. 45. 7776 / 8140
V 45. 45. 4806 / 5075
VI 45. 45. 3486/ 3740

 

 

 

 

 

 

Solutions for each imperfect mesh were generated for database construction.

54

4.3 3D Test Problem I: SO Flow In a Cubic Chamber with Four Curved Inlet

Ducts and One Exit Duct

In this section, the description and formulation for this 30 test problem will be first
presented. Follows are the details of how to obtain grid-independent solution

and generate imperfect solutions.

4.3.1 Problem Description and Formulate

A schematic diagram of this problem is shown in Fig. 13. It is a SD cubic
chamber with four curved inlet ducts and one outlet duct. The edge length of the
cubic chamber is equal to 150mm. The inlets are four curved square ducts
(25mm x 25mm). The outlet is a straight square duct (25mm x 25mm) with

length equal to 100mm. This geometry is symmetric about the X-axis and Y-axis.

The water is used as the flow fluid with constant 11 = 1.003E-Skg/m-s and
constant p = 998.2kg/m3. In Fig. 14, the inlet velocity is uniform. The flow speed
(0.01m/s) is determined so that the flow remains laminar. The standard
atmosphere pressure and the no—slip condition are still employed at the outlet

and all walls, respectively.

55

 

 

 

 

Figure 13. Geometry description of 3D test problem I

  
   
  
 

No-slip
W . 11

/
Velocny
Inlet

21‘

Pressure
0 tlet

Figure 14. Flow conditions of SD test problem I

As noted, this flow problem was selected because it possesses complicated flow
features such as secondary flow, flow separation, jets, and recirculation zones.
Because of symmetry in this problem, we consider only one quarter of the
geometry. Similarly as the previous 20 problem, we are only interested in the
steady-state solution of this incompressible laminar flow and solution to the

governing equations was obtained by using the Fluent code.
4.3.2 Grid-Independent Solution Generation

The procedure used to obtain grid-independent solution is the same as the one
described in Section 4.2.2. Table 6 gives the series of uniform meshes used to
generate this grid-independent solution. A total of 11 grids/solutions were
generated/computed and the cell number of the final grid-independent mesh is

535,680. Table 7 lists the average relative errors (also see Fig. 15) and absolute

57

error between successive levels. For this 3D problem, because of its difficulties
in computation, we define another error to further compute the difference

between solutions. This error describes the relative error in velocity magnitude
and can be expressed by

N (I VXJ I —Ivk+l,i/
2I IV... |

_k _ i=1
E (4.3)

re,ma —
g N

 

In this study, the solution on level 11 grid is taken to be the grid-independent
solution. This level was selected because

1) s1 is taken to be 2.5%. E}: in Table 7 is well below 81;

2) 52 is taken to be 1.5%. E12“ in Table 7 is well below 22;

3) The velocity profiles at selected planes between level 10 and level 11 varies

by less than 82;

Table 6. Uniform grids used to get grid-independent solution for 3D test problem I

 

Level Mesh Spacing A Inlet duct Chamber Exit duct Total

 

 

 

 

Cell no. Cell no. Cell no. Cell no.
0.00625(1/1 60) 608 1 ,728 128 2,464
0.004167(1/240) 2,016 5,832 432 8,280
0.003125(1/320) 4,864 13,824 1,024 19,712
0.002500(1/400) 9,800 27,000 2,000 38,800

 

0.002083(1/480 17,568 46,656 3,456 67,680

0.001785(1/560) 26,656 74,088 5,488 106,362
0.001563(1/640) 39,424 1 10,592 8,192 158,208
0.001 389(1 /7 20) 57,672 157,464 1 1 ,664 226,800
0.001250(1/800) 78,400 216,000 16,000 310,400
0.001136(1/880) 110,352 287,496 21,296 419,144
0.001042(1/960) 134,784 373,248 27,648 535,680

 

 

 

 

 

 

 

 

 

 

 

 

 

figmenxioamrswro-e

 

58

Table 7. Average error between successive levels for SD test problem I

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Levels 13}, (m/s) 3,, (°/.) Es... (%)
1-2 4.896x10“t 17.6 13.6
2-3 3.070x10“T 11.3 8.11
3-4 2.012x104 6.94 5.02
4-5 1.470x10‘4 5.15 3.63
5-6 1.167x10" 4.17 2.87
6-7 9.018x10’5 3.52 2.45
7-8 7.945x10'5 2.97 1.98
8-9 6.510x10’5 2.47 1.64
9-10 6.791x10'5 2.70 1.81

10-11 4.920x10'5 1.91 1.24

02

 

 

 

045

005

 

. l 1 A l I m _L 1 n I 1 1
0,612.00 2.01:+05 4.0E+05
Cell number

 

 

Figure 15. Average error level between successive meshes

Figure 16 shows the flow patterns at several slices of the grid-independent result.
In Fig. 16(a) (z = 55mm slice), the strong jet-like flow comes out from the curved
inlet duct into the cubic chamber. When it impinges on the opposite wall, the flow
is split to form two recirculating flows. Figure 16(b) shows the complex flow
structures at z = 12.5mm slice, which is located very close to the symmetry
plane. One big recirculating zone exists on the left part and one small

recirculating flow appears at the top corner. The pair of secondary flows formed

59

in the curved inlet duct can be seen to persist at the chamber inlet (Fig. 16(0)).
These secondary flows embedded in the jet make the flow in the chamber more
complicated. Figure 16(d) (x = 140m slice) depicts how the main flow behaves

when hitting the opposite wall.

 

(c) Slice X = 5mm (d) Slice X = 140mm
Figure 16. Flow pattern of grid-independent solution for SD test problem | (see

definition of X, Z in Fig. 13)

60

4.3.3 Imperfect Grids

Once the grid-independent solution was obtained, a series of imperfect grid are
generated. For this SD problem, a total of 9 ‘imperfect’ grids are generated with
different aspect ratio, smoothness, and cell shape. Solutions for each imperfect
mesh are generated for database construction. Tables 8 to 10 give the detail

information of these meshes.

Table 8. Meshes with different aspect ratio

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Aspect ratio N0. of Cells No. of Cells No. of Cells Total No. of
Level in Inlet duct in Chamber in Exit duct Cells
1 21,312 21,600 2,304 45,216
2 20,736 14,400 1 ,728 36,864
3 3,200 11,520 1,200 15,920
4 2,400 9,216 1,000 12,616
Table 9. Meshes with different smoothness
Smoothness No. of Cells No. of Cells No. of Cells Total No. of
Level in Inlet duct in Chamber in Exit duct Cells
1 22,176 21,600 4,608 48,384
2 15,400 1 1,520 3,200 30,120
3 6,400 21,952 3,200 31,552
4 6,400 13,824 3,200 23,424
Table 10. Meshes with different cell shape
Alignment No. of Cells No. of Cells N0. of Cells Total No. of
Level in Inlet duct in Chamber in Exit duct Cells
1 23,616 39,600 9,216 72,432

 

 

 

 

 

 

61

 

 

 

4.4 3D Test Problem II: 3D Flow In 8 Cubic Chamber with Two Opposing

Inlet Ducts and Two Exlt Ducts

In this section, the description and formulation for this 3D problem with opposing
flow will be first presented. Afterwards, the details of how to obtain grid-

independent solution and generate imperfect solutions are described.

4.4.1 Problem Description and Formulate

Fig. 17 presents a schematic diagram of this problem. This is a 3D cubic
chamber with two opposing inlet ducts and two outlet ducts. The edge length of
cubic is equal to 150mm. The inlets are two square ducts (25mm x 25mm). The
exits are also square ducts (50mm x 50mm). The length of all inlets and exits
are equal to 100 mm. This geometry is symmetric about the X-axis, Y-axis, and

Z-axis.

The fluid is still water with constant I1 = 1.003E-3 kg/m-s and constant p = 998.2
Kg/ma. Fig. 18 shows its flow conditions. The velocity profile is uniform with u =
0.02 m/s at two opposing inlets so that the flow remains laminar. The conditions

imposed at outlet and all walls are exactly same as previous SD problem.

62

 

Figure 17. Geometry description of 3D test problem ll

63

 
 

,, Pressure
Outlet

  

Pressure

Figure 18. Flow conditions of SD test problem ll

As noted, this flow problem was selected because it possesses complicated flow
features such as opposing flows, flow separation, jets, and recirculation zones.
Because of symmetry in this problem, we consider only one eighth of the
geometry. Similarly as the previous problems, we are only interested in the
steady-state solution of this incompressible laminar flow and solution to the

governing equations was obtained by using the Fluent code.

4.4.2 Grid-Independent Solution Generation

By using the same method described in section 4.2.2, we obtained the grid-

independent solution for this 3D problem. Table 11 shows the series of uniform

mesh. A total of 11 grids/solutions were generated/computed and the final level

contains 442,368 cells. Table 12 lists the average relative errors (also see Fig.

19) and absolute error between successive levels. The value of 81 and 62 is still

taken to be 2.5% and 1.5%, respectively. From Table 12, E}: is well below 61

and E” is well below 82. The velocity profiles at selected planes between level

running
10 and level 11 also varies by less than 82. Thus, in this study, the solution on

level 11 grid is considered as the grid-independent solution.

Table 11. Uniform grids used to get grid-independent solution for SD problem Il

 

Level Mesh Spacing A Inlet duct Chamber Exit duct Total

 

 

 

 

 

 

Cell no. Cell no. Cell no. Cell no.
0.00625(1/160) 64 1 ,728 256 2,048
0.004167(1/240) 216 5,832 864 6,912
0.003125(1/320) 512 13,824 2,048 16,384
0.002500(1/400) 1 ,000 27,000 4,000 32,000
0.002083(1 I480) 1 ,728 46,656 6,912 55,296
0.001 785(1/560) 2,744 74,088 1 0,976 87,808

 

0.001563(1/640) 4,096 1 10,592 16,384 131,072
0.001 389(1 /7 20) 5,832 157,464 23,328 186,624
0.001250(1/800) 8,000 216,000 32,000 256,000
0.001 136(1/880) 10,648 287,496 42,592 340,736
0.001042(1/960) 13,824 373,248 55,296 442,368

 

 

 

 

 

 

 

 

 

 

 

jaomwmmth-s

 

Table 12. Average errors between successive levels for SD problem ll

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Levels E... (m/SI E... (%l is... ("/o)
1-2 6.301x10‘4 16.9 10.3
2-3 4.292x10'4 13.5 7.82
3-4 2.844x10‘4 8.95 5.24
4-5 1.949x10'4 5.96 3.68
5-6 1 .488x10" 4.72 2.99
6-7 1.231 x1 0" 4.08 2.52
7-8 1.036x10'4 3.49 2.13
8-9 8.723x10'5 2.94 1.80

9-10 7.406x10'5 2.48 1.53

10-11 6.423><10°5 2.14 1.33

 

 

65

 

 

 

 

0.15

 

 

0.19E+00 1.0E+05 2.0E+05 3.0E+05
Cell Number

Figure 19. Average error level between successive meshes

Figure 20 shows the basic flow features of this grid-independent result. Because
of symmetry, only one-eighth of the entire flow domain is presented. Figure
20(a) shows the flow pattern at the slice y = 0mm. It is clear that when the strong
main flow enters the cubic chamber and hits the flow from the opposing direction,
it changes direction immediately. One big vortex is created at the top of chamber
and some complex separation flows also show up. In Fig. 20(b) (x = 0 slice), it
can be seen that the main flow turns its direction to the exit duct after hitting the
opposing flow. Figures 20(c) and (d) shows the z = 0 and z = 32.5 mm slices. In
these plots, some complex flow structures (e.g., recirculating regions, jets, flow

separation, etc.) can be found in the main chamber.

 

 

 

 

 

 

 

 

 

 

(a) Slice y = 0 mm (b) Slice x = 0 mm

66

 

 

 

(0) Slice 2 = 0 mm (d) Slice 2 = 32.5 mm

Figure 20. Flow pattern of grid-independent solution for SD problem lI

4.4.3 Imperfect Grids

Same as before, a series of imperfect grid is generated for the database
construction. A total of 14 simulations are computed here with different aspect
ratio, smoothness, and cell shape. The detail features of these meshes are
summarized in Tables 13 to 15. Solutions on each of imperfect mesh are

computed for database construction.

Table 13. Meshes with different aspect ratio

 

67

Table 14. Meshes with different smoothness

 

 

 

 

 

 

 

 

 

 

 

 

Smoothness No. of Cells No. of Cells No. of Cells Total No. of
Level in Inlet duct in Chamber in Exit duct Cells
1 512 23,040 3,072 26,624
2 512 24,192 2,560 27,264
3 512 19,008 2,048 21,568
4 640 28,800 3,840 33,280
5 640 25,200 3,200 29,040
6 512 24,960 3,072 28,544

 

 

Table 15. Meshes with different cell shape

 

Alignment No. of Cells No. of Cells No. of Cells Total No. of
Level in Inlet duct in Chamber in Exit duct Cells

 

1 1 ,536 24,576 4,608 30,720

 

 

 

 

 

 

2 1 ,536 20,480 4,608 26,624

 

 

68

Chapter 5

Evaluation of Grid-Quality Measures

As noted in Chapter 4, three test problems with complicate flow features were
selected to evaluate grid-quality measures. For each problem, a grid-
independent solution was generated. Also, for each problem, a series of
imperfect grids and their corresponding solutions were generated. Solutions
obtained on these imperfect grids are not grid independent, and hence have
errors. All solutions generated, solution errors, and newly developed grid-quality
measures can be used to construct a database. In this chapter, we describe the
details of this database. Aftenivards, we show how measures are evaluated by

interrogating this database.

5.1 Database: Construction and Organization

In this study, the procedure used to construct the database is as follows:

1. Generate a grid-independent solution for each problem. In this study, grid-
independent solutions were generated for three test problems (see Section
4.2, 4.3, and 4.4).

2. Generate a series of imperfect solutions for each problem and obtain
solutions on these grids. As noted in Sections 4.2 to 4.4, 13 different

solutions were generated for the 2D test problem, 9 different solutions were

69

generated for the 3D test problem I, and 14 different solutions were generated
for the SD test problem II.

3. Calculate the absolute error, the relative error, and the change in relative
error at each cell of every solution by using Eq. (3.1).

4. Calculate the value of measures at each cell of every solution by using Eqs.

(2.1) to (2.8).

At the conclusion of these four steps, we have a database made of three sub-
databases, one for each of the three test problems. Each sub-database has the
following information:

DB, = (Eab, Era, A(E..), o1, o2, , on)” (5.1)
where, j = 1, 2, 3, represents each sub-database. n is the total number of
measures. i = 1, 2, N and N is the total number of cell considered in all

solutions for a given test problem.

The goal in evaluating the data given by Eq. (5.1) is to determine how well the
errors -- Eab, Em, A(E..) — correlate with the measures -- o1, o2, , o... This
database can be interrogated by sorting and then grouping the data into a finite
number of partitions. The correlation between errors and measures can then be

examined by using statistical analysis.

When data is partitioned, two types of organizations were used. The first type

considers only one measure at a time; i.e.,

70

E”: F(¢,), i= 1, 2, 11 (5.2a)
E..= F(l>.).i= 1. 2. n (5.2b)

A(Ere) = FI‘I’i ). I = 1, 2, ..., n (5.20)

Since Eq. (5.2) is one-dimensional, the data can be partitioned according to the
value of the error or according to the value of the measure (i.e., the independent
coordinate can be either error or measure). The reason for using Eq. (5.2) is that
we want to examine which measure can exert influence on the local errors. The
goal is to find a set of measures that have strong correlation with errors. The
strong correlation between error and measure means that in partitions with low
error, the measure should have values close to its ideal value. Conversely, the
measure should show very bad values in partitions with the high error.
Furthermore, this correlation should be universal at least within all the

simulations in each test problem we studied.

Figure 21 summarizes what we did in this study to find all the possible correlation
between measures and errors. As seen from the figure, for each measure, we
make the examination by accessing data from one simulation or from all
simulations of each test problem. In each individual examination, data is
partitioned according to different types of error. When Eq. (5.2) is used, there is
a distribution of data within each partition. In each partition, all the following
statistics parameters are computed: mean, median value, mode (most-probable

value), variance, skewness, and kurtosis. In this study, the data descriptors

71

(e.g., mean, median, and mode) within each partition are used to correlate with
error. After correlation between measures and errors have been found, other
parameters (e.g., variance, skewness, and kurtosis) are calculated for further
examination of the statistical features in each partition. If the variance is low,
then the measure being studied is considered dominant in determining the local
error. Otherwise, it means that there are other factors affecting local error, or the

error cannot be determined by using only one local measure.

Therefore, for those measures that have the strongest correlation with error, the
second type of partitioning will be performed. It considers two or more measures
at a time, i.e.,

E = F(¢i. <12. 011) (53)
where M is the number of measures. If each of the M measures in Eq. (5.3) are
divided into N monotonically increasing intervals, then the number of partitions is
equal to N”. For Eq. (5.3), linear least square is used to represent the data in
each partition. If there are M measures, then Eq. (5.3) has the following form in
each of the partitions:

E=an+a1o1+a2o2+~+aMoM (5.4)
where a, a2, at). are the coefficients. By using linear least square for each
partition, the function E = F(o1, o2, o...) is discontinuous at the boundaries of
the partitions. The statistical parameters are then calculated for checking the

statistical features in each partition. If the variance of the data in each partition in

72

this step is low, it means that the local error can be determined by the local

measures.

In order to determine grid-quality measures that have correlation with grid-
induced errors, extensive numerical experiments were conducted based on the
database described above. In these experiments, many grid-quality measures
that combine information about the geometry and the solution were developed.
For each of these measures, a detailed statistical analysis was performed to
check its correlation with error. Of all the measures developed and evaluated,
only six measures described in Chapter 2 were found have some correlation with
error. In the following sections, we will describe the details of the evaluation of

these measures.

73

 

 

Grid-Quality Measure
Study
I

 

 

 

 

 

 

 

 

r
Types of Error

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

l
Solution-Based Test Problems Data
Measures Statistics
Directional _ 2D L Absolute Error _ Average
Resolution
Gradient , , Relative Error Median
Resolution _ 13 Simulations F
_ Deformation , , i _ Change in Mode
Tensor Resolution .. COWb'Mt'P" 0' Relative Error
- all Simulations -
Directional ’ Standard
Smoothness 1. 30 Variance
Problem I ,
Streamline Skewness I
Smmothness _ 9 simulations I g
I. Relative Velocity Kurtosis l

 

Difference

 

 

 

 

 

 

 

 

 

 

 

 

 

 

_ Combination of
all simulations

 

 

 

 

 

 

 

e SD
Problem II

 

 

14 simulations]

 

 

 

_ all simulations

 

Combination ofl

 

 

 

 

 

 

 

 

 

 

74

Figure 21 . Summary of grid-quality study

 

 

 

 

 

 

 

 

5.2 Evaluation of Error with One Measure at a Time

In this section, the correlation between error with one measure at a time is
evaluated. one is used to introduce the detailed evaluation process. Afterwards,

the results about how other measures correlate with error are presented.

5.2.1 Evaluation of Directional Resolution

In this subsection, we describe the procedure of correlating one with different
types of errors, present the correlation curves in different 2D/3D simulations, and

also analyze the statistical features in different partitions.

5.2.1.1 Correlation with Different Errors

Figure 22 shows how the directional resolution (one) correlates with different
types of errors. In this figure, each curve represents one sub-database. As
mentioned, when using Eq. (5.2), each sub-database can be partitioned
according to the error value or according to the measure value. Here, we
partition according to the value of errors first. An equal-population approach was
employed to partition each sub-database. With this approach, each partition is
an interval on the error coordinate line. The width of each interval, however, is

variable because it depends on how the error is described in the sub-database.

75

For each interval of the partition, we assign a single value. Here, we choose the

mean value to represent the central tendency.

For the correlation between one and the absolute error, which is shown in Fig.
22(a), all three curves changes randomly and do not show any trends. Figure
22(b) gives the relationship between one and the relative error. For the 2D test
problem, the mean value of one is less than 10 when E... is small (less than 5%).
As E... grows, one was found to increase monotonically. For the 3D test problem I
and II, similar results can be seen. The trend of the correlation curve between
Ere and one is almost monotonic and linear except the last partition for 3D
problem I. Figure 22(0) shows the result between one and the change in relative
error, A(E,e). Looking at the curve of SD problem I, the curve remains almost flat
when one is smaller than 10 and increases linearly and monotonically until one
reaches 28. Then the slope of curve becomes very sharp when one continues to

increase. For the results of 2D problem and SD problem II, the similar results

can be found. Compared with the figures of Eat, and En, the correlation between
ME.) and measure one seems not so strong as that for E... but much better than

Eab does.

76

 

   
      
    

2D Problem
30 Problem 1
30 Problem 11

—9— 2D Problem
—*— 3D Problem 1
l—— 30 Problem 11

  

 

 

3.5E-03
305-03
2.5E-03

,9 2013-03
1 .5E-03

1 .05-03
5.0E-04
0.0E+00 0

1.0 :-

 

Ere

 

 

DR

 

(a)

 

‘
I r rrr WT]

1—-—a— 2D Problem
._,e,_ 3D Problem 1
,_.._... 3D Problem 11

 

 

.0
at

 

.0
O)
T 7r

A(Ere)
P
M
°.'T7T'1j1fiirIT Y

.0
h

 

P
o

(C)

Figure 22. Correlation between one and errors, Partition by Error

(3) 5.1, = F(¢i )- (b) E... = F(¢i )- (C) ME.) = FOP.)

Next, we continue to seek out how the measure one is sensitive to the different
types of error by partitioning the data according to the value of measure one.
Figure 23(a) shows that there is still no clear relationship between one and the
absolute error in all problems. For the correlation between one and E... (Fig.
.23(b)), when one is small, the mean value of the relative error is small. When one

increases, E... goes up monotonically except some partitions in 2D problem.

77

Similarly, in Figure 23(0), A(Ere) increases in an oscillatory manner as one grows

and the trend of correlation curve is much better than that of Eab.

   
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1.0E-03 0.6 ~
L
8.08.04 _
i 0.4 P
60504 I. 1—8— 2D Problem
.0 . 8 _ :- 3D Problem I
13 ; m : 30 Problem 11
4.0E-04 :- ’
. 0.2 ’.,
: ._9_ 2D Problem
2.05-04 - —é— 3D Probleml ,
I ._._ 3D Problem 11 . .r
W1 1 l l l 0 - i L i
0.0E+000 20 40 60 80 100 120 0 20 40 60 80 100 120
DR ¢ DR
(3) E86 (b) Ere
0.3

 

0.2
2'

0.1
21) Problem
30 Problem 1
31) Problem 11

O 0 20 40 60 80 100 120
9 on
(C) A( Ere)

Figure 23. Correlation between one and different errors, Partition by one

By comparing the results in Figs. 22 and 23, it can be seen that whether the error
or the measure is used for partitioning, the relative error has the strongest
correlation with measure one. Also, the relative error can be problem
independent and is more meaningful in CFD analysis. So, in the following

sections, the relative error is considered as the default solution error for measure

78

evaluation. For the convenience of comparison for different measures’ behavior,

the data will be partitioned by the relative error in later section.

5.2.1.2 Evaluation of Generality for all ZD/SD Simulations

Though Figs. 22 and 23 show a strong correlation between the relative error and
one, we want to further examine these correlation. In constructing Figs. 22 and
23, the entire sub-database was used for each test problem. In this sub-section,

we examine each individual simulation.

Figure 24(a) shows the result for the 20 test problem. For each of the 13
simulations performed, a E = f(one) curve was constructed. Each “thin” curve
represents one simulation with different aspect ratios, smoothness, or alignment
characteristics. As seen from the curve, the correlation between E... and one in
all the simulations looks similar to the combined result, which is represented by
the “thick” solid line. But data spread in these curves is still big. Figure 24(b)
presents the result for SD problem I and there are a total of 9 simulations. The
generality among simulations shows that the 3D cases are much better than 20
case. Almost all the curves match together except two simulations and the data
spread in curves is very small, especially in the low error intervals. The result of
3D problem II also has the similar trend. A total of 14 simulations is included in
the Fig. 24(0). The correlation curves of all these simulations almost match

together and the spread in curves is very small.

79

 

 

 

 

 

 

 

 

 

 

 

 

   

 

 

 

 

 

 

 

:1—8-— All chases I ._g_ All9cases
: --+-'_ Alignment .\ : ————— Aumt '\\ \\ '
0.8 I: ------ Aspect Ratio \ 08 '_' ------ Aspect Ratio I \\:"e 1
: —----~-~~—---- Smoothness I : --~-~-o-------- Smoothness ,’ . I
0.6 :- I
D e l 4 A
1.13 In ; '
0.4 :-
0.2 I-
0 ’ .
0
(a) 2D Problem (b) SD Problem |
1 .-
i 1—3— All ”cases
i -—+—— Alignment
0.8 :' " " ”*- ‘ - ASWRadO 4/2 ‘- [1”
; -------°----~-- Smoothness ”(g/)7" 7’
0.6 :- '
o .
ii ;
0.4 :-
0.2 :-
041414%.1L111.l..411
0 10 20 30 40

DR

(0) SD Problem ll

Figure 24. Generality in all simulations

Above results indicate the general trend of correlation between one and the
relative error in all simulations. Though certain data spread still exists in the
correlation curves, it is expected since only one measure was examined at a
time. For example, though one may be near ideal, the local grid quality is still bad
if the mesh spacing is too coarse or the mesh smoothness is too high. So the
detected spread in data implies the need to examine two or more grid-quality

measures at a time in the form of Eq. (5.3). This is described in Section 5.3.

80

5.2.1.3 Evaluation of Different Data Descriptor - Mean, Median, Mode

According to statistical concepts, mean value is one of the most popular central
tendency measures. It is known as the most important and commonly used
descriptor of data. For n points in each partition, if all observations are given

equal weights, the mean value is given by:

(5.5)

><|
ll
3 I 1—
M
'34

where xi, = a sample point, and k = 1, 2, n.

The other two types of central tendency descriptors are median and mode. The
median xm is defined as the point that divides each division into two equal parts,
that is 50% of the data are above xm and 50% are below X... Then median can
be determined by ranking the n values in the division in decreasing order, 1 to n.
If n is an odd number, the median is the value with a rank of (n+1)/2. If n is an
even number, the median equals the average of the two middle value, which are
with ranks N2 and (n/2)+1. The mode xd is defined as the point of the highest
percent for the frequency of occurrence. Thus, mode represents the most
probable value. In order to determine this point, each division will be sub-divided

and the subdivision with highest frequency of point will be sought out.

In previous analysis, the database is partitioned into small intervals according to
the relative error and only the mean value of the measures is used to represent

each partition. In order to evaluate the performance of measure, the mean,

81

median, and mode were calculated to represent the data tendency in each

partition.

Figure 25(a) shows the results of the 2D problem. The data is still partitioned by
the relative error but three data descriptors are computed instead. When Er9 is
small, the median value of (boa ranges from 0 to 80. When the error is high
(greater than 10%), $03 is almost fixed at 85. From the mode value curve, om
does not change a lot when the error increases. When the mean value is
compared with these two descriptors, it is clear that the mean value is the best
descriptor for data in each partition and has the strongest the correlation with
solution error. Of course, the other two descriptors also provide some useful
information about the actual data distribution in each partition. The mode value
curve locates at the left side, so it means that the most probable value is very
close to the ideal value of ¢on- The mean value and median value are far from
each other in each partition, it indicates that the data distribution is highly

skewed.

The result of the 3D problem I is given by Fig. 25(b). The mean value curve still
shows the best correlation with the error. The median value curve locates at
right and changes in an oscillatory manner. The mode value curve remains at
left and does not correlate well with the error. Similar results can be found in Fig.

25(c) for the 3D problem ll database.

82

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1 1
—e— Averagevalue A
---b-- Medianvalue .
0.8 -—-o--- Modevaluc :
4.
l
I
g l—a— Averagevalue?
LL] ---¢5--- Medianvalue A
---o--- Modevalue is
:3»
.Q'
0:25“
A:--%
l 7* . A 1
4O 60 80
¢mz
(b) 30 Probleml
1-
E o f] ?
0.8 f : I
. n l—e— Average value :
I p ---A--- Medianvalue A
0.6: I: ---o--- Modevalue ,'
2
m P f ",.A
0.4- ° A---_
: 2 f
L' A----
l- A1
0.2 - 14%---3
Ann-L"-
E A,,__A.~:-:-3'S°
0 Ljéfjgl....1... I
0 20 4O 60 80
¢DR
(0) 30 Problem ll

Figure 25. Comparison of different data descriptors

According to above tests, it is obvious that when different descriptors are used to
represent data partition, their behaviors are much different from each other and
the mean value correlates best with the solution error. The difference among
these descriptors also implies that though the mean value can convey certain
information about the data sample, it cannot completely characterize the feature
of each partition. So, it is important to further study the detailed statistical

characteristics of data in each division.

83

5.2.1.4 Analysis of Statistics in Data Partition — Variance, Skewness, Kurtosis

Standard deviation is the square root of the variance. It has the same unit as the
variable and the mean. Therefore, it is a good indicator of the dispersion or

spread of the data distribution. The standard deviation of the data points is

denoted by a and expressed as:

 

 

a=‘[ 1120,35): (5.6)

’1‘ k=l

where x.. = a sample point, and l = 1, 2, n., if is the average value in data

sample.

Additionally, according to statistical theorem, the moments are very useful
descriptors of data. A moment can be referenced to any point on the
measurement axis. The origin and mean value are the most common reference
points. Generally, the mean value can be considered as the first moment about
origin and the variance can be defined as the second moment about the mean.
Although most data analysis only considers two moments, it is important for
some probabilistic and statistical studies to examine the third moment about the

mean, which is Skewness. Skewness can be mathematically defined by

n

 

n ZCX,‘ - Y)3
C = "=‘ 5.7
5 (n—l)(n-2) 0'3 ( )
where xk = a sample point, and l = 1, 2, n., )7 , o is the average value and

standard deviation, respectively. Equation (5.7) and has unit of the cube of the

84

random variable. It is a measure of the lack of symmetry. A symmetric
distribution has skewness. of zero. Furthermore, depending on the direction of
skewness, a non-symmetric distribution can be positive or negative. When the
peak value is at the right side, the skewness is positive. The skewness is

negative when more data is locates to the left of mean.

A further type of departure from normality is called kurtosis. In a data population,

it can be expressed as

 

n _"X' 4 n _ —' 2 2
n2 -2n+3 Em ) 3(2n-3) lg“ X) 1
Ir = 4 - 4 (5'8)
(n -1)(n — 2)(n — 3) 0' n(n — l)(n — 2)(n — 3) 0'
where xk = a sample point, and l = 1, 2, n., 5(- , o is the average value and

standard deviation, respectively. For the normal distribution, this ratio has the
value 3. If ratio is less than 3, the curve has a flatter top than normal. If the

value exceeds 3, there is usually a sharper top than normal.

The mean value in each partition is shown in Fig. 26(a) for reference. Figure
26(b) shows the data variance in each partition. It can be seen that the standard
deviation has similar order of magnitude as the mean values. This means that
the data variance in each partition is very large. Figure 26(c) presents the
behavior of skewness. As observed, the skewness remains positive in all
partitions. So, the distribution of skewness in all intervals is non-symmetric and
at the right side. When error is low, the skewness changes randomly. When

error is high, the skewness value is closer to zero, which indicates symmetric

85

distribution. Figure 26(d) gives the results of kurtosis. The kurtosis values are
higher than 3 in most of data partitions. This means that the data distributions in
those divisions highly departure from the normal distribution and has a sharper
top than normal. As the error increases, the kurtosis value reduces. In other
word, the data distribution in higher error partitions should be flatter than that of

the lower ones.

 

 

—e-— 2D Problem —9— 2D Problem
—a—— 3D Problem I —— 3D Problem I
—— 3D Problem 11 ——0— 3D Problem II

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 
     
   

 

 

 

 

 

 

 

 

 

 

(b) Variance
1 - 1 ..
E 21) Problem i —e— 20mm
I —a.— 301:“,ka I —-t— 3DProbleml
0-8 : »——.— 39mm“ 0.8 _- l—o— sommn
0.6 :- - :
. 0 .
5 C 5 :
0.4 :- . _
0.2 E- . : k
. g, o * . ' 7
o 6 o 25

75 100

2 A L 4 so
503(¢oa) KO ’(4. on)
(c) Skewness (d) Kurtosis

Figure 26. Comparison of different statistical parameters

86

5.2.1.5 Detailed Data Distribution in Each Partition

To further examine the statistical characteristics in different data partitions, Fig.
27 gives the detailed data distribution within selected partitions. In this figure, X

axis represents (pm and Y axis represents the percentage of cells that falls into

certain ranges of (Don.

To illustrate the meaning of Fig. 27, consider the 3D test problem II (see Fig.
27(0)). In this figure, it can be seen that in the partition with the lowest error
(denoted as #1), the data distribution has sharp kurtosis and is highly skewed to
the left. Also, the data variance is small. In partitions with higher error (denoted
as #7 and #13), the distribution curve moves to the right. The sharpness and
skewness in data sample are reduced, whereas the data variance increases. For
partition with the highest relative error (denoted as #20), the distribution curve
becomes much flatter and less skewed compared to other partitions. It also has

the largest data variance among all partitions.

Figures 27(a) and (b) show the results for the 2D test problem and 30 test
problem I, respectively. From these figures, it can be seen that in partitions with
low error, the data generally has low (pm, and the data variance is small. In
partitions with high error, the variance of (Don is large. It means that in some flow
regions, though the local (pm is perfect, the local relative error still can be very

high. Why? One possible explanation is that the local error is affected by

87

multiple measures instead of one measure. For example, though the local 4);»; is
perfect, the local grid quality can still be poor if other measures such as

resolution and smoothness are not good enough.

 

 
  
  

 

 

 

 

 

 

O 8 0.3
.. -- —--o--- Division”
——-v— glivrsion :23 i. ——o— Division #13
4 ---o--- Vision 0.2 ---o--- Divisionm

  

 

9
OJ
TTYf'Tj—rfififY'

 

 

 

 

 

 

 

   
 

59 0.4 89
0.1 .
0.2 L-
E I .
0 0 ‘ 0 o A _ 40
(a) 20 problem (b) 30 problem I
0.3 - Division #1
E - - -A- - - Division #7
_ l—v— Division #13
- - -o- - - Division #20

 

 

  

 

 

(c) 30 problem lI

Figure 27. Detailed data distribution in selected partitions

As discussed, cog combines the geometry and solution information, and its
definition enable it to be applied to structured and unstructured meshes in 2D
and 3D. (Don is also non-dimensional. Finally, it correlates well with local

solution error. All of these properties make it an ideal grid-quality measure. We

88

just finished the evaluation of (ion. In the next section, we present results of all

other measures.

5.2.2 Evaluation of Remaining Measures

In Chapter 2, we presented six solution-based grid-quality measures, which are

. Local, solution-based measure: (boa. «pom

. Macro, solution-based measure: (boa. ¢os

. Physical, solution-based measure: (1)33, cgvo

These measures possess these characteristics important to grid-quality — cell
resolution, cell smoothness, and cell shape. than and com are related to the local
cell resolution. (two describes the cell resolution perpendicular to streamline. (boa
presents the macro cell resolution feature in recirculating zones. ¢os can be
considered as the description of local cell smoothness and (1)35 is specified for the

smoothness feature along streamline. Here, we classify measures (those

developed here and those geometric ones developed by others) into three
groups, which are resolution-type (¢on. (boa. (born. and (pm/o), smoothness-type
(geometric smoothness, (bus, and ¢ss). and cell-shape type (geometric
skewness). Each of these measures will be evaluated by using the same
statistical process described for ¢on and their correlation with the local solution

error will be shown.

89

5.2.2.1 Resolution-Type Measures

¢on accounts for the direction change of velocity vector and belongs to the
resolution-type measure. It has already been addressed in Section 5.2.1.
Besides the regions with rapid change in flow direction, the regions with high
solution gradient also needs more grid nodes. In Chapter 2, L33 is called the
effective grid spacing and defined as the largest distance spanned by a cell in a
direction perpendicular to the velocity since most velocity gradients are typically
largest in that direction. Figure 28(a) shows the correlation between Lea and the
relative error. Each database is still partitioned according to error and there has
nearly identical population of cells in each partition. The mean value of Lea is
computed to represent its tendency of each partition. As seen from this figure,
Lea has certain nonlinear and monotonic correlation with the error for all
problems. The only problem is that Lea is a dimensional measure and needs a
length scale to become dimensionless. Then, than is defined as Lea over a
characteristic length and referred as the non-dimensional gradient resolution.
Figure 28(b) presents the correlation curve between (pea and errors. For ZD test
problem, the correlation is good except for some small oscillations. For SD test
problems, the correlation is only good for low error regions and behaves

randomly when the relative error is greater than 10%.

(Porn is another measure that gauges resolution. Its mathematic definition and

physical meaning are defined in Chapter 2. Figure 28(c) shows its evaluation

result. It can be seen that the correlation of 2D problem shows very good
monotonic behavior. But (born does not correlate so well with error in 3D
problems. In these figures, (ppm correlates with Er9 only when error is higher than

20% and is almost inverse proportional to Ere when Em (less than 20%) Is small.

<l>nvo can also be considered a resolution type measure and describe the local
gradient behavior by computing the relative velocity difference perpendicular to
streamline. Figure 28(d) shows the relationship between (have and error. In this
figure, the correlation is very similar to that of (pom. Only the curve of 2D problem

shows the monotonic feature, while the curves of 3D problems behaves in an

oscillatory manner.

 

 

 
 
  
  
 

 

 

 

 

 

 

 

 

 

 

   

 

P —9—— ZD Problem
t —e— 21) Problem . ._._ 31) Problem]
F —&— 3D Problem] 1 T ._._ 3D Problem II
0.8 F .—o— 3D Problem 11 ;
0.8 :-
0.6 :- I
o l. L
La : 00.6 :
0.4 _- til :
- 0.4 :-
02 T 0.2 L-
i , :1,“ _
or '.;::""—'-I'.“'l°‘**4“:l 0" 1 ...1....1Dr..|
0 0.002 0.004 0.006 0.008 0 25.03 45.03 35.05 BEE-06 15.05
LGR (m) GR
(3) Ere ~ LGR (b) Ere ‘7 ¢GR

91

 

 

 

 

 

 

 

 

 

    

 

 

 

 

w—e— 21) Problem I ——8— 20 Problem
”—4—— 30 Problem] : 4‘. 3D Probleml
I ._._ 31) Problem 11 0.8 _- ——._ 31) Problem 11
0.6 f
5.3 ;
0.4 -
0.2 -
0' __-.»..;1g.111.. 31
o 02 0.4 0.6 0.8
MM)
(C) Ere ~ (born (d) Ere ~ ¢Rvo

Figure 28. Evaluation of resolution-type measures

5.2.2.2 Smoothness-Type Measures

The grid smoothness is another important factor that can affect the local grid
quality. The simplest one is the traditional grid smoothness, which is defined as
the maximum volume ratio of the current cell and all its neighbors. This one only
includes the geometry information instead of any local solution behavior. Figure
29(a) shows its correlation curves in all cases. Obviously, it does not show any
relationship with the solution error, which indicates that this geometric measure is

not very helpful to estimate the grid-induced error.

Figure 29(b) is the result of ¢Ds. which is defined as the changes in effective grid
spacing along the direction of the velocity. As seen from this figure, ¢os performs
much better than the traditional one since ops increases monotonically as the

error goes up and this trend can be found in all curves.

92

¢ss is another smoothness type measure, which describes the smoothness
feature along streamline. Figure 29(c) shows the relationship between ties and
error. This correlation is similar to that of (has. but it is not so monotonic since

some small oscillations exist in the curves of 3D problems.

 

      
   

 

  

a 2D Problem 2]) Problem
*—‘— 3D Pmblcm l 31) Problem 1
1 :- .__.__ 3D Problem 11 1 3D Problem [1

 

 

 

0.8

0.6

0.4

0.2

 

  

lllllALllll

1.5 2 2.5““3 1 1.2 1.4 1.6

 

Smoothness Ilibs
(a) Ere ~ Geometry Smoothness (b) Er.a ~ ops

20 Problem
3D Problem I
30 Problem [I

Ere

 

0.1 0.2 0.3 0.4 0.5 0.6
¢ss

(C) Ere ~ ¢ss

Figure 29. Evaluation of smoothness-type measures

93

5.2.2.3 Cell Shape-Type Measures

In order to find factors that can fully determine the local grid quality, the
traditional skewness, which only accounts for the geometric characteristics of
local cell shape, is also evaluated. Figure 30 is the result of the equiangle
skewness. Same as other geometry only measures, skewness does not show
any correlation with solution error in all cases. It indicates that the cell shape is
not so critical for local grid quality as the resolution-type and smoothness-type

measures.

 

—e——— 2D Problem
—&— 3D Problem]
—0— 3D Problem 11

 

 

 

0.8 :-

0.6 }

Ere

0.4 :-

02 L

  

 

l 1 1 1 1 l 1 1 1
0 0.05 0.1 0.15 0.2
Skewness

Figure 30. Evaluation of skewness

5.2.3 Summary of Evaluation One Measure at a Time

In this section, a set of grid-quality measures, which describe the characteristics
of grid resolution, grid smoothness, and cell shape, is evaluated. These
measures include the newly developed solution-based measures and some
geometry only measures. According to the results by correlating the relative

error with one measure at a time, the directional resolution rim and the

94

directional smoothness (Dos have the strongest correlation with solution error. We
regard them strongest because the correlation curves of (boa and thus with error
are monotonic, and this feature is nearly universal for all three test problems.
The streamline smoothness ¢ss behaves similarly to (tips but its correlation curve
still has small oscillation. Cell shape-type measure, like skewness, does not
seem to be a dominant factor for error estimation when compared with
resolution-type and smoothness-type measures. Several measures that relate to
the local gradient feature, such aspen, (born. and ¢nvo. were examined. But they

do not show so strong relationship with the local relative error as do (Don and cos.

Even for $011 and (Dos. which has the best relationship with solution error, there still
exist large variance in each data partition. This means that the correlation
between error and (Don or ¢os is only strong on the average basis. The error still
cannot be accurately estimated by (boa or ops. One possibility is that the local
error is affected by multiple measures instead of just one measure. So, in next
section, we will examine the correlation between solution error and multiple

measures.

95

5.3 Evaluation of Error with Mum-Measures

5.3.1 Evaluation of Errors with Two Measures at a Time

Since (Don and thus are the best measures based on the analysis in section 5.2, in
this section, we will first examine Ere = f((bon, $05) by using Eq. (5.4) and see if the
large variance in each partition can be reduced. Later, we include more

measures in the function.

5.3.1.1 Contour Pattern

Figure 31 shows the results of Em = f(oon, ops). In this figure, lion and «fins
represent the X-axis and Y-axis, respectively. The relative error is represented
by the gray scale contour. Each of the two measures is divided into ten intervals
by using the equal population approach and the total number of partitions is
equal to 102. Then the linear least square method is used to represent error level

in each partition.

Figure 31(a) shows the results of the 20 test problem. It can be seen that the
error is very small (less than 10 to 15%) when lien is less than 4° and cos is less
than 1.15. When lion increases and cos is fixed, the error level rises
monotonically. On the other hand, if ops increases and 4199 is fixed, error

increases with some oscillation. The roughness and jaggedness of figures is a

result of using linear least square for each partition and not having enough

96

partitions from a lack of data to span the spectrum of the grid-quality values.
Figure 31 (b) shows the evaluation result for SD test problem I. When (pm is very
small (less than 4°), the error is normally lower than 5% except for (has higher
than 1.3. When ¢on changes from 4° to 12°, the error increases from 5% to 20%.
When «lion continues to grow, the error becomes even higher. In other word, (boa

seems more dominant in determining error level than (lbs. The similar result can

also be seen in Fig. 31 (c) for result of SD problem ll.

 

(a) 2D problem (b) 3D problem |

 

(c) SD problem II

Figure S1. Er = f(¢oR. (lbs)

97

5.3.1.2 Mean Value in Partition

In this subsection, the mean value of relative error in each partition is used to

describe the correlation Em: f(lpon, ops). As seen from Fig. 32, the X-axis

represents the partition level of (boa & ops and the Y-axis represents the mean
value of error in each partition. Each of the two measures is divided into ten
intervals, which are same as the data partitions used in the previous section. In
this figure, the legend (1-1) represents the partition with ideal lion and $05 value

while (10-10) represents the partition with the worst (ion and $03 value.

Figure 32(a) is the result of the 2D test problem. As can be observed, when (Don
and cos are both ideal, the mean value of error in this partition remains very low
level (approximately 1%). When ¢on and cos change from level (2-2) to (10—10),
the error increases monotonically. Figures 32(b) and (c) show the result of 3D
problems, respectively. The monotonic trends of correlation between error and

measures are also very clear.

 

 

0.25: 05:.
1 10—1 '
I: R 05: 104:
0.2 '_- ' ;
l C
» 0.4:-
0.15E : .2
0 O I' H
a E ‘5 0.3.? 545“ 7A7 A
0.1 I.- 9.9 : 22 34 H A
. A 02— I A A
~ M : "1
005’- 55 H 71 A : ‘
' ~ ' 0.1-
1 2- 3 .
4411111 .
WA ..WWM
0012345678910 0012345678910
to... 110,. Level om, ops, Level
(a) 2D problem (b) 3D problem I

98

 

012345678910
om, (DDS, Level

(c) 3D problem II

Figure 32. En = f(¢oa.¢os)

5.3.1.3 Distribution Feature in Each Partition

In above evaluations, the function Ere = f(¢oa, 005) is examined by using linear
least-square approach and the mean value of error in each partition, respectively.
Strong correlation between the error and measures can be found. Next, we will

check if the data variance in data partition gets reduced as we expected.

Figure 33 shows the detailed distribution in selected partitions for all cases. In
this figure, each curve represents the distribution of error value in one partition,
and the partition level (1-1), (4-4), (7-7), and (10-10) are selected. The X-axis
represents the actual value of Er.a and the Y-axis represents the percentage of
cells that fall into certain ranges of error values. To illustrate the meaning of Fig.
33, consider the SD test problem I (see Fig. 33(b)). As observed, when (Don and
003 both have the ideal values (denoted as 1-1), the distribution of error has

sharp kurtosis and is highly skewed to the left. The data variance in this partition

99

is small. When the value of fine and 005 get worse (denoted as 4-4 and 7-7), the
distribution curve moves to the right. The sharpness and skewness of
distribution are reduced but the data variance is still very large. When then and
095 both have the worst values (denoted as 10-10), the curve becomes very flat
and least skewed, whereas the data variance is largest among all partitions.
Similar results can be found in figures of 2D problem (see Fig. 33(a)) and SD

problem ll (see Fig. SS(c)).

 

 

 

 

 

 

 

 

 

 

 

 

 

 

0.8 - 0.4 .-
. 4 i
. a Division #1 '1 i —-8— Division #1 -l
0.6 - “‘0'" DIVPIWM4 0.3 . ---o--- DivisionM-«f
~ —v— Division #7 -7 - —9— Division #7 -7
L - - '0 - - Division ”IO-10 : .. - .0. - .. Division #1040
89 0.4 :- 52 0.2 P
r .
0.2 l 0.1
0 0 A 0.051 0.1 0.15 02 ° 7
Ere
(a) 2D problem (b) SD problem I

0.4 ~

 

—a— Division #1 -I
--+-- Division#4-4
—9— Division #7 -7
- - -0- - - Division 1110-10

 
 
 

 

 

  

 

 

(c) SD problem II

Figure 33. Detailed data distribution in selected partitions

100

5.3.2 Analysis of Correlation Coefficient

In the results of evaluating error with two measures at a time, the data variance is
pretty large in many data partitions. Compared with the results of evaluating
error with only one measure at a time, there is no significant improvement for the
variance feature and the error still can not be exactly determined. According to
the statistical theorem, the examination of correlation coefficient between the
independent variables is very important to evaluate the results. For this [5,., =
f((lbn. 003) function, if the correlation coefficient between 1pm and (bus is very high,
it means that they are highly correlated to each other and only one measure is
necessary to predict the error. The mathematical expression of the correlation

coefficient can be expressed by

:xuxu
r = n i=1 n
(Z 3‘12; )(Z xix)

(5.9)

 

Table 16 shows all the computed correlation coefficients, which can be denoted
as r-values. The first row of this table shows the r-values between the error and
each of these measures. As can be seen, the r-values between E“, and (pun/ops
are 0543/0503, 0728/0819, and 0.785/0.816 for 2D problem, SD problem I,
and SD problem II, respectively. Compared with the r-values between the
relative error and other measures, it is very clear that one and (hos do have the

strongest correlation with the relative error. This result further confirms our

101

results in Section 5.2, that is, one and 093 are the best error indicators among all

measures.

Next, we examine whether $011 and (1)03 are independent to each other. As seen
from Table 16, the r-values between one and 003 are 0.474, 0.752, and 0.710 for
the 2D problem, SD problem I, and SD problem II, respectively. Compared with
the r-values between other measures, they are the most highly correlated
measures except the correlation between the (Dora and 0mm. In other word,
though one and 093 describe different features of local quality, they are not
independent to each other and only one measure is needed. Perhaps, that is the
reason that the magnitude of variance in E“, = f(¢0n. (bus) is similar to that in E“, =
f(¢oa). If that is the case, the measures used to correlate error should be as

independent as possible to each other.

To determine a set of independent measures, the r-values in Table 16 are again
examined. First, since one. has the strongest correlation with error, it should be
one of the measures considered. By examining the r-values between one and all
other measures, we find that than and (born have the least correlation with (Ice in
all three test problems. The r-value between (pea and (ppm is also very low. On
the other hand, ¢nvo has very strong correlation with (born and 033 can also be
considered correlated with 009. Thus, we can select $09. (lion. and (born as a set

of independent measures to correlate with error. In next section, we build the

102

function E... = f(¢0n, (ban. 00m) and check if the data variance in data partition gets

reduced as we expected.

Table 16. Correlation coefficient between the relative error and new measures

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

(a) 2D problem
Ere one we ¢orn (Dave 00;. ¢S§
Ere 1 .000 0.543 0.278 0.309 0.387 0.503 0.407
(Don - 1 .000 0.217 0.316 0.463 0.474 0.467
¢GR ' - 1 .000 0.095 0.202 0.352 0.147
¢DTR ' ' - 1 .000 0.652 0.297 0.496
(ano ' ' - - 1 .000 0.482 0.331
(1)08 ' " ' - - 1 .000 0.340
¢ss ' ‘ - - - - 1.000
(b) SD problem I
Ere (1)011 (Don (Pom (ano 111th these
Era 1 .000 0.728 0.356 0.512 0.579 0.819 0.536
(pop. - 1 .000 0.339 0.581 0.620 0.752 0.625
069 - - 1.000 0.249 0.321 0.474 0.272
(Porn ' ' - 1 .000 0.856 0.609 0.632
¢Rvo ' ' - - 1 .000 0.696 0.570
lbs ' - - - - 1.000 0.565
¢ss ' ' - - - - 1.000
(0) 3D problem ll
Ere (boa ¢GR ¢orn (Dave <l>os ¢ss
Ere 1 .000 0.785 0.487 0.366 0.530 0.816 0.488
(pm - 1 .000 0.446 0.429 0.563 0.710 0.564
(1)99 - - 1 .000 0.285 0.467 0.620 0.263
(ppm - - - 1 .000 0.813 0.500 0.492
¢Rvo " ‘ - - 1 .000 0.672 0.459
files " ' - - - 1 .000 0.455
058 ' ' ' - - - 1.000

 

 

 

 

 

 

 

 

 

103

 

 

 

5.3.3 Evaluation of Errors with Three Measures at a Time

In this section, three independent measures are considered and E... = f(¢on. 063,
com) is examined by using Eq. (5.4) and see if the large variance in each
partition can be reduced. Each of the three measures is divided into ten intervals
by using the equal population approach and the total number of partitions is
equal to 103. Then the linear least square method is used to represent error level

in each partition.

Figure 34 shows the detailed distribution in selected partitions for all cases. In
this figure, each curve represents the distribution of error value in one partition,
and the partition level (1-1-1), (4-4-4), (7-7-7), and (10-10-10) are selected. (1-1-
1) represents the partition with ideal (boa. than. and 00m value while (10-10-10)
represents the partition with the worst 003, (Don. and (born value. The X-axis
represents the actual value of Er.a and the Y-axis represents the percentage of

cells that fall into certain ranges of error values in that partition.

As seen from the results of evaluating error with three independent measures at
a time, the data variance is still large in many data partitions. Compared with the
results of evaluating error with only one measure at a time or with than and 005 at

a time, there is no significant improvement for the variance feature.

104

0.4 .-

 

 

 

 

0.8 - i _
‘ 1
r l r
t “ —e— Division #1 -1 -1 t —e— Division #1 -1 -1
0,5- 1 ---A--- Division#4-4-4 0.3- ---or-- Division#4-4-4
“ ——v-— Division #7 -7 -7 ; —v—— Division in -7 -7
i - - o - - Division #1010 -10 . - - -°- - - Division #1040 ~10

 

 

 

 

580.2 -

1 fl fi‘Tfi Y

 

 

   

 

 

 

 

Ere
(a) 2D problem (b) 3D problem I
0.4 —
i —a—— Division #1 -1 -1
0.31 --.n.-- Division#4-4-4
' —v— Division #7 -7 -7
; - - -o- — - Division #1040 -10
as 0.2

0.1

 

 

(0) SD problem II

Figure 34. Detailed data distribution in selected partitions

105

5.4 Summary

In this chapter, a set of grid-quality measures that can account for both local grid
features and local solution information are evaluated using statistical analysis
based on three sub-databases. The evaluation results show that these
measures have certain correlation with the solution error and perform much
better than those traditional measures that only include geometry information and

almost have no clear correlation with actual grid-induced error.

Also, the generality of correlation curves between solution-based measures and
error is checked based on different 2D/SD simulations using different imperfect
meshes. Among the six solution-based measures, lion, which represents the grid
resolution for velocity direction change, not only shows the strongest relationship
with solution error but also has the best generality feature. The directional
smoothness, 093, also correlates well with error and shows good generality in
different simulations. The streamline smoothness, 033, is similar to (bus but not so
good as ops. Measure (Ion. 00m, and ¢nvo describe the local flow gradient

information and only shows very weak correlation with solution error.

However, the correlation between error and measures is still not good enough to
make error estimation. As described, in the evaluations of error with one
measure at a time, though the mean value of these new measures in each

partition shows certain monotonic and universal correlation with error, the data

106

variance in some partitions is still very large. It means that the local error can not
be accurately determined by single measure. By checking the detailed data
distribution in these partitions, it can be noted that in the regions with high
relative error, many cells still have very good measure values. That is why the

data variance in these partitions is so large.

One possibility is that the local error may be determined by multiple measures
instead of by only one measure. Thus, the function E... = 1(01, 02, ...) was then
established to analyze the correlation between error and multiple measures.
Since 0m and fins are the measures correlating best with local error, the function
E... = f(¢on. $05) is first evaluated based on statistical methods. Though the
strong correlation can be detected between error and measures 00381093 by
using the linear least-square and mean value in each data partition, the data
variance in some partitions is still very large. It indicates that some other factors
exert influences on the local error besides the local measure then and ¢0s- By
using the r-values check, we find that con and (1)03 are not independent variables,
which means only one measure is necessary to predict the error. So, a set of
independent measures -- (Don. (1)011. and 00m - were selected to evaluate the
function Ere = f(opn, than. 00m). However, the data variance in some partitions is
still very large. This means that the error cannot be accurately predicted even
multiple local measures are used. There may have other factors affecting the

error other than the local grid-quality measures.

107

In additional, by summarizing the physical meanings of these solution-based
measures, we find that (Don gives that the local mesh resolution should be fine
enough to resolve the complex flow structures, like recirculating zones, flow
separation, and stagnation flow, etc. ¢Ds shows that the smoothness of local
effective spacing should not exceed certain ratio. Besides the resolution
requirement for regions with sharp change in flow direction and the smoothness
requirement for complex flow situation, we all have experiences that the local cell
resolution is extremely critical if the local flow gradient is high. In order to find out
such a measure, which can account for the local gradient information and has
strong correlation with the local error, we developed several gradient-related
measures, like than, 0079, and (have. But none of them shows strong correlation
with the local relative error. Many other evaluation ways were tried without any

success.

All these numerical results and analysis show that the correlation between the
grid quality and grid-induced error is extremely complicated. In other word,
besides the local grid-quality measures, there do have some other factors that
can affect the error level. At least, some gradient-related measures should exert
influence on error in some way that we still did not figure out. In Part II of this
dissertation, we will further study the possible correlation between error and grid-
quality measures by applying these measures in practical problems and using

grid adaptation approach.

108

PART II
Application of Grid-Quality Measures as Error lndictors and for

Grid-Refinement

109

Chapter 6

Application of Grid-Quality Measures as Error Indicators

6.1 Introduction

As noted in Part I, we developed a set of solution-based grid-quality measures
and then evaluated these measures by interrogating a database, which was
made of three sub-databases. The evaluation showed that when error is
correlated with one measure at a time, (lion 81 1003 have very strong correlation
with solution error but the variance in data partitions is very large, and other
measures only have weak correlation with solution error. When error is
correlated with multiple independent measures, the variance in data partitions
still cannot be reduced. All these results mean that there may have other factors

that affect the errors besides the local grid-quality measures.

Thus, to further study the correlation between error and grid-quality measures, it

is very necessary to make clear following basic questions:

1. How does error pattern change when grid quality gets improved?

2. Since some measures, such as con and 003, already show strong correlation
with solution error, how good they are when they are used in practice?

3. What are the other factors that can influence the error besides the local grid

quality?

4. Why gradient-related measures do not have strong correlation with error?

110

In this and following chapters, we give possible answers to above questions.
Investigations are committed by applying new measures to practical problems

and using grid adaptation approach.

6.2 The Error Pattern of Previous 20 Test Problem

Obviously, the best and most direct way to answer the first question in Section
6.1 is to review the previous 2D test problem and compute the detailed error
pattern on the series of uniform meshes. Then, the change in error pattern
obtained on grids with different quality can be checked. We want to know that
whether the error level can get successively reduced when the grid quality gets
improved? If no, it means that other factors may exert more influence on the
error than the grid quality even when same physical model and numerical
schemes are utilized. If that is the case, probably it makes no sense to set up
the correlation between error and grid quality. Additionally, we like to check how

error region changes when the uniform grid spacing decreases.

0 Average error level check

In Chapter 3, a set of uniform grids is generated to compute the grid-independent

solution for the 2D test problem. The solution obtained on the level 8 was taken

to be the grid-independent solution. The average error level on each coarse grid

111

can be computed according to this grid-independent solution. From Table 17
and Fig. 35, it is clear that the average relative error converges monotonically to
zero as the total cell number increases. The average absolute errors are also
calculated and have the same trend as that of the relative error. This result
indicates that the grid quality is the most critical factor to determine the error level

when same physical model and numerical schemes are used.

Table 17. Average error level in uniform meshes

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Levels Absolute error (m/s) Relative Error (%)
1 3.268x10“ 17.21
2 2.066x104 10.98
3 9.576x10'5 4.85
4 5.755x10'5 2.68
5 2.626x10’5 1.1 1
6 1.170x10'5 0.511
7 4.175x10'6 0.224

0.2 L-
0.151
"2 0.1
0.05:
00 i ‘ ' 10000 i ‘ 20000 1 30000 40000

Cell Number

Figure 35. Average error level in uniform meshes

112

0 Detailed Error Pattern Check

Figure 36 shows the exact relative error pattern on each coarse mesh. In these
plots, the relative error is represented by the gray scale contour and the darker
regions indicate higher errors. In the plot of level1, the regions with high error
mainly locate along the side of the main flow and at the two upper corners.
When the mesh is refined to level2, the errors still appear at the same regions as
plot of level1, but the high error area and error level reduces a little bit. This
trend keeps going on when the mesh is successively refined from levelS to

level7. At level7, the relative error cannot be easily detected anymore.

According to above results, the answer to the first question in Section 6.1 can be

summarized as follows:

1. The average error level can be successively reduced as the grid quality gets
improved. So, the grid quality is truly dominant in determining the relative
error when same model and numerical schemes are used.

2. The errors always appear in some fixed regions with certain flow
characteristics instead of showing up randomly.

Of course, this summary is made only from this 2D test problem. It still needs to

be verified in other flow cases.

113

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Ere
03
0.7
0.6 ..,,... (a) Level1
05
0.4 1
0.3 1
0.3
02
0.1 ”g
°°° .. (b) Level2
. ..., ..,,..-” . ._ j
(c) LevelS
l
f
(d) Level4
-.I F .- . _ m- war—— 1
(e) Level5
‘ . ”'" l
‘l
(f) Level6
.., l
7.1
(g) Level7
‘7 l

 

 

Figure 36. Error pattern on uniform meshes of 2D case

114

6.3 Application of Grid-Quality Measures as Error Indicator

In this section, the goal is to answer the second question in Section 6.1, i.e.,
since these newly developed measures have been proved to have strong
correlation with solution error based on statistical analysis, how good they are

when they are applied as error indicators in practical problems?

As we know, among these measures, (boa has the best performance and thus is
also pretty good. For measures 093, (ion. 00m and thaw. they only have certain

correlation with relative error and behave not so well as 4an and 005.

In following subsections, two 20 problems are selected to assess the quality of
these measures in practical use. One is a 2D moving wall cavity flow and the
other is a 20 chamber flow with one inlet and one exit. First, for each problem,
the grid-independent solution is generated using the same approach as we
introduced in Chapter 4. Second, some coarse meshes and their consequent
‘imperfect’ solutions are selected as examples. Based on the selected mesh and
solution, measures can be calculated on each cell and used to indicate the
possible high error regions. Third, the exact error pattern on these meshes can
be calculated according to the difference between the baseline solution and the
‘imperfect solution’. By comparing the real error pattern with the predicted error

patterns, we can easily assess whether these measures are good error indicators

115

for the relative error. The performance of these measures will also be discussed

and compared with the previous results got from statistical analysis.

6.3.1 2D Moving Wall Cavity Flow

6.3.1.1 Problem Description and Formulate

Figure 37 is a schematic diagram of this 20 moving wall cavity flow. The cavity
has an aspect ratio of 2 to 1 (L = 200 mm and h = 100 mm). The fluid in the
cavity is water with density p = 998.2. Kg/ms and dynamic viscosity p. = 1.003E-3
kg/m-s. The velocity of moving wall is V = 0.01 m/s so that the flow remains

laminar.
V

L” v

L

 

Figure 37. Schematic of cavity-driven problem studied.
L = 200 mm, h = 100 mm, V = 0.01m/s.

This problem, though seemingly simple, was selected because there exist many
complex flow structures in the flow field such as flow separation, recirculation

zones, and high shear layers. Also, we are only interested in the steady-state

116

solution of this incompressible laminar flow and solution to the governing

equations was obtained by using the Fluent code.

6.3.1.2 Grid-Independent Solution Generation

By using the same method described in Section 4.2.2, we obtained the grid-
independent solution for this 2D cavity flow problem. Table 18 shows the series
of uniform mesh. A total of 9 grids/solutions were generated/computed and the
final level contains 32,768 cells. Table 19 lists the average relative error and

absolute error between successive levels. Figure 38 shows the convergence

rate when grid gets refined. As can be seen, E: in Table 19 is well below 0.5%.

Thus, in this study, the solution on level 9 grid is considered as the grid-

independent solution.

Table 18. Uniform grids used to get grid-independent solution

 

 

 

 

 

 

 

 

 

 

Level Grid Spacing (m) Total Cell No.

1 6.250e-3 (1/160) 512

2 4.1676-3 (1/240) 1,152

3 3.125e-3 (1/320) 2,048

4 20836-3 (1/480L 4,608

5 1 .563e-3 (1 /640) 8,192

6 1 .2506-3 (1 /800) 12,800

7 1 .042e-3 (1/960) 18,432

8 8.929e-4 (1/1 120) 25,088

9 78136-4 (1/1280) 32,768

 

 

 

 

 

 

117

Table 19. Convergence process to get grid-independent solution

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Levels En (m/s) E... (%)
1-2 2.668x10" 17.32
23 1.308x10" 8.86
3-4 9.482x10'5 6.19
4-5 3.448x10'5 2.30
5-6 1 .603x10'5 1 .10
6-7 8.978x10j 0.650
7-8 5.642x10'6 0.438
8-9 3.831x10'6 0.309

0.2

0.15

0.05 P

 

 

 

1 1 1 1 l 1 1 1
0 10000 20000 30000
Call number

Figure 38. Convergence process to get grid-independent solution

Figure 39 shows the streamline pattern of the grid-independent solution. The
direction of driven wall is from left to right with velocity magnitude equals 0.01
m/s. There is one strong clockwise vertex at the right side of cavity and one
weak anti-clockwise vertex at the left part of cavity. Some other complex flow
structures, such as flow separation and small circulation flow, exist at the left side
of wall and corners. Thus, this 2D cavity flow does posses a number of

complicate flow characteristics.

118

 

 

 

Figure 39. Flow pattern of grid-independent solution,
obtained by using a mesh with 32,768 aspect ratio unity cells.

6.3.1.3 Error Generation Check

0 Average error level check

Similarly as before, the grid-independent solution is projected on each of the
coarse uniform mesh. Then the average error level is computed according to the
difference between the ‘imperfect' solution on coarse mesh and final baseline
solution. This result is given in Table 20 and Figure 40. As seen from the plot,
with the increasing of cell number, the average error level drops down
monotonically. At the beginning levels, the error reduces very quickly as mesh
gets refined. When the grid spacing continues to decrease, the error level

reduces very slowly and becomes convergent to zero.

119

Table 20. Convergence process to get grid-independent solution

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Leve's 75.. (m/SI E. (%)

1 5.352x10" 28.82
2 2.855x10" 16.31
3 1.643x10‘4 9.76
4 6.855x10‘5 4.16
5 3.434x10'5 2.08
6 1.727x1075 1.11

7 8.703x10'6 0.613
8 3.831x10“6 0.309

03

0.25 :
0.2 :

g E
u‘ 0.15:

0.1 I

 

0.05 :

 

 

1 1 1 1 l 1 1 1
0 10000 20000 30000
Call number

Figure 40. Convergence process of the 20 moving wall cavity flow

0 Detailed Error Pattern Check

The exact error pattern of each coarse level is shown in Fig. 41. As observed,
the high errors always happen in certain regions. When the global grid quality
becomes better and better, the area of the high error region becomes smaller
and smaller. When the grid is fine enough, the error can hardly be detected.
This result is consistent with the summary made from the previous 20 test

problem.

120

 

 

 

 

 

(c) Leve|3 (d) Level4

 

 

 

 

   

 

 

 

 

 

 

 

 

 

 

(e) Level5 (f) Level6
.l i
(9) Level? (h) Level8

Figure 41. Error pattern on uniform meshes of 2D cavity case

121

6.3.1.4 Example #1: Uniform Coarse Grid

As shown in Fig. 42, a uniform coarse mesh with 24x48 cells is taken as the first
example test. An “imperfect" solution is then computed based on this coarse

mesh.

 

Figure 42. Uniform example mesh

Suppose that the grid-independent solution is not available for this problem now.
The newly developed measures will be used to indicate the possible high error
regions. In order to compute measure patterns on this example, the geometry
information needed includes the x and y coordinates of each node and the
connectivity information of each cell. The solution information needed Includes
the x and y velocity components at each cell, which are used to compute 41m,
003, 035, and (bravo. The first and second derivatives of velocity are also needed

and used to calculate lien and 09m.

From the statistical analysis in Chapter 5, (lion and 1|le show the strongest
relationship with the exact solution error. Their correlation with error is not only

monotonic but also very universal in different simulations. Thus, 00R and 003 are

122

taken as the main error indicators for this problem. 035, ion. 410m, and ¢nvo are
only used for comparison and reference. Figure 43 shows the measure patterns,
which are represented by the gray scale contour. The regions covered by dark
gray mean that those cells have bad measure value and thus high relative error

is predicted.

Figure 43 (a) presents the result of (ion with range from 0 to 120. As seen from
the plot, the region under the left side of moving wall indicates the largest area
with bad 0m. By comparing with the flow pattern, it can be seen that regions with
high con value include the rotating centers of the two large vortical structures, the
middle of bottom wall where the two large vortical structures meet, and the right
bottom corner. Some other regions, like the left and the right upper comers, also
are indicated with bad (ion value. The (jigs pattern is shown as Fig. 43 (b). The
regions with bad tics is a little similar to that of ion. which include vortex centers,
stagnation flow, and corners. Figure 43 (0) shows the result of 035, which is very
similar to the (ins pattern. Figure 43 (d), (e) and (f) show the $06. 00m, and 0WD
patterns. Compared with (lion and $08 pattern, they look much different. The
regions with high (lien include the corners, the zone where the two rotating flows
meet, and some positions surrounding the right vortex. The regions with high

00m and 0mm mainly locate on the boundary walls and the two rotating centers.

123

 

(c) 033 Pattern

   

(6) 0073 Pattern (f) 0WD Pattern

Figure 43. Measure patterns for example #1

Some traditional measures, such as cell volume, cell aspect ratio, and cell
smoothness, are popularly used as criteria to check the grid quality in practical
problems. For this example case, the cell volume is constant and the aspect
ratio is equal to one everywhere because all the cells have uniform size. As for

cell smoothness, which is defined as the maximum ratio between the local cell

124

and all its neighbors, it is also equal to one for all cells in this case. 80,
according to the check of cell aspect ratio and cell smoothness, which only
include geometry information, the grid quality in this case is absolutely perfect.
As for cell volume, it also cannot provide us any information for error prediction.

Thus, all these traditional measures seem not useful for this problem.

a The real grid-induced error pattern

Figure 44 shows the exact relative error pattern over the flow field. It is
computed based on the grid-independent solution. The error level is expressed
as the gray scale contour ranging from 0 to 80%. Regions with high relative error
are labeled in figure for convenience of discussion on their location and cause.
Region #1, which locates under the left side of moving wall, is the largest error
area showing up in the flow field. By comparing with Figure 39, it can be seen
that region #1 contains a strong recirculating flow that impinges on the moving
wall. In region #2, the high error is due to the stagnation flow where the flow
separates into upward and downward direction. Region #3 and #7 are around
the centers of two large rotating flows, respectively. Region #5 locates at where
the two rotating flows meet. Region #4 and #6 are located at the left and the

right corners, respectively, and contain small circulation flows.

125

#1 #7

#4

 

Figure 44. Exact relative error pattern for example#1

By comparing this real error pattern with the measure patterns in Fig. 43, it can
be seen that for can, the largest region with bad ¢Dn value exactly matches the
region #1 in the actual error pattern. The positions of actual error region #2, #3,
#5, #6, and #7 are also predicted accurately on the (ion pattern. The only region
that cannot be indicated is region #2, which is located on the left comer of cavity.
So, almost all the major high error regions are successfully indicated by the
measure ¢on- As for the prediction by ¢os1 it can accurately indicate the region
#6 in the actual error pattern. Region #1 is only partly indicated by (1)03. Region
#2, #3, #4, and #5 are also predicted, but not so accurate as the prediction by
¢0n- There are also some areas that are indicated bad 003 value but the actual
error is very low. The behavior of 055 is very similar to that of 1105- For the
estimation based on (baa, 00m, and ¢nv0. they are not good error indicators.
Many regions with bad values do not exactly match the regions with high relative

error.

126

6.3.1.5 Example #2: Non-Uniform Coarse Grid

Figure 45 is another example to test the performance of measures as error
indicators. It is a non-uniform 32x64 grid, with the vertical smoothness ratio
equal to 1.04 for double sides and the horizontal ratio equal to 1.08 for double

sides.

 

Figure 45. Non-uniform examples grid

Similarly, the measure patterns of 00a, tics. 053, 11166. 00m, and (bravo are computed
and shown in Fig. 46. By comparing Fig. 46 (a) with Fig. 43 (a), the regions
indicated with bad ¢DR value are very similar. The only differences are that the
region in the left bottom corner of this example is indicated with high 113011- For
other measures, their patterns are also very similar to the measure patterns in

Fig. 43, but the predicted value level is a little bit lower.

127

 

 

(6) 1ion; Pattern (f) ¢fivo Pattern

Figure 46. Measure patterns for example #2

The traditional geometry-only measures are also calculated for discussion. For
cell volume, those cells near comers have very small size and those cells in the
center of cavity have big value. Smoothness is equal to 1.08 everywhere. The
pattern of cell aspect ratio is shown in Fig. 47, which is perfect along two

diagonals and has high value in cells far from diagonals.

128

 

Figure 47. Cell aspect ratio pattern for example #2

o The real grid-induced error pattern

Figure 48 is the exact error pattern on this non-uniform coarse mesh. This grid is
much different from the example #1 grid. First, they have different cell resolution.
Second, this grid is non-uniform. By comparing Fig. 48 with Fig. 44, which is the
exact error of example #1, it can be noted that the regions indicated with high
relative error in both figures look very similar. Region #1, #2, #4, #5, and #6 in
both plots show the highest relative error. However, the area of these regions is
a little smaller than that in previous uniform example. This means that the error
probably shows up at certain regions and the global error can be reduced when

the global grid quality gets improved.

When comparing the real error pattern with the measure patterns for this
example, it can be observed that then again successfully indicates the actual high
error region #1, #2, #3, #5, #6, and #7. It also indicates the region #4, which

contains a small vertex. Similarly, <l>0s can also predict some high error regions

129

but not so accurate as ¢DR~ 035 performs very similar to clips. The prediction by
(pen, 00m, and (have are the worst among these new measures. As for the
traditional ones, like the smoothness, cell volume, and aspect ratio, still cannot
provide any useful information for error estimates.

#1 #7

 

 

Figure 48. Exact relative error pattern for example#2

As can be seen, similar results are obtained in above two examples, which uses
uniform and non-uniform mesh, respectively. Basically, eon, (Dos. and 053 show
their advantages over (lien, 09m, and 0mm as well as the traditional measures.
Furthermore, 00a is so far the best error indicator to indicate highest relative
errors. 005 and 053 are highly correlated with one. and work also well. than, 00m,
and 0WD, which are originally developed to account for the high gradient
information, do not show any strong relationship with the relative error. These
results are consistent with the previous results we got from the statistical

analysis.

130

6.3.2 20 Chamber Flow

In this subsection, another 2D problem is selected to further assess the quality of
these measures as error indicators. The performance of these measures will be

discussed and compared.

6.3.2.1 Problem Description and Formulate

Figure 49 is a schematic diagram of this 20 chamber flow. As seen from the
figure, the main chamber is a 125mm x 100mm rectangle. One inlet and one exit
locate on the top of chamber with duct width equal to 25mm and 12.5mm,
respectively. The fluid in the cavity is water with density p = 998.2 Kg/m3 and
dynamic viscosity p = 1.003E-3 kg/m-s. The velocity at inlet is uniform at V =

0.01 m/S.

This problem was selected because many complex flow structures exist in the
flow field. Also, we are only interested in the steady-state solution of this
incompressible laminar flow and solution to the governing equations was

obtained by using the Fluent code.

131

Uniform
Velocity
V Inlet

L“ III II

 

 

Figure 49. Schematic of 2D chamber flow problem

6.3.2.2 Grid-Independent Solution Generation

Similarly, the grid-independent solution was obtained by performing a series of
simulations, all of which used meshes with uniformly distributed grid points with
aspect ratio of unity cells so as to minimize the grid-induced errors. Table 21
summarizes the grid spacing and total number of cells used. A total of 9

grids/solutions were generated/computed and the final level contains 23,552

132

cells. Table 22 lists the average relative error and absolute error between

successive levels. Figure 50 shows the convergence rate when grid gets refined.
As can be seen, E: in Table 22 is well below 0.5%. Thus, in this study, the

solution on level 9 grid is considered as the grid-independent solution.

Table 21. Uniform grids used to get grid-independent solution

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Coarse Grid Spacing (m) Total Cell no.
Level

1 62506-3 (1/160) 368

2 4.167e-3 (1/240) 828

3 3.125e-3 (1 /320) 1 ,472
4 20836-3 (1 I480) 3,312
5 1.563e-3 (1/640) 5,888
6 12506-3 (1/800L 9,200
7 1 .042e-3 (1/960) 13,248
8 89296-4 (1/1 120) 18,032
9 7.8136-4(1/1280) 23,552

 

Table 22. Average relative errors between successive levels

 

 

 

 

 

 

 

 

 

 

Levels Absolute error Relative error

(W8) (‘70
1-2 3.695x10‘ 23.69
23 1 .672x104 7.44
3-4 1.702x10" 8.14
4-5 6.403x10'5 3.23
5-6 3.137x10'5 1.66
6-7 1 .932x10'5 0.963
7-8 1 .369x10'5 0.624
8-9 1.024x10'5 0.431

 

 

 

 

133

Q4-

0.3 L

 

 

0‘ 10000 20000 30000
C011 Numbor

Figure 50. Average error level in uniform meshes

Figure 51 shows the flow pattern of this grid-independent solution. On the top
part, one main flow with high velocity goes from the inlet to the exit duct. There
is one large clockwise recirculation flow created in the middle of rectangular
chamber. This flow merges into the main flow at the left upper corner. Then, the
two flows separate at the right upper corner. Several small size vortex flows are
induced at the middle-top of chamber and the two bottom comers. Thus, this 2D

flow posses many complicate flow characteristics.

 

Figure 51. Flow pattern of grid-independent solution

134

6.3.2.3 Error Generation Check

0 Average error level check

Similarly, the average error level of each coarse mesh is computed according to
the grid-independent solution. This result is shown in Table 23 and Fig. 52. It
can be seen that the average error level drops down monotonically as the mesh

gets refined.

Table 23. Average relative errors on each coarse mesh

 

 

 

 

 

 

 

 

 

Levels Absolute error (m/s) Relative Error (%)
1 8.434x10“ 37.10
2 5.138x10‘r 20.09
3 3.393x10" 13.85
4 1.514x10‘ 6.14
5 8.502x10'5 3.17
6 4.183x10'5 1.64
7 2.127x10'r 0.864
8 1.024x10'5 0.431

 

 

 

 

 

 

 

 

1 1 l 1
10000 20000 30000
C011 Number

Figure 52. Average error level on coarse meshes

135

. Detailed error pattern check

The detailed error pattern of each coarse level is also computed to check how
these error distributions change according to the increasing of cell number.
These patterns are shown in Fig. 53. Same as before, the error always happens
in certain regions instead of showing up randomly. When the grid gets refined
step by step, the global errors decrease little by little. When the grid spacing is

fine enough, the error can be finally removed.

So far we have examined three sets of error patterns, including one 20 test
problem and two 20 example problems. In all these problems, the highest
relative errors always show up at regions with complex flow structures, like the
strong reCirculation flow, stagnation flow, corners with small vortex, and zones
where two flows merge or separate, etc. It means that the local grid quality in
regions with complex flow situations is much more important than the grid quality
in other regions with benign flow situations, and once it is not good enough, high

relative error will easily be activated.

136

 

    

1 :76..— s.—

 

(0) Levels

   

 

 

 

 

 

(e) Level5 (f) Level6

 

 

Figure 53. Error pattern on
uniform meshes of 2D
chamber case

 

 

 

 

 

 

(9) Level? (h) Level8

137

6.3.2.4 Example #1: Uniform Coarse Grid

As shown in Fig. 54, a uniform coarse mesh with A = 1/320 and totally 1,472 cells
was created. Then, the preliminary solution is computed based on this coarse

mesh.

 

Figure 54. Uniform mesh for example #1

Same as before, measures are computed on each cell by using both the
geometry and solution information. Figure 55 (a) shows the pattern of 1106- As
seen from the plot, the regions on top of the chamber and at the two bottom
corners are indicated with bad (ion value, which predict the possible high relative
error. Comparing with the grid-independent flow pattern, all these regions have
small recirculation flows. Other high (lion regions include the center of the large
vertex, the zones where the main flow and the large vertex flow merge and

separate. So, almost all the zones with complex flow structures have high than

138

value. Figure 55 (b) and (c) present the (lbs and 053 patterns, respectively. It is

clear that both 003 and 035 patterns are very similar to that of (ion-

Figure 55 (d), (e), and (f) show the ¢GR1 lime. and (have patterns, respectively. As
noted, the regions with high (ion value mainly include the complex flow structures
at bottom comers and the left-up part of chamber where the two flows merge
together. Some regions inside the inlet and exit ducts are also indicated with
high (pen value. As for them and 0mm pattern, the highlighted regions mainly
concentrate on the boundary walls. By comparing them with the con and <l>0s
patterns, we can see that the patterns of than, 0973, and ¢nvo look much different,

especially for cells in the inlet and exit duct.

139

 

 

 

 

 

 

 

(e) 11073 (f) ¢RVD

Figure 55. Measure patterns for example #1

o The Exact Error Pattern

Figure 56 shows the exact relative error distribution over the flow field. Same as
before, regions with high relative error are labeled for discussion on their location
and cause. Regions #1, #4, and #5 are the zones containing small vortical
structures. Region #3 locates at the center of the large rotation flow. Region #2
indicates where the main flow and the large rotation flow merge together. In

region #6, the large rotation flow separates from the main flow.

Comparing with the exact error pattern, ¢DR successfully indicates all the labeled
regions. $03 and 035 are similar to ion and also can predict the main high error
regions. As for 4km. 09m, and (have, they can only partly estimate the error
regions, like region #2, #4, #5, and many regions with bad 09m and own values

do not show any high error.

 

 

Figure 56. Exact error pattern of example #1 mesh

141

6.3.2.5 Example #2: Non-Uniform Coarse Grid

0 Mesh and solution of example case

Figure 57 is another example to test the quality of new measures for error
indication. It is a non-uniform grid with totally 3,168 cells. Compared with the
previous uniform mesh, grid nodes In this mesh are clustered on the boundary

walls.

 

Figure 57. Non-uniform example #2 mesh

Similarly, the patterns of (jinn, tips, 035, 063, $0111. and (bravo are computed and
shown in Fig. 58. By comparing Fig. 58 with Fig. 55, the regions indicated with
bad measure value are very similar to each other. The only difference is that the
measure values along boundary walls in this example are much lower than that
of the uniform example. This is because more points are put in the boundary

regions for this non-uniform mesh.

142

(e) ¢DTR

 

 

 

 

 

. ltlulllmullmlililuu

 

 

(f) ¢Rvo

Figure 58. Measure patterns for example #2

143

 

 

 

 

 

 

o The Exact Error Pattern

Figure 59 is the exact error pattern on this non-uniform coarse mesh. By
comparing Fig. 59 with Fig. 56, it can be noted that the regions indicated with
high relative error in both figures look very similar. This again confirms that the

relative error always happens at certain regions no matter how grid changes.

When comparing the actual error pattern with the measure patterns of this
example, it can be observed that fine again successfully indicates almost all the
actual high error regions. 003 can also predict some high error region but not so
good as ice. and 033 performs very similar to 093. The prediction by (lien, 00m,

and <l>nvo are still the worst in these new measures.

 

 

 

 

 

Figure 59. Exact error patterns for example #2

6.4 Summary

In this chapter, the correlation between error and grid-quality measures is further
studied by applying new measures to practical problems. In order to answer the
first basic question: How does error pattern change when grid quality gets
improved? Three sets of error patterns, including the 2D test problem and two
20 example problems, are examined. In all these problems, we can find that the
grid quality is truly the dominant factor and the global error can be monotonically
reduced with the improvement of the global grid quality. Furthermore, the
highest relative errors always show up at regions with complex flow structures.
This means that the solution is more sensitive to local grid quality in regions with

complex flow situations than regions with benign flow situations.

Next, in order to answer the second question: Since some measures, such as
(Don and <l>0s. already show strong correlation with solution error, how good they
are when they are used in practice? All measures are computed on the cells on
the two 2D example problems and then compared with the real errors, which are
calculated according to their own grid-independent solution. The results show
that (lion is the best error indicator. It can successfully indicate almost all the
actual high error regions. tips can also predict some high error region but not so
good as (lion, and 053 performs very similar to tips. The prediction by (pen, 00m,
and ¢Rvo are still the worst in these new measures. These results are consistent

with our previous statistical analysis.

145

Chapter 7
Differentiating Between Source and Location of Error for

Solution-Adaptive Mesh Refinement

7.1 Introduction

In chapter 6, a set of numerical experiments were performed and clearly showed
that the measures 00R, tins, and 053 are good error indicators in practical use.
But the gradient-related measures, like (ion. 00m, and envo, still do not show
strong correlation with errors. In this chapter, we will continue to study the
correlation between the local grid quality and local error by using grid adaptation
approach, and try to solve the third and fourth question in section 6.1: What are
the other factors that can influence the error besides the local grid quality? Why

gradient-related measures do not have strong correlation with error?

According to the previous statistical analysis and numerical experiments, some
measures have proved to be very good error indicators, so it is natural to use
them to refine the local grid and consequently reduce the errors. This is the
basic idea of grid-adaptation. Generally, two issues should be concerned in the
adaptation process. One is what kind of adaptation procedure is the best choice.
The other is which error estimator should be used. In next section, the widely
used h-refinement method is adopted to commit adaptation and (lion is selected

as the first adaptation criteria. In the adaptation process, the value of adaptation

146

criteria in each cell is compared with a threshold value. Once the criteria value is
greater than the threshold value, the local cell will be marked for refinement.
Finally, all the marked cells are isotropically subdivided. For example, a

quadrilateral cell will be split into four quadrilaterals like Fig. 60.

 

 

f

1:) it ii

0

 

 

 

 

 

 

 

 

Figure 60. Schematic for isotropic subdivision

For convenience, we directly use the adaptation function in Fluent5 by adding the
new criteria into Fluent5 variable list. Two different approaches were used to
determine the number of cells that get refined in each adaptation stage: fixed-
percentage and fixed criteria. For the fixed-percentage approach, the
percentage was set at 25% of the current number of cells. Thus, if we start with
1000 cells, then the first refinement will allow h-refinement to be applied to the
250 cells with the highest criteria values. Since each h-refinement produces four
cells from one cell in two dimensions, the total cell number after one h-refinement
is 1,750. In the fixed-criteria approach, the user specifies a bound on the error
indicator. All cells having criteria values greater than the specified bound will be
refined. With this approach, the number of cells that get refined becomes smaller

as the refinement process proceeds.

147

7.2 Fix-Percentage Adaptation by con and Exact Solution Error

7.2.1 2D Moving Wall Problem

In this subsection, we will commit adaptations for the 20 moving wall problem by
using different criteria. For all mesh refinements described in this section, the
starting “coarse” mesh and solution are taken from the level 2 in Table 18, which

has a grid spacing of 4.167x10’3 m and 1,152 cells.

0 Adaptation committed by ton

Figure 61 (0) shows the starting coarse grid and its consequent exact relative
error pattern with respect to its grid-independent solution. As listed in Table 24,
the average error of starting mesh is 16.31%. 11011 is computed for all cells and
the pattern is also shown in Fig. 61 (0). The legend of fine for all plots is from 0.
~ 120. Then, the computed (lion data are exported into Fluent5 and used for the

lso—Value adaptation.

Table 24. Fixed-percentage adaptation process by using (lion

 

 

 

 

 

 

 

 

 

 

 

 

 

Refined (pm Marked Total Ear, Ere
Stage (%) Cells Cells (m/s) (%)
Original 1,152 2.855x10'4 16.31
1 14.5 289 2,019 2.228x10" 12.75
2 12.9 506 3,603 1.949x10" 10.97
3 12.0 897 6,456 1.847x10“ 10.91
4 11.8 1,624 11,493 1.815x10" 10.78
5 11.2 2,805 19,992 1 .807x10’4 10.76

 

148

 

Mesh Exact Relative error ¢DR Pattern

(0) Starting case

(1) Stage 1 refinement

 

(2) Stage 2 refinement

 

(3) Stage 3 refinemenm

    

'7 ‘ug—nsfi : "s—.__

5:- ‘
....... . L 315??» a \. 2
(4) Stage 4 refinement

1;; ,- w" z —

2 14..., 0 g

 

() Stage rinement

Figure 61. Adaptation by ¢DR using fixed-percentage method

149

 
  
 

As seen from Table 24 and Fig. 61 (1), the first refinement will allow h-refinement
to be applied to the 289 cells with the highest (>03 (greater than 14.5), and thus
total cell number becomes 2,019. Then, an improved solution is computed on
this refined mesh and the updated error pattern can be generated, which is
shown in Fig. 61 (1). Comparing with the error pattern on the starting grid, the
solution after one refinement gets some improved. As noted, the errors in
regions labeled #2, #3, and #6 are greatly reduced (the region label was defined
in Fig. 44). The errors in regions #1 and #4 drop down a little bit. As for regions
#5 and #7, the local errors do not change. Table 24 also presents that the
average error level reduces to 12.75%. From this result, it indicates that the
adaptation process based on (Don can reduce the solution error but the error
reduction is not so efficient as expected. Of course, only one stage of adaptation
may be not enough for solution improvement, further adaptation should be

committed by using the same criteria.

Similarly, Table 24 shows that 506 cells got marked for the stage 2 and the total
cell number becomes 3,603 after refinement. As seen from the Fig. 61 (2), the
refinements are mostly embedded within previous refinements. This means that
though the regions with bad 4m get improved in stage1 refinement, the local grid
quality in these regions still needs to be further refined. Figure 61 (2) also
presents the consequent error pattern on this refined stage 2 grid. To our

surprise, the actual error pattern does not make big improvement compared with

150

the previous refinement level. Only errors in regions #1 and #4 get a little bit

reduced. No big improvement can be detected for other high error regions.

Furthermore, the same process continues for further adaptation and Table 24
lists all the detailed procedure. Totally five stages of adaptation were committed
and the final cell number is 19,992. According to the exact error patterns for
Stage 3, 4, and 5 in Fig. 61, it can be seen that the solution does not get
improved even further refinements were committed. Table 24 also lists the
average relative errors of each stage. As observed, when the first-stage
refinement is committed, the error level comes down about 3.56% with the total
cell number increasing 25%. In the second-adaptation stage, the error level only
decreases another 1.78%. For further stages, the error level almost keeps at the
same level. it means that further refinement by rpm; can not reduce the error so

efficiently as the first stage of refinement does.

Also, it indicates that, though (Don has proven to be a very good error indicator,
the adaptation process based on than is not so efficient for solution improvement
as expected. One possibility is that (pop. is a good error indicator but it still can not
perfectly match the exact error. In order to confirm this, we compute ¢on pattern
on each stage of refined mesh, which is shown in Fig. 61. As noted, in stage 1,
the (>03 pattern can indicate most of the high error regions, and the area of
regions indicated with bad rpm; is also smaller than that of starting case. This is

consistent with the change of relative error from starting case to stage 1. In

151

stage 2, measure (pun can still predicts the high error regions and the level of «pm
continues to reduce. So far, ¢on performs very well. When the refinement
proceeds, the exact error pattern almost keeps unchanged, but the value of $03
in all regions continues to decrease. This means that although (boa can still

indicate the highest error regions, it is not a perfect quantitative estimator for the

exact relative error under some situations.

80, it is natural for us to take the exact solution error as the ‘perfect’ error
indicator. This is supposed to reduce error most efficiently. If the efficiency of
adaptation using exact error is much higher than using 4m, it means that rpm is
still far away from the ideal error indicator and more efforts are needed to find out
better error indicators. On the contrary, if they refine the similar regions and their
efficiencies for error reduction are almost in the same level, it further confirms

that rpm is truly a good indicator for the exact relative error.

. Adaptation committed by Ere

In this section, we show the results of using the exact relative error as the error
indicator. The goal is to study what “perfect" error indicator can do to reduce
relative error by using as few grid points as possible. By comparing the effects of
using ‘perfect’ error indicator with the results of using rpm, it can further evaluate
the quality of (lion- The exact relative error is computed aocording to the grid-

independent solution and then imported into Fluent5. The same iso-value

152

adaptation approach is used to commit refinement by selecting the exact relative

error from the variable list in Fluent.

Figure 62 shows the refined mesh and its consequent error pattern of each
stage. Similarly as adaptation by 099, the first refinement is applied to the 290
cells with the highest E“a (greater than 23.3%) and total cell number reaches
2,022. By comparing its consequent error pattern with the ‘starting’ error pattern,
it can be seen that the solution truly gets some improved, especially for regions
labeled #1, #2, #3, and #6 (the region label was defined in Fig. 44). When the
second stage of adaptation is committed, the refinement can only reduce the
error in region #1. For further stages, no obvious reduction for errors can be
detected. Table 25 shows the detailed adaptation process by using Ere. For the
average error level, the error drops down 3.19% after the first-stage of refinement
and reduces another 1.81% on the second-stage. For further stages, the error

level does not change anymore.

Table 25. Fixed-percentage adaptation process by using Em

 

 

 

 

 

 

 

 

 

 

 

 

 

Refined (pm Marked Total Eab Em
Stage (%) Cells Cells (m/s) (%)
Original 1,152 2.855x10" 16.31
1 23.3 290 2,022 2.427x10" 13.12
2 22.9 507 3,714 2.127x10‘4 11.31
3 27.9 888 6,546 2.054x10“ 10.78
4 34.9 1,598 1 1,670 2.046x104 10.35
5 49.0 2,818 20,559 2.044x10“ 10.84

 

 

153

Mesh Actual Relative error $012 pattern

 

(0) Starting case

  

(1) Stage 1 refinement

 
   

{
111
“‘1' - .111. men".

 

(2) Stage 2 refinement

 
 

 
  
 

in»
., i
l

 
 
   

(3) Stage 3 refinement

 

(4) Stage 4 refinement

minus.
I

 

(5) Stage 5 refinement

Figure 62. Adaptation by Ere using fixed-percentage method

154

. Comparison of Error Reduction Efficiency between 00;. and Ere

As observed, there are many common features in Fig. 61 and Fig. 62. First, in
both plots, the refinements in each stage are mostly embedded within previous
refined area. Second, the refined regions by using 0m and En, on the same
stage look similar. Third, the efficiencies of error reduction by both indicators are
almost same. This efficiency can be described as how much the error can be
reduced when the same number of cells is used. Table 26 summarizes the error
reduction of each stage by using 00p. and Em. It clearly shows that their error
reduction efficiencies are in the same level. All above comparisons prove that

0011 is a very good indicator for the exact relative error.

Table 26. Comparison of error reduction

 

 

 

 

 

 

Refined Adaptation By (Don Adaptation By Em
Stag:
1 3.56 3.1 9
2 1.78 1.81
3 0.06 0.53
4 0.23 -0.07
5 0.02 0.01

 

 

 

 

 

For the adaptation by using the exact error, which can be considered as the
“perfect” error indicator, its efficiency of error reduction is expected to be the
best. However, similar as that of 0m, its error reduction after the first-stage of
refinement is only 3.19% and the error reductions in further refinements are even
lower. In other word, though the perfect error indicator is used, the refinement

process is still not so efficient as we imagined. Since the exact solution error has

155

been used as error indicator, we can no longer put the blame on the error

indicators.

Perhaps, the problem is that we are not refining in the right place. Perhaps, the
regions requiring refinement are not necessarily located in regions where the
relative error is the highest. On this, we note that when we refine a mesh in a
certain region, the relative error reduces not only just in that region but also in

other regions where mesh spacing have not been reduced.

Let us review the three types of error defined in chapter 3. The first one is the
relative error, which is the most meaningful error and is used by default in all the
examples to indicate the exact error level. The second one is the absolute error.
It is recognized that the absolute error is not very useful since it is problem
dependent. The final one is the change of relative error. It is defined as the
maximum change of relative error between the local cell and all its neighboring
cells. The original idea to define it is to explain regions in which the local grid
quality is good but the local error is high. Such situations occur because the
errors were created elsewhere where grid qualities are poor. So the change in
relative error would reveal information about how the local grid quality affect the
change of solution error in neighboring cells. In the previous statistical analysis,
the change of relative error also shows certain relationship with local grid quality,
though not so strong as the relative error. It indicates that if the local error is high

and the local grid quality is perfect, the source for this error may exist in other

156

regions and transport into the local cell. Perhaps the error source and actual
error region are neighbors, just as the idea for defining the change of relative
error. It is also possible that these two regions are far from each other and only

have certain correlation based on hidden physics.

So, among the three errors, the relative error is meaningful to indicate the actual
error level. The analysis on the change of relative error suggests us that the
local error may be remotely affected by grid quality in other regions. As for the
absolute error, though it is problem dependent and not a meaningful measure of
solution error, it is a measure relating to truncation error in the CFD analysis. If
the error source does exist, perhaps, it has some kind of relationship with the

high truncation errors.

0 Adaptation committed by Em,

First, let us examine the actual absolute error regions on the starting mesh,
which is shown in Fig. 63. The error level is expressed by the gray scale contour
ranging from 0 to 1x10'3m/s. Regions with high absolute error are labeled for
discussion on their location and cause. Region #1 locates just under the left side
of moving wall. Region #2 presents the largest error area in the flow field, which
surrounds the largest vortical structures locating at the right part of cavity.

Comparing Fig. 63 with the exact relative error pattern, which is given by Fig. 44,

157

it can be seen that the regions with high relative error and with high absolute

error have very different locations.

 

 

Figure 63. Exact absolute error pattern for solution on starting mesh

Next, we show the results of using the exact absolute error as the error indicator.
The goal is to study the efficiency of error reduction by this error indicator.
Similarly, the exact absolute error on each cell is computed according to the grid-
independent solution and then imported into Fluent5. Figure 64 shows the
refined mesh and its consequent error pattern of each stage. The first refinement
is applied to the 290 cells with the highest Eab (greater than 4.40x104 m/s) and
total cell number reaches 2,022. Figure 64 (1) presents the exact error pattern
after the first-stage refinement, it is clear that almost all the regions labeled with
high relative error gets much improved. Especially for regions #2, #3, #5, and #6
(the region labels were defined in Fig. 44), the local grid quality keeps exactly
same as before, but the error level reduces greatly. The only explanation is that
these errors are affected remotely by the grid quality in the refined regions. in

other word, the refined region, which has high absolute error, may be the

158

possible global error source region. When further adaptations are committed, the
refinements in each stage are mostly embedded within previous refined area and

the error level keeps decreasing.

Table 27 shows the detailed adaptation process by using Eab. For the average
error level, the error drops down 6.74% after the first-stage of refinement. This is
much more efficient than what 003 and Ere did, which are only 3.56% and 3.19%,
respectively. When the refinement continues, the average error keeps going
down. This is another difference from the results by using 001:. and Era, in which
the error can not be further reduced even the further stages of adaptation were

committed.

Table 27. Fixed-percentage adaptation process by using Ed,

 

 

 

 

 

 

 

 

 

 

 

 

 

Refined 0m Marked Total Eat, Em
Stage (m/s) Cells Cells (m/s) (%)
Original 1,152 2.855x10“ 16.31
1 4.40x10'4 290 2,022 1.631x10‘4 9.57
2 3.33x10'4 504 3,645 1.390x10:‘ 9.50
3 3.63x10'4 336 6,516 1.124x10’4 3.23
4 3.41 x1 0“ 1 ,596 1 1,533 8.809x10'5 6.63
5 3.14x10" 2,322 20,202 7.931x10'5 5.94

 

 

Mesh Actual Relative error than Pattern

   

      
    
    
  

(0) Starting case

(1) Stage 1 refinement

..
1

E111, "r‘ll..mmw

I"

k

(2) Stage 2 refinement

  

"w 7’ Fm— (’

g f:

5._ = ‘ i t_ ‘
............................... . K .. x A
(3) Stage 3 refinement

W ~ gnaw—a, ’

"3- ’ 3 5 Q
~ 3
’ ..—

F is: a. R. .e‘. E ‘5 J
(4) Stage 4 refinement

: ' E

(5) Stage 5 refinement

Figure 64. Adaptation by Ear, using fixed-percentage method

7.2.2 2D Chamber Problem

In this subsection, we will commit adaptations for the 2D chamber problem by
using different criteria. The goal is to further confirm the effects of using $011. the
exact relative error, and absolute error as the error indicator. For all h-type mesh
refinements described in this section, the starting “coarse” mesh and solution are
taken from the level 3 case in Table 21, which has a grid spacing of 3.125x10'3 m

and 1,472 cells.

0 Adaptation committed by (Don

First, we still commit the refinement by using (>03. It can be seen from Table 28

that it is not helpful to reduce the average error if only the regions with high than
get refined. Figure 65 (a) shows the mesh and consequent exact relative error

after 4 refinements.

Table 28. Fixed-percentage adaptation process by using 003

 

 

 

Refined rpm Marked Total Eat, Er.
Stag 4%) Cells Cells (m/s) (%)
(flinal 1,472 3.393x10“ 13.35

1 9.3 371 2,535 3.207x10'4 14.72

 

10.4 642 4,565 3.335x10'4 16.36

 

 

2
3 10.5 1,151 3,153 3.571x10“ 17.97
4 10.4 2,054 14,463 3.545x10‘4 13.16

 

 

 

 

 

 

 

 

161

Mesh Actual Relative error ¢on pattern

(a) Starting
Mesh

(b) Adaptation
BY ¢DR

 

 

 

(c) Adaptation
By Ere

(d) Adaptation
By Eab

 

Figure 65. Adaptation process in 2D chamber test case

162

. Comparison of E,,, and E31, pattern

Next, let us compare the exact relative error and absolute error pattern on the
starting “coarse” mesh, which is shown in Fig. 66. All the regions with high error
are labeled for the convenience of discussion. The location and cause for high
relative error has been analyzed in previous chapter. For the regions with high
absolute error, it can be seen from Fig. 66 (b) that Region #1 contains small
vortical structures. Regions #2 and #3 locate on both sides of the main flow.
There are also some regions with high absolute error inside the inlet and exit
ducts. Comparing the two patterns, it is clear that they look very different except
the overlap in region #1. This means that different regions will be refined when

the adaptation is committed by Eli and Eab, respectively.

 

 

 

 

 

(a) Ere (D) Em

Figure 66. Exact relative and absolute error

163

. Adaptation committed by Ere

1n the adaptation process committed by Era, totally 4 stages of refinement are
performed and in each stage approximately 25% cells with highest Em get
refined. The total cell number in final stage is 14,882. Figure 65 (b) shows the
mesh and consequent exact relative error after 4 refinements. For each stage,
the refinements are mostly embedded within previous refinements. By
comparing the final error with the starting error, regions #2, #3, and #6 get much
improved, and the errors in regions #1, #4, and #5, which contains small
vortexes, only reduce a little bit. Table 29 lists the detail adaptation process by
using Era. The average relative error drops down 1.77% after the first stage of
refinement and only 1.19% in the second stage adaptation. Then, the error level

almost keeps unchanged when further refinements are committed.

Table 29. Fixed-percentage adaptation process by using Em

 

 

 

 

 

 

 

 

 

 

 

 

Refined q)", Marked Total Eat, Er.
Stage (%) Cells Cells (m/s) (%)
Origijal 1,472 3.393x10" 13.35
1 19.2 363 2,576 3.095x10“ 12.03
2 22 641 4,676 3.026x104 10.39
3 27.3 1,149 3,331 3.050x104 10.36
4 33 2,053 14,332 3.020x10'4 10.72

 

 

164

. Adaptation committed by Eab

For the adaptation committed by Eab, four successive refinements are performed
and the total cell number in final stage is 15,248. It can be seen from Table 30
that the average relative error reduces 5.7% after the first stage of refinement
and 4.05% in second stage adaptation. Compared with the error reduction by
using Ere, the efficiency by Eat, is much higher. Figure 65 (c) shows the mesh
and consequent exact relative error after 4 refinements. As observed, the grid is
mainly refined in the inlet/exit ducts and the top part of chamber. These regions
have the highest absolution error. Afterwards, the high relative errors in the
lower part of chamber are almost eliminated even the grid quality in that part is

notchanged.

Table 30. Fixed-percentage adaptation process by using Eat,

 

 

 

 

 

 

 

 

 

 

 

 

Refined q)", Marked Total Eab Er.
Stag: (%) Cells Cells (m/s) (%)
Original 1,472 3.393x10‘4 13.35
1 4.40x10'4 367 2,573 1.309x10" 3.15
2 3.31 x1 0'4 645 4,530 7.284x10'5 4.1 1
3 3.01x10'4 1,150 3,291 5.022x10'5 3.97
4 2.42x10'4 2,093 15,243 5.072x10'5 4.02

 

165

 

7.3 Error Source Region and Error Apparent Region Model

As observed from the two examples in previous section, when a mesh is refined
in one region, the relative error can be reduced not only in the refined region but
also in some other unrefined regions. This is especially true when the absolute
error is used as the error indicator. On the other hand, when we exactly refine
the high relative error regions, it seems not the most efficient way to reduce the
local error. In this section, an ESR/EAR model is described based on the
previous results. This model is trying to figure out the possible correlation
between the global grid quality and solution errors. New adaptation strategies

will be developed based on this model.

1. Definition of ESR & EAR

We hypothesize that for a given flow field, there exist two regions that are related
to the actual solution error, and thus need mesh refinement. One region,
referred to as “Error-Source Region” (ESR), is where relative error initiates but is
not apparent. The other region, referred to as “Error-Apparent Region” (EAR), is

where the relative error presents high value.

11. Nature of error generation
There are two necessary factors to determine the final error level in EAR. First,
error initiates in the ESR when the grid in ESR is not fine enough to minimize the

local truncation errors. Second, the local grid quality in EAR is not good enough

166

to resolve the local flow pattern. In other word, if there is no error source
initiated, even the local grid quality in EAR is bad, probably no high error will
show up in EAR. On the other hand, if error source has initiated, the error will

expand into the flow field and easily present high relative error in EAR.

lll. Correlation between error generation and grid-quality
The initiation of error source in ESR and the appearance of actual error in EAR
are determined by the local grid quality. According to previous discussion,
measure 003 is a good indicator for relative error. So it is more reasonable to call
¢0n as the Error-Apparent Indicator (EAI). The features for EAl can be
summarized as:

- Combination of the geometry and solution information based physics;

- Good correlation with the exact relative error;

- Dimensionless, general for different ZD/3D cells;
The position of ESR and EAR may be overlapped but generally they are
separated. So, another indicator, which can accurately indicate the ESR, is
defined as “Error-Source Indicator" (ESI). A good ESI should have following
features:

- Combination of the geometry and solution information based on physics;

- Good correlation with possible error source, like absolute error related

regions;

- General for different ZD/3D cells;

167

IV. Grid adaptation

When the adaptation is committed, grid quality in ESR needs refined first to
minimize the possible error source. Next, the grid quality in EAR should be
improved so as to reduce the possible apparent error. When the grid quality in

both regions is improved, the error level should be reduced most efficiently.

This model is also a good explanation to some difficulties we once met in
previous analysis. First, in the previous effort of establishing a perfect function
between the local solution error and local grid-quality measure, we found that the
local error can not be totally determined by the local grid quality. Large data
variances were found in many partitions. In other word, we noted that when the
local grid-quality measure is pretty good, sometimes the local error is still very
high. On the other hand, when the local grid quality is very bad, sometimes the
local error is not so bad. Based on the ESR/EAR model, we hypothesize that the
generation of error depends on grid features in both ESR and EAR. The grid-
quality in ESR controls the initiation of error source and the grid quality in EAR
controls the apparent actual error. 80 the grid quality in EAR can not totally
determine the local error. This is also the answer to the third question in section

6.1.

Second, in order to correlate the local grid quality with the local solution error, we

make all the efforts to look for a complete set of measures that can fully describe

the local grid quality. The important grid feature should at least include the cell

168

resolution, cell smoothness, and cell skewness information. Based on the
statistical analysis and numerical experiments, 11109, (bps, and 053 prove to be
correlated well with the solution error. These measures can account for the
resolution of recirculation flow and represent the local cell smoothness feature.
For cell skewness, it seems not so critical to the solution error as resolution and
smoothness type measures. As for gradient-related measures, we developed
measures like 069, com, and 0WD. The previous analysis shows that these
measures do not locally correlate with relative error. However, as we know, the

poor grid quality in high gradient regions will definitely increase the solution error.

The ESR/EAR model emphasizes that the error should be determined by grid
quality in both ESR and EAR. The absolute error is a measure of truncation error
and has certain relationship with ESR. So it is probably that the regions with high
solution gradient can correlate well with the absolute error and is a good indicator
for ESR. Then, the grid quality in high gradient regions will remotely influence
the local error in EAR, but does not show any local correlation with it. This is the

possible answer to the fourth question in section 6.1.

Third, in previous grid-adaptation study, we noted that when $03 or even the
exact relative error is used to make refinement, the error reduction is not very
efficient. The ESR/EAR model presents that EAR normally locates in the regions
with complex flow structures, like recirculating flow, stagnation flow, and the

zones where two flow merge and separate, etc., while ESR possibly relates to

169

the high gradient regions, like boundary wall. 30 EAR and ESR are probably in
different locations. $03 or the exact relative error can only indicate the EAR. The
error can not be completely removed until the grid quality in ESR also gets

refined and the error source is eliminated.

7.4 Error. Source Indicator (ESI)

7.4.1 Two new developed ESI

According to the previous results of adaptation by using the exact absolute error,
it shows that the absolute error is a possible indicator of ESR. But it is
impossible to know the absolute error before the final grid independent solution is
available. Furthermore, probably the absolute error is only related to ESR, but
not the best error source indicator (ESI). In this section, the goal is to find out

some good solution-based ESIs.

The discussion in previous section gives us a direction to look for ESI. First, a
good ESI should combine the geometry and solution information together based
on certain physical meaning. Second, ESI should have certain correlation with
the exact absolute error. Third, ESI may be gradient-related measures. After
extensive numerical experiments, two measures are found to have certain

correlation with the absolute error and also account for the gradient information.

170

Measure I: Velocity Gradient 11>.
This measure is defined as the undivided Laplacian of the selected solution

variable multiples a characteristic length scale, which is given by:
¢,-%¢MW%V1$ (In

The length scale is the square root of cell area and the solution variable is using
the velocity magnitude. This measure has been popularly used in many

commercial codes for gradient adaptation.

Measure 11: Velocity Difference 0..
This measure is defined as:
(1)" = max (AUi), i=1,2,..,N (7.2)

where (see Fig. 67)

Au=fiqfifil (73

EC; = (xi — x0); '1' (yr — yo);t (74)
170 = 140? + v0] (7.5)
17‘. = u; + vi; , (7.6)

Where, EC, is a vector pointing from the center of a neighboring cell i (i=1,,..2, n)
to the center of the given cell. 170 is the velocity vector on cell Co. 17”. is called

effective velocity, which is the projection of velocity vector l7,on {I}, direction. N in

Eq. (7.2) is the number of neighboring cells, whose angle 3., defined by

171

Q = COS-1(uo(xi-x(:)+t:o(yi— yo))
leCJIIVOI

3; is between 45 and 135 degrees. This measure indicates the gradient of the

(7.7)

 

velocity vector in the direction normal to the streamline.

 

 

 

 

Figure 67. Definition of terms in 0.. (Eqs. (7.2)-(7.7)).

According to above definitions, both measures combine the local geometry
feature with the solution information and physically account for high gradient

regions.

Next, we will evaluate their relationship with the absolute error. Similarly as what
we did to evaluate the previous measures based on three sub-databases, each
database is partitioned by the absolute error and the average value of 0. and 0.. in
each partition is computed. The results are shown in Fig. 68. As seen from Fig.
68 (a), for partitions with small E», the average value of 41. keeps in low level.
For partitions with high absolute error, 0. value is also high. The correlation curve

is almost monotonic for 2D test problem and 3D test problem I, and has some

172

oscillation for 3D test Problem ll. Figure 68 (b) presents the correlation between
0.. and the absolute error. Similarly as 0., 0.. also shows strong correlation with
the absolute error for 2D problem and 3D problem I. As for 3D problem II, it is

not so good as other cases but still presents certain correlation with the absolute

  

  
      
   

 

   

 

error.
2D Problem
3D Problem I o_oo3 _
3D Problem 11 —
0.003 11
0.002 -
0.002 .0 t
a 13
0.001 0001 ,_9_ 2DPmblcm
+ 3D Probleml
—+—— 3D Probleml]

 

 

 

 

 

0 515-07 1E-061.5E-06 25-06 0 l t ‘05 A A 1 10.004‘
010/1118) %2(m/5)
(8)9. (bliu

Figure 68. Evaluation of 0. and 0..

7.4.2 Test by 2D Moving Wall Problem

. Adaptation committed by 0.

Figure 69 shows the 0. pattern on the starting coarse mesh. It can be observed
that the regions with highest 0. mainly locate under the moving wall and

surrounding the large vortex at the right part of cavity. It looks very similar to the

pattern of exact absolute error, which is given by Fig. 63.

173

 

 

Figure 69. 0. pattern

As what we did in the previous section, 5 stages of refinement are committed by
refining fixed-percentage cells with the highest 0. value each time. Figure 70
shows the refined mesh and the consequent exact relative error pattern of each
stage. In the first-stage, 288 cells are marked with highest (11. (greater than
4.05x10'7) and total cell number becomes 2,016 after refinement. As can be
seen from the figure, the mesh spacing is only refined in regions under the
moving wall and surrounding the large recirculating flow. However, the error
level in the whole cavity is greatly reduced. This behavior is similar to what we
observed in the adaptation committed by the absolute error. From first to third
stage of adaptation, the refinement is still mainly embedded within the previous
refinements. When the refinement continues to the fourth and fifth stage, the
refined region expands to other areas. The global error keeps going down and

finally decreases to a very low level.

174

Mesh Actual Relative error ¢Dn Pattern

    
 

(0) Starting case

 

(1) Stage 1 refinement

                
    
    
    

«hi
111‘
\l
‘1.
.1

(2) Stage 2 refinement

r» _ 1 r r:
- a 2 -f‘ _ a;
3* I a,
nun-nun“ I ‘- Ax, a". h 5‘. j
(3) Stage 3 refinement
'7 -" “-12 ' \ M
v: 5: ,
a 3 . '1: 4e
3’ j
-- . . : c- X. .r .. L k
(4) Stage 4 refinement
r» M “E 1 - - __\___
e "t, 2k
. .1
$— ‘ «a .r r L :\ {"7

(5) Stage 5 refinement

Figure 70. Adaptation by 6. using fixed-percentage method

175

Table 31 lists the details of adaptation process by using 0.. Compared with Table
30, the efficiency of error reduction by 11>. is even higher than that by the absolute
error. As noted, the average error of final stage by 0. is 2.94% with total cell
number equal to 19,908, while the final error level by Eu, is only 5.94% with total
cell number equal to 20,202. It means that 0. is a better error source indicator

than the exact absolute error for this case.

Table 31. Fixed-percentage adaptation process by using it).

<|>m

 

. Adaptation committed by 0..
Figure 71 shows the 0.. pattern on the starting mesh. As noted, this pattern is

very similar to 0. pattern.

 

 

Figure 71. 0.. pattern

176

Table 32 summarizes the detailed adaptation process by using 0... In the first
stage, the average error level drops from the starting 16.31% down to 6.30%.
This efficiency is higher than that of using Em and 0., which are only 9.57% and
7.71%, respectively. When the refinement continues, the error level keeps going
down in the second stage and then remains the same level in the further
refinements. The final average error is 4.14% with total cell number equal to
19,896. Its final efficiency is better than using Eab (5.94%), but not so good as
using 0. (2.94%). Figure 72 shows the refined mesh and the exact error on each
stage, it can be seen that all the refinements are embedded in the previous
refinements instead of expanding to other regions. The evolution of error is

consistent with Table 32.

Table 32. Fixed-percentage adaptation process by using 0..

 

 

 

 

 

 

 

 

 

 

 

 

 

Refined q)", Marked Total Eab E...
Stage Cells Cells (m/s) (%)
Original 1,152 2.855x10" 16.31
1 6.5x10'4 233 2,016 1.073x10" 6.30
2 8.2x10" 507 3,576 6.823x10'5 4.65
3 8.9x104 900 6,313 5.951x10’5 4.30
4 8.1)(10‘4 1,641 11,274 5.662x10'5 4.21
5 6.25x104 2,355 19,396 5.920x10’5 4.14

 

 

177

Mesh

  
   

(0) Starting case

 
 

(1) Stage 1 refinement

 

  

(2) Stage 2 refinement

 

(4) Stage 4 refinement

(5) Stage 5 refinement

    

a. a... I". 5

 

Actual Relative error

 
 

0D... pattern

55. I
i. ' ' z;- -
e . a?
i- i i. J

E?“
..___‘ i E. “—1
.. - t 3
" a, '3 ((1 B i J
.. ”3" , .‘_ .- g

x J .1—

Figure 72. Adaptation by 11).. using fixed-percentage method

178

. Fixed-percentage adaptation summary

In this subsection, we will make a summary of the previous adaptation processes
committed by 110..., 15..., 15...... 1b., and 11>... Figure 73 shows the error reduction
efficiency of these error indicators. In the figure, X-axis represents the total cell

number in each simulation and Y-axis is the average relative error.

0.)... and Em indicate the regions with high relative error, which are defined as
“EAR”. As observed, their curves almost match together and have the lowest
efficiencies among all error indicators. For the adaptation performed by E», the
refined regions correlate with the ESR and its error reduction efficiency is much
better than that of (Don and Er... As for the refinements by using 0 and 11>... they
behave very similar and have the best efficiency. Furthermore, comparing with
11>", 11>. is less efficient at the initial stages and much more efficient in the final
stages. So 11>. and 1).. can be considered as better ESI than the exact absolute
error. Mostly important, 11>. and 11>" can be directly computed based on the
geometry and preliminary solution information on the coarse mesh without the

grid-independent result.

179

 

 

 

 

 

 

 

 

Ere (%)

- —B— ¢=¢DR ‘9
""A'" ¢=9re
l—_ ¢=eab
""°'" 11:41.
-—0—- ¢=¢ll

 

 

 

 

 

1 1 1 1 t1 1 1 1 § 1 1 1 1 1 1 1 1 1 J
0 5000 1 0000 1 5000 20000
Cell No.

Figure 73. Summary of fixed-percentage adaptation

According to above results, it can be clearly seen that the refinement committed
only on the ESR has higher efficiency in error reduction than the adaptations only
on the EAR. In Fig. 73, we also found that 0. behaves a little special compared to
other measures. For the curves of other measures, like 00.1, E19, E1... and 0.., the
error reduction in the initial stages is pretty high but not so efficient in latter
stages. However, for measure 0., even in its final stage, the error level still has a
big drop. What makes 0. so special? By comparing the refined mesh on final
stage of 0. and 0.., we noted that for 0.., all the refinements are embedded in
previous refinements. This also happens to adaptation processes by using 00.1,
E,.,, and E31,. As for 0., the final refinements expand to other unrefined regions

instead of being embedded in previous refinements.

180

7.5 Fixed-Criterion Adaptation

As mentioned, two different approaches were popularly used in grid adaptation to
determine the number of refined cells during each stage: fixed-percentage and
fixed-criterion. In previous sections, all the refinement is based on the fixed-
percentage. In the fixed-criterion approach, the user specifies a bound on the
error indicator. All cells having error indicators greater than the specified bound

will be refined.

In this section, the 20 moving wall case is taken as the example and the fixed-
criterion approach will be used to commit adaptation by taking <l>on. E19, E31,, 0.,
and 0.. as error indicator. The goal is to analyze their behaviors in error
reduction, and then seek new strategy to further improve adaptation efficiency

based on our ESR/EAR model.

. Adaptation committed by 003

Table 33 lists the detailed information of the adaptation committed by error
indicator 0.3.1.. The threshold value is set to 5, which is approximately compatible
to the average 5% relative error. Totally three stages of refinement are
performed. Figure 75 (1) shows the final mesh and its exact relative error
pattern. Comparing with the starting error pattern, the global error level on the

final stage is greatly reduced. As seen from Table 33, the error decreases to

181

10.19% after the first stage and finally reduces to 8.80%. It also can be noted
that in the first stage approximately 70% of current cells get refined and then the
number of cells that get refined becomes smaller and smaller as the refinement

process proceeds.

Table 33. Fixed-criteria adaptation process by using 0.3.1, 01.. = 5

 

 

 

 

 

 

Refined Marked Marked Total E81, E...
Stage Cells % Cells (m/s) (%)
Original 1,152 2.355x10‘4 16.31
1 794 69 3,534 1.793x10" 10.19
2 1 ,627 46 3,466 1 .546x10‘4 3.60
3 2,513 30 16,053 1.543x10‘ 3.30

 

 

 

 

 

 

 

. Adaptation committed by Era

Table 34 summarizes the adaptation by using Era. The bound value is set to 5%.
In the first stage, almost all the cells get refined and the error rapidly reduces to
4.77%. Figure 75 (2) gives the final mesh and error pattern of this adaptation.
As noted, all the possible EAR and ESR get refined in the first stage. When
further refinements are committed, only the cells in EAR get refined. Thus, the

errors do not reduce as efficiently as that of the first stage.

Table 34. Fixed-criteria adaptation process by using E11,, 0... = 5.0%

 

 

 

 

 

Refined Marked Marked Total E... E...
Stage Cells % Cells (m/s) (%)
Original 1,152 2.855x10’4 16.31
1 1,043 91 4,231 8.105x10'5 4.77
2 1,532 37 9,034 6.846x10'5 4.45
3 4221 46 22,093 5.521x10'5 3.75

 

 

 

 

 

 

 

 

182

. Adaptation committed by E1...

Table 35 shows the adaptation process by using E31,. In this process, an average
absolute error of 8.0 x 10‘5 is approximately equal to an average relative error of
5%, and is set to the bound value. As observed from table, the average error
reduces to 3.54% after the first refinement stages. Comparing with the result of
using fixed-percentage approach, in which the average error only decreases to
9.57% for the first stage, the efficiency by using fix-criterion is much higher.
Figure 75 (3) gives the final mesh and error pattern of this adaptation. By
comparing the refined meshes based on two approaches, we found that for the
fixed-percentage approach, the mesh spacing is only reduced in the highest
absolute error region, which is the possible ESR. When using the fixed-criterion
approach, it improves the grid quality not only in possible ESR but also in some
EAR. This can greatly increase the efficiency in error reduction. Further
refinements mainly concentrate on the highest absolute error regions and the
error continues to go down. The final average error reduces to 1.98% with the

total cell number equal to 22,437.

Table 35. Fixed-criteria adaptation process by using E31,, 0... = 8.0><10‘5

 

 

 

 

 

 

Refined Marked Marked Total E51, E...
Stage Cells % Cells (m/s) (%)
Ofl'ginal 1,152 2.855x10“ 16.31
1 339 77 3,319 5.930x10'5 3.54
2 1 ,215 32 7,479 3.702x10's 2.60
3 1,977 26 13,512 2.619><10'5 2.06
4 2,917 22 22,437 2.158x10'5 1.93

 

 

 

 

 

 

 

 

183

. Adaptation committed by 0. and 0..

Table 36 and 37 present the refinements performed by 0. and 0.., respectively. In
these processes, an average 0., 0.. of 5.0 x 10'8 and 2.0 x 10“ are approximately
equal to an average relative error of 5%, and are set as their bound values. As
seen from the table, the final average error level reduces to 2.05%, 3.77% with
the cell number equal to 14,526 and 21 ,945, respectively. Figure 75 (3) and (4)
give the final mesh and error pattern of these two adaptations. As seen from
their first refinement stage, the grid quality in the ESR and part of the EAR gets

improved so the error reduction efficiency is much higher than only refining the

ESR.

Table 36. Fixed-criteria adaptation process by using 0., 0... = 5.><10'8

 

 

 

 

 

 

 

 

 

 

 

 

Refined Marked Marked Total E1... E1.
Stage Cells % Cells (m/s) (%)
Original 1,152 2.855x10’4 16.31
1 353 74 3,711 9.545x10'5 5.16
2 1,563 42 3,451 3.133x10'5 2.34
3 1,327 16 12,495 2.530x10'5 2.05
4 665 5 14,526 2.603x10'5 2.06

 

Table 37. Fixed-criteria adaptation process by using 0.., 0... = 2.x10"

 

 

 

 

 

 

 

 

 

 

 

Refined Marked Marked Total E1... E...
Stage Cells % Cells (m/s) (%)
Original 1,152 2.3551104 16.31
1 732 63 3,493 9.141x10‘5 5.44
2 2192 63 10,137 5.189x10'5 3.33
3 3920 33 21,945 4.863x10'5 3.77

 

184

 

 

. Fixed-criterion adaptation summary

For the convenience of comparison and discussion, Figure 74 shows the error
reduction efficiency of 001., E19, E11,, 0., and 0.. by using the fixed-criteria approach.
As seen from the figure, refinement by 0D... is still the least efficient in error
reduction compared with other error indicators. The curves of E... and 0.. are
almost same to each other and have the medium efficiency. The error reduction
behavior of E1... and 0. are similar and have the highest efficiency. Additionally,
for all curves, the error reduces very quickly in the initial refinement stages and

almost levels off when the adaptation proceeds.

 

 

—E— ¢=¢Dn
--«5-- ¢=ew
—‘V_ ¢=Oab
---o--- ¢=¢.
—'O-— ¢=¢u

Ere (%)

 

 

 

 

 

 

 

01 1 1501001 110000 1150001 12100001 J
Cell No.
Figure 74. Summary of fixed-criterion approach
According to the adaptation results by using fixed-criterion approach, it further
confirms the possible existence of Error Source Region and Error Apparent
Region. So, this ESR/EAR model can reasonably explain the different behaviors
in error reduction by using fixed-percentage and fixed-criterion approaches. In

following sections, our goal is to develop new strategy of adaptation that can

further improve the efficiency based on this model.

185

Mesh

  

(0) Starting case

 

(1) Ad ptation by 0.3.1

    

 

 

 

    
  

(2) Adaptation by Ere
:: 5‘25.
(3) Adaptation by Eab
:33: - é .,
(4) Adaptation by 0.

 

(5) Adaptation by 0..

Actual Relative error

09;. pattern

 

. a: ‘
.’ 1. Z. A

a //

l w...

11*"

:R

a. .1

"as

.1113.

Figure 75. Adaptation by using fixed-criterion method

186

7.6 New Adaptation Strategy Using ESR/EAR Concept

As what we learned in previous analysis, if the refinements are completely
committed on EAR, the efficiency is worst. If the refinements only concentrate on
ESR, the efficiency is better but still not the best. By using the fixed—criterion
approach, if both ESR and EAR can get effectively refined when setting a very
low bound value, the efficiency is the highest. This result are encouraging in the
sense that if we can accurately predict the ESR/EAR using error indicators and
then refine these regions simultaneously or alternatively, maybe that is the

optimal approach for mesh refinements.

7.6.1 2D moving wall case

Based on the new ESR/EAR model, the error is determined by both ESR and
EAR. The high gradient region is the possible ESR, which can be indicated by
measure 0. or 0... When the grid quality in ESR is poor, the error source is easily
initiated and then expanded into flow field. Thus, this flow field is highly polluted
by the potential errors. Next, than can accurately indicate the regions with
complex flow structures and is a good error apparent indicator. If the local grid
quality in EAR is also bad, the potential errors will be activated and present as
high relative error. According to this ESR/EAR concept, only reducing the grid
spacing in EAR can not really eliminate the apparent error in EAR because the

error source still exists. On the other hand, though improving the grid quality in

187

ESR can be very helpful to reduce the error source, it does not directly work on
the apparent error. So, if the ESR/EAR concept is reasonable, a new adaptation
strategy that alternately improves the grid quality in both ESR and EAR should

be more efficient for error reduction.

. Alternate adaptation approach A

First, we use (Pu and 003 as the error indicators for ESR and EAR, respectively.
Table 38 gives the detailed alternate adaptation process. In the first stage, 302
(around 25% of the current cells) cells are marked with highest 0.. (greater than
6.3x10‘) and total cell number becomes 2,058 after refinement. In second
stage, we continue to reduce the mesh spacing in ESR. After two stages of
refinement by 0.., the error reduces to 4.47%. Then, 0.3.; is used to refine the
EAR. Finally, another ¢u refinement and one more 0.3.. refinement are committed
alternately to further improve the grid quality in ESR and EAR. In the final stage,
the error decreases to 2.25% with cell number equal to 19,521. Figure 76 gives

the final mesh and exact error pattern of this adaptation.

Table 38. Alternate adaptation A by using 0.. and 0.3.1, fixed-percentage=25%

 

 

 

 

 

 

 

 

Refined 11> 11>“1 Marked Total E1... E...
Stgge Cells Cells (m/s) (%)
OLginal 1,152 2.355x10‘ 16.31
1 0.. 6.3x10‘4 302 2,053 1.064x10" 6.21
2 0.. 3.51104 501 3,594 6.690x10'5 4.47
3 00,. 7° 944 6,474 5.172x10'5 3.20
4 0.. 6. x10‘4 1,660 1 1,499 4.518x10’5 2.99
5 1.101. 5° 2641 19,521 3.715x10’5 2.25

 

 

 

 

 

 

 

 

188

(0) Starting case

 

(1) Adaptation by 0..

N

    

(2) Adaptation by 0..

 

‘5.._____.
11,
J j
- r. n. \
(3) Adaptation by 00;.
_ _ a , . 7
e E3
- ‘1 7%

    
  

(4) Adaptation by 0..

I'M"
l
r
it
I
\

 

(5) Adaptation by 0.3;.

Figure 76. Adaptation by alternate approach A

189

. Alternate adaptation B

Table 39 shows another type of adaptation process by alternately using $11 and
0D... In this approach, the 0.. and 0.3.. are set to a low bound. In the first stage, all
the cells with 0.. larger than 1.x10" get refined. Then one more stage is
committed by 0.. to further reduce the error source. Finally, the grid quality in
EAR is improved by refining cells with highest 0.3.1 (larger than 5). The final

average error decreases to 1.64%. Figure 77 gives the final mesh and exact

error pattern of this adaptation.

Table 39. Alternate adaptation B by using 0.. and 00.1, set 0..:1. x104, 0nn=5.

 

 

 

 

 

 

 

 

 

 

 

 

 

# 11> Marked Marked Total E1... E19
Cells % Cells (m/s) (%)
Origllal 1,152 2.355110" 16.31
1 .1... 395 77 3,337 7.343><10'5 4.50

2 .1... 1,752 45 9,193 3.520><10'5 2.97

3 11.0.. 1,667 13 14,214 2.253x10'5 1.64

 

190

Mesh Actual Relative error 0.... pattern

If.“ :1:
., *1; J

1 h a

   

(0) Starting case

 

(1) Adaptation by 0..

 

(2) Adaptation by 0..

 

.1‘....' 2.. .11 ' ' cg" ' ' : . ,. . .
1 U V. 1" t';.fl,..°.‘.; '
' .1 o . ..' .; . .' '.
' N . I: .' "' 9

1:: -- ., ‘

:3 .: , t ‘

1. , _ . ' .
.:E:‘ 3' -
. 1: 4' 5; :::;’ E? .-
"I 1 ' 4 : :1“ " '
__ :.5 .. .. .‘ .

:5 3 \1 . g
H ‘0. III-III1§.~ .. ' ‘ _ - ..“ 3;; ‘ ’ ‘3
- m: .. ...'.-..." 3: 1....... in . . 9'. _. ~ 0 ' k . .

(3) Adaptation by 0.3..

Figure 77. Adaptation by alternate approach B

191

. Efficiency comparison

Figure 78 summarizes the efficiency of the adaptation committed by E11,, 0., 0..,
and two alternate approaches. It can be seen from the figure that even the exact
absolute error regions get refined, the efficiency is not perfect. It is even lower
than the efficiencies of the new developed ESI 0. and 0... Based on the ESR/EAR
model, we developed a new adaptation strategy by alternately refining the ESR
and EAR. The result shows that the alternate approaches can further improve

the efficiency of error reduction.

 

 

is?

E.”
--«A—-- ¢=9ab
--5b--¢=¢.
"0"' ¢=¢ll

—-B— 0 = Alternate A
—A-— 0 = Alternate B

 

 

 

 

1111111111111Ll1111l

0 5000 1 0000 1 5000 20000
C611 No.

Figure 78. Efficiency comparison

7.6.2 20 Chamber Problem

In order to further confirm the results we got in previous subsection, we perform

same numerical tests on the 2D chamber problem.

192

o Adaptation committed by 0.

Figure 79 shows the 0. pattern and Table 40 presents the details of adaptation
process by using 0. based on the fixed-percentage approach. Figure 80 (b) gives

the mesh and exact error pattern on the final refined mesh.

 

 

 

Figure 79. 0. pattern

Table 40. Fixed-percentage adaptation process by using 0.

0m

 

193

Mesh Actual Relative error 0.... pattern

(a) Starting
Mesh

(b) Adaptation
By 91

(c) Adaptation
BY ¢ll

 

Figure 80. Adaptation by 0. , 0..

194

o Adaptation committed by 0..

Figure 81 shows the 0.. Pattern. Table 41 summarizes the details of adaptation
process by using ¢ll based on the fixed-percentage approach. Figure 80 (c)

gives the mesh and exact error pattern on the final refined mesh.

 

 

 

Figure 81. 0.. pattern

Table 41. Fixed-percentage adaptation process by using 0..

¢th

 

o Efficiency Comparison

Figure 82 shows the error reduction efficiency of these error indicators. It can be

seen that the efficiencies of 0.3.1. and E... are the lowest among all error indicators.

195

The curves of 0.. and Eab almost match together and show the highest efficiency.
The efficiency in initial stage of 0. is not very good, but finally it achieves the
same level as that of 0.. and E1... This result is consistent with what we got in the

previous 20 moving wall problem.

 

155-
10:-

 

 

Ere (%)

----- ¢=em
_ -—"— ¢=9ab

--«>-- ¢=¢1
-—'0— ¢=¢ll

 

 

 

J

 

l I L I

 

0” 5000 510000' '15000
Cell No.
Figure 82. Efficiency comparison

. Alternate Adaptation
Table 42 and 43 gives the details of two alternate adaptation approaches. Figure

83 gives their mesh and exact error pattern on the final refined mesh.

Table 42. Alternate adaptation A by using 0.. and 0.3a, fixed-percentage=25%

 

 

 

 

 

 

 

 

Refined 11>... Marked Total E... E...
Stage Cells Cells (m/s) (%)
@inal 1,472 3.3931104 13.35
1 1.2x10'3 374 2,594 9.929x10'5 4.45
2 1 .3x10’3 337 5,117 4.457><10'5 3.61
3 3.3 1,231 3,960 3.895x10’5 3.06
4 9.0x10" 2,323 16,049 3.920x10'5 2.59

 

 

 

 

 

 

196

Table 43. Alternate adaptation B by using 0.. and 0.31:1, set 0..:1. x104, ¢0a=5.

# ¢

%

 

Mesh Actual Relative error 0.... pattern

Alternate
Adaptation A

Alternate
Adaptation B

 

Figure 83. Adaptation by alternate adaptation A, B

197

. Efficiency comparison

Figure 84 presents the efficiency of the adaptations committed by E31,, 0., 0.., and
two alternate approaches. This result is almost same as what we obtained in
previous 2D moving wall problem. As noted, even the exact absolute error is
used as the error indicator, the efficiency is not perfect. The efficiencies of using
0. and 0.. is similar to that of E8... Both alternate adaptations have much higher

efficiencies than using only ESI or EAI.

 

 

 

 

 

 

 

20 g‘
15 g—
10 :—
53 5 —
[11.3 r
~ --A—-- ¢=°ab
—-v—-- ¢=¢.
""0"" ¢=¢ll
l—B— 0=AlternateA
l—A— 0=AlternateB
1.111111111111111
0 5000 1 0000 1 5000

Cell No.
Figure 84. Efficiency comparison

198

Chapter 8

Summary and Conclusion

This paper presents a framework to study the correlation between grid quality
and grid-induced errors in computed CFD solutions of complex flows. Based on
the hypothesis that the grid-induced errors defined with respect to the grid-
independent solutions is a function of a set of grid-quality measures and that this
function may be fairly universal, six grid-quality measures were developed so as
to provide a posteriori estimates of the solution errors. These measures can
account for both local grid features and the scalar, vector, and tensor nature of
the computed solution. Among them, the directional resolution (0.3.1.) gauges the
resolution of the directional change in the velocity vector. The gradient resolution
(061.), the deformation tensor resolution (0.3.3), and the relative velocity difference
(0.11m) are measures relating to the local flow gradient information. The
directional smoothness (005) and the streamline smoothness (035) describe the

local smoothness feature of a given mesh.

These measures are evaluated using statistical approaches based on three test
problems. One problem involves a 2-D flow in a square chamber with two inlet
ducts and one exit duct. The second one describes a 3-D flow in a cubic
chamber with four curved inlet ducts and one exit duct. The third is a 3-D flow in
a cubic chamber with two opposing inlet ducts and two exit ducts. The result of

evaluation shows that all these solution-based measures have certain correlation

199

with the relative error. Among these measures, 0.3.. has the strongest correlation

with the relative error in all simulations. 003 also correlates well with error but not

as good as (Don. 093 performs very similar to 0.33. 09.1, 09..., and (lava only

present very weak correlation with solution error.

However, by checking the detailed data distribution in all partitions, it can be
noted that the data variance in many partitions is very large. One possibility is
that the local error may be determined by multiple measures instead of by only
one measure. Thus, the function E", = f(0., 02, ...) was established to further
analyze the correlation between error and local grid quality. For the tests of
functions Ere = «0.3.1, 0.3...) and E... = f(0on, 03.1, 00...), the result shows that the
data variance in some partitions is still very large. This means that the error
cannot be accurately predicted even multiple measures are used, and the
correlation between the grid quality and grid-induced error is extremely
complicated. Though there has correlation between local grid-quality measures
and solution error, probably the error can not be quantitatively estimated by only

using comprehensive local geometry and local solution information.

In order to check how useful these measures are when applied as error
indicators in practical problems, a set of numerical test was committed based on
two 20 examples. Of these two examples, one is a 20 Moving wall cavity flow
and the other is a 2D chamber flow with one inlet and one exit. The result of

these applications shows that the highest relative errors always show up at

200

certain regions with complex flow structures like recirculating regions, vortex, and
separation zones, etc. 0.3.. can successfully indicate almost all the high error
regions and is the best error indicator among all measures. 005 can also predict
some high error region but not as good as 0.3.1, and 053 performs very similar to
0.33. The ability of predicting errors by 0G9, 0.3711, and (have are still the worst
among all measures. These results are very consistent with our previous result

of evaluation.

In order to study how to use these measures in mesh refinement to reduce
solution errors, a series of adaptations are committed based on above 2D
examples by using different criteria including real errors Ere and E31,, and three
measures. The three measures are 0.3.1, which can indicate high relative error
regions, and two newly developed measures 0. and 011. which are gradient-related
and closer to the absolute error. The result of experiments showed that the error
reduction by refining high 0.3.. and E... regions is not so efficient as refining high
Ea... regions. We also found that the local errors can be significantly affected by
grid quality in remote regions. Thus, this further confirms that the error can not
be completely determined by the local grid quality. Instead, the grid quality in
some specific regions can remotely control the error level. Further numerical
experiments showed that the error reduction by refining high 0. and 0.. regions is
even more efficient than refining high Eab regions. However, that is still not the
best way to reduce error. A new adaptation strategy, in which the grid quality in

high 0.. regions and high 00.. regions were alternately improved, was shown to be

201

superior to other indicators.

Thus, a new ESR/EAR concept is finally induced as a possible explanation of
error generation and the correlation between the grid quality and solution errors.
In this model, we hypothesize that in a given flow field, there exist two regions
that are related to the actual solution errors, and thus need mesh refinement.
One region, referred to as “Error-Source Region” (ESR), is where relative error
initiates but is not apparent. The other region, referred to as “Error-Apparent
Region” (EAR), is where the relative error presents high value. ESR normally
has very high gradient but the local grid resolution is not high enough. Thus, high
truncation errors are easily induced as error source and then expanded into the
whole flow field. Measures 0. and 0.. are good indicator to this region. EAR
normally has complex flow structures, like recirculating flow, separation flow, and
stagnation flow, etc., but the local resolution or smoothness is not good enough.
Thus, high error will be easily presented here it the whole flow field has been
heavily polluted by truncation error. Measure 0.3.. is a good indicator to EAR.
The position of ESR and EAR may be overlapped but generally they are
separated. When estimating the local error, the grid quality in both ESR and
EAR should be considered. When committing the adaptation to reduce error, the
grid quality in ESR needs refined first to minimize the possible error pollution.
Furthermore, the grid quality in EAR should be improved so as to reduce the
possible apparent error. When the grid quality in both regions is improved, the

error level should be reduced most efficiently.

202

Appendix A

FLUENT Code Validation

A.1 Introduction

The driven cavity flow is a validation problem, used as a benchmark by many
researchers. Accordingly, many accurate numerical results are available for
comparison (Ho, OJ. and Lin PH [1], Karadag, A. [2]). In this appendix, the
Fluent code is assessed by computing a 2D laminar driven cavity flow and

compared with the calculations presented in reference [1] [2].

A.2 Problem Description

A schematic diagram of this test problem is shown in Fig. 1. The cavity is a
square whose edges are 0.1 m long. The viscosity and density of the fluid are
0.001 kg/m-s and 1 mm", respectively. The lower and side walls are stationary
while the upper wall is moving at a constant rate of 1.0 We in the horizontal
direction. Reynolds number based on the edge length and upper wall velocity is

100. The no-slip condition is employed at all boundaries.

As noted, we are only interested in the steady-state solution of this

incompressible laminar flow. The equations govern the flow are the continuity

203

and momentum equations (full Navior-Stokes equation). Since p. at f(T), the
energy equation is decoupled from the continuity and momentum equations and

was not used in this study.

1%

 

—h

0.1m

\\\\\\\j
\\\\\V\

 

\
‘3,

//////////
1— 0.1m -->

 

 

Figure 1. Problem description of cavity flow

A.3 Grid

By using the same method described in Section 4.2.2, we obtained the grid-
independent solution for this 2D cavity flow problem. Table 1 shows the series of
uniform mesh. A total of 8 grids/solutions were generated/computed and the final
level contains 16,384 cells. Table 2 lists the average relative error and absolute

error between successive levels. Figure 2 shows the convergence rate when

grid gets refined. E], in Table 2 is only 0.169% and the solution on level 8 grid is

considered as the grid-independent solution.

204

Table 1. Uniform Grids Used to Get Grid-Independent Solution

 

 

 

 

 

 

 

 

 

 

 

 

Level Grid Spacing (m) Total Cell No.

1 12509-2 (1/80) 64

2 62509-3 (1/160) 256

3 31259-3 (1/320) 1 ,024
4 20839-3 (1/480) 2,304
5 15639-3 (1/640) 4,096
6 1.250e-3 (1/800) 6,400
7 89299-4 (1/1120) 12,544
8 78139-4 (1/1280) 16,384

 

Table 2. Convergence Process to Get Grid Independent Solution

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Levels E... (m/s) E1. (%)
1-2 3.406x10'2 19.67
23 1.342x10'2 7.17
3-4 2.657x10'3 1.52
4-5 9.575x10“‘ 0.730
5-6 4.663x10" 0.429
6-7 4.167x10“ 0.321
7-3 1 .320x10'4 0.169

0.2.;
015
l
30.. ..
II .
2
0.05:-
00 I I 1 5000 10000 A' 15000

C011 number

Figure 2. Convergence Process to Get Grid-Independent Solution

205

A.4 Results

Figure A.3 shows the streamline pattern of the grid-independent solution. The
direction of driven wall is from left to right with velocity magnitude equals 1 m/s.
There is one strong clockwise vertex at the center of cavity. Two small

circulation flows exist at the corners.

 

 

 

Figure 3. Flow pattern of grid-independent solution,

obtained by using a mesh with 128x128 aspect ratio unity cells.

Figures 4 and 5 show the horizontal velocity profile at the mid-plane in the x-
direction and vertical velocity profile at the mid-plane in the y-direction,
respectively. Computations by Ho, OJ. and Lin F.H. [1] and Karadag, A. [2] are
also included. As can be seen from these plots, the u, v profile in selected plane

generated by FLUENT are in good match with the benchmark data.

206

 

 

 

Fluent
» a He at al.
0-75 r . FLOTRAN

 

 

-0.25 '1

 

Figure 4. Horizontal velocity profile along the vertical centerline

 

0.2 . -——- Fluent
: ‘ ' - , o Hoetal.
0'15 ' . . FLOTRAN

 

 

 

0.1
0.05 ,
o .'
> -0.05
-0.1
-0.15
-0.2

-.21111111111111111111.
050 0.2 0.4 0.6 0.8 1

le
Figure 5. Vertical velocity profile along the horizontal centerline

 

207

A.6 Conclusion

The Fluent code has been validated against benchmark solutions of flows in a
driven cavity. The results predicted by FLUENT using the second-order

discretization scheme are in good agreement with the benchmark data.

Reference:

[1] Ho, OJ. and Lin, F. H, “Numerical Simulation of 3D lncompressible Flow,” lntemational
Journal for Numerical Methods in Fluids, Vol. 23, No. 10, pp. 1073-1084, 1996

[2] Karadag, A., Ganjoo, D. K., Ostergaard, D. F., and Shi, T. l-P., “ A Density-Based Method for
Computing lncompressible Multiphase Flows Part 2: Verification and Application,” AIAA paper

2000-0460, Jan., 2000

208

Appendix B

Publications Associated with this Work

Gu, X. and Shih, T.l-P., “Error Source and Error Appearance Locations for Error
Estimation and Solution-Adaptive Refinement in Complex Flows,” AIAA
Journal, 2002, manuscript in preparation.

Gu, X., Wu, H.-W., Schock, H.J., and Shih, T. I-P., “T wo-Equation Versus
Reynolds-Stress Modeling in Predicting Flow and Heat Transfer in a
Smooth U-Duct with and without Rotation,” 47th ASME Gas Turbine and
Aeroengine Technical Congress, 3-6 June 2002, Amsterdam, The
Netherlands.

Gu, X. and Shih, T.l-P., “Differentiating between Source and Location of Error for

Solution-Adaptive Mesh Refinement,” 15th AIAA Computational Fluid
Dynamics Conference, Anaheim, California, June 2001.

Gu, X., Hernandez, E., Chu, D., Sun, R., Keller, P., and Shih, T.l-P., “Grid-Quality

Measures for Error Estimation of CFD Solutions,” Computational Fluid and
Solid Mechanics, Vol. 1, Edited by K.J. Bathe, Elsevier, 2001, pp.843-846.

Gu, X., Shih, T.l-P., Schock, H.J., Hernandez, E., Sun, R., and Keller, F., “Grid-

Quality Measures for Structured and Unstructured Meshes,” 39‘“ AIAA
Aerospace Sciences Meeting and Exhibit, Reno, Nevada, January 2001.

Shih, T.l-P., Gu, X., and Chu, D., “Grid-Quality Measures and Error Estimates,”

Numerical Grid Generation in Computational Field Simulations, edited by
B.K. Soni, et al., Mississippi State University, 2000, pp. 799-808.

209

BIBLIOGRAPHY

210

BIBLIOGRAPHY

Adjerid, 8., J. Flaherty, and l. Babuska, "A posteriori error estimation for the finite
element method-of-Iines solution of parabolic problems,“ Mathematical
Models and Methods in Applied Sciences, Vol. 9, No. 2, 1999, pp. 261-
286.

Adjerid, S., Belguendouz, B., and Flaherty, J. E., “A Posteriori Finite Element
Error Estimation for Diffusion Problems,” SIAM Journal on Scientific
Computing, Volume 21, 1999, pp. 728-746.

Aftosmis, M.J., "Upwind Method for Simulation of Viscous Flow on Adaptively
Refined Meshes,” AIAA Journal, Vol. 32, No. 2, p. 268, February 1994.

Aftosmis, M., Gaitonde D., and Tavares T. 8., “Behavior of Linear
Reconstruction Techniques on Unstructured Meshes,” AIAA Journal, Vol.
33, No. 11, Nov. 1995, pp. 2038-2049.

Ainsworth, M. and Oden, J. T., “A Procedure for a-Posteriori Error
Estimation for h-p Finite Element Methods,” Computer Methods in Applied
Mechanics and Engineering, Vol. 101, 1992, pp. 73-96.

Ait-AIi-Yahia, D., Habashi, W. G., Tam, A., Valet, M-G., and Fortin, M., “A
Directionally Adapative Methodology Using an Edge Based Error Estimate
on Quadrilateral Grids,” International Journal for Numerical Methods in
Fluids, Vol. 23, 1996, pp. 673-690.

Anderson, D. A., Tannehill, J. C., and Pletcher, R. H., “Computational fluid
mechanics and heat transfer”, McGraw-Hill, c1984.

Babuska, l. and Aziz, A. K., “On the Angle Condition in the Finite Element
Method,” SIAM Journal Numerical Analysis, Vol. 13, No.2, April 1976, pp.
215-226.

Babuska, I. and Rheinboldt, W. C., “A-Posteriori Error Estimates for the
Finite Element Method,” Int. Journal Comp. Mesh. Engineering, Vol. 12,
1978. PP. 1597-1615.

Babuska, I. and Rheinboldt, W. 0., “Error Estimates for Adaptive Finite Element
Calculations,” SIAM Journal Numerical Analysis, Vol. 15(4), 1978, pp.
736-753.

Babuska, l. and Dorr, M. R., “Error Estimates for the Combined h and p versions

of the Finite Element Method,” Numerical Mathematics, Vol. 37, 1981, pp.
257-277.

211

Babuska, l. and Szabo, B., “On the rates of convergence of the finite
element method,” lntemational Journal Numerical Methods in Engineering,
Vol. 18, 1982, pp. 323-341.

Babuska, l., Zienkiewicz, O. C., Gago, J., and Oliveira, E.R.A., “Accuracy
Estimates and Adaptive Refinements in Finite Element Computations,”
AEARFEC, Wiley, 1986.

Babuska, l., Strouboulis, T., Upadhyay, C. S., and Gangaraj, S. K., “A Posteriori
Estimation and Adaptive Control of the Pollution Error in the h-Version of
the Finite Element Method,” International Journal for Numerical Methods
in Engineering, Vol.38, 1995, pp. 4207-4235.

Baker, T. J., “Prospects and Expectations for Unstructured Methods,” Workshop
on Surface Modeling, Grid Generation, and Related Issues in CFD
Solutions, NASA Lewis, May 1995.

Bank, R. E. and Weiser, A., “Some a-Posteriori Error Estimators for elliptic Partial
Differential Equations,” Math. Comput., Vol.44, 1985, pp. 283-301.

Bank, R. E. and Welfert, D. B., “Some A Posteriori Error Estimates for the Stokes
Equations: A Comparison,” Computer Methods in Applied Mechanics and
Engineerlng, Vol. 87, 1990, pp. 323-340.

Bank, R. E., Cao J., Mantel, B., and Periaux, J. , “Mesh Quality Control for
Industrial Navier—Stokes Problems via a Posteriori Error Estimates,”
Proceedings of ENUMATH-95 (A. Bespalov and Y. Kuznetsov, eds.),
VSP, Netherlands

Bansch, E. and Siebert, K. G., “A Posteriori Error Estimation for Nonlinear
Problems by Duality Techniques,”

Barth, T. J., “Numerical Aspects of Computing Viscous High Reynolds Flows on
Unstructured Meshes,” AIAA Paper 91 -0721, 29th Aerospace Sciences
Meeting, Jan. 1991, Reno Nevada

Barth, T.J. and MG. Larson, “ A Posteriori Error Estimation for Discontinuous
Galerkin Approximations of Hyperbolic Systems”, 1st International
Conference on DG Methods, May 24-26, 1999, NASA Report NAS-99-
010.

Berger, M. J. and Oliger, J., “Adaptive Mesh Refinement for Hyperbolic Partial
Differential Equations,” Journal of Computational Physics, Vol. 53, 1984,
pp. 484-512.

Berger, M. J. and Colella, P., “Local Adaptive Mesh Refinement for Shock

212

Hydrodynamics,” Journal of Computational Physics, Vol. 82, 1989, pp. 64-
84.

Berzins, M., “A Solution-Based Triangular and Tetrahedral Mesh Quality
Indicator,” SIAM J. of Scientific Computing, 1998; Vol. 19, No. 6, 1998, pp.
2051-2060.

Berzins, M., “Mesh Quality - Geometry, Error Estimates or Both?, ” Engineering
and Computers, Vol. 15, 1999, pp. 236-247.

Berzins, M., “Solution-Based Mesh Quality Indicators For Triangular and
Tetrahedral Meshes,” International Journal of Computational Geometry 81
Applications, Vol. 10(3), 2000, pp. 333-346.

Bhatia, R. P. and Lawrence, K. L., “T wo-Dimensional Finite Element Mesh
Generation Based on Stripwise Automatic Triangulation,” Computers 8
Structures, Vol. 36(2), 1990, pp. 209-319.

Biswas, R. and Strawn, R., 'A New Procedure for Dynamic Adaption of Three-
Dimensional Unstructured Grids,” AIAA-93-0672, 31st AIAA Aerospace
Sciences Meeting, Reno, NV, January 1993.

Bowyer A., “Computing Dirichlet tessalations,” The Computer Journal, Vol. 24,
1981, pp. 162-166.

Brackbill, J.U., “An Adaptive Grid with Directional Control,“ Journal of
Computational Physics, Vol. 108, 1993, pp. 38.

Brandt, A., “Multi-Ievel adaptive solutions to boundary-value problems,” Math.
Comp., 31, 1977, pp.333-390.

Carey, G. F., “Computational Grids Generation, Adaptation, and Solution
Strategies,” Taylor & Francis, Washington, DC, 1997.

Caruso, 8., “Adaptive Grid Techniques for Fluid Flow Problems.” Ph.D. thesis,
Thermosciences Division, Department of Mechanical Engineering,
Standford University, 1985

Celik, I. and Zhang, W. M., “Application of Richardson Extrapolation to Some
Simple Turbulent Flow Calculations,” Quantification of Uncertainty in
Computational Fluid Dynamics, FED-Vol.158, ASME, 1993, pp.29-38.

Celik, I. and Zhang, W.M., “Calculation of Numerical Uncertainty Using
Richardson Extrapolation: Application to Some Simple Turbulent Flow
Calculations,” ASME, Journal of Fluids Engineering, Vol. 117, September,
1995, pp. 439-445.

213

 

 

Chang, S. M. and Haworth, D. C., “Kinetic-Energy-Balance Based Solution-
Adaptive Mesh Refinement,” Quantification of Uncertainty in
Computational Fluid Dynamics, FED-Vol. 213, ASME, 1995, pp. 7-12.

Chen, W. L., Lien, F. S., and Leschziner, M. A., “ Local Mesh Refinement within a
multi-block structured-grid scheme for general flows,” Computer Methods
in Applied Mechanics and Engineering, Vol. 144, 1997, pp. 327-369

Coelho, P., Pereira, J. C. F., and Carvalho, M. G., “Calculation of Laminar
Recirculating Flows using a Local Non-Staggered Grid Refinement
System,” lntemational Journal of Numerical Methods in Fluids, Vol. 12(6),
1991 , pp. 535-557.

Coleman, H. W. and Stern, F., “Uncertainties and CFD Code Validation,” ASME,
Journal of Fluids Engineering, Vol. 119, December 1997, pp. 795-803.

Connell, S. D. and Holmes, D. G., “A 3D Unstructured Adaptive Multigrid Scheme
for the Euler Equations,” AIAA Journal, Vol. 32, 1994, pp. 1626-1632

Connell, S.D., Sober, J.S. and Lamson, S.H., “Grid Generation and Surface
Modeling for CFD,” Proceedings of the Surface Modeling, Grid
Generation, and Related Issues in Computational Fluid Dynamics
Workshop, NASA Conference Publication 3291, p. 29, NASA Lewis
Research Center, Cleveland, OH, May 1995.

Dandekar, H. W., Hlavacek, V., and Degreve, J., “An Explicit SD Finite-Volume
Method for Simulation of Reactive Flows Using a Hybrid Moving Adaptive-
Grid,” Numerical Heat Transfer, part B, Vol. 24(1), 1993, pp. 1-29.

Dannelongue, H. H. AND Tanguy, P.A., “ Three-Dimensional Adaptive Finite
Element Computations and Applications to Non-Newtonian Fluids,”
International Journal for Numerical Methods in Fluids, Vol. 13, 1991, pp.
145-165.

Davis, S. F. and Flaherty, J. E., “An Adaptive Finite Element Method for Initial-
Boundary Value Problems for Partial Differential Equations,” SIAM journal
on scientific and statistical computing, Vol. 3, 1982, pp. 6-27.

Demuren, A. O. and Wilson, R. V., “Estimating Uncertainty in Computations of
Two —Dimensional Separated Flows,” ASME, Journal of Fluids
Engineering, Vol. 116, June, 1994, pp. 216-220.

Diaz, A. R., Kikuchi, N., and Taylor, J. E., “A Method of Grid Optimization for
Finite Element Methods,” Computer Methods in Applied Mechanics and
Engineering, Vol. 41, 1983, pp. 29-45.

Dolce, J. S., Silva D., and Shih, T. l-P, “Applications of Variational Techniques for

214

Constructing Streching Functions,” 10th Brazilian Congress of Mechanical
Engineering, 1989.

Dwyer, H. A., “Grid Adaptation for Problems in Fluid Dynamics,” AIAA Journal,
Vol. 22(12), Dec. 1984, pp. 1705-1712.

Ekeland, l. and Teman, R., “Convex analysis and variational problems,”
Amsterdam : North-Holland Pub. Co., 1976

Eriksson, K., “ Improved Accuracy by Adapted Mesh-Refinements in the Finite
Element Method,” Mathematics of Computation, Vol. 44, 1985, pp. 321-
343.

Eriksson, K. and Johnson, C., “An Adaptive Finite Element Method for Linear
Elliptic Problems,” Mathematics of Computation, Vol. 50, 1988, pp. 361-
383.

Eriksson, K.; Johnson, 0.: Adaptive finite element methods for parabolic
problems I: A linear model problem, SIAM Journal on Numerical Analysis,
Vol. 28, 1991, pp. 43-47.

Ferziger, J. H., “Estimation and Reduction of Numerical Error,” Quantification of
Uncertainty in Computational Fluid Dynamics, FED-Vol.158, ASME, 1993,

pp.1-7.

Field, D. A., “Qualitative Measures for Initial Meshes,” International Journal for
Numerical Methods in Engineering, Vol.47, 2000, pp. 887-906.

Fluent User's Guide Volume, Fluent Inc.

Freitag, L. A. and Patrick M. K., “T etrahedral Element Shape Optimization via the
Jacobian Determinant and Condition Number,” Proceedings, 8th
International Meshing Roundtable, South Lake Tahoe, CA, USA,
pp.247-258, October 1999.

Frink, N. T., Pirzadeh, S., and Parikh, P., “An Unstructured-Grid Software System
for Solving Complex Aerodynamic Problems,” NASA CP-3291, pp. 289-
308, May 9-11, 1995.

Gitlin, CS. and Johnson, C.R., “Meshview: A tool for exploring 3d unstructured
tetrahedral meshes,” in 5th International Meshing Roundtable, 1996, pp.
333--345.

“Guide for the Verification and Validation of Computational Fluid Dynamics
Simulations,” AIAA G-077-1998.

215

Habashi, W. G., Dompierre, J., Bourgault, Y., Fortin, M., and Vallet, M. G.,
“Certifiable Computational Fluid Dynamics Through Mesh Optimization,”
AIAA Journal, Vol. 36, No. 5, 1998, pp.703-711.

Hawken, D. F., Gottlieb, J. J., and Hansen, J. S., “Review of some Adaptive
Node-Movement Techniques in Finite-Element and Finite-Difference
Solutions of Partial Differential Equation,” Journal of Computational
Physics, Vol. 95, 1991, pp. 254-302.

Haworth, D. 0., El Tahry, S. H., and Huebler, M. S., “A Global Approach to Error
Estimation and Physical Diagnostics in Multidimensinal Fluid Dynamics,”
lntemational Journal of Numerical Methods in Fluids, Vol. 17(1), 1993, pp.
75-97.

Hoffman, J. 0., “Relationship between the Truncation Errors of Centered Finite-
Difference Approximations on Uniform and Nonuniform Meshes,” Journal
of Computational Physics, Vol.46, 1982, pp. 469-474.

Ilinca, A., Camarero, R., Trepanier, J. Y., and Reggio, M., “Error estimator and
adaptive moving grids for finite volumes schemes,” AIAA Journal, vol. 33,
No. 11, Nov. 1995. PP. 2058-2065.

Jacquotte, O.-P., “Grid Optimization Methods for Quality Improvement and
Adaptation,” in Handbook of Grid Generation, edited by JP. Thompson,
B.K. Soni, and NR Weatherill, CRC Press, Boca Raton, 1999.

Jasak, H. and Gosman, A.D., “Local Problem Error Estimate in Finite Volume
Discretisation,” Submitted to Comp. Meth. Appl. Mech. Engineering, 1998.

Jasak, H. and Gosman, A.D., “Residual error estimate for the Finite Volume
Method,” Submitted to International Journal for Numerical Methods in
Fluids 1999.

Jeng, Y. M, Lin S. Y., Lee Z-S, and kung B-Y, “On the Smoothing of the
Multiple One-Dimensional Adaptive Grid Procedure,” 37th AIAA Aerospace
Sciences Meeting and Exhibit, Jan. 11-14, 1999/Reno, NV.

Johnson, C., Rannacher, R. and Boman, M., “Numerics and Hydrodynamics
Stability Theory: Towards Error Control in CFD,” SIAM Journal on
Numerical Analysis, Vol. 32, 1995, pp. 1058-1079.

Johnson, C. and Szepessy, A., “Adaptive Finite Element Methods for
Conservation Laws Based on a Posteriori Error Estimates,”
Communications on Pure and Applied Mathematics, Vol. 48, No. 3, 1995,
pp. 199-234.

Kallinderis, Y. and Vijayan, P., ”Adaptive Refinement-Coarsening Scheme for

216

Three-Dimensional Unstructured Meshes,” AIAA Journal, Vol. 31, No. 8, p.
1440, August 1993.

Kelly, D. W., “The self-equilibration of residuals and “upper-bound” error
estimates for the finite element method,” Chapter 8, Accuracy Estimates
and Adaptive refinements in Finite element Computations, Edit by
Babuska, 1., et al., J Wiley 8. Sons, Chichester, 1986.

Khayat, R. E. and Genouvrier, D., “An Adaptive Lagrangian Boundary Element
Approach for Three-Dimensional Transient Free-Surface Stokes Flow as
Applied to Extrusion, Therrnoforrning, and Rheometry,” lntemational
Journal for Numerical Methods in Fluids, Vol. 36,2001, pp. 1-33.

Knupp, P., ”Mesh Generation Using Vector-Fields,” Journal of Computational
Physics, Vol. 119, 1995, pp. 142.

Krizek, M., “On the Maximum Angle Condition for Linear Tetrahedral Elements,”
SIAM Journal of Numerical Analysis, Vol. 29, 1992, pp. 513-520.

Lee, D. and Tsuei, Y. M., “A Formula for Estimation of Truncation Errors of
Convection Terms in a Curvilinear Coordinate System,” Journal of
Computational Physics, Vol. 98, 1992, pp. 90-100.

Lepape, M.-C. and Jacquotte, O.-P., ”Acceleration Techniques for a 3D
Structured Multiblock Mesh Optimization Method,“ Numerical Grid
Generation in Computational Field Simulation and Related Fields, N.P.
Weatherill, P.R. Eiseman, J. Hauser and J.F. Thompson (Eds.), p. 95,
Proceedings of the 4th lntemational Grid Conference, Pineridge Press
Limited: Swansea Wales (UK), 1994.

Lisa, D., “Evaporation: a Technique for Visualizing Mesh Quality,” Proceedings,

8th International Meshing Roundtable, South Lake Tahoe, CA, USA,
pp.259-265, October 1999.

Liu, A. and Joe, B., “Relationship between Tetrahedron Shape Measures,” BIT,
Vol. 34, 1994, pp. 268-287.

Liseikin, V. D., “Grid Generation Methods,” Springer, 1999.

Lo, S. H., “Optimization of Tetrahedral Meshes based on Element Shape
Measures,” Computers 8. Structures, Vol. 63, No. 5, 1997, pp. 951-961.

Lohner, R. and Baum J. D., “Adaptive H-refinement on 3D Unstructured Grids for

Transient Problems,” International Journal for Numerical Methods in
Fluids, Vol. 14, 1992, PP. 1407-1419.

217

Marcum, D. L., “Generation of unstructured grids for viscous flow applications,”
AIAA Paper No. 95-0212, 1995

Mavriplis, D. J., “Three-dimensional multigrid for the Euler Equations,” AIAA
Journal, Vol. 30, 1992, pp. 1753-1761

Mavriplis, D. J., “Unstructured Mesh Generation and Adaptivity,” NASA
Constractor Report 195069, ICASE Report No. 95-26.

Mavriplis, D. J., “Adaptive Meshing Techniques for Viscous Flow Calculations on
Mixed Element Unstructured Meshes,” International Journal for Numerical
Methods in Fluids, Vol. 34, 2000, pp. 93-111.

McGuirk, J.J. and Rodi, W., “A depth-averaged mathematical model for the near
field of side discharges into open channel flow,” Journal of Fluid
Mechanics, Vol. 86, 1978, pp.761.

McRae, D. S. and Laflin, K.R., “Dynamic Grid Adaptation and Grid Quality,” in
Handbook of Grid Generation, edited by J.F. Thompson, B.K. Soni, and
NP. Weatherill, CRC Press, Boca Raton, 1999.

Mehta, U. B., “Some Aspects of Uncertainty in Computational Fluid Dynamics
Results,” ASME, Journal of Fluids Engineering, Vol. 113, December,
1991 , pp. 538-543.

Mehta, U.B., ”Guide to Credible Computer Simulations of Fluid Flows,“ AIAA
Journal of Propulsion and Power, Vol. 12, No. 5, September-October
1996, pp. 940-948. (Also AIAA Paper 95-2225).

Merriam, M. L., “An Efficient Advancing Front Algorithm for Delaunay
Triangulation,” AIAA, 29th Aerospace Science Meeting, Jan.7—10,
1991/Reno, NV.

Mitchell, A. R., Phillips, G., and Wachspress, E., “Forbidden Shapes in the Finite
Element Method,” Journal of the Institute for Mathematics and
Applications, Vol. 8, 1971, pp. 260-269.

Muzaferija, S. and Gosman, D., “Finite-volume CFD Procedure and Adaptive
Error Control Strategy for Grids of Arbitrary Topology,” Journal of Comp.
Physics, Vol. 138(2), 1997, pp. 766-787.

Oberkampf, W. L. and Blottner, F. 6., “Issues in Computational Fluid Dynamics
Code Verification and Validation,” AIAA Journal, Vol. 36, No. 5, 1998, pp.
687-695.

Oberkampf, W. L. and Trucano, T. G., “Validation Methodology in Computational

218

Fluid Dynamics,” AIAA Paper 00-2549, Fluids 2000, Jun. 2000, Denver,
CO.

Oden, J. T., Demkowicz, L., Strouboulis, T., and Devloo, P., “Adaptive Methods
for Problems in Solid and Fluid Mechanics,” Chapter 14, Accuracy
Estimates and Adaptive refinements in Finite element Computations, Edit
by Babuska, I., et al., J Wiley 8 Sons, Chichester, 1986.

Oden, J. T., Demkowicz, L., Rachowicz, L., and Hardy, O., “Toward a Universal
h-p Adaptive Finite Element Strategy. Part I: Constrained Approximation
and Data Structure,” Computer Methods in Applied Mechanics and
Engineering, Vol. 77, 1989, pp. 79-112.

Oden, J. T., Demkowicz, L., Rachowicz, W., and Westermann, T., “Toward a
Universal h-p Adaptive Finite Element Strategy. Part II: A Posteriori Error
Estimation,” Computer Methods in Applied Mechanics and Engineering,
Vol.77, 1989, pp. 113-180.

Oden, J. T., Wu, W., and Ainsworth, M., “An a-Posteriori Error Estimate for Finite
Element Approximations of the Navier Stokes Equations,” Computer
Methods in Applied Mechanics and Engineering, Vol. 111, 1994, pp. 185-
202.

Owen, S. J., “Constrained Triangulation: Application to Hex-Dominant Mesh
Generation,” Proceedings, 8th International Meshing Roundtable, South
Lake Tahoe, CA, U.S.A., pp.31-41, October 1999.

Parthasarathy, V.N. and Kodiyalam S., “A Constrained Optimization Approach to
Finite Element Mesh Smoothing,” Finite Elements in Analysis and Design,
Vol. 9, 1991, pp. 309-320.

Parthasarathy, V. N., et al, “A comparison of tetrahedron quality measures,”
Finite Elements in Analysis and Design, Vol. 15, 1993, pp. 255-261.

Pelletier, D. and lgnat, L., “On the Accuracy of the Grid Convergence Index and
the Zhu-Zienkiewicz Error Estimator,” Quantification of Uncertainty in
Computational Fluid Dynamics, FED-Vol.213, ASME, 1995, pp.1-7.

Peraire, J., Vahdati, M., Morgan, K., and Zienkiewicz, O.C., “Adaptive
Remeshing for Compressible Flow Calculations,” Journal of
Computational Physics, Vol.22, 1976, pp. 131-149.

Pirzadeh, S., “Structured Background Grids for Generation of Unstructured Grids
by Advancing-Front Method,” AIAA Journal, Vol. 31, No. 2, Feb. 1993,
pp.257-265.

Pirzadeh, 8., “Progress Toward a User-Oriented Unstructured Viscous Grid

219

Generator,” 34‘“ Aerospace Sciences Meeting & Exhibit, Jan.15—18,
1996/Reno, NV.

Pirzadeh, S., “Three-Dimensinal Unstructured Viscous Grids by the Advancing-
Layers Method,” AIAA Journal, Vol.34, No. 1, Jan. 1996, pp. 43-49.

Pyun S. B. and Y00 H. S., “An a-Posteriori Error Estimate and Mesh Generation
Algorithm using Boundary Node Control,” Numerical Grid Generation in
Computational Field Simulations, edited by B.K. Soni, J.F. Thompson, J.
Hauser, and P. R. Eiseman, ISGG, Mississippi State University, 2000, pp.
819-825.

Richardson, L. F. and Gaunt, J. A., “T he deferred approach to the limit, “Trans.
Roy. Soc. London, Ser. A, vol. 226,229-361, 1927

Richter, R., 'Discontinuity Capturing Using Auto-Adaptive Finite Elements,“
Numerical Grid Generation in Computational Field Simulation and Related
Fields, N.P. Weatherill, P.R. Eiseman, J. Hauser and J.F. Thompson
(Eds.), p. 627, Proceedings of the 4th lntemational Grid Conference,
Pineridge Press Limited: Swansea Wales (UK), 1994.

Rippa, 8., “Long and Thin Triangles Can Be Good for Linear Interpolation,” SIAM
Journal of Numerical Analysis, Vol. 29, 1992, pp. 257-270.

Roache, P. J., “Perspective: A Method for Uniform Reporting of Grid Refinement
Studies,” Journal of Fluids Engineering, Vol. 116, 1994, pp. 405-413.

Roache, P. J., “Verification and Validation in Computational Science and
Engineering,” Hermosa Publishers, Albuquerque, New Mexico, 1998.

Robinson, J., “Some New Distortion Measures for Quadrilaterals,” Finite
Elements in Analysis and Design, Vol. 3, 1987, pp. 183-197.

Robinson, J., “Distortion Measures for Quadrilaterals with curved boundaries,”
Finite Elements in Analysis and Design, Vol. 4, 1988, pp. 115-131.

Robinson, J., “Quadrilateral and Hexahedron Shape Parameters,” Finite
Elements in Analysis and Design, Vol. 16, 1994, pp. 43-52.

Schonfeld, T., "Automatic Local Grid Refinement for Vortical Flow Computations
Around Delta Wings,” Numerical Grid Generation in Computational Field
Simulation and Related Fields, N.P. Weatherill, P.R. Eiseman, J. Hauser
and J.F. Thompson (Eds.), p. 589, Proceedings of the 4th lntemational
Grid Conference, Pineridge Press Limited: Swansea Wales (UK), 1994.

Shih, T. l-P., Gu, X., and Chu, D., “Grid-Quality Measures and Error Estimates,”

220

Numerical Grid Generation in Computational Field Simulations, edited by
B.K. Soni, J.F. Thompson, J. Hauser, and P. R. Eiseman, ISGG,
Mississippi State University, 2000, pp. 799-808.

Simpson, R. B., “Anisotropic Mesh Transformations and Optimal Error Control,”
Applied Numerical Mathematics, Vol. 14, 1994, pp. 183-198.

Sonar, T., Hannemann, V., and Hempel, D., “Dynamic Adaptivity and Residual
Control in Unsteady Compressible Flow Computation,” Mathematical and
Computer Modelling, Vol. 20(10-1 1), 1994, pp. 201-213.

Sorokin, A.M., ”Structured and Unstructured Grid Generation: An Evolutionary
Approach," Numerical Grid Generation in Computational Field Simulation
and Related Fields, N.P. Weatherill, P.R. Eiseman, J. Hauser and J.F.
Thompson (Eds.), p. 287 Proceedings of the 4th lntemational Grid
Conference, Pineridge Press Limited: Swansea Wales (UK), 1994.

Spragle, G. and Smith W. A., “Hanging Node Solution Adaptation on Hybrid
Grids,” Proceedings of the 5 lntemational Conference on Numerical Grid
Generation in Computational Field Simulations, Starkville, M. 3., Sony, B.
K. et al., Mississippi State University, 1996, pp. 1221-1230.

Steinbrenner, JP. and Noack, R.W., “Three-Dimensional Hybrid Grid Generation
Using Advancing Front Techniques,” NASA CP-3291, Surface Modeling,
Grid Generation, and Related Issues in CFD Solutions, NASA Lewis
Research Center, Cleveland, OH, May 1995.

Straalen, B. V., “A Posteriori Error Estimation for Finite Volume Simulations of
Fluid Flows,” University of Waterloo, Canada, 1995, Master Thesis.

Strouboulis, T. and Hague, K. A., “Recent Experiences with Error Estimation and
Adaptivity IzReview of Error Estimators for Scalar Elliptic Problems,”
Computer Methods in Applied Mechanics and Engineering, Vol. 97, 1992,
pp. 399-436.

Suli, E., “A Posteriori Error Analysis and Adaptivity for Finite Element
Approximations of Hyperbolic Problems,” In An Introduction to Recent
Developments in Theory and Numerics for Conservation Laws, Vol. 5 of
LNCSE, pp. 112-194, Springer-Venag, Heidelberg, 1998.

Synge, J. L., “The hypercircle in mathematical physics; a method for the
approximate solution of boundary value problems,” Cambridge [Eng.]
University Press, 1957.

Szmelter, J., Marchant, M., Evans, A. and Weatherill, N.P., “Solution of the Two-

221

Dimensional Compressible Navier-Stokes Equations on Embedded
Structured Multiblock Meshes,“ Numerical Grid Generation in
Computational Field Simulation and Related Fields, A.S. Arcilla, J.
Hauser, P.R. Eiseman, and J.F. Thompson (Eds.), p. 287, Proceedings of
the 3rd International Grid Conference, North Holland, Barcelona, Spain,
June 1991.

Thacker, W. C., “A brief Review of Techniques for Generating Irregular
Computational Grids,” International Journal for Numerical Methods in
Engineering, Vol. 15, 1980, pp. 1335-1341.

Teigland, R. and Eliassen, I. K., “A Multiblock/Multilevel Mesh Refinement
Procedure for CFD computations,” International Journal for Numerical
Methods in Fluids, Vol. 36, 2001, pp. 519-538.

Thompson, M. C. and Ferziger, J. H., “An Adaptive Multigrid Technique for the
lncompressible Navier-Stokes Equations,” Journal of Computational
Physics, Vol. 82, 1989, pp. 94-121.

Thompson, J. F., Soni, B. K., and Weatherill, N. P., “Handbook of Grid
Generation,” CRC Press, 1998.

Venditti, D. A., and Darmofal, D. L., “A Grid Adaptive Methodology for Functional
Outputs of Compressible Flow Simulations,” AIAA Paper 01-2659, 15th
Computational Fluid Dynamics Conference, Jun. 2001, Anaheim, CA.

Verfurth, R., “A Posteriori Error Estimates for Nonlinear Problems. Finite Element
Discretization of Elliptic Equations,” Mathematics of Computation, Vol. 62,
No. 206, 1994, pp. 445-475.

Verfurth, R., “A Posteriori Error Estimation and Adaptive Mesh-Refinement
Techniques,” Journal of Computational and Applied Mathematics, Vol. 50,
1994. PP. 67-83.

Vilsmeier, R. and Hanel, D., “Adaptive Methods on Unstructured Grids for Euler
and Navier-Stokes Equations,” Computers and Fluids, Vol. 22(4-5), 1993,
pp. 485-499.

Vinokur, M., “On One-Dimensional Stretching Functions for Finite-Difference
Calculations,” Journal of Computational Physics, Vol. 50, 1983, pp. 215-
234.

Wa, K. and Chen, Z., “A Simple and Effective Mesh Quality Metric for Hexahedral

and Wedge Elements,” Proceedings, 9th International Meshing
Roundtable, Sandia National Laboratories, pp.325-333, October 2000.

222

Warsi, Z.U.A. and Thompson, J.F., "Application of Variational Methods in The
Fixed and Adaptive Grid Generation," Computers 81 Mathematical
Applications, Vol. 19, 1990, pp. 31.

Webster, B.E., Shephard, M.S., Rusak, Z. and Flaherty, J.E., “Automated
Adaptive Time-Discontinuous Finite Element Method for Unsteady
Compressible Airfoil Aerodynamics,“ AIAA-93-0339, 3tst AIAA Aerospace
Sciences Meeting, Reno, NV, January 1993, AIAA Journal, Vol. 32, No. 4,
p. 748, April 1994.

Wilcox, D. C., “Turbulence Modeling for CFD,” DCW Industries, Inc., La Canada,
California, 1993.

Wilson, F. W., Goodrich, R. K., and Spratte, W., “Lawson‘s Triangulation is
Nearly Optimal for Controlling Error Bounds,” SIAM Journal of Numerical
Analysis, Vol. 27(1), 1990, pp. 190-197.

Yuan, K-Y., Huang, Y-S., and Plan, H. H.,“Inverse Mapping and Distortion
Measures for Quadrilaterals with curved boundaries,” lntemational Journal
for Numerical Methods in Engineering, Vol. 37, 1994, pp. 861-875.

Zhu, J. 2. and Zienkiewicz, O. 0., “Adaptive Techniques in the Finite Element
Method,” Communications in Applied Numerical Methods, vol.4, 1988, pp.
197-204.

Zienkiewicz, O. C. and Taylor, R. L., “The finite element method,” New York,
McGraw-Hill, 1989.

Zienkiewicz, DC. and Zhu, J.Z., “A Simple Error Estimator and Adaptive
Procedure for Practical Engineering Analysis,” lntemational Journal for
Numerical Methods in Engineering, Vol.24, 1987, pp.337-357.

Zienkiewicz, O. C. and Zhu, J. 2., “The Superconvengent Patch Recovery and a

Posteriori Error Estimates. Part 1 : The Recovery Technique,”
International Journal for Numerical Methods in Engineering, Vol.33, 1992,
pp. 1331-1364.

Zienkiewicz, O. C. and Zhu, J. 2., “The Superconvergent Patch Recovery and a
Posteriori Error Estimates. Part 2 : Error Estimates and Adaptivity,”
International Journal for Numerical Methods in Engineering, Vol.33, 1992,
pp. 1365-1382.

Zlamal, M., “On the Finite Element Method,” Numerical Math., Vol. 12, 1968, pp.
394-409.

223