THÈSE DE DOCTORAT

DE
L’UNIVERSITÉ D’ÉVRY-VAL D’ESSONNE

Spécialité : Mécanique

Présentée par
Zhiming ZHANG

pour obtenir le grade de

DOCTEUR DE L’UNIVERSITÉ D’ÉVRY-VAL D’ESSONNE

Sujet de la thèse :

Modélisation mécanique des interfaces multi-contacts dans une pile à combustible

Soutenance prévue le 30 novembre 2010 devant le jury composé de :

Q.-C. HE Professeur, Université Paris-Est, Rapporteur

M. HJIAJ Professeur, INSA de Rennes, Rapporteur
C. VALLÉE Professeur, Université de Poitiers, Examinateur
Z.-Q. FENG Professeur, Université d’Évry - Val d’Essonne, Directeur de thèse
C. RENAUD Maître de Conférences, Université d’Évry - Val d’Essonne, Co-directeur de thèse

—————————————————————————————————–
L ABORATOIRE DE M ÉCANIQUE ET D ’É NERGÉTIQUE D ’É VRY (LMEE, EA 3332)
40, rue du Pelvoux, 91020 Évry Cedex
 i

Résumé :

La pile à combustible transforme l’énergie chimique en énergie électrique grâce à un

empilement de différentes structures lamellaires. Des phénomènes mécaniques présents aux
interfaces agissent de manière plus ou moins directe sur les performances et la durée de vie
de la pile. Nous avons montré que la précompression par boulons tirants pendant l’assem-
blage, la déformation de la couche de diffusion gazeuse (GDL), le contact et la configuration
de la plaque bipolaire (BPP) avaient des influences importantes sur la résistance de contact,
la porosité et la perméabilité de la pile. La résistance de contact est déterminée par la zone
de contact et la pression de contact. La porosité et la perméabilité sont liées à la déformation
aux interfaces. Le contact entre les différentes structures a un rôle majeur dans le fonction-
nement de la pile. Ce problème a été modélisé par la méthode des éléments finis. Différents
paramètres de la pile comme la précompression, la structure géométrique des dents de la
plaque BPP ou encore la porosité de la GDL ont été étudiés et ont permis de connaître l’état
de contact ainsi que les déformations dans les structures. L’influence de certains paramètres
sur les résultats mécaniques de la pile a ensuite été abordée. L’objectif de cette étude est
de fournir des valeurs optimales de ces paramètres pour obtenir la meilleure performance
possible de la pile.
Mots clés : pile à combustible - phénomènes mécaniques - contact - éléments finis

Abstract :
The fuel cell transforms chemical energy to electrical power sources through a stack of dif-
ferent planar structures. Mechanical phenomena presented on the multi-contact interface acts
more or less the fuel cell’s performance and lifetime. We have shown that the pre-load by
stacking bolts, the deformation of the gas diffusion layer (GDL), the contact and the confi-
guration of the bipolar plate (BPP) had influences on the contact resistance, the porosity and
the permeability of the fuel cell. The contact resistance is determined by the contact area and
contact pressure. The porosity and permeability are related to the interfacial deformation.
The contact between the different structures has a major role in the fuel cell operation. This
problem is solved by the finite element method. Various parameters of the fuel cell as the
pre-load, the geometric structure of teeth of BPP as well as the porosity of the GDL were
studied and allowed to know the contact behavior and the deformation. The influence of
some parameters on the mechanical results of fuel cell stack was then tackled. The purpose
of this study is to provide optimum values of these parameters to obtain the best performance
of fuel cells.
Key words : fuel cell - mechanical phenomena - contact - finite element
ii
 iii

Table des matières

Avant-propos i

Résumé i

Abstract i

Table des matières iii

Introduction générale 1
Cadre du problème . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
Objectif de l’étude . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
Plan du document . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1

1 Présentation de la pile à combustible 7

1.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
1.2 Problématique énergétique et hydrogène . . . . . . . . . . . . . . . . . . . 7
1.2.1 Problématique énergétique . . . . . . . . . . . . . . . . . . . . . . 7
1.2.2 Technologie d’hydrogène . . . . . . . . . . . . . . . . . . . . . . . 8
1.3 Description générale d’une pile à combustible . . . . . . . . . . . . . . . . 9
1.3.1 Aperçu historique de la pile . . . . . . . . . . . . . . . . . . . . . 9
1.3.2 Types de piles à combustible . . . . . . . . . . . . . . . . . . . . . 12
1.4 Piles à combustible de type PEM . . . . . . . . . . . . . . . . . . . . . . . 13
1.4.1 Principe de fonctionnement de la pile . . . . . . . . . . . . . . . . 13
1.4.2 Performance actuelle de la pile . . . . . . . . . . . . . . . . . . . . 15
1.4.3 Structure de la pile à combustible . . . . . . . . . . . . . . . . . . 17
1.4.4 Assemblage de la pile à combustible . . . . . . . . . . . . . . . . . 18
1.4.5 Le processus opératoire de la PAC . . . . . . . . . . . . . . . . . . 18
1.4.6 Composants de la PAC, leurs rôles et leurs caractéristiques . . . . . 21
1.4.6.1 Assemblage membrane des électrodes (MEA) . . . . . . 21
1.4.6.2 Couches de diffusion gazeuse (GDL) . . . . . . . . . . . 22
1.4.6.3 Plaques bipolaires (BPP) . . . . . . . . . . . . . . . . . 23
iv TABLE DES MATIÈRES

1.4.6.4 Joints d’étanchéité . . . . . . . . . . . . . . . . . . . . . 25

1.5 Étude mécanique de la structure avec des interfaces multi-contacts . . . . . 26
1.5.1 Force de serrage . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
1.5.2 Porosité de la GDL . . . . . . . . . . . . . . . . . . . . . . . . . . 28
1.5.3 Optimisation des composants . . . . . . . . . . . . . . . . . . . . 32
1.6 Interactions des effets mécaniques aux interfaces d’une PAC . . . . . . . . 33
1.6.1 Forces de serrage . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
1.6.2 Déformation de la GDL . . . . . . . . . . . . . . . . . . . . . . . 34
1.6.3 Contact . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
1.6.4 Endommagement mécanique . . . . . . . . . . . . . . . . . . . . . 35
1.6.5 Pression fluidique . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
1.7 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36

2 Méthodes numériques et résolution du problème de contact 37

2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37
2.2 Méthode des éléments finis (MEF) . . . . . . . . . . . . . . . . . . . . . . 37
2.2.1 Principe . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
2.2.2 Maillage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
2.2.3 Fonctions d’interpolation . . . . . . . . . . . . . . . . . . . . . . . 39
2.2.4 Assemblage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
2.2.4.1 Élasticité linéaire . . . . . . . . . . . . . . . . . . . . . 41
2.2.4.2 Élasticité non linéaire . . . . . . . . . . . . . . . . . . . 44
2.3 Cinématique du contact . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
2.3.1 Notion du problème de contact . . . . . . . . . . . . . . . . . . . . 46
2.3.2 Conditions de contact de Signorini . . . . . . . . . . . . . . . . . . 47
2.3.3 Loi de frottement - modèle de Coulomb . . . . . . . . . . . . . . . 49
2.3.4 Loi de contact avec frottement . . . . . . . . . . . . . . . . . . . . 51
2.4 Méthodes de résolution des problèmes de contact . . . . . . . . . . . . . . 52
2.4.1 Formulation par des fonctions de pénalité . . . . . . . . . . . . . . 53
2.4.1.1 Formulation implicite . . . . . . . . . . . . . . . . . . . 53
2.4.1.2 Formulation explicite . . . . . . . . . . . . . . . . . . . 54
2.4.2 Formulations par des multiplicateurs de Lagrange . . . . . . . . . . 55
2.4.2.1 Formulation par adhérence . . . . . . . . . . . . . . . . 55
2.4.2.2 Formulation par glissement sans frottement . . . . . . . . 55
2.4.3 Formulations lagrangiennes augmentées . . . . . . . . . . . . . . . 56
2.4.3.1 Pseudo-potentiels . . . . . . . . . . . . . . . . . . . . . 56
2.4.3.2 Bi-potentiel . . . . . . . . . . . . . . . . . . . . . . . . 57
2.5 Application numérique . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
 v

2.6 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64

3 Étude des effets mécaniques dans une cellule de la pile 65

3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
3.2 Simplification et propriétés du modèle numérique . . . . . . . . . . . . . . 65
3.2.1 Description et simplification du modèle . . . . . . . . . . . . . . . 65
3.2.2 Propriétés mécaniques du modèle . . . . . . . . . . . . . . . . . . 66
3.3 Modèle des éléments finis de la PAC . . . . . . . . . . . . . . . . . . . . . 67
3.3.1 Description du modèle . . . . . . . . . . . . . . . . . . . . . . . . 67
3.3.2 Maillage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
3.3.3 Définition de l’objet “zones de contact” . . . . . . . . . . . . . . . 70
3.3.4 Conditions aux limites . . . . . . . . . . . . . . . . . . . . . . . . 72
3.3.5 Chargement de précompression de la PAC . . . . . . . . . . . . . . 72
3.4 Études des contacts dans la PAC . . . . . . . . . . . . . . . . . . . . . . . 73
3.4.1 Étude de deux zones de contacts (GDL/BPP et GDL/MEA) . . . . 73
3.4.2 Étude d’une zone de contact (GDL/BPP) avec l’angle droit . . . . . 75
3.4.3 Modélisation avec congé au niveau de la dent . . . . . . . . . . . . 76
3.4.3.1 Contrainte dans les plaques . . . . . . . . . . . . . . . . 77
3.4.3.2 Pression de contact . . . . . . . . . . . . . . . . . . . . 80
3.4.3.3 Aire de contact (longueur de contact) . . . . . . . . . . . 80
3.4.3.4 Déformation de la plaque GDL . . . . . . . . . . . . . . 82
3.4.3.5 Déformation de la plaque GDL avec pores . . . . . . . . 84
3.5 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88

4 Influence des principaux paramètres sur la performance des piles 89

4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89
4.2 Précompression de la PAC . . . . . . . . . . . . . . . . . . . . . . . . . . 89
4.2.1 Pression de contact . . . . . . . . . . . . . . . . . . . . . . . . . . 90
4.2.2 Aire de contact (longueur de contact) . . . . . . . . . . . . . . . . 92
4.2.3 Porosité de la plaque GDL . . . . . . . . . . . . . . . . . . . . . . 94
4.3 Rayon du congé des dents de la plaque BPP . . . . . . . . . . . . . . . . . 97
4.3.1 Pression de contact . . . . . . . . . . . . . . . . . . . . . . . . . . 97
4.3.2 Aire de contact (longueur de contact) . . . . . . . . . . . . . . . . 100
4.3.3 Porosité de la plaque GDL . . . . . . . . . . . . . . . . . . . . . . 101
4.4 Largeur des dents de la plaque BPP . . . . . . . . . . . . . . . . . . . . . 103
4.5 Épaisseur de la plaque GDL . . . . . . . . . . . . . . . . . . . . . . . . . 107
4.6 Épaisseur de la membrane MEA . . . . . . . . . . . . . . . . . . . . . . . 108
4.7 Étude de la géométrie à l’interface . . . . . . . . . . . . . . . . . . . . . . 110
vi TABLE DES MATIÈRES

4.7.1 Influence des aspérités de l’interface . . . . . . . . . . . . . . . . . 112

4.7.2 Influence des hauteurs de dents de la BPP . . . . . . . . . . . . . . 116
4.8 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119

Conclusions et perspectives 119

Références bibliographiques 125

Liste des figures 135

Liste des tableaux 138

 Introduction générale 1

Introduction générale

Cadre du problème :
De nos jours, avec le développement de l’industrie, nous consommons de plus en plus de
ressources en énergies, comme le charbon, le pétrole et le gaz naturel, etc.. Les limites de
ces ressources deviennent un problème auquel le public s’intéresse. Les questions énergé-
tiques recouvrent à l’heure actuelle deux domaines. L’un est lié au risque d’épuisement des
ressources fossiles et fissiles, l’autre est le problème de la pollution environnementale, tels
que le réchauffement par effet de serre et les pluie acides.
Face à la diminution des ressources énergétiques et à l’augmentation de la pollution des
énergies, il est impératif de diminuer la pollution dans le monde (un premier pas a été le
Protocole de Kyoto signé en 1997), et de rechercher des alternatives énergétiques possédant
les mêmes propriétés que les hydrocarbures. Dans ce contexte, l’hydrogène est une solution
sérieuse et intéressante en tant que source d’énergie durable et renouvelable.
L’hydrogène, qui n’existe pas à l’état naturel, peut en effet être synthétisé à partir des
énergies renouvelables. Il est particulièrement intéressant car son stockage est possible sous
différents états (gaz, liquide et solide). D’autre part, il est léger et "propre". Mais il ne peut
pas être utilisé directement, et nécessite la transformation d’une énergie chimique en énergie
électrique. Cette production d’énergie a lieu dans une pile à combustible.
La pile à combustible est un générateur efficace d’énergie électrique à partir de l’éner-
gie chimique, elle utilise la réaction électrochimique d’oxydoréduction entre l’hydrogène et
l’oxygène pour donner de l’électricité. En fait, des recherches sur les piles à combustibles
ont été faites depuis une centaine d’années. Les avantages offerts par la technologie des piles
à combustibles comme l’efficacité, la protection de l’environnement et la construction mo-
dulaire, sont tels qu’elles ont suscité un intérêt permanent pour leur développement depuis
1893.
Après 170 ans de recherches, différentes piles ont déjà été développées en particulier, la
pile à combustible de type PEM (Membrane Échangeuse des Protons), qui peut être em-
ployée dans une position fixe ou embarquée (dans le transport comme par exemple le véhi-
cule électrique). La stabilité des composants comme la nouvelle membrane, la plaque bipo-
laire, le flux et la température nominaux ont été étudiés, mais ce n’est pas encore suffisant
2 Introduction générale

pour la commercialisation à grande échelle de la pile. En effet, les influences mécaniques sur
la performance des piles sont reconnues et les interactions entre les effets mécanique et la
performance des piles doivent être étudiés. Dans ces travaux, nous nous sommes intéressés
aux interactions des effets mécaniques au sein de la pile à combustible sur la performance et
la durée de vie des piles.
Dans cette thèse, nous entendrons par :
– Empilement de la pile à combustible (PAC) : un empilement de cellules électrochi-
miques sans les auxiliaires (système de refroidissement, d’alimentation en gaz).
– Système PAC : une pile à combustible avec ses auxiliaires (système de refroidissement,
d’alimentation en gaz) et ses organes de contrôle et de commande.
Nous considérerons principalement le cas d’un empilement de la pile. Nous ne nous inté-
resserons jamais au système entier de la pile à combustible.

Objectifs de l’étude :
Ce travail a été réalisé au sein du laboratoire LMÉÉ (Laboratoire de Mécanique et d’Éner-
gétique d’Évry). La thèse est consacrée aux travaux de modélisation des effets mécaniques
sur les interfaces de contact en vue d’étudier ses interactions sur la performance des piles à
combustibles de type de PEM et sur leur durée de vie.
Le principe de fonctionnement de la pile est une réaction électrochimique inverse de
l’électrolyse à partir d’hydrogène et d’oxygène. En fait, la pile à combustible qui n’est pas
parfaitement utilisée, a des pertes de performance et des dégradations irréversibles du cœur
de la pile qui font chuter la durée de vie des piles. La perte de la performance est surtout dûe
à la résistance de contact de la pile à combustible causée par la structure lamellaire. La pile
à combustible est constituée d’un empilement d’éléments comme le MEA (Assemblage des
membranes électrodes) et la GDL (Couche de Diffusion Gazeuse), qui sont pris en sandwich
entre des plaques bipolaires (BPP) formant ainsi une “cellule”. Ces cellules sont reliés en
série électrique pour former un “empilement”. Pour un fonctionnement normal de la pile, la
force de serrage par boulons tirants constitue le chargement qui maintient les éléments en
contact. Ainsi la pile à combustible peut être regardée comme une structure lamellaire que
nous appellerons structure multi-contact.
Le contact est un phénomène non linéaire, complexe et s’appliquant sur des zones pas
toujours faciles à déterminer, sa prise en compte entre les éléments de la PAC est pourtant
importante pour étudier le comportement mécanique aux interfaces. Il est évident que les
effets mécaniques sur les interfaces comme la force de contact joue un rôle fondamental
dans le comportement de la structure : non seulement la surface de contact et l’état de contact,
mais aussi sa déformation, son déplacement, la distribution des efforts doivent être connus
afin d’améliorer la performance de la pile.
 Introduction générale 3

Lorsqu’un empilement de la pile est serré par boulons tirants, la force de serrage induit
des efforts normaux et tangentiels de contact aux différentes interfaces. La force de contact,
la déformation et la zone de contact sont liées entre elles. La déformation des éléments est
liée à leur porosité et à leur perméabilité aux gaz réactifs, la zone de contact détermine donc
la voie de transfert des électrons. Par ailleurs, une force de contact très élevée peut provoquer
la fissure mécanique. Ainsi une étude de ces interactions est pertinente afin de caractériser
leur influence sur la performance de la pile et sa durée de vie. Sur les interfaces, plus la
force de contact est élevée, et l’aire de contact grande, plus la résistance de contact est petite
et meilleure est la performance de la PAC. En revanche, plus la porosité et la perméabi-
lité de GDL sont petites (sous de fortes déformations) et plus le risque d’endommagement
mécanique sur MEA et GDL est réel car ils sont très minces et fragiles. Comme la pile à
combustible est un empilement multi-contacts en série électrique, l’endommagement d’un
élément de la pile peut conduire au dysfonctionnement de tout l’empilement. Ainsi le com-
promis entre une perméabilité élevée et une petite résistance de contact est un problème clé
dans la recherche d’une puissance optimale de la pile.
De par la configuration rectangulaire des petits dents de la BPP, des concentrations de
contraintes peuvent apparaître aux angles droits des dents, la distribution de la force de
contact sur les interfaces entre GDL et les petits dents de BPP n’est en effet pas uniforme. Par
conséquent, la distribution de la perméabilité et de la résistance de contact sur les interfaces
ne sont pas uniformes, ce qui peut provoquer des problèmes : d’une part, une distribution
non uniforme de l’alimentation des gaz réactifs, d’autre part, une distribution non uniforme
de la résistance de contact entre des éléments de la PAC quand il faut obtenir une réaction
électrochimique la plus régulière possible dans la PAC.
Ainsi, l’influence de ces paramètres interagissant au niveau des interfaces multi-contacts
sur la performance de la pile à combustible doit être étudiée.
L’étude expérimentale des phénomènes de contact étant complexe, une simulation numé-
rique est plus simple et appropriée. Une modélisation est proposée dans nos études. Elle est
basée sur la méthode des éléments finis et prend en considération les problèmes de contact
qui ne peuvent pas être résolus analytiquement. Les objectifs de la thèse sont les suivants :
– Identifier les effets mécaniques et leurs interactions sur les interfaces multi-contacts
liées à la performance des piles à combustibles, déterminer des paramètres importants
qui peuvent influencer la performance de la pile ;
– Modéliser une cellule de la pile avec ses différents contact en utilisant la méthode des
éléments finis afin d’obtenir les résultats précis ;
– Développer un modèle paramétrique permettant d’optimiser l’état de contact méca-
nique sur les interfaces (distribution des contraintes, distribution de la pression de
contact, zone de contact et déformation) afin d’étudier leurs influences sur la perfor-
mance de la pile.
4 Introduction générale

Avant d’aborder les paramètres importants dans un empilement de pile, un travail sur la
modélisation et sur la caractérisation des éléments mécaniques de la pile à combustible s’est
avéré nécessaire.
La pile à combustible est un système multi-physique : mécanique, électrochimique, ther-
mique et fluidique. Notre approche se limite à l’étude des phénomènes mécaniques au sein
d’un empilement de cellules soumis à un chargement purement mécanique, sans considérer
les aspects thermique et fluidique. Nous avons retenu une approche simple sans considérer
les conditions opératoires comme la pression des gaz réactifs, la température et l’hydratation
des membranes. Dans cette approche, nous abordons la modélisation des états de contact
mécaniques au sein de la pile à combustibles de type PEM pour connaître leur impact sur la
performance et la durée de vie.

Plan du document :
Dans le cadre de la modélisation d’un empilement de la pile, les modèles développés sont
bidimensionnels et les matériaux élastiques isotropes dans l’étude.
Le premier chapitre présente les objectifs et les enjeux de la thèse. Nous commençons
par la description d’un état problématique de nature énergétique et environnementale dans le
monde. L’énergie d’hydrogène est une solution prometteuse et favorable à l’environnement
qui permet d’obtenir de l’électricité grâce à la pile à combustible. Nous étudierons ensuite le
principe, la structure et la performance d’une pile à combustible. Parce que la pile à combus-
tible de type PEM est une structure lamellaire, l’état de contact mécanique entre éléments
de la pile a un rôle important sur sa performance. Les effets mécaniques comme le contact,
la déformation sont décrits en détail, leurs interactions sur la performance de la pile sont
discutées. La modélisation numérique est un moyen d’approfondir les connaissances sur les
phénomènes mécaniques aux interfaces multi-contacts des piles à combustibles.
Le second chapitre présente l’approximation d’un problème de contact avec frottement
par la méthode des éléments finis. On commencera par présenter la méthode des éléments
finis, puis nous aborderons la formulation des problèmes de contact unilatéral avec frotte-
ment de Coulomb. Nous présenterons ensuite différentes méthodes de résolution numérique
de ces problèmes, en donnant, pour terminer, une application numérique d’un cas test.
Dans le troisième chapitre, les phénomènes mécaniques au niveau des interfaces entre
éléments de la pile ont été étudiés. Afin de simplifier la modélisation, nous nous sommes
limités à regarder les effets mécaniques pour une seule cellule de l’empilement soumise à un
chargement statique. Après avoir décrit le modèle éléments finis utilisé (géométrie, maillage,
conditions aux limites, etc.), nous avons présenté les résultats obtenus sur les états de contact
(zones de contact, pression de contact) et sur les déformations. Un modèle prenant en compte
la porosité de la GDL est également abordé.
 Introduction générale 5

La quatrième et dernière partie porte sur l’influence de certains paramètres sur la perfor-
mance de la pile. Les résultats obtenus en faisant varier des paramètres comme la charge
de précompression, le rayon du congé et la largeur des dents de la BPP, les épaisseurs des
éléments GDL et MEA ou encore des aspérités de l’interface sont présentés. L’objectif est
ici de dégager des valeurs optimales de ces paramètres pour la performance de la pile.
Nous terminerons la thèse par une conclusion générale qui synthétisera les résultats obte-
nus et présentera des perspectives de recherche.
6 Introduction générale
 7

Chapitre 1

Présentation de la pile à combustible

1.1 Introduction
Nous débutons le chapitre par une description de la problématique énergétique et environ-
nementale. L’énergie provenant de l’hydrogène est un des meilleurs choix en matière d’éner-
gie propre à condition de savoir le transformer en énergie électrique, ce que permet la pile à
combustible. Nous retracerons brièvement l’historique de la pile à combustible. Par la suite,
nous présenterons les caractéristiques intéressantes de la pile à combustible de type PEM
(membrane échangeuse des protons) : basse température, légèreté et facilité de construc-
tion. Elle est prometteuse dans le cadre d’une application au transport (véhicule électrique).
Des problèmes existent dans la structure de la PAC comme une température non uniforme
sur les interfaces, une résistance de contact élevée et un endommagement des composants,
ce qui peut conduire à une fissure mécanique et à un manque de performance de la PAC.
L’exploration des interactions des états mécaniques concernant la performance de la pile et
sa durabilité apparaît comme une phase essentielle avant la commercialisation de la pile à
combustible. Ces interactions sont présentées en détail.

1.2 Problématique énergétique et hydrogène

1.2.1 Problématique énergétique

Avec le développement de l’industrie depuis plusieurs décennies, nous consommons de
plus en plus de ressources d’énergies pour pouvoir satisfaire les besoins croissants en énergie
de la population. Presque 80% de la production mondiale d’énergie provient du charbon, du
pétrole, du gaz naturel et le reste du nucléaire, etc. [Hoo03]. L’épuisement de ces ressources
devient être un problème qui préoccupe chacun d’entre nous, les questions énergétiques re-
couvrent à l’heure actuelle deux domaines. L’un est lié au risque d’épuisement des ressources
8 Chapitre 1 : Présentation de la pile à combustible

fossiles et fissiles, notamment celle du pétrole qui au rythme actuel de consommation, est
prévue pour dans un siècle ou moins. L’autre est le problème de la pollution environnemen-
tale, notamment l’émission de CO2 (gaz à effet de serre) et de gaz polluants (SO2, NOx, CO,
CH4, photofluorographies, particules solides, etc.), provenant des véhicules de transport ou
de l’entreprise industrielle (figure 1.1).

F IG . 1.1: Pollution environnementale

Depuis le 19ème siècle, de plus en plus de CO2 est libéré (figure 1.2), il est désormais
important de stopper cette progression et de la réduire le plus possible.

F IG . 1.2: Accroissement de CO2 depuis le 19ème siècle

1.2.2 Technologie d’hydrogène

Face à l’accroissement incessant de la consommation énergétique et aux problèmes envi-
ronnementaux qu’il soulève, il est urgent d’opérer des choix de société. Deux solutions sont
proposées. Un est de réduire la consommation globale d’énergie. Plusieurs gouvernements,
conscients de cette situation, se sont mis d’accord au Japon en 1997 et ont signé le Protocole
de Kyoto qui vise à réduire et à stabiliser les émissions du CO2 pour la période 2008-2012,
à des valeurs correspondant à 1990. Malheureusement, on connaît les limites de ces accords.
La seconde solution consiste à développer des technologies nouvelles plus favorables à
Description générale d’une pile à combustible 9

l’environnement. L’énergie d’hydrogène apparaît comme l’une des solutions les plus pro-
metteuses car elle présente plusieurs caractéristiques intéressantes [Hoo03] :
– énergie efficace : par rapport au pétrole ou au charbon, l’hydrogène de poids équivalent
libère environ 3 fois plus d’énergie que le pétrole, environs 6 fois plus d’énergie que le
charbon.
– énergie propre : production d’eau sans pollution.
– source fiable : l’hydrogène est très abondant et très accessible dans la nature.
Grâce aux nombreux avantages, l’hydrogène est une énergie très respectueuse de l’envi-
ronnement. La recherche de nouvelles technologies d’utilisation a été encouragée et entre-
prise, afin de développer des systèmes de conversion ou de production d’énergie électrique.
Les piles à combustible apparaissent comme l’une des meilleures solution pour la transfor-
mation de l’hydrogène en énergie électrique.

1.3 Description générale d’une pile à combustible

Une pile à combustible est un générateur d’énergie électrique. Elle transforme directe-
ment l’énergie chimique du combustible en énergie électrique. C’est un système qui ne pro-
duit pratiquement pas de nuisances sonores, puisqu’il ne comporte pas de composants mé-
caniques en mouvement, comme les turbines et les moteurs. De plus, le courant électrique
est produit tant que la pile est alimentée conjointement en combustible (hydrocarbures, al-
cools, biomasse, gaz naturel, hydrogène) et en comburant (oxygène de l’air). C’est ce qui
la différencie des batteries, accumulateurs et autres piles, où se trouve stockée sous forme
chimique une quantité limitée d’énergie électrique et qui doivent soit être rechargés lorsque
c’est possible (batterie de véhicule), soit être remplacés (piles pour poste de radio).

1.3.1 Aperçu historique de la pile

En 1806, Sir Humphry Davy réalisait l’électrolyse de l’eau distillée et obtenait de l’hydro-
gène et de l’oxygène en consommant de l’électricité. C’est Christian Friedrich Schoenbein
le premier qui, en 1838, observe le principe des piles à combustible. Dans son expérience,
il utilise un tube en U avec deux électrodes de platine. Grâce à un courant électrique, il par-
vient à obtenir de l’hydrogène et de l’oxygène. En coupant l’alimentation, il constate que
les gaz produisent un courant électrique en sens inverse du premier. William Robert Grove
a rencontré Schoenbein lors d’une conférence à Birmingham en 1839. Les deux hommes
sympathisèrent et se mirent au courant de leurs recherches. En 1839, Grove a démontré le
principe de la pile à combustible [Gro39].
Dans sa célèbre expérience, il utilisait des tubes en U (figure 1.3), il s’agissait d’une cel-
lule hydrogène-oxygène avec des électrodes de platine poreux et de l’acide sulfurique comme
10 Chapitre 1 : Présentation de la pile à combustible

F IG . 1.3: Pile inventée par Grove [Gro39]

électrolyte. Un courant constant circule entre les électrodes en platine dont une extrémité est
immergée dans un récipient d’acide sulfurique et l’autre dans des récipients scellés permet-
tant de recueillir l’oxygène et de l’hydrogène. Les récipients scellés ont tenu l’eau comme
les gaz, et il a noté que le niveau d’eau était monté dans les deux tubes tant que le courant
circulait. En combinant plusieurs ensembles de ces électrodes dans un circuit en série, il a
créé ce qu’il a appelé une “ Pile de gaz ” , nom donné en 1889 par Ludwig Mond et Carl
Langer [Mos03] qui ont introduit les catalyseurs (platine) et ont perfectionné l’électrolyte.
À la fin du 19ème siècle, Friedrich Wilhelm Ostwald, le physico-chimiste qui a fourni
une grande partie de la compréhension théorique du fonctionnement des piles, en étudiant la
relation entre la propriété physique et les réactions chimiques, a amélioré la pile de Grove.
En 1893, il a déterminé expérimentalement les rôles interactifs entre différentes composants
de la pile : des électrodes, des électrolytes, oxydants et ions, etc. Ses travaux en chimie ont
donné les bases de la pile à combustible aux chercheurs. Plus tard, il a également démontré
que les piles étaient plus efficaces que le moteur à combustion interne [Sri06].
Dans les années 1930, un autre scientifique célèbre, Francis Thomas Bacon (figure1.4) a
entrepris de développer un dispositif opérationnel à partir de l’expérience de Grove. Il a amé-
lioré le catalyseur platine onéreux d’une PAC hydrogène-oxygène en utilisant un électrolyte
alcaline moins corrosif et une électrode moins coûteuse en nickel, appelé pile de combus-
tible alcaline (AFC, Alkaline Fuel Cell). En 1959, Bacon et ses collègues ont abouti à la
réalisation d’une pile à combustible hydrogène-oxygène d’une puissance de 5 KWs capable
d’alimenter une machine. La température de fonctionnement allait de 40 ◦ C à 200 ◦ C. Plus
tard, la pile de Bacon a été modifiée pour le projet spatial de Gemini de 1960 afin d’obtenir
de l’électricité et de l’eau pour les voyages dans l’espace [SNCH00].
Pour que la pile à combustible soit aussi efficace que la pile de Bacon, l’électrolyte doit
être tel que les ions le traversent aussi facilement qu’un flux d’électrons. On découvre plus
tard que plusieurs cellules empilées permettent d’obtenir une puissance plus élevée.
La pile à combustible moderne avec une membrane (Électrolyte solide à membrane échan-
geuse de protons [PEM]) a permis l’application dans un milieu immobile ou mobile.
Description générale d’une pile à combustible 11

F IG . 1.4: Francis Thomas Bacon et la pile à combustible [Sri06]

La pile à combustible de type PEM est aussi appelé polymère solide pile à combustible
en raison de l’utilisation de l’électrolyte solide. La première pile à combustible de type PEM
représentée sur la figure 1.5, a été conçue par la compagnie General Electric au début des
années 1960 grâce au travail de Thomas Grubb et Leonard Niedrach [GN60].

F IG . 1.5: Pile à combustible de PEM utilisée dans le programme Gemini [GHV69]

La technologie de la pile de PEM a servi dans le cadre du Projet-Gemini début 1960. En

1970, l’entreprise Dupont a mis au point, pour un autre usage, la membrane Nafion, qui a
permis de relancer les piles à combustible de PEM. Un grand nombre de produits dérivés des
modifications de Nafion et d’autres polymères (par exemple sulfonées polyetherketones ou
SPEK) font leur apparition dans un large éventail de prototypes de la pile. Le développement
des piles à combustible a largement suivi l’histoire des membranes.
Le très fort développement des recherches sur les piles à combustible de PEM dans les
années 1970 résulte de la première crise du pétrole. En 1972, plus de 30 projets de recherche
sont menés aux États-Unis sur le stockage de l’hydrogène embarqué ou sur sa fabrication
embarquée.
12 Chapitre 1 : Présentation de la pile à combustible

