FR 3064311 A1 20180928 FR 1752337 A 20170322 1 L'invention concerne le domaine des dispositifs pour convertir l'énergie des vagues en énergie électrique ou hydraulique. Les ressources énergétiques renouvelables connaissent un fort intérêt depuis quelques années. Propres, gratuites et inépuisables, autant d'atouts majeurs dans un monde rattrapé par la diminution inexorable des ressources fossiles disponibles et prenant conscience de la nécessité de préserver la planète. Parmi ces ressources, l'énergie houlomotrice, source relativement méconnue au milieu de celles largement médiatisées comme l'éolien ou le solaire, contribue à la diversification indispensable de l'exploitation des énergies renouvelables. Les dispositifs, couramment appelés appareils "houlomoteurs", sont particulièrement intéressants, car ils permettent de produire de l'électricité à partir de cette source d'énergie renouvelable (l'énergie potentielle et cinétique des vagues), sans émission de gaz à effet de serre. Ils sont particulièrement bien adaptés pour fournir de l'électricité à des sites insulaires isolés. Par exemple, les demandes de brevet FR 2876751, FR 2973448 et WO 2009/081042 décrivent des appareils pour capter l'énergie produite par le flot marin. Ces dispositifs sont composés d'un support flottant sur lequel est disposé un pendule monté mobile par rapport au support flottant. Le mouvement relatif du pendule par rapport au support flottant est utilisé pour produire de l'énergie électrique au moyen d'une machine de conversion de l'énergie (par exemple une machine électrique). La machine de conversion fonctionne en tant que générateur et en tant que moteur. En effet, pour fournir un couple ou une force qui entraîne le mobile, on fournit une puissance à la machine de conversion afin de le mettre en résonance avec les vagues (mode moteur). En revanche, pour produire un couple ou force qui résiste au mouvement du mobile, on récupère une puissance via la machine de conversion (mode générateur). Le mouvement des moyens mobiles est donc contrôlé par la machine de conversion de l'énergie pour favoriser la récupération d'énergie. Afin d'optimiser l'énergie électrique récupérée par les systèmes houlomoteurs, différents procédés de commande de la machine de conversion ont été envisagés. Certains de ces procédés ne sont pas optimaux, car la prédiction de la houle n'est pas considérée. De plus, ces procédés ne prennent pas en compte les pertes énergétiques lors de la conversion de l'énergie dans le système houlomoteur. Par exemple, la demande de brevet FR 2973448 (WO 2012/131186) décrit un tel procédé. La commande « PI » est une approche bien connue pour la commande des systèmes houlomoteurs. Par exemple, le document Ringwood, J.V.. Control optimisation and 3064311 2 parametric design, dans Numerical Modelling of Wave Energy Converters: State-of-the art techniques for single WEC and converter arrays (M. Folly Ed.), Elsevier, 2016 décrit une telle approche. Des versions adaptatives par rapport à l'état de mer (c.-à-d. dont les paramètres varient par rapport à l'état de mer) ont été présentées dans la littérature, par exemple dans 5 les documents suivants : - Whittaker, T., Collier, D., Folley, M., Osterried, M., Henry, A., & Crowley, M. (2007, September). The development of Oyster-a shallow water surging wave energy converter. In Proceedings of the 7th European Wave and Tidal Energy Conference, Porto, Portugal (pp. 11-14). 10 - Hansen, Rico H and Kramer, Morten M., "Modelling and control of the Wavestar prototype". In: Proceedings of 2011 European Wave and Tidal Energy Conference (2011) - Hals, J., Faines, J. and Moan, T., "A Comparison of Selected Strategies for Adaptive Control of Wave Energy Converters". In: J. Offshore Mech. Arct. Eng. 133.3 (2011) Pour toutes ces méthodes, l'approche utilisée est le gain scheduling, ou pré 15 programmation des gains : un jeu de gains ou paramètres optimaux (un gain pour la commande P, deux gains pour la commande PI) est calculé hors ligne, analytiquement ou numériquement, pour un ensemble d'états de mer, afin de remplir des cartographies (tableaux) de gains en fonction de l'état de mer. Pour les méthodes décrites dans ces documents, la mise à jour des gains PI s'effectue toujours à partir d'estimations moyennées, 20 sur des fenêtres temporelles de plusieurs minutes (par exemple entre 10 et 30 minutes). Ces méthodes ne permettent donc pas d'adapter les gains « de vague à vague », c'est-à-dire à une fréquence qui est celle de la commande en temps réel (par exemple comprise entre 10 et 100 Hz, soit une période d'échantillonnage comprise entre 10 et 100 ms). Ainsi, le temps de réaction élevé de ces méthodes ne permet une commande optimale du système 25 houlomoteur' la puissance récupérée n'est donc pas optimale. De plus, ces méthodes ne permettent d'optimiser la puissance produite par la machine de conversion, car elles ne prennent pas en compte le rendement de la machine de conversion. Pour pallier ces inconvénients, la présente invention propose un procédé de 30 commande d'un système houlomoteur, dans lequel on estime la force exercée par la houle sur le moyen mobile du système houlomoteur, puis on détermine au moins une fréquence dominante de la force exercée par la houle sur le moyen mobile au moyen d'un filtre de Kalman sans parfum (de l'anglais « Unscented Kalman Filter » UKF), et on détermine la commande du système houlomoteur par une loi de commande PI à gains variables, dont les 35 coefficients sont fonction de la fréquence dominante. Le procédé de commande selon l'invention permet d'adapter les coefficients (gains) des actions P, et I, automatiquement et en continu, à l'état de mer courant de manière à maximiser la puissance récupérable avec la 3064311 3 structure de commande PI. Ceci est réalisé à travers une estimation en ligne, par filtrage UKF (Unscented Kalman Filter), de la fréquence dominante (ou des fréquences dominantes) du spectre de la houle. 5 Le procédé selon l'invention L'invention concerne un procédé de commande d'un système houlomoteur qui convertit l'énergie de la houle en énergie électrique ou hydraulique, ledit système houlomoteur comportant au moins un moyen mobile en liaison avec une machine de conversion de l'énergie, et ledit moyen mobile effectuant un mouvement oscillatoire par 10 rapport à ladite machine de conversion. Pour ce procédé, on réalise les étapes suivantes : a) on mesure la position et/ou la vitesse et/ou l'accélération dudit moyen mobile ; b) on estime la force exercée par la houle sur ledit moyen mobile au moyen de ladite mesure de ladite position et/ou de ladite vitesse dudit moyen mobile ; c) on détermine au moins une fréquence dominante de la force exercée par ladite houle 15 sur ledit moyen mobile au moyen d'un filtre de Kalman non parfumé ; d) on détermine une valeur de commande de ladite force exercée par ladite machine de conversion sur ledit moyen mobile pour maximiser la puissance générée par ladite machine de conversion, au moyen d'une loi de commande proportionnelle intégrale PI à gains variables, dont les coefficients sont déterminés au moyen de ladite 20 fréquence dominante de ladite force exercée par ladite houle sur ledit moyen mobile ; et e) on commande ladite machine de conversion au moyen de ladite valeur de commande. 25 Selon un mode de réalisation, on détermine le coefficient de proportionnalité Kp de ladite loi de commande PI au moyen d'une courbe de résistance de charge optimale à partir de ladite fréquence dominante de ladite force exercée par ladite houle sur ledit moyen mobile. Avantageusement, on détermine le coefficient intégral Ki de ladite loi de commande PI 30 au moyen d'une courbe de réactance de charge optimale à partir de ladite fréquence dominante de ladite force exercée par ladite houle sur ledit moyen mobile. De manière avantageuse, on détermine ladite résistance Rc de la charge optimale et ladite réactance Xc de charge optimale par résolution du problème d'optimisation suivant : MinRc,Xc (Xc+Xi)2 +(Rc+Ri)- 1 71n-nP (Xc arctan (x4))) , (71 p - ir Rc \Rc/ avec Rc 3064311 4 Ri (w)= Bpa (w) Xi(co) = + + Mpa (w) - Kt:2a ) co la fréquence d'excitation, Bpa(w) la résistance à la radiation dudit moyen mobile, M la masse dudit moyen mobile, Mpa(w) la masse ajoutée et M, masse ajoutée à l'infini dudit moyen mobile, Kpa la raideur hydrostatique dudit moyen mobile, Tip le rendement moteur de ladite machine de conversion, et Tipn le rendement générateur de ladite machine de 5 conversion. De prééfrence, lesdites courbes de résistance de charge optimale et/ou de réactance de charge optimale sont déterminées préalablement à la mise en oeuvre des étapes du procédé. Conformément à une mise en oeuvre, on maximise ladite puissance générée par ladite 10 machine de conversion en prenant en compte le rendement de ladite machine de conversion. Selon une caractéristique, ladite loi de commande PI s'écrit par une équation du type : = Kpv(t) + Kip(t) avec fu(t) la commande de la force exercée par ladite machine de conversion sur ledit moyen mobile, v(t) la vitesse dudit moyen mobile, p(t) la position dudit 15 moyen mobile, Kp le coefficient de proportionnalité de ladite loi de commande PI, et Ki le coefficient intégral de ladite loi de commande PI. Conformément à un mode de réalisation, ladite loi de commande est une loi de commande proportionnelle intégrale dérivée PID à gains variables. Avantageusement, ladite loi de commande proportionnelle intégrale dérivée PID s'écrit 20 sous la forme : fu(t) = Kpv(t) + Kip(t) + Kda(t) avec MO la commande de la force exercée par la machine de conversion sur ledit moyen mobile, v(t) la vitesse dudit moyen mobile, p(t) la position dudit moyen mobile, a(t) l'accélération dudit moyen mobile, Kp le coefficient de proportionnalité de ladite loi de commande PID, Ki le coefficient intégral de ladite loi de 25 commande PID, et Kd le coefficient dérivé de ladite loi de commande PID. Selon une option de réalisation, on détermine la fréquence dominante de ladite force exercée par ladite houle sur ledit moyen mobile au moyen d'une modélisation de ladite force exercée par la houle sur ledit moyen mobile comme un signal sinusoïdal, ou comme la somme de deux signaux sinusoïdaux. 30 Selon une mise en oeuvre de l'invention, on estime la position et/ou la vitesse dudit moyen mobile en utilisant un modèle dynamique qui représente l'évolution de la position et de la vitesse du moyen mobile. 3064311 5 De manière avantageuse, ledit modèle dynamique comporte un modèle de la force de radiation. De préférence, ladite machine de conversion est une machine électrique ou une machine hydraulique. 5 Présentation succincte des figures D'autres caractéristiques et avantages du procédé selon l'invention, apparaîtront à la lecture de la description ci-après d'exemples non limitatifs de réalisations, en se référant aux figures annexées et décrites ci-après. 10 La figure 1 illustre les étapes du procédé selon un mode de réalisation de l'invention. La figure 2 illustre la détermination de la valeur de commande selon un mode de réalisation de l'invention. La figure 3 illustre une modélisation du système houlomoteur et de sa commande. La figure 4 illustre un système houlomoteur selon un exemple de réalisation. 15 La figure 5a illustre une courbe de résistance de charge optimale pour un premier exemple d'application. La figure 5b illustre une courbe de réactance de charge optimale pour le premier exemple d'application. La figure 6 illustre un spectre de l'une des vagues irrégulières pour le premier exemple 20 d'application. La figure 7 illustre une estimation de la fréquence dominante de la force de la vague au moyen du procédé de commande selon l'invention pour le premier exemple d'application. La figure 8 illustre une estimation de l'amplitude de la force de la vague au moyen du procédé de commande selon l'invention pour le premier exemple d'application. 25 La figure 9 illustre une courbe comparative entre l'énergie récupérée par le procédé de commande selon l'invention et par un procédé de commande selon l'art antérieur pour le premier exemple d'application. La figure 10 illustre un spectre d'une vague irrégulière pour un deuxième exemple d'application. 30 La figure 11 illustre une courbe comparative entre l'énergie récupérée par le procédé de commande selon l'invention et par un procédé de commande selon l'art antérieur pour le deuxième exemple d'application. 3064311 6 Description détaillée de l'invention L'invention concerne un procédé de commande d'un système houlomoteur qui comprend au moins un moyen mobile (par exemple un flotteur) qui coopère avec au moins une machine de conversion de l'énergie (également appelée PTO de l'anglais « Power Take- 5 Off »). Le moyen mobile effectue un mouvement oscillatoire par rapport à la machine de conversion, sous l'action des vagues (ou houle) et de la machine de conversion. La machine de conversion convertit l'énergie mécanique du mouvement du moyen mobile en énergie électrique. Dans ce but, la machine de conversion peut être une simple machine électrique ou un dispositif plus complexe qui inclut d'autres machines, comme une machine 10 hydraulique. La machine de conversion peut être considérée comme l'actionneur par lequel le système de commande pilote le fonctionnement du système houlomoteur. Notations Au cours de la description, les notations suivantes sont utilisées : 15 fu : force exercée par la machine de conversion sur le moyen mobile, f' : force exercée par la houle sur le moyen mobile, p : position du moyen mobile par rapport à son point d'équilibre. : vitesse du moyen mobile. - a : accélération du moyen mobile. 20 : fréquence du mouvement du moyen mobile. M : masse du moyen mobile. Zpa : impédance de radiation, il s'agit d'un paramètre, fonction de la fréquence et déterminé expérimentalement ou issu du calcul des coefficients hydrodynamiques du moyen mobile, qui permet de rendre compte du phénomène de radiation, selon 25 lequel le mouvement du moyen mobile dans l'eau crée une vague radiée qui l'amortit (c'est-à-dire qui amortit le mouvement) Kpa : coefficient de raideur hydrostatique. Bpa : résistance de radiation, soit la partie réelle de l'impédance de radiation Mpa : masse ajoutée, il s'agit d'un paramètre, fonction de la fréquence et déterminé 30 expérimentalement, ou issu du calcul des coefficients hydrodynamiques du moyen mobile, qui permet de rendre compte du phénomène qui augmente la masse équivalente du moyen mobile, à cause des particules d'eau qui sont entrainées par son mouvement : masse ajoutée à fréquence infiniment élevée, 35 Pa : puissance moyenne générée par le système houlomoteur. t : temps continu s : variable de Laplace 3064311 7 k : temps discret : rendement de la conversion de l'énergie, avec ^ tip : rendement moteur de la machine de conversion, il s'agit d'une donnée constructeur ou d'une donnée déterminée expérimentalement. 5 ^ ï : rendement générateur de la machine de conversion, il s'agit d'une donnée constructeur ou d'une donnée déterminée expérimentalement. - Zi : impédance intrinsèque du moyen mobile du système houlomoteur, paramètre connu, fonction de la fréquence, issu de la modélisation du moyen mobile sur la base de la théorie linéaire des vagues, déterminé à partir des coefficients 10 hydrodynamiques du moyen mobile et éventuellement de mesures expérimentales. - Ri : résistance intrinsèque du moyen mobile du système houlomoteur, soit la partie réelle de l'impédance intrinsèque - Xi : réactance intrinsèque du moyen mobile du système houlomoteur, soit la partie imaginaire de l'impédance intrinsèque 15 - Kp : coefficient de proportionnalité de la loi de commande PI. - Ki : coefficient intégral de la loi de commande PI. - Kd : coefficient dérivé de la loi de commande PID. - Z, : impédance de la loi de commande (ou de charge) - R, : résistance de la loi de commande (ou de charge), avec 20 ^ R° résistance optimale pour un état de mer donné - X, : réactance de la loi de commande, avec ^ e. résistance optimale pour un état de mer donné - A : amplitude de la houle - rex : estimation de la fréquence dominante de la vague 25 f' : estimation de la force exercée par la vague sur le moyen mobile - TS : période d'échantillonnage : déphasage du signal - A, : matrice du modèle d'état - C : matrice du modèle d'état 30 - v(k) : bruit gaussien avec matrice de covariance Q - 1,1(k) : bruit gaussien avec matrice de covariance R - P, : matrice de covariance de x(k) Pour ces notations, les valeurs estimées sont généralement marquées d'un accent circonflexe. Le temps est noté t (variable continue), ou k (variable discrète). 35 Dans la suite de la description et pour les revendications, les termes vagues, flots marins et houle sont considérés comme équivalents. 3064311 8 L'invention concerne un procédé de commande d'un système houlomoteur. La figure 1 représente les différentes étapes du procédé selon l'invention : 1. Mesure(s) de la position et/ou la vitesse du moyen mobile (p, v) 5 2. Estimation de la force exercée par les vagues (EST) 3. Détermination de la fréquence dominante (UKF) 4. Détermination de la valeur de commande (COEFF) 5. Commande de la machine de conversion (COM) 10 Les étapes 1 à 5 sont réalisées en temps réel, dans une boucle temps réel. Toutefois, selon un mode de réalisation de l'invention, la détermination de la valeur de commande peut comprendre un calcul préalable de courbes de résistance et de réactance de charge optimale. 15 De manière avantageuse, le procédé de commande selon l'invention peut être mis en oeuvre par des moyens informatiques, par exemple un ordinateur. ou un processeur, en particulier un processeur embarqué. Etape 1) Mesure(s) de la position et/ou la vitesse du moyen mobile (p,v) 20 Lors de cette étape, la position et/ou la vitesse du moyen mobile sont mesurées. La position correspond au mouvement (par exemple distance ou angle) par rapport à la position d'équilibre du moyen mobile. Ces mesures peuvent être réalisées au moyen de capteurs, généralement implémentés sur un système houlomoteur pour son contrôle et/ou sa supervision. 25 Selon une mise en oeuvre de l'invention, lors de cette étape, on peut également mesurer ou estimer l'accélération du moyen mobile, car cette mesure ou estimation peut être utilisée dans les étapes suivantes du procédé selon l'invention. Par exemple, l'accélération peut être mesurée au moyen d'un accéléromètre placé sur le moyen mobile. 30 Etape 2) Estimation de la force exercée par la houle (EST) Pour cette étape, on estime en temps réel la force exercée par la houle sur le moyen mobile. L'estimation de la houle est réalisée à partir des mesures disponibles (position et/ou vitesse et/ou accélération) obtenues à l'étape précédente. L'estimation de la force exercée par la houle sur le moyen mobile est réalisée en ligne et en temps réel, afin de pouvoir 35 réaliser une commande en temps réel. Dans le but d'une commande ayant un temps de réponse optimal, on peut choisir une méthode d'estimation rapide. 3064311 9 Pour cette étape du procédé, tout type d'estimation de la force exercée par la houle sur le moyen mobile peut être envisagé. Selon un mode de réalisation de l'invention, l'estimation de la force exercée par la 5 houle sur le moyen mobile, peut être mise en oeuvre au moyen d'un estimateur basé sur un modèle dynamique du système houlomoteur. Dans ce cas, on peut construire un modèle dynamique du système houlomoteur. Le modèle dynamique représente le comportement dynamique, traduisant le mouvement des éléments constituant le système houlomoteur sous l'action de la houle et sous l'action de la commande de force donnée à la machine de 10 conversion. Le modèle dynamique est un modèle qui relie la vitesse du moyen mobile à la force exercée par les vagues sur le moyen mobile à la commande de force donnée à la machine de conversion qui est à son tour traduite en une force exercée par la dite machine de conversion sur le moyen mobile. Conformément à un mode de réalisation de l'invention, le modèle de dynamique peut 15 être obtenu par application du principe fondamental de la dynamique sur le moyen mobile du système houlomoteur. Pour cette application, on prend en compte notamment la force exercée par la houle sur le moyen mobile et la force exercée par la machine de conversion sur le moyen mobile. Selon une mise en oeuvre de l'invention, on peut considérer un système houlomoteur 20 avec une partie flottante (moyen mobile) dont le mouvement oscillant de translation ou de rotation est contraint dans une seule dimension. On peut alors supposer que le mouvement de translation ou de rotation puisse être décrit par un modèle linéaire sous forme d'état qui inclut la dynamique du flotteur avec son interaction avec la houle et la dynamique du système de prise de puissance (PTO), ou machine de conversion, qui constitue l'actionneur 25 du système. Dans la suite de la description, seul un mouvement unidirectionnel est considéré pour le modèle dynamique. Toutefois, le modèle dynamique peut être développé pour un mouvement multidirectionnel. Selon un exemple de réalisation, on peut estimer en temps-réel la force exercée par 30 les vagues sur le moyen mobile en mettant en oeuvre un procédé de détermination de la force d'excitation exercée par la houle incidente sur un moyen mobile d'un système houlomoteur au moyen d'un modèle de la force de radiation, tel que décrit dans la demande de brevet dont le numéro de dépôt est FR 16/53109. Pour rappel, la force de radiation est la force appliquée sur le moyen mobile et générée par le mouvement même du moyen mobile, 35 contrairement à la force d'excitation de la houle qui est générée par la houle uniquement. Etape 3) Détermination de la fréquence dominante 3064311 10 Lors de cette étape, on détermine au moins une fréquence dominante de la force exercée par la houle sur le moyen mobile (déterminée à l'étape précédente). On appelle fréquence dominante, la fréquence correspondante au pic (maximum) de son spectre. La fréquence dominante est déterminée au moyen d'un filtre de Kalman non parfumé (de 5 l'anglais « unscented Kalman filter » UKF). Un filtre UKF est basé sur la théorie de la transformation « non parfumée », qui permet d'obtenir un estimateur pour un système non linéaire sans avoir besoin de le linéariser au préalable pour l'application du filtre. Le filtre UKF utilise une distribution statistique de l'état qui est propagée à travers les équations non linéaires. Un tel filtre présente l'avantage de permettre une stabilité, donc une robustesse de 10 l'estimation. Selon un mode de réalisation de l'invention, la force d'excitation de la vague est modélisée comme un signal sinusoïdal variant dans le temps : lex(t) = A(t) sin(co(Ot + 0(0) où A(t), (0(0 et 0(0 sont respectivement l'amplitude, la fréquence et le déphasage 15 du signal. Il s'agit d'une approximation, car cette force n'est pas en réalité une sinusoïde à paramètres variants. Dans la théorie linéaire des vagues, elle est plutôt modélisée comme une superposition de sinusoïdes à paramètres constants. Alternativement, la force d'excitation de la vague peut être modélisée comme la somme de deux signaux sinusoïdaux variant dans le temps : fex(t) = A1(t) sin(coi(t) + 01m) + A2 (t) sin(w2(t) + (t)) 20 Ces paramètres de la modélisation de la force d'excitation de la vague (variant dans le temps) sont à estimer. Puisqu'ils rentrent de manière non linéaire dans l'équation ci-dessus, il s'agit d'un problème d'estimation non linéaire. La méthode du filtre de Kalman non parfumé (UKF, Unscented Kalman Filter, en 25 anglais) est utilisée pour estimer A(t), (0(0 et o(t). On note que d'autres méthodes d'estimations non linéaires pourraient, en principe, être utilisées pour réaliser cette étape tels que le filtre de Kalman étendu (EKF, Extended Kalman Filter) ou le filtre particulaire (Particle Filter), mais l'UKF donne des résultats particulièrement bons. L'utilisation du filtre EKF en particulier a été évoquée dans la 30 littérature dans ce but. Cependant, le filtre EKF étant basé sur la linéarisation locale d'un modèle non linéaire, il ne donne pas les garanties de stabilité, et donc de robustesse, de l'estimation que donne l'UKF. Pour appliquer l'UKF, on met d'abord l'équation modélisant la force d'excitation sous forme d'état et à temps discret. 3064311 11 Soit n la période d'échantillonnage à laquelle le filtre devra tourner. La force d'excitation de la vague sera alors estimée au temps discret t = kTs, k = 0, 1, 2, ..., que l'on notera simplement k. En définissant {xi(k) = A sin(kTsco + tp) , x2(k) = A cos(kTsco + tp) x3(k) = co 5 et dans l'hypothèse que co(k) change lentement dans le temps (par rapport à la période d'échantillonnage), on peut obtenir le modèle sous forme d'état suivant : x1(k) cos(Tsx3(k - 1)) sir(Tsx3 (k - 1)) 0 x1(k - 1) v1(k - 1) x2 (k) = x3 (k - 1)) cos(Tsx3(k - 1)) 0 x2(k - 1) + v2(k - 1) [ x3(k) - sin(Ts.0 x3(k - 1) v3(k - 1) 01 x1(k) fex(k) = [1 0 0] x2(k) +µ(k) x3(k) Dans ce modèle, des incertitudes ont été rajoutées pour prendre en compte les erreurs de modélisation. Plus particulièrement, v1(k - 1) et v2(k - 1) servent à pallier la nature variant dans le temps de A(t), co(t) et (p(t) qui n'est pas pris en compte dans la 10 définition des états x1(k) et x2(k). v3(k - 1) est une incertitude sur le troisième état (la fréquence à estimer), qui est naturellement corrélée avec v1(k - 1) et v2(k - 1). 1,1(k) est une incertitude que l'on peut assimiler à une erreur de mesure sur fex(k), qui sert à prendre en compte le fait que fex(k) n'est pas exactement une sinusoïde. 15 En notant x(k) = [x1(k) x2(k) x3(k)r, v(k) = [v1(k) v2(k) v3 (MT Y (k) = fex(k) et {Ax = cos(Tsx3(k - 1)) sin(Ts.x3 (k - 1)) 0 C = [1 - sin(Tsx3(k - 1)) cos(Tsx3(k - 1)) 0 0 0 1 o o] les équations du modèle d'état s'écrivent : {y(k)= Axx(k - 1) + v(k), ly(k) = Cx(k) + ,u(k) 20 On adopte les hypothèses suivantes : - l'état initial x(0) est un vecteur aléatoire de moyenne m(0) = E[x(0)] et covariance P(0) = E[(x(0) - m(0))(x(0) - m(0))T] ; 3064311 12 - v(k) et 1,1(k) sont des bruits gaussiens avec matrices de covariance Q et R, respectivement ; et les notations suivantes : - x(kIk - 1) est l'estimation de x(k) à partir des mesures jusqu'au temps 5 k - 1, soit y(k - 1), y(k - 2), ... - x(klk) est l'estimation de x(k) à partir des mesures jusqu'au temps k, soit y(k), y(k - 1), ... - Px(kik - 1) est la matrice de covariance de x(k) à partir des mesures jusqu'au temps k - 1, soit y(k - 1), y(k - 2), ... 10 - Px(kik) est la matrice de covariance de x(k) à partir des mesures jusqu'au temps k, soit y(k), y(k - 1), ... Il y a trois étapes dans la méthode UKF : 1. Calcul des points sigma 15 Dans cette première étape, on calcule un ensemble d'échantillons dans l'espace d'état, appelés points sigma, qui représentent bien la distribution probabiliste de l'état selon ses paramètres de moyenne et covariance. Soient Wô= n+Â.' 20 W,,c = n+ (1 - a2 + p), -FÂ. TAjm wc vv1 2(n-Fa), I 1,2, ..., 2n où À = (a2 - 1)n est un paramètre de mise à l'échelle, a est un paramètre qui détermine la diffusion des points sigma autour de x(k - 11k - 1) et auquel on donne généralement une petite valeur positive, par exemple 10-3, fi est un paramètre utilisé 25 pour incorporer la connaissance a priori de la distribution de x : pour une distribution gaussienne, fl = 2 est optimal. Au temps k - 1, on considère le choix suivant de points sigma (ensemble de points qui encodent de manière exacte les informations de moyenne et variance) : xo(k - 1) = x(k - 11k - 1), xi(k - 1) = x(k - 11k - 1) + À Si(k - 1), i = 1,2, ... , n 30 xi+n(k - 1) = x(k - 11k - 1) -nE À Si(k - 1), i = 1,2, ..., n où Si(k - 1) est la i-ème colonne de la racine carrée matricielle de Px(k 11k - 1), soit Px(kik - 1) = S(k - 1)T S(k - 1) 2. Mise à jour des prédictions 3064311 13 Chaque point sigma est propagé à travers le modèle non linéaire qui représente l'évolution des états : xj(k lk - 1) = Axxi (k - 1), j = 0,1, ..., 2n La moyenne et la covariance de 2(klk - 1), la prédiction de x(klk - 1) 2(klk - 1) = e1,3V1r2j(kik - 1), 5 Px(klk - 1) = Elno Wic - 1) - X(IC lk - 1)) - 1) - X(IC lk - 1))T+Q. Les états prédits .5ei(klic - 1) sont utilisés dans l'équation d'état de sortie, ce qui donne : 3y/clic - 1) = C.5ei(klic - 1) 10 La moyenne et la covariance de 3î(kik - 1) sont calculées comme - 1) lno 2, WirnMICIIC - 1), Py(klk - 1) _ Elnowic (Mkik - 1) - y(klk - 1)) (j771(kik - 1) - y(klk - T 1)) R tandis que la covariance croisée entre 2(klk - 1) et j(kik - 1) est : 2n Pxy(klk - 1) =IWic (21(klk - 1) - x(klk - 1)) (3îj(klk - 1) - y(klk - 1)) T i=0 15 3. Mise à jour à partir des mesures Comme dans le filtre de Kalman, l'estimation finale de l'état est obtenue en corrigeant la prédiction avec une rétroaction sur l'erreur de prédiction sur la sortie (mesurée) : 20 2(k) = 2(klk - 1) + K(3-(k) - 3î(kik - 1)) où le gain K est donné par K = Pxy(klk - 1)Py(klk - 1)-1 La covariance a posteriori de l'estimation est mise à jour avec la formule suivante : Px(klk) = Px(klk - 1) - KPy(klk - 1)KT 25 Etape 4) Détermination de la valeur de commande Lors de cette étape, on détermine la valeur de commande de la force exercée par la machine de conversion sur le moyen mobile pour maximiser la puissance générée par la machine de conversion. Cette détermination est mise en oeuvre au moyen d'une loi de 3064311 14 commande proportionnelle intégrale à gains variables, dont les coefficients (les gains variables) sont déterminés en fonction de la fréquence dominante de la force exercée par la houle sur le moyen mobile. Selon une mise en oeuvre de l'invention, la loi de commande peut être une loi de 5 commande de type proportionnelle intégrale dérivée à gains variables. Conformément à un mode de réalisation, on peut écrire un modèle dynamique force-vitesse de ce type de système houlomoteur sous forme fréquentielle K a\ (i(x)A1 Zpa@co) 3)1)(ico) = fex(ico) Mco) co où : 10 - v(j(o) est la vitesse du moyen mobile ; - fex(JW) et fit@ co) sont la force d'excitation de la vague incidente et la force appliquée par la machine de conversion au moyen mobile, respectivement ; -M est la masse du moyen mobile (flotteur, par exemple) et de toutes les autres parties du système houlomoteur solidaires de ce moyen mobile ; 15 -Zpa(j(0) est l'impédance de radiation ; - Kpa est la raideur hydrostatique. Ce modèle est issu de l'équation intégro-différentielle de Cummins, et ses coefficients Kpa et Zpa(j(0) (et ceux qui sont issus de sa décomposition, ci-après) peuvent être calculés à l'aide de codes d'hydrodynamique basés sur la méthode des éléments frontières (BEM), tels 20 que WAMIT, Diodore ou NEMOH. L'impédance de radiation Zpa(j(0), qui, dans la théorie linéaire des vagues, décrit l'effet du mouvement libre du flotteur dans l'eau, est le résultat d'une approximation de la réponse impulsionnelle de radiation par un filtre à réponse impulsionnelle infinie. Elle peut( être décomposée comme il suit : Zpa (jc0) = Bpa(jco) + jcoMa(jco) + M Co) = Hpa(1(0) +IWMoo où Bpa(j(o) est la résistance de radiation, Mpa(j(0) est la masse ajoutée après que la 25 singularité à l'infini M', a été retirée et Hpa(j(o) = Bpa(j(o) + jcoMpa(j(o). La vitesse du flotteur, en fonction des forces qui lui sont appliquées, peut se récrire : 1 v@co) = Vex @co) fu @co)) Zj (jco ) où l'impédance intrinsèque Zi(j(0) est définie comme Zi (kt)) = Bpa((o) + co + Mo, + Mpa(co) - KP2a) = Ri(w) +./Xi(CO) où 3064311 15 {Ri (w) = Bpa (w) (w) = (M+ Mo, + Mpa(w) - Kt:2a) sont respectivement la résistance et la réactance (partie réelle et partie imaginaire de l'impédance) intrinsèques du système. Lors de cette étape, on peut chercher à optimiser les performances de la loi de 5 commande, en termes de maximisation de la puissance électrique produite en moyenne Pa : T 13,7 1 = fu(t)v(t)dt - T t=o où i est un coefficient qui représente le rendement de la machine de conversion. Si = 1, on considère la machine de conversion parfaite, sans pertes de conversion d'énergie. Bien que cette hypothèse soit peu réaliste, elle est souvent adoptée dans la littérature car elle simplifie grandement les calculs, en l'occurrence ceux des paramètres optimaux d'une 10 loi de commande de type PI. Mais elle correspond à considérer que prélever de la puissance du réseau par la machine de conversion a le même coût que d'en délivrer, ce qui est en général faux et peut conduire à obtenir une production d'énergie électrique très inférieure à ce qui est attendu, voire à un gaspillage d'énergie du réseau (Pa négative). Selon une mise en oeuvre, le rendement i peut être considéré comme une fonction de 15 la puissance instantanée fuv, définie comme suit 17p si fuv 0 tiin si fuv 0 Rc,J(c (Xc. Xi)2 p n (R' + Ri)2 z R, R, 25 on obtient les deux courbes de résistance et réactance de charge optimales des figures 5a et 5b, qui sont utilisés dans le calcul en ligne des paramètres de la loi de commande PI. La figure 5a illustre, en fonction de la fréquence w, la résistance optimale R, issue de l'optimisation précédente. Sur cette courbe est représentée en traits pointillés, la résistance 30 intrinsèque Ri, et en trait continu la résistance optimale R, 3064311 19 La figure 5b illustre, en fonction de la fréquence w, la réactance optimale X, issue de l'optimisation précédente. Sur cette courbe est représentée en traits pointillés, la réactance intrinsèque Xi, et en trait continu la réactance optimale X0. Pour valider le procédé proposé, le système houlomoteur, piloté par le procédé de 5 commande PI adaptatif selon l'invention, a été soumis à une série de vagues irrégulières de test. Le spectre S (Nm2s/rad) de l'une de ces vagues, longue d'un peu plus que 1000 s et avec une fréquence dominante d'environ 5 rad/s, est représenté en figure 6. La bonne performance de l'étape d'estimation en ligne de la fréquence dominante de la force de la vague w (rad/s) en fonction du temps t(s) par le procédé de commande selon 10 l'invention est montrée par la figure 7. La figure 8 illustre l'amplitude A de la vague wa, et l'amplitude estimée Aest (en trait foncé) par le procédé de commande selon l'invention. Cette figure montre que l'algorithme UKF utilisé dans le procédé de commande selon l'invention pour estimer les paramètres du modèle sinusoïdal à paramètres variables de la force d'excitation de la vague est capable de 15 gérer de fortes transitions (p.ex. vers t =1040 s dans la figure) tout en permettant d'obtenir des estimations stables de l'amplitude A(t) et de la fréquence dominante w(t)), grâce à sa robustesse vis-à-vis des non linéarités. D'autres approches d'estimation, et en particulier l'EKF, ne seraient pas capables d'une pareille performance. La Figure 9 illustre l'énergie récupérée par la commande PI adaptative selon l'invention 20 INV, et par une commande PI non adaptative selon l'art antérieur AA. Il s'agit de l'énergie récupérée par la commande PI selon l'invention ou selon l'art antérieur à paramètres fixes sur une vague irrégulière. Ainsi, cette courbe compare les performances en récupération d'énergie sur la même vague de la loi de commande PI adaptative selon l'invention et d'un PI à paramètres fixes optimisés à partir des caractéristiques du spectre de l'état de mer, 25 selon un état de l'art. On observe que la loi de commande PI selon l'invention permet même de récupérer plus d'énergie que dans une situation où des paramètres fixes selon l'art antérieur devraient suffire. Un deuxième test a été effectué avec une vague irrégulière dont la fréquence 30 dominante varie de manière linéaire entre la fréquence de la première vague testée (5 rad/s) et celle d'une vague de fréquence 1.8 rad/s. La figure 10 représente le spectre S (Nm2s/rad) de la vague irrégulière en fonction de la fréquence w (rad). La figure 11 est une courbe similaire à la courbe de la figure 9 pour la vague irrégulière dont le spectre est illustré en figure 10. Comme le montre la figure 11, le procédé selon 35 l'invention INV permet dans ce cas de récupérer beaucoup plus d'énergie qu'une commande PI à paramètres fixes calibrée sur le premier spectre, selon un art antérieur. 3064311 20 Par conséquent, le procédé de commande selon l'invention, permet une commande en ligne qui optimise l'énergie récupérée. REVENDICATIONS1) Procédé de commande d'un système houlomoteur qui convertit l'énergie de la houle en énergie électrique ou hydraulique, ledit système houlomoteur comportant au moins un moyen mobile en liaison avec une machine de conversion de l'énergie, et ledit moyen mobile effectuant un mouvement oscillatoire par rapport à ladite machine de conversion, caractérisé en ce qu'on réalise les étapes suivantes : a) on mesure la position et/ou la vitesse et/ou l'accélération dudit moyen mobile ; b) on estime la force exercée par la houle sur ledit moyen mobile au moyen de ladite mesure de ladite position et/ou de ladite vitesse dudit moyen mobile ; c) on détermine au moins une fréquence dominante de la force exercée par ladite houle sur ledit moyen mobile au moyen d'un filtre de Kalman non parfumé ; d) on détermine une valeur de commande de ladite force exercée par ladite machine de conversion sur ledit moyen mobile pour maximiser la puissance générée par ladite machine de conversion, au moyen d'une loi de commande proportionnelle intégrale PI à gains variables, dont les coefficients sont déterminés au moyen de ladite fréquence dominante de ladite force exercée par ladite houle sur ledit moyen mobile ; et e) on commande ladite machine de conversion au moyen de ladite valeur de commande. 20 2) Procédé selon la revendication 1, dans lequel on détermine le coefficient de proportionnalité Kp de ladite loi de commande PI au moyen d'une courbe de résistance de charge optimale à partir de ladite fréquence dominante de ladite force exercée par ladite houle sur ledit moyen mobile. 3) Procédé selon l'une des revendications précédentes, dans lequel on détermine le coefficient intégral Ki de ladite loi de commande PI au moyen d'une courbe de réactance de charge optimale à partir de ladite fréquence dominante de ladite force exercée par ladite houle sur ledit moyen mobile. 4) Procédé selon les revendications 2 et 3, dans lequel on détermine ladite résistance Rc de la charge optimale et ladite réactance Xc de charge optimale par résolution du problème d'optimisation suivant : MinRc)(c (Xc+Xi)2 +(Rc+Ri)- 1 7171-71P (Xc arctan (x4))) , (11 p - ir Rc \Rc/ avec 25 30 Rc 3064311 22 {Ri(co) = Bpa (w) (w) = + + Mpa (w) - Kt:2a) co la fréquence d'excitation, Bpa(w) la résistance à la radiation dudit moyen mobile, M la masse dudit moyen mobile, Mpa(w) la masse ajoutée et M, masse ajoutée à l'infini dudit moyen mobile, Kpa la raideur hydrostatique dudit moyen mobile, Tip le rendement moteur de ladite machine de conversion, et Tipn le rendement générateur de ladite machine de conversion. 5) Procédé selon l'une des revendications 2 à 4, dans lequel lesdites courbes de résistance de charge optimale et/ou de réactance de charge optimale sont déterminées préalablement à la mise en oeuvre des étapes du procédé. 6) Procédé selon l'une des revendications précédentes, dans lequel on maximise ladite puissance générée par ladite machine de conversion en prenant en compte le rendement de ladite machine de conversion. 15 7) Procédé selon l'une des revendications précédentes, dans lequel ladite loi de commande PI s'écrit par une équation du type : = Kpv(t) + Kip(t) avec MO la commande de la force exercée par ladite machine de conversion sur ledit moyen mobile, v(t) la vitesse dudit moyen mobile, p(t) la position dudit moyen mobile, Kp le coefficient de proportionnalité de ladite loi de commande PI, et Ki le coefficient intégral de ladite loi de commande PI. 8) Procédé selon l'une des revendications 1 à 6, dans lequel ladite loi de commande est une loi de commande proportionnelle intégrale dérivée PID à gains variables. 9) Procédé selon la revendication 8, dans lequel ladite loi de commande proportionnelle intégrale dérivée PID s'écrit sous la forme : fu(t) = Kpv(t) + Kip(t) + Kda(t) avec MO la commande de la force exercée par la machine de conversion sur ledit moyen mobile, v(t) la vitesse dudit moyen mobile, p(t) la position dudit moyen mobile, a(t) l'accélération dudit moyen mobile, Kp le coefficient de proportionnalité de ladite loi de commande PID, Ki le coefficient intégral de ladite loi de commande PID, et Kd le coefficient dérivé de ladite loi de commande PID. 5 10 3064311 23 10) Procédé selon l'une des revendications précédentes, dans lequel on détermine la fréquence dominante de ladite force exercée par ladite houle sur ledit moyen mobile au moyen d'une modélisation de ladite force exercée par la houle sur ledit moyen mobile comme un signal sinusoïdal, ou comme la somme de deux signaux sinusoïdaux. 5 11) Procédé selon l'une des revendications précédentes, dans lequel on estime la position et/ou la vitesse dudit moyen mobile en utilisant un modèle dynamique qui représente l'évolution de la position et de la vitesse du moyen mobile. 10 12) Procédé selon la revendication 11, dans lequel ledit modèle dynamique comporte un modèle de la force de radiation. 13) Procédé selon l'une des revendications précédentes, dans lequel ladite machine de conversion est une machine électrique ou une machine hydraulique. 15