Keywords: architecture [34], architecture [34], 3, 5 Distributed MPC architecture [34], 3, 6 Classification of dynamic optimization solution methods [18]

Le travail présenté dans ce mémoire a été mené à CentraleSupélec/Laboratoire des Signaux et Systèmes (L2S).

Je tiens à exprimer ma profonde gratitude au Professeur Didier Dumur pour m'avoir accueilli dans son équipe au sein du L2S. Sa gentillesse, sa disponibilité sur le plan professionnel et humain, ses compétences scientifiques m'ont permis d'effectuer ce travail dans les meilleures conditions. Je remercie également Madame Sihem Tebbani qui a co-encadré cette thèse, pour tous les conseils qu'elle m'a prodigués, la confiance et le suivi régulier qu'elle m'a accordés qui ont largement contribué à rendre ces années de thèse très agréables.

Je remercie le Professeur Hugues Mounier pour avoir présidé ce jury. Je suis également honoré que les professeurs Mohammed M'Saad et Alain Vande Wouwer aient accepté d'être mes rapporteurs. Malgré leurs emplois du temps chargés, ils ont pris le temps de juger mon travail. Je les en remercie vivement. Je tiens également à adresser mes sincères remerciements à Madame Estelle Courtial qui a accepté de faire partie du jury.

Je remercie également tout le personnel du

Dédicaces

Je dédie ce modeste travail :

A mes chers parents Yasmina et Djamel Eddine. Loin de vous, votre sacrifice et votre amour m'ont toujours donné de la force pour prospérer dans la vie. J'ose espérer que ma mère verra à travers cette thèse de doctorat une sorte de concrétisation de tous les efforts qu'elle a fournis pour m'éduquer tout au long de sa vie.

A tonton Mostefa paix à son âme. J'adresse une pensée toute particulière à celui qui restera pour moi une inspiration pour donner le meilleur de moi même.

A ma deuxième mère tata Farah qui m'a accueilli à bras ouverts. En témoignage de l'attachement et de l'affection que je te porte.

A ma tendre femme Kaouter pour son éternel soutien, sa patience et surtout sa précieuse aide chaque fois que nécessaire.

A mes deux frères Abdelhamid et Mohamed Amine ainsi qu'à leurs femmes et ma nièce Manessa.

A ma belle famille et spécialement Djaoued, Ilyes, Asma, Alya et Mami.

i

Résumé

Les dernières années ont permis des développements très rapides, tant au niveau de l'élaboration que de l'application, d'algorithmes de commande prédictive non linéaire (CPNL), avec une gamme relativement large de réalisations industrielles. Un des obstacles les plus significatifs rencontré lors du développement de cette commande est lié aux incertitudes sur le modèle du système. Dans ce contexte, l'objectif principal de cette thèse est la conception de lois de commande prédictives non linéaires robustes vis-à-vis des incertitudes sur le modèle. Classiquement, cette synthèse peut s'obtenir via la résolution d'un problème d'optimisation min-max. L'idée est alors de minimiser l'erreur de suivi de la trajectoire optimale pour la pire réalisation d'incertitudes possible. Cependant, cette formulation de la commande prédictive robuste induit une complexité qui peut être élevée ainsi qu'une charge de calcul importante, notamment dans le cas de systèmes multivariables, avec un nombre de paramètres incertains élevé. Pour y remédier, la principale approche proposée dans ces travaux consiste à simplifier le problème d'optimisation min-max, via l'analyse de sensibilité du modèle vis-à-vis de ses paramètres afin d'en réduire le temps de calcul.

Dans un premier temps, le critère est linéarisé autour des valeurs nominales des paramètres du modèle. Les variables d'optimisation sont soit les commandes du système soit l'incrément de commande sur l'horizon temporel. Le problème d'optimisation initial est alors transformé soit en un problème convexe, soit en un problème de minimisation unidimensionnel, en fonction des contraintes imposées sur les états et les commandes. Une analyse de la stabilité du système en boucle fermée est également proposée.

En dernier lieu, une structure de commande hiérarchisée combinant la commande prédictive robuste linéarisée et une commande par mode glissant intégral est développée afin d'éliminer toute erreur statique en suivi de trajectoire de référence. L'ensemble des stratégies proposées est appliqué à deux cas d'études de commande de bioréacteurs de culture de microorganismes.

v vi

Executive Summary

The last few years have led to very rapid developments, both in the formulation and the application of Nonlinear Model Predictive Control (NMPC) algorithms, with a relatively wide range of industrial achievements. One of the most significant challenges encountered during the development of this control law is due to uncertainties in the model of the system.

In this context, the thesis addresses the design of NMPC control laws robust towards model uncertainties. Usually, the above design can be achieved through solving a min-max optimization problem. In this case, the idea is to minimize the tracking error for the worst possible uncertainty realization. However, this robust approach tends to become too complex to be solved numerically online, especially in the case of multivariable systems with a large number of uncertain parameters. To address this shortfall, the main proposed approach consists in simplifying the min-max optimization problem through a sensitivity analysis of the model with respect to its parameters, in order to reduce the calculation time.

First, the criterion is linearized around the model parameters nominal values. The optimization variables are either the system control inputs or the control increments over the prediction horizon. The initial optimization problem is then converted either into a convex optimization problem, or a one-dimensional minimization problem, depending on the nature of the constraints on the states and commands. The stability analysis of the closedloop system is also addressed.

Finally, a hierarchical control strategy is developed, that combines a robust model predictive control law with an integral sliding mode controller, in order to cancel any tracking error. The proposed approaches are applied through two case studies to the control of microorganisms culture in bioreactors.

Preliminaries

Notations, basic definitions and properties are introduced and will be further used in the next chapters.

Notations

Notation 1. Let N, R, R ≥0 , Z and Z ≥0 denote natural, real, non-negative real, integer and non-negative integer number sets, respectively.

Notation 2. 0 n×m ∈ R n×m is the zero matrix of dimension n × m and I n ∈ R n×n is the identity matrix of dimension n × n. A symmetric n×n real matrix A is said to be positive semidefinite if the scalar z Az is non-negative for every non-zero column vector z of n real numbers. It is denoted A 0 Definition 2. A symmetric n×n real matrix A is said to be positive definite if the scalar z Az is positive for every non-zero column vector z of n real numbers. It is denoted A 0 Definition 3. A matrix A ∈ R m×n is full column rank if and only if A A is invertible.

Definition 4. Given an affine space A, a set B ⊆ A is said to be convex if ∀x, y ∈ B, ∀t ∈ [0, 1] : (1t)x + ty ∈ B. Definition 5. Let A be a convex set and let f :

A -→ R. f is said to be convex if ∀a, b ∈ A, ∀t ∈ [0, 1] : f (ta + (1 -t)b) ≤ tf (a) + (1 -t)f (b).
is called convex if the objective function f : R n -→ R is convex on R n and the constraints g 1 , . . . , g m : R n -→ R are convex. Definition 7. A function f is said to be positive definite if f (x) > 0 for all x > 0. Definition 8. A function f : R -→ R is said to be C n function with n ∈ Z ≥0 , if the first n derivatives f (.), f (.), . . . , f (n) (.) all exist and are continuous with respect to their argument. Definition 9. A function f (x) is said to be locally lipschitz with respect to its argument x if there exists a positive scalar L f x (so-called Lipschitz constant) such that |f (x 1 )f (x 2 )| ≤ L f x |x 1x 2 | for all x 1 and x 2 in a given region of x. Definition 10. A function α(.): R ≥0 -→ R ≥0 is a K-function (or of class K) if it is continuous, positive definite, strictly increasing and α(0) = 0. Definition 11. A function β(.): R ≥0 -→ R ≥0 is a K ∞ -function if it is a K-function and also β(s) -→ ∞ as s -→ ∞. xxvi Definition 12. A function γ(., .): R ≥0 × Z ≥0 -→ R ≥0 is of class KL if, for each fixed t ≥ 0, γ(., t) is of class K, for each fixed s ≥ 0, is decreasing and γ(s, t) -→ 0 as t -→ ∞.

Definition 13. Let us consider an autonomous system

x k+1 = f (x k , w k ), k ≥ 0, x 0 = x (3)
where x k ∈ X is the state of the system, w k ∈ W is the disturbance vector (s.t. X and W are compact set that contain the origin).

A function V (.): R n -→ R ≥0 is called a Lyapunov function for system (3), if there exist sets X and K ∞ -functions τ 1 , τ 2 , τ 3 s.t.

τ 1 (|x k |) ≤ V (x k ) ≤ τ 2 (|x k |) ∆V (x k , w k ) = V (f (x k , w k )) -V (x k ) ≤ -τ 3 (|x k |) (4)
with x k ∈ X and w k ∈ W Definition 14. Let us consider the following discrete-time nonlinear system given by:

x k+1 = f (x k , w k ), k ≥ 0, x 0 = x (5)
with

|w k | ≤ γ(|x k |) + µ (6)
where γ(.) is a K-function and µ is a modelled bound of uncertainties.

A function V (.): R n -→ R ≥0 is called a robust Lyapunov function if ∃ K ∞functions α 1 , α 2 , α 3 and σ, and a K-function ζ s.t.

α 1 (|x k |) ≤ V (x k ) ≤ α 2 (|x k |) + σ(η) V (f (x k , w k )) -V (x k ) ≤ -α 3 (|x k |) + ζ(η) (7)
with x k ∈ X and w k ∈ W Definition 15. A set Φ ⊂ R n is a robust positively invariant set for the system (5), if f (x k , w k ) ∈ Φ, ∀x k ∈ Φ and ∀w k ∈ W.

Definition 16. Every pair (x e , u e ) satisfying F (x e , u e ) = 0 is called a steadystate of the ODE ẋ(t) = F (x, u). This means that a process is at steady state (x e , u e ) if and only if it remains at the same state when the input u e is applied.

Definition 17. Consider two vector fields f (x) and g(x) in R n . The Lie bracket operation generates a new vector field:

[f, g] = ∂g ∂x f - ∂f ∂x g (8)
xxvii Then, higher order Lie brackets can be defined as follows:

         ad f g [f, g] ad 2 f g [f, [f, g]] . . . ad k f g [f, ad k-1 f g] (9)
xxviii Chapitre 1

Résumé

La commande des systèmes non linéaires soumis à des contraintes physiques sur l'état et l'entrée prend de plus en plus d'importance aux yeux de la communauté de l'automatique. Les méthodes classiques de commande non linéaire telles que le retour d'état linéarisant et les approches fondées sur le formalisme de Lyapunov offrent des solutions élégantes. Malheureusement, le développement de telles lois de commande devient de plus en plus difficile de par la complexité croissante du modèle mathématique nécessaire à leur mise en oeuvre. C'est pour cette raison que la commande prédictive non linéaire (CPNL, nonlinear Model Predictive Control, NMPC en anglais) est une excellente alternative dont la méthodologie de synthèse reste assez simple tout en prenant en considération les contraintes et les incertitudes du modèle. Par ailleurs, la commande des systèmes complexes, fortement non linéaires, incertains tels que les bioprocédés, s'avère être une tâche très délicate. En effet, dans ce cas, les paramètres du modèle sont généralement connus uniquement avec un intervalle de confiance associé (déterminé à partir d'une procédure d'identification par exemple). Par conséquent, toute stratégie de commande et en particulier la commande prédictive doit être robustifiée pour compenser le manque d'information et/ou inexactitudes paramètres. Il existe deux alternatives principales pour la robustification de la commande prédictive : l'approche stochastique (théorie probabiliste) et l'approche déterministe (formulation d'un problème d'optimisation de type min-max fondé sur la théorie des jeux), qui apparaissent comme des approches coûteuses en termes de temps de calcul. Ainsi, le défi qui doit être relevé dans ce travail consiste à synthétiser une loi de commande prédictive robuste vis-à-vis des incertitudes paramétriques avec un temps de calcul raisonnable, pour permettre une implémentation en temps réel. En partant de la formulation min-max de la commande prédictive, de nouvelles structures sont développées pour satisfaire le compromis entre la robustesse et la charge de calcul, 1 Résumé avec une application dédiée aux bioprocédés.

Cette thèse est structurée comme suit. Le Chapitre 2, non repris et résumé ci-après, propose une introduction portant sur le contexte, les motivations, les contributions et les publications issues des résultats obtenus pendant ces travaux. Le Chapitre 3 présente un état de l'art sur les techniques de commande prédictive existantes. Par la suite, la commande prédictive non linéaire permettant la poursuite d'une trajectoire de référence pour un système non linéaire avec une formulation continue/discrète est développée dans le Chapitre 4. De plus, une variante de la CPNL y est proposée pour permettre de tenir compte de l'écart entre la prédiction du modèle et la sortie du système. Dans le Chapitre 5, la commande prédictive robuste vis-à-vis d'incertitudes paramétriques est formulée en un problème d'optimisation du type min-max. Cependant, cette formulation induit une complexité qui peut être élevée ainsi qu'une charge de calcul importante, notamment dans le cas de systèmes multivariables, avec un nombre de paramètres incertains élevé. Pour y remédier, les deux approches proposées dans ces travaux consistent à simplifier le problème d'optimisation min-max, via l'analyse de sensibilité du modèle vis-à-vis de ses paramètres et la linéarisation des sorties prédites autour des valeurs nominales des paramètres du modèle et de la séquence de commande de référence, respectivement. La deuxième approche débouche sur la résolution d'un problème d'optimisation scalaire. Une analyse de la stabilité du système en boucle fermée est également proposée. La prise en compte de contraintes de type inégalité sur les variables d'optimisation transforme le problème initial en un problème convexe à deux niveaux : un niveau haut avec une optimisation scalaire et un niveau bas avec un problème de programmation quadratique. Le Chapitre 6 détaille deux améliorations de la loi de commande proposée. En premier lieu, la linéarisation des sorties prédites autour des valeurs nominales des paramètres du modèle et de la séquence de commandes optimales obtenue à l'itération précédente. Cette modification permet de rendre la solution moins sensible aux bruits de mesure. En second lieu, une structure de commande hiérarchisée combinant la commande prédictive robuste linéarisée et une commande auxiliaire (proportionnel intégrateur ou mode glissant intégral) est développée afin d'éliminer toute erreur statique en suivi de trajectoire de référence. Le Chapitre 7 est dédié à l'application des lois de commande développées à un cas d'étude de commande de bioréacteur de culture de microorganismes. Enfin, le Chapitre 8 reprend et résume le contenu du manuscrit ainsi que les perspectives proposées dans la continuité de ce travail. Chacun de ces chapitres se trouve résumé ci-dessous.

1.1 Chapitre 3 : Commande prédictive -état de l'art et principales stratégies

Le but de ce chapitre est de donner une introduction sur l'une des stratégies les plus puissantes permettant de résoudre les problèmes de commande non linéaires. Nous présentons ci-dessous la formulation du problème de commande prédictive ainsi qu'un état de l'art non exhaustif des principales stratégies prédictives rencontrées dans la littérature.

Contexte

La théorie de la commande optimale non linéaire prend son essor dans les années 1950-1960 résultant du principe du maximum de Pontryagin [START_REF] Pontryagin | The Mathematical Theory of Optimal Processes[END_REF] et de la programmation dynamique développée par Bellman [START_REF] Bellman | The theory of dynamic programming[END_REF]. Le principe de l'horizon fuyant, concept clé de la commande prédictive, a été proposé par Propoi [START_REF] Propoi | [END_REF]. La philosophie de l'horizon fuyant est relativement proche du problème de la commande optimale en temps minimal et de la programmation linéaire [START_REF] Zadeh | On optimal contro and linear programming[END_REF]. Dans les années 70, la commande prédictive, grâce au progrès des méthodes de résolution, est devenue populaire dans l'ingénierie. Richalet et al. [START_REF] Richalet | Model predictive heuristic control: Application to industrial processes[END_REF] et Cutler et al. [START_REF] Cutler | Dynamic matrix control-a computer control algorithm[END_REF] ont été les premiers à proposer l'application de la commande prédictive linéaire dans l'industrie. En général, beaucoup de systèmes sont fortement non linéaires. Dans ce cas, l'utilisation d'un modèle linéaire est souvent insuffisante pour bien décrire la dynamique du système de sorte qu'un modèle non linéaire doit être utilisé. Divers travaux sur la commande prédictive non linéaire ont été proposés dans la littérature [START_REF] Chen | On receding horizon feedback control[END_REF][START_REF] Keerthi | Optimal infinite-horizon feedback laws for a general class of constrained discrete-time systems: Stability and moving-horizon approximations[END_REF][START_REF] Michalska | Robust receding horizon control of constrained nonlinear systems[END_REF][START_REF] Parisini | A receding-horizon regulator for nonlinear systems and a neural approximation[END_REF][START_REF] Alamir | Stability of a truncated infinite constrained receding horizon scheme: the general discrete nonlinear case[END_REF][START_REF] Magni | Stability margins of nonlinear recedinghorizon control via inverse optimality[END_REF][START_REF] Jadbabaie | On the stability of receding horizon control with a general terminal cost[END_REF]. La commande prédictive a impacté significativement le monde des techniques de régulation industrielle. Par conséquent, de nombreuse applications se retrouvent dans différents domaines : robotique [START_REF] Lelić | Generalized pole-placement self-tuning controller part 1, basic algorithm[END_REF], anesthésie clinique [START_REF] Linkers | Advances in Model-Based Predictive Control, chapter Generalized Predictive Control in Clinical Anesthesia[END_REF], industrie du ciment et usines de pâte à papier [START_REF] Sánchez | Adaptive Predictive Control: From the Concepts to Plant Optimization[END_REF], tours de séchage et bras de robot [START_REF] Clarke | Application of generalized predictive control to industrial processes[END_REF], colonnes de distillation [START_REF] Huang | Advanced step nonlinear model predictive control for air separation units[END_REF][START_REF] Richalet | Industrial applications of model based predictive control[END_REF], systèmes biochimiques [START_REF] Ashoori | Optimal control of a nonlinear fed-batch fermentation process using model predictive approach[END_REF], secteur pétrochimique [START_REF] García | Model predictive control: Theory and practice-a survey[END_REF], aérospatial [START_REF] Hartley | Predictive control using an FPGA with application to aircraft control[END_REF][START_REF] Bourgeois | Launcher atmospheric guidance based on nonlinear model predictive control[END_REF][START_REF] Crassidis | Nonlinear predictive control of spacecraft[END_REF], automobile [START_REF] Re | Automotive Model Predictive Control: Models, Methods and Applications[END_REF], etc.

Principes de la commande prédictive

Soit le système discret non linéaire suivant :

x k+1 = f (x k , u k ) (1.1) oú x k ∈ R nx est l'état du système et u k ∈ R nu est l'entrée de commande.
Étant donné le système (1.1), une formulation générale du problème de com-mande prédictive non linéaire est donnée par min Tous les algorithmes de commande prédictive ont les éléments communs suivants.

