Axial forces in capillary liquid bridges of polymer solutions

Abstract

Liquid bridges form between particles during wet mixing with binders or by condensation due to ambient humidity. The consequences of capillary bridges can be quite drastic, creating macroscopic cohesion, as seen in sandcastles and in the formation of particulate agglomerates. Bulk effects in cohesive particles arise from forces generated by capillary bridges, so particle-scale measurements are needed to develop predictive models. Most existing studies at the particle scale assume Newtonian liquids. Yet many binders in industry and in the environment can exhibit viscoelastic behavior. In this study, we measure the axial force generated by liquid bridges of viscoelastic polymer solutions between two spherical beads during controlled uniaxial separation. We vary the polymer concentration, separation velocity, and particle size, and track the force as the bridge thins and ruptures. At quasi-static rates, the axial force remains dominated by capillarity and is not significantly affected by polymer rheology. However, increasing the stretching rate increases the peak force through viscous dissipation and promotes the formation of a viscoelastic filament, thereby delaying rupture. The peak axial forces collapse when rescaled by a capillary number and particle size, while the effective rupture distance collapses with a Weissenberg number. These results provide a simple first-order particle-scale force law for polymeric binders.

---

1. Introduction

Liquid bridges between solid particles are ubiquitous, appearing in contexts ranging from soil moisture^1–3 to industrial granules and wet agglomeration processes.^4–7 Many natural and industrial solids are heterogeneous, for instance, due to extracellular polymeric substances (EPS) in moisture trapped between soil particles⁸ or the use of polymeric binders in industrial slurries.^4,9–11 In high-intensity erosion and debris flows, rapid deformations with large localized strain rates at particle scales can trigger microscale polymeric cohesion.^12–15 Similarly, precise modeling of particle-scale interactions with binders is critical for industrial granulation, fluidization, and asphalt pavement design.^16–18 However, a direct measurement of polymeric axial forces between a pair of individual particles is lacking in the literature.

Particle-scale axial forces in Newtonian liquid bridges, shown in Fig. 1(a), are well-established.^19–21 Early investigations focused on theoretical descriptions of capillary cohesion in soil^1,22 and submerged, gravity-free measurements of inter-particle adhesion.¹⁹ Lian et al.²³ later provided numerical solutions for capillary force and rupture distance, which were validated experimentally using silicone oil by Willett et al.²¹ and Pitois et al.²⁰ Pitois et al. and Ennis et al.²⁴ extended these measurements to fully wetting viscous fluids. Force measurements in industrially relevant water-based solutions add challenges involving evaporation, contact angle hysteresis,²⁵ gravitational effects,^26–28 and three-body interactions.²⁹ More recently, more complex expressions have been developed for capillary bridges between spheres, including models for perfectly and partially wetting particles.^30–33 The present work extends these investigations to include the impact of polymeric viscoelasticity. Indeed, although some studies have considered the breakup of viscoelastic liquid bridges and the liquid transfer between separating solid surfaces with free or pinned contact lines, for instance, in geometries relevant to printing applications, direct measurements of the axial force induced by a polymeric bridge between two spherical particles remain elusive.^34–37

Fig. 1 Thinning of a liquid bridge of (a) Water, and (b) 4 M PEO with 1% by mass concentration prepared in water. The liquid bridge is positioned between two beads of radius R = 2 mm and separated with a velocity v = 10 mm s⁻¹. (c) Schematic of the experimental setup.

Inter-particle axial forces are primarily controlled by capillarity. The capillary thinning dynamics of Newtonian and non-Newtonian liquid bridges are well established.^38–48 For Newtonian liquid bridges, the regime in which the breakup occurs (viscous, inertial, and viscous-inertial) is controlled by the Ohnesorge number, so the relevant pre-elastic thinning dynamics depend on both the liquid properties and the bridge size.⁴⁹ Previous investigations have shown that polymers modify capillary flow in uniaxial extension, where the thinning dynamics are demarcated by distinct Newtonian and viscoelastic regimes, introducing new length scales to the flow.^44,50–53 Related studies on viscoelastic particulate suspensions further showed that coupling polymer elasticity with microstructural heterogeneity modifies the transition to the viscoelastic regime.⁵⁴ Techniques such as FISER,^40,41 CaBER,^42,43 droplet pinch-off,^44,45,55 Dripping-onto-Substrate (DoS),^46,56,57 and the Slow Retraction Method (SRM)^47,48 have established that polymeric axial flows exhibit an initial Newtonian regime, followed by a sharp coil-stretch transition to a viscoelastic regime characterized by cylindrical ligaments. We show this for a capillary liquid bridge between spherical particles in Fig. 1(b). Recent studies have reported a smoother transition at semi-dilute entangled concentrations due to lower critical coil-stretch strain rates.^44,47 More generally, recent work has shown that the onset of the elasto-capillary regime is not determined only by a single relaxation time, but also depends on the initial stretching history, the bridge geometry, and finite extensibility of the polymers.^{47,48,58–61} In the Newtonian regime, coiled polymers contribute to viscous dissipation that scales with concentration, generating larger axial forces than the solvent. Consequently, the combined effects of increased force and delayed rupture distinguish polymeric bridges from Newtonian ones, requiring modified bulk descriptions for granular–polymer mixtures.

This study investigates the axial forces in polymeric liquid bridges across a range of polymer concentrations and separation velocities, v, spanning four orders of magnitude. We specifically use a constant separation velocity, contrasting with the exponential or gravity-driven profiles in previous studies,^43,44 due to the relevance of grain kinematics in Discrete Element Method (DEM) models. Experimental methods and polymer rheology characterization are detailed in Section 2. Section 3 examines these forces under quasi-static conditions (v = 0.01 mm s⁻¹), where thinning is strictly controlled by capillarity, as well as under dynamic separations, where viscoelastic effects become prominent. We then present a first-order analytical framework to describe the evolution of the force across both regimes. In Section 4, we develop dimensionless scaling laws: we collapse the peak force—which governs the bridge's adhesive strength—using the Capillary number, confirm its linear scaling with particle size, and rescale the extended rupture distance using the Weissenberg number. Finally, Section 5 provides our concluding remarks.

2. Experimental methods

2.1. Setup

The custom-built experimental setup [Fig. 1(c)] features a Futek LSB 200 force sensor with an accuracy of ±5 μN. Data acquisition is performed at 10–100 Hz. This rate is well below the sensor's 140 Hz resonant frequency, yet sufficient to resolve the viscoelastic dynamics governed by the polymer relaxation time (≈100 ms). To minimize mechanical noise, the lower bead is mounted on the force sensor, which rests on an active vibration-isolation system (Nexus, Thorlabs). Meanwhile, the upper bead is attached to a linear translation stage (Thorlabs NRT150) that provides 0.01 mm positional precision and a velocity resolution of 0.01 mm s⁻¹. The stage's maximum travel (∼30 mm) exceeds the maximum rupture distance observed in this study (see Section 4.3). Both beads (ruby-doped sapphire, Edmund Optics) have a radius R = 2 mm and a surface roughness on the order of a few microns. Additionally, a subset of experiments utilizes particles with radii ranging from R = 1 to 3 mm. A custom environmental chamber encloses the particle pair to maintain a relative humidity >80%, effectively suppressing evaporation over the experimental timescale.⁶² Finally, a high-speed camera (Phantom VEO 710) records the liquid bridge dynamics during separation at 100–1000 FPS. A validation test for the setup with AP100 silicone oil is described in Section S1 in the SI.

