1. Introduction
The high efficiency and manoeuvrability of aquatic animals (Fish Reference Fish2020) have driven the design of numerous underwater vehicles using bio-inspired propulsive mechanisms (Lauder et al. Reference Lauder, Madden, Tangorra, Anderson and Baker2011; Roper et al. Reference Roper, Sharma, Sutton and Culverhouse2011; Siddall & Kovač Reference Siddall and Kovač2014; Katzschmann et al. Reference Katzschmann, DelPreto, MacCurdy and Rus2018; Zhu et al. Reference Zhu, White, Wainwright, Di, Valentina, George and Bart-Smith2019; White et al. Reference White, Lauder and Bart-Smith2021; Zhong et al. Reference Zhong, Zhu, Fish, Kerr, Downs, Bart-Smith and Quinn2021b ). Now, attention is turning to operating these bio-inspired vehicles as multi-agent collectives or schools (Berlinger et al. Reference Berlinger, Gauci and Nagpal2021) to allow them to perform distributed tasks and to harness the hydrodynamic benefits of schooling, such as the reduced energy expenditure that fish in a school experience (Zhang & Lauder Reference Zhang and Lauder2024). Previous hydrodynamic studies on pairs of swimmers have found that canonical formations such as in-line formations (Boschitsch et al. Reference Boschitsch, Dewey and Smits2014; Ramananarivo et al. Reference Ramananarivo, Fang, Oza, Zhang and Ristroph2016a ; Muscutt et al. Reference Muscutt, Weymouth and Ganapathisubramani2017; Kurt & Moored Reference Kurt and Moored2018; Heydari & Kanso Reference Heydari and Kanso2020; Kurt et al. Reference Kurt, Mivehchi and Moored2021) and side-by-side formations (Verma et al. Reference Verma, Novati, Koumoutsakos and Sethian2018; Kurt et al. Reference Kurt, Panah and Moored2020; Zhong et al. Reference Zhong, Han, Moored and Quinn2021a ; Ormonde et al. Reference Ormonde, Kurt, Mivehchi and Moored2021; Han et al. Reference Han, Zhong, Mivehchi, Quinn and Moored2023) have the potential to substantially increase the thrust and/or efficiency when there is proper synchronisation between swimmers.
Several mechanisms have been proposed to be behind these schooling benefits. Weihs (Reference Weihs1973, Reference Weihs1975) proposed an ideal model of schooling where fish would minimise their energy expenditure by arranging themselves into two-dimensional layers with every other row of fish staggered, giving rise to the now classical diamond formation used as a model arrangement for several studies (Stöcker Reference Stöcker1999; Hemelrijk et al. Reference Hemelrijk, Reid, Hildenbrandt and Padding2015; Pan & Dong Reference Pan and Dong2022; Kelly et al. Reference Kelly, Pan, Menzer and Dong2023; Wei et al. Reference Wei, Hu, Li and Shi2023; Kelly & Dong Reference Kelly and Dong2024). Weihs argued that this formation would exploit several mechanisms at once, leading to reduced energy expenditure. One benefit would occur when fish in each row would ideally synchronise with an out-of-phase relationship to their nearest neighbour. This would lead to direct body – body interactions akin to ground effect that would increase their thrust (Dong & Lu Reference Dong and Lu2007; Dewey et al. Reference Dewey, Quinn, Boschitsch and Smits2014; Quinn et al. Reference Quinn, Moored, Dewey and Smits2014; Zhong et al. Reference Zhong, Han, Moored and Quinn2021a ). A second benefit would be derived when a row of fish downstream of a leading row would swim in between the reverse von Kármán vortex streets shed by the leaders. This would cause the follower fish to experience a reduced flow region in the time average and, consequently, reduced drag – a form of a drafting mechanism. Later, it was appreciated that unsteady vortex–body interactions could also lead to substantial thrust and efficiency enhancements for a follower when its tail motion is synchronised with the induced time-varying flow from an impinging or nearby vortex street (Gopalkrishnan et al. Reference Gopalkrishnan, Triantafyllou, Triantafyllou and Barrett1994; Streitlien et al. Reference Streitlien, Triantafyllou and Triantafyllou1996; Liao et al. Reference Liao, Beal, Lauder and Triantafyllou2003; Liao Reference Liao2004; Beal et al. Reference Beal, Hover, Triantafyllou, Liao and Lauder2006; Akhtar et al. Reference Akhtar, Mittal, Lauder and Drucker2007; Boschitsch et al. Reference Boschitsch, Dewey and Smits2014; Portugal et al. Reference Portugal, Hubel, Fritz, Heese, Trobe, Voelkl, Hailes, Wilson and Usherwood2014; Kurt & Moored Reference Kurt and Moored2018; Li et al. Reference Li, Nagy, Graving, Bak-Coleman, Xie and Couzin2020; Alben Reference Alben2021; Han et al. Reference Han, Pan, Liu and Dong2022) thereby maximising the follower’s effective angle of attack (Muscutt et al. Reference Muscutt, Weymouth and Ganapathisubramani2017). Contrary to these two-dimensional studies of a pair of swimmers, studies examining three-dimensional swimmers in infinite schools show that wake breakdown eliminates coherent vortex–body interactions and leaves only direct body–body interactions to provide a schooling benefit (Daghooghi & Borazjani Reference Daghooghi and Borazjani2015). Yet, three-dimensional fish-like swimmers can experience performance boosting vortex–body interactions when swimming in a pair (Verma et al. Reference Verma, Novati, Koumoutsakos and Sethian2018) or not (Li et al. Reference Li, Kolomenskiy, Liu, Thiria and Godoy-Diana2019) depending upon their kinematics and/or geometry. Understanding the interplay and applicability of schooling mechanisms from body–body interactions and drafting to vortex–body interactions and wake breakdown for three-dimensional swimmers in school sizes beyond a pair of swimmers has yet to be fully understood.
Beyond the direct performance benefits that occur from various schooling mechanisms it is also known that these interactions can give rise to forces that push and pull swimmers to settle into hydrodynamically stable formations, a phenomenon known as the Lighthill conjecture (Ramananarivo et al. Reference Ramananarivo, Fang, Oza, Zhang and Ristroph2016b ). Previous work on unconstrained foils has discovered streamwise stable (Becker et al. Reference Becker, Masoud, Newbolt, Shelley and Ristroph2015; Peng et al. Reference Peng, Huang and Lu2018; Newbolt et al. Reference Newbolt, Zhang and Ristroph2019; Heydari & Kanso Reference Heydari and Kanso2021; Wei et al. Reference Wei, Hu, Li and Shi2023), and both streamwise and cross-stream stable formations of two-dimensional (Lin et al. Reference Lin, Wu, Zhang and Yang2021, Reference Lin, Wu, Yang and Dong2022) and three-dimensional (Ormonde et al. Reference Ormonde, Kurt, Mivehchi and Moored2021) swimmers. Constrained or tethered swimmers have also been used to approximate free-swimming states and to characterise the stability of formations in order to qualify their performance benefits (Li et al. Reference Li, Kolomenskiy, Liu, Thiria and Godoy-Diana2019). To date, our knowledge of the stability of formations as the school size increases beyond a pair of swimmers is quite limited to only a few studies (Wei et al. Reference Wei, Hu, Li and Shi2023; Newbolt et al. Reference Newbolt, Lewis, Bleu, Wu, Mavroyiakoumou, Ramananarivo and Ristroph2024).
In this study, our aim is to detail schooling interactions of three-dimensional swimmers are moderately high (
$\textit{Re}\gt 5$
$000$
), biologically relevant Reynolds numbers, extending the current body of experimental works that only address schools of two swimmers (Ligman et al. Reference Ligman, Lund and Fürth2023). To achieve that, we impose pitching kinematics and a flow speed to the foils that yield a positive thrust force and reverse von Kármán wakes of momentum surplus, as typically observed for oscillatory swimmers (Smits Reference Smits2019). In contrast, foils oscillating at sufficiently low Strouhal numbers produce zero or negative thrust (drag), and an associated von Kármán wake characterised by a momentum deficit, which is not of interest for the study of fish-like swimmers. Inspired by Weihs’ schooling hypothesis (Weihs Reference Weihs1973, Reference Weihs1975) we study a school of two leaders and one follower where we can directly examine the ‘back half’ of the classic diamond formation and probe the associated schooling mechanisms of wake breakdown, drafting, vortex–body interactions and body–body interactions. We also further examine the stability characteristics of the three-foil school.
The paper is organised as follows. Section 2 describes the experimental set-up and methodology used throughout the study. Section 3 details the three-dimensional flow fields generated by the small school. Section 4 probe drafting, vortex–body, and body–body mechanisms for the downstream foil. Section 5 presents the performance of the leaders. Section 6 presents the collective performance of the school and examines its stability. Section 7 provides a discussion of the results and summarises the conclusions of the study.
2. Methods
To study the flow interactions, and their implications, that occur in fish schools we focus on an idealised model of a fish school as three pitching hydrofoils, as adopted previously (Moored & Quinn Reference Moored and Quinn2019). Pitching hydrofoils capture the salient features of unsteady flows of circulatory forces, added mass forces and vortex shedding that are key characteristics underlying the hydrodynamics of fish swimming. These hydrofoils can be viewed as a representation of the caudal fins of fish, which, for many species, are the primary propulsive elements responsible for generating the bulk of a swimmer’s thrust (Smits Reference Smits2019).
2.1. Water-channel experiments
Water-channel experimental set-up. The follower positioning in the horizontal
$x{-}y$
plane is automated using a traverse system. Stereo particle image velocimetry (PIV) is performed at a fixed horizontal plane and the multiple spanwise PIV layers are obtained by traversing all three foils up and down using a vertical traverse system. Force/torque (F/T) sensors are connected to each foil.