Dans les années 1980, la technologie de la pile a commencé à être testée par les services
publics et les constructeurs automobiles. Des avancés techniques ont été faite avec le déve-
loppement du premier véhicule muni d’une pile à combustible de PEM en 1993 par la société
canadienne Ballard [Can08].
Depuis, l’objet de la recherche sur les cellules de la pile est au niveau de l’application
pratique : l’amélioration de la performance et de la durée de vie en utilisation.

1.3.2 Types de piles à combustible

La classification des piles à combustible se fait généralement selon la nature de l’électro-
lyte car celui-ci détermine, d’une part, la température à laquelle la pile fonctionne et, d’autre
part, le type d’ion assurant la conduction ionique. La classification repose essentiellement
sur :
– La nature de la membrane : liquide ou solide
– La température de fonctionnement de la pile : basse (60 ◦ C-250 ◦ C) ou haute (600 ◦ C-
1000 ◦ C)

TAB . 1.1: Types des piles à combustibles

Il existe 6 types de piles à combustible qui se différencient entre elles par leur température
de fonctionnement et par l’électrolyte donnant son nom à la pile.
Les piles à basse température fonctionnent entre 60 ◦ C-250 ◦ C :
– la pile à combustible alcaline ou alkaline fuel cell (AFC)
– la pile à combustible à membrane échangeuse de proton (polymer electrolyte membrane
fuel cell : PMEFC)
– la pile à combustible à méthanol direct (direct methanol fuel cell : DMFC)
– la pile à combustible à acide phosphorique (phosphoric acid fuel cell : PAFC)
Les piles à haute température fonctionnent entre 600 ◦ C-1000 ◦ C :
Piles à combustible de type PEM 13

– la pile à combustible à oxyde solide (solid oxide fuel cell : SOFC)

– la pile à combustible à carbonates fondus (molten carbonate fuel cell : MCFC)
Dans le tableau 1.1, sont présentées les principales caractéristiques de ces différentes
familles de piles ainsi que leurs domaines d’application.

1.4 Piles à combustible de type PEM

Parmi toutes les familles existantes, la pile à combustible de type PEM suscite de nom-
breux travaux de recherche et développement à travers le monde. La technologie évolue vite,
d’autant plus qu’elle est poussée par la volonté des constructeurs de piles de proposer le plus
rapidement possible des produits économiquement viables et fiables.
Nous ne nous intéresserons qu’à la pile à combustible de type PEM, parce que cette tech-
nologie semble effectivement être la plus proche de la commercialisation dans le domaine
du transport. En comparaison avec d’autres types des piles à combustibles, la pile de PEM a
plusieurs caractéristiques attractives [HCCP05] comme :
– Densité de puissance plus élevée ;
– Fort rendement de transformation d’énergie chimique en énergie électrique ;
– Basses températures, rapide et facile à démarrer ;
– Polymère solide ce qui réduit les soucis liés à la construction, au transport et à la fiabi-
lité ;
– Plus compacte et légère : une meilleure densité de puissance volumique ;
– Modulaire donc facile à installer.

1.4.1 Principe de fonctionnement de la pile

Une pile à combustible fonctionne sur le principe inverse de l’électrolyse de l’eau. Des
réactifs chimiques permettent la production d’énergie électrique. De manière générale, le
combustible est de l’hydrogène et il se combine avec de l’oxygène pour former de l’eau
selon une réaction chimique globale universellement connue :

1
H2 + O2 → H2 O (1.1)
2

La réaction électrochimique met en œuvre une membrane d’électrolyte solide conductrice

de protons, et théoriquement imperméable aux gaz réactifs (Hydrogène et Oxygène). L’élec-
trolyte est prise en sandwich entre deux électrodes (l’anode et la cathode). L’assemblage
membrane d’électrodes (MEA, Membrane Electrode Assembly), est représentée schémati-
quement sur la figure 1.6.
Les électrodes mettent en jeu des catalyseurs pour activer la réaction d’un côté. À l’anode,
14 Chapitre 1 : Présentation de la pile à combustible

F IG . 1.6: Schéma d’une MEA relatif au principe de fonctionnement d’une pile de PEM

les ions et les électrons sont produits par l’oxydation du combustible (hydrogène), selon
l’équation 1.2. Les électrons sont toujours transférés de l’anode à la cathode. Et puis, à la
cathode, les électrons permettent la réduction de l’oxydant (oxygène de l’air) et la produc-
tion des ions d’oxygènes comme l’équation 1.3. Les réactions électrochimiques peuvent être
représentées de façon schématique par :
À l’anode :

H2 → 2H + + 2e− (1.2)

À la cathode :

1
O2 + 2e− → O2− (1.3)
2

Comme l’électrolyte est un conducteur protonique, le matériau ne permet que le passage

des ions. Les ions produits sur l’anode sont toujours transférés de l’anode à la cathode. Fi-
nalement, à la cathode, les ions se combinent avec les ions oxygènes pour donner des molé-
cules d’eau. Les réactions électrochimiques peuvent être représentées de façon schématique
par l’équation 1.4.
À la cathode :

O2− + 2H + → H2 O (1.4)

Cette relation correspond à celle de la combustion de l’hydrogène et l’oxygène. L’énergie

de la réaction est libérée sous forme de chaleur et sous forme d’électricité. La combustion se
produit à l’aide d’un catalyseur en platine (Pt).
Piles à combustible de type PEM 15

La réaction chimique globale de la pile est la suivante :

1
O2 + 2H + + 2e− → H2 O+ électricité + chaleur (1.5)
2

1.4.2 Performance actuelle de la pile

La pile à combustible convertit directement l’énergie chimique en énergie électrique.
L’énergie chimique libérée peut être calculée par la variation de l’énergie libre de Gibbs
(∆gf ), qui est la différence entre l’énergie des produits et l’énergie des réactives. L’énergie
libre de Gibbs est utilisée pour calculer l’énergie disponible pour effectuer le travail externe.
La variation de l’énergie libre de Gibbs pour la pile à combustible est [LD03] :

1
∆gf = (gf ) des produits − (gf ) des réactifs = (gf )H2 O − (gf )H2 − (gf )O2 (1.6)
2

La variation de l’énergie libre de Gibbs dépend de la température et de la pression comme :

1
pH2 pO2 2
∆gf = ∆gf0 − RT ln[ ] (1.7)
pH2 O

où ∆gf0 est la variation de l’énergie libre de Gibbs à la pression standard (1 bar) qui dépend
de la température T exprimée en Kelvin. pH2 , pO2 et pH2 O sont les pressions d’hydrogène,
d’oxygène et de vapeur d’eau. R est la constante universelle des gaz (8,31451 J.kg −1 .K −1 ).
Notons que la valeur de ∆gf0 est négative (-237.2 kJ.mol−1 ) car l’énergie est libérée par la
réaction. S’il n’y avait pas de pertes dans la pile à combustible, toute l’énergie libre de Gibbs
serait convertie en énergie électrique. Pour chaque môle d’hydrogène, deux électrons passent
par le circuit électrique externe et le travail électrique est égal à la variation de l’énergie libre
de Gibbs si le système est sans pertes, le travail électrique effectué est :

∆gf = −nF E (1.8)

où F est la constante de Faraday (96,485 Coulombs/mole) qui représente la charge élec-

trique d’une môle d’électrons. n correspond au nombre de moles d’électrons de la réaction :
pour l’hydrogène n = 2, mais pour d’autres carburant, cette valeur peut être plus élevée (n =
6 pour le méthanol, par exemple). E est la tension de circuit ouvert de la pile à combustible.
La tension ouverte de circuit de pile à combustible peut donc être écrite comme :
1
−∆gf −∆gf0 RT pH2 pO2 2
E= = + ln[ ] (1.9)
2F 2F 2F pH2 O

En pratique, le fonctionnement des piles à combustible s’accompagnent de pertes, une

16 Chapitre 1 : Présentation de la pile à combustible

partie de l’énergie chimique est convertie en chaleur. Le terme −∆gf /2F varie en fonction
du point de fonctionnement. Il est égal à 1,229 volt à l’état standard (25 ◦ C) et 1 bar. On peut
exprimer la tension E sous la forme [BL05] :

1
E = 1.299 − 0.85 × 10−3 ¦ (T − 298.15) + 4.3085 × 10−5 T ¦ [ln(pH2 ) + (pO2 )] (1.10)
2

La tension dépend de la température, de la pression de l’hydrogène et de l’oxygène

(Eq.1.10).
La courbe de polarisation de la pile de PEM est décrite comme un modèle statique ; pro-
posé par de nombreux travaux ([CFC04] [Iqb03] et [LLA98]). En fonctionnement, lorsque
la pile produit du courant, la tension est toujours inférieure à la tension à vide. En fait, des
pertes de tension existent dans la pile à combustible à cause de différentes polarisations. Ces
chutes de tension sont parfois appelées polarisations ou surtensions en référence à l’électro-
lyse (overpotential en anglais). Elles sont liées :
– aux vitesse de transferts de charges (on parle de polarisations d’activation),
– aux vitesse de transfert (d’alimentation) des réactifs (on parle de polarisations de concen-
tration),
– aux résistances (on parle de polarisations d’ohm)

F IG . 1.7: Courbe de polarisation de la pile PEM

La tension VF C est la somme de quatre termes : la tension à vide E, la surtension d’ac-

tivation Vact (ou chute d’activation de la région 1), la surtension ohmique Vohm (ou chute
ohmique de la région 2) et la surtension de concentration Vconc (ou chute de concentration de
la région 3) :

VF C = E − Vact − Vohm − Vconc (1.11)

Piles à combustible de type PEM 17

Les pertes ohmiques sont dûes à la résistance qu’opposent les électrodes et les plaques
bipolaires à la circulation des électrons et l’électrolyte au passage des protons. La chute de
tension ohmique correspondante s’écrit :

Vohm = Rm (IF C + in ) (1.12)

où Rm est la résistance totale de la pile à combustible, IF C est le courant délivré par la pile
à combustible, in est le courant interne permettant de tenir compte d’une éventuelle traversée
de gaz et/ou d’électrons à travers l’électrolyte. Sur la figure 1.7, nous pouvons constater que
la perte ohmique correspond à un grand domaine dans la courbe, la résistance a donc un
rôle important dans la performance des piles. Ces pertes sont proportionnelles au courant
électrique. La résistance totale comprend la résistance ohmique et la résistance de contact.
La résistance ohmique est dûe à la résistance électrique du matériau par les composants eux-
mêmes de la pile comme MEA, GDL et BPP. La résistance de contact est dûe à la résistance
de la conductivité sur les interfaces entre les composants de la pile telles que MEA et GDL
ou GDL et BPP. La résistance de contact peut conduire à une performance médiocre de la
pile à combustible et produire de la chaleur susceptible de provoquer une fissure thermique.
Les résistances existantes sont dues à la structure de la pile à combustible et au processus
opératoire que nous allons présenter par la suite.

1.4.3 Structure de la pile à combustible

Comme mentionné dans la performance actuelle de la pile, la tension d’une seule cellule
à l’état standard (25 ◦ C et 1 bar) est d’environ 1.299 V. Par ailleurs, la tension d’une cellule
sous le chargement est environ 0.6 V - 0.7 V. Par conséquent, il est nécessaire de disposer des
cellules en série électrique, formant un “empilement” afin d’atteindre une tension suffisante.
La tension d’un empilement est alors la somme de toutes les tensions des cellules. L’empile-
ment de la pile peut satisfaire les besoins en tension de systèmes alimentés électriquement.

La structure de la pile de PEM est simple par rapport aux autres types de piles à combus-
tibles parce que les composants de la pile sont solides et modulaires. Aujourd’hui, beaucoup
de développeurs ont proposé leur conception des piles. Une structure typique de la pile de
PEM a été introduite, la figure 1.8 montre un schéma de la structure d’un empilement de la
pile PEM.
Un empilement de la pile est constitué de cellules unitaires, serrées par des plaques de
serrage (PS) aux deux extrémités pour la conductivité électrique entre ces plaques.
Une cellule unitaire comprend 5 éléments : Un assemblage membrane des électrodes
(MEA), deux couches de diffusion gazeuses (GDL) et deux plaques bipolaires (BPP). La
plaque MEA est insérée entre deux GDLs en carbone poreuses enduites de PTFE [Pou63]
18 Chapitre 1 : Présentation de la pile à combustible

F IG . 1.8: Schéma de la structure d’un empilement de la pile de PEM

pour distribuer les combustibles sur les électrodes. Les GDLs et la plaque MEA sont en-
core pris en sandwich entre deux BPPs en graphite, sur lesquelles la distribution des canaux
gravés et des dents est alternée.
La réalisation complète d’une telle structure comprend 4 étapes principales comme :
– l’alimentation des gaz réactifs ;
– la réaction électrochimique ;
– la conductivité ionique par l’électrolyte et conductivité électrique de l’anode à la ca-
thode ;
– l’évacuation de l’eau produite.

1.4.4 Assemblage de la pile à combustible

L’assemblage de la pile à combustible est indiqué dans la figure 1.9 et la figure 1.10.

1.4.5 Le processus opératoire de la PAC

Le schéma 1.11 montre le processus opératoire dans un empilement de la pile avec deux
cellules comme exemple. Des gaz réactifs (hydrogène et oxygène) alimentent la pile à com-
bustible parallèlement, ils circulent dans les canaux gravés de la BPP, puis traversent la GDL
et atteignent le MEA où la réaction électrochimique se passe conformément à ce qui est
décrit dans le principe de la pile à combustible de PEM (figure 1.6).
Comme dans l’équation 1.2, l’oxydation de l’hydrogène produit les ions et les électrons à
l’anode du MEA à l’aide du catalyseur. Le transfert ionique est en sens inverse par rapport à
Piles à combustible de type PEM 19

F IG . 1.9: Assemblage d’une cellule de la PAC de PEM avec tous les composants

F IG . 1.10: Processus de l’assemblage d’une cellule de la PAC

20 Chapitre 1 : Présentation de la pile à combustible

F IG . 1.11: Schéma du processus opératoire de la réaction de la pile

la conductivité électrique dans le schéma 1.11. Les ions passent par la membrane échangeuse
de l’anode à la cathode. En revanche, les électrons doivent parcourir une longue distance,
traverser les GDLs, BPPs de l’anode d’un MEA à l’autre cathode du MEA suivant. Ainsi,
la conductivité électrique se passe entre les plaques de l’empilement de la PAC d’un côté à
l’autre.
La résistance électrique est inévitable pendant la conductivité électrique, ce qui conduit à
une perte de performance comme on le présente dans la courbe de polarisation (figure 1.7).
La résistance des plaques et la résistance entre les plaques ont des influences importantes
sur la performance de la PAC. Pour les ions, la résistance est déterminée par l’épaisseur
du MEA. Une épaisseur mince facilite la conductivité des ions et réduit donc la résistance.
La conductivité électrique est dûe non seulement aux résistances des plaques eux-même,
mais aussi aux résistances de contact aux interfaces entre GDL et MEA, GDL et BPP. Les
caractéristiques mécaniques aux interfaces de la PAC sont donc importantes.
Tout comme une bonne conductivité des électrons, l’alimentation efficace des gaz est
aussi importante. Il est nécessaire d’alimenter les gaz réactifs de manière continue pour une
production continue des électrons de la pile. Par ailleurs, la demande des gaz réactifs est forte
pour un niveau élevé en courant électrique. L’alimentation uniforme des gaz offre un temps
suffisant de la réaction. L’alimentation des gaz réactifs est réalisée par canaux gravés sur la
BPP et par porosité de la GDL. D’une part, la configuration, la forme et la taille des canaux
des BPPs sont fondamentales dans l’alimentation des gaz réactifs ; et d’autre part, la porosité
Piles à combustible de type PEM 21

des GDLs est importante dans le transfert des gaz réactifs vers les MEAs pour la réaction
électrochimique.

1.4.6 Composants de la PAC, leurs rôles et leurs caractéristiques

La réalisation de la PAC met en œuvre un grand nombre de composants et de couches
multiples dans la structure lamellaire des piles à combustibles. Ainsi une connaissance ap-
profondie des composants, leurs rôles et leurs caractéristiques mécaniques est fortement de-
mandée pour l’amélioration de la performance de la pile.

1.4.6.1 Assemblage membrane des électrodes (MEA)

Sur la figure 1.12, un profil de MEA est représenté. Il se compose d’une membrane échan-
geuse (ou membrane d’électrolyte de polymère) solide serrée entre deux électrodes, l’anode
et la cathode. Les électrodes fournissent l’interface entre les gaz réactifs et l’électrolyte. Elles
doivent permettre le passage des gaz réactifs humides, fournir une surface de réaction où les
gaz entrent en contact avec l’électrolyte, être conductrices aux électrons libres et assurer le
passage des ions vers la cathode.

F IG . 1.12: Profil de la MEA (électrolyte, anode et cathode) et GDL [CEA06]

Un catalyseur à base de platine est intégré entre membrane et électrodes. Le catalyseur

22 Chapitre 1 : Présentation de la pile à combustible

améliore la réaction chimique en fournissant les sites de réaction mais n’est pas consommé
dans le procédé. Le platine est typiquement utilisé en raison de sa haute activité, sa stabilité
et sa conductivité électrique.
La membrane d’électrolyte est très mince, son épaisseur varie de 50 à 200 microns, dans
le but de conduire les protons de l’anode vers la cathode tout en étant un excellent isolant
des électrons. Son épaisseur est issue d’un compromis entre le gain en performance et une
certaine épaisseur permettant de séparer les gaz hydrogène et oxygène, sans quoi, les deux
demi-réactions d’oxydoréduction seraient alors localisées au même endroit et aucun électron
ne circulerait par le circuit externe.
Enfin, la membrane très mince et fragile doit tenir mécaniquement car étant au cœur de
la pile, son endommagement peut diminuer la durée de vie des piles. D’ailleurs dans un
empilement de la pile, assemblé en série électrique, un dysfonctionnement d’un composant
d’une cellule peut conduire aux dysfonctionnements de tout l’empilement de la pile.
Un certain nombre de membranes commerciales sont disponibles, dont la série de N af ionr
[DuP06] de perfluorosulfonic acid (PFSA) d’utilisation très large, produit par Dupont, dont
les caractéristiques sont indiquées dans le tableau 1.2.

TAB . 1.2: Caractéristiques des MEAs de la série de Nafion

Propriétés Unit Nafion-112 Nafion-113 Nafion-115 Nafion-117
Épaisseur mm 0.051 0.091 0.127 0.183
3
Densité g/cm 1.0 1.9 2.5 3.6
Conductivité S/cm 0.083 0.083 0.083 0.083
Limite tensile maximale
50%RH Mpa 43 43 43 43
23 Mpa 34 34 34 34
100 Mpa 25 25 25 25

1.4.6.2 Couches de diffusion gazeuse (GDL)

La figure 1.12 montre que la GDL s’est imprégnée du MEA, pour fournir une plus grande
interface possible et jouer un triple transport des gaz réactifs, des électrons et de l’eau pro-
duite sur toute la surface de l’électrode.
Cette couche conductrice et poreuse est classiquement en tissu, mousse, hydrophobe (non-
mouillable) et non-corrosive. Les couches poreuses sont là d’une part pour assurer une répar-
tition uniforme des gaz de canaux sur l’ensemble des électrodes, d’autre part, pour évacuer
l’eau produite à la cathode durant la réaction électrochimique. Des particules hydrophobes de
polytétrafluoroéthylène (PTFE) sont ajoutées aux couches poreuses, afin d’éviter l’accumu-
lation d’eau durant la réaction électrochimique. Afin de mieux connaître la structure d’une
Piles à combustible de type PEM 23

GDL, des images par microscope électron scanners (SEM) ont été capturées sur la couche
de fibre carbone de la GDL (figure 1.13).

F IG . 1.13: Image SEM des couches fibres carbones de la GDL [Mac09]

L’épaisseur des couches est de quelques centaines de microns pour faciliter la conductivité
des électrons de l’anode à la cathode. En fait, l’utilisation des GDLs est venue de TOARY
au Japon, SGL en Allemagne et BALLARD au Canada, etc.. Etant donné le coût de la GDL,
la série de TGP-H de TORAY est prise en compte dans nos études. Les caractéristiques des
GDLs de la série de TGP-H de TORAY sont indiqués dans le tableau 1.3 [Jap08].

TAB . 1.3: Tableau des caractéristiques de GDL de la série de TGP-H

Propriétés Unit TGP-H-030 TGP-H-060 TGP-H-090 TGP-H-120
Épaisseur mm 0.11 0.192 0.28 0.377
3
Densité g/cm 0.40 0.44 0.44 0.45
Porosité % 80 78 78 78
Rugosité µm 8 8 8 8
Perméabilité du gaz ml·mm/ 2500 1900 1700 1500
(cm2 ·hr·mmAq)
Résistivité électrique
Par plat mΩcm 80 80 80 80
Dans plat mΩcm 5.8 5.8 5.6 4.7
Limite en flexion Mpa 40 40 40 40
Module en flexion Gpa 8 10 10 10

1.4.6.3 Plaques bipolaires (BPP)

Les plaques bipolaires sont comme une frontière entre deux cellules de la pile et inter-
viennent dans leur tenue mécanique. Elles constituent 80% du poids d’un empilement de
24 Chapitre 1 : Présentation de la pile à combustible

la pile. L’utilisation des plaques nécessite un usinage de petits canaux et de petits dents al-
ternés à la surface des plaques bipolaires afin d’approvisionner les gaz et de transférer les
électrons. Les plaques peuvent prendre des formes diverses et variées. On trouve ainsi les
configurations parallèles, serpentins, parallèle-serpentins, grille, etc. comme indiqué dans la
figure 1.14.

F IG . 1.14: Profils des plaques bipolaires commerciales [Men00]

Elles ont plusieurs rôles comme :

– l’approvisionnement des gaz (hydrogène et oxygène) à travers les microcanaux sur les
interfaces réactives,
– la tenue mécanique des cellules élémentaires,
– la conductivité des électrons à travers les microdents de l’anode d’une cellule à la ca-
thode de la cellule suivante,
– évacuer les gaz non consommés et l’eau produite par la réaction électrochimique.
La plupart des plaques bipolaires sont généralement réalisées dans un matériau à base
de carbone (graphite). La structure graphite a une masse volumique de 2g/cm3 . Les plaques
bipolaires doivent être minces avec des zones de contacts uniformes avec les GDLs afin
de réduire la résistance de contact électrique et thermique. Sans diminuer la performance
des piles, il est intéressant de minimiser l’épaisseur de la BPP pour abaisser le volume, le
poids et le coût de la procédure de fabrication [Hoo03]. Dans une cellule unitaire, la tenue
mécanique est réalisée en serrant les MEA, GDL et BPP ensemble, en créant des zones de
contact permettant la conductivité du courant électrique par les interfaces, en particulier,
Piles à combustible de type PEM 25

entre la GDL et les petits dents de la BPP. La configuration des dents est donc un paramètre
important pour la performance de la PAC. Les canaux sont significatifs pour l’alimentation
des gaz réactifs et les dents de la BPP en contact avec la GDL déterminent la conductivité
des électrons produits.
Des caractéristiques des bipolaires et des paramètres sont indiqués dans le tableau suivant
[All09].

TAB . 1.4: Caractéristiques des BPPs de la SGL

Propriétés Unit BPP4 PPG86
3
Densité g/cm 1.98 1.85
Rugosité µm 4 4
Conductivité
Par plat s/cm 210 500
Dans plat mΩcm 80 180
Limite flexion Mpa 53.1 40
Limite compression Mpa 76 76
Limite traction Mpa 25.1 25.1
Température d’opération ◦ C 180 80

1.4.6.4 Joints d’étanchéité

La parfaite étanchéité est une condition indispensable pour le bon fonctionnement de la

pile à combustible. Des joints d’étanchéité doivent être interposés périphériquement entre
les BPPs pour éviter toute perte de fluides possible. En effet, la perte des gaz peut réduire
la puissance de la pile et la rencontre directe de l’oxygène et l’hydrogène conduit à une
explosion dans la pile [AUW05][HSK07]. Les joints d’étanchéité correspondent à la taille
de la BPP et sont conçus avec une dimension intérieure de la taille de la GDL. L’épaisseur
des joints doit correspondre aux épaisseurs de MEA et GDL comme indiqué dans la figure
1.15.

F IG . 1.15: Schéma de la localisation des joints dans une cellule de la PAC

Une compression suffisante des joints est nécessaire pour assurer une parfaite étanchéité
26 Chapitre 1 : Présentation de la pile à combustible

de la PAC. En effet, le niveau minimum de compression des joints est donné par :

∆d
ε= × 100% (1.13)
d0

où ε est le ratio de compression des joints, d0 l’épaisseur de la section et ∆d est le niveau

de compression des joints. Dans une étude sur l’assemblage d’un empilement de pile à com-
bustible, la condition indispensable à soumettre est la parfaite étanchéité afin d’éviter toutes
pertes des gaz au sein de chaque cellule. En fait, le niveau de compression des joints est de
10% à 30% [Fri01].

1.5 Étude mécanique de la structure avec des interfaces

multi-contacts
Comme le procédé de la réaction dans une structure lamellaire, le contact sur les inter-
faces entre les couches multiples des piles est une condition indispensable pour garantir
la conductivité électrique. La structure de la pile peut être regardée comme une structure
multi-contacts. La figure 1.16 montre un schéma permettant de localiser les interfaces multi-
contacts dans une cellule unitaire. Des contacts multiples existent sur les interfaces comme
(1) : contact entre MEA et GDL, (2) : contact entre BPP et GDL ou plus précisément entre
GDL et les petits dents de BPP.
Comme l’indique la courbe de la polarisation, la pile n’est pas parfaitement utilisée. Elle a
des pertes de performance et des dégradations irréversibles du cœur de la pile qui font chuter
la durée de vie. La perte de la performance est principalement dûe à la résistance de contact
dans la PAC et à la chute de l’alimentation gazeuse dans la pile [ZLS+ 06] [BSP+ 08] [ZL06]
[KFT08].
Le transfert des électrons est influencé par les effets mécaniques aux interfaces. Les re-
cherches pour améliorer la performance de la pile à combustible ont porté sur la préparation
de la membrane, l’efficacité du catalyseur, le matériel de la BPP, etc.. Mais des études ré-
centes montrent que des effets mécaniques jouent un rôle important sur la performance de la
pile [NHM08] [MABHAJ07].

1.5.1 Force de serrage

Des études récentes ont constaté que presque 59% de la perte de puissance est dû à la
résistance de contact entre la BPP et la GDL [RA09]. En effet, la résistance de contact
dépend fortement d’un parfait état du contact entre les composants de la PAC. Pour cette
raison, une étude de contact de la structure multi-contacts de la PAC paraît nécessaire.
Étude mécanique de la structure avec des interfaces multi-contacts 27

F IG . 1.16: Schéma d’une cellule unitaire représente des interfaces multi-contacts

Un paramètre principal qui influence l’état de contact est la force de serrage par boulons
tirants. Cette force permet de maintenir les éléments de l’empilement en contact les uns avec
les autres pendant l’assemblage. Une force élevée conduit à l’augmentation de la conductivité
et de la densité de puissance de la PAC [CHWC07]. Lee et al. [LHVZM99] ont fait des
expériences avec trois types de GDLs : TORAY, ELAT, une combinaison de TORAY et
CARBEL, sous trois types de forces de serrage, 100, 125 and 150 in.lbf (1 in.lbf= 0.113
N.M). Avec ces trois forces de serrage, ils ont constaté que chaque GDL a sa propre force
de serrage optimale pour acquérir une puissance maximale de la pile. Les résultats montrent
que la performance de la pile est influencée par la force de serrage. Ge et al. [GHL06] ont
trouvé qu’il existe une force de serrage optimale donnant la meilleure performance des piles
possible : si la force de serrage est inférieure à la force optimale et qu’elle augmente, la
performance augmente. Chang et al. [CHWC07] ont trouvé que l’influence de la force de
serrage sur la performance des piles est principalement dû au changement de la résistance
de l’interface entre GDL et BPP. Une force élevée peut diminuer la résistance de contact
aux interfaces de la PAC. De plus, Wang et al. [WSZ08] ont étudié l’influence de la force
de serrage uniforme sur la performance de la PAC. Dans son expérience, ils ont trouvé que
la performance de la PAC avec une force répartie uniformément est meilleure qu’avec une
force conventionnelle. Mishra et al. [MYP05] et Wang et al. [WST03] ont développé une
relation entre la résistance de contact et la force de serrage, ils ont constaté une diminution
de la résistance de contact lorsqu’on augmente la force de serrage.
28 Chapitre 1 : Présentation de la pile à combustible

Les modélisations théoriques étant complexes, des expériences ont permis de découvrir
certains liens entre les efforts appliqués à la pile et l’efficacité de celle-ci. Cependant les ex-
périences sont souvent coûteuses et n’offrent qu’un résultat global sur la pile, ne permettant
pas d’identifier clairement les interactions des phénomènes mécaniques aux interfaces sur la
performance de la PAC. Des modélisations numériques ont donc été abordées.
D’après les expériences de Mishra et al. [MYP04], Mishra et al. [MYP05] et Wang et al.
[WST03], Zhou et al. [ZWM06] ont évalué la résistance de contact par rapport à l’aire de
contact et à la pression de contact entre la BPP et la GDL par une équation semi-empirique.
Dans leur modèle 2D, la méthode des éléments finis est utilisée afin d’étudier l’influence
de différents configurations de la BPP sur la résistance aux interfaces entre la GDL et la
BPP. Zhang et al. [ZLS+ 06] ont étudié la résistance de contact par une expérience, ils ont
obtenu une équation précise pour évaluer la résistance de contact. Ils ont ensuite développé
un modèle numérique 3D afin d’étudier l’influence de différentes forces de serrages sur la
résistance de contact aux interfaces. La simulation montre que la résistance de contact est
liée à l’aire de contact par l’équation suivante :

A B
RBP P /GDL = + (1.14)
Acontact C

où, RBP P /GDL est la résistance de contact entre la BPP et la GDL, Acontact est l’aire de
contact, A, B sont les facteurs constants obtenus par l’expérience, et C est le chargement
extérieur. L’équation de Zhang et al. [ZLS+ 06] montre que l’aire de contact est un paramètre
important dans l’évaluation de la résistance de contact dans la PAC. Ainsi toute étude per-
mettant de connaître précisément l’aire de contact permettra de bien évaluer la résistance de
contact.

1.5.2 Porosité de la GDL

Dans la pile à combustible, la GDL en papier carbone est classiquement en tissu et mousse
avec une porosité initiale élevée afin de faciliter l’alimentation du gaz pour la réaction élec-
trochimique au sein de la PAC. La capacité des GDLs à transporter les gaz réactifs est me-
surée par la perméabilité. Des études sur la perméabilité de la plaque GDL [Sch97] [IML04]
[WBB+ 04] [KH98] [APT02] [UFS+ 96] [JSB+ 00] ont montré l’importance de son influence
sur la performance de la PAC. Feser et al. [FPA06] et Gostick et al. [GFP+ 06] ont mesuré
la perméabilité de la GDL comme une fonction de la compression et ont constaté qu’elle
diminuait avec l’augmentation de la force de serrage. Un calcul général de la perméabilité
des gaz est basé sur la loi de Darcy avec un faible nombre de Reynolds.

q·δ
k = K·υ = (1.15)
A · ∆p
Étude mécanique de la structure avec des interfaces multi-contacts 29

où q est le flux de fuite gazeuse ; υ est l’épaisseur de la couche ; A est l’aire de la couche ;
∆p est la différence de pression ; δ est la vélocité gazeuse, K est la perméabilité gazeuse
de la couche et k est la perméabilité gazeuse spécifique. Mais la loi de Darcy est une loi
idéale, une relation entre la perméabilité et la propriétés des pores fut proposée par Kozeny
et modifiée par Carman. L’équation résultante est connue sous le nom de Kozeny-Carman
(KC), elle est une prédiction plus précise, fonction de la porosité et de la taille du filament
des fibres. Depuis sa première apparition, cette équation a pris plusieurs formes, celle qui est
couramment utilisée [FPA06] est :

ϕ3
k=C (1.16)
(1 − ϕ)2

où k est le cœfficient de perméabilité ; C est un cœfficient constant ; ϕ est la porosité.

C’est l’une des relations très largement utilisée dans l’étude des milieux poreux. Mais pour
atteindre une estimation précise, l’équation de KC demande une évaluation d’un facteur
constant C, qui doit être mesuré de façon expérimental. Un modèle plus précis et robuste a
été créé par Tomadakis et al. [TR05] dans une étude des milieux poreux, il offre une équation
pour mesurer la perméabilité :

ϕ (ϕ − εp )α+2
k = r2 (1.17)
8ln2 ϕ (1 − εp )α [(α + 1)ϕ − εp ]2