2.2. Methodology

The force measurement procedure is as follows: a solution droplet (V = 0.5–1 µL) is deposited onto the lower bead using a micropipette (0.1–2.5 µL, Eppendorf). Since polymer viscosity and setup geometry make precise dispensing challenging, the actual bridge volume is measured using a custom image-processing routine (see SI Section S6). The upper bead is then lowered onto the droplet and oscillated vertically by a few micrometers at low velocity to ensure axial symmetry without disrupting the initial polymer microstructure. The inter-particle gap is subsequently zeroed by monitoring force readings, which become negative upon contact. Following this, data acquisition begins by synchronizing the sensor, camera, and translation stage to track the force F as a function of the gap S. Separation velocity is varied within the range v ∈ [0.01,10] mm s⁻¹. Finally, force data and high-speed videos are post-processed using custom Python and ImageJ routines. Each experiment was repeated at least three times. Continuous force and thinning curves show representative trials, while scalar quantities reported in the figures, including F_peak, S_peak, and S_rup, are averages over n = 3 independent trials. The corresponding standard deviations are generally smaller than the symbol size or line thickness.

2.3. Rheology of the polymers

Polymer solutions are prepared by dissolving polyethylene oxide (PEO, Sigma-Aldrich) powder with molecular weight M_w = 4000 kg mol⁻¹ in deionized (DI) water; these are hereafter referred to as 4 M PEO. Mixtures are homogenized on a roller (Thermo Fisher ScientificTM) at 40 rpm for 24 hours prior to use. The intrinsic viscosity is [η] = 0.072M_w^0.65 = 1.41 m³ kg⁻¹, yielding an overlap concentration c* = 0.77/[η] = 0.055%.⁵² We test concentrations c = 0, 0.1, 0.5, 1 and 1.5%, spanning the semi-dilute to dense regimes (2 ≤ c/c* ≤ 30). To minimize chain degradation, solutions of each concentration are prepared individually rather than by diluting a stock. The shear rheology for the solutions, summarized in Fig. S2(a) in SI, exhibits strong shear-thinning behavior fitted by the Carreau–Yasuda model. Fig. S2(b) details the corresponding viscous and elastic moduli.

The extensional rheology is summarized in Fig. S3. The relaxation time λ_R, determined via droplet pinch-off, scales as c^0.7, consistent with previous observations.^44,52 The measurement protocol is outlined in Section S3, and further details are available in the literature.^44,46,51,52 We use this relaxation time, λ_R, as a reference timescale for analyzing the axial force and rupture dynamics of the liquid bridge. However, recent studies have reported that relaxation times inferred from capillary-thinning experiments can depend on the deformation history, bridge size, finite extensibility, and the molecular-weight distribution sampled during stretching. As a result, such thinning timescales should be interpreted as protocol-dependent measures rather than as unique geometry-independent material constants.^{47,48,58–60,63}

2.4. Interfacial properties

Interfacial properties are summarized in SI Fig. S4. Fig. S4(a) shows that adding 0.1% polymer reduces the surface tension from γ = 72.8 mN m⁻¹ for water to ≈63 mN m⁻¹. Subsequent concentration increases yield negligible changes, as measured by an Attension Tensiometer. While direct measurement at c = 1.5% proved difficult, we assume γ = 60 mN m⁻¹ consistent with high-molecular-weight solutions.⁶⁴ The solid–liquid contact angle θ in Fig. S4(b) is extracted via custom Python image processing as a function of separation distance S. At small S, θ is initially large because the bridge remains convex at very small gaps, and it then decreases as the meniscus evolves toward a concave shape with increasing separation, before stabilizing at larger gaps.⁶⁵ Although θ varies slightly with concentration and volume, it stabilizes beyond initial gaps to θ ∈ [30°, 50°] for all c.

3. Results

3.1. Quasi-static axial forces

In this subsection, we present axial forces measured under quasi-static conditions. For the bridge volumes used here, a truly static configuration is difficult to maintain because the liquid slowly drains under gravity. We therefore use the lowest accessible separation velocity, v = 0.01 mm s⁻¹, as a practical quasi-static limit. Fig. 2(a) shows that, at a given gap, the liquid bridge contours remain very similar across 4 M PEO solutions with c = 0–1.5%. Consistently, the rescaled force profiles in Fig. 2(b) also vary only weakly with concentration. This indicates that the quasi-static response is controlled mainly by capillarity and bridge geometry, while polymer rheology plays only a minor role in this regime. This interpretation is consistent with the interfacial measurements reported in Fig. S4. Under these conditions, the axial force can therefore be estimated from the liquid bridge shape, as for water.³⁰ Although the measured liquid bridge volume varies slightly across concentrations (V = 0.6–0.9 µL), it does not significantly influence the observed curvature. As a result, capillarity dominates, and we expect similar axial forces for c = 0–1.5%. Therefore, under quasi-static conditions, capillary forces can be estimated by quantifying the liquid bridge shape, since the primary difference between solutions is surface tension, as in water.

Fig. 2 (a) Contours of a quasi-static liquid bridge of c = 0–1.5% 4 M PEO in water for S = 0.2 to S = 1.2 mm at separation velocity v = 0.01 mm s⁻¹ for various gaps S between the particles. (b) Quasi-static axial force as a function of the gap S. The dashed lines show the prediction of the reduced capillary model in eqn (3).

Deriving an exact analytical solution for the capillary force between particles requires solving the Young–Laplace equation, which is generally complex.⁶⁶ Therefore, many studies often use numerical solutions or fitted expressions based on numerical integration of the Young–Laplace equation.^20,21,23 More recently, analytical or near-closed-form expressions have been developed for perfectly wetting bridges, finite contact angles, suction-controlled bridges, and unequal contact angles between the two particles.^30–33,67 In the present work, however, we aim to use a reduced description that can be used in both the quasi-static regime and the later dynamic viscoelastic regime. We therefore express the force in terms of the measured minimum neck radius of the liquid meniscus, R_min, which is directly accessible from the bridge images throughout stretching. Such neck-based descriptions have been widely used for extensional flows of inviscid, viscous, and viscoelastic liquids,^40,43,68 and provide a convenient framework between the static force problem and the thinning dynamics.

