# Modélisation numérique de l'étendard d'une plume, une membrane fortement anisotrope

Jean Jouve, Florence Bertails-Descoubes, Victor Romero, Rahul Narain,
Theodore Kim

► To cite this version:

Jean Jouve, Florence Bertails-Descoubes, Victor Romero, Rahul Narain, Theodore Kim. Modélisation numérique de l'étendard d'une plume, une membrane fortement anisotrope. 16ème Colloque National en Calcul de Structures (CSMA 2024), CNRS; CSMA; ENS Paris-Saclay; CentraleSupélec, May 2024, Hyères, France. pp.1-8. http://hal-04610994

HAL Id: hal-04610994

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

Submitted on 3 Dec 2024

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 Attribution 4.0 International License
# Modélisation numérique de l'étendard d'une plume, une membrane fortement anisotrope

Jean Jouve<fnref n="1" />, Victor Romero<fnref n="1" />, Rahul Narain<fnref n="2" />, Theodore Kim<fnref n="3" />, Florence Bertails-Descoubes<fnref n="1" />

**Résumé** — À partir de mesures expérimentales, nous proposons un modèle mécanique simple pour capturer le comportement macroscopique de la membrane d'une plume d'oiseau. Cette membrane étant fortement anisotrope, résoudre le problème résultant par des éléments finis standards conduit à un phénomène de verrouillage numérique ainsi qu'à un mauvais conditionnement du problème discrétisé. Pour limiter le verrouillage et améliorer le conditionnement, nous alignons les éléments du maillage sur la direction de forte raideur et traitons la direction de forte raideur par une contrainte d'inextensibilité.

**Mots clés** — Élasticité, Anisotropie, Verrouillage numérique.

## 1 Introduction

Les pennes (ou plumes de contour) sont des plumes qui parent les ailes et queux des oiseaux. En raison de leur structure multi-échelle hiérarchique, elles ont des propriétés hors du commun. Cette structure est constituée d'une tige centrale nommé *rachis*, de cette tige sortent des centaines d'autres appelées *barbes*, et sur chacune d'elles se trouve une centaine de branches nommées *barbules*. Ces barbules ont une structure qui diffère en fonction de leur orientation par rapport à la barbe. Celles dirigées vers la pointe extérieure du rachis ont des crochets à leur extrémité, celles dirigées vers la base du rachis ont la forme d'une gouttière dans laquelle les crochets peuvent venir se loger. Cette particularité permet aux barbes de s'accrocher entre elles, formant ainsi une surface le long du rachis dénommée *étendard*. Cette structure permet à l'étendard de dissiper les chocs en se fracturant de façon réversible. L'oiseau peut ensuite reformer l'étendard abîmé en accrochant les barbules entre elles à travers le passage de son bec dans ses plumes [11].

<figure><img src="image_2.png" /><figcaption>FIGURE 1 — Schéma de la structure hiérarchique d'une plume. De chaque côté de la tige centrale, nommée *rachis*, se trouvent les *barbes*. Comme le rachis, les barbes sont bordées d'autres tiges, *barbules*. Dans le cas spécifique des pennes, les barbules pointant vers le haut de la plume ont des crochets qui leur permettent de s'accrocher à d'autres barbules. Les barbes attachées entre elles par leurs barbules forment l'étendard.</figcaption></figure>

En plus d'avoir ce comportement inélastique particulier, l'étendard se comporte comme une coque fortement anisotrope, similaire aux coques renforcées par des tiges inextensibles. En raison de la grande rigidité des barbes comparée à celle des barbules, toute déformation de l'étendard qui introduit une extension des barbes a une raideur supérieure de plusieurs ordres de grandeur à celle des déformations évitant cette extension.

<footnote n="1">¹ Univ. Grenoble Alpes, Inria, CNRS, Grenoble INP, LJK, {jean.jouve, florence descoubes} @ inria.fr</footnote>

<footnote n="2">² Indian Institute of Technology Delhi</footnote>