The three pitching hydrofoils are machined aluminium with a NACA 0012 cross-section, a chord length of
$c=$
0.0488 m and an aspect ratio of
. The hydrofoils are immersed in a recirculating water channel (
$W \times H \times L= 0.38 \times 0.45 \times 1.52$
m) at the University of Virginia. Figure 1 shows a schematic of the water-channel set-up. The flow speed is maintained at
$U=0.174$
m s−1 using an ultrasonic flow meter (Dynasonics Series TFXB). The foils are constrained at fixed positions during each trial. The foil shaft is connected to a six-axis F/T transducer (ATI Mini40) to measure the hydrodynamic loads experienced by the foil, and an optical encoder (US Digital E2-5000) measures the pitching angle. Transducer and encoder data are acquired through a National Instruments data acquisition card at a sampling rate of 1000 Hz. The foil is driven by a servo motor (Teknic CPM-MCPV-2341S-ELN) coupled to a 5 : 1 gearbox (SureGear PGCN23-0525). The positioning of the three foils is enabled by a four-axis motorised traverse system, in which the follower is traversed in the
$x$
- and
$y$
-directions independently. A sinusoidal oscillatory pitching motion about the leading-edge axis is prescribed to all three swimmers.
We use multi-layer phase-averaged stereoscopic particle image velocimetry to measure the three-dimensional flow field around the pitching foils. Similar techniques have been used previously (King et al. Reference King, Kumar and Green2018; Zhong et al. Reference Zhong, Han, Moored and Quinn2021a
; Zhu & Breuer Reference Zhu and Breuer2023; Zhu et al. Reference Zhu, Liu, Han, Feng, Moored, Zhong and Quinn2025). We seed the water flow using neutrally buoyant 50
$\unicode{x03BC}$
m silver-coated hollow ceramic particles (Potters Industries), and illuminate the particles using two 5 mm-thickness laser sheets firing from each side of the water channel to minimise foil shadows. The laser sheets are created using continuous wave lasers (532 nm, 5W Raypower MGL-W-532 and 10W CNI MGL-W-532A) with laser guiding arms and sheet optics. Two high-speed cameras (Phantom SpeedSense M341,
$2956\times 1877$
pixels) with Scheimpflug adaptors (Dantec Dynamics) and 50 mm lenses (Zeiss) beneath the water channel are used for recording the PIV image pairs. The raw images are processed using Dantec Dynamic Studio 6.9 by an adaptive PIV algorithm (minimum interrogation window
$32\times 32$
pixels, maximum interrogation window
$64\times 64$
pixels). To obtain the three-dimensional (3-D) velocity field, we fix the laser sheet position and traverse the foils in the
$z$
-axis with a step size of 5 mm. For each spanwise layer, we phase average 750 instantaneous two-dimensional, three-component (2D3C) velocity fields over 25 cycles (i.e. 30 vector fields per cycle). Due to symmetry, we measure the flow field of the bottom half of the foils plus 2 cm below the foil tip (20 layers) and then mirror the flow field about the mid-span plane to construct the three-dimensional, three-component (3D3C) velocity field (39 layers). The 3D3C velocity field is reconstructed using MATLAB (R2023a), attaching all of the three-component (3C) velocity fields (U,V,W) extracted from each two-dimensional (2-D) PIV plane obtained for each position along the spanwise direction.
To visualise the 3-D flow structures, we use the
$Q$
-criterion (Jeong & Hussain Reference Jeong and Hussain1995),
$Q=0.5(\Vert \boldsymbol{\varOmega }\Vert ^2-\Vert \boldsymbol{S}\Vert ^2)$
, where
$Q$
is the second invariant of the velocity gradient tensor,
$\boldsymbol{\varOmega }$
is the vorticity tensor and
$\boldsymbol{S}$
is the strain-rate tensor. The flow structures are plotted as iso-
$Q=10$
surfaces coloured by the non-dimensional spanwise vorticity
$\omega _zc/U$
, where
$\omega_z = \partial v/\partial x - \partial u/\partial y$
and
$c$
is the chord length.
Schematics of the arrangement of the hydrofoils in the experimental domain. Dots represent the locations of the follower’s leading edge. The domain is
$1.2\leqslant X^* \leqslant 4.0$
,
$-1.6 \leqslant Y^* \leqslant 1.6$
. Leader 1 and the follower are always in phase, and leader 2 is either in phase with leader 1 (
$\phi = 0$
), or out of phase (
$\phi = \pi$
).

