.. _failure-anisotropic:

Anisotropic Failure
===================

Background
----------

A common feature of many of the modular failures model available in LAMÉ is damage accumulation through integration by a scalar equivalent plastic strain rate variable, e.g. :eq:`mat:eq:modularfailure_tp`, :eq:`mat:eq:modularfailure_johnsoncook`, :eq:`mat:eq:modularfailure_wilkins`, :eq:`mat:eq:modularfailure_modfail`, :eq:`mat:eq:modularfailure_BCJ_nucleation_HG`, :eq:`mat:eq:modularfailure_BCJ_growth_CA`. This anisotropic failure approach leverages that common form to introduce a heuristic anisotropic dependence.

The plastic strain rate tensor :math:`\dot{\varepsilon}^p_{ij}` is defined through the flow rule deriving from the normality principle of associated plasticity, by scalar plasticity rate variable :math:`\dot{\gamma}`, effective stress :math:`\sigma_{\star}`, and stress tensor :math:`\sigma_{ij}`:

.. math::

   \dot{\varepsilon}^p_{ij} = \dot{\gamma} \frac{\partial \sigma_{\star}}{\partial \sigma_{ij}}

This relation can be expressed equivalently as an equivalent plastic strain rate scalar :math:`\dot{\bar{\varepsilon}}^p` and a flow direction tensor :math:`N`. The LAMÉ implementation defines the flow tensor as the effective stress derivative with respect to the stress tensor, without any further normalization:

.. math::

   \dot{\varepsilon}^p_{ij} = \dot{\bar{\varepsilon}}^p N_{ij}

.. math::
   :label: q:flow_rule

   \dot{\bar{\varepsilon}}^p = \dot{\gamma} ,
   \quad
   N_{ij} = \frac{\partial \sigma_{\star}}{\partial \sigma_{ij}}

.. math::
   :label: eq:eqpsdot

   \dot{\bar{\varepsilon}}^p = \frac{\sqrt{\dot{\varepsilon}^p_{ij}\dot{\varepsilon}^p_{ij}}}{\sqrt{N_{ij}N_{ij}}}

For example, in :math:`J_2` plasticity, :math:`N_{ij} = \frac32 \frac{s_{ij}}{\sigma_{\star}} = \sqrt{\frac32} \frac{s_{ij}}{s_{ij}s_{ij}}`, so :math:`N_{ij}N_{ij} = \frac32`. The plastic strain rate scalar is then the familiar  :math:`\dot{\bar{\varepsilon}}^p =\sqrt{ \frac23 \dot{\varepsilon}^p_{ij} \dot{\varepsilon}^p_{ij}}`.

Anisotropic Weighting Concept
-----------------------------

In order to create anisotropic failure behavior, the plastic strain rate scalar used for failure model integration is weighted using an *anisotropic weight tensor*; the weighted equivalent plastic strain scalar is denoted as :math:`\hat{\dot{\bar{\varepsilon}}}^p`.
Crucially, the weighting is *not* intended to modify the material constitutive (i.e. plasticity, non-failure) model in any way. 

Second-order and fourth-order anisotropic weighting is proposed, using symmetric second-order :math:`A_{ij}` or fourth-order :math:`B_{ijkl}` symmetric weighting tensors [:footcite:`Stershic_AnisotropicFailure_Memo`, :footcite:`Stershic_AnisotropicFailure_USNCTAM`]: 

2nd-order formulation:

.. math::
   :label: mat:eq:anisotropic_Aij

   \hat{\dot{\bar{\varepsilon}}}^p
   = \sqrt{ \frac{  \dot{\varepsilon}^p_{ik} A_{ij} \dot{\varepsilon}^p_{jl} \delta_{kl} } {N_{mn}N_{mn}} }
   = \dot{\bar{\varepsilon}}^p \sqrt{ \frac{ N_{ik} A_{ij} N_{jl} \delta_{kl} } {N_{mn}N_{mn}} }

