2.1 Problem set-up

The two-dimensional computational model in the current study employs a rigid NACA0015 airfoil with a slightly rounded trailing edge, immersed in an incompressible fluid with free-stream velocity $U_{\infty }$ . The slight rounding of the trailing edge ensures that the flow is well resolved around it, and we have verified that this has no significant effect on the aerodynamic characteristics and forces on the airfoil. The governing incompressible Navier–Stokes equations expressed in the dimensionless form as follows:

(2.1) $$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\boldsymbol{u}}{\unicode[STIX]{x2202}t}+\boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{u}=-\unicode[STIX]{x1D735}p+\frac{1}{Re}\unicode[STIX]{x1D6FB}^{2}u,\end{eqnarray}$$

where $Re=\unicode[STIX]{x1D70C}U_{\infty }C/\unicode[STIX]{x1D707}$ is the chord-based Reynolds number. The torsional elasticity of the airfoil is modelled using a linear torsional spring with spring constant $k$ and damping coefficient $b$ , attached to the rigid airfoil at a prescribed chordwise location ( $X_{e}$ ), which we refer to as the elastic axis (normalized location denoted as $X_{e}^{\ast }=X_{e}/C$ ). The equilibrium angular position of the spring is denoted by angle $\unicode[STIX]{x1D703}_{0}$ . A schematic of this set-up is shown in figure 1.

Figure 1.

(a) Schematic of the aeroelastic system used in this study; (b) computational domain and close-up of the Cartesian computational grid.

The structural system is governed by the equation for a forced spring–mass–damper oscillator. The equation is scaled by the characteristic variables used for the flow (length: $C$ ; time: $c/U_{\infty }$ ), and the resulting non-dimensionalized equation governing the dynamics in pitch of the airfoil is given by

(2.2) $$\begin{eqnarray}I^{\ast }\ddot{\unicode[STIX]{x1D703}}+b^{\ast }\dot{\unicode[STIX]{x1D703}}+k^{\ast }(\unicode[STIX]{x1D703}-\unicode[STIX]{x1D703}_{0})=C_{M},\end{eqnarray}$$

where $C_{M}=M/(0.5\unicode[STIX]{x1D70C}U_{\infty }^{2}C^{3})$ is the coefficient of aerodynamic pitching moment, and $I^{\ast }=2I/(\unicode[STIX]{x1D70C}C^{4})$ , $b^{\ast }=2b/(\unicode[STIX]{x1D70C}U_{\infty }C^{3})$ , $k^{\ast }=2k/(\unicode[STIX]{x1D70C}U_{\infty }^{2}C^{2})$ are the dimensionless moment of inertia, damping and spring stiffness respectively. Here, $C_{M}$ and $I^{\ast }$ are calculated with respect to $X_{e}^{\ast }$ . The spring stiffness is commonly expressed in terms of a reduced velocity, $U^{\ast }=U_{\infty }/f_{s}C$ where $f_{s}=(1/2\unicode[STIX]{x03C0})\sqrt{k/I}$ is the natural frequency of the spring. In essence, $U^{\ast }$ is a ratio of the natural time scale of the elastic system to the convective time scale of the flow. However, as we will show, the choice of fluid time scale has significant implications for the dynamics of the system.

Thus, even this highly simplified configuration has six governing parameters ( $Re$ , $X_{e}^{\ast }$ , $I^{\ast }$ , $b^{\ast }$ , $U^{\ast }$ , $\unicode[STIX]{x1D703}_{0}$ ) and this is indicative of the inherent complexity of this configuration. In the current study, we fix the values of $Re$ and $I^{\ast }$ to 1000 and 4.1 respectively and explore the effect of the other four parameters. The Reynolds number chosen is high enough so as to generate robust vortex shedding phenomenon, which drives flutter, but is also low enough so as to allow resolved simulations at a reasonable computational expense. The choice of the moment of inertia about the elastic axis corresponds to a solid-to-fluid density ratio of ${\approx}120$ . While the effect of varying $I^{\ast }$ in this system could be of significant engineering interest, we choose to use a fixed value throughout this study to reduce the number of independent parameters. However, it must be noted that varying $I^{\ast }$ is equivalent to changing the natural frequency as $f_{s}\sim 1/\sqrt{I^{\ast }}$ , which is an important parameter in this study.

