# Shear-thickening suspensions down inclines: from Kapitza to Oobleck waves

Baptiste Darbois Texier, Henri Lhuissier, Bloen Metzger, Yoel Forterre

► To cite this version:

Baptiste Darbois Texier, Henri Lhuissier, Bloen Metzger, Yoel Forterre. Shear-thickening suspensions down inclines: from Kapitza to Oobleck waves. Journal of Fluid Mechanics, 2023, 959, pp.A27. 10.1017/jfm.2023.162. hal-03763665v2

HAL Id: hal-03763665

https://hal.science/hal-03763665v2

Submitted on 21 Nov 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Shear-thickening suspensions down inclines: from Kapitza to Oobleck waves

Baptiste Darbois Texier<fnref n="1" />,<fnref n="2" />, Henri Lhuissier<fnref n="1" />, Bloen Metzger<fnref n="1" />
Yoël Forterre <fnref n="1" />†

(Received xx; revised xx; accepted xx)

We investigate experimentally and theoretically the stability of a shear-thickening suspension flowing down an inclined plane. In a previous paper (Darbois Texier et al. 2020), we have shown that for particle volume fractions $\phi$ above the discontinuous shear-thickening fraction $\phi_{DST}$, long surface waves spontaneously grow at a flow Reynolds number much below one. This motivated a simplified analysis based on a purely inertialess mechanism, called ‘Oobleck waves’ mechanism, which couples the negatively-sloped rheology of the suspension with the free-surface deflection and captures well the experimental instability threshold and the wave speed, for $\phi > \phi_{DST}$. However, neglecting inertia does not allow to describe the inertial Kapitza regime observed for $\phi < \phi_{DST}$, nor does it allow to discriminate between Oobleck waves and other inertial instabilities expected above $\phi_{DST}$. This paper fills this gap by extending our previous analysis, based on a depth-averaged approach and the Wyart-Cates constitutive shear-thickening rheology, to account for inertia. The extended analysis recovers quantitatively the experimental instability threshold in the Kapitza regime, below $\phi_{DST}$, and in the Oobleck waves regime, above $\phi_{DST}$. By providing additional measurements of the wave growth rate and investigating theoretically the effect of a strain delay in the rheology, it also confirms that the instability observed above $\phi_{DST}$ stem from the non-inertial Oobleck wave mechanism, which is specific to free-surface flows and dominates modes of inertial origin. These results emphasize the variety of instability mechanisms for shear-thickening suspensions and might be relevant to free-surface flows of other complex fluids displaying velocity-weakening rheology.

**Key words:** Shear-thickening suspension; free-surface flow; linear stability analysis;

## 1. Introduction

The resistance to flow of a shear-thickening suspension, such as an aqueous suspension of starch particles, increases steeply with increasing strain rate. Though it is not thicker than milk when it is gently stirred, the suspension may suddenly become rock-solid under high stresses or upon impact. This intriguing behavior has been puzzling scientists for more than eighty years since the first study by Freundlich & Röder (1938). It is also an important question in industry (LaFarge 2013; Abdesselam et al. 2017; Blanco et al. 2019; Zarei & Aalaie 2020), where sudden thickening or jamming of the suspension can

<footnote n="1">¹Aix-Marseille University, CNRS, IUSTI, 13453 Marseille, France.</footnote>

<footnote n="2">²Université Paris-Saclay, CNRS, FAST, 91405, Orsay, France.</footnote>

<footnote>† Email address for correspondence: yoel.forterre@univ-amu.fr.</footnote>
damage mixers or clog pipes, but it can also be harnessed to design new impact-resistant materials.

Shear-thickening arises when the suspension particles interact through a short-range repulsive force, which can stem from surface physical-chemistry effects or Brownian motion. The repulsive force implies that the contacts between the particles transition from frictionless, under a small shear stress, to frictional, when the stress is large enough. This results in a large variation in the suspension viscosity at constant volume fraction, because the jamming volume fraction of the suspension depends on the frictional state between the particles (Guazzelli & Pouliquen 2018). This frictional transition scenario, first reported by Seto et al. (2013), has been supported by discrete numerical simulations (Mari et al. 2014; Dong & Trulsson 2017; Singh et al. 2018) and experiments performed at both the contact and the flow scales (Guy et al. 2015; Lin et al. 2015; Clavaud et al. 2017; Comtet et al. 2017; Hsu et al. 2018; Clavaud et al. 2020). It has been rationalized by Wyart & Cates (2014) through a simple constitutive law assuming a stress-dependent jamming volume fraction, which successfully reproduces the different continuous shear-thickening (CST), discontinuous shear-thickening (DST) and shear-jamming (SJ) regimes observed experimentally (Guy et al. 2015; Mari et al. 2015a; Rathee et al. 2017; Morris 2018; Richards et al. 2019; Guy et al. 2020).

In particular, Wyart-Cates rheology and its later refinements (Singh et al. 2018; Richards et al. 2019; Ramaswamy et al. 2021) have a remarkable feature. Above a critical volume fraction, called $\phi_{DST}$, the flow curve becomes S-shaped, with a negatively-sloped region where the shear rate decreases with increasing stress. Such a non-monotonicity is known to promote unstable flow conditions (Yerushalmi et al. 1970; Spenley et al. 1996; Olmsted 1999; Goddard 2003; Olmsted 2008; Nakanishi & Mitarai 2011; Divoux et al. 2016), and, indeed, shear-thickening suspension flows often destabilize and grow highly unsteady and inhomogeneous structures (Boersma et al. 1991; Lootens et al. 2003; Nagahiro et al. 2013; von Kann et al. 2011; Von Kann et al. 2013; Mari et al. 2015b; Hermes et al. 2016; Rathee et al. 2017; Chacko et al. 2018; Saint-Michel et al. 2018; Richards et al. 2019; Ovarlez et al. 2020; Sedes et al. 2020; Gauthier et al. 2021). In most models, these instabilities are understood as an immediate consequence of the coupling between the S-shape rheology and inertia. Indeed, it can be shown that, if the shear rate and shear stress are instantaneously related through a decreasing flow curve, a simple shear flow is unstable along the flow direction only if inertia is taken into account (Spenley et al. 1996; Nakanishi & Mitarai 2011; Mari et al. 2015b).

Interestingly, we have reported recently an instability in the flow of a shear-thickening suspension down an inclined plane which does not rely on inertia (Darbois Texier et al. 2020). For a volume fraction above $\phi_{DST}$, long surface waves spontaneously grow, in spite of a flow Reynolds number much smaller than one. This instability was first observed by Balmforth et al. (2005) but could not be modeled at the time due to the lack of appropriate flow rule for shear-thickening suspensions. We have proposed that these waves originate from the coupling between the free-surface deformation and the negatively-sloped rheology of the suspension, when the latter is forced into the discontinuously shear-thickening region. The mechanism, which we coined ‘Oobleck waves’ instability, is specific to surface flows and does not require inertia. It actually stems from the amplification of kinematic surface waves, by a mismatch between hydrostatics and the basal stress rheology. It has been supported by a depth-averaged analysis of the flow, neglecting inertia and using Wyart-Cates rheology, which has provided predictions in fair agreement with the instability threshold and wave speed measured above $\phi_{DST}$ (Darbois Texier et al. 2020).

However, this study leaves important open questions regarding the actual role of inertia
on the formation of the surface waves. First, the non-inertial instability mechanism applies only for a volume fraction above $\phi_{\text{DST}}$, when the flow curve is negatively-sloped. Yet, growing surface waves were also observed below $\phi_{\text{DST}}$, where the flow curve is monotonous, though it was at a much larger Reynolds number than above $\phi_{\text{DST}}$ (Darbois Texier et al. 2020). These finite Reynolds number waves are certainly reminiscent of the Kapitza, or roll-waves, instability, which is observed for Newtonian (Jeffreys 1925; Kapitza & Kapitza 1948) and complex fluids, such as power law fluids (Ng & Mei 1994; Hwang et al. 1994; Allouche et al. 2017), mud (Trowbridge 1987; Liu & Mei 1994; Balmforth & Liu 2004) or granular materials (Forterre & Pouliquen 2003; Forterre 2006). In all these cases, the Kapitza instability is inertia-driven and emerges above a critical Reynolds number (or Froude number), whose value depends on the precise rheology of the fluid. Therefore, inertia must be considered to obtain a complete description of the instability, including below $\phi_{\text{DST}}$, and to understand the transition between the inertial Kapitza regime and the overdamped Oobleck wave regime. Addressing these questions represents a non-trivial test for the constitutive law of shear-thickening suspensions, which to date have been primarily confronted with steady rheological measurements.

A second issue about inertia, not addressed in our previous study, concerns its influence on the Oobleck wave instability itself, i.e., above $\phi_{\text{DST}}$. A stability analysis neglecting inertia has proven sufficient to predict the correct behavior for the instability threshold, suggesting that inertia is not involved in the instability mechanism. However, as mentioned above, even a small inertial component is known to give unstable modes for a negatively-sloped flow curve, regardless of whether the flow has a free-surface or not (Mari et al. 2015b). This raises an important fundamental question: Is the instability observed above $\phi_{\text{DST}}$ a purely non-inertial instability, resulting from the novel Oobleck wave mechanism specific to free-surface flows, or does it belong to the same class of inertial instabilities that have been reported so far for rheometric or confined shear-thickening flows (Richards et al. 2019)?

This paper addresses these questions by considering in details the role of inertia in the surface destabilization of a shear-thickening suspensions flow down an incline. §2 details the experimental setup, already used in Darbois Texier et al. (2020), and provides additional measurements of the instability growth rate, in both the dilute and the concentrated regimes. §3 presents a linear stability analysis of the flow, using depth-averaged equations, assuming homogeneous volume fraction and accounting for hydrostatic contribution, Wyart-Cates rheology and the flow inertia. The predictions of the analysis are compared to the experimental observations in §4. Finally, in §5, the results and the competition between inertial and non-inertial modes are discussed in light of a refinement of the Wyart-Cates law introducing a strain delay in the rheology (Mari et al. 2015b; Han et al. 2018; Chacko et al. 2018; Richards et al. 2019). The conclusion (§6) confirms the novelty of the instability reported in the discontinuous shear-thickening regime. Although inertial unstable modes also exist, the instability that actually emerges stem from the intrinsically non-inertial Oobleck mechanism, which is specific to free-surface flows.

## 2. Experiments

The same set-up as in Darbois Texier et al. (2020) was used to obtain complementary measurements of the instability growth rate, both below and above $\phi_{\text{DST}}$. We provide, below, more details about the set-up and the different protocols used to characterize the instability onset in the two regimes.
## 2.1. Shear-thickening suspension: composition and rheology

We use an aqueous suspension of commercial organic cornstarch (Maisita®), www.agrana.com) prepared at a volume fraction $\phi$, which is determined from the dry mass and density of the starch, $\rho_p = 1550 \text{ kg/m}^3$. The starch particles (shown in figure 1a) are polydisperse angular grains, with an average size of about 15 µm. The rheology of the suspension was characterized in Darbois Texier et al. (2020) using a cylindrical-Couette rheometer. The shear stress, $\tau$, was imposed and the shear rate, $\dot{\gamma}$, was measured to obtain, for different volume fractions, the flow curves, $\tau(\dot{\gamma})$, which are reproduced in figure 1b. The measurements are fitted with the Wyart-Cates constitutive laws (Wyart & Cates 2014). The latter assume that the effective viscosity of the suspension diverges at a critical volume fraction $\phi_J$, according to $\eta(\phi, f) = \eta_s(\phi_J(f) - \phi)^{-2}$, with $\eta_s$ a prefactor proportional to the solvent viscosity. The jamming fraction itself depends on the fraction of frictional contacts, $f$, according to $\phi_J(f) = (1-f)\phi_0 + f\phi_1$, where $\phi_0$ and $\phi_1$ are the jamming fractions for a suspension of frictionless and frictional particles, respectively. The fraction of frictional contacts is assumed to follow $f = e^{-\tau^*/\tau}$, with $\tau^*$ the critical stress scale above which frictional contacts are activated. We follow the fitting procedure of Guy et al. (2015) to fit our measurements with the model and obtain $\eta_s = 0.91 \pm 0.01 \text{ mPa s}$, $\phi_0 = 0.52 \pm 0.005$, $\phi_1 = 0.43 \pm 0.005$, and $\tau^* = 12 \pm 2 \text{ Pa}$. With these parameters, the Wyart-Cates model captures fairly well (i) the low-stress part (frictionless regime) of the rheogram for all $\phi$, (ii) the continuously shear-thickening part observed for moderate $\phi$, and (iii) the onset of discontinuous shear-thickening, i.e., the lowest stress at which the curve presents a negative slope ($d\tau/d\dot{\gamma} < 0$), for volume fractions above $\phi_{\text{DST}} \equiv \phi_0 - 2e^{-1/2}(\phi_0 - \phi_1) \simeq 0.41$. Above the threshold stress of discontinuity, the flow inside the rheometer is highly unsteady and inhomogeneous (Guy et al. 2015; Saint-Michel et al. 2018; Richards et al. 2019; Ovarlez et al. 2020; Gauthier et al. 2021), and the experimental rheogram can no longer be fitted with the model rheology.

## 2.2. Determination of the instability threshold

### 2.2.1. Experiments at low $\phi$ (Kapitza waves)

