Leading-edge noise is a complex phenomenon that occurs when a turbulent fluid encounters a solid object, and is a notable concern in various engineering applications. This study enhances a mathematical leading-edge noise model (Hales et al., J. Fluid Mech., vol. 970, 2023, A29) for anisotropic flow and porous boundaries. The model has two key components. First, we adjust the velocity spectrum to account for the possibility of anisotropy in the flow. This paper rigorously introduces a third dimension for the turbulence spectrum that preserves the turbulence kinetic energy and mathematical definitions for integral length scales. Second, we adapt the fully analytical acoustic transfer function to account for different boundaries by implementing convective impedance boundary conditions when formulating the gust-diffraction problem. This problem is then solved using the Wiener–Hopf technique. We discuss important aspects of this method, including the factorisation of a non-trivial scalar kernel function and the application of suitable edge conditions for the problem. Each modification is inspired by experimental leading-edge noise data using a series of different porous leading edges and anisotropic turbulence generated by a cylinder upstream of the edge. Experimental data demonstrate the interplay between anisotropy and leading-edge modifications while achieving the characteristic mid-frequency noise reduction expected from porous leading edges. Our model is adapted to best fit the trends of the data via a tailored impedance function, leading to good agreement with all datasets across an extended frequency range. This tailored function is used to successfully validate the model against other datasets from a different set of experiments.

Ice-crystal icing (ICI) in aircraft engines is a major threat to flight safety. Due to the complex thermodynamic and phase-change conditions involved in ICI, rigorous modelling of the accretion process remains limited. The present study proposes a novel modelling approach based on the physically observed mixed-phase nature of the accretion layers. The mathematical model, which is derived from the enthalpy change after accretion (the enthalpy model), is compared with an existing pure-phase layer model (the three-layer model). Scaling laws and asymptotic solutions are developed for both models. The onset of ice accretion, the icing layer thickness and solid ice fraction within the layer are determined by a set of non-dimensional parameters including the Péclet number, the Stefan number, the Biot number, the melt ratio and the evaporative rate. Thresholds for freezing and non-freezing conditions are developed. The asymptotic solutions present good agreement with numerical solutions at low Péclet numbers. Both the asymptotic and numerical solutions show that, when compared with the three-layer model, the enthalpy model presents a thicker icing layer and a thicker water layer above the substrate due to mixed-phased features and modified Stefan conditions. Modelling in terms of the enthalpy poses significant advantages in the development of numerical methods to complex three-dimensional geometrical and flow configurations. These results improve understanding of the accretion process and provide a novel, rigorous mathematical framework for accurate modelling of ICI.

A long-wave asymptotic model is developed for a viscoelastic falling film along the inside of a tube; viscoelasticity is incorporated using an upper convected Maxwell model. The dynamics of the resulting model in the inertialess limit is determined by three parameters: Bond number Bo, Weissenberg number We and a film thickness parameter $a$. The free surface is unstable to long waves due to the Plateau–Rayleigh instability; linear stability analysis of the model equation quantifies the degree to which viscoelasticity increases both the rate and wavenumber of maximum growth of instability. Elasticity also affects the classification of instabilities as absolute or convective, with elasticity promoting absolute instability. Numerical solutions of the nonlinear evolution equation demonstrate that elasticity promotes plug formation by reducing the critical film thickness required for plugs to form. Turning points in travelling wave solution families may be used as a proxy for this critical thickness in the model. By continuation of these turning points, it is demonstrated that in contrast to Newtonian films in the inertialess limit, in which plug formation may be suppressed for a film of any thickness so long as the base flow is strong enough relative to surface tension, elasticity introduces a maximum critical thickness past which plug formation occurs regardless of the base flow strength. Attention is also paid to the trade-off of the competing effects introduced by increasing We (which increases growth rate and promotes plug formation) and increasing Bo (which decreases growth rate and inhibits plug formation) simultaneously.

Describing the evolution of a wind turbine's wake from a top-hat profile near the turbine to a Gaussian profile in the far wake is a central feature of many engineering wake models. Existing approaches, such as super-Gaussian wake models, rely on a set of tuning parameters that are typically obtained from fitting high-fidelity data. In the current study, we present a new engineering wake model that leverages the similarity between the shape of a turbine's wake normal to the streamwise direction and the diffusion of a passive scalar from a disk source. This new wake model provides an analytical expression for a streamwise scaling function that ensures the conservation of linear momentum in the wake region downstream of a turbine. The model also considers the different rates of wake expansion that are known to occur in the near- and far-wake regions. Validation is presented against high-fidelity numerical data and experimental measurements from the literature, confirming a consistent good agreement across a wide range of turbine operating conditions. A comparison is also drawn with several existing engineering wake models, indicating that the diffusion-based model consistently provides more accurate wake predictions. This new unified framework allows for extensions to more complex wake profiles by making adjustments to the diffusion equation. The derivation of the proposed model included the evaluation of analytical solutions to several mathematical integrals that can be useful for other physical applications.

