# Analyses multi-échelles du premier et du second ordre des instabilités dans les matériaux hétérogènes

Saeid Nezamabadi, Hamid Zahrouni, Michel Potier-Ferry

► To cite this version:

Saeid Nezamabadi, Hamid Zahrouni, Michel Potier-Ferry. Analyses multi-échelles du premier et du second ordre des instabilités dans les matériaux hétérogènes. 11e colloque national en calcul des structures, CSMA, May 2013, Giens, France. hal-01722098

HAL Id: hal-01722098

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

Submitted on 2 Mar 2018

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 CC0 1.0 - Universal - International License
# Analyses multi-échelles du premier et du second ordre des instabilités dans les matériaux hétérogènes

Saeid NEZAMABADI<fnref n="1" />, Hamid ZAHROUNI<fnref n="2" />*, Michel POTIER-FERRY<fnref n="2" />

**Résumé** — On présente dans cet article une modélisation des instabilités dans les matériaux hétérogènes en considérant les approches numériques multi-échelles du premier et du second ordre. Les techniques proposées combinent la méthode des éléments finis multi-échelle (EF<fnref n="2" />) et la méthode asymptotique numérique (MAN). La méthode EF<fnref n="2" /> nous permet une grande souplesse quant au choix de la morphologie des microstructures et la MAN nous permet de suivre avec une grande efficacité les branches solutions même complexes et avec un temps de calcul réduit. L'apport de l'homogénéisation du second ordre est mis en lumière au travers un exemple numérique présentant des instabilités au niveau de la microstructure.

**Mots clés** — Homogénéisation numérique, Méthode asymptotique numérique, Instabilité, Flambage

## 1 Introduction

Les matériaux composites et matériaux hétérogènes sont généralement considérés comme «multi-échelle», dans lesquels les dimensions caractéristiques de leurs constituants sont très différentes de la dimension de la structure. L'étude des phénomènes d'instabilité est un objectif clé de l'analyse multi-échelle de ces structures. Ces phénomènes peuvent survenir à la fois aux niveaux macroscopique et microscopique et s'influencer mutuellement. Dans nos précédents travaux ([1, 2]), nous avons proposé une méthode numérique, appelée «MAN-multi-échelle», pour résoudre des problèmes nonlinéaires dans le cadre des matériaux hétérogènes. Elle combine le modèle des éléments finis multi-échelle du premier ordre et une technique de perturbation basée sur la Méthode Asymptotique Numérique (MAN). Dans ces travaux, la précision et l'efficacité de la méthode proposée ont été testées et illustrées par plusieurs exemples numériques impliquant des instabilités locales et globales. Néanmoins, quelques restrictions liées à l'homogénéisation du premier ordre limitent le domaine d'utilisation de cette technique : (i) cette méthode n'est pas sensible à la taille de la microstructure malgré la prise en compte de la fraction volumique, de la distribution et de la morphologie des constituants. Dans le contexte de l'homogénéisation du premier ordre, on suppose que la taille de la microstructure est infiniment petite par rapport à celle de la macrostructure. Toutefois, il arrive que la taille de la microstructure soit comparable à la taille caractéristique de la structure macroscopique ; (ii) l'hypothèse de l'uniformité des champs macroscopiques (contrainte/déformation) associés à chaque VER, exigée dans le concept de la séparation des échelles, n'est pas appropriée dans les régions de forts gradients où les champs macroscopiques peuvent varier considérablement. En outre, si un adoucissement se produit en un point du matériau au niveau macroscopique, la solution obtenue à l'aide de l'homogénéisation du premier ordre va dépendre du maillage macroscopique. Cette dépendance vient du fait que le problème aux limites macroscopique est mal posé (voir [3, 4]).

De multiples auteurs ont travaillé sur ce sujet et ont proposé différentes techniques comme par exemple la contribution sur l'homogénéisation du second ordre proposée par Kousnetzova et ses collaborateurs [3, 4], dans laquelle le double gradient de déformation macroscopique est considéré dans la transition cinématique micro-macro. Le tenseur des contraintes macroscopiques et le tenseur des contraintes d'ordre élevé sont extraits en utilisant une extension de la condition de Hill-Mandel. Un milieu continu d'ordre supérieur est obtenu au niveau macroscopique, tandis que le problème aux limites microscopique

<footnote n="1">¹ LMGC, Université Montpellier 2, Montpellier, France, saeid.nezamabadi@univ-montp2.fr</footnote>

<footnote n="2">² LEM3, Université de Lorraine, Metz, France, hamid.zahrouni@univ-lorraine.fr & michel.potier-ferry@univ-lorraine.fr</footnote>

<footnote>* Auteur correspondant</footnote>
reste classique. Ainsi, le calcul d'homogénéisation du second ordre établit directement une réponse de comportement du second ordre au niveau microscopique sans hypothèses supplémentaires au niveau macroscopique.

