Dynamic wetting by concentrated granular suspensions

Abstract

Many functional materials, such as paints and inks used in applications like coating and 3D printing, are concentrated granular suspensions. In such systems, the contact line dynamics and the internal structure of the suspension interact through shear rate dependent viscosity and microstructural rearrangements. The local shear rate increases sharply near moving contact lines, leading to the non-Newtonian rheology of dense suspensions in this region. While hydrodynamic solutions can describe dilute suspensions, their applicability near advancing contact lines in dense suspensions remains unclear. This study quantifies the deviation from the Newtonian solution by systematically varying interparticle interactions through the choice of dispersion medium. We use silica particles suspended in two refractive index-matched fluids: (i) aqueous 2,2′-thiodiethanol (weak interactions) and (ii) aqueous sodium thiocyanate solution (strong interactions). These systems exhibit substantially different rheological responses, shear-thickening and yield-stress behaviour, respectively. Using astigmatism particle tracking velocimetry (APTV), we resolve the three-dimensional trajectories of tracer particles within a drop driven over a substrate, in an arrangement enabling tracking of the internal flows over a long travel distance of the drop. We observe distinct flow behaviours depending on the particle interactions and the resulting suspension rheology. The more the particle interactions play a role, i.e., the more pronounced the non-Newtonian effects, the more strongly the measured flow profiles differ from the Newtonian solution of the hydrodynamic equations. In the case of the shear-thickening suspension, a notable deviation from Newtonian behaviour is observed. Conversely, the yield-stress suspension exhibits plug flow over the substrate, with Newtonian-like behaviour restricted to the yielded region near the substrate.

---

1 Introduction

Understanding the dynamic wetting process of complex liquids requires a solid understanding of their rheology, surface tension, and dynamic wetting properties of substrates.^1,2 ‘Dynamic wetting’ refers to the physical process in which the three-phase contact line (where the solid, liquid, and gas phases meet) moves across a substrate due to spontaneous or forced mechanisms. This motion causes the advancing/receding contact angles to deviate from their static values, resulting in advancing/receding dynamic contact angles (θ_adv/rec).^3–5 Granular suspensions are one class of complex fluids of an especially high practical relevance. In this context, the term ‘granular’ pertains to particles that are significantly smaller than the characteristic fluid volume, yet sufficiently large, such that the effects of Brownian motion can be disregarded. Consequently, the Peclet number becomes large, indicating that particle motion is primarily controlled by advection due to external forces rather than diffusion. Such suspensions are widely employed in both industrial and scientific settings, including slot-die coating,⁶ 3D-printing,⁷ cosmetics,⁸ painting,⁹ dip-coating,¹⁰ and drop impact on solid surfaces.¹¹ For instance, in the paint industry, pigments are added to the paint. The rheology must be optimized to prevent pigment sedimentation and enable easy re-dispersion through stirring. Additionally, it must display lower viscosity under shear and higher viscosity at rest, i.e., shear-thinning behaviour, facilitating smooth application and adherence to vertical surfaces without excessive dripping.

Hydrodynamic models describe fluid flow near moving contact lines within their range of validity. When extrapolated beyond this, solutions to these models are physically unrealistic due to the contradiction between contact line motion and the boundary condition, characterized by a singularity in the velocity gradient and/or pressure.^12–14 To address this issue, a small slip length was introduced at the substrate, allowing finite velocity of the liquid on the solid substrate and resolving the velocity gradient singularity by permitting slip.¹⁵ However, the divergence in the pressure persists.¹⁴ A detailed mathematical analysis of this subject can be found in ref. 16. The velocity gradient diverges proportionally to the inverse distance to the contact line.³ Consequently, a high shear rate zone is created. Fluids having non-Newtonian properties, such as granular suspensions or polymer solutions, can exhibit even more complex behaviour near the contact line, affecting the distribution of viscous dissipation.^17–19

The classical hydrodynamic description of moving contact line dynamics involves the relationship between the dynamic contact angle and the velocity of the moving contact line. The Cox–Voinov relationship, eqn (1), is generally used to predict the relationship between the dynamic advancing or receding contact angle θ_adv/rec and capillary number Ca near the contact line.^20–22

The capillary number Ca is defined as the ratio of viscous and capillary forces Ca = ηU/γ, where η, γ, and U are the viscosity, surface tension, and contact line velocity, respectively. θ_adv/rec is the dynamic advancing or receding contact angle, and θ _{0, adv/rec} corresponds to the static advancing or receding contact angle.²³ Additionally, ℓ_M is established as the macroscopic length, encompassing the capillary length or the dimensions of a spreading droplet, whereas ℓ_m pertains to a microscopic length. The numerical constant α assumes a non-universal nature and depends on the intricacies within the microscopic and macroscopic boundary conditions. However, most substrates are not ideally smooth. The hydrodynamic model also summarizes microscopic effects, such as pinning on small length scales, using the logarithm of the microscopic term (αℓ_M/ℓ_m) in eqn (1).^16,24 The additional length scale arising from suspended particles lies between the microscopic ℓ_m and macroscopic ℓ_M length scales. A schematic representation of the Cox–Voinov solution near the advancing contact line is provided in Fig. 1.

Fig. 1 Schematic representation of a droplet's advancing contact line and the lubrication flow near it. The Stokes equation is solved within the lubrication approximation, following the individual works of Cox²⁰ and Voinov.²¹

Beyond the continuum-level hydrodynamic description, the behaviour of dense suspensions is strongly influenced by the microscopic contact forces between individual particles. Granular particles experience various types of friction during suspension deformation, similar to their millimeter-scale counterparts.^25,26 These include sliding, rolling, and twisting motions, each associated with distinct frictional behaviours that have been documented at the single-particle level.^27,28 Surface roughness contributes to resistance upon contact. Additional resistance can result from deformation of particle surfaces during interaction.²⁹ In addition to interparticle friction, normal forces arise through the mechanisms described by DLVO theory, eqn (2): van der Waals (vdW) forces and electrostatic forces.

V_Total = V_vdW + V_Elect (2)

The magnitude of vdW interactions is strongly influenced by the properties of the dispersing medium, which can either diminish or augment their effects. The van der Waals attraction between two identical spheres of radius R, separated by a distance h, is proportional to the Hamaker constant (A_H). A_H can be approximated by considering the differences in dielectric constants (ε) and refractive indices (n) between the particles and the dispersing medium, ignoring their frequency dependence.^30,31 When the refractive indices of the particles and the surrounding medium are identical, A_H becomes very small and so does the vdW interaction.^32–35

Rheological measurements provide critical insight into how suspensions respond to applied shear forces. At low particle concentrations, suspensions exhibit Newtonian-like behaviour, with viscosity increasing as a function of volume fraction (ϕ), consistent with the classical expressions of Einstein^36,37 and Batchelor & Green.^38,39 However, as the particle concentration approaches the dense regime, hydrodynamic interactions and interparticle contact forces become increasingly important.⁴⁰ These interactions give rise to pronounced non-Newtonian behaviours. Comprehensive reviews^41,42 outline the current understanding of these complex flow behaviours.

Plotting the steady-shear viscosity versus particle volume fraction reveals a characteristic divergence: as ϕ approaches the critical volume fraction ϕ_c, the viscosity diverges as η ∼ (1 − ϕ/ϕ_c)^−β. Numerous phenomenological expressions exist to describe this divergence. For example, the classical empirical model by Maron–Pierce assumes a fixed β = 2.⁴³ Microscopically, this divergence corresponds to a transition from lubrication-mediated interactions to a frictional contact-dominated network that arrests particle motion.^44,45 The exact value of ϕ_c depends on particle shape, size distribution, and surface roughness.⁴¹ Near jamming, momentum transfer is dominated by particle contacts rather than by the viscous fluid, making the mechanics of these dense suspensions closely analogous to dry granular media.^45,46

At a fixed particle concentration below ϕ_c, as the shear rate increases, dense suspensions can undergo continuous or discontinuous shear-thickening, marked by a transition from fluid-like to solid-like behaviour. It typically begins with a first regime of Newtonian-like or slightly shear-thinning behaviour. This is followed by a second regime of shear-thickening, where viscosity rises sharply and stabilizes at a substantially elevated value.^47–50