Figure 1c shows a sketch of the experimental set-up used to characterize the stability threshold at low volume fractions ($\phi < \phi_{\text{DST}}$). The set-up consists of a 1 m long and 10 cm wide plane, which can be tilted with an angle $\theta$, varied between 2° and 22°. The inclined plane is covered with a diamond lapping film of typical roughness 45 µm to prevent wall-slip. The flow is controlled by the gravity-driven drainage of a reservoir of suspension through a gate located at the top of the plane. Two low-incidence laser sheets and two cameras are used to measure the mean film thickness, $h_0 \sim 2-10$ mm, and the crest-to-crest amplitude of the waves, $\Delta h_1$ and $\Delta h_2$, at distances $x_1 = 10$ cm and $x_2 = 70$ cm from the gate. The calibration of the lasers incidence yields a precision in the local measurement of $h_0$, $\Delta h_1$ and $\Delta h_2$ of $\sim 10$ µm. The current flow rate $q$ of the suspension is measured with a scale placed at the bottom end of the incline. The current Reynolds number of the flow is computed from the current flow rate $q$ and mean film thickness $h_0$, using the relation $Re = 3q^2/(gh_0^3 \sin\theta)$, with $g$ the gravitational acceleration. This definition, which does not explicitly depend on the suspension viscosity $\eta$, is convenient since it can be used whatever the rheology of the fluid. The factor 3 is chosen so as to recover $Re = \rho \bar{u}_0 h_0 / \eta_0$, with $\bar{u}_0 = q/h_0$ the depth-averaged velocity, for a steady Newtonian flow (Landau & Lifshitz 2013).

To determine the instability threshold a small perturbation is imposed to the flow, while the flow rate decreases quasi-steadily because of the slow drainage of the reservoir
<figure><img src="image_2.png" /><figcaption>FIGURE 1. Rheograms and experiments at low $\phi$ to characterize the Kapitza instability. (a) Image of the cornstarch grains. (b) Rheograms of the aqueous cornstarch suspension for various volume fractions. Solid lines: Wyart-Cates rheology with $\eta_s = 0.91$ mPa s, $\phi_0 = 0.52$, $\phi_1 = 0.43$ and $\tau^* = 12$ Pa. The region where $d\tau/d\dot{\gamma} < 0$ is highlighted in blue. (c) Sketch of the set-up used to characterize the instability below $\phi_{\text{DST}}$ and typical picture of the Kapitza waves ($\phi = 0.33$, $\theta = 2^\circ$, $Re \approx 37$). (d) Spatio-temporal plots showing the transverse displacement of the intersection between the laser sheet and the flow surface, at the top and at the bottom of the incline ($\phi = 0.33$, $\theta = 2^\circ$, $Re \approx 37$). (e) Reynolds number of the flow, and amplitude of the perturbation at the top, $\Delta h_1$, and at the bottom, $\Delta h_2$, of the incline ($\phi = 0.36$, $\theta = 3^\circ$, $Re \approx 28$). The black dashed line indicates the instability threshold $Re_c$. (d) and (e) are reproduced from Darbois Texier et al. (2020).</figcaption></figure>

(the variation is sufficiently slow to ensure a uniform flow rate along the incline). The perturbation is forced by sinusoidally modulating the aperture of the gate (at 3 Hz, with an amplitude of $\pm 100$ µm), with the help of a translating stage. The perturbation is convected and its amplification or damping is monitored by measuring the amplitude at $x_1$ and $x_2$ (see figure 1d).

Figure 1e shows a typical evolution of the Reynolds number, $Re$, together with the wave amplitude at the top ($x_1$) and at the bottom ($x_2$) of the incline, starting from an
<figure><img src="image_3.png" /><figcaption>FIGURE 2. Experiments at high $\phi$ to characterize the Oobleck waves. (a) Sketch of the set-up used for $\phi > \phi_{DST}$ and typical image of Oobleck waves ($\phi = 0.45, \theta = 10^\circ, Re \approx 1.14 \approx 0.2Re_{Kap}$). (b) Image of the flow surface intersected by the laser sheet (same condition as in (a)). (c) Reynolds number of the flow and amplitude of the perturbation $\Delta h_1$ and $\Delta h_2$ (same condition as in (a)). The black dashed line indicates the instability threshold $Re_c$. (e) Normalized wave amplitude $\Delta h/h_0$ as a function of $x$ (same $\phi$ and $\theta$, $Re/Re_c = 1.05$). The growth rate is measured over the region highlighted in blue.</figcaption></figure>

unstable situation where $\Delta h_2 > \Delta h_1$. The instability threshold is determined from the current flow rate at the time $\Delta h_2 = \Delta h_1$, which sets the critical Reynolds number, $Re_c$ (dashed line in figure 1e), the critical flow thickness, $h_c$, the critical mean flow velocity, $u_c$, and the critical basal shear stress, $\tau_c = \rho gh_c \sin \theta$.

### 2.2.2. Experiments at high $\phi$ (Oobleck waves)

For a volume fraction above $\phi_{DST}$, the instability changes qualitatively. The perturbation is either dampened or amplified and saturated over a very short distance (~1 cm), which compares with the flow thickness, instead of increasing or decreasing gently all along the inclined plane, as for $\phi < \phi_{DST}$. Forcing the instability is not useful anymore because the most unstable modes of the perturbative noise background dominates wave formation. Moreover, for $\phi > \phi_{DST}$, it is not possible to set the flow with the draining reservoir because the jamming of the suspension at the gate creates large perturbations, which prevent studying the stability over the incline. To circumvent these issues, a modified injection system is used above $\phi_{DST}$. The suspension is discharged from a large funnel into an upper pool, which lets the discharge perturbations decay before feeding the incline by a gentle overflow (see figure 2a). To increase the suspension flow rate quasi-steadily, the funnel's aperture is opened slowly with the help of a translating stage. In this case, the wave amplitude grows over a short distance (see figure 2b), which allows characterizing the wave growth rate with a single laser-sheet and camera. Figure 2c presents the simultaneous evolution of $Re$, $\Delta h_1$ and $\Delta h_2$, as obtained with this protocole, starting from a stable situation where $\Delta h_2 < \Delta h_1$. As previously, the stability threshold is reached when $\Delta h_2 = \Delta h_1$, providing $Re_c$, $h_c$, $u_c$ and $\tau_c$.
For both protocols (above and below $\phi_{DST}$) we have verified that the same instability criteria is obtained from successive steady state measurements at various constant flow rates. For each volume fraction investigated, experiments are repeated at least four times and for each repetition a new, freshly prepared, suspension is used to avoid starch aging or evaporation issues.

### 2.2.3. Wave speed and growth rate measurements

Beside the instability threshold, three important properties characterizing the surface wave propagation are extracted from these experiments. From the measured steady state relation $q(h_0)$ between the average flow rate and the mean layer thickness, we obtain an experimental determination of the kinematic wave speed, $c_{kin} = dq/dh_0$, which will turn to be important to discriminate between the different instability mechanisms.

From the evolution of the amplitude and phase of the wave along the plane, we measure the growth rate and wave speed. It was not possible to obtain experimentally the complete dispersion relation as a function of the wave frequency, because for $\phi > \phi_{DST}$, the waves are most often dominated, within a very short distance, by the non-linear growth of the most unstable mode of the background noise, regardless of the forcing frequency. Therefore, to characterize the strength of the instability we focused the growth rate measurements on the most unstable mode just above the instability threshold, i.e., at an arbitrary distance above the threshold ($Re - Re_c/Re_c = 0.05$). For low volume fractions $\phi < \phi_{DST}$, the wave grows exponentially all along the incline. The spatial growth rate $\sigma$ is obtained from the amplitude measurements at $x_1$ and $x_2$, according to $\sigma = \ln(\Delta h_2/\Delta h_1)/(x_2 - x_1)$. For large volume fractions ($\phi > \phi_{DST}$), the amplitude of $\Delta h(x)/h_0$ saturates within a shorter distance, as shown in figures 2b & 2d. In this case, the growth rate is measured by fitting the short initial exponential regime, which is highlighted in blue. In both cases, the reported wave speed is that of the most unstable mode.

## 3. Linear stability analysis

To rationalize the instability observed experimentally, we perform a linear stability analysis of the flow. The depth-averaged approach and the approximation of homogeneous volume fraction used in Darbois Texier et al. (2020) is extended to include inertial terms. This approach has the advantage of embedding the complex rheology of the suspension in a single term, the basal stress, while not limiting significantly the scope of the analysis, since the most unstable modes will turn to be slender-sloped. The rheology of the shear-thickening suspension is modeled by the Wyart-Cates flow rule introduced in §2.1.

### 3.1. Base flow

We compute, first, the base flow, i.e., the steady uniform flow of a shear-thickening suspension, with volume fraction $\phi$, density $\rho$, and thickness $h_0$, down a plane with slope $\theta$. The base state will be denoted by the subscript 0. For a layer with a stress-free-surface, the momentum balance imposes that the shear stress $\tau_0$ increases linearly with the depth, $h_0 - z$, according to

$$ \tau_0(z) = \rho g (h_0 - z) \sin \theta. \qquad (3.1) $$

On the other hand, the shear stress is related to the shear rate by

$$ \tau_0(z) = \eta(\phi, z) \frac{du_0(z)}{dz}, \qquad (3.2) $$
<figure><img src="image_4.png" /><figcaption>FIGURE 3. (a) Sketch of the notations. (b) Velocity profiles of the base flow for Wyart-Cates' rheology with the parameters obtained from figure 1b and for various volume fractions ($\tau_{b,0}/\tau^* = 2$).</figcaption></figure>

with $u_0(z)$ the suspension velocity parallel to the plane and $\eta(\phi, z)$ the suspension viscosity, which is generally not uniform. Combining (3.1) with (3.2) yields the velocity profile in terms of the reduced variable $\tau_0$

$$u_0(\tau_0) = \frac{1}{\rho g \sin \theta} \int_{\tau_0}^{\tau_{b,0}} \frac{\tau'}{\eta(\phi, \tau')} d\tau', \quad (3.3)$$

where $\tau_{b,0} \equiv \tau_0(z=0) = \rho g h_0 \sin \theta$ is the basal shear stress. Finally, the viscosity is given by Wyart-Cates' expression

$$\eta(\phi, \tau) = \eta_s \left[ \phi_0 (1 - e^{-\tau^*/\tau}) + \phi_1 e^{-\tau^*/\tau} - \phi \right]^{-2}, \quad (3.4)$$

where $\eta_s, \tau^*, \phi_0$ and $\phi_1$ are the rheological parameters introduced in §2.1.

Figure 3 shows the base flow velocity profile obtained by numerically integrating (3.3) using (3.4), for volume fractions between 0.30 and 0.48. For low $\phi$, the velocity profile is semi-parabolic, as expected for a Newtonian fluid. For increasing $\phi$, the concavity of the profile reverses, which reflects the increase in the suspension viscosity at the bottom of the layer where the stress is the largest. The flowing region even localizes close to the free-surface when the suspension jams beneath, i.e., when the basal stress reaches

$$\tau_{b,SJ} = \frac{\tau^*}{\ln\left(\frac{\phi_0 - \phi_1}{\phi_0 - \phi}\right)}. \quad (3.5)$$

Note that the vertical gradient of shear is expected to drive particle migration from the bottom to the top of the layer, which in turn should slightly modify the velocity profile (Carpen & Brady 2002; Dhas & Roy 2022). For simplicity, we do not consider this coupling between flow and volume fraction variation, which is not essential to account for the instabilities studied here.

In the following of the analysis, we will use depth-averaged quantities and restrict the calculations to $\tau_{b,0} < \tau_{b,SJ}$, which does not affect the flow stability prediction (see figures 7a and 9c). From (3.3), the depth-averaged velocity of the base flow, $\bar{u}_0 = \frac{1}{h_0} \int_0^{h_0} u_0(z) dz = \frac{1}{\tau_{b,0}} \int_0^{\tau_{b,0}} u_0(\tau) d\tau$, is given by

$$\bar{u}_0 = \frac{h_0}{\tau_{b,0}^2} \int_0^{\tau_{b,0}} \int_{\tau}^{\tau_{b,0}} \frac{\tau'}{\eta(\phi, \tau')} d\tau' d\tau. \quad (3.6)$$
This expression can be recast into a formal effective rheological law relating the basal stress, $\tau_{b,0}$, with the effective shear rate, $\bar{u}_0/h_0$, as follows

$$ \frac{\bar{u}_0}{h_0} = \frac{\tau_{b,0}}{3\eta(\phi, \tau_{b,0})} \mathcal{G}(\phi, \tau_{b,0}), \qquad (3.7) $$

where the function $\mathcal{G}$ is defined as

$$ \mathcal{G}(\phi, \tau_b) = \frac{3\eta(\phi, \tau_b)}{\tau_b^3} \int_0^{\tau_b} \int_{\tau}^{\tau_b} \frac{\tau'}{\eta(\phi, \tau')} d\tau' d\tau. \qquad (3.8) $$

For a Newtonian fluid with a uniform viscosity, $\mathcal{G} = 1$. One recovers the basal stress relation for a steady uniform Newtonian flow, $\tau_{b,0} = 3\eta\bar{u}_0/h_0$, where the factor 3 is a signature of the semi-parabolic velocity profile. For Wyart-Cates' shear-thickening law (3.4), $\mathcal{G}$ is no longer constant and depends on both the volume fraction, $\phi$, and the relative basal stress, $\tau_{b,0}/\tau^*$.

Finally, the Reynolds number of the base flow is given by

