# Identification sans loi de la réponse mécanique de matériaux à partir de mesures de champs hétérogènes via une approche " Data Driven Computational Mechanics "

Adrien Leygue, Michel Coret, Marie Dalemat, Julien Réthoré, Laurent Stainier, Erwan Verron

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

## To cite this version:

Adrien Leygue, Michel Coret, Marie Dalemat, Julien Réthoré, Laurent Stainier, et al. Identification sans loi de la réponse mécanique de matériaux à partir de mesures de champs hétérogènes via une approche "Data Driven Computational Mechanics". CFM 2017 - 23ème Congrès Français de Mécanique, Aug 2017, Lille, France. hal-03465375

HAL Id: hal-03465375

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

Submitted on 3 Dec 2021

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published 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
# Identification sans loi de la réponse mécanique de matériaux à partir de mesures de champs hétérogènes via une approche « Data Driven Computational Mechanics »

Adrien Leygueª, Michel Coretª, Marie Dalematª, Julien Réthoréᵇ, Laurent Stainierᵇ et Erwan Verronª

a. Ecole Centrale de Nantes, GeM/PMM, [prénom].[nom]@ec-nantes.fr

b. Ecole Centrale de Nantes, GeM/MS, [prénom].[nom]@ec-nantes.fr

## Résumé :

Kirchdoerfer & Ortiz [1] ont récemment proposé un nouveau paradigme pour la simulation numérique de structures, dans lequel la loi de comportement est remplacée par une base de données constituée d'un ensemble fini de couples contrainte/déformation admissibles pour le matériau. Cette approche originale, baptisée « Data Driven Computational Mechanics » repose toutefois sur la disponibilité de telles bases de données matériau pour une large gamme de déformations et de contraintes.

Dans ce travail, nous proposons une méthode, basée sur l'approche de Kirchdoerfer et Ortiz qui permet de construire de telles bases de données pour des matériaux élastiques, sans postuler de loi de comportement. À partir d'un ensemble de mesures de champs [2, 3] hétérogènes de déformation pour différents cas de chargement appliqués à une éprouvette, on identifie simultanément la base de données et un champ de contraintes équilibrées au sein de l'éprouvette pour chaque cas de chargement.

## Abstract :

Kirchdoerfer & Ortiz [1] recently introduced a new paradigm in computational mechanics, where the constitutive law is replaced by a database of discrete strain-stress couples representing the mechanical behavior of the material. This original and promising approach nevertheless requires the availability of such a database for a wide range of strain and stress.

In this work we propose a method, based on the work of Kirchdoerfer and Ortiz which produces such databases for elastic materials, without any assumption on an underlying constitutive law. From a set of heterogeneous DIC displacement fields [2, 3] for different loading conditions of a sample, we simultaneously identify the database and a statically admissible stress field in the sample for each loading condition.

**Mots clefs :** Data-Driven Computational Mechanics, corrélation d'image, loi de comportement
# 1 Introduction

Le concept de loi de comportement est central en ingénierie mécanique et en calcul de structure car il permet de relier directement, dans le cas d'un comportement élastique, les déformations et les contraintes au sein d'un matériau. Une loi de comportement comporte en général un jeu de paramètres qui sont ajustés à partir d'un ensemble de données expérimentales [3]. Au delà de la simple description de la réponse mécanique, la loi en elle même a plusieurs rôles :

- Elle fourni une relation déformation contrainte lissée et ainsi débarrassée d'un certain bruit expérimental.

- Pour un type donné de chargement (e.g. traction uniaxiale) elle permet d'interpoler et d'extrapoler les contraintes parmi les quelques mesures effectuées.

- Sa forme tensorielle permet d'étendre la relation déformation-contrainte à des cas de chargement multiaxiaux possiblement difficile à produire et à instrumenter expérimentalement.

Récemment, Kirchdoerfer & Ortiz [1] ont introduit le concept de "Data-Driven Computational Mechanics" (DDCM) pour le calcul de structure impliquant des matériaux élastiques. Dans cette approche, le concept de loi de comportement s'efface au profit d'une base de données de couples déformation-contrainte qui échantillonnent la réponse mécanique du matériau. Ainsi, le lissage et l'interpolation des données expérimentales traditionnellement réalisés par la loi de comportement sont effectués au sein du calcul de structure. Les premiers résultats de la méthode sont encourageants et on entrevoit de nombreuses perspectives du point de vue de la modélisation et de la simulation.

