2.1 Phase-reduction analysis of periodic flows

Let us consider the dynamics of incompressible flow, represented as

(2.1) $$\begin{eqnarray}\dot{\boldsymbol{q}}=\boldsymbol{F}(\boldsymbol{q}),\end{eqnarray}$$

where $\boldsymbol{q}$ is the flow variable which has a stable periodic flow (limit cycle) $\boldsymbol{q}_{0}$ satisfying

(2.2) $$\begin{eqnarray}\boldsymbol{q}_{0}(\boldsymbol{x},t+T)=\boldsymbol{q}_{0}(\boldsymbol{x},t).\end{eqnarray}$$

For this limit cycle, the periodicity is denoted by $T$ which gives the natural frequency of $\unicode[STIX]{x1D714}_{n}=2\unicode[STIX]{x03C0}/T$ . As expressed in (2.2), flows are generally described in terms of spatial and temporal variables. For incompressible flow, we can take $\boldsymbol{q}$ to be the velocity vector. In this work, we instead describe the state of the oscillatory flow using the phase $\unicode[STIX]{x1D703}$ of the limit cycle, such that

(2.3) $$\begin{eqnarray}\dot{\unicode[STIX]{x1D703}}(t)=\unicode[STIX]{x1D714}_{n},\text{ where }\unicode[STIX]{x1D703}(t)\in [-\unicode[STIX]{x03C0},\unicode[STIX]{x03C0}].\end{eqnarray}$$

The main idea here is to reduce the description of the flow from the flow state $\boldsymbol{q}(\boldsymbol{x},t)$ to a single scalar phase variable $\unicode[STIX]{x1D703}$ (Winfree Reference Winfree1967; Guckenheimer Reference Guckenheimer1975; Kuramoto Reference Kuramoto1984). If we have a periodic flow, the phase can point back to the full periodic flow field $\boldsymbol{q}_{0}(\unicode[STIX]{x1D703})$ . The description of the phase is given in figure 1 for an example of periodic cylinder flow, which will be examined in detail later in § 3. In this example, the phase dynamics of the flow can be understood by examining the lift $C_{L}$ on the cylinder, which holds the same synchronization properties as the flow. The limit cycle of periodic cylinder flow is shown in blue over the ${\dot{C}}_{L}$ – $C_{L}$ phase space, with the vorticity fields visualized at representative phases.

Figure 1. Definition of phase $\unicode[STIX]{x1D703}$ and phase function $\unicode[STIX]{x1D6E9}$ with $C_{L}$ and ${\dot{C}}_{L}$ for cylinder flow at a Reynolds number of $Re=100$ . Vorticity fields are shown at representative phases.

The above characterization of the phase is limited to flow states that reside on the limited cycle $\boldsymbol{q}_{0}$ . We can extend the concept of phase in the basin of the limit cycle and define a phase function $\unicode[STIX]{x1D6E9}(\boldsymbol{q}(\boldsymbol{x},t))$ that returns the phase, i.e. $\unicode[STIX]{x1D703}=\unicode[STIX]{x1D6E9}(\boldsymbol{q})$ . The phase function can capture the phase dynamics off the limit cycle, as illustrated in red in figure 1. By combining (2.1) and (2.3), we observe that the phase function satisfies

(2.4) $$\begin{eqnarray}\dot{\unicode[STIX]{x1D703}}(t)=\dot{\unicode[STIX]{x1D6E9}}(\boldsymbol{q})=\unicode[STIX]{x1D735}_{\boldsymbol{q}}\unicode[STIX]{x1D6E9}(\boldsymbol{q})\boldsymbol{\cdot }\dot{\boldsymbol{q}}=\unicode[STIX]{x1D735}_{\boldsymbol{q}}\unicode[STIX]{x1D6E9}(\boldsymbol{q})\boldsymbol{\cdot }\boldsymbol{F}(\boldsymbol{q})=\unicode[STIX]{x1D714}_{n}.\end{eqnarray}$$

This phase function enables the phase to be defined over the entire phase space. A surface of constant phase is referred to as an isochron and is directed perpendicular to $\unicode[STIX]{x1D735}\unicode[STIX]{x1D6E9}(\boldsymbol{q})$ (Guckenheimer Reference Guckenheimer1975). The notion of isochron is closely related to the eigenfunction of the Koopman operator (Mauroy & Mezić Reference Mauroy and Mezić2012).