$$ Re = \frac{3\bar{u}_0^2}{gh_0 \sin \theta} = \frac{\tau_{b,0}^3 \mathcal{G}(\phi, \tau_{b,0})^2}{3\eta(\phi, \tau_{b,0})^2 \rho g^2 \sin^2 \theta}. \qquad (3.9) $$

The latter depends on three of the four main dimensionless parameters of the problem, namely, the Reynolds number based on the suspending liquid viscosity, $Re_s = \tau_{b,0}^3/(3\eta_s^2\rho g^2 \sin^2 \theta)$, the volume fraction, $\phi$, and the magnitude of the basal shear stress relative to the repulsive stress, $\tau_{b,0}/\tau^*$, two of which being controlled by the flow thickness $h_0$. The fourth parameter is the inclination angle $\theta$.

The base state flow rule (3.7-3.8) summarizes the rheological behavior of the suspension flow. It will be used in the following to study stability.

## 3.2. Depth-averaged equations

To study flow stability, we take advantage of the longness of the observed waves, whose wavelength ($\sim 10$ cm) is much larger than the layer thickness ($h_0 \lesssim 1$ cm). In this long wave limit, the vertical momentum balance implies that the pressure distribution is hydrostatic to the lowest order and that horizontal viscous stress gradients can be neglected. Integrating the mass and horizontal momentum equations across the flow, for an incompressible medium, yields the depth-averaged, or Saint-Venant, equations

$$ \frac{\partial h}{\partial t} + \frac{\partial h \bar{u}}{\partial x} = 0, \qquad (3.10) $$

$$ \rho \left( \frac{\partial h \bar{u}}{\partial t} + \frac{\partial h \bar{u}^2}{\partial x} \right) = \rho g h \sin \theta - \tau_b - \rho g h \cos \theta \frac{\partial h}{\partial x}, \qquad (3.11) $$

with $h(x, t)$ the flow thickness, $\bar{u}(x, t) = (1/h) \int_0^h u(x, z, t) dz$ the depth-averaged velocity, $\bar{u}^2 = (1/h) \int_0^h u^2(x, z, t) dz$ the averaged square velocity, and $u(x, z, t)$ the parallel velocity component. The right-hand terms in (3.11) correspond to the gravity term, the basal shear stress, and the resultant of the horizontal gradient of hydrostatic pressure, from left to right. To derive this last term, the normal stress tensor of the fluid is assumed isotropic at the lowest order.

To solve the system, closure relations are required for the basal stress $\tau_b$ and momentum flux term $\bar{u}^2$. Following a common approach in roll wave studies (Kapitza & Kapitza 1948; Trowbridge 1987; Ng & Mei 1994; Forterre & Pouliquen 2003), we assume that the base state flow rule (3.7-3.8), derived for a steady uniform flow, remains valid for an unsteady,
non-uniform flow in the long wavelength limit, which implies

$$ \frac{\bar{u}}{\bar{h}} = \frac{\tau_b}{3\eta(\phi, \tau_b)} \mathcal{G}(\phi, \tau_b) \equiv \dot{\gamma}(\tau_b). \qquad (3.12) $$

Similarly, we rewrite the momentum flux term as $\bar{u}^2 = \alpha \bar{u}^2$, and assume that the factor $\alpha$, which is set by the shape of the velocity profile, is constant and equal to the base state value. From equation (3.3), we obtain

$$ \alpha = \frac{\tau_{b,0} \int_{0}^{\tau_{b,0}} \left( \int_{\tau}^{\tau_{b,0}} \frac{\tau'}{\eta(\phi, \tau')} d\tau' \right)^2 d\tau}{\left( \int_{0}^{\tau_{b,0}} \int_{\tau}^{\tau_{b,0}} \frac{\tau'}{\eta(\phi, \tau')} d\tau' d\tau \right)^2}. \qquad (3.13) $$

For a Newtonian fluid (uniform viscosity), $\alpha = 6/5$. This value increases as the flow localizes closer and closer beneath the surface. We will see that the instability threshold can be significantly shifted by the value of $\alpha$ at large volume fractions.

### 3.3. Linearization

To analyse the linear stability of the base state flow we non-dimensionalize equations using $\tilde{h} = h/h_0$, $\tilde{x} = x/h_0$, $\tilde{u} = \bar{u}/\bar{u}_0$, $\tilde{t} = t \bar{u}_0/h_0$, $\tilde{\tau}_b = \tau_b/\tau_{b,0}$ and $\tilde{\gamma} = \dot{\gamma} h_0/\bar{u}_0$. The conservation equations and flow rule (3.10-3.12) become

$$ \frac{\partial \tilde{h}}{\partial \tilde{t}} + \frac{\partial \tilde{h} \tilde{u}}{\partial \tilde{x}} = 0, \qquad (3.14) $$

$$ \frac{Re}{3} \left( \frac{\partial \tilde{h} \tilde{u}}{\partial \tilde{t}} + \alpha \frac{\partial \tilde{h} \tilde{u}^2}{\partial \tilde{x}} \right) = \tilde{h} - \tilde{\tau}_b - \frac{\tilde{h}}{\tan \theta} \frac{\partial \tilde{h}}{\partial \tilde{x}}, \qquad (3.15) $$

$$ \frac{\tilde{u}}{\tilde{h}} = \tilde{\gamma}(\tilde{\tau}_b), \qquad (3.16) $$

where $Re$ is given by (3.9).

Considering a small perturbation of the base flow, $\tilde{h} = 1 + h_1$, $\tilde{u} = 1 + u_1$, $\tilde{\tau}_b = 1 + \tau_1$, with $|h_1|, |u_1|, |\tau_1| \ll 1$, equations (3.14-3.16) become, at the lowest order,

$$ \frac{\partial h_1}{\partial \tilde{t}} + \frac{\partial h_1}{\partial \tilde{x}} + \frac{\partial u_1}{\partial \tilde{x}} = 0, \qquad (3.17) $$

$$ \frac{Re}{3} \left( \frac{\partial u_1}{\partial \tilde{t}} + (\alpha - 1) \frac{\partial h_1}{\partial \tilde{x}} + (2\alpha - 1) \frac{\partial u_1}{\partial \tilde{x}} \right) = h_1 - \tau_1 - \frac{1}{\tan \theta} \frac{\partial h_1}{\partial \tilde{x}}, \qquad (3.18) $$

$$ u_1 - h_1 = A\tau_1, \qquad (3.19) $$

where $A$ is defined as

$$ A \equiv \left( \frac{d\tilde{\gamma}}{d\tilde{\tau}_b} \right)_{\tilde{\tau}_b=1} = \frac{\tau_{b,0}h_0}{\bar{u}_0} \left( \frac{d\tilde{\gamma}}{d\tilde{\tau}_b} \right)_{\tilde{\tau}_b=\tau_{b,0}} = \frac{3}{\mathcal{G}(\phi, \tau_{b,0})} - 2, \qquad (3.20) $$