Dans cette étude, le problème non linéaire résultant de l'homogénéisation du second ordre est résolu à l'aide de la MAN. Notons que la technique MAN employée pour résoudre le problème d'homogénéisation du premier ordre, reste assez similaire sur le plan algorithmique pour le problème d'homogénéisation du second ordre. La seule différence concerne les conditions aux limites qui sont un peu plus complexes dans le cadre de l'homogénéisation du second ordre car elles concernent à la fois les variables cinématiques ou mécaniques du premier ordre ainsi que celles du second ordre. Un grand avantage de la MAN est qu'elle nous permet de transformer le problème nonlinéaire micro-macro en une suite de problèmes linéaires. Sur chaque problème linéaire, on peut utiliser les outils de l'homogénéisation linéaire comme par exemple l'application du principe de superposition pour le calcul du module tangent macroscopique. Les techniques de localisation et d'homogénéisation sont utilisées une seule fois par pas de calcul ce qui permet de réduire de façon significative le temps de calcul. Rappelons que nous utilisons la technique MAN comme prédicteur d'ordre élevé sans aucune itération de correction ; les modules tangents sont donc construits une seule fois par pas.

## 2 Problème multi-échelle

Nous présentons dans cette section les formulations du problème multi-échelle du premier et du second ordre et la différence entre les deux. La résolution de ces problèmes est effectuée à l'aide de la méthode asymptotique numérique.

### 2.1 Problème microscopique

Au niveau microscopique, on considère la présence de différentes phases avec différentes caractéristiques mécaniques connues. On considère également une microstructure périodique caractérisée par un volume élémentaire représentatif (VER) qui occupe un domaine $\omega$ avec le bord extérieur $\partial\omega$ dans sa configuration de référence. Les problèmes au niveau microscopique pour les deux formulations sont similaires. La forme faible associée aux équations de ces problèmes en l'absence des forces volumiques est donnée par (voir [1]):

Trouver $\mathbf{u}$ satisfaisant les conditions aux limites microscopiques, tel que :

$$ \int_{\omega} t \mathbf{P} : \delta \mathbf{F} d\omega = 0, \qquad (1) $$

où $\mathbf{P}$ représente le premier tenseur des contraintes de Piola-Kirchhoff et $\mathbf{F}$ est le tenseur gradient de déformation défini par $\mathbf{F} = \mathbf{I} + \nabla \mathbf{u}$, où $\mathbf{I}$ est le tenseur identité d'ordre deux et $\mathbf{u}$ est le champ de déplacement. La contrainte $\mathbf{P}$ peut être définie par :

$$ \mathbf{P} = \frac{\partial \Psi^{(r)}(\mathbf{F})}{\partial \mathbf{F}}, \qquad (2) $$

où $\Psi^{(r)}$ représente la fonction d'énergie de déformation décrivant le comportement matériel de la phase (r). Selon la forme de cette fonction, diverses lois de comportement, linéaires ou non-linéaires, peuvent être prises en compte dans le cadre de notre procédure d'homogénéisation [5].

La seule différence au niveau microscopique concerne les conditions aux limites. En considérant les conditions de périodicité, les conditions aux limites microscopiques peuvent être définies pour les deux cas comme suit (voir la figure 1) :

Homogénéisation du premier ordre :

$$ \mathbf{u}^{+} - \mathbf{u}^{-} = (\bar{\mathbf{F}} - \mathbf{I}) \cdot (\mathbf{X}^{+} - \mathbf{X}^{-}) \quad \text{sur } \partial\omega, \qquad (3) $$

Homogénéisation du second ordre :

$$ \mathbf{u}^{+} - \mathbf{u}^{-} = (\bar{\mathbf{F}} - \mathbf{I}) \cdot (\mathbf{X}^{+} - \mathbf{X}^{-}) + \frac{1}{2} \bar{\mathbf{G}} : (\mathbf{X}^{+} \mathbf{X}^{-} - \mathbf{X}^{-} \mathbf{X}^{+}) \quad \text{sur } \partial\omega, \qquad (4) $$
$$ \int_{\partial\omega^{-}} \mathbf{u}^{-} d\Gamma = (\bar{\mathbf{F}} - \mathbf{I}) \cdot \int_{\partial\omega^{-}} \mathbf{X}^{-} d\Gamma + \frac{1}{2} \bar{\mathbf{G}} : \int_{\partial\omega^{-}} \mathbf{X}^{-} \mathbf{X}^{-} d\Gamma, \quad (5) $$

où $\bar{\mathbf{G}}$ est le tenseur, d'ordre trois, gradient de déformation du second ordre macroscopique défini par : $\bar{\mathbf{G}} = \nabla \cdot \bar{\mathbf{F}}$; $\mathbf{X}$ représente les coordonnées d'un point dans la configuration de référence. Notons que la notation (·) sera utilisée pour indiquer les quantités macroscopiques. Pour plus de détails sur les conditions aux limites au niveau microscopique, voir, par exemple, [1, 3].

<figure><img src="image_2.png" /><figcaption>Fig. 1 – Schéma d'un VER avec des conditions de périodicité.</figcaption></figure>

