ASSESSMENT OF TAYLOR-SERIES MODEL FOR ESTIMATION OF VORTEX
FLOWS USING SURFACE STRESS MEASUREMENTS
By
Gaurang Shrikhande

A THESIS
Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of
MASTER OF SCIENCE
Mechanical Engineering
2012

ABSTRACT
ASSESSMENT OF TAYLOR-SERIES MODEL FOR ESTIMATION OF
VORTEX FLOWS USING SURFACE STRESS MEASUREMENTS
By
Gaurang Shrikhande
This research is motivated by the development of ground-based sensor arrays for detection
of wake vortices behind an aircraft that would enable take off and landing at faster rates at
commercial airports in the Unites States. It is envisioned that multiple sensor-array units
installed along the sides of runways would accurately detect the strength and the location
of wake vortices in real time, enabling efficient and safe management of airplane take-offs
and landings. The focus of this research is on examining the feasibility of a method that can
estimate the details of the velocity field of wake vortices from distributed measurements of
the wall shear stress and pressure. The method relies on obtaining a Taylor-series expansion
of the velocity field starting from the ground. This idea is based on the ability to derive the
coefficients of the Taylor-series expansion up to an arbitrary order using only wall-shear-stress
and -pressure information coupled with recursive relations that are derived from the NavierStokes equations without any knowledge of the far-field boundary condition. In the current
research, we examined the possibility that this concept can be extended to estimate the wakevortex velocity field from near-wall-velocity and -pressure information. The analysis show
that though the Taylor-series model works in principle, the convergence of the series is very
slow, enabling accurate estimation of only the flow within the boundary layer. Furthermore,
the increase in size of the accurately estimated domain above the wall with increase in the
series order reaches a ’break even’ point at which further increase in the order is offset by
inaccuracies in calculating high-order derivatives.

Copyright by
GAURANG SHRIKHANDE
2012

To
My Family

iv

ACKNOWLEDGMENTS

I would like to express my sincere gratitude to Dr. Naguib for giving me the opportunity
to work at Flow Physics and Control Laboratory (FPaCL), for his patience and continuous
guidance during my gradute studies, and for helping me in some of my difficult times both
academically and personally. I am thankful to Dr. Dominique Fourguette for her guidance
and expertise. I am thankful to Michigan Aerospace Corporation for generously supporting
the initial phase of the project under contract # F1373-02012010. I would like to express my
special thanks to Dr. M. M. Koochesfahani for giving us access to the experimental data,
providing a discussion forum through group meetings and giving feedback on research, as
and when required.
I owe my deepest gratitude to my grand-parents late Shalini, late Digambar, late Arvind,
Indira, and Madhuri, my parents Swati and Sanjay, my sister Gayatri, and aunt Kalpana for
their unconditional love, support, encouragement, and for being my strength.
Special thanks to Rohit Nehe, Kyle Bade, and Malek Al-Aweni for making my time
at ERC eventful. I would like to thank my closest friends from elementary school and
beyond - Sonali, Suchita, Abdul Husein, Bhaargav, and Ankur. Sincere thanks to my friend
Krupesh for his support and encouragement at all times. Lastly, thanks to International
Students Association (ISA) and Office for International Students and Scholars (OISS) for
truely internationalizing my life.

v

TABLE OF CONTENTS

List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
List of Figures

ix

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

x

1 Introduction
1.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.2 Approach . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.2.1 Step I - Construction and Validation of the Taylor-series model. . . .
1.2.2 Step II - Numerical Computation to evaluate the applicability of the
Taylor-series model to estimate the flow field of a vortex-pair above a
wall. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.2.3 Step III - Study of Cartesian Taylor Series model using simulations. .

1
1
4
4

2 Taylor-series Model For Estimation of Incompressible Flow Fields Using
Wall-Stress Information
2.1 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2 Concept in 2D . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2.1 Step I. Expanding the stream function (ψ) in a two-dimensional, infinite Taylor-series form . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2.2 Step II. Expressing the velocity components in series form by differentiating the series for the stream function given by Equation 2.1 . . . .
2.2.3 Step III. Obtaining series expressions for partial derivatives of the
velocity components in the momentum equation (the components of
which in the x and y directions are given by Equations 2.20 and 2.24
respectively) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2.4 Step IV. Obtaining expression relating the lowest-order Taylor series
coefficients a1j to the wall shear stress using the constitutive equation
for Newtonian fluids . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2.5 Step V. Obtaining expression relating the next-order Taylor series coefficients a2j to the wall pressure gradient using the momentum equation
in the streamwise direction . . . . . . . . . . . . . . . . . . . . . . . .
2.2.6 Step VI. Obtaining expressions for higher-order Taylor series coefficients in terms of a1j and a2j and their time derivatives by successive
differentiation (in the wall-normal direction) of the 2D momentum
equations and evaluating the result at the wall . . . . . . . . . . . . .
2.3 Demonstration In Laminar Flows With Known Analytical Solution . . . . .
2.3.1 Blasius Boundary Layer . . . . . . . . . . . . . . . . . . . . . . . . .
vi

5
6

7
7
8
8
9

10

11

12

14
17
19

2.4

2.3.2 Falkner-Skan Boundary Layer (m = 1) . . . . . . . . . . . . . . . . .
2.3.3 Stokes Oscillatting Stream Problem . . . . . . . . . . . . . . . . . . .
Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3 Validation of The Computational Approach
3.1 Background . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2 Description of The Experimental Data . . . . . . . . . . . .
3.2.1 Credits . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2.2 Introduction . . . . . . . . . . . . . . . . . . . . . . .
3.2.3 Experimental Facility . . . . . . . . . . . . . . . . . .
3.2.4 Observations . . . . . . . . . . . . . . . . . . . . . .
3.3 Initial Conditions For The Simulation Model . . . . . . . . .
3.3.1 Gaussian Vorticity Distribution . . . . . . . . . . . .
3.3.1.1 Two-Parameter Fit . . . . . . . . . . . . . .
3.3.1.2 Three-Parameter Fit . . . . . . . . . . . . .
3.3.2 Comparison of Numerical Results With Experiments
3.3.2.1 Effect of computational time step . . . . . .
3.3.2.2 Effect of number of iterations per time step
3.3.2.3 Effect of grid size . . . . . . . . . . . . . . .
3.4 Summary of Fluent Parameters . . . . . . . . . . . . . . . .
3.4.0.4 General Attributes: . . . . . . . . . . . . . .
3.4.0.5 Solution Methods: . . . . . . . . . . . . . .
3.4.0.6 Solution Controls: . . . . . . . . . . . . . .

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

4 Simulation of A Pair of Line Vortices Impinging On A Solid Wall And
Estimation of Velocities Using The Cartesian Taylor Series Model
4.1 Geometry of the Computational Model . . . . . . . . . . . . . . . . . . . . .
4.2 Simulations Using ANSYS R -Fluent . . . . . . . . . . . . . . . . . . . . . . .
4.3 Simulation Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.4 Estimations Using Cartesian Taylor-series Model . . . . . . . . . . . . . . . .
4.5 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.5.1 Comparison of estimations based on different grid resolutions . . . . .
4.5.2 The same grid resolution (600 × 600) but two different estimation
locations: x = 25 mm and x = 30 mm . . . . . . . . . . . . . . . . .

22
26
29
32
33
34
34
34
35
35
39
39
43
43
44
54
58
58
60
60
62
62

63
64
65
68
76
77
77
83

5 Summary, Conslusions and Recommendations
5.1 Summary and Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5.2 Recommendations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

86
86
88

APPENDICES

90
vii

A Taylor-series Model for Axisymmetric Flow
91
A.1 Derivation of the Taylor-series Model for Axisymmetric Flow . . . . . . . . . 91
A.1.1 Step I. Expanding the stream function (ψ) in a two-dimensional, infinite Taylor-series form . . . . . . . . . . . . . . . . . . . . . . . . . . 91
A.1.2 Step II. Expressing the velocity components in a series form by differentiating the series for the stream function given by Equation A.1 . . 92
A.1.3 Step III. Obtaining series expressions for partial derivatives of the velocity components in the momentum equations A.20 . . . . . . . . . . 93
A.1.4 Step IV. Obtaining expression relating the lowest-order Taylor series
coefficients a0j to the wall shear stress using constitutive equation for
Newtonian fluids . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94
A.1.5 Step V. Obtaining expression relating the next-order Taylor series coefficients a1j to the wall pressure gradient using the momentum equation
in the streamwise direction . . . . . . . . . . . . . . . . . . . . . . . . 95
A.1.6 Step VI. Obtaining expressions for higher-order Taylor series coefficients in terms of a0j and a1j and their time derivatives by successive
differentiation (in the wall-normal direction) of the 2D momentum
equations and evaluating the result at the wall . . . . . . . . . . . . . 97
A.2 Demonstration In Axisymmetric Laminar Flows With Known Analytical Solution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99
A.2.1 Axisymmetric Flow . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100
BIBLIOGRAPHY

. . . . . . . . . . . . . . . . . . . . . . . . . . .

viii

105

LIST OF TABLES

Analytical flows employed in validation of MATLAB R code used to
derive Taylor series coefficients . . . . . . . . . . . . . . . . . . . . .

18

Table 3.1

Initial Vortex parameters obtained from the Gaussian fit

. . . . . .

44

Table 3.2

Computational cases and associated initial vortex parameters. Listed
experimental values are those given by Gendrich et al. [3] . . . . . .

54

Comparison of the y extent from the wall over which the u-velocity
estimation has an error less than ±0.1% at a location directly beneath
the center of the primary vortex (x = 25 mm) for different grid sizes

82

Comparison of the y extent from the wall over which the v-velocity
estimation has an error less than ±0.1% at a location directly beneath
the center of the primary vortex (x = 25 mm) for different grid sizes

82

Comparison of the y extent from the wall over which the u-velocity
estimation has an error less than ±0.1% at a location directly beneath
the center of the primary vortex (x = 25 mm) and x = 30 mm . . .

83

Comparison of the y extent from the wall over which the u-velocity
estimation has an error less than ±0.1% at a location directly beneath
the center of the primary vortex (x = 25 mm) and x = 30 mm . . .

85

Table 2.1

Table 4.1

Table 4.2

Table 4.3

Table 4.4

Table A.1

Analytical flows employed in validation of MATLAB R code used to
derive Taylor series coefficients . . . . . . . . . . . . . . . . . . . . . 100

ix

LIST OF FIGURES

Figure 1.1

Figure 2.1

Figure 2.2

Figure 2.3

Figure 3.1

Schematic of wake vortices behind a finite-span airfoil (top) and an
aircraft (bottom). (Web source [13]) . . . . . . . . . . . . . . . . . .

2

Series expansions up to the 45th order of the streamwise velocity
compared to the analytical solution for Blasius boundary layer. (For
interpretation of the references to color in this and all other figures,
the reader is referred to the electronic version of this thesis.) . . . .

23

Series expansions up to the 45th order of the streamwise velocity compared to the analytical solution for Falkner-Skan (m = 1) boundary
layer . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

27

Series expansions up to the 45th order of the streamwise velocity
compared to the analytical solution for Stokes oscillating stream flow.
Different plots represent different phases of the oscillation cycle . . .

30

Experimental Set up (units are in cm) - Axisymmetric vortex impenging on a wall [3] . . . . . . . . . . . . . . . . . . . . . . . . . .
Coordinate axes, D0 = 0.0364 m and the core diameter 2Rc = 1.06
cm (used in this thesis as Rc = 0.0053 m[3]) . . . . . . . . . . . .

36

Figure 3.3

Flow visualization using LIF [3] . . . . . . . . . . . . . . . . . . . .

37

Figure 3.4

Near-wall velocity vectors, time in seconds [3] . . . . . . . . . . . . .

38

Figure 3.5

Geometry of the computational model to be used in ANSYS R -Fluent
based on experimental geometry given in Figure 3.2 [3] . . . . . . .

40

Experimental vorticity at t = 1.4 s along lines passing through the
vortex core center while parallel to r and z axis . . . . . . . . . . . .

42

Two-parameter fit, with rpeak = 0.018 m and zpeak = 0.02 m is taken
from from the experiments at t = 1.4 s, compared to experimental
data for the r-profile . . . . . . . . . . . . . . . . . . . . . . . . . . .

45

Figure 3.2

Figure 3.6

Figure 3.7

x

36

Figure 3.8

Figure 3.9

Figure 3.10

Figure 3.11

Figure 3.12

Figure 3.13
Figure 3.14

Figure 3.15
Figure 3.16

Figure 3.17
Figure 3.18

Figure 3.19
Figure 3.20

Two-parameter fit, with rpeak = 0.018 m and zpeak = 0.02 m is taken
from from the experiments at t = 1.4 s, compared to experimental
data for the z-profile . . . . . . . . . . . . . . . . . . . . . . . . . . .
Three-parameter fit compared to experimental data for the r-profile.
Initial position of the center of vortex rpeak and zpeak is obtained
from the fit . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Three-parameter fit compared to experimental data for the z-profile.
Initial position of the center of vortex rpeak and zpeak is obtained
from the fit . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

45

46

46

Experimental, r and z vorticity profiles compared with two-parameter
and three-parameter Gaussian fits . . . . . . . . . . . . . . . . . . .

47

Vorticity (ωθ ) contour of the experimental flow-field at t = 0.3 s (measured relative to the time of occurrence of the velocity field used to
set the initial Gaussian vortex parameters for the computation: t =
1.4 s relative to the solenoid opening) . . . . . . . . . . . . . . . . .
Vorticity (ωθ ) contour of the simulated flow-field at t = 0.3 s . . . .

50
50

Vorticity (ωθ ) contour of the experimental flow-field at t = 0.5 s (measured relative to the time of occurrence of the velocity field used to
set the initial Gaussian vortex parameters for the computation: t =
1.4 s relative to the solenoid opening) . . . . . . . . . . . . . . . . .
Vorticity (ωθ ) contour of the simulated flow-field at t = 0.5 s . . . .

51
51

Vorticity (ωθ ) contour of the experimental flow-field at t = 0.8 s (measured relative to the time of occurrence of the velocity field used to
set the initial Gaussian vortex parameters for the computation: t =
1.4 s relative to the solenoid opening) . . . . . . . . . . . . . . . . .
Vorticity (ωθ ) contour of the simulated flow-field at t = 0.8 s . . . .

52
52

Vorticity (ωθ ) contour of the experimental flow-field at t = 1 s (measured relative to the time of occurrence of the velocity field used to
set the initial Gaussian vortex parameters for the computation: t =
1.4 s relative to the solenoid opening) . . . . . . . . . . . . . . . . .
Vorticity (ωθ ) contour of the simulated flow-field at t = 1 s . . . . .

53
53

Comparison of the temporal evolution of the vorticity at the center of
the primary vortex. Results from simulations based on three-parameter
r,z and avg(r,z) Gaussian profiles are compared with those from MTV
experiments. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

55

xi

Figure 3.21

Figure 3.22

Figure 3.23

Figure 3.24

Axial- (z), and radial- (r) trajectories of the center of the primary
vortex. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

56

Effect of the computational time step on maximum vorticity of the
primary-vortex (a) and boundary layer (b) vorticity.(Γ0 = 46.11 ×
10−4 m2 /s, Rc = 0.005012 m, 300 iterations/time step, 600 × 600
grid) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

57

Effect of no iterations per time step on maximum vorticity of the
primary-vortex (a) and boundary layer (b) vorticity.(Γ0 = 46.11 ×
10−4 m2 /s, Rc = 0.005012 m, Time step = 0.005s, 600 × 600 grid)

59

Effect of different grid resolution on maximum vorticity of the primaryvortex (a) and boundary layer (b) vorticity. (Γ0 = 46.11 × 10−4
m2 /s, Rc = 0.005012 m, 300 iterations/time step, Time step = 0.005
s) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

61

Figure 4.1

Illustration of the numerical domain and initial/boundary conditions
used for simulation of a pair of counter rotating vortices having initial circulation Γ0 = 46.11 × 10−4 m2 /s and initial core radius
Rc = 0.005012 m introduced in a 0.06 × 0.06 m2 square-domain
at (0.018337,0.019602) m at t = 0. Because of problem symmetry,
only the upper half domain containing one of the vortices is computed. 65

Figure 4.2

Comparison of the Vorticity at the core center of the primary vortex
for an axisymmetric vortex ring (ωθ ) and 2D vortex (ωz ) . . . . . .

67

