1 Introduction
Successful control of vortex-induced vibration (VIV) can lead to safer and more cost-effective structures in offshore, aeronautical and civil engineering. In the past few decades, various passive control techniques (Owen & Bearman Reference Owen and Bearman2001; Choi, Jeon & Kim Reference Choi, Jeon and Kim2008; Baek & Karniadakis Reference Baek and Karniadakis2009; Yu et al.
Reference Yu, Xie, Yan, Constantinides, Oakley and Karniadakis2015; Law & Jaiman Reference Law and Jaiman2017) have been explored via geometry modification and by adding auxiliary surfaces to alter the flow dynamics without any energy input. While passive VIV control methods offer some simplicity, they do not have the ability to work on demand and may not be effective from the perspective of wake stabilization, drag reduction and VIV suppression in a wide range of operational conditions. More importantly, passive devices, e.g. strakes, splitter plates or fairings, are not easy to implement for certain situations such as square-shaped multicolumn offshore platforms (Chakrabarti Reference Chakrabarti2005) and subsea pipelines undergoing VIV in the proximity of the seabed (Sumer & Fredsoe Reference Sumer and Fredsoe1997). The motivation of this study is to develop a feedback active control algorithm based on a reduced-order model and to demonstrate the performance of the algorithm in stabilizing the wake flow and the VIV for a canonical two-dimensional circular cylinder problem. To investigate the proposed control scheme, we consider the flow past a stationary circular cylinder at low Reynolds number,
$Re=UD/\unicode[STIX]{x1D708}$
, and the vibrating cylinder as a function of the reduced velocity,
$U_{r}=U/(f_{N}D)$
, and mass ratio,
$m^{\ast }=4m/(\unicode[STIX]{x1D70C}\unicode[STIX]{x03C0}D^{2})$
, with zero mass-damping parameter,
$m^{\ast }\unicode[STIX]{x1D701}=0.0$
. Here,
$U$
,
$\unicode[STIX]{x1D70C}$
,
$m$
,
$D$
,
$f_{N}$
,
$\unicode[STIX]{x1D701}$
and
$\unicode[STIX]{x1D708}$
are the free-stream velocity, the density of fluid, the mass of cylinder, the diameter of the cylinder, the structural natural frequency, the damping ratio and the kinematic viscosity respectively.
Through external input of a small amount of tunable energy into the surrounding flow, active VIV control techniques offer a better alternative due to their adaptive and efficient performance. Kim & Choi (Reference Kim and Choi2005) utilized blowing/suction slots placed on the top and bottom of the circular cylinder to stabilize the unstable wake. The applied forcing in the slots was sinusoidal along the spanwise direction but kept steady in time. In another interesting study by Dong, Triantafyllou & Karniadakis (Reference Dong, Triantafyllou and Karniadakis2008), combined steady windward suction and leeward blowing (WSLB) was found to be an effective strategy to eliminate the vortex street and to suppress VIV in the cross-flow direction. This WSLB method in Dong et al. (Reference Dong, Triantafyllou and Karniadakis2008) is essentially an active flow control (i.e. external energy is required to maintain steady suction and blowing), but can be customized as a passive technique by deploying porous surfaces to form connecting channels between the windward and leeward stagnation points of a circular cylinder. The aforementioned blowing/suction control schemes provide only open-loop alteration of unsteady flow, whereby the control input is prescribed and is independent of the flow states. In another recent study, Wang et al. (Reference Wang, Tang, Yu, Duan, Wang, Tang, Yu and Duan2016) utilized windward-suction and leeward-blowing feedback control to suppress VIV in both the inline and cross-flow directions at a specific condition
$(m^{\ast },U_{r})=(2,5)$
. However, it is not certain whether the system is free of VIV for a range of reduced velocity
$U_{r}$
and mass ratio
$m^{\ast }$
. Moreover, a simplified closed-loop control strategy based on the proportional (P), integral (I) and proportional–integral (PI) schemes was applied to manipulate the WSLB velocities through the standard deviation of the surrounding flow velocity. The PI control was found to outperform the P and I control schemes with respect to the effectiveness of VIV suppression. To the best of our knowledge, feedback control of VIV via vertical suction and blowing has not been studied previously. Furthermore, feedback control of VIV based on model reduction, which can be important in both numerical and experimental settings, has not been explored in earlier studies.
Various types of active flow control strategies have been explored for the flow past a circular cylinder, such as full-state feedback control, neural networks and proportional closed-loop feedback control. In particular, active feedback or closed-loop control of unsteady flow over a bluff body has recently been investigated via numerical simulations (Ahuja & Rowley Reference Ahuja and Rowley2010; Parkin, Thompson & Sheridan Reference Parkin, Thompson and Sheridan2014; Flinois & Morgans Reference Flinois and Morgans2016). To implement the optimal linear control in an efficient manner, it is imperative to develop a low-order linear model by retaining the significant dynamics of the original system (Ahuja & Rowley Reference Ahuja and Rowley2010; Flinois & Morgans Reference Flinois and Morgans2016). A linear reduced-order model (ROM) provides a way to trace the eigenspectrum of a dynamical system while maintaining approximately an order of magnitude efficiency improvement to construct the essential dynamics of the system. Instead of deriving the system model directly from the linearized governing equations of the fluid flow, the system identification method is particularly desirable, because it is non-intrusive and can be used directly in an existing Navier–Stokes (NS) solver to generate approximate system matrices using only the input–output data sequences. The eigensystem realization algorithm (ERA) (Juang & Pappa Reference Juang and Pappa1985) is a well-established system identification method for the construction of a linear ROM from the stable system linear impulse response. Ma, Ahuja & Rowley (Reference Ma, Ahuja and Rowley2011) proved that the ROM constructed by the ERA is mathematically equivalent to balance truncation. Recently, Flinois, Morgans & Schmid (Reference Flinois, Morgans and Schmid2015) provided a mathematically rigorous proof that unmodified balance truncation (designed for a stable system) is applicable for an unstable system. Based on the work of Ma et al. (Reference Ma, Ahuja and Rowley2011) and Flinois et al. (Reference Flinois, Morgans and Schmid2015), the ERA is used for the stabilization of unstable wake flow (Flinois & Morgans Reference Flinois and Morgans2016).
Although the suction/blowing control strategy has been extensively studied for unstable flow and VIV past a circular cylinder both numerically (Kim & Choi Reference Kim and Choi2005; Dong et al. Reference Dong, Triantafyllou and Karniadakis2008; Mao, Blackburn & Sherwin Reference Mao, Blackburn and Sherwin2015; Wang et al. Reference Wang, Tang, Yu, Duan, Wang, Tang, Yu and Duan2016) and experimentally (Fransson, Konieczny & Alfredsson Reference Fransson, Konieczny and Alfredsson2004; Chen et al. Reference Chen, Xin, Xu, Li, Ou and Hu2013, Reference Chen, Gao, Yuan, Li and Hu2015), the earlier research relies on open-loop control strategies, except for the recent study of Wang et al. (Reference Wang, Tang, Yu, Duan, Wang, Tang, Yu and Duan2016) via a simplified non-adaptive control procedure. Ahuja & Rowley (Reference Ahuja and Rowley2010) and Flinois & Morgans (Reference Flinois and Morgans2016) designed a feedback control law to stabilize unstable wake flow over a plate and bluff body using a model reduction method; however, the actuation was modelled as a body force and could not be implemented as a practical actuator. The recognition of a linear mechanism during the self-sustaining behaviour of VIV in our recent study (Yao & Jaiman Reference Yao and Jaiman2017) has inspired us to develop a feedback suppression strategy based on a linear control theory. Motivated by the insight about the frequency lock-in process during VIV, it is possible to design a linear controller to minimize the unsteady vortex shedding and the VIV effects. For that purpose, a control input based on vertical blowing and suction at the surface can be utilized (Kim & Choi Reference Kim and Choi2005).
The key contribution of this paper is to develop an active feedback blowing and suction (AFBS) procedure based on an ERA-based ROM to control the wake instability and the VIV. Using the ERA method, a low-order fluid model is constructed, and is coupled with the structure via the linear quadratic regulator (LQR) optimal control scheme. The proposed AFBS procedure ensures a VIV free system within a large parameter space of reduced velocity (
$U_{r}$
) and mass ratio (
$m^{\ast }$
), and it can handle both one degree-of-freedom (1-DOF) and two degree-of-freedom (2-DOF) VIV. The results from the ROM solver are compared with those of full-order simulations based on the incompressible NS equations. We employ the ROM model to predict the performance of the AFBS procedure through the eigenvalue distribution in the complex plane. The present study is based on two questions pertaining to the VIV physics and control. (i) Can we suppress self-sustaining VIV by assuming a linear lock-in mechanism? (ii) How much additional energy is required for the suppression of VIV compared with the stationary nonlinear vortex shedding state?
The problem set-up for feedback control of an unsteady wake and VIVs: (a) computational domain and boundary conditions for the flow past a freely vibrating cylinder in uniform flow; (b) the proposed new actuator configuration
$BS0$
with blowing/suction over the porous surface of the circular cylinder. The positive control input is defined as suction from the bottom and blowing at the top of the cylinder.
The organization of the article is as follows. The numerical details of the full-order model based on the NS equations and the ERA-based ROM are presented in § 2. Section 3 introduces the active feedback control based on the vertical suction and blowing for the unstable wake flow of a stationary cylinder and VIV. A sensitivity study on the number of suction/blowing actuators, the angular arrangement of actuators and the combined versus independent control architecture is also carried out in § 3. Section 4 summarizes the main conclusions.
2 Numerical methodology
2.1 Problem set-up
Figure 1(a) shows a schematic of the problem set-up used in our simulation study for a flexibly mounted circular cylinder in a flowing stream. At the inlet boundary,
$\unicode[STIX]{x1D6E4}_{in}$
, a stream of incompressible fluid enters the domain at a horizontal velocity
$(u,v)=(U,0)$
, where
$u$
and
$v$
denote the streamwise and transverse velocities in the
$x$
and
$y$
directions respectively. For the VIV configuration, a circular cylinder with mass
$m$
is elastically mounted on a linear spring and is allowed to vibrate only in the transverse direction. A no-slip wall condition is implemented on the surfaces of the bluff body and a traction-free boundary condition is implemented along the outlet,
$\unicode[STIX]{x1D6E4}_{out}$
, while a slip wall condition is implemented on the top,
$\unicode[STIX]{x1D6E4}_{top}$
, and bottom,
$\unicode[STIX]{x1D6E4}_{bottom}$
, boundaries. Except where stated otherwise, all length scales are normalized by the cylinder diameter
$D$
and velocities by the free-stream velocity
$U$
. The numerical domain extends from
$-10D$
at the inlet to
$30D$
at the outlet and from
$-15D$
to
$15D$
in the transverse direction. For the flow control, we consider blowing and suction on the porous cylinder surface through fluidic actuators. It is known that moderate levels of suction/blowing into the surrounding flow can have a great impact on the boundary layer, the separation point and the wake characteristics (Fransson et al.
Reference Fransson, Konieczny and Alfredsson2004; Chen et al.
Reference Chen, Xin, Xu, Li, Ou and Hu2013, Reference Chen, Gao, Yuan, Li and Hu2015). Through the active feedback control, the suction mechanism can delay the separation of the boundary layer (i.e. narrower wake width and reduced drag), whereas blowing can have the opposite effect. On the other hand, continuous blowing generally tends to decrease the Strouhal number for the flow around a porous cylinder, while continuous suction has the opposite influence on the vortex shedding frequency (Fransson et al.
Reference Fransson, Konieczny and Alfredsson2004).
In the present study, we propose a feedback control based on a configuration with three pairs of suction/blowing actuators, as depicted in figure 1(b). In this proposed configuration, termed as
$BS0$
, there are a total of six suction and blowing slots distributed over the cylinder surface, with pairs of slots at the windward,
$\unicode[STIX]{x1D703}=(135^{\circ },225^{\circ })$
, midward,
$\unicode[STIX]{x1D703}=(90^{\circ },270^{\circ })$
, and leeward,
$\unicode[STIX]{x1D703}=(45^{\circ },315^{\circ })$
, sides. Here,
$\unicode[STIX]{x1D703}$
is the deviation angle between the centreline of each suction/blowing actuator and the base suction point. As shown in figure 1(b), the positive control input is defined as suction from the bottom of the cylinder and blowing at the top surface. Similarly to Pastoor et al. (Reference Pastoor, Henning, Noack, King and Tadmor2008), we consider the actuation slot width to be
$\unicode[STIX]{x1D70E}_{c}=\unicode[STIX]{x03C0}D/72$
, whereby the energy supply from the actuation is characterized by the momentum coefficient as
$C_{\unicode[STIX]{x1D707}}=2N\unicode[STIX]{x1D70C}V_{c}^{2}\unicode[STIX]{x1D70E}_{c}/(\unicode[STIX]{x1D70C}U^{2}D)$
, where
$N$
is the number of slots and
$V_{c}$
denotes the time-dependent suction and blowing velocity. Owing to the body-conforming Lagrangian–Eulerian coupling for fluid–structure interaction, the actuation conditions for the blowing and suction are accurately enforced by the Dirichlet boundary condition. The present full-order fluid–structure model relies on a variational finite-element formulation and a semi-discrete time stepping. While the NS equations are discretized in space using
$\mathbb{P}_{n}/\mathbb{P}_{n-1}$
isoparametric finite elements for the fluid velocity and pressure, the second-order backward scheme is used for the time discretization, where
$\mathbb{P}_{n}$
denotes the standard
$n\text{th}$
order Lagrange finite-element space on the discretized fluid domain. Details of the numerical techniques, the verification and the mesh convergence study based on
$\mathbb{P}_{2}/\mathbb{P}_{1}$
isoparametric elements are documented in Yao & Jaiman (Reference Yao and Jaiman2017).
2.2 Full-order model
For the sake of completeness, we first present the full-order model (FOM) based on the NS equations for the moving incompressible viscous fluid domain
$\unicode[STIX]{x1D6FA}^{f}(t)$
as
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D70C}\left(\frac{\unicode[STIX]{x2202}\boldsymbol{u}}{\unicode[STIX]{x2202}t}\bigg|_{\unicode[STIX]{x1D74C}}+(\boldsymbol{u}-\boldsymbol{w})\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{u}\right)=\unicode[STIX]{x1D735}\boldsymbol{\cdot }\unicode[STIX]{x1D748}\quad \text{on }\unicode[STIX]{x1D6FA}^{f}(t), & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{u}=0\quad \text{on }\unicode[STIX]{x1D6FA}^{f}(t), & \displaystyle\end{eqnarray}$$
where the time derivative is taken with the referential coordinate
$\unicode[STIX]{x1D74C}$
held fixed and the Cauchy stress tensor for a Newtonian fluid is
$\unicode[STIX]{x1D748}=-p\unicode[STIX]{x1D644}+\unicode[STIX]{x1D707}(\unicode[STIX]{x1D735}\boldsymbol{u}+(\unicode[STIX]{x1D735}\boldsymbol{u})^{\text{T}})$
. Here,
$p$
,
$\boldsymbol{u}$
,
$\boldsymbol{w}$
,
$\unicode[STIX]{x1D707}$
and
$\unicode[STIX]{x1D644}$
denote the fluid pressure, the fluid velocity, the mesh velocity, the dynamic viscosity and the identity tensor respectively. The mesh nodes on the fluid domain
$\unicode[STIX]{x1D6FA}^{f}(t)$
are updated by solving a linear steady pseudo-elastic material model
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D735}\boldsymbol{\cdot }\unicode[STIX]{x1D748}^{m}=\mathbf{0}, & & \displaystyle\end{eqnarray}$$
where
$\unicode[STIX]{x1D748}^{m}$
is the stress experienced by the arbitrary Lagrangian–Eulerian (ALE) mesh due to the strain induced by the rigid-body movement, which is defined as
$$\begin{eqnarray}\displaystyle \qquad \unicode[STIX]{x1D748}^{m}=(1+k_{m})[(\unicode[STIX]{x1D735}\boldsymbol{\unicode[STIX]{x1D702}}^{\,f}+(\unicode[STIX]{x1D735}\boldsymbol{\unicode[STIX]{x1D702}}^{\,f})^{\text{T}})+(\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{\unicode[STIX]{x1D702}}^{\,f})\unicode[STIX]{x1D644}], & & \displaystyle\end{eqnarray}$$
where
$\boldsymbol{\unicode[STIX]{x1D702}}^{\,f}$
represents the ALE mesh node displacement and
$k_{m}$
is a mesh stiffness variable chosen as a function of the element area to limit the distortion of small elements located in the immediate vicinity of the fluid–body interface Liu, Jaiman & Gurugubelli (Reference Liu, Jaiman and Gurugubelli2014).
Given a base flow
$\boldsymbol{u}_{0}$
, the corresponding linearized NS equations can be written in a semi-discrete form as
$$\begin{eqnarray}\unicode[STIX]{x1D640}\frac{\text{d}\unicode[STIX]{x1D64C}}{\text{d}t}=\unicode[STIX]{x1D641}\unicode[STIX]{x1D64C}+\unicode[STIX]{x1D642}\unicode[STIX]{x1D651}_{m},\end{eqnarray}$$
where the matrices and vectors in (2.5) are
$$\begin{eqnarray}\unicode[STIX]{x1D641}=\left(\begin{array}{@{}cc@{}}-(\,)\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{u}_{0}-\boldsymbol{u}_{0}\boldsymbol{\cdot }\unicode[STIX]{x1D735}(\,)+\unicode[STIX]{x1D707}(\unicode[STIX]{x1D735}(\,)+\unicode[STIX]{x1D735}^{\text{T}}(\,)) & -\unicode[STIX]{x1D735}(\,)\\ \unicode[STIX]{x1D735}\boldsymbol{\cdot }(\,) & \mathbf{0}\end{array}\right),\end{eqnarray}$$
$$\begin{eqnarray}\unicode[STIX]{x1D642}=\left(\begin{array}{@{}cc@{}}(\,)\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{u}_{0} & \mathbf{0}\\ \boldsymbol{ 0} & \mathbf{0}\end{array}\right),\quad \unicode[STIX]{x1D640}=\left(\begin{array}{@{}cc@{}}\unicode[STIX]{x1D644} & \mathbf{0}\\ \mathbf{0} & \mathbf{0}\end{array}\right),\quad \unicode[STIX]{x1D64C}=\left(\begin{array}{@{}c@{}}\boldsymbol{u}\\ p\end{array}\right),\quad \unicode[STIX]{x1D651}_{m}=\left(\begin{array}{@{}c@{}}\boldsymbol{w}\\ \boldsymbol{ 0}\end{array}\right).\end{eqnarray}$$
The linearized NS equation for a stationary cylinder is obtained by setting the mesh velocity
$\boldsymbol{w}=\mathbf{0}$
. After the discretization, the generalized eigenvalue problem of the linearized NS equation can be written as
$(\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D65B}}+\unicode[STIX]{x1D70E}\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D65B}})\unicode[STIX]{x1D665}=\mathbf{0}$
, where the non-symmetric matrices
$\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D65B}}$
and
$\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D65B}}$
result from the spatial and temporal discretizations,
$\unicode[STIX]{x1D70E}$
denotes the eigenvalue of the discretized system and
$\unicode[STIX]{x1D665}$
is the right eigenvector (forward modes). The corresponding discrete adjoint problem can be obtained as
$\unicode[STIX]{x1D666}(\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D65B}}+\unicode[STIX]{x1D70E}\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D65B}})=\mathbf{0}$
, where
$\unicode[STIX]{x1D666}$
is the left eigenvector of the discrete system and represents the approximation of the adjoint modes (Giannetti & Luchini Reference Giannetti and Luchini2007). The linearized NS equations are solved by the semi-discrete variational procedure, which is employed for the nonlinear fluid–structure equations in Liu et al. (Reference Liu, Jaiman and Gurugubelli2014) and Jaiman, Sen & Gurugubelli (Reference Jaiman, Sen and Gurugubelli2015). The fluid–structure coupling is achieved through a partitioned staggered procedure (Jaiman et al.
Reference Jaiman, Geubelle, Loth and Jiao2011). We next present the ERA-based ROM via the input–output dynamics of the FOM.
2.3 Feedback control via the ROM
The ERA-based ROM is constructed by the linear input/output dynamics of the NS equations given in (2.1) and (2.2). The input vector
$\unicode[STIX]{x1D66A}=[Y,V_{c}]^{\text{T}}$
for the fluid system is the transverse displacement
$Y$
and the suction and blowing vertical velocity
$V_{c}$
, while the output is the total lift coefficient
$C_{l}$
over the structural body. The fluid ERA-based ROM formulated in state-space form at discrete times
$t=k\unicode[STIX]{x0394}t$
,
$k=0,1,2,\ldots ,$
with a constant sampling time
$\unicode[STIX]{x0394}t$
reads
$$\begin{eqnarray}\left.\begin{array}{@{}c@{}}\displaystyle \unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D65B}}(k+1)=\unicode[STIX]{x1D63C}\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D65B}}(k)+\unicode[STIX]{x1D63D}\unicode[STIX]{x1D66A}(k),\\ \displaystyle C_{l}(k)=\unicode[STIX]{x1D63E}\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D65B}}(k)+\unicode[STIX]{x1D63F}\unicode[STIX]{x1D66A}(k),\end{array}\right\}\end{eqnarray}$$
where
$\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D65B}}$
is an
$n_{r}$
-dimensional state vector and the integer
$k$
is a sample index for the time stepping. The system matrices are
$(\unicode[STIX]{x1D63C},\unicode[STIX]{x1D63D},\unicode[STIX]{x1D63E},\unicode[STIX]{x1D63F})$
, which are obtained by the ERA method. To construct the ERA-based ROM, the impulse response of the NS equation is first defined as
$\unicode[STIX]{x1D66E}$
, based on which the generalized block Hankel matrix
$r\times s$
can be constructed as
$$\begin{eqnarray}\unicode[STIX]{x1D643}(k-1)=\left[\begin{array}{@{}cccc@{}}\unicode[STIX]{x1D66E}_{k+1} & \unicode[STIX]{x1D66E}_{k+2} & \cdots \, & \unicode[STIX]{x1D66E}_{k+s}\\ \unicode[STIX]{x1D66E}_{k+2} & \unicode[STIX]{x1D66E}_{k+3} & \cdots \, & \unicode[STIX]{x1D66E}_{k+s+1}\\ \vdots & \vdots & \ddots & \vdots \\ \unicode[STIX]{x1D66E}_{k+r} & \unicode[STIX]{x1D66E}_{k+r+1} & \cdots \, & \unicode[STIX]{x1D66E}_{k+(s+r-1)}\end{array}\right],\end{eqnarray}$$
and by applying singular value decomposition of the Hankel matrix
$\unicode[STIX]{x1D643}(0)$
as
$$\begin{eqnarray}\unicode[STIX]{x1D643}(0)=\unicode[STIX]{x1D650}\unicode[STIX]{x1D72E}\unicode[STIX]{x1D651}^{\ast }=\left[\unicode[STIX]{x1D650}_{1}\quad \unicode[STIX]{x1D650}_{2}\right]\left[\begin{array}{@{}cc@{}}\unicode[STIX]{x1D72E}_{1} & \mathbf{0}\\ \mathbf{0} & \unicode[STIX]{x1D72E}_{2}\end{array}\right]\left[\begin{array}{@{}c@{}}\unicode[STIX]{x1D651}_{1}^{\ast }\\ \unicode[STIX]{x1D651}_{2}^{\ast }\end{array}\!\right],\end{eqnarray}$$
where the diagonal matrix
$\unicode[STIX]{x1D72E}$
is the Hankel singular values (HSVs). The block matrix
$\unicode[STIX]{x1D72E}_{2}$
contains the zeros or negligible elements. By truncating the dynamically less significant states, we estimate
$\unicode[STIX]{x1D643}(0)\approx \unicode[STIX]{x1D650}_{1}\unicode[STIX]{x1D72E}_{1}\unicode[STIX]{x1D651}_{1}^{\ast }$
. The state-space matrices
$(\unicode[STIX]{x1D63C},\unicode[STIX]{x1D63D},\unicode[STIX]{x1D63E},\unicode[STIX]{x1D63F})$
are obtained by
$$\begin{eqnarray}\left.\begin{array}{@{}c@{}}\displaystyle \unicode[STIX]{x1D63C}=\unicode[STIX]{x1D72E}_{1}^{-1/2}\unicode[STIX]{x1D650}_{1}^{\ast }\unicode[STIX]{x1D643}(1)\unicode[STIX]{x1D651}_{1}\unicode[STIX]{x1D72E}_{1}^{-1/2},\\ \displaystyle \unicode[STIX]{x1D63D}=\unicode[STIX]{x1D72E}_{1}^{1/2}\unicode[STIX]{x1D651}_{1}^{\ast }\unicode[STIX]{x1D640}_{\unicode[STIX]{x1D662}},\\ \displaystyle \unicode[STIX]{x1D63E}=\unicode[STIX]{x1D640}_{\unicode[STIX]{x1D669}}^{\ast }\unicode[STIX]{x1D650}_{1}\unicode[STIX]{x1D72E}_{1}^{1/2},\\ \displaystyle \unicode[STIX]{x1D63F}=\unicode[STIX]{x1D66E}_{1}.\end{array}\!\right\}\end{eqnarray}$$
Here,
$\unicode[STIX]{x1D640}_{\unicode[STIX]{x1D662}}^{\ast }=[\unicode[STIX]{x1D644}_{q}\mathbf{0}]_{q\times N}$
and
$\unicode[STIX]{x1D640}_{\unicode[STIX]{x1D669}}^{\ast }=[\unicode[STIX]{x1D644}_{p}\mathbf{0}]_{p\times M}$
, where
$N=s\times q$
,
$M=r\times p$
and the
$\unicode[STIX]{x1D644}_{p,q}$
are the identity matrices. The matrices
$\unicode[STIX]{x1D63D}$
and
$\unicode[STIX]{x1D63F}$
can be rewritten as
$\unicode[STIX]{x1D63D}=[\unicode[STIX]{x1D63D}_{Y},\unicode[STIX]{x1D63D}_{V_{c}}]$
and
$\unicode[STIX]{x1D63F}=[D_{Y},D_{V_{c}}]$
, where the subscripts denote the input components defined in the vector
$\unicode[STIX]{x1D66A}$
. The ERA-based ROM is constructed in the vicinity of a given base flow at
$t=0$
(
$k=0$
) and the impulse signal starts from
$t=\unicode[STIX]{x0394}t$
(
$k=1$
).
The VIV system is simplified to a transversely vibrating circular cylinder with one degree of freedom (1-DOF), and the non-dimensional structural equation in the state-space form is written as (Yao & Jaiman Reference Yao and Jaiman2016)
$$\begin{eqnarray}\dot{\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}}}=\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}}\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}}+\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}}C_{l}.\end{eqnarray}$$
The state matrices and vectors are
$$\begin{eqnarray}\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}}=\left[\begin{array}{@{}cc@{}}0 & 1\\ -(2\unicode[STIX]{x03C0}F_{s})^{2} & -4\unicode[STIX]{x1D701}\unicode[STIX]{x03C0}F_{s}\end{array}\right],\quad \unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}}=\left[\begin{array}{@{}c@{}}0\\ {\displaystyle \frac{2}{\unicode[STIX]{x03C0}m^{\ast }}}\end{array}\right],\quad \unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}}=\left[\begin{array}{@{}c@{}}Y\\ {\dot{Y}}\end{array}\right],\end{eqnarray}$$
where
$Y$
is the transverse displacement,
$C_{l}$
is the lift coefficient and
$F_{s}$
is the non-dimensional reduced structural frequency defined as
$F_{s}=f_{N}D/U=1/U_{r}$
. Based on the ERA-based ROM, the resulting closed-loop system for VIV can be formulated by coupling the structure (2.11) and the fluid system (2.7) as
$$\begin{eqnarray}\left.\begin{array}{@{}l@{}}\displaystyle \unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k+1)=\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k)+\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}\unicode[STIX]{x1D66A}_{c}(k),\\ \displaystyle \unicode[STIX]{x1D66E}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k+1)=\unicode[STIX]{x1D63E}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k)+\unicode[STIX]{x1D63F}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}\unicode[STIX]{x1D66A}_{c}(k),\end{array}\!\!\right\}\end{eqnarray}$$
where
$$\begin{eqnarray}\left.\begin{array}{@{}c@{}}\underset{(n_{r}+2)\times (n_{r}+2)}{\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}}=\left[\begin{array}{@{}cc@{}}\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D659}}+\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D659}}D_{Y}\unicode[STIX]{x1D63E}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D659}} & \unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D659}}\unicode[STIX]{x1D63E}\\ \unicode[STIX]{x1D63D}_{Y}\unicode[STIX]{x1D63E}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D659}} & \unicode[STIX]{x1D63C}\end{array}\right],\quad \underset{(n_{r}+2)\times 2}{\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}}=\left[\begin{array}{@{}cc@{}}\mathbf{0} & \unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D659}}D_{V_{c}}\\ \mathbf{0} & \unicode[STIX]{x1D63D}_{V_{c}}\end{array}\right],\\ \underset{3\times (n_{r}+2)}{\unicode[STIX]{x1D63E}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}}=\left[\begin{array}{@{}cc@{}}\unicode[STIX]{x1D644} & \mathbf{0}\\ D_{Y}\unicode[STIX]{x1D63E}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D659}} & \unicode[STIX]{x1D63E}\end{array}\right],\quad \underset{3\times 2}{\unicode[STIX]{x1D63F}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}}=\left[\begin{array}{@{}cc@{}}\mathbf{0} & \mathbf{0}\\ 0 & D_{V_{c}}\end{array}\right].\end{array}\right\}\end{eqnarray}$$
The structural time discrete matrices are defined as
$\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D659}}=\text{e}^{\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}}\unicode[STIX]{x0394}t}$
,
$\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D659}}={\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}}}^{-1}(\text{e}^{\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}}\unicode[STIX]{x0394}t}-\unicode[STIX]{x1D644})\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}}$
and
$\unicode[STIX]{x1D63E}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D659}}=[1,0]$
, where the state vector
$\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}=[\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}},\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D65B}}]^{\text{T}}$
is an
$(n_{r}+2)$
-dimensional vector and the output vector is defined as
$\unicode[STIX]{x1D66E}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}=[\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}},C_{l}]^{\text{T}}$
. The present ERA-based ROM allows access to the most controllable and observable modes of the coupled fluid and structure system. With regard to control design, the ROM is valid in the neighbourhood of the unstable steady state.
As shown in figure 2, the optimal feedback control based on the LQR is utilized to calculate the optimal gain matrix
$\unicode[STIX]{x1D646}$
such that the state-feedback law
$\unicode[STIX]{x1D66A}_{c}(k)=-\unicode[STIX]{x1D646}\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k)$
minimizes the quadratic cost function for the discrete system as
$$\begin{eqnarray}J=\mathop{\sum }_{k=1}^{\infty }\{\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}^{\ast }\unicode[STIX]{x1D64C}\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}+\unicode[STIX]{x1D66A}_{c}^{\ast }\unicode[STIX]{x1D64D}\unicode[STIX]{x1D66A}_{c}\},\end{eqnarray}$$
where
$\unicode[STIX]{x1D64C}$
and
$\unicode[STIX]{x1D64D}$
set the relative weights of state deviation and input usage respectively and the asterisk
$^{\ast }$
denotes the transpose of a matrix. The
$\unicode[STIX]{x1D64C}$
is chosen as the identity matrix
$\unicode[STIX]{x1D644}$
for simplicity, whereas
$\unicode[STIX]{x1D64D}=c\unicode[STIX]{x1D644}>0$
provides the input to the cost function
$J$
. The coefficient
$c>0$
gives a relative weighting of output and input norms and can be tuned for an optimal tradeoff between the efficiency of VIV regulation and the energy input effort. The control input is defined as blowing and suction with vertical velocity magnitude
$\unicode[STIX]{x1D66A}_{c}=(0,V_{c})$
over the surface of a two-dimensional cylinder. For the LQR algorithm, the full-state vector
$\unicode[STIX]{x1D66D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k)$
should be observable to calculate the control input
$\unicode[STIX]{x1D66A}_{c}$
, while this is not always possible in practice. By assuming the measurements of force, velocity or acceleration on the body-mounted sensors, we can estimate the full-state vector via the Kalman filter as
$$\begin{eqnarray}\hat{\unicode[STIX]{x1D66D}}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k+1)=\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}\hat{\unicode[STIX]{x1D66D}}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k)+\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}\unicode[STIX]{x1D66A}_{c}(k)+\unicode[STIX]{x1D647}[\unicode[STIX]{x1D66E}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k)-\unicode[STIX]{x1D63E}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}\hat{\unicode[STIX]{x1D66D}}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k)-\unicode[STIX]{x1D63F}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}\unicode[STIX]{x1D66A}_{c}(k)],\end{eqnarray}$$
where
$\hat{\unicode[STIX]{x1D66D}}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}$
is the minimum mean-square estimator of the state vector and
$\unicode[STIX]{x1D647}$
is the filter gain matrix. After this, the control input vector
$\unicode[STIX]{x1D66A}_{c}$
can be determined as
$\unicode[STIX]{x1D66A}_{c}(k)=-\unicode[STIX]{x1D646}\hat{\unicode[STIX]{x1D66D}}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}(k)$
, resulting in a closed-loop VIV system as given in (2.13). Although the formulation presented here is general for any vibrating structure, results are presented for a canonical bluff body of a circular cylinder. We next demonstrate the designed controller for the unstable flow past a stationary circular cylinder and the feedback control of a freely vibrating cylinder.
Feedback control of VIV using the ROM: schematic of closed-loop control
$G_{CL}$
with the ERA-based ROM, where
$\unicode[STIX]{x1D6FF}$
and
$\unicode[STIX]{x1D646}$
represent the impulse input and gain matrix respectively; Kalman denotes the filter defined in (2.16).
The leading POD mode at
$Re=60$
: (a) streamwise velocity and (b) cross-stream velocity. The contour levels are from
$-0.01$
to
$0.01$
in increments of
$0.0025$
.
3 Results
3.1 Active feedback control of unsteady wake flow
We first demonstrate the feedback control scheme for the flow past a stationary circular cylinder at
$Re=60$
with
$BS0$
. The system input is the blowing and suction vertical velocity
$V_{c}$
and the output is the fluctuating lift coefficient
$C_{l}$
. The corresponding fluid ROM can be obtained by setting
$\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D654}}=\mathbf{0}$
and
$D_{Y}=0$
in (2.7). To start with the ERA-based ROM construction, the unstable steady state (i.e. the base flow) is first computed by a fixed point iteration without the time-dependent term in (2.1). For this purpose, 800 impulse response outputs (
$C_{l}$
) are stacked at each time step
$\unicode[STIX]{x0394}t=0.05$
by imposing an impulse of
$\unicode[STIX]{x1D6FF}(t)=10^{-4}$
on the blowing and suction vertical velocity
$V_{c}$
. Subsequently, the ERA-based ROM is obtained by performing a singular value decomposition of a
$500\times 200$
Hankel matrix. The order of the ERA-based ROM is determined by examining the singular values of the Hankel matrix. The linearity of the impulse response outputs (
$C_{l}$
) is confirmed by comparing the impulse responses of two different values,
$\unicode[STIX]{x1D6FF}(t)=10^{-4}$
and
$10^{-3}$
. To visualize the most excited flow structures, the leading proper orthogonal decomposition (POD) mode shown in figure 3 is extracted via the POD method. For this purpose, the snapshots of the flow field are stacked at each time step during the impulse response modelling. The aforementioned ERA-based ROM construction procedure can be considered as an open-loop identification process. Based on the wavemaker region (Giannetti & Luchini Reference Giannetti and Luchini2007), we determine the locations with high sensitivity and strong response, which are computed by taking the pointwise product of the forward and adjoint global modes, as shown in figure 4(a,b). The modes are obtained directly by solving a generalized eigenvalue problem of the linearized NS equation (2.5) in the neighbourhood of the base flow. As shown in figure 4(c), the body-mounted force sensor and the suction and blowing actuators are within the wavemaker region, which allows the feedback control to produce the largest drift of unstable observable and controllable modes.
The spatial distribution of the flow field: (a) forward and (b) adjoint velocity amplitudes, and (c) wavemaker region for a circular cylinder at
$Re=60$
. In (c), actuation slots are shown as triangles and the body-mounted force sensor as a red line over the cylinder for the
$BS0$
configuration. In (a) and (b), the contour levels are from
$0.002$
to
$0.018$
in increments of
$0.002$
. In (c), the contour levels are from
$0.02$
to
$0.2$
in increments of
$0.01$
.
The impulse response of a stationary cylinder in the
$BS0$
configuration: temporal variation of (a) the lift coefficient
$C_{l}$
and (b) the control input
$V_{c}$
predicted by the ROM and compared with the FOM at
$Re=60$
and
$c=10^{2}$
. The controller is switched on after
$tU/D=50$
convective time units.
Once the ERA-based ROM has been constructed, the LQR method is employed to design the optimal feedback gain
$\unicode[STIX]{x1D646}$
for different values of the weighting parameter
$c=10^{2},10^{3},10^{4}$
. As shown in figure 2, the closed-loop model can be directly constructed through the ERA-based ROM. The effectiveness of the closed-loop ROM with
$c=10^{2}$
is examined by comparing the temporal variations of
$C_{l}$
and
$V_{c}$
against the FOM counterparts, as illustrated in figure 5. An excellent match is found between the ROM with a number of modes
$n_{r}=33$
and the FOM, thereby confirming the accuracy and convergence of the present ERA-based ROM. The impulse response of the open-loop system is rapidly attenuated once the controller is switched on at
$t\geqslant 50$
.
The eigenvalues are calculated by
$\unicode[STIX]{x1D706}=\log (\text{eig}(\unicode[STIX]{x1D63C}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}-\unicode[STIX]{x1D63D}_{\unicode[STIX]{x1D668}\unicode[STIX]{x1D65B}}\unicode[STIX]{x1D646}))/\unicode[STIX]{x0394}t$
, where
$\unicode[STIX]{x0394}t=0.05$
is used for all of the computations. Figure 6(a) shows a comparison of eigenvalue distributions between the open-loop (OL) and closed-loop (CL) systems. As expected, a pair of eigenvalues of the OL system is in the unstable right half-plane as
$Re=60>Re_{cr}$
(where the critical Reynolds number is
$Re_{cr}\approx 46.8$
). On the other hand, the eigenvalues of the CL system are all in the stable left half-plane, which further confirms that the CL system is stable via the process of placing the poles in the stably damped locations in the complex plane. As
$c$
decreases, figure 6(a) also indicates that the eigenvalue moves further leftward and the CL system becomes more stable. However, a smaller value of
$c$
introduces a larger control input
$V_{c}$
; thus, there is a tradeoff between aggressive control input and rapid suppression of wake instability.
Feedback control of a stationary cylinder at
$Re=60$
with
$BS0$
: (a) eigenspectrum for OL and CL with different values of
$c$
; time variation of (b)
$C_{l}$
and (c)
$C_{d}$
(base flow drag subtraction); (d) control input
$V_{c}$
corresponding to CL response of the full nonlinear system. While the system has an impulse at
$t=0$
, the controller is switched on after
$tU/D=175.5$
.
The stabilization effect on the vortex shedding of a stationary circular cylinder at
$Re=60$
with
$BS0$
and
$c=10^{2}$
. Snapshots of spanwise vorticity contours at
$tU/D=$
(a) 175, (b) 225, (c) 250 and (d) 325. The contour levels are from
$-1$
to
$1$
in increments of
$0.1$
.
To demonstrate the AFBS controller to stabilize the saturated vortex street, figure 6(b,c) shows the effective attenuation of the fluctuating lift
$C_{l}$
and drag
$C_{d}$
with the control turned on at
$tU/D=175.5$
. It is expected that
$c=10^{4}$
requires the smallest control input
$V_{c}$
and the longest time scale to reach the target state, followed by
$c=10^{3}$
and
$c=10^{2}$
. The corresponding values of
$C_{\unicode[STIX]{x1D707}}$
are
$0.086$
,
$0.22$
and
$0.52$
for
$c=10^{4}$
,
$10^{3}$
and
$10^{2}$
respectively. The flow-field evolution during the suppression process is illustrated in figure 7.
3.2 Sensitivity study for unsteady wake flow control
Before proceeding to the application of the proposed ERA-based active jet control for a VIV system, we present a parametric study of a set of representative actuator configurations based on suction/blowing pairs, the effect of an angular arrangement of actuators and the combined versus independent control system architectures.
3.2.1 Effect of actuator configurations
To compare the effectiveness of the reference case
$BS0$
, we consider three additional configurations, as depicted in figure 8, to analyse the active feedback control based on the ERA-based ROM. The configuration
$BS1$
has only midward actuation slots, whereas the configuration
$BS2$
has the slots on the leeward side. By removing the windward slots from
$BS0$
at
$\unicode[STIX]{x1D703}=(135^{\circ },225^{\circ })$
, the new configuration of
$BS3$
is recovered. Figure 9(a) shows the comparison of the eigenvalues of CL between
$BS0$
and the other three actuation configurations with
$c=10^{2}$
. While
$BS0$
,
$BS1$
and
$BS3$
provide similar damped eigenvalues, the configuration
$BS2$
is the least effective actuator, which is further confirmed by the CL response of the full nonlinear system in figure 9(b,c). The figures also show that the fastest suppression of the vortex street is achieved by
$BS0$
, while the
$BS2$
actuator takes the longest time to eliminate vortex shedding. On the other hand, the actuation configurations of
$BS1$
and
$BS3$
behave similarly with respect to the suppression of the vortex street. The baseline
$BS0$
with
$\unicode[STIX]{x1D703}=45^{\circ }$
is considered to be the most effective actuator configuration. On comparing the
$BS1$
and
$BS2$
configurations, the midward
$BS1$
actuation slots have a better control performance than the leeward
$BS2$
configuration. This implies that the active suction/blowing control in the boundary layer (before the separation) is relatively efficient. Based on the above study, the configuration
$BS0$
is found to be more effective in reducing mean drag and suppressing fluctuating forces. Next, we investigate the effect of the suction/blowing angle
$\unicode[STIX]{x1D703}$
in the baseline
$BS0$
configuration.
An additional three (
$BS1$
,
$BS2$
,
$BS3$
) configurations of actuators as blowing and suction slots over the surface of the circular cylinder. The positive control input is defined as suction at the bottom and blowing at the top of the cylinder. While the
$BS1$
configuration forms a symmetric configuration of suction–blowing pairs,
$BS2$
and
$BS3$
are asymmetric with respect to the quadrants of the cylinder.
The effect of the placement of the fluidic actuator-based suction–blowing pairs on the cylinder surface: (a) distribution of the eigenspectrum of the CL system for different actuation slot placements (
$BS0{-}BS3$
); temporal variations of fluctuating (b)
$C_{l}$
and (c)
$C_{d}$
; (d) control input
$V_{c}$
for the CL response of the full nonlinear system at
$(Re,c)=(60,10^{2})$
.
3.2.2 Effect of angle for suction–blowing pairs
In this section, the sensitivity of the actuation angle for
$BS0$
is investigated by varying it from the baseline
$\unicode[STIX]{x1D703}=45^{\circ }$
to
$\unicode[STIX]{x1D703}=(30^{\circ },60^{\circ })$
. We keep the midward pair of suction/blowing actuation at
$\unicode[STIX]{x1D703}=(90^{\circ },270^{\circ })$
, but we change the other two diametrically opposed pairs of suction and blowing actuations. The distribution of the eigenspectrum, as shown in figure 10(a), suggests that the
$BS0$
configuration with
$\unicode[STIX]{x1D703}=(30^{\circ },45^{\circ },60^{\circ })$
provides a similar effectiveness with regard to the least damped eigenvalues. Next, the configuration
$BS0$
with different angle
$\unicode[STIX]{x1D703}$
is applied to stabilize the vortex shedding in a nonlinear saturated state. The results in figure 10(b,c,d) show similar suppression trends for the lift,
$C_{l}$
, and the drag,
$C_{d}$
, signals and the control input,
$V_{c}$
, obtained for the three angles. Figure 10(b,c) also shows that
$BS0$
with
$\unicode[STIX]{x1D703}=60^{\circ }$
has a slightly larger overshoot when the controller is switched on, followed by
$\unicode[STIX]{x1D703}=30^{\circ }$
and
$\unicode[STIX]{x1D703}=45^{\circ }$
. In all of our previous studies, the combined controller has been designed as the macro-manipulator for controlling the global dynamics. In other words, the actuators are not allowed to vary their control input speed
$V_{c}$
independently; therefore, the same control input
$V_{c}$
is generated in all suction–blowing pairs. It is interesting to study the effect of the control system architecture, whereby there is no coupling between the suction/blowing subsystems and the actuators can generate independent control inputs.
The effect of the angle of the fluidic actuators on the cylinder surface: (a) distribution of the eigenspectrum of the CL system for different actuation slot placements; time variations of fluctuating (b)
$C_{l}$
and (c)
$C_{d}$
; (d) control input
$V_{c}$
for the CL response of the full nonlinear system at
$(Re,c)=(60,10^{2})$
.
Schematics of combined and independent (decoupled) suction/blowing control system architectures. In contrast to the combined controller
$BS0$
(a),
$BS0D$
(b) partitions the controller pairs into different pieces with different blowing and suction velocities denoted by
$V_{c1}$
,
$V_{c2}$
and
$V_{c3}$
respectively.
Comparison between combined and decoupled control system architectures. (a) Distribution of the eigenspectrum:
$c=10^{2}$
for
$BS0$
and
$c=(10,10^{2},10^{3})$
for
$BS0D$
. Time variations of fluctuating (b)
$C_{l}$
and (c)
$C_{d}$
. (d) Control input
$V_{c}$
for the CL response of the full nonlinear system at
$c=10^{2}$
and
$Re=60$
. In (c), the vector (
$V_{c1},V_{c2},V_{c3}$
) is the control input of the
$BS0D$
configuration.
3.2.3 Combined versus independent controller
To understand the effect of the control architecture, we decouple the pairs of suction/blowing actuators and compare the performance against the combined controller counterpart. The decoupled controller configuration with different blowing and suction velocities, termed as
$BS0D$
, is shown in figure 11. As compared with the combined controller configuration
$BS0$
with 1-DOF control input velocity
$V_{c}$
, the controller of
$BS0D$
has three independent control inputs
$V_{c1}$
,
$V_{c2}$
and
$V_{c3}$
. As demonstrated in figure 12(a), the combined configuration
$BS0$
is more effective in damping the unstable eigenvalues than the decoupled
$BS0D$
counterpart. Moreover, the configuration
$BS0D$
is found to be less sensitive for the least damped eigenvalues for an identical range of
$c$
values. To further assess the performance of the
$BS0$
and
$BS0D$
controllers, figure 12(b,c) illustrates the nonlinear saturated state suppression through the force time histories. After removing the same matching control input constraint in the LQR algorithm for the
$BS0D$
controller, the force trends suggest that the controller becomes less effective when using different blowing and suction velocities at the same value of
$c$
. This is further confirmed in figure 12(d), which shows that the control velocity inputs (
$V_{c1}$
,
$V_{c2}$
,
$V_{c3}$
) become smaller at the same
$c=10^{2}$
, resulting in a larger time for suppression of the vortex street.
As shown in the previous section, the baseline
$BS0$
configuration with
$\unicode[STIX]{x1D703}=45^{\circ }$
is the most effective suction/blowing controller configuration. Therefore, we will employ this configuration for the active feedback control of VIV. As reported in Owen & Bearman (Reference Owen and Bearman2001), passive control techniques that work well in suppressing loads for fixed cylinders may not be effective for elastically mounted configurations. We next demonstrate our proposed active feedback control scheme for elastically mounted circular cylinders, which are free to vibrate in transverse only (1-DOF) and coupled streamwise/transverse directions (2-DOF).
The impulse response of a CL VIV system based on AFBS control: temporal evolution of (a) the lift coefficient
$C_{l}$
, (b) the transverse displacement
$Y$
and (c) the control input
$V_{c}$
predicted by the ROM and compared with the FOM for
$(Re,m^{\ast })=(60,10)$
,
$F_{s}=0.176$
and
$c=10^{2}$
with the controller switched on at
$tU/D=0$
.
3.3 Feedback control of VIV
In Yao & Jaiman (Reference Yao and Jaiman2017), we discussed how the onset of VIV lock-in is related to the instability exchange between the structural mode (SM) and the fluid mode (WM), where the SM and WM represent two distinct eigenvalue branches of the ERA-based VIV ROM system. The critical reduced natural frequency
$F_{s}$
or onset reduced velocity
$U_{r}=1/F_{s}$
can be pinpointed by the SM instability branch. The objective of the present active feedback control is to drive all unstable eigenvalues of the VIV system to the stable left half-complex-plane. While the OL VIV system (without feedback control) can be obtained by simply setting
$\unicode[STIX]{x1D66A}_{c}=0$
in (2.13), the optimal gain is computed for the VIV system at
$(Re,m^{\ast })=(60,10)$
using a similar procedure to that discussed for the stationary case.
As shown in figure 13, the CL VIV system response subject to an impulse is completely stabilized and the ROM results nearly overlap with the FOM counterparts. The linear stability is then determined by the eigenvalue analysis for both OL and CL systems, as illustrated in figure 14(a). Both the WM and part of the SM
$(0.147<F_{s}\leqslant 0.179)$
are unstable
$(Re(\unicode[STIX]{x1D706})>0)$
for the OL system, while the eigenvalues of the CL system are all in the stable left half-plane for
$0.005\leqslant F_{s}\leqslant 0.5$
(
$2\leqslant U_{r}\leqslant 200$
) with
$c=10^{2},10^{3}$
. The results show that the unstable eigenvalues are driven to the lower stable left half-plane as the coefficient
$c$
decreases, as indicated by the dash-dotted arrow in figure 14(a). It is important to note that the CL VIV system with
$c=10^{4}$
only delays the VIV onset and remains unstable for
$F_{s}\leqslant 0.117$
(
$U_{r}\geqslant 8.55$
), which indicates that the amplitude may grow continually as
$F_{s}$
decreases.
Root loci of the OL and CL VIV systems (
$0.005\leqslant F_{s}\leqslant 0.5$
): (a) sensitivity of the
$c$
parameter to the eigenspectrum at
$(Re,m^{\ast })=(60,10)$
; (b) effect of the mass ratio
$m^{\ast }$
on the control performance at
$(Re,c)=(60,10^{3})$
. The unstable right half-plane is shaded in grey. The parameter
$c\in [10^{2},10^{3}]$
ensures a VIV free system for the whole range of reduced frequency.
Suppression of the VIV response of a circular cylinder using the AFBS procedure: time variation of (a)
$C_{l}$
for saturated (thin dashed) and suppressed (thick solid) states and (b) transverse displacement
$Y$
and actuator input
$V_{c}$
with the controller switched on at
$tU/D=150$
,
$c=10^{3}$
and
$F_{s}=0.176$
; snapshots of spanwise vorticity contours at
$tU/D=$
(c) 155, (d) 400, (e) 600 and (f) 800. The contour levels are from
$-1$
to
$1$
in increments of
$0.1$
.
In figure 14(b), the root loci of different mass ratios are plotted, indicating that the unstable eigenvalues of lower mass ratio are damped more effectively with
$c=10^{3}$
. The controller with
$c=10^{3}$
is also able to suppress the saturated VIV response, as shown in figure 15. The lift coefficient
$C_{l}$
and the transverse displacement
$Y$
are suppressed once the controller is switched on after
$tU/D=150$
, when the VIV response reaches the saturated state. As compared with the saturated vortex shedding state of the stationary cylinder, which requires the maximum
$V_{c}\approx 0.65$
, the VIV case requires approximately
$4.5$
times larger maximum control input,
$V_{c}\approx 2.91$
with
$C_{\unicode[STIX]{x1D707}}\approx 4.40$
. Figure 15(c,d) shows snapshots of the spanwise vorticity contours to illustrate the attenuation of vortex shedding when the controller is switched on at
$t\geqslant 150$
. We can infer from the success of the proposed linear AFBS control for VIV that there exists a key linear mechanism during the self-sustained VIV oscillation that can be controlled effectively at low Reynolds number. Nonlinear effects of the fluid flow, which attempt to saturate vortex shedding and to form a limit cycle in VIV, become dominant in the later stage. It is worth pointing out that identification and control of the linear mechanism in nonlinear VIV is not equivalent to the prediction of nonlinear VIV with a linear model.
The OL saturated state of a 2-DOF VIV system at
$(Re,m^{\ast })=(60,10)$
: (a) figure-of-eight trajectory of the
$X$
and
$Y$
displacements; (b) the normalized power spectrum
$P$
versus the frequency
$f^{\ast }$
of the lift coefficient
$C_{l}$
and the transverse displacement
$Y$
, where
$f^{\ast }=f/F_{s}$
. Suppression of the 2-DOF VIV response: (c) transverse
$Y$
and inline
$X$
response trends; (d) lift coefficient
$C_{l}$
and actuator input
$V_{c}$
. The controller is switched on after
$tU/D=150$
with
$c=10^{3}$
and
$F_{s}=0.176$
.
Feedback control of the 2-DOF VIV of a circular cylinder at
$(Re,m^{\ast })=(60,10)$
. Snapshots of the spanwise vorticity contours at
$tU/D=$
(a) 175, (b) 200, (c) 400 and (d) 500. The contour levels are from
$-1$
to
$1$
in increments of
$0.1$
.
Comparison of the control inputs between a stationary cylinder, a 1-DOF VIV system and a 2-DOF VIV with
$BS0$
at
$c=10^{3}$
. The root mean square (r.m.s.) and maximum (
$\max$
) are computed for
$t\in [t_{0},t_{0}+200]$
;
$t_{0}$
is when the controller is switched on after
$t_{0}=175.5$
and
$150$
for the stationary and VIV systems respectively.
To demonstrate whether the designed controller for the 1-DOF VIV system can be utilized for 2-DOF VIV system feedback control, we introduce free vibration in the streamwise (
$X$
) direction along with the transverse (
$Y$
) direction. A typical figure-of-eight VIV saturated trajectory is shown in figure 16(a) at
$F_{s}=0.176$
. The frequency plot shown in figure 16(b) confirms the occurrence of frequency lock-in. As shown in figure 16(c,d), similar suppression is achieved for the 2-DOF VIV system with the same controller as designed for the 1-DOF VIV system (figure 15). As compared with 1-DOF, which requires maximum
$V_{c}\approx 2.91$
, the same controller requires maximum
$V_{c}\approx 3.0$
for 2-DOF VIV saturated response suppression, or approximately
$3\,\%$
larger than its 1-DOF VIV counterpart. Figure 17 further illustrates the similar flow-field evolution during the suppression process for the 2-DOF system. Overall, the results suggest that the proposed feedback control scheme based on vertical blowing/suction is able to suppress the 2-DOF VIV system effectively with approximately the same amount of control energy as the 1-DOF system.
The control input is summarized in table 1 for a stationary cylinder and VIV configurations with both 1-DOF and 2-DOF response at
$c=10^{3}$
. As shown in the table, the 2-DOF VIV system only requires a slightly larger control energy to suppress the saturated VIV state than its 1-DOF counterpart. The results suggest that the transverse response is the most dynamically significant for the isolated circular cylinder VIV system. Therefore, the proposed controller designed for a 1-DOF VIV system can be applied to reduce 2-DOF VIV response.
4 Concluding remarks
Through ERA-based dimensionality reduction, an AFBS concept is proposed to suppress the vortex street and VIV for flexibly mounted structures. A variational finite-element formulation has been considered for the full-order fluid–structure model and the generalized eigenvalue problem of the linearized NS system. The control scheme relies on the optimal feedback gain by the LQR synthesis and state estimation by the Kalman filter. Based on combined vertical blowing and suction, we employed the feedback control to obtain suitable gains that stabilize the wake instability and the VIV. We found that the essential elements in the self-sustaining VIV process are linear, and are subject to the active feedback control. From the ERA-based computations and feedback control of unstable eigenmodes, it can be deduced that the increase in the VIV amplitude of a cylinder occurs primarily through the linear instability process.
By applying the AFBS procedure to the vortex shedding of a stationary cylinder, we observed a remarkable reduction in the time-dependent fluctuating components of both the lift and drag forces. By means of eigenspectrum distributions via ERA-based simulations, we have explored several configurations to confirm the generality and sensitivity of the AFBS procedure with respect to a range of parameters. While the configuration
$BS0$
with six actuation slots performed efficiently, followed by
$BS1$
and
$BS3$
, the configuration
$BS2$
exhibited poor performance as a result of having the actuators on the leeward side of cylinder. In general, the performance of the
$BS0$
configuration is not very sensitive to the jet angle
$\unicode[STIX]{x1D703}$
. When it is applied for a nonlinear saturated state, similar suppression is achieved with
$\unicode[STIX]{x1D703}=(30^{\circ },45^{\circ },60^{\circ })$
, while a slightly larger overshoot is found when the controller is switched on for the angle
$\unicode[STIX]{x1D703}=60^{\circ }$
. The performance of
$BS0$
becomes less effective when the same control input velocity constraint is removed, thereby suggesting an improved performance of the combined control architecture compared with the decoupled independently designed control system.
With respect to VIV suppression, the designed controller performed well for a range of the mass ratio
$m^{\ast }\in [5,100]$
and the reduced natural structural frequency
$F_{s}\in [0.005,0.5]$
at
$Re=60$
. In contrast to stationary vortex shedding, suppression of VIV requires approximately four times larger control input for the same Reynolds number
$Re=60$
to suppress both the fluctuating lift and the VIV amplitude. The controller designed for a transverse VIV system can also be adapted to 2-DOF VIV suppression with only approximately
$3\,\%$
larger control input than the transverse 1-DOF VIV system. Since the present ROM does not require any adjoint solver, the proposed ERA-based feedback control can be directly used for actual physical problems and experimental settings. The ERA-based ROM approach is relatively straightforward and computationally efficient. Furthermore, the model offers a reasonable accuracy for the development of stabilized feedback controller design for unstable flows and VIV. Three-dimensional and high-
$Re$
flows will be of interest for future study.
Acknowledgement
The first author would like to thank Singapore Maritime Institute Grant (SMI-2014-OF-04) for the financial support.