Not all the information in a turbulent field is relevant for understanding particular regions or variables in the flow. Here, we present a method for decomposing a source field into its informative $\boldsymbol {\varPhi }_{I}(\boldsymbol {x},t)$ and residual $\boldsymbol {\varPhi }_{R}(\boldsymbol {x},t)$ components relative to another target field. The method is referred to as informative and non-informative decomposition (IND). All the necessary information for physical understanding, reduced-order modelling and control of the target variable is contained in $\boldsymbol {\varPhi }_{I}(\boldsymbol {x},t)$, whereas $\boldsymbol {\varPhi }_{R}(\boldsymbol {x},t)$ offers no substantial utility in these contexts. The decomposition is formulated as an optimisation problem that seeks to maximise the time-lagged mutual information of the informative component with the target variable while minimising the mutual information with the residual component. The method is applied to extract the informative and residual components of the velocity field in a turbulent channel flow, using the wall shear stress as the target variable. We demonstrate the utility of IND in three scenarios: (i) physical insight into the effect of the velocity fluctuations on the wall shear stress; (ii) prediction of the wall shear stress using velocities far from the wall; and (iii) development of control strategies for drag reduction in a turbulent channel flow using opposition control. In case (i), IND reveals that the informative velocity related to wall shear stress consists of wall-attached high- and low-velocity streaks, collocated with regions of vertical motions and weak spanwise velocity. This informative structure is embedded within a larger-scale streak–roll structure of residual velocity, which bears no information about the wall shear stress. In case (ii), the best-performing model for predicting wall shear stress is a convolutional neural network that uses the informative component of the velocity as input, while the residual velocity component provides no predictive capabilities. Finally, in case (iii), we demonstrate that the informative component of the wall-normal velocity is closely linked to the observability of the target variable and holds the essential information needed to develop successful control strategies.

We analyse the motion of a flagellated bacterium in a two-fluid medium using slender body theory. The two-fluid model is useful for describing a body moving through a complex fluid with a microstructure whose length scale is comparable to the characteristic scale of the body. This is true for bacterial motion in biological fluids (entangled polymer solutions), where the entanglement results in a porous microstructure with typical pore diameters comparable to or larger than the flagellar bundle diameter, but smaller than the diameter of the bacterial head. Thus, the polymer and solvent satisfy different boundary conditions on the flagellar bundle and move with different velocities close to it. This gives rise to a screening length $L_B$ within which the fluids exchange momentum and the relative velocity between the two fluids decays. In this work, both the solvent and polymer of the two-fluid medium are modelled as Newtonian fluids with different viscosities $\mu _s$ and $\mu _p$ (viscosity ratio $\lambda = \mu _p/\mu _s$), thereby capturing the effects solely introduced by the microstructure of the complex fluid. From our calculations, we observe an increased drag anisotropy for a rigid, slender flagellar bundle moving through this two-fluid medium, resulting in an enhanced swimming velocity of the organism. The results are sensitive to the interaction between the bundle and the polymer, and we discuss two physical scenarios corresponding to two types of interaction. Our model provides an explanation for the experimentally observed enhancement of swimming velocity of bacteria in entangled polymer solutions and motivates further experimental investigations.