..    %=  \sqrt{ \frac{ \mbox{tr}(\dot{\varepsilon}^p \cdot A \cdot  \dot{\varepsilon}^p)} {N:N} }
..    %=   \dot{\bar{\varepsilon}}^p  \sqrt{ \frac{ \mbox{tr}(N \cdot A \cdot N)} {N:N} }

4th-order formulation:

.. math::
   :label: mat:eq:anisotropic_Bijkl

   \hat{\dot{\bar{\varepsilon}}}^p
   = \sqrt{ \frac{ \dot{\varepsilon}^p_{ij} B_{ijkl} \dot{\varepsilon}^p_{kl}}  {N_{mn}N_{mn}} }
   = \dot{\bar{\varepsilon}}^p \sqrt{ \frac{ N_{ij} B_{ijkl} N_{kl}}  {N_{mn}N_{mn}} }

..    %= \sqrt{ \frac{ \dot{\varepsilon}^p : B : \dot{\varepsilon}^p } {N:N} }
..    %= \dot{\bar{\varepsilon}}^p \sqrt{ \frac{ N : B : N } {N:N} }

Symmetry
--------

As stress tensors and plastic flow direction tensors are assumed to be symmetric, the anisotropic weight tensors are constrained to have major (:math:`A_{ij}=A_{ji}`; :math:`B_{ijkl}=B_{klij}`) and minor (:math:`B_{ijkl}=B_{jikl}=B_{ijlk}=B_{jilk}`) symmetries. This reduces the number of unique parameters to six for the 2nd-order tensor :math:`A_{ij}` and 21 for the 4th-order tensor :math:`B_{ijkl}`.

Mandel Notation
---------------

It is convenient to express the fourth-order tensor :math:`B_{ijkl}` in Mandel notation as a 6x6 matrix :math:`B_{IJ}`. The conversion between the full forms of the flow tensor (:math:`N_{ij}`) and weight tensor (:math:`B_{ijkl}`) to Mandel notation :math:`\{N_I,B_{IJ}\}` are as follows:

.. math::

   N_I = \left[
   \begin{array}{c}
   N_{11} \\
   N_{22} \\
   N_{33} \\
   \sqrt{2} N_{12} \\
   \sqrt{2} N_{23} \\
   \sqrt{2} N_{31}
   \end{array}
   \right]

.. math::

   B_{IJ} = \left[
   \begin{array}{c c c c c c}
   B_{1111} & B_{1122} & B_{1133} & \sqrt{2} B_{1112} & \sqrt{2} B_{1123} & \sqrt{2} B_{1131} \\
   B_{1122} & B_{2222} & B_{2233} & \sqrt{2} B_{2212} & \sqrt{2} B_{2223} & \sqrt{2} B_{2231} \\
   B_{1133} & B_{2233} & B_{3333} & \sqrt{2} B_{3312} & \sqrt{2} B_{3323} & \sqrt{2} B_{3331} \\
   \sqrt{2} B_{1112} & \sqrt{2} B_{2212} & \sqrt{2} B_{3312} & 2 B_{1212} & 2 B_{1223} & 2 B_{1231} \\
   \sqrt{2} B_{1123} & \sqrt{2} B_{2223} & \sqrt{2} B_{3323} & 2 B_{1223} & 2 B_{2323} & 2 B_{2331} \\
   \sqrt{2} B_{1131} & \sqrt{2} B_{2231} & \sqrt{2} B_{3331} & 2 B_{1231} & 2 B_{2331} & 2 B_{3131} 
   \end{array}
   \right]

The anisotropic weighting can now be expressed in terms of the Mandel-converted tensors.

.. math::
   :label: mat:eq:anisotropic_BIJ

   \hat{\dot{\bar{\varepsilon}}}^p
   = \dot{\bar{\varepsilon}}^p \sqrt{ \frac{ N_{ij} B_{ijkl} N_{kl}}  {N_{mn}N_{mn}} }
   = \dot{\bar{\varepsilon}}^p \sqrt{ \frac{ N_{I} B_{IJ} N_{J}}  {N_{K}N_{K}} }

