6. Resolution algorithm#

The variational problem (-) can be written in the form:

(6.16)#\[ F (U) = {L} ^ {i} (U) - {L} ^ {\ mathrm {meca}} =0\]

where \(U\) refers to widespread travel.

The associated tangent operator is noted \(DF=\frac{\partial F}{\partial U}\).

6.1. Generalized stresses and deformations#

We note \({U}^{\mathrm{el}}\) the vector of nodal unknowns on the element \(\mathrm{el}\) , \({E}_{\mathrm{pi}}^{\mathrm{el}}\) the vector of generalized deformations at the point of integration \(\mathrm{pi}\) of the element \(\mathrm{el}\) , \({\Sigma }_{\mathrm{pi}}^{\mathrm{el}}\) the vector of generalized constraints

(6.16)#\[\begin{split} {\ Sigma} _ {\ mathrm {pi}}} ^ {\ mathrm {el}} = {\ left\ {\ begin {array} {c}\ tau\ text {'}\\ -p\\ w\\\\\\ W\\\\ W\\\\ W\\\\ W\\\ W\\\ W\\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\\ W\ p\\ {p} ^ {\ text {-}}} -p\ end {array}\ right\}} ^ {\ mathrm {el}}\end{split}\]
(6.16)#\[\begin{split} {E} _ {\ mathrm {pi}}} ^ {\ mathrm {el}}} = {\ left\ {\ begin {array} {c} u\\ p\\\nabla p\\ {p}} ^ {\ text {p}} ^ {\ text {+}} ^ {\ text {-}}\\ {q}} ^ {\ text {-}} ^ {\ text {-}} ^ {\ text {+}}}\\ {q} ^ {q} ^ {q} ^ {q} ^ {q} ^ {\ text {-}} ^ {\ text {+}}}\\ {q} ^ {q} ^ {q} ^ {q} ^ {q} ^ {q} ^ {q} ^ {q} ^ {\ text {+}}\ text {-}}\ end {array}\ right\}}} ^ {\ mathrm {el}}\end{split}\]

Generalized deformations are obtained by the relationship

(6.16)#\[ {E} _ {\ mathrm {pi}}} ^ {\ mathrm {el}} = {Q} _ {\ mathrm {pi}}} ^ {\ mathrm {el}}}\ mathrm {el}}\ mathrm {.} {U} ^ {\ mathrm {el}}\]

where \({Q}_{\mathrm{pi}}^{\mathrm{el}}\) is the transition matrix from nodal degrees of freedom to generalized deformations at integration points \(\mathrm{pi}\). It is calculated from the shape functions and the orientation of the element. In fact, the mechanical movements of each node are defined in the global coordinate system. In order to extract the components normal and tangential to the joint element, a rotation matrix \({\Theta }_{\mathit{pi}}^{\mathit{el}}\) is applied to the nodal displacement vectors.

6.2. Integration#

To integrate the terms of the vector onto the element, a selective integration method is chosen. This makes it possible to avoid numerical oscillations for problems where the structure is subject to a shock and where mechanical phenomena predominate. This method consists in integrating at the vertices of the element the terms involving a derivative with respect to time and in integrating the permanent terms into the Gauss points.

The internal force vector is written

(6.16)#\[ {L} ^ {i} (U) = {\ sum} _ {\ mathrm {el}}} ({\ sum}} _ {g} {R} _ {g} _ {g} ^ {\ mathrm {el}} ^ {\ mathrm {el}}}} ({U}} ^ {\ mathrm {el}}}) + {\ sum} _ {s} {\ omega} _ {g} _ {g} ^ {\ mathrm {el}}} ^ {\ mathrm {el}}}} ({U} ^ {\ mathrm {el}}}) + {\ sum} _ {s} {\ omega} _ {g} _ {g}} ^ {\ mathrm {el}}} ^ {\ mathrm {el}}} {s} ^ {\ mathrm {el}}} ({U} ^ {\ mathrm {el}}))\]

