TH' ESE Docteur L' Ecole Nationale Sup Erieure D'arts Et M Etiers

N˚ : 2007 ENAM 0032
École doctorale n˚432 : Science des Métiers de l’Ingénieur
THÈSE
pour obtenir le grade de
Docteur
de
l’École Nationale Supérieure d’Arts et Métiers

Spécialité : MÉCANIQUE
présentée et soutenue publiquement par
Munif JOURIEH
Développement d’un modèle représentatif d’une

éolienne afin d’étudier l’implantation de plusieurs
machines sur un parc éolien.
à soutenir le 20 décembre 2007 à L’ENSAM, devant le Jury composé de :
F. BAKIR Professeur, LEMFI, ENSAM, Paris Président

G. DESCOMBES Professeur, CNAM, Paris Rapporteur
I. AMJEH Professeur, Université de Damas Rapporteur
F. MASSOUH Maı̂tre de Conférences HDR, LMF, ENSAM, Paris Examinateur
P. KUSZLA Professeur agrégé, SINUMEF, ENSAM, Paris Examinateur
G. NOTTON Maı̂tre de Conférences HDR, Université de Corse, Ajaccio Examinateur
M. RAPIN Ingénieur de recherche, ONERA, Chatillon Invité
Laboratoire de Simulation Numérique en Mécaniques des Fluides
ENSAM CER de PARIS
L’ENSAM est un Grand Établissement dépendant du Ministère de l’Éducation Nationale, composé de huit centres :
AIX-EN-PROVENCE ANGERS BORDEAUX CHÂLONS-EN-CHAMPAGNE CLUNY LILLE METZ PARIS
Table des matières
1 Introduction générale 7
1.1 Evolution de l’énergie éolienne . . . . . . . . . . . . . . . . . . . . . . . . . 8
1.2 Sillage éolien et interaction entre machines . . . . . . . . . . . . . . . . . . 12
1.3 Prédiction du comportement d’un rotor aérodynamique : état de l’art . . . 14
1.4 Modèles hybrides . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
1.5 Conclusion : . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
2 Exploration en soufflerie du champ de vitesse autour d’une éolienne 21

2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.2 Description des moyens d’essais . . . . . . . . . . . . . . . . . . . . . . . . 22
2.2.1 Soufflerie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.2.2 Éolienne testée . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
2.3 Mesures PIV (Particle Image Velocimetry) . . . . . . . . . . . . . . . . . . 25
2.3.1 Mise en place des mesures PIV . . . . . . . . . . . . . . . . . . . . 27
2.3.2 Résultats . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
2.3.3 Champs de vitesse . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.3.4 Analyse des résultats . . . . . . . . . . . . . . . . . . . . . . . . . . 32
2.4 Mesures par fil chaud . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
2.4.1 Principe de l’anémométrie à fil chaud . . . . . . . . . . . . . . . . . 35
2.4.2 Méthodologie suivie dans l’exploration au fil chaud . . . . . . . . . 37
2.4.3 Plan d’exploration . . . . . . . . . . . . . . . . . . . . . . . . . . . 37
2.4.4 Représentation du sillage . . . . . . . . . . . . . . . . . . . . . . . . 39
2.4.5 Étude de trajectoire des tourbillons marginaux . . . . . . . . . . . . 40
2.4.6 Effet de l’éolienne sur l’écoulement à l’amont du rotor . . . . . . . . 44
2.4.7 Champs de vitesse tangentielle . . . . . . . . . . . . . . . . . . . . . 46
2.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
2
3 Développement d’un modèle hybride basé sur le concept du disque actif 51
3.1 Introduction : . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
3.2 Méthode de l’élément de pale . . . . . . . . . . . . . . . . . . . . . . . . . 52
3.2.1 Théorie de l’élément de pale . . . . . . . . . . . . . . . . . . . . . . 53
3.2.2 Bilan de quantité de mouvement appliqué à un élément de pale . . 55
3.3 Application de la méthode de l’élément de pale . . . . . . . . . . . . . . . . 57
3.3.1 Algorithme de calcul . . . . . . . . . . . . . . . . . . . . . . . . . . 58
3.3.2 Discrétisation de la pale et critère de convergence . . . . . . . . . . 58
3.4 Modèle de disque actif . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
3.4.1 Equation de continuité . . . . . . . . . . . . . . . . . . . . . . . . . 62
3.4.2 Bilan de quantité de mouvement . . . . . . . . . . . . . . . . . . . 63
3.4.3 Coefficient de puissance . . . . . . . . . . . . . . . . . . . . . . . . 63
3.4.4 Limite de Betz . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
3.5 Description du modèle hybride . . . . . . . . . . . . . . . . . . . . . . . . . 64
3.6 Etude de la discrétisation du disque actif . . . . . . . . . . . . . . . . . . . 66
3.7 Interfaçage avec le logiciel Fluent . . . . . . . . . . . . . . . . . . . . . . . 68
3.8 Répartition volumique des forces appliquées sur le disque actif . . . . . . . 68
3.9 Algorithme de calcul . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
3.10 Application du modèle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
3.11 Résultats et analyse . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
3.11.1 Eoliennes comparées . . . . . . . . . . . . . . . . . . . . . . . . . . 73
3.11.2 Le domaine de calcul . . . . . . . . . . . . . . . . . . . . . . . . . . 75
3.11.3 Comparaison des résultats numériques et expérimentaux . . . . . . 77
3.12 Comparaison entre le sillage mesuré d’une éolienne et le sillage calculé par
le modèle hybride . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78
3.13 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80
4 Développement d’un modèle hybride basé sur le concept de la ligne

active 83
4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84
4.2 Développement d’un modèle de cylindre actif . . . . . . . . . . . . . . . . . 84
4.3 Répartition des efforts dans un cylindre actif . . . . . . . . . . . . . . . . . 87
4.3.1 Choix du diamètre du cylindre actif . . . . . . . . . . . . . . . . . . 88
4.4 Application du modèle de cylindre actif aux cas NREL s809 phase II et VI 96
4.4.1 Algorithme de calcul . . . . . . . . . . . . . . . . . . . . . . . . . . 96
4.4.2 Surface de référence . . . . . . . . . . . . . . . . . . . . . . . . . . . 99
3
4.4.3 Domaine de calcul . . . . . . . . . . . . . . . . . . . . . . . . . . . 100
4.5 Résultats et analyse . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101
4.6 Comparaison simulation/ expérience du sillage d’une éolienne Rutland 503
modifiée . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103
4.7 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 105
5 Interactions entre des éoliennes 108

5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109
5.2 Simulation complète de l’interaction entre deux éoliennes à l’aide d’un
modèle hybride . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110
5.3 Effet de l’éloignement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112
5.4 Effet de la vitesse à l’infini amont . . . . . . . . . . . . . . . . . . . . . . . 118
5.4.1 Effet de la vitesse en amont avec un modèle hybride CA . . . . . . 118
5.4.2 Effet de la vitesse en amont avec un modèle hybride DA . . . . . . 121
5.5 Interaction entre six machines en cascade . . . . . . . . . . . . . . . . . . . 122
5.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123
6 Conclusion générale 125
7 Annexe 130
4
Liste des symboles
V1 : Vitesse du vent à l’infini amont. [m/s]

Vd : Vitesse du vent dans le plan du rotor . [m/s]
V2 : Vitesse du vent à l’infini aval. [m/s]
Vl : vitesse locale. [m/s]
W : Vitesse relative du vent. [m/s]
D : Diamètre du rotor. [m]
d : Diamètre du cylindre actif (CA). [m]
R : Rayon de l’extrémité de la pale. [m]
r : Rayon local. [m]
B : Nombre des pales. [−]
Cx : Coefficient de traı̂née. [−]
Cz : Coefficient de portance. [−]
Cn : Coefficient de la force normale au plan de rotation. [−]
Ct : Coefficient de la force tangentielle au plan de rotation. [−]
Cp : Coefficient de puissance. [−]
β : Angle de calage. [rad]
α : Angle d’incidence. [rad]
φ : Angle d’écoulement relatif. [rad]
θ : Angle azimutal. [rad]
ρ : Masse volumique. [kg/m3 ]
ν : Viscosité cinématique. [m2 /s]
Ω : Vitesse angulaire du rotor. [rad/s]
N : Vitesse de rotation du rotor. [tr/min]
c : corde de pale.[m]
a : Facteur d’induction axial. [−]
a′ : Facteur d’induction tangentiel. [−]
P + : Pression en amont du disque actif. [P a]
P − : Pression en aval du disque actif. [P a]
P1 : Pression atmosphérique. [P a]
A1 : Surface du tube de courant à l’infini amont. [m2 ]
Ad : Surface du disque actif. [m2 ]
A2 : Surface du tube de courant à l’infini aval. [m2 ]
σ : Solidité de la pale. [−]
λ : Rapidité spécifique. [−]
5
P : Puissance extraite par une éolienne. [W ]
6
Chapitre 1
Introduction générale
7
1.1 Evolution de l’énergie éolienne
L’énergie éolienne est probablement une des plus anciennes sources d’énergie. Cette
énergie propre et renouvelable existe depuis toujours, mais jusqu’à présent son exploitation
reste difficile. L’utilisation de l’énergie éolienne a commencé en 1700 avant Jésus-Christ
environ, Hammourabi, roi de Babylone, actuellement l’Irak, a pensé bénéficier de cette
énergie propre pour l’irrigation. Il a utilisé la puissance du vent pour le pompage de l’eau
avec des éoliennes à axe vertical [1, 2].
Au neuvième siècle avant Jésus-Christ, la Perse, royaume qui correspond maintenant
à l’est de l’Iran et à l’Afghanistan, a utilisé l’énergie du vent dans le but de moudre le
grain et de pomper l’eau. La première conception documentée connue était celle d’un
moulin à vent persan, le panémone, utilisant des voiles verticales faites de roseaux ou de
bois qui étaient attachées à un axe central vertical. La figure 1.1 montre le principe de
fonctionnement d’un tel moulin [3, 4].
Vent
Voile
Fig. 1.1 – Principe du premier moulin à vent persan à axe vertical.
En Inde, l’exploitation de cette énergie renouvelable a pris son essor au quatrième

siècle avant Jésus-Christ [5].
8
Ensuite, trois siècles avant Jésus-Christ, les Egyptiens ont commencé à bénéficier de
cette énergie propre. L’inventeur égyptien Héron d’Alexandrie utilisa l’énergie éolienne
grâce à un moulin à vent à axe horizontal. Ce dispositif permettait, contrairement à
ceux utilisés pour les tâches agricoles, d’utiliser l’énergie éolienne pour alimenter en air
comprimé un orgue [1].
Au Moyen-Age, en Europe, l’exploitation de cette énergie a commencé par l’apparition
des moulins de vent, principalement en France, en Italie, en Espagne et au Portugal. Les
premiers travaux écrits ont ainsi vu le jour en Normandie en 1180. Certains auteurs
pensent que la technologie du moulin à vent à axe vertical développée en Perse s’étendit
en Europe suite aux Croisades au Moyen-Orient en évoluant sous la forme de moulins
à vent à axe horizontal. On les rencontre, un peu plus tard, en Grande-Bretagne, en
Hollande et en Allemagne sous le type de machine à axe horizontal comportant quatre
ailes placées en croix. Elles servaient principalement à moudre du grain [6].
Les Hollandais participèrent activement au développement des moulins en Europe,
grâce aux nombreuses améliorations dans la conception et à l’invention de différents types
de moulins. En effet, pendant le treizième siècle, la Hollande utilisa les moulins à vent
pour pomper l’eau et ainsi assécher les polders. Accouplées à des roues à godet ou à des
vis d’Archimède, ces machines pouvaient élever l’eau jusqu’à cinq mètres. L’utilisation
des moulins à vent a connu un grand succès jusqu’au dix neuvième siècle. Ils ont fourni à
l’homme l’énergie mécanique nécessaire pour certains travaux agricoles.
Mais avec la révolution industrielle et l’apparition de la machine à vapeur, du moteur à
combustion et plus tard le développement de l’électricité, le développement des éoliennes et
l’exploitation des moulins à vent sont délaissés. Leurs utilisations ont décliné jusqu’après
la deuxième guerre mondiale. Face à la présence de nouveaux moyens de production
d’énergie, les aérogénérateurs n’arrivent pas à s’imposer [2].
Pendant le vingtième siècle, plusieurs projets d’éoliennes voient le jour. L’éolienne
lente multipale est développée en Amérique par la Rural Electrification Administration
[7]. L’éolienne rapide, inventée en France par l’Académicien français Darrieus, entraı̂nait
des générateurs électriques [1]. En 1950, Johannes Juul [6] développa un modèle éolien
avec trois pales, utilisant des dispositifs de réglages aérodynamiques de la puissance dans
le cas de décrochage et du contrôle de dérapage.
Cependant dans les années 1970 et après la crise du pétrole causée par la guerre au
Moyen-Orient, l’administration du président Carter aux USA et d’autres dans plusieurs
pays comme la Suède, le Canada et la Grande-Bretagne, ont pris la décision de développer
la recherche dans le domaine de l’énergie renouvelable [8, 9]. En effet, avec la diminution du
stock mondial d’hydrocarbures, la demande énergétique sans cesse croissante, et la crainte
9
d’une pollution de plus en plus envahissante, l’énergie propre et renouvelable attire les
puissances mondiales. L’énergie éolienne revient au premier plan de l’actualité car son
exploitation peut s’avérer très rentable dans les nombreuses régions ventées du globe.
Ainsi, au début des années 1980, dans des pays tels que le Danemark ou l’Allemagne,
où il n’existe pas des ressources d’énergie aussi importantes que l’énergie éolienne, le
développement d’une industrie éolienne performante a été privilégié [10, 11, 7].
Les éoliennes ont ainsi continué à évoluer au cours des 20 dernières années, et le coût
global de l’énergie nécessaire à la production d’électricité à partir du vent est maintenant
concurrentiel avec les sources d’énergie traditionnelles comme les combustibles fossiles
[12, 13]. Cette réduction du coût de l’électricité est le résultat de progrès importants de
la technologie utilisée par cette industrie (amélioration des conceptions aérodynamiques,
amélioration des matériaux utilisés, . . .) [14, 15].
Actuellement, l’énergie éolienne est bien implantée parmi les autres sources d’énergie
avec une croissance très forte [16, 17]. La figure 1.2 présente la puissance issue de l’énergie
éolienne des dix premiers marchés éoliens du monde en 2005 [18]. A partir de celle-ci, on
observe que cette énergie est bien développée aux USA, Allemagne, Espagne et en Inde
par rapport aux autres pays.
Fig. 1.2 – Production d’énergie éolienne dans les dix premiers marchés éoliens en 2005.
10
Fig. 1.3 – Puissance et croissance du secteur de l’énergie en France (MW, MW/an).
La croissance de l’exploitation de l’énergie éolienne dans le monde est influencée par un

facteur très important, c’est la pollution, car les autres types d’énergie (Pétrole, Nucléaire)
sont forcément polluants, ce qui amène à prendre des initiatives dans le monde pour
ajouter des conditions d’utilisation afin de diminuer cette pollution (protocole de Kyoto
[19, 20]). Ce qui donne des avantages supplémentaires à l’utilisation de l’énergie éolienne.
La figure 1.3 représente l’évolution de l’énergie éolienne en France de 1996 à 2005.
On remarque une forte croissance de l’exploitation de l’énergie éolienne au cours des cinq
dernières années.
Étant donnée cette forte croissance et le fait que les zones où l’énergie éolienne est
techniquement exploitable sont limités, il est important de bien gérer les parcs éoliens.
L’installation d’un parc éolien doit donc être optimisée de manière à tirer parti au mieux
de la ressource en vent. Cette recherche de l’installation optimale fait l’objet de cette thèse
où nous allons tenter de mettre en place des modèles de rotors d’éoliennes qui permettront
de représenter un parc de machines et d’en évaluer les performances.
11
1.2 Sillage éolien et interaction entre machines
Dans les parcs éoliens, il est nécessaire d’espacer les éoliennes afin d’éviter que le sillage
et le déficit de vitesse existant derrière chaque machine n’affectent trop la production
énergétique et l’intégrité mécanique des éoliennes situées plus en aval [21, 22]. En règle
générale, les éoliennes dans les parcs sont toutes identiques et la distance entre les machines
est définie d’après des règles simples données par les constructeurs, à savoir : cinq à neuf
fois le diamètre dans la direction des vents dominants et de trois à cinq fois le diamètre dans
la direction perpendiculaire [23]. Ces règles ne constituent évidemment pas une garantie
d’optimisation de l’implantation des machines et ne permettent pas d’évaluer la prise en
compte de la topographie du terrain. Actuellement, la simulation numérique en mécanique
des fluides a fait d’énormes progrès et permet de donner des informations intéressantes
sur l’écoulement d’air autour des éoliennes. Les codes industriels modernes permettent
d’obtenir des résultats quantitatifs satisfaisants [24]. La simulation pourrait donc être
un outil appréciable pour analyser le sillage éolien [25] et optimiser l’implantation des
éoliennes dans un parc en fonction des caractéristiques des machines et de la topographie
du site éolien [26, 27, 28].
Cependant, la simulation de l’écoulement à travers le rotor éolien et dans son sillage
est toujours une opération lourde nécessitant une préparation et un temps de calcul im-
portant. Les moyens de calcul modernes permettent de simuler finement le comportement
aérodynamique d’une pale d’éolienne et éventuellement d’un rotor complet [29]. Il est
toutefois hors de question, dans l’état actuel du matériel informatique disponible, de si-
muler l’ensemble des machines composant un parc éolien. D’un point de vue expérimental,
la non uniformité du champ de vitesse sur un site éolien, de même que l’encombrement
des machines rend pratiquement impossible l’étude en soufflerie de l’ensemble d’un parc
éolien. Ainsi, la voie expérimentale se limite à des études partielles des parcs éoliens
avec un nombre très réduit de machines. Il s’agit donc dans le travail réalisé dans cette
thèse de proposer un modèle capable de dépasser les limitations imposées tant par le
dispositif expérimental que par les outils de simulation. Cette thèse partira de l’idée qui
consiste à rechercher un modèle équivalent du rotor éolien permettant de bien représenter
le sillage et son développement tout en évitant de calculer en détail l’écoulement autour
des différentes pales. Le temps de calcul ainsi épargné, grâce à l’usage d’un tel modèle,
doit permettre de simuler et d’optimiser l’implantation des machines d’un parc éolien
dans un environnement donné.
L’approche habituellement utilisée pour modéliser les rotors éoliens et les rotors d’hélicoptères
aussi, est basée sur un disque actif [30, 31, 32]. Le rotor complet est ainsi remplacé par un
12
disque infiniment mince à nombre infini de pales et qui génère un écoulement équivalent
à l’écoulement réel moyenné. Concrètement, c’est une condition aux limites particulière
(discontinuité de pression) appliquée sur toute la surface balayée par les pales, condition
qui produit la même variation de quantité de mouvement et d’énergie que le rotor en
fonctionnement.
L’objectif recherché consiste à trouver un meilleur modèle que le disque actif et qui
rendra possible la prise en compte du nombre fini des pales du rotor. Ainsi, il est question
de développer une modélisation capable de reproduire chaque pale individuellement.
13
1.3 Prédiction du comportement d’un rotor aérodynamique :
état de l’art
Nous allons ici faire un bref tour d’horizon des méthodes utilisées dans la simulation des
écoulements autour des rotors éoliens. La plupart des méthodes appliquées pour prévoir
le comportement d’une éolienne dans le cas axisymétrique, utilisent la théorie du disque
actif, où le rotor est remplacé par des forces distribuées sur un disque d’épaisseur nulle et
perméable [33, 4]. Cette théorie a été introduite par Froude [34] et Rankine [35], elle est
développée ensuite par Lanchester en 1915 [36] puis par Betz en 1920 [5]. Ils ont démontré
16
que la puissance maximale extraite par une éolienne est égale à 27 de l’énergie cinétique
de l’écoulement incident, soit environ 59.3%.
En 1935, la théorie de l’élément de pale est proposée par Glauert [37], cette méthode
est basée sur la division de l’écoulement en volumes de contrôles annulaires, auxquels
on applique le bilan de quantité de mouvement et d’énergie. Ces anneaux s’étendent de
l’infini amont jusqu’à l’infini aval par rapport au rotor, figure 1.4.
Disque actif
V1 R
X
dr
Fig. 1.4 – Tube de courant au rayon local r à travers une éolienne tripales remplacée par
un disque actif.
Du fait de sa simplicité, cette méthode est largement utilisée dans le domaine industriel
[38]. L’application de la théorie de l’élément de pale sur le disque actif donne la base la
plus utilisable pour construire des modèles hybrides [39, 40, 41].
14
La limitation principale de la méthode du disque actif est qu’elle n’est valide que pour
les cas des éoliennes axisymétriques, et que cette méthode distribue les forces de façon
égale dans la direction azimutale du disque actif et ne représente pas individuellement
chaque pale. Ainsi, l’influence des pales est prise comme une quantité intégrée dans la
direction azimutale.
Pour surmonter ce problème et prendre en compte la présence des pales, Sørensen et
Shen [42], ont introduit la technique de la ligne active, où les forces créées par une pale
sont distribuées le long d’une ligne qui représente cette pale, figure 1.5. Dans le concept
de la ligne active, chaque pale est représentée par une ligne discrétisée par n points. Les
forces sont distribuées autour de chaque point dans un plan perpendiculaire à la ligne
active selon une distribution spécifique. Par exemple, une distribution Gaussienne dans le
travail de Mikkelesen [43], ou une distribution uniforme dans le modèle de cylindre actif
présenté dans cette thèse au chapitre quatre.
Ligne active
x
Fig. 1.5 – Concept de ligne active
La première formulation du concept de ligne active a été proposée par Sørensen et Shen
[42], leur analyse montre l’accord avec la courbe de puissance mesurée pour l’éolienne
15
Nordtank tri pale 500kW. Plus tard, Mikkelsen [43] a proposé un modèle hybride en
combinant cette méthode avec le code EllipSys-3D. Les résultats obtenus sont comparés
avec les résultats expérimentaux de la turbine de Tjaereborg 1 pour plusieurs cas de
dérapage, et donnent des valeurs assez proches et comparables [39].
Ivanell [44] a utilisé la même méthode que Mikkelsen et Sørensen, et a montré que le
contrôle de la distribution Gaussienne des forces autour de la ligne active est influencé
par quelques paramètres. Il a mis en évidence l’influence du nombre de Reynolds et de
la résolution du maillage. Ivanell a aussi observé l’amélioration des résultats en utilisant
la correction de Prandtl afin de mieux prendre en compte les phénomènes d’extrémité de
pale.
D’autre part, Massouh et Dobrev [45] ont proposé un nouveau modèle hybride, dans
lequel chaque pale est remplacée par une surface de discontinuité de pression. De cette
manière, on surmonte les difficultés de la distribution des forces autour de la ligne active.
Cette distribution de pression varie le long de la corde. La représentation de la pale et des
conditions initiales pour le développement du sillage est ainsi améliorée par rapport aux
méthodes de disque actif ou de ligne active. La validation de ce modèle hybride est faite
avec le cas de l’éolienne NREL2 s8093 phase IV et aboutit à des résultats satisfaisants.
On observe que les méthodes basées sur le concept du disque actif ou sur le concept de
la ligne active et la théorie de l’élément de pale BEM (Blade Element Momentum), ont
bien validé leur capacité à représenter différents types de rotors(pales droites ou vrillées)
dans différentes conditions de fonctionnement (dérapage par exemple) [46, 47, 48].
Toutes les méthodes trouvées dans la littérature sont des représentations plus ou moins
approchées des rotors réels. De ce fait on observe une dispersion importante des résultats
fournis par chacun des auteurs. Ceci est particulièrement visible avec la comparaison
simulation-expérience établie par le NREL en Décembre 2000 [49]. Dans cette étude, 19
méthodes de simulation numérique d’une éolienne NREL phase VI sont comparées. Elles
utilisent différents modèles, de la méthode la plus simple, la méthode de l’élément de
pale (BEM) jusqu’aux simulations tridimensionnelles turbulentes complètes. La figure 1.6,
représente les modèles participants à la comparaison pour une vitesse de rotation faible
et sans dérapage, la ligne noire avec les points représente les mesures expérimentales pour
1
La turbine à axe horizontal de Tjaereborg, d’un diamètre de 61 m, a trois pales. Le profil des pales
est un NACA 4412-43. La longueur de la corde vaut 0.9 m à l’extrémité de pale et augmente jusqu’à 3.3
m au rayon extérieur de 6 m. L’axe de rotation est situé à 60 m du sol. Le vrillage de pale est de 1° tous
les 3 m.
2
NREL : National Renewable Energy Laboratory.
3
Pour plus des information autour le profil s809, voir l’annexe A.
16
l’éolienne NREL phase VI dans la soufflerie Ames de la NASA4 [50]. On voit nettement
sur cette figure qu’il y a plusieurs méthodes qui ont donné des résultats assez loin des
données expérimentales. Même d’un simple point de vue qualitatif, certaines approches
semblent mises en echec. Ces résultats mettent en évidence la difficulté de produire des
simulation prédictives.
Fig. 1.6 – Présentation des méthodes de prédiction de couple de l’éolienne NREL phase
VI en fonction de la vitesse du vent.
1.4 Modèles hybrides