In oscillatory measurements, the storage modulus (G′) represents the in-phase (elastic) component of the material response, corresponding to the energy reversibly stored during cyclic deformation. Conversely, the loss modulus (G″) corresponds to the out-of-phase (viscous) component, quantifying the energy dissipated as heat through internal friction as particles and fluid elements experience relative motion.⁵¹ In hard-sphere suspensions, G′ arises from transient particle cages or force-bearing contacts. However, such elastic behaviour is significant only at high volume fractions and within small strain amplitudes that preserve the particle structure. At larger strains or lower concentrations, cage structures yield, and the rheology is dominated by viscous dissipation rather than elastic storage. Although yield-stress is commonly defined as the threshold below which a material does not flow, this transition is often complex and not captured by a single stress value.^52–54 The nature of these contacts, including their strength and geometry, can affect the transmission of forces between these force-carrying particles. The dynamic yield-stress (σ^d_y) refers to the constant steady-state stress observed at low shear rates; below this, some slow creeping motion might still be possible, but the corresponding viscosities can become extremely large.

In summary, considering suspensions as a combination of particles and carrier liquid rather than a single phase provides us with the complete picture to explain the non-Newtonian phenomena observed.⁵⁵ Other studies examining the correlation between rheology and suspension microstructure suggest that merely considering viscosity as a function of volume fraction is insufficient to fully understand suspension behaviour.^56–58 Particularly at higher particle volume fractions, it's crucial to account for shear-induced structures within the suspension.

Addition of particles to the droplet introduces the new length scale of the particle diameter into the dynamic wetting problem. The confined geometry due to the contact angle limits the accessibility of this region for the suspension particles.⁵⁹ The visco-capillary region denotes the microscale zone adjacent to the advancing contact line, where viscous and capillary stresses govern the dynamics, and gravitational effects are negligible. Within a range of a couple of particle diameters, there is an interplay between the wetting dynamics and the non-Newtonian rheological behaviour of suspensions.⁶⁰ Recent studies demonstrated that for particles entering the visco-capillary region at moderate concentrations (volume fraction ϕ = 0.4), the advancing contact line dynamics of the suspension is altered in response to variations in the viscosity.⁵⁹ Still, the Cox–Voinov equation can be used to predict θ³_adv in terms of capillary number by taking the effective viscosity of the suspension into account.⁶⁰

Most studies of receding suspension contact lines focus on dip-coating, where a substrate is immersed in or withdrawn from a particle-laden pool. If the particles are smaller than the fluid film thickness, they become entrained, with entrainment governed by withdrawal velocity, particle concentration, and the receding contact angle. As velocity increases, the system transitions from no entrainment to a particle monolayer and eventually to thick films.^10,61,62 There is a stagnation point on the liquid-air interface near the receding contact line, separating shear flow into the film from recirculating flow into the drop bulk, thereby favoring particle accumulation at the meniscus.⁶³ When particles cluster together, their collective hydrodynamic drag and interfacial deformation enable them to overcome capillary barriers and penetrate into the film.⁶³

Despite extensive research efforts, the complex wetting dynamics of dense granular suspensions remains poorly understood. The local rheology of these suspensions interacts with the shear stress distribution inside the drop, particularly within the high-shear region near the contact line. In this study, we investigate the dynamic wetting behaviour of two carefully designed suspension systems with contrasting rheological properties, achieved by varying the properties of the dispersion medium. One system exhibits viscous and shear-thickening behaviour, while the other displays yield-stress and shear-thinning characteristics. We systematically compare the resulting flow profiles in moving drops to classical hydrodynamic solutions. While the Cox–Voinov relation is a standard reference in the literature, our analysis of measured flow profiles highlights the deviations from such hydrodynamic models when applied to dense granular suspensions. We probe the coupling between contact line motion and rheology by tracking the three-dimensional trajectories of tracer particles near the advancing contact line. The trajectories are resolved using astigmatism particle tracking velocimetry (APTV). To support our interpretation of these flow fields, we perform a detailed rheological characterization of the samples.

2 Materials, methods and sample characterization

2.1 Sample preparation and characterization

For the APTV measurements, it is necessary to ensure that the refractive index of the dispersion medium is matched to that of the dispersed particles. This optical matching renders the silica particles invisible, allowing the analysis to focus on the fluorescently labelled tracer particles introduced into the system. The suspensions were prepared using 5 µm diameter silica particles (Orbit 100 Sil 5 µm, MZ0930-BULK, MZ-Analysentechnik GmbH, Germany), Fig. 2. Monodisperse polystyrene 4.47 µm diameter particles, red-fluorescent (PS-FluoRed): Ex/Em 530 nm/607 nm (microParticles GmbH, Germany) have been used as tracer particles. The tracer particle concentration was around 0.0001 wt%, a negligible value at which any influence on the suspension's hydrodynamics can be safely excluded. Two distinct methods were employed to achieve refractive index (RI) matching. Either a solution of 2,2′-thiodiethanol, 99% (TDE, Sigma-Aldrich) or a solution of sodium thiocyanate, 99% (NaSCN, Sigma-Aldrich) is used. The specifications of the dispersion media are summarized in Table 1.

Fig. 2 Scanning electron microscopy image of spherical silica particles with a diameter of 5 µm, used in the preparation of the suspensions.

Table 1 Specifications of high refractive index solutions used as dispersion media

Notation of the model fluid	Composition	ρ (g cm^−3)	γ (mN m^−1)	η (mPa s)	RI
2,2′-Thiodiethanol (TDE) solution	69.75 wt% in water	1.14	56	8.3 ± 0.1	1.457
Sodium thiocyanate (NaSCN) solution	11.65 mol L^−1 in water	1.27	79	2.8 ± 0.1	1.457

---

The density of each solution is determined using a 25 mL Gay-Lussac pattern pycnometer made of borosilicate glass (BRAND BLAUBRAND). The surface tensions measured using OCA 35L (DataPhysics Instruments GmbH, Germany) and all the refractive index measurements were conducted at the D-Line wavelength of 589.29 nm using a spectral-refractometer (SCHMIDT + HAENSCH GmbH, Germany). The comprehensive details of the refractive index matching approach can be found in other research papers.^28,64–66 Microscopy cover slides (thickness 170 µm) were cleaned by sonicating in water, acetone, and ethanol, then left in a vacuum oven to dry for 12 hours. There was no additional coating applied to the substrates to keep the (dynamic) contact angles below 90°, the ideal range for APTV measurements near the contact line.

The suspensions were prepared by weighing the fluid into glass bottles, then adding the required particle mass to achieve the target weight fraction. Using hand mixing followed by gentle mixing on a rolling device ensures thorough mixing while minimizing air bubble entrapment. The porous silica particles used in this study exhibit a lower effective density (∼1.9 g cm⁻³) than bulk fused silica (2.2 g cm⁻³) due to internal voids. Because of uncertainties associated with this estimation, particle concentrations are reported in terms of weight fraction rather than volume fraction (see the SI, Section S8.1).

2.2 Measurement procedure and data analysis

2.2.1 DropSlider. The DropSlider setup⁶⁷ ensures locally stable moving contact lines, enabling their visualization within the microscope's field of view, for an extended duration. The droplet was pinned while the substrate underwent controlled motion using a piezoelectric motor (LPS-45, Physik Instrumente, Germany). To fix the droplet, the geometry was selected based on the suspension rheology: a disk geometry for TDE-based and a prism geometry for NaSCN-based suspensions (Fig. 3). The disk geometry has a diameter of 5 mm and a gap of 2 mm between the disk and the substrate. The prism geometry uses a right-angle prism with a face length and width of 25 mm. For all experiments, a droplet volume of 50 µL is used.

Fig. 3 Schematic representation of the DropSlider configuration mounted on an inverted microscope for APTV measurements. Two geometries are used: (a) disk geometry for TDE-based suspensions, and (b) prism geometry for NaSCN-based suspensions.

