Prediction of energy harvesting efficiency through a wake–foil interaction model for oscillating foil arrays

Abstract

This research investigates the wake–foil interactions between two oscillating foils in a tandem configuration undergoing energy harvesting kinematics. Oscillating foils have been shown to extract hydrokinetic energy from free-stream flows through a combination of periodic heave and pitch motions, at relatively higher amplitudes and lower reduced frequency than thrust generating foils. When placed in tandem, the wake–foil interactions can govern the energy harvesting efficiency of the system due to a reduced relative flow velocity in combination with a structured and coherent wake of vortices shed from the high amplitude flapping of upstream foils. This work utilizes simulations of two tandem foils to parameterize and model the energy harvesting performance as a function of array configuration and foil kinematics. Once the wake of the leading foil has been fully parameterized, the placement, phase angle and kinematic stroke of the second foil is utilized to estimate the time-dependent power curve. The algorithm predicts the power of the second foil through the mean and unsteady wake characteristics, including the direct impingement of a vortex with the trailing foil.

1. Introduction

In this paper we investigate wake–foil interactions and the impact on power production within arrays of oscillating foils. An oscillating foil can generate power through a periodic pitch and heave motion. In contrast to a thrust producing foil, oscillating foils for energy harvesting are drag producing, and operate at lower reduced frequencies and higher pitch and heave amplitudes (Kinsey & Dumas 2008). As a result of these high amplitudes, the kinematic motion produces a sequence of vortices that form a structured wake, in which the vortex pattern and wake topology are a complex function of the foil kinematics (Ribeiro & Franck 2022b). In array configurations, these unsteady vortices convect downstream and can significantly impact performance of trailing foils. The magnitude of the impact depends on foil kinematics, spatial configuration and inter-foil phase angle, and has been shown to both increase and decrease power generation (Ashraf et al. 2011; Xu & Xu 2017). The goal of this research is to predict the time-dependent power coefficient of a trailing foil in a two-foil array configuration utilizing the wake profile of a single foil.

Power generation in an oscillating foil can be produced through heave and through pitch motions. The heave power is composed of the lift force times the heave velocity, whereas the pitch power is defined by the torque about the pitching axis times the pitch velocity. Both heave and pitch power are a function of the oscillation frequency, heave and pitch amplitude (Xiao & Zhu 2014; Young, Lai & Platzer 2014; Laws & Epps 2016; Wu et al. 2020). For motions with high energy harvesting efficiency, heave power dominates since the average pitch power is close to zero (Kinsey & Dumas 2008; Zhu 2011; Ribeiro & Franck 2022a). During the heave stroke, power is augmented by the formation and shedding of a coherent leading edge vortex (LEV) as the associated low pressure region causes an increase in the lift force (Baik et al. 2012; Ribeiro & Franck 2019). Furthermore, the LEV strength is directly associated with the foil's relative angle of attack, and in a power generation regime, a stronger LEV is desired. Depending on the foil parameters, a trailing edge vortex (TEV) can also form, and/or more than one LEV, forming a two-dimensional structured wake intricately linked to the underlying foil kinematics (Ribeiro & Franck 2022b).

Due to the multiple degrees of freedom in the oscillating foil motion, there is a wide range of kinematics that yield high efficiency power conversions. Thus introducing an array of two foils, each with their own oscillation kinematics, whose relative spacing and timing must be determined, is an enormous parameter space only partially explored. Many researchers have considered a tandem array configuration with the same kinematic motion for both foils, varying only the inter-foil phase, $\psi$, and inter-foil spacing, $S_x$. Numerical (Ashraf et al. 2011; Broering & Lian 2012; Broering, Lian & Henshaw 2012; Xu, Sun & Tan 2016; Xu & Xu 2017; Ma et al. 2019) and experimental (Platzer et al. 2009; Kinsey et al. 2011; Karakas & Fenercioglu 2017; Oshkai et al. 2022) work show these two parameters greatly affect array performance due to the timing of wake–foil interactions. To establish a relationship between the trailing foil motion and the oncoming wake, Kinsey & Dumas (2012) defined a global phase parameter, $\varPhi$, combining inter-foil phase with the wake trajectory assuming mean convection at the free-stream velocity,

(1.1)\begin{equation} \varPhi = 2{\rm \pi} \frac{S_x \,f}{U_\infty} + \psi, \end{equation}