To describe the capillarity-dominated axial force on the lower bead, we consider the liquid volume V_s below the minimum meniscus radius R_min, as shown in Fig. 1(c):

F = 2πγR_min − πR_min² ΔP − V_sρg (1)

The first term represents the axial contribution of surface tension along the perimeter 2πR_min. The second term accounts for the Laplace pressure difference across the curved interface, while the third term corrects for gravitational distortions.²¹ The pressure difference ΔP is given by the Young–Laplace equation:

where R₁ and R₂ are the principal radii of curvature of the liquid–air interface. We approximate the first principal radius (inside the liquid) as R₁ ≈ R_min. Image analysis indicates that the local curvature near the neck is dominated by the azimuthal component, such that 1/R₁ ≫ 1/R₂. This assumption is also consistent with the contact-angle measurements for these water-based solutions [see Fig. S4(b) and discussion in Section 5 in the SI]. Finally, for the reduced description adopted here, we neglect the gravitational term. These approximations simplify eqn (1) to:

F ≃ πγR_min (3)

We should emphasize that eqn (3) is a reduced approximation and not a general replacement for full Young–Laplace or closed-form capillary-bridge expressions. Its role here is to capture the measured quasi-static force decay using a single geometric quantity, R_min(S), that can be tracked during the thinning. Eqn (3) is plotted against the experimental data as dashed lines in Fig. 2(b) and captures the experimental results reasonably well. Therefore, this approximation captures the accessible quasi-static force curves sufficiently well for the present first-order model.

Experimentally, we observe that the force initially increases to a peak value before decreasing as S increases.^{28,29,65,69,70} This non-monotonic behavior is consistent with the evolution of the liquid bridge from a convex meniscus at small gaps after the deposition of the liquid to a concave meniscus at larger gaps, together with the associated evolution of the contact line and contact angle.⁶⁵ At small gaps, a contact angle of θ ≳ 90° can produce a repulsive contribution associated with the convex meniscus. As the bridge is stretched, the attractive force reaches a maximum and then decreases as the bridge thins further.⁶⁵ Small differences in the initial placement of the droplet between the beads can also slightly modify the initial wetting state and bridge volume, because forming a reproducible micron-scale contact line on two curved particles is experimentally challenging. The wetting history, initial preparation, and instantaneous bridge geometry may therefore contribute to the observed variations in the peak axial force.⁶⁹ The peak axial force is F ≃ 0.18 ± 0.03 mN across all concentrations (c = 0–1.5%). Furthermore, the F(S) profile is similar for all solutions. Minor variations, particularly for c = 1%, result from a smaller bridge volume (V = 0.6 µL versus 0.7–0.9 µL for the other concentrations). A smaller volume leads to faster force decay, as evidenced by the rupture distance scaling S ∝ V^1/3.²³ Consequently, the measured forces are essentially independent of c. For c ≥ 0.1%, weak viscoelastic effects appear only milliseconds before rupture and are negligible.

3.2. Dynamic liquid bridges: axial force in the pre-elastic dynamic regime

When the separation velocity v > 0.01 mm s⁻¹, the viscoelastic regime becomes increasingly relevant. The transition between regimes is triggered by the imposed strain rate, which, when sufficiently large, leads to the onset of viscoelasticity.⁴⁴ In velocity-controlled setups, such as FISER or the present study, we use the global kinematic measure ε_v = S/R, resulting in an imposed strain rate [small epsi, Greek, dot above]_v = v/R. If [small epsi, Greek, dot above]_v is sufficiently large, it triggers the abrupt uncoiling and extension of polymer chains,⁷¹ resulting in a viscoelastic regime characterized by cylindrical ligaments.

Fig. 3(a) shows the time evolution of the rescaled minimum neck radius, R_min(t)/R, for 4 M PEO solutions at v = 1 mm s⁻¹. As expected, the initial thinning regime (t < t_c) is dominated by Newtonian effects, resulting in similar profiles for all concentrations (c = 0–1.5%).^42,44 At t ∼ t_c, the dynamics become viscoelastic [dashed lines in Fig. 3(a)]. Although t_c varies with concentration, we indicate a representative transition time t_c ≃ 1.3 s for visual guidance. For t > t_c, the bridge undergoes exponential thinning:

R_min ∝ R_ce^−t/3λ_e (4)

where R_c is the critical radius at the onset of the viscoelastic regime and λ_e is a characteristic timescale for relaxation of the stretched polymers in the liquid bridge. Unlike interpretations of CaBER or droplet pinch-off in which the exponential thinning rate directly returns a unique longest relaxation time, recent works have shown that the apparent elasto-capillary timescale can depend on the pre-stretch history, bridge size, and finite extensibility of the polymers, and can also vary with concentration.^{47,48,59,60,63} We therefore interpret λ_e here as an effective thinning timescale for the present constant-velocity liquid bridge protocol, rather than as a direct measurement of a unique longest material relaxation time. The observation that λ_e < λ_R is thus consistent with incomplete chain stretching prior to the formation of the cylindrical filament.⁵⁸ In practice, a more rigorous method to define the transition point S_c is to take the intersection of the asymptotic fits of the Newtonian and viscoelastic regimes, as shown in droplet pinch-off configurations.^44,54 Physically, for the range of polymer concentrations investigated, the peak force F_peak occurs at a gap smaller than the physical transition point S_c (S_peak < S_c). This is because F_peak is associated with the initial dynamics of the contact angle as the liquid bridge moves from a convex to a concave shape, occurring before the onset of any elastic dynamics. However, for a more straightforward analysis and to simplify the modeling of the force profiles, S_c is chosen here as a constant. Under this simplified representation, the gap at which the peak force is measured, S_peak, is consistently smaller than S_c. The same S_c is then used in eqn (7) and (8) to capture the subsequent force evolution.

Fig. 3 (a) Time evolution of the rescaled minimum neck radius R_min(t)/R for liquid bridges of 4 M PEO in water (c = 0 − 1.5%) separating at v = 1 mm s⁻¹. The thinning dynamics transition from Newtonian (t < t_c) to viscoelastic (t ≥ t_c). (b) Corresponding rescaled axial force F/πγR versus separation gap S = vt. Dashed lines correspond to the model predictions in the Newtonian (eqn (6)) and viscoelastic (eqn (7)) regimes. The inset compares the fitted critical force F_c from eqn (7) with the prediction from eqn (8).

The corresponding axial forces are shown in Fig. 3(b). In the Newtonian regime, the axial force measured at v = 1 mm s⁻¹ increases with concentration. Furthermore, an increase in concentration sustains the force for a longer duration following the transition to the viscoelastic regime. In summary, both breakup time and rupture distance increase with v and c, distinguishing viscoelastic bridges from their Newtonian counterparts under dynamic conditions.