2.2.2 Astigmatism particle tracking velocimetry (APTV).
2.2.2.1 Setup configuration. APTV stands out in scenarios with limited access to the sample, as it can measure the 3D flow field using a single objective and only 2D images. APTV operates with not just one, but two primary focal planes. For APTV measurements, the DropSlider configuration is positioned on the microscope table,⁶⁷ above the objective lens of an inverted fluorescence microscope (Leica DMI6000 B). A combination of a 20× microscope objective and a 150 mm cylindrical lens was employed for the measurements. The cylindrical lens was positioned between the microscope output and the camera sensor, close to the microscope output. The images were acquired using the FASTCAM Mini AX200 camera at a frame rate of 50 Hz, Fig. 4. The substrate travels at a constant velocity of 200 µm s⁻¹ covering a maximum range of 22 mm in 110 seconds.

Fig. 4 Schematic representation of the APTV setup utilized for three-dimensional tracking of fluorescently labelled tracer particles by integrating a cylindrical lens in front of the camera sensor. The framed region highlights the deformation of particle shapes based on their relative distance to the focal planes.

2.2.2.2 Measurement procedure. A steady state and continuous flow was attained approximately 20 s after initiating the movement of the substrate. This steady state was maintained for at least 60 s. The initial 20 seconds period of transient dynamics of the droplet is excluded from the subsequent analysis. In our experiments, the droplet was allowed to rest for 5 minutes after deposition to minimize disturbances from particle sedimentation before starting the experiment. Consequently, there will be a denser layer of particles close to the substrate. In addition, each experiment was repeated at least three times to account for small, unavoidable variations introduced during droplet deposition. This ensured the reproducibility of the reported results.
2.2.2.3 Analysis procedure. The coordinate system is defined such that the x-axis aligns with the direction of droplet motion relative to the substrate, the y-axis lies parallel to the contact line, and the z-axis is normal to the x–y plane, representing the out-of-plane direction. The deformation of particle images from circular to elliptical shapes is correlated with their axial (z) position relative to the focal planes. This relationship is quantified using the aspect ratio, based on calibration procedures described in ref. 67–69. A particle image shape can be reliably determined in the region between z_min and z_max, shown in Fig. 4. This range is determined by the particle intensity being high enough for a reliable analysis. For all the data provided, the plane at z = 0 corresponds to the substrate. The minor vertical displacement of the substrate was corrected by analyzing the trajectories of particles located on the substrate or confined at the contact line. To project data in a coordinate system that is co-moving with the contact line, the position of the contact line was determined in each image. This adjustment ensures that the resulting x coordinate accurately represents the distance of the particles from the contact line, with x = 0 being the contact line position. Further details of the experimental setup and analysis are published elsewhere.⁶⁷ On average, out of 14 360 particle velocities determined from tracking 58 particles, 0.6% of these velocities exceeded the substrate velocity of 200 µm s⁻¹. These are considered to be outliers due to the misdetection of the particle trajectories and are discarded from further analysis.

2.2.3. Side view imaging of droplet shape. Side-view images of the droplet were recorded using a FASTCAM Mini AX200 camera equipped with a Nativar 12X Zoom lens. A SCHOTT KL2500 LED light source with a diffuser was positioned behind the experimental setup to ensure uniform illumination. The droplet shapes were captured at a frame rate of 50 Hz. Droplet surface profiles were analyzed using edge-detection algorithms implemented in MATLAB.

2.3 Particle–particle interactions

Diverse methodologies are employed to elucidate particle–particle interactions and the resulting flow properties. AFM measurements enable direct quantification of nanoscale normal forces.^70,71 Macroscopically, angle of repose measurements primarily reflect rolling friction,⁷² while rheological measurements capture bulk flow behaviour. The following sections detail these measurements.

2.3.1. Rheology measurements.
2.3.1.1 Protorheology, dripping onto substrate (DoS). Jamming in dense granular suspensions denotes the onset of a rigid, solid-like state, marked by a loss of flowability and a divergence in viscosity or the emergence of a finite yield stress. While this transition is typically associated with exceeding a critical particle volume fraction (ϕ_c), it can also be dynamically induced at lower concentrations under sufficient shear stress.^41,42,73 To evaluate the proximity of each sample system to the jamming threshold, we employed dripping onto substrate (DoS) rheometry across a concentration series. DoS rheometry quantifies capillary driven pinch-off dynamics to extract rheological properties and is broadly applicable to Newtonian and complex fluids.⁷⁴

The results of these measurements, shown in Fig. S1a and b of the SI, indicate that the maximum flowable concentration beyond which a transition to a jammed state occurs is system specific. For silica suspensions in TDE, the threshold is 35 wt%, while in NaSCN, the transition occurs at a lower concentration of 30 wt%. We attribute this shift to enhanced interparticle friction in NaSCN, as discussed in Section 2.3.2, which promotes the formation of stress bearing contact networks at lower solid loadings. These findings indicate that the concentration ranges used in subsequent rheometry and wetting experiments are close to the jamming threshold. Complementary rheological characterization using conventional rheometry was limited to 31 wt% in TDE and 28 wt% in NaSCN, the highest concentrations that could be reliably measured in the rheometer due to challenges at higher loadings. Our dynamic wetting experiments focus on suspensions containing 30 and 33 wt% silica in TDE and 30 wt% in NaSCN solutions. These concentrations represent the maximum number of particles within the droplets that still allow their movement. An overview of all concentrations used throughout the study is provided in Table S1 of the SI.

2.3.1.2 Rheology measurements. Shear rate dependent viscosity and oscillatory strain sweep measurements were conducted to characterize the rheological response of the suspensions. The former is presented here, while the latter is detailed in the SI, Section S2. The suspensions consisted of silica particles dispersed in either a TDE or NaSCN solution, with concentrations of 31 wt% and 28 wt%, respectively. During shear rate dependent viscosity experiments, the samples were subjected to continuous shear by gradually increasing the shear rate from 0.1 s⁻¹ to 100 s⁻¹, while the corresponding stress response was monitored. The results of these measurements are presented in Fig. 5. Additional details regarding the measurement procedure and data reproducibility can be found in Fig. S3 of the SI.

Fig. 5 Shear rate dependent viscosity of dense granular suspensions composed of 5 µm silica particles dispersed in TDE and NaSCN solutions, measured at 31 wt% and 28 wt% respectively.

The dispersion of silica particles in TDE solution exhibits non-Newtonian behaviour, characterized by slight shear-thinning at low to intermediate shear rates followed by a pronounced shear-thickening at higher shear rates (Fig. 5). These rheological transitions arise from stress dependent interparticle interactions.⁴⁹ While two-body interactions are important in dilute suspensions,⁷⁵ multi-body interactions dominate dense suspensions.⁴⁹ Prior to shear, particles are randomly distributed, and their inherent excluded volume interactions hinder flow.⁷⁵ At low shear rates, weak hydrodynamic interactions dominate, allowing the particles to gradually rearrange and adopt more anisotropic configurations aligned with the flow direction. This restructuring reduces flow resistance and results in mild shear-thinning.⁷⁶ As the shear rate increases, particles are driven into closer proximity, narrowing the interstitial gaps.⁷⁵ Once the applied stress surpasses a critical threshold, frictional contacts generate transient load-bearing networks that dominate the suspension's rheology.^47,77,78 Rheological measurements were conducted at concentrations slightly below those used in wetting experiments. Nonetheless, for the TDE suspension, the concentration lies well within the regime where shear-thickening is clearly observed (Fig. 5). Notably, in non-Brownian suspensions, shear-thickening occurs at an approximately fixed critical stress across all investigated volume fractions, as established both experimentally and numerically.^47–50 Consequently, this trend is anticipated to hold even at higher concentrations.

On the other hand, the dispersion of silica particles in NaSCN salt solution exhibits a yield-stress fluid characteristic with a plateau at low frequencies, followed by considerable shear-thinning as it flows, Fig. 5. In granular suspensions, bulk stress arises from particle-level interactions, including adhesion and the frictional network of contacts between particles, which collectively govern the complex macroscopic behaviour.⁷⁹ For yield-stress granular suspensions where particle adhesion is weak (discussed further in Section 2.3.3: AFM measurement), yield-stress is mainly due to frictional and geometrical interactions between particles.^80,81 Particles form a network based on their size, shape, and arrangement, with the friction creating resistance to deformation. Therefore, the yield stress emerges from the necessity to break and rearrange the interlocking particles. Shearing disrupts the network formed by interparticle interactions, allowing the suspension to flow, which results in shear-thinning.