where $U_\infty$ is the free-stream velocity and $f$ is the oscillation frequency. More recently, by quantifying the mean wake velocity from various leading foil kinematics, Ribeiro et al. (2021) replaced $U_\infty$ with a measured mean wake velocity, generating the wake phase parameter. With this update, they noted that foil performance has a strong relationship with wake phase over a wide range of operating kinematics. A wake phase of $0^{\circ }$ corresponds to the trailing foil oscillating in sync with the wake, directly impinging with vortex structures, whereas a wake phase of $180^{\circ }$ corresponds to high trailing foil efficiency since the motion is out of phase with the wake, avoiding destructive vortex–foil interactions.

The interaction between vortex gusts and foils has been previously analysed in literature in the context of how vortex–body interactions may affect the onset of vortex formation and body loading (Rockwell 1998). More recently, by analysing vortex gusts at different positions with respect to a stationary foil, Peng & Gregory (2015) classified the vortex–foil interactions into three categories: close interaction, very close interaction and collision. Within this classification, they identified changes in the vortex dynamics as it interacts with the leading edge and boundary layer, both of which are functions of vortex–foil proximity, Reynolds number and vortex rotation. Also investigating the effect of vortex rotation, Barnes & Visbal (2018a,b) found that the downwash from a clockwise vortex causes separation and transition to turbulence to be partially suppressed on the upper foil surface, which delays the LEV formation. In contrast, an early interaction between a counterclockwise vortex and the foil is manifested by the rapid flow separation at the leading edge due to the increased angle of attack caused by vortex-induced upwash.

To predict the effects of vortex–foil interactions, Biler et al. (2021) experimentally investigated gusts on a stationary foil and noted similar trends between the gust-induced angle-of-attack profile over time and the transient lift force. Similarly, Turhan, Wang & Gursul (2022) analysed a vortex wake interacting with a stationary foil, and found a directly proportional relationship between the effective angle of attack and lift force.

While most work considers a stationary foil interacting with a vortex gust performance, Xu, Duan & Xu (2017) considers vortex–foil interactions of an oscillating propulsive foil, finding that each interaction translates to an instantaneous change in the lift and effective angle-of-attack profiles. Using a gust-induced angle of attack, Muscutt, Weymouth & Ganapathisubramani (2017) predicted forces on an oscillating virtual foil through the steady-state aerodynamic theory. Although their methodology captures the effects of destructive vortex–foil interactions, whenever there is vortex–foil avoidance, the lift prediction is not as accurate.

This paper focuses on the wake–foil interactions within a two-foil array undergoing high amplitude and high heave oscillations for the purpose of energy harvesting. In this configuration the oscillation kinematics are such that a large coherent wake pattern is introduced, which highly influences the energy conversion efficiency of downstream oscillating foils depending on their distance and phase angle with respect to the lead foil. Unlike prior work on stationary or propulsive foils, the wake interactions are impacting a downstream foil that is also undergoing a high amplitude heave and pitch motion. Thus, the baseline state of the foil (without wake interactions) relies on massive leading edge separation and LEV formation to generate maximum power, a process which can be extenuated, accelerated or diminished due to wake interactions. To shed light on this process, this paper presents a methodology to extract velocity profiles from the wake of a single foil and utilize it to predict the energy efficiency of downstream foils in various configurations (spacing, phase angle or kinematic stroke). Thus, a physics-based approach is developed from mean wake and unsteady vortex–foil interactions to estimate the time-dependent power coefficient in a two-foil turbine array. The power coefficient predictions are then compared against two-foil simulations from Ribeiro et al. (2021), and the model limitations are discussed.

The paper is organized as follows. In § 2 we introduce the numerical methods utilized in this research. In § 3 we develop the correlation between power generation and the wake kinematics. In § 4 we evaluate the model at different conditions and discuss its limitations, and finally, summarize the findings in § 5.

2. Numerical methods

This section introduces the computational data utilized in the analysis, describes the extraction of wake velocity data, and defines the kinematic motion and power generation in tandem two-foil arrays.

2.1. Definition of the foil kinematic motion

To generate the oscillatory motion of the foil, an active kinematic stroke is applied as

(2.1)\begin{equation} h(t)=-h_{o} \cos(2{\rm \pi} f t),\end{equation}

and

(2.2)\begin{equation} \theta(t)=-\theta_{o} \sin(2{\rm \pi} f t),\end{equation}

where $h(t)$ and $\theta (t)$ are heave and pitch, respectively, with the pitching motion about the midpoint of the chord. The reduced frequency of oscillation, $f$, heave amplitude, $h_o$, and pitch amplitude, $\theta _o$, are the parameters that control the foil motion and are non-dimensionalized by the chord length, $c$, and free-stream velocity, $U_\infty$.