One foil, designated leader 1, is used as the reference with kinematics
$\theta _{1}(t) = \theta _0 \sin (2\pi f t)$
, where
$\theta _{1}(t)$
is the instantaneous pitching angle,
$\theta _0$
is the pitching amplitude and
$f$
is the oscillation frequency. The second foil, denoted as leader 2, is positioned to the side of leader 1 with prescribed kinematics
$\theta _{2}(t) = \theta _0 \sin (2\pi f t + \phi )$
, where
$\phi$
is the synchrony, or phase difference of leader 2 relative to leader 1. Experiments were performed for the leaders swimming in-phase (
$\phi = 0$
) and out-of-phase (
$\phi = \pi$
). The third foil, designated as the follower, oscillates in phase with leader 1 with identical kinematics, i.e.
$\theta _3(t) = \theta _1(t)$
. All three foils have a dimensionless peak-to-peak amplitude of motion of
$A^* = A/c = 0.24$
, with
$A=2 \sin {\theta _0}$
, Strouhal number
$\textit{St} = \textit{fA}/U = 0.29$
, reduced frequency
$k = f c/U = 1.14$
and chord-based Reynolds number
$\textit{Re} =$
8 450. The relative positions of the leaders remain constant throughout the entire study in a side-by-side formation with a dimensionless cross-stream distance of
$D_{2,1}^* = D_{2,1}/c=1.1$
, which is the closest spacing possible for the experimental apparatus. The position of the follower is varied throughout a 2-D grid presented in figure 2, in the ranges
$1.2 \leqslant X^* \leqslant 4.0$
,
$-1.6 \leqslant Y^* \leqslant 1.6$
. A constant cross-stream grid spacing of
$\Delta Y^*=0.2c$
is employed. The streamwise spacing
$\Delta X^* = 0.2$
c is used in the region
$1.2 \leqslant X^* \leqslant 3.0$
. Further downstream data are collected at points
$X^*=[3.25, \ 3.5, \ 4.0]$
. Each blue dot represents a position of the leading edge of the follower where measurements were collected.
Measurements of the thrust (
$T$
), lift (
$L$
) and pitching moment (
$M_z$
), allows us to calculate the thrust (
$C_T$
), lift (
$C_{\!L}$
) and power (
$C_{\!P}$
) coefficients as well as the efficiency (
$\eta$
) for each foil
where
$q_{\textit{dyn}} = ({1}/{2}) \rho U^2$
is the dynamic pressure,
$\rho$
is the fluid density,
$s$
is the foil span and
$U$
is the free-stream speed. The
$\overline {(\boldsymbol{\cdot })}$
symbol indicates the mean values computed from a total of five independent trials, each comprising 80 oscillating cycles. The four quantities defined in (2.1) can also be defined for the collective denoted with a
$C$
subscript
In this study we compare the performance of a foil that is in a school with that of a single, isolated foil. The forces and power coefficients are thus normalised by their equivalent isolated foil quantities, noted with the
$(\boldsymbol{\cdot })^{\textit{iso}}$
superscript. This is carried out to mitigate any systematic bias in the measurements. The normalised quantities are noted with the
$(\boldsymbol{\cdot })^*$
superscript and are defined as
The measured thrust of the isolated foils, are
$C_{T,1}^{\textit{iso}} = 0.18 \pm 0.03$
,
$C_{T,2}^{\textit{iso}} = 0.17 \pm 0.02$
and
$C_{T,3}^{\textit{iso}} = 0.12 \pm 0.01$
, with the intervals consisting of one standard deviation. The root mean square of the fluctuations of the resultant force coefficients produced by the isolated foils are
$F^{\textit{rms}}_{1,\textit{iso}} = 4.8$
,
$F^{\textit{rms}}_{2,{iso}} = 4.4$
and
$F^{\textit{rms}}_{3,\textit{iso}} = 4.0$
, where the coefficient of the resultant force is
$F = \sqrt {T^2 + L^2}/(q_{\textit{dyn}}cs)$
. The values for the static drag coefficients, measured with the foils maintained at a zero angle of attack and averaged over a total of three independent trials of 68 convective time units
$c/U$
are
$C_{D,1}^{\textit{static}} = 0.16 \pm 0.02$
,
$C_{D,2}^{\textit{static}} = 0.12 \pm 0.02$
and
$C_{D,3}^{\textit{static}} = 0.12 \pm 0.01$
. The power coefficients are
$C_{\!P,1}^{\textit{iso}} = 0.989\pm 0.006$
,
$C_{\!P,2}^{\textit{iso}} = 0.836 \pm 0.001$
and
$C_{\!P,3}^{\textit{iso}} = 0.739\pm 0.001$
, and the efficiencies are
$\eta _1^{\textit{iso}} = 0.18\pm 0.03$
,
$\eta _2^{\textit{iso}} = 0.21\pm 0.02$
and
$\eta _3^{\textit{iso}} = 0.17\pm 0.01$
. The coefficients for thrust, drag, power and efficiency reported for the isolated foils were computed from a total of three independent trials of 80 oscillating cycles each. The force and torque transducers are tared before each trial for static foils positioned at their respective location in the domain with the water channel running at the set free-stream speed, yielding net thrust and lift forces, i.e. the forces produced by each foil subtracted from their mean values when held stationary at their respective position within the school. We note that the differences in the measured performance of each foil in isolation are due to experimental errors such as manufacturing imperfections, sensor background noise and stiffness and the effective pitching amplitude of the follower being up to 3 % lower than prescribed for that specific servo motor, causing a reduction in the thrust and power output of the follower compared with the two leaders, for which their respective motors achieved the desired prescribed amplitude of motion.
The average performance metrics for the three-foil school can also be compared with the collective performance of three isolated foils performance to define normalised collective performance metrics as
3. Three-dimensional flow fields
Isosurfaces of
$Q=10$
from phase-averaged PIV data of the two leaders in the absence of the follower. (a,d) Side view. Black dashed lines highlight the spanwise compression of the wake structures and their wake wavelength
$\lambda$
is marked. (b,e) Top view. For the out-of-phase case (e) the iso-surfaces from the two wakes merge together, starting around
$X^*=3.8$
, as the two horizontally oriented portions of the surfaces (light grey colour) generated by the wing tips make contact in region
$C_1$
, eventually merging completely as highlighted by
$C_2$
. (c, f) Perspective view.

We start by presenting the 3-D flow fields generated by the two side-by-side leaders in the absence of the downstream follower. The leaders are positioned side by side with leader 1 located at
$(X^*_1,\,Y^*_1) = (0,\,0.55)$
and leader 2 at
$(X^*_2,\,Y^*_2) = (0,\,-0.55)$
. Figure 3 presents vortex structures that form the wake shed by the side-by-side leaders. All plots shown are constructed using either time-averaged or phase-averaged PIV data for a zero pitching angle (
$\theta = 0$
), which is at the beginning of the pitching cycle
$t^* = t/T = 0$
(
$T = 1/f$
is the cycle period). For the in-phase leaders (figure 3
a–c), both foils are pitching counter-clockwise (
$-\hat {z}$
), and for the out-of-phase leaders, leader 1 is pitching counter-clockwise while leader 2 is pitching clockwise (
$+\hat {z}$
). The 3-D structure of the wake shows the shedding of alternating vortex rings that are intertwined amongst themselves, generated from the tips and trailing edge of the foils once every half-cycle of the pitching motion (Buchholz & Smits Reference Buchholz and Smits2008; King et al. Reference King, Kumar and Green2018). The wake structures narrow along the spanwise dimension
$z$
as they advect downstream, a phenomenon observed for low-aspect-ratio rectangular panels (Buchholz & Smits Reference Buchholz and Smits2008) and trapezoidal panels of aspect ratio 4.17 (King et al. Reference King, Kumar and Green2018). The dashed lines on figures 3(a) and 3(d) highlight this spanwise compression.
Flow field at the midspan plane
$z=0$
for two leaders. (a) and (c) Cycle-average vorticity for the in-phase (a) and out-of-phase (c) cases. Here,
$A_1$
and
$A_2$
highlight regions of more pronounced breakdown of the vortices shed by leader 1 and leader 2, respectively. Contours presented for beginning of pitching cycle
$t/T = 0$
. (b) and (d) Time-average vorticity contours for in-phase (b) and out-of-phase (d) cases. Dashed lines highlight the deflection of the positive (red) and negative (blue) vortices due to the interaction between the two wakes. (e) and ( f) Mean streamwise velocity fields for in-phase (e) and out-of-phase ( f) cases. The dark blue region between the out-of-phase leaders in ( f) show a reduction of up to 17 % in the mean streamwise speed due to the synchrony between the wakes for
$1.1\leqslant X^* \leqslant 2.2$
.

Figures 4(a) and 4(b) show the
$z$
-component of vorticity,
$\omega _z$
, at the midspan of the foils (
$z = 0$
) and at
$t^* = 0$
. A typical reverse von Kármán vortex street is observed up to 2.5 wake wavelengths downstream of the trailing edge (at approximately
$X^*=4$
). This can be observed at any given spanwise location within
$-1\lt Z^*\lt 1$
from the 3-D structures shown in figure 3, with well-defined regions of spanwise vorticity of alternating sign. The wakes thus retain features characteristic of 2-D flows for a relatively long convective time due to the aspect ratio (
$\textit{AR}=3$
) and Strouhal number (
$\textit{St}=0.29$
) chosen for this study. Further downstream, at
$X^*\gt 4$
, the wakes start breaking down and the spanwise-oriented sections of the vortex rings from both hydrofoils unravel, interact and eventually merge, leading to a disruption of the vortical structures observed upstream. The breakdown of coherent vortex structures is even more pronounced for the out-of-phase case, as can be seen in figure 3(d–f). Figure 3(e) shows the coherent vortex structures from both leaders making contact as highlighted by region
$C_1$
, and completely merging downstream in region
$C_2$
, which enhances their breakdown. The mean vorticity contours depicted in figures 4(b) and 4(d) also aid in visualising how the interaction of the two wakes alters their development and advection trajectories. For the in-phase case the clockwise vortices (blue) from leader 1 and counter-clockwise vortices (red) from leader 2 curve towards each other starting at
$X^*\approx 2.8$
, eventually meeting at the symmetry line
$y=0$
at
$X^*=4.5$
, spreading out the wakes of each one of the leaders along the cross-stream dimension. Figures 4(e) and 4( f) present the time-average streamwise velocity
$\overline {u}^*$
in the wake of the leaders at the midspan (
$z=0$
). The flow fields show two distinct jets (coloured red) with mean flow speeds up to 20 % higher than the free stream. The region between the two jets (blue) shows a reduction of the streamwise flow speed, with maximum reductions of 11 % and 17 % for the in-phase and out-of-phase cases, respectively. According to the drafting mechanism proposed by Weihs (Reference Weihs1973), a third swimmer (follower) positioned within this region would benefit from a reduction in the viscous drag force, thus improving its thrust and efficiency. In § 4.1 it is shown that the performance of the follower decreases in this region, contrary to the predictions of the drafting mechanism. In fact, the best thrust benefits for the follower occur in the accelerated-flow region of the jets where there are direct vortex–body interactions that are regulated by the phase synchrony between the pitching motion of the follower and the impinging wake vortices–something not considered in the drafting hypothesis proposed in Weihs (Reference Weihs1973). It will be shown that the hypothesised performance gains originating from viscous drag reduction in the reduced flow region are not observed for a school of three pitching hydrofoils.
Isosurfaces of
$Q=10$
from phase-averaged PIV data of the three-foil school. (a) and (d) Side view of the flow field structures located behind the vertical plane
$y=0$
. The wake of leader 1 impinges directly onto the follower. Black dashed lines show the spanwise compression of the wake structures generated by the follower. (b) and (e) Top view. The structures highlighted by
$S_{L1}$
are generated by leader 1 and greatly deformed by the presence of the follower for both the in-phase (b) and out-of-phase (e) cases. The iso-surfaces from the two wakes merge, starting around
$X^*=4$
for the in-phase case (b) and
$X^*=3$
for the out-of-phase case (e), as the two horizontally oriented portions of the surfaces (light grey colour) generated by the wing tips make contact in region
$C_1$
, eventually merging completely as highlighted by
$C_2$
. (c) and ( f) Perspective view. Only structures behind the vertical plane
$y=0$
are depicted for clarity.