## 2.2 Problème macroscopique

On considère un domaine $\Omega$ dans $\mathbb{R}^D$, $D$ étant la dimension du domaine, avec un bord extérieur $\partial\Omega$ séparé en deux parties disjointes et complémentaires $\partial\Omega_u$ et $\partial\Omega_t$ ($\partial\Omega = \partial\Omega_t \cup \partial\Omega_u$ et $\partial\Omega_t \cap \partial\Omega_u = 0$). Sur ces deux parties sont imposés les déplacements et les efforts, respectivement. Ce domaine désigne la structure au niveau macroscopique dans sa configuration de référence. La forme faible associée aux équations du problème macroscopique en l'absence des forces volumiques pour l'homogénéisation du premier ordre est donnée par (voir [1]):

Trouver $\bar{\mathbf{u}}$, $\mathbf{u} = \hat{\mathbf{u}}$ sur $\partial\Omega_u$, tel que :

$$ \int_{\Omega} {}^{t}\bar{\mathbf{P}} : \delta\bar{\mathbf{F}} \, d\Omega = \lambda \int_{\partial\Omega_t} \mathbf{f} \cdot \delta\bar{\mathbf{u}} \, d\Gamma, \quad (6) $$

où $\hat{\mathbf{u}}$ et $\mathbf{f}$ sont les déplacements et les efforts imposés sur les bords de la macrostructure. $\lambda$ est un paramètre qui donne l'intensité de la charge appliquée.

Dans le cadre de l'homogénéisation du second ordre, la forme faible associée au problème macroscopique en l'absence des forces volumiques est définie par (voir [3]):

Trouver $\bar{\mathbf{u}}$, $\mathbf{u} = \hat{\mathbf{u}}$ sur $\partial\Omega_u$, tel que :

$$ \int_{\Omega} \{ {}^{t}\bar{\mathbf{P}} : \delta\bar{\mathbf{F}} + \bar{\mathbf{Q}} : \delta\bar{\mathbf{G}} \} \, d\Omega = \lambda \int_{\partial\Omega_t} (\mathbf{f} \cdot \delta\bar{\mathbf{u}} + {}^{t}\mathbf{r} : \delta\bar{\mathbf{F}}) \, d\Gamma, \quad (7) $$

où $\mathbf{r}$ est un tenseur de traction de second ordre défini par rapport à la contrainte du second ordre macroscopique, $\bar{\mathbf{Q}}$. Les relations de comportement entre les contraintes et les déformations sont inconnues au niveau macroscopique. Ces contraintes macroscopiques sont obtenues directement en résolvant un problème non linéaire local d'éléments finis lié à la microstructure périodique.

## 2.3 Couplage entre les échelles

Le couplage entre les échelles microscopique et macroscopique est basé sur des théorèmes de moyenne. La première relation de couplage propose que le gradient de déformation macroscopique, $\bar{\mathbf{F}}$, est égal à la valeur moyenne du gradient de déformation microscopique, $\mathbf{F}$ :

$$ \bar{\mathbf{F}} = \frac{1}{|\omega|} \int_{\omega} \mathbf{F} d\omega, \quad (8) $$

où $|\omega|$ est le volume du VER.
La deuxième relation de couplage permet de déterminer la contrainte effective $\bar{\mathbf{P}}$ à un point particulier $\bar{\mathbf{X}}$ du domaine macroscopique $\Omega$, où l'état de contrainte est supposé homogène dans le voisinage de $\bar{\mathbf{X}}$ :

$$ \bar{\mathbf{P}} = \frac{1}{|\omega|} \int_{\omega} \mathbf{P} d\omega . \qquad (9) $$

Finalement, la contrainte du second ordre macroscopique, $\bar{\mathbf{Q}}$, est liée aux variables microscopiques comme suit :

$$ \bar{\mathbf{Q}} = \frac{1}{2|\omega|} \int_{\omega} ({}^t\mathbf{P}\mathbf{X} + \mathbf{X}\mathbf{P}) d\omega . \qquad (10) $$

### 3 Résolution des problèmes à l'aide de la MAN

La résolution des systèmes d'équations non linéaires, présentés ci-dessus, dans le cadre de l'homo-généisations du premier et du second ordre est effectuée en utilisant la MAN. Dans ce cas, le problème non linéaire est transformé en une série de problèmes linéaires, ce qui permet d'utiliser le principe de superposition. Dans cet article, nous présentons seulement la résolution du problème d'homogénéisation du second ordre en sachant que celle du premier ordre est similaire.

Pour les problèmes multi-échelles non linéaires du second ordre, toutes les variables aux niveaux micro et macro sont présentées dans le vecteur mixte $\mathbf{U} = \{\bar{\mathbf{u}}, \bar{\mathbf{F}}, \bar{\mathbf{G}}, \bar{\mathbf{P}}, \bar{\mathbf{Q}}, \mathbf{u}, \mathbf{F}, \mathbf{P}, \dots\}$. Dans le cadre de la MAN, ces variables sont développées sous forme de séries entières par rapport à un paramètre de chemin "a" au voisinage d'une solution de départ connue, $(\mathbf{U}_0 \text{ et } \lambda_0)$ :

