# Réduction de modèles par CPOD et Krigeage
## Application à l'optimisation de forme d'un conduit d'admission

Manyu Xiao, Rajan Filomeno Coelho, Piotr Breitkopf, Catherine Knopf-Lenoir, Pierre Villon, Maryan Sidorkiewicz

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

To cite this version:

Manyu Xiao, Rajan Filomeno Coelho, Piotr Breitkopf, Catherine Knopf-Lenoir, Pierre Villon, et al. Réduction de modèles par CPOD et Krigeage Application à l'optimisation de forme d'un conduit d'admission. 9e Colloque national en calcul des structures, CSMA, May 2009, Giens, France. hal-01408405

HAL Id: hal-01408405

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

Submitted on 4 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_2.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Réduction de modèles par CPOD et Krigeage
## Application à l'optimisation de forme d'un conduit d'admission

Manyu Xiao¹, Rajan Filomeno Coelho¹, Piotr Breitkopf¹, Catherine Knopf-Lenoir¹, Pierre Villon¹, Maryan Sidorkiewicz²

¹ Laboratoire Roberval, UTC – CNRS (UMR 6253)

Université de Technologie de Compiègne

Centre de Recherches de Royallieu, BP 20259, Compiègne 60205 Cedex

{manyu.xiao, rajan.filomeno-coelho, piotr.breitkopf, catherine.vayssade, pierre.villon}@utc.fr

² Renault, Technocentre, Avenue du Golf 1, 78288 Guyancourt

maryan sidorkiewicz@renault.com

**Résumé** – Une technique de réduction de modèles combinant une variante contrainte de la méthode de décomposition aux valeurs propres avec un métamodèle par krigeage est proposée pour représenter de manière économique et précise des écoulements au sein de la chambre de combustion d'un moteur thermique. Le modèle généré est ensuite directement exploitable pour l'optimisation multicritère de la forme du conduit afin de trouver un front de Pareto pour deux objectifs antagonistes : la puissance du moteur et le respect de l'environnement.

**Mots clés** – Optimisation multicritère, réduction de modèles, krigeage, POD, POD contrainte, conduit d'admission.

# 1 Introduction

Dans cette étude, on s'intéresse au développement et à l'utilisation de modèles réduits pour l'optimisation. Le calcul de la fonction objectif est coûteux en temps CPU si une simulation numérique haute fidélité (éléments finis, volumes finis, etc.) est utilisée. Même si avec l'aide de plates-formes de calcul parallèle on peut diminuer le volume de calcul, celui-ci reste trop élevé à cause du nombre de simulations et du coût de chaque simulation. De plus, des contraintes informatiques peuvent également entrer en considération (gestion des licences et environnements de calcul, utilisation efficace de plates-formes de calcul à haute performance).

En conséquence, une approche largement utilisée consiste à construire des métamodèles associant les variables de conception v (par exemple les paramètres définissant la géométrie d'une pièce mécanique) aux réponses y (masse, tenue mécanique, etc.). Une large gamme de méthodes d'approximation a été proposée à cet effet dans la littérature, parmi lesquelles on peut citer les surfaces de réponse polynomiales, le krigeage, l'approximation diffuse, les réseaux de neurones ou encore les réseaux de fonctions à base radiale [1].

Néanmoins, cette méthodologie comporte une limitation importante : pour des phénomènes physiques complexes, ou lorsque le nombre de variables de conception est élevé, des surfaces de réponse liant directement les variables de conception v aux réponses y peuvent lisser
abusivement le comportement réel. Ce phénomène a été observé notamment dans le cas de l'optimisation d'une aile flexible 2D soumise à un couplage fluide-structure [2].

Pour construire des modèles réduits de meilleure qualité, au lieu de représenter les grandeurs globales, l'approximation peut être réalisée directement sur les champs vectoriels résultant de la simulation numérique. Cependant, comme ces grandeurs comportent généralement un grand nombre de composantes (e.g. le champ des vitesses nodales pour un maillage tridimensionnel), une technique de décomposition est envisagée.

Dans cet article, une stratégie de réduction de modèles est proposée pour répondre à ces problématiques, combinant une stratégie originale de décomposition aux valeurs propres et une méthode d'approximation par krigeage des coefficients. Ensuite, l'application de la méthodologie à l'écoulement du mélange air-combustible dans une chambre de combustion d'un moteur Renault est présentée.

## 2 Réduction de modèles par POD contrainte

La première étape consiste donc à réduire un champ $\mathbf{p}=(p_1, p_2, ..., p_n)^T$, que l'on supposera de taille $n$ (correspondant par exemple au nombre de nœuds d'un maillage). L'objectif est de caractériser ce champ $\mathbf{p}$ à l'aide d'un nombre restreint de paramètres scalaires.