and use has been made of the identity $\frac{d}{d\tau} \left( \int_0^\tau \int_{\tau'}^\tau \frac{\tau''}{\eta(\phi, \tau'')} d\tau'' d\tau' \right) = \frac{\tau^2}{\eta(\phi, \tau)}$. The parameter $A$ represents the dimensionless inverse slope of the flow rule between the effective shear rate $\tilde{u}/h$ and the basal stress $\tilde{\tau}_b$. For a shear-thickening suspension following the Wyart-Cates' flow rule, $A$ depends on $\phi$ and $\tau_{b,0}/\tau^*$. It is equal to 1, for a Newtonian flow, and negative for discontinuous shear-thickening.

Overall, the linearized system of equations (3.17-3.19) involves four dimensionless parameters, $\theta$, $Re$, $\alpha$ and $A$ (which are alternative to those listed above $\theta$, $Re_s$, $\tau_{b,0}/\tau^*$ and $\phi$).
### 3.4. Modes and stability diagram

The system (3.17-3.19) is solved for a normal mode, $h_1 = He^{i(k\tilde{x}-\tilde{\omega}\tilde{t})}$, $u_1 = Ue^{i(k\tilde{x}-\tilde{\omega}\tilde{t})}$, with dimensionless wave number $\tilde{k}$ and dimensionless pulsation $\tilde{\omega}$. A non-trivial solution exists only if

$$ \det \left( i(\tilde{k} - \tilde{\omega}) \frac{i}{\tan \theta} \tilde{k} + \frac{Re}{3}(\alpha - 1)i\tilde{k} - \left(1 + \frac{1}{A}\right) \frac{Re}{3}(i\tilde{k}(2\alpha - 1) - i\tilde{\omega}) + \frac{1}{A} \right) = 0, \quad (3.21) $$

which provides the dispersion relation

$$ -\frac{Re}{3}\tilde{\omega}^2 + \left(\frac{2Re}{3}\alpha\tilde{k} - \frac{i}{A}\right)\tilde{\omega} + \left(\frac{1}{\tan\theta} - \frac{Re\alpha}{3}\right)\tilde{k}^2 + \left(1 + \frac{2}{A}\right)i\tilde{k} = 0. \quad (3.22) $$

We conduct the temporal stability analysis with $\tilde{k}$ real and $\tilde{\omega}$ complex. Equation (3.22) is of order 2 in $\omega$ and has two branches. Each of them may actually embed different instabilities depending on point of the phase space considered. To get insight into the physical meaning and stability of the branches, it is instructive to study their behavior at low $\tilde{k}$, before giving the exact solutions. The structure of the dispersion relation ensures that the growth rate $\tilde{\sigma} = \Im[\tilde{\omega}(\tilde{k})]$ is monotonous and does not change sign with $\tilde{k}$, which means that the stability criterion at low $\tilde{k}$ is valid for all wavenumbers. Expanding the pulsation as $\tilde{\omega} = ia_0 + ck + ia_2\tilde{k}^2$ in the dispersion relation (3.22) gives the two following solutions at the lowest order in $\tilde{k}$

$$ \tilde{\omega}_1 \approx (2+A)\tilde{k} + iA \left[ \frac{Re}{3} \left[ (2+A)(2+A-2\alpha) + \alpha \right] - \frac{1}{\tan\theta} \right] \tilde{k}^2, \quad (3.23) $$

$$ \tilde{\omega}_2 \approx (2\alpha - 2 - A)\tilde{k} - i \frac{3}{ARe}. \quad (3.24) $$

The first branch, $\tilde{\omega}_1(\tilde{k})$, is the 'kinematic' branch, since its wave speed in the long wave limit ($\tilde{k} \to 0$), $\tilde{c}_1 = \Re(\tilde{\omega}_1)/\tilde{k} = 2+A$, is that of kinematic waves, i.e., the slender small amplitude waves, which propagate at the speed $c_{kin} = (dq/dh)_0 = \bar{u}_0+h_0(\bar{d}\bar{u}/dh)_0$, obtained by combining the steady flow rule $\bar{u}(h)$ with the mass equation (3.10) (Whitham 2011). Indeed, $c_{kin}$ can be expressed in terms of $A$ by noting that $\bar{u}(h)$ satisfies the force balance $\tau_b[\bar{u}(h)/h] = \rho gh \sin\theta$, in the base state. Differentiating with respect to $h$ and making use of the definition of $A$ (3.20), one recovers $\tilde{c}_{kin} = 1+(h_0/\bar{u}_0)(\bar{d}\bar{u}/dh)_0 = 2+A$.

The kinematic branch $\tilde{\omega}_1(\tilde{k})$ is unstable when the growth rate $\tilde{\sigma}_1 = \Im(\tilde{\omega}_1)$ is positive. Depending on the sign of $A$, two cases must be considered, which will be shown to concern two different instabilities. For $A > 0$, i.e., when the effective rheology (3.12) is monotonous, the kinematic branch is unstable for large Reynolds numbers

$$ Re > Re_{\text{Kap}} = \frac{3}{[(2+A)(2+A-2\alpha) + \alpha] \tan \theta}, \quad (3.25) $$

which extends the classical inertial Kapitza instability criteria to the shear-thickening rheology. In the Kapitza regime ($A > 0$), inertia introduces a lag, which tends to amplify kinematic waves, while gravity tends to spread and stabilize them. The instability arises when the speed of kinematic waves is larger than the speed of gravity waves (Whitham 2011). For a Newtonian fluid ($A = 1$ and $\alpha = 6/5$), the threshold of the Kapitza instability predicted by (3.25) is $1/\tan\theta$, which slightly overestimates the exact prediction $Re_{Kap,Newt} = (5/6)\tan\theta$ obtained from a rigorous long wave expansion of the Navier-Stokes equations (Benjamin 1957; Yih 1963). This well-documented discrepancy stem from assuming a fixed shape of the velocity profile. For a continuously shear-thickening suspension ($0 < A < 1$), (3.25) predicts an increase in the critical Reynolds number
relative to the Newtonian case. This is consistent with previous studies on power-law rheology fluids, which have shown that shear-thickening has a stabilizing effect on the flow (Ng & Mei 1994; Hwang et al. 1994).

For $A < 0$, i.e, when the effective flow rule (3.12) becomes negatively-sloped, the stability condition is reversed. The kinematic branch is unstable for

$$Re < Re_{\text{Kap}} = \frac{3}{[(2 + A)(2 + A - 2\alpha) + \alpha] \tan \theta}, \quad (3.26)$$

which means, surprisingly, that the kinematic branch is unstable at low Reynolds number, while inertia has now a stabilizing effect. This low Reynolds number instability, appearing for a negatively-sloped flow rule ($A < 0$), corresponds to the mechanism of formation of the Oobleck waves proposed by Darbois Texier et al. (2020). Indeed, in the limit of vanishing inertia ($Re = 0$), the dispersion relation (3.22) reduces to

$$\tilde{\omega} = (2 + A)\tilde{k} - \frac{A}{\tan \theta}i\tilde{k}^2, \quad (3.27)$$

or, equivalently, in the spatio-temporal domain

$$\frac{\partial h_1}{\partial t} + (2 + A) \frac{\partial h_1}{\partial x} = \frac{A}{\tan \theta} \frac{\partial^2 h_1}{\partial x^2}. \quad (3.28)$$

One recognizes an advection-diffusion equation for the perturbative wave $h_1$, which predicts that waves propagate at the speed of kinematic waves $\tilde{c}_{kin} = 2 + A$, while diffusing with an effective diffusion coefficient $A/\tan \theta$. For $A < 0$, waves anti-diffuse, i.e., grow during propagation. As discussed in Darbois Texier et al. (2020), this instability can be understood, physically, as follows. In the absence of inertia, the balance of forces (3.11) between the gravity term, the basal stress and the pressure term implies that a locally positive (resp. negative) slope of the free-surface causes a decrease (resp. an increase) in the basal stress. Because of the negative slope of the flow rule ($A < 0$), the basal stress variation induces anti-correlated velocity variations (positive upstream of a bump and negative downstream), which amplify the initial perturbation.

The analysis above confirms that, although they are both kinematic modes, the extended Kapitza instability ($A > 0$) and Oobleck waves ($A < 0$) are fundamentally different. For the latter, the destabilizing mechanism is non-inertial and inertia has only a stabilizing effect, which stabilizes high Reynolds number flow.

The second branch $\tilde{\omega}_2(\tilde{k})$, with growth rate $\tilde{\sigma}_2 = \Im(\tilde{\omega}_2) = -3i/(ARe)$, is unstable only if $A < 0$, regardless of Reynolds number. The condition on $A$ is the same as for Oobleck waves. However, the instability mechanism is, once again, fundamentally different. For the second branch, any perturbation is amplified when $A < 0$, independently of whether a free-surface is present or not, because inertia introduces a mismatch between the basal stress and the driving gravity force. The branch is not specific to free-surface flows and disappears in the strict absence of inertia ($Re = 0$). For this reason, we call it ‘inertial branch’.

The two critical curves, $A = 0$ and $Re = Re_{Kap}$, lead to the stability diagram shown in figure 4, for an arbitrary plane inclination $\theta = 10^\circ$. For the sake of simplicity, the predictions are plotted for a fixed value of $\alpha (= 1$, corresponding to a plug velocity profile). This simplification permits a two-dimensional representation, without altering the stability diagram, qualitatively. Note that the assumption of $\alpha = 1$ is only made in figure 4 while the rest of the analysis considers the exact value of $\alpha$ obtained from (3.13). In this case, the critical Reynolds number of the kinematic branch reduces to $Re_{Kap} = 3/[(1+A)^2 \tan \theta]$ (dark solid line in figure 4a). As discussed above, the extended Kapitza instability develops for $A > 0$ and $Re > Re_{Kap}$, Oobleck waves for $A < 0$ and
<figure><img src="image_5.png" /><figcaption>FIGURE 4. Stability diagram ($Re, A$) for (a) the kinematic branch and (b) the inertial branch ($\theta = 10^{\circ}$ and a plug flow profile ($\alpha = 1$) is assumed for simplicity, see text). (a) Black line: Kapitza instability threshold ($Re = Re_{Kap}$). Red line: Oobleck waves instability threshold ($A = 0$). (b) Dashed-blue line: inertial branch instability threshold ($A = 0$). (a-b) The green line indicates the Newtonian case ($A = 1$). The coloured trajectories indicate the evolution of $Re$ and $A$ for various volume fractions and increasing flow rates (or basal stress $\tau_{b,0}$) ($\theta = 10^{\circ}$ and the rheological parameters are those obtained from figure 1b). For most volume fractions above $\phi_{DST}$, the discontinuous shear-thickening condition $Re = Re_{A=0}$ (i.e., $A = 0$) is expected to be reached before (lower flow rate) the Kapitza instability onset ($Re = Re_{Kap}$).</figcaption></figure>

$Re < Re_{Kap}$, whereas the inertial branch, shown in figure 4b, is unstable for $A < 0$ and $Re > 0$.

To determine which criterion is reached first and what instability is expected to be observed in practice, it is crucial to understand how $Re$ and $A$ vary in experiments given their coupled dependance on $\tau_{b,0}/\tau^*$, $\theta$ and $\phi$. To this end, we display in figure 4 the trajectories followed by $A$ and $Re$ for an increasing flow rate (i.e., increasing $\tau_{b,0}/\tau^*$ or flow thickness) and a fixed angle ($\theta = 10^{\circ}$), which mimics the experimental protocol. The different trajectories correspond to different volume fractions and the rheological parameters are those measured for the cornstarch suspensions (see §2.1). Below $\phi_{DST}$, the trajectories only cross the $Re = Re_{Kap}$ critical line, since $A$ remains strictly positive for all flow rates. This means that the Kapitza instability is expected, provided the flow rate is increased sufficiently. By contrast, above $\phi_{DST}$, one can, a priori, expect either the Kapitza instability or one of the two other instabilities (Oobleck wave and inertial branch), depending on which criterion ($A = 0$ or $Re = Re_{Kap}$) is reached first when the flow rate is increased. This condition is given by the respective value of the two Reynolds numbers defined by

$$ Re_{A=0} = \frac{3}{4} \frac{\tau_{b,A=0}^3}{\rho[g\eta(\phi, \tau_{b,A=0}) \sin\theta]^2}, \quad \text{with} \quad \mathcal{G}(\phi, \tau_{b,A=0}) = \frac{3}{2}, \quad (3.29) $$

corresponding to the intersection of the iso-$\phi$ trajectory with the vertical axis $A = 0$ (purple circle in figure 4), and

$$ Re_{Kap,A=0} = \frac{3}{[(4 - 3\alpha) \tan \theta]}, \quad (3.30) $$

corresponding to the intersection between the Kapitza threshold and the vertical axis $A = 0$ (black circle in figure 4). If $Re_{A=0} < Re_{Kap,A=0}$, as in figure 4, the trajectory intersects the $A = 0$ criteria first, meaning that Oobleck waves and inertial branches are expected to be observed first, for an increasing flow rate. In the opposite case ($Re_{A=0} > Re_{Kap,A=0}$) the trajectory first encounters the Kapitza threshold (with $A$ still positive) and the
<figure><img src="image_6.png" /><figcaption>FIGURE 5. Dispersion relation of the two instabilities of the kinematic branch ($\theta = 10^\circ$). Left: Kapitza instability ($\phi = 0.33$). Right: Oobleck waves instability ($\phi = 0.45$). (a) Temporal growth rate $\Im(\tilde{\omega}_1)$. Note the highly different scales of the y-axis. (b) Wave speed $\tilde{c}_1$. (a-b) The yellow circles indicate the theoretical growth rates, for $Re = 1.05Re_c$ and at the wavelength observed experimentally, which are compared with measurements in figure 8b.</figcaption></figure>

Kapitza instability is expected to develop first. The above condition between $Re_{A=0}$ and $Re_{Kap,A=0}$ non-trivially involves $\phi$, the rheological parameters and the inclination angle $\theta$. However, as figure 4 shows, for cornstarch and provided the plane remains far from the vertical ($\theta \ll 90^\circ$), the onset of discontinuous shear-thickening ($A = 0$) is reached before the Kapitza threshold ($Re_{Kap}$) for almost all volume fractions above $\phi_{DST}$.

In the following, the value of the Reynolds number when the first instability criterion is met for increasing flow rate and a fixed angle (i.e. following the iso-$\phi$ trajectories in figure 4) will be denoted by $Re_c$, in order to match the experimental definition. In practice, for our range of parameters (cornstarch rheology, plane far from vertical), $Re_c = Re_{Kap}$ for $\phi < \phi_{DST}$ and $Re_c = Re_{A=0}$ for $\phi > \phi_{DST}$.

### 3.5. Dispersion relation

The previous analysis has focused on the limit of vanishing $k$. We solve, now, the dispersion relation exactly for an arbitrary wavenumber. The two solutions of (3.22) are

$$ \tilde{\omega}_{1,2} = \frac{\frac{2Re}{3} \alpha \tilde{k} \pm \left( \frac{D+\sqrt{D^2+C^2}}{2} \right)^{1/2}}{\frac{2Re}{3}} + i \frac{-\frac{1}{A} \pm \frac{C}{2} \left( \frac{2}{D+\sqrt{D^2+C^2}} \right)^{1/2}}{\frac{2Re}{3}}, \quad (3.31) $$

where $\pm$ stands for + for the kinematic branch $\tilde{\omega}_1$ and - for the inertial branch $\tilde{\omega}_2$, with

$$ C = \frac{4Re}{3} \left( 1 + \frac{2-\alpha}{A} \right) \tilde{k}, \quad \text{and} \quad D = -\frac{1}{A^2} + \frac{4Re}{3} \left( \frac{Re}{3} \alpha(\alpha-1) + \frac{1}{\tan\theta} \right) \tilde{k}^2. \quad (3.32) $$

Figure 5a presents the growth rate for the two instabilities of the kinematic branch, namely, the Kapitza instability ($\phi = 0.33 < \phi_{DST}$, left panel) and the Oobleck wave instability ($\phi = 0.45 > \phi_{DST}$, right panel). The different colors stand for increasing values of $Re$ close to $Re_c$ (i.e., either $Re_{Kap}$ or $Re_{A=0}$). For both instabilities the growth
<figure><img src="image_7.png" /><figcaption>FIGURE 6. Dispersion relation of the inertial branch ($\theta = 10°$). Left: $\phi = 0.33$. Right: $\phi = 0.45$. (a) Temporal growth rate $\Im(\tilde{\omega}_1)$. (b) Wave speed $c_2$. The yellow circles indicate the theoretical growth rates, for $Re = 1.05Re_c$ and at the wavelength observed experimentally, which are compared with measurements in figure 8b.</figcaption></figure>

rate is null for $\tilde{k} = 0$ and increases monotonically up to a plateau value at large $\tilde{k}$, which is the signature of a zero wavenumber instability. Interestingly, for a given wavenumber and distance to the threshold, the growth rate is several orders of magnitude larger for Oobleck waves than for the Kapitza instability (the vertical scale between the two panels differs by a factor $10^3$). Figure 5b presents the wave speed of the two instabilities for the same parameters. Close to the instability threshold the wave speed depends only weakly on $\tilde{k}$. It drops from about 3 for the Kapitza instability to 2 for Oobleck waves, in agreement with the long wave limit $2+A$, since $A \approx 1$ at the Kapitza threshold and $A=0$ at the Oobleck waves threshold.

The growth rate of the inertial branch is shown in figure 6a for the same set of $\phi$ and $Re/Re_c$ as previously. We recover that the branch is unconditionally stable below $\phi_{DST}$ and unstable above $\phi_{DST}$ for Reynolds numbers larger than $Re_{A=0}$. By contrast with the kinematic branch, the growth rate is non-null at $\tilde{k} = 0$ and actually diverges at the instability threshold ($A = 0$), where it changes sign, before decreasing with increasing $Re$ above the threshold. The corresponding wave speeds are presented in figure 6b. For $\phi > \phi_{DST}$, the wave speed at threshold is lower than for the kinematic branch. Although the kinematic and inertial branches share the same criterion of stability for concentrated suspensions ($A = 0$), the comparison between figures 5 and 6 suggests that they strongly differ in terms of growth and propagation speed. We address this point in the next section, where we compare predictions with experiments.

# 4. Comparison with experiments

## 4.1. Stability threshold

We compare, first, the stability criteria derived above with the wave onset conditions observed in experiments. Figure 7a presents $Re_c/Re_{Kap,Newt}$, i.e., the critical Reynolds number normalized by that for a Newtonian liquid ($Re_{Kap,Newt} = (5/6)/\tan\theta$, see §3.4),
<figure><img src="image_8.png" /></figure>

FIGURE 7. Destabilization threshold: comparison with experiments. (a) Critical Reynolds number $Re_c$ normalized by the Newtonian Kapitza threshold $Re_{Kap,Newt} = (5/6)/\tan\theta$ versus volume fraction $\phi$. (b) Critical basal shear stress $\tau_c$ versus $\phi$. The symbols indicate measurements. Their shape encodes the inclination $\theta$: ◇ 2°, ▽ 3°, ▷ 6°, ◁ 9°, ◦ 10°, △ 22°. The error bars indicate the standard deviation of the measurements at a given $\phi$. The lines represent the predictions for the different modes (as labeled above the graph), i.e., the critical Reynolds number at which each iso-$\phi$ trajectory (see figure 4) reaches the $Re = Re_{Kap}$ or $Re = Re_{A=0}$ condition (the black line is interrupted above $\phi \approx 0.42$ since no steady flow verifies $Re = Re_{Kap}$). The shear-jamming limit corresponds to the value of the basal shear stress when the flow first jams at $z=0$ (see eq. 3.5). Inset: basal shear stress $\tau_{b,0} = \rho gh_0 \sin\theta$ versus mean shear rate $u_0/h_0$ for different volume fractions obtained from Wyart-Cates rheological laws (dashed lines) and from direct measurements on the inclined plane (plain lines). The red line and black crosses highlight the condition $d\tau_b/d\dot{\gamma} = 0$ for both set of curves.

as a function of $\phi$, whereas figure 7b reports the critical basal shear stress $\tau_{b,c}$, also versus $\phi$. The symbols represent the experimental measurements for various volume fractions. Each one is obtained by varying the flow rate at a fixed inclination of the plane (encoded by the shape of the symbol), as detailed in §2. The solid lines represent the theoretical predictions for the same protocol, i.e., the critical Reynolds number $Re_c$ at which each iso-$\phi$ trajectory (see figure 4) reaches the $Re = Re_{Kap}$ or $Re = Re_{A=0}$ condition. The relative threshold $Re_c/Re_{Kap,Newt}$ is close to 1 at low volume fraction, where shear-thickening is mild, and increases with increasing $\phi$ to reach $\approx 6$ at $\phi = 0.41 \approx \phi_{DST}$. This illustrates the significant stabilization effect of continuous shear-thickening in the Kapitza regime, in fair agreement with the evolution of $Re_{Kap}$ predicted by equation (3.25) (solid black line in figure 7a). Similarly, the steep increase in the critical stress $\tau_{b,c}$, which is observed experimentally close to $\phi = 0.41$ (figure 7b), agrees with the expected divergence of the effective viscosity coefficient of the flow, $d\tau_{b,0}/d\dot{\gamma} \propto A^{-1}$, in $\phi_{DST} \approx 0.41$. This confirms the Kapitza-like nature of the instability (destabilizing inertia vs stabilizing gravity) below $\phi_{DST}$.

Above $\phi_{DST}$, both the critical Reynolds number and the critical shear stress observed experimentally drop drastically, by up to two orders of magnitude, for $Re_c$, at $\phi = 0.47$ (figure 7a). The drop in $Re_c$ is correctly captured by the theoretical prediction $Re_c = Re_{A=0}$ of (3.29), which, once again, applies to both Oobleck waves and the inertial waves instability, and does not permit to distinguish between them. The agreement is also reasonable when the instability threshold is expressed in terms of $\tau_{b,c}$, although the prediction underestimates the measured value by about a factor 2 (figure 7b). To
<figure><img src="image_9.png" /><figcaption>FIGURE 8. Wave speed and growth rate: comparison with experiments. (a) Wave speed $c_c$ normalized by the mean fluid velocity $u_c$ at the instability onset versus $\phi$. Inset: normalized wave velocity $c_c/u_c$ versus normalized speed of kinematic waves $c_{kin}/u_c$ at the onset. (b) Normalized spatial growth rate $\tilde{\sigma}$ versus $\phi$ for $Re/Re_c = 1.05$.</figcaption></figure>

clarify this discrepancy, we note that the theoretical prediction relies on the value of the rheological parameters ($\eta_s$, $\phi_0$, $\phi_1$ and $\tau^*$) as obtained from the cylindrical Couette rheometry (see §2.1). In the inset of figure 7b, we test these parameters more directly versus the incline flow configuration by comparing the depth-averaged flow rule, $\tau_{b,0}$ vs $u_0/h_0$, they predict (3.7) with the one measured directly in the experiments on the incline. A significant difference is observed between the two flow rules, showing that the steady uniform flow down the incline is not well predicted from the rheological parameters obtained with the Couette rheometer. Such a difference has been reported previously in the case of non-shear-thickening suspensions (Bonoit et al. 2010) and could arise from the modification of the velocity profile due to particle migration effects, which are not considered here (see §3.1). Remarkably however, when the expression for the critical shear stress $\tau_{b,c}$ is computed from the flow rule measured with the inclined plane (i.e., from the points $d\dot{\gamma}/d\tau_{b,0} = 0$ highlighted by the black crosses in the inset of figure 7b), the agreement between the theoretical predictions and measurements becomes quantitative (red blue line versus crosses in figure 7b). This suggests that the mild quantitative discrepancy between theory and experiments for $\tau_c$ does not stem from a limitation of the linear stability analysis but rather from the calibration of the base flow itself.

The comparison above confirms that the onset of a negatively-sloped flow rule ($A = 0$) is, experimentally, the condition for flow stability above $\phi_{DST}$. However, it does not allow to determine which of the kinematic branch or inertial branch is observed. To do so, the predictions for the growth rate and celerity of the waves, which differ between the two instability mechanisms, have to be compared with experiments.

## 4.2. *Wave speed and growth rate*

Figure 8a reports the waves speed, at the instability threshold, measured for various volume fractions. The ratio $c_c/u_c$ is almost constant around 3 at low volume fractions, decreases to about 2 between $\phi \approx 0.37$ and $\phi = 0.41 \approx \phi_{DST}$, and decreases further, a little below 2, for higher $\phi$. This behavior agrees well with the prediction of the shear-thickening Kapitza regime expected below $\phi_{DST}$ (black solid line). Above $\phi_{DST}$, the
measurements are found to match better the prediction $\tilde{c}_c = 2$ for the kinematic branch than the prediction $\tilde{c}_c = 2(\alpha - 1) \simeq 0.6$ for the inertial branch, which suggests that the mechanism of the instability observed in experiments is that of Oobleck waves, rather than that of the inertial branch instability. This result is confirmed by a direct comparison between the measured wave speed and the speed of the kinematic waves $c_{kin} \equiv dq/dh_0$, as deduced from the experimental base flow measurements. As shown in the inset of figure 8a, $c_c$ is fairly close to $c_{kin}$ over the whole range of volume fraction studied.

The growth rates of the instability are compared in figure 8b. The measurements are performed when the Reynolds number of the flow is 5% above the observed critical value ($Re/Re_c = 1.05$). The theory is computed for the wavelength observed experimentally at each $\phi$ (see figures 5 and 6). Here again, below $\phi_{DST}$, the growth rate is correctly predicted by the Kapitza instability accounting for continuous thickening (black solid line). Above $\phi_{DST}$, the prediction for the kinematic branch (red solid line) matches the measurements fairly well, whereas that for the inertial branch (blue dashed line) overestimates the observed growth rate by about two orders of magnitude.

The previous results indicate that including inertia to the depth-averaged analysis provides a fair description for both the shear-thickening Kapitza regime observed below $\phi_{DST}$ and the low Reynolds number Oobleck wave regime observed above $\phi_{DST}$. Nonetheless, one important question remains. Above $\phi_{DST}$ the inertial branch has the same instability condition as Oobleck waves, but since the former is expected to amplify two orders of magnitude faster (see figure 8b), why don't we observe, experimentally, the inertial mode rather than Oobleck waves?

## 5. Role of a delay in the rheology

To explain the apparent paradox of the sub-dominance of the inertial branch above $\phi_{DST}$, it is important to realize that the inertial mode in the previous analysis has a singular behavior. It disappears in the strict absence of inertia, but its growth rate diverges as the Reynolds number tends to zero (see 3.24). This singularity at $Re = 0$ results from the assumption of a steady flow rule, which implies that viscosity adapts to change in stress, instantaneously. In reality, when a shear-thickening suspension flow is perturbed, a finite strain $\gamma_0$ is required to relax the fraction of frictional contacts $f$ (hence the viscosity) to the new steady state value (Mari et al. 2015b; Han et al. 2018; Chacko et al. 2018; Richards et al. 2019). As shown in previous studies, this delay in the rheology may strongly modify the stability of the flow in the negatively-sloped region (Nakanishi & Mitarai 2011; Chacko et al. 2018; Richards et al. 2019). For instance, in the case of a simple shear flow driven by a heavy rheometric tool, the addition of a delay stabilizes the flow, by shifting the instability condition on $d\dot{\gamma}/d\tau$ towards negative values (Richards et al. 2019).

To study the influence of a strain delay on our predictions, we extend the linear stability analysis of §3, by adding an evolution equation for the fraction $f$ of frictional contacts. Following Mari et al. (2015b); Han et al. (2018); Chacko et al. (2018); Richards et al. (2019), we assume

$$ \frac{df}{dt} = -\frac{\dot{\gamma}}{\gamma_0}(f - f_{eq}), \qquad (5.1) $$

where $f_{eq} = e^{-\tau^*/\tau}$ is the equilibrium value of $f$, obtained for a steady flow, and $\gamma_0$ is a relaxation strain, whose typical value for frictional spheres is $\sim 10^{-2} - 10^{-1}$ (Mari et al. 2015b; Chacko et al. 2018; Richards et al. 2019). The details of the stability analysis are given in appendix B. We only discuss here the main predictions. The first consequence
of the introduction of a delay is to lead to three modes instead of two. One of them is a modified kinematic mode and the other two are derived from the former inertial branch.

The stability diagram of the modified kinematic mode is shown in figure 9a. In the long wavelength limit ($k \to 0$) and for a small delay ($\gamma_0 \tan \theta \ll 1$), the instability thresholds are set by the following Reynolds number and critical value of A,

$$ Re_{\text{Kap}}^{\gamma_0} \approx \left( 1 + \frac{\gamma_0 \tan \theta}{2A} (2 - A - A^2) \right) Re_{\text{Kap}} \quad \text{and} \quad A_{\perp} \approx -\gamma_0 \tan \theta. \quad (5.2) $$

For $A > 0$, we recover a Kapitza regime above the modified critical Reynolds number $Re_{\text{Kap}}^{\gamma_0}$, which is shifted relative to the zero-delay threshold $Re_{\text{Kap}}$. The delay also shifts slightly the criterion for Oobleck waves, which is now given by $Re < Re_{\text{Kap}}^{\gamma_0}$ and $A < A_{\perp}$ (instead of $A < 0$ without delay). As a result, there exists a narrow range of negative value of A ($A_{\perp} < A < 0$), between the Kapitza and the Oobleck waves instabilities, where the flow is always stable. Therefore, the addition of a small strain delay in the rheology slightly stabilises Oobleck waves but it does not change the properties of the kinematic mode relative to the case without delay, significantly. This is confirmed by verifying that the wave speed in the long wave limit has still the expression for kinematic waves, $\tilde{c}_{\text{kin}} = 2 + A$ (see appendix B), as was the case without delay.

The conclusion is different for the inertial mode, which is now twofold in the presence of a delay. In the long wave ($\tilde{k} \to 0$) and small delay ($\gamma_0 \tan \theta \ll 1$) limit, these two modes are unstable for

$$ A < A_c = -\frac{3\gamma_0}{2Re}. \quad (5.3) $$

This instability threshold is similar to that obtained by Mari et al. (2015b) and Richards et al. (2019) (and previously by Nakanishi & Mitarai (2011) using a different S-shape rheology) when considering the influence of a delay on the stability of a simple shear flow in the DST regime. However, it differs strongly from the simple $A < 0$ condition without delay, as illustrated in the stability diagram of the inertial mode presented in figure 9b. In the negatively-sloped region ($A < 0$), there is now a large domain of stability, which is all the more extended for a low inertia and long delay. Importantly, the addition of a delay also removes the singularity of the inertial mode at $Re = 0$ discussed above. The growth rate now vanishes at the stability threshold instead of diverging (see B18) in appendix B). This means that, slightly above the threshold, the inertial mode grows slowly and not several orders of magnitude faster than the kinematic mode, as predicted in the case without delay. Finally, the addition of a delay also qualitatively changes the frequency of the inertial mode in the long wave limit. Instead of vanishing as $\Re(\tilde{\omega}) \propto \tilde{k}$, the frequency at threshold with delay remains finite at $\tilde{k} = 0$ and is given by $\Re(\tilde{\omega}) \approx \sqrt{6/(\gamma_0 Re)} \approx 2\sqrt{A_c}$. Thus, the long wave limit of the inertial modes with delay is an oscillatory instability, with a wave speed diverging in $\tilde{k} \to 0$. Overall, these results for the inertial modes at $\tilde{k} \to 0$ recover the purely temporal analysis performed by Richards et al. (2019) in the case of a confined Couette flow.

