Mehdi Ghommem : A Low-Dimensional Tool For Predicting Force Decomposition Coefficients For Varying Inflow Conditions

Download as pdf or txt
Download as pdf or txt
You are on page 1of 14

368 Progress in Computational Fluid Dynamics, Vol. 13, No.

6, 2013

A low-dimensional tool for predicting force


decomposition coefficients for varying inflow
conditions

Mehdi Ghommem*
Division of Physical Sciences and Engineering (PSE),
King Abdullah University of Science and Technology (KAUST),
Thuwal 23955-6900, Saudi Arabia
E-mail: mehdi.ghommem@kaust.edu.sa
*Corresponding author

Imran Akhtar
Department of Mechanical Engineering,
NUST College of Electrical and Mechanical Engineering,
National University of Sciences and Technology (NUST),
Islamabad, 44000 Pakistan
E-mail: imran.akhtar@ceme.nust.edu.pk

Muhammad R. Hajj
Department of Engineering Science and Mechanics,
Virginia Tech,
Blacksburg, VA 24061, USA
E-mail: mhajj@vt.edu

Abstract: We develop a low-dimensional tool to predict the effects of unsteadiness in the inflow
on force coefficients acting on a circular cylinder using proper orthogonal decomposition (POD)
modes from steady flow simulations. The approach is based on combining POD and linear
stochastic estimator (LSE) techniques. We use POD to derive a reduced-order model (ROM) to
reconstruct the velocity field. To overcome the difficulty of developing a ROM using Poisson’s
equation, we relate the pressure field to the velocity field through a mapping function based on
LSE. The use of this approach to derive force decomposition coefficients (FDCs) under unsteady
mean flow from basis functions of the steady flow is illustrated. For both steady and unsteady
cases, the final outcome is a representation of the lift and drag coefficients in terms of velocity
and pressure temporal coefficients. Such a representation could serve as the basis for
implementing control strategies or conducting uncertainty quantification.

Keywords: reduced-order modelling; proper orthogonal decomposition; POD; linear stochastic


estimator; LSE; force decomposition coefficient; FDC; unsteady inflow.

Reference to this paper should be made as follows: Ghommem, M., Akhtar, I. and Hajj, M.R.
(2013) ‘A low-dimensional tool for predicting force decomposition coefficients for varying
inflow conditions’, Progress in Computational Fluid Dynamics, Vol. 13, No. 6, pp.368–381.

Biographical notes: Mehdi Ghommem is currently a Postdoctoral Fellow at the Division of


Physical Sciences and Engineering, King Abdullah University of Science and Technology
(KAUST). He obtained his BEng degree (with distinction) and MS in Computational Mechanics
from Tunisia Polytechnic School in 2007 and 2008, respectively. He completed his PhD in
Engineering Mechanics from Virginia Tech in 2011. His research interests include reduced-order
modeling, analysis and identification of non-linear systems and fluid-structure interactions.

Imran Akhtar is an Assistant Professor at the Department of Mechanical Engineering, the College
of E&ME, National University of Sciences and Technology (NUST), Islamabad, Pakistan since
2011. He received his Bachelors degree in Aerospace Engineering in 1996. He received his MS
in Aerospace from The George Washington University, Washington DC in 2003. He then moved
to Virginia Tech and received his PhD in Engineering Mechanics in 2008. He later joined
Interdisciplinary Center for Applied Mathematics at Virginia Tech as a Postdoctoral Associate.
His research interests include high performance computing, reduced-order modelling, non-linear
systems, flow control and energy efficient buildings.

Copyright © 2013 Inderscience Enterprises Ltd.


A low-dimensional tool for predicting force decomposition coefficients for varying inflow conditions 369

Muhammad R. Hajj is a Professor of Engineering Science and Mechanics (ESM) at Virginia


Tech. He received his BEng degree (with distinction) from the American University of Beirut in
1983. He received his MS and PhD degrees from the University of Texas at Austin in 1985 and
1990. He joined the ESM Department in 1992 after spending two years as a Research Associate
and Lecturer at the University of Texas at Austin. His research interests include reduced-order
modelling, analysis, and identification of non-linear systems with applications in fluid dynamics,
fluid-structure interactions, and structures vibrations. He is currently the Director of Research
and Graduate Studies of the Department of Engineering Science and Mechanics at Virginia
Tech.

1 Introduction (2009b) for more details]. Flow quantities such as force