The prescribed sinusoidal heave and pitch of the single foil generates a time-varying effective angle of attack, $\alpha (t)$, with respect to the free-stream flow. A representative effective angle of attack is evaluated when the foil is at maximum heave velocity, which occurs at one quarter of the cycle period, $T$, or

(2.3)\begin{equation} \alpha_{T/4} = \alpha(t=0.25T) = \theta_o - \tan^{-1} \left(\frac{2{\rm \pi} fh_o}{U_{\infty}}\right) , \end{equation}

with $\alpha (t)$ assumed to be in radians and the term $2{\rm \pi} f h_o$ is obtained through the time derivative of the heave motion $\dot {h}(t)$ at $t=0.25T$ (Kim et al. 2017).

2.2. Computational data utilized in analysis

With the foil kinematic motion defined, the flow over the oscillating foils is simulated at $Re=1000$ using a second-order accurate finite volume, pressure-implicit split-operator algorithm in ${OpenFOAM}$ (Weller et al. 1998). The foil shape is a $10\,\%$ thick ellipse, which is convenient for tidal energy due to its fore-aft symmetry. A two-dimensional unstructured dynamic mesh is utilized and the refinement analysis along with the validation of the dynamic mesh against a stationary mesh are presented in Ribeiro et al. (2021).

Two data sets are considered in this paper. The first set of simulations are a sweep of kinematics of a single oscillating foil with a fully resolved wake (Ribeiro & Franck 2022b). These simulations are used to extract velocity profiles from the wake under various oscillation kinematics, varying the parameters in (2.1) and (2.2).

The second set of simulations have two foils operating in a tandem array configuration at a fixed distance of 6 chords separation. For each array configuration, the inter-foil phase angle, $\psi$, is varied. The kinematic motion for the leading foil (foil 1) is given as

(2.4a,b)\begin{equation} h_1(t)=-h_{o,1} \cos(2{\rm \pi} f t) \quad \theta_1(t)=-\theta_{o,1} \sin(2{\rm \pi} f t),\end{equation}

and for the trailing foil (foil 2),

(2.5a,b)\begin{equation} h_2(t)=-h_{o,2} \cos(2{\rm \pi} f t + \psi) \quad \theta_2(t)=-\theta_{o,2} \sin(2{\rm \pi} f t + \psi),\end{equation}

where the frequency, $f$, remains constant to maintain the same relative phase separation in each stroke. In addition, to reduce the parameter space, the same amplitudes are applied to both foils ($h_{o,1}=h_{o,2}=h_o$; $\theta _{o,1}=\theta _{o,2}=\theta _o$), although this could be varied in future simulations. The inter-foil phase angle ranges from $-180^{\circ }$ to $+180^{\circ }$ with an increment of $30^{\circ }$. The range of foil kinematic parameters investigated in this paper includes $f = 0.10\unicode{x2013}0.15$, $h_o=0.75\unicode{x2013}1.5$ and $\theta _o=55^{\circ }\unicode{x2013}75^{\circ }$ for a total of $16$ sets of kinematics. The parameter range selected is ideal for oscillating foils in energy harvesting mode (Xiao & Zhu 2014). For more details on the foil kinematics, see Appendix A.

2.3. Quantifying wake velocity

To visualize the parameter space, the schematic in figure 1 displays the foil parameters, the inter-foil spacing, $S_x$, and the swept area, $Y_p$, of a tandem two-foil array where a trailing foil is placed in the leading foil wake. The parameter $Y_p$ is defined as the distance from tip to tip of the foil for a full cycle motion multiplied by the foil span. The wake velocity profile at a streamwise distance $x_w$ downstream from the leading foil is given as

(2.6) \begin{equation} \boldsymbol{u_w}(y,t : x_w) = u(y,t)\hat{\imath} + v(y,t)\hat{\jmath}, \end{equation}

where $u$ and $v$ are the streamwise and cross-flow velocity components, respectively. In this paper the wake velocity extraction is performed at $x_w=5c$, one chord length upstream from the second foil at $x=6c$. The choice of $x_w$ is to provide an accurate representation of the oncoming energy flux in the immediate vicinity of the trailing foil. A recommendation for other configurations (including staggered configurations and in-line of other separation distances) would be to sample the wake at approximately 1 chord upstream of the second foil.

Figure 1. Leading foil (foil 1) parameters and placement of a trailing foil (foil 2) to form a two tandem foil array. Vorticity flow field at time $t$ is illustrated along with vortex window $\Delta y$ at wake probe location $x_w$.