As previously, the onset of instability in the experiments (where the flow rate, hence $\tau_b$, are slowly varied) is expected at the intersection between the critical stability curves and the system trajectories in the $A-Re$ plan for fixed $\phi$ and $\theta$. For $\phi < \phi_{\text{DST}}$, the Kapitza regime is expected, with an onset at $Re = Re_{\text{Kap}}^{\gamma_0}(A)$ (5.2). For $\phi > \phi_{\text{DST}}$, the onset of instability is expected at $A = A_{\perp}(Re)$, for the kinematic branch, and $A = -3\gamma_0/(2Re)$, for the inertial branches. These predictions are plotted in terms of the critical shear stress in figure 9c for four values of the strain delay $\gamma_0 = 0.05, 0.1, 0.2$ and $0.5$, which bound the experimental value determined by Richards et al. (2019). The color intensity of the lines corresponds to increasing values of the delay and the darkest color recalls the stability
<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /></figure>

<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /></figure>

FIGURE 9. Role of a strain delay on stability. (a-b) Stability diagrams (Re, A) for the kinematic branch (a) and the inertial branch (b) ($\theta = 10°$, $\alpha = 1°$). (a) The gray and red domains indicate the unstable regions for $\gamma_0 = 0.1$. (b) The blue dashed line is the instability threshold (5.3) for $\gamma_0 = 0.1$. The coloured lines indicate the system trajectories for increasing flow rate at fixed angle (same as in figure 4). The vertical dark lines indicate the Newtonian case ($A = 1°$). (c-d) Critical shear stress $\tau_c$ (c) and normalized wave velocity $c_c/u_c$ at the threshold (d) versus $\phi$. The symbols report the measurements presented in figures 7(b) and 8(c). The dark and red lines indicate the predictions for the kinematic branch below and above $\phi_{DST}$, respectively. The blue lines is the prediction for the inertial branch. The color intensity of the lines correspond to delay strains of $\gamma_0 = 0.05, 0.1, 0.2$ and $0.5$.