où la perméabilité k est comme une fonction de la porosité ϕ, du diamètre du filament de

la GDL r, de l’orientation du filament α et d’un cœfficient constant εp . Dans les couches, des
filaments sont positionnés aléatoirement (figure 1.17).

F IG . 1.17: Profile de la plaque de la GDL [TR05]

30 Chapitre 1 : Présentation de la pile à combustible

Le tableau 1.5 donne les paramètres concernant les filaments de la plaque GDL utilisés
dans le modèle de Tomadakis.

TAB . 1.5: Paramètres du modèle de Tomadakis

r (µm) εp α
3.5 0.11 0.785

Selon l’équation 1.17, la proportion entre la perméabilité et la porosité de la plaque GDL

peut être décrite par une courbe (figure 1.18). qui montre que la perméabilité de la plaque
GDL augmente avec la porosité de la plaque.

F IG . 1.18: Relation "perméabilité-porosité"

Actuellement, la plupart des études dans ce domaine fait l’hypothèse simplificatrice d’une
porosité de la GDL constante. Cependant, ceci ne peut pas refléter l’importance de la porosité
de la GDL sur la performance des piles. Tout changement de la morphologie de la GDL peut
conduire à une influence importante sur les performances des piles à travers le changement
de la porosité de la GDL.
La porosité est un facteur permettant d’évaluer l’espace au transport des fluides. C’est le
pourcentage de l’espace poreux par unité de volume total. La porosité est définie comme :

Vp
ϕ= (1.18)
Vt

où Vp est le volume poreux et Vt est le volume total. Le volume solide comme le volume
des fibres de la GDL est Vs = Vt − Vp . Si la plaque GDL est déformée sous une force
extérieure, la porosité de la GDL ϕ après la déformation par rapport à la porosité initiale ϕ0
Étude mécanique de la structure avec des interfaces multi-contacts 31

devient :

(Vt0 − ∆Vt ) − (Vs0 − ∆Vs )

ϕ= (1.19)
Vt0 − ∆Vt

où Vt0 est la valeur initiale du volume total, ∆Vt est sa variation ; Vs0 est la valeur initiale
du volume solide, ∆Vs est sa variation.

ϕ0 − εv + ∆Vs /Vt0 ϕ 0 − εv + εs
ϕ= = (1.20)
1 − εv 1 − εv

avec,

Vt0 − Vs0 ∆Vt εx + εy + εz ∆Vs

ϕ0 = , εv = = , εs = (1.21)
Vt0 Vt0 Vt0 Vt0

Si le changement du volume solide des fibres peut être négligée, c’est-à-dire, ∆Vs = 0,
l’équation 1.21 devient :

ϕ 0 − εv ϕ0 − 1
ϕ= =1+ (1.22)
1 − εv 1 − εv

La porosité initiale ϕ0 de la GDL, est déterminée par les caractéristiques du matériau.

L’équation 1.22 montre que la porosité de la GDL décroît lors de l’augmentation de la dé-
formation de la GDL. Une déformation élevée de la GDL conduit à une baisse de la porosité
et de la perméabilité de la GDL qui doit être étudiée.
Ainsi, la force de serrage influence non seulement la résistance de contact aux interfaces
mais aussi la porosité de la GDL de deux façons contraires. Une force de serrage forte peut
conduire à une diminution de la résistance de contact ce qui facilite le transfert des électrons
produits, mais conduit à une baisse de la porosité de la GDL, ce qui réduit l’alimentation du
gaz réactifs. Ainsi, l’aire de contact et la déformation de la GDL doivent être considérées en-
semble. La configuration de dents de la BPP de la PAC est aussi un des paramètres importants
influençant l’aspect mécanique aux interfaces, comme l’aire de contact et la déformation de
la GDL.
Afin d’améliorer la durabilité de la PAC, il s’avère aussi nécessaire de renforcer la connais-
sance des effets mécaniques de la pile comme la contrainte, la déformation découlant pendant
l’assemblage de la structure de la pile. Des effets mécaniques de la PAC ont été étudiés par
Mehta et al. [MC03] qui ont considéré que la PAC était une structure dans laquelle des phé-
nomènes multi-physiques étaient couplés. La possibilité de la fissure mécanique existe et a
été considérée par Weber et al. [WN04]. Lee et al. [LHH05] ont étudié la distribution de
la contrainte d’une seule cellule soumise à la force de serrage. La contrainte obtenue sur le
MEA a été comparée avec des résultats expérimentaux. Kusoglu et al. [KKS+ 06][KKS+ 07]
ont mis en évidence l’évolution de la contrainte du MEA sous différentes conditions de
32 Chapitre 1 : Présentation de la pile à combustible

l’humidité et la température. Borgrachev et al. [BGGM08] ont développé un modèle 2D sy-

métrique avec une seule cellule. Les résultats montrent que la plasticité peut apparaître sur
le MEA si la force de serrage est trop importante.

1.5.3 Optimisation des composants

Un grand nombre d’études portent sur l’optimisation de la conception des composants
individuels comme la BPP, la MEA et la GDL de la PAC pour une meilleure performance.
La BPP est un composant essentiel de la PAC, car elle fournit l’alimentation du gaz et
la collecte des électrons. Elle permet aussi la tenue mécanique des cellules de la pile. La
configuration des dents et des canaux de la BPP est importante dans la performance de la
PAC. À l’origine, la conception de la BPP avait des canaux droits parallèles comme cela a
été présenté par Pellegri et al. [PS80]. Cette conception avec canaux droits avait un certain
nombre de désavantages, conduisant à l’accumulation d’eau et à la chute de pression insuf-
fisante dans des canaux, résultant de la distribution de gaz réactifs pauvres. Pour y remédier,
Watkins et al. [WDE91][WDE92] ont développé un canal continu, en serpentin. Ce canal
a une meilleure gestion de l’eau grâce à l’amélioration de la chute de pression. Aussi cette
distribution des canaux a permis une augmentation en puissance de la PAC de 50%. Cette dis-
position a aussi été modifiée pour inclure des circuits parallèles de refroidissement [FCP+ 96]
et les circuits en serpentin segmenté permettant une meilleure distribution des gaz réactifs
[CHW97]. Récemment, Shimpalee et al. ont établi un modèle fluidique [SLVZNN01] afin
d’étudier l’influence de différents longueurs du canal et de différentes largeurs de canal/dent
sur la performance de la PAC [SLVZ06] [SGVZ07]. Ils ont trouvé qu’un facteur de canal/dent
optimal est plus sensible sur la performance de la PAC que la longueur de canal. La configu-
ration de la BPP est importante dans la distribution de la pression aux interfaces. Un canal
droit de la BPP peut conduire à une fissure aux interfaces [SLC+ 08].
La conception de la MEA se base sur les transferts de gaz et électrons et le transfert
efficace de protons vers les réactions électrochimiques [LM04]. C’est difficile d’optimiser la
MEA tout seule parce que les transferts dépendent aussi des effets mécaniques. Paganin et al.
[PTG95] ont présenté une étude sur l’analyse de trois épaisseurs de la membrane Nafion. Les
trois membranes analysées étaient Nafion 112, 115 et 117 avec des épaisseurs respectives
de 50, 125 et 175 um. Ils ont constaté que la plus mince membrane Nafion 112 fournit la
meilleure performance de la PAC, offrant non seulement une plus faible résistance et mais
aussi une augmentation de la densité de courant.
La GDL joue aussi un grand rôle dans la pile, puisqu’elle alimente la MEA en gaz réac-
tifs par les canaux et qu’elle assure le transfert des électrons. La capacité de la GDL pour
transporter efficacement des gaz est généralement quantifiée par la perméabilité. Bien qu’il
existe un certain nombre d’études qui ont mesuré la perméabilité GDL par la simulation nu-
Interactions des effets mécaniques aux interfaces d’une PAC 33

mérique, tel que Gostick et al. [GFP+ 06] et Feser et al. [FPA06], une simulation simple a été
présentée pour observer la performance de la PAC en fonction de la perméabilité. Williams
et al. [WBB+ 04] ont étudié certains effets de la variation de la perméabilité de la GDL et
ont constaté une faible proportionnalité entre la perméabilité et la performance de la PAC.
Cependant, ces études n’ont pas évoqué d’interactions entre des effets mécaniques et la per-
méabilité de la GDL.
En effets, la pile à combustible est une structure multi-contacts dans laquelle des phé-
nomènes multi-physiques sont couplés. Des interactions de ces phénomènes influencent la
performance de la PAC et sa durée de vie. Il est donc important d’étudier en détails les effets
mécaniques pour mieux connaître les interactions au niveau de la pile.

1.6 Interactions des effets mécaniques aux interfaces d’une

PAC
La pile à combustible de PEM doit être étudiée à partir des conditions opératoires et de la
composition de la pile (forme des plaques bipolaires, des canaux, des dents et des épaisseur
de membrane...). L’amélioration des performances des piles à combustibles doit passer par
une meilleure compréhension des différentes interactions des composants mécaniques dans
un empilement des piles.
La figure 1.19 représente un schéma d’une cellule unitaire qui permet de localiser les dif-
férentes interactions des effets mécaniques principaux au sein d’une pile à combustible. Ces
effets mécaniques sont décrits en allant de l’extérieur vers l’intérieur : force de serrage par
boulons tirants, contact sur les interfaces entre composants, déformation de GDL, endom-
magement mécanique de MEA et pression fluidique intérieure des canaux de BPP.

1.6.1 Forces de serrage

Pour le fonctionnement normal de la pile à combustible, l’assemblage des composants
empilés est réalisé en vissant l’empilement par des boulons tirants. Ces forces de serrage
s’exercent en allant de deux extrémités d’un empilement vers l’intérieur de l’empilement de
la PAC. Des composants de la PAC sont maintenus en contact aux interfaces traversées par le
courant électrique. En même temps, ces forces de serrage compressent les joints d’étanchéité
entre plaques bipolaires et évitent les fuites de gaz réactifs. Ainsi la force de serrage permet
d’avoir :
– une étanchéité de la pile à combustible ;
– une tenue mécanique sur les composants de la PAC ;
– une production du contact aux interfaces de la PAC.
34 Chapitre 1 : Présentation de la pile à combustible

F IG . 1.19: Schéma des phénomènes mécaniques existants dans une pile de PEM

La force de serrage permet tout d’abord l’étanchéité des piles. Une parfaite étanchéité est
une condition indispensable pour le fonctionnement normal de la pile. En effet, la perte des
gaz peut réduire la puissance de la pile et la rencontre directe de l’oxygène et l’hydrogène
conduit à une explosion dans la pile [AUW05][HSK07].
Sous l’action de la force de serrage par boulons tirants, les MEAs, GDLs et BPPs sont
assemblés pour former une cellule. Les composants rentrent en contact aux interfaces. À
l’interface BPP/GDL, ce sont plus exactement les dents de la BPP qui sont en contact avec la
GDL. Dans notre étude, le contact aux interfaces est une condition essentielle pour garantir
le transfert des électrons.
Au niveau des efforts, la cellule d’une pile est en équilibre, la force de contact sur les
interfaces est égale à la précharge extérieure. Dans nos études, les matériaux sont élastiques,
l’augmentation de la précharge extérieure peut conduire à l’augmentation de la force de
contact et modifier l’état de contact sur les interfaces. Dans certains cas, la force de contact
peut provoquer des fissures mécaniques sur les interfaces.

1.6.2 Déformation de la GDL

Les composants de la pile sont élastiques, l’augmentation de la force de serrage conduit
à l’augmentation des déformations des composants. À contrainte équivalente, la GDL va se
Interactions des effets mécaniques aux interfaces d’une PAC 35

déformer plus que la BPP et la MEA car son module de Young est petit [JL09]. La GDL
est en papier carbone, avec une porosité élevée afin de faciliter le transfert des gaz. Pour
une GDL de type de TGP-H120, la porosité initiale est de ϕ0 = 78% [All09]. Pendant
l’assemblage de la pile, la plaque GDL est déformée, ce qui peut nuire au passage des gaz
réactifs et diminuer la performance de la PAC. Chang et al. [CHWC07] ont étudié qu’une
force de serrage élevée peut abaisser l’épaisseur, la porosité et la perméabilité de la plaque
GDL et bloquer le transfert des gaz et des électrons, diminuant ainsi la performance de la
pile.

1.6.3 Contact
Le contact est un phénomène non linéaire, complexe et s’appliquant sur des zones pas
toujours évidentes à déterminer, sa prise en compte entre composants de la PAC est impor-
tante afin d’étudier le comportement mécanique aux interfaces. Il est évident que les effets
mécaniques sur les interfaces comme la force de contact joue un rôle fondamental dans le
comportement de la structure.
La résistance de contact permet à un courant électrique de traverser deux composants
dissociables. Une petite résistance peut favoriser le transfert des électrons au sein des piles.
En effet, la résistance de contact dépend fortement d’un parfait état du contact entre les
composants de la PAC. Une séparation de contact peut fortement diminuer la densité de
puissance de la PAC. Il est donc très important d’étudier le contact au sein d’une pile afin de
pouvoir diminuer la résistance de contact dans la PAC.

1.6.4 Endommagement mécanique

L’accroissement de la force de contact conduit à l’augmentation des contraintes méca-
niques au sein de la pile. Une force de contact élevée conduit de façon certaine à l’augmen-
tation des risques d’endommagement mécanique des composants. Les MEAs et GDLs sont
des membranes très minces et fragiles, situées au cœur de la PAC, leur endommagement peut
diminuer la vie des PACs. Ainsi l’endommagement comme la fracture et la fissure des MEAs
ou GDLs doivent être considérées et évitées. Ainsi, les contraintes au sein de la MEA et de
la GDL doivent être étudiées pour éviter tout endommagement mécanique.
Su et al. [SLC+ 08] ont constaté qu’une fissure était apparue sur GDL à cause d’une
concentration de contraintes (Figure 1.20) aux coins des dents de la BPP. Stanic et al. [SH04]
ont étudié l’endommagement mécanique d’un MEA et ont trouvé que l’endommagement est
souvent dû au “Pinhole” (figure 1.21). Ceci apparaît comme une pointe où la contrainte est
maximale.
Il paraît donc intéressant d’étudier et d’optimiser la distribution des contraintes dans les
composants MEAs et GDLs en vue d’obtenir une performance maximale.
36 Chapitre 1 : Présentation de la pile à combustible

F IG . 1.20: Fissure sur GDL [SLC+ 08] F IG . 1.21: ’Pinhole’ sur MEA [SH04]

1.6.5 Pression fluidique

L’alimentation de gaz réactifs dans une pile se fait de la façon suivante : les gaz circulent
dans les canaux, qui sont gravés sur les plaques bipolaires. La conception des canaux est
particulièrement importante pour assurer une répartition homogène des gaz sur toute la sur-
face des électrodes. Les gaz réactifs sont mis sous pression afin d’atteindre le cœur d’une
cellule pour la réaction électrochimique. À cause de la mécanique des fluides, la pression
partielle du gaz est généralement inférieure à la pression totale des gaz à cause de la perte
des pressions. Dans ce travail, ces pressions mécaniques des fluides n’ont pas été pris en
compte parce que les pressions totales et le gradient de pression sont supposées faibles par
rapport aux charges mécaniques issues des boulons. Nous avons intentionnellement négligé
les effets des fluides dans notre étude.

1.7 Conclusions
Au début de ce premier chapitre, nous avons présenté le problème énergétique dans le
monde, l’énergie d’hydrogène est une énergie renouvelable pour le futur, la pile à combus-
tible permet de transformer efficacement de l’énergie chimique en énergie électrique. De
toute la famille des piles à combustibles, la pile à combustible avec la membrane échan-
geuse des protons (PEM) est la plus adaptée à une application dans le transport. Les effets
mécaniques existants dans la pile à combustible influencent la performance de la PAC et sa
durée de vie. Il est donc important d’étudier ces effets. La modélisation numérique est une
approche efficace afin d’étudier les multi-contacts aux interfaces de la PAC.
 37

Chapitre 2

Méthodes numériques et résolution du

problème de contact

2.1 Introduction
Dans le cadre de cette thèse, la modélisation de la pile nécessite une phase de calcul des
forces de contact dans le cas de contacts multiples et simultanés avec frottement. Dans ce
chapitre, nous allons nous intéresser à la résolution du problème du contact. Pour cela, nous
procéderons comme suit :
– dans la première partie, nous effectuerons un rappel de la méthode des éléments finis,
– dans deuxième partie, nous effectuerons un rappel des notions de cinématique du contact,
– dans la troisième partie, nous aborderons la mise en œuvre de la résolution numérique
des méthodes les plus utilisées pour intégrer la résolution du contact avec la méthode
des éléments finis.

2.2 Méthode des éléments finis (MEF)

Le problème de contact est non linéaire, à cause des zones de contact inconnues, et de la
loi de frottement. Depuis longtemps, les chercheurs ont consacré des travaux aux problèmes
de contact. Les notions du contact élastique ont été établies à la fin du 19e siècle par Hertz
et d’autres chercheurs. La théorie de Hertz est utilisée pour résoudre le problème de contact
élastique, dans un cas géométrique simple. Le problème de contact est non linéaire de par la
non connaissance des zones de contact mais aussi par les lois de frottement aux frontières
en contact. De plus, l’évolution du contact dépend de l’histoire du changement. Il est donc
impossible de résoudre les problèmes de configuration complexe par une solution analytique.
L’apparition de la méthode des éléments finis a permis de nouveaux développements et
de nouvelles applications en mécanique du contact. Les études effectuées à partir d’élé-
38 Chapitre 2 : Méthodes numériques et résolution du problème de contact

ments finis offrent certaines voies de résolution des phénomènes mécaniques aux interfaces
de contact.
La méthode des éléments finis s’applique pour l’analyse des déformations et des contraintes
dans des structures mécaniques. Elle est aussi utilisée dans des problèmes de conduction ther-
mique, d’écoulement de fluides et de flux gazeux ou magnétiques, etc.. Plus généralement,
cette technique permet de résoudre des équations différentielles en temps et aux dérivées
partielles en espace avec conditions aux limites [Nes04].
Pour une étude plus approfondie de la méthode des éléments finis, nous conseillons les
cours de Garrigues [Gar99] et de Hunter [HP98] et les ouvrages de références de Zienkiewicz
[ZWHT84], de Bathe [Bat82] et de Hughes [Hug87].

2.2.1 Principe
La méthode des éléments finis est un outil de discrétisation. Elle permet de diviser un
milieu continu complexe (un objet) en un certain nombre fini d’éléments géométriques (un
maillage) relativement simples appelés éléments finis, comme des tétraèdres ou des hexa-
èdres. Les champs physiques sont interpolés sur chaque élément en fonction de leur valeur
en certains points nommés nœuds. Dans les cas les plus simples ces nœuds sont situés aux
sommets de l’élément. A chaque élément est défini un système d’équations obtenu à partir
d’une formulation intégrale variationnelle des équations. Puis un processus d’assemblage
prenant en compte la connectivité entre tous les éléments mène à la construction d’un en-
semble d’équations à résoudre. La résolution de ce système d’équations, en tenant compte
des conditions imposées aux limites du système, permet d’obtenir une approximation conti-
nue des problèmes [CB97].

2.2.2 Maillage
Les éléments finis peuvent être de différentes formes : en dimension 2 les formes les
plus couramment utilisées sont les triangles ou les quadrilatères. En dimension 3, ceux sont
les tétraèdres ou les hexaèdres. Il est possible d’utiliser différents types d’éléments dans
un même maillage mais cela complique significativement l’écriture du système. Dans la
majorité des cas, un seul type d’éléments est utilisé pour un maillage donné, (figure 2.1).

F IG . 2.1: Exemples d’éléments finis surfaciques et volumiques.

Méthode des éléments finis (MEF) 39

Plus le maillage est fin, plus le calcul est précis, mais coûteux. Un compromis s’impose
donc entre la finesse de la discrétisation des objets et le coût du calcul. Nous pouvons aussi
réaliser un maillage grossier dans certaines régions et plus fin dans d’autres, en fonction de
l’importance des déformations et afin d’augmenter la précision là où le calcul s’en fait res-
sentir. Ainsi nous pouvons adapter la densité d’éléments selon les nécessités de la simulation.
Un avantage supplémentaire est apporté par cette représentation locale des propriétés mé-
caniques. Comme les paramètres mécaniques sont associés à chaque élément d’un maillage,
il est possible de construire des modèles hétérogènes composés de différentes structures im-
briquées et de propriétés mécaniques différentes.
À chaque élément fini est associé un certain nombre de nœuds. Les nœuds forment l’en-
semble discret des points du système en lesquels les propriétés physiques sont calculées.
L’ensemble des nœuds ne se limite pas forcément aux sommets des éléments finis. Il est
possible de définir des nœuds additionnels, par exemple un quadrilatère à 8 nœuds ou un
tétraèdre à 10 nœuds (figure 2.2) afin d’assurer une meilleure continuité à la frontière entre
les sous domaines et d’approcher au mieux la solution cherchée.

F IG . 2.2: Interpolation linéaire et quadratique

2.2.3 Fonctions d’interpolation

Le choix du nombre de nœuds de chaque élément est lié au choix d’un schéma d’interpo-
lation. Le principe de la méthode des éléments finis étant d’interpoler les champs physiques
en tout point du solide en fonction de leur valeur aux nœuds, diverses fonctions d’interpola-
tion peuvent être choisies. La plus simple est une fonction d’interpolation linéaire. Dans ce
cas les seuls nœuds utilisés sont les sommets des éléments finis. En choisissant une fonction
d’interpolation quadratique, un nœud additionnel est nécessaire sur chacune des arêtes. Ce
cas correspond au cas du quadrilatère à 8 nœuds ou au tétraèdre à 10 nœuds (figure 2.2). La
précision du résultat est améliorée lorsqu’une fonction d’interpolation de degré plus élevé est
utilisée. Il est possible d’augmenter librement le degré de la fonction d’interpolation mais en
pratique des fonctions de degré supérieur à 3 sont rarement utilisées, étant donné la lourdeur
des calculs que cela implique.
40 Chapitre 2 : Méthodes numériques et résolution du problème de contact

Il est aussi possible de construire un espace de fonctions d’interpolation directement sur

chaque maille réelle, mais généralement, nous procédons en deux temps : nous construisons
d’abord un espace de fonctions d’interpolation Ni sur une maille de référence standard topo-
logiquement équivalente à la maille réelle, puis nous le transformons pour qu’il devienne un
espace de fonctions d’interpolation sur les mailles réelles. J(3 × 3) est la matrice jacobienne
qui correspond à cette transformation. Ce procédé a l’avantage de faire gagner du temps de
calcul (figure 2.3).

F IG . 2.3: La transformation iso-paramétrique d’un élément cubique.

La formule suivante représente un exemple de fonctions d’interpolations utilisées dans le

cas d’un élément cubique de 8 sommets défini dans un espace local (ξ, η, ζ) :
 

 N1 = 1
(1 − ξ)(1 − η)(1 − ζ) 


 8 


 N2 = 1
(1 + ξ)(1 − η)(1 − ζ) 


 8 


 


 N3 = 1
(1 + ξ)(1 + η)(1 − ζ) 


 8 

 N = 1
(1 − ξ)(1 + η)(1 − ζ) 
4 8

 N5 = 1
(1 − ξ)(1 − η)(1 + ζ) 


 8 


 N6 = 1
(1 + ξ)(1 − η)(1 + ζ) 


 8 


 


 N7 = 1
(1 + ξ)(1 + η)(1 + ζ) 


 8 

 N = 1
(1 − ξ)(1 + η)(1 + ζ) 
8 8

Pour un élément fini e composé de s nœuds, on en déduit la relation :

u = Ne u(e) (2.1)
 ¯ ¯ ¯ 
N1 0 0 ¯ N2 0 0 ¯ ... ¯ Ns 0 0
¯ ¯ ¯
 ¯ ¯ ¯ 
avec : Ne =  0 N1 0 ¯ 0 N2 0 ¯ ... ¯ 0 Ns 0 
¯ ¯ ¯
0 0 N1 ¯ 0 0 N2 ¯ ... ¯ 0 0 Ns

n oT
(e) (e) (e)
u (e) (e)
= u1 u2 . . . us , ui = {ui vi wi }T
Méthode des éléments finis (MEF) 41

(e)
ui représente le vecteur déplacement du nœud i appartenant à l’élément e, de même on
(e)
peut lui associer son vecteur position xi = {xi , yi , zi }T .

2.2.4 Assemblage
Dans un maillage d’éléments finis, un nœud donné appartient en général à plusieurs
éléments et la déformation de chacun d’entre eux induit une force en ce nœud. Le calcul
des forces se fait ainsi élément par élément avec des systèmes locaux. Les sommets com-
muns sont des nœuds de connexion qui permettent l’assemblage de tous les systèmes lo-
caux en un grand système global contenant les équations de tous les nœuds. Cela nécessite
une renumérotation des nœuds à travers une table de connectivité. Pour l’ensemble des n
nœuds qui discrétisent l’objet déformable S, on définit le vecteur de position global X =
{x1 , x2 , x3 , ........, xn }T et le vecteur de déplacement globale U = [u1 , u2 , u3 , ........, un ]T .
Le processus d’assemblage est relativement coûteux numériquement. Dans ce qui suit nous
aborderons uniquement la formulation de modèles physiques simples qui peuvent être réso-
lus de manière rapide par la méthode des éléments finis.

2.2.4.1 Élasticité linéaire

Un matériau élastique se déforme sous l’action de certaines forces et regagne son état de
repos, en passant par le même chemin, une fois que ces forces disparaissent. Il restitue ainsi
toute l’énergie qu’il avait absorbée. L’élasticité est de plus linéaire lorsqu’il y a proportion-
nalité entre déformations et efforts imposés. Ce cas particulier n’est constaté qu’en petites
déformations (voir figure 2.4). En général, on considère qu’un matériau élastique n’est plus
linéaire lorsque les déplacements élastiques dépassent 10% de la taille de l’objet considéré.
Au delà, il faut utiliser la théorie de la MMC (Mécanique des Milieux Continus) en grandes
déformations [PDA01]. Pour plus de détails sur les concepts abordés dans ce paragraphe,
nous conseillons les ouvrages de François et al. [FPZ90], de Fung [Fun65] et de Spencer
[Spe84] ; de Ciarlet [Cia88], [Cia90], [Cia97], de Oden [Ode67] et de Ogden [Ogd84].

F IG . 2.4: Approximation linéaire d’une loi de comportement

42 Chapitre 2 : Méthodes numériques et résolution du problème de contact

Nous allons aborder dans la suite quelques éléments de la théorie de la MMC en élasticité
linéaire. Les petites déformations peuvent être écrites d’une manière condensée comme :
 

 εxx 


 
 ³ ´ 

 εyy 
 

 
  ε = ∂u , εxy = 1 ∂u ∂v
+ ∂x 

 ε 
 
 xx ∂x 2 ∂y 

zz ∂v
¡
1 ∂u ∂w
¢
ε= avec εyy = ∂y , εxz = 2 ³∂z
+ ∂x ´

 2εxy 
 
 


 
  ε = ∂w , εyz = 1 ∂v
+ ∂w 

 2εxz 
 zz ∂z 2 ∂z ∂y

 


 2ε 

yz

On définit ainsi un opérateur de dérivation S entre les petites déformations ε et les dépla-
cements u = {u, v, w}T tel que :

ε = Su (2.2)

 
∂
∂x
0 0
 
 0 ∂
0 
 ∂y 
 0 0 ∂ 
 ∂z 
avec S =  ∂ 
 ∂y ∂
0 
 ∂x 
 ∂ 0 ∂ 
 ∂z ∂x 
∂ ∂
0 ∂z ∂y

Loi de Hooke : La loi de Hooke est une loi de comportement du matériau qui lie les
petites déformations ε aux contraintes σ et qui s’applique aux matériaux isotropes (propriétés
du matériau identiques dans toutes les directions), homogènes. Cette loi est linéaire et peut
se mettre sous la forme suivante :

σ = Dε (2.3)

avec
Méthode des éléments finis (MEF) 43

 
λ + 2µ λ λ 0 0 0
 
 λ λ + 2µ λ 0 0 0 
 
 λ λ λ + 2µ 0 0 0 
 
D= 
 0 0 0 µ 0 0 
 
 0 0 0 0 µ 0 
 
0 0 0 0 0 µ

D est le tenseur d’élasticité avec µ le coefficient de résistance au cisaillement et λ le

cœfficient de compression, appelés également cœfficients de Lamé et qui sont définis par :

E νE
µ= λ=
2 (1 + ν) (1 − 2ν) (1 + ν)

où E et ν sont deux paramètres caractéristiques du matériau appelés respectivement mo-

dule de Young (exprimé en MPa) et coefficient de Poisson (sans unité). E caractérise la
raideur du matériau, il peut varier de 0.02 Mpa pour des tissus mous à plus de 200 000 Mpa
pour des aciers. ν caractérise la compressibilité du matériau, sa valeur est comprise entre
-1.0 et 0.5. ν=0.5 pour un matériau dit «incompressible».

Dans le cas d’un modèle statique ou quasi-statique, on ne prend pas en compte les forces
inertielles liées aux accélérations. La méthode des éléments finis basée sur la théorie de la
MMC en élasticité linéaire permet de construire une matrice de rigidité constante.

Les équations du modèle statique sont obtenues par application du principe des travaux
virtuels appliqué à chaque élément :

(e) (e)
δWint + δWext = 0 (2.4)

avec
Z
(e)
δWint =− δεT σ dv
Ωe

(e) T (e)
δWext = δu(e) fext

En utilisant les relations 2.1, 2.2, 2.3, on obtient :

(e) T T (e)
δWint = −δu(e) K(e) u(e) = δu(e) fint (2.5)
44 Chapitre 2 : Méthodes numériques et résolution du problème de contact

avec
Z
(e)
K = BT DB dv, B = SN(e)
Ωe

En utilisant l’équation 2.4, la relation 2.5 et par application du principe des travaux vir-
tuels, on obtient les équations d’équilibre :

(e)
K(e) u(e) = fext (2.6)

Le calcul de la matrice de raideur K(e) est élémentaire lorsque les fonctions de forme sont
linéaires car la matrice B est constante et indépendante des variables d’intégration. Dans le
cas contraire, le calcul de K(e) est discrétisé en des points particuliers appelés «points de
Gauss» [Zie77]. Les matrice K(e) sont indépendantes de u(e) et donc des déformations.
La table de connectivité associée au maillage, nous permet de construire une matrice de
raideur globale K(3n × 3n) où n représente le nombre de nœud du maillage. Le système
linéaire global s’écrit alors sous la forme matricielle :

KU = Fext (2.7)

Cette matrice K née du regroupement des K(e) est creuse puisque un nœud ne peut appar-
tenir qu’à un nombre très limité d’éléments. Elle est par contre symétrique, définie positive et
il est possible de la mettre sous la forme d’une matrice bande afin d’appliquer des méthodes
optimisées de résolution du système d’équations 2.7.

2.2.4.2 Élasticité non linéaire

Lorsque les déplacements élastiques dépassent 10% de la taille de l’objet considéré, on

est en élasticité non linéaire. Dans ce cas les déformations ne peuvent plus être représentées
par le tenseur de Cauchy ε. Il faut alors considérer le tenseur de Green-Lagrange E :

1¡ T ¢
E = ε+ ∇u ∇u (2.8)
2

∇ représente l’opérateur gradient.

Le tenseur de Green-Lagrange est invariant en rotation, il permet donc de prendre en
compte les grandes rotations sans qu’il soit nécessaire de faire appel à des approches coro-
tationnelles. Dans le cadre de la méthode des éléments finis, le tenseur de Green-Lagrange
inclut formellement les termes linéaires(BL u) et non linéaires par rapport aux déplacements
Méthode des éléments finis (MEF) 45

nodaux (BN L u).

³
1 ´
E = BL + BN L (u) u (2.9)
2

Dans le cas des lois élastiques ou hyper-élastiques, il existe un potentiel élastique W (ou
densité d’énergie) qui est une fonction scalaire du tenseur des déformations. Les contraintes
représentent la dérivée par rapport aux déformations de cette fonction scalaire. Dans le cas
particulier des modèles isotropes de Saint-Venant-Kirchhoff, on obtient une relation entre
contraintes et déformations similaire à la loi de Hooke :

S = DE (2.10)

D est le tenseur d’élasticité et S le second tenseur des contraintes de Piola-Kirchhoff.

Comme dans le cas linéaire, par application du principe des travaux virtuels, on obtient le
vecteur des forces internes :
Z
(e) ¡ ¢T
fint = − BL + BN L (u) S dV (2.11)
Ωe

On peut montrer qu’après calculs de dérivation [FPH06, FVFP06], on obtient :

