# Méthodes de décomposition de domaine et problèmes semi-définis positifs

Jean-Michel Cros

<figure><img src="image_1.png" /></figure>

## To cite this version:

Jean-Michel Cros. Méthodes de décomposition de domaine et problèmes semi-définis positifs. 9e Colloque national en calcul des structures, CSMA, May 2009, Giens, France. hal-01421849

HAL Id: hal-01421849

https://hal.science/hal-01421849v1

Submitted on 23 Dec 2016

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>
# Méthodes de décomposition de domaine et problèmes semi-définis positifs

Jean-Michel CROS

Laboratoire de Mécanique et d'Énergétique d'ÉVRY

Université d'Évry-Val d'Essonne

40 rue du Pelvoux, CE1455 Courcouronnes

91020 ÉVRY Cedex

jean-michel.cros@iup.univ-evry.fr

**Résumé** — Les méthodes de décomposition de domaine FETI-DP et BDDC se caractérisent par un problème auxiliaire, indispensable pour assurer l'extensibilité numérique, qui repose, notamment, sur des contraintes imposées aux coins des sous-domaines. Pour les problèmes d'élasticité linéaire semi-définis positifs, la singularité du problème global induit celle de ce problème auxiliaire. L'extraction de ces singularités permet de construire un filtre des modes rigides du problème global afin de résoudre ce dernier par un algorithme du gradient conjugué préconditionné.

**Mots clés** — décomposition de domaine, mode rigide, méthode itérative.

## 1 Introduction

Les problèmes semi-définis positifs se rencontrent fréquemment en calcul des structures pour des analyses statiques ou dynamiques. Dans le cadre de l'élasticité linéaire, les codes de calculs par éléments finis, grâce à différentes techniques [3], peuvent traiter une déficience liée aux conditions aux limites en déplacement. Toutefois, si on utilise des méthodes de décomposition de domaine, il convient de pouvoir également gérer ce problème alors que la matrice de rigidité assoxiée au problème global n'est pas assemblée et/ou factorisée. Ce problème a été abordé [2, 10] pour la méthode FETI (Finite Element Tearing and Interconnecting). Dans cette communication, on considère les méthodes de décomposition de domaine, FETI-DP (Finite Element Tearing and Interconnecting - Dual Primal) [4] et BDDC (Balancing Domain Decomposition by Constraints) [1, 9] dont le problème auxiliaire ou grossier repose, notamment, sur des contraintes imposées aux coins des sous-domaines. Le problème grossier qui couple l'ensemble des sous-domaines est indispensable pour limiter la dépendance du nombre des itérations pour converger via-à-vis du nombre des sous-domaines et donc assurer l'extensibilité numérique de la méthode. Par ailleurs, il s'avère, sous certaines conditions, que la singularité du problème global induit la singularité du problème grossier. L'extraction des singularités de ce dernier permet de reconstituer par combination linéaire les modes rigides du problème global. Il est alors possible de résoudre le problème global par une méthode de décomposition de domaine qui utilise un algorithme de type gradient conjugué.
## 2 Méthodes de décomposition de domaine

On rappelle quelques éléments des méthodes de décomposition de domaine. On s'intéresse en particulier à la méthode du complément de Schur sous forme primale. Dans le cadre de la théorie de l'élasticité linéaire et d'une discrétisation par éléments finis, on considère le système linéaire suivant :

$$K u = f, \tag{1}$$

où $K$ est une matrice carrée symétrique, supposée définie positive, provenant de la discrétisation par éléments finis d'un problème d'élasticité posé sur un domaine $\Omega$. Le domaine $\Omega$ est découpé en $n_s$ sous-domaines $\Omega^s$ qui ne se recouvrent pas. Notons $K^s$ la matrice de rigidité correspondant au sous-domaine $\Omega^s$, $u^s$ est le vecteur des degrés de liberté de tous les éléments de $\Omega^s$ et $N^s$ une matrice booléenne qui assure le passage entre la numérotation locale et la numérotation globale des degrés de liberté. Par ailleurs, l'interface entre les sous-domaines est notée $\Gamma$ telle que $\Gamma = \cup_{s=1}^{n_s} \Gamma^s$ avec $\Gamma^s = \partial\Omega^s \setminus \partial\Omega$. On distingue alors les degrés de liberté (notés $u_b^s$) de l'interface $\Gamma^s$ d'un sous-domaine des autres degrés de liberté (notés $u_i^s$) du sous-domaine. Sans nuire à la généralité de la méthode, les degrés de liberté sont ordonnés comme suit :

$$K^s = \begin{bmatrix} K_{ii}^s & K_{ib}^s \\ K_{ib}^{sT} & K_{bb}^s \end{bmatrix}, \quad \{u^s\} = \begin{Bmatrix} u_i^s \\ u_b^s \end{Bmatrix}, \quad \{f^s\} = \begin{Bmatrix} f_i^s \\ f_b^s \end{Bmatrix}, \quad \text{et } N^s = [N_i^s N_b^s]. \tag{2}$$