2.2 Phase-reduction analysis of perturbed flows

Now, let us consider perturbing the governing dynamics by

(2.5) $$\begin{eqnarray}\dot{\boldsymbol{q}}=\boldsymbol{F}(\boldsymbol{q})+\unicode[STIX]{x1D700}\boldsymbol{f}(t),\end{eqnarray}$$

where we have added a small perturbation $\unicode[STIX]{x1D700}\boldsymbol{f}(t)$ with $\unicode[STIX]{x1D700}\ll 1$ and $\Vert \,\boldsymbol{f}\Vert =1$ . The dynamics of the phase under perturbation can be determined through

(2.6) $$\begin{eqnarray}\dot{\unicode[STIX]{x1D703}}(t)=\dot{\unicode[STIX]{x1D6E9}}(\boldsymbol{q}(t))=\unicode[STIX]{x1D735}_{\boldsymbol{q}}\unicode[STIX]{x1D6E9}(\boldsymbol{q})\boldsymbol{\cdot }\dot{\boldsymbol{q}}(t)=\unicode[STIX]{x1D735}_{\boldsymbol{q}}\unicode[STIX]{x1D6E9}(\boldsymbol{q})\boldsymbol{\cdot }[\boldsymbol{F}(\boldsymbol{q})+\unicode[STIX]{x1D700}\boldsymbol{f}(t)].\end{eqnarray}$$

We can evaluate $\unicode[STIX]{x1D735}_{\boldsymbol{q}}\unicode[STIX]{x1D6E9}$ along the limit cycle with $\boldsymbol{q}_{0}$ , assuming that the higher-order terms are negligible since the added perturbation is weak. Recalling (2.4), we then arrive at

(2.7) $$\begin{eqnarray}\dot{\unicode[STIX]{x1D703}}(t)=\unicode[STIX]{x1D714}_{n}+\unicode[STIX]{x1D700}\boldsymbol{Z}(\unicode[STIX]{x1D703})\boldsymbol{\cdot }\boldsymbol{f}(t),\end{eqnarray}$$

where $\boldsymbol{Z}(\unicode[STIX]{x1D703})\equiv \unicode[STIX]{x1D735}_{\boldsymbol{q}}\unicode[STIX]{x1D6E9}(\boldsymbol{q})|_{\boldsymbol{q}=\boldsymbol{q}_{0}}$ is called the phase-sensitivity function. Given the perturbation function, we can determine the influence of the perturbation function on the phase dynamics through (2.7) if $\boldsymbol{Z}(\unicode[STIX]{x1D703})$ is known.

This phase-sensitivity function can be determined in two ways (Nakao Reference Nakao2016). We can find $\boldsymbol{Z}(\unicode[STIX]{x1D703})$ by examining the phase response to impulse perturbations (direct method) or by considering the adjoint dynamics (adjoint method; Ermentrout & Terman Reference Ermentrout and Terman2010). In this work, we consider the former impulse-response approach which can be used in both computational and experimental studies. We consider here the influence of a weak impulse perturbation $\boldsymbol{I}=I\unicode[STIX]{x1D6FF}(t-t_{0})h(\boldsymbol{x})\hat{\boldsymbol{e}}_{f}$ on the state variable $\boldsymbol{q}$ , where $\hat{\boldsymbol{e}}_{f}$ is the unit vector in the forcing direction and $h(\boldsymbol{x})$ describes the spatial profile of the impulse. This impulse can be added to the governing dynamics, such that

(2.8) $$\begin{eqnarray}\dot{\boldsymbol{q}}=\boldsymbol{F}(\boldsymbol{q})+\boldsymbol{I}\end{eqnarray}$$

at some phase $\unicode[STIX]{x1D703}$ (at a chosen $t_{0}$ ). With the introduction of this impulse, the phase dynamics becomes affected, causing a shift in the phase along the limit cycle state at steady oscillating conditions.

The phase change at asymptotic state in response to the impulse is known as the phase-response function $g(\unicode[STIX]{x1D703};I\hat{\boldsymbol{e}}_{f})$ and is dependent on the phase at which the impulse is added as well as its location and strength, as illustrated in figure 2(a). Given a weak perturbation, the phase-sensitivity function can be evaluated as