threshold without delay ($\gamma_0 = 0°$). We recover the weak influence of the delay on the onset of the kinematic branch mentioned in the previous discussion, as well as the strong shift in the predicted critical shear stress $\tau_c$ for the inertial branches, which is observed to increase with $\gamma_0$.

Remarkably, for volume fractions sufficiently above $\phi_{DST}$ the inertial mode is stable whatever the basal stress (see the shift of the right-hand boundary of the blue domains in figure 9b), although the slope of the flow rule reaches strongly negative values. This is because the Reynolds number decreases at large basal stress for these high values of $\phi$, such that the trajectory never reaches the stability threshold $A_c = -3\gamma_0/(2Re)$ (for the inclination $\theta = 10°$ considered, which is representative of the experimental range $2°-22°$).
Clearly, this predicted stabilization of the inertial mode at high volume fraction does not match experiments, where waves are all the more unstable as the volume fraction is high. By contrast, no band of unconditionally stable volume fractions is expected for Oobleck waves (see the red curves in figure 9b), which match the observations. Another important difference concerns the wave velocity (figure 9d), which is predicted to diverge for the inertial branches with delay (blue dashed lines), in stark contrast with the prediction for the kinematic mode, which is not far from experimental observation.

These considerations shed some light on the paradoxe of the sub-dominance of the inertial modes above $\phi_{\text{DST}}$. While a strain delay $\gamma_0$ in the rheology only marginally modifies the kinematic mode, it may turn the inertial modes stable even for largely negatively-sloped flow rules, provided the Reynolds number is low enough (which is all the more true as $\phi$ exceeds $\phi_{\text{DST}}$ largely). More precisely, comparing (5.2) with (5.3), indicates that inertial waves are expected before Oobleck waves $(3\gamma_0/(2Re) < \gamma_0 \tan\theta)$, only for $Re/Re_{\text{Kap,Newt}}$ above 9/5, and regardless of the value assumed for $\gamma_0$. Given the smallness of the critical Reynolds numbers actually observed (see figure 7), this suggests that Oobleck waves dominates the inertial modes for most of the volume fractions above $\phi_{\text{DST}}$, which would explain the experimental observations. Another possibility is that the two modes actually coexist but that the inertial modes are not detected in our measurements. Indeed, since Oobleck waves are a zero wavenumber instability with finite wave velocity, while the inertial modes are an oscillatory instability, the two instability are presumably decoupled both spatially and temporally. Therefore, it is possible that the inertial mode be related with the high frequency ‘jittering’ reported by Balmforth et al. (2005), which is not characterized in our experiments. However, such a statement would require more investigations that could be the topic of future studies.

# 6. Conclusion

This study has addressed the stability of a free-surface layer flow of a shear-thickening suspension down an incline. It has shown that the onset of instability and the main characteristics of the waves observed experimentally close to the onset can be rationalized, on the basis of a depth-averaged analysis considering Wyart-Cates effective rheology and inertia, in the different regimes where the flow is continuously shear-thickening, discontinuously shear thickening, or shear jamming.

Below the onset particle volume fraction for discontinuous shear-thickening, $\phi_{\text{DST}}$, the analysis predicts a modified Kapitza (roll wave) instability of inertial origin, which develops above a critical Reynolds number, $Re_{\text{Kap}}$, increasing with increasing $\phi$, as a result of the continuous shear-thickening of the suspension, and which propagates at the velocity of kinematic waves, in good agreement with experimental observations.

Above $\phi_{\text{DST}}$, two other unstable branches are identified, which both stem from the negative slope of the suspension flow rule ($A \equiv d\dot{\gamma}/d\tau < 0$) when the bottom of the flow shear-thickens discontinuously. They have, consequently, the same critical Reynolds number, $Re_{A=0}$, or critical shear-stress, $\tau_{b,A=0}$, set by $A = 0$ (which are lower than the Kapitza threshold values for most of the volume fractions above $\phi_{\text{DST}}$, for cornstarch and provided the plane is far from vertical ($\theta \ll 90^\circ$)). Nonetheless, the mechanisms behind these two branches are fundamentally different. The ‘kinematic branch’ is an amplification of kinematic surface waves due to a mismatch between the free-surface deformation and the basal stress rheology. For a negatively-sloped rheology, this branch is unstable when $Re < Re_{\text{Kap}}$, does not require inertia, and corresponds to the Oobleck waves instability identified in Darbois Texier et al. (2020) from the analysis conducted at $Re = 0$. By contrast the ‘inertial branch’, results from the acceleration of the flowing
layer when the basal stress is velocity-weakening. It requires a negatively-sloped rheology and inertia, but no coupling with the free-surface deformation. These different instability mechanisms yield very different predictions for the growth and propagation of the waves. Those for the Oobleck wave mechanism match the measured wave speed and growth rate much better than those for the inertial branch instability, which supports the mechanism of Oobleck wave formation proposed in Darbois Texier et al. (2020). This conclusion is confirmed by extending the stability analysis to a modified shear-thickening law including a strain delay, which has been shown to have important consequences on the stability of other shear-thickening flows (Nakanishi & Mitarai 2011; Mari et al. 2015b; Han et al. 2018; Chacko et al. 2018; Richards et al. 2019). The addition of a delay is found to modify only slightly the predictions for Oobleck waves, which still agree with measurements, whereas it predicts a strong stabilization of the inertial branch at large volume fraction, where the instability is still observed experimentally.

Overall, this study confirms that the non-inertial Oobleck wave mechanism proposed in Darbois Texier et al. (2020) is at the origin of the wave formation above $\phi_{DST}$, and not inertial modes studied in previous works (Mari et al. 2015b; Richards et al. 2019). To our knowledge, the only other inertialess mechanism reported so far in shear-thickening suspensions is a dynamic vorticity banding instability observed in overdamped discrete numerical simulations (Chacko et al. 2018). However, unlike the Oobleck wave mechanism, this instability (i) develops in the direction transverse to the flow, (ii) requires a velocity-driven configuration (as opposed to the stress-driven configuration of the inclined plane) and (iii) requires that an extra order parameter be added to the flow rule, such that $\dot{\gamma}$ and $\tau$ are not instantaneously related.

Developments remain to be done to obtain a complete description of shear-thickening waves. Extending the (one-dimensional) depth-averaged linear analysis to two dimensions would be needed to capture the correct dissipation at short wavelength. A further extension would be to relax the assumption of a constant volume fraction and account for particle migration induced by the inhomogeneous stress profile in the stability analysis, in order to seek for possible new stabilizing or destabilizing mechanisms (Carpen & Brady 2002; Chacko et al. 2018; Dhas & Roy 2022). Also, waves quickly become highly nonlinear, especially in the discontinuous shear-thickening regime (see figure 2a), and investigating experimentally and theoretically the properties of large amplitude waves, or solitons, would be relevant. The theoretical part of this extension would require considering the effect of higher order terms in the wave dynamics equations. More generally, this study shows the relevance of S-shape constitutive rheological laws for predicting novel hydrodynamic instabilities in shear-thickening suspensions. It also emphasizes the variety of the instability mechanisms resulting from such a rheology, which could be increased by addressing the case of capillary flows, where subtile stabilizing effects can be expected. Finally, beyond shear-thickening suspensions, the theoretical framework adopted in this study could be applied to other complex fluids showing a non-monotonous effective flow curve down slopes, which may be found in geophysical and industrial contexts.

## Appendix A. Experimental data

This appendix compiles the experimental data obtained for the destabilization of the cornstarch suspension flow down the inclined plane. Table 1 lists for each experiment the following quantities: the suspension volume fraction $\phi$, the slope angle $\theta$, the Kapitza's criterion $Re_{Kap}$, the critical flow thickness $h_c$, the critical mean flow velocity $u_c$, the critical Reynolds number $Re_c$, the ratio $Re_c/Re_{Kap}$, the critical basal shear stress $\tau_{b,c}$,
<table><caption>TABLE 1. Experimental data for surface destabilization of a flow of cornstarch suspension down an inclined plane.</caption><thead><tr><th>φ</th><th>θ (°)</th><th>Re<sub>Kap</sub></th><th>h<sub>c</sub> (mm)</th><th>u<sub>c</sub> (cm/s)</th><th>Re<sub>c</sub></th><th>Re<sub>c</sub>/Re<sub>Kap</sub></th><th>τ<sub>b,c</sub> (Pa)</th><th>c<sub>c</sub> (cm/s)</th><th>c<sub>c</sub>/u<sub>c</sub></th><th>c<sub>kin</sub>/u<sub>c</sub></th><th>σ̃</th></tr></thead><tbody><tr><td>0.30</td><td>2</td><td>23.9</td><td>3.91</td><td>12.9</td><td>37.0</td><td>1.55</td><td>1.56</td><td>39.9</td><td>3.10</td><td>2.65</td><td>4.2 × 10<sup>-4</sup></td></tr><tr><td>0.33</td><td>3</td><td>15.9</td><td>2.45</td><td>10.2</td><td>26.4</td><td>1.66</td><td>1.46</td><td>32.9</td><td>3.23</td><td>2.76</td><td>3.8 × 10<sup>-4</sup></td></tr><tr><td>0.34</td><td>9.3</td><td>5.1</td><td>1.78</td><td>11.5</td><td>13.6</td><td>2.67</td><td>3.35</td><td>36.1</td><td>3.14</td><td>2.49</td><td>-</td></tr><tr><td>0.35</td><td>6</td><td>7.9</td><td>2.26</td><td>11.1</td><td>16.0</td><td>2.02</td><td>2.76</td><td>35</td><td>3.16</td><td>2.94</td><td>7.2 × 10<sup>-4</sup></td></tr><tr><td>0.366</td><td>9.3</td><td>5.1</td><td>2.55</td><td>11.6</td><td>9.9</td><td>1.96</td><td>4.85</td><td>36.8</td><td>3.24</td><td>3.03</td><td>-</td></tr><tr><td>0.37</td><td>6</td><td>7.9</td><td>3.40</td><td>14.3</td><td>17.8</td><td>2.24</td><td>4.14</td><td>43.3</td><td>3.03</td><td>2.49</td><td>5.3 × 10<sup>-4</sup></td></tr><tr><td>0.38</td><td>9.3</td><td>5.1</td><td>2.40</td><td>11.2</td><td>11.2</td><td>2.21</td><td>4.60</td><td>33.3</td><td>2.98</td><td>2.51</td><td>-</td></tr><tr><td>0.39</td><td>6</td><td>7.9</td><td>4.75</td><td>18.4</td><td>20.9</td><td>2.63</td><td>5.91</td><td>43.0</td><td>2.35</td><td>3.10</td><td>8.7 × 10<sup>-4</sup></td></tr><tr><td>0.40</td><td>10</td><td>4.7</td><td>5.00</td><td>22.2</td><td>16.5</td><td>3.49</td><td>10.4</td><td>40.2</td><td>1.81</td><td>2.05</td><td>7.7 × 10<sup>-4</sup></td></tr><tr><td>0.41</td><td>10</td><td>4.7</td><td>7.82</td><td>30.7</td><td>21.6</td><td>4.57</td><td>16.3</td><td>-</td><td>-</td><td>-</td><td>5.4 × 10<sup>-4</sup></td></tr><tr><td>0.414</td><td>22</td><td>2.06</td><td>8.60</td><td>37.7</td><td>13.5</td><td>6.56</td><td>38.5</td><td>71.7</td><td>1.90</td><td>2.01</td><td>-</td></tr><tr><td>0.43</td><td>10</td><td>4.7</td><td>11.55</td><td>22.4</td><td>7.66</td><td>1.62</td><td>24.3</td><td>46.4</td><td>2.07</td><td>1.75</td><td>7.7 × 10<sup>-2</sup></td></tr><tr><td>0.44</td><td>10</td><td>4.7</td><td>11.51</td><td>13.6</td><td>2.84</td><td>0.60</td><td>24.3</td><td>20.6</td><td>1.51</td><td>1.52</td><td>9.6 × 10<sup>-2</sup></td></tr><tr><td>0.45</td><td>10</td><td>4.7</td><td>10.79</td><td>8.9</td><td>1.29</td><td>0.27</td><td>22.9</td><td>14.8</td><td>1.67</td><td>1.53</td><td>5.5 × 10<sup>-2</sup></td></tr><tr><td>0.46</td><td>10</td><td>4.7</td><td>8.70</td><td>5.7</td><td>0.66</td><td>0.14</td><td>18.6</td><td>9.0</td><td>1.58</td><td>2.05</td><td>6.7 × 10<sup>-2</sup></td></tr><tr><td>0.47</td><td>10</td><td>4.7</td><td>6.79</td><td>3.9</td><td>0.39</td><td>0.08</td><td>14.6</td><td>5.9</td><td>1.53</td><td>1.58</td><td>9.5 × 10<sup>-2</sup></td></tr><tr><td>0.48</td><td>10</td><td>4.7</td><td>5.42</td><td>2.1</td><td>0.13</td><td>0.03</td><td>11.7</td><td>2.8</td><td>1.37</td><td>1.35</td><td>1.3 × 10<sup>-1</sup></td></tr></tbody></table>

