HVNM Protocol Equations & Derivations

This page provides detailed derivations of HVNM predictions for each rheological protocol. All results assume the four-subnetwork model (P + E + D + I) in simple shear with constant kinetic rates unless otherwise stated.

For the governing equations, see HVNMLocal — Full Model Reference. For HVM protocol derivations (P + E + D), see HVM Protocol Equations & Derivations. For the underlying VLB single-network results, see VLB Protocol Equations & Derivations.

Simple Shear Geometry

ODE state vector (17 components, 18 with \(D_{int}\)):

[mu_E_xx, mu_E_yy, mu_E_xy,              # E-network distribution (3)
 mu_E_nat_xx, mu_E_nat_yy, mu_E_nat_xy,  # E-network natural state (3)
 mu_D_xx, mu_D_yy, mu_D_xy,              # D-network distribution (3)
 mu_I_xx, mu_I_yy, mu_I_xy,              # I-network distribution (3)
 mu_I_nat_xx, mu_I_nat_yy, mu_I_nat_xy,  # I-network natural state (3)
 gamma,                                    # accumulated strain (1)
 D]                                        # permanent damage (1)
 # D_int                                   # interfacial damage (if enabled)

The P-network is tracked analytically: \(\sigma_{P,xy} = (1-D) G_P X(\phi) \gamma(t)\).

I-network component equations (parallel to E-network, with amplification):

\[\begin{split}\dot{\mu}^I_{xx} &= 2 X_I \dot{\gamma}\,\mu^I_{xy} + k_{BER}^{int}(\mu^{I,nat}_{xx} - \mu^I_{xx}) \\ \dot{\mu}^I_{yy} &= k_{BER}^{int}(\mu^{I,nat}_{yy} - \mu^I_{yy}) \\ \dot{\mu}^I_{xy} &= X_I \dot{\gamma}\,\mu^I_{yy} + k_{BER}^{int}(\mu^{I,nat}_{xy} - \mu^I_{xy})\end{split}\]

I-network natural-state evolution:

\[\dot{\mu}^{I,nat}_{ij} = k_{BER}^{int}(\mu^I_{ij} - \mu^{I,nat}_{ij})\]

These have the identical mathematical form as the E-network in HVM, with \(k_{BER}^{mat} \to k_{BER}^{int}\) and \(\dot{\gamma} \to X_I \dot{\gamma}\).

Shorthand: \(k_m \equiv k_{BER}^{mat,0}(T)\), \(k_I \equiv k_{BER}^{int,0}(T)\), \(k_D \equiv k_d^D(T)\), \(\hat{\tau}_m = 1/(2k_m)\), \(\hat{\tau}_I = 1/(2k_I)\), \(\tau_D = 1/k_D\), \(X = X(\phi)\), \(X_I = X(\phi_{eff})\).

Flow Curve Derivation

Objective: Steady-state \(\sigma_{12}(\dot{\gamma})\).

Step 1 – Permanent network:

\(\sigma_{P,xy} = (1-D)\tilde{G}_P \dot{\gamma} t\) (unbounded, no steady state).

Step 2 – Exchangeable and interphase at steady state:

Both have evolving natural states, so at true steady state \(\mu^E = \mu^{E,nat}\) and \(\mu^I = \mu^{I,nat}\):

\[\sigma_E^{ss} = 0, \quad \sigma_I^{ss} = 0\]

Step 3 – Dissociative network (unchanged from HVM):

\[\sigma_D^{ss} = \eta_D \dot{\gamma}\]

Total steady-state stress (transient networks only):

\[\boxed{\sigma_{12}^{ss} = \eta_D \dot{\gamma}}\]

Key result

Both the exchangeable and interphase networks carry zero stress at true steady state – their natural states fully track the deformation. Only the dissociative network contributes viscous stress.

Intermediate quasi-steady state: On timescales \(t \gg \hat{\tau}_m\) but \(t \ll \hat{\tau}_I\), the matrix BER has relaxed but the interphase has not equilibrated:

\[\sigma_{12}^{qs} \approx (1-D)\tilde{G}_P \dot{\gamma} t + (1-D_{int}) G_{I,eff} X_I \dot{\gamma} \hat{\tau}_I + \eta_D \dot{\gamma}\]