<footnote n="3">³ Yale University</footnote>
Nous considérons que ces propriétés hors du commun de la plume en font un objet d'étude intéressant pour augmenter les connaissances concernant l'ornithologie, permettre la création future de matériaux par biomimétisme, ainsi que produire des outils aidant à la création d'images animés par ordinateur de plumes. À notre connaissance, aucune modélisation de ce comportement de l'étendard n'existe à ce jour. Les études de la mécanique des pennes se sont concentrées sur des aspects différents. L'interaction entre deux barbes due à l'imbrication des barbules a été étudiée à plusieurs reprises [4, 3, 11]. Certaines études ont aussi analysé le comportement mécanique de la penne complète, rachis compris [1, 12]. Ainsi, nous proposons un modèle de membrane caractérisant la forte anisotropie de la plume.

La modélisation numérique de matériaux fortement anisotropes voire quasiment inextensibles dans une direction ou deux directions, est confrontée à un phénomène de verrouillage lorsque des éléments finis standards sont utilisés [9, 10]. C'est-à-dire que la vitesse de convergence de ces méthodes est considérablement réduite, ou que leur convergence vers la solution du problème continue est impossible [2]. Plusieurs travaux ont proposé des solutions pour limiter ces effets. Par exemple, dans le cadre de matériaux ayant deux directions quasi inextensibles, le verrouillage est soit résolu par un alignement des éléments avec ces directions, soit par une intégration réduite de la partie raide de l'énergie élastique, ou encore par l'introduction d'éléments mixtes [10, 7]. Pour les cas où l'extension dans une direction seulement est raide, par exemple pour les matériaux renforcés par des fibres, les solutions proposées se reposent sur l'utilisation d'éléments finis mixtes [9, 6].

Notre modèle de plume présentant une grande raideur dans une direction (celle des barbes), nous avons choisi de le résoudre numériquement à travers des éléments finis linéaires alignés selon la direction de plus grande raideur. De plus, afin d'améliorer le conditionnement du problème, nous remplaçons la partie raide de l'énergie par une contrainte.

Pour résumer, les contributions qui seront présentées sont :

- Un modèle de membrane pour l'étendard des pennes ;

- Une comparaison de ce modèle avec des tests de traction sur des échantillons réels de plumes de cygne ;

- Une analyse du phénomène de verrouillage et du mauvais conditionnement du problème qui émerge de la forte anisotropie présente dans notre modèle. Cette étude confirme les analyses précédentes de l'état de l'art.

## 2 Méthode expérimentale

Afin de modéliser les étendards de pennes et de valider le modèle résultant, nous avons analysé plusieurs échantillons de plumes de cygne à travers des tests de traction. Comme la plume est clairement anisotrope due à la présence des barbes, les échantillons sont découpés à différents angles θ entre barbes et direction d'extension. Nous avons choisi les directions transverse et longitudinale pour leur simplicité d'analyse et une direction en biais à 45° dont le but est d'activer les déformations de cisaillement (voir Figure 2). Nous avons testé 10 échantillons différents pour la direction transverse, 9 échantillons pour la direction longitudinale et 8 échantillons pour la direction en biais.

De ces tests de traction, nous tirons des mesures de contraintes de Green-Lagrange E et du second tenseur de déformation de Piola-Kirchhoff S (voir Figure 7). On peut noter que la raideur du test longitudinal est quatre ordres de grandeur plus grande que celle des autres tests, montrant bien le comportement fortement anisotrope de l'étendard. De plus, dans le test transverse le comportement semble linéaire malgré les grandes déformations, cela nous amène à choisir une relation linéaire entre déformation et contrainte pour notre modèle. Les résultats du scénario en biais présentent de grandes différences : cette disparité est discutée dans la section 5.

## 3 Modèle

Afin de décrire notre modèle de membrane pour l'étendard, nous fixons le référentiel de l'espace de référence de sorte que l'axe x soit colinéaire aux barbes et que la base soit orthonormée. Dans la section précédente, nous avons identifié une relation linéaire entre contrainte et déformation dans les résultats expérimentaux. Ainsi, en prenant les mêmes notations que dans la section précédente, cette relation
<figure><img src="image_3.png" /><figcaption>FIGURE 2 – Haut : Plume de cygne et exemples de découpes tracées en pointillé pour chaque orientation d'échantillon testé. Bas gauche : Photo du montage expérimental. Bas droite : Gros plan sur le test de traction.</figcaption></figure>