We consider the initial ‘slumping phase’ of a lock-release gravity current (GC) on a down slope with focus on particle-driven (turbidity) flows, in the inertia–buoyancy (large Reynolds number) and Boussinesq regime. We use a two-layer shallow-water (SW) model for the depth-averaged variables, and compare the predictions with previously published experimental data. In particular, we analyse the empirical conclusion of Gadal et al. (J. Fluid Mech., vol. 974, 2023, A4) that the slumping displays a constant speed for a significant range of slopes and particle-sedimentation speeds. We emphasize the physical definition of the slumping phase (stage): the adjustment process during which (a) the fluid in the lock is set into motion by the dam break, then (b) forms a tail from the backwall to the nose. We focus on the question of if and when the propagation speed $u_N$ of the nose (front) of the GC is constant during this process (there is consensus that a significant deceleration of $u_N$ appears in the post-slumping stage.) The SW theory predicts correctly the adjustment of the flow field during the slumping stage, but indicates that a constant $u_N$ appears only for the classical case ($\gamma =E=c_D=\beta =0$) where $\gamma, E, c_D, \beta$ are the slope, entrainment and drag coefficients, and the scaled particle settling speed for a particle-driven GC. However, since $\gamma, E, c_D, \beta$ are typically small, the change of $u_N$ during the slumping phase is also small in many cases of interest. The interaction between the various driving and hindering mechanisms is elucidated. We show that, in a system with a horizontal (open) top (typical laboratory experiments), the height of the ambient increases along the slope, and this compensates for buoyancy loss due to particle sedimentation. We point out the need for further experimental and simulation studies for a better understanding of the slumping phase and transition to the next phases, and further assessment/improvement of the SW predictions.

Metastructures composed of a closely spaced plate array have been widely used in bespoke manipulation of waves in contexts of acoustics, electromagnetics, elasticity and water waves. This paper focuses on wave scattering by discrete plate array metastructures of arbitrary cross-sections, including isolated vertical metacylinders, periodic arrays and horizontal surface-piercing metacylinders. A suitable transform-based method has been applied to each problem to reduce the influence of barriers in a two-dimensional problem to a set of points in a one-dimensional wave equation wherein the solution is constructed using a corresponding Green's function. A key difference from the existing work is the use of an exact description of the plate array rather than an effective medium approximation, enabling the exploration of wave frequencies above resonance where homogenisation models fail but where the most intriguing physical findings are unravelled. The new findings are particularly notable for graded plate array metastructures that produce a dense spectrum of resonant frequencies, leading to broadband ‘rainbow reflection’ effects. This study provides new ideas for the design of structures for the bespoke control of waves with the potential for innovative solutions to coastal protection schemes or wave energy converters.

A gas bubble sitting at a liquid–gas interface can burst following the rupture of the thin liquid film separating it from the ambient, owing to the large surface energy of the resultant cavity. This bursting bubble forms capillary waves, a Worthington jet and subsequent droplets for a Newtonian liquid medium. However, rheological properties of the liquid medium like elastoviscoplasticity can greatly affect these dynamics. Using direct numerical simulations, this study exemplifies how the complex interplay between elasticity (in terms of elastic stress relaxation) and yield stress influences the transient interfacial phenomenon of bursting bubbles. We investigate how bursting dynamics depends on capillary, elastic and yield stresses by exploring the parameter space of the Deborah number ${{\textit {De}}}$ (dimensionless relaxation time of elastic stresses) and the plastocapillary number $\mathcal {J}$ (dimensionless yield-stress of the medium), delineating four distinct characteristic behaviours. Overall, we observe a non-monotonic effect of elastic stress relaxation on the jet development while plasticity of the elastoviscoplastic (EVP) medium is shown to affect primarily the jet evolution only at faster relaxation times (low ${{\textit {De}}}$). The role of elastic stresses on jet development is elucidated with the support of energy budgets identifying different modes of energy transfer within the EVP medium. The effects of elasticity on the initial progression of capillary waves and droplet formation are also studied. In passing, we study the effects of solvent–polymer viscosity ratio on bursting dynamics and show that polymer viscosity can increase the jet thickness apart from reducing the maximum height of the jet.

The spherical Couette system consists of two differentially rotating concentric spheres with the space in between filled with fluid. We study a regime where the outer sphere is rotating rapidly enough so that the Coriolis force is important and the inner sphere is rotating either slower or in the opposite direction with respect to the outer sphere. We numerically study the sudden transition to turbulence at a critical differential rotation seen in experiments at BTU Cottbus-Senftenberg, Germany, and investigate its cause. We find that the source of turbulence is the boundary layer on the inner sphere, which becomes centrifugally unstable. We show that this instability leads to generation of small-scale structures which lead to turbulence in the bulk, dominated by inertial waves, a change in the force balance near the inner boundary, the formation of a mean flow through Reynolds stresses and, consequently, to an efficient angular momentum transport. We compare our findings with axisymmetric simulations and show that there are significant similarities in the nature of the flow in the turbulent regimes of full three-dimensional and axisymmetric simulations but differences in the evolution of the instability that leads to this transition. We find that a heuristic argument based on a Reynolds number defined using the thickness of the boundary layer as a length scale helps explain the scaling law of the variation of critical differential rotation for transition to turbulence with rotation rate observed in the experiments.

