Eigenmodes of averaged small-amplitude perturbations to a turbulent channel flow – which is one of the most fundamental canonical flows – are identified for the first time via an extensive set of high-fidelity graphics processing unit accelerated direct numerical simulations. While the system governing averaged small-amplitude perturbations to turbulent channel flow remains unknown, the fact such eigenmodes can be identified constitutes direct evidence that it is linear. Moreover, while the eigenvalue associated with the slowest-decaying anti-symmetric eigenmode mode is found to be real, the eigenvalue associated with the slowest-decaying symmetric eigenmode mode is found to be complex. This indicates that the unknown linear system governing the evolution of averaged small-amplitude perturbations cannot be self-adjoint, even for the case of a uni-directional flow. In addition to elucidating aspects of the flow physics, the findings provide guidance for development of new unsteady Reynolds-averaged Navier–Stokes turbulence models, and constitute a new and accessible benchmark problem for assessing the performance of existing models, which are used widely throughout industry.

We consider a three-dimensional acoustic field of an ideal gas in which all entropy production is confined to weak shocks and show that similar scaling relations hold for such a field as for forced Burgers turbulence, where the shock amplitude scales as $(\unicode[STIX]{x1D716}d)^{1/3}$ and the $p$th-order structure function scales as $(\unicode[STIX]{x1D716}d)^{p/3}r/d$, $\unicode[STIX]{x1D716}$ being the mean energy dissipation per unit mass, $d$ the mean distance between the shocks and $r$ the separation distance. However, for the acoustic field, $\unicode[STIX]{x1D716}$ should be replaced by $\unicode[STIX]{x1D716}+\unicode[STIX]{x1D712}$, where $\unicode[STIX]{x1D712}$ is associated with entropy production due to heat conduction. In particular, the third-order longitudinal structure function scales as $\langle \unicode[STIX]{x1D6FF}u_{r}^{3}\rangle =-C(\unicode[STIX]{x1D716}+\unicode[STIX]{x1D712})r$, where $C$ takes the value $12/5(\unicode[STIX]{x1D6FE}+1)$ in the weak shock limit, $\unicode[STIX]{x1D6FE}=c_{p}/c_{v}$ being the ratio between the specific heats at constant pressure and constant volume.

Granular flows occur in a wide range of situations of practical interest to industry, in our natural environment and in our everyday lives. This paper focuses on granular flow in the so-called inertial regime, when the rheology is independent of the very large particle stiffness. Such flows have been modelled with the $\unicode[STIX]{x1D707}(I),\unicode[STIX]{x1D6F7}(I)$-rheology, which postulates that the bulk friction coefficient $\unicode[STIX]{x1D707}$ (i.e. the ratio of the shear stress to the pressure) and the solids volume fraction $\unicode[STIX]{x1D719}$ are functions of the inertial number $I$ only. Although the $\unicode[STIX]{x1D707}(I),\unicode[STIX]{x1D6F7}(I)$-rheology has been validated in steady state against both experiments and discrete particle simulations in several different geometries, it has recently been shown that this theory is mathematically ill-posed in time-dependent problems. As a direct result, computations using this rheology may blow up exponentially, with a growth rate that tends to infinity as the discretization length tends to zero, as explicitly demonstrated in this paper for the first time. Such catastrophic instability due to ill-posedness is a common issue when developing new mathematical models and implies that either some important physics is missing or the model has not been properly formulated. In this paper an alternative to the $\unicode[STIX]{x1D707}(I),\unicode[STIX]{x1D6F7}(I)$-rheology that does not suffer from such defects is proposed. In the framework of compressible $I$-dependent rheology (CIDR), new constitutive laws for the inertial regime are introduced; these match the well-established $\unicode[STIX]{x1D707}(I)$ and $\unicode[STIX]{x1D6F7}(I)$ relations in the steady-state limit and at the same time are well-posed for all deformations and all packing densities. Time-dependent numerical solutions of the resultant equations are performed to demonstrate that the new inertial CIDR model leads to numerical convergence towards physically realistic solutions that are supported by discrete element method simulations.