.. %= \dot{\bar{\varepsilon}}^p \sqrt{ \frac{ N : B : N } {N:N} }
 
Note that the anisotropic failure model parameters are specified as individual entries of the fourth-order tensor - not as entries of the equivalent Mandel matrix.

Weighting Tensor Forms
----------------------

Common forms for the anisotropic weighting tensors are as follows:

*Isotropic, Second-Order*:

.. math::

   A_{ij}^{iso} = \left[
   \begin{array}{c c c}
   1 & 0 & 0 \\
   0 & 1 & 0 \\
   0 & 0 & 1 
   \end{array}
   \right]

*Coordinate-Axis Aligned, Second-Order*:

.. math::

   A_{ij}^{diag} = \left[
   \begin{array}{c c c}
   A_{11} & 0 & 0 \\
   0 & A_{22} & 0 \\
   0 & 0 & A_{33} 
   \end{array}
   \right]

*General (Rotated), Second-Order*:

.. math::

   A_{ij}^{gen} = \left[
   \begin{array}{c c c}
   A_{11} & A_{12} & A_{13} \\
   A_{12} & A_{22} & A_{23} \\
   A_{13} & A_{23} & A_{33} 
   \end{array}
   \right]
   = A_{kl}^{diag} R_{ik} R_{lj}

Note that the rotation matrix can equivalently be applied to the weighting tensor (shown here) or to the flow direction tensor.

*Fourth-Order Equivalence to Second-Order*

.. math::

   B_{IJ} = \left[
   \begin{array}{c c c c c c}
   A_{11} & 0 & 0 & \frac{1}{\sqrt{2}} A_{12} & 0 & \frac{1}{\sqrt{2}} A_{13} \\
   0 & A_{22} & 0 & \frac{1}{\sqrt{2}} A_{12} & \frac{1}{\sqrt{2}} A_{23} & 0 \\
   0 & 0 & A_{33} & 0 & \frac{1}{\sqrt{2}} A_{23} & \frac{1}{\sqrt{2}} A_{13} \\
   \frac{1}{\sqrt{2}} A_{12} & \frac{1}{\sqrt{2}} A_{12} & 0 & \frac12 (A_{11}+A_{22}) & \frac12 A_{13} & \frac12 A_{23} \\
   0 & \frac{1}{\sqrt{2}} A_{23} & \frac{1}{\sqrt{2}} A_{23} & \frac12 A_{13} & \frac12 (A_{22}+A_{33}) & \frac12 A_{12} \\
   \frac{1}{\sqrt{2}} A_{13} & 0 & \frac{1}{\sqrt{2}} A_{13} & \frac12 A_{23} & \frac12 A_{12} & \frac12 (A_{33}+A_{11})
   \end{array}
   \right]

This tensor form allows for a second-order weighting while using the fourth-order tensor. Any deviation from this tensor form necessarily induce a fundamentally fourth-order weighting.

*Isotropic, Fourth-Order*:

.. math::

   B_{IJ}^{iso} = \left[
   \begin{array}{c c c c c c}
   1 & 0 & 0 & 0 & 0 & 0 \\
   0 & 1 & 0 & 0 & 0 & 0 \\
   0 & 0 & 1 & 0 & 0 & 0 \\
   0 & 0 & 0 & 1 & 0 & 0 \\
   0 & 0 & 0 & 0 & 1 & 0 \\
   0 & 0 & 0 & 0 & 0 & 1
   \end{array}
   \right]

*Coordinate-Axis Aligned, Fourth-Order*:

.. math::

   B_{IJ}^{diag} = \left[
   \begin{array}{c c c c c c}
   A_{11} & 0 & 0 & 0 & 0 & 0 \\
   0 & A_{22} & 0 & 0 & 0 & 0 \\
   0 & 0 & A_{33} & 0 & 0 & 0 \\
   0 & 0 & 0 & \frac12 (A_{11}+A_{22}) & 0 & 0 \\
   0 & 0 & 0 & 0 & \frac12 (A_{22}+A_{33}) & 0 \\
   0 & 0 & 0 & 0 & 0 & \frac12 (A_{33}+A_{11})
   \end{array}
   \right]