Figure 5 shows the vortex structures measured for the three-foil school. Here the follower is located almost directly behind leader 1 at
$(X^*,\,Y^*) = (2.4,\,0.6)$
for both the in-phase and out-of-phase cases. Figures 5(a), 5(c), 5(d) and 5( f) show the flow structures behind the vertical plane
$y=0$
, i.e. only the structures at
$Y^*\gt 0$
are displayed and the structures in front of the plane (
$Y^*\lt 0$
) are omitted for clarity. Figures 5(a) and 5(d) show leader 1 and the follower from a side view, as seen from their left-hand side. The wake structures from leader 1 can be seen impinging onto the follower, and partially covering it. The 3-D structure of the wake generated by the follower is similar to those generated by the leaders, with alternating vortex rings that show spanwise compression as highlighted by the black dashed lines. A perspective view of such flow features are also shown for the in-phase (figure 5
c) and out-of-phase (figure 5
f) cases. Figures 6(a) and 6(c) show this direct vortex impingement at the midspan (
$z=0$
). We can observe how the positive (red) vorticity within region A forms at the leading edge of the follower, remains attached to its surface, and connects to the downstream vorticity generated from the wake of leader 1. The importance of the synchronisation between the follower kinematics and the wake is discussed in § 4 where the performance of the follower is investigated. Figures 6(b) and 6(d) show the time-average vorticity field for the school, and figures 6(e) and 6( f) show the time-average streamwise velocity field for the in-phase and out-of-phase cases, respectively. The black and white lines downstream of the follower highlight the mean jets generated by the follower. Region B in figure 6(e) highlights a region of around 15 % reduction to the time-average streamwise velocity for the in-phase case along the
$y=0$
vertical plane. Such a reduced-speed region was predicted by Weihs (Weihs Reference Weihs1973) and appears here, indicating that some flow field features predicted by 2-D potential flow theory are still present for the 3-D schooling foils.
4. Hydrodynamic performance of the follower
Flow field at the midspan plane
$z=0$
for the minimal school. (a) and (c) Cycle-average vorticity for the in-phase (a) and out-of-phase (c) cases. Region
$A$
highlights the positive (red) vorticity generated at the left-side surface of the follower merging to the wake vortex (WV) shed by leader 1. Contours presented for beginning of pitching cycle
$t/T = 0$
. (b) and (d) Time-average vorticity contours for in-phase (b) and out-of-phase (d) cases. (e) and ( f) Mean streamwise velocity fields for in-phase (e) and out-of-phase ( f) cases. Black and white lines downstream of the follower highlight the follower mean streamwise jets of increased flow speed. Region B in (e) highlights an area of reduced flow speed along the symmetry plane
$y=0$
for the in-phase case.

Normalised performance of the follower for the in-phase and out-of-phase cases. (a,b) Thrust coefficient
$C^*_{T,3}$
. Black lines are isolines of
$\omega _z c/U = 2.5$
and 10; white lines are isolines of
$\omega _z c/U = -2.5$
and −10. (c,d) Power coefficient
$C^*_{P,3}$
. (e, f): Efficiency
$\eta ^*_{3}$
.

