lidiffuseinterfacemarangoniinstability2025.md 39 KB


category: literaturenote citekey: lidiffuseinterfacemarangoniinstability2025 title: A Diffuse-Interface Marangoni Instability authors: "Li, Xiangwei; Wan, Dongdong; Hao, Haohao; Diddens, Christian; Zhang, Mengqi; Tan, Huanshu" year: 2025 date: 2025-00-00 2025 doi: 10.1017/jfm.2025.10895 url: "https://arxiv.org/abs/2506.09945v1" zotero_key: ABSSP5EU zotero_storage: 2M9XD4XS collections: transpordi_planeerimine folder: 001_artiklid firstAuthor: "Li, Xiangwei"

status: converted

A Diffuse-Interface Marangoni Instability

– 1

Xiangwei Li1 , Dongdong Wan2 , Haohao Hao1 , Christian Diddens3 , Mengqi Zhang2 , and Huanshu Tan1 †

1Multicomponent Fluids group, Center for Complex Flows and Soft Matter Research, Department of Mechanics and Aerospace Engineering, Southern University of Science and Technology, Shenzhen, 518055, Guangdong, P.R. China

2Department of Mechanical Engineering, National University of Singapore, 9 Engineering Drive 1, 117575, Republic of Singapore

3Physics of Fluids Department, Max-Planck Center Twente for Complex Fluid Dynamics and J. M. Burgers Centre for Fluid Dynamics, University of Twente, P. O. Box 217, 7500AE Enschede, The Netherlands

We investigate a novel Marangoni-induced instability that arises exclusively in diffuse fluid interfaces, absent in classical sharp-interface models. Using a validated phase-field Navier–Stokes–Allen–Cahn framework, we linearize the governing equations to analyze the onset and development of interfacial instability driven by solute-induced surface tension gradients. A critical interfacial thickness scaling inversely with the Marangoni number, δcr ∼ M a−1 , emerges from the balance between advective and diffusive transport. Unlike sharp-interface scenarios where matched viscosity and diffusivity stabilize the interface, finite thickness induces asymmetric solute distributions and tangential velocity shifts that destabilize the system. We identify universal power-law scalings of velocity and concentration offsets with a modified Marangoni number M aδ , independent of capillary number and interfacial mobility. A critical crossover at M aδ ≈ 590 distinguishes diffusiondominated stabilization from advection-driven destabilization. These findings highlight the importance of diffuse-interface effects in multiphase flows, with implications for miscible fluids, soft matter, and microfluidics where interfacial thickness and coupled transport phenomena are non-negligible.

Key words: Marangoni instability, Diffuse interface, Multicomponent fluids

1. Introduction

Marangoni instability, arising from the transverse mass transfer of surface-active species across a fluid interface (Sternling & Scriven 1959; Kovalchuk & Vollhardt 2006; Schwarzenberger et al. 2014), plays a fundamental role in multicomponent flow systems, such as evaporating droplets (Tan et al. 2016; Lohse & Zhang 2020; Diddens et al. 2024), self-motived droplets (Levich & Krylov 1969; Michelin 2023), microextraction solutions (Eckert et al. 2012; Mitra et al. 2020), and liquid films (Darhuber & Troian 2003; Sultan et al. 2005). Its onset and evolution are governed by fluid property contrasts, particularly viscosity and mass diffusivity ratios between the interacting layers. Traditional linear stability analysis, based on an idealized sharp and non-deformable interface, establishes critical instability criteria (Sternling & Scriven 1959; Reichenbach & Linde 1981). However, direct numerical simulations with a finite-thickness, deformable

Figure 1. (a) Schematic of a two-liquid layer system with a constant gradient of solute ∇c in y-direction. (b) Phase-field variable ϕeq distribution in the diffuse-interface method.

interface validate these criteria while also revealing that excessive interface thickness can itself induce Marangoni instability (Verschueren et al. 2001; Li et al. 2025). This study seeks to further clarify the role of interface thickness in Marangoni instability and the underlying mechanisms.

A perfectly sharp fluid interface – one with zero thickness – does not exist. Even in immiscible liquids, the phase transition region spans nanometer scales (Mitrinovi´c et al. 2000; Senapati & Berkowitz 2001). Jump conditions are commonly employed to approximate interface by neglecting the continuous variation of fluid properties across it (Scardovelli & Zaleski 1999; Reder et al. 2024), provided that the characteristic system size is much larger than the interfacial thickness ˆδ. However, this assumption breaks down when the interfacial thickness becomes comparable to the length scale of the phenomena being examined (Anderson et al. 1998) – either due to interface thickening in miscible or near-critical fluids, or the shrinking of the system in sub-micron liquid systems (Cahn & Hilliard 1958; Stanley 1987; Joseph et al. 1993). In such cases, the transition occurs over a finite region still large compared to the molecular scale, motivating the diffuse-interface concept originally introduced by Rayleigh (1892) and van der Waals (1893), as illustrated in Figure 1a. In such cases, the variation of properties within the interface can no longer be ignored (Cahn & Hilliard 1958).

