It is becoming increasingly clear that the strong spatial and temporal fluctuations observed in a narrow Reynolds number regime around the laminar–turbulent transition in shear flows can best be understood using the concepts and techniques from a seemingly unrelated discipline – statistical mechanics. During the last few years, a consensus has begun to emerge that these phenomena reflect an underlying non-equilibrium phase transition exhibited by a model of interacting particles on a crystalline lattice, directed percolation, that seems very far from fluid mechanics. Now, Chantry et al. (J. Fluid Mech., vol. 824, 2017, R1) have developed a truncated-mode computation of a model shear flow, capable of simulating systems far larger and longer than any previous study and have for the first time generated enough statistical data that a high-precision test of theory is feasible. The results broadly confirm the theory, extending the class of flows for which the directed percolation scenario holds and removing any remaining doubts that non-equilibrium statistical mechanical critical phenomena can be exhibited by the Navier–Stokes equations.

The influence of surface contamination upon the mass transfer rate of a low diffusivity gas across a flat surface is studied using direct numerical simulations. The interfacial mass transfer is driven by isotropic turbulence diffusing from below. Similar to Shen et al. (J. Fluid Mech., vol. 506, 2004, pp. 79–115) the surface contamination is modelled by relating the normal gradient of the horizontal velocities at the top to the horizontal gradients of the surfactant concentrations. A broad range of contamination levels is considered, including clean to severely contaminated conditions. The time-averaged results show a strong correlation between the gas transfer velocity and the clean surface fraction of the surface area. In the presence of surface contamination the mass transfer velocity $K_{L}$ is found to scale as a power of the Schmidt number, i.e. $Sc^{-q}$, where $q$ smoothly transitions from $q=1/2$ for clean surfaces to $q=2/3$ for very dirty interfaces. A power law $K_{L}\propto Sc^{-q}$ is proposed in which both the exponent $q$ and the constant of proportionality become functions of the clean surface fraction.

The use of static or dynamic roughness elements has been shown in the past to delay the separation of a laminar boundary layer from a solid surface. Here, we examine analytically the effect of such elements on the local and breakaway separation points, corresponding respectively to the position of zero skin friction and presence of a singularity in the roughness region, for flow over a hump embedded within the boundary layer. Two types of roughness elements are studied: the first is small and placed near the point of vanishing skin friction; the second is larger and extends downstream. The forced flow solution is found as a sum of Fourier modes, reflecting the fixed frequency forcing of the dynamic roughness. Solutions for both the static and dynamic roughness show that the presence of the roughness element is able to move the separation points downstream, given an appropriate choice of roughness frequency, height, position and width. This choice is found to be qualitatively similar to that observed for leading-edge separation. Furthermore, for a negative static roughness a small region of separated flow forms at high roughness depth, although there is a critical depth above which boundary-layer breakaway moves suddenly upstream.

In this paper we investigate, using theory and direct numerical simulations (DNS), the forward in time (FIT) and backward in time (BIT) probability density functions (PDFs) of the separation of inertial particle pairs in isotropic turbulence. In agreement with our earlier study (Bragg et al., Phys. Fluids, vol. 28, 2016, 013305), where we compared the FIT and BIT mean-square separations, we find that inertial particles separate much faster BIT than FIT, with the strength of the irreversibility depending upon the final/initial separation of the particle pair and their Stokes number $St$. However, we also find that the irreversibility shows up in subtle ways in the behaviour of the full PDF that it does not in the mean-square separation. In the theory, we derive new predictions, including a prediction for the BIT/FIT PDF for $St\geqslant O(1)$, and for final/initial separations in the dissipation regime. The prediction shows how caustics in the particle relative velocities in the dissipation range affect the scaling of the pair-separation PDF, leading to a PDF with an algebraically decaying tail. The predicted functional behaviour of the PDFs is universal, in that it does not depend upon the level of intermittency in the underlying turbulence. We also analyse the pair-separation PDFs for fluid particles at short times, and construct theoretical predictions using the multifractal formalism to describe the fluid relative velocity distributions. The theoretical and numerical results both suggest that the extreme events in the inertial particle-pair dispersion at the small scales are dominated by their non-local interaction with the turbulent velocity field, rather than due to the strong dissipation range intermittency of the turbulence itself. In fact, our theoretical results predict that for final/initial separations in the dissipation range, when $St\gtrsim 1$, the tails of the pair-separation PDFs decay faster as the Taylor Reynolds number $Re_{\unicode[STIX]{x1D706}}$ is increased, the opposite of what would be expected for fluid particles.