Startup Shear Derivation

Objective: \(\sigma_{12}(t)\) for constant \(\dot{\gamma}\) from rest.

Initial conditions: All tensors at \(\mathbf{I}\), \(D = D_{int} = 0\).

Permanent Network

\[\sigma_{P,xy}(t) = \tilde{G}_P X \dot{\gamma} t\]

Exchangeable Network

Identical to HVM (see Startup Shear Derivation):

\[\sigma_{E,xy}^+(t) = \frac{G_E^{eff} \dot{\gamma}}{2 k_m} \left(1 - e^{-2 k_m t}\right) = \frac{\eta_E^{eff}}{2} \dot{\gamma} \left(1 - e^{-t/\hat{\tau}_m}\right)\]

Dissociative Network

\[\sigma_{D,xy}^+(t) = \eta_D \dot{\gamma}\left(1 - e^{-t/\tau_D}\right)\]

Interphase Network (New)

The I-network follows the same mathematics as the E-network with \(k_m \to k_I\) and \(\dot{\gamma} \to X_I \dot{\gamma}\):

\[\boxed{\sigma_{I,xy}^+(t) = \frac{(1-D_{int}) G_{I,eff} X_I \dot{\gamma}}{2 k_I} \left(1 - e^{-2 k_I t}\right) = \frac{(1-D_{int}) \eta_I^{eff}}{2} \dot{\gamma} \left(1 - e^{-t/\hat{\tau}_I}\right)}\]

where \(\eta_I^{eff} = G_{I,eff} X_I / k_I\).

Total Startup Stress

\[\sigma_{xy}^+(t) = \tilde{G}_P X \dot{\gamma} t + \frac{\eta_E^{eff}}{2} \dot{\gamma} (1 - e^{-t/\hat{\tau}_m}) + \eta_D \dot{\gamma} (1 - e^{-t/\tau_D}) + \frac{(1-D_{int}) \eta_I^{eff}}{2} \dot{\gamma} (1 - e^{-t/\hat{\tau}_I})\]

Short-time (elastic) limit (\(t \ll \min(\hat{\tau}_m, \tau_D, \hat{\tau}_I)\)):

\[\sigma_{xy}^+(t) \approx G_{tot}^{NC} \dot{\gamma} t\]

where \(G_{tot}^{NC} = \tilde{G}_P X + G_E^{eff} + G_D^{eff} + (1-D_{int}) G_{I,eff} X_I\) is the total instantaneous nanocomposite modulus.

Double-overshoot signature: With TST kinetics, the BER rates accelerate as stress builds. The interphase, with its higher activation energy and strain amplification \(X_I\), shows a later and larger overshoot than the matrix exchangeable network. This produces a distinctive double-overshoot in startup flow at certain temperatures – a key experimental fingerprint of dual-kinetics nanocomposites.

Stress Relaxation Derivation

Objective: \(G(t) = \sigma_{12}(t)/\gamma_0\) after step strain \(\gamma_0\).

Kinematics: For \(t > 0\): \(\dot{\gamma} = 0\).

Initial conditions (immediately after step): All distribution tensors at \(\mu_{xy} = \gamma_0\) (affine), all natural-state tensors at \(\mu^{nat}_{xy} = 0\).

Each subnetwork relaxes independently:

Permanent: \(\sigma_P(t) = (1-D) \tilde{G}_P X \gamma_0\) (plateau).

Exchangeable: \(\sigma_E(t) = G_E^{eff} \gamma_0 e^{-t/\hat{\tau}_m}\).

Dissociative: \(\sigma_D(t) = G_D^{eff} \gamma_0 e^{-t/\tau_D}\).

Interphase: \(\sigma_I(t) = (1-D_{int}) G_{I,eff} X_I \gamma_0 e^{-t/\hat{\tau}_I}\).

Total relaxation modulus:

\[\boxed{G(t) = (1-D) \tilde{G}_P X + G_E^{eff} e^{-t/\hat{\tau}_m} + G_D^{eff} e^{-t/\tau_D} + (1-D_{int}) G_{I,eff} X_I e^{-t/\hat{\tau}_I}}\]