In order to better understand the hydrodynamic mechanisms that occur within the school, we turn our attention to the interactions between the leaders’ wakes and the follower. Figures 7(a) and 7(b) present the normalised thrust coefficient of the follower for all spatial arrangements of the school as described in figure 2. The highest thrust improvement occurs when the follower is within one chord length downstream of the leaders, with values of
$C_{T,3}^* = 2.6$
for the in-phase case and
$C_{T,3}^* = 2.9$
for the out-of-phase case occurring at the same location, downstream of leader 1 at
$(X^*,\, Y^*) = (1.2,\, 0.6)$
. The power coefficient is directly related to the behaviour of the thrust production, with the local minima and maxima of thrust and power nearly coinciding with each other for the two cases studied. Overall, areas of increased thrust also display increased efficiency despite the power expenditure increases (figures 7
c and 7
d). Figures 7(e) and 7( f) show that the follower achieves the highest measured efficiencies in close proximity to an upstream leader with phase synchrony
$\phi =0$
. For the in-phase leaders, the peak follower efficiency is
$\eta ^*_{3,\textit{max}}=2$
and occurs downstream of both leader 1 and leader 2. For the out-of-phase leaders the peak efficiency achieved is the same, but only occurs downstream of leader 1 which is in-phase with the follower. For this schooling case (out-of-phase leaders) the second highest peak in efficiency is shifted approximately
$0.4$
chords downstream of leader 2 since the synchrony between leader 2 and the follower is
$\phi = \pi$
.
4.1. Reduced flow region and the drafting mechanism
With the follower performance presented, we now revisit the mean flow field of the leaders in figures 4(e) and 4( f), where we see a reduced flow region between the vortical wakes of the leaders, as predicted by Weihs (Reference Weihs1973). For example, at
$(X^*,\, Y^*) = (1.6,\, 0)$
the mean flow speed is 10 % lower than the free-stream speed for the in-phase case and 15 % lower for the out-of-phase case. If placed in a steady flow at that reduced speed (in the absence of the time-varying wake structures) the foil is expected to have a 19 % reduction in its viscous drag force for the in-phase case, and a 28 % reduction for the out-of-phase case (assuming a constant drag coefficient and a
$U^2$
drag law), which would lead to a higher net thrust relative to our baseline case–a thrust benefit. However, the follower foil at that location has a 15 % and 4 % thrust penalty compared with an isolated foil for the in-phase and out-of-phase cases, respectively. Figure 7 shows that the performance of the follower is not improved along the centre line (
$Y^* = 0$
) in precisely the locations where the follower is postulated to gain a benefit from drafting (Weihs Reference Weihs1973). In fact, the data from the three-foil school show that the performance benefits are observed in the accelerated jet flow regions directly behind the leaders where the drafting hypothesis would predict the greatest performance penalties to occur.
While drafting is a widely accepted mechanism for drag reduction, our results indicate that the vortex interactions outweigh the viscous effects from drafting in this case. Previous studies of in-line formations show that vortex interactions can increase the thrust of a downstream follower even when placed in the leader’s accelerated-flow wake that increases its time-average viscous drag (Boschitsch et al. Reference Boschitsch, Dewey and Smits2014; Kurt & Moored Reference Kurt and Moored2018; Ormonde et al. Reference Ormonde, Kurt, Mivehchi and Moored2024). Now, our data show that the time-varying effects of the vortex interactions also eclipse the expected thrust benefits due to reduced viscous drag in a reduced-flow wake. The drafting mechanism thus cannot be assumed to improve the performance of oscillating propulsors without explicitly accounting for the significant impacts of the unsteady vortex interactions, as will be shown in the next section. It is worth noting that the observed behaviour may change for leaders placed at different lateral distances or for different pitching amplitudes
$A^*$
and Strouhal numbers
$\textit{St}$
since these parameters will directly affect the width of the wakes produced by the leaders such that the interaction of the follower with the wakes at the centre line
$Y^*=0$
may change significantly. Additionally, pitching hydrofoils are simple models of only a propulsive fin of a fish, whereas fish have a propulsive fin and a body that mostly generates drag (Smits Reference Smits2019). With a larger viscous drag source, it might still be possible that fish-like swimmers extract performance benefits from drafting within this reduced flow region. While the present results indicate that the drafting mechanism may not be as important as previously hypothesised, this is still an open question that deserves further investigation. Finally, we also stress that our experimental methodology does not allow for a direct measurement of the skin-friction drag force of a pitching foil, which should be possible to quantify through viscous simulations, instead.
4.2. Vortex–body mechanism
The thrust coefficient map (figure 7) shows a very distinctive spatial arrangement with regions of high thrust located within the accelerated jet flow behind the leaders. Specifically, these regions show a subtle double peak structure with the peaks located along the paths of the wake vortices shed by the leaders, represented by the black (white) isolines of positive (negative) time-average vorticity from measurements of just the two leaders alone. The thrust production is thus maximally enhanced by a direct impingement of the oncoming vortices onto the follower as observed in previous work (Ormonde et al. Reference Ormonde, Kurt, Mivehchi and Moored2024). The spatial periodicity of these high-thrust regions coincides with the wake wavelength
$\lambda$
of the leader’s wake shown in figures 4(a) and 4(b), in agreement with previous findings for 2-D and 3-D foils located directly downstream of a single in-line fore foil (Kurt & Moored Reference Kurt and Moored2018). The direct vortex impingement changes the instantaneous effective angle of attack of the follower, which in turn changes its cycle-averaged performance. Figure 7 indicates that the same mechanism exists for the wakes shed by the two interacting leaders. Following a direction perpendicular to the vortex trajectories we also see a decay in the thrust of the follower and a decay as the follower is placed further downstream of the leaders.
The synchrony between the vortex impingement timing and the foil motion dictates where the peaks are located in the 2-D plane. Figure 8 shows the follower interacting with the wake of leader 2 at the spatial location
$(X^*,\, Y^*) = (1.8,\, -0.6)$
, for the in-phase case (figure 8
a) and the out-of-phase case (figure 8
b). The change in the phase-averaged thrust caused by the vortex–foil interaction is plotted as
$\Delta C^*_T$
, being the difference between the schooling foil and an isolated foil, normalised by the isolated foil as
Figure 8(a) shows that for the in-phase case at instant (1) the WV is at the leading edge, moving towards the lower surface at the beginning of the cycle (downstroke). It then interacts with the leading-edge vortex (LEV) at instant (2) and induces a downwards velocity that detaches the LEV from the foil surface during the upstroke at instants (3) and (4). The time-varying thrust is modestly affected by the vortex interaction, and the follower’s mean thrust is somewhat elevated from an isolated foil with
$C_{T,3} = 1.15$
, however, this is lower than the maximum achievable mean thrust at that spatial location. This interaction also leads to a 14 % reduction in power.
Wake–foil interaction and its impact on the phase-average thrust production for the follower foil at
$(X^*,Y^*) = (1.8,-0.6)$
.:(a) In-phase case
$\varPhi _{\textit{spatial}}=3.81$
rad. The interaction of the WV with the LEV causes it to prematurely detach from the foil surface during the cycle; (b) Out-of-phase case
$\varPhi _{\textit{spatial}}=0.67$
rad. The WV interacts with the LEV such that it remains attached to the foil surface up until it merges with the vorticity shedding from the trailing edge, increasing the thrust generation. Non-dimensional time at snapshots: (1)
$t/T = 1/30$
; (2)
$t/T = 8/30$
; (3)
$t/T = 18/30$
; (4)
$t/T = 24/30$
.

In figure 8(b) the follower is at the same spatial location, but is out-of-phase with leader 2, which, consequently, changes the synchrony between the foil motion and the impinging vortex. At instant (1) the WV impinges onto the follower at the beginning of the cycle (downstroke). The LEV remains attached to the upper surface throughout all of the downstroke, and remains attached during the earlier portion of the upstroke at instant (2) and later in the upstroke at instant (3). At instant (4) the LEV merges with the trailing-edge vortex and detaches from the upper surface. Snapshots (3,4) show the next WV (blue) interacting with the next LEV (red) on the lower surface during the upstroke. The time-varying thrust is observed to be significantly enhanced throughout most of the cycle and the mean thrust,
$C^*_{T,3} = 2.0$
, is twice the thrust of an isolated foil while the mean power also increases by 20 %.
Descriptions of the synchrony between an impinging wake and a downstream foil have been presented previously (Portugal et al. Reference Portugal, Hubel, Fritz, Heese, Trobe, Voelkl, Hailes, Wilson and Usherwood2014; Kurt & Moored Reference Kurt and Moored2018; Lee et al. Reference Lee, Hrynuk and Moored2023) and described in terms of a spatial phase defined as
$\varPhi _U=2\pi f\,X/U + \phi$
, taking the free-stream (or swimming) speed as the advection velocity of the wake vortices. We propose an improved definition for the spatial phase utilising the measured wake wavelength
$\lambda$
instead of the estimated value of
$\lambda _U = U/f$
based on the free-stream speed. An average wake wavelength is measured for the two leaders in the absence of a follower. The measured wavelength is defined as the mean distance along the streamwise direction between the first and third centroids of vorticity within successive triplets of shed vortices, as shown in figure 9(a).
Relationship between the wake dynamics of leader 1 and the thrust of a follower foil directly downstream of leader 1, at
$(X^*,Y^*)=([1.1 ,\ 4],0.6)$
. Results for in-phase leaders. (a) Measured wake wavelength of the leaders in the absence of a follower, defined as the average distance spanning each set of three consecutive centroids in the streamwise direction. The vortex boundary is defined by the iso-line of constant vorticity (black line)
$|\omega |=0.02\max (\omega )$
. Snapshot shown for
$t/T = 0$
. (b) Thrust map of the follower foil highlighting streamwise locations of local maxima, approximately spaced by
$\lambda ^*$
. (c) Streamwise decay of vorticity within each identifiable vortex for all instants of the flapping cycle
$t/T=[0,1]$
.
$\|\omega /\max (\omega ) \|= C_1/(X^*+C_2)$
, where
$C_1=1$
and
$C_2=-0.18$
. (d) Predictive model of the follower’s thrust when directly downstream of the leader
$C^*_{T,3} = 1 + C_1/(X^*+C_2)(1+ \sin {\varPhi })$
, where
$\varPhi =2\pi X^*/\lambda ^*+\phi$
. The vertical dotted lines illustrate how the measured wake wavelength
$\lambda ^*$
in (a) matches the streamwise spacing of local peaks of thrust in (b) and (d).