Notre contribution à la simulation des écoulements autour des éoliennes consistera à
développer, dans la lignée des méthodes existantes un nouveau modèle hybride. C’est à
dire que nous tenterons de combiner des modèles de rotor (éléments de pale, disque actif,
cylindre actif) avec une simulation numérique basée sur la résolution des équations de
Navier-Stokes tridimensionnelles turbulentes.
Dans ce travail nous présentons deux types de modèles hybrides, le premier est basé
sur le concept du Disque Actif (DA), et le deuxième est basé sur le concept de la ligne
active et sera un modèle de Cylindre Actif (CA).
4
NASA : National Aeronautics and Spase Administration.
17
– Avec le modèle hybride DA, nous représentons le rotor par un disque de diamètre
égal à celui du rotor et d’épaisseur liée à la corde de la pale5 . La distribution des
forces exercées par les pales est faite de manière uniforme sur le disque. La présence
des pales est donc moyennée. Ce type de modèle hybride n’est utilisable que dans
les cas axisymétriques et ne peut pas calculer un rotor en dérapage.
– Pour surmonter cette limitation, nous utilisons ensuite un modèle hybride dit modèle
de cylindre actif (CA) qui prend en compte le nombre de pales. Ce modèle est basé
sur le concept de la ligne active, chaque pale est remplacée par un cylindre de
longueur égale à la longueur de la pale et de diamètre lié à la corde de la pale.
La validation de ces modèles hybrides est faite de deux manières différentes, par rap-
port au sillage éolien et par rapport à l’évolution de la puissance en fonction de la vitesse
du vent.
Pour le sillage éolien, nous avons réalisé des travaux expérimentaux pour une éolienne
tri pales (Rutland 503 modifiée) dans la soufflerie du laboratoire de mécanique des fluides
(LMF) à l’École Nationale Supérieure d’Arts et Métiers (ENSAM) à Paris. Ce travail
présente l’exploration du sillage éolien dans la veine d’essai de la soufflerie, celui-ci est
exploré par deux méthodes différentes. La première méthode utilisée est l’anémométrie
à fil chaud à deux composantes, cela donne l’exploration du champ des vitesses axiale
et tangentielle du sillage. La deuxième méthode est la méthode PIV (Particle Image
Velocity). Avec celle-ci, on explore le champ des vitesses axiales et radiales. Pour avoir
une présentation assez détaillée sur le sillage derrière l’éolienne, on a fait des explorations
à plusieurs angles azimutaux.
Pour l’évolution de puissance en fonction de la vitesse du vent, la validation de ces
modèles est faite avec une éolienne NREL phase II [51] et NREL phase VI [52].
Enfin, nous avons étudié l’interaction entre des machines en utilisant ces modèles
hybrides.
1.5 Conclusion :
Ce chapitre donne une idée générale sur de l’énergie éolienne et son intérêt. Cette
énergie renouvelable a suivi son chemin depuis plusieurs années avec une croissance an-
nuelle d’utilisation très importante dans le monde. Au début du vingtième siècle, plusieurs
hypothèses sont apparues dans le domaine de l’aérodynamique concernant les éoliennes.
Les théories du disque actif et la ligne active sont les théories les plus utilisées pour étudier
5
Par abus de langage, nous retiendrons ici l’expression ”disque actif” bien que ce disque possède une
épaisseur non nulle et soit en fait un cylindre de hauteur très petite devant son rayon.
18
la performance d’une éolienne avec les équations de Navier-Stokes. La progression dans
le domaine de la simulation numérique permet de représenter une géométrie reélle d’une
éolienne par un modèle hybride beaucoup plus simple et facile à réaliser.
19
20
Chapitre 2
Exploration en soufflerie du champ

de vitesse autour d’une éolienne
21
2.1 Introduction
Cette partie présente des explorations de l’écoulement en amont et en aval d’une
éolienne tri pales à axe horizontal. Ces travaux permettent d’avoir les données expérimentales
nécessaires pour la validation de nos modèles hybrides, validation qui sera présentée dans
les chapitres suivants.
Dans ce travail, nous avons utilisé deux méthodes différentes, la première est la méthode
de vélocimétrie par image de particules PIV (Particle Image Velocimetry). Par cette
méthode, nous avons exploré le sillage proche derrière la machine testée. Les mesures
ont été synchronisées avec la position d’une pale de référence. Pour avoir une idée claire
de la structure du sillage en trois dimensions, d’autres mesures ont été réalisées avec un
retard contrôlé afin d’explorer plusieurs plans azimutaux derrières la pale, l’angle azimutal
entre deux plans successif est égal à 30°. Avec cette méthode, nous n’avons pas pu explo-
rer les deux champs de vitesse axiale et tangentielle simultanément. Pour atteindre cet
objectif, nous avons utilisé une autre méthode, la méthode de l’anémométrie à fil chaud.
Avec l’anémométrie à fil chaud, nous avons fait une exploration pour les deux champs de
vitesse (axiale et tangentielle). Cette exploration est réalisée en amont du plan de rotation
jusqu’à 0.4D1 , et en aval jusqu’à 1.5D environ du celui-ci.
Ces expériences sont réalisées dans le veine d’essai de la soufflerie de l’ENSAM-Paris.
A cause des dimensions limitées de notre veine d’essai, il est possible d’explorer le sillage
proche uniquement. Mais, l’intérêt ici est d’avoir des données expérimentales concernant le
sillage proche, et que ces données permettent de valider des modèles de rotors simplifiés,
couplés avec la simulation numérique. De plus, le sillage proche donne les conditions
initiales pour le développement du sillage lointain [53, 54].
2.2 Description des moyens d’essais

2.2.1 Soufflerie
La soufflerie du laboratoire de mécanique des fluides de l’ENSAM-Paris est de type
Prandtl avec une veine semi guidée sans parois latérales et ouverte à la pression at-
mosphérique. Elle comporte les parties suivantes : veine d’essai, diffuseur, deux séries
de coudes avant la veine de retour, ventilateur, et enfin deux séries de coudes avant la
chambre de tranquillisation puis le convergent d’alimentation de la veine, voir la figure
2.1.
1
D : Diamètre du rotor testé.
22
La veine d’essai a une section rectangulaire 1.35 m × 1.65 m avec une longueur de 2 m.
Cette partie est ouverte à la pression atmosphérique et la vitesse maximale de l’écoulement
est 40 m/s. La veine d’essai est dotée d’un plafond transparent en plexiglas qui facilite la
visualisation des écoulements. L’équipement de la veine d’essai comporte un explorateur
commandé par ordinateur et capable de se déplacer dans les trois directions de l’espace.
La veine est suivie par un diffuseur avec un angle d’ouverture voisin de 7 degrés pour
éviter les décollements de l’écoulement sur les parois.
Chambre de Explorateur. Veine d’essai

tranquillisation. 1.35m x 1.65m x 2m.
8.4 m.
7 m.
22.5m.
Ventilateur électrique Explorateur. Veine de retour
75kW. 3 pales. 3m x 3m x 6m.
Fig. 2.1 – Soufflerie de l’ENSAM-Paris.
La veine de retour possède une section carrée de 3 m × 3 m sur une longueur de 6 m

et la vitesse maximale d’écoulement est d’environ 10 m/s. Cette veine est intéressante par
ses dimensions, mais elle n’offre pas de bonnes caractéristiques d’uniformité du champ de
vitesses et pour ces raisons, l’exploration du champ de vitesse dans la veine d’essai est
préférable.
Le ventilateur comporte une hélice tripale de trois mètres de diamètre pouvant fonc-
tionner jusqu’à 750 tour/min. Il est entraı̂né par un moteur électrique de 75 kW . La
variation de vitesse de rotation du moteur permet d’assurer la régulation de la vitesse
d’écoulement. Après le ventilateur, l’air est orienté vers la chambre de tranquillisation qui
comporte une grille permettant d’aligner les filets fluides et de diminuer la turbulence. La
chambre de tranquillisation alimente la veine d’essai par l’intermédiaire d’un convergent
qui permet une répartition uniforme des vitesses dans sa section contractée. Le convergent
23
entre la chambre de tranquillisation et la veine d’essai possède un rapport de contraction
de 1 :12.5, ce qui permet de rendre l’écoulement globalement uniforme. Ainsi, le vent
atteint dans la veine d’essai une vitesse maximale de 40 m/s (144 km/h) avec une très
bonne qualité d’écoulement.
Enfin, il y a quatre coudes à 90° munis d’aubages qui imposent au fluide un écoulement
en plans parallèles. Les coudes ont pour effet de guider au mieux l’écoulement pour éviter
la formation de structures tourbillonnaires lors des changements de direction.
2.2.2 Éolienne testée

La maquette utilisée dans ce travail est une éolienne commerciale (Rutland 503 de la
Société Marlec) qui a été modifiée pour les essais, figure 2.2. En effet, l’éolienne d’origine
avait 6 pales, mais trois pales ont été ôtées afin d’obtenir une éolienne tripales semblable
aux éoliennes utilisées dans les fermes éoliennes. Le diamètre de cette éolienne est de 0.5
m avec un moyeu de 0.135 m de diamètre. Les pales ne sont pas vrillées et sont montées
avec un angle de calage constant de 10°. La corde au pied de la pale est 0.065 m et 0.045
m à l’extrémité. L’axe de rotation de cette éolienne est placé à une hauteur de 0.7 m au
milieu de la veine d’essai, à l’aide d’un mat de diamètre 0.037 m pour éviter l’influence
du plafond et du plancher sur le comportement du sillage et pour permettre aux lasers
fixés sur le plafond transparent de la veine d’essais d’éclairer le plan d’exploration avec
une intensité suffisante.
L’exploration du sillage de cette maquette est réalisée avec une vitesse axiale de 9.4
m/s et une vitesse de rotation de 110 rad/sec. Le contrôle de la vitesse de rotation a été
assuré par un rhéostat connecté à la sortie de la génératrice électrique équipant l’éolienne.
24
Fig. 2.2 – L’éolienne testée dans la veine d’essai.
2.3 Mesures PIV (Particle Image Velocimetry)

Depuis une vingtaine d’années, la technique PIV est utilisée pour mesurer la vitesse
d’un fluide et déterminer le champs de vorticité instantanée. Cette technique permet de
mesurer le champ instantané de vitesse au lieu de mesurer la vitesse en un point en
fonction du temps.
Le principe de cette méthode est basé sur l’ensemenssement de l’écoulement et la
mesure, en fonction du temps, de la vitesse des particules et sur la détermination de leurs
trajectoires dans un plan illuminé par une nappe laser [55]. L’illumination intensive par
le laser d’un plan de l’écoulement, permet d’avoir une image instantanée, et de détecter
les positions instantanées des particules d’ensemencement.
On enregistre les images des particules à deux instants successifs. La synchronisation
entre les tirs laser et la prise d’images permet d’obtenir un champ de vitesse instantanée.
Le principe de mesure est expliqué sur la figure 2.4 qui représente un plan rectangulaire
détecté par une caméra.
25
Fig. 2.3 – Éolienne testée avec une nappe Laser.
Pour simplifier, on suppose la présence de 6 particules qui passeront des positions

indiquées A,B,C,D,E et F dans la première image aux positions A1,B1,C1,D1,E1 et F1
dans la deuxième image pendant un temps égal à △t.
A1 B1
A
B
D1
C1
C D
E1
E F F1
Fig. 2.4 – Déplacement des particules.
A partir de cette figure 2.4, les deux composantes de la vitesse instantanée selon X et
Y sont calculables.
VX = △X △ t.
VY = △Y △ t.
△t : Temps entre deux images successives (Le temps entre deux tirs laser).
26
△X : Distance entre la première et la deuxième position pour une particule selon l’axe
X.
△Y : Distance entre la première et la deuxième position pour une particule selon l’axe
Y.
V : Vitesse instantanée pour une particule dans la tranche éclairée.
En réalité, l’image PIV contient un grand nombre de particules. Un calcul d’inter-
corrélation est utilisé pour identifier la position des particules dans les deux images suc-
cessives. La répétition des mesures permet d’obtenir le champ de vitesses moyennées
[56, 57].
2.3.1 Mise en place des mesures PIV

Dans la soufflerie de l’ENSAM-Paris, la partie en plexiglas transparente dans le pla-
fond de la veine d’essai a permis d’installer les lasers Nd-Yag au-dessus de la veine.
Comme la veine est semi guidée et sans parois latérales, la caméra a été installée à coté
de la veine pour assurer la prise d’images dans le sens normal au plan d’éclairage laser.
Ainsi, le matériel PIV est entièrement à l’extérieur de la veine d’essais et ne perturbe pas
l’écoulement, voir la figure 2.5.
L’ensemencement est réalisé grâce à l’utilisation de gouttelettes d’huile d’olive. Pour
avoir un écoulement avec une répartition homogène des particules d’ensemencement dans
la tranche éclairée par le laser, qui aide à bien représenter les vecteurs du champ de
vitesse, un générateur est utilisé pour pulvériser des gouttelettes d’huile dans le diffuseur
de la soufflerie. Le temps nécessaire pour parcourir la veine de retour et la chambre de
tranquilisation a permis une excellente répartition des gouttelettes dans l’écoulement en
arrivant à la veine d’essai.
Pour synchroniser les mesures PIV avec la position angulaire du rotor, un codeur
optique situé sur le mat de l’éolienne a été utilisé pour viser une cible tournante collée sur
le moyeu et émettre un signal de référence temporelle qui indique la position angulaire du
rotor. Ce signal permet de contrôler le déclenchement des tirs laser et la prise d’images.
Le temps entre deux tirs laser a été réglé à : △ T = 150 µs.
Pour explorer le volume occupé par le sillage et les trajectoires des tourbillons mar-
ginaux, nous avons mesuré le champ de vitesse dans plusieurs plans azimutaux avec les
angles 0, 30, 60 et 90º ; le zéro correspondant à la position verticale de pale, voir la figure
2.5.
A cause de la limite de la capacité de notre matériel, il est impossible d’explorer le
champ de vitesse derrière la machine par une seule image. Pour surmonter ce problème,
27
Fig. 2.5 – Banc d’éssai.
nous avons divisé le champ de vitesses en six fenêtres (3 horizontales × 2 verticales), voir
la figure 2.6.
L’échelle des fenêtres et leurs positions par rapport au rotor ont été définies à l’aide de
mires placées dans le plan de l’exploration après chaque série d’essais. Ainsi, quatre séries
de 6 fois 95 paires d’images ont été réalisées afin d’obtenir les différents plans d’exploration
azimutaux définis avec les angles 0, 30, 60 et 90º ; figure 2.5.
En raison de la vitesse de rotation élevée de l’éolienne (autour de 1050 tour/min),
la caméra prenait une paire d’images tous les deux tours. Pour chaque fenêtre, la prise
d’images était répétée 95 fois afin de reproduire une séquence temporelle d’environ 12
secondes et d’améliorer le calcul de la vitesse moyenne.
28
h2 h3
h1
m2 m3
m1
Fig. 2.6 – Six fenêtres d’exploration h1, h2, h3 (haut), et m1, m2, m3 (bas) et champs
des vitesses moyennées.
2.3.2 Résultats
La figure 2.7 présente une image brute prise dans la fenêtre h1 (1600 × 1200 pixels)
directement derrière le rotor en position verticale de la pale. Cette image montre la pale
(à gauche) et les noyaux des tubes tourbillonnaires issus des extrémités des autres pales.
Noyaux
R=0.25 m
Pale
Fig. 2.7 – Image brute, fenêtre h1.
29
En raison des vitesses induites par les tubes tourbillonnaires et qui écartent les par-
ticules d’ensemencement du centre du noyau, ce dernier apparaı̂t sur l’image comme un
cercle noir par manque de particules éclairées. Le traitement des images brutes permet de
constituer les champs des vitesses instantanées et moyennées.
Après avoir obtenu un nombre suffisant d’images par la technique PIV, le traitement de
celles-ci aura lieu par un algorithme itératif multi-résolutions, Susset [58]. Nous présentons
les résultats dans les parties suivantes.
2.3.3 Champs de vitesse

Dans ce paragraphe, nous allons représenter les champs des vitesses instantanées et
moyennées. La figure 2.8 montre le champ de vitesse instantanée correspondant à l’image
brute dans la fenêtre h1, déjà présentée sur la figure 2.7. Ce champ montre clairement l’in-
tersection des tubes tourbillonnaires émanant des extrémités des pales et leurs interactions
avec le plan d’exploration.
Par ailleurs, on peut voir facilement l’augmentation du diamètre du tube de courant
comme résultat du ralentissement de l’écoulement créé par l’éolienne, comme le prévoit
la théorie de Froude-Rankine [59].
Fig. 2.8 – Champ de la vitesse instantanée, fenêtre h1.
30
Pour obtenir le champ de vitesse moyennée pour chaque fenêtre, on a moyenné les
champs instantanés résultant des séries d’images prises (95 images). La figure 2.9 montre
le champ de la vitesse moyennée pour la fenêtre h1.
Fig. 2.9 – Champ de la vitesse moyennée, fenêtre h1.
Il est à signaler que l’examen des champs successifs dans le temps montre une fluc-
tuation de la position des noyaux des tourbillons marginaux au delà d’une distance de
l’ordre 0.75 D par rapport au plan de rotation. Cette fluctuation est due à la variation
temporelle de la puissance du rotor. En effet, lors d’une augmentation de la puissance
absorbée, la vitesse moyenne dans le sillage diminue conduisant ainsi à une diminution du
pas des tourbillons marginaux [56]. En conséquence, le respect de l’équation de continuité
fait que le rayon du tube de courant augmente et déplace les tourbillons vers l’extérieur.
On peut noter que la fluctuation radiale des noyaux dépend faiblement de leurs positions
axiales. Par contre, leur battement axial augmente rapidement, parce ce que ce battement
est amplifié en fonction du nombre de pas parcourus. En conséquence du battement du
noyau, une réduction de l’intensité du tourbillon apparaı̂t lors du calcul de la moyenne
des champs instantanés. Ainsi, le champ de vitesse stationnaire (champ moyenné) n’est
pas tout à fait représentatif du développement du sillage.
A l’aide des mesures des champs de la vitesse axiale moyennée dans la fenêtre h1 et m1
dans les quatre positions azimutales, 0, 30, 60 et 90°, une représentation tridimensionnelle
31
est obtenue pour le sillage proche du rotor, voir la figure 2.10. Nous observons à l’extrémité
de chaque pale un tube tourbillonnaire qui sort en trajectoire hélicoı̈dale. Ainsi, cette figure
2.10, donne une idée de la forme de sillage et le volume occupé par celui-ci.
0°
30°
60°
90°
Fig. 2.10 – Représentation 3D du sillage proche.
2.3.4 Analyse des résultats

L’analyse des vitesses dans le sillage du rotor, figure 2.11, montre que le moyeu est
une source de perturbations et d’instationnarités aérodynamiques fortes. Ceci vient de la
forme encombrante du moyeu qui comporte la génératrice.
Par ailleurs, le faible angle de calage au pied de la pale conduit à un blocage (avec
le moyeu) important dans le chemin de l’écoulement. Par conséquence, un ralentissement
important aura lieu derrière le moyeu, qui donne une raison de l’existence d’un courant
de retour.
L’incertitude sur les mesures de vitesse au voisinage de la pale vient des difficultés
spécifiques de la technique PIV pour explorer l’écoulement près des parois, par exemple,
32
Courant
de retour
Fig. 2.11 – Lignes de courant dans les plans h1 et m1, angle azimutal 0°.
V=0
Fig. 2.12 – Champ de la vitesse axiale induite par les tubes tourbillonnaires.
la réflexion du laser sur la surface de la pale, la coupe du trajet suivi par un particule au
voisinage de la pale et la valeur élevée de la vitesse tangentielle à quelques millimètres de
la pale, voir figure 2.9.
Pour étudier l’influence des tourbillons marginaux, nous avons présenté le champ de
la vitesse induite par les tubes tourbillonnaires (vitesse axiale locale moins la vitesse en
amont), le résultat est représenté sur la figure 2.12. Ainsi, la zone tracée par une ligne
en pointillés au niveau des noyaux, avec une vitesse nulle, divise l’écoulement en deux
33
parties : une partie interne où la vitesse axiale est ralentie par les tourbillons marginaux
et une partie externe où la vitesse est accélérée par les tourbillons. Cette image peut
rappeler à un certain point un modèle tourbillonnaire simple (modèle de Prandtl) où le
rotor est représenté par un disque de tourbillons annulaires et le sillage par un cylindre
avec des tourbillons annulaires espacés d’une distance égale au pas de l’hélicoı̈de divisée
par le nombre de pales.
Le développement du sillage tourbillonnaire (variation du diamètre) présenté sur la
figure 2.11 et figure 2.12 montre clairement pourquoi les modèles tourbillonnaires qui
supposent un diamètre de sillage constant, ne permettent pas d’obtenir de bons résultats.
34
2.4 Mesures par fil chaud
Les mesures PIV, permettent d’obtenir les composantes axiale et radiale de vitesse
dans les plans azimutaux. Pour obtenir la composante tangentielle, les explorations sont
faites à l’aide de l’anémométrie à fil chaud.
L’anémométrie à fil chaud est une technique de mesure qui permet de déterminer la
vitesse instantanée d’un fluide s’écoulant autour d’une sonde à fil chaud placée au sein de
l’écoulement.
Depuis son invention par L.V.King en 1914 [60], cet instrument de mesure est très
utilisé dans le domaine de la mécanique des fluides. Il y a deux méthodes de mesure par
fil chaud : l’anémomètre à température constante (CTA) [61] et l’anémomètre à courant
constant (CCA) [62]. Dans nos travaux, nous avons utilisé l’anémomètre à température
constante (CTA), car celui-ci a l’avantage d’avoir une plus grande bande passante ainsi
qu’une fréquence de coupure plus élevée allant jusqu’à 1 MHz.
2.4.1 Principe de l’anémométrie à fil chaud