Vorticity (ωz in s−1 ) contour of Cartesian flow at t = 0.005 s
showing evolution of the primary vortex. . . . . . . . . . . . . . . .
Vorticity (ωz in s−1 ) contour of Cartesian flow at t = 0.25 s showing
evolution of the primary vortex. . . . . . . . . . . . . . . . . . . . .

69

Vorticity (ωz in s−1 ) contour of Cartesian flow at t = 0.50 s showing
evolution of the primary vortex. . . . . . . . . . . . . . . . . . . . .
Vorticity (ωz in s−1 ) contour of Cartesian flow showing at t = 0.605
s evolution of the primary vortex. . . . . . . . . . . . . . . . . . . .

70

Figure 4.3
Figure 4.4

Figure 4.5
Figure 4.6

Figure 4.7
Figure 4.8

Vorticity (ωz in s−1 ) contour of Cartesian flow at t = 0.75 s showing
evolution of the primary vortex. . . . . . . . . . . . . . . . . . . . .
Vorticity (ωz in s−1 ) contour of Cartesian flow at t = 1 s showing
evolution of the primary vortex. . . . . . . . . . . . . . . . . . . . .
xii

69

70

71
71

Figure 4.9

τwall profile at time t = 0.605 s . . . . . . . . . . . . . . . . . . . .

72

Figure 4.10

pwall profile at time t = 0.605 s . . . . . . . . . . . . . . . . . . . .

72

Figure 4.11

∂pwall
profile at time t = 0.605 s . . . . . . . . . . . . . . . . . . .
∂x

73

Figure 4.12

Figure 4.13

Figure 4.14
Figure 4.15

Figure 4.16

Figure 4.17

Streamwise (u-) and wall-normal (v-) profile in a plane at x = 30
mm (near outer edge of the primary vortex) at time t = 0.605 s as
determined by Fluent simulations for the case of Domain-600 . . . .

73

Streamwise (u-) and wall-normal (v-) profile in a plane at x = 25 mm
(nearest pixel to the maximum vorticity, at the center of the primary
vortex) at time t = 0.605 s as determined by Fluent simulations for
the case of Domain-600 . . . . . . . . . . . . . . . . . . . . . . . . .

74

Comparison of the wall-shear stress profile along the wall, at t = 0.605
s for Domain-500, -600, -715 . . . . . . . . . . . . . . . . . . . . . .
Comparison of the wall-pressure gradient along the wall, at t = 0.605
s for Domain-500, -600, -715 . . . . . . . . . . . . . . . . . . . . . .
u-Velocity estimation using Taylor-series expansions up to 15th order for the case of Domain-500 in a plane through the center of the
primary vortex (x = 25 mm) . . . . . . . . . . . . . . . . . . . . . .
v-Velocity estimation using Taylor-series expansions up to 15th order for the case of Domain-500 in a plane through the center of the
primary vortex (x = 25 mm) . . . . . . . . . . . . . . . . . . . . . .

Figure 4.18

Figure 4.19

u-Velocity estimation using Taylor-series expansions up to 15th order for the case of Domain-600 in a plane through the center of the
primary vortex (x = 25 mm) . . . . . . . . . . . . . . . . . . . . . .
v-Velocity estimation using Taylor-series expansions up to 15th order for the case of Domain-600 in a plane through the center of the
primary vortex (x = 25 mm) . . . . . . . . . . . . . . . . . . . . . .

Figure 4.20

Figure 4.21

u-Velocity estimation using Taylor-series expansions up to 15th order for the case of Domain-715 in a plane through the center of the
primary vortex (x = 25 mm) . . . . . . . . . . . . . . . . . . . . . .
v-Velocity estimation using Taylor-series expansions up to 15th order for the case of Domain-715 in a plane through the center of the
primary vortex (x = 25 mm) . . . . . . . . . . . . . . . . . . . . . .
xiii

75
75

78

78

79

79

80

80

Figure 4.22

Figure 4.23

Figure A.1
Figure A.2

u-Velocity estimation using Taylor-series expansions up to 15th order
for the case of Domain-600 in a plane through outer edge (away from
the axis) of the Primary vortex (x = 30 mm) . . . . . . . . . . . . .
v-Velocity estimation using Taylor-series expansions up to 15th order
for the case of Domain-600 in a plane through outer edge (away from
the axis) of the Primary vortex (x = 30 mm) . . . . . . . . . . . . .

84

84

Series expansions up to the 40th order of the radial velocity compared
to the analytical solution for axisymmetric stagnation flow . . . . . 103
Series expansions up to the 40th order of the wall-normal velocity
compared to the analytical solution for axisymmetric stagnation flow 103

xiv

Chapter 1
Introduction

1.1

Motivation

This research is motivated by the development of ground-based sensor arrays for detection
of wake vortices behind an aircraft that would enable take off and landings at faster rates
at commercial airports in the Unites States. This research, undertaken at the Flow Physics
and Control Laboratory (FPaCL) at Michigan State University, coupled with resources and
expertise from Michigan Aerospace Corporation (MAC), the sponsor of the current work, is
intended to lead to designing such a sensor system to be developed commercially (by MAC)
and installed at commercial airports in the country.
Current Federal Aviation Administration (FAA) regulations limit the rate at which aircrafts landings and take-offs happen due to mandatory spacing requirements between two
successive commercial aircrafts. It is required for a medium-sized aircraft to maintain 4 -6
NM (Nautical Mile) distance behind a heavy aircraft during landings and 2 - 3 NM during
take-offs [1]. This limit is imposed to avoid accidents due to wake vortices forming in the
1

wake of the leading aircraft interfering with the trailing aircraft during landing or take-off
(see Figure 1.1 [13]). Wake vortices are the vortical structures that are formed at the wingtip of an aircraft (or any lift producing system). These vortices could be quite strong and
they persist over a long time period depending upon the type of aircraft. Heavy aircrafts
generate stronger vortices that last longer compared to light-weight ones. These wake vortices interact with a trailing airplane and may cause it to roll and crash due to the associated
swirling flow or cause loss of altitude due to downwash (see Figure 1.1).

Figure 1.1: Schematic of wake vortices behind a finite-span airfoil (top) and an aircraft
(bottom). (Web source [13])

2

The ground-based sensor system will allow landings and take-offs at faster rate (15% [1])
through accurate detection of the vortex strength and the location of the wake vortices. It
is envisioned that multiple sensor-array units will be installed along the sides of runways
to capture near-ground-velocity (i.e. wall shear stress) and -pressure information near the
runway. The assessment of a method for utilization of this information to estimate the
above-ground flow field of wake vortices is the subject of the current study. Unlike other
surface-based flow estimation methods, such as used for feedback flow control, the uniqueness
of the current approach is that it does not require empirical information and/or approximate
models.
Current wake-vortex measurement systems at aiports include wind anemometers, pulse
and continuous-wave Light Detection and Ranging (LIDAR) and Radar Acoustic Sounding
Sensor (RASS). These technologies capture integrated data related to overall vortex circulation and not the detailed flow-field of the vortex structure.
The method examined here relies on obtaining a Taylor-series expansion of the velocity
field starting from the wall. The basic idea of this technique was published by Dallman et al.
[2] who demonstrated the ability to derive the coefficients of the Taylor-series expansion up
to an arbitrary order using only wall-shear-stress and -pressure information coupled with
recursive relations that are derived from the Navier-Stokes equations without any knowledge
of the far-field boundary condition. In the current research, we examine the possibility
that this concept can be extended to estimate the wake-vortex velocity field where the
far-field boundary is not known but knowledge of wall stresses can be obtained from the
near-wall-velocity and -pressure signatures captured by sensors installed near the runways.
If successful, the estimated velocity field would be analyzed to locate the cores of wake
3

vortices and find their strengths in real time.
The specific objectives of this research are to demonstrate the usability of the Taylorseries model for velocity estimation, understand the strengths and limitations of the method,
study parameters influencing its performance, and test its applicability for the wake vortex
problem.

1.2

Approach

In order to address the objectives outlined above, we followed a three-step approach-

1.2.1

Step I - Construction and Validation of the Taylor-series
model.

The expressions for the Taylor-series coefficients were initially hand-derived to understand
the underlying mathematical process and to identify recursive patterns that may be used
to construct a computer program to derive the coefficients. Using the identified patterns,
a MATLAB R program is developed to automate the process of deriving the expressions in
terms of the wall-shear and -pressure information for series coefficients up to any arbitrary
order. Additionally, a separate MATLAB R code is developed that is used in conjunction
with the program that generates the series-coefficients expressions to compute the values
of the coefficients based on the time history of the wall-shear stress and the wall-pressure
fields. This program, after computing the coefficients, is used to expand the series for the
streamwise u- and the wall-normal v-velocities up to a specified order.
The constructed programs required validation to ensure that the expressions for the
coefficients are computed accurately. We identified three fundamental flows with known
4

analytical solutions to compare against that estimated using the Taylor-series model for
validation purpose 1. Blasius Boundary Layer
2. Falkner Skan (m = 1; i.e. two-dimensional stagnation-point flow)
3. Stokes Oscillating Stream Problem

1.2.2

Step II - Numerical Computation to evaluate the applicability of the Taylor-series model to estimate the flow field of
a vortex-pair above a wall.

As discussed earlier, the motivation behind this study was to analyze the Taylor-series model
to see if it can be applied to estimate vortex flows. Given the infancy state of the estimation
model, it is sensible to test the method first under simple, controlled, yet relevant, flow
conditions. The selected flow problem is that of a counter-rotating vortex pair above a wall.
This problem can be viewed as analogous to the wake-vortex problem where the counterrotating vortices are generated at the tip of an aircraft wing which is symmetric about the
aircraft axis (parallel to the centerline of the runway). The vortices then convect down
towards the ground under the effect of their mutually-induced velocity.
The database of the vortex-pair-above-a-wall flow field was obtained from numerical simulation using ANSYS R -Fluent CFD solver. In order to validate the accuracy of the computational approach/parameters, it is necessary to compare against experiments. Though no
such data were accessible for the specific problem of interest, measurements in the closelyrelated problem of an axisymmetric ring above a wall by Gendrich et al. [3] could be found.
Thus, prior to computation of the two-dimensional vortex-pair simulation, a calculation of
5

the axisymmetric problem was undertaken. Upon validation, the same computational parameters (e.g. initial and boundary conditions, time step size, grid resolution, etc.) were
then used to generate the numerical data for the counter-rotating linear vortex pair above a
wall to study the Cartesian Taylor-series model.

1.2.3

Step III - Study of Cartesian Taylor Series model using simulations.

The numerical data for the Cartesian simulations yield the velocity-vector and the pressure
field. This can be used to compute the wall-stresses (and their time derivatives) by using the
central difference scheme. The computed stresses are then used to obtain the Taylor series
estimation of the velocity field and the result is compared against the actual flow field to
assess the accuracy and speed of convergence of the series. It is obvious that the accuracy
of the Taylor-series model will depend on how accurately the derivatives used to obtain the
series coefficients are computed. Finer grid should give better accuracy of the derivatives.
Hence, three spatial resolutions were used to study their effect on the estimation.

6

Chapter 2
Taylor-series Model For Estimation of
Incompressible Flow Fields Using
Wall-Stress Information
It is demonstrated in this chapter that the velocity field for an incompressible flow can be estimated using a Taylor-series-based mathematical model with coefficients that can be determined solely
from knowledge of the time history of the wall-shear stress and wall-pressure gradient.

2.1

Background

Wu et al. [4] developed a generalized theory for the interaction between a three-diamensional
moving body with compressible viscous flows. Their analysis showed that for incompressible
flows, it is possible to estimate the entire flow field with knowledge of initial conditions
as well as, the vorticity and its flux at the wall without having any knowledge of the farfield boundary conditions. The authors were able to reduce the Navier-Stokes equations
7

into a form that could be expressed in terms of generalized recursive relation between the
wall-vorticity and the wall-vorticity flux. These expressions can be used for computing wallnormal derivatives of the velocity and pressure field at the wall up to any arbitrary order.
Meaning, the entire flow field can be obtained from an expansion in terms of a continuous
(and infinite) series whose coefficients are essentially functions of the wall vorticity and
wall-vorticity flux, assuming that the series does converge. Dallman et al. [2] elaborated
on the same concept for two-dimentional incompressible flows by analyzing the fourth-order
stream-function (ψ) form of Navier-Stokes equations.

2.2

Concept in 2D

It is known that the two-dimensional (2D), velocity field (u, v in x, y directions respectively
for Cartesian coordinate system) for an incompressible flow can be expressed in terms of
spatial derivatives of a scalar field given by ψ(x, y, t), the stream function.
Typically, the following steps are involved in deriving the series coefficients in a Cartesian
system (similar details in cylindrical coordinate system are outlined in Appendix A for the
interested reader) :

2.2.1

Step I. Expanding the stream function (ψ) in a two-dimensional,
infinite Taylor-series form
∞ ∞
ψ=

aij xj y i+1

(2.1)

i=1 j=0
Where aij are the series coefficients (which are time dependent), and x and y are stream8

wise and wall-normal coordinates respectively. Note that the lowest power of y in the series
is 2. This is so because of the wall boundary conditions:

ψy=0 = 0

(2.2)

∂ψ
uy=0 =
=0
∂y y=0

(2.3)

This also guarantees that

vy=0 = −

2.2.2

∂ψ
=0
∂x y=0

(2.4)

Step II. Expressing the velocity components in series form
by differentiating the series for the stream function given by
Equation 2.1

Since the series defined for the stream function (ψ) is continuous and differentiable in space,
we differentiate it to obtain the series form for the velocity components. Considering the
streamwise velocity component (in x direction):

u=

∞ ∞
u=

∂ψ
∂y

aij (i + 1)xj y i

i=1 j=0
and the wall-normal velocity component (in y direction):
9

(2.5)

(2.6)

v=−

∞ ∞

∂ψ
∂x

aij (j)xj−1 y i+1

v=−

(2.7)

(2.8)

i=1 j=1

2.2.3

Step III. Obtaining series expressions for partial derivatives
of the velocity components in the momentum equation (the
components of which in the x and y directions are given by
Equations 2.20 and 2.24 respectively)

Starting with Equation 2.6 and differentiating with respect to x:

∂u
=
∂x

∂ 2u
=
∂x2

∞ ∞

aij (i + 1)(j)xj−1 y i

(2.9)

aij (i + 1)(j)(j − 1)xj−2 y i

(2.10)

i=1 j=1

∞ ∞
i=1 j=2

Similarly, starting with Equation 2.6 but differentiating with respect to y:

∂u
=
∂y

∂ 2u
=
∂y 2

∞ ∞

aij (i + 1)(i)xj y i−1

(2.11)

aij (i + 1)(i)(i − 1)xj y i−2

(2.12)

i=1 j=0

∞ ∞
i=2 j=0

10

Starting with Equation 2.8 and differentiating with respect to x:

∂v
=−
∂x

∞ ∞

aij (j)(j − 1)xj−2 y i+1

(2.13)

i=1 j=2

Differentiating a second time with respect to x:
∂ 2v
=−
∂x2

∞ ∞

aij (j)(j − 1)(j − 2)xj−3 y i+1

(2.14)

i=1 j=3

Similarly, Equation 2.8 is used to obtain the first and second derivatives of v with respect
to y:
∂v
=−
∂y

∂ 2v
=−
∂y 2

2.2.4

∞ ∞

aij (j)(i + 1)xj−1 y i

(2.15)

aij (i + 1)(i)(j)xj−1 y i−1

(2.16)

i=1 j=1

∞ ∞
i=1 j=1

Step IV. Obtaining expression relating the lowest-order Taylor series coefficients a1j to the wall shear stress using the
constitutive equation for Newtonian fluids

Starting with the constitutive equation for Newtonian fluids applied at the wall:

τwall = µ

∂u
∂y y=0

(2.17)

Substituting the series expression for ∂u (Equation 2.11), and evaluationg the expression
∂y
11

at the wall i.e. y = 0:
∞ ∞
τwall = µ

aij (i + 1)(i)xj y i−1

i=1 j=0

(2.18)
y=0

We observe that, only terms with y  will survive (i.e. i = 1) and the expression reduces
to the form below τwall
=
2µ

2.2.5

∞

a1j xj

(2.19)

j=0

Step V. Obtaining expression relating the next-order Taylor
series coefficients a2j to the wall pressure gradient using the
momentum equation in the streamwise direction

We use the momentum equation in the streamwise direction (x) and replace all derivatives
by their corresponding series forms derived in step 3 above:

ρ

∂u
∂u
∂u
+u
+v
∂t
∂x
∂y