The dynamics of irrotational shallow water wave turbulence forced at large scales and dissipated at small scales is investigated. First, we derive the shallow water analogue of the ‘four-fifths law’ of Kolmogorov turbulence for a third-order structure function involving velocity and displacement increments. Using this relation and assuming that the flow is dominated by shocks, we develop a simple model predicting that the shock amplitude scales as $(\unicode[STIX]{x1D716}d)^{1/3}$, where $\unicode[STIX]{x1D716}$ is the mean dissipation rate and $d$ the mean distance between the shocks, and that the $p$th-order displacement and velocity structure functions scale as $(\unicode[STIX]{x1D716}d)^{p/3}r/d$, where $r$ is the separation. Then we carry out a series of forced simulations with resolutions up to $7680^{2}$, varying the Froude number, $F_{f}=(\unicode[STIX]{x1D716}L_{f})^{1/3}/c$, where $L_{f}$ is the forcing length scale and $c$ is the wave speed. In all simulations a stationary state is reached in which there is a constant spectral energy flux and equipartition between kinetic and potential energy in the constant flux range. The third-order structure function relation is satisfied with a high degree of accuracy. Mean energy is found to scale approximately as $E\sim \sqrt{\unicode[STIX]{x1D716}L_{f}c}$, and is also dependent on resolution, indicating that shallow water wave turbulence does not fit into the paradigm of a Richardson–Kolmogorov cascade. In all simulations shocks develop, displayed as long thin bands of negative divergence in flow visualisations. The mean distance between the shocks is found to scale as $d\sim F_{f}^{1/2}L_{f}$. Structure functions of second and higher order are found to scale in good agreement with the model. We conclude that in the weak limit, $F_{f}\rightarrow 0$, shocks will become denser and weaker and finally disappear for a finite Reynolds number. On the other hand, for a given $F_{f}$, no matter how small, shocks will prevail if the Reynolds number is sufficiently large.

The squirmer model of Lighthill and Blake has been widely used to analyse swimming ciliates. However, real ciliates are covered by hair-like organelles, called cilia; the differences between the squirmer model and real ciliates remain unclear. Here, we developed a ciliate model incorporating the distinct ciliary apparatus, and analysed motion using a boundary element–slender-body coupling method. This methodology allows us to accurately calculate hydrodynamic interactions between cilia and the cell body under free-swimming conditions. Results showed that an antiplectic metachronal wave was optimal in the swimming speed with various cell-body aspect ratios, which is consistent with former theoretical studies. Exploiting oblique wave propagation, we reproduced a helical trajectory, like Paramecium, although the cell body was spherical. We confirmed that the swimming velocity of model ciliates was well represented by the squirmer model. However, squirmer modelling outside the envelope failed to estimate the energy costs of swimming; over 90 % of energy was dissipated inside the ciliary envelope. The optimal swimming efficiency was given by the antiplectic wave; the value was 6.7 times larger than in-phase beating. Our findings provide a fundamental basis for modelling swimming micro-organisms.

This paper reports on a theoretical analysis of convection in an inclined layer of mercury, a common low-Prandtl-number fluid ($Pr=0.025$). The investigation is based on the standard Oberbeck–Boussinesq equations, which are explored as a function of the inclination angle $\unicode[STIX]{x1D6FE}$ and for Rayleigh numbers $R$ in the vicinity of the convection onset. Along with the conventional Galerkin methods to study convection rolls and their secondary instabilities, we employ direct numerical simulations for fluid layers with quite large aspect ratios. It turns out that, even for small inclination angles $\unicode[STIX]{x1D6FE}\lesssim 6^{\circ }$, the secondary instabilities of the basic rolls lead either to oscillatory three-dimensional patterns or to stationary ones, which appear alternately with increasing $\unicode[STIX]{x1D6FE}$. Due to the competition of these instabilities the patterns may show a complex dynamics.

The effectiveness of streak modes in controlling the oblique-type breakdown in a supersonic boundary-layer at Mach 2.0 is investigated using direct numerical simulations. Investigations in the literature have shown the effectiveness of streak modes in delaying the onset of transition dominated by two-dimensional waves, but in oblique breakdown, three-dimensional waves and a strong streak mode dominate the transition process. Paredes et al. (J. Fluid Mech., vol. 831, 2017, pp. 524–553) discussed the possible stabilization of supersonic boundary layers by optimally growing streaks using parabolized stability equations. However, no study has as yet been reported regarding direct nonlinear control of oblique breakdown. This study deals with the effects of large-amplitude decaying streak modes generated by a blowing–suction strip at the wall to control full breakdown in a reference case. Modes with four to five times the fundamental wavenumber are found to be beneficial for controlling the transition. In the first region after the control-mode forcing, the beneficial mean-flow distortion (MFD), generated by inducing the control mode, is solely responsible for hampering the growth of the fundamental-mode. On the whole, the MFD and the three-dimensional part of the control contribute equally towards controlling the oblique breakdown. The results show significant suppression of transition, and substantial improvements have been observed in the levels of the skin-friction coefficient and wall-temperature in comparison to the uncontrolled case. Moreover, refreshing the control using an additional downstream control strip increases the gain. However, the forcing amplitude must be carefully chosen in order not to introduce a generalized inflection point in the spanwise averaged mean flow invoking enhanced disturbance growth.