Le principe de cet appareil CTA, est basé sur la loi de transfert de chaleur entre un
écoulement fluide et un élément très sensible et très fin, de l’ordre de 5 µm de diamètre et
de longueur 1 mm, chauffé par un courant électrique, et placé dans cet écoulement, figure
(Fig.2.13).
V1
E
Q
Pont de Wheatstone Servoamplificateur
Fil chaud
1
V
Q
Fig. 2.13 – Principe du fil chaud.
Il s’agit de maintenir à température constante un fil en tungstène traversé par un

courant électrique. Si l’on met une sonde à fil chaud dans un écoulement fluide, le fil est
refroidi par convection. Pour maintenir la température du fil constante, on doit augmenter
35
l’intensité électrique le traversant de façon à apporter une énergie (RΩ × I 2 ) égale à celle
perdue par convection. RΩ : Résistance de fil chaud. I : Intensité du courant électrique.
Le fil chaud est connecté avec un pont de Wheatstone et chauffé par un courant
électrique. Une variation de la résistance du fil chaud, qui est fonction de la température,
rend le pont déséquilibré. Un servoamplificateur est connecté au pont, celui-ci, maintient
toujours le pont en équilibre en contrôlant le courant qui passe dans le fil chaud, par
conséquence, la température du fil reste constante. Sans écoulement, on s’arrange pour
avoir un courant de chauffage I. Avec un écoulement, une quantité de chaleur Q est
transférée à l’écoulement et on mesure une tension E (la tension de pont de Wheatstone).
Afin de maintenir la température du fil chaud constante, et donc le pont en équilibre, le
courant I, renforcé par la variation la tension E, chauffe le fil chaud. Donc, une variation
du courant I, et de la tension E, permet de capter les variations de la vitesse d’écoulement
selon la loi de King :
E 2 = G + jVlϕ (2.1)
ϕ1
E2 − G

Vl = (2.2)
j
- ϕ ≈ 0.5
- G, j : constantes de calibrage de la pont de Wheatstone.
- Vl : vitesse locale.
A partir du principe du fil chaud, plusieurs catégories de fil chaud sont réalisables.
Avec l’utilisation d’un fil chaud, on peu mesurer une seule composante de la vitesse. La
construction de deux ou trois fils chauds ensemble, donne la possibilité de mesurer deux ou
trois composantes de vitesse simultanément. La figure 2.14, représente les cas précédents.
Fig. 2.14 – Trois catégories de sondes à fil chaud.
36
Un avantage d’utiliser le fil chaud est qu’il permet la mesure du taux de turbulence.
De plus, les sondes à deux ou trois fils permettent la mesure des écoulements turbulents
car ils ont un faible encombrement et ont un temps de réponse très rapide [63]. Ils ne
perturbent presque pas l’écoulement et permettent l’analyse de fluctuations très rapides.
2.4.2 Méthodologie suivie dans l’exploration au fil chaud

Le but de cette exploration est d’étudier les deux champs des vitesses, axiale et tan-
gentielle, en amont et en aval d’une éolienne.
Des explorations du sillage de la même éolienne testée précédemment par la méthode
PIV, sont réalisées dans la veine d’essai, figure (Fig. 2.3), avec une vitesse axiale d’écoulement
d’environ 9.4 m/s.
Par contre, la vitesse de rotation de l’éolienne était approximativement de 1050 tr/min.
Le fil chaud utilisé pour ce travail est de type 55 P11 (Dantec Mesurement Tchnology), il
est fixé à l’aide d’un support porté par un robot d’exploration et branché avec un câble
de (5 m) à la chaı̂ne de mesure.
2.4.3 Plan d’exploration

Dans la veine d’essai de notre soufflerie, nous avons essayé d’explorer le plan le plus
grand possible. Malgré les dimensions limitées de la veine(1.35m × 1.65m × 2m), nous
avons pu explorer l’écoulement, jusqu’à 0.4D en amont du plan de rotation et 1.5D en
aval environ.
Le plan d’exploration a été défini de la façon suivante : on commence la mesure au
niveau de l’axe de rotation, le pas selon l’axe Y (axe vertical) est toujours égal à 0.05 m,
et les pas horizontaux sont variables. En amont, on a commencé à la distance 0.015 m du
plan de rotation, et chaque pas est égal à 0.03 m. Par contre, en aval, on a commencé à
partir de 0.015 m du plan de rotation avec un pas de 0.05 m jusqu’à la huitième mesure,
puis avec un pas égal à 0.1 m pour les trois dernières mesures, figure 2.15.
Pour chaque rayon de mesure choisi, le fil chaud enregistre sur chaque tour de rotation
800 valeurs pour la vitesse, ces valeurs sont synchronisées azimutalement par rapport aux
pales grâce à un capteur optique qui envoie un signal de synchronisation à chaque passage
d’une cible collée sur le moyeu, figure 2.16. L’acquisition est répétée des centaines de fois
pour obtenir la vitesse moyennée.
37
0.03 m
0.05 m
0.1 m
0.05 m
Fig. 2.15 – Plan d’exploration du sillage autour du rotor.
Fig. 2.16 – Points de mesures sur des cercles complets.
38
2.4.4 Représentation du sillage
A partir de la figure 2.16, une interpolation de vitesse à partir des points obtenus est
réalisée afin d’avoir un champ de vitesse continu sur tout le domaine exploré comme le
présente la figure 2.17.
Fig. 2.17 – Représentation du sillage.
La coupe longitudinale de ce volume est représentée sur la figure 2.18, sur cette figure
on observe que le déficit de la vitesse axiale est très clair derrière le rotor. On observe
aussi l’intersection entre cette coupe et les tubes tourbillonnaires émis par chaque pale.
39
Fig. 2.18 – Champ de vitesse axiale dans le sillage.
2.4.5 Étude de trajectoire des tourbillons marginaux

Pour visualiser les phénomènes tourbillonnaires dans le sillage étudié, nous pouvons
examiner des coupes cylindriques et coupes axiales.
Les coupes cylindriques sont réalisées à trois différents diamètres. D’après la figure
2.19, où le diamètre est de 540 mm, les tourbillons sont visibles seulement à coté du rotor.
Parce que les tourbillons, en s’éloignant du plan de rotation, s’éloignent aussi de l’axe de
rotation et leurs positions deviennent dehors du diamètre 540 mm. Avec un diamètre de
sillage égal à 560 mm, on trouve que des tourbillons sont visibles en s’éloignant du rotor
, figure 2.20. Par contre, pour la figure 2.21, R = 290 mm, on trouve que les tourbillons
lointains sont plus clairs. On constate d’après les résultats des trois figures 2.19, 2.20 et
2.21 que les tourbillons se détachent des extrémités des pales, prennent des trajectoires
hélicoı̈dales par rapport à l’axe de rotation, et s’éloignent du celui-ci en fonction de leur
distance au plan de rotation.
40
Fig. 2.19 – Champ de la vitesse axiale, R = 270 mm.
41
1 2
Fig. 2.22 – Vitesse axiale, Z/D = 0.32. Fig. 2.23 – Vitesse axiale, Z/D = 0.48.
3 4
5 6
Fig. 2.26 – Vitesse axiale, Z/D = 1. Fig. 2.27 – Vitesse axiale, Z/D = 1.2.
42
7 8
D’autre part, nous avons présenté plusieurs sections axiales (normales à l’axe de ro-
tation) successives du sillage en fonction de leurs distances au plan de rotation pour
montrer aussi les trajectoires des tourbillons au fur et à mesure de l’éloignement du plan
de rotation.
A partir de la figure 2.22, à la distance Z/D = 0.32, on trouve les tourbillons qui
quittent les extrémités des pales. Nous avons indiqué un tourbillon par une flèche, celui-ci
vient de se détacher de l’extrémité d’une pale. Sur la figure 2.23, à la distance Z/D = 0.48,
les tourbillons s’éloignent de l’axe de rotation avec des positions azimutales différentes, le
tourbillon indiqué a pris un angle azimutal plus grand que dans le cas précédent. Dans la
figure 2.24, à la distance Z/D = 0.62, on trouve que les tourbillons s’éloignent encore de
l’axe de rotation.
Par contre, la figure 2.25, à la distance Z/D = 0.8, montre bien que les tourbillons sont
clairement d’un rayon plus grand que R. Ainsi, nous observons que le tourbillon indiqué
a tourné autour l’axe de rotation de plus que 120°.
Les figures 2.26, 2.27 et 2.28 montrent bien l’éloignement des tourbillons de l’axe
de rotation de l’éolienne et les mouvements azimutaux. A la distance égale à Z/D =
1.4, les tourbillons ont complètement disparu parce qu’ils quittent le domaine de notre
exploration, voir la figure 2.29.
43
2.4.6 Effet de l’éolienne sur l’écoulement à l’amont du rotor
Pour présenter l’influence du rotor sur l’écoulement en amont, nous avons réalisé quatre
sections successives à différentes distances par rapport au plan de rotation.
Fig. 2.30 – Vitesse axiale à Z/D = 0.4 du plan de rotation en amont.
Pour la première section, figure 2.30, à une distance égale à Z/D = 0.4, on observe
que qu’il y a quatre régions de vitesses différentes à cause de l’influence du rotor éolien,
la diminution de la vitesse axiale est plus grande à l’approche de l’axe de rotation de
l’éolienne, parce que la corde de la pale augmente de 0.047 m à l’extrémité de pale jusqu’à
0.067 m au pied de pale où se trouve le moyeu qui fait une obstacle devant le l’écoulement.
La vitesse la plus faible dans cette section est d’environ 3 m/s.
En nous rapprochant du plan de rotation à une distance de Z/D = 0.25, figure 2.31,
on observe l’apparition de l’influence des trois pales du rotor, et la surface où la vitesse
est inférieure à 3 m/s est plus grande par rapport à la section précédente.
A une distance du plan de rotation égale à Z/D = 0.125, figure 2.32, la présence de
l’influence des pales devient plus nette, la vitesse est plus faible que dans le cas précédent,
environ 1 m/s.
44
45
A la distance de Z/D = 0.02 du plan de rotation, figure 2.33, la présence des pales
et du moyeu est encore plus nette, l’influence des tourbillons détachés de l’extrémité de
chaque pale est présente avec des régions de vitesse élevées (environ 11 m/s ).
2.4.7 Champs de vitesse tangentielle

Nous avons utilisé l’anémométrie à fil chaud, pour avoir la possibilité de mesurer le
champ de vitesse tangentielle qu’on n’a pas pu mesurer avec la PIV. Dans ce paragraphe,
nous allons présenter le champ de vitesse tangentielle en amont et en aval de l’éolienne
testée. Le vent arrive face à l’éolienne sans rotation, mais par l’influence de la rotation des
pales, une vitesse tangentielle va apparaı̂tre. La figure 2.34, présente la vitesse tangentielle
à distance de Z/D = 0.048 en amont du plan de rotation ; a partir de cette figure 2.34,
nous observons que, la vitesse tangentielle prend des valeurs très faibles (très proche de
zéro) sur la majorité de cette surface, mais dans les zones très proches des pales, la vitesse
augmente jusqu’à 1.5 m/s.
Après le passage du vent à travers les pales de rotor, la vitesse tangentielle va augmen-
ter par rapport à la vitesse tangentielle en amont du plan de rotation, voir la figure 2.35.
Cette augmentation est due aux influences des tourbillons attachés, des tourbillons margi-
46
Fig. 2.34 – Vitesse tangentielle à Z/D = 0.048 en amont du plan de rotation.
naux et des tourbillons des pieds des pales qui forment une vitesse induite très importante
[37].
Fig. 2.35 – Vitesse tangentielle à Z/D = 0.048 en aval du plan de rotation.
La figure 2.35, représente la vitesse tangentielle à une distance de Z/D = 0.048 en

aval du plan de rotation, nous observons que, au pied des pales, la vitesse tangentielle est
élevée à cause de l’influence des tourbillons en pieds des pales. Par contre, sur cette figure
47
2.35, la vitesse tangentielle est aussi augment à l’extrémités des pales par l’influence des
tourbillons marginaux.
En s’éloignant en aval du plan de rotation, la vitesse tangentielle diminue car l’influence
des tourbillons est diminuée. Sur la figure 2.36 à la distance de Z/D = 0.1 du plan de
rotation, la vitesse tangentielle est moins forte que sur la figure 2.35 à Z/D = 0.048 en
aval du plan de rotation.
Fig. 2.36 – Vitesse tangentielle à Z/D = 0.1 en aval du plan de rotation.
Enfin, la figure 2.37 représente un champ de vitesse tangentielle le long du plan d’explo-
ration. Nous observons que la vitesse tangentielle est très importante sur les frontières du
sillage, dans la zone d’intersection entre les tubes tourbillonnaires et le plan d’exploration.
2.5 Conclusion
Dans ce travail expérimental, nous avons mis en œuvre la technique PIV et la technique
d’anémomètrie à fil chaud pour l’exploration en soufflerie de l’écoulement et l’obtention
de données quantitatives sur le champ de vitesse dans le sillage proche d’une éolienne à
axe horizontal.
Avec la technique PIV, la synchronisation des tirs lasers et des prises d’images par
rapport à la position azimutale de la pale, a permis de visualiser l’écoulement dans un
repère lié au rotor et de réaliser une reconstitution tridimensionnelle du champ de vitesse.
Ainsi, la position des tourbillons marginaux a pu être localisée. Les résultats montrent que
les tourbillons marginaux issus des extrémités des pales ne sont pas situés sur une surface
48
Fig. 2.37 – Champ de vitesse tangentielle.
cylindrique comme le suppose la théorie tourbillonnaire linéaire [64]. Ils se déplacent vers
l’extérieur en augmentant le diamètre du tube de courant comme le prévoit la théorie de
Froude-Rankine [59]. Pour élargir le champ mesuré, les explorations PIV, ont été réalisés
sur 6 fenêtres voisines présentant un certain chevauchement.
Le traitement effectué après les mesures a permis de caler l’ensemble des résultats de
ces fenêtres et d’obtenir le champ de vitesse jusqu’à 1.5D en aval le rotor. L’analyse des
résultats montre l’interaction entre l’écoulement en aval du rotor et les vitesses induites par
les tourbillons marginaux. Les mesures ont révélé la présence de structures tourbillonnaires
importantes en aval du moyeu et près du pied des pales.
La technique de l’anémométrie à fil chaud a permis d’explorer le champ de vitesse
axiale et tangentielle en amont et en aval du rotor. Les essais par fil chaud, montrent
bien le volume occupé par le sillage proche, les détachements des tourbillons marginaux
de l’extrémité de chaque pale et leurs trajectoires qui ne sont pas situées sur une surface
cylindrique. L’exploration de l’éolienne dans la veine de retour permet d’étudier le sillage
lointain jusqu’à 6D, celui-ci est assez perturbé.
49
50
Chapitre 3
Développement d’un modèle hybride

basé sur le concept du disque actif
51
3.1 Introduction :
Dans ce chapitre nous présentons les éléments théoriques que nous avons utilisés pour
la construction d’un modèle hybride basé sur le concept du disque actif de cette thèse.
Nous commençons par la présentation de la méthode de l’élément de pale BEM (Blade
Element Momentum), qui est utilisée pour déterminer les performances d’une éolienne
et les forces appliquées par celle-ci sur l’écoulement. Ensuite, on présente le modèle du
disque actif, où le rotor est remplacé par un disque d’épaisseur nulle et de diamètre égal
à celui-ci du rotor. Dans la pratique, ce disque aura une épaisseur E.
Des efforts aérodynamiques volumiques sont calculés à l’aide de la théorie de l’élément
de pale et sont répartis sur un ensemble d’anneaux coaxiaux d’épaisseur radiale dr.
L’écoulement global est calculé par la résolution numérique des équations de Navier-Stokes
turbulentes moyennées. Cet écoulement fournit en particulier la vitesse locale incidente à
chaque pale de l’éolienne et permet d’évaluer les efforts aérodynamiques appliqués sur le
disque actif. Comme nous l’avons déjà fait remarquer, cette simplification est faite pour
diminuer le coût du calcul en diminuant le nombre de cellules du maillage. De part sa
simplicité géométrique, il suffit de mailler un disque actif et son environnement, le temps
consacré à la réalisation des maillages est bien plus court.
Nous validons notre modèle à l’aide de données expérimentales issue du NREL et des
essais effectués dans la soufflerie de l’ENSAM Paris. Trois éoliennes de référence sont
choisies ici. Les puissances en fonction de la vitesse du vent sont validées par les données
des éoliennes NREL phase II [51] et NREL phase VI [65]. Le champ de vitesse calculé est
comparé aux mesures expérimentales relevées sur une éolienne Rutland 503 testée dans
la soufflerie à l’ENSAM Paris.
3.2 Méthode de l’élément de pale

La méthode la plus utilisée pour calculer la charge appliquée par l’écoulement sur une
éolienne et le comportement d’une éolienne, consiste à utiliser le bilan de quantité de
mouvement appliqué à un élément de pale (Blade Element Momentum). Cette méthode
est basée sur la division de l’écoulement en volumes de contrôle annulaires d’épaisseur
dr, auxquels on applique le bilan de quantité de mouvement et d’énergie. Ces anneaux
s’étendent de l’infini amont jusqu’à l’infini aval par rapport au rotor, figure 1.4, page 14.
Les hypothèses principales de cette méthode sont :
– que la vitesse induite dans le plan de rotation est égale à la moitié de la vitesse
induite à l’infini aval.
52
– que l’on peut analyser l’écoulement par la division de la pale en nombre d’éléments
indépendants et que la force d’un élément de pale est seule responsable de la variation
de quantité de mouvement de l’air qui passe dans l’anneau balayé par cet élément.
Cette hypothèse est acceptable si on suppose qu’aucune interaction radiale n’existe
entre les écoulements qui passent dans des anneaux voisins balayés par les éléments
de pales. Pratiquement, le facteur d’induction axial de l’écoulement est rarement
uniforme, mais les résultats expérimentaux de Lock (1924) sur une hélice traversée
par un écoulement prouvent que l’hypothèse de l’indépendance radiale est acceptable
[37].
3.2.1 Théorie de l’élément de pale