by noting \({R}_{g}^{\mathrm{el}}\) and \({R}_{s}^{\mathrm{el}}\) the values respectively at the Gauss point \(g\) and at the top \(s\) of the nodal forces and \({\omega }_{g}\) and \({\omega }_{s}\) the integration weights of \(g\) and \(s\) respectively.

The tangent operator is written

(6.16)#\[ DF (U) = {\ sum} _ {\ mathrm {el}} ({\ mathrm {el}}}} ({\ sum}} _ {g} _ {g} ^ {\ mathrm {el}} (U) + {\ mathrm {el}}} (U) + {\ mathrm {el}}} (U) + {\ sum}}} (U) + {\ sum}} _ {\ mathrm {el}} (U))\]

by noting \(D{F}_{g}^{\mathrm{el}}\) and \(D{F}_{s}^{\mathrm{el}}\) the values respectively at the Gauss point \(g\) and at the vertex \(s\) of the tangent operator.

The tangent operator and the nodal forces are therefore calculated differently at Gauss points and at the vertices. In contrast, all the components of the generalized stresses and all the internal variables are calculated at both the Gauss points and the integration points.

6.3. Internal forces vector: options RAPH_MECA and FULL_MECA#

The following decomposition is adopted

(6.16)#\[ {{\ stackrel {} {E}}} _ {\ mathrm {pi}}} ^ {\ mathrm {el}}}} ^ {T}\ mathrm {.} {\ stackrel {} {\ Sigma}}} _ {\ mathrm {pi}}} ^ {\ mathrm {el}} = {\ stackrel {} {\ Sigma}} _ {u} {\ stackrel {}} {\ stackrel {}}} _ {p} {\ stackrel {}}} _ {p} {\ stackrel {}}} _ {p} {\ stackrel {}} _ {p} {\ stackrel {}}} _ {p} {\ stackrel {}} p} + {\ stackrel {} {\ Sigma}} {\ Sigma}} _ {\nabla p}}\nabla\ stackrel {} {\ Sigma}} _ {{p}} _ {{p}} ^ {\ text {+}} ^ {\ text {+}}} {\ stackrel {+}}} _ {{p}}} _ {{p}}} _ {{p}} ^ {\ stackrel {}}} _ {{p}}} _ {{p}} ^ {\ stackrel {}}} _ {{p}}} _ {{p}} ^ {\ stackrel {}}} _ {{p}}} _ {{p}} ^ {\ stackrel {}}} _ {Sigma}} _ {{p} ^ {\ text {-}}}} {\ stackrel {-}}} {\ stackrel {-}} ^ {\ stackrel {} {\ Sigma}}} _ {{q}}} _ {{q}} ^ {\ text {-}}} {\ stackrel {+}}} + {{q}} _ {{q}}} _ {{q}}} _ {{q}} _ {q}} _ {{q}} _ {q}} _ {{q}} ^ {\ stackrel {}}} _ {{q}}} _ {{q}} _ {q}} _ {{q}}} _ {{q}}} _ {{q}}} _ {{q}}}\ Sigma}} _ {{q} ^ {\ text {-}}}} {\ stackrel {} {q}}} ^ {\ text {-}} ^ {\ text {-}}\]

where \({\stackrel{ˆ}{E}}_{\mathrm{pi}}^{\mathrm{el}}=(〚{\stackrel{ˆ}{u}}_{g}〛,\stackrel{ˆ}{p},\nabla \stackrel{ˆ}{p},{\stackrel{ˆ}{p}}^{\text{+}},{\stackrel{ˆ}{p}}^{\text{-}},{\stackrel{ˆ}{q}}^{\text{+}},{\stackrel{ˆ}{q}}^{\text{-}})\) is a virtual generalized deformation calculated from a virtual generalized displacement vector.

From the discrete variational formulations and by distributing the stationary terms at the Gauss points and the unsteady terms at the vertices of the element, we obtain:

  • at Gauss points