Solids dispersion is an important part of hydraulic fracturing, both in helping to understand phenomena such as tip screen-out and spreading of the pad, and in new process variations such as cyclic pumping of proppant. Whereas many frac fluids have low viscosity, e.g. slickwater, others transport proppant through increased viscosity. In this context, one method for influencing both dispersion and solids-carrying capacity is to use a yield stress fluid as the frac fluid. We propose a model framework for this scenario and analyse one of the simplifications. A key effect of including a yield stress is to focus high shear rates near the fracture walls. In typical fracturing flows this results in a large variation in shear rates across the fracture. In using shear-thinning viscous frac fluids, flows may vary significantly on the particle scale, from Stokesian behaviour to inertial behaviour across the width of the fracture. Equally, according to the flow rates, Hele-Shaw style models give way at higher Reynolds number to those in which inertia must be considered. We develop a model framework able to include this range of flows, while still representing a significant simplification over fully three-dimensional computations. In relatively straight fractures and for fluids of moderate rheology, this simplifies into a one-dimensional model that predicts the solids concentration along a streamline within the fracture. We use this model to make estimates of the streamwise dispersion in various relevant scenarios. This model framework also predicts the transverse distributions of the solid volume fraction and velocity profiles as well as their evolutions along the flow part.

Investigations about the role of nuclei and nucleation for the inception and formation of cavitation have been part of cavitation research since Harvey et al. (J. Cell. Physiol., vol. 24 (1), 1944, pp. 1–22) postulated the existence of gas filled crevices on surfaces and particles in liquids. In a supersaturated liquid, surface nuclei produce small gas bubbles due to mass transfer of gas or themselves work as weak spots in the liquid that are necessary for a phase change under technically relevant static pressures. Although various theories and models about nuclei and nucleation have found their way into standard literature, there is a lack of experimentally validated theories that describe the process of diffusion-driven nucleation in hydrodynamic cavitation. In order to close this gap we give new theoretical insights into the physics of this nucleation mechanism at technically relevant low supersaturations validated with extensive experimental results. The nucleation rate, the number of produced bubbles per second, is proportional to the supersaturation of the liquid and shows a nonlinear dependence on the shear rate at the surface nucleus. A model for the Strouhal number as dimensionless nucleation rate is derived allowing the estimation of nucleation rates from surface nuclei in hydrodynamic cavitation. The model provides three asymptotes, being a function of Péclet number, Weber number, the supersaturation of the liquid $\unicode[STIX]{x1D701}$ and gas solubility $\unicode[STIX]{x1D6EC}$ for three different detachment mechanisms, $Sr\propto \unicode[STIX]{x1D701}\unicode[STIX]{x1D6EC}We^{n}Pe^{1/3}$ with $n=1/3,3/4,1$. The theoretical findings are in good agreement with experimental results, leading to a new assessment of the role of diffusion in cavitating flows.