Progress in roughness research, mapping any given roughness geometry to its fluid dynamic behaviour, has been hampered by the lack of accurate and direct measurements of skin-friction drag, especially in open systems. The Taylor–Couette (TC) system has the benefit of being a closed system, but its potential for characterizing irregular, realistic, three-dimensional (3-D) roughness has not been previously considered in depth. Here, we present direct numerical simulations (DNSs) of TC turbulence with sand grain roughness mounted on the inner cylinder. The model proposed by Scotti (Phys. Fluids, vol. 18, 031701, 2006) has been modified to simulate a random rough surface of monodisperse sand grains. Taylor numbers range from $Ta=1.0\times 10^{7}$(corresponding to $Re_{\unicode[STIX]{x1D70F}}=82$) to $Ta=1.0\times 10^{9}$ ($Re_{\unicode[STIX]{x1D70F}}=635$). We focus on the influence of the roughness height $k_{s}^{+}$ in the transitionally rough regime, through simulations of TC with rough surfaces, ranging from $k_{s}^{+}=5$ up to $k_{s}^{+}=92$. We analyse the global response of the system, expressed both by the dimensionless angular velocity transport $Nu_{\unicode[STIX]{x1D714}}$ and by the friction factor $C_{f}$. An increase in friction with increasing roughness height is accompanied with enhanced plume ejection from the inner cylinder. Subsequently, we investigate the local response of the fluid flow over the rough surface. The equivalent sand grain roughness $k_{s}^{+}$ is calculated to be $1.33k$, where $k$ is the size of the sand grains. We find that the downwards shift of the logarithmic layer, due to transitionally rough sand grains exhibits remarkably similar behaviour to that of the Nikuradse (VDI-Forsch., vol. 361, 1933) data of sand grain roughness in pipe flow, regardless of the Taylor number dependent constants of the logarithmic layer. Furthermore, we find that the dynamical effects of the sand grains are contained to the roughness sublayer $h_{r}$ with $h_{r}=2.78k_{s}$.

We present large-eddy simulation (LES) of flow past different airfoils with $Re_{c}$, based on the free-stream velocity and airfoil chord length, ranging from $10^{4}$ to $2.1\times 10^{6}$. To avoid the challenging resolution requirements of the near-wall region, we develop a virtual wall model in generalized curvilinear coordinates and incorporate the non-equilibrium effects via proper treatment of the momentum equations. It is demonstrated that the wall model dynamically captures the instantaneous skin-friction vector field on arbitrary curved surfaces at the resolved scale. By combining the present wall model with the stretched-vortex subgrid-scale model, we apply the wall-modelled LES approach to three different airfoil cases, spanning different geometrical parameters, different attack angles and low to high $Re_{c}$. The numerical results are verified with direct numerical simulation (DNS) at low $Re_{c}$, and validated with experiment data at higher $Re_{c}$, including typical aerodynamic properties such as pressure coefficient distributions, velocity components and also more challenging measurements such as skin-friction coefficient and Reynolds stresses. All comparisons show reasonable agreement, providing a measure of validity that enables us to further probe simulation results into aspects of flow physics that are not available from experiments. Two techniques to quantify hitherto unexplored physics of flows past airfoils are employed: one is the construction of the anisotropy invariant map, and the second is skin-friction portraits with emphasis on flow transition and unsteady separation along the airfoil surface. The anisotropy maps for all three $Re_{c}$ cases, show clearly that a portion of the flow field is aligned along the axisymmetric expansion line, corresponding to the turbulent boundary layer log-law behaviour and the appearance of turbulent transition. The instantaneous skin-friction portraits reveal a monotonic shrinking of the near wall structure scale. At $Re_{c}=10^{4}$, the interaction between the primary separation bubble and the secondary separation bubble contributes to turbulent transition, similar to the case of flow past a cylinder. At higher $Re_{c}=10^{5}$, the primary separation breaks into several small separation bubbles. At even higher $Re_{c}=2.1\times 10^{6}$, near the turbulent separation, the skin-friction lines show small-scale reversal flows that are similar to those observed in DNS of the flat plate turbulent separation. A notable feature of turbulent separation in flow past an airfoil is the appearance of turbulence structures and small-scale reversal flows in the spanwise direction due to the vortex shedding behaviour.