Dans les cas particulier des treillis de barres, l'algorithme de Kirchdoerfer & Ortiz cherche à attribuer à chaque barre $e$ de la structure deux états, ou couples déformation/contrainte. Le premier état ($\epsilon_e, \sigma_e$), que nous appelons état mécanique satisfait, à l'échelle de la structure, à la fois la compatibilité des déformations $\epsilon_e$ et l'équilibre mécanique des contraintes $\sigma_e$ :

$$ \epsilon_e = \sum_j \mathbf{B}_{ej} \mathbf{u}_j, \qquad (1) $$

et :

$$ \sum_e w_e \mathbf{B}_{ej} \sigma_e = f_j. \qquad (2) $$

Dans les équations ci-dessus, $\mathbf{u}_j$ représente les valeurs nodales de déplacement, $\mathbf{f}_j$ représente les forces nodales, $w_e$ est le volume de la barre $e$ et $\mathbf{B}_{ej}$ traduit la géométrie et la connectivité du treillis. Le second état associé à $e$ et noté $(\epsilon_{ie}^*, \sigma_{ie}^*)$ est appelé état matériau. Cet état est choisi au sein d'une base de données d'états matériau admissibles pour le matériau considéré : $(\epsilon_i^*, \sigma_i^*)$, avec $i \in 1:N^*$. L'indice $ie \in 1:N^*$ donne l'état matériau particulier de l'élément $e$. L'algorithme proposé cherche, pour chaque élément un état mécanique et un état matériau aussi similaires que possible, mais tels que le premier satisfait la compatibilité des déformations (1) et l'équilibre des contraintes (2). On obtient le problème de minimisation suivant :

$$ \text{solution} = \arg \min_{\epsilon_e, \sigma_e, ie} \frac{1}{2} \sum_e w_e \left\| \left( \epsilon_e - \epsilon_{ie}^*, \sigma_e - \sigma_{ie}^* \right) \right\|_C^2, \qquad (3) $$

sous les contraintes (1) et (2). Dans l'équation ci-dessus, $\|\| (\epsilon, \sigma) \|\|_C$ est une norme énergétique caractérisant le désaccord entre l'états mécanique et matériau. La figure 1 illustre le concept de double état associé à chaque élément du treillis.
<figure><img src="image_3.png" /><figcaption>FIGURE 1 – Illustration des deux états $(\epsilon_e, \sigma_e)$ et $(\epsilon_{ie}^*, \sigma_{ie}^*)$ associés à un élément barre $e$. La ligne interrompue représente la disparité entre les deux états. Sur la gauche, on observe que le mapping $ie$ entre éléments et états matériau peut associer le même état matériau à deux éléments différents.</figcaption></figure>

Pour la norme :

$$ ||(\epsilon_e, \sigma_e)||_C^2 = \frac{1}{2}(C_e\epsilon_e^2 + C_e^{-1}\sigma_e^2), \quad (4) $$

les auteurs proposent un algorithme efficace pour la résolution de ce problème de complexité combinatoire. La constante $C_e$ est le seul paramètre de la méthode et peut-être interprété comme une raideur virtuelle associée à la disparité entre états mécaniques et matériaux.

Cette approche, sans loi de comportement, nécessite toutefois de disposer d'une base de données d'états matériau. Pour des approches purement numériques la construction de ces données, par exemple via une approche FEM² peut être numériquement couteuse. D'un point de vue expérimental il n'est pas trivial de générer et de mesurer déformations et contraintes sur une large gamme afin de peupler la base de données.

Dans ce travail on propose une méthode basée sur l'approche DDCM pour extraire des états matériaux à partir d'une collection de mesures de champ non homogènes correspondant à différents cas de chargement d'une même éprouvette. Ces mesures sont typiquement obtenues à l'aide de méthodes de corrélation d'image[2]. Ces travaux préliminaires ne considèrent que le cas des petites déformations et d'une réponse purement élastique mais non-linéaire.

## 2 Identification "Data-driven" pour un comportement élastique

Dans cette section, on développe une méthode d'identification du comportement élastique d'un matériau basée sur l'approche DDCM présentée dans la section précédente. Dans un soucis de simplicité on se limite à des déformations infinitésimales, même si l'extension à des déformations finies ne pose pas de problèmes fondamentaux. On défini, un état (mécanique ou matériau) comme un couple constitué d'un tenseur de déformation $\epsilon$ et d'un tenseur des contraintes de Cauchy $\sigma$, chacun appartenant à un espace à six dimensions (après prise en compte des symétries).

