# Une nouvelle procédure d'identification des paramètres de lois cohésives

Benoît Blaysat, Johan P M Hoefnagels, Gilles Lubineau, Marc Geers

► To cite this version:

Benoît Blaysat, Johan P M Hoefnagels, Gilles Lubineau, Marc Geers. Une nouvelle procédure d'identification des paramètres de lois cohésives. 11e colloque national en calcul des structures, CSMA, May 2013, Giens, France. hal-01717061

HAL Id: hal-01717061

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

Submitted on 28 Feb 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
# Une nouvelle procédure d'identification des paramètres de lois cohésives

Benoît BLAYSAT<fnref n="1" />, Johan P.M. HOEFNAGELS<fnref n="1" />, Gilles LUBINEAU<fnref n="2" />, Marc G.D. GEERS<fnref n="1" />

**Résumé** — L'objectif est ici d'introduire un nouvel outil de caractérisation de lois cohésives, modèle le plus couramment utilisé pour décrire et simuler le phénomène de délaminage. Cet outil est basé sur une méthode de corrélation d'image globale, où la définition de l'espace de recherche cinématique est adéquatement choisie. Après avoir présenté le principe de cette méthode, on s'intéresse à sa validation, ainsi qu'à son comportement au bruit. Cette validation est réalisée numériquement, (i) en construisant un ensemble d'images à l'aide d'un outil de simulation évolué du délaminage et (ii) en identifiant à partir de ces images les caractéristiques du comportement de l'interface.

**Mots clés** — G-DIC, Zone Cohésive, Délaminage, Loi cohésive

## 1 Introduction

Le phénomène de délaminage, qui s'initie aux interfaces à une échelle micro, entraîne la plupart du temps la rupture, macro, du système. Au sein des composites, matériaux hétérogènes par définition et présentant une haute densité d'interfaces entre différents matériaux, ce mécanisme de rupture est de loin le plus courant. Une illustration des conséquences nuisibles du délaminage sur une puce microélectronique est présentée figure 1. Ces puces, constituées d'un circuit de cuivre enveloppé par une résine d'époxy, présentent bien les caractéristiques d'un matériau composite. Lors de leur utilisation, le courant traversant le circuit de cuivre a tendance à chauffer celui-ci, soumettant alors l'interface à des contraintes singulières qui initient le délaminage.

<figure><img src="image_2.png" /><figcaption>Fig. 1 — (a,b) images de microscope à balayage électronique (MEB) montrant (cf flèches) le délaminage d'une interface entre du cuivre et un diélectrique (Courtesy of NXP Semiconductors). (c) Schéma structurel d'un ensemble électronique noyé dans une résine époxy. Toutes les interfaces où du délaminage est susceptible de s'initier et de se propager sont représentées en pointillés</figcaption></figure>

Afin d'optimiser le design de ces composites et donc de minimiser le risque d'initiation mais aussi de propagation du délaminage interfacial, les concepteurs ont de plus en plus recours à des simulations éléments finis, où l'interface est souvent représentée à l'aide d'une loi cohésive. Introduit initialement dans [1, 2], le modèle de zone cohésive permet de représenter le comportement d'interface à l'aide d'une loi d'ouverture, liant l'effort d'interface à son ouverture, comme illustré figure 2. Dans [3, 4], les auteurs ont introduit les éléments cohésifs, qui permettent d'implémenter un modèle cohésif dans un problème

<footnote n="1">¹ Eindhoven University of Technology (TU/e), Department of Mechanical Engineering,
Den Dolech 2, 5612 AZ, Eindhoven, the Netherlands, b.blaysat/j.p.m.hoefnagels/m.g.d.geers@tue.nl</footnote>

<footnote n="2">² King Abdullah University of Science and Technology (KAUST), Physical Science and Engineering Division, Cohmas Laboratory,
Thuwal 23955-6900, Saudi Arabia, gilles.lubineau@kaust.edu.sa</footnote>
éléments finis. Aujourd'hui, l'utilisation conjointe des éléments finis et des éléments cohésifs est de loin l'approche prédictive du délaminage la plus commune.

Enfin, de sorte à enrichir au mieux le modèle continu sous-jacent au modèle éléments finis, il est nécessaire d'avoir une vision fidèle et précise du comportement d'une interface. Nous proposons, dans ce papier, une méthode qui permet de caractériser un modèle de cohésif d'une fissure interfacial.