The current study addresses structurally damped as well as undamped systems with damping being defined by $\unicode[STIX]{x1D701}^{\ast }=b/b_{cr}$ where $b_{cr}=2\sqrt{kI}$ , is the critical damping for the harmonic oscillator. For the undamped system $\unicode[STIX]{x1D701}^{\ast }=0$ and for the damped system we prescribe $\unicode[STIX]{x1D701}^{\ast }=0.15$ where this lightly damped structure enables us to understand the effect of damping without overly damping the pitch oscillations. The remaining two parameters, the elastic axis ( $X_{e}^{\ast }$ ) and the spring stiffness ( $U^{\ast }$ ), are a particular focus of the current study. For $X_{e}^{\ast }$ , we use three values – $X_{e}^{\ast }=0.50$ , 0.33 and 0.25 – and we will show that these generate very significant variations in the pitch response. The final parameter $U^{\ast }$ represents the spring stiffness (a larger $U^{\ast }$ corresponds to a softer spring) and this parameter is varied from about 2 to 13.

2.2 Numerical method

We simulate the coupled fluid–structure system using the sharp-interface immersed boundary method based solver ‘ViCar3D’ described in Mittal et al. (Reference Mittal, Dong, Bozkurttas, Najjar, Vargas and von Loebbecke2008) and Seo & Mittal (Reference Seo and Mittal2011). This allows us to preserve sharp interfaces along the surface of our geometry using a body-non-conformal Cartesian grid. This is particularly useful in fluid–structure interaction studies as it allows us to simulate a variety of shapes as well as large-amplitude motions on simple, non-adaptive Cartesian grids. In particular, this approach has the advantage of not being constrained by highly deformed grid cells or the need to periodically recompute the grid during high-amplitude motions (Mittal & Iaccarino Reference Mittal and Iaccarino2005). The fluid equations are solved using a second-order fractional-step method, and the pressure Poisson equation is solved using the geometric multigrid method. Second-order finite differences are used for all spatial derivatives, and the time integration is performed using a second-order Adams–Bashforth method. ViCar3D employs a boundary imposition scheme that is second-order accurate and consequently, models surface quantities such as pressure and shear with high fidelity (Mittal et al. Reference Mittal, Dong, Bozkurttas, Najjar, Vargas and von Loebbecke2008). ViCar3D has been validated and verified extensively for a variety of stationary and moving boundary problems in earlier papers (Ghias, Mittal & Dong Reference Ghias, Mittal and Dong2007; Mittal et al. Reference Mittal, Dong, Bozkurttas, Najjar, Vargas and von Loebbecke2008; Seo & Mittal Reference Seo and Mittal2011), as well as during the current study (see appendix B).

The fluid–structure coupling employs a loosely coupled approach wherein the flow and dynamical equation are solved sequentially. The aerodynamic forces/moments are calculated on the Lagrangian marker points on the surface of the solid body, and passed on to the solid dynamical equation (2.2). The angular velocity predicted from this equation is subsequently imposed on the marker points. The aeroelastic system is immersed in a large $21C\times 23C$ computational domain, where the airfoil is placed 10 chord-lengths from the downstream boundary, and the isotropic grid resolution around the solid body corresponds to approximately 125 points along the chord. The grid is stretched in all directions away from the rectangular region that surrounds the foil and the near wake resulting in a baseline grid of $384\times 320$ points. A Dirichlet velocity boundary condition is used at the inlet boundary of the domain, and zero-gradient Neumann conditions are specified at all other boundaries. Grid refinement studies described in appendix B confirm that the results on this grid are well converged.

For all the cases simulated here, we initialize the airfoil at its equilibrium angle $\unicode[STIX]{x1D703}_{0}$ , with zero initial angular velocity. The constant free-stream flow is then imposed and the dynamics allowed to evolve until the system reaches a stationary state. For the cases simulated here, it takes between $O(20)$ and $O(200)$ oscillation cycles before this stationary state is achieved. Average quantities are computed after the stationary state is reached.