$$ \left\{ \begin{matrix} \mathbf{U}(a) \\ \lambda(a) \end{matrix} \right\} = \left\{ \begin{matrix} \mathbf{U}_0 \\ \lambda_0 \end{matrix} \right\} + \sum_{i=1}^{N} a^i \left\{ \begin{matrix} \mathbf{U}_i \\ \lambda_i \end{matrix} \right\}. \qquad (11) $$

Dans ce cadre, l'introduction du développement asymptotique (11) dans le problème non linéaire (1) et les conditions aux limites microscopiques (4) et (5) permet d'obtenir les problèmes linéaires suivants :

$$ \begin{array}{l} \text{Ordre 1 :} \\ \mathcal{L}(\mathbf{u}_1, \delta\mathbf{u}) = 0 \\ +C.L.1: \mathbf{u}_1^+ - \mathbf{u}_1^- = \bar{\mathbf{F}}_1 \cdot (\mathbf{X}^+ - \mathbf{X}^-) + \frac{1}{2}\bar{\mathbf{G}}_1 : (\mathbf{X}^+ \mathbf{X}^+ - \mathbf{X}^- \mathbf{X}^-) \\ +C.L.2: \int_{\partial\omega^-} \mathbf{u}_1^- d\Gamma = \bar{\mathbf{F}}_1 \cdot \int_{\partial\omega^-} \mathbf{X}^- d\Gamma + \frac{1}{2}\bar{\mathbf{G}}_1 : \int_{\partial\omega^-} \mathbf{X}^- \mathbf{X}^- d\Gamma \end{array} \qquad (12) $$

$$ \begin{array}{l} \text{Ordre } p > 1 : \\ \mathcal{L}(\mathbf{u}_p, \delta\mathbf{u}) = \mathcal{F}_p^{nl}(\delta\mathbf{u}) \\ +C.L.1: \mathbf{u}_p^+ - \mathbf{u}_p^- = \bar{\mathbf{F}}_p \cdot (\mathbf{X}^+ - \mathbf{X}^-) + \frac{1}{2}\bar{\mathbf{G}}_p : (\mathbf{X}^+ \mathbf{X}^+ - \mathbf{X}^- \mathbf{X}^-) \\ +C.L.2: \int_{\partial\omega^-} \mathbf{u}_p^- d\Gamma = \bar{\mathbf{F}}_p \cdot \int_{\partial\omega^-} \mathbf{X}^- d\Gamma + \frac{1}{2}\bar{\mathbf{G}}_p : \int_{\partial\omega^-} \mathbf{X}^- \mathbf{X}^- d\Gamma \end{array} \qquad (13) $$

où $\mathcal{L}$ est un opérateur tangent de rigidité qui dépend de la solution initiale :

$$ \mathcal{L}(\mathbf{u}_p, \delta\mathbf{u}) = \int_{\Omega} {}^t\nabla \mathbf{u}_p : \mathcal{H}^{(r)} : \nabla \delta\mathbf{u} \, d\Omega \quad p = 1 \cdots N. \qquad (14) $$

avec $\mathcal{H}^{(r)} = \frac{\partial^2 \Psi^{(r)}(\mathbf{F})}{\partial \mathbf{F}^2}$, et $\mathcal{F}_p^{nl}$ est un vecteur résidu défini par :

$$ \mathcal{F}_p^{nl}(\delta\mathbf{u}) = - \int_{\Omega} {}^t\mathbf{P}_p^{nl} : \nabla \delta\mathbf{u} \, d\Omega . \qquad (15) $$

Dans cette équation, $\mathbf{P}_p^{nl}$ est un terme résiduel qui contient des termes dépendant des ordres inférieures à $p$ ($p > 1$).
## 3.1 Localisation et homogénéisation

On résout partiellement les problèmes linéaires au niveau microscopique à chaque ordre en considérant les conditions de couplage entre les problèmes macroscopiques et microscopiques. Dans le cadre de l'homogénéisation du second ordre, la solution obtenue permet de construire numériquement deux tenseurs de localisation : le tenseur de localisation du premier ordre et le tenseur de localisation du second ordre.

Le déplacement microscopique à l'ordre $p$, $\mathbf{u}_p$, peut être décomposé sous forme d'une combinaison linéaire des modes solutions des problèmes aux limites élémentaires résolus indépendamment en impitant les conditions aux limites microscopiques du second ordre (4) et (5). Par exemple en 2D, on peut écrire en considérant la symétrie mineure de $\tilde{\mathbf{G}}(\tilde{\mathbf{G}}_{ijk} = \tilde{\mathbf{G}}_{kji})$ :

$$ \mathbf{u}_p = \mathbf{u}_p^l + \mathbf{u}_p^{nl} \text{ dans } \omega, \qquad (16) $$

