Roll patterns on floating ice shelves have been suggested to arise from viscous buckling under compressive stresses. A model of this process is explored, allowing for a power-law fluid rheology for ice. Linear stability theory of uniformly compressing base flows confirms that buckling modes can be unstable over a range of intermediate wavelengths when gravity does not play a dominant role. The rate of compression of the base flow, however, ensures that linear perturbations have wavelengths that continually shorten with time. As a consequence, linear instability only ever arises over a certain window of time $t$, and its strength can be characterised by finding the net amplification factor a buckling mode acquires for $t\to \infty$, beginning from a given initial wavenumber. Bi-axial compression, in which sideways straining flow is introduced to prevent the thickening of the base flow, is found to be more unstable than purely two-dimensional (or uni-axial) compression. Shear-thinning enhances the degree of instability in both uni-axial and bi-axial flow. The implications of the theoretical results for the glaciological problem are discussed.

This study experimentally investigates passive drag reduction on a sphere using azimuthally spaced surface protrusions under subcritical Reynolds numbers, focusing on the effects of the protrusion number at fixed surface coverage. The proposed surface modification strategy, termed partial protrusions, maintains a constant total protruded area while varying the number of protrusions $N$, thereby adjusting their azimuthal spacing. The objective is to determine whether such configurations can outperform the conventional full protrusion, in which protrusions continuously surround the azimuthal direction, and to elucidate the flow mechanisms behind any observed enhancement. Drag and flow field measurements reveal that increasing $N$ significantly improves aerodynamic performance. When $N$ exceeds a certain threshold, the partial-protrusion configuration achieves a greater drag reduction than the full-protrusion case, despite using only half the surface coverage. For low $N$, asymmetric pressure distributions across the protruded and smoothed sides induce unsteady separation delay, leading to shear-layer oscillations and elevated turbulent kinetic energy. As $N$ increases, the azimuthal spacing between protrusions decreases, promoting stable interaction between the two sides and leading to separation delay farther downstream than in the full-protrusion case, along with suppression of flow unsteadiness. These results demonstrate that a well-designed partial-protrusion configuration can outperform the full-protrusion configuration in drag reduction and unsteadiness control, offering new insights into effective passive flow control strategies for bluff body flows.

We study the combined effects of natural convection and rotation on the dissolution of a solute in a solvent-filled circular cylinder. The density of the fluid increases with increasing concentration of the dissolved solute, and we model this using the Oberbeck–Boussinesq approximation. The underlying moving-boundary problem has been modelled by combining the Navier–Stokes equations with the advection–diffusion equation and a Stefan condition for the evolving solute–fluid interface. We use highly resolved numerical simulations to investigate the flow regimes, dissolution rates and mixing of the dissolved solute for $Sc = 1$, $Ra \in [10^5, 10^8]$ and $\varOmega \in [0, 2.5]$. In the absence of rotation and buoyancy, the distance of the interface from its initial position follows a square root relationship with time ($r_d \propto \sqrt {t}$), which ceases to exist at a later time due to the finite-size effect of the liquid domain. We then explore the rotation parameter, considering a range of rotation frequency – from smaller to larger, relative to the inverse of the buoyancy-induced time scale – and Rayleigh number. We show that the area of the dissolved solute varies nonlinearly with time depending on $Ra$ and $\varOmega$. The symmetry breaking of the interface is best described in terms of $Ra/\varOmega ^2$.

We investigate the unsteady lift response of compliant membrane wings in hovering kinematics by combining analytical inviscid theory with experimental results. An unsteady aerodynamic model is derived for a compliant thin aerofoil immersed in incompressible inviscid flow of variable free-stream velocity at high angles of attack. The model, representing a spanwise section of a hovering membrane wing, assumes small membrane deformation and attached flow. These assumptions are supported by experiments showing that passive membrane deformation suppresses flow separation when hovering at angles of attack up to $55^\circ$. An analytically derived expression is obtained for the unsteady lift response, incorporating the classical Wagner and Theodorsen functions and the membrane dynamic response. This theoretical expression is validated against experimental water-tank measurements that are performed on hovering membrane wings at angles of attack of $35^\circ$ and $55^\circ$. Data from membrane deformation measurements is applied to the theoretical lift expression, providing the theoretical lift response prediction for each of the available experimental scenarios. Results of the comparison show that the proposed theory accurately predicts unsteady lift contributions from membrane deformation at high angles of attack, provided the deformation remains small and the flow is attached. This agreement between inviscid theory and experimental measurements suggests that when flow separation is suppressed, the unsteady aerodynamic theory is valid well beyond the typical low-angle-of-attack regime.