A theoretically based relationship for the Darcy–Weisbach friction factor $f$ for rough-bed open-channel flows is derived and discussed. The derivation procedure is based on the double averaging (in time and space) of the Navier–Stokes equation followed by repeated integration across the flow. The obtained relationship explicitly shows that the friction factor can be split into at least five additive components, due to: (i) viscous stress; (ii) turbulent stress; (iii) dispersive stress (which in turn can be subdivided into two parts, due to bed roughness and secondary currents); (iv) flow unsteadiness and non-uniformity; and (v) spatial heterogeneity of fluid stresses in a bed-parallel plane. These constitutive components account for the roughness geometry effect and highlight the significance of the turbulent and dispersive stresses in the near-bed region where their values are largest. To explore the potential of the proposed relationship, an extensive data set has been assembled by employing specially designed large-eddy simulations and laboratory experiments for a wide range of Reynolds numbers. Flows over self-affine rough boundaries, which are representative of natural and man-made surfaces, are considered. The data analysis focuses on the effects of roughness geometry (i.e. spectral slope in the bed elevation spectra), relative submergence of roughness elements and flow and roughness Reynolds numbers, all of which are found to be substantial. It is revealed that at sufficiently high Reynolds numbers the roughness-induced and secondary-currents-induced dispersive stresses may play significant roles in generating bed friction, complementing the dominant turbulent stress contribution.

We present a mesoscale kinetic model for multicomponent flows, augmented with a short range forcing term, aimed at describing the combined effect of surface tension and near-contact interactions operating at the fluid interface level. Such a mesoscale approach is shown to (i) accurately capture the complex dynamics of bouncing colliding droplets for different values of the main governing parameters, (ii) predict quantitatively the effective viscosity of dense emulsions in micro-channels and (iii) simulate the formation of the so-called soft flowing crystals in microfluidic focusers.

Conceptually, tidal rivers are seen as narrow channels along which the cross-section geometry remains constant and the bed is horizontal. As tidal waves propagate along such a channel, they decrease exponentially in height. The more rapid the decrease, the stronger the river flow. Near the coast, the tidally averaged width and depth change little throughout the year, even if the river discharge varies strongly between the seasons. However, further upstream, the water depth varies considerably with the river discharge. Recent observations from the Kapuas River, Indonesia, show that the water surface forms a backwater profile when the river flow is low. In this case, the depth converges, i.e. it gradually decreases between the river mouth and the point where the bed reaches sea level. This effect distinctly influences how tidal waves propagate up river so that their wave height does not decrease exponentially any more. We present a theoretical analysis of this phenomenon, which reveals several so far overlooked aspects of river tides. These aspects are particularly relevant to low river flow. Along the downstream part of the tidal river, depth convergence counteracts frictional damping so that the tidal range is higher than expected. Along the upstream parts of the tidal river, the low depth increases the damping so that the tide more rapidly attenuates. The point where the bed reaches sea level effectively limits the tidal intrusion, which carries over to the overtide and the subtidal water level set-up.

We investigate the hydrodynamic stability of compressible boundary layers over adiabatic walls with fluids at supercritical pressure in the proximity of the Widom line (also known as the pseudo-critical line). Depending on the free-stream temperature and the Eckert number that determines the viscous heating, the boundary-layer temperature profile can be either sub-, trans- or supercritical with respect to the pseudo-critical temperature, $T_{pc}$. When transitioning from sub- to supercritical temperatures, a seemingly continuous phase change from a compressible liquid to a dense vapour occurs, accompanied by highly non-ideal changes in thermophysical properties. Using linear stability theory (LST) and direct numerical simulations (DNS), several key features are observed. In the sub- and supercritical temperature regimes, the boundary layer is substantially stabilized the closer the free-stream temperature is to $T_{pc}$ and the higher the Eckert number. In the transcritical case, when the temperature profile crosses $T_{pc}$, the flow is significantly destabilized and a co-existence of dual unstable modes (Mode II in addition to Mode I) is found. For high Eckert numbers, the growth rate of Mode II is one order of magnitude larger than Mode I. An inviscid analysis shows that the newly observed Mode II cannot be attributed to Mack’s second mode (trapped acoustic waves), which is characteristic in high-speed boundary-layer flows with ideal gases. Furthermore, the generalized Rayleigh criterion (also applicable for non-ideal gases) unveils that, in contrast to the trans- and supercritical regimes, the subcritical regime does not contain an inviscid instability mechanism.