<figure><img src="image_3.png" /><figcaption>Fig. 2 – Détails du modèle de zone cohésive. (a) description de l'interface à l'aide d'un effort *t* et de l'ouverture associée $\Delta$, (b) évolution des quantités d'interface en fonction de la coordonnée d'interface et (c) loi d'ouverture</figcaption></figure>

Notre outil repose sur trois points clés : (i) le MMMB, qui permet d'extraire des images de la propagation d'une fissure interfaciale à des échelles inédites, (ii) une description de la cinématique lors d'un essai DCB, rendant possible le lien fort entre champs de déplacement et loi d'ouverture de l'interface, et enfin, (iii) un algorithme de corrélation d'image permettant d'identifier les paramètres matériaux sous-jacent à un ensemble d'images. Ce papier s'intéresse aux seules étapes numériques, i.e. (ii) et (iii).

# 2 Principe

## 2.1 Essai

De nombreuses techniques expérimentales, dédiées à l'étude des propriétés d'interfaces, ont été développées (par exemple, cf [5, 6, 7, 8, 9]). La plupart de ces approches cherchent à identifier des grandeurs issues de la mécanique de la fracture, la grande polyvalence de cette description ayant déjà largement été prouvée. Un essai, dédié à la caractérisation de ces quantités “mécanique de la fracture” se doit d’offrir un pilotage précis du processus de délaminage. Ensuite, la courbe globale effort-déplacement est utilisée pour d’obtenir le taux de relaxation d’énergie (CERR) ou la ténacité d’interface. De nombreux auteurs ([10, 11, 12]) ont montré que le CERR dépend grandement du mode d’ouverture. De plus, le fait d’avoir deux matériaux différents à une interface augmente la complexité du post-traitement de la courbe globale effort-déplacement, et ainsi l’identification du CERR. Par conséquent, on cherche une méthode expérimentales qui offre la possibilité de mesurer le CERR, pour tous les modes possibles de I au mode II.

Les techniques actuelles permettant de mesurer le CERR sont par exemple l'essai dit de "double cantilever beam" (DCB) (pure mode I, cf ASTM D 5528-01, 2001, [10, 11, 12]), ou l'essai dit de "end notch flexure" (ENF) (pure mode II, cf [13]) ou enfin l'essai de mode mixte, dit "mixed mode bending" (MMB) (cf [14, 15, 16]). Ces trois types d'essais sont illustrés sur la figure 3. De tous ces appareils d'essais communs, seul le MMB permet de faire des mesures dans une large plage du pure mode I au pure mode II.

Un nouvel appareil d'essai, le micro mixed mode bending set-up (MMMB) a été récemment présenté (cf [17]). Celui-ci permet d'effectuer un essai de délaminage en mode mixte au sein d'un microscope, optique mais aussi électronique à balayage. Ainsi, il est maintenant possible d'obtenir des images de qualité d'un délaminage interfacial. Mais, extraire de ces images les caractéristiques relevant du délaminage reste un challenge. Notablement car les images issues de microscopes présentent un bruit de mesure particulièrement élevée.
<figure><img src="image_4.png" /><figcaption>Fig. 3 – Représentation schématique de trois types classiques d'essai dédiée à la caractérisation d'interface : (a) DCB, (b) ENF et (c) MMB</figcaption></figure>

## 2.2 Procédure inverse

Le maillon central de notre méthode est une procédure de corrélation d'images globale (cf [18]). Soient deux images, (f) (image de référence) et (g) (image déformée). On note la fonction coût (I), du déplacement $\underline{u}$, comme la mesure du résidu sur une région d'intérêt (RdI) entre les deux images (f) et (g), lorsque (g) est corrigée de $\underline{u}$ :

$$I(\underline{u}) = \int_{RdI} \{f(\underline{x}) - g(\underline{x} + \underline{u})\}^2 d\Omega \quad (1)$$

Cette méthode est dite globale dès lors que le déplacement est définit en tout point de la région d'intérêt. Cet espace de déplacement étant de dimension infinie, on le réduit à un espace de dimension finie, que l'on appelle espace d'admissibilité cinématique $\mathcal{U}$. Le déplacement solution est alors celui minimisant la fonction coût $I(\underline{u})$ sur $\mathcal{U}$ :

$$\underline{u}^{\text{opti}} = \underset{\underline{u} \in \mathcal{U}}{\operatorname{Argmin}} \{I(\underline{u})\} \quad (2)$$

Cette méthode est illustrée sur la figure 4

