# Sur l'approximation, par éléments finis d'ordre un, et la résolution, par pénalisation-dualité, d'une classe de problèmes de Dirichlet non linéaires

Roland Glowinski, Americo Marrocco

► To cite this version:

Roland Glowinski, Americo Marrocco. Sur l'approximation, par éléments finis d'ordre un, et la ré-solution, par pénalisation-dualité, d'une classe de problèmes de Dirichlet non linéaires. [Rapport de recherche] IRIA-RR-115, IRIA. 1975, pp.57. hal-04716085

HAL Id: hal-04716085

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

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_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
RIA
laboria

Institut de Recherche
d'Informatique
et d'Automatique

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

laboratoire de recherche
en informatique
et automatique

ARCHIVES

Sur l'approximation,
par éléments finis d'ordre un,
et la résolution,
par pénalisation-dualité,
d'une classe de problèmes
de Dirichlet non linéaires

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

R.Glowinski, A.Marrocco.

Rapport de Recherche n° 115.

Avril 1975.
# SUR L'APPROXIMATION, PAR ÉLÉMENTS FINIS D'ORDRE UN, ET LA RÉSOLUTION, PAR PÉNALISATION-DUALITÉ, D'UNE CLASSE DE PROBLEMES DE DIRICHLET NON LINÉAIRES

R. Glowinski, A. Marrocco

## RESUME

On étudie dans cet article l'approximation par éléments finis et la résolution, par pénalisation-dualité, du problème de Dirichlet non linéaire

$$ - \nabla \cdot (|\nabla u|^{p-2} \nabla u) = f \quad \text{dans } \Omega, u = 0 \quad \text{sur } \partial\Omega \\ \text{avec } 1 < p < +\infty. $$

On montre la possibilité de généraliser les méthodes utilisées à d'autres problèmes non linéaires ; des résultats numériques sont présentés.

## ABSTRACT

We study in this paper the numerical analysis of the following non linear Dirichlet problem

$$ - \nabla \cdot (|\nabla u|^{p-2} \nabla u) = f \quad \text{in } \Omega, u = 0 \quad \text{ou } \partial\Omega \\ \text{with } 1 < p < +\infty. $$

The approximation is achieved by a finite element method and error estimates are given ; the approximated problem is solved by an iterative algorithm combining S.O.R., penalty and duality. Extension to other non linear problems is given and numerical results are presented.

Ce travail est une version étendue d'un article à paraître en 1975 dans la Revue Française d'Automatique, Informatique et Recherche Opérationnelle.
# I. INTRODUCTION

On étudie (1) l'approximation et la résolution numérique de certains problèmes aux limites non linéaires de type monotone; on montre également que la méthode itérative utilisée pour résoudre les problèmes approchés permet de traiter d'autres problèmes non linéaires de la physique et de la mécanique.

## REMARQUE 1.1.

Parmi les travaux consacrés à l'analyse numérique des problèmes non linéaires monotones on citera CIARLET -SCHULTZ- VARGA [1], BREZIS- SIBONY [2], SIBONY [3], LEBAUD [4] et la bibliographie de ces travaux. On trouvera en particulier dans [1] un cadre et des résultats généraux pour l'approximation des problèmes non linéaires monotones; il semble cependant que les résultats des N°6 et 8 améliorent ceux que l'on obtiendraient en appliquant directement le théorème 2.1. de [1,N°2].

# 2. LE PROBLEME CONTINU