The wake profile at $x_w$ is dramatically influenced by the periodic structure of the lead foil's wake. The presence of a strong vortex will increase the velocity magnitude, add rotation to the flow and strongly impact the relative velocity seen by the trailing foil. To capture how these transient and periodic structures affect the trailing foil, a vortex windowing scheme is implemented. Figure 1 demonstrates this concept of a window, length $\Delta y$, that corresponds to the size of the wake disturbance relative to the trailing foil. As the trailing foil oscillates the vortex window is centred at its leading edge, $y_{LE}$, and translates along $x=x_w$. To determine the optimal window size, the vorticity flow field is visually inspected and the $\Delta y$ distance is selected to encompass the maximum vortex diameter. For the configurations investigated in this paper, a size of $\Delta y = 1.2c$ is sufficient to capture the induced velocity of the primary wake vortex.

With the vortex window defined, the mean wake velocity $u_w$, measured at $x_w$, is the spatial and time-averaged magnitude of $\boldsymbol {u_w}$,

(2.7)\begin{equation} u_w = \frac{1}{T(Y_p+\Delta y)} \int^T_0 \int_{-(Y_p+\Delta y)/2}^{+(Y_p+\Delta y)/2} \sqrt{u^2(y,t)+v^2(y,t)} \, {{\rm d} y}\,{\rm d} t, \end{equation}

where the limits of integration expand beyond $Y_p$ to encompass the energy of vortices that surpass the trailing foil swept area.

2.4. Definition of normalized power generation

The unsteady power, P, generation in oscillating foil arrays is defined as the sum of heave and pitch power. Although both components contribute to energy extraction, heave power dominates as the time-averaged pitch power is close to zero (Kinsey & Dumas 2008; Zhu 2011; Ribeiro & Franck 2022a). Thus, the power coefficient for the leading foil, $C_{p,1}$, is approximated as

(2.8)\begin{equation} C_{p,1}(t) = \frac{P_{extraction}}{P_{available}} = \frac{\dot{h}_1 L_1(t)}{\tfrac{1}{2}\rho U_\infty^3 c}, \end{equation}

where $L$ is the lift force on the foil and $\frac {1}{2}\rho U_\infty ^3 c$ is the total power available from the free-stream velocity per planform area of the foil.

For the trailing foil, two modifications are made. The first is that the power extracted is a function of the phase angle $\psi$ between the operating kinematics, as this determines the nature of the wake–foil interactions. To account for the horizontal spacing, frequency and phase angle, the wake phase parameter (Ribeiro et al. 2021) is utilized, defined as

(2.9)\begin{equation} \varPhi = 2{\rm \pi} \frac{S_x \,f}{u_w} + \psi,\end{equation}

where the lift force and power extraction are both functions of $\varPhi$. This parameter defines a non-dimensional wake wavelength and adjusts the phase angle appropriately. Secondly, the average power available to the trailing foil is defined by the mean wake velocity in (2.7). Thus, the power coefficient for the trailing foil, $C_{p,2}$, is given by

(2.10)\begin{equation} C_{p,2}(\varPhi,t) = \frac{\dot{h}_2 L_2(\varPhi,t)}{\tfrac{1}{2}\rho u_w^3 c}.\end{equation}

The advantages of this definition are displayed in the example kinematics of figure 2. At a wake phase of $\varPhi =180^{\circ }$, the trailing foil is out of phase with the wake wavelength, minimizing interactions. As a result, the power coefficients presented in figure 2(a) align well over the upstroke within $t/T=0\unicode{x2013}0.5$ (and is symmetric on the downstroke for these kinematics). Figure 2(c) confirms that there is no direct vortex–foil impingement at this wake phase. Thus, the normalization of the power curves has taken into consideration the decrease in mean flow due to the average wake deficit generated from the leading foil.

Figure 2. Analysis of vortex–foil interactions at wake phases $\varPhi =180^{\circ }$ and $\varPhi =0^{\circ }$ and foil parameters $f=0.10, h_o=1.0, \theta _o=55^{\circ }$. These phases illustrate the power contributions of the vortex gusts on the trailing foil when, for this case, there is a vortex-foil impingement ($\varPhi =0^{\circ }$) or vortex avoidance ($\varPhi =180^{\circ }$). The instantaneous vorticity flow fields are plotted at $t/T=0.33$ (black markers). (a) Similar power coefficient between foils, (b) power coefficient is affected by vortex–foil interactions, (c) vortex–foil avoidance and (d) vortex–foil impingement.