A cet effet, une décomposition POD (Proper Orthogonal Decomposition) [3] est investiguée à partir d'un plan d'expériences {$\mathbf{p}^{(k)}$; $k = 1, ..., M$} déterminé par exécution d'un modèle « haute fidélité» pour un échantillon représentatif de «clichés» correspondant aux points dans l'espace de conception. Le but est de déterminer une base $\mathbf{\Phi}$ telle que :

$$ \mathbf{p}^{(k)} = \bar{\mathbf{p}} + \sum_{i=1}^{M} \alpha_i^{(k)} \varphi_i = \bar{\mathbf{p}} + \mathbf{\Phi} \alpha^{(k)}, \quad (1) $$

expression dans laquelle $\bar{\mathbf{p}}$ est la contribution moyenne sur les $M$ clichés, $\mathbf{\Phi}$ est la base POD constituée de $M$ vecteurs $\varphi_i$ appelés modes et $\alpha^{(k)}$ sont les coefficients de la combinaison linéaire. Le calcul de $\mathbf{\Phi}$ est décrit dans [4] ; on observe qu'en général, l'essentiel de la contribution énergétique est localisé dans les premiers modes [5]. L'approximation POD consiste alors à ne conserver qu'un nombre restreint $m$ de modes :

$$ \mathbf{p}^{\text{POD}, (k)} = \bar{\mathbf{p}} + \sum_{i=1}^{m} \alpha_i^{(k)} \varphi_i = \bar{\mathbf{p}} + \mathbf{\Phi}_{,m} \alpha^{(k)}. \quad (2) $$

Néanmoins, cette troncature de la base entraîne inévitablement une non-conservation de grandeurs intégrales, notamment celles qui seraient définies comme des applications linéaires $\ell$ sur $\mathbf{p}$ :

$$ q^{\text{POD}, (k)} = \ell(\mathbf{p}^{\text{POD}, (k)}) \neq q^{(k)} = \ell(\mathbf{p}^{(k)}). \quad (3) $$

Pour obtenir une approximation du vecteur $\mathbf{p}$ qui garantit la préservation des applications linéaires (3), une méthode originale consiste à remplacer les coefficients $\alpha^{(k)}$ par les coefficients $\gamma^{(k)}$, en gardant la même base $\mathbf{\Phi}$. Les coefficients $\gamma^{(k)} = (\gamma_1^{(k)} ... \gamma_m^{(k)})^T$, $k= 1, ..., m$ sont déterminés de manière à minimiser l'erreur $J(\gamma^{(k)})$ commise par l'approximation, tout en imposant la conservation des $r$ applications linéaires $\ell_j, j = 1, ..., r$ :

$$ J(\gamma^{(k)}) = \frac{1}{2} \| \mathbf{p}^{(k)} - (\bar{\mathbf{p}} + \mathbf{\Phi}_{,m} \gamma^{(k)}) \|^2 \quad (4) $$

$$ \text{soumis à : } \ell_j(\bar{\mathbf{p}} + \mathbf{\Phi}_{,m} \gamma^{(k)}) = \ell_j(\mathbf{p}^{(k)}) = \mathbf{c}_j; j = 1, \dots, r. \quad (5) $$