linéaire entre contrainte et déformation en notation de Voigt s'écrit

$$ \begin{bmatrix} S_{xx} \\ S_{yy} \\ S_{xy} \end{bmatrix} = \underbrace{\begin{bmatrix} C_{11} & C_{12} & C_{13} \\ C_{12} & C_{22} & C_{23} \\ C_{13} & C_{23} & C_{33} \end{bmatrix}}_{c} \begin{bmatrix} E_{xx} \\ E_{yy} \\ E_{xy} \end{bmatrix}. $$

Ici, le tenseur d'élasticité C est symétrique pour des raisons de conservation d'énergie. Les paramètres de ce tenseur peuvent être réduits en utilisant des arguments géométriques. L'étendard a un axe de symétrie local en tous points aligné sur les barbes (voir Figure 3), ce qui induit un modèle orthotrope, c'est-à-dire $C_{13} = C_{23} = 0$. De plus, nous avons pu voir que les échantillons ne subissent aucun changement de largeur dans les tests transversaux et longitudinaux, ainsi il n'y a aucun phénomène de Poisson et on a $C_{12} = 0$. Avec cela, le modèle est réduit à la relation

$$ \begin{bmatrix} S_{xx} \\ S_{yy} \\ S_{xy} \end{bmatrix} = \begin{bmatrix} C_{xx} & 0 & 0 \\ 0 & C_{yy} & 0 \\ 0 & 0 & C_{xy} \end{bmatrix} \begin{bmatrix} E_{xx} \\ E_{yy} \\ E_{xy} \end{bmatrix}, $$

où les paramètre du tenseur d'élasticité ont été renommé pour mettre en valeur les quantités qu'ils relient. Aussi, nous dénoterons parfois $C_{xx}$ par module longitudinal, $C_{yy}$ par module transverse et $C_{xy}$ par module de cisaillement. Grâce aux résultats expérimentaux on peut déduire que le rapport entre le module longitudinal et le module transverse $\frac{C_{xx}}{C_{yy}}$ est de l'ordre de $10^4$, ce qui indique bien l'aspect fortement anisotrope des pennes. Dans la suite du résumé, nous appelons le rapport $\frac{C_{xx}}{C_{yy}}$ rapport de raideur.

## 4 Verrouillage

Nous implémentons notre modèle numériquement grâce au simulateur libre de plaques minces élastiques qui nous est le plus familier ARCSim [5] que nous modifions pour utiliser Ipopt comme solveur non-linéaire [8]. Ce simulateur ne supporte actuellement que des éléments triangulaires linéaires, l'impact du type d'éléments utilisé est gardé pour des travaux futurs. La modélisation de notre modèle dans le cas fortement anisotrope présente deux désagréments. Un de ces ennuis est la présence d'un phénomène de verrouillage numérique, c'est-à-dire que plus le rapport de raideur est élevé, plus la vitesse de convergence du modèle d'élément fini est faible, voire la convergence devient impossible. L'autre souci est le mauvais conditionnement du problème, qui impacte considérablement la vitesse de convergence du solveur non-linéaire.
<figure><img src="image_4.png" /><figcaption>FIGURE 3 – En tous points, l'étendard admet un axe de symétrie colinéaire aux barbes (-----).</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>FIGURE 4 – Maillage utilisé pour une simulation du scénario expérimental en biais. L'alignement des éléments sur l'orientation des barbes réduit considérablement le phénomène de verrouillage numérique.</figcaption></figure>