Interfacial surface tension arises from tangential forces across the finite-thickness transition layer between two fluids, where anisotropic intermolecular forces act parallel to the interface (Berry 1971; Marchand et al. 2011). Surface-active solutes locally reduce tension according to a first-order approximation, γ = γ0 − βc, where c is the solute concentration and β is a proportionality constant. While classical theories assume a sharp interface, this concept extends to finite-thickness interfaces, where surface force SF may develop heterogeneously across sublayers (Fig. 1a inset), driving localized solutal Marangoni flows (Anderson et al. 1998; Hao et al. 2025). These flows typically form symmetric convection rolls centered at stagnation points (Fig. 1b), and their multi-scale interactions lead to interfacial turbulence (Ruckenstein & Berbente 1964). The linear theory by Sternling & Scriven (1959), based on the sharp-interface model, highlights the roles of viscosity and diffusivity ratios between source (ϕs) and receiver (ϕr) fluids, but is intrinsically incapable of capturing the coupled process (Manikantan & Squires 2020) within diffuse interfaces (Fig. 1b inset). The interfacial thickness impact on the onset and development of Marangoni instability warrants further study.

Phase-field methods model fluid interfaces with finite thickness and have been successfully applied to both diffuse-interface problems and sharp-interface two-phase Navier-Stokes flows (Anderson et al. 1998; Jacqmin 1999). Derived from fluid free energy formulations (Cahn & Hilliard 1958), the phase field recasts the moving-boundary problem into a continuous framework. In simulations, interface distortion due to flow is regulated by the mobility parameter Mϕ, balancing anti-diffusion and convection of the phase field ϕ (Demont et al. 2023). As the interface thickness ˆδ decreases, phase-field solutions converge toward sharp-interface results (Magaletti et al. 2013). The physicbased formulation and its asymptotic consistency make the phase-field method well-suited for capturing interfacial flow dynamics with finite-thickness effects.

The paper is organized as follows. Section 2 introduces the problem formulation using a phase-field model and its linearization. Section 3 identifies the instability unique to diffuse interfaces, establishes a scaling law for the critical thickness, and reveals the underlying mechanisms through interface offset analysis. Section 4 concludes with key findings and their physical implications.

2. Problem Formulation and Linearization

We use a validated phase-field Navier–Stokes–Allen–Cahn (NSAC) model to study solute-driven interfacial instability (Li et al. 2025). A brief description follows.

As depicted in Figure 1b, the system is characterized by a solute concentration gradient Gc oriented normal to a diffuse interface separating two fluid phases. These fluids are denoted by ϕs and ϕr, corresponding to the solute-supplying and solute-receiving phases, respectively. Both phases are confined between two flat, no-flow-penetration (NFP) boundaries, while periodic boundary conditions apply laterally. The characteristic length is defined as one quarter of the spacing distance, i.e., L = H/4. The interface is captured using the conservative Allen-Cahn equation in dimensionless form, which accounts for advective transport,

$$\frac{\partial \phi_j}{\partial t} + \nabla \cdot (\phij \boldsymbol{u}) = \frac{1}{Pe{\phi}} \left{ \nabla \cdot \left[ \nabla \phi_j - B\phi_j \left( 1 - \phi_j \right) \frac{\nabla \phi_j}{|\nabla \phi_j|} \right] + \alpha_j \right}, \tag{2.1}$$

where j = s,r (no j summation), and B = 4/δ artanh(1−2δn) controls the dimensionless interfacial thickness δ = ˆδ/L, and δn defines the numerically resolvable transition width. The equilibrium phase-field profile takes the form ϕeq = {1 − tanh[(y − 2L)B/2]}/2 (Fig. 1a inset). The Lagrange multipliers αj enforce the constraint of ϕs + ϕr = 1.

The governing equations for fluid momentum and solute transport in dimensionless form are given by,

$$\frac{\partial \boldsymbol{u}}{\partial t} + \boldsymbol{u} \cdot \nabla \boldsymbol{u} = \frac{Sc}{Ma} \left{ -\nabla p + \nabla \cdot \left[ \eta \left( \nabla \boldsymbol{u} + (\nabla \boldsymbol{u})^T \right) \right] + \frac{SF}{Ca} \right}, \quad 0 = \nabla \cdot \boldsymbol{u}, \quad (2.2a)$$

