Equation of State

Energy Equation and Hugoniot Curve

The equation for conservation of energy equates the increase in internal energy per unit mass, $E_{m}$ , to the rate at which work is being done by the stresses and the rate at which heat is being added. In the absence of heat conduction the energy equation can be written as

ρ \frac{\partial E_{m}}{\partial t} = (p + p_{b v}) \frac{1}{ρ} \frac{\partial ρ}{\partial t} + S : \dot{e} + ρ \dot{Q},

where p is the pressure stress defined as positive in compression, $p_{b v}$ is the pressure stress due to the bulk viscosity, $S$ is the deviatoric stress tensor, $\dot{e}$ is the deviatoric part of strain rate, and $\dot{Q}$ is the heat rate per unit mass.

The equation of state is assumed for the pressure as a function of the current density, $ρ$ , and the internal energy per unit mass, $E_{m}$ :

p = f (ρ, E_{m}),

which defines all the equilibrium states that can exist in a material. The internal energy can be eliminated from the above equation to obtain a p versus V relationship (where V is the current volume) or, equivalently, a p versus $1 / ρ$ relationship that is unique to the material described by the equation of state model. This unique relationship is called the Hugoniot curve and is the locus of p–V states achievable behind a shock (see Figure 1).

The Hugoniot pressure, $p_{H}$ , is a function of density only and can be defined, in general, from fitting experimental data.

An equation of state is said to be linear in energy when it can be written in the form

p = f + g E_{m},

where $f (ρ)$ and $g (ρ)$ are functions of density only and depend on the particular equation of state model.

Mie-Grüneisen Equations of State

A Mie-Grüneisen equation of state is linear in energy. The most common form is

p - p_{H} = Γ ρ (E_{m} - E_{H}),

where $p_{H}$ and $E_{H}$ are the Hugoniot pressure and specific energy (per unit mass) and are functions of density only, and $Γ$ is the Grüneisen ratio defined as

Γ = Γ_{0} \frac{ρ_{0}}{ρ},

where $Γ_{0}$ is a material constant and $ρ_{0}$ is the reference density.

The Hugoniot energy, $E_{H}$ , is related to the Hugoniot pressure by

E_{H} = \frac{p_{H} η}{2 ρ_{0}},

where $η = 1 - ρ_{0} / ρ$ is the nominal volumetric compressive strain. Elimination of $Γ$ and $E_{H}$ from the above equations yields

p = p_{H} (1 - \frac{Γ_{0} η}{2}) + Γ_{0} ρ_{0} E_{m} .

The equation of state and the energy equation represent coupled equations for pressure and internal energy. Abaqus/Explicit solves these equations simultaneously at each material point.

Linear U_s − U_p Hugoniot Form

A common fit to the Hugoniot data is given by

p_{H} = \frac{ρ_{0} c_{0}^{2} η}{{(1 - s η)}^{2}},

where $c_{0}$ is the bulk speed of sound and s is the slope of the linear $U_{s} - U_{p}$ Hugoniot form of the equation of state. Material parameters $c_{0}$ and s define the linear relationship between the shock velocity, $U_{s}$ , and the particle velocity, $U_{p}$ , as follows:

U_{s} = c_{0} + s U_{p} .

With the above assumptions the linear $U_{s} - U_{p}$ Hugoniot form is written as

p = \frac{ρ_{0} c_{0}^{2} η}{{(1 - s η)}^{2}} (1 - \frac{Γ_{0} η}{2}) + Γ_{0} ρ_{0} E_{m},

where $ρ_{0} c_{0}^{2}$ is equivalent to the elastic bulk modulus at small nominal strains.

There is a limiting compression given by the denominator of this form of the equation of state

η_{l i m} = \frac{1}{s}

or

ρ_{l i m} = \frac{s ρ_{0}}{s - 1} .

At this limit there is a tensile minimum; thereafter, negative sound speeds are calculated for the material.

Input File Usage

Use both of the following options:

DENSITY (to specify the reference density  $ρ_{0}$ )
EOS, TYPE=USUP (to specify the variables  $c_{0}$ , s, and  $Γ_{0}$ )

Abaqus/CAE Usage

Property module: material editor:
GeneralDensity (to specify the reference density  $ρ_{0}$ )
MechanicalEos: Type: Us - Up (to specify the variables  $c_{0}$ , s, and  $Γ_{0}$ )

Initial State

The initial state of the material is determined by the initial values of specific energy, $E_{m}$ , and pressure stress, p. Abaqus/Explicit will automatically compute the initial density, $ρ$ , that satisfies the equation of state, $p = f (ρ, E_{m})$ . You can define the initial specific energy and initial stress state (see Initial Conditions). The initial pressure used by the equation of state is inferred from the specified stress states. If no initial conditions are specified, Abaqus/Explicit will assume that the material is at its reference state:

\begin{matrix} E_{m} & ​ = 0, \\ p & ​ = 0, \\ ρ & ​ = ρ_{0} . \end{matrix}

Input File Usage

Use either or both of the following options, as required:

INITIAL CONDITIONS, TYPE=SPECIFIC ENERGY
INITIAL CONDITIONS, TYPE=STRESS

Abaqus/CAE Usage

Load module: Create Predefined Field: Step: Initial: choose Mechanical for the Category and Stress for the Types for Selected Step

Initial specific energy is not supported in Abaqus/CAE.

Tabulated Equation of State

The tabulated equation of state provides flexibility in modeling the hydrodynamic response of materials that exhibit sharp transitions in the pressure-density relationship, such as those induced by phase transformations. The tabulated equation of state is linear in energy and assumes the form

p = f_{1} (ε_{vol}) + ρ_{0} f_{2} (ε_{vol}) E_{m},

where $f_{1} (ε_{vol})$ and $f_{2} (ε_{vol})$ are functions of the logarithmic volumetric strain $ε_{vol}$ only, with $ε_{vol} = \ln (ρ_{0} / ρ)$ , and $ρ_{0}$ is the reference density.

You can specify the functions $f_{1} (ε_{vol})$ and $f_{2} (ε_{vol})$ directly in tabular form. The tabular entries must be given in descending values of the volumetric strain (that is, from the most tensile to the most compressive states). Abaqus/Explicit will use a piecewise linear relationship between data points. Outside the range of specified values of volumetric strains, the functions are extrapolated based on the last slope computed from the data.