On peut réécrire le système (4-5) de la manière suivante :
$$
\left\{
\begin{aligned}
J(\gamma^{(k)}) &= \frac{1}{2} \int_{\Omega} (\mathbf{p}^{(k)}(\mathbf{x}) - p^{\text{CPOD}, (k)}(\mathbf{x})) d\Omega && (6) \\
\text{soumis à : } \ell_j(\bar{\mathbf{p}} + \mathbf{\Phi}_{,m} \gamma^{(k)}) &= \ell_j(\mathbf{p}^{(k)}) = \mathbf{c}_j; \quad j = 1, \dots, r. && (7)
\end{aligned}
\right.
$$

Cette stratégie porte le nom de POD contrainte, ou CPOD (Constrained Proper Orthogonal Decomposition).

En utilisant en un point $\mathbf{x} = (x, y, z)^T$ l'approximation par éléments finis $p^{(k)}(\mathbf{x}) = \mathbf{N}(\mathbf{x})^T \cdot \mathbf{p}^{(k)}$, $p^{\text{CPOD}, (k)}(\mathbf{x}) = \mathbf{N}(\mathbf{x})^T \cdot \mathbf{p}^{\text{CPOD}, (k)} = \mathbf{N}(\mathbf{x})^T \cdot (\bar{\mathbf{p}} + \mathbf{\Phi}_{,m} \gamma^{(k)})$, où $p^{(k)}(\mathbf{x})$ et $p^{\text{CPOD}, (k)}(\mathbf{x})$ sont les champs continus, $\mathbf{N}(\mathbf{x})$ les fonctions de forme et $\mathbf{p}^{(k)}, \mathbf{p}^{\text{CPOD}, (k)}$ les vecteurs de valeurs nodales, on peut réécrire le système (6-7) avec la matrice $\mathbf{W}$ définissant les $r$ applications linéaires :

$$
\left\{
\begin{array}{l}
J(\gamma^{(k)}) = \frac{1}{2} \gamma^{(k)\mathrm{T}} \mathbf{\Phi}_{,m} \mathrm{T} \left[ \underbrace{\int_{\Omega} \mathbf{N}\mathbf{N}^{\mathrm{T}} d\Omega}_{= \mathbf{M} \text{ (matrice de masse)}} \right] \mathbf{\Phi}_{,m} \gamma^{(k)} - (\mathbf{p}^{(k)} - \bar{\mathbf{p}})^{\mathrm{T}} \left[ \int_{\Omega} \mathbf{N}\mathbf{N}^{\mathrm{T}} d\Omega \right] \mathbf{\Phi}_{,m} \gamma^{(k)} + \frac{1}{2} (\mathbf{p}^{(k)} - \bar{\mathbf{p}})^{\mathrm{T}} \left[ \int_{\Omega} \mathbf{N}\mathbf{N}^{\mathrm{T}} d\Omega \right] (\mathbf{p}^{(k)} - \bar{\mathbf{p}}) \\
\hfill (8)
\end{array}
\right.
$$

$$
\text{soumis à : } (\mathbf{W}^{\mathrm{T}} \cdot (\bar{\mathbf{p}} + \mathbf{\Phi}_{,m} \gamma^{(k)})) = \mathbf{c}. \qquad (9)
$$

En définissant le lagrangien $\mathcal{L}(\gamma^{(k)}, \boldsymbol{\lambda}) = J(\gamma^{(k)}) + \boldsymbol{\lambda}^{\mathrm{T}} \cdot (\mathbf{W}^{\mathrm{T}} \cdot (\bar{\mathbf{p}} + \boldsymbol{\Phi}_{,m}\gamma^{(k)}) - \mathbf{c})$ et en écrivant les conditions d'optimalité du premier ordre par rapport à $\gamma^{(k)}$ et à $\boldsymbol{\lambda}$, on aboutit au système suivant :

$$
\left\{
\begin{array}{l}
\displaystyle \frac{\partial \mathcal{L}}{\partial \gamma^{(k)}} = \Phi_{,m}^{\mathrm{T}} \mathbf{M} \Phi_{,m} \gamma^{(k)} - \Phi_{,m}^{\mathrm{T}} \mathbf{M} (\mathbf{p}^{(k)} - \bar{\mathbf{p}}) + \Phi_{,m}^{\mathrm{T}} \mathbf{W} \boldsymbol{\lambda} = 0 & (10) \\[1em]
\displaystyle \frac{\partial \mathcal{L}}{\partial \boldsymbol{\lambda}} = \mathbf{W}^{\mathrm{T}} \cdot (\bar{\mathbf{p}} + \Phi_{,m} \gamma^{(k)}) - \mathbf{c} = 0. & (11)
\end{array}
\right.
$$

$$
\frac{\partial \mathcal{L}}{\partial \boldsymbol{\lambda}} = \mathbf{W}^{\mathrm{T}} \cdot (\bar{\mathbf{p}} + \Phi_{,m} \gamma^{(k)}) - \mathbf{c} = 0. \quad (11)
$$

Sous forme matricielle, on aboutit finalement à l'expression :

$$
\left[
\begin{array}{cc}
\Phi_{,m}^{\mathrm{T}} \mathbf{M} \Phi_{,m} & \Phi_{,m}^{\mathrm{T}} \mathbf{W} \\
\mathbf{W}^{\mathrm{T}} \Phi_{,m} & \mathbf{O}
\end{array}
\right]
\left\{
\begin{array}{c}
\gamma^{(k)} \\
\lambda
\end{array}
\right\}
=
\left\{
\begin{array}{c}
\Phi_{,m}^{\mathrm{T}} \mathbf{M} (\mathbf{p}^{(k)} - \bar{\mathbf{p}}) \\
\mathbf{c} - \mathbf{W}^{\mathrm{T}} \bar{\mathbf{p}}
\end{array}
\right\} .
\quad (12)
$$

La dimension du système est égale à la somme du nombre de modes *m* et du nombre de restrictions *r*.

3 Métamodèles par POD contrainte et krigeage

En pratique, la stratégie complète de réduction de modèles, initialement proposée dans le cas
de problèmes couplés fluide-structure [6] et utilisée ici en combinaison avec la CPOD, se divise
en quatre étapes :

1. exécution d'un plan d'expériences de simulations numériques « haute fidélité » pour obtenir *M* clichés ;

2. calcul de la contribution moyenne $\bar{\mathbf{p}}$ et de la base POD $\mathbf{\Phi}$ ;

3. calcul des coefficients $\gamma^{(k)}$ par CPOD (12) ;

4. à partir des données calculées aux points 2 et 3, un métamodèle généraliste permet de déterminer pour n'importe quel jeu de variables de conception **v** une approximation $\tilde{\gamma}^{(k)}(\mathbf{v})$ des coefficients CPOD correspondants, et donc de reconstruire le
champ **p**, et d'en déduire les grandeurs post-traitées. Dans ce travail, c'est un métamodèle basé sur le krigeage qui est utilisé.

La section suivante illustre l'utilisation de ces modèles réduits sur dans le cas du dimensionnement d'un conduit d'admission équipant un moteur Renault.

## 4 Application à l'optimisation multiobjectif d'un conduit d'admission

L'application industrielle étudiée est un conduit d'admission guidant le mélange dans la chambre de combustion. Le but est de maximiser deux critères antagonistes : le *tumble t* (ou *vorticité*), caractérisant la qualité du mélange de l'air et du combustible dans la chambre de combustion (réduction de la pollution), et le *débit q* (performance). Ces deux grandeurs dépendent de la géométrie du conduit pilotée par six variables CAO $\mathbf{v} = (v_1, ..., v_6)^T$, et sont calculées en fonction du champ de vitesse dans une zone fixe située dans le cylindre (voir figure 1, à gauche).

<figure><img src="image_3.png" /><figcaption>Figure 1 - lignes d'écoulement obtenues par simulation CFD « haute fidélité» (à gauche) et front de Pareto obtenu par métamodèle de krigeage sur les coefficients CPOD (à droite).</figcaption></figure>

<table><caption>Tableau 1- Paramètres de MOGA et Krigeage</caption><tr><td rowspan="2">MOGA</td><td>Nombre d'itérations</td><td>50</td></tr><tr><td>La taille de population initiale</td><td>250 (146 points du plan d'expériences + 104 points aléatoires)</td></tr><tr><td rowspan="2">Krigage</td><td>Corrélation</td><td>Gaussienne</td></tr><tr><td>Régression</td><td>polynôme d'ordre 0</td></tr></table>

L'optimisation multicritère a été réalisée à l'aide de l'algorithme génétique MOGA [7] disponible dans l'environnement DAKOTA développé par Sandia Labs en 1994, et dont l'implémentation se classe dans la catégorie des algorithmes parallèles à gros grain [8]. Elle utilise une population initiale étant constituée notamment des points du plan d'expériences, pour progressivement converger vers les solutions de compromis entre les deux critères (également dénommées *front de Pareto*).

Les algorithmes de krigeage sont ceux de la boîte à outils DACE [9]. Les paramètres de MOGA et du krigeage sont donnés sur le tableau 1.

Les résultats de l'optimisation multiobjectif sont présentés sur la figure 1 (à droite). On constate d'emblée le caractère interpolant du métamodèle : les valeurs obtenues pour les points du plan d'expérience coïncident exactement avec les points initiaux calculés à partir de la
simulation « haute fidélité ». Le résultat obtenu domine le plan d'expériences, on trouve donc un front de Pareto. Des résultats plus détaillés, incluant notamment une analyse fine des paramètres de la méthode (nombre de modes choisis pour la CPOD, l'influence du choix des paramètres du krigeage) seront présentés et discutés lors de l'exposé oral.

5 Conclusions

Dans ce travail, une stratégie originale de décomposition aux valeurs propres contrainte (CPOD) a été proposée, combinée à un métamodèle (krigeage). Cette approche permet de décomposer une grandeur physique (champs de vitesses, de pression) tout en assurant la préservation de réponses définies par applications linéaires (débit, portance, etc.). La méthodologie de réduction de modèles ainsi proposée est ensuite appliquée à l'optimisation multiobjectif d'un conduit d'admission développé par Renault, illustrant l'intérêt de la CPOD pour représenter des champs complexes préservant des grandeurs intégrales comme le débit ou le tumble.

6 Remerciements

Ce travail a été soutenu par le Réseau National en Technologies Logicielles (RNTL) de
l'Agence Nationale de la Recherche (ANR), dans le cadre du projet OMD (« Optimisation
Multi-Disciplinaire », référence ANR-05-RNTL-009).