We calculate the leading-order correction to the time period of rotation of a neutrally buoyant spheroid of arbitrary aspect ratio, in a simple shear flow ($\boldsymbol{u}^{\infty }=\dot{\unicode[STIX]{x1D6FE}}y\mathbf{1}_{1}$; $\mathbf{1}_{1}$ is the unit vector in the flow direction, $y$ being the coordinate along the gradient direction), in its long-time orbit set up by the weak fluid inertial drift at $O(Re)$. Here, $Re$ is the microscale Reynolds number, a dimensionless measure of the fluid inertial effects on the length scale of the spheroid, and is defined as $Re=\dot{\unicode[STIX]{x1D6FE}}L^{2}\unicode[STIX]{x1D70C}/\unicode[STIX]{x1D707}$, where $L$ is the semimajor axis of the spheroid, $\unicode[STIX]{x1D707}$ and $\unicode[STIX]{x1D70C}$ are respectively the viscosity and density of the fluid, and $\dot{\unicode[STIX]{x1D6FE}}$ is the shear rate. This long-time orbit is the tumbling orbit for prolate spheroids; for oblate spheroids, it is the spinning orbit for aspect ratios greater than $0.137$, and can be either the tumbling or the spinning orbit for oblate spheroids of aspect ratios less than $0.137$. We also calculate the leading-order correction to the time period of rotation of a neutrally buoyant triaxial ellipsoid in a simple shear flow, rotating with its intermediate principal axis aligned along the vorticity of the flow; the latter calculation is in light of recent evidence, by way of numerical simulations (Rosen, PhD dissertation, 2016, Stockholm), of the aforementioned rotation being stabilized by weak inertia. The correction to the time period for arbitrary $Re$ is expressed as a volume integral using a generalized reciprocal theorem formulation. For $Re\ll 1$, it is shown that the correction at $O(Re)$ is zero for spheroids (with aspect ratios of order unity) as well as triaxial ellipsoids in their long-time orbits. The first correction to the time period therefore occurs at $O(Re^{3/2})$, and has a singular origin, arising from fluid inertial effects in the outer region (distances from the spheroid or triaxial ellipsoid of the order of the inertial screening length of $O(LRe^{-1/2})$), where the leading-order Stokes approximation ceases to be valid. Since the correction comes from the effects of inertia in the far field, the rotating spheroid (triaxial ellipsoid) is approximated as a time-dependent point-force-dipole singularity, allowing for the reciprocal theorem integral to be evaluated in Fourier space. It is shown for all relevant cases that fluid inertia at $O(Re^{3/2})$ leads to an increase in the time period of rotation compared with that in the Stokes limit, consistent with the results of recent numerical simulations at finite $Re$. Finally, combination of the $O(Re^{3/2})$ correction derived here with the $O(Re)$ correction derived earlier by Dabade et al. (J. Fluid Mech., vol. 791, 2016, 631703) yields a uniformly valid description of the first effects of inertia for spheroids of all aspect ratios, including prediction of the arrest of rotation for extreme-aspect-ratio spheroids.

Rayleigh–Taylor instability induced turbulence between two compressible miscible Newtonian fluids is studied in a strongly stratified configuration at a moderate Atwood number. A direct numerical simulation has been carried out with an auto-adaptive multidomain Chebyshev–Fourier–Fourier numerical method. The spatial resolution is increased up to $(9\times 100)\times 1000^{2}=900M$ collocation points. These numerical data are compared with those obtained from a simulation carried out at a lower Reynolds number and at the same Atwood number, and those obtained from a simulation carried out within the Boussinesq approximation at the same Reynolds number. A comprehensive data analysis is reported. Physical-variable mean profiles – density, concentration, temperature, entropy, velocity, vorticity, helicity and palinstrophy – are provided. Anisotropy is studied in the spectral space. The intermediate-scale isotropy and the small-scale anisotropy are exhibited for the scalars, i.e. concentration and temperature. Velocity is anisotropic at all scales but this anisotropy is more marked at small scales. The data are also analysed with the Favre-averaged equations. Sources of the turbulent kinetic energy, mass flux, root-mean-square density and energy equations are analysed. Compressibility effects are discussed in particular with the Kovàsznay-mode decomposition. A statistical study is reported where skewnesses, flatnesses and probability density functions (PDFs) are displayed and commented. A flow visualization is also given. Finally, the temperature field appears to be the slave of the mixing. This conclusion is drawn from the comparison of power spectra, anisotropy spectra, skewnesses, flatnesses, PDFs and correlation coefficients. There is however a significant time lag between the density and temperature evolution.