(2.9) $$\begin{eqnarray}Z_{j}(\unicode[STIX]{x1D703})=\lim _{I\rightarrow 0}\frac{g(\unicode[STIX]{x1D703};I\hat{\boldsymbol{e}}_{f})}{I}\approx \frac{g(\unicode[STIX]{x1D703};I\hat{\boldsymbol{e}}_{f})}{I}.\end{eqnarray}$$

By introducing the impulse perturbations over a range of $\unicode[STIX]{x1D703}$ , the phase-sensitivity function can be determined, as shown in figure 2(b). If the phase dynamics of the state variable $\boldsymbol{q}$ is followed by an observable (sensor output), the observable, such as lift in figure 2(a), can be used to determine the phase response. This enables experimental efforts with limited measurements to determine the phase-sensitivity function without difficulty.

Figure 2. (a) Determination of the phase-response function $g(\unicode[STIX]{x1D703};\boldsymbol{I})$ from the impulse input. (b) Phase-sensitivity function. (c) Phase-coupling function and stability condition for the phase difference $\unicode[STIX]{x1D719}$ between the system dynamics and the oscillatory actuation. For graphical clarity, $I=0.025$ is used here. All subsequent phase-reduction analysis uses $I=0.01$ .

2.3 Synchronization of phase dynamics to external oscillation

With the oscillatory dynamics of fluid flow reduced to its phase dynamics, we can further analyse the condition for which the flow exhibits synchronization to an external periodic perturbation $\boldsymbol{f}(t)$ with a frequency $\unicode[STIX]{x1D714}_{f}$ (and period $T_{f}=2\unicode[STIX]{x03C0}/\unicode[STIX]{x1D714}_{f}$ ) that is close to $\unicode[STIX]{x1D714}_{n}$ . To facilitate our analysis, we consider the relative phase $\unicode[STIX]{x1D719}(t)$ between the phase $\unicode[STIX]{x1D703}(t)$ of the flow and the phase of the perturbation $\unicode[STIX]{x1D714}_{f}t$ ,

(2.10) $$\begin{eqnarray}\unicode[STIX]{x1D719}(t)\equiv \unicode[STIX]{x1D703}(t)-\unicode[STIX]{x1D714}_{f}t.\end{eqnarray}$$

Using (2.7), we find that the time rate of change of the relative phase is

(2.11) $$\begin{eqnarray}\dot{\unicode[STIX]{x1D719}}(t)=\unicode[STIX]{x1D714}_{n}-\unicode[STIX]{x1D714}_{f}+\unicode[STIX]{x1D700}\boldsymbol{Z}(\unicode[STIX]{x1D719}(t)+\unicode[STIX]{x1D714}_{f}t)\boldsymbol{\cdot }\boldsymbol{f}(t)=\unicode[STIX]{x1D700}[\unicode[STIX]{x1D6E5}+\boldsymbol{Z}(\unicode[STIX]{x1D719}(t)+\unicode[STIX]{x1D714}_{f}t)\boldsymbol{\cdot }\boldsymbol{f}(t)],\end{eqnarray}$$

where the difference in the natural frequency $\unicode[STIX]{x1D714}_{n}$ and the perturbation frequency $\unicode[STIX]{x1D714}_{f}$ is expressed as $\unicode[STIX]{x1D700}\unicode[STIX]{x1D6E5}=\unicode[STIX]{x1D714}_{n}-\unicode[STIX]{x1D714}_{f}$ , with $\unicode[STIX]{x1D6E5}=\mathit{O}(1)$ .

By examining the stability characteristics, we can determine when the oscillatory dynamics of the flow can synchronize to external oscillation. That is, when $|\dot{\unicode[STIX]{x1D719}}|\rightarrow 0$ , the oscillation frequency of the flow converges towards the actuator frequency. At this moment, equation (2.11) describes the dynamics of the phase difference $\unicode[STIX]{x1D719}(t)$ in a non-autonomous manner. To perform a stability analysis, it is desirable to have (2.11) expressed in autonomous form. By realizing that $\unicode[STIX]{x1D719}(t)$ varies slowly compared with $\unicode[STIX]{x1D714}_{f}t$ , we can approximate the right-hand side of (2.11) by its average over a period of actuation (Kuramoto Reference Kuramoto1984; Ermentrout & Terman Reference Ermentrout and Terman2010) to arrive at

(2.12) $$\begin{eqnarray}\dot{\unicode[STIX]{x1D719}}(t)=\unicode[STIX]{x1D700}[\unicode[STIX]{x1D6E5}+\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})],\end{eqnarray}$$

