# Décomposition spectrale généralisée et séparation de variables : réduction de modèle pour la résolution de problèmes stochastiques

Anthony Nouy

► To cite this version:

Anthony Nouy. Décomposition spectrale généralisée et séparation de variables : réduction de modèle pour la résolution de problèmes stochastiques. 9e Colloque national en calcul des structures, CSMA, May 2009, Giens, France. hal-01413977

HAL Id: hal-01413977

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

Submitted on 12 Dec 2016

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 CC BY 4.0 - Attribution - International License
# Décomposition spectrale généralisée et séparation de variables : réduction de modèle pour la résolution de problèmes stochastiques

A. Nouy¹

¹ GeM - Institut de Recherche en Génie Civil et Mécanique, UMR CNRS 6183 (Université de Nantes, Centrale Nantes)

2 rue de la Houssinière, 44322 NANTES Cedex 3

anthony.nouy@univ-nantes.fr

**Résumé** — On propose une stratégie de résolution d'EDP stochastiques définies sur des espaces de probabilité de grande dimension. Cette stratégie repose sur la construction *a priori* d'une représentation en variables séparées de la solution. Une méthode de décomposition spectrale généralisée est utilisée pour découpler la résolution de problèmes déterministes et de simples équations stochastiques, ces dernières étant résolues par un algorithme de séparation de variables en grande dimension exploitant la structure produit tensoriel de l'espace des fonctions stochastiques.

**Mots clés** — Méthodes spectrales stochastiques, Décomposition spectrale généralisée, Séparation de variables, Réduction de modèle.

## 1 Introduction

La quantification et la propagation des incertitudes dans les modèles physiques apparaissent comme des voies essentielles vers l'amélioration de la prédiction de leur réponse. Dans ce contexte, les méthodes spectrales stochastiques de type Galerkin [1, 2, 3] offrent des outils robustes conduisant à des prédictions numériques précises. La réponse est recherchée dans un espace d'approximation qui est le produit tensoriel d'espaces d'approximation déterministe et stochastique. Le calcul d'une solution approchée demande alors la résolution de problèmes de très grande taille, conduisant à des coûts de calcul généralement prohibitifs. De plus, elles nécessitent une bonne maîtrise de la structure mathématique du modèle afin d'étendre les solveurs déterministes usuels au cadre stochastique (algorithmes non-linéaires, préconditionneurs dédiés, ...). Elles nécessitent donc souvent des développements logiciels spécifiques pour une classe particulière de problèmes.

Récemment, une méthode de réduction de modèle, baptisée décomposition spectrale généralisée (GSD), a été proposée afin de pallier les inconvénients mentionnés ci-dessus [4, 5, 6]. Cette méthode exploite la structure produit tensoriel de l'espace de travail et conduit à une représentation optimale de la solution sous la forme d'une décomposition en série de produits de fonctions déterministes et de variables aléatoires, cette série possédant les mêmes propriétés de convergence que des décompositions spectrales classiques (de type Karhunen-Loève). Les algorithmes proposés [5] demandent simplement la résolution de problèmes déterministes et de problèmes stochastiques sur bases réduites (équations stochastiques algébriques), ces problèmes étant découplés.
La méthode GSD ne répond cependant que partiellement à la “malédiction de la dimensionalité”, attribuée à l’augmentation dramatique de la dimension des espaces d’approximation stochastique, lorsqu’on travaille en grande dimension stochastique. Dans ce papier, nous proposons un mariage de la méthode GSD avec une technique de représentation en variables séparées [7, 8, 9], qui exploite la structure produit tensoriel de l’espace fonctionnel au niveau stochastique. La séparation de variables est utilisée pour la représentation des fonctions de carré intégrable à valeurs vectorielles (ou vecteur aléatoire du second ordre) définies sur un espace de probabilité de grande dimension. Dans la méthode GSD, ces fonctions sont solutions des problèmes stochastiques sur bases réduites. La stratégie proposée permet de calculer une solution approchée du modèle, vivant dans des espaces d’approximation de très grande taille (10²⁰, 10⁵⁰...).