4.1 Energy maps

Our analysis in the preceding sections indicates a need to better understand the amplitude response of the flow-induced pitching airfoil for a set of given system parameters. We were able to demonstrate the mechanism governing the onset and frequency of oscillation as a competition of time scales. However, the variations in amplitude with $U^{\ast }$ and $X_{e}^{\ast }$ are complex and non-intuitive. Here, we propose the use of ‘energy maps’ to gain an understanding of the complex amplitude response of the system. For an airfoil undergoing sinusoidal pitching oscillations in a free stream, we can derive (see appendix A) the following equation for the growth of the amplitude over one cycle

(4.1) $$\begin{eqnarray}A_{\unicode[STIX]{x1D703}}(t_{n}+T)-A_{\unicode[STIX]{x1D703}}(t_{n})=\frac{E^{\ast }}{k^{\ast }\bar{A}_{\unicode[STIX]{x1D703}}},\end{eqnarray}$$

where $\bar{A}_{\unicode[STIX]{x1D703}}$ is the average amplitude over a cycle. In the above equation $E^{\ast }=E_{f}^{\ast }-E_{d}^{\ast }$ , where $E_{f}^{\ast }$ is the energy extracted by the airfoil from the free stream over one oscillation cycle and $E_{d}^{\ast }$ is the energy lost to structural damping; $E_{f}^{\ast }$ is given by

(4.2) $$\begin{eqnarray}E_{f}^{\ast }=\int _{t_{n}}^{t_{n}+T}C_{M}\dot{\unicode[STIX]{x1D703}}\,\text{d}t,\end{eqnarray}$$

where $\dot{\unicode[STIX]{x1D703}}$ is the angular velocity and $T$ is the period of oscillation.

In the case of purely sinusoidal oscillations where $\dot{\unicode[STIX]{x1D703}}=\dot{\unicode[STIX]{x1D6E9}}\sin (2\unicode[STIX]{x03C0}t/T)$ and $C_{M}=K\sin (2\unicode[STIX]{x03C0}t/T+\unicode[STIX]{x1D719})$ , the energy transfer is given by $E_{f}^{\ast }=(1/2)\dot{\unicode[STIX]{x1D6E9}}KT\cos \unicode[STIX]{x1D719}$ and this simple model shows that the sign of energy extraction, which determines whether the airfoil extracts energy from the flow, is governed by the phase difference between the pitch velocity and the pitch moment. The physical interpretation of this phase difference, in the context of flow-induced oscillations, is related to the timing between the pitching and the resulting response of the flow. In § 3.3, our description of the mechanism determining the frequency response showed that this timing determines the critical $U^{\ast }$ for the bifurcation to large amplitude oscillations. The phase difference between pitch and pitch moment is also a strong factor when considering the changes in the location of the elastic axis. This is demonstrated in figure 12 where we show time series plots of the pitch deflection and moment coefficient for an airfoil forced to oscillate at $A_{\unicode[STIX]{x1D703}}=25^{\circ }$ , about two different locations of the elastic axis. These sample cases show that changing the elastic axis, while retaining the kinematics, produces very different forcing on the airfoil. This is a result of the timing associated with the shedding as well as convection of vortices past the hinge location. These effects were also demonstrated in figure 1. The difference in timing leads to the case with $X_{e}^{\ast }=0.33$ being damped, while $X_{e}^{\ast }=0.50$ extracts energy. These examples suggest that the phase difference, and hence energy extraction, is indeed an important factor in analysing the response of flow-induced flutter.

Figure 12.

Comparison of $\unicode[STIX]{x0394}\unicode[STIX]{x1D703}$ , $C_{M}$ and extracted energy $E_{f}^{\ast }$ , for airfoils forced to oscillate with $A_{\unicode[STIX]{x1D703}}=25^{\circ }$ and $f_{p}^{\ast }=0.25$ about two different hinge locations, (a) $X_{e}^{\ast }=0.33$ ; (b) $X_{e}^{\ast }=0.50$ .