$$\frac{\partial c}{\partial t} = \frac{1}{Ma} \nabla \cdot [D(\nabla c)] - \nabla \cdot (c\mathbf{u}), \qquad (2.2b)$$

where the surface force term SF = 2 [γ(κsns + κrnr) + (Is − nsns + Ir − nrnr)∇γ] W follows the formulation by Kim (2012), with nj = ∇ϕj/|∇ϕj |, κj = −∇ · nj , and (Ij − njnj )·∇ representing the surface gradient operator in terms of phase field ϕj . The surface tension decreases linearly with solute concentration under dilute conditions (Picardo et al. 2016), modeled as γ = 1 − Ca c, where Ca is the capillary number. The interface localization function is W = Aϕsϕr|∇ϕs||∇ϕr|, with A = 30/B. Viscosity and diffusivity are set as η = ϕs + ϕr/ζη and D = ϕs + ϕr/ζD. We choose ζη = ζD = 1 to isolate Marangoni instability driven by property ratios (Sternling & Scriven 1959). Gravity is neglected and density ratio is unity.

The dimensionless groups appearing in the equations are Schmidt number Sc = η/(ρD), capillary number Ca = Uη/γ0, and Marangoni number M a = UL/D. The characteristic velocity is defined as the Marangoni flow U = β|Gc|L/η that prevails in the stationary flow, with Gc = ∆c/(4L), β = dγ/dc. The P´eclet number for the phase field is given by P eϕ = UL/Mϕ, where the interfacial mobility Mϕ ∝ δ α, with an optimal

Figure 2. (a) Growth rate ωi versus wavenumber k for different δ. The gray-shaded area indicates ωi > 0, i.e., instability. Stability is restored when Marangoni force is off. Parameters: Ma = 104 , Ca = 0.01, α = 1.7. (b–e) Instability intervals (enclosed regions) under varying Sc, Ca, P eϕ, and Ma, respectively.

scaling parameter of α ≈ 1.7, describes the rate of convergence of the diffuse-interface solution to the sharp-interface solution (Demont et al. 2023).

The governing equations (Eqs. 2.1, 2.2) are linearized around a base state Q b = (U b , Pb , Cb , ϕb s ) T with zero velocity U b = 0, uniform pressure P b = const, an equilibrium phase field ϕ b s = ϕeq, and a steady solute gradient C b = Gc(H − y). The perturbations are assumed to take the form q ′ = q˜(y)eikx−iωt + c.c., where k is the wavenumber, and ω = ωr + iωi is the complex frequency. A positive growth rate (ωi > 0) indicates linear instability, while ωi < 0 implies stability. The real part ωr represents the oscillation frequency. This leads to a generalized eigenvalue problem of the form, ωMfq˜ = Leq˜, where Mf and Le are linear operators derived from the linearized equations as following,

$$-i\omega\tilde{\phi}{s} = -\frac{\mathrm{d}\phi{s}^{b}}{\mathrm{d}y}\tilde{v} + \frac{1}{Pe{\phi}}\frac{\mathrm{d}^{2}}{\mathrm{d}y^{2}}\tilde{\phi}{s} - \frac{B}{Pe{\phi}}\left[ (2\phi{s}^{b} - 1)\frac{\partial\tilde{\phi}{s}}{\partial y} + 2\frac{\mathrm{d}\phi{s}^{b}}{\mathrm{d}y}\tilde{\phi}_{s} \right],\tag{2.3a}$$

$$-i\omega\tilde{c} = -\frac{\mathrm{d}C^b}{\mathrm{d}y}\tilde{v} + \frac{1}{Ma}(\frac{\mathrm{d}^2}{\mathrm{d}y^2} - k^2)\tilde{c}, \tag{2.3b}$$

$$-\mathrm{i}\omega\tilde{u} = \frac{Sc}{Ma} \left[ -\mathrm{i}k\tilde{p} + (\frac{\mathrm{d}^2}{\mathrm{d}y^2} - k^2)\tilde{u} + \frac{W^b}{Ca} \frac{1}{|\nabla \phi\mathrm{s}^b|} \frac{\mathrm{d}\Gamma^b}{\mathrm{d}y} \mathrm{i}k\tilde{\phi}\mathrm{s} - \mathrm{i}k\tilde{c}W^b \right], \tag{2.3c}$$

