KTH Engineering Sciences

Drag reduction using plasma actuators

Romain Futrzynski

Licentiate Thesis
Stockholm, Sweden
2015
Academic thesis with permission by KTH Royal Institute of Technology,
Stockholm, to be submitted for public examination for the degree of
Licentiate in Engineering Mechanics, Friday the 27th of March, 2015
at 10.00, in the Vehicle Engineering Lab, Teknikringen 8, KTH - Royal
Institute of Technology, Stockholm, Sweden.
TRITA-AVE 2015:10
ISSN 1651-7660
ISBN 978-91-7595-479-0

c Romain Futrzynski, 2015

Postal address: Visiting address: Contact:

KTH Teknikringen 8 [email protected]
Farkost och Flyg
SE-100 44 Stockholm
Abstract
This thesis is motivated by the application of active flow control on the
cabin of trucks, thereby providing a new means of drag reduction. Partic-
ularly, the work presented strives to identify how plasma actuators can
be used to reduce the drag caused by the detachment of the flow around
the A-pillars. This is achieved by conducting numerical simulations, and
is part of a larger project that also includes experimental.
The effect of plasma actuators is modeled through a body force, which
adds very little computational cost and is suitable for implementation
in most CFD solvers. The spatial distribution of this force is described
by coefficients which have been optimized against experimental data,
and the model was shown to be able to accurately reproduce the wall jet
created by a single plasma actuator in a no-flow condition.
A half cylinder geometry—a simplified geometry for the A-pillar of a
truck—was used in a preliminary Large Eddy Simulation (LES) study
that showed that the actuator alone, operated continuously, was not
sufficient to achieve a significant reduction of the drag. Nevertheless, a
significant drag reduction was obtained by simply increasing the strength
of the body force to a higher value, showing that this type of actuation
remains relevant for the reduction of drag.
In the course of finding ways to improve the efficiency of the actuator,
dynamic mode decomposition was investigated as a post-processing tool
to extract structures in the flow. Such structures are identified by their
spatial location and frequency, and might help to understand how the
actuator should be used to maximize drag reduction. Thus a parallel code
for dynamic mode decomposition was developed in order to facilitate the
treatment of the large amounts of data obtained by LES. This code and
LES itself were thereafter evaluated in the case of a pulsating channel flow.
By using the dynamic mode decomposition it was possible to accurately
extract oscillating profiles at the forcing frequency, although harmonics
with lower amplitude compared to the turbulence intensity could not be
obtained.

Keywords: flow control, drag reduction, plasma actuator, DMD, LES,

optimization, pulsating flow

iii
Sammanfattning
Denna avhandling behandlar tillämpningen av aktiv strömningskontroll
för lastbilshytter, vilket är en ny metod för minskning av luftmotståndet.
Mer i detalj är det övergripande målet att visa på hur plasmaaktuatorer
kan användas för att minska luftmotståndet orsakat av avlösningen runt
A-stolparna. In denna avhandling studeras detta genom numeriska
simuleringar. Arbetet är en del av ett projekt där även experimentella
försök görs.
Effekten av plasmaaktuatorer modelleras genom en masskraft, vilket inte
ger nämnvärd ökning av beräkningstiden och är lämplig för implementer-
ing i de flesta CFD-lösare. Den rumsliga fördelningen av kraften bestäms
av koefficienter vilka i detta arbete beräknades utifrån experimentella
data. Modellen har visat sig kunna återskapa en stråle nära väggen med
god noggrannhet av en enskild plasmaaktuator för en halvcylinder utan
strömning.
Samma geometri—en halvcylinder som här används som förenklad geo-
metri av A-stolpen på en lastbil—användes i en preliminär LES studie
som visade att enbart aktuatorn vid kontinuerlig drift inte var tillräckligt
för att uppnå en signifikant minskning av luftmotståndet. En signifikant
minskning av luftmotståndet erhölls genom att helt enkelt öka styrkan
på kraften, vilket visats att denna typ av strömningskontroll är relevant
för minskning av luftmotståndet.
I syfte att förbättra effektiviteten hos aktuatorn, studerades dynamic
mode decomposition, som ett verktyg för efterbehandling för att få fram
flödesstrukturer. Dessa strukturer identifieras genom deras rumsupplös-
ning och frekvens och kan hjälpa till att förstå hur aktuatorerna bör an-
vändas för att minska luftmotståndet. En parallelliserad kod för dynamic
mode decomposition utvecklades för att underlätta efterbehandlingen
av de stora datamängder som fås från LES-beräkningarna. Slutligen,
utvärderades denna kod och LES-beräkningar på ett strömningsfall med
pulserande kanalflöde. Metoden, dynamic mode decomposition, visade
sig kunna extrahera de oscillerande flödesprofilerna med hög noggran-
nhet för den påtvingade frekvensen. Övertoner med lägre amplitud
jämfört med turbulensintensiteten kunde dock inte erhållas.

Nyckelord: strömningskontroll, motståndsminskning, plasma ställdon,

DMD, LES, optimering, pulserande flöde

iv
Acknowledgements
The work was financially supported by the Swedish Energy Agency
within the project Flow Research on Active and Novel Control Efficiency
(FRANCE), project number 34186-1.
The computations were performed on resources provided by the Swedish
National Infrastructure for Computing (SNIC) at the High Performance
Computing Center North (HPC2N), National Supercomputer Centre
at Linköping University (NSC) and PDC Center for High Performance
Computing.

Per Elofsson and Guillaume Mercier from Scania AB are gratefully thanked
for their contribution to the project. Thanks are also owed to Julie Vernet,
Ramis Örlü and P. Henrik Alfredsson for making this project more than
isolated simulations and for providing valuable ideas and comments for
my work.

I would of course like to thank Gunilla Efraimsson for offering me to

work on this fascinating project and for supervising my research. I really
appreciate your enthusiasm, ideas and availability. You have always been
helpful in situations ranging from dealing with “the hierarchy” to getting
out of research dead ends. Thank you too for taking care of making such
a nice environment to work in, it is a pleasure to come to work every day!
Naturally having such a good workplace is only possible because of the
people. A heartfelt thank you goes to everyone at the department, office-
mates, neighbors, PhDs, Masters, Post-docs, Professors, lunch-group
members or simply acquaintances, past and present, for providing stimu-
lating discussions during weekdays, and fine company outside of work.

Thanks also to my friends in Stockholm, and those who have since moved,
for making it a cool place to live in.

Romain Futrzynski
Stockholm, 4th March 2015

v
Dissertation
This thesis consists of two parts: The first part gives an overview of
the research area and work performed. The second part contains the
appended research papers 1–3. The work was divided between authors
as follows:

Paper 1
Effect of a SDBD on the Drag of a Half-Submerged Cylinder in Crossflow.
R. Futrzynski, G. Efraimsson. In ASME 2014 4th Joint US-European Flu-
ids Engineering Division Summer Meeting collocated with the ASME
2014 12th International Conference on Nanochannels, Microchannels,
and Minichannels (peer-reviewed).

R. Futrzynski created the computational meshes, performed the com-

putations, discussed the results and wrote the paper together with G.
Efraimsson.

Paper 2
Dymode, A parallel dynamic mode decomposition software.
R. Futrzynski, G. Efraimsson. Internal report, KTH 2015. TRITA-AVE
2014:78, ISSN 1651-7660, ISBN: 978-91-7595-386-1

R. Futrzynski wrote the dymode code, discussed various computational

methods for DMD, described the particular implementation of dymode
and wrote the report together with G. Efraimsson.

Paper 3
Numerical study of the Stokes layer in oscillating channel flow.
R. Futrzynski, C. Weng, S. Boij, G. Efraimsson, A. Hanifi. To be submitted.

R. Futrzynski created the computational mesh, performed the computa-

tion of the LES and the DMD. C. Weng performed the DNS and provided
the literature review. R. Futrzynski discussed results and wrote the paper
together with C. Weng, G. Efraimsson, S. Boij, and A. Hanifi

vii
Publications not included in this thesis
Conference papers
Numerical simulation of a plasma actuator on a half- submerged cylinder.
R. Futrzynski, G. Efraimsson, P. H. Alfredsson. In 4th International
Conference on Jets, Wakes and Separated Flows, September 17-21, 2013,
Nagoya, Japan.

Aero-acoustic source analysis of landing gear noise via dynamic mode decompos-
ition.
J. Dahan, R. Futrzynski, C. O’Reilly, G. Efraimsson. In Proceedings of the
21st International Congress on Sound and Vibration: In Depth Sound and
Vibration Research, 2014 / [ed] Malcolm J. Crocker, Marek Pawelczyk,
Jing Tian, 2014 (reviewed).

viii
Contents

I OVERVIEW 1
1 Introduction 3
1.1 Truck aerodynamics . . . . . . . . . . . . . . . . . . . . . . 3
1.2 Active flow control . . . . . . . . . . . . . . . . . . . . . . . 4
1.3 Contents of the thesis . . . . . . . . . . . . . . . . . . . . . . 6

2 Numerical models 7
2.1 Flow simulations . . . . . . . . . . . . . . . . . . . . . . . . 7
2.1.1 Large Eddy Simulations . . . . . . . . . . . . . . . . 7
2.1.2 Numerical solver . . . . . . . . . . . . . . . . . . . . 9
2.2 Plasma actuators . . . . . . . . . . . . . . . . . . . . . . . . 10
2.2.1 The dielectric barrier discharge . . . . . . . . . . . . 10
2.2.2 Overview of plasma models . . . . . . . . . . . . . . 12
2.2.3 The exponential model . . . . . . . . . . . . . . . . . 13

3 Drag reduction study 15

3.1 The model problem . . . . . . . . . . . . . . . . . . . . . . . 15
3.2 Model coefficients . . . . . . . . . . . . . . . . . . . . . . . . 17
3.2.1 The modified test case . . . . . . . . . . . . . . . . . 17
3.2.2 Optimization of the coefficients . . . . . . . . . . . . 20
3.3 Effect of the actuator . . . . . . . . . . . . . . . . . . . . . . 21
3.3.1 Numerical setup . . . . . . . . . . . . . . . . . . . . 21
3.3.2 Turbulence resolution of the un-actuated flow case 24
3.3.3 Drag reduction of the actuated flow case . . . . . . 26
3.4 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28

ix
CONTENTS

4 Tools for flow improvement 29

4.1 Dynamic mode decomposition . . . . . . . . . . . . . . . . 30
4.2 LES of pulsating channel flow . . . . . . . . . . . . . . . . . 32
4.2.1 Numerical setup . . . . . . . . . . . . . . . . . . . . 33
4.2.2 Mean profiles . . . . . . . . . . . . . . . . . . . . . . 36
4.2.3 Oscillating profiles . . . . . . . . . . . . . . . . . . . 38

5 Conclusions and outlook 43