(e)
∂fint
K(e) = = K(e) (e) (e)
e + Kσ + Ku (2.12)
∂u

avec
Z
Ke(e) = BTL DBL dV (2.13)
V0

Z
∂BN L
Kσ(e) = S dV (2.14)
V0 ∂u

Z
¡ ¢
Ku(e) = BTL DBN L + BTN L DBL + BTN L DBN L dV (2.15)
V0

On remarque que la matrice de raideur tangente K(e) associée à un élément e se décom-

(e)
pose en trois matrices, la matrice de rigidité élastique Ke qui est constante et identique au
(e) (e)
cas linéaire, la matrice de rigidité Kσ des contraintes initiales et la matrice de rigidité Ku
des déplacements initiaux.
46 Chapitre 2 : Méthodes numériques et résolution du problème de contact

2.3 Cinématique du contact

Dans ce paragraphe, nous introduisons les notions nécessaires pour prendre en compte le
contact dans la formulation d’un problème de mécanique. Pour une question de simplicité et
de clarté, nous considérons uniquement deux corps déformables en contact. La formulation
de problèmes mettant en jeu un plus grand nombre de solides se fait par simple extension de
la démarche présentée dans la suite.

2.3.1 Notion du problème de contact

F IG . 2.5: Repère local de contact

Deux corps déformables B α (figure 2.5), α = 1, 2 sont considérées. Chacun des deux
occupe, à l’instant initial, les ensembles fermés Ωα ⊂ R3 , constitués de points X α . Les
deux corps sont élastiques et peuvent subir de grandes déformations. La frontière Γα de
chaque corps est supposée suffisamment régulière pour qu’un vecteur normal unitaire vers
l’extérieur, noté par nα , puisse être défini à tout point P sur Γα . À chaque fois t ⊂ I,
où I = [0, T ] donne l’intervalle de temps correspondant à la procédure de chargement, la
fontière Γα du corps B α peut, en fait, être partitionnée en trois. Γαu où le déplacement uα
est donné, Γασ où la charge sur la frontière σ α est donnée, et Γαc est défini comme la zone
potentielle de contact sur laquelle B 1 et B 2 peuvent éventuellement venir en contact à un
instant tc . Cette nouvelle décomposition de la frontière de chacun des deux solides doit
vérifier :

Γα = Γαu ∪ Γασ ∪ Γcα ,

∅ = Γαu ∩ Γασ = Γασ ∩ Γcα = Γαu ∪ Γcα . (2.16)

Dans la suite, les deux corps B α , α = 1, 2 (l’un étant l’impacteur et l’autre étant l’obs-
Cinématique du contact 47

tacle) doivent subir une condition de non pénétration que l’on peut traduire par :

B1 ∩ B2 = ∅ (2.17)

Pour la partie sur la frontière, on peut écrire :