<figure><img src="image_5.png" /><figcaption>Fig. 4 – Schéma simplifié de la procédure inverse dite de corrélation d'image</figcaption></figure>

L'originalité de la méthode présentée est de définir cet espace d'admissibilité cinématique en faisant intervenir des quantités propres au comportement de l'interface. On parle alors de I-DIC pour méthode de corrélation d'image intégrée, puisque des informations sur le comportement du matériau ont donc été intégrées à la procédure d'identification. De plus une unique minimisation est réalisé, puisque l'étape corrélation fournit directement les grandeurs recherchées.

# 3 Validation

La validation de l'outil consiste à sa mise en oeuvre à partir d'images numériquement construites. Illustrée ici en pure mode I, des essais DCB sont donc simulés à l'aide d'une analyse éléments finis. A partir de la solution obtenue, un ensemble d'images est construit. On connaît alors dans ce cas exactement (a) les paramètres définissant le comportement de l'interface ainsi que (b) le bruit ajouté lors de la construction des images. Afin d'avoir une solution FE de la meilleure qualité possible, on utilise un outil de simulation FE où l'interface est modélisée à l'aide d'éléments cohésifs auto-adaptatifs (cf [19]).

## 3.1 Images virtuelles

Avant tout, on simule un essai DCB. On considère ainsi une éprouvette constituée de deux couches d'acier (module d'Young de 210 GPa, et coefficient de poisson de 0.35), de dimension 20 [mm] x 3
[mm] × 0.1 [mm]. On suppose la contrainte $\sigma_{zz}$ nulle et considère ainsi le problème plan. Un problème discret 2D est construit, ou chaque couche est discrétisée par 1000 × 10 éléments. Raccordant les deux couches, des éléments cohésifs sont utilisés pour modéliser l'interface. Pour ce premier exemple, on considère la loi d'ouverture bilinéaire, dont les paramètres sont les suivants :

$$ \varphi = 14.4 \, [\text{MPa}\cdot\mu\text{m}], \text{ et } t_c = 0.137 \, [\text{MPa}] \qquad (3) $$

où $\varphi$ est la travail dissipé lors de l'ouverture, et $t_c$ la contrainte normal maximale. Enfin, on impose un déplacement de +0.25 [mm] et de -0.25 [mm] aux extrémités supérieure et inférieure gauche de l'éprouvette.

La figure 3.1 présente le champs de déplacement de la solution éléments obtenue.

<figure><img src="image_6.png" /><figcaption>Fig. 5 – Carte de déplacement (haut : composante -x, bas : composante -y) de la solution éléments finis</figcaption></figure>

Une fois le problème éléments finis résolu, la construction des images virtuelles se fait en deux étapes. La première consiste en la définition d'un mouchetis réaliste, qui est associé à la géométrie de l'éprouvette. Ensuite, on applique sur ce mouchetis le déplacement solution du problème éléments finis simulant l'essai DCB. Avec toujours l'objectif d'avoir un ensemble de photos réaliste, on ajoute en même temps un bruit blanc.

Illustrée en haut de la figure 6, l'image de référence (f) est la réduction d'un mouchetis à la géométrie de l'éprouvette. Ce mouchetis est construit en sommant trois bruits Gaussien, de même moyenne (1/6) et de même variance (10%). Ces trois champs sont ensuite interpolés sur trois échelles différentes afin de construire le mouchetis. On utilise les longueurs caractéristiques suivantes : 1, 8 et 32 [pixels].

<figure><img src="image_7.png" /><figcaption>Fig. 6 – Illustration d'un mouchetis virtuel et de sa réduction à une éprouvette (cf (f)) ainsi que de sa déformée (g), calculée en y appliquant une solution éléments finis et un bruit blanc. La résolution de ces images et de 1024 par 80 pixels, avec une dynamique codée sur 8 bits</figcaption></figure>

Afin d'appliquer le déplacement éléments finis à l'image de référence, un problème non linéaire peut être résolu pour trouver quel point de (f) se retrouve aux coordonnées des pixels dans (g) (cf [20]). Un bruit blanc Gaussien est ensuite ajouté, toujours dans un soucis de travailler avec des images réalistes. Ce bruit blanc est définit avec une variance de 3%. En effet, les appareils photos classiques présentent un bruit blanc qui peut facilement être modélisé par un bruit Gaussien, mais de variance alors plus faible puisque d'ordre de grandeur de 1%. La figure 6 présente l'image (g), qui est la déformée de (f) par application de la solution FE et du bruit blanc.
## 3.2 Application