*Coordinate-Axis Aligned Rotated, Fourth-Order*:

.. math::

   B_{IJ}^{diag} = B_{KL} M_{IK} M_{LJ}
   = B_{ijkl}^{diag,rot} = B_{mnpq}^{diag} R_{im} R_{jn} R_{kp} R_{lq}

where :math:`M` is a Mandel-notation rotation tensor [:footcite:`Auld1973`, :footcite:`Mehrabadi1990`].

*General, Fourth-Order*:

.. math::

   B_{IJ} = \left[
   \begin{array}{c c c c c c}
   B_{1111} & B_{1122} & B_{1133} & \sqrt{2} B_{1112} & \sqrt{2} B_{1123} & \sqrt{2} B_{1131} \\
   B_{1122} & B_{2222} & B_{2233} & \sqrt{2} B_{2212} & \sqrt{2} B_{2223} & \sqrt{2} B_{2231} \\
   B_{1133} & B_{2233} & B_{3333} & \sqrt{2} B_{3312} & \sqrt{2} B_{3323} & \sqrt{2} B_{3331} \\
   \sqrt{2} B_{1112} & \sqrt{2} B_{2212} & \sqrt{2} B_{3312} & 2 B_{1212} & 2 B_{1223} & 2 B_{1231} \\
   \sqrt{2} B_{1123} & \sqrt{2} B_{2223} & \sqrt{2} B_{3323} & 2 B_{1223} & 2 B_{2323} & 2 B_{2331} \\
   \sqrt{2} B_{1131} & \sqrt{2} B_{2231} & \sqrt{2} B_{3331} & 2 B_{1231} & 2 B_{2331} & 2 B_{3131} 
   \end{array}
   \right]

Calibration
-----------

For simple loading states, it may be possible to determine the anisotropic weighting tensor entries analytically, using experimentally observed failure strains and presumed plastic strain flow directions for several different loading directions. Examples of this are included in references [:footcite:`Stershic_AnisotropicFailure_Memo`, :footcite:`Stershic_AnisotropicFailure_USNCTAM`]. For more complex load paths or rate-dependence and temperature effects included, analytical calibration may be intractable, and a numerical optimization approach should be used instead.

Permissive Option
-----------------

Note that the numerator of the anisotropic weighting relations :eq:`mat:eq:anisotropic_Aij` and :eq:`mat:eq:anisotropic_Bijkl` include a square root of the tensor products. With no constraints imposed on :math:`A_{ij}` and :math:`B_{ijkl}`, this creates the possibility of having a negative argument to the square root, implying a weighting with imaginary terms. By default, LAMÉ will error out if a negative argument is passed to the square root. As an alternative, the ``ANISOTROPIC_EQPS_WEIGHT_PERMISSIVE = 1`` option allows the failure model to continue in these cases, by setting the argument (and weighting) to zero.

Choosing a positive-semi-definite weighting tensor would ensure that this issue is prevented. Further, since the weighting tensors are contracted with traceless flow direction *tensors* rather than vectors, the positive-semi-definite requirement can be weakened. For second-order tensors, it is only required for tensor invariant :math:`I_1(A)\geq 0` and :math:`I_2(A) \geq 0`, but :math:`I_3(A)` may be negative.
LAMÉ returns a *warning* when the minimum eigenvalue of :math:`A_{ij}` is negative and returns an *error* if :math:`I_2(A)` is negative.

This approach is difficult to generalize for the fourth-order weighting tensor, so no warnings or errors based on the structure of :math:`B_{ijkl}` are provided.

Material Coordinate Rotation
----------------------------

