VLBNonlocal — Shear Banding PDE

Quick Reference

Model Summary

Model Class

VLBNonlocal

Physics

Nonlocal VLB with tensor diffusion for shear banding

Key Parameters

\(G_0, k_d^0, \eta_s, D_\mu\) (+ \(\nu\) for Bell, \(L_{\max}\) for FENE)

Protocols

FLOW_CURVE (steady shear PDE), STARTUP, CREEP

Key Features

Spatially-resolved 1D PDE, shear banding, cooperativity length

Reference

Dhont (1999). PRE 60, 4534; Vernerey et al. (2017). JMPS 107, 1-20

Import:

from rheojax.models import VLBNonlocal

Overview

VLBNonlocal solves the VLB constitutive equation in one spatial dimension, modeling the gap of a Couette rheometer. The distribution tensor \(\boldsymbol{\mu}(y, t)\) varies across the gap, coupled by a diffusion term that represents cooperative rearrangements in the network. For the underlying nonlocal PDE theory, see the nonlocal formulation section in VLB Advanced Theory & Numerical Methods.

The PDE is:

\[\frac{\partial \boldsymbol{\mu}}{\partial t} = k_d(\mathbf{I} - \boldsymbol{\mu}) + \mathbf{L} \cdot \boldsymbol{\mu} + \boldsymbol{\mu} \cdot \mathbf{L}^T + D_\mu \nabla^2 \boldsymbol{\mu}\]

where \(D_\mu\) (m2/s) is the distribution tensor diffusivity.

Shear banding arises when the Bell breakage rate creates a non-monotonic (S-shaped) constitutive curve \(\sigma(\dot\gamma)\). In the unstable region, the flow spontaneously separates into coexisting high- and low-shear-rate bands. The diffusion term regularizes the interface, setting its width to the cooperativity length:

\[\xi = \sqrt{D_\mu / k_d^0}\]

Constructor

model = VLBNonlocal(
    breakage="constant",     # or "bell" for banding
    stress_type="linear",    # or "fene"
    n_points=51,             # spatial grid resolution
    gap_width=1e-3,          # gap width in meters
)

Parameters

Parameter

When Present

Units

Description

\(G_0\)

Always

Pa

Network modulus

\(k_d^0\)

Always

1/s

Unstressed dissociation rate

\(\eta_s\)

Always

Pa·s

Solvent viscosity (regularization)

\(D_\mu\)

Always

m2/s

Distribution tensor diffusivity

\(\nu\)

breakage="bell"

Force sensitivity (Bell model)

\(L_{\max}\)

stress_type="fene"

Maximum chain extensibility (FENE-P)

Simulation Methods

simulate_steady_shear

Approach to steady state under imposed average shear rate. Uses a stress feedback controller to enforce the velocity constraint.

result = model.simulate_steady_shear(
    gamma_dot_avg=1.0,    # imposed average shear rate
    t_end=100.0,          # simulation time
    dt=1.0,               # output time step
    perturbation=0.05,    # initial symmetry-breaking noise
)

# Returns dict with keys:
# 't': time array
# 'y': spatial grid
# 'mu_xy': mu_xy profiles (N_t, N_y)
# 'gamma_dot': shear rate profiles (N_t, N_y)
# 'stress': wall stress Sigma(t)

simulate_startup

Startup from rest (equivalent to steady shear with small perturbation).

result = model.simulate_startup(
    gamma_dot_avg=1.0, t_end=50.0, dt=0.5
)

simulate_creep

Stress-controlled creep with spatial resolution.

result = model.simulate_creep(
    sigma_0=100.0, t_end=100.0, dt=0.1
)

Banding Detection

banding = model.detect_banding(result, threshold=0.1)

if banding["is_banding"]:
    print(f"Band contrast: {banding['band_contrast']:.1f}")
    print(f"Band width: {banding['band_width']*1e3:.2f} mm")
    print(f"Band location: {banding['band_location']*1e3:.2f} mm")

The detection checks the spatial variation of the final shear rate profile. threshold is the relative standard deviation cutoff.

Velocity Profile

v = model.get_velocity_profile(result)
# v(0) = 0, v(H) = gamma_dot_avg * H