Références

[1] D. R. Jones. A Taxonomy of Global Optimization Methods Based on Response Surfaces, *Journal of Global Optimization*, 21: 345–383, 2001.

[2] R. Filomeno Coelho, P. Breitkopf, C. Knopf-Lenoir. Model reduction for multidisciplinary optimization – Application to a 2D wing, *Structural and Multidisciplinary Optimization*, 10.1007/s00158-007-0212-5, 2007.

[3] G. Berkooz, P. Holmes, J. L. Lumley. The Proper Orthogonal Decomposition in the Analysis of Turbulent Flows, *Annual Reviews of Fluid Mechanics*, vol. 25, pp. 593–575, 1993.

[4] A. Newman. Model reduction via the Karhunen-Loève expansion, Part I: An exposition, Technical Report TR96-32, *Inst. Systems Research*, avril 1996.

[5] P. A. LeGresley, J. J. Alonso. Improving the Performance of Design Decomposition Methods with POD, *10th AIAA/ISSMO Multidisciplinary Analysis and Optimization Conference*, 30 août – 1er septembre 2004, Albany, New York, USA, AIAA 2004-4465.

[6] R. Filomeno Coelho, P. Breitkopf, C. Knopf-Lenoir, P. Villon. Bi-level Model Reduction for Coupled Problems – Application to a 3D wing, *Structural and Multidisciplinary Optimization*, DOI 10.1007/s00158-008-0335-3, 2008.

[7] J. Eddy, K. Lewis. Effective generation of Pareto sets using genetic programming, *Proc. of DETC'01 ASME 2001 Design Engineering Technical Conferences and Computers and Information in Engineering Conference Pittsburgh*, PA, September 9-12, 2001.
[8] M. S. Eldred, B. D. Schimel. Extended parallelism models for optimization on massively parallel computers, *In Proc. 3rd World Congress of Structural and Multidisciplinary Optimization (WCSMO-3)*, number 16-POM-2, Amherst, NY, May 17–21, 1999.

[9] S.N. Lophaven, H.B. Nielsen, J. Søndergaard, Aspects of the Matlab Toolbox DACE. Report IMM-REP-2002-13, Informatics and Mathematical Modelling, DTU, 2002.