Soit $\Omega$ un ouvert borné de $\mathbb{R}^N$, de frontière $\Gamma$ régulière et $p$ avec $1<p<+\infty$. On note $V$ l'espace $W_0^1, P(\Omega)$, Banach réflexif pour la norme $\|v\|_1 = (\int_{\Omega} |\nabla v|^p dx)^{1/p}$. Soit $p'$ l'exposant conjugué de $p$, d'où $\frac{1}{p} + \frac{1}{p'} = 1$ ; on note $V'$ le dual $W^{-1}, P'(\Omega)$ de $V$ et $\|\cdot\|_\star$ la norme duale de $\|\cdot\|_1$.

Soit $f \in V'$, on considère le problème aux limites non linéaire

(2.1) $$ \left\{ \begin{array}{ll} -\nabla \cdot (|\nabla u|^{p-2} \nabla u) = f & \text{dans } \Omega \\ u=0 & \text{sur } \Gamma \cdot \end{array} \right. $$

On notera dans ce qui suit A l'opérateur de $V \to V'$ défini par

$$ (2.2) \qquad A(v) = -\nabla \cdot (|\nabla v|^{p-2} \nabla v) ; $$

<footnote>(1) Le présent travail constitue un développement de GLOWINSKI-MAROCCO [5],[6].</footnote>
on remarque que pour $p = 2$ on a $A = -\Delta$.

Soit $J : V \to R$ définie par

$$ (2.3) \qquad J(v) = \frac{1}{p} \int_{\Omega} |\nabla v|^p dx - \langle f, v \rangle $$

où, dans (2.3), $\langle \cdot, \cdot \rangle$ dénote la forme bilinéaire de la dualité entre $V'$ et $V$; le gradient $J'$ de $J$ en $v$ est alors défini par

$$ (2.4) \qquad J'(v) = A(v) - f $$

et il en résulte qu'il y a équivalence entre résoudre (2.1) dans $V$ et résoudre le problème de minimisation

$$ (2.5) \qquad \left\{ \begin{array}{ll} J(u) \le J(v) & \forall v \in V \\ u \in V. & \end{array} \right. $$

On a $J$ **strictement convexe, différentiable**, avec $\lim_{\|v\|_1 \to +\infty} J(v) = +\infty$, il en résulte (Cf, par exemple LIONS [7,Ch.1]) que (2.5) admet une solution et une seule caractérisée par

$$ (2.6) \qquad \left\{ \begin{array}{ll} \int_{\Omega} |\nabla u|^p \chi_{\{|\nabla u| > \Delta v\}} dx = \langle f, v \rangle & \forall v \in V \\ u \in V & \end{array} \right. $$

ce qui est une formulation **variationnelle** de (2.1).

### REMARQUE 2.1.

L'opérateur $A$ et une variante$^{(2)}$ du problème (2.1) interviennent (avec $1<p<2$) dans certains modèles mathématiques décrivant les déformations mécaniques de la glace (Cf M.C. PELISSIER [8] et la bibliographie de ce travail).
### 3. LE PROBLEME APPROCHE

#### 3.1. DEFINITION DU PROBLEME APPROCHE

On suppose (3) (pour simplifier l'exposé) que $\Omega$ est un polygone de $R^2$ ; soit $\mathcal{L}_h$ une triangulation finie de $\Omega$ vérifiant

$$ (3.1) \qquad T \subset \bar{\Omega} \quad \forall T \in \mathcal{L}_h, \cup_{T \in \mathcal{L}_h} T = \bar{\Omega} $$

$$ (3.2) \qquad \left\{ \begin{array}{l} T, T' \in \mathcal{L}_h \rightarrow T \cap T' = \emptyset \text{ ou } T \text{ et } T' \\ \text{n'ont qu'un sommet commun ou n'ont qu'un côté commun.} \end{array} \right. $$

On prend h égal à la longueur du plus grand côté des $T \in \mathcal{L}_h$ et on approche $V$ par

$$ (3.3) \qquad V_h = \{v_h | v_h \in C^\circ(\bar{\Omega}), v_h|_T = 0, v_h|_T \in P_l, \forall T \in \mathcal{L}_h\} $$

avec $P_l$ = espace des polynômes de degré $\le 1$ ; on a $V_h \subset V$ et on approche (2.1) (et (2.6)) par le problème en dimension finie

$$ (3.4) \qquad \left\{ \begin{array}{l} \int_{\Omega} |\nabla u_h|^{p-2} \nabla u_h \cdot \nabla v_h \, dx = \langle f, v_h \rangle \quad \forall v_h \in V_h \\ u_h \in V_h. \end{array} \right. $$

### REMARQUE 3.1.

On s'est limité dans ce travail à une approximation par éléments finis d'ordre un, dans la mesure où, pour $p \ne 2$, on peut avoir $u \notin V \cap W^2, P(\Omega)$ même avec des données **très régulières** (on renvoie au N°8.1 ci-après pour de tels exemples) ; dans ces conditions l'intérêt d'utiliser des éléments d'ordre supérieur à **un** semble limité.

<footnote>(3) On adapterait sans peine ce qui suit au cas où $\Omega$ est un polyèdre de $R^3$.</footnote>
## 3.2. RESOLUBILITE DU PROBLEME APPROCHE

$V_h$ étant un sous-espace fermé de $V$ on a la

### PROPOSITION 3.1.

Le problème approché (3.4) admet une solution et une seule.

## 4. CONVERGENCE DES SOLUTIONS APPROCHEES

On va démontrer, dans ce numéro, un théorème de convergence pour les solutions approchées, les estimations de l'erreur d'approximation feront l'objet du N°6 ci-après.

L'espace $\mathcal{H}(\Omega)$ étant dense dans $V$ on en déduit le

### THEOREME 4.1.

Si les angles de $t_h$ sont bornés inférieurement, uniformément en $h$, par $\theta_0 > 0$, on a, u étant la solution de (2.1)

$$ (4.1) \qquad \lim_{h \to 0} \|u_h - u\|_1 = 0. $$

### DEMONSTRATION

On fait $v_h = u_h$ dans (3.4) d'où

$$ (4.2) \qquad \|u_h\|_1^{p-1} \le \|f\|_{\star} = \|u\|_1^{p-1}. $$

Soit $\varphi \in \mathcal{H}(\Omega)$, on note $\Pi_h \varphi$ l'interpolé "linéaire" de $\varphi$ sur $t_h$, soit

$$ \left\{ \begin{array}{ll} \Pi_h \varphi \in V_h & \\ \\ \Pi_h \varphi(P) = \varphi(P) & \forall P \text{ sommet de } \partial t_h ; \end{array} \right. $$

sous les hypothèses sur $t_h$ de l'énoncé, il est standard (Cf. par exemple CIARLET -RAVIART [9]) que
$$ \| \Pi_h \varphi - \varphi \|_1 \to 0 \quad \text{lorsque} \quad h \to 0. $$

Remarquant que (3.4) équivaut à

$$ \begin{cases} J(u_h) \le J(v_h) & \forall v_h \in V_h \\ u_h \in V_h & \end{cases} $$

on en déduit

$$ (4.3) \qquad J(u_h) \le J(\Pi_h \varphi) \quad \forall \varphi \in \mathcal{H}(\Omega). $$

Compte tenu de (4.2) on peut extraire de $(u_h)_h$ une sous-suite, encore notée $(u_h)_h$, telle que

$$ (4.4) \qquad u_h \to u^* \quad \text{dans } V \text{ faible ;} $$

par *semi-continuité inférieure faible* on en déduit à la limite dans (4.3)

$$ J(u^*) \le J(\varphi) \quad \forall \varphi \in \mathcal{H}(\Omega), $$

d'où par densité de $\mathcal{H}(\Omega)$ dans $V$

$$ J(u^*) \le J(v) \quad \forall v \in V; $$

on a donc $u^* = u$ et *toute* la suite $u_h$ converge vers $u$.

La convergence forte résulte de la convergence faible et de

$$ (4.5) \qquad \lim_{h \to 0} J(u_h) = J(u); $$

en effet (4.5) implique $\lim_{h \to 0} \|u_h\|_1 = \|u\|_1$ et $V$ étant *uniformément convexe* pour $\|\cdot\|_1$, la convergence faible et la convergence de la norme entraînent la convergence forte. ■

### REMARQUE 4.1.

On donnera au N°7 une démonstration du Théorème 4.1 utilisant les relations d'ellipticité et de continuité des N°5.1 et 5.2. ■
# 5. PROPRIETES D'ELLIPTICITE ET DE CONTINUITE

## 5.1. RESULTATS D'ELLIPTICITE

Soit $z = (z_1, z_2) \in \mathbb{R}^2$, $|z| = \sqrt{z_1^2 + z_2^2}$; on a le

### LEMME 5.1.

Si $p \ge 2$ on a

$$ (5.1) \qquad (|z|^{p-2} z - |y|^{p-2} y, z-y)_{\mathbb{R}^2} \ge \alpha |z-y|^p $$

avec $\alpha > 0$, indépendant de $y$ et $z$.

### DEMONSTRATION

Si $y = z$, (5.1) est vérifiée pour tout $\alpha$; de même si $y$ et/ou $z$ est nul, (5.1) est vérifiée pour $0 < \alpha \le 1$. Supposant $y \ne z$, $y \ne 0$ on va étudier la fonction $\phi_1$ définie par

$$ \phi_1 (y,z) = \frac{(|z|^{p-2} z - |y|^{p-2} y, z-y)_{\mathbb{R}^2}}{|z-y|^p} . $$

La fonction $\phi_1$ étant homogène (d'ordre zéro), on peut toujours supposer $|y|=1$ et même $y = (1,0)$. Remarquant que $\lim_{|z| \to +\infty} \phi_1(y,z) = 1$ et que $\phi_1(y,z) > 0 \quad \forall y,z,y \ne z$, il suffit, pour démontrer (5.1), de vérifier que $\lim_{z \to y} \phi_1(y,z) > 0$.

On pose

$$ z_1 = 1 + \rho \cos \theta, z_2 = \rho \sin \theta $$

d'où

$$ |z-y| = \rho, \quad |z| = \sqrt{1+2\rho\cos\theta + \rho^2}, $$

$$ (|z|^{p-2} z - |y|^{p-2} y, z-y)_{\mathbb{R}^2} = 1 + (1+2\rho\cos\theta+\rho^2)^2 - (1+\rho\cos\theta)(1+(1+2\rho\cos\theta+\rho^2)^2)^{\frac{p}{2}-1} ; $$

un développement limité à l'ordre 2 montre que si $y = (1,0)$ on a
$$
(5.2) \quad \phi_1(y,z) = \frac{1+(p-2)\cos^2\theta + \varepsilon(\rho,\theta)}{\rho^{p-2}}
$$

avec $\lim_{\rho \to 0} \varepsilon(\rho, \theta) = 0$ uniformément en $\theta$. Il en résulte que si $y \neq 0$,

$$
\lim_{z \to y} \phi_1(y,z) = +\infty \text{ si } p>2 \text{ et } 1 \text{ si } p=2; \text{ ceci achève de démontrer (5.1).}
$$

On en déduit immédiatement la propriété d'ellipticité de la

PROPOSITION 5.1.

$$
\text{On a } \forall u, v \in V \text{ et si } p \ge 2
$$

$$
(5.3) \quad \langle A(v) - A(u), v - u \rangle \ge \alpha \|v - u\|_1^p
$$

avec $\alpha > 0$, indépendant de $u$ et $v$.

LEMME 5.2.

$$
\text{Si } 1 < p \le 2 \text{ on a } \forall y, z \in \mathbb{R}^2
$$

$$
(5.4) \quad \left\{
\begin{array}{@{}l@{}}
\text{avec } \alpha > 0 \text{ indépendant de } y \text{ et } z. \\
\frac{(|z|+|y|)^{2-p} (|z|^{p-2} z - |y|^{p-2} y, z-y)}{|z-y|^2} \ge \alpha |z-y|^2
\end{array}
\right.
$$

DEMONSTRATION

\text{On pose}
$$
\phi_2(y,z) = \frac{(|z|+|y|)^{2-p} (|z|^{p-2}z - |y|^{p-2}y, z-y)}{|z-y|^2}
$$

il suffit là encore de vérifier que $\lim_{z \to y} \inf \phi_2(y,z) > 0$ pour $y = (1,0)$.

Avec les notations du lemme précédent on montre

$$
(5.5) \quad \phi_2(y,z) = 2^{2-p} (1 - (2-p)\cos^2\theta + \varepsilon(\rho,\theta))
$$

avec $\lim_{\rho \to 0} \varepsilon(\rho, \theta) = 0$ uniformément en $\theta$ ; ceci démontre le lemme puisque

$1 < p \le 2.$

On en déduit la
## PROPOSITION 5.2.

On a $\forall u,v \in V \text{ et si } 1 < p \le 2$

$$ (5.6) \qquad (\|v\|_1 + \|u\|_1)^{2-p} \langle A(v) - A(u), v-u \rangle \ge \alpha \|v-u\|_1^2 $$

avec $\alpha > 0$, indépendant de $u$ et $v$.

## DEMONSTRATION

Compte tenu du lemme 5.2 on a p.p.

$$ (5.7) \qquad (|\nabla v| + |\nabla u|)^{2-p} (|\nabla v|^{p-2} \nabla v - |\nabla u|^{p-2} \nabla u) \cdot \nabla (v-u) \ge \alpha |\nabla (v-u)|^2 \quad \forall u, v \in V. $$

Procédant comme dans TARTAR [10, 2ème partie] à propos d'un problème voisin on élève les deux membres de (5.7) à la puissance $\frac{p}{2}$, d'où après intégration

$$ (5.8) \qquad \alpha \frac{p}{2} \int_{\Omega} |\nabla (v-u)|^p dx \le \int_{\Omega} (|\nabla v| + |\nabla u|)^{(2-p)\frac{p}{2}} ((|\nabla v|^{p-2} \nabla v - |\nabla u|^{p-2} \nabla u) \cdot \nabla (v-u))^{\frac{p}{2}} dx ; $$

on a

$$ (|\nabla v| + |\nabla u|)^{(2-p)\frac{p}{2}} \in L^{\frac{2}{p-2}}(\Omega) $$

et

$$ ((|\nabla v|^{p-2} \nabla v - |\nabla u|^{p-2} \nabla u) \cdot \nabla (v-u))^{\frac{p}{2}} \in L^p(\Omega), $$

on déduit donc de (5.8) et de l'inégalité de Hölder

$$ \alpha^{\frac{p}{2} \|v-u\|_1^p} \le \langle A(v) - A(u), v-u \rangle^{\frac{p}{2}} \| |\nabla v| + |\nabla u| \|^{\frac{p}{2}(2-p)} \chi_{L^p(\Omega)} $$

soit encore

$$ \alpha \|v-u\|_1^2 \le \langle A(v) - A(u), v-u \rangle \|v-u\|^{2-p} \chi_{L^p(\Omega)} \le \langle A(v) - A(u), v-u \rangle (\|v\|_1 + \|u\|_1)^{2-p}, $$

c.q.f.d.
## 5.2. RESULTATS DE CONTINUITE

### LEMME 5.3.

$$ 
\begin{gather*}
\text{Si } p \ge 2 \text{ on a } \forall y, z \in \mathbb{R}^2 \\
\qquad | |z|^{p-2}z - |y|^{p-2}y | \le \beta |z-y| (|z| + |y|)^{p-2} \\
(5.9) \qquad \text{avec } \beta > 0, \text{ indépendant de } y \text{ et } z.
\end{gather*}
$$

### DEMONSTRATION

On pose

$$ \phi_3(y,z) = \frac{| |z|^{p-2}z - |y|^{p-2}y |}{|z-y| (|z| + |y|)^{p-2}} $$

Il suffit pour démontrer (5.9) de vérifier que

$$ \limsup_{z \to y} \phi_3(y,z) < +\infty, \text{ avec } y = (1,0). $$

En procédant comme dans les lemmes précédents on démontre

$$ (5.10) \qquad \phi_3(y,z) = \frac{1}{2^{p-2}} (1+(p-2)\cos^2\theta + \epsilon(\rho,\theta)) $$

$$ \text{avec } \lim_{\rho \to 0} \epsilon(\rho, \theta) = 0 \text{ uniformément en } \theta, \text{ d'où le résultat.} $$

### PROPOSITION 5.3.

$$ 
\begin{gather*}
\text{Si } p \ge 2 \text{ on a, } \forall u, v \in V \\
(5.11) \qquad \|A(v)-A(u)\|_{\star} \le \beta \|v-u\|_{1} (\|u\|_{1} + \|v\|_{1})^{p-2}
\end{gather*}
$$

$$ \text{avec } \beta > 0, \text{ indépendant de } u,v. $$

### DEMONSTRATION

$$ \text{On a} \\
\|A(v)-A(u)\|_{\star} = \sup_{\substack{w \in V \\ w \neq 0}} \frac{|\langle A(v)-A(u), w \rangle|}{\|w\|_{1}} $$

avec
$$
\langle A(v) - A(u), w \rangle = \int_{\Omega} (|\nabla v|^{p-2} \nabla v - |\nabla u|^{p-2} \nabla u) \cdot \nabla w \, dx,
$$

d'où, d'après successivement la proposition 5.3, l'inégalité de Hölder
et l'inégalité triangulaire

$$
\begin{align*}
|\langle A(v) - A(u), w \rangle| &\le \int_{\Omega} |\nabla v|^{p-2} \nabla v - |\nabla u|^{p-2} \nabla u| \, |\nabla w| \, dx \\
&\le \beta \int_{\Omega} |\nabla (v-u)| (|\nabla v| + |\nabla u|)^{p-2} |\nabla w| \, dx \\
&\le \beta \|v-u\|_1^p |\nabla v| + |\nabla u| \|w\|_1^{p-2} \|w\|_1 \\
&\le \beta \|v-u\|_1 (\|v\|_1^p + \|u\|_1^p)^{p-2} \|w\|_1
\end{align*}
$$

c.q.f.d.

■

LEMME 5.4.

$$
\text{Si } 1 < p \le 2 \text{ on a, } \forall y, z \in \mathbb{R}^2 \\
(5.12) \qquad \left\{
\begin{array}{l}
| |z|^{p-2}z - |y|^{p-2}y | \le \beta |z-y|^{p-1}, \\
\text{avec } \beta > 0, \text{ indépendant de } y \text{ et } z.
\end{array}
\right.
$$

DEMONSTRATION

On pourrait utiliser les lemmes précédents en posant

$$
z' = |z|^{p-2} z \quad \text{d'où} \quad z = |z'|^{p'-2} z'
$$

avec p' exposant conjugué de p. On a préféré raisonner directement ;
on pose

$$
\phi_4 (y, z) = \frac{\left| |z|^{p-2} z - |y|^{p-2} y \right|}{|z-y|^{p-1}}
$$

il suffit de montrer que $\lim_{z \to y} \sup \phi_4(y,z) < +\infty$ pour $y = (1,0)$ ; on a,

avec les notations des lemmes précédents,

$$
\phi_4 (y, z) = \rho^{2-p} (1 - (2-p) \cos^2 \theta + \epsilon (\rho, \theta))
$$

avec $\lim_{\rho \to 0} \epsilon(\rho, \theta) = 0$ uniformément en $\theta$, d'où le résultat.

■
### PROPOSITION 5.4.

Si $1 < p \le 2$ on a $\forall u,v \in V$

$$ (5.13) \qquad \|A(v)-A(u)\|_{\star} \le \beta \|v-u\|_{1}^{p-1} $$

avec $\beta > 0$, indépendant de $u,v$.

### DEMONSTRATION

Variante immédiate de celle de la proposition 5.3.

### REMARQUE 5.1.

Tous les résultats précédents valent avec les mêmes $\alpha$ et $\beta$ si $\Omega \subset \mathbb{R}^N$, $N > 2$.

## 6 - ESTIMATIONS DE L'ERREUR SI $u \in V \cap W^2, P(\Omega)$.

On notera $\|v\|_2$ la norme de $V$ dans $V \cap W^2, P(\Omega)$ ; u et $u_h$ étant respectivement solutions de (2.6) et (3.4) on a le

### LEMME 6.1.

$$ \|u_h\|_1 \le \|u\|_1. $$

### DEMONSTRATION

On a

$$ \langle A(u_h), u_h \rangle = \langle f, u_h \rangle = \langle A(u), u_h \rangle, $$

d'où

$$ \|u_h\|_1^{p-1} \le \|A(u)\|_{\star} $$

or

$$ \|A(u)\|_{\star} = \|u\|_1^{p-1} $$

d'où le résultat.

Relativement à l'erreur d'approximation on a le
## THEOREME 6.1.

Sous les hypothèses sur $\ell_h$ du théorème 4.1, on a

$$ (6.1) \quad \|u_h - u\|_1 \le C \|u\|_2^{\frac{p-2}{p-1}} \|u\|_1^{\frac{p-1}{p-2}} \frac{1}{p-1} \quad \forall u \in V \cap W^2, P(\Omega), p \ge 2 $$

$$ (6.2) \quad \|u_h - u\|_1 \le C \|u\|_2^{\frac{3-p}{3}} \|u\|_1^{\frac{2-p}{3-p}} \frac{1}{h^{3-p}} \quad \forall u \in V \cap W^2, P(\Omega), 1 < p \le 2 $$

avec $c$ indépendant de $h$ et $u$.

## DEMONSTRATION

On a

$$ (6.3) \quad \langle A(u), u_h - u \rangle = \langle f, u_h - u \rangle $$

$$ (6.4) \quad \langle A(u_h), v_h - u_h \rangle = \langle f, v_h - u_h \rangle \quad \forall v_h \in V_h. $$

Par addition on déduit de (6.3), (6.4)

$$ \text{soit encore } \langle A(u_h) - A(u), u_h - u \rangle = \langle A(u_h), v_h - u \rangle - \langle f, v_h - u \rangle, \\ \langle A(u_h) - A(u), u_h - u \rangle = \langle A(u_h) - A(u), v_h - u \rangle + \langle A(u) - f, v_h - u \rangle $$

$$ \text{mais } A(u) = f, \quad \text{d'où} $$

$$ (6.5) \quad \langle A(u_h) - A(u), u_h - u \rangle = \langle A(u_h) - A(u), v_h - u \rangle \quad \forall v_h \in V_h. $$

L'ouvert $\Omega$ étant bidimensionnel borné à frontière $\Gamma$ lipschitzienne, on a

$$ W^2, P(\Omega) \subset C^\circ(\bar{\Omega}) \quad \forall p > 1 $$

avec injection continue.

Si $u \in V \cap W^2, P(\Omega)$ on peut donc prendre pour $v_h$ l'interpolé linéaire $\Pi_h u$ de $u$ sur $\ell_h$ i.e.

$$ (6.6) \quad \left\{ \begin{array}{ll} \Pi_h u \in V_h & \\ \Pi_h u(P) = u(P) & \forall P \text{ sommet de } \ell_h; \end{array} \right. $$
d'après [9] ou STRANG-FIX [11], on a sous les hypothèses sur $\mathcal{L}_h$ faites ci-dessus

$$ (6.7) \qquad \|\Pi_h u - u\|_1 \le C_o h \|u\|_2 $$

avec $C_o$ indépendant de $h$ et $u$ (mais non de $\theta_o$).

Si $p \ge 2$ on déduit des propositions 5.1, 5.3, du lemme 6.1, de (6.5) et de (6.7)

$$ \begin{aligned} \alpha \|u_h - u\|_1^p &\le \|A(u_h) - A(u)\|_{\star} \|\Pi_h u - u\|_1 \\ &\le \beta C_o (\|u\|_1 + \|u_h\|_1)^{p-2} \|u_h - u\|_1 \|u\|_2^h \\ &\le 2^{p-2} \beta C_o \|u\|_1^{p-2} \|u_h - u\|_1 \|u\|_2^h \end{aligned} $$

$$ \text{d'où} \qquad \|u_h - u\|_1 \le (2^{p-2} \frac{\beta}{\alpha} C_o) \frac{1}{p-1} \|u\|_1^{\frac{p-2}{p-1}} \|u\|_2^{\frac{1}{p-1}} \frac{1}{h^{p-1}}, $$

d'où le résultat (6.1) en posant

$$ C = (2^{p-2} \frac{\beta}{\alpha} C_o)^{\frac{1}{p-1}} . \qquad \blacksquare $$

Si $1 < p \le 2$, on déduit des propositions 5.2, 5.4, du lemme 6.1, de (6.5) et de (6.7)

$$ \begin{aligned} \alpha \|u_h - u\|_1^2 &\le 2^{2-p} \|u\|_1^{2-p} \|A(u_h) - A(u)\|_{\star} \|\Pi_h u - u\|_1 \\ &\le 2^{2-p} \beta C_o \|u\|_1^{2-p} \|u_h - u\|_1^{p-1} \|u\|_2^h, \end{aligned} $$

$$ \text{d'où} \qquad \|u_h - u\|_1 \le (2^{2-p} \frac{\beta}{\alpha} C_o)^{\frac{1}{3-p}} \|u\|_1^{\frac{2-p}{3-p}} \|u\|_2^{\frac{1}{3-p}} \frac{1}{h^{3-p}}, $$

d'où (6.2) en posant

$$ C = (2^{2-p} \frac{\beta}{\alpha} C_o)^{\frac{1}{3-p}} . \qquad \blacksquare $$
7. UNE DEUXIEME DEMONSTRATION DU THEOREME 4.1.

On va compléter la remarque 4.1 ; on a montré au N°6 que

$$
\langle A(u_h) - A(u), u_h - u \rangle = \langle A(u_h) - A(u), v_h - u \rangle \quad \forall v_h \in V_h .
$$

Si $p \ge 2$ on a

$$
(7.1) \quad \alpha \|u_h - u\|_1^{p-1} \le 2^{p-2} \beta \|u\|_1^{p-2} \|v_h - u\|_1 \quad \forall v_h \in V_h .
$$

Soit $\varphi \in L^1(\Omega)$, on a $\|\Pi_h \varphi - \varphi\|_1 \le C_o h^2 \|\varphi\|_2$ et en posant $v_h = \Pi_h \varphi$ dans (7.1)
on en déduit

$$
(7.2) \quad \alpha \|u_h - u\|_1^{p-1} \le 2^{p-2} \beta \|u\|_1^{p-2} (C_o h^2 \|\varphi\|_2 + \|\varphi - u\|_1) \quad \forall \varphi \in L^1(\Omega);
$$

compte tenu de $\overline{\Omega}(\Omega) = V$, (7.2) implique la convergence forte. On raisonnerait de même pour $1 < p \le 2$. $\quad \blacksquare$

8 - ESTIMATIONS DE L'ERREUR SI $u \in V \cap W^{s,p}(\Omega')$, $1 \le s \le 2$.

8.1. MOTIVATION DE CETTE ETUDE.

On peut avoir $u \notin W^2, P(\Omega)$, même avec des données très régulières comme le montre l'exemple ci-après ; on prend $\Omega = ]-1,+1[$, $<f,v> = C \int_{\Omega} v \, dx$, d'où le problème

$$
(8.1) \quad \left\{
\begin{array}{l}
\displaystyle \int_{-1}^{+1} \left| \frac{du}{dx} \right|^{p-2} \frac{du}{dx} \frac{dv}{dx} \, dx = C \int_{-1}^{+1} v \, dx \quad \forall v \in V \\
\\
u \in V.
\end{array}
\right.
$$

La solution de (8.1), également solution de

$$
\begin{cases}
- \frac{d}{dx} \left( \left| \frac{du}{dx} \right|^{p-2} \frac{du}{dx} \right) = C \\
u(-1) = u(1) = 0,
\end{cases}
$$

est donnée par

$$
(8.2) \quad u(x) = \frac{p-1}{p} C \frac{1}{p-1} (1 - |x|^{\frac{p}{p-1}}) ;
$$
on a $u \in V \cap W^2, P(\Omega)$ si $1 < p < \frac{3+\sqrt{5}}{2}$, cette propriété de régularité ne vaut plus si $p \ge \frac{3+\sqrt{5}}{2}$, auquel cas on vérifiera en utilisant la définition d'Aronszajn des espaces $W^s, P(\Omega)$ pour s quelconque (Cf NEČAS [12], Ch.2) que

$$u \in V \cap W^s, P(\Omega) \quad \text{pour tout } s < 1 + \frac{2p-1}{p-1}$$

(on a $1+\frac{2p-1}{p-1} \le 2$ si $p \ge \frac{3+\sqrt{5}}{2}$).

## 8.2. ESTIMATIONS DE L'ERREUR D'APPROXIMATION.

Pour obtenir des estimations de l'erreur d'approximation lorsque $u \in V \cap W^s, P(\Omega)$, $1 \le s \le 2$, on va procéder par *interpolation non linéaire* en utilisant le théorème 2 de TARTAR [13]; on va d'abord démontrer le

### LEMME 8.1.

Soit $u, w \in V$, quelconques, $u_h, w_h$ les solutions approchées correspondant respectivement à $f = A(u)$ et $f = A(w)$; on a alors

$$(8.3) \quad \|w_h - u_h\|_1 \le \left(\frac{\beta}{\alpha}\right)^{\frac{1}{p-1}} \|w-u\|_1^{\frac{1}{p-1}} (\|w\|_1 + \|u\|_1)^{\frac{p-2}{p-1}} \quad \text{si } p \ge 2$$

$$(8.4) \quad \|w_h - u_h\|_1 \le \frac{\beta}{\alpha} \|w-u\|_1^{p-1} (\|w\|_1 + \|u\|_1)^{2-p} \quad \text{si } 1 < p \le 2.$$

### DEMONSTRATION

On a

$$\langle A(u_h), w_h - u_h \rangle = \langle A(u), w_h - u_h \rangle$$

$$\langle A(w_h), u_h - w_h \rangle = \langle A(w), u_h - w_h \rangle$$

d'où par addition

$$(8.5) \quad \langle A(w_h) - A(u_h), w_h - u_h \rangle = \langle A(w) - A(u), w_h - u_h \rangle.$$

Si $p \ge 2$ on déduit de (8.5)

$$\alpha \|w_h - u_h\|_1^{p-1} \le \beta (\|w\|_1 + \|u\|_1)^{p-2} \|w-u\|_1$$
soit

$$
\| w_h - u_h \|_1 \le \left( \frac{\beta}{\alpha} \right)^{\frac{1}{p-1}} \left( \| w \|_1 + \| u \|_1 \right)^{\frac{p-2}{p-1}} \| w - u \|_1^{\frac{1}{p-1}}
$$

■

Si $1 < p \le 2$ on déduit de (8.5)

$$
\alpha \|w_h - u_h\|_1 \le \beta (\|w\|_1 + \|u\|_1)^{2-p} \|w-u\|_1^{p-1}
$$

d'où (8.4).

■

Notant $\|v\|_s$ la norme de $V$ dans $V \cap W^{s,P}(\Omega)$, on a les estimations du

THEOREME 8.1.

Sous les hypothèses sur $\tilde{t}_h$ du théorème 4.1., on a pour

$1 \leqq s \leqq 2$

$$
(8.6) \quad \| u_h - u \|_1 \le C \| u \|_1^{\frac{p-2}{p-1}} \| u \|_s^{\frac{1}{p-1}} \frac{s-1}{h^{p-1}} \quad \forall u \in V \cap W^{s,p}(\Omega), p \ge 2
$$

$$
(8.7) \quad \|u_h - u\|_1 \le C \|u\|_1^{\alpha'} \|u\|_s^{\beta'} h^{\gamma'} \quad \forall u \in V \cap W^{s,p}(\Omega), \ 1 < p \le 2
$$

avec C indépendant de h et u ; dans (8.7) on a

$$
\alpha' = \frac{(2-p)[(2-s)+(s-1)(p-1)]}{(2-s)+(s-1)(p-1)(3-p)}, \quad \beta' = \frac{p-1}{(2-s)+(s-1)(p-1)(3-p)}, \quad \gamma' = (s-1)\beta'.
$$

DEMONSTRATION

Soit $u \in V$, $u_h$ la solution du problème approché correspondant à $f = A(u)$ ; on définit l'application $T : V \to V$ par

$$
(8.8) \qquad T(u) = u_h - u.
$$

Si $p \ge 2$, (6.1) équivaut à

$$
(8.9) \quad \| T(u) \|_1 \le C \| u \|_2^{\frac{p-2}{p-1}} \| u \|_1^{\frac{p-1}{p-2}} \frac{1}{h^{p-1}} \quad \forall u \in V \cap W^{2,p}(\Omega);
$$

on a par ailleurs

$$
\|w-u\|_1 = \|w-u\|_1^{\frac{p-2}{p-1}} \|w-u\|_1^{\frac{1}{p-1}} \le \|w-u\|_1^{\frac{p-2}{p-1}} (\|w\|_1 + \|u\|_1)^{\frac{1}{p-1}}
$$
d'où

$$
(8.10) \quad \|w-u\|_1 \le \|w-u\|_1^{p-1} (\|w\|_1 + \|u\|_1)^{\frac{p-2}{p-1}} .
$$

Compte tenu de (8.10), on déduit de (8.3)

$$
(8.11) \quad \|T(w)-T(u)\| \le (1+\left(\frac{\beta}{\alpha}\right)^{\frac{1}{p-1}}) \|w-u\|_1^{\frac{1}{p-1}} (\|w\|_1+\|u\|_1)^{\frac{p-2}{p-1}}.
$$

Par application du théorème 2 de [13], on déduit de (8.9) et (8.11),
*par interpolation non linéaire* entre V et V ∩ W²,ₚ(Ω) que

$$
\|T(u)\|_1 \le C \|u\|_1^{\frac{p-2}{p-1}} \|u\|_s^{\frac{1}{p-1}} h^{\frac{s-1}{p-1}} \quad \forall u \in V \cap W^{s,p}(\Omega)
$$

avec C indépendant de h et u (et distinct du C de (8.9)) ;
on a donc démontré (8.6).

Si $1 < p \le 2$, (6.2) équivaut à

$$
(8.12) \quad \|T(u)\|_1 \le C \|u\|_2^{\frac{1}{3-p}} \|u\|_1^{\frac{2-p}{3-p}} h^{\frac{1}{3-p}} \quad \forall u \in V \cap W^{2,p}(\Omega);
$$

par ailleurs

$$
\|w-u\|_1 = \|w-u\|_1^{p-1} \|w-u\|_1^{2-p} \le \|w-u\|_1^{p-1} (\|w\|_1 + \|u\|_1)^{2-p}
$$

d'où

$$
(8.13) \quad \|w-u\|_1 \le \|w-u\|_1^{p-1} (\|w\|_1 + \|u\|_1)^{2-p} .
$$

Compte tenu de (8.13) on déduit de (8.4)

$$
(8.14) \quad \|T(w)-T(u)\|_1 \le (1+\frac{\beta}{\alpha}) \|w-u\|_1^{\alpha-1} (\|w\|_1+\|u\|_1)^{2-\alpha}.
$$

Par interpolation non linéaire entre V et V ∩ W²,ₚ(Ω) on déduit de (8.12)
et (8.14)

$$
\|T(u)\|_1 \le C \|u\|_1^{\alpha'} \|u\|_s^{\beta'} h^{\gamma'} \quad \forall u \in V \cap W^{s,p}(\Omega)
$$

avec C indépendant de h et u et α', β', γ' comme dans l'énoncé du théorème.
## REMARQUE 8.1.

On peut mettre $\alpha', \beta', \gamma'$ sous les formes (plus parlantes) suivantes :

$$ \alpha' = \frac{\frac{2-p}{3-p} \frac{1 + \frac{(2-s)}{(s-1)(p-1)}}{1 + \frac{(2-s)}{(s-1)(p-1)(3-p)}} $$

$$ \beta' = \frac{1}{(3-p)(s-1)} \frac{1}{1 + \frac{(2-s)}{(s-1)(p-1)(3-p)}} $$

$$ \gamma' = \frac{\frac{1}{3-p}}{1 + \frac{(2-s)}{(s-1)(p-1)(3-p)}} $$

# 9 - UNE METHODE ITERATIVE DE RESOLUTION DU PROBLEME APPROCHE

## 9.1. GENERALITES

Il y a équivalence entre le problème approché (3.4) et le problème de minimisation

$$ (9.1) \qquad \begin{cases} J(u_h) \le J(v_h) & \forall v_h \in V_h \\ u_h \in V_h & \end{cases} $$

avec

$$ (9.2) \qquad J(v_h) = \frac{1}{p} \int_{\Omega} |\nabla v_h|^p dx - \langle f, v_h \rangle. $$
Les différentes méthodes de surrelaxation non linéaires décrites dans GLOWINSKI - MARROCCO [14] sont pratiquement inopérantes (4), appliquées à la résolution de (9.1), si $p < 1.5$ et $p > 10$ ; de même la méthode d'opérateur auxilliaire décrite dans [2] ne fonctionne en toute généralité appliquée à (9.1), que pour $p$ égal à 2. La méthode de Newton-Raphson ne peut être appliquée pour $p < 2$ puisque dans ce cas l'application $v_h \to J(v_h)$ n'est pas de classe $C^2$ ; si $p > 2$ la matrice hessienne de $v_h \to J(v_h)$ n'est en général que semi-définie positive, donc, là également, la méthode de Newton-Raphson ne peut être appliquée. Quant aux méthodes de pas factionnaires ou de direction alternées (A.D.I.) l'absence de directions de dérivation privilégiées, ou de décomposition particulière de l'opérateur A, font que ces méthodes semblent non applicables, directement tout au moins.

Pour pallier ces difficultés on va augmenter le nombre de variables, tout en simplifiant la structure non linéaire de (9.1), (9.2), en posant

$$ (9.3) \qquad z_h = \nabla v_h ; $$

on découple ensuite $z_h$ et $v_h$ par **pénalisation et dualisation simultanées** de $z_h - \nabla v_h = 0$ (suivant un principe dû à HESTENES [15]). En effet la **pénalisation utilisée seule conduit à un problème dont la solution est distincte** de celle du problème initial, et en est d'autant plus voisine que le paramètre de pénalisation est petit ; cette approche présente un inconvénient majeur qui est que le problème pénalisé est d'autant plus mal conditionné -toutes choses égales par ailleurs - que le paramètre de pénalisation est petit. Par ailleurs la dualité utilisée seule conduit à un Lagrangien, soit $\mathcal{L}$ dans ce qui suit (Cf. N°9.2), coercif en $z_h$ et linéaire, donc non coercif, en $v_h$ ; ceci a pour effet de rendre impossible la résolution de (9.1) par une méthode de type Uzawa (Cf. ARROW - HURWICZ-UZAWA [16]) appliquée à $\mathcal{L}$.

## 9.2. UNE PROBLEME D'OPTIMISATION EQUIVALENT. RESULTATS DE POINT-SELLE.

Remarquant que $\nabla v_h$ est constant sur tout $T \in \mathcal{L}_h$, on note $\chi_T$ la fonction caractéristique de T et on définit

$$ (9.4) \qquad L_h = \{z_h \mid z_h = \sum_{T \in \mathcal{T}_h} z_T \chi_T, z_T \in \mathbb{R}^2\} $$

(4) même lorsqu'elles sont théoriquement convergentes.
puis

$$
j : V \times L^p(\Omega) \times L^p(\Omega) \to \mathbb{R} \quad \text{par}
$$

$$
(9.5) \quad j(v,z) = \frac{1}{p} \int_{\Omega} |z|^p dx - \langle f, v \rangle .
$$

Le problème (9.1) équivaut à

$$
(9.6) \quad \left\{
\begin{array}{ll}
j(u_h, y_h) \le j(v_h, z_h) & \forall (v_h, z_h) \in V_h \times L_h, \nabla v_h \cdot z_h = 0 \\
(u_h, y_h) \in V_h \times L_h & \nabla u_h \cdot y_h = 0
\end{array}
\right.
$$

et à (9.6) on peut associer le langrangien $\mathcal{L}$ défini par

$$
\mathcal{L}(v,z;\mu) = j(v,z) - \int_{\Omega} \mu \cdot (z - \nabla v) dx
$$

qui est bien du type mentionné au N° 9.1.

Par pénalisation et dualisation de $\nabla v - z = 0$ on est conduit à introduire
(pénalisation)

$$
j_{\varepsilon} = j + \frac{1}{2\varepsilon} \|z - \nabla v\|_{L^2(\Omega)}^2, \quad \varepsilon > 0
$$

puis (dualisation) le langrangien $\mathcal{L}_{\varepsilon}$ défini par

$$
\mathcal{L}_{\varepsilon}(v,z;\mu) = j_{\varepsilon}(v,z) - \int_{\Omega} \mu \cdot (z - \nabla v) dx .
$$

On démontre alors (soit directement, soit en utilisant ROCKAFELLAR [ 17, N°28])
la

PROPOSITION 9.1.

$$
\mathcal{L}_{\varepsilon} \text{ admet sur } V_h \times L_h \times L_h \text{ un point-selle de la forme } \\
(u_h, \nabla u_h, \lambda_h) ; \quad \text{ solution de (9.1) et } \lambda_h = |\nabla u_h|^{p-2} \nabla u_h.
$$

9.3. DESCRIPTION DE LA METHODE ITERATIVE

Compte tenu de la proposition 9.1 on va appliquer
l'algorithme d' Uzawa à $\mathcal{L}$, sur $V_h \times L_h \times L_h$, d'où

$$
(9.7) \qquad \lambda_h^0 \in L_h, \text{ donné}
$$
$\lambda_h^n$ connu on détermine $u_h^n, y_h^n, \lambda_h^{n+1}$ par

$$
(9.8) \quad \mathcal{L}_{\varepsilon}(u_h^n, y_h^n; \lambda_h^n) \le \mathcal{L}_{\varepsilon}(v_h^n, z_h^n; \lambda_h^n) \quad \forall v_h^n, z_h^n \in L_h^n ; u_h^n \in V_h^n, y_h^n \in L_h^n
$$

$$
(9.9) \quad \lambda_h^{n+1} = \lambda_h^n - \rho_n (y_h^n - \nabla u_h^n), \quad \rho_n > 0.
$$

9.4. CONVERGENCE DE LA METHODE ITERATIVE

On considère la fonctionnelle $z \to \frac{1}{p} \int_{\Omega} |z|^p dx$ de $L^p(\Omega) \times L^p(\Omega) \to \mathbb{R}$ ; cette fonctionnelle est *différentiable* et son gradient $B$ est défini par

$$
(9.10) \quad B(z) = |z|^{p-2} z ,
$$

on remarquera que $B(z) \in L^p'(\Omega) \times L^p'(\Omega)$, $p'$ exposant conjugué de $p$.

On note

$$
(B(z)-B(y),z-y) = \int_{\Omega} (B(z)-B(y)) \cdot (z-y) dx \\
|z|_{p} = \|z\|_{L^{p}(\Omega) \times L^{p}(\Omega)} = \left( \int_{\Omega} (z_{1}^{2}+z_{2}^{2})^{\frac{p}{2}} dx \right)^{\frac{1}{p}}.
$$

On vérifiera les propriétés suivantes

$$
(9.11) \quad (B(z)-B(y),z-y) \ge \alpha |z-y|_p^p \quad \text{si } p \ge 2
$$

$$
(9.12) \quad (B(z)-B(y),z-y) \ge \alpha \frac{|z-y|_p^2}{(|z|_p^2 + |y|_p^2)^{2-p}} \quad \text{si } 1 < p \le 2,
$$

$\alpha$ ayant les mêmes valeurs que dans les propositions 5.1 ($p \ge 2$) et 5.2 ($1<p<2$).

Pour démontrer la convergence de l'algorithme (9.7)-(9.9) on aura besoin
des
## PROPOSITIONS 9.2.

$$ (3.13) \qquad \left\{ \begin{array}{l} \text{Le couple } (u_h, y_h) \text{ solution de (9.6) est caractérisé par} \\[1em] \qquad (B(y_h), z_h) + \frac{1}{\epsilon}(y_h - \nabla u_h, z_h - \nabla v_h) = \langle f, v_h \rangle + (\lambda_h, z_h - \nabla v_h) \\[1em] \qquad \forall v_h \in V_h, \forall z_h \in L_h \\[1em] \qquad y_h - \nabla u_h = 0 \\[1em] \qquad (u_h, y_h; \lambda_h) \in V_h \times L_h \times L_h. \end{array} \right. $$

### DEMONSTRATION

Immédiate compte tenu de la proposition 9.1.

## PROPOSITION 9.3.

$$ (9.14) \qquad \left\{ \begin{array}{l} \text{Le couple } (u_h^n, y_h^n) \text{ solution de (9.8) est caractérisé par} \\[1em] \qquad (B(y_h^n), z_h^n) + \frac{1}{\epsilon}(y_h^n - \nabla u_h^n, z_h^n - \nabla v_h^n) = \langle f, v_h^n \rangle + (\lambda_h^n, z_h^n - \nabla v_h^n) \\[1em] \qquad \forall v_h^n \in V_h^n, \forall z_h^n \in L_h^n, \\[1em] \qquad u_h^n \in V_h^n, y_h^n \in L_h^n. \end{array} \right. $$

### DEMONSTRATION

Immédiate.

Compte tenu de ces résultats on a le

## THEOREME 9.1.

$$ \begin{array}{l} \underline{Si} \ 0 < r_0 \le \rho_n \le r_1 < \frac{2}{\epsilon}, \ \underline{on a, lorsque} \ n \to +\infty \ \underline{et} \\[1em] \forall \lambda_h^o \in L_h \\[1em] \qquad u_h^n \to u_h \\[1em] \qquad y_h^n \to \nabla u_h \\[1em] \qquad \lambda_h^n \to |\nabla u_h|^{p-2} \nabla u_h \end{array} $$

$u_h$ étant la solution du problème (9.1).
## DEMONSTRATION

On note $\bar{y}_h^n = y_h^n - y_h$, $\bar{u}_h^n = u_h^n - u_h$, $\bar{\lambda}_h^n = \lambda_h^n - \lambda_h$.

Par soustraction entre (9.13) et (9.14) on déduit

$$ (9.15) \qquad (B(y_h^n) - B(y_h), \bar{y}_h^n) + \frac{1}{\epsilon} |\bar{y}_h^n - \nabla \bar{u}_h^n|_2^2 = (\bar{\lambda}_h^n, \bar{y}_h^n - \nabla \bar{u}_h^n); $$

par ailleurs $y_h - \nabla u_h = 0$ implique

$$ (9.16) \qquad \lambda_h = \lambda_h - \rho_n (y_h - \nabla u_h) \quad \forall n $$

d'où par soustraction de (9.16) à (9.9)

$$ (9.17) \qquad \bar{\lambda}_h^{n+1} = \bar{\lambda}_h^n - \rho_n (\bar{y}_h^n - \nabla \bar{u}_h^n). $$

De (9.17) on déduit

$$ (9.18) \qquad |\bar{\lambda}_h^n|_2^2 - |\bar{\lambda}_h^{n+1}|_2^2 = 2\rho_n (\bar{\lambda}_h^n, \bar{y}_h^n - \nabla \bar{u}_h^n) - \rho_n^2 |\bar{y}_h^n - \nabla \bar{u}_h^n|_2^2 $$

et (9.15), (9.18) impliquent

$$ (9.19) \qquad |\bar{\lambda}_h^n|_2^2 - |\bar{\lambda}_h^{n+1}|_2^2 = 2\rho_n (B(y_h^n) - B(y_h), \bar{y}_h^n) + \rho_n \left(\frac{2}{\epsilon}\rho_n\right) |\bar{y}_h^n - \nabla \bar{u}_h^n|_2^2. $$

De (9.19) on déduit que si

$$ (9.20) \qquad 0 < r_o \le \rho_n \le r_l < \frac{2}{\epsilon} $$

la suite $|\bar{\lambda}_h^n|_2$ est décroissante; étant bornée inférieurement par zéro cette suite est convergente d'où

$$ (9.21) \qquad \lim_{n \to +\infty} (|\bar{\lambda}_h^n|_2^2 - |\bar{\lambda}_h^{n+1}|_2^2) = 0. $$

Compte tenu de (9.19), (9.20) et (9.21) impliquent

$$ (9.22) \qquad \lim_{n \to +\infty} (B(y_h^n) - B(y_h), y_h^n - y_h) = 0 $$

$$ (9.23) \qquad \lim_{n \to +\infty} |\bar{y}_h^n - \nabla \bar{u}_h^n|_2 = 0. $$
Si $p \ge 2$ on a, d'après (9.11)

$$ (B(y_h^n) - B(y_h), y_h^n - y_h) \ge \alpha |y_h^n - y_h|^p, \quad \alpha > 0, $$

(9.22) implique donc $y_h^n \to y_h$ et, compte tenu de ce résultat,

(9.23) implique que $u_h^n \to u_h$.

La convergence de $\lambda_h^n$ résulte alors de celle des suites $u_h^n, y_h^n$ et de (9.14).■

Si $1 < p \le 2$, on a , d'après (9.12)

$$ (9.24) \qquad (B(y_h^n) - B(y_h), y_h^n - y_h) \ge \frac{|y_h^n - y_h|^2}{(|y_h^n|^p + |y_h|^p)^{2-p}}, \quad \alpha > 0. $$

Compte tenu de (9.22), (9.24) la suite $y_h^n$ est bornée et ces deux relations impliquent alors que $y_h^n \to y_h$ ; la convergence de $u_h^n$ vers $u_h$ résulte alors de la convergence de $y_h^n$ et de (9.23). ■

### REMARQUE 9.1.

Dans certains cas l'algorithme (9.7)-(9.9) appliqué directement au problème continu est convergent.

### REMARQUE 9.2.

Pour $p = 2$, la méthode ci-dessus, appliquée à la résolution du problème (2.1), a peu d'intérêt ; ceci étant pour $\rho_n = \frac{1}{\varepsilon}$ il y a convergence de $u_h^n$ en deux itérations, que (9.7)-(9.9) soit appliqué au problème approché ou au problème continu ; cela va résulter de la

### PROPOSITION 9.4.

Pour $p = 2$ et $\rho_n = \frac{1}{\varepsilon}$ il y a convergence en deux itérations de la suite $u_h^n$ définie par l'algorithme (9.7)-(9.9).

### DEMONSTRATION

La démonstration ci-après vaut encore si (9.7)-(9.9) est appliquée au problème continu.

Puisque $p = 2$ on a $B(z) = z$ et (9.9), (9.13), (9.14) impliquent pour $n \ge 0$
$$
(9.25) \qquad \frac{1}{\epsilon} (\nabla \bar{u}_h^n - \bar{y}_h^n, \nabla v_h) = -(\bar{\lambda}_h^n, \nabla v_h) \qquad \forall v_h \in V_h
$$

$$
(9.26) \qquad \bar{y}_h^n + \frac{1}{\epsilon} (\bar{y}_h^n - \nabla \bar{u}_h^n) = \bar{\lambda}_h^n
$$

$$
(9.27) \qquad \bar{\lambda}_h^{n+1} = \bar{\lambda}_h^n - \frac{1}{\epsilon} (\bar{y}_h^n - \nabla \bar{u}_h^n).
$$

De (9.25), (9.27) on déduit

$$
(9.28) \qquad (\bar{\lambda}_h^{n+1}, \nabla v_h) = 0 \qquad \forall v_h \in V_h, n \ge 0.
$$

De (9.25), (9.28) on déduit

$$
(9.29) \qquad (\nabla \bar{u}_h^n - \bar{y}_h^n, \nabla v_h) = 0 \qquad \forall v_h \in V_h, \forall n \ge 1,
$$

et (9.26), (9.28), (9.29) impliquent

$$
(9.30) \qquad (\bar{y}_h^n, \nabla v_h) = 0 \qquad \forall v_h \in V_h, \forall n \ge 1;
$$

compte tenu de (9.29), (9.30), on a alors

$$
(9.31) \qquad \int_{\Omega} \nabla \bar{u}_h^n \cdot \nabla v_h \, dx = 0 \qquad \forall v_h \in V_h, \forall n \ge 1,
$$

d'où en prenant $v_h = \bar{u}_h^n$ dans (9.31), $\bar{u}_h^n = 0$ pour $n \ge 1$.

REMARQUE 9.3.

A µ *fixé*, L_ est strictement convexe en (v,z) et quadra-
tique en v; il en résulte que le problème (9.8) peut être résolu par une
variante de type *relaxation* en z_h, *surrelaxation* en v_h, de la méthode de
surrelaxation usuelle du problème de Dirichlet linéaire; les résultats
de CEA -GLOWINSKI [18] s'appliquent à cette variante, très facile à mettre
en oeuvre puisqu'à v_h, µ_h *donnés* la minimisation en z_h de L_ se décompose
en Card(C_h) problèmes à deux variables aisément résolubles; ceci est une
des justifications de l'algorithme (9.7)-(9.9) et on reviendra sur cet
aspect au N°12 relatif aux expériences numériques.
## 10 - REMARQUE SUR LA DUALITE

On va expliciter un problème dual du problème de minimisation (2.5) équivalent à (2.1). On utilise l'approche par le min-max (Cf. EKELAND TEMAM [19]); compte tenu de l'équivalence entre (2.5) et

$$ (10.1) \qquad \left\{ \begin{array}{l} \displaystyle \text{Min} \left[ \frac{1}{p} \int_{\Omega} |z|^p dx - \langle f, v \rangle \right], \\[1em] (z, v) \in L^p(\Omega) \times L^p(\Omega) \times V, \nabla v - z = 0 \end{array} \right. $$

on introduit le lagrangien $\mathcal{L}$ défini par

$$ (10.2) \qquad \mathcal{L}(v,z;\mu) = \frac{1}{p} \int_{\Omega} |z|^p dx + \int_{\Omega} \mu \cdot (\nabla v - z) dx - \langle f, v \rangle . $$

Le problème dual est alors défini par

$$ (10.3) \qquad \left\{ \begin{array}{ll} \displaystyle \text{Max} & \text{Min } \mathcal{L}(v,z;\mu) \\ & \mu \in (v,z) \\ & \mu \in L^{p'}(\Omega) \times L^{p'}(\Omega), (v,z) \in V \times L^p(\Omega) \times L^p(\Omega). \end{array} \right. $$

La forme explicite de (10.3) est donnée par

$$ (10.4) \qquad \left\{ \begin{array}{l} \displaystyle \text{Min}_{\mu \in H} \frac{1}{p'} \int_{\Omega} |\mu|^{p'} dx \\[1em] H = \{\mu | \mu \in L^{p'}(\Omega) \times L^{p'}(\Omega), \nabla \cdot \mu + f = 0\}. \end{array} \right. $$

Soit u la solution de (2.1), on démontrerait aisément la

### PROPOSITION 10.1

Le problème dual (10.4) admet une solution et une seule, soit $\lambda$, et on a

$$ (10.5) \qquad \lambda = |\nabla u|^{p-2} \nabla u. $$

■

Pour résoudre (10.4) on peut utiliser une méthode d'approximation, de type éléments finis, par des fonctions vérifiant (de façon exacte ou approchées) la condition $\nabla \cdot \mu + f = 0$; en particulier il conviendrait de voir dans quelle mesure les résultats de CROUZEIX-RAVIART [20] relatifs à $\nabla \cdot \mu = 0$ peuvent
être adaptés à la présente situation. Il convient de remarquer qu'avec de telles approximations on approchera l'opérateur non linéaire $|\mu|^{p'}-2\mu$, gradient de la fonctionnelle $\frac{1}{p'}\int_{\Omega}|\mu|^{p'}dx$, par un opérateur en dimension fini, *non diagonal* par rapport aux degrés de liberté définissant l'approximation. ■

On peut également traiter la contrainte $\nabla \cdot \mu + f = 0$ par dualisation et/ou pénalisation; on fera l'étude sur le cas continu (certaines intégrations, totalement formelles dans le cas continu, seraient justifiées pour le problème approché).

On associe à (10.4) le lagrangien $\mathcal{L}^* : L^p'(\Omega) \times L^p'(\Omega) \times V \to R$ défini par

$$ (10.6) \qquad \mathcal{L}^*(\mu, v) = \frac{l}{p'} \int_{\Omega} |\mu|^{p'} dx + \langle v, \nabla \cdot \mu + f \rangle; $$

on démontre facilement que $\mathcal{L}^*$ admet comme point-selle unique, sur $L^p'(\Omega) \times L^p'(\Omega) \times V$, le couple $(|\nabla u|^{p'-2}\nabla u, u)$, u solution de (2.1).

Un algorithme de dualité associé à $\mathcal{L}^*$ est alors

$$ (10.7) \qquad u^\circ \in V, \text{ donné} $$

$u^n$ donné on détermine $\lambda^n$ et $u^{n+1}$ par

$$ (10.8) \qquad \mathcal{L}^*(\lambda^n, u^n) \le \mathcal{L}^*(\mu, u^n) \quad \forall \mu \in L^p'(\Omega) \times L^p'(\Omega), \lambda^n \in L^p'(\Omega) \times L^p'(\Omega) $$

$$ (10.9) \qquad u^{n+1} = u^n + \rho_n (-\Delta)^{-1} (\nabla \cdot \lambda^n + f) , \quad \rho_n > 0. $$

Compte tenu de (10.8) on a $\lambda^n = |\nabla u^n|^{p'-2}\nabla u^n$, d'où $\nabla \cdot \lambda^n = -A(u^n)$ et

$$ (10.10) \qquad u^{n+1} = u^n - \rho_n (-\Delta)^{-1} (A(u^n) - f). $$

L'algorithme (10.7)-(10.10) n'est rien d'autres que la méthode du gradient avec *opérateur auxilliaire*, appliquée à la résolution de (2.5). Cet algorithme étudié et utilisé, entre autres, dans [2], n'est réellement opérationnel, appliqué à (2.5), que pour p égal à 2 comme on l'a déjà signalé. ■

On peut espérer améliorer la situation en procédant comme au N°9 par pénalisation et dualisation simultanées.

Soit $\mathcal{L}^* : L^p'(\Omega) \times L^p'(\Omega) \times V \to \mathbb{R}$ défini par
$$ (10.11) \qquad \mathcal{L}_{\varepsilon}^{*}(\mu, \nu) = \mathcal{L}^{*}(\mu, \nu) + \frac{1}{2\varepsilon} \langle (-\Delta)^{-1} (\nabla \cdot \mu + f), \nabla \cdot \mu + f \rangle ; $$

on montrerait que $(|\nabla u|^{p-2}\nabla u, u)$ est point-selle de $\mathcal{L}_{\varepsilon}^*$ sur $L^p'(\Omega) \times L^p'(\Omega) \times V$, d'où l'algorithme

$$ (10.12) \qquad u^{\circ} \in V, \text{ donné} $$

$u^n$ donné on détermine $\lambda^n$ et $u^{n+1}$ par

$$ (10.13) \qquad \mathcal{L}_{\varepsilon}^{*}(\lambda^{n}, u^{n}) \le \mathcal{L}_{\varepsilon}^{*}(\mu, u^{n}) \quad \forall \mu \in L^{p'}(\Omega) \times L^{p'}(\Omega), \lambda^{n} \in L^{p'}(\Omega) \times L^{p'}(\Omega) $$

$$ (10.14) \qquad u^{n+1} = u^n + \rho_n (-\Delta)^{-1} (\nabla \cdot \lambda^n + f), \quad \rho_n > 0. $$

On démontrerait, en procédant comme au N°9.4, que si $0 < r_o \le \rho_n \le r_l < \frac{2}{\varepsilon}$ on a

$$ (10.15) \qquad \lim_{n \to +\infty} |\lambda^n|^{p'-2} \lambda^n = \nabla u $$

$$ (10.16) \qquad \lim_{n \to +\infty} u^n = u. $$

L'utilisation de (10.12)-(10.14) sera justifiée si le calcul de $\lambda^n$ à partir de $u^n$, i.e. la résolution de (10.13), est un problème simple. Posons $G = (-\Delta)^{-1}$, $\lambda^n$ est alors solution du problème *non linéaire*

$$ (10.17) \qquad |\lambda^n|^{p'-2} \lambda^n - \frac{1}{\varepsilon} \nabla G (\nabla \cdot \lambda^n) = \frac{1}{\varepsilon} \nabla G f + \nabla u^n. $$

L'opérateur linéaire $\mu \to -\nabla G(\nabla \cdot \mu)$ de $L^2(\Omega) \times L^2(\Omega) \to L^2(\Omega) \times L^2(\Omega)$ a une interprétation très simple; on a en effet

$$ (10.18) \qquad L^2(\Omega) \times L^2(\Omega) = L_1 \oplus L_2 $$

avec

$$ (10.19) \qquad L_1 = \nabla(H_o^1(\Omega)) $$

$$ (10.20) \qquad L_2 = L_1^\perp = \{\mu | \mu \in L^2(\Omega) \times L^2(\Omega), \nabla \cdot \mu = 0\}, $$

l'opérateur $\mu \to -\nabla G(\nabla \cdot \mu)$ n'est rien d'autre que le projecteur de $L^2(\Omega) \times L^2(\Omega) \to L_1$. Cet opérateur a l'inconvénient d'être *non local* comme le montre l'exemple ci-après; si $\Omega = ]0,1[$ on a alors
$$ (10.21) \qquad -\nabla G(\nabla \cdot \mu) = \mu - \int_0^1 \mu(x) dx \qquad \forall \mu \in L^2(0,1). $$

Ce caractère *non local* de l'opérateur du premier membre de (10.17), ajouté à la non linéarité, a pour effet de rendre compliquée la détermination de $\lambda^n$, donc la mise en oeuvre de l'algorithme (10.12)-(10.14).

### REMARQUE 10.1

Il serait intéressant d'étudier l'application des méthodes d'éléments finis hybrides, développées, entre autres, par J.M. THOMAS [21], aux différents problèmes duals considérés dans ce N°10.

## 11. APPLICATION A D'AUTRES PROBLEMES NON LINEAIRES.

On peut appliquer la méthode de pénalité-dualité du N°9, avec quelques modifications mineures, aux problèmes suivants (en continu ou sur le problème approché).

### TORSION ELASTO-PLASTIQUE D'UNE BARRE CYLINDRIQUE

$$ (11.1) \qquad \Min_{v \in K} \left[ \frac{1}{2} \int_{\Omega} |\nabla v|^2 dx - C \int_{\Omega} v dx \right] $$

avec

$$ (11.2) \qquad K = \{v | v \in H_0^1(\Omega), |\nabla v| \le 1 \text{ p.p.}\}. $$

### ECOULEMENT D'UN FLUIDE VISCO-PLASTIQUE DANS UNE CONDUITE CYLINDRIQUE.

$$ (11.3) \qquad \Min \left[ \frac{1}{2} \int_{\Omega} |\nabla v|^2 dx + g \int_{\Omega} |\nabla v| dx - \int_{\Omega} fv dx \right], \quad v \in H_0^1(\Omega). $$

### SURFACES MINIMA

$$ (11.4) \qquad \Min \int_{\Omega} \sqrt{1 + |\nabla v|^2} dx, \quad v = g \text{ sur } \Gamma. $$

D'une façon générale la méthode ci-dessus est très bien adaptée aux problèmes elliptiques d'ordre 2, non linéaires, lorsque la non-linéarité porte sur $\nabla v$ ; c'est le cas des problèmes précédents et du problème de l'alternateur en magnéto-statique traité dans [14] par surrelaxation non linéaire, de l'écoulement subsonique d'un fluide parfait *compressible* autour d'un profil de $\mathbb{R}^2$, etc... L'application de l'algorithme de pénalité-dualité à ces différents problèmes fera l'objet d'un travail ultérieur de FORTIN-GLOWINSKI-MARROCCO [22].
### REMARQUE 11.1

Les considérations de dualité du N°10 s'appliquent bien entendu aux problèmes ci-dessus.

### REMARQUE 11.2

D'un point de vue plus informatique il convient de souligner que tous les problèmes mentionnés ci-dessus peuvent être résolus par le même programme d'ordinateur (pour un même type d'approximation évidemment).

## 12. APPLICATIONS NUMÉRIQUES

### 12.1 Triangulation de $\Omega$

On a pris pour $\Omega$ le disque $\{x \mid x \in \mathbb{R}, x_1^2 + x_2^2 < R\}$ sur lequel on a utilisé les triangulations données par les figures 1,2,3.

<figure><img src="image_3.png" /><figcaption>Fig.1</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Fig.2</figcaption></figure>
<figure><img src="image_5.png" /></figure>

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

12.2 Exemple 1

On a pris $f = C^{te} = C$, dans ce cas la solution de (2-1)
est donnée, en fonction de $r = \sqrt{x_2^2 + x_1^2}$ par
$$ (12-1) \quad u(r) = \frac{P-1}{P} \left(\frac{CR}{2}\right)^{\frac{P}{P-1}} R \left(1-\left(\frac{r}{R}\right)^{\frac{P}{P-1}}\right) $$

Les calculs ont été effectués avec R= 1. On a indiqué
sur les figures 4,5,6 les variations de diverses erreurs relatives en fonction
de p et pour les trois triangulations du domaine Ω

<figure><img src="image_7.png" /></figure>
<figure><img src="image_8.png" /><figcaption>ERREUR RELATIVE MOYENNE AUX C.D.G.</figcaption></figure>

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

La figure 4 donne les variations de l'erreur relative maximale aux sommets de $\mathcal{C}_h$, la figure 5 les variations de l'erreur relative moyenne aux sommets de $\mathcal{C}_h'$, i.e.

$$ (12-2) \qquad \frac{1}{N_h} \sum_{i=1}^{N_h} \left| \frac{u_h(m_i) - u(m_i)}{u(m_i)} \right| $$
où $N_h$ est le nombre de sommets de $\mathcal{L}_h$ intérieurs à $\Omega$, $(m_i)$ est l'ensemble, supposé ordonné, des noeuds intérieurs, $u(m_i)$ est la valeur de la solution exacte calculée au noeud $m_i$, $u_h(m_i)$ est la valeur de la solution calculée, évaluée au même point. La figure 6 donne l'erreur relative moyenne, calculée cette fois-ci au centre de gravité des triangles. Les divers essais numériques effectués nous ont conduit à faire les remarques suivantes.

### REMARQUE 12.1

Les différentes erreurs relatives ci-dessus apparaissent comme étant indépendantes de C.

### REMARQUE 12.2

Les erreurs aux centres de gravité des triangles sont beaucoup plus importantes que les erreurs aux sommets de ces mêmes triangles. Un éventuel phénomène de *super-convergence* aux sommets de $\mathcal{L}_h$ serait donc à étudier. A titre indicatif nous avons donné sur les figures 7,8, 9, 10 des coupes suivant un diamètre de la solution exacte pour les valeurs suivantes de $p$: 1,1 ; 1,75 ; 10 ; 50 ; ... Les valeurs de C retenues sont indiquées sur les figures.

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

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

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

Pour la résolution numérique du problème (3-4), nous avons tout d'abord
utilisé l'algorithme EGSN (décrit et utilisé dans GLOWINSKI - MARROCCO
{14}), cela à cause des formulations très voisines du problème, et aussi
parce que EGSN avait donné d'excellents résultats pour le problème de
magnétostatique.

Les résultats numériques présentés montrent que cet algorithme reste très efficace lorsque p reste "voisin de 2", mais devient inopérant en dehors de cet intervalle.

12.2.1 UTILISATION DE L'ALGORITHMÉ EGSN.

Nous allons nous ramener à une formulation semblable à celle
utilisée dans (GLOWINSKI - MARROCCO {14}). La fonctionnelle d'énergie (2-3)
peut aussi s'écrire

$$
(12-3) \quad \mathfrak{F}(u) = \frac{1}{P} \int_{\Omega} \left( \left( / \text{gradu} / \right)^2 \right)^{\frac{P}{2}} dx \quad - \langle f, u \rangle \qquad (*)
$$

La fonction $\psi(\eta)$ est ici $\eta = \frac{P}{2}$ et la fonction $v$
(inverse de la perméabilité magnétique dans le problème de magnétostatique)

$$
est tout simplement \quad v(\eta) = \eta \quad (où \quad \eta = \sqrt{\frac{P}{\lambda \left(\frac{1}{2}\right)})}
$$

<footnote>(*) L'équation en magnétostatique s'écrivait
$$
\begin{cases}
- \sum_{i=1}^{2} \frac{\partial}{\partial x_i} \left( v \left( |\text{gradu}|^2 \right) \frac{\partial u}{\partial x_i} \right) = j \\
u|_{\eta=0} = 0
\end{cases}
$$</footnote>
Nous donnons sur la figure 11 quelques courbes $v(\eta)$ pour différentes valeurs de $p$. Quelle que soit la valeur de $p$, la courbe $v(\eta)$ reste tout de même différente de la courbe utilisée dans { 14 } .

Au vu des courbes représentées sur la figure 11, il faut s'attendre à certaines difficultés de résolution. Le caractère $v$ strictement croissante, positive et bornée n'est jamais constaté ici.

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

Rappelons maintenant l'algorithme EGSN.

Si $u_h$ est la solution approchée calculée, elle peut être représentée par les valeurs qu'elle prend aux noeuds intérieurs de $\mathcal{C}_h$ (supposés numérotés de 1 à $N_h$)

$$u_h = (u_1, u_2, \dots, u_{N_h})$$

La fonctionnelle (12.3) s'écrit donc explicitement comme fonction des variables scalaires $u_i$ ($i=1$ à $N_h$), $\mathfrak{F}(u_1, u_2, \dots, u_{N_h})$, si bien que EGSN s'exprime par

$$ (12.4) \quad \left\{ \begin{array}{ll} u_h^\circ = (u_1^\circ, u_2^\circ, \dots, u_{N_h}^\circ) & \text{donné dans } v_h \\ u_i^{n+1} = u_i^n - \omega \frac{\frac{\partial \mathfrak{F}}{\partial v_i}(u_1^{n+1}, \dots, u_{i-1}^{n+1}, u_i^n, u_{i+1}^n, \dots, u_{N_h}^n)}{\frac{\partial^2 \mathfrak{F}}{\partial v_i^2}(u_1^{n+1}, \dots, u_{i-1}^{n+1}, u_i^n, u_{i+1}^n, \dots, u_{N_h}^n)} & \text{et } 0 < \omega < \omega_m \le 2 \\ & i=1,2, \dots, N_h \end{array} \right. $$
On ne peut avec EGSN initialiser l'algorithme avec
$u_h = 0$ a cause de la valeur de $v'(0)$ et $v''(0)$.

La fonction $u_h$ choisie pour initialiser EGSN a été
l'interpolée linéaire de $U(r) = (R - r)(R + r)$ (avec $r = \sqrt{x_1^2 + x_2^2}$ et $R = 1$)

### 12.2.2 ANALYSE DES RESULTATS OBTENUS AVEC EGSN.

Nous avons représenté sur les figures 12,13, 14, 15, le
nombre d'itérations (valeur de n dans 12-4) nécessaires à l'obtention
du test

$$ (12.5) \quad \text{Res}(n) = \frac{\sum_{i=1}^{n} \left| u_i^{n+1} - u_i^n \right|}{\sum_{i=1}^{n} \left| u_i^{n+1} \right|} \le 10^{-5} $$

en fonction du paramètre de relaxation $\omega$, cela pour $p = 1.75$, $p = 3$, $p = 4$ et $p = 7$ ; les valeurs extrêmes de $p$ délimitent l'intervalle (sur l'axe des $p$) pour lequel l'algorithme EGSN reste efficace.

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

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

<figure><img src="image_17.png" /></figure>
Sur chaque figure (12 à 15) on trouvera
plusieurs courbes ; la courbe 1 correspond
à C = 0,1 ; la courbe 2 à C = 1, la
courbe 3 à C = 4, la courbe 4 à C = 10
et la courbe 5 (si elle est présente) à
C = 20.

Le choix du ω optimal (minimisant le nombre
d'itérations nécessaires à 1 obtention de
12.5) dépend apparemment plus du second
membre lorsque p < 2 ; on pourra aussi con-
sulter la figure 16 qui donne la valeur du
paramètre de relaxation optimal en fonction
du second membre C. Chaque courbe est numé-
rotée et la valeur correspondante de p est
indiquée sur la figure. Sur la figure 17, on trouvera les courbes donnant
le nombre minimum d'itérations nécessaires (i.e. avec ω opt) en fonction du
second membre C. Comme pour la figure 16 les valeurs de p sont indiquées
sur la figure.

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

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

<figure><img src="image_20.png" /></figure>
Sur les figures 16 et 17 on trouvera aussi des courbes (ou plutôt des morceaux de courbes) relatives à $p = 1,5$.

La détermination du $\omega$ opt est très difficile. Le nombre d'itérations nécessaires à l'obtention de (12.5) est très élevé et de plus les courbes analogues aux courbes des figures 12,13,14,15 sont "très pointues", ce qui exprime que le nombre d'itérations est très sensible à la variation de $\omega$.

De plus pour $p = 1.5$ des ennuis de résolution numérique sont survenus, ennuis liés au fait que $|\text{gradu}|^2$ calculé peut être nul sur un triangle et alors $\nabla$ et $\nabla'$ prennent des valeurs infinies, il devient donc nécessaire de faire des tests supplémentaires alourdissants l'algorithme.

Nous avons aussi tenté d'appliquer EGSN lorsque $p = 1.25$, là l'échec a été total même en faisant de la sous-relaxation.

Une autre constatation que l'on peut faire au vu des courbes de la figure 17, est que pour $p > 4$, le nombre d'itérations nécessaires à l'obtention de (12.5) augmente avec $p$ et avec le second membre $C$.

Sur la figure 18 nous donnons la variation de $\text{Res}(n)$ donné en (12.5) au cours des itérations, pour différentes valeurs de $p$ (allant de 1,5 à 10), le second membre $C$ étant toujours égal à 1. Sur les courbes 4,5,6 de la figure 18, correspondant respectivement aux valeurs 5.5, 7 et 10 pour $p > 4$, on constate que la quantité $\text{Res}(n)$ reste relativement constante pendant quelques itérations, le "palier" allant grandissant (et très rapidement) avec la valeur de $p$.

<figure><img src="image_21.png" /><figcaption>Fig. 18</figcaption></figure>
Ainsi par exemple pour $p=10$ (courbe 6 de la figure 18) Res (n) garde une valeur apparemment constante de l'itération 30 à l'itération 380, la valeur $10^{-5}$ étant ensuite atteinte rapidement (20 itérations environ).

Les résultats présentés ci-dessus ont été fait avec la triangulation 2 (figure 2) comportant 256 triangles et 113 noeuds intérieurs. Le temps calcul pour effectuer une itération est, en simple précision, sur ORDINATEUR CII 10 070, de 0.8 s

### 12.2.3 REMARQUES SUR LA MISE EN OEUVRE NUMERIQUE DE L'ALGORITHME (9.7) - (9.9)

Soit $\mu_h \in L_h$ ; le couple $(\hat{u}_h, \hat{y}_h)$ solution de

$$ (12.6) \quad \left\{ \begin{array}{l} \mathcal{L}_{\varepsilon} ( (\hat{u}_h, \hat{y}_h ; \mu_h) ) \le \mathcal{L}_{\varepsilon} ( v_h, z_h ; \mu_h ) \quad \forall v_h \in V_h, z_h \in L_h \\ \hat{u}_h \in V_h, \hat{y}_h \in L_h \end{array} \right. $$

est caractérisé par

$$ (12.7) \quad \left\{ \begin{array}{l} \displaystyle \int_{\Omega} \nabla \hat{u}_h \cdot \nabla v_h \, dx - \int_{\Omega} \hat{y}_h \cdot \nabla v_h \, dx = \varepsilon \langle f, v_h \rangle - \varepsilon \int_{\Omega} \mu_h \cdot \nabla v_h \, dx \quad \forall v_h \in V_h \\ \\ \hat{u}_h \in V_h \\ \\ \displaystyle \int_{\Omega} |\hat{y}_h|^{p-2} \hat{y}_h \cdot z_h \, dx + \frac{1}{\varepsilon} \int_{\Omega} (\hat{y}_h - \nabla \hat{u}_h) \cdot z_h \, dx = \int_{\Omega} \mu_h \cdot z_h \, dx \quad \forall z_h \in L_h \\ \\ \hat{y}_h \in L_h ; \end{array} \right. $$

(12.7) n'est autre qu'un système non linéaire en $\hat{u}_h, \hat{y}_h$ écrit sous forme variationnelle.
A $\mu_h$ et $\hat{y}_h$ donnés la première équation variationnelle (12.7) est un problème de Dirichlet, linéaire et homogène, approché. A $\mu_h$ et $\hat{u}_h$ donnés la seconde équation (12.7) équivaut à un système de 2 Card ($\mathcal{L}_h$) équations non linéaires, mais $\hat{y}_h, z_h, \nabla \hat{u}_h, \mu_h$ étant constants sur chaque triangle, ce système se décompose en Card ($\mathcal{L}_h$) systèmes, indépendants les uns des autres, de 2 équations non linéaires à 2 inconnues, admettant chacun une solution unique dans $\mathbb{R}^2$ ; de façon plus précise, on détermine $\hat{y}_h$ en résolvant pour tout $T \in \mathcal{L}_h$ le système

$$ (12.8) \qquad \left\{ \begin{array}{ll} \left|q_T\right|^2_{p-2} q_{1T} + \frac{1}{\varepsilon} q_{1T} = \frac{1}{\varepsilon} \left(\frac{\partial \hat{u}_h}{\partial x_1}\right)_T & + \mu_{1T} \\ \\ \left|q_T\right|^2_{p-2} q_{2T} + \frac{1}{\varepsilon} q_{2T} = \frac{1}{\varepsilon} \left(\frac{\partial \hat{u}_h}{\partial x_2}\right)_T & + \mu_{2T} \end{array} \right. $$

avec $q_T = (q_{1T}, q_{2T})$, $\mu_T = (\mu_{1T}, \mu_{2T}) - (z_h)_T$ désignant de façon générale la restriction de la fonction $z_h$ à $T$-ou

En combinant les 2 équations (12.8) on en déduit que $|q_T|$ est solution de l'équation

$$ (12.9) \qquad \left|q_T\right|^2_{p-1} + \frac{1}{\varepsilon} \left|q_T\right| - \left|\frac{1}{\varepsilon} (\nabla \hat{u}_h)_T + \mu_T\right| = 0, $$

et $|q_T|$ calculé on en déduit facilement $q_{1T}$ et $q_{2T}$ à partir de (12.8)

La fonction $g: \mathbb{R} \to \mathbb{R}$ définie par le premier membre de (12.9), fonction en la variable $|q_T|$ est convexe et croissante si $\ge p 2$, concave et croissante si $1<p<2$. La méthode utilisée pour résoudre (12.9) est la méthode de Newton à une variable, soit

$$ (12.10) \qquad \xi^\circ \mathbb{R}_+ \text{ donné} $$

$$ (12.11) \qquad \xi^{m+1} = \xi^m - \frac{g(\xi^m)}{g'(\xi^m)} $$
ξ° étant choisi de telle façon que $g(ξ°)$ $g''(ξ°) \ge 0$ ; ce type d'initialisation a une grande importance, surtout pour les "grandes valeurs" de $p$ et au voisinage de $p=1$

En ce qui concerne la résolution proprement dite du système (12.7), on peut
utiliser des méthodes de type *relaxation* (cf [18]) ; étant donné la structure
de (12.7), on peut envisager les deux méthodes ci-dessous

■ Première méthode

Partant du couple initial $(\hat{u}_h^o, \hat{y}_h^o)$ on passe de $(\hat{u}_h^m, \hat{y}_h^m)$ à
$\hat{u}_h^{m+1}$ puis $\hat{y}_h^{m+1}$ en deux étapes.

1ère Etape - Tout élément $v_h$ de $V_R$ étant complètement défini par les valeurs prises sur les noeuds de $\mathcal{C}_h$, on calcule $\hat{u}_h^{m+1}$ par un seul balayage de type surrelaxation ponctuelle linéaire des noeuds de $\mathcal{C}_h$.

2ème Etape - La fonction $\hat{u}_h^{m+1}$ étant déterminée on calcule ensuite $\hat{y}_h^{m+1}$ par résolution des Card $(\mathcal{C}_h)$ systèmes non linéaires du type (12.8).

Par une généralisation facile des méthodes de [18] on démontrerait la conver-
gence de l'algorithme ci-dessus si le paramètre $\omega$ de surrelaxation vérifie
$1 < \omega < 2$ ; dans les applications qui suivent on a pris pour $\omega$ la valeur optimale,
dans la résolution par surrelaxation ponctuelle linéaire, du problème de
Dirichlet approché correspondant à la première équation (12.7).

Test d'arrêt des itérations.

Dans la résolution de (12.9) par la méthode de Newton (12.10), (12.11),
on utilise le test suivant (si $\xi^{m+1} \neq 0$)

$$
(12.12) \qquad \frac{|\xi^{m+1} - \xi^m|}{|\xi^{m+1}|} \le n_1;
$$

Dans la résolution du système (12.7) par la méthode décrite ci-dessus,
on a retenu un test d'arrêt ne portant que sur la suite $(\hat{u}_h^m)_m$, soit
$$
(12.13) \qquad \frac{\sum_i |\hat{u}_h^{m+1}(m_i) - \hat{u}_h^m(m_i)|}{\sum_i |\hat{u}_h^{m+1}(m_i)|} \le \eta_2
$$

ou $(m_i)_{i \in I}$ a déjà été défini en (12.2)

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

Deuxième méthode

On utilise la méthode par *blocs* suivante

$$
(12.14) \qquad (\hat{u}_h^c, \hat{y}_h^o) \text{ donné,} \\
\qquad \qquad (\hat{u}_h^m, \hat{y}_h^m) \text{ connu, on calcule } \hat{u}_h^{m+1} \text{ en résolvant le système linéaire}
$$

$$
(12.15) \quad \left\{
\begin{array}{@{}l@{}}
  \displaystyle \int_{\Omega} \nabla \hat{u}_h^{m+1} \cdot \nabla v_h \, dx - \int_{\Omega} \hat{y}_h^m \cdot \nabla v_h \, dx = \varepsilon \langle f, v_h \rangle - \varepsilon \int_{\Omega} \mu_h \cdot \nabla v_h \, dx \\
  \qquad \wedge v_h \in V_h ; \hat{u}_h^{m+1} \in V_h ;
\end{array}
\right.
$$

puis $\hat{y}_h^{m+1}$ en résolvant, comme indiqué précédemment, le système non linéaire

$$
(12.16) \quad \left\{
\begin{array}{@{}l@{}}
  \displaystyle \int_{\Omega} |\hat{y}_h^{m+1}|^{p-2} \hat{y}_h^{m+1} \cdot z_h \, dx + \frac{1}{\varepsilon} \int_{\Omega} (\hat{y}_h^{m+1} - \nabla \hat{u}_h^{m+1}) \cdot z_h \, dx = \int_{\Omega} \mu_h z_h \, dx \\
  \qquad \forall z_h \in L_h ; \hat{y}_h^{m+1} \in L_h
\end{array}
\right.
$$

La convergence de (12.14) - (12,16) peut se démontrer à partir des résultats
de [18] .

Le problème (9.8) étant un cas particulier du problème (12.7) (avec $\mu_h = \lambda_h^n$);
on peut utiliser, pour déterminer $(u_h^n, y_h^n)$, les deux méthodes décrites

ci-dessus (initialisées par $(u_h^{n-1}, y_h^{n-1})$). D'un point de vue pratique,
on a utilisé la première méthode, telle qu'elle a été décrite, et la variante de la deuxième méthode obtenue en se limitant à *une seule itération* par bloc du type (12.15) - (12.16) (la convergence globale de la variante de l'algorithme (9.7) - (9.9) ainsi obtenue est démontrée dans B.Mercier [23] et D.Gabay - B.Mercier [24]).

Comme test d'arrêt de l'algorithme (9.7) - (9.9) on a pris

$$ (12.17) \quad \text{Res}(n) = \frac{\sum_i \left| u_h^{n+1}(m_i) - u_h^n(m_i) \right|}{\sum_i \left| u_h^{n+1}(m_i) \right|} \le \eta_3 $$

### 12.2.4 ANALYSE DES RESULTATS OBTENUS

Les résultats numériques présentés dans cette section sont relatifs à la triangulation de la figure 2 (256 triangles et 113 noeuds intérieurs). On a représenté sur les figures 19 à 40, la dépendance en $\varepsilon$ du comportement des algorithmes décrits au n°12.2.3 ; les expériences numériques correspondantes ont été effectuées en *simple précision* (7 chiffres significatifs) sur ordinateur CII 10070.

Détaillons les figures ci-dessous :

Les courbes 1,2, 3, 4 correspondent respectivement à $C = 0.1, 1, 4, 10$. Les figures 19 à 28 sont relatives au *premier algorithme* décrit au n° 12.2.3 et indiquent le nombre total de *balayages* (i.e le nombre total d'itérations de type relaxation utilisées dans la résolution du problème (9.8) nécessaires à l'obtention de la convergence, au sens du test (12.17), lorsque $\rho_n = \rho = \frac{1}{\varepsilon}$ dans (9.9), $\eta_1 = 10^{-4}$ dans (12.12), $\eta_2 = 10^{-5}$ dans (12.13), $\eta_3 = 10^{-5}$ dans (12.17). Les figures 29 à 40 donnent des résultats de même nature, pour la variante de la *deuxième méthode* obtenue en se limitant à *une seule itération* par bloc du type (12.15) - (12.16) ; les itérations considérées sont celles indexées par $n$ et $\rho_n$, $\eta_1, \eta_2$ étant choisis comme ci-dessus on a pris $\eta_3 = 10^{-5}$ (courbes (a)) et $\eta_3 = 5 \cdot 10^{-5}$ (courbes (b)) dans (12.17).
NOMÎBRE DE IBALAYAGES

P=1.1

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

NOMÎBRE DE IBALAYAGES

P=1.1

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

NOMÎBRE DE IBALAYAGES

P=1.25

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

NOMÎBRE DE IBALAYAGES

P=1.25

<figure><img src="image_26.png" /></figure>
<figure><img src="image_27.png" /></figure>

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

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

<figure><img src="image_30.png" /></figure>
<figure><img src="image_31.png" /></figure>

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

### REMARQUE 12.3 :

La valeur du paramètre $ε$ indiqué sur les figures 19 à 45 est en réalité 2 fois la valeur du paramètre de même nom utilisé dans les n° 9 et N° 12.2.3.
<figure><img src="image_33.png" /></figure>

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

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

<figure><img src="image_36.png" /></figure>
<figure><img src="image_37.png" /></figure>

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

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

<figure><img src="image_40.png" /></figure>
<figure><img src="image_41.png" /></figure>

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

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

<figure><img src="image_44.png" /></figure>
## REMARQUE 12.4

Au vu des figures 19 à 40 il apparaît que pour $\epsilon$ "grand", le nombre d'itération augmente considérablement (ou même on ne peut plus réaliser le test 12.17 avec $\eta = 10^{-5}$), ce phénomène étant d'autant plus sensible que p est voisin de 1 ou grand (disons $p \ge 7$). Ceci est lié aux *erreurs d'arrondi*, en effet par un calcul en *double précision* (11 chiffres significatifs), on obtient – toutes choses égales par ailleurs – les résultats présentés sur les figures 41 et 42 ; les courbes montrent que par rapport à la simple précision, le choix de $\epsilon$ n'est pas critique et la remontée brusque du nombre d'itérations (ou la non convergence) pour $\epsilon$ grand se produit pour des valeurs considérablement plus élevées. Les courbes tracées sur ces figures correspondent respectivement à $p = 1.1$ et $p = 1.5$ pour la valeur $C = 0.1$, nous avons indiqué en pointillés les courbes correspondant aux calculs simple précision (fig. 19 et 23).

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

<figure><img src="image_46.png" /></figure>
### REMARQUE 12.5

Les figures 29 à 40 montrent l'importance du choix de $η_3$ pour l'obtention de la convergence au sens du test (12.17) ; les phénomènes de non convergence qui peuvent se produire lorsque $η_3$ est trop petit sont - là encore - conséquence des erreurs d'arrondis.

### REMARQUE 12.6

Pour l'algorithme du N° 12.2.3 (1ère méthode), le nombre d'itérations en $λ$ (indice n de l'algorithme (9.8) - (9.10)) ne semble ni dépendre de p, ni du second membre C cela pour un intervalle assez large de choix du $ε$ de pénalisation (la valeur commune constatée est 4). Une autre constatation qui vaut aussi bien pour cet algorithme que pour la variante de la méthode par bloc (12.15)-(12.16) limitée à 1 seul cycle est : le nombre de balayages pour la première méthode où le nombre d'itérations pour la variante dépendent de $ε$, (bien entendu), mais apparemment ne dépendent pas du second membre C. Plus précisement des essais numériques effectués nous avons extrait la propriété suivante :

Si N est le nombre d'itérations (ou de balayages) pour $ε_0, C = Co$ et p donnés, on obtient le même nombre N pour $C = C_1$ si on choisit $ε_1$ de la façon suivante :

$$ (12.18) \qquad \epsilon_1 = \epsilon_0 \left( \frac{C_1}{C_0} \right)^{\frac{2-p}{p-1}} $$

Cette formule est de vérification facile pour $1 < p < 1.5$, au vu des résultats numériques ; pour $p = 2$ on trouve que le $ε$ ne dépend pas du second membre, et pour $p > 2$, on constate bien sur les figures 25 à 28 et 36 à 40 que la variation de $ε$ par rapport à C se trouve inversée.

### REMARQUE 12.7

A titre indicatif, nous donnons sur les figures 43,44, 45, les variations de la quantité (12.17) pour la variante de la deuxième méthode exposée dans le N°12.2.3. Les valeurs de p retenues sont respectivement 1.1, 4 et 50, le second membre $C = 0.1$ et les différentes valeurs de $ε$ représentées sont indiquées.
<figure><img src="image_47.png" /></figure>

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

<figure><img src="image_49.png" /></figure>
Pour différentes valeurs de p on a indiqué dans le tableau 12.1, les temps de calcul en secondes, en simple précision (sur CII 10070) obtenus avec $\rho_n = \frac{1}{\varepsilon}$ $\varepsilon$ ayant sa valeur optimale. Les colonnes T1 et T2 sont respectivement relatives à la première méthode et à la variante (à une itération par blocs) de la deuxième méthode décrites ci-dessus.

### 12.3 Exemple 2

Si $p > 2$, on a $-\Omega$ étant bidimensionnel $- V = W_o^{1,p}(\Omega) \subset C^0(\bar{\Omega})$ avec injection continue ; il en résulte que $\forall x_o \in \Omega$, la mesure de Dirac en $x_o$, $\delta(x_o)$ appartient à $V'$.

On a représenté sur les figures 46,47, 48 les équipotentialles de la solution du problème approché (3.4) ; calculée sur la triangulation de la figure 3 pour $p = 3$ avec

$$f = \delta_{(0,0)} - \delta_{(0,-0.875)} \quad \text{pour figure 46}$$

$$f = \delta_{(0,0)} - \delta_{(0, 0.25)} \quad \text{pour figure 47}$$

$$f = \delta_{(0,0)} - \frac{1}{2}\delta_{(0,5,0)} - \frac{1}{2}\delta_{(0,0.5)} - \frac{1}{2}\delta_{(-0.5,0)} - \frac{1}{2}\delta_{(0,-0.5)}$$

pour la figure 48

<table><caption>Tableau 12.1</caption><thead><tr><th>P</th><th>T1</th><th>T2</th></tr></thead><tbody><tr><td>1,1</td><td>100</td><td>110</td></tr><tr><td>1,25</td><td>75</td><td>55</td></tr><tr><td>1.5</td><td>70</td><td>45</td></tr><tr><td>1.75</td><td>65</td><td>27</td></tr><tr><td>2.</td><td>45</td><td>7</td></tr><tr><td>3.</td><td>150</td><td>30</td></tr><tr><td>4.</td><td>150</td><td>40</td></tr><tr><td>10</td><td>150</td><td>80</td></tr><tr><td>50</td><td>150</td><td>210</td></tr></tbody></table>
<figure><img src="image_50.png" /><figcaption>Fig. 46</figcaption></figure>

<figure><img src="image_51.png" /><figcaption>Fig. 47</figcaption></figure>

## REMARQUE 12.8

Si $\Omega$ est le disque de rayon $R$ et si $f = C\delta^{(o,o)}$, la solution de (2-1) est alors donnée par

$$ (12.19) \quad u(r) = \frac{P-1}{P-2} \left(\frac{C}{2\pi}\right)^{\frac{1}{P-1}} \left( R^{\frac{P-2}{P-1}} - r^{\frac{P-2}{P-1}} \right) \quad \text{si } p > 1, p \neq 2 $$

$$ (12.20) \quad u(r) = \frac{C}{2\pi} \log \frac{r}{R} \quad \text{si } p = 2 $$

Il résulte de (12.19), (12.20) que $u \notin W_o^{1,p}(\Omega)$ si $p \le 2$ ; ceci étant, l'utilisation des méthodes d'approximation et de résolution précédentes, conduit pour $p \le 2$ à des résultats très convenables - sauf au voisinage de l'origine - puisque par exemple, pour $p = 2$, sur la triangulation de la figure 3, l'erreur relative moyenne aux noeuds est $3.7 \cdot 10^{-3}$
<figure><img src="image_52.png" /><figcaption>Fig. 48</figcaption></figure>

Remerciements.

Nous voudrions remercier M.L.Tartar pour les discussions très utiles que nous avons eues à propos de ce travail.
BIBLIOGRAPHIE

[1] CIARLET P.G., SCHULTZ M.H., VARGA R.S.,
Numerical Methods of High-Order Accuracy for Nonlinear
Boundary Value Problems, V. Monotone Operator Theory,
Numer. Math. 13, 1969, 51-77.

[2] BREZIS H., SIBONY M.,
Methodes d'Approximation et d'Itération pour les Opérateurs
Monotones, Arch. of Rat. Mech. and Analysis 28, 1968,59-82.

[3] SIBONY H.,
Méthodes Itératives pour les Equations et Inéquations
aux Dérivées Partielles Non Linéaires de Type Monotone,
Calcolo 7, Fasc. 1-2, 1970, 65-183.

[4] LEBEAUD G.,
Sur l'Approximation de l'Equation Au = - ∑ (D_i.u |P-2|D_i.u) = T
Renconti di Matematica (3-4), vol. 2, Série VI, 1969, 443-471.

[5] GLOWINSKI R., MARROCCO A.,
Sur l'Approximation par Éléments Finis d'Ordre 1, et la
Résolution par Pénalisation-Dualité, d'une classe de problèmes
de Dirichlet non linéaires, C.R.A.S. Paris, T.278, 26 juin 1974,
Série A, 1649-1652.

[6] GLOWINSKI R., MARROCCO A.,
On the solution of a class of non linear Dirichlet problems
by a penalty-duality method and finite element of order one,
à paraitre dans les Proceedings du 6ème Congrès IFIP-TC7 sur
l'optimisation Lecture Notes in Comp. Sciences, Springer.

[7] LIONS J.L.,
Contrôle optimal des systèmes gouvernés par des équations
aux dérivées partielles, Dunod et Gauthier-Villars. Paris, 1968.

[8] PELISSIER M.C., Thèse, Paris Orsay, à paraitre.

[9] CIARLET P.G., RAVIART P.A.,
General Lagrange and Hermite interpolation in Rⁿ with applications
to finite element methods, Arch. Rat. Mech. Anal. 46, 1972, 177-199.
[10] TARTAR L., Interpolation non linéaire et Applications, Thèse Université Paris VI, 1971.

[11] STRANG G., FIX G., An Analysis of the Finite Element Method, Prentice-Hall, Englewood Cliffs, N.J. 1973.

[12] NECAS J., Les Methodes Directes en Théorie des Equations Elliptiques. Masson, Paris 1967.

[13] TARTAR L., Interpolation non Linéaire et Régularité, Journal of Functional Analysis, vol.9, N°4, 1972,469-489.

[14] GLOWINSKI R., MARROCCO A., Analyse Nûmérique du champ magnétique d'un alternateur par éléments finis et surrelaxation ponctuelle non linéaire, Comp. Meth. Appl. Mech. Engineering 3, 1974, 55-85.

[15] HESTENES M., Multiplier and Gradient Methods, J. of Optimization Theory and Applications, vol.4, N°5, 1969, 303-320.

[16] ARROW K.J, HURWICZ L., UZAWA H. Studies in Linear and non linear programming, Stanford University Press, 1958.

[17] ROCKAFELLAR, T., Convex Analysis, Princeton University Press, 1970.

[18] CEA J., GLOWINSKI R., Sur des méthodes d'optimisation par relaxation, Revue Française d'Automatique, Informatique et Rech. Opérationnelle Dec. 1973, R-3, 5-32.

[19] EKELAND I., TEMAM R., Analyse Convexe et Problèmes Variationnels. Dunod et Gauthier-Villars, Paris, 1974.

[20] CROUZEIX M., RAVIART P.A., Conforming and non conforming finite element methods for solving stationary Stokes equations, Revue Française d'Automatique, Informatique et Rech. Operationnelle, Dec. 1973, R-3, 33-76.

[21] THOMAS J.M., A paraitre.

[22] FORTIN M., GLOWINSKI R., MARROCCO A., A paraitre.

[23] MERCIER B, C.R.A.S., A paraitre.

[24] GABAY D., MERCIER B., A paraitre.