u k ,...,u k+Np-1 Np i=1 φ(x k+i , u k+i-1 ) (1.2) s.c xk+i+1 = f (x k+i , u k+i ), pour i = 0, N p -1 xk = x k (1.

Modèle de prédiction

Un modèle du système est une représentation mathématique de la réalité. Il vise à prévoir l'évolution temporelle de chaque variable d'état et chaque sortie du système. En général, les modèles sont liés à la fois aux paramètres, les signaux d'entrée et de sortie. L'utilisation explicite du modèle (1.3) est due à la nécessité de calculer les états prédits xk+i sur un horizon de prédiction déterminé de longueur N p . Ces états prédits dépendent des valeurs connues jusqu'à l'instant k (i.e. x k ) et de la séquence de commandes futures, u k+i , i = 0, N p -1, qui sont les variables d'optimisation.

Fonction de coût

Le critère qui doit être optimisé est une fonction quadratique qui peut prendre par exemple la forme suivante :

Np i=1 φ(x k+i , u k+i-1 ) = Np i=1 ||x k+i -x r k+i || 2 α i + Nc i=1 ||∆u k+i-1 || 2 β i (1.4) avec ∆u k+j = 0, j = N c , N p -1 et ∆u k+i-1 = u k+i-1 -u k+i-2
, et où x est la prédiction du modèle, x r la référence sur la fenêtre de prédiction et ∆u l'effort de commande nécessaire pour atteindre l'objectif spécifié. Le critère considère généralement les paramètres d'optimisation suivants :

• N p et N c sont les horizons de prédiction et de commande respectivement, qui ne sont pas nécessairement égaux (N c ≤ N p ).

• x r est l'état de référence utilisé pour spécifier le comportement en boucle fermée et les performances de suivi attendues.

• α i ∈ R nx×nx et β i ∈ R nu×nu sont les matrices de pondération sur l'erreur de suivi pour l'état et sur la commande respectivement (ici l'incrément de commande).

Les paramètres de réglage doivent être ajustés afin d'éviter l'instabilité du système et d'avoir de bonnes performances en boucle fermée. De plus, des contraintes peuvent être prises en compte dans le problème d'optimisation. Par exemple, les actionneurs peuvent avoir un champ d'action limité. En ajoutant ces contraintes au problème d'optimisation, ce dernier devient plus complexe, de sorte que la solution ne peut pas être obtenue de manière explicite comme dans le cas sans contrainte. En outre, dans ce cas, le paradigme lié à l'obtention d'une solution faisable au problème d'optimisation doit être abordé.

Calcul de la loi de commande

L'optimiseur est un élément fondamental de la stratégie CPNL car il fournit les actions de commande en minimisant la fonction de coût (1.4). Une solution analytique peut être obtenue dans le cas d'un critère quadratique si le modèle est linéaire et qu'il n'y a pas de contraintes. Dans le cas contraire, la séquence de commande est déterminée à partir d'une stratégie d'optimisation en temps réel, impliquant la résolution d'un problème de programmation non linéaire. La taille du problème d'optimisation dépend du nombre de variables et de l'horizon de prédiction considéré. La démarche de résolution du problème d'optimisation peut être formalisée par les étapes suivantes :

1. Obtention des mesures 2. Calcul des prédictions via le modèle du système sur un certain horizon de prédiction 3. Calcul de la séquence de commandes optimales par minimisation de la fonction de coût 4. Implémentation de la première commande au système 5. Retour à l'étape 1.

État de l'art

En général, la commande prédictive peut être divisée en deux classes méthodologiques : linéaire et non linéaire. Maciejowski et Camacho donnent un bon aperçu de la théorie linéaire [START_REF] Maciejowski | Predictive Control: With Constraints[END_REF], tandis que Allgöwer et Zheng donnent un aperçu des méthodes non linéaires [START_REF] Allgöwer | Nonlinear Model Predictive Control. Progress in Systems and Control Theory[END_REF]. Les différents algorithmes de la commande prédictive résultent de la façon dont la fonction de coût à minimiser et les contraintes du système sont spécifiées. Selon les différentes stratégies de modélisation et de formulation du problème, de nombreuses variantes de la CPL existent dans la littérature : Dynamic Matrix Control (DMC), Model Algorithmic Control (MAC), Predictive Functional Control (PFC), Extended Prediction Self Adaptive Control (EP-SAC), Extended Horizon Adaptive Control (EHAC), Generalized Predictive Control (GPC).

Cas non linéaire

Du fait des limitations de la commande prédictive linéaire vis-à-vis des procédés ayant une dynamique fortement non linéaire, soumis à des contraintes et/ou régis par un changement fréquent de régimes de fonctionnement, l'application de la commande prédictive non linéaire (CPNL) est à privilégier. L'utilisation d'un modèle non linéaire pour la prédiction transforme le problème quadratique convexe en un problème non convexe, plus difficile à résoudre. En conséquence, la convergence vers l'optimum global va fortement dépendre de l'étape d'initialisation. La CPNL est une méthode fondée sur l'optimisation pour la commande par retour d'état des systèmes non linéaires. Ses principales applications sont la stabilisation et les problèmes de suivi de trajectoires de référence. La CPNL est fondée comme la CPL sur le principe de l'horizon fuyant, où un problème de commande optimale en boucle ouverte sur un horizon fini est résolu à chaque instant d'échantillonnage et la séquence de commande optimale est appliquée jusqu'à ce qu'une nouvelle séquence de commande optimale soit disponible à l'instant suivant d'échantillonnage. Sa philosophie est donc similaire à celle de la CPL. Les méthodes considérées comme les plus représentatives sont les suivantes : CPNL à horizon infini, CPNL à horizon fini avec contrainte terminale d'égalité, CPNL à horizon fini avec contrainte terminale d'inégalité, CPNL à horizon fini avec coût terminal, CPNL à horizon quasi-infini, CPNL à retour d'état/sortie, CPNL économique. Plusieurs approches ont été proposées dans la littérature. Elles diffèrent les unes des autres par les hypothèses faites sur le type d'interactions entre les différentes composantes de l'ensemble du système (algorithmes totalement ou partiellement connectés), le modèle d'échange d'informations entre les sous-systèmes, et la technique de conception du régulateur utilisé localement.

Architectures de commande prédictive

Robustesse

La CPNL est souvent utilisée dans des problèmes de commande. Cependant, les hypothèses de nominalité réduisent son utilisation ou ses performances. Par conséquent, d'autres types de commande prédictive sont développés pour tenir compte des incertitudes liées au modèle et des perturbations agissant sur ce système. Une large gamme de formulations qui incluent la robustesse dans la formulation du problème d'optimisation, existe dans la littérature. On distingue quatre types d'approches principales :

• Formulations de type LPV (Linear Parameter Varying) [START_REF] Angeli | Constrained predictive control of nonlinear plants via polytopic linear system embedding[END_REF][START_REF] Casavola | Predictive control of constrained nonlinear systems via LPV linear embedding[END_REF][START_REF] Wan | Efficient scheduled stabilizing output feedback model predictive control for constrained nonlinear systems[END_REF].

• Min-max NMPC [START_REF] Limon | Input-to-state stable MPC for constrained discrete-time nonlinear systems with bounded additive uncertainties[END_REF][START_REF] Magni | Assessment and Future Directions of Nonlinear Model Predictive Control, chapter Robustness and Robust Design of MPC for Nonlinear Discrete-Time Systems[END_REF][START_REF] Michalska | Robust receding horizon control of constrained nonlinear systems[END_REF][START_REF] Limon | Input to state stability of min-max MPC controllers for nonlinear systems with bounded uncertainties[END_REF][START_REF] Magni | Robust model predictive control for nonlinear discrete time systems[END_REF][START_REF] Magni | Assessment and Future Directions of Nonlinear Model Predictive Control, chapter Robustness and Robust Design of MPC for Nonlinear Discrete-Time Systems[END_REF] • H ∞ -NMPC [START_REF] Gautam | Robust H ∞ receding horizon control for a class of coordinated control problems involving dynamically decoupled subsystems[END_REF][START_REF] Wang | h ∞ control for networked predictive control systems based on the switched lyapunov function method[END_REF][START_REF] Wang | h ∞ controller design for networked predictive control systems based on the average dwell-time approach[END_REF].

• Commande prédictive stochastique [START_REF] Cannon | Nonlinear Model Predictive Control: Towards New Challenging Applications, chapter Successive Linearization NMPC for a Class of Stochastic Nonlinear Systems[END_REF][START_REF] Goodwin | Nonlinear Model Predictive Control: Towards New Challenging Applications, chapter A Vector Quantization Approach to Scenario Generation for Stochastic NMPC[END_REF].

Stabilité

L'une des questions clés de la commande prédictive est certainement la stabilité du système en boucle fermée. En effet, après avoir défini les principes théoriques de base de la commande prédictive, des sujets plus avancés tels que la stabilité robuste en présence de perturbations, l'estimation des performances et l'efficacité des algorithmes numériques ont été abordés dans la littérature. La stabilité de la commande prédictive pour les systèmes non linéaires contraints nécessite le recours à la théorie de la stabilité au sens de Lyapunov, qui peut être exprimée commodément par le biais des fonctions dites fonctions de comparaison, introduites dans la théorie de la commande non linéaire par Sontag [START_REF] Sontag | Smooth stabilization implies coprime factorization[END_REF]. La stabilité du système sans perturbation est appelée stabilité nominale, tandis que la stabilité du système en présence de perturbations est appelée stabilité robuste.

Stabilité nominale

La stabilité nominale d'un système en boucle fermée est prouvée par l'existence d'une fonction de Lyapunov pour le système. Les stratégies qui garantissent la stabilité nominale rencontrées dans la littérature sont : CPNL à horizon infini, CPNL à horizon fini avec contrainte terminale d'égalité, CPNL à horizon fini avec contrainte terminale d'inégalité, CPNL à horizon fini avec coût terminal, CPNL à horizon quasi-infini.

Stabilité robuste

Comme indiqué précédemment, la conception de lois de commande prédictive robuste se fait en prenant en compte les incertitudes d'une manière explicite afin d'optimiser la fonction objectif pour la pire configuration des incertitudes. Il existe de nombreuses approches permettant d'analyser la stabilité robuste, comme la structure entrée-état (ISS), la marge de stabilité robuste et la théorie des ensembles invariants couplée au cadre de la stabilité ISS-Lyapunov [START_REF] Limon | Input-to-state stable MPC for constrained discrete-time nonlinear systems with bounded additive uncertainties[END_REF].

Avantages et inconvénients

La commande prédictive présente une série d'avantages en comparaison avec les autres lois de commandes :

• Formulation simple, basée sur des concepts bien compris et maitrisés.

• Méthodologie permettant des extensions futures.

• Applicabilité à une grande variété de systèmes, y compris les systèmes non linéaires et des systèmes à retard.

• Preuve de la stabilité pour les systèmes linéaires et non linéaires avec des contraintes d'entrée et d'état, sous certaines conditions spécifiques.

• Très utile quand les références futures sont connues a priori.

• Cas multivariable pouvant être facilement traité.

• Extension pour la prise en compte des contraintes (limitation des actionneurs).

• Utilisation explicite du modèle.

• Facilité de maintenance et de mise en oeuvre en cas de changement du modèle ou de ses paramètres.

Cependant, la commande prédictive présente certains inconvénients :

• la commande prédictive nécessite généralement un temps de calcul plus important par rapport aux stratégies de commande classiques. La prise en compte des contraintes de type inégalité augmente encore plus la charge de calcul.

• En général, elle induit un manque de preuve de stabilité (sauf cas particuliers).

Cependant, la plus grande exigence est la nécessité d'un modèle approprié. En effet, la stratégie de commande est fondée sur une connaissance préalable du modèle, mais ne dépend pas d'une structure de modèle spécifique.

Il est cependant évident que les avantages obtenus seront affectés par l'écart existant entre le processus réel et le modèle de prédiction.

Chapitre 4 : Commande prédictive non linéaire

Dans ce chapitre, on propose la formulation et la mise en oeuvre d'une stratégie de commande prédictive non linéaire permettant la poursuite d'une trajectoire de référence.

Formulation du problème d'optimisation

Soit le système non linéaire à temps continu suivant :

ẋ(t) = F (x(t), u(t), θ), x 0 = x y(t) = Hx(t) (1.5) où x ∈ X ⊆ R nx est
θ nom = θ + + θ - 2 (1.7)
et δθ le vecteur des incertitudes paramétriques. L'application

F : R nx × R nu × R n θ -→ R nx ,
x k+1 = f (x k , u k , θ), k ≥ 0, x 0 = x y k = Hx k (1.
u k u(τ ) = u(t k ), τ ∈ [t k , t k+1 [ (1.9)
Soit la trajectoire d'état discrète g, définie comme étant la solution, à l'instant t k+1 , du système (1.5) :

x k+1 = g(t 0 , t k+1 , x, u t k t 0 , θ) y k = Hx k (1.10)
avec l'état initial x 0 , et u t k t 0 la séquence de commande à partir de l'instant initial t 0 jusqu'à l'instant t k . La méthode utilisée pour obtenir g est supposée être la même que celle utilisée pour obtenir (1.8). Il apparait clairement que :

f (x k , u k , θ) ≡ g(t k , t k+1 , x k , u k , θ) (1.11)
Par la suite, le modèle (1.10) sera utilisé dans la stratégie CPNL afin de prédire le comportement futur du système. Dans cette étude, le problème CPNL est formulé à des fins de suivi de trajectoire. L'objectif principal est de contraindre la sortie y à suivre la trajectoire de référence y r tandis que l'entrée de commande est contrainte à suivre la consigne u r . En outre, les saturations sur les vecteurs d'états et d'entrée de commande avec des seuils minimal et maximal x min , x max , u min et u max , respectivement peuvent être prises en considération (i.e. X = [x min , x max ] et U = [u min , u max ]). Ces contraintes inégalité peuvent résulter de contraintes physiques et opérationnelles sur le système commandé.

Le contrôleur prédictif prédit l'évolution future du processus ŷk+Np k+1 sur un horizon temporel fini de longueur N p T s , en utilisant un modèle dynamique non linéaire. À chaque instant t k , la séquence de commande optimale sur l'horizon de prédiction est obtenue en minimisant une fonction coût quadratique fondée sur les erreurs de suivi tout en assurant que les contraintes sont respectées. Seule la première valeur de la séquence optimale est appliquée au système, et toute la procédure est répétée à l'instant d'échantillonnage suivant selon le principe de l'horizon fuyant [START_REF] Maciejowski | Predictive Control: With Constraints[END_REF].

La fonction coût à minimiser est définie comme étant la somme de deux fonctions quadratiques, en fonction des erreurs de suivi sur l'horizon fuyant :

Π NMPC (x k , u k+Np-1 k ) ||u k+Np-1 k -u r,k+Np-1 k || 2 V + ||ŷ k+Np k+1 -y r,k+Np k+1 || 2 W (1.12)
Le modèle de prédiction est donné par :

xk+j = g(t k , t k+j , x k , u k+j-1 k , θ) ŷk+j = H xk+j , ∀j = 1, N p (1.13)
sous les contraintes suivantes reformulées matriciellement :

    I nx 0 nu -I nx 0 nu 0 nx I nu 0 nx -I nu     xk+j u k+j-1 ≤     x max -x min u max -u min     , ∀j = 1, N p (1.14) avec x k le vecteur d'état à l'instant t k , u k+Np-1 k =    u k . . . u k+Np-1    le vecteur d'optimisation, u r,k+Np-1 k =    u r k . . . u r k+Np-1  
 la séquence de commande de référence, N p est la longueur de l'horizon de prédiction. V 0 et W 0 sont les matrices de pondération.

ŷk+Np k+1 =      Hg(t k , t k+1 , x k , u k , θ nom ) Hg(t k , t k+2 , x k , u k+1 k , θ nom ) . . . Hg(t k , t k+Np , x k , u k+Np-1 k , θ nom )      la séquence
En supposant une parfaite connaissance du vecteur des paramètres θ (i.e. θ = θ nom ), la formulation du problème d'optimisation est traduite en un problème de programmation non linéaire sur un horizon de prédiction fini N p T s 13 Résumé à chaque instant d'échantillonnage. La séquence de commande optimale est obtenue en minimisant le critère de performance (1.12) avec prise en compte des contraintes (1.14) et de la dynamique du système comme suit : On constate cependant que l'amélioration en termes de robustesse induite par cette structure reste limitée, nécessitant donc l'élaboration de véritables stratégies robustes.

NMPC : u k+Np-1 k = arg min u k+Np-1 k Π NMPC (x k , u k+Np-1 k ) (1.
RNMPC : u k+Np-1 k = arg min u k+Np-1 k max δθ Π RNMPC (x k , u k+Np-1 k , δθ) (1.18) soumis à                xk+j = g(t k , t k+j , x k , u k+j-1 k , θ = θ nom + δθ), j = 1, N p θ ∈ Θ     I nx 0 nu -I nx 0 nu 0 nx I nu 0 nx -I nu     xk+j u k+j-1 ≤     x max -x min u max -u min     , ∀j = 1, N p (1.19)
où θ et θ nom sont donnés par (1.6) et (1.7). La fonction coût est donnée par :

Π RNMPC (x k , u k+Np-1 k , δθ) ||u k+Np-1 k -u r,k+Np-1 k || 2 V + ||ŷ k+Np k+1 -y r,k+Np k+1 || 2 W (1.20)
Les sorties prédites ŷk+Np k+1 ont la mêmes formulation que dans (1.15) :

ŷk+Np k+1 =      Hg(t k , t k+1 , x k , u k , θ) Hg(t k , t k+2 , x k , u k+1 k , θ) . . . Hg(t k , t k+Np , x k , u k+Np-1 k , θ)      (1.21)

Résumé

La séquence de commande optimale u k+Np-1 k est déterminée en minimisant l'erreur de suivi en considérant toutes les trajectoires possibles [START_REF] Du | Minimax and Applications[END_REF][START_REF] Kasperski | Discrete Optimization with Interval Data: Minmax Regret and Fuzzy Approach[END_REF].

Il apparait clairement que le temps de calcul augmente en fonction de la taille du vecteur des paramètres, du nombre d'entrée de commande et de la taille de l'horizon de prédiction pouvant vite devenir prohibitif. Ainsi, le défi principal est de réduire le temps de calcul tout en maintenant de bonnes performances en termes de précision de poursuite.

Commande prédictive robuste réduite

Étant donné que le problème d'optimisation min-max (1.18)-(1.20) est coûteux en termes de temps de calcul, il peut se simplifier en réduisant le nombre de paramètres θ incertains pris en compte dans l'optimisation à partir d'une analyse de sensibilité du modèle vis-à-vis de ses paramètres.

Analyse de sensibilité

La fonction de sensibilité S x i θ j représente la sensibilité de chaque état x i aux faibles variations de chaque paramètre θ j . Elle est donnée par l'expression suivante :

S x i θ j ∂x i ∂θ j , i = 1, n x et j = 1, n θ (1.22)
La dynamique des fonctions de sensibilité pour le système (1.5) est calculée comme suit :

Ṡx i θ j = d dt ∂x i ∂θ j = ∂ ∂θ j dx i dt = ∂F i ∂θ j + nx k=1 ∂F i ∂x k ∂x k ∂θ j (1.23)
avec la condition initiale suivante : ∂x i ∂θ j = 0. L'analyse de l'évolution temporelle des fonctions de sensibilité selon l'ordre de grandeur de leurs amplitudes, permet de sélectionner les paramètres κ qui ont le plus d'influence.

Reformulation du problème

Par la suite, seuls les paramètres les plus influents vont être considérés dans le problème d'optimisation min-max, au lieu de tous les paramètres (avec θ [κ, ζ]). Les paramètres influents κ sont définis par :

   κ = κ nom + δκ (1.24)
κ nom = κ + + κ - 2 (1.
rRNMPC : u k+Np-1 k = arg min u k+Np-1 k max δκ Π rRNMPC (x k , u k+Np-1 k , δκ) (1.26) soumis à                    xk+j = g(t k , t k+j , x k , u k+j-1 k , κ nom + δκ ζ nom ), j = 1, N p κ ∈ Θ κ     I nx 0 nu -I nx 0 nu 0 nx I nu 0 nx -I nu     xk+j u k+j-1 ≤     x max -x min u max -u min     , ∀j = 1, N p (1.27)
La nouvelle fonction coût est définie par :

Π rRNMPC (u k+Np-1 k , δκ) ||u k+Np-1 k -u r,k+Np-1 k || 2 V + ||ŷ k+Np k+1 -y r,k+Np k+1 || 2 W (1.28)
Les sorties prédites ŷk+Np k+1 sont données par :

ŷk+Np k+1 =      Hg(t k , t k+1 , x k , u k , [κ , ζ nom ] ) Hg(t k , t k+2 , x k , u k+1 k , [κ , ζ nom ] ) . . . Hg(t k , t k+Np , x k , u k+Np-1 k , [κ , ζ nom ] )      (1.29)
Néanmoins, cette approche ne peut pas être appliquée pour les systèmes dont les paramètres ont une influence similaire, aucune réduction des paramètres n'est alors possible. Par conséquent, cette approche ne peut pas être utilisée dans tous les cas. , pour j = 1, N p :

xk+j ≈ xnom,k+j + ∇ θ g(t k+j )δθ + ∇ u g(t k+j ) u k+j-1 k -u r,k+j-1 k (1.30) avec                              xnom,k+j = g(t k , t k+j , x k , u r,k+j-1 k , θ nom ) (1.31) ∇ θ g(t k+j ) = ∂g(t k , t k+j , x k , u k+j-1 k , θ) ∂θ u k+j-1 k = u r,k+j-1 k θ = θnom (1.32) ∇ u g(t k+j ) = ∂g(t k , t k+j , x k , u k+j-1 k , θ) ∂u k+j-1 k u k+j-1 k = u r,k+j-1 k θ = θnom (1.33)
Ainsi, les sorties prédites sont données par l'expression suivante :

ŷk+Np k+1 = G k+Np nom,k+1 + G k+Np θ,k+1 δθ + G k+Np-1 u,k (u k+Np-1 k -u r,k+Np-1 k ) (1.34) avec G k+Np nom,k+1 =         H xnom,k+1 . . . H xnom,k+j . . . H xnom,k+Np         , G k+Np θ,k+1 =         H∇ θ g(t k+1 ) . . . H∇ θ g(t k+j ) . . . H∇ θ g(t k+Np )         , G k+Np-1 u,k =         H∇ u g(t k+1 ) . . . H∇ u g(t k+j ) . . . H∇ u g(t k+N p )        
En utilisant l'équation (1.34), la fonction coût initiale Π RNMPC (1.20) est approchée par :

Π(x k , u k+Np-1 k , δθ) ≈ ||u k+Np-1 k -u r,k+Np-1 k || 2 V + ||G k+Np nom,k+1 -y r,k+Np k+1 + G k+Np θ,k+1 δθ + G k+Np-1 u,k (u k+Np-1 k -u r,k+Np-1 k )|| 2 W (1.35)
Finalement, le problème d'optimisation (1.18) devient :

u k+Np-1 k = arg min u k+Np-1 k max δθ Π(x k , u k+Np-1 k , δθ) (1.36) soumis à θ ∈ Θ, x ∈ X, u ∈ U (1.
x k+1 = l(x k , w k ), k ≥ 0, x 0 = x (1.
x k+1 = f (x k , u k , θ nom + δθ s ) f (x k , u r k , θ nom ) + ∇ u f (u r k , θ nom ).(u k -u r k )+ ∇ θ f (u r k , θ nom ).δθ s + ϑ(|u k -u r k | 2 ) + ϑ(|δθ s | 2 ) (1.39) où ϑ(|.| 2 )
est le reste du développement en série de Taylor. Soit le modèle de prédiction à l'instant t k+1 :

xk+1|k = f p (x k , u k , θ nom + δθ) f (x k , u r k , θ nom ) + ∇ u f (u r k , θ nom ).(u k -u r k ) + ∇ θ f (u r k , θ nom ).δθ (1.40)
Après développement, on obtient une borne sur l'erreur de prédiction :

∃α ∈ R + , ∀x, ∀u, ∀θ, |∇ θ f | ≤ α (1.41) et ∃η ∈ R + , tel que η = max(|ϑ(|u k -u r k | 2 )+ϑ(|δθ s | 2 )|, max(|δθ|, |δθ s |)) (1.42) alors |x k+1 -xk+1|k | ≤ (2α + 1)η Λ(η) (1.43) où Λ est une fonction K ∞ .
La fonction coût optimale Π peut être bornée comme suit :

α ψ (|x k |) ≤ Π(x k , u, δ θ) ≤ β t (|x k |) + ϕ(η) + N p χ(η) (1.44) où α ψ , β t , et χ sont des fonctions K ∞ et ϕ est une fonction K.
De plus, la variation de la fonction coût ∆Π est donnée par l'expression suivante : Le critère J(λ) est défini par l'expression suivante :

∆Π Π(x k+1 , ȗ, δ θ) -Π(x k , u, δ θ) ≤ -α ψ (|x|) + χ(η) (1.45) avec χ(η) = χ(η) + L tx L Np-1 f x + L ψx Np-2 j=0 L j f x Λ(η) (1.
J(λ) = ||z(λ)|| 2 V + λ||δθ max || 2 + ||G k+Np-1 u,k z(λ) -y r,k+Np k+1 + G k+Np nom,k+1 || 2 W (λ) (1.48) avec z(λ) = (V + G k+Np-1 u,k W (λ)G k+Np-1 u,k ) † (G k+Np-1 u,k W (λ)(y r,k+Np k+1 -G k+Np nom,k+1 )) (1.49) et W (λ) = W + W G k+Np θ,k+1 (λI -G k+N p θ,k+1 W G k+Np θ,k+1 ) † G k+N p θ,k+1 W (1.50)
Étape 2. La séquence de commande u k+Np-1 k est donnée par :

u k+Np-1 k =u r,k+Np-1 k + (V + G k+Np-1 u,k W (λ )G k+Np-1 u,k ) † × (G k+Np-1 u,k W (λ )(y r,k+Np k+1 -G k+Np nom,k+1 )) (1.
z(λ) = arg min z≤z≤z z E(λ)z -2B(λ) z (1.53) où              E(λ) = V + G k+Np-1 u,k W (λ)G k+Np-1 u,k (1.54)
B(λ) = G k+Np-1 u,k W (λ) y r,k+Np k+1 -G k+Np nom,k+1 (1.55)
z = u min I nu -u r,k+Np-1 k (1.56) z = u max I nu -u r,k+Np-1 k (1.
u k+Np-1 k = u r,k+Np-1 k + z(λ ) (1.
δ u k+Np-1 k = arg min δu k+Np-1 k ∈U δ max δθ∈Θ δ Π δ (x k , δu k+Np-1 k , δθ) (1.59) soumis à xk+j = g(t k , t k+j , x k , u k+j-1 k , θ = θ nom + δθ), j = 1, N p (1.60) où Π δ (x k , δu k+Np-1 k , δθ) ||δu k+Np-1 k || 2 V + ||ŷ k+Np k -y r,k+Np k || 2 W (1.61)
Les variables d'optimisation δu k+Np-1 k sont définies comme suit :

δu k+Np-1 k =         u k -u k-1 . . . u k+j -u k+j-1 . . . u k+Np-1 -u k+Np-2                δu k . . . δu k+j . . . δu k+Np-1        (1.62) avec u k-1 l'entrée de commande appliquée à l'instant k -1 (solution du problème d'optimisation (1.59) à l'instant k -1).
Les sorties prédites sont données par l'expression suivante :

ŷk+Np k+1 = G k+Np nom,k+1 + G k+Np θ,k+1 δθ + G k+Np-1 u,k (Ξ k+Np-1 k + T Np δu k+Np-1 k ) (1.63) avec Ξ k+Np-1 k =    u k-1 . . . u k-1    -    u k-1 . . . u k+Np-2    (
z(λ) = V + T Np G k+Np-1 u,k W (λ)G k+Np-1 u,k T Np † T Np G k+Np-1 u,k W (λ) × ȳk+Np k+1 -G k+Np nom,k+1 -G k+Np-1 u,k Ξ k+Np-1 k (1.66)
Étape 2. La séquence de commande u k+Np-1 k est donnée par : La boucle de régulation est formée par une loi de commande prédictive robuste (entrée de commande notée u ) combinée à un correcteur auxiliaire (entrée de commande notée ũ). Le régulateur prédictif permet de suivre la trajectoire de référence, alors que le régulateur supplémentaire annule toute erreur de suivi résiduelle. Deux types de régulation additionnelle sont étudiés : une structure proportionnelle-intégrale (PI) et une stratégie de commande par mode glissant intégrateur (ISM).

u k+Np-1 k =    u k-1 . . . u k-1    + T Np z(λ ) (1.67
En premier lieu, la loi de commande proportionnel intégrateur est choisie pour sa simplicité. Elle est donnée par l'expression suivante :

ũ(t) = K p (y(t) -ŷ(t)) + K p T i t t k-1 (y(τ ) -ŷ(τ ))dτ (1.68)
où K p et T i sont les deux paramètres de réglage du PI et ŷ la prédiction du modèle.

En second lieu, on considère une loi de commande comme régulateur auxiliaire. La loi de commande par mode glissant intégrateur présente l'avantage d'être robuste vis-à-vis des perturbations. La surface de glissement est choisie comme suit :

φ(x, t) = Z 2 (t) + ξ 3 Z 1 (t) (1.69) avec      Z 1 (t) = t t k-1 (y(τ ) -ŷ(τ )) dτ Z 2 (t) = 1 ξ 1 (y(t) -ŷ(t) -ξ 2 Z 1 (t)) (1.70)
L'attractivité de la surface de glissement est donnée ci-dessous :

φ(x, t) = -K s sign(φ(x, t)), K s ∈ R + (1.71)
Ainsi, l'expression de la loi de commande par mode glissant intégrateur (ISM) est la suivante :

ũ(t) = 1 η (-ϕ -χu k-1 -ξ 1 K s sign(φ(x, t)) + (ξ 2 -ξ 1 ξ 3 )(ξ 1 Z 2 (t) + ξ 2 Z 1 (t))) (1.72) avec    ϕ = H(f x (x(t), θ) -f x (x(t), θ nom )), χ = H(f u (x(t), θ) -f u (x(t), θ nom )), η = Hf u (x(t), θ) (1.73) et ẋ(t) = f x (x(t), θ) + f u (x(t), θ)u(t), ∀t > t 0 , x(t 0 ) = x 0 y(t) = Hx(t) (1

Modélisation du système

Le modèle dynamique considéré suppose que le photobioréacteur est en mode continu (le débit moyen de soutirage est égal à celui d'alimentation, conduisant à un volume constant), sans aucune biomasse supplémentaire dans l'alimentation et en négligeant l'effet des échanges gazeux. Les bilans de matières conduisent aux équations différentielles suivantes :

   Ẋ(t) = µ(Q(t), I(t)) X(t) -DX(t) Q(t) = ρ(S(t)) -µ(Q(t), I(t)) Q(t) Ṡ(t) = (S in -S(t))D -ρ(S(t))X(t) (1.75) avec                  µ(Q, I) = μ(1 - K Q Q )µ I (Droop) ρ(S) = ρ m S S + K s (Michaelis-Menten) µ I = I I + K sI + I 2 K iI (Haldane) (1.76)
où X est la concentration en biomasse (µm 3 L -1 ), Q le quota interne (µmolµm -3 ) et S la concentration en substrat (µmol L -1 ). D (d -1 ) est le taux de dilution défini comme le ratio entre le débit et le volume (d -1 ), µ(.) est la vitesse de croissance (d -1 ), μ le taux de croissance maximal, K Q la quantité de quota interne en-dessous de laquelle les bactéries ne croissent plus (µmol µm -3 ), ρ(.) la vitesse d'absorption (µmol µm -3 d -1 ), ρ m le taux d'absorption maximal (µmol µm -3 d -1 ) et K s la constante de saturation par le substrat (µmol L -1 ). La fonction µ I représente l'effet de la lumière sur le taux de croissance (-), I est l'intensité lumineuse (µE m -2 s -1 ), K sI la constante de saturation par la lumière (µE m -2 s -1 ) et K iI la constante d'inhibition par la lumière (µE m -2 s -1 ). Les paramètres du modèle (1.75) sont donnés dans le tableau suivant :

Tableau 1.1 -Paramètres du modèle. Paramètre Valeur Unité μ 2 d -1 ρ m 9,3 µmol µm -3 d -1 K Q 1,8 µmol µm -3 K s 0,105 µmol L -1 K sI 150 µE m -2 s -1 K iI 2000 µE m -2 s -1

Stratégie de commande

L'objectif de commande est de réguler la concentration en biomasse à une valeur de consigne en agissant sur le taux de dilution tout en respectant les contraintes biologiques.

Résultats en simulation

Dans cette partie, les performances des stratégies de commande prédictives sont évaluées dans le cas de paramètres incertains et en présence de bruit de mesure. Les conditions de simulation sont résumées dans le tableau 1.2. Les quatre lois de commande prédictives sont tout d'abord comparées : la commande prédictive non linéaire (NMPC), la formulation min-max robuste (RNMPC), la formulation min-max réduite (rRNMPC) et la commande prédictive robuste linéarisée (LRMPC) sans prise en compte des contraintes explicitement et saturation de la commande a posteriori. Grâce à l'analyse de sensibilité du modèle vis-à-vis de ses paramètres, il apparait que la constante de saturation par le substrat K s et la quantité de quota interne minimale K Q sont les paramètres qui ont le plus d'influence sur l'évolution de la concentration en biomasse (κ Les mesures de la concentration en biomasse sont entachées d'un bruit blanc gaussien de moyenne nulle et d'écart type égale à 0,1 et l'intensité lumineuse I est égale à I opt = √ K sI K iI . De plus, afin de tester la robustesse des algorithmes de commande par rapport aux incertitudes du modèle, les paramètres du système réel sont incertains de 30% au maximum. Les figures 1.7 et 1.8 présentent l'évolution de la concentration en biomasse et du taux de dilution. Les résultats obtenus montrent que la diminution du taux de dilution entraine une augmentation de la concentration en biomasse et vice versa. D'autre part, on peut remarquer aussi que les lois de commande (r)RNMPC et LRMPC ont de meilleures performances que la commande prédictive NMPC en termes de précision de poursuite. Dans le cas de la NMPC, la présence d'une grande erreur de poursuite est due au fait que la dispersion entre le modèle et le système n'est pas prise en compte dans l'étape de prédiction pendant l'optimisation. Par ailleurs, la RNMPC régule mieux la concentration en biomasse que la rRNMPC et LRMPC mais induit une grande charge de calcul comme indiqué dans le tableau 1.3. En revanche, la méthode rRNMPC permet de réduire le temps de calcul mais on constate une nette diminution de la qualité de poursuite. Enfin, la stratégie LRMPC permet d'avoir le meilleur compromis entre la précision de poursuite et le temps de calcul. • Deux améliorations de la loi de commande prédictive robuste linéarisée ont été proposées dans la troisième partie. Tout d'abord, la linéarisation des sorties prédites autour des valeurs nominales des paramètres du modèle et de la séquence de commande optimale obtenue à l'itération précédente. Cette modification permet de rendre la solution moins sensible aux bruits de mesure. Ensuite, une structure de commande hiérarchisée combinant la commande prédictive robuste linéarisée et une commande auxiliaire (proportionnel intégrateur ou mode glissant intégrateur) est développée afin d'éliminer toute erreur statique en suivi de trajectoire de référence.

= [K s , K Q ] ).
• Le dernier chapitre est dédié à l'application des lois de commande développées à un cas d'étude de commande de photobioréacteur de culture de microalgues pour la régulation de la concentration en biomasse en présence d'incertitudes paramétriques, bruit de mesure et fluctuations de l'intensité lumineuse.

Perspectives

Dans la continuité de ces travaux, les perspectives suivantes peuvent être envisagées :

• Application des lois de commande proposées au cas d'un système multivariable comprenant un nombre de paramètres incertains plus important.

• Extension de la formulation du problème d'optimisation aux cas des perturbations et des cinétiques incertaines.

• Synthèse d'un observateur afin d'estimer le vecteur d'état à partir des mesures.

• Amélioration de la qualité du modèle linéarisé en utilisant un développement en série de Taylor limité au second ordre au lieu du premier ordre.

• Étude de l'influence de la résolution du problème d'optimisation en utilisant la dualité lagrangienne sur la stabilité du système en boucle fermée.

• Analyse de la stabilité du système bouclé avec la stratégie de commande hiérarchisée.

• Détermination des conditions à satisfaire pour assurer la stabilité du système avec le couplage régulateur/estimateur.

Chapter 2 Introduction

Context and motivations

The control of nonlinear systems subject to physical constraints on the input and state is undoubtly a challenging and important problem. The well known systematic nonlinear control methods such as feedback linearization and constructive Lyapunov-based methods lead to very elegant solutions, but they depend on a complicated design procedure that does not scale well to large systems and they cannot handle constraints easily or in a systematic manner. Based on this, the concept of optimal control and in particular the Model Predictive Control (MPC) approach and its nonlinear version NMPC appears to be an attractive alternative since the complexity of the control design increases moderately with the size and complexity of the system. The NMPC strategy is therefore put forward because of its ability to deal with uncertainties and constraints. The essence of NMPC is to optimize, over an open-loop time sequence of controls, the process response using a model of the system to forecast the future process behavior over a prediction horizon. Moreover, even if NMPC strategies have proved to be efficient in many industrial applications due to their ability to operate the process safely under physical constraints, the monitoring of complex, highly nonlinear, uncertain systems becomes more and more a delicate task. This is the case for example of bioprocesses, that will be the application field considered in this work. Indeed, in such cases, model parameters are generally known with a confidence interval only (determined from an identification procedure for example). The main objective is thus to elaborate an adequate robust control strategy in order to guarantee that the process will yield the reference trajectory under model parameter uncertainties.

Therefore, NMPC strategies must be extended to provide robustness fea-tures, developing robust control strategies, which can compensate for the lack of parameters information and/or accuracy. There are two popular alternatives for making decision with incomplete knowledge: the stochastic solution (probabilistic theory) and the min-max solution (game theory), which appear as expensive approaches. Indeed, the robust NMPC (RNMPC) law can be formulated as a nonlinear min-max optimization problem where the effect of the uncertainties must be taken into account in the design procedure. However, this approach tends to become too complex to be solved numerically online. Consequently, the total calculation time is an important factor that must be reduced as much as possible.

As a consequence, the challenge that must be taken up in this work is to design a MPC controller that would be robust against unknown but bounded model parameters uncertainties and computationally more tractable in calculating the optimal control, which makes it suitable for online implementation. For that purpose, starting from the standard RNMPC formulations, new structures will be developed to match this compromise between robustness and computational load, with a dedicated application to bioprocesses.

Outline

We briefly summarize the main results of this thesis, chapter by chapter, as the following:

• Chapter 3: In this chapter, a state of the art of existing Model Predictive Control (MPC) techniques is given with a focus on nonlinear control problems. It starts with the basic MPC problem formulation and provides a review of linear and nonlinear MPC schemes, addresses robustness issues and presents briefly nominal and robust stability concepts.

• Chapter 4: This chapter investigates the Nonlinear MPC formulation in order to solve a trajectory tracking problem for continuous/discretetime nonlinear systems. A continuous/discrete formulation is considered in this thesis, since models of bioprocesses are time-continuous but with discrete-time measurements. A variant of NMPC is presented in section 4.3. It consists in considering the accumulation of the gap between the system and the model at each sampling time over prediction intervals. An illustrative example (for a bioprocess) is presented in section 4.4 and numerical results are provided in order to demonstrate the effectiveness of the NMPC strategy and highlight its drawbacks and limitations, which motivate the development of RNMPC strategies.

• Chapter 5: This chapter deals with a robust formulation of the NMPC in order to take into account the model parameters mismatch. The robust NMPC (RNMPC) law is formulated as a nonlinear min-max optimization problem. However, this approach tends to become too complex to be solved numerically online. To overcome this problem, two major solutions have been proposed. The first one consists in reducing the number of uncertain parameters through a sensitivity analysis.

Only the main influencial parameters on the model are considered in the min-max problem whereas the other parameters are kept constant and equal to their nominal values. Secondly, a new approach has been proposed, based on the linearization of the predicted trajectory to turn the original optimization problem into a more tractable one. It consists in simplifying the min-max optimization problem through a sensitivity analysis of the model with respect to its parameters, in order to reduce the calculation time, which makes it suitable for online implementation. The criterion is linearized around the model parameters nominal values and the reference control sequence (control inputs over the prediction horizon). The derived optimization problem is an unidimensionnal optimization problem (referred to as LRMPC). Moreover, stability analysis of the closed-loop system while applying the LRMPC law is established under some assumptions. Furthermore, the inequality constraints on the optimization variables are taken into account in the LRMPC problem formulation. This leads to a bilevel problem (referred to as CLRMPC), with a scalar optimization problem in the upper level, and a quadratic programming one in the lower level. The performances of the proposed control strategies will be illustrated in simulations (the same system as in chapter 4).

• Chapter 6: This chapter proposes two improvements for the (C)LRMPC law. The first one consists in linearizing the model at each sampling time along the nominal trajectory which is defined based on the nominal parameter values and the optimal control sequence obtained at the previous iteration, instead of reference control sequence as done in the previous chapter. This modification is motivated by the fact that the optimal control sequence obtained at the previous iteration is non-model-based. As a consequence, the proposed approach appears to be less sensitive against measurement noise thanks to the penality term on the control evolution. The derived optimization prob-lem is a bilevel one or a scalar optimization problem depending on the considered constraints. However, the proposed control strategies have accuracy issue as they deal with a first-order approximation of the nonlinear system. A hierarchical control strategy is proposed in order to cancel the residual tracking error due to the linearization drawback and possible disturbances. The proposed control strategy combines a robust model predictive control law with a Proportional Integral (PI) law or Integral Sliding Mode (ISM) controller. The predictive controller guarantees the tracking of the reference trajectory, whereas the other regulator ensures a good tracking accuracy. These two improvements are applied to the illustrative example.

• Chapter 7: In this chapter, all proposed approaches in the previous chapters are applied to a case study to control of microalgae culture in a continuous photobioreactor. The aim is to regulate the microalgae concentrations at a given reference value, despite disturbances (light intensity fluctuation) and model parameters uncertainties. First, the model of the bioprocess is described and analysed. Secondly, the determination of the equilibrium is addressed. Finally, simulation results are given for reference trajectory tracking and disturbance rejection.

• Chapter 8: The last chapter summarizes the developed work in this PhD thesis and proposes some recommendations for future directions.

List of Publications

The thesis work has resulted in several accepted/submitted publications to international journal and conferences. They are given hereafter:

Journal paper

Chapter 3

Predictive control: State of the art and main strategies

Introduction

The purpose of this chapter is to give an introduction to one of the most powerful strategies that can be used to address nonlinear control problems. We present the basic MPC problem formulation and provide a review on previous works that appear in the literature. The key advantages/drawbacks of MPC are outlined and some of the theoretical, computational and implementation aspects are discussed. It also serves as a centralized literature study for the following chapters. In any case, this introduction does not represent an exhaustive overview of predictive control, it however provides the reader with the sufficient background for the next chapters.

Background

The nonlinear optimal control theory was developed in the 1950's and 1960's, resulting from the maximum principle of Pontryagin [START_REF] Pontryagin | The Mathematical Theory of Optimal Processes[END_REF] and dynamic programming method developed by Bellman [START_REF] Bellman | The theory of dynamic programming[END_REF]. The receding horizon principle was proposed by Propoi [START_REF] Propoi | [END_REF], within the frame of open loop optimal feedback. The philosophy of the receding horizon problem is closely related to the minimum time optimal control problem and to linear programming [START_REF] Zadeh | On optimal contro and linear programming[END_REF].

In the late 1970s, due to the progress in algorithms for solving constrained linear and quadratic optimization problems, linear MPC, also referred to as receding or moving horizon control, became popular in control engineering. As a consequence, the interest in MPC for process control grew quickly. Many works appeared principally Richalet et al. [START_REF] Richalet | Model predictive heuristic control: Application to industrial processes[END_REF] and Cutler et al. [START_REF] Cutler | Dynamic matrix control-a computer control algorithm[END_REF] were the first in the control community to propose industrial applications of MPC strategy. MPC has become the de-facto standard advanced control method in the process industries due to its ability to handle large scale multivariable processes. Many systems are however in general inherently nonlinear. In this case, linear models are often inadequate to describe the process dynamics and nonlinear models have to be used.

The practical interest is driven by the fact that process nonlinearities and contraints are explicitly considered. This motivates the use of nonlinear model predictive control (henceforth abbreviated as NMPC), the extension of well established linear MPC. MPC involves the repetitive solution of an optimal control problem at each sampling instant in a receding horizon window. The key challenge that has to be met is to ensure that the implementation of a sequence of open-loop optimal control values will be stable, when considering the closed loop system. Subsequently, NMPC problems have witnessed a steadily increasing attention from academic control theoretists. This leads to increasing advances with focus on design that guarantee stability and robustness. A lot of pioneer works have emerged in this topic: Chen and Shaw [START_REF] Chen | On receding horizon feedback control[END_REF], Keerthi and Gilbert [START_REF] Keerthi | Optimal infinite-horizon feedback laws for a general class of constrained discrete-time systems: Stability and moving-horizon approximations[END_REF], Michalska and Mayne [START_REF] Michalska | Robust receding horizon control of constrained nonlinear systems[END_REF], Parisini and Zoppoli [START_REF] Parisini | A receding-horizon regulator for nonlinear systems and a neural approximation[END_REF], Alamir and Bornard [START_REF] Alamir | Stability of a truncated infinite constrained receding horizon scheme: the general discrete nonlinear case[END_REF], Magni and Sepulchre [START_REF] Magni | Stability margins of nonlinear recedinghorizon control via inverse optimality[END_REF], Chen and Allgöwer [START_REF] Chen | A quasi-infinite horizon nonlinear model predictive control scheme with guaranteed stability[END_REF], De Nicolao et al. [START_REF] De Nicolao | Stabilizing recedin-horizon control of nonlinear time-varying systems[END_REF], Mayne et al. [START_REF] Mayne | Constrained model predictive control: Stability and optimality[END_REF], and Jadbabaie and Hauser [START_REF] Jadbabaie | On the stability of receding horizon control with a general terminal cost[END_REF]. In addition to the research on NMPC, researchers have developed socalled Real Time Optimization (RTO) [START_REF] Chachuat | Adaptation strategies for real-time optimization[END_REF] approaches which are similar to NMPC, as they are generally based on nonlinear programming, and use static nonlinear models, while NMPC uses dynamic nonlinear models.

Receding Horizon Control (RHC) is cited as one of the most popular advanced techniques for industrial process applications and has been widely adopted in the field of process control, due to the simplicity of the algorithm. MPC has induced a significant impact on industrial control engineering world, as a consequence, there are many applications of predictive control strategy in very various domains, e.g. robot manipulators [START_REF] Lelić | Generalized pole-placement self-tuning controller part 1, basic algorithm[END_REF], clinical anesthesia [START_REF] Linkers | Advances in Model-Based Predictive Control, chapter Generalized Predictive Control in Clinical Anesthesia[END_REF], cement industry and pulp factories [START_REF] Sánchez | Adaptive Predictive Control: From the Concepts to Plant Optimization[END_REF], drying towers and robot arms [START_REF] Clarke | Application of generalized predictive control to industrial processes[END_REF], distillation columns [START_REF] Huang | Advanced step nonlinear model predictive control for air separation units[END_REF][START_REF] Richalet | Industrial applications of model based predictive control[END_REF], Polyvinyl chloride (PVC) plants, steam generators and oils refining [START_REF] Richalet | Model predictive heuristic control: Application to industrial processes[END_REF], solar power plant [START_REF] Camacho | Model Predictive Control[END_REF], thermomechanical pulping [START_REF] Harinath | Control and optimization strategies for thermo-mechanical pulping processes: Nonlinear model predictive control[END_REF][START_REF] Nagy | Robust nonlinear model predictive contro of fedbatches processes[END_REF], biochemical systems [START_REF] Ashoori | Optimal control of a nonlinear fed-batch fermentation process using model predictive approach[END_REF][START_REF] Tebbani | Nonlinear predictive control for maximization of CO 2 bio-fixation by microalgae in a photobioreactor[END_REF][START_REF] Tebbani | Nonlinear Predictive Control of fed-batch Cultures of E. coli[END_REF], motor control and food extruder process [START_REF] Wang | Model Predictive Control System Design and Implementation Using MATLAB R[END_REF], petrochemical sector [START_REF] García | Model predictive control: Theory and practice-a survey[END_REF], aerospace [START_REF] Hartley | Predictive control using an FPGA with application to aircraft control[END_REF][START_REF] Bourgeois | Launcher atmospheric guidance based on nonlinear model predictive control[END_REF][START_REF] Crassidis | Nonlinear predictive control of spacecraft[END_REF], automotive [START_REF] Re | Automotive Model Predictive Control: Models, Methods and Applications[END_REF], mining metallurgy [START_REF] Balchen | State space model predictive control of a multi stage electro-metallurgical process[END_REF][START_REF] Nguyen | Nonlinear model predictive control of steel slab walkingbeam reheating furnace based on a numerical model[END_REF], etc.

MPC is a powerful approach of great promise that has proven itself in several applications for many classes of systems. On the other side, nonlinear MPC is limited in its industrial impact due to the challenges of guaranteeing a global solution to the resulting nonlinear optimization problem within the real time requirements. For additional information about the use of linear and nonlinear MPC in practice, applications can be found in many articles and books, e.g. in [START_REF] Gyurkovics | Assessment and Future Directions of Nonlinear Model Predictive Control, chapter Conditions for MPC Based Stabilization of Sampled-Data Nonlinear Systems Via Discrete-Time Approximations[END_REF][START_REF] Allgöwer | Nonlinear Model Predictive Control. Progress in Systems and Control Theory[END_REF][START_REF] Kouvaritakis | Non-linear Predictive Control: Theory and Practice[END_REF]. It is not so easy to list all referenced MPC applications, more can be found in [START_REF] Qin | A survey of industrial model predictive control technology[END_REF].

The basic principles of MPC

We consider a general discrete-time, state-space nonlinear model of the plant, in the form

x k+1 = f (x k , u k ) (3.1)
where x k ∈ R nx and u k ∈ R nu are the plant state and control action vectors at time step k, respectively. Given the plant (3.1), a general formulation of MPC problem can be described as the following min

u k ,...,u k+Np-1 Np i=1 φ(x k+i , u k+i-1 ) (3.2) s.t xk+i+1 = f (x k+i , u k+i ), for i = 0, N p -1 xk = x k (3.3)
where

• N p is the prediction horizon.

• xk+i is the predicted state value.

• x k is the initial condition which is the plant state at time step k.

• φ(., .) is the cost per stage that is generally a quadratic term that minimizes the difference between the predicted state and the setpoint. However, other norms can be considered such as the L 1 or infinity norms.

The MPC controller is implemented in a moving horizon framework. At current time step k, the plant state x k is used as the initial condition and the optimization problem

Prediction model

A system model is a mathematical representation of the reality. It aims at predicting the time evolution of each state variable. Generally, models are related to both parameters, input and output signals. The explicit use of the model (3.3) is determined by the necessity to calculate the predicted states xk+i for a determined horizon of length N p . These predicted states depend on the known values up to time instant k (i.e. x k ) and on the future control sequence, u k+i , i = 0, N p -1, which are the optimization variables.

The process model plays, as a consequence, a decisive role in the controller performance. Moreover, the design of the model should contain the necessary mechanisms for obtaining the best possible model while being simple to be implemented and to be understood. Various models can be considered, the following being the most commonly used [START_REF] Camacho | Model Predictive Control[END_REF]: impulse response, step response, transfer function, state space, neural networks and fuzzy models.

Cost function

For setpoint tracking MPC, the index to be optimized is a quadratic (nonnegative) function measuring the distance between the predicted model state x and a determined known reference sequence x r over the horizon window, plus a quadratic function weighting the control effort necessary to achieve the specified goal.

If the quadratic norm is used (which is usually the case), the general expression for such an objective function can be:

Np i=1 φ(x k+i , u k+i-1 ) = Np i=1 ||x k+i -x r k+i || 2 α i + Nc i=1 ||∆u k+i-1 || 2 β i (3.4) with ∆u k+j = 0, j = N c , N p -1 and ∆u k+i-1 = u k+i-1 -u k+i-2 .
The criterion usually considers the following tuning parameters:

• N p and N c are the length of the prediction horizon and the control horizon respectively, which are not necessarily equal (i.e. N p ≥ N c ≥ 1) in order to reduce the computational load.

• x r denotes the reference state trajectory which is used to specify closedloop behavior and tracking performances.

• α i ∈ R nx×nx and β i ∈ R nu×nu are the weighting matrices for the state tracking error and the control, respectively. The form of the cost function in (3.4) implies that the deviations of the predicted controlled states x from a reference x r and the control moves are penalized at every point in the prediction horizon.

A conceptual picture of the prediction horizon is shown in Figure 3.2. These tuning parameters should be adjusted to avoid unstability and to give satisfactory dynamic performance, since all affect the behavior of the closedloop combination of system and predictive controller. It should be noted that it is possible to use a non-quadratic cost, so that e.g. absolute values of errors can be penalized, rather than square values as detailed in [START_REF] Maciejowski | Predictive Control: With Constraints[END_REF].

Future/Prediction Past

In addition, constraints can be taken into account in the optimization problem. For example, the actuators may have a limited field of action. By adding these constraints to the optimization problem, the latter becomes more complex, so that the solution can not be obtained explicitly as in the unconstrained case. Moreover, in this case, the paradigm of finding a feasible solution to the optimization problem has to be addressed.

Control law calculation

The optimizer is a fundamental part of the MPC strategy as it provides the control actions by minimizing the cost function (3.4). An analytic solution can be obtained for the quadratic criterion case if the model is linear and there are no constraints. Otherwise, the control sequence results from a real-time optimization strategy, involving solving a nonlinear programming problem (NLP) (for more details see section 3.5). The size of the optimization problem depends on the number of variables and on the prediction horizon considered. Summarizing, the basic MPC scheme works as follows:

1. Obtain measurements 2. Compute the model prediction over a certain prediction horizon 3. Get an optimal input signal by minimizing a given cost function 4. Implement the first value of the optimal input signal until new measurements are available.

Continue with 1.

While it can be desirable for computational and performances reasons to choose unequal lengths of the prediction and control horizon, we assume in the following that N p = N c for the study.

Literature review

In general, MPC can be split into two different classes of methodologies: linear and nonlinear. Maciejowski and Camacho provide a good overview of linear theory [START_REF] Maciejowski | Predictive Control: With Constraints[END_REF][START_REF] Camacho | Model Predictive Control[END_REF], while Allgöwer and Zheng provide an overview of nonlinear methods [START_REF] Allgöwer | Nonlinear Model Predictive Control. Progress in Systems and Control Theory[END_REF]. The various MPC algorithms result from the way the cost function to be minimized and the constraints of the system are specified. This section aims at briefly presenting the main MPC strategies found in the literature, for each of them the reader is referred to appropriate references for more details.

Linear case

Generally, linear MPC uses a low order linear model to represent the system. In case of linear constraints, it has the advantage of solving small optimization problem, which can be done fast enough at each sampling time in order to be implemented in the moving horizon framework. In this case, the solution of the optimization problem is the solution of a quadratic or linear program (QP/LP), which is known to be convex and for which there exists a variety of numerical methods and softwares.

Remark 3.1. An analytic solution can be obtained for the quadratic criterion if the model is linear and there are no constraints.

On the other hand, linear MPC also suffers from drawbacks such as plantmodel mismatch because the model may be only valid around an operating point. Depending on different modelling and solution strategies, there are many variations of linear MPC, methods considered to be the most representative will be briefly cited [START_REF] Camacho | Model Predictive Control[END_REF]:

• Dynamic Matrix Control (DMC) [START_REF] Cutler | Dynamic matrix control-a computer control algorithm[END_REF]: the DMC uses a step response model of the process. The output of the plant is then expressed through the time response parameters. The setpoint trajectory can be approached by a smoother trajectory (called projected desire trajectory). Shell engineers developed QDMC by posing the DMC problem as a quadratic programming (QP) problem in which constraints appear explicitly. Formulations based on linear impulse response models can be found in [START_REF] Normey-Rico | A smithpredictor-based generalised predictive controller for mobile robot pathtracking[END_REF].

• Model Algorithmic Control (MAC) [START_REF] Richalet | Algorithmic control of industrial processes[END_REF][START_REF] Richalet | Model predictive heuristic control: Application to industrial processes[END_REF][START_REF] Ramine | Model algorithmic control (MAC); basic theoretical properties[END_REF]: this approach is also known as model predictive heuristic control. It uses an impulse response model (only for stable processes) and introduces a reference trajectory as the output of a first order system which evolves from the actual output to the setpoint according to a determined time constant.

• Predictive Functional Control (PFC) [START_REF] Richalet | Predictive functional control application to fast and accurate robots[END_REF][START_REF] Richalet | Pratique de la commande prédictive[END_REF]: this strategy is dedicated to fast processes although it can accomodate all kinds of systems and uses a state space model representation of the process. It has two distinctive characteristics: the use of coincidence points by considering only a subset of points in the prediction horizon and the parametrization of the control signal by means of a set of polynomial basis functions.

• Extended Prediction Self Adaptive Control (EPSAC) [START_REF] De Keyser | Extended prediction self-adaptive control[END_REF]: in this method, the process is modelled by a discrete transfer function for prediction. The control signal is set constant over the prediction horizon. The feedforward effect can be included by extending the model with a measurable disturbance term. The control signal can be calculated analytically leading to a very simple control law structure.

• Extended Horizon Adaptive Control (EHAC) [START_REF] Ydstie | Extended horizon adaptive control[END_REF]: the implementation of EHAC is similar to the previous method but it allows a longer time to drive the plant output to its reference value. For prediction, the process is modelled by a transfer function without taking a model of the disturbances into account. It aims at minimizing the discrepancy between the model and the reference at final instant of the prediction window.

• Generalized Predictive Control (GPC) [START_REF] Clarke | Generalized predictive control-Part I. the basic algorithm[END_REF]: the GPC calculates a control sequence in a such way that it results from the minimization of a multi-stage cost function defined over a prediction horizon. An analytic solution is derived in the unconstrained case. Unstable and non-minimum phase plants could be dealt with.

Nonlinear case

There are some processes for which the nonlinearities are so crucial or so important that a linear model is not sufficient. For these systems, a linear MPC will not be very effective, so some solutions have been proposed to cope with this problem. Using a nonlinear model changes the problem from an online convex quadratic problem to a repeated online possibly non-convex nonlinear problem, which is more difficult to solve. As a consequence, the convergence and success of a given optimization problem depend largely on the initial guess provided for the solution. Nonlinear model predictive control (NMPC) is an optimization-based method for the feedback control of nonlinear systems. Its primary applications are stabilization and tracking problems. NMPC is based on the receding horizon principle, where a finite horizon open-loop optimal control problem is solved at each sampling instant and the optimized control trajectory is implemented until a new optimized control trajectory is available at the next sampling instant. Its philosophy is therefore similar to linear MPC. The key characteristics of the NMPC are as follows:

• it allows the use of an explicit nonlinear model for prediction.

• it allows the explicit consideration of state and input constraints.

• specified performance criteria can be minimized online.

• it requires solving online an open-loop optimal control problem.

• system states must be measured or estimated in order to perform the prediction.

Remark 3.2. It should be noted that except the first point, all others are present in the linear MPC law.

This subsection aims at drawing up a picture of some classical NMPC algorithms:

• Explicit schemes [START_REF] Bemporad | The explicit linear quadratic regulator for constrained systems[END_REF][START_REF] Grancharova | Explicit Nonlinear Model Predictive Control[END_REF][START_REF] Olaru | A parameterized polyhedra approach for explicit constrained predictive control[END_REF]: they compute the control sets offline by enumerating all the possible states or their approximations. Then, online control actions are chosen from these sets based on where the state lies.

• Infinite-horizon NMPC [START_REF] Keerthi | Optimal infinite-horizon feedback laws for a general class of constrained discrete-time systems: Stability and moving-horizon approximations[END_REF]: this problem has only theoretical significance, because an infinite horizon NMPC formulation is that at every sampling instant an infinite dimensional optimization problem must be solved as follows:

min u k ,...,u∞ ∞ i=1 ||x k+i -x r k+i || 2 α i + ||∆u k+i-1 || 2 β i (3.5) s.t xk+i+1 = f (x k+i , u k+i ), for i ≥ 0 xk = x k (3.6)
Moreover, infinite-horizon NMPC has the property that the open-loop predictions are identical to the closed-loop response in nominal application. Ideally, we would like to use infinite-horizon NMPC formulations due to stability properties. Unfortunately, the main problem is that infinite-horizon schemes can often not be applied in practice.

• Finite-horizon NMPC with terminal equality constraint [START_REF] Keerthi | Optimal infinite-horizon feedback laws for a general class of constrained discrete-time systems: Stability and moving-horizon approximations[END_REF]: the infinite horizon NMPC can be approximated by a finite horizon formulation with terminal equality constraints as follows:

min u k ,...,u k+Np-1 Np i=1 ||x k+i -x r k+i || 2 α i + ||∆u k+i-1 || 2 β i (3.7) s.t xk+i+1 = f (x k+i , u k+i ), for i = 0, N p -1 xk = x k , xk+Np = x r k+Np (3.8)
The terminal equality constraint, xk+Np = x r k+Np , basically ensures that at the end of the finite horizon, the closed-loop system approaches the same steady state as in the infinite horizon.

• Finite-horizon NMPC with terminal inequality constraint [START_REF] Morari | Model predictive control: past, present and future[END_REF][START_REF] Scokaert | Suboptimal model predictive control (feasibility implies stability)[END_REF]: in this method, a terminal constraint set at the end of the horizon is introduced in the problem formulation. The finite-horizon NMPC goal is to steer the plant to the constraint set X Np in order to guarantee the stability of the closed-loop system. The optimization problem is then given by: min

u k ,...,u k+Np-1 Np i=1 ||x k+i -x r k+i || 2 α i + ||∆u k+i-1 || 2 β i (3.9) s.t xk+i+1 = f (x k+i , u k+i ), for i = 0, N p -1 xk = x k , xk+Np ∈ X k+Np (3.10)
• Finite-horizon NMPC with terminal cost [START_REF] Bitmead | Adaptive Optimal Control: The Thinking Man's GPC[END_REF]: the proposed idea to ensure the closed-loop stability is to add a terminal cost function at the end of the finite-horizon as follows:

min u k ,...,u k+Np-1 Np i=1 ||x k+i -x r k+i || 2 α i + ||∆u k+i-1 || 2 β i + T (x k+Np ) (3.11) s.t xk+i+1 = f (x k+i , u k+i ), for i = 0, N p -1 xk = x k (3.12)
where T (x k+Np ) is the terminal cost function.

• Quasi-infinite horizon NMPC (QIH-NMPC) [START_REF] Chen | A quasi-infinite horizon nonlinear model predictive control scheme with guaranteed stability[END_REF]: this strategy consists in combining both the terminal cost function and terminal inequality constraint into the finite-horizon formulation. min

u k ,...,u k+Np-1 Np i=1 ||x k+i -x r k+i || 2 α i + ||∆u k+i-1 || 2 β i + T (x k+Np ) (3.13) s.t xk+i+1 = f (x k+i , u k+i ), for i = 0, N p -1 xk = x k , xk+Np ∈ X k+Np (3.14)
The terminal state penalty term and the terminal inequality constraint have to be chosen in order to guarantee the asymptotic stability of the closed-loop system [START_REF] Chen | A quasi-infinite horizon nonlinear model predictive control scheme with guaranteed stability[END_REF].

Remark 3.3. Thanks to the principle of optimality, we can recover infinite-horizon NMPC from finite-horizon NMPC with an appropriate choice of the terminal penality T (.) and terminal region X Np [START_REF] Chen | A quasi-infinite horizon nonlinear model predictive control scheme with guaranteed stability[END_REF].

• State/Output-Feedback NMPC [START_REF] Findeisen | State and output feedback nonlinear model predictive control: An overview[END_REF]: one of the key obstacles of NMPC is that it is inherently a state feedback control scheme using the current state and system model for prediction. Thus, for an application of predictive control in general, the full state information is necessary and must be reconstructed from the available output information. In many applications, however, the system state can not be fully measured, i.e. only some outputs are directly available for feedback.

Various researchers have addressed the question of output feedback NMPC using observers for state recovery. To achieve non-local stability results of the observer-based output feedback NMPC controller, two possibilities seem to be attractive:

-Certainty equivalence approach: the state observer is used as the real system state following the "certainty equivalence" principle.

The stability of the closed-loop system is established thanks to the separation of the observer error from the state feedback by time scale separation. Another possibility is to use observers for which the speed of convergence of the observation error can be made sufficiently fast and the absolute achieved observation error can be made sufficiently small. Semi-regional stability results for the closed-loop can be established in this case [START_REF] Findeisen | State and output feedback nonlinear model predictive control: An overview[END_REF].

-Taking into account the observation error in the NMPC controller by using some bounds on the estimation error. This solution is closely related to the design of robust stabilizing NMPC schemes and typically requires observers that deliver an estimate of the observation error.

• Economic NMPC (E-NMPC) [START_REF] Diehl | A lyapunov function for economic optimizing model predictive control[END_REF]: the E-NMPC has essentially the same characteristics of the formulation and implementation of a NMPC, the only difference is that the objective function is a general economically-oriented cost function of the state and manipulated variables instead of a tracking quadratic function. For example, the objective function for an extractive distillation column [START_REF] Santamaria | Economic oriented nmpc for an extractive distillation column using an index hybrid dae model based on fundamental principles[END_REF] is a trade-off between the energy cost, reboiler duty, the production revenue and distillate flow rate.

Predictive Control architecture

There are large-scale systems for which control problems to be solved are too complex by using a unique controller. The reason is that these systems to be controlled are often composed by many interacting subsystems with an increasing complexity, e.g. water systems, process plants, interconnected power systems, manufacturing systems and traffic networks. They necessitate new ideas for dividing the analysis and synthesis of the overall system into independent or almost independent subproblems, for dealing with the incomplete information about the system. To overcome this difficulty, many control structures (i.e. distributed and decentralized) have been developped [START_REF] Magni | Stabilizing decentralized model predictive control of nonlinear systems[END_REF][START_REF] Christofides | Distributed model predictive control: A tutorial review and future research directions[END_REF][START_REF] Scattolini | Architectures for distributed and hierarchical model predictive control-A review[END_REF]. In the sequel, three types of implementation will be presented in order to control the system.

• Centralized control: in a centralized control architecture, the controller manages all the subsystems. Figure 3.3 is a schematic representation of a centralized MPC architecture (with y r the reference, y the output, u the control input and x the state vector). All the MPC laws presented previously in sections 3.4.1-3.4.2 are centralized control strategies.

• Decentralized control: this kind of strategy is based on control architectures where the control input and the controlled variables are gathered into disjoint sets. These sets are then coupled to produce non-overlapping pairs for which local regulators are designed to operate in a completely independent fashion, like the simple example shown in Figure 3.4 (with same notations as in Figure 3. 3). An illustrative application to formation control with collision avoidance for a multi-UAV system can be found in [START_REF] Kuriki | Formation control with collision avoidance for a multi-UAV system using decentralized MPC and consensusbased control[END_REF]. • Distributed control: in distributed control systems, the achievement of a global control task is obtained by the cooperation of many controllers, each one computing a subset of control commands individualy under a possibly limited exchange of information with the other controllers. In fact, the local regulators are designed with MPC, the information transmitted typically consists in the future predicted control computed locally, so that any regulator can predict the interaction effects over the considered prediction horizon. An example of distributed strategy is reported in Figure 3.5 (with same notations as in Figure 3.3). Several approaches have been proposed in the literature. They differ from each other in the assumptions made on the kind of interactions between different components of the overall system (fully or partially connected algorithms), the model of information exchange between subsystems, and the control design technique used for each subsystem. An illustrative application to Hydro-Power Valleys with French main electricity provider EDF can be found in [START_REF] Florez | Explicit coordination for mpc-based distributed control with application to hydropower valleys[END_REF].

Methods for dynamic optimization

The optimization problems in which the model and constraints are linear and the objective is quadratic are a well studied class of problems, and the dynamic optimization problem is solved using well established QP algorithms. In this section, we discuss solution strategies to efficiently solve nonlinear dynamic optimization problems. In general the solution of a nonlinear (in general non-convex) optimization problem can be computationally expensive. Methods to solve the open-loop optimal control problem (3.2-3.3) can be classified into three different approaches (Figure 3.6):

• Hamilton-Jacobi-Bellmann partial differential equations/dynamic programming: this method is based on the direct solution of so-called Hamilton-Jacobi-Bellmann partial differential equations. This approach is valid only for small systems due to the curse of dimensionality, i.e. since the complete solution is considered at once, it is in general computationally intractable. Moreover, the inequality constraints on the state variables, as well as dynamical systems with switching points, lead to discontinuous partial derivatives and cannot be easily included. The application of this methodology is restricted to the case of continuous state systems.

• Euler-Lagrange differential equations/calculus of variations/maximum principle: classical calculus of variations is used in order to obtain an explicit solution of the input as a function of time.

• Direct solution using a finite parametrization of the control signal: The input is parametrized by a finite number of values allowing an approximation of the original open-loop control problem.

The basic idea behind the direct solution using a finite parametrization of the controls is to approximate the original infinite dimensional problem into a finite dimensional nonlinear programming problem. Then, the input control signals are parametrized, e.g. using a piecewise constant function over a sampling time. As classified in [START_REF] Biegler | An overview of simultaneous strategies for dynamic optimization[END_REF][START_REF] Diehl | Nonlinear Model Predictive Control: Towards New Challenging Applications[END_REF] there are two main strategies:

• Sequential approach: in every iteration step of the optimization strategy, the differential equations are solved exactly by a numerical integration, i.e. the solution of the system dynamics is implicitly determined during the calculation of the cost function and only the input vector u k+i-1 appears directly in the optimization problem while the intermediate states xk+i disappear from the problem formulation. This leads to a Direct Single Shooting strategy (DSS) [START_REF] Kraft | Computational Mathematical Programming, chapter On Converting Optimal Control Problems into Nonlinear Programming Problems[END_REF]. This approach is also known as Control Vector Parametrization (CVP) [START_REF] Vassiliadis | Computational solution of dynamic optimization problems with general differential algebric constraints[END_REF];

• Simultaneous approach: the simultaneous approach requires not only an initial control trajectory guess, but also one for the state trajectory. The Ordinary Differential Equations (ODE) are discretized in time and the resulting finite set of nonlinear algebraic equations are treated as nonlinear equality constraints. The intermediate states xk+i are treated as unkown variables together with the parametrized control signal. This leads generally to a Direct Multiple Shooting strategy (DMS) [START_REF] Deuflhard | A modified newton method for the solution of illconditioned systems of nonlinear equations with application to multiple shooting[END_REF] and Collocation methods [START_REF] Biegler | Solution of dynamic optimization problems by successive quadratic programming and orthogonal collocation[END_REF].

For both approaches, the resulting optimization problem is often solved using sequential quadratic programming algorithm (SQP). Remark 3.4. For an online solution of the NMPC problem, only the sequential approach is used in general.

In the case of NMPC, some algorithms were proposed to simplify the optimization procedure:

• Neighboring extremals [START_REF] Würth | Neighboring-extremal updates for nonlinear model-predictive control and dynamic real-time optimization[END_REF]: this method was applied to optimal control problems in [START_REF] Bryson | Applied Optimal Control: Optimization, Estimation and Control[END_REF]. It consists in solving a full solution of the optimal control problem offline using indirect methods and perform a full iteration online to find the approximated optimal solution. Moreover, a neighboring extremal update approach was proposed in [START_REF] Würth | Neighboring-extremal updates for nonlinear model-predictive control and dynamic real-time optimization[END_REF], and more specifically applied to NMPC design. It consists in solving an optimal control problem over a long time horizon. Then, during each sampling time, a fast update, determined from a QP problem, is performed for the control variable.

• Newton-type controller [START_REF] Li | Process control strategies for constrained nonlinear systems[END_REF][START_REF] Dehaan | A new real-time perspective on non-linear model predictive control[END_REF][START_REF] Diehl | A real-time iteration scheme for nonlinear optimization in optimal feedback control[END_REF][START_REF] De Oliveira | An extension of Newton-type algorithms for nonlinear process control[END_REF]: it performs a single full Newton step of NMPC online to allow a quick return of control action. The main advantage is that it rejects the fast disturbance quickly. The Newton's method is based on linearization and it is used in nonlinear programming to solve algebraic equations closely related to the first order optimality conditions of the optimization problem, known as the Karush-Kuhn-Tucker (KKT) conditions. Obviously, the KKT conditions also involve inequalities which means that Newton's method cannot be applied directly.

• NLP sensitivity-based controller [START_REF] Büskens | Online Optimization of Large Scale Systems, chapter Sensitivity Analysis and Real-Time Control of Parametric Optimal Control Problems Using Nonlinear Programming Methods[END_REF][START_REF] Kadam | Sensitivity-based solution updates in closed-loop dynamic optimization[END_REF]: it aims at overcoming the drawbacks posed by the Newton-type controllers. NLP sensitivitybased controllers use a direct method to solve the optimal control problem. The sensitivity analysis for NLP provides information on regularity and curvature conditions at KKT points, evaluates which variables play dominant roles in the optimization, and provides first order estimates for parametric nonlinear programs. In closed-loop, this approach is very effective to handle uncertainty while requiring only a minimum number of full real-time optimizations, reducing the online computational load.

• Advanced step NMPC (asNMPC) [START_REF] Zavala | The advanced-step NMPC controller: Optimality, stability and robustness[END_REF]: a kind of NLP sensitivitybased NMPC which uses the current control action to predict the future plant state in order to solve the future optimization problem in advance, while the current sampling period evolves.

Robust NMPC schemes

Nonlinear MPC is often used in control problems. However, assumptions of nominality reduce its use or performances. This is why other kinds of predictive control methods are developed to take account of uncertainties related to the system model and disturbances acting on it. Hereafter, the main robust NMPC strategies will be briefly presented. There exists a wide range of NMPC formulations that include robustness into the formulation of the optimization problem. We distinguish four main types of approaches:

• Linear parameter varying formulations: these methods are based on the reformulation of nonlinear models as linear parameter-varying (LPV) models which allow the use of linear and bilinear matrix inequality (LMI and BMI) formulations [START_REF] Angeli | Constrained predictive control of nonlinear plants via polytopic linear system embedding[END_REF][START_REF] Casavola | Predictive control of constrained nonlinear systems via LPV linear embedding[END_REF][START_REF] Wan | Efficient scheduled stabilizing output feedback model predictive control for constrained nonlinear systems[END_REF].

• Min-max NMPC: in this formulation, the standard NMPC setup is kept, however, now, the cost function to be optimized considers the worst disturbance sequence occuring. The disturbance/uncertainty in the system is considered as a player working against the input. Thus, the problem is formulated as a classical min-max problem over a finitehorizon. There are two formulations of min-max NMPC:

-Open-loop formulation: this approach guarantees the robust stability of the system and the robust feasibility of the optimization problem, but it may be conservative since the control sequence has to ensure constraints fulfilment for all possible uncertainty scenarii without considering the fact that future measurements of the state contain information about past uncertainty values [START_REF] Limon | Input-to-state stable MPC for constrained discrete-time nonlinear systems with bounded additive uncertainties[END_REF][START_REF] Magni | Assessment and Future Directions of Nonlinear Model Predictive Control, chapter Robustness and Robust Design of MPC for Nonlinear Discrete-Time Systems[END_REF][START_REF] Michalska | Robust receding horizon control of constrained nonlinear systems[END_REF]. This method may have a small feasible set and a poor performance due to the fact that it does not include the effect of feedback provided by the receding horizon framework.

-Closed-loop formulation: the min-max problem represents a differential game where the controller is the minimizing player and the disturbance is the output of the maximizing player [START_REF] Limon | Input to state stability of min-max MPC controllers for nonlinear systems with bounded uncertainties[END_REF][START_REF] Magni | Robust model predictive control for nonlinear discrete time systems[END_REF][START_REF] Magni | Assessment and Future Directions of Nonlinear Model Predictive Control, chapter Robustness and Robust Design of MPC for Nonlinear Discrete-Time Systems[END_REF].

The controller chooses the control input as a function of the current state so as to ensure that the effect of the disturbance on the system output is sufficiently small for any choice made by the maximizing player. This method would guarantee a larger feasible set and a higher level of performances in comparison with the open-loop min-max NMPC. By doing this, the reaction of the controller to the uncertainty is incorporated in the prediction and the conservativeness is mitigated.

The robust min-max strategy and more specifically the closed-loop formulation, will be detailed in Chapter 5.

• H ∞ -NMPC : the standard H ∞ problem is implemented in a receding horizon framework by considering a particular choice of the cost function. More details concerning this approach can be found in [START_REF] Gautam | Robust H ∞ receding horizon control for a class of coordinated control problems involving dynamically decoupled subsystems[END_REF][START_REF] Wang | h ∞ control for networked predictive control systems based on the switched lyapunov function method[END_REF][START_REF] Wang | h ∞ controller design for networked predictive control systems based on the average dwell-time approach[END_REF].

• Stochastic Model Predictive Control (SMPC): all previous strategies consider a deterministic framework by assuming that pertubations and uncertainties belong to compact bounded sets. However, it can be sometimes more realistic to solve problems based on a probabilistic description of the uncertainties and perturbations. This corresponds to a stochastic approach. Thus, formulations based on a probabilistic description of uncertainty can also be characterized as problems involving random model uncertainty with known probability distribution. Furthermore, stochastic NMPC can also be characterized similarly to minmax NMPC as described above, either open-loop formulation [START_REF] Cannon | Nonlinear Model Predictive Control: Towards New Challenging Applications, chapter Successive Linearization NMPC for a Class of Stochastic Nonlinear Systems[END_REF] or closed-loop formulation [START_REF] Goodwin | Nonlinear Model Predictive Control: Towards New Challenging Applications, chapter A Vector Quantization Approach to Scenario Generation for Stochastic NMPC[END_REF].

Stability

One of the key questions in the MPC strategy is certainly the stability of the closed-loop system. After basic theoretical principles of MPC had been clarified, more advanced topics like robust stability under perturbations, performance estimates and efficiency of numerical algorithm were addressed. The stability of MPC for constrained nonlinear systems necessitates the use of Lyapunov stability theory that can be expressed conveniently via so-called comparison functions, which were introduced in nonlinear control theory by Sontag [START_REF] Sontag | Smooth stabilization implies coprime factorization[END_REF].

The stability of the system without disturbances is called nominal stability, while the stability of the system in the presence of disturbances is termed robust stability.

Nominal stability

Robust stability

As stated before, the robust control design of NMPC is done by taking into account uncertainties in an explicit manner in order to optimize the objective function for the worst situation of the uncertainties. There are many approaches to analyze the robust stability, such as Input-to-State (ISS) framework, robust stability margin and invariant set theory coupled with ISS-Lyapunov stability framework. Further details and a more comprehensive treatment of this topic can be found in [START_REF] Magni | Robustness and robust design of MPC for nonlinear discrete-time systems[END_REF][START_REF] Magni | Assessment and Future Directions of Nonlinear Model Predictive Control, chapter Robustness and Robust Design of MPC for Nonlinear Discrete-Time Systems[END_REF][START_REF] Pin | Robust model predictive control of nonlinear systems with bounded and state-dependent uncertainties[END_REF][START_REF] Raimondo | Robust Model Predictive Control Algorithms for Nonlinear Systems: An Input-to-State Stability Approach[END_REF][START_REF] Limon | Input-to-state stable MPC for constrained discrete-time nonlinear systems with bounded additive uncertainties[END_REF][START_REF] Limon | Input to state stability of min-max MPC controllers for nonlinear systems with bounded uncertainties[END_REF][START_REF] Allwright | Advances in Model-Based Predictive Control, chapter on min-max Model-Based Predictive Control[END_REF][START_REF] Morari | Robust model predictive control[END_REF]. Fundamental results on stability for constrained and nonlinear predictive control for state-space models are well summarized and categorized in [START_REF] Limon | Nonlinear Model Predictive Control: Towards New Challenging Applications, chapter Input-to-State Stability: A Unifying Framework for Robust Model Predictive Control[END_REF][START_REF] Raimondo | Min-max model predictive control of nonlinear systems: A unifying overview on stability[END_REF].

Advantages and drawbacks

The RHC presents a series of advantages over the other existing control strategies, summarized as follows:

• Straightforward formulation, based on well understood concepts,

• Open methodology which allows future extensions,

• Compensation for dead times,

• Applicability to a great variety of systems, including nonlinear systems and time-delayed systems,

• Proof of stability for linear and nonlinear systems with input and state constraints, under certain specific conditions,

• Very useful when future references are known a priori, the system can react before the change has effectively been made inducing an anticipative effect, thus minimizing the effects of delay in the process response,

• Multivariable case can easily be dealt with,

• Extension to the treatment of the constraints (take into account actuators limitation),

• Development time relatively shorter than for competing advanced control methods.

• Easier to maintain and to implement in case of changing the model or its parameters. Indeed, model modification does not require complete redesign.

It also has its drawbacks:

• One of these is that the control law requires a longer computation time compared with conventional control strategies, e.g. Lyapunov-based design. When inequality constraints are considered, the amount of computation required is even higher,

• Generally, lack of stability proof except for special cases (as mentioned in the point above in the advantage section).

However, the greatest requirement is the need for an appropriate model of the process to be available. In fact, the design algorithm is based on a prior knowledge of the model but is not dependant on a specific model structure. It is however obvious that benefits obtained will be affected by the gap existing between the real process and the prediction model.

Concluding remarks

The receding horizon control is one of the most popular advanced control strategies due to its simplicity of implementation and interesting properties mentioned above. This chapter presented a general formulation of the model predictive control that will be used in the context of trajectory tracking problem. The brief review focused on the area of linear and nonlinear MPC schemes, addressing robustness issues and nominal and robust stabilities.

The specific task we will focus on in the next chapter is the setpoint tracking NMPC for discrete-time nonlinear systems.

Chapter 4

Nonlinear Model Predictive Control

Introduction

In this chapter, we investigate the strategy that will be further considered and modified to solve a trajectory tracking problem for continuous/discrete-time nonlinear systems. This leads to explore the NMPC formulation discussed in Chapter 3. The chapter will be organized as follows. Section 4.2 is devoted to the NMPC problem formulation and the class of nonlinear systems that will be considered. A variant of NMPC is presented in section 4.3 and an illustrative example (a class of bioprocesses) is presented in section 4.4.

Numerical results are provided in order to demonstrate the effectiveness and limitations of the NMPC strategy. Finally, conclusions are stated in section 4.5.

Problem formulation

This section will introduce the formulation of a basic standard nonlinear finite-horizon optimal control problem starting from a continuous-time model.

Continuous/discrete formulation

As mentionned in Chapter 3, the first step to implement the NMPC strategy is achieved with the use of a prediction model. In this context, it is important to obtain a mathematical representation as properly as possible reproducing the behavior of the system to be controlled. In our case, we will use a continuous/discrete formulation: time-continuous model for prediction and discrete model for control. Consider a system described by an uncertain continuous-time nonlinear model:

ẋ(t) = F (x(t), u(t), θ), x 0 = x y(t) = Hx(t) (4.1)
where

• x ∈ X ⊆ R nx is the state vector with X the compact set of admissible states;

• x is the initial state vector;

• y ∈ Y ⊆ R ny is the measured output with Y the compact set of admissible outputs;

• u ∈ U ⊆ R nu represents the control input with U the compact set of admissible controls;

• θ ∈ R n θ is the vector of uncertain parameters that are assumed to lie in the compact set Θ = [θ -, θ + ] defined as follows:

θ = θ nom + δθ (4.2)
where θ nom is the nominal parameters vector defined as the average value (centroid) of the compact set:

θ nom = θ + + θ - 2 (4.3)
and δθ the parameters uncertainties vector;

• The mapping

F : R nx × R nu × R n θ -→ R nx , of class C 1
with respect to all its arguments, represents the nonlinear process dynamics;

• The measurement matrix is given by H ∈ R ny×nx ;

Remark 4.1. The sets X and U are generally polyhedral convex sets, taking into account e.g. physical constraints acting on the system.

Remark 4.2. The measurement could be nonlinear with respect to the state. It is assumed linear to simplifiy mathematical developments.

Exogenous inputs can act on system (4.1). They are omitted to simplify notation (but are applied to the system).

Most models of real life processes are given as continuous-time models, usually in the form of differential equations (4.1).

In order to convert these models for the design purposes, a continuous/discrete formulation is used.

The discrete-time outputs are obtained at each constant sampling time T s by the integration of the continuous-time state space model (4.1) using for example the Runge-Kutta method with an integration time step T d (lower than the sampling time T s ). The control input u(t) is parametrized using a piecewise-constant approximation over a time interval [t k , t k+1 ] [kT s , (k + 1)T s )]:

u k u(τ ) = u(t k ), τ ∈ [t k , t k+1 [ (4.4)
Let us define the discrete state trajectory g (see Figure 4.1) as the solution, at time t k+1 , of system (4.1):

x k+1 = g(t 0 , t k+1 , x, u t k t 0 , θ) y k = Hx k (4.5)
with initial state x 0 , and u t k t 0 the control sequence from the initial time instant t 0 to the time instant t k . It can be obtained with a Runge-Kutta method for example. Thus, the prediction model could be defined by the following recursive equations:

x k+1 = f (x k , u k , θ), k ≥ 0, x 0 = x y k = Hx k (4.6)
where x k+1 is the state at time t k+1 , k is the time index, x k and y k are the discrete state vector and the sampled output at time t k , respectively.

It can be easily shown that:

f (x k , u k , θ) ≡ g(t k , t k+1 , x k , u k , θ) (4.7)
In the sequel, the model (4.5) will be used in the NMPC strategy to predict the future behavior of the system.

Control objectives

In this study, the NMPC problem is formulated for trajectory tracking purposes. The main objective is to force the output signal y to follow a given reference trajectory y r , while the control input u is constrained to track a reference u r . In addition, saturations on the state vector and control input signal with minimum and maximum thresholds x min , x max , u min and u max , respectively can be included (i.e. X = [x min , x max ] and U = [u min , u max ]). These inequality constraints may result from both physical and operational constraints on the controlled system. The predictive controller predicts the plant future evolution ŷk+Np k+1 over a finite time receding horizon of length N p T s , using a nonlinear dynamic model. At each time instant t k , the optimal control sequence over the prediction horizon is computed by minimizing a cost function expressed as a quadratic criterion based on the future tracking errors, while ensuring that all constraints are respected. The first control in the optimal sequence is applied to the system until the next time step, when the measurement becomes available. The optimization problem is solved again at the next sampling time according to the well-known receding horizon principle [START_REF] Camacho | Model Predictive Control[END_REF][START_REF] Maciejowski | Predictive Control: With Constraints[END_REF] as shown on Figure 4.2.

Remark 4.3. This chapter is dedicated to the description of NMPC formalism, without consideration of robustness issues directly in the cost function. Therefore, even if uncertainties are included for simulation purposes, the system model will be considered as the nominal one with θ = θ nom , for the derivation of the optimization problem.

Derivation of the control law

Following the general idea of the previous section, the cost function that will be minimized in this chapter is expressed as the sum of two quadratic functions based on tracking errors over the receding horizon and defined as (at time t k ): The discrete-time prediction model is chosen as:

Π NMPC (x k , u k+Np-1 k ) ||u k+Np-1 k -u r,k+Np-1 k || 2 V + ||ŷ k+Np k+1 -y r,k+Np k+1 || 2 W (4.8) Future Past k k + 1 k + 2 k + Np -1 ŷk ŷk+1 ŷk+2 y k y k+1 y k+2 ŷk+Np-1 y k+Np-1 y k-1 y k-2 u k u k+1 u k+Np-1 u k-1 u k-2
xk+j = g(t k , t k+j , x k , u k+j-1 k , θ) ŷk+j = H xk+j , ∀j = 1, N p (4.9)
subject to the constraints in a matrix form:

    I nx 0 nu -I nx 0 nu 0 nx I nu 0 nx -I nu     xk+j u k+j-1 ≤     x max -x min u max -u min     , ∀j = 1, N p (4.10)
with x k the state vector at time t k , N p is the length of the prediction horizon. V 0 and W 0 are tuning weighting matrices.

u k+Np-1 k =    u k . . . u k+Np-1    the optimization variable, u r,k+Np-1 k =    u r k . . . u r k+Np-1    the reference control sequence, ŷk+Np k+1 =      Hg(t k , t k+1 , x k , u k , θ nom ) Hg(t k , t k+2 , x k , u k+1 k , θ nom ) . . . Hg(t k , t k+Np , x k , u k+Np-1 k , θ nom )     
As mentionned before, assuming a perfect knowledge of the parameter vector θ (i.e. θ = θ nom , determined from an identification procedure for example), the formulation of the optimization problem is moved into NLP problem over a finite prediction horizon N p T s at each sampling time t k . The optimal control sequence is obtained by minimizing the performance criterion (4.8) with the constraints (4.10) as follows:

u k+Np-1 k = arg min u k+Np-1 k Π NMPC (x k , u k+Np-1 k ) (4.12) s.t. (4.9)-(4.10)
This problem is in general non-convex, it is solved classicaly using algorithms for constrained optimization problems. However, in the case where the constraints are only bounds on the optimization variables, the optimization problem is more tractable.

Remark 4.4. This formulation of the optimization problem belongs in fact to the sequential family of dynamic optimization methods (as mentioned in section 3.5).

Remark 4.5. For the NMPC problem (4.12), at each sampling time (k+1)T s the optimization variable u k+Np k+1 is initialized by the optimal control sequence u k+Np-1 k+1 obtained from the optimization (4.12) at time instant k as follows:

u k+Np k+1 ini = [ u k+Np-1 k+1 , u k+Np-1 ], ∀k ≥ 0
The NMPC problem is implemented as stated in Figure 4.3.

NMPC tuning

The computational burden of a NMPC algorithm is high and there are numbers of factors that affect it. These factors include the prediction horizon N p , the optimization method, the number of states n x and control variables n u , the sampling time step T s , the weighting matrices V and W , the nonlinearity of the system (4.1), the calculation time of the model related to its complexity, and the speed of the computer processor. The main factors are discussed hereafter:

• Prediction horizon: the meaning of N p is rather intuitive. It must be chosen in such a way that it gives a sufficient vision of the system behaviour. It marks the limit of the instants in which it is desirable for the output to follow the reference.

Remark 4.6. Increasing the length of the prediction window N p can lead to a loss of tracking accuracy. In fact, a large time window penalizes prediction when the model is far from the real one.

• Weighting matrices: the weights may be dictated by the economic objectives for the control system. Increasing the weighting matrix on the control action, V , relatively to the weighting matrix W has the effect of reducing the control effort and viceversa.

• Sampling time: in general, a major issue when solving an optimization problem is related to the discretization of the model (when using a discrete/discrete formulation). Indeed, the validity of the discretized model requires choosing an adequate sampling time in agreement with the dynamics of the system. However, when solving a NLP problem, this sampling time must be large enough to permit online implementation. This is a trade-off to consider. In our case, we considered a continuous/discrete formulation. The choice helps reducing the impact of the sampling time value on the discrete model and thus on the controller stability. The choice of the sampling time will mainly condition the controller performance, leading to the same trade-off as mentionned previously (trade-off between computation time and controller performance). it should be mentionned that considering a small sampling time in general induces a large value of N p , for a given prediction horizon.

• Optimization method: strategies for solving the NMPC problem are typically based on direct optimization methods using a finite parametrization of the input to find an open-loop solution to problem (4.12) that is implemented in a receding horizon framework. More details can be found in section 3.5.

A variant of NMPC

Model uncertainties can be taken into account by assuming that the gap between the system and the model at time instant k, denoted ε s/m , is propagated over the prediction horizon, i.e. considering the accumulation of errors over j prediction intervals [START_REF] Tebbani | Nonlinear Predictive Control of fed-batch Cultures of E. coli[END_REF]. Thus, the predicted output of the system ŷk+j|k is related to the predicted output through the model ŷm k+j|k after j prediction intervals through the following relation, where y k denotes the measured output of the system at time instant kT s :

ŷk+j|k = ŷm k+j|k + jε s/m k , j = 1, N p ε s/m k = y k -ŷm k|k (4.13)
The term jε s/m k represents the integration of the modelling error up to time j. As a consequence, it could improve the control law performance with respect to model mismatch in comparison to including no error at all as presented before, i.e. ŷk+j|k = ŷm k+j|k , or even including a constant error signal, i.e. we assume that the difference between the system and the model remains constant over the prediction horizon as in [START_REF] Richalet | Predictive functional control application to fast and accurate robots[END_REF]:

ŷk+j|k = ŷm k+j|k + ε s/m k (4.14)
q out q in X S The mass balances lead to the following system of ODE's that describes the microorganisms growth rate and the substrate consumption, respectively:

               Ẋ = µ(S)X 1 -DX 2 Ṡ = - µ(S) Y X 3 + DS i 4 -DS 5 (4.15)
with 1. The growth of biomass;

2. The dilution effect;

3. The consumption of substrate for the biomass growth;

4. The amount of incoming substrate into the reactor through the feed;

5. The dilution impact.

The dilution rate, D (h -1 ), is defined as

D = q in V = qout V .
V is the volume of the bioreactor, S i the inlet substrate concentration (g.L -1 ), Y the yield production coefficient and µ(S) is the specific growth rate (h -1 ). The specific growth rate µ(S) is given by a Monod kinetic:

µ(S) = µ m S S + K m (4.16)
where µ m and K m are the maximal specific growth rate and semi-saturation constant, respectively. The yield production coefficient, Y , represents the efficiency of conversion of substrate into biomass. Following the state space formalism (4.1), the state vector x, the control input u and the parameter vector θ are given by:

x = X S , u = D, θ =   µ m K m Y   (4.17)
The control objective for the bioreactor given by (4.15) is to regulate the biomass concentration X (output variable) at a predetermined setpoint X r , by manipulating the dilution rate D (control input) which remains closely to a specified reference, D r nom , while satisfying some accuracy requirements. In order to illustrate the performance of the closed-loop system, numerical simulations were carried out by considering the values given in Table 4.1 for the parameters involved in the bioreactor model.

     S r nom = S i - X r Y D r nom = µ m S r nom S r nom + K m (4.18)
In this nominal case, the biomass concentration X will follow its reference value X r , while constraining the dilution rate D to track the reference D r nom . It appears clearly that if the model is uncertain, the equilibrium for the considered values of X r differs from the one given by (4.18). Indeed, the equilibrium corresponding to system under parameter uncertainties is determined as follows:

     S r real = S i - X r Y + δY D r real = (µ m + δµ m ) S r real S r real + K m + δK m (4.19)
where δµ m , δK m and δY are the parameters uncertainties.

In this case, the biomass concentration X should follow its reference value X r , while constraining the dilution rate D to track the reference D r real . In section 4.4.2, the references S r real and D r real will be provided in order to better understand the behavior of the controller in the presence of parameters uncertainties.

Remark 4.7. If there are no uncertainties, the equilibrium given by (4. [START_REF] Biegler | Advances in sensitivitybased nonlinear model predictive control and dynamic real time optimization[END_REF]) is equivalent to the equilibrium given by (4.18).

The maximum admissible dilution rate D max is set to 0.3 h -1 . The inlet substrate concentration S i is assumed to be perfectly known (35 g.L -1 ). The initial conditions of the concentrations in the bioreactor are:

S(0) = S r real (0) g.L -1 X(0) = X r (0) g.L -1 (4.20)
The initial substrate concentration is set to the setpoint S r real in order to cancel the transient effect and to focus the study only on the behavior during the setpoint changes. The time critical code including ODE's was written in C language using Matlab CMEX-functions. The optimization was run on Microsoft PC (Intel(R) Core(TM) i7 -3770, 3.40 GHz, 8GB Ram). The optimization problem was solved using lsqnonlin from Matlab optimization toolbox. Hereafter, the performance of the predictive controller is first studied in the nominal case (section 4.4.1) and then in the case of model parameters mismatch (section 4.4.2).

Influence of the tuning parameters

In this section, we shortly outline the influence of the prediction window N p T s and weighting matrices V and W on the closed-loop trajectory tracking behavior. For this purpose, we consider the NMPC as described above with the assumption that all states can be accessed directly by measurements. Simulations have been conducted considering the nominal case (i.e. without model mismatch).

First, the length of the prediction horizon and the weighting matrices are fixed at the values N p = 6, V = I Np and W = 0.01I Np respectively. The results are compared for three case studies concerning the choice of sampling time T s (case 1: T s = 5 min, case 2: T s = 10 min and case 3: T s = 20 min). The integration time step T d is chosen T d = T s /50. It can be observed that the dilution rate D is decreased at time around 1h (Figure 4.8), this leads to increase the biomass concentration X in order to follow the desired setpoint X r (Figure 4.6). Indeed, the dilution rate is decreased, so that it is lower than the biomass growth rate, leading to the increase of the micro-organisms concentration. When X reaches X r , the controller sets D to D r nom , so that the biomass concentration is maintained at X r . In this case, the substrate concentration decreases since it is consumed by the micro-organisms (Figure 4.7) and then it reaches its reference value. Furthermore, it can be noticed the anticipation behavior to a setpoint change, due to the prediction of the setpoint trajectory future evolution over the moving horizon. In the case of falling edge at time t = 15 h, the dilution rate is set at its maximal value (Figure 4.8) to dilute the culture, X decreases as depicted by Figure 4.6 while S increases as shown in Figure 4.7. Then, when the setpoint is reached, the dilution rate is kept equal to its reference value, as well as the state variables X and S. The results as depicted in Figures 4. 6-4.8 show that the larger sampling time leads to the larger prediction window. The model should be sampled sufficiently fast in order to guarantee its accuracy as much as possible but the sampling time T s should be large enough to take into account the computation burden due to potential state estimator and/or online determination of the optimal trajectory. On the other hand, it can be noticed that the choice of T s is linked to the size of the prediction window for a given value of the prediction horizon N p . The sampling time T s is hereafter chosen equal to 5 min, which is reasonable according to the dynamics of the system. Secondly, the influence of the prediction horizon N p is studied for the same conditions cited previously with T s = 5 min, V = I Np and W = 0.01I Np as shown in Figures 4.9-4.11. The results are compared for three case studies (case 1: N p = 3, case 2: N p = 6 and case 3: N p = 12). The integration time step will be given by T d = T s /50 and thus hereafter is set to 0.1 min. It appears clearly that increasing the prediction horizon N p leads to a larger anticipation time to a setpoint change as shown in Figure 4.9 (black box). The prediction horizon N p must be chosen to satisfy a compromise between the computation time and a sufficient vision of the system behavior in the future. In fact, the dimension of the optimization problem increases with N p leading to a greater computation time as shown in Table 4.2 (it provides computation times of iterations). Even if this time remains negligible compared to T s in our case study, it can became critical for more complex systems with large number of parameters. For this example, the prediction horizon is set to N p = 6, which ensures a good trade-off between the computation time and a sufficient vision of the system behavior in the future.

Time (h)

Finally, the length of the prediction horizon, the sampling time and the integration time step are fixed at the values N p = 6, T s =5 min and T d =0.1 min respectively. The results are compared for two case studies concerning the choice of weighting matrices (case 1: V = W = I Np and case 2: V = I Np W = 0.01I Np ). The matrix W is chosen less than one in order to emphasize the control tracking error in the cost function, and to induce a slight correction with the term related to the output tracking error (to smoothen the control input evolution). Thus, the choice of the weighting matrices directly conditions the state and control behaviors. They can be for example chosen based on the variables magnitude. The obtained results are depicted in Figures 4.12-4.14. The chosen values of the weights V and W in case 2 lead to an anticipation time of 30 min approximatively as shown in Figure 4.12 (green box) (this is also the case for case 1 even if their effect is here less visible). This result is in agreement with the size of the desired prediction window: N p T s =30 min. In the case 2, the control evolution is smoother and smaller, leading to a slightly longer response time. In order to ensure a quicker response time, the weighting matrices V and W will be chosen equal to I Np hereafter.

Influence of model parameters mismatch

In this section, the NMPC tuning parameters are chosen as follows: N p = 6, T s =5 min, T d =0.1 min and V = W = I Np . Simulations were performed considering that the prediction model is erroneous. The maximum admissible dilution rate D max is set to 0.5h -1 (it was set before at 0.3h -1 in order to emphasize the effects of the studied NMPC parameters and make them more visible). The parameters values of the system are chosen on the parameter subspace border referred to

θ real = [µ + m K - m Y + ].
The parameters θ are assumed to be ±20% uncertain maximum. Figures 4.15-4.17 compare the results obtained when running simulation with and without addition of the error signal jε s/m (see (4.13)) during prediction in the NMPC algorithm, as detailed in section 4.3. The addition of signal ε s/m (see (4.14)) is also considered. The controller that includes the jε s/m signal is denoted NMPC-jε s/m whereas the one that includes ε s/m is denoted NMPC-ε s/m . The dilution rate moves away from its reference value D r nom (calculated with the nominal model (4.18)) and converges to D r real (calculated with (4.19)) as shown in Figure 4.17. It can be seen that NMPC-jε s/m has better performances than the classical NMPC under parameters uncertainties. It can be noticed that without the (j)ε s/m term, the output is not able to track the specified setpoint in the presence of parameters uncertainties, due to the fact that the gap between the system and the model is not considered during the prediction step inside the minimization procedure. Taking the signal jε s/m into account during the optimization strategy helps to improve significantly the tracking accuracy in comparison with when considering a constant signal ε s/m (Figure 4.15). This is due to the fact that the control input of NMPC-jε s/m converges exactly to D r real which is not the case for NMPC-ε s/m . Thus, the predictive controller regulates the biomass to its reference value with good performances in the nominal case. It however lacks accuracy in the case of model mismatch. This example will be considered in the next chapters to evaluate the performance of the proposed predictive algorithms.

Concluding remarks

In this chapter, we presented the NMPC problem formulation for setpoint tracking purpose. This strategy induces solving, online, an optimization problem which is expressed as a nonlinear programming problem. The output signal is not able to track accurately the specified setpoint in the presence of parameters uncertainties, due to the fact that the mismatch between the system and the model is not considered during the prediction step inside the minimization procedure. Consequently, there is a need for improvement in the case of parameters uncertainties. One improvement could be to include, during the prediction step, the difference between the system and the model outputs in order to increase the control law robustness against model mismatch. Several simulations were performed in order to compare the predictive strategy performances for an illustrative biological application. They show the efficiency of the proposed control law, but also highlighted its limits. The assumption on including the prediction error is quite restrictive due to the fact that we do not exactly know, in a realistic application, how the error will evolve over the prediction horizon. If the assumption is false then the performances will be degraded. The drawbacks of this method are a motivation to investigate a new approach in order to satisfy the objectives of robustness. An optimization problem taking into account the model pa-rameters uncertainties must be formulated in order to ensure the robustness against model mismatch. The next chapter will be then devoted to the robust model predictive control.

Chapter 5

Robust Nonlinear Model Predictive Control

Introduction

In the previous chapter, a predictive controller was designed to ensure a stable real time operation of the plant, close to a certain state or desired profile. However, the performances of the NMPC law usually decrease when the true plant evolution deviates significantly from the one predicted by the model. Robust variants of NMPC [START_REF] Kerrigan | Feedback min-max model predictive control using a single linear program: Robust stability and the explicit solution[END_REF][START_REF] Limon | Robust stability of minmax MPC controllers for nonlinear systems with bounded uncertainties[END_REF] as detailed in section 3.6 are able to take into account set bounded disturbance and/or constraints. In addition, when considering uncertain systems, or when only a limited amount of data is available, it is necessary to use a robust formulation of the controller where the effect of the uncertainties can be taken into account in the design procedure. The robust NMPC (RNMPC) law can be formulated as a nonlinear min-max optimization problem (where the objective function is minimized for the worst possible uncertainty realization). However, this approach tends to become too complex to be solved numerically online. Consequently, the total calculation time is an important factor that must be reduced as much as possible. We propose mainly two solutions to address this problem. The first one consists in reducing the number of uncertain parameters through a sensitivity analysis. Only the main influencial parameters on the model are considered in the min-max problem whereas the other parameters are kept constant and equal to their nominal values. This approach will be referred to as reduced RNMPC (rRNMPC). The second approach aims at transforming the min-max problem into a more tractable optimization problem using a model linearization technique. The objective is then to propose a new formulation that is computationally more tractable in calculating the opti-mal control compared to a classical min-max robust approach, which makes it suitable for online implementation. This approach will be referred to as Linearized RMPC (LRMPC).

Stability properties of the robust model predictive control strategy taking into account bounded uncertainties has been analyzed in [START_REF] Magni | Robustness and robust design of MPC for nonlinear discrete-time systems[END_REF][START_REF] Mayne | Control of constrained dynamic systems[END_REF][START_REF] Raimondo | Min-max model predictive control of nonlinear systems: A unifying overview on stability[END_REF][START_REF] Raimondo | Robust Model Predictive Control Algorithms for Nonlinear Systems: An Input-to-State Stability Approach[END_REF][START_REF] Mayne | Constrained model predictive control: Stability and optimality[END_REF][START_REF] Biegler | Advances in sensitivitybased nonlinear model predictive control and dynamic real time optimization[END_REF][START_REF] Scokaert | Suboptimal model predictive control (feasibility implies stability)[END_REF]. Based on work developped by [START_REF] Limon | Robust stability of minmax MPC controllers for nonlinear systems with bounded uncertainties[END_REF] and [START_REF] Raimondo | Min-max model predictive control of nonlinear systems: A unifying overview on stability[END_REF], and taking the objective function as the Lyapunov candidate function, the robust stability of the closed-loop system while applying the LRMPC law is established under some assumptions.

The LRMPC controller represents our main contribution. Thus, the stability and performances analysis will be focused on this approach.

The chapter is organized as follows. The next section is devoted to the min-max problem RNMPC formulation. A variant of RNMPC, based on a reduction of the number of the parameters θ that will be optimized, is presented in section 5.3. Section 5.4 is dedicated to the use of a model linearization technique in order to approach the min-max problem. Stability analysis of the closed-loop system is also discussed. The unconstrained formulation (LRMPC) is detailed in section 5.5. In section 5.6, a robust predictive controller (CLRMPC), similar to the above one with consideration of inequality constraints on the control input, is presented. Numerical results are provided in section 5.7 and conclusions are stated in section 5.8.

Min-max strategy

Since the predictive controller is model-based, it is very sensitive to model uncertainties, and more specifically to the model parameters values. In this context, we will assume that the parameter vector θ is uncertain and belongs to a known set Θ as stated in section 4.2.1. Thus, robustification in the presence of model uncertainties naturally leads to the formulation of a nonlinear min-max optimization problem [START_REF] Allwright | Advances in Model-Based Predictive Control, chapter on min-max Model-Based Predictive Control[END_REF][START_REF] Magni | Assessment and Future Directions of Nonlinear Model Predictive Control, chapter Robustness and Robust Design of MPC for Nonlinear Discrete-Time Systems[END_REF][START_REF] Michalska | Robust receding horizon control of constrained nonlinear systems[END_REF]. The control sequence that minimizes a worst case cost function is derived from the following optimization problem (at time index k and using the same notation as in section 4.2):

u k+Np-1 k = arg min u k+Np-1 k max δθ Π RNMPC (x k , u k+Np-1 k , δθ) (5.1) subject to                xk+j = g(t k , t k+j , x k , u k+j-1 k , θ = θ nom + δθ), j = 1, N p θ ∈ Θ     I nx 0 nu -I nx 0 nu 0 nx I nu 0 nx -I nu     xk+j u k+j-1 ≤     x max -x min u max -u min     , ∀j = 1, N p (5.2)
where θ and θ nom are given in (4.2) and (4.3).

The cost function associated to the future evolution of the system depends on the future control actions and the uncertainties as follows:

Π RNMPC (x k , u k+Np-1 k , δθ) ||u k+Np-1 k -u r,k+Np-1 k || 2 V + ||ŷ k+Np k+1 -y r,k+Np k+1 || 2 W (5.
3)

The predicted outputs ŷk+Np k+1 having the same formulation as in (4.11):

ŷk+Np k+1 =      Hg(t k , t k+1 , x k , u k , θ) Hg(t k , t k+2 , x k , u k+1 k , θ) . . . Hg(t k , t k+Np , x k , u k+Np-1 k , θ)      (5.4)
Relations (5.2) and (5.3) are thus similar to the classical NMPC approach, but by considering a value for θ given by the optimization algorithm instead of θ nom . The optimal control sequence u k+Np-1 k is determined to minimize the tracking error by considering all trajectories over all possible data scenarii [START_REF] Du | Minimax and Applications[END_REF][START_REF] Kasperski | Discrete Optimization with Interval Data: Minmax Regret and Fuzzy Approach[END_REF]. The RNMPC problem is implemented as depicted in Figure 5.1.

Remark 5.1. For the min-max problem (5.1), at sampling time (k + 1)T s the optimization variables u k+Np k+1 and δθ are initialized by the optimal control sequence u k+Np-1 k+1 and the optimal parameter vector δ θ obtained from the optimization (5.1) at time instant k as follows:

u k+Np k+1 ini = [ u k+Np-1 k+1 , u k+Np-1 ], ∀k ≥ 0 δθ ini = δ θ = arg max δθ Π RNMPC (x k , u k+Np-1 k , δθ)
Remark 5.2. It appears clearly that the computation load grows with the size of the parameters vector, the number of control inputs and the prediction horizon, while the control strategy has to be implemented online.

The challenge is to reduce the computation burden while maintaining good performances in term of accuracy.

Reduced robust predictive controller

Since the initial min-max optimization problem (5.1)-(5.3) is time consuming, it will be simplified by reducing the number of the parameters θ that will be optimized from a sensitivity analysis of the model with respect to its parameters.

Sensitivity analysis

The sensitivity functions S x i θ j represent the sensitivity of each state x i to (small) variations in each model parameter θ j . It is defined as the partial derivatives of the state variable vector x i with respect to θ j :

S x i θ j ∂x i ∂θ j , i = 1, n x and j = 1, n θ (5.5)
Different approaches are possible to determine the sensitivity functions. The most precise method involves analytical derivation [START_REF] Dochain | Automatic control of bioprocesses[END_REF]. In this case, the dynamics of sensitivities are calculated as follows for system (4.1):

Ṡx i θ j = d dt ∂x i ∂θ j = ∂ ∂θ j dx i dt = ∂F i ∂θ j + nx k=1 ∂F i ∂x k ∂x k ∂θ j (5.6)
with as an initial condition: ∂x i ∂θ j = 0. From the analysis of the sensitivity functions temporal evolution, and according to their magnitude order, one can select the parameters which are significantly the most influential on the model [START_REF] Dochain | Automatic control of bioprocesses[END_REF]. These parameters will be denoted κ.

Problem reformulation

In the sequel, only the most influential parameters, κ ∈ Θ κ ⊆ Θ, are considered in the min-max optimization, instead of the full model parameters (with θ [κ, ζ]). The influential parameters κ are defined as follows:

   κ = κ nom + δκ
(5.7)

κ nom = κ + + κ - 2 (5.8)
where κ nom is the nominal influential parameters and δκ the parameters uncertainties vector. The other parameters, ζ, are set to their nominal values with

ζ nom = (ζ + + ζ -)/2.
Thanks to the sensitivity analysis, the min-max optimization (5.1)-(5.3) is approached by the following optimization problem (at time index k):

u k+Np-1 k = arg min u k+Np-1 k max δκ Π rRNMPC (x k , u k+Np-1 k , δκ) (5.9) subject to                    xk+j = g(t k , t k+j , x k , u k+j-1 k , κ nom + δκ ζ nom ), j = 1, N p κ ∈ Θ κ     I nx 0 nu -I nx 0 nu 0 nx I nu 0 nx -I nu     xk+j u k+j-1 ≤     x max -x min u max -u min     , ∀j = 1, N p (5.10)
where κ nom and δκ are given above.

The new cost function is defined as follows:

Π rRNMPC (u k+Np-1 k , δκ) ||u k+Np-1 k -u r,k+Np-1 k || 2 V + ||ŷ k+Np k+1 -y r,k+Np k+1 || 2 W (5.11)
The predicted outputs ŷk+Np k+1 are given by:

ŷk+Np k+1 =      Hg(t k , t k+1 , x k , u k , [κ , ζ nom ] ) Hg(t k , t k+2 , x k , u k+1 k , [κ , ζ nom ] ) . . . Hg(t k , t k+Np , x k , u k+Np-1 k , [κ , ζ nom ] )      (5.12)
This approach presents however some drawbacks:

• A sensitivity analysis of the model with respect to it parameters is required, leading to a more complex and time demanding developments.

• For some systems, all parameters can have high influence on the model, and thus no parameter reduction is possible (e.g. for systems with a small number of parameters).

Consequently, this approach cannot be used in all cases (it is then equivalent to the RNMPC). Hereafter, this approach is referred to as reduced RNMPC (rRNMPC).

Linearized robust predictive controller

Since the min-max optimization problem (5.1)-(5.4) even when reducing the number of uncertain parameters as in (5.9)-(5.12), is time consuming, it will be approached further by a more tractable optimization problem. This step focuses on reducing the computational burden of the initial problem. The key idea is to approach the predicted outputs based on the nonlinear model through linearization technique.

As a direct result of the linearization, the non-convex problem will be approached by a convex one. So, more flexible tools will be accessible to better handle the optimization problem.

Main principle

In the following proposed approach, the outputs in the moving time frame are predicted by Taylor series expansion. A similar dual problem for robust state estimation, consisting in the design of receding-horizon observer, was presented in [START_REF] Goffaux | Design of a robust nonlinear receding-horizon observer-Application to a biological system[END_REF]. In this section, we propose to apply the same approach to the case of NMPC law design.

Based on a first order Taylor series expansion of (4.7), the prediction model g for time t k+1 , starting from state x k , is linearized around the reference trajectory given by the reference control u r k and for the nominal parameters θ nom :

xk+1 ≈ xnom,k+1 + ∇ θ g(t k+1 )δθ + ∇ u g(t k+1 )(u k -u r k ) f p (x k , u k , θ nom + δθ) (5.13) with                                δθ = θ -θ nom (5.14) xnom,k+1 = g(t k , t k+1 , x k , u r k , θ nom ) (5.15) ∇ θ g(t k+1 ) = ∂g(t k , t k+1 , x k , u k , θ) ∂θ u k = u r k θ = θnom (5.16) ∇ u g(t k+1 ) = ∂g(t k , t k+1 , x k , u k , θ) ∂u k u k = u r k θ = θnom
(5.17)

Generalizing, the predicted state for time t k+j , starting from state at t k , is linearized around the reference trajectory given by the reference control sequence u r,k+j-1 k and for θ nom . Using the same approach as in (5.13) for j = 1, N p , it comes:

xk+j ≈ xnom,k+j + ∇ θ g(t k+j )δθ + ∇ u g(t k+j ) u k+j-1 k -u r,k+j-1 k (5.18)
In order to simplify the calculation of the gradients ∇ θ g and ∇ u g, finite differences are used to approximate numerically the derivatives ∇ θ g(t k+j ) and ∇ u g(t k+i ). From (4.11) and (5.18), the predicted ouputs over the moving horizon are expressed as follows:

ŷk+Np k+1 = G k+Np nom,k+1 + G k+Np θ,k+1 δθ + G k+Np-1 u,k (u k+Np-1 k -u r,k+Np-1 k ) (5.25)
where

G k+Np nom,k+1 =        H xnom,k+1 . . . H xnom,k+j . . . H xnom,k+Np       
, is the column vector containing the predicted output for the nominal case.

G k+Np θ,k+1 =        H∇ θ g(t k+1 ) . . . H∇ θ g(t k+j ) . . . H∇ θ g(t k+Np )       
, regroups the Jacobian matrices related to the parameters.

G k+Np-1 u,k =       
H∇ u g(t k+1 ) . . . H∇ u g(t k+j ) . . .

H∇ u g(t k+N p )       
, regroups the Jacobian matrices related to the control sequence. Assuming that the uncertain parameters are uncorrelated and recalling that

θ -≤ θ ≤ θ + (5.26)
and

θ nom = θ + + θ - 2 (5.27) Thus, θ --θ + 2 ≤ θ -θ nom ≤ θ + -θ - 2 (5.28)
Then, the bounded parametric error δθ can be expressed by: δθ = γδθ max (5.29) with δθ max = (θ +θ -)/2 (5.30) and ||γ|| ≤ 1 (5.31)

The initial objective function Π RNMPC (5.3) is substituted by a cost function using the equation (5.25). The result is given by the following expression (with the same notations as in (4.11)):

Π RNMPC (x k , u k+Np-1 k , δθ) ≈ ||u k+Np-1 k -u r,k+Np-1 k || 2 V + ||G k+Np nom,k+1 -y r,k+Np k+1 + G k+Np θ,k+1 δθ + G k+Np-1 u,k (u k+Np-1 k -u r,k+Np-1 k )|| 2 W Π(x k , u k+Np-1 k , δθ) (5.32)
The new optimization problem is given by:

u k+Np-1 k = arg min u k+Np-1 k max δθ Π(x k , u k+Np-1 k , δθ) (5.33) subject to θ ∈ Θ, x ∈ X, u ∈ U δθ = γδθ max (5.34)

Stability analysis

In this section, the robust stability of the closed-loop system (4.6) with (5.33)-(5.34) is analysed by exploiting the results obtained in [START_REF] Biegler | Advances in sensitivitybased nonlinear model predictive control and dynamic real time optimization[END_REF][START_REF] Raimondo | Min-max model predictive control of nonlinear systems: A unifying overview on stability[END_REF][START_REF] Limon | Robust stability of minmax MPC controllers for nonlinear systems with bounded uncertainties[END_REF].

More specifically, we will use the following theorem:

Theorem 5.1. Consider a discrete-time nonlinear system given by:

x k+1 = l(x k , w k ), k ≥ 0, x 0 = x (5.35)
where x k ∈ X is the state of the system, w k ∈ W is the disturbance vector (s.t. W is a compact set that contains the origin).

If system (5.35) admits a robust Lyapunov function, then it is robustly stable.

Proof. see [START_REF] Limon | Robust stability of minmax MPC controllers for nonlinear systems with bounded uncertainties[END_REF].

We adapt the results in the above cited references to address the stability analysis of the proposed control strategy. First, a bound on the prediction error will be determined. Secondly, we give an upper and lower bounds on the optimal cost. Finally, we establish the robust stability of the closed-loop system.

In the sequel, for the stability analysis of the closed-loop system (4.6) with (5.1)-(5.32), we need to consider the following assumptions: Assumption 5.1. The state of the plant x k is measured at each sampling time.

Assumption 5.2. The state x and the control u of the plant must fulfill the following constraints: x k ∈ X and u k ∈ U. where X and U are compact sets, both of them contain the origin.

Bound on prediction error

In this subsection, an upper bound on the prediction error provided by the linearization step is derived. Consider the real system for time t k+1 , starting from state x k at time t k :

x k+1 = f (x k , u k , θ nom + δθ s ) (5.36)
where θ nom is the nominal parameters vector and δθ s the real parameter values mismatch.

Assumption 5.3. The function f is of class C 2 with respect to all its arguments.

From Assumption 5.3, and using Taylor developments (around θ nom and u r k ), the system dynamics can be rewritten as follows:

x k+1 = f (x k , u r k , θ nom ) + ∇ u f (u r k , θ nom ).(u k -u r k ) + ∇ θ f (u r k , θ nom ).δθ s + ϑ(|u k -u r k | 2 ) + ϑ(|δθ s | 2 ) (5.37)
where ϑ(|.| 2 ) is the remainder term of the Taylor series expansion limited to the first order.

Assumption 5.4. The error with respect to the first Taylor expansion, w p , defined as:

w p ϑ(|u k -u r k | 2 ) + ϑ(|δθ s | 2 ) (5.38)
is assumed to be bounded as follows:

∃η 1 ∈ R + , such that |w p | ≤ η 1 (5.39)
Consider the prediction model (Taylor series expansion) for time t k+1 , starting from state x k at time t k :

xk+1|k = f p (x k , u k , θ nom + δθ) f (x k , u r k , θ nom ) + ∇ u f (u r k , θ nom ).(u k -u r k ) + ∇ θ f (u r k , θ nom ).δθ (
Let us define η ∈ R + by:

η = max(η 1 , η 2 ) (5.43)
From (5.37) and (5.40), the prediction error at time index k + 1 for u k = u k is given by:

x k+1 -xk+1|k = ∇ θ f (u r k , θ nom ).(δθ s -δθ) + w p (5.44)
Thanks to the triangle inequality and using the Assumptions 5.3 and 5.5, we obtain an upper bound on the prediction error as follows:

|x k+1 -xk+1|k | ≤ |∇ θ f (u r k , θ nom )|.|(δθ s -δθ)| + |w p | ≤ |∇ θ f (u r k , θ nom )|. (|δθ s | + |δθ|) + |w p | (5.45)
From Assumptions 5.2, 5.3 and 5.5, we have:

∃α ∈ R + , ∀x, ∀u, ∀θ, |∇ θ f | ≤ α (5.46)
where α is an upper bound of |∇ θ f |. Thanks to Assumptions 5.4-5.5, (5.43) and (5.46), it comes

|x k+1 -xk+1|k | ≤ 2α.η 2 + η 1 ≤ (2α + 1)η Λ(η) (5.47)
where Λ is a K ∞ -function.

Upper and lower bounds on the optimal cost

In this subsection, we would like to find bounds on the optimal cost in order to satisfy condition [START_REF] Ashoori | Optimal control of a nonlinear fed-batch fermentation process using model predictive approach[END_REF] in Definition 14. In the sequel, for manipulation purposes, the optimal cost function (5.32) is rewritten as follows:

Π(x k , u k+Np-1 k , δ θ) k+Np-1 t=k ψ(x t|k , u t ) + T f (x k+Np|k ) (5.48) with          xt|k = f p (x t-1|k , u t-1 , θ nom + δ θ), t = k + 1, k + N p xk|k = x k ψ(x t , u t ) = ũ t vũ t + ỹ t wỹ t , t = k, k + N p -1 T f (x k+Np|k ) = ỹ k+Np wỹ k+Np (5.49) and ũt = u t -u r t ỹt = H xt -y r t (5.50)
The stage cost ψ(x, u) is definite positive, while the terminal cost is denoted by T f (x) : R nx -→ R ≥0 . Without any lack of generality, the weighting matrices V and W are chosen in diagonal form V = vI and W = wI to simplify mathematical developments hereafter.

Remark 5.4. Formulation in (5.48) and (5.32) are the same, T f being clearly split in (5.48).

Remark 5.5. The terminal stage T f is a K ∞ -function.

Assumption 5.6. Let assume the existence of a terminal set Φ, an admissible robust positively invariant set for the system (5.36) which is controlled by the control law u k = π(x k ) ∈ U s.t. the origin is in its interior. Let assume that T f is an associated robust Lyapunov function s.t. for all x k ∈ Φ and for all δθ satisfying (5.42), we have that:

α t (|x k |) ≤ T f (x k ) ≤ β t (|x k |) + ϕ(η) T f (f (x k , u k , θ nom + δθ)) -T f (x k ) ≤ -ψ(x k , u k ) + χ(η) (5.51)
where α t , β t , and χ are K ∞ -functions and ϕ is a K-function.

Lemma 1. Let us consider the system (5.36) and suppose that the uncertainty on θ is modelled by |δθ| ≤ γ(|x|) + η (γ is a K-function). Let Φ and T f (x) satisfy Assumption 5.6, then ∀x ∈ Φ we have that

Π(x k , u, δ θ) ≤ T f (x k ) + N p χ(η) (5.52)
Proof. see Lemma 3 of section 5 in [START_REF] Limon | Robust stability of minmax MPC controllers for nonlinear systems with bounded uncertainties[END_REF].

From (5.52) and Assumption 5.6, we get

Π(x k , u, δ θ) ≤ β t (|x k |) + ϕ(η) + N p χ(η) (5.53)
Assumption 5.7. The stage cost (non-negative) is such that

ψ(x, u) ≥ α ψ (|x|) (5.54)
where α ψ is a K ∞ -function.

From (5.48) and Assumption 5.7, and since ψ and T f are positive functions, then

Π(x k , u, δ θ) ≥ ψ(x k|k , u k ) ≥ α ψ (|x k|k |) = α ψ (|x k |) (5.55)
Thanks to (5.53) and (5.55), it comes:

α ψ (|x k |) ≤ Π(x k , u, δ θ) ≤ β t (|x k |) + ϕ(η) + N p χ(η) (5.56)
Thus, the optimal cost is bounded as given by (5.56).

Robust stability

Theorem 5.2. Consider system (5.35) and suppose that uncertainties are modelled by

|w k | ≤ γ(|x k |) + η.
Then, the uncertain system controlled by the controller u k = π(x k ) is robust stable for any initial x 0 ∈ X Np (Φ). X Np (Φ) is the set of admissible states at time k + N p . Furthermore, the optimal cost is a robust Lyapunov function.

Proof. see [START_REF] Limon | Robust stability of minmax MPC controllers for nonlinear systems with bounded uncertainties[END_REF].

Remark 5.6. The i-step robust stabilizable set X i (Φ) is the set of admissible states which can be steered to the target set Φ in i steps or less by a sequence of admissible control law π(x k ) for all possible realizations of the uncertainty. This set satisfies that X i-1 (Φ) ⊆ X i (Φ) for i ≥ 0 with X 0 (Φ) = Φ. Thanks to the invariance of the terminal set, the feasible region of the controller X Np (Φ) is a robust invariant set for the closed loop system [START_REF] Mayne | Control of constrained dynamic systems[END_REF]. Now, we consider Π(x k , u, δ θ) as our candidate robust Lyapunov function. Then, the optimal cost function at time index k + 1 is defined as follows:

Π(x k+1 , ȗk+Np k+1 , δ θ) k+Np t=k+1 ψ(x t|k+1 , ȗt ) + T f (x k+Np+1|k+1 ) (5.57) with xt|k+1 = f p (x t-1|k+1 , ȗt-1 , θ nom + δ θ), t = k + 2, k + N p + 1 xk+1|k+1 = x k+1
(5.58) where xt|k+1 denotes the state obtained applying the input sequence ȗt-1

k+1 to the prediction model with the initial condition x k+1 . ȗk+Np k+1 denote an admissible solution of the optimization problem at time index k + 1. In the proposed algorithm, it is based on the optimal solution at time index k:

ȗk+Np k+1 = [ u k+Np-1 k+1 , u k+Np-1 ] (5.59)
Notations used in this subsection are detailed in Figure 5.2.

Assumption 5.8. The parameters uncertainties are assumed to be constant throughout the prediction horizon.

Assumption 5.9. The function f p is Lipschitz with respect to x with Lipschitz constant L f x .

Proposition 5.1. Let us define the following residual at time index l:

x (l) xl|k+1 -xl|k , l = k + 1, k + N p -1 (5.60)
Then, with Assumption 5.9, Proof. Using the result obtained in (5.47) and Assumption 5.9, we get that (reminding that ȗl = u l for l ∈

| x (l)| ≤ L l-k-1 f x Λ(η) (
[k + 1, k + N p -1] and xk+1|k+1 = x k+1 ) | x (k + 1)| = |x k+1|k+1 -xk+1|k | ≤ Λ(η) | x (k + 2)| = |f p (x k+1|k+1 , ȗk+1 , θ nom + δ θ) -f p (x k+1|k , u k+1 , θ nom + δ θ)| ≤ L f x |x k+1|k+1 -xk+1|k | ≤ L f x Λ(η) . . . | x (l)| ≤ L l-k-1 f x Λ(η)
(5.62) Now, it is easy to check the result by recurrence. We assume that the proposition holds for l, let us prove that it is also the case for l + 1 (from (5.62)):

| x (l + 1)| = |f p (x l|k+1 , ȗl , θ nom + δ θ) -f p (x l|k , u l , θ nom + δ θ)| ≤ L f x | x (l)| ≤ L l-k f x Λ(η) (5.63)
Which completes the proof by recurrence.

Let define the difference ∆Π as

∆Π Π(x k+1 , ȗ, δ θ) -Π(x k , u, δ θ) = k+Np-1 t=k+1 (ψ(x t|k+1 , ȗt ) -ψ(x t|k , u t )) + ψ(x k+Np|k+1 , ȗk+Np ) -ψ(x k , u k ) + T f (x k+Np+1|k+1 ) -T f (x k+Np|k ) (5.64)
Assumption 5.10. The stage cost ψ is Lipschitz with respect to x with Lipschitz constant L ψx .

Using the Assumptions 5.9-5.10, and from (5.62) we have that

|ψ(x t|k+1 , ȗt ) -ψ(x t|k , u t )| ≤ L ψx |x t|k+1 -xt|k | ≤ L ψx L t-k-1 f x Λ(η) (5.65)
and

| k+Np-1 t=k+1 (ψ(x t|k+1 , ȗt ) -ψ(x t|k , u t ))| ≤ L ψx Np-2 j=0 L j f x Λ(η) (5.66)
Assumption 5.11. The terminal function T f is Lipschitz with respect to x with Lipschitz constant L tx .

Under Assumption 5.11, and from Proposition 5.1 we have that

|T f (x k+Np|k+1 ) -T f (x k+Np|k )| ≤ L tx | x (k + N p )| ≤ L tx L Np-1 f x Λ(η) (5.67)
The last term of (5.64) is rewritten as follows:

T f (x k+Np+1|k+1 ) -T f (x k+Np|k ) =T f (x k+Np+1|k+1 ) -T f (x k+Np|k+1 ) + T f (x k+Np|k+1 ) -T f (x k+Np|k ) (5.68)
If the Assumption 5.6 holds, it comes

T f (x k+Np+1|k+1 ) -T f (x k+Np|k+1 ) ≤ -ψ(x k+Np|k+1 , ȗk+Np ) + χ(η) (5.69)
Thanks to (5.67) and (5.69) and considering the fact that for any scalar ν ∈ R, ν ≤ |ν|, we have that

T f (x k+Np+1|k+1 ) -T f (x k+Np|k ) ≤ -ψ(x k+Np|k+1 , ȗk+Np ) + χ(η) + L tx L Np-1 f x Λ(η) (5.70)
By substituting the equations (5.66) and (5.70) in (5.64), and using the Assumption 5.7, we get that

∆Π ≤ k+Np-1 t=k+1 (ψ(x t|k+1 , ȗt ) -ψ(x t|k , u t )) -ψ(x k , u k ) + χ(η) + L tx L Np-1 f x Λ(η) ≤ L ψx Np-2 j=0 L j f x Λ(η) -ψ(x k , u k ) + χ(η) + L tx L Np-1 f x Λ(η) ≤ -α ψ (|x|) + χ(η) (5.71) with χ(η) = χ(η) + L tx L Np-1 f x + L ψx Np-2 j=0 L j f x Λ(η) (5.72)
where χ is a K ∞ -function.

According to the results obtained in (5.56) (bounds on the optimal cost) and (5.71) (evolution of the optimal cost), the optimal cost (5.48) is a robust Lyapunov function (according to Definition 14). Finally, based on the Theorem 5.2, the system (5.36) controlled by π(x) = u k is robustly stable in Φ for any uncertainty θ ∈ Θ and for the considered linearized prediction model.

Unconstrained Linearized Robust Model Predictive Controller

The optimization problem (5.33)-(5.34) is solved by means of a robust regularized least squares approach in the presence of uncertain data, following an approach developed by Sayed et al. [START_REF] Sayed | A regularized robust design criterion for uncertain data[END_REF]. Since it is an unconstrained formulation, the control signal will be a posteriori saturated for real-time implementation. This approach is thus dedicated to the case where the constraints on the control inputs are bounds on it and there are no constraints on the state vector.

Robust regularized least squares problem

The different steps to solve a robust regularized least squares problem are presented below. • z is an unknown n-dimensional column vector;

• b is a known m × 1 vector;

• W 0 is a positive-definite weighting matrix.

Assume that the matrix A is tall (m ≥ n) and of full column rank (see Definition 3). Then, the OLS problem (5.73) has a unique solution, given by:

z = (A W A) -1 A W b (5.74)
However, this formulation can not be used for all problems (e.g. in the case when the matrix A is not full column rank). In order to overcome the drawbacks, a regularized least squares (RLS) can be used instead of OLS approach (5.73).

Let us consider the following RLS problem:

min z J (z) (5.75) with J (z) ||z|| 2 V + ||Az -b|| 2 W (5.76)
V 0 is a positive-definite weighting matrix.

The solution of the optimization problem (5.75) is given by

z = (V + A W A) -1 A W b (5.
||z|| 2 V + ||(A + δA)z -(b + δb)|| 2 W (5.78)
where uncertainties δA ∈ R m×n and δb ∈ R m can be structured under the following factored form:

δA = C∆E a (5.79) δb = C∆E b (5.80)
where ∆ denotes an arbitrary contraction term with ||∆|| ≤ 1, with a known matrix C ∈ R m×n ξ not identically null and where E a and E b are known quantities of appropriate dimensions.

In the sequel, the uncertainties δA and δb are replaced by a perturbation vector ξ ∈ R n ξ which is assumed to satisfy the following factored form:

Cξ = δAz -δb = C∆ (E a z -E b ) (5.81)
Since ||∆|| ≤ 1, ξ is therefore constrained as follows:

||ξ|| ≤ ||E a z -E b || Γ(z) (5.82)
The nonnegative function Γ(z) is assumed to be a known bound on the perturbation ξ and is a function of z only. Thanks to (5.81) and (5.82), the optimization problem (5.78) can be expressed as follows:

min z max ||ξ||≤Γ(z) ||z|| 2 V + ||Az -b + Cξ|| 2 W (5.83)
The maximization subproblem is transformed to a standard form, which will enable further to define the corresponding Lagrange dual problem: The constrained subproblem on ξ is solved by considering the Lagrangian duality [START_REF] Boyd | Convex Optimization[END_REF]. We define the Lagrangian L : R n × R n ξ × R + -→ R associated with the optimization problem (5.84) as

L(z, ξ, λ) -||z|| 2 V -||Az -b + Cξ|| 2 W + λ(||ξ|| 2 -Γ(z) 2 ) (5. 85
)
where λ is the Lagrange multiplier associated to the inequality constraint (5.82) on ξ. Consequently, the problem (5. Since L(z, ξ, λ) is a convex quadratic function of ξ, we can find an explicit solution of ξ which depends on the two variables z and λ by cancelling the gradient of the Lagrangian with respect to ξ, leading to:

ξ (z, λ) = (λI -C W C) † C W (Az -b) (5. 87
)
where I is the identity matrix with appropriate dimension. Due to the fact that the Hessian of the Lagrangian function (5.85) with respect to ξ must be non-negative at the optimum: with

∂ 2 L ∂ξ 2 = -C W C + λI 0 (5.
L(z, λ) = -||z|| 2 V -||Az -b|| 2 W (λ) -λΓ(z) 2 (5.91)
in which the modified weighting matrix W (λ) is derived from W via:

W (λ) = W + W C(λI -C W C) † C W (5.92)

Linearized Robust Model Predictive Control

Based on this formalism, the approximated min-max optimization problem which is defined by problem (5.1) with criterion (5.32), is written in the form (5.78)-(5.82) with:

min z max ||ξ||≤||Eaz-E b || ||z|| 2 V + ||Az -b + Cξ|| 2 W (5.102) and              z = u k+Np-1 k -u r,k+Np-1 k A = G k+Np-1 u,k b = y r,k+Np k+1 -G k+Np nom,k+1 C = G k+Np θ,k+1 ∆ = γ, E a = 0, E b = -δθ max (5.103)
The application of (5.97-5.101) provides the solution of (5.1) with criterion (5.32) as follows: step 1. The scalar λ is computed from the following minimization problem:

λ = arg min λ≥||G k+N p θ,k+1 W G k+Np θ,k+1 || J(λ) (5.104)
where the function J(λ) is defined by (see (5.101)):

J(λ) = ||z(λ)|| 2 V + λ||δθ max || 2 + ||G k+Np-1 u,k z(λ) -y r,k+Np k+1 + G k+Np nom,k+1 || 2 W (λ)
(5.105) with (from (5.97))

z(λ) = (V + G k+Np-1 u,k W (λ)G k+Np-1 u,k ) † (G k+Np-1 u,k W (λ)(y r,k+Np k+1 -G k+Np nom,k+1 )) (5.106)
and (from (5.92))

W (λ) = W + W G k+Np θ,k+1 (λI -G k+N p θ,k+1 W G k+Np θ,k+1 ) † G k+N p θ,k+1 W (5.107) step 2. The control sequence u k+Np-1 k
is then given by (from (5.106) for λ = λ ):

u k+Np-1 k =u r,k+Np-1 k + (V + G k+Np-1 u,k W (λ )G k+Np-1 u,k ) † × (G k+Np-1 u,k W (λ )(y r,k+Np k+1 -G k+Np nom,k+1 )) (5.108)
with W (λ ) given by (5.107) for λ = λ .

To summarize, the predictive controller consists in solving online an unidimensional optimization problem (5.104)-( 5.105) at each sampling time, instead of solving the min-max problem (5.1)-(5.3). In the sequel, this predictive control law will be referred to as Linearized Robust Model Predictive Controller (LRMPC). The LRMPC algorithm is summarized hereafter: LRMPC Algorithm Inputs: T s : sampling time; T d : integration time step; y r , u r : reference outputs and control inputs, respectively; x 0 : initial state vector; θ nom : nominal parameters; δθ max : maximum parameters uncertainties; N p : length of the prediction horizon; W , V : weighting matrices on the outputs and the control inputs, respectively;

Outputs: 1. Initialisation: k = 1 2. Update x k , y r,k+Np k+1 , u r,k+Np-1 k 3. Compute G k+Np-1 nom,k , G k+Np-1 θ,k , G k+Np-1 u,k
according to equations (5.19-5.21)

4.

Optimize λ by solving the unidimensional problem (5.104)

5.

Compute u

k+Np-1 k according to the equation (5.108) 6.

Apply u k to the system (4.6) after saturation between u min and u max 7.

Save

x k+1 8. k ←-k + 1 9. return to 2
The LRMPC problem is implemented as stated in Figure 5.3. Optimal control sequence

[u k (λ ), . . . , u k+Np-1 (λ )]
Application of u k to the system

Constrained linearized robust predictive controller

As mentionned before, the control problem that determines the control sequence has been formulated till now in an unconstrained case, with a posteriori saturation of the control signal for real-time implementation. In the sequel, inequality constraints on the control input u are taken into account in the optimization problem (e.g. corresponding to physical constraints on the actuators).

In this section, we will therefore consider the robust nonlinear predictive problem which is defined by (5.33)-(5.34) subject to the following inequality constraint on the optimization variable u k+Np-1 k :

u min ≤ u k+j-1 ≤ u max , j = 1, N p (5.109)
where u min and u max are the lower and upper bounds on the control signal, respectively. No constraints on the state are considered.

Based on the formalism presented in section 5.5.1, the optimization problem (5.1-5.3) with constraints (5.109) can be approached by the two player game problem (5.78). In the sequel, we extend results presented in [START_REF] Sayed | A regularized robust design criterion for uncertain data[END_REF] to the case of constrained robust regularized least squares. The main contribution is to solve the robust regularized optimization problem (5.78) taking into account constraints on the optimization variable z. Then, the unconstrained optimization problem defined in (5.93) with (5.94) becomes:

min λ≥||C W C|| min z≤z≤z J(z, λ) (5.110)
where the cost function J(z, λ) is still defined by (5.94). From (5.82), the optimization problem (5.110) is rewritten as:

min λ≥||C W C|| min z≤z≤z ||z|| 2 V + λ||E a z -E b || 2 + ||Az -b|| 2 W (λ)
(5.111)

Bilevel optimization problem

The corresponding formulation of the optimization problem (5.111) into a bilevel optimization problem can be written as follows:

λ = arg min λ≥||C W C|| J(z(λ), λ) s.t. z(λ) = arg min z≤z≤z J(z, λ) (5.112)
The problem (5.112) is made of two levels of optimization problems:

• Lower-level The minimum z(λ) is the solution of the following quadratic programming problem:

min z 1 2 z Hz + F z subject to I 0 0 -I z ≤ z -z (5.113) with H = 2 V (λ) + A W (λ)A F = -2 A W (λ)b + λE a E b (5.114)
where V (λ) and W (λ) are given by (5.100) and (5.92), respectively.

• Upper-level

The nonnegative scalar parameter λ ∈ R + is computed from the following unidimensional minimization problem:

λ = arg min λ≥||C W C|| ||z(λ)|| 2 V + λ||E a z(λ) -E b || 2 + ||Az(λ) -b|| 2 W (λ) (5.115)
Finally, the bilevel problem (5.112) has a unique global minimum z given by (5.113) for λ = λ (i.e. z = z(λ )).

Remark 5.7. Of course, if there are no constraints on z, the solution of problem (5.112) is equivalent to the one of the problem (5.101).

Constrained Linearized Robust Model Predictive Control

The application of (5.113)-( 5.115) provides the solution of (5.1) with criterion (5.32) and inequality constraints on u k+Np-1 k as follows:

Step 1. The scalar λ is computed from the following minimization problem:

λ = arg min λ≥||G k+N p θ,k+1 W G k+Np θ,k+1 || J(z(λ), λ) (5.116)
where the function J(z(λ), λ) is defined by (see (5.105)):

J(z(λ), λ) = ||z(λ)|| 2 V + λ||δθ max || 2 + ||G k+Np-1 u,k z(λ) -y r,k+Np k+1 + G k+Np nom,k+1 || 2 W (λ)
(5.117) and z(λ) is given by

z(λ ) = arg min z≤z≤z z E(λ)z -2B(λ) z (5.118)
where

             E(λ) = V + G k+Np-1 u,k W (λ)G k+Np-1 u,k
(5.119)

B(λ) = G k+Np-1 u,k W (λ) y r,k+Np k+1 -G k+Np nom,k+1
(5.120)

z = u min I nu -u r,k+Np-1 k (5.121) z = u max I nu -u r,k+Np-1 k (5.122)
with W (λ) given in (5.107).

Step 2. The control sequence u k+Np-1 k is derived from (5.118) for λ = λ :

u k+Np-1 k = u r,k+Np-1 k + z(λ ) (5.123)
To summarize, the predictive controller consists in solving online a bilevel optimization problem (5.112) instead of solving a min-max problem (5.1-5.3): a quadratic programming problem (5.113) in the lower level, and an unidimensional optimization problem (5.115) in the upper level. Since there are very efficient algorithms for this kind of optimization problems and that the two problems are convex, the obtained optimization problem remains more tractable than the min-max one. In the sequel, this predictive control law will be referred to as Constrained Linearized Robust Model Predictive Controller (CLRMPC). The CLRMPC algorithm is summarized hereafter.

CLRMPC Algorithm

Inputs: T s : sampling time; T d : integration time step; y r , u r : reference outputs and control inputs, respectively; x 0 : initial state vector; θ nom : nominal parameters; δθ max : maximum parameters uncertainties; N p : length of the prediction horizon; W , V : weighting matrices on the outputs and the control inputs, respectively;

Outputs: 1. Initialisation: k=1 2. Update x k , y r,k+Np k+1 , u k+Np-2 k-1 ; 3. Compute G k+Np-1 nom,k , G k+Np-1 θ,k , G k+Np-1 u,k
according to equations (5.19-5.21); 4.

Solve the bilevel optimization problem (5.112); 4.1.

Optimize λ by solving the unidimensional problem (5.116), by solving QP problem (5.118) for each λ; Apply u k to the system (4.6); 6.

Save x k+1 ; 7.

k ←k + 1 8. return to 2

The CLRMPC problem is implemented as stated in Figure 5.4.

It should be mentionned that the CLRMPC approach can be used in the case of more complex constraints on the control input (not only bounds). In this case, the problem (5.118) becomes NLP problem with quadratic cost function and nonlinear constraints. It can be solved with an SQP algorithm for example.

Numerical results and discussion

The following section illustrates the properties of both the RNMPC and the linearized formulation, using the bioprocess described in section 4.4. In fact, for this bioprocess, the rRNMPC cannot be used since the size of the parameters vector is small. This approach will be discussed in Chapter 7 on a more complex system. First, the tuning of the LRMPC is discussed in section 5.7.1. Secondly, the section 5.7.2 will compare unconstrained and constrained formulations of the LRMPC. Since the control is bounded (0 ≤ D ≤ D max = 0.5 h -1 ), an a posteriori saturation is applied to control inputs determined by LRMPC. Simulations have been carried out considering the uncertain parameter case discussed in the previous chapter (see section 4.4.2).

LRMPC tuning

In the sequel, the performances of the LRMPC under parameter uncertainties are analyzed, and more specifically the impact of its tuning parameters (T s and N p ) on the tracking accuracy with V = W = I Np . The same scenario of reference trajectory as in section 4.4 is considered. First, the results as depicted in Figures 5. 5-5.8 show the influence of the sampling time T s on the concentrations of biomass and substrate, tracking error and dilution rate evolutions with a prediction horizon N p = 6. The integration time step, T d , is chosen equal to T s /50. It can be observed that the smallest sampling time provides the best tracking accuracy for the closed-loop system (Figure 5.7). In fact, T s must be small to guarantee that the first order Taylor series expansion is accurate. It should be reminded that a compromise is required to properly select an appropriate sampling time taking into account the computation burden due to potential state estimator and/or online determination of the optimal trajectory. It appears clearly that for this example T s =5 min allows satisfying a good trade-off between linearization accuracy and computational burden.

Secondly, the choice of the prediction horizon N p is studied for the same conditions cited previously, with T s = 5 min as shown in Figures 5.9-5.12. The prediction horizon N p is chosen to satisfy a compromise between the computation time and a sufficient vision of the system behaviour. In fact, increasing N p leads to a loss of accuracy (quite small in this case) which is due to the prediction over the moving horizon using an approximated model (linearization) in the presence of model parameters uncertainties. It will nevertheless lead to a better anticipation to the future variation of reference trajectory. Hereafter, the prediction horizon is set to N p = 6, which represents the best trade-off between these two criteria for this example. In the case of the LRMPC strategy, the output converges to the setpoint value with an additional time d equal to 1.5 T s (as shown in the green box of Figure 5.13) in comparison with the CLRMPC algorithm. In fact, for CLRMPC the control is set quickly to its maximal value D max while the LRMPC control input takes longer time to converge to D max (see Figure 5.15). The LRMPC controller leads to suboptimal solution since it consists in solving online a scalar optimization problem (5.104-5.105) with an a posteriori saturation of the control input, whereas the constrained LRMPC law takes into account the constraints on the control inputs. On the other hand, the average computation time of the LRMPC is about fifty times less than the CLRMPC one as shown in Table 5.1. Both LRMPC and CLRMPC are equivalent outside the saturated zone (i.e. in case of inactive constraints).

Comparison of (C)LRMPC algorithms

In the case of bioprocess control, the operating points are usually outside the saturated zone in particular when addressing the issue of disturbance rejection. This is a motivation for using the LRMPC strategy in the sequel.

(s) X X X X X X X X X X X X Algo.

Comparison of predictive controllers

In this section, three predictive control laws will be tested (Figures 5. 16-5.18): a classical Nonlinear Model Predictive Control (denoted as NMPC), a robust one using criterion (5.1) (denoted as RNMPC) and the proposed one (LRMPC). The tuning parameters are the same for all strategies (N p = 6, T s = 5 min, T d = 0.1 min and V = W = I Np ). Both RNMPC and LRMPC are compared to the NMPC in order to quantify the impact of tracking error due to the parametric uncertainty. It can be noticed the anticipation behavior to a setpoint change (Figure 5.16), due to the prediction of the output future evolution over the moving horizon and the knowledge of the setpoint trajectory in the future. The obtained results show that both RN-MPC and LRMPC have better performances than NMPC under parameter uncertainties. In the NMPC law, the output is not able to track the specified setpoint in the presence of parameters uncertainties, due to the fact that the mismatch between the system and the model is not considered during the prediction step inside the minimization. The robust formulation tackles this drawback. In addition, the obtained results show that RNMPC has better performances than the LRMPC controller under parameter uncertainties in term of tracking accuracy. In the case of LRMPC, the static error is due to the approximation of the model through linearization. On the other hand, the LRMPC offers a very significant computational load reduction in comparison with the RNMPC as shown in Table 5.2. In fact, this can be explained by the fact that RNMPC is an optimization problem of dimension N p × n u × n θ while LRMPC is a unidimensional optimization problem. Consequently, when considering a more complex model with a greater number of state variables and parameters, the computation time increases quickly in the RNMPC strategy.

(s) X X X X X X X X X X X X Algo.
Perf. indices min mean max NMPC < 10 -5 0.02 0.32 RNMPC 0.06 0.31 6.77 LRMPC < 10 -5 0.01 0.1

Therefore LRMPC appears to be a good compromise between calculation time and accuracy. Moreover, a possible improvement to reduce the static error which is due to the approximation of the model through linearization, will be presented in the next chapter section 6.3.

Remark 5.9. From Tables 5.1-5.2, it appears that the mean computation time for CLRMPC at each sampling time when considering the worst uncertain parameter case is greater than in the RNMPC one. This result is mainly due to the fact that the RNMPC considers the worst case in the problem formulation (pessimistic approach) which is the case for the consid-ered model mismatch. This leads to a reduction of the mean computation time to solve the problem. In addition, the max subproblem of RNMPC is initialized from one step to another with the optimal parameters values of the previous optimization (see Remark 5.1), thus saving time, which is not the case of CLRMPC due to a different structure of the problem. On the contrary, in the nominal case (without parameters uncertainties) the mean computation time of CLRMPC at each time instant is less than the RNMPC one as shown in Table 5.3.

Concluding remarks

In this chapter, the robust MPC was investigated. Considering a system with parameters that are within a given confidence interval, the robust NMPC is designed in order to take into account these parameters uncertainties. It leads to a min-max optimization problem where the optimal control sequence is determined so that the maximum deviation for all trajectories overall possible data scenarii is minimized. To reduce the computation burden, we propose two alternative approaches to solve the min-max problem. In the first one, the computation time is reduced by considering only the most influencial parameters (so-called reduced RNMPC). It uses the sensitivity analysis of the model with respect to its parameters to determine the parameters to be considered in the min-max problem. In the second approach, we propose to linearize the predicted trajectory to turn the original optimization problem into a more tractable one. The derived optimization problem can be an unidimensionnal optimization problem (LRMPC) by applying the saturation on the control input a posteriori. Whereas taking into account the inequality constraints on the optimization variable in the design controller leads to a bilevel problem (CLRMPC), with a scalar optimization problem in the upper level, and a quadratic programming one in the lower level. The main advantage of LRMPC is to be computationally tractable in calculating the optimal control compared to both the bilevel problem and a classical min-max robust approach, which makes it suitable for online implementation (even if in our case, all calculation times remain much smaller than the sampling period, this should not be the same for a more complex system, with more uncertain parameters). Several simulations were performed in order to compare the LRMPC strategy to the classical RNMPC law in the case of model parameters uncertainties to control bacteria growth in a bioreactor.

The LRMPC was shown to ensure a good trade-off between computational load and tracking trajectory accuracy. However, there are several issues that deserve further investigation (e.g. non-zero tracking error). Some improvements of this control law are detailed in the next chapter.

Chapter 6

Some improvements of LRMPC

Introduction

In the previous chapter, in the LRMPC strategy, the min-max optimization problem was approached by a more tractable formulation, derived from the linearization of the system around nominal values of the parameters and reference control inputs. In this chapter, some improvements are proposed for the LRMPC controller.

A first proposed idea consists in linearizing the model around the nominal parameter values and the optimal control sequence obtained at the previous iteration instead of a reference control sequence determined at the steady state behavior. This modification is motivated by the fact that the optimal control sequence obtained at the previous iteration is non-model-based, and thus the controller robustness with respect to model uncertainties is increased. As in the case of LRMPC presented in Chapter 5, the derived optimization problem is here a bilevel one or a scalar optimization problem depending on the considered constraints. This controller will be referred as incremental LRMPC.

Furthermore, as highlighted in section 5.7.1, the LRMPC accuracy can be an issue (depending on the choice of the sampling period) as it deals with a first-order approximation of nonlinear functions. To go further with model uncertainties and linearization drawbacks, the idea here is to develop a hierarchical control scheme. The proposed control strategy combines a robust model predictive control law with a Proportional Integral (PI) law or Integral Sliding Mode (ISM) controller regulator. The predictive controller guarantees the tracking of the reference trajectory, whereas the other regulator (PI or ISM) ensures cancellation of any residual error.

The chapter is organized as follows. The section 6.2 deals with the de-sign of the incremental LRMPC. Section 6.3 presents the hierarchical control scheme that combines a robust predictive law with an auxiliary controller. Performances of the two controllers mentionned above will be illustrated on the illustrative example studied in the previous chapters. Finally, conclusions are summarized in section 6.4.

Variant of LRMPC

This section is motivated by the problem of designing a predictive controller in the presence of bounded uncertainties and inputs. In comparison with the previous chapter, the robust cost that will be minimized is expressed as a quadratic function measuring the control efforts instead of the difference between the control inputs and their reference values which are model-based. The idea is to use a model linearization along the nominal parameters value and the optimal control sequence obtained at the previous iteration. In addition to be more robust to errors in the determination of the reference trajectory than the previous approach, it presents the advantage to lead to a solution less sensitive to noise measurements, due to the fact that the control increments are explicitly considered in the cost function [START_REF] Tebbani | Nonlinear Predictive Control of fed-batch Cultures of E. coli[END_REF].

Problem formulation

The control sequence that minimizes a worst case cost function is derived from the following optimization problem (at time index k) instead of the optimization problem (5.1):

δ u k+Np-1 k = arg min δu k+Np-1 k ∈U δ max δθ∈Θ δ Π δ (x k , δu k+Np-1 k , δθ) (6.1)
Subject to

xk+j = g(t k , t k+j , x k , u k+j-1 k , θ = θ nom + δθ), j = 1, N p (6.2)
where the robust cost function is defined as

Π δ (x k , δu k+Np-1 k , δθ) ||δu k+Np-1 k || 2 V + ||ŷ k+Np k -y r,k+Np k || 2 W (6.
3)

The optimization variables δu k+Np-1 k are defined as the control increments as follows:

δu k+Np-1 k =         u k -u k-1 . . . u k+j -u k+j-1 . . . u k+Np-1 -u k+Np-2                δu k . . . δu k+j . . . δu k+Np-1        (6.4)
with u k-1 the control input applied at time index k -1 (i.e. solution of (6.1) at time index k -1).

As in Chapter 5, the predicted outputs over the receding horizon are given by: U δ and Θ δ are compact sets that contains the origin, defined as follows:

ŷk+Np k+1 =         Hg(t k , t k+1 , x k , u k , θ) . . . Hg(t k , t k+j , x k , u k+j-1 k , θ) . . . Hg(t k , t k+Np , x k , u k+Np-1 k , θ)         (6.
Θ δ = [-δθ max , δθ max ] with δθ max = (θ + -θ -)/2 (6.6) U δ = [δu, δu] (6.7)
V 0 and W 0 are tuning weighting matrices.

The predicted state for time t k+j , starting from state at t k , is linearized around the trajectory given for the nominal parameters θ nom and a control sequence ūk+Np-1 k defined as the optimal control sequence of the optimization problem (6.1) obtained at the previous sampling time (at time index k -1):

ūj = u j-1 , j = k, k + N p -1 (6.8)
As for LRMPC, a first order Taylor series (local) expansion of (4.5) for j = 1, N p is used:

xk+j ≈ xnom,k+j + ∇ θ g(t k+j )δθ + ∇ u g(t k+j ) u k+j-1 k -ūk+j-1 k (6.9) with                              xnom,k+j = g(t k , t k+j , x k , ūk+j-1 k , θ nom ) (6.10) ∇ θ g(t k+j ) = ∂g(t k , t k+j , x k , u k+j-1 k , θ) ∂θ u k+j-1 k = ūk+j-1 k θ = θnom (6.11) ∇ u g(t k+j ) = ∂g(t k , t k+i , x k , u k+j-1 k , θ) ∂u k+j-1 k u k+j-1 k = ūk+j-1 k θ = θnom (6.
12)

The dynamics of sensitivity function with respect to θ and u, defined in (6.11-6.12), can be computed for time t ∈ [t k , t k+Np ] as detailed in section 5.4.1. From (6.5) and (6.9), it comes:

                 xk+1 ≈ xnom,k+1 + ∇ θ g(t k+1 )δθ + ∇ u g(t k+1 )(u k -ūk ) . . . xk+j ≈ xnom,k+j + ∇ θ g(t k+j )δθ + ∇ u g(t k+j )(u k+j-1 k -ūk+j-1 k ) (6.13) . . . xk+Np ≈ xnom,k+Np + ∇ θ g(t k+Np )δθ + ∇ u g(t k+Np )(u k+Np-1 k - ūk+Np-1 k )
Using (6.4) and (6.8), the terms u k+j -ūk+j , j = 1, N p -1 are written as functions of ∆u j , j = 1, N p -1 as follows:

                                 u k -ūk = u k -u k-1 = δu k u k+1 -ūk+1 = δu k+1 u k+1 0 -u k + u k -u k-1 δu k + u k-1 -u k = δu k+1 + δu k + u k-1 -u k . . . u k+j -ūk+j = δu k+j + . . . + δu k + u k-1 -u k+j-1 (6.14) . . . u k+Np-1 -ūk+Np-1 = δu k+Np-1 + . . . + δu k + u k-1 -u k+Np-2
And thus:

u k+Np-1 k - ūk+Np-1 k =      δu k δu k+1 + δu k + u k-1 -u k . . . δu k+Np-1 + . . . + δu k + u k-1 -u k+Np-2      (6.15)
Finally, it comes

u k+Np-1 k - ūk+Np-1 k =    u k-1 . . . u k-1    -    u k-1 . . . u k+Np-2    +      I nu 0 . . . 0 . . . . . . . . . . . . . . . . . . 0 I nu . . . . . . I nu      δu k+Np-1 k Ξ k+Np-1 k + T Np δu k+Np-1 k (6.16)
where

Ξ k+Np-1 k =    u k-1 . . . u k-1    -    u k-1 . . . u k+Np-2    (6.17)
and T Np ∈ R Np×Np : unitary lower triangular matrix (6.18)

From (4.11), (6.9) and (6.16), the predicted ouputs ŷk+Np k+1 over the moving horizon are expressed as follows:

ŷk+Np k+1 = G k+Np nom,k+1 + G k+Np θ,k+1 δθ + G k+Np-1 u,k (Ξ k+Np-1 k + T Np δu k+Np-1 k ) (6.19) with G k+Np nom,k+1 , G k+Np θ,k+1 and G k+Np-1 u,k
defined as in section 5.4.1 (see (5.25)). Consequently, the linearized prediction outputs have similar structure as in the LRMPC (5.25), only the last term differs. This controller will be referred to as incremental LRMPC (LRMPC-δu).

Stability analysis

In the previous chapter, the stability analysis has been addressed for the closed-loop with the LRMPC law. The same methodology could be used for the LRMPC-δu under the same assumptions as in section 5.4.2.

The difference is that the model is linearized around θ nom and ū instead of θ nom and u r .

Consequently, the bound on the prediction error given by (5.47) as well as upper and lower bounds on the optimal cost given by (5.56) as still valid.

For the robust stability (as in section 5.4.2.3), the admissible control input is still chosen as:

ȗk+Np k+1 = [ u k+Np-1 k+1 , u k+Np-1 ] (6.20)
and the relation (5.71) is still valid. Consequently, considering the control sequence ū instead of u r do not change the stability of the controller (which is in this case mainly linked to the considered prediction model).

Derivation of the control law

Problem formulation

The constrained min-max optimization problem (6.1)-( 6.3) is approached by a robust RLS problem (6.21)-(6.22) when applying (5.29) and (6.19) in the presence of uncertain data, as follows:

min z≤z≤z max ||ξ||≤||Eaz-E b || ||z|| 2 V + ||Az -b + Cξ|| 2 W (6.21) with                  z = δu k+Np-1 k A = G k+Np-1 u,k T Np b = y r,k+Np k+1 -G k+Np nom,k+1 -G k+Np-1 u,k Ξ k+Np-1 k C = G k+Np θ,k+1 ∆ = γ, E a = 0, E b = -δθ max z = δu, z = δu (6.22)
where the perturbation vector ξ is assumed to satisfy the following form:

ξ = ∆(E a z -E b ) (6.23)
with ||∆|| ≤ 1 (6.24)

In this case, only bounds on the control input are considered (but the control strategy can be generalized to a more general constraints, i.e. z ∈ Z c ).

The obtained optimization problem (6.21)-(6.22) is similar to the one obtained when considering reference control inputs in the cost function (see (5.103)) if the correspondances summarized in Table 6.1 are considered.

Bilevel optimization problem

The constrained optimization problem (6.21)-(6.22) is solved by considering the Lagrangian dual problem of the maximization subproblem (i.e. the maximization of the error over all possible values of model parameters), following a similar approach as in section 5.6.2. Indeed, problem (5.102) has the same structure as problem (6.21) as mentionned previously. The solution of (6.21) with (6.22) is then given by: Lower-level: The optimal value of z = δu k+Np-1 k is derived from the following QP problem:

z(λ) = arg min z≤z≤z z E(λ)z -2B(λ) z (6.25)
u k+Np-1 k -u r,k+Np-1 k δu k+Np-1 k A G k+Np-1 u,k G k+Np-1 u,k T Np b y r,k+Np k+1 -G k+Np nom,k+1 y r,k+Np k+1 -G k+Np nom,k+1 -G k+Np-1 u,k Ξ k+Np-1 k C G k+Np θ,k+1 G k+Np θ,k+1 ∆ γ γ E a 0 0 E b -δθ max -δθ max where E(λ) = V + T Np G k+Np-1 u,k W (λ)G k+Np-1 u,k T Np (6.26)
and

B(λ) = T Np G k+Np-1 u,k W (λ) y r,k+Np k+1 -G k+Np nom,k+1 -G k+Np-1 u,k Ξ k+Np-1 k (6.27)
Upper-level: The scalar λ is computed from the following minimization problem: λ = arg min

λ≥||G k+N p θ,k+1 W G k+Np θ,k+1 || J(z(λ), λ) (6.28)
where the function J(z(λ), λ) is defined by (see (5.105)):

J(z(λ), λ) =||z(λ)|| 2 V + λ||δθ max || 2 + ||G k+Np-1 u,k T Np z(λ) -y r,k+Np k+1 + G k+Np nom,k+1 + G k+Np-1 u,k Ξ k+Np-1 k || 2 W (λ) (6.29)
with W (λ) given by (5.107).

The control sequence u k+Np-1 k is then given by:

u k+Np-1 k =    u k-1 . . . u k-1    + T Np z (6.30)
where z is obtained from (6.25) for λ = λ (i.e. z = z(λ )).

Remark 6.1. The same methodology could be applied for the unconstrained problem. Then, the resulting optimization problem will be an unidimensional problem as presented in section 5.5.2: step 1. The scalar λ is computed from the following minimization problem:

λ = arg min λ≥||G k+N p θ,k+1 W G k+Np θ,k+1 || J(z(λ), λ) (6.31)
where J(z(λ), λ) and W (λ) are given by (6.29) and (5.107), respectively, with (from (5.97))

z(λ) = V + T Np G k+Np-1 u,k W (λ)G k+Np-1 u,k T Np † T Np G k+Np-1 u,k W (λ) × ȳk+Np k+1 -G k+Np nom,k+1 -G k+Np-1 u,k Ξ k+Np-1 k (6.32) step 2. The control sequence u k+Np-1 k
is then given by (from (6.32) for λ = λ ):

u k+Np-1 k = [ u k-1 , . . . , u k-1 ] + T Np z(λ ) (6.33)
The saturation is then applied a posteriori to the first control input before its application to the plant.

Numerical results

The following section illustrates the properties of the modified LRMPC formulation using the process described in section 4.4. To illustrate its benefits with respect to measurement noise rejection, a uniform random noise of magnitude 0.2 g.L -1 at the maximum is added to the biomass measurements. The predictive tuning parameters are chosen as follows: N p =6, T s =5 min, T d = 0.1 min, V = I Np and W = 0.01I Np . In this example, more weight is given to V in order to promote the smoothness of the control input. The chosen values of V and W take into account the magnitude order of the control evolution in the cost function (smaller than the tracking error on the output). First, simulations have been carried out considering the nominal case (no uncertainties) in order to focus only on the effect of the measurement noise. The achieved results are shown in Figures 6.1-6.3. The obtained results in Figure 6.3 show that the control evolution is smoother when considering a penality on the control evolution (LRMPC vs LRMPCδu). The tracking performances are quite similar. Secondly, the achieved results shown in Figures 6.4-6.6 compare the performances of LRMPC and LRMPC-δu, with simulation conditions similar to those considered previously except we consider the uncertain parameter case discussed in the chapter 4 (see section 4.4.2). It can be observed that the control evolution is smoother for LRMPC-δu in comparison with the LRMPC (as shown in Figure 6.6) but the tracking accuracy on the biomass concentration is decreased (Figure 6.4).

Indeed, there is a trade-off between tracking accuracy and control smoothness. This decrease of accuracy observed with LRMPC-δu can be overcome if another tuning of V and W is considered.

Hierarchical control strategy

As mentionned previously, the proposed robust predictive approaches can have accuracy issues especially for large sampling period (since they deal with a first-order approximation of nonlinear functions). In order to better handle this approximation, an interesting solution to increase the quality of the linearized model, may be considering a second order expansion rather than the first order approximation, as proposed in [START_REF] Goffaux | Design of a robust nonlinear receding-horizon observer-Application to a biological system[END_REF] for robust state estimation.

Hereafter, we propose to use as an alternative strategy a hierarchical control scheme as depicted in Figure 6.7. This structure is similar to the Internal Model Control (IMC) (even if their principles are quite different) [START_REF] Garcia | Internal model control[END_REF]. The controller is formed by a robust predictive law (control input denoted u ) coupled to an auxiliary controller (control input denoted ũ). The predictive controller allows tracking the reference trajectory, whereas the additional controller is added to cancel any residual tracking error, also taking into account that this additional controller should be as simple as possible in order not to increase the calculation time too much. In the sequel, two controllers will be investigated: the Proportional-Integral (PI) law and the Integral Sliding Mode (ISM) controller.

Remark 6.2. The proposed hierarchical scheme can be used for all the predictive controllers (NMPC, RNMPC, rRNMPC, (C)LRMPC and LRMPCδu) but in this section the hierarchical approach will be only presented in the case of LRMPC.

PI controller design

The proposed control strategy consists in a hierarchical control scheme (see Figure 6.8). The controller is formed by a robust predictive law (LRMPC) coupled to a Proportional-Integral (PI) controller. The predictive controller allows tracking the reference trajectory, whereas the PI is added to reduce any residual tracking error. The motivation for choosing the PI controller is the simplicity of the control design in comparison with other advanced control strategies. In fact, the PI is a widely used control strategy that has its dynamics depending on the values of two constant gains: proportional and integral. At each time instant t k = kT s , the optimal control law u k obtained from the predictive controller is completed by a Proportional-Integral law ũ = u P I . The idea is to drive the system to track the predicted output ŷ in order to cancel the difference between the model prediction output and the system output. The PI design is performed in continuous-time, then discretized for implementation. Let us define the model prediction resulting from the application of the previous control input for t ∈ [t k-1 , t k ], as follows:

ŷ(t) = Hg(t k-1 , t, x k-1 , u k-1 , θ nom ) (6.34)
In the case of the considered control problem, the control law, ũ, derived from the PI strategy is given by (with K p and T i the controller parameters and are its tuning parameters):

ũ(t) = K p (y(t) -ŷ(t)) + K p T i t t k-1 (y(τ ) -ŷ(τ ))dτ (6.35)
Finally, with (6.35) evaluated at t = t k and a discretization via an Euler scheme, the control input to be applied to the plant is obtained as the sum of two parts, given by:

u(t k ) = u k + ũ(t k ) (6.36)
The component u k (the first value of the optimal control sequence) is generated by the LRMPC controller, while ũ(t k ) is generated by the PI controller.

The main drawback of this controller is the tuning of the PI controller (especially to ensure the closed-loop plant stability).

ISM controller design

An alternative to the implementation of a PI controller is to use the hierarchical control scheme formed by a Sliding Mode Control (SMC) law and a robust predictive controller as shown in Figure 6.9 similarly to the one proposed in [START_REF] Rugabotti | Robust model predictive control of continuous-time sampled data nonlinear systems with interval sliding mode[END_REF], which may in some cases improve performances compared to LRMPC-PI. The reasons for the choice of the SMC are:

• it is a robust control technique which guarantees the complete elimination of the tracking error due to the approximation of the model through the linearization;

• it offers a strong disturbance rejection;

• it guarantees the finite-time convergence to the sliding manifold.

Details related to SMC can be found in [START_REF] Utkin | Integral sliding mode in systems operating under uncertainty conditions[END_REF]. The SMC considered in this study is designed according to the so-called Integral Sliding Mode (ISM) approach. In the latter case, an integral action is included in the SMC. Let us assume that the considered system has a single output y and single input u and is control-affine, which is a special case of (4.1). We will then consider for the ISM development, the following model:

ẋ(t) = f x (x(t), θ) + f u (x(t), θ)u(t), ∀t > t 0 , x(t 0 ) = x 0 y(t) = Hx(t) (6.37)
At each time instant t = kT s , the goal is to complete the optimal control law u obtained from the predictive controller by an Integral Sliding mode control law ũ(t) in order to cancel the error between the predicted output and the system output (same philosophy as LRMPC-PI). In this study, the ISM design is performed in continuous-time, then discretized for implementation. The sliding mode control design consists in choosing the control input in such a way that the system is driven to reach a sliding manifold and be maintained there for all future time. The goal is to track the predicted output ŷ in order to cancel the difference between the model prediction output and the system output.

Let us define the modelling error variables for t

∈ [t k-1 , t k ] [145]:      Z 1 (t) = t t k-1 (y(τ ) -ŷ(τ )) dτ Z 2 (t) = 1 ξ 1 (y(t) -ŷ(t) -ξ 2 Z 1 (t)) (6.38)
where the model prediction ŷ resulting from the application of the previous control input u k-1 is given by (6.34). ξ 1 ∈ R and ξ 2 ∈ R are the ISM tuning parameters. Differentiating (6.38) with respect to time and rearanging (6.38), we obtain:

   Ż1 (t) = ξ 1 Z 2 + ξ 2 Z 1 Ż2 (t) = 1 ξ 1 ( ẏ(t) -ẏ(t) -ξ 2 Ż1 (t)) (6.39)
A time-varying sliding surface φ(x, t) is defined in the state space R nx as φ(x, t) = Z 2 (t) + ξ 3 Z 1 (t) (6.40)

The ISM control law needs to be designed so that the invariance of the sliding manifold is satisfied, i.e.:

∀x ∈ R nx , φ(x, t) = 0 (6.41)

Furthermore, the local attractivity of the sliding surface φ can be expressed by the condition: ∀x ∈ R nx φ(x, t)φ(x, t) < 0 (6.42)

Hence, the sliding surface (6.40) is made attractive by choosing [START_REF] Utkin | Integral sliding mode in systems operating under uncertainty conditions[END_REF]:

φ(x, t) = -K s sign(φ(x, t)) (6.43)
where the switching gain K s is a strictly positive constant. Then, the attractive equation which implies that the distance to the sliding surface decreases along all system trajectories is satisfied since (from (6.42)-(6.43)): φ(x, t)φ(x, t) = -K s |φ(x, t)| < 0 (6.44)

From (6.40) and (6.41), it comes:

Z 2 (t) + ξ 3 Z 1 (t) = 0 (6.45)
Thanks to (6.39), equation (6.45) becomes:

1 ξ 1 Ż1 (t) - ξ 2 ξ 1 Z 1 (t) + ξ 3 Z 1 (t) = 0 (6.46)
Then, we have that

Ż1 (t) = (ξ 2 -ξ 1 ξ 3 )Z 1 (t) (6.47)
In order to ensure the convergence of Z 1 , the following condition must be satisfied:

ξ 2 -ξ 1 ξ 3 < 0 (6.48)
Consequently, differentiating the sliding surface vector (6.40), we obtain:

φ(x, t) = 1 ξ 1 ( ẏ(t) -ẏ(t) -(ξ 2 -ξ 1 ξ 3 )(ξ 1 Z 2 (t) + ξ 2 Z 1 (t))) (6.49)
The system output is obtained by the application of the previous control input u k-1 (reminding that a ZOH is applied, i.e.

u(t) = u k-1 , ∀t ∈ [t k-1 , t k ])
combined with the sliding mode control law ũ(t) and the predicted output is calculated by applying only the previous input u k-1 as follows:

ẏ(t) = H (f x (x(t), θ) + f u (x(t), θ)(u k-1 + ũ(t))) ẏ(t) = H (f x (x(t), θ nom ) + f u (x(t), θ nom )u k-1 ) (6.50)
where the difference between the system output and the nominal model prediction output is due only to parameters uncertainties:

ẏ(t) -ẏ(t) = ϕ + χu k-1 + ηũ(t) (6.51) with    ϕ = H(f x (x(t), θ) -f x (x(t), θ nom )), χ = H(f u (x(t), θ) -f u (x(t), θ nom )), η = Hf u (x(t), θ) (6.52)
From (6.51), it is deduced that the control law can be accordingly found as:

ũ(t) = 1 η ẏ(t) -ẏ(t) -ϕ -χu k-1 (6.53)
Note that the term η must be regular. And from (6.49),

ẏ(t) -ẏ(t) = ξ 1 φ + (ξ 2 -ξ 1 ξ 3 )(ξ 1 Z 2 (t) + ξ 2 Z 1 (t)) (6.54)
By replacing φ by (6.43), the control law can then be found as

ũ(t) = 1 η (-ϕ -χu k-1 -ξ 1 K s sign(φ(x, t)) + (ξ 2 -ξ 1 ξ 3 )(ξ 1 Z 2 (t) + ξ 2 Z 1 (t)
)) (6.55) In (6.52), θ is replaced by θ nom since θ is uncertain. It could be replaced by an estimate of θ (with the use of an estimation algorithm). The use of θ in (6.55) is the main drawback of this approach. It should be recalled that the proposed controller can only be used to control-affine system which is another limitation of this strategy.

Another drawback of the ISM is the presence of a chattering phenomenon. In order to reduce it, a hyperbolic function can be used instead of the switching function sign(φ(x, t)).

Finally, with (6.55) evaluated at t = t k , the control input is obtained as the sum of two parts, given by:

u(t k ) = u k + ũ(t k ) (6.56)
The component u k (the first value of the optimal control sequence) is generated by the LRMPC controller, while ũ(t k ) is generated by the Integral sliding mode controller (6.55).

Numerical results

In this section, the efficiency of the proposed hierarchical control strategies is evaluated in simulation. The performances of both the LRMPC-PI and LRMPC-ISM are compared for the worst uncertain parameters case discussed in the previous chapter (see section 4.4.2). The sampling time T s is chosen equal to 5 min and T d the integration time step to 0.1 min. The controllers tuning parameters are given in Table 6.2. They were obtained by a trial-anderror technique (the best setting found). The obtained results as depicted in Figures 6.10-6.12 compare the LRMPC, the LRMPC-PI and the LRMPC-ISM controllers performances.

In the case of LRMPC law, it can be noticed the presence of a static error (Figure 6.10) due to the approximation of the model through linearization. On the other side, both the LRMPC-PI and LRMPC-ISM help to improve signifiantly the tracking accuracy in comparison with the LRMPC law. In fact, the dilution rate given by the hierarchical controllers moves away from its reference value D r nom , and converges to D r real as shown in Figure 6.12. In addition, it can be observed that the LRMPC-PI law leads to a larger overshoot in biomass concentration evolution during setpoint decrease in comparison with the LRMPC-ISM performances (Figure 6.10). The PI loop has a slower dynamic than the ISM one (Figure 6.12). Still, the overshoot is small and the LRMPC-PI presents a good global behavior. Thus, the hierarchical structure allows obtaining good closed-loop performances in terms of transient response and tracking accuracy.

Concluding remarks

In this chapter, two improvements have been proposed to the LRMPC controller. First, a new robust model predictive controller (denoted LRMPC-δu) is presented. The min max problem is turned into a more tractable optimization problem through the linearization of the predicted trajectory over the nominal parameters and the optimal contol sequence obtained at the previous sampling time. This leads to a control law less sensitive against noise measurements thanks to the inclusion of the penality term on the control evolution.

In a second step, a hierarchical control strategy is proposed. It consists in combining a robust predictive law (LRMPC) with an auxiliary controller (PI or ISM). The LRMPC allows tracking the reference trajectory, whereas the additional controller is added to cancel any residual tracking error. For the auxiliary controller, the choice of the PI presents the advantage of design simplicity in comparison with ISM controller which is a model-based control strategy. Nevertheless, the PI controller presents a slight loss of performances (response time) in comparison with the ISM regulator. To conclude, the PI law ensures the best trade-off between complexity of implementation and performances.

All control strategies presented till now will be tested in the case of a complex system (in comparison with the previous illustrative example). In the next chapter, the predictive controllers will be applied to control microalgae culture in a continuous photobioreactor.

Chapter 7

Illustrative example: Microalgae cultivation system

Introduction

One of the most important challenges of this century is to satisfy our energetic needs while addressing environemental issues. In the mean time, the influence of climate change, and the decline of fossil fuels have motivated a growing interest in the search of renewable sources of energy. The biotechnology industry is an important sector due to its economical, societal and environmental potential. Microalgae culture are among the most interesting types of biocultures. A microalga is a microscopic plant living mainly in aquatic environment, representing the basic level of the food chain in the ocean. These organisms have an increasing interest due to their ability to fix CO 2 [START_REF] Olaizola | Commercial development of microalgal biotechnology: from the test tube to the marketplace[END_REF] and to produce lipid or hydrogen that can be used as biofuel [START_REF] Chisti | Biodiesel from microalgae[END_REF][START_REF] Huntley | CO 2 mitigation and renewable oil from photosynthetic microbes: A new appraisal[END_REF][START_REF] Metting | Biodiversity and application of microalgae[END_REF]. The industrial success of the microalgae cultivation is due to its biochemical characteristics. The microalga is of particular interest for the growing demand of organic products intended to a large number of industrial applications: feed, food, pharmacology, chemistry, cosmetics production [START_REF] Enzing | Microalgae-based products for the food and feed sector: an outlook for Europe[END_REF][START_REF] Spolaore | Commercial applications of microalgae[END_REF], and has recently emerged as an interesting source for sustainable energy production at large scale, such as wastewater treatment [START_REF] Filipe | Biological Wastwater Treatment[END_REF] or decomposition of different classes of toxic compounds [START_REF] Boussiba | Microalgal biotechnology for environmental remediation[END_REF]. For these reasons, microalgae cultivation attracted the interest from large companies. In this context, the control strategy of microalgae culture is becoming a key research topic and consequently it received high attention from scientific community leading to many studies. This requires first the design of mathematical models to describe more accurately the system behavior. It is however quite delicate to derive a simple model which can be exploited during online optimization phases because the microalgae are living organisms and their characteristics can evolve in time (mutation, biofilm formation, stress, etc.). Biochemical processes are systems where nonlinear effects are significant enough to justify the use of nonlinear model to give a sufficiently adequate representation of the system behavior. In addition, generally, the process model is identified and uncertain parameters are estimated with evaluated confidence intervals, which motivates the development of robust control laws in the presence of modelling uncertainties.

In industrial applications, the control of bioreactors is most often limited to simple regulation loops like: partial pressure of dissolved oxygen and dissolved carbon dioxide [START_REF] Diaz | Adaptive predictive control of dissolved oxygen concentration in a laboratoryscale bioreactor[END_REF], pH and temperature. Hence, to increase significantly the process performances, new challenges emerge related to the control of the biological variables. In the literature of microalgae cultivations, several nonlinear control strategies have been developed: predictive [START_REF] Abdollahi | Lipid production optimization and optimal control of heterotrophic microalgae fed-batch bioreactor[END_REF][START_REF] Becerra-Celis | Nonlinear predictive control for continuous microalgae cultivation process in a photobioreactor[END_REF][START_REF] Tebbani | Nonlinear predictive control for maximization of CO 2 bio-fixation by microalgae in a photobioreactor[END_REF][START_REF] Santos | Nonlinear model predictive control of fed-batch cultures of micro-organisms exhibiting overflow metabolism: Assessment and robustness[END_REF][START_REF] Dewasme | State estimation and predictive control of fed-batch cultures of hybridoma cells[END_REF], adaptive [START_REF] Marcos | Output feedback adaptive extremum seeking control of a continuous stirred tank bioreactor with monod's kinetics[END_REF][START_REF] Cougnon | Online optimization of fedbatch bioreactors by adaptive extremum seeking control[END_REF][START_REF] Mailleret | Nonlinear adaptive control for bioreactors with unknown kinetics[END_REF], sliding mode [START_REF] Selisteanu | Sliding mode and adaptive sliding mode control of a class of nonlinear bioprocesses[END_REF], feedback linearization [START_REF] Becerra-Celis | Control strategy for continuous microalgae cultivation process in a photobioreactor[END_REF][START_REF] Ifrim | Multivariable feedback linearizing control of chlamydomonas reinhardtii photoautotrophic growth process in a torus photobioreactor[END_REF][START_REF] Tebbani | Nonlinear control of continuous cultures of Porphyridium purpureum in a photobioreactor[END_REF], and backstepping [START_REF] Dochain | Adaptive backstepping nonlinear control of bioprocesses[END_REF][START_REF] Toroghi | Observer based backstepping controller for microalgae cultivation[END_REF] approaches. They however do not specifically focus on robustness features. Our aim is therefore to design a robust predictive controller which would be able to find an optimal feeding strategy in order to guarantee that the process will yield the desired amount of biomass along the cultivation period under model parameter uncertainties. Here, the challenge is to lay down a stable real-time operation, insensitive to various disturbances, then, close to a certain state or desired profile. This requires the application of advanced optimal control strategies to ensure the bioprocess efficiency.

In this chapter, the proposed approaches developped in the previous chapters are applied to control microalgae culture in photobioreactor through a specific case study: the cultivation of I. galbana in a continuous photobioreactor [START_REF] Masci | Microalgal biomass surface productivity optimization based on a photobioreactor model[END_REF]. A Droop model is used to describe the internal nutrient quantity per unit of biomass evolution. The aim here is to control biomass concentration in the photobioreactor through the dilution rate.

System modelling

Model Description

The specificity of microalgae is that inorganic substrate uptake and growth are decoupled thanks to an intracellular storage of nutrients [START_REF] Bernard | Hurdles and challenges for modelling and control of microalgae for CO 2 mitigation and biofuel production[END_REF]. In order to take into account this phenomenon, the growth of microalgae is represented by a Droop model [START_REF] Bernard | Transient behavior of biological loop models, with application to the droop model[END_REF][START_REF] Bernard | Hurdles and challenges for modelling and control of microalgae for CO 2 mitigation and biofuel production[END_REF] which uncouples growth from substrate uptake, leading to the definition of an internal cell quota (i.e., the internal nutrient quantity per unit of biomass), and describes the growth rate as a function of the internal quota only. The mass balance model (see Appendix B) involves three state variables: the biomass concentration (denoted X, in µm 3 L -1 ), the internal quota (denoted Q, in µmol µm -3 ), and the substrate concentration (denoted S, in µmol L -1 ). The considered dynamic model assumes that the photobioreactor is in continuous mode (medium withdrawal flow rate equals its supply one, leading to a constant effective volume), without any additional biomass in the feed and neglecting the effect of gas exchanges. The time varying equations resulting from mass balances are given by [START_REF] Masci | Microalgal biomass surface productivity optimization based on a photobioreactor model[END_REF]: • S in the input substrate concentration (µmol L -1 ).

   Ẋ(t) = µ(Q(t), I(t)) X(t) -DX(t) Q(t) = ρ(S(t)) -µ(Q(t), I(t)) Q(t) Ṡ(t) = (S in -S(t))D -ρ(S(t))X(t)
The specific uptake rate ρ(S) is given by a Monod kinetic:

ρ(S) = ρ m S S + K s (7.
2)

The parameters K s and ρ m represent respectively the substrate half saturation constant and the maximal specific uptake rate.

The specific growth rate µ(Q, I), on the other hand, can be defined as a function of the intercellular quota Q as follows (Droop model):

µ(Q, I) = μ 1 - K Q Q µ I (I) (7.3)
The theoretical specific maximal growth rate (at infinite internal quota) is denoted μ and K Q represents the minimal cell quota, for which no algal growth can take place.

The light intensity has a direct effect on growth (photosynthesis), while uptake can continue in the dark. The modelling of light effect consists in including the term µ I in (7.3) which is represented by a Haldane type kinetics to model the photoinhibition [START_REF] Peeters | The relationship between light intensity and photosynthesis[END_REF]:

µ I (I) = I I + K sI + I 2 K iI (7.4)
where I is the light intensity (µE m -2 s -1 ) and K sI and K iI are light saturation and inhibition constants respectively. The optimal light intensity that maximises the function µ I is given by I opt = √ K sI K iI . In the sequel, the light intensity is either set at this optimal value I opt or is time varying. The parameters of the model used in this study are given in Table 7.1 [START_REF] Goffaux | Design of a robust nonlinear receding-horizon observer-Application to a biological system[END_REF], [START_REF] Munoz-Tamayo | Getting the most out of it: optimal experiments for parameter estimation of microalgae growth models[END_REF].

ρ m 9.3 µmol µm -3 d -1 K Q 1.8 µmol µm -3 K s 0.105 µmol L -1 K sI 150 µE m -2 s -1 K iI 2000 µE m -2 s -1

Model analysis

The nonlinear model (7.1) is represented thereafter in the state-space formalism (4.1) as follows:

ẋ(t) = F (x, u, w, θ), x(t 0 ) = x 0 y = X (7.5) with                    x =   X Q S   , w = S in I , u = D F =   µ(Q, I)X -DX ρ(S) -µ(Q, I)Q (S in -S)D -ρ(S)X   θ = ρ m K s μ K Q K sI K iI T (7.6)
where the state variables are assembled in a vector denoted x and x 0 its initial value. The nonlinear process dynamics is denoted F . The measurements are related to the vector y whereas the inputs are represented by the vector u. The other exogenous inputs are denoted w. Finally, the parameters refer to the vector θ. The maximum growth rate denoted µ m is reached for the maximum internal quota Q l . From (7.3), it comes:

µ m = µ(Q l , I) = μ 1 - K Q Q l µ I (7.7)
From the dynamics of Q, the maximal specific uptake is related to µ m by:

ρ m = µ m Q l (7.8)
From (7.7)-(7.8), the maximum cell quota Q l obtained in conditions of nonlimiting nutrient is given by:

Q l = ρ m μµ I + K Q with μµ I = 0 (7.9)
The condition μµ I = 0 holds ∀I = 0 which is the case here.

In addition, the cell quota is constrained to remain between two bounds [START_REF] Bernard | Hurdles and challenges for modelling and control of microalgae for CO 2 mitigation and biofuel production[END_REF]:

K Q ≤ Q ≤ Q l (7.10)
The growth rate is also bounded:

0 ≤ µ(Q) ≤ µ m (7.11)
where (by replacing (7.9) in (7.8))

µ m = ρ m μµ I ρ m + K Q μµ I (7.12)
Finally, the state and control variables are restricted to fulfill the following constraints [START_REF] Bernard | Transient behavior of biological loop models, with application to the droop model[END_REF]:

         X > 0 K Q ≤ Q ≤ ρ m µ I μ + K Q 0 ≤ S < S in D ≥ 0 (7.13)

Determination of equilibrium

The steady states of the system are derived from three nonlinear equations, given in order to cancel out the model's dynamic equations, i.e. F (x e , u e , w, θ) = 0 (7.14) where w is assumed constant and known). For a given value of X, the goal here is to characterize the corresponding values for Q, S and D for a constant light intensity I. Then, the equilibrium points are defined as follows (from (7. Then, rearranging (7.15), the following system of equations has to be solved algebraically :

                       D e = μ 1 - K Q Q e µ I S e = µ(Q e , I)Q e K s ρ m -µ(Q e , I)Q e S in - µ(Q e , I)Q e K s ρ m -µ(Q e , I)Q e -Q e X e = 0 (7.16)
Taking Q e as an unknown variable, the following quadratic equation must be solved:

[μµ I X e ]Q 2 e -[(S in + K s )μµ I + (ρ m + μµ I K Q )X e ]Q e + [(ρ m + μµ I K Q )S in + μµ I K Q K s ] = 0 (7.17)
The discriminant of equation (7.17) is given by:

∆ = aX 2 e + bX e + c = ( √ aX e - √ c) 2 + (2 √ a √ c + b)X e (7.18)
where

   a = (ρ m + μµ I K Q ) 2 > 0 b = 2μµ I ((K s -S in )(ρ m + μµ I K Q ) -2μµ I K Q K s ) c = μ2 µ 2 I (S in + K s ) 2 > 0 (7.19)
We have that

2 √ a √ c + b =2(ρ m + μµ I K Q )μµ I (S in + K s ) + 2μµ I ((K s -S in )(ρ m + μµ I K Q ) -2μµ I K Q K s ) 2 √ a √ c + b =4μµ I K s ρ m > 0 (7.20) Finally, it holds ∆ = ((ρ m + μµ I K Q )X e -μµ I (S in + K s )) 2 + 4μµ I K s ρ m X e (7.21)
Since the discriminant ∆ is strictly positive (and reminding that X e > 0, see (7.13)), there are two real solutions for (7.17):

Q * 1,2 = S in + K s 2X e + 1 2
ρ m μµ I + K Q ∓ √ ∆ 2μµ I X e (7.22)
By using equations (7.9), (7.18) and (7.19), we obtain:

√ ∆ 2μµ I X e = 1 2 Q l 1 + α 1 X e + β 1 X 2 e (7.23) with        α = b a β = c a > 0 (7.24)
Then (7.22) becomes:

Q * 1,2 = S in + K s 2X e + 1 2 Q l 1 ∓ 1 + α 1 X e + β 1 X 2 e (7 .25)
From (7.13) the choice of the equilibrium X e must be chosen to fulfill the following conditions:

• Condition Q * 1,2 ≥ 0:      Q * 1 ≥ 0 iff X e -X 2 e + αX e + β ≥ - S in + K s Q l (7.26) Q * 2 ≥ 0, ∀X e (7.27) • Condition Q * 1,2 ≤ Q l :          Q * 1 iff X e + X 2 e + αX e + β ≥ S in + K s Q l (7.28) Q * 2 iff X e -X 2 e + αX e + β ≥ S in + K s Q l (7.29) • Condition 0 ≤ S e < S in : ρ m ≥ µ(Q e , I)Q e (7.30)
which is equivalent to have Q e ≤ Q l (and hence Q * 1,2 ≤ Q l which is the previous condition).

Control strategy

The main objective of the controller is to regulate the biomass concentration X to a reference value X r in the presence of parameters uncertainties and noise measurement, while the dilution rate D is constrained to track the reference D r where 0 ≤ D ≤ D max (D max is the maximal dilution rate).

The controllability of this model was checked (See Appendix C). The NMPC problem is then formulated as:

min 0≤D k+Np-1 k ≤Dmax ||D k+Np-1 k -D r,k+Np-1 k || 2 V + ||X k+Np k+1 -X r,k+Np k+1 || 2 W (7.31)
and the RNMPC by: min

0≤D k+Np-1 k ≤Dmax max θ∈[θ -,θ + ] ||D k+Np-1 k -D r,k+Np-1 k || 2 V + ||X k+Np k+1 -X r,k+Np k+1
|| 2 W (7.32) where the uncertain parameters subspace [θ -, θ + ] is given by [0.7θ nom , 1.3θ nom ] with θ nom given in Table 7.1. This 30% mismatch has been chosen as a rather classical percenttage. A more rigourous approach could be to proceed with an identification procedure as in [START_REF] Breitenecker | Parameter identification of the droop model using optimal experiment design[END_REF] to determine the confidence intervals for all parameters.

The performance of the controllers in case of disturbances (here the light intensity fluctations) will be also studied and evaluated.

Simulation results

The efficiency of the proposed control strategies is validated in simulation. The initial biomass concentration value is set close to the setpoint in order to cancel the transient effect and focus only on the behavior during setpoint changes (rising and falling edge respectively). The light intensity is assumed to be non-measured, non-corrupted with noise. In sections 7.4.2 and 7.4.3, it is constant, equal to I opt that maximizes µ(Q, I) defined in (7.3). In section 7.4.4, it is time-varying. The performances of the predictive algorithms presented in chapters 4, 5 and 6 are compared in a worst uncertain parameters case. Thanks to the monotonocity properties of the Droop model as discussed in [START_REF] Goffaux | Design of a robust nonlinear receding-horizon observer-Application to a biological system[END_REF], the worst-case prediction can be generated using parameters bounds {θ -, θ + } only, rather than by exploring the full parameter space [θ -, θ + ]. The parameters values of the system are chosen on the parameter subspace border (θ

real = [ρ + m , K - s , μ+ , K - Q , K - sI , K + iI ]
). All the simulation conditions are summarized in Table 7.2.

Determination of the reference

Two configurations may be considered for the determination of the pair (D r , X r ) of reference signals:

• Case 1: D r -→ X r
The biomass reference trajectory X r is obtained by applying the dilution rate reference trajectory D r to the model.

• Case 2: X r -→ D r In case of constant X r (with an assumed constant light intensity), the dilution rate reference trajectory D r is computed from the knowledge of the target setpoint X r using relations at equilibrium (7.16):

D r = μ 1 - K Q Q r µ I (7.33)
with Q r the reference internal quota given by (7.25). From (7.25), two admissible solutions are possible for a given X r , depending on the chosen setpoint.

In case of time-varying X r or X r constant with a time-varying light intensity, the dilution rate reference trajectory D r could be determined by solving the following constrained open-loop optimization problem:

D r (t) = arg min 0≤D(t)≤Dmax |X(t) -X r (t)| 2 (

Setpoint tracking

First, the performances of setpoint tracking is studied. The goal is to track biomass concentration setpoint (rising and falling step changes as illustrated in Figure 7.4). Four predictive control laws are tested (Figure 7. It can be noticed the anticipation behavior to a setpoint change (Figure 7.4) for all controllers, due to the prediction of the setpoint trajectory future evolution over the moving horizon. The dilution rate decrease leads to an increase of the cell concentration (e.g. between 0.2 d and 0.3 d in Figure 7.5) and vice versa (e.g. between 0.7 d and 0.8 d), which agrees with the biological aspect. The obtained results show that both (r)RNMPC and LRMPC have better performances than the classical NMPC under parameter uncertainties and measurement noise. In the NMPC law, the biomass concentration is not able to track the specified setpoint in the presence of parameters uncertainties, due to the fact that the mismatch between the system and the model is not considered during the prediction step inside the minimization procedure. Furthermore, the RN-MPC has better performances than the rRNMPC and LRMPC controllers under parameter uncertainties in term of tracking accuracy. The LRMPC performs better than the rRNMPC since it considers all model parameters. Furthermore, the classical RNMPC has better results than the reduced one but the computational burden of the former is much higher. The LRMPC algorithm performs well and offers a very significant computational load reduction comparing with (r)RNMPC as shown in Table 7.3.

Remark 7.1. A comparison between NMPC-ε s/m and Generic Model Control (GMC [START_REF] Lee | Generic model control (GMC)[END_REF]) is provided in Appendix C. The latter is studied since it is widely used for bioprocess control. The predictive controller presents better performances than the GMC law, highlighting the benefits of the proposed controller. For the same simulation conditions as previously, the robust predictive controllers with penality term on the control variation (6.3) (denoted LRMPCδu) and the initial formulation (5.33) are compared. The obtained results in Figure 7.9 (green box) show that the control evolution is slightly smoother when considering control increments (6.3) in comparison with the LRMPC. Thus, the LRMPC insensitivity against measurement noise is improved thanks to the penality term on the control evolution in the cost function. As mentionned previously, the LRMPC or LRMPC-δu accuracy in tracking the setpoint is lower than the one obtained by RNMPC (Figure 7.4). Thus, a hierarchical strategy with a PI and with ISM controllers as detailed in section 6.3.1, are applied.

The PI control input in this case is given by (from (6.35)):

D(t) = K p (X(t) -X(t)) + K p T i t t k-1 (X(τ ) -X(τ ))dτ (7.36)
where X denotes the biomass model prediction.

The control law derived from the ISM strategy is given by (6.55) (as detailed in section 6.3.2):

D(t) = - 1 X(t) ((u k-1 -µ(Q, I))(X(t) -X(t)) -ξ 1 K s tanh(φ(x, t)) + (ξ 2 -ξ 1 ξ 3 )(ξ 1 Z 2 (t) + ξ 2 Z 1 (t))) (7.37)
Both PI and ISM tuning parameters are determined by a trial-and-error technique (see Table 7.4). The use of a hierarchical scheme allows cancelling the static error in comparison with LRMPC law as shown in Figure 7.12 (the obtained accuracy for LRMPC-PI is equivalent to the one obtained with RNMPC in Figure 7.4).

The obtained results show that the LRMPC-ISM has slightly better performances than the LRMPC-PI (lower overshoot and lower time response). Nevertheless, in the sequel, the performances analysis will be focused on the LRMPC-PI since it presents the best trade-off between performances and simplicity of implementation as mentionned in section 6.3.1. Figure 7.16 illustrates more spcifically the results obtained by the LRMPC-PI (biomass concentration, and predictive and PI control inputs). It can be observed that the LRMPC-PI control input is non-zero during setpoint change (rising edge) due to the fact that the PI control is not cancelled. The LRMPC, thanks to it predictive property, cancels the dilution rate so that the growth is maximized, whereas, the PI does not take into account future reference evolution. As a consequence, the time response of the LRMPC-PI is longer than the one achieved by the LRMPC. A statistical analysis of the robustness is considered in order to highlight the advantage of the hierarchical approach. To emphasize this aspect, simulation for a large number of independent tests are performed with parameters variations for the same conditions as previously. Based on Monte-Carlo procedure, 100 tests have been conducted with a simultaneous random non-correlated variation of 30% at the maximum in all the parameters θ. Figure 7.17 compares NMPC, LRMPC and LRMPC-PI performances. It clearly appears that using the hierarchical algorithm (solid line in blue) allows reducing the standard deviation of tracking error in comparison with a classical NMPC (dotted line in red) and LRMPC (dashed line in green). 7.5 gives the mean and standard deviation of tracking error for the whole simulation.

These results confirm those obtained in the worst case (Figure 7.12). Furthermore, the results show that the LRMPC has better results than the NMPC in term of accuracy. It can be seen that including the difference between the system and the model during the prediction in the NMPC law does not improve the tracking accuracy (Figure 7.19). This is due to the fact that the hypothesis on the dynamics of the error between the system and the model, modeled with jε s/m , is not realistic in the case of time-varying reference trajectory.

On the other hand, both RNMPC and the hierarchical approaches (LRMPC-PI and LRMPC-ISM) have similar performances and lead to a slight dynamic error. The RNMPC leads to the most accurate tracking but presents strong computation time. The above case study highlights the limits to the applicability of the NMPC with jε s/m signal.

Disturbance rejection

The light intensity was set constant in the previous simulations (equal to I opt ). In this section, the behavior of the controller in case of time-varying light intensity is discussed. Hereafter, a day/light like variation is considered [START_REF] Tebbani | Model-based versus model-free control designs for improving microalgae growth in a closed photobioreactor: Some preliminary comparisons[END_REF]:

I real (t) = Ī(max{0, sin(2πt)}) 2 + Iopt (7.38)
where the time t is in days and Ī is set to 280 µE m -2 s -1 . I nominal = I opt represents the light energy provided by panels to the bioreactor and I real the perturbation modelled as a day/night cycle (i.e. non perfectly isolated culture) as shown in Figure 7.25. The biomass concentration setpoint is constant (X r = 25µm For both the controllers, when the light I varies between times t = 0.3 and t = 0.7 d (Figure 7.25) leading to values higher than the optimal light intensity I opt , and thus, the specific growth rate µ(Q, I) decreases (inhibition phenomenon as shown in Figure 7.26). As a consequence, the dilution rate D also decreases accordingly in order to maintain the biomass concentration X approximately constant (Figure 7.26) (reminding that D must equals µ to cancel Ẋ). The reference tracking presents an offset in the case of LRMPC as shown in Figure 7.23. On the other side, the hierarchical approach (LRMPC-PI) maintains the biomass concentration at its reference value in the presence of parameters uncertainties and counters the effect of fluctuations in light intensity. In the case of LRMPC-ISM, the non-zero tracking error is due mainly to the use of µ(Q, I), evaluated at θ nom , in (7.37). This is the main drawback of the LRMPC-ISM. The specific growth rate, µ, could be estimated online, and its estimated value used in the controller.

Conclusion

In this chapter, predictive controllers were applied to control microalgae culture in a continuous photobioreactor (Droop model). First, the obtained results show that both RNMPC and LRMPC have better than the classical NMPC performances under parameters uncertainties and measurements affected by noise. The LRMPC accuracy in tracking the setpoint is lower than the one obtained by RNMPC due to the linearization of the prediction model. However, the RNMPC suffers from the difficulty to derive a solution as quickly as possible. The LRMPC offers a strong computational load reduction. Secondly, the use of a hierarchical scheme allows cancelling the static error in comparison with LRMPC law. For a time-varying reference trajectory, the hierarchical control scheme ensures a good tracking accuracy which is not the case for the NMPC-jε s/m . Thus, it highlights the limits to the applicability of the latter approach. Finally, for time-varying light intensity, the hierarchical approach (LRMPC-PI) maintains the output at its reference value in the presence of parameters uncertainties and counters the effect of fluctuations in light intensity. This example illustrated and highlighted the advantages of the proposed predictive controllers. More specifically, the hierarchical structure that combines a linearized robust predictive controller with a PI controller proposes a good solution for microalgae cultivation control.

The applicability of this control strategy to other bioprocesses and its possible theoretical and numerical improvements will be discussed in the next chapter.

Chapter 8

General conclusions and future directions

This thesis has investigated new structures of robust nonlinear predictive control, which are able to deal with the compromise between robustness features, tracking accuracy and dynamic performances, on the one hand, and reduced computational load on the other hand. This aspect of computational burden is indeed of major concern, when applying the control strategy in real-time. Even if the bioprocess considered in this work for the application framework has limited real-time constraints, this may not be the case in particular for systems with higher number of uncertain parameters. Therefore the aspect is crucial for RNMPC.

Based on this, the proposed structure, combining a linearized MPC controller with an additional one to deal with residual tracking accuracy, appears to be a good candidate to the paradigm mentioned above.

This development has been elaborated in a progressive way, through the contribution of successive chapters, and assessed through simulations results in the case of bioprocesses.

In this chapter, we summarize the thesis contributions and present recommendations for future work.

Thesis summary and contributions

Chapter 3 reviewed previous works in the area of predictive control strategies and presents their key advantages and drawbacks. Theoretical, computational and implementation aspect were briefly discussed. It served as a centralized literature study in the thesis.

Chapter 4 made use of the Nonlinear Model Predictive Control (NMPC) discussed in Chapter 3 and appplied it in the context of trajectory tracking problem for continuous/discrete-time nonlinear systems. This strategy induces solving, online, an optimization problem which is expressed as a nonlinear programming problem. This formulation highlighted its limits in the case of plant-model mismatch. One studied improvement consisted to include, during the prediction step, the difference between the system and the model outputs. Unfortunately, the assumption on including the prediction error is quite restrictive due to the fact that we do not know, in a realistic application, how the error will evolve over the prediction horizon.

Chapter 5 focused on the Robust Model Predictive Control. Robust variants of NMPC are able to take into account set bounded uncertainties in the design procedure. In this study, the robust NMPC law was formulated as a min-max optimization problem where the objective function is minimized for the worst possible uncertainty realization. In order to reduce the computation burden, two solutions have been proposed. On one side, thanks to a sensitivity analysis, a reduction of the number of uncertain parameters is performed. Only the main influencial parameters on the model are considered in the min-max problem whereas the other parameters are kept constant and equal to their nominal values. This approach cannot be always used due to two major drawbacks: the necessity to consider a study of the respective impact of the parameters through a sensitivity analysis, leading to a more complex and time demanding study and its applicability only in some cases (when all parameters have high influence on the model, no parameter reduction is possible). On th other side, a new approach based on a model linearization technique (first-order Taylor series expansion) has been developed. It aims to turn the min-max problem into a more tractable optimization problem in order to reduce the computation time as much as possible and to make it suitable for online implementation. The derived optimization problem is an unidimensional optimization (so-called LRMPC). Taking into account the inequality constraints on the optimization variable leads to a bilevel problem (so-called CLRMPC), with a quadratic programming problem in the lower level, and a scalar one in the upper level. This approach leads to a loss of accuracy (linearizaton drawback) but it is computationally tractable in calculating the optimal control compared to the min-max optimization, which makes it suitable for online implementation. Moreover, stability properties of the proposed LRMPC approach has been analyzed.

Chapter 6 proposed two improvements for the developed method. The first idea deals with the linearization of the model around the nominal parameter values and the optimal control sequence obtained at the previous iteration. This modification led to an optimal control sequence less sensitive ato noise measurements thanks to the inclusion of the penality term on the control evolution. The derived optimization problem is a bilevel one or a scalar minimization, depending on the considered constraints. In a second step, to go further with linearization drawback, a hiearchical control scheme has been proposed. It combines a robust model predictive control law with a proportional Integral (PI) law or Integral Sliding Mode (ISM) controller. The predictive controller guarantees the tracking of the reference trajectory, whereas the additional regulator ensures cancelling any residual tracking error. For the auxiliary controller, the PI law ensures the best trade-off between design complexity and performances, without requiring specific knowledge on the system parameters.

Chapter 7 was devoted to the case study of the control of microalgae culture in continuous photobioreactor. Simulations have been carried out in the case of constant and time-varying reference trajectories to highlight the advantage of the proposed control approach under parameters uncertainities, external perturbation and measurement affected by noise.

Recommendations for future directions

In this work, the obtained results have highlighted many issues that have to be considered as perspectives of this work. Some recommendations for future work in which further research may be carried out are proposed below, with different time scales.

Short term issues are more related to the application field considered in this work:

• In this dissertation, the proposed control strategies were applied in the case of Single Input Single Output (SISO) system with a small number of uncertain parameters. The application of the proposed methods to a highly nonlinear multivariable system with much more uncertain parameters is under progress (based on the experimental work of [START_REF] De | Optimization of biomass and lipid production in heterotrophic microalgal cultures[END_REF]). The studied MIMO sytem is a heterotrophic microalgae fed-batch bioreactor. The aim is to control biomass and product concentrations through two substrate feed rates. This objevtive has to be finalized.

• Based on results on the robustness of the LRMPC with respect to model parameters uncertainties, future work can consider extending the study to the case of disturbances and unmodeled uncertainties (e.g. model mismatch of the growth kinetics). In fact, these unmodeled mismatches in the kinetic part could come from the difficulty to select an accurate model structure, a possible evolution of the model parameters and/or a metabolic change over time.

• The biological variables are sometimes not accessible to be measured online and usually measured offline using expensive sensors. In this context, it is important to design an estimation algorithm or so-called soft sensors to rebuild time evolution of the state and to develop a constructive procedure for designing controllers robust against additional estimation errors. Contributions in this area are particularly important.

On top of that, the experimental validation of the proposed robust hierarchical structure on a real bioreactor, either lab-scale or industrial, would be the best validation to fully assess the theoretical developments. This is however a complex task, starting first, for a dedicated benchmark, with the model definition, the parameters identification with their confidence intervals, then the implementation of the control law (coupled to soft sensors) in a dedicated software. This perspective is a really important one to demonstrate the efficiency and advantages of the proposed control strategies.

Long term issues are more related to theoretical developments:

• In order to improve the performances of the LRMPC-δu in the presence of model mismatch and measurements affected by noise (section 6.2.4), an alternative would be to include the penalty term on the control variation (defined in (6.3)) as an additional term in the criterion (5.3) (i.e. to consider simultaneously penalty on the control increments and on the difference between u and its reference u r ). We already applied this strategy in the case of NMPC law, leading to good performances (see example in Appendix C). This new criterion should be considered for LRMPC and the control input calculation extended to this new optimization problem.

• In order to increase the quality of the linearized model, there are several issues that deserve further investigation. An interesting perspective may be considering a second order expansion rather than the first order approximation to improve the robustness and accuracy of the LRMPC as in [START_REF] Goffaux | Design of a Robust Nonlinear Receding-Horizon Observer -First-Order and Second-Order Approximations[END_REF].

• In particular, theoretical questions like stability need to be addressed. In fact, stability analysis of the closed-loop system as presented in section 5.4.2 leads to an open question related to the impact of solving the optimization problem by using the Lagrangian duality on the stability of the closed-loop system. Nevertheless, the derived optimization problem being a convex one, leading to better convergence properties than the original min-max problem, an interesting idea is to focus the analysis on the impact of the convergence of the optimization algorithm on the stability of the closed-loop system.

• On the other side, an interesting perspective may be the determination of sufficient conditions ensuring robust stability of the overall hierarchical control scheme, presented in Chapter 6, in case of bounded uncertainties and/or trajectory constraints, starting from the work already done in [START_REF] Rugabotti | Robust model predictive control of continuous-time sampled data nonlinear systems with interval sliding mode[END_REF].

• It would be also interesting to analyse the robust stability of the closed loop system (LRMPC strategy coupled with the observer,e.g., Extended Kalman Filter) as done in [START_REF] Huang | Robust stability of nonlinear model predictive control with extended Kalman filter and target setting[END_REF].

• Haldane kinetics [START_REF] Andrews | A mathematical model for continuous culture of microorganism utilizing inhibitory substrate[END_REF]:

µ = µ m S S + K s + S 2 K i (B.5)
where K i is the inhibition constant (g.L -1 ). The growth is inhibited by the substrate S.

• Cantois kinetics [START_REF] Contois | Kinetics of bacterial growth: Relationship between population density and specific growth rate of continuous cultures[END_REF]:

µ = µ m S S + XK c (B.6)
where K c is the Contois saturation constant . The growth depends upon the concentrations of both substrate S and biomass X with growth being inhibited at high concentrations of biomass.

These kinetics are also used to model substrate uptake.

In the case of microalgae, the above laws that describe reaction kinetics are used to model both the limitation and the inhibition by the light intensity.

Remark C.1. It should be mentionned that there exist several controllability approaches such as controllability, local controllability, weak controllability and local weak controllability. For further details we refer the reader to [START_REF] Hermann | Nonlinear controllability and observability[END_REF].

Hereafter, condition (C.3) is checked for the studied system. After developments, we have that

ad fx f u (x) =       0 K s μρ m (S -S in ) (K s + S) 2 - K s ρ m (s -s in )X (K s + S) 2       (C.

C.2 Dynamics of sensitivity functions

This part details the development of the dynamics of sensitivity functions for the Droop model (section 7.2). The dynamics of sensitivities are calculated based on analytical derivation [START_REF] Dochain | Automatic control of bioprocesses[END_REF]

=                 0 0 1 - K Q Q I I+K sI + I 2 K iI -μ Q I I+K sI + I 2 K iI -μ 1 - K Q Q I (I+K sI + I 2 K iI ) 2 μ 1 - K Q Q I 3 K 2 iI (I+K sI + I 2 K iI ) 2 0                 , ∂ρ(S) ∂θ =             S S+Ks 1 -Q Q l -ρ m 1 -Q Q l S (S+Ks) 2 0 0 0 0 ρ m SQ S+Ks 1 Q 2 l             and ∂µ(Q, I) ∂x =       0 μ I I + K sI + I 2 K iI K Q Q 2 0       , ∂ρ(S) ∂x =       0 - ρ m Q l S S + K s ρ m 1 - Q Q l K s (S + K s ) 2      

C.3 Sensitivity analysis

Hereafter, a sensitivity analysis is carried out to determine the most influencial parameters. To simplify the study of the sensitivity functions, each function is rescaled by multiplying it by the parameter value under study. The state x i is said to be sensitive to a given parameter if a change in the parameter's value significantly affects the predictive quality of the model. The sensitivity of one state with respect to all parameters is to be considered separately and the different orders of magnitude have to be compared. For the considered system (i.e. 3 state variables and 6 parameters), 3 × 6 sensitivity functions are then computed for biomass concentration X = 25µm 3 L -1 > X l (satisfying (7.35)). The corresponding dilution rate is calculated at the equilibrium using (7.33). The sensitivity analysis according to the magnitude order of the sensitivity functions allows to select the parameters which are significantly the most influential.

According to results obtained in on the biomass concentration evolution (100 times higher than the other parameters).

Table C.1: The ranking of parameters according to their influence on the model (from more to less).

X Q S K Q K Q K Q K s K s K s μ μ ρ m K sI K sI μ K iI K iI K sI ρ m ρ m K iI
Thereafter, only these two parameters (K s ,K Q ) are considered in the minmax optimization problem, instead of the 6 model parameters (the 4 other parameters are set to their nominal values).

C.4 Additional simulation results

In this section, the NMPC-ε s/m and Generic Model Control (GMC) performances for reference tracking are compared in the case of model mismatch and measurement noise.

C.4.1 Predictive controller

The formulation of the NMPC-ε s/m optimization problem at each sampling instant kT s is as follows: where λ is the control weighting factor ensuring that D remains close to its reference D r , β is a penalty weighting factor on the control variation term and This penalization term was chosen to allow a good trade-off between noise rejection and trajectory tracking. The GMC tuning parameters are chosen as follows: ξ =1 and ω 0 =5 rad d -1 . For the tuning parameters ξ =1 and ω 0 =5 rad d -1 , it can be observed (Figure C.2) that the peak time is well respected. The GMC law leads to a larger overshoot in biomass concentration evolution during a setpoint change, with a longer transient and response time in comparison with the NMPC performances (mainly due to the lack of anticipation). In the case of the NMPC, it can be noticed the anticipation of a setpoint change due to the prediction of the model behavior in the moving horizon. More details on this example can be found in [START_REF] Benattia | Nonlinear model predictive control for regulation of microalgae culture in a continuous photobioreactor[END_REF]. Abstract : The last few years have led to very rapid developments, both in the formulation and the application of Nonlinear Model Predictive Control (NMPC) algorithms, with a relatively wide range of industrial achievements. One of the most significant challenges encountered during the development of this control law is due to uncertainties in the model of the system. In this context, the thesis addresses the design of NMPC control laws robust towards model uncertainties. Usually, the above design can be achieved through solving a min-max optimization problem. In this case, the idea is to minimize the tracking error for the worst possible uncertainty realization. However, this robust approach tends to become too complex to be solved numerically online, especially in the case of multivariable systems with a large number of uncertain parameters. To address this shortfall, the main proposed approach consists in simplifying the min-max optimization problem through a sensitivity analysis of the model with respect to its parameters, in order to reduce the calculation time. First, the criterion is linearized around the model parameters nominal values. The optimization variables are either the system control inputs or the control increments over the prediction horizon. The initial optimization problem is then converted either into a convex optimization problem, or a onedimensional minimization problem, depending on the nature of the constraints on the states and commands. The stability analysis of the closedloop system is also addressed. Finally, a hierarchical control strategy is developed, that combines a robust model predictive control law with an integral sliding mode controller, in order to cancel any tracking error. The proposed approaches are applied through two case studies to the control of microorganisms culture in bioreactors.

5. 1

 1 Comparison of (C)LRMPC algorithms in terms of computation time at each sampling time (mismatched parameters). . . 5.2 Comparison of the proposed algorithms in terms of computation time at each sampling time (mismatched parameters). . . 5.3 Comparison of the proposed algorithms in terms of computation time at each sampling time (nominal parameters). . . . . 6.1 Comparison of LRMPC and LRMPC-δu algorithms in term of parameters of the optimization problem. . . . . . . . . . . . 6.2 Controllers tuning parameters. . . . . . . . . . . . . . . . . . . 7.1 Droop model parameters. . . . . . . . . . . . . . . . . . . . . . 7.2 Simulation conditions for the Droop model. . . . . . . . . . . . 7.3 Comparison of the predictive algorithms in terms of computation time at each sampling time. . . . . . . . . . . . . . . . . . 7.4 PI & ISM tuning parameters. . . . . . . . . . . . . . . . . . . 7.5 Comparison of the proposed algorithms in terms of tracking error distribution features. . . . . . . . . . . . . . . . . . . . . C.1 The ranking of parameters according to their influence on the model (from more to less). . . . . . . . . . . . . . . .

Notation 3 .Notation 4 .Notation 5 .Notation 7 .Notation 8 .- 1

 345781 The notation A * denotes the conjugate transpose of the matrix A. The notation A † denotes the pseudo inverse of the matrix A such thatA † lim δ→0 (A * A + δI) -1 A * Given n ∈ Z ≥0 ,an arbitrary norm of a vector x ∈ R n is denoted as |x|. Notation 6. ||z|| 2 P = z P z is the Euclidean norm weighted by the matrix P . Matrix norm ||A|| is given by ||A|| = σ(A * A) with σ(A) the maximum eigenvalue of A. The signum function of a real number x (denoted sign) is defined as follows: sign(x) :=   if x < 0 0 if x = 0 1 if x > 0

Definition 6 .

 6 An optimization problem of the form min x f (x) s.t. g i (x) ≤ 0, i = 1, . . . , m

Figure 1 . 1 -

 11 Figure 1.1 -Schéma de la structure générale de la commande prédictive.

Figure 1 . 4 -

 14 Figure 1.4 -Architecture de la commande prédictive distribuée [34].

Figure 1 . 5 -

 15 Figure 1.5 -CPNL incluant le signal ε s/m .

Étape 2 . 1 k

 21 57) avec W (λ) donnée par (1.50). La séquence de commande u k+Npest calculée à partir du problème d'optimisation (1.53) pour λ = λ :

Figure 1 . 6 -

 16 Figure 1.6 -Schéma de la stratégie de commande hiérarchisée.

)Figure 1 . 9 -

 19 Figure 1.9 -Évolution temporelle de l'erreur de poursuite pour une variation aléatoire non-corrélée de tous les paramètres (Monte-Carlo).

Figure 3 . 1 :

 31 Figure 3.1: Structure of the Model Predictive Control strategy.

Figure 3 . 2 :

 32 Figure 3.2: Receding horizon principle: the basic idea (SISO case).

2 Figure 3 . 3 :

 233 Figure 3.3: Centralized MPC architecture [34].

2 Figure 3 . 4 :

 234 Figure 3.4: Decentralized MPC architecture [34].

1 Figure 3 . 5 :

 135 Figure 3.5: Distributed MPC architecture [34].

Figure 3 . 6 :

 36 Figure 3.6: Classification of dynamic optimization solution methods [18].

Figure 4 . 1 :

 41 Figure 4.1: Discrete state path (SISO case).

Figure 4 . 2 :

 42 Figure 4.2: Principle of the NMPC (SISO case).

  values, where the subscript is related to the time instant.

Figure 4 . 3 :

 43 Figure 4.3: Diagram of the NMPC algorithm.

Figure 4 . 5 :

 45 Figure 4.5: Schematic representation of a bioreactor.

Figure 4 . 6 :

 46 Figure 4.6: Influence of the sampling time T s on the behavior of the biomass concentration for NMPC law.

Figure 4 . 7 :

 47 Figure 4.7: Influence of the sampling time T s on the behavior of the substrate concentration for NMPC law.

Figure 4 . 8 :

 48 Figure 4.8: Influence of the sampling time T s on the behavior of the control input for NMPC law.

Figure 4 . 9 :

 49 Figure 4.9: Influence of the prediction horizon N p on the behavior of the biomass concentration for NMPC law.

Figure 4 . 10 :

 410 Figure 4.10: Influence of the prediction horizon N p on the behavior of the substrate concentration for NMPC law.

Figure 4 . 11 :Table 4 . 2 :

 41142 Figure 4.11: Influence of the prediction horizon N p on the behavior of the control input for NMPC law.

Figure 4 . 12 :

 412 Figure 4.12: Influence of the weighting matrices on the biomass concentration evolution with time for NMPC law.

Figure 4 . 13 :

 413 Figure 4.13: Influence of the weighting matrices on the substrate concentration evolution with time for NMPC law.

Figure 4 . 14 :

 414 Figure 4.14: Influence of the weighting matrices on the control input evolution with time for NMPC law.

Figure 4 . 16 :Figure 4 . 17 :

 416417 Figure 4.16: Substrate concentration evolution with time for NMPC-(j)ε s/m laws in the case of parameter uncertainties.

Figure 5 . 1 :

 51 Figure 5.1: Diagram of RNMPC based on min-max optimization problem.

Figure 5 . 2 :

 52 Figure 5.2: Predictions considered for the stability proof.

step 1 .||Az -b|| 2 W•

 12 Solution of the OLS problemLet us consider first the following ordinary least squares (OLS) problem: min z A is a known m × n matrix;

77 ) step 2 .

 772 Introduction of uncertainties in the RLS problem (5.75) The cost function (5.76) and the optimization problem (5.75) can be modified to take into account uncertainties on A and b as follows: min z max δA,δb

2 V 2 W( 5 . 84 ) 3 .

 225843 -||Azb + Cξ|| step Inclusion of the previous constraint (5.82) in the robust RLS problem (5.84)

step 4 .

 4 Solving the problem(5.86)

88 )

 88 It turns out that the dual variable λ must satisfy the following inequality constraint λ ≥ ||C W C|| (5.89)Thanks to the(5.87) and (5.89), problem (5

Figure 5 . 3 :

 53 Figure 5.3: Diagram of the LRMPC algorithm.

k

  by solving the QP problem (5.118) for λ = λ ; 5.

Figure 5 . 4 :

 54 Figure 5.4: Diagram of the CLRMPC algorithm.

Figure 5 . 5 :Figure 5 . 6 :Figure 5 . 7 :Figure 5 . 8 :Figure 5 . 9 :Figure 5 . 10 :

 5556575859510 Figure 5.5: Biomass concentration evolution with time for LRMPC strategy for several values of the sampling time T s .

Figure 5 . 11 :Figure 5 . 12 :

 511512 Figure 5.11: Tracking error evolution with time for LRMPC strategy for several values of the prediction horizon N p .

Figures 5 .dFigure 5 . 13 :

 5513 Figures 5.13-5.15 compares both the constrained and unconstrained LRMPC (as presented in sections 5.5.2 and 5.6.2) performances obtained with the following tuning parameters: N p = 6, T s = 5 min, T d = 0.1 min, V = W = I Np .

Figure 5 . 14 :Figure 5 . 15 :

 514515 Figure 5.14: Substrate concentration evolution with time for (C)LRMPC strategies.

Remark 5 . 8 .

 58 The same phenomenon can be observed for a setpoint change in falling edge as shown in Figure5.15.

Figure 5 . 17 :

 517 Figure 5.17: Substrate concentration evolution with time for NMPC, RN-MPC and LRMPC strategies.

Figure 5 . 18 :

 518 Figure 5.18: Control input evolution with time for NMPC, RNMPC and LRMPC strategies.

Figure 6 . 3 :

 63 Figure 6.3: Control input evolution with time for LRMPC and LRMPC-δu strategies (nominal case).

Figure 6 . 5 :Figure 6 . 6 :

 6566 Figure 6.5: Substrate concentration evolution with time for LRMPC and LRMPC-δu strategies (mismatched parameters).

Figure 6 . 7 :

 67 Figure 6.7: Scheme of the hierarchical control strategy [131].

Figure 6 . 8 :

 68 Figure 6.8: Scheme of the hierarchical control strategy based on robust MPC and PI controller.

Figure 6 . 9 :

 69 Figure 6.9: Scheme of the hierarchical control strategy based on robust MPC and ISM controller.

Figure 6 . 10 :

 610 Figure 6.10: Biomass concentration evolution with time for LRMPC, LRMPC-PI and LRMPC-ISM strategies.

Figure 6 . 11 :Figure 6 . 12 :

 611612 Figure 6.11: Substrate concentration evolution with time for LRMPC, LRMPC-PI and LRMPC-ISM strategies.

•

  D represents the dilution rate (d -1 , d: day) which is the ratio of the inlet flow rate to the volume of the culture.

  e , I) -D e = 0 ρ(S e )µ(Q e , I)Q e = 0 (S in -S e )D eρ(S e )X e = 0 (7.15)

Figure 7 . 2 :Figure 7 . 3 :

 7273 Figure 7.2: Graphical representation of the condition (7.28).

)Figure 7 . 4 :

 74 Figure 7.4: Biomass concentration evolution with time for NMPC, RNMPC, rRNMPC and LRMPC strategies.

)Figure 7 . 8 :

 78 Figure 7.8: Biomass concentration evolution with time for LRMPC-(δu) strategies.

Figure 7 . 11 :

 711 Figure 7.11: Substrate concentration evolution with time for LRMPC-(δu) strategies.

Figure 7 . 12 :

 712 Figure 7.12: Temporal evolution of biomass concentration for LRMPC, LRMPC-PI and LRMPC-ISM strategies.

Figure 7 . 16 :

 716 Figure 7.16: Output and control input evolution with time for LRMPC-PI strategy.

Figure 7 . 17 :

 717 Figure 7.17: Time evolution of tracking error for simultaneous random non correlated variation in all the parameters (Monte-Carlo).

Figure 7 .

 7 Figure 7.18 shows the resulting histogram of the distribution of the tracking error, and Table7.5 gives the mean and standard deviation of tracking error for the whole simulation.

EFigure 7 . 18 :

 718 Figure 7.18: Histogram of the average tracking error for simultaneous random non correlated variation in all the parameters (Monte-Carlo).

Figure 7 . 19 :

 719 Figure 7.19: Biomass concentration evolution with time for NMPC-jε s/m , RNMPC, LRMPC-PI and LRMPC-ISM strategies in case of time varying reference trajectory.

)Figure 7 . 23 :

 723 Figure 7.23: Biomass concentration evolution with time for LRMPC, LRMPC-PI and LRMPC-ISM strategies (disturbance rejection).

Figure 7 . 25 :

 725 Figure 7.25: Light intensity evolution with time for LRMPC, LRMPC-PI and LRMPC-ISM strategies (disturbance rejection).

- 2 s 4 + 4 +

 244 Ks μK Q ρmX(S-S in ) Q 2 (Ks+S) 2 K μ2 X(S-S in )-Ksρ 2 m SX(Ks-S+2S in ) (Ks+S) Ks μρm(S-S in ) (Ks+S) 2 Ksρ 2 m SX 2 (Ks-S+2S in )-Ks 2 ρ 2 m X 2 (S-S in ) (Ks+S)Ks determinent of the controllability matrix C o is given by:det(C o ) = K 2 s μK Q ρ 2 m X(S -S in ) 2 (S -S in + QX) Q 2 (K s + S) 4 (C.7)Thus, we obtaindetC o = 0 iff S -S in = 0 S -S in + QX = 0 (C.8)Based on Theorem C.1, the affine-system (C.1)-(C.2) is locally controllable if the condition (C.8) holds. These conditions are verified for the considered simulation tests.

Figure C. 1 :

 1 Figure C.1: Evolution of the scaled sensitivity functions (for biomass concentration).

  Hg(t k , t k+j , x k , D k+j-1 k , θ nom ), j = 1, N p D k ≥ 0, X k ≥ 0, ∀k ∈ N (C.12)

Figure C. 2 :

 2 Figure C.2: Evolution of biomass concentration and dilution rate with measurement affected by noise, GMC and NMPC-ε s/m laws.

Université

  Paris-Saclay Espace Technologique / Immeuble Discovery Route de l'Orme aux Merisiers RD 128 / 91190 Saint-Aubin, France Titre : Robustification de la commande prédictive non linéaire -Application à des procédés pour le développement durable Mots clés : Commande prédictive non linéaire, Robustesse, Systèmes incertains, Bioprocédé Résumé : Les dernières années ont permis des développements très rapides, tant au niveau de l'élaboration que de l'application, d'algorithmes de commande prédictive non linéaire (CPNL), avec une gamme relativement large de réalisations industrielles. Un des obstacles les plus significatifs rencontré lors du développement de cette commande est lié aux incertitudes sur le modèle du système. Dans ce contexte, l'objectif principal de cette thèse est la conception de lois de commande prédictives non linéaires robustes visà-vis des incertitudes sur le modèle. Classiquement, cette synthèse peut s'obtenir via la résolution d'un problème d'optimisation min-max. L'idée est alors de minimiser l'erreur de suivi de la trajectoire optimale pour la pire réalisation d'incertitudes possible. Cependant, cette formulation de la commande prédictive robuste induit une complexité qui peut être élevée ainsi qu'une charge de calcul importante, notamment dans le cas de systèmes multivariables, avec un nombre de paramètres incertains élevé. Pour y remédier, la principale approche proposée dans ces travaux consiste à simplifier le problème d'optimisation min-max, via l'analyse de sensibilité du modèle vis-à-vis de ses paramètres afin d'en réduire le temps de calcul. Dans un premier temps, le critère est linéarisé autour des valeurs nominales des paramètres du modèle. Les variables d'optimisation sont soit les commandes du système soit l'incrément de commande sur l'horizon temporel. Le problème d'optimisation initial est alors transformé soit en un problème convexe, soit en un problème de minimisation unidimensionnel, en fonction des contraintes imposées sur les états et les commandes. Une analyse de la stabilité du système en boucle fermée est également proposée. En dernier lieu, une structure de commande hiérarchisée combinant la commande prédictive robuste linéarisée et une commande par mode glissant intégral est développée afin d'éliminer toute erreur statique en suivi de trajectoire de référence. L'ensemble des stratégies proposées est appliqué à deux cas d'études de commande de bioréacteurs de culture de microorganismes. Title : Robustification of Nonlinear Model Predictive Control -Application to sustainable development processes Keywords : Nonlinear model predictive control, Robustness, Uncertain Systems, Bioprocess

  Architectures de commande prédictive . . . . . . . . . 1.1.5 Robustesse . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.6 Stabilité . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.6.1 Stabilité nominale . . . . . . . . . . . . . . . 1.1.7 Stabilité robuste . . . . . . . . . . . . . . . . . . . . . 1.1.8 Avantages et inconvénients . . . . . . . . . . . . . . . . 1.2 Chapitre 4 : Commande prédictive non linéaire . . . . . . . . . 1.2.1 Formulation du problème d'optimisation . . . . . . . . 1.2.2 Une variante de la CPNL . . . . . . . . . . . . . . . . 1.3 Chapitre 5 : Commande prédictive non linéaire robuste . . . . 1.3.1 Stratégie min-max . . . . . . . . . . . . . . . . . . . . 1.3.2 Commande prédictive robuste réduite . . . . . . . . . . 1.3.2.1 Analyse de sensibilité . . . . . . . . . . . . . 1.3.2.2 Reformulation du problème . . . . . . . . . . 1.3.3 Commande prédictive robuste linéarisée . . . . . . . . 1.3.3.1 Principe général . . . . . . . . . . . . . . . . 1.3.3.2 Analyse de la stabilité . . . . . . . . . . . . . 1.3.4 Implémentation de la commande prédictive robuste linéarisée en absence de contraintes (LRMPC) . . . . . . . . . . Commande prédictive robuste linéarisée avec prise en compte des contraintes . . . . . . . . . . . . . . . . . . 1.4 Chapitre 6 : Améliorations de la commande prédictive robuste linéarisée . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4.1 Commande prédictive robuste linéarisée avec incréments de commande . . . . . . . . . . . . . . . . . . . . . . . 1.4.2 Stratégie de commande hiérarchisée . . . . . . . . . . . 1.5 Chapitre 7 : Application à un système de culture de microalgues 1.5.1 Modélisation du système . . . . . . . . . . . . . . . . . 1.5.2 Stratégie de commande . . . . . . . . . . . . . . . . . . 1.5.3 Résultats en simulation . . . . . . . . . . . . . . . . . . 1.6 Conclusions et perspectives . . . . . . . . . . . . . . . . . . . . Literature review . . . . . . . . . . . . . . . . . . . . . . . . . 3.4.1 Linear case . . . . . . . . . . . . . . . . . . . . . . . . 3.4.2 Nonlinear case . . . . . . . . . . . . . . . . . . . . . . . 3.4.3 Predictive Control architecture . . . . . . . . . . . . . 3.5 Methods for dynamic optimization . . . . . . . . . . . . . . . 3.6 Robust NMPC schemes . . . . . . . . . . . . . . . . . . . . . . 3.7 Stability . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.7.1 Nominal stability . . . . . . . . . . . . . . . . . . . . . 3.7.2 Robust stability . . . . . . . . . . . . . . . . . . . . . . 3.8 Advantages and drawbacks . . . . . . . . . . . . . . . . . . . . 3.9 Concluding remarks . . . . . . . . . . . . . . . . . . . . . . . . x 7.23 Biomass concentration evolution with time for LRMPC, LRMPC-PI and LRMPC-ISM strategies (disturbance rejection). . . . . 173 7.24 Control input evolution with time for LRMPC, LRMPC-PI and LRMPC-ISM strategies (disturbance rejection). . . . . . . 174 7.25 Light intensity evolution with time for LRMPC, LRMPC-PI and LRMPC-ISM strategies (disturbance rejection). . . . . . . 174 7.26 Growth rate and dilution rate evolution for LRMPC-PI strategy (disturbance rejection). . . . . . . . . . . . . . . . . . . . 175 C.1 Evolution of the scaled sensitivity functions (for biomass concentration). . . . . . . . . . . . . . . . . . . . . . . . . . . . . 193 C.2 Evolution of biomass concentration and dilution rate with measurement affected by noise, GMC and NMPC-ε s/m laws. . 196 xix List of Tables 1.1 Paramètres du modèle. . . . . . . . . . . . . . . . . . . . . . . 1.2 Conditions de simulation. . . . . . . . . . . . . . . . . . . . . 1.3 Comparaison des algorithmes prédictifs en termes de temps de calcul à chaque instant d'échantillonnage. . . . . . . . . . . . . 1.4 Comparaison des lois de commande en termes de distribution de l'erreur de poursuite. . . . . . . . . . . . . . . . . . . . . .

	CONTENTS
	1.3.5
	ix

viii Contents 1 Résumé 1.1 Chapitre 3 : Commande prédictive -état de l'art et principales stratégies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.1 Contexte . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.2 Principes de la commande prédictive . . . . . . . . . . 1.1.2.1 Modèle de prédiction . . . . . . . . . . . . . . 1.1.2.2 Fonction de coût . . . . . . . . . . . . . . . . 1.1.2.3 Calcul de la loi de commande . . . . . . . . . 1.1.3 État de l'art . . . . . . . . . . . . . . . . . . . . . . . . 1.1.3.1 Cas linéaire . . . . . . . . . . . . . . . . . . . 1.1.3.2 Cas non linéaire . . . . . . . . . . . . . . . . 1.1.4 1.6.1 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . 1.6.2 Perspectives . . . . . . . . . . . . . . . . . . . . . . . . 2 Introduction 2.1 Context and motivations . . . . . . . . . . . . . . . . . . . . . 2.2 Outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.3 List of Publications . . . . . . . . . . . . . . . . . . . . . . . . 3 Predictive control: State of the art and main strategies 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.3 The basic principles of MPC . . . . . . . . . . . . . . . . . . . 3.3.1 Prediction model . . . . . . . . . . . . . . . . . . . . . 3.3.2 Cost function . . . . . . . . . . . . . . . . . . . . . . . 3.3.3 Control law calculation . . . . . . . . . . . . . . . . . . 3.4 4.1 Model parameters for system (4.15). . . . . . . . . . . . . . . . 4.2 Influence of the prediction horizon N p on the computation time.

  . . . . . Nonlinear Model Predictive Control with ε s/m signal NMPC-jε s/m Nonlinear Model Predictive Control with jε s/m signal ODE

	Acronyms
		NLP	NonLinear Programming
		NMPC	Nonlinear Model Predictive Control
	asNMPC NMPC-ε s/m advanced step Nonlinear Model Predictive Control
	BMI		Bilinear Matrix Inequality
	CLRMPC	Constrained Linearized Robust Model Predictive Controller Ordinary Differential Equation
	CSTR	OLS	Continuous Stirred Tank Reactor Ordinary Least Squares
	CVP	PDE	Control Vector Parametrization Partial Differential Equation
	DMC	PFC	Dynamic Matrix Control Predictive Functional Control
	DMS	PI	Direct Multiple Shooting Proportional Integral
	DPC	PVC	Distributed Predictive Control PolyVinyl Chloride
	DSS	Direct Simple Shooting QIH-NMPC Quasi-Infinite Horizon Nonlinear Model Predictive Control
	D-RTO	QP	Dynamic Real Time Optimization Quadratic Programming
	EHAC	RHC	Extended Horizon Adaptive Control Receding Horizon Control
	EPSAC	RLS	Extended Prediction Self Adaptive Control Robust Least Squares
	E-NMPC RNMPC Economic Nonlinear Model Predictive Control Robust Nonlinear Model Predictive Control
	GMC	Generic Model Control rRNMPC reduced Robust Nonlinear Model Predictive Control
	GPC	RTO	Generalized Predictive Control Real Time Optimization
	IMC	SISO	Internal Model Control Single Input Single output
	ISM	SMC	Integral Sliding Mode Sliding Mode Control
	ISS	SMPC	Input-to-State Stability Stochastic Model Predictive Control
	LMI	SQP	Linear Matrix Inequality Sequential Quadratic Programming
	LP	UAV	Linear Programming Unmanned Aerial Vehicule
	LPV	ZOH	Linear Parameter Varying Zero Order Hold
	LRMPC	Linearized Robust Model Predictive Controller
	LRMPC-δu	incremental Linearized Robust Model Predictive Control
	LRMPC-ISM Linearized Robust Model Predictive Control-Integral Sliding Mode
	LRMPC-PI	Linearized Robust Model Predictive Control-Proportional Integral
	MAC		Model Algorithmic Control
	MIMO		Multiple Input Multiple output
	MPC		Model Predictive Control
			xxi xxiii	xxiv

  3) où N p est l'horizon de prédiction, xk+i est l'état prédit, x k est la condition initiale correspondant à l'état du système à l'instant k et φ(., .) est le coût instantané (généralement un terme quadratique). La commande prédictive non linéaire est implémentée selon le principe de l'horizon fuyant. Á l'instant k, l'état du système x k est utilisé comme condition initiale et le problème d'optimisation (1.2)-(1.3) est résolu sur un horizon de prédiction de longueur N p . Néanmoins, seule la première valeur de la séquence de commande est appliquée au système, i.e. u k . La procédure est répétée à la période d'échantillonnage suivante. La structure de base globale d'une boucle de régulation prédictive est représentée sur la figure 1.1.

	Fonction de coût				entrée de commande	Système	état du système
	Optimiseur	futures entrées	prédiction Modèle de	prédites sorties	consigne de réfèrence
			future errors
	Commande prédictive

  de classe C 2 vis-à-vis de tous ces arguments, représente la dynamique du système non linéaire. La matrice d'observation est donnée par H ∈ R ny×nx . Les sorties discrétisées sont obtenues à chaque période d'échantillonnage T s constante par l'intégration du modèle d'état continu (1.5) en utilisant par exemple la méthode de Runge-Kutta avec un pas d'intégration T d inférieur Résumé à T s . Ainsi le modèle de prédiction est défini par les équations récursives suivantes :

  8) où x k+1 est l'état à l'instant t k+1 , k est l'indice de temps, x k et y k sont le vecteur d'état et des mesures discrétisées à l'instant t k , respectivement. L'entrée de commande u(t) est paramétrée en utilisant une approximation constante par morceaux sur un intervalle de temps [t k , t k+1 ] [kT s , (k + 1)T s )] :

  1.3 Chapitre 5 : Commande prédictive non linéaire robuste Ce chapitre s'attache à proposer une nouvelle stratégie de commande prédictive robuste. Elle est fondée sur l'approximation du problème min-max via la linéarisation du modèle de prédiction et l'utilisation de la dualité Lagrangienne. Dans ce contexte, deux méthodes ont été développées. Dans la première méthode, la séquence de commande optimale est calculée à partir d'un problème d'optimisation scalaire. La seconde méthode permet de prendre en considération explicitement les contraintes de type inégalité.

1.3.1 Stratégie min-max

On suppose que le vecteur des paramètres est incertain et appartient à l'ensemble Θ défini par (1.6)-(1.7). Ainsi, la robustification de la commande prédictive non linéaire induit la formulation d'un problème d'optimisation de type min-max

[START_REF] Allwright | Advances in Model-Based Predictive Control, chapter on min-max Model-Based Predictive Control[END_REF][START_REF] Magni | Assessment and Future Directions of Nonlinear Model Predictive Control, chapter Robustness and Robust Design of MPC for Nonlinear Discrete-Time Systems[END_REF][START_REF] Michalska | Robust receding horizon control of constrained nonlinear systems[END_REF]

. La séquence de commande qui minimise la fonction coût est la solution du problème d'optimisation suivant :

1.22 1.23 1.24 1.25 1.26 1.27

														Résumé
	35 mande prédictive non linéaire permettant la poursuite d'une trajectoire 25.3 850
	30 de référence pour un système non linéaire à temps discret. 25.25 800	NMPC LRMPC LRMPC-PI I nominal I real
	25 • Dans la deuxième partie, la commande prédictive robuste vis-à-vis des 25.2 X r LRMPC LRMPC-PI 750 incertitudes paramétriques est formulée en un problème d'optimisa-
	LRMPC-ISM tion de type min-max. Afin de réduire la charge de calcul induite ) Number of scores 15 20 X(µm 3 .L -1 ) 25.1 .s -1 700 par la formulation min-max, l'approche proposée consiste à simpli-25.15 650 fier le problème d'optimisation initial, via l'analyse de sensibilité du modèle vis-à-vis de ses paramètres en linéarisant les sorties prédites I (µE.m -2 autour des valeurs nominales des paramètres du modèle et de la sé-
	600 quence de commande de référence. Deux méthodes d'implantation ont 10 25.05 été développées. Dans la première méthode ne tenant pas compte des
	contraintes, la séquence de commande optimale est calculée à partir 550 5 25 d'un problème d'optimisation scalaire. La seconde méthode permet de
	500 prendre en compte explicitement les contraintes inégalité sur l'entrée
	-0.6 0 0.1 de commande. 0 24.95	0	0.2	-0.4 0.3 0.1	0.2	0.4	Time (d) E[X -X r ] (µm 3 .L -1 ) -0.2 0 0.5 0.6 0.7 Time (d) 0.3 0.4 0.5 0.6	0.2 0.8 0.7	0.9 0.8	0.4	1	0.9	0.6 1
	Figure 1.10 -Histogramme de la moyenne de l'erreur de poursuite pour Figure 1.11 -Évolution temporelle de la concentration en biomasse. Figure 1.13 -Perturbation sur l'intensité lumineuse.
		une variation aléatoire non-corrélée de tous les paramètres (Monte-Carlo). La figure 1.11 montre que les trois lois de commande arrivent à maintenir la 1.4 concentration en biomasse constante en présence de perturbations sur l'in-
		On constate que la stratégie hiérarchisée LRMPC-PI permet d'avoir la pour-suite la plus précise en termes d'écart type sur l'erreur de poursuite comme tensité lumineuse. Dans le cas de la LRMPC, on constate la présence d'une 1.3 importante erreur statique due à l'étape de linéarisation durant la prédiction.
	) D (d -1	indiqué dans le tableau 1.4. En revanche, la stratégie LRMPC-PI arrive à maintenir la sortie autour de 1.2 la valeur de consigne grâce à l'action du régulateur PI. En ce qui concerne Tableau 1.4 -Comparaison des lois de commande en termes de distribution la seconde commande hiérarchisée (LRMPC-ISM), on note la présence d'une de l'erreur de poursuite. 1.1 certaine erreur statique non négligeable due à l'évaluation du taux de crois-sance aux valeurs nominales des paramètres dans l'expression de la loi de 1 commande par mode glissant. En effet, ce dernier point constitue l'inconvé-
		nient le plus important de l'utilisation de la commande par mode glissant.
		0.9
		0.8 1.6 Conclusions et perspectives	LRMPC D r
			0.3			0.4	0.5	0.6	0.7			LRMPC-PI
													LRMPC-ISM
		0.7 1.6.1 Conclusions
		0	0.1	0.2	0.3	0.4	0.5	0.6	0.7	0.8		0.9	1
		Time (d) Dans cette thèse, des algorithmes de commande prédictive robustes pour des
		systèmes non linéaires à temps discret soumis à des incertitudes paramé-
		Figure 1.12 -Évolution temporelle du taux de dilution. triques ont été proposés. Les principales contributions de cette thèse peuvent
		être résumées comme suit :
		• La première partie propose la mise en oeuvre d'une stratégie de com-

  Tebbani and D. Dumur. Robust nonlinear model predictive control of microalgae culture in a continuous photobioreactor. 1st IFAC Conference on Modelling identification and Control of Nonlinear Systems, Saint-Petersburg, Russia , pp. 192-197, June 24-26, 2015.

• S.E. Benattia, S.

Tebbani

and D. Dumur. Multivariable robust MPC of heterotrophic microalgae fed-batch bioreactor. To be submitted to Journal of Process Control, 2016. Proceedings • S.E. Benattia, S. Tebbani and D. Dumur. Nonlinear model predictive control for regulation of microalgae culture in a continuous photobioreactor. Proceedings of the 22nd Mediterranean Conference on Control and Automation, Palermo, Italy, pp. 469-474, June 16-19, 2014. • S.E. Benattia, S. Tebbani, D. Dumur and D. Selisteanu. Robust nonlinear model predictive controller based on sensitivity analysis -application to a continuous photobioreactor. Multi-Conference on Systems and Control, Antibes/Nice, France, pp. 1705-1710, October 8-10, 2014. • S.E. Benattia, S. • S.E. Benattia, S. Tebbani and D. Dumur. Hierarchical control strategy based on robust MPC and integral sliding mode -application to a continuous photobioreactor. 5th IFAC Conference on Nonlinear Model Predictive Control, Seville, Spain, pp. 212-217, September 17-20, 2015. • S.E. Benattia, S. Tebbani and D. Dumur. Hierarchical predictive control strategy of microalgae culture in a photobioreactor. 19th International Conference on System Theory, Control and Computing, Cheile Gradistei-Fundata, Romania , pp. 231-236, October 14-16, 2015. • S.E. Benattia, S. Tebbani and D. Dumur. Robust model predictive control -application to microalgae cultivation. Submitted to Conference on Decision and Control, Las Vegas, USA, 2016. Conferences with abstract • S.E. Benattia, S. Tebbani, D. Dumur. Nonlinear model predictive control of microalgae culture. Young Algaeneers Symposium, Montpellier & Narbonne, France, April 3-5, 2014. • S.E. Benattia, S. Tebbani, D. Dumur. Regulation of microalgae culture in a continuous photobioreactor. ESBES-IFIBIOP Symposium, Lille, France, September 7-10, 2014. Others (oral presentation) • S.E. Benattia, S. Tebbani and D. Dumur. Robust nonlinear model predictive control for regulation of microalgae culture in a continuous photobioreactor. Journée du Groupe de Travail de Commande Prédictive Non Linéaire, Paris, June 11, 2015. • S.E. Benattia, S. Tebbani and D. Dumur. Hierarchical robust predictive control strategy of microalgae cultivation process. Journée du Groupe de Travail de Commande Prédictive Non Linéaire, Paris, March 24, 2016.

Table 4 .

 4 1: Model parameters for system(4.15). The simulation time is set to T f = 20 hours. The reference trajectory X r has two changes: a rising and falling edge respectively at two different instants as shown in Figure4.6. The reference dilution rate D r nom and the reference substrate concentration S r nom are determined at the equilibrium for the considered values of X r and model parameters of Table4.1, as follows (from (4.15) and (4.16)):

	Parameter Value Unit
	µ m	0.3	h -1
	K m	1.75 g.L -1
	Y	0.9	-

Table 5 .

 5 1: Comparison of (C)LRMPC algorithms in terms of computation time at each sampling time (mismatched parameters).Computation time

Table 5 . 2

 52

: Comparison of the proposed algorithms in terms of computation time at each sampling time (mismatched parameters).

Computation time

Table 5 .

 5 3: Comparison of the proposed algorithms in terms of computation time at each sampling time (nominal parameters).

		Computation time (s)
	X Algo. X X X CLRMPC X X X X X X Perf. indices min mean X X 0.2 0.47	max 0.62
	RNMPC	0.04 1.55	2.92

Table 6 .

 6

	1: Comparison of LRMPC and LRMPC-δu algorithms in term of
	parameters of the optimization problem.
	LRMPC	LRMPC-δu
	z

Table 6 .

 6

		2: Controllers tuning parameters.
		LRMPC	PI	ISM
	Parameter N p V	W K p T i ξ 1	ξ 2 ξ 3 K s
	Value	6 I Np I Np 0.1 0.1 0.1 -0.1 2 15

Table 7 .

 7

		1: Droop model parameters.
	Parameter Value	Unit
	μ	2	d -1

Table 7 .

 7 2: Simulation conditions for the Droop model.

		Variable Value Unit
	sampling time	T s	10	min
	integration time step	T d	12	sec
	simulation time	T f	1	d
	inlet substrate concentration	S in	100	µmol L -1
	optimal light intensity	I opt	547	µE m -2 s -1
	maximum cell quota	Q l	9	µmol µm -3
	maximal admissible dilution rate D max	1.6	d -1
	prediction horizon	N p	5	-
	weighting matrix on control	V	I Np	-
	weighting matrix on state	W	I Np	-
	initial biomass concentration	X(0)	24.95 µm 3 L -1
	initial internal quota	Q(0)	4	µmol µm -3
	initial substrate concentration	S(0)	0.05	µmol L -1

  Biomass concentration measurements, y k , are assumed to be corrupted by a centred Gaussian white noise with 0.1 standard deviation. Thanks to the sensitivity analysis of the model with respect to its parameters (see Appendix C), only the most influential parameters κ are considered in the min-max problem instead of the full model parameters as presented in section 5.3 (rRNMPC). In this case, the substrate half saturation constant K s and the minimal cell quota K Q are the most influential parameters on the biomass concentration evolution (i.e. κ = [K s , K Q ] ). The other parameters, ζ = [ρ m , μ, K sI , K iI ] , are set to their nominal values (Table7.1).

						5): a classi-
	cal Nonlinear Model Predictive Control (4.12) (denoted as NMPC), a robust
	NMPC using criterion (5.3) (denoted as RNMPC), a reduced robust one
	using criterion (5.11) (denoted as rRNMPC) and the linearized one (5.104)-
	0 (5.105) (LRMPC). Time [d] 0.1 0.2 0.3 0.4 0.5 X(µm 3 .L -1	0.6	0.7	0.8	0.9	1

Table 7 .

 7

	3: Comparison of the predictive algorithms in terms of computation
	time at each sampling time.
		Computation time (s)
	X Algo. X X X NMPC X X X X X Perf. indices X X X	min < 10 -5 0.014 0.29 mean max
	RNMPC	0.55	1.85 19.31
	rRNMPC	0.047 0.078 0.36
	LRMPC	< 10 -5 0.01	0.09

  Figure 7.9: Dilution rate evolution with time for LRMPC-(δu) strategies.

		0.08		1.6
		0.07		1.4
		0.06		1.2
	)	0.05		1
	(µmol.L -1	0.04	D (d -1 )	0.8
	S	0.03		0.6				1.4
								1.35
		0.02		0.4				1.3
								1.25
		0.01		0.2				1.2				S r	D r
												LRMPC	LRMPC
		0	0	0.1 0	0	0.1 0.2	0.2 0.3	0.3 Time [d] 0.4 0.4 0.5 0.38	Time [d] 0.5 0.6 0.7 0.4	0.6 0.42 0.8	0.7 LRMPC-δu 0.8 0.9	0.9 LRMPC-δu 1	1
				4.05								Q r
													LRMPC
				4								LRMPC-δu
				3.95
			(µmol.µm -3 )	3.85 3.9
			Q	3.8
				3.75
				3.7
				3.65
					0	0.1	0.2	0.3	0.4	0.5	0.6	0.7	0.8	0.9	1
										Time [d]
	Figure 7.10: Internal quota evolution with time for LRMPC-(δu) strategies.

Table 7 .

 7 4: PI & ISM tuning parameters.

	Parameter K p	T i	ξ 1 ξ 2 ξ 3 K s
	Value	1 0.01 1 -1 0.1 1
	The obtained results as depicted in Figures 7.12-7.15 compare the LRMPC
	and hierarchical controller (LRMPC-PI and LRMPC-ISM) performances con-
	sidering the uncertain parameter worst case cited previously with biomass
	affected by noise measurement.

Table 7 . 5 :

 75 Comparison of the proposed algorithms in terms of tracking error distribution features.In section 4.4.2, simulations have shown the NMPC performances with jε s/m signal for setpoint reference trajectory. In the sequel, a time-varying reference trajectory is selected as shown in Figure7.19 instead of step-setpoint trajectory. The purpose of this choice is to investigate the NMPC-jε s/m (section 4.3) performances in the case of variable reference trajectory to highlight the advantage of the proposed strategy (both LRMPC-PI and LRMPC-ISM). Four controllers will be tested: a classical Nonlinear Model Predictive Control with addition of the error signal jε s/m during prediction step (denoted as NMPC-jε s/m ), a robust one using criterion (5.1) (denoted as RNMPC) and two hierarchical control strategies (LRMPC-PI and LRMPC-ISM). Biomass concentration measurements are assumed to be non-corrupted by noise. The light intensity is assumed to be measured and equal to I opt .

	X X Algo. X X X Perf. indices mean standard deviation X X X X X X X NMPC -0.07 0.219
	LRMPC	-0.06	0.098
	LRMPC-PI	-0.05	0.014

.16 25.62 25.64 25.66 25.68 25.7

  Figure 7.21: Internal quota evolution with time for NMPC-jε s/m , RNMPC, LRMPC-PI and LRMPC-ISM strategies in case of time varying reference trajectory.

				4.05									NMPC-jε s/m
													RNMPC
													LRMPC-PI
				4									LRMPC-ISM Q r
				3.95
			)
		1.4 1.6	Q (µmol.µm -3	3.9
				3.85
		1.2
				3.8
	)
	D (d -1	1		3.75
					0	0.1	0.2	0.3	0.4	0.5	0.6	0.7	0.8	0.9	1
		0.8							Time [d]
		0.6										D r NMPC-jε s/m
												RNMPC
												LRMPC-PI
												LRMPC-ISM
		0.4
		0	0.1	0.08	0.2	0.3	0.4	0.5	0.6	0.7	0.8	0.9	1
								Time [d]
	Figure 7.20: Dilution rate evolution with time for NMPC-jε s/m , RNMPC, 0.07
	LRMPC-PI and LRMPC-ISM strategies in case of time varying reference trajectory. 0.06
			)
			S (µmol.L -1	0.04 0.05
				0.03
													NMPC-jε s/m
				0.02									RNMPC LRMPC-PI
													LRMPC-ISM
													S r
				0.01
				0	0.1	0.2	0.3	0.4	0.5	0.6	0.7	0.8	0.9	1
									Time [d]
	Figure 7.22: Substrate concentration evolution with time for NMPC-jε s/m ,
	RNMPC, LRMPC-PI and LRMPC-ISM strategies in case of time varying
	reference trajectory.

  3 /L). The goal is thus to maintain the biomass concentration as close as possible to this value, despite the light intensity fluctuations. Figures 7.23-7.26 compare the LRMPC, LRMPC-PI and LRMPC-ISM controllers performances. The light intensity considered in the model for prediction is I nominal .

	.L -1
	X(µm 3
	0	0.1	0.2	0.3	0.4	0.5	0.6	0.7	0.8	0.9	1
						Time (d)

  as follows:

	It follows that
	                      	d dt d dt d dt	∂X ∂θ j ∂Q ∂θ j ∂S ∂θ j	= = = -∂µ ∂θ j ∂ρ(S) X + (µ(Q, I) -D) ∂θ j -∂µ(Q, I) ∂θ j Q -µ(Q, I) ∂X ∂θ j + ∂µ(Q, I) ∂Q ∂Q ∂θ j + ∂ρ(S) ∂Q X ∂θ j ∂S ∂S ∂θ j ∂S ∂θ j D -∂ρ(S) ∂θ j X -ρ(S) ∂X ∂θ j -∂ρ(S) ∂S ∂S ∂θ j X	-	∂µ(Q, I) ∂Q	∂Q ∂θ j	Q
													(C.10)
	with
	∂µ(Q, I)
		∂θ
				d dt	∂x i ∂θ j	=	∂ ∂θ j	dx i dt	=	∂f i ∂θ j	+	nx k=1	∂f i ∂x k	∂x k ∂θ j	(C.9)

X X X X X X X X X X X XAlgo. Ind. perf. moyenne écart type NMPC -0,07 0,219 LRMPC -0,06 0,098 LRMPC-PI -0,05 0,014 Ainsi, à l'aide de cet outil statistique on peut conclure que la structure hiérarchisée permet d'annuler l'erreur de poursuite due à l'étape de linéarisation dans la partie synthèse prédictive. En dernier lieu, on compare la commande prédictive robuste linéarisée (LRMPC) avec deux stratégies hiérarchisées (LRMPC-PI et LRMPC-ISM) en absence de bruit de mesure. L'objectif de commande est de maintenir la concentration en biomasse à une valeur de consigne fixe (X r = 25µm 3 /L) en dépit des perturbations sur l'intensité lumineuse (figure 1.13).

Remerciements

The structure of the NMPC law with consideration of the model-system error signal is given in Although this strategy remains simple to implement and improves performances, even better results will be obtained in the following chapters by considering robust NMPC structures.

Numerical illustrative example

In the sequel, to illustrate the strategies presented previously, the control of a bacteria growth in a bioreactor is considered, in order to evaluate the performance of the setpoint tracking NMPC controller. Schematic representation of a bioreactor is shown in Figure 4.5. Generally, components are introduced into the bioreactor at the rate q in and others are removed from the bioreactor at the rate q out . Classically, three operating modes are used in practice:

• Batch mode: components are neither introduced nor removed i.e.

q in = q out = 0.

• Fed-batch mode: only components are introduced (q in > 0) and q out = 0.

• Continuous mode: components are introduced and removed at the same rate i.e. q in = q out > 0.

Let us consider a biological process in continuous mode with a single biomass (bacteria) and a single substrate (nutrient). The dynamics of biomass and substrate concentrations, denoted X and S respectively, are obtained using mass balances for a continuous stirred tank reactor (CSTR) [START_REF] Bastin | On-line estimation and adaptive control of bioreactors[END_REF]. More details concerning the mass balance modelling can be found in Appendix B.

with

Different approaches are possible for determining the sensitivity functions with respect to the parameters vector and the control sequence, defined in (5.20) and (5.21), respectively as presented in section 5.3.1.

The dynamics of the sensitivity function with respect to θ can be computed for time t ∈ [t k , t k+Np ] by solving numerically the following differential equation (from (5.6)):

with as an initial condition:

Indeed, ∇ θ g is the sensitivity function of the state x with respect to the parameter θ.

An alternative procedure is to use the finite differences in order to approximate numerically the derivatives ∇ θ g for each parameter θ l , l ∈ [1, n θ ] and ∇ u g(t k+j ) for each control

The finite difference method approximates the (i, j)-th element of the jacobian of a vector function g(z) as

for some small δ > 0. A too large δ will induce inaccuracies due to the nonlinearity of g i , since the method computes the average slope between two points.

Remark 5.3. The most accurate result and computationally most efficient approach is to calculate gradients analytically (by symbolic differentiation). Doing this by hand, or even using symbolic computations in Maple or Mathematica, may quickly become intractable for MPC problems that may contain a large number of variables and parameters.

Details concerning how to obtain (5.91) and (5.92) can be found in Appendix A.

The optimization problem (5.90) is further replaced by: min

where the cost function J(z, λ) is defined as follows

The optimal value of z for every fixed value of λ is determined by cancelling the derivative of J(z, λ) with respect to z (i.e. ∇ z J(z, λ) = 0). Consequently, the minimum z must satisfy the following equation

), the gradient of Γ(z) 2 with respect to z is given by:

From (5.96), the solution of the equation (5.95), which dependents on λ, is given by:

with

The modified weighting matrix V (λ) is obtained from V via:

The invertibility of E(λ) is guaranteed by the positive definiteness of V . By replacing z(λ) by its optimal value given by (5.97) in the right side of (5.93), the nonnegative scalar parameter λ ∈ R solution of (5.93) is determined from the following unidimensional minimization problem:

(5.101) Finally, the problem has a unique global minimum z given by (5.97) for λ = λ , solution of (5.101) (i.e. z = z(λ )). More details related to this development can be found in [START_REF] Sayed | A regularized robust design criterion for uncertain data[END_REF].

In the following simulations, only case 2 will be considered, starting with the knowledge of the reference biomass concentration X r . In one side, it can be observed that Q * 1 is always non-negative (Figure 7.1) and lower than Q l (condition 7.2 holds) for values of X e > X l ≈ 12.5 µm 3 L -1 as depicted in Figure 7.2 (green box). On the other side, it appears clearly that Q * 2 is a non-admissible solution since it does not respects the condition (7.29) as shown in Figure 7.3.

Then, to have a feasible solution (here only one, which is Q * 1 ), the reference biomass concentration X r should satisfy the following constraint:

In our case, if the reference is constant, X r must be chosen greater than X l ≈ 12.5 µm 3 L -1 . In the following, during constant tracking of a setpoint, X r will be chosen between 25 and 27 µm 3 L -1 .

Appendix A

A.1 Robust regularized Least squares problem

This part details the development from (5.85) to (5.91) in section 5.5.

Let us consider the Lagrangian given by (5.85):

The above equation (A.1) can be expressed as:

The optimal solution ξ that minimizes L is given by:

2) and after some mathematical manipulationsn, it comes:

Then, we have that:

Finally, the Lagrangian (A.1) becomes a function of z and λ as follows:

Then, (A.6) becomes:

which is the expression (5.91).

Appendix B B.1 Biological systems modelling

This section summarizes general principles in bioprocess modelling. More details can be found in [START_REF] Bastin | On-line estimation and adaptive control of bioreactors[END_REF]. Biological systems are described by a set of macroscopic components such as biomass, substrates and metabolic products. The mechanistic approach is the most popular bioprocess modelling technique. For optimization, monitoring and control purposes, this approach is usually macroscopic in essence, i.e., it makes use of the concept of macroscopic reaction scheme involving a few reactants, products and catalysts considered as macroscopic entities. The material flow set that defines the transformation of reactants into products is represented in the following form [START_REF] Bastin | On-line estimation and adaptive control of bioreactors[END_REF]:

where • Π .,k is the k th pseudo stoichiometric coefficient or yields coefficients. They are negative when they relate to a reactant and positive when they relate to a product.

• r k (.) is the k th reaction rate.

• The constants n r , n s and n p are respectively the number of reactions, substrates and products.

• S i and P j represent respectively the i th component consumed and j th component produced.

Mass balance modelling

Assuming that the reactor is homogeneous, the macroscopic mass balances for each component involved in the reaction network (B.1) lead to a general differential state-space model [START_REF] Bastin | On-line estimation and adaptive control of bioreactors[END_REF]:

where

• r(.) is the reaction rate vector, (made of r k (.)).

• K the matrix of yields coefficients, (made of Π .,k ).

• Λ the vector of components concentrations.

• Λ in the vector of feed concentrations.

• D the dilution rate which is defined as the medium feed rate over the effective reactor volume.

• G(Λ) describes the exchange between the gas and the liquid phases.

B.2 Reaction kinetics modelling

The reaction rates describes how the compenents interact in the reactor based on activation, limitation and inhibition phenomena. The expression of the reaction rates is typically a nonlinear function depending on the components concentrations and some kinetic parameters. As an example, for the reactions associated with a biomass X, the growth reaction is expressed as follows [START_REF] Bastin | On-line estimation and adaptive control of bioreactors[END_REF]:

where µ k denotes the specific growth rate, which is used to describe the bacterial growth. In the literature, commonly used laws, that describe the reaction kinetics, are determined by the biologists. The most cited one are presented below:

• Monod model [START_REF] Monod | Recherches sur la croissance des cultures bactériennes[END_REF]:

where µ m is the maximum specifi growth rate (h -1 ) and K s is the half saturation constant (g.L -1 ). It models the limitation of the growth by the substrate S.

Appendix C

In this section, additional results for the Droop model studied in Chapter 7 are proposed. First, the controllability of model is discussed. Secondly, the dynamics of the sensitivity functions is detailed. Finally, the NMPC-ε s/m performances are compared to those obtained with the GMC law.

C.1 Controllability

Controllability is concerned with whether one can design control input to steer the state to arbitrary values. Then, in the next, the controllability of the studied system is checked.

The nonlinear system (7.5) is expressed as a control-affine system as follows:

with

Theorem C.1. The nonlinear system (C.1) is locally controllable at x if and only if the following Lie algebra rank condition is satisfied:

Proof. see [START_REF] Isidori | Nonlinear Control Systems: An Introduction[END_REF].

Finally, it comes:

Xk+j|k is the predicted biomass concentration. Model uncertainties are taken into account by assuming that the gap between the system and the model at sample time k, denoted ε s/m , is included in the prediction as detailed in section 4.3.

C.4.2 Generic model control

The GMC is a kind of Nonlinear Proportional Integral regulator. It includes proportional and integral actions, cancelling static errors, and performs an input/output linearization of the system. Details related to the GMC law can be found in [START_REF] Lee | Generic model control (GMC)[END_REF]. In the case of the considered control problem, the dilution rate derived from the GMC law is given by:

The GMC tuning parameters (G 1 = 2ξω 0 and G 2 = ω 2 0 ) are chosen based on the fact that a second order behavior is specified for the closed-loop system. The damping factor ξ gives a desired shape of response, and the natural frequency ω 0 an appropriate settling time.

C.4.3 Simulation results

The initial biomass concentration value is set close to the setpoint in order to cancel the transient effect and to focus only on the behavior during setpoint changes. The simulation time is fixed at 3 d, the sample time T s is chosen equal to 30 min. Only uncertainties on the first four parameters of vector θ are considered. In this case, the four first parameters of the mismatched model are chosen on the parameter subspace border (θ real = [θ + 1 , θ - 2 , θ + 3 , θ - 4 ]), where the uncertain parameters subspace [θ -, θ + ] is given by [ θnom 2 , 3θnom 2 ]. The two last parameters K sI and K iI are set to their nominal values. Biomass concentration measurements, y k , are assumed to be corrupted by a Gaussian white noise of zero mean and 0.05 standard deviation. The light intensity is assumed to be measured online, non-corrupted with noise.

The tuning parameters of the NMPC are determined by a trial-and-error technique. The prediction horizon of the NMPC law is set to N p =10, chosen to satisfy a compromise between the computation time and a sufficient vision of the system behaviour in the future. The weighting factor λ is selected at 1 to promote the tracking error on biomass concentration relatively to the tracking error on dilution rate in the objective function. In addition, the penalty weighting factor on the control variation β is chosen equal to 10λ.