An emerging long obstacle placed in a boundary layer developing under a free surface generates a complex horseshoe vortex (HSV) system, which is composed of a set of vortices exhibiting a rich variety of dynamics. The present experimental study examines such flow structure and characterizes precisely, using particle image velocimetry (PIV) measurements, the evolution of the HSV geometrical and dynamical properties over a wide range of dimensionless parameters (Reynolds number $Re_{h}\in [750,8300]$ , boundary layer development ratio $h/\unicode[STIX]{x1D6FF}\in [1.25,4.25]$ and obstacle aspect ratio $W/h\in [0.67,2.33]$ ). The dynamical study of the HSV is based on the categorization of the motions of HSV vortices that result in an enhanced specific bi-dimensional typology, separating a coherent (due to vortex–vortex interactions) and an irregular evolution (due to the appearance of small-scale instabilities). This precise categorization is made possible thanks to the use of vortex tracking methods applied to PIV measurements; a semi-empirical model for the motion of the HSV vortices is then proposed to highlight some important mechanisms of the HSV dynamics, such as (i) the influence of the surrounding vortices on vortex motion and (ii) the presence of a phase shift between the motion of all vortices. Finally, the study of the HSV’s geometrical properties (vortex position and characteristic lengths and frequencies) evolution with the flow parameters shows that strong dependencies exist between the streamwise extension of the HSV and the obstacle width, and between the HSV vortex number and its elongation. Comparison of these data with prior studies for immersed obstacles reveals that emerging obstacles lead to greater adverse pressure gradients and down-flows in front of the obstacle. This implies a precocious separation of the boundary layer, leading to a larger HSV streamwise extension, and a lower vertical extension of the HSV, leading to smaller HSV vortices.

Integration across a fully developed turbulent channel flow of the transport equations for the mean and turbulent parts of the scalar dissipation rate yields relatively simple relations for the bulk mean scalar and wall heat transfer coefficient. These relations are tested using direct numerical simulation datasets obtained with two isothermal boundary conditions (constant heat flux and constant heating source) and a molecular Prandtl number Pr of 0.71. A logarithmic dependence on the Kármán number $h^{+}$ is established for the integrated mean scalar in the range $h^{+}\geqslant 400$ where the mean part of the total scalar dissipation exhibits near constancy, whilst the integral of the turbulent scalar dissipation rate $\overline{\unicode[STIX]{x1D700}_{\unicode[STIX]{x1D703}}}$ increases logarithmically with $h^{+}$. This logarithmic dependence is similar to that established in a previous paper (Abe & Antonia, J. Fluid Mech., vol. 798, 2016, pp. 140–164) for the bulk mean velocity. However, the slope (2.18) for the integrated mean scalar is smaller than that (2.54) for the bulk mean velocity. The ratio of these two slopes is 0.85, which can be identified with the value of the turbulent Prandtl number in the overlap region. It is shown that the logarithmic $h^{+}$ increase of the integrated mean scalar is intrinsically associated with the overlap region of $\overline{\unicode[STIX]{x1D700}_{\unicode[STIX]{x1D703}}}$, established for $h^{+}$ (${\geqslant}400$). The resulting heat transfer law also holds at a smaller $h^{+}$ (${\geqslant}200$) than that derived by assuming a log law for the mean temperature.