The observed differences in flow behaviour, characterized by viscous response with shear-thickening in the TDE-based suspension and yield-stress behaviour in the NaSCN-based system, are further substantiated by oscillatory measurements presented in the SI, Section S2.2. The TDE suspension consistently exhibits a dominant viscous response (G″ > G′) and a pronounced increase in both moduli beyond intermediate strain amplitudes, indicative of strain-hardening driven by the formation of stress-bearing microstructures. In contrast, the NaSCN suspension displays an early crossover between G′ and G″, followed by a monotonic decline in both, indicative of yielding and structural breakdown. These distinct nonlinear viscoelastic responses reflect a mechanistic difference in the underlying stress-bearing network that cannot be attributed solely to variations in effective packing fraction.

2.3.2 Angle of repose measurements. We examined how changing the medium between particles affects the friction coefficient by measuring the angle of repose, which indicates the rolling friction of the particles, with further details provided in Fig. S5 of the SI.^72,82,83 In the TDE solution, the angle of repose was (8 ± 2)°, while in the NaSCN solution it was much higher at (52 ± 5)°. These results indicate that particle friction, particularly rolling friction, is significantly higher in the NaSCN solution. The electrostatic double layer functions as a lubricating interface, mitigating direct particle contact and thereby reducing friction.⁸⁴ In a high ionic strength solution such as the NaSCN solutions, the electrostatic double layer undergoes compression due to the ion screening,³⁰ resulting in increased direct contacts and consequently higher friction. Refer to the SI, Section S4, for a more detailed discussion.

2.3.3 Atomic force microscopy measurements. Particle–particle interactions in the normal direction were measured using Colloidal Probe Atomic Force Microscopy (CP-AFM), employing a silica colloidal probe from the same batch of silica particles used in the other experiments. Experimental details are given in the SI, Section S4. The two index-matched liquids, TDE solution and NaSCN solution, were used as dispersion media for the measurements. The results of adhesion force measurements are presented in Fig. S6 and S7 of the SI. The adhesion forces measured by CP-AFM experiments are slightly higher in TDE solution (−5.1 ± 4.8 nN) compared to NaSCN solution (−1.5 ± 1.3 nN).

Modeling and simulation studies consistently demonstrate that contact forces, particularly the frictional forces arising from sliding and rolling, alongside normal forces, are pivotal in determining the suspension's response to applied shear stress.^85–87 The magnitude of the sliding friction is directly proportional to the normal force and is determined by the static and dynamic friction coefficients. Rolling resistance typically exerts a weaker influence compared to tangential sliding forces.⁸⁷ However, the dominance of one type of force over another is contingent upon particle characteristics, including size, surface roughness, and electrostatic forces, which collectively determine the relative strength of these contact forces.

Based on the characterization experiments, the rheological analysis reveals a significant difference in interparticle dynamics between the two granular suspension systems. Particles in the NaSCN solution exhibit notably enhanced interactions, including frictional contacts and adhesive forces. AFM measurements show no substantial difference in adhesive forces between the systems, with slightly higher adhesion in the TDE medium. Therefore, the increased resistance to particle motion in the NaSCN suspension is mainly due to higher sliding and rolling friction coefficients, as indicated by the angle of repose measurements. We categorize the suspension of silica particles in TDE as weakly interacting particles. Conversely, the suspension of silica particles in a NaSCN solution is described as strongly interacting particles.

3 Results and discussion

3.1 Transient dynamics of receding contact line

To elucidate the interplay between the particulate phase and contact line dynamics, we analyze two representative cases using side-view imaging: (i) a reference case comprising the pure dispersion medium (TDE solution without particles), and (ii) a dense suspension containing 33 wt% silica particles dispersed in the same medium (see Fig. 6). Receding and advancing surface profiles of the moving droplet are extracted from side view images and superimposed to visualize the temporal evolution of the droplet shape. These measurements distinguish the role of the suspending fluid from that of the particle-laden phase, particularly near the receding contact line, Fig. 7.

Fig. 6 Side-view images of a suspension droplet containing 33 wt% silica particles with a diameter of 5 µm dispersed in TDE solution, translating relative to the substrate at a constant velocity of 200 µm s⁻¹.

Fig. 7 Time evolution of droplet surface profiles near the contact lines. (a) and (b) Pure TDE droplet used as a Newtonian reference near the advancing and receding fronts. (c) and (d) Suspension droplet containing 33 wt% silica particles dispersed in TDE. The color gradient from blue to yellow encodes the temporal progression. While the Newtonian case exhibits monotonic interface evolution, the suspension reveals a significant difference at the receding front, where a non-monotonic evolution of curvature emerges.

Side view experiments with the pure dispersion medium indicate clear dewetting close to the receding side. The contact line recedes without leaving a residual film on the substrate. The fluid alone does not form a film during recession. See also the receding contact angle measurements at very low velocities, provided in Fig. S8 of the SI.

In contrast, when employing a 33 wt% concentrated suspension in TDE solution, the dynamics changes fundamentally. At a substrate velocity of 200 µm s⁻¹, we observe no dewetting. The receding contact angle diminishes steadily and eventually decays to near-zero. There is no classical signature of contact line recession. Instead, we observe a continuous particle-laden layer remaining on the substrate in the wake of the contact line. This layer is substantially thicker than the particle diameter (see Fig. S9 of the SI). As noted in the introduction, particle entrainment into thin films requires collective hydrodynamic assembly to overcome capillary resistance.⁶³ In our case, the suspension concentration lies well above the dilute limit. Moreover, interactions between particles and the substrate, as well as interparticle interactions at high concentration, likely contribute to the formation of a thick deposit. We did not further investigate the influence of individual mechanisms.

For the concentrated suspension, we observe a transient evolution of the droplet interface profile near the receding side (Fig. 7d). The slope of the interface initially decreases, indicating a flattening or reduction in the steepness of the curve. This is followed by a subsequent increase in the slope, indicating a re-steepening or rise in the curvature. The droplet typically travels 1 to 1.5 times its initial footprint diameter before attaining a steady-state shape. The comparison of the advancing and receding surface profiles of the suspension in Fig. 7c and d shows that while the advancing side reaches a stationary state after a brief dynamic phase of approximately 10 seconds, the receding side requires over 30 seconds to stabilize. This difference arises from the droplet depositing a suspension layer as it recedes.

The deposition of a thin suspension layer as the droplet recedes is a well-documented phenomenon in similar systems.⁶¹ This behaviour is consistent with the Landau–Levich–Derjaguin relationship, which predicts the thickness of the liquid film entrained during the withdrawal of a plate from a bath, accounting for the dependence of the suspension bulk viscosity on ϕ.⁶² At the liquid–air interface near the receding contact line, the net interfacial stresses acting on fully wetted and neutrally buoyant particles are similar to the Laplace pressure in a pure liquid, thereby justifying the applicability of such classical capillary–hydrodynamic scaling to dense suspensions.⁸⁸

3.2 Advancing contact line of weakly interacting particles

In what follows, we examine the wetting dynamics of suspensions in TDE solution. We first describe the measured tracer trajectories and the resulting average flow fields, and then compare them with theory. Our analysis focuses on how the local flow structure near the advancing contact line interacts with the non-Newtonian rheology of granular suspensions, rather than on testing the contact angle-capillary number relationship.