The aim in this section is to use our knowledge of the energy extraction, as described above, to understand the amplitude response of the dynamical system as a function of the various parameters explored in this study. To do this, we first create a map of energy transfer as a function of the parameters of interest. This is done by performing simulations of airfoils that are forced to undergo sinusoidal pitching oscillations over a range of amplitudes and frequencies, while holding the other parameters such as $X_{e}^{\ast }$ and $Re$ fixed. The pitching and pitch-moment data from each of these simulations is then used to calculate the energy extracted by the airfoil. The validity of these energy maps, which are obtained from forced oscillations, for the case of free, flow-induced oscillations rests primarily on the assumption that the free oscillations are strictly sinusoidal. We find that this condition is indeed satisfied to a very large degree for all the cases simulated here (see figures 2 and 5) and a sinusoidal fit to the pitching time series yields an R-square value of roughly 0.99 for nearly all cases studied. Indeed, energy maps have been used quite successfully to understand flow-induced vibration of bluff bodies (Morse & Williamson Reference Morse and Williamson2009) and Bhat & Govardhan (Reference Bhat and Govardhan2013) also demonstrated the use of an energy map to study flow-induced flutter of an airfoil.

Figure 13.

(a) Contours of energy extraction, $E_{f}^{\ast }$ , as a function of oscillation amplitude and frequency for $X_{e}^{\ast }=0.50$ , $\unicode[STIX]{x1D703}_{0}=15^{\circ }$ . Symbols represent stationary state $f_{p}^{\ast }$ and $A_{\unicode[STIX]{x1D703}}$ for undamped (○) and damped (♢) flow-induced oscillation cases. Dashed and dotted lines represent undamped ( $E_{f}^{\ast }=0$ ) and damped ( $E_{f}^{\ast }=E_{d}^{\ast }$ ) equilibrium curves respectively. $E_{d}^{\ast }$ is calculated using $\unicode[STIX]{x1D701}=0.15$ for the darkest damped equilibrium curve, and increasing values of $\unicode[STIX]{x1D701}$ for subsequently lighter equilibrium curves ( $\unicode[STIX]{x1D701}=0.30$ , 0.45); (b) reproduction of figure 10(a) with inset box showing range of $U^{\ast }$ and $A_{\unicode[STIX]{x1D703}}$ computed in the energy map; (c) zoom-in of flow-induced oscillation amplitudes for cases compared with the energy map.

The first map shown in figure 13 is for the case with $X_{e}^{\ast }=0.50$ and $Re=1000$ . This map is the result of a total of ${\sim}300$ simulations over a range of non-dimensional pitch frequencies $f_{p}^{\ast }=f_{p}C/U_{\infty }:[0.10,0.75]$ and pitch amplitudes $A_{\unicode[STIX]{x1D703}}:[0,50^{\circ }]$ . For each case, we simulate over 50 total cycles to ensure that a stationary state is achieved, and obtain $E_{f}^{\ast }$ by integrating over the last 20 cycles. The dashed line in the energy map corresponds to $E_{f}^{\ast }=0$ and the dotted line corresponds to $E_{f}^{\ast }=E_{d}^{\ast }$ , where $E_{d}^{\ast }$ is calculated using the value of structural damping previously specified in § 2.1. These represent potential equilibrium conditions in the case of flow-induced flutter.

The primary features in the map are the nearly vertical demarcation between regions of positive and negative energy transfer at approximately $f_{p}^{\ast }\approx 0.25$ and the small region of positive energy transfer around $f_{p}^{\ast }\approx 0.65$ . From this energy map, we expect that sustained flow-induced oscillations for cases with structural damping would occur in the region with $E_{f}^{\ast }>0$ , whereas cases without structural damping should saturate along the $E_{f}^{\ast }=0$ contour line. To verify this, we superimpose on this map, the stationary-state amplitude and frequencies from the flow-induced flutter simulations for $X_{e}^{\ast }=0.50$ . We find that for all these cases, the system saturates very close to the contour line corresponding to net energy equilibrium, i.e. $E_{f}^{\ast }=0$ for undamped oscillations and $E_{f}^{\ast }=E_{d}^{\ast }$ for damped oscillations. These observations provide clear confirmation of the validity of these energy maps as a tool for understanding the flow-induced flutter response.