A model of imbibition dynamics in a channel of flattened triangular cross-section is presented, taking into account the liquid film flow in the corners of the channel. The quasi-analytical solutions are derived on the basis of a lubrication approximation. The analysis encompasses two imbibition scenarios corresponding to a constant flow rate or constant pressure imposed in the wetting fluid at the inlet of the channel. In the former case, the process starts with a liquid film flow regime in the corners that is followed by a bulk and corner film flow regime characterised by a triple point advancing (far) ahead of the bulk meniscus after its entrance in the channel. In the latter case, the occurrence of the bulk and corner film flow regime is conditioned by an imposed pressure yielding a capillary pressure at the inlet smaller than a threshold capillary pressure. Above this threshold, the liquid film regime remains. For both imbibition scenarios under concern, important features are highlighted, including (i) the time scalings of the dynamics of both the triple point and apex of the bulk meniscus (when it exists), (ii) the contrast in the positions of these two points showing that the classical Washburn approach, which neglects the effect of the corner films, overpredicts the dynamics of the bulk meniscus. The important consequence is an early wetting fluid breakthrough at the channel outlet much before the bulk meniscus arrival. Comparisons with experimental data available in the literature are provided, validating the approach proposed in this work.

Capsules, which are potentially active fluid droplets enclosed in a thin elastic membrane, experience large deformations when placed in suspension. The induced fluid–structure interaction stresses can potentially lead to rupture of the capsule membrane. While numerous experimental studies have focused on the rheological behaviour of capsules until rupture, there remains a gap in understanding the evolution of their mechanical properties and the underlying mechanisms of damage and breakup under flow. We here investigate the damage and rupture of bioartificial microcapsules made of ovalbumin reticulated with terephthaloyl chloride and placed in simple shear flow. We characterize damage by identifying how the surface shear modulus of the capsule membrane changes over time. Rupture is then characterized by comparing the number and size distribution of capsules before and after exposure to shear, while varying the shear rates and time during which capsules are sheared. Our findings reveal how the percentage of ruptured capsules increases with their size, exposure time to shear, and the ratio of viscous to elastic forces at rupture.

The coherent vortical structures in turbulent flow through a strong 16 : 1 three-dimensional contraction are studied using time-resolved volumetric measurements. Visualization using the vorticity magnitude criterion shows the emergence of long, stretched cylindrical vortices aligned with the mean flow. This alignment is quantified by probability density functions (p.d.f.s) of the direction cosines. We propose two measures to quantify the alignment, the peak height in the probability and a coefficient from the moment of the p.d.f., both of which reaffirm the strong streamwise alignment. The root mean square streamwise vorticity grows within the contraction to become 4.5 times larger than the transverse component, at the downstream location where the contraction ratio is $C=11$. The characteristic vortices become as long as the measurement volume, or more than 4 times the integral scale at the entrance to the contraction. We also characterize the vorticity enhancement along individual vortices, measuring 65 % strengthening over the distance where $C$ goes from 4 to 11. The prevalence of these coherent structures is estimated from 700 000 measured volumes, showing that near the outlet, it is more likely to have one or two of these structures present than none.

We explore the application of the reference map technique, originally developed for Eulerian simulation of solid mechanics, in Lagrangian kinematics of turbulent flows. Unlike traditional methods based on explicit particle tracking, the reference map facilitates the calculation of flow maps and gradients without the need for particles. This is achieved through an Eulerian update of the reference map, which records the take-off positions of fluid particles. This approach is found to be mathematically equivalent to the work of Leung (J. Comput. Phys., vol. 230, issue 9, 2011, pp. 3500–3524), who computed the flow map of simple two-dimensional flows using an Eulerian approach. We discuss important modifications necessary for its first application to complex three-dimensional turbulent flows, including the conservative, low-dissipation update of the flow map and the treatment of periodic boundary conditions. We first demonstrate the accuracy of finite-time Lyapunov exponent (FTLE) calculations based on the reference map against the standard particle-based approach in a two-dimensional Taylor–Green vortex. Then we apply it to turbulent channel flow at $Re_\tau =180$, where Lagrangian coherent structures identified as ridges of the backward-time FTLE are found to bound vortical regions of flow, consistent with Eulerian coherent structures from the $Q$-criterion. The reference map also proves suitable for material surface tracking despite not explicitly tracking particles. This capability can provide valuable insights into the Lagrangian landscape of turbulent momentum transport, complementing Eulerian velocity field analysis. The evolution of initially wall-normal material surfaces in the viscous sublayer, buffer layer and log layer sheds light on the Reynolds stress-generating events from a Lagrangian perspective.