=−

∂p
+µ
∂x

∂ 2u ∂ 2u
+
∂x2 ∂y 2

(2.20)

Rearranging terms in Equation 2.20 to obtain series expression for the pressure gradient
12

in the streamwise (x) direction:

−

1 ∂p
(2.21)
ρ ∂x
∞ ∞
=
aij (i + 1)xj y i
˙
i=1 j=0
∞ ∞ ∞ ∞
+
aij apq (i + 1)(p + 1)(q)xj+q−1 y i+p
i=1 j=0 p=1 q=1
∞ ∞ ∞ ∞
−
aij apq (j)(p)(p + 1)xj+q−1 y i+p
i=1 j=1 p=1 q=0
∞ ∞
∞ ∞
j−2 y i +
aij (i + 1)(i)(i − 1)xj y i−2 ]
aij (i + 1)(j)(j − 1)x
− v[
i=2 j=0
i=1 j=2

where, a indicates differentiation of a with respect to time. If we evaluate Equation 2.21 at
˙
the wall,

−

1 ∂p
(2.22)
ρ ∂x W all
∞ ∞
=
aij (i + 1)xj y i
˙
i=1 j=0
∞ ∞ ∞ ∞
+
aij apq (i + 1)(p + 1)(q)xj+q−1 y i+p
i=1 j=0 p=1 q=1
∞ ∞ ∞ ∞
−
aij apq (j)(p)(p + 1)xj+q−1 y i+p
i=1 j=1 p=1 q=0
∞ ∞
∞ ∞
j−2 y i +
− v[
aij (i + 1)(j)(j − 1)x
aij (i + 1)(i)(i − 1)xj y i−2 ]
y=0
i=1 j=2
i=2 j=0

13

In equation 2.22, only terms where the power of y is 0 will survive. Considering each
of the terms on the right hand side, the lowest power of y is 1 in the first term, 2 in the
second term, 2 in the third term, 1 in the fourth term and 0 in the last term. Thus, only
the last term has components that do not vanish at the wall, giving an expression for the
wall-pressure gradient in the streamwise direction in terms of a2j :
∞
1 ∂p
=
a2j xj
6µ ∂x wall
j=0

2.2.6

(2.23)

Step VI. Obtaining expressions for higher-order Taylor series coefficients in terms of a1j and a2j and their time derivatives by successive differentiation (in the wall-normal direction) of the 2D momentum equations and evaluating the
result at the wall

The series expression of the streamwise momentum equation is given by Equation 2.21 above.
The momentum equation in the wall-normal (y) direction is given by-

ρ

∂v
∂v
∂v
+u
+v
∂t
∂x
∂y

=−

∂p
+µ
∂y

∂ 2v ∂ 2v
+
∂x2 ∂y 2

(2.24)

Substituting for series expressions of the velocity components and their derivatives (obtained in step 3):
14

−

1 ∂p
(2.25)
ρ ∂y
∞ ∞
=
aij (i + 1)(j)xj−1 y i+1
˙
i=1 j=1
∞ ∞ ∞ ∞
−
aij apq (i + 1)(q)(q − 1)xj+q−2 y i+p+1
i=1 j=0 p=1 q=2
∞ ∞ ∞ ∞
+
aij apq (j)(q)(p + 1)xj+q−2 y i+p+1
i=1 j=1 p=1 q=1
∞ ∞
∞ ∞
j−3 y i+1 +
aij (i + 1)(i)(j)xj−1 y i−1 ]
aij (j)(j − 1)(j − 2)x
+ v[
i=1 j=1
i=1 j=3

We combine the expressions for both momentum equations (2.21 and 2.25) and rearrange
the terms in such a way that the highest starting value of index i, relating to the series order
with respect to y, in the summation will be on the left hand side and other lower-order y
terms appear on the right hand side. To do so, we first take the derivative with respect to
y of Equation 2.21, and exchange the order of the x and y derivatives of the pressure on the
left hand side. Subsequently, Equation 2.25 is used to express
in terms of series expansion. This leads to:

15

∂p
in the resulting equation
∂y

∞ ∞
v

aij (i + 1)(i)(i − 1)(i − 2)xj y i−3

(2.26)

i=3 j=0
∞ ∞
∞ ∞
=
aij (i)(i + 1)xj y i−1 +
˙
aij (j)(j − 1)xj−2 y i+1
˙
i=1 j=0
i=1 j=2
∞ ∞ ∞ ∞
+
aij apq (i + 1)(p + 1)(q)(i + p)xj+q−1 y i+p−1
i=1 j=0 p=1 q=1
∞ ∞ ∞ ∞
−
aij apq (j)(p + 1)(p)(i + p)xj+q−1 y i+p−1
i=1 j=1 p=1 q=0
∞ ∞ ∞ ∞
aij apq (i + 1)(q)(q − 1)(j + q − 2)xj+q−3 y i+p+1
+
i=1 j=0 p=1 q=0
∞ ∞ ∞ ∞
aij apq (p + 1)(q)(j)(j + q − 2)xj+q−3 y i+p+1
−
i=1 j=1 p=1 q=1
∞ ∞
− v[2
aij (i + 1)(i)(j)(j − 1)xj−2 y i−1 ]
i=1 j=2
∞ ∞
aij (j)(j − 1)(j − 2)(j − 3)xj−4 y i+1 ]
− v[
i=1 j=4

By applying Equation 2.26 at the wall, i.e. setting y = 0, we obtain equations for a3j
in terms of a1j , a2j and their time derivatives. To obtain expressions for a4j . Equation
2.26 is differented with respect to y, and the result is evaluated at the wall. The process
can continue indefinitely in this recursive manner to derive expressions for any aij (i > )
in terms of a1j , a2j and their time derivatives.
Some low-order coefficients were hand-derived initially to understand the process of obtaining these coefficients prior to developing a MATLAB R program to automate the deriva16

tion process. Given below are sample hand-derived expressions for some coefficients of the
series. In addition to clarifying the derivation process, the hand-derived expressions were
also helpful in validating the those obtained from the MATLAB R program.

a30 =

1
[2a − 8νa12 ]
˙
24ν 10

(2.27)

a31 =

1
[2a − 24νa13 ]
˙
24ν 11

(2.28)

a37 =

1
[2a − 288νa19 ]
˙
24ν 17

(2.29)

1
[6a + 4a10 a11 − 24νa22 ]
˙
120ν 20

(2.30)

1
[6a21 + 8a10 a12 + 4a11 2 − 72νa23 ]
˙
120ν

(2.31)

1
[6a + 12a10 a13 + 12a11 a12 − 144νa24 ]
˙
120ν 22

(2.32)

a40 =
a41 =
a42 =

1
a50 =
[24a30 + 4a12 + 24a10 a21 − 48νa14 − 96νa32 ]
˙
˙
720ν

2.3

(2.33)

Demonstration In Laminar Flows With Known Analytical Solution

In order to accurately estimate the velocity profile above the wall, it may be necessary to
employ high-order Taylor series expansion. As seen from Equation 2.26, the expressions
for the series coefficients become more complex at higher order. The underlying process
to derive these coefficients also becomes too time consuming and tedious. Therefore, it is
important to automate the process of computing the series coeficients up to any arbitrary
order. To this end, MATLAB R programs are designed and built based on the process
17

mentioned in section 2.2 to compute coefficients up to any arbitrary oder. The resulting
expressions are saved and used directly in expansion of series for different applications, to
avoid re-deriving the expressions for each flow situation, and hence reduce computational
time. However, prior to using the MATLAB-derived coefficient expressions, it is important to
validate these coefficients. This is done in two ways: (1) for low-order coefficients up to order
2, the MATLAB R expressions are compared to their analytically-derived counterpart; (2) by
comparing estimations of simple flows with their analytical solutions. Table 2.1 summarizes
the different types of flow analyzed for the purpose of validation.

Table 2.1: Analytical flows employed in validation of MATLAB R code used to derive Taylor
series coefficients
Case Flow
No.
1.
Blasius
Layer

2.

3.

Nature of flow

Boundary 1. 2D Steady flow.
2. Only wall-shear stress information is required
to obtain series coefficients.
3. Zero pressure gradient.
Falkner-Skan Bound- 1. 2D Steady flow.
ary Layer
2. Stagnation line flow.
(m = 1)
3. Wall-shear stress and wall-pressure gradient
information is required to obtain series coefficients.
4. Flow with favorable pressure gradient.
5. The streamwise wall-shear stress and the
wall-pressure gradient are linear.
Stokes
Oscillating 1. Unsteady flow.
Stream Problem
2. Wall-shear stress and wall pressure gradient
is sinusoidal in time.
3. Series coefficients are computed dynamically.

It is important to note that when validating the series predictions against analytical solutions, both the wall-shear stress and wall-pressure are available in analytical forms that can
be differentiated up to an arbitrary order without loss of accuracy. Therefore these valida18

tion cases enable verifications of the complex analytical forms obtained from the MATLAB R
program to relate the high-order series coefficients to the wall-shear and wall-pressure series
coefficients. They also help to provide a sense of the speed of convergence of the series. On
the other hand, the validation process does not provide an assessment of the effect of the
accuracy of the wall shear stress and pressure, and their derivatives up to the required order,
on the estimation accuracy.

2.3.1

Blasius Boundary Layer

The following derivation may be found in most standard textbooks on fluid mechancis (e.g.
Panton [9], pages 498:501).
Blasius boundary layer is a 2D steady boundary layer formed on a semi-infinite plate
oriented parallel to the flow direction. The solution is expressed in terms of a similarity
variable in the wall-normal direction:

η=y

U0
2xν

(2.34)

Where, y is the the wall-normal coordinate measured from the plate, x is the distance
along the plate, U0 is the free stream velocity and ν is the kinematic viscosity.
For a Blasius boundary layer problem, the wall-shear stress is given by:

τwall = µU0 f (0)

U0
2xν

u
where U = f (η). Equation 2.35 may be simplified to:

19

(2.35)

Ω
τwall = √
x

(2.36)

where,
Ω = µU0 f (0)

U0
2ν

(2.37)

The Taylor series coefficients in this case can be derived from knowledge of the wall-shear
stress only due to the absence of pressure-gradient. To do so, the coefficients of the infinite
series for τwall (on the left hand side of the Equation 2.19) are first computed as follows -

∞

a1j xj

(2.38)

j=0
1
x0
1 ∂
x1
1 ∂2
x2
(τwall )
+
(τwall )
+
(τwall )
2µ
0!
2µ ∂x
1!
2µ ∂x
2!
3
3
4
4
1 ∂
x
1 ∂
x
+
(τwall )
+
(τwall )
+ ...
2µ ∂x
3!
2µ ∂x
4!
=

Equation 2.38 can be rewritten by substituting for τwall from Equation 2.36 to get:

∞

a1j xj

(2.39)

j=0
1 (−1)Ω x1
1 Ω
1 (−1)(−3)Ω x2
√ )+
=
(
+
3 1!
5
2µ x
2µ
2µ
2!
2x 2
22 x 2
1 (−1)(−3)(−5)(−7)Ω x4
1 (−1)(−3)(−5)Ω x3
+
+
+ ...
7
9
2µ
3!
2µ
4!
23 x 2
24 x 2

20

Comparing the LHS and the RHS in Equation 2.39, expressions for and any a1j may be
obtained, such as:
1 Ω 1
0! 2µ 1
x2

(2.40)

1 Ω (−1)
1! 2µ 3
2x 2

(2.41)

1 Ω (−1)(−3)
5
2! 2µ
22 x 2

(2.42)

1 Ω (−1)(−3)(−5)
7
3! 2µ
23 x 2

(2.43)

1 Ω (−1)(−3)(−5)(−7)
9
4! 2µ
24 x 2

(2.44)

a10 =

a11 =

a12 =

a13 =

a14 =

Where x is the location at the wall where the series expansion is obtained. Given the above
coefficients, and similar higher order ones, the Taylor-series expansion for u(y) (Equation
2.6) may be written as:

u(y) = 2a10 y + 3a20 y 2 + 4a30 y 3 + 5a40 y 4 + 6a50 y 5 + ...

(2.45)

Or, using the expressions 2.27 through 2.33 and similar higher-order ones:

4
[2a10 − 8νa12 ]y 3
˙
24ν
5
+
[6a + 4a10 a11 − 24νa22 ]y 4
˙
120ν 20

u(y) = 2a10 y + 3a20 y 2 +

+

(2.46)

6
[24a30 + 4a12 + 24a10 a21 − 48νa14 − 96νa32 ]y 5 + ...
˙
˙
720ν

Eliminating terms containing the time derivatives (since the flow is steady) and a2j (since
21

∂p
= )from Equation 2.46:
∂x

a a
2νa14 5
4
y + ...
u(y) = 2a10 y − a12 y 3 + 10 11 y 4 −
3
6ν
5ν

(2.47)

Figure 2.1 shows different-order expansions of Taylor-series to estimate the streamwise
velocity profile in the Blasius boundary layer, plotted in similarity variables. Expansions
up to the 45th order where the order represents the order of y in the stream function, are
depicted in the figure and compared to the exact solution.
Two main observations can be made from Figure 2.1. First, the series expansion coincides
exactly with the solution over a certain η range, beyond which the series diverges. The η
range over which the agreement is accomplished becomes progressively bigger with increasing
series order. This provides validation of the MATLAB R code used to derive the series
coefficients. Specifically, the results show the appropriate behavior of a series expression of
the solution that increases in accuracy with increasing number of terms in the series. This
complements the validation done by comparing the code’s output to hand-derived expressions
(which could only be done for very low-order expressions).
A second observation is that the series convergence appears to be very slow. Though the
lowest order series expansion shown (5th order) accurately captures the velocity profile over
more than third of the 99% boundary layer thickness (η = 4.9), accurate representation over
the full boundary layer thickness is not possible even with as high of an order as 45th .

2.3.2

Falkner-Skan Boundary Layer (m = 1)

The following derivation may be found in most standard textbooks on fluid mechanics (e.g.
Panton [9], pages 508:512).
22

6

Order 5
Order 10
Order 15
Order 20
Order 25
Order 30
Order 35
Order 40
Order 45
Analytical

5

η

4
3
2
1
0
0

0.5
u/Uo

1

Figure 2.1: Series expansions up to the 45th order of the streamwise velocity compared to
the analytical solution for Blasius boundary layer. (For interpretation of the references to
color in this and all other figures, the reader is referred to the electronic version of this
thesis.)

The Flakner-skan boundary layer problem with favorable pressure gradient (m = 1) is a
steady state stagnation flow whose self-similar solution is expressed in terms of the following
similarity variable in the wall-normal direction:

η=

y√
Re
L

(2.48)

Where, y is the wall-normal coordinate, L is a characteristic length, Re is Reynolds
x
number based on L and u0 , and ue is the external velocity (ue (x) = u ( L )); x being the
streamwise coordinate.
The Falkner-Skan solution is given by-

x
u(x, y) = u0 ( )f (η)
L
23

(2.49)

Differentiating Equation 2.49 with respect to y:
√
u0 Re
∂u
x
=
f (η)
∂y
L
L

(2.50)

u Re
∂ 2u
x
= 0 f (η)
2y
2
L
∂
L

(2.51)

and the shear-stress on the wall:
√
u0 Re
x
τwall = µ
f (0)
L
L

(2.52)

We substitute Re = u L in Equation 2.52 to get a function for τwall :
ν
µ u 3
τwall = √ ( 0 ) 2 f (0)x
ν L

(2.53)

Now, the x-momemtum equation (Equation 2.20) if evaluated at the wall reduces to:
∂p
∂ 2u
=µ
∂x wall
∂ 2 y y=0

(2.54)

Combining equations 2.54 and 2.51:

u
∂p
= ρ( 0 )2 f (0)x
∂x wall
L

(2.55)

From Equations 2.53 and 2.55, it is evident that, the wall-shear stress and wall-pressure
gradient are linear in x. We now set the characteristic length L to 1 m and free stream
velocity u0 to 1 m/s for simplicity of calculations. If we compare 2.53 with the Equation
24

2.19, we get the coefficients a1j ’s as below-

a10 = 0

a11 =

f (0)
√
2 ν

a12 , a13 , a14 , ...a1∞ = 0

(2.56)

(2.57)
(2.58)

Similarly, if we compare 2.55 with the Equation 2.23, we get coefficients a2j ’s as below-

a20 = 0

a21 =

f (0)
6ν

a22 , a23 , a24 , ...a2∞ = 0

(2.59)

(2.60)
(2.61)

Given the above coefficients (Equations 2.56 through 2.61), and similar higher order ones
(anj ; n > ), the Taylor-series expansion for u(y) (Equation 2.6) can be written as:

u(y) = (2a10 y + 3a20 y 2 + 4a30 y 3 + ...)

(2.62)

+ (2a11 y + 3a21 y 2 + 4a31 y 3 + 5a41 y 4 + 6a51 y 5 + ...)x
+ (2a12 y + 3a22 y 2 + 4a32 y 3 + ...)x2 + ...

Combining Equations 2.62, 2.56 through 2.61 and 2.27 through 2.33:
a 2
u(y) = (2a11 y + 3a21 y 2 + 11 y 4 + ...)x
6ν
25

(2.63)

x
Or (recalling, ue (x) = u ( L ), u = L = ),
a 2
u
= 2a11 y + 3a21 y 2 + 11 y 4 + ...
ue
6ν

(2.64)

Figure 2.2 shows the expansions of Taylor-series to estimate the streamwise velocity
profile of (m = 1) boundary layer. Similar to the Blasius boundary layer case, comparison
with the exact solution gives confidence in the MATLAB R code used to calculate the series
coefficients. Unlike the Blasius boundary layer case, here the derived coefficients depend
also on the wall pressure gradient. Thus, the results give added validation of terms in the
coefficient expressions that are related to the wall pressure gradient. It is also noted here
that the series convergence seen in Figure 2.2 is also found to be very slow, with the 45th
series not able to accurately estimate the velocity over the entire boundary layer thickness.

2.3.3

Stokes Oscillatting Stream Problem

The following derivation may be found in most standard textbooks on fluid mechancis (e.g.
Panton [9], pages 221:224).
Stokes Oscillating stream problem is an unsteady flow (sinusoidal in time) that helps
us validate Taylor-series-coefficient expressions that contain high-order time derivatives of
the wall-shear stress and the wall-pressure gradient. This also helps us validate the Taylorseries-based flow-estimation due to dynamically changing coefficients.
The oscillating free-stream problem is analyzed based on the assumption that the fluid
has a single velocity component u(y, t). The x-momentum equation in this case reduces to:
∂P
∂ 2U
∂U
=
−
2Y
∂X
∂T
∂
26

(2.65)

4

Order 5
Order 10
Order 15
Order 20
Order 25
Order 30
Order 35
Order 40
Order 45
Analytical

3.5
3

η

2.5
2
1.5
1
0.5
0
0

0.2

0.4

0.6
u/Uo

0.8

1

Figure 2.2: Series expansions up to the 45th order of the streamwise velocity compared to
the analytical solution for Falkner-Skan (m = 1) boundary layer

p
x
u
Where - P = ρu ΩL , X = L , U = u , Ω is the frequency of oscillation (s−1 ), L is

a characteristic length, u0 is the free stream velocity, Y is similarity variable defined as
Y =y

Ω , and T is the time similarity variable defined as T = Ωt.
ν

The solution to the above problem subject to the boundary conditions associated with a
uniform stream oscillating above a fixed wall is given by:

Y
Y
U (y, t) = − sin(T − √ ) exp(− √ ) + sin T
2
2

(2.66)

Differentiating 2.66 with respect to T:

∂U
Y
Y
= − cos(T − √ ) exp(− √ ) + cos T
∂T
2
2
27

(2.67)

Differentiating 2.66 with respect to Y twice:
Y
exp(− √ )
Y
Y
√ 2 [sin(T − √ ) + cos(T − √ )]
2
2
2

(2.68)

∂ 2U
Y
Y
= − cos(T − √ ) exp(− √ )
2Y
2
2
∂

∂U
=
∂Y

(2.69)

Combining equations 2.65, 2.67, 2.69 and evaluating at the wall (Y = 0):

∂P
= − cos T
∂X wall

(2.70)

Evaluating 2.68 at the wall (Y = 0) and combining with 2.19 to obtain the series expression for τwall , we get:

∞
j=0

1
a1j xj = √ (sin T + cos T )
2 2

(2.71)

It is evident from 2.71 that the wall-shear stress is independent of X. Therefore, we get:

1
a10 = √ (sin T + cos T )
2 2

(2.72)

a11 , a12 , a13 , a14 , ...a1∞ = 0

(2.73)

Similarly, comparing Equation 2.70 and 2.23, we get:

a20 = −

cos T
6µ

a21 , a22 , a23 , a24 , ...a2∞ = 0

(2.74)

(2.75)

Given the above coefficients, and similar higher order ones, the Taylor-series expansion
28

for U (Equation 2.6) can be written as:

U = 2a10 Y + 3a20 Y 2 + 4a30 Y 3 + 5a30 Y 4 + 6a30 Y 5 + ...

(2.76)

Combining Equations 2.76, 2.72 through 2.75:

a
˙
a
˙
a
¨
U = 2a10 Y + 3a20 Y 2 + 10 Y 3 + 20 Y 4 + 10 Y 5 + ...
3ν
4ν
60ν 2

(2.77)

It is important to note that, for this problem all the higher-order spatial derivatives
are zero but time-derivatives (up to any order) of a10 and a20 exist, thereby making the
coefficients of the Cartesian Taylor-series completely independent of spatial derivatives.
Figures 2.3a, 2.3b, 2.3c, 2.3d, 2.3e show expansions to capture the boundary layer profile
using series up to the order of 45 for different time periods. It is observed that the near-wall
region, where the largest variation in the velocity takes place, was captured by the 10th
order series but even the 45th order series is not sufficient to capture the entire flat zone
of the profile. The convergence rate is relatively faster compared to Blasius boundary layer
and Falkner-Skan.

2.4

Summary

Some of the key findings of the present chapter are:
1. The constructed MATLAB R programs are capable of generating Taylor-series coefficients
up to any order accurately.
2. The agreement of series prediction with the theoretical solution gets progressively better
with increasing order of the expansion, but, with a very slow rate of convergence.
29

Order 5
Order 10
Order 15
Order 20
Order 25
Order 30
Order 35
Order 40
Order 45
Analytical

1

u/u0

0.5
0
−0.5
−1
−1.5
0

5

10
Y

15

1.5

Order 5
Order 10
Order 15
Order 20
Order 25
Order 30
Order 35
Order 40
Order 45
Analytical

1
0.5
u/u0

1.5

0
−0.5
−1
−1.5
0

20

5

(a) T = - π
2

1

u/u0

0.5
0
−0.5
−1
10
Y

15

Order 5
Order 10
Order 15
Order 20
Order 25
Order 30
Order 35
Order 40
Order 45
Analytical

1
0.5
0
−0.5
−1
−1.5
0

20

5

10
Y

15

20

(d) T = π
4

(c) T = 0
1.5

Order 5
Order 10
Order 15
Order 20
Order 25
Order 30
Order 35
Order 40
Order 45
Analytical

1
0.5
u/u0

20

1.5

u/u0

Order 5
Order 10
Order 15
Order 20
Order 25
Order 30
Order 35
Order 40
Order 45
Analytical

5

15

(b) T = - π
4

1.5

−1.5
0

10
Y

0
−0.5
−1
−1.5
0

5

10
Y

15

20

(e) T = π
2
Figure 2.3: Series expansions up to the 45th order of the streamwise velocity compared to
the analytical solution for Stokes oscillating stream flow. Different plots represent different
30
phases of the oscillation cycle

3. The domain of convergence of the series expansion (for up to the order of 45) is less
than the boundary layer thickness (except in the case of the unsteady Stokes flow, where the
convergence seems faster).

31

Chapter 3
Validation of The Computational
Approach
Having validated the concept to determine the flow field from knowledge of wall-stress information (without knowing the far-field boundary condition) using the Taylor-series based
expansion model, it is now required to examine the feasibility of this concept to estimate
a flow field that is relevant to the wake-vortex problem. As discussed in Chapter 1, the
selected problem is that of a counter rotating line-vortex pair impinging on a wall. For the
purpose of assessing the Taylor-series model it is necessary to have space-time information
of the flow field and associated wall stresses. The latter can then be used in conjunction
with the Taylor-series model to estimate the velocity field and compare the outcome with
the actual field.
Numerical calculation, employing Fluent, is used to produce the required database of the
vortex pair impinging on a wall. However, prior to using the database to validate the Taylor
series, it is necessary to ensure the accuracy of the computation. Ideally, this would be done
32

by comparing against experimental data of the same flow problem. Unfortunately, no such
data were available. Instead, it was possible to obtain experimental data in the closelyrelated problem of an axisymmetric vortex ring impinging on a wall. Thus, simulation of the
axisymmetric problem was undertaken first to validate the computational approach and the
choice of numerical parameters (grid resolution, domain size, time step, etc). Once confidence
was established in the latter, the same numerical approach and parameters were used for
computing the Cartesian counterpart problem of line vortex pair. In this chapter, details
of the computation of the axisymmetric configuration and validation against experimental
data are described. Information regarding the Cartesian-configuration computation and its
use to assess the Taylor-series model are left to the following chapter.

3.1

Background

The experimental data used for validation were generated in a study done at Turbulent
Mixing and Unsteady Aerodynamics Laboratory (TMUAL) at Michigan State University
(MSU) by Gendrich et al. [3]. This study was done for a vortex ring approaching a solid
wall (i.e. axisymmetric co-ordinate system).
This flow can be simulated using the proprietary code - ANSYS R -Fluent. Confidence
that such a computational model will provide results that match quantitatively with the
known experimental data is based on a study by Fabris et al. [5] who developed in-house
code to simulate a vortex ring impinging on a wall and did quantitative comparisons with
experiments.
The following steps are followed to generate the required numerical data set for final
estimations 33

1. Analyze the experiments performed by Gendrich et al. [3] to determine the initial conditions of the simulation that are consistent with the experiment.
2. Compare the simulation results with those from the experiments.
3. Optimize the simulation parameters, such as the time step, grid resolution, etc. by quantitative comparison with the experimental results.
4. Finalize parameters for the simulated model of the vortex ring problem.
5. Use the parameters fixed for the axisymmetric simulations and simulate the 2D planar
counterpart problem of the counter-rotating vortex pair.
6. The results from step 5 are then employed to generate the numerical data used to assess
the feasibility of using the Taylor series model to estimate vortex-dominated flows.

3.2
3.2.1

Description of The Experimental Data
Credits

The description in this entire section is based on the paper by Gendrich et al. [3]. Kindly
refer to this paper for more detailed illustration.

3.2.2

Introduction

The flow-field due to vortex ring interaction with a solid wall was studied using 2-color LaserInduced Fluorescence (LIF) and Molecular Tagging Velocimetry (MTV). A high Reynolds
number, based on circulation, (ReΓ = 4500) vortex ring was introduced in a flow-domain
and allowed to traverse under the action of self-induced velocity to interact with the wall.
34

This section gives a description of the experimental set up, MTV data, wall signature study
due to impingement of a vortex ring.

3.2.3

Experimental Facility

As seen in 3.1, a gravity-driven vortex ring generator was used to introduce a 0.06 m long
slug of water from a vertical tube for a short duration of time (0.5 seconds) controlled by a
solenoid valve. During the solenoid actuation, the head change was very small, less than 3
× 10−4 m. The wall was placed at a distance of 0.07 m normal to the vortex ring generator
axis. The initial vortex ring diameter was approximated to be D0 = 0.0364 m (Figure 3.1)
and the ring convection speed was U0 = 0.051 m/s. The resulting Reynolds number based
on the initial ring diameter was 1860 and 4500 based on the initial circulation (Γ0 = 45 ×
10−4 m2 /s).
The velocity field was quantified using Molecular Tagging Velocimetry (MTV). This
method is a non-intrusive flow diagnostic which measures two components of the velocity
simultaneously at many points in a 2D plane.

3.2.4

Observations

The vortex ring interaction with the wall exhibits different stages as it convects down from
the time it was introduced. Some of the prominent observations as seen in Figure 3.3 and
3.4 are listed below. These observations are consistent with those published by Fabris et al.
[5], Chu et al. [6], Orlandi et al. [7], Naguib et al. [8].
1. Vortex stretching effect is observed as the ring convects towards the wall. This is discerned
from the apparent reduction in the core size of the vortex ring.
35

reservoir tube
D = 3.81

Hayward SO120V19W
3/4” solenoid valve
generator tube
D = 2.54

93.5

7.0

(a)

Figure 3.1: Experimental Set up (units are in cm) - Axisymmetric vortex impenging on a
wall [3]

D0

2Rc

z
r
Figure 3.2: Coordinate axes, D0 = 0.0364 m and the core diameter 2Rc = 1.06 cm (used in
this thesis as Rc = 0.0053 m[3])

36

2. Visocus boundary layer growth is observed on the solid wall (light-colored dye).
3. The vortex ring causes boundary layer separation (t = 1.7 s in Figure 3.3).
4. The separated boundary layer forms a secondary vortex having opposite vorticity to
that of the primary vortex (marked with light-colored dye in Figure 3.3 at t = 1.8s). The
secondary vortex is seen to orbit around the primary one before traveling upwards under
self-induction effects (t = 2.0 to 3.8 s in Figure 3.3).
5. A tertiary vortex is also seen to form from boundary layer separation (see the outer edges
of the image at t = 2.67 s in Figure 3.3).

Figure 3.3: Flow visualization using LIF [3]

37

Figure 3.4: Near-wall velocity vectors, time in seconds [3]

38

3.3

Initial Conditions For The Simulation Model

In order to simulate the vortex ring interaction with the wall, an axisymmetric grid was
generated in Gambit. The numerical simulation was done using ANSYS R -Fluent solver.
The solver restricts the horizontal axis to be the symmetry axis for axisymmetric problems.
Hence, the geometry of the simulated problem is as shown in Figure 3.5 with the impingement
wall oriented vertically and the vortex ring’s axis is horizontal. Different approaches may
be used to initialize the computation. One method is to directly use experimental data at a
time instant that is early during the vortex ring evolution. In such an approach, one would
need to interpolate the measurements in order to map the experimental velocity data onto
the much finer computational grid. In the present work, such interpolation is avoided by
modeling the initial velocity field as that associated with a finite viscous core vortex having
Gaussian vorticity distribution. Gendrich et al.[3] showed that such a model represents the
actual vorticity distribution with good accuracy.
A user-defined function (UDF) was written in C and ’hooked’ to the fluent solver to compute the radial and wall-normal velocity components induced by the ’Gaussian-core’ vortex
and initialize the calculation. The specific characteristics of the vortex model (i.e. initial
circulation, core-radius, and core-center location) were determined from the experimental
data as described in the following section.

3.3.1

Gaussian Vorticity Distribution

The Gaussian vorticity distribution underlying the velocity field used to initialize the computation is given by ωθ = ω0 e

(r−rpeak )2 +(z−zpeak )2
−(
)
Rc 2
39

(3.1)

r
Streamwise direction

Circulation = 0.0045 m2/s
RC = 0.0053 m
Core (0.018, 0.020) m

Axis of symmetry
0.06 m

Wall

0.06 m
Wall

Wall

z
Wall-normal direction

Figure 3.5: Geometry of the computational model to be used in ANSYS R -Fluent based on
experimental geometry given in Figure 3.2 [3]

where,
ωθ : Out-of-plane (azimuthal) vorticity at a particular location in the flow domain.
ω0 : This is the initial peak vorticity (at the vortex core center) of the primary vortex ring
when introduced in the flow-domain (at t = 1.4 s, in experiments).
Rc : Initial core radius of the primary vortex.
rpeak , zpeak : The initial radial and wall-normal coordinate, respectively, of the vortex core
center. It is also the location of the peak vorticity.
r and z: The radial and wall-normal coordinates, respectively, of locations within the flow
domain.

The parameters in Equation 3.1 (ω0 , rpeak , zpeak and Rc ) are determined by fitting the
40

vorticity distribution given by Equation 3.1 to that extracted from the experimental velocity
vector field at time t = 1.4 s. As seen from the LIF visualization images in Figure 3.3,
the selected time corresponds to an instant where the vortex ring is approaching the wall
and there is very little activity within the boundary layer flow (light-colored dye) induced
by the ring. Because the Gaussian model provided good, but not perfect, representation of
the experimental vorticity distribution, the characteristics of the experimental distribution
exhibited some variation depending on the direction along which the distribution is extracted.
This can be seen more clearly in Figure 3.6 where the experimental vorticity distribution is
shown along two different lines that pass through the vortex core center (i.e. where the peak
vorticity is found) at t = 1.4 s: one line parallel to the r and the other to the z axis. As seen
from the figure, though the two vorticity distributions are qualitatively similar, they have
quantiative differences that will lead to different vortex parameters when fitting Equation 3.1
to each. Thus, for the purpose of the computation, three different sets of vortex parameters
are employed: one based on the vorticity distribution along r (referred to as r-profile), the
other based on that along z (or z-profile) and the last based on an average of the parameters
extracted from the r- and z- profiles.

In order to do the Gaussian fit, we are required to know the following information from
the experimental data: ω(r) and ω(z) through the vortex core center, and location of the
core center (rpeak , zpeak ). The former is known from the vorticity values at different
locations in the flow-field from the MTV data set [3]. Additionally, Gendrich et al. [3] give
the approximate location of the center of the votex ring at t = 1.4 s. For fitting the r-profile
(similar procedure is followed for the z-profile equation (3.1) can be rewritten as)41

ωθ(r)

60

ωθ(z)

ωθ (s−1)

50
40
30
20
10
0
0.005

0.01

0.015

0.02
0.025
r or z (m)

0.03

0.035

Figure 3.6: Experimental vorticity at t = 1.4 s along lines passing through the vortex core
center while parallel to r and z axis

ωθ = ω0 e

−(

r−rpeak 2
Rc )

(3.2)

With rpeak known from the experiments, this equation is fitted to determine the values
of ω0 and Rc . A second alternative fitting approach is also used in the present study. Since
the location of the center of the vortex ring is approximated to the nearest data grid point
in paper [3], better accuracy may be gained by choosing to make the initial position (rpeak ,
zpeak ) a parameter of the fit. These two different ways of doing the Gaussian fit will be
referred to as: two-parameter fit and three-parameter Fit
42

3.3.1.1

Two-Parameter Fit

The initial core center position, (rpeak , zpeak ) = (0.018,0.02) m, is known from the MTV
data [3][8]. A linear fit can be done after taking natural-log of both sides of Equation (3.2)-

1
ln(ωθ ) = ln(ω0 ) −
(r − rpeak )2
2
(Rc )

(3.3)

This expression is of the form y = mx + c, with x = (r − rpeak ) and y = ln(ωθ ). From
the linear fit (refer to Figures 3.7 and 3.8), we can determine theThe y-intercept, c = ln(ω ) which gives the initial vorticity (ω0 ),and the slope m = −

Rc
which gives the initial vortex core radius (Rc ). It is seen in the Figures 3.7 and 3.8 that the
experimental data points fall above and below the fit, with some scatter, in an asymmetric
manner. One of the two data branches represent data at r > Rc and the other at r < Rc .
For a perfect Gaussian behavior, the two branches should yield identical vorticity values (i.e.
they should collapse on one another). The discrepancy seen in Figures 3.7 and 3.8 is believed
to be due to inaccurate initial core center as well as deviation from true Gaussian behavior.

3.3.1.2

Three-Parameter Fit

Here, the vortex core center location (rpeak , zpeak ) is one of the fit parameters. Taking
natural-log of both sides of the Equation 3.2 -

ln(ωθ ) =

2rpeak
rpeak
−1 2
r +(
)r + (ln(ω0 ) − (
)2 )
2
2
Rc
Rc
Rc

(3.4)

This expression is of the form y = ax + bx + c, and a quadratic polynomial fit (refer to
Figures 3.9 and 3.10) will yield the constants43

a = − which gives the initial vortex core radius (Rc )

Rc
rpeak
b=
 that gives the initial location of the vortex center, and
Rc
rpeak
c = ln(ω ) − ( R ) that gives the initial peak vorticity (ω0 ).
c
Figure 3.11 shows comparison of the Gaussian fits compared to the data extracted along
a line parallel to r and z axis using two-parameter and three-parameter fits. Table 3.1
summarizes the values of three parameters rpeak or zpeak , Rc , and ω0 . The values of initial
maximum vorticity (ω0 ) and initial vortex core radius (Rc ) are used to compute the initial
circulation (Γ = πRc  ω : a relation that can be derived given the Guassian shape of the
vorticity distribution).
Table 3.1: Initial Vortex parameters obtained from the Gaussian fit
Profile Type
Maximum
of fit (n- Vorticity ω0
Parameter) (1/s)

Initial
Circulation
Γ = πRc  ω
×10−4
(m2 /s)

r
z
r
z

3.3.2

Position of Initial core
maximum radius Rc
vorticity
×10−3 (m)
rpeak or
zpeak
×10−3 (m)
18.000
20.000
18.337
19.602

35.7345
48.5570
33.5664
46.1137

2
2
3
3

60.2509
56.7836
62.2480
58.4330

4.4345
5.217
4.143
5.012

Comparison of Numerical Results With Experiments

The initial condition of the computation is specified using the initial vortex parameters listed
in Table 3.1. Although the parameters were obtained using two- and three- parameter fits,
the latter is more accurate as it extracts the initial location of the primary vortex center with
an accuracy better than the spacing of the measurement grid. Thus, for the purpose of the
simulations, the initial vortex location is set to - (rpeak , zpeak ) = (0.018337,0.019602) m.
44

4.5

Exparimental Data
Two−Parameter Fit

3.5

θ

ln(ω (r)) (s−1)

4

3
2.5
0

0.5

1
1.5
2 2
(r−rpeak) (m )

2

2.5

−5

x 10

Figure 3.7: Two-parameter fit, with rpeak = 0.018 m and zpeak = 0.02 m is taken from
from the experiments at t = 1.4 s, compared to experimental data for the r-profile

4.5

Exparimental Data
Two Parameter Fit

ln(ωθ(z)) (s−1)

4
3.5
3
2.5
0

1

2

2

2

(z−zpeak) (m )

3

4

−5

x 10

Figure 3.8: Two-parameter fit, with rpeak = 0.018 m and zpeak = 0.02 m is taken from
from the experiments at t = 1.4 s, compared to experimental data for the z-profile

45

4.2

ln(ωθ(r)) (s−1)

4
3.8
3.6
3.4
3.2
3
2.8

Exparimental Data
Three−parameter Fit
0.014 0.016 0.018 0.02 0.022 0.024 0.026
r (m)

Figure 3.9: Three-parameter fit compared to experimental data for the r-profile. Initial
position of the center of vortex rpeak and zpeak is obtained from the fit

4.2

ln(ωθ(z)) (s−1)

4
3.8
3.6
3.4
3.2
3
2.8

Exparimental Data
Three Parameter Fit
0.014 0.016 0.018 0.02 0.022 0.024 0.026
z (m)

Figure 3.10: Three-parameter fit compared to experimental data for the z-profile. Initial
position of the center of vortex rpeak and zpeak is obtained from the fit

46

Experimental r−Profile
Three Parameter r−Fit
Two Parameter r−Fit
Experimental z−Profile
Three Parameter z−Fit
Two Parameter z−Fit

ωθ(z) or ωθ(r) (1/s)

60
50
40
30
20
10
0
0

0.01

0.02

0.03 0.04
r or z (m)

0.05

0.06

Figure 3.11: Experimental, r and z vorticity profiles compared with two-parameter and
three-parameter Gaussian fits

To initialize the simulation, the velocity field induced by the Gaussian vorticity distribution
must be specified. The radial and wall-normal components of the velocity field are computed
from:
u=−

Γ(r)
y
2π x2 + y 2

(3.5)

v=+

Γ(r)
x
2 + y2
2π x

(3.6)

where, Γ(r) is the circulation given byr
Γ(r) =

ω(r)dr
0

contained in a circle of radius, r =

x + y 
47

(3.7)

The simulation employed 0.06 × 0.06 m square domain having 360,000 grid points (600
× 600 uniform grid resolution, corresponding to a resolution of 0.1 mm in either direction:
about one-fifth of the initial vortex core radius) and a time step of 0.005 s which is more
than two orders of magnitude bigger than the initial vortex core radius in each direction.
The fluid used was water with kinematic viscosity, ν = 1.004 × 10−6 m2 /s. The choices of
the spatial and the temporal resolution were confirmed to be appropriate as will be discussed
in the following sections. The boundary condition was set to the no-slip condition on three
of the four boundaries, and to symmetry axis on the fourth boundary. This agrees with the
experiment except for the boundary opposite to the impingement wall, which is a free surface
in the experiment. However in both the experiment and computation, the boundaries of the
domain are made sufficiently far away such that the boundary condition at the boundaries
opposite to the impingement wall and axis of symmetry should not influence the evolution
of the flow during the time of interest.

Figures 3.12 through 3.19 show a comparison of the experimental and computational
vorticity contour plots for different time instants during the flow evolution. It is interesting to
find that the vorticity contours of the simulated model match closely with their experimental
counterpart at the same time instant. This matching is particularly true for the evolution of
the primary vortex. The secondary vortex development also agrees closely, but not exactly.
Specifically, the development of the secondary vortex captured in the simulation seems to
lag by a small amount relative to the one seen in the experiments. This is likely because
the initial condition in the simulation matches the vorticity distribution associated with
the primary vortex while ignoring the boundary layer on the impingement wall (i.e. no
boundary layer exists at t = 0 in the simulation). As a result, it is probably not surprising
48

for the boundary layer development in the computation to lag behind the experiments. This
would cause a lag in the development of the secondary vortex as well, which forms from
the boundary layer. It is also worth noting that the boundary-layer vorticity level in the
simulation is larger than that seen in the experiment; e.g. compare Figures 3.14 and 3.15.
This is believed to be caused by the much finer resolution of the computation, yielding
vorticity values much closer to the wall than in the experiment.

The above discussion demonstrates the consistency of the flow physics extracted from the
computation in comparison to the experiments. To give a more quantitative measure of the
accuracy of the computation, the evolution of the peak vorticity at the center of the primary
vortex is obtained from the three computational cases and compared to the experimental
counter part in Figure 3.20. Note that since the inital peak vorticity for the computation is
obtained from the three-parameter curve fitting rather than taken directly from the experimental data, the plot in Figure 3.20 gives the vorticity at any time instant relative to the
initial vorticity. Inspection of Figure 3.20 shows that out of the three computational cases,
case 2 (with the initial condition based on the three-parameter fit to the z-profile of vorticity)
agrees the best with the experimental results. This may not be too surprising since ReΓ
0
(see Table 3.2)value computed from the vortex parameters is extracted from the z-profile
curve fit which is the closest to that extracted from the experiments. Additional quantitative measures that gives confidence in the computation may be obtained by comparing the
axial- and the radial- trajectory of the center of the primary vortex core from the Fluent
computation and the experiments (Figure 3.21). The comparison is shown in Figure 3.21.
Note that the initial coordinates of the vortex core (i.e. at zero time) do not agree exactly
with their counterpart from the experiment. This is because the latter are located based
49

250

0.03

200

0.025
z (m)

0.035

150

0.02

100

0.015

50

0.01

0

0.005
0

−50
0.01

0.02 0.03
r (m)

0.04

0.05

Figure 3.12: Vorticity (ωθ ) contour of the experimental flow-field at t = 0.3 s (measured
relative to the time of occurrence of the velocity field used to set the initial Gaussian vortex
parameters for the computation: t = 1.4 s relative to the solenoid opening)

0.04
300

0.035

250

0.03

200

z (m)

0.025

150

0.02

100

0.015

50

0.01

0

0.005

−50
0.01

0.02 0.03
r (m)

0.04

0.05

Figure 3.13: Vorticity (ωθ ) contour of the simulated flow-field at t = 0.3 s

50

0.035

150

0.03

100

0.025

50

0.015

z (m)

0.02

0

0.01

−50

0.005
0

0.01

0.02 0.03
r (m)

0.04

0.05

Figure 3.14: Vorticity (ωθ ) contour of the experimental flow-field at t = 0.5 s (measured
relative to the time of occurrence of the velocity field used to set the initial Gaussian vortex
parameters for the computation: t = 1.4 s relative to the solenoid opening)

0.04

250

0.035

200

0.03

150

z (m)

0.025

100

0.02
0.015

50

0.01

0

0.005

−50
0.01

0.02 0.03
r (m)

0.04

0.05

Figure 3.15: Vorticity (ωθ ) contour of the simulated flow-field at t = 0.5 s

51

0.035

100

0.03
50

z (m)

0.025
0.02

0

0.015
0.01

−50

0.005
0

0.01

0.02 0.03
r (m)

0.04

0.05

Figure 3.16: Vorticity (ωθ ) contour of the experimental flow-field at t = 0.8 s (measured
relative to the time of occurrence of the velocity field used to set the initial Gaussian vortex
parameters for the computation: t = 1.4 s relative to the solenoid opening)

0.04
100

0.035
0.03

50

z (m)

0.025
0.02

0

0.015
0.01

−50

0.005
0.01

0.02 0.03
r (m)

0.04

0.05

Figure 3.17: Vorticity (ωθ ) contour of the simulated flow-field at t = 0.8 s

52

150

0.035
0.03

100

z (m)

0.025
50

0.02
0.015

0

0.01
−50

0.005
0

0.01

0.02 0.03
r (m)

0.04

0.05

Figure 3.18: Vorticity (ωθ ) contour of the experimental flow-field at t = 1 s (measured
relative to the time of occurrence of the velocity field used to set the initial Gaussian vortex
parameters for the computation: t = 1.4 s relative to the solenoid opening)

0.04
0.035

150

0.03
z (m)

0.025

100

0.02

50

0.015

0

0.01
0.005

−50
0.01

0.02 0.03
r (m)

0.04

0.05

Figure 3.19: Vorticity (ωθ ) contour of the simulated flow-field at t = 1 s

53

on the location of peak vorticity at t = 0: a location that coincides with an experimental
data grid point. In comparison, as described earlier, the comutational results are based on
curve-fitting of a Gaussian vortex model to the experimental data, which yields a center
location which may fall in-between the experimental data grid points.

Table 3.2: Computational cases and associated initial vortex parameters. Listed experimental values are those given by Gendrich et al. [3]
Case

Profile

1
2
3
4

r
z
avg(r,z)
Experimental

3.3.2.1

Initial core Circulation Γ0
radius Rc ×10−4 (m2 /s)
×10−3
(m)
4.143
33.566
5.012
46.114
4.578
39.840
5.300
45.000

Reynolds No based
on initial circulation ReΓ = Γ
ν
3347
4593
3968
4500

Effect of computational time step

To physically justify the choice of the time step of 0.005 s, we consider the two different
time scales associated with the vortex-wall interaction: a diffusive time scale (τDif f usive ),
associated with the diffusion of viscous effects from the wall and a convective time scale
(τConvective ), associated with the advection of the vortex ring. If we compare the diffusive
time and the convective time, we find that the convective time scale is much smaller, and
therefore, it sets the more stringent constraint on the computational time step.This can be
explained as follows τConvective =

2
2
4Rc
Do
=
Γ0
Γ0

(3.8)

τDif f usive =

2
2
Do
4Rc
=
ν
ν

(3.9)

54

40

ωθ − ωθ0 (s−1)

30
20
10
R−Profile
Z−Profile
Average
Experimental

0
−10
0

0.2

0.4

t (s)

0.6

0.8

1

Figure 3.20: Comparison of the temporal evolution of the vorticity at the center of the
primary vortex. Results from simulations based on three-parameter r,z and avg(r,z) Gaussian
profiles are compared with those from MTV experiments.

From equations 3.8 and 3.9 -

τConvective
ν
1
=
=
τDif f usive
Γ0
ReΓ
0

(3.10)

ReΓ >> 1
0

(3.11)

τConvective << τDif f usive

(3.12)

Since-

then-

The expression in Equation 3.12 implies that the computational time step must be sufficiently small in comparison to the convective time scale. Our choice of 0.005 s is about four
55

0.035
0.03

z or r (m)

0.025
Axial−trajectory (Simulations)
Axial−trajectory (Experiments)
Radial−trajectory (Simulations)
Radial−trajectory (Experiments)

0.02
0.015
0.01
0.005
0

0.2

0.4

0.6
t (s)

0.8

1

Figure 3.21: Axial- (z), and radial- (r) trajectories of the center of the primary vortex.

times smaller than the convective time scale (0.0217 s). To ascertain that this is sufficiently
small, three other computations similar to that reported above except for the time step size
were undertaken. The evolution of the vorticity at the primary vortex center as well as of
the boundary layer for all four cases is depicted in Figure 3.22a and 3.22b respectively. Note
that the boundary layer vorticity is determined by the slope of the streamwise velocity at
the wall (ωBoundaryLayer = ∂u
).
∂z z=
From Figures 3.22a and 3.22b, we see that the vorticity profiles converge for all time
steps smaller or equal to 0.005 s. It is also observed that time steps greater than 0.005
s underestimate the magnitude of the vorticity of the primary vortex and boundary layer
vorticity. Based on this observation, the computation time step was set to 0.005 s.
56

100

ωθ (s−1)

90
80
70

0.01 s
0.005 s
0.0025 s
0.00125 s

60
0

0.2

0.4

t (s)

0.6

0.8

1

(a) Vorticity of Primary Vortex

−100

ωθ (s−1)

−150
−200
−250
−300

0.01 s
0.005 s
0.0025 s
0.00125 s

−350
−400
0

0.2

0.4

t (s)

0.6

0.8

1

(b) Boundary Layer Vorticity
Figure 3.22: Effect of the computational time step on maximum vorticity of the primaryvortex (a) and boundary layer (b) vorticity.(Γ0 = 46.11 × 10−4 m2 /s, Rc = 0.005012 m,
300 iterations/time step, 600 × 600 grid)
57

3.3.2.2

Effect of number of iterations per time step

Number of iterations per time step is one of the important parameters affecting the accuracy
of the computed velocities in Fluent. The utilization of an integral global-error criterion to
determine convergence of the solution at each time step did not yield accurate results. We
believe this stemmed from the fact that the computational domain is much larger than the
vortex size, and hence the velocity at most nodes was close to zero. To remedy this problem,
we systematically varied the number of iterations to see its effect on the accuracy of the
computed vorticity (primary vortex and boundary layer vorticity) for the time step of 0.005
s. The number of iterations was set to 200, 300, 400 and 500 per time step to verify the
convergence. The results are shown in Figures 3.23a and 3.23b.
It is evident from Figures 3.23a and 3.23b that 200 iterations underestimate vorticities
but larger number of iterations, 300 through 500, produce converging profiles. Thus, 300
iterations per time step is selected in the final computation. Larger number of iterations
does not give more accuracy and would increase the computational time considerably.

3.3.2.3

Effect of grid size

Spatial resolution is important for accurate computation. The MSU research license provided
by Fluent allow a maximum of 715 × 715 grid points. Hence, we decided to analyze the
accuracy of the computations for four grid sizes, viz. 1. 300 × 300
2. 500 × 500
3. 600 × 600
4. 715 × 715
58

100

ωθ (s−1)

90
80
70

200
300
400
500

60
0

0.2

0.4

t (s)

0.6

0.8

1

(a) Vorticity of Primary Vortex

−100
−150

ωθ (s−1)

−200
200
300
400
500

−250
−300
−350
−400
0

0.2

0.4

t (s)

0.6

0.8

1

(b) Boundary Layer Vorticity
Figure 3.23: Effect of no iterations per time step on maximum vorticity of the primary-vortex
(a) and boundary layer (b) vorticity.(Γ0 = 46.11 × 10−4 m2 /s, Rc = 0.005012 m, Time
step = 0.005s, 600 × 600 grid)
59

The resulting evolution of primary-vortex and boundary-layer vorticity for all four cases
can be observed in Figures 3.24a and 3.24b respectively. The figure demonstrates that a grid
size of 500 × 500 or bigger would resolve the primary vortex evolution properly. On the other
hand, the boundary-layer vorticity does not completely converge with an approximately 5%
difference seen between computation employing the largest two grid sizes. Of course, it is
possible that convergence is achieved with the 715 × 715 grid, but this can’t be ascertained
without using an even larger grid. Since it is not possible to increase the number of grid
points (given the Fluent license limitation), the reader is advised that computations based
on the finest possible grid are uncertain to within 5%. Additionally, to gauge the influence
of this error on the Taylor-series-based estimation, we will carry out flow estimations using
the Taylor-series model for data obtained from the three different grid sizes.

3.4
3.4.0.4

Summary of Fluent Parameters
General Attributes:

Solver: Fluent (12.1 version), 2D- Planar, Laminar without gravity effect
Solver Type: Pressure-Based
Velocity Formulation: Absolute
Time: Transient
Model used: Viscous - Laminar
Material: Water (Liquid) [Density (ρ): 998.2 kg/m3 , Dynamic Viscosity (µ): 0.001003
kg/m − s]

60

100

ωθ (s−1)

90
80
70

300 x 300
500 x 500
600 x 600
715 x 715

60
0

0.2

0.4

t (s)

0.6

0.8

1

(a) Vorticity of Primary Vortex

0

ωθ (s−1)

−100
−200
300 x 300
500 x 500
600 x 600
715 x 715

−300
−400
0

0.2

0.4

0.6
t (s)

0.8

1

(b) Boundary Layer Vorticity
Figure 3.24: Effect of different grid resolution on maximum vorticity of the primary-vortex
(a) and boundary layer (b) vorticity. (Γ0 = 46.11 × 10−4 m2 /s, Rc = 0.005012 m, 300
iterations/time step, Time step = 0.005 s)

61

3.4.0.5

Solution Methods:

Pressure-Velocity Coupling:
Scheme: PISO (Pressure Implicit solution by Split Operator method)
Skewness Correction: 0
Neighbor Correction: 1
Skewness-Neighbor Coupling: No
Spatial Discretization:
Gradient: Least Squares Cell Based
Pressure: Standard
Momentum: Third-Order MUSCL
Transient Formulation:
Second-Order Implicit with no non-iterative Time Advancement and no Forzen Flux formulation.

3.4.0.6

Solution Controls:

Under-Relaxation Factors:
Pressure: 0.3
Density: 1
Body Forces: 1
Momentum: 0.7

62

Chapter 4
Simulation of A Pair of Line Vortices
Impinging On A Solid Wall And
Estimation of Velocities Using The
Cartesian Taylor Series Model
In Chapter 3, we discussed how an axisymmetric vortex ring impinging on a wall can be
simulated accurately in ANSYS R -Fluent. The parameters used to simulate the vortex
ring problem are now used to simulate a pair of linear vortices impinging on a wall; i.e.
the Cartesian geometry counterpart. The numerical database will contain time-resolved
velocity-field information, which will be used to calculate the wall stresses (wall-shear stress
and wall-pressure gradient) that, in turn, enable us to apply the Cartesian Taylor-series
model to estimate the velocity components. The estimated field is then compared to the
true field to assess the viability of the model to capture the characteristics of the vortex pair
63

(primarily the center location and the strength of the vortex).

4.1

Geometry of the Computational Model

A pair of counter rotating vortices is introduced in a 0.06 × 0.06 m2 flow-domain (Figure
4.1). As determined in the previous chapter, based on the three-parameter Gaussian fit,
a vortex is introduced at (xpeak , ypeak ) = (0.018337,0.019602) m. The vortex has initial
circulation (Γ0 ) of 46.11 × 10−4 m2 /s with the core radius Rc = 0.005012 m. Only the
vortex in the first-quadrant is shown in the geometry given in Figure 4.1 since symmetry is
imposed relative to the y-axis in Fluent and hence there is no need to explicitly solve for the
line vortex on the other side of the y-axis. Each of the vortices in the pair induce a velocity
on the other vortex causing the vortex to approach the vertical wall placed along the x-axis.
Subsequently, the image (virtual) vortex in the wall (required to maintain zero wall-normal
velocity at the wall) induce a velocity on the real vortex, causing it to travel in the positive x
direction, hereafter, referred to as the streamwise direction (while the y axis will be referred
to as the wall-normal direction). The reference point is at the top-left is where the pressure
value is used as a reference in calculating the surface-pressure coefficient.
The computational grid is prepared in Gambit. Three different grids are used in order
to study the effect of spatial resolution 1. 500 × 500 - Called Domain-500 - The square domain is divided into 500 equal intervals
along both directions.
2. 600 × 600 - Called Domain-600 - The square domain is divided into 600 equal intervals
along both directions.
3. 715 × 715 - Called Domain-715 - The square domain is divided into 715 equal intervals
64

x
Streamwise direction
Wall

Circulation = 0.004611 m2/s
RC = 0.00501 m
Core (18.337e-3, 19.602e-3) m

Axis of symmetry
0.06 m

Wall

0.06 m
Wall

Reference Point

y
Wall-normal direction

Figure 4.1: Illustration of the numerical domain and initial/boundary conditions used for
simulation of a pair of counter rotating vortices having initial circulation Γ0 = 46.11 ×
10−4 m2 /s and initial core radius Rc = 0.005012 m introduced in a 0.06 × 0.06 m2 squaredomain at (0.018337,0.019602) m at t = 0. Because of problem symmetry, only the upper
half domain containing one of the vortices is computed.

along both directions. As noted in Chapter 3, this is the maximum allowable number of grid
points under the current MSU ANSYS R -Fluent license.

4.2

Simulations Using ANSYS R -Fluent

ANSYS R -Fluent solver is used to carry out the simulation described in the previous section.
The grid geometry from Gambit is imported in Fluent and the solution is initialized (with
the Gaussian vortex of initial circulation, initial core radius and center of the vortex core
as given in section 4.1) using a User Defined Function (UDF) in Fluent. The vortex profile
65

is estimated to be Gaussian and the solver is run to obtain the solution over a flow time
period of 1.11 seconds. This duration is long enough for the primary vortex to approach the
wall and induce the formation of a boundary layer and subsequent formation of a secondary
vortex.
At each time steps (of 0.005 s) 300 iterations were carried out and the data files are
saved and the computed Fluent data files are converted into text files using TecPlot-10
built-in capability of importing Fluent data files. These text files are then processed using
MATLAB R to obtain the wall stresses. Once read in MATLAB R the wall-shear stress
and the wall-pressure gradient and their time derivatives are computed using the central
difference scheme to obtain the higher order coefficients in the Taylor-series model. The
central difference formula used to compute the wall-shear stress, wall-pressure gradient, and
their time derivatives is given below-

f (x) =

f (x + h) − f (x − h)
2h

(4.1)

where h is a distance between two adjacent nodes.
Before beginning to analyze the Cartesian simulations we compared time evolution of
the vorticity at the core of the primary vortex against that obtained in the axisymmetric
simulation. Figure 4.2 shows the comparison of the two cases. It is interesting to note
how the early stages of development of the primary vortex are affected by the coordinate
geometry. The primary vortex of the same initial conditions evolves in different ways in
cylindrical and Cartesian system. The cylindrical - axisymmetric vortex ring - problem show
initial increase in the vorticity of the primary vortex as it approaches the wall. This is caused
by stretching of the vortex ring diameter, which results in proportional decrease in the core
66

size to conserve the ring’s volume. Subsequently, conservation of angular momentum (over a
time scale during which viscosity has no significant influence) results in amplification of the
vorticity, as seen in Figure 4.2. This causes increase in the vorticity, hence we see the peak
on a vorticity profile in Figure 4.2. This observation is consistent with results published by
Fabris et al. [5], and Chu et al. [6]. In case of a Cartesian - counter rotating line pair of
vortex - system, no vortex stretching is possible and the vorticity continues decrease.

100

ωz or ωθ (s−1)

90
80
70

Axisymmetric Vortex Ring
2D vortex

60
50
40
0

0.2

0.4

0.6
t (s)

0.8

1

Figure 4.2: Comparison of the Vorticity at the core center of the primary vortex for an
axisymmetric vortex ring (ωθ ) and 2D vortex (ωz )

Figures 4.3, 4.4, 4.5, 4.6, 4.7, 4.8 depict how the primary vortex evolves. Figure 4.3,
shows the vorticity contour at the initial time t = 0.005 s one time step after the primary
vortex is introduced. As the computations continue beyond this time, we see the vortex
moving away from the axis of symmetry and approaching the wall causing formation of the
boundary layer. The boundary layer thickens and at one point separates (see Figure 4.7)
and rolls up into a secondary vortex with the opposite-sign vorticity to that of the primary
67

vortex.

4.3

Simulation Results

The time instant t = 0.605 s (file number 121, see Figure 4.6) was selected for checking
the viability of the Taylor series estimation model. This time instant was deemed suitable
since both the primary vortex and boundary layer are prominent (a situation akin to that
of the wake vortex above ground, though at a much lower Reynolds number) but without
the added complexity of a well formed secondary vortex.
Figures 4.9, 4.10 and 4.11 depict the wall-shear stress, the wall-pressure and the wallpressure gradient profiles at t = 0.605 s, respectively. These profiles are important as their
nodal information is used to compute spatial derivatives. Figures 4.12 and 4.13 show the
u- and v- velocity profiles in a plane through x = 0.025 m (location of the primary vortex
core) and x = 0.03 m (outer boundary of the primary vortex). These are the velocity profiles
we will be estimating using our Cartesian Taylor-series model. Readers should note that
it is difficult to identify the exact vortex core location in the computational grid and it is
approximated to the grid point nearest to the actual center of the core by finding the location
of maximum vorticity within the vortex. The almost zero v-velocity profile in Figure 4.13 is
an indication that we are near the core.
Also, in order to expand the Taylor-series model up to a high order, high-order time
derivatives of the wall-shear stress and the wall-pressure gradient are required. The stencil
size of the central-finite-difference scheme used to compute these derivatives increases with
the increase in derivative order. Thus, more time snap shots of the flow must be available
before and after the instance at which the derivatives are computed. The selected time for
68

0.06

200

0.05

150

y (m)

0.04
100

0.03

50

0.02

0

0.01
0
0

0.02

x (m)

0.04

0.06

−50

Figure 4.3: Vorticity (ωz in s−1 ) contour of Cartesian flow at t = 0.005 s showing evolution
of the primary vortex.

0.06

100

0.05

y (m)

0.04

50

0.03
0

0.02
0.01
0
0

0.02

x (m)

0.04

0.06

−50

Figure 4.4: Vorticity (ωz in s−1 ) contour of Cartesian flow at t = 0.25 s showing evolution
of the primary vortex.

69

0.06

100

0.05

y (m)

0.04

50

0.03
0.02

0

0.01
0
0

0.02

x (m)

0.04

0.06

−50

Figure 4.5: Vorticity (ωz in s−1 ) contour of Cartesian flow at t = 0.50 s showing evolution
of the primary vortex.

0.06

100

0.05

y (m)

0.04

50

0.03
0.02

0

0.01
0
0

0.02

x (m)

0.04

0.06

−50

Figure 4.6: Vorticity (ωz in s−1 ) contour of Cartesian flow showing at t = 0.605 s evolution
of the primary vortex.

70

0.06

100

0.05

y (m)

0.04

50

0.03
0

0.02
0.01
0
0

0.02

x (m)

0.04

0.06

−50

Figure 4.7: Vorticity (ωz in s−1 ) contour of Cartesian flow at t = 0.75 s showing evolution
of the primary vortex.

80

0.05

60

0.04
y (m)

0.06

40

0.03

20
0

0.02

−20

0.01
0
0

−40
0.02

x (m)

0.04

0.06

Figure 4.8: Vorticity (ωz in s−1 ) contour of Cartesian flow at t = 1 s showing evolution of
the primary vortex.

71

0.15

τwall (Pa)

0.1

0.05

0

−0.05
0

0.02

x (m)

0.04

0.06

Figure 4.9: τwall profile at time t = 0.605 s

1

pwall (Pa)

0.5
0
−0.5
−1
−1.5
0

0.02

x (m)

0.04

Figure 4.10: pwall profile at time t = 0.605 s

72

0.06

150
100

(dp/dx)

wall

(Pa/m)

50
0
−50
−100
−150
−200
0

0.02

Figure 4.11:

x (m)

0.04

0.06

∂pwall
profile at time t = 0.605 s
∂x

0.06

u−Velocity
v−Velocity

0.05

y (m)

0.04
0.03
0.02
0.01
0
−0.1

−0.05

0
0.05
Velocity (m/s)

0.1

Figure 4.12: Streamwise (u-) and wall-normal (v-) profile in a plane at x = 30 mm (near
outer edge of the primary vortex) at time t = 0.605 s as determined by Fluent simulations
for the case of Domain-600

73

0.06

u−Velocity
v−Velocity

0.05

y (m)

0.04
0.03
0.02
0.01
0
−0.1

−0.05

0
0.05
Velocity (m/s)

0.1

0.15

Figure 4.13: Streamwise (u-) and wall-normal (v-) profile in a plane at x = 25 mm (nearest
pixel to the maximum vorticity, at the center of the primary vortex) at time t = 0.605 s as
determined by Fluent simulations for the case of Domain-600

the estimation is roughly half-way through the computed flow duration. This makes sufficient
data available before/after the time of estimation, enabling calculation of time derivatives
up to 15th order (and Taylor series order up to 45).
Figures 4.14 and 4.15 compare the effect of grid size on the wall-shear stress and the wallpressure gradient. All three curves, corresponding to different grid sizes (Domain-500, -600,
-715) converge, giving confidence that the grid size (spatial resolution) has no detectable
effect on the accuracy of the computed wall-shear stress and the wall-pressure gradient.
However, as we employ high order derivatives, minute differences may have significant influence on the accuracy of the derivatives at high order. Hence, we decided to estimate the
velocity components for all three grid sizes and compare the convergence of Taylor-series to
the actual velocity profiles (from Fluent simulations) for all three cases.
74

0.15

Domain−500
Domain−600
Domain−715

τwall (Pa)

0.1

0.05

0

−0.05
0

0.02

x (m)

0.04

0.06

Figure 4.14: Comparison of the wall-shear stress profile along the wall, at t = 0.605 s for
Domain-500, -600, -715

150

(dp/dx)

wall

(Pa/m)

100
50
0
−50
−100

Domain−500
Domain−600
Domain−715

−150
−200
0

0.02

x (m)

0.04

0.06

Figure 4.15: Comparison of the wall-pressure gradient along the wall, at t = 0.605 s for
Domain-500, -600, -715

75

4.4

Estimations Using Cartesian Taylor-series Model

The data set generated by Fluent simulations for three different grid sizes (Domain-500,
Domain-600, and Domain-715) is further processed to estimate the u- and v-velocity (streamwise and wall-normal, respectively) components. The data set provides the velocity components and the pressure within the flow. Wall-pressure is extracted from the Fluent data
file and the central difference scheme is applied to compute the wall-pressure gradient along
the streamwise direction. The time derivatives of the wall-pressure gradient at each node
on the wall is computed by applying the central difference method on the data available at
the same node from adjacent data files, from the previous (t-0.005 s) and the next time step
(t+0.005 s) (see Equation 4.1).
For a given estimation, the velocity profiles are estimated at the same x location as
the wall point at which the wall-shear and wall-pressure-gradient derivatives are obtained.
Thus, the series expansion is only done in y, but not in x (i.e. using a one-dimensional
series). That is, estimates at a different x location are accomplished by moving the point of
wall-stresses observation to the new x location rather than fixing the latter and using a twodimensional series expansion. We believe this approach provides better accuracy by avoiding
the expansion in the second dimension. In the following, there are two sets of comparisons
done for the final velocity estimations based on I. The same x location (x = 25 mm; directly beneath the vortex core center) but different
grid resolutions (three cases, Domain-500, Domain-600, and Domain-715)
II. The same grid resolution (600 × 600) but two different estimation locations: x = 25 mm
(similar to part 1 above) and x = 30 mm (below outer periphery of the primary vortex where
the wall-stress signatures are weaker and velocity variation is less).
76

4.5
4.5.1

Results and Discussion
Comparison of estimations based on different grid resolutions

u- and v-Velocity estimations are done for the three cases of different grid sizes (Domain-500,
Domain-600, and Domain-715). The idea is to study the effect of spatial resolution on
the final estimations. As observed before (Figures 4.14 and 4.15), grid size doesn’t have
significant effect on the wall stresses at a given node on the wall, but high order expansions of
the Taylor-series may still exhibit sensitivity to the grid resolution due to error in computing
high order derivatives using the central difference scheme with different spatial step size
(0.12, 0.10 and 0.084 mm for Domain-500, Domain-600, and Domain-715 respectively).
The results are shown in Figure 4.16, 4.18 and 4.20 for u estimations and 4.17, 4.19
and 4.21 for the v estimations. In each figure, the true velocity profile, obtained from
Fluent, is shown using red diamonds. Other symbols/line colors represent the estimation
done using different series orders up to 15. Note that the series order is the order of y in
the stream function. Also, the domain extent in the wall-normal direction is truncated to
only encompass the y range for which accurate estimation could be obtained. This range
is approximately less than one tenth the height of the vortex core center above the wall for
all Taylor-series orders considered. Outside the y range shown, the series estimate deviates
considerably from the true velocity value such that it is not meaningful to show results for
a wider domain.
Inspection of Figures 4.16 through 4.21 reveals a pattern in which increasing the Taylorseries estimation order, starting from the first order, systematically increases the range of
77

1.5

−3

x 10

Fluent
y1
y2
y3
y4
y5
y6
y7
y8
y9
y10
y11

y (m)

1

0.5

0
0

0.02

0.04
u (m/s)

y12
y13
y14
y15

0.06

Figure 4.16: u-Velocity estimation using Taylor-series expansions up to 15th order for the
case of Domain-500 in a plane through the center of the primary vortex (x = 25 mm)
2

Fluent

−3

x 10

y1
y2
y3

y4
y5
y6
y7
y8
y9
y10
y11

y (m)

1.5

1

0.5

y12
y13
0
0

0.005
v (m/s)

0.01

y14
y15

Figure 4.17: v-Velocity estimation using Taylor-series expansions up to 15th order for the
case of Domain-500 in a plane through the center of the primary vortex (x = 25 mm)
78

1.5

−3

Fluent
y1

x 10

y2
y3
y4

1

y (m)

y5
y6
y7
y8
y9
y10
y11

0.5

0
0

0.02

0.04

u (m/s)

y12
y13
y14
y15

0.06

Figure 4.18: u-Velocity estimation using Taylor-series expansions up to 15th order for the
case of Domain-600 in a plane through the center of the primary vortex (x = 25 mm)

2

Fluent
y1

−3

x 10

y2
y3
y4
y5
y6
y7

y (m)

1.5

1

y8
y9
y10
y11

0.5

0
0

0.005
v (m/s)

0.01

y12
y13
y14
y15

Figure 4.19: v-Velocity estimation using Taylor-series expansions up to 15th order for the
case of Domain-600 in a plane through the center of the primary vortex (x = 25 mm)
79

1.5

Fluent
y1

−3

x 10

y2
y3
y4
y5
y6
y7
y8

y (m)

1

y9
y10
y11

0.5

0
0

0.02

0.04
u (m/s)

y12
y13
y14
y15

0.06

Figure 4.20: u-Velocity estimation using Taylor-series expansions up to 15th order for the
case of Domain-715 in a plane through the center of the primary vortex (x = 25 mm)

2

Fluent
y1

−3

x 10

y2
y3
y4
y5
y6
y7
y8
y9
y10
y11

y (m)

1.5

1

0.5

0
0

0.005
v (m/s)

0.01

y12
y13
y14
y15

Figure 4.21: v-Velocity estimation using Taylor-series expansions up to 15th order for the
case of Domain-715 in a plane through the center of the primary vortex (x = 25 mm)
80

y over which good agreement is obtained between the series and the true velocity. This
trend, however, continues only up to a certain order beyond which the range of agreement
becomes smaller rather than bigger. This is unlike what one would expect from both a
theoretical point of view (i.e. increasing the order should increase the accuracy of the series)
as well as from results presented in Chapter 2 concerning the series expansion for simple
flows with known analytical solution. The key difference, however, between the latter and
the present case is that for the simple flows, the wall stresses are known analytically, and
hence they can be differentiated up to an arbitrary order without loss of accuracy. In the
present case, the stresses are known from a numerical solution and one would expect that,
even with very high spatial and temporal resolution and minimum ’numerical noise’, the
accuracy of the derivatives will start to deteriorate beyond a certain order. Thus, the order
of the series beyond which the results in Figures 4.16 through 4.21 become less accurate may
be considered as a break-even point; at which further increase in the estimation accuracy
by increasing the series order is more than offset by the deterioration in the accuracy of the
computed derivative.

In attempt to quantify the assessment of the Taylor series effectiveness in capturing the
true velocity field as well as to compare the influence of the grid resolution systematically,
the series order that gives agreement between the true and estimated velocity (to within
0.1%) over the largest y range is identified. A list of this order for the three grid resolutions,
along with the y range of agreement is given in Tables 4.1 and 4.2 for estimations of u and
v respectively. One of the key observation for the table is that the y range over which the
series convergence is satisfactory is very small compared to the core center of the primary
vortex. This shows that the Taylor-series-based estimation would be unusable for capturing
81

the vortex velocity field and associated characteristics. Specifically, the convergence of the
series is too slow to be a useful tool for wake vortex detection.
Finally, Tables 4.1 and 4.2 show that although some differences are seen in the most
accurate series order and corresponding y range for different grid resolutions, these differences
are not substantial. Essentially for all cases the most accurate estimate is obtained for a
series order of 8 - 12 over a domain range of approximately 1 mm (or 5.7057% of the core
vortex center).

Table 4.1: Comparison of the y extent from the wall over which the u-velocity estimation
has an error less than ±0.1% at a location directly beneath the center of the primary vortex
(x = 25 mm) for different grid sizes
Case

Series
Order

Location % of the
above
core vortex
the wall center
(mm)

Domain-500
Domain-600
Domain-715

12
12
12

0.9000
0.9500
0.7972

5.4152
5.7057
4.7860

Actual
velocity
(Fluent)
×10−3
m/s
58.44
60.26
53.78

Estimated
velocity
×10−3 (m/s)

58.39
60.28
53.79

Table 4.2: Comparison of the y extent from the wall over which the v-velocity estimation
has an error less than ±0.1% at a location directly beneath the center of the primary vortex
(x = 25 mm) for different grid sizes
Case

Series
Order

Location % of the
above
core vortex
the wall center
(mm)

Domain-500
Domain-600
Domain-715

8
12
10

0.780
0.950
0.797

4.6931
5.7057
4.7848

82

Actual
velocity
(Fluent)
×10−3
m/s
2.839
3.537
2.910

Estimated
velocity
×10−3 (m/s)

2.837
3.540
2.913

4.5.2

The same grid resolution (600 × 600) but two different estimation locations: x = 25 mm and x = 30 mm

Similar comparison is done for the case of Domain-600 in which we expanded series in two
different planes. The two planes selected are 1. x = 25 mm : A plane passing through the center of the primary vortex.
2. x = 30 mm : A plane passing through the outer edge of the Primary vortex.
The Tables (4.3 and 4.4) are self-explanatory. It shows that u-velocity can be estimated
with low order series in a plane through x = 30 mm as compared to the one passing through
the center of the vortex but the domain of estimation above the wall is very small as compared
to that with 12th order series for the expansions in plane x = 25 mm. v-velocity, on the
other hand, has the lower order expansion capturing more distance above the wall in the
plane x = 30 mm as compared to the case of plane x = 25 mm.
Table 4.3: Comparison of the y extent from the wall over which the u-velocity estimation
has an error less than ±0.1% at a location directly beneath the center of the primary vortex
(x = 25 mm) and x = 30 mm
Case

Series
Order

Location % of the
above
core vortex
the wall center
(mm)

x = 25 mm
x = 30 mm

12
6

0.950
0.350

5.7057
2.1021

83

Actual
velocity
(Fluent)
×10−3
m/s
60.26
55.19

Estimated
velocity
×10−3 (m/s)

60.28
55.18

2.5

−3

Fluent
y1

x 10

y2
y3
y4
y5
y6
y7
y8
y9
y10
y11

2

y (m)

1.5
1
0.5
0
0

0.02

0.04
u (m/s)

0.06

y12
y13
y14
y15

Figure 4.22: u-Velocity estimation using Taylor-series expansions up to 15th order for the
case of Domain-600 in a plane through outer edge (away from the axis) of the Primary vortex
(x = 30 mm)
2

Fluent
y1
y2
y3
y4

−3

x 10

y (m)

1.5

1

0.5

0
0

0.005
0.01
v (m/s)

y5
y6
y7
y8
y9
y10
y11
y12
y13
y14
y15

Figure 4.23: v-Velocity estimation using Taylor-series expansions up to 15th order for the
case of Domain-600 in a plane through outer edge (away from the axis) of the Primary vortex
(x = 30 mm)

84

Table 4.4: Comparison of the y extent from the wall over which the u-velocity estimation
has an error less than ±0.1% at a location directly beneath the center of the primary vortex
(x = 25 mm) and x = 30 mm
Case

Series
Order

Location % of the
above
core vortex
the wall center
(mm)

x = 25 mm
x = 30 mm

12
8

0.950
1.050

5.7057
6.3036

85

Actual
velocity
(Fluent)
×10−3
m/s
3.537
4.748

Estimated
velocity
×10−3 (m/s)

3.540
4.747

Chapter 5
Summary, Conslusions and
Recommendations
As outlined in Chapter I, we followed a three-step approach during this research. This
chapter captures key results and conclusions from each step.

5.1

Summary and Conclusions

We successfully constructed a Taylor-series model, based on the ideas discussed by Dallman
et al. [2], that enables the estimation of the velocity field from wall-stress information in
Cartesian geometry. MATLAB R programs were constructed to derive and compute Taylorseries coefficients up to any arbitrary order. The model and the MATLAB R program
were validated using known analytical flows, namely Blasius boundary layer, 2D steady
stagnation-point Falkner Skan (m = 1) flow, and Stokes oscillating stream problem. These
validation tests show that the estimations do agree very closely with the analytical solution
but one of the important findings was that the series converges very slowly. Steady flows such
86

as Blasius boundary layer and Falkner Skan (m = 1) show very slow convergence pattern.
Specifically, though the solution over more than one-third of the boundary layer thickness is
captured accurately with a 5th order series, a 45th order series was insufficient to represent
the solution over the entire boundary layer extent. The solution of Stokes oscillating stream
problem showed faster rate of convergence. Though the reason for this is unclear, it is notable
that in this case the wall shear stress and pressure gradient are invariant in space, and hence
the series coefficients are independent of spatial derivatives of the wall stresses.
To address the motivation of potentially using the Taylor-series model in sensing aircraft
wake vortices, the model was used to estimate the velocity field of a counter-rotating vortex pair above a wall. The flow field was generated using numerical simulation employing
ANSYS R -Fluent. In order to perform the simulation accurately we needed some baseline
flow to validate the simulated flow. We identified the experiments - of the closely-related
problem of an axisymmetric vortex ring impinging on a solid wall [2]- conducted at Michigan
State University - as a good bench mark to validate the computational approach against.
The computation was initialized using Gaussian vorticity distribution, the parameters of
which were set by matching the distribution to the experimental data, in the least-squares
sense. Grid and time-step sensitivity tests were undertaken, showing that the employed
2
resolutions of 0.918 (normalized based on the factor Rc ) and 0.02 (normalized based on
Γ0
the initial vortex core radius (Rc )) are sufficient to accurately resolve the evolution of the
primary vortex ring.
Having identified the computational parameters for accurate simulation of a vortex above
a wall, the flow of a two-dimensional counter-rotating vortex pair was computed. We performed three simulations for three different grid sizes (500 × 500, 600 × 600, 715 × 715).
87

We observed that there is no significant effect of grid size on computation of the wall-shear
stress and the wall-pressure gradient but it did influence the accuracy of estimated velocities
with different orders of expansions. Most significantly, similar to the simple analytical flow,
the series convergence was extremely slow, enabling the capture of the velocity field of only a
portion of the boundary layer (approximately 5% of the vortex core center height). However,
unlike the simple analytical flows, the accuracy of the series did not increase monotonically
with increasing series order. The trend only held up to 8th to 10th order, then further
increase in series order resulted in less accurate estimation. This is attributed to the loss of
accuracy in computing high-order derivatives from numerical data.
Overall, as much as succeeded to validate the Taylor-series model for estimation of velocities in the flow-field, the accuracy of the estimation and its practical usage is severely
limited by the spatial and temporal resolution required to accurately compute high order
derivatives. The model application is also limited by several factors such as knowledge of
flow stresses at large number of surface points, data captured over a large time range, and the
slow convergence rate. Therefore, regarding the original research motivation, we conclude
that the Taylor-series model is incapable of estimating the vortex flows behind an aircraft.

5.2

Recommendations

There is further scope for analysis of the Taylor-series model as well as to refine the computation of the vortex-ring/counter-rotating-vortex-pair impinging on a wall. Some of the
associated recommendations are listed below1. As we determined that the high-order expansions do not provide good estimation over
the complete boundary layer, we can cover the near-wall region with the Taylor-series model
88

and use this information as a boundary condition to estimate the remaining domain using
efficient function bases. This approach will reduce the number of coefficients considerably
and requires fewer boundary conditions. These function bases can be obtained, for instance,
from Proper Orthogonal Decomposition of wake-vortex data.
2. The model needs to be tested for different temporal resolutions to examine the effect of
the latter on the accuracy and quality of the estimation.
3. Some study could be done to see if the Taylor-series model can be used for Boundary
layer control (since the latter application may only require very-near-wall information).
4. It will be advantageous to perform simulations with non-uniform grid in the near-wall
region. This will improve the accuracy of the wall-shear stress and the wall-pressure, thereby
reducing the estimation error propagating while computing high order derivatives.
5. Current approach is based on introduction of a primary vortex in flow domain at t =
0. This is equivalent to the experimental time t = 1.4 s, when the vortex was initially
introduced. MTV data set provides us with the velocity vector array captured at t = 1.4
s. It would be interesting to simulate the flow by supplying velocities at each point in the
flow domain as initial condition in Fluent, instead of introducing the primary vortex. This
approach will simulate the initial boundary layer (if any) in the experiments and give us
better near wall information, which will ultimately affect the accuracy of the Taylor-series
model.

89

APPENDICES

90

APPENDIX A
Taylor-series Model for Axisymmetric
Flow
We derived expressions for Taylor-series coefficients in Cartesian system in section 2.2. Similar exercise was done to derive coefficients for the Cylindrical Taylor-series model.

A.1

Derivation of the Taylor-series Model for Axisymmetric Flow

A.1.1

Step I. Expanding the stream function (ψ) in a two-dimensional,
infinite Taylor-series form
∞ ∞

aij rj+1 z i+2

ψ=

(A.1)

i=0 j=0
Where aij are the series coefficients (which are time dependent), and r and z are streamwise and wall-normal coordinates respectively. Note that the lowest power of z in the series
91

is 2. This is so because of the wall boundary conditions:

ψz=0 = 0

Vry=0 =

−1 ∂ψ
=0
r ∂z z=0

(A.2)

(A.3)

This also guarantees-

1 ∂ψ
=0
vzy=0 =
r ∂r z=0

A.1.2

(A.4)

Step II. Expressing the velocity components in a series form
by differentiating the series for the stream function given
by Equation A.1

Since the series defined for the stream function (ψ) is continuous and differentiable in space,
we differentiate it to obtain the series form for the velocity components. Considering the
streamwise velocity component (in r direction):

vr =

−1 ∂ψ
r ∂z

∞ ∞

aij (i + 2)rj z i+1

vr = −
i=0 j=0

92

(A.5)

(A.6)

Wall-normal velocity component (in z direction):

1 ∂ψ
r ∂r

vz =
∞ ∞

aij (j + 1)rj−1 z i+2

vz = −

(A.7)

(A.8)

i=0 j=1

A.1.3

Step III. Obtaining series expressions for partial derivatives
of the velocity components in the momentum equations A.20

Starting with Equation A.6 and differentiating with respect to r:

∂vr
=−
∂r

∞ ∞

aij (i + 2)(j)rj−1 z i+1

(A.9)

aij (i + 2)(j)(j − 1)rj−2 z i+1

(A.10)

i=0 j=1

Differentiating A.9 with respect to r :
∂ 2 vr
=−
∂r2

∞ ∞
i=0 j=2

Differentiating equation A.5 with respect to z:

∂vr
=−
∂z

∞ ∞

aij (i + 2)(i + 1)rj z i

(A.11)

aij (i + 2)(i + 1)(i)rj z i−1

(A.12)

i=0 j=0

Differentiating A.11 with respect to z:
∂ 2 vr
=−
∂z 2

∞ ∞
i=1 j=0

93

Similarly, differentiating equation with respect to z:

∂vz
=
∂z

∞ ∞

aij (i + 2)(j + 1)rj−1 z i+1

(A.13)

i=0 j=1

Differentiating Equation A.13 with respect to z:

∂vz
=
∂z 2

∞ ∞

aij (i + 2)(i + 1)(j + 1)rj−1 z i

(A.14)

i=0 j=1

Differentiating Equation A.7 with respect to r:

∂vz
=
∂r

∞ ∞

aij (j − 1)(j + 1)rj−2 z i+2

(A.15)

i=0 j=2

Differentiating Equation A.15 with respect to r:
∂ 2 vz
=
∂r2

A.1.4

∞ ∞

aij (j + 1)(j − 1)(j − 2)rj−3 z i+2

(A.16)

i=0 j=3

Step IV. Obtaining expression relating the lowest-order Taylor series coefficients a0j to the wall shear stress using constitutive equation for Newtonian fluids

Starting with the constitutive equation for Newtonian fluids:

τwall = µ

∂vr
∂z z=0

(A.17)

Substituting series expression for ∂vr (Equation A.11), and evaluationg the expression
∂z
94

at wall i.e. z = 0:
∞ ∞
τwall = −µ

aij (i + 2)(i + 1)rj z i

i=0 j=0

(A.18)
z=0

We observe that, only terms with z  will servive (i.e. i = 0) and the expression reduces
to the form belowτwall
=
−2µ

A.1.5

∞

a0j rj

(A.19)

j=0

Step V. Obtaining expression relating the next-order Taylor
series coefficients a1j to the wall pressure gradient using the
momentum equation in the streamwise direction

We use momemtum equation in streamwise direction (r) and replace all differential terms by
their correponsding series forms derived in step 4 above.

ρ

∂vr
∂vr
∂vr
+ vr
+ vz
∂t
∂r
∂z

=−

∂p
∂
+µ
∂r
∂r

1 ∂
∂ 2 vr
(rvr ) +
r ∂r
∂z 2

(A.20)

Rearranging terms in Equation A.20 to obtain series expression for pressure gradient in
streamwise (r) direction. We get95

−

1 ∂p
ρ ∂r

(A.21)
∞ ∞

=−

aij (i + 2)rj z i+1
˙

i=0 j=0
∞ ∞ ∞ ∞
+
i=0 j=0 p=0 q=1
∞ ∞ ∞ ∞
−

aij apq (i + 2)(p + 2)(q)rj+q−1 z i+p+2
aij apq (j + 1)(p + 1)(p + 2)rj+q−1 z i+p+2

i=0 j=1 p=0 q=0
∞ ∞
∞ ∞
j−2 z i+1 +
aij (i + 1)(i + 2)(i)rj z i−1 ]
aij (i + 2)(j + 1)(j − 1)r
+ v[
i=1 j=0
i=0 j=2

where, a indicates differentiation with respect to time. If we evaluate the Equation A.21
˙
at the wall,

−

1 ∂p
(A.22)
ρ ∂r W all
∞ ∞
=−
aij (i + 2)rj z i+1
˙
i=0 j=0
∞ ∞ ∞ ∞
+
aij apq (i + 2)(p + 2)(q)rj+q−1 z i+p+2
i=0 j=0 p=0 q=1
∞ ∞ ∞ ∞
−
aij apq (j + 1)(p + 1)(p + 2)rj+q−1 z i+p+2
i=0 j=1 p=0 q=0
∞ ∞
∞ ∞
j−2 z i+1 +
+ v[
aij (i + 2)(j + 1)(j − 1)r
aij (i + 1)(i + 2)(i)rj z i−1 ]
y=0
i=0 j=2
i=1 j=0

96

In equation A.22, only terms where the power of z is 0 will survive. Considering each of
the terms on the right hand side, the lowest power of z is 1 in the first term, 2 in the second
term, 2 in the third term, 1 in the fourth term and 0 in the last term. Thus, only the last
term contributes and has components that do not vanish at the wall, giving an expression
for the wall-pressure gradient in the streamwise direction in terms of a2j :

∞
1 ∂p
=−
a1j rj
6µ ∂r wall
j=0

A.1.6

(A.23)

Step VI. Obtaining expressions for higher-order Taylor series coefficients in terms of a0j and a1j and their time derivatives by successive differentiation (in the wall-normal direction) of the 2D momentum equations and evaluating the
result at the wall

The series expression of the streamwise momentum equation is given by Equation A.21
above. The momentum equation in the wall-normal (z) direction is given by-

ρ

∂vz
∂vz
∂vz
+ vr
+ vz
∂t
∂r
∂z

=−

∂p
1 ∂
+µ
∂z
r ∂r

97

r

∂vz
∂r

+

∂ 2 vz
∂z 2

(A.24)

−

1 ∂p
ρ ∂z
∞ ∞
=
aij (j + 1)rj−1 z i+2
˙
i=0 j=1
∞ ∞ ∞ ∞
−
aij apq (i + 2)(q + 1)(q − 1)rj+q−2 z i+p+3
i=0 j=0 p=0 q=2
∞ ∞ ∞ ∞
+
aij apq (j + 1)(q + 1)(p + 2)rj+q−2 z i+p+3
i=0 j=1 p=0 q=1
∞ ∞
aij (i + 2)(j + 1)(j − 1)(j − 1)rj−3 z i+2 ]
− v[
i=0 j=3
∞ ∞
− v[
aij (j + 1)(i + 2)(i + 1)rj−1 z i ]
i=0 j=1

(A.25)

We combine expressions for both the momentum equations (A.21 and A.25 and rearrange
the terms in such a way that highest order expansion coefficients (indicated by index i) will

98

be on the left hand side and other lower order terms appear on the right hand side. We get∞ ∞
v

aij (i + 2)(i + 1)(i)(i − 1)rj z i−2

(A.26)

i=2 j=0
∞ ∞
∞ ∞
=
aij (i + 2)(i + 1)rj z i +
˙
aij (j + 1)(j − 1)rj−2 z i+2
˙
i=0 j=0
i=0 j=2
∞ ∞ ∞ ∞
−
aij apq (i + 2)(p + 2)(q)(i + p + 2)rj+q−1 z i+p+1
˙ ˙
i=0 j=0 p=0 q=1
∞ ∞ ∞ ∞
+
aij apq (j + 1)(p + 2)(p + 1)(i + p + 2)rj+q−1 z i+p+1
˙ ˙
i=0 j=1 p=0 q=0
∞ ∞ ∞ ∞
aij apq (i + 2)(q + 1)(q − 1)(j + q − 2)rj+q−3 z i+p+3
˙ ˙
−
i=0 j=0 p=0 q=2
∞ ∞ ∞ ∞
aij apq (p + 2)(q + 1)(j + 1)(j + q − 2)rj+q−3 z i+p+3
˙ ˙
−
i=0 j=1 p=0 q=1
∞ ∞
− v[2
aij (i + 2)(i + 1)(j + 1)(j − 1)rj−2 z i ]
˙
i=0 j=2
∞ ∞
aij (j + 1)(j − 1)2 (j − 3)rj−4 z i+2 ]
˙
− v[
i=0 j=4

A.2

Demonstration In Axisymmetric Laminar Flows
With Known Analytical Solution

In this section we will analyzing the axisymmetric flow (see Table A.1) to validate the
coefficients derived by using MATLAB R code.
99

Table A.1: Analytical flows employed in validation of MATLAB R code used to derive Taylor
series coefficients
Case Flow
No.
1.
Axisymmetric Flow

A.2.1

Nature of flow
1. Steady flow.
2. Axisymmetric stagnation point flow.
3. The radial wall-shear stress and the wall
pressure gradient are linear.

Axisymmetric Flow

The following derivation may be found in most standard textbooks on fluid mechancis (e.g.
White [10], Section 3.8).
White [10] defines the streamwise velocity component in terms of similarity variable (η)
as:
vr =

Br
f (η)
n−1

(A.27)

And the wall normal velocity component in terms of similarity variables as:
√
vz = − Bνf (η)

(A.28)

And, the similarity variable is defines as:

η = −z

B
ν

(A.29)

Where z is wall-normal coordinate, r is streamwise coordinate, B is a constant, and ν is
kinematic viscosity.
100

Differentiating A.27 with respect to z:

∂vr
Br
=
∂z
n−1

B
f (η)
ν

(A.30)

∂ 2 vr
B2r
=
f (η)
(n − 1)ν
∂ 2z

(A.31)

Differentiating A.30 with respect to z:

Shear stress at the wall can be obtained by multiplying Equation A.30, and evaluating
it at the wall:
B
f (0)
ν

Br
τwall = µ
n−1

(A.32)

For axisymmetric flow (n = 3), τwall is given by:
B
τwall = µ
2

B
f (0)r
ν

(A.33)

Comparing Equations A.33 and A.34 to get the series for τwall :

∞

a0j rj =

j=0

B
f (0)r;
ν

B
4

(A.34)

Equation A.34 implies that the wall-shear stress is linear in ’r’ and gives a0j coefficients.

a00 = 0

a01 =

B
4
101

B
f (0)
ν

(A.35)

(A.36)

a02 , a03 , a04 , ...a0∞ = 0

(A.37)

Now, Equation A.20 when evaluated at the wall (z = ) gives:
∂pwall
∂ 2 vr
=µ
∂r
∂z 2

(A.38)

Combining Equations A.31 and A.38 (for n = 3, and at wall - y = 0):
∂pwall
µB 2
=
f (0)r
∂r
(2)ν

(A.39)

Combining Equations A.39 and A.23
∞

a1j rj = −

j=0

B2
f (0)r
12ν

(A.40)

Equation A.40 implies that the wall-prssure gradient is linear in ’r’ and gives a1j coefficients, as belowa10 = 0
a11 = −

B2
f (0)
12ν

a12 , a13 , a14 , ...a1∞ = 0

(A.41)
(A.42)
(A.43)

Having computed the coefficients using MATLAB R , we now exstimate series for vr and
vz . Figure A.1 and A.2 depict the series expansions up to 40th order of the stream function.
It is evident from the figure that the convergence is very slow. It implies that we need to a
polynomial series higher than 40th order to estimate the complete analytical profile.

102

3.5

Order 5
Order 10
Order 15
Order 20
Order 25
Order 30
Order 35
Order 40
Analytical

3

η

2.5
2
1.5
1
0.5
0
0

0.5
u/Uo

1

Figure A.1: Series expansions up to the 40th order of the radial velocity compared to the
analytical solution for axisymmetric stagnation flow

3.5

Order 5
Order 10
Order 15
Order 20
Order 25
Order 30
Order 35
Order 40
Analytical

3

η

2.5
2
1.5
1
0.5
0
−0.06

−0.04 −0.02
u/U0

0

Figure A.2: Series expansions up to the 40th order of the wall-normal velocity compared to
the analytical solution for axisymmetric stagnation flow
103

BIBLIOGRAPHY

104

BIBLIOGRAPHY

[1] Speijker, N., Vidal A., Barbaresco F., Gerz, T., Barny, H., Winckelmans, G. (2005)
ATC-Wake: Integrated Wake Vortex Safety and Capacity System (ATC Wake D6 2),
IST2001-34729.
[2] Dallmann, U. and Gebing, H. (1994) Flow attachment at flow separation lines: on uniqueness problems between wall-flows and off-wall fields, Acta Mechanica [Suppl] 4, 47.
[3] Gendrich, C. P., Bohl, D., and Koochesfahani, M. M. [1997] Whole field measurements of
unsteady separation in a vortex ring/wall interaction, AIAA Paper No. AIAA-97-1780.
[4] Wu, J.Z., Wu, J.M. and Wu, C.J. (1987) A general three-dimensional viscous compressible theory on the interaction of solid body and vorticity-dilatation field, University of
Tennessee Space Institute Report 87/031.
[5] D. Fabris, D. Liepmann, and D. Marcus, Quantitative experimental and numerical investigation of a vortex ring impinging on a wall, Phys. Fluids 8, 2640 1996.
[6] Chu, C.-C.,Wang, C.-T. and Chang, C.-C. [1995] A vortex ring impinging on a solid
plane surface vortex structure and surface force, Phys. Fluids, 7(6), 1391
[7] Orlandi, P. and Verzicco, R. [1993], Vortex rings impinging on walls: axisymmetric and
threedimensional simulations, J. Fluid Mech., 256, 615
[8] Naguib, A. and Koochesfahani, M. M. (2004) On wall-pressure sources associated with
the unsteady separation in a vortex-ring/wall interaction, Phys. Fluids, Vol. 16, No. 7,
2613-2622.
[9] Panton, R. L. (2005) Incompressible Flow Third Edition. New Delhi, India: Wiley India
(P) Ltd.
105

[10] White F.M. (1991), Viscous Fluid Flow, McGraw Hill International.
[11] Lele SK (1990) Compact finite difference schemes with spectral-like resolution, CTR
Manuscript 107, NASA Technical Reports.
[12] Etebari A, Vlachos PP (2005) Improvements on the accuracy of derivative estimation
from DPIV velocity measurements. Exp Fluids 39: 1040-1050.
[13] http://www.aerospaceweb.org/question/nature/q0237.shtml

106