## 2 Méthodes spectrales stochastiques

Les méthodes spectrales stochastiques fournissent une démarche générale pour la prédiction de la réponse d'une grande classe de modèles régis par des équations aux dérivées partielles stochastiques (EDPS). On considère ici que les incertitudes peuvent être modélisées par un ensemble fini de m variables aléatoires $\xi = (\xi_1, \dots, \xi_m)$, de loi de probabilité $P_\xi$, définissant un espace probabilisé $(\Xi, \mathcal{B}_\Xi, P_\xi)$, avec $\Xi \subset \mathbb{R}^m$. Les EDPS peuvent généralement se ramener à la recherche d'une variable aléatoire $u$ à valeurs dans un espace fonctionnel $\mathcal{V}$, et vérifiant presque sûrement

$$u(\xi) \in \mathcal{V}, \quad a(u(\xi), v; \xi) = b(v; \xi) \quad \forall v \in \mathcal{V} \quad (1)$$

où $a(\cdot, \cdot; \xi)$ est une forme bilinéaire (voire semi-linéaire) sur $\mathcal{V}$ et $b(\cdot; \xi)$ est une forme linéaire sur $\mathcal{V}$. Les approches spectrales stochastiques consistent à rechercher une approximation de la solution sous la forme $u(\xi) = \sum_{\alpha=1}^{P} u_\alpha H_\alpha(\xi)$, où les $H_\alpha$ forment une base d'approximation de $S_P \subset S := L^2(\Xi, dP_\xi)$, généralement choisie orthonormée pour le produit scalaire naturel sur $L^2(\Xi, dP_\xi)$. Pour la définition de cette approximation, les approches de type Galerkin s'appuient sur une formulation faible de (1) au niveau stochastique :

$$u \in \mathcal{V} \otimes S, \quad A(u, v) = B(v) \quad \forall v \in \mathcal{V} \otimes S \quad (2)$$

avec

$$A(u, v) := \int_{\Xi} a(u(\mathbf{y}), v(\mathbf{y}); \mathbf{y}) dP_{\xi}(\mathbf{y}), \quad B(v) := \int_{\Xi} b(v(\mathbf{y}); \mathbf{y}) dP_{\xi}(\mathbf{y}) \quad (3)$$

## 3 Décomposition spectrale généralisée

### 3.1 Principe

Le principe de la méthode de décomposition spectrale généralisée (GSD) [4, 5] consiste à rechercher une représentation optimale de la solution sous la forme $u_M = \sum_{i=1}^{M} w_i \lambda_i$, avec $w_i \in \mathcal{V}$ et $\lambda_i \in S$, et où ni les fonctions déterministes $w_i$ ni les fonctions stochastiques $\lambda_i$ ne sont choisies a priori. La méthode peut être interprétée comme une méthode de réduction de modèle a priori, conduisant à la construction d'une base réduite de fonctions déterministes et stochastiques, sans connaître a priori la solution du problème (2). Rappelons brièvement les bases de cette approche.

On note $W = (w_i)_{i=1}^M \in (\mathcal{V})^M$ et $\Lambda = (\lambda_i)_{i=1}^M \in (S)^M$ les bases réduites déterministes et stochastiques recherchées, et on note par convention $u_M = W\Lambda$ la décomposition d'ordre M associée. On introduit alors les applications $F : (S)^M \to (\mathcal{V})^M$ et $f : (\mathcal{V})^M \to (S)^M$ respectivement définies par

$$A(F(\Lambda)\Lambda, W^*\Lambda) = B(W^*\Lambda) \quad \forall W^* \in (\mathcal{V})^M, \quad (4)$$

$$A(Wf(W), W\Lambda^*) = B(W\Lambda^*) \quad \forall \Lambda^* \in (S)^M \quad (5)$$
Pour $\Lambda$ fixé, $W = F(\Lambda)$ est solution d'une équation qui s'apparente à une approximation de Galerkin du problème (2) dans l'espace réduit $\mathcal{V} \otimes S_M$, avec $S_M := \text{span}(\Lambda)$. De même, pour $W$ fixé, $\Lambda = f(W)$ est solution d'une équation qui s'apparente à une approximation de Galerkin du problème (2) dans l'espace réduit $\mathcal{V}_M \otimes S$, avec $\mathcal{V}_M := \text{span}(W)$.