We theoretically and experimentally study gravity currents of a Newtonian fluid advancing in a two-dimensional, infinite and saturated porous domain over a horizontal impermeable bed. The driving force is due to the density difference between the denser flowing fluid and the lighter, immobile ambient fluid. The current is taken to be in the Darcy–Forchheimer regime, where a term quadratic in the seepage velocity accounts for inertial contributions to the resistance. The volume of fluid of the current varies as a function of time as $\sim T^{\gamma }$, where the exponent parameterizes the case of constant volume subject to dam break ($\gamma =0$), of constant ($\gamma =1$), waning ($\gamma <1$) and waxing inflow rate ($\gamma >1$). The nonlinear governing equations, developed within the lubrication theory, admit self-similar solutions for some combinations of the parameters involved and for two limiting conditions of low and high local Forchheimer number, a dimensionless quantity involving the local slope of the current profile. Another parameter $N$ expresses the relative importance of the nonlinear term in Darcy–Forchheimer's law; values of $N$ in practical applications may vary in a large interval around unity, e.g. $N\in [10^{-5},10^{2}]$; in our experiments, $N\in [2.8,64]$. Sixteen experiments with three different grain sizes of the porous medium and different inflow rates corroborate the theory: the experimental nose speed and current profiles are in good agreement with the theory. Moreover, the asymptotic behaviour of the self-similar solutions is in excellent agreement with the numerical results of the direct integration of the full problem, confirming the validity of a relatively simple one-dimensional model.

The impact of a chemical reaction, $A+B \rightarrow C$, on the stability of a miscible radial displacement in a porous medium is established. Our study involves a comprehensive analysis employing both linear stability analysis and nonlinear simulations. Through linear stability analysis, the onset of instability for monotonic as well as non-monotonic viscosity profiles corresponding to the same end-point viscosity are discussed and compared. We establish a $(R_b,R_c)$ phase plane for a wide range of Damköhler number ($Da$) and Péclet number ($Pe$) into stable and unstable regions. Here, $R_b=\ln (\mu _B/ \mu _A)$ and $R_c=\ln (\mu _C/ \mu _A)$ and $\mu _{i}$ is the viscosity of fluid $i$ $\in \{A$, $B$, $C$}. The stable zone in the $(R_b, R_c)$ phase plane contracts with increased $Da$ and $Pe$ but never vanishes. It exists even for $Da \rightarrow \infty$. Interestingly, we obtain a $Da$ independent stable region in the neighbourhood of $R_c=R_b$ where no transition occurs in stability despite changes in reaction rate. The study allows us to acquire knowledge about the transition of the stability for varying $Da, Pe$ and different reactions classified using $R_b, R_c$.

We study numerically the flow around a spherical droplet set fixed in a linear shear flow with moderate shear rates ($Sr\leq 0.5$, $Sr$ being the ratio between the velocity difference across the drop and the relative velocity) over a wide range of external Reynolds numbers ($0.1<{{Re}}\leq 250$, ${{Re}}$ based on the slip velocity and the viscosity of the external fluid) and drop-to-fluid viscosity ratios ($0.01\leq \mu ^\ast \leq 100$). The flow structure, the vorticity field and their intrinsic connection with the lift force are analysed. Specifically, the results on lift force are compared with the low-${{Re}}$ solution derived for droplets of arbitrary $\mu ^\ast$, as well as prior data at finite ${{Re}}$ available in both the clean-bubble limit ($\mu ^\ast \to 0$) and the solid-sphere limit ($\mu ^\ast \to \infty$). Notably, at ${{Re}}=O(100)$, the lift force exhibits a non-monotonic transition from $\mu ^\ast \to 0$ to $\mu ^\ast \to \infty$, peaking at $\mu ^\ast \approx 1$. This behaviour is related to an internal three-dimensional flow bifurcation also occurring under uniform-flow conditions, which makes the flow to evolve from axisymmetric to biplanar symmetric. This flow bifurcation occurs at low-but-finite $\mu ^\ast$ when the internal Reynolds number (${{Re}}^i$, based on the viscosity of the internal fluid) exceeds approximately 300. In the presence of shear, the corresponding imperfect bifurcation enhances the extensional rate of the flow in the wake. Consequently, the streamwise vortices generated behind the droplet can be more intense compared with those behind a clean bubble. Given the close relation between the lift and these vortices, a droplet with ${{Re}}=O(100)$ and $\mu ^\ast \approx 1$ typically experiences a greater lift force than that in the inviscid limit.