Examining the energy maps further, we note that the nearly vertical $E_{f}^{\ast }=0$ contour line lies very close to the frequency corresponding to the reciprocal of the critical velocity, which is the frequency at which the bifurcation to large-amplitude oscillations occurs for the flow-induced oscillation with $X_{e}^{\ast }=0.50$ . This feature of the current energy map agrees qualitatively with that generated via experiments of Bhat & Govardhan (Reference Bhat and Govardhan2013), despite the fact that their Reynolds number is an order of magnitude larger than the current one and their pitch axis is at quarter chord. Also, the small region of positive energy transfer on the right side of the map has a frequency that corresponds to the Kármán oscillation mode identified earlier (see figure 4).

The topology of the energy map and the equilibrium curve provide useful insights regarding the dynamics of this configuration. For instance, the nearly vertical demarcation between the regions of energy growth and decay at $f_{p}^{\ast }\approx 0.25$ suggests that an undamped system would experience rapid growth in flutter amplitude with a relatively small decrease (increase) in pitch frequency (reduced velocity). Indeed, figure 10(a) shows that the undamped system exhibits what seems to be a subcritical bifurcation at $U^{\ast }\approx 4.0$ with a jump in amplitude of nearly $40^{\circ }$ . As the structural damping is increased, the equilibrium curve bends to the left (i.e. increasingly negative slope) indicating a smoother onset of flutter. Furthermore, increased damping significantly limits the maximum amplitude attained by the system. Other implications of the topology of the energy map for the system dynamics are explored in the next section.

4.2 Amplitude response for $X_{e}^{\ast }=0.33$

The case with $X_{e}^{\ast }=0.33$ exhibits a very complex, non-monotonic response with $U^{\ast }$ and we now demonstrate the utility of the energy maps in providing insight into this response. As shown in (4.1) and appendix A, oscillations in the absence of structural damping are expected to grow (decay) in regions of positive (negative) $E_{f}^{\ast }$ . The requirements for a stable equilibrium are

(4.3a,b ) $$\begin{eqnarray}E_{f}^{\ast }=0;\quad \text{d}E_{f}^{\ast }/\text{d}A<0.\end{eqnarray}$$

Figure 14 shows the energy map for this case which is generated as before by conducting ${\sim}200$ simulations of prescribed flutter over a large range of pitch frequencies and amplitudes. On this energy map, we have identified the curves corresponding to the equilibrium condition ( $E_{f}^{\ast }=0$ ) as well as the regions of stability on these equilibrium curves ( $\text{d}E_{f}^{\ast }/\text{d}A<0$ ). A number of qualitative observations can immediately be made regarding this energy map. First, the topology is very different from that for $X_{e}^{\ast }=0.50$ . Since free oscillations can only exist within the region of $E_{f}^{\ast }\geqslant 0$ , the energy map suggests an upper limit of $A_{\unicode[STIX]{x1D703}}\approx 40^{\circ }$ on the flow-induced oscillation amplitude in this case. The topology of the map is also extremely complicated: two main regions of energy growth and decay are separated by a complex equilibrium curve that extends over the entire range of frequencies considered here. There also exist multiple ‘islands’ with $E_{f}^{\ast }<0$ within the large region of energy growth in the low-amplitude ( ${<}10^{\circ }$ ) portion of the map, and this points to the possibility of multiple equilibrium states for a given conditions as well as abrupt changes in oscillation amplitude with changes in the underlying parameters.

Figure 14.

