Giesekus Model — Handbook¶

Quick Reference¶

Use when: Polymer melts/solutions with shear-thinning, normal stress differences, stress overshoot
Parameters: 4 (\(\eta_p\), \(\lambda\), \(\alpha\), \(\eta_s\))
Key equation: \(\boldsymbol{\tau} + \lambda \overset{\nabla}{\boldsymbol{\tau}} + \frac{\alpha \lambda}{\eta_p} \boldsymbol{\tau} \cdot \boldsymbol{\tau} = 2 \eta_p \mathbf{D}\)
Diagnostic: \(N_2/N_1 = -\alpha/2\) (direct experimental route to \(\alpha\))
Test modes: Flow curve, oscillation, startup, relaxation, creep, LAOS
Material examples: Polymer melts, concentrated solutions, wormlike micelles

Notation Guide¶

Symbol	Meaning
\(\boldsymbol{\tau}\)	Polymer extra stress tensor (Pa)
\(\eta_p\)	Polymer viscosity (Pa·s). Zero-shear polymer contribution.
\(\lambda\)	Relaxation time (s). Characteristic stress decay time.
\(\alpha\)	Mobility factor (dimensionless, \(0 \leq \alpha \leq 0.5\)). Controls shear-thinning.
\(\eta_s\)	Solvent viscosity (Pa·s). Newtonian background contribution.
\(\eta_0\)	Zero-shear viscosity, \(\eta_0 = \eta_p + \eta_s\)
\(G\)	Elastic modulus, \(G = \eta_p / \lambda\)
\(\text{Wi}\)	Weissenberg number, \(\text{Wi} = \lambda \dot{\gamma}\)
\(\text{De}\)	Deborah number, \(\text{De} = \lambda / t_{\text{obs}}\)
\(N_1\)	First normal stress difference, \(N_1 = \tau_{xx} - \tau_{yy}\)
\(N_2\)	Second normal stress difference, \(N_2 = \tau_{yy} - \tau_{zz}\)
\(\Psi_1\)	First normal stress coefficient, \(\Psi_1 = N_1 / \dot{\gamma}^2\)
\(\Psi_2\)	Second normal stress coefficient, \(\Psi_2 = N_2 / \dot{\gamma}^2\)
\(\eta^*\)	Complex viscosity, \(\eta^* = \eta' - i\eta''\)
\(J(t)\)	Creep compliance, \(J(t) = \gamma(t) / \sigma_0\)
\(\overset{\nabla}{\boldsymbol{\tau}}\)	Upper-convected derivative (frame-invariant time derivative)
\(\mathbf{c}\)	Conformation tensor (average molecular conformation)

Overview¶

The Giesekus model (1982) is a nonlinear differential constitutive equation that extends the Upper-Convected Maxwell (UCM) model with a quadratic stress term representing anisotropic molecular mobility. It provides a physically motivated description of:

Shear-thinning viscosity: Viscosity decreases with increasing shear rate
Normal stress differences: Both \(N_1 > 0\) and \(N_2 < 0\)
Stress overshoot: Peak stress in startup flow at constant rate
Faster-than-exponential relaxation: Due to the quadratic stress term

The model is particularly valuable because it predicts both first and second normal stress differences with a fixed ratio \(N_2/N_1 = -\alpha/2\), providing a direct experimental route to determine the mobility parameter \(\alpha\).

Historical Context¶

Hanswalter Giesekus introduced this model in 1982 [1] as a “simple constitutive equation based on the concept of deformation-dependent tensorial mobility.” The key insight was that molecular mobility in polymer melts is not isotropic—molecules aligned by flow experience different friction in different directions.

The model became widely adopted because:

It uses only one additional parameter (\(\alpha\)) beyond the Maxwell model
It captures essential nonlinear features with simple mathematics
The parameter \(\alpha\) has clear physical interpretation
Predictions agree well with experimental data for many polymeric systems

Physical Foundations¶

Molecular Picture: Anisotropic Drag¶

The Giesekus model arises from considering how polymer chains experience drag in a flowing medium. When chains are stretched and aligned by flow:

Isotropic drag (UCM model):: Chains experience the same friction regardless of orientation. Result: No shear-thinning, \(N_2 = 0\)
Anisotropic drag (Giesekus model):: Aligned chains slip more easily along their backbone than perpendicular to it. Result: Shear-thinning, \(N_2 < 0\)

The mobility parameter \(\alpha\) quantifies this anisotropy:

\(\alpha = 0\): Isotropic drag; recovers UCM model
\(\alpha = 0.5\): Maximum anisotropy; strongest thinning
Typical values: 0.1–0.4 for most polymer melts and solutions

Network Interpretation¶

Alternatively, the Giesekus model can be derived from a temporary network theory where:

Polymer chains form a transient network of entanglements
Network junctions break and reform with rate dependent on local stress
Higher stress leads to faster junction breakage and lower effective viscosity

The quadratic \(\boldsymbol{\tau} \cdot \boldsymbol{\tau}\) term represents the stress-induced acceleration of network relaxation.

Stress Decomposition¶

The total Cauchy stress for an incompressible Giesekus fluid is split into solvent and polymeric contributions:

\[\boldsymbol{\sigma} = -p\mathbf{I} + \boldsymbol{\sigma}_s + \boldsymbol{\tau}\]

where:

\(\boldsymbol{\sigma}_s = 2\eta_s \mathbf{D}\) is the Newtonian solvent stress
\(\boldsymbol{\tau}\) is the polymer extra stress evolving via the Giesekus law
\(p\) is the isotropic pressure

Some texts use \(\boldsymbol{\sigma}_p\) in place of \(\boldsymbol{\tau}\) for the polymeric stress. Throughout this handbook we use \(\boldsymbol{\tau}\) to denote the polymer contribution, consistent with the rest of RheoJAX documentation.

Conformation Tensor Form¶

An alternative and often numerically preferred formulation uses the conformation tensor \(\mathbf{c}\) representing the average molecular conformation. The stress–configuration relation is:

\[\boldsymbol{\tau} = \frac{\eta_p}{\lambda}(\mathbf{c} - \mathbf{I})\]

The evolution of \(\mathbf{c}\) follows:

\[\lambda \overset{\nabla}{\mathbf{c}} + (\mathbf{c} - \mathbf{I}) + \alpha (\mathbf{c} - \mathbf{I})^2 = 0\]

This form is preferred in CFD applications because it guarantees positive-definiteness of \(\mathbf{c}\) when combined with appropriate numerical methods [14].

Material Functions¶

The Giesekus model defines the following material functions, which are measurable experimentally:

Shear viscosity (from steady shear):

\[\eta(\dot{\gamma}) = \frac{\sigma_{xy}}{\dot{\gamma}} = \frac{\tau_{xy}}{\dot{\gamma}} + \eta_s\]

Complex viscosity (from oscillatory shear):

\[\eta^*(\omega) = \eta'(\omega) - i\eta''(\omega)\]

Normal stress coefficients (from steady shear):

\[ \begin{align}\begin{aligned}\Psi_1(\dot{\gamma}) = \frac{N_1}{\dot{\gamma}^2} = \frac{\tau_{xx} - \tau_{yy}}{\dot{\gamma}^2}\\\Psi_2(\dot{\gamma}) = \frac{N_2}{\dot{\gamma}^2} = \frac{\tau_{yy} - \tau_{zz}}{\dot{\gamma}^2}\end{aligned}\end{align} \]

Crossover frequency (from SAOS):

\[\omega_c = \frac{1}{\lambda} \quad \text{where } G'(\omega_c) = G''(\omega_c) - \eta_s \omega_c\]

Governing Equations¶

Kinematics and Notation¶

The velocity field \(\mathbf{v}(\mathbf{x}, t)\) defines the velocity gradient tensor \(\nabla\mathbf{v}\) and the rate-of-deformation tensor:

\[\mathbf{D} = \frac{1}{2}\bigl(\nabla\mathbf{v} + (\nabla\mathbf{v})^T\bigr)\]

The upper-convected derivative of a tensor \(\mathbf{A}\) is the frame-invariant time derivative:

\[\overset{\nabla}{\mathbf{A}} = \frac{D\mathbf{A}}{Dt} - (\nabla\mathbf{v})^T \cdot \mathbf{A} - \mathbf{A} \cdot (\nabla\mathbf{v})\]

where \(D/Dt = \partial_t + \mathbf{v} \cdot \nabla\) is the material derivative. For homogeneous flows (spatially uniform stress), the convective term \(\mathbf{v} \cdot \nabla\boldsymbol{\tau}\) vanishes and the material derivative reduces to the ordinary time derivative.

Simple shear geometry:

The velocity field \(\mathbf{v} = (\dot{\gamma} y, 0, 0)\) defines:

\[\begin{split}\nabla\mathbf{v} = \begin{pmatrix} 0 & \dot{\gamma} & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{pmatrix}, \qquad \mathbf{D} = \frac{1}{2}\begin{pmatrix} 0 & \dot{\gamma} & 0 \\ \dot{\gamma} & 0 & 0 \\ 0 & 0 & 0 \end{pmatrix}\end{split}\]

The polymer stress tensor in simple shear has the structure:

\[\begin{split}\boldsymbol{\tau} = \begin{pmatrix} \tau_{xx} & \tau_{xy} & 0 \\ \tau_{xy} & \tau_{yy} & 0 \\ 0 & 0 & \tau_{zz} \end{pmatrix}\end{split}\]

Constitutive Equation¶

The polymer stress \(\boldsymbol{\tau}\) satisfies the Giesekus constitutive equation [1]:

\[\boldsymbol{\tau} + \lambda \overset{\nabla}{\boldsymbol{\tau}} + \frac{\alpha \lambda}{\eta_p} \boldsymbol{\tau} \cdot \boldsymbol{\tau} = 2 \eta_p \mathbf{D}\]

The measurable total shear stress is:

\[\sigma_{xy} = \eta_s \dot{\gamma} + \tau_{xy}\]

The normal stress differences are:

\[N_1 = \tau_{xx} - \tau_{yy}, \qquad N_2 = \tau_{yy} - \tau_{zz}\]

Component-wise ODE System¶

In simple shear, the constitutive equation reduces to four coupled ODEs for the stress components [3]:

\[\frac{d\tau_{xx}}{dt} = -\frac{\tau_{xx}}{\lambda} + 2\dot{\gamma}\,\tau_{xy} - \frac{\alpha}{\eta_p}(\tau_{xx}^2 + \tau_{xy}^2)\]

\[\frac{d\tau_{yy}}{dt} = -\frac{\tau_{yy}}{\lambda} - \frac{\alpha}{\eta_p}(\tau_{xy}^2 + \tau_{yy}^2)\]

\[\frac{d\tau_{xy}}{dt} = -\frac{\tau_{xy}}{\lambda} + \dot{\gamma}\,\tau_{yy} - \frac{\alpha}{\eta_p}\tau_{xy}(\tau_{xx} + \tau_{yy}) + \frac{\eta_p}{\lambda}\dot{\gamma}\]

\[\frac{d\tau_{zz}}{dt} = -\frac{\tau_{zz}}{\lambda} - \frac{\alpha}{\eta_p}\tau_{zz}^2\]

Each equation has three contributions:

Linear relaxation: \(-\tau_{ij}/\lambda\) (exponential decay toward equilibrium)
Convective coupling: terms involving \(\dot{\gamma}\tau_{ij}\) (flow-induced stress transfer)
Quadratic nonlinearity: terms involving \(\alpha \tau^2/\eta_p\) (anisotropic drag)

The \(\tau_{zz}\) component decouples from the other three and relaxes to zero from any initial condition.

Dimensionless Formulation¶

Define dimensionless variables:

Weissenberg number: \(\text{Wi} = \lambda \dot{\gamma}\)
Dimensionless stress: \(\tau_{ij}^* = \tau_{ij} \lambda / \eta_p\)
Dimensionless time: \(t^* = t / \lambda\)

The ODEs become:

\[\frac{d\tau_{xx}^*}{dt^*} = -\tau_{xx}^* + 2\,\text{Wi}\;\tau_{xy}^* - \alpha(\tau_{xx}^{*2} + \tau_{xy}^{*2})\]

\[\frac{d\tau_{yy}^*}{dt^*} = -\tau_{yy}^* - \alpha(\tau_{xy}^{*2} + \tau_{yy}^{*2})\]

\[\frac{d\tau_{xy}^*}{dt^*} = -\tau_{xy}^* + \text{Wi}\;\tau_{yy}^* - \alpha\,\tau_{xy}^*(\tau_{xx}^* + \tau_{yy}^*) + \text{Wi}\]

This formulation is useful because:

All behavior is parameterized by just two numbers: \(\text{Wi}\) and \(\alpha\)
Universal behavior curves collapse data at different rates and relaxation times
Improved numerical conditioning when \(\eta_p/\lambda\) spans many orders of magnitude

Analytical Steady-State Solutions¶

At steady state (\(d/dt = 0\)), the ODE system reduces to a nonlinear algebraic system that admits closed-form solutions [1] [7].

Define the auxiliary discriminant:

\[\Lambda = \sqrt{1 + 16\,\alpha(1-\alpha)\,\text{Wi}^2}\]

and the auxiliary function:

\[f(\text{Wi}) = \frac{1 - \Lambda}{8\,\alpha(1-\alpha)\,\text{Wi}^2}\Bigl[1 + \Lambda + 2(1-2\alpha)\,\text{Wi}^2\Bigr]\]

The steady-state polymer stress components are:

\[\tau_{xy,\text{ss}} = \frac{\eta_p}{\lambda} \cdot \frac{(1 - f)\,\text{Wi}}{1 + (1-2\alpha)\,f}\]

\[\tau_{xx,\text{ss}} = \frac{\eta_p}{\lambda} \cdot \frac{2\,(1-f)^2\,\text{Wi}^2}{[1 + (1-2\alpha)\,f]\,[1 - \alpha\,f]}\]

\[\tau_{yy,\text{ss}} = \frac{\eta_p}{\lambda} \cdot \frac{-2\,\alpha\,f\,(1-f)\,\text{Wi}^2}{[1 + (1-2\alpha)\,f]\,[1 - \alpha\,f]}\]

\[\tau_{zz,\text{ss}} = 0\]

The steady-state shear viscosity is:

\[\eta(\dot{\gamma}) = \eta_s + \eta_p \cdot \frac{1 - f}{1 + (1 - 2\alpha)\,f}\]

where the term \((1-f)/[1 + (1-2\alpha)\,f]\) is the polymeric viscosity reduction factor.

Normal stress coefficients at steady state:

\[\Psi_1(\dot{\gamma}) = \frac{\tau_{xx,\text{ss}} - \tau_{yy,\text{ss}}}{\dot{\gamma}^2}, \qquad \Psi_2(\dot{\gamma}) = \frac{\tau_{yy,\text{ss}}}{\dot{\gamma}^2}\]

Limiting Behaviors¶

Quantity	Low Wi (\(\text{Wi} \ll 1\))	High Wi (\(\text{Wi} \gg 1\))	Notes
\(\eta\)	\(\eta_0 = \eta_p + \eta_s\)	\(\sim \text{Wi}^{-1}\)	Shear-thinning
\(\Psi_1\)	\(\Psi_{1,0} = 2\eta_p\lambda\)	\(\sim \text{Wi}^{-2}\)	Decreases
\(\Psi_2\)	\(\Psi_{2,0} = -\alpha\,\eta_p\lambda\)	\(\sim \text{Wi}^{-2}\)	Negative
\(N_2/N_1\)	\(-\alpha\)	\(-\alpha/2\)	Rate-independent at high Wi

Protocol-Specific Equations¶

This section presents the complete equations for each experimental protocol supported by the Giesekus model. Each protocol specifies the imposed kinematic or stress condition, the resulting ODE system (or algebraic system), initial conditions, and characteristic output observables.

Steady Shear (Flow Curve)¶

Protocol: Constant shear rate \(\dot{\gamma} = \text{const}\), solve at steady state (\(\partial_t \boldsymbol{\tau} = 0\)).

Governing system: Setting all time derivatives to zero in the component ODEs yields the nonlinear algebraic system:

\[\frac{\tau_{xx}}{\lambda} - 2\dot{\gamma}\,\tau_{xy} + \frac{\alpha}{\eta_p}(\tau_{xx}^2 + \tau_{xy}^2) = 0\]

\[\frac{\tau_{yy}}{\lambda} + \frac{\alpha}{\eta_p}(\tau_{xy}^2 + \tau_{yy}^2) = 0\]

\[\frac{\tau_{xy}}{\lambda} - \dot{\gamma}\,\tau_{yy} + \frac{\alpha}{\eta_p}\tau_{xy}(\tau_{xx} + \tau_{yy}) = \frac{\eta_p}{\lambda}\dot{\gamma}\]

\[\frac{\tau_{zz}}{\lambda} + \frac{\alpha}{\eta_p}\tau_{zz}^2 = 0 \quad \Rightarrow \quad \tau_{zz} = 0\]

Solution method: Use the analytical formulas from the previous section or Newton–Raphson iteration.

Output observables:

Flow curve: \(\sigma_{xy}(\dot{\gamma}) = \eta_s \dot{\gamma} + \tau_{xy,\text{ss}}(\dot{\gamma})\)
Viscosity: \(\eta(\dot{\gamma}) = \sigma_{xy}/\dot{\gamma}\)
Normal stresses: \(N_1 = \tau_{xx} - \tau_{yy} > 0\), \(N_2 = \tau_{yy} < 0\)

Shear-thinning behavior:

Wi range	\(\eta\) behavior	Physics
\(\text{Wi} \ll 1\)	\(\eta \approx \eta_0\)	Newtonian plateau
\(\text{Wi} \sim 1\)	Onset of thinning	Nonlinear drag effects begin
\(\text{Wi} \gg 1\)	\(\eta \sim \text{Wi}^{-1}\)	Power-law region

Normal stress ratio:

\[\frac{N_2}{N_1} \approx -\frac{\alpha}{2} \quad (\text{at high Wi})\]

This ratio is approximately independent of shear rate, making it the primary experimental route to determine \(\alpha\) [10] [11].

Startup of Steady Shear¶

Protocol: Apply constant shear rate from rest: \(\dot{\gamma}(t) = \dot{\gamma}_0\,H(t)\), where \(H(t)\) is the Heaviside step function.

Initial conditions: \(\tau_{xx}(0) = \tau_{yy}(0) = \tau_{xy}(0) = \tau_{zz}(0) = 0\)

ODE system:

\[\frac{d\tau_{xx}}{dt} = -\frac{\tau_{xx}}{\lambda} + 2\dot{\gamma}_0\,\tau_{xy} - \frac{\alpha}{\eta_p}(\tau_{xx}^2 + \tau_{xy}^2)\]

\[\frac{d\tau_{yy}}{dt} = -\frac{\tau_{yy}}{\lambda} - \frac{\alpha}{\eta_p}(\tau_{xy}^2 + \tau_{yy}^2)\]

\[\frac{d\tau_{xy}}{dt} = -\frac{\tau_{xy}}{\lambda} + \dot{\gamma}_0\,\tau_{yy} - \frac{\alpha}{\eta_p}\tau_{xy}(\tau_{xx} + \tau_{yy}) + \frac{\eta_p}{\lambda}\dot{\gamma}_0\]

Output: \(\sigma_{xy}(t) = \eta_s \dot{\gamma}_0 + \tau_{xy}(t)\) and \(N_1(t) = \tau_{xx}(t) - \tau_{yy}(t)\).

Characteristic features:

Time/strain regime	Behavior
\(t \ll \lambda\) (linear elastic)	\(\tau_{xy} \approx G\,\dot{\gamma}_0\,t\) (affine, slope = \(G\))
\(\gamma \sim O(1)\) (overshoot)	Stress peaks above steady state, \(N_1\) also overshoots
\(t \gg \lambda\) (steady state)	Stress relaxes to \(\tau_{xy,\text{ss}}\)

Overshoot characteristics:

Peak strain: \(\gamma_{\text{peak}} \sim 2\text{–}3\) strain units (depends on Wi and \(\alpha\))
Overshoot ratio: \(\sigma_{\text{peak}}/\sigma_{\text{ss}}\) increases with Wi
Higher \(\alpha\) gives smaller overshoot (stronger nonlinear damping)
High-Wi scaling: \(\gamma_{\text{peak}} \sim \text{const}\) (2–3), \(\sigma_{\text{peak}}/\sigma_{\text{ss}} \sim \text{Wi}^{1/2}\)

Stress Relaxation¶

Protocol: Apply instantaneous step strain \(\gamma_0\) at \(t = 0\), then \(\dot{\gamma}(t > 0) = 0\).

Initial conditions (from instantaneous elastic response):

\[\tau_{xy}(0^+) = G\,\gamma_0 = \frac{\eta_p}{\lambda}\,\gamma_0\]

\[\tau_{xx}(0^+) = 2\,G\,\gamma_0^2\]

\[\tau_{yy}(0^+) = 0, \qquad \tau_{zz}(0^+) = 0\]

Relaxation ODEs (with \(\dot{\gamma} = 0\)):

\[\frac{d\tau_{xx}}{dt} = -\frac{\tau_{xx}}{\lambda} - \frac{\alpha}{\eta_p}(\tau_{xx}^2 + \tau_{xy}^2)\]

\[\frac{d\tau_{yy}}{dt} = -\frac{\tau_{yy}}{\lambda} - \frac{\alpha}{\eta_p}(\tau_{xy}^2 + \tau_{yy}^2)\]

\[\frac{d\tau_{xy}}{dt} = -\frac{\tau_{xy}}{\lambda} - \frac{\alpha}{\eta_p}\tau_{xy}(\tau_{xx} + \tau_{yy})\]

Linear regime (small \(\gamma_0\), quadratic terms negligible):

\[G(t) = \frac{\tau_{xy}(t)}{\gamma_0} = G\,e^{-t/\lambda}\]

Nonlinear regime (finite \(\gamma_0\)):

The quadratic \(\alpha\)-terms accelerate relaxation when stress is high, giving faster-than-exponential initial decay:

\[\sigma(t) < \sigma_0 \exp(-t/\lambda)\]

Damping function (quantifies strain-dependent relaxation):

\[h(\gamma_0) = \frac{G(t, \gamma_0)}{G(t)} \quad \text{at early times}\]

For the Giesekus model, the instantaneous response obeys the Lodge–Meissner rule (\(h(\gamma) = 1\) at \(t = 0^+\)), but nonlinear effects emerge during the relaxation process.

Time-strain separability (approximate):

\[G(t, \gamma_0) \approx G(t) \cdot h(\gamma_0)\]

where \(G(t) = G\,e^{-t/\lambda}\) and \(h(\gamma_0)\) is the damping function.

Creep (Step Stress)¶

Protocol: Apply constant total shear stress \(\sigma_{xy}(t) = \sigma_0\,H(t)\).

Stress-control closure: The applied stress constraint gives:

\[\dot{\gamma}(t) = \frac{\sigma_0 - \tau_{xy}(t)}{\eta_s}\]

This makes the shear rate a dependent variable computed from the evolving polymer stress.

Coupled ODE system (5 equations: 4 stress + strain):

\[\frac{d\tau_{xx}}{dt} = -\frac{\tau_{xx}}{\lambda} + 2\dot{\gamma}\,\tau_{xy} - \frac{\alpha}{\eta_p}(\tau_{xx}^2 + \tau_{xy}^2)\]

\[\frac{d\tau_{yy}}{dt} = -\frac{\tau_{yy}}{\lambda} - \frac{\alpha}{\eta_p}(\tau_{xy}^2 + \tau_{yy}^2)\]

\[\frac{d\tau_{xy}}{dt} = -\frac{\tau_{xy}}{\lambda} + \dot{\gamma}\,\tau_{yy} - \frac{\alpha}{\eta_p}\tau_{xy}(\tau_{xx} + \tau_{yy}) + \frac{\eta_p}{\lambda}\dot{\gamma}\]

\[\frac{d\gamma}{dt} = \dot{\gamma}(t) = \frac{\sigma_0 - \tau_{xy}(t)}{\eta_s}\]

where \(\dot{\gamma}\) in the stress equations is evaluated from the closure at each time step.

Initial conditions: \(\tau_{xx} = \tau_{yy} = \tau_{xy} = \tau_{zz} = \gamma = 0\)

Creep compliance:

\[J(t) = \frac{\gamma(t)}{\sigma_0}\]

Limiting behaviors:

Time	\(J(t)\)	Physics
\(t \to 0^+\)	\(J_0 = 1/G = \lambda/\eta_p\)	Instantaneous elastic compliance
\(t \to \infty\)	\(J(t) \sim t/\eta_0\)	Steady-state viscous flow

Recovery after unloading (stress removed at \(t = t_1\)):

Elastic strain recovered: \(\Delta\gamma_{\text{rec}} \approx \sigma_0/G\)
Permanent (viscous) strain: \(\gamma_{\text{perm}} = \gamma(t_1) - \Delta\gamma_{\text{rec}}\)

Note

When \(\eta_s = 0\) (no solvent), the stress-control closure becomes singular. This case requires a DAE (differential-algebraic equation) solver or reformulation.

Small-Amplitude Oscillatory Shear (SAOS)¶

Protocol: \(\gamma(t) = \gamma_0 \sin(\omega t)\) with \(\gamma_0 \ll 1\).

In the linear limit, the quadratic \(\alpha\)-term is negligible and the Giesekus model reduces to the Oldroyd-B/Maxwell response. The SAOS moduli are therefore independent of \(\alpha\):

Storage modulus:

\[G'(\omega) = G \frac{(\omega\lambda)^2}{1 + (\omega\lambda)^2} = \frac{\eta_p \omega^2 \lambda}{1 + (\omega\lambda)^2}\]

Loss modulus:

\[G''(\omega) = G \frac{\omega\lambda}{1 + (\omega\lambda)^2} + \eta_s \omega = \frac{\eta_p \omega}{1 + (\omega\lambda)^2} + \eta_s \omega\]

where \(G = \eta_p/\lambda\) is the elastic modulus.

Complex viscosity:

\[\eta'(\omega) = \frac{\eta_p}{1 + (\omega\lambda)^2} + \eta_s\]

\[\eta''(\omega) = \frac{\eta_p \omega\lambda}{1 + (\omega\lambda)^2}\]

Limiting behaviors:

Frequency	\(G'\)	\(G''\)
\(\omega \to 0\)	\(G' \sim \omega^2\)	\(G'' \sim \omega\)
\(\omega \to \infty\)	\(G' \to G = \eta_p/\lambda\)	\(G'' \sim \eta_s \omega\) (solvent)

Crossover frequency:

\[\omega_c = 1/\lambda \quad \text{where } G'(\omega_c) = G''(\omega_c) - \eta_s\omega_c\]

Large-Amplitude Oscillatory Shear (LAOS)¶

Protocol: \(\gamma(t) = \gamma_0 \sin(\omega t)\) with \(\gamma_0\) finite.

The shear rate is \(\dot{\gamma}(t) = \gamma_0 \omega \cos(\omega t)\).

Full ODE system: The component equations are the same as the general simple shear system with the time-dependent \(\dot{\gamma}(t)\) inserted:

\[\frac{d\tau_{xx}}{dt} = -\frac{\tau_{xx}}{\lambda} + 2\dot{\gamma}(t)\,\tau_{xy} - \frac{\alpha}{\eta_p}(\tau_{xx}^2 + \tau_{xy}^2)\]

\[\frac{d\tau_{yy}}{dt} = -\frac{\tau_{yy}}{\lambda} - \frac{\alpha}{\eta_p}(\tau_{xy}^2 + \tau_{yy}^2)\]

\[\frac{d\tau_{xy}}{dt} = -\frac{\tau_{xy}}{\lambda} + \dot{\gamma}(t)\,\tau_{yy} - \frac{\alpha}{\eta_p}\tau_{xy}(\tau_{xx} + \tau_{yy}) + \frac{\eta_p}{\lambda}\dot{\gamma}(t)\]

The total shear stress is \(\sigma(t) = \tau_{xy}(t) + \eta_s \dot{\gamma}(t)\).

Fourier decomposition of the periodic stress response [16]:

\[\sigma(t) = \sum_{n=1,3,5,\ldots} \bigl[\sigma_n' \sin(n\omega t) + \sigma_n'' \cos(n\omega t)\bigr]\]

Only odd harmonics appear due to the symmetry of shear flow.

First harmonic moduli (strain-amplitude dependent):

\[G_1'(\omega, \gamma_0) = \frac{\sigma_1'}{\gamma_0}, \qquad G_1''(\omega, \gamma_0) = \frac{\sigma_1''}{\gamma_0}\]

Third harmonic ratio (primary nonlinearity measure):

\[I_{3/1} = \frac{\sqrt{\sigma_3'^2 + \sigma_3''^2}}{\sqrt{\sigma_1'^2 + \sigma_1''^2}}\]

MAOS scaling (medium-amplitude regime):

\[I_{3/1} \sim \gamma_0^2 \quad \text{as } \gamma_0 \to 0\]

Chebyshev decomposition:

\[\sigma(\gamma, \dot{\gamma}) = \gamma_0 \sum_{n \text{ odd}} \bigl[e_n\,T_n(x) + v_n\,T_n(y)\bigr]\]

where \(x = \gamma/\gamma_0\), \(y = \dot{\gamma}/(\gamma_0\omega)\), and \(T_n\) are Chebyshev polynomials.

Pipkin diagram regimes:

\(\text{De} = \omega\lambda\)	\(\text{Wi} = \gamma_0 \omega \lambda\)	Regime
Any	\(\ll 1\)	Linear viscoelastic (SAOS)
\(\ll 1\)	Any	Quasi-steady nonlinear
\(\gg 1\)	\(\gg 1\)	Highly nonlinear viscoelastic

Giesekus LAOS signatures:

Feature	Giesekus behavior
Strain softening	\(G_1'\) decreases with \(\gamma_0\) (from \(\alpha > 0\))
Higher harmonics	Present due to quadratic stress term
Lissajous shape	Ellipse, tilted/distorted at high \(\gamma_0\)
\(I_{3/1}\) scaling	\(I_{3/1} \sim \gamma_0^2\) in MAOS regime

Multi-Mode Giesekus¶

Motivation¶

Real polymer systems have a broad spectrum of relaxation times arising from polydispersity and the range of molecular conformations. A single-mode Giesekus model cannot capture the broad frequency dependence typically observed in \(G'(\omega)\) and \(G''(\omega)\) data. The multi-mode extension addresses this by superposing \(N\) independent Giesekus modes.

Constitutive Equation¶

The total stress is:

\[\boldsymbol{\sigma} = -p\mathbf{I} + 2\eta_s\mathbf{D} + \sum_{k=1}^{N} \boldsymbol{\tau}_k\]

where each mode \(k\) evolves independently:

\[\boldsymbol{\tau}_k + \lambda_k \overset{\nabla}{\boldsymbol{\tau}_k} + \frac{\alpha_k \lambda_k}{\eta_{p,k}} \boldsymbol{\tau}_k \cdot \boldsymbol{\tau}_k = 2\eta_{p,k}\,\mathbf{D}\]

Each mode has its own relaxation time \(\lambda_k\), polymer viscosity \(\eta_{p,k}\), and mobility factor \(\alpha_k\).

Linear Viscoelastic Spectra (SAOS)¶

For multi-mode SAOS, the moduli superpose linearly:

\[G'(\omega) = \sum_{k=1}^{N} \frac{\eta_{p,k}\,\omega^2\,\lambda_k}{1 + (\omega\lambda_k)^2}\]

\[G''(\omega) = \sum_{k=1}^{N} \frac{\eta_{p,k}\,\omega}{1 + (\omega\lambda_k)^2} + \eta_s\,\omega\]

Zero-Shear Properties¶

\[\eta_0 = \sum_{k=1}^{N} \eta_{p,k} + \eta_s\]

\[\Psi_{1,0} = 2 \sum_{k=1}^{N} \eta_{p,k}\,\lambda_k\]

\[\Psi_{2,0} = -\sum_{k=1}^{N} \alpha_k\,\eta_{p,k}\,\lambda_k\]

Multi-Mode ODE State Vector¶

For transient simulations, the state vector has \(4N\) components (4 stress components per mode). Each mode evolves independently with its own parameters but shares the same velocity field \(\dot{\gamma}(t)\):

\[\mathbf{y} = [\tau_{xx}^{(1)}, \tau_{yy}^{(1)}, \tau_{xy}^{(1)}, \tau_{zz}^{(1)}, \ldots, \tau_{xx}^{(N)}, \tau_{yy}^{(N)}, \tau_{xy}^{(N)}, \tau_{zz}^{(N)}]\]

Fitting Strategy¶

Discrete spectrum from SAOS: Fit \(G_k = \eta_{p,k}/\lambda_k\) and \(\lambda_k\) to SAOS data
Logarithmic spacing: Place \(\lambda_k\) at logarithmically spaced points across the frequency window
Regularization: Use non-negative least squares (NNLS) or Tikhonov regularization to avoid overfitting
Typical mode count: 5–10 modes cover 4–6 decades in frequency
Fix \(\alpha_k\) from nonlinear data: The linear spectrum determines \(\eta_{p,k}\) and \(\lambda_k\); fit \(\alpha_k\) to flow curve or normal stress data

Parameters¶

Giesekus Model Parameters¶
Parameter	Symbol	Units	Bounds	Physical Meaning
\(\eta_p\)	\(\eta_p\)	Pa·s	(1e-3, 1e6)	Polymer zero-shear viscosity
\(\lambda\)	\(\lambda\)	s	(1e-6, 1e4)	Characteristic relaxation time
\(\alpha\)	\(\alpha\)	—	[0, 0.5]	Mobility anisotropy factor
\(\eta_s\)	\(\eta_s\)	Pa·s	[0, 1e4)	Solvent/Newtonian viscosity

Parameter Interpretation¶

Polymer viscosity \(\eta_p\):

Dominant contribution to zero-shear viscosity
Scales with molecular weight: \(\eta_p \sim M_w^{3.4}\) above entanglement
Temperature dependent via Arrhenius/WLF

Relaxation time \(\lambda\):

Time for stress to decay to 1/e of initial value
Scales with molecular weight: \(\lambda \sim M_w^{3.4}\)
Defines crossover frequency: \(\omega_c = 1/\lambda\)

Mobility factor \(\alpha\):

\(\alpha = 0\): Isotropic mobility (UCM limit)
\(\alpha = 0.5\): Maximum anisotropy
Directly measurable: \(\alpha = -2 N_2/N_1\)
Typical values: - Polymer melts: 0.1–0.3 - Concentrated solutions: 0.2–0.4 - Wormlike micelles: 0.3–0.5

Solvent viscosity \(\eta_s\):

Newtonian background contribution
Important for dilute/semi-dilute solutions
Often negligible for melts (\(\eta_s \ll \eta_p\))

Physical Constraints¶

\(0 \leq \alpha \leq 0.5\) for most physical systems
\(\alpha > 0.5\) can produce unphysical behavior at high Wi (non-monotonic flow curves)
\(\eta_s \geq 0\), \(\lambda > 0\), \(\eta_p > 0\)

Typical Parameter Ranges by Material¶

Material	\(\eta_p\) (Pa·s)	\(\lambda\) (s)	\(\alpha\)	\(\eta_s\) (Pa·s)
Polymer solutions	0.1–1000	0.001–10	0.1–0.5	0.001–1
Polymer melts	100–10⁶	0.1–1000	0.1–0.5	~0
Wormlike micelles	1–100	0.1–10	0.3–0.5	0.001–0.1

Derived Quantities¶

Zero-shear viscosity: \(\eta_0 = \eta_p + \eta_s\)
Elastic modulus: \(G = \eta_p/\lambda\)
Weissenberg number: \(\text{Wi} = \lambda \dot{\gamma}\)
Deborah number: \(\text{De} = \lambda/t_{\text{obs}}\)

Validity and Assumptions¶

Model Assumptions¶

Incompressibility: Constant density during deformation
Homogeneous deformation: No spatial gradients in material properties
Isothermal conditions: Temperature held constant
Upper-convected derivative: Frame-invariant stress transport
Single relaxation time: Monodisperse or narrow distribution

Validity Range¶

Condition	Range	Notes
Weissenberg number	\(\text{Wi} \lesssim 100\)	Numerical stability limit
Shear rate	\(\dot{\gamma} < 1/\lambda\) to \(100/\lambda\)	Power-law region
Strain (startup)	\(\gamma \lesssim 10\)	Overshoot captured
Temperature	Near reference T	Use TTS for other temperatures

Limitations¶

Single relaxation time: Real polymers have spectra (use multi-mode)
No extensional hardening: Underpredicts extensional viscosity
Fixed \(N_2/N_1\) ratio: Cannot vary independently
Numerical stiffness: High Wi may require adaptive solvers

When NOT to Use¶

Extensional flows: Use FENE-P or PTT for extensional hardening
Broad relaxation spectra: Use multi-mode Giesekus
Thixotropic materials: Use fluidity models
Yield stress fluids: Use EVP models (Saramito)

Regimes and Behavior¶

Weissenberg Number Regimes¶

Regime	Wi Range	Viscosity	Physics
Newtonian	\(\text{Wi} \ll 1\)	\(\eta \approx \eta_0\)	Linear response, no thinning
Transition	\(\text{Wi} \sim 1\)	Onset of thinning	Nonlinear effects begin
Power-law	\(\text{Wi} \gg 1\)	\(\eta \sim \text{Wi}^{n-1}\)	Strong shear-thinning

Effect of \(\alpha\) on Behavior¶

\(\alpha\) value	Shear-thinning	\(N_2/N_1\)	Example materials
0	None (UCM)	0	Ideal elastic liquid
0.1	Weak	−0.05	Some polymer melts
0.3	Moderate	−0.15	Typical polymers
0.5	Maximum	−0.25	Wormlike micelles

What You Can Learn¶

From SAOS Data¶

Extractable parameters: \(\eta_p\), \(\lambda\), \(\eta_s\), \(G\)

Observable	Extracted quantity
\(G''/\omega\) as \(\omega \to 0\)	\(\eta_0 = \eta_p + \eta_s\)
\(G''/\omega\) as \(\omega \to \infty\)	\(\eta_s\) (high-frequency limit)
Crossover \(G' = G'' - \eta_s\omega\)	\(\lambda = 1/\omega_c\)
\(G'\) plateau (\(\omega \to \infty\))	\(G = \eta_p/\lambda\)

What this reveals:

Elastic modulus \(G\): Network strength / entanglement density
Relaxation time \(\lambda\): Molecular weight, longest relaxation mode
Relaxation spectrum width: Single-mode fit quality indicates how narrow/broad the spectrum is

Note

The mobility parameter \(\alpha\) is not determinable from SAOS data because SAOS is \(\alpha\)-independent (linear regime).

From Steady Shear (Flow Curve)¶

Extractable parameters: \(\eta_0\), \(\lambda\) (onset), \(\alpha\) (shape)

Observable	Extracted quantity
Low-rate plateau \(\eta(\dot{\gamma} \to 0)\)	\(\eta_0 = \eta_p + \eta_s\)
Onset shear rate for thinning	\(\lambda \approx 1/\dot{\gamma}_{\text{onset}}\)
Shape of thinning curve	\(\alpha\) (controls power-law slope)

What this reveals:

Molecular weight (via \(\eta_0\) and \(\lambda\) scaling laws)
Entanglement density: \(\eta_0 \sim c^{3.4}\) for entangled systems
Cross-validation: Compare \(\eta_0\) and \(\lambda\) from SAOS

From Normal Stress Measurements¶

Primary output: Direct \(\alpha\) determination.

\[\alpha = -\frac{2\,N_2}{N_1}\]

What this reveals:

Degree of molecular anisotropy: Higher \(\alpha\) indicates more anisotropic drag
Material classification: Polymer melts (\(\alpha \sim 0.1\text{–}0.3\)), wormlike micelles (\(\alpha \sim 0.3\text{–}0.5\))
Experimental techniques: Cone-and-plate (for \(N_1\)), parallel-plate edge measurements or cone-partitioned plate (for \(N_2\))

From Startup Flow¶

Primary outputs:

Overshoot ratio \(\sigma_{\text{max}}/\sigma_{\text{ss}}\): Increases with Wi, quantifies nonlinear viscoelastic character
Strain at peak: \(\gamma_{\text{peak}} \sim 2\text{–}3\) — network deformation scale
Time to steady state: \(\sim 3\text{–}5\lambda\) — validates relaxation time
Initial slope: \(d\sigma/d\gamma|_{t \to 0} = G\) — instantaneous elastic modulus

From Stress Relaxation¶

Primary outputs:

Exponential vs. faster-than-exponential decay: Faster initial decay confirms \(\alpha > 0\)
Relaxation modulus: \(G(t) = \tau_{xy}(t)/\gamma_0\) — full time-dependent response
Damping function \(h(\gamma_0)\): Quantifies nonlinear strain effects (strain thinning)
Time-strain separability: Whether \(G(t, \gamma_0) \approx G(t) \cdot h(\gamma_0)\) holds

From Creep¶

Primary outputs:

Instantaneous compliance: \(J_0 = 1/G = \lambda/\eta_p\)
Steady-state viscosity: Long-time slope \(dJ/dt \to 1/\eta_0\)
Elastic recovery (after unloading): Recoverable strain \(\approx \sigma_0/G\)
Retardation spectrum: Transition from elastic to viscous response

From LAOS¶

Primary outputs [16]:

Strain softening onset: Critical \(\gamma_0\) where \(G_1'\) begins decreasing — identifies linear-to-nonlinear transition
Third harmonic ratio \(I_{3/1}\): Quantifies strength of nonlinearity
MAOS scaling \(I_{3/1} \sim \gamma_0^2\): Material time exponent from intrinsic nonlinearity
Lissajous shapes: Visual nonlinear fingerprint — ellipse distortion at high amplitude
Chebyshev coefficients \(e_n, v_n\): Decompose intracycle elastic and viscous nonlinearity

Combined Multi-Protocol Analysis¶

Recommended fitting sequence:

SAOS: \(\eta_p, \lambda, \eta_s\) (linear parameters, \(\alpha\)-independent)
Flow curve: refine \(\eta_p, \lambda\); determine \(\alpha\) from thinning shape
Normal stresses: fix \(\alpha = -2N_2/N_1\) (most direct route)
Startup: validate overshoot predictions, refine \(\alpha\)
Relaxation/Creep: confirm time constants, validate nonlinear response

Parameter-to-data mapping:

Data type	\(\eta_p\)	\(\lambda\)	\(\alpha\)	\(\eta_s\)	Strength
SAOS	✓	✓	—	✓	Best for linear params
Flow curve	✓	✓	✓	✓	Thinning shape gives \(\alpha\)
\(N_1, N_2\)	—	—	✓✓	—	Most direct \(\alpha\)
Startup	✓	✓	✓	—	Overshoot validates model
Relaxation	✓	✓	(✓)	—	Decay rate confirms \(\lambda\)
Creep	✓	✓	(✓)	✓	Compliance confirms \(G\)

(✓✓ = primary route; ✓ = determinable; (✓) = weakly sensitive; — = not accessible)

Experimental Design¶

When to Use Giesekus¶

Use the Giesekus model when your material exhibits:

Shear-thinning viscosity
Measurable \(N_2\) (negative second normal stress difference)
Stress overshoot in startup flow
SAOS that fits Maxwell/Generalized Maxwell
Single or narrow relaxation time distribution

Decision Tree¶

Is N_2 measurable (negative)?
├── YES → Giesekus captures N_2/N_1 = -α/2
│
└── NO → Is only shear-thinning needed?
    ├── YES → Consider simpler Carreau/Cross
    └── NO → Consider PTT or FENE-P for extensional

Recommended Protocol Sequence¶

SAOS first: Determine \(\eta_p\), \(\lambda\), \(\eta_s\) from linear regime
Flow curve: Confirm thinning, refine parameters
Normal stresses: Measure \(N_2/N_1\) to determine \(\alpha\)
Startup flow: Validate overshoot predictions
Relaxation: Confirm faster-than-exponential decay

Material-Specific Recommendations¶

Material	Typical \(\alpha\)	n_modes	Key protocols
Polymer melts	0.1–0.3	3–5	Flow curve + SAOS + \(N_2\)
Polymer solutions	0.2–0.4	1–3	Startup + SAOS
Wormlike micelles	0.3–0.5	1	Startup overshoot + relaxation
Biological fluids	0.2–0.4	2–3	SAOS + low-Wi flow curve

Computational Implementation¶

RheoJAX Implementation¶

The Giesekus model in RheoJAX uses:

JAX acceleration: JIT-compiled kernels for fast predictions
diffrax integration: Adaptive ODE solvers (Tsit5) for transients
Analytical solutions: Where available (steady shear, SAOS)
Float64 precision: Essential for accurate stress calculations

Architecture¶

GiesekusBase (ABC)
├── GiesekusSingleMode
│   ├── Analytical: flow_curve, SAOS
│   └── ODE: startup, relaxation, creep, LAOS
│
└── GiesekusMultiMode
    ├── SAOS superposition (analytical)
    └── Extended state vector ODE

Numerical Considerations¶

Steady-state solver:

Newton iteration for auxiliary function f(Wi)
Converges in 5–10 iterations typically
May need damping at very high Wi

ODE integration:

Tsit5 (Runge-Kutta 5(4)) for accuracy
Adaptive step size with PIDController
rtol=1e-6, atol=1e-8 default tolerances

Numerical stability:

High Wi (>100) may require reduced tolerances
Very small \(\alpha\) (<0.01) approaches UCM singularities
Use log-residuals for fitting flow curves

Fitting Guidance¶

Initial Parameter Estimates¶

From SAOS data:

# At crossover (G' = G'')
lambda_1 = 1 / omega_crossover
G = G_prime_at_crossover * 2  # G' = G'' = G/2 at crossover
eta_p = G * lambda_1

From flow curve:

# Zero-shear plateau
eta_0 = stress[0] / gamma_dot[0]  # At lowest rate

# Onset of thinning
lambda_1 = 1 / gamma_dot_onset  # Where η starts dropping

\(\alpha\) estimation:

# From normal stresses (if available)
alpha = -2 * N2 / N1

# From thinning slope (rough estimate)
# High-Wi slope of η vs γ̇ in log-log ≈ (n-1)
# For Giesekus: n ≈ 0.5 at alpha = 0.5

From transient data:

Observable	Estimated parameter
Time to steady state	\(\lambda \approx t_{\text{ss}} / (3\text{–}5)\)
Overshoot magnitude	Higher \(\alpha\) gives smaller overshoot
Initial slope in startup	\(G = \eta_p/\lambda\) from \(d\sigma/d\gamma\|_0\)

Parameter Estimation Summary¶

Parameter	From SAOS	From steady shear	From transient
\(\eta_p\)	\(\eta_0 - \eta_s\)	Low-rate plateau − \(\eta_s\)	Initial slope / \(\lambda\)
\(\lambda\)	\(1/\omega_c\)	\(1/\dot{\gamma}_{\text{onset}}\)	\(t_{\text{ss}}/(3\text{–}5)\)
\(\alpha\)	— (not accessible)	Thinning shape; \(-2N_2/N_1\)	Overshoot ratio
\(\eta_s\)	High-\(\omega\) \(G''/\omega\)	High-rate plateau	—

Fitting Strategy¶

Fix \(\eta_s\) if known (pure solvent viscosity)
Fit SAOS first for \(\eta_p\), \(\lambda\) (\(\alpha\)-independent)
Fit flow curve to refine and get \(\alpha\)
Validate with startup for dynamic behavior

Multi-Mode Fitting¶

Discrete spectrum from SAOS: Fit \(G_k, \lambda_k\) pairs to \(G'(\omega), G''(\omega)\) using logarithmically spaced relaxation times
Non-negative least squares (NNLS): Ensures \(\eta_{p,k} \geq 0\)
Tikhonov regularization: Prevents overfitting when the number of modes exceeds data quality
Fix \(\alpha_k\) after linear fit: Determine mobility factors from nonlinear data (flow curve, normal stresses)
Typical: 5–10 modes for 4–6 decades in frequency

Troubleshooting¶

Problem	Likely Cause	Solution
Poor flow curve fit	Wrong \(\alpha\)	Use \(N_2/N_1\) to fix \(\alpha\), then fit others
Overshoot too small	\(\alpha\) too low	Increase \(\alpha\) toward 0.5
No convergence at high Wi	Numerical stiffness	Reduce max Wi, use adaptive solver
Relaxation too slow	\(\lambda\) too long	Fit SAOS crossover more carefully
SAOS mismatch	Single mode inadequate	Use multi-mode Giesekus

Usage Examples¶

Basic Single-Mode¶

from rheojax.models.giesekus import GiesekusSingleMode
import numpy as np

# Create model with parameters
model = GiesekusSingleMode()
model.parameters.set_value("eta_p", 100.0)  # Pa·s
model.parameters.set_value("lambda_1", 1.0)  # s
model.parameters.set_value("alpha", 0.3)     # dimensionless
model.parameters.set_value("eta_s", 10.0)    # Pa·s

# Predict flow curve
gamma_dot = np.logspace(-2, 2, 50)
sigma = model.predict(gamma_dot, test_mode='flow_curve')

# Get viscosity
_, eta, _ = model.predict_flow_curve(gamma_dot, return_components=True)

Predict SAOS¶

# SAOS is alpha-independent (linear regime)
omega = np.logspace(-2, 3, 50)
G_prime, G_double_prime = model.predict_saos(omega)

# Complex modulus
G_star = np.sqrt(G_prime**2 + G_double_prime**2)

Normal Stress Prediction¶

# Normal stress differences
gamma_dot = np.logspace(-1, 2, 30)
N1, N2 = model.predict_normal_stresses(gamma_dot)

# Verify diagnostic ratio
ratio = N2 / N1  # Should equal -alpha/2 = -0.15 (for alpha=0.3)

Startup with Overshoot¶

# Startup flow at constant rate
t = np.linspace(0, 10, 500)
sigma_t = model.simulate_startup(t, gamma_dot=10.0)

# Find overshoot
sigma_max = np.max(sigma_t)
sigma_ss = sigma_t[-1]
overshoot_ratio = sigma_max / sigma_ss  # > 1 indicates overshoot

# Get full stress tensor evolution
result = model.simulate_startup(t, gamma_dot=10.0, return_full=True)

Multi-Mode Giesekus¶

from rheojax.models.giesekus import GiesekusMultiMode

# Create 3-mode model
model = GiesekusMultiMode(n_modes=3)

# Set per-mode parameters
model.set_mode_params(0, eta_p=100.0, lambda_1=10.0, alpha=0.3)
model.set_mode_params(1, eta_p=50.0, lambda_1=1.0, alpha=0.25)
model.set_mode_params(2, eta_p=20.0, lambda_1=0.1, alpha=0.2)
model.parameters.set_value("eta_s", 5.0)

# SAOS captures broad spectrum
omega = np.logspace(-3, 3, 100)
G_prime, G_double_prime = model.predict_saos(omega)

Bayesian Fitting¶

from rheojax.core.data import RheoData

# Create data object
data = RheoData(x=omega, y=G_star, test_mode='oscillation')

# NLSQ warm-start
model.fit(data)

# Bayesian inference
result = model.fit_bayesian(
    data,
    num_warmup=1000,
    num_samples=2000,
    num_chains=4,
    seed=42
)

# Get credible intervals
intervals = model.get_credible_intervals(result.posterior_samples)

Model Comparison¶

vs. Upper-Convected Maxwell (UCM)¶

Feature	UCM (\(\alpha = 0\))	Giesekus (\(\alpha > 0\))
Viscosity	Constant	Shear-thinning
\(N_1\)	Positive	Positive
\(N_2\)	Zero	Negative
Startup	Overshoot (weak)	Overshoot (strong)
Relaxation	Exponential	Faster than exponential

vs. Phan-Thien–Tanner (PTT)¶

Feature	Giesekus	PTT
Thinning mechanism	Anisotropic drag	Network destruction
\(N_2/N_1\)	Fixed = \(-\alpha/2\)	Adjustable
Extensional	Bounded	Bounded (stronger)
Parameters	4	4-5
Best for	Shear flows	Mixed flows

vs. FENE-P¶

Feature	Giesekus	FENE-P
Mechanism	Anisotropic drag	Finite extensibility
Extensional	Moderate	Strong hardening
Shear thinning	Strong	Moderate
\(N_2\)	Nonzero	Zero
Best for	Shear + \(N_2\)	Extensional flows