This is a four-mode relaxation spectrum (compared to two for HVM):

Mode

Timescale

Origin

Fast

\(\tau_D\) (ms–s)

Physical bond turnover

Intermediate

\(\hat{\tau}_m = \tau_m/2\) (min–hr)

Matrix BER

Slow

\(\hat{\tau}_I = \tau_I/2\) (hr–days)

Interfacial BER (new NP mode)

Plateau

\(\infty\)

\((1-D) \tilde{G}_P X\) (amplified permanent)

The appearance of a distinct slow relaxation mode associated with the interphase is a key experimental signature that distinguishes vitrimer nanocomposites from unfilled vitrimers.

Verification conditions:

  • \(G(0^+) = G_{tot}^{NC}\) (instantaneous modulus)

  • \(G(\infty) = (1-D) \tilde{G}_P X\) (amplified equilibrium modulus)

With diffusion-limited slow mode (Karim et al. 2025):

When include_diffusion=True, each BER mode acquires a long-time tail:

\[G(t)_{diff} = G_E^{eff} e^{-t/\hat{\tau}_m} e^{-k_{diff}^{mat} t} + (1-D_{int}) G_{I,eff} X_I e^{-t/\hat{\tau}_I} e^{-k_{diff}^{int} t} + \ldots\]

The \(e^{-k_{diff} t}\) factors produce slow exponential decay beyond BER relaxation, capturing experimentally observed incomplete relaxation at intermediate times.

Creep Derivation

Objective: \(J(t) = \gamma(t)/\sigma_0\) under constant stress.

Instantaneous response: \(\gamma(0^+) = \sigma_0 / G_{tot}^{NC}\).

Long-time saturation: \(\gamma(\infty) = \sigma_0 / [(1-D)\tilde{G}_P X]\).

Full creep compliance (Prony series form):

\[J(t) = \frac{1}{(1-D)\tilde{G}_P X} + \sum_{\alpha \in \{E,D,I\}} \frac{G_\alpha^{eff}} {(1-D)\tilde{G}_P X \cdot ((1-D)\tilde{G}_P X + G_\alpha^{eff})} \left(1 - e^{-t/t_c^\alpha}\right)\]

where the three retardation times are:

\[\begin{split}t_c^E &= \hat{\tau}_m \cdot \frac{(1-D)\tilde{G}_P X + G_D^{eff} + (1-D_{int}) G_{I,eff} X_I}{G_{tot}^{NC}} \\ t_c^D &= \tau_D \cdot \frac{(1-D)\tilde{G}_P X + G_E^{eff} + (1-D_{int}) G_{I,eff} X_I}{G_{tot}^{NC}} \\ t_c^I &= \hat{\tau}_I \cdot \frac{(1-D)\tilde{G}_P X + G_E^{eff} + G_D^{eff}}{G_{tot}^{NC}}\end{split}\]

NP effect on creep: The nanocomposite creeps less than the unfilled vitrimer because:

  1. Hydrodynamic amplification raises the permanent modulus by \(X(\phi)\)

  2. The slow interphase adds a long retardation time

  3. Interphase modulus contributes to shielding the permanent network

The long-time compliance \(J(\infty) = 1/[(1-D)\tilde{G}_P X]\) is reduced by the factor \(1/X(\phi)\) relative to the unfilled system.

SAOS Derivation