We compare the average flow fields obtained from experiments across three representative cases: dispersion medium (TDE solution), moderately concentrated suspension (30 wt% silica in TDE), and dense suspension (33 wt% silica in TDE) with analytical predictions for Newtonian fluids near the advancing contact line. We interpret our experimental data using the Moffatt analytical solution for flow near a sharp corner.¹³ Moffatt determined the solution to the Stokes equation in a wedge geometry and associated boundary conditions.¹³ Both the Moffatt and Cox–Voinov solutions describe flow near a moving contact line, differing mainly in their treatment of the innermost region very close to the line, which typically spans from nanometers to a few micrometers.³ Nevertheless, both descriptions are consistent within the outer, Stokes-flow region. Given that our measurements are performed using 5 µm tracer particles, geometric constraints near the contact line set a lower bound on the accessible region. As the presence of the tracer particles perturbs the local flow, reliable velocity measurements can only be obtained at distances l ≫ D, where D is the tracer particle diameter, practically at least 20 µm from the contact line. Our measurements therefore lie within the outer region, where the Moffatt solution provides an appropriate framework for interpreting the flow fields. Physical processes occurring closer to the contact line, from nanometers to a few micrometers, such as van der Waals forces, viscous bending of the interface, and particle-induced interfacial deformations are not resolved in our measurements and lie beyond the scope of this analysis. The analytical expression of the Moffatt solution, along with the corresponding theoretical flow field used for comparison, is provided in the SI, Section S6.1.

Since the analysis involves comparison with a two-dimensional theoretical model, all flow fields are evaluated within the x–z plane. Details of the mean flow field reconstruction from APTV data are provided in the SI, Section S6.2. The resulting two-dimensional velocity fields are shown in the first column of Fig. 8. Each vector's color corresponds to its absolute velocity value, computed as , where v_x,exp and v_z,exp are the experimentally measured velocity components in the horizontal (x) and vertical (z) directions, and the third component parallel to the contact line v_y,exp is excluded.

Fig. 8 Comparison of experimental mean velocity fields obtained from APTV measurements (first column), flow fields predicted by the Moffatt solution (second column) and the deviation fields (third column), which quantify the discrepancy between experiment and theory at the same advancing contact angle. Shown are results for (a)–(c) the TDE dispersion medium, (d)–(f) the 30 wt% moderately concentrated suspension in TDE, and (g)–(i) the 33 wt% dense suspension in TDE, with corresponding contact angles of 38°, 53°, and 80°. For the experimental flow fields, the averaging is performed in the x–z plane. Empty grid points indicate regions with insufficient tracer statistics for reliable averaging. For the theoretical flow fields, velocity vectors are plotted at the center of each 20 µm × 20 µm cell and are color-coded by their magnitude, calculated as . The solid black line marks the advancing contact angle and defines the upper boundary of the wedge-shaped flow domain.

To quantify the difference between experimental measurements and theoretical predictions, we compute a normalized deviation field. For each grid cell, the local deviation is defined as follows:

v_x,dev = v_x,exp − v_x,theo, v_z,dev = v_z,exp − v_z,theo (3)

where v_x,theo and v_z,theo denote the theoretical velocity components in the horizontal (x) and vertical (z) directions predicted by the Moffatt solution. v_x,dev and v_z,dev represent the corresponding deviations between experiment and theory.

The magnitude of deviation is then expressed as a normalized scalar quantity S_dev,

This formulation captures both magnitude and directional discrepancies. The resulting spatial map of v_dev and S_dev is shown in the third column of Fig. 8.

Fig. 8(a–f) illustrate the experimental flow fields, the theoretical predictions, and normalized deviation maps for the dispersion medium and the moderately concentrated suspension (30 wt%). In both cases, there is strong agreement between the measured velocity fields and the analytical Moffatt solution. The deviation magnitude remains generally below 25%, confirming the applicability of classical hydrodynamic theory to these systems. Minor deviations arise from physical and experimental constraints. These include the finite spatial resolution imposed by the tracer particle size, the mismatch between the three-dimensional droplet curvature and the idealized two-dimensional wedge geometry of the Moffatt solution, and the local breakdown of the continuum assumption when approaching particle-scale dimensions (within approximately 50 µm of the contact line). Hence, three-dimensional trajectory visualizations are provided (see Fig. 9a) to capture the spatial complexity of granular suspension flows. Despite the increase in dynamic contact angle from 38° to 53° observed in the 30 wt% suspension, attributable to the enhanced effective viscosity at moderate particle loading, the measured flow field remains consistent with the hydrodynamic solutions.

Fig. 9 Three-dimensional particle trajectories near the advancing contact line: (a) 30 wt% silica particles in TDE solution, hydrodynamic solutions can predict the flow. (b) 33 wt% silica particles in TDE solution, substantial deviation from hydrodynamic flow patterns. The data is shifted along the flow direction such that the y axis coincides with the advancing contact line. The color bar indicates the particle velocity magnitude.

As the particle concentration increases to 33 wt%, entering the dense regime, we observe a clear alteration in the internal dynamics of the suspension near the contact line, a systematic departure from predictions of classical hydrodynamic models. Fig. 8(g–i) present the average velocity field and corresponding deviation map of the dense suspension. The advancing contact angle increases to 80°. Two distinct types of tracer particle trajectories appear at approximately 250 µm from the contact line: one moving with vanishing velocity toward the contact line, while the other moves toward the bulk. Within approximately 250 µm of the contact line, tracer particles exhibit minimal displacement, with velocities remaining extremely low and approaching zero. This behaviour is indicative of a localized stagnation zone that persists throughout the duration of the experiment, which we interpret as the formation of a mechanically arrested structure. Over time, particles located within the stagnant region near the advancing contact line exhibit a gradual migration toward the edges of the droplet. Ultimately, the particles are either re-entrained into the flow or deposited onto the surface. This behaviour is clearly resolved in the three-dimensional trajectory data shown in Fig. 9b, where particles diverge away from the centerline, defined along the direction of relative motion. Additional evidence is provided by the top-view of the average flow fields presented in Fig. S6.3 of the SI. Beyond this region, the flow resumes and particle trajectories show a circulation pattern. The particles initially move toward the advancing contact line. When they reach the proximity of the stagnation zone, they reorient parallel to the substrate and subsequently return to the bulk flow. The deviation field shown in Fig. 8(i) emphasizes the extent of mismatch between experimental data and continuum model predictions. 3D trajectories for another run of the same experiment are provided in Fig. S12 of the SI, reaffirming its reproducibility.

Another characteristic of Fig. 9b is the slower particle movement on the substrate in comparison to the substrate velocity. Specifically, the particle's highest velocity within the TDE suspension droplet, measured close to the substrate, was 150 µm s⁻¹. This is noticeably slower than the substrate's velocity, which was set to 200 µm s⁻¹. This behaviour is consistent with particles rolling and sliding relative to one another and along the substrate, rather than moving affinely with the imposed substrate motion.

Near the advancing contact line, the local shear rates are sufficiently high to trigger non-Newtonian effects. If we consider a distance d in the range of 5 µm (approximately one particle diameter) to 1 mm from the contact line, and assume a characteristic contact line velocity U = 0.2 mm s⁻¹, the local shear rate γ = U/d is estimated to be in the range ∼0.2 to 40. The altered flow pattern reflects the strong coupling between local shear and the strong shear-thickening response observed in shear rate dependent viscosity measurements, Fig. 5. This interpretation is further supported by oscillatory strain sweep measurements, which reveal pronounced strain-induced stiffening in the TDE-based suspension, as evidenced by a sharp increase in both G′ and G″ beyond approximately 30% strain amplitude (see the SI, Fig. S4a). Furthermore, considering the viscosity-volume fraction relationship,⁴³ the suspension concentration lies close to ϕ_c, so that even modest external stresses can drive the system toward a shear-jammed, solid-like state.⁸⁹

When subjected to shear, granular suspensions respond differently depending on solid fraction. Dilute suspensions compact, whereas dense suspensions dilate. This dilation originates from the steric and mechanical constraints inherent to densely packed suspensions, where frictional contacts and geometric exclusion inhibit affine particle motion, thereby necessitating a local expansion of the network to accommodate shear-induced rearrangements.^90,91 In our system with 33 wt% silica particles in TDE, the solid fraction is above the dilation threshold, placing the suspension in the dense regime (see the SI, Section S1). Within the proximity of the high shear zone near the advancing contact line, as the system attempts to dilate under shear, particles press against one another more strongly. This dilation is associated with an increase in particle pressure and the emergence of normal stress differences.⁹² At the particle scale, shear is expected to develop an anisotropic contact and force network, and in the present dense regime this anisotropy is accompanied by shear-induced dilation.⁴⁰ Similarly, we expect such anisotropic stress states to influence local flow resistance and particle dynamics near the contact line, as observed in our APTV measurements.

