\[\newcommand{\vector}[1]{\underline{#1}}
\newcommand{\vectorZero}{\vector{0}}
\newcommand{\tensTwo}[1]{\boldsymbol{#1}}
\newcommand{\tensTwoZero}{\tensTwo{0}}
\newcommand{\tensTwoUnit}{\tensTwo{I}}
\newcommand{\tensFour}[1]{\mathbb{#1}}
\newcommand{\inverse}[1]{{#1}^{-1}}
\newcommand{\transpose}[1]{{#1}^{T}}
\newcommand{\inverseTranspose}[1]{ {#1}^{-T}}
\newcommand{\vectorCmpCO}[2]{#1_{\left(#2\right)}}
\newcommand{\vectorCmpCT}[2]{#1^{\left(#2\right)}}
\newcommand{\tensTwoCmpCO}[3]{#1_{\left(#2 #3\right)}}
\newcommand{\tensTwoCmpCT}[3]{#1^{\left(#2 #3\right)}}
\newcommand{\tensTwoInva}[2]{I^{#1}_{#2}}
\newcommand{\tensTwoDevia}[1]{\boldsymbol{\tilde{#1}}}
\newcommand{\tensFourCmpCO}[5]{#1_{\left(#2 #3 #4 #5\right)}}
\newcommand{\tensFourCmpCT}[5]{#1^{\left(#2 #3 #4 #5\right)}}
\newcommand{\domain}{\Omega}
\newcommand{\domainRefe}{\domain^{0}}
\newcommand{\domainCurr}{\domain^{t}}
\newcommand{\bound}{\partial\domain}
\newcommand{\boundRefe}{\bound^{0}}
\newcommand{\boundCurr}{\bound^{t}}
\newcommand{\boundN}{\bound_{\textrm{N}}}
\newcommand{\boundNRefe}{\boundN^{0}}
\newcommand{\boundNCurr}{\boundN^{t}}
\newcommand{\boundD}{\bound_{\textrm{D}}}
\newcommand{\boundDRefe}{\boundD^{0}}
\newcommand{\boundDCurr}{\boundD^{t}}
\newcommand{\normal}{\vector{n}}
\newcommand{\normalRefe}{\normal^{0}}
\newcommand{\normalCurr}{\normal^{t}}
\newcommand{\posi}{\vector{x}}
\newcommand{\posiRefe}{\posi^{0}}
\newcommand{\posiCurr}{\posi^{t}}
\newcommand{\disp}{\vector{u}}
\newcommand{\dispVirt}{\delta\vector{v}}
\newcommand{\funcTransfor}{\vector{\varphi}^{t}}
\newcommand{\gradTransfor}{\tensTwo{F}}
\newcommand{\jacobTransfor}{J}
\newcommand{\posiIncr}{\vector{dx}}
\newcommand{\posiIncrRefe}{\posiIncr^{0}}
\newcommand{\posiIncrCurr}{\posiIncr^{t}}
\newcommand{\strainCmp}{\varepsilon}
\newcommand{\strain}{\tensTwo{\strainCmp}}
\newcommand{\EGLCmp}{E}
\newcommand{\ECGDroiteCmp}{C}
\newcommand{\ECGGaucheCmp}{B}
\newcommand{\ELOGCmp}{{E}_{ln}}
\newcommand{\EGL}{\tensTwo{\EGLCmp}}
\newcommand{\ECGDroite}{\tensTwo{\ECGDroiteCmp}}
\newcommand{\ECGGauche}{\tensTwo{\ECGGaucheCmp}}
\newcommand{\ELOG}{\tensTwo{\ELOGCmp}}
\newcommand{\EDilCmp}{C^{\star}}
\newcommand{\EDil}{\tensTwo{\EDilCmp}}
\newcommand{\divTensTwo}[1]{\vector{\nabla} {\cdot} #1 }
\newcommand{\divVector}[1]{\vector{\nabla} {\cdot} #1 }
\newcommand{\gradScal}[1]{\vector{\nabla} \times #1 }
\newcommand{\gradVector}[1]{\vector{\nabla} \times #1 }
\DeclareMathOperator{\rand}{rand}
\DeclareMathOperator{\round}{round}
\DeclareMathOperator{\trace}{trace}
\newcommand{\metric} {\tensTwo{g}}
\newcommand{\vectorBaseCV}[1]{\vector{g_{#1}}}
\newcommand{\vectorBaseCT}[1]{\vector{g^{#1}}}
\newcommand{\metricRefe} {\metric^{t}}
\newcommand{\metricCurr} {\metric^{0}}
\newcommand{\discVect}[1]{\lbrace #1 \rbrace}
\newcommand{\discVectLigne}[1]{\langle #1 \rangle}
\newcommand{\discVectZero}{\discVect{0}}
\newcommand{\discVectCmp}[2]{\lbrace #1 \rbrace_{(#2)}}
\newcommand{\discVectLigneCmp}[2]{\langle #1 \rangle_{(#2)}}
\newcommand{\discMatr}[1]{\left[ #1 \right]}
\newcommand{\discMatrZero}{\discMatr{0}}
\newcommand{\discMatrCmp}[3]{\left[ #1 \right]_{(#2#3)}}
\newcommand{\onQuadPoint}[2]{{#1}_{#2}}
\newcommand{\onNode}[2]{{#1}^{#2}}
\newcommand{\quadOrder}{k_{Q}}
\newcommand{\nbQuadPoint}{n_Q}
\newcommand{\quadWeight}[1]{\omega_{#1}}
\newcommand{\quadPointIndex}{i_{\textrm{pg}}}
\newcommand{\nodeIndex}{i_{\textrm{no}}}
\newcommand{\pres}{p}
\newcommand{\presRefe}{\pres^{0}}
\newcommand{\presCurr}{\pres^{t}}
\newcommand{\stressCmp}{\sigma}
\newcommand{\stress}{\tensTwo{\stressCmp}}
\newcommand{\stressPKTwoCmp}{S}
\newcommand{\stressPKTwo}{\tensTwo{\stressPKTwoCmp}}
\newcommand{\work}[1]{W^{#1}}
\newcommand{\enerInterne}{\Psi_i}
\newcommand{\yieldStress}{\sigmaCmp_{Y}}
\newcommand{\youngModulus}{E}
\newcommand{\youngModulusCplx}{E^{\star}}
\newcommand{\poissonCoef}{\nu}
\newcommand{\poissonCoefCplx}{\nu^{\star}}
\newcommand{\shearModulus}{G}
\newcommand{\bulkModulus}{K}
\newcommand{\modulusTangent}{\tensFour{K}}
\newcommand{\anglDila}{\psi}
\newcommand{\anglFric}{\varphi}
\newcommand{\cohesion}{c}
\newcommand{\presCapi}{\pres_{c}}
\newcommand{\presEau}{\pres_{w}}
\newcommand{\presAir}{\pres_{a}}
\newcommand{\shapeFunc}{\Phi}
\newcommand{\shapeDFunc}{B}
\newcommand{\measLine}{\, \mathrm{d} l}
\newcommand{\measDomain}{\, \mathrm{d} \domain}
\newcommand{\measBound}{\, \mathrm{d} \Gamma}
\newcommand{\normEucl}[1]{\Vert{#1} \Vert }
\newcommand{\derivee}[2]{#1_{, #2}}
\newcommand{\soundSpeed}{c}
\newcommand{\soundSpeedComp}{c_{P}}
\newcommand{\soundSpeedCisa}{c_{S}}
\newcommand{\posiTang}{\posi_{\tau}}
\newcommand{\posiNorm}{{x}_{n}}
\newcommand{\dispTang}{\disp_{\tau}}
\newcommand{\dispNorm}{{u}_{n}}
\newcommand{\basisVector}[1]{\vector{e}_{#1}}
\newcommand{\sigmVector}{\vector{t}}
\newcommand{\sigmVectorTang}{\sigmVector_{\tau}}
\newcommand{\sigmVectorNorm}{t_{n}}
\newcommand{\lameLambda}{\lambda}
\newcommand{\lameMu}{\mu}\]
3. Discretization
3.1. Calculated options
Finite skin elements (surface elements immersed in three-dimensional space) are used to discretize the real and virtual movements involved in surface expressions such as eq-5 and (2.14). The latter make it possible to express respectively the second member vector and the stiffness matrix due to pressure, the use of which by the algorithm of STAT_NON_LINE is specified in [R5.03.01]
We therefore develop four options:
RIGI_MECA_PRSU_R: stiffness matrix for a follower pressure as a real constant
RIGI_MECA_PRSU_F: stiffness matrix for a follower pressure as a real function
CHAR_MECA_PRSU_R: second-member vector for a follower pressure as a real constant
CHAR_MECA_PRSU_F: second-member vector for a follower pressure as a real function
These options are developed for 3D elements, D_ PLAN and AXIS. It is possible to apply normal follower pressure but no follower tangent shear. Pressure can be a real function of time or a real constant.
3.2. Discretization of differential geometry terms
We note \(\discVect{G_{\alpha}}\) the discretized version of the vectors of the covariant basis with \({\alpha={1,2}}\). `
The first two vectors of the covariant base \(\discVect{G_{\alpha={1,2}}}\) are calculated from the vector of discretized displacements at the nodes \(\discVect{U}\) and from the derivatives of the form functions \(\discMatr{B_{\alpha}}\):
\[\]
: label: eq-20
discVect {G_ {alpha}}
=
discMatr {B_ {alpha}}},discVect {U}
Normal \(\normal\) is calculated as the cross product of these first two \(\vector{g_{\alpha}}\) vectors:
\[\]
: label: eq-21
nnormal
=
frac {{vector {g}} _ {1} {wedge} {wedge} {vector {g}} {vector {g}} _ {1} {wedge} {wedge} {wedge} {vector {g}}} _ {2}Green}
We can also calculate the discretized metric tensor in covariant components \(\discMatr{G_{\alpha \beta}}\):
(3.1)\[ \ discMatr {G_ {\ alpha\ beta}}
=
\ discVect {G_\ alpha}\ discVectLine {{G} _ {\ beta}}\]
And Jacobian sound:
(3.2)\[ \ JacobTransfor
=
\ det\ discMatr {G_ {\ alpha\ beta}}\]
We can calculate the contravariant metric matrix that will be noted \(\discMatr{G^{\delta \gamma}}\)
(3.3)\[ \ discMatr {G^ {\ delta\ gamma}}
=
\ inverse {\ discMatr {G_ {\ alpha\ beta}}}\]
And finally extract the contravariant database \(\discVectCmp{G}{\delta}\):
(3.4)\[ \ discVect {G^\ delta}
=
\ discMatr {G^ {\ delta\ gamma}}}\,\ discVect {G_\ alpha}\]
3.3. Vector of the following forces
The calculation of the virtual work of pressure forces eq-5 is in fact identical to that carried out in small movements, provided that the geometry is updated beforehand. We start from the expression virtual works:
(3.5)\[ \ work {\ pres} (\ disp)
=
{\ int} _ {\ bound NCurr (\ disp)} {-\ pres\,\normalCurr\ cdot\ dispVirt\ measBound}\]
In discretized form:
\[\]
: label: eq-27
work {pres} (disp)
=
discVectLine {delta V}discVect {{F} ^ {pres}left (dispright)}
The variation in displacements is written from the form functions:
\[\]
: label: eq-28
DispVirt
=
discMatr {shapeFunc},left{delta Vright}
We have discretized all the terms of differential geometry in the previous paragraph, all we have to do is discretize the integral using a Gauss diagram with weights \(\quadWeight{\quadPointIndex}\). On any quantity \(A\), we will have
(3.6)\[ {\ int} _ {\ bound NCurr (\ disp)} {A\ measBound}
=
\ sum_ {\ quadpointIndex} {\ onQuadpoint {A} {\ quadpointIndex}\,\ quadWeight {\ quadpointIndex}}\]
The \(\pres\) pressure given by AFFE_CHAR_MECA is located at the nodes,
We must therefore interpolate the pressure \(\pres\) from its values on the nodes that we note \(\onNode{\pres}{\nodeIndex}\)
to the Gauss point with coordinates \(\onQuadPoint{\xi}{\quadPointIndex}\):
\[\]
: label: eq-30
onQuadpoint {pres} {quadpointIndex}
=
sum_ {nodeIndex} {onNode {ShapeFunc (onQuadPoint {xi} {quadPointIndex})} {nodeIndex}\,onNode {pres} {pres} {nodeIndex}}
We denote \(\onNode{\shapeFunc(\onQuadPoint{\xi}{\quadPointIndex})}{\nodeIndex}\) the value of the shape function of the node \(\nodeIndex\) at the point of gauss with coordinates \(\onQuadPoint{\xi}{\quadPointIndex}\). Finally:
\[\]
: label: eq-31
discVect {{F} ^ {pres}left (dispright)}
= -
sum_ {{quadPointIndex}} {
onQuadpoint {pres} {quadpointIndex}
,
quadWeight {quadPointIndex}
,
onQuadpoint {JacobTransfor} {quadpointIndex} {quadpointIndex}
,
transpose {discMatr {onQuadpoint {ShapeFunc} {quadPointIndex}}}}
,
discVect {onQuadpoint {normal} {quadpointIndex}}
}
3.4. Matrix of the following forces
The variation in the virtual work of pressure efforts is expressed by the expression (2.14)
(3.7)\[ \ frac {\ partial\ work {\ pres} (\ pres} (\ disp)} {\ partial\ disp}\ cdot\ delta\ disp\ cdot\ dispVirt
=
{\ int} _ {\ bound NCurr (\ disp)} {
-\ presCurr\,
\ left [
\ left (\ frac {\ partial\ delta\ disp} {\ partial {\ theta} ^ {\ alpha}}\ cdot\ VectorBasect {\ alpha}\ disp} {\ right)
\,
\ left (\ DispVirt\ cdot\normal\ right)
-
\ left (\ frac {\ partial\ delta\ disp} {\ partial {\ theta} ^ {\ alpha}}\ cdot\normal\ disp}} {\ partial\ delta\ disp} {\ partial\ disp} {\ partial {\ theta} ^ {\ alpha}}\ cdot\normal\
\,
\ left (\ dispVirt\ cdot\ VectorBasect {\ alpha}\ right)
\ right]\ measBound}\]
It allows you to extract the tangent matrix \(\discMatr{{K}^{\pres}\left(\disp\right)}\) from the following forces:
(3.8)\[ \ frac {\ partial\ work {\ pres} (\ pres} (\ disp)} {\ partial\ disp}\ cdot\ delta\ disp\ cdot\ dispVirt
=
\ discVectLine {\ delta V}\ discMatr {{-K} ^ {\ pres}\ left (\ disp\ right)}\ discVect {\ delta U}\]
The minus sign comes from the fact that the matrix’s contribution is to the first member.
Finally:
(3.9)\[ \ discMatr {{K} ^ {\ pres}\ left (\ disp\ right)}
=
\ sum_ {{\ quadPointIndex}} {
\ onQuadpoint {\ pres} {\ quadpointIndex}
\,
\ quadWeight {\ quadPointIndex}
\,
\ onQuadpoint {\ JacobTransfor} {\ quadpointIndex} {\ quadpointIndex}
\,
\ transpose {\ discMatr {\ onQuadpoint {\ shape DFunc} {\ quadpointIndex}}}
\ left (
\ discVect {\ onQuadpoint {G^\ delta} {\ quadpointIndex}}
\,
\ discVectLine {\ onQuadpoint {\normal} {\ quadpointIndex}}
-
\ discVect {\ onQuadpoint {\ normal} {\ quadpointIndex}}
\,
\ discVectLine {\ onQuadpoint {G^\ delta} {\ quadpointIndex}}
\ right)
\,
{\ discMatr {\ onQuadpoint {\ shapeFunc} {\ quadpointIndex}}}
}\]
with \(\discMatr{\onQuadPoint{ \shapeDFunc} {\quadPointIndex} }\) the matrix of derivatives of form functions.
3.5. Choice of the matrix
In general, the matrix \(\discMatr{{K}^{\pres}\left(\disp\right)}\) is not symmetric (except in the specific case of a structure subjected to constant internal or external pressure). We also note in practice that for strong variations in geometry (using hyper-elastic behavior in large deformations like ELAS_HYPER), the fact of symmetrizing this matrix is not a good strategy (convergence failure). We therefore decide to keep this non-symmetric matrix, despite the (slight) additional cost induced by the factorization of such a matrix.