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

HAL
open science

# Recherche numérique d'un champ de déplacements en élastoplasticité

Bertrand Mercier

► To cite this version:

Bertrand Mercier. Recherche numérique d'un champ de déplacements en élastoplasticité. [Rapport de recherche] IRIA-RR-071, IRIA. 1974, pp.23. hal-04715990

HAL Id: hal-04715990

https://inria.hal.science/hal-04715990v1

Submitted on 1 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished 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>

Distributed under a Creative Commons Attribution 4.0 International License
RIA
laboria

laboratoire de recherche
en informatique
et automatique

Recherche numérique
d'un champ de déplacements
en élastoplasticité

B.Mercier

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

Institut de Recherche
d'Informatique
et d'Automatique

Domaine de Voluceau
Rocquencourt
78150 - Le Chesnay
Tél. 954.90.20

Rapport de Recherche n° 71
mai 1974
RECHERCHE NUMÉRIQUE
D'UN CHAMP DE DÉPLACEMENTS EN ÉLASTOPLASTICITÉ

Bertrand Mercier

(présenté par R. Glowinski)

RÉSUMÉ

On recherche théoriquement et numériquement le champ de déplacements dans un milieu élasto-plastique soumis à la loi de Hencky.

On donne une démonstration de son existence et des essais numériques dans le cas de la torsion d'une barre cylindrique.

ABSTRACT

In this paper, we shall study, both theoretically and numerically, the displacement field of an elastoplastic medium subject to the law of Hencky.

A proof of existence and numerical results are given for the torsion of a cylindrical beam.
# RECHERCHE NUMERIQUE
D'UN CHAMP DE DEPLACEMENTS EN ELASTOPLASTICITE.

(*) (**) Bertrand MERCIER

## INTRODUCTION.

Pour les problèmes élastique-parfaitement-plastiques (cf. DUVAUT-LIONS [3]), quasi-statiques ou stationnaires, le problème qui est mathématiquement bien posé est celui de la recherche du champ des contraintes (on adopte la loi de Hencky).

Mais pour la résolution numérique, il vaut mieux introduire un **problème aux déplacements** qui apparaît comme le dual du précédent.
C'est un problème non linéaire pour lequel on peut définir des solutions généralisées (cf. EKELAND-TEMAM [4]), mais comme il est sans contraintes, il se prête bien à une résolution par des méthodes de S.O.R. non linéaire (cf. par exemple GLOWINSKI-MARROCCO [5]).

Un problème technique particulier vient de ce que, de par sa provenance, la fonctionnelle à minimiser n'est qu'une seule fois continûment différentiable et très faiblement coercive, ce qui nécessite une adaptation de la méthode de Newton à une dimension.

Dans le cas particulier de la torsion élasto-plastique, on peut démontrer, à l'aide d'un résultat de BREZIS [12], l'existence du champ de déplacement solution dans $H^1(\Omega)$. C'est dans le cadre de cet exemple que l'on a traité les essais numériques. Ils peuvent être comparés à ceux relatifs au Problème aux contraintes, qui est facile à résoudre dans ce cas. (Cf. GLOWINSKI-LIONS-TREMOLIERES [6]).

* IRIA-LABORIA - 78150 ROCQUENCOURT

** Travail effectué alors que l'auteur était stagiaire au Laboratorio di Analisi Numerica del C.N.R. 27100 PAVIA (Italie)

<footnote>(1) S.O.R. : successive over-relaxation (cf. Varga [11]).</footnote>
# I. - DUALITE EN PLASTICITE A LA HENCKY.

## 1.1. Problème initial aux contraintes.

DUVAUT-LIONS (loc. cit.) montrent que pour un milieu élasto-plastique soumis à la loi de Hencky, le problème aux contraintes peut s'écrire

$$ (\mathcal{P}) \qquad \inf_{K \cap M} \frac{1}{2} \mathcal{A}(\varphi, \varphi) - \int_{\Gamma_U} U_i \varphi_{ij} n_j \, d\Gamma $$

où $\Omega$ est l'ouvert de frontière $\Gamma = \Gamma_U \cup \Gamma_F$ limitant le milieu continu,

* $K = \{\varphi_{ij} \in L^2(\Omega) \mid \varphi_{ij} = \varphi_{ji} \text{ et } \mathfrak{F}(\varphi) \le 0 \text{ p.p.}\}$ est le convexe des champs de contraintes plastiquement admissibles.

* $M = \{\varphi_{ij} \in L^2(\Omega) \mid \varphi_{ij} = \varphi_{ji}, \varphi_{ij,j} + f_i = 0 \text{ (1) dans } \Omega \text{ et } \varphi_{ij} n_j = F_i \text{ sur } \Gamma_F\}$ est le convexe des champs de contraintes statiquement admissibles,

* $\mathcal{A}(\varphi, \psi) = \frac{1}{9K_o} \int_{\Omega} \varphi_{kk} \psi_{jj} + \frac{1}{2\mu} \int_{\Omega} \varphi_{ij}^D \psi_{ij}^D$ est la partie bilinéaire de l'énergie élastique, $\mu$ étant le module de cisaillement et $K_o$ celui de rigidité à la compression.