The proposed spatial phase
$\varPhi$
is thus defined as
where
$\lambda ^*=\lambda /c$
is the dimensionless measured wake wavelength, and
$\phi$
is the temporal synchrony between the leader whose wake is impinging on the follower. Our data show that the measured wavelength is 29 % larger than the estimated
$\lambda _U$
based on the free-stream speed. This discrepancy is expected due to the fact that the pitching foils produce a net thrust force, thus creating a jet of downstream accelerated flow, as shown in figures 4(e) and 4( f). This region of accelerated flow directly impacts the advection speed of the wake vortices, stretching it along the streamwise direction compared with the free-stream-based wavelength prediction.
Our goal now is to connect the properties of the reverse von Kármán wake shed by the leaders, without the follower present, to the performance of the follower foil when it is located directly downstream of them. Figure 9(b) shows the streamwise distance between thrust peaks of the follower, which are approximately separated by a distance of
$\lambda ^*$
. A simple model of the follower’s thrust is postulated to be proportional to a periodic function of the measured spatial phase
$\varPhi$
such as
$C^*_{T,3} \propto \sin { (\varPhi )}$
and to have a streamwise decay directly tied to the decay of vorticity in the wake. The streamwise decay of vorticity in the wake of the leaders (in the absence of a downstream follower) is due to the spreading of vorticity as the vortices advect downstream and is shown in figure 9(c). The average vorticity contained within each vortex is measured for all instants
$t/T=[0,1]$
as a function of the vortex centroid position. The measured magnitude of positive (red) and negative (blue) vortices decay downstream, and is well described by an inverse function of the form
$\|\omega /\max (\omega )\|=C_1/(X^*+C_2)$
, where
$C_1=1$
and
$C_2=-0.18$
. The hyperbolic decay used in our model is inspired by the decay law predicted by Lamb’s vortex model as described by Ponta (Reference Ponta2010), which was observed in their measurements of a von Kármán cylinder wake. Our proposed model for the thrust of the follower is given by
with
$C_1 = 1$
,
$C_2=-0.18$
and
$\varPhi = 2\pi f/1.29U X^* + \phi$
for the follower directly behind one of the leaders. The term
$2\pi f/1.29 U$
accounts for the measured wake wavelength downstream of leader 1. The model successfully predicts the streamwise behaviour of the thrust coefficient of the follower as shown in figure 9(d). The spatial oscillations are primarily due to the spatial phase
$\varPhi$
, which governs how the WV interacts with the follower’s LEV as depicted in figure 8, while the streamwise wake strength decay affects the effective angle of attack of the follower foil and is a function of the distance
$X^*$
from the follower to the trailing edge of the leader.
The spatial phase has been previously defined throughout literature using the estimated (approximate) wake wavelength involving the free-stream speed
$U$
and the downstream distance
$x$
(Portugal et al. Reference Portugal, Hubel, Fritz, Heese, Trobe, Voelkl, Hailes, Wilson and Usherwood2014; Newbolt et al. Reference Newbolt, Zhang and Ristroph2019; Li et al. Reference Li, Nagy, Graving, Bak-Coleman, Xie and Couzin2020; Heydari & Kanso Reference Heydari and Kanso2021). This is a suitable approximation as long as the wake follows a straight path in the
$x$
-direction and its advection speed is close to
$U$
in magnitude and remains constant as it advects downstream. For the wakes presented here, however, these assumptions do not hold true. The spatial phase parameter has also been used to study arrays of tandem oscillating foils for energy harvesting. A similar correction of the wake advection speed has been shown to improve the prediction of the optimal spacing for the performance of the downstream foil (Ribeiro et al. Reference Ribeiro, Su, Guillaumin, Breuer and Franck2021; Handy-Cardenas et al. Reference Handy-Cardenas, Zhu and Breuer2025) compared with predictions based on an uncorrected wake advection velocity equal to the free-stream speed (Kinsey & Dumas Reference Kinsey and Dumas2012). Unlike the wakes discussed in this study, which produce regions of mean accelerated flow (see figures 4
c and 4
f) and a reverse von Kármán vortex street, the wakes generated by drag-producing energy-harvesting foils exhibit a mean flow deficit and different vortical structures. Despite this difference, the spatial phase parameter encapsulates the combined effects of temporal synchronisation and spatial separation between two tandem foils into a single parameter. The spatial phase describes the timing between the downstream foil kinematics and the impingement of the wake shed by the upstream foil.
The model’s insight is that only the wake information from the leaders alone is needed to accurately predict the optimal vortex–body phase synchronisation of downstream followers given a von Kármán or reverse von Kármán wake of well-defined vortices sufficiently distinct from one another. Therefore a school with optimal phase synchronisation of any size can be designed by progressively marching from leader to follower in the following way. First, the wake wavelength for a leader’s wake can be calculated/measured and then the optimal phase of a follower placed in this flow can be chosen. Then the wake wavelength from the leader–follower pair can be calculated/measured and another follower can be placed downstream with an optimal synchronisation, and so forth. For significantly different wake topologies other than the von Kármán and reverse von Kármán types, future work is required to determine the validity of the model proposed in (4.3). The vortex–body interactions observed in figure 8 are not to be expected for substantially different wakes, which should also directly affect the impact of the interaction onto the foil’s performance.
5. Hydrodynamic performance of the leaders
Performance landscape of leader 1 and leader 2 as a function of the position of the follower. The black-coloured foil indicates the foil to which the figure corresponds to. The small grey dots indicate the grid points where data were collected. Thrust and power coefficients are normalised by the performance of each respective leader in the side-by-side arrangement in the absence of a follower. (a,b) Thrust of leader 1 and leader 2, respectively, for an in-phase synchrony. (c,d) Thrust of leader 1 and leader 2, respectively, for an out-of-phase synchrony. (e, f) Power of leader 1 and leader 2, respectively, for an in-phase synchrony. In ( f), points A and B are equidistant from leader 2. (g,h) Power of leader 1 and leader 2, respectively, for an out-of-phase synchrony. In (h), points A and B are equidistant from leader 2. Measurement uncertainty is found in figure 13 in the Appendix A.

The performance of the leaders is affected by the follower when it is within approximately one chord length away from either leader 1 or leader 2. Figure 10 shows the thrust and power maps for leader 1 and leader 2 as a function of the follower’s position. The effect of a downstream follower on a single upstream foil has been extensively documented in the literature (Boschitsch et al. Reference Boschitsch, Dewey and Smits2014; Kurt & Moored Reference Kurt and Moored2018) and is also observed in our results for two side-by-side leaders. Figure 10(a) shows, for the in-phase case, that leader 1 produces higher thrust when the follower is directly interacting with its wake. The same is true for leader 2 for the in-phase case shown in figure 10(b), denoted by the deep red region directly behind leader 2 in the figure. For the out-of-phase case, leader 1 produces higher thrust when the follower is in close proximity (figure 10 c). However, for the out-of-phase case, figure 10(d) shows no changes in the thrust of leader 2 when the follower is in its close vicinity.
The power coefficient of the leaders as a function of the position of the follower is presented in figure 10(e–f). For the in-phase case, figure 10(e) shows that leader 1 experiences a 10 % increase in its power expenditure when the follower is directly behind it (deep red region at
$1.2\leqslant X^*\leqslant 1.8$
and
$Y^*=D^*/2$
). Conversely, when the follower is behind of leader 2, leader 1 experiences a 5 %–10 % decrease in its power expenditure, as shown by the deep blue region at
$1.2\leqslant X^* \leqslant 1.4$
and
$Y^*=-D^*/2$
. Figure 10( f) shows a symmetric pattern in the power coefficient map of leader 2: the power expenditure of leader 2 increases when the follower is closely behind it (deep red region at
$1.2 \leqslant X^* \leqslant 1.8$
and
$Y^*=-D^*/2$
), and decreases when the follower is closely behind of leader 1 (deep blue region at
$1.2\leqslant X^* \leqslant 1.4$
and
$Y^*= D^*/2$
).
The performance maps for the thrust and power of the two leaders shown in figure 10 reveal how performance changes can propagate upstream along the school for sufficiently compact formations. The standard deviation of the scalar fields shown in figure 10 is presented in the Appendix A. They also show that the temporal synchronisation between foils is an important parameter dictating the impact that the downstream foil has on its upstream neighbours, something that has been previously established for in-line formations (Ramananarivo et al. Reference Ramananarivo, Fang, Oza, Zhang and Ristroph2016a ; Kurt & Moored Reference Kurt and Moored2018) and has now been experimentally confirmed for a larger sized school of three foils and several different spatial formations. The data also reveal indirect, or ‘cascading’ effects where one foil affects their nearest upstream neighbour (second foil), which in turn affects a third foil. Figures 10( f) and 10(h) show this phenomenon for the power coefficient of leader 2. Points A and B are equidistant from the trailing edge of leader 2. When the follower is positioned at point B (pink circle), below leader 2, there is no significant change in the power coefficient of leader 2 (black foil). However, when the follower is placed at point A, above leader 2, the data show a 5 % decrease in the power coefficient of leader 2 for the in-phase case (figure 10 f) and a 10 % increase for the out-of-phase case (figure 10 h). We hypothesise that this is caused by both a weak direct effect from the follower on leader 2 and an indirect effect where the proximity of the follower to leader 1 alters the bound circulation of leader 1, which in turn affects its interaction with leader 2. So, in essence, the upstream effect of the follower onto leader 1 is propagated to leader 2. Further studies investigating this phenomenon, especially with more foils interacting, are needed to clarify this hypothesis.
6. Collective performance, fluid mediated forces, and formation stability
We now shift our focus to the collective performance and stability of the school. On the previous sections, we established how each foil is affected by the school’s spatial formation and the temporal synchronisation of the three foils. Here, we present maps of the collective performance of the school and also analyse the stability of those formations. We then discuss the relationship between high-performance formations and school stability.
6.1. Collective performance
Normalised collective performance for the in-phase and out-of-phase cases. (a,b) Thrust coefficient
$C^*_{T,C}$
. (c,d) Power coefficient
$C^*_{P,C}$
. (e, f) Efficiency
$\eta ^*_{C}$
.