In the Newtonian regime (t < t_c), polymer chains in solution remain coiled, and the hydrodynamic interactions between these coils primarily drive viscous dissipation.⁷¹ The relative magnitude of this viscous contribution can be estimated using the Ohnesorge number, . Assuming a characteristic length scale L = 1 mm, density ρ = 1000 kg m⁻³, and viscosity η_eff evaluated at [small gamma, Greek, dot above] ∼ v/L for v = 1 mm s⁻¹ gives Oh = 3.7 × 10⁻³ for water, Oh ≃ 0.1 for the 0.1% PEO solution, and Oh ≃ 1.9–45 for c = 0.5–1.5%. Thus, water and the most dilute PEO solution remain mainly capillary-dominated, whereas viscous stresses become comparable to or larger than capillarity for the higher-concentration PEO solutions. We therefore retain both capillary and viscous contributions in the pre-elastic force balance, without assuming that all fluids are in an asymptotically viscous-breakup regime. With the effective viscosity specified below, the viscous contribution to the axial force can be expressed as:^20,72

We use a prefactor of 3π/2 instead of 6π, as the former derives from the lubrication approximation,²⁰ whereas the latter is associated with Stokes drag. The critical parameter here is the increased solution viscosity η due to the presence of polymers. As a first estimate of the viscous contribution, we use an effective shear-rate-dependent viscosity extracted from the Carreau–Yasuda fit to the shear rheology [see SI Fig. S2(a)], evaluated at [small gamma, Greek, dot above]_eff ∼ v/2R_min. For the force curves, this effective viscosity is evaluated using the measured R_min(S); for the peak-force collapse, the representative value is evaluated at S = S_peak, which occurs very early in the thinning dynamics [see Fig. 3(b)] where the corresponding minimum neck radius is R_min ≈ 0.5 mm (or neck diameter d_min ≈ 1 mm, see Fig. 3(a)). For the separation velocity v = 1 mm s⁻¹ shown in Fig. 3, this yields an estimated shear rate of [small gamma, Greek, dot above]_eff ∼ 1 s⁻¹. Similar estimates for [small gamma, Greek, dot above]_eff are used for the model fits presented in Fig. 5. We emphasize that this is an effective estimate for the pre-elastic regime, not a direct equivalence between shear and extensional viscosities. A more complete treatment would require simultaneous measurements or modeling of the extensional viscosity, polymer stretch, and dynamic contact angle within the bridge, which is beyond the scope of the present measurements. Additionally, a correction factor K = 0.5 accounts for the finite volume of the liquid bridge.²⁰

Fig. 4 Temporal evolution of the liquid bridge profile for 1% 4 M PEO solution at (a) v = 0.01 mm s⁻¹, (b) v = 0.1 mm s⁻¹, (c) v = 1 mm s⁻¹, and (d) v = 10 mm s⁻¹.

Fig. 5 Axial forces in a 1% 4 M polyethylene oxide solution prepared in water as a function of velocity in the range v = 0.01 − 10 mm s⁻¹.

The total axial force in the Newtonian regime is the sum of capillary and viscous contributions:

Eqn (6), plotted as dashed lines in Fig. 3(b), accurately describes the experimental measurements in the Newtonian regime, except as S → 0, where the modeled force diverges. As noted in Section 3.1, we attribute this discrepancy to large contact angles that result in particle–particle repulsion. Altogether, we present a first-order approximation that captures the axial force for a given polymer concentration and separation velocity using independently measured material properties and the literature correction factor K = 0.5. In the Newtonian regime, this description does not require any additional fitting to the force data.

3.3. Axial force in the viscoelastic regime

As noted with eqn (4), the viscoelastic regime (t ≥ t_c) exhibits exponential thinning characterized by cylindrical ligaments. While previous studies using FISER,^40,41 CaBER,^73,74 and droplet thinning,³⁹ have estimated axial stresses in the ligaments, our focus is on the total force exerted on the particles, relevant to bulk cohesive granular flows.^75,76 The end geometry of the liquid bridge likely modifies this total force compared to the localized stresses in the cylindrical region.

From Fig. 3(b), we observe that for t > t_c, the axial force decays exponentially, similar to the thinning dynamics. Motivated by the Oldroyd-B elasto-capillary thinning law and by the measured force decay, we therefore use the phenomenological form

where S_c is the critical gap at the onset of the viscoelastic transition, λ_e is the effective thinning timescale measured from the ligament evolution, and F_c is the transition force obtained from the fit. As shown in Fig. 3, eqn (7) provides a good approximation for the axial forces in this regime.

The prefactor F_c is the critical elastic force at this transition (t = t_c). At its onset, the critical minimum bridge radius is R_c (see eqn (4))—the radius at which polymers are sufficiently stretched to dominate the dynamics. To obtain an order-of-magnitude estimate for F_c, we equate the Newtonian force (eqn (6)) to the polymeric force F ≈ πR_c²τ_zz at S = S_c, and assume that the elastic stress balances the capillary pressure (τ_zz,c ≈ γ/R_c). This scaling argument gives:

Previous studies have shown that the critical radius scales with concentration as R_c ∝ c^0.15.^44,45 Substituting this scaling into eqn (8) allows us to compare the model's predictions with the fitted values of Ḟ_c, as shown in the inset of Fig. 3(b). As Rajesh et al.⁴⁴ noted, R_c increases more rapidly in the semi-dilute entangled regime, which explains the slight overestimation by eqn (8). Nevertheless, this approach establishes a framework for describing the viscoelastic axial force using parameters derived entirely from the thinning dynamics.

In summary, combining eqn (6) and (7) yields a first-order model that captures the axial forces across both the Newtonian and viscoelastic regimes in polymeric liquid bridges. The piecewise framework presented here is a simplified model, neglecting early elastic contributions at higher concentrations and approximating the elasto-capillary transition as discrete regimes. In the following sections, we further investigate these dynamic liquid bridges by varying the separation velocity v, thereby controlling the imposed strain rate [small epsi, Greek, dot above].

3.4. Axial force dependence on stretching velocities

To trigger the viscoelastic transition, the local meniscus strain rate [small epsi, Greek, dot above] must exceed the critical unwinding strain rate, [small epsi, Greek, dot above]_c, of the polymers.⁷¹ This critical rate [small epsi, Greek, dot above]_c depends on polymer properties such as concentration, molecular weight, and solvent viscosity.^44,52 Meanwhile, the imposed strain rate, [small epsi, Greek, dot above]_v = v/R, is velocity-controlled and distinct from the local strain rate [small epsi, Greek, dot above]. Because the concentration changes both the viscosity and the relaxation time, it is useful to separate material-based and process-based measures of elasticity. For the present bridge geometry, we can define Deborah numbers as De_γ = λ_R/t_γ, with t_γ = (ρR³/γ)^1/2, or as De_η = λ_R/t_η, with t_η = η_effR/γ. These quantities compare the intrinsic elasticity of different polymer solutions. However, in the present velocity-controlled experiments, the relevant observation time is t_v = R/v, giving De_v = λ_R/t_v = λ_Rv/R, identical to the imposed Weissenberg number used below. We therefore use Wi to describe the rate-dependent transition and rupture dynamics, while De_γ and De_η remain useful parameters when comparing different polymer solutions or geometries. Fig. 4 illustrates the effect of [small epsi, Greek, dot above]_v on the capillary flow dynamics of a 1% 4 M PEO solution, where similar bridge volume is maintained across different velocities. All else being equal, we observe that the duration of the viscoelastic regime increases with v. In this subsection, we quantify the evolution of the axial force with the separation gap across velocities spanning four orders of magnitude (v = 0.01–10 mm s⁻¹), highlighting the maximum force, F_peak, which controls the adhesive strength of the liquid bridge in bulk granular materials.²⁴