Smooth surface features were recently found to stabilise stationary cross-flow instability (CFI) of swept-wing boundary layers, thus holding potential for passive laminar flow control. Notably, the effect of surface features on the transition location exhibited a significant dependence on the CFI amplitude. In this work, numerical solutions of the harmonic Navier–Stokes (HNS) equations are used to explore the impact of a smooth surface hump on the linear and nonlinear development of stationary CFI under various perturbation amplitudes. Linear simulations identify regions of successive inhibited and enhanced perturbation growth. Despite the recovery of the base flow and perturbation kinetic energy to the reference (i.e. no-hump) state, significantly reduced perturbation growth is observed. The distorted perturbation profile due to the interaction with the hump is postulated to be responsible for this. Increasing the perturbation amplitude results in a response of the flow that is qualitatively similar to the linear case, albeit with increasing local destabilisation of new fundamental (i.e. primary wavelength) structures and higher-order harmonics near the wall. An energy budget analysis reveals that the growth of the fundamental incoming CFI is inhibited through the reduced effectiveness of the lift-up mechanism downstream of the hump. This is preceded by a spatial perturbation shape deformation, governed by (spanwise) transport terms. The results suggest that stabilisation of incoming stationary CFI via smooth surface humps is most effective at low incoming perturbation amplitudes. At higher perturbation amplitudes, newly formed near-wall structures, pre-conditioned by the incoming CFI, overtake the incoming CFI and could anticipate the transition process.

The modal and non-modal stability of laminar flow in a rectangular microchannel is investigated by incorporating the effects of Coriolis forces due to rotation, cross-sectional aspect ratio and superhydrophobic wall slip. The full Navier–Stokes equations are linearised into modified Orr–Sommerfeld and Squire equations, which are then formulated as an eigenvalue problem using small disturbances of the Tollmien–Schlichting type. These equations are subsequently solved by the spectral collocation method. The transition to instability in rotating microchannel flows, influenced by aspect ratio and slip conditions, is analysed through eigenvalue spectra and neutral stability curves. For non-modal analysis, we express the solution in matrix exponential form and then, using the singular value decomposition method, calculate the maximum energy growth. The study reveals that the flow becomes unstable in the presence of rotation at a critical Reynolds number of $ Re_c \approx 40$ for a low aspect ratio and $ Re_c \approx 50.4$ for a high aspect ratio. We find that instability is more pronounced in spanwise-rotating flows at higher aspect ratios compared with those at lower aspect ratios. Rotation induces disturbances from both walls along the spanwise direction, forming secondary flow structures near the centreline. Furthermore, we examine the influence of anisotropic slip by separately considering streamwise and spanwise slip as limiting cases. The numerical results demonstrate that while streamwise slip has a stabilising effect on rotating flows at small scales, a sufficiently large spanwise slip length can trigger instability at Reynolds numbers lower than those observed in the no-slip case. Rotation has the potential to enhance non-modal transient energy growth, while streamwise slip can effectively suppress this instability. These findings suggest that the onset of instability and transient energy growth can be effectively regulated by adjusting the aspect ratio and spanwise slip of the channel walls.