$$-i\omega\tilde{v} = \frac{Sc}{Ma} \left[ -\frac{\partial \tilde{p}}{\partial y} + (\frac{d^2}{dy^2} - k^2)\tilde{v} + \frac{W^b}{Ca} \frac{1}{|\nabla \phi_s^b|} \Gamma^b(-k^2)\tilde{\phi}_s \right], \tag{2.3d}$$

$$0 = ik\tilde{u} + \frac{\partial \tilde{v}}{\partial y},\tag{2.3e}$$

where the base-state of surface tension Γ b = 1 − CaCb and weight function Wb = Aϕb s (1 − ϕ b s ) ∇ϕ b s . We refer the reader to Li et al. (2025) for the detailed derivation.

3. Results and Discussions

3.1. Unveiling instabilities in diffuse interfaces

We identify an instability that arises solely in diffuse interfaces and is absent in the sharp-interface limit, driven by Marangoni stresses. We term this the diffuse-interface Marangoni instability.

As shown in Figure 2a, short-wavelength instabilities emerge when the interface has finite thickness, even under conditions where sharp-interface theory predicts stability (equal viscosity and diffusivity in both phases) (Sternling & Scriven 1959; Schwarzenberger et al. 2014). For M a = 104 , Ca = 0.01, α = 1.7, and Sc = 1000, an unstable band appears at δ = 0.2 and narrows as δ decreases, disappearing entirely as δ approaches 0.02—an estimated sharp-interface limit. Disabling Marangoni stresses by omitting the −ikcW˜ b term in Eq. 2.3c suppresses the instability, yielding a single δ-independent growth curve (blue triangles in Fig.2a). This confirms that the observed instability stems from solutal Marangoni effects rather than numerical artifacts.