6 Summary of appended papers 45

Bibliography 47

II APPENDED PAPERS 1-3 53

x
 Part I

OVERVIEW
 1 Introduction

This chapter offers basic understanding of the driving forces behind the
research performed in this thesis. The incentive for this work is to solve
practical development problems in industry applications, however, the
research performed is general and independent of application.

1.1 Truck aerodynamics

A lower fuel consumption can lead to significant savings for owners of
fleets of long-haul vehicles. Furthermore, it will also reduce the amount
of particulates and greenhouse gas for which public concern has greatly
grown since the invention of automobiles.
For a 40 tonnes tractor plus semi-trailer driving at 90 km/h on a
flat road, it can be estimated that one third of the fuel consumption is
due to rolling resistance, one third is due to aerodynamic drag, and
the remaining is due to factors such as engine efficiency or equipment
powering as illustrated in fig. 1.1. It can also be estimated that about
90% of the aerodynamic losses are caused by form drag, with air friction
causing only 10% of the losses. As such, aerodynamics represent a key
aspect of low fuel consumption, and receives high attention from truck
manufacturers.
Indeed, a large number of techniques to reduce aerodynamic drag
have been developed and improved between each generation of trucks.
The most natural one being to streamline the shape of the truck. For
instance, on a bus-type vehicle, replacing sharp leading corners with
rounded ones may decrease the aerodynamic drag coefficient as mush
as 60% [1]. Other techniques include the use of roof and side deflectors

3
CHAPTER 1. INTRODUCTION

3% Friction drag

33% Other losses

30% Form drag

33% Rolling resistance

Figure 1.1: Causes of fuel consumption for a 40 tonnes truck driving at 90 km/h.