Fig. 5 quantifies the evolution of F(S) across separation velocities v ∈ [0.01,10] mm s⁻¹. The measured forces align well with our qualitative observations: an increase in v prolongs the persistence of the force signal. In addition to the increase in viscous dissipation, rate-dependent changes in the contact angle may also contribute to the early-stage force response, since dynamic contact angle hysteresis in liquid bridges becomes more pronounced as the loading rate increases, thereby modifying the capillary-force response.²⁵ Quantifying this phenomenon is highly relevant to the bulk transport and advection of cohesive grains—such as in fluidization or silo flows^77,78—where large inter-particle distances are typical, and the presence of polymers can significantly modify particle dynamics. To model the force across these velocities, we apply eqn (6) in the Newtonian regime and eqn (7) in the viscoelastic regime. Both equations describe F accurately using parameters extracted directly from shear rheology (η) and thinning dynamics (R_min, λ_e).

The inset of Fig. 5 shows the peak force, F_peak, required to separate the particles, which increases with velocity. This increase in F_peak, coupled with the extended duration of the thinning dynamics, is attributed to enhanced viscous dissipation caused by the uncoiling of a larger fraction of polymers at higher imposed strain rates ([small epsi, Greek, dot above]_v). In the following sections, we further explore the relationship between v and F_peak.

4. Discussion

4.1. Strength of a liquid bridge

Macroscopic cohesive granular flows are typically characterized by their bulk yield stress, τ_z, which is linked to the particle-scale bridge strength via the Rumpf scaling:⁷⁹ τ_z ∝ F_peak/R². Although macroscopic techniques such as powder shear cells are commonly used to measure static yield stress, they are often expensive and prone to high variability. Predicting τ_z directly from particle-scale measurements of F_peak provides a robust alternative and establishes a framework for developing structure–property relationships.^75,80,81 Beyond static strength measurements, particle-scale cohesion is also known to control wet aggregate growth and the flow of cohesive granular assemblies, underscoring the need for particle-level force laws.^6,82,83 In this section, we characterize the evolution of F_peak in polymeric liquid bridges across varying polymer concentrations and separation velocities spanning four orders of magnitude. By rescaling this relationship with the dimensionless Capillary number, Ca = ηv/γ, we demonstrate that the scaling framework originally developed for viscous Newtonian liquids²⁴ can be successfully extended to polymer solutions.

The inset of Fig. 6 shows F_peak for 4 M PEO solutions (c = 0–1.5%) across separation velocities v ∈ [0.01,10] mm s⁻¹. As expected, F_peak scales with both velocity and polymer concentration. Because F_peak is measured within the Newtonian regime, a Newtonian fluid-based description provides a good starting point, despite polymeric contributions to the viscosity. Motivated by Ennis et al.,²⁴ we evaluate eqn (6) at the gap S = S_peak, where the peak force, F_peak = F (S = S_peak), occurs. We then rescale this expression using the Capillary number, yielding:

Fig. 6 Bridge strength (F_peak) measured for polymer solutions across various concentrations c = 0 − 1.5% as a function of the capillary number, Ca. Inset: The bridge strength dependence on separation velocities. All data represent the average of n = 3 replicates.

To test this model, we plot the experimentally measured F_peak rescaled against Ca in Fig. 6. The data collapses onto a master curve bounded by eqn (9) for the experimentally observed gap limits S_peak ∈ [0.1,1] mm (see Fig. S7). The data correspond to averages over n = 3 independent replicates and the standard deviations are generally smaller than the symbol size. We observe a slight deviation in the rescaled data from S_peak = 1 mm, which primarily corresponds to the highest separation velocities. This likely reflects the increasing limitations of shear-rheology-based viscosity estimation under these strongly stretched polymer-bridge conditions. Overall, the collapse shows that, despite complex polymer–solvent interactions that modify the bulk viscosity, at first order the liquid bridge strength F_peak is set mainly by capillary forces and viscous dissipation.

4.2. Role of particle-size on peak force

To quantify the geometric dependence of the axial force, we measure the force for a fixed polymer concentration (1% 4 M PEO) across particle radii R = 1.0 to 3.0 mm. With increasing particle size, we also scale the bridge volume such that R ∝ V^1/3, resulting in a volume range V ∈ [0.12,3.4] µL. The evolution of the measured force as a function of the separation gap is shown in Fig. S8 of the SI, where we note distinct Newtonian and viscoelastic regimes across the various particle sizes. We also observe that the rescaled axial force, F/πγR, collapses onto a single curve.

The peak axial force for various R is plotted in Fig. 7 over a range of separation velocities, increasing with velocity according to a power-law behavior. Because the exponent remains consistent across different particle sizes, we can directly rescale the peak force. As shown in the inset of Fig. 7, normalizing the peak force by the particle size (F_peak/πγR) successfully collapses the data across the full range of radii. This recovers behavior similar to that of Newtonian liquids, in which capillary-bridge-driven axial forces scale with particle size.⁷⁵ This collapse should be interpreted together with the geometric constraint used in the experiment. In capillary-bridge descriptions, the dimensionless force can be written as F/(πγR) = f(S/R, V/R³,θ,…). By scaling the bridge volume approximately as V ∝ R³, we keep the dimensionless volume fixed to first order, so a collapse of F_peak/πγR is expected if the wetting state does not change strongly with particle size. The measured steady contact angle for R = 1–3 mm lies in the range 37°–45°, consistent with the steady value for c = 1% shown in Fig. S4. Thus, the leading order expression remains F ∝ πγR_min, while small changes in volumes or contact-line position contribute to the residual scatter. This interpretation is consistent with recent work showing that wettability can influence early viscoelastic bridge dynamics, even when the later elasto-capillary thinning remains more robust.⁵⁷

Fig. 7 Peak axial force, F_peak for a liquid bridge of 1% 4 M PEO solution measured for various particle radius R ∈ [1,3] mm. Inset: Rescaled peak force, F_peak/πγR, which is independent of particle size. All data represent the average of n = 3 replicates.

4.3. Rupture distance for the liquid bridge

