\[\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{\vectorCmp}[2]{#1_{\left(#2\right)}} \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} #1 } \newcommand{\gradVector}[1]{\vector{\nabla} \times #1 } \DeclareMathOperator{\rand}{rand} \DeclareMathOperator{\round}{round} \DeclareMathOperator{\trace}{trace} \DeclareMathOperator{\rank}{rang} \DeclareMathOperator{\dev}{dev} \DeclareMathOperator{\cot}{cot} \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{\valeAbs}[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} \newcommand{\temp}{T} \newcommand{\therCp}{C_{p}} \newcommand{\therFlux}{\vector{q}} \newcommand{\therCondTher}{\lambda} \newcommand{\therSource}{s} \newcommand{\tempVirt}{\delta T} \newcommand{\beamArea}{S} \newcommand{\beamLength}{L} \newcommand{\setC}{\mathbb{C}} \newcommand{\setR}{\mathbb{R}}\]

4. Principe des « doubles Lagrange »

4.1. Méthodologie

La méthode proposée ici est celle mise en œuvre dans le code CASTEM (communication personnelle de Thierry Charras et Pierre Verpeaux). Une présentation intuitive pourrait en être la suivante:

On voit que le problème dualisé a des termes nuls sur la diagonale : ceux correspondants aux degrés de liberté de Lagrange. Cette propriété se remarque aussi sur le Lagrangien : il n’y a pas de termes quadratiques en \(\vector{\lambda}\) .

Cette nullité des termes diagonaux empêche certaines permutations de lignes et colonnes : on ne peut placer un Lagrange avant les degrés de liberté physiques qu’il contraint.

L’idée est alors de décomposer chaque coefficient de Lagrange \(\vector{\lambda}\) en deux parties égales \({\vector{\lambda}}^{1}\) et \({\vector{\lambda}}^{2}\). L’équation \(\tensTwo{B} \vector{u}=\vector{d}\) est alors remplacée par:

\[\begin{split}\left \lbrace \begin{array}{l} \tensTwo{B} \vector{u}-\alpha \left({\vector{\lambda}}^{1}-{\vector{\lambda}}^{2} \right) = \vector{d}\\ \tensTwo{B} \vector{u}+\alpha \left({\vector{\lambda}}^{1}-{\vector{\lambda}}^{2} \right) = \vector{d} \end{array} \right .\end{split}\]

où \(\alpha\) est une constante non nulle.

Montrons l’équivalence de l’ancien problème et du nouveau:

Problème 1 : « simple Lagrange »

Trouver:

\[\begin{split}\left \lbrace \begin{array}{l} \vector{u}\in {\setR}^{n} \\ \vector{\lambda} \in {\setR}^{p} \end{array} \right .\end{split}\]

tels que :

\[\begin{split}(S): \left \lbrace \begin{array}{l} \tensTwo{A} \vector{u} + \tensTwo{B}^{T} \vector{\lambda} = \vector{b} \\ \tensTwo{B} \vector{u} = \vector{d} \end{array} \right .\end{split}\]

On a les équivalences suivantes

\[\begin{split}(S) \iff \left \lbrace \begin{array}{l} {\vector{\lambda}}^{1} = {\vector{\lambda}}^{2}\\ \vector{\lambda} = {\vector{\lambda}}^{1} + {\vector{\lambda}}^{2}\\ \tensTwo{A} \vector{u}+{\tensTwo{B}}^{T}\vector{\lambda} = \vector{b}\\ \tensTwo{B} \vector{u} = \vector{d} \end{array} \right . \iff \left \lbrace \begin{array}{l} \tensTwo{A} \vector{u}+{\tensTwo{B}}^{T} {\vector{\lambda}}^{1}+{\tensTwo{B}}^{T}{\vector{\lambda}}^{2}=\vector{b}\\ \tensTwo{B} \vector{u}-\alpha \vector{\lambda}^{1}+\alpha \vector{\lambda}^{2}=\vector{d}\\ \tensTwo{B} \vector{u}+\alpha \vector{\lambda}^{1}-\alpha \vector{\lambda}^{2}=\vector{d} \end{array} \right .\end{split}\]

\(\alpha\) est une constante \(\ne 0\)

D’où le problème équivalent au précédent:

Problème 2 : « double Lagrange »

Trouver:

\[\begin{split}\left \lbrace \begin{array}{l} \vector{u}\in {\setR}^{n}\\ {\vector{\lambda}}^{1},{\vector{\lambda}}^{2}\in {\setR}^{p}\times {\setR}^{p} \end{array} \right .\end{split}\]

tels que

\[\begin{split}(S^\prime)=: \left \lbrace \begin{array}{l} \tensTwo{A} \vector{u}+\tensTwo{B}^{T}{\vector{\lambda}}^{1}+\tensTwo{B}^{T}{\vector{\lambda}}^{2}=\vector{b}\\ \tensTwo{B} \vector{u}-\alpha \vector{\lambda}^{1}+\alpha \vector{\lambda}^{2} = \vector{d}\\ \tensTwo{B} \vector{u}+\alpha \vector{\lambda}^{1}-\alpha \vector{\lambda}^{2} = \vector{d} \end{array} \right .\end{split}\]

Le nouveau problème peut s’écrire:

\[\tensTwo{K^{\prime}}\vector{X^{\prime}} = \vector{F^{\prime}}\]

avec:

\[\begin{split}\left \lbrace \begin{array}{l} \vector{X^{\prime}} = ( \vector{u},{\vector{\lambda}}_{1},{\vector{\lambda}}_{2})\\ \vector{F^{\prime}} = ( \vector{b} , \vector{d}, \vector{d}) \end{array} \right .\end{split}\]

et

\[\begin{split}\discMatr{K^{\prime}} = \left[\begin{array}{ccc} \tensTwo{A}& {\tensTwo{B}}^{T}& {\tensTwo{B}}^{T}\\ \tensTwo{B}& -{\alpha \tensTwoUnit}& {\alpha \tensTwoUnit}\\ \tensTwo{B}& \alpha \tensTwoUnit & - \alpha \tensTwoUnit \end{array}\right]\end{split}\]

Le problème correspond à rendre extrémale la fonctionnelle:

\[\mathcal{L^{\prime}} \left( \vector{u},{\vector{\lambda}}^{1},{\vector{\lambda}}^{2} \right) = \frac{1}{2} \left( \tensTwo{A} \vector{u} , \vector{u} \right) - \left( \vector{b} , \vector{u} \right) + \left({\vector{\lambda}}^{1}, \tensTwo{B} \vector{u}-\vector{d} \right) + \left({\vector{\lambda}}^{2}, \tensTwo{B} \vector{u}-\vector{d} \right) - \frac{\alpha}{2}\left( {\vector{\lambda}}^{1}-{\vector{\lambda}}^{2},{\vector{\lambda}}^{1}-{\vector{\lambda}}^{2} \right)\]

On peut montrer (Cf Annexe 2) que si on observe une certaine règle de numérotation des inconnues, et en choisissant la constante \(\alpha >0\) , la matrice \(\tensTwo{K^{\prime}}\) vérifie la condition cond1.

Cette règle est la suivante:

Soit une relation de blocage \(\tensTwo{B} \vector{u}-\vector{d}=0\) , il lui correspond 2 multiplicateurs de Lagrange \({\vector{\lambda}}^{1}\) et \({\vector{\lambda}}^{2}\) . Cette relation fait intervenir un certain nombre de degrés de liberté physiques.

Règle \({R}_{0}\) :

Pour chaque relation de blocage, il faut placer \({\vector{\lambda}}^{1}\) avant le premier degré de liberté physique contraint et \({\vector{\lambda}}^{2}\) après le dernier degré de liberté physique contraint. Pour diminuer l’occupation mémoire de la matrice \(\tensTwo{K}\) , il faut chercher à minimiser la largeur de bande. C’est ce que l’on fait dans en « encadrant » les relations « au plus près » : \({\vector{\lambda}}^{1}\) est placé juste avant le premier degré de liberté contraint, \({\vector{\lambda}}^{2}\) est placé juste après le dernier ddl contraint.

4.2. Illustration

Soit un problème à 4 degrés de liberté physiques \(\left( {u}_{1}, {u}_{2}, {u}_{3}, {u}_{4} \right)\) .

Ce système est soumis à 2 conditions:

\[\begin{split}\left \lbrace \begin{array}{l} {a}_{11}{u}_{1} + {a}_{13}{u}_{3} = {b}_{1}\\ {a}_{22}{u}_{2} + {a}_{24}{u}_{4} = {b}_{2} \end{array} \right .\end{split}\]

appelons \(\left({\lambda}_{1}^{1},{\lambda}_{1}^{2} \right)\) les deux degrés de liberté de Lagrange associés à la première condition et \(\left({\lambda}_{2}^{1},{\lambda}_{2}^{2} \right)\) ceux associés à la deuxième condition.

En supposant que les degrés de liberté physiques aient été numérotés dans l’ordre \(\left( {u}_{1}, {u}_{2}, {u}_{3}, {u}_{4} \right)\), la numérotation totale des degrés de liberté retenue est alors \(\left( {\lambda}_{1}^{1},{u}_{1},{\lambda}_{2}^{1},{u}_{2}, {u}_{3} , {\lambda}_{1}^{2},{u}_{4},{\lambda}_{2}^{2} \right)\).

\({\lambda}_{1}^{1}\) et \({\lambda}_{1}^{2}\) encadrent « au plus près » les degrés de liberté contraints (\({u}_{1}\) et \({u}_{3}\)), tandis que \({\lambda}_{2}^{1}\) et \({\lambda}_{1}^{2}\) encadrent « au plus près » les degrés de liberté contraints (\({u}_{2}\) et \({u}_{4}\)).

La technique des « doubles Lagrange » associée à la règle \({R}_{0}\) permet donc de résoudre tout système linéaire physiquement bien posé avec l’algorithme de \(LD{L}^{T}\) sans permutation. La démonstration suppose néanmoins que la matrice \(\tensTwo{A}\) soit symétrique et positive (pas forcément définie).

Remarque

Les hypothèses \(\tensTwo{A}\) symétrique et \(\tensTwo{A}\) positive sont nécessaires pour utiliser \({{LDL}}^{T}\) (ou \({LU}\)) sans permutation comme le montrent les deux contre-exemples suivants:

• \({A}_{1}=\left[\begin{array}{cc}0& 1\\ 1& 0\end{array}\right]\) est symétrique mais non positive

• \({A}_{2}=\left[\begin{array}{cc}0& 1\\ -1& 1\end{array}\right]\) est positive mais non symétrique