The velocity profile is computed by integrating the shear rate: \(v(y) = \int_0^y \dot\gamma(y')\,dy'\).

Cooperativity Length

The cooperativity length \(\xi\) sets the shear band interface width:

xi = model.get_cooperativity_length()
print(f"Cooperativity length: {xi*1e6:.1f} um")

Grid Resolution

The number of spatial grid points (n_points) should resolve the band interface. A rule of thumb is \(\Delta y < \xi / 3\):

# Estimate required grid points
xi = model.get_cooperativity_length()
n_min = int(3 * model.gap_width / xi) + 1

# Recreate with sufficient resolution
model = VLBNonlocal(breakage="bell", n_points=max(n_min, 51))

For convergence studies, compare results with 25, 51, and 101 points.

Boundary Conditions

The model uses Neumann (zero-flux) boundary conditions:

\[\left.\frac{\partial \boldsymbol{\mu}}{\partial y}\right|_{y=0} = \left.\frac{\partial \boldsymbol{\mu}}{\partial y}\right|_{y=H} = 0\]

This means the structure tensor gradient vanishes at the walls.

Limitations

  • VLBNonlocal does not support fit() or fit_bayesian(). Use the simulate_* methods directly for forward predictions.

  • Computational cost scales as \(O(N_y \times N_t)\) where \(N_y\) is the number of spatial points.

  • Only 1D (planar Couette) geometry is supported.

See Also

  • VLB Advanced Theory & Numerical Methods — Nonlocal PDE theory, cooperativity length derivation

  • DMTNonlocal — DMT thixotropic model with similar nonlocal PDE structure (structure parameter diffusion instead of tensor diffusion)

  • FluidityNonlocal — Fluidity model with spatial diffusion of fluidity variable

API Reference

class rheojax.models.vlb.VLBNonlocal(breakage='constant', stress_type='linear', n_points=51, gap_width=0.001)[source]

Bases: VLBBase

Nonlocal VLB with tensor diffusion for shear banding.

Solves a 1D PDE across the gap of a Couette cell. The state at each spatial point is (mu_xx, mu_yy, mu_zz, mu_xy), plus a single wall stress Sigma (spatially uniform at low Reynolds number).

Shear banding occurs when the Bell breakage rate creates a non-monotonic flow curve. The diffusion term D_mu * nabla^2(mu) regularizes the interface with width ~ xi = sqrt(D_mu / k_d_0).

Parameters:
  • breakage (Literal['constant', 'bell']) – “constant” or “bell”

  • stress_type (Literal['linear', 'fene']) – “linear” or “fene”

  • n_points (int) – Spatial grid resolution

  • gap_width (float) – Gap width (m)

__init__(breakage='constant', stress_type='linear', n_points=51, gap_width=0.001)[source]

Initialize VLBNonlocal model.

property G0: float
property k_d_0: float
get_cooperativity_length()[source]

Cooperativity length xi = sqrt(D_mu / k_d_0).

This sets the shear band interface width.

Returns:

Cooperativity length (m)

Return type:

float

simulate_steady_shear(gamma_dot_avg, t_end=100.0, dt=0.1, perturbation=0.01)[source]

Simulate approach to steady state under imposed average shear rate.

Parameters:
  • gamma_dot_avg (float) – Imposed average shear rate (1/s)

  • t_end (float) – Simulation end time (s)

  • dt (float) – Output time step (s)

  • perturbation (float) – Initial spatial noise amplitude

Returns:

‘t’: time array ‘y’: spatial grid ‘mu_xy’: mu_xy profiles (N_t, N_y) ‘gamma_dot’: shear rate profiles (N_t, N_y) ‘stress’: wall stress Sigma(t)

Return type:

dict

simulate_startup(gamma_dot_avg, t_end=50.0, dt=0.05)[source]

Simulate startup from rest with banding evolution.

Parameters:
  • gamma_dot_avg (float) – Imposed average shear rate (1/s)

  • t_end (float) – End time (s)

  • dt (float) – Output time step (s)

Returns:

Same format as simulate_steady_shear

Return type:

dict

simulate_creep(sigma_0, t_end=100.0, dt=0.1)[source]

Simulate stress-controlled creep with spatial resolution.

In creep, the stress Sigma is held fixed (no feedback).

Parameters:
  • sigma_0 (float) – Applied stress (Pa)

  • t_end (float) – End time (s)

  • dt (float) – Output time step (s)

Returns:

‘t’, ‘y’, ‘gamma_dot’, ‘mu_xy’, ‘velocity’

Return type:

dict

detect_banding(result, threshold=0.1)[source]

Detect shear banding from steady-state profiles.

Parameters:
  • result (dict) – Result from simulate_steady_shear()

  • threshold (float) – Relative variation threshold for banding detection

Returns:

‘is_banding’: bool ‘band_contrast’: max/min shear rate ratio ‘band_width’: approximate band width (m) ‘band_location’: center of high-shear band (m)

Return type:

dict

get_velocity_profile(result)[source]

Compute velocity profile from final shear rate profile.

v(y) = integral_0^y gamma_dot(y’) dy’

Parameters:

result (dict) – Result from simulate_steady_shear()

Returns:

Velocity profile v(y)

Return type:

ndarray

plot_profiles(result, ax=None)[source]

Plot spatial profiles (shear rate and mu_xy).

Parameters:
  • result (dict) – Result from simulate_steady_shear()

  • ax (matplotlib axes, optional) – If None, creates new figure

Return type:

matplotlib figure