We investigate the effect of high wind speeds on the breakup mechanisms that govern the formation of a spray from nozzles that form liquid sheets, which subsequently break up. The fragmentation mechanism of liquid sheets from spray nozzles has recently been described in detail under quiescent conditions. With high wind speeds, measurements of the droplet size distribution reveal two rather than one characteristic drop sizes, suggesting the existence of two distinct breakup mechanisms. High-speed images of the spray are used to identify these two mechanisms. We show that the smaller droplets result from the breakup of ‘bags’ formed in the spray sheet by the wind, while the larger droplets result from the breakup of the remaining perforated sheet. Based on the two mechanisms, a probability density function is constructed and fitted to the measured droplet size distributions. We show that the spray sheet destabilises due to the Rayleigh–Taylor instability induced by the airflow, and that the experimentally observable breakup length and size of the holes blown in the sheet are predicted by the fastest growing wavenumber. From this, a theoretical prediction for the droplet size from bag breakup and remaining sheet breakup is derived.

The study presents observations on the interaction of double-blade propeller tip vortices with a smooth-wall turbulent boundary layer (TBL). The wall-bounded helicoidal vortices from the propeller modify the velocity profiles and turbulence statistics. The effects of two different tip clearances, $\epsilon = 0.1\delta _0$ and $0.5\delta _0$, at a matched thrust, are explored with particle image velocimetry to understand the dynamics of tip-vortex formation within the logarithmic and wake regions of the boundary layer. The measurements are performed with $\lambda =U_{tip}/U_{\infty }$ in the range 5.3–5.9, and a blade passing frequency ($\,f_{prop}$) of the same order of the boundary-layer time scale ($\,f_{TBL}$). Observations indicate a reduction in the extent of the log region and an enhancement of the wake parameter $\varPi$, mirroring the behaviour seen in TBLs under adverse pressure gradient conditions. Notably, the slipstream most contracted region exhibits a significant reduction in the skin friction coefficient $C_f$ and an amplification of the velocity fluctuation statistics across the entire boundary layer. At a clearance of $\epsilon = 0.1\delta _0$, there is evidence of the formation of paired coherent wall-bounded structures. The presence of the wall decreases the amplitude of both periodic and stochastic fluctuations obtained with a phase-locked triple decomposition. An exception is observed behind the propeller for the stochastic fluctuations of the wall-normal component of the flow, which become amplified as the blades move away from the wall. This leads to the creation of a more intense phase-locked two-point spatial coherence than that observed in fluctuations aligned with the streamwise direction. Furthermore, results reveal that reduced tip clearances lead to higher viscous dissipation and more active energy exchange between the mean flow and organized motions.

We present new constrained and free-swimming experiments and simulations in the inertial regime, with Reynolds number $\mbox{Re} = O(10^4)$, of a pair of two-dimensional and three-dimensional pitching hydrofoils interacting in a minimal school. The hydrofoils have an out-of-phase synchronisation, and they are varied through in-line, staggered and side-by-side formations within the two-dimensional interaction plane. It is discovered that there is a two-dimensionally stable equilibrium point for a side-by-side formation. This formation is super-stable, meaning that hydrodynamic forces will passively maintain this formation even under external perturbations, and the school as a whole has no net forces acting on it that cause it to drift to one side or the other. Previously discovered one-dimensionally stable equilibria driven by wake vortex interactions are shown to be, in fact, two-dimensionally unstable, at least for an out-of-phase synchronisation. Additionally, it is discovered that a trailing-edge vortex mechanism provides the restorative force to stabilise a side-by-side formation. The stable equilibrium is further verified by experiments and simulations for freely swimming foils where dynamic recoil motions are present. When constrained, swimmers in compact side-by-side formations experience collective efficiency and thrust increases up to 40 % and 100 %, respectively, whereas slightly staggered formations output an even higher efficiency improvement of 84 %, with an 87 % increase in thrust. Freely swimming foils in a stable side-by-side formation show efficiency and speed enhancements of up to 9 % and 15 %, respectively. These newfound schooling performance and stability characteristics suggest that fluid-mediated equilibria may play a role in the control strategies of schooling fish and fish-inspired robots.