Dans cette théorie, les forces agissant sur un élément de pale sont calculées en uti-
lisant les caractéristiques de portance et de traı̂née d’un profil bidimensionnel. Ces ca-
ractéristiques sont obtenues en utilisant l’angle d’attaque déterminé à partir de la vitesse
relative de l’écoulement dans le plan perpendiculaire à l’élément de pale (la vitesse W ,
figure 3.2). Ainsi, on néglige l’interaction entre les tubes d’écoulement correspondant aux
éléments voisins [66].
On considère une éolienne avec B pales, leur rayon est R , leur corde est c, et l’angle
de calage est β mesuré entre la ligne de portance nulle du profil et le plan de rotation.
La vitesse de rotation de la pale est Ω, la vitesse du vent à l’infini amont est V1 . Ωr est
la vitesse tangentielle de l’élément de pale tandis que l’on utilise les facteurs d’induction
axial a et tangentiel a′ pour obtenir la vitesse relative de l’écoulement 3.1. D’après la 3.2
on peut écrire :
q
W = V12 (1 − a)2 + Ω2 r2 (1 + a′ )2 (3.1)
On peut également tirer les relations suivantes :
V1 (1 − a) Ωr(1 + a′ )
sin φ = , cos φ = (3.2)
W W
L’angle d’incidence :
α=φ−β (3.3)
53
Ωrá
V(1-a)
1
dr
Ωr r
r
Ωr
Fig. 3.1 – Un élément de pale au rayon local r et un anneau balayé par cet élément.
W`i
W`i `
W
i
v
Fn
W
Fz W
V1 (1-a)
V1
φ
Fx α
Ft β -θ
- Ωr Vθ
Fig. 3.2 – Vitesses et forces agissant sur un élément de pale.
−
→
La vitesse induite Wi a pour composantes :
−
→
Wi = −aV1~z − a′ Ωr−
→
eθ (3.4)
54
La portance générée par un élément de pale d’épaisseur dr, normale à la direction de
W est :
dFz = 0.5ρW 2 Cz c dr (3.5)

où Cz : Coefficient de portance.
−
→
Et la traı̂née parallèle à W est :
dFx = 0.5ρW 2 Cx c dr (3.6)

où Cx est le coefficient de traı̂née.
On peut reprojeter ces efforts afin d’obtenir les composantes axiales ou normales et
tangentielles. La force axiale est donnée par :
dFn = dFz cosφ + dFx sinφ (3.7)
La force tangentielle est donnée par :
dFt = dFz sinφ − dFx cosφ (3.8)
En adimensionnant ces relations on peut obtenir le coefficient de force normale au plan

de rotation :
Cn = Cz cos φ + Cx sin φ (3.9)

De même, le coefficient de force parallèle au plan de rotation, pour un élément de pale
est :
Ct = Cz sin φ − Cx cos φ (3.10)
3.2.2 Bilan de quantité de mouvement appliqué à un élément

de pale
L’étude que nous venons de faire et plus particulièrement l’équation (3.9) nous permet
d’évaluer la résultante des forces aérodynamiques de B éléments de pales dans la direction
axiale au rayon r :
1
dN = ρW 2 Bc(Cx sin φ + Cz cos φ)dr (3.11)
2
D’autre part, cette même résultante est égale à la variation de la quantité de mouve-
ment axiale de l’air passant à travers un anneau de rayon r, d’où :
55
1
dN = 4πρV12 a(1 − a)rdr + ρ(2a′ Ωr)2 2πrdr (3.12)
2
où on a tenu compte de la contribution axiale de la variation de quantité de mouvement
produite par le facteur d’induction tangentielle.
On peut alors égaler les deux écritures de dN :
1
ρW 2 Bc(Cz cos φ + Cx sin φ)dr = 4πρ[V12 a(1 − a) + (a′ Ωr)2 ]rdr (3.13)
2
On note que µ = Rr , ceci permet de simplifier l’expression ci-dessus :
W2 c
2
B (Cz cos φ + Cx sin φ) = 8π(a(1 − a) + (a′ λµ)2 )µ (3.14)
V1 R
Nous pouvons mener la même étude sur les couples et sur la variation de quantité de
mouvement angulaire. Le couple exercé sur B éléments de pale au rayon r s’écrit :
1
dC = ρcB(Cz sin φ − Cx cos φ)W 2 rdr (3.15)
2
D’autre part, la variation de quantité de mouvement angulaire de l’air passant à travers
un anneau est :
dC = ρV1 (1 − a)Ωr2a′ r2πrdr = 4πρV1 (Ωr)a′ (1 − a)r2 dr (3.16)

L’égalité des deux expression du couple dC permet d’écrire :
1
ρW 2 Bc(Cz sin φ − Cx cos φ)rdr = 4πρV1 (Ωr)a′ (1 − a)r2 dr (3.17)
2
Une méthode BEM classique utilise le système de deux équations (3.13) et (3.17)
à deux inconnues a et a′ . Ce système, généralement résolu à l’aide d’une méthode de
Newton, permet d’obtenir les caractéristiques locales de l’éolienne.
A partir de ces équations, on trouve que :
1
a= 2 (3.18)
( 4FPσCsinn φ ) + 1
1
a′ = (3.19)
( 4FP sin
σCt
φ cos φ
) −1
où la correction de Prandtl FP et la solidité de l’élément de pale σ, sont déterminables
selon les équations suivantes [67, 68] :
2 B R−r
FP = arccos(e(− 2 r sin φ ) ) (3.20)
π
56
Bc
σ= (3.21)
2πr
Il suffit ensuite d’intégrer ces efforts et couples élémentaires le long du rayon de la pale
pour obtenir les performances globales de l’éolienne.
Ainsi, le couple total a pour expression :
Rmax Rmax
B B
Z Z
Cmeca = ρcCt W 2 rdr = ρc Ct W 2 rdr (3.22)
Rmin 2 2 Rmin
Ce couple permet aussi de calculer la puissance mécanique sur l’arbre du rotor, en

intégrant la puissance élémentaire :
B B
dPmeca = dCmeca Ω = ( ρcCt W 2 rdr)Ω = ρcCt W 2 rΩdr
2 2
Rmax Rmax
B B
Z Z
Pmeca = ρcCt W 2 rΩdr ⇒ Pmeca = ρcΩ Ct W 2 rdr (3.23)
Rmin 2 2 Rmin
La puissance électrique est donnée par une équation caractéristique de l’éolienne NREL
phase II [51] :
Pelec = 0.9036 · Pmeca − 0.847 (3.24)

Nous venons ici de rappeler les principaux éléments de la méthode des éléments de
pale. Nous allons maintenant l’appliquer au cas du calcul des performances du rotor de
l’éolienne NREL Phase II.
3.3 Application de la méthode de l’élément de pale

Nous avons appliqué la méthode de l’élément de pale sur une éolienne NREL phase
II, cette éolienne a trois pales de diamètre égal à 10.1 m, le diamètre du moyeu est égal
à 0.75 m, la corde des pales est constante c = 0.458 m, les pales ne sont pas vrillées, et
l’angle de calage est fixé à 12 degrés. L’essai est fait à une vitesse de rotation égale à 72
tr/min, le rotor est fixé sur une nacelle, et cet ensemble est situé à une hauteur de 16.8
m [51].
57
3.3.1 Algorithme de calcul
La figure 3.3 présente l’algorithme de calcul utilisé pour calculer les forces appliquées
par le rotor sur l’écoulement à l’aide de la théorie de l’élément de pale. A partir de celle-
ci, on observe que le calcul a commencé en connaissant trois valeurs, la vitesse du vent
à l’infini amont V1 , le facteur d’induction axial a et tangentiel a′ de l’écoulement. Ici on
suppose que a = a′ =0, puis on prend un élément de pale situé au rayon local r, on calcule
les caractéristiques du profil sur celle-ci, ensuite on détermine les nouvelles valeurs de a
et a′ , donc la différence entre a au départ de notre calcul et la nouvelle valeur de a est
déterminable, si cette différence est supérieure de la différence souhaitée, on recommence
le calcul jusqu’à convergence. La même méthode permet d’obtenir a′ .
Ensuite on applique cette opération sur n éléments de pale à partir du rayon inférieur
jusqu’au rayon supérieur de la pale. Enfin, nous avons les caractéristiques propres des n
éléments de pales qui nous permettent de calculer la puissance du rotor.
3.3.2 Discrétisation de la pale et critère de convergence

Le calcul par la méthode de l’élément de pale impose une discrétisation du domaine de
calcul. Afin de choisir une discrétisation satisfaisante, nous allons évaluer la puissance de
l’éolienne obtenue en fonction de la finesse de la discrétisation choisie. On effectue cette
étude sur le cas de l’éolienne NREL s809 phase II. Son rayon d’extrémité est égal à 5.05
m, le rayon du moyeu (le pied de la pale) est égal à 0.75 m, donc la longueur de la pale
est égale à 4.3 m. Nous l’avons divisée en n parties, avec n = 20, 50, 100, 250, 500, 700, 800
et 1000. Nous testons ici uniquement la méthode de calcul représentée sur la 3.3 sans
aucun couplage avec un logiciel externe, la vitesse incidente est imposée directement par
l’utilisateur.
Nous observons qu’à partir de 700 éléments de pale, dr = 0.006143 m, l’évolution de
la puissance en fonction de la vitesse du vent ne varie plus. Dans la suite de ce travail,
nous choisissons n = 1000 afin de s’assurer que la discrétisation spatiale est suffisante.
L’algorithme itératif permettant d’obtenir les coefficients d’induction a et a′ doit
également être étudié de manière à évaluer le critère de convergence requis pour arrêter
l’algorithme. On quitte la boucle itérative sur a et a′ lorsque la différence de ces valeurs
entre deux itérations successives |an+1 − an | et |a′ n+1 − a′ n | est inférieure à une valeur
prescrite notée ǫ. On teste notre algorithme sur l’éolienne NREL s809 Phase II pour des
valeurs de ǫ comprises entre 10−1 et 10−5 . Les résultats de cette étude sont représentés
sur la figure 3.5. A partir de ǫ < 10−3 , le calcul de la puissance en fonction de la vitesse
58
Initialisation des V1, a et a′
Pour R min ≤ r ≤ R max
Calculer les angles α et φ
Lecture Cx et Cz
Calculer C t et C n
Calculer FP
n+1
Calculer an+1et a′
n
| a n+1- a |≤ ε ? Non
n+1 n
| a′ - a′ |≤ ε ?
Oui
Non
r = R max ?
Oui
Calculer la puissance
Afficher les résultats
Fin
Fig. 3.3 – Algorithme de calcul selon la théorie de l’élément de pale.
59
Fig. 3.4 – Influence du nombre de divisions de la pale sur l’évolution de la puissance,
NREL phase II.
du vent ne change plus. Le critère de convergence appliqué dans tous nos calculs est fixé
à ǫ = 10−5 .
Les deux paramètres évoqués ci-dessus (discrétisation spatiale et critère de conver-
gence) ayant été validés. Nous présentons sur la figure 3.6, les résultats produits par la
méthode de l’élément de pale en comparant ses résultats à la puissance électrice mesurée
expérimentalement sur l’éolienne NREL s809 Phase II. L’accord entre les deux résultats
est satisfaisant étant donné que nous avons ici utilisé la méthode de l’élément de pale avec
un champ incident imposé par l’utilisateur.
60
Fig. 3.5 – Influence de la précision de calcul des coefficients d’induction sur l’évolution
de la puissance, éolienne NREL phase II.
Fig. 3.6 – Comparaison des puissances calculées par la méthode de l’élément de pale et
les données expérimentales pour l’éolienne NREL phase II.
61
3.4 Modèle de disque actif
En mécanique des fluides, le disque actif est défini comme une surface de discontinuité
où des forces de surface agissent sur l’écoulement. Ce modèle est extrêmement simplifié
et repose sur les hypothèses suivantes [59, 44] :
V1
Vd
P+ D V2
P-
Fig. 3.7 – Modélisation du rotor éolien par un disque actif.
- La géométrie du rotor est effacée et ce dernier n’est représenté que par un disque
d’épaisseur nulle de diamètre D.
- Le fluide est incompressible, non visqueux et non pesant.
- Les vitesses V1 à l’infini amont, Vd dans le plan du disque et V2 dans la veine à l’infini
aval sont uniformes et axiales.
- L’énergie spécifique de l’écoulement comporte deux parties : cinétique et potentielle
de pression.
La vitesse axiale dans le plan de rotation est définie en fonction de la vitesse à l’infini
amont par l’introduction du facteur d’induction axial a , soit :
Vd = (1 − a)V1 (3.25)
3.4.1 Equation de continuité

L’application de l’équation de continuité permet d’écrire :
ρA1 V1 = ρAd Vd = ρA2 V2 (3.26)
62
3.4.2 Bilan de quantité de mouvement
Quand le vent passe dans le tube de courant comportant le disque actif, il y a un
changement de vitesse égal à (V1 − V2 ), et le taux de variation de quantité de mouve-
ment est égal à la somme des efforts extérieurs appliqués. Comme le tube de courant est
complètement entouré par le vent à la pression atmosphérique, les forces à l’origine du
changement de quantité de mouvement viennent uniquement de la différence de pression
créée par le disque actif :
(Pd+ − Pd− )Ad = (V1 − V2 )ρAd Vd

Comme Vd = V1 (1 − a)
Il vient alors :
(Pd+ − Pd− )Ad = (V1 − V2 )ρAd V1 (1 − a) (3.27)

Pour obtenir la différence de pression (Pd+ −Pd− ), nous utilisons l’équation de Bernoulli
entre l’infini amont et le disque et entre le disque et l’infini aval :
A l’amont,
1 1
P1 + ρV12 = Pd+ + ρVd2 (3.28)
2 2
A l’aval,
1 1
P1 + ρV22 = Pd− + ρVd2 (3.29)
2 2
D’où :
1
∆P = Pd+ − Pd− = ρ(V12 − V22 ) (3.30)
2
A l’aide des équations (3.25), (3.27) et (3.30) on trouve que :
V2 = (1 − 2a)V1 (3.31)
En comparant les équations (3.25) et (3.31) on trouve que la vitesse induite dans le
plan du rotor est égale à la moitié de la vitesse induite à l’infini aval.
3.4.3 Coefficient de puissance

D’après les équations (3.27) et (3.31), on exprime la force appliquée sur le disque actif :
Ff orce = (Pd+ − Pd− )Ad = 2ρAd V12 a(1 − a) (3.32)

La puissance transmise au disque est :
63
P = Ff orce Vd = 2ρAd V13 a(1 − a)2 (3.33)
Le coefficient de puissance est défini par le rapport entre la puissance transmise au
disque actif et une valeur de référence correspondant à la puissance du vent amont tra-
versant une surface égale à celle du disque actif :
2ρAd V13 a(1 − a)2

Cp = 1
2
ρV13 Ad
Donc :
Cp = 4a(1 − a)2 (3.34)
3.4.4 Limite de Betz

L’équation (3.34) montre que le coefficient de puissance dépend du facteur d’induction
axial a. La valeur maximale de Cp est déterminée par :
dCp
= 4a(1 − a)(1 − 3a) = 0
da
soit
1
a=
3
Ce qui correspond à,
16
Cpmax = = 0.593 (3.35)
27
Cette valeur est appelée la limite de Betz et montre la limite supérieure théorique de
la puissance que l’on peut extraire du vent incident avec une éolienne.
3.5 Description du modèle hybride

L’épaisseur d’un disque actif est généralement considérée nulle [33, 69]. Toutefois,
l’introduction de notre modèle hybride dans un code industriel permet un meilleur contrôle
des efforts appliqués localement si l’on considère des efforts répartis sur un volume et non
sur une surface. En conséquence, nous utilisons ici un disque actif d’épaisseur E, comme
dans la figure 3.8.
Pour choisir l’épaisseur E, nous avons calculé par ce modèle hybride la puissance de
l’éolienne NREL phase II en fonction de l’épaisseur E pour une vitesse à l’infini amont de
10 m/s. A partir de la figure 3.9 on observe que, l’épaisseur du modèle hybride qui donne
64
Y
R
Moyeu
Z
θ
dr
∆θ X
r
Fig. 3.8 – Le modèle du rotor.
la puissance la plus proche de la puissance expérimentale (6.3 kW) pour l’éolienne NREL
phase II est de manière assez prévisible E = c × sin(β). Pour l’éolienne en question, on
a c × sin(β) = 0.208 × c = 0.09522 m. Cette épaisseur du modèle hybride est très proche
de l’épaisseur du volume cylindrique réel E´balayé par les pales qui est représenté sur la
figure 3.10.
Dans la suite de notre travail, on choisira donc l’épaisseur du disque actif comme suit :
E = c × sin β où β est l’angle de calage et c est la corde de la pale.
Ainsi, on divise le disque en un nombre d’anneaux n, chacun de longueur E et de
largeur radiale de dr, figure 3.8.
65
Fig. 3.9 – Évolution de la puissance en fonction de l’épaisseur du disque E/c, la vitesse
à l’infini amont est V1 = 10m/s, éolienne NREL phase II.
β
E΄ E
Fig. 3.10 – Définition de l’épaisseur du disque.
3.6 Etude de la discrétisation du disque actif

Nous étudions ici l’influence de la finesse du maillage hexaédrique utilisé pour discrétiser
le disque actif. A cette fin, nous allons calculer la puissance restituée par l’éolienne pour
un vent à l’infini V1 = 10 m/s avec des maillages de finesse croissante. Nous présentons
l’influence de ce facteur dans le cas NREL phase II.
Le rayon de disque est maillé avec 1000 mailles, le pas dans la direction azimutale
est égal à ∆θ = 3.6 degrés. Le calcul de la puissance de l’éolienne, avec cette densité du
maillage, la vitesse à l’infini amont étant égale à 10 m/s, donne une puissance égale à
7.51597 kW (la puissance expérimentale correspondante est égale à 6.3 kW ).
Ensuite, nous augmentons la densité du maillage progressivement avec la diminution
du pas dans la direction azimutal. Nous avons trouvé qu’après un pas azimutal égal à
∆θ = 0.90657 degré, la puissance obtenue par ce modèle hybride ne change pas et est
66
égale à 7.38852 kW , figure 3.11. Dans notre travail, nous avons maillé le modèle hybride
DA avec une densité du maillage correspondante à ∆θ = 0.90657 degré.
Fig. 3.11 – Influence de la finesse du maillage sur la puissance calculée, V1 = 10 m/s,

NREL phase II.
67
3.7 Interfaçage avec le logiciel Fluent
Une fonction définie par l’utilisateur UDF1 , permet de paramétrer le logiciel Fluent et
de le coupler avec un calcul spécial. Dans le cas de cette étude, l’intérêt va se porter sur la
possibilité de calculer une distribution volumique de force différente pour chaque maille
du système et non plus une distribution commune à toutes les mailles du même système.
Une UDF est une fonction programmée par l’utilisateur qui peut être liée avec le
solveur Fluent pour améliorer les capacités standard du code de calcul. Les UDF sont
programmées dans le langage informatique C. Elles utilisent à la fois les librairies de
fonctions de C et les macros prédéfinies avec le logiciel Fluent qui permettent d’accéder
aux données du solveur.
Les données du solveur auxquelles accède l’UDF sont :
* Coordonnées du centre de la maille.
* Vecteur vitesse associé à la maille.
Les paramètres utilisés comme données d’entrée par le solveur sont :
*Force volumique dans la direction x.
*Force volumique dans la direction y.
*Force volumique dans la direction z.
Ces paramètres ne peuvent être calculés que pour une maille volumique ; en conséquence,
il faut adapter le modèle surfacique de disque actif en un modèle volumique. Ceci fait l’ob-
jet du paragraphe suivant, qui reprend la théorie des éléments de pale.
3.8 Répartition volumique des forces appliquées sur

le disque actif
A partir de la théorie de l’élément de pale, on peut démontrer qu’un élément de pale
reçoit une force élémentaire par l’écoulement égale à :
1
dF~ = ρcCF W 2 dr~e (3.36)
2
3
ρ : Masse volumique de l’air. [kg/m ]
c : Corde du profil. [m]
CF : Coefficient caractéristique du profil.
~e : Vecteur unitaire définissant la direction dF~ .
1
UDF : User-Defined Function.
68
Cette force élémentaire multipliée par le nombre de pales peut être répartie uni-
formément sur un anneau de rayon r et d’épaisseur radiale infiniment petite dr. Or on
peut aussi ramener cette distribution à une distribution volumique de force en choisissant
un élément de volume annulaire.
Cet élément de volume est situé au rayon r, et a pour épaisseur radiale dr et pour
longueur axiale E. E représente l’épaisseur du cylindre 3.29, page 63 qui remplace le rotor.
L’élément de volume s’écrit donc comme :
dv = 2πrEdr (3.37)
En divisant les deux équations (3.36) et (3.37), on obtient :
dF~ 1
2
ρcCF W 2 dr
= ~e (3.38)
dv 2πrEdr
Soit l’intensité de la force égale à :
dF~ BρcCF W 2
f~ = B = ~e (3.39)
dv 4πrE
En projetant cette équation sur les cordonnées x, y et z, voir la figure 3.12, et à l’aide
de la figure 3.2 page 54 on trouve que :
La composante de l’intensité de la force appliquée par un élément de pale sur l’écoulement
dans la direction x est :
BρcCt W 2
fx = − sin θ (3.40)
4πrE
dans la direction y est :
BρcCt W 2
fy = − cos θ (3.41)
4πrE
θ : Angle azimutal. [deg]
dans la direction z est :
BρcCn W 2
fz = (3.42)
4πrE
Le coefficient de la force normale dans le plan de rotation est :
Cn = Cx sin φ + Cz cos φ (3.43)
69
Le coefficient de la force tangentielle dans le plan de rotation est :
Ct = Cz sin φ − Cx cos φ (3.44)
dr fZ
fX
Ω
r
fy
θ X
Fig. 3.12 – Position angulaire d’une maille et forces agissant sur elle-ci.
3.9 Algorithme de calcul

L’UDF programmée pour cette étude, calcule la force volumique au centre d’une maille
(point de rayon r et d’angle d’azimut θ ), appartenant au disque modélisant le rotor, puis
l’affecte à la maille dans son ensemble.
Cependant, puisque l’intensité de la force volumique f est appliquée uniformément
sur un volume- en l’occurrence une maille- la distribution de f sur le cylindre rotor est
discrétisée et la qualité de la discrétisation dépend de celle du maillage : autrement dit,
dr apparaı̂t implicitement et est fonction de la finesse du maillage du cylindre rotor.
Pour chaque maille appartenant au cylindre rotor, on va effectuer la boucle qui est
décrite dans la figure 3.13.
70
Lecture des coordonnées du centre de la maille. Mcentre = ( x centre ,ycentre ,z centre )
Calcul du rayon du centre de la maille. r=√x²centre + y²centre
Calcul de l’azimut du centre de la maille. tanθ = ycentre / x centre
Lecture de la vitesse axiale du vent sur la maille. Vair = |V z |
Calcul de la vitesse tangentielle. V tan = Vy cos θ - Vx sin θ
Calcul de la vitesse relative. W = √ V²air+( Ω r + V tan)²
Calcul de l’angle relatif φ et tanφ = Vair / ( Ω r + Vtan )

de l’angle d’incidence α. α=φ−β
Calcul des coefficients de traînée et de C x = f (α, Re), Cz = f (α, Re)

portance du profil, C x et Cz
Calcul des coefficients des forces dans et Cn = Cz cos φ + Cx sin φ

normale au plan de rotation, Ct et C n C t = Cz sin φ − C x cos φ
Calcul de l’intensité des forces sur les mailles. B ρ c Ct W²

fx= - sin θ
4лrE
B ρ c Ct W²
fy = - cos θ
4лrE
B ρ c Cn W²
fz =
4лrE
Fig. 3.13 – UDF algorithme.
Pour la simulation du sillage par le modèle hybride basé sur le concept du disque
actif, cette méthode donne accès à toutes les variables permettant de calculer plusieurs
grandeurs mécaniques : 71
* Effort de traı̂née de l’éolienne.
* Effort tangentiel qui s’exerce sur l’éolienne.
*Couple mécanique créé sur l’arbre de l’éolienne.
Ces calculs sont réalisés dans une UDF qui s’exécute au terme de chaque itération de
la simulation numérique. Pour calculer les efforts, on procède en multipliant l’intensité de
l’effort f par le volume de la maille considérée, puis en sommant les produits pour toutes
les mailles du cylindre rotor.
Pour l’effort axial ou effort de traı̂née, le calcul est :
nmailles
X
N= (fz,i V olumemailles,i ) (3.45)
i=1
Avec
nmailles : le nombre de mailles constituant le cylindre rotor.
fz,i : la force volumique dans la direction z pour la maille i.
V olumemaille,i : le volume de la maille i.
L’effort tangentiel est calculé comme suit :
nmailles
X q
2 2
T = ( fx,i + fy,i · V olumemaille,i ) (3.46)
i=1
Et le couple a pour valeur :

nmailles
X q
2 2
Couplemeca = ( fx,i + fy,i · ri · V olumemaille,i ) (3.47)
i=1
Avec fx,i et fy,i les forces volumiques dans les directions x et y pour la maille i.
La puissance mécanique se calcule en multipliant le couple mécanique par la vitesse
de rotation du rotor :
Pmeca = Couplemeca · Ω (3.48)

La puissance électrique est donnée par une équation caractéristique de l’éolienne
(NREL phase II) [51] :
Pelec = 0.9036 · Pmeca − 0.847 (3.49)
Avec la puissance électrique, on dispose d’une valeur que l’on va pouvoir comparer aux
résultats expérimentaux du NREL sur la même éolienne, et donc discuter de la validité
du modèle étudié. Par contre, la validation pour l’éolienne NREL phase VI a été faite par
rapport à la puissance mécanique.
72
3.10 Application du modèle
La validation de ce modèle hybride, le modèle hybride basé sur le concept du disque
actif, a été faite de deux manières différentes par rapport à trois types d’éoliennes. Tout
d’abord, nous avons comparé l’évolution des puissances obtenues par ce modèle hybride en
fonction de la vitesse du vent à l’infini amont, avec les données expérimentales concernant
les éoliennes de type NREL phase II et VI [51, 70]. Le profil de ces éoliennes est le s809,
ce profil est bien étudié à l’université de TU-Delft2 qui’a déterminé ses caractéristiques
aérodynamiques [71, 72, 73], nous présentons ces caractéristiques aérodynamiques dans
l’annexe A.
D’autre part, nous avons comparé le sillage simulé par ce modèle hybride qui représente
l’éolienne Rutland 503 avec le sillage que l’on a mesuré derrière cette éolienne dans la
soufflerie de l’ENSAM-Paris. Les propriétés aérodynamiques utilisées, du profil éolien
(Rutland 503), ont été obtenues numériquement.
La méthode de l’élément de pale (BEM) utilisée dans ce calcul, est basée sur les
polaires de profil en cas bidimensionnelle et donne des résultats réalistes dans le cas d’un
angle d’attaque avant le décrochage [74, 75]. Quand cette condition n’est pas respectée,
l’écoulement est tri-dimensionnel. Ici, l’écoulement décroché qui tourne avec la pale, crée
un écoulement dans la couche limite le long de la pale à cause des forces centrifuges, et le
long de la corde à cause des forces de Coriolis. Par conséquence, le décrochage est diminué
par rapport au cas bidimensionnelle, et donc le coefficient de portance est plus élevé que
dans le cas bidimensionnelle, spécialement près du pied de pale [76, 77, 78, 79].
Tous les résultats sont obtenus dans le cas d’un écoulement incompressible et station-
naire. Le modèle de turbulence utilisé est k − w − SST .
3.11 Résultats et analyse

Dans cette partie de notre travail, on va présenter les résultats du modèle hybride basé
sur le concept du disque actif sur trois types d’éoliennes.
3.11.1 Eoliennes comparées

Les tests pour les deux éoliennes NREL phase II et VI sont faits dans une grande souf-
flerie située à la NASA Ames Research Center Moffett Field, California. Cette soufflerie
est de grande dimensions, sa section est de 24.4 m × 36.6 m [80].
2
TU-Delft : Delft University of Technology.
73
Fig. 3.14 – Eolienne NREL phase VI dans la soufflerie de la NASA, au centre de recherche
d’Ames à Moffett Field, California, USA. La section de la soufflerie est 24.4 m × 36.6 m.
- L’éolienne NREL phase VI : cette éolienne a deux pales de diamètre de 10.05 m, la

corde n’est pas constante, et les pales sont vrillées. La vitesse de rotation pendant l’essai
est fixée à 71.63 tours/min, le rotor est fixé sur une nacelle, et cet ensemble est situé à
une hauteur de 12.192 m [65].
- L’éolienne Rutland 503 : cette éolienne a trois pales de diamètre égal à 0.5 m,
l’angle de calage est constant le longe de la pale est égale à 10°, la corde de la pale varie
linéairement entre 0.065 m au pied de la pale jusqu’à 0.045 m à l’extrémité, le rayon du
moyeu est égale à 0.065 m, la vitesse de rotation appliquée pendant tous les essais est
égale à 1050 tours/min. Les essais pour cette éolienne sont faits dans la veine d’essai de
la soufflerie de l’ENSAM-paris. Cette veine est semi-guidée avec une section de 1.65 m ×
1.35 m, et une longueur de 2 m.
74
3.11.2 Le domaine de calcul
Tout d’abord, le modèle hybride basé sur le concept de disque actif, qui représente
l’éolienne NREL phase II, est de forme cylindrique de diamètre de 10.1 m, et d’épaisseur
de 0.09522 m ( c × sinβ, figure 3.29, page 63. Ce modèle hybride contient 360000 mailles,
situé à 5 D de l’entrée principale et représente le rotor, figure 3.15.
Le domaine de calcul est constitué de six blocs cylindriques de plus en plus grands
autour du modèle hybride, ceci afin d’obtenir un maillage assez fin autour de celui-ci et
de plus en plus grand en s’éloignant, avec un nombre de mailles optimale.
Le bloc 1 à un diamètre de 1.5D, et une longueur de 0.25D, ce bloc englobe le modèle
hybride avec un nombre de mailles égal à 7900 mailles. Le bloc 2 englobe le bloc 1, son
diamètre extérieur est égal à 2D et de longueur de 0.5D, celui-ci a 11526 mailles. Le bloc
3, englobe le bloc 2, son diamètre extérieur est égal à 3D et de longueur de 1D, celui-ci
a 39250 mailles. Le bloc 4, englobe le bloc 3, son diamètre extérieur est égal à 5D et
de longueur de 2D, celui-ci a 217594 mailles. Le bloc 5, englobe le bloc 4, son diamètre
extérieur est égal à 7D et de longueur de 7.5D, celui-ci a 232768 mailles. Le bloc 6, englobe
le bloc5, son diamètre extérieur est égal à 10D et de longueur de 20, celui-ci a 314606
mailles. Donc le total des mailles du domaine de calcul est de 1180000.
Bloc 1
Entrée
Bloc 2
Bloc 3
Modèle Bloc 4 Sortie
10 D
hybride Bloc 5
Bloc 6
20 D
Fig. 3.15 – Domaine de calcul pour le modèle hybride basé sur le concept du disque actif.
75
Pour l’éolienne NREL phase VI, le maillage est très semblable au modèle précédent,
l’épaisseur du volume cylindrique qui remplace le rotor est égal 0.0364 m, le nombre total
des mailles dans ce cas est égal à 360000 mailles. Dans ce cas, le bloc 1 reste avec le même
diamètre extérieur de 1.5D mais avec un nombre de mailles égal à 8956. Le bloc 2, reste
avec les mêmes dimensions mais avec un nombre de mailles de 14890. Par contre, pour
les autres blocs, 3, 4, 5 et 6 aucun changement, et le nombre total de mailles est égal à
1120000.
76
3.11.3 Comparaison des résultats numériques et expérimentaux
On va commencer nos comparaisons avec l’éolienne NREL phase II puis avec l’éolienne
NREL phase VI.
Fig. 3.16 – Comparaison des puissances entre le modèle hybride DA et les données
expérimentales, NREL phase II.
Pour l’éolienne NREL phase II, les données expérimentales obtenues par le TU-Delft
[71] sont comparées aux résultats du modèle hybride, voir la figure 3.16.
Les résultats de puissance en fonction de la vitesse de vent du modèle hybride sont
proches des données expérimentales. Le résultat de notre modèle hybride est très semblable
aux résultats trouvés par YawDyn qui a utilisé la méthode BEM dans son calcul [81].
On utilise les données expérimentales de l’éolienne NREL phase VI [70], et on les
compare aux résultats obtenus par le modèle hybride basé sur le concept du disque actif,
figure 3.17. On remarque à partir de cette figure que, le modèle hybride étudié donne
des résultats proches des données expérimentales jusqu’à la vitesse de vent environ 10
m/s, après cette valeur de la vitesse de vent, le modèle hybride donne des résultats qui
s’éloignent des données expérimentales car on arrive dans la zone de décrochage.
77
Fig. 3.17 – Comparaison des puissances entre le modèle hybride DA et les données
expérimentales, NREL phase VI .
3.12 Comparaison entre le sillage mesuré d’une éolienne

et le sillage calculé par le modèle hybride
Nous allons comparer deux champs de vitesse axiale, le première a été mesuré en aval
de l’éolienne Rutland 503. Cette exploration a été faite dans la veine d’essai de la soufflerie
de l’ENSAM-Paris. Le deuxième champ comparé a été calculé par le modèle hybride DA
pour la même machine.
La figure 3.18, représente le champ de la vitesse axiale moyenne derrière le rotor, et
sur une surface de longueur allant jusqu’à Z/D=2.1 à partir de Z/D=0.25 du plan de
rotation, et de hauteur allant jusqu’à Z/D=0.9 de l’axe de rotation. On a utilisé la même
surface correspondante à celle-ci dans le domaine de calcul derrière le modèle hybride,
figure 3.19. La différence entre les deux champs de vitesse est donnée par la figure 3.20.
78
Fig. 3.18 – Champ de vitesse axiale moyennée expérimentale mesurée dans la soufflerie
de l’ENSAM-Paris.
Fig. 3.19 – Champ de vitesse axiale obtenu par le modèle hybride basé sur le concept du
disque actif.
79
Fig. 3.20 – Différence entre les deux champs des vitesses axiales, modèle hybride et Essai.
A partir de la figure 3.20, on observe que la différence entre les deux champs de vitesse
axiale est faible et que la majorité de cette erreur est assez proche de 0. En revanche, il y
a des zones où la différence est assez nette, spécialement près du rotor, où cette différence
est due à l’absence de la géométrie réelle des pales dans le cas du modèle hybride.
3.13 Conclusion
Dans ce travail, nous avons présenté un modèle hybride qui combine un rotor simplifié
avec un solveur des équations de Navier Stokes tridimensionnel. On l’a validé sur trois
types d’éoliennes : éoliennes Rutland 503, NREL phase II et NREL phase VI.
Le rotor est représenté par un disque de diamètre égal au diamètre du rotor et
d’épaisseur liée à la corde de la pale de l’éolienne. La validation de ce modèle hybride
est faite de deux façons différentes.
-Pour les performances de l’éolienne en fonction de la vitesse du vent : cette validation
est faite pour l’éolienne NREL phase II et VI. Pour l’éolienne NREL phase II, on trouve
qu’il y a une corrélation forte entre les données expérimentales et les résultats obtenus par
le modèle hybride pour toutes les valeurs des vitesses de vent utilisées. Par contre, pour
l’éolienne NREL phase VI, les résultats obtenus par le modèle hybride sont assez proches
80
des données expérimentales jusqu’à V1 = 10 m/s, au delà de cette valeur, le modèle
hybride donne des résultats différents des données expérimentales à cause du décrochage
mal pris en compte par le modèle DA.
-Pour le champ de vitesse axiale dans le sillage : la validation est faite entre le champ
de vitesse axiale obtenu par le modèle hybride et le champ mesuré dans la soufflerie de
l’ENSAM-Paris en aval de l’éolienne Rutland 503. Le résultat montre que les deux champs
sont assez proches. La figure 3.20 montre que la majorité des différences entre les deux
champs est proche de zéro. Les différences notables sont près du rotor, car la géométrie
réelle de la pale est complètement absente.
En conclusion générale, les résultats obtenus par le modèle hybride basé sur le concept
du disque actif ont été présentés dans ce chapitre, ils sont comparables avec les données
expérimentales concernant les trois cas d’éoliennes présentées dans ce travail. La validation
de ce modèle, avec les restrictions qui s’y appliquent, est considérée comme acquise.
81
82
Chapitre 4
Développement d’un modèle hybride

basé sur le concept de la ligne active
83
4.1 Introduction
Le chapitre précédent a présenté un modèle hybride basé sur le concept du disque actif.
Les forces axiales et tangentielles dans le disque actif sont distribuées radialement avec une
intensité variable mais la distribution dans la direction azimutale reste uniforme et cette
distribution ne permet donc pas la représentation individuelle des pales. Pour surmonter
cette limitation, et obtenir un sillage plus proche de la réalité, Sørensen et Shen [42] ont
développé une représentation tenant compte du nombre des pales. Dans leur modèle, la
géométrie réelle des pales est remplacée dans le solver RANS1 par des forces volumiques
distribuées le long de lignes actives représentant les pales. Leurs résultats ont montré un
bon accord avec la courbe de puissance mesurée pour l’éolienne Nordtank tri pales 500
kW [42]. Par la suite, Mikkelsen [43] a utilisé cette approche de ligne active avec une
distribution des forces dans des plans normaux à chaque ligne active selon une répartition
Gaussienne. Ensuite, Ivanell [44] a utilisé la même approche que Mikkelsen et Sørensen et a
montré l’influence des paramètres qui contrôlent la distribution de Gauss sur les résultats
du calcul. Il a aussi montré qu’on peut améliorer les résultats avec l’utilisation de la
correction de Prandtl. Le modèle hybride présenté dans ce chapitre est basé sur le concept
de la ligne active, mais chaque pale est remplacée par un volume cylindrique (Cylindre
Actif- CA). Les forces exercées par un élément de pale sur le fluide sont distribuées de
façon uniforme à l’intérieur de l’élément de cylindre correspondant.
La validation de ce modèle hybride est faite de deux façons différentes. D’abord par
rapport à l’évolution de la puissance en fonction de la vitesse du vent où des comparaisons
sont faites avec l’éolienne NREL phase II [51] et l’éolienne NREL phase VI [65]. Ensuite,
une autre comparaison est faite entre le sillage simulé par ce modèle hybride et le sillage
de l’éolienne Rutland 503 testée dans la soufflerie de l’ENSAM-Paris.
4.2 Développement d’un modèle de cylindre actif

Dans ce modèle hybride, une combinaison est faite entre la méthode de l’élément de
pale et les équations de Navier-Stokes. Chaque pale est remplacée par un cylindre actif
CA de longueur égale à la longueur de la pale, figure 4.1, et d’un diamètre qui sera défini
suite à une étude d’optimisation spécifique.
L’approche de l’élément de pale sera utilisée pour calculer les efforts exercés par la
pale sur le fluide à l’aide d’un sous programme intégré au sein d’un logiciel de simu-
lation numérique. Pour définir ces efforts, les caractéristiques aérodynamiques bidimen-
1
RANS : Reynolds-averaged Navier-Stokes.
84
y n
dr
θ x
Ω
z t
Fig. 4.1 – Modèle de Cylindre Actif représentant une éolienne tripale.
sionnelle des profils des pales sont utilisées conjointement avec le vent incident obtenu
numériquement.
Les efforts appliqués par l’élément de pale seront répartis uniformément à l’intérieur
du segment de cylindre correspondant à cet élément. Il s’agit donc de forces volumiques
réparties sur des cellules de calcul et représentant un terme source de quantité de mou-
vement.
85
n
dF
dFz
φ dFx
Plan de rotation φ t
β
φ α
W
Fig. 4.2 – Forces exercées par l’écoulement sur un élément de pale.
La figure 4.2 représente les forces aérodynamiques agissant sur un élément de pale dr
situé à la distance r de l’axe de rotation.
Le vecteur W , situé dans le plan du profil, représente la vitesse de l’air relative au
profil après prise en compte de la vitesse induite.
A partir de la figure 4.2 on trouve que les forces élémentaires de traı̂née et de portance
dFx et dFz respectivement, appliquées sur l’élément de pale, sont :
dFx = 0.5 ρ W 2 c dr Cx

(4.1)
dFz = 0.5 ρ W 2 c dr Cz

(4.2)
Où Cx et Cz sont les coefficients de traı̂née et de portance du profil.
Pour la suite, nous définissons un référentiel t, n, r lié à la pale tournante. r représente
l’axe de la pale qui tourne avec l’axe t autour l’axe fixe n.
La projection des forces élémentaires dFx et dFz sur l’axe t (dans le plan de rotation),
donne la force élémentaire tangentielle dFt appliquée sur l’élément de pale dr.
dFt = dFx cosφ − dFz sinφ = 0.5 ρ W 2 c dr (Cx cosφ − Cz sinφ)

(4.3)
86
La projection des forces élémentaires dFx et dFz dans la direction normale au plan de
rotation selon l’axe n, donne la force élémentaire normale dFn appliquée sur l’élément de
pale dr.
dFn = dFx sinφ + dFz cosφ = 0.5 ρ W 2 c dr (Cx sinφ + Cz cosφ)

(4.4)
Pour dFt et dFn , on définit les coefficients Ct et Cn tels que :
dFt = 0.5 ρ W 2 c dr Ct

(4.5)
dFn = 0.5 ρ W 2 c dr Cn

(4.6)
A partir des équations (4.3) et (4.5) on trouve que le coefficient de la force tangentielle
dans le plan de rotation, est :
Ct = Cx cosφ − Cz sinφ (4.7)

A partir des équations (4.4) et (4.6) on trouve que le coefficient de la force normale
au plan de rotation, est :
Cn = Cx sinφ + Cz cosφ (4.8)

Précisons que, dFt et dFn sont des forces appliquées sur le profil par l’écoulement. Par
contre, les forces appliquées sur le fluide par le profil sont −dFt et −dFn .
4.3 Répartition des efforts dans un cylindre actif

Dans le modèle étudié ici, on représente chaque pale par un cylindre actif, figure 4.1.
Les forces créées par la pale sont distribuées comme des forces volumiques à l’intérieur
d’un cylindre actif qui remplace la pale. Selon la méthode de l’élément de pale, on divise
chaque pale et chaque cylindre actif en n éléments.
D’après les propriétés aérodynamiques du profil de pale, on trouve qu’un élément de
pale d’épaisseur radiale dr situé au rayon r reçoit une force élémentaire de la part de
l’écoulement égale à :
1
dF~ = ρ c CF W 2 dr~e (4.9)
2
CF : Coefficient de force appliquée au profil.
~e : Vecteur unitaire définissant la direction de F~ .
87
Cette force sera supposée repartie uniformément à l’intérieur de l’élément dv du cy-
lindre actif :
πd2
dv = dr (4.10)
4
d : Diamètre du cylindre actif.
En divisant l’équation (4.9) par l’équation (4.10), on trouve que l’intensité de force
par unité du volume est :
dF~ 2 ρ c CF W 2
f~ = = ~e (4.11)
dv π d2
Pour faciliter la communication avec le logiciel de simulation de l’écoulement, on ana-
lyse dF~ dans le référentiel x, y, z utilisé dans ce dernier, voir la figure 4.1.
Par projection de f sur x, y, z, on trouve que l’intensité de la force appliquée par
l’écoulement sur un volume élémentaire dv, par rapport aux axes x, y et z est :
Dans la direction x :
2 ρ c Ct W 2
fx = sin θ (4.12)
π d2
Dans la direction y :
2 ρ c Ct W 2
fy = cos θ (4.13)
π d2
Dans la direction z :
2 ρ c Cn W 2
fz = − (4.14)
π d2
Ces forces sont appliquées par l’écoulement sur un élément de pale.
4.3.1 Choix du diamètre du cylindre actif

Pour utiliser l’approche du cylindre actif proposée dans cette thèse, la première ques-
tion à régler est celle concernant le choix du diamètre du cylindre actif qui est capable de
bien représenter une pale réelle. Pour aider à faire ce choix, une étude préparatoire est faite.
Elle consiste à comparer les résultats de simulation numérique pour le champ de vitesse
créé par le cylindre actif et celui créé par la pale réelle. Pour y arriver, une étude spécifique
par simulation numérique a été faite en 2D entre le profil s809 de l’éolienne NREL phase
II, figure 4.3, et le profil du cylindre actif constitué par un cercle de diamètre d, figure
4.4.
88
Domaine de calcul
La figure 4.3 représente le domaine de calcul autour du profil s809, celui-ci a une
longueur égale à 120 c (c : la corde de la pale) et une largeur de 80 c. Le profil s809 est
situé à 40 c de l’entrée, ce profil est divisé en deux parties, extrados et intrados, chaque
partie est maillée avec 125 mailles resserrées au bord de fuite et au bord d’attaque. Un
maillage de couche limite est plaqué sur l’intrados et l’extrados du profil, l’épaisseur de
la première maille de cette couche limite est égale à 0.00065 c de manière à obtenir des
valeurs de y + de l’ordre de 30. La taille des mailles augmente en s’éloignant de l’extrados
et l’intrados selon une progression géométrique de raison 1.15. Les lignes AB et BF sont
maillées avec 125 mailles, les lignes AJ, BC et FD sont maillées avec 250 mailles chacune.
Les trois surfaces BFDC, BCJA et FOA sont maillées avec 31250 mailles chacune. Le
maillage utilisé ici est structuré et quadrangulaire.
Fig. 4.3 – Domaine de calcul autour du profil s809.
89
A
B J
D
Fig. 4.4 – Domaine de calcul autour d’un cercle équivalent au profil s809.
La figure 4.4, représente le domaine de calcul autour d’un cercle destiné à remplacer
l’influence du profil s809. Afin de pouvoir comparer les solutions numériques produites
dans ces deux situations, on génère des maillages ayant le même nombre de mailles et la
même taille caractéristique. Le domaine maillé est un cercle de diamètre de 40 c. Les lignes
AB et BC sont maillées par 125 mailles. Les lignes ADC et AJC sont maillées avec 250
mailles chacune. Les deux surfaces identiques ABCD et ABCJ sont maillées avec 31250
mailles chaque avec un maillage structuré quadrangulaire.
Validation de maillage
Pour justifier la construction des maillages présenter précédemment sur le domaine

de calcul, nous allons calculer les coefficients de traı̂née et portance du profil s809 en
fonction de l’angle d’incidence, et comparer les résultats obtenus avec des autre obtenus
par l’université de Dlft2 et l’université de OSU3 [82].
2
Delft University of Technology.
3
Ohio State University.
90
Tous les calculs sont réalisés à l’aide du logiciel Fluent version 6.2.18, le modèle de tur-
bulence choisi dans tous les cas est k − w − SST , car ce modèle de turbulence permet une
prévision précise des écoulements aérodynamiques avec de forts gradients opposés de pres-
sion [83], ce modèle est plus robuste que les autres modèles, même pour des applications
complexes [84, 85].
Fig. 4.5 – Coefficients de portance pour le profil s809.
Le nombre de Reynolds utilisé dans cette étude est égal à 106 , ce choix correspond à
la valeur utilisée par le NREL, l’université de Delft et l’université d’Ohio pour calculer
les coefficients de traı̂née et de portance du profil s809 [86, 52, 82]. Nous avons calculé les
coefficients de portance et de traı̂née du profil s809 par un code de simulation numérique,
les résultats sont présentés sur les figures (4.5 et 4.6).
A partir de ces deux figures, on observe que les résultats de notre calcul sont assez
proches des résultats données par OSU et par Delft pour les deux coefficients, portance
et traı̂née. Nous pouvons donc valider nos paramètres de calcul.
Calcul des champs de vitesse axiale autour du profil s809 et du cercle équivalent
Les simulations effectuées pour trouver le cylindre actif équivalent à la pale réelle sont
faites pour quatre angles de calages différents, 0, 6, 12 et 18 degrés, et plusieurs valeurs
d/c du diamètre du cercle entre 0.05 et l.
Dans notre calcul, nous avons calculé les forces agissant sur le profil s809, ensuite,
nous avons appliqué ces forces sur le disque qui remplace le profil. Nous présentons les
91
Fig. 4.6 – Coefficients de traı̂née pour le profil s809.
comparaisons sur des domaines de calcul de forme carrée de longueur égale à six cordes
pour bien montrer le comportement du sillage.
Fig. 4.7 – Champ de vitesse axiale autour du profil s809, l’angle de calage est β = 6°.
92
La (Fig.4.7) présente le champ de vitesse axiale autour du profil s809 avec un angle
de calage égal à 6°. On observe bien l’augmentation de la vitesse axiale sur l’extrados du
profil et la diminution de cette vitesse sur l’intrados. Les variations importantes du champ
de vitesse se concentrent clairement dans une région proche du profil et dans son sillage.
Des comparaisons ont été faites pour 13 valeurs de d/c et quatre angles d’incidence,
soit 52 cas de calcul (voir les annexes B, C, D et E). Nous allons présenter trois cas,
d = 0.05 c, 0.4 c et 1 c.
Fig. 4.8 – Champ de vitesse axiale autour d’un cercle de diamètre de 0.05 c, β = 6°.
A partir des figures 4.8, 4.9 et 4.10, on observe qualitativement les mêmes variations
de vitesse autour des cercles qu’autour du profil s809. Afin de mieux quantifier les écarts
entre les vitesses produites par le profil s809 et ses cercles équivalents, nous allons tracer
des profils de vitesse à une corde en amont du centre de chaque cercle. C’est en effet dans
cette zone que sont prises les vitesses de référence utilisées dans les différents modèles de
rotors équivalents.
La figure 4.11 présente la vitesse axiale extraite du calcul complet autour des cercles
équivalents comparée à la vitesse axiale du profil s809, sur la figure 4.12 nous avons tracé la
vitesse tangentielle. Ces résultats montrent que les profils de vitesse à une corde à l’amont
du centre des cercles sont très peu sensibles au diamètre du cercle équivalent utilisé si
bien que la puissance globale de l’éolienne sera la même pour tous les cas considérés. En
revanche, le sillage et le champ proche sont fortement dépendants du diamètre utilisé.
93
Fig. 4.9 – Champ de vitesse axiale autour d’un cercle de diamètre de 0.4 c, β = 6°.
Fig. 4.10 – Champ de vitesse axiale autour d’un cercle de diamètre de 1 c, β = 6°.
Des valeurs de d = 0.3 c à d = 0.6 c donnent des résultats corrects. Le choix d’un
diamètre permettant d’obtenir la même surface que celle du profil s809 nous semble donc
94
Fig. 4.11 – Comparaison entre des profils de vitesse axiale.
Fig. 4.12 – Comparaison entre des profils de vitesse tangentielle.
une bonne solution pour le choix du diamètre du cercle équivalent. Ce critère nous donne
ici un diamètre d = 0.4 c qui sera utilisé par la suite dans toutes nos simulations.
95
4.4 Application du modèle de cylindre actif aux cas
NREL s809 phase II et VI
La validation du modèle hybride CA est faite de deux manières différentes, la première
concerne l’évolution de la puissance en fonction de la vitesse à l’infini amont, cette valida-
tion est faite pour deux types d’éoliennes, NREL phase II et VI. La deuxième validation
a été faite sur le champ de vitesse axiale de l’éolienne Rutland 503. Tout d’abord, nous
allons présenter l’algorithme de calcul utilisé dans ce modèle hybride.
4.4.1 Algorithme de calcul

Dans la construction de cet algorithme, les équations (4.7), (4.8), (4.12), (4.13) et
(4.14) sont principalement utilisées. Pour chaque maille appartenant au cylindre actif, on
effectue la boucle présentée sur la figure (Fig.4.13).
Ainsi, à partir de cet algorithme, on peut calculer :
- l’effort de traı̂née dans chaque cellule.
- l’effort tangentiel qui s’exerce par chaque cellule.
- le couple mécanique créé autour de l’axe de rotation du rotor par chaque cellule.
Le calcul est fait à l’aide d’un programme écrit en langage C et intégré dans le logiciel
Fluent. Ce programme calcule la force volumique au centre de chaque maille. Chaque
maille dans le cylindre actif donne une force volumique différente en fonction du volume
de la maille et de la distance de cette maille par rapport à l’axe de rotation r.
Ces calculs sont réalisés à chaque itération de la simulation numérique. Pour calculer les
efforts appliquées par une maille du cylindre actif sur l’écoulement, on multiple l’intensité
de l’effort f par le volume de la maille considérée, puis la somme des forces volumiques
de toutes les mailles, donne les forces globales.
L’effort dans la direction axiale est :
nmailles
X
Fn = (fz,i · V olumemaille,i ) (4.15)
i=1
L’effort tangentiel est calculé comme suit :

nmailles
X q
2 2
Ft = ( fx,i + fy,i · V olumemaille,i ) (4.16)
i=1
Et le couple vaut :
96
Lecteur des coordonnées du centre de la maille. xcentre ycentre zcentre
Calcul de rayon du centre de la maille. r=√x²centre + y²centre
Calcul de l’azimut du centre de la maille. tanθ = ycentre / x centre
Lecture de la vitesse dans le plan Vx , V y , Vz .

de réference.
Calcul de la vitesse tangentielle. V tan = Vy cos θ - Vx sin θ
Calcul de la vitesse relative. W = √ V²z +( Ω r + V tan)²
Calcul de l’angle d’écoulement relatif φ et tanφ = Vz / ( Ω r + Vtan )

l’angle d’incidence α . α=φ−β
Calcul de la traînée et de la portance du profil. Cx = f (α, Re) , Cz = f (α, Re)

Cx et Cz
Calcul des coefficients C t et C n Cn = Cz cos φ + Cx sin φ

C t = Cx cos φ − C z sin φ
Calcul l’intensité des forces volumiques 2 ρ c Ct W²

f x= - sin θ
dans les mailles du cylindre actif. л d2
2 ρ c Ct W²
fy = - cos θ
л d2
2 ρ c Cn W²
fz =
л d2
Fig. 4.13 – Algorithme de calcul pour le modèle hybride CA.
97
nmailles
X q
2 2
Couplemeca = ( fx,i + fy,i · ri · V olumemaille,i ) (4.17)
i=1
La puissance mécanique se calcule par la multiplication du couple mécanique par la

vitesse de rotation du rotor Ω :
Pmeca = Couplemeca · Ω (4.18)

La puissance électrique est donnée par une équation caractéristique de l’éolienne NREL
phase II [51] :
Pelec = 0.9036 · Pmeca − 0.847 (4.19)

Avec la puissance électrique, on dispose d’une valeur que l’on va pouvoir comparer
aux résultats expérimentaux du NREL phase II et donc discuter de la validité du modèle
proposé.
Pour l’éolienne NREL phase VI, on compare directement la puissance mécanique.
98
4.4.2 Surface de référence
Cette surface de référence nous permet de calculer les vitesses relatives et incidentes
en évitant l’effet local des forces volumiques créées par le cylindre actif sur l’écoulement.
Après plusieurs essais, et plusieurs visualisations des champs des vitesses, nous avons
trouvé que la meilleure position pour cette surface est à une distance égale à une corde
de pale en amont du plan de rotation.
En raison de la périodicité du calcul (éolienne possédant trois pales), nous avons pris
un tiers de cylindre correspondant à un tiers du domaine global d’étude. On divise la
surface de référence en n tiers d’anneaux de largeur radiale dr, on calcule les vitesses
moyennées sur chaque tiers d’anneau. La vitesse axiale moyennée calculée sur un tiers
d’anneau situé au rayon r, est utilisée par la théorie de l’élément de pale pour calculer les
forces agissant sur un élément de CA à partir de l’élément de pale correspondant, figure
4.14.
Surface de référence Cylindre actif
dr
V1
c
Z
Fig. 4.14 – Surface de référence dans le cas d’une éolienne tripales.
99
4.4.3 Domaine de calcul
Le domaine de calcul est périodique, nous avons divisé un domaine global de forme
cylindrique de longueur 10D et de diamètre de 10D.
Dans le cas de l’éolienne NREL phase II, on prend un tiers du domaine global de
diamètre égal à 10D et de longueur 10D, figure 4.15. Le domaine est divisé en cinq parties,
figure 4.16, la première de forme cylindrique de diamètre égal à la corde de la pale, l’axe
de ce cylindre est le même que celui du cylindre actif (CA) qui remplace la pale, cette
partie englobe le CA, le maillage de cette partie est très fin pour bien détecter le champ
de vitesse autour du CA, avec un nombre de mailles égal à 39709.
Fig. 4.15 – Domaine de calcul.
La deuxième partie du domaine a une forme cylindrique de longueur égale à 2 cordes,

parallèle à la direction de l’axe de rotation, de diamètre égal à 1.2D, la densité du maillage
de cette partie est plus faible que pour la première, le nombre de mailles est de 42055.
La troisième partie est de la même forme que la deuxième, mais avec une longueur de
10 cordes et un diamètre égal à 1.6D avec un nombre de mailles de 60663. La quatrième
partie a la forme d’un grand cylindre de diamètre 2D, et de longueur 1.8D, le nombre
de mailles de ce bloc est égal à 22078. La cinquième partie a une longueur égale à 10D
et un diamètre égal à 10D, cette partie est maillée avec une faible densité de maillage,
le nombre de mailles est égal à 300164. Le cylindre actif qui remplace la pale est placé à
100
Bloc 1 5D
Bloc 2
Bloc 3
Bloc 4
Bloc 5
10D
Fig. 4.16 – Section longitudinale du domaine de calcul montrant les cinq blocs autour du
cylindre actif.
trois diamètres de l’entrée du domaine, le nombre de mailles du cylindre actif est égal à
20695 mailles, figure 4.16. Donc, le nombre total de mailles est égal à 482504 mailles. La
densité du maillage diminue au fur et à mesure que l’on s’éloigne du cylindre actif.
Dans le cas de l’éolienne NREL phase VI, le domaine de calcul a les mêmes dimensions
que dans le cas précédent avec la même structure des blocs, 10D de longueur et 10D de
diamètre et 5 blocs. L’éolienne NREL phase VI a deux pales, donc on va prendre dans le
calcul périodique un demi-domaine, le diamètre du cylindre actif qui remplace la pale est
égal à 0.4 fois la corde moyenne de la pale, parce que la corde n’est pas constante le long
de la pale.
Le maillage est analogue au cas précédent et compte 700000 cellules.
4.5 Résultats et analyse

Le modèle hybride CA est introduit dans le logiciel Fluent 6.2.16. Tous les résultats
sont obtenus dans le cas d’un écoulement incompressible et stationnaire, utilisant le modèle
de turbulence k − w − SST .
Le modèle hybride est appliqué au cas des éoliennes NREL phase II et VI. La figure
4.17, représente les données expérimentales et les résultats obtenus par le modèle hybride
pour la puissance de l’éolienne NREL phase II. On trouve qu’il y a une bonne corrélation
pour les vitesses de vent faibles, jusqu’à 12m/s, au-delà de cette valeur, les résultats du
101
modèle hybride s’éloignent des données expérimentales jusqu’à la valeur de 20m/s. Les
résultats de notre modèle hybride est très semblable aux résultats du modèle de YawDyn
[81].
Fig. 4.17 – Comparaison des puissances, modèle hybride CA, NREL phase II.
Par contre, la figure 4.18 représente la comparaison entre les résultats obtenus pour
le modèle hybride et les données expérimentales de l’éolienne NREL phase VI. Sur cette
figure, on trouve qu’il y a une bonne corrélation des résultats jusqu’à la vitesse du vent
de 9 m/s, après cette valeur, on arrive dans la zone de décrochage et les résultats sont
nettement plus éloignés.
102
Fig. 4.18 – Comparaison des puissances, modèle hybride CA, NREL phase VI.
4.6 Comparaison simulation/ expérience du sillage

d’une éolienne Rutland 503 modifiée
Dans ce paragraphe, nous allons présenter la comparaison entre le sillage mesuré en
aval d’une éolienne Rutland 503 tripale et le sillage calculé derrière un modèle hybride
CA qui remplace celle-ci. La machine a été testée dans la veine d’essai de la soufflerie de
l’ENSAM-Paris. Le calcul du modèle hybride a été induit dans le Fluent version 6.2.16 et
le modèle de turbulence est k − w − SST . Dans les deux cas, la vitesse à l’infini amont
égale à 9.4 m/s.
La figure 4.19, représente le champ de vitesse axiale moyenne derrière le rotor. On
représente la vitesse axiale produite par le modèle CA dans le même domaine géométrique,
figure 4.20.
103
Fig. 4.19 – Champ de la vitesse axiale moyenne expérimentale en aval d’un éolien Rutland
503 modifiée.
Fig. 4.20 – Champ de vitesse axiale obtenu par le modèle hybride CA.
Sur la figure 4.21, nous présentons la différence entre les deux champs de vitesse axiale,
mesurée et calculée, nous observons qu’il des zones où la différence est proche de zéro.
Ainsi, il y a des zones où la différence arrive jusqu’à une valeur de 5 m/s à côté du rotor,
où cette différence est due à l’absence de la géométrie réelle des pales dans le cas du
modèle hybride.
104
Fig. 4.21 – Différence entre les deux champs des vitesses axiales, modèle hybride et Essai.
4.7 Conclusion
Dans le travail présenté dans ce chapitre, un modèle hybride de cylindre actif est décrit
pour évaluer les performances aérodynamiques d’une éolienne. Une combinaison entre la
théorie de l’élément de pale et un code de Navier-Stokes tridimensionnelle est utilisée
par ce modèle hybride. L’utilisation de la théorie d’élément de pale permet de calculer
d’une part, les forces aérodynamiques appliquées sur l’écoulement par le cylindre actif
qui remplace la pale. D’autre part, la résolution itérative des équations de Navier-Stokes
tridimensionnelle donne le champ de vitesses en amont et en aval du rotor.
Pour représenter le rotor, chaque pale est remplacée par un cylindre actif de longueur
égale à celle de la pale et de diamètre égal à 0.4 fois la corde de pale.
La validation de ce modèle hybride a été faite de deux façons différentes :
- L’évolution de la puissance en fonction de la vitesse à l’infini amont a été comparée
sur deux l’éoliennes. Pour l’éolienne NREL phase II, on observe que le modèle hybride
proposé donne des résultats raisonnables par rapport aux résultats expérimentaux. Le cas
de NREL phase IV, est plus difficile à représenter mais les résultats restent satisfaisant.
- Concernant les sillages en aval le modèle CA et l’éolienne remplacée, une comparaison
a été faite entre le champ de la vitesse axiale calculé en aval le modèle hybride et celui
mesure en aval de la machine testée. La résultat est présenté sur la figure 4.21, on observe
qu’il y a un bon accord entre les deux champs malgré la différence près du pied de la pale.
105
Le modèle hybride basé sur le concépt de la ligne active (CA) fournit des résultats
qualitativement et quantitativement correct sur le cas de validation proposés. L’effort
pour produire de tels résultats numériques est plus important que dans le cas du modèle
CA : le maillage est en effet plus difficile à réaliser.
106
107
Chapitre 5
Interactions entre des éoliennes
108
5.1 Introduction
Le comportement du sillage d’une éolienne est une question très importante qu’il faut
prendre en compte lors de la construction d’un parc éolien. En général, dans le sillage
d’une éolienne, la vitesse du vent diminue et le niveau de turbulence augmente. Dans ces
conditions, si une éolienne se trouve dans le sillage d’une autre, elle recevra un vent qui
a déjà cédé une bonne partie de son énergie et sa production d’énergie sera affectée. De
plus, les phénomènes instationnaires existants dans le sillage sont à l’origine de charges
dynamiques qui peuvent à long terme accélérer la fatigue des matériaux et réduire la durée
de vie des pales et des autres composants sollicités [87, 88, 89].
Pour éviter ces problèmes, il faut prévoir une distance minimale entre les éoliennes.
Cependant, la valeur financière du terrain et le besoin des infrastructures dans un parc
éolien tendent à rapprocher les machines le plus possible. Alors, pour déterminer cette
distance minimale entre les éoliennes et ainsi augmenter la rentabilité d’un parc éolien, il
est important d’étudier et de prendre en compte les effets des éoliennes les unes sur les
autres et leur contribution comme éléments de rugosité dans le parc éolien [90, 91, 92].
Les premières fermes offshore furent construites à Helgoland en Allemagne en 1989, à
Blekinge en Suède en 1990 et à Vindeby au Danemark en 1991 [93]. Dans ces cas, les
effets du sillage (diminution de la vitesse et augmentation de la turbulence en aval des
machines) se propagent sur de plus grandes distances en aval de l’éolienne par rapport au
cas d’un parc terrestre, ceci étant dû à la plus faible rugosité de la surface de la mer par
rapport à celle de la surface de la terre [94].
Les paramètres les plus importants dans l’étude de l’interaction entre plusieurs ma-
chines sont : la disposition relative de chaque machine et la vitesse du vent à l’infini
amont. A un degré moindre, on pourrait citer également l’intensité de la turbulence de
l’écoulement incident ainsi que la stabilité atmosphérique (naturel, stable, instable)1 [95].
Dans cette partie de la thèse, une étude de l’interaction entre deux éoliennes identiques
situées l’une derrière l’autre va être présentée. Nous montrerons l’évolution de la puissance
de deux machines en fonction de la distance entre celles-ci et en fonction de la vitesse du
vent amont. Nous avons vu que le modèle DA fournissait de bon résultats de port sur
simplicité, le modèle DA est donc à privilégier. Cependant, il nous semble intéressant
d’effectuer l’étude des deux modèles. Cette étude utilisera donc les deux modèles hybrides
1
Naturel : les particules d’air restent à leurs nouvelles places après la disparition des forces agissant
sur les particules. Stable : chaque particule va revenir à sa position initiale après la disparition des forces
agissant sur les particules. Instable : chaque particule va continuer son chemin même après la disparition
des forces agissent sur les particules.
109
développés dans les chapitres précédents. Enfin, nous présenterons le comportement de
plusieurs machines situées l’une en aval de l’autre. Les calculs seront effectués en supposant
la présence d’éoliennes de type NREL phase II.
5.2 Simulation complète de l’interaction entre deux

éoliennes à l’aide d’un modèle hybride
La construction d’un parc éolien est souvent basée sur des résultats de calcul d’un
sillage seul. Normalement, pour prendre en compte les effets combinés des différents
sillages, l’hypothèse de superposition qui permet de cumuler l’effet des différents sillages
est la plus utilisée [96, 97]. Lissaman [98] est le premier à avoir utilisé la ”superposition
linéaire des perturbations” créées par le sillage des différentes machines présentes sur le
parc. Cependant, cette méthode échoue à la détermination de grandes perturbations car
elle surestime les déficits de vitesse, ce qui pourrait aboutir à un résultat absurde i.e.
l’apparition de vitesses négatives lors de la superposition de nombreux sillages. Katic et
al [99] ont utilisé le principe de superposition linéaire des zones de déficit de vitesse. Dans
ce cas, l’effet cumulatif calculé du sillage de plusieurs machines s’avère être plus faible par
rapport à la méthode utilisée par Lissaman. La méthode de Katic permet ainsi d’obtenir
des résultats satisfaisants en comparaison avec des données expérimentales.
Dans notre travail, nous adoptons une méthode plus directe que les méthodes précédentes,
nous utilisons un domaine du calcul qui comporte deux modèles hybrides représentant
les deux machines en même temps. Ceci donne une simulation numérique complète de
l’interaction éolienne-éolienne à l’aide des modèles hybrides en évitant les problèmes de
superposition.
Ce travail a été fait avec le logiciel Fluent version 6.2.16 avec les modèles DA et CA.
Tous les résultats sont obtenus dans le cas d’un écoulement stationnaire et incompressible
avec un modèle de turbulence k − w − SST .
Domaine de calcul
Le domaine de calcul utilisé dans le cas du modèle hybride basé sur le concept de ligne
active est le tiers d’un domaine cylindrique de longueur 20D et de diamètre 10D (D étant
le diamètre de l’éolienne), figure 5.1, le rotor étudié est tripale.
Chaque pale (ici 3) est modélisée par un volume cylindrique (Cylindre Actif, CA) de
diamètre 0.4c, et de longueur égale à la longueur de la pale ; le maillage utilisé pour ce
110
volume compte 20695 mailles. Chaque cylindre actif est entouré par un volume cylindrique,
appelé première partie du domaine, de diamètre égal à la corde et d’axe confondu avec
celui du CA ; ce volume est maillé avec 39709 éléments. La deuxième partie du domaine a
une forme cylindrique de longueur égale à 2 c dans la direction parallèle à l’axe de rotation
de l’éolienne, et de diamètre égal à 1.2D ; le maillage de cette partie est plus grossier que
la première partie du domaine avec 42055 mailles. La troisième partie du domaine est de
la même forme que la deuxième mais avec une longueur de 10c, un diamètre de 1.6D et un
nombre de mailles égal à 60663. La quatrième partie a une forme cylindrique de diamètre
2D, de longueur 1.8D et un nombre de mailles égal à 22078. Ces parties sont réalisées
identiquement pour le deuxième modèle hybride qui remplace la deuxième machine. A
la fin, une grande partie de longueur 20D et de diamètre 10D entoure toutes les parties
précédentes et est maillée avec 724766 mailles. Donc le nombre total des mailles, qui prend
en compte les deux machines modélisées, est égal à 973840 mailles.
20 D
V1 Bloc 1
Bloc 2
Sortie
5D
Entrée
Bloc 3
Bloc 4 Modèle
Bloc 5 hybride
Premièr groupe Deuxième groupe
Fig. 5.1 – Domaine de calcul contenant deux modèles hybrides basés sur le concept de la
ligne active.
Le modèle hybride qui remplace la première machine est situé à une distance 5D de
l’entrée du domaine. Le deuxième modèle hybride remplaçant la deuxième machine placée
en aval de la première est situé à différentes distances du premier modèle permettant ainsi
d’observer l’influence de la distance entre les deux machines.
Le domaine de calcul utilisé dans le cas du modèle hybride basé sur le concept de disque
actif est un domaine cylindrique d’axe confondu avec l’axe de rotation de la machine de
longueur 20D et de diamètre 10D, figure 5.2. La première éolienne est modélisée par un
disque d’épaisseur non nul, ici 0.09522 m (0.09522 = c × sinβ, figure 3.29, page 63 et de
diamètre égal au diamètre de l’éolienne. Ce disque est maillé avec 360000 éléments et est
situé à une distance égale à 5D de l’entrée principale du domaine. La deuxième éolienne
est modélisée par un disque de même géométrie et ayant le même maillage. Différentes
111
positions par rapport à la première éolienne seront prises pour étudier l’influence de la
distance entre les deux machines. Chaque disque est ensuite entouré par un cylindre,
désigné par bloc 1, d’épaisseur c et de diamètre 1.5D, maillé avec 45700 éléments. Le
deuxième cylindre, bloc 2, qui englobe le bloc 1, a une longueur de 0.5D, un diamètre de
2D et un maillage utilisant 11526 mailles. Le bloc 3, englobant le bloc 2, a un diamètre
égal à 3D, une longueur égale à 1D et 69250 mailles. Le bloc 4, entourant le bloc 3, possède
un diamètre égal à 5D, une longueur de 2D et 217594 mailles. Enfin, le bloc 5, entourant
le bloc 4, est un cylindre de diamètre 7D, de longueur 7.5D et utilisant 232768 mailles.
Ces parties sont réalisées identiquement pour le deuxième modèle hybride qui remplace
la deuxième machine. A la fin, une grande partie de longueur de 20D et diamètre 10D
entoure toutes les parties précédentes et est maillée avec 485694 mailles. Donc le nombre
total des mailles, en prend en compte les deux machines remplacées, est égal à 1700000
mailles environ.
Entrée
Bloc 1
Bloc 2 Modèle
Bloc 3 hybride
V1 Bloc 4
10 D
Bloc 5
Sortie
Bloc 6
Premièr groupe Deuxième groupe
20 D
Fig. 5.2 – Domaine de calcul contenant deux modèles hybrides basés sur le concept du
disque actif.
5.3 Effet de l’éloignement

La distance entre deux machines est un élément très important qu’il faut prendre en
compte lors de la construction un parc éolien. Généralement, cette distance est prise entre
112
5 et 9 D dans la direction dominante et entre 3 et 5 D dans la direction perpendiculaire
à la direction dominante.
Dans ce travail, nous présentons l’effet de l’éloignement sur plusieurs distances entre
les deux machines (5 D, 6 D, 7 D, 8 D, 10 D et 15 D), dans la direction dominante, en
utilisant les deux types de modèles hybrides étudiés dans cette thèse. Les éoliennes testées
sont de type NREL phase II avec un profil des pales de type s809.
1- Étude de l’éloignement par le modèle hybride de Cylindre Actif (CA) : dans cette
partie du travail, nous avons remplacé chaque pale du rotor par un Cylindre Actif de
diamètre égale à 0.4 de la corde de la pale et de même longueur que celle-ci. Cette étude
a été faite pour les distances entre éoliennes indiquées précédemment. La vitesse à l’infini
amont est de 10 m/s et l’intensité de la turbulence est égale à 10 %.
Les résultats obtenus sont présentés sur la figure (Fig.5.3) :
113
Distance entre les deux machines P1 P2 ((P1 − P2 )/P1 ) × 100
5D 7729.48 7304.46 5.4986
6D 7730.34 7378.62 4.5501
7D 7730.69 7452.21 3.6023
8D 7730.61 7499.77 2.9859
10D 7731.02 7599.84 1.6968
15D 7730.38 7702.43 0.3615
Fig. 5.3 – Influence de la distance entre deux éoliennes sur la puissance d’une éolienne
située dans le sillage de l’autre.
P1 : Puissance de la première machine.[W ]

P2 : Puissance de la deuxième machine.[W ]
A partir de la figure 5.3, on observe que la valeur de perte de puissance [((P1 −
P2 )/P1 ) × 100] varie fortement lorsque l’espacement des machines est compris entre 5D
et 10D. Au-delà d’un espacement de 10D, la variation est moins forte. On peut donc
dire que l’influence d’une éolienne sur une autre située dans son sillage est fortement
diminuée après d’un espacement est supérieur ou égal à 10D. La figure 5.4 montre bien
cette évolution.
Fig. 5.4 – Influence de la distance entre éoliennes sur la puissance, modèle hybride CA.
114
5D 7500.43 6675.86 10.99
6D 7518.11 6769.02 9.96
7D 7518.75 6845.85 8.94
8D 7519.67 6915.1 8.04
10D 7519.24 7017.44 6.673
15D 7518.76 7246.46 3.62
Fig. 5.5 – Influence de la distance sur la puissance d’une éolienne située dans le sillage
de l’autre en utilisant un modèle hybride DA.
2- Étude de l’effet de l’éloignement par le modèle hybride de Disque Actif (DA) :

dans ce cas, chaque éolienne est remplacée par un disque d’épaisseur liée à la corde de la
pale et de diamètre égal au diamètre du rotor. Nous avons étudié l’interaction entre deux
éoliennes dans les mêmes conditions que l’étude avec le modèle hybride de CA présentée
précédemment. Les résultats sont sur la figure 5.5 :
L’évolution de la perte de la puissance est bien représentée par la figure (Fig.5.6).
de l’autre, modèle DA.
115
Pour comparer les résultats obtenus par le modèle hybride de DA et les résultats
obtenus par le modèle hybride de CA, nous avons présenté touts les résultats sur la même
figure 5.7.
Fig. 5.7 – Comparaison entre les deux modèles hybrides DA et CA.
A partir de la figure 5.7, on trouve qu’il y a des différences importantes entre les deux
modèles hybrides. Pour une distance entre les deux machines égale à 5D, la perte de
puissance est égale à 11 % en utilisant le modèle hybride de DA, par contre, à la même
espacement entre les deux machines, la perte de puissance est diminuée jusqu’à 5.5 % en
utilisant le modèle hybride CA.
Cette différence entre les deux modèles hybrides est due au fait que le sillage derrière
un modèle hybride de DA est beaucoup plus long que le sillage derrière un modèle hybride
de CA.
Pour bien démontrer la différence entre les sillages selon le modèle hybride utilisé, nous
avons présenté une comparaison entre les deux sillages obtenus (selon le modèle hybride
DA et le modèle hybride CA) dans le cas d’interaction entre deux éoliennes. Le cas choisi
d’interaction est celui où l’espacement entre les deux machines est égal à 10D, la vitesse
à l’infini amont est égale à 10 m/s, les deux éoliennes remplacées sont du type NREL
phase II. Après convergence, nous obtenons les résultats présentés dans les deux figures
suivantes 5.8 et 5.9. On observe bien à partir de ces figures que, le sillage derrière le
modèle hybride de CA est plus court que le sillage derrière le modèle hybrides DA, et que
116
le déficit de vitesse en aval des modèles hybride DA est plus fort que celui en aval des
modèles hybrides CA.
Fig. 5.8 – Champ de vitesse axiale dans un cas d’interaction, modèle hybride CA.
L’explication de cette différence est que, dans le cas d’un modèle hybride de DA,
les forces calculées à partir des pales, sont distribuées sur un disque de diamètre D et
d’épaisseur E. Donc, dans ce cas, il y a une surface égale à [π(R2 − r2 )] qui creé donc
un déficit de vitesse plus long que le CA, R est le rayon d’extrémité de la pale, r est le
rayon de moyeu. Par contre, dans le cas d’un modèle hybride de CA, les forces obtenues à
partir des pales ; sont distribuées dans des cylindres remplaçant les pales, dans notre cas
nous avons 3 cylindres. La longueur de chaque cylindre est égal à la longueur de la pale
remplacée par celui-ci, et son diamètre égal est à 0.4 c (c est la corde de la pale). Donc,
la surface qui fait obstacle au vent est égale à 3(0.4c(R − r)).
117
Fig. 5.9 – Champ de vitesse axiale dans un cas d’interaction, modèle hybride DA.
5.4 Effet de la vitesse à l’infini amont

Dans cette partie, nous avons étudié l’influence de la vitesse à l’infini amont. Cette
étude a été faite pour deux machines situées l’une derrière l’autre avec un espacement de
10D.
5.4.1 Effet de la vitesse en amont avec un modèle hybride CA

Pour montrer l’effet de la vitesse du vent en amont sur l’interaction entre deux
éoliennes, nous avons choisi trois valeurs pour cette vitesse, 8 m/s, 10 m/s et 12 m/s.
Nous avons réalisé ce travail sur les deux modèles présentés précédemment, le modèle de
DA et le modèle de CA.
La figure 5.10 représente les résultats obtenus en utilisant le modèle hybride CA pour
une vitesse à l’infini amont égale à 8 m/s. On observe que la valeur perte de puissance
en fonction de la distance entre les deux éoliennes est plus grande que dans le cas où la
vitesse à l’infini amont est égale à 10 m/s.
Par contre, la figure 5.11, représente l’effet d’une éolienne sur l’autre située dans son
sillage pour plusieurs distances entre les deux avec une vitesse à l’infini amont égale à 12
118
Distance entre les deux machines P1 P2 ((P1 − P2 )/P1 )%
5D 3460.99 3130.51 9.5487
6D 3460.31 3192.84 7.7294
7D 3460.65 3254.69 5.9514
8D 3460.54 3296.94 4.7276
10D 3461.36 3378.14 2.4043
15D 3461.04 3446.66 0.4154
de l’autre, V1 = 8 m/s, modèle hybride CA.
m/s. On observe que la perte de puissance dans ce cas, est plus faible que dans les cas à
8 et 10 m/s.

5D 11198.7 10944.9 2.2665
6D 11211.8 10981.1 1.9432
7D 11212.3 11036.7 1.5608
8D 11212.9 11062.7 1.3402
10D 11212.6 11124.1 0.7896
15D 11212.3 11193.7 0.1663
d’une autre, V1 = 12 m/s, modèle hybride CA.
119
Pour comparer les résultats des trois vitesses étudiées, 8 m/s, 10 m/s et 12 m/s, nous
présentons ces résultats dans la figure (Fig.5.12).
Fig. 5.12 – Influence de la vitesse à l’infini amont sur l’interaction entre deux éoliennes,
modèle hybride CA.
La figure 5.12, récapitule tous les résultats d’interaction pour le modèle CA. On y
observe la diminution de la perte de puissance avec l’éloignement ainsi que la diminution
de cette même perte lors que V1 augmente.
120
5.4.2 Effet de la vitesse en amont avec un modèle hybride DA
Dans cette partie, nous allons présenter l’interaction entre deux éoliennes en utilisant
le modèle hybride de DA. Nous présenterons ici, l’effet de l’espacement pour deux valeurs
de vitesses : 8 m/s et 12 m/s ( le cas à 10 m/s a déjà été présenté).

5D 3217.91 2568.65 20.176
6D 3233.07 2646.67 18.137
7D 3233.55 32708.8 16.228
8D 3234.41 2764.52 14.52
10D 3233.91 2845.53 12.00
15D 3233.45 3016.8 6.69
Fig. 5.13 – Influence de l’espacement entre deux machines sur la puissance, V1 = 8m/s,
modèle hybride DA.

5D 10880.2 10276.3 5.55
6D 10878.6 10324.9 5.09
7D 10878.7 10373.9 4.64
8D 10878.7 10420.7 4.21
10D 10878.8 10494.8 3.53
15D 10878.5 10669.7 1.92
Fig. 5.14 – Influence de l’espacement entre deux machines sur la puissance, V1 = 12m/s,
modèle hybride DA
Les figures 5.5, Fig.5.13, Fig.5.14, représentent les résultats de l’interaction pour une
vitesse à l’infini amont égale à 8 m/s, 10 m/s et 12 m/s, l’intensité de la turbulence est
égale à 10%. A partir de ces figures, on trouve que la perte de puissance diminue avec
l’augmentation de la vitesse à l’infini amont. La figure 5.15 représente l’évolution de la
perte de puissance en fonction de la vitesse à l’infini amont et la distance entre les deux
machines.
121
Fig. 5.15 – Influence de la vitesse à l’infini amont sur l’interaction entre deux éoliennes,
modèle hybride DA.
5.5 Interaction entre six machines en cascade

Pour finir, nous présentons les champs de vitesse axiale pour six éoliennes représentées
par des modèles hybrides DA, la distance entre deux machines consécutives est égale à
5D, la vitesse à l’infini amont est égale à 10 m/s et l’éolienne est de type NREL phase II.
Sur la figure 5.16 on observe que chaque machine est située dans le sillage d’une autre, et
le déficit de la vitesse axiale derrière chaque machine, augmente différemment au fur et à
mesure qu’on s’éloigne de l’entrée du domaine. Le déficit entre la première et la deuxième
et plus grand que celui-ci entre la deuxième et la troisième et ainsi de suite.
Les puissances résultantes de ces éoliennes sont sur la figure 5.17, on observe que la
perte de puissance entre deux machines consécutives diminue au fur et à mesure que l’on
s’éloigne de la première machine. La différence de puissance entre la première machine et
la deuxième est de 10 %, la différence entre la troisième machine et la quatrième est de 5
%, la différence entre la quatrième et la cinquième est de 1.5 % et la différence entre la
cinquième machine et la sixième est inférieure à 1 %.
122
Fig. 5.16 – Champ de vitesse axiale pour 6 rotors consécutifs.
Fig. 5.17 – Evolution de la puissance des éoliennes situées sur le même axe de rotation,
Ri = P i
P.
5.6 Conclusion
Dans cette partie, nous avons étudié les effets induits par la présence de plusieurs
éoliennes coaxiales en cascade. Le paramètre le plus important semble être l’espacement
123
entre les éoliennes. En effet, on observe naturellement que la puissance restituée par une
éolienne située dans le sillage d’une autre augmente avec la distance inter-éolienne. D’autre
part, la vitesse à l’infini amont est également un facteur important. On observe que la
perte de puissance d’une éolienne dans le sillage d’une autre est moins sensible lorsque la
vitesse augmente. Il est à noter que nous avons également effectué des tests sur le niveau
de turbulence en entrée du domaine de calcul, cette influence est très négligeable ainsi que
l’ont montré Sicot et al. [100] ainsi que Chohran [101]. Ces résultats prouvent la faisa-
bilité de l’étude proposée, nous regrettons toutefois l’absence de résultats expérimentaux
accessibles pour nous fournir un point de comparaison. Une des perspectives de ce travail
consisterait donc à étudier expérimentalement l’interaction entre deux éoliennes de taille
réduite.
124
Chapitre 6
Conclusion générale
125
Dans ce travail, nous nous sommes intéressés aux écoulements autour des éoliennes et
plus particulièrement à leur sillage, avec comme objectif le développement d’une méthode
de simulation numérique capable de représenter le sillage éolien sans nécessiter un temps
excessif de calcul. Ainsi, avec une telle méthode, il sera possible d’étudier les perfor-
mances de plusieurs éoliennes dans un parc. Pour cela, notre travail a été orienté vers le
développement de modèles hybrides couplant des modèles de rotor simplifiés, comme le
disque ou le cylindre actif, à des simulations numériques de l’écoulement autour du rotor.
Pour analyser les phénomènes existants dans le sillage, avoir une idée claire de la struc-
ture du sillage proche d’une éolienne, et obtenir des données pouvant servir de référence
à la validation de la simulation, un travail expérimental en soufflerie a été réalisé sur une
maquette d’éolienne Rutland 503 modifiée. Des explorations, synchronisées avec une pale
de référence, du sillage ont été faites par la technique PIV dans différents plans azimu-
taux. Les mesures PIV donnent les composantes de vitesse axiale et radiale. Pour explorer
le champ de vitesse tangentielle, la technique de l’anémométrie à fil chaud a été utilisée
pour obtenir les champs de vitesses (tangentielle et axiale) en amont et en aval du ro-
tor. La synchronisation des mesures, par rapport à la position azimutale de la pale de
référence a permis de faire une reconstruction tridimensionnelle du champ de vitesse dans
le sillage. Ces explorations ont permis de visualiser les tourbillons marginaux émis des
extrémités des pales. L’intersection entre les plans d’exploration azimutaux et les tubes
tourbillonnaires hélicoı̈daux a permis de visualiser les foyers tourbillonnaires et le champ
de vitesse induite par chaque foyer. Ainsi, la position des tourbillons marginaux a pu être
localisée et leur pas déterminé. Le diamètre du tube de courant augmente en aval du ro-
tor à cause du ralentissement de l’écoulement créé par l’éolienne. Nos résultats montrent
que les tourbillons marginaux issus des extrémités des pales ne sont pas situés sur une
surface cylindrique comme le suppose la théorie tourbillonnaire linéaire ; ils se déplacent
vers l’extérieur en augmentant le diamètre du tube de courant comme le prévoit la théorie
de Froud-Rankine. L’examen des champs successifs dans le temps montre une fluctuation
et une instationnarité de la position des noyaux des tubes tourbillonnaires dans le sillage.
Cette fluctuation est due à la variation temporelle de la puissance du rotor. En effet, lors
d’une augmentation de la puissance absorbée, la vitesse moyenne dans le sillage diminue
conduisant ainsi à une diminution du pas des tubes tourbillonnaires. En conséquence,
l’équation de continuité fait que le rayon du tube de courant augmente et déplace les
tourbillons vers l’extérieur. On peut noter que la fluctuation radiale des noyaux dépend
faiblement de leurs positions axiales. Par contre, leur battement axial augmente vite parce
que ce battement est amplifié en fonction du nombre de pas parcourus. Ainsi, il est pos-
sible à l’aide de coupes cylindriques et des coupes normales à l’axe du sillage, d’observer
126
les tourbillons marginaux et de suivre leurs trajectoires hélicoı̈dales. Pour les travaux
de simulation numérique du sillage éolien, deux modèles hybrides ont été étudiés : un
modèle de disque actif (DA) et un modèle de cylindre actif (CA) couplés avec une simu-
lation numérique basée sur la résolution des équations de Navier-Stokes. La résolution
itérative des équations de Navier-Stokes donne le champ de vitesse et le calcul des forces
appliquées sur le DA ou CA est faite par la méthode de l’élément de pale. La validation
de ces modèles hybrides a été faite de deux manières différentes : calcul de la puissance en
fonction de la vitesse du vent et calcul du champ de vitesse dans le sillage. Les puissances
calculées par les modèles hybrides ont été comparées avec les données expérimentales
concernant les éoliennes, NREL phase II et phase VI. Les résultats montrent qu’il y a de
bonnes corrélations entre les valeurs calculées et les données expérimentales notamment
par faible vent. Avec l’augmentation de la vitesse du vent, les puissances fournies par les
modèles hybrides sont moins prédictives, plus particulièrement dans le cas NREL phase
VI. Cette divergence est probablement due à l’entrée dans une zone de fonctionnement
en décrochage.
En ce qui concerne le sillage, nous avons comparé le champ de vitesse axiale mesuré en
aval l’éolienne Rutland 503 avec le champ de vitesse axiale calculé en aval par les modèles
hybrides. La comparaison entre les deux champs, montre que la différence est assez faible
et permet de valider nos modèles. Après cette validation, une étude d’interaction entre
des éoliennes a été faite à l’aide des deux modèles hybrides présentés précédemment. Les
résultats montrent que le phénomène d’interaction est influencé par deux facteurs impor-
tants : l’espacement entre les machines et la vitesse du vent à l’infini amont. En effet,
lorsqu’une éolienne est située dans le sillage d’une autre machine, l’augmentation de la
distance entre les deux éoliennes, conduit à une augmentation de la puissance fournie par
la machine à l’aval. Ce résultat s’explique facilement par le fait que le sillage de la première
machine a plus de temps pour récupérer l’énergie cinétique du vent. En revanche, avec une
vitesse à l’infini amont plus élevée, le phénomène de mélange entre le sillage d’une machine
et du vent environnant devient plus importante. Ce mélange participe à la récupération
de l’énergie cinétique par le sillage et augmente la vitesse de ce dernier. Dans ce cas, la
production d’énergie de la deuxième machine est améliorée. L’étude de l’interaction de
plusieurs machines en cascade avec un espacement constant entre deux machines succes-
sives, montre que la différence des puissances entre deux éoliennes successives diminue au
fur à mesure que l’on s’éloigne de la première machine.
La suite de ce travail peut consister à affiner les conditions de simulation pour mieux
l’approcher de la réalité physique. En effet, pour développer les deux modèles hybrides
étudiés dans cette thèse, on a considéré une distribution uniforme de la vitesse à l’infini
127
amont et il serait souhaitable de prendre en compte la présence du gradient de vitesse
crée par la couche limite terrestre. De même, la présence du mat de l’éolienne doit être
prise en compte dans la simulation.
Par manque de données expérimentales sur l’interaction entre éoliennes, il serait intéressant
d’étudier en soufflerie le développement du sillage et la production d’énergie en cas d’inter-
action entre deux éoliennes. Avec ces développements, on pourra améliorer la simulation
d’un ensemble d’éoliennes placées dans parc éolien.
128
Chapitre 7
Annexe
130
Annexe A
Profil s809
Le profil s809, a été développé par le NREL et a été optimisé pour améliorer la produc-
tion d’énergie éolienne. Celui-ci a une base de données expérimentales bien documentée
qui inclut des distributions de pression, des endroits de séparation de la couche limite,
des données de traı̂née, des valeurs de traı̂née et de portance (Somers 1997, Butterfield et
al 1992) [65, 73].
Fig. 1 – Coefficients de portance et de traı̂née pour l’éolienne NREL phase II, profil s809.
Les cordonnées du profil s809 sont sur la figure (Fig.2), L’allrle du profil est donnée
par la figure (Fig.3).
Extrados Intrados
x y z x y z
0.00037 0.00275 0 0.00140 -0.00498 0
0.00575 0.01166 0 0.00933 -0.01272 0
0.01626 0.02133 0 0.02321 -0.02162 0
0.03158 0.03136 0 0.04223 -0.03144 0
0.05147 0.04143 0 0.06579 -0.04199 0
0.07568 0.05132 0 0.09325 -0.05301 0
0.1039 0.06082 0 0.12397 -0.06408 0
0.1358 0.06972 0 0.15752 -0.07467 0
0.17103 0.07786 0 0.19362 -0.08447 0
0.2092 0.08505 0 0.23175 -0.09326 0
0.24987 0.09113 0 0.27129 -0.1006 0
0.29259 0.09594 0 0.31188 -0.10589 0
0.33689 0.09933 0 0.35328 -0.10866 0
0.38223 0.10109 0 0.39541 -0.10842 0
0.42809 0.10101 0 0.43832 -0.10484 0
0.47384 0.09843 0 0.48234 -0.09756 0
0.52005 0.09237 0 0.52837 -0.08697 0
0.56801 0.08356 0 0.57663 -0.07442 0
0.61747 0.07379 0 0.62649 -0.06112 0
0.66718 0.06403 0 0.67710 -0.04792 0
0.71606 0.05462 0 0.72752 -0.03558 0
0.76314 0.04578 0 0.77668 -0.02466 0
0.80756 0.03761 0 0.82348 -0.01559 0
0.84854 0.03017 0 0.86677 -0.00859 0
0.88537 0.02335 0 0.90545 -0.00370 0
0.91763 0.01694 0 0.93852 -0.00075 0
0.94523 0.01101 0 0.96509 0.00054 0
0.96799 0.00600 0 0.98446 0.00065 0
0.98528 0.00245 0 0.99612 0.00024 0
0.99623 0.00054 0 1 0 0
1 0 0 0 0 0
Fig. 2 – Cordonnées du profil s809.

Fig. 3 – Profil s809.
Annexes B, C, D, E
Les annexes B, C, D et E présentent les résultats de simulation du profil s809 et des
cercles équivalents (d = 0.05c − 1c), qui ont permis de choisir le diamètre du cylindre
actif optimal. Les forces aérodynamiques appliquées sur le profil s809 sont uniformement
distribuées à l’intérieur du cercle qui remplace le profil.
Calages du profil s809 étudies : 0°, 6°, 12°, 18°.
Annexe B
Champs de vitesse axiales autour du profil s809 et des cercles équivalentes
(d = 0.05c − 1c).
L’angle de calage du profil est β = 0°.
Fig. 4 – Champ de vitesse axiale autour le profil s809, β = 0°.

Fig. 5 – Champ de vitesse axiale au- Fig. 6 – Champ de vitesse axiale au-
tour d’un cercle de diamètre 0.05 c, β = tour d’un cercle de diamètre 0.1 c, β =
0°. 0°.
0°. 0°.
0°. 0°.
0°. 0°.
tour d’un cercle de diamètre 0.8 c, β = tour d’un cercle de diamètre 1 c, β =
0°. 0°.
Annexe C
Présentation des champs de vitesse axiales autour du profil s809 et des cercles
équivalentes.

tour d’un cercle diamètre de 0.05 c, β = tour d’un cercle de diamètre 0.1 c, β =
6°. 6°.
6°. 6°.
6°. 6°.
6°. 6°.
6°. 6°.
Annexe D
équivalentes.

tour d’un cercle de diamètre de 0.05 c, tour d’un cercle de diamètre de 0.1 c,
β = 12°. β = 12°.
β = 12°. β = 12°.
β = 12°. β = 12°.
β = 12°. β = 12°.
tour d’un cercle de diamètre de 0.8 c, tour d’un cercle de diamètre de 1 c, β =
β = 12°. 12°.
Annexe E
équivalentes.

18°. 18°.
18°. 18°.
18°. 18°.
18°. 18°.
18°. 18°.
Bibliographie
[1] D. Gouriéres. Energie éolienne. Université de Dakar, 1982.

[2] R. Ouziaux, P. P. Jean, and J. Perrier. Mécanique des fluides appliquée. Dunod,
1966.
[3] Dodge. D. M. Illustrated history of wind power development. Darrell Dodge and
TelosNet Web Development, Septembre, 2001.
[4] B. Sørensen. History of, and Recent Progress in, Wind-Energy Utilization. Annual
Review of Energy and the Environment, 20(1) :387–424, 1995.
[5] A. Mirecki. Étude comparative de chaı̂nes de conversion d’énergie dédiées à une
éolienne de petite puissance. Thèse de doctorat à l’institut polytechnique de Toulouse,
Avril, 2005.
[6] Joris PEETERS. Simulation of dynamic drive train loads in a wind turbine. Thèse
de doctorat de l’Université de Katholieke, Belgium, 2006.
[7] J. Vestergaard, L. Brandstrup, and III R. D. Goddard. A Brief History of the Wind
Turbine Industries in Denmark and the United States. Energy Policy, pp. 322-327,
November 2004.
[8] E. Golding. The Generation of Electricity by Wind Power. E.et F. N. Spon, London,
1977.
[9] A. Reeves. Wind Energy for Electric Power : A REPP Issue Brief. Renewable
Energy, July, 2003.
[10] V. Lauber. REFIT and RPS : Options for a harmonised Community framework.
Energy Policy, 32(12) :1405–1414, 2004.
[11] N. I. Meyer. Danish wind power development. Energy for Sustainable Development,
1 :18–25, 1995.
[12] M. Eric Buchet. L’énergie eolienne. 5e Journée de l’Environnement et du
Développement Durable - DESS Gestion de la Planète / Association GAÏA - Nice,
28 septembre 2002.
151
[13] A. R. Henderson, C. Morgan, B. Smith, H. C. Sørensen, R. J. Barthelmie, and
B. Boesmans. Offshore wind energy in europe a review of the state of the art. Wind
Energ. 6 :35–52 (DOI : 10.1002/we.82), 2003.
[14] J. Berger. Charging Ahead : The business of Renewable Energy and What it Means
for America. University of California Press, Berkley, CA, Canyon, 1997.
[15] R. Harrison, E. Hau, and H. Snel. Large Wind Turbines : Design and Economics.
Wiley, 2000.
[16] A. Benlemilim. Promotion de l’electricité Produite à Partir des Eoliennes et Prepa-
ration d’un Cadre Reglementaire Approprie. FIER 2002, Tétouan-Maroc, 2002.
[17] G. Ste-Marie. Le Développement de la Filière Eolienne au Québec et Ses Coûts.
Hydro-Québec, Sections 957, 1500, 2000 et 4250 du SCFP, Juin 2005.
[18] B. Chabot and L. Buquet. Le Développement de l’énergie Eolienne en France en
2005. DEWI Magazin Nr. 29, August 2006.
[19] F. Jésus. Les mécanismes de projet du protocole de Kyoto . Ministère de l’économie
des finances et de l’industrie, Paris, 04 juillet 2005.
[20] B. Mandat. Protocole de Kyoto à la Convention-Cadre des Nations Unies sur les
Changements Climatiques . FCCC/INFORMAL/83, Nations Unies 1998.
[21] H. G. Beyer, B. Lange, and H. P. Waldl. Modelling Tools for Wind Farm Upgrading.
Proc. European Union Wind Energy Conf, Sweden, pages 1069–1072, 20-24 may,
1996.
[22] B. Lange, H. P. Waldl, A. G. Guerrero, D. Heinemann, and R. J. Barthelmie. Mo-
delling of Offshore Wind Turbine Wakes with the Wind Farm Program FLaP. Wind
Energy, 6(1) :87–104, 2003.
[23] Davide Medici. Wind turbine wakes - control and vortex shedding. Technical report,
KTH Mechanics, Royal Institute of Technology, Stockholm, Sweden, 2004.
[24] M. Magnusson, K. G. Rados, and S. G. Voutsinas. A Study of the Flow Downstream
of a Wind Turbine Using Measurments and Simulations. Wind Engineering Vol.20
N.6, 1996.
[25] S. G. Voutsinas, K. G. Rados, and A. Zervos. On the analysis of wake effects in
wind parks. Wind Engng, 14 :204±219, 1990.
[26] K. Badreddinne, H. Ali, and A. David. Optimum project for horizontal axis wind
turbines ‘OPHWT’. Renewable Energy, 30(13) :2019–2043, 2005.
[27] C. Leclerc. Simulation numérique de l’écoulement tridimensionnel turbulent dans
un parc éolien. National Library of Canada= Bibliothque nationale du Canada,
2000.
[28] H. Mattsson and E. Djerf. Evaluation of the software program windfarm and com-
parisons with mesure data from Alsvik. The Aeronautical Research Institute of
Sweden, FFA TN 30, 2000.
[29] J. Tangler and G. Bir. Evaluation of RCAS inflow models for wind turbine analysis.
NREL/TP-500-35109, February 2004.
[30] F. Le Chuiton. Actuator disc modelling for helicopter rotors. Aerospace Science
and Technology, 8(4) :285–297, 2004.
[31] D. M. O’Brien and M. J. Smith. Analysis of rotor-fuselage interactions using va-
rious rotor models. 43rd AIAA Aerospace Sciences Meeting and Exhibit, Reno, NV,
January 10-13, 2005.
[32] H. Snel. Review of the present status of rotor aerodynamics. Wind Energy, 1(s
1) :46–69, 1998.
[33] D. J. Sharpe. A general momentum theory applied to an energy-extracting actuator
disc. Wind Energy, 7(3) :177–188, 2004.
[34] R. E. Froude. On the part played in propulsion by differences of fluid pressure.
Transactions of the Institute of Naval Architects, 30 :390–405, 1889.
[35] W. J. M. Rankine. On the Mechanical Principles of the Action of Propellers. Ins-
titution of Naval Architects, Sweden, 1865.
[36] F. W. Lanchester. A Contribution to the Theory of Propulsion and the Screw
Propeller’. Transactions Institution of Naval Architects, pages 57–98, 1915.
[37] T. Burton, D. Sharpe, N. Jenkins, and E. Bossanyi. Wind Energy Handbook. John
Wiley & Sons, 2001.
[38] J. Martinez, L. Bernabini, O. Probst, and C. Rodriguez. An Improved BEM Model
for the Power Curve Prediction of Stall-regulated Wind Turbines. Wind Energy-
Chichester, 8(4) :385, 2005.
[39] A. Ahlström. Aeroelastic simulation of wind turbine dynamics. Thèse de doctorat
de Department of Mechanics, Royal Institute of Technology, SE-100 44 Stockholm,
Sweden, April, 2005.
[40] C. Alinot and C. Masson. Aerodynamics of wind turbines in thermally stratified
turbulent atmospheric boundary layer. Proceedings of the 10th Annual Conference
of the CFDSC, Ottawa, Canada, pages 553–559, 9-11 Juine 2002.
[41] W. Z. Shen, Sørensen J. N., and R. Mikkelsen. Tip Loss Correction for
Actuator/Navier–Stokes Computations. Journal of Solar Energy Engineering,
127 :209, 2005.
[42] J. N. Sørensen and W. Z. Shen. Numerical modeling of wind turbine wakes. Journal
of Fluids Engineering(Transactions of the ASME), 124(2) :393–399, 2002.
[43] R. Mikkelsen. Actuator disc methods applied to wind turbines. Thèse de doctorat
à la Technical University de Denmark, 2003.
[44] S. S. A. Ivanell. Numerical computations of wind turbine wakes. Royal Institute of
Technology, Stockholm, Sweden, April 2005.
[45] F. Massouh, I. Dobrev, and M. Rapin. Numerical Simulation of Wind Turbine
Performance Using a Hybrid Model. 44th AIAA Aerospace Sciences Meeting and
Exhibit, pages 1–10, 2006.
[46] C. Masson, A. Sma, and C. Leclerc. Aerodynamic Analysis of HAWTs Operating
in Unsteady Conditions. Wind Energ, 2001.
[47] C. Masson and A. Smaı̈li. Numerical Study of Turbulent Flow around a Wind
Turbine Nacelle. Wind Energy-Chichester, 9(3) :281, 2006.
[48] L. Sons. Re-examining the Precepts of the Blade Element Momentum Theory for
Coning Rotors. Wind Energ. 9 :457–478, June 2006.
[49] D. Simms, S. Schreck, M. Hand, and L. J. Fingersh. NREL Unsteady Aerodynamics
Experiment in the NASA-Ames Wind Tunnel : A Comparison of Predictions to
Measurements. National Renewable Energy Laboratory, USA, June 2001.
[50] N. N. Sørensen and J. Michelsen. Navier-Stokes predictions of the NREL Phase VI
rotor in the NASA Ames 80-by-120 wind tunnel. ASME Wind Energy Symposium,
pages 94–105, 2002.
[51] E. P. N. Duque and W. Johnson. Numerical Predictions Of Wind Turbine Power
And Aerodynamic Loads For The NREL phase II Combined Experiment Rotor.
AIAA Paper, 38, 2000.
[52] D. A. Simms. Unsteady aerodynamics experiment phases II-VI test configurations
and available data campaigns. National Renewable Energy Laboratory, USA, July,
1999.
[53] M. Magnusson, K. G. Rados, and S. G. Voutsinas. Study of the flow downstream of
a wind turbine using measurements and simulations. Wind Engineering, 20(6) :389–
403, 1996.
[54] S. G. Voutsinas, K. G. Rados, and A. Zervos. Wake effects in wind parks. A
new modelling approach. Proc. European Community Wind Energy Conf, Lübeck-
Travemünde, Germany, page 444±447, 8—12 March, 1993.
[55] V. P. A. R. I. DE and P. I. V OU. “Vélocimétrie et Granulométrie Laser en
Mécanique des Fluides”. Institut von Karman de Dynamique des Fluides, Belgique,
2005.
[56] F. Massouh, I. Dobrev1, and M. Jourieh. Etude par PIV du sillage proche d’une
éolienne. Actes du Colloque FLUVISU11, 7-9 juin, ECL, Ecully, FRANCE, 2005.
[57] R. Theunissen, P. Corieri, and M. L. Riethmuller. Application de la technique PTV
à la modélisation d’écoulements d’aérosols dans les voies aériennes pulmonaires. 9e
Congrès Francophone de Vélocimétrie Laser, Bruxelles, 14-17 septembre, 2004.
[58] Most J.M. D. Honoré Susset, A. Développement et validation d’un algorithme de
super-résolution par corrélation directe pour la PIV : méthodes d’accélération du
traitement. 8e Congrès Francophone de Vélocimétrie Laser, Orsy, France, 17-20
septembre, 2002.
[59] D. Dimitrov and V. Lazarov. Source d’energie renouvlables. université technique de
sofia, 1999.
[60] L. V. King. On the Convection of Heat from Small Cylinders in a Stream of Fluid :
Determination of the Convection Constants of Small Platinum Wires with Appli-
cations to Hot-Wire Anemometry. Philosophical Transactions of the Royal Society
of London. Series A, Containing Papers of a Mathematical or Physical Character,
214 :373–432, 1914.
[61] A. M. Savill. Fundamentals of Hot Wire Anemometry . Journal of Fluid Mechanics
Digital Archive, 237 :692–693, 2006.
[62] M. A. Kegerise and E. F. Spina. A comparative study of constant-voltage and
constant-temperature hot-wire anemometers. Experiments in Fluids, 29(2) :154–
164, 2000.
[63] S. Blidi and H. Miton. Use of the Hot Wire Anemometry for Velocity and Tempera-
ture Measurements in a Turbomachine. Journal de Physique III, 5(10) :1513–1535,
1995.
[64] F. Nabah, A. Ettaouil, and N. Guennoun. Calcul du sillage d’un rotor éolien à axe
horizontal calculation of wake of an horizontal axis wind rotor. Forum International
sur les Énergies Renouvelables (FIER), Tétouan-Maroc, 2002.
[65] M. Hand, D. Simms, L. J. Fingersh, D. Jager, J. Cotrell, M. Robinson, S. Schreck,
and S. Larwood. Unsteady aerodynamics experiment phase VI : wind tunnel test
configurations and available data campaigns. DRAFT to be published as and NREL
Technical Report, 2001.
[66] J. M. Jonkman. Modeling of the UAE Wind Turbine for Refinement of FAST-AD.
Technical report, NREL National Reneweable Energy Laboratory, NREL/TP-500-
34755, 2003.
[67] P. J. Moriarty. ACH, AeroDyn Theory Manual. National Renewable Energy Labo-
ratory, 2005. Technical report, NREL/TP-500-36881.
[68] J. N. Sørensen and C. W. Kock. Model for unsteady rotor aerodynamics. Journal
of Wind Engineering and Industrial Aerodynamics, 58(3) :259–275, 1995.
[69] J. N. Sørensen, W. Z. Shen, and X. Munduate. Analysis of wake states by a full-field
actuator disc model. Wind Energy, 1(2) :73–88, 1998.
[70] C. Lindenburg. Modelling of Rotational Augmentation Based on Engineering Consi-
drations and Measurments. European Wind Energy Conference, London, 22-25 No-
vember, 2004.
[71] D. M. Somers. Design and experimental results for the s809 airfoil. NREL/SR–440-
6918, National Renewable Energy Lab, Golden, CO (USA), 1997.
[72] J. L. van Ingen, L. M. M. Boermans, and J. J. H. Blom. Low-speed airfoil section
research at Delft University of Technology. International Council of the Aeronautical
Sciences, Congress, 12th, Munich, West Germany, Proceedings, pages 401–416, Oct
1980.
[73] W. P. Wolfe and S. S. Ochs. CFD Calculations of s809 Aerodynamic Characteristics.
16. American Society of Mechanical Engineers wind energy symposium, Reno, NV
(United States), 6-9 Jan, 1997.
[74] J. Tangler and J. D. Kocurek. Wind turbine post-stall airfoil performance charac-
teristics guidelines for blade-element momentum methods. 43rd AIAA Aerospace
Sciences Meeting and Exhibit, Nivada, Reno, October 2005.
[75] R. Van Rooij, A. Bruining, and J. G. Schepers. Validation of some rotor stall models
by analysis of the IEA Annex XVIII field data. Proceedings of EWEC European
Wind Energy Conference Madrid, June, 2003.
[76] Dr. P. K. C. et al. Investigating Three-Dimensional and Rotational Effects on Wind
Turbine Blades by Means of a Quasi-3D Navier-Stokes Solver. Journal of Fluids
Engineering, 122 :330, 2000.
[77] P. K. Chaviaropoulos, I. G. Nikolaou, K. A. Aggelis, N. N. Sørensen, J. Johansen,
M. O. L. Hansen, M. Gaunaa, T. Hambraus, et al. Viscous and Aeroelastic Effects
on Wind Turbine Blades. The VISCEL project. Part I : 3D Navier-Stokes Rotor
simulations. Wind Energy, 6(4) :365–385, 2003.
[78] Z. Du and M. S. Selig. The effect of rotation on the boundary layer of a wind turbine
blade. Renewable Energy, 20(2) :167–181, 2000.
[79] H. Dumitrescu and V. Cardos. Rotational effects on the boundary-layer flow in
wind turbines. AIAA Journal, 42(2) :408–411, 2004.
[80] D. Simms, S. Schreck, M. Hand, L. Fingersh, J. Cotrell, K. Pierce, and M. Robinson.
Plans for Testing the NREL Unsteady Aerodynamics Experiment 10-m Diameter
HAWT in the NASA Ames Wind Tunnel. National Renewable Energy Laboratory,
2000.
[81] G. Xu and L. N. Sankar. Computational study of horizontal axis wind turbines.
Journal of Solar Energy Engineering, 122 :35, 2000.
[82] JL Tangler. The Nebulous Art of Using Wind-Tunnel Airfoil Data for Predicting
Rotor Performance. 21st ASME Wind Energy, Reno, Nevada, 5(2/3), 14-17 January
2002.
[83] F. R. Menter, M. Kuntz R., and Langtry. Ten Years of Industrial Experience with
the SST Turbulence Model, Turbulence. Heat and Mass Transfer, 4 :625–632, 2003.
[84] A. G. Kanaris, Mouza K. A., and Paras S. V. Designing Novel Compact Heat
Exchagers For Improved Efficiency Using a CFD Code. 1st International Conference
“From Scientific Computing to Computational Engineering”, 1st IC-SCCE Athens,
8-10 September, 2004.
[85] L.Davidson. An Introduction to Turbulence Models. Department of Thermo and
Fluid Dynamics,Chalmers University of Technology, Götemberg, Sweden, 2001.
[86] RF Ramsay, MJ Hoffman, and GM Gregorek. Effects of grit roughness and pitch
oscillations on the S809 airfoil. Technical report, NREL/TP–442-7817, National
Renewable Energy Lab., Golden, CO (USA), 1995.
[87] A. Crespo, J. Hernandez, and S. Frandsen. Survey of modelling methods for wind
turbine wakes and wind farms. Wind Energy, 2(1) :1–24, 1999.
[88] A. E. Feijoo and J. Cidras. Modeling of wind farms in the load flow analysis. Power
Systems, IEEE Transactions on, 15(1) :110–115, 2000.
[89] N. D. Kelley and H. J. Sutherland. Damage Estimates from Long-term Structural
Analysis of a Wind Turbine in a US Wind Farm Environment. Wind Energy, 1997.
[90] C. Crafoord. An estimate of the interaction of a limited array of windmills. N-77-
13539 ; DM-16, Stockholm University.(Sweden). Meteorologiska Institutionen, 1975.
[91] C. Crafoord. Interaction in limited arrays of windmills. Report DM-26, Department
of Meteorology, University of Stockholm, 1979.
[92] R. J. Templin. An estimation of the interaction of windmills in widespread ar-
rays. National Aeronautical Establisshment, Laboratory Report LTR-LA-171, Ot-
tawa, Décembre, 1974.
[93] D. P. Rodmell and M. L. Johnson. The Development of Marine Based Wind Energy
Generation and Inshore Fisheries in UK Waters : Are They Compatible ? Scar-
borough Centre for Coastal Studies, Department of Biology, University of Hull at
Scarborough, Filey Road, Scarborough, YO11 3AZ, UK, 2004.
[94] W. Bergström, Schlez, K. Rados, B. Lange, P. Vølund, S. Neckelmann, S. Mogen-
sen, and L. Folkerts. ENDOW (Efficient Development of Offshore Wind Farms) :
Modelling Wake and Boundary Layer Interactions. Wind Energy, 7(3) :225–245,
2004.
[95] K. Rados, G. Larsen, R. Barthelmie, W. Schlez, U. Hassan, B. Lange, I. Waldl,
G. Schepers, T. Hegberg, and M. Magnusson. A Comparison of Wake Model Per-
formances in an Offshore Environment. Ewec-Conference, pages 781–784, 2001.
[96] S. Frandsen and M. L. Thøgersen. Integrated fatigue loading for wind turbines
in wind farms by combining ambient turbulence and wakes. Wind Engineering,
23(6) :327–340, 1999.
[97] M. Magnusson and A. S. Smedman. Air flow behind wind turbines. Journal of
Wind Engineering & Industrial Aerodynamics, 80(1-2) :169–189, 1999.
[98] P. B. S. Lissaman. Energy effectiveness of arbitrary arrays of wind turbines. Ame-
rican Institute of Aeronautics and Astronautics, Aerospace Sciences Meeting, 17th,
New Orleans, 15-17 Jan, 1979.
[99] I. Katic, J. Højstrup, and N. O. Jensen. A simple model for cluster efficiency.
Proceedings of EWEC, 86 :407–410, Rome, Italy, 1986.
[100] C. Sicot, P. Devinant, T. Laverne, S. Loyer, and J. Hureau. Experimental study of
the effect of turbulence on horizontal axis wind turbine aerodynamics. Wind Energ.
9 :361–370, 2006.
[101] C. Cochran. The Influence of Atmospheric Turbulence on the Kinetic Energy Avai-
lable During Small Wind Turbine Power Performance Testing . IEA Expert Meeting
on : Power Performance of Small Wind Turbines Not Connected to The Grid, Spain,
April 2002.

TH' ESE Docteur L' Ecole Nationale Sup Erieure D'arts Et M Etiers

Transféré par

Droits d'auteur :

Formats disponibles

TH' ESE Docteur L' Ecole Nationale Sup Erieure D'arts Et M Etiers

Transféré par

Informations du document

Titre original

Copyright

Formats disponibles

Partager ce document

Partager ou intégrer le document

Options de partage

Avez-vous trouvé ce document utile ?

Ce contenu est-il inapproprié ?

Droits d'auteur :

Formats disponibles

TH' ESE Docteur L' Ecole Nationale Sup Erieure D'arts Et M Etiers

Transféré par

Droits d'auteur :

Formats disponibles

N˚ : 2007 ENAM 0032

École doctorale n˚432 : Science des Métiers de l’Ingénieur

l’École Nationale Supérieure d’Arts et Métiers

Développement d’un modèle représentatif d’une

à soutenir le 20 décembre 2007 à L’ENSAM, devant le Jury composé de :

F. BAKIR Professeur, LEMFI, ENSAM, Paris Président

2 Exploration en soufflerie du champ de vitesse autour d’une éolienne 21

4 Développement d’un modèle hybride basé sur le concept de la ligne

5 Interactions entre des éoliennes 108

6 Conclusion générale 125

V1 : Vitesse du vent à l’infini amont. [m/s]

Fig. 1.1 – Principe du premier moulin à vent persan à axe vertical.

En Inde, l’exploitation de cette énergie renouvelable a pris son essor au quatrième

La croissance de l’exploitation de l’énergie éolienne dans le monde est influencée par un

Fig. 1.5 – Concept de ligne active

1.4 Modèles hybrides

Exploration en soufflerie du champ

2.2 Description des moyens d’essais

Chambre de Explorateur. Veine d’essai

Fig. 2.1 – Soufflerie de l’ENSAM-Paris.

La veine de retour possède une section carrée de 3 m × 3 m sur une longueur de 6 m

2.2.2 Éolienne testée

2.3 Mesures PIV (Particle Image Velocimetry)

Pour simplifier, on suppose la présence de 6 particules qui passeront des positions

Fig. 2.4 – Déplacement des particules.

2.3.1 Mise en place des mesures PIV

Fig. 2.7 – Image brute, fenêtre h1.

2.3.3 Champs de vitesse

Fig. 2.8 – Champ de la vitesse instantanée, fenêtre h1.

Fig. 2.9 – Champ de la vitesse moyennée, fenêtre h1.

Fig. 2.10 – Représentation 3D du sillage proche.

2.3.4 Analyse des résultats

2.4.1 Principe de l’anémométrie à fil chaud

Pont de Wheatstone Servoamplificateur

Fig. 2.13 – Principe du fil chaud.

Il s’agit de maintenir à température constante un fil en tungstène traversé par un

Fig. 2.14 – Trois catégories de sondes à fil chaud.

2.4.2 Méthodologie suivie dans l’exploration au fil chaud

2.4.3 Plan d’exploration

Fig. 2.15 – Plan d’exploration du sillage autour du rotor.

Fig. 2.16 – Points de mesures sur des cercles complets.

Fig. 2.17 – Représentation du sillage.

2.4.5 Étude de trajectoire des tourbillons marginaux

Fig. 2.20 – Champ de la vitesse axiale, R = 280 mm.

Fig. 2.21 – Champ de la vitesse axiale, R = 290 mm.

Fig. 2.30 – Vitesse axiale à Z/D = 0.4 du plan de rotation en amont.

Fig. 2.32 – Vitesse axiale à Z/D = 0.125 du plan de rotation en amont.

2.4.7 Champs de vitesse tangentielle

Fig. 2.35 – Vitesse tangentielle à Z/D = 0.048 en aval du plan de rotation.

La figure 2.35, représente la vitesse tangentielle à une distance de Z/D = 0.048 en

Fig. 2.36 – Vitesse tangentielle à Z/D = 0.1 en aval du plan de rotation.

Développement d’un modèle hybride

3.2 Méthode de l’élément de pale

3.2.1 Théorie de l’élément de pale

Fig. 3.2 – Vitesses et forces agissant sur un élément de pale.

dFz = 0.5ρW 2 Cz c dr (3.5)

dFx = 0.5ρW 2 Cx c dr (3.6)

dFn = dFz cosφ + dFx sinφ (3.7)

La force tangentielle est donnée par :