Calculation of tangent stiffness
==============================

In order to allow the global problem (equilibrium equations) to be solved by a Newton method, it is necessary to determine the coherent tangent matrix of the incremental problem. To do this, we once again adopt the convention for writing symmetric tensors of order 2 in the form of vectors with 6 components. So, for a :math:`a` tensor:

:math:`a\text{=}{}^{t}\text{}\left[\begin{array}{ccc}{a}_{\text{xx}}& {a}_{\text{yy}}& {a}_{\text{zz}}\end{array}\begin{array}{ccc}\sqrt{2}{a}_{\text{xy}}& \sqrt{2}{a}_{\text{xz}}& \sqrt{2}{a}_{\text{yz}}\end{array}\right]` [:ref:`éq4-1 <éq4-1>`]

If we also introduce the hydrostatic vector :math:`1` and the deviatoric projection matrix :math:`P`:

:math:`1\text{=}{}^{t}\text{}\left[\begin{array}{cccccc}1& 1& 1& 0& 0& 0\end{array}\right]` [:ref:`éq4-2 <éq4-2>`]

:math:`P\text{=}\text{Id}\text{-}\frac{1}{3}1\otimes 1` [:ref:`éq4-3 <éq4-3>`]

So the coherent tangent stiffness matrix is written for elastic behavior:

:math:`\frac{\partial s}{\partial \Delta \varepsilon }\text{=}K1\otimes 1\text{+}2\mup` [:ref:`éq4-4 <éq4-4>`]

and for plastic behavior:

:math:`\frac{\partial s}{\partial \Delta \varepsilon }\text{=}K1\otimes 1\text{+}2\mu (1\text{-}\frac{3\mu \Deltap }{{s}_{\text{eq}}^{e}})P\text{+}9{\mu }^{2}(\frac{\Deltap }{{s}_{\text{eq}}^{e}}\text{-}\frac{1}{{R}^{\text{'}}(p)+\frac{3}{2}(2\mu \text{+}C)})(\frac{{\tilde{s}}^{e}}{{s}_{\text{eq}}^{e}}\otimes \frac{{\tilde{s}}^{e}}{{s}_{\text{eq}}^{e}})` [:ref:`éq4-5 <éq4-5>`]

The initial tangent matrix, used by option RIGI_MECA_TANG, is obtained by adopting the behavior of the previous step (elastic or plastic, signified by an internal variable :math:`x` equal to 0 or 1) and by taking :math:`\Deltap \text{=}0` into the equation [eq].

**Note:**

*RIGI_MECA_TANG is the operator linearized with respect to* **time (cf.* [:external:ref:`R5.03.01 <R5.03.01>`] *,* [:external:ref:`R5.03.05 <R5.03.05>`] *) and corresponds to what is called the speed problem; in this case, linearization with respect to* :math:`\Deltau` *, in* :math:`\Deltau \text{=}0` *, provides the same expression.*

It is now proposed to demonstrate the expression [eq]. By differentiating between [eq] and [eq] at a fixed temperature, we immediately obtain:

:math:`\delta \sigma \text{=}\left[K1\otimes 1\text{+}2\mup \right]\delta \varepsilon \text{-}2\mu \delta {\varepsilon }^{p}` [:ref:`éq4-6 <éq4-6>`]

If the behavior regime is plastic, the incremental flow law [eq] then provides:

:math:`\delta {\varepsilon }^{p}=\frac{3}{2}\deltap \frac{{\tilde{s}}^{e}}{{s}_{\text{eq}}^{e}}\text{+}\frac{3}{2}\Deltap \delta (\frac{{\tilde{s}}^{e}}{{s}_{\text{eq}}^{e}})` [:ref:`éq4-7 <éq4-7>`]

As for :math:`\mathrm{dp}`, it is obtained by differentiating the implicit equation [eq]:

:math:`\left[\frac{3}{2}(2\mu \text{+}C)\text{+}{R}^{\text{'}}(p)\right]\deltap \text{=}\delta {s}_{\text{eq}}^{e}` [:ref:`éq4-8 <éq4-8>`]

Finally, all that's left to do is provide the :math:`{\tilde{s}}^{e}` variations:

:math:`\delta {\tilde{s}}^{e}\text{=}2\mu \delta \tilde{\varepsilon }{\mathrm{ds}}_{\text{eq}}^{e}\text{=}3\mu \frac{{\tilde{s}}^{e}}{{s}_{\text{eq}}^{e}}\cdot \delta \tilde{\varepsilon }\delta (\frac{{\tilde{s}}^{e}}{{s}_{\text{eq}}^{e}})\text{=}\frac{1}{{s}_{\text{eq}}^{e}}(2\mu \text{-}3\mu \frac{{\tilde{s}}^{e}}{{s}_{\text{eq}}^{e}}\otimes \frac{{\tilde{s}}^{e}}{{s}_{\text{eq}}^{e}})\cdot \delta \tilde{\varepsilon }` [:ref:`éq4-9 <éq4-9>`]

By then replacing [eq], [eq] and [eq] in [eq], we get the expression [eq].

This expression is formally identical to that defined in R5.03.02: [eq] and is written as:

:math:`\frac{\partial \sigma }{\partial \Delta \varepsilon }\text{=}K1\otimes 1\text{+}2\mu (1\text{-}\frac{3\mu \xi \Deltap }{{s}_{\text{eq}}^{e}})(\text{Id}-\frac{1}{3}1\otimes 1)\text{+}9{\mu }^{2}\xi (\frac{\Deltap }{{s}_{\text{eq}}^{e}}\text{-}\frac{1}{{R}^{\text{'}}\text{+}\frac{3}{2}(2\mu \text{+}C)})(\frac{{\tilde{s}}^{e}}{{s}_{\text{eq}}^{e}}\otimes \frac{{\tilde{s}}^{e}}{{s}_{\text{eq}}^{e}})`

with :math:`\xi \text{=}1` if :math:`\Delta \varepsilon` leads to lamination, and :math:`\xi \text{=}0` otherwise.

Using [eq], we find:

:math:`\frac{\partial s}{\partial \Delta \varepsilon }\text{=}{\lambda }^{\text{*}}\overrightarrow{1}\otimes \overrightarrow{1}\text{+}2{\mu }^{\text{*}}\text{Id}\text{-}\xi \frac{9{\mu }^{2}}{H(p)}(1\text{-}\frac{{R}^{\text{'}}\text{.}\Delta p}{R(p)})\frac{1}{{R}^{\text{'}}+\frac{3}{2}(2\mu \text{+}C)}(\frac{{\sigma }^{\text{dev}}}{R(p)}\otimes \frac{{\sigma }^{\text{dev}}}{R(p)})`

with :math:`{\lambda }^{\text{*}}\text{=}K\text{-}\frac{2\mu }{3}\frac{G(\Delta p)}{H(\Delta p)}2{\mu }^{\text{*}}\text{=}2\mu \frac{G(\Delta p)}{H(\Delta p)}`

for option FULL_MECA: :math:`{\sigma }^{\text{dev}}\text{=}\tilde{\sigma }-X`

for option RIGI_MECA_TANG: :math:`{\sigma }^{\text{dev}}\text{=}{\tilde{\sigma }}^{-}\text{-}{X}^{-}`

with :math:`H(\Delta p)\text{=}1+\frac{\frac{3}{2}(2\mu \text{+}C)\xi \Delta p}{R(p)}`

and :math:`G(\Delta p)=1+\frac{3}{2}C\xi \frac{\Delta p}{R(p)}`