Compared to Moffatt's solution, the Cox–Voinov model is inherently more complex, as it employs matched asymptotic expansions to bridge the slip-dominated inner region near the contact line with the outer viscous flow. Previous studies have demonstrated that Cox–Voinov solution can predict the advancing contact line behaviour of granular suspensions for concentrations up to 75% of their critical jamming concentration ϕ_c when the suspension's effective viscosity is taken into account.^59,60,93 In our experiments, however, we find that although the apparent contact angle increases with an increase in particle concentration, the average flow field near the advancing contact line of a dense suspension deviates significantly from the predictions of hydrodynamic models, for the outer region defined in the Cox–Voinov framework.

The Cox–Voinov model was developed under the Stokes approximation, assuming small capillary number and viscous (Newtonian) flow where stress scales linearly with strain rate. Cox himself noted that non-Newtonian effects could also regularize the moving contact line singularity.²⁰ While subsequent work has extended the framework to complex fluids by incorporating shear rate dependent viscosities, such as shear-thinning,⁹⁴ no corresponding theory exists for concentrated granular suspensions. In such suspensions, momentum transfer is governed by frictional particle contacts rather than viscous stresses. These contact networks generate strong non-linearities that dominate the flow near the moving contact line, rendering the classical Cox–Voinov model inadequate for describing dense suspensions.

The Cox–Voinov solution assumes a quasi-steady state, with viscous dissipation as the sole energy dissipation mechanism. In dense granular suspensions, however, energy is dissipated not only via flow through the particle network but also at the interparticle contacts within the force-bearing structure.⁴⁰ This introduces localized flow variations and stress fluctuations, rendering the quasi-steady assumption inadequate for describing such systems.

Furthermore, the Cox–Voinov solution relies on the continuum hypothesis, which treats the liquid as a uniform medium across all scales. In our system with 5 µm particles, the discrete nature of the suspended phase generates localized flow disturbances. At low to moderate particle concentrations, the Cox–Voinov description remains applicable, as the disturbances induced by individual particles are effectively averaged out. Our spatial resolution (20 µm), approximately four times larger than the particle size, effectively smooths over the particle scale heterogeneities. Thus, the suspension can still be considered as an effective continuum. As the suspension approaches the jamming threshold, particle motion becomes increasingly collective, and the associated microstructural length scales grow substantially. In this regime, the continuum approximation breaks down, and our flow measurements exhibit pronounced deviations from Cox–Voinov predictions near the advancing contact line.

3.3 Advancing contact line of strongly interacting particles

In what follows, we first present the droplet shapes and particle trajectories for the NaSCN-based suspension, and then discuss how these observations relate to its yield-stress behaviour. Unlike Newtonian fluids, yield-stress droplets spread spontaneously only when capillary stresses exceed the yield threshold.⁹⁵ Capillary stresses may induce limited flow near the substrate, but the upper region remains solid-like, and spreading eventually halts at a finite contact angle.⁹⁵ In our experiments we consider forced wetting, where substrate motion imposes shear stresses exceeding the yield stress, thereby sustaining localized shear-driven flow near the advancing contact line.

For silica suspensions in NaSCN solution, the droplet retains a solid-like character in the absence of substrate motion. This behaviour is consistent with strong interparticle rolling and sliding friction that maintains a persistent contact network, effectively inhibiting flow. This is consistent with the high yield-stress observed in shear rate dependent viscosity measurements (Fig. 5) and further supported by oscillatory strain sweep measurements (SI, Fig. S4b). As a result, the droplet retains its deposited shape and resists deformation, Fig. 10. The circular disk cannot exert sufficient force to hold the droplet in place during substrate motion. To overcome this limitation, a prism geometry is used to hold the droplet while the substrate moves. A side-view image of the yield-stress droplet in motion is provided in Fig. 10. The shape of the droplet tip following deposition with a pipette reflects a balance between yield-stress and surface tension, where elevated yield-stress suppresses capillary-driven thinning.^96,97

Fig. 10 Side view of a droplet containing 30 wt% silica particles in NaSCN solution pushed over a substrate. The blunted tip of the droplet, along with the entrapped air bubbles, indicates a high yield stress. The red circle marks the distribution of shear rate within the droplet near the advancing contact line.

Once the substrate starts moving, particle layers near the substrate experience a high shear stress. When the local stress exceeds the yield threshold, the suspension in proximity to the substrate undergoes yielding and initiates flow. This is consistent with the shear-thinning observed in our rheological measurements (Fig. 5). The resulting particle velocity distribution indicates the presence of three different regions, as shown in the schematic representation of particle trajectories in Fig. 11b. The suspension near the substrate transitions into a liquid-like state, whereas the rest of the droplet experiences creep deformation, with a transitional region occurring between these two regions. Here, creep deformation refers to the very slow, gradual rearrangement of particles under stresses below the yield threshold⁹⁸ and should not be confused with creeping flow (Stokes flow), which describes inertialess, viscosity-dominated motion in Newtonian fluids at low Reynolds numbers. Creep deformation happens in regions where localized stresses are not high enough to cause a strong flow but slightly exceed the local frictional thresholds between particles. The particle structure rearranges itself, allowing gradual deformation through minute movements of particles relative to each other. Particles situated within the bulk of the droplet, where the prevailing shear stress has yet to exceed the yielding threshold, demonstrate a tendency to gradually migrate towards the transitional zone. Upon entry into this transitional domain, a change in their type of motion takes place, causing them to adopt a faster motion parallel to the substrate. The system reaches a dynamic equilibrium, where particle rearrangement balances the applied shear, resulting in steady flow near the high-shear region adjacent to the substrate. Consequently, when observed from the outside, the motion seems akin to a mobile plug traversing the surface. This coexistence of yielded and non-yielded regions within the droplet is consistent with earlier observations,⁹⁵ which also reported a non-flowing region far from the substrate.

Fig. 11 Three-dimensional trajectories of tracer particles near the advancing contact line of a suspension containing 30 wt% silica particles in NaSCN solution. (a) Experimental measurements, where the color bar indicates particle velocity magnitude. A yield-stress suspension is formed; upon initiation of relative motion, the droplet bulk undergoes creep deformation. Particles migrate toward the contact line and, upon entering the high-shear region close to the substrate, they move parallel to it. (b) Schematic illustration highlighting the distinct flow regimes: creeping flow (blue), transitional domain (green), and yielded region (red). The imposed substrate velocity is 0.2 mm s⁻¹.

The measured trajectories of tracer particles near the advancing contact line of the yield-stress suspension are shown in Fig. 11a. Tracer particles move at velocities consistently below that of the substrate. The maximum tracer velocity measured is a factor of two smaller than the substrate velocity, suggesting persistent rolling and sliding along the substrate. Additionally, the average particle velocity in the NaSCN solution is notably lower than in the TDE solution. Specifically, the maximum particle velocity in the TDE suspension is (150 ± 10) µm s⁻¹, compared to 94 ± 10 µm s⁻¹ in the NaSCN suspension. The reduced velocity in the NaSCN solution highlights the profound influence of solvent properties on the rolling and sliding of particles relative to each other and the substrate.

A quantitative estimate of the locally fluidized layer thickness is obtained from bulk rheology data and compared with the flow-field measurements (SI, Section S7). Fitting the NaSCN-based suspension rheology to a Herschel–Bulkley constitutive relationship yields a characteristic yield-onset shear rate [small gamma, Greek, dot above]₀ ≈ 1.9 s⁻¹; together with the imposed substrate velocity U = 0.2 mm s⁻¹, this implies a characteristic yielded-layer thickness δ ∼ U/[small gamma, Greek, dot above]₀ ≈ 100 µm. Consistently, the height–velocity histogram of the tracer velocity magnitude (SI, Fig. S14) shows a narrow band of elevated velocities confined to the first ≈50 µm above the substrate, beyond which the velocity rapidly decays into a low-value plateau.