In contrast, figures 2(b) and 2(d) display the vorticity flow field and the power curve for the same foil kinematics but with a wake phase of $\varPhi =0^{\circ }$. In this configuration the foil intercepts a strong clockwise vortex on its upstroke, causing a large decrease in the instantaneous power coefficient at $t/T=0.33$. In contrast, if the vortex direction was counterclockwise, it may encourage the vortex formation over the downstream foil, and thus, generate a constructive vortex–foil interaction. The latter is typically found in staggered foil arrays (Kinsey & Dumas 2012) while the former is observed in tandem arrangements such as those investigated in this paper. Thus, the differences between these two curves, shaded yellow (negative) and green (positive), represent portions of the cycle where unsteady wake–foil interactions are influencing the power extraction.

3. Prediction model

The unsteady vortex–foil interactions presented in § 2 correspond to the power difference between foils normalized by their respective oncoming flow velocity. In this section, these vortex disturbances are associated with the change in the effective angle of attack of a trailing foil. This relationship will be used to predict the power generation from the trailing foil at different wake phases.

3.1. Effective angle of attack in the presence of vortex disturbances

Using the wake velocity profiles described in § 2, an effective angle of attack is computed within the moving vortex window upstream of the trailing foil. First, an instantaneous velocity vector, $\boldsymbol {u^v}$, is computed by spatially averaging the velocity profile over $\Delta y$,

(3.1)\begin{equation} \boldsymbol{u^v}(\varPhi, t) = \frac{1}{\Delta y} \int_{y_{LE}(t+t^v) - \Delta y/2}^{y_{LE}(t+t^v) + \Delta y/2} \boldsymbol{u_w}(y, \varPhi, t) \, {{\rm d} y}.\end{equation}

Time is shifted by $t^v$,

(3.2)\begin{equation} t^v = \frac{S_x - x_w}{u_w}\end{equation}

to account for the convection time between the measured wake at $x=x_w$ and the trailing foil at $x=S_x$. The value of $y_{LE}(t+t^v)$ corresponds to the trailing foil's leading edge position at $t+t^v$, which is when the gust–foil interactions occur. Next, using the spatially averaged velocity within the vortex window, the effective angle of attack of the trailing foil, $\alpha ^v$, is calculated as

(3.3)\begin{equation} \alpha^v(\varPhi, t) = \theta(t) - \tan^{-1} \left(\frac{-\dot{h}(t) - v^v(\varPhi, t)}{u^v(\varPhi, t)}\right),\end{equation}

as illustrated by the velocity triangle and foil heave velocity in figure 3. Although the effective angle of attack is computed at $x=x_w$ it is assumed constant as the vortex convects from $x=x_w$ to $x=S_x$.

Figure 3. Effective angle of attack of the trailing foil in the presence of wake disturbances.

3.2. Introduction of a reference foil with an equivalent mean flow

To quantify the effects of the vortex gust on the trailing foil, an equivalent foil operating in a uniform flow with velocity $u_w$ is introduced as a reference foil illustrated in figure 4.

Figure 4. Effective angle of attack of the reference foil (no wake disturbances).

Thus, the kinematic parameters are now normalized by $u_w$ and are given as

(3.4a–c)\begin{equation} f^* = \frac{f}{u_w}, \quad h^*(t)=-h_{o} \cos(2{\rm \pi} f^* t), \quad \text{and} \quad \theta^*(t)=-\theta_{o} \sin(2{\rm \pi} f^* t),\end{equation}

where the $^*$ superscript denotes the reference foil. The effective angle of attack of the reference foil is defined by $\alpha ^*$,

(3.5)\begin{equation} \alpha^*(t) = \theta^*(t) - \tan^{-1} \left(\frac{-\dot{h}^*(t)}{u_w}\right).\end{equation}

The difference in effective angle of attack between the reference foil and the trailing foil is thus

(3.6)\begin{equation} \Delta \alpha (\varPhi, t) = \alpha^v(\varPhi, t) - \alpha^*(t),\end{equation}

which is an indication of the strength and direction of the unsteady flow in the near vicinity of the trailing foil.

The power coefficient for the reference foil is then calculated from $L^*$ and $\dot {h}^*$ as

(3.7)\begin{equation} C_p^*(t) = \frac{\dot{h}^*L^*(t)}{\tfrac{1}{2}\rho u_w^3 c}.\end{equation}

The power difference, $\Delta C_p$, between the gust–foil interaction and the equivalent reference foil is

(3.8)\begin{equation} \Delta C_p (\varPhi,t) = C_{p,2}(\varPhi,t) - C_p^*(t) = \frac{\dot{h}_2L_2(\varPhi, t)}{\tfrac{1}{2}\rho u_w^3 c} - \frac{\dot{h}^*L^*(t)}{\tfrac{1}{2}\rho u_w^3 c}.\end{equation}