### 2.1 Procédure d'identification

On considère un ensemble de mesures de champ de déplacement réalisées, par exemple par corrélation d'image, sur une éprouvette pour différents cas de chargement. Ces mesures sont fournies sur un maillage de type Elements Finis pour lequel chaque point de quadrature $e$ a un poids d'intégration associé $w_e$. Pour chaque cas de chargement ou "snapshot" indexé par $X$ on suppose qu'on a accès aux quantités suivantes :

— Les déplacements nodaux $u_j^X$,
— La géométrie et la connectivité du maillage, encodée à travers les matrices $\mathbf{B}_{ej}^X$, permettant de
calculer les déformations mécaniques associées aux observations : $\epsilon_e^X = \sum_j \mathbf{B}_{ej}^X \cdot \mathbf{u}_j^X$

— Les forces nodales $\mathbf{f}_j^X$ correspondant au chargement,

— Les conditions limites de Dirichlet.

On cherche à calculer un ensemble $N^*$ d'états matériau $(\epsilon_i^*, \sigma_i^*)$ tels que :

— pour chaque "snapshot" $X$ et chaque point de quadrature $e$ on peut calculer une contrainte mécanique $\boldsymbol{\sigma}_e^X$ qui satisfait l'équilibre mécanique du "snapshot" $X$,

— pour chaque "snapshot", on peut assigner un état matériau $(\epsilon_{ie}^*, \sigma_{ie}^*)$ à chaque point de quadrature $e$. Cet état matériau doit être aussi proche que possible de l'état mécanique au sens d'une certaine norme $|| \cdot ||_{\mathcal{C}_e^2}$.

En s'inspirant de Kirchdoerfer et Ortiz, on choisi :

$$ ||(\epsilon_e, \sigma_e)||_{\mathcal{C}_e}^2 = \frac{1}{2}(\epsilon_e : \mathcal{C}_e : \epsilon_e + \sigma_e : \mathcal{C}_e^{-1} : \sigma_e), \quad (5) $$

où $\mathcal{C}_e$ est un tenseur du quatrième ordre symétrique et défini positif, assimilable à un tenseur d'élasticité.
Le problème de minimisation décrit ci-dessus se traduit comme :

$$ \text{solution} = \arg \min_{\boldsymbol{\sigma}_e^X, \epsilon_i^*, \sigma_i^*, ie^X} \mathcal{E}(\boldsymbol{\sigma}_e^X, \epsilon_i^*, \sigma_i^*, ie), \quad (6) $$

avec

$$ \mathcal{E}(\boldsymbol{\sigma}_e^X, \epsilon_i^*, \sigma_i^*, ie) = \sum_X \sum_e w_e^X ||(\epsilon_e^X - \epsilon_{ie}^*{}^X, \boldsymbol{\sigma}_e^X - \boldsymbol{\sigma}_{ie}^*{}^X)||_{\mathcal{C}_e}^2, \quad (7) $$

et les contraintes d'équilibre mécanique de chaque "snapshot" :

$$ \sum_e w_e^X \mathbf{B}_{ej}^{X T} \cdot \boldsymbol{\sigma}_e^X = \mathbf{f}_j^X \quad \forall X, j. \quad (8) $$

Dans le problème ci-dessus, toutes les inconnues sont continues, à l'exception du mapping $ie^X$ qui est discret. En supposant que $ie^X$ est donné, on impose l'équilibre mécanique à l'aide de multiplicateurs de Lagrange $\eta_j^X$. Les variables continues $(\boldsymbol{\sigma}_e^X, \epsilon_i^*, \sigma_i^*)$ sont alors solutions du problème suivant :

$$ \delta \left( \sum_X \sum_e \left( w_e^X (\epsilon_e^X - \epsilon_{ie}^*{}^X) : \mathcal{C}_e^X : (\epsilon_e^X - \epsilon_{ie}^*{}^X) + w_e^X (\boldsymbol{\sigma}_e^X - \boldsymbol{\sigma}_{ie}^*{}^X) : \mathcal{C}_e^{X-1} : (\boldsymbol{\sigma}_e^X - \boldsymbol{\sigma}_{ie}^*{}^X) \right. \right. \\ \left. \left. - \sum_j (w_e^X \mathbf{B}_{ej}^{X T} \cdot \boldsymbol{\sigma}_e^X - \mathbf{f}_j^X) \cdot \boldsymbol{\eta}_j^X \right) \right) = 0 \quad (9) $$