(a) Energy map for $X_{e}^{\ast }=0.33$ and $\unicode[STIX]{x1D703}_{0}=15^{\circ }$ , with the equilibrium curve ( $E_{f}^{\ast }=0$ ) shown as a dashed line, and stable equilibrium ( $E_{f}^{\ast }=0$ , $\text{d}E^{\ast }/\text{d}A<0$ ) highlighted along the equilibrium curve (thick line); (b) maximum pitch deflections versus $U^{\ast }$ for cases with $\unicode[STIX]{x1D701}=0$ corresponding to figure 10(b); the stationary-state oscillation amplitude and frequency of all cases shown in (b) are plotted on the energy map using circles (○). Also shown on the energy map are bars indicating the local resolution in $f_{p}^{\ast }$ and $A_{\unicode[STIX]{x1D703}}$ in various regions of the energy map. Local regions of the map that contain complex topological features have been provided with higher resolution.

To compare our prediction for the equilibrium positions on this energy map, we superimpose the stationary-state locations of our flow-induced data for $X_{e}^{\ast }=0.33$ on the energy map using circle (○) symbols. This data set was described earlier in § 3.5 (see figure 10 b) and is also shown as maximum pitch deflection versus $U^{\ast }$ in figure 14(b). These cases correspond to zero structural damping, and are hence expected to lie along stable equilibria identified on the $E_{f}^{\ast }=0$ curve. In figure 14 we see that there is again very good agreement between the flow-induced oscillations data and the energy map. The slight discrepancies in this comparison are a result of the finite resolution of the energy map, which is not able to fully resolve the complicated structure of the map. This is especially apparent at very small oscillation amplitudes, where multiple equilibria exist in close proximity to each other. Hence on figure 14 we have also plotted error bars at various locations indicating the local step size in $f_{p}^{\ast }$ and $A_{\unicode[STIX]{x1D703}}$ used in computing generating the energy map.