The emergence of both partially and fully yielded regions within the droplet aligns with the two step yielding framework.⁵⁴ In this framework, jammed suspensions yield sequentially through an initial rupture of weak interparticle bonds (such as adhesive contacts), followed by a second yielding event associated with larger scale structural rearrangements, such as cage escape or aggregate breakdown. In our system, both the droplet bulk and the region near the substrate undergo the first yielding stage. However, only the high-shear region adjacent to the substrate experiences the second yielding stage, wherein the constraints imposed by densely packed particle cages are overcome and the material transitions into a fully fluidized state.

4 Conclusion

This study elucidates the complex interplay between contact line dynamics and the rheology of dense granular suspensions. By employing APTV to resolve three-dimensional flow fields near the advancing contact line and complementary rheological measurements, we directly probed the shear rate dependent internal flow fields of suspensions comprising 5 µm silica particles dispersed in two distinct refractive index-matched media. This allowed us to systematically explore how particle interactions, due to the different dispersion media, influence the wetting dynamics and the corresponding suspension flow behaviour. Near the advancing contact line, there is a localized high-shear zone where the non-Newtonian nature of these granular suspensions emerges. For instance, at extremely high concentrations of particles close to the jamming concentration, silica suspensions in TDE solutions exhibited a transition from hydrodynamic lubrication-dominated flow to frictional contact-dominated behaviour as particle concentration increased, with strong deviations from the flow pattern predicted by conventional solutions of hydrodynamic models near the contact line. In contrast, silica suspensions in NaSCN solution, characterized by yield-stress and shear-thinning behaviour, showed distinct flow states; the droplet showed a plug flow over the substrate with a Newtonian-like behaviour within the boundary layer adjacent to the substrate. These differences underscore the importance of particle-scale interactions, which were validated through complementary CP-AFM and angle of repose measurements. These results highlight the limitations of traditional hydrodynamic models in capturing the complex dynamics of dense suspensions. The assumptions underlying classical relationships, such as the Cox–Voinov relationship, become invalid under conditions of high particle concentration or significant interparticle friction. Instead, the discrete nature of particles, shear-induced structural transitions, and non-Newtonian rheological responses must be integrated into theoretical frameworks to accurately predict contact line behaviour. This understanding has far-reaching implications for industrial applications, where controlling wetting dynamics and suspension flow behaviour is essential for achieving uniform deposition, structural stability, and reliable performance in processes such as 3D printing, coatings, and cosmetics.

Author contributions

R. Azizmalayeri (lead), P. Rostami (equal), T. Witzmann (equal), C. O. Klein (equal), and G. K. Auernhammer (equal) all contributed to: data curation, formal analysis, investigation, methodology, software, validation, visualization, writing – original draft, and writing – review & editing. G. K. Auernhammer additionally contributed to: conceptualization, funding acquisition, project administration, supervision, and resources.

Conflicts of interest

There are no conflicts of interest to declare.

Data availability

Supplementary information is available. See DOI: https://doi.org/10.1039/d5sm00987a.

Further data supporting the findings of this study are available on Zenodo research repository at https://doi.org/10.5281/zenodo.14587081.

Acknowledgements

This research was funded by the Deutsche Forschungsgemeinschaft (DFG), Germany, under Project No. 265191195–SFB 1194, “Interaction between Transport and Wetting Processes”, particularly through projects A06 and A02.

References