Previous studies show that at the small scales of stably stratified turbulence, the scale-dependent buoyancy flux reverses sign, corresponding to a conversion of turbulent potential energy (TPE) back into turbulent kinetic energy (TKE). Moreover, the magnitude of the reverse flux becomes stronger with increasing Prandtl number $\textit{Pr}$. Using a filtering analysis we demonstrate analytically how this flux reversal is connected to the mechanism identified in Bragg & de Bruyn Kops (2024 J. Fluid Mech. vol. 991, A10) that is responsible for the surprising observation that the TKE dissipation rate increases while the TPE dissipation rate decreases with increasing $\textit{Pr}$ in stratified turbulence. The mechanism identified by Bragg & de Bruyn Kops, which is connected to the formation of ramp–cliff structures in the density field, is shown to give the scale-local contribution to the buoyancy flux. At the smallest scales this local contribution dominates and explains the flux reversal, while at larger scales a non-local contribution is important. Direct numerical simulations of three-dimensional statistically stationary, strongly stably stratified turbulence confirm the theoretical analysis, and indicate that, while on average the local contribution only dominates the buoyancy flux at the smallest scales, it remains strongly correlated with the buoyancy flux at all scales. The results show that ramp cliffs are not only connected to the reversal of the local buoyancy flux but also the non-local part. At the small scales (approximately below the Ozmidov scale), ramp structures contribute exclusively to reverse buoyancy flux events, whereas cliff structures contribute to both forward and reverse buoyancy flux events.

Direct numerical simulation is performed to study the effects of spanwise curvature on transitioning and turbulent boundary layers. Turbulent transition is induced with an array of resolved cuboids. Spanwise curvature is prescribed using a novel approach with a body force that is applied orthogonally to the bulk flow to curve the mean free-stream streamlines at a set radius. The flows are analysed in a streamline-aligned coordinate system. Although the radius of curvature is large compared with the size of the boundary layer, its effects on the development of the boundary layer are appreciable. The results indicate that spanwise curvature induces a non-uniform mean secondary flow and alters the structure of turbulence within the boundary layer. Analytical expressions for the crossflow are derived in the viscous sublayer and log layer. These alterations are visible as changes in the distribution of the turbulent stresses and alignment of the vortical structures with the mean flow. These modifications are responsible for a misalignment between the Reynolds stress tensor and the velocity gradient tensor, which has important consequences for the validity of the widely used Boussinesq turbulent viscosity hypothesis in Reynolds-averaged Navier–Stokes models. Spanwise curvature was observed to decrease turbulent kinetic energy. These results have important implications on the development of turbulence in general applications, such as the flow over a prolate spheroid.

Aeroacoustic noise generated by wings under the wing-in-ground (WIG) effect is prevalent in various industrial applications, such as WIG vehicles, tower–blade interactions, and slat device noise issues. At chord-based Reynolds number 50 000 and freestream Mach number 0.3, the introduction of an engine jet transforms the separated stall noise of an NACA 4412 aerofoil under the influence of the WIG effect into laminar boundary layer vortex shedding (LBL-VS) noise. This study investigates the underlying mechanisms of this LBL-VS noise. Instead of relying on acoustic analogy, the unapproximated acoustic field is captured with high fidelity using direct numerical simulations. We identify the vorticity transfer process around the trailing edge as a key mechanism in the generation of LBL-VS noise. The results show that as the ground clearance decreases, the overall noise intensity is reduced. When the clearance becomes sufficiently small (10 % chord length), the well-organised vortex structures above the aerofoil break down under high adverse pressure, transitioning into a turbulent state that disrupts the vorticity transfer process. At this clearance, the dominant noise frequency drops from the vortex shedding frequency to an intermittent bursting frequency. This intermittent behaviour arises because only when certain vortices are amplified by acoustic feedback can they shed from the trailing edge, triggering the vorticity transfer process and generating pressure fluctuations. These findings provide new insights into the LBL-VS noise mechanisms under WIG conditions, and can inform strategies for noise reduction in relevant applications.

We integrate a discrete vortex method (DVM) with complex network analysis to strategise dynamic stall mitigation over aerofoils with active flow control. The objective is to inform the actuator placement and the timing to introduce control inputs during the highly transient process of dynamic stall. To this end, we treat a massively separated flow as a network of discrete vortical elements and quantify the interactions among the vortical nodes by tracking the spread of displacement perturbations between each pair of vortical elements using a DVM. This allows us to perform network broadcast mode analysis to identify an optimal set of discrete vortices, the critical timing and the direction to seed perturbations as control inputs. Motivated by the objective of dynamic stall mitigation, the optimality is defined as maximising the reduction of total circulation of the free vortices generated from the leading edge over a prescribed time horizon. We demonstrate the use of the analysis on a two-dimensional flow over a flat plate aerofoil and a three-dimensional turbulent flow over an SD$7003$ aerofoil. The results from the network analysis reveal that the optimal timing for introducing disturbances occurs slightly after the onset of flow separation, before the shear layer rolls up and forms the core of the dynamic stall vortex. The broadcast modes also show that the vortical nodes along the shear layer are optimal for introducing disturbances, hence providing guidance to actuator placement. Leveraging these insights, we perform nonlinear simulations of controlled flows by introducing flow actuation that targets the shear layer slightly after the separation onset. We observe that the network-guided control results in a $21 \,\%$ and $14\,\%$ reduction in peak lift for flows over the flat plate and SD$7003$ aerofoil, respectively. A corresponding decrease in vorticity injection from the aerofoil surface under the influence of control is observed from simulations, which aligns with the objective of the network broadcast analysis. The study highlights the potential of integrating the DVMs with the network analysis to design an effective active flow control strategy for unsteady aerodynamics.