Objective: \(G'(\omega)\) and \(G''(\omega)\) in the linear regime.

Using the relaxation modulus from the relaxation derivation:

Total storage and loss moduli:

\[\boxed{G'(\omega) = (1-D)\tilde{G}_P X + G_E^{eff} \frac{\omega^2 \hat{\tau}_m^2}{1 + \omega^2 \hat{\tau}_m^2} + G_D^{eff} \frac{\omega^2 \tau_D^2}{1 + \omega^2 \tau_D^2} + (1-D_{int}) G_{I,eff} X_I \frac{\omega^2 \hat{\tau}_I^2}{1 + \omega^2 \hat{\tau}_I^2}}\]
\[\boxed{G''(\omega) = G_E^{eff} \frac{\omega \hat{\tau}_m}{1 + \omega^2 \hat{\tau}_m^2} + G_D^{eff} \frac{\omega \tau_D}{1 + \omega^2 \tau_D^2} + (1-D_{int}) G_{I,eff} X_I \frac{\omega \hat{\tau}_I}{1 + \omega^2 \hat{\tau}_I^2}}\]

Key features vs unfilled HVM:

  • Three loss peaks in \(G''(\omega)\) at \(\omega \sim 1/\tau_D\), \(1/\hat{\tau}_m\), \(1/\hat{\tau}_I\) (instead of two)

  • Elevated high-frequency plateau: \(G'(\omega \to \infty) = G_{tot}^{NC} > G_{tot}^{HVM}\)

  • Elevated low-frequency plateau: \(G'(\omega \to 0) = (1-D)\tilde{G}_P X > G_P\)

  • Shifted terminal crossover: Additional slow mode pushes it to lower frequency

Zero-shear viscosity:

\[\eta_0^{NC} = G_E^{eff} \hat{\tau}_m + G_D^{eff} \tau_D + (1-D_{int}) G_{I,eff} X_I \hat{\tau}_I\]

The interphase BER contributes significantly because \(\hat{\tau}_I \gg \hat{\tau}_m\).

Temperature master curves: Below \(T_v^{int}\), the interphase freezes and \(G'(\omega)\) shows a secondary plateau that does not shift with temperature – the failure of simple TTS (thermorheological complexity) is a diagnostic for dual-kinetics.

LAOS Derivation

Objective: Stress response under \(\gamma(t) = \gamma_0 \sin(\omega t)\).

Linear Regime (Constant Rates)

With constant \(k_m\), \(k_I\), \(k_D\), each subnetwork’s shear stress is exactly sinusoidal (linear ODE with sinusoidal forcing). The total \(\sigma_{12}\) is sinusoidal – the model predicts purely linear LAOS.

Nonlinearity Sources

Four mechanisms generate higher harmonics in the HVNM:

  1. TST stress-activated BER (both \(k_{BER}^{mat}\) and \(k_{BER}^{int}\)): Rate increases at peak stress, producing odd harmonics (\(3\omega, 5\omega, \ldots\)).

  2. Strain-amplified interphase (\(X_I > 1\)): The interphase experiences amplified strain, reaching nonlinear territory before the bulk matrix. For \(\gamma_0 X_I > \gamma_{NL}\) but \(\gamma_0 < \gamma_{NL}\), the interphase is nonlinear while the matrix is still linear – intracycle strain softening.

  3. Interfacial damage evolution (\(\dot{D}_{int} > 0\)): At large \(\gamma_0\), interfacial debonding occurs cyclically. If healing rate \(h_{int}\) is slow compared to \(\omega\), \(D_{int}\) ratchets up – progressive intercycle softening (Payne + Mullins in LAOS).

  4. Force-dependent dissociation (\(k_d^D(\boldsymbol{\sigma}^D)\)): Same as HVM.

Critical LAOS strain for nonlinearity onset:

\[\boxed{\gamma_c^{NC} = \frac{\gamma_c^{mat}}{X_I}}\]

The nanocomposite becomes nonlinear at \(1/X_I\) times the strain of the unfilled system.

LAOS descriptors: Fourier decomposition \(\sigma_{12}(t) = \sum_n [G_n' \sin(n\omega t) + G_n'' \cos(n\omega t)]\) for odd \(n\). Third harmonic ratio \(I_{3/1} = |G_3|/|G_1|\) quantifies nonlinearity. The HVNM predicts \(I_{3/1}\) onset at lower \(\gamma_0\) than HVM due to strain amplification.

SAOS Validity Criterion

The analytical SAOS uses constant \(k_{BER}\) rates. For HVNM, both matrix and interphase networks have TST coupling, so the validity depends on two mechanochemical coupling numbers:

\[\Lambda_{mat} = \frac{V_{act} \, G_E}{RT}, \qquad \Lambda_{int} = \frac{V_{act}^{int} \, G_{I,eff} \, X_I}{RT}\]

The critical strain amplitude is set by the more restrictive network:

\[\gamma_c = \min\!\left( \frac{0.14}{\Lambda_{mat}},\; \frac{0.14}{\Lambda_{int}} \right)\]

Nanocomposites typically have \(\Lambda_{int} > \Lambda_{mat}\) (the interphase is stiffer and more strongly coupled), so the interphase network usually limits the SAOS validity. Combined with strain amplification \(X_I\), the effective SAOS window narrows as \(\phi\) increases.

Note

Use model.check_saos_validity(gamma_0) to verify the analytical SAOS is appropriate. The method returns \(\Lambda_{mat}\), \(\Lambda_{int}\), both critical strains, and a validity flag.

Cyclic Loading & Mullins Effect

Objective: Stress-strain response under cyclic loading to maximum strain \(\gamma_{max}\) followed by unloading to zero stress.

First Loading

All four subnetworks respond in parallel. The stress-strain curve follows the startup response (Startup Shear Derivation) truncated at \(\gamma_{max}\).

First Unloading (Hysteresis)

Elastic unloading from all networks. The unloading path lies below the loading path, enclosing a hysteresis loop. The dissipated energy per cycle:

\[W_{diss} = \oint \sigma\,d\gamma = W_{diss}^E + W_{diss}^D + W_{diss}^I + W_{diss}^{dam}\]

where \(W_{diss}^E, W_{diss}^I\) come from BER exchange, \(W_{diss}^D\) from physical bond turnover, and \(W_{diss}^{dam}\) from damage.

Second Loading (Mullins Effect)

The reloading curve lies below the first because:

  1. \(D\) may have increased (permanent chain scission), reducing \(\tilde{G}_P X\)

  2. \(D_{int}\) may have increased (interfacial debonding), reducing \((1-D_{int}) G_{I,eff} X_I\)

  3. \(\boldsymbol{\mu}^E_{nat}\) and \(\boldsymbol{\mu}^I_{nat}\) have evolved toward the deformed state, reducing stored elastic energy

Recovery Between Cycles

If a rest period \(t_{rest}\) is allowed at zero strain:

  • Matrix BER: \(\boldsymbol{\mu}^E, \boldsymbol{\mu}^E_{nat} \to \mathbf{I}\) on timescale \(\hat{\tau}_m\) (full recovery)

  • Interfacial BER: \(\boldsymbol{\mu}^I, \boldsymbol{\mu}^I_{nat} \to \mathbf{I}\) on timescale \(\hat{\tau}_I\) (slower recovery)

  • Interfacial healing: \(D_{int}\) decreases on timescale \(1/h_{int}\) (if \(T > T_v^{int}\))

  • Permanent damage: \(D\) is irreversible, never recovers

Key prediction

The Mullins effect in vitrimer nanocomposites is partially recoverable with time and temperature, unlike conventional filled rubbers where it is permanent. The recovery fraction increases with temperature and rest time, with interphase healing as the rate-limiting step.

Temperature dependence:

  • Above \(T_v^{int}\): softening partially recovered between cycles (self-healing active)

  • Below \(T_v^{int}\): softening is permanent (no interfacial healing)

  • Below \(T_v^{mat}\): all exchange frozen, behavior like filled thermoset

Protocol Comparison: HVNM vs HVM

Protocol

HVM (P + E + D)

HVNM addition (I-network)

Flow curve

\(\sigma_E = 0\) at SS

\(\sigma_I = 0\) at SS; intermediate quasi-steady

SAOS

Two Maxwell modes + \(G_P\)

Third mode + amplified \(G_P X\); three \(G''\) peaks

Relaxation

Bi-exponential + plateau

Quad-exponential + amplified plateau

Startup

TST overshoot

Double overshoot (matrix + interphase)

Creep

Two retardation modes

Three retardation modes; reduced \(J(\infty)\) by \(1/X\)

LAOS

\(\gamma_c^{mat}\)

\(\gamma_c^{NC} = \gamma_c^{mat}/X_I\) (earlier onset)

Cyclic

Not applicable

Partially recoverable Mullins + Payne

For HVM protocol derivations, see HVM Protocol Equations & Derivations.