1. G. Lu, X.-D. Wang and Y.-Y. Duan, Adv. Colloid Interface Sci., 2016, 236, 43–62.
2. M. Sudheer, P. Yadav, B. Thomas and U. U. Ghosh, Eur. Phys. J.:Spec. Top., 2023, 232, 769–780.
3. J. H. Snoeijer and B. Andreotti, Annu. Rev. Fluid Mech., 2013, 45, 269–292.
4. M. Roché, L. Talini and E. Verneuil, Langmuir, 2024, 40, 2830–2848.
5. J. W. Drelich, L. Boinovich, E. Chibowski, C. Della Volpe, L. Holysz, A. Marmur and S. Siboni, Surf. Innovations, 2020, 8, 3–27.
6. P. Fuaad, A. Swerin, F. Lundell and M. Toivakka, J. Coat. Technol. Res., 2021, 1–10.
7. A. M'barki, L. Bocquet and A. Stevenson, Sci. Rep., 2017, 7, 6017.
8. K. Timm, C. Myant, H. Nuguid, H. Spikes and M. Grunze, Int. J. Cosmet. Sci., 2012, 34, 458–465.
9. M. C. Karlsson, R. Álvarez-Asencio, R. Bordes, A. Larsson, P. Taylor and B.-M. Steenari, J. Coat. Technol. Res., 2019, 16, 607–614.
10. D.-H. Jeong, M. K. H. Lee, V. Thiévenaz, M. Z. Bazant and A. Sauret, J. Fluid Mech., 2022, 936, A36.
11. L. Jørgensen, Y. Forterre and H. Lhuissier, J. Fluid Mech., 2020, 896, R2.
12. C. Huh and L. E. Scriven, J. Colloid Interface Sci., 1971, 35, 85–101.
13. H. K. Moffatt, J. Fluid Mech., 1964, 18, 1–18.
14. C. Huh and S. G. Mason, J. Fluid Mech., 1977, 81, 401–419.
15. E. Dussan, Annu. Rev. Fluid Mech., 1979, 11, 371–400.
16. M. Fricke and D. Bothe, Eur. Phys. J.:Spec. Top., 2020, 229, 1849–1865.
17. X. Li, F. Bodziony, M. Yin, H. Marschall, R. Berger and H.-J. Butt, Nat. Commun., 2023, 14, 4571.
18. H. Grewal, H. Nam Kim, I.-J. Cho and E.-S. Yoon, Sci. Rep., 2015, 5, 14159.
19. P. Rostami, M. Fricke, S. Schubotz, H. Patel, R. Azizmalayeri and G. K. Auernhammer, J. Fluid Mech., 2024, 988, A51.
20. R. Cox, J. Fluid Mech., 1986, 168, 169–194.
21. O. Voinov, Fluid Dyn., 1976, 11, 714–721.
22. L. Hocking, Q. J. Mech. Appl. Math., 1983, 36, 55–69.
23. H.-J. Butt, J. Liu, K. Koynov, B. Straub, C. Hinduja, I. Roismann, R. Berger, X. Li, D. Vollmer and W. Steffen, et al., Curr. Opin. Colloid Interface Sci., 2022, 59, 101574.
24. F. Henrich, D. Fell, D. Truszkowska, M. Weirich, M. Anyfantakis, T.-H. Nguyen, M. Wagner, G. K. Auernhammer and H.-J. Butt, Soft Matter, 2016, 12, 7782–7791.
25. C. Zhao, C. Li and L. Hu, Phys. A, 2018, 492, 181–191.
26. A. P. Santos, D. S. Bolintineanu, G. S. Grest, J. B. Lechman, S. J. Plimpton, I. Srivastava and L. E. Silbert, Phys. Rev. E, 2020, 102, 032903.
27. J. Tomas, Granular Matter, 2004, 6, 75–86.
28. J. Wenzl, R. Seto, M. Roth, H.-J. Butt and G. Auernhammer, Granular Matter, 2013, 15, 391–400.
29. R. Fuchs, T. Weinhart, J. Meyer, H. Zhuang, T. Staedler, X. Jiang and S. Luding, Granular Matter, 2014, 16, 281–297.
30. H.-J. Butt and M. Kappl, Surface and interfacial forces, John Wiley & Sons, 2018.
31. C. Randall and J. Van Tassel, in Encyclopedia of Materials: Science and Technology, ed. K. J. Buschow, R. W. Cahn, M. C. Flemings, B. Ilschner, E. J. Kramer, S. Mahajan and P. Veyssière, Elsevier, Oxford, 2001, pp. 2733–2738.
32. P. S. Yadav, D. Dupre, R. Tadmor, J. S. Park and D. Katoshevski, Surf. Sci., 2007, 601, 4582–4585.
33. M. Wozniak, T. Graule, Y. de Hazan, D. Kata and J. Lis, J. Eur. Ceram. Soc., 2009, 29, 2259–2265.
34. M. Wozniak, Y. de Hazan, T. Graule and D. Kata, J. Eur. Ceram. Soc., 2011, 31, 2221–2229.
35. J. C. Hopkins, D. M. Dryden, W.-Y. Ching, R. H. French, V. A. Parsegian and R. Podgornik, J. Colloid Interface Sci., 2014, 417, 278–284.
36. A. Einstein, Ann. Phys., 1906, 324, 289–306.
37. A. Einstein, Ann. Phys., 1911, 339, 591–592.
38. G. K. Batchelor and J. T. Green, J. Fluid Mech., 1972a, 56, 401–427.
39. G. K. Batchelor and J. T. Green, J. Fluid Mech., 1972b, 56, 375–400.
40. R. Seto and G. G. Giusteri, J. Fluid Mech., 2018, 857, 200–215.
41. É. Guazzelli and O. Pouliquen, J. Fluid Mech., 2018, 852, P1.
42. J. F. Morris, Phys. Rev. Fluids, 2020, 5, 110519.
43. S. H. Maron and P. E. Pierce, J. Colloid Sci., 1956, 11, 80–95.
44. R. Seto, A. Singh, B. Chakraborty, M. M. Denn and J. F. Morris, Granular Matter, 2019, 21, 82.
45. A. L. Liu and S. R. Nagel, Nature, 1998, 396, 21–22.
46. V. Trappe, V. Prasad, L. Cipelletti, P. Segre and D. A. Weitz, Nature, 2001, 411, 772–775.
47. R. Seto, R. Mari, J. F. Morris and M. M. Denn, Phys. Rev. Lett., 2013, 111, 218301.
48. B. M. Guy, M. Hermes and W. C. K. Poon, Phys. Rev. Lett., 2015, 115, 088304.
49. J. F. Morris, Annu. Rev. Fluid Mech., 2020, 52, 121–144.
50. N. Malbranche, A. Santra, B. Chakraborty and J. F. Morris, Front. Phys., 2022, 10, 946221.
51. T. Mezger, The rheology handbook: for users of rotational and oscillatory rheometers, European Coatings, 2020.
52. A. Otsuki, Electrophoresis, 2018, 39, 690–701.
53. M.-G. Li, F. Feng, W.-T. Wu and M. Massoudi, Energies, 2020, 13, 6635.
54. A. Ahuja, A. Potanin and Y. M. Joshi, Adv. Colloid Interface Sci., 2020, 282, 102179.
55. J. F. Morris, Rheol. Acta, 2009, 48, 909–923.
56. L. E. Edens, E. G. Alvarado, A. Singh, J. F. Morris, G. K. Schenter, J. Chun and A. E. Clark, Soft Matter, 2021, 17, 7476–7486.
57. A. Singh, G. L. Jackson, M. Van Der Naald, J. J. De Pablo and H. M. Jaeger, Phys. Rev. Fluids, 2022, 7, 054302.
58. Y.-F. Lee, K. Whitcomb and N. J. Wagner, Phys. Fluids, 2023, 35, 073308.
59. M. Zhao, M. Oléron, A. Pelosse, L. Limat, E. Guazzelli and M. Roché, Phys. Rev. Res., 2020, 2, 022031.
60. A. Pelosse, É. Guazzelli and M. Roché, J. Fluid Mech., 2023, 955, A7.
61. A. Gans, E. Dressaire, B. Colnet, G. Saingier, M. Z. Bazant and A. Sauret, Soft Matter, 2019, 15, 252–261.
62. S. Palma and H. Lhuissier, J. Fluid Mech., 2019, 869, R3.
63. C. E. Colosqui, J. F. Morris and H. A. Stone, Phys. Rev. Lett., 2013, 110, 188302.
64. S. Wiederseiner, N. Andreini, G. Epely-Chauvin and C. Ancey, Exp. Fluids, 2011, 50, 1183–1206.
65. S. F. Wright, I. Zadrazil and C. N. Markides, Exp. Fluids, 2017, 58, 1–39.
66. G. K. Auernhammer, S. Fataei, M. A. Haustein, H. P. Patel, R. Schwarze, E. Secrieru and V. Mechtcherine, Mater. Des., 2020, 193, 108673.
67. B. B. Straub, H. Schmidt, P. Rostami, F. Henrich, M. Rossi, C. J. Kähler, H.-J. Butt and G. K. Auernhammer, Soft Matter, 2021, 17, 10090–10100.
68. C. Cierpka, M. Rossi, R. Segura and C. J. Kähler, Meas. Sci. Technol., 2010, 22, 015401.
69. M. Rossi and C. J. Kähler, Exp. Fluids, 2014, 55, 1–13.
70. H. J. Butt, B. Cappella and M. Kappl, Surf. Sci. Rep., 2005, 59, 1–152.
71. J. Comtet, G. Chatté, A. Nigues, L. Bocquet, A. Siria and A. Colin, Nat. Commun., 2017, 8, 15633.
72. H. M. B. Al-Hashemi and O. S. B. Al-Amoudi, Powder Technol., 2018, 330, 397–417.
73. J. F. Morris, Rheol. Acta, 2023, 62, 617–629.
74. J. Dinic, L. N. Jimenez and V. Sharma, Lab Chip, 2017, 17, 460–473.
75. J. R. Melrose and R. C. Ball, J. Rheol., 2004, 48, 937–960.
76. N. J. Wagner and J. F. Brady, Phys. Today, 2009, 62, 27–32.
77. S. Gallier, E. Lemaire, F. Peters and L. Lobry, J. Fluid Mech., 2014, 757, 514–549.
78. M. Trulsson, E. DeGiuli and M. Wyart, Phys. Rev. E, 2017, 95, 012605.
79. C. Ancey, P. Coussot and P. Evesque, J. Rheol., 1999, 43, 1673–1699.
80. A. Fall, F. Bertrand, G. Ovarlez and D. Bonn, Phys. Rev. Lett., 2009, 103, 178301.
81. A. Singh, MRS Commun., 2023, 13, 971–979.
82. M. Y. Matsuo, D. Nishiura and H. Sakaguchi, Granular Matter, 2014, 16, 441–447.
83. T. Doan, B. Indraratna, T. T. Nguyen and C. Rujikiatkamjorn, Int. J. Geomech., 2023, 23, 04022250.
84. C. Clavaud, A. Bérut, B. Metzger and Y. Forterre, Proc. Natl. Acad. Sci. U. S. A., 2017, 114, 5147–5152.
85. J. Tomas, Part. Sci. Technol., 2009, 27, 337–351.
86. R. Seto, R. Botet, G. K. Auernhammer and H. Briesen, Eur. Phys. J. E, 2012, 35, 1–12.
87. S. Luding, Granular Matter, 2008, 10, 235–246.
88. J. Château, É. Guazzelli and H. Lhuissier, J. Fluid Mech., 2018, 852, 178–198.
89. I. R. Peters, S. Majumdar and H. M. Jaeger, Nature, 2016, 532, 214–217.
90. C. Ness, R. Seto and R. Mari, Annu. Rev. Condensed Matter Phys., 2022, 13, 97–117.
91. I. Srivastava, L. E. Silbert, G. S. Grest and J. B. Lechman, J. Fluid Mech., 2021, 907, A18.
92. A. Deboeuf, G. Gauthier, J. Martin, Y. Yurkovetsky and J. F. Morris, Phys. Rev. Lett., 2009, 102(10), 108301.
93. J. Han and C. Kim, Langmuir, 2012, 28, 2680–2689.
94. D. Bonn, J. Eggers, J. Indekeu, J. Meunier and E. Rolley, Rev. Mod. Phys., 2009, 81, 739–805.
95. G. Martouzet, L. Jørgensen, Y. Pelet, A.-L. Biance and C. Barentin, Phys. Rev. Fluids, 2021, 6, 044006.
96. L.-H. Luu and Y. Forterre, J. Fluid Mech., 2009, 632, 301–327.
97. M. Aytouna, J. Paredes, N. Shahidzadeh-Bonn, S. Moulinet, C. Wagner, Y. Amarouchene, J. Eggers and D. Bonn, Phys. Rev. Lett., 2013, 110(3), 034501.
98. M. Houssais, C. Maldarelli and J. F. Morris, Phys. Rev. Fluids, 2021, 6(1), L012301.

---