In this study, we revisit the validity of eddy viscosity models for predicting wave-induced airflow disturbances over ocean surface waves. We first derive a turbulence curvilinear model for the phase-averaged Navier–Stokes equations, extending the work of Cao, Deng & Shen (2020 J. Fluid Mech. 901, A27), by incorporating turbulence stress terms previously neglected in the linearised viscous curvilinear model. To verify our formulation, we perform a priori tests by numerically solving the model using mean wind and turbulence stress profiles from large-eddy simulations (LES) of airflow over waves across various wave ages. Results show that including turbulence stress terms improves wave-induced airflow predictions compared with the previous viscous curvilinear model. We further show that using a standard mixing-length eddy viscosity yields inaccurate predictions at certain wave ages, as it fails to capture wave-induced turbulence, which fundamentally differs from mean shear-driven turbulence. The LES data show that accurate representations of wave-induced stresses require a complex-valued eddy viscosity. The maximum magnitude of this eddy viscosity scales as $\sim \!u_\tau \zeta _{\textit{inner}}$, where $u_\tau$ is the friction velocity and $\zeta _{\textit{inner}}$ is the inner-layer thickness, the height at which the eddy-turnover time matches the wave advection time scale. This scaling aligns with the prediction by Belcher & Hunt (1993 J. Fluid Mech. 251, 109–148). Overall, the findings demonstrate that traditional eddy viscosity models are inadequate for capturing wave-induced turbulence. More sophisticated turbulence models are essential for the accurate prediction of airflow disturbances and form drag in wind–wave interaction models.

This work presents wavepacket models for supersonic round twin jets operating at perfectly expanded conditions, computed via plane-marching parabolised stability equations based on mean flows obtained from the compressible Reynolds-averaged Navier–Stokes (RANS) equations. High-speed schlieren visualisations and non-time-resolved PIV measurements are performed to obtain experimental datasets for validating the modelling strategy. The RANS solutions are found to be in good quantitative agreement with the particle image velocimetry (PIV) mean-flow measurements, confirming the ability of the approach to capture the interaction between jets at the mean-flow level. The obtained wavepackets consist of toroidal and flapping fluctuations of the twin-jet system, and show similarities with those of single axisymmetric jets. However, for the case of closely spaced jets, they exhibit deviations in the phase speed of structures travelling in the outer mixing layer and those travelling in the inner one, leading to different non-axisymmetric behaviours. In particular, toroidal twin-jet wavepackets feature tilted ring-like structures with respect to the jet axis, while flapping twin-jet wavepackets are distorted and lose the clean chequerboard pattern typically observed in $m = 1$ modes in axisymmetric jets. A quantitative comparison of the modelled wavepackets with experimentally educed coherent structures is performed in terms of their structural agreement measured through an alignment coefficient, providing a first validation of the modelling strategy. Alignment coefficients are found to be particularly high in the intermediate range of studied frequencies.