Static surface shapes of a magnetic fluid volume between two plates in a non-uniform magnetic field are investigated theoretically and experimentally. Abrupt changes and hysteresis of the magnetic fluid surface shape are observed in the experiments when the current in the coil increases and decreases quasi-statically. The necessary and sufficient conditions for a local minimum of the energy functional are derived theoretically. A method to find stable/unstable surface shapes is developed. The ambiguity in the determination of the magnetic fluid surface shape at the same value of the current is shown. It is found that the experimentally observed surface shapes of the given magnetic fluid volume coincide with the shapes obtained numerically, and practically all of them satisfy the derived necessary and sufficient conditions of the minimum energy. The stability curves of the magnetic fluid bridge between the plates are determined experimentally and theoretically.

We examine the problem of prescribing the macroscale boundary condition to the solute convective–diffusive mass transport equation at a heterogeneous surface consisting of reactive circular disks distributed uniformly on a non-reactive surface. The reaction rate at the disks is characterized by a first-order kinetics. This problem was examined by Shah & Shaqfeh (J. Fluid Mech., vol. 782, 2015, pp. 260–299) who obtained the boundary condition in terms of an effective first-order rate constant, which they determined as a function of the Péclet number $Pe=\dot{\unicode[STIX]{x1D6FE}}a^{2}/D$, the fraction $\unicode[STIX]{x1D719}$ of the surface area occupied by the reactive disks and the non-dimensional reaction rate constant $K=ka/D$. Here, $a$ is the radius of the disks, $D$ is the solute diffusivity, $\dot{\unicode[STIX]{x1D6FE}}$ is the wall shear rate and $k$ is the first-order surface-reaction rate constant. Their analysis assumed that $Pe$ and $K$ are $O(1)$ while the ratio of the microscale $a$ to the macroscale $H$ is small. The macroscale transport process is convection–diffusion dominated under these conditions. We examine here the case when the non-dimensional numbers based on the macroscale $H$ are $O(1)$. In this limit the microscale transport problem is reaction rate dominated. We find that the boundary condition can be expressed in terms of an effective rate constant only up to $O(\unicode[STIX]{x1D716})$, where $\unicode[STIX]{x1D716}=a/H$. Higher-order expressions for the mass flux involve both the macroscopic concentration and its surface gradient. The $O(\unicode[STIX]{x1D716})$ microscale problem is relatively easy to solve as the convective effects are unimportant and it is possible to obtain analytical expressions for the effective rate constant as a function of $\unicode[STIX]{x1D719}$ for both periodic and random arrangement of the disks without having to solve the boundary integral equation as was done by Shah and Shaqfeh. The results thus obtained are shown to be in good agreement with those obtained numerically by Shah and Shaqfeh for $Pe=0$. In a separate study, Shah et al. (J. Fluid Mech., vol. 811, 2017, pp. 372–399) examined the inverse-geometry problem in which the disks are inert and the rest of the surface surrounding them is reactive. We show that the two problems are related when $Pe=0$ and $kH/D=O(1)$. Finally, a related problem of determining the current density at a surface consisting of an array of microelectrodes is also examined and the analytical results obtained for the current density are found to agree well with the computed values obtained by solving the integral equation numerically by Lucas et al. (SIAM J. Appl. Maths, vol. 57(6), 1997, pp. 1615–1638) over a wide range of parameters characterizing this problem.

We perform a static analysis of a circular cylinder that forms a barrier between surfactant-laden and surfactant-free portions of a liquid–gas interface. In addition to determining the general implications of the balances for forces and torques, we quantify how the imbalance $\unicode[STIX]{x0394}\unicode[STIX]{x1D6FE}=\unicode[STIX]{x1D6FE}_{a}-\unicode[STIX]{x1D6FE}_{b}$ between the uniform surface tension $\unicode[STIX]{x1D6FE}_{a}$ of the surfactant-free portion of the interface and the uniform surface tension $\unicode[STIX]{x1D6FE}_{b}$ of the surfactant-laden portion of the interface influences the load-bearing capacity of a hydrophobic cylinder. Moreover, we demonstrate that the difference between surface tensions on either side of a cylinder with a cross-section of arbitrary shape induces a horizontal force component $f^{h}$ equal to $\unicode[STIX]{x0394}\unicode[STIX]{x1D6FE}$ in magnitude, when measured per unit length of the cylinder. With an energetic argument, we show that this relation also applies to a rod-like barrier with cross-sections of variable shape. In addition, we apply our analysis to amphiphilic Janus cylinders and we discuss practical implications of our findings for Marangoni propulsion and surface pressure measurements.