For Newtonian fluids, the rupture distance S_rup is the maximum gap between the particles at which the liquid bridge remains stable.^23,30 However, as noted in previous sections, the viscoelastic regime in polymer solutions modifies it. In these solutions, exponential thinning is typically followed by either satellite drop formation⁸⁴ or a blistering instability.^85–87 The magnitude of the axial force during this instability regime, however, is negligible (F/πγR < 0.01) and unlikely to influence bulk granular kinematics. Consequently, we adopt this force limit as the threshold for defining the effective rupture distance for the polymeric liquid bridges investigated here.

The approximate rupture distance for a quasi-static Newtonian liquid bridge is a function of bridge volume, V, and contact angle, θ (with θ expressed in radians), and may be written as:²³

For the bridge volumes used in the present work (V ∈ [0.5,1] µL), this yields S_rup ∈ [1.14,1.44] mm, which agrees well with experimental observations of S_rup ∈ [1.25,1.89] mm for quasi-static liquid bridges of 4 M PEO (c = 0–1.5%) shown in Fig. 2(b). We note a slight increase in S_rup with c. In this section, we define this quasi-static Newtonian rupture distance as S_rup = S_rup,Newt.

Eqn (10) does not hold for dynamic liquid bridges of polymer solutions, where the rupture distance is significantly larger (S_rup,elast ≫ S_rup,Newt). In Fig. 8, we plot the non-dimensionalized rupture distance S* = S_rup,elast/S_rup,Newt, rescaled with respect to the Weissenberg number imposed on the system by the stretching velocity, Wi = λ_R[small epsi, Greek, dot above]_v. Here, any dynamic effects introduced by the polymer are accounted for in the Weissenberg number. Physically, S_rup,Newt represents the distance at which a liquid bridge without any polymers, separated quasi-statically, would rupture, as described using eqn (10). For the collapse in Fig. 8, S_rup,Newt is estimated using a representative bridge volume V = 0.75 µL, corresponding to the approximate average volume across the measurements. This choice likely contributes part of the residual scatter because the individual bridge volumes span V ∈ [0.5,1] µL. We prefer this quasi-static normalization to a dynamic Newtonian rupture distance, since a reliable model for the rate dependence of Newtonian rupture is not yet available. The relaxation time, λ_R, is obtained from droplet pinch-off measurements (see SI S3). Over the parameter range explored here, rescaling S_rup,elast for polymer solutions collapses the data onto the empirical relation:

S* = 1 + 2Wi (11)

Fig. 8 The rescaled rupture distance S* = S_rup,elast/S_rup,Newt with respect to the Weissenberg number, Wi, of the polymers. The rescaled S* collapses onto a curve described by eqn (11). Inset shows the experimentally measured rupture distance, S_rup. All data represent the average of n = 3 replicates.

Eqn (11) fits the evolution of S* (Wi) very well. The numerical prefactor 2 should be interpreted as an empirical value for the present bead geometry, bridge-volume normalization, PEO molecular weight, concentration range, and wetting conditions, rather than as a universal constant. Establishing whether a similar prefactor applies to other polymer molecular weights, bead radii, contact angles, or bridge volumes will require a broader parametric study. We highlight that rescaling by Wi rather than Ca is appropriate here, as the elasto-capillary balance dominates the rupture distance in the viscoelastic regime. At small [small epsi, Greek, dot above]_v (Wi ∈ [10⁻⁴,10⁻²]), we find S* ∼ 1, aligning with the Newtonian critical rupture gap S_rup,Newt. At larger imposed rates (Wi > 10⁻²), the scaling shifts to S* ∼ Wi, confirming that elastic effects control liquid bridge rupture.

5. Conclusions

Liquid bridges between solid particles are ubiquitous in contexts ranging from soils^1–3 to industrial granules.^4,9 In many of these systems, the liquid phase exhibits heterogeneity due to the presence of dissolved polymers.⁸ However, the influence of such heterogeneities on liquid bridge dynamics remains poorly understood, as previous studies have predominantly focused on ideal, fully wetting Newtonian fluids like silicone oil^20,21 or water.^29,30 In the present work, we investigate polymeric liquid bridges suspended between two solid particles. We quantify the evolution of the axial force as a function of the separation gap across a range of polymer concentrations, separation velocities, and particle sizes.

In the quasi-static regime (v = 0.01 mm s⁻¹), we demonstrated—both qualitatively through the shape of the liquid bridge meniscus and quantitatively via force measurements and scaling arguments—that the axial forces are dominated by capillarity. Within the range of polymer concentrations investigated (2 < c/c* < 30), contributions from viscosity and elastic stresses remain negligible. Under dynamic conditions (v ≥ 0.02 mm s⁻¹), however, we observe and quantify the growing relevance of viscoelastic effects with increasing velocity. To capture this behavior, we divided the axial force evolution into two discrete regimes—Newtonian and viscoelastic—to derive a first-order approximation using parameters obtained solely from shear rheology and thinning dynamics. We then described the evolution of the peak axial force, F_peak, across various velocities and concentrations through a Capillary-number-based rescaling, and confirmed that F_peak scales linearly with the particle radius R. Finally, to account for the significantly extended rupture distance in polymer solutions, we introduced a rescaling based on the dimensionless Weissenberg number. This approach allows us to describe the modified rupture distance across a wide range of v and c using a single equation.

This study extends capillary-force measurements to polymeric liquid bridges, beyond earlier work focused on fully wetting liquids^20,21,24 or water.^28,29 We also note a few limitations of the present work. First, our study focused on purely extensional flows between perfectly spherical particles. In real granular systems, however, liquid bridges experience a combination of extension and shear.^88,89 While shear contributions can sometimes be approximated using a geometric correction factor,⁷² the shear-thinning nature of polymer solutions introduces additional complexity that warrants further investigation. Additionally, real granular media consist of irregularly shaped particles, which present another necessary avenue for follow-up work.^90,91

Despite these limitations, the results provide a practical first-order framework for modeling bulk cohesive granular flows, such as in discrete element method (DEM) simulations.^92,93 When the bridge volume and contact angle are known accurately, more detailed closed-form or approximate capillary-force expressions are also available for efficient DEM implementations in the quasi-static limit.^94,95 By dividing the axial force into distinct Newtonian and viscoelastic regimes, the attractive normal force in a DEM implementation may be written as

where the negative sign denotes attraction and S_rup,elast is defined by eqn (11). Such reduced particle-scale laws are particularly useful for connecting bridge-scale physics to macroscopic behaviors observed in cohesive granular systems, including, for instance, wet accretion and cohesive granular flows.^6,82,83

Future work could refine this model by developing a continuous, rather than piecewise, description of the transition from the Newtonian to the viscoelastic regime, with a specific focus on the critical force F_c. While eqn (8) offers a preliminary model for this boundary, a rigorous study of the transition dynamics is needed. Ultimately, by capturing the influence of polymer additives on bridge rupture and cohesive strength, this framework takes a first step towards bridging the gap between micro-scale rheology and the macroscopic mechanics of polymer-based cohesive granular materials.