(6.16)#\[\begin{split} \ begin {array} {ccc} {\ stackrel {} {\ Sigma}}} _ {u} ^ {\ mathrm {el}} &\ tau\ text {=} &\ tau\ text {'} - {p} - {p} _ {p} _ {p}} - {p}} - {p}} - {p}} - {p} - {p}} - {p} _ {p}} - {p} _ {p}} _ {\ mathrm {el}}} - {p} _ {p}} - {p} _ {p}} _ {\ mathrm {el}}} - {p} _ {p}} - {p} _ {p}} _ {\ mathrm {el}}} - {p}} - {p}\ Delta t ({q} _ {g} ^ {\ text {+}}} + {\ text {+}}} + {\ text {-}})\\ {\ stackrel {} {\ Sigma}}} _ {\nabla p}} _ {\nabla p}} ^ {\ mathrm {el}}} ^ {\ mathrm {el}}} + {q}} _ {\ mathrm {el}}} + {g} ^ {g} ^ {g} ^ {\ mathrm {el}}} &\ text {=}} &\ Delta tW\\ {\ stackrel {} {\ Sigma}}} _ {{p} ^ {\ text {+}}}}} ^ {\ mathrm {el}}} ^ {\ mathrm {el}}} ^ {\ text {+}}\\ {\ stackrel {} {} {\ Sigma}} {\ Sigma}}} _ {{p}}} _ {{p}} ^ {\ text {-}} ^ {\ mathrm {el}} ^ {\ mathrm {el}}} ^ {\ mathrm {el}}}\\ text {+}}}\\ {\ stackrel {}} {\ stackrel {}} {\ stackrel {} {\}} {\ stackrel {}} {\ Sigma}} {\ Sigma}}} _ {{p}} q} ^ {\ text {-}}\\ {\ stackrel {} {\ stackrel {}} {\ sigma}} _ {{q} ^ {\ text {+}}} ^ {\ mathrm {el}}} &\ text {=}} & {=}} & {=}} & {p}} & {p}} ^ {p} ^ {p} ^ {{p}} ^ {\ text {+}} -p\\ {\ stackrel {}} {\ sigma}}} _ {{q}} ^ {{q} ^ { \ text {-}}}} ^ {\ mathrm {el}}} &\ text {=} & {p} ^ {\ text {-}}} -p\ end {array}\end{split}\]
  • at the summits

(6.16)#\[\begin{split} \ begin {array} {ccc} {\ stackrel {} {\ Sigma}}} _ {u} ^ {\ mathrm {el}} &\ text {=} & 0\\ {\ stackrel {} {} {\ Sigma}} {\ Sigma}}} _ {\ Sigma}}} _ {n+1} _ {n} - {w}} {\ Sigma}}} _ {n+1} _ {n+1}\\ {\ stackrel {} {\ Sigma}}} _ {\nabla p}} _ {\nabla p}} ^ {\ mathrm {el}} & 0\\ {\ stackrel {} {\ sigma}} {\ sigma}}} _ {{p}}} _ {{p}}} _ {{p}}} _ {{p}} ^ {\ text {+}}}} ^ {\ mathrm {el}}} &\\ {\ stackrel {=} & 0\\ {\ stackrel}} &\\ {\ stackrel {=} & 0\\ {\ stackrel}} & 0\\ {\ stackrel}} & 0\\ {\ stackrel}} {} {\ Sigma}} _ {{p} ^ {\ text {-}}}} ^ {\ mathrm {el}} &\ text {=} & 0\\ {\ stackrel {} {\ Sigma}} {\ sigma}}} _ {{q}}} _ {{q}}} _ {{q}} ^ {\ text {+}}}} ^ {\ mathrm {el}}} &\\ {\ stackrel {=} & 0\\ {\ stackrel {=} & 0\\ {\ stackrel}} & 0\\ {\ stackrel}} & 0\\ {\ stackrel {\ stackrel}} & 0\\ {\ stackrel {\ stackrel}} & 0\\ {\ stackrel} {\ Sigma}} _ {{q}} ^ {\ text {-}}}} ^ {\ mathrm {el}} &\ text {=}} & 0\ end {array}\end{split}\]

In addition, we have:

(6.16)#\[ {{\ stackrel {} {U}}}} ^ {\ mathrm {el}}} ^ {T}\ mathrm {.} {R} _ {\ mathrm {pi}}} ^ {\ mathrm {el}}} = {{\ stackrel {} {E}}} _ {\ mathrm {pi}}} ^ {\ mathrm {pi}}} ^ {\ mathrm {el}}} ^ {T}\ mathrm {.} {\ stackrel {} {\ Sigma}}} _ {\ mathrm {pi}} ^ {\ mathrm {el}}\]

What gives

(6.16)#\[ {R} _ {\ mathrm {pi}}} ^ {\ mathrm {el}}} = {{Q} _ {\ mathrm {pi}} ^ {\ mathrm {el}}}} ^ {el}}} ^ {el}}} ^ {T}\ mathrm {.} {\ stackrel {} {\ Sigma}}} _ {\ mathrm {pi}} ^ {\ mathrm {el}}\]

6.4. Tangent operator: options FULL_MECA and RIGI_MECA_TANG#

The tangent operator at the Gauss point \(g\) is given by:

(6.16)#\[\begin{split} \ frac {\ partial {\ Sigma} _ {g}} ^ {\ mathrm {el}}} {\ partial {E} _ {g} ^ {\ mathrm {el}}}} =\ left [\ begin {array} {\ array} {cccccc} {cccccc} {\ cccccc}\ frac {\ mathrm {el}}}} =\ left [\ begin {array} {{array} {{array} {cccccc}}\ frac {\ partial {el}}}} =\ left [\ begin {array} {{array} {{array} {cccccc}}\ frac {\ partial {el}}} {\ partial {el}}} =\ left [\ begin {array} {{array} {{array} {cccccc}\ frac {\ partial {\ tau} _ {n}\ text {'}}} {\ partial {u} _ {t}} & -1& 0& 0& 0\\ frac {\ partial {\ tau} 0\\ frac {\ partial {\ tau} _ {\ tau} _ {t}} {t}} {t}} {t}} {t}} {t}} {t}} {t}} {t}} {t}} {t}} {\ partial {u} _ {n}}} &\ frac {\ partial {\ tau} _ {t}} {\ partial} {u} _ {t}} & 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0&\ Delta t&\\ Delta t\ frac {\ partial W} {\ partial W} {\ partial W} {\ partial W} & -\ delta t\ frac {\ partial p} & -\ delta t\ frac {\ partial W} & -\ delta t\ frac {\ partial W}\ frac {\ partial W} & -\ delta t\ frac {\ partial W}\ frac {\ partial W} & -\ delta t\ frac {\ partial W} & -\ delta t\ frac {\ partial W} {\ partial W} & -\ delta t\ frac {\ partial W}\ frac {\ partial W} 0& 0\\ 0& 0& 0& 0& 0& 0& 0& -\ Delta t& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& -\ Delta t \\ 0& 0& -1& 0& 0& 1& 0& 0& 0& 0& 0\\ 0& 0& -1& 0& 0& 0\ end {array}\ right]\end{split}\]

The tangent operator at vertex \(s\) is given by:

(6.16)#\[\begin{split} \ frac {\ partial {\ Sigma} _ {s}} ^ {\ mathrm {el}}} {\ partial {E} _ {s} ^ {\ mathrm {el}}}} =\ left [\ begin {array} {\ begin {array} {cccccc} {cccccc} 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0\\ -\ frac {\ partial w} {\ partial {u} _ {n}} _ {n}} & 0& -\ frac {\ partial w} {\ partial p} & 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 00& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\ & 0\ end {array}\ right]\end{split}\]

The operator tangent to integration point \(\mathrm{pi}\) is finally given by:

(6.16)#\[ D {F} _ {\ mathrm {pi}}} ^ {\ mathrm {el}}} = {{Q} _ {\ mathrm {pi}}} ^ {\ mathrm {el}}}} ^ {el}}}} ^ {T}\ mathrm {.} \ frac {\ partial {\ Sigma} _ {\ mathrm {pi}}} ^ {\ mathrm {el}}} {\ partial {E} _ {\ mathrm {pi}}} ^ {\ mathrm {pi}}} ^ {\ mathrm {el}}}}\ mathrm {.} {Q} _ {\ mathrm {pi}}} ^ {\ mathrm {el}}\]

