# Analyse isogéométrique et corrélation d'image

Julien Réthoré, Thomas Elguedj

► To cite this version:

Julien Réthoré, Thomas Elguedj. Analyse isogéométrique et corrélation d'image. CFM 2009 - 19ème Congrès Français de Mécanique, Aug 2009, Marseille, France. hal-00517508

HAL Id: hal-00517508

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

Submitted on 26 Mar 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

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

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Analyse isogéométrique et corrélation d'images

J. RÉTHORÉ, T. ELGUEDJ

Université de Lyon, CNRS
INSA-Lyon, LaMCoS UMR5259, F-69621, Villeurbanne, France

## Résumé :

Depuis quelques années, certains groupes de recherche se sont intéressés à l'utilisation des fonctions utilisées dans les modèles de CAO (NURBS) pour les simulations numériques en mécanique. D'un point de vue analytique, ces fonctions présentent plusieurs avantages en vue d'un enrichissement de type k. On montrera comment et pourquoi utiliser ces fonctions dans le cadre de la mesure de champs de déformation par corrélation d'images numériques en s'appuyant sur des cas artificiels.

## Abstract :

During the past few years, research groups have been interested in using Computer Aided Design (CAD) functions (NURBS) for the numerical analysis of mechanical problem. From a mathematical point of view, these functions have several advantages for k refinement. Based on an artificial example, we show how and why one may use NURBS functions in the context of strain measurement by digital image correlation.

Mots clefs : corrélation d'images, déformation, NURBS

## 1 Introduction

La corrélation d'images numériques ([1] par exemple ou [2] pour une caractérisation des performances de cette technique) est aujourd'hui une des techniques de mesure de champs cinématiques les plus utilisées que ce soit dans un cadre académique ou industriel. Le plus souvent, les champs de déplacement qui sont obtenus sont dérivés pour obtenir des champs de déformation qui sont ensuite utilisés à des fins d'identification. L'effet du bruit de mesure est inévitablement amplifié et affecte la qualité de l'étape finale d'identification. Il peut donc être judicieux d'essayer de diminuer l'influence du bruit, soit a posteriori ([3] pour l'extraction des facteurs d'intensité des contraintes en mécanique de la rupture) ou alors lors de la mesure elle-même. C'est cette deuxième solution que nous exploitons ici.

## 2 Corrélation d'images

L'équation de conservation des niveaux de gris pour la détermination du flot optique que nous cherchons à résoudre dans le cadre de la corrélation d'images

$$f(\mathbf{x}) = g(\mathbf{x} - \mathbf{u}(\mathbf{x})) \quad (1)$$

est un problème mal posé qui nécessite donc une régularisation. Plus on intensifie cette régularisation plus le conditionnement du système s'améliore et plus la sensibilité au bruit et les incertitudes de mesures chutent. Les approches habituelles consistent à travailler sur des zones d'études indépendantes (qui se recouvrent) dont la taille est importante et à l'intérieure desquelles la description cinématique est raisonnablement pauvre.

Une alternative proposée par [4] consiste à résoudre l'équation (1) en moyenne et à assurer la continuité du champ de déplacement sur l'intégralité du domaine étudié en utilisant des fonctions éléments finis. Cette stratégie permet d'obtenir des niveaux de performance supérieure aux approches usuelles mais le champ de déformation qui en résulte reste assez pauvre en terme de description (linéaire et continu par morceau). En vue d'une utilisation de ce champ de déformation pour l'identification de lois de comportement il peut être intéressant d'opérer un raffinement dit p des fonctions éléments finis (augmentation du degré). A titre d'exemple, pour un élément quadrangle, le passage d'une interpolation de degré 1 à une interpolation de degrés 2 nécessite l'ajout de 4 degrés de liberté par composante et surtout ne permet d'assurer qu'une continuité $C^0$ de la solution. La discontinuité du champ de déformation à la frontière de éléments est donc maintenue et l'augmentation du coût (en terme de nombre de degrés de liberté) est prohibitive. L'augmentation du degré d'interpolation dans le cadre des éléments finis permet d'augmenter la continuité mais seulement à l'intérieur des éléments. L'augmentation du degré de continuité est donc local et ne permet pas de mieux contraindre le
problème de conservation du flot optique de manière globale. De plus l'accroissement du nombre des inconnues à déterminer rend l'utilisation d'une telle stratégie inappropriée dans le cadre de la corrélation d'image. Dans la suite, on propose d'utiliser des fonctions connues sous l'appellation de NURBS et dont l'augmentation du degré génére également une augmentation de la classe de continuité et ce à moindre coût (faible augmentation du nombre d'inconnues).