where

(2.13) $$\begin{eqnarray}\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})=\frac{1}{T_{f}}\int _{t}^{t+T_{f}}\boldsymbol{Z}(\unicode[STIX]{x1D719}+\unicode[STIX]{x1D714}_{f}\unicode[STIX]{x1D70F})\boldsymbol{\cdot }\boldsymbol{f}(\unicode[STIX]{x1D70F})\,\text{d}\unicode[STIX]{x1D70F}=\frac{1}{2\unicode[STIX]{x03C0}}\int _{-\unicode[STIX]{x03C0}}^{\unicode[STIX]{x03C0}}\boldsymbol{Z}(\unicode[STIX]{x1D719}+\unicode[STIX]{x1D713})\boldsymbol{\cdot }\boldsymbol{f}(\unicode[STIX]{x1D713}/\unicode[STIX]{x1D714}_{f})\,\text{d}\unicode[STIX]{x1D713}.\end{eqnarray}$$

The function $\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})$ is $2\unicode[STIX]{x03C0}$ -periodic and is known as the phase-coupling function. With (2.12), we now have a one-dimensional autonomous evolution equation for  $\unicode[STIX]{x1D719}$ .

Through inspection of (2.12), we find that the dynamics exhibits stable behaviour if $\min _{\unicode[STIX]{x1D719}}\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})<-\unicode[STIX]{x1D6E5}<\max _{\unicode[STIX]{x1D719}}\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})$ , which translates to

(2.14) $$\begin{eqnarray}\unicode[STIX]{x1D700}\min _{\unicode[STIX]{x1D719}}\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})<\unicode[STIX]{x1D714}_{f}-\unicode[STIX]{x1D714}_{n}<\unicode[STIX]{x1D700}\max _{\unicode[STIX]{x1D719}}\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719}).\end{eqnarray}$$

The phase-coupling function $\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})$ is shown for the example of cylinder flow in figure 2(c) with an illustration of the stability condition. The value of the right-hand side of (2.12) describes the directions of the solution with respect to the equilibrium points, which are indicated by the arrows in figure 2(c). These conditions determine the region of stability over the $\unicode[STIX]{x1D714}_{f}$ – $\unicode[STIX]{x1D700}$ space, known as the Arnold tongue (Arnold Reference Arnold1997). This will be shown later in § 3.3 (see figure 5). The decay of the frequency difference between the natural limit cycle and the forcing input constitutes the synchronization of the dynamics (they can have a constant offset in phase under synchronization). Strictly speaking, the above condition for synchronization holds for small forcing input $\unicode[STIX]{x1D700}\ll 1$ . While the above analysis assumes that $\unicode[STIX]{x1D714}_{f}$ and $\unicode[STIX]{x1D714}_{n}$ are close to one another for the averaging operation to be a valid approximation, the analysis can be easily extended to consider subharmonic and superharmonic frequencies.

For fluid flows, such synchronization amounts to the so-called lock on. The above result enables us to determine the combination of actuation frequency and amplitude for the flow to synchronize with periodic excitations. As long as we can examine the system response against impulse inputs, the overall formulation can be performed through numerical simulations or experiments, which is attractive for fluid mechanics.

3.1 Model problem

We apply phase-reduction analysis to examine the synchronization of a circular cylinder wake to periodic forcing. The base periodic flow is taken to be the two-dimensional incompressible flow over a circular cylinder at a diameter-based Reynolds number of $Re\equiv U_{\infty }d/\unicode[STIX]{x1D708}=100$ , where $U_{\infty }$ , $d$ and $\unicode[STIX]{x1D708}$ denote the characteristic free-stream velocity, cylinder diameter and kinematic viscosity of the flow respectively. The dynamics of the flow (2.1) is governed by the non-dimensional incompressible Navier–Stokes equations,

(3.1a,b ) $$\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+Re^{-1}\unicode[STIX]{x1D6FB}^{2}\boldsymbol{u}+\unicode[STIX]{x1D700}\boldsymbol{f}\quad \text{and}\quad \unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{u}=0,\end{eqnarray}$$

where $\boldsymbol{u}$ is the velocity and $p$ is the pressure. In this equation, the forcing (actuation) input $\unicode[STIX]{x1D700}\boldsymbol{f}$ is taken to be a localized force (momentum injection)