The collective thrust, power and efficiency of the school are presented in figure 11 for the in-phase and out-of-phase cases. High collective performance improvements are achieved for both cases, with the peak collective normalised thrust of
$C^*_{T,C} = 1.34$
at
$(X^*,\,Y^*) = (1.2,\,0.6)$
for the in-phase case and
$C^*_{T,C} = 1.58$
at
$(X^*,\,Y^*) = (1.2,\,0.6)$
for the out-of-phase case. The out-of-phase case shows thrust enhancement for almost all the formations investigated since the two leaders by themselves already experience a thrust improvement. For the in-phase case this is not observed, and the minimum thrust value achieved is
$C^*_{T,C} = 0.8$
. This number is largely due to the baseline thrust produced by the two in-phase leaders. In the absence of a follower the thrust coefficient of two in-phase leaders is 17 % lower than the baseline case for an isolated swimmer.
The streamwise spatial distribution of power and efficiency peaks coincide with the streamwise locations of peak thrust. While the power requirements increase for higher-thrust formations, figure 11(e–f) shows that the normalised efficiency is still greater than one, indicating an efficiency enhancement over an isolated foil. Overall, high-thrust formations yield higher propulsive efficiencies. Peak efficiencies of
$\eta ^*_C = 1.24$
are achieved at
$(X^*,\,Y^*) = (1.4,\,0.4)$
for the in-phase case (figure 11
e) and
$\eta ^*_C = 1.21$
at
$(X^*,\,Y^*) = (1.8,\,-0.8)$
for the out-of-phase case (figure 11
f). The highest propulsive performances are found for compact arrangements, where the follower is positioned within one chord length of the trailing edge of either of the leaders.
6.2. Lateral forces and stability considerations
Figures 12(a) and 12(b) show the time-average lift coefficient of the follower as a function of its position for the in-phase and out-of-phase cases, respectively. For the in-phase case, the follower produces small lift forces of less than
$\left |C_{\!L} \right | \lt 0.05$
within a relatively narrow region, highlighted by the white isolines in figure 12(a). Higher lift forces are produced outside of this region. At
$Y^* \gt D^*/2$
, regions of positive (upwards) lift up to
$C_{\!L} \approx 0.45$
are observed, and regions of negative (downwards) lift of the same magnitude can be found at
$Y^* \lt -D^*/2$
. Based on the follower’s lift coefficient map, an unconstrained oscillating foil located anywhere in the domain at
$Y^*\gt D/2$
or
$Y^*\lt -D/2$
will experience a net lateral force that pushes it away from the two upstream leaders. As a result, the follower is required to alter its kinematics in order to maintain its lateral position relative to the leaders, counterbalancing the flow field effects caused by its upstream neighbours. On the other hand, if the follower stays within the region bounded by the white iso-lines, it will experience very small lift forces. As a consequence, significantly less control input is required from a follower to maintain its lateral position within that region. Figure 12(b) shows the same map for the out-of-phase case. Similarly to the in-phase case, high-magnitude lift forces of up to
$\left |C_{\!L} \right | \approx 0.45$
are produced by the follower in the two regions of
$\left |Y^*\right |\gt D^*/2$
. However, in-between the two leaders at
$\left |Y^*\right |\lt D^*/2$
, the lift forces only stay below the 0.05 threshold within a very small, sinuous region that varies quickly along the streamwise direction, as highlighted by the white isolines. As a result, figures 12(a) and 12(b) suggest that, compared with the in-phase leaders that produce a region of low-magnitude lift forces, the two out-of-phase leaders require significantly more input control/corrections from the follower to maintain a given spatial formation downstream of its upstream neighbours.
(a) and (b) Time-average lift coefficient of the follower
$C_{\!L,3}$
for the in-phase and out-of-phase cases, respectively. The white lines are iso-lines of
$|C_{\!L,3} | = 0.05$
. (c) and (d) Time-average vector fields of a schooling follower compared with the isolated case
$\boldsymbol{F}_{\Delta }(x,y)$
for the in-phase and out-of-phase cases, respectively. Vectors represent the direction and magnitude of the time-average force produced by the follower at each location
$\boldsymbol{F}_{\Delta } = ( C_{T,3} - C_{T,\textit{iso}})\hat {x} + ( C_{\!L,3} - C_{\!L,\textit{iso}})\hat {y}$
. The black solid lines are iso-lines of
$|F_{\Delta }| = 0.05$
. (e) and ( f) Root mean square fluctuations of the vector force of the follower normalised by that of an isolated follower
$F_{3}^{\textit{rms}}/F_{3,\textit{iso}}^{\textit{rms}}$
for the in-phase and out-of-phase cases, respectively.