Pour atténuer le phénomène de verrouillage, nous nous inspirons des travaux de YU, CARTWRIGHT et al. [10] et TEN THIJE et AKKERMAN [7], dans lesquels ils et elles considèrent un matériau ayant deux directions de très grande raideur comparée au module de cisaillement associé à ces directions. Une simulation naïve d'un tel matériau présente un phénomène de verrouillage. Une de leurs solutions pour pallier cela est d'aligner les éléments selon les directions de forte raideur. Nous la reproduisons dans notre contexte en faisant en sorte que chaque élément triangulaire ait au moins une arête alignée selon la direction des barbes (voir Figure 4). Comme attendu, cet alignement limite considérablement le phénomène de verrouillage numérique (voir Figure 8).

Malgré son atténuation du phénomène de verrouillage, la simulation avec le maillage aligné présente un désagrément. Pour de grands rapports de raideur, le problème est mal conditionné, demandant donc un grand nombre d'itérations du solveur non-linéaire pour être résolu correctement (voir Figure 9). Nous remédions à cela en remplaçant la partie de l'énergie liée à la direction de grande raideur par une contrainte prévenant l'extension dans cette direction. Cette modification ne dégrade pas la convergence des éléments finis (voir Figure 8) et permet de résoudre le problème en un temps réduit comparé à la méthode non contrainte (voir Figure 9).

# 5 Validation

Notre modèle a trois paramètres indépendants $C_{xx}$, $C_{yy}$ et $C_{xy}$. Les deux premiers paramètres peuvent être mesurés à travers les tests de tractions transverses et longitudinaux<fnref n="1" />. Le troisième $C_{xy}$ nécessite un scénario en biais. Ainsi, la calibration seule du modèle demanderait d'utiliser nos trois scénarios de laboratoire, ce qui nous contraindrait à produire une nouvelle expérience afin de valider le modèle.

Pour éviter l'exécution laborieuse de nouvelles expériences en laboratoires, nous choisissons de ne pas valider notre modèle directement sur les expériences. Au lieu de cela, nous introduisons un modèle de membrane plus basse échelle pouvant être validé avec seulement nos trois scénarios. Dans le reste du document, nous appellerons ce modèle intermédiaire *modèle micro* et notre modèle d'étendard présenté plus tôt *modèle macro*. Ainsi après avoir validé le modèle micro, nous l'utilisons pour produire en simulation les scénarios nous manquant pour la validation du modèle macro.

Le modèle micro pour l'étendard se trouve à l'échelle des barbes. C'est un modèle de membrane hétérogène fait de bandes raides représentant les barbes connectées ensemble par des bandes souples représentant les barbules. Les matériaux de ces bandes sont décrits par le modèle de Saint Venant-Kirchhoff isotrope dont le coefficient de Poisson est nul. Ainsi, le modèle a seulement deux paramètres physiques, le module de Young des barbes $E_B$ et le module de Young des barbules $E_b$ (voir Figure 10).

Ce modèle n'ayant que deux paramètres physiques, seulement deux des scénarios sont nécessaires pour les calibrer. Les autres paramètres étant géométriques, ils peuvent être mesurés sur les images des échantillons prises au début de chaque expérience. Les mesures des modules $E_B$ et $E_b$ sont données dans le tableau 1.

<footnote n="1">1. En pratique, comme le paramètre $C_{xx}$ est remplacé par une contrainte (voir Section 4), il n'a pas besoin d'être mesuré.</footnote>
<table><caption>TABLE 1 – Mesure des paramètres physiques du modèle micro et macro.</caption><thead><tr><th>E<sub>b</sub> [Nm<sup>-2</sup>]</th><th>E<sub>B</sub> [Nm<sup>-2</sup>]</th><th>C<sub>yy</sub> [Nm<sup>-2</sup>]</th><th>C<sub>xy</sub> [Nm<sup>-2</sup>]</th></tr></thead><tbody><tr><td>3.7 &times; 10<sup>4</sup></td><td>2.3 &times; 10<sup>9</sup></td><td>4.5 &times; 10<sup>4</sup></td><td>3.1 &times; 10<sup>4</sup></td></tr></tbody></table>

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

FIGURE 5 – Photographies de deux tests de traction. L'expérience montrée dans la sous-figure 5a présente un comportement non élastique, à l'état au repos (gauche) les barbules du bord sont alignée, après extension (droite) ces barbules se réarrange. Ce réarrangement des barbules permet aux barbes de glisser sans induire de cisaillement important et ainsi réduire la contrainte sur le matériau. La sous-figure 5b montre une expérience où ce comportement non élastique n'est pas présent.