(as seen on the R640 shown in fig. 1.2, which ease the transition from the
tractor unit to the trailer of the truck. Guiding vanes may also be utilized.
The use of skirts, wheel fairings or even a boat tail shape at the end of the
trailer could also contribute to drag reduction, although they are seldom
used for practical or legal reasons. However, these techniques could be
called passive, since they are fixed during the design phase and cannot
adapt to real-time flow conditions. For instance, the radius of the leading
corners on a Scania truck is optimized to produce low drag over a range
of side wind angles, each angle being assumed to happen for a certain
fraction of the time as described in the SAE standard [2]. Furthermore,
passive techniques may conflict with other aspects of truck design, such
as clearances or structural integrity. In the case of the rounded leading
corners, for example, space is lost that could have been used for the
underhood systems.

1.2 Active flow control

A different set of techniques has gained more and more interest over
the last decades, generally referred to as active flow control. This time,
various devices are used at key locations of a flow in order to change
its proprieties. These devices are active only when powered up by an

4
 1.2. ACTIVE FLOW CONTROL

Figure 1.2: Scania trucks LBS 140 (1969) and R640 (2010). A number of techniques have
been used on the more recent model in order to reduce drag.

energy source, and thus may be modulated to provided the best output
relatively to the instantaneous flow.

Commonly researched devices for active flow control include suction

and/or blowing orifices or slots [3–6], synthetic jets, sometimes referred
to as zero net mass flux jets [7, 8], or microjet arrays [9]. A compilation of
flow control research is presented in ref. [10]. Applications range from
suppression of vortex shedding, separation control on turbine blades
or turbulence transition control. A large part of the ongoing research
is dedicated to applications for airfoils, such as lift increase and stall
delay, as well as to ground vehicles. Another type of devices are plasma
actuators, which have in the recent years received a lot of attention since
they present extremely fast reaction time, nonexistent risk of mechanical
failure and low power consumption. Such actuators are discussed in
more details in Chapter 2.

Usage of active flow control is also pushed forward by the devel-

opments in control theory for fluid mechanics [11, 12], numerical simu-
lations and reduced order modeling [13–15], as well as sensing techno-
logy [16]. Together, these disciplines enable the application of closed-loop
actuation to adapt control in real time and obtain the best performance
for any instantaneous flow.

5
CHAPTER 1. INTRODUCTION

1.3 Contents of the thesis

The thesis begins with a description of the equations used for flow mod-
eling and the numerical schemes employed to solve them. Further, a
description of plasma actuators is given. Several models of plasma actu-
ators are first discussed, and thereafter the plasma actuator model used
in this work is presented in more details.
In Chapter 3, this model is used in a set-up designed to resemble the
A-pillars of a truck. The model of actuator is first used in a no-flow
configuration and compared to experimental data in order to optimize
the value of its coefficients. Then, it is used in a crossflow configuration,
where its effect on the drag of a half-submerged cylinder is measured.
The next chapter describes dynamic mode decomposition and how it
can be used to help make the actuation more efficient. The scientific
code written to compute such mode decompositions for large datasets is
thereafter used in a pulsating channel flow case, and the results compared
with profiles phase averaged at known frequencies.
Conclusions from the work performed are made in Chapter 5, where
directions for the continuing of the research are given as well. Finally, a
summary of the main papers produced is given in Chapter 6.

6
 2 Numerical models

Simulating the physical world on a computer requires mathematical

equations describing its phenomena, and a way to solve these equations
using binary logic. This chapter describes how the turbulent fluid flows
encountered in the scope of this thesis were solved. In addition, a de-
scription of plasma actuators is provided together with a practical way
to incorporate their effect into the flow solution.

2.1 Flow simulations

2.1.1 Large Eddy Simulations
In the simulations presented in this thesis, the incompressible Navier-
Stokes equations are solved:

 ∂ui = 0

 x

i
(2.1)
 ∂ui ∂ui u j 1 ∂p ∂Sij

 + = − + 2ν + f i
 ∂t xj ρ ∂xi xj
 ∂u j

where Sij = 12 ∂u i
x j + xi is the strain-rate tensor and f i is the body force
per mass unit acting on the flow.
Like the RANS equations [17], the equations solved in LES are ob-
tained by applying a filtering operation to the Navier-Stokes equations.
However it is different in that instead of being a temporal filtering, it
is a spatial filtering that transforms quantities at one position into its
weighted average over a volume surrounding that point. In that sense,

7
CHAPTER 2. NUMERICAL MODELS

turbulent length scales that are smaller than the characteristic length of
the averaging volume are lost, while larger length scales remain. In prac-
tice, a top hat filtering function in combination with the computational
grid itself is often used. This is very convenient since averaged values
over each cell are already what finite volume solvers use, and therefore
the computation of a 3-dimensional convolution product can be avoided.
This type of “implicit filtering” was used for the LES simulations presen-
ted here.
Again, similarly to the RANS filtering, the filtered term ui u j cannot
be expressed solely using the filtered quantities ui . It is therefore split
into two parts as ui u j = ui u j + τij , where τij is a residual stress tensor
relative to subgrid-scales. Hence (2.1) can be rewritten as

∂u
 i =0

 xi

. (2.2)
 ∂ui ∂ui u j 1 ∂p ∂Sij ∂τij

 + = − + 2ν − + f i
xj xj

∂t ρ ∂xi ∂x j

An eddy-viscosity model is used to close the equations, where the trace-

less part of τij is proportionally related to the strain rate via the turbulent
viscosity νt such that
δij
τij − τ = −2νt Sij . (2.3)
3 kk
The coefficient νt can be related to a characteristic length and velocity of
the subgrid-scale eddies. However, these characteristic scales have to be
modeled since such eddies are not resolved. In the LES performed in this
thesis, two models have been used.
The Smagorinsky model [18] was used in Paper 3. It defines νt as
  2  
νt = f v min κd, Cs V 1/3 S , (2.4)

where V is the cell volume, Cs = 0.1, d is the wall distance, κ = 0.41 is

the Von Karman constant, and f v is the Van Driest damping function [19]
computed as
 +
y
f v = 1 − exp − . (2.5)
A
Here, y+ is the dimensionless wall distance normalized with the viscous
length scale at the nearest wall face, and the constant A equals 25.

8
 2.1. FLOW SIMULATIONS

The WALE model (Wall-Adapting Local-Eddy Viscosity) [20] was used

in Paper 1. This time, νt is computed as
  2
νt = min κd, Cw V 1/3 Sw , (2.6)

where Cw = 0.544, Sw is a deformation parameter and other symbols are

identical to those in the Smagorinsky-Lilly model. Introducing the tensor
 !2  
2
∂uk 2
 
d 1 ∂u i ∂u j 1
Sij =  +  − δij , (2.7)
2 ∂x j ∂xi 3 ∂xk

the deformation parameter is defined as

 3/2
Sijd Sijd
Sw =  5/4 . (2.8)

5/2
Sijd Sijd + Sij Sij

Unlike the Smagorinsky-Lilly model, the WALE model has the advantage
of accounting for the effect of small scale strain rates as well as rotation
rates, and producing a value of νt that inherently goes to zero at the wall
following a more physically accurate cubic law.

2.1.2 Numerical solver

All the simulations performed in this thesis used the finite volume solver
Star-CCM+, versions 8 to 9. The only exception are the DNS results in
Paper 3 which were obtained using the spectral solver SIMSON [21] and
was performed by co-author C. Weng.
An implicit scheme with second-order accuracy was used to advance
all the simulations in time throughout the thesis. The time steps were
fixed, and were set such that they at least satisfy that the convective
Courant number remains below one throughout the domain and over
the entire simulation time, that is,

kuk ∆t
< 1. (2.9)
V 1/3
Here kuk is the norm of the velocity vector and V is the cell volume. Even
though an implicit solver is used, this Courant number criteria helps to

9
CHAPTER 2. NUMERICAL MODELS

make sure that no time-filtering of the turbulence occurs, which would

reduce the amount of the resolved turbulence [22].
A central differencing scheme was used in space in Paper 3 and a
bounded central differencing scheme was used in Paper 1. The bounded
central differencing scheme blends the central differencing and the second
order upwind schemes according to a local flow values in order to bal-
ance accuracy and robustness. It can also revert to a first order upwind
scheme when the convection boundedness criterion is not satisfied. As a
result, it may be more dissipative than the central differencing scheme.
All the flows were solved in a segregated manner using the SIMPLE
algorithm, and a Gauss-Seidel relaxation scheme.

2.2 Plasma actuators

2.2.1 The dielectric barrier discharge
Creating a plasma at atmospheric pressure is a remarkably easy task.
Bring two electrodes connected to a high voltage generator close enough
to each other, and voila, plasma will form. All that is required is that
the electric field be high enough to cause the electric breakdown of the
ambient air atoms; that is, to dissociate them into a globally neutral
volume of freely moving ions and electrons. The electrodes may come in
any shape or form as long as the electric field created is strong enough.
However, there is a particular configuration in which the electrodes form
a so-called dielectric barrier discharge (DBD), that has been particularly
widely used. In this configuration, two flat electrodes, one larger than
the other, are separated by an insulating dielectric layer as illustrated
in fig. 2.1. This is what is referred to as the “plasma actuator” hereafter.
The actuator will generally be embedded within a surface, leaving only
one exposed electrode. By applying an alternative high voltage between
them, a plasma can be sustained. Reasons for operating in an alternative
way rather than in a direct way include a lower breakdown voltage
requirement and the lack of a strong continuous current that would lead
to higher power consumption and electrode corrosion [23]. The size of
the electrodes is typically of the order of the centimeter for flow actuation
purposes, and the voltage used is in the range of 1–10 kV at an alternative
frequency of 1–20 kHz.
For the time when the plasma exists, momentum is added to the
surrounding air. While the physics of the phenomenon is not clearly

10
 2.2. PLASMA ACTUATORS

Actuator location reference

Induced flow
plasma

dielectric

Figure 2.1: Sketch of a single dielectric barrier discharge (SDBD) embedded in a solid wall.

understood, it is generally admitted that the charged particles subjected

to the electric field transfer their momentum to the adjacent non-ionized
fluid by collision. Interestingly, this momentum seems to always be
directed from the exposed electrode to the embedded one, even though an
alternative voltage is applied. This is usually admitted to be caused by the
asymmetry between the two electrodes [24, 25]. As a result, the actuator
can be utilized to create a wall jet downstream of the exposed electrode,
making it a serious candidate for active flow control applications.

Indeed, unlike zero-net mass flux devices, plasma actuators have no

moving parts, making them less prone to failure. They do not require any
network of pipes to provide pressurized air like suction or blowing slots
do. In fact, they are quite small and may simply be laid on the surface
where actuation is required with minimal effort. Furthermore, their
reaction time is orders of magnitude faster than the time scale of the flow.
However, a high voltage source is required and sustaining the plasma
produces ozone. Moreover, the exposed electrode can be potentially
dangerous due to the high voltage involved. Finally, the actuators as
built in ref. [26, 27] have been observed to burn out after an hour or so of
operation, so that better materials of manufacturing techniques should be
employed for commercial applications. Nevertheless, plasma actuators
are being investigated thoroughly for applications as wide as separation
control for wings [28–31] and other bluff bodies [32] or transition to
turbulence control [33].

11
CHAPTER 2. NUMERICAL MODELS

2.2.2 Overview of plasma models

Detailed simulations of plasma actuators are particularly computation-
ally expensive, in that they have to solve, in addition to the Navier-Stokes
equations, Maxwell’s equations to take into account the electrodes and
the time-dependent space distribution of charged particles within the
plasma. A complication is that the time scales for these processes are or-
ders of magnitude smaller than the time scales of the ambient flow alone.
Even when just three generic species are considered—ions, electrons and
neutrals, as opposed to a declination of these for each gas making up
common air—within a 2-dimensional domain, the computational cost
remain out of reach for a use in optimization of more practical cases
[34]. Boltzman-based simulations have also been employed bringing
considerable insight to the understanding of the fundamental physics,
but again the cost is too high for practical applications [35].
Oftentimes the interest of a simulation is not on the plasma itself,
but rather on how it can be used to achieve flow control. Thus, it is
sufficient to have a model of the effect of the actuator on the flow, so that
the underlying electrohydrodynamic problem can be omitted and the
problem complexity falls back to solving the Navier-Stokes equations
for the flow. An approach for modeling this effect is the use of a body
force acting on the fluid. This spawned a very popular branch of models
since it is faithful to the idea that momentum is transmitted by collision
between the charged particles of the plasma and the ambient fluid. A
first simplification consists in computing the electric field E(x, t) between
the two electrodes assuming a constant charge density q(x).
A body force to be applied on the fluid can then be formulated as

f (x, t) = q(x) · E(x, t). (2.10)

In this way the force acts as a model for the charged particles following
the electric field lines and colliding with neutral species, but the particles
themselves and their interactions does not have to be simulated. How-
ever, E(x, t) created by the electrodes still needs to be computed using
electro-magnetic equations, increasing the problem complexity. There-
fore, additional models for E have been proposed. A model of plasma
actuator proposed in ref [36] proposes to split the volume directly above
the actuator into N parallel slices. Each slice is then considered in par-
allel as an electrical circuit having certain dissipative and capacitance
properties. The electrical analogy permits to compute an approximation
of E(x, t) using less resources.

12
 2.2. PLASMA ACTUATORS

A further simplification can be made by considering that the time-

dependent phenomena occurring within the plasma repeat at a frequency
of the order of the kilo-Hertz which is much higher than any time scale
of the flow. Hence the body force experienced by the flow can reasonably
well be modeled as a constant related to the mean force exerted over a
period. In a model given in ref [24], a geometric simplification of the
electric field lines are used to approximate the average direction of the
force as shown in fig. 2.2. The strength of the force is then assumed to
decrease linearly from a maximum at the edge of the exposed electrode,
so that

f (x) = ( f (0) − k1 x − k2 y) E, (2.11)

where k1 and k2 are two constants defining how fast the body force
decreases away from the actuator in each direction. Yet a further simpli-

Figure 2.2: Sketch of the electric field lines between the electrodes (left), approximated as
oblique parallel lines of slope a/b (right).

fication can be made by only keeping the component of the body force
that is tangential to the surface [37].

2.2.3 The exponential model

The purpose of the research presented in this thesis is to numerically
study how plasma actuators can be used to decrease the drag of road
vehicles. This requires detailed simulations of the flow field at high
Reynolds numbers, which are in themselves computationally costly. Sim-
ulating the physics of the plasma, on the other hand, is not required as
long as the effect of the actuator can be reproduced accurately. Therefore,
numerical simulations have been conducted using a body force model

13
CHAPTER 2. NUMERICAL MODELS

which includes the simplifications of time independence and single non-

zero component tangent to the surface. The force is considered to be
maximum at the reference location of the actuator, which is the edge of
the exposed electrode. However, an exponential decay of the force is
used in both directions away from the actuators, as it was observed in ex-
perimental work. The force is null upstream and below the actuator. Such
a model has previously been successfully compared with experimental
results [31].
Since in the present thesis, the actuator is always located on a half
cylinder, we define the Cartesian coordinate system (Ω, τ, r) for any angle
θ where the actuator might be placed, as illustrated in fig. 2.3. The angle
θ identifies the location of the trailing edge of the exposed electrode. By

τ
Ω

Figure 2.3: Coordinate system for the exponential model.

calling (x; y) the coordinates of a point in this coordinate system, the body
force used to model the effect of the plasma actuator is expressed as
(
Fmax e−ξx e−ζy · τ if (x; y) ≥ (0; 0),
f = (2.12)
0 otherwise.

The coefficient Fmax is the strength of the body force at the edge of the
electrode, and ξ and ζ are constant defining how fast the force decays in
the tangential and wall-normal directions, respectively. Here, the values
of ξ and ζ will not be derived from electro-magnetic considerations, but
an optimization method will be used to determine them empirically, see
Chapter 3.

14
 3 Drag reduction study

In this chapter the flow and plasma models are used to quantify the drag
reduction achievable on a simplified geometry. The optimization method
used to obtain empirical model coefficients is first described. Results for a
half-cylinder in crossflow using the calibrated model are then discussed.

3.1 The model problem

Even though plasma actuators may be utilized for many applications,

they are here studied in the context of truck cabins. Sources of aerody-
namic drag caused by the cabin include side mirrors, wipers, air intake
and leading edges when they cause flow detachment. Rounding the
leading edges is a quite efficient way of reducing drag [1], and in fact
edges and corners of current truck models have curvature radii designed
to minimize drag in a variety of wind conditions. However, there are two
main issues with rounded edges: The first issue is that they still relate to
passive techniques, and so the drag reduction obtained may not be the
absolute maximum for any wind configuration. The second issue is that
rounding the edges reduces the space available in-hood, which might
otherwise have been used to improve other areas of the cabin design.
Therefore the leading edges, and the A-pillars in particular (see fig. 3.1),
have been chosen as the target for the actuators.
Thus the goal is set to determine a way to use plasma actuators located
on the A-pillars of the truck in order to reduce the aerodynamic drag.
Using a full-scale truck cabin driving at 90 km/h would give a Reynolds

15
CHAPTER 3. DRAG REDUCTION STUDY

Figure 3.1: Truck cabin with left A-pillar marked in red. It is the target location selected for
the use of plasma actuators.

number of approximately

25 × 2
Re = ≈ 3 × 106 , (3.1)
1.5 × 10−5
at which performing either simulations or experiments would present
a challenge in itself—even if removing the fine geometry details, and
may not be suitable to the study of the plasma actuators. Furthermore,
separation on the roof and floor edges may occur and even interfere with
the separation around the A-pillars, causing the effect of the actuators to
be difficult to measure.
A simpler case that would still resemble a truck geometry is a forward
facing step. The step could be thought of as a half truck, infinitely high,
with the A-pillar being the edge of the step. However, this case requires
a number of parameters to be set, such as the height and length of the
step. For instance, should the step be as long as the truck cabin plus
trailer? And if so, should the flow at the back of the trailer be part of the
study? Moreover, one can expect that by using a sharp step edge, flow
separation will always occur and that actuators would have very low
impact. What, then, should the curvature radius be? All these questions
would be better answered at a later stage of the truck design, when the
actuators and their effects are better understood.
Further simplification of the forward facing step brings the case of the

16
 3.2. MODEL COEFFICIENTS

flow over a hill. In order to limit the number of parameters necessary for
the study, the shape of the hill is taken to be half a circle. With such a
shape, separation is certain to occur downstream of the hill, so that the
effect of the actuator can be measured.

3.2 Model coefficients

3.2.1 The modified test case
As discussed in Chapter 2, the numerical simulations require a model for
the actuators. The exponential model has the advantage of a very light
computational cost, and it will thus be used for the present study. As it
can be seen in the expression of the model (2.12), the model requires val-
ues for a few coefficients defining the strength and reach of the actuator.
The value of these coefficients was obtained by comparing simulations
and experiments on a no-flow set-up of the half-cylinder case. In this set-
up, the half-cylinder is placed in a closed test-section, ensuring that the
only source of momentum within the test section is the plasma actuator.
The actuator itself is placed at the top of the half-cylinder, occupying the
whole span of the test-section. When the actuator is turned on, a wall
jet is created along the surface of the half-cylinder downstream of the
actuator, see fig. 3.2.
Velocity profiles of the wall jet were measured experimentally using
laser Doppler velocimetry (LDV) at several locations downstream of the
actuator, and have been used as reference for the determination of the
model’s coefficients. They correspond to the jet from an actuator run at
12 kV, 1 kHz as described in [26]. These profiles reveal that the maximum
of the velocity within the jet is of the order of 4 m/s. This velocity is
low enough so that the flow within the test section is considered to be
laminar. By neglecting the effect of the sides of the test section, this flow
is also considered to be 2-dimensional.
Therefore, the simulations for this particular case will use RANS
modeling and a 2-dimensional domain. The unstructured mesh used
was constituted of square cells that were trimmed at the boundaries
to fit the geometry. The surface of the cylinder was covered with a
prismatic layer, which provides a hyperbolic tangent growth of the cells
in the wall-normal direction. In addition, refinement was used at and
downstream from the position of the actuator, see fig. 3.3. The reason
for this refinement is that both the body force distribution (having an

17
CHAPTER 3. DRAG REDUCTION STUDY

Figure 3.2: Wall jet and recirculation driven by the plasma actuator, inside the closed
test-section.

exponential decay in space) and the high velocity region of the jet present
large gradients in space. The refinement is thus intended to provide a
better resolution for both of them.
The same approach was used in [31] to calibrate the coefficients of a
similar, albeit simplified, exponential model. The calibration was done
for a plasma actuator on a flat plate and comparing velocity profiles at a
single position downstream of the actuator. The model was simplified
by observing that the plasma region is smaller in the wall-normal than
in the tangential direction, and thus setting ζ = 2ξ. In this study, it was
reported that values of Fmax = 70 × 103 N/m3 and ξ = 3.25 × 103 m-1
(summarized as Original in table 3.1) provided a good agreement between
the experimental and numerical velocity profiles. However, when used
in the half-cylinder case, using the same values produced a wall jet much
stronger than the one measured experimentally. This is not an alarming
result, since there is no standard way to build a plasma actuator, and

18
 3.2. MODEL COEFFICIENTS

Figure 3.3: 2D mesh close to the actuator’s position at the top of the cylinder.

it is likely that the actuators from the two studies would have different
characteristics. Perhaps even more importantly, the operating conditions
were not the same in the two studies, as the actuator in ref. [31] was run
at 10 kV, 4 kHz. A naive solution to this issue is to simply halve the Fmax
coefficient in order to reduce the velocity of the jet (corresponding to
Modified in table 3.1). Unfortunately, it can be seen in fig. 3.4 that the
resulting jet profiles display conflicting features: On the one hand, the
peak velocity of the profile closest to the actuator is comparable with
the experimental data, while the peak of the profile furthest from the
actuator is about twice as large as in the experiment. On the other hand,
the location of the peak in the closest profile is about twice as high as
in the experiment, while further from the actuator the heights of the
peaks seem to coincide. Hence it becomes clear that in order to reproduce
the correct jet profiles, all of the Fmax , ξ and ζ values have to be tuned
together.
Moreover, it was observed during the experiments that the plasma
region does not start immediately at the edge of the exposed electrode,
but that some distance separates this edge to the most intense plasma
location. To reflect this behavior, even though the location of the actuator
is θ = 90◦ , which is identified by the edge of the exposed electrode, the
most upstream location where the body force is applied is θ = 90◦ + δθ,
where δθ is an additional unknown to our model. This makes four
tangled coefficients that have to be determined. Proceeding with a trial

19
CHAPTER 3. DRAG REDUCTION STUDY

5
4
3
2
1
0
−1
−2
−3
Experimental data
−4
Exponential model
−5
5 10 15 20

Figure 3.4: Experimental and numerical velocity profiles created by the actuator in the
no-flow test case. One unit equals one tenth of the diameter for lengths, and two meters
per second for velocities.

and error method by hand would prove tedious at best. Nevertheless,

the 2-dimensional, laminar simulations are not particularly costly so that
a trial and error method can be undertaken if proper automating tools
are used.

3.2.2 Optimization of the coefficients

The body force model can be re-written as
(
Fmax e−ξxθ+δθ e−ζyθ+δθ · τ θ+δθ if (xθ+δθ ; yθ+δθ ) ≥ (0; 0),
f = (3.2)
0 otherwise.

While θ is known and corresponds to the angle where the actuator is

located, the four parameters Fmax , ξ, ζ and δθ were iteratively guessed
and tested by running a new simulation. The agreement between experi-
mental data and one of these simulations was evaluated by calculating
the quantity
s
g = ∑ uexp (xi ) − usim (xi )

2
(3.3)
xi

which is the root mean square of the velocity difference with the sim-
ulation at every point xi for which there was experimental data of the

20
 3.3. EFFECT OF THE ACTUATOR

velocity available. With starting values corresponding to the Modified

value set, the SHERPA (Simultaneous Hybrid Exploration that is Robust,
Progressive, and Adaptive) algorithm [38] from the HEEDS software
package was used to iteratively find better and better values that min-
imize g. The range of allowed values was set to be extremely large, and
simulation where normalized residuals did not decrease below 10−3 were
discarded. The algorithm was allowed to run 400 simulations to find the
best set of coefficients. Figure 3.5 shows how the performance in min-
imizing g increase when one of the parameters, δθ, converges towards a
slightly negative value. Finally, the coefficients used in the simulation

Figure 3.5: Performance of the simulations for each coefficient set, plotted as a function of
the δθ used.

that rated the highest performance are shown in table 3.1 as Optimized.
Although these coefficients seem to be very similar to the Original ones,
the velocity profiles obtained for the wall jet are now much closer to the
experimental ones, as it can be seen in fig. 3.6.

3.3 Effect of the actuator

3.3.1 Numerical setup
The flow over the half cylinder is hereafter studied without and with
actuation, as acted by the body force model. A half cylinder of dia-
meter D = 0.1 m is subject to a cross-flow having free stream velocity

21
CHAPTER 3. DRAG REDUCTION STUDY

Figure 3.6: Comparison between experimental (symbols) and numerical (lines) velocity
profiles using the optimized coefficients at several locations downstream of the actuator
(bottom). The velocity field obtained from the simulation as well as the sampling locations
are also shown (top).

22
 3.3. EFFECT OF THE ACTUATOR

Table 3.1: Various values used for the model coefficients.

Fmax (N/m3 ) ξ (m-1 ) ζ (m-1 ) δθ (◦ )

Original 70000 3250 2ξ = 3250 0
Modified 35000 3250 2ξ = 3250 0
Optimized 21900 2050 3375 -5.7
Effective 70000 2050 3375 -5.7

U∞ = 10 m/s directed in the x-direction; the origin of the Cartesian co-

ordinate system is taken at the upstream edge of the half cylinder, as
depicted in fig. 3.7. Thus the Reynolds number for the case is
U∞ D
Re = ≈ 63 × 103 , (3.4)
ν
making the flow fully turbulent. LES with the WALE sub-grid scale
model using the bounded central difference scheme were used to simu-
late the flow. The dimensions of the domain have been chosen to match
the size of the the test section in a parallel wind tunnel experiment. The
top and bottom boundaries are thus defined as no-slip walls, and the
downstream boundary is set as a pressure outlet. The side boundaries,
however, were set as periodic boundaries. This was done to simulate
an infinitely long cylinder. In the experiments the cylinder was taken
long enough such that side effects would be negligible for measurements

Figure 3.7: Geometry and size of the computational domain.

23
CHAPTER 3. DRAG REDUCTION STUDY

made in the middle of the test section. Finally, the inlet was set with a
uniform velocity of U∞ , similarly to what has been done in ref. [39], and
turbulence was left to develop within the domain. This was motivated
by the ability of the wind tunnel used in the experiments to deliver a
very homogeneous flow with low turbulent intensity at the beginning of
the test section. Numerically, this is also very straightforward since an
independent turbulent flow does not have to be computed to serve as a
time dependent inlet boundary condition. Moreover, the final interest of
this study concerns the drag of a bluff body moving through quiescent
air where free stream turbulence levels are expected to be low.
The mesh contained 3 × 106 cells and was made of cube cells in
the volume that were trimmed on the boundaries to fit the geometry.
The bottom boundary was covered with seven layers of prismatic cells
growing in the wall normal direction according to an hyperbolic tangent
law. Several levels of refinements were used around the cylinder and in
the wake in order to provide a good resolution of the turbulent eddies. An
additional level of refinement was used around and slightly downstream
of the actuation location, so that the body force would be well resolved in
space as well. This refinement was present even for the un-actuated case,
so that the mesh is strictly the same for both simulations. The mesh has
been refined by successively doubling the number of cells to check that
no significant changes in terms of drag coefficient value and behavior
could be identified.

3.3.2 Turbulence resolution of the un-actuated flow case

A simulation was first run with no body force applied to the flow. The
results were used to check that the flow was properly resolved, and
served as a baseline to gauge the effect of activating the body force.
In order to estimate the accuracy of the LES simulation, the ratio of
the time-averaged subgrid scale viscosity with the fluid viscosity, νt /ν,
was checked. Furthermore, the time average of the resolved turbulent
kinetic energy was estimated as

1 3 2
2∑
Kres = ui − ui 2 . (3.5)
i=1

This quantity is compared to the time averaged subgrid scale turbulent

kinetic energy, defined as

Ksgs = Ct νt |S| , (3.6)

24
 3.3. EFFECT OF THE ACTUATOR

where Ct = 3.5 is a time scale coefficient of the WALE model.

(a)

(b)

Figure 3.8: (a) Mean viscosity ratio and (b) comparison of the mean resolved and modeled
turbulent kinetic energies in a plane of constant z-coordinate .

All the aforementioned quantities are shown in fig. 3.8. Recalling

from (2.6) that νt is proportional to the local cell characteristic length,
sudden jumps are to be expected from a given level of mesh refinement
to the next one. Indeed, while flow properties vary continuously in space,
such jumps are visible for νt /ν and Ksgs . Nevertheless, the turbulent
viscosity ratio never exceeds 25. A rule of thumb for the turbulent kinetic

25
CHAPTER 3. DRAG REDUCTION STUDY

energy [40] is that

Kres
> 0.8, (3.7)
Kres + Ksgs

which is true in the present simulation, suggesting that a sufficient num-

ber of eddies are being resolved.

3.3.3 Drag reduction of the actuated flow case

Figure 3.9: Time- and spanwise-averaged pressure coefficient on the surface of the half
cylinder.

Figure 3.9 shows the time- and spanwise-average of the pressure

coefficient,
p
Cp = 1
, (3.8)
2
2 ρU∞

for the non-actuated case in solid red line. It can be seen from the non-
actuated case that the pressure coefficient gradient reverses at a position
x = 0.3D, which corresponds to an angle of θ = 78.5◦ . The adverse
pressure gradient after that point will eventually cause the flow to detach.
Thus the position of the actuation in the actuated case was set to θ = 78.5◦ ,
with the aim to delay the reversal of the pressure gradient. This is in

26
 3.3. EFFECT OF THE ACTUATOR

turn expected to delay separation and finally reduce the drag of the half
cylinder. It should be noted that since the electrodes are not modeled in
the simulations, the offset δθ between the edge of the exposed electrode
and the beginning of the body force is omitted for simplicity, and the
location identifies hereafter the beginning of the body force rather than
the edge of the exposed electrode.
A simulation was thus run with actuation activated, with the body
force model using the Optimized set coefficients indicated in table 3.1.
Unfortunately, this resulted in the drag of the half cylinder being virtually
unchanged compared to the baseline without actuation. Notwithstanding
the lack of effect of the body force, the hypothesis of a bad modeling
of the actuator can be reasonably discarded since no significant drag
reduction could be achieved experimentally either [27]. The hypothesis
of a poor choice of location for the actuation does not hold much better,
as various locations θ were tried experimentally without more success.
Therefore it was concluded that the present actuator, used in a steady
manner, was not good enough to influence the drag in a significant
manner. But this does not imply that such a body force distribution, or
that any plasma actuator, cannot achieve drag reduction objectives. In an
attempt to determine just how strong an actuator should need to be to
satisfy such a goal, another simulation was run with an increased value
of Fmax , which was then set to the Original value of 70000 N/m-3 . Even
though this triples the Optimized Fmax value, such a value was obtained
with another actuator which still makes it a realistic value [31] (although
the total amount of force would differ slightly due to different spatial
coefficients used). The set of coefficients used for this new simulations is
labeled Effective in table 3.1, and was used to obtain the results shown as
the actuated case.
The mean pressure coefficient, averaged in the spanwise direction, is
shown in fig. 3.9. Small variations with the baseline case can be seen that,
summed over the whole surface, create a 4.6% relative decrease in the
drag coefficient. This decrease can be seen in fig. 3.10: A simulation was
run without actuation until the flow became statistically steady, then the
time history of the drag coefficient was recorded. The state of the flow at
t∗ = 10 was used as an initial condition for a simulation with actuation,
for which the drag was also recorded. Despite the fluctuations due to
turbulence, the mean drag is smaller in the actuated case. Here,

U∞ t
t∗ = (3.9)
D

27
CHAPTER 3. DRAG REDUCTION STUDY

denotes the dimensionless time based of the cylinder diameter, and one
throughflow of the entire domain corresponds to t∗ = 23.65.

Figure 3.10: Drag coefficient with and without actuation as a function of dimensionless
time.

3.4 Summary
Rather than modeling all the physical phenomena occurring within the
plasma, a body force model was used to model directly the effect of the
actuator on the ambient fluid. This model uses coefficients to describe the
strength and spatial distribution of the body force, and the resulting effect
was seen to be dependent on the particular combination of coefficients.
Thus an optimization algorithm was used to determine the value of
each coefficient such that the velocity profiles of a wall jet created by the
actuator would match experimental data. This model was subsequently
used in an attempt to reduce the drag of a half submerged cylinder by
adding actuation where the pressure gradient reverses. However, it was
found that the strength of the body force calibrated this way was not
large enough to produce any significant change, but using a higher value
could decrease the drag coefficient by about 4%.

28
 4 Tools for flow
improvement

As presented in the previous chapter, by applying a tangential body force

on top of the half-cylinder the drag coefficient could be reduced by a
significant 4%. However, the magnitude of the force necessary to do so
was about three times larger than than the one obtained by optimization
with respect to an actual actuator. This is coherent with experiments,
where a single actuator did not have a significant effect on the flow,
and two actuators “in series” have subsequently been used. Yet in the
numerical simulation, the actuator was set at a location—the pressure
gradient reversal line—that was merely an educated guess. Therefore
it is possible that significant drag reduction may be obtained with the
force amplitude found in the experimental study, e.g. if the location of
the force were shifted. One may then undertake a parametric study of
the drag reduction with respect to the angular position of the actuator.
A second parameter of this study could, of course, be the amplitude
of the body force. It has also been widely observed that the efficiency
of plasma actuators increases when they are used in an intermittent
way [41]. This gain in efficiency is explained by two factors: first, because
an intermittent wall jet may maintain the same effect on a flow as a steady
jet, or even have a greater effect due to interactions with incoming flow
structures; second, by the energy saved while the actuator is off. Thus
in order to find the best regime for the actuator, one may add two new
parameters to the parametric study: a frequency of actuation and a duty
cycle—the fraction of a period during which the actuator is on.
On the other hand, the cost of running LES simulations for a study
of four independent parameters can be expected to be quite high, and in
order to get physical insight to guide the work, one may look instead for
flow analysis tools. By using dynamic mode decomposition (DMD) it is

29
CHAPTER 4. TOOLS FOR FLOW IMPROVEMENT

possible to identify, within a turbulent flow, coherent structures. The res-

ult of a DMD typically consists of the mean flow (from which the location
of e.g. the separation line can be deduced), as well as secondary struc-
tures (such as vortex shedding or periodic instabilities). The advantage
of DMD over other decomposition techniques, like the more extens-
ively used proper orthogonal decomposition, is that each 3-dimensional
structure is associated with a single frequency. Thus, should particular
structures be found to cause a great part of the drag, a straightforward
approach would be to set the actuator at the location of this structure,
and to operate it at a related frequency in order to cause great interaction
with the natural structure.
This chapter introduces the foundation of DMD and describes the
method that was used to compute the modes from snapshot matrices.
Further, the accuracy of the modes obtained by DMD and LES is assessed
against DNS results in the case of a channel flow oscillating at a fixed,
forced frequency.

4.1 Dynamic mode decomposition

The objective of mode decomposition is to separate time and space de-
pendent flow information, which may be particularly difficult to appre-
hend for complex turbulent flows, into several structures that may be
more familiar or easier to analyze, and may be looked at individually. In
other words, if s(t) is a vector describing the state of a quantity at time t,
then one is looking to find a number of pairs {α(t), ϕ} such that

J
s(t) = ∑ α(t) ϕj . (4.1)
j=1

Several methods have been proposed to find such pairs, and perhaps
the most used is currently the proper orthogonal decomposition (POD),
which was introduced in 1967 [42] and has been used extensively ever
since. The principle of POD is to create an orthogonal base that spans
the space of s(t), choose the vectors ϕ as the base vectors, then find the
associated α(t) that minimize energy distribution across the modes. The
advantage of this method is that the energy of the flow is distributed in as
few modes as possible. Detailed derivations are provided in the literature
[43, 44]. The DMD on the other hand, has been introduced relatively
recently to fluid mechanics [45, 46], although it is based upon long lived

30
 4.1. DYNAMIC MODE DECOMPOSITION

ideas of linear algebra and dynamical systems. The modes obtained from
DMD, rather than being orthogonal in space between them, present the
particularity of being periodic in time. The attention to DMD has been
rapidly increasing due to its ability to extract coherent spatial structures
presenting a frequency-based temporal evolution which is both intuitive
and common in signal processing. As a result DMD seems particularly
well suited to analyzing flows with periodic phenomena.
The constraint of periodicity that we put on α(t) can be expressed by
writing
 t
α(t) = reiω , (4.2)

where r is a growth or decay rate that modes are also allowed to have.
Thus, with λ = reiω , and assuming that the relation 4.1 holds we then
have (in discrete time)

 

s(tn+1 ) = Φ λn+1 

 

= ΦΛ λnj 
(4.3)
 

= ΦΛΦ−1 M λn 

= ΦΛΦ−1 s(tn )

With Λ a diagonal matrix whose elements are the complex values λ. Since
ΦΛΦ−1 is the diagonal form of a matrix A, it appears that the modes
present in the columns of Φ are the eigenvectors of A, while the values
of λ are its eigenvalues. Although the matrix A is not known a priori, a
similar matrix B can be computed from the Krylov subspace {s(tn )}n∈1..J .
A numerically robust way to do so is to create a snapshots matrix S whose
columns are made of N vectors s(tn ) taken at equidistant times. Let the
notation Š designate the first N − 1 columns of S, while Ŝ designate the
last N − 1 columns. The matrix Š is then replaced by its singular value
decomposition

Š = UΣV t . (4.4)

31
CHAPTER 4. TOOLS FOR FLOW IMPROVEMENT

From 4.3 and from properties of the SVD matrices, we have

Ŝ = AŠ
Ŝ = AUΣV t (4.5)
U t ŜVΣ+ = U t AU.

Since U is unitary, U t AU is indeed a matrix similar to A, and it can be

computed using the known matrices U t , Ŝ, V and Σ+ , the pseudo-inverse
of Σ. By definition of a similar matrix, A and B = U t AU will share the
same eigenvalues. Moreover, the eigenvectors of A can be obtained from
the eigenvectors of B by multiplying by U from the left. Therefore the
modes ϕ and their temporal coefficients λ (containing their frequency ω
and growth rate r) can be found by solving the eigen-problem for B.
In principle the turbulent part of the flow would require, by nature,
an infinite number of modes to be decomposed.
The number of positions where flow signals can be sampled can
be quite high, particularly in the case of numerical simulations where
it may easily match the number of grid points. In LES or DNS, this
gives a number of rows of the order m of a million. In order to achieve
convergence of the modes, a large number n of snapshot vectors should
also be used. Depending on the time of the simulation and the time
resolution, the number of snapshots can reach several thousands. This
results in possibly very large snapshot matrices, and the steps required
to compute the DMD can become extremely costly in CPU time and
memory. In particular, the time necessary to compute the SVD alone is of
the order of mn2 .
Therefore a C++ computer code named dymode was developed to per-
form DMD computations. Dymode was developped using the Eigen C++
template library for linear algebra [47] and has the advantage of using
the parallel library ScaLAPACK [48] for its computations, which allows
to reduce both the time and the memory required per computing node.
A user manual for the code, which provides more detailed explanations
about the DMD, a description of the numerical implementation and a
compilation and usage help, is present in Paper 2.

4.2 LES of pulsating channel flow

For the continuity of this work, it is interesting to assess how well the
results from LES compare to DNS in the case of flows where small amp-

32
 4.2. LES OF PULSATING CHANNEL FLOW

litude oscillations are forced, and how accurate the modes obtained by
DMD are. The pulsating channel flow is a good flow case for this purpose
since coherent structures can be reduced to 1-dimensional profiles and
the forcing frequency is known, making it easy to use other filtering
techniques such as phase averaging in order to obtain reference data.
Although this case has been investigated experimentally over several
decades [49–53], it is still a source of many interrogations. For instance,
at certain frequencies, the amplitude of the wall shear stress becomes
even smaller than when the pulsating Reynolds stress is null. The reason
remains unclear, even after LES were conducted to shed light on the
phenomenon [54]. Thus the flow is complicated enough (and has eluded
good modeling for long enough) that there are non-trivial structures to
be found and used to assess DMD modes.

4.2.1 Numerical setup

The channel is formed by two parallel planes specified as non-permeable
walls with no slip. It is simulated as being infinite in the streamwise
and spanwise dimensions by the use of periodic boundary conditions. A
sketch of the computational domain is shown in fig. 4.1. Simulations have

y
mean flow direction

Ly = 2h

Lz

Lx
z

Figure 4.1: Geometry of the channel. The computational domain is marked with plain lines.

been performed with two codes. A DNS has been performed using the
pseudo-spectral solver SIMSON [21]. In addition, an LES was performed

33
CHAPTER 4. TOOLS FOR FLOW IMPROVEMENT

with the finite volume software STAR-CCM+ 9.06.009. The dimensions of

the computational domain in the streamwise, wall-normal and spanwise
directions, respectively L x , Ly and Lz , illustrated in fig. 4.1 are reported
for both simulations in table 4.1 together with the meshes characteristic
lengths. The meshes were made of regular hexahedral cells with constant
sizes in the streamwise and spanwise directions. In the wall-normal
direction, cells were stretched using a hyperbolic tangent law for the LES,
and a Gauss-Lobatto law for the DNS.
The flow is driven by a time dependent pressure gradient

dp dp
(t) = [1 + β sin(ωt)] (4.6)
dx dx
in the form of a homogeneous volume force in the streamwise direction
dp
for the DNS; and by applying a pressure jump L x dx (t) between the
streamwise boundaries for the LES. In order to obtain a mean Reynolds
number Reτ = uτ h/ν = 350 based on the mean friction velocity uτ , the
mean pressure gradient is chosen to be

dp Re2 µ2
= − 3τ . (4.7)
dx h ρ

The corresponding Reynolds number based on the mean velocity at the

center of the channel is approximately Rec = uc h/ν ≈ 7000.
The frequency and amplitude of the flow oscillations are determined
by the values of β and ω. The non-dimensional angular frequency
ων
ω+ = = 0.01 (4.8)
u2τ

was used since it lies at the heart of a transition region between two well
modeled types of flow. When ω +  0.01 or ω +  0.01, respectively the
quasi-static or quasi-laminar models provide good accuracy. Further-
more, the oscillating wall shear stress has the peculiarity of dropping
around this frequency, and, while useful for engineering purposes such
as acoustic damping, the reasons for this phenomenon are still unclear.
The data was post-processed using phase averaging at the forcing fre-
quency for both simulations, and DMD was also used to post-process the
LES data. Since several DMD modes may have the same frequency (but
different growth or decay rates), the modes were selected based on their
norm. Indeed, this norm is proportional to the energy contribution of the

34
 4.2. LES OF PULSATING CHANNEL FLOW

Table 4.1: Discretization of the domains.

DNS LES
Length Cells ∆+ Length Cells ∆+
x 4πh 384 12 10h 80 43.75
y 2h 193 0.19 - 6 2h 64 0.37 - 32.7
z 2πh 384 6 3h 64 16.4

mode to the flow. For modes representing coherent structures, this norm
should have a fixed value. However the modes resulting from the inco-
herent fluctuations should see their energy decrease when more modes
are added. Therefore if the DMD is computed from a sufficiently large
dataset, the physical modes can be identified by looking for relatively
higher mode norms. Particularly, with the data obtained by LES, several
modes systematically appeared at zero frequency. In that case the mode
plotted and compared with mean profiles were the mode with the largest
norm. The same criteria was applied to select the modes compared with
the oscillating profiles when several modes appeared with frequencies
close to the forcing frequency.
By computing the phase average h f i of a fluctuating quantity f , its
time average f = h f i and coherent part fe = h f i − f can be obtained. The
remainder f 0 is assumed to be purely turbulent fluctuations, which leads
to the following triple decomposition of the fluctuating quantities [55]

f (x, t) = f (x) + fe(x, t) + f 0 (x, t). (4.9)

While data is available for f by sampling the simulation result, and can
be post-processed directly using either phase averaging or DMD, it is not
f 0 f 0 . First f 0 was
as straightforward to obtain quantities such as f 0 f 0 or g
obtained as

f 0 (x, t) = f (x, t) − f (x) − fe(x, t), (4.10)

where f (x) and fe(x, t) have been obtained by phase averaging. Then f 0
(or derived quantities) was processed using either phase averaging or
DMD. This was done to ensure that both the phase averaging and the
DMD had the same input so that the results could be compared fairly.

35
CHAPTER 4. TOOLS FOR FLOW IMPROVEMENT

4.2.2 Mean profiles

Streamwise velocity
First, the mean streamwise velocity profiles obtainable by LES are com-
pared with the ones from DNS. Two velocity profiles corresponding to
the LES data are plotted in fig. 4.2: one is obtained from the traditional
time averaging of the snapshots; the other corresponds to the DMD mode
with zero frequency and highest norm. Both are virtually identical, and
come very close to the DNS profile.

22
DNS
20 LES
LES - DM D
18 DNS (stea dy)
Scotti (DNS, 7 0% )
16 y+ ; ln (y+ )/0 .39+ 5

14

12
u+

10

0
−2 −1 0 1 2 3
10 10 10 10 10 10
y+

Figure 4.2: Mean velocity u vs. wall distance; normalized in wall units.

It has been assumed in the past that adding pulsations in a turbulent

channel flow would not modify mean profiles as long as the amplitude of
the pulsations remains below 100% of the mean velocity [53]. Since in the
present case the amplitude of the oscillations at the center of the channel
is only 15%, the mean velocity profiles should be identical to the profiles
of a steady channel flow at Reτ = 350. This was verified by running a
DNS simulation where β = 0. All the profiles come within a few percents
of each other, and even fit a logarithmic profile characteristic of steady
channel flows [40].
A profile of the mean velocity in a case where the amplitude of the
oscillations was 70% and obtained by Scotti et al. using LES [56] is also
shown for comparison. This velocity seems to be slightly overestimated,
which suggests that the threshold given in ref. [53] should be decreased
lest mean profiles deviate from the steady case.

36
 4.2. LES OF PULSATING CHANNEL FLOW

Reynolds stresses

10 DNS
LES
LES - DM D 1
+

M o ser 3 9 5
u′ u′

5 Sco tti DNS Total stress

DNS −u ′ ′+
0.9 1 v 2+
0 DNS µ du dy
0.8
LES Resolved stress
LES −u ′+
1 v 2+
′
1 0.7
LES µ du dy
0.6 LES Scotti
+

+
v′ v′

0.5 LES - DMD −u′ v′

τ12
+
0.5

0 0.4

2 0.3

1.5
0.2
+
w′ w′

1
0.1
0.5

0 0
0 50 100 150 200 250 300 350 0 50 100 150 200 250 300 350
y+ y+

Figure 4.3: Mean profiles of various stresses, normalized in wall units.

When looking at the turbulent stresses, once again DNS results from
the pulsating case match the results from a steady flow DNS. Profiles
from Moser et al. [57] are also included here as a reference for the steady
channel flow. Although the DNS by Moser et al. was performed at
the slightly higher Reynolds number of Reτ = 395, the profiles of the
mean symmetric stresses, scaled in wall units, come close the the profiles
obtained by the current DNS as illustrated in figure 4.3. Only a small
variation that puts it consistently below can be seen. This could be
a result of having oscillations within the flow, or simulation variations.
Compared to any of the DNS, the LES results underestimate the spanwise
and wall-normal Reynolds stresses. This can be expected for an LES
simulation since some of the turbulent energy is modeled rather than
resolved in the flow. However, the peak in the streamwise Reynolds
stress profile is overestimated, while the remaining of the profile is close
to DNS values. This behavior of the symmetric Reynolds stresses is
usual in LES simulations, and has been displayed with several models
and codes [58]. The LES from ref. [56] exhibits similar profiles for the
spanwise and wall-normal stresses, although no overshoot can be seen
in the profile of the streamwise stress. Since that LES and the present one
use similar spatial resolutions, it is likely that the variation in the profiles

37
CHAPTER 4. TOOLS FOR FLOW IMPROVEMENT

come from the larger amplitude of oscillations used.

A balance of the mean shear stress is also provided in fig. 4.3. It can
be seen that that the total stress obtained by DNS follows a linear profile
in the wall normal direction, which is what would be expected from a
steady channel flow. LES profiles of the turbulent stress u0 v0 obtained
by averaging and by DMD present extremely little differences. Both are
slightly lower lower than the DNS profile, as a result of some turbulence
being modeled. This difference can be seen in the total resolved stress
profile which deviates from a linear profile close to the wall.
The mean profiles obtained by DNS for this pulsating channel flow are
virtually identical with profiles in a steady channel flow. This agrees with
the low amplitude of oscillation being used and the statement in ref. [53].
The LES profiles are also similar to profiles that would be expected from
a steady flow, with the mean velocity relatively well captured and a
deficit in resolved turbulent stresses due to sub-grid scale modeling, with
the exception of u0 u0 whose peak is overestimated. The method used to
obtain mean profiles, phase averaging or DMD, has no influence on the
results.

4.2.3 Oscillating profiles

The coherent part of the fluctuations need not have a sinusoidal wave-
form [53], and in fact they were found to carry several frequency com-
ponents. This is illustrated in fig. 4.4 where the oscillating part of u0 v0
obtained by phase averaging is shown, normalized in wall units. The
same remark can be made for the oscillating part of u, although frequen-
cies other than the forcing frequency are weaker than for ug 0 v0 . By com-
puting the Fourier transform of oscillating quantities at all y-locations, it
is possible to plot profiles of amplitude and phase for every frequency
component.
On the other hand, when computing the DMD of flow snapshots, a
mode should be obtained for each coherent frequency. The modes them-
selves contain the amplitude and phase profiles, such that no additional
step is required after computing the DMD. However, for all the quantities
assessed, the only DMD modes standing out from the background noise
were found at zero and at the forcing frequency. In other words, none
of the harmonic components present in the flow and obtained by phase
averaging could be identified by DMD.
An explanation can be found by looking at the Fourier transform of
the flow history. As an example, the Fourier transform of u(y+ = 8.15, t) is

38
 4.2. LES OF PULSATING CHANNEL FLOW

0.15

0.1

0.05

0
+
g
u′ v′

−0.05

−0.1

y+ = 4 .7
−0.15 y+ = 7 .9
y+ = 1 6.8
y+ = 7 2.3
−0.2
0 0.2 0.4 0.6 0.8 1
φ

Figure 4.4: Oscillating part of u0 v0 at various y+ -locations, vs. pulsation phase.

2
10

1
10

0
10

−1
10
Am plitude

−2
10

−3
10

−4
10

−5
10

−6
10
−0.04 −0.02 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16
ω+

Figure 4.5: Fourier transform of u(t) sampled at y+ = 8.15.

plotted in fig. 4.5. The location y+ = 8.15 was chosen since the harmonic
component is particularly strong there. It can be seen that while the mean
and forcing frequency components stand out very clearly from the noise
level, the component at twice the forcing frequency only stands out from
the noise in neighbor frequencies, but has comparable amplitude with
the low frequency noise. Furthermore, the amplitude of this component

39
CHAPTER 4. TOOLS FOR FLOW IMPROVEMENT

is location dependent, and is null over a large part of the channel. On the
other hand, DMD does not produce modes at fixed and equally spaced
frequencies, but determines the frequencies based on frequency contri-
butions over the whole domain. Therefore, the component at twice the
forcing frequency was too weak to produce a mode at that particular fre-
quency, and more modes were produced instead at low frequencies since
this is where most of the energy is found. This is illustrated in fig. 4.6,
which shows the norm of the modes as a function of there frequency.

0
10
M o de no rm

0 0.02 0.04 0.06 0.08

ω+

Figure 4.6: Spectrum of the DMD modes for u.

Even though the harmonic components were not found by the DMD,
profiles of the coherent fluctuations at the forcing frequency can still be
compared with those obtained by phase averaging and Fourier transform.
As it can be seen in fig. 4.6, the mode having the second largest norm
is found at an angular frequency of 0.00986, only 1.4% away from the
forcing frequency. The modulus and argument of this mode is shown
as a function of y+ in fig. 4.7, where it is compared with the profiles
obtained by phase averaging the both the LES and DNS data. While the
profiles obtained from LES data differ from the DNS profile, both DMD
and phase averaging yield the same profiles of amplitude and phase.
0 v0 are also shown in fig. 4.7. Again, a single mode from
Profiles for ug
the DMD of u0 v0 stood out around the forcing frequency, at 0.00979. The
amplitude and particularly the phase of this mode, however, differ from
the phase averaged profiles.

40
 4.2. LES OF PULSATING CHANNEL FLOW

1.1

0.2
DNS
1.05 LES
0.15 Quasi−static
u|∗ω +

Non−equilibrium
|e

-u′ v′ - +
-+
ω
-

1 0.1
-g
-

0.05
0.95

0.8
DNS
LES 400
0.6 LES - DM D
Lam inar m o del
u)

Q uasi-Sta tic m o del 300

a rg ω + (e

0.4
Non-equilibrium m o del
2
g
′ v′

Scotti (LES 70% am plitude) 200

a rgω + u

0.2
1

100
0

0
−0.2
0 0.2 0.4 0.6 0.8 1
y∗ −100
0 50 100 150 200 250 300 350

Figure 4.7: Profiles in amplitude and relative phase for u (left) and u0 v0 (right) of the forcing
frequency components.

41
 5 Conclusions and
outlook

Passive flow control is reaching its limits, and is either as optimized as it

will ever be, or it is greatly competing with other areas of vehicle design,
such as under-hood integration, structural integrity or legal regulations.
This has lead to the dawn of active control as means to push even further
flow control in many areas of aeronautics and vehicle design. This rise is
supported by advances in other engineering domains such as material
sciences, electronics, MEMS, control theory, metrology and so forth.
In particular, plasma actuators are receiving a lot of attention from
researchers and the industry since they combine mechanical simplicity,
high reaction times and low power consumption. Although the mech-
anism by which they operate is extremely costly to simulate and not
entirely understood, the plasma physics is often replaced by the effect
of this physics in simulations of practical cases. In the present work, a
model of this effect using a body force approach was adapted from the
literature, and an empirical optimization method has been employed
to calibrate the model against experimental data. The result of using
this model matched very well with measurements at several locations
downstream of the actuator in a case where the actuator was the only
source of momentum.
This model was furthermore utilized in a highly resolved LES of the
flow over a half-submerged cylinder, chosen for its similarity with the
A-pillar of a truck. The actuated turbulent flow showed no significant
change compared to the baseline case. However, a similar conclusion
was drawn experimentally, so that there is no reason to dismiss the model
used. In fact, arbitrarily increasing the strength of the model yielded a
drag reduction of more than 4%. The amount of force used to achieve
this result was the same as reported in a literature paper, so that it is

43
CHAPTER 5. CONCLUSIONS AND OUTLOOK

reasonable to assume that actuators can be built to offer the same strength.
Although it is possible that the actuator used in the experiments
and simulated in the LES was simply not strong enough to achieve
drag reduction, other possibilities include that the location of actuation
was not chosen appropriately, or that a steady use of the plasma is not
appropriate for the incoming flow. Rather than trying all possibilities, a
parallel code to compute the DMD of large datasets has been written. The
ambition is that the DMD will allow to identify locations and frequencies
within the natural flow where actuation will have the greatest impact.
The DMD code was evaluated in a pulsating channel flow case, for
which an LES of the case necessitating reasonable computing power was
run. It was then post-processed in two ways: with a usual phase aver-
aging, and with the DMD code. It was shown that the DMD is able to
identify oscillating structures in the flow having high energy, although
quieter harmonics do not appear as easily and remain buried in noise.
The LES data itself compares rather well compared to DNS when it comes
to mean and oscillating profiles of first order quantities. Stresses how-
ever are somewhat less accurate, but they may not always be necessary,
depending on the applications.

The model used in the simulations could only be compared with ex-
perimental data in the no-flow case. Indeed, in the case with flow, no
effect was measured, either experimentally or numerically. While this is
not a reason to dismiss the model, it is not a reason to validate it either.
Therefore it is desirable to compare the model in a flow case where an
effect can be measured. Experimental data is now available where two
successive actuators are used and a change in drag was recorded for
several positions of the pair. A good step will have been made if these
measurements can be reproduced in a future LES using the model at
the two same locations. Furthermore, DMD should be applied to the
natural flow in order to extract modes. A link between these modes and
actuation efficiency then needs to be made.

44
6 Summary of
appended papers

Paper 1
Effect of a SDBD on the Drag of a Half-Submerged Cylinder in Crossflow.
R. Futrzynski, G. Efraimsson. In ASME 2014 4th Joint US-European Flu-
ids Engineering Division Summer Meeting collocated with the ASME
2014 12th International Conference on Nanochannels, Microchannels,
and Minichannels (peer-reviewed).

In this paper the effect of a SDBD-type plasma actuator on the flow

over a half-submerged cylinder is investigated numerically. The actuator
is modeled via a body force, which is steady in time and where an expo-
nential decay in space is assumed. First, the parameters in the numerical
actuator model are determined for the case of no flow by optimization
relative to experimental data. Thereafter, numerical solutions for the
case with flow are studied numerically with and without actuation. A
grid study is performed to check that the flow structures are resolved in
both space and time. The effect of the actuator is examined. Although
no significant change is observed when using the optimized parameters,
using a stronger body force yields a reduction in drag of the order of 5%.

Paper 2
Dymode, A parallel dynamic mode decomposition software.
R. Futrzynski, G. Efraimsson. Internal report, KTH 2015. TRITA-AVE
2014:78, ISSN 1651-7660, ISBN: 978-91-7595-386-1

Computing the dynamic mode decomposition may necessitate a lot of

45
CHAPTER 6. SUMMARY OF APPENDED PAPERS

computing resources when it is applied to large datasets. Dymode is C++

program that uses parallel computing in order to accelerate the extraction
of DMD modes and lighten the memory requirement per CPU. Several
parameters can be specified in order to control the computational aspects
of the program as well as the input and output of the decomposition,
particularly how the modes are sorted. The dymode package also in-
cludes dymodem, a Matlab implementation of the code which accepts
the same arguments as dymode, when they are relevant, and produces
the same output. It can be useful to use dymodem when dealing with
smaller datasets, or to validate the output from dymode. The principle
of DMD and numerical algorithms are discussed, and a description of
the implementation is provided.

Paper 3
Numerical study of the Stokes layer in oscillating channel flow.
R. Futrzynski, C. Weng, S. Boij, A. Hanifi, G. Efraimsson. To be submitted.

Here we propose to compare methodologies relevant for the numer-

ical study of pulsating ows. DNS and LES of the same case are compared,
and the data is processed by using Reynold’s triple decomposition. In
addition the LES data is also processed with dynamic mode decomposi-
tion, a relatively new post-processing method suitable to the extraction
of frequency-based coherent structures. The case studied was chosen
in continuity of past simulations by other authors, and is a pulsating
channel flow. The frequency of the pulsation is chosen to be in the trans-
ition range between the quasi-laminar and quasi-static states since it
proves particularly challenging for models; and the amplitude of the
oscillations is chosen relatively low in order to ease the comparison with
linear models. The quality of each method is evaluated by comparing
several quantities of interest for such a flow case, and, when applicable,
by highlighting similarities and differences with other flow cases.

46
Bibliography

[1] R. H. Barnard, Road Vehicle Aerodynamic Design. MechAero Publish-

ing, third ed., 2009.
[2] J1252, “SAE Wind Tunnel Test Procedure for Trucks and Buses,”
1979.
[3] S. Kang and H. Choi, “Suboptimal feedback control of turbulent flow
over a backward-facing step,” Journal of Fluid Mechanics, vol. 463,
pp. 201–227, July 2002.
[4] S. Dong, G. Triantafyllou, and G. Karniadakis, “Elimination of Vor-
tex Streets in Bluff-Body Flows,” Physical Review Letters, vol. 100,
p. 204501, May 2008.
[5] E. Wassen, S. Eichinger, and F. Thiele, “Simulation of Active Drag
Reduction for a Square-Back Vehicle,” Active Flow Control II, pp. 241–
255, 2010.
[6] M. El-alti, V. Chernoray, M. Jahanmiri, and L. Davidson, “Experi-
mental and computational studies of active flow control on a model
truck-trailer,” in EPJ Web of Conferences, vol. 25, 2012.
[7] R. Duvigneau and M. Visonneau, “Optimization of a synthetic jet
actuator for aerodynamic stall control,” Computers & Fluids, vol. 35,
pp. 624–638, July 2006.
[8] S. H. Kim and C. Kim, “Separation control on NACA23012 using
synthetic jet,” Aerospace Science and Technology, vol. 13, pp. 172–182,
June 2009.
[9] S. Aubrun, J. McNally, F. Alvi, and A. Kourta, “Separation flow
control on a generic ground vehicle using steady microjet arrays,”
Experiments in Fluids, vol. 51, pp. 1177–1187, June 2011.

47
BIBLIOGRAPHY

[10] R. E. King, Active Flow Control Papers contributed to the Conference

"Active Flow Control 2006", Berlin, Germany, September 27 to 29, 2006.
2006.
[11] T. R. Bewley, “Flow control: new challenges for a new Renaissance,”
Progress in Aerospace Sciences, vol. 37, pp. 21–58, Jan. 2001.
[12] S. S. Collis, R. Joslin, and A. Seifert, “Issues in active flow control:
theory, control, simulation, and experiment,” Progress in Aerospace,
vol. 972, no. 3, 2004.
[13] I. Mezić, “Spectral Properties of Dynamical Systems, Model Re-
duction and Decompositions,” Nonlinear Dynamics, vol. 41, no. 1-3,
pp. 309–325, 2005.
[14] M. Morzynski, W. Stankiewicz, B. R. Noack, R. King, F. Thiele, and
G. Tadmor, “Continuous Mode Interpolation for Control-Oriented
Models of Fluid Flow,” in Active Flow Control (R. E. King, ed.),
pp. 260–278, Springerlink, 2007.
[15] C. L. Rumsey, “Successes and challenges for flow control simula-
tions,” tech. rep., American Institute of Aeronautics and Astraunics,
2008.
[16] A. Berns, U. Buder, E. Obermeier, A. Wolter, and A. Leder,
“AeroMEMS sensor array for high-resolution wall pressure measure-
ments,” Sensors and Actuators A: Physical, vol. 132, pp. 104–111, Nov.
2006.
[17] O. Reynolds, “On the Dynamical Theory of Incompressible Viscous
Fluids and the Determination of the Criterion,” Philosophical Transac-
tions of the Royal Society of London. A, vol. 186, no. 1895, pp. 123–164,
1895.
[18] J. Smagorinsky, “General circulation experiments with the primitive
equations. Part I: The basic experiment,” Monthly Weather Review,
vol. 91, no. 3, pp. 99–164, 1963.
[19] E. R. Van Driest, “On Turbulent Flow Near a Wall,” Journal of the
Aeronautical Sciences (Institute of the Aeronautical Sciences), vol. 23,
pp. 1007–1011, 1956.
[20] F. Nicoud and F. Ducros, “Subgrid-Scale Stress Modelling Based on
the Square of the Velocity Gradient Tensor,” Flow, Turbulence and
Combustion, pp. 183–200, 1999.

48
 BIBLIOGRAPHY

[21] M. Chevalier, P. Schlatter, A. Lundbladh, and D. S. Henningson,

“SIMSON A Pseudo-Spectral Solver for IncoMpreSsible BOuNdary
Layer Flows,” tech. rep., KTH Mechanics, Stockholm, 2007.
[22] C. Mockett, A comprehensive study of detached-eddy simulation. PhD
thesis, Technischen Universität Berlin, 2009.
[23] T. C. Corke, C. L. Enloe, and S. P. Wilkinson, “Dielectric Barrier
Discharge Plasma Actuators for Flow Control,” Annual Review of
Fluid Mechanics, vol. 42, pp. 505–529, Jan. 2010.
[24] W. Shyy, B. Jayaraman, and a. Andersson, “Modeling of glow
discharge-induced fluid dynamics,” Journal of Applied Physics, vol. 92,
no. 2002, pp. 6434–6443, 2002.
[25] C. L. Enloe, T. E. McLaughlin, R. D. Van Dyken, K. D. Kachner, E. J.
Jumper, and T. C. Corke, “Mechanisms and Responses of a Single
Dielectric Barrier Plasma Actuator: Plasma Morphology,” 2004.
[26] J. Vernet, R. Örlü, G. Efraimsson, and P. H. Alfredsson, “Experi-
mental Study of the Electric Wind Induced by a Dielectric Barrier
Discharge Plasma,” in 4th International Conference on Jets, Wakes and
Separated Flows, ICJWSF2013 September 17-21, 2013, Nagoya, Japan,
pp. 1–6, 2013.
[27] J. A. Vernet, R. Örlü, and P. H. Alfredsson, “Flow separation delay on
trucks A-pillars by means of Dielectric Barrier Discharge actuation,”
submitted, 2014.
[28] J. R. Roth, “Aerodynamic flow acceleration using paraelectric and
peristaltic electrohydrodynamic effects of a one atmosphere uniform
glow discharge plasma,” in Physics of Plasmas, vol. 10, pp. 2117–2126,
2003.
[29] M. L. Post and T. C. Corke, “Separation Control on HIgh Angle
of Attack Airfoil Using Plasma Actuators,” AIAA Journal, vol. 42,
pp. 2177–2184, Nov. 2004.
[30] Y.-C. Cho and W. Shyy, “Adaptive flow control of low-Reynolds
number aerodynamics using dielectric barrier discharge actuator,”
Progress in Aerospace Sciences, vol. 47, pp. 495–521, Oct. 2011.
[31] D. Greenblatt, T. Schneider, and C. Y. SchuÌĹle, “Mechanism of flow
separation control using plasma actuation,” Physics of Fluids, vol. 24,
no. 7, p. 077102, 2012.

49
BIBLIOGRAPHY

[32] T. N. Jukes and K.-S. Choi, “Control of unsteady flow separation

over a circular cylinder using dielectric-barrier-discharge surface
plasma,” Physics of Fluids, vol. 21, no. 9, p. 094106, 2009.

[33] R. Dadfar, A. Hanifi, and D. S. Henningson, “Feedback Control for

Laminarization of flow over Wings,” Flow, Turbulence and Combustion,
vol. 94, no. 4, pp. 43–62, 2014.

[34] K. P. Singh and S. Roy, “Simulation of an asymmetric single dielectric

barrier plasma actuator,” Journal of Applied Physics, vol. 98, no. 2005,
2005.

[35] B. Jayaraman and W. Shyy, “Modeling of dielectric barrier discharge-

induced fluid dynamics and heat transfer,” Progress in Aerospace
Sciences, vol. 44, pp. 139–191, Apr. 2008.

[36] D. M. Orlov, T. C. Corke, and M. P. Patel, “Electric circuit model

for aerodynamic plasma actuator,” in 44th AIAA Aerospace Sciences
Meeting and Exhibit 9 - 12 January 2006, Reno, Nevada, 2006.

[37] D. P. Rizzetta, M. R. Visbal, and P. E. Morgan, “A high-order compact

finite-difference scheme for large-eddy simulation of active flow
control,” Progress in Aerospace Sciences, vol. 44, pp. 397–426, Aug.
2008.

[38] Red Cedar Technology, “SHERPA - An Efficient and Robust Optim-

ization/Search Algorithm.”.

[39] C. Mockett, R. Perrin, T. Reimann, M. Braza, and F. Thiele, “Analysis

of detached-eddy simulation for the flow around a circular cylinder
with reference to PIV data,” Flow turbulence and combustion, vol. 85,
no. 2, pp. 167–180, 2010.

[40] S. B. Pope, Turbulent flows. Cambridge University Press, 2000.

[41] B. Göksel, D. Greenblatt, I. Rechenberg, Y. Kastantin, C. N. Nayeri,

and C. O. Paschereit, “Pulsed Plasma Actuators for Active Flow
Control at MAV Reynolds Numbers,” in Active Flow Control (R. E.
King, ed.), pp. 42–55, Springer, 2007.

[42] J. L. Lumley, “The structure of inhomogeneous turbulent flows,”

Atmospheric turbulence and radio wave propagation, pp. 166–178, 1967.

50
 BIBLIOGRAPHY

[43] W. Cazemier, R. W. C. P. Verstappen, and a. E. P. Veldman, “Proper

orthogonal decomposition and low-dimensional models for driven
cavity flows,” Physics of Fluids, vol. 10, no. 7, p. 1685, 1998.

[44] M. Manhart and H. Wengle, “A spatiotemporal decomposition of a

fully inhomogeneous turbulent flow field,” Theoretical and Computa-
tional Fluid Dynamics, vol. 5, pp. 223–242, Nov. 1993.

[45] C. W. Rowley, I. Mezić, S. Bagheri, P. Schlatter, and D. S. Henningson,

“Spectral analysis of nonlinear flows,” Journal of Fluid Mechanics,
vol. 641, p. 115, Nov. 2009.

[46] P. J. Schmid, L. Li, M. P. Juniper, and O. Pust, “Applications of the

dynamic mode decomposition,” Theoretical and Computational Fluid
Dynamics, vol. 25, pp. 249–259, Aug. 2010.

[47] G. Guennebaud, B. Jacob, and Others, “Eigen v3,” 2010.

[48] L. S. Blackford, J. Choi, A. Cleary, E. D’Azevedo, J. Demmel, I. Dhil-

lon, J. Dongarra, S. Hammarling, G. Henry, A. Petitet, K. Stanley,
D. Walker, and R. C. Whaley, ScaLAPACK User’s Guide. Philadelphia,
PA: Society for Industrial and Applied Mathematics, 1997.

[49] A. K. M. F. Hussain and W. C. Reynolds, “The mechanics of an

organized wave in turbulent shear flow. Part 2. Experimental results,”
Journal of Fluid Mechanics, vol. 54, pp. 241–261, Mar. 1972.

[50] D. Ronneberger and C. D. Ahrens, “Wall shear stress caused by

small amplitude perturbations of turbulent boundary-layer flow: an
experimental investigation,” Journal of Fluid Mechanics, vol. 83, no. 3,
pp. 433–464, 1977.

[51] M. S. Howe, “On the absorption of sound by turbulence and other

hydrodynamic flows,” IMA Journal of Applied Mathematics (Institute
of Mathematics and Its Applications), vol. 32, pp. 187–209, 1984.

[52] Z.-X. Mao and T. J. Hanratty, “Studies of the wall shear stress in a
turbulent pulsating pipe flow,” Journal of Fluid Mechanics, vol. 170,
pp. 545–564, 1986.

[53] S. F. Tardu, G. Binder, and R. F. Blackwelder, “Turbulent channel

flow with large-amplitude velocity oscillations,” Journal of Fluid
Mechanics, vol. 267, pp. 109–151, Apr. 1994.

51
BIBLIOGRAPHY

[54] P. Comte, M. Haberkorn, G. Bouchet, V. Pagneux, and Y. Aurégan,

“Large-Eddy Simulation of acoustic propagation in a turbulent
channel flow,” in Direct and Large-Eddy Simulation VI (E. Lambal-
lais, R. Friedrich, B. J. Geurts, and O. Métais, eds.), pp. 521–528,
Dordrecht: Springer Netherlands, 2006.

[55] A. K. M. F. Hussain and W. C. Reynolds, “The mechanics of an

organized wave in turbulent shear flow,” Journal of Fluid Mechanics,
vol. 41, pp. 241–258, Mar. 1970.
[56] A. Scotti and U. Piomelli, “Numerical simulation of pulsating tur-
bulent channel flow,” Physics of Fluids, vol. 13, no. 5, pp. 1367–1384,
2001.

[57] R. D. Moser, J. Kim, and N. N. Mansour, “Direct numerical simu-

lation of turbulent channel flow up to Reτ = 590,” Physics of Fluids,
vol. 11, no. 4, pp. 943–945, 1999.
[58] A. Rasam, S. Wallin, G. Brethouwer, and A. V. Johansson, “Large
eddy simulation of channel flow with and without periodic constric-
tions using the explicit algebraic subgrid-scale model,” Journal of
Turbulence, vol. 15, no. 11, pp. 752–775, 2014.

52
 Part II

APPENDED PAPERS 1-3