En considérant toutes les variations possibles, on obtient les équations suivantes :

$$ \delta \epsilon_i^* \Rightarrow \sum_X \sum_{ie=i} w_e^X \mathcal{C}_e^X : (\epsilon_e^X - \epsilon_{ie}^*{}^X) = 0 \quad \forall i \quad (10) $$

$$ \delta \boldsymbol{\sigma}_i^* \Rightarrow \sum_X \sum_{ie=i} w_e^X \mathcal{C}_e^{X-1} : (\boldsymbol{\sigma}_e^X - \boldsymbol{\sigma}_{ie}^*{}^X) = 0 \quad \forall i \quad (11) $$

$$ \delta \boldsymbol{\sigma}_e^X \Rightarrow w_e^X \mathcal{C}_e^{X-1} : (\boldsymbol{\sigma}_e^X - \boldsymbol{\sigma}_{ie}^*{}^X) - \sum_j w_e^X \mathbf{B}_{ej}^X \cdot \boldsymbol{\eta}_j^X = 0 \quad \forall e, X \quad (12) $$

$$ \delta \boldsymbol{\eta}_j^X \Rightarrow \sum_e (w_e^X \mathbf{B}_{ej}^{X T} \cdot \boldsymbol{\sigma}_e^X - \mathbf{f}_j^X) = 0 \quad \forall j, X \quad (13) $$
En combinant les équations (11), (12) et (13) on obtient le système suivant pour déterminer simultanément $\sigma_i^*$ et $\sigma_e^X$ (à travers $\eta_j^X$ et (12)) :

$$ \sum_k \sum_e w_e^X \mathbf{B}_{ej}^{X,T} : \mathbb{C}_e^X : \mathbf{B}_{ek}^X \eta_k^X + \sum_e w_e^X \mathbf{B}_{ej}^{X,T} \sigma_{ieX}^* = \mathbf{f}_j^X \quad \forall j, X, \quad (14) $$

$$ \sum_{ie=i} \sum_{X} \sum_{j} w_e^X \mathbf{B}_{ej}^X \eta_j^X = 0 \quad \forall i. \quad (15) $$

Ce système linéaire peut s'interpréter comme la recherche des champs de déplacement virtuels $\eta_j^X$ satisfaisant, pour la raideur $\mathbb{C}_e$ et après prise en compte des contraintes matériau, l'équilibre mécanique de tous les "snapshots" ; sous la contrainte que la moyenne des déformations virtuelles associées à chaque état matériau est nulle.

Dans le cas où le mapping ie, n'est pas connu, on propose l'algorithme de résolution suivant pour le calcul de $\epsilon_i^*, \sigma_i^*, \sigma_e^X$ et ie :

1. Initialisation de $\epsilon_i^*$ et ie par un algorithme de type kmeans [4, 5] appliqué à $\epsilon_e^X$.

2. Calcul de $\sigma_i^*$ et $\eta_j^X$ via (14) et (15).

3. Mise à jour des valeurs de $\sigma_e^X$ à l'aide de l'équation (12).

4. Calcul d'un nouveau mapping ie$^X$ :

$$ ie^X = \arg \min_{ie^X} \sum_X \sum_e w_e^X \left\| \left( \epsilon_e^X - \epsilon_{ieX}^*, \sigma_e^X - \sigma_{ieX}^* \right) \right\|_{\mathbb{C}_e}^2. \quad (16) $$

5. Mise à jour de $\epsilon_i^*$ via (10).

6. Itération des étapes 2,3,4 et 5 jusqu'à convergence de ie$^X$.

## 2.2 Résultats numériques

Pour des données manufacturées, on applique la méthode décrite ci-dessus à l'identification de la réponse mécanique d'un matériau élastique non linéaire incompressible, dans l'hypothèse des contraintes planes. On considère un maillage EF 2D de 1340 noeuds et 2416 éléments triangulaires, soumis à différents cas de chargement ($N^X = 40$). Le maillage et la déformée de deux cas de chargement est représentée en figure 2. La loi de comportement utilisée pour générer les données est de la forme :

$$ \boldsymbol{\sigma} = G(\boldsymbol{\epsilon} + \alpha\boldsymbol{\epsilon}^3) - p\mathbf{I}, \quad (17) $$

$$ p = -(\epsilon_{xx} + \epsilon_{yy}) - \alpha(\epsilon_{xx} + \epsilon_{yy})^3. \quad (18) $$