An experimental and computational analysis of a wing tip at moderate angle of attack highlights the leading role of the wing-tip vortex wandering along the direction grazing the wing-tip corner in generating far-field noise. The cases of Reynolds numbers $ \textit{Re}_c=0.6\times 10^6$ and $1.0\times 10^6$ at angle of attack $\alpha =10^\circ$ are presented. The vorticity field shows the existence of a system of three wing-tip vortices that co-rotate to form a helical structure. The vortices have wandering motions that develop as they travel downstream. Surface pressure measurements indicate the unsteadiness in the primary vortex to be coherent at a chord-based Strouhal number $ \textit{St}_c\approx 9$. The coherence between the surface pressure fluctuations and the far-field noise is the highest at the primary vortex crossover from the tip surface to the suction surface, which also occurs at $ \textit{St}_c\approx 9$. This is supported by computational results, where the crossover position on the wing surface experiences local maxima of pressure fluctuations at $ \textit{St}_c=9$, and the dilatation shows a wavefront emanating from the vortex crossover location. Given the downstream convection of the unsteadiness along the primary vortex, the crossover is suggested to be converting the pressure fluctuations in the vortex to acoustic waves rather than being a source of a new spectral feature. The causality correlation calculated between the surface pressure and the proper orthogonal decomposition modes of the flow field identifies the vortex kinematic modes that contribute the most to the surface pressure fluctuations at the vortex crossover.

Interfacial interactions between gas bubbles and the free surface are a hallmark of flows involving aqueous foams. In practice, bubble foams commonly arise from processes such as breaking waves at the ocean–atmosphere interface, plunging liquid jets and the effervescence of carbonated liquids. Once generated, bubbles within foam layers remain afloat at the free surface for finite durations before finally bursting into a fine spray of droplets. While the birth and bursting of bubble foams have received considerable attention, the understanding of floating bubbles is limited mainly to a single bubble. To build on this, in this article, we undertake numerical simulations of two or more floating bubbles in various canonical settings to examine their geometry and self-organising nature, with implications for real-world phenomena such as ocean spray production. Under lateral confinement, floating bubbles are prone to form vertically stacked layers. To this end, we analyse the geometry of coaxial pairs of floating bubbles and link geometrical differences between single and coaxial bubbles to various aspects of the ensuing bursting stage. Furthermore, we extend the existing theory of isolated floating bubbles to obtain unified analytical expressions for the shape parameters of single and coaxial bubbles of small sizes. Next, we investigate a pair of side-by-side floating bubbles, which serves as a minimal configuration to understand the formation of bubble rafts through self-organisation. We discover that Bond numbers in the range $10\leqslant \textit{Bo}\leqslant 50$ are more favourable for raft formation due to pronounced capillary attraction. The time required for two floating bubbles to assemble through capillary attraction grows exponentially with their initial separation. We also develop a linear model to capture the evolution of bubble spacing during capillary migration at low Bond numbers. Lastly, we extend the two-bubble configuration and showcase the emergent dynamics of a swarm of floating bubbles in mono- and bilayer configurations.

Pulsed gravity currents are generated by the sequential release of dense material into a lighter ambient. We investigate the dynamics of pulsed gravity currents using physical scale experiments, two-dimensional depth-averaged shallow water equation (SWE) based models and three-dimensional lattice Boltzmann method (LBM) simulations. Integrating these results we show for the first time that short duration pulsed releases generate intrusive layers, which accelerate front propagation relative to an instantaneously released current of the same total volume. Conversely, a long delay time between pulses produces a current that propagates slower than an equivalent instantaneous release. This finding is supported by physical experiments and depth-resolving LBM simulations. The depth-resolving simulations show that intrusions in pulsed flows experience less drag resistance than those generated by instantaneous releases. The depth-averaged model considered in the present study does not accurately capture the intrusive flow dynamics of pulsed currents. However, the limitations of the finite-depth SWE model may be mitigated by extensions to incorporate entrainment and density stratification. The results also motivate further research into the impact of buoyancy Reynolds number and channel slope on the propagation of pulsed currents.