The fundamental mechanisms of a hitherto unstudied approach to control the crossflow-induced transition in a three-dimensional boundary layer employing unsteady control vortices are investigated by means of direct numerical simulations. Using a spanwise row of blowing/suction or volume-force actuators, subcritical travelling crossflow vortex modes are excited to impose a stabilizing (upstream) flow deformation (UFD). Volume forcing mimics the effects of alternating current plasma actuators driven by a low-frequency sinusoidal signal. In this case the axes of the actuators are aligned with the wave crests of the desired travelling mode to maximize receptivity and abate the influence of other unwanted, misaligned modes. The resulting travelling crossflow vortices generate a beneficial mean-flow distortion reducing the amplification rate of naturally occurring steady or unsteady crossflow modes without invoking significant secondary instabilities. It is found that the stabilizing effect achieved by travelling control modes is somewhat weaker than that achieved by the steady modes in the classical UFD method. However, the energy requirements for unsteady-UFD plasma actuators would be significantly lower than for steady UFD because the approach makes full use of the inherent unsteadiness of the plasma-induced volume force with alternating-current-driven actuators. Also, the input control amplitude can be lower since unsteady crossflow vortex modes grow stronger in the flow.

We investigate the growth of gas bubbles in a water solution at rest with a supersaturation level that is generally associated with diffusive mass transfer. For $\text{CO}_{2}$ bubbles, it has been previously observed that, after some time of growing in a diffusive regime, a density-driven convective flow enhances the mass transfer rate into the bubble. This is due to the lower density of the gas-depleted liquid which surrounds the bubble. In this work, we report on experiments with different supersaturation values, measuring the time $t_{conv}$ it takes for convection to dominate over the diffusion-driven growth. We demonstrate that by considering buoyancy and drag forces on the depleted liquid around the bubble, we can satisfactorily predict the transition time. In fact, our analysis shows that this onset does not only depend on the supersaturation, but also on the absolute pressure, which we corroborate in experiments. Subsequently, we study how the depletion caused by the growth of successive single bubbles influences the onset of convection. Finally, we study the convection onset around diffusively growing nitrogen $\text{N}_{2}$ bubbles. As $\text{N}_{2}$ is much less soluble in water, the growth takes much longer. However, after waiting long enough and consistent with our theory, convection still occurs as for any gas–liquid combination, provided that the density of the solution sufficiently changes with the gas concentration.

The interaction of disturbance modes behind an isolated cylindrical roughness element in a laminar boundary layer is investigated by means of hot-film anemometry and particle image velocimetry in a low-turbulence laminar water channel. Both sinuous and varicose disturbance modes are found in the wake of a roughness with unit aspect ratio (diameter/height $=$ 1). Interestingly, the frequency of the varicose mode synchronizes with the first harmonic of the sinuous mode when the critical Reynolds number from three-dimensional global linear stability theory is exceeded. The coupled motion of sinuous and varicose modes is explained by frequency lock-in. This mechanism is of great importance in many aspects of nature, but has not yet received sufficient attention in the field of boundary-layer theory. A Fourier mode decomposition provides detailed analyses of sinuous and varicose modes. The observation is confirmed by a second experiment with the same aspect ratio at a different position in the laminar boundary layer. When the aspect ratio is increased, the flow is fully governed by the varicose mode. Thus, no frequency lock-in can be observed in this case. The significance of this work is to explain how sinuous and varicose modes can co-exist behind a roughness and to propose a mechanism which is well established in physics but not encountered often in boundary-layer theory.

Multicomponent liquid drops in a host liquid are very relevant in various technological applications. Their dissolution or growth dynamics is complex. Differences in solubility between the drop components combined with the solutal Marangoni effect and natural convection contribute to this complexity, which can be even further increased in combination with the ouzo effect, i.e. the spontaneous nucleation of microdroplets due to composition-dependent miscibilities in a ternary system. The quantitative understanding of this combined process is important for applications in industry, particularly for modern liquid–liquid microextraction processes. In this work, as a model system, we experimentally and theoretically explore water–ethanol drops dissolving in anethole oil. During the dissolution, we observed two types of microdroplet nucleation, namely water microdroplet nucleation in the surrounding oil at drop mid-height, and oil microdroplet nucleation in the aqueous drop, again at mid-height. The nucleated oil microdroplets are driven by Marangoni flows inside the aqueous drop and evolve into microdroplet rings. A one-dimensional multiphase and multicomponent diffusion model in combination with thermodynamic equilibrium theory is proposed to predict the behaviour of spontaneous emulsification, i.e. microdroplet nucleation, that is triggered by diffusion. A scale analysis together with experimental investigations of the fluid dynamics of the system reveals that both the solutal Marangoni flow inside the drop and the buoyancy-driven flow in the host liquid influence the diffusion-triggered emulsification process. Our work provides a physical understanding of the microdroplet nucleation by dissolution of a multicomponent drop in a host liquid.