Ainsi, avec ces notations, le système linéaire (1) s'écrit :

$$\begin{bmatrix} K_{ii}^1 & \cdots & 0 & K_{ib}^1 N_b^{1T} \\ \vdots & \ddots & \vdots & \vdots \\ 0 & \cdots & K_{ii}^{n_s} & K_{ib}^{n_s} N_b^{n_s T} \\ N_b^1 K_{ib}^{1T} & \cdots & N_b^{n_s} K_{ib}^{n_s T} & \sum_{s=1}^{n_s} N_b^s K_{bb}^s N_b^{sT} \end{bmatrix} \begin{Bmatrix} u_i^1 \\ \vdots \\ u_i^{n_s} \\ u_b \end{Bmatrix} = \begin{Bmatrix} f_i^1 \\ \vdots \\ f_i^{n_s} \\ \sum_{s=1}^{n_s} N_b^s f_b^s \end{Bmatrix}. \tag{3}$$

Par élimination des degrés de liberté $u_i^s$, on obtient le système linéaire réduit suivant :

$$\left( \sum_{s=1}^{n_s} N_b^s (K_{bb}^s - K_{ib}^{sT} K_{ii}^{s-1} K_{ib}^s) N_b^{sT} \right) u_b = \sum_{s=1}^{n_s} N_b^s \left( f_b^s - K_{ib}^{sT} K_{ii}^{s-1} f_i^s \right). \tag{4}$$

La matrice $S = \sum_{s=1}^{n_s} N_b^s S^s N_b^{sT}$, avec $S^s = K_{ii}^s - K_{ib}^{sT} K_{bb}^{s-1} K_{ib}^s$ est appelée matrice du complément de Schur. Le système (4) est résolu à l'aide d'une méthode itérative de type gradient conjugué préconditionné, puisque la pertinence des méthodes de décomposition de domaine par rapport aux méthodes de sous-structuration est de ne pas assembler la matrice $S$. L'efficacité de la méthode repose en grande partie sur le préconditionneur. Ce dernier doit être peu sensible à la finesse de la discrétisation ainsi qu'à l'augmentation du nombre des sous-domaines et facilement calculable en parallèle. On s'intéresse, en particulier, à la méthode BDDC [1, 9] qui est un préconditionneur de la matrice du complément de Schur. Cette méthode a pour origine la résolution de problèmes discrétisés à l'aide d'éléments de plaque ou de coque. En effet, dans ce cas la sous-structuration engendre des coins qui sont définis comme des nœuds appartenant à plus de deux sous-domaines et/ou comme étant susceptibles d'appartenir à l'intersection de l'interface entre les sous-domaines et de la frontière du domaine d'étude. Or, ces coins induisent des discontinuités qui ralentissent considérablement la convergence. La solution fût de traiter à part et explicitement ces coins en considérant les sous-domaines connectés uniquement par ces derniers. Il s'avère que cette méthode est également très performante pour la résolution de problème 2D et 3D pour peu que dans ce dernier cas, on prenne en compte d'autres contraintes comme la moyenne sur les faces des
sous-domaines.

La description qui va suivre est générale. En effet, on peut utiliser un changement de base [8] pour transformer n'importe quelle contrainte primale (moyenne sur les faces des sous-domaines ou sur les arêtes) de manière explicite de sorte qu'elle soit associée à un degré de liberté qui sera alors identifié comme un coin. On introduit quelques notations, puisqu'il convient de distinguer les degrés de liberté associés aux coins ($u_c$) des autres degrés de liberté ($u_r$) au sein du vecteur des déplacements :

$$ \{u\} = \left\{ \begin{matrix} u_r \\ u_c \end{matrix} \right\} = \left\{ \begin{matrix} u_r^1 \\ \vdots \\ u_r^{n_s} \\ u_c \end{matrix} \right\}. \qquad (5) $$

On introduit alors une matrice booléenne $N_c^s$ qui permet d'obtenir la restriction du déplacement correspondant aux coins dans un sous-domaine :

$$ u_c^s = N_c^{sT} u_c, \qquad (6) $$

et une matrice de pondération [7] définie par :

$$ \sum_{s=1}^{n_s} N_b^s D^s N_b^{sT} = I_{\Gamma}. \qquad (7) $$

Enfin, il faut assurer le passage entre les inconnues $u_i^s$ et $u_b^s$ et les inconnues introduites au niveau de l'expression (5) :

$$ \{u^s\} = \left\{ \begin{matrix} u_r^s = [R_{ri}^s R_{rb}^s] \left\{ \begin{matrix} u_i^s \\ u_b^s \end{matrix} \right\} \\ u_c^s = R_c^s u_b^s \end{matrix} \right\}. \qquad (8) $$