Conflicts of interest

There are no conflicts to declare.

Data availability

The data supporting this study are available upon request from the authors.

The supplementary information (SI) contains setup validation, shear and extensional rheology characterization, interfacial measurements, bridge-volume estimation, image-processing details, and additional force, particle-size, and rupture-distance analyses. See DOI: https://doi.org/10.1039/d6sm00314a.

Acknowledgements

This work was supported by the National Science Foundation under NSF CBET PMP Grant No. 2533460, by the U.S. Army Research Office under Grant No. W911NF-23-2-0046, by the International Fine Particle Research Institute, and by the Gordon and Betty Moore Foundation, Grant DOI 10.37807/GBMF13831. The authors thank Anuj Acharya for their help in preliminary studies. SR would like to thank Kirk Fields and Dustin Goetz for initial discussions on the experimental setup.

Notes and references

1. W. B. Haines, J. Agric. Sci., 1925, 15, 529–535.
2. A. Carminati, P. Benard, M. A. Ahmed and M. Zarebanadkouki, Plant Soil, 2017, 417, 1–15.
3. P. Benard, J. R. Schepers, M. Crosta, M. Zarebanadkouki and A. Carminati, Water Resour. Res., 2021, 57, e2021WR030052.
4. C. Reynolds, J. Lam, L. Yang and E. Kendrick, Mater. Des., 2022, 222, 111104.
5. V. Calabrese, ACS Macro Lett., 2026, 15, 368–379.
6. G. Saingier, A. Sauret and P. Jop, Phys. Rev. Lett., 2017, 118, 208001.
7. P. Jop, G. Saingier and A. Sauret, EPJ Web Conf., 2021, 249, 09012.
8. J. Malarkey, J. H. Baas, J. A. Hope, R. J. Aspden, D. R. Parsons, J. Peakall, D. M. Paterson, R. J. Schindler, L. Ye and I. D. Lichtman, et al., Nat. Commun., 2015, 6, 6257.
9. W. J. Lee, N. Park, J. I. Park, J. Nam, K. H. Ahn and J. M. Kim, J. Colloid Interface Sci., 2024, 663, 508–517.
10. X. Zhang, M. Du, H. Fang, M. Shi, C. Zhang and F. Wang, Constr. Build. Mater., 2021, 299, 124290.
11. B. Bitsch, J. Dittmann, M. Schmitt, P. Scharfer, W. Schabel and N. Willenbacher, J. Power Sources, 2014, 265, 81–90.
12. J. Kameda and H. Yohei, Sci. Rep., 2021, 11, 1493.
13. R. Kostynick, H. Matinpour, S. Pradeep, S. Haber, A. Sauret, E. Meiburg, T. Dunne, P. Arratia and D. Jerolmack, Proc. Natl. Acad. Sci. U. S. A., 2022, 119, e2209109119.
14. W. J. Orts, A. Roa-Espinosa, R. E. Sojka, G. M. Glenn, S. H. Imam, K. Erlacher and J. S. Pedersen, J. Mater. Civ. Eng., 2007, 19, 58–66.
15. S. M. Zeinali and S. L. Abdelaziz, Using Biopolymers for Beach Erosion Mitigation, Geotechnical Frontiers, 2025, pp. 21–30.
16. B. J. Ennis, G. Tardos and R. Pfeffer, Powder Technol., 1991, 65, 257–272.
17. G. Zhou, Y. Xu, Q. Wang, B. Jiang, B. Ren, X. Zhang and L. Yi, J. Environ. Manage., 2023, 330, 117097.
18. F. Beainy, S. Commuri, M. Zaman and I. Syed, International Journal of Geomechanics, 2013, 13, 581–594.
19. G. Mason and W. Clark, Chem. Eng. Sci., 1965, 20, 859–866.
20. O. Pitois, P. Moucheront and X. Chateau, J. Colloid Interface Sci., 2000, 231, 26–31.
21. C. D. Willett, M. J. Adams, S. A. Johnson and J. P. Seville, Langmuir, 2000, 16, 9396–9405.
22. R. Fisher, J. Agric. Sci., 1926, 16, 492–505.
23. G. Lian, C. Thornton and M. J. Adams, J. Colloid Interface Sci., 1993, 161, 138–147.
24. B. J. Ennis, J. Li and P. Robert, et al., Chem. Eng. Sci., 1990, 45, 3071–3088.
25. Z. Shi, Y. Zhang, M. Liu, D. A. Hanaor and Y. Gan, Colloids Surf., A, 2018, 555, 365–371.
26. D. N. Mazzone, G. I. Tardos and R. Pfeffer, J. Colloid Interface Sci., 1986, 113, 544–556.
27. S. Wang, F. Liu, C. Pu, J. Cui and Z. Zeng, Powder Technol., 2022, 407, 117662.
28. B. Mielniczuk, O. Millet, G. Gagneux and M. S. El Youssoufi, Granular Matter, 2018, 20, 14.
29. J.-P. Wang, E. Gallo, B. François, F. Gabrieli and P. Lambert, Powder Technol., 2017, 305, 89–98.
30. N. P. Kruyt and O. Millet, J. Fluid Mech., 2017, 812, 129–151.
31. G. Lian and J. Seville, Adv. Colloid Interface Sci., 2016, 227, 53–62.
32. C.-F. Zhao, N. P. Kruyt and O. Millet, Powder Technol., 2018, 339, 827–837.
33. C.-F. Zhao, N. P. Kruyt and O. Millet, Powder Technol., 2020, 360, 622–634.
34. H. Chen, A. Ponce-Torres, J. M. Montanero and A. Amirfazli, J. Colloid Interface Sci., 2021, 582, 1251–1256.
35. H. Pingulkar, J. Peixinho and O. Crumeyrolle, Langmuir, 2021, 37, 10348–10353.
36. A. K. Sankaran and J. P. Rothstein, J. Non-Newtonian Fluid Mech., 2012, 175–176, 64–75.
37. J. A. Lee, J. P. Rothstein and M. Pasquali, J. Non-Newtonian Fluid Mech., 2013, 199, 1–11.
38. V. Entov and E. Hinch, J. Non-Newtonian Fluid Mech., 1997, 72, 31–53.
39. A. B. Bazilevskii and A. N. Rozhkov, Fluid Dyn, 2015, 50, 800–811.
40. G. H. McKinley and T. Sridhar, Annu. Rev. Fluid Mech., 2002, 34, 375–415.
41. S. L. Anna, G. H. McKinley, D. A. Nguyen, T. Sridhar, S. J. Muller, J. Huang and D. F. James, J. Rheol., 2001, 45, 83–114.
42. L. E. Rodd, T. P. Scott, J. J. Cooper-White and G. H. McKinley, Appl. Rheol., 2005, 15, 12–27.
43. G. H. McKinley and A. Tripathi, J. Rheol., 2000, 44, 653–670.
44. S. Rajesh, V. Thiévenaz and A. Sauret, Soft Matter, 2022, 18, 3147–3156.
45. P. J. Dekker, M. A. Hack, W. Tewes, C. Datt, A. Bouillant and J. H. Snoeijer, Phys. Rev. Lett., 2022, 128, 028004.
46. J. Dinic, L. N. Jimenez and V. Sharma, Lab Chip, 2017, 17, 460–473.
47. A. Gaillard, M. A. Herrada, A. Deblais, C. van Poelgeest, L. Laruelle, J. Eggers and D. Bonn, J. Fluid Mech., 2025, 1005, A10.
48. A. Gaillard, M. A. Herrada, A. Deblais, J. G. Eggers and D. Bonn, Phys. Rev. Fluids, 2024, 9, 073302.
49. Y. Li and J. E. Sprittles, J. Fluid Mech., 2016, 797, 29–59.
50. V. M. Entov and A. L. Yarin, Fluid Dyn, 1984, 19, 21–29.
51. Y. Amarouchene, D. Bonn, J. Meunier and H. Kellay, Phys. Rev. Lett., 2001, 86, 3558.
52. V. Tirtaatmadja, G. H. McKinley and J. J. Cooper-White, Phys. Fluids, 2006, 18, 043101.
53. V. Calabrese, A. Q. Shen and S. J. Haward, Macromolecules, 2024, 57, 9668–9676.
54. V. Thiévenaz and A. Sauret, Phys. Rev. Fluids, 2021, 6, L062301.
55. S. Rajesh, S. S. Peddada, V. Thiévenaz and A. Sauret, J. Non-Newtonian Fluid Mech., 2022, 310, 104921.
56. J. Dinic and V. Sharma, Proc. Natl. Acad. Sci. U. S. A., 2019, 116, 8766–8774.
57. K. Zinelis, T. Abadie, G. H. McKinley and O. K. Matar, Soft Matter, 2024, 20, 8198–8214.
58. C. Clasen, J. P. Plog, W.-M. Kulicke, M. Owens, C. Macosko, L. E. Scriven, M. Verani and G. H. McKinley, J. Rheol., 2006, 50, 849–881.
59. L. Campo-Deaño and C. Clasen, J. Non-Newtonian Fluid Mech., 2010, 165, 1688–1699.
60. A. Aisling, R. Saraka and N. J. Alvarez, J. Non-Newtonian Fluid Mech., 2024, 333, 105321.
61. K. Zinelis, T. Abadie, G. H. McKinley and O. K. Matar, J. Fluid Mech., 2024, 998, A4.
62. M. Huisman, P. Digard, W. C. Poon and S. Titmuss, Phys. Rev. Lett., 2023, 131, 248102.
63. V. Calabrese, A. Q. Shen and S. J. Haward, Phys. Rev. X, 2025, 15, 021025.
64. W. H. Cao and M. Won Kim, Faraday Discuss., 1994, 98, 245–252.
65. F. Xiao, J. Jing, S. Kuang, L. Yang and A. Yu, Powder Technol., 2020, 363, 59–73.
66. F. Orr, L. Scriven and A. P. Rivas, J. Fluid Mech., 1975, 67, 723–742.
67. D. Wu, P. Zhou, G. Wang, B. Zhao, T. Howes and W. Chen, Powder Technol., 2020, 376, 390–397.
68. J. Eggers, Rev. Mod. Phys., 1997, 69, 865–930.
69. C. D. Willett, M. J. Adams, S. A. Johnson and J. P. K. Seville, Powder Technol., 2003, 130, 63–69.
70. D. Megias-Alguacil and L. J. Gauckler, AIChE J., 2009, 55, 1103–1109.
71. P. G. D. Gennes, J. Chem. Phys., 1974, 5030–5042.
72. S. T. Nase, W. L. Vargas, A. A. Abatan and J. McCarthy, Powder Technol., 2001, 116, 214–223.
73. J. Van Aeken, L. Passaro and C. Clasen, Rheol. Acta, 2022, 61, 191–206.
74. C. O. Klein, I. F. Naue, J. Nijman and M. Wilhelm, Soft Mater., 2009, 7, 242–257.
75. R. S. Sharma and A. Sauret, Soft Matter, 2025, 21, 2193.
76. V. Richefeu, M. S. El Youssoufi and F. Radjai, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2006, 73, 051304.
77. R. S. Sharma, A. Leonelli, K. Zhao, E. Meiburg and A. Sauret, Phys. Rev. Lett., 2026, 136, 068204.
78. P. Philippe and M. Badiane, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2013, 87, 042206.
79. H. Rumpf, Chem. Ing. Tech., 1974, 46, 1–11.
80. N. Mitarai and F. Nori, Adv. Phys., 2006, 55, 1–45.
81. S. Herminghaus, Adv. Phys., 2005, 54, 221–261.
82. A. Gans, A. Abramian, P.-Y. Lagrée, M. Gong, A. Sauret, O. Pouliquen and M. Nicolas, J. Fluid Mech., 2023, 959, A41.
83. R. S. Sharma, W. Sarlin, L. Xing, C. Morize, P. Gondret and A. Sauret, Phys. Rev. Fluids, 2024, 9, 074301.
84. C. Wagner, Y. Amarouchene, D. Bonn and J. Eggers, Phys. Rev. Lett., 2005, 95, 164504.
85. R. Sattler, S. Gier, J. Eggers and C. Wagner, Phys. Fluids, 2012, 24, 023101.
86. R. Sattler, C. Wagner and J. Eggers, Phys. Rev. Lett., 2008, 100, 164502.
87. A. Deblais, K. Velikov and D. Bonn, Phys. Rev. Lett., 2018, 120, 194501.
88. Q. Song, K. Liu, W. Sun, R. Chen, J. Ji, Y. Jiao, T. Gao and J. Ye, Langmuir, 2021, 37, 11737–11749.
89. Q. Song, J. Ye, S. Zhang, J. Yuan, R. Chen and K. Liu, Colloids Surf., A, 2025, 709, 136098.
90. H.-J. Butt and M. Kappl, Adv. Colloid Interface Sci., 2009, 146, 48–60.
91. N. Chatterjee and M. Flury, Langmuir, 2013, 29, 7903–7911.
92. Y. Guo and J. S. Curtis, Annu. Rev. Fluid Mech., 2015, 47, 21–46.
93. F. Radjai, J.-N. Roux and A. Daouadji, J. Eng. Mech., 2017, 143, 04017002.
94. A. Argilaga and C. Zhao, Powder Technol., 2023, 427, 118702.
95. M. Bagheri, S. Roy and T. Pöschel, Comput. Part. Mech., 2024, 11, 2179–2190.

---