* $U_i \in H^2(\Gamma)$ représente les déplacements imposés sur la partie de la frontière $\Gamma_U$ et $u_{ad} = \{v_i \in H^1(\Omega) \mid v_i = u_i \text{ sur } \Gamma_U\}$.

Comme le problème $(\mathcal{P})$ revient à la minimisation d'une fonctionnelle quadratique coercive sur un convexe : sous l'hypothèse $K \cap M \neq \emptyset$, il existe $\sigma$ solution unique du problème $(\mathcal{P})$.

(1) On fait la convention de sommation de l'indice répété :

La virgule indique une dérivation spatiale partielle :

$$ \varphi_{ij,j} \equiv \sum_{j=1}^{n} \frac{\partial \varphi_{ij}}{\partial x_j} . $$

L'opérateur des déformations $\varepsilon(v)$ que l'on introduira par la suite vaut :

$$ \varepsilon_{ij}(v) = \frac{1}{2}(v_i + v_j). $$
**Remarque :** $(\rho)$ peut aussi s'écrire :

$$
(1.1) \quad \inf_{\varphi \in K} \frac{1}{2} \sigma(\varphi, \varphi) - \int_{\Gamma_U} u_i \varphi_{ij} n_j d\Gamma + \psi_M(\varphi)
$$
où $\psi_M(\varphi)$ est la fonction indicatrice du convexe $M$ qui vaut $\begin{cases} 0 & \text{si } \varphi \in M \\ +\infty & \text{sinon.} \end{cases}$

1.2. **Dualité.**

On part, pour introduire la dualité, de la remarque suivante :

$$
(1.2) \quad \psi_M(\varphi) = \sup_{v \in u_{ad}} \int_{\Omega} (\varphi_{ij,j} + f_i) v \, dx + \int_{\Gamma_F} (F_i - \varphi_{ij} n_j) d\Gamma \\
\text{en effet si :} \quad \varphi_{ij,j} + f_i \neq 0 \text{ sur } \Omega \\
\text{ou si :} \quad F_i - \varphi_{ij} n_j \neq 0 \text{ sur } \Gamma_F,
$$

alors le membre de droite peut être rendu infini par un choix convenable
de v dans u_ad et le "sup" vaut donc + ∞,

et sinon $\varphi \in M$, et le membre de droite vaut toujours zéro,
donc le "sup" vaut zéro, ce qui montre que (1.2) est vrai.

Transformons (1.2) par la formule de Green, on obtient :

$$
(1.3) \quad \psi_M(\varphi) = \int_{\Gamma_U} U_i \varphi_{ij} n_j \, d\Gamma + \sup_{v \in u_{ad}} \left[ (f,v) + \int_{\Gamma_F} F_i v_i \, d\Gamma - (\varepsilon(v), \varphi) \right]
$$

Donc le problème $(\rho)$ peut s'écrire (grâce à (1.1), (1.3))

$$
(1.4) \quad \inf_{\varphi \in K} \sup_{v \in u_{ad}} \left[ \frac{1}{2} \sigma(\varphi, \varphi) - (\varepsilon(v), \varphi) + (f, v) + \int_{\Gamma_F} F_i v_i \, d\Gamma \right].
$$

On peut alors associer à ce problème un problème dual, en permutant
le "Inf" et le "Sup" (cf. [4]).

$$
(\rho*) \quad \sup_{v \in u_{ad}} \left[ \inf_{\varphi \in K} \left[ \frac{1}{2} \sigma(\varphi, \varphi) - (\varepsilon(v), \varphi) \right] + (f,v) + \int_{\Gamma_F} F_i v_i \, d\Gamma \right].
$$

On recherche donc le champ de déplacements comme solution de ce
problème dual, mais malheureusement, on ne peut rien dire de son existence,
en général.
### 1.3. Cas particulier du critère de plasticité de Von Mises.

C'est le cas où :

$$ (1.5) \qquad \mathcal{F}(\sigma_{ij}) = \frac{1}{2} \sigma_{ij}^D \sigma_{ij}^D - k^2 $$

où $\sigma_{ij}^D$ représente le déviateur du tenseur des contraintes, soit :

$$ (1.6) \qquad \sigma_{ij}^D \equiv \sigma_{ij} - \frac{\sigma_{kk}}{3} \delta_{ij} $$

On évalue alors, comme dans DUVAUT-LIONS [3]

$$ (1.7) \qquad \sup_{\tau \in K} \left[ (\varepsilon(v), \nu) - \frac{1}{2} \alpha(\tau, \nu) \right] = \int_{\Omega} \left[ \frac{1}{2} K_o (\mathrm{Div} \, v)^2 + \Phi(\varepsilon^D(v)) \right] dx $$

avec