On s'intéresse dans cette partie à l'étude des premiers résultats de la procédure inverse proposée, appliquée aux images virtuelles. On étudie tout d'abord en détails une première application où une loi bilinéaire a été utilisée pour représenter le comportement d'interface. Ensuite on s'intéressera à valider la méthode sur une large gamme de loi d'ouverture.

### 3.2.1 Cas d'une loi bilinéaire

La loi considérée est la même que celle utilisée dans la partie 3.1. La figure 7 présente le résidu final après corrélation. Cette carte des résidus met en avant la faible résolution spatiale des images. En effet, seuls 8 pixels représentent chacun des deux substrats. Ensuite, on peut remarquer deux points propres à l'étude d'une carte des résidus :

- la dynamique de la carte des résidus correspond au bruit blanc ajouté lors de la construction de l'image (g);

- l'amplitude des résidus est homogène sur l'ensemble de la carte.

De ces deux points, on peut déduire que la définition de l'espace cinématique est suffisamment riche pour décrire un essai DCB.

<figure><img src="image_8.png" /><figcaption>Fig. 7 – Résidu final. On peut relever deux remarques notables : (1) le dynamique du résidu est du même ordre de grandeur que le bruit ajouté lors de la construction des images et (2) la carte de résidu est homogène</figcaption></figure>

Cette étape de validation, basée sur une simulation éléments finis, permet de comparer directement les champs déplacement. La figure 8 illustre ainsi l'erreur en déplacement (composante-x en haut, composante-y en bas) entre la solution EF et la solution retrouvée par la procédure inverse. Encore, les ordres de grandeurs des deux champs sont faibles (3% [pixels] pour les déplacement verticaux, 0.3% [pixels] pour les déplacement horizontaux).