coefficients are then obtained by integrating the pressure
The governing equations for fluid flows consist of coefficients. As pointed by Akhtar et al. (2009b), this model
partial-differential equations (PDEs) which involve an works well only for low Reynolds numbers. Furthermore, it
infinite number of degrees of freedom. Recent advances in cannot be implemented to extend the use of the ROM from
computer hardware and software have extended the use of one flow condition to another. An example would be the
direct numerical simulations whereby, the aforementioned need to derive force coefficients under unsteady mean flow
equations are discretised and integrated with numerical conditions (e.g., gust) using parameters of a ROM derived
algorithms. However, the resulting large number of degrees under steady mean flow.
of freedom and associated computational costs could inhibit POD-based ROMs are often limited in their application
the capability to perform additional analyses such us due to
implementing a control strategy or conducting uncertainty
quantification. This presents the need to develop tools to a its construction being computationally expensive for
significantly reduce the number of degrees of freedom, large-scale problems, especially for turbulent flows
simplify the complexity of fluid dynamical systems, capture b lack of robustness with varying flow or geometric
the main physics, and predict accurately flow quantities parameters in the ROMs.
such as force coefficients.
Most widely used reduced-order model (ROM) Lie and Farhat (2007) addressed the former issue and
techniques for systems involving fluid flows are based presented an interpolation method based on the Grassmann
on the proper orthogonal decomposition (POD) (Bakewell manifold and its tangent space at a point that is applicable
and Lumley, 1967; Sirovich, 1987; Deane et al., 1991; to structural, aerodynamic, aeroelastic and many other
Berkooz et al., 1993; Holmes et al., 1996). The most ROMs based on projection schemes. For the latter issue,
common implementation of this approach is to represent Amsallem and Farhat (2008) developed a strategy to
the flow properties by a collection of flow variables pre-compute bases with discrete flight parameters for
(velocity and pressure fields) of the system (commonly ROMs to predict the onset of flutter. However, it required
referred to as snapshots). These snapshots are originally training and interpolation of parameter set according to
obtained through numerical simulations or experiments. flight conditions; subsonic, transonic, or supersonic. As for
The POD is then used to find a low-dimensional set the incorporation of the effect of turbulence in POD-based
of basis functions by processing snapshot data (Holmes ROMs, various closure models and discretisation techniques
et al., 1996; Sirovich, 1987). These functions are then have been proposed and successfully implemented (Wang
used to derive a low-dimensional dynamical system that is et al., 2011, 2012; Akhtar et al., 2012).
typically obtained by Galerkin projection. For instance, a To cover small variations in geometry, viscosity or
ROM for the velocity field can be obtained by projecting mean flow, sensitivity analysis has been proposed and
the Navier-Stokes equations onto the space formed by used to increase the spectrum of POD-based ROMs to
the velocity POD modes. This projection yields a set varying flow and shape parameters. Hay et al. (2009)
of ordinary-differential equations (ODEs) that can be included sensitivity of flow parameter (viscosity) within
directly integrated to reconstruct the velocity field. The key the POD basis (extrapolated basis) and Galerkin projection
advantage of these models is their suitability to be applied (expanded basis) and showed that the ROMs developed for
for flow control design (Graham et al., 1999; Singh et al., the baseline flow, i.e., at a particular Reynolds number, also
2001; Arian et al., 2000; Bergmann et al., 2005; Siegel works for the variation in the Reynolds number (viscosity).
et al., 2006; Akhtar and Nayfeh, 2010). Later, Hay et al. (2010, 2012) used a similar analysis for
In a similar fashion, to develop a ROM for the pressure variation in geometry and angle of attack. Akhtar et al.
field, one could project, in incompressible flows, the (2010) developed a ROM for shape variation from circular
Poisson equation onto the pressure POD modes. The result to elliptic cylinders using sensitivity analysis.
is a set of algebraic equations that express the pressure
temporal coefficients as a function of the corresponding
temporal coefficients of the velocity field [see Akhtar et al.
370 M. Ghommem et al.

In the literature, the presence of multiple spurious is presented in Section 5. At the end, the effect of flow
limit cycles in reduced-order has been observed for unsteadiness is incorporated in the model and corresponding
long time integration. For instance, Foias et al. (1991) results are presented.
and Aubry et al. (1993) investigated the existence of
multiple spurious steady states in the Galerkin expansion
of the Kuramoto-Sivashinsky equation. For the flow 2 Numerical methodology
past a circular cylinder, Sirisup and Karniadakis (2004)
showed that the onset of divergence from the correct A parallel DNS code is used to simulate the incompressible
limit cycle depends on the number of modes in the flow past a two-dimensional circular cylinder (Akhtar,
Galerkin expansion, the Reynolds number, and the flow 2008; Akhtar et al., 2009b). The governing equations are
geometry. They used a spectral vanishing viscosity (SVV) the Navier-Stokes equations in curvilinear coordinates (ξm )
method (Tadmor, 1989), which adds a small amount of and can be written in strong conservation form as follows:
mode-dependent dissipation satisfying the entropy condition
while retaining the spectral accuracy. Akhtar et al. (2009b) ∂Um
= 0, (1)
used a shooting method (Nayfeh and Balachandran, 1995) ∂ξm
to home on the limit cycle by adjusting the initial
∂(J −1 ui ) ∂Fim
conditions so that the system does not drift to a spurious + = 0, (2)
limit cycle. ∂t ∂ξm
POD-based ROM technique has been also applied in where the flux is defined as
the field of aeroelasticity (Dowell et al., 1998; Hall et al.,
2000; Tang et al., 2001; Thomas et al., 2003). For instance, ∂ξm 1 mn ∂ui
Fim = Um ui + J −1 p− G . (3)
Dowell et al. (1998) developed a ROM for unsteady ∂xi Re ∂ξn
aerodynamic problems using eigen analysis. Thomas et al.
where i = 1, 2 for two-dimensional
( i) simulations. Here
(2003) applied this technique to a three-dimensional
Re = U∞ν D ; J −1 = det ∂x
∂ξj is the inverse of the Jacobian
wing configuration in transonic regime. They reduced
the direct numerical simulation (DNS) model with over or the volume of the cell; Um = J −1 ∂ξ ∂xj uj is the volume
m

seven million degrees of freedom to a few dozen degrees flux (contravariant velocity multiplied by J −1 ) normal to
of freedom while retaining the accuracy of unsteady the surface of constant ξm ; and Gmn = J −1 ∂ξ m ∂ξn
∂xj ∂xj is the
aerodynamics. They showed that POD vectors based on ‘mesh skewness tensor’.
one set of structural mode shapes can be used for different A non-staggered-grid layout is employed to solve
structural modes provided the solution snapshots at the the transformed Navier-Stokes equations. The Cartesian
endpoints of the desired frequency range are present in the velocity components (u, v) and pressure (p) are defined at
snapshot data. the centre of the control volume in the computational space
In our work, we focus on the pressure POD and and the volume fluxes (U, V ) are defined at the mid points
force coefficients in time domain in contrast to frequency of the corresponding faces. All of the spatial derivatives are
domain to develop a ROM. We follow a different approximated with second-order accurate central differences
approach and develop a low-dimensional tool to predict except for the convective terms, for which QUICK is
the effects of unsteadiness in the inflow on the force employed (Leonard, 1979).
coefficients acting on a circular cylinder. Particularly, A semi-implicit scheme is employed to advance the
the linear stochastic estimator (LSE) is implemented to solution in time. The diagonal viscous terms are advanced
directly relate temporal coefficients of the pressure on implicitly using the second-order accurate Crank-Nicolson
the surface of the cylinder to the coefficients of the method, whereas all of the other terms are advanced
velocity field. The outcome is a representation of the using the second-order accurate Adams-Bashforth method.
lift and drag coefficients in terms of the velocity and The Adams-Bashforth scheme was chosen because of its
pressure temporal coefficients that are obtained from a computational efficiency when coupled with the fractional
ROM. We also decompose the pressure POD modes on step method. For more details of the numerical algorithm,
the surface to determine force decomposition coefficients validation and verification, the reader is referred to Akhtar
(FDCs) for each mode, both for lift and drag coefficients. (2008) and Akhtar et al. (2009a, 2009b).
The implementation of the proposed tool is performed under In this study, we perform numerical simulation of the
both steady and unsteady mean flow conditions to show flow past a circular cylinder at Re = 180. We chose this
how mean flow unsteadiness can be represented within a Reynolds number since the effect of third-dimensionality
model reduction framework that makes use of steady flow initiates around a Reynolds number of 183 (Williamson,
parameters. 1996; Ma and Karniadakis, 2002). Mode A instability
The manuscript is organised as follows. Section 2 briefly causes spanwise disturbances of the order of four times
describes the numerical methodology to record flow data the cylinder diameter. As the Reynolds number crosses the
for the velocity and pressure fields. The technique of value of 260, Mode B instability causes the disturbance in
the low-dimensional modelling is outlined in Section 3. the spanwise direction of the order of the cylinder diameter.
Section 4 discusses in detail the development of ROM for Beyond Reynolds number of 500, the wake starts to become
the pressure while the procedure and significance of FDC turbulent and three dimensional results are dominant (Mittal
A low-dimensional tool for predicting force decomposition coefficients for varying inflow conditions 371

and Balachandar, 1996). Since our analysis is performed yield the pressure coefficients. The final outcome is a
on a two-dimensional flow, we keep the Reynolds number representation of the lift and drag coefficients in terms of
below the Mode A transition range. the velocity and pressure temporal coefficients.

3 Low-dimensional representation of flow quantities 3.1 The POD

A schematic of different components of the proposed tool A classical way to compute the POD modes is to
to derive a ROM for the force coefficients is shown in perform a singular value decomposition (SVD), however,
Figure 1. The approach combines POD-based models to this approach may have a limited application, especially
predict the flow field and LSE for approximating the when dealing with a large mesh size as in DNS cases
pressure on the cylinder. First, we generate instantaneous which require fine mesh and high computational cost.
solutions or snapshots of the velocity and pressure fields Alternatively, one could use the method of snapshots
from flow simulations. Then, we apply the method of (Sirovich and Kirby, 1987) which allows for a significant
snapshots to compute the velocity and pressure POD modes. reduction of the large datasets. In this method, sets
By projecting the Navier-Stokes equations onto the space of instantaneous solutions (or snapshots) of the flow
formed by the POD modes, we derive a ROM for the parameters (velocity and pressure fields) obtained from the
velocity field. The result is a set of ODEs in the form DNS are first generated and stored in N × n matrices,
of q̇ = F (q) that governs the variations of the temporal denoted by S where N and n denote, respectively, the
coefficients of the velocity field. To develop a ROM for the number of grid points and snapshots. Then, an eigen
pressure field, one approach is to project Poisson equation analysis of the matrix S T S is performed; that is (Sirovich
onto the space formed by the pressure POD modes (Akhtar and Kirby, 1987),
et al., 2008). The result is a set of algebraic equations which
relates the temporal velocity and pressure coefficients. As 1
S T S ∈ Rn×n : S T SVi = σi2 Vi and Ui = SVi
noted by Akhtar et al. (2008), this model works well σi
only for low Reynolds numbers and cannot be extended to
unsteady mean flow cases. In our approach, we generate The selection of POD modes is optimal in the sense that
the temporal coefficients of the pressure fluctuations using a the error between each snapshot and its projection on the
mapping function M (q) that is based on LSE technique. We space spanned by those modes is minimised (Burkardt et al.,
use these coefficients with pressure POD modes to compute 2002). Besides, the square of the singular values represent
the pressure distribution, which is integrated to determine a measure of the energy content of each POD mode and, as
the fluid loads. Again the idea here is to relate pressure and such, provide guidance for the number of modes that should
velocity coefficients through the LSE instead of using the be considered in order to capture the main physics of the
Poisson equation and then build a ROM that can directly system.

Figure 1 A descriptive schematic showing the implementation of the low-dimensional tool (see online version for colours)

Velocity field
CFD
Fluid Loads

Snapshots

SVD/Method of Pressure POD Fluid Loads


Snapshots modes
a
Validation

Velocity POD Galerkin Projection Reduced


modes -Order Model

Linear Stochastic
q Estimator

Velocity field
372 M. Ghommem et al.

The POD basis vectors Φi and their corresponding a pattern of pair similarities. Looking at the first six
singular values σi are obtained by applying the method of POD modes of the streamwise and crossflow velocity
snapshots as described above to data from flow simulations components in Figures 3 and 4, respectively, it is observed
over a circular cylinder. We note that the singular values that the POD modes ϕui (i = 1, 2, 5, 6) of the streamwise
obtained here are in fact the eigenvalues of the correlation component are antisymmetric with respect to the x-axis;
matrix obtained S T S. They are related to the singular that is,
values obtained in the SVD of the snapshot matrix. In
Figure 2(a), we plot the singular values of the first ϕui (x, −y) = −ϕui (x, y), (4a)
20 modes obtained from a steady mean flow condition. All while the ϕui (i = 3, 4) are symmetric; that is,
singular values are normalised with the corresponding value
of the first mode (i.e., σi /σ1 ). We observe that the two ϕui (x, −y) = ϕui (x, y). (4b)
modes of each pair have values of the same order. These
values decrease from one pair to the next in approximately On the other hand, the POD modes ϕvi of the crossflow
a geometric progression. component are symmetric for i = 1, 2, 5, 6 with respect to
the x-axis; that is,
Figure 2 Normalised singular values and cumulative
energy content of velocity POD modes,
ϕvi (x, −y) = ϕvi (x, y), (5a)
(a) normalised singular values of velocity POD and antisymmetric for i = 3, 4; that is,
modes (b) normalised cumulative energy of velocity
POD modes (see online version for colours) ϕvi (x, −y) = −ϕvi (x, y). (5b)

0
10
3.2 ROM of the velocity field

10
−2
To obtain the ROM, we expand the fluctuations in a
Galerkin representation in terms of the POD modes (Φi );
−4 i.e., we let
σ

10


M
10
−6
u(x, t) ≈ u(x) + qi (t)Φi (x) (6)
i=1
−8
10
0 5 10 15 20 where u(x, t) = (u, v)T , u(x) is the mean flow,
Number of modes
Φi (x) = (ϕui (x), ϕvi (x))T , and M is the number of POD
(a)
modes used in the projection. Then, we project the
1
Navier-Stokes equations onto the space formed by the POD
modes as
0.9 ∂u 1 2
< Φk , + (u.∇)u + ∇p − ∇ u >= 0, (7)
Cumulative energy

∂t Re
0.8 ∫
where < F, G >= Ω (F · G)dΩ represents the inner
0.7 product between F and G and k = 1, 2, ..., M . Substituting
equation (6) into equation (7) and based on the
0.6 orthogonality of the POD basis functions, one obtains a set
of M ordinary differential equations (Akhtar et al., 2009b)
0.5
0 5 10
Number of modes
15 20

M ∑
M ∑
M
q̇k (t) = Ak + Bki qi (t) + Ckij qi (t)qj (t), (8)
(b)
i=1 i=1 j=1

To gain insight into the contribution of each mode to the where


total energy
∑M of the system, we define the cumulative energy 1
as E = i=1 σi and the cumulative contribution of the first
∑j Ak = < Φk , ∇2 u > − < Φk , u.∇u >,
Re
j modes as cj = ( i=1 σi )/E. We remark that most of the
Bki = − < Φk , u.∇Φi > − < Φk , Φi .∇u >
energy is contained in the first few modes [see Figure 2(b)].
1
More specifically, the first ten modes contain more than + < Φk , ∇2 Φi >
99.9% of the total flow energy. Re
Ckij = − < Φk , Φi ∇Φj > .
The streamwise and crossflow velocity modes in the
steady flow are plotted in Figures 3 and 4. These modes
constitute the dominant spatial structures of the flow field.
Similarly to the singular values, the POD modes show
A low-dimensional tool for predicting force decomposition coefficients for varying inflow conditions 373

Figure 3 Contours of the streamwise velocity POD modes (upper frame: ϕui , i = 1, 2, centre frame: ϕui , i = 3, 4,
and lower frame: ϕui , i = 5, 6) for the case of steady mean flow (see online version for colours)

Figure 4 Contours of the crossflow velocity POD modes (upper frame: ϕvi , i = 1, 2, centre frame: ϕvi , i = 3, 4,
and lower frame: ϕvi , i = 5, 6) for the case of steady mean flow (see online version for colours)
374 M. Ghommem et al.

In equation (8), A is an M × 1 vector resulting from the 4 ROM for the pressure field
average flow field, B is the linear part of the dynamical
system, which originates from the interaction of the average In incompressible flows, the Poisson equation relates the
field with the eigenfunctions and the linear dissipative velocity and the pressure fields. In the fractional-step
operator, and C is a tensor representing the quadratic method used in the current numerical methodology, most of
non-linearity of the Navier-Stokes equations (Noack et al., the computation time at every time iteration is consumed in
2005; Akhtar et al., 2009b). It should be noted here that solving the pressure-Poisson equation. POD-based models
using Green’s theorem and the divergence-free property and do not lend themselves easily for a direct computation
for the case of p = 0 on the outer flow boundary Ωso , of the loads. Thus, a modest approach is to apply the
the pressure term is not part of equation (8) (Ma and POD technique on the pressure field and later integrate
Karniadakis, 2002). The POD eigenfunctions are identically the pressure field onto surface to compute the loads.
zero on the inflow boundary because the average flow Hence, similar to the velocity field, the pressure field
is subtracted from the total flow. However, in case of is decomposed into mean and fluctuating components.
Neumann boundary conditions on Ωso , the contribution of Expanding the latter in terms of the pressure POD modes
the pressure term is not exactly zero for the cylinder wake. Ψi , we write
At 10D downstream from the cylinder (D is the cylinder
diameter), i.e., the outflow boundary, the pressure effects of ∑
M

the wake are minimal. As such, their representative terms p(x, t) ≈ p̄(x) + ai (t)Ψi (x), (10)
can be assumed to vanish in the ROM (Noack et al., 2005). i=1

To assess the validity of the ROM, we project each snapshot


where p̄(x) is the average pressure. Using equations (10)
onto the POD modes and compute the following inner
and (6) in the pressure-Poisson equation leads to an
product
algebraic ROM for the pressure (Akhtar et al., 2009b).
q̃i (tj ) =< Φi (x), (u(x, tj ) − u(x)) > . (9) The model works well for low Reynolds number steady
cases. However, explicit non-linear dependence on the
We use the first ten POD modes to develop a ROM for velocity field may result in errors at the higher Reynolds
the case of steady mean flow. We perform the Galerkin numbers where more modes may be needed in the model.
projection and compute Ak , Bki , and Ckij , where k, m, n = Furthermore, we were unsuccessful in our initial effort to
1, 2, ..., 10 in equation (8). Thus, the DNS with 193 × 257 extend the use of the ROM developed under a particular
degrees of freedom is reduced to a dynamical system flow condition (e.g., steady mean flow) to other flow and
with ten dimensions. Then, we integrate equation (8) and geometric conditions. To overcome these issues, we use
perform a two-dimensional projection of the phase portrait stochastic estimation techniques to estimate the pressure
on the plane (q1 , q2 ) and compare it with the projection field as discussed below.
obtained on the (q̃1 , q̃2 ) plane in Figure 5. Clearly, the two
projections of the ROM and the computed modes from the
numerical simulations are in good agreement. 4.1 Pressure POD modes

Figure 5 A two-dimensional projection of the POD phase We compute the pressure POD mode from the pressure
portrait on (q1 , q2 )-plane (solid) and snapshot snapshot data obtained in a similar manner as in the
portrait (stars) (see online version for colours) case of the velocity field. However, we are interested
only in the distribution of the POD pressure modes on
the surface of the cylinder rather than the whole domain.
3 Using the method of snapshots, we also obtain the singular
values for the pressure. Unlike the velocity POD modes,
2 the lower singular values of the pressure do not occur
in pairs while the higher ones do occur in pairs, as
1 shown in Figure 6(a). From Figure 6(b), we can clearly
observe that the first pressure POD mode is the most
dominant and contains more than 85% of the cumulative
q2

energy contribution. Figure 7 shows the distribution


−1 of the pressure POD modes on the cylinder surface. In
Figure 7(a), we plot the average pressure over the cylinder
−2 surface. We note that θ = 0◦ corresponds to the base point
and θ = 180◦ corresponds to the stagnation point on the
−3 cylinder. Figures 7(b) to 7(f) show the first five pressure
−3 −2 −1 0
q1
1 2 3 POD modes on the cylinder surface. From Figure 7(b), we
observe that the absolute maxima of the first POD mode
occur at θ ≈ 90◦ and θ ≈ 270◦ . For the higher modes,
Note: 10 POD modes are kept in the truncated expansion given the absolute maxima occur close to the base point. Based
by equation (6). on the singular values, the contributions of the higher
A low-dimensional tool for predicting force decomposition coefficients for varying inflow conditions 375

modes are much less than that of the first mode. As such, estimator (QSE), and so on (Taylor and Glauser, 2004;
Akhtar and Nayfeh (2010) used this result and showed that Ausseur et al., 2006). We note that the use of linear or QSE
placing fluidic suction actuators at the location of maximum would depend on the level of non-linearity involved in the
pressure of the first POD mode could produce significant dynamical process. In our study, considering only the linear
reduction in the lift and drag coefficients. terms in equation (11) was found enough to accurately
estimate the fluid forces as will be shown below.
Figure 6 Normalised singular values and cumulative energy Denoting by ãi the temporal coefficients obtained by
content of pressure POD modes, projecting each pressure snapshots onto the POD mode Ψi ;
(a) normalised singular values of pressure POD that is,
modes (b) normalised cumulative energy of pressure
POD modes (see online version for colours) ãi (t) =< Ψi (x), (p(x, t) − p̄(x)) >, (12)
we compute the mapping matrix L by minimising the mean
0
10 square error defined as (Taylor and Glauser, 2004; Ausseur
et al., 2006)
10
−2
⟨( )2 ⟩
Ean = an (t) − ãn (t) , (13)
T
−4 ⟨ ⟩ ∫T
σ

10
where (.), (.) = 0 (.)(.)dt. The solution of this
T
−6 minimisation problem is obtained by solving the following
10
linear system:
⟨ ⟩ ⟨ ⟩ ⟨ ⟩ 
10
−8
q 1 , q1 q 1 , q2 . . . q1 , q M  
0 5 10 15 20
 T T T 
Ln1
Number of modes
 .. .. ..   .. 
 . ⟩ ⟨ . ⟩ . ⟩   . 
(a) ⟨ ⟨ 
qM , q 1 qM , q2 . . . qM , qM LnM
1
⟨ ⟩ 
T T T

0.98 ãn , q1
 T 
 .. 
= 
Cumulative energy

0.96
. (14)
⟨ ⟩ 
0.94
ãn , qM
T
0.92

0.9 We consider the first eight POD modes (which contain more
than 99% of the total energy) and compare in Figure 8
0.88
the estimated pressure POD temporal coefficients obtained
0.86 from LSE with those obtained by projecting snapshots data
0 5 10 15 20
Number of modes for the pressure as given by equation (12). A good match
(b) between the two sets of data shows the efficiency of LSE
to estimate the pressure POD temporal coefficients.
4.2 Linear stochastic estimator

For each mode n, the estimated pressure POD temporal 5 Force decomposition coefficients
coefficients are described as a series expansion of velocity
POD temporal coefficients which can be determined from a The forces (lift and drag) applied on the cylinder can be
ROM or a snapshot projection as shown above; that is, written as the addition of the contributions of the pressure
and shear stresses acting on the surface of the cylinder; that
an (t) = Lni qi (t) + Qnij qi (t)qj (t)
is (Akhtar et al., 2009b),
+ Cnijk qi (t)qj (t)qk (t) + · · · (11)
∫ 2π ( 1 ∂v(x, t)
Akhtar et al. (2009b) developed the pressure POD ROM by CL (t) = − p(x, t) sin(θ) − (
0 Re ∂x
projecting the pressure Poisson equation onto the pressure )
∂u(x, t)
POD modes. The algebraic model obtained comprised − ) cos(θ) dθ
linear and non-linear components. Time histories and phase ∂y
∫ 2π (
portraits of POD-basis coefficients of the velocity qi and
CD (t) = − p(x, t) cos(θ) (15)
pressure ai fields showed that the linear dependence is 0
dominant than the quadratic component. Thus, equation (11) 1 ∂v(x, t) ∂u(x, t) )
can be truncated to include just the linear term and + ( − ) sin(θ) dθ
Re ∂x ∂y
neglecting higher-order terms. This technique is called
LSE. Including quadratic term leads to quadratic stochastic
376 M. Ghommem et al.

Figure 7 Pressure POD modes on the surface of the cylinder for steady mean flow condition

0.1
0.6

0.4
0.05
0.2

ψ
0
ψ

-0.2
-0.05
-0.4

-0.6 -0.1

90 180 270 90 180 270


θ θ
(a) (b)

0.2
0.2

0.15
0.15

0.1 0.1

0.05 0.05
ψ

0
ψ

-0.05 -0.05

-0.1 -0.1

-0.15 -0.15

-0.2 -0.2
90 180 270 90 180 270
θ θ
(c) (d)

0.3
0.3

0.2
0.2

0.1
0.1
ψ

0
ψ

-0.1 -0.1

-0.2 -0.2

-0.3 -0.3
90 180 270 90 180 270
θ θ
(e) (f)
A low-dimensional tool for predicting force decomposition coefficients for varying inflow conditions 377

Figure 8 Pressure coefficients for steady mean flow conditions where i = 1, 2, · · · , M.


obtained from snapshot projection and LSE An important feature of FDC is that it reduces the
(see online version for colours) spatial representation of a POD mode, contributing in
pressure and shear stresses, to scalar quantities in lift and
drag directions. Thus, ROM is developed for both lift and
Snapshot projection drag coefficients. Equation (16) provides a representation of
2 1
LSE
the lift and drag forces in terms the velocity and pressure
temporal coefficients. The time-varying coefficients q ′ s are
1

2
0 0
a

a
obtained through a ROM as given by equation (8), and
−2
400 402 404 406
−1
400 402 404 406 the coefficients a′ s are obtained using the LSE technique.
0.5 0.5 Table 1 provides values for the lift and drag FDC,
(Lpi and Dpi ), for the first six modes. We note that Lp0 and
a4
3

0 0
Dp0 are the mean values of the lift and drag coefficients,
a

−0.5 −0.5
respectively. It is interesting to note that Lp0 = 0 which
400 402 404 406 400 402 404 406 suggests that there is zero mean lift on the cylinder.
0.05 0.05 Similarly, Dp0 = 1.0146 pertains to the mean pressure drag
of the circular cylinder. Clearly, the parameters used in
5

6
a

0 0
a

reduced-order have physical attributes and signify flow


−0.05 −0.05 physics. Furthermore, from the values given in Table 1
400 402 404 406 400 402 404 406
and the amplitudes of the a′ s as shown in Figure 8, one
0.05 0.05
could identify the relative contribution of the pressure POD
modes to the FDC. Clearly, the first few pressure POD
a8
7

0 0
a

modes contribute the most to the lift and drag coefficients.


−0.05
400 402 404 406
−0.05
400 402 404 406
In Figure 9, we plot the lift and drag coefficients
Time Time obtained from DNS and reconstructed from the
low-dimensional tool by combining POD and LSE. We
note that ten POD modes are used in equation (16) to
compute the lift and drag coefficients. A good agreement
Substituting equations (6) and (10) into equation (15),
between these data is clearly observed. This demonstrates
and integrating over the cylinder surface, we obtain a
the capability of the developed tool to accurately predict
representation of the lift and drag forces in terms of FDC
the loads.
given by
Table 1 Lift and drag FDC

M ∑
M
CL (t) = Lpi ai (t) + Lvi qi (t) Lift FDC
i=0 i=0 Lp0 0

M ∑
M Lp1 –0.314
CD (t) = Dpi ai (t) + Dvi qi (t) (16) Lp2 –0.210
i=0 i=0 Lp3 –0.0007
Lp4 0.0002
where
∫ 2π ( )
Lp5 –0.0451
Lp6 0.0517
Lp0 = − p̄(x) sin(θ) dθ,
0 Drag FDC
∫ 2π ( )
Dp0 1.0146
Dp 0 = − p̄(x) cos(θ) dθ,
0 Dp1 –0.0001
∫ 2π ( ) Dp2 0.0002
Lpi = − Ψi (x) sin(θ) dθ, Dp3 0.0750
0
∫ 2π ( ) Dp4 0.1315
Dp i = − Ψi (x) cos(θ) dθ Dp5 –0.0006
0 Dp6 0.0003
∫ 2π ( )
1 ∂v̄(x) ∂ ū(x)
Lv0 = ( − ) cos(θ) dθ,
0 ReD ∂x ∂y
∫ 2π ( )
1 ∂v̄(x) ∂ ū(x)
Dv 0 =− ( − ) sin(θ) dθ, 6 Modelling of unsteady flows
0 ReD ∂x ∂y
∫ 2π ( )
1 ∂ϕv (x) ∂ϕu (x) Flow quantities such as aerodynamic or hydrodynamic
Lvi = ( − ) sin(θ) dθ,
0 ReD ∂x ∂y loads, heat transfer and noise radiation can be significantly
∫ 2π ( )
1 ∂ϕ (x) ∂ϕu (x)
v affected by mean flow unsteadiness which is inevitably a
Dvi =− ( − ) cos(θ) dθ, part of many natural flows and is encountered in many
0 ReD ∂x ∂y
378 M. Ghommem et al.

technological applications. In this section, we expand our We emphasise the use of the ΦSi since it constitutes
approach to develop ROMs that cover unsteady flow cases. a key contribution of this effort, at least for the velocity
field. Computing new POD modes for the modified
Figure 9 Lift and drag coefficients computed from the ROM flow field, i.e., unsteady inflow, and developing a ROM
and DNSs, (a) lift coefficient (b) drag coefficient based on modified POD modes is not the objective of
(see online version for colours) this work. In fact, it is one of the drawbacks of the
POD-based ROMs which has been addressed by developing
0.8 ROM the ROM for unsteady inflow using the POD basis of
DNS the steady flow field. In our opinion, this approach will
0.6
add robustness to the POD-based ROMs and increase the
0.4 spectrum of its application. As a first step in modelling
the inflow unsteadiness, we have kept the frequency of
0.2
CL

disturbance (fm ) close to the shedding frequency (fs ).


0 Adding unsteadiness to the mean flow that has multiple
frequencies will render the flow field complex and is
−0.2
beyond the scope of current paper and can be addressed in
−0.4 a future effort using a similar approach.
Substituting equation (17) into the Navier-Stokes
160 165 170 175 180 185 190
T ime equations, expanding each term, and projecting these
equations on the ΦSi , one obtains
(a)


M M ∑
∑ M

ROM q̇k (t) = Ak + Bki qi (t) + Ckij qi (t)qj (t)


1.06
DNS i=1 i=1 j=1


M
1.04 + Hk γ̇(t) + Iki qi (t)γ(t) (18)
i=1
CD

1.02 + Kk γ(t) + Lk γ 2 (t)

1 where

0.98 Hk = − < ΦSk , uS >,


160 165 170 175 180 185 190
Iki = − < ΦSk , uS .∇ΦSi > − < ΦSk , ΦSi .∇uS >
T ime 1
Kk = < ΦSk , ∇2 uS > −2 < ΦSk , uS .∇uS >,
(b) ReD
Lk = − < ΦSk , uS .∇uS > .

In the unsteady case, we add to the mean flow a sinusoidal It is important to note that equation (17) contains two
component with a specific amplitude A/U∞ = 0.0532 time-dependent states; qi and γ. Therefore, we obtain two
(approximately 5% disturbance in the freestream direction) time derivatives in the modified reduced order model. The
and non-dimensional frequency fm = 0.2052 that is convection and diffusion terms of Navier-Stokes equations
considered close to the shedding frequency of the steady give rise to linear and quadratic terms in qi , γ, and their
flow fs = 0.1862; that is, combinations in the model. In the modified reduced order
model given by equation (18), the first part is the same
U∞ (t) = U∞ + γ(t) as the usual ROM for steady flow and the second part is
A an explicit function of time that models the unsteadiness
= U∞ (1 + sin(2πfm t)) in the mean flow. We use 10 POD modes (i.e., M = 10)
U∞
and integrate equation (18) simultaneously and compute
To enable modelling of the unsteady mean flow in terms of the qk (t). Then, we recreate the velocity field using
the steady flow POD modes, we expand the velocity field equation (17). In Figure 10, we plot the streamwise velocity
as at the probe in the vortex shedding region of the cylinder as
obtained from DNS and reconstructed from the POD-based

M reduced order model given by equation (18). The fairly
u(x, t) ≈ uS (x)(1 + γ(t)) + qi (t)ΦSi (x) (17) good match between the two sets of data shows that the
i=1 reduced order model captures the main features of the
velocity field.
where ΦSi are the POD modes of the steady mean flow uS .
A low-dimensional tool for predicting force decomposition coefficients for varying inflow conditions 379

Figure 10 Streamwise velocity at the probe obtained from relatively low Reynolds number flow regimes, higher-order
DNS and POD-based ROM (unsteady mean stochastic estimators might be needed to perform similar
flow case) (see online version for colours) analysis for high Reynolds number cases where non-linear
flow aspects are more pronounced.

Figure 11 Lift and drag coefficients for the unsteady inflow


ROM obtained from the ROM and DNS,
DNS (a) lift coefficient (b) drag coefficient (see online
0.6
version for colours)
0.4

0.2
UP

0.8 ROM
0 DNS
0.6
−0.2
0.4

−0.4
0.2

CL
0
300 310 320 330 340 350 360 370 380 390 400
T ime
−0.2

−0.4
Note: Ten POD modes are kept in the truncated expansion given
by equation (17). −0.6

80 90 100 110 120 130 140


T ime
We attempted to follow a similar approach to use the steady
POD modes to compute the mapping function which relates (a)
the pressure and velocity fields. However, we found that the
use of the steady pressure POD modes is not appropriate 1.3 ROM
to capture the unsteadiness effects. Alternatively, we used DNS
the unsteady pressure POD modes and equations (16) and 1.2
(17) to reconstruct the lift and drag coefficients for the
unsteady mean flow case. We plot, in Figure 11, the lift 1.1
CD

and drag coefficients for the unsteady case obtained from


DNS and reconstructed from the low-dimensional tool using 1
the unsteady POD modes. A good match between the two
sets of data is observed and shows the capability of the 0.9
developed tool to reconstruct the force coefficients even for
the unsteady mean flow case. Noise in the drag coefficient 0.8
80 90 100 110 120 130 140
is attributed to numerical errors which may be present T ime
due to grid or numerical scheme. However, the mean and (b)
fluctuating amplitudes of drag coefficient match the results
in the literature (Akhtar, 2008; Akhtar et al., 2009b).
References
Akhtar, I. (2008) Parallel Simulations, Reduced-Order Modeling,
7 Conclusions
and Feedback Control of Vortex Shedding using Fluidic
Actuators, PhD thesis, Virginia Tech, Blacksburg, VA.
In this effort, we have developed a low-dimensional tool to
Akhtar, I., Borggaard, J. and Hay, A. (2010) ‘Shape sensitivity
predict the effects of unsteadiness in the inflow on force analysis in flow models using a finite-difference approach’,
coefficients acting on a circular cylinder. The approach Mathematical Problems in Engineering, Article ID 209780,
is based on combining POD and LSE techniques. We doi:10.1155/2010/209780, pp.1–22.
used POD to derive a ROM of the velocity field. Then, Akhtar, I., Marzouk, O.A. and Nayfeh, A.H. (2009a) ‘A van
to overcome the difficulty of developing a ROM using der Pol-Duffing oscillator model of hydrodynamic forces
Poisson’s equation, we have related the pressure field to the on canonical structures’, Journal of Computational and
velocity field through a mapping function based on LSE. Nonlinear Dynamics, Vol. 4, No. 4, p.041006.
Results for force coefficients obtained from this approach Akhtar, I. and Nayfeh, A.H. (2010) ‘Model based control
were verified by comparison with results from DNS. This of vortex shedding using fluidic actuators’, Journal of
finding shows the capability of the developed tool to Computational and Nonlinear Dynamics, Vol. 5, No. 4,
accurately predict the loads and then, could serve as the p.041015.
basis to apply additional analyses such as implementing
control strategies. Although LSE has been successfully used
to produce accurate estimation of force coefficients for
380 M. Ghommem et al.

Akhtar, I., Nayfeh, A.H. and Ribbens, C.J. (2008) ‘A Galerkin Hay, A., Akhtar, I. and Borggaard, J.T. (2012) ‘On the use
model of the pressure field in incompressible flows’, of sensitivity analysis in model reduction to predict flows
in Proceedings of the 46th Aerospace Sciences Meeting and for varying inflow conditions’, International Journal for
Exhibit. AIAA Paper No. 2008-611. Numerical Methods in Fluids, Vol. 68, No. 1, pp.122–134.
Akhtar, I., Nayfeh, A.H. and Ribbens, C.J. (2009b) ‘On the Hay, A., Borggaard, J., Akhtar, I. and Pelletier, D. (2010)
stability and extension of reduced-order Galerkin models in ‘Reduced-order models for parameter dependent geometries
incompressible flows: a numerical study of vortex shedding’, based on shape sensitivity analysis’, Journal of
Theoretical and Computational Fluid Dynamics, Vol. 23, Computational Physics, Vol. 229, No. 4, pp.1327–1352.
No. 3, pp.213–237. Hay, A., Borggaard, J.T. and Pelletier, D. (2009) ‘Local
Akhtar, I., Wang, Z., Borggaard, J. and Iliescu, T. (2012) improvements to reduced-order models using sensitivity
‘A new closure strategy for proper orthogonal decomposition analysis of the proper orthogonal decomposition’, Journal of
reduced-order models’, Journal of Computational and Fluid Mechanics, Vol. 629, pp.41–72.
Nonlinear Dynamics, Vol. 7, No. 3, p.034503. Holmes, P., Lumley, J.L. and Berkooz, G. (1996) Turbulence,
Amsallem, D. and Farhat, C. (2008) ‘Interpolation method Coherent Structures, Dynamical Systems and Symmetry,
for adapting reduced-order models and application to Cambridge University Press, Cambridge, UK.
aeroelasticity’, AIAA Journal, Vol. 46, No. 7, pp.1803–1813. Leonard, B.P. (1979) ‘A stable and accurate convective modeling
Arian, E., Fahl, M. and Sachs, E. (2000) ‘Trust-region procedure based on quadratic upstream interpolation’,
proper orthogonal decomposition for optimal flow control’, Computational Methods in Applied Mechanical Engineering,
Technical Report ICASE. Vol. 19, No. 1, pp.59–98.
Aubry, N., Lian, W.Y. and Titi, E.S. (1993) ‘Preserving Lie, T. and Farhat, C. (2007) ‘Adaptation of aeroelastic
symmetries in the proper orthogonal decomposition’, SIAM reduced-order models and application to an f-16
Journal of Scientific Computing, Vol. 14, No. 2, pp.483–505. configuration’, AIAA Journal, Vol. 45, No. 6, pp.1244–1257.
Ausseur, J.M., Pinier, J.T., Glausser, M.N. and Higuchi, H. Ma, X. and Karniadakis, G. (2002) ‘A low-dimensional model
(2006) ‘Experimental development of a reduced-order model for simulating three-dimensional cylinder flow’, Journal of
for flow separation control’, in Proceedings of the AIAA Fluid Mechanics, Vol. 458, pp.181–190.
30th Aerospace Sciences Meeting and Exhibit, AIAA Paper Mittal, R. and Balachandar, S. (1996) ‘Direct numerical
No. 2006-1251. simulation of flow past elliptic cylinders’, Journal of
Bakewell, H.P. and Lumley, J.L. (1967) ‘Viscous sublayer and Computational Physics, Vol. 124, No. 2, pp.351–367.
adjacent wall region in turbulent pipe flow’, The Physics of Nayfeh, A.H. and Balachandran, B. (1995) Applied Nonlinear
Fluids, Vol. 10, No. 9, pp.1880–1889. Dynamics: Analytical, Computational, and Experimental
Bergmann, M., Cordier, L. and Brancher, J.P. (2005) ‘Optimal Methods, Wiley, New York, NY, pp. 449–454.
rotary control of the cylinder wake using POD reduced-order Noack, B.R., Papas, P. and Monkewitz, P.A. (2005) ‘The need for
model’, Physics of Fluids, Vol. 17, No. 9, pp.097101–20. a pressure-term representation in empirical Galerkin models
Berkooz, G., Holmes, P. and Lumley, J.L. (1993) ‘The proper of incompressible shear flows’, Journal of Fluid Mechanics,
orthogonal decomposition in the analysis of turbulent flows’, Vol. 523, pp.339–365.
Annual Review of Fluid Mechanics, Vol. 25, pp.539–575. Siegel, S., Cohen, K. and McLaughlin, T. (2006) ‘Numerical
Burkardt, J., Du, Q., Gunzburger, M.D. and Lee, H.C. (2002) simulations of a feedback-controlled circular cylinder wake’,
‘Reduced order modeling of complex systems’, in Biennal AIAA Journal, Vol. 44, No. 6, pp.1266–1276.
Conferences on Numerical Analysis, University of Dundee, Singh, S.N., Myatt, J.H., Addington, G.A., Banda, S. and
Scotland, UK. Hall, J.K. (2001) ‘Optimal feedback control of vortex
Deane, A.E., Kevrekidis, I.G., Karniadakis, G.E. and Orsag, S.A. shedding using proper orthogonal decomposition models’,
(1991) ‘Low-dimensional models for complex geometry Journal of Fluids Engineering, Vol. 123, No. 3, pp.612–618.
flows: application to grooved channels and circular cylinder’, Sirisup, S. and Karniadakis, G.E. (2004) ‘A spectral viscosity
Physics of Fluids A, Vol. 3, No. 10, pp.2337–2354. method for correcting the long-term behavior of pod
Dowell, E.H., Hall, K.C. and Romanowski, M.C. (1998) models’, Journal of Computational Physics, Vol. 194, No. 1,
‘Eigenmode analysis in unsteady aerodynamics: pp.92–116.
reduced-order models’, Applied Mechanics Review, Vol. 50, Sirovich, L. (1987) ‘Turbulence and the dynamics of coherent
No. 6, pp.371–385. structures’, Quarterly of Applied Mathematics, Vol. 45,
Foias, C., Jolly, M.S., Kevrekidis, I.G. and Titi, E.S. (1991) pp.561–590.
‘Dissipativity of the numerical schemes’, Nonlinearity, Sirovich, L. and Kirby, M. (1987) ‘Low-dimensional procedure
Vol. 4, No. 3, pp.591–613. for the characterization of human faces’, Journal of Optical
Graham, W.R., Peraire, J. and Tang, K.Y. (1999) ‘Optimal Society of America A, Vol. 4, No. 3, pp.529–524.
control of vortex shedding using low-order models. Tadmor, E. (1989) ‘Convergence of spectral methods for
Part I – open-loop model development’, International nonlinear conservation laws’, SIAM Journal of Numerical
Journal for Numerical Methods in Engineering, Vol. 44, Analysis, Vol. 26, No. 1, p.30.
No. 7, pp.945–972. Tang, D., Kholodar, D., Junag, J-N. and Dowell, E.H. (2001)
Hall, K.C., Thomas, J.P. and Dowell, E.H. (2000) ‘Proper ‘System identification and proper orthogonal decomposition
orthogonal decomposition technique for transonic unsteady method applied to unsteady aerodynamics’, AIAA Journal,
aerodynamic flows’, AIAA Journal, Vol. 38, No. 10, Vol. 39, No. 8, pp.1569–1576.
pp.1853–1862.
A low-dimensional tool for predicting force decomposition coefficients for varying inflow conditions 381

Taylor, J.A. and Glauser (2004) ‘Towards practical flow sensing Wang, Z., Akhtar, I., Borggaard, J. and Iliescu, T. (2012) ‘Proper
and control via pod and LSE based low-dimensional tools’, orthogonal decomposition closure models for turbulent flows:
Journal of Fluids Engineering, Vol. 126, No. 3, pp.337–345. a numerical comparison’, Computer Methods in Applied
Thomas, J.P., Dowell, E.H. and Hall, K.C. (2003) ‘Three Mechanics and Engineering, Vol. 237–240, No. 1, pp.10–26.
dimensional transonic aeroelasticity using proper orthogonal Williamson, C.H.K. (1996) ‘Vortex dynamics in the cylinder
decomposition-based reduced-order models’, Journal of wake’, Annual Review of Fluid Mechanics, Vol. 28,
Aircraft, Vol. 40, No. 3, pp.544–551. pp.477–539.
Wang, Z., Akhtar, I., Borggaard, J. and Iliescu, T. (2011)
‘Two-level discretizations of nonlinear closure models for
proper orthogonal decomposition’, Journal of Computational
Physics, Vol. 230, No. 1, pp.126–146.

You might also like