Input File Usage

Use both of the following options:

DENSITY (to specify the reference density  $ρ_{0}$ )
EOS, TYPE=TABULAR (to specify  $f_{1}$  and  $f_{2}$  as functions of  $ε_{vol}$ )

Abaqus/CAE Usage

Property module: material editor:
GeneralDensity (to specify the reference density  $ρ_{0}$ )
MechanicalEos: Type: Tabular (to specify  $f_{1}$  and  $f_{2}$  as functions of  $ε_{vol}$ )

Initial State

The initial state of the material is determined by the initial values of specific energy, $E_{m}$ , and pressure stress, p. Abaqus/Explicit automatically computes the initial density, $ρ$ , that satisfies the equation of state. You can define the initial specific energy and initial stress state (see Initial Conditions). The initial pressure used by the equation of state is inferred from the specified stress states. If no initial conditions are specified, Abaqus/Explicit assumes that the material is at its reference state:

\begin{matrix} E_{m} & ​ = 0, \\ p & ​ = 0, \\ ρ & = ρ_{0} (ε_{vol} = 0) . \end{matrix}

Input File Usage

Use either or both of the following options, as required:

INITIAL CONDITIONS, TYPE=SPECIFIC ENERGY
INITIAL CONDITIONS, TYPE=STRESS

Abaqus/CAE Usage

Load module: Create Predefined Field: Step: Initial: choose Mechanical for the Category and Stress for the Types for Selected Step

Initial specific energy is not supported in Abaqus/CAE.

User-Defined Equation of State

The user-defined equation of state provides a general capability for modeling the volumetric response of materials through user subroutine VUEOS (see VUEOS). The equation of state defines the pressure as a function of the current density, $ρ$ , and the internal energy per unit mass, $E_{m}$ : $p = f (ρ, E_{m})$ . Abaqus/Explicit solves the energy equation together with the equation of state using an iterative method. The pressure stress, $p$ , and the derivatives of the pressure with respect to the internal energy and to the density, $\partial p / \partial E_{m}$ and $\partial p / \partial ρ$ , must be provided by user subroutine VUEOS. The latter is needed for the evaluation of the effective bulk modulus of the material, which is necessary for the stable time increment calculation.

Optionally, you can also specify the number of property values needed as data in the user subroutine as well as the number of solution-dependent variables (see About User Subroutines and Utilities).

Input File Usage

Use the following option:

EOS, TYPE=USER, PROPERTIES=n

Abaqus/CAE Usage

The user-defined equation of state is not supported in Abaqus/CAE.

Initial State

You need to make sure that the initial specific energy, the initial stress, and the initial density satisfy the equation of state. If you do not specify the initial conditions, Abaqus/Explicit assumes that the material is at its reference state:

\begin{matrix} E_{m} & ​ = 0, \\ p & ​ = 0, \\ ρ & ​ = ρ_{0} . \end{matrix}

Input File Usage

Use either or both of the following options to define the initial specific energy and/or initial pressure stress:

INITIAL CONDITIONS, TYPE=SPECIFIC ENERGY
INITIAL CONDITIONS, TYPE=STRESS

Use the following option to define the initial density:

DENSITY

Abaqus/CAE Usage

Load module: Create Predefined Field: Step: Initial: choose Mechanical for the Category and Stress for the Types for Selected Step

Initial specific energy is not supported in Abaqus/CAE.

P–alpha Equation of State

The $P - α$ equation of state is designed for modeling the compaction of ductile porous materials. The implementation in Abaqus/Explicit is based on the model proposed by Hermann (1968) and Carroll and Holt (1972). The constitutive model provides a detailed description of the irreversible compaction behavior at low stresses and predicts the correct thermodynamic behavior at high pressures for the fully compacted solid material. In Abaqus/Explicit the solid phase is assumed to be governed by either the Mie-Grüneisen equation of state or the tabulated equation of state. The relevant properties of the porous material in the virgin state, to be discussed later, and the material properties of the solid phase are specified separately.

The porosity of the material, n, is defined as the ratio of pore volume, $V_{p}$ , to total volume, $V = V_{s} + V_{p}$ , where $V_{s}$ is the solid volume. The porosity remains in the range $0 \leq n < 1$ , with 0 indicating full compaction. It is convenient to introduce a scalar variable $α$ , sometimes referred to as “distension,” defined as the ratio of the density of the solid material, $ρ_{s}$ , to the density of the porous material, $ρ$ , both evaluated at the same temperature and pressure:

α = \frac{ρ_{s}}{ρ} \geq 1.

For a fully compacted material $α = 1$ ; otherwise, $α$ is greater than 1. Assuming that the density of the pores is negligible compared to that of the solid phase, $α$ can be expressed in terms of the porosity n as

α = \frac{ρ_{s}}{ρ} = \frac{V}{V_{s}} = \frac{V}{V - V_{p}} = \frac{1}{1 - V_{p} / V} = \frac{1}{1 - n} .

An equation of state is assumed for the pressure of the porous material as a function of $α$ ; current density, $ρ$ ; and internal energy per unit mass, $E_{m}$ , in the form

p = p (α, ρ, E_{m}) .

Assuming that the pores carry no pressure, it follows from equilibrium considerations that when a pressure p is applied to the porous material, it gives rise to a volume-average pressure in the solid phase equal to $p_{s} = α p$ . Assuming that the specific internal energies of the porous material and the solid matrix are the same (that is, neglecting the surface energy of the pores), the equation of state of the porous material can be expressed as

p (α, ρ, E_{m}) = \frac{1}{α} p_{s} (α ρ, E_{m}) = \frac{1}{α} p_{s} (ρ_{s}, E_{m}),

where

p_{s} (ρ_{s}, E_{m})

is the equation of state of the solid material. For the fully compacted material (that is, when

α = 1

), the

P - α

equation of state reduces to that of the solid phase, therefore predicting the correct thermodynamic behavior at high pressures.