We present numerical and theoretical results concerning the properties of turbulent flows with strong multi-scale helical injection. We perform direct numerical simulations of the Navier–Stokes equations under a random helical stirring with power-law spectrum and with different intensities of energy and helicity injections. We show that there exists three different regimes where the forward energy and helicity inertial transfers are: (i) both leading with respect to the external injections, (ii) energy transfer is leading and helicity transfer is sub-leading and (iii) both are sub-leading and helicity is maximal at all scales. As a result, the cases (ii)–(iii) give flows with Kolmogorov-like inertial energy cascade and tuneable helicity transfers/contents. We further explore regime (iii) by studying its effect on the kinetics of point-like isotropic helicoids, particles whose dynamics is isotropic but breaks parity invariance. We investigate small-scale fractal clustering and preferential sampling of intense helical flow structures. Depending on their structural parameters, the isotropic helicoids either preferentially sample co-chiral or anti-chiral flow structures. We explain these findings in limiting cases in terms of what is known for spherical particles of different densities and degrees of inertia. Furthermore, we present theoretical and numerical results for a stochastic model where dynamical properties can be calculated using analytical perturbation theory. Our study shows that a suitable tuning of the stirring mechanism can strongly modify the small-scale turbulent helical properties and demonstrates that isotropic helicoids are the simplest particles able to preferentially sense helical properties in turbulence.

When a layer of static grains on a sufficiently steep slope is disturbed, an upslope-propagating erosion wave, or retrogressive failure, may form that separates the initially static material from a downslope region of flowing grains. This paper shows that a relatively simple depth-averaged avalanche model with frictional hysteresis is sufficient to capture a planar retrogressive failure that is independent of the cross-slope coordinate. The hysteresis is modelled with a non-monotonic effective basal friction law that has static, intermediate (velocity decreasing) and dynamic (velocity increasing) regimes. Both experiments and time-dependent numerical simulations show that steadily travelling retrogressive waves rapidly form in this system and a travelling wave ansatz is therefore used to derive a one-dimensional depth-averaged exact solution. The speed of the wave is determined by a critical point in the ordinary differential equation for the thickness. The critical point lies in the intermediate frictional regime, at the point where the friction exactly balances the downslope component of gravity. The retrogressive wave is therefore a sensitive test of the functional form of the friction law in this regime, where steady uniform flows are unstable and so cannot be used to determine the friction law directly. Upper and lower bounds for the existence of retrogressive waves in terms of the initial layer depth and the slope inclination are found and shown to be in good agreement with the experimentally determined phase diagram. For the friction law proposed by Edwards et al. (J. Fluid. Mech., vol. 823, 2017, pp. 278–315, J. Fluid. Mech., 2019, (submitted)) the magnitude of the wave speed is slightly under-predicted, but, for a given initial layer thickness, the exact solution accurately predicts an increase in the wave speed with higher inclinations. The model also captures the finite wave speed at the onset of retrogressive failure observed in experiments.

The effect of permeability heterogeneities and viscosity variations on miscible displacement processes in porous media is examined using high-resolution numerical simulations and reduced theoretical modelling. The planar injection of one fluid into a fluid-saturated, two-dimensional porous medium with a permeability that varies perpendicular to the flow direction is studied. Three cases are considered, in which the injected fluid is equally viscous, more viscous or less viscous than the ambient fluid. In general it is found that the flow in each case evolves through three regimes. At early times, the flow exhibits the concentration evolves diffusively, independent of both the permeability structure and the viscosity ratio. At intermediate times, the flow exhibits different dynamics including channelling and fingering, depending on whether the injected fluid is more or less viscous than the ambient fluid, and depending on the relative magnitude of the viscosity and permeability variations. Finally, at late times, the flow becomes independent of the viscosity ratio and dominated by shear-enhanced (Taylor) dispersion. For each of the regimes identified above, we develop reduced-order models for the evolution of the transversely averaged concentration and compare them to the full numerical simulations.