$$ (1.8) \qquad \Phi(\varepsilon^D) = \begin{cases} \mu \varepsilon_{ij}^D \cdot \varepsilon_{ij}^D & \text{si } \varepsilon_{ij}^D \cdot \varepsilon_{ij}^D \le \frac{k^2}{2\mu^2} \\ k \left( (2\varepsilon_{ij}^D \cdot \varepsilon_{ij}^D)^{\frac{1}{2}} - \frac{k}{2\mu} \right) & \text{si } \varepsilon_{ij}^D \cdot \varepsilon_{ij}^D > \frac{k^2}{2\mu^2} \end{cases} $$

Le problème dual ($\mathcal{P}^*$) devient donc, en changeant de signe :

$$ (1.9) \qquad \inf_{v \in \mathcal{U}_{\text{ad}}} \left\{ \int_{\Omega} \left[ \frac{1}{2} K_o (\mathrm{Div} \, v)^2 + \Phi(\varepsilon^D(\sigma)) \right] dx - (f,v) - \int_{\Gamma_F} F_i v_i d\Gamma \right\} $$

cette fonctionnelle n'est pas coercive sur $\mathcal{U}_{\text{ad}}$, par contre, elle sera infinie à l'infini dans $W^{1,1}(\Omega)$ moyennant certaines hypothèses et on renvoie à DUVAUT-LIONS [3] et à EKELAND-TEMAM [4] pour l'introduction de solutions généralisées.
## II. - APPLICATIONS A LA TORSION D'UNE BARRE CYLINDRIQUE.

### 2.1. **Problème initial.**

C'est un cas particulier où le domaine est un cylindre de section droite $\Omega \subset \mathbb{R}^2$ et de longueur h, dont la surface latérale est libre, et les faces terminales soumises à un déplacement imposé.

Les symétries du problème permettent de montrer (LANCHON [7]) que les seules composantes non nulles du tenseur $\sigma_{ij}$ des contraintes sont :

$$ \sigma_{13} = \sigma_{31}(x_1, x_2) \quad \text{et} \quad \sigma_{23} = \sigma_{32}(x_1, x_2), $$

ce qui permet de se ramener à un problème bidimensionnel en posant :

$$ (2.1) \qquad P_1 = \sigma_{13} \quad \text{et} \quad P_2 = \sigma_{23}. $$

Le problème aux contraintes s'écrit donc :

$$ (2.2) \qquad \inf_{p \in K \cap M} \frac{1}{2C} \int_{\Omega} |p|^2 - \int_{\Omega} \beta \cdot p $$

où

$$ K = \{ p \in (L^\infty(\Omega))^2 \mid |p| \le k \ p \cdot p \} $$

$$ (2.3) \qquad M = \{ p \in (L^\infty(\Omega))^2 \mid \operatorname{div} p = 0 \text{ sur } \Omega \text{ et } p \cdot n = 0 \text{ sur } \Gamma \} $$

avec $C = 2\mu\alpha$, où $\alpha$ est l'angle de torsion unitaire, et $\beta$ est une fonction de $(L^\infty(\Omega))^2$ donnée, qui vaut

$$ (2.4) \qquad \beta = \begin{cases} -x_2 & (\Omega \text{ est supposé borné}) \\ x_1 & \end{cases} $$

(2.2) étant un problème de minimisation d'une fonctionnelle convexe coercive sur un convexe non vide, il admet une solution $\bar{p}^*$.

### 2.2. **Dualité.**

En appliquant les résultats du § 1 (nous sommes dans un cas particulier, avec $f_i = F_i = 0$), on obtiendrait un problème aux déplacements du type de (1.9), coercif seulement dans le quotient de $W^{1,1}(\Omega)$ par les constantes additives $(W^{1,1}(\Omega)/C)$.
On va utiliser les résultats et les notations d'EKELAND-TEMAM [4] (Ch.III § 4), et on introduit :

$$ \bullet \quad V = W^{1,1}(\Omega)/C \qquad \bullet \quad Y = (L^1(\Omega))^2 $$

$$ Y^* = (L^\infty(\Omega))^2 \text{ est le dual de } Y $$

$$ \bullet \quad L : V \to Y \quad \text{application linéaire continue de } V \text{ dans } Y \\ v \mapsto L(v) = \text{grad } v $$

$$ \bullet \quad F(v) = 0 $$

$$ \bullet \quad G(p) = \frac{1}{2} \int_{\Omega} \Phi(p) dX \quad \text{où} \quad \Phi: (L^{\infty}(\Omega))^2 \to L^{\infty}(\Omega) \text{ est définie} $$

presque partout par :

$$ (2.5) \quad \Phi(p) = \begin{cases} C \mid p + \beta \mid^2 & \text{si } C \mid p + \beta \mid \le k \\ 2k \mid p + \beta \mid - \frac{k^2}{C} & \text{sinon} \end{cases} $$

On vérifie que $F$ et $G$ sont deux fonctionnelles convexes s.c.i. et $G$ est même une fois continûment Gateaux-différentiable ; on peut en effet calculer son gradient $G'(p) \in Y^*$ :