Les prédictions de la simulation sont comparées avec le scénario en biais. Dans ce scénario le comportement dépend fortement de l'aspect de forme de l'échantillon utilisé. Pour rendre compte de cette dépendance, des échantillons ayant des rapport de forme variant de 1.0 à 3.82 ont été utilisés (voir Figure 11). Les calculs en simulation prédissent une contrainte légèrement plus élevée que la mesure expérimentale. Cela peut être expliqué en partie par des événements non-élastiques de réarrangement des barbules ayant lieu dans ces expériences (voir Figure 5). Notre modèle étant complètement élastique, il ne peut rendre compte de ces événements dissipateurs d'énergie. Dans certains scénarios en biais exécutés en laboratoire, aucun réarrangement des barbules n'a pu être identifié visuellement : pour ces scénarios la simulation prédit un comportement bien plus proche de l'expérience.

Le modèle micro étant validé, nous considérons qu'il peut être utilisé pour valider le modèle macro. Nous reproduisons en simulation avec le modèle micro les tests de traction sur un échantillon de rapport de forme 1.38 avec différentes orientations de barbes allant de 0° à 0.53°. Les paramètres $C_{yy}$ et $C_{xy}$ du modèle macro sont calibrés sur les orientations 0° et 45° (voir Tableau 1). Les résultats de la validation montrent un bon accord entre le modèle micro et macro (voir Figure 6).

## 6 Conclusion

Dans ce résumé étendu, nous avons décrit notre montage expérimental qui nous a permis de suggérer un modèle que nous pensons être approprié pour l'étendard des pennes. Nous implémentons ce modèle grâce à des éléments finis linéaires. Pour pallier le phénomène de verrouillage numérique et de mauvais

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

FIGURE 6 – Comparaison entre les pentes des courbes de contraintes-déformations du modèle micro et du modèle macro dans des scénarios de traction sur un échantillon de rapport de forme 1.38 avec des orientations de barbes à différents angles. Les pentes sont adimentionnées en les divisant par le module de Young des barbules dans le modèle micro $E_b$.
conditionnement dus à la forte anisotropie de notre modèle, nous proposons d'aligner les éléments le long des barbes et de remplacer la partie raide de l'énergie élastique par une contrainte d'inextensibilité. De plus, nous avons validé notre modèle indirectement, en utilisant un modèle intermédiaire plus bas niveau et validé expérimentalement, et ce pour réduire le nombre de scénarios différent à réaliser en laboratoire.

Références

[1] K. BHAR, B. CHANG et al., *How localized force spreads on elastic contour feathers*, Journal of The Royal Society Interface, p. 20190267, 2019.

[2] F. BREZZI et M. FORTIN, *Mixed and Hybrid Finite Element Methods*, 1991.

[3] Q. CHEN, S. GORB et al., *An analytical hierarchical model explaining the robustness and flaw-tolerance of the interlocking barb-barbule structure of bird feathers*, EPL (Europhysics Letters), p. 24001, 2016.

[4] A. KOVALEV, A. E. FILIPPOV et S. N. GORB, *Unzipping bird feathers*, Journal of The Royal Society Interface, 2014.

[5] T. PFAFF, R. NARAIN et al., *Adaptive Tearing and Cracking of Thin Sheets*, ACM Transactions on Graphics, xx:1-9, juil. 2014.

[6] J. SCHRÖDER, N. VIEBAHN et al., *A novel mixed finite element for finite anisotropic elasticity; the SKA-element Simplified Kinematics for Anisotropy*, Computer Methods in Applied Mechanics and Engineering, p. 475-494, 2016.

[7] R. TEN THIJE et R. AKKERMAN, *Solutions to intra-ply shear locking in finite element analyses of fibre reinforced materials*, Composites Part A: Applied Science and Manufacturing, p. 1167-1176, 2008.