In spite of the difficulty in fully resolving the topology of the energy map, it allows us to explain the complex, non-monotonic response of the flow-induced flutter for the undamped $X_{e}^{\ast }=0.33$ case. The map clearly shows that for $f_{p}^{\ast }>0.25$ $(U^{\ast }<4)$ the system is trapped in a lower branch of the energy map which has very small amplitude. However, as $f_{p}^{\ast }$ reduces below $0.25$ $(U^{\ast }>4)$ the lower branch all but disappears. Thus for the case with $f_{p}^{\ast }>0.237\;(U^{\ast }=4.22)$ the system finds itself in a region of positive energy growth and the amplitude grows and saturates at the first stable equilibrium along the vertical line representing its natural frequency in the map, which is approximately $A_{\unicode[STIX]{x1D703}}=22^{\circ }$ . Between $f_{p}^{\ast }=0.22$ and 0.25, the system continues to grow to the large amplitudes associated with the upper equilibrium curve, however, for $f_{p}^{\ast }<0.22\;(U^{\ast }>4.5)$ , a lower, stable equilibrium branch appears and the system locks on this lower branch. Hence, after the bifurcation to large-amplitude oscillations at approximately $U^{\ast }=4$ , the amplitude of the system drops significantly ( ${\approx}3^{\circ }$ for $U^{\ast }>4.54$ . In contrast, for the case of $X_{e}^{\ast }=0.50$ , there is no such lower branch in the energy map (see figure 13). Hence the oscillation amplitude increases monotonically along the $E_{f}^{\ast }=0$ curve for that case.

Figure 15.

Prediction of $A_{\unicode[STIX]{x1D703}}$ using (A 5), compared with calculations from flow-induced oscillation simulations. Here, $X_{e}^{\ast }=0.33$ and $\unicode[STIX]{x1D703}_{0}=15^{\circ }$ ; (a) stationary-state amplitudes and frequencies from the model (▫), and flow-induced simulations corresponding to figure 14 (○) plotted on the energy map; (b) corresponding amplitudes plotted against $U^{\ast }$ using the same symbols as in (a).

Having established that the energy map indeed describes the stationary states of the flow-induced flutter system, we now demonstrate its use in making a priori predictions of the stationary-state amplitude without the need for flow simulations. This serves as a demonstration of the energy map as a predictive tool as well as a validation of our model for amplitude growth (shown in appendix A). From (A 5), we can calculate the amplitude growth of an oscillator at every cycle, with the knowledge of its spring stiffness and cycle-wise energy extraction. As mentioned before, this assumes that the oscillations are purely sinusoidal and occur at the natural frequency of the system. With these assumptions, we can pick the oscillation frequency of the system and predict the growth of amplitude based on the energy map and an initial condition. We iterate (A 5) using a semi-implicit Crank–Nicholson scheme, and the iteration is performed until a stationary state is achieved. We demonstrate this for some selected cases in figure 15, where we compare the amplitude response observed in flow-induced oscillations (○), with that predicted by our model (▫). We see that we are able to predict the amplitude response accurately for most cases by purely using the energy map. One case however, with $U^{\ast }=4.62$ , does not agree with the model’s predictions. As is apparent in figure 15(a), this is because the energy map computed here predicts the existence of a small-amplitude stable equilibrium at this oscillation frequency. This is a result of finite resolution in computing the map, as was mentioned earlier. Furthermore, the fact that this region of the map is especially sensitive to small changes in frequency and amplitude magnifies the finite resolution effects. This sensitivity of the energy prediction to input conditions was also observed by Leontini et al. (Reference Leontini, Stewart, Thompson and Hourigan2006).

4.3 Transient response

The energy maps also allow us to understand the transient response of the various flow-induced flutter cases. In figure 16 we plot the pitching amplitude time series and trajectories of two cases of flow-induced flutter, where the circles represent the location of the oscillator on the frequency–amplitude space during each cycle. This is calculated by dividing the pitch response time series into individual cycles, and calculating the amplitude and time scale associated with each cycle. In the case of $U^{\ast }=4.50$ , we observe that the trajectory of the system follows a nearly vertical line on the energy map. This is because this case corresponds to the post-bifurcation regime, and the frequency is selected based on the natural frequency of the system. The system starts at rest, close to $\unicode[STIX]{x0394}\unicode[STIX]{x1D703}=0$ , and the amplitude grows due to the positive energy extracted from the flow. We see in figure 16(c) that this amplitude growth occurs until it reaches a stable equilibrium position corresponding to its oscillation frequency. However, before reaching this stationary state, the oscillator spends a long time close to a small-amplitude equilibrium region, as shown by the density of circles in the zoomed-in trajectory shown in figure 16(d). The effect of this is also seen in the pitch amplitude time series of this case, shown in figure 16(a). The system exhibits small amplitude oscillations for an extended period of time, before a sudden growth in amplitude caused by its escape from the small-amplitude stable region.

Figure 16.

Comparison of two flow-induced oscillation cases with $U^{\ast }=4.50$ and $U^{\ast }=4.95$ . Top panels show time series of pitch oscillations for (a) $U^{\ast }=4.50$ and (b) $U^{\ast }=4.95$ ; (c) amplitude–frequency trajectories of these cases on the energy map. Each circle represents the amplitude and frequency during one oscillation cycle. The colour intensity represents time; (d) zoom-in of trajectories in the inset box shown in (c).

On slightly increasing $U^{\ast }$ to 4.95, the system gets trapped in the low-amplitude stable branch, as seen in figure 16(c). Further, we see from the pitch response time series in figure 16(b) that the oscillation amplitude in this case shows an interesting beating behaviour. This suggests that the oscillator is trapped in a bi-stable region that allows it to move between two stable equilibrium positions. While the exact structure of this stable region is difficult to resolve, figure 16(d) shows that the system gets trapped in a nearly vertical branch of $E_{f}^{\ast }=0$ that would allow the existence of multiple stable regions. Further, the vertical spread of instantaneous amplitudes after the system has reached a stationary state (shown by the spread of circles) suggests that the system does indeed move between two equilibrium states. The analysis of these cases shows the utility of the energy map in analysing the stationary-state response, as well as the trajectory of the oscillator from the initial condition to the stationary state. We also see that the complicated structure of the energy map leads to interesting behaviour in the pitching response of this system, which would be very difficult to interpret without knowledge of the corresponding energy map.

Figure 17.

(a) Energy map for $X_{e}^{\ast }=0.33$ , with equilibrium curve, stable equilibria and stationary-state frequency and amplitudes as in figure 14. Also shown are cases demonstrating hysteresis, represented by diamonds (♢), by continuing along the directions indicated by the arrows; (b) the stationary-state oscillation amplitude of all cases shown in the energy map, plotted against $U^{\ast }$ using the same symbols.

Figure 18.

Trajectories of flow-induced oscillation cases showing hysteresis. Top panels show time series of pitch oscillations for (a) $U^{\ast }=3.30$ and (b) $U^{\ast }=4.95$ , where the lighter time series represents the approach to the initial state (before hysteresis) and the darker time series shows the evolution after the change of $U^{\ast }$ ; (c) amplitude–frequency trajectories of the case with $U^{\ast }=3.30$ on the energy map, where the approach to the initial condition as well as the subsequent evolution at $U^{\ast }=3.30$ is shown; (d) amplitude–frequency trajectories of the case with $U^{\ast }=4.95$ .

The energy map also allows us to predict hysteretic behaviour due to the presence of multiple stable equilibria at particular oscillation frequencies. In figure 17 we show the energy map, as in figure 14, with additional cases that demonstrate hysteresis shown using diamond (♢) symbols. These simulations were carried out using a previously computed stationary state along the equilibrium curve as an initial condition, changing the value of $U^{\ast }$ , and allowing the dynamics to evolve to a new stationary state. We use arrows pointing from the initial to the final stationary state to denote the change in $U^{\ast }$ in figure 17. We see that these cases agree with the equilibrium and stability conditions predicted by the energy map too.

The transient response of these cases is especially interesting to analyse as it also demonstrates the ability to use these energy maps in flutter amplitude control. In figure 18 we plot the pitching time series and trajectories of two such cases showing hysteretic behaviour with $U^{\ast }=3.30$ and $U^{\ast }=4.95$ . Figures 18(a) and 18(b) show time series plots of pitching amplitude, where the grey portion represents the last few cycles of the case used as an initial condition, and the black time series shows the subsequent amplitude response leading to the new stationary state. In figures 18(c) and 18(d) we show the trajectories of these cases where the circles show the state of the system on the frequency–amplitude space during each oscillation cycle. Here we show the trajectory leading up to the stationary state that is used as an initial condition, as well as the trajectory from the previous stationary state to the new one. The change in $U^{\ast }$ is indicated by the arrow.

For $U^{\ast }=3.30$ , the case using the static initial condition shows very small-amplitude oscillations ( $A_{\unicode[STIX]{x1D703}}<1^{\circ }$ ) as seen in figure 10(b) and is described by the pre-bifurcation regime on § 3.3. However, when initialized with an oscillation amplitude ${\sim}25^{\circ }$ , we see that it reaches a stationary-state amplitude of approximately $10^{\circ }$ . This new stationary state is achieved through a loss of energy, which is different from all the cases analysed in this paper so far. The loss of energy is a consequence of the fact that the initial amplitude of oscillation, which corresponds to a stable equilibrium at a different natural frequency, is in a region of negative energy transfer at the natural frequency corresponding to $U^{\ast }=3.30$ . This results in a reduction in amplitude of oscillation as described by (A 6), until it reaches a new stable equilibrium. This is an example of the ability to reduce the amplitude of flutter using a knowledge of the energy map.

The case of $U^{\ast }=4.95$ shows a case that is initialized using a stationary state on one portion of the equilibrium curve, and ends up on higher-amplitude portion of the curve on changing its natural frequency. This occurs as the change in natural frequency pushes the system into a region of positive energy transfer, as indicated by the arrow in figure 18(d). Due to continuous extraction of energy from the flow, the amplitude grows until it reaches a new stable equilibrium. It must be noted that this is the same $U^{\ast }$ value shown in figures 16(b) and 16(d), where the static initial condition was used.

With these above examples we have been able to demonstrate the use of the idea of an energy map to analyse the dynamics of the flow-induced oscillator a priori, based on our knowledge of the frequency selection mechanism. We have also been able to explain the bizarre non-monotonic amplitude response seen in the case of $X_{e}^{\ast }=0.33$ using this tool.