On montre [1] que le préconditionneur de la méthode BDDC, noté $M$, se formule sous forme additive comme suit :

$$ \begin{aligned} M &= M_{\text{local}} + M_{\text{auxiliaire}} \\ &= \sum_{s=1}^{n_s} N_b^s D^s [R_{rb}^{sT} K_{rr}^{s-1} R_{rb}^s] D^s N_b^{sT} \\ &\quad + \sum_{s=1}^{n_s} N_b^s D^s (N_c^s (R_c^s - K_{rc}^{sT} K_{rr}^{s-1} R_{rb}^s))^T S_c^{-1} \sum_{s=1}^{n_s} (N_c^s (R_c^s - K_{rc}^{sT} K_{rr}^{s-1} R_{rb}^s)) D^s N_b^{sT}. \end{aligned} \qquad (9) $$

avec

$$ S_c = \sum_{s=1}^{n_s} N_c^s [K_{cc}^s - K_{rc}^{sT} K_{rr}^{s-1} K_{rc}^s] N_c^{sT} = \sum_{s=1}^{n_s} N_c^s S_c^s N_c^{sT}. \qquad (10) $$

Le premier terme de (9) correspond à la somme pondérée des solutions de problèmes (indépendants) de Neumann posés dans chaque sous-domaine. Le second terme de (9) comprend un problème auxiliaire ou grossier (car ne faisant intervenir qu’un faible nombre d’inconnues) qui couple l’ensemble des sous-domaines par les coins. Il est possible de donner une forme plus synthétique au problème auxiliaire :

$$ M_{\text{auxiliaire}} = \sum_{s=1}^{n_s} N_b^s D^s \Psi^s \left( \sum_{s=1}^{n_s} \Psi^{sT} K^s \Psi^s \right)^{-1} \sum_{s=1}^{n_s} \Psi^{sT} D^s N_b^{sT} \qquad (11) $$
où par analogie avec la condensation statique, $\Psi^s$ correspond à un mode de coin :

$$ \Psi^s = \left\{ -K_{rr}^{s-1} K_{rc}^s N_c^{sT} \middle| N_c^{sT} \right\}. \qquad (12) $$

Dans l'expression (11) les degrés de liberté au sein de la matrice $K^s$ sont ordonnés comme suit : en premier les degrés de liberté ($u_r$) puis les degrés de liberté associés aux coins ($u_c$).

On ne présente pas la méthode FETI-DP [4] car le problème d'interface est semblable [5, 11] au préconditionneur de la méthode BDDC. Le préconditionneur le plus efficace pour la méthode FETI-DP étant la matrice du complément de Schur.

## 3 Détermination des modes rigides

L'exposé du paragraphe précédent suppose que chaque sous-domaine possède un nombre suf-fisant de coins, afin que la matrice $K_{rr}^{s-1}$, qui apparaît (9) dans le terme $M_{\text{local}}$ soit inversible. Dans le cas contraire, il est impératif d'ajouter artificiellement des coins de sorte à lever la singularité. Cela signifie que lorsque l'on traite un problème global mal posé, toute la singularité du problème se concentre dans le problème auxiliaire. Ainsi, en calculant le noyau de la matrice $S_c$, définie en (10), on peut obtenir par combinaison linéaire des modes de coin, les modes rigides de la structure globale. À partir de là, on utilise la technique classique de projection orthogonale pour filtrer les modes rigides.

En ce qui concerne la validation numérique, la méthode proposée est implantée dans un code de calcul parallèle. Le parallélisme est réalisé par échanges de messages via la librairie MPI selon le modèle SPMD. Chaque sous-domaine est affecté à un processeur. La mise en œuvre peut se résumer en trois étapes :

1. Il s'agit de calculer les singularités de la matrice associée au problème auxiliaire. Cette étape est la plus délicate et réclame une attention toute particulière. Deux cas peuvent se présenter, en fonction de la manière dont l'implantation est réalisée : la matrice est assemblée ou non. Dans le premier cas, on peut utiliser différentes techniques [3] : algébrique, géométrique, ... Pour le second cas, on peut utiliser un algorithme de gradient conjugué [10] pour calculer les modes rigides du problème auxiliaire.

2. L'étape suivante consiste simplement à réaliser une projection orthogonale du second mem-bre du système linéaire global par rapport au noyau de la matrice de rigidité globale. À partir de la connaissance des modes rigides, cette étape ne pose aucune difficulté.

3. Enfin, au cours des itérations de l'algorithme du gradient conjugué préconditionné, il faut veiller à filtrer les modes rigides avant et après la résolution du problème auxiliaire. Ceci est nécessaire pour se prémunir de toute dérive numérique (arrondi, ...).