Un choix optimal des bases réduites d'approximation, au sens de Galerkin, consiste à choisir $W$ et $\Lambda$ tels qu'on ait simultanément $W = F(\Lambda)$ et $\Lambda = f(W)$. Ces deux conditions peuvent s'écrire sous la forme de deux problèmes équivalents :

$$W = T(W) \quad \text{ou} \quad \Lambda = T^{\diamond}(\Lambda) \qquad (6)$$

avec $T = F \circ f$ et $T^{\diamond} = f \circ F$. Ces deux problèmes peuvent être interprétés comme des problèmes aux valeurs propres [5]. On montre dans de nombreux cas que l'espace propre dominant de dimension $M$ de $T$ (resp. $T^{\diamond}$), noté $\mathcal{V}_M$ (resp. $S_M$), conduit à une décomposition spectrale généralisée optimale $Wf(W)$ (resp. $F(\Lambda)\Lambda$).

## 3.2 Algorithmes

Les algorithmes GSD proposés dans [4, 5] s'inspirent des techniques de résolution de problèmes aux valeurs propres classiques (méthode des puissances, Arnoldi, ...) et permettent d'obtenir une bonne approximation de l'espace propre dominant de l'opérateur $T$ (ou $T^{\diamond}$). Certaines variantes de ces algorithmes ne demandent que la résolution de problèmes déterministes (calcul du type $w = F(\lambda)$ pour un $\lambda$ donné), résolu de manière classique, ainsi que la résolution d'équations stochastiques algébriques (calcul du type $\lambda = f(w)$ pour $w$ donné, voire calcul de $\Lambda = f(W)$ pour $W$ donné). En assimilant $\Lambda \in (S)^M$ à un vecteur $\Lambda \in \mathbb{R}^M \otimes S$, le problème définissant $\Lambda = f(W)$ pour $W$ donné peut s'écrire :

$$\int_{\Xi} \Lambda(\mathbf{y})^{*T} \mathbf{A}(\mathbf{y}) \Lambda(\mathbf{y}) dP_{\xi}(\mathbf{y}) = \int_{\Xi} \Lambda(\mathbf{y})^{*T} \mathbf{b}(\mathbf{y}) dP_{\xi}(\mathbf{y}) \quad \forall \Lambda^* \in \mathbb{R}^M \otimes S \qquad (7)$$

où $(\mathbf{A}(\xi))_{ij} = a(w_j, w_i; \xi)$ et $(\mathbf{b}(\xi))_i = b(w_i; \xi)$. L'équation (7) est une formulation faible d'un système de $M$ équations algébriques stochastiques.

L'utilisation d'une approche spectrale stochastique classique consiste alors à injecter dans (7) un espace d'approximation $S_P \subset S$. Dans le cas d'une grande dimension stochastique $m$, la dimension du problème (7) devient prohibitive, en raison de la grande dimension de $S_P$.

# 4 Décomposition en variables séparées pour la résolution d'équations algébriques stochastiques

Dans le cas de variables aléatoires indépendantes, l'espace fonctionnel au niveau stochastique a une structure produit tensoriel : $L^2(\Xi, dP_\xi) \simeq L^2(\Xi_1, dP_{\xi_1}) \otimes \dots \otimes L^2(\Xi_m, dP_{\xi_m})$. On propose ici une extension naturelle de la méthode de décomposition spectrale généralisée pour l'approximation de $\Lambda$ sous la forme d'une série de fonctions à variables séparées :

$$\Lambda_Z(\xi) = \sum_{k=1}^{Z} \varphi_k^0 \varphi_k^1(\xi_1) \dots \varphi_k^m(\xi_m), \quad \varphi_k^0 \in \mathbb{R}^M, \quad \varphi_k^i \in L^2(\Xi_i, dP_{\xi_i}) \qquad (8)$$