Compact solutions are presented for planetary, non-divergent, barotropic Rossby waves generated by (i) an impulsive point source and (ii) a sustained point source of curl of wind stress. Previously, only cumbersome integral expressions were known, rendering them practically useless. Our simple expressions allow for immediate numerical visualization/animation and further mathematical analysis.

The goal of the present work is to derive the closed-form expressions of coherence and admittances to describe the spatial distribution of lift on rectangular cylinders in turbulent flow, which can be used to investigate the three-dimensional effects of turbulence. The coherence of the three-dimensional aerodynamic admittance (3D AAF), which takes into full account the spanwise variations in the vertical velocity fluctuations, is introduced to assess the validity of the strip assumption. A theoretical coherence model expressed in a double-exponential form is derived starting from the two-wavenumber spectral tensor of the lift on a thin aerofoil in Fourier space, providing us with explicit insight into the coherence of the lift force. Notably, it is an intrinsic property that the lift force on the structure is more strongly correlated than the oncoming flow and 3D AAF. This coherence model is extended to rectangular cylinders by the introduction of three floating parameters into the decay parameters of the 3D AAF. Based on theoretical and experimental investigations, it is shown that the three-dimensional effects of turbulence grow more prominent as the difference between the decay parameters of the 3D AAF and vertical velocity fluctuations decreases. A generalized approach for rapidly deriving the closed-form expressions of the admittances is proposed to study the unsteady behaviour of the lift force and the distortion of the free stream passing through the rectangular cylinders.

Unsteady inviscid flow models of wings and airfoils have been developed to study the aerodynamics of natural and man-made flyers. Vortex methods have been extensively applied to reduce the dimensionality of these aerodynamic models, based on the proper estimation of the strength and distribution of the vortices in the wake. In such modelling approaches, one of the most fundamental questions is how the vortex sheets are generated and released from sharp edges. To determine the formation of the trailing-edge vortex sheet, the classical steady Kutta condition can be extended to unsteady situations by realizing that a flow cannot turn abruptly around a sharp edge. This condition can be readily applied to a flat plate or an airfoil with cusped trailing edge since the direction of the forming vortex sheet is known to be tangential to the trailing edge. However, for a finite-angle trailing edge, or in the case of flow separation away from a sharp corner, the direction of the forming vortex sheet is ambiguous. To remove any ad hoc implementation, the unsteady Kutta condition, the conservation of circulation as well as the conservation laws of mass and momentum are coupled to analytically solve for the angle, strength and relative velocity of the trailing-edge vortex sheet. The two-dimensional aerodynamic model together with the proposed vortex-sheet formation condition is verified by comparing flow structures and force calculations with experimental results for several airfoil motions in steady and unsteady background flows.

A theory for the low-Reynolds-number gravity-driven flow of two Newtonian fluids separated by a density interface in a two-dimensional annular geometry is developed. Solutions for the governing time-dependent equations of motion, in the limit that the radius of the inner and outer boundaries are similar, and in the case that the interface is initially inclined to the horizontal, are analysed numerically. We focus on the case in which the fluid is arranged symmetrically about a vertical line through the centre of the annulus. These solutions are successfully compared with asymptotic solutions in the limits that (i) a thin film of dense fluid drains down the outer boundary of the annulus, and (ii) a thin layer of less dense fluid is squeezed out of the narrow gap between the base of the inner annulus and dense fluid. Application of the results to the problem of mud displacement by cement in a horizontal well is briefly discussed.