$$ \begin{aligned} \mathbf{u}_p^l = & \bar{F}_p^{(11)} \tilde{\mathbf{u}}^{(110)} + \bar{F}_p^{(12)} \tilde{\mathbf{u}}^{(120)} + \bar{F}_p^{(21)} \tilde{\mathbf{u}}^{(210)} + \bar{F}_p^{(22)} \tilde{\mathbf{u}}^{(220)} + \bar{G}_p^{(111)} \tilde{\mathbf{u}}^{(111)} + 2 \bar{G}_p^{(112)} \tilde{\mathbf{u}}^{(112)} \\ & + \bar{G}_p^{(121)} \tilde{\mathbf{u}}^{(121)} + 2 \bar{G}_p^{(122)} \tilde{\mathbf{u}}^{(122)} + \bar{G}_p^{(212)} \tilde{\mathbf{u}}^{(212)} + \bar{G}_p^{(222)} \tilde{\mathbf{u}}^{(222)}, \end{aligned} \qquad (17) $$

où $\tilde{\mathbf{u}}^{(ijk)}$ ($i, j = 1 \cdots 2$ et $k = 0 \cdots 2$) sont les solutions des problèmes suivants :

$$ \left\{ \begin{array}{ll} \mathcal{L}(\tilde{\mathbf{u}}^{(ijk)}, \delta\mathbf{u}) = 0 & \text{dans } \omega \\ \tilde{\mathbf{u}}^{(ijk)+} - \tilde{\mathbf{u}}^{(ijk)-} = \mathbf{X}^{(ijk)+} - \mathbf{X}^{(ijk)-} & \text{sur } \partial\omega \\ \int_{\partial\omega^{-}} \tilde{\mathbf{u}}^{(ijk)-} d\Gamma = \int_{\partial\omega^{-}} \mathbf{X}^{(ijk)-} d\Gamma & \end{array} \right. \qquad (18) $$

avec

$$ \mathbf{X}^{(110)} = \begin{bmatrix} 1 & 0 \\ 0 & 0 \end{bmatrix} \mathbf{X}, \quad \mathbf{X}^{(120)} = \begin{bmatrix} 0 & 1 \\ 0 & 0 \end{bmatrix} \mathbf{X}, $$

$$ \mathbf{X}^{(210)} = \begin{bmatrix} 0 & 0 \\ 1 & 0 \end{bmatrix} \mathbf{X}, \quad \mathbf{X}^{(220)} = \begin{bmatrix} 0 & 0 \\ 0 & 1 \end{bmatrix} \mathbf{X}, $$

$$ \mathbf{X}^{(111)} = \frac{1}{2} \begin{bmatrix} 1 & 0 \\ 0 & 0 \end{bmatrix} \mathbf{X} \cdot \begin{bmatrix} 1 & 0 \\ 0 & 0 \end{bmatrix} \mathbf{X}, \quad \mathbf{X}^{(112)} = \frac{1}{2} \begin{bmatrix} 1 & 0 \\ 0 & 0 \end{bmatrix} \mathbf{X} \cdot \begin{bmatrix} 0 & 1 \\ 0 & 0 \end{bmatrix} \mathbf{X}, $$

$$ \mathbf{X}^{(121)} = \frac{1}{2} \begin{bmatrix} 0 & 0 \\ 1 & 0 \end{bmatrix} \mathbf{X} \cdot \begin{bmatrix} 0 & 0 \\ 1 & 0 \end{bmatrix} \mathbf{X}, \quad \mathbf{X}^{(122)} = \frac{1}{2} \begin{bmatrix} 0 & 0 \\ 1 & 0 \end{bmatrix} \mathbf{X} \cdot \begin{bmatrix} 0 & 0 \\ 0 & 1 \end{bmatrix} \mathbf{X}, $$

$$ \mathbf{X}^{(212)} = \frac{1}{2} \begin{bmatrix} 0 & 1 \\ 0 & 0 \end{bmatrix} \mathbf{X} \cdot \begin{bmatrix} 0 & 1 \\ 0 & 0 \end{bmatrix} \mathbf{X}, \quad \mathbf{X}^{(222)} = \frac{1}{2} \begin{bmatrix} 0 & 0 \\ 0 & 1 \end{bmatrix} \mathbf{X} \cdot \begin{bmatrix} 0 & 0 \\ 0 & 1 \end{bmatrix} \mathbf{X}, \qquad (19) $$

et $\mathbf{u}_p^{nl}$ est la solution du problème :

$$ \mathcal{L}(\mathbf{u}_p^{nl}, \delta\mathbf{u}) = \mathcal{F}_p^{nl}(\delta\mathbf{u}). \qquad (20) $$

Notons que la méthodologie proposée pour 2D peut être facilement étendue au cadre 3D. Enfin, la solution $\mathbf{u}_p$ peut s'écrire sous la forme suivante :

$$ \mathbf{u}_p = \mathbb{A} : \bar{\mathbf{F}}_p + \tilde{\mathbb{A}} : \bar{\mathbf{G}}_p + \mathbf{u}_p^{nl}, \qquad (21) $$

où le tenseur $\mathbb{A}$ représente un tenseur d'ordre trois donné par $A_{ijk} = \tilde{\mathbf{u}}_i^{(jk0)}$ et $\tilde{\mathbb{A}}$ est un tenseur d'ordre quatre défini par $\tilde{A}_{ijkl} = \tilde{\mathbf{u}}_i^{(jkl)}$ ($i, j, k, l = 1 \cdots 3$). En dérivant (21) par rapport à $\mathbf{X}$ et en utilisant $\nabla \mathbf{u}_p = \mathbf{F}_p$, on obtient :

$$ \mathbf{F}_p = \mathbb{A}_{,X} : \bar{\mathbf{F}}_p + \tilde{\mathbb{A}}_{,X} : \bar{\mathbf{G}}_p + \mathbf{u}_p^{nl}, \qquad (22) $$

où $\mathbb{A}_{,X}$ est un tenseur d'ordre quatre de localisation, $\tilde{\mathbb{A}}_{,X}$ est un tenseur d'ordre cinq qui est identifié comme un tenseur de localisation du second ordre, et $\mathbf{u}_p^{nl}$ est un tenseur d'ordre deux.

Pour obtenir la relation de comportement au niveau macroscopique à l'ordre $p$, on peut écrire la relation suivante :

$$ \mathbf{P}_p = \mathbb{H}^{(r)} : \mathbf{F}_p + \mathbf{P}_p^{nl}. \qquad (23) $$
En substituant l'équation (22) dans (23), la moyenne de l'équation (23) donne la contrainte macroscopique à chaque ordre $p$, $\bar{\mathbf{P}}_p$, sous la forme suivante :

$$ \bar{\mathbf{P}}_p = \bar{\mathbb{L}}_1 : \bar{\mathbf{F}}_p + \bar{\mathbb{L}}_2 : \bar{\mathbf{G}}_p + \bar{\mathbf{P}}_p^{nl}. \quad (24) $$

où

$$ \bar{\mathbb{L}}_1 = \frac{1}{|\omega|} \int_{\omega} (\mathbb{H}^{(r)} : \mathbb{A}, \mathbf{x}) d\omega, \qquad \bar{\mathbb{L}}_2 = \frac{1}{|\omega|} \int_{\omega} (\mathbb{H}^{(r)} : \tilde{\mathbb{A}}, \mathbf{x}) d\omega \\ \bar{\mathbf{P}}_p^{nl} = \frac{1}{|\omega|} \int_{\omega} (\mathbb{H}^{(r)} : \mathbf{u}_p^{nl}, \mathbf{P}_p^{nl}) d\omega. $$

Pour obtenir la contrainte macroscopique du second ordre, on développe d'abord l'équation de couplage (10). Donc, en utilisant la relation (24), on obtient à l'ordre $p$ :

$$ \bar{\mathbf{Q}}_p = \bar{\mathbb{L}}_3 : \bar{\mathbf{F}}_p + \bar{\mathbb{L}}_4 : \bar{\mathbf{G}}_p + \bar{\mathbf{Q}}_p^{nl}. \quad (25) $$

où :

$$ \bar{\mathbb{L}}_3 = \frac{1}{2|\omega|} \int_{\omega} [{}^t(\mathbb{H}^{(r)} : \mathbb{A}, \mathbf{x})\mathbf{X} + \mathbf{X}(\mathbb{H}^{(r)} : \mathbb{A}, \mathbf{x})] d\omega, \\ \bar{\mathbb{L}}_4 = \frac{1}{2|\omega|} \int_{\omega} [{}^t(\mathbb{H}^{(r)} : \tilde{\mathbb{A}}, \mathbf{x})\mathbf{X} + \mathbf{X}(\mathbb{H}^{(r)} : \tilde{\mathbb{A}}, \mathbf{x})] d\omega, \\ \bar{\mathbf{Q}}_p^{nl} = \frac{1}{2|\omega|} \int_{\omega} [{}^t(\mathbb{H}^{(r)} : \mathbf{u}_p^{nl}, \mathbf{P}_p^{nl})\mathbf{X} + \mathbf{X}(\mathbb{H}^{(r)} : \mathbf{u}_p^{nl}, \mathbf{P}_p^{nl})] d\omega. $$

## 3.2 Problème macroscopique à résoudre

Pour utiliser directement les éléments finis classiques basés sur le déplacement au niveau macroscopique (équation (7)), on a besoin de fonctions d'interpolation au moins continues de classe $C^1$ ; c'est-à-dire, la fonction et sa dérivée première sont continues à travers les éléments. Toutefois, le développement des éléments continus $C^1$ adaptés pour les problèmes non linéaires multidimensionnels n'est pas simple car ces éléments nécessitent généralement plus d'hypothèses (par exemple géométrique). Pour ces raisons, les éléments continus de classe $C^0$ sont très utilisés dans la littérature. Le développement de ce type d'éléments pour un problème du second ordre est basé généralement sur une formulation mixte. La stratégie consiste à introduire un autre champ tensoriel inconnu qui présente un gradient de déformation, noté $\hat{\mathbf{F}}$, en plus du champ de déplacement inconnu $\bar{\mathbf{u}}$. La condition cinématique entre ces deux champs est imposée au sens faible par les multiplicateurs de Lagrange, $\mu$ :

$$ \int_{\Omega} \{ {}^t \bar{\mathbf{P}} : \delta \bar{\mathbf{F}} + \bar{\mathbf{Q}} : \delta \hat{\mathbf{G}} + \delta \{ {}^t \mu : (\hat{\mathbf{F}} - \bar{\mathbf{F}}) \} d\Omega = \lambda \int_{\partial \Omega_t} (\mathbf{f} \cdot \delta \bar{\mathbf{u}} + {}^t \mathbf{r} : \delta \hat{\mathbf{F}}) d\Gamma. \quad (26) $$

On décompose l'équation (26) en trois parties par rapport à $\delta \bar{\mathbf{u}}$, $\delta \hat{\mathbf{F}}$ et $\delta \mu$ :

$$ \left\{ \begin{array}{l} \int_{\Omega} \{ {}^t \bar{\mathbf{P}} : \delta \bar{\mathbf{F}} - {}^t \mu : \delta \bar{\mathbf{F}} \} d\Omega = \lambda \int_{\partial \Omega_t} \mathbf{f} \cdot \delta \bar{\mathbf{u}} d\Gamma \\ \int_{\Omega} \{ \bar{\mathbf{Q}} : \delta \hat{\mathbf{G}} + {}^t \mu : \delta \hat{\mathbf{F}} \} d\Omega = \lambda \int_{\partial \Omega_t} {}^t \mathbf{r} : \delta \hat{\mathbf{F}} d\Gamma \\ \int_{\Omega} \{ {}^t \delta \mu : (\hat{\mathbf{F}} - \bar{\mathbf{F}}) \} d\Omega = 0 \end{array} \right. \quad (27) $$

Au niveau macroscopique, l'introduction du développement asymptotique dans le problème non linéaire (27) permet d'obtenir une suite de problèmes linéaires qui s'écrivent à l'ordre $p$ comme suit :

$$ \left\{ \begin{array}{l} \int_{\Omega} \{ {}^t \bar{\mathbf{P}}_p : \delta \bar{\mathbf{F}} - {}^t \mu_p : \delta \bar{\mathbf{F}} \} d\Omega = \lambda_p \int_{\partial \Omega_t} \mathbf{f} \cdot \delta \bar{\mathbf{u}} d\Gamma \\ \int_{\Omega} \{ \bar{\mathbf{Q}}_p : \delta \hat{\mathbf{G}} + {}^t \mu_p : \delta \hat{\mathbf{F}} \} d\Omega = \lambda_p \int_{\partial \Omega_t} {}^t \mathbf{r} : \delta \hat{\mathbf{F}} d\Gamma \\ \int_{\Omega} \{ {}^t \delta \mu : (\hat{\mathbf{F}}_p - \bar{\mathbf{F}}_p) \} d\Omega = 0 \end{array} \right. \quad (28) $$
En substituant les relations de comportement (24) et (25) dans les équations (28), la forme finale du problème linéaire macroscopique du second ordre peut s'écrire à l'ordre p sous la forme suivante :

$$
\left\{
\begin{array}{l}
\int_{\Omega} \left\{ \left[ \mathbf{t} \bar{\mathbf{F}}_p : \bar{\mathbb{L}}_1 + \hat{\mathbf{G}}_p : \bar{\mathbb{L}}_2 \right] : \delta \bar{\mathbf{F}} - {}^t \mu_p : \delta \bar{\mathbf{F}} \right\} d\Omega = \lambda_p \int_{\partial \Omega_t} \mathbf{f} \cdot \delta \bar{\mathbf{u}} d\Gamma - \int_{\Omega} {}^t \bar{\mathbf{P}}_p^{nl} : \delta \bar{\mathbf{F}} d\Omega \\[2ex]
\int_{\Omega} \left\{ \left[ \mathbf{t} \bar{\mathbf{F}}_p : \bar{\mathbb{L}}_3 + \hat{\mathbf{G}}_p : \bar{\mathbb{L}}_4 \right] : \delta \hat{\mathbf{G}} + {}^t \mu_p : \delta \hat{\mathbf{F}} \right\} d\Omega = \lambda_p \int_{\partial \Omega_F} {}^t \mathbf{r} : \delta \hat{\mathbf{F}} d\Gamma - \int_{\Omega} \bar{\mathbf{Q}}_p^{nl} : \delta \hat{\mathbf{G}} d\Omega \\[2ex]
\int_{\Omega} \left\{ {}^t \delta \mu : (\hat{\mathbf{F}}_p - \bar{\mathbf{F}}_p) \right\} d\Omega = 0
\end{array}
\right.
\qquad (29)
$$

La résolution par éléments finis du problème (29) donne les déplacements ū ainsi que les déformations macroscopiques ˆF et les multiplicateurs de Lagrange µ au niveau de la structure macroscopique. Ensuite, le calcul des tenseurs ˆF et ˆG en chaque point d'intégration permet de finir les calculs des différentes variables au niveau microscopique.

4 Exemple numérique

Pour évaluer la précision et l'efficacité des approches proposées, nous avons étudié la compression
d'un matériau composite à fibres longues. La structure macroscopique et la microstructure sont présen-
tées par les figures 2a et 2b. La microstructure est maillée avec des éléments quadrangulaires à huit
nœuds. En revanche, la structure macroscopique est maillée avec des éléments quadrangulaires à huit
nœuds pour l'homogénéisation du premier ordre et avec des éléments mixtes QU8F4L1 pour l'homogé-
néisation du second ordre [3] (voir figure 2c). Grâce à la symétrie, seule la moitié de la structure a été
utilisée pour résoudre ce problème. Le problème proposé a été discrétisé à l'aide d'éléments finis à deux
dimensions sous l'hypothèse des contraintes planes.

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

Fig. 2 – (a) Structure macroscopique du problème de compression d’un matériau composite à fibres longues, P = 1 N/mm ; (b) microstructure associée au matériau composite à fibres longues (l=2mm et l/h=20) ; (c) élément QU8F4L1 pour discrétisation de la macrostructure dans le cadre de l’homogénéisation du second ordre (• : déplacement **u**, □ : gradient de déformation **F**, × : multiplicateur de Lagrange **μ**).

Différents maillages macroscopiques sont considérés afin de montrer l'influence du maillage sur les résultats des deux approches. La réponse charge-déplacement des deux techniques pour différents maillages est reportée sur la figure 3 (le maillage *m* × *n* désigne un maillage avec *m* éléments dans la longueur et *n* éléments dans la largeur). Cela montre une convergence assez rapide de la charge critique de flambage pour l'homogénéisation du second ordre contrairement de l'homogénéisation du premier ordre. Cette convergence confirme que la solution obtenue à l'aide de l'homogénéisation du second ordre ne
dépend pas du maillage macroscopique lorsque la taille de la microstructure n'est pas négligeable par rapport à la longueur d'onde du champ de déformation macroscopique (maillages 6 × 6, 6 × 9 et 6 × 12). Dans le cas contraire (maillage 6 × 3), les résultats obtenus par la modélisation du second ordre coïncident approximativement avec ceux de l'analyse du premier ordre.

<figure><img src="image_4.png" /><figcaption>Fig. 3 – Le diagramme charge-déplacement pour les différents maillages macroscopiques pour les homogénéisations du premier ordre et du second ordre.</figcaption></figure>

## 5 Conclusion

Dans cette étude, nous avons présenté une procédure MAN-multi-échelle dans le cadre de l'homogénéisation de second ordre. Cette technique permet un meilleur conditionnement du problème macroscopique en présence d'instabilités au niveau de la microstructure. Le problème macroscopique est classique et reste inchangé par rapport à la technique MAN-multi-échelle appliquée dans le cadre de l'homogénéisation du premier ordre. En revanche, le problème macroscopique fait intervenir des variables de premier et second ordre. Un élément fini adapté est donc utilisé pour la discrétisation de ce problème. Nous montrons au travers un exemple numérique d'instabilité que l'homogénéisation de second ordre permet d'obtenir la solution du problème sans dépendance par rapport au maillage ce qui n'est pas le cas lorsque la théorie du premier ordre est utilisée.

## Références

[1] S. Nezamabadi, J. Yvonnet, H. Zahrouni, M. Potier-Ferry. A multilevel computational strategy for microscopic and macroscopic instabilities.Computer Methods in Applied Mechanics and Engineering, 198, 2099–2110, 2009.

[2] S. Nezamabadi, H. Zahrouni, J. Yvonnet, M. Potier-Ferry. A multiscale finite element approach for buckling analysis of elastoplastic long fiber composites. International Journal for Multiscale Computational Engineering, 8, 287–301, 2010.

[3] V.G. Kouznetsova, M.G.D. Geers and W.A.M. Brekelmans, Multi-scale constitutive modelling of heterogeneous materials with a gradient-enhanced computational homogenization scheme. International Journal for Numerical Methods in Engineering, 54, 1235–1260, 2002.

[4] V.G. Kouznetsova, M.G.D. Geers and W.A.M. Brekelmans, Multi-scale second-order computational homogenization of multi-phase materials: a nested finite element solution strategy. Computer Methods in Applied Mechanics and Engineering, 193, 5525–5550, 2004.

[5] S. Nezamabadi, H. Zahrouni, J. Yvonnet. Solving hyperelastic material problems by Asymptotic Numerical Method, Computational Mechanics, Springer, 47: 77-92, 2011.