Différents algorithmes ont été proposés pour la construction de cette décomposition en variables séparées [7, 8, 9], qui s'interprète comme une décomposition spectrale généralisée (ou décomposition en valeurs singulières généralisée). Si la solution est définie par un problème de minimisation (e.g. si (2) est un problème variationnel ou si on utilise une reformulation de (7) sous la forme d'une minimisation du résidu), les algorithmes classiques s'interprètent comme des procédures de
minimisation alternées. Pour $i \in \{0, ..., m\}$ donné, la minimisation est opérées sur les fonctions $\varphi_k^i$, les fonctions $\varphi_k^j$ étant fixées pour $j \neq i$. Ces algorithmes permettent dans notre cas de découpler la résolution d'équations stochastiques algébriques définies sur des espaces de probabilité de dimension 1. Les Z ensembles de $\{\varphi_k^0, \varphi_k^1, ..., \varphi_k^m\}$ de la décomposition (8) peuvent être recherchés les uns après les autres ou bien simultanément, cette dernière option conduisant à une convergence en Z souvent sensiblement meilleure.

En introduisant des espaces d'approximation $S_p^{(i)} \subset L^2(\Xi_i, dP_{\xi_i})$ de dimension $p$, l'algorithme avec construction simultanée des Z ensembles de fonctions ne nécessite que la résolution de systèmes d'équations de taille $MZ$ et $pZ$. Avec une construction successive des ensembles de fonctions, il nécessite la résolution de systèmes de taille $M$ ou $p$. Une approche conventionnelle, tentant de résoudre directement l'équation (7) dans l'espace d'approximation $S_P = S_p^{(1)} \otimes \dots \otimes S_p^{(m)}$, nécessiterait la résolution (impossible en pratique) d'un système de taille $Mp^m$ (par exemple $10^{31}$ pour $(M, p, m) = (10, 10, 30)$). Par ces techniques classiques, on introduit de préférence des espaces d'approximation obtenus par tensorisation creuse des espaces $S_p^{(i)}$. Cela permet de réduire la dimension de l'espace d'approximation à $M(p+m)!/p!m!$ (soit $8.10^9$ pour l'exemple ci-dessus), ce qui est souvent insuffisant en pratique.

Pour les approches traditionnelles, le coût de stockage mémoire de la solution est également limitant. Avec la stratégie proposée, seulement $Mn + Z(M + pm)$ réels doivent être stockés, où $n$ désigne la dimension de l'espace d'approximation au niveau déterministe.

# 5 Illustration

On considère un problème d'advection-diffusion défini sur le domaine $\Omega = (0, 1) \times (0, 1)$. On impose des conditions de Dirichlet homogènes sur $\partial\Omega$ et considère un terme source volumique uniforme. On pose donc pour ce problème $\mathcal{V} = H_0^1(\Omega)$ et

$$a(u, v; \xi) = \int_{\Omega} \mu(\xi) \nabla u \cdot \nabla v + \int_{\Omega} c(\xi) \cdot \nabla u \cdot \nabla v \quad (9)$$

$$l(v) = \int_{\Omega} v \quad (10)$$

où $\mu(\xi)$ et $c(\xi) = (c_x(\xi), c_y(\xi))$ dépendant de $3m'$ variables aléatoires indépendantes $\xi = (\xi_i)_{i=1}^{3m'}$ de distribution uniforme sur $(-1/m', 1/m')$. On choisit la paramétrisation suivante : $\mu(\xi) = 1 + 0.3 \sum_{i=1}^{m'} \xi_i$, $c_x(\xi) = 1 + 0.3 \sum_{i=m'+1}^{2m'} \xi_i$ et $c_y(\xi) = 1 + 0.3 \sum_{i=2m'+1}^{3m'} \xi_i$. On choisit pour cet exemple $m' = 6$, si bien que l'espace de probabilité $(\Xi, \mathcal{B}_{\Xi}, P_{\xi})$, avec $\Xi = (-1, 1)^m$, est de dimension $m = 18$. On choisit une approximation éléments finis au niveau spatial en introduisant un espace $\mathcal{V}_n \subset \mathcal{V}$ de dimension 1240 et on choisit une approximation par chaos polynomial de degré 4 suivant chaque dimension stochastique, ce qui correspond à des espaces $S_p^{(i)}$ de dimension $p = 5$. L'espace d'approximation $S_P = S_p^{(1)} \otimes \dots \otimes S_p^{(m)}$ est de dimension $5^{18} \approx 3.10^{12}$. L'approximation de Galerkin dans $\mathcal{V}_n \otimes S_P$ n'est donc bien sûr pas calculable !

La méthode proposée parvient à capturer une bonne approximation de la solution sur une base réduite déterministe de faible dimension. Pour $M = 6$ fonctions déterministes, et $Z = 4$ fonctions à variables séparées pour l'approximation de $\Lambda_Z(\xi)$ (voir équation (8)), l'erreur, calculée en évaluant le résidu de l'équation (2), est inférieure à $10^{-2}$. Les 6 fonctions déterministes $w_i$ sont montrées sur la figure 1. Sur la figure 2, on trace la densité de probabilité de la solution au point $(1/2, 1/2)$ pour une décomposition en $M = 1, 3$ et 6 modes. Ces densités sont comparées avec celle obtenue par une approche de Monte-Carlo classique avec 20000 tirages. On observe une convergence très rapide de cette quantité d'intérêt en fonction de $M$.
<figure><img src="image_2.png" /><figcaption>Figure 1 – Base réduite déterministe calculée : fonctions $\{w_i\}_{i=1}^6$</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>Figure 2 – Densité de probabilité de la solution au point (1/2, 1/2). Comparaison Monte-Carlo et décompositions GSD d'ordre M = 1, 3 et 6.</figcaption></figure>

# 6 Conclusion

On a proposé une méthode de résolution d'équations aux dérivées partielles stochastiques rendant possible le traitement de problèmes à grande dimension stochastique. La méthode proposée conjugue une approche par décomposition spectrale généralisée (GSD) et une approche par séparation de variables en grande dimension, qui exploite la structure produit tensoriel de l'espace des fonctions stochastiques. Elle permet d'envisager la résolution de problèmes ayant une dimension telle que des approches standards deviennent inutilisables.

# Remerciements

Ces travaux sont soutenus par l'Agence National de la Recherche (projet ANR-06-JCJC-0064).
Références

[1] R. Ghanem and P. Spanos. *Stochastic finite elements : a spectral approach*. Springer, Berlin, 1991.

[2] H. G. Matthies and A. Keese. Galerkin methods for linear and nonlinear elliptic stochastic partial differential equations. *Computer Methods in Applied Mechanics and Engineering*, 194(12-16):1295-1331, 2005.

[3] A. Nouy. Recent developments in spectral stochastic methods for the numerical solution of stochastic partial differential equations. *Archives of Computational Methods in Engineering*, 2009. To appear.

[4] A. Nouy. A generalized spectral decomposition technique to solve a class of linear stochastic partial differential equations. *Computer Methods in Applied Mechanics and Engineering*, 196(45-48):4521-4537, 2007.

[5] A. Nouy. Generalized spectral decomposition method for solving stochastic finite element equations : invariant subspace problem and dedicated algorithms. *Computer Methods in Applied Mechanics and Engineering*, 197:4718-4736, 2008.

[6] A. Nouy and O. Le Maître. Generalized spectral decomposition method for stochastic non linear problems. *Journal of Computational Physics*, 228(1):202-235, 2009.

[7] G. Beylkin and M. J. Mohlenkamp. Algorithms for numerical analysis in high dimensions. *SIAM J. Sci. Comput.*, 26(6):2133-2159, 2005.

[8] A. Doostan, G. Iaccarino, and N. Etemadi. A least-squares approximation of high-dimensional uncertain systems. Technical report, Center for Turbulence Research, Annual Research Briefs, 2007.

[9] A. Ammar, B. Mokdad, F. Chinesta, and R. Keunings. A new family of solvers for some classes of multidimensional partial differential equations encountered in kinetic theory modelling of complex fluids. *Journal of Non-Newtonian Fluid Mechanics*, 139(3):153-176, 2006.