# Bergstrom-Boyce (BB) Model

## Introduction

The Bergstrom-Boyce (BB) model is an advanced material model that I developed for predicting the time-dependent, large-strain behavior of elastomer-like materials. The model has been shown to be accurate for both traditional engineering rubbers, and soft biomaterials. The version of the BB-model that is discussed here is the PolyUMod version of the BB-model. This model is part of the PolyUMod library, and has Material Model (MM) id=5.

The model uses the following material parameters:

Parameter NameDescription
muShear modulus of network A
lambdaLLocking stretch
kappaBulk modulus
sRelative stiffness of network B
CStrain exponential
tauBaseFlow resistance
mStress exponential
tauCutNormalized cut-off stress for flow

General notes about the material parameters:

• There is no need to search for lambdaL unless the experimental data contains large strains.
• The parameter m should not be larger than 20 in order to get good convergence in a FE simulation.
• Having a very small tauBase compared to stiffness can cause convergence problems.

Here is an exemplar set of stress-strain predictions from the BB-model created by MCalibration. ## BB-model in MCalibration

The material model can be solved using the internal MCalibration solver, or using any of the supported FE solvers. It is recommended to use the MCalibration native solver when possible since it runs significantly faster. If kappa=0, the applied load is uniaxial, and MCalibration is set as the solver, then a fast incompressible solver will be used to calculate the stress-strain response.

## Model Theory

In the Bergstrom-Boyce (BB) model the applied deformation gradient is acting on two parallel macromolecular networks: $$\mathbf{F} = \mathbf{F}_A = \mathbf{F}_B$$, see the following rheological representation. The deformation gradient acting on network B is further decomposed into elastic and viscoelastic components: $$\mathbf{F}_B = \mathbf{F}_B^e \mathbf{F}_B^v$$. The response of network A is given by the eight-chain model:
$\boldsymbol{\sigma}_A = \frac{\,\mu\,}{J \overline{\lambda^*}} \, \frac{\mathcal{L}^{-1}\!\left(\overline{\lambda^*}/\lambda_L\right)}{\mathcal{L}^{-1}\!\left(1/\lambda_L\right)} \, dev[\mathbf{b}^*] + \kappa (J-1) \mathbf{I},$
The stress on network B is also given by the eight-chain model, but with a different effective shear modulus:
$\boldsymbol{\sigma}_B = \frac{s\,\mu\,}{J_B^e \overline{\lambda^{e*}_B}} \, \frac{\mathcal{L}^{-1}\!\left(\overline{\lambda^{e*}_B}/\lambda_L\right)}{\mathcal{L}^{-1}\!\left(1/\lambda_L\right)} \, dev[\mathbf{b}^{e*}_B] + \kappa (J_B^e-1) \mathbf{I},$
where s is a dimensionless material parameter specifying the shear modulus of network B relative to network A, and $$\overline{\lambda_B^{e*}}$$ is the chain stretch in the elastic part of Network B. Using this representation the total Cauchy stress is given by
$\boldsymbol{\sigma} = \boldsymbol{\sigma}_A + \boldsymbol{\sigma}_B.$

The velocity gradient on network B, $$\mathbf{L}_B = \dot{\mathbf{F}}_B \mathbf{F}_{B}^{-1}$$, can be decomposed into elastic and viscous components:
$\mathbf{L}_B = \left[\frac{d}{dt}\left(\mathbf{F}_B^e \mathbf{F}_B^v \right)\right] \left(\mathbf{F}_B^e \mathbf{F}_B^v \right)^{-1} \\ = \left[ \dot{\mathbf{F}}_B^e \mathbf{F}_B^v + \mathbf{F}_B^e \dot{\mathbf{F}}_B^v \right] \left(\mathbf{F}_B^v\right)^{-1} \left(\mathbf{F}_B^e\right)^{-1} \\ = \dot{\mathbf{F}}_B^e \left(\mathbf{F}_B^e\right)^{-1} + \mathbf{F}_B^e \dot{\mathbf{F}}_B^v \left(\mathbf{F}_B^v\right)^{-1} \left(\mathbf{F}_B^e\right)^{-1} \\ = \mathbf{L}_B^e + \mathbf{F}_B^e \mathbf{L}_B^v (\mathbf{F}_B^e)^{-1} \\ = \mathbf{L}_B^e + \tilde{\mathbf{L}}_B^v,$
where
$\mathbf{L}_B^v = \dot{\mathbf{F}}_B^v \left(\mathbf{F}_B^v\right)^{-1} = \mathbf{D}_B^v + \mathbf{W}_B^v, \\ \tilde{\mathbf{L}}_B^v = \tilde{\mathbf{D}}_B^v + \tilde{\mathbf{W}}_B^v.$
To make the unloading unique, prescribe $$\tilde{\mathbf{W}}_B^v \equiv 0$$. The rate of viscous deformation of network B is constitutively prescribed by:
$\tilde{\mathbf{D}}_B^v = \dot{\gamma}_B(\boldsymbol{\sigma}_B,\mathbf{b}_B^{e*}) \, \mathbf{N}_B^v,$
where
$\mathbf{N}_B^v = \frac{dev[\boldsymbol{\sigma}_B]}{\tau} = \frac{dev[\boldsymbol{\sigma}_B]}{|| dev[\boldsymbol{\sigma}]_B ||_F}.$
and $$\tau$$ is the effective stress driving the viscous flow. The time derivative of $$\mathbf{F}_B^v$$ can be derived as follows:
$\tilde{\mathbf{L}}_B^v = \dot{\gamma}_B^v \mathbf{N}_B^v,\\ \Rightarrow\qquad \mathbf{F}_B^e \dot{\mathbf{F}}_B^v \left(\mathbf{F}_B^v\right)^{-1} \left(\mathbf{F}_B^e\right)^{-1} = \dot{\gamma}_B^v \mathbf{N}_B^v,\\ \Rightarrow\qquad \dot{\mathbf{F}}_B^v = \dot{\gamma}_B^v \left(\mathbf{F}_B^e\right)^{-1} \frac{dev[\boldsymbol{\sigma}_B]}{|| dev[\boldsymbol{\sigma}]_B ||_F} \mathbf{F}_B^e \mathbf{F}_B^v.$
The rate-equation for viscous flow is given by:
$\dot{\gamma}_B^v = \dot{\gamma}_0 \left(\overline{\lambda_B^v} – 1 + \xi \right)^C \, \left[ R\left( \frac{\tau}{\tau_{\mathit{base}}} – \hat{\tau}_{\mathit{cut}} \right) \right]^m,$
$$R(x) = (x + |x|) / 2$$ is the ramp function, $$\hat{\tau}_{\mathit{cut}}$$ is a cut-off stress below which no flow will occur, and where $$\dot{\gamma}_0 \equiv 1$$/s is a constant introduced to ensure dimensional consistency,
$\overline{\lambda_B^v} = \sqrt{\frac{tr[\mathbf{b}_B^v]}{3}}.$
is the viscoelastic chain stretch. The effective stress driving the viscous flow is:
$\tau = || dev[\boldsymbol{\sigma}_B] ||_F = \sqrt{ tr \left[ \boldsymbol{\sigma}_B’ \boldsymbol{\sigma}_B’ \right] }.$