The resonance modes in Mach 0.94 turbulent flow over a cavity having a length-to-depth ratio of five were explored using time-resolved particle image velocimetry (TR-PIV) and time-resolved pressure sensitive paint (TR-PSP). Mode switching was quantified in the velocity field simultaneous with the pressure field. As the mode number increased from one through three, the resonance activity moved from a region downstream within the recirculation region to areas further upstream in the shear layer, an observation consistent with linear stability analysis. The second and third modes contained organized structures associated with shear layer vortices. Coherent structures occurring in the velocity field during modes two and three exhibited a clear modulation in size with streamwise distance. The streamwise periodicity was attributable to the interference of downstream-propagating vortical disturbances with upstream-travelling acoustic waves. The coherent structure oscillations were approximately $180^{\circ }$ out of phase with the modal surface pressure fluctuations, analogous to a standing wave. Modal propagation (or phase) velocities, based on cross-correlations of bandpass-filtered velocity fields were found for each mode. The phase velocities also showed streamwise periodicity and were greatest at regions of maximum constructive interference where coherent structures were the largest. Overall, the phase velocities increased with modal frequency, which coincided with the modal activity residing at higher portions of the cavity where the local mean flow velocity was elevated. Together, the TR-PIV and TR-PSP provide unique details not only on the distribution of modal activity throughout the cavity, but also new understanding of the resonance mechanism as observed in the velocity field.

This article describes an unexplored transport phenomenon where a mildly viscoelastic medium encroaches a narrow capillary channel under the action of surface-tension force. The ultimate goal of the study is to provide the penetration length and the intrusion rate of the liquid as functions of time. The resulting analysis would be instrumental in building an inexpensive and convenient rheometric device which can measure the temporal scale for viscoelastic relaxation from the stored data of the aforementioned quantities. The key step in the formulation is a transient eigenfunction expansion of the instantaneous velocity profile. The time-dependent amplitude of the expansion as well as the intruded length are governed by a system of integro-differential relations which are derived by exploiting the mass and momentum conservation principles. The obtained integro-differential equations are simultaneously solved by using a fourth-order Runge–Kutta method assuming a start-up problem from rest. The resulting numerical solution properly represents the predominantly one-dimensional flow which gradually slows down after an initial acceleration and subsequent oscillation. The computational findings are independently verified by two separate perturbation theories. The first of these is based on a Weissenberg number expansion revealing the departure in the unsteady imbibition due to small but finite viscoelasticity. In contrast, the second one explains the long-time behaviour of the system by analytically predicting the decay features of the dynamics. These asymptotic results unequivocally corroborate the simulation inferring the accuracy of the numerics as well as the utility of the simplified mathematical models.

The $\unicode[STIX]{x1D707}(I)$-rheology was recently proposed as a potential candidate to model the incompressible flow of frictional grains in the dense inertial regime. However, this rheology was shown to be ill-posed in the mathematical sense for a large range of parameters, notably in the low and large inertial number limits (Barker et al., J. Fluid Mech., vol. 779, 2015, pp. 794–818). In this rapid communication, we extend the stability analysis of Barker et al. (J. Fluid Mech., vol. 779, 2015, pp. 794–818) to compressible flows. We show that compressibility regularizes the equations, making the problem well-posed for all parameters, with the condition that sufficient dissipation be associated with volume changes. In addition to the usual Coulomb shear friction coefficient $\unicode[STIX]{x1D707}$, we introduce a bulk friction coefficient $\unicode[STIX]{x1D707}_{b}$, associated with volume changes and show that the problem is well-posed if $\unicode[STIX]{x1D707}_{b}>1-7\unicode[STIX]{x1D707}/6$. Moreover, we show that the ill-posed domain defined by Barker et al. (J. Fluid Mech., vol. 779, 2015, pp. 794–818) transforms into a domain where the flow is unstable but remains well-posed when compressibility is taken into account. These results suggest the importance of taking into account dynamic compressibility for the modelling of dense granular flows and open new perspectives to investigate the emission and propagation of acoustic waves inside these flows.