(3.2) $$\begin{eqnarray}\unicode[STIX]{x1D700}\boldsymbol{f}=\unicode[STIX]{x1D6FF}(\boldsymbol{x}-\boldsymbol{x}_{f})\unicode[STIX]{x1D702}(t)\hat{\boldsymbol{e}}_{f}\end{eqnarray}$$

having a magnitude of $\unicode[STIX]{x1D700}\ll 1$ , position of $\boldsymbol{x}_{f}$ and forcing direction specified by a unit vector $\hat{\boldsymbol{e}}_{f}$ . The Dirac delta in space is approximated using the three-cell discrete delta function (Roma, Peskin & Berger Reference Roma, Peskin and Berger1999). The overall set-ups of the forcing input for the two cases considered below are illustrated in figure 3.

Figure 3. The cylinder flow forcing set-up to analyse (a) the impulse response with a top actuator only and (b) the effect of periodic forcing with top and bottom actuators. The selected forcing directions are the $x$ direction ( $x$ ), the $y$ direction ( $y$ ), tangential ( $t$ ) and normal ( $n$ ).

The function $\unicode[STIX]{x1D702}(t)$ specifies the temporal form of the forcing input as

(3.3) $$\begin{eqnarray}\unicode[STIX]{x1D702}(t)=\left\{\begin{array}{@{}ll@{}}I\unicode[STIX]{x1D6FF}(t)\approx {\displaystyle \frac{I}{\sqrt{2\unicode[STIX]{x03C0}}\unicode[STIX]{x1D70E}}}\exp \left[-{\displaystyle \frac{1}{2}}\left({\displaystyle \frac{t-t_{0}}{\unicode[STIX]{x1D70E}}}\right)^{2}\right]\quad & \text{for impulse perturbation},\\ {\displaystyle \frac{\unicode[STIX]{x1D700}}{2}}[1+\sin (\unicode[STIX]{x1D714}_{f}t+\unicode[STIX]{x1D713})]\quad & \text{for periodic perturbation.}\end{array}\right.\end{eqnarray}$$

The impulse perturbation, approximated as a narrow Gaussian in time (with $\unicode[STIX]{x1D70E}=0.005$ ), is introduced to the flow to determine the phase-response, phase-sensitivity and phase-coupling functions. The sinusoidal forcing function is chosen to examine the lock-on behaviour of the flow to periodic flow control input. While the impulse perturbation is introduced only at the top separation point ( $58.4^{\circ }$ from the aft separation point, positioned $3\unicode[STIX]{x0394}x$ away from the surface), we add the sinusoidal forcing input at the top and bottom separation points with a phase difference of $\unicode[STIX]{x1D719}$ (top $\unicode[STIX]{x1D713}=0$ and bottom $\unicode[STIX]{x1D713}=\unicode[STIX]{x03C0}$ ). The actuation is introduced near the separation points to effectively modify the wake dynamics. However, the perturbations can be added to any locations in the flow field for this analysis. The forcing amplitudes are chosen such that $I=0.01$ and $\unicode[STIX]{x1D700}\in [0,0.1]$ . The corresponding steady and oscillatory momentum coefficients for the periodically forced cases are less than 0.02 and 0.007 respectively for $\unicode[STIX]{x1D700}\leqslant 0.1$ (Munday & Taira Reference Munday and Taira2013).

3.2 Computational set-up

The immersed boundary projection method simulates the unsteady flow on a Cartesian grid with the cylinder represented through the boundary forcing (Taira & Colonius Reference Taira and Colonius2007). For the present study, the fast version of the immersed boundary projection method is chosen (Colonius & Taira Reference Colonius and Taira2008), with the computational domain extending over $(x/d,y/d)\in [-31,33]\times [-32,32]$ . The grid near the cylinder has $\unicode[STIX]{x0394}x_{min}/d=\unicode[STIX]{x0394}y_{min}/d=0.02$ and the time step is chosen to meet $U_{\infty }\unicode[STIX]{x0394}t/\unicode[STIX]{x0394}x_{min}=0.5$ . The computational algorithm and set-up have been validated extensively in past studies, including a similar study on control of cylinder flow (Munday & Taira Reference Munday and Taira2013). The present study captures the baseline flow accurately and agrees with findings from three independent numerical codes, as summarized in table 1.

Table 1. Comparison of lift, drag and Strouhal number ( $St$ ) for circular cylinder flow at $Re=100$ .

3.3 Analysis of the cylinder wake

In this section, we determine the phase-response, phase-sensitivity and phase-coupling functions to characterize the cylinder flow dynamics in terms of phase on the ${\dot{C}}_{L}$ – $C_{L}$ space. Once these functions are determined, we can find the conditions for synchronization of the cylinder wake to periodic perturbations.

We first characterize the phase response of the cylinder wake to weak impulse perturbations added near the top separation point. By evaluating the change in phase, we can determine the phase-response function $g(\unicode[STIX]{x1D703};I\hat{\boldsymbol{e}}_{f})$ . The magnitude of the impulse is set to $I=0.01$ , which is small in magnitude compared with the force on the cylinder. Here, we consider adding the impulse in four different directions ( $x$ , $y$ , surface-tangential and surface-normal directions) from the top separation point, as illustrated in figure 3. These directions are chosen to form two pairs of orthogonal directions, such that they can serve as bases for subsequent analyses.

Based on the impulse responses, we can compile in figure 4(a) the phase-sensitivity function $Z(\unicode[STIX]{x1D703})$ evaluated with (2.9) over the phase space for the four directions, shown by the solid lines. Symbols represent where actual simulations were performed; 10 simulations were performed for each direction. Impulse added in the tangential and $x$ directions leads to phase advancement. The tangential forcing exhibits the largest phase sensitivity among the four directions examined. Forcing in the $y$ direction causes phase delay, and normal forcing does not result in much of a phase change. We also show with dashed lines in figure 4(a) the phase-sensitivity functions found by rotating the coordinate systems. That is, the $x$ and $y$ components are evaluated based on the tangential and normal components, and vice versa. The close agreement between the solid and dashed lines suggests that the chosen impulse magnitude is adequately small for the linear approximations made in (2.7) and (2.9) to hold. Hence, it suffices to perform the impulse-response analysis in two independent forcing directions to map out the phase-sensitivity function over the phase space. In this example, the phase-sensitivity function quantifies the sensitivity of the shear layer dynamics aft of the separation point (as well as the overall wake vortex dynamics) to the added impulse perturbation.

Figure 4. (a) The phase-sensitivity function $Z(\unicode[STIX]{x1D703})$ and (b) phase-coupling function $\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})$ for the cylinder wake and periodic forcing input. The dashed lines for the phase-sensitivity function show approximations from complementary bases. Symbols represent where actual simulations were performed.