the critical wave speed $c_c$, the ratio $c_c/u_c$, the normalized critical kinematic speed $c_{kin}/u_c$ and the normalized growth rate $\tilde{\sigma}$ for $Re = 1.05 Re_c$.

## Appendix B. Linear stability analysis with a strain delay

This section details the linear stability analysis for Wyart-Cates' rheology with a strain delay $\gamma_0$.

First, the effective basal rheology (3.12) is written in terms of the fraction of frictional contact $f_b$ at the wall ($z = 0$)

$$ \frac{\bar{u}}{h} = \frac{\tau_b}{3\eta(f_b)} \mathcal{G}_2(f_b) \equiv \dot{\Gamma}(\tau_b, f_b), \qquad (\text{B 1}) $$

where $\eta(f_b) = \eta_s (\phi_J(f_b) - \phi)^{-2}$ is the viscosity, $\phi_J(f_b) = \phi_0(1 - f_b) + \phi_1 f_b$ is the critical volume fraction and $\mathcal{G}_2(f_b)$ is a dimensionless corrective factor accounting for the shape of the velocity profile, such that for a steady flow $f_b = f_{eq}(\tau_b) = e^{-\tau^*/\tau_b}$ and $\mathcal{G}_2(f_b) = \mathcal{G}(\tau_b)$. Then, we express the evolution equation (5.1) for the fraction of frictional contact at the bottom of the flowing layer

$$ \frac{\partial f_b}{\partial t} + u_b \frac{\partial f_b}{\partial x} = - \frac{\dot{\gamma}_b}{\gamma_0} [f_b - f_{\text{eq}}(\tau_b)], \qquad (\text{B 2}) $$

where $u_b = 0$ is the velocity at the wall (no-slip condition) and $\dot{\gamma}_b$ is the shear rate at the wall. Writing the flow rule at the bottom $\tau_b = \eta(f_b)\dot{\gamma}_b$ and identifying with (B 1) gives $\dot{\gamma}_b = (3/\mathcal{G}_2(f_b))\bar{u}/h$.

Normalizing and linearizing (B 1) and (B 2) together with the mass and momentum equations (3.10-3.12) around the base state, using $\tilde{x} = x/h_0$, $\tilde{t} = t\bar{u}_0/h_0$, $\tilde{\Gamma} = \dot{\Gamma}h_0/\bar{u}_0$, $\tilde{h} = h/h_0 = 1 + h_1$, $\tilde{u} = \bar{u}/\bar{u}_0 = 1 + u_1$, $\tilde{\tau}_b = \tau_b/\tau_{b,0} = 1 + \tau_1$ and $f_b = f_{\text{eq}}(\tau_{b,0}) + f_1$,
with $|h_1|, |u_1|, |\tau_1|, |f_1| \ll 1$, give

$$ \frac{\partial h_1}{\partial \tilde{t}} + \frac{\partial h_1}{\partial \tilde{x}} + \frac{\partial u_1}{\partial \tilde{x}} = 0, \qquad (B3) $$

$$ \frac{Re}{3} \left[ \frac{\partial u_1}{\partial \tilde{t}} + (\alpha - 1) \frac{\partial h_1}{\partial \tilde{x}} + (2\alpha - 1) \frac{\partial u_1}{\partial \tilde{x}} \right] = h_1 - \tau_1 - \frac{1}{\tan \theta} \frac{\partial h_1}{\partial \tilde{x}}, \qquad (B4) $$

$$ u_1 - h_1 = \tau_1 + \left( \frac{\partial \tilde{\Gamma}}{\partial f_b} \right)_0 f_1, \qquad (B5) $$

$$ \frac{\partial f_1}{\partial \tilde{t}} = -\frac{1}{\gamma_0^*} \left[ f_1 - \left( \frac{d f_{eq}}{d \tilde{\tau}_b} \right)_0 \tau_1 \right], \qquad (B6) $$

where $\gamma_0^* = \gamma_0 \mathcal{G}(\tau_{b,0})/3$. Making use of $\dot{\Gamma} = u/h$, the equations can be recast into the following linear system of $h_1$, $\dot{\Gamma}_1 = u_1 - h_1$ and $f_1$

$$ \frac{\partial h_1}{\partial t} + 2 \frac{\partial h_1}{\partial x} + \frac{\partial \dot{\Gamma}_1}{\partial x} = 0, \qquad (B7) $$

$$ \frac{Re}{3} \left( \frac{\partial \dot{\Gamma}_1}{\partial \tilde{t}} + 2(\alpha - 1) \frac{\partial \dot{\Gamma}_1}{\partial \tilde{x}} + (3\alpha - 4) \frac{\partial h_1}{\partial \tilde{x}} \right) = h_1 - \dot{\Gamma}_1 - \frac{1-A}{B} f_1 - \frac{1}{\tan \theta} \frac{\partial h_1}{\partial \tilde{x}}, \quad (B8) $$

$$ \frac{\partial f_1}{\partial \tilde{t}} = -\frac{1}{\gamma_0^*} \left( A f_1 - B \dot{\Gamma}_1 \right), \qquad (B9) $$

with

$$ A \equiv \left( \frac{d \tilde{\gamma}}{d \tilde{\tau}_b} \right)_0 = 1 + \left( \frac{\partial \tilde{\Gamma}}{\partial f_b} \right)_0 B, \quad \text{and} \quad B \equiv \left( \frac{d f_{eq}}{d \tilde{\tau}} \right)_0. \qquad (B10) $$

The system has non-trivial solutions of the form $h_1 = He^{i(k\tilde{x} - \tilde{\omega}t)}$, $\Gamma_1 = \Gamma e^{i(k\tilde{x} - \tilde{\omega}t)}$ and $f_1 = Fe^{i(k\tilde{x} - \tilde{\omega}t)}$, with $\tilde{k}$ real and $\tilde{\omega}$ complex, only if

$$ \det \left( \begin{array}{ccc} 2i\tilde{k} - i\tilde{\omega} & i\tilde{k} & 0 \\ \frac{Re}{3}(3\alpha - 4)i\tilde{k} - 1 + \frac{i\tilde{k}}{\tan\theta} & \frac{Re}{3}(-i\tilde{\omega} + 2(\alpha - 1)i\tilde{k}) + 1 & \frac{1-A}{B} \\ 0 & -\frac{B}{\gamma_0^*} & -i\tilde{\omega} + \frac{A}{\gamma_0^*} \end{array} \right) = 0, \qquad (B11) $$

i.e.,

$$ \begin{aligned} & -i \frac{Re}{3} \tilde{\omega}^3 + \left( 1 + \frac{Re}{3} \frac{A}{\gamma_0^*} + \frac{2Re}{3} \alpha i \tilde{k} \right) \tilde{\omega}^2 + \left( \frac{i}{\gamma_0^*} - 3 \tilde{k} - \frac{2Re}{3} \frac{A}{\gamma_0^*} \alpha \tilde{k} - \left( \frac{Re}{3} \alpha - \frac{1}{\tan \theta} \right) i \tilde{k}^2 \right) \tilde{\omega} \\ & - \frac{i}{\gamma_0^*} (2 + A) \tilde{k} + \frac{A}{\gamma_0^*} \left( \frac{Re}{3} \alpha - \frac{1}{\tan \theta} \right) \tilde{k}^2 = 0. \end{aligned} \qquad (B12) $$

The relation dispersion (B 12) has three solutions, which are solved numerically to obtain the behavior of the branches for arbitrary values of $\tilde{k}$. The long wave asymptotics can also be obtained analytically. At order $O(\tilde{k}^2)$, the first mode is

$$ \tilde{\omega}_1 = (2 + A)\tilde{k} + i \left[ \frac{A}{\tan \theta} \left( \frac{Re}{Re_{\text{Kap}}} - 1 \right) + \gamma_0^* (A^2 + A - 2) \right] \tilde{k}^2, \qquad (B13) $$

where $Re_{\text{Kap}}$ is the Kapitza threshold without delay depending on $A$, $\alpha$ and $\theta$ (see 3.25). This is the kinematic mode, with a wave speed $\tilde{c} = \Re(\tilde{\omega}_1)/\tilde{k} = 2+A$ that is unchanged by the addition of the delay. However, the stability threshold given by $\Im(\tilde{\omega}_1) = 0$ is modified by the delay.

We consider, first, the case without inertia ($Re = 0$). Equation (B 13) shows that the
kinematic mode is unstable for $A < A_-$ and for $A > A_+$, where

$$A_{\mp} = \frac{1 - \gamma_0^* \tan \theta \mp \sqrt{(1 - \gamma_0^* \tan \theta)^2 + 8\gamma_0^{*2} \tan^2 \theta}}{2\gamma_0^* \tan \theta}, \quad (B 14)$$

are the roots of the polynomial $-A + \gamma_0^* \tan \theta (A^2 + A - 2)$. In practice $\gamma_0^* \tan \theta \ll 1$, such that $A_- \simeq -2\gamma_0^* \tan \theta$ and $A_+ \simeq \frac{1}{2\gamma_0^* \tan \theta}$. Therefore, the delay tends to slightly stabilize the kinematic mode in the case $A < 0$ (S-shape flow rule), whereas it destabilizes the flow for $A \gg 1$ (highly shear-thinning fluid). Since $\gamma_0^* = \gamma_0 \mathcal{G}(\tau_{b,0})/3$ and $\mathcal{G}(\tau_{b,0}) \simeq 3/2$ for $A \simeq 0$ (see eq. 3.20), we have $\gamma_0^* \simeq \gamma_0/2$ close to the $A_-$ threshold. Therefore, for a rheology with a small strain delay and in the absence of inertia, the kinematic mode is unstable for

$$A < A_- \simeq -\gamma_0 \tan \theta, \quad (B 15)$$

where the wave velocity $\tilde{c}_{kin} = 2 + A_- \simeq 2 - \gamma_0 \tan \theta$ is slightly below the value of 2 obtained when there is no delay (see figure 9d).

For a finite inertia ($Re > 0$), the stability threshold of the kinematic mode is set by the modified Kapitza Reynolds number

$$Re_{\text{Kap}}^{\gamma_0} = Re_{\text{Kap}} \left( 1 + \frac{\gamma_0^* \tan \theta}{A} (2 - A - A^2) \right), \quad (B 16)$$

and the growth rate can be rewritten as

$$\Im(\tilde{\omega}_1) = \frac{1}{\tan \theta} \left( \frac{Re}{Re_{\text{Kap}}^{\gamma_0}} - 1 \right) (A + \gamma_0^* \tan \theta (2 - A - A^2)) \tilde{k}^2. \quad (B 17)$$

The stability of the kinematic branch depends on the sign of $(Re/Re_{Kap}^{\gamma_0}) - 1$ and the sign of the polynomial in A. For $A > A_+$ the branch is always unstable. This situation corresponds to a strongly shear-thinning fluid, not considered here. For $0 < A < A_+$, the mode is unstable for $Re > Re_{Kap}^{\gamma_0}$, which corresponds to the inertial Kapitza waves modified by the presence of a delay. For $A < A_-$ (negatively-sloped flow rule), the mode is unstable for low Reynolds numbers $Re < Re_{Kap}^{\gamma_0}$, which corresponds to Oobleck waves. Finally, a small gap exists, $A_- < A < 0$, where the kinematic mode $\tilde{\omega}_1$ is always stable.

The second and third solutions of (B 12) at the lowest order in $\tilde{k}$ are

$$\tilde{\omega}_{2,3} = \frac{-i \left( 1 + \frac{ReA}{3\gamma_0^*} \right) \mp i \sqrt{\Delta}}{2Re/3}, \quad (B 18)$$

with $\Delta = \left(1 + \frac{ReA}{3\gamma_0^*}\right)^2 - \frac{4Re}{3\gamma_0^*}$. They correspond to the inertial branch, which has split in two.

For $\Delta < 0$, the two inertial branches are unstable ($\Im(\tilde{\omega}_{2,3}) > 0$) when $1 + \frac{ReA}{3\gamma_0^*} < 0$. For $A > 0$, this condition is never fulfilled and the two inertial branches are always stable. However, for $A < 0$, the onset of instability is given by

$$A < A_c = -\frac{3\gamma_0^*}{Re} \simeq -\frac{3\gamma_0}{2Re}, \quad (B 19)$$

using as previously $\gamma_0^* \simeq \gamma_0/2$ for small A. In this case, the frequency of the wave at the onset in the long wave limit is $\Re(\tilde{\omega}_{2,3}) = \pm\sqrt{3/(Re\gamma_0^*)} \simeq \pm\sqrt{6/(Re\gamma_0)} \simeq \pm 2\sqrt{A_c}$. The wave velocity $\tilde{c} = \Re(\tilde{\omega}_{2,3})/\tilde{k}$ thus diverges in $\tilde{k} \to 0$.

When $\Delta > 0$, the frequencies $\tilde{\omega}_{2,3}$ are purely imaginary. The imaginary parts are
positive (the flow is unstable) if $1 + \frac{ReA}{3\gamma_0^2} < 0$, which is the same condition as (B 19).
However, in practice the condition $\Delta > 0$ is reached after that given by (5.3) and the
onset of instability is ruled by (5.3) (for $A < 0$).

Declaration of Interests: The authors report no conflict of interest.

Funding: This work was supported by the European Research Council under the Euro-
pean Union Horizon 2020 Research and Innovation program (ERC grant agreement no.
647384) and by ANR ScienceFriction (ANR-18-CE30-0024).