To illustrate the relationship between these quantities, figure 5 shows $\Delta \alpha$ and $\Delta C_p$ profiles within the upstroke foil motion ($t/T=0\unicode{x2013}0.5$) for the wake phase $\varPhi =0^{\circ }$ with foil parameters $f=0.10$, $h_o=1.0$ and $\theta _o=55^{\circ }$. There is a significant $\Delta \alpha$ drop at approximately $t/T=0.33$, which translates to a vortex–foil interaction that is detrimental to the formation of vortices over the foil. This destructive interaction is also observed by a $\Delta C_p$ drop at approximately the same time.

Figure 5. Comparison between the difference in effective angle of attack and power coefficient at $\varPhi =0^{\circ }$ for foil parameters of $f=0.10$, $h_o=1.0$ and $\theta _o=55^{\circ }$.

3.3. Power prediction based on the effective angle of attack

Next, a relationship between $\Delta \alpha$ and $\Delta C_p$ is derived in order to complete the model for the energy harvesting efficiency of the trailing foil. Following from (3.8), which is the measured difference in power between the trailing foil and an equivalent reference foil, a modelled power difference is constructed utilizing input data from the upstream foil wake kinematics and the kinematics of the trailing foil. Since the effective angle of attack is proportional to the lift force on the foil (Biler et al. 2021), the difference in power, $\Delta \widetilde {C_p}$, is modelled as

(3.9)\begin{equation} \Delta \widetilde{C_p}(\varPhi, t) = \beta (\varPhi)\left( \frac{\dot{h}_2 \alpha^v(\varPhi,t) - \dot{h}^*\alpha^*(t)}{\tfrac{1}{2}\rho u_w^3 c} \right ), \end{equation}

where the parameter $\beta$ represents a coefficient of proportionality between the power coefficient and effective angle of attack. Thus, a model of the time-dependent power coefficient of the trailing foil, $C_p^v$, as a function of wake phase $\phi$, can be constructed as

(3.10)\begin{equation} C_p^v(\varPhi, t) = C_p^* (t)+ \Delta \widetilde{C_p} (\varPhi, t).\end{equation}

The optimal value of $\beta$ (for each phase difference $\varPhi$) can be computed using data from two-foil simulations by minimizing the root-mean-square difference between the model and the instantaneous power profiles, demonstrated by the black lines in figure 6. For this subset of two-foil kinematics, $\beta$ peaks in the vicinity of $\varPhi =0^{\circ }$, coinciding with strong wake–foil interactions. In contrast, at phases closer to $180^{\circ }$, $\beta$ is generally smaller when the vortex–foil interactions are weaker. To maintain periodicity, a sinusoidal equation is fit to the data with a nonlinear least squares regression algorithm, generating a $\beta$ profile (blue curve in figure 6) given by

(3.11)\begin{equation} \beta(\varPhi) = 0.75\cos \left(\frac{\varPhi{\rm \pi}}{180} - 0.12{\rm \pi}\right)+0.41,\end{equation}

with $\varPhi$ given in degrees. A physical interpretation of the parameter $\beta$ is that it accounts for changes in the non-circulatory lift forces as the surrounding flow around the airfoil is modified by the presence of the vortex.

Figure 6. Optimal coefficient of proportionality, $\beta$, determined from two-foil data (black lines) and general $\beta$ profile implemented in the model (blue line).

4. Model evaluation

In this section the model performance is evaluated at different foil kinematics with respect to time and wake phase, and the limitations of the model are discussed.

4.1. Time-dependent power prediction

In figure 7 the instantaneous power coefficient developed by the model is compared against the computed power from a two-foil simulation, and compared against a reference foil of equivalent mean free-stream velocity. Whereas the simulation is the exact value, the reference foil represents the baseline power curve for a steady inlet flow. Three representative oscillation kinematics (cases A, B and C) are shown that span a range of $\alpha _{T/4}$ from $18^{\circ }$ to $28^{\circ }$. For each case, three phase angles are shown at $\varPhi =-30^{\circ }$, $0^{\circ }$ and $30^{\circ }$, representing the regime where vortex–foil interactions have the highest probability of occurring.

Figure 7. Comparison between power prediction ($C_p^v$), trailing foil ($C_{p,2}$) and reference foil ($C_p^*$) at three wake phases and configurations. Kinematics: Case A: $f=0.12$, $h_o=1.0$, $\theta _o=55^{\circ }$; Case B: $f=0.10$, $h_o=1.0$, $\theta _o=55^{\circ }$; Case C: $f=0.12$, $h_o=1.0$, $\theta _o=65^{\circ }$.