Les méthodes de décomposition de domaine peuvent également s'utiliser en association avec une méthode de recherche de valeurs propres (Lanczos, sous-espace, ...). Pour les problèmes d'analyse modale, on peut s'affranchir du caractère semi-défini positif par un décalage spectral avec la matrice de masse. Toutefois, il n'est pas forcément simple de choisir la valeur adéquate du décalage. Il est possible, à partir de la connaissance des modes rigides de la structure, de modifier [6] légèrement la séquence standard qui engendre le sous-espace formé par les vecteurs $y_r$ :

$$ 1) \quad g_r = \hat{M} y_r \qquad (13) $$

$$ 2) \quad y_{r+1} = K^{-1} g_r \qquad (14) $$
où $\hat{M}$ représente la matrice de masse de la structure globale. Il faut s'assurer qu'à l'issue de l'étape 1, la charge $g_r$ auto-équilibre la structure et qu'à l'issue de l'étape 2 le nouveau vecteur $y_{r+1}$ soit M-orthogonalisé par rapport aux modes rigides de la structure.

## 4 Exemple d'illustration

On considère un problème d'élasticité linéaire 2d. Un bloc est simplement posé, il supporte sur une partie de sa face supérieure une force répartie. Le maillage comprend 1 148 nœuds et 1 080 éléments finis Q1. La figure 1 représente le maillage de la structure ainsi que la décomposition en cinq sous-domaines. Lors de la résolution, l'unique mode rigide est parfaitement identifié. La déformée (amplifiée) est représentée sur la figure 2. Les deux méthodes (BDDC et FETI-DP) ré-solvent le problème en 15 itérations pour une précision donnée de $10^{-6}$. Le coût de mise en œuvre est relativement faible comparativement au coût des autres opérations (deux résolutions de sys-tème linéaire dans chaque sous-domaine par itération). Cet exemple simple démontre l'efficacité de la méthode, il en est de même sur des cas plus complexes.

<figure><img src="image_3.png" /><figcaption>Figure 1 – maillage décomposé en cinq sous-domaines</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 2 – Contrainte équivalente de von Mises</figcaption></figure>

*Remerciements :* les calculs ont été effectués sur la machine parallèle IBM Power 6 de l'Institut du Développement et des Ressources en Informatique Scientifique (IDRIS).
Références

[1] Jean-Michel Cros. A preconditioner for the Schur complement domain decomposition method. *Domain Decomposition Methods in Science and Engineering*, I. Herrera, D. E. Keyes, and O. B. Widlund, 373–380, 2002.

[2] Jean-Michel Cros, Françoise Léné. Parallel iterative methods to solve large-scale eigenvalue problems in structural dynamics. *Domain Decomposition Methods in Sciences and Engineering*, Petter E. Bjørstad and Magne Espedal and David Keyes, 318–324, 1997.

[3] Charbel Farhat, Michel Géradin. On the general solution by a direct method of a large-scale singular system of linear equations : application to the analysis of floating structures. *International Journal for Numerical Methods in Engineering*, John Wiley & Sons, 41(4):675–696, 1998.

[4] Charbel Farhat, Michel Lesoinne, Patrick LeTallec, Kendall Pierson, Daniel Rixen. FETI-DP : a dual-primal unified FETI method - part I : A faster alternative to the two-level FETI method. *International Journal for Numerical Methods in Engineering*, John Wiley & Sons, 50(7):1523–1544, 2001.

[5] Yannis Fragakis. A study on the lumped preconditioner and memory requirements of FETI and related primal domain decomposition methods. *International Journal for Numerical Methods in Engineering*, John Wiley & Sons, 73:1865–1884, 2008.

[6] Michel Géradin, Daniel Rixen. *Théorie des vibrations : Application à la dynamique des structures*, Masson, 1995.

[7] Patrick Le Tallec. Domain decomposition methods in computational mechanics. *Computational Mechanics Advances*, J. Oden, 121–220, 1994.

[8] Jing Li, Olof Widlund. FETI-DP, BDDC, and block Cholesky methods. *International Journal for Numerical Methods in Engineering*, Wiley, 66:250–271, 2006.

[9] Jan Mandel, Clark R. Dohrmann. Convergence of a balancing domain decomposition by constraints and energy minimization. *Numerical Linear Algebra with Applications*, John Wiley & Sons, 10:639–659, 2003.

[10] Daniel Rixen. Dual Schur Complement Method for Semi-Definite Problems. *Contemporary Mathematics*, AMS, 218:341–348, 1998.

[11] Bedrich Sousedik, Jan Mandel. On the Equivalence of Primal and Dual Substructuring Preconditioners. *Computational Linear Algebra with Applications*, ETNA, http://www.citebase.org/abstract?id=oai:arXiv.org:0802.4328, 2008.