REFERENCES

ABDESSELAM, Y., AGASSANT, J.-F., CASTELLANI, R., VALETTE, R., DEMAY, Y., GOURDIN, D. & PERES, R. 2017 Rheology of plastisol formulations for coating applications. *Polymer Engineering & Science* **57** (9), 982–988.

ALLOUCHE, MOHAMED HATEM, BOTTON, VALÉRY, MILLET, SÉVERINE, HENRY, DANIEL, DAGOIS-BOHY, SIMON, GÜZEL, BÜLENT & HADID, H BEN 2017 Primary instability of a shear-thinning film flow down an incline: experimental study. *Journal of Fluid Mechanics* **821**.

BALMFORTH, NJ & LIU, JJ 2004 Roll waves in mud. *Journal of Fluid Mechanics* **519**, 33–54.

BALMFORTH, N. J., BUSH, J. W. M. & CRASTER, R. V. 2005 Roll waves on flowing cornstarch suspensions. *Physics Letters A* **338** (6), 479–484.

BENJAMIN, T BROOKE 1957 Wave formation in laminar flow down an inclined plane. *Journal of Fluid Mechanics* **2** (6), 554–573.

BLANCO, E., HODGSON, D. J. M., HERMES, M., BESSELING, R., HUNTER, G. L., CHAIKIN, P. M., CATES, M. E., VAN DAMME, I. & POON, W. C. K. 2019 Conching chocolate is a prototypical transition from frictionally jammed solid to flowable suspension with maximal solid content. *Proceedings of the National Academy of Sciences* **116** (21), 10303–10308.

BOERSMA, WH, BAETS, PJM, LAVÈN, JOZUA & STEIN, HNJ 1991 Time-dependent behavior and wall slip in concentrated shear thickening dispersions. *Journal of Rheology* **35** (6), 1093–1120.

BONNOIT, C., DARNIGE, T., CLEMENT, E. & LINDNER, A. 2010 Inclined plane rheometry of a dense granular suspension. *Journal of Rheology* **54** (1), 65–79.

CARPEN, ILEANA C & BRADY, JOHN F 2002 Gravitational instability in suspension flow. *Journal of Fluid Mechanics* **472**, 201–210.

CHACKO, RAHUL N, MARI, ROMAIN, CATES, MICHAEL E & FIELDING, SUZANNE M 2018 Dynamic vorticity banding in discontinuously shear thickening suspensions. *Physical Review Letters* **121** (10), 108003.

CLAVAUD, C., BÉRUT, A., METZGER, B. & FORTERRE, Y. 2017 Revealing the frictional transition in shear-thickening suspensions. *Proceedings of the National Academy of Sciences* **114** (20), 5147–5152.

CLAVAUD, C., METZGER, B. & FORTERRE, Y. 2020 The darcytron: a pressure-imposed device to probe the frictional transition in shear-thickening suspensions. *Journal of Rheology* **64**, 395.

COMTET, J., CHATTÉ, G., NIGUÈS, A., BOCQUET, L., SIRIA, A. & COLIN, A. 2017 Pairwise frictional profile between particles determines discontinuous shear thickening transition in non-colloidal suspensions. *Nature Communications* **8**, 15633.

DARBOIS TEXIER, BAPTISTE, LHUISSIER, HENRI, FORTERRE, YOËL & METZGER, BLOEN 2020 Surface-wave instability without inertia in shear-thickening suspensions. *Communications Physics* **3**.

DHAS, DARISH JESWIN & ROY, ANUBHAB 2022 Stability of gravity-driven particle-laden flows—roles of shear-induced migration and normal stresses. *Journal of Fluid Mechanics* **938**.

DIVOUX, THIBAUT, FARDIN, MARC A, MANNEVILLE, SEBASTIEN & LEROUGE, SANDRA 2016 Shear banding of complex fluids. *Annual Review of Fluid Mechanics* **48**, 81–103.
DONG, JUNHAO & TRULSSON, MARTIN 2017 Analog of discontinuous shear thickening flows under confining pressure. *Physical Review Fluids* 2 (8), 081301.

FORTERRE, YOEL 2006 Kapiza waves as a test for three-dimensional granular flow rheology. *Journal of Fluid Mechanics* 563, 123–132.

FORTERRE, YOËL & POULIQUEN, OLIVIER 2003 Long-surface-wave instability in dense granular flows. *Journal of Fluid Mechanics* 486, 21–50.

FREUNDLICH, H. & RÖDER, H. L. 1938 Dilatancy and its relation to thixotropy. *Transactions of the Faraday Society* 34, 308–316.

GAUTHIER, ANAÏS, PRUVOST, MICKAËL, GAMACHE, OLIVIER & COLIN, ANNIE 2021 A new pressure sensor array for normal stress measurement in complex fluids. *Journal of Rheology* 65 (4), 583–594.

GODDARD, JD 2003 Material instability in complex fluids. *Annual Review of Fluid Mechanics* 35 (1), 113–133.

GUAZZELLI, ÉLISABETH & POULIQUEN, OLIVIER 2018 Rheology of dense granular suspensions. *Journal of Fluid Mechanics* 852.

GUY, B. M., HERMES, M. & POON, W. C. K. 2015 Towards a unified description of the rheology of hard-particle suspensions. *Physical Review Letters* 115 (8), 088304.

GUY, BEN M, NESS, CHRISTOPHER, HERMES, MICHIEL, SAWIAK, LAURA J, SUN, JIN & POON, WILSON CK 2020 Testing the wyart–cates model for non-brownian shear thickening using bidisperse suspensions. *Soft matter* 16 (1), 229–237.

HAN, ENDAO, WYART, MATTHIEU, PETERS, IVO R & JAEGER, HEINRICH M 2018 Shear fronts in shear-thickening suspensions. *Physical Review Fluids* 3 (7), 073301.

HERMES, MICHIEL, GUY, BEN M, POON, WILSON CK, POY, GUILHEM, CATES, MICHAEL E & WYART, MATTHIEU 2016 Unsteady flow and particle migration in dense, non-brownian suspensions. *Journal of Rheology* 60 (5), 905–916.

HSU, CHIAO-PENG, RAMAKRISHNA, SHIVAPRAKASH N, ZANINI, MICHELE, SPENCER, NICHOLAS D & ISA, LUCIO 2018 Roughness-dependent tribology effects on discontinuous shear thickening. *Proceedings of the National Academy of Sciences* 115 (20), 5117–5122.

HWANG, CHI-CHUAN, CHEN, JUN-LIANG, WANG, JAW-SHI & LIN, JENN-SEN 1994 Linear stability of power law liquid film flows down an inclined plane. *Journal of physics D: applied physics* 27 (11), 2297.

JEFFREYS, HAROLD 1925 The flow of water in an inclined channel of rectangular section. *Philosophical Magazine* 49 (293), 793–807.

VON KANN, STEFAN, SNOEIJER, JACCO H, LOHSE, DETLEF & VAN DER MEER, DEVARAJ 2011 Nonmonotonic settling of a sphere in a cornstarch suspension. *Physical Review E* 84 (6), 060401.

KAPITZA, P. L. & KAPITZA, S. P. 1948 Wave flow of thin viscous fluid layers. *Zh. Eksp. Teor. Fiz* 18 (1), 3–28.

LAFARGE 2013 Superplasticizers: the wonder of fluid concrete. https://www.youtube.com/watch?v=CSZxjQwDKF0 38.

LANDAU, LEV DAVIDOVICH & LIFSHITZ, EVGENII MIKHAILEVICH 2013 *Fluid Mechanics: Landau and Lifshitz: Course of Theoretical Physics, Volume 6*, vol. 6. Elsevier.

LIN, N. Y. C., GUY, B. M, HERMES, M., NESS, C., SUN, J., POON, W. C. K. & COHEN, I. 2015 Hydrodynamic and contact contributions to continuous shear thickening in colloidal suspensions. *Physical Review Letters* 115 (22), 228304.

LIU, K. & MEI, C. C. 1994 Roll waves on a layer of a muddy fluid flowing down a gentle slope—a bingham model. *Physics of Fluids* 6 (8), 2577–2590.

LOOTENS, DIDIER, VAN DAMME, HENRI & HÉBRAUD, PASCAL 2003 Giant stress fluctuations at the jamming transition. *Physical review letters* 90 (17), 178301.

MARI, R., SETO, R., MORRIS, J. F & DENN, M. M 2014 Shear thickening, frictionless and frictional rheologies in non-brownian suspensions. *Journal of Rheology* 58 (6), 1693–1724.

MARI, R., SETO, R., MORRIS, J. F. & DENN, M. M. 2015a Discontinuous shear thickening in brownian suspensions by dynamic simulation. *Proceedings of the National Academy of Sciences* 112 (50), 15326–15330.

MARI, ROMAIN, SETO, RYOHEI, MORRIS, JEFFREY F & DENN, MORTON M 2015b Nonmonotonic flow curves of shear thickening suspensions. *Physical Review E* 91 (5), 052302.
MORRIS, JEFFREY F 2018 Lubricated-to-frictional shear thickening scenario in dense suspensions. *Physical Review Fluids* **3** (11), 110508.

NAGAHIRO, SHIN-ICHIRO, NAKANISHI, HIIZU & MITARAI, NAMIKO 2013 Experimental observation of shear thickening oscillation. *Europhysics Letters* **104** (2), 28002.

NAKANISHI, HIIZU & MITARAI, NAMIKO 2011 Shear thickening oscillation in a dilatant fluid. *Journal of the Physical Society of Japan* **80** (3), 033801.

NG, C.-O. & MEI, C. C. 1994 Roll waves on a shallow layer of mud modelled as a power-law fluid. *Journal of Fluid Mechanics* **263**, 151–184.

OLMSTED, PETER D 1999 Two-state shear diagrams for complex fluids in shear flow. *EPL (Europhysics Letters)* **48** (3), 339.

OLMSTED, PETER D 2008 Perspectives on shear banding in complex fluids. *Rheologica Acta* **47** (3), 283–300.

OVARLEZ, GUILLAUME, LE, ANH VU NGUYEN, SMIT, WILBERT J, FALL, ABDOULAYE, MARI, ROMAIN, CHATTÉ, GUILLAUME & COLIN, ANNIE 2020 Density waves in shear-thickening suspensions. *Science advances* **6** (16), eaay5589.

RAMASWAMY, MEERA, GRINIASTY, ITAY, LIARTE, DANILO B, SHETTY, ABHISHEK, KATIFORI, ELENI, DEL GADO, EMANUELA, SETHNA, JAMES P, CHAKRABORTY, BULBUL & COHEN, ITAI 2021 Universal scaling of shear thickening transitions. *arXiv preprint arXiv:2107.13338*.

RATHEE, VIKRAM, BLAIR, DANIEL L & URBACH, JEFFREY S 2017 Localized stress fluctuations drive shear thickening in dense suspensions. *Proceedings of the National Academy of Sciences* **114** (33), 8740–8745.

RICHARDS, J. A., ROYER, J. R., LIEBCHEN, B., GUY, B. M. & POON, W. C. K. 2019 Competing timescales lead to oscillations in shear-thickening suspensions. *Physical Review Letters* **123** (3), 038004.

SAINT-MICHEL, B., GIBAUD, T. & MANNEVILLE, S. 2018 Uncovering instabilities in the spatiotemporal dynamics of a shear-thickening cornstarch suspension. *Physical Review X* **8** (3), 031006.

SEDES, OMER, SINGH, ABHINENDRA & MORRIS, JEFFREY F 2020 Fluctuations at the onset of discontinuous shear thickening in a suspension. *Journal of Rheology* **64** (2), 309–319.

SETO, R., MARI, R., MORRIS, J. F. & DENN, M. M. 2013 Discontinuous shear thickening of frictional hard-sphere suspensions. *Physical Review Letters* **111** (21), 218301.

SINGH, ABHINENDRA, MARI, ROMAIN, DENN, MORTON M & MORRIS, JEFFREY F 2018 A constitutive model for simple shear of dense frictional suspensions. *Journal of Rheology* **62** (2), 457–468.

SPENLEY, NA, YUAN, XF & CATES, ME 1996 Nonmonotonic constitutive laws and the formation of shear-banded flows. *Journal de Physique II* **6** (4), 551–571.

TROWBRIDGE, JH 1987 Instability of concentrated free surface flows. *Journal of Geophysical Research: Oceans* **92** (C9), 9523–9530.

VON KANN, STEFAN, SNOEIJER, JACCO H & VAN DER MEER, DEVARAJ 2013 Velocity oscillations and stop-go cycles: The trajectory of an object settling in a cornstarch suspension. *Physical Review E* **87** (4), 042301.

WHITHAM, G. B. 2011 *Linear and nonlinear waves*,, vol. 42. John Wiley & Sons.

WYART, M. & CATES, M. E. 2014 Discontinuous shear thickening without inertia in dense non-brownian suspensions. *Physical Review Letters* **112** (9), 098302.

YERUSHALMI, JOSEPH, KATZ, STANLEY & SHINNAR, REUEL 1970 The stability of steady shear flows of some viscoelastic fluids. *Chemical Engineering Science* **25** (12), 1891–1902.

YIH, CHIA-SHUN 1963 Stability of liquid flow down an inclined plane. *The physics of Fluids* **6** (3), 321–334.

ZAREI, MOHAMMAD & AALAIE, JAMAL 2020 Application of shear thickening fluids in material development. *Journal of Materials Research and Technology* **9** (5), 10411–10433.