To probe the instability mechanism, we vary key parameters individually: Marangoni number M a ∈ [103 , 105 ], capillary number Ca ∈ [10−3 , 10−1 ], Schmidt number Sc ∈ [102 , 104 ], and mobility index α ∈ [1.3, 1.7]. These span Re ∈ [1, 100] and P eϕ ∈ [0.3, 1.6 with Re = M a/Sc and P eϕ = UL/δα. Default values are M a = 104 , Ca = 0.01, α = 1.7, and Sc = 1000.

Figure 2b-e show the neutral stability curves kNS versus interfacial thickness δ. Sc has little influence, consistent with prior studies (Reichenbach & Linde 1981; Li et al. 2025). In contrast, M a, Ca, and P eϕ ∼ δ −α significantly alter the instability. Lower Ca or higher α suppress short-wavelength modes but leave the critical thickness δcr unchanged. In contrast, reducing M a narrows the unstable wavenumber band and increases δcr, shrinking the unstable regime. These results indicate that while δ, M a, Ca, and P eϕ govern instability onset, δcr is uniquely determined by M a.

3.2. Scaling laws for critical diffuse-interface thickness

In the laminar regime, the critical thickness δcr depends uniquely on M a, not Ca or P eϕ, indicating that δcr is primarily governed by solutal Marangoni advection–diffusion dynamics near the interface (Eqns. 2.2b and 2.3b). The circulation (Fig. 1b) transports solute toward the interface, steepening concentration gradients, while diffusion tends to smooth them. This competition features the Marangoni interfacial flows (Manikantan & Squires 2020). Balancing diffusion (timescale tD ∼ (δL) 2/D) and advection (tA ∼ δL/U) yields the scaling,

$$\delta_{\rm cr} \sim M a^{-1}$$ . (3.1)

As M a increases, stronger advection requires a thinner interface to maintain stability. This M a - δcr scaling is confirmed as δcr → 0 in Figure 3a, and remains robust under changes in Ca or P eϕ, as shown in Figures 2cd. This relation provides a physical criterion for selecting δ in diffuse-interface simulations with Marangoni effects, ensuring convergence to the sharp-interface limit Jacqmin (1999). While thin-interface limit criteria exist for passive two-phase flows (Magaletti et al. 2013; Demont et al. 2023), the present finding extends them to Marangoni-driven instabilities, consistent with thermocapillary flow observations by Verschueren et al. (2001).

To explain the deviation from the scaling law (Eqn. 3.1) at larger interfacial thickness (i.e., for M a ≲ 4000), we perform an energy-based analysis to quantify the competing mechanisms. Multiplying the linearized perturbation equations (Eqns. 2.3c and 2.3d) by the velocity perturbations yields an energy budget,

$$\frac{\partial E{\text{kin}}}{\partial t} = \underbrace{\frac{1}{Re} \left[ - (\nabla u')^2 - (\nabla v')^2 \right]}{\mathcal{E}_{\eta}: \text{ Viscous dissipation}} + \underbrace{\frac{1}{Re} \frac{1}{Ca} \left[ \frac{\Gamma^b}{|\nabla \phi_s^b|} \frac{\partial^2 \phis'}{\partial x^2} v' \right] W^b}{\mathcal{E}_{\gamma}^{\text{cap}}: \text{ Capillary effect}} + \underbrace{\frac{1}{Re} \frac{1}{Ca} \left[ \frac{1}{|\nabla \phi_s^b|} \frac{\partial \phis'}{\partial x} \frac{d\Gamma^b}{dy} u' \right] W^b}{\mathcal{E}{\gamma}^{\text{talt}}: \text{ Tilt--induced Marangoni}} + \underbrace{\frac{1}{Re} \left[ -\frac{\partial c'}{\partial x} u' \right] W^b}{\mathcal{E}_{\gamma}^{\text{mag}}: \text{ Solutal Marangoni effect}} . \tag{3.2}$$

Figure 3. (a) Scaling of the critical interface thickness δcr with Ma, showing δcr = 590/Ma for Ma ≳ 4000, independent of Ca and α. Deviations from the scaling emerge with increasing α or decreasing Ca. (b) Energy budget showing the relative contributions of viscous dissipation, capillary forces, and Marangoni stresses as δcr increases. (c,d) Base-state phase field ϕ b s and solute perturbation c for the stable case (δ = 0.05) and unstable case (δ = 0.8). Velocity vectors are overlaid on ϕ b s ; streamlines are shown in c ′ . (e) Cross-sectional profiles of u (at x = π) and c (at x = π/2) along the y-direction for the stable (solid lines) and unstable (dashed lines) cases.

Here, the kinetic energy density is defined as Ekin = ρ b (u ′2 + v ′2 )/2, since the base flow U b is zero and total velocity reduces to the perturbation field (u = u ′ ). The prefactor M a/Sc is replaced by Re ≲ 4.

The viscous dissipation Eη reflects energy loss due to shear. The solutal Marangoni term E mag γ , originates from tangential solute gradients −∂c′/∂x. The capillary term E cap γ arises from interface curvature ∂ 2ϕ ′ u/∂x2 and base-state tension Γ b . The tilt-induced Marangoni term E tilt γ stems from interfacial tilt ∂ϕ′ u/∂x, which projects normal interfacial tension gradients dΓ b/dy across the diffuse layer (Berry 1971; Marchand et al. 2011) onto the tangential direction – an effect inherently tied to the finite thickness of diffuse interfaces and absent in sharp-interface models.

Along the M a - δcr paths (Fig. 3a) , integrating over the domain [2π/k × 4L] gives the power balance ⟨Eη⟩ + ⟨Ecap γ ⟩ + ⟨Etilt γ ⟩ + ⟨Emag γ ⟩ = 0 (ωi = 0). Divergence terms vanish under periodic, homogeneous boundary conditions (Zhang et al. 2013). Figure 3b shows the normalized contributions of each term by ⟨Emag γ ⟩ as δcr increases. For small δcr (≲ 0.1), the balance is primarily between viscous dissipation Eη and solutal Marangoni forcing E mag γ , consistent with the scaling law. As δcr increases, the influence of capillary E cap γ and tilt-induced Marangoni effects E tilt γ becomes significant at a certain δcr, where the deviation from the ideal δcr ∼ M a−1 relation starts exactly. Thus, interface deformation and non-uniform tension projection across the diffuse layer play key roles in setting the limits of the scaling law's validity.

Figures 3cd show the base-state phase filed ϕ b s (blue), velocity field u = u ′ (arrows), and perturbation concentration field c ′ (red) for two representative cases: a stable system (δ = 0.05) and an unstable one (δ = 0.8), corresponding to the solid and open star symbols in Figure 3a, respectively. In the stable case (Fig. 3c), vortices are symmetrically

Figure 4. (a) Schematic of interfacial solute and velocity fields in the sharp-interface limit, where symmetric advection yields no net concentration change. (b) Profiles of solute concentration C(y) and tangential velocity u(y) across the diffuse interface.(c) Mechanism of velocity offset ℓ : asymmetric solute advection enhances Marangoni flow on the source side, displacing the velocity peak away from the centerline. (d) Mechanism of solute offset ℓ δ : coupled vertical flow and interfacial diffusion shift the position of zero net solute change. (e, f) Normalized offsets ℓ δ /δcr and ℓ δ /δcr versus the modified Marangoni number Maδ = Ma δ, showing collapse across varying Ca and α. Instability sets in at Maδ cr ≈ 590. (g) Phase diagram in the ℓ u δ – ℓ plane. Stable (solid) and unstable (open) regimes are separated by the scaling ℓ u δ ∼ (ℓ ).

distributed about the interface, as evidenced by the velocity component u profiles at x = π (Fig. 3e, black solid line). This balanced Marangoni recirculation gives rise to an anti-symmetric c ′ field, confirmed by the profile at x = π/2 (Fig. 3e, red solid line). In contrast, the unstable case (Fig. 3d) exhibits anti-symmetry breaking, with flow field shifting downward into source phase (Fig. 3e, black dashed line). The zero-crossing point of c ′ shifts upward from the interface centerline (y = 2), indicating solute enrichment in the receiver phase and depletion in the source phase (Fig. 3e red dashed line). These symmetry-breaking shifts in u and c ′ motivate the next section, in which we investigate their physical origin and impact on the instability.

3.3. Mechanism of instability in diffuse interfaces: Anti-flux vs. Pro-flux shifts

In classical sharp-interface systems with equal viscosity and diffusivity between two fluids, interfacial Marangoni instability is absent (Sternling & Scriven 1959). As illustrated in Figure 4a, symmetric convection brings high (C H) and low (C L) concentrations toward the interfacial stagnation point. Equal diffusivity and viscosity ensure that net interfacial concentration remains unchanged, yielding no amplification of the Marangoni flow um. In contrast, a diffuse interface introduces local asymmetry. Figure 4b shows that inflow of low-concentration fluid reduces interfacial concentration C − on the receiver phase ϕr side, while inflow of high-concentration fluid increases C + on the source ϕs side. Within the interfacial layer of thickness δ, the resulting imbalance cannot be instantly homogenized, producing a tangential concentration gradient C(y) and velocity profile u(y).

As shown in figure 4c, the local concentration decrease on the receiver side (ϕr) reduces the interfacial tension gradient (∇sγ) −, thereby weakening the Marangoni flow (Anderson et al. 1998). Conversely, the increased concentration on the source side (ϕs) enhances (∇sγ) +, strengthening the flow. This asymmetry in surface tension gradients induces a shear layer within the interface, with suppressed velocity near ϕr and enhanced velocity near ϕs. The resulting shear promotes an interfacial instability reminiscent of the Kelvin–Helmholtz type (Chandrasekhar 2013). We quantify the shift of the peak velocity (um) away from the interface centerline (y = 2L) by the anti-flux offset, ℓ u δ = 2L − yum, which characterizes the destabilizing effect of advection-driven shear. Diffusion counters this destabilization by smoothing concentration gradients (C − ↑ and C + ↓), as depicted in Figure 4d. The position where the net concentration change vanishes (C 0 ) shifts along the solute flux into the receiver side. This stabilizing effect is captured by the pro-flux offset ℓ c δ = yC0 − 2L.

To quantify the competition between advection and diffusion across the interface, we introduce a modified Marangoni number M aδ = M a δ, using the interfacial thickness ˆδ as the characteristic length scale. As shown in Figures 4ef, both the velocity offsets ℓ u δ and the concentration offset ℓ ′ δ increase monotonically with M aδ , following distinct powerlaw trends that collapse across different Ca ∈ [10−3 , 10−1 ] and α ∈ [1.3, 1.7], confirming their universality. A critical transition emerges at M aδ ≈ 590: below this threshold, the system remain stable, with ℓ u δ growing sublinearly and ℓ c ′ δ superlinearly; Beyond it, the system becomes unstable, and both offsets grow with fixed exponents, approximately 1.13 for ℓ u δ and 0.75 for ℓ c δ .

To further capture this interplay, Figure 4g plots ℓ u δ against ℓ c δ , revealing a clear separation between stable (bluish solid) and unstable (reddish open) regimes by the power-law relation ℓ u δ ∼ (ℓ c δ ) m. The critical exponent m = 1 at M aδ = 590 (black solid line with star symbols) marks a crossover from diffusion-dominated stabilization (m < 1) to advective-driven destabilization (m > 1), highlighting the offset competition as the underlying mechanism governing interfacial instability.

4. Conclusion

We have uncovered a Marangoni-induced diffuse-interface instability that is absent in the sharp-interface limit. Through linearized analysis of a validated phase-field NSAC model, we show that finite interfacial thickness introduces asymmetric solute distributions and tangential velocity shifts that destabilize the interface, even in fluids with matched viscosity and diffusivity. A critical interfacial thickness δcr ∼ M a−1 emerges from the competition between diffusive and advective timescales. Energy analysis reveals that deviations at large δcr arise from interface deformation and asymmetric tension projection. By introducing a modified Marangoni number M aδ = M aδ, we reveal universal power-law behavior in velocity and concentration offsets, independent of Ca and α. The crossover at M aδ ≈ 590, where the exponent ℓ u δ ∼ (ℓ c δ ) m with m = 1, marks the transition from diffusive stabilization to advective destabilization.

These findings not only clarify the physical origin of instability in diffuse interfaces but also provide a predictive framework for multicomponent multiphase systems where sharp-interface assumptions fail. This work opens new directions for studying interfacial transport in miscible fluids, soft materials, and microfluidic systems (Shim 2022; Michelin 2023; Ault & Shin 2024), where interfacial thickness and multi-physics coupling critically influence dynamics.

H.T. acknowledges support from the National Natural Science Foundation of China (Grant Nos. 12472271, 12102171) and the Guangdong Basic and Applied Basic Research Foundation (Grant No. 2024A1515010614). M.Z. acknowledges support from the Ministry of Education, Singapore (WBS No. A-8001172-00-00).

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

REFERENCES

  • Anderson, Daniel M, McFadden, Geoffrey B & Wheeler, Adam A 1998 Diffuse-interface methods in fluid mechanics. Annual review of fluid mechanics 30 (1), 139–165.
  • Ault, Jesse T & Shin, Sangwoo 2024 Physicochemical hydrodynamics of particle diffusiophoresis driven by chemical gradients. Annual Review of Fluid Mechanics 57.
  • Berry, Michael V 1971 The molecular mechanism of surface tension. Physics Education 6 (2), 79.
  • Cahn, John W & Hilliard, John E 1958 Free energy of a nonuniform system. i. interfacial free energy. The Journal of chemical physics 28 (2), 258–267.
  • Chandrasekhar, Subrahmanyan 2013 Hydrodynamic and hydromagnetic stability. Courier Corporation.
  • Darhuber, Anton A & Troian, Sandra M 2003 Marangoni driven structures in thin film flows. In Gallery of Fluid Motion, p. S9.
  • Demont, T.H.B., Stoter, S.K.F. & Van Brummelen, E.H. 2023 Numerical investigation of the sharp-interface limit of the Navier-Stokes-Cahn-Hilliard equations. Journal of Fluid Mechanics 970, A24.
  • Diddens, Christian, Dekker, Pim J & Lohse, Detlef 2024 Non-monotonic surface tension leads to spontaneous symmetry breaking in a binary evaporating drop. arXiv preprint arXiv:2402.17452 .
  • Eckert, K., Acker, M., Tadmouri, R. & Pimienta, V. 2012 Chemo-Marangoni convection driven by an interfacial reaction: Pattern formation and kinetics. Chaos: An Interdisciplinary Journal of Nonlinear Science 22 (3), 037112.
  • Hao, Haohao, Li, Xiangwei, Liu, Tian & Tan, Huanshu 2025 Enhanced profile-preserving phase-field model of two-phase flow with surfactant interfacial transport and marangoni effects. Journal of Computational Physics p. 114058.
  • Jacqmin, David 1999 Calculation of Two-Phase Navier Stokes Flows Using Phase-Field Modeling. Journal of Computational Physics 155 (1), 96–127.
  • Joseph, Daniel D, Renardy, Yuriko Y, Joseph, Daniel D & Renardy, Yuriko Y 1993 Fluid dynamics of two miscible liquids with diffusion and gradient stresses. Fundamentals of Two-Fluid Dynamics: Part II: Lubricated Transport, Drops and Miscible Liquids pp. 324–395.
  • Kim, Junseok 2012 Phase-Field Models for Multi-Component Fluid Flows. Communications in Computational Physics 12 (3), 613–661.
  • Kovalchuk, Nina & Vollhardt, D. 2006 Marangoni instability and spontaneous non-linear oscillations produced at liquid interfaces by surfactant transfer. Advances in Colloid and Interface Science .
  • Levich, VG & Krylov, VS 1969 Surface-tension-driven phenomena. Annual review of fluid mechanics 1 (1), 293–316.
  • Li, Xiangwei, Wan, Dongdong, Zhang, Mengqi & Tan, Huanshu 2025 Marangoni interfacial instability induced by solute transfer across liquid–liquid interfaces. International Journal of Multiphase Flow 184, 105073.
  • Lohse, Detlef & Zhang, Xuehua 2020 Physicochemical hydrodynamics of droplets out of equilibrium. Nature Reviews Physics 2 (8), 426–443.
  • Magaletti, F., Picano, F., Chinappi, M., Marino, L. & Casciola, C. M. 2013 The sharpinterface limit of the Cahn–Hilliard/Navier–Stokes model for binary fluids. Journal of Fluid Mechanics 714, 95–126.

  • Manikantan, Harishankar & Squires, Todd M. 2020 Surfactant dynamics: hidden variables controlling fluid flows. Journal of Fluid Mechanics 892, P1.

  • Marchand, Antonin, Weijs, Joost H, Snoeijer, Jacco H & Andreotti, Bruno 2011 Why is surface tension a force parallel to the interface? American Journal of Physics 79 (10), 999–1008.

  • Michelin, Sebastien ´ 2023 Self-propulsion of chemically active droplets. Annual Review of Fluid Mechanics 55 (1), 77–101.

  • Mitra, Shirsendu, Singh, Sunil Kumar, Shevchenko, Ekaterina, Sachan, Mohit, Ghosh, Abir, Basak, Mitali & Gooh Pattader, Partho Sarathi 2020 Efficient microextraction process exploiting spontaneous interfacial convection driven by marangoni and electric field induced instability: A computational fluid dynamics study. Physics of Fluids 32 (1).

  • Mitrinovic, Dragoslav M, Tikhonov, Aleksey M, Li, Ming, Huang, Zhengqing & ´ Schlossman, Mark L 2000 Noncapillary-wave structure at the water-alkane interface. Physical review letters 85 (3), 582.

  • Picardo, Jason R., Radhakrishna, T. G. & Pushpavanam, S. 2016 Solutal marangoni instability in layered two-phase flows. Journal of Fluid Mechanics 793, 280–315.

  • Rayleigh, Lord 1892 Xx. on the theory of surface forces.—ii. compressible fluids. The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science 33 (201), 209–220.

  • Reder, Martin, Prahs, Andreas, Schneider, Daniel & Nestler, Britta 2024 Viscous stress approximations in diffuse interface methods for two-phase flow based on mechanical jump conditions. Computer Methods in Applied Mechanics and Engineering 432, 117341.

  • Reichenbach, J & Linde, H 1981 Linear perturbation analysis of surface-tension-driven convection at a plane interface (marangoni instability). Journal of Colloid and Interface Science 84 (2), 433–443.

  • Ruckenstein, Eli & Berbente, C 1964 The occurrence of interfacial turbulence in the case of diffusion accompanied by chemical reaction. Chemical Engineering Science 19 (5), 329– 347.

  • Scardovelli, Ruben & Zaleski, Stephane ´ 1999 Direct numerical simulation of free-surface and interfacial flow. Annual review of fluid mechanics 31 (1), 567–603.

  • Schwarzenberger, Karin, Kollner, Thomas, Linde, Hartmut, Boeck, Thomas, ¨ Odenbach, Stefan & Eckert, Kerstin 2014 Pattern formation and mass transfer under stationary solutal Marangoni instability. Advances in Colloid and Interface Science 206, 344–371.

  • Senapati, Sanjib & Berkowitz, Max L 2001 Computer simulation study of the interface width of the liquid/liquid interface. Physical review letters 87 (17), 176101.

  • Shim, Suin 2022 Diffusiophoresis, diffusioosmosis, and microfluidics: surface-flow-driven phenomena in the presence of flow. Chemical Reviews 122 (7), 6986–7009.

  • Stanley, H Eugene 1987 Introduction to phase transitions and critical phenomena.

  • Sternling, C. V. & Scriven, L. E. 1959 Interfacial turbulence: Hydrodynamic instability and the marangoni effect. AIChE Journal 5 (4), 514–523.

  • Sultan, Eric, Boudaoud, Arezki & Amar, Martine Ben 2005 Evaporation of a thin film: diffusion of the vapour and marangoni instabilities. Journal of Fluid Mechanics 543, 183–202.

  • Tan, Huanshu, Diddens, Christian, Lv, Pengyu, Kuerten, Johannes GM, Zhang, Xuehua & Lohse, Detlef 2016 Evaporation-triggered microdroplet nucleation and the four life phases of an evaporating ouzo drop. Proceedings of the National Academy of Sciences 113 (31), 8642–8647.

  • Verschueren, M., Van De Vosse, F. N. & Meijer, H. E. H. 2001 Diffuse-interface modelling of thermocapillary flow instabilities in a Hele-Shaw cell. Journal of Fluid Mechanics 434, 153–166.

  • van der Waals, J 1893 Verhandel. konink. akad. weten. amsterdam (sect1). Philos. Mag 1, 56.

  • Zhang, Mengqi, Lashgari, Iman, Zaki, Tamer A. & Brandt, Luca 2013 Linear stability analysis of channel flow of viscoelastic Oldroyd-B and FENE-P fluids. Journal of Fluid Mechanics 737, 249–279.