Boundary-layer instability and transition control have drawn extensive attention from the hypersonic community. The acoustic metasurface has become a promising passive control method, owing to its straightforward implementation and lack of requirement for external energy input. Currently, the effects of the acoustic metasurface on the early and late transitional stages remain evidently less understood than the linear instability stage. In this study, the transitional stage of a flat-plate boundary layer at Mach 6 is investigated, with a particular emphasis on the nonlinear mode–mode interaction. The acoustic metasurface is modelled by the well-validated time-domain impedance boundary condition. First, the resolvent analysis is performed to obtain the optimal disturbances, which reports two peaks corresponding to the oblique first mode and the planar Mack second mode. These two most amplified responses are regarded as the dominant primary instabilities that trigger the transition. Subsequently, both optimal forcings are introduced upstream in the direct numerical simulation, which leads to pronounced detuned modes before breakdown. The takeaway is that the location of the acoustic metasurface is significant in minimising skin friction and delaying transition onset simultaneously. The bispectral mode decomposition results reveal the dominant energy-transfer routine along the streamwise direction – from primary modes to low-frequency detuned modes. By employing the acoustic metasurface, the nonlinear triadic interaction between high- and low-frequency primary modes is effectively suppressed, ultimately delaying transition onset, whereas the late interaction related to lower-frequency detuned modes is reinforced, promoting the late skin friction. The placement of the metasurface in the linearly unstable region of the second mode delays the transition, which is due to the suppressed streak in the oblique breakdown scenario. However, in the late stage of the transition, the acoustic metasurface induces an undesirable increment of skin friction overshoot due to the augmented shear-induced dissipation work, which mainly arises from reinforced detuned modes related to the combination resonance. Meanwhile, by restricting the location of the metasurface upstream of the overshoot region, this undesirable augmentation of skin friction can be eliminated. As a result, the reasonable placement of the metasurface is crucial to damping the early instability while causing less negative impacts on the late transitional stage.

This direct numerical simulation study analyses the transition-to-turbulence in plane Couette flow (PCF) with three-dimensional (3-D) roughness. Square ribs of height $k=0.2h$ (where $h$ is the half-channel height) and a streamwise pitch separation of $\lambda =10k$, classified as k-type roughness, are mounted on the stationary wall. This configuration features alternating rough and smooth zones in the spanwise direction and is referred to as 3-D k-type roughness. We compare the behaviour of 3-D k-type roughness in the transitional regime with that of two-dimensional (2-D) k-type roughness (Gokul & Narasimhamurthy 2024 J. Fluid Mech. vol. 1000, p. A40). The route-to-transition in 3-D k-type roughness confirms the formation of laminar–turbulent patterns, which exist in the transitional Reynolds number range $\textit{Re} \in [325, 350]$. This range interestingly overlaps with those for smooth PCF ($\textit{Re} \in [325, 400]$) and 2-D k-type roughness ($\textit{Re} \in [300, 325]$), thereby indicating that 3-D k-type roughness amalgamates the characteristics of both the rough and smooth zones. In striking contrast to the 2-D k-type roughness, the laminar–turbulent bands in the 3-D k-type configuration are of non-uniform bandwidth. The 3-D k-type roughness surprisingly introduces continuous competition among bands of opposite orientations, a characteristic unique in this case. Due to this competition, the large-scale flow associated with the oblique bands is never fully aligned with the diagonal direction. The smooth zones in the 3-D k-type roughness exhibit complex interactions, evident from oscillatory velocity signals and multiple high-energy peaks in the frequency spectra, which likely contribute to the competing patterns with opposite orientations.

We study the surfing motion of an active particle along a planar interface, separating a semi-infinite layer of gas from a deep layer of liquid. The interface-trapped particle self-propels, thanks to an uneven distribution of surface tension in its immediate vicinity, which itself results from a non-uniform release of an active agent from the particle’s surface. We use the reciprocal theorem in conjunction with singular perturbation expansions to calculate the leading-order contributions to the propulsion speed of the surfer due to the advective transport of mass and momentum when the Péclet and Reynolds numbers (denoted by $\textit{Pe}$ and $\textit{Re}$, respectively) are small but finite. Assuming that the surface tension varies linearly with the concentration of the agent with a slope of negative $\alpha$, we show, perhaps unexpectedly, that the normalised speed for a purely translating (but otherwise arbitrarily shaped) particle, independent of the agent discharge mechanism, can be expressed as $\mathscr{U} = 1 + \mathscr{A} ( 2 \textit{Pe} \ln \textit{Pe} + \textit{Re} \ln \textit{Re} ) + \mathscr{O}(\textit{Pe}) + \mathscr{O}(\textit{Re})$, where the prefactor $\mathscr{A}$ is positive for negative $\alpha$ and vice versa. For reference, the self-propulsion speed of autophoretic Janus spheres varies with $\textit{Pe}$ as $\mathscr{U} = 1 + \mathscr{B} \, \textit{Pe} + {\cdots}$, where $\mathscr{B}$ is positive when the mobility coefficient of the particle is negative and vice versa. Also, the speed of spherical squirmers changes with $\textit{Re}$ as $\mathscr{U} = 1 + \mathscr{C} \, \textit{Re} + \mathscr{O}(\textit{Re})^2$, with $\mathscr{C}$ being positive for pushers and negative for pullers. Our asymptotic formula reveals that the speed of a Marangoni surfer is a non-monotonic function of the Péclet and Reynolds numbers, hinting at the existence of optimal values for both $\textit{Pe}$ and $\textit{Re}$. The information contained within the multiplier $\mathscr{A}$ also offers guidance for customising the shape of the surfer, as well as the release rate and configuration of the agent, to enhance the self-surfing performance. Our general theoretical analysis is complemented by detailed numerical simulations for a representative spherical surfer. These simulations confirm our theoretical predictions and shed light on the effects of intermediate and large values of $\textit{Pe}$ and $\textit{Re}$ on the performance of Marangoni surfers.