Dans la procédure d'identification, on a utilisé un tenseur $\mathbb{C}_e$ correspondant à l'équation ci-dessus avec $\alpha = 0$ et $G = 1$. Dans ce cas, l'initialisation de ie$^X$ fournie par l'algorithme kmeans est si efficace que quelques itérations seulement sont nécessaires à la convergence pour $N^* = 500$.

Afin de visualiser les résultats, on montre dans le figure 3 la première valeur propre des contraintes matériau ($\sigma_I^*$) en fonction des valeurs propres de la déformation matériau correspondante ($\epsilon_I^*, \epsilon_{II}^*$). On remarque que tous les points sont sur la surface qu'on peut construire à partir de la loi de comportement utilisée pour générer les données :

$$ \frac{\sigma_I}{G} = \epsilon_I + \alpha * \epsilon_I^3 + (\epsilon_I + \epsilon_{II}) + \alpha * (\epsilon_I + \epsilon_{II})^3. \quad (19) $$
<figure><img src="image_4.png" /><figcaption>FIGURE 2 – Maillage 2D et déformées correspondant à deux cas de chargement. La couleur traduit l'amplitude du déplacement.</figcaption></figure>

On a ainsi construit un ensemble d'états matériau qui fourni une bonne description discrète de la "vraie" réponse du matériau, pour la gamme de déformations présente dans les données.

<figure><img src="image_5.png" /><figcaption>FIGURE 3 – Première valeur propre du tenseur de contraintes en fonction des valeurs propres du tenseur de déformation. Les symboles sont calculés à partir des états matériau $(\epsilon^*, \sigma^*)$ et la surface est celle décrite par l'équation (17).</figcaption></figure>

Outre les états matériau la méthode proposée calcule aussi, pour chaque "snapshot", un champ de contraintes mécaniques équilibrées et en accord avec le comportement identifié. Dans la figure 4 on montre la contrainte de Von Mises prédite pour un des "snapshots" ainsi que l'erreur relative. La prédiction a une erreur essentiellement inférieure à 10% pour la majorité des éléments et atteint 50% dans quelques uns des 2416 éléments, confirmant ainsi la précision de la méthode.
<figure><img src="image_6.png" /><figcaption>FIGURE 4 – Gauche : Contrainte de Von Mises prédite pour un "snapshot" particulier. Droite : erreur relative sur la contrainte de Von Mises par rapport au modèle utilisé pour générer les données.</figcaption></figure>

## 3 Résumé et conclusion

On a proposé une méthode inspirée du formalisme DDCM qui, pour un matériau élastique, permet de calculer un ensemble de couples déformations contraintes décrivant la réponse mécanique du matériau. La méthode d'identification requiert uniquement un ensemble de mesures de champ de déformation, par exemple fourni par corrélation d'image, et les forces appliquées pour chaque cas de chargement. Ces états matériau peuvent ensuite être utilisés directement en remplacement d'une loi de comportement dans l'approche DDCM, ou bien servir de données pour l'ajustement d'un modèle classique de comportement. La méthode permet en outre de calculer, pour chaque donnée, un champ de contraintes équilibrées satisfaisant le comportement identifié au sens de la DDCM. Cette approche est pour le moment uniquement applicable à des comportements élastiques non linéaires, pour lesquels la contrainte est une fonction de la déformation uniquement. Des extensions à des classes de comportement plus complexes comme la viscoélasticité et l'élastoplasticité sont en cours de développement.

## Références

[1] Kirchdoerfer, T., & Ortiz, M. (2016). Data-driven computational mechanics. Comput. Methods Appl. Mech. Engrg., 304, 81–101. http://doi.org/10.1016/j.cma.2016.02.001

[2] Sutton, M. A., Orteu, J. J., & Schreier, H. (2009). Image Correlation for Shape, Motion and Deformation Measurements. Springer Science & Business Media.

[3] Réthoré, J., Muhibulah, Elguedj, T., Coret, M., Chaudet, P., & Combescure, A. (2013). Robust identification of elasto-plastic constitutive law parameters from digital images using 3D kinematics. International Journal of Solids and Structures, 50(1), 73-85.

[4] MacQueen, J. (1967). Some methods for classification and analysis of multivariate observations (Vol. 1, pp. 281–297). Presented at the Proceedings of the fifth Berkeley symposium on Mathematical Statistics and Probability.

[5] Lloyd, S. (1982). Least squares quantization in PCM. IEEE Transactions on Information Theory, 28(2), 129–137. http://doi.org/10.1109/tit.1982.1056489