[8] A. WÄCHTER et L. T. BIEGLER, *On the implementation of an interior-point filter line-search algorithm for large-scale nonlinear programming*, Mathematical Programming, p. 25-57, 2006.

[9] P. WRIGGERS, B. HUDOBIVNIK et J. SCHRÖDER, *Finite and Virtual Element Formulations for Large Strain Anisotropic Material with Inextensive Fibers*, 2018.

[10] X. YU, B. CARTWRIGHT et al., *Intra-ply shear locking in finite element analyses of woven fabric forming processes*, Composites Part A: Applied Science and Manufacturing, p. 790-803, 2006.

[11] F. ZHANG, L. JIANG et S. WANG, *Repairable cascaded slide-lock system endows bird feathers with tear-resistance and superdurability*, Proceedings of the National Academy of Sciences, p. 10046-10051, 2018.

[12] J.-S. ZHAO, J. ZHANG et al., *Shaking the wings and preening feathers with the beak help a bird to recover its ruffled feather vane*, Materials & Design, p. 108410, 2020.
<figure><img src="image_8.png" /><figcaption>FIGURE 7 – Résultats expérimentaux pour chaque orientation (haut) et dimensions des échantillons utilisés (bas) $l_x$ et $l_y$ sont données en millimètre et $h$ en micromètre. Pour la direction transverse, les pentes des courbes varient entre 37 kPa et 61 kPa, et entre 263 MPa et 457 MPa pour la direction longitudinale.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>FIGURE 8 – Contrainte $S$ pour une déformation de 5 % en fonction du nombre d'éléments dans une simulation du scénario en biais. Dans ces simulations le rapport de raideur est fixé) $C_{xx} = 10^4$. Les courbes de convergence sont données pour une simulation utilisant un maillage arbitraire (●) et pour une simulation utilisant un maillage aligné. Les courbes pour la simulation avec le maillage aligné sont données pour le cas d'une simulation sans contrainte (■) et d'une simulation sous contrainte (■). La valeur de référence affichée (---) est celle obtenue avec la simulation sans contrainte utilisant le maillage aligné constituée de $128 \times 10^3$ éléments. La simulation utilisant un maillage non-aligné a une vitesse de convergence plus faible que les simulations utilisant un maillage aligné. Cela met en valeur le phénomène de verrouillage dans le cadre de cette simulation et montre que l'utilisation d'un maillage aligné pallie ce phénomène.</figcaption></figure>
<figure><img src="image_10.png" /></figure>

FIGURE 9 – Courbes contraintes-déformations pour le scénario en biais calculé avec un maillage aligné. Les calculs sont faits avec et sans contrainte. Les calculs non contraints sont faits avec 50 (—) et 100 (---) itérations du solveur non-linéaire. Ceux contraints utilisent au plus 50 itérations du solveur non-linéaire (---). À cause du mauvais conditionnement du problème, les résultats obtenus sans utiliser la contrainte oscillent. Utiliser une contrainte permet de calculer une solution sans oscillation à bien moindre coût.

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

FIGURE 10 – Schéma du modèle intermédiaire micro utilisé pour la validation de notre modèle d'étendard. Ce modèle de membrane est formé d'une alternance de bandes élastiques raides et isotropes représentant les barbes (■) et de bandes élastiques souples et isotropes représentant les barbules (■) accrochant les barbes entre eux. Le modèle a quatre paramètres : deux paramètres physiques, le module d'Young des bandes représentant les barbes $E_B$, et celui des bandes représentant les barbules $E_b$ ; deux paramètres géométriques, la largeur des barbes $w_B$ et la distance entre deux barbes centre à centre $d_B$.

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

FIGURE 11 – Comparaison de la simulation et de l'expérience sur le scénario en biais. Dans certaines expériences les barbules se réarrangent, dissipant ainsi l'énergie (gauche) (voir Figure 5). Notre modèle d'étendard étant strictement élastique, il prédit une contrainte légèrement plus haute que celle observée dans l'expérience. D'autres expériences ne présentent aucun réarrangement (droite), en conséquence la simulation fournit de meilleures prédicitions.