(
∅ si non contact
Γ1c ∩ Γ2c =
Γ1c = Γ2c si contact

Les deux corps sont en contact à l’instant tc si les éléments de frontières ont une partie
commune, sinon, aucun contact n’existe entre les deux corps.
En cette partie commune aux deux frontières, nous pouvons décrire à chaque instant t les
α
domaines de déplacement uα défini sur la Ω . Sur la surface de contact, une normale unique
n orientée vers B 1 (n ≡ n2 ) est définie, et le plan tangent, orthogonal à n dans le R3 , est
donné par T. Afin de construire une base locale orthogonale, deux vecteurs unitaires tx et ty
sont définis dans le plan T. Pour décrire l’interaction du contact frottant qui peut se produire
sur Γc , nous introduisons la vélocité relative en ce qui concerne B 2 .

u̇ = u̇1 − u̇2 (2.18)

où u̇1 et u̇2 sont les vitesses des déplacements respectivement de B 1 et B 2 . Soit r la force
de contact exercée par B 1 sur B 2 en P . Le principe d’action-réaction nous permet d’écrire
que le solide B 1 subit la réaction de contact −r. Dans le repère local du contact, u̇ et r peut
être décomposés de manière unique sous la forme :

u̇ = u̇t + u̇n n, u̇t ∈ T, u̇n ∈ R (2.19)

r = rt + rn n, rt ∈ T, rn ∈ R (2.20)

– u est le vecteur déplacement.

– r est le vecteur force de contact qu’exerce le corps B 2 au point P .

2.3.2 Conditions de contact de Signorini

Dans le repère local de contact, le loi de contact se présente sous la forme de trois re-
lations : une condition géométrique de la non-pénétration, une condition de non-adhérence
en statique et une condition complémentaire mécanique. Ces trois conditions sont appelées
conditions de Signorini.

1. condition de non-pénétration : elle permet de s’assurer qu’aucune particule de B 1 ne

pénètre dans le solide B 2 . Cette condition est définie dans la configuration d’équilibre
48 Chapitre 2 : Méthodes numériques et résolution du problème de contact

c’est à dire à l’instant (t + ∆t). On peut en déduire :

g(X) = min[(X1 − X2 ) · n] ≥ 0 (2.21)

Où g représente le «gap» défini à l’instant t et n la normale à l’instant (t + ∆t), et

Xα (t) = Xα (t = 0) + uα (2.22)

La position du vecteur X2 est trouvé comme la projection la plus proche du point

X1 ∈ Γ1c sur la surface Γ2c . Considérant h comme l’intervalle initial obtenu au début
de chaque pas de temps, on a :

h = (X1 − X2 ) · n ≥ 0 (2.23)

2. condition de non-adhésion : cette seconde condition traduit le fait que la réaction de

contact exercée par B 2 sur B 1 est dirigée vers l’extérieur de B 2 et à donc pour effet
de repousser B 1 . Avec la convention d’orientation adoptée pour n, la non-adhérence
s’écrit :

rn (x1t+∆t ) ≥ 0, (2.24)

3. condition complémentaire : cette dernière condition assure qu’il ne peut y avoir une
force de réaction entre B 1 et B 2 que s’il y a contact :

g(x1t+∆t ) · rn (x1t+∆t ) = 0. (2.25)

rn

contact

non-contact g

F IG . 2.6: Graphe des conditions de Signorini.

Cinématique du contact 49

Notons g0 la distance initiale entre X10 et X20 :

g0 = (X10 − X20 ) · n ≥ 0. (2.26)

Les conditions de Signorini peuvent alors être réécrites sous la forme :

g0 + un ≥ 0 , rn ≥ 0 , (g0 + un ) · rn = 0. (2.27)

Considérons maintenant le cas de deux solides initialement en contact sur une partie de
Γc . Sur cette partie de Γc , les conditions de Signorini se réduisent à :

un ≥ 0 , rn ≥ 0 , un · rn = 0. (2.28)

Ainsi, pour tout t ∈ I, les surfaces potentielles de contact Γαc (α = 1, 2) peuvent être
découpées en deux parties distinctes : + Γαc sur laquelle les deux corps sont effectivement en
contact et − Γαc pour laquelle il n’y a pas, à ce moment précis, de contact. Nous avons donc :

Γαc =+ Γαc ∪ − Γαc et + α − α

Γc ∩ Γc = ∅. (2.29)

En revanche Γαc , + Γαc et − Γαc changent avec le temps t et peuvent être vides pour certains
t ∈ I. Nous devons souligner que avec l’équation (2.28) traduit le fait que lorsque les deux
solides sont en contact, seule une séparation est autorisée. Elle ne permet pas d’exprimer la
mise en contact de nouvelles zones.
Dans le contexte de la dynamique, il est possible d’exprimer les conditions de Signorini
en + Γαc en terme de vitesse :

+ α
u̇n ≥ 0 , rn ≥ 0 , u̇n rn = 0 sur Γc . (2.30)

Il est important de noter que cette dernière formulation (équation 2.30) n’est valable que
si les deux solides sont initialement en contact.

2.3.3 Loi de frottement - modèle de Coulomb

Il est maintenant nécessaire de définir la loi qui permettra, lors des calculs, de déterminer
la partie tangentielle des réactions de contact. Dans ce travail, nous ne considérons que les
cas de frottement sec. En effet, le problème de frottement lubrifié fait intervenir un troisième
corps, en l’occurrence un liquide, entre les deux solides et sort de notre domaine d’étude.
Dans cette loi de frottement, la force tangentielle appliquée au solide en contact doit
dépasser un certain seuil pour qu’un glissement se produise. Le seuil à dépasser pour obtenir
un glissement relatif entre les deux solides en contact est une constante, fixée à priori et qui
50 Chapitre 2 : Méthodes numériques et résolution du problème de contact

dépend des matériaux. Afin de refléter ce qu’il a observé par expérimentation, Coulomb a
imposé une dépendance du seuil de glissement vis à vis de la réaction normale. Ainsi, la loi
de frottement de Coulomb s’écrit :

 Si krt k < µrn alors u̇t = 0 (adhérence)
rt (2.31)
 Si krt k = µrn alors ∃λ > 0 tel que u̇t = −λ . (glissement)
krt k

où µ est le cœfficient de frottement. La forme incrémentale est :

(
krt k < µrn si k∆xt k = 0
Coulomb(∆xt , rt ) ⇔ ∆xt
(2.32)
rt = −µrn k∆xt k si k∆xt k 6= 0

∆xt = xt (t + ∆t) − xt (t) mais du fait que le point Q soit le point projeté orthogonalement
de P, on a xt (t) = 0 ce qui implique que ∆xt = xt (t + ∆t) = xt . On a donc :
(
krt k < µrn si kxt k = 0
Coulomb(xt , rt ) ⇔ (2.33)
rt = −µrn kxxtt k si kxt k 6= 0

Comme les conditions de Signorini pour le problème de contact, le modèle de Coulomb

définissant le frottement constitue une difficulté majeure pour la résolution. Là aussi, la re-
lation entre la force tangentielle et le glissement relatif n’est pas définie par une fonction
régulière comme le montre la figure 2.7.

rt
glissement
µrn

adhérence
−u̇t

−µrn
glissement

F IG . 2.7: Loi de frottement de Coulomb

En couplant cette loi avec les conditions de Signorini, nous pouvons définir l’ensemble
fermé et convexe K µ des forces de réactions admissibles par :
© ª
K µ = r ∈ R3 tel que krt k − µrn ≤ 0 . (2.34)

La représentation graphique de l’ensemble K µ (appelé cône de Coulomb) est proposée sur

Cinématique du contact 51

la figure 2.8. Nous y retrouvons les trois statuts de contact possibles : non-contact, contact
adhérant et contact avec glissement.

rn

adhérence
glissement

rt2

non-contact rt1

F IG . 2.8: Cône de Coulomb

Il est possible de différencier le cœfficient de frottement qui s’applique lors de la phase de

glissement de celui qui sert de seuil lors de la phase d’adhérence. On parle alors de cœfficient
d’adhérence et de cœfficient de frottement. En effet, les expérimentations montrent que la
force nécessaire pour maintenir un glissement relatif entre deux solides est généralement
plus faible que la « force seuil » qui a été nécessaire pour le provoquer.
Par ailleurs, la loi de Coulomb ne constitue pas la seule alternative pour introduire le phé-
nomène du frottement. La tribologie est la science qui vise à étudier ce phénomène méca-
nique et de nombreux modèles peuvent être adaptés selon le cas étudié. Dans certains travaux
s’inspirant du domaine de l’élastoplasticité, nous trouvons de nouveaux modèles de frotte-
ment qui permettent un mouvement relatif réversible des surfaces de contact. Ce type de loi
prend en compte les déformations élastiques, au niveau microscopique, des aspérités de cha-
cune des surfaces de contact et des frottements anisotropes [Cur84, HC93, Wri02, HFdM04].

2.3.4 Loi de contact avec frottement

Nous considérons à présent la loi de frottement de Coulomb présentée précédemment que
nous allons coupler avec les conditions de Signorini afin d’obtenir une loi de contact avec
frottement complète. Cette loi fera intervenir les trois statuts de contact possibles : l’absence
de contact, l’adhérence et le glissement.
52 Chapitre 2 : Méthodes numériques et résolution du problème de contact

Une manière élégante de mettre en avant les différents statuts de contact consiste à écrire
la loi complète de contact avec frottement sous la forme d’une double boucle de type « Si ...
alors ... sinon.

si rn = 0 alors u̇n ≥ 0 !non-contact

sinon r ∈ K µ alors u̇n = 0 et u̇t = 0 ! adhérence
sinon (r ∈ ∂K µ et rn > 0)
½ ¾ (2.35)
rt
u̇n = 0 et ∃ λ > 0 tel que − u̇t = λ ! glissement
krt k
fin

où K µ et ∂K µ sont respectivement l’intérieur et le bord du cone de Coulomb.

La loi que nous venons de présenter permet de déterminer la vitesse relative u̇ lorsque l’on
connaît la force de réaction r. La loi de contact inverse, i.e. la relation r(−u̇) peut également
être mise sous la forme d’une double boucle.

si u̇n > 0 alors r = 0 !non-contact

sinon u̇ = 0 alors r ∈ K µ ! adhérence
sinon (nu̇ ∈ T − 0 alors o (2.36)
−u̇t
u̇n = 0 et rt = µrn k− u̇t k
! glissement
fin

2.4 Méthodes de résolution des problèmes de contact

L’analyse des problèmes de contact avec frottement est d’une grande importance dans
de multiples domaines d’ingénierie. Dans la littérature, de nombreuses méthodes ont été
développées pour traiter ce genre de problèmes :
– la méthode des multiplicateurs de Lagrange : dans cette méthode, les contraintes liées
au contact et au frottement sont introduites dans la formulation par le biais d’inconnues
supplémentaires (les multiplicateurs de Lagrange). Cette méthode permet de vérifier
exactement les conditions de contact mais elle a le désavantage d’augmenter la taille
du système à résoudre. Nous trouvons quelques exemples d’utilisation de cette méthode
dans [HTS+ 76, BN91, CTK91, PS98, SP98],
– la méthode de pénalisation : dans cette autre méthode, la loi de contact ainsi que la
loi de frottement sont régularisées par l’intermédiaire de deux cœfficients de pénalité
[CT71, TY73, KS81, CB86, Par89]. Cette méthode a pour principal avantage sa sim-
plicité d’implantation. Cependant, il faut souligner que cette méthode ne permet pas de
vérifier avec exactitude les conditions de contact. L’importance de la pénétration auto-
risée dépend du choix (difficile) des cœfficients de pénalité [Hun92] qui n’ont pas de
Méthodes de résolution des problèmes de contact 53

signification physique et qui doivent être déterminés pour chaque cas étudié,
– la méthode du Lagrangien augmenté : dans cette dernière méthode, les déplacements et
les réactions de contact sont déterminés simultanément [AC91, SL92, dSF91, Kla92,
dSF98, CKPS98]. Cette méthode présente l’avantage de fournir un résultat qui vérifie
exactement les conditions de contact. De plus, le paramètre utilisé dans la méthode
n’influe pas sur la précision du résultat mais seulement sur la vitesse de convergence
de l’algorithme.
La formulation du problème de contact avec la méthode de fonction de pénalisation
conduit à une seule résolution globale en déplacements car les forces de contact sont explici-
tement définies en fonction des déplacements. Par contre, dans les deux autres méthodes, les
forces de contact sont définies d’une manière implicite et constituent des inconnues supplé-
mentaires qui s’ajoutent aux inconnus en déplacement dans le système d’équations d’équi-
libre global des corps en contact. Deux approches sont alors possibles pour déterminer ces
deux inconnus. La première consiste à les résoudre simultanément d’une manière globale.
La deuxième approche dite «mixte» consiste à les résoudre séparément de manière succes-
sive. Les forces de contact sont tout d’abord déterminées d’une manière locale par réduction
des équations d’équilibre aux points de contact. Les déplacements sont alors obtenus par une
résolution globale du système d’équations d’équilibre en considérant les forces de contact
comme des forces externes [Fen95]. La méthode du Lagrangien augmenté permet de séparer
la non linéarité du comportement mécanique due aux grands déplacements, de la non linéa-
rité due aux lois de contact frottant. Elle permet aussi de créer une approche modulaire en
introduisant dans le processus de simulation un solveur des forces de contact indépendant du
solveur global des déplacements.

2.4.1 Formulation par des fonctions de pénalité

Le principe de la méthode consiste à régulariser les fonctions associées aux lois de contact
de manière à obtenir des graphes d’application.

2.4.1.1 Formulation implicite

Les conditions de Signorini et les conditions de la loi de Coulomb peuvent ainsi être
reformulées par [RF03] :

Signorini(un , rn ) ⇒ rn = ProjR+ (−ρn un ) avec ρn > 0 (2.37)

Coulomb(ut , rt ) ⇒ rt = ProjDµ (−ρt ut ) avec ρt > 0 (2.38)

54 Chapitre 2 : Méthodes numériques et résolution du problème de contact

Avec ProjR+ et ProjDµ des opérateurs de projection définis par :

( (
z
0 si z < 0 h kzk si z ∈
/ Dµ
ProjR+ (z ∈ R) = , ProjDµ (z ∈ R2 ) =
z si z ≥ 0 z si z ∈ int(Dµ )

Dµ représente le disque de Coulomb défini dans R2 et de rayon h (h = µrn ).

Cette formulation constitue une approximation des lois de contact dans laquelle les para-
mètres ρn , ρt représentent des raideurs ajoutées respectivement suivant la direction normale
et tangentielle au point de contact considéré. En effet, pour avoir une force de contact, il faut
nécessairement que un < 0 c’est-à-dire qu’il y ait interpénétration. Dans cette approche les
forces de contact agissent comme des forces de rappel pour compenser l’interpénétration et
le glissement. En augmentant la valeur de ces raideurs, on augmente la pente des courbes
représentées sur les graphes de 2.9 et on se rapproche ainsi des lois de contact. Cependant,
on ne peut pas excéder une certaine limite qui entraînerait un mauvais conditionnement nu-
mérique de la matrice tangente. En effet dans cette formulation les forces de contact sont
fonctions des déplacements contraints et contribuent au calcul des forces internes modifiant
ainsi le calcul de la matrice tangente KT du système.
Cette formulation se prête bien à la programmation mais elle nécessite la modification de
la matrice tangente qui change de forme suivant le statut de chaque point de contact (Non
contact, adhérence ou glissement) ce qui dans une résolution numérique itérative est très
coûteux.

F IG . 2.9: Régularisation des lois de contact

2.4.1.2 Formulation explicite

Les conditions de Signorini et la loi de Coulomb peuvent être aussi reformulées par :

Signorini(ũn , rn ) ⇒ rn = ProjR+ (−ρn ũn ) avec ρn > 0 (2.39)

Coulomb(ũt , rt ) ⇒ rt = ProjDµ (−ρt ũt ) avec ρt > 0 (2.40)

Méthodes de résolution des problèmes de contact 55

Cette formulation diffère de la précédente car elle ne prend pas en compte les déplace-
ments contraints mais les déplacements libres uniquement. Dans cette formulation les forces
de contact s’ajoutent aux calculs des forces externes, la matrice tangente reste inchangée.
Pour calculer les déplacements libres, on peut procéder de manière locale par projection des
équations d’équilibre au point de contact sans qu’il soit nécessaire de résoudre globalement
le système. Les paramètres ρn , ρt doivent être choisis avec précaution en relation avec les
raideurs des deux solides au point de contact. Des valeurs trop élevées entraînent une ré-
activité trop forte du point de contact et une instabilité numérique. Des valeurs trop faibles
entraînent des violations de contraintes (interpénétration, adhérence) trop importantes.

2.4.2 Formulations par des multiplicateurs de Lagrange

2.4.2.1 Formulation par adhérence

On bloque le point de contact dans les 3 directions normales et tangentielles (un =

0, ut = 0), on en déduit une force de blocage ~λ = λn ~n + ~λt dans le repère local (calculée
en prenant en compte les équations d’équilibre). On projette ensuite cette force de manière à
satisfaire les conditions de contact. Les conditions de Signorini et la loi de Coulomb peuvent
ainsi être reformulées par :

Signorini(un , rn ) ⇒ rn = ProjR+ (λn ) (2.41)

Coulomb(ut , rt ) ⇒ rt = ProjDµ (λt ) (2.42)

2.4.2.2 Formulation par glissement sans frottement

On bloque le point de contact uniquement dans la direction normale (un = 0 ), on en

déduit une force de blocage ~λ = λn ~n dans la direction normale (calculée en prenant en
compte les équations d’équilibre) et ũt le déplacement libre suivant les directions tangentes.
Les conditions de Signorini et la loi de Coulomb peuvent ainsi être reformulées par :

Signorini(un , rn ) ⇒ rn = ProjR+ (λn ) (2.43)

ũt
Coulomb(∆ut , rt ) ⇒ rt = −µrn (2.44)
kũt k

Le choix de l’une ou l’autre formulation dépend du problème considéré. Si c’est un pro-

56 Chapitre 2 : Méthodes numériques et résolution du problème de contact

blème dans lequel il y a plus d’adhérence que de glissement alors c’est la première méthode
qui est préférable sinon c’est la seconde. D’une manière générale la seconde méthode est
plus stable que la première. En effet dans la seconde méthode, l’effet dissipatif lié aux frot-
tements est garanti car la force tangentielle s’oppose aux mouvements tangentiels relatifs,
tandis que la première méthode relaxe la force tangentielle au risque de libérer de l’énergie
de déformation de manière trop brutale et d’entraîner des problèmes de stabilité numérique.
Ces deux méthodes ne simulent pas correctement les changements de statuts (glissement,
adhérence).

2.4.3 Formulations lagrangiennes augmentées

Nous allons voir maintenant deux méthodes qui considèrent la force de contact et les
déplacements contraints comme étant deux extrêmes de potentiels convexes non différen-
tiables. La première méthode introduite par Moreau [Mor88, Mor94], définie deux potentiels
non différentiables appelés «pseudo-potentiels». Le premier est associé aux conditions de Si-
gnorini et le deuxième aux conditions de la loi de Coulomb. La deuxième méthode introduite
par De Saxcé et Feng [dSF91], [dSF98], défini un unique potentiel non différentiable appelé
«bi-potentiel» associé à la fois aux conditions de Signorini et à la loi de Coulomb. L’avan-
tage de cette dernière méthode permet de mieux prendre en compte le couplage entre la force
normale et la force tangentielle au niveau des lois du contact lors de la résolution numérique.
Ces deux méthodes conduisent à l’élaboration d’opérateurs de projection à partir de consi-
dérations issues de l’analyse convexe et de la technique du lagrangien augmenté.

2.4.3.1 Pseudo-potentiels

Conditions de Signorini

0
³ 0
´
Signorini(xn , rn ) ⇔ rn ≥ 0 et ∀rn ≥ 0 rn − rn x n ≥ 0 (2.45)

Montrons qu’il y a bien équivalence :

0 0
Si rn > 0 alors il existe au moins un rn > 0 tel que rn < rn ce qui implique xn = 0 et
on vérifie donc bien Signorini.
0
Si rn = 0 alors rn .xn ≥ 0 ce qui implique xn ≥ 0 et on vérifie donc bien de nouveau
Signorini.
On peut montrer que l’équation 2.45 definit −xn comme un sous gradient d’une fonction
potentielle non différentiable construite à partir d’une fonction indicatrice.
Méthodes de résolution des problèmes de contact 57

On peut remplacer l’inégalité dans 2.45 par :

¡ 0 ¢ ¡ 0 ¢
rn − rn xn ≥ 0 ⇔ rn − rn ρxn ≥ 0 avec ρ > 0
¡ 0 ¢
⇔ rn − rn (rn − rn + ρxn ) ≥ 0 avec ρ > 0 (2.46)
¡ 0 ¢
⇔ rn − rn (rn − rn∗ ) ≥ 0 avec rn∗ = rn − ρxn

rn∗ constitue la force de contact normale augmentée.

0
³ 0
´
rn ≥ 0 et ∀rn ≥ 0 rn − rn (rn − rn∗ ) ≥ 0 ⇔ rn = ProjR+ (rn∗ ) (2.47)

¡ 0 ¢
En effet rn > 0 alors rn − rn peut être de signe quelconque ce qui implique rn = rn∗ par
0
contre si rn = 0 alors −rn∗ rn ≥ 0 ce qui implique 0 ≥ rn∗ .
Loi de Coulomb

0
³ 0
´
Coulomb(xt , rt ) ⇔ rt ∈ Dµ , ∀rt ∈ Dµ , rt −rt xt ≥ 0 (2.48)

Montrons qu’il y a bien équivalence :

Si rt est à l’intérieur du cône de Coulomb alors il existe un disque centré sur rt inclus
¡ 0 ¢
dans ce cône et on peut avoir toutes les directions rt −rt possibles ce qui implique que
xt = 0, ce qui définit bien l’adhérence. Si rt est au bord du disque de Coulomb (de rayon
non nul) alors l’inégalité traduit que −xt est dirigé suivant la normale extérieure impliquant
∃λ > 0, xt = −λ krrtt k , ce qui définit bien le glissement. On peut montrer que 2.48 définit
−xt comme un sous gradient d’un fonction potentielle non différentiable construite à partir
d’une fonction indicatrice.
On modifie l’inégalité dans 2.48 par :
¡ 0 ¢ ¡ 0 ¢
rt −rt xt ≥ 0 ⇔ − rt −rt ρxt ≤ 0
¡ 0 ¢ (2.49)
⇔ rt −rt (r∗t −rt ) ≤ 0 avec r∗t = rt − ρxt

En raisonnant de la même manière que précédemment, il est évident que lorsque rt est à
l’intérieur du disque de Coulomb alors nécessairement on a r∗t = rt . Lorsque rt est sur le
¡ 0 ¢
bord du disque cela veut dire que rt −rt est orienté suivant la normale extérieure du disque.
rt est donc le projeté de r∗t sur le bord du disque lorsque r∗t est à l’extérieur (voir figure 2.10).

2.4.3.2 Bi-potentiel

Nous avons vu dans les paragraphes précédents que pour chaque nœud de contact, la déter-
mination des réactions de contact consiste en la recherche d’un couple (u̇; r) qui, d’une part,
est solution de l’équation du mouvement et qui, d’autre part, vérifie l’équation de contact
avec frottement.
58 Chapitre 2 : Méthodes numériques et résolution du problème de contact

F IG . 2.10: Projection sur le disque de Coulomb


Signorini(xn , rn ) 
 0 ¡ 0 ¢
r ∈ K µ , ∀rt ∈ K µ , r −r v ≥ 0
+ ⇔ (2.50)

 avec v = (µ kxt k + xn )n + xt
Coulomb(xt , rt )

Montrons qu’il y a bien équivalence :

¡ 0 ¢
Si r ∈ (K µ ) alors r −r peut prendre n’importe quelle direction ce qui implique que :
( (
µ kxt k + xn = 0 xn = 0
v=0⇔ ⇔ (2.51)
xt = 0 xt = 0

Cela définit bien les conditions d’adhérence. Si r ∈ Bd(K µ ) − {0} alors −xt est ortho-
gonal au cône de Coulomb impliquant :
(
µ kxt k + xn = λµ
avec λ > 0 (2.52)
−xt = λ krrtt k

λ représente la norme de xt ce qui entraîne xn =0. On retrouve les conditions de contact

0 0
avec glissement. Si r = {0} alors on a ∀rt ∈ K µ , r (−v) ≤ 0 ce qui entraîne −v ∈ Kµ∗ où
Kµ∗ représente le cône polaire (ou cône dual) dont les bords sont orthogonaux aux bords du
cône de Coulomb (voir Figure 2.11) et qui est défini par Kµ∗ = {z ∈ R3 , µ kzt k + zn ≤ 0}.

−v ∈ Kµ∗ ⇔ µ k−xt k − xn − µ kxt k ≤ 0

(2.53)
⇔ 0 ≤ xn

Ceci est conforme avec la perte de contact. On modifie l’inégalité dans 2.50 par la même
Application numérique 59

adhérence
glissement
non contact

F IG . 2.11: Projection sur le cône de Coulomb

technique du lagrangien augmenté :

¡ 0 ¢ ¡ 0 ¢
r −r .v ≥ 0 ⇔ − r −rt ρv ≤ 0 avecρ> 0
¡ 0 ¢ (2.54)
⇔ r −r . (r∗ −r) ≤ 0 avec r∗ = r − ρv

Si r ∈ Bd(K µ ) alors r = r∗ , si r ∈ Bd(K µ ) − {0} alors (r∗ −r) est orienté suivant la
normale extérieure du cône de Coulomb. rt est donc le projeté orthogonal de r∗t sur le bord
du cône de Coulomb lorsque celui-ci est à l’extérieur de K µ et de Kµ∗ (voir Figure 2.11). On
a donc :

Signorini(xn , rn ) 

+ ⇔ r = ProjKµ (r∗ ) avec r∗ = r − ρv (2.55)


Coulomb(xt , rt )

2.5 Application numérique

Cette partie est consacré à la présentation d’une application numérique permettant de
mettre en évidence la robustesse et la précision du code de calcul FER (basé sur la méthode
du bi-potentiel) en comparaison avec ANSYS. Pour cela, nous proposons un problème 2D
connu de contact avec frottement.
Cet exemple est souvent étudié pour la validation de codes de calcul. Il s’agit d’un barreau
élastique de section rectangulaire comprimé sous l’action de forces verticales et latérales
contre un plan rigide (figure 2.12).
Cet exemple est intéressant car nous avons trois zones de contact différents dans un même
temps : la zone de non-contact (AB), la zone de glissement (BC) et la zone d’adhérence
(CD). En raison de la symétrie de la structure, seule la moitié du barreau est modélisée. Les
caractéristiques de cet exemple sont les suivants :
1. module de YOUNG E = 130000 Mpa ; cœfficient de Poisson ν = 0.2 ; cœfficient de
frottement µ=1
60 Chapitre 2 : Méthodes numériques et résolution du problème de contact

F IG . 2.12: Un barreau avec des conditions aux limites

2. conditions aux limites : ux = 0 sur ED (condition symétrique)

3. Charge : f=50 Mpa sur GE, F= 100 Mpa sur GA.

La moitié de la barre est modélisée avec 561 nœuds et 512 éléments quadrilatéraux li-
néaires. La figure 2.13 montre le maillage initial et déformé (avec un facteur d’amplification
de 500).
La figure 2.14 montre la distribution de la contrainte de cisaillement sur laquelle une
concentration est observée à la frontière de la zone de glissement (BC) et de la zone d’adhé-
rence (CD).
La pression normale et la pression tangentielle sur la zone AD sont données par ANSYS
et le code FER dans la figure 2.15, les résultats numériques correspondent bien.
La pression normale calculée par les deux codes est aussi donnée dans la figure 2.16.
Nous pouvons voir que les résultats de ANSYS dépendent du cœfficient de pénalité Kn
(inconnues à priori). Les résultats sont meilleurs avec un cœfficient de Kn = 106 comparé
avec Kn = 105 pour ANSYS.
Les figures 2.17 et 2.18 montrent le déplacement tangentiel et le déplacement normal (pé-
nétration) sur la zone AD entre ANSYS et le code FER. La remarque concernant l’influence
du cœfficient “Kn ” est toujours valable. Il faut noter que la pénétration est inévitable dans
ANSYS. Tous ces résultats illustrent que la stabilité numérique et la précision du code FER
sont meilleures qu’ANSYS. Un autre avantage est qu’il n’y a pas de paramètres numériques
à choisir (tels que les facteurs de pénalité).
Afin d’étudier l’influence du frottement sur le comportement de contact, des cœfficients
différents ont été considérés : µ = 0, 0.4, 0.8 et 1.3. Le tableau 2.1 montre des résultats ob-
tenus par les deux codes : FER et ANSYS. UA et VA sont respectivement le déplacement
horizontal et le déplacement vertical du point A. Le temps de CPU augmente avec le cœffi-
Application numérique 61

F IG . 2.13: Maillage initial et déformé

F IG . 2.14: Distribution de la contrainte de cisaillement

62 Chapitre 2 : Méthodes numériques et résolution du problème de contact

F IG . 2.15: Pression normale et pression tangentielle par ANSYS et FER

F IG . 2.16: Pression normale par ANSYS et FER

Application numérique 63

F IG . 2.17: Déplacement tangentiel sur la zone AD

F IG . 2.18: Déplacement normal sur la zone AD

64 Chapitre 2 : Méthodes numériques et résolution du problème de contact

cient de frottement, le code de FER est supérieur à ANSYS par rapport au temps de calcul
(tableau 2.1).

TAB . 2.1: Temps de CPU

µ Code AB BC CD UA (mm) VA (mm) CPU (s)
0.0 FER 0 40 0 0.025846 0 4
AN SY S 1 0 40 0 0.025838 -0.5E-04 16
AN SY S2 0 40 0 0.019072 -0.38823E-05 19
0.4 FER 1.25 35 3.75 0.019122 0.0000178 4
AN SY S1 0 34.9816 4.99963 0.018760 -0.000032 17
AN SY S2 0 36.231 3.74995 0.019072 0.0000039 20
0.8 FER 2.5 22.5 15 0.015574 0.0003721 5
AN SY S1 2.498 22.487 14.9996 0.01521 0.0003039 19
AN SY S2 2.498 23.7364 13.75 0.015548 0.0003389 22
1 FER 3.75 17.5 18.75 0.014597 0.0005576 6
AN SY S1 2.498 18.7381 18.7496 0.014197 0.0005188 20
AN SY S 2 3.747 18.7384 17.5 0.014571 0.000524 24
1.3 FER 3.75 15 21.25 0.013642 0.0008009 10
AN SY S1 3.747 14.99 21.2497 0.013193 0.0007985 22
AN SY S2 3.747 16.2393 20 0.013598 0.0007858 26
AN SY S : facteur de pénalité Kn = 10 ; AN SY S : facteur de pénalité Kn = 106
1 5 2

2.6 Conclusions
Au cours de ce deuxième chapitre, nous avons abordé la résolution des problèmes de
contact avec frottement. Dans un premier temps, nous avons rappelé le principe de la mé-
thode des éléments finis. Puis nous avons présenté la cinématique du contact, la loi de contact
unilatéral et la loi de frottement. Ensuite nous avons effectué une brève revue des méthodes
fréquemment utilisées pour traiter le contact dans la résolution. Enfin, nous avons présenté
les résultats numériques donnés par deux algorithmes différents (méthode de pénalité et mé-
thode du bi-potentiel) sur un cas simple de contact avec frottement.
 65

Chapitre 3

Étude des effets mécaniques dans une

cellule de la pile

3.1 Introduction
La force de serrage joue un rôle important sur les effets mécaniques au sein de l’assem-
blage de la PAC, en particulier la force de contact, l’état de contact et la déformation de la
GDL qui influent sur la performance de la PAC et sur sa durée de vie. La force de serrage
peut être considérée comme un chargement statique par rapport au chargement dynamique
en situation de transport. La modélisation numérique de la précompression peut identifier la
distribution de la force de contact, la zone de contact et la déformation aux interfaces de la
PAC. Le maillage, les conditions aux limites et la précompression représentent le cas réel.

3.2 Simplification et propriétés du modèle numérique

3.2.1 Description et simplification du modèle

Pendant l’assemblage d’un empilement de la pile de PEM, des composants de la structure
sont mis en contact. La force de serrage par boulons tirants se transmet au sein de la pile
par des points de contact aux interfaces. Le contact dans la structure multi-contact de la PAC
existe principalement aux interfaces entre la plaque GDL et la plaque BPP et entre la plaque
GDL et la plaque MEA. La distribution de la pression de contact n’est pas uniforme aux
interfaces à cause de petits dents présentes sur les BPPs. La force appliquée extérieure et la
configuration des plaques aux contacts sont importantes pour la zone de contact. Le problème
de contact dans la PAC est impossible à résoudre de manière analytique mais il est réalisable
avec la méthode des éléments finis. Dans nos études, les matériaux des composants des piles
sont linéaires et élastiques, le modèle de la PAC a été développé en utilisant des hypothèses
66 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

suivantes :
1. Toutes les propriétés des matériaux sont isotropes,
2. La force de précompression appliquée est imposée par un déplacement en réalisant la
parfaite condition d’étanchéité des joints d’étanchéité dans la PAC,
3. Considérant la symétrie de la structure et des conditions aux limites, la modélisation
de la structure peut être simplifiée (figure 3.1),
4. La contrainte dans le modèle n’atteint pas la limite d’élasticité, c’est donc une étude
purement élastique,
5. L’étude se limite aux comportements mécaniques des éléments solides constituant la
pile. On ne considère jamais la réaction chimique ni les phénomènes fluidiques.

3.2.2 Propriétés mécaniques du modèle

Divers types des piles à combustible existent dans différents domaines : stationnaire, mo-
bile, véhiculaire, aéronautique, etc. Les propriétés des piles dépendent du domaine d’utilisa-
tion. Dans nos études, le type des piles de PEM est constitué par la BPP de SGL [All09], la
GDL de TGP-H-120 [Jap08] et la MEA de Nafion 117 de Dupont [DuP06] intervenant dans
la domaine du transport.

TAB . 3.1: Propriétés mécaniques de la plaque BPP, la plaque GDL et la plaque MEA

Propriétés Unités BPP GDL MEA

Densité Kg/m3 1980 450 1980
Module de Young Mpa 13000 6.13 64
Cœfficient de Poisson 0.21 0.09 0.09
Porosité initiale 78%
2
Surface mm 140 × 140 100 × 100 100 × 100

La plaque BPP en graphite vient de la plaque bipolaire de SGL Inc, qui est un matériaux
spécialement développé pour la PAC avec super-conducteur du courant électrique. Elle est
rigide par rapport à la plaque GDL et à la plaque MEA.
La plaque GDL du type TGP-H-120 vient de TORAY Inc. Elle est en papier carbone avec
une porosité initiale élevée, un bon passage gazeux et une résistance électrique basse. De
plus, elle est anti-corrosive et son coût est économique.
La plaque MEA est une membrane des électrodes comme la Nafion 117 de Dupont Inc,
qui est une membrane commerciale utilisée dans les PACs. La membrane est mince et très
conductrice afin de réduire la résistance ionique.
Modèle des éléments finis de la PAC 67

3.3 Modèle des éléments finis de la PAC

Le modèle de la PAC a été réalisé par la méthode des éléments finis, en particulier la
création du modèle géométrique, la définition des paramètres, le maillage, la précompression
et les conditions aux limites.

3.3.1 Description du modèle

La structure lamellaire d’un empilement des piles est constitué par plusieurs cellules uni-
taires pour offrir une puissance suffisante. Chaque cellule est une structure symétrique (fig
3.1). De plus, la structure de la demi-cellule est aussi symétrique verticalement (figure 3.1).
Parce que la symétrie de la condition aux limites et en particulier le déplacement appliqué
sur la plaque BPP est aussi symétrique, un modèle symétrique a donc été pris en compte.
Enfin, un modèle 2D est utilisé par rapport au modèle 3D car avec la symétrie de la structure
et la symétrie de la charge, aucune différence se trouve dans le profondeur de la structure de
la PAC.

F IG . 3.1: Modèle symétrique 2D de la PAC

Un modèle paramétrique est considéré afin d’éviter des travaux répétitifs de la modélisa-
tion, de facilement modifier des paramètres dans le but d’étudier l’influence des paramètres
mécaniques sur la performance de la pile. De plus le modèle paramétrique peut servir de
base pour de futures études.
68 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

La figure 3.2 montre les paramètres définis dans le modèle de contact de la PAC. Des
paramètres sont définis comme l’épaisseur de la plaque BPP (TBPP ), la largeur de la plaque
BPP (WBPP ), la largeur des canaux de la plaque BPP (WCN ), le profondeur des canaux
(TCN ), la largeur des dents (WRIB ) , le rayon du congé des dents (RRIB ) de la plaque BPP,
l’épaisseur de la plaque GDL (TGDL ), la largeur de la plaque GDL (WGDL ), l’épaisseur de la
plaque MEA (TMEA ) et la largeur de la plaque MEA (WMEA ). La largeur de la plaque BPP
(WBPP ) est donc égale à la largeur de la plaque GDL (WGDL ) et à la largeur de la plaque
MEA (WMEA ). La largeur de la plaque BPP (WBPP ) est le somme de la largeur des canaux
de la plaque BPP (WCN ) et de la largeur des dents (WRIB ).

F IG . 3.2: Modèle paramétrique de la PAC

Les propriétés matérielles et géométriques [JL09] sont définies dans le tableau 3.2.
Modèle des éléments finis de la PAC 69

TAB . 3.2: Propriétés et géométries du modèle de contact paramétrique de la PAC

Paramètres Définitions Valeurs Sources

TM EA Épaisseur de MEA (mm) 0.05 [DuP06]
WM EA Largeur de MEA (mm) 1.6 [DuP06]
EM EA Module de Young de MEA (Mpa) 64 [DuP06]
PM EA Cœfficient de Poisson de MEA 0.09 [DuP06]
TGDL Épaisseur de GDL (mm) 0.375 [Jap08]
WGDL Largeur de GDL (mm) 1.6 [Jap08]
EGDL Module de Young de GDL (Mpa) 6.13 [Jap08]
PGDL Cœfficient de Poisson de GDL 0.09 [Jap08]
TBP P Épaisseur de BPP (mm) 1.5 [All09]
WBP P Épaisseur de BPP (mm) 1.6 [All09]
WCN Largeur des canaux de BPP (mm) 0.8 [All09]
TCN Profondeur des canaux de BPP (mm) 1 [All09]
WRIB Largeur des dents de BPP (mm) 0.8 [All09]
RRIB Rayon des congés des dents de BPP (mm) 0.05 [All09]
EBP P Module de Young de BPP (Mpa) 13000 [All09]
PBP P Cœfficient de Poisson de BPP 0.21 [All09]

3.3.2 Maillage

Dans une modélisation par éléments finis, le maillage est une étape importante. En effet,
les résultats d’une simulation numérique dépendent du maillage et en particulier du type et
de la taille des éléments. Un maillage très fin permet d’obtenir des résultats précis mais le
coût du calcul peut être important. Il faut donc optimiser le maillage et trouver un compromis
entre la précision des résultats et le coût de la ressource du calcul numérique.

Dans le modèle de la PAC, l’élément quadrilatéral à 8 nœuds est utilisé. Le maillage est
plus fin aux interfaces dans le but de simuler correctement le contact. De plus, un maillage
fin est également utilisé au niveau des congés des dents de la BPP car la zone de contact
peut potentiellement s’étendre sous la précompression. Pour cette raison, la fine taille des
éléments à l’interface est imposée à 0.025 mm et la taille des éléments des autres zones à
0.1 mm. Ainsi, le total nombre des éléments dans le modèle est 964, le nombre des nœuds
est 3122.
70 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

F IG . 3.3: Maillage du modèle éléments finis

3.3.3 Définition de l’objet “zones de contact”

Pour tout problème de contact, nous définissons dans le fichier de données sources une
liste de zones potentielles de contact. La schéma 3.4 montre quels sont les éléments consti-
tutifs de chacune de ces zones de contact. Cette façon de déclarer au préalable les zones
potentielles de contact a un désavantage : il faut avoir une idée assez précise de ce qui va se
passer au cours de la simulation. Dans le cas contraire, il est nécessaire de déclarer toutes les
zones de contact possibles.

Zone de contact

Liste des nœuds de contact « impacteurs »

Liste des éléments cibles

Coefficient de frottement µ

Distance d’alerte

F IG . 3.4: Composantes d’une zone de contact

Par ailleurs, la désignation de l’impacteur et de l’obstacle peut avoir son importance dans
le cadre des géométries discrétisées utilisées dans la méthode des éléments finis. En effet,
dans le cas de la figure 3.5, si le solide B 1 est désigné comme impacteur, l’algorithme de
détection de contact renseignera que le nœud X est à l’intérieur de B 1 . Dans le cas contraire,
Modèle des éléments finis de la PAC 71

où B 2 est désigné comme étant l’impacteur, l’algorithme ne détecte aucun nœud de B 2 à

l’intérieur de B 1 et renseigne qu’il n’y a pas de contact. Cette difficulté peut aisément être
traitée en effectuant deux passages de l’algorithme de détection de contact et en inversant
les rôles des deux solides entre les deux passages. Toutefois, il est bon de noter qu’il est
préférable de choisir l’objet ayant la surface la plus régulière comme obstacle.

B1

B2

F IG . 3.5: Erreur de détection du contact

Il existe plusieurs autres traitements du contact. En particulier, il est possible (et avan-
tageux) d’adopter un formalisme de type nœud-contre-nœud sous l’hypothèse des petites
perturbations. D’autre part, l’utilisation d’une méthode segment-contre-segment en 2D et
segment-contre-facette en 3D permet d’éviter l’utilisation d’un algorithme à double passage
pour traiter le problème défini sur la figure 3.5 [PL04].

Dans le modèle de la PAC (figure 3.1), le contact existe aux interfaces entre la plaque
GDL et la plaque BPP et entre la plaque GDL et la plaque MEA. Le contact entre la plaque
GDL et la plaque MEA est de surface plane à surface plane. Le contact entre la plaque GDL
et la plaque BPP est différent. La zone de contact à l’interface entre la plaque GDL et la
plaque BPP est constituée par deux parties comme indiqué dans la figure 3.6 : une partie est
le contact de plan à plan entre la plaque GDL et la partie plane inférieure de la dent de la
plaque, l’autre partie est le contact de courbe à plan entre la plaque GDL et le congé de la
dent de la plaque BPP. La plaque GDL est comme l’obstacle car sa surface est régulière par
rapport à la plaque BPP avec un congé. Ainsi, la surface inférieure de la dent de la plaque
BPP est définie comme l’impacteur.
72 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

F IG . 3.6: La zone de contact entre les plaques GDL et BPP

3.3.4 Conditions aux limites

Le modèle numérique des éléments finis doit être assortie de conditions aux limites afin
de d’empêcher tout mouvement d’ensemble et refléter le cas réel. La figure 3.7 montre les
conditions aux limites utilisées. Le chargement de la précompression est imposé sur la ligne
(A) à travers un déplacement vertical dont la valeur est obtenue par le niveau de compression
des joints soumis à une parfaite condition d’étanchéité des piles. Utilisant la symétrie de la
structure, une condition de symétrie est imposée sur la ligne (B), c’est-à-dire, on impose un
déplacement horizontal nul à tous les nœuds de (B). La symétrie horizontale de la cellule
nous permet d’appliquer un déplacement vertical nul sur les nœuds de la ligne (C) corres-
pondant à la ligne de symétrie au milieu de la plaque MEA. On pourrait encore considérer
seulement la demi-structure de celle présentée sur la figure 3.7, ce que l’on fera plus tard.

3.3.5 Chargement de précompression de la PAC

Comme nous l’avons déjà énoncé précédemment dans le chapitre 1, la parfaite étanchéité
est une condition indispensable pour le fonctionnement normal de la PAC. Un chargement
de précompression est appliqué afin de garantir cette condition. Le chargement est effectué
par une force de serrage par boulons tirants dans la PAC. Un niveau de compression suffisant
des joints d’étanchéité dépend du chargement appliqué. Un niveau de 25% de compression
des joints d’étanchéité est pris en compte dans nos études d’après les expériences de Parsons
et al. [Par09]. L’épaisseur des joints étant identique à la sommes des épaisseurs de la plaque
MEA et de la plaque GDL (figure 1.15), l’épaisseur de joints d’étanchéité est de 0.8 mm
(2 ∗ TGDL + TM EA ). Le déplacement imposé reflétant la précompression dans le modèle est
0.2 mm. Pour le demi-modèle symétrique de la PAC, le déplacement appliqué est UY = 0.1
Études des contacts dans la PAC 73

F IG . 3.7: Conditions aux limites et la charge

mm sur la plaque BPP (figure 3.7). Parce que le joint d’étanchéité n’influence ni l’état de
contact, ni la déformation, les joints d’étanchéité ne sont pas modélisés dans le modèle de la
PAC.
Par ailleurs, l’application d’un déplacement “de serrage”, plutôt qu’une force de serrage
présente l’avantage de garantir l’uniformité de l’étanchéité et d’éviter toutes les pertes flui-
diques dans la pile.

3.4 Études des contacts dans la PAC

3.4.1 Étude de deux zones de contacts (GDL/BPP et GDL/MEA)

Dans un premier temps, un modèle de contact de la PAC avec deux zones de contact sous
précompression est présenté, ce sont les zones de contact entre les plaques GDL et MEA, et
entre les plaques GDL et BPP.
La figure 3.8 montre une distribution de la déformation du modèle. Sur la figure 3.9, sur
laquelle sont représentées uniquement les déformations de la plaque GDL, nous pouvons
trouver que le résultat est tout à fait symétrique. La déformation de la plaque GDL est plus
importante que les plaques MEA et BPP car son module d’élasticité est faible par rapport
aux autres. La zone de la GDL sous la dent de la plaque BPP est évidemment compressée
et la déformation y est négative. La zone initialement en contact reste en contact sur toute la
74 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

partie inférieure de la dent. En dehors de cette zone, on remarque que la surface supérieure
de la GDL remonte légèrement, le gonflement de la GDL influence l’état de contact de deux
interfaces (GDL et BPP, GDL et MEA), ce qui conduit la plaque GDL à s’approcher la
plaque BPP et à s’éloigner de la plaque MEA.

F IG . 3.8: Déformation totale du modèle F IG . 3.9: Déformation de la GDL

Sur la figure 3.10 du déplacement normal, une perte de contact apparaît aux extrémités de
l’interface entre la plaque GDL et la plaque MEA. La figure 3.11 montre qu’un décollement
apparaît à l’interface de la plaque GDL et de la plaque MEA. On retrouve ce décollement
lorsqu’on trace les déplacements normaux de la plaque GDL et de la plaque MEA à l’inter-
face. De plus, les figures 3.12 et 3.13 représentant la distribution de la pression sur l’interface
des plaques GDL et MEA confirment qu’il y a une zone de pression nulle où le décollement
apparaît.

F IG . 3.10: Une séparation de contact F IG . 3.11: déplacement normal

Le décollement entre la plaque MEA et la plaque GDL conduit à l’augmentation de la

résistance de contact ce qui est très dommageable pour la performance de la pile. En plus,
le gonflement de la GDL dans le canal de la plaque BPP, même s’il est faible, diminue le
Études des contacts dans la PAC 75

F IG . 3.12: Pression nulle F IG . 3.13: Distribution de la pression

volume du canal, peut gêner la circulation fluidique et diminuer l’efficacité de la réaction de

la PAC. Ainsi une technologie comme le recuit thermique (Thermal annealing) [HWDA07]
doit être appliquée à l’interface entre la plaque GDL et la plaque MEA afin d’éliminer la sé-
paration de contact et diminuer le gonflement de la GDL dans le canal de la BPP. Le recuit est
réalisé sous une température optimale avec une compression suffisante. Des catalyseurs sont
insérés entre la plaque GDL et la plaque MEA afin de coller les deux plaques et diminuer la
résistance de contact. Pour cette raison, un collage des deux plaques (GDL et MEA) a éga-
lement été toujours réalisé dans le modèle numérique de la PAC. Le collage permet d’éviter
le mouvement relatif des nœuds à l’interface entre la plaque GDL et la plaque MEA, et de
réaliser un contact parfait entre la plaque GDL et la plaque MEA. Dans les études suivantes,
on prend en compte ce collage.

3.4.2 Étude d’une zone de contact (GDL/BPP) avec l’angle droit

On considère maintenant une seule zone de contact entre la plaque GDL et la plaque BPP.
Ce contact a lieu sur l’interface entre la surface supérieure de la plaque GDL et la surface
inférieure de la dent formant le canal qui est gravé sur la plaque bipolaire. Ainsi la topologie
de la BPP a une rôle dans ce contact interfacial. En effet une dent rectangulaire de la BPP
peut conduire au problème mécanique comme la fissure sur la plaque GDL observée par Su
et al. [SLC+ 08]. La fissure de la plaque GDL peut être provoquée par la concentration de
contraintes sur la surface supérieure en contact avec les angles droits des dents de la plaque
BPP.
La figure 3.14a montre la distribution de la contrainte σV onmises dans les plaques GDL,
BPP et MEA. Sur la figure 3.14a, nous retrouvons bien un résultat symétrique, le modèle et
le chargement étant symétriques. La contrainte maximale apparaît évidemment aux coins des
dents rectangulaires qui sont en contact. Avec une concentration de contrainte de 10.1 Mpa,
cela peut conduire à la fissure mécanique à l’interface de la GDL [BGGM08]. La figure 3.14b
76 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

donne la distribution de la pression de contact à l’interface. Nous constatons que la pression

est la plus élevée aux angles des dents rectangulaires, elle est de 9.02 Mpa. En revanche, en
dehors de cette zone, la pression peut être considérée comme uniforme le long de la zone de
contact. Par rapport à la pression maximale, la plus petite pression de 1.5 Mpa sur l’interface
apparaît sur l’axe de symétrie (c’est-à-dire, X = 0 sur la figure 3.14b).

(a) (b)

F IG . 3.14: (a) Concentration de contrainte et (b) pression de contact

Le modèle présenté ici met en évidence une concentration de contraintes due à la forte
pression de contact aux angles des dents. Cette concentration peut être la cause de fissure
mécanique dans la plaque GDL comme celle étudiée par Su et al. [SLC+ 08]. Stanic et al.
[SH04] ont étudié que l’endommagement au sein de la pile intervient souvent dans la zone
de la concentration de contraintes, une pression et contrainte élevée à l’interface causent des
dommages sur la plaque GDL et raccourcissent la durée de vie de la PAC.
De plus, la non uniformité de la pression de contact sur la plaque GDL est aussi un para-
mètre important pour la performance de la PAC. Comme le montre de l’experience de Wang
et al [WSZ08], plus la distribution de la pression est uniforme, plus l’état de contact et la
résistance de contact sont uniformes, plus la performance des piles est bonne. Il faut donc
chercher un moyen permettant de réduire cette concentration de contraintes. L’angle droit
représente une singularité géométrique responsable de cette pression élevée. En remplaçant
cet angle par un congé, la concentration devrait s’atténuer.

3.4.3 Modélisation avec congé au niveau de la dent

Pour éviter les problèmes de concentration de contraintes avec la dent rectangulaire de la
BPP évoqués précédemment, on remplace donc l’angle droit par un congé de rayon RRIB .
Études des contacts dans la PAC 77

Le but de ce congé au niveau de la dent de la plaque BPP est de diminuer la concentration

de contrainte et d’améliorer l’uniformité de la pression à l’interface entre les plaques GDL
et BPP. Les résultats obtenus seront comparés avec les résultats d’une zone de contact avec
une dent rectangulaire. L’aire de contact à l’interface entre les plaques GDL et BPP et la
déformation de la plaque GDL vont aussi être présentées dans cette partie.
Pour affiner la modélisation à l’interface dent de la BPP/GDL, on utilise toutes les symé-
tries pour réduire au minimum la géométrie (figure 3.15a).

(a) (b)

F IG . 3.15: (a) Modèle du contact actuel et (b) modèle symétrique

Sur la figure 3.15b, les conditions aux limites n’ont pas changé, sauf une condition de
symétrie (Ux=0) verticalement sur B, le déplacement sur la surface supérieure de la plaque
BPP n’influence pas le résultat numérique car la déformation de la plaque BPP peut être
négligée grâce à sa rigidité élevée. Aussi pouvons nous appliquer ce même déplacement sur
la partie supérieure de la dent (A). Les propriétés et géométries du modèle de contact de la
PAC restent les mêmes (tableau 3.2). Le maillage contient 892 nœuds et 266 éléments (figure
3.16a).

3.4.3.1 Contrainte dans les plaques

Afin de garantir la durée de vie des composants de la PAC, une étude la contrainte maxi-
male sur les composants fragiles, comme les plaques MEA et GDL est aussi réalisée.
La figure 3.17a montre la distribution de la contrainte de Von-mises dans les plaques.
Nous constatons que la contrainte maximale (zone rouge : 3.352 Mpa) est présente au niveau
des congés des dents de la plaque BPP. Par rapport aux contraintes (10.1 Mpa) avec une
78 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

dent rectangulaire (figure 3.14a), la concentration de contraintes est divisée par 3. De plus, la
contrainte maximale de la plaque MEA (figure 3.17b) est de 2.87 Mpa. Un congé est permet
donc d’éviter la concentration de contraintes dans la plaque GDL mais aussi la concentration
de contraintes dans la plaque MEA (Stanic et al. [SH04]).
Pour une température de 23◦ C avec une humidité de 50%, la limite élastique de la plaque
MEA est de 43 Mpa. Dans le cas de l’eau infiltré, la limite élastique de la plaque MEA est
de 34 Mpa ; dans le cas de l’eau infiltré et une température de 100◦ C, la limite élastique
s’abaisse à 25 Mpa. Par contre, la limite élastique de la plaque GDL de la série de TGP-H
de TORAY est 40 Mpa. Les contraintes dans les plaques GDL et MEA de notre modèle sont
encore loin de la limite élastique des matériaux. C’est donc bien une étude élastique.

(a) (b)

F IG . 3.16: (a) Maillage et (b) Maillage fin

Pour étudier l’influence de la taille du maillage sur la précision de la contrainte, nous

avons comparé la contrainte avec un maillage différent. Sans changer la géométrie et les
conditions aux limites, on affine le maillage. La taille des éléments passe de 0.0125 mm dans
la zone de contact contre 0.025 mm avant. Dans les autres zones, la taille des éléments est
diminuée : 0.05 mm contre 0.1 mm avant, le nouveau nombre d’éléments dans le modèle
numérique est 623 et le nombre des nœuds est 2039.
La distribution des contraintes avec le maillage fin (figure 3.18a) montre que la contrainte
maximale de la plaque GDL est de 3.359 Mpa et la contrainte maximale de la plaque MEA
est de 2.88 Mpa. Si on compare avec les valeurs précédentes, la différence des contraintes est
de 0.01 Mpa, l’ancien maillage était donc déjà suffisant pour la précision attendue. Pour les
études suivantes, on conserve donc le maillage initial avec les tailles d’éléments suivantes :
0.025mm à l’interface entre la GDL et la BPP et 0.1 mm dans les autres zones du modèle de
Études des contacts dans la PAC 79

(a) (b)

F IG . 3.17: (a) Contraintes de trois plaques et (b) contraintes de la MEA

(a) (b)

F IG . 3.18: (a) Contraintes de trois plaques et (b) contraintes de la MEA avec maillage fin
80 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

contact de la PAC.

3.4.3.2 Pression de contact

La distribution et l’uniformité de la pression de contact à l’interface vont être étudiées.

La figure 3.19a montre la distribution de la pression de contact entre la plaque GDL et
la plaque BPP. Cette distribution dépend de la zone de contact (figure 3.6), qui est plan à
plan puis courbe à plan. La pression le long de la zone de contact plan à plan est presque
uniforme, alors qu’elle varie dans la zone de contact correspondant au congé.
La pression maximale est toujours localisée au niveau du congé de la dent de la plaque
BPP (figure 3.19b) car le congé ne permet pas d’annuler complètement l’effet de l’angle
droit. Toutefois, la valeur maximale est inférieure (3.44 Mpa), contre 9.02 Mpa (figure 3.14b)
en présence de l’angle droit.
Sur la zone de contact coïncidant avec l’axe de symétrie (x=0), la pression est de 1.86
Mpa ce qui représente la valeur de la pression uniforme sur la zone de contact surface à
surface. Par rapport à la pression maximale obtenue dans le modèle avec congé, la différence
de pression à l’interface est de 1.58 Mpa. L’uniformité de la pression à l’interface a été
améliorée par rapport au modèle avec l’angle droit pour lequel la différence de pression était
de 7.5 Mpa. La présence d’un congé à la place de l’angle droit sur la dent de la plaque BPP
permet donc de diminuer la concentration de contraintes et d’améliorer l’uniformité de la
pression à l’interface.
De plus sur la zone de contact courbe à plan, correspondant au congé de la dent (figure
3.19b), la pression augmente jusqu’à la zone jaune et rouge puis diminue régulièrement
comme la théorie de Hertz peut le prévoir. En plus, le changement de la pression au congé
aussi montre que l’état de contact peut aussi changer au congé de la dent de la plaque BPP.

3.4.3.3 Aire de contact (longueur de contact)

Nous avons vu au chapitre 1 que la performance de la pile est liée à la conductivité élec-
trique. Or, la capacité de la conductivité électrique entre la GDL et la BPP dépend de la
résistance de contact et l’aire de contact à l’interface. Un contact parfait entre la plaque GDL
et la plaque BPP permet d’obtenir une petite résistance et une bonne performance de la PAC.
Si l’on regarde maintenant l’étendue de la zone de contact : avec un angle droit sur la dent
de la plaque BPP, le contact a lieu seulement sur la surface rectiligne inférieure et s’arrête au
niveau de l’angle droit. Avec un congé à la place de l’angle droit, le contact va s’étendre le
long du congé, au delà de la surface rectiligne. L’étendue de la zone de contact est donc plus
grande.
L’aire de contact peut être obtenue en additionnant toutes les longueurs des éléments
0
en contact. La demi longueur maximale de la zone de contact avec le congé est WRIB =
Études des contacts dans la PAC 81

(a) (b)

F IG . 3.19: (a) Pression à l’interface et (b) pression agrandie aux fonds des congés

WRIB +(π−2)×RRIB
2
= 0.42854 mm. Le contact est insuffisant si la zone de contact n’atteint
pas cette dimension.
Les éléments en contact sont détectés selon la condition Signorini. Sous une précompres-
sion de 0.1 mm, l’aire de contact de tous les éléments de contact à l’interface est indiquée
dans le tableau 3.3 :

TAB . 3.3: Aire de contact à l’interface entre la plaque GDL et la plaque BPP

No. éléments Longueur de contact (mm) No. éléments Longueur de contact (mm)
291 0.29163E-01 299 0.29166E-01
292 0.29166E-01 300 0.29167E-01
293 0.29167E-01 301 0.29167E-01
294 0.29167E-01 302 0.29167E-01
295 0.26177E-01 303 0.29167E-01
296 0.29165E-01 304 0.26174E-01
297 0.29164E-01 305 0.26175E-01
298 0.29167E-01

La somme des longueurs des éléments de contact est de 0.42852 mm ce qui est très proche
de la valeur théorique. Avec une précompression de 0.1 mm, la plaque GDL et la plaque
BPP sont donc parfaitement en contact à l’interface, la structure lamellaire de la PAC a une
résistance de contact minimale et une meilleure performance.
82 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

3.4.3.4 Déformation de la plaque GDL

Le congé peut améliorer la concentration de contrainte à l’interface, et la précompression

élevée permet d’obtenir une petite résistance de contact dans la structure de la PAC. Une phé-
nomène mécanique inévitable est la déformation de la GDL amenée par la précompression
de la PAC. La déformation élevée peut conduire à une faible perméabilité de la plaque GDL
ce qui influence la performance de la PAC. Parce que la capacité des GDLs à transporter
des gaz réactifs est mesurée par la perméabilité, une petite perméabilité de la plaque GDL
conduit à une alimentation des gaz insuffisante et une perte de performance de la PAC. Ainsi
une petite déformation de la plaque GDL est préférable pour la PAC.
La non uniformité de la distribution de la pression conduit à une distribution non uni-
forme de la déformation à l’interface ce qui conduit à la non uniformité de la porosité et de
la perméabilité de la GDL existent à l’interface entre les plaques GDL et BPP. La déforma-
tion élevée conduit à une petite porosité et une petite perméabilité de la GDL (d’après les
équations 1.17 et 1.22), qui conduisent à une difficulté d’alimenter les gaz réactifs et à une
perte de performance de la PAC. Ainsi la distribution de la déformation de la GDL locale est
étudiée afin de connaître son interaction sur la performance des piles.

F IG . 3.20: Déformation de la plaque GDL, dont zone AB : zone compressée sous la dent de
la BPP ; zone BC : zone proche de la dent et sous le canal de la BPP ; zone CD : zone loin de
la dent et sous le canal de la BPP.

La figure 3.20 montre la distribution de la déformation du modèle de la PAC. Parmi les

trois composants (MEA, GDL et BPP), la déformation principale du modèle concerne la
plaque GDL, la déformation maximale de la plaque GDL est de 0.548. En effet les plaques
Études des contacts dans la PAC 83

MEA et la BPP sont plus rigides que la plaque GDL. La déformation maximale parmi des
plaques MEA et BPP est 9.45E-5, ce qui peut être négligé par rapport à la déformation de la
GDL.

Sur la figure 3.20, la distribution de la déformation correspond à la distribution de la

contrainte (figure 3.18a). Pour bien étudier la distribution de la déformation de la plaque
GDL, la plaque GDL est divisé par trois zones. La zone AB est la zone compresée sous la
dent de la BPP ; la zone BC est la zone proche de la dent et sous le canal de la BPP ; et la zone
CD est la zone loin de la dent et sous le canal de la BPP. La déformation de la plaque GDL
a principalement lieu sur les zones AB et BC, la zone CD n’est quasiment pas déformée. La
zone AB de la GDL située sous la dent de la plaque BPP est évidemment compressée. De
plus, la zone supérieure de BC et CD de la GDL est gonflée.

La zone compressée et la zone gonflée peuvent individuellement influencer la perfor-

mance de la PAC. En effet, la déformation importante de la zone AB de la plaque GDL peut
entraîner une perte de porosité. Quand au “gonflement” de la partie supérieure de la zone
BC, il va provoquer une diminution de la surface du canal de la plaque BPP entraînant, ainsi
une accélération du flux des gaz et donc une diminution du temps réactifs des gaz.

Dans la zone compressée AB de la GDL, la déformation maximale de la plaque GDL de

0.548 apparaît aux points inférieurs des dents de la BPP ; vers l’axe de symétrie, la distribu-
tion de la déformation est plus uniforme. Sous la précompression de la PAC, le changement
de la déformation est donc à surveiller au niveau des congés des dents de la BPP.

La GDL est en papier carbone, son rôle est de fournir les gaz réactifs et d’évacuer l’eau. La
couche GDL a une porosité initiale élevée qui lui permet d’acheminer les gaz efficacement.
La porosité initiale de la GDL de TORAY est par exemple de 78%. Pendant l’assemblage de
l’empilement de la PAC, la déformation de la plaque GDL augmente, alors que la porosité
ou la perméabilité de la GDL diminue. La porosité minimale de la GDL apparaît au niveau
des congés où la déformation de la GDL est la plus élevée. Par contre autour la zone de
l’axe de symétrie, la déformation est plus uniforme et inférieure à la déformation maximale.
Sur la figure 3.20, la déformation maximale au fond du congé εv = εx + εy = 0.73 ; la
0 0 0
déformation sur l’axe de symétrie εv = εx + εy = 0.3. Selon l’équation 1.22, la porosité
de la GDL au fond du congé est presque 18%, la porosité de la GDL sur l’axe de symétrie
est de 69%. Par rapport à la porosité initiale de la GDL (78%), la porosité de la GDL a
donc effectivement diminué, la plus importante baisse se situant au niveau du congé. Il faut
noter que la diminution de la porosité est inévitable sous une précompression de la PAC qui
garantie la parfaite étanchéité. Cependant la précompression peut être contrôlée à un niveau
acceptable permettant d’obtenir une aire de contact élevée pour une résistance de contact
minimale une porosité correcte, tout en garantissant l’étanchéité de l’empilement des piles.
84 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

3.4.3.5 Déformation de la plaque GDL avec pores

Comme nous l’avons vu, la porosité de la plaque GDL est décrite par l’équation 1.22,
notamment la porosité aux zones locales. Bien que l’équation 1.17 puisse donner une per-
méabilité précise de la plaque GDL, elle ne permet pas facilement d’étudier le changement
de la porosité de la plaque GDL. Une description globale de la porosité est pourtant néces-
saire pour connaître le changement de la porosité. Un modèle numérique tenant compte de
petits pores qui existent dans la plaque GDL va permettre d’étudier ce phénomène. La taille
du rayon ri et le nombre de petits pores np sont importants et sont déterminés par la porosité
initiale de la plaque GDL ϕ0 et par le volume de la plaque GDL Vt . Le volume poreux Vp
(Vp = Vt × ϕ0 ) est défini par l’équation :
np
X
Vp = πri 2 (3.1)
i=1

Un modèle éléments finis pour étudier la porosité globale est présenté sur la figure 3.21a,
ce modèle est soumis aux mêmes conditions aux limites que le modèle précédent (figure
3.15b), la porosité initiale de la plaque GDL est de 78%.

(a) (b)

F IG . 3.21: (a) Plaque GDL avec des pores avec np = 24 et porosité de 78% de la GDL (b)
porosité globale de la GDL après déformation, dont zone AB : zone compressée sous la dent
de la BPP ; zone BC : zone proche de la dent et sous le canal de la BPP ; zone CD : zone loin
de la dent et sous le canal de la BPP.

La porosité globale de la GDL après déformation est montrée dans la figure 3.21b. La
porosité minimale de la GDL apparaît sous le congé de la dent de la BPP. Ce résultat corres-
Études des contacts dans la PAC 85

pond au résultat obtenu par l’équation 1.22, parce que le pore sous le congé de la dent de la
BPP a la déformation maximale comme nous avons étudié dans la figure 3.20. Pour les trois
autres pores autour de l’axe de symétrie dans la zone AB, la porosité est plus uniforme que
pour le pore sous le congé ce qui correspond à la conclusion ce que nous avons obtenu par
l’équation 1.22.
Par ailleurs, la figure 3.22 montre la porosité de toute la plaque GDL après déformation
et la porosité initiale. Sur la figure 3.22, nous pouvons voir que le changement de la porosité
des pores en zone AB est supérieur à la zone BC, qui lui-même est supérieur à la zone CD. La
zone CD n’a pas changé, sa porosité correspond à la porosité initiale. En plus, le changement
de la porosité des pores supérieurs de la GDL est plus important que les pores intermédiaires
qui est lui-même plus important que les pores inférieurs de la GDL. C’est donc un modèle
intéressant permettant de décrire la porosité totale de la GDL après déformation.

F IG . 3.22: Porosité de la GDL déformée et porosité initiale de la GDL

La figure 3.23 montre le même exemple mais avec une porosité initiale de 50% de la
plaque GDL.
Sur la figure 3.23, nous pouvons voir que le changement de la porosité dans la plaque
GDL est identique à la description de la figure 3.22, mais à cause d’une porosité initiale plus
faible, le pore au fond du congé est plus compressée que sur la figure 3.22.
D’après le résultat numérique de la distribution de la porosité totale de la GDL, on peut
noter que la taille des pores aux zones importantes où la déformation est maximale doit être
optimisée pendant la fabrication afin de garantir une porosité suffisante.
En réalité, la distribution des pores dans la plaque GDL n’est sûrement pas uniforme
comme dans le modèle très simple précédemment décrit (figure 3.21a). Un modèle éléments
finis a donc été développé avec une distribution aléatoire des pores dans la plaque GDL
(figure 3.25a). Les espaces blancs sont les pores dans la plaque GDL. Ce modèle est toujours
soumis aux mêmes conditions aux limites que le modèle de la figure 3.15b.
La porosité globale de la GDL avec des pores distribués aléatoirement est indiquée dans
la figure 3.25b. La distribution des pores dans la GDL étant aléatoire, la taille et le nombre
des pores sont incertains, mais l’espace des pores est soumise à la condition d’une porosité
86 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

F IG . 3.23: (a) Plaque GDL avec des pores avec np = 24 et porosité de 50% de la GDL (b)
porosité globale de la GDL après déformation, dont zone AB : zone compressée sous la dent
de la BPP ; zone BC : zone proche de la dent et sous le canal de la BPP ; zone CD : zone loin
de la dent et sous le canal de la BPP.

F IG . 3.24: Porosité de la GDL déformée et porosité initiale de la GDL

Études des contacts dans la PAC 87

(a) (b)

F IG . 3.25: (a) Modèle de distribution aléatoire des pores et (b) porosité globale de la GDL
après déformation, zone AB : zone compresée sous la dent de la BPP ; zone BC : zone proche
de la dent et sous le canal de la BPP ; zone CD : zone loin de la dent et sous le canal de la
BPP.

F IG . 3.26: Porosité de la GDL déformée et porosité initiale de la GDL

88 Chapitre 3 : Étude des effets mécaniques dans une cellule de la pile

initiale de 78%, la déformation d’un petit pore dans la GDL n’est pas aussi évidente que sur
la figure 3.21b. Cependant nous pouvons encore contrôler que la plus grande déformation
apparaît au niveau des congés des dents de la BPP, mais on ne peut pas dire que la porosité
minimale est encore localisée ici parce que la distribution des pores étant aléatoire, la taille
des trous n’est pas uniforme dans la plaque GDL. La porosité dans la partie inférieure du
congé est encore grande ce qui est mieux pour alimenter des gaz réactifs. La diminution de
la porosité dans la zone AB est supérieure à la diminution de la porosité dans la zone BC qui
elle même est supérieure à la diminution de la porosité dans la zone CD. Cette situation est
montrée dans la figure 3.26.
Par rapport à la structure uniforme des pores dans la plaque GDL (figure 3.22), une dis-
tribution aléatoire des pores paraît meilleure pour conserver une porosité suffisante.

3.5 Conclusions
Dans ce chapitre, nous avons présenté un modèle éléments finis d’une cellule de la PAC
(maillage, conditions aux limites, définition des zones de contact). Nous avons vu qu’un
décollement apparaissait entre la plaque GDL et la plaque MEA. Un recuit permet de s’af-
franchir de ce problème. En considérant ces deux éléments collés, nous avons ensuite vu que
l’angle droit de la dent de la plaque BPP engendrait une concentration de contraintes. L’uti-
lisation d’un congé permet de réduire ce problème. Les pressions de contact et les aires de
contact ont ici été étudiés dans les deux cas : avec l’angle droit et avec le congé. Différents
maillages ont également été utilisés pour voir leur impact sur les résultats. Enfin, la défor-
mation de la plaque GDL est étudiée pour décrire la porosité de la plaque GDL en utilisant
un modèle éléments finis avec des pores.
 89

Chapitre 4

Influence des principaux paramètres sur

la performance des piles

4.1 Introduction
D’après l’étude menée dans le chapitre précédent, nous avons vu que certains effets mé-
caniques ont des influences importantes sur la performance de la PAC et sur la durée de vie.
Notamment la contrainte, la déformation, la pression de contact et l’aire de contact peuvent
être responsables d’une diminution de la performance de la PAC. Ces effets mécaniques dé-
pendent des paramètres de la PAC comme la précompression, les propriétés des composants,
et les paramètres géométriques de l’interface. Une étude de ces paramètres s’avère nécessaire
pour obtenir la meilleure performance possible de la PAC et optimiser sa structure.
Dans cette partie, une étude sera effectuée sur les influences des paramètres modifiant
les effets mécaniques aux interfaces de la structure de la PAC. Les paramètres tels que la
précompression, la taille du congé des dents de la plaque BPP, la largeur des dents de la
plaque BPP, l’épaisseur de la plaque GDL et de la membrane MEA, etc. vont être étudiés.
Une étude sur l’influence de l’aspect de surface de la GDL et de la BPP sera également
présentée. Ces études mécaniques pourront aider à déterminer des paramètres optimaux pour
une meilleure performance et une conception optimale de la PAC.

4.2 Précompression de la PAC

L’empilement de la PAC est serré par des boulons tirants pour le fonctionnement normal,
cette force produit une précompression dans tous les composants de la PAC. L’état de contact
aux interfaces entre les plaques GDL et BPP est principalement influencé par la précompres-
sion. En effet, une précompression insuffisante conduit à un contact imparfait à l’interface,
car il existe alors une zone où le contact ne s’établit pas encore, la PAC n’atteint pas la résis-
90 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

tance de contact minimale. À l’opposé, une précompression trop élevée diminue la porosité
de la plaque GDL, ce qui peut rétrécir les pores dans la plaque GDL et réduire l’alimenta-
tion des gaz réactifs. L’objectif est donc de trouver une précompression optimale offrant un
compromis entre ces deux cas extrêmes.
Dans le modèle de la PAC, la configuration géométrique du modèle est inchangée, seule
la précompression appliquée dans le modèle varie. Nous étudions ici l’influence de la pré-
compression sur la distribution de la pression de contact, l’aire de contact, la déformation et
la porosité de la plaque GDL.

4.2.1 Pression de contact

Dans un premier temps, la variation de la pression de contact (en Mpa) entre la plaque
GDL et la plaque BPP pour différentes précompressions (en mm) est indiquée dans les fi-
gures ci-dessous.

F IG . 4.1: Précompression = 0.01 mm F IG . 4.2: Précompression = 0.02 mm

F IG . 4.3: Précompression = 0.04 mm F IG . 4.4: Précompression = 0.06 mm

Sur les figures (figure 4.1 à figure 4.8), la distribution de la pression de contact est simi-
laire avec différentes précompressions : la pression est relativement uniforme à l’interface
sauf dans la zone du congé de la dent de la plaque BPP où la variation de la pression est
importante. Lorsqu’on observe l’évolution de la répartition de la pression de contact lorsque
la précompression augmente, la pression maximale localisée au niveau du congé comme la
pression au niveau de l’axe de symétrie augmentent également. Toutefois la répartition de
pression au niveau du congé se modifie au fur et à mesure du l’augmentation du chargement.
Précompression de la PAC 91

F IG . 4.5: Précompression = 0.08 mm F IG . 4.6: Précompression = 0.1 mm

F IG . 4.7: Précompression = 0.12 mm F IG . 4.8: Précompression = 0.14 mm

La pression sur l’axe de symétrie est considérée comme la pression moyenne uniforme
à l’interface et la pression maximale est donnée par la pression au congé de la dent de la
plaque BPP. La variation de ces deux pressions (courbe 1 et courbe 2) en fonction de la
précompression est indiqué dans la figure 4.9. La pression de contact augmente évidemment
avec la précompression, non seulement la pression à l’axe de symétrie (courbe 1), mais aussi
la pression maximale au congé (courbe 2).
En force équivalente, la somme des produits de la pression et de l’aire de contact à l’in-
terface est la précompression extérieure appliquée dans la PAC, ce qui peut être représentée
par :

n
X
F = p i × Ai , (4.1)
i=1

où F est la précompression de la PAC, pi est la pression de iime élément de contact, Ai

est l’aire de contact de iime élément de contact, n est la somme des éléments de contact à
l’interface. D’après l’équation 4.1, en cas de contact imparfait, l’augmentation de la pré-
compression conduit à l’augmentation de la pression de contact et de l’aire de contact à
l’interface. Par contre en cas de contact parfait, l’aire de contact est constante, la pression de
contact peut toujours augmenter avec la précompression de la PAC comme les courbes 1 et
2 indiqués dans la figure 4.9.
Pmax −Paxe
Le pourcentage Puni = Pmax
× 100% peut être utilisée pour représenter l’uniformité
92 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

F IG . 4.9: Variation de la pression avec la précompression de la PAC

de la pression de contact, c’est-à-dire, le pourcentage de la différence de pression de contact

entre la pression maximale et la pression à l’axe de symétrie par rapport à la pression maxi-
male Pmax . Sur la figure 4.9, le pourcentage Puni (courbe 3) diminue avec la précompression
de la PAC, c’est-à-dire, la variation de la pression à l’interface est de plus en plus faible
avec l’augmentation de la précompression de la PAC. On peut donc conclure que dans une
certaine mesure, l’augmentation de la précompression permet d’améliorer l’uniformité de la
pression normale de contact. Une distribution de la pression uniforme à l’interface facilite
la production de courant électrique et permet d’obtenir une température uniforme à l’inter-
face, l’uniformité de la pression de contact joue donc un rôle important dans la performance
de la PAC comme mentionné dans l’expérience effectuée par Wang et al. [WSZ08]. Ainsi
une précompression de la PAC assez élevée est nécessaire pour réaliser un contact parfait à
l’interface de la PAC et de ce fait, améliorer la performance de la PAC.

4.2.2 Aire de contact (longueur de contact)

L’équation 4.1 montre que si la précompression F croît ce qui est le cas pendant l’assem-
blage de la PAC, non seulement la pression de contact augmente avec la précompression,
mais aussi l’aire de contact entre la plaque GDL et la dent de la plaque BPP grâce au chan-
gement de la zone de contact au niveau des congés de la plaque BPP. L’aire de contact est
un paramètre principal dans la résistance de contact. Un contact parfait avec une large aire
de contact permet d’avoir une faible résistance de contact à l’interface, obtenant ainsi une
Précompression de la PAC 93

meilleure performance. Ainsi une étude de la variation de la longueur de la zone de contact

avec la précompression de la PAC est représentée dans la figure 4.10 pour trouver un état de
contact parfait à l’interface avec la plus grande aire de contact possible.

F IG . 4.10: Variation de l’aire de contact avec la précompression de la PAC

La figure 4.10 montre que la zone de contact s’élargit avec la précompression de la PAC
jusqu’à 0.1 mm. En revanche, pour une précompression supérieure à 0.1 mm, la zone de
contact ne varie plus, sa longueur se stabilise à 0.42852 mm. On constate donc que la pré-
compression de 0.1 mm est la précompression minimale pour réaliser un contact parfait à
l’interface avec une longueur de contact de 0.42852 mm. Toute la zone de contact potentielle
à l’interface entre les plaques GDL et BPP est alors en contact. Bien sur l’augmentation de
la précompression de la PAC mène également à l’augmentation de la pression à l’interface
comme indiqué dans la figure 4.9.
Comme nous l’avons étudié dans le chapitre 3 (figure 3.6), la zone de contact à l’interface
entre la plaque GDL et la plaque BPP est constituée de deux parties : une partie est le contact
de surface plane à surface plane entre la plaque GDL et la surface basse de la dent de la plaque
BPP, et l’autre partie est le contact surface plane à surface courbée entre la plaque GDL et
le congé de la dent de la plaque BPP. Dans les figures (figure 4.1 à figure 4.8), les interfaces
rentrent en contact initialement au niveau des surfaces planes puis entre le congé et la plaque
plane de la GDL avec l’augmentation de la précompression. Quand la précompression est
assez élevée, tout le côté du congé de la dent est en contact avec la plaque GDL. Dans le
modèle, la précompression de 0.1 mm est la force minimale qui permet de réaliser "le contact
parfait" entre le congé de la BPP et la plaque GDL. En revanche, quand la précompression
dépasse 0.1 mm, il n’y a plus de zones de contact potentielles, c’est-à-dire, l’aire de contact
n’augmente plus à l’interface. Ainsi pour une précompression de 0.1 mm, un état de contact
94 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

idéal permet d’obtenir une résistance de contact minimale à l’interface de la PAC.

Quand la précompression est suffisamment élevée, la résistance de contact est constante,
ce qui correspond aux expériences effectuées par Lee et al.[LHVZM99], Ge et al. [GHL06]
et Chang et al. [CHWC07], etc.

4.2.3 Porosité de la plaque GDL

Nous venons de voir qu’une certaine valeur de la précompression permet d’établir une
zone de contact de taille maximale et donc une résistance de contact minimale ce qui est bon
pour la performance de la pile. Cependant la précompression provoque aussi la déformation
des trois plaques BPP, MEA et GDL au niveau du contact. Comme un matériau en tissu,
mousse déformante, la plaque GDL a une grande déformation par rapport aux plaques BPP
et MEA à cause de sa faible rigidité. La déformation est importante dans la zone compressée,
sous la dent de la plaque BPP et dans la zone gonflée sous la canal de la plaque BPP. La
déformation est liée à la porosité et à la perméabilité de la GDL et peut agir indirectement
sur l’alimentation des gaz réactifs dans la PAC. Une grande déformation pouvant réduire la
porosité de la plaque GDL, une étude de la variation de la déformation de la GDL en fonction
de la précompression de la PAC a été réalisée dans cette partie.
Dans un premier temps, la variation de la déformation du modèle avec différentes pré-
compressions est représentée dans les figures ci-dessous :

F IG . 4.11: Précompression = 0.01 mm F IG . 4.12: Précompression = 0.02 mm

F IG . 4.13: Précompression = 0.04 mm F IG . 4.14: Précompression = 0.06 mm

Précompression de la PAC 95

F IG . 4.15: Précompression = 0.08 mm F IG . 4.16: Précompression = 0.1 mm

F IG . 4.17: Précompression = 0.12 mm F IG . 4.18: Précompression = 0.14 mm

D’après les figures 4.11 à 4.18, la déformation de la plaque GDL augmente avec la pré-
compression de la PAC, non seulement la déformation au congé (point 1, figure 3.6) mais
aussi la déformation sur l’axe de symétrie (point 3, figure 3.6). L’augmentation de la pré-
compression accroît l’aire de contact à l’interface, accroît la zone compressée de la plaque
GDL et diminue la zone non déformée jusqu’à la précompression de 0.1 mm.
La figure 4.19 montre la variation de la déformation au point 1 et au point 3 (figure 3.6)
et la porosité de la plaque GDL avec la précompression de la PAC. La déformation de la
plaque GDL est quasi uniforme sauf au niveau du congé de la BPP. Elle peut donc être
représentée par la déformation sur l’axe de symétrie. La déformation maximale apparaît au
point d’intersection entre le congé et la surface plane représente (point 1, figure 3.6). Ainsi
la déformation maximale (courbe 1), la déformation sur l’axe de symétrie (courbe 2), les
porosités correspondantes (courbe 3 et courbe 4) sont présentées dans la figure 4.19.
Sur la figure 4.19, la variation de la déformation maximale de la plaque GDL au niveau du
congé (courbe 1) et de la déformation uniforme sur l’axe de symétrie (courbe 2) augmentent
avec la précompression. La porosité étant contraire à la déformation (équation 1.22), la po-
rosité de la plaque GDL au congé (courbe 3) et la porosité sur l’axe de symétrie (courbe 4)
diminuent donc avec la précompression de la PAC comme indiqué dans la figure 4.19. Une
porosité importante ne sera garantie que par une faible précompression. Or il est nécessaire
de conserver une porosité assez élevée pour la performance de la pile.
En effet, une faible porosité conduit à une circulation insuffisante des gaz réactifs : l’eau
produite par la réaction électrochimique ne s’évacue plus, ce qui conduit à un phénomène
96 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

d’inondation (gouttelette de l’eau) dans la plaque GDL et bloque la circulation des gaz dans
la plaque GDL. Ainsi, il est important dans la conception de la plaque GDL d’éviter des
zones de porosité nulle et il faut garder une grande porosité relative dans la plaque GDL.
Sur la figure 4.19 pour une précompression de 0.1 mm, la porosité minimale de la GDL au
point 1 (figure 3.6) est environ de 18%. Quand la précompression est de 0.12 mm, on a une
déformation de 0.84 et la porosité au point 1 est nulle. En ce qui concerne la porosité de la
GDL sur l’axe de symétrie (point 3, figure 3.6), la porosité avec une précompression de 0.1
mm est de 69%, qui n’abaisse que de 9 % la porosité initiale de valeur 78%. La diminution de
la porosité sur l’axe de symétrie est plus douce que celle au point 1. Toutefois si la précom-
pression dépasse 0.1 mm, cette porosité diminuera plus vite. Parce que la précompression de
0.1 mm peut réaliser un contact parfait à l’interface avec une aire de contact maximale dans
la PAC, la précompression de 0.1 mm est favorable non seulement à l’aire de contact mais
aussi à la porosité de la plaque GDL.

F IG . 4.19: Variation de la déformation et de la porosité avec la précompression

D’après cette étude, nous pouvons voir que la précompression de la PAC influence non
seulement l’aire de contact qui est liée à la résistance de contact de la PAC, mais aussi la dé-
formation de la GDL qui est liée à la porosité de la plaque GDL. L’aire de contact augmente
avec la précompression la PAC, une précompression élevée peut réaliser un contact parfait à
l’interface, ce qui conduit à une petite résistance de contact, favorable à la performance de
la PAC. En revanche, la déformation de la plaque GDL augmente avec la précompression et
conduit à une baisse de la porosité de la plaque GDL, ce qui provoque une alimentation insuf-
fisante des gaz réactifs et abaisse la performance de la PAC. Ainsi un compromis concernant
Rayon du congé des dents de la plaque BPP 97

la valeur de la précompression devra être fait pour permettre d’établir la plus grande zone de
contact possible tout en conservant une porosité maximale.
Notre étude numérique de la PAC est en accord avec des expériences effectuées pour me-
surer la distribution de la performance de la PAC pour différentes précompressions [LLA98]
[WSZ08] [GHL06] : avant la précompression optimale, l’augmentation de la précompres-
sion conduit à l’augmentation de l’aire de contact à l’interface et à une performance élevée
de la PAC ; si la précompression dépasse cette valeur optimale, l’aire de contact n’augmente
plus, mais la diminution de la porosité continue ce qui abaisse la performance de la PAC.
Ces études indiquent que pour une meilleure performance de la PAC, une précompres-
sion optimale existe et garantit la condition de l’étanchéité de la PAC. Cette précompression
produit une aire de contact maximale et une résistance de contact minimale à l’interface,
mais la précompression doit permettre de conserver une porosité élevée dans la plaque GDL.
D’après nos résultats, une précompression de 0.1 mm paraît être un choix optimal pour l’aire
de contact dans ce cas là .

4.3 Rayon du congé des dents de la plaque BPP

La variation de la pression de contact et de l’aire de contact avec la précompression
montrent que la précompression de la PAC peut influencer l’état de contact à l’interface,
le changement de l’état de contact intervient principalement au niveau des congés des dents
de la plaque BPP. La présence du congé améliore l’effet de la singularité de l’angle droit
de la dent rectangulaire de la plaque BPP. Une variation du rayon du congé sur les effets
mécaniques à l’interface dans la PAC va donc être étudiée. Dans cette partie, on regarde plus
précisément l’influence du rayon du congé sur la pression, l’aire de contact et sur la porosité
de la GDL.

4.3.1 Pression de contact

Nous avons vu dans le chapitre 3 que le remplacement de l’angle droit de la dent de
la plaque BPP qui apporte une singularité au niveau du contact par un congé permettait
d’obtenir une distribution de pression de contact plus régulière avec une valeur maximale
moins élevée. Dans un premier temps, la variation de la pression de contact (en Mpa) entre
la plaque GDL et la plaque BPP sous différentes rayons du congé (en mm) est indiquée sur
les figures de 4.20 à 4.31 :
Sur les figures 4.20 à 4.31, l’importance du rayon du congé sur la distribution de la pres-
sion à l’interface est évidente. On observe deux choses sur ces figures : la pression maximale
de contact diminue à mesure que le rayon augmente et les pressions maximales se déplacent
petit à petit vers l’axe de symétrie. La figure 4.32 montre la variation de la pression de contact
98 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

F IG . 4.20: Rayon RRIB = 0.01 mm F IG . 4.21: RRIB = 0.03 mm

F IG . 4.22: RRIB = 0.05 mm F IG . 4.23: RRIB = 0.07 mm

F IG . 4.24: RRIB = 0.09 mm F IG . 4.25: RRIB = 0.12 mm

F IG . 4.26: RRIB = 0.16 mm F IG . 4.27: RRIB = 0.20 mm

Rayon du congé des dents de la plaque BPP 99

F IG . 4.28: RRIB = 0.25 mm F IG . 4.29: RRIB = 0.30 mm

F IG . 4.30: RRIB = 0.35 mm F IG . 4.31: RRIB = 0.40 mm

au niveau du congé (point 1, figure 3.6) et la pression sur l’axe de symétrie (point 3, figure
3.6) avec le rayon. Nous pouvons voir que l’augmentation du rayon du congé conduit à la di-
minution de la pression de contact, notamment la pression au point 1. La pression maximale
au point 1 est très sensible au rayon du congé. La pression chute rapidement lorsqu’on passe
d’un rayon de 0.01 mm à 0.03 mm.
Avec un rayon RRIB de 0.01 mm, la forme de la dent est encore quasi-rectangulaire, on
constate donc l’effet singulier de l’angle droit sur la pression. La précompression de la PAC
produit une pression de contact élevée de 11.7 Mpa au point 1. Mais avec un rayon RRIB de
0.05 mm, la pression de contact au point 1 diminue jusqu’à 3.44 Mpa, ce qui diminue de 72
% la pression obtenue avec un rayon de 0.01 mm, l’effet de l’angle droit à l’interface a été
clairement amélioré. Enfin quand RRIB est supérieur à 0.05 mm, la diminution de la pression
de contact au point 1 est de plus en plus douce. Cependant, la distribution de la pression sur
l’axe de symétrie (figure 4.32) augmente très légèrement avec l’augmentation du rayon du
congé parce que le contact est plan à plan autour de la zone proche de l’axe de symétrie, le
rayon du congé influence donc principalement la pression au niveau du congé à l’interface.
De plus, dans la figure 4.32, quand le rayon est supérieur à 0.3 mm, la pression sur l’axe
de symétrie augmente légèrement avec le rayon du congé jusqu’au rayon de 0.4 mm. Parce
que la largeur de la dent de la plaque BPP est de 0.8 mm dans le modèle de la PAC, le largeur
de la demi dent est de 0.4 mm dans le modèle symétrique. Quand le rayon est supérieur à 0.3
mm, la fin du congé arrive presque au milieu de la dent de la BPP, les pressions maximales
se rapprochent de l’axe de symétrie comme indiqué dans les figures 4.29, 4.30 et 4.31.
100 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

F IG . 4.32: Variation de la pression de contact avec le rayon du congé

Quand le rayon RRIB est 0.4 mm, la dent de la plaque BPP est tout à fait semi-circulaire,
la pression maximale est au point le plus bas de la dent, et la pression maximale est sur l’axe
de symétrie (figure 4.31). C’est pour cette raison que la courbe de la pression maximale et
la courbe de la pression sur l’axe de symétrie, dans la figure 4.32, se rejoignent pour un
rayon de 0.4 mm. À ce moment, la pression maximale à l’interface est 2.52 Mpa, ce qui est
inférieure de la pression maximale (3.44 Mpa) avec un rayon de 0.05 mm.

4.3.2 Aire de contact (longueur de contact)

Dans un premier temps, l’influence du rayon du congé sur la pression a été étudiée. La
présence du congé permet d’améliorer la concentration de la pression de contact (figures 4.20
à 4.31) et peut aussi influencer la longueur de contact à l’interface. La figure 4.33 montre la
variation de la longueur de la zone de contact avec le rayon du congé. La courbe n’est pas
linéaire car elle augmente au début et puis diminue jusqu’à la fin avec l’augmentation du
rayon du congé.
L’augmentation du rayon du congé dans la figure 4.33 va de 0 mm (dent rectangulaire) à
0.4 mm. Avec l’augmentation de la taille du rayon, la demi-longueur de la zone de contact
augmente jusqu’à la valeur maximale de 0.428 mm pour RRIB = 0.05 mm. D’après la relation
entre l’aire de contact et la résistance de contact dans la PAC, quand RRIB = 0.05 mm, la
PAC a une résistance minimale et une performance élevée.
Quand la taille du rayon du congé est supérieure à 0.05 mm, l’aire de contact entre la
plaque GDL et la plaque BPP commence à diminuer. Plus large du rayon du congé, moins de
Rayon du congé des dents de la plaque BPP 101

l’aire de contact à l’interface. L’aire de contact minimale est de 0.42 mm avec un rayon du
congé de 0.4 mm quand la dent de la plaque BPP est tout à fait semi-circulaire à ce moment.
L’aire de contact avec RRIB = 0.4 mm est réduite de 50% par rapport à l’aire maximale
obtenue pour RRIB = 0.05 mm. Même si la pression obtenue avec un rayon de 0.4 mm est
inférieure à la pression pour un rayon de 0.05 mm (figure 4.32), un rayon de 0.05 mm est une
meilleure alternative pour la performance de la PAC. Autrement dit, une dent rectangulaire
avec un petit congé donnera une meilleure performance que une dent circulaire de la plaque
BPP.

F IG . 4.33: Variation de l’aire de contact avec le rayon du congé

Enfin, il est bon de noter que dans la figure 4.33, quand le rayon du congé est inférieure
à 0.09 mm, l’aire de contact est plus grande avec une dent qu’avec l’angle droit de la plaque
BPP. La conception de la dent de la PAC avec un congé est donc meilleure que la conception
de la dent avec un coin droit, non seulement pour l’amélioration de la concentration de
pression de contact, mais aussi pour une aire de contact plus étendue.
Dans l’étude de la variation du rayon du congé sur l’étendue du contact, une valeur maxi-
male de l’aire de contact est donnée pour un rayon (RRIB = 0.05 mm) qui optimise la per-
formance des piles.

4.3.3 Porosité de la plaque GDL

Tout comme l’influence du rayon RRIB sur l’aire de contact et la pression de contact,
nous allons étudier l’influence du rayon du congé sur la porosité de la plaque GDL. En effet,
le rayon RRIB modifie la distribution de la déformation de la plaque GDL. La figure 4.34
montre la variation de la déformation au point 1 (figure 3.6) et la déformation sur l’axe de
102 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

symétrie (point 3, figure 3.6) avec le rayon RRIB et ses porosités correspondantes. Sur la
figure 4.34, on constate que la déformation au point 1 (courbe 1), diminue avec le rayon
du congé alors que la porosité associée augmente rapidement (courbe 3). En revanche, la
déformation sur l’axe de symétrie (courbe 2) et sa porosité (courbe 4) varient très lentement
mais toujours inversement.
La diminution de la déformation au point 1 s’explique comme dans le chapitre précédent
par la diminution de l’effet de l’angle droit. Le congé améliore la circulation des gaz réactifs
dans le cas d’une dent rectangulaire de la plaque BPP. Avec l’augmentation du rayon du
congé, la concentration de contrainte a été améliorée et la déformation de la plaque GDL
diminuée. Quand la dent de la BPP est juste circulaire avec un rayon de 0.4 mm, le point
1 coïncide alors avec l’axe de symétrie. Ainsi, dans la figure 4.34, la déformation au pint 1
et la déformation sur l’axe de symétrie sont identiques pour un rayon de 0.4 mm. Mais la
porosité de la GDL varie inversement à la déformation, l’augmentation du rayon du congé
conduit donc à l’accroissement de la porosité de la plaque GDL au point 1. La porosité sur
l’axe de symétrie varie très doucement car le contact est régulier (plan à plan ) autour l’axe
de symétrie.

F IG . 4.34: Variation de la déformation et la porosité avec le rayon du congé

Sur la figure 4.34, pour RRIB = 0.01 mm, la porosité de la plaque GDL sur l’axe de
symétrie est de 68.8 % ; au même endroit pour RRIB = 0.4 mm, la porosité est de 63.2 %, la
variation de la porosité sur l’axe de symétrie est d’environ 6 % ce qui est faible. En revanche
au point 1, la déformation diminue rapidement avec le rayon du congé et la porosité de la
GDL augmente logiquement . La déformation au point 1 pour 0.01 mm est assez élevée et
sa porosité est nulle. Avec l’augmentation du rayon du congé, la porosité atteint 12 % pour
un rayon de 0.05 mm et 42 % pour 0.09 mm. Pour un rayon supérieur à 0.15 mm, la courbe
Largeur des dents de la plaque BPP 103

de la porosité au point 1 (courbe 3) est presque proportionnelle au rayon du congé jusqu’à

ce qu’elle atteigne 63 % pour un rayon de 0.4 mm.
Considérant l’influence du rayon du congé sur l’aire de contact, la pression à l’interface
et la porosité de la plaque GDL, un rayon de 0.05 mm semble être le meilleur compromis
pour la performance de la PAC même s’il désavantage la porosité de la pile. Cependant on
pourra considérer un rayon pour le congé compris entre 0.05 mm et 0.09 mm si une baisse
du coût de fabrication le nécessite.

4.4 Largeur des dents de la plaque BPP

Comme indiquée dans la figure 1.14, la plaque BPP de la PAC comprend des canaux et des
dents afin de réaliser le fonctionnement de la conductivité des électrons et de la distribution
des gaz réactifs dans la PAC. La largeur des dents dans la PAC détermine largement l’aire de
contact initiale entre les plaques GDL et BPP, ce qui représente un chemin pour le passage
des électrons produits. La largeur des canaux influence la circulation des gaz réactifs comme
l’hydrogène et l’oxygène, c’est-à-dire la capacité de flux des gaz réactifs pour la réaction
électrochimique. Ces deux largeurs entrent en jeux pour la performance des piles. Un petit
canal conduit à une alimentation insuffisante des gaz réactifs et une petite dent conduit à
une résistance de contact élevée à l’interface. La dimension de la plaque BPP étant fixe,
il faut donc essayer de trouver un bon rapport entre les dimensions du canal et de la dent.
Pour l’instant, les différentes études présentées ont été réalisées en considérant une taille
identique pour la dent et le canal. Mais l’augmentation de la largeur des dents permet peut
être d’augmenter l’aire de contact à l’interface. Dans ce modèle, la largeur des dents WRIB
et celle des canaux WCN varient l’une et l’autre en conservant WRIB + WCN = 1.6 mm.

F IG . 4.35: Variation de l’aire de contact avec la largeur des dents de la BPP

104 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

La figure 4.35 montre la variation de l’aire de contact avec différentes largeurs des dents
de la BPP (WRIB =0.4 mm, 0.8 mm et 1.2 mm) sous différentes précompressions de la PAC.
L’aire de contact augmente comme prévu avec la largeur des dents. Plus les dents de la plaque
BPP sont larges, plus l’aire de contact est grande, moins la résistance de contact est élevée.
En revanche, un canal large est plus favorable à la circulation des gaz réactifs. Or un canal
plus large implique une dent plus étroite car WRIB + WCN doît rester constant. Il faudra là
encore trouver un compromis.
Les figures 4.36, 4.37, 4.38 représentent la distribution de la différence entre la pression
au point 1 et la pression sur l’axe de symétrie (point 3), de la porosité de la GDL au point 1
et de la porosité de la GDL sur l’axe de symétrie (point 3). Sur la figure 4.36, nous pouvons
voir qu’avec l’augmentation de la précompression, l’uniformité de la pression à l’interface
pour différentes largeurs de dents est très proche. Sur les figures 4.37, 4.38, on constate que
la largeur des dents de la BPP n’influence ni la porosité au point 1, ni la porosité sur l’axe de
symétrie. La largeur des dents de la plaque BPP n’influence donc que principalement l’aire
de contact.

F IG . 4.36: Différence de la pression de contact au poin 1 et 3 avec différentes largeurs

Pour une plaque BPP de longueur LBP P , de largeur WBP P , si la largeur des dents est
WRIB , la largeur des canaux est WCN , la longueur des dents est LBP P , la longueur de contact
0
avec des congés est WRIB , l’aire de contact totale aux interfaces multi-contacts AT entre les
petits dents de la plaque BPP et la plaque GDL est le somme des aires de contact individuelles
Largeur des dents de la plaque BPP 105

F IG . 4.37: Porosité au point 1 pour différentes largeurs

F IG . 4.38: Porosité sur l’axé de symétrie avec différentes largeurs des dents
106 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

d’une dent avec la plaque GDL.

N
X 0
AT = WRIB (i) × LBP P (4.2)
i=1

N est le nombre de dent de la plaque BPP, i est la iime zone de contact.

WBP P
N= (4.3)
WRIB + WCN

Parce que le déplacement appliqué s’applique de la même façon sur toute la largeur de la
BPP, toutes les zones de contact ont la même aire de contact. On peut donc écrire :
0
0 W × LBP P × WBP P
AT = N × WRIB × LBP P = RIB (4.4)
WRIB + WCN

si le rapport entre la largeur des canaux et la largeur des dents est RWCN /WRIB , c’est-à-
dire, le WCN = RWCN /WRIB × WRIB , on a
0
W LBP P × WBP P
AT = RIB ( ) (4.5)
WRIB 1 + RWCN /WRIB

0
puisque WRIB = WRIB + (π − 2) × RRIB ,

(π − 2) × RRIB LBP P × WBP P

AT = (1 + ) (4.6)
WRIB 1 + RWCN /WRIB

Parce que la dimension de la plaque BPP est constante et carrée, LBP P et WBP P dans
l’équation 4.6 sont constants, l’aire de contact totale est déterminée par les paramètres RRIB ,
WRIB et le rapport RWCN /WRIB . Le RRIB est un paramètre positive à l’aire de contact totale
aux interfaces multi-contacts en cas d’un contact parfait, par contre le WRIB et le rapporteur
RWCN /WRIB sont les paramètres négatives à l’aire de contact totale aux interfaces multi-
contacts. Au cas du contact parfait, l’aire de contact augmente avec le rayon RRIB avant d’un
rayon de 0.05 mm comme être étudié dans la figure 4.33. Le RRIB de 0.05 mm est le rayon
qui réalise l’aire de contact maximale aux interfaces. Par contre, pour les deux paramètres
négatives à l’aire de contact, le WRIB et le rapporteur RWCN /WRIB , ses petits valeurs sont
favorables à l’aire de contact totale aux interfaces. C’est-à-dire, la largeur des dents plus
fin conduit à une plus large de l’aire de contact totale aux interfaces. Afin de trouver un bon
compromis entre une largeur de dent importante (mais un canal plus étroit) qui favorise l’aire
de contact et une largeur de canal assez grande qui favorise la circulation des gaz réactifs,
nous garderons une largeur 0.8 mm identique pour les deux largeurs.
Épaisseur de la plaque GDL 107

4.5 Épaisseur de la plaque GDL

Nous avons vu que le contact parfait était déterminé par la précompression de la PAC, la
pression de contact influencée par le rayon du congé, l’aire de contact déterminée en grande
partie par la largeur des dents de la plaque BPP. L’épaisseur de la plaque GDL est aussi un
paramètre important pour la performance de la PAC, qui intervient dans la porosité de la
plaque GDL. La plaque GDL est poreuse et permet le passage des gaz réactifs. Comme nous
avons étudié dans l’équation 1.22, la porosité et la perméabilité de la GDL sont déterminées
par la déformation volumique de la GDL. Le changement de l’épaisseur de la plaque GDL
peut influencer la distribution de la déformation volumique et la porosité. Il apparaît donc
intéressant de considérer l’influence de l’épaisseur de la plaque GDL sur la performance de
la pile.
Actuellement quatre types de plaque GDL sont disponibles dans la commerce dans la
série TGP-H, ce sont TGP-H-120, TGP-H-90, TGP-H-60 et TGP-H-30. Leurs épaisseurs
sont 0.375 mm, 0.285 mm, 0.195 mm et 0.11 mm. Dans un premier temps, la porosité de la
plaque GDL au point 1 et la porosité sur l’axe de symétrie pour différentes épaisseurs de la
plaque GDL sont représentées dans les figures 4.39 et 4.40.

F IG . 4.39: Variation de la porosité au point 1 pour différentes épaisseurs

Les figures 4.39 et 4.40 montrent que l’épaisseur de la plaque GDL a une influence im-
portante sur la porosité de la plaque GDL. Plus l’épaisseur de la plaque augmente, plus la
porosité est importante. Afin d’optimiser la performance de la PAC, l’épaisseur de la GDL
recommandée sera de 0.375 mm.
108 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

F IG . 4.40: Variation de la porosité sur l’axe de symétrie pour différentes épaisseurs

4.6 Épaisseur de la membrane MEA

Nous étudions maintenant l’influence de l’épaisseur de la membrane MEA sur la perfor-
mance de la PAC. Afin de diminuer la résistance propre de la membrane MEA et d’augmenter
la conductivité des protons, on considère une membrane MEA de plus en plus fine dans la
PAC. Actuellement quatre types de la membrane MEA commerciales sont disponibles dans
la série de Nafion. Leur épaisseurs sont 0.185 mm de Nafion 117, 0.125 mm de Nafion 115,
0.09 mm de Nafion 113 et 0.05 mm de Nafion 112. L’influence de l’épaisseur de la membrane
MEA sur les effets mécaniques sera effectuée dans cette partie.
Dans un premier temps, la différence des pressions entre la valeur maximale au point
1 et la valeur sur l’axe de symétrie est indiquée dans la figure 4.41 pour différentes épais-
seurs des plaques MEAs. Cette différence de pression à l’interface s’améliore lorsque la pré-
compression de la PAC augmente. Nous pouvons constater que l’épaisseur de la membrane
n’influence pas la différence de pression.
Les figures 4.42 et 4.43 montrent la variation de la porosité au point 1 et de la porosité
sur l’axe de symétrie pour différentes épaisseurs des membranes. La porosité de la plaque
GDL diminue avec la précompression de la PAC, mais nous constatons que l’épaisseur de la
membrane n’agit pas sur la porosité. En fait la précompression de la PAC influence la défor-
mation de la plaque GDL qui n’est pas liée à l’épaisseur des membranes MEAs. L’épaisseur
de la plaque GDL peut donc influencer la déformation et la porosité comme nous l’avons
étudié, mais pas l’épaisseur des membranes MEAs.
Épaisseur de la membrane MEA 109

F IG . 4.41: Uniformité de pression de contact pour différentes épaisseurs des MEAs

F IG . 4.42: Porosité au point 1 pour différentes épaisseurs des membranes

110 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

F IG . 4.43: Porosité sur l’axe de symétrie pour différentes épaisseurs des membranes

En revanche l’épaisseur de la membrane MEA peut influencer la distribution de la contrainte

au sein de la membrane MEA. Les figures 4.44 et 4.45 montrent la distribution de la contrainte
dans une membrane MEA fine (0.05 mm) et dans la membrane MEA épaisse (0.183 mm).
La distribution de la contrainte avec les deux épaisseurs est similaire, la contrainte maximale
de deux épaisseurs est tout au centre de la membrane MEA. Mais la contrainte maximale est
de 0.5 Mpa (figure 4.44) avec la membrane MEA fine, ce qui est 35 % plus élevée que la
contrainte de 0.37 Mpa (figure 4.45) avec la membrane MEA épaisse. La membrane la plus
fine, soumise à des plus fortes contraintes sera donc probablement moins résistante dans le
temps (fatigue, endommagement) ce qui aura un impact sur la durée de vie de la pile. Toute-
fois dans nos études, les contraintes des membranes MEAs sont encore relativement petites
par rapport à la contrainte limite élastique de la membrane MEA.
De plus, la résistance propre de la membrane est plus petite avec une membrane MEA
fine qu’avec une membrane MEA épaisse, la plaque MEA de 0.05 mm sera donc privilégiée
dans le modèle de contact de la PAC.

4.7 Étude de la géométrie à l’interface

Comme nous l’avons vu, la géométrie comme la largeur des dents et des canaux de la BPP,
le rayon du congé, l’épaisseur de la GDL, l’épaisseur de la MEA influencent la performance
Étude de la géométrie à l’interface 111

F IG . 4.44: Contrainte dans la MEA (0.05 mm)

F IG . 4.45: Contrainte dans la MEA (0.183 mm)

112 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

des piles à combustibles. La géométrie à l’interface joue un rôle direct sur l’état de contact
entre les plaques BPP et GDL. On considère maintenant la tolérance géométrique issue de
la procédure de fabrication industrielle. Dans cette partie, deux études géométriques sont
considérées : la première porte sur la géométrie de l’interface, c’est-à-dire, la prise en compte
d’une interface rugueuse, la seconde sur la tolérance des hauteurs des dents de la plaque BPP.
En effet, sous une précompression uniforme, différentes hauteurs de dents peuvent conduire
à différents états de contact entre la plaque BPP et la plaque GDL.

4.7.1 Influence des aspérités de l’interface

D’un point de vue microstructure, la plaque GDL est en papier carbone, la plaque BPP
est gravée en graphite. Á cause des imperfections de la fabrication industrielle, leurs sur-
faces sont inévitablement rugueuses au niveau microscopique. La surface de la plaque GDL
est rugueuse car son matériau est en tissu avec des pores. La série de TGP-H a une rugosité
unique de 8 % comme indiquée dans le tableau 1.3 du chapitre 1. La rugosité de la plaque
BPP est plus difficile à contrôler que celle de la plaque GDL parce que la plaque BPP, en
graphite, est rigide. Besmann et al.[BKH03] ont mesuré l’aire de contact entre une plaque
GDL et une plaque BPP rugueuse avec trois types de rugosité différente. Dans ces expé-
riences, ils ont étudié l’influence des différentes aspérités de la plaque BPP sur la résistance
à l’interface. Une relation existe entre l’aspérité et la résistance comme indiqué dans la fi-
gure 4.46[BKH03]. Sur la figure 4.46, la résistance augmente avec l’aspérité de la surface,
c’est-à-dire, plus de l’interface est rugueuse, plus la résistance est élevée. Mais Besmann et
al. n’ont pas donné d’explication pour cette expérience. Aussi nous avons trouvé intéressant
de prendre en compte une interface rugueuse dans notre modèle afin d’étudier l’influence des
aspérités sur l’aire de contact entre les plaques GDL et BPP.
En fait, il existe deux façons de produire une surface rugueuse des plaques GDL et BPP.
Une façon est le VEECO (Vacuum Electronic Equipment Company) [Vee07] et l’autre façon
est la distribution de Gauss pour simuler l’aspérité de l’interface rugueuse. VEECO prend un
profilomètre opticien pour capturer les données des points de contact à la surface rugueuse.
Puis, les données de la surface rugueuse sont transférées dans le code des éléments finis pour
simuler la surface.
L’autre façon est de créer des aspérités avec une distribution de GAUSS pour simuler
la procédure de la fabrication industrielle. La taille de l’aspérité est supposée vérifier une
distribution Gaussienne normale avec une moyenne nulle et une déviation standard.
Dans le modèle de contact de la PAC, nous ferons l’hypothèse que la distribution des
aspérités à l’interface rugueuse des plaques GDL et BPP est soumise à une distribution
Gaussienne normale. Deux paramètres sont importants : l’un est la déviation standard de
la hauteur de l’aspérité, l’autre est le nombre de sommets sur l’interface rugueuse, c’est-à-
Étude de la géométrie à l’interface 113

F IG . 4.46: La tendance entre la résistance de contact et la rugosité différente

dire, la densité des aspérités. Les différentes hauteurs des aspérités ainsi que leur nombre
peuvent influencer l’état de contact à l’interface et la résistance. Les hauteurs des aspérités
à l’interface rugueuse sont indiquées dans le tableau 4.1 selon l’expérience réalisée par Bes-
mann et al [BKH03]. L’aspérité de 8 mm de la surface de la plaque GDL est donnée dans le
tableau 1.3 dans le chapitre 1.

TAB . 4.1: Aspérités des plaques BPP et GDL

Aspérité de la BPP (µm) 2.9 8.83 10.2 [BKH03]

Aspérité de la GDL (µm) 8 8 8 [Jap08]

Le modèle des éléments finis avec une interface rugueuse est représenté dans la figure
4.47. Ce modèle n’est pas symétrique car la rugosité ne l’est pas. La surface haute de la
plaque GDL est rugueuse avec une rugosité de 8 µm, et la surface basse de la plaque BPP
a une rugosité de : 2.9 µm, 8.83 µm ou 10.2 µm. Le maillage du modèle est indiqué dans
la figure 4.48. Le maillage est plus fin que le modèle utilisé avec l’interface lisse pour tenir
compte des petits aspérités sur les surfaces rugueuses. Les propriétés et les conditions aux
limites du modèle sont identiques à celles des études précédents.
Dans un premier temps, l’influence des différentes aspérités de la surface de la plaque
BPP sur l’aire de contact est étudiée. La figure 4.49 montre la variation de l’aire de contact
avec trois types d’aspérité (2.9 µm, 8.3 µm et 10.2 µm). L’aire de contact avec une surface
lisse est différente de celle avec une surface rugueuse, car quand la précompression est faible,
114 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

F IG . 4.47: Modélisation des rugosités

F IG . 4.48: Maillage de l’interface rugueuse

Étude de la géométrie à l’interface 115

l’aspérité produit des espaces non en contact entre les plaques BPP et GDL. Mais quand la
précompression est assez élevée comme 0.1 mm ou plus, l’aire de contact avec interface
rugueuse est proche de l’aire de contact avec interface lisse. Une précompression élevée
diminue donc les aspérités à l’interface et conduit à l’augmentation de l’aire de contact.
Plus les aspérités sont faibles, plus de l’aire de contact se rapproche de celle obtenue avec
une surface lisse. En revanche, les aspérités importantes défavorisent le taille de la zone de
contact, ceci d’autant plus que la pré-compression est faible.
Quand la précompression de la PAC est assez élevée, l’influence des aspérités sur l’aire
de contact se fait moindre, parce que la précompression élevée réduit les aspérités et les met
en contact avec la plaque GDL.

F IG . 4.49: Influence des aspérités sur l’aire de contact

La figure 4.50 montre que l’influence de la densité des aspérités de la plaque BPP sur
l’aire de contact avec une taille des aspérités de 2.9 mm. La densité des aspérités est le
nombre des aspérités par unité de surface. Les aires de contact avec les densités d’aspérités
de 20/mm, 40/mm et 80/mm sont représentées dans la figure 4.50. Nous pouvons voir que
plus il y a d’aspérités à l’interface, moins l’aire de contact est importante, plus d’aspérités
produit plus d’espaces à l’interface et ces espaces ne sont pas en contact.
Enfin, dans la figure 4.50, une précompression de la PAC élevée augmente toujours l’aire
de contact. Ainsi une précompression élevée et une petite densité d’aspérités à l’interface
permettent de réaliser une aire de contact plus large. La précision de la fabrication de la
plaque BPP est donc importante dans la performance de la PAC. Il faut avoir la plus faible
densité d’aspérités possible.
116 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

F IG . 4.50: Influence de la densité des aspérités

Pour conclure, on peut dire que l’interface entre la plaque GDL et la plaque BPP doit
être la plus lisse possible, pour avoir une aire de contact large et une performance de la
PAC supérieure. La précision dans la fabrication de la surface de la plaque BPP est donc
très importante. Toutefois, si les conditions de la fabrication ne le permettent pas (coût de
fabrication élevé), une précompression assez élevée devra être appliquée afin de diminuer les
effets de la rugosité à l’interface.

4.7.2 Influence des hauteurs de dents de la BPP

Les études précédentes ont été réalisées avec un modèle de contact singulier, c’est-à-
dire, un modèle d’un contact entre une seule dent de la plaque BPP et la plaque GDL. Ce
modèle permet de représenter les contacts multiples entre les dents de la plaque BPP et la
plaque GDL à condition que chaque dent soit identique à celle du modèle. Or, la plaque BPP,
en graphite, est un matériau rigide avec de nombreuse dents et canaux, qui sont difficiles
à graver. La précision de la configuration de la plaque est donc difficile à contrôler. Les
hauteurs des dents de la plaque BPP doivent notamment être bien contrôlées, car elle peuvent
avoir une influence sur l’état de contact.
Normalement la plaque BPP a une hauteur unique DCN et une épaisseur unique TBP ,
ainsi les surfaces basses de toutes les dents de la plaque BPP sont à la même position hori-
zontale. Si une différence existe parmi les hauteurs de dent, la précompression uniforme de
la PAC peut produire un contact non-uniforme aux interfaces multiples. Ainsi l’influence de
la tolérance que l’on peut avoir sur les hauteurs de dents de la BPP est intéressante à étudier.
Étude de la géométrie à l’interface 117

Sur la figure 4.51, les hauteurs des dents de la BPP ont été modifiées et suivent une dis-
tribution Gaussienne normale. Le maillage et les propriétés du modèle de contact multiples
sont identiques au modèle avec une seule dent.

F IG . 4.51: Modèle avec différentes hauteurs de dent

Dans un premier temps, la distribution de la pression de contact entre la GDL et la BPP

est présentée dans la figure 4.52. Sous une précompression uniforme, la distribution de la
pression dans chaque zone de contact est évidemment différente car les différentes hauteurs
de dents de la plaque BPP conduisent à différents états de contact. Certaines zones de contact
sont totales pour des hauteurs élevées, certaines zones de contact sont incomplètres quand
la hauteur est courte. Sous un chargement la précompression uniforme, c’est la dent avec la
hauteur la plus élevée parmi les dents de la plaque BPP qui est la première en contact avec
la plaque GDL. Il est évident que les différents états de contact conduisent à différentes aires
de contact et à différentes déformations de la plaque GDL.

F IG . 4.52: Pression avec différentes dents

La figure 4.53 montre la distribution de la déformation de la plaque GDL pour diffé-

rentes hauteurs des dents de la plaque BPP. La déformation de la plaque GDL est différente
pour chaque dent, ce qui conduit à une porosité non uniforme. Pour des hauteurs élevées, la
contrainte élevée conduit à une déformation élevée et à une faible porosité de la plaque GDL.
En revanche, les zones tout juste en contact obtenues pour de faibles hauteurs des dents ne
sont pas fortement compressées et la porosité y est donc meilleure.
L’aire de contact, la pression de contact au point 1 et la pression de contact sur l’axe de
symétrie de chaque dent sont indiquées dans le tableau 4.2.
118 Chapitre 4 : Influence des principaux paramètres sur la performance des piles

F IG . 4.53: Déformation de la plaque GDL avec différentes hauteurs des dents

F IG . 4.54: Numéro des dents

TAB . 4.2: Variation des états de contact sur les dents

No. des Aire de Pression au Pression sur l’axe Différence de Hauteurs

dents contact (mm) point 1 (Mpa) de symétrie (Mpa) pression (Mpa) des dents (mm)
1 0.85704 4.36 2.26 2.10 1.046
2 0.80469 2.62 1.35 1.27 0.971
3 0.85704 4.14 2.13 2.01 1.017
4 0.85704 4.05 2.08 1.97 1.005
5 0.85704 3.79 1.93 1.86 1.002
6 0.80469 2.43 1.22 1.21 0.97
7 0.80469 3.24 1.68 1.56 0.979
8 0.80469 2.19 1.03 1.16 0.966
9 0.83087 3.82 1.78 2.04 1.003
10 0.85704 4.36 2.27 2.09 1.046
11 0.80469 3.34 1.63 1.71 0.993
12 0.85705 3.95 1.94 2.01 1.01

La différence des hauteurs de dents de la plaque BPP conduit à différents états de contact
entre toutes les petites dents de la plaque BPP et la plaque GDL : la pression et l’aire de
contact sont différentes, ce qui a une influence sur la performance des piles. En effet, les
résistances de contact sont différentes, les températures aussi. Tout ceci n’est pas favorable à
une performance optimale de la PAC et peut agir sur la durée de vie des piles.
 Influence des principaux paramètres sur la performance des piles 119

4.8 Conclusions
Dans ce chapitre, une étude est effectuée sur les influences des paramètres variables at-
tribués aux effets mécaniques aux interfaces de la structure de la PAC. La variation des
paramètres comme la précompression, le rayon du congé et la largeur des dents de la plaque
BPP, les épaisseurs de GDL et MEA est étudiée. On a trouvé que la précompression a une
influence importante pour le contact à l’interface, la précompression de 0.1 mm permet de
réaliser un contact complet entre BPP/GDL, d’obtenir une bonne uniformité de la pression
de contact à l’interface, et la porosité la plus grande possible pour la plaque GDL. Le rayon
du congé de dents change la pression de contact, la longueur de contact entre BPP/GDL et
la déformation de la plaque GDL. On a trouvé qu’un rayon de congé de 0.05 mm est un bon
compromis pour la performance de la PAC. En ce qui concerne l’influence des épaisseurs de
GDL et MEA, une plaque GDL épaisse de 0.375 mm permet de garder une grande porosité.
De plus, une membrane de 0.05 mm est considérée grâce à sa petite résistance propre pour
la performance de la PAC. L’étude de l’influence de la géométrie à l’interface a montré que
les aspérités réduisent la surface de contact et diminuent donc la performance de la pile. Il
en est de même avec les irrégularité des hauteurs de dents de la plaque BPP.
Ces défauts de fabrication devrons donc être bien contrôlés au moment de la mise en
œuvre de la pile.
120 Influence des principaux paramètres sur la performance des piles
 121

Conclusions et perspectives

Conclusions
Le but de ce travail était de développer un modèle de la pile basé sur la méthode des
éléments finis et permettant d’étudier le comportement mécanique au niveau des interfaces
multi-contacts dans la pile à combustible. La pile de type PEM ayant une structure lamel-
laire, tous les éléments sont en contact les uns avec les autres. Il est évident que le contact
joue un rôle important dans le fonctionnement et la performance des piles. Aussi nous a-t-il
paru pertinent d’approfondir les connaissances sur ces différents contacts et de proposer des
valeurs optimales pour certains paramètres qui pourraient être prises en comptes au moment
de la conception.

Tout au long de ce mémoire nous avons présenté des points importants qui ont permis
de mieux connaître les phénomènes de contact au sein d’un pile. En résumé, nous pouvons
définir trois points pour lesquels les études dans ce travail sont originales :

– le traitement du contact : pour ce premier point, nous avons utilisé la méthode du bi-
potentiel développée par De Saxcé et Feng pour la résolution des problèmes de contact
avec frottement en statique dans la structure de la pile.

– le rôle du contact dans la pile à combustible : pour ce second point, nous avons étu-
dié le rôle du contact dans la performance de la pile à travers : la pression de contact,
la zone de contact et la déformation sur la zone de contact.

– L’influence de certains paramètres sur la performance : pour ce troisième point,

nous avons étudié l’influence de différents paramètres sur les interfaces de la pile à
combustible : précompression, rayon du congé et largeur des dents de la BPP, épais-
seurs de la plaque GDL et MEA, etc. afin de proposer une structure optimale pour la
performance de la pile à combustible.

Comme nous l’avons présenté dans la chapitre 2, la résolution des problèmes de contact
122

par la méthode du bi-potentiel, développée par De Saxcé et Feng [Fen95], est précise et
rapide.
Dans le chapitre 3, la modélisation de la plaque bipolaire (BPP), de la couche de diffusion
gazeuse (GDL) et de la membrane échangeuse de protons (MEA) a permis de montrer que
des zones de décollements existaient entre la MEA et la GDL. Considérant ensuite que ces
deux éléments étaient collés (physiquement possible par un recuit), nous avons constaté que
la dent rectangulaire avec angle droit engendrait une répartition de la pression normale de
contact d’allure hyperbolique à l’approche du coin. L’introduction d’un congé à la place de
l’angle droit a permis d’atténuer de manière conséquente cette concentration de contraintes.
Une modélisation prenant en compte la porosité des matériaux de la GDL a ensuite montré
qu’une déformation importante diminuait la porosité et génait la circulation des gaz réactifs.
Le chapitre 4 a été consacré à l’influence de certains paramètres de chargement ou d’ordre
géométrique sur le fonctionnement de la pile. Des valeurs optimales ont été proposées pour
les paramètres suivants : précompression, largeur et rayon du congé des dents de la BPP,
épaisseur de la GDL et MEA. Nous avons également vu que les aspérités de surface pro-
venant d’une fabrication non parfaite des surface lisses diminuaient indirectement la perfor-
mance de la pile. Ces valeurs pouront par la suite être utilisées pour la conception des piles
et pour de futures études sur les piles à combustibles.

Il est à noter que ce travail a donné lieu à une publication dans un journal scientifique
[ZRFY10] et deux communications dans des congrès [ZRF09b][ZRF09a].

Perspectives
Les différents résultats obtenus dans ce travail nous encouragent fortement à poursuivre
dans cette voie. De nombreuses modifications peuvent être envisagées à plus ou moins court
terme dans le but d’étendre nos connaissances sur la façon d’améliorer le fonctionnement
des piles à combustible. Ces recherches à venir peuvent être regroupées en deux parties in-
dépendantes : l’optimisation de la structure multi-contact de la PAC et la modélisation des
phénomènes couplés (mécano-électro-chimiques) au sein de la pile.

La première partie apparaît comme une suite logique du dernier chapitre de cette thèse.
En effet nous avons vu que les effets inter-agissaient parfois de manière contraire, et qu’il
fallait trouver des compromis. Il serait donc très intéressant de suivre un plan d’expériences,
appliqué au numérique, pour dégager les paramètres les plus influents sur la performance de
la pile et ainsi proposer une optimisation globale de la pile. Nous pourrions peut être ajouter
de nouveaux paramètres à ceux déjà étudiés comme par exemple, le profondeur DCN des
dents de la plaque BPP afin d’étudier l’influence de DCN sur le flux des gaz réactifs.
 123

Dans le second axe de recherche, nous envisageons d’étendre le champ de l’étude mé-
canique aux phénomènes chimiques et électroniques présents dans la pile à combustible.
Différents couplages peuvent être étudiés :

– mécanique-électronique : nous avons vu dans le chapitre 3 que la pression de contact

et la zone de contact sont liées à la résistance de contact et au courant électrique de
la pile. De nombreux développements peuvent être apportés au niveau du choix d’un
modèle entre les paramètres mécaniques et le phénomène électrique. En effet, chaque
domaine d’application propose ses spécificités et il est nécessaire d’adapter au mieux
notre modélisation aux cas étudiés.

– mécanique-chimique : dans notre travail dans le chapitre 4, nous avons observé l’in-
fluence des paramètres mécaniques sur la porosité de la plaque GDL. Nous avons vu
que l’alimentation des gaz réactifs était influencée par la déformation mécanique de
la plaque GDL et par la largeur des dents de la plaque BPP. Cependant, il existe bien
d’autres phénomènes chimiques pour lesquels l’étude mécanique-chimique serait es-
sentielle.

Nous espérons que les travaux présentés dans ce mémoire pourront constituer une base
à de futures recherches dans le domaine de la modélisation des phénomènes couplés méca-
niques, électriques, chimiques présents au sein d’une pile à combustible.
124
 125

Références bibliographiques

[AC91] P. Alart and A. Curnier. A mixed formulation for frictional contact problems
prone to Newton like solution methods. Comp. Meth. Appl. Mech. Engng.,
92, :353–375, 1991.
[All09] SGL Co. Allemagne. http ://www.sglcarbon.com, 2009.
[APT02] E. Antolini, R. R. Passos, and E. A. Ticianelli. Effects of the carbon powder
characteristics in the cathode gas diffusion layer on the performance of poly-
mer electrolyte fuel cells. Journal of Power Sources, 109 :477–482, 2002.
[AUW05] M. Aoki, H. Uchida, and M. Watanabe. Novel evaluation method for de-
gradation rate of polymer electrolytes in fuel cells. Electrochem Commun,
12 :1434–1438, 2005.
[Bat82] K. J. Bathe. Finite element procedures in engineering analysis. Prentice-
Hall, 1982.
[BGGM08] D. Bograchev, M. Gueguen, J.-C. Grandidier, and S. Martemianov. Stress
and plastic deformation of MEA in fuel cells-stress generated during cell
assembly. Journal of Power Sources, 180 :393–401, 2008.
[BKH03] T. M. Besmann, J. W. Klett, and J. J. Henry. Carbon composite bipolar plates.
Progress Report, 1-4, 2003.
[BL05] J. J. Baschuk and X. Li. A general formulation for a mathematical PEM fuel
cell model. Journal of Power Sources, 134-153, 2005.
[BN91] T. Belytschko and M. O. Neal. Contact-impact by the pinball algorithm with
penalty and lagrangian methods. Int. J. Num. Meth. Engng., 31 :547–572,
1991.
[BSP+ 08] M. Barber, T. S. Sun, E. Petrach, X. Wang, and Q. Zou. Contact mecha-
nics approach to determine contact surface area between bipolar plates and
current collector in proton exchange membrane fuel cells. Journal of Power
Sources, 185 :1252–1256, 2008.
[Can08] Ballard Co. Canada. http ://www.ballard.com, 2008.
126

[CB86] A. B. Chaudhary and K. J. Bathe. A solution method for static and dynamic
analysis of three dimensional contact problems with friction. Computers &
Structures, 24 :855–873, 1986.
[CB97] T. R. Chandrupatla and A. D. Belegundu. Introduction to finite elements in
engineering. Prentice-Hall, Upper Saddle Rive, 1997.
[CEA06] http ://www.cea.fr, 2006.
[CFC04] J. M. Corréa, F. A. Farret, and L. N. Canha. An electrochemical based fuel
cell model suitable for electrical engineering automation approach. IEEE
Trans. on Industrial Electronics, 1103 - 1112, 2004.
[CHW97] C. Cavalca, S. T. Homeyer, and E. Walsworth. Flow field plate for use in a
proton exchange membrane fuel cell. US Patent, No. 5,686,199, 1997.
[CHWC07] W. R. Chang, J. J. Hwang, F. B. Weng, and S. H. Chan. Effect of clamping
pressure on the performance of a PEM fuel cell. Journal of Power Sources,
115 :149–154, 2007.
[Cia88] P. G. Ciarlet. Mathematical Elasticity, Vol. 1 : Three-dimensional elasticity.
North-Holland, 1988.
[Cia90] P. G. Ciarlet. Plates and Junctions in Elastic multi-structures. SpringerVer-
lag, 1990.
[Cia97] P. G. Ciarlet. Mathematical elasticity, Vol.2 : Theory of plates. Elsevier
Science Publishers, 1997.
[CKPS98] P. W. Chritensen, A. Klarbring, J. S. Pang, and N. N. Strömberg. Formulation
and comparison of algorithms for frictional contact problems. Int. J. Numer.
Meth. Engng., 42 :145–173, 1998.
[CT71] S. H. Chan and I. S. Tuba. A finite element method for contact problems of
solid bodies. Int. J. Mech. Sci., 13 :615–639, 1971.
[CTK91] N. J. Carpenter, R. L. Taylor, and M. G. Katona. Lagrange constraints for
transient finite element surface contact. Int. J. Num. Meth. Engng., 32 :103–
128, 1991.
[Cur84] A. Curnier. A theory of friction. Int. J. Solids Structures, 20(7) :637–647,
1984.
[dSF91] G. de Saxcé and Z.-Q. Feng. New inequality and functional for contact
with friction : The implicit standard material approach. Mech. Struct. Mach.,
19 :301–325, 1991.
[dSF98] G. de Saxcé and Z.-Q. Feng. The bi-potential method : a constructive ap-
proach to design the complete contact law with friction and improved nume-
 127

rical algorithms. Mathematical and Computer Modeling, 28(4-8) :225–245,

1998.
[DuP06] DuPont. http ://www2.dupont.com, 2006.
[FCP+ 96] N. J. Fletcher, C. Y. Chow, E. G. Pow, H. H. Wozniczka, G. H. Voss, and
D. P. Wilkinson. Canadian Patent, No. 2,192,170, 1996.
[Fen95] Z.-Q. Feng. 2D or 3D frictional contact algorithms and applications in a large
deformation context. Comm. Numer. Meth. Engng., 11 :409–416, 1995.
[FPA06] J. Feser, A. Prasad, and S. Advani. Experimental characterization of in-plane
permeability of gas diffusion layers. Journal of Power Sources, 162 :1226–
1231, 2006.
[FPH06] Z.-Q. Feng, F. Peyraut, and Q.-C. He. Finite deformations of Ogden’s mate-
rials under impact loading. Int. J. Non-linear Mech., 41 :575–585, 2006.
[FPZ90] D. François, A. Pineau, and A. Zaoui. Comportement mécanique des mate-
riaux. Hermes, 1990.
[Fri01] F. Frisch. PEM fuel cell stack sealing using silicone elastomers. Sealing
Technology, 93 :7–9, 2001.
[Fun65] Y. C. Fung. Foundations of solid mechanics. Prentice-Hall, 1965.
[FVFP06] Z.-Q. Feng, C. Vallée, D. Fortuné, and F. Peyraut. The 3é hyperelastic model
applied to the modeling of 3D impact problems. Finite Elements in Analysis
and Design, 43 :51–58, 2006.
[Gar99] J. Garrigues. Eléments finis, http ://esm2.imt-mrs.fr/gar/ef.html, 1999.
[GFP+ 06] J. T. Gostick, M. W. Fowler, M. D. Pritzker, M. A. Ioannidis, and L. M.
Behra. Experimental characterization of in-plane permeability of gas diffu-
sion layers. Journal of Power Sources, 162 :228–238, 2006.
[GHL06] J. Ge, A. Higier, and H. Liu. Effect of gas diffusion layer compression on
PEM fuel cell performance. Journal of Power Sources, 922-927, 2006.
[GHV69] J. M. Grimwood, B. C. Hacker, and P. J. Vorzimmer. Project gemini techno-
logy and operations :. A Chronology NASA, SP-4002, 1969.
[GN60] W. T. Grubb and L. W. Niedrach. Batteries with solid ion-exchange mem-
brane electrolytes. Jounal of The Electrochimical Society, 107 :131–135,
1960.
[Gro39] W. R. Grove. On voltaic series and the combination of gases by platinum.
The London and Edimburgh philosophical magazine and journal of science,
14 :127–130, 1839.
128

[HC93] Q.-C. He and A. Curnier. Anisotropic dry friction between two orthotropic
surfaces undergoing large displacements. European Journal of Mechanics
-A/Solids, 12 :631–666, 1993.
[HCCP05] R. Hayre, S. W. Cha, W. Colella, and F. B. Prinz. Fuel Cell Fundamentals.
Wiley, 2005.
[HFdM04] M. Hjiaj, Z.-Q. Feng, G. de Saxcé, and Z. Mróz. There-dimensional finite
element computations for frictional contact problems with non-associated
sliding rule. Int. J. Numer. Meth. Engng., 60 :2045–2076, 2004.
[Hoo03] G. Hoogers. Fuel cell Technology Handbook. CRC Press, 2003.
[HP98] P. Hunter and A. Pullan. FEM/BEM notes. Bioengineering Research Group,
Departement of Engineering Science, University of Auckland(New Jersey),
1998.
[HSK07] A. Husar, M. Serra, and C. Kunusch. Description of gasket failure in a 7 cell
PEMFC stack. Journal of Power Sources, 169 :85–91, 2007.
[HTS+ 76] T. J. R. Hughes, R. L. Taylor, J. L. Sackman, A. Curnier, and W. Kanok-
nukulchai. A finite element method for a class of contact-impact problems.
Comp. Meth. Appl. Mech. Engng., 8 :249–276, 1976.
[Hug87] T. J. R. Hughes. The finite element method : Linear static and dynamic finite
element engng. Pretic-Hall,Inc, Enlewood cliffs :New Jersey, 1987.
[Hun92] I. Huněk. On a penalty formulation for contact-impact problems. Computers
& Structures, 48(2) :192–203, 1992.
[HWDA07] J. E. Hensleya, J. D. Waya, S. F. Decb, and K. D. Abneyc. The effects of ther-
mal annealing on commercial nafion membranes. J. Membr. Sci., 298 :190–
201, 2007.
[IML04] J. Ihonen, M. Mikkola, and G. Lindbergh. Flooding of gas diffusion backing
in PEMFCs physical and electrochimical characterization. J. Electrochemi-
cal Society, 151 :1152–1161, 2004.
[Iqb03] M. T. Iqbal. Modeling and control of a wind fuel cell hybrid energy system.
Renewable Energy, 223 - 237, 2003.
[Jap08] Toray Co. Japon. http ://www.toray.com, 2008.
[JL09] K. Jiao and X. G. Li. Effects of various operating and initial conditions on
cold start performance of polymer electrolyte membrane fuel cells. Int. J.
Hydr. Energ., 34 :8171–8184, 2009.
[JSB+ 00] L. R. Jordan, A. K. Shukla, T. Behrsing, N. R. Avery, B. C. Muddle, and
M. Forsyth. Diffusion layer parameters influencing optimal fuel cell perfor-
mance. Journal of Power Sources, 86 :250–254, 2000.
 129

[KFT08] J. Kleemanna, F. Finsterwaldera, and W. Tillmetz. Characterisation of me-

chanical behaviour and coupled electrical properties of polymer electrolyte
membrane fuel cell gas diffusion. Journal of Power Sources, 153 :130–135,
2008.
[KH98] T. H. Ko and L. C Huang. The influence of cobaltous chloride modifica-
tion on physical properties and microstructure of modified PAN fiber during
carbonization. J. Appl. Polymer Sci., 70 :2409–2415, 1998.
[KKS+ 06] A. Kusoglu, A. M. Karlsson, M. H. Santare, S. Cleghorn, and W. B. John-
son. Mechanical response of fuel cell membranes subjected to hygro-thermal
loading. Journal of Power Sources, 161 :987–996, 2006.
[KKS+ 07] A. Kusoglu, A. M. Karlsson, M. H. Santare, S. Cleghorn, and W. B. Johnson.
Mechanical behavior of fuel cell membranes under humidity cycles and ef-
fect of swelling anisotropy on the fatigue stresses,. Journal of Power Sources,
170 :345–358, 2007.
[Kla92] A. Klarbring. Mathematical programming and augmented Lagrangian me-
thods for frictional contact problems. In A. Curnier, editor, Contact Me-
chanics International Symposium, pages 409–422. Presses Polytechniques et
Universtaires Romandes, Lausanne, Switzerland, 1992.
[KS81] N. Kikuchi and Y. J. Song. Penalty/finite element approximations of a class
of unilateral problems in linear elasticity. Quaterly of Applied Mathematics,
39 :1–22, 1981.
[LD03] J. Larminie and A. Dicks. Fuel cell systems explained. John Wiley and Sons,
Inc., New York, 2003.
[LHH05] S. J. Lee, C. D. Hsu, and C. H. Huang. Analyses of the fuel cell stack assem-
bly pressure. Journal of Power Sources, 145 :353–361, 2005.
[LHVZM99] W. K. Lee, C. H. Ho, J. W. Van Zee, and M. Murthy. The effects of compres-
sion and gas diffusion layers on the performance of a PEM fuel cell. Journal
of Power Sources, 84 :45–51, 1999.
[LLA98] J. H. Lee, T. R. Lalk, and A. J. Appleby. Modeling electrochemical perfor-
mance in large scale proton exchange membrane fuel cell stacks. Journal of
Power Sources, 258 - 268, 1998.
[LM04] S. Litster and G. McLean. PEM fuel cell electrodes. Journal of Power
Sources, 130 :61–76, 2004.
[MABHAJ07] A. R. Maher, S. Al-Baghdadi, A. K. Haroun, and S. Al-Janabi. Influence of
the design parameters in a proton exchange membrane (PEM) fuel cell on the
130

mechanical behavior of the polymer membrane. Journal of Power Sources,

21 :2258–2267, 2007.
[Mac09] C. S. MacDonald. Effect of Compressive Force on PEM Fuel Cell. PhD
thesis, University of Waterloo, 2009.
[MC03] V. Mehta and J. S. Cooper. Review and analysis of PEM fuel cell design and
manufacturing. Journal of Power Sources, 114 :32–53, 2003.
[Men00] T. Mennola. Design and Experimental Characterization of Polymer Electro-
lyte Membrane Fuel Cells. PhD thesis, Helsinki University of Technology,
2000.
[Mor88] J. J. Moreau. Unilateral contact and dry friction in finite freedom dyna-
mics. In J.J. Moreau and P.D. Panagiotopoulos, editors, Non Smooth Me-
chanics and Applications, CISM Courses and Lectures, volume 302, pages
1–82. Springer-Verlag, Wien, New York, 1988.
[Mor94] J. J. Moreau. Some numerical methods in multibody dynamics : application
to granular materials. Eur. J. Mech., A/Solids, 13 :93–114, 1994.
[Mos03] R. Mosdale. Transport électrique routier, véhicules électriques à piles à com-
bustible. Techniques de l’ingénieur, 12 :1–16, 2003.
[MYP04] V. Mishra, F. Yang, and R. Pitchumani. Measurement and prediction of elec-
trical contact resistance between gas diffusion layers and bipolar plate for
applications to PEM fuel cells. Journal of Fuel Cell Sci. Technol., 1 :2–9,
2004.
[MYP05] V. Mishra, F. Yang, and R. Pitchumani. Analysis and design of PEM fuel
cells. Journal of Power Sources, 141 :47–64, 2005.
[Nes04] M. Nesme. Modèle déformable par éléments finis élasticité linéaire et grands
déplacements, pour application à la simulation chirurgicale. PhD thesis, Ins-
titut National Polytechnique de Grenoble, 2004.
[NHM08] I. Nitta, O. Himanen, and M. Mikkola. Contact resistance between gas dif-
fusion layer and catalyst layer of PEM fuel cell. Electrochemistry Commu-
nications, 10 :47–51, 2008.
[Ode67] J. T. Oden. Mechanics of elastic structures. Mc Graw-Hill, 1967.
[Ogd84] R. W. Ogden. Non-linear elastic deformations. Ellis Horwood, 1984.
[Par89] H. A. Parisch. Consistent tangent stiffness matrix for three-dimensional non-
linear contact analysis. Int. J. Numer. Meth. Engng., 28 :1803–1812, 1989.
[Par09] J. Parsons. Low cost, durable seals for PEM fuel cells,
http ://www.hydrogen.energy.gov/, 2009.
 131

[PDA01] G. Picinbono, H. Delingette, and N. Ayache. Non-linear and anisotropic

elastic soft tissue models for medical simulation. ICRA 2001 (IEEE Interna-
tional Conference on Robotics and Automation, Seoul, Korea, 21-26 May),
pages 1370–1375, 2001.
[PL04] M. A. Puso and T. A. Laursen. A mortar segment-to-segment contact method
for large deformation solid mechanics. Comp. Meth. Appl. Mech. Engng.,
193 :601–629, 2004.
[Pou63] M. Pourbaix. Atalas d’equilibres electrochimiques. In Engineering and Tech-
nology, Gauthier-Villars, 1963.
[PS80] A. Pellegri and P. M. Spaziante. Bipolar separator for electrochemical cells
and method of preparation thereof. US Patent, No.4,197,178, 1980.
[PS98] P. Papadopoulos and J. M. Solberg. A Lagrange multiplier method for the
finite element solution of frictionless contact problems. Mathematical and
Computer Modelling, 28 :373–384, 1998.
[PTG95] V. Paganin, E. A. Ticianelli, and E. R. Gonzalez. Development and elec-
trochemical studies of gas diffusion electrodes for polymer electrolyte fuel
cells. Journal of Applied Electrochemistry, 26 :297–304, 1995.
[RA09] http ://bama.ua.edu/∼rreddy/projects/fuel-cells.htm, 2009.
[RF03] C. Renaud and Z.-Q. Feng. BEM and FEM analysis of Signorini contact
problems with friction. Comput. Mech., 31 :390–399, 2003.
[Sch97] G. G. Scherer. Interfacial aspects in the development of polymer electrolyte
fuel cells. Solid State Ionics, 94 :249–257, 1997.
[SGVZ07] S. Shimpalee, S. Greenway, and J. W. Van-Zee. Numerical studies on rib
channel dimension of flow-field on PEMFC performance. Journal of Power
Sources, 32 :842–856, 2007.
[SH04] V. Stanic and M. Hoberecht. Mechanism of pin-hole formation in membrane
electrode assemblies for PEM fuel cells. Proceedings of 4th International
Symposium on Proton Conducting Membrane Fuel Cells, 1891, 2004.
[SL92] J. C. Simo and T. A. Laursen. An augmented Lagrangian treatment of contact
problems involving friction. Computers & Structures, 42 :97–116, 1992.
[SLC+ 08] Z. Y. Su, C. T. Liu, H. P. Chang, C. H. Li, K. J. Huang, and P. C. Sui. A
numerical invertigation of the effects of compression force on PEM fuel cell
performance. Journal of Power Sources, 183 :182–192, 2008.
[SLVZ06] S. Shimpalee, W. K. Lee, and J. W. Van-Zee. The impact of channel path
length on PEMFC flow field design. Journal of Power Sources, 160 :398–
406, 2006.
132

[SLVZNN01] S. Shimpalee, W. K. Lee, J. W. Van-Zee, and H. Naseri-Neshat. Advances in

computational fluid dynamics modeling for PEM fuel cells. Proceedings of
ASME IECEC IECEC-2001-10, Savannah :GA, 2001.
[SNCH00] P. Stevens, F. Novel-Cattin, and A. Hammou. Pile à combustible, Techniques
de l’Ingénieur. Trainté Génie Electrique, 2000.
[SP98] A. M. Solberg and P. Papadopoulos. A finite element method for
contact/impact. Finite Elements in Analysis and Design, 30 :297–311, 1998.
[Spe84] A. J. M. Spencer. Continum Thoery of Fiber-Reinforced Composites. Spring-
Verlag, New York, 1984.
[Sri06] S. Srinivasan. Fuel Cells : From Fundamentals to Applications. Springer
Science, 2006.
[TR05] M. M. Tomadakis and T. J. Robertson. The viscous permeability of random
fiber structure : comparison of elctrical and diffusional estimates with expe-
rimental and analytical results. J. composite materials, 163-187, 2005.
[TY73] T. Tsuta and S. Yamaji. Finite element analysis of contact problem. In Tokyo
seminar on element analysis, editor, Theory and practice in finite element
structural analysis, pages 178–194, 1973.
[UFS+ 96] M. Uchida, Y. Fukuoka, Y. Sugawaya, N. Eda, and A. Ohta. Effects of mi-
crostructure of carbon support in the catalyst layer on the performance of
polymer-electrolyte fuel cells. J. Electrochemical Society, 143 :2245–2252,
1996.
[Vee07] http ://www.veeco.com, 2007.
[WBB+ 04] M. Williams, E. Begg, L. Bonville, H. R. Kunz, and J. M. Fenton. Characte-
rization of gas diffusion layers for PEMFC. Journal of the Electrochemical
Society, 8 :1173–1180, 2004.
[WDE91] D. S. Watkins, K. W. Dircks, and D. G. Epp. Novel fuel cell fluid flow field
plate. US Patent, No. 4,988,583, 1991.
[WDE92] D. S. Watkins, K. W. Dircks, and D. G. Epp. Fuel cell fluid flow field plate.
US Patent, No. 5,108,849, 1992.
[WN04] A. Z. Weber and J. Newman. A theoretical study of membrane constraint in
polymer-electrolyte fuel cells. AIChE Journal, 50 :3215–3226, 2004.
[Wri02] P. Wriggers. Computational contact mechanics. John Wiley & Sons, 2002.
[WST03] H. L. Wang, M. A. Sweikart, and J. A. Turner. Stainless steel as bipolar
plate material for polymer electrolyte membrane fuel cells. Journal of Power
Sources, 115 :243–251, 2003.
 133

[WSZ08] X. T. Wang, Y. Song, and B. Zhang. Experimental study on clamping pres-

sure distribution in PEM fuel cells. Journal of Power Sources, 179 :305–309,
2008.
[Zie77] O. C. Zienkiewicz. The Finite Element Method. Mc Graw-Hill, London,
1977.
[ZL06] T. Zhou and H. Liu. Effects of the electrical resistances of the gdl in a PEM
fuel cell. Journal of Power Sources, 161 :444–453, 2006.
[ZLS+ 06] L. H. Zhang, Y. Liu, H. M. Song, S. X. Wang, Y. Y. Zhou, and S. J. Hu.
Estimation of contact resistance in proton exchange membrane fuel cells.
Journal of Power Sources, 162 :1165–1171, 2006.
[ZRF09a] Z. M. Zhang, C. Renaud, and Z.-Q. Feng. Numerical analysis of mechanical
multi-contacts on the interfaces in a PEM fuel cell stack. In Computational
Structure Engineering, Shanghai-Chine, 2009.
[ZRF09b] Z. M. Zhang, C. Renaud, and Z.-Q. Feng. Simulation numérique des in-
terfaces multi-contacts dans une pile à combustible de pem. In Colloque
National en calcul des structures, Giens-France, 2009.
[ZRFY10] Z. M. Zhang, C. Renaud, Z.-Q. Feng, and H. P. Yin. Mechanical contact ana-
lysis on the interfaces in a proton exchange membrane fuel cell. Mécanique
et Industries, In Press, 2010.
[ZWHT84] O. C. Zienkiewicz, W. L. Wood, L. W. Hine, and R. L. Taylor. A unified set
of single step algorithms. part 1. general formulation and application. Int. J.
Numer. Meth. Engng., 20 :1529–1552, 1984.
[ZWM06] P. Zhou, C. W. Wu, and G. J. Ma. Contact resistance prediction and structure
optimization of bipolar plates. Journal of Power Sources, 159 :1115–1122,
2006.
134
 135

Table des figures

1.1 Pollution environnementale . . . . . . . . . . . . . . . . . . . . . . . . . . 8

1.2 Accroissement de CO2 depuis le 19ème siècle . . . . . . . . . . . . . . . . 8
1.3 Pile inventée par Grove [Gro39] . . . . . . . . . . . . . . . . . . . . . . . 10
1.4 Francis Thomas Bacon et la pile à combustible [Sri06] . . . . . . . . . . . 11
1.5 Pile à combustible de PEM utilisée dans le programme Gemini [GHV69] . 11
1.6 Schéma d’une MEA relatif au principe de fonctionnement d’une pile de PEM 14
1.7 Courbe de polarisation de la pile PEM . . . . . . . . . . . . . . . . . . . . 16
1.8 Schéma de la structure d’un empilement de la pile de PEM . . . . . . . . . 18
1.9 Assemblage d’une cellule de la PAC de PEM avec tous les composants . . . 19
1.10 Processus de l’assemblage d’une cellule de la PAC . . . . . . . . . . . . . 19
1.11 Schéma du processus opératoire de la réaction de la pile . . . . . . . . . . . 20
1.12 Profil de la MEA (électrolyte, anode et cathode) et GDL [CEA06] . . . . . 21
1.13 Image SEM des couches fibres carbones de la GDL [Mac09] . . . . . . . . 23
1.14 Profils des plaques bipolaires commerciales [Men00] . . . . . . . . . . . . 24
1.15 Schéma de la localisation des joints dans une cellule de la PAC . . . . . . . 25
1.16 Schéma d’une cellule unitaire représente des interfaces multi-contacts . . . 27
1.17 Profile de la plaque de la GDL [TR05] . . . . . . . . . . . . . . . . . . . . 29
1.18 Relation "perméabilité-porosité" . . . . . . . . . . . . . . . . . . . . . . . 30
1.19 Schéma des phénomènes mécaniques existants dans une pile de PEM . . . 34
+
1.20 Fissure sur GDL [SLC 08] . . . . . . . . . . . . . . . . . . . . . . . . . . 36
1.21 ’Pinhole’ sur MEA [SH04] . . . . . . . . . . . . . . . . . . . . . . . . . . 36

2.1 Exemples d’éléments finis surfaciques et volumiques. . . . . . . . . . . . . 38

2.2 Interpolation linéaire et quadratique . . . . . . . . . . . . . . . . . . . . . 39
2.3 La transformation iso-paramétrique d’un élément cubique. . . . . . . . . . 40
2.4 Approximation linéaire d’une loi de comportement . . . . . . . . . . . . . 41
2.5 Repère local de contact . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
2.6 Graphe des conditions de Signorini. . . . . . . . . . . . . . . . . . . . . . 48
2.7 Loi de frottement de Coulomb . . . . . . . . . . . . . . . . . . . . . . . . 50
2.8 Cône de Coulomb . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
136

2.9 Régularisation des lois de contact . . . . . . . . . . . . . . . . . . . . . . . 54

2.10 Projection sur le disque de Coulomb . . . . . . . . . . . . . . . . . . . . . 58
2.11 Projection sur le cône de Coulomb . . . . . . . . . . . . . . . . . . . . . . 59
2.12 Un barreau avec des conditions aux limites . . . . . . . . . . . . . . . . . . 60
2.13 Maillage initial et déformé . . . . . . . . . . . . . . . . . . . . . . . . . . 61
2.14 Distribution de la contrainte de cisaillement . . . . . . . . . . . . . . . . . 61
2.15 Pression normale et pression tangentielle par ANSYS et FER . . . . . . . . 62
2.16 Pression normale par ANSYS et FER . . . . . . . . . . . . . . . . . . . . 62
2.17 Déplacement tangentiel sur la zone AD . . . . . . . . . . . . . . . . . . . 63
2.18 Déplacement normal sur la zone AD . . . . . . . . . . . . . . . . . . . . . 63

3.1 Modèle symétrique 2D de la PAC . . . . . . . . . . . . . . . . . . . . . . 67

3.2 Modèle paramétrique de la PAC . . . . . . . . . . . . . . . . . . . . . . . 68
3.3 Maillage du modèle éléments finis . . . . . . . . . . . . . . . . . . . . . . 70
3.4 Composantes d’une zone de contact . . . . . . . . . . . . . . . . . . . . . 70
3.5 Erreur de détection du contact . . . . . . . . . . . . . . . . . . . . . . . . 71
3.6 La zone de contact entre les plaques GDL et BPP . . . . . . . . . . . . . . 72
3.7 Conditions aux limites et la charge . . . . . . . . . . . . . . . . . . . . . . 73
3.8 Déformation totale du modèle . . . . . . . . . . . . . . . . . . . . . . . . 74
3.9 Déformation de la GDL . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
3.10 Une séparation de contact . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
3.11 déplacement normal . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
3.12 Pression nulle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
3.13 Distribution de la pression . . . . . . . . . . . . . . . . . . . . . . . . . . 75
3.14 (a) Concentration de contrainte et (b) pression de contact . . . . . . . . . . 76
3.15 (a) Modèle du contact actuel et (b) modèle symétrique . . . . . . . . . . . 77
3.16 (a) Maillage et (b) Maillage fin . . . . . . . . . . . . . . . . . . . . . . . . 78
3.17 (a) Contraintes de trois plaques et (b) contraintes de la MEA . . . . . . . . 79
3.18 (a) Contraintes de trois plaques et (b) contraintes de la MEA avec maillage fin 79
3.19 (a) Pression à l’interface et (b) pression agrandie aux fonds des congés . . . 81
3.20 Déformation de la plaque GDL . . . . . . . . . . . . . . . . . . . . . . . . 82
3.21 Déformation de la plaque GDL avec pores de la porosité de 78% . . . . . . 84
3.22 Porosité de la GDL déformée et porosité initiale de la GDL . . . . . . . . . 85
3.23 Déformation de la plaque GDL avec pores de la porosité de 50% . . . . . . 86
3.24 Porosité de la GDL déformée et porosité initiale de la GDL . . . . . . . . . 86
3.25 Déformation de la plaque GDL avec pores aléatoires . . . . . . . . . . . . 87
3.26 Porosité de la GDL déformée et porosité initiale de la GDL . . . . . . . . . 87

4.1 Précompression = 0.01 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 90

 137

4.2 Précompression = 0.02 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 90

4.3 Précompression = 0.04 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 90
4.4 Précompression = 0.06 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 90
4.5 Précompression = 0.08 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 91
4.6 Précompression = 0.1 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 91
4.7 Précompression = 0.12 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 91
4.8 Précompression = 0.14 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 91
4.9 Variation de la pression avec la précompression de la PAC . . . . . . . . . 92
4.10 Variation de l’aire de contact avec la précompression de la PAC . . . . . . . 93
4.11 Précompression = 0.01 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 94
4.12 Précompression = 0.02 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 94
4.13 Précompression = 0.04 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 94
4.14 Précompression = 0.06 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 94
4.15 Précompression = 0.08 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 95
4.16 Précompression = 0.1 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 95
4.17 Précompression = 0.12 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 95
4.18 Précompression = 0.14 mm . . . . . . . . . . . . . . . . . . . . . . . . . . 95
4.19 Variation de la déformation et de la porosité avec la précompression . . . . 96
4.20 Rayon RRIB = 0.01 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
4.21 RRIB = 0.03 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
4.22 RRIB = 0.05 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
4.23 RRIB = 0.07 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
4.24 RRIB = 0.09 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
4.25 RRIB = 0.12 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
4.26 RRIB = 0.16 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
4.27 RRIB = 0.20 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
4.28 RRIB = 0.25 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99
4.29 RRIB = 0.30 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99
4.30 RRIB = 0.35 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99
4.31 RRIB = 0.40 mm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99
4.32 Variation de la pression de contact avec le rayon du congé . . . . . . . . . . 100
4.33 Variation de l’aire de contact avec le rayon du congé . . . . . . . . . . . . 101
4.34 Variation de la déformation et la porosité avec le rayon du congé . . . . . . 102
4.35 Variation de l’aire de contact avec la largeur des dents de la BPP . . . . . . 103
4.36 Différence de la pression de contact au poin 1 et 3 avec différentes largeurs 104
4.37 Porosité au point 1 pour différentes largeurs . . . . . . . . . . . . . . . . . 105
4.38 Porosité sur l’axé de symétrie avec différentes largeurs des dents . . . . . . 105
4.39 Variation de la porosité au point 1 pour différentes épaisseurs . . . . . . . . 107
138

4.40 Variation de la porosité sur l’axe de symétrie pour différentes épaisseurs . . 108
4.41 Uniformité de pression de contact pour différentes épaisseurs des MEAs . . 109
4.42 Porosité au point 1 pour différentes épaisseurs des membranes . . . . . . . 109
4.43 Porosité sur l’axe de symétrie pour différentes épaisseurs des membranes . 110
4.44 Contrainte dans la MEA (0.05 mm) . . . . . . . . . . . . . . . . . . . . . 111
4.45 Contrainte dans la MEA (0.183 mm) . . . . . . . . . . . . . . . . . . . . . 111
4.46 La tendance entre la résistance de contact et la rugosité différente . . . . . . 113
4.47 Modélisation des rugosités . . . . . . . . . . . . . . . . . . . . . . . . . . 114
4.48 Maillage de l’interface rugueuse . . . . . . . . . . . . . . . . . . . . . . . 114
4.49 Influence des aspérités sur l’aire de contact . . . . . . . . . . . . . . . . . 115
4.50 Influence de la densité des aspérités . . . . . . . . . . . . . . . . . . . . . 116
4.51 Modèle avec différentes hauteurs de dent . . . . . . . . . . . . . . . . . . . 117
4.52 Pression avec différentes dents . . . . . . . . . . . . . . . . . . . . . . . . 117
4.53 Déformation de la plaque GDL avec différentes hauteurs des dents . . . . . 118
4.54 Numéro des dents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118
 139

Liste des tableaux

1.1 Types des piles à combustibles . . . . . . . . . . . . . . . . . . . . . . . . 12

1.2 Caractéristiques des MEAs de la série de Nafion . . . . . . . . . . . . . . . 22
1.3 Tableau des caractéristiques de GDL de la série de TGP-H . . . . . . . . . 23
1.4 Caractéristiques des BPPs de la SGL . . . . . . . . . . . . . . . . . . . . . 25
1.5 Paramètres du modèle de Tomadakis . . . . . . . . . . . . . . . . . . . . . 30

2.1 Temps de CPU . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64

3.1 Propriétés mécaniques de la plaque BPP, la plaque GDL et la plaque MEA . 66

3.2 Propriétés et géométries du modèle de contact paramétrique de la PAC . . . 69
3.3 Aire de contact à l’interface entre la plaque GDL et la plaque BPP . . . . . 81

4.1 Aspérités des plaques BPP et GDL . . . . . . . . . . . . . . . . . . . . . . 113

4.2 Variation des états de contact sur les dents . . . . . . . . . . . . . . . . . . 118