$$ (2.6) \quad (G'(p),q) = \int_{\Omega} \Phi'(p) \cdot q \, dX \quad \forall q \in Y $$

avec

$$ (2.7) \quad \Phi'(p) = \begin{cases} C(p+\beta) & \text{si } C|p+\beta| \le k \\ k \frac{p+\beta}{|p+\beta|} & \text{sinon} \end{cases} $$

On peut alors introduire le problème aux déplacements, analogue de (1.9), qui sera cette fois-ci le problème primal :

$$ (2.8) \quad \begin{array}{ll} \text{Trouver } \bar{u} & \text{minimisant} \\ \text{sur } V = W^{1,1}(\Omega)/C & F(u) + G(Lu) = \frac{1}{2} \int_{\Omega} \Phi(\text{gradu}) \, dX \end{array} $$

On va montrer que le problème dual :

$$ (2.9) \quad \text{Inf}_{p^* \in Y^*} F^*(-L^*p^*) + G^*(p^*) $$

n'est autre que le problème initial (2.2).
Il faut pour cela évaluer les fonctions $F^*$ et $G^*$ conjuguées de $F$ et $G$.

Calculons d'abord :

$$
(2.10) \quad F^*(-L^*p^*) = \sup_{v \in V} -\langle L^*p^*, v \rangle = \sup_{v \in Y} -\langle p^*, Lv \rangle
$$

puisque $L^*$ est l'application adjointe de $L$. Donc, en appliquant formellement la formule de Green :

$$
(2.11) \quad F^*(-L^*p^*) = \sup_{v \in V} \int_{\Omega} p^* \operatorname{grad} v = \sup_{v \in V} \int_{\Gamma} v p^* \cdot n d\Gamma + \int_{\Omega} v \operatorname{div} p^* dX \\
\qquad \text{quantit\e qui vaut} \quad \begin{cases} 0 & \text{si } p^* \cdot n = 0 \text{ sur } \Gamma \text{ et } \operatorname{div} p^* = 0 \text{ sur } \Omega \\ +\infty & \text{sinon.} \end{cases}
$$

$$
\text{Evaluons maintenant } G^*(p^*) = \sup_{p \in Y} \langle p^*, p \rangle - G(p).
$$

On a le :

LEMME 2.1 :

$$
(2.12) \quad G^*(p^*) = \begin{cases} \displaystyle \frac{1}{2C} \int_{\Omega} |p^*|^2 - \int_{\Omega} \beta p^* & \text{si } |p^*| \le k \quad \text{p.p.} \\ +\infty & \text{sinon} \end{cases}
$$

et lorsque $|p^*| \le k$, le "sup" est atteint pour :
$p \in Y$

$$
(2.13) \quad p = (1+\lambda) \frac{p^*}{C} - \beta \quad \text{où} \quad \lambda \in L^1(\Omega)
$$

est défini presque partout par :

$$
(2.14) \quad \left\{
\begin{array}{ll}
\lambda = 0 & \text{si } |p^*| < k \\
\lambda \ge 0 & \text{qu}\ell\text{elconque si } |p^*| = k.
\end{array}
\right.
$$

Démonstration.

On a donc $G^*(p^*) = \sup_p \int_{\Omega} (p^* \cdot p - \frac{1}{2} \Phi(p)) dx$.
1°) Montrons tout d'abord que si $|p^*| \ge k_1 > k$ p.p. sur $\Omega_p$
de mesure non nulle, alors $G^*(p^*)$ est infini.

En effet, sur $\Omega_p$, pour $p$ colinéaire à $p^*$ assez grand, on aura :

$$
(2.15) \quad p^* \cdot p - \frac{1}{2} \Phi(p) = |p^*| \cdot |p| - k|p| + \beta|p| + \frac{k^2}{2C} \ge (k_1 - k)|p| - k|\beta| + \frac{k^2}{2C}
$$

et comme $\beta \in (L^\infty(\Omega))^2$, en choisissant

$$
p = \lambda \cdot p^* \text{ p.p. avec } \lambda \begin{cases} \ge 0 \\ = 0 \end{cases} \text{ sur } \begin{cases} \Omega_p \\ \Omega - \Omega_p \end{cases} \quad \lambda \in L^1(\Omega)
$$

on a donc

$$
\int_{\Omega} (p^*p - \frac{1}{2} \Phi(p)) > c_1 \|\lambda\|_L^1 + c_2
$$

qui peut être rendu aussi grand que l'on veut par un choix convenable
de $\lambda$ dans $\frac{L^1}{L^\infty}$.

2°) Lorsque $|p^*| \le k$ p.p. on va montrer que le "sup" est atteint.

Pour cela, on peut essayer de trouver localement $p$ tel que

$$
(2.16) \quad p^* - \Phi'(p) = 0 \quad (\text{dans } \mathbb{R}^2)
$$