With the phase-sensitivity function $Z(\unicode[STIX]{x1D703})$ available, the phase-coupling function $\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})$ can be found through (2.13). For the present periodic forcing case, we need to sum the influences of the top and bottom actuators, which essentially doubles the magnitude of $\unicode[STIX]{x1D702}(t)$ in the computation of $\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})$ . The phase-coupling functions for the four representative directions are shown in figure 4(b). As expected, large variations in $\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D719})$ over the relative phase $\unicode[STIX]{x1D719}$ are given by the periodic perturbations in the tangential and $x$ directions. Forcing in the $y$ direction has opposing effects with reduced variation. Normal forcing seems to be the least influential among these four directions. It is interesting to note that $Z(\unicode[STIX]{x1D703})$ has a zero average over $\unicode[STIX]{x1D703}$ for the tangential forcing case, which removes the influence of non-zero mean of the forcing input.

Figure 5. Synchronization of the cylinder wake to periodic actuation for (a) tangential and (b) normal forcing inputs. Synchronized and unsynchronized cases from DNS (Munday & Taira Reference Munday and Taira2013) are shown by ○ and $\times$ respectively. The criteria for synchronization (2.14) correspond to the region bounded by the red lines.

It should be recalled that greater variations in the phase-coupling function widen the regions for achieving synchronization of the wake dynamics to actuation fluctuation, as stated in condition (2.14). For example, tangential and normal direction forcing yield $(\min \unicode[STIX]{x1D6E4},\max \unicode[STIX]{x1D6E4})=(-3.52,4.13)$ and $(-1.10,0.33)$ respectively. We can take these minimum and maximum values for the synchronization condition and compare them with the results from direct numerical simulations (DNS) (Munday & Taira Reference Munday and Taira2013), as summarized in figure 5. In the work of Munday & Taira (Reference Munday and Taira2013), the cylinder wake is considered to be synchronized (locked on) to the actuation frequency when the lift power spectrum exhibits a single peak at the actuation frequency.