Next, in figures 12(c) and 12(d) we combine the lift and thrust force coefficients of the schooling follower into a resultant vector and compare it with a follower in isolation. The resulting vector field describes the change in the resultant force coefficient produced by the schooling follower relative to the isolated case. The resultant force coefficient vector
$ \boldsymbol{F}_{\Delta }(x,y)$
is non-dimensionalised by
$1/2\rho U^2cs$
and is defined as
where
$C_{T,\textit{iso}}$
is the isolated thrust coefficient and
$C_{\!L,\textit{iso}}$
is in principle zero, given its symmetric motion. A small lift bias
$C_{\!L,3}^{\textit{iso}} = -0.04$
was measured for the isolated follower, however, due to small foil misalignments and other unknown experimental errors.
Figures 12(c) and 12(d) show the force coefficient vector fields for the in-phase and out-of-phase cases, respectively, where the background contour is the magnitude of the non-dimensional force vector
$|\boldsymbol{F}_{\Delta }(x,y)|$
, used to help visualise the magnitude of the magenta-coloured vectors. Isolines of value
$\left | F_{\Delta } \right | = 0.05$
are also used to indicate regions of low change in the total force coefficient. It becomes apparent, from the vector fields, that the schooling follower experiences substantially larger changes in its force production in the lateral direction (lift) compared with the streamwise direction (thrust/drag). Figure 12(a–d) suggests that, in order to stay close to the two side-by-side leaders at the regions
$\left | Y^* \right | \gt D^*/2$
, a follower is required to alter its kinematics compared with an isolated foil. Otherwise, the downstream foil will be laterally pushed away from the leaders until it is effectively no longer schooling, i.e. it is no longer hydrodynamically interacting with its neighbours. In order to produce the required time-averaged zero lift required to maintain group cohesion, a self-propelling follower could possibly alter its mean angle of attack (or yaw angle), adopt an asymmetric oscillatory motion, or in the case of a fish or a fish-like robot, employ their control surfaces such as pectoral/dorsal fins to counteract the lift force produced by the caudal fin. In addition, Figures 12(c) and 12(d) clearly highlight that formations of high collective performance (see figure 11) are not stable, since the higher thrust production of the follower is associated with large lift forces. As a consequence, in order for a freely swimming follower to extract large performance benefits from the wake of its upstream neighbours, it is required to counter-balance these lift forces at the expense of some hydrodynamic performance, which has not been accounted for in this study.
The root mean square (r.m.s.) of the resultant force fluctuations for the follower normalised by that of an isolated foil
$F^{\textit{rms}}_{3}/F^{\textit{rms}}_{3,\textit{iso}}$
are presented in figures 12(e) and 12( f). The thrust and lift fluctuations are calculated by subtracting the mean values from their corresponding time-varying signals
$\tilde {T} = T(t) - \overline {T}$
and lift
$\tilde {L} = L(t) - \overline {L}$
. The resultant force is the vector
$\tilde {\mathbf{F}} = \tilde {T}\hat {x} + \tilde {L}\hat {y}$
, so the r.m.s. of the resultant fluctuations is calculated as
$F^{\textit{rms}} = \sqrt ({{1}/{N}) \sum ( \tilde {T}^2 + \tilde {L}^2 ) }$
, where
$N$
is the number of samples of the time series and the summation is performed over the entire length of the measured signal of 80 oscillating cycles. The r.m.s. values are then normalised by the r.m.s. of the follower foil in isolation, measured to be
$F_{3,\textit{iso}}^{\textit{rms}}/(0.5 \rho U^2cs)$
= 4.0 in non-dimensional form, or
$F_{\textit{iso}}^{\textit{rms}} = 0.43$
N. There is a direct correlation between the fluctuations of figures 12(e) and 12( f), and the follower’s coefficients of thrust and power as shown in figure 7(a–d). The similarities with the power coefficient are particularly evident, since regions of high and low
$C_{\!P}$
coincide almost exactly with the spatial distribution of the fluctuations for both synchrony cases. The same overlap cannot be said about the spatial distribution of the lateral forces, since the peaks of lift forces from figures 12(a) and 12(d) do not overlap with the peaks of force fluctuations. The data indicate that the force fluctuations of the follower directly impact its power expenditure, but are not directly correlated to large mean lateral forces. High mean thrust formations, on the other hand, seem to correlate with increased force fluctuations.
Consider an unconstrained follower swimming downstream of the two leaders. In order for it to gain energetic benefits from the wake of the leaders, the constrained performance data presented suggest it will have to exert some effort to remain at a high-performance position. For example, in the out-of-phase case the follower must counteract a mean lateral force that will push it away from the
$Y^* = 0$
line in between the two leaders. The large force fluctuations detected in the high-performance regions also mean that the follower would likely experience large amplitude recoil motions. Some regions provide improved thrust and efficiency with increased fluctuations, while others are more stable and should require little control of a follower to maintain its position. The force maps from figure 12, together with the performance maps from figure 7, thus provide a framework to navigate through different spatial formations downstream of the leaders that generally offer a trade-off between efficiency and stability.
7. Conclusions
The present study examines the performance and flow structures of three, 3-D pitching hydrofoils across a range of formations that examine the ‘back half’ of the classical diamond formation with both in-phase and out-of-phase synchronisation between the leaders. The mechanisms of drafting, vortex–body interactions, body-to-body interactions and wake breakdown were probed for their role in driving the performance of the collective.
The drafting hypothesis proposed by Weihs (Weihs Reference Weihs1973, Reference Weihs1975) was shown to play no role in the performance benefits of the three-foil school. Reduced flow regions of up to a 17 % reduction in the free-stream speed were identified between the rows of vortices shed from the leaders, as predicted by Weihs. However, the reduced flow regions showed no accompanying regions of increased follower net thrust (reduced drag) when the follower was placed there. In fact, the highest thrust benefits of the follower (190 % more thrust than an isolated foil) occurred directly in the accelerated jet flow regions within the rows of vortices shed by the leaders, where the drafting theory would predict a drag increase and thus a thrust penalty.
The follower foil experienced the largest performance variations of all three foils, achieving high performance for formations where the wake of a leader directly impinges onto it. A new spatial phase definition was proposed to identify the optimal synchrony between the foil oscillation and the oncoming wake, as well as the downstream decay in thrust with increasing distance between the follower and the upstream leader. The key features of the model are the use of the actual wake wavelength instead of the estimated wavelength previously used in literature, as well as the wake’s streamwise vorticity decay. The optimal spatial phase that maximises the efficiency represents the synchronisation between the foil’s motion and an impinging vortex such that the oncoming vortex acts to maximise the angle of attack of the follower. A LEV is also produced by the follower, which pairs with the impinging vortex to propagate along the foil and merge with the forming trailing-edge vortex.
The follower foil also affects the upstream leader foils when it is sufficiently close to them. It is observed that not only can the follower affect its nearest neighbour, but it triggers a ‘cascading’ effect where the nearest leader affects the other leader as well, a phenomenon not typically reported for experiments with oscillatory foils. Additionally, the thrust data for the leaders suggest that the side-by-side formation is more sensitive to a downstream follower than a single upstream swimmer. We hypothesise that this occurs due to the fact that, when the follower disrupts the wake from one leader, the wake from the second leader is also affected, slightly altering the thrust production of both leaders. Future studies performed at different Reynolds and Strouhal numbers, with different foil profiles and lateral spacings between the two upstream leaders, are still required to better evaluate the importance of each of the hydrodynamic mechanisms identified in this study as contributing to the performance of schooling foils. Extending these parameter spaces will thus allow for a better understanding of the relative significance of the drafting mechanism, body–body interactions and vortex–body interactions. Additionally, it will be valuable to test the validity of the follower’s thrust model proposed in (4.3) for different wake wavelengths relative to the foil chord length.
The collective performance of the school shows high maximum thrust and efficiency gains for formations where the follower is in close proximity to the leaders and when the follower is in the wake of either of the two leaders. Collective thrust and efficiency enhancements of up to 58 % and 24 % higher than an isolated foil, respectively, are observed.
The lateral forces of the follower were also examined. The magnitude of the changes in the lift forces are significantly larger than the changes in thrust, and the regions of highest thrust coincide with high net lift forces. A narrow region of reduced total lift and thrust changes were found for both in-phase and out-of-phase leaders. This suggests that schooling fish that want to remain in a certain formation may be faced with a dilemma: regions of high performance may require strong lateral control measures and constant trajectory corrections, whereas regions of relative ‘calm’, or low lateral forces, are not as energetically beneficial.
In summary, it was discovered that (a) drafting provided no performance benefits despite the presence of reduced flow regions, (b) wake breakdown did not eliminate performance benefits from vortex–body mechanisms within three chord lengths downstream of the leaders, (c) vortex–body interactions drove the performance benefits of the follower and collective, (d) body-to-body interactions played a role for compact formations with the follower less than one chord length downstream of the leader and (e) the highest collective performance benefits occurred for compact formations.
Funding
This work was funded by the Office of Naval Research under Program Director Dr R. Brizzolara on MURI grant number N00014-22-1-2616.
Declaration of interests
The authors report no conflict of interest.
Data availability statement
The data that support the findings of this study are openly available at https://doi.org/10.18130/V3/UL6CJE.
Appendix A. Measurement uncertainty of the performance of the leaders
Figure 13 shows the measurement uncertainty of the performance of the leaders. The 95 % confidence interval standard deviation is plotted for the thrust coefficient
$\sigma _{C_T}$
of the in-phase leaders 1 and 2 are shown in figures 13(a) and 13(b), respectively. Figures 13(c) and 13(d) show the uncertainty for the out-of-phase leaders 1 and 2, respectively. Figures 13(e) and 13( f) present the standard deviation for the power coefficient
$\sigma _{C_{\!P}}$
for the in-phase leaders 1 and 2, respectively, and figures 13(g) and 13(h) present the same data for the out-of-phase leaders 1 and 2, respectively. It can be observed that the power coefficient data show very small uncertainty levels, whereas the thrust data of the leaders are substantially more noisy.
Measurement uncertainty of the performance of the leaders.


x−y
1.2⩽X∗⩽4.0
−1.6⩽Y∗⩽1.6
ϕ=0
ϕ=π
Q=10
λ
X∗=3.8
C1
C2
z=0
A1
A2
t/T=0
1.1⩽X∗⩽2.2
Q=10
y=0
SL1
X∗=4
X∗=3
C1
C2
y=0
z=0
A
t/T=0
y=0
CT,3∗
ωzc/U=2.5
ωzc/U=−2.5
CP,3∗
η3∗
(X∗,Y∗)=(1.8,−0.6)
Φspatial=3.81
Φspatial=0.67
t/T=1/30
t/T=8/30
t/T=18/30
t/T=24/30
(X∗,Y∗)=([1.1, 4],0.6)
|ω|=0.02max(ω)
t/T=0
λ∗
t/T=[0,1]
‖ω/max(ω)‖=C1/(X∗+C2)
C1=1
C2=−0.18
CT,3∗=1+C1/(X∗+C2)(1+sinΦ)
Φ=2πX∗/λ∗+ϕ
λ∗
CT,C∗
CP,C∗
ηC∗
CL,3
|CL,3|=0.05
FΔ(x,y)
FΔ=(CT,3−CT,iso)x^+(CL,3−CL,iso)y^
|FΔ|=0.05
F3rms/F3,isorms