When performing the anisotropic weighting, the failure model first rotates the flow direction tensors from the global coordinate system to a material coordinate system, if one is provided. This is the same material coordinate system used by anisotropic plasticity models like Hill Plasticity. This is performed by applying the rotation matrices to the plastic flow direction tensor: :math:`N^{\text{matl}} = R^T N^{\text{global}} R`. The weighting is then done using :math:`N^{\text{matl}}` and anisotropic weight tensors, which are presumed to be in the material coordinates (i.e. :math:`A = A^{\text{matl}}`). 

Equivalently, the rotation matrices could be applied to the weighting tensor, e.g. :math:`A^{\text{global}} = R A^{\text{matl}} R^T`, which could be done either every timestep or even during user parameter specification. While the latter would be more efficient, it limits the use cases to those where the coordinate transformation is invariant in time and space (e.g. not a cylindrical coordinate system).

User Guide
----------

.. code-block:: sierrainput

   #
   # Anisotropic Failure model definitions
   #
   ANISOTROPIC_EQPS_WEIGHT_ORDER      = <int> 0 | 2 | 4 (0)
   ANISOTROPIC_EQPS_WEIGHT_PERMISSIVE = <bool> 0 | 1 (0)
   #
   # Definitons for second-order anisotropic weighting
   # (symmetric, 6 unique)
   ANISOTROPIC_EQPS_WEIGHT_XX = <real> (1.0)
   ANISOTROPIC_EQPS_WEIGHT_YY = <real> (1.0)
   ANISOTROPIC_EQPS_WEIGHT_ZZ = <real> (1.0)
   ANISOTROPIC_EQPS_WEIGHT_XY = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_YZ = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_ZX = <real> (0.0)
   #
   # Definitons for fourth-order anisotropic weighting
   # (major/minor symmetric, 21 unique)
   ANISOTROPIC_EQPS_WEIGHT_XXXX = <real> (1.0)
   ANISOTROPIC_EQPS_WEIGHT_YYYY = <real> (1.0)
   ANISOTROPIC_EQPS_WEIGHT_ZZZZ = <real> (1.0)
   ANISOTROPIC_EQPS_WEIGHT_XXYY = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_YYZZ = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_ZZXX = <real> (0.0)
   #
   ANISOTROPIC_EQPS_WEIGHT_XXXY = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_XXYZ = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_XXZX = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_YYXY = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_YYYZ = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_YYZX = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_ZZXY = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_ZZYZ = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_ZZZX = <real> (0.0)
   #
   ANISOTROPIC_EQPS_WEIGHT_XYXY = <real> (0.5)
   ANISOTROPIC_EQPS_WEIGHT_XYYZ = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_XYZX = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_YZYZ = <real> (0.5)
   ANISOTROPIC_EQPS_WEIGHT_YZZX = <real> (0.0)
   ANISOTROPIC_EQPS_WEIGHT_ZXZX = <real> (0.5)

Output variables available for this model are listed in :ref:`anisofail_stvar`.

.. _anisofail_stvar:

.. csv-table:: State Variables for Anisotropic Failure Model
   :align: center
   :delim: &
   :header: Name, Description

   ``FAILURE_ROT_XX`` & material rotation tensor, :math:`R_{11}`
   ``FAILURE_ROT_YY`` & material rotation tensor, :math:`R_{22}`
   ``FAILURE_ROT_ZZ`` & material rotation tensor, :math:`R_{33}`
   ``FAILURE_ROT_XY`` & material rotation tensor, :math:`R_{12}`
   ``FAILURE_ROT_YZ`` & material rotation tensor, :math:`R_{23}`
   ``FAILURE_ROT_ZX`` & material rotation tensor, :math:`R_{31}`
   ``FAILURE_ROT_YX`` & material rotation tensor, :math:`R_{21}`
   ``FAILURE_ROT_ZY`` & material rotation tensor, :math:`R_{32}`
   ``FAILURE_ROT_XZ`` & material rotation tensor, :math:`R_{13}`

.. raw::
   html

   <hr>

.. footbibliography::