The $P - α$ equation of state must be supplemented by an equation that describes the behavior of $α$ as a function of the thermodynamic state. This equation takes the form

α = A (p, α_{m i n}),

where

α_{m i n}

is a state variable corresponding to the minimum value attained by

α

during plastic (irreversible) compaction of the material. The state variable is initialized to the elastic limit

α_{e}

for a material that is at its virgin state. The specific form of the function

A (p, α_{m i n})

used by Abaqus/Explicit is illustrated in Figure 2 and is discussed next.

The function $A (p, α_{m i n})$ captures the general behavior to be expected in a ductile porous material. The unloaded virgin state corresponds to the value $α_{0} = 1 / (1 - n_{0})$ , where $n_{0}$ is the reference porosity of the material. Initial compression of the porous material is assumed to be elastic. Recall that decreasing porosity corresponds to a reduction in $α$ . As the pressure increases beyond the elastic limit, $p_{e}$ , the pores in the material start to crush, leading to irreversible compaction and permanent (plastic) volume change. Unloading from a partially compacted state follows a new elastic curve that depends on the maximum compaction (or, alternatively, the minimum value of $α$ ) ever attained during the deformation history of the material. The absolute value of the slope of the elastic curve decreases as $α_{m i n}$ decreases, as will be quantified later. The material becomes fully compacted when the pressure reaches the compaction pressure $p_{S}$ ; at that point $α = α_{m i n} = 1$ , a value that is retained forever. The function $A (p, α_{m i n})$ therefore has multiple branches: a plastic branch, $A_{p l} (p)$ , and multiple elastic branches, $A_{e l} (p, α_{m i n})$ , corresponding to elastic unloading from partially compacted states. The appropriate branch of A is selected according to the following rule:

α = A (p, α_{m i n}) = {\begin{matrix} A_{p l} (p) & if & A_{p l} (p) \leq α_{m i n} \\ A_{e l} (p, α_{m i n}) & if & A_{p l} (p) > α_{m i n} \end{matrix}

These expressions can be inverted to solve for p:

p = P (α, α_{m i n}) = {\begin{matrix} P_{p l} (α) & if & α \leq α_{m i n} \\ P_{e l} (α, α_{m i n}) & if & α > α_{m i n} \end{matrix}

The equation for the plastic curve takes the form

A_{p l} (p) = 1 + (α_{e} - 1) {(\frac{p_{S} - p}{p_{S} - p_{e}})}^{2}

or, alternatively,

P_{p l} (α) = p_{S} - (p_{S} - p_{e}) {(\frac{α - 1}{α_{e} - 1})}^{\frac{1}{2}} .

The elastic curve originally proposed by Hermann (1968) is given by the differential equation

\frac{d {\bar{A}}_{e l}}{d p} (α) = \frac{α^{2}}{K_{0}} (1 - \frac{1}{h^{2} (α)}), h (α) = 1 + \frac{(c_{e} - c_{s}) (α - 1)}{c_{s} (α_{0} - 1)},

where $K_{0} = ρ_{s 0} c_{s}^{2}$ is the elastic bulk modulus of the solid material at small nominal strains; $ρ_{s 0}$ is the reference density of the solid; and $c_{s}$ and $c_{e}$ are the reference sound speeds in the solid and virgin (porous) materials, respectively.

If the solid phase is modeled using the Mie-Grüneisen equation of state, $c_{s}$ is given directly by the reference sound speed, $c_{0}$ . On the other hand, if the solid phase is modeled using the tabulated equation of state, $c_{s}$ is computed from the initial bulk modulus and reference density of the solid material, $c_{s} = \sqrt{K_{0} / ρ_{s 0}}$ . In this case the reference density is required to be constant; it cannot be a function of temperature or field variables.

Following Wardlaw et al. (1996), the above equation for the elastic curve in Abaqus/Explicit is simplified and replaced by the linear relations

A_{e l} (p, α_{m i n}) = α_{m i n} + {(p - P_{p l} (α_{m i n})) \frac{d {\bar{A}}_{e l}}{d p} |}_{α = α_{m i n}}

and

P_{e l} (α, α_{m i n}) = P_{p l} (α_{m i n}) + \frac{(α - α_{m i n})}{{\frac{d {\bar{A}}_{e l}}{d p} |}_{α = α_{m i n}}} .

Input File Usage

Use the following option to specify the reference density of the solid phase, $ρ_{s 0}$ :

DENSITY

Use one of the following options to specify additional material properties for the solid phase:

EOS, TYPE=USUP (if the solid phase is modeled using the Mie-Grüneisen equation of state)
EOS, TYPE=TABULAR (if the solid phase is modeled using the tabulated equation of state)

Use the following option to specify the properties of the porous material (the reference sound speed, $c_{e}$ ; the reference porosity, $n_{0}$ ; the elastic limit, $p_{e}$ ; and the compaction pressure, $p_{S}$ ):

EOS COMPACTION

Abaqus/CAE Usage

Property module: material editor:
GeneralDensity (to specify the reference density  $ρ_{0}$ )

Use one of the following options to specify additional material properties for the solid phase:

MechanicalEos: Type: Us - Up (if the solid phase is modeled using the Mie-Grüneisen equation of state)
MechanicalEos: Type: Tabular (if the solid phase is modeled using the tabulated equation of state)

Use the following option to specify the properties of the porous material:

MechanicalEos: Suboptions Eos Compaction (to specify the reference sound speed,  $c_{e}$ ; the porosity of the unloaded material,  $n_{0}$ ; the pressure required to initialize plastic behavior,  $p_{e}$ ; and the pressure at which all pores are crushed,  $p_{S}$ )

Initial State

The initial state of the porous material is determined from the initial values of porosity, $n = (α - 1) / α$ ; specific energy, $E_{m}$ ; and pressure stress, p. Abaqus/Explicit automatically computes the initial density, $ρ$ , that satisfies the equation of state, $p = f (α, ρ, E_{m})$ . You can define the initial porosity, initial specific energy, and initial stress state (see Initial Conditions). If no initial conditions are given, Abaqus/Explicit assumes that the material is at its virgin state:

\begin{matrix} E_{m} & ​ = 0, \\ p & ​ = 0, \\ α & = α_{0}, (n = n_{0}), \\ ρ & ​ = ρ_{s 0} / α_{0} . \end{matrix}

Abaqus/Explicit will issue an error message if the initial $p - α$ state lies outside the region of allowed states (see Figure 2). When initial conditions are specified only for p (or for $α$ ), Abaqus/Explicit will compute $α$ (or p) assuming that the $p - α$ state lies on the primary (monotonic loading) curve.

Input File Usage

Use some or all of the following options, as required:

INITIAL CONDITIONS, TYPE=SPECIFIC ENERGY
INITIAL CONDITIONS, TYPE=STRESS
INITIAL CONDITIONS, TYPE=POROSITY

Abaqus/CAE Usage

Load module: Create Predefined Field: Step: Initial: choose Mechanical for the Category and Stress for the Types for Selected Step

Initial specific energy and initial porosity are not supported in Abaqus/CAE.

JWL High Explosive Equation of State

The Jones-Wilkins-Lee (or JWL) equation of state models the pressure generated by the release of chemical energy in an explosive. This model is implemented in a form referred to as a programmed burn, which means that the reaction and initiation of the explosive is not determined by shock in the material. Instead, the initiation time is determined by a geometric construction using the detonation wave speed and the distance of the material point from the detonation points.

The JWL equation of state can be written in terms of the internal energy per unit mass, $E_{m}$ , as

p = A (1 - \frac{ω ρ}{R_{1} ρ_{0}}) \exp (- R_{1} \frac{ρ_{0}}{ρ}) + B (1 - \frac{ω ρ}{R_{2} ρ_{0}}) \exp (- R_{2} \frac{ρ_{0}}{ρ}) + ω ρ E_{m},

where $A, B, R_{1}, R_{2},$ and $ω$ are user-defined material constants; $ρ_{0}$ is the user-defined density of the explosive; and $ρ$ is the density of the detonation products.

Input File Usage

Use both of the following options:

DENSITY (to specify the density of the explosive  $ρ_{0}$ )
EOS, TYPE=JWL (to specify the material constants  $A, B, R_{1}, R_{2},$  and  $ω$ )

Abaqus/CAE Usage

Property module: material editor:
GeneralDensity (to specify the density of the explosive  $ρ_{0}$ )
MechanicalEos: Type: JWL (to specify the material constants  $A, B, R_{1},$ 
 $R_{2},$  and  $ω$ )

Arrival Time of Detonation Wave

Abaqus/Explicit calculates the arrival time of the detonation wave at a material point $(t_{d}^{m p})$ as the distance from the material point to the nearest detonation point divided by the detonation wave speed:

t_{d}^{m p} = min [t_{d}^{N} + \sqrt{(x^{m p} - x_{d}^{N}) \cdot (x^{m p} - x_{d}^{N})} / C_{d}],

where $x^{m p}$ is the position of the material point, $x_{d}^{N}$ is the position of the Nth detonation point, $t_{d}^{N}$ is the detonation delay time of the Nth detonation point, and $C_{d}$ is the detonation wave speed of the explosive material. The minimum in the above formula is over the N detonation points that apply to the material point.

Burn Fraction

To spread the burn wave over several elements, a burn fraction, $F_{b}$ , is computed as

F_{b} = min [1, \frac{(t - t_{d}^{m p}) C_{d}}{B_{s} l_{e}}],

where $B_{s}$ is a constant that controls the width of the burn wave (set to a value of 2.5) and $l_{e}$ is the characteristic length of the element. If the time is less than $t_{d}^{m p}$ , the pressure is zero in the explosive; otherwise, the pressure is given by the product of $F_{b}$ and the pressure determined from the JWL equation above.

Defining Detonation Points

You can define any number of detonation points for the explosive material. Coordinates of the points must be defined along with a detonation delay time. Each material point responds to the first detonation point that it sees. The detonation arrival time at a material point is based on the time that it takes a detonation wave (traveling at the detonation wave speed $C_{d}$ ) to reach the material point plus the detonation delay time for the detonation point. If there are multiple detonation points, the arrival time is based on the minimum arrival time for all the detonation points. In a body with curved surfaces care should be taken that the detonation arrival times are meaningful. The detonation arrival times are based on the straight line of sight from the material point to the detonation point. In a curved body the line of sight might pass outside of the body.

Input File Usage

Use both of the following options to define the detonation points:

EOS, TYPE=JWL
DETONATION POINT

Abaqus/CAE Usage

Property module: material editor: MechanicalEos: Type: JWL: 
SuboptionsDetonation Point

Initial State

Explosive materials generally have some nominal volumetric stiffness before detonation. It might be useful to incorporate this stiffness when elements modeled with a JWL equation of state are subjected to stress before initiation of detonation by the arriving detonation wave. You can define the pre-detonation bulk modulus, $K_{p d}$ . The pressure will be computed from the volumetric strain and $K_{p d}$ until detonation, at which time the pressure will be determined by the procedure outlined above. The initial relative density ( $ρ / ρ_{0}$ ) used in the JWL equation is assumed to be unity. The initial specific energy $E_{m_{0}}$ is assumed to be equal to the user-defined detonation energy $E_{0}$ .

If you specify a nonzero value of $K_{p d}$ , you can also define an initial stress state for the explosive materials.

Input File Usage

Use the following option to define the initial stress:

INITIAL CONDITIONS, TYPE=STRESS

Optionally, you can also define the initial specific energy directly:

INITIAL CONDITIONS, TYPE=SPECIFIC ENERGY

Abaqus/CAE Usage

Load module: Create Predefined Field: Step: Initial: choose Mechanical for the Category and Stress for the Types for Selected Step

Initial specific energy is not supported in Abaqus/CAE.

Ignition and Growth Equation of State

The ignition and growth equation of state models shock initiation and detonation wave propagation of solid high explosives. The heterogeneous explosive is modeled as a homogeneous mixture of two phases: the unreacted solid explosive and the reacted gas products. Separate JWL equations of state are prescribed for each phase:

p_{s} = {\tilde{F}}_{1 s} (ρ_{s}) - {\tilde{F}}_{1 s} (ρ_{0}) + {\tilde{F}}_{2 s} (ρ_{s}) E_{m s},

p_{g} = {\tilde{F}}_{1 g} (ρ_{g}) + {\tilde{F}}_{2 g} (ρ_{g}) (E_{m g} + E_{d}),

where

{\tilde{F}}_{1 i} (ρ_{i}) = A_{i} (1 - \frac{ω_{i} ρ_{i}}{R_{1 i} ρ_{0}}) \exp (- R_{1 i} \frac{ρ_{0}}{ρ_{i}}) + B_{i} (1 - \frac{ω_{i} ρ_{i}}{R_{2 i} ρ_{0}}) \exp (- R_{2 i} \frac{ρ_{0}}{ρ_{i}})

and

{\tilde{F}}_{2 i} (ρ_{i}) = ω_{i} ρ_{i}, (i = s, g) .

The subscript s refers to the unreacted solid explosive, and g refers to the reacted gas products. $A_{i}, B_{i}, R_{1 i}, R_{2 i},$ and $ω_{i}$ are user-defined material constants used in the JWL equations; $E_{d}$ is the detonation energy; $ρ_{0}$ is the user-defined reference density of the explosive, and $ρ_{i}$ is the density of the unreacted explosive or the reacted products.

Input File Usage

Use both of the following options:

DENSITY(to specify the density of the explosive  $ρ_{0}$ )
EOS, TYPE=IGNITION AND GROWTH, DETONATION ENERGY= $E_{d}$ 
(to specify the material constants  $A, B, R_{1},$   $R_{2},$  and  $ω$ 
 of the unreacted solid explosive and the reacted gas product)

Abaqus/CAE Usage

Property module: material editor:
GeneralDensity (to specify the density of the explosive  $ρ_{0}$ )
MechanicalEos: Type: Ignition and growth: Detonation energy:  $E_{d}$ ;
Solid Phase tabbed page and Gas Phase tabbed page
(to specify the material constants  $A, B, R_{1},$   $R_{2},$  and  $ω$ 
 of the unreacted solid explosive and the reacted gas product)

Mass Fraction

The mixture of unreacted solid explosive and reacted gas products is defined by the mass fraction

F_{i} = \frac{m_{i}}{(m_{s} + m_{g})},

where $m_{s}$ is the mass of the unreacted explosive, and $m_{g}$ is the mass of the reacted products. It is assumed that the two phases are in thermomechanical equilibrium:

p_{s} = p_{g} and T_{s} = T_{g} .

It is also assumed that the volumes are additive:

V = V_{s} + V_{g} or \frac{1}{ρ} = (1 - F) \frac{1}{ρ_{s}} + F \frac{1}{ρ_{g}} .

Similarly, the internal energy is assumed to be additive:

(m_{s} + m_{g}) E_{m} = m_{s} E_{m s} + m_{g} E_{m g},

where

E_{m} = E_{m_{0}} + \int_{θ_{0} - θ^{Z}}^{θ - θ^{Z}} c_{v} (T) d T .

Hence, the specific heat of the mixture is given by

c_{v} = (1 - F) c_{v s} + F c_{v g} .

Input File Usage

Use the following options to define the specific heat of the unreacted solid explosive:

EOS, TYPE=IGNITION AND GROWTH
SPECIFIC HEAT, DEPENDENCIES=n

Use the following options to define the specific heat of the reacted gas product:

EOS, TYPE=IGNITION AND GROWTH
GAS SPECIFIC HEAT, DEPENDENCIES=n

Abaqus/CAE Usage

Use the following options to define the specific heat of the unreacted solid explosive:

Property module: material editor:
MechanicalEos: Type: Ignition and GrowthThermalSpecific Heat

Use the following options to define the specific heat of the reacted gas product:

Property module: material editor:
MechanicalEos: Type: Ignition and growth:
Gas Specific tabbed page: Specific Heat

You can toggle on Use temperature-dependent data to define the specific heat as a function of temperature or select the Number of field variables to define the specific heat as a function of field variables.

Reaction Rate

The conversion of unreacted solid explosive to reacted gas products is governed by the reaction rate. The reaction rate equation in the ignition and growth model is a pressure-driven rule, which includes three terms:

\frac{d F}{d t} = {\dot{F}}_{i g} + {\dot{F}}_{G_{1}} + {\dot{F}}_{G_{2}} .

These three terms are defined as follows:

\begin{matrix} {\dot{F}}_{i g} & , = I {(1 - F)}^{b} {(\frac{ρ_{s}}{ρ_{0}} - 1 - a)}^{x}, \\ {\dot{F}}_{G_{1}} & , = G_{1} {(1 - F)}^{c} F^{d} {(\frac{p}{p_{r e f}})}^{y}, \\ {\dot{F}}_{G_{2}} & , = G_{2} {(1 - F)}^{e} F^{g} {(\frac{p}{p_{r e f}})}^{z}, \end{matrix}

where $I, G_{1}, G_{2}, a, b, c, d, e, g, x, y$ , and $z$ are reaction rate constants; and $p_{r e f}$ is a reference pressure.

The first term, ${\dot{F}}_{i g}$ , describes hot spot ignition by igniting some of the material relatively quickly but limiting it to a small proportion of the total solid $F_{i g}^{m a x}$ . The second term, ${\dot{F}}_{G_{1}}$ , represents the growth of reaction from the hot spot sites into the material and describes the inward and outward grain burning phenomena; this term is limited to a proportion of the total solid $F_{G_{1}}^{m a x}$ . The third term, ${\dot{F}}_{G_{2}}$ , is used to describe the rapid transition to detonation observed in some energetic materials.

\begin{matrix} {\dot{F}}_{i g} = 0 & if F \geq F_{i g}^{m a x} \\ {\dot{F}}_{G_{1}} = 0 & if F \geq F_{G_{1}}^{m a x} \\ {\dot{F}}_{G_{2}} = 0 & if F \leq F_{G_{2}}^{m i n} \end{matrix}

Input File Usage

Use both of the following options to define the reaction rate:

EOS, TYPE=IGNITION AND GROWTH
REACTION RATE

Abaqus/CAE Usage

Property module: material editor:
MechanicalEos: Type: Ignition and growth:
Reaction Rate tabbed page

Initial State

The initial mass fraction of the unreacted solid explosive is assumed to be one. The initial relative density ( $ρ / ρ_{0}$ ) used in the ignition and growth equation is assumed to be unity. The initial specific energy can be defined for the unreacted explosive.

Input File Usage

Use the following option to define the initial specific energy:

INITIAL CONDITIONS, TYPE=SPECIFIC ENERGY

Abaqus/CAE Usage

Initial specific energy is not supported in Abaqus/CAE.

Ideal Gas Equation of State

An ideal gas equation of state can be written in the form of

p + p_{A} = ρ R (θ - θ^{Z}),

where $p_{A}$ is the ambient pressure, R is the gas constant, $θ$ is the current temperature, and $θ^{Z}$ is the absolute zero on the temperature scale being used. It is an idealization to real gas behavior and can be used to model any gases approximately under appropriate conditions (for example, low pressure and high temperature).

One of the important features of an ideal gas is that its specific energy depends only on its temperature; therefore, the specific energy can be integrated numerically as

E_{m} = E_{m_{0}} + \int_{θ_{0} - θ^{Z}}^{θ - θ^{Z}} c_{v} (T) d T,

where $E_{m_{0}}$ is the initial specific energy at the initial temperature $θ_{0}$ and $c_{v}$ is the specific heat at constant volume (or the constant volume heat capacity), which depends only on temperature for an ideal gas.

Modeling with an ideal gas equation of state is typically performed adiabatically; the temperature increase is calculated directly at the material integration points according to the adiabatic thermal energy increase caused by the work $p d v$ , where v is the specific volume (the volume per unit mass, $v = 1 / ρ$ ). Therefore, unless a fully coupled temperature-displacement analysis is performed, an adiabatic condition is always assumed in Abaqus/Explicit.

When performing a fully coupled temperature-displacement analysis, the pressure stress and specific energy are updated based on the evolving temperature field. The energy increase due to the change in state will be accounted for in the heat equation and will be subject to heat conduction.

For the ideal gas model in Abaqus/Explicit you define the gas constant, R, and the ambient pressure, $p_{A}$ . For an ideal gas R can be determined from the universal gas constant, $\tilde{R}$ , and the molecular weight, $M W$ , as follows:

R = \frac{\tilde{R}}{M W} .

In general, the value R for any gas can be estimated by plotting $p v / (θ - θ^{Z})$ as a function of state (for example, pressure or temperature). The ideal gas approximation is adequate in any region where this value is constant. You must specify the specific heat at constant volume, $c_{v}$ . For an ideal gas $c_{v}$ is related to the specific heat at constant pressure, $c_{p}$ , by

R = c_{p} - c_{v} .

Input File Usage

Use both of the following options:

EOS, TYPE=IDEAL GAS
SPECIFIC HEAT, DEPENDENCIES=n

Abaqus/CAE Usage

Property module: material editor:
MechanicalEos: Type: Ideal Gas
ThermalSpecific Heat

Initial State

There are different methods to define the initial state of the gas. You can specify the initial density, $ρ$ , and either the initial pressure stress, $p_{0}$ , or the initial temperature, $θ_{0}$ . The initial value of the unspecified field (temperature or pressure) is determined from the equation of state. Alternatively, you can specify both the initial pressure stress and the initial temperature. In this case the user-specified initial density is replaced by that derived from the equation of state in terms of initial pressure and temperature.

By default, Abaqus/Explicit automatically computes the initial specific energy, $E_{m_{0}}$ , from the initial temperature by numerically integrating the equation

E_{m_{0}} = \int_{0}^{θ_{0} - θ^{Z}} c_{v} (T) d T .

Optionally, you can override this default behavior by defining the initial specific energy for the ideal gas directly.

Input File Usage

Use some or all of the following options, as required:

DENSITY, DEPENDENCIES=n
INITIAL CONDITIONS, TYPE=STRESS
INITIAL CONDITIONS, TYPE=TEMPERATURE

Use the following option to specify the initial specific energy directly:

INITIAL CONDITIONS, TYPE=SPECIFIC ENERGY

Abaqus/CAE Usage

Property module: material editor: General Density

Load module: Create Predefined Field: Step: Initial: choose Other for the Category and Temperature for the Types for Selected Step

Load module: Create Predefined Field: Step: Initial: choose Mechanical for the Category and Stress for the Types for Selected Step

Initial specific energy is not supported in Abaqus/CAE.

Absolute Zero Value

When a non-absolute temperature scale is used, you must specify the value of absolute zero temperature.

Input File Usage

PHYSICAL CONSTANTS, ABSOLUTE ZERO= $θ^{Z}$

Abaqus/CAE Usage

Any module: ModelEdit Attributesmodel_name: Absolute zero temperature

Special Case

In the case of an adiabatic analysis with constant specific heat (both $c_{v}$ and $c_{p}$ are constant), the specific energy is linear in temperature

E_{m} = c_{v} (θ - θ^{Z}) .

The pressure stress can, therefore, be recast in the common form of

p + p_{A} = (γ - 1) ρ E_{m},

where $γ = c_{p} / c_{v}$ is the ratio of specific heats and can be defined as

γ = \frac{n + 2}{n},

where

$n = 3$: for a monatomic;
$n = 5$: for a diatomic; and
$n = 6$: for a polyatomic gas.

Comparison with the Hydrostatic Fluid Model

The ideal gas equation of state can be used to model wave propagation effects and the dynamics of a spatially varying state of a gaseous region. For cases in which the inertial effects of the gas are not important and the state of the gas can be assumed to be uniform throughout a region, the hydrostatic fluid model (About Surface-Based Fluid Cavities) is a simpler, more computationally efficient alternative.

Deviatoric Behavior

The equation of state defines only the material's hydrostatic behavior. It can be used by itself, in which case the material has only volumetric strength (the material is assumed to have no shear strength). Alternatively, Abaqus/Explicit allows you to define deviatoric behavior, assuming that the deviatoric and volumetric responses are uncoupled. Three options are available for the deviatoric response: a linear isotropic elastic model, a viscous model, and a nonlinear viscoelastic model. The material's volumetric response is governed then by the equation of state model, while its deviatoric response is governed by either the linear isotropic elastic model, the viscous fluid model or the nonlinear viscoelastic model.

Elastic Shear Behavior

For the elastic shear behavior the deviatoric stress is related to the deviatoric strain as

S = 2 μ e^{e l},

where $S$ is the deviatoric stress and $e^{e l}$ is the deviatoric elastic strain. See Defining Isotropic Shear Elasticity for Equations of State in Abaqus/Explicit for more details.

Input File Usage

Use both of the following options to define elastic shear behavior:

EOS
ELASTIC, TYPE=SHEAR

Abaqus/CAE Usage

Property module: material editor: Mechanical Elasticity Elastic ; Type: Shear; Shear Modulus

Viscous Shear Behavior

For the viscous shear behavior the deviatoric stress is related to the deviatoric strain rate as

S = 2 η \dot{e} = η \dot{γ},

where $S$ is the deviatoric stress, $\dot{e}$ is the deviatoric part of the strain rate, $η$ is the viscosity, and $\dot{γ} = 2 \dot{e}$ is the engineering shear strain rate.

Abaqus/Explicit provides a wide range of viscosity models to describe both Newtonian and non-Newtonian fluids. These are described in Viscosity.

Input File Usage

Use both of the following options to define viscous shear behavior:

EOS
VISCOSITY

Abaqus/CAE Usage

Property module: material editor: Mechanical Viscosity

Nonlinear Viscoelastic Shear Behavior

The combination of the elastic and viscous shear behaviors results in a nonlinear viscoelastic shear definition (based on a Maxwell viscoelastic model), such that the deviatoric stress and the deviatoric strain satisfy the following equation:

\frac{\dot{S}}{2 μ} + \frac{S}{2 η} = \dot{e},

where

S

is the deviatoric stress,

\dot{e}

is the deviatoric part of the strain rate,

η

is the viscosity, and

μ

is the elastic shear modulus.

In addition to modeling nonlinear viscoelastic shear behavior, you can also use this model to simulate purely viscous fluids with very large viscosity by carefully choosing a penalty value for the elastic shear modulus, $μ$ , to achieve a reasonable stable time increment. You can combine elastic shear behavior with viscous shear behavior only; you cannot add any additional plastic shear behavior to the material definition.

Input File Usage

Use both of the following options to define the combined elastic and viscous shear behavior:

EOS,
ELASTIC, TYPE=SHEAR
VISCOSITY

Abaqus/CAE Usage

Property module: material editor:

Mechanical Elasticity Elastic ; Type: Shear; Shear Modulus

Mechanical Viscosity

Use with the Mises or the Johnson-Cook Plasticity Models

An equation of state model, with the exception of the ignition and growth model, can be used with the Mises (Classical Metal Plasticity) or the Johnson-Cook (Johnson-Cook Plasticity) plasticity models to model elastic-plastic behavior. In this case you must define the elastic part of the shear behavior. The material's volumetric response is governed by the equation of state model, while the deviatoric response is governed by the linear elastic shear and the plasticity model.

Input File Usage

Use the following options:

EOS
ELASTIC, TYPE=SHEAR
PLASTIC

Abaqus/CAE Usage

Property module: material editor:
MechanicalElasticityElastic; Type: Shear
MechanicalPlasticityPlastic

Initial Conditions

You can specify initial conditions for the equivalent plastic strain, ${\bar{ε}}^{p l}$ (Initial Conditions).

Input File Usage

INITIAL CONDITIONS, TYPE=HARDENING

Abaqus/CAE Usage

Load module: Create Predefined Field: Step: Initial, choose Mechanical for the Category and Hardening for the Types for Selected Step

Use with the Extended Drucker-Prager Plasticity Models

An equation of state model, with the exception of the ignition and growth model, can be used in conjunction with the extended Drucker-Prager (Extended Drucker-Prager Models) plasticity models to model pressure-dependent plasticity behavior. This approach can be appropriate for modeling the response of ceramics and other brittle materials under high velocity impact conditions. In this case you must define the elastic part of the shear behavior. The material's deviatoric response is governed by the linear elastic shear and the pressure-dependent plasticity model, while the volumetric response is governed by the equation of state model. In particular, no plastic dilation effects are taken into account (if you specify a dilation angle other than zero, the value is ignored and Abaqus/Explicit issues a warning message).

High-velocity impact of a ceramic target illustrates the use of an equation of state model with the extended Drucker-Prager plasticity models.

Input File Usage

Use the following options:

EOS
ELASTIC, TYPE=SHEAR
DRUCKER PRAGER
DRUCKER PRAGER HARDENING

Abaqus/CAE Usage

Property module: material editor:
MechanicalElasticityElastic; Type: Shear
MechanicalPlasticityDrucker Prager: SuboptionsDrucker Prager Hardening

Initial Conditions

You can specify initial conditions for the equivalent plastic strain, ${\bar{ε}}^{p l}$ (Initial Conditions).

Input File Usage

INITIAL CONDITIONS, TYPE=HARDENING

Abaqus/CAE Usage

Load module: Create Predefined Field: Step: Initial, choose Mechanical for the Category and Hardening for the Types for Selected Step

Use with the Tensile Failure Model

An equation of state model (except the ideal gas equation of state) can also be used with the tensile failure model (Dynamic Failure Models) to model dynamic spall or a pressure cutoff. The tensile failure model uses the hydrostatic pressure stress as a failure measure and offers a number of failure choices. You must provide the hydrostatic cutoff stress.

You can specify that the deviatoric stresses should fail when the tensile failure criterion is met. In the case where the material's deviatoric behavior is not defined, this specification is meaningless and is, therefore, ignored.

The tensile failure model in Abaqus/Explicit is designed for high-strain-rate dynamic problems in which inertia effects are important. Therefore, it should be used only for such situations. Improper use of the tensile failure model may result in an incorrect simulation.

Input File Usage

Use the following options:

EOS
TENSILE FAILURE

Abaqus/CAE Usage

Property module: material editor: Mechanical Eos : Suboptions Tensile Failure

Adiabatic Assumption

An adiabatic condition is always assumed for materials modeled with an equation of state unless a dynamic coupled temperature-displacement procedure is used. The adiabatic condition is assumed irrespective of whether an adiabatic dynamic stress analysis step has been specified. The temperature increase is calculated directly at the material integration points according to the adiabatic thermal energy increase caused by the mechanical work

ρ c_{v} (θ) \frac{\partial θ}{\partial t} = (p - p_{b v}) \frac{1}{ρ} \frac{\partial ρ}{\partial t} + S : \dot{e},

where $c_{v}$ is the specific heat at constant volume. Specifying temperature as a predefined field has no effect on the behavior of this model.

When performing a fully coupled temperature-displacement analysis, the specific energy is updated based on the evolving temperature field using

ρ \frac{\partial E_{m}}{\partial t} = ρ c_{v} (θ) \frac{\partial θ}{\partial t} .

Modeling Fluids

A linear $U_{s} - U_{p}$ equation of state model can be used to model incompressible viscous and inviscid laminar flow governed by the Navier-Stokes equation of motion. The volumetric response is governed by the equations of state, where the bulk modulus acts as a penalty parameter for the incompressible constraint.

To model a viscous laminar flow that follows the Navier-Poisson law of a Newtonian fluid, use the Newtonian viscous deviatoric model and define the viscosity as the real linear viscosity of the fluid. To model non-Newtonian viscous flow, use one of the nonlinear viscosity models available in Abaqus/Explicit. Appropriate initial conditions for velocity and stress are essential to get an accurate solution for this class of problems.

To model an incompressible inviscid fluid such as water in Abaqus/Explicit, it is useful to define a small amount of shear resistance to suppress shear modes that can otherwise tangle the mesh. Here the shear stiffness or shear viscosity acts as a penalty parameter. The shear modulus or viscosity should be small because flow is inviscid; a high shear modulus or viscosity will result in an overly stiff response. To avoid an overly stiff response, the internal forces arising due to the deviatoric response of the material should be kept several orders of magnitude below the forces arising due to the volumetric response. This can be done by choosing an elastic shear modulus that is several orders of magnitude lower than the bulk modulus. If the viscous model is used, the shear viscosity specified should be on the order of the shear modulus, calculated as above, scaled by the stable time increment. The expected stable time increment can be obtained from a data check analysis of the model. This method is a convenient way to approximate a shear resistance that will not introduce excessive viscosity in the material.

If a shear model is defined, the hourglass control forces are calculated based on the shear resistance of the material. Thus, in materials with extremely low or zero shear strengths such as inviscid fluids, the hourglass forces calculated based on the default parameters are insufficient to prevent spurious hourglass modes. Therefore, a sufficiently high hourglass scaling factor is recommended to increase the resistance to such modes.

Elements

In general, the equations of state can be used with any solid (continuum) elements in Abaqus/Explicit except plane stress elements. The ignition and growth equation of state cannot be used with plane strain, plane stress, and axisymmetric elements. For three-dimensional applications exhibiting high confinement, the default kinematic formulation is recommended with reduced-integration solid elements (see Section Controls).

Output

In addition to the standard output identifiers available in Abaqus (Abaqus/Explicit Output Variable Identifiers), the following variables have special meaning for the equation of state models:

PALPH: Distension, $α$ , of the $P - α$ porous material. The current porosity is equal to one minus the inverse of $α$ : $n = 1 - α^{- 1} .$
PALPHMIN: Minimum value, $α_{m i n}$ , of the distension attained during plastic compaction of the $P - α$ porous material.
PEEQ: Equivalent plastic strain, ${\bar{ε}}^{p l} = {{\bar{ε}}^{p l} |}_{0} + \int_{0}^{t} \sqrt{\frac{2}{3} {\dot{ε}}^{p l} : {\dot{ε}}^{p l}} d t,$ where ${{\bar{ε}}^{p l} |}_{0}$ is the initial equivalent plastic strain (zero or user-specified; see Initial Conditions). This is relevant only if the equation of state model is used in combination with the Mises, Johnson-Cook, or extended Drucker-Prager plasticity models.

References

Carroll, M., and A. C. Holt, “Suggested Modification of the $P - α$ Model for Porous Materials,” Journal of Applied Physics, vol. 43, no. 2, pp. 759–761, 1972.
Dobratz, B. M., “LLNL Explosives Handbook, Properties of Chemical Explosives and Explosive Simulants,” UCRL-52997, Lawrence Livermore National Laboratory, Livermore, California, January 1981.
Herrmann, W., “Constitutive Equation for the Dynamic Compaction of Ductile Porous Materials,” Journal of Applied Physics, vol. 40, no. 6, pp. 2490–2499, 1968.
Lee, E., M. Finger, and W. Collins, “JWL Equation of State Coefficients for High Explosives,” UCID-16189, Lawrence Livermore National Laboratory, Livermore, California, January 1973.
Wardlaw, A. B., R. McKeown, and H. Chen, “Implementation and Application of the $P - α$ Equation of State in the DYSMAS Code,” Naval Surface Warfare Center, Dahlgren Division, Report Number: NSWCDD/TR-95/107, May 1996.

Equation of State

Energy Equation and Hugoniot Curve

Mie-Grüneisen Equations of State

Linear Us − Up Hugoniot Form

Initial State

Tabulated Equation of State

Initial State

User-Defined Equation of State

Initial State

P–alpha Equation of State

Initial State

JWL High Explosive Equation of State

Arrival Time of Detonation Wave

Burn Fraction

Defining Detonation Points

Initial State

Ignition and Growth Equation of State

Mass Fraction

Reaction Rate

Initial State

Ideal Gas Equation of State

Initial State

Absolute Zero Value

Special Case

Comparison with the Hydrostatic Fluid Model

Deviatoric Behavior

Elastic Shear Behavior

Viscous Shear Behavior

Nonlinear Viscoelastic Shear Behavior

Use with the Mises or the Johnson-Cook Plasticity Models

Initial Conditions

Use with the Extended Drucker-Prager Plasticity Models

Initial Conditions

Use with the Tensile Failure Model

Adiabatic Assumption

Modeling Fluids

Elements

Output

References

Linear U_s − U_p Hugoniot Form