6.5. Vector nodal forces: option FORC_NODA#

Option FORC_NODA is used by STAT_NON_LINE during the prediction phase.

At integration point \(\mathrm{pi}\), the vector nodal forces at time \(n+1\) is defined by:

(6.16)#\[ {F} _ {\ mathrm {pi}, n+1}} ^ {\ mathrm {el}} = {{Q} _ {\ mathrm {pi}}} ^ {\ mathrm {el}}}} ^ {\ mathrm {el}}}} ^ {T}\ mathrm {.} {\ stackrel {} {S}}} _ {\ mathrm {pi}} ^ {\ mathrm {el}}\]

where \({\stackrel{ˉ}{S}}_{\mathrm{pi}}^{\mathrm{el}}\) is zero at the vertices and is given at the Gauss points by:

(6.16)#\[\begin{split} \ begin {array} {ccc} {\ stackrel {} {S}}} _ {u, n+1}} ^ {\ mathrm {el}} & {\ tau} _ {g, n} - {g, n} - {p} - {p} - {p}} - {p}} - {p, n}} - {p, n}} - {p, n} - {p}} - {p}, n} - {p}} - {p}} - {p}, n} - {p}} - {p}} - {p}, n} - {p}} - {p}} - {p}, n} - {p}} - {p}} - {p}, n} - {p}} - {p}, n} - {p}}} &\ text {=} &\ Delta t ({q} _ {g, n} _ {g, n} _ {q}} _ {g, n} ^ {\ text {-}})\\ {\ stackrel {} {S}} _ {\ text {-}})\\ {\ stackrel {} {S}} _ {\ text {-}})\\ {\ stackrel {} {S}}} _ {\ text {-}})\\ {\ stackrel {} {S}}} _ {\nabla p, n+1}} ^ {\ mathrm {el}}} &\ text {=}})\\ {\ stackrel {} {\ stackrel {} {S}} {S}} _ {} _ {g, n}\\ {\ stackrel {} {S}}} _ {{p}} _ {{p}} ^ {\ text {+}}, n+1} ^ {\ mathrm {el}}} &\ text {=}} &\\ text {=}} &\\ text {=}} &\ text {=}} &\ text {=}} &\ text {=}} & -\ text {=}} & -\ Delta t {q}} & -\ Delta t {q} {S}} & -\ Delta t {q}} & -\ Delta t {q} {S}} & -\ Delta t {q} {S}} & -\ Delta t {q}} _ {q} _ {p} _ {p}} ^ {\ text {-}}, n+1}} ^ {\ mathrm {el}} ^ {\ mathrm {el}}} ^ {\ text {-}} ^ {\ text {-}}\\ {\ stackrel {} {el}}} ^ {\ mathrm {el}}}}\\ {\ stackrel {} {S}}} &\ text {S}}} _ {{q}} &\ text {=} & {p} _ {g, n} ^ {\ text {+}} - {p} _ {g, n}\\ {\ stackrel {} {S}} _ {{q}} ^ {\ text {-}}, n+1} ^ {\ mathrm {el}} &\ text {=} & {p} _ {g, n} _ {g, n}} ^ {\ text {-}} ^ {\ text {-}}} - {p} _ {g, n}\ end {array}} _ {array}\end{split}\]

6.6. Internal variables#

As for the classical THM models, the variables from \(1\) to \(\mathrm{NVIM}\) (\(\mathrm{NVIM}\) depending on the law of mechanical behavior used) relate to the law of mechanics used. The following internal variables are:

  • \({V}_{\mathrm{NVIM}+1}\): \(\rho -{\rho }_{0}\), density variation,

  • \({V}_{\mathrm{NVIM}+2}\): \(\varepsilon\), opening the crack.