The coalescence and breakup of drops are classic examples of flows that feature singularities. The behaviour of viscoelastic fluids near these singularities is particularly intriguing – not only because of their added complexity, but also due to the unexpected responses they often exhibit. In particular, experiments have shown that the coalescence of viscoelastic sessile drops can differ significantly from that of their Newtonian counterparts, sometimes resulting in a sharply distorted interface. However, the mechanisms driving these differences in dynamics, as well as the potential influence of the contact angle are not fully known. Here, we study two different flow regimes effectively induced by varying the contact angle and demonstrate how that leads to markedly different coalescence behaviours. We show that the coalescence dynamics is effectively unaltered by viscoelasticity at small contact angles. The Deborah number, which is the ratio of the relaxation time of the polymer to the time scale of the background flow, scales as $\theta ^3$ for $\theta \ll 1$, thus rationalising the near-Newtonian response. On the other hand, it has been shown previously that viscoelasticity dramatically alters the shape of the interface during coalescence at large contact angles. We study this large contact angle limit using two-dimensional numerical simulations of the equation of motion. We show that the departure of the coalescence dynamics from the Newtonian case is a function of the Deborah number and the elastocapillary number, which is the ratio between the shear modulus of the polymer solution and the characteristic stress in the fluid.

We study natural convection in porous media using a lattice Boltzmann method that recovers the incompressible Navier–Stokes–Fourier dynamics. The porous structure consists of a staggered two-dimensional cylinder array with half-cylinders at the walls, forming a Darcy continuum at the domain scale. Hydrodynamic reference simulations reveal distinct flow regimes: laminar (Darcy), steady inertial (Forchheimer) and vortex shedding. We then analyse the effects of porosity and solid-to-fluid conductivity ratio ($k_s/k_{\!f}$) on natural convection. At low porosity ($\varphi = 33\,\%$), convection is highly sensitive to thermal coupling, particularly for insulating solids, whereas conductive matrices buffer this effect through lateral diffusion. Increasing porosity ($\varphi = 43\,\%$) smooths the transition as solid and fluid phases become more balanced. Across the explored range, two inertial regimes emerge governed by plume-scale confinement. The transition from Darcy to inertia-driven convection begins once the dynamics resembles the Forchheimer regime of the reference simulations. Based on our data, the system is governed by the confinement parameter $\varLambda$, which relates the plume-neck width, equivalent to the thermal boundary-layer thickness, to the pore scale: for $\varLambda \gtrsim 1$, the dynamics follows Forchheimer scaling, while for $\varLambda \lt 1/2$ it shifts toward Rayleigh–Bénard behaviour. Comparison with experimental data shows the same trend: the nominal Darcy–Rayleigh-to-porous-Prandtl ratio, $Ra^*/\textit{Pr}_{\!p} \approx 1$, holds for $\varLambda \gt 10$, but weaker confinement causes earlier departure. Finally, we revise benchmark Nusselt numbers for a cavity with square obstacles, showing that the reference by Merrikh & Lage (2005 Intl J. Heat Transfer 48(7), 1361–1372) misrepresents trends due to improper normalisation.