Across the data displayed, the occurrence of vortex–foil interactions in the trailing foil are well captured by the model, as indicated by the peaks and troughs in the power coefficient. Each wake–foil interaction may be classified as constructive or destructive depending on the vortex sign and positioning with respect to the foil (Ribeiro et al. 2021). The reference foil (dark red line) is the baseline power for equivalent stroke kinematics without any unsteady vortex interaction. Thus, power achieved above the reference foil is described as a constructive vortex–foil interaction, whereas a power coefficient lower than the reference foil represents a destructive interaction. For the tandem configuration explored in this paper, the majority of the vortex–foil interactions are destructive, resulting in the computed and modelled power coefficient less than the reference foil. Most all of these vortex–foil interactions are captured by the model. However, there are instances where the timing and/or magnitude of the power is either over- or under-predicted. In the data shown, there are also two moments of constructive interaction in cases B and C at $\varPhi =-30^{\circ }$, identified by the instantaneous power coefficient surpassing that of the reference foil. The model captures this increase although the amplitude is under-predicted.

In general, when analysing the time-dependent power in each wake phase, the $\varPhi =0^{\circ }$ power prediction profiles are the closest to the simulation, representing where the maximum vortex–foil interaction is expected by the model. At wake phases just before ($\varPhi =-30^{\circ }$) or just after ($\varPhi =30^{\circ }$) the predicted power may have a time shift compared with the simulation as illustrated in case A at $\varPhi =-30^{\circ }$. This is due to the trailing foil interacting with the wake at an earlier time than the estimated vortex convection time.

These discrepancies between model and simulation tend to increase with higher values of $\alpha _{T/4}$, which are known to produce more chaotic wake structures (Ribeiro et al. 2021). The amplitude deviation between model and simulation can be partially explained by the model considering a uniform $\beta$ profile for all kinematics. With stronger wake vortices, the effects of vortex–foil interactions are more apparent and, thus, the coefficient of proportionality between power and angle of attack may be higher than the $\beta$ given in the universal profile.

4.2. Time-averaged power coefficient as a function of wake phase

The time-averaged power coefficient is computed for the three representative cases and displayed in figure 8 as a function of wake phase. The dark red lines are the time-averaged power coefficient from the reference foil that is independent of wake phase since it assumes a steady flow. The $\alpha _{T/4}$ value increases from cases A to C and, thus, so does the strength of the vortex–foil interactions.

Figure 8. Time-averaged trailing foil power comparison between simulation and prediction model with respect to wake phase for three sets of foil kinematics. The dark red line corresponds to the mean power from the reference foil. Results are shown for (a) $\alpha _{T/4}=18^{\circ }$, (b) $\alpha _{T/4}=23^{\circ }$ and (c) $\alpha _{T/4}=28^{\circ }$.

Case A shows a sinusoidal power trend with respect to the wake phase, which occurs from the continuous and smooth interaction typically found in cases with similar values of relative angle of attack (Ribeiro & Franck 2022b). Overall, the model is able to capture the mean power trend at this low angle of attack with a slight under-prediction in magnitude. Similarly, a sinusoidal power variation is shown by the simulation data for case B. However, the model starts showing a more localized power variation around $\varPhi =0^{\circ }$ and a roughly constant power prediction at phase angles farther away from the main vortex interaction at $\varPhi =0^{\circ }$. Case C corresponds to a higher power of the reference foil since a stronger LEV and higher lift are found at the foil's mid-stroke position (Ribeiro, Frank & Franck 2020). Furthermore, a stronger and localized mean power variation is observed due to the more coherent wake vortices and stronger vortex–foil interactions compared with the other cases (Ribeiro et al. 2021). For instance, as discussed by Ribeiro et al. (2021), the circulation of the primary vortex in case C is approximately $\varGamma /(U_\infty c) \approx 1.1$, which is higher than case B ($\varGamma /(U_\infty c) \approx 0.9$), and case A ($\varGamma /(U_\infty c) \approx 0.5$). For cases with leading foil kinematics at high relative angles of attack ($\alpha _{T/4} \geq 28^ {\circ }$), a TEV starts to form that also contributes to stronger wake–foil interactions. The strength of such vortices were previously correlated with $\alpha _{T/4}$ in the literature (Ribeiro et al. 2021; Lee et al. 2022). The more localized power variation explains the non-sinusoidal trend as the vortex–foil interactions are stronger and more localized compared with lower $\alpha _{T/4}$. The prediction model, however, is still able to capture this power trend shift as a function of wake phase and foil kinematics.

4.3. Model capabilities and limitations