<figure><img src="image_9.png" /><figcaption>Fig. 8 – Cartes d'erreur en déplacement entre la solution éléments finis et la solution retrouvée par la méthode inverse proposée (en haut : composante-x, en bas : composante-y). Les unités sont en [mm], (1 pixel = 8.73 [µm])</figcaption></figure>
La carte d'erreur en déplacement met aussi en avant la non-homogénéité de l'erreur. Ainsi, on re-
marque que proche de conditions aux limites particulières (comme, par exemple, l'interface) l'erreur est
la plus élevée.

Enfin, la figure 9 présente la loi d'ouverture qu'il est possible d'extraire de ces images. Sur ce gra-
phique, plusieurs lois sont représentées. Tout d'abord, en bleu, la loi d'entrée, i.e. la loi utilisée pour la
simulation EF est tracée. Ensuite, 5 courbes sont représentées, toutes retraçant les lois d'ouverture que la
procédure inverse permet de retrouver à partir de 5 couples d'images. Chacun de ces couples est construit
en ajoutant un bruit blanc Gaussien différent, allant d'un bruit nul à un bruit d'une variance de 5%.

<figure><img src="image_10.png" /><figcaption>Fig. 9 – Comparaison de la loi d'ouverture choisie pour décrire l'interface dans l'approche EF et des lois d'oubverture retrouvées par la procédure inverse (pour différentes variances du bruit blanc Gaussien ajouté)</figcaption></figure>

Dans l'étude présentée, la loi d'entrée est bilinéaire. Or, dans notre procédure inverse, on n'introduit pas une telle discontinuité. Il est donc impossible de retrouver exactement la même loi. Toutefois, on illustre la faisabilité de la méthode proposé sur cet exemple précis pour montrer que même si on ne peut capturer l'évolution exacte de la loi d'ouverture, les principales caractéristiques sont estimées avec une grande précision. En effet, l'effort d'interface critique, l'énergie dissipée ainsi que l'ouverture critique sont très bien identifiés.

Enfin, on peut retenir de cette dernière figure la robustesse de la méthode introduite. En effet, même avec un bruit élevé, la loi d'ouverture retrouvée permet toujours de caractériser avec précision les caractéristiques les plus importantes.

3.2.2 Applications à plusieurs lois d'ouverture

On s'intéresse maintenant à valider la flexibilité de la méthode inverse introduite. Pour cela, plusieurs
simulations EF sont réalisées, utilisant la même géométrie pour l'éprouvette que celle utilisée dans 3.2.1.

<table>
  <thead>
    <tr>
      <th rowspan="2">Application</th>
      <th colspan="4">Modèle de zone cohésive</th>
      <th>Chargement</th>
      <th colspan="2">Propriétés des images</th>
    </tr>
    <tr>
      <th>TSL</th>
      <th>φ<br>[MPa.µm]</th>
      <th>t<sub>c</sub><br>[MPa]</th>
      <th>Δ<sub>c</sub><br>[µm]</th>
      <th>U<sub>p</sub><br>[mm]</th>
      <th>Résolution<br>[pixels<sup>2</sup>]</th>
      <th>Taille des pixels<br>[µm]</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <td>#1</td>
      <td>PPR</td>
      <td>33</td>
      <td>1</td>
      <td>11.0</td>
      <td>0.25</td>
      <td>1024×80</td>
      <td>8.77</td>
    </tr>
    <tr>
      <td>#2</td>
      <td>PPR</td>
      <td>33</td>
      <td>0.3</td>
      <td>68.9</td>
      <td>0.25</td>
      <td>1024×80</td>
      <td>8.77</td>
    </tr>
    <tr>
      <td>#3</td>
      <td>bi-linear</td>
      <td>14.4</td>
      <td>0.137</td>
      <td>30.0</td>
      <td>0.25</td>
      <td>1024×80</td>
      <td>8.77</td>
    </tr>
    <tr>
      <td>#4</td>
      <td>PPR</td>
      <td>3.3</td>
      <td>0.3</td>
      <td>6.9</td>
      <td>0.05</td>
      <td>1024×62</td>
      <td>4.85</td>
    </tr>
    <tr>
      <td>#5</td>
      <td>PPR</td>
      <td>33</td>
      <td>0.3</td>
      <td>23.0</td>
      <td>0.30</td>
      <td>1024×91</td>
      <td>8.79</td>
    </tr>
  </tbody>
</table>

Tableau 1 – Propriétés utilisées pour chaque application. Le modèle PPR a besoin de quatre paramètres, qui sont φ et t_c mais aussi un paramètre de forme ainsi qu'un paramètre de pente initiale. Pour simplifier le tableau, ces deux paramètres n'y sont pas présentés
Par contre, le comportement de l'interface est maintenant décrit en utilisant le modèle dit PPR (cf [21]).
Ce modèle de zone cohésif, polynomial, offre une grande flexibilité pour la forme de la loi d'ouverture.
Le tableau 1 récapitule les paramètres utilisés, ainsi que les propriétés des images.

<figure><img src="image_11.png" /><figcaption>Fig. 10 – Comparaison entre différentes lois d'ouverture choisies pour décrire l'interface dans l'approche EF (traits pointillés) et les lois retrouvées par la procédure inverse proposée (traits continus). Pour tous ces essais virtuels, un bruit blanc Gaussien, d'une variance de 3%, a été ajouté aux images</figcaption></figure>

La figure 10 illustre les lois d'entrée, utilisées pour la simulations EF (lignes pointillées), ainsi que les lois que la procédure introduite est capable de retrouver (lignes continues), pour chacune des applications. De la même manière que pour l'étude d'une unique loi d'ouverture (cf Section.3.2.1), les caractéristiques importantes du phénomène de délaminage sont très bien retrouvées, quelque soit la loi d'entrée.

# 4 Conclusion et perspectives

Appliquée à la caractérisation d'une interface cohésive, les premiers résultats de la procédure inverse proposée dans ce papier sont prometteurs. Cette procédure, basée sur technique inverse dite de corrélation d'image intégrée, utilise une écriture fine des champs de déplacement dans le cas d'un essai type DCB et relie ceux-ci aux paramètres d'interface.

Dans ce papier, on s'est intéressé à la validation de la méthode. Celle-ci a été rendu possible par la fabrication d'images virtuelles, à l'aide d'un outils de simulations éléments finis d'essais DCB. À partir de ces images, les paramètres caractéristiques du modèle cohésif du délaminage, que sont l'effort d'interface critique, l'énergie dissipé lors de son ouverture ou la taille de la zone cohésive sont clairement bien retrouvés. De plus, on a monté que l'ajout d'un bruit blanc réaliste aux images ne change pas la qualité des résultats. La méthode proposée est donc robuste.

Aujourd'hui, des essais sont en cours de réalisation afin d'appliquer la méthode présentée sur des images "réelles" et de caractériser ainsi, d'une manière précise, le comportement d'interfaces. De plus, la méthode introduite n'est validée que dans des cas d'ouverture en pure mode I, sa validation dans des cas de mode II ainsi que de modes mixtes est aussi en cours.

# Références

[1] D.S. Dugdale. Yielding of steel sheets containing slits. *Journal of the Mechanics and Physics of Solids*, 8(2):100–104, 1960.

[2] GI Barenblatt. The mathematical theory of equilibrium cracks in brittle fracture. *Advances in applied mechanics*, 7(55-129):104, 1962.

[3] X. Xu and A. Needleman. Void nucleation by inclusion debonding in a crystal matrix. *Modelling and Simulation in Materials Science and Engineering*, 1:111, 1993.
[4] X.P. Xu and A. Needleman. Numerical simulations of fast crack growth in brittle solids. *Journal of the Mechanics and Physics of Solids*, 42(9):1397-1434, 1994.

[5] AG Evans, M. Ruhle, BJ Dalgleish, and PG Charalambides. The fracture energy of bimaterial interfaces. *Materials Science and Engineering: A*, 126(1-2):53-64, 1990.

[6] W. Li and T. Siegmund. An analysis of the indentation test to determine the interface toughness in a weakly bonded thin film coating-substrate system. *Acta materialia*, 52(10):2989-2999, 2004.

[7] L. Banks-Sills, V. Boniface, and R. Eliasi. Development of a methodology for determination of interface fracture toughness of laminate composites—the 0°/90° pair. *International journal of solids and structures*, 42(2):663–680, 2005.

[8] J.A.J. Wang, I.G. Wright, M.J. Lance, and K.C. Liu. A new approach for evaluating thin film interface fracture toughness. *Materials Science and Engineering: A*, 426(1-2):332-345, 2006.

[9] I. Ocaña, JM Molina-Aldareguia, D. Gonzalez, MR Elizalde, JM Sánchez, JM Martínez-Esnaola, J. Gil Sevillano, T. Scherban, D. Pantuso, B. Sun, et al. Fracture characterization in patterned thin films by cross-sectional nanoindentation. *Acta materialia*, 54(13):3453-3462, 2006.

[10] JW Hutchinson and Z. Suo. Mixed mode cracking in layered materials. *Advances in applied mechanics*, 29(63):191, 1992.

[11] A. Kuhl and J. Qu. A technique to measure interfacial toughness over a range of phase angles. *Journal of Electronic Packaging*, 122:147, 2000.

[12] XQ Shi, XR Zhang, and JHL Pang. Determination of interface fracture toughness of adhesive joint subjected to mixed-mode loading using finite element method. *International journal of adhesion and adhesives*, 26(4):249-260, 2006.

[13] LA Carlsson, JW Gillespie, and RB Pipes. On the analysis and design of the end notched flexure (enf) specimen for mode ii testing. *Journal of Composite Materials*, 20(6):594, 1986.

[14] J.R. Reeder and J. Crews. Mixed-mode bending method for delamination testing. *AIAA journal*, 28(7):1270-1276, 1990.

[15] WO Soboyejo, G.Y. Lu, S. Chengalva, J. Zhang, and V. Kenner. A modified mixed-mode bending specimen for the interfacial fracture testing of dissimilar materials. *Fatigue & Fracture of Engineering Materials & Structures*, 22(9):799-810, 1999.

[16] N. Blanco, E.K. Gamstedt, J. Costa, and D. Trias. Analysis of the mixed-mode end load split delamination test. *Composite structures*, 76(1-2):14-20, 2006.

[17] M. Kolluri, JPM Hoefnagels, JAW van Dommelen, and MGD Geers. An improved miniature mixed-mode delamination setup for in situ microscopic interface failure analyses. *Journal of Physics D: Applied Physics*, 44:034005, 2011.

[18] G. Besnard, F. Hild, and S. Roux. Finite-element displacement fields analysis from digital images : application to portevin-le Châtelier bands. *Experimental Mechanics*, 46(6):789-803, 2006.

[19] M. Samimi. *A self-adaptative cohesive zone model for interfacial delamination*. PhD thesis, TU/e, 2011.

[20] Julien Réthoré, Muhibullah, Thomas Elguedj, Michel Coret, Philippe Chaudet, and Alain Combescure. Robust identification of elasto-plastic constitutive law parameters from digital images using 3d kinematics. *International Journal of Solids and Structures*, 50(1):73-85, 2013.

[21] K. Park, G.H. Paulino, and J.R. Roesler. A unified potential-based cohesive model of mixed-mode fracture. *Journal of the Mechanics and Physics of Solids*, 57(6):891-908, 2009.