on peut donc selon que l'on choisit $C \mid p+\beta \mid \le k$ , choisir $p$ tel que
(d'après (2.7)) :

$$
(2.17) \quad p^* - C(p + \beta) = 0 \quad \Leftrightarrow \quad p = \frac{p^*}{C} - \beta
$$

ou selon que $C \mid(p + \beta)\mid \ge k$ choisir $p$ tel que :

$$
(2.18) \quad p^* - k \frac{p+\beta}{|p+\beta|} = 0
$$

ce qui ne peut arriver que si $|p^*| = k$ et

$$
p = (1+\lambda) \frac{p^*}{C} - \beta \text{ avec } \lambda \ge 0 \text{ quelconque, car } C|p+\beta| \ge k.
$$
Donc, pour résumer, le sup est atteint en tout point si :

$$ 
\begin{array}{l@{\qquad}l}
p = (1+\lambda) \bar{p}_c^* - \beta & \text{avec } \lambda \ge 0 \text{ qqc si } |p^*| = k \\
& \text{et } \lambda = 0 \text{ sinon}
\end{array}
$$

Le sup. sera atteint dans $L^1$ si, de plus, $\lambda \in L^1$, ce qui achève la démonstration du Lemme. ■

Avec ce Lemme, on vérifie que le problème (2.9), dual de (2.8) est bien le problème initial aux contraintes (2.2).

On va maintenant relier l'existence du déplacement v solution du problème (2.8) à celle d'un multiplicateur pour le convexe K (v quant à lui apparaît comme un multiplicateur pour le convexe M) à l'aide de la Proposition suivante : (+)

**PROPOSITION 2.1.**

Les affirmations suivantes sont équivalentes :

(i) $\exists \bar{u} \in W^1, \infty(\Omega)/C$ ($\subset V$) solution du Problème aux déplacements (2.8)

(ii) $\exists \bar{\lambda}$ tel que $(\bar{p}^*, \bar{\lambda})$ soit point-selle de :

$$ \frac{1}{2} \int_{\Omega} |p^*|^2 - \int_{\Omega} p^* \beta + \frac{1}{2} \int_{\Omega} \lambda (|p^*|^2 - 1) \quad \text{sur } M \times (L^\infty(\Omega))^+ $$

On a de plus la relation suivante entre $\bar{u}$ et $\bar{\lambda}$ :

$$ (2.19) \qquad \text{grad} \bar{u} = (1+\bar{\lambda}) \frac{\bar{p}^*}{c} - \beta \quad \text{où } \bar{\lambda} \text{ satisfait à (2.14)} $$

La relation (2.19) est d'ailleurs la loi de Hencky qui relie les déformations et les contraintes (cf. [3]).

**Démonstration.**

D'après EKELAND-TEMAM [4] (Proposition III-4.1), l'existence de $\bar{u} \in V$ solution de (2.8) est équivalente aux relations d'extrêmalité.

$$ (2.20) \qquad \begin{cases} F(\bar{u}) + F^*(-L^*\bar{p}^*) + \langle \bar{p}^*, \bar{L} \rangle = 0 \\ G(\bar{L}u) + G^*(\bar{p}^*) = \langle \bar{p}^*, \bar{L}u \rangle \end{cases} $$

(+) H. LANCHON a montré un résultat analogue [7]
La première est trivialement vérifiée car $\bar{p}^* \in M$.

La seconde équivaut à dire que $\bar{L}$ réalise le maximum de $\sup_{p \in Y} \langle p^*, p \rangle - G(p)$ et donc, d'après le Lemme 2.1, il existe $\bar{\lambda} \in L^1(\Omega)$ vérifiant (2.14) tel que :

$$ (2.21) \qquad \operatorname{grad} \bar{u} = (1 + \bar{\lambda}) \frac{\bar{p}^*}{C} - \beta $$

ce qui montre déjà (2.19).

On voit de plus que $\bar{u} \in W^1, \infty(\Omega) \Leftrightarrow \bar{\lambda} \in L^\infty(\Omega)$

Mais (2.21) peut aussi s'écrire :

$$ (2.22) \qquad \int_{\Omega} \operatorname{grad} \bar{u} \cdot p^* = \frac{1}{C} \int_{\Omega} (1 + \bar{\lambda}) \bar{p}^* \cdot p^* - \int_{\Omega} \beta \cdot p^* \quad \forall p^* \in Y^* $$

équivalent à dire que :

$$ (2.23) \qquad \bar{p}^* \text{ minimise } \frac{1}{2C} \int_{\Omega} |p^*|^2 - \int_{\Omega} \beta \cdot p^* + \frac{1}{2C} \int_{\Omega} \bar{\lambda} - (|p^*|^2 - 1) - \int_{\Omega} \operatorname{grad} \bar{u} \cdot p^* \\ \qquad \text{ sur } Y^*, $$

qui implique à son tour que

$$ (2.24) \qquad \bar{p}^* \text{ minimise } \frac{1}{2C} \int_{\Omega} |p^*|^2 - \int_{\Omega} \beta \cdot p^* + \frac{1}{2C} \int_{\Omega} \bar{\lambda} (|p^*|^2 - 1) \\ \qquad \text{ sur } M $$

(car $p^* \in M \Rightarrow (\operatorname{grad} \bar{u}, p^*) = (0, 0)$, cf. par exemple (2.11)) et (2.24) est bien équivalent à (ii).

Pour montrer la réciproque, il suffit de démontrer que (2.24) $\Rightarrow$ (2.23) car toutes les autres étapes sont des équivalences.

Or, le problème (2.24) $\bar{\lambda}$ étant donné correspond à :

$$ (2.25) \qquad \inf_{p^* \in Y^*} \left[ \frac{1}{2C} \int_{\Omega} (1+\bar{\lambda}) |p^*|^2 - \int_{\Omega} \beta \cdot p^* - \frac{1}{2C} \int_{\Omega} \bar{\lambda} + \sup_{V \in \mathcal{V}} \left[ -\int_{\Omega} \operatorname{grad} v, p^* \right] \right] $$

qui admet $\bar{p}^*$ comme solution.

Son problème dual (obtenu par inversion de l' "Inf" et du "Sup") s'écrit
$$
(2.26) \quad \sup_{v \in V} \left[ \text{Inf}_{p^* \in Y^*} \left[ \frac{1}{2C} \int_{\Omega} (1+\bar{\lambda}) |p^*|^2 - \int_{\Omega} \beta \cdot p^* - \int_{\Omega} \text{grad} \cdot p^* \right] - \frac{1}{2C} \int_{\Omega} \bar{\lambda} \right]
$$

1'"Inf" est atteint pour
$p^*$

$$
|p^*| = C \frac{\beta + \operatorname{grad} v}{1 + \bar{\lambda}}
$$

et vaut $-\frac{C}{2} \int_{\Omega} \frac{|\beta + \text{grad} v|^2}{1 + \bar{\lambda}}$ qui est bien coercif sur $H^1(\Omega)/C$ (réflexif)
$\quad \text{car } \bar{\lambda} \in L^\infty(\Omega)$

$$
\mathrm{car} \quad \bar{\lambda} \in L^{\infty}(\Omega)
$$

et (2.26) admet donc une solution $\bar{u} \in H^1(\Omega)/C \subset W^{1,1}(\Omega)/C$ qui satisfera
justement à (2.23).

On va maintenant montrer que l'existence de $\bar{\lambda}$ résulte de BREZIS [12]

PROPOSITION 2.2.

Si $\Omega$ est sans trôus :

il existe $\bar{u}$ solution du problème aux déplacements (2.8).

Démonstration.

D'après la Proposition 2.1, il suffit de montrer l'existence de $\bar{\lambda}$.
Or, $p^* \in M \Leftrightarrow \exists \psi \in W_0^{1,\infty}(\Omega)$ appelée fonction de contraintes,
(stress function)

telle que :

$$
(2.27) \qquad \left\{
\begin{array}{l}
P_1^* = - \displaystyle\frac{\partial \psi}{\partial x_2} \\
P_2^* = \displaystyle\frac{\partial \psi}{\partial x_1}
\end{array}
\right.
$$

En remarquant que,

$$
(2.28) \quad \int_{\Omega} p^* \beta = \int_{\Omega} x \cdot \text{grad} \, \psi = \int_{\Omega} \psi \quad (\text{en intégrant par parties}).
$$

On en déduit que la condition (ii) est donc équivalente à trouver
($\bar{\psi}$, $\bar{\lambda}$) point-selle de

$$ \frac{1}{2C} \int_{\Omega} |\mathrm{grad} \, \psi|^2 \, dX - \int_{\Omega} \psi \, dX + \frac{1}{2C} \int_{\Omega} \lambda (\mathrm{grad} \, \psi)^2 - 1 \, dX \\ \text{sur } W_o^{1;\infty}(\Omega) \times L^\infty(\Omega)^+ . $$

Or, BREZIS [12] montre que ce point-selle existe sur $H_o^1(\Omega) \times L^\infty(\Omega)^+$, ce qui est plus fort. ■

### Remarque 2.1

On peut étendre sans difficulté au cas général le Lemme et la Proposition 2.1. Mais, pour avoir l'existence, il faudrait obtenir un résultat analogue à celui de BREZIS dans un cas a priori plus compliqué du fait de l'inexistence d'une fonction de contraintes.

## III. - APPROXIMATION ET CONVERGENCE DE LA SOLUTION DU PROBLEME AUX DEPLACEMENTS DE LA TORSION ELASTO-PLASTIQUE.

On va discrétiser le problème aux déplacements (2.8). En dimension finie, les fonctions intégrables étant également de carré intégrable, on approximera $V = W^{1,1}(\Omega)/C$ de la même façon que l'on approximerait $H^1(\Omega)/C$.

On le fera comme dans [5] en prenant des éléments finis linéaires par morceaux sur une triangulation "régulière" du domaine $\Omega$ à N noeuds. De façon à approximer $H^1(\Omega)/C$, on peut rechercher des solutions nulles en un noeud, par exemple central, et ainsi, on définit :

$$ \begin{gathered} v_h \subset V \\ \text{et } (\rho_h) \qquad \begin{array}{l} \operatorname{Inf} J(\psi) \\ \psi \in v_h \end{array} \end{gathered} $$

qui admet une solution unique.

Evidemment, la nature mal posée du problème sur $H^1(\Omega)/C$ rend ardue la question de la convergence qui est traitée par VANENDE [10] dans le cas voisin des hypersurfaces minima.
IV. - RESOLUTION NUMERIQUE DU PROBLEME APPROCHE PAR S.O.R. NON LINEAIRE.

D'après l'analyse des méthodes faite dans GLOWINSKI-MARROCCO [5]
on a envisagé plutôt des méthodes du type surrelaxation ponctuelle non
linéaire. Il faut signaler que la"méthode de Newton"globale (cf. [5])diverge
apparemment mais cela peut fort bien s'expliquer par le fait que la
fonctionnelle à minimiser n'est pas deux fois continûment différentiable.

4.1. <u>Description de l'algorithmè.</u>

On a repris l'algorithmme SNL1 de [8].

On part de $\psi_h^o$ donné et par récurrence :

$\psi_h^n$ connu, on calcule $\psi_h^{n+1}$ coordonnée par coordonnée, par :

$$
\begin{equation} \tag{4.1}
\left\{
\begin{array}{ll}
\bar{\psi}_i^{n+1} & \text{minimise} \\[1ex]
\psi_i^{n+1} = \psi_i^n + \omega (\bar{\psi}_i^{n+1} - \psi_i^n), & 0 < \omega < \omega_M \le 2 \\[1ex]
& i = 1, 2, \dots, N,
\end{array}
\right.
\quad \text{en } \xi \text{ sur } \mathbb{R}
\end{equation}
$$

et ainsi de suite.

La convergence pour $\omega_M > 1$ est montrée dans SCHECHTER [8] lorsque
la fonctionnelle est $C^\infty$ ce qui n'est pas le cas ici.

En revanche, pour $\omega = 1$, on peut appliquer un résultat de CEA-GLOWINSKI [2].

4.2. <u>Minimisation à une variable.</u>

On voit que, à chaque pas élémentaire de l'algorithmme, on devra
minimiser une fonction réelle de classe C¹.

En fait, on résoudra donc l'équation non linéaire du type :

$$
g(\xi) = 0 \quad \text{avec} \quad g(\xi) = \frac{\partial}{\partial \xi} \left( J(\psi_1^{n+1}, \dots, \psi_{i-1}^{n+1}, \xi, \psi_{i+1}^n, \dots, \psi_N^n) \right)
$$

où $g'(\xi)$ est seulement continu par morceaux.
On utilisera donc la méthode de Newton à une dimension en partant de $\xi^0 = \psi_i^n$, quitte à corriger dans les cas défavorables (cf. Fig. 1).

La méthode de correction choisie lors de l'apparition d'une divergence (cf. Fig. 1) caractérisée par :

$$ (\xi^{n+2} - \xi^n)(\xi^{n+1} - \xi^n) \le 0 $$

consiste à prendre alors :

$$ \xi^{n+2} = \frac{\xi^{n+1} + \xi^n}{2} $$

et de recommencer l'algorithme de Newton en partant de cette nouvelle valeur initiale.

<figure><img src="image_4.png" /><figcaption>Figure 1. : Divergence de la méthode de Newton à une dimension.</figcaption></figure>

**Remarque 4.1.** : On a essayé la variante de l'algorithme SNL1 signalée dans SCHECHTER [8] et qui consiste à ne faire à chaque pas qu'une seule itération de la méthode de Newton, mais on a obtenu une divergence.

**Remarque 4.2.** : On constate une accélération de la convergence si l'on part de la solution élastique plutôt que de 0. Mais celle-ci étant presque aussi longue à calculer que la solution plastique, en pratique, on commence l'algorithme par un balayage élastique.
### 4.3. Résultats numériques.

On a choisi $\Omega = ]0,1[^2$.

$$k = 1 \qquad C = 2\mu\alpha = 10$$

On a pris une triangulation uniforme sur un maillage régulier de pas $h = \frac{1}{m}$ (Fig. 2).

Les résultats ont tous été obtenus avec l'algorithme SNL1, le point de départ choisi étant 0. (Figures 3., 4., 5.).

#### Remarque 4.3 : Test d'arrêts.

Pour la méthode de Newton à une dimension, on a choisi :

$$\frac{|\xi^{n+1} - \xi^n|}{\xi^n} \le \varepsilon$$

et pour les itérations de surrelaxation :

$$\frac{\sum_{i=1}^{N} |\psi_{i}^{n+1} - \psi_{i}^{n}|}{\sum_{i=1}^{N} |\psi_{i}^{n+1}|} \le \varepsilon$$

en choisissant comme dans GLOWINSKI-MARROCCO (loc.cit.) $\varepsilon = 10^{-5}$.

#### Remarque 4.4. :

**Forme des équipotentielles** trouvées pour le champ de déplacements : elles ont à peu près la même forme que dans le cas élastique (Fig. 2), mais les déformations étant évidemment plus grandes dans le cas plastique, les valeurs seront plus grandes.

<figure><img src="image_5.png" /><figcaption>Fig. 2 : Forme des équipotentielles.</figcaption></figure>

On voit qu'il y a donc 4 axes d'antisymétrie : les deux axes de coordonnées et les deux bissectrices. (Cf. Sokolnikoff [9] solution exacte de la torsion élastique d'un cylindre carré).
Figure 3 : Résultats pour $h = \frac{1}{6e}$ . Triangulation : zones plastiques (en grisé)

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

Les valeurs du déplacement normal $\psi$ sont inscrites sur le bord.
On peut comparer avec les résultats du problème élastique :

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

On voit donc que les résultats en élasto-plasticité sont un peu moins symétriques, les déplacements sont plus grands et les équipotentialles, quoique gardant une forme analogue, se déforment, les déplacements les plus grands étant enregistrés maintenant au coeur des zones plastiques.
**Figure 4** : Variation du nombre d'itérations en fonction de $\omega$ when $h = \frac{1}{6e}$

Comparaison entre les problèmes élastiques et élastoplastiques

<figure><img src="image_8.png" /></figure>
**Figure 5.** : Résultats pour $h = \frac{1}{12^e}$ ;

Triangulation : zones élastiques et plastiques

169 variables

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

en grisé : zones plastiques

Résultats obtenus par S.O.R. non linéaires avec :

$\omega_{op}t = 1.78$ 49 itérations
CPU time : 74" d'IBM 360/44 (1)

(~ 1' environ de CII 10070)

<footnote>(1) Les calculs ont été effectués sur l'IBM 360/44 du Centre de Calcul de l'Université de Pavie.</footnote>
<figure><img src="image_10.png" /><figcaption>Figure 6. Variation du nombre d'itérations en fonction de $\omega$.</figcaption></figure>

$$h = \frac{1}{12^e} \qquad \omega_{opt} = 1.78 (49 \text{ itérations}).$$
**Figure 7.** Résultats pour $h = \frac{1}{24^e}$.

Triangulation, zones élastiques et plastiques.

Quelques valeurs du déplacement normalisé $\psi$.

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

## 625 variables

Résultat obtenu par S.O.R. non linéaire avec $\omega = 1.88$ (optimal)
en 87 itérations
en 7' 15" de 360/44

On avait obtenu pour
$\omega = 1.78$ 190 itérations
$\omega = 1.83$ 146 itérations
$\omega = 1.88$ 87 itérations
$\omega = 1.93$ 143 itérations
$\omega = 1.98$ divergence.
## CONCLUSION :

Les résultats peuvent être comparés avec, par exemple, ceux de BOURGAT [1], qui résout le problème aux contraintes : on obtient des zones plastiques et des temps de calcul comparables.

On voit donc que la méthode proposée aboutit à un moyen efficace (rapide et peu encombrant) de calcul des déplacements dans un problème élasto-plastique. Pour la torsion élasto-plastique, la méthode présente l'avantage de s'étendre sans modifications au cas des ouverts avec trous.

Elle présente surtout l'avantage de pouvoir se généraliser au cas de l'élasto-plasticité à deux ou trois dimensions, car, comme on l'a vu, le problème aux déplacements reste d'ordre 2 (*) alors que le problème aux contraintes devient, quant à lui, d'ordre 4, ce qui rend sa résolution très difficile.

---- § ----

<footnote>(*) Problème d'ordre 2 : ne fait intervenir dans la fonctionnelle à minimiser que des produits de dérivées partielles premières.</footnote>
BIBLIOGRAPHIE

[1] J.F. BOURGAT Thèse de 3ème cycle, Paris (1972).

[2] J. CEA-R. GLOWINSKI RAIRO (1974).

[3] G. DUVAUT-J.L. LIONS Les Inéquations en Mécanique et en Physique
Dunod (1972).

[4] I. EKELAND-R. TEMAM Analyse convexe et Problèmes variationnels,
Dunod (1973)

[5] R. GLOWINSKI-A. MARROCCO Analyse numérique du champ magnétique d'un
alternateur par éléments finis.

[6] R. GLOWINSKI-J.L. LIONS-R. TREMOLIERES : Approximation numérique des
solutions des inéquations en Mécanique et en Physique
Dunod (1975)

[7] H. LANCHON Thèse Paris (1972).

[8] S. SCHECHTER Relaxation methods for convex problems. SIAM J. on
Numerical Analysis 5. pp. 601-612.

[9] SOKOLNIKOFF Mathematical theory of elasticity. McGRAWHILL (1956)

[10] VANENDE Thèse de 3ème cycle, Paris.

[11] VARGA Matrix iterative analysis. Prentice Hall.

[12] H. BREZIS Multiplicateur de Lagrange en torsion élasto-plastique. Arch. Rat. Mech. Analysis (1972)49-32

---- § ----