The proposed wake–foil interaction model is able to predict the effects on instantaneous power coefficient as a result of constructive and destructive vortex–foil interactions. A single-foil simulation is utilized for input data, extracting the power coefficient of the baseline flow and the time-dependent wake data. Using this information, several proposed configurations and kinematics of the second foil can be easily modelled. These can be expanded beyond the tandem configurations currently proposed to include staggered configurations and various combinations of kinematic strokes.

Although the model predicts events when there is a direct impingement or weak interactions, there is a magnitude mismatch at wake phases close to $|\varPhi | \sim 90^{\circ }$, as seen in figure 8(b). A reason for the larger error at these phases is the failure of the vortex window to capture the entire wake disturbance that affects the power distribution. A potential solution may be to increase the vortex window size to include wake disturbances that have secondary effects on power generation. The vortex window size may become a limitation especially when considering cases with high $\alpha _{T/4}$ ($\alpha _{T/4} > 30^{\circ }$) since the vortex wakes in these cases contain not only LEVs but also TEVs and stronger secondary vortices that can affect power generation (Ribeiro et al. 2021). For more details, Appendix B presents an error quantification between model and simulations based on the mean power coefficient.

Another potential improvement to the model is to account for the vortex trajectory and variation in vortex convection speed within the wake. It is assumed the primary vortex moves at a constant speed, however, our measurements have shown it does vary within the wake region. For this reason, the model relies on sampling the velocity field in close proximity (approximately 1 chord length upstream) to the second foil for the best estimate of the impact of the vortex on the local flow field. If the model incorporated a vortex trajectory and better convection speed, the sampled position could be moved upstream and provide more flexibility to the model.

A final consideration is that the proposed model only considers two-dimensional flows, which is often a good assumption for the high aspect ratio wings deployed for energy harvesting. Prior experimental work has shown that adding end plates maintains an approximately two-dimensional wake and improves the efficiency of the oscillating foils by mitigating tip losses (Kim et al. 2017). Furthermore, these results have good agreement with simulations at a matching Reynolds number (Ribeiro et al. 2020). However, an extension to this model could incorporate a spanwise profile that accounts for the tip vortex and associated loss of lift and power.

5. Conclusion

In this paper we develop a power prediction model for oscillating foil turbine arrays where the foil arrangement is governed by the spacing between the two foils and the wake phase. The goal is to develop an estimation of the time-dependent power curve of a trailing foil using only the wake velocity data from single-foil simulations and the proposed position, phase angle and stroke kinematics of the trailing foil. Such a prediction eliminates the need for costly simulations exploring the wide configuration and kinematic parameter space for an array of two oscillating foils.

First, the model introduces a reference foil, which normalizes the expected power coefficient of the trailing foil with respect to the reduced wake velocity. As the first foil extracts a high percentage of energy from the free-stream flow, the mean wake velocity available to the trailing foil is decreased. Normalizing by this new reference velocity provides a baseline power coefficient for the trailing foil. This rescaling, however, does not capture the unsteady interactions of the trailing foil with incoming vortices in the wake.

The vortex–foil interactions are modelled as deviations from the reference foil's power curve. It is assumed that these time-dependent deviations in power production are proportional to the difference in relative angle of attack in the vicinity of the trailing foil. Using the kinematics of the trailing foil, a moving window is constructed to quantify the local velocity magnitude and relative angle of attack with respect to the heaving and pitching trailing foil. The instantaneous velocity vectors are extracted from the unsteady wake data of a single-foil simulation. To complete the model, a coefficient of proportionality is computed from available two-foil simulations and found to be a function of wake phase.

The results show that the prediction model captures both power trends and magnitudes across the range of wake phases and foil kinematics explored. Depending on the wake phase, the model prediction can be remarkably close to the simulation, especially at wake phases close to a direct vortex–foil impingement ($\varPhi \sim 0^{\circ }$) and regions of minimal vortex interaction ($\varPhi \sim 180^{\circ }$). At wake phases in between, typically around $|\varPhi | \sim 90^{\circ }$, the differences between model and simulation are more apparent. This is likely due to the shortcoming of the vortex window in capturing secondary wake disturbances and, thus, correctly matching with the trailing foil power variation. Additionally, in cases with $\alpha _{T/4} > 30^{\circ }$ the wake vortices are stronger and multiple vortices are interacting with the foil, which makes the prediction more challenging.

The advantage of this model is the ability to predict the time-dependent power over a range of potential two-foil configurations based solely on single-foil simulations. Although a limited set of kinematics and configurations are explored in this paper, the model can be applied to two-foil systems operating with different kinematic parameters and staggered configurations, both of which can improve the overall efficiency of the system.