For the tangential forcing cases shown in figure 5(a), the condition for synchronization from the phase-response analysis agrees well with the full nonlinear wake simulations capturing the region of the Arnold tongue, especially for low forcing amplitude $\unicode[STIX]{x1D700}$ until nonlinear effects become apparent at higher $\unicode[STIX]{x1D700}\gtrsim 0.07$ . As mentioned earlier, the phase-sensitivity function for tangential forcing removes the influence of the non-zero forcing input for synchronization. In fact, the region of lock on is not altered much for a zero-mean forcing case, as presented in the DNS findings (Munday & Taira Reference Munday and Taira2013). For $\unicode[STIX]{x1D700}\gtrsim 0.097$ , the wake becomes synchronized to the oscillatory forcing for all considered actuation frequencies, which is beyond the validity of the present phase-reduction analysis.

For the normal forcing case, we observe that the Arnold tongue is much narrower than that of the tangential forcing case, as is evident from figure 5(b). The phase-reduction analysis predicts the narrow width of the synchronization region well, but slightly biases the tongue towards lower frequency, which becomes obvious for higher forcing amplitudes. This deviating trend from normal forcing is probably caused by the inaccuracy in the evaluation of the small-amplitude phase-sensitivity function as well as by the nonlinear effects.

The identification of the Arnold tongue from the phase-response analysis is attractive since we only require a limited number of simulations to determine the phase-response function over the phase $\unicode[STIX]{x1D703}\in [-\unicode[STIX]{x03C0},\unicode[STIX]{x03C0}]$ . Instead of having to perform a large number of full DNS (or experiments) as a parametric study to locate the Arnold tongue, we can analytically determine the boundary of synchronization from the phase-coupling function. Moreover, for two-dimensional flow, we can consider two independent forcing directions to determine the impulse response for any forcing direction, as depicted in figure 6. By constructing the phase-sensitivity function over the forcing direction $\unicode[STIX]{x1D6FC}$ based on the tangential and normal components or the $x$ and $y$ components, we can determine the minimum and maximum values of $\unicode[STIX]{x1D6E4}$ as a function of $\unicode[STIX]{x1D6FC}$ , as shown in figure 6. In this figure, we shade out the directions that would be directed below the tangential direction (into the body) as well as the upstream direction (acting against the free stream). Restricting our interest to downstream forcing in the range $\unicode[STIX]{x1D6FC}\in [0,\unicode[STIX]{x03C0}/2]$ , we find that synchronization can be achieved for the widest range of $\unicode[STIX]{x1D714}_{f}$ with tangential forcing ( $\unicode[STIX]{x1D6FC}\approx 2^{\circ }$ ). On the other hand, normal forcing would be the most difficult to achieve synchronization with, due to its narrowest bound over $\unicode[STIX]{x1D6FC}$ . We emphasize that the bounds for synchronization have been determined with only a limited number of simulations and do not require an extensive parameter sweep, which significantly reduces the need for computational (and experimental) resources.

Figure 6. Minimum and maximum values of the phase-coupling function $\unicode[STIX]{x1D6E4}$ over the forcing angle $\unicode[STIX]{x1D6FC}$ measured from the tangential direction. Extrema based on the $x$ and $y$ directions as well as the tangential and normal directions are shown. The forcing directions and angle $\unicode[STIX]{x1D6FC}$ are depicted on the left.

While synchronization in this particular example leads to drag increase for the majority of the synchronized cases, phase-response analysis shows tremendous potential to identify the condition of lock on for a range of fluid flow problems. In this example, we know from our previous work (Munday & Taira Reference Munday and Taira2013) that the addition of non-zero mean in tangential forcing can achieve drag reduction and with minimal effect of altering the lock-on region. In fact, the zero-mean case also shows good agreement with the predicted condition for synchronization. There are other cylinder flow examples that can achieve drag reduction, for example the case of synchronizing the wake to the sinusoidal rotation of the cylinder (Tokumaru & Dimotakis Reference Tokumaru and Dimotakis1991).

For the presently considered periodically forced cylinder wake, we do not observe subharmonic or superharmonic synchronizations from DNS. There are of course other fluid flows that exhibit such behaviour. Phase-response analysis can be extended to determine the condition for synchronization near harmonic frequencies.