## 3 Analyse isogéométrique

<figure><img src="image_2.png" /><figcaption>FIG. 1 – Fonctions B-Spline pour un vecteur de noeuds périodique. En haut : cas impair p=3 ; en bas, cas pair p=2.</figcaption></figure>

L'analyse isogéométrique est une nouvelle méthode numérique développée par Hughes et al. depuis quelques années (cf. [5]). Elle introduit un nouveau paradigme qui tente d'améliorer la description géométrique, la représentation de solutions complexes et le raffinement de maillage tels qu'ils sont connus dans les éléments finis traditionnels. L'analyse isogéométrique est fondée sur les bases de la représentation géométrique utilisée en Conception Assistée par Ordinateur et en infographie. La première implémentation de l'analyse isogéométrique fut basée sur les fonctions NURBS (Non-Uniform Rational B-Splines), une technologie standard dans la majorité des modeleurs CAO du marché [6], et a été récemment étendue aux T-Splines. Comme il a été démontré dans de nombreux travaux sur le sujet, la méthode possède de nombreux avantages par rapport à la méthode des éléments finis traditionnels.

### 3.1 B-Splines et NURBS

Les B-Splines sont des fonctions polynomiales par morceaux possédant une classe de continuité donnée. Elles sont construites à partir d'un vecteur de noeuds, qui est lui-même constitué d'un ensemble de coordonnées dans l'espace paramétrique $\Xi = \{\xi_1, \xi_2, \dots, \xi_{n+p+1}\}$, où $\xi_i \in \mathbb{R}$ est le $i^{\text{ème}}$ noeud, $i$ est l'index, $i = 1, 2, \dots, n + p + 1$, $p$ est le degré polynomial et $n$ est le nombre de fonctions considérées. Les fonctions d'un ordre $p$ donné quelconque sont définies récursivement à l'aide de la formule de Cox-de Boor à partir du vecteur de noeud $\Xi$. On commence par définir les fonctions constantes par morceaux ($p=0$) :

$$ N_{i,0}(\xi) = \begin{cases} 1 & \text{if } \xi_i \le \xi < \xi_{i+1}, \\ 0 & \text{otherwise.} \end{cases} \qquad (2) $$

Pour $p \ge 1$, la relation de récurrence est :

$$ N_{i,p}(\xi) = \frac{\xi - \xi_i}{\xi_{i+p} - \xi_i} N_{i,p-1}(\xi) + \frac{\xi_{i+p+1} - \xi}{\xi_{i+p+1} - \xi_{i+1}} N_{i+1,p-1}(\xi). \qquad (3) $$

Une courbe B-Spline dans $\mathbb{R}^{d_s}$ estdéfinie comme suit :

$$ \mathbf{C}(\xi) = \sum_{i=1}^{n} \mathbf{P}_{i} N_{i,p}(\xi), \qquad (4) $$
où $\mathbf{P}_i \in \mathbb{R}^{d_s}$ est le i-ème point de contrôle. Les techniques de raffinement en espace, dit “h”, et en degré polynomial, dit “p”, telles qu’elles existent en éléments finis classiques sont tranposables ici. La différence étant qu’elle ne sont pas commutatives. On effectuant le raffinement en degré polynomial avant le raffinement en espace, on obtient un nombre assez faible de fonctions B-Spline ayant une continuité maximale sur la majeure partie du domaine. Cette technique, appelée raffinement $k$, est un avantage majeur sur les éléments finis traditionnels. Les concepts précédant sont étendus en 2D et 3D par une structure de type produit tensoriel comme c’est le cas en éléments finis. Le fait de travailler sur un maillage rectangulaire régulier dans l’approche présentée ici nous permet de travailler directement avec des fonctions B-Spline (et non des NURBS) et de faire partiellement abstraction de la notion de points de controles puisque la géométrie n’est pas représentée explicitement. De plus, contrairement à la forme appliquée à la simulation numérique de l’analyse isogéométrique, il n’est pas nécessaire d’imposer des conditions aux limites de type Dirichlet au problème considéré. Ceci nous permet de travailler avec des vecteurs de noeuds périodiques, et non ouverts, et de maximiser la continuité des fonctions ($p-1$ pour des fonctions de degrés $p$) sur l’ensemble du domaine y compris sur les bords (Figure 1). De plus on constate que l’augmentation de degré $p$ à $p+1$ de ces fonctions ne “coûte” qu’un seul degré de liberté par dimension ce qui rend le raffinement de type $k$ envisageable dans le cadre de corrélation d’images numériques.

## 4 Résolution

<figure><img src="image_3.png" /><figcaption>FIG. 2 – Déplacement imposé et déformation correspondante</figcaption></figure>

Après avoir choisi une base de discrétisation

$$ \mathbf{u}(\mathbf{x}) = \sum_{i \in \mathcal{N}} u_i \mathbf{N}_{i,p}(\mathbf{x}) = [\mathbf{N}(\mathbf{x})] \{\mathbf{U}\}, \quad (5) $$

on minimise l’écart quadratique moyen à la conservation du flot optique :

$$ \min_{u_i} \iint_{\Omega} [f(\mathbf{x}) - g(\mathbf{x} - u_i N_{i,p}(\mathbf{x}))]^2(\mathbf{x}) d\mathbf{x}. \quad (6) $$

En utilisant un algorithme de Newton, on résout une succession de systèmes linéaires

$$ [\mathbf{M}]^i \{d\mathbf{U}\} = \{\mathbf{b}\}^i, \quad (7) $$

avec

$$ M_{nm}^i = \iint_{\Omega} (\mathbf{N}_n \cdot \nabla g(\mathbf{x} - [\mathbf{N}]\{\mathbf{U}\}^i))(\mathbf{N}_m \cdot \nabla g(\mathbf{x} - [\mathbf{N}]\{\mathbf{U}\}^i))d\mathbf{x}, \quad (8) $$

et

$$ b_n^i = \iint_{\Omega} (\mathbf{N}_n \cdot \nabla g(\mathbf{x} - [\mathbf{N}]\{\mathbf{U}\}^i))(f(\mathbf{x}) - g(\mathbf{x} - [\mathbf{N}]\{\mathbf{U}\}^i))d\mathbf{x}. \quad (9) $$

Cette technique de résolution globale proposée dans [4] permet d’utiliser des fonctions définies non pas par morceaux comme dans les approches usuelles mais sur l’intégralité du domaine d’étude. Notons que des fonctions B-Spline ont déjà été utilisées dans le cadre d’une technique de corrélation usuelle [7]. Ici les gradients
des images sont calculés par différences finies et l'interpolation subpixel des niveaux de gris est de type spline cubique. Comme dans [4] et comme décrit en détail dans [8], un algorithme multi-échelle est mis en place pour restaurer progressivement les plus petites longueurs d'onde de la solution. Par soucis de simplicité, les échelles grossières sont décrites avec des fonctions de degrés 1. Notons que des alternatives sont possibles.

## 5 Exemple

<figure><img src="image_4.png" /><figcaption>FIG. 3 – Carte de déformation pour 16P1 (à gauche) et 16P3 (à droite)</figcaption></figure>

En vue de tester les performances des fonctions NURBS dans le cadre de la corrélation d'images, on impose artificiellement à une image de 512 × 512 pixels, un déplacement vertical de type sinus à fréquence variable :

$$U_x(X) = U_0 \sin \left[ \frac{\omega_x}{2\pi} \left( 1 - \left| 2 \frac{X}{\max(X)} - 1 \right| \right)^2 \right], \quad (10)$$

avec $\omega_x = 9\pi^2$ et $U_0 = 1$ pixel. Ceci produit une déformation à fréquence et amplitude variable telle que présentée sur la Figure 2. Ainsi tous les pixels suivant la direction y ont le même déplacement et une analyse systématique de l'erreur commise peut être menée.

Sur la Figure 4, on peut observer le champ de déplacement expérimental et théorique pour des éléments de 16 pixels et des fonctions B-Spline de degré 3 ou des fonctions de degrés 1. La Figure 3 montre quant à elle les cartes de déformation pour ces deux types de description du champ de déplacement. On observe qu'une description P1 du déplacement ne suffit à capturer efficacement le déplacement imposé à l'image. De plus, dans la direction dans laquelle l'image a été déformée, la déformation avec des fonction de degrés 1 est constante par morceaux et bien que les niveaux obtenus aux centres des éléments soient satisfaisants les discontinuités de déformation aux frontières des éléments peuvent s'avérer pénalisantes. Au contraire, à taille d'élément équivalente (16 pixels), des fonctions de degrés 3 se comportent beaucoup mieux. On peut constater que la déformation est également décrite de manière précise.

La Figure 5 montre les résultats de l'analyse systématique de l'erreur. On peut observer sur cette figure l'évolution de l'incertitude (écart type de la différence entre le champ imposé et le champ mesuré) sur le déplacement et la déformation en fonction de la taille d'élément pour $p = 1$ à $p = 6$. Pour $p = 1$, l'erreur de discretisation l'emporte et l'incertitude en déplacement augmente quand la taille des éléments augmente. Pour la déformation on obtient une valeur quasiment constante de l'ordre de $1.5 \times 10^{-3}$. Dès lors que l'on augmente le degré d'interpolation en maintenant une continuité maximale, on peut distinguer deux régimes : pour les petites tailles d'éléments c'est l'erreur de corrélation qui domine (l'erreur due au fait que le problème est mal posé) puis quand on augment la taille des éléments c'est l'erreur de discrétisation qui devient prépondérante. Il existe donc une taille d'élément optimale, par exemple 16 pixels pour $p = 3$ pour laquelle on a une incertitude de $2.5 \times 10^{-3}$ pixel pour le déplacement et $3.0 \times 10^{-4}$ pour la déformation. La Figure 6 montre l'évolution de l'erreur moyenne estimée dans les mêmes conditions. Les niveaux sont inférieurs à ceux obtenus pour l'incertitude et c'est l'erruer de discrétisation qui domine puisque l'on obtient des niveaux d'erreur qui augmentent quand la taille des éléments augmente. Pour 16 pixels et $p = 3$, on a une erreur moyenne de $1.0 \times 10^{-5}$ pixel pour le déplacement et $2.0 \times 10^{-6}$ pour la déformation.
<figure><img src="image_5.png" /><figcaption>FIG. 4 – Comparaison entre les champs mesurés et la solution pour 16P1 (en haut) et 16P3 (en bas) : déplacement (à gauche) et déformation (à droite). L’épaisseur du trait bleu correspond aux différentes valeurs obtenues dans la direction *y* orthogonale au déplacement imposé.</figcaption></figure>

## 6 Conclusion

On montre dans cette contribution que les fonctions NURBS pour lesquelles on procède à un raffinement de type *k* peuvent être utilisées avantageusement dans le cadre de la mesure de champs de déformation par corrélation d’images numériques. En effet, une telle stratégie de raffinement permet d’augmenter le degré d’interpolation tout en conservant une continuité maximale et à un coût très faible en terme d’augmentation du nombre de degrés de liberté. Ainsi le conditionement du problème de corrélation ne s’en trouve que peu dégradé et les perfomances sont accrues. On montre par un exemple où une image est déformée artificiellement les bénéfices apportés par la stratégies proposée qui permet *in fine* d’obtenir des champs de déformation continus avec une grande précision et une grande robustesse.

## Références

[1] Sutton M., Wolters W., Peters W., Ranson W., and McNeill S. Determination of displacements using an improved digital correlation method. Image Vision Computing, 1(3), 133–139, 1983.

[2] Bornert M., Brémand F., Doumalin P., Dupré J., Fazzini M., Grédiac M., Hild F., Mistou S., Molimard J., Orteu J., Robert L., Surrel Y., Vacher P., and Wattrisse B. Assessment of digital image correlation measurement errors : Methodology and results. Experimental Mechanics, 2008.

[3] Réthoré J., Roux S., and Hild F. Noise-robust stress intensity factor determination from kinematic field measurements. Engineering Fracture Mechanics, 75(13), 3763–3781, 2006.

[4] Besnard G., Hild F., and Roux S. ‘finite-element’ displacement fields analysis from digital images : Application to portevin-le châtelier bands. Experimental Mechanics, 46(6), 789–803, 2006.
<figure><img src="image_6.png" /><figcaption>FIG. 5 – Incertitude versus taille d'élément pour différentes valeurs de p : déplacement (à gauche) et déformation (à droite)</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>FIG. 6 – Erreur moyenne versus taille d'élément pour différentes valeurs de p : déplacement (à gauche) et déformation (à droite)</figcaption></figure>

[5] Hughes T. J., Cottrell J., and Bazilevs Y. Isogeometric analysis : CAD, finite elements, NURBS, exact geometry and mesh refinement. Computer Methods in Applied Mechanics and Engineering, 194, 4135–4195, 2005.

[6] Piegl L. and Tiller W. The NURBS book. Springer, 1997.

[7] Cheng P., Sutton M. A., Schreier H. W., and McNeill S. R. Full-field speckle pattern image correlation with b-spline deformation function. Experimental Mechanics, 42(3), 344–352, 2002.

[8] Réthoré J., Hild F., and Roux S. Shear-band capturing using a multiscale extended digital image correlation technique. Computer Methods in Applied Mechanics and Engineering, 196(49-52), 5016–5030, 2007.