TNT FENE-P (Finite Extensibility) — Handbook¶

Quick Reference ¶

TNT FENE-P at a Glance¶
Use when	Networks where chains have finite extensibility (strain stiffening at large deformations, rubber-like materials, elastomers with physical crosslinks)
Parameters	4 parameters: \(G\) (modulus, Pa), \(\tau_b\) (bond lifetime, s), \(L_{\text{max}}\) (maximum extensibility, dimensionless), \(\eta_s\) (solvent viscosity, Pa·s)
Key equation	FENE-P stress: \(\boldsymbol{\sigma} = G \, f(\mathbf{S}) \, (\mathbf{S} - \mathbf{I}) + 2\eta_s \mathbf{D}\) where \(f(\mathbf{S}) = L_{\text{max}}^2 / (L_{\text{max}}^2 - \text{tr}(\mathbf{S}))\)
Test modes	All 6 protocols: FLOW_CURVE, OSCILLATION, STARTUP, CREEP, RELAXATION, LAOS
Material examples	Rubber networks, highly extensible hydrogels, FENE polymers, silicone elastomers, polyacrylamide hydrogels, PDMS networks
Key characteristics	Strain stiffening, bounded chain extension, stress saturation at high strain rates, non-Gaussian chain statistics

When to choose FENE-P:

Strain stiffening observed in startup or LAOS experiments
Materials with finite maximum chain extensibility
Rubber-like networks with physical or chemical crosslinks
Need to prevent unphysical infinite extension of chains

When NOT to use:

Small-strain linear viscoelasticity only (use base Tanaka-Edwards)
Rigid or semi-flexible polymers (assumptions of Gaussian chains break down)
Materials showing only strain softening (consider non-affine or Bell variants)

Notation Guide ¶

Mathematical Symbols¶
Symbol	Type	Meaning
\(\mathbf{S}\)	Tensor	Conformation tensor (dimensionless, measures chain end-to-end vector)
\(G\)	Scalar	Network modulus (Pa)
\(\tau_b\)	Scalar	Bond lifetime / relaxation time (s)
\(L_{\text{max}}\)	Scalar	Maximum chain extensibility (dimensionless, \(L_{\text{max}}^2 = N_K\))
\(f(\mathbf{S})\)	Scalar	FENE-P spring function (dimensionless)
\(\text{tr}(\mathbf{S})\)	Scalar	Trace of conformation tensor \(S_{xx} + S_{yy} + S_{zz}\)
\(\eta_s\)	Scalar	Solvent viscosity (Pa·s)
\(R_{\text{max}}\)	Scalar	Maximum chain end-to-end distance (m)
\(R_0\)	Scalar	Equilibrium chain end-to-end distance (m)
\(N_K\)	Scalar	Number of Kuhn segments per chain (dimensionless)
\(\boldsymbol{\kappa}\)	Tensor	Velocity gradient tensor \(\nabla \mathbf{v}\)
\(\mathbf{D}\)	Tensor	Rate-of-strain tensor \((\boldsymbol{\kappa} + \boldsymbol{\kappa}^T)/2\)
\(\boldsymbol{\sigma}\)	Tensor	Cauchy stress tensor (Pa)
\(\gamma_0\)	Scalar	Strain amplitude (dimensionless)
\(\omega\)	Scalar	Angular frequency (rad/s)
\(\text{Wi}\)	Scalar	Weissenberg number \(\text{Wi} = \tau_b \dot{\gamma}\)

The TNT FENE-P model extends the classical Tanaka-Edwards transient network theory by incorporating finite chain extensibility through the FENE-P (Finitely Extensible Nonlinear Elastic with Peterlin closure) spring law. Real polymer chains have a finite contour length and cannot extend indefinitely — as the chain approaches its maximum extension, the restoring force diverges, preventing unphysical behavior.

Historical development:

Warner (1972): Introduced the FENE spring model to describe finite extensibility in polymer kinetic theory.
Peterlin (1966): Developed the Peterlin closure approximation for treating non-Gaussian chain statistics.
Bird, Dotson, Johnson (1980): Formulated the FENE-P constitutive equation for dilute polymer solutions.
Tanaka & Edwards (1992): Transient network theory for physical gels and associating polymers.
Integration: FENE-P closure combined with transient network framework to model extensible networks.

Key Distinctions from Base Model ¶

The base Tanaka-Edwards model assumes Hookean (Gaussian) springs, valid only at small chain extensions. The FENE-P variant introduces:

Nonlinear spring function \(f(\mathbf{S})\) that diverges as \(\text{tr}(\mathbf{S}) \to L_{\text{max}}^2\)
Strain stiffening at large deformations (stress increases faster than linear)
Bounded chain extension preventing unphysical infinite stretch
Stress saturation at high strain rates due to maximum chain extension

Comparison:

Gaussian vs FENE-P Spring Laws¶
Property	Gaussian (Base)	FENE-P
Spring function	\(f = 1\)	\(f = L_{\text{max}}^2 / (L_{\text{max}}^2 - \text{tr}(\mathbf{S}))\)
Small extension	Linear	Linear (recovers Gaussian)
Large extension	Unbounded	Diverges (stiffens)
Stress	\(\sigma \propto S\)	\(\sigma \propto f(S) \cdot S\)
Material	Weak gels	Rubber networks, elastomers

Physical Foundations ¶

Gaussian vs Non-Gaussian Chain Statistics ¶

A freely-jointed chain with \(N_K\) Kuhn segments of length \(b_K\) has:

Equilibrium end-to-end distance: \(R_0 = b_K \sqrt{N_K}\)
Maximum contour length: \(R_{\text{max}} = N_K b_K\)
Extensibility parameter: \(L_{\text{max}}^2 = R_{\text{max}}^2 / R_0^2 = N_K\)

At small extensions \(R \ll R_{\text{max}}\), the chain behaves as a Hookean spring (Gaussian statistics):

\[F \approx \frac{3 k_B T}{R_0^2} R\]

As \(R \to R_{\text{max}}\), the exact force (from inverse Langevin function \(\mathcal{L}^{-1}\)) diverges:

\[F = \frac{k_B T}{b_K} \mathcal{L}^{-1}\left(\frac{R}{R_{\text{max}}}\right)\]

The FENE-P approximation replaces this with a tractable form:

\[F \approx \frac{3 k_B T}{R_0^2} \frac{L_{\text{max}}^2}{L_{\text{max}}^2 - R^2/R_0^2} R\]

This is the Peterlin closure — it captures the divergence at maximum extension while remaining analytically tractable.

Interpretation of \(L_{\text{max}}\)¶

The dimensionless extensibility parameter \(L_{\text{max}}\) has direct physical meaning:

\[L_{\text{max}}^2 = \frac{R_{\text{max}}^2}{R_0^2} = N_K\]

where \(N_K\) is the number of Kuhn segments in the chain.

Examples:

Short chains (\(N_K \sim 10\), \(L_{\text{max}} \sim 3\)): Strongly non-Gaussian, early strain stiffening
Moderate flexibility (\(N_K \sim 100\), \(L_{\text{max}} \sim 10\)): Typical for flexible polymers
Long chains (\(N_K \sim 1000\), \(L_{\text{max}} \sim 30\)): Nearly Gaussian until very high strains

From molecular weight:

\[L_{\text{max}} \approx \sqrt{\frac{M_w}{M_K}}\]

where \(M_w\) is the chain molecular weight and \(M_K\) is the Kuhn monomer mass.

Mechanical Analog ¶

The FENE-P spring can be visualized as a nonlinear spring with hardening:

At rest, the spring constant is \(k_0 = 3 k_B T / R_0^2\)
As extension \(R\) increases, the effective spring constant increases: \(k_{\text{eff}} = k_0 f(R)\)
At maximum extension, \(k_{\text{eff}} \to \infty\) (infinite stiffness prevents further extension)

This is analogous to a “rubber band” that becomes progressively harder to stretch.

Representative Materials ¶

Materials Exhibiting FENE-P Behavior¶
Material Class	Examples
Natural rubber	Vulcanized rubber, latex networks
Silicone elastomers	PDMS networks, silicone gels
Hydrogels	Polyacrylamide (PAAm), alginate gels with physical crosslinks
Thermoplastic elastomers	Styrene-butadiene-styrene (SBS), triblock copolymers
Associating polymers	PEO with hydrophobic associating groups, telechelic polymers

Governing Equations ¶

FENE-P Spring Function ¶

The central feature of the FENE-P model is the nonlinear spring function:

\[f(\mathbf{S}) = \frac{L_{\text{max}}^2}{L_{\text{max}}^2 - \text{tr}(\mathbf{S})}\]

where \(\text{tr}(\mathbf{S}) = S_{xx} + S_{yy} + S_{zz}\) is the trace of the conformation tensor.

Properties:

\(f(\mathbf{I}) = L_{\text{max}}^2 / (L_{\text{max}}^2 - 3)\) at equilibrium (isotropic state \(\mathbf{S} = \mathbf{I}\))
\(f(\mathbf{S}) > 1\) for \(\text{tr}(\mathbf{S}) > 3\) (extension)
\(f(\mathbf{S}) \to \infty\) as \(\text{tr}(\mathbf{S}) \to L_{\text{max}}^2\) (singularity at maximum extension)
\(f(\mathbf{S}) \to 1\) as \(L_{\text{max}} \to \infty\) (recovers Gaussian limit)

Conformation Tensor Evolution ¶

The conformation tensor \(\mathbf{S}\) evolves according to:

\[\frac{D\mathbf{S}}{Dt} = \boldsymbol{\kappa} \cdot \mathbf{S} + \mathbf{S} \cdot \boldsymbol{\kappa}^T - \frac{1}{\tau_b} \left( f(\mathbf{S}) \mathbf{S} - \mathbf{I} \right)\]

where:

\(D/Dt\) is the material derivative
\(\boldsymbol{\kappa} = \nabla \mathbf{v}\) is the velocity gradient tensor
The term \(\boldsymbol{\kappa} \cdot \mathbf{S} + \mathbf{S} \cdot \boldsymbol{\kappa}^T\) represents affine deformation
The term \(-(1/\tau_b)(f(\mathbf{S})\mathbf{S} - \mathbf{I})\) represents bond breakage and reformation

Key difference from Gaussian model: The breakage term is multiplied by \(f(\mathbf{S})\), making the effective relaxation time strain-dependent:

\[\tau_{\text{eff}}(\mathbf{S}) = \frac{\tau_b}{f(\mathbf{S})}\]

At high extension, \(f(\mathbf{S})\) increases, so \(\tau_{\text{eff}}\) decreases — chains relax faster when highly stretched.

Stress Tensor ¶

The Cauchy stress tensor is:

\[\boldsymbol{\sigma} = G \, f(\mathbf{S}) \, (\mathbf{S} - \mathbf{I}) + 2\eta_s \mathbf{D} - p \mathbf{I}\]

where:

\(G\) is the network modulus (Pa)
\(f(\mathbf{S})(\mathbf{S} - \mathbf{I})\) is the elastic stress from chain extension
\(2\eta_s \mathbf{D}\) is the solvent viscosity contribution
\(p\) is the hydrostatic pressure (determined by incompressibility)

For simple shear, the shear stress is:

\[\sigma_{xy} = G \, f(\mathbf{S}) \, S_{xy} + \eta_s \dot{\gamma}\]

For oscillatory shear, only the elastic contribution matters (assuming \(\eta_s = 0\)):

\[\sigma_{xy}(t) = G \, f(\mathbf{S}(t)) \, S_{xy}(t)\]

Steady Shear Flow ¶

In steady simple shear at rate \(\dot{\gamma}\), the components of \(\mathbf{S}\) satisfy:

\[\begin{split}\dot{\gamma} S_{yy} - \frac{1}{\tau_b} (f S_{xy} - 0) &= 0 \\ \dot{\gamma} (S_{xx} + S_{yy}) - \frac{1}{\tau_b} (f S_{xx} - 1) &= 0 \\ -\frac{1}{\tau_b} (f S_{yy} - 1) &= 0 \\ -\frac{1}{\tau_b} (f S_{zz} - 1) &= 0\end{split}\]

where \(f = L_{\text{max}}^2 / (L_{\text{max}}^2 - (S_{xx} + S_{yy} + S_{zz}))\).

Unlike the Gaussian case, these equations are implicit and must be solved numerically (no closed-form solution exists).

Limiting behavior:

Low Wi (\(\text{Wi} = \tau_b \dot{\gamma} \ll 1\)): \(f \approx f_{\text{eq}} = L_{\text{max}}^2/(L_{\text{max}}^2 - 3)\), recovers Newtonian \(\eta \approx G \tau_b f_{\text{eq}}\)
High Wi (\(\text{Wi} \gg 1\)): \(S_{xx}\) grows, \(f\) increases, stress saturates near \(G L_{\text{max}}^2\)

Small Amplitude Oscillatory Shear (SAOS)¶

For infinitesimal strain \(\gamma_0 \to 0\), linearize around equilibrium \(\mathbf{S} = \mathbf{I}\):

\[f_{\text{eq}} = \frac{L_{\text{max}}^2}{L_{\text{max}}^2 - 3}\]

The linearized ODE gives:

\[\begin{split}G'(\omega) &= G \, f_{\text{eq}} \, \frac{(\omega \tau_b)^2}{1 + (\omega \tau_b)^2} \\ G''(\omega) &= G \, f_{\text{eq}} \, \frac{\omega \tau_b}{1 + (\omega \tau_b)^2}\end{split}\]

Effect of finite extensibility:

Moduli are scaled by \(f_{\text{eq}} > 1\)
For \(L_{\text{max}} = 10\): \(f_{\text{eq}} = 100/97 \approx 1.03\) (3% increase)
For \(L_{\text{max}} = 3\): \(f_{\text{eq}} = 9/6 = 1.5\) (50% increase)
Relaxation time remains \(\tau_b\) (to leading order)

Startup of Shear Flow ¶

Following a step imposition of shear rate \(\dot{\gamma}\) at \(t = 0\), the stress evolves as:

\[\sigma_{xy}(t) = G \, f(\mathbf{S}(t)) \, S_{xy}(t) + \eta_s \dot{\gamma}\]

where \(\mathbf{S}(t)\) is obtained by numerically integrating the ODE.

Characteristic features:

Initial linear growth: \(\sigma \approx G f_{\text{eq}} \dot{\gamma} t\) for \(t \ll \tau_b\)
Strain stiffening: As \(\text{tr}(\mathbf{S})\) increases, \(f\) increases, leading to an upturn in stress
Overshoot: For \(\text{Wi} > 1\), stress peaks at \(t \sim \tau_b\), then relaxes
Steady state: Approaches steady-shear solution with \(f > f_{\text{eq}}\)

The FENE-P overshoot can be more pronounced than Gaussian due to stiffening before the peak.

Stress Relaxation ¶

After cessation of flow from steady state at \(\text{Wi} = \text{Wi}_0\), the stress relaxes:

\[\sigma_{xy}(t) = \sigma_0 \exp\left(-\int_0^t \frac{f(\mathbf{S}(s))}{\tau_b} ds\right)\]

where \(\mathbf{S}(t)\) evolves under \(\boldsymbol{\kappa} = 0\):

\[\frac{d\mathbf{S}}{dt} = -\frac{1}{\tau_b} (f(\mathbf{S}) \mathbf{S} - \mathbf{I})\]

Key feature: The relaxation is faster than exponential initially (while \(f > f_{\text{eq}}\)), then slows down as \(\mathbf{S} \to \mathbf{I}\).

Effective relaxation time:

\[\tau_{\text{eff}}(t) = \frac{\tau_b}{f(\mathbf{S}(t))}\]

FENE Damping Function ¶

The step-strain damping function for the FENE-P TNT model is:

\[h(\gamma_0) \approx \frac{L_{\max}^2 - 3}{L_{\max}^2 - (2 + \gamma_0^2)}\]

For small strains (\(\gamma_0 \ll L_{\max}\)), \(h \to 1\) (linear regime). As \(\gamma_0 \to \sqrt{L_{\max}^2 - 2}\), the damping function diverges — reflecting the stress stiffening as chains approach maximum extension.

Initial FENE factor after step strain:

\[f_0 = \frac{L_{\max}^2}{L_{\max}^2 - (2 + \gamma_0^2)}\]

This produces a very fast initial relaxation for large \(\gamma_0\) approaching \(L_{\max}\), because the high initial FENE force drives rapid stress decay.

Creep and Recovery ¶

Under constant stress \(\sigma_0\), the system evolves via 5 coupled ODEs (for \(S_{xx}, S_{yy}, S_{zz}, S_{xy}\), plus \(\gamma\) from kinematics).

Elastic contribution (instantaneous jump):

\[\gamma_e(0^+) = \frac{\sigma_0}{G f_{\text{eq}}}\]

Viscous flow: As \(\mathbf{S}\) evolves, the strain rate \(\dot{\gamma}(t) = \sigma_0 / \eta(t)\) depends on the instantaneous viscosity:

\[\eta(t) = G \, f(\mathbf{S}(t)) \, \tau_b\]

Unlike the Gaussian case (constant \(\eta = G\tau_b\)), the FENE-P viscosity increases with extension.

Creep Compliance Saturation ¶

The FENE stress function prevents infinite chain extension, which manifests in creep as compliance saturation: the creep compliance \(J(t)\) levels off sharply as chains approach \(L_{\max}\). Unlike the Hookean model where \(\gamma \to \infty\) under constant stress, FENE chains have a maximum recoverable strain set by \(L_{\max}\).

At high applied stress, the creep rate \(\dot{\gamma}(t)\) shows a transient acceleration followed by saturation — the FENE spring hardens and resists further extension.

Large Amplitude Oscillatory Shear (LAOS)¶

For \(\gamma(t) = \gamma_0 \sin(\omega t)\) with \(\gamma_0 \gg 1\), the FENE-P nonlinearity generates rich harmonic content:

\[\sigma(t) = \sum_{n=1,3,5,\ldots} [G'_n(\gamma_0, \omega) \sin(n\omega t) + G''_n(\gamma_0, \omega) \cos(n\omega t)]\]

Sources of nonlinearity:

\(f(\mathbf{S})\) varies cyclically with \(\mathbf{S}(t)\)
Asymmetry between extension and compression (chains stiffen more in extension)

Lissajous curves (\(\sigma\) vs \(\gamma\)) show characteristic stiffening loops at large \(\gamma_0\).

Parameter Table ¶

TNT FENE-P Parameters¶
Symbol	Default	Bounds	Units	Description
\(G\)	1000	\((10^0, 10^8)\)	Pa	Network modulus (plateau modulus)
\(\tau_b\)	1.0	\((10^{-6}, 10^4)\)	s	Bond lifetime / relaxation time
\(L_{\text{max}}\)	10.0	\((2.0, 100)\)	dimensionless	Maximum chain extensibility (\(\sqrt{N_K}\))
\(\eta_s\)	0.0	\((0, 10^4)\)	Pa·s	Solvent viscosity (optional)

Notes:

Lower bound \(L_{\text{max}} \ge 2\) ensures FENE-P closure validity (minimum 4 Kuhn segments)
For \(L_{\text{max}} > 50\), behavior is nearly Gaussian (little practical difference)
\(\eta_s\) is often negligible for dense networks, but important for dilute solutions

Parameter Interpretation ¶

Effect of \(L_{\text{max}}\)¶

Extensibility Regimes¶
\(L_{\text{max}}\)	Regime	Behavior
\(\to \infty\)	Gaussian limit	Recovers base Tanaka-Edwards model, no strain stiffening
\(50 - 100\)	Nearly Gaussian	Stiffening only at very high strains (\(\gamma > 10\))
\(10 - 30\)	Moderate extensibility	Typical for flexible polymers, stiffening at \(\gamma \sim 1-5\)
\(3 - 10\)	Strongly non-Gaussian	Early onset of stiffening, significant FENE effects
\(2 - 3\)	Highly extensible limit	Dramatic stiffening, low chain flexibility

Connection to Chain Structure ¶

From polymer physics:

\[L_{\text{max}}^2 = N_K = \frac{M_w}{M_K}\]

where:

\(N_K\) = number of Kuhn segments per chain
\(M_w\) = chain molecular weight (g/mol)
\(M_K\) = Kuhn monomer mass (g/mol)

Example (polyacrylamide):

Kuhn length \(b_K \approx 2\) nm
Molecular weight \(M_w = 10^5\) g/mol
Monomer mass \(M_0 = 71\) g/mol
Kuhn monomer \(M_K \approx 2 \times M_0 = 142\) g/mol
\(N_K \approx 10^5 / 142 \approx 700\)
\(L_{\text{max}} \approx \sqrt{700} \approx 26\)

Relationship Between Parameters ¶

Plateau modulus and network density:

\[G = \nu k_B T = \frac{\rho R T}{M_c}\]

where \(\nu\) is the network strand density (chains/m³) and \(M_c\) is the molecular weight between crosslinks.

Zero-shear viscosity (Gaussian limit):

\[\eta_0 = G \tau_b f_{\text{eq}} + \eta_s \approx G \tau_b \left(1 + \frac{3}{L_{\text{max}}^2 - 3}\right)\]

For large \(L_{\text{max}}\), \(f_{\text{eq}} \to 1\), giving \(\eta_0 \approx G\tau_b\).

Strain at onset of stiffening:

Heuristic estimate from \(f(\mathbf{S}) \approx 1.5 f_{\text{eq}}\):

\[\gamma_{\text{stiff}} \sim \frac{L_{\text{max}}}{\sqrt{3}}\]

For \(L_{\text{max}} = 10\), \(\gamma_{\text{stiff}} \sim 6\).

Dimensionless Groups ¶

Key Dimensionless Parameters¶
Group	Definition	Physical Meaning
Weissenberg number	\(\text{Wi} = \tau_b \dot{\gamma}\)	Ratio of relaxation to deformation timescales
Deborah number	\(\text{De} = \tau_b \omega\)	Ratio of relaxation to observation timescales
Extensibility	\(L_{\text{max}}^2 = N_K\)	Number of Kuhn segments
Normalized trace	\(\text{tr}(\mathbf{S}) / L_{\text{max}}^2\)	Fraction of maximum extension (must be < 1)

Validity and Assumptions ¶

Peterlin Closure Approximation ¶

The FENE-P model uses the Peterlin approximation:

\[\langle \mathbf{R} \otimes \mathbf{R} \rangle \approx \frac{\langle \mathbf{R} \rangle \otimes \langle \mathbf{R} \rangle}{\langle R^2 \rangle / R_{\text{max}}^2}\]

This decouples the orientation and extension fluctuations. Limitations:

Underestimates stress at very high extensions (compared to exact FENE or bead-spring simulations)
Error increases as \(L_{\text{max}}\) decreases
For \(L_{\text{max}} \ge 10\), approximation is typically accurate to within 10-20%

Affine Deformation Assumption ¶

The model assumes affine deformation: the conformation tensor \(\mathbf{S}\) deforms with the macroscopic flow field.

Breaks down when:

Network heterogeneity (non-affine rearrangements)
Slip at crosslinks (see non-affine variant)
Highly entangled systems

Mean-Field Averaging ¶

All chains are assumed to have the same conformation tensor \(\mathbf{S}\). In reality, there is a distribution of chain extensions.

Consequences:

Overestimates the rate of stress increase in startup (real networks have dispersion)
Cannot capture chain-length polydispersity effects

Valid for Flexible Polymers ¶

The FENE-P model assumes flexible, Gaussian-like chains with many Kuhn segments.

Not valid for:

Rigid rods (use wormlike chain models)
Semi-flexible polymers with persistence length \(\ell_p \sim L\) (use WLC-based models)
Colloidal particles (use Brownian dynamics)

Constant Breakage Rate ¶

Bond breakage rate \(1/\tau_b\) is assumed independent of force.

To incorporate force-dependent kinetics, combine with the Bell variant:

\[\tau_b \to \tau_b(\mathbf{S}) = \tau_{b0} \exp\left(\frac{F(\mathbf{S})}{F_b}\right)\]

where \(F(\mathbf{S}) \propto f(\mathbf{S}) \sqrt{\text{tr}(\mathbf{S})}\).

Regimes and Behavior ¶

Small Strain Regime ¶

Condition: \(\text{tr}(\mathbf{S}) - 3 \ll L_{\text{max}}^2\)

In this regime:

\[f(\mathbf{S}) \approx f_{\text{eq}} = \frac{L_{\text{max}}^2}{L_{\text{max}}^2 - 3}\]

The stress response is nearly Gaussian:

\[\sigma_{xy} \approx G f_{\text{eq}} S_{xy}\]

Observable: SAOS moduli are \(G' \approx G f_{\text{eq}} (\omega\tau_b)^2 / (1 + (\omega\tau_b)^2)\).

Moderate Strain Regime ¶

Condition: \(\text{tr}(\mathbf{S}) \sim L_{\text{max}}^2 / 2\)

The FENE function begins to increase noticeably:

\[f(\mathbf{S}) \approx 2 f_{\text{eq}}\]

Observable: Onset of strain stiffening in startup experiments. The stress grows faster than \(t\) (upward curvature in \(\sigma(t)\)).

Near Maximum Extension Regime ¶

Condition: \(L_{\text{max}}^2 - \text{tr}(\mathbf{S}) \sim O(1)\)

The FENE function diverges:

\[f(\mathbf{S}) \sim \frac{L_{\text{max}}^2}{\epsilon}, \quad \epsilon \ll 1\]

Stress saturation:

\[\sigma_{\text{max}} \sim G L_{\text{max}}^2\]

Observable: Stress plateaus at high \(\text{Wi}\) in steady shear. This is a signature prediction of finite extensibility.

Warning: In practice, sample failure (chain scission, network rupture) often occurs before reaching this limit.

Extensional Viscosity Bound ¶

Unlike the upper-convected Maxwell (UCM) model which predicts divergent extensional viscosity at \(Wi = 1/2\), the FENE-P model gives a bounded extensional viscosity:

\[\eta_E \leq 3\eta_0 \frac{L_{\max}^2}{L_{\max}^2 - 3}\]

This bound is approached as \(\dot{\epsilon} \to \infty\). For typical values \(L_{\max} = 10\), the maximum extensional viscosity is approximately \(3.1 \eta_0\) — a physically realistic prediction.

Steady Shear Regimes ¶

Flow Curve Regimes¶
\(\text{Wi}\)	Conformation	Viscosity
\(\ll 1\)	Isotropic (\(\mathbf{S} \approx \mathbf{I}\))	Newtonian \(\eta \approx G\tau_b f_{\text{eq}}\)
\(\sim 1\)	Moderate extension	Shear thinning begins
\(\gg 1\)	High extension (\(S_{xx} \gg 1\))	Bounded \(\eta \sim G\tau_b L_{\text{max}}^2 / \text{Wi}\)

Key difference from Gaussian: The FENE-P viscosity saturates at high \(\text{Wi}\) rather than decreasing indefinitely.

Startup Regimes ¶

For \(\text{Wi} > 1\), startup shows:

Linear growth (\(t \ll \tau_b\)): \(\sigma \approx G f_{\text{eq}} \dot{\gamma} t\)
Strain stiffening (\(t \sim \tau_b\)): Stress curves upward as \(f\) increases
Overshoot (\(t \sim 2\tau_b\)): Stress peaks at \(\sigma_{\text{max}} > \sigma_{\text{ss}}\)
Relaxation to steady state (\(t \gg \tau_b\)): Stress decays to steady-shear value

FENE-P overshoot can be larger than Gaussian due to stiffening before the peak.

LAOS Regimes ¶

LAOS Behavior vs Strain Amplitude¶
\(\gamma_0\)	Response	Harmonics
\(\ll 1\)	Linear (SAOS)	Only \(n=1\)
\(\sim 1\)	Weakly nonlinear	Weak \(n=3,5\)
\(\gg 1\)	Strongly nonlinear	Strong higher harmonics, stiffening loops

Lissajous signature: Clockwise loops in \(\sigma\) vs \(\gamma\) due to strain stiffening (elastic dominance).

LAOS Strain-Stiffening Signature ¶

In LAOS, FENE stress function produces a distinctive strain-stiffening signature:

Lissajous curves become rectangular (box-like) at large \(\gamma_0\) — high stress at peak strain due to the FENE divergence as \(\text{tr}(\mathbf{S}) \to L_{\max}^2\)
Intra-cycle stiffening: The tangent modulus \(G'_L\) (large-strain modulus from Ewoldt decomposition) exceeds \(G'_M\) (minimum-strain modulus)
Higher harmonics grow with \(\gamma_0\), but the dominant signature is stiffening rather than the softening seen in Bell variants

What You Can Learn ¶

From SAOS Data ¶

Information from SAOS¶
Observable	Extracted Parameter
Plateau modulus \(G_0\)	\(G \times f_{\text{eq}}\) (need to decouple \(L_{\text{max}}\) from nonlinear data)
Crossover frequency \(\omega_c\)	Relaxation time \(\tau_b \approx 1/\omega_c\)
Slope of \(\tan\delta\) at low \(\omega\)	No direct \(L_{\text{max}}\) information (linearized model)

Limitation: \(L_{\text{max}}\) is not identifiable from SAOS alone (only appears via \(f_{\text{eq}} \approx 1\) correction).

From Startup and Flow Curves ¶

Information from Nonlinear Shear¶
Observable	Extracted Parameter
Strain at stiffening onset	\(L_{\text{max}}\) (via \(\gamma_{\text{stiff}} \sim L_{\text{max}}/\sqrt{3}\))
Stress overshoot magnitude	Combined \(G, \tau_b, L_{\text{max}}\) (stronger overshoot for low \(L_{\text{max}}\))
Steady-state stress saturation	\(G \times L_{\text{max}}^2\) (maximum network stress)
High-Wi viscosity	Bounded value confirms FENE-P vs Gaussian

From LAOS Data ¶

Information from LAOS¶
Observable	Extracted Parameter
Strain amplitude for nonlinearity	\(L_{\text{max}}\) (via \(\gamma_0 \sim L_{\text{max}}\))
Lissajous loop shape	Stiffening vs softening (FENE-P gives stiffening)
Higher harmonics \(I_{3/1}, I_{5/1}\)	Degree of non-Gaussian behavior (larger for low \(L_{\text{max}}\))

Chain Architecture from \(L_{\text{max}}\)¶

Once \(L_{\text{max}}\) is known, infer chain structure:

\[N_K = L_{\text{max}}^2\]

If Kuhn length \(b_K\) is known (from literature or scattering):

\[R_{\text{max}} = L_{\text{max}} \times b_K \sqrt{N_K} = L_{\text{max}}^2 b_K\]

If molecular weight is known:

\[M_K \approx \frac{M_w}{L_{\text{max}}^2}\]

Example: \(L_{\text{max}} = 10\) from startup → \(N_K = 100\) Kuhn segments.

Stress Saturation Level ¶

The maximum stress the network can sustain:

\[\sigma_{\text{max}} \approx G L_{\text{max}}^2\]

Practical use:

Predict failure stress (if chains break before full extension)
Estimate safe operating strain limits
Design formulations to avoid rupture

Experimental Design ¶

Best Protocols for Determining \(L_{\text{max}}\)¶

Protocol Recommendations¶
Protocol	\(\text{Wi}\) or \(\gamma_0\)	Observable
Startup	\(\text{Wi} > 5\)	Strain stiffening upturn before overshoot
LAOS	\(\gamma_0 > 1\)	Stiffening loops, higher harmonics
Flow curve	\(\text{Wi} > 10\)	Stress saturation at high \(\dot{\gamma}\)

Optimal strategy: Combine SAOS (for \(G, \tau_b\)) + Startup (for \(L_{\text{max}}\)) + LAOS (for validation).

SAOS First ¶

Purpose: Obtain \(G\) and \(\tau_b\) from linear response.

Protocol:

Frequency sweep: \(\omega = 10^{-2}\) to \(10^2\) rad/s
Strain amplitude: \(\gamma_0 = 0.01\) (linear regime)
Extract \(G_0 \approx G f_{\text{eq}}\) and \(\tau_b\) from \(\omega_c\)

Startup at Multiple Wi ¶

Purpose: Capture strain stiffening and determine \(L_{\text{max}}\).

Protocol:

Step shear rates: \(\dot{\gamma} = 0.1, 1, 10, 100\) s⁻¹ (adjust to span \(\text{Wi} = 0.1\) to \(100\))
Duration: \(t = 10 \tau_b\) (long enough for steady state)
Look for upward curvature before overshoot

Analysis:

Fit early-time slope to get \(G f_{\text{eq}}\)
Fit stiffening region to extract \(L_{\text{max}}\)

LAOS for Validation ¶

Purpose: Confirm FENE-P nonlinearity, check for additional physics.

Protocol:

Strain amplitudes: \(\gamma_0 = 0.1, 0.5, 1, 5, 10\)
Frequency: \(\omega = 1/\tau_b\) (Deborah number \(\sim 1\))
Extract Lissajous curves and higher harmonics

Check:

Stiffening (clockwise) loops confirm FENE-P
Softening (counterclockwise) loops suggest additional physics (non-affine, yielding)

Experimental Artifacts ¶

Common Artifacts¶
Artifact	Mitigation
Edge fracture	Use cone-plate or parallel-plate with edge trim
Sample expulsion	Lower strain amplitude, use sandblasted plates
Chain scission	Limit maximum \(\gamma\) or \(\dot{\gamma}\), check for hysteresis
Wall slip	Use roughened surfaces, check with gap-dependent tests
Inertia	Keep \(\rho \omega R^2 / \eta \ll 1\) (low Reynolds number)

Complementary Techniques ¶

Complementary Measurements¶
Technique	Information
SANS/SAXS	Direct measurement of \(R_0, R_{\text{max}}\); validation of \(L_{\text{max}}\)
Birefringence	Chain orientation, validates \(\mathbf{S}\) predictions
Microscopy	Network heterogeneity, defects
DLS	Relaxation time \(\tau_b\), diffusion coefficient

Computational Implementation ¶

JAX Implementation ¶

The FENE-P stress kernel is implemented in rheojax/models/tnt/_kernels.py as:

def stress_fene(S: Array, params: Array) -> Array:
    """FENE-P stress function.

    Args:
        S: Conformation tensor (3x3)
        params: [G, L_max]

    Returns:
        Stress tensor (3x3)
    """
    G, L_max = params
    tr_S = jnp.trace(S)

    # Clamp trace to prevent singularity
    epsilon = 1e-6
    tr_S_safe = jnp.minimum(tr_S, L_max**2 - epsilon)

    f = L_max**2 / (L_max**2 - tr_S_safe)
    sigma = G * f * (S - jnp.eye(3))
    return sigma

FENE Singularity Handling ¶

The key numerical challenge is the singularity \(f \to \infty\) as \(\text{tr}(\mathbf{S}) \to L_{\text{max}}^2\).

Clamping strategy:

epsilon = 1e-6  # Safety margin
tr_S_safe = jnp.minimum(tr_S, L_max**2 - epsilon)

This prevents division by zero while allowing \(f\) to become very large (up to \(\sim L_{\text{max}}^2 / \epsilon\)).

Alternative: Use a smoothed FENE function:

\[f_{\text{smooth}}(\mathbf{S}) = \frac{L_{\text{max}}^2}{L_{\text{max}}^2 - \text{tr}(\mathbf{S}) + \epsilon}\]

This avoids sharp cutoffs but changes the physics slightly at high extension.

ODE Stiffness ¶

Near the FENE singularity, the ODE becomes stiff due to rapid changes in \(f(\mathbf{S})\).

Solution: Use adaptive step-size control in Diffrax:

from diffrax import diffeqsolve, Dopri5, PIDController

solver = Dopri5()
stepsize_controller = PIDController(rtol=1e-6, atol=1e-8)

solution = diffeqsolve(
    terms, solver, t0, t1, dt0, y0,
    stepsize_controller=stepsize_controller
)

The PID controller automatically reduces step size when \(f\) increases rapidly.

JIT Compilation ¶

The stress kernel is JIT-compiled for performance:

from jax import jit

@jit
def residuals_fene(params, x, y, test_mode):
    # ... compute predictions
    return y_pred - y

Dispatch at Python level: The stress type (“gaussian”, “fene”, “bell”) is selected before JIT compilation, avoiding dynamic branches inside the JIT region.

Memory Considerations ¶

For time-series predictions (startup, LAOS), the full trajectory \(\mathbf{S}(t)\) is stored:

Memory usage: \(O(N_{\text{time}} \times 9)\) (9 components of \(\mathbf{S}\))
Optimization: Use jax.lax.scan instead of jax.lax.fori_loop to reduce memory

Example:

def step(S, t):
    S_new = integrate_one_step(S, t, dt)
    return S_new, S_new

_, S_trajectory = jax.lax.scan(step, S_initial, t_array)

This stores only the final result in memory during compilation.

Fitting Guidance ¶

Two-Stage Fitting Strategy ¶

Recommended workflow:

Stage 1: SAOS → Fit \(G, \tau_b\) using linearized model (fix \(L_{\text{max}} = \infty\) or large value)
Stage 2: Nonlinear → Fit startup/LAOS data with \(G, \tau_b\) fixed, optimize \(L_{\text{max}}\)

Rationale: \(L_{\text{max}}\) is weakly constrained by SAOS (only via \(f_{\text{eq}} \approx 1\) correction), but strongly constrained by large-strain data.

# Stage 1: SAOS
model_saos = TNTSingleMode(stress_type="gaussian")  # Gaussian approximation
model_saos.fit(omega, G_star, test_mode='oscillation')
G_fit, tau_b_fit = model_saos.params['G'].value, model_saos.params['tau_b'].value

# Stage 2: Startup
model_fene = TNTSingleMode(stress_type="fene")
model_fene.params['G'].value = G_fit
model_fene.params['tau_b'].value = tau_b_fit
model_fene.params['G'].fixed = True
model_fene.params['tau_b'].fixed = True
model_fene.fit(t, sigma_startup, test_mode='startup', gamma_dot=10.0)
L_max_fit = model_fene.params['L_max'].value

Bayesian Inference ¶

Recommended priors:

import numpyro.distributions as dist

priors = {
    'G': dist.LogNormal(jnp.log(1e3), 1.0),
    'tau_b': dist.LogNormal(jnp.log(1.0), 1.0),
    'L_max': dist.LogNormal(jnp.log(10.0), 0.5),  # Mode at 10
    'eta_s': dist.HalfNormal(1e2)
}

Why LogNormal for :math:`L_{text{max}}`:

Ensures \(L_{\text{max}} > 0\)
Allows wide range (e.g., 5-20) while concentrating prior mass near typical values
Reflects multiplicative uncertainty (factor of 2-3 uncertainty is common)

Posterior checks:

Verify \(L_{\text{max}}\) posterior is unimodal (multimodality indicates non-identifiability)
Check \(\hat{R} < 1.01\) for convergence
Ensure ESS > 400 for reliable credible intervals

Identifiability of \(L_{\text{max}}\)¶

\(L_{\text{max}}\) is only identifiable from nonlinear data where \(\text{tr}(\mathbf{S})\) becomes large.

Symptoms of non-identifiability:

Flat posterior over wide range of \(L_{\text{max}}\)
High correlation between \(G\) and \(L_{\text{max}}\) (both scale stress)
Poor convergence (low ESS, high \(\hat{R}\))

Solution: Include data at high \(\text{Wi}\) or \(\gamma_0\) where FENE effects are pronounced.

NLSQ Convergence Issues ¶

Common failure modes:

Singularity divergence: Optimizer drives \(\text{tr}(\mathbf{S}) \to L_{\text{max}}^2\), causing \(f \to \infty\)
Poor initialization: Starting \(L_{\text{max}}\) too small leads to early saturation

Mitigations:

# Initialize L_max larger than observed max trace
tr_S_max = compute_max_trace_from_data(t, sigma)
L_max_init = jnp.sqrt(tr_S_max) * 1.5  # 50% safety margin

model.params['L_max'].value = L_max_init

# Use bounded optimization
from rheojax.utils.optimization import nlsq_curve_fit
result = nlsq_curve_fit(
    model_fn, x, y, params,
    bounds={'L_max': (2.0, 100.0)}  # Hard bounds
)

Multi-Start Optimization ¶

For complex data (multiple overshoot peaks, noise), use multi-start:

from rheojax.utils.optimization import nlsq_optimize

L_max_guesses = [3.0, 10.0, 30.0]
results = []
for L_max_init in L_max_guesses:
    model.params['L_max'].value = L_max_init
    result = nlsq_optimize(residuals, params)
    results.append(result)

# Select best fit
best = min(results, key=lambda r: r.cost)

Usage Examples ¶

Basic Fitting to Startup Data ¶

from rheojax.models.tnt import TNTSingleMode
from rheojax.core import RheoData
import jax.numpy as jnp

# Generate synthetic startup data
t = jnp.linspace(0, 10, 200)
gamma_dot = 10.0

# True parameters
model_true = TNTSingleMode(stress_type="fene")
model_true.params['G'].value = 1000.0
model_true.params['tau_b'].value = 1.0
model_true.params['L_max'].value = 10.0
stress_true = model_true.predict(t, test_mode='startup', gamma_dot=gamma_dot)

# Add noise
stress_noisy = stress_true + jax.random.normal(jax.random.PRNGKey(0), stress_true.shape) * 10

# Fit model
model = TNTSingleMode(stress_type="fene")
rheo_data = RheoData(x=t, y=stress_noisy, test_mode='startup', gamma_dot=gamma_dot)
model.fit(rheo_data)

print(f"Fitted G: {model.params['G'].value:.1f} Pa")
print(f"Fitted tau_b: {model.params['tau_b'].value:.3f} s")
print(f"Fitted L_max: {model.params['L_max'].value:.2f}")

# Predict and plot
stress_pred = model.predict(t, test_mode='startup', gamma_dot=gamma_dot)
import matplotlib.pyplot as plt
plt.plot(t, stress_noisy, 'o', label='Data')
plt.plot(t, stress_pred, '-', label='FENE-P fit')
plt.xlabel('Time (s)')
plt.ylabel('Stress (Pa)')
plt.legend()
plt.show()

Bayesian Inference with Warm-Start ¶

from rheojax.models.tnt import TNTSingleMode
import numpyro.distributions as dist

# Step 1: NLSQ point estimate
model = TNTSingleMode(stress_type="fene")
model.fit(t, stress, test_mode='startup', gamma_dot=10.0)

# Step 2: Bayesian with informative priors from NLSQ
G_nlsq = model.params['G'].value
tau_b_nlsq = model.params['tau_b'].value
L_max_nlsq = model.params['L_max'].value

priors = {
    'G': dist.LogNormal(jnp.log(G_nlsq), 0.2),
    'tau_b': dist.LogNormal(jnp.log(tau_b_nlsq), 0.2),
    'L_max': dist.LogNormal(jnp.log(L_max_nlsq), 0.3),
    'eta_s': dist.HalfNormal(10.0)
}

result = model.fit_bayesian(
    t, stress, test_mode='startup', gamma_dot=10.0,
    priors=priors,
    num_warmup=1000, num_samples=2000, num_chains=4
)

# Extract credible intervals
intervals = model.get_credible_intervals(result.posterior_samples, credibility=0.95)
print("95% Credible Intervals:")
for param, (lower, upper) in intervals.items():
    print(f"  {param}: [{lower:.3f}, {upper:.3f}]")

# Diagnostic plots
import arviz as az
idata = az.from_dict(result.posterior_samples)
az.plot_trace(idata)
plt.tight_layout()
plt.show()

Simulating Startup Response ¶

from rheojax.models.tnt import TNTSingleMode
import jax.numpy as jnp

model = TNTSingleMode(stress_type="fene")
model.params['G'].value = 1000.0
model.params['tau_b'].value = 1.0
model.params['L_max'].value = 10.0

# Startup at multiple Wi
gamma_dots = [0.1, 1.0, 10.0, 100.0]
t = jnp.linspace(0, 10, 500)

import matplotlib.pyplot as plt
for gamma_dot in gamma_dots:
    stress = model.predict(t, test_mode='startup', gamma_dot=gamma_dot)
    Wi = model.params['tau_b'].value * gamma_dot
    plt.plot(t, stress, label=f'Wi = {Wi:.1f}')

plt.xlabel('Time (s)')
plt.ylabel('Stress (Pa)')
plt.legend()
plt.title('FENE-P Startup: Strain Stiffening at High Wi')
plt.show()

LAOS Simulation and Harmonic Analysis ¶

from rheojax.models.tnt import TNTSingleMode
import jax.numpy as jnp
import numpy as np
from scipy.fft import fft

model = TNTSingleMode(stress_type="fene")
model.params['G'].value = 1000.0
model.params['tau_b'].value = 1.0
model.params['L_max'].value = 10.0

# LAOS parameters
gamma_0 = 2.0  # Large amplitude
omega = 1.0    # De ~ 1
n_cycles = 10
t = jnp.linspace(0, 2 * jnp.pi * n_cycles / omega, 2000)

# Simulate
stress = model.predict(t, test_mode='laos', gamma_0=gamma_0, omega=omega)

# Extract last cycle
t_cycle = t[-200:]
stress_cycle = stress[-200:]
gamma_cycle = gamma_0 * jnp.sin(omega * t_cycle)

# Fourier analysis
stress_fft = fft(np.array(stress_cycle))
n_harmonics = 5
harmonics = np.abs(stress_fft[:n_harmonics])

print("Harmonic intensities:")
for n in range(1, n_harmonics):
    ratio = harmonics[n] / harmonics[1]
    print(f"  I_{n}/I_1 = {ratio:.4f}")

# Lissajous plot
plt.figure(figsize=(10, 4))
plt.subplot(1, 2, 1)
plt.plot(gamma_cycle, stress_cycle)
plt.xlabel('Strain')
plt.ylabel('Stress (Pa)')
plt.title('Lissajous Curve (Stiffening)')

plt.subplot(1, 2, 2)
plt.stem(range(1, n_harmonics), harmonics[1:n_harmonics] / harmonics[1])
plt.xlabel('Harmonic n')
plt.ylabel('I_n / I_1')
plt.title('Harmonic Spectrum')
plt.tight_layout()
plt.show()

Comparing Gaussian vs FENE-P ¶

from rheojax.models.tnt import TNTSingleMode
import jax.numpy as jnp

# Shared parameters
G = 1000.0
tau_b = 1.0
gamma_dot = 10.0
t = jnp.linspace(0, 5, 300)

# Gaussian model
model_gaussian = TNTSingleMode(stress_type="gaussian")
model_gaussian.params['G'].value = G
model_gaussian.params['tau_b'].value = tau_b
stress_gaussian = model_gaussian.predict(t, test_mode='startup', gamma_dot=gamma_dot)

# FENE-P model
model_fene = TNTSingleMode(stress_type="fene")
model_fene.params['G'].value = G
model_fene.params['tau_b'].value = tau_b
model_fene.params['L_max'].value = 10.0
stress_fene = model_fene.predict(t, test_mode='startup', gamma_dot=gamma_dot)

# Plot comparison
import matplotlib.pyplot as plt
plt.plot(t, stress_gaussian, '--', label='Gaussian (base)')
plt.plot(t, stress_fene, '-', label='FENE-P (L_max=10)')
plt.xlabel('Time (s)')
plt.ylabel('Stress (Pa)')
plt.legend()
plt.title('Strain Stiffening in FENE-P vs Gaussian')
plt.show()

Composition with Other Variants ¶

FENE + Bell (Force-Dependent Breakage)¶

Combine finite extensibility with force-dependent bond kinetics:

from rheojax.models.tnt import TNTSingleMode

model = TNTSingleMode(breakage="bell", stress_type="fene")
# Parameters: G, tau_b0, L_max, F_b, eta_s

Physical scenario: Rubber networks where bonds break faster under high tension (e.g., physical crosslinks with thermal activation).

Effect: Stress saturation occurs at lower levels than pure FENE-P (bonds break before reaching maximum extension).

FENE + Non-Affine (Finite Extension with Slip)¶

model = TNTSingleMode(non_affine=True, stress_type="fene")
# Parameters: G, tau_b, L_max, xi, eta_s

Physical scenario: Elastomers with imperfect crosslinks, allowing partial slip while still exhibiting strain stiffening.

Effect: Reduces stress overshoot magnitude (slip dissipates energy) but retains stiffening at high strain.

FENE + Stretch-Dependent Creation ¶

model = TNTSingleMode(creation="stretch", stress_type="fene")
# Parameters: G, tau_b, L_max, tau_c, alpha_c, eta_s

Physical scenario: Self-healing networks where new bonds form more readily under extension.

Effect: Enhances stress hardening (reformation counteracts breakage).

Multi-Variant Composition ¶

Combine all three variants:

model = TNTSingleMode(
    breakage="bell",
    stress_type="fene",
    non_affine=True,
    creation="stretch"
)

Use case: Complex materials (e.g., dynamic covalent networks, vitrimers) requiring multiple physics.

Failure Mode: Chain Snap ¶

The FENE-P stress function diverges as \(\text{tr}(\mathbf{S}) \to L_{\max}^2\):

\[f(\mathbf{S}) = \frac{L_{\max}^2}{L_{\max}^2 - \text{tr}(\mathbf{S})} \to \infty\]

This represents catastrophic chain extension near the maximum extensibility limit. Physically, this corresponds to chains being pulled taut between crosslinks, with the entropic spring force diverging as the end-to-end distance approaches the contour length.

Practical implications:

Numerical integration requires clipping: \(\text{tr}(\mathbf{S}) < 0.99 L_{\max}^2\)
Very high \(Wi\) or \(\gamma_0\) near \(L_{\max}\) causes ODE stiffness
Physical chains would rupture or detach before reaching this limit (combine with Bell breakage for realistic behavior)

API Reference ¶

class rheojax.models.tnt.TNTSingleMode(breakage='constant', stress_type='linear', xi=0.0)[source]

Single-mode Transient Network Theory model.

Implements the Green-Tobolsky / Tanaka-Edwards transient network model with composable physics variants. The conformation tensor S tracks the average chain configuration between reversible crosslinks.

The constitutive equation is:

dS/dt = L·S + S·L^T + g₀·I - β(S)·S

Stress is computed from S via σ = G·f(S)·(S - I) + η_s·γ̇.

Parameters:

breakage (Literal['constant', 'bell', 'power_law', 'stretch_creation']) – Breakage rate function: - “constant”: β = 1/τ_b (Tanaka-Edwards, UCM-like) - “bell”: β = (1/τ_b)·exp(ν·(stretch-1)) (force-dependent) - “power_law”: β = (1/τ_b)·stretch^m - “stretch_creation”: β = (1/τ_b), g₀ = (1+κ·stretch)/τ_b
stress_type (Literal['linear', 'fene']) – Stress formula: - “linear”: σ = G·(S - I) (Gaussian chains) - “fene”: σ = G·f(tr(S))·(S - I) (finitely extensible)
xi (float) – Gordon-Schowalter slip parameter (0=upper-convected, 1=corotational)

parameters

Model parameters

Type:: ParameterSet

fitted_

Whether the model has been fitted

Type:: bool

Examples

Basic Tanaka-Edwards model:

>>> model = TNTSingleMode()
>>> gamma_dot = np.logspace(-2, 2, 50)
>>> sigma = model.predict(gamma_dot, test_mode='flow_curve')

Bell force-dependent breakage:

>>> model = TNTSingleMode(breakage="bell")
>>> # Now has additional parameter 'nu' (force sensitivity)

References ¶

[Warner1972]

Warner, H. R. (1972). “Kinetic Theory and Rheology of Dilute Suspensions of Finitely Extendible Dumbbells”. Industrial & Engineering Chemistry Fundamentals, 11(3), 379-387. https://doi.org/10.1021/i160043a017

[Peterlin1966]

Peterlin, A. (1966). “Hydrodynamics of Linear Macromolecules”. Pure and Applied Chemistry, 12(1-4), 563-586. https://doi.org/10.1351/pac196612010563

[Bird1980]

Bird, R. B., Dotson, P. J., & Johnson, N. L. (1980). “Polymer Solution Rheology Based on a Finitely Extensible Bead-Spring Chain Model”. Journal of Non-Newtonian Fluid Mechanics, 7(2-3), 213-235. https://doi.org/10.1016/0377-0257(80)85007-5

[TanakaEdwards1992]

Tanaka, F., & Edwards, S. F. (1992). “Viscoelastic Properties of Physically Crosslinked Networks”. Macromolecules, 25(5), 1516-1523. https://doi.org/10.1021/ma00031a024

[Keunings1997]

Keunings, R. (1997). “On the Peterlin Approximation for Finitely Extensible Dumbbells”. Journal of Non-Newtonian Fluid Mechanics, 68(1), 85-100. https://doi.org/10.1016/S0377-0257(96)01497-8

[Wiest1989]

Wiest, J. M. (1989). “A Differential Constitutive Equation for Polymer Melts”. Rheologica Acta, 28(1), 4-12. https://doi.org/10.1007/BF01354763

[DoiEdwards1986]

Doi, M., & Edwards, S. F. (1986). The Theory of Polymer Dynamics. Oxford University Press. ISBN: 978-0198519768

[Bird1987]

Bird, R. B., Armstrong, R. C., & Hassager, O. (1987). Dynamics of Polymeric Liquids, Volume 1: Fluid Mechanics (2nd ed.). Wiley. ISBN: 978-0471802457

Recommended Reading Order:

Warner (1972) — FENE spring foundation
Peterlin (1966) — Closure approximation
Bird et al. (1980) — FENE-P constitutive equation
Tanaka & Edwards (1992) — Transient network base theory
Keunings (1997) — FENE-P accuracy and limitations

Note

Practical Tips:

Always fit SAOS data first to get \(G\) and \(\tau_b\)
\(L_{\text{max}}\) requires nonlinear data (startup or LAOS)
Clamp \(\text{tr}(\mathbf{S}) < L_{\text{max}}^2 - \epsilon\) to avoid singularity
Use multi-start optimization if fits fail to converge
Validate FENE-P assumptions with complementary scattering data

Warning

Numerical Stability:

Near \(\text{tr}(\mathbf{S}) \to L_{\text{max}}^2\), ODEs become stiff
Use adaptive step size control (Diffrax PIDController)
Initialize \(L_{\text{max}} > \sqrt{\text{max observed } \text{tr}(\mathbf{S})}\)
Check for optimizer divergence (cost → ∞) indicating singularity issues