Learning-based model predictive control for constrained nonlinear systems
Abstract
Esta tesis está dedicada al control de sistemas, bajo la hipótesis de que no se conoce nada sobre la dinámica del sistema a controlar. En vez de eso, solamente las entradas y salidas son accesibles, y por tanto se puede tener acceso a un histórico de datos. El objetivo principal es el control de la planta en condiciones eficientes y seguras usando únicamente dichas medidas. Con este fin, se usará un conjunto de métodos de aprendizaje automático conocido como link inferencia, para modelar sistemas no lineales desconocidos, usando controladores predictivos basados en modelo. Por ello, esta tesis presenta contribuciones en dos campos distintos. En primer lugar, se extenderán las técnicas de link inferencia, proponiendo métodos con dos objetivos: reducir tanto el tiempo de cálculo de los algoritmos como el error de predicción cometido por ellos. En segundo lugar, se desarrollarán controladores predictivos robustos y con la habilidad de aprender basándose en datos. Estos controladores serán estables por diseño, capaces de satisfacer restricciones robustamente y de mejorar su actuación beneficiándose de nuevas medidas recogidas en línea.
Full text
DOCTORAL THESIS Learning-based model predictive control for constrained nonlinear systems Jos´e Mar´ıa Manzano Crespo Seville, 2020 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
ESCUELA T´ ECNICA SUPERIOR DE INGENIER´ IA DEPARTAMENTO DE INGENIER´ IA DE SISTEMAS Y AUTOM´ ATICA Doctoral Thesis Learning-based model predictive control for constrained nonlinear systems Jos´e Mar´ıa Manzano Crespo Thesis submitted for the degree of Doctor of Philosophy in the University of Seville Supervised by: Daniel Lim´on Marruedo David Mu˜noz de la Pe˜na Sequedo Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
A mi familia, mis amigos y mis profesores I Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
Acknowledgements Agradecimientos A Daniel Lim´on, por todo lo que he aprendido contigo estos a˜nos, de lo cual lo menos imporante que me llevo es todo lo que me has ense˜nado de control predictivo. Espero que el fin de esta tesis no suponga el fin de nuestros almuerzos, charlas y por supuesto del trabajo juntos. Muchas gracias, Dani. A David Mu˜noz, por transmitirme su punto de vista y eficacia, que atesoro como ´utiles cualidades de lo que requiere ser buen investigador y docente. Muchas gracias por tantas horas de dedicaci´on a todos los resultados de esta tesis. Huge thanks to Jan Calliess, not only for receiving me in Oxford during my research stay, but also for all the work done together, and all the teachings and enlightening discussions that arose these years. A Teo ´ Alamo, por todas las ideas y debates de pizarra que han conseguido que mejore mi visi´on cr´ıtica de la investigaci´on. A los dem´as profesores del Grupo de Estimaci´on, Predicci´on, Optimizaci´on y Control: Daniel Rodr´ıguez, Ignacio Alvarado y Mario Pereira, por todos los consejos y ayuda prestada. A Jos´e Ram´on Salvador y Pablo Krupa, por todos los a˜nos compartiendo despacho y etapa de doctorando, as´ı como a los dem´as miembros del grupo y de despacho, Juan, Joaqu´ın, Daniel, V´ıctor y Jose, y dem´as profesores de la Escuela T´ecnica Superior de Ingenier´ıa. Muchas gracias a Bernardo P´erez, por hacer que me sintiera en Oxford como en casa. A todos mis profesores durante todas las etapas de mi educaci´on, como Ana Montero, que me ense˜n´o el primer problema de matem´aticas o Antonio Almansa, que hizo que las amase. Especialmente gracias a Manolo Fern´andez Arias: sin tus ense˜nanzas no estar´ıa escribiendo hoy esta tesis. II Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
Acknowledgements Muchas gracias a toda mi familia y amigos, no s´olo por vuestro apoyo durante estos a˜nos, sino tambi´en por vuestra contribuci´on en mi educaci´on. Todos hab´eis aportado para convertirme en quien soy hoy, y por tanto este trabajo es en parte vuestro; espero que os haga sentir orgullosos. Finalmente, me gustar´ıa dedicar esta tesis a Marina, por su paciencia y ayuda los a˜nos que ha durado la elaboraci´on. Sevilla, 2020. III Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
Abstract This thesis stands upon the assumption that nothing is known about the dynamics of the system to be controlled. Instead, only inputs and outputs are measurable, and hence an historical set of these measurements may be available. The main objective is to control this system under efficient and safe conditions using only these measurements. To this aim, the class of learning rules known as kinky inference is used to model nonlinear output-feedback systems, in a model predictive control framework. Hence, the contributions of this thesis are two-fold. First, the kinky inference methods are extended, proposing methods both to reduce the computational burden and to decrease the prediction errors. Second, robust learning-based predictive controllers are derived. These controllers are stable by design, and are able to ensure robust constraint satisfaction and to improve their performance by benefiting from new information gathered online. IV Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
Resumen Esta tesis est´a dedicada al control de sistemas, bajo la hip´otesis de que no se conoce nada sobre la din´amica del sistema a controlar. En vez de eso, solamente las entradas y salidas son accesibles, y por tanto se puede tener acceso a un hist´orico de datos. El objetivo principal es el control de la planta en condiciones eficientes y seguras usando ´unicamente dichas medidas. Con este fin, se usar´a un conjunto de m´etodos de aprendizaje autom´atico conocido como kinky inference, para modelar sistemas no lineales desconocidos, usando controladores predictivos basados en modelo. Por ello, esta tesis presenta contribuciones en dos campos distintos. En primer lugar, se extender´an las t´ecnicas de kinky inference, proponiendo m´etodos con dos objetivos: reducir tanto el tiempo de c´alculo de los algoritmos como el error de predicci´on cometido por ellos. En segundo lugar, se desarrollar´an controladores predictivos robustos y con la habilidad de aprender bas´andose en datos. Estos controladores ser´an estables por dise˜no, capaces de satisfacer restricciones robustamente y de mejorar su actuaci´on benefici´andose de nuevas medidas recogidas en l´ınea. V Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
Contents 1 Introduction 1 1.1 Motivation and objectives . . . . . . . . . . . . . . . . . . . . 1 1.2 Context .............................. 4 1.2.1 Identification and machine learning . . . . . . . . . . . 5 1.2.2 Efficient and safe data-driven control . . . . . . . . . . 7 1.3 Data-based models . . . . . . . . . . . . . . . . . . . . . . . . 8 1.3.1 Nonlinear autoregressive models . . . . . . . . . . . . . 8 1.3.2 Machine learning . . . . . . . . . . . . . . . . . . . . . 10 1.4 Kinkyinference.......................... 13 1.5 Model predictive control . . . . . . . . . . . . . . . . . . . . . 17 1.5.1 Learning-based MPC . . . . . . . . . . . . . . . . . . . 20 1.6 Mainresults............................ 22 1.7 Thesisoverview.......................... 23 1.8 Publications............................ 25 2 Output feedback MPC based on smoothed projected kinky inference 27 2.1 Introduction............................ 27 2.2 Smoothed projected kinky inference . . . . . . . . . . . . . . . 28 2.2.1 Smoothed kinky inference . . . . . . . . . . . . . . . . 28 2.2.2 Projected kinky inference . . . . . . . . . . . . . . . . 29 2.3 Data-based predictive control . . . . . . . . . . . . . . . . . . 33 VI Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
CONTENTS 2.3.1 Stability analysis of the controller . . . . . . . . . . . . 35 2.4 Casestudy............................. 42 3 Robust learning-based MPC for nonlinear constrained systems 54 3.1 Introduction............................ 54 3.2 Problemsetting.......................... 55 3.3 Stabilizing data-based NMPC . . . . . . . . . . . . . . . . . . 56 3.4 Casestudy............................. 64 4 Componentwise H¨older inference for robust learning-based MPC 68 4.1 Introduction............................ 68 4.2 Extension of H¨older properties . . . . . . . . . . . . . . . . . . 69 4.2.1 Componentwise H¨older continuity . . . . . . . . . . . . 69 4.2.2 Componentwise H¨older kinky inference . . . . . . . . . 72 4.3 Reachability sets . . . . . . . . . . . . . . . . . . . . . . . . . 80 4.4 Stabilizing conditions for CHoKI-MPC . . . . . . . . . . . . . 83 4.5 Casestudies............................ 84 4.5.1 Quadruple tank process . . . . . . . . . . . . . . . . . 85 4.5.2 Isothermal chemical reactor . . . . . . . . . . . . . . . 91 5 Online learning constrained MPC based on double prediction 93 5.1 Introduction............................ 93 5.2 Double prediction framework . . . . . . . . . . . . . . . . . . . 94 5.3 Online learning-based MPC . . . . . . . . . . . . . . . . . . . 95 5.4 Casestudy.............................108 6 Online learning robust MPC: an exploration-exploitation approach 113 6.1 Introduction............................113 VII Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.1. MOTIVATION AND OBJECTIVES The process of learning, when referred to machines instead of humans, reminds of the artificial intelligence (AI) field. By “artificial intelligence” we describe all machines able to make decision or to carry out tasks in a “smart” way, like beating chess players or recognizing speech or patterns in images. Among all the available technologies within the AI field, the machine learning (ML) sub-field refers to those algorithms that are able to learn (and to improve this knowledge) as the exposure to data increases, without being explicitly programmed to do so. The development of computer capacity has favoured the continuous appearance of more powerful AI and ML applications. In this thesis, ML techniques will be used to learn the behaviour of dynamical systems and to compute predictions in different scenarios. The second objective of this thesis characterises the control itself: what kind of decision-making process do we want to apply? Once the response of the system to a given input is modelled, a typical closed-loop controller penalises the deviation of the output to the desired reference, so control inputs that steer the system to the given reference are applied. Considering more advanced control techniques, we aim to design controllers that orbit around two crucial aspects: efficiency and safety. By “efficient”, we address the objective of implementing controllers that provide the best control actions, according to some performance index that measures what “best” is. With respect to safety, we account for the guarantee that the closed-loop system is not going to violate the security limits, by means of admissible control actions, while asymptotically converging to its reference. Safety is a key issue in most control applications. For instance, consider the automatic control (autopilot) of an aircraft. This is considered a critical system, in the sense that if something fails, the potential damage is devastating. The same applies to the control systems of drones or large chemical plants. Achieving safety and stability by design of the controlled systems is one of the main challenges of this thesis, due to the black-box nature of the plants to be controlled. In addition, note that because of the dynamical nature of the systems in consideration and their constraints, control algorithms require the decision making process to be solved in real time. For this reason, it is of critical importance to develop computationally efficient implementations of them. Summing up, the main objective of this thesis is the design of efficient and safe control strategies for constrained nonlinear systems, that only use both historical data and real-time measurements in the decision making process. A schematic representation of the framework is represented in Figure 1.1. To achieve these objectives, the model predictive control (MPC) frame- 2 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.1. MOTIVATION AND OBJECTIVES Plant MPC D Model uˆ y Cost Constraints u∗˜y Figure 1.1: Scheme of the overall framework. work will be used, since it is an advanced control technique that, based on a prediction model, minimizes a performance index in the operation, while being able to guarantee safety by taking the constraints into account. In this thesis, a learning-based approach is taken, in which machine learning techniques that offer accurate predictions will be tailored to identify the unknown dynamics in a computationally tractable form. It has been proven that in their standard form, predictive controllers can ensure robust stability of the closed-loop if (among other requirements) the model’s transition function is Lipschitz continuous and if the uncertainty is bounded with a known worst-case value. Because of this relevant property, the machine learning method considered in this thesis is a class of Lipschitz interpolation techniques known as kinky inference, whose resulting predictor is Lipschitz continuous, and it provides guaranteed uncertainty bounds and an accurate estimation of the Lipschitz constant, which will be used for the controller design. Consequently, this thesis contributes to both the fields of machine learning and model predictive control. Regarding the first, several improvements to the class of kinky inference methods are proposed, aimed to reduce computational times and to improve the quality of the predictions, decreasing their estimation error. Second, new formulations of robust learning-based model predictive controllers are developed, able to ensure stability in the presence of both hard and soft constraints in inputs and outputs, as well as to include new measurements of the plant in an online fashion, thus learning and improving during the operation. This chapter presents an introductory analysis of the work proposed in 3 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.2. CONTEXT this thesis, its context in the existing literature, its motivation and objectives. Section 1.2 presents the context and the state-of-the-art. The problem setting is addressed in Section 1.3. Section 1.4 introduces the kinky inference technique, and Section 1.5 the model predictive control framework. The main contributions of the thesis are sketched in Section 1.6. Finally, Section 1.7 offers an overview of the contents of each chapter and appendices, and Section 1.8 lists the publications resulting from their elaboration. 1.2 Context In the light of growing successes of machine learning algorithms, data-based control techniques have become increasingly popular in the control community in recent years. Although the conditions under which a model of the system can be obtained from input-output data have been studied since the late 1980’s [41, 76], the recent development of advanced machine learning techniques has favoured the appearance of data-based controllers, using different approaches, like dual [8], adaptive [17], or reinforcement learning [150] controllers. Another control methodology in which ML finds application consists in deriving data-based models of the system, and to use model predictive controllers based on those, yielding the learning-based MPC field [13]. In contrast to the standard identification approach [86], these new techniques offer a wider class of functions they are able to learn, and have proven to be useful when learning dynamical systems, e.g. [73]. In these cases, a machine learning technique can be used to learn and to predict the evolution of the plant. Among the different techniques applied to learning-based MPC, one of the most popular consists in using a machine learning algorithm to obtain a predictor from input-output data, such as direct weight optimization [137, 136], Gaussian processes [52, 94] or neural networks [159], among many others. Such techniques may handle data in a deterministic or stochastic way, e.g., Lipschitz interpolation [40] for the former and Gaussian processes [2, 18] for the latter. These approaches have been applied in different model predictive controllers [63]. One of the most interesting properties of using ML techniques for control is its natural capability to learn from fresh data, that is, to enhance its predictions by adding new information, like the one collected online from feedback. Control techniques that incorporate this ability are referred to as online learning controllers. While being a popular topic in reinforcement 4 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.2. CONTEXT learning (RL) [78, 115] and adaptive control [1, 153], the matter of online learning still seems to be somewhat under-researched in data-based model predictive control [43, 121]. The different results found in the literature can be classified according to the use they make of new data, like improving an initial feasible solution in an optimization problem [138], designing safety filters [156], or improving the model of the system [62, 87]. In conclusion, there exists many data-driven controllers based on machine learning. However, when the systems in consideration are complex and have constraints, these techniques must be able able to guarantee adequate properties over the closed-loop, like safety of operation and stability. This is currently an open and challenging field, in which this thesis is centered. Further motivation and state-of-the-art of such fields is presented next. 1.2.1 Identification and machine learning One of the main objectives in this thesis is the study and design of novel data-based methods tailored to real-time decision making processes. Over the last two decades, different techniques such as statistical learning [53], pattern recognition [30], and kernel-based methods [140] have been applied to different fields, like process monitoring [128], fault detection [65] and nonlinear functions approximation [151]. In statistical learning, in order to manage the trade-off between model complexity and performance, regularization (or penalization) parameters to enhance robustness of the results are included in the learning process [58]. Another source of uncertainty in machine learning methods is the lack of sufficient historical data [6]. Techniques to minimize the amount of data required have been studied, using probabilistic methods based on synthetic generation of scenarios [33]. Besides, the machine learning field is not only based on statistics, but also on optimization. This is the case, for example, of random and kernel methods [5]. Specific optimization techniques must be used to handle the resulting high dimensional problems, often with non differentiable objective functions. There exists a manifold of optimization techniques used in machine learning, like accelerated methods, proximal methods, dual decomposition or coordinate descent [147]. With the appropriate data analysis methods, different inference problems in nonlinear systems have been addressed, like state estimation, prediction and reachability analysis. State estimation is a key aspect in many fields 5 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.2. CONTEXT related to dynamical systems. In this context, there exists many strategies in the machine learning field, like parametric [126] or non-parametric methods [29], or multivariate statistical distributions. Recently, seeking to handle uncertainty, nonlinearities and lack of data, new works have satisfactorily applied the tools of Bayesian inference [50], which benefit from the ever increasing computational power available. This thesis will focus on combining classic identification scopes with databased inference techniques, aiming to improve the results available in the literature. Existing inference methods will be extended to state estimation, based on the information provide by historical data. In order to design processes for optimal decision making in real-time, it is necessary to obtain guaranteed predictions of the future behaviour of the plant. The reachability analysis of a dynamical system usually leads to highly complex problems that require different approximation methods, like those based on ellipsoids, zonotopes [3], or interval algebra [28]. One of the main challenges on reachability analysis is the scaling to higher dimensional problems [48], as well as to combine available methods with techniques based on data to bound the reachable sets of the systems in a more precise and efficient way. This problem will be studied in this thesis. Beyond data-based inference, for an efficient real-time operation in the presence of uncertainty, it is necessary to asses the quality of the models used online, being able to adapt such models and to quantify the uncertainty [2]. Data validation and online model adaptation (both structural and parametric) is a key aspect to guarantee safe operation of a system [138]. Therefore, the derivation of data-based methods able to learn and to adapt themselves online, recursively, as more data become available during the operation will also be studied in this thesis. Finally, it is worth mentioning methods such as data completion [152], that have been applied to various fields and which use stored data to find connecting patterns and to recover lost information. Similarly, data-based virtual sensors [148] usually offer a method to obtain estimated values of physical variables that cannot be measured by a physical sensor. All these technologies motivate the study of black-box systems carried out in this thesis. 6 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.2. CONTEXT 1.2.2 Efficient and safe data-driven control Efficiency is a term that measures the behaviour of a given system, taking into account, for example, energy consumption, use of actuators, evolution according to a prescribed setting, etc. This term can encompass a collection of measurements, and in a control engineering terminology, this is comprised in a general function named performance index or economic cost. The controller must choose the action that each time step provides the least economic cost, associated to the dynamical evolution of the system [44]. The decision-making process must also ensure that the system, during its evolution, guarantees the imposed operation limits, derived for example from environmental conditions or physical constraints [60]. This problem becomes much more complex if the inherent uncertainty of the unknown plant, as well as the modelling errors due to the learning methods applied are taken into account [82]. Designing control laws to optimize the economic cost of the closed-loop operation of the plant is challenging, due to the general metric used to measure efficiency, which does not necessarily penalize the tracking error w.r.t. a set-point, as it is customary in classic control approaches. Recently, model predictive control techniques have been proposed to deal with this problem, leading to the so-called economic MPC. Following this line, stable designs of economic predictive controllers for large processes have been proposed in the last decade [60]. One of the most important properties of model predictive control [38] is that the resulting controller may, by design, guarantee stability of the closed-loop system, which together with the guarantee of constraints satisfaction, yields safe operation of the system. Model predictive control has a strong dependence on the quality of the predictions. In order to improve closed-loop performance, predictive controllers able to adapt their models online have been developed, both linear [13] and nonlinear [40]. Currently, there is ongoing research on this trending topic, as shown for example in the excellent survey paper by Hewing et al. [63]. In contrast to the majority of approaches presented in that review, the work proposed in this thesis will address the problem from a deterministic point of view, instead of stochastically, in order to give ensured guarantees of robust performance, which is a crucial aspect in critical systems. 7 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.3. DATA-BASED MODELS 1.3 Data-based models This section explains how a regressive model of a system can be derived using only input-output data. Then, the machine learning field is introduced to provide a framework able to learn such regressive models. 1.3.1 Nonlinear autoregressive models Consider a nonlinear, discrete time, output-feedback, multiple inputs multiple outputs (MIMO) dynamical system in which y(k)∈Rnyis the vector of outputs, x(k)∈Rnxthe state, and u(k)∈Rnuthe manipulable inputs at time step k∈N0. Note that by system, one typically refers to the real plant to be controlled. However, it is often not possible to have a perfect knowledge of the real system behaviour, so an approximated mathematical model is used instead. In order to design a controller, the first objective of every control engineer is to obtain such model of the system, to simulate the plant’s behaviour, and thus taking into account the possible response of the closed-loop system in the design. A model is even more important in predictive control strategies, where the model is explicitly used by the controller, in order to predict future evolution of the plant. A first attempt to construct a model may be based on first principles. This approach consists in describing the set of physical, ordinary differential equations that govern the plant. However, as it was stated before, in this thesis we assume that the inner signals of the system are not accessible by the user, and only input-output samples are available. For this reason, a regressive form of the state must be constructed in order to describe the system effectively. For the case of linear dynamics, autoregressive exogenous (ARX) models provide an appropriate way to accurately describe the behaviour from inputoutput information, such that y(k+ 1) + a1y(k) + . . . +anay(k−na) = b1u(k) + . . . +bnbu(k−nb).(1.1) Here, the model is said to be autoregressive since y(k+1) depends on previous outputs, up to the memory horizon na, and since it depends on the external previous inputs u(k−1) up to u(k−nb), it is said to be exogenous. The extension to nonlinear systems is not straightforward (see [67]). An overview of nonlinear black-box identification methods by the end of last 8 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.3. DATA-BASED MODELS century was done by Sj¨oberg et al. in [141]. Nonlinear ARX models were introduced in the late 1980’s [41, 76], such that y(k+ 1) = fy(k), . . . , y(k−na), u(k), . . . , u(k−nb).(1.2) Previous outputs and inputs have an effect on the following output, so they can be grouped in the state x(k) = (y(k), . . . , y(k−na), u(k−1), . . . , u(k−nb)),(1.3) where na∈N0represents the memory horizon for the outputs, and nb∈N0 the memory horizon for the inputs. Although the state is fully defined for nb≥1, with a slight abuse of notation nb= 0 indicates that no past term of the inputs appears in the state. Note that x∈Rnx, with nx= (na+ 1)ny+nbnu. Taking into account (1.3), (1.2) can be rewritten as y(k+ 1) = f(x(k), u(k)).(1.4) The conditions under which a system can be described as a NARX model were presented in [105]. These conditions have been widely studied during the last 40 years. One of the first results on this topic was given by Sontag [144], relating the existence of this model to the observability property of the system. Later, Chen and Billings [41] proved that local NARX models can be obtained if the system is locally observable. A comprehensive study of this problem, for local and global estimators, was presented by Levin and Narendra [77]. Further extensions can be found in [158], where global inputoutput models are proven valid via transversality theory. In [77, Theorem 3] it was proven that if a linearised model at an equilibrium point (xs, us) is observable, then the dynamics of the system can be locally described by a NARX model. When a global NARX model is required, the strong observability property must be ensured, from the linearised system, for instance. However, in [77, Theorem 6] the authors proved that only generic observability is necessary, which is a property that most systems enjoy in practice. From this result, the conditions required to ensure the existence of this model for a sampled continuous-time system can be derived: Theorem 1.1. Consider that each component of the system function, fi, with i∈ I ny 1, is continuous and such that for all ˜xwhere the gradient is null, i.e., ∇xfi(˜x, u) = 0, the Hessian matrix ∇xxfi(˜x, u)is non-singular. Then, a NARX model can be determined as valid for almost every input sequence. Besides, the model function ˆ fiis continuous and the horizons can be taken as na=nb= 2n, where nis the minimum number of states of the system. 9 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.3. DATA-BASED MODELS Note that nothing is mentioned here regarding the quality of the inputoutput data available to derive the model. This fact has a crucial effect on the determination of the NARX model function using, for example, ML techniques, and it is analysed in Appendix A. Methods to study the model order and structure from input-output data are presented in [26], [56] and [59]. Based on the available input-output data, and assuming that they are informative enough, the structure of the NARX model and its (hyper-)parameters must be determined. The field of machine learning offers a wide range and well-established methods to cope with this problem, as it will be summarized next. 1.3.2 Machine learning Machine learning [113] is a whole area of computer science, often said to be a subset of artificial intelligence. The objective of machine learning is to develop algorithms that automatically improve through experience. Here, automatically means without being explicitly programmed to do so; and to improve means that they should be able to perform better according to some metric/objective. Typically, this objective translates into being able to predict function values of given queries, from functions that are unknown a priori, and performance implies accuracy in the predictions. This thesis will be centered in the kinky inference method, which is a supervised regression technique. The static function to be learnt is defined as follows. Definition 1.1 (Ground truth function).Given two Hilbert spaces, W ⊆ Rnw and Y ⊆ Rny, referred to as input and output space, respectively, the ground truth function is defined as the mapping f:W → Y such that y=f(w).(1.5) Unless otherwise indicated, we will consider Euclidean metrics for both spaces, although J. Calliess proposed to use pseudo-metrics in the general kinky inference method [34], in order to account, e.g., for periodic signals. As it was stated before, supervised regression algorithms are going to be considered. Hence, it is assumed that a set of NDinput-output observations of the ground truth function is available. They may be also be referred to as samples,measurements, or data points. Note that in a real scenario, the measurements of the output are subject to possible noise ε(k)∈Rny, and hence the noise-corrupted observation is denoted ˜y, such that ˜y(k+ 1) = f(w(k)) + ε(k).(1.6) 10 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.3. DATA-BASED MODELS Then, the observations are gathered in the data set D={(˜yi, wi), i ∈ I ND 1}.(1.7) The data set containing only inputs will be denoted WD. From now onwards, it is assumed that the possible noise is bounded. Assumption 1.1. The output noise ε∈Rnyis upper-bounded by some constant ¯ e∈Rny, that is, |ε| ≤¯ e.(1.8) In other words, the noise lies in a set defined as E ⊂ Rny={ε:|ε| ≤ ¯ e}, that is, ε∈ E =B(¯ e). This is considered a sensible assumption and besides, notice that the discrete nature of every real system forbids the measurements to become infinite. In addition, the research in this thesis will often be applied to compact spaces, where every signal is bounded. Then, a machine learning algorithm will be used to build the predictor ˆy=ˆ f(w;θ, D),(1.9) where ˆystands for predicted outputs and θ∈Rnθstands for the set of (hyper-)parameters of the predictor. This predictor may use the set of training data Dto make predictions, or only use it during the learning (tuning) process of θ. This thesis will frequently make use of the term prediction error, which is defined as follows: Definition 1.2 (Prediction error).Given the ground truth function fand a predictor ˆ ffor this function, the prediction error d∈Rnyfor a given query w is computed as d(w) = |f(w)−ˆ f(w;θ, D)|. Besides, we will denote the maximum prediction error as follows: Definition 1.3 (Worst-case prediction error).The maximum prediction error using the predictor ˆ f, given the parameters θand data set Dis denoted µ, such that |ˆ f(q;θ, D)−f(q)|=d(q)≤µ∈Rny,∀q∈ W.(1.10) Then, the last step in order to fully identify a dynamical NARX model consists in fixing a structure and tuning the model’s parameters accordingly. Hence, the dynamical system is transformed into a static function, suitable to be processed by machine learning algorithms, via regression. To this end, we define the regressor as follows: 11 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.5. MODEL PREDICTIVE CONTROL to more basic approaches, where an explicit solution can be computed (like PID [12]). In addition, the plant may be subject to operation constraints, restricting both the control actions to u(k)∈ U and the states to x(k)∈ X, for all k. The resulting optimization problem in optimal control schemes is denoted P∞(x(k)): min κ∞(x)V∞(x(k)) (1.20a) s.t. x(j|k) = F(x(j|k), κ∞(x(j|k))), j ∈ I ∞ 0,(1.20b) κ∞(x(j|k)) ∈ U, j ∈ I ∞ 0,(1.20c) x(j|k)∈ X, j ∈ I ∞ 0.(1.20d) This yields the best possible controller, considering the stage cost, and it guarantees stability of the closed-loop system and constraints satisfaction for every initial state for which (1.20) is feasible, as proven in Appendix C. However, the solution to this problem may be impossible to compute, given the potential infinite number of decision variables. Note that equation (1.20b) implies a differential constraint defined by the systems’ transition function, so it leads to a dynamical optimization problem, which is known as a variational problem. The solution to this optimization problem is known for especial cases, such as unconstrained linear systems with quadratic cost functions (which leads to the so-called linear quadratic regulator (LQR) [31]). In general, this variational problem can be solved by determining the solution of certain set of equations derived from the Hamilton-Jacobi-Bellman equations [16] or from Pontryagin’s minimum principle [127], which are very difficult to solved. This impossibility to solve optimal control problems led to model-based predictive controllers (MPC), which are based on two premises. First, a fix, finite prediction horizon is used, typically denoted N, allowing numerical solving of the OCP. Second, a receding horizon policy is used, such that given the obtained control sequence of Nterms u(k) = (u(k), u(k+ 1), . . . , u(k+N)),(1.21) only the first one is applied to the system, yielding a new state and computing the problem again every time step. This receding horizon technique provides feedback and continuous updating of the cost. The standard MPC problem, 18 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.5. MODEL PREDICTIVE CONTROL denoted PN(x(k)), is the following: min u(k)VN(x(k),u(k)) (1.22a) s.t. x(j|k) = F(x(j|k), u(k+j)), j ∈ I N−1 0,(1.22b) u(k+j)∈ U, j ∈ I N−1 0,(1.22c) x(j|k)∈ X, j ∈ I N−1 0,(1.22d) x(N|k)∈ Xf,(1.22e) where the optimization cost is VN(x(k), u(k)) = N−1 X j=0 `(x(j|k), u(k+j)) + Vf(x(N|k)),(1.23) which includes a terminal cost Vf(·). Besides, the terminal constraint (1.22e) is considered, forcing the system to be at the terminal region Xfat time step k+N. MPC has had an enormous success in industry applications [129], given its optimal nature and its ability to manage constraints. However, taking a finite horizon may loose the well-established properties of optimal controllers. For this reason, stability of MPC has been thoroughly studied since the 1980’s. Authors that have laid the foundations of model predictive control are D.Q. Mayne and J.B. Rawlings [133], E.F. Camacho and C. Bordons [38], M. Morari [117], and J. Maciejowski [92], among others. The extension to nonlinear MPC has been studied, together with the previous authors, by F. Allg¨ower and R. Findeisen [7], L. Gr¨une [55] or D. Limon [79], among others. Mayne et al. [107] proved that both stability and constraints satisfaction can be ensured adding an appropriate terminal cost and terminal constraint, as shown in Appendix C. Recently, it has been proven that an MPC without terminal conditions can be stable for large memory horizons [55]. In contrast to more basic control strategies, the method requires an accurate prediction model of the system and to solve an optimization problem each sampling time, so it demands higher computational capacity and an advanced knowledge of the optimization field. Some optimization techniques that have been studied for solving the MPC formulation have been proposed by D.P. Bertsekas [20], L.T. Biegler [22], S. Boyd [27] or M. Diehl [64]. Several sub-fields in MPC have also been of interest to the research community, like distributed MPC, which has been studied in depth by E. Camponogara [39], P.D. Christofides, R. Scattolini and D. Mu˜noz de la Pe˜na [42], among others; or economic MPC, which can be found in the 19 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.5. MODEL PREDICTIVE CONTROL works by M. Diehl [44], R. Amrit [10], D. Angeli [11, 132], M. M¨uller [118] , T. Faulwasser [47] or S. Lucia [89]. The problem of tracking MPC, for changing references has been studied by D. Limon, A. Ferramosca, A.H. Gonz´alez or M. M¨uller [84, 85]. Further literature review, together with an introduction to the stability and robustness analysis of model predictive control is presented in Appendix C. 1.5.1 Learning-based MPC As it was introduced in Section 1.2, recently, machine learning techniques are being used in the control field to enhance the performance of the resulting closed-loop. Many control methodologies that use data in their processes have gained attention in the last few years. Among others, it is worth mentioning adaptive control, for which the data-driven approach was recently reviewed by M. Benosman in [17]. Some adaptive controllers that merge data and MPC techniques can be found in the works by Adetola et al. [1] or Tanaskovic et al. [153]; the former for systems linear in their parameters and the latter for linear time-varying systems. Likewise, dual control seeks to identify and to control the plant at the same time. The derivation of data-based approaches for dual control can be found in the recent survey by A. Mesbah [110]. Apart from that, reinforcement learning is another useful technique that benefits from the experience through data observation. As we will discuss below, the learning-based approach usually enhances the performance of the controllers, although few results address the safety and stability analysis. Regarding reinforcement learning, Berkenkamp et al [19] recently proposed a safe approach with stability guarantees, able to improve the performance and to explore the safe region of the state-space, under some assumptions. Gros & Zanon have recently proposed a RL scope tailored to economic MPCs [54]. Other techniques do not require a model of the system, for instance, the direct approach is a data-driven technique in which a set of past trajectories is stored, so the control law is learnt, as shown in [125]. Later, Salvador et al. [137] proposed an offset-free controller based on these historian databases. Among other advantages, there is no tuning step, and new data may easily be incorporated into the process. Identification via machine learning provides tools that favour the determination and improvement of prediction models. Model predictive control 20 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.5. MODEL PREDICTIVE CONTROL provides the methodologies and theoretical support for the practical implementation of optimal control strategies. Therefore, apart from the techniques mentioned above, the encounter of both fields has recently gained increasing attention in the research community. MPC can benefit from the use of machine learning algorithms, since data-driven techniques have proven useful to various purposes, like learning (or improving) the prediction model or learning the control laws, as reviewed in [63]. Learning the system model versus learning the control law does not necessarily have to be opposite approaches. Piga et al. recently proposed a method that takes into account the control objective in the model learning step [124]. In the last few years, machine learning techniques have proven valid to describe wide ranges of classes. Consider for example the autonomous driving applications that are currently being developed by the groups of F. Borrelli [135], S. Di Cairano [43] or M. Zeilinger [70], among others. One of the most challenging issues of using learning-based MPC controllers is how to explicitly recover the safe-by-design properties of nominal formulations. Most off-the-shelf results do not provide stabilizing designs or do not consider modelling errors or constraints, forbidding the deployment of a stability analysis of the closed-loop system. The inherent uncertainty in the ML-based prediction model brings one of the most challenging problems in learning-based control: safety. When the system to be controlled is subject to limits in the manipulable inputs or in the measurable outputs, the resulting controller must make a decision coping with this uncertainty in order to provide a safe evolution of the system. This requires the ML method to provide a description of the uncertainty or prediction error between the real evolution of the plant and the estimated one. Within this topic, there are two approaches that arise when working with data-based schemes: deterministic and stochastic. Then, many applications could benefit from the stability analysis, as for example the recent works on robotic systems by Fisac et al. [52] using the former, or Ostafew et al. [122] using the latter. Dealing with the uncertainty following a stochastic approach is an open problem in the MPC community. During the last few years several research groups have proposed excellent techniques to address this problem, either studying safe control methodologies (see for example [88]) or using identification methods that provide a probabilistic distribution on the uncertainty of the predictions, like Gaussian processes regression [131]. The estimation on the prediction error offers a probabilistic bound, used to derive safe controllers, typically according to a given confidence interval. For example, the group of M. Zeilinger uses Gaussian processes to design stochastic robust bounds for MPC, via designing safety filters [156] or cautious controllers [62]. 21 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.6. MAIN RESULTS Other couple of research lines that consider the stability analysis of Gaussian process-based MPCs are those recently proposed by Koller et al. [74] and Maiworm et al. [93]. Besides, both controllers are able to adapt themselves, enhancing the performance by taking into account new measurements. On the other hand, the problem can be addressed from a deterministic approach. As it was mentioned before, the safe and robust design of controllers requires bounded errors, but there are few methods that guarantee an ensured bound on the uncertainty. Lipschitz interpolation methods [15] provides a bound on the prediction error, as the robust design requires. They are also referred to as nonlinear set membership (NSM), used by Canale et al. to derive nonlinear data-based MPCs [40]. Lorenzen et al. have recently proposed an online framework that ensures stability of linear systems with parametric uncertainty and additive disturbance [87]. Besides, Terzi et al. also proposed recently a general framework for state-estimation using NSM [154], for linear systems as well. Another robust predictive controller for linear systems was recently proposed by Soloperto et al [143], considering state-dependent uncertainties. Lucia & Karg presented in [90] the power of deep-learning approaches to approximate robust NMPC laws. Summing up, on the one hand, stochastic approaches yield less rigorous analysis of the stability results [109], so deterministic scopes are needed when ensured robust guarantees are sought. On the other hand, there are few results that address the problem of controlling nonlinear black-box systems while taking into account the uncertainty, in order to consider a robust stability scenario. This thesis aims to contribute to this topic, offering the results that are sketched in the following section. 1.6 Main results This section presents the main contributions of this thesis, whose structure can be found in the following section. As it was stated before, these contributions belong to two fields: machine learning and model predictive control. An improved kinky inference algorithm, named projected KI, has been developed in order to reduce computational times of the prediction method. A new approach on H¨older continuity of functions has been proposed, which explicitly takes into account the contribution of each input onto each output, named componentwise H¨older continuity. 22 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.7. THESIS OVERVIEW Based on the latter, an novel KI predictor has been proposed: the componentwise H¨older kinky inference (CHoKI) method, which significantly reduces estimation errors. A new method to estimate the H¨older constant from the available data has been proposed. It is an hybrid method that merges the two most commonly used methods: LACKI and POKI. The online version of the KI algorithms has been studied. Several policies to update the data set recursively have been proposed. Two new control laws have been derived for learning-based model predictive controllers subject to uncertainty: a controller for systems with constraint on the inputs and its extension to systems with constraints also in the outputs. These novel control laws do not require a terminal constraint, avoiding the need for the calculation of any kind of invariant set, which is cumbersome in learning-based control. The procedure and analysis to achieve stability by design are provided. Based on the CHoKI predictor, predictive controllers with larger domain of attractions have been presented, benefiting from two facts: lower and more accurate uncertainty bounds and larger feasibility regions. An online learning MPC based on two different prediction models has been proposed. A safe model guarantees robust stability, while an updated model improves the performance by means of including new measurements during the operation. An online learning MPC for tracking based on an exploration-exploita- tion criteria has been proposed. Even in the case of low-dense data sets, the closed-loop system is able to explore unseen regions of the space, driving the system to the best reachable reference while robustly satisfying constraints. 1.7 Thesis overview The rest of the thesis is structured as follows: Chapter 2 proposes a robust predictive controller that uses a model derived with the class of kinky inference. The controller formulation considers 23 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.7. THESIS OVERVIEW soft constraints in the outputs, hard constraints in the inputs and guarantees closed-loop robust stability as well as performance by means of the use of different control and prediction horizons and a weighted terminal cost, without the need of a terminal region. With respect to the learning method, a modification of the standard KI algorithm is proposed: the projected KI, which aims to reduce computational times, in order to allow the application to real-time systems. Chapter 3 addresses the problem of robust and stable control of systems subject to hard constraints in the outputs. A learning-based MPC and a design procedure to guarantee that the closed-loop system satisfies the constraints while steering the system to the reference without the calculation of a terminal invariant set are proposed. This is one of the key contributions of the thesis. Chapter 4 presents an improvement of kinky inference, which is able to explicitly take into account the different effects of the inputs of the ground truth function into different outputs, in the so-called componentwise H¨older approach. The proposed new predictor yields smaller prediction errors than the standard KI. The controller presented in Chapter 3 is extended to the proposed CHoKI prediction model, enhancing the performance and with larger domains of attraction than the previous ones. Chapter 5 exploits the fact that during the operation of the system, new measurements become available. Hence, incorporating these new measurements into the model, in an online fashion, could improve the prediction and therefore the performance. To this end, in this chapter an online learning MPC is proposed, based on a double prediction model: an initial safe model is in charge of the constraints satisfaction and feasibility of the controller, while an updated model is used to improve the performance of the closed-loop. Chapter 6 the online learning approach presented in the previous chapter, but taking an exploration-exploitation approach that cautiously forces the system to stay in a safe region in which stability is guaranteed, while preventing the prediction model from becoming computationally intractable. A model predictive control for tracking is proposed which, together with the exploration policy, allows the system to move to unexplored regions while satisfying constraints, even in the case of low-dense initial data sets. In addition, three appendices are included. Appendix A analyses several facts that must be taken into account when dealing with data-based methods. Appendix B presents further results of the kinky inference class of learning rules, and Appendix C addresses the stability analysis of standard MPC. 24 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.8. PUBLICATIONS 1.8 Publications Part of the research developed during the elaboration of this thesis has been published in the following papers: Journal papers Manzano, J. M., Limon, D., Mu˜ noz de la Pe˜ na, D., and Calliess, J.-P. Robust learning-based MPC for nonlinear constrained systems. Automatica 117 (2020), 108948 Manzano, J. M., Limon, D., Mu˜ noz de la Pe˜ na, D., and Calliess, J.-P. Output feedback MPC based on smoothed projected kinky inference. IET Control Theory & Applications 13, 6 (2019), 795–805 Manzano, J. M., Limon, D., Mu˜ noz de la Pe˜ na, D., and Calliess, J.-P. Online learning constrained MPC based on double prediction. Accepted in the International Journal of Robust and Nonlinear Control (2020) Manzano, J. M., Mu˜ noz de la Pe˜ na, D., Calliess, J.-P., and Limon, D. Componentwise H¨older inference for robust learning-based MPC. Submitted to the IEEE Transactions on Automatic Control (2020) Manzano, J. M., Mu˜ noz de la Pe˜ na, D., and Limon, D. Oraclebased economic predictive control. Submitted to Computer Applications in Chemical Engineering (2020) Book chapters Alamo, T., Manzano, J. M., and Camacho, E. Robust design through probabilistic maximization. In Uncertainty in Complex Networked Systems. Springer, 2018, pp. 247–274 Conference papers Manzano, J. M., Limon, D., Mu˜ noz de la Pe˜ na, D., and Calliess, J.-P. Data-based robust MPC with componentwise H¨older 25 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
1.8. PUBLICATIONS kinky inference. In 2019 IEEE 58th Conference on Decision and Control (CDC) (2019), IEEE, pp. 6449–6454 Manzano, J. M., Nadales, J. M., Limon, D., and Mu˜ noz de la Pe˜ na, D. Oracle-based economic predictive control. In 2019 IEEE 58th Conference on Decision and Control (CDC) (2019), IEEE, pp. 4246– 4251 Blaas, A., Manzano, J. M., Limon, D., and Calliess, J.-P. Localised kinky inference. In 2019 18th European Control Conference (ECC) (2019), IEEE, pp. 985–992 Manzano, J. M., Limon, D., and Mu˜ noz de la Pe˜ na, D. T´ecnicas de aprendizaje automatizado para la operaci´on econ´omica basada en datos de sistemas ciberf´ısicos. XVII Simposio CEA de Ingenier´ıa de Control y V Seminario de Innovaci´on Docente en Autom´atica (2019) Manzano, J. M., Limon, D., Mu˜ noz de la Pe˜ na, D., and Calliess, J.-P. Robust data-based model predictive control for nonlinear constrained systems. IFAC-PapersOnLine 51, 20 (2018), 505–510 Maiworm, M., Limon, D., Manzano, J. M., and Findeisen, R. Stability of Gaussian process learning based output feedback model predictive control. IFAC-PapersOnLine 51, 20 (2018), 455–461 Manzano, J. M., Limon, D., Alamo, T., and Callies, J.-P. Control predictivo basado en datos. Actas de las XXXVIII Jornadas de Autom´atica (2017) Nadales, J. M., Manzano, J. M., Barriga, A., and Limon, D. FPGA parallel implementation of Lipschitz interpolation techniques. In Submitted to the 59th IEEE Conference on Decision and Control (2020), IEEE Manzano, J. M., Calliess, J.-P., Mu˜ noz de la Pe˜ na, D., and Limon, D. Online learning robust MPC: an exploration-exploitation approach. In Accepted for publication in the 21st IFAC World Congress (2020), IFAC 26 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
Chapter 2 Output feedback MPC based on smoothed projected kinky inference 2.1 Introduction As a non-parametric regression method, the computational effort when using KI for making predictions grows with the sample size. In MPC, the prediction model is frequently evoked to compute large numbers of predictions. Hence, the computational complexity of non-parametric methods often forbid their applicability in MPC frameworks. Therefore, reducing computational effort for prediction is crucial. To this end, in this chapter we propose a variation of kinky inference that is designed to speed up predictions. Instead of making predictions on the basis of all observed data, our approach is based on using a different reduced subset of data for each region of a partition of the state space, to implement the smoothed kinky inference method presented in [34]. The use of different reduced data sets depending on the current state decreases the computational effort of the inference method. The smoothed predictions yield smooth control actions that result in better performance. The whole prediction method is called smoothed projected kinky inference (SPKI). In addition, a new controller formulation that guarantees closed-loop robust stability by means of the use of different control and prediction horizons and a weighted terminal cost is also presented. This controller is an extension of the controller proposed in [83], which provides a stabilizing design without 27 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.3. DATA-BASED PREDICTIVE CONTROL To derive the predictive controller, a positive definite stage cost function `(y, u) is defined as: `(y, u) = `t(y, u) + `b(y), where the term `t(y, u) penalizes the tracking error of inputs and outputs w.r.t. the reference set-point (assumed to be given by (0,0)) and `b(y) is a barrier function that implements soft constraints on the outputs. The barrier cost function `b(·) must satisfy that `b(y) = 0 for all y∈ Ysoft and that there exists two K-functions αband βbsuch that αb(d(y, Ysoft)) ≥`b(y)≥βb(d(y, Ysoft)) ∀y, where d(y, Ysoft) is a measure of the distance from yto Ysoft. Assumption 2.1. There exist two K-functions αyand αusuch that `(y, u)≥αy(kyk) + αu(kuk). Besides, `(y, u)is assumed to be continuous. The optimization problem PNc,Np(x(k),D) to be solved is the following: min u VNc,Np(x(k),u) = Nc−1 X i=0 `(ˆy(i|k), u(k+i)) + Np−1 X i=Nc `(ˆy(i|k), κf(ˆx(i|k)) + λVf(ˆx(Np|k)) (2.9a) s.t. ˆx(0|k) = x(k) (2.9b) ˆx(j+ 1|k) = ˆ F(ˆx(j|k), u(k+j); θ, D), j ∈ I Np−1 0(2.9c) ˆy(j|k) = Mˆx(j|k), j ∈ I Np−1 0(2.9d) u(k+j)∈ U, j ∈ I Np−1 0,(2.9e) where λ≥1 is a design parameter of the controller, Np∈Nis the prediction horizon, Nc≤Npis the control horizon and Uis the compact set defined by the constraints on the inputs. Asymptotic stability is ensured by designing a suitable terminal cost function Vf(·) and a terminal control law κf(·). Notice that no terminal constraint is considered, since the terminal cost function Vf is weighted by a factor λ≥1, following the procedure proposed in [83], and extended for a prediction horizon larger than the control horizon. The MPC control law is u(k) = κMPC(x(k); θ, D) = u∗(0), which depends on the NARX state of the system and the collected historical data set D. 34 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.3. DATA-BASED PREDICTIVE CONTROL 2.3.1 Stability analysis of the controller In order to demonstrate the stability of the proposed controller, first the asymptotic stability of the nominal case will be proven, that is, when the proposed controller is used to regulate a system whose model is equal to the prediction one. Theorem 2.1. Suppose there exists a terminal control law κf(·) : Rnx→Rnu, a terminal cost function Vf(·) : Rnx→Rand a region Ωγ={x:Vf(x)≤γ}, where γis a positive constant such that for all x∈Ωγ:2 α1(kxk)≤Vf(x)≤α2(kxk), Vf(x+)−Vf(x)≤ −`t(y, κf(x)), Mx ∈ Ysoft, κf(x)∈ U. Consider the so-called domain of attraction, defined as XNc,Np(λ) = nx:V∗ Nc,Np(x(k)) ≤(Np−n)φ+λγo, where φis a positive constant and n= max(na, nb). Then, given Np≥n, for all λ≥1and all x(0) ∈XNc,Np(λ), the closedloop system x+=ˆ F(x, κMPC(x;θ, D); θ, D)is asymptotically stable. Proof. Notice that if x∈Ωγ,y∈ Ysoft so `b= 0, and thus `(y, u) = `t(y, u). For the sake of clarity it may be used xjor x(j|k) indistinctly. The proof uses the following three lemmas: Lemma 2.1. If the assumptions of Theorem 2.1 hold, then V∗ Nc,Np(x(k)) ≤λVf(x(k)),∀x(k)∈Ωγ. Proof.For any system controlled by u(k+j) = κf(x(j|k)), Vf(x(j|k)) −Vf(x(j+ 1|k)) ≥`(y(j|k), κf(x(j|k))). 2The notation x+implies x(k+ 1), given x(k). 35 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.3. DATA-BASED PREDICTIVE CONTROL Summing for the whole sequence: Np−1 X j=0 Vf(x(j|k)) −Vf(x(j+ 1|k)) = Vf(x(k)) −Vf(x(Np|k)) ≥ Np−1 X j=0 `(y(j|k), κf(x(j|k))) and by optimality, provided that λ≥1, ∀Np≥Nc≥1 λVf(x(k)) ≥ Np−1 X j=0 `(y(j|k), κf(x(j|k))) + λVf(x(Np|k)) ≥V∗ Nc,Np(x(k)).(2.10) Lemma 2.2. If x∗(j|k)/∈Ωγ,∀j= 0,··· , Np−1, then V∗ Nc,Np(x(k)) ≥(Np−n)φ+λγ. Proof.Let define `ext(x, u) as the stage cost extended to consider a positive definite function of the whole state: `ext(x(k), u(k)) = 1 n n X i=0 `(y(k−n+i), u(k−n+i)), where n= max(na, nb), and with a slight abuse of notation, the terms given by y(k−n+i) and u(k−n+i) are set to 0 if y(k−n+i) or u(k−n+i) are not included in x(k), respectively. Then, φis defined as a positive constant such that ∀x /∈Ωγ, `ext(x, u)≥ φ. Hence, V∗ Nc,Np(x(k)) = Np−1 X j=0 `(y∗(j), u∗(j)) + λVf(x∗(Np)) ≥ jNp nk−1 X k=0 n X i=0 `(y∗(k+n−i), u∗(k+n−i)) + λVf(x∗(Np)) ≥Np nnφ +λγ ≥(Np−n)φ+λγ. 36 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.3. DATA-BASED PREDICTIVE CONTROL Lemma 2.3. If the assumptions of Theorem 2.1 hold, and if x∗(Np|k)/∈Ωγ, then x∗(j|k)/∈Ωγ,∀j= 0, . . . , Np−1. Proof.This lemma is proven by contradiction: (i) Assume that there exists any instant i<Ncin which x∗ i∈Ωγ. From Lemma 2.1 it is inferred that for a problem of horizons Nc−i, Np−i, λVf(x∗ i)≥V∗ Nc−i,Np−i(x∗ i). From Bellman’s optimality principle, and provided that `(y, u) is positive definite, it is derived that V∗ Nc−i,Np−i(x∗ i)≥λVf(x∗ Np), and given that x∗ Np/∈Ωγ, then λVf(x∗ Np)> λγ. Hence, λVf(x∗ i)> λγ, which is a contradiction, so x∗ i/∈Ωγ. (ii) Contrary, assume there exists an instant i≥Ncin which x∗ i∈Ωγ. Then, x∗ j+1 =ˆ F(x∗ j, kf(x∗ j)) for j=i, . . . , Np. Since Ωγis a positive invariant set for this system, we have that x∗ j∈Ωγfor j=i, . . . , Np, leading again to a contradiction. Given these three lemmas, we have that ∀x∈XNc,Np(λ), if x∗ Np/∈Ωγ, then x∗ j/∈Ωγ∀j= 0, . . . , Np−1, and therefore V∗ Nc,Np(x(k)) ≥(Np−n)φ+λγ. So on the contrary, if V∗ Nc,Np(x(k)) ≤(Np−n)φ+λγ, then x∗ Np∈Ωγ. Next, we are ready to prove the statement of the theorem. First, recursive feasibility is proven: Assume that x(k)∈XNc,Np(λ), then x∗ Np∈Ωγ. Define the sequence of future control inputs for x(k+ 1) as ˜u(k+j+ 1) = u∗(k+j)j= 0, . . . , Nc−1 κf(˜x(j|k+ 1)) j=Nc, . . . , Np−1, and define ˜x(j|k+ 1) as the predicted trajectory. Then, from standard arguments of MPC proofs, this sequence is feasible for the optimization problem, 37 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.3. DATA-BASED PREDICTIVE CONTROL since u∗(k) is feasible. Besides, the difference between the costs associated to these trajectories is: ˜ VNc,Np(x(k+ 1)) −V∗ Nc,Np(x(k)) = −`(y(k), u∗(k)) +h`(˜y(Np|k+ 1), κf(˜x(Np|k+ 1))) + λVf(˜x(Np|k+ 1)) −λVf(x∗(Np|k))i. The sum among brackets is negative, since x∗ Np∈Ωγ. Thus, ˜ VNc,Np(x(k+ 1)) −V∗ Nc,Np(x(k)) ≤ −`(y(k), u∗(k)), and by optimality, V∗ Nc,Np(x(k+ 1)) ≤˜ VNc,Np(x(k+ 1)), V∗ Nc,Np(x(k+ 1)) −V∗ Nc,Np(x(k)) ≤ −`(y(k), u∗(k)).(2.11) So since x(k)∈XNc,Np(λ) implies that V∗ Nc,Np(x(k)) ≤(Np−n)φ+λγ, then V∗ Nc,Np(x(k+ 1)) ≤(Np−n)φ+λγ, implying that x(k+ 1) ∈XNc,Np(λ). Therefore the closed-loop system is recursively feasible. Then, asymptotic stability is proven. Given the definition of x(k), k X j=k−n `(y(j), u(j)) ≥ n X j=0 αy(ky(k−j)k) + αu(ku(k−j)k) ≥αx(kx(k)k) + αu(ku(k)k) ≥αx(kx(k)k),(2.12) where n= max(na, nb) and αy, αu,and αxare K-functions. Defining W(xk) = n X j=0 V∗ Nc,Np(x(k−j)), we can infer that W(x(k+ 1)) −W(x(k)) = V∗ Nc,Np(x(k+ 1)) −V∗ Nc,Np(x(k−n)) ≤ − k X j=k−n `(y(j), u(j)) ≤ −αx(kx(k)k), where the first inequality is derived by recursion of (2.11), and the second inequality comes from (2.12). Besides, W(x)≥ k X j=k−n `(y(j), u(j)) ≥αx(kxk), 38 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.3. DATA-BASED PREDICTIVE CONTROL and, by Prop. B.25 in [133], there exists a K-function α2(·) such that W(x)≤ α2(kxk). Hence, W(x(k)) is a Lyapunov function of the system, which is asymptotically stable. Remark 2.1. Whereas Npand λare design parameters, φand γare problemdependent parameters. They require solving three global optimization problems, whose ingredients are sketched below. According to the definition of γin Theorem 2.1, let us define γ1as: γ1= arg min xVf(x) (2.13) s.t. λVf(ˆ F(x, κf(x); θ, D)) −λVf(x)≥ −`(Mx, κf(x)). It must also be guaranteed that for all Vf(x)≤γ,κf(x)∈ U, so we define γ2= max γγ(2.14) s.t. κf(x)∈ U Mx ∈ Ysoft,∀x∈ {x:Vf(x)≤γ}. Then, γ= min(γ1, γ2). In addition, φis obtained as: φ= min `ext(x, κf(x)) (2.15) s.t. Vf(x)≥γ, where the extended stage cost `ext is defined in Lemma 2.2 (on page 36). Note that these calculations are done offline. Corollary 2.1. Under the assumptions of Theorem 2.1, it can be proven that: (i)XNc,Np(λ1)⊆XNc,Np(λ2)∀λ1≤λ2(2.16) (ii)XNc,N1(λ)⊆XNc,N2(λ)∀N1≤N2(2.17) Proof. First, it is proven that if λ1≤λ2and considering that x∈XNc,Np(λ1), then x∈XNc,Np(λ2). Define ˆu∗,ˆx∗and ˆ V∗ Nc,Npas the optimal solution for x∈ XNc,Np(λ1) and ˜u∗,˜x∗and ˜ V∗ Nc,Npas the optimal solution for x∈XNc,Np(λ2). 39 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.3. DATA-BASED PREDICTIVE CONTROL Taking into account that ˆu∗is a feasible solution of the optimization problem for λ2, ˜ V∗ Nc,Np(x) = Np−1 X i=0 `(˜x∗ i,˜u∗ i) + λ2Vf(˜x∗ Np) ≤ Np−1 X i=0 `(ˆx∗ i,ˆu∗ i) + λ2Vf(ˆx∗ Np) = Np−1 X i=0 `(ˆx∗ i,ˆu∗ i) + λ1Vf(ˆx∗ Np) + λ2Vf(ˆx∗ Np)−λ1Vf(ˆx∗ Np) =ˆ V∗ Nc,Np(x)+(λ2−λ1)Vf(ˆx∗ Np) ≤ˆ V∗ Nc,Np(x)+(λ2−λ1)γ ≤(Np−n)φ+λ1γ+ (λ2−λ1)γ = (Np−n)φ+λ2γ, and hence x∈XNc,Np(λ2). Now it is proven that if N1≤N2and x∈XNc,N1(λ), then x∈XNc,N2(λ). Provided that V∗ Nc,N1(x(k)) ≤(N1−n)φ+λγ, ∀x∈XNc,N1(λ) and N1≤N2, from (2.10) it is derived that Vf(x∗ N1)≥ N2−1 X j=N1 `(y∗ j, κf(x∗ k)) + Vf(x∗ N2), which implies that V∗ Nc,N1≥V∗ Nc,N2, so V∗ Nc,N2≤V∗ Nc,N1≤(N1−n)φ+λγ ≤(N2−n)φ+λγ, and hence XNc,N1⊆XNc,N2∀N1≤N2. Remark 2.2. From standard theory of MPC, it is well known that increasing the control horizon Ncincreases the domain of attraction XNc,Np(λ), with the respective increase in the number of decision variables. Corollary 2.1 provides other two ways of increasing the domain of attraction of the MPC: (i) increasing the weighting factor λand (ii) increasing the prediction horizon Np. A larger weighting factor provides a larger domain of attraction, but the performance might be worse. On the other hand, a larger prediction horizon does not only increase the domain of attraction, but it also improves the performance, since V∗ Nc,N2≤V∗ Nc,N1. The number of decision variables does not increase with Np, but it increases the number of constraints in the optimization problem, which might lead to a mild increase of the computational cost. 40 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.3. DATA-BASED PREDICTIVE CONTROL Once the nominal asymptotic stability is proven, then it is shown that the real plant, controlled with the proposed predictive control law, is input- to-state stable w.r.t. the model mismatches under the following assumption. Assumption 2.2. It is assumed that the radius Riin equation (2.4) is large enough to ensure that the function ˆ fSP (·,·)is H¨older continuous. Definition 2.1. A system x+=f(x) + dis input-to-state stable (ISS) w.r.t. dif there exists a KL function βand a Kfunction αsuch that kx(k)k ≤ β(kx(0)k, k) + supj∈[0,k]α(kd(j)k). Now, the stability result is stated in the following theorem. Theorem 2.2. Assume that the stage cost function `(·,·), the terminal control law κf(·)and the terminal cost function Vf(·)are uniformly continuous functions in the feasibility region XNc,Np(λ). Assume that the hypotheses of Theorem 2.1 hold. Then, the real plant, controlled by u(k) = κMPC(x(k); θ, D) is input-to-state stable w.r.t. the model mismatch signal: d(k) = |f(x(k), u(k)) −ˆ f(x(k), u(k); θ, D)|. Proof. First, notice that the model of real plant can be posed as f(x(k), u(k)) = ˆ f(x(k), u(k)) + d(k) and the signal d(k) is bounded by µ. Given that the function ˆ f(·,·) is H¨older continuous, the model function ˆ F(·,·) is also H¨older continuous. Since the stage cost function `(·,·) and the terminal cost function Vf(·) are continuous functions and XNc,Np(λ) is a compact set, then they are uniformly continuous. From Proposition 1, case C1, in [81], it is inferred that the closed-loop system is ISS w.r.t. the signal that describes the model mismatch, d(k). Remark 2.3. This theorem resorts on Assumption 2.2, which may be difficult to demonstrate. Notice that the assumption can always be satisfied since for Di=Dthe model function is H¨older continuous. In case one cannot ensure H¨older continuity, stability is not lost. Instead, according to the definition of input-to-state practical stability (ISpS, see Definition 6 in [81]), it can be proven that if the discontinuity jump is upper bounded by c, then kˆ fSP (w+ ∆w, u)−ˆ fSP (w, u)k ≤ σ(c) + ˆ Lk∆wkp, where σis a K-function. Using similar arguments to the proof of Theorem 2.2 and Definition 6 in [81], it can be proven that the closed-loop system is ISpS, since there exists a KL-function βand K-functions αand ωsuch that: kx(k)k ≤ β(kx(0)k, k) + sup j∈[0,k] α(kd(j)k) + ω(c). 41 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.4. CASE STUDY Remark 2.4. The calculation of the control action requires the solution of a mathematical programming problem on line. This is typically a gradientbased algorithm that may be sensitive to abrupt variations of the gradients, which may be derived from discontinuity of the derivative of the models. As it is illustrated in the case study, this may produce convergence issues leading to spikes in the evolution of the control action. To avoid this effect, the tracking stage cost function can be modified, adding a term penalizing the variation of the control input, i.e. ku(k+i)−u(k+i−1)k2 S. The resulting tracking stage cost can be recasted, w.l.o.g., as a function of the state (that is, the regressor), `t(x, u), provided that nb≥1. The stability analysis presented can be extended to this case. 2.4 Case study In this section we consider the control of the continuously-stirred tank reactor presented in [139]. The experiments carried out to generate the data sets and the procedure to tune the proposed controller are described in detail. The input of the reactor is the reference of the coolant temperature Tr(K) and the output is the concentration of the reactant, CA(mol l−1), in the reaction A→B. The temperature of the tank and the coolant are given by Tand Tc(K), respectively. It is assumed that the evolution of the plant is given by the following set of differential equations: dCA(t) dt =q0 V·(CAf −CA(t)) −k0·exp −E R·T(t)·CA(t) (2.18a) dT(t) dt =q0 V·(Tf−T(t)) + (−∆Hr)·k0 %·Cp·exp −E R·T(t)·CA(t) +U·A V·%·Cp·(Tc(t)−T(t)) (2.18b) dTc(t) dt =Tr(t)−Tc(t) τc (2.18c) The parameters of the model are given in Table 2.1. Note that the model is only used to carry out simulations, no information is used to design the controller. It is also assumed that the concentration sensor adds an error of 2.5% of the measurement. The error is generated randomly for each measurement using an uniform distribution. The constraints in the input are Tmin r= 335 K and Tmax r= 372 K. 42 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.4. CASE STUDY Param. Definition Value Units q0Input flow of the reactive 10 l min−1 VLiquid volume in the tank 150 l k0Frequency constant 6 ×1010 min−1 E/R Arrhenius constant 9750 K −∆HrEnthalpy of the reaction 10000 J mol−1 UA Heat transfer coefficient 70000 J min−1K−1 %Density 1100 g l−1 CpSpecific heat 0.3 J g−1K−1 τTime constant 1.5 min CAf CAin the input flow 1 mol l−1 TfTemperature (input flow) 370 K Table 2.1: Parameters of the system First, the set of equilibrium points of the system is estimated using a sequence of steps in the input, from Tmin rto Tmax rwith 0.5 K increments, and each step long enough to reach a steady sate. The result is shown in Figure 2.1. Apart from obtaining the equilibrium points of the system, this test is used to adjust the sampling time, which is set to τs/20 = 30 s, where τs stands for the mean settling time of the sequence of steps applied. This test is also used to bound the output variable CA, which is defined by Cmax A= 0.8053 mol l−1and Cmin A= 0.2 mol l−1(note that the system has an inverse static characteristic). This restriction is treated by means of a soft constraint, using a barrier cost. In addition, an equilibrium point is chosen to be the reference operating point, given by Cref A= 0.439 mol l−1,Tref r= 356 K. After defining the set of equilibrium points, a set of experiments are carried out to obtain the data for the predictor. The experiments are designed using the methodologies presented in [134] to identify the dynamics of a system within a workspace: a relay test is carried out to estimate the crossover frequency of the system around the operating point [12], with limits Tmax r and Tmin r. The crossover frequency obtained is 3.7 mHz. Figure 2.2 shows the simulation. Then, a sequence of chirp signals covering the workspace are applied to generate the raw data set containing the trajectories of concentrations and temperatures, Draw. The parameters of the chirp signals (length of the signal τf, initial and final frequency f0, ff, center Tr0and amplitude A) are shown in Table 2.2. The simulation is shown in Figure 2.3. The set of 43 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.4. CASE STUDY P= 183.48 30.828 −0.2224 −1.1638 30.828 14.994 −0.1079 −0.5647 −0.2224 −0.1079 0.0008 0.0041 −1.1638 −0.5647 0.0041 11.021 , K=0.4856 0.2278 −0.0016 −0.0086 . Regarding Theorem 2.1, the parameters φand γhave been computed. Setting the control horizon to Nc= 4 and the prediction horizon to Np= 14 results in φ= 1.2046 and γ= 0.0294. If x(0) = [0.6 300 300]T, then the optimal cost is V∗ Nc,Np(x(0)) = 13.721, so x(0) ∈XNc,Np(λ). Regarding Theorem 2.2, as far as we analysed the outcomes of the predictor, Assumption 2.2 holds for every query point observed, making `(·,·), κf(·) and Vf(·) uniformly continuous. Once the controller is designed, a hundred simulations are carried out, with the same initial state, but a hundred different realizations of the 2.5% random noise. The following figures represent the output and the control input. Within each subplot, the constraints appear in black, the reference in red, and the 100 trajectories are represented by the gray band, whose mean value is plotted in blue. A first experiment is carried out in which the model used by the MPC is the real plant (defined by its set of ODEs (2.18)). The simulations provide a measure of the best closed-loop behaviour, achievable assuming perfect knowledge of the model. This result is shown in Figure 2.9. Secondly, an MPC controller based on the standard KI prediction is applied, using the data set D. Notice that neither the data set was divided into partitions, nor the filtered KI was implemented. The result of this simulation is shown in Figure 2.10. These results are not satisfactory, since the solver is failing to converge due to the fact that the KI prediction may be non-differentiable or due to the large number of iterations. Finally, the SPKI-MPC is applied, to address the issues mentioned above. The result is shown in Figure 2.11. The behaviour is smoother, with a smaller band than the previous experiment, which indicates greater robustness and a performance closer to the ideal one (Figure 2.9) than the simpler KI set. Apart from that, in order to compare the performance of the presented results, a hundred different simulations of the MPC (with random initial state) are carried out, whose results are compared in Figure 2.12. The boxplots represent maximum, minimum, median and deviation of the performance costs and computational times of the experiments. The performance index of each 50 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.4. CASE STUDY Model: ODEs. Plant: Noisy ODEs. 0 2 4 6 8 10 12 14 16 18 20 Time (min) 0.3 0.4 0.5 0.6 0.7 CA (mol/l) Reference Band Mean 0 2 4 6 8 10 12 14 16 18 20 Time (min) 340 350 360 370 Tr (K) Reference Constraints Band Mean Figure 2.9: Ideal ODEs-MPC. 0 2 4 6 8 10 12 14 16 18 20 Time (min) 340 350 360 370 Tr (K) Reference Constraints Band Mean Model: KI. Plant: Noisy ODEs. 0 2 4 6 8 10 12 14 16 18 20 Time (min) 0.3 0.4 0.5 0.6 0.7 CA (mol/l) Reference Band Mean Figure 2.10: Primary KI-MPC. simulation is defined as: Φ = tsim X i=1 `t(x(i), u(i)) + `b(y(i)). 51 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.4. CASE STUDY Model: SPKI. Plant: Noisy ODEs. 0 2 4 6 8 10 12 14 16 18 20 Time (min) 0.3 0.4 0.5 0.6 0.7 CA (mol/l) Reference Band Mean 0 2 4 6 8 10 12 14 16 18 20 Time (min) 340 350 360 370 Tr (K) Reference Constraints Band Mean Figure 2.11: Data-based SPKI predictive control applied to a CSTR. 52 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
2.4. CASE STUDY ODEs KI PKI SPKI 60 80 100 120 140 160 180 200 220 Performance Index Performance Index (a) Performance ODEs KI PKI SPKI 0 1 2 3 4 5 Time (s) Computational Time per Iteration (b) Computational time Figure 2.12: Boxplot comparison of the different prediction models for 100 simulations of the MPC. 53 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
Chapter 3 Robust learning-based MPC for nonlinear constrained systems 3.1 Introduction Robust constraint satisfaction is particularly important for safe operation of systems, and challenging for learning-based controllers, as they have to consider not only the possible perturbations of the system, but also the uncertainty introduced by the inference method. The main limitation of the controller presented in the previous chapter is that it cannot guarantee robust satisfaction of output constraints. In this chapter, the main contribution is a new robust predictive controller formulation for systems also subject to hard constraints on the outputs. A design method that takes into account the prediction error bounds in an explicit way to tighten the problem constraints and that guarantees closedloop constraint satisfaction and input-to-state stability (ISS) [81] is provided. One of the main characteristics of this design is that it is not based on a terminal region constraint. In general, this terminal constraint is based on robust invariant sets which, for the class of systems considered (that is, unknown systems, possibly nonlinear, for which a priori only input-output data are available), are difficult to obtain. Given the hard feature of the constraints, inference methods that provide a guaranteed bound on the prediction error must be used. Besides, knowledge of the Lipschitz constant allows the estimation of the reachable sets, in order to account for a robust design of the controller. For these two reasons, kinky inference is an useful method to address this problem. The procedure 54 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
3.2. PROBLEM SETTING to define the tightened constraints of the MPC optimization problem is provided. This procedure is specifically tailored to the KI model used, in order to obtain the least conservative possible bounds. The results in this chapter were published in [98, 102]. 3.2 Problem setting The problem setting is inherited from Chapter 2, with the additional hard constraint in the output, that is, y(k)∈ Y,∀k, (3.1) where Y ⊂ Rnyis a compact set. MPC requires repeated optimization of the predicted control inputs subject to constraints. Therefore, in order to give guarantees on the controller’s closed-loop performance, recursive feasibility and constraint satisfaction must be ensured. That is, it is necessary to ensure that all constraints remain satisfiable during runtime, or equivalently, to guarantee that the controlled system will not leave the feasibility region. However, since the controller will not be based on the ground-truth dynamics f, but on the learnt model ˆ f inferred from samples of the ground-truth, recursive feasibility can only be guaranteed if a bound on the discrepancy between fand ˆ fis known a priori, and taken into account by the controller. The LACKI method explained in Appendix B ensures that if the model function is H¨older and the noise is bounded, then the estimation error is bounded, which is required to design a deterministic robust controller to regulate the plant. Any worst-case guarantee inevitably requires a priori knowledge. Hence, in the following hypothesis, it is assumed that the maximum prediction error µis available for the design of the controller. Assumption 3.1. It is assumed that for L,p, and D, a bound on the error between the estimated output and the real output is known, denoted µ∈Rny, such that |ˆ fi(x, u;θ, D)−fi(x, u)| ≤ µi,(3.2) for all i∈ I ny 1,x∈ Yna+1 ×Unb, and u∈ U. Remark 3.1. From a practical point of view, the problem of how to calculate the error bound must be addressed. Kinky inference methods enjoy the property of providing a deterministic error bound if the Lipschitz constant and an 55 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
3.3. STABILIZING DATA-BASED NMPC upper bound of the noise are known (see Appendix B). Moreover, if a bound on the second derivative is known, it is also possible to derive an estimation error bound. If these parameters are not known a priori, which is usual in practice, then they must be estimated from experimental data. Consequently, the validity of the results presented in this section is conditioned to the validity of the estimated error bound. This is the reason why this is considered as a standing assumption. 3.3 Stabilizing data-based NMPC In this section, a model predictive controller is derived based on a prediction model learnt from data of the plant. Since the prediction model is not accurate, the effect of the estimation error on the predictions must be analysed, to be taken into account in the design of the controller. For this analysis, it is convenient to define the NARX model of the plant in a state–space form, which is done as in Chapter 2, in equations (2.8) (on page 33). Then, the proposed robust MPC is based on nominal predictions and tightened constraints. To guarantee robustness, a bound on the propagation of the prediction error (see Figure 3.1) is calculated from the following lemma: Lemma 3.1. Assume that at sampling time k, the state of the plant is x(k) and a sequence of future control inputs u(k+j)for j∈ I N−1 0is given. Let ˆx(j|k)and ˆy(j|k)be the predicted states and outputs, respectively, derived from (2.8a) for the given sequence of future control inputs and the current state x(k), i.e. ˆx(0|k) = x(k). Assume that at sampling time k+ 1, the current output y(k+ 1) is measured, and hence the current state x(k+1) is known. Based on these new measurements, an updated sequence of states and outputs ˆx(j|k+1) and ˆy(j|k+1) is predicted based on (2.8a) with ˆx(0|k+ 1) = x(k+ 1) and the remaining sequence of the given future control inputs. Let c1∈Rnybe a vector such that |y(k+ 1) −ˆy(1|k)| ≤ c1.(3.3) Then, the mismatch between the predictions satisfies 1 |ˆy(j−1|k+ 1) −ˆy(j|k)| ≤ cj, j ∈ I N 1,(3.4a) kˆx(j−1|k+ 1) −ˆx(j|k)kX≤rj, j ∈ I N 1,(3.4b) 1k·kXis a norm for the state–space such that kxkX=k(x, 0)kW. 56 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
3.3. STABILIZING DATA-BASED NMPC y(k) y(k+ 1) y(k+j) ˆy(1|k) ˆy(j|k) ˆy(j−1|k+ 1) cj c1 Figure 3.1: Propagation of the prediction error. where cj∈Rnyand rj∈Rare obtained from the recursion cj+1,i =Lirpi j,(3.5) and rj=kΞjkX,j∈ I N−1 1,i∈ I ny 1, where Ξj=B(cj)×···×B(cσ(j))×{0}×···×{0} | {z } nb+1−σ(j−1) times ⊆Rnx, with σ(j) = max(1, j −na). Proof. Provided that ˆy(j−1|k+ 1) = ˆ f(x(j−2|k+ 1, u(k+j−1)) and ˆy(j|k)) = ˆ f(ˆx(j−1|k), u(k+j−1)), it can be derived that ∀i∈ I ny 1 |ˆyi(j−1|k+ 1) −ˆyi(j|k))| ≤ Likˆx(j−2|k+ 1) −ˆx(j−1|k)kpi X. Note that ˆx(j−2|k+ 1) −ˆx(j−1|k) = ˆy(j−2|k+ 1) −ˆy(j−1|k),..., ˆy(σ(j−1) −1|k+ 1) −ˆy(σ(j−1)|k), 0,...,0. Then, ˆx(j−2|k+ 1) −ˆx(j−1|k)∈Ξj−1. Assuming that cj−1is known, kˆx(j−2|k+ 1) −ˆx(j−1|k)kX≤ kΞj−1kX=rj−1, which implies the stated result. 57 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
3.3. STABILIZING DATA-BASED NMPC Remark 3.2. Note that the formulation of the previous Lemma and in what is to follow treats each of the nyoutputs of fas different functions, considering different Li, pifor each i∈ I ny 1. Remark 3.3. If the infinity norm is chosen as the norm of the input space, then rj= max s∈ I j σ(j)kcsk∞. Based on the derived bounds on the prediction error, the problem of robust constraint satisfaction is addressed, by means of a set of tightened constraints on the outputs [133], computed offline for the maximum possible prediction error, i.e. taking c1=µ. These sets are defined as follows: Yj=Y B(dj),(3.6) where dj, known as back-off in MPC literature, is obtained as dj= j X s=1 cs.(3.7) These constraints sets will be used to prove recursive feasibility of the controller, following standard procedures. Lemma 3.2. The sets Yjare such that for all y∈ Yjand for all deviation ∆y∈ B(cj),y+ ∆y∈ Yj−1. Proof. Since for j≥1, dj=dj−1+cj, it follows that B(dj) = B(dj−1)⊕B(cj). By definition, y+ ∆y∈ Yj⊕B(cj) = Y B(dj)⊕B(cj), and hence Yj=Y B(dj) = Y B(dj−1)B(cj), so y+ ∆y∈ Yj⊕B(cj) =Y B(dj−1)B(cj)⊕B(cj) ⊆ Y B(dj−1) = Yj−1. 58 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
3.3. STABILIZING DATA-BASED NMPC In order to ensure that the proposed controller is feasible, the tightened set of constraints must be non-empty along the prediction horizon, as stated in the following assumption. Assumption 3.2. The prediction horizon Nand estimation error bound µ are such that the set YNis non-empty. Based on the previous definitions, the optimization problem PN(x(k); D) of the proposed predictive controller is: min u VN(x(k),u) = N−1 X i=0 `(ˆx(i|k), u(k+i)) + λVf(ˆx(N|k)) (3.8a) s.t. ˆx(0|k) = x(k) (3.8b) ˆx(j+ 1|k) = ˆ F(ˆx(j|k), u(k+j)), j ∈ I N−1 0(3.8c) ˆy(j|k) = Mˆx(j|k), j ∈ I N−1 0(3.8d) u(k+j)∈ U, j ∈ I N−1 0(3.8e) ˆy(j|k)∈ Yj, j ∈ I N−1 0,(3.8f) where λ≥1 is a weighting parameter and ˆ Fand Mdefine the state-space (see equations (2.8)). Note that this problem is nonlinear, non-convex and nondifferentiable. Its ingredients are required to meet the following assumption, which is similar to the standard MPC ones [133]: Assumption 3.3. 1. The stage cost function `(·,·)is a H¨older continuous positive definite function such that `(·,·)≥αy(kxkX)for a certain K-function αy, and its H¨older parameters are Lxand px. 2. There exists a control law u=κf(x), a function Vfand a level set Ωγ={x:Vf(x)≤γ} ⊆ Rnx for some γ > 0such that for all x∈Ωγthe following conditions hold: (a) Vfis a H¨older continuous positive definite function, with H¨older constants LVf, pVf. Given two K-functions αf, βf, αf(kxkX)≤Vf(x)≤βf(kxkX), Vf(ˆ F(x, κf(x))) −Vf(x)≤ −`(x, κf(x)). 59 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
3.4. CASE STUDY 50 100 150 200 250 0.4 0.6 0.8 1 50 100 150 200 250 280 290 300 310 Figure 3.2: Simulation obtained applying a sequence of chirp signals. from (3.8), but with the set of ODEs (2.18) as the state-feedback prediction model, for which µis the maximum noise, 0.025 mol l−1. This aims to resemble the ideal case of perfect knowledge of the plant. The result is shown in the first row of Figure 3.3. As expected, the data-based controller performs slower than the ideal since, unlike the latter, an output-feedback uncertain framework is considered. Second, the controller proposed in Chapter 2 is applied. It is based on KI and guarantees closed-loop stability, but does not take into account output constraints, so as shown in the second row of Fig. 3.3, the minimum CA limit is violated. To sum up, the controller proposed in this chapter is able to robustly satisfy hard constraints in the outputs, learning the model from input-output data with a closed-loop performance similar to the ideal case. 66 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
3.4. CASE STUDY 0 2 4 6 8 10 12 0.4 0.6 0.8 1 0 2 4 6 8 10 12 0.4 0.6 0.8 1 0 2 4 6 8 10 12 0.4 0.6 0.8 1 Figure 3.3: Closed-loop output for 100 simulations of the ideal MPC (top), the robust learning MPC for systems without output constraints (middle) and the proposed robust constrained learning MPC (bottom). The grey band groups the trajectories, the blue dashed line represents its mean, the green dotted one the reference and the black dash-dotted one the constraints. 67 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
Chapter 4 Componentwise H¨older inference for robust learning-based MPC 4.1 Introduction The objective of this chapter is to introduce a novel learning methodology that extends and improves KI, able to significantly reduce the prediction error. To this end, an extended version of continuity of functions is used, namely componentwise H¨older continuity. Then, leveraging on this notion, the KI is extended to exploit this property, yielding the componentwise H¨older kinky inference (CHoKI) predictor. This improved inference method could be used in the two MPCs proposed in chapters 2 and 3. In H¨older continuity, the H¨older constant describes the aggregated contribution of the variation of the inputs on the outputs. The idea of componentwise H¨older continuity is to disaggregate this contribution, considering an individual constant and exponent for the effect of each input on each output. These parameters can be described by a gain matrix and an exponent matrix. Hence, componentwise H¨older continuity can be regarded as if the scalar H¨older constant and exponent were extended to matrices. It is first proven that the proposed learning method can estimate any H¨older continuous function, by proving that componentwise H¨older continuity over a compact domain is equivalent to H¨older continuity. Besides, a non-parametric estimation method and the procedure to calculate its hyper-parameters from measured (possibly noisy) input-output data 68 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.2. EXTENSION OF H ¨ OLDER PROPERTIES is presented. It is also shown that the proposed method encompasses and extends existing methods, such as kinky inference or NSM. It is rigorously proven that the proposed learning method ensures a bounded estimation of the prediction error and that the method is a learning algorithm. This is a key aspect, since it allows one to relate the accuracy of the predictions to the density of the available observations. Besides, it renders the online version of such algorithm applicable, that is, makes it prone to decrease the prediction error as more data points are taking into consideration by the predictor. Based on the properties of the proposed estimation method, the predictive control law presented in Chapter 3 has been extended to prediction models based on CHoKI. All the results presented in this chapter encompass the standard KI approach, which could be considered as a sub-case of the more general CHoKI framework. Hence, the method proposed here is in general more precise and yields smaller prediction errors without loss of generality. The properties of the proposed estimation method could also be used to enhance existing controllers, such as [13] or [40]. The main advantage of using CHoKI as prediction method in the design of safe LBMPC, in contrast to standard KI, is two-fold. First, the prediction error bounds obtained are decreased. This is particularly important for the design of robust MPC controllers, since this reduces the size of the sets used to tighten the constraints to guarantee safety. Second, the accuracy in the estimation of the reachable sets is improved, by reducing the expansion rate of the reachable sets, leading also to a significant reduction of the conservativeness of the tightened constraints. A procedure is proposed to calculate the reachable sets and it is proven that the result yields the tightest sets that can be computed using the componentwise H¨older condition. Various case studies illustrate that CHoKI notably reduces the conservatism of the design when compared to standard KI methods, thus enhancing the performance and enlarging the domain of attraction of the proposed controller. The results in this chapter were presented in [99, 103]. 4.2 Extension of H¨older properties 4.2.1 Componentwise H¨older continuity The standard H¨older continuity condition allows one to estimate the effect on the outputs of a variation on the inputs, bounding the worst case. In this section, a new approach on the H¨older property is considered. It exploits the 69 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.2. EXTENSION OF H ¨ OLDER PROPERTIES fact that a function may vary sharply along one dimension of the input, while varying smoothly with respect to other input dimensions. Besides, if the function maps the inputs to several outputs, it seems intuitive to distinguish different parameters of the H¨older continuity for each output, as if they were different functions. Hence, we propose to use matrices as the parameters of the H¨older property: L,P ∈ Rny×nw, yielding the following definition: Definition 4.1 (Componentwise H¨older continuity).Given the matrices L and P, a function f:W → Y is componentwise L-P-H¨older continuous if ∀w1, w2∈ W and ∀i∈ I ny 1 |fi(w1)−fi(w2)| ≤ nw X j=1 Lij|w1,j −w2,j|Pij .(4.1) In H¨older continuity, Laggregates the effects of the inputs on the outputs into a single constant. On the contrary, the proposed componentwise approach uses each Li,j to take into account separately the effect of each input on each output. This componentwise H¨older continuity condition may be rewritten in a more compact form. To this end, new notation is introduced. Given a vector w∈Rnwand two matrices L,P ∈ Rny×nw, we define dP L(w) := a:ai= nw X j=1 Li,jwPi,j j,∀i∈ I ny 1!.(4.2) Then, componentwise H¨older continuity (cf. (4.1)) can be written as |f(w1)−f(w2)| ≤ dP L(|w1−w2|),∀w1, w2.(4.3) The following theorem states under which conditions H¨older continuity and componentwise H¨older continuity are equivalent. Theorem 4.1. Let f:W ⊆ Rnw→ Y ⊆ Rny. 1. If fis H¨older continuous in W, then fis componentwise H¨older continuous in W. 2. If Wis compact and fis componentwise H¨older continuous in W, then fis H¨older continuous in W. 70 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.2. EXTENSION OF H ¨ OLDER PROPERTIES Proof (point 1). Consider two parameters Land psuch that f:W → Y is H¨older continuous. Provided that for a given w∈Rnw,kwk ≤ √nwkwk∞, define L∞=Ln(p/2) w. Then, we have that ky1−y2k ≤ Lkw1−w2kp(4.4a) ≤L∞kw1−w2kp ∞(4.4b) =L∞max j(|w1,j −w2,j|)p(4.4c) ≤L∞ nw X j=1 |w1,j −w2,j|p,(4.4d) where inequality (4.4d) holds provided that |w1,j −w2,j|p≥0,∀j. Then, defining L=L∞1ny×nwand P=p1ny×nw, the function fis componentwise H¨older continuous. Proof (point 2). Given that for any w1, w2, j ∈ I nw 1,|w1,j −w2,j|≤kw1−w2k∞, |yi,1−yi,2| ≤ nw X j=1 Li,j|w1,j −w2,j|Pi,j ≤ nw X j=1 Li,jkw1−w2kPi,j ∞. If Wis a compact space, there exists a csuch that kwk∞≤c, ∀w∈ W, so kw1−w2k ≤ 2c, ∀w1, w2∈ W. Besides, if x∈[0,1] and p1≥p2, with p1, p2∈(0,1], then xp1≤xp2. Let p= min i,j Pi,j and let L= max i ny P j=1 Li,j(2c)Pi,j −p. Then Li,jkw1−w2kPi,j ∞≤ Li,j(2c)Pi,j kw1−w2k 2cPi,j ≤ Li,j(2c)Pi,j kw1−w2k∞ 2cp =Li,j(2c)Pi,j −pkw1−w2kp ∞. 71 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.2. EXTENSION OF H ¨ OLDER PROPERTIES Hence, |fi(w1)−fi(w2)| ≤ ny X j=1 Li,jkw1−w2kPi,j ∞ ≤ ny X j=1 Li,j(2c)Pi,j −pkw1−w2kp ∞ ≤Lkw1−w2kp ∞. And since ∃i∗:|fi∗(w1)−fi∗(w2)|=kf(w1)−f(w2)k∞, we have that kf(w1)−f(w2)k∞≤Lkw1−w2kp ∞. Corollary 4.1. If P=p1ny×nw, the equivalence in Theorem 4.1 holds for any space, even if it is not compact. With a slight abuse of notation, in what is to follow the mapping in (4.2) may also be used for sets, such that for a given set A⊆Rnw, dP L(A) := {dP L(x)|x∈A}.(4.5) Corollary 4.2. In the Lipschitz case and for every output component i∈ I ny 1, if the Lipschitz constant is such that L=kLik∞, then Lkwik∞≥dLi(|wi|). Proof. Using the Cauchy-Schwarz inequality, dL(|w|) = hL,|w|i ≤ kLk∞kwk∞=Lkwk∞, where h·,·i denotes the canonical inner product. 4.2.2 Componentwise H¨older kinky inference Proceeding as in the standard KI approach, assuming that fis H¨older continuous, and given a data set of input-output observations, this section presents the componentwise H¨older kinky inference (CHoKI) predictor. 72 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.2. EXTENSION OF H ¨ OLDER PROPERTIES Provided that the componentwise H¨older condition (4.1) holds in virtue of Theorem 4.1, the proposed estimation method is: ˆ fC(q; Θ,D) = 1 2min i=1,...,ND (˜yi+dP L(|q−wi|)) + max i=1,...,ND (˜yi−dP L(|q−wi|)), (4.6) where Θ = ( L,P). The resulting prediction function is componentwise H¨older continuous, with the same parameters, as stated in the following lemma. Lemma 4.1. For a given Θ = ( L,P), the CHoKI predictor ˆ fCgiven by (4.6) is componentwise L-P-H¨older continuous. Proof. Given two scalars x≥0 and 0 < p ≤1, xpis concave, so kx1kp−kx2kp≤ kx1−x2kp. Then, given L,P,dP L(|w|) is H¨older continuous, since for each i-th component nw X j=1 Li,j|w1,j|Pi,j − nw X j=1 Li,j|w2,j|Pi,j ≤ nw X j=1 Li,j(|w1,j −w2,j|)Pi,j . The sum of two componentwise H¨older continuous functions f, g with H¨older parameters Land Pis also componentwise H¨older, since |f(w1) + g(w1)−f(w2)−g(w2)| ≤ |f(w1)−f(w2)|+|g(w1)−g(w2)| ≤dP L(|w1−w2|) + dP L(|w1−w2|). Finally, the minimum (or, equivalently the maximum) of componentwise H¨older functions is also componentwise H¨older: Let denote h(w) = min(f, g) and assume, w.o.l.g., that f(w1)> g(w2). Then, if f(w1)≤g(w1), |h(w1)−h(w2)|=|f(w1)−g(w2)|=f(w1)−g(w2) ≤g(w1)−g(w2)≤dP L(|w1−w2|). Hence, ˆ fC(cf. eq. (4.6)) is componentwise H¨older continuous, with H¨older parameters Land P. In case that the parameters Land Pare unknown a priori, they must be estimated from the available input-output data. In the LACKI method, the 73 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.2. EXTENSION OF H ¨ OLDER PROPERTIES H¨older constant can be derived from Dby a procedure based on sample consistency (see Appendix B). However, in the componentwise case, the LACKI method cannot be used, since kw1−w2kprovides an aggregated measurement of the effect of the inputs on the outputs, and there is no direct information of the contribution of each input on a particular output. In order to infer this contribution, an optimization method is adopted, extending the results of POKI to obtain the matrices Land P. The estimation of the parameters Θ = ( L,P) is based on solving an optimization problem offline, which depends on a regularization parameter η∈Rnyand two data sets: the data set D, used for estimation, and a data set Dtest used for validation. Thus, and hybrid method that merges LACKI and POKI is proposed, to compute the parameters L,Pas solutions of the following optimization problem, where a measure of the performance of the prediction over the data set Dtest,g(Θ,D,Dtest), and a regularization term are minimized, subject to a constraint that ensures their consistency with the samples of the data set: Θ = arg min Θg(Θ,D,Dtest) + τLkL− L0k1,(4.7a) s.t.|˜yi−˜yj|−η≤dP L(|wi−wj|) ∀wi, wj∈ WD, wi6=wj,(4.7b) 0<Pij ≤1, i ∈ I ny 1, j ∈ I nw 1,(4.7c) where τLis a design regularization hyper-parameter used to ensure boundedness of L, and L0stands for any possible prior guess of L. Besides, the cost function gmust be positive and bounded, for any size of D. A possible choice of the performance measure is the mean squared prediction error: g(Θ,D,Dtest) = 1 NDX wi∈WDtest kˆ fC(wi; Θ,D)−˜yik2.(4.8) Remark 4.1. For an analysis on the effect of the regularization hyperparameter η, the reader is referred to [36], where η≥2¯ eis taken for each component to effectively smooth out the effect of the noise in the prediction. Remark 4.2. The regularization term of the cost function, τLkL − L0k1, prevents the problem from overfitting the noise, while ensuring boundedness of L. If Dtest is separate from Dand such that gis bounded for all L, then it can be removed, setting τL= 0. Remark 4.3. In practice, it may be easier to fix Pa priori, and to optimize over L, provided that the assumptions hold for the chosen P. 74 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.2. EXTENSION OF H ¨ OLDER PROPERTIES Next, based on T. Mitchell’s definition of a learning algorithm [113], it will be proven that CHoKI is a learning method. For this, the following lemma, which proves sample consistency, is needed. Lemma 4.2 (Sample consistency of CHoKI).If Land Pare obtained as in (4.7), the CHoKI predictor (4.6) is sample-consistent (up to η/2+¯ e), that is, |ˆ fC(wk)−f(wk)| ≤ η 2+¯ e,∀wk∈ WD.(4.9) Proof.Let us denote i= arg min n˜yn+dP L(|wk−wn|), j= arg max n˜yn−dP L(|wk−wn|). Then, ˆ fC(wk) = 1 2˜yi+dP L(|wk−wi|) | {z } B +1 2˜yj−dP L(|wk−wj|) | {z } A =1 2B+1 2A. It is first proven that A≥˜yk. If j=kthis is immediate. Otherwise, note that A= ˜yj−dP L(|wk−wj|) ≥˜yk−dP L(|wk−wk|) = ˜yk. Then, it is proven that A≤˜yk+η. From (4.7b) we have that |˜yk−˜yj| ≤ dP L(|wk−wj|) + η. Provided that ˜yj≥A≥˜yk, A= ˜yj−dP L(|wk−wj|) ≤˜yk+|˜yj−˜yk|−dP L(|wk−wj|)≤˜yk+η. The same procedure is applied to prove that ˜yk−η≤B≤˜yk. Then, ˆ fC(wk) is such that 1 2(˜yk−η) + 1 2˜yk≤ˆ fC(wk)≤1 2(˜yk+η) + 1 2˜yk, or equivalently, ˜yk−η 2≤ˆ fC(wk)≤˜yk+η 2,|ˆ fC(wk)−˜yk| ≤ η 2. 75 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.3. REACHABILITY SETS Proof. Provided that w(k)=(x(k), u(k)), and that the sequence of future inputs u(k+i) for i∈ I N−1 0is given, |ˆw(j|k)−ˆw(j−1|k+ 1)|= |ˆy(j|k)−ˆy(j−1|k+ 1)| |ˆy(j−1|k)−ˆy(j−2|k+ 1)| . . . |ˆy(j−na|k)−ˆy(j−na−1|k+ 1)| |u(k+j−1) −u(k+j−1)| . . . |u(k+j−nb)−u(k+j−nb)| |u(k+j)−u(k+j)| .(4.20) These ˆyare predicted values if j > na, and real measurements if not. Equation (4.20) directly translates into (4.15b), provided that if ˆy(j|k) is a real measurement, |y(j|k)−y(j−1|k+ 1)|equals 0 ny, and belongs to Mj otherwise. To obtain Mjwe make use of the componentwise H¨older continuity of the predictor ˆ fC. Given equation (4.13), we have that M1= B (c1). Note that ˆy(j+ 1|k) = ˆ fC(ˆx(j|k), u(k+j); Θ,D), and that |ˆ fC,i −ˆ fC,j| ≤ dP L(|wi−wj|). Hence, |∆y|=|ˆy(j+ 1|k)−ˆy(j|k+ 1)| =|ˆ fC(ˆx(j|k), u(k+j)) −ˆ fC(ˆx(j−1|k), u(k+j))| ≤dP L |ˆy(j|k)−ˆy(j−1|k+ 1)| |ˆy(j−1|k)−ˆy(j−2|k+ 1)| . . . |ˆy(j−na|k)−ˆy(j−na−1|k+ 1)| |u(k+j−1) −u(k+j−1)| . . . |u(k+j−nb)−u(k+j−nb)| |u(k+j)−u(k+j)|) . From this inequality follows that |∆y| ≤ dP L(Gj), so the minimum set that contains all possible values of |∆y|is the Cartesian closed topological hull, that is, Mj+1 =B(dP L(Gj)). 82 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.4. STABILIZING CONDITIONS FOR CHOKI-MPC If the set Yis a polytope, as it is customary, then the resulting tightened constraints are also polytopes. Notice that these calculations are done only once, offline. 4.4 Stabilizing conditions for CHoKI-MPC In order to recover the safe-by-design properties of the controller (3.8) based on KI, some of the assumptions are extended to the componentwise framework as follows: Assumption 4.1. The stage cost function `(x, u)is a componentwise H¨older continuous, positive definite function such that `(x, u)≥α(kxk), where αis aK-function. Given P, its H¨older constant is L`. Assumption 4.2. Vfis a componentwise H¨older continuous, positive definite function, with H¨older parameters Lf,P, such that Assumption 3.3 holds. Assumption 4.3. Assumption 3.4 holds redefining νas ν(µ) = N X j=1 kdP L`(Gj)k∞+λkdP Lf(GN+1)k∞.(4.21) Next, we present the stability result for the CHoKI-MPC, which follows the same line of reasoning as the stability proof presented in Chapter 3 for KI predictors. Theorem 4.3 (ISS stability).Suppose that assumptions 3.1, 4.1-4.3 hold for the optimization problem PN(·)obtained with the CHoKI predictor. Let κN(x) be the control law derived from the solution of PN(x)applied using a receding horizon policy. Then, for any x(0) ∈Γ, the system to be controlled by the control law u(k) = κN(x(k)) is input-to-state stable with respect to the estimation error, x(k)∈Γ, and the constraints are always satisfied, i.e. y(k)∈ Y,∀k. Proof. The proof of this theorem follows the same steps as the proof of Theorem 3.1 (see page 61). Lemma 3.2 is modified as follows: Lemma 4.6. For all y∈ Yjand for all ∆ysuch that |∆y|∈Mj, the sets Yjare such that y+ ∆y∈ Yj−1. 83 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.5. CASE STUDIES Proof By definition, we have that if |∆y| ∈ Mjthen ∆y∈ Rj. Thus, since the origin is contained in Rj, y+ ∆y∈ Yj⊕Rj=Yj−1Rj⊕Rj⊆ Yj−1. Corollary 3.1 can be easily extended to analyse this theorem with respect to the validity of Assumption 3.1. As stated in Remark 3.6, if the prediction error is larger than the estimated µfor a period of time, Theorem 4.3 still holds if x(k)∈Γ for that period of time. Furthermore, if Assumption 3.1 does not hold in general, we propose to have a backup controller, similar in spirit to the approach proposed by L. Hewing et al. [61]. Specifically, the controller proposed in Chapter 2 admits the extension to the CHoKI-based model, and it guarantees stability relaxing the hard constraints in the outputs to soft constraints included in the cost of the MPC. It is important to remark that the benefits of the robust predictive controller based on the componentwise H¨older approach are two-fold, compared to the KI-based MPC of Chapter 3. First, the enhanced learning method potentially provides less conservative estimation errors (see Fig. 4.1 on page 79), yielding more accurate predictions. This leads to an improvement of the closed-loop performance of the controlled system. Second, note that the recursion used for the set of tightened constraints (cf. equation (4.12)) is obtained using the componentwise H¨older metric, in contrast to the standard version presented in Section 3.3. Hence, recalling Corollary 4.2 (on page 72), even if the maximum prediction errors were the same (which in general they are not), the back-off of the set of tightened constraints are less conservative (RCHoKI N⊆ RKI N), thus yielding larger regions of feasibility of the proposed controller. This double benefit of the proposed method will be illustrated in the following section, using two case studies. 4.5 Case studies This section considers a couple of case studies to illustrate the performance of the proposed method. 84 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.5. CASE STUDIES h1h2 h3 h4 Tank 1 Tank 2 Tank 3 Tank 4 γaγb qaqb q1 q2 q3q4 Figure 4.3: Quadruple-tank process scheme. 4.5.1 Quadruple tank process The first system considered is the quadruple-tank process described in [9, 68]. The process consists of four tanks, where the two on top discharge on the inferior ones. The tanks are fed with two pumps, which send two flows represented by qaand qb. These flows enter two three-ways valves, which divide each flow into two branches, determined by the fractions γaand γb. Thus, qaγagoes into the tank number one, qa(1−γa) into tank four, qbγbinto the tank number two and qb(1 −γb) into tank number three. Tank number three discharges into tank number one, and the fourth one into the second one, as represented in Figure 4.3. There are two control inputs, the flows qaand qb(m3h−1). The heights of the tanks are denoted hi(m), i ∈ I 4 1. The outputs of the system are the heights of the two lower tanks, i.e., h1and h2. The dynamics of the plant are nonlinear, and are well-modelled by the following set of differential equations: A1 dh1(t) dt =−a1p2gh1(t) + a3p2gh3(t) + γa qa(t) 3600,(4.22a) A2 dh2(t) dt =−a2p2gh2(t) + a4p2gh4(t) + γb qb(t) 3600,(4.22b) A3 dh3(t) dt =−a3p2gh3(t) + (1 −γb)qb(t) 3600,(4.22c) 85 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.5. CASE STUDIES Param. Definition Value Units AArea of the four tanks 0.03 m2 a1Eq. area of the hole of tank 1 1.31 ×10−4m2 a2Eq. area of the hole of tank 2 1.51 ×10−4m2 a3Eq. area of the hole of tank 3 9.27 ×10−5m2 a4Eq. area of the hole of tank 4 8.82 ×10−5m2 γaFraction of three-ways valve a0.3 - γbFraction of three-ways valve b0.4 - gGravity acceleration 9.8 m s−2 Table 4.1: Parameters of the quadruple tank process. A4 dh4(t) dt =−a4p2gh4(t) + (1 −γa)qa(t) 3600,(4.22d) where Ai(m2) denotes the area of tank iand ai(m2) the equivalent area of the hole of tank i. The parameters of the model are given in Table 4.1. It is also assumed that the height sensors have a 2% measuring error. The error is generated randomly for each measurement using an uniform distribution. The constraints in the inputs are 0 ≤qa,b ≤2.6 m3h−1, and the constraints in the heights are given by 0 ≤h1≤1.25 m and 0 ≤h2≤1.42 m. Note that this is a really interesting problem from the academical control point of view, since it is a MIMO system with high coupling, which induces non-minimal phase zeros. First, the set of equilibrium points is estimated using a set of steps in the inputs from qmin to qmax with increments of 0.1 m3h−1, where each step is long enough to reach a steady sate. The result is shown in Figure 4.4. In addition to obtaining the equilibrium points of the system, this test is used to adjust the sampling time, which is set to τs/5 =5 s, where τsstands for the mean settling time of the sequence of steps applied. An equilibrium point is chosen to be the reference operating point, i.e. href = [0.65 0.65](m), and qref = [1.63 1.99](m3h−1). After defining the set of equilibrium points, a set of experiments is carried out to obtain the data sets. The experiments are designed using the methodologies presented in Chapter 2 to identify the dynamics of a system within a workspace: a sequence of chirp signals covering the workspace are applied to generate the raw data set containing the trajectories of heights and flows, Draw, with ND= 36000. In addition, several tests with random input 86 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.5. CASE STUDIES Figure 4.4: Set of equilibrium points of the quadruple tank process. 0 50 100 0 0.5 1 0 50 100 0 0.5 1 0 50 100 0 1 2 0 50 100 0 1 2 Figure 4.5: Part of the PRS applied to obtain data sets for cross-validation. signals are carried out in order to obtain data sets for cross-validation, Dtest, (see Figure 4.5). The values of all the signals are scaled between 0 and 1. The regression vector w(k) is constructed for different values of naand nb. The prediction horizon is set to N= 5. Cross validation tests are used to estimate the parameter L, such that the sum of the areas of Yjis maximized, i.e., setting g(·,·,·) in equation (4.7) (on page 74) to g(Θ,D,Dtest) = − N X j=1 kYjk,(4.23) with τL= 0 and η= 0.04 12. The minimum is obtained for na= 3 87 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.5. CASE STUDIES and nb= 1, yielding µ= [0.031 0.032]. 1The set of tightened constraints is obtained considering the Lipschitz case. The resulting Lis L=2 0.8 5.3 1.3 3.7 1.2 7.9 2.4 5.8 1.2 2.5 4.3 0.4 1.9 1.561.3 4.4 1.9 7.6 1.9 6.2 2.9 3.5. The solver chosen for the optimization problem is MATLAB ® ’s fmincon. The stage and terminal costs are defined as `(x, u) = kx−xrefk2 Q+ku−urefk2 R Vf(x) = kx−xrefk2 P, where Q= 100I10,R=I2, and λ= 10. The terminal matrix Pis obtained solving an LQR for a linearisation of the CHoKI model around the reference point. The initial state is set to hi= 0.45 m, i ∈ I 4 1, and the proposed CHoKI- based controller is applied to the system. In order to compare the results, the same setup is applied to other two controllers, whose models are (i) the ideal state-feedback set of ODEs and (ii) the KI used in Chapter 3. The results are shown in Figure 4.6, for 100 simulations subject to random noise. Besides, the performance of these simulations is measured according to Φ = tsim X i=1 `(x(i), u(i)),(4.24) which is compared in Figure 4.7. Note that the data-based control problem is able to perform in a similar way to the ideal MPC, whereas the standard KI exhibits a worse performance, which illustrates the main properties of CHoKI. Note that the main advantage of CHoKI with respect to standard KI is not only the improvement of the prediction, which in general leads to better closed-loop performance results, but also the enlargement of the tightened constraints presented in Section 4.3, which implies that the multivariate bound is less conservative. If we used the obtained Din the standard KI approach presented in Section 3.3, the prediction error would be increased up to µKI = [0.087 0.088] 1(recall that µCHoKI = [0.031 0.032]). Then, the maximum prediction horizon such that Assumption 3.2 holds, i.e., that YNis not empty, is NKI = 2. The sets Yjand Rjare represented in Figure 4.8, 1 both for the CHoKI and the KI approaches. Note that kRCHoKI 5k= 0.014 while kRKI 5k= 54.96. 1 1Recall that the signals are scaled between 0 and 1. 88 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.5. CASE STUDIES 01234 0.4 0.6 0.8 01234 0.4 0.6 0.8 01234 0.4 0.6 0.8 01234 0.4 0.6 0.8 01234 0.4 0.6 0.8 01234 0.4 0.6 0.8 Figure 4.6: A hundred closed-loop simulations of the quadruple-tank process, with three MPCs whose models are (top row) the ideal state-feedback set of ODEs, (middle row) the standard KI and (bottom row) the proposed CHoKI- based MPC. 2000 4000 6000 8000 Figure 4.7: Performance index comparison of the three controllers of Fig. 4.6. 89 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.5. CASE STUDIES (a) Yjwith KI (b) Yjwith CHoKI (c) Rjwith KI (d) Rjwith CHoKI Figure 4.8: Sets Yjand Rj, both for the KI and the CHoKI approaches 1in the quadruple-tank case. Note that Yjis empty for j≥3 in the KI setup, and notice the different scales between the sets Rjin figures (c) and (d). 90 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
4.5. CASE STUDIES 4.5.2 Isothermal chemical reactor An isothermal chemical reactor (ICR) is considered next. It was studied in [132]. The following consecutive-competitive reactions take place: P0+B→P1 P1+B→P2. The inputs considered are the flows of P0and B, denoted qaand qb, respectively; and the outputs are the measured concentrations of P1and P2. The dimensionless equations that model this system are dP0(t) dt =qa−P0−k1P0B(4.25a) dB(t) dt =qb−B−k1P0B−k2BP1(4.25b) dP1(t) dt =−P1+k1P0B−k2BP1(4.25c) dP2(t) dt =−P2+k2BP1,(4.25d) in which the parameters k1and k2are set to 1 and 0.4, respectively. The constraints in the inputs are given by 5 ≤qa, qb≤10, and in the outputs by 1.59 ≤P1≤3.82 and 0.71 ≤P2≤2.50. The reference point is given by yref = [3.09 1.66] and uref = [8.28 7.76]. We proceed the same way as in the previous case study, obtaining the data sets, the regressors, and the prediction error, which is minimized for na= 2 and nb= 0. The costs are also defined the same way, with Q= 100I6,R= 10I2and λ= 100. The resulting MPC, applied to 100 simulations, is shown in Figure 4.9(a). Noticing the overshoot in the concentration of P1, a second set of experiments is carried out, in order to demonstrate the ability of the proposed controller to maintain robust satisfaction of the hard constraints. To this end, the upper constraint of P1is lowered to from 3.82 to 3.3, and the same controller is applied. Note that the MPC is still able to satisfy the constraint. 91 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
5.3. ONLINE LEARNING-BASED MPC Remark 5.2. Equation (5.10b) in Assumption 5.2 requires that the Lyapunov condition holds for the online prediction model. Since the terminal ingredients must be computed offline, the online model is not available for this design. Taking into account that the safe model can be regarded as a sampled subset of the online model, since D(0) ⊆ D(k), the online model is contained in a difference inclusion centered in the safe model and given by the estimation error of the LACKI method, which is explicitly known. Therefore, robust design methods can be used to calculate the terminal ingredients, satisfying Assumption 5.2. Definition 5.1. Making use of the definition of cj(c1)and rj(c1)in (5.6) and of Assumption 5.2, the following function, which is a K-function for every input, is defined: θa(c1) = N−1 X j=0 θ`(kcj+1(c1)k) + θf(rN+1(c1)).(5.11) Assumption 5.3. The set Υ = {x:`(Mx, 0) ≤θa(2µ)}(5.12) is contained in Ωγ. Definition 5.2. φis a positive constant such that `(Mx, 0) > φ, ∀x /∈Ωγ.(5.13) Similarly to Chapter 3, it can be proven that φ≥θa(2µ). The region of feasibility of the problem (5.7) is denoted XN, and a level set of its optimal cost is defined as: Γ = {x∈XN:V∗ N(x)≤Nφ +λγ}.(5.14) The control law is given in Algorithm 1, where the updating policy of the online data set guarantees closed-loop robust stability. At each time step k, the online data set may be updated with the current data point, that is, D(k) = D(k−1)∪(y(k), w(k−1)), subject to the following triple updating criteria: 1. The current data point must not be close to any point already contained in D(k−1), for a given threshold τt≥0. This rule prevents the 98 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
5.3. ONLINE LEARNING-BASED MPC cardinality of Dfrom becoming large, yielding the method computationally expensive. The policy is given by the Boolean function: πexp(w, D) = min wi∈WDkw−wik> τt. A reasonable estimate for τtis the noise level ¯ e. 2. Adding the new data point must not increase the cost, calculated with the shifted input sequence which, together with the control policy and update rule, is defined in Algorithm 1. πV(y, w, D) = ¯ VN(x, ¯ u;¯ D)≤¯ VN(x, ¯ u;D). 3. The new candidate (y(k), w(k−1)) must be consistent with the Lipschitz constant. This rule is optional, included so that the online model maintains the properties of the LACKI method (see Appendix B). πL(y, w, L, D) = max i∈ I ND 1 |y−yi|−2¯ e kw−wik≤L!. Hence, the aggregated update rule is the Boolean function π(y, w, L, D) = πexp ∧πV∧πL.(5.15) In the following theorem it is proven that the proposed control algorithm guarantees that the closed-loop system is input-to-state stable (ISS). Theorem 5.1. Consider that Assumptions 5.1-5.3 hold, and let κMPC(x)be the control law derived from the solution of PN(k)applied using Algorithm 1. Then, for any feasible state x(0) ∈Γ, the system controlled by the control law u(k) = κMPC(x(k)) is input-to-state stable w.r.t. the estimation error don(k), and the constraints are fulfilled along the operation, i.e., y(k)∈ Y and u(k)∈ U,∀k. Proof. In this proof, different predicted trajectories will be considered, as shown in Figure 5.1. These trajectories are obtained predicting from x(k) with u∗(k), or from x(k+ 1) with ¯ u(k+ 1), and using the safe or the online models. It is first proven that Γ is an invariant set of the closed-loop system, i.e., if x(k)∈Γ, then x(k+ 1) ∈Γ. To this end, as it is standard in MPC proofs, 99 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
5.3. ONLINE LEARNING-BASED MPC Algorithm 1 Online update and control law. while automatic control is on do Read y(k) x(k)←(y(k), . . . , y(k−na), u(k−1), . . . , u(k−nb)) .State function Shifted Sequence(¯ u).Shifted Sequence ¯u(j)←u∗(j), j ∈ I N−2 0 ¯u(N)←κf(ˆxon(N|k−1)) return ¯ u end function ¯ D(k)← D(k−1) ∪(y(k), w(k−1)) .Updated data set if π(y(k), w(k−1), L, D(k−1)) then D(k)←¯ D(k).Update else D(k)← D(k−1) .Do not update end if u∗(k)←κMPC(x(k); L, D(0),D(k)) .Solution of PN(k) u(k)←u∗(0) .Receding horizon Apply u(k) to the system end while it is proven that the shifted trajectory presented in Algorithm 1 is feasible for x(k+ 1). Provided that x(k)∈Γ, due to its definition in equation (5.14), it follows that V∗ N(x(k)) ≤Nφ +λγ, (5.16) and therefore it can be proven [83] that ˆxon(N|k)∈Ωγ.(5.17) Feasibility of the shifted sequence ¯ u(k+ 1) is proven for any x(k)∈Γ, showing that (i) ¯u(k+j+ 1) ∈ U,∀j∈ I N−1 0. (ii) ˆys(j|k+ 1) ∈ Yj,∀j∈ I N−1 0. Constraint (i) holds because u∗(k)∈ U, and since ˆxon(N|k)∈Ωγ, because of equation (5.9a), the shifted sequence defined in Algorithm 1 is feasible, i.e., ¯u(k+j+ 1) ∈ U,∀j∈ I N−1 0. To address (ii), the following lemmas are used: 100 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
5.3. ONLINE LEARNING-BASED MPC y(k)y(k+ 1) y(k+j) ˆyon(1|k) ˆyon(j|k) ˆys(1|k) ˆys(j|k) ˆys(j−1|k+ 1) aj cj Figure 5.1: Propagation of the predictions. Lemma 5.1. The following inequalities hold: |ˆyon(j−1|k+ 1) −ˆyon(j|k)| ≤ cj(don(k)),(5.18a) kˆxon(j−1|k+ 1) −ˆxon(j|k)k ≤ rj(don(k)),(5.18b) where cj(c1)and rj(c1)are obtained from the recursion (5.6). The bounds (5.18) hold for the safe model, taking c1=µinstead of don(k). Proof.Given the triangular inequality that applies to the definition of the state vector x, kx(k+1)−ˆxon(1|k)k≤ky(k+1)−ˆyon(1|k)k=kdon(k)k=kc1k=r1.(5.19) Assuming that rj−1(c1) and cj−1(c1) are known, we can derive, from the Lipschitz continuity of ˆ fon, that |ˆyon(j−1|k+ 1) −ˆyon(j|k)| =|ˆ fon(ˆxon(j−2|k+ 1), u(k+j−1)) −ˆ fon(ˆxon(j−1|k), u(k+j−1))| ≤Lkˆxon(j−2|k+ 1) −ˆxon(j−1|k)k ≤Lrj−1=cj.(5.20) The error of the estimated state can be bounded as follows, given that ˆxon(j|k) 101 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
5.3. ONLINE LEARNING-BASED MPC contains min(na+1, j) estimated outputs (and real measurements otherwise) kˆxon(j−1|k+ 1) −ˆxon(j|k)k ≤ kˆyon(j−1|k+ 1) −ˆyon(j|k)k+. . . +ky(k+ 1) −ˆyon(j|k)k ≤ j X i=σjkcik=rj.(5.21) Note that the same holds if the safe model is used, with ˆxsand ˆys, being c1=µ. Lemma 5.2. For all y∈ Yjand all ∆y∈ B(cj(µ)), the sets Yjare such that y+ ∆y∈ Yj−1. Proof.First, it is proven that B(δ1)⊕ B(δ2)⊆ B(δ1+δ2). Indeed, for all y=y1+y2, with y1∈ B(δ1) and y2∈ B(δ2), we have that kyk ≤ ky1k+ky2k ≤ δ1+δ2,(5.22) and thus y∈ B(δ1+δ2). Since (for any c1)dj=cj+dj−1, it is inferred that B(dj−1)⊕B(cj)⊆ B(dj), and then, Yj=Y B(dj(µ)) ⊆ Y B(dj−1(µ)) B(cj(µ)) = Yj−1B(cj(µ)).(5.23) Given that ∆y∈ B(cj(µ)) and y∈ Yj, it is derived that y+ ∆y∈ Yj⊕B(cj(µ)) ⊆ Yj−1B(cj(µ)) ⊕B(cj(µ)) ⊆ Yj−1.(5.24) It is known that ˆys(j|k)∈ Yjup to N−1, and that |ˆys(j|k)−ˆys(j−1|k+ 1)| ≤ cj(µ).(5.25) Hence, because of Lemma 5.2, ˆys(j|k+ 1) ∈ˆys(j+ 1|k)⊕B(cj+1(µ)) ⊆ Yj+1 ⊕B(cj+1(µ)) ⊆ Yj. For the last prediction, the fact that ˆxon(N|k)∈Ωγimplies that the last predicted output ˆyon(N|k)∈ YN(see equation (5.9b)). Another lemma is introduced: 102 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
5.3. ONLINE LEARNING-BASED MPC Lemma 5.3. The following inequalities hold: |ˆyon(j|k)−ˆys(j|k)| ≤ aj,(5.26) where ajis given by the recursion aj+1 =Lkajk+µ, (5.27) with a1=µ. Proof.This lemma is proven by recursion. For j= 1 the condition holds, given the definition of ˆ fon in (5.4): |ˆyon(1|k)−ˆys(1|k)| ≤ a1=µ. (5.28) Then, assuming that aj−1is known, we can obtain. |ˆys(j|k)−ˆyon(j|k)| =|ˆ fs(ˆxs(j−1|k), u(k+j−1)) −ˆ fon(ˆxon(j−1|k), u(k+j−1))| (5.29a) ≤ |ˆ fs(ˆxs(j−1|k), u(k+j−1)) −ˆ fs(ˆxon(j−1|k), u(k+j−1))| +|ˆ fs(ˆxon(j−1|k), u(k+j−1)) −ˆ fon(ˆxon(j−1|k), u(k+j−1))| (5.29b) ≤Lkˆys(j−1|k)−ˆyon(j−1|k)k+µ(5.29c) =Lkaj−1k+µ=aj.(5.29d) Given the definition of cj(µ) in Lemma 5.1 and ajin Lemma 5.3, the following bound is obtained: |ˆyon(j|k)−ˆys(j−1|k+ 1)| ≤ |ˆyon(j|k)−ˆys(j|k)|+|ˆys(j|k)−ˆys(j−1|k+ 1)| ≤aj+cj(µ).(5.30) Using Lemma 5.3 and equation (5.30) it follows that the error between the predictions of the safe and the online models is bounded, so |ˆyon(N|k)−ˆys(N−1|k+ 1)| ≤ aN+cN(µ). 103 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
5.3. ONLINE LEARNING-BASED MPC Then, because of equation (5.9b), ˆys(N−1|k+1) ∈Mˆxon(N|k)⊕B(aN)⊕B(cN(µ)) ⊆ YN⊕B(cN(µ)) ⊆ YN−1, (5.31) which completes the proof of (ii). Given the definition of the cost in (5.7a), the measurements of y(k) and y(k+ 1), and the optimal and shifted sequences u∗(k) and ¯ u(k+ 1) (see Algorithm 1), the following equality holds: VN(ˆxon(1|k),¯ u(k+ 1); D(k)) −V∗ N(x(k),u∗(k); D(k)) =`(ˆyon(N|k), κf(ˆxon(N|k))) + Vfˆ Fon (ˆxon(N|k), κf(ˆxon(N|k))) −`(y(k), u∗(k)) −Vf(ˆxon(N|k)).(5.32) Since ˆxon(N|k)∈Ωγand taking into account (5.10b) the following inequality is obtained: Vfˆ Fon (ˆxon(N|k), κf(ˆxon(N|k)))−Vf(ˆxon(N|k)) +`(ˆyon(N|k), κf(ˆxon(N|k))) ≤0,(5.33) which implies that VN(ˆxon(1|k),¯ u(k+ 1); D(k)) −V∗ N(x(k),u∗(k); D(k)) ≤ −`(y(k), u∗(k)). (5.34) In addition, VN(x(k+ 1),¯ u(k+ 1); D(k)) −VN(ˆxon(1|k),¯ u(k+ 1); D(k)) = N−1 X j=0 `(ˆyon(j|k+ 1),¯u(k+j+ 1)) − N−1 X j=0 `(ˆyon(j+ 1|k),¯u(k+j+ 1)) +Vf(ˆxon(N|k)) −Vfˆ Fon (ˆxon(N|k), κf(ˆxon(N|k))).(5.35) Using Lemma 5.1 and taking into account continuity of `and Vf, there exists certain K-functions θ`and θfand a function θasuch that VN(x(k+ 1),¯ u(k+ 1); D(k)) −VN(ˆxon(1|k),¯ u(k+ 1); D(k)) ≤ N−1 X j=0 θ`(kcj+1(don(k))k) + θf(rN+1(don(k))) = θa(don(k)), (5.36) 104 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
5.3. ONLINE LEARNING-BASED MPC Now, merging (5.34) and (5.36) yields VN(x(k+ 1),¯ u(k+ 1); D(k)) −V∗ N(x(k),u∗(k); D(k)) ≤ −`(y(k), u∗(k)) + θa(don(k)).(5.37) Moreover, because of the updating policy of Algorithm 1, VN(x(k+ 1),¯ u(k+ 1); D(k+ 1)) ≤VN(x(k+ 1),¯ u(k+ 1); D(k)) ,(5.38) and by optimality, V∗ N(x(k+ 1); D(k+ 1)) ≤¯ VN(x(k+ 1),¯ u(k+ 1); D(k+ 1)).(5.39) Summing up, V∗ N(x(k+ 1),u∗(k+ 1); D(k+ 1))−V∗ N(x(k),u∗(k); D(k)) ≤ −`(y(k), u∗(k)) + θa(don(k)).(5.40) Next, recursive feasibility is proven ensuring that x(k+1) ∈Γ, for which it is necessary that V∗ N(x(k+1); D(k+1)) ≤Nφ+λγ. Consider that x(k)/∈Υ. Then, `(y(k), u∗(k)) > θa(2µ), provided that the maximum online prediction error satisfies don ≤2µ, given the definition of ˆ fon in (5.4). Therefore, using (5.16) and (5.40) it is derived that x(k+ 1) ∈Γ. If x(k)∈Υ, then x(k)∈Ωγ(Assumption 5.3). From standard arguments of MPC it follows that V∗ N(x(k)) ≤λVf(x(k)) ≤λγ, (5.41) which taking into account θa(don(k)) ≤φimplies that V∗ N(x(k+1),u∗(k+1); D(k+1)) ≤φ+λγ −`(y(k), u(k)) ≤Nφ+λγ, (5.42) and hence x(k+ 1) ∈Γ. To prove input-to-state stability, (5.40) and the continuity of the stage cost is used to derive that V∗ N(x(k+ 1); D(k+1)) −V∗ N(x(k); D(k)) ≤ −αy(ky(k)k)−αu(ku(k)k) + θa(don(k)).(5.43) Defining W(x(k); D(k)) = n X j=0 V∗ N(x(k−j); D(k−j)),(5.44) 105 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
5.3. ONLINE LEARNING-BASED MPC with n= max(na, nb+ 1), it follows that αx(kx(k)k)≤ n X j=0 αy(ky(k−j)k) + αu(ku(k−j−1)k), as in Chapter 2, and hence W(x(k+ 1);D(k+ 1)) −W(x(k); D(k)) ≤ −αx(kx(k)k)+(n+ 1)θamax j=0,...,n don(k−j),(5.45a) α1(kxk)≤W(x(k)) ≤α2(kxk).(5.45b) Hence, W(x(k)) is an ISS Lyapunov function. The proposed MPC has been proven to be stable and to robustly satisfy the output constraints under certain assumptions. Note that it lacks a terminal constraint, avoiding the calculation of a robust invariant set, unlike what is common in the design of robust MPCs [133]. The performance of the controller is enhanced by the inclusion of fresh data from the operation of the plant, as it will be illustrated in the case study of the following section. Not only does the behaviour improve, but also the convergence rates are decreased as the updated estimation error decreases, because the closed-loop system is proven to be input-to-state stable with respect to the estimation error of the online model. This error decreases in average with time, for increasing data sets. If the terminal cost is designed offline with the safe model, the following corollary proves that input-to-state practical stability can still be derived: Corollary 5.1. Consider that Assumptions 5.1 and 5.3 hold, and that Assumption 5.2 holds replacing equation (5.10b) by Vf(ˆ Fs(x, κf(x))) −Vf(x)≤ −`(Mx, κf(x)),(5.46) that is, satisfying the Lyapunov condition for the safe model. Let κMPC(x)be the control law derived from the solution of PN(k)applied using Algorithm 1. Then, for any feasible state x(0) ∈Γ, the system controlled by the control law u(k) = κMPC(x(k)) is input-to-state practically stable w.r.t. don(k)and µ. Besides the constraints are fulfilled along the operation, i.e., y(k)∈ Y. 106 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
5.3. ONLINE LEARNING-BASED MPC Proof. The proof is similar to the one of Theorem 5.1, but replacing (5.33) as follows. First, see that from the uniform continuity of Vfand Lemma 5.3 (equation (5.28)), we have that Vfˆ Fon (ˆxon(N|k), κf(ˆxon(N|k)))−Vfˆ Fs(ˆxon(N|k), κf(ˆxon(N|k))) is upper bounded by θf(kˆ Fon (ˆxon(N|k), κf(ˆxon(N|k))) −ˆ Fs(ˆxon(N|k), κf(ˆxon(N|k)) k), which is equal to θf(kˆ fon (ˆxon(N|k), κf(ˆxon(N|k))) −ˆ fs(ˆxon(N|k), κf(ˆxon(N|k)) k)≤θf(kµk). Taking this into account, (5.33) can be rewritten in this case as Vfˆ Fon (ˆxon(N|k), κf(ˆxon(N|k)))−Vf(ˆxon(N|k)) +`(ˆyon(N|k), κf(ˆxon(N|k))) ≤θf(kµk),(5.47) Following the subsequent steps, (5.45a) would be rewritten as follows, given n= max(na, nb+ 1): W(x(k+ 1); D(k+ 1)) −W(x(k); D(k)) ≤ −αx(kx(k)k)+(n+ 1)θamax j=0,...,n don(k−j)+θf(kµk).(5.48) Consequently, the controlled system would be ISS w.r.t. don and µ. Remark 5.3. The domain of attraction of the controller is defined by the feasibility region Γin (5.14). This set increases as the weighting factor λ increases, as proven in Chapter 2. Hence, λcan be chosen arbitrarily big in order to enlarge the domain of attraction of the MPC. Remark 5.4. Soft constraints in the outputs, y(k)∈ Ysoft, may also be considered, so that the performance will be penalized if the outputs go beyond Ysoft. To account for the soft constraint, we will make use of barrier functions, adding a penalizing term to the stage cost. Hence, we can describe the stage cost `(y, u)as the standard cost to track the reference, `t, plus the barrier function `b, as in Chapter 2: `(y, u) = `t(y, u) + `b(y),(5.49) where `b(y)=Λd(y, Ysoft).Λis a large constant, and d(y, Ysoft)is a measurement of the distance of yto the set Ysoft. 107 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
6.2. EXPLORATION-EXPLOITATION APPROACH The results presented in this chapter have been published in [95]. 6.2 Exploration-exploitation approach 6.2.1 Exploration The objective is the design of predictive controllers able to control the system in regions with low data density. In such areas, the prediction error increases rapidly, probably exceeding the robustness bound that the controller can afford. This problem can be tamed considering an exploration technique in which the control strategy forces the system to stay close to observed data, in a safe region where we have enough information to guarantee a worst case upper bound of the prediction error. This safe region is defined as Wr=w: min(kw−wik)≤τr,∀i∈ I ND 1,(6.1) for certain threshold τr≥0. Using KI with the Lipschitz case, the following property allows the method to relate the exploration distance τrwith the estimation error bound. Property 6.1. The prediction error µis bounded by µ=Lfτr+ 2¯ e.(6.2) Note that this bound is based on the true Lipschitz constant Lf, which in general is unknown. In this chapter we assume that this constant is known, equal to the estimated constant L. Assumption 6.1. L=Lf. The veracity of Assumption 6.1 depends on the density of the data set, and it conditions the validity of the results presented in this chapter. Previous works on the estimation of the Lipschitz constant [35] provided a Pareto probability distribution on the Bayesian estimation of Lf. Besides, if the probability of an underestimation is bounded by ρ, i.e., P(Lf> L)≤ρ, the applicability of this chapter is extended to a confidence level (1 −ρ), similar in spirit to Corollary 3.1 (see page 63). 114 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
6.2. EXPLORATION-EXPLOITATION APPROACH 6.2.2 Exploitation It has been proven that the prediction error d(k) vanishes when the density of the data set becomes infinite, up to 2¯ e. In practice, adding data points to the model increases computation times, and hence, it may not be the ideal procedure to include every new data point observed. Instead, in this work we propose to implement an exploitation policy, adding only informative data points, that is, those that are not close to data points already seen. We characterize the term close by another threshold of the distance, such that a new data point qis not informative if it belongs to the well-known region, defined as Wt=w: min(kw−wik)≤τt,∀i∈ I ND 1,(6.3) for a threshold 0 ≤τt≤τr. This threshold hyper-parameter has to be appropriately chosen, according to the information added by the inclusion of qin D. The dependence of new data points w.r.t. the data stored in Dwas studied in [72], in the context of Gaussian processes for learning-based control. A procedure to prune uninformative sample points given a Lipschitz constant estimation is given in [34]. An example of the exploration-exploitation algorithm for a two-dimensional input space is shown in Figure 6.1. In this figure, an initial data set is considered with ND= 3. New data points are drown randomly within the safe region Wr, but only added if they do not belong to Wt, with τr= 0.1 and τt= 0.05. After 25 iterations ND= 68. The data set is updated every time step, denoting D(k) the data set at time instant k, and D(0) the initial one. Note that the safe and the wellknown regions are also time-dependent (i.e., Wr(k),Wt(k)). This update policy is such that D(k+ 1) = D(k) if w(k)∈ Wt(k) D(k)∪(y(k+ 1), w(k)) if w(k)/∈ Wt(k).(6.4) The same occurs with the estimation of the Lipschitz constant. The recalculation of Lis done recursively [36]. It can be proven that this estimation tends to the real Lfwhen the data set becomes infinitely dense. Computationally, this recursion is linear w.r.t. the cardinality of the data set, O(ND), in contrast to other existing methods (such as Gaussian processes, which are quadratic, O(N2 D)). 115 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
6.2. EXPLORATION-EXPLOITATION APPROACH 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 (a) Initial data set 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 (b) Safe region 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 (c) New data points 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 (d) Well-known region 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 (e) Updated data set 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 (f) Dafter 25 iterations Figure 6.1: Exploration-exploitation algorithm. 116 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
6.3. ONLINE LEARNING CONTROLLER Notice that from the computation point of view, the exploration-exploita- tion algorithms barely increase the calculation times. It is not necessary to calculate a closed-form of the sets Wr(k) and Wt(k). Instead, it is only necessary to check whether a given query point qbelongs to them. This is carried out evaluating the minimum distance to the data set WD; that is, q∈Wt Wriff min(kq−wik)≤τt τr,∀i∈ I ND 1.(6.5) It is also important to remark that these distances kq−wikare used to make predictions for that regressor, so they were already computed in the prediction step. 6.3 Online learning controller In this section, a predictive controller that makes use of the data-based prediction model and the exploration-exploitation approach of Section 6.2 is presented. In order to be able to follow references (possibly outside of the initial data set), we propose to use a robust MPC for tracking [85] that explicitly takes into account at each time step the current safe region, in order to guarantee a given uncertainty bound on the predictions. MPC for tracking is designed to guarantee stability in presence of sudden reference changes, even if they are not reachable, which may be the case in the exploration scenario. To this end, the MPC optimization problem considers an artificial reference (us, ys) as additional decision variables. The deviation of the system to this reference is penalised along the prediction horizon, by means of a stage cost of the form `(y−ys, u −us). A term VO(ys−yt) is added to the cost function, in order to penalise the deviation of the artificial reference to the true reference, denoted (ut, yt). Therefore, we propose the following controller, capable of exploring unseen regions, by forcing the system to stay in the explored area Wr. Combined with the exploitation algorithm of adding data points only if q /∈ Wt, one can step by step move onto unexplored areas, while maintaining a prediction error bound suitable to ensure robust stability. The resulting optimization 117 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
6.3. ONLINE LEARNING CONTROLLER problem PNis: min u,ys,us VN(x(k),u, us, ys;yt) = N−1 X i=0 `(ˆy(i|k), u(k+i); ys, us) + VO(ys−yt) (6.6a) s.t. ˆx(0|k) = x(k) (6.6b) ˆx(j+ 1|k) = ˆ F(ˆx(j|k), u(k+j)), j ∈ I N−1 0(6.6c) ˆy(j|k) = Mˆx(j|k), j ∈ I N 1(6.6d) u(k+j)∈ U, j ∈ I N−1 0(6.6e) ˆy(j|k)∈ Yj, j ∈ I N 1(6.6f) us∈U(6.6g) ys=ˆ f(xs, us) (6.6h) ys∈ YN(6.6i) ˆy(N|k) = ys(6.6j) ˆw(j|k)∈ Wr(k), j ∈ I N 1(6.6k) ws∈ Wr(k),(6.6l) where xs= (ys, . . . , ys, us, . . . , us), ws= (xs, us), and is a design parameter lower than (and close to) 1. For simplicity we may group the constraints u∈ Uand y∈ Y under the notation (y, u)∈ Z. Note that a terminal equality constraint is included. Besides, a set of tightened constraints Yjis considered, as in Chapter 3. They are defined as Yj=Y B(dj(µ)),(6.7) where djis a function of µthat can be calculated as per equation (3.7) (on page 58). Notice that the updating policy guarantees recursive feasibility of the controller: since the regions Wrare dynamical, they evolve with time, and hence if (y(k), u(k)) ∈ Wr(k) then (y(k+ 1), u(k+ 1)) ∈ Wr(k+ 1). The set of tightened constraints counteracts the effect of the error between the real plant and the data-based model (which is bounded by µ) on the constraints. Hence, this set must not be empty for any j∈ I N 1. This is stated as an assumption in chapters 3, 4 and 5, conditioning the feasibility of the controllers. In the proposed approach, the bound on the prediction error within the safe regions is a design parameter (see Property 6.1). This implies that in 118 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
6.3. ONLINE LEARNING CONTROLLER general, τris chosen so that YNis not empty, which is an important property from the implementation point of view. If the error bound for the whole state space were considered, it could be too large to obtain non-empty tightened constraints. Given the procedure to calculate such sets in Chapter 3, and the definition of µin (6.2), the maximum admissible value of the exploration radius τmax rcan be explicitly obtained via reverse engineering, provided that Assumption 6.1 holds. The closed-loop properties of the proposed controller are presented next. 6.3.1 Stability analysis The ingredients of the optimization problem are required to satisfy the following assumption: Assumption 6.2. 1. The stage cost function `(y, u;ys, us)is a positive definite function, and `(y, u)≤αy(ky−ysk) + αu(ku−usk), for two K-functions αy, αu. 2. The offset cost function VO(ys−yt)is a subdifferentiable convex positive definite function such that the best reachable reference y0 s= arg min ys∈YN VO(ys−yt) is unique; and VO(ys−yt)−VO(y0 s−yt)≥αO(|ys−y0 s|), for a given K∞-function αO. Define Yeas a convex set of reachable equilibrium points, Ye⊆ {y:∃us∈U:ˆ f(xs, us) = ys}.(6.8) Then, the following assumptions are needed to derive the recursive feasibility of the optimization problem. Assumption 6.3. YN−1⊆ Ye. Assumption 6.4. For all xsuch that y∈ YN−1there exists a continuous control law uF=κF(y)∈ U such that ˆ f(x, uF)∈ YN. 119 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
6.3. ONLINE LEARNING CONTROLLER Theorem 6.1. Suppose that Assumptions 6.1-6.4 hold for the optimization problem PN. Let κN(x)be the control law derived from the solution of PN applied using a receding horizon policy. Then, for any x(0) ∈ Z, the system controlled by the control law u(k) = κN(x(k)) is recursively feasible, stable, and the constraints are always satisfied, i.e. u(k)∈ U, y(k)∈ Y,∀k. Proof. This proof is sketched in two parts. First, recursive feasibility is proven. Denote the optimal solution of the optimization problem PNat time step kand state x(k) as u∗(k). Let define the shifted sequence candidate for the optimization problem for x(k+ 1) as ¯ u: ¯ u(k+ 1) = ¯u(j|k+ 1) = u∗(j−1|k), j ∈ I N−2 0 ¯u(N−1|k+ 1) = uF ¯us(k+ 1) = ¯us, (6.9) where ¯usis such that y(N|k+ 1) = ˆ f(x(N|k+ 1),¯us). From (6.6e) it is known that ¯u(j|k+ 1) ∈ U for j∈ I N−2 0; and from (5.5) and Chapter 3 it is derived that ¯y(j|k+ 1) ∈y∗(j−1|k)⊕B(dj−1)∈ Yj. Then, for j=N−1, given Assumption 6.4 it is known that uF∈ U and that ¯y(N−1|k+ 1) ∈y∗(N|k)⊕B(dN)∈ YN⊕B(dN)⊆ YN−1. Moreover, taking Assumption 6.4, and equations (6.6h) and (6.6i) it is derived that for ¯ys=ˆ f(¯x(N−1|k+ 1), uF)∈ YN⊆ Ye there exists a ¯us∈U:ˆ f(¯xs,¯us) = ys. Hence given Assumption 6.3, and equation (6.6j), y(N|k+1) ∈ YN, which proves that the problem is recursively feasible. In second place, stability is proven. Let ¯ uin (6.9) denote the candidate shifted sequence, and define and intermediate sequence (solution of the nominal problem) as ˜ u(k+ 1) = ˜u(j|k+ 1) = u∗(j−1|k)j∈ I N−2 0 ˜u(N−1|k+ 1) = u∗ s(k) ˜us(k+ 1) = u∗ s(k). (6.10) 120 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
6.3. ONLINE LEARNING CONTROLLER Given the definition of VNin (2.9a)1and its continuity ensured by Assumption 6.2, we know that the difference VN(x(k+ 1),¯ u(k+ 1),¯ys(k+ 1),¯us(k+ 1)) −VN(x(k+ 1|k),¯ u(k+ 1),¯ys(k+ 1),¯us(k+ 1)) is upper bounded by some K-function θ1(kx(k+ 1) −x(k+ 1|k)k)≤θ1(µ). Analogously, VN(x(k+ 1|k),¯ u(k+ 1),¯ys(k+ 1),¯us(k+ 1)) −VN(x(k+ 1|k),˜ u(k+ 1),˜ys(k+ 1),˜us(k+ 1)) is equal to "N−2 P j=0 `(y(j|k)−¯ys,¯u(k+j)−¯us)−`(y(j|k)−˜ys,˜u(k+j)−˜us)# +h`(y(N|k)−¯ys, uF−¯us)−`(y(N|k)−˜ys, u∗ s−˜us)i +hVO(¯ys−yt)−VO(˜ys−yt)i.(6.11) Provided that ¯u(k+j) = ˜u(k+j)∀j∈ I N−2 0, and given the definition of the stage cost in Assumption 6.2, then the first bracket is upper bounded by (N−1)σ1(k¯ys−˜ysk) + σ2(k¯us−˜usk), for some K-functions σ1, σ2. Let gs(ys) denote the function that yields equilibrium inputs of the model, given its equilibrium output. Since the model is Lipschitz continuous and Z is compact, then gsis uniformly continuous, and then σ2(k¯us−˜usk) = σ2(kgs(¯ys)−gs(˜ys)k)≤σ3(k¯ys−˜ysk). Recall that for a given K-function σ, σ(k¯ys−˜ysk) = σ(ky(N|k)−y(N|k+ 1)k) =σ(kˆ f(x(N|k),˜us)−ˆ f(x(N−1|k+ 1), uF)k) =σa(µ) + σb(k˜us−uFk), 1The dependence of VNwith ytis omitted, for the sake of conciseness. 121 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
6.3. ONLINE LEARNING CONTROLLER where given the continuity of the model, the derivation into σais explained in Chapter 3, and taking κF(y) = gs(ys) if y∈ YN, σb(kuF−˜usk) = σb(kκF(y(N−1|k+ 1)) −κF(˜ys)k) ≤σc(ky(N−1|k+ 1) −˜ysk)≤σd(µ). The second bracket of (6.11) is upper bounded by σ4(k¯ys−˜ysk) + σ5(kuF−¯usk), for some K-functions σ4, σ5, provided that u∗ s= ˜us. Notice that, for the K-functions σe,f,g,h, σe(kuF−¯usk)≤σf(kuF−˜usk) + σg(k˜us−¯usk)≤σh(µ). The third bracket is also upper bounded by some K-function σ6(k¯ys−˜ysk), provided that Assumption 6.2 holds. Apart from that, the nominal tracking of changing piece-wise constant reference signals [85] upper bounds VN(x(k+ 1|k),¯ u(k+ 1),¯ys(k+ 1),¯us(k+ 1)) −V∗ N(x(k)) by some K-function −α(ky(k)−˜ys(k)k). Hence, VN(x(k+ 1),¯u(k+ 1),¯xs(k+ 1),¯us(k+ 1)) −V∗ N(x(k)) ≤θ(µ)−α(ky(k)−y∗ s(k)k). Thus, by optimality V∗ N(x(k+ 1), u∗(k+ 1),x∗ s(k+ 1), u∗ s(k+ 1)) ≤VN(x(k+ 1),¯u(k+ 1),¯xs(k+ 1),¯us(k+ 1)), which yields a bound on the increment of V∗ N. From this theorem, the following corollary can also be proven. Corollary 6.1 (Convergence).In case that the prediction error of the KI model µ(k)tends to 0 , the system converges to the best reachable reference y0 s. 122 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
6.4. CASE STUDY Figure 6.2: Set of equilibrium points of the CSTR, depending on the parameters; and set of constraints considered in each chapter. Remark 6.1. Note that the online algorithm presented here decreases the prediction error while operating the system. However, even in the absence of the exploitation policy (i.e. τt= 0) and infinitely dense data sets, the maximum prediction error µvanishes up to a factor of the noise level [36]. Illustrations of the operation of the proposed controller are presented in the following section. 6.4 Case study The continuously stirred tank reactor (CSTR) presented in Chapter 2 (on page 42) is considered. The constraints in the input are 300 K ≤Tr≤400 K, and in the output 0 ≤CA≤0.88 mol l−1. The different configurations of the CSTR used in this thesis are represented in Figure 6.2. Online learning We first consider a case in which, at the beginning of the simulation, very few data points are known, just 300 corresponding to some equilibrium points of the system, obtained from the experiment shown in Figure 6.3. 123 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
7.1. SUMMARY OF CONTRIBUTIONS To overcome such problems, several methods to enhance the standard kinky inference learning class have been proposed. Specifically, the projected KI presented in Chapter 2 aims to reduce the computational times of the inference, which allowed for the derivation of output-feedback model predictive controllers. Second, the componentwise H¨older predictor presented in Chapter 4 is proven to improve the standard KI method, yielding smaller worst-case prediction errors. Besides, it was proven that this method enlarges the feasibility region of controllers based on it. With respect to the controllers proposed, the objective has been to develop stabilizing and robust designs of MPC without terminal constraint, which is itself an important contribution to the field of predictive control. This was motivated by the avoidance of the calculation of invariant robust sets, which is undoubtedly the most difficult ingredient in the design of MPCs, even more in the case of prediction models based on data. In this context, new controllers have been proposed. First, an output-feedback MPC based on the smoothed projected version of the KI algorithm was derived in Chapter 2. It was proven able to robustly satisfy soft constraints in the outputs and hard constraints in the inputs, while converging to the reference. Next, an MPC able to satisfy hard constraints in the output, which is a much more complex case, was proposed in Chapter 3. Instead of adding a terminal constraint based on invariant sets, the terminal cost was weighted by a factor λ. Based on this approach, another controller with larger domains of attraction was derived using the componentwise H¨older predictor of Chapter 4. Chapter 5 opened the door to the online approach, in which the controllers proposed sought to improve their performance by the inclusion of new measurements that became available during the operation of the plant. An online learning MPC was proposed, based on a double prediction model: a safe model was derived offline with the initial data set, in charge of ensuring robust stability and constraints satisfaction. Then, an updated model was used to improve the predictions, yielding smaller costs in the optimization problem, and hence more efficient performance of the closed-loop. Finally, an online learning MPC for tracking, based on an explorationexploitation criteria, was presented in Chapter 6. The design proposed in this chapter is able to uncouple the uncertainty and the domain of attraction of the controller, by means of combining the exploration approach with the tracking formulation. Hence, the controller was proven able to maintain robust stability of the closed-loop even in the case of low-dense data sets. 130 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
7.2. FUTURE WORK Besides, the exploitation criteria prevented the model to become computationally intractable. These facts allowed the system to be steered to the best reachable reference, even if it was in unexplored areas of the workspace, while robustly satisfying the constraints. 7.2 Future work The previous paragraphs have briefly summarized the contributions made to the field of data and learning-based predictive control, which has gained increasing attention during the years of elaboration of this thesis. As a consequence of the research presented, further works will certainly arise motivated by those already mentioned. Among other lines, we will be interested in the validation of the proposed methodologies in experimental plants, which will include the design of the experiments and processing of the data sets, as well as the case-dependant design of the models and of the control law, in order to illustrate the proven behaviour of the proposed strategies. In addition, note that several methods were proposed to enhance the predictions. Combining these methods is a short or medium term future objective, since it may result in a broader and more general inference algorithm suitable for many machine learning applications, as well as to model predictive control. Hence, novel and more efficient implementations of the KI class will be proposed in the future. In line with the previous point, the parallelization of the kinky inference algorithm is a project that we have already begun [119], and that will undoubtedly allow for efficient and fast implementations of the prediction models. This will consequently increase the list of case studies where the methodologies presented in this thesis may find application, like embedded and real-time systems, which could be combined with the controllers proposed. On other matters, the stochastic study of the controllers derived must be addressed. For example, the scenario approach, order statistics or probabilistic optimization [4] may open useful frameworks for the data-based approaches considered in this thesis. The combination of other inference methods with the ones proposed in these chapters may lead to enhancing the overall performance, so deriving the stability analysis for models based on, e.g., Gaussian processes [94] will 131 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
7.2. FUTURE WORK be an interesting research topic. In addition, the application of these methods to changing operation conditions or to systems varying in time (either their structure, parameters, etc.) can be studied, using the online configurability of the techniques proposed in this thesis. To this end, dynamical data sets, able to forget past measurements in detriment to new ones can be considered. Finally, another current line of research [105] is the extension of the proposed control laws to take into account economic criteria in the optimization step, so further investigation on how the learning-based approach could be adapted to such controllers is a promising field. 132 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
Appendix A Working with data In general, when identifying a model based on input-output data, these data must be informative and rich enough to correctly capture the dynamics of the system within the considered workspace. To this end, this appendix presents some aspects that must be taken into account when working with data-based algorithms. A.1 Input design The quality of the available data is a key aspect in system identification, even more in the nonlinear case (as this thesis considers), since the nonlinear behaviour forbids the model to extrapolate well enough when estimating over unseen areas of the workspace of the plant. The input-output data must capture the dynamics and modes of the system correctly. To this aim, input design is a crucial step when obtaining valid data sets for input-output identification. Consequently, this thesis in general followed the procedures explained by ˚ Astr¨om & H¨agglund [12], Billings & Fadzil [23] and Rivera & Jun [134] for the design of input signals. Therefore, to obtain data trajectories in order to learn the case studies considered in this thesis, we avoided the use of grids in the inputs and/or in the initial states. Although grids in the considered spaces may result a quick and easily implementable approach in simulated case studies, they are unlikely to represent applicable methods in real plants, and they would probably fail in capturing most of dynamics in the nonlinear system. Instead, there are two other types of input trajectories that are considered: 133 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
A.1. INPUT DESIGN Pseudorandom binary sequences (PRBS) PRBS are widely used in the identification and telecommunications fields. The signals have a square shape, alternating between two values at apparently random times. It is said to be pseudorandom since it is generated via deterministic computer algorithms. Following this spirit, this thesis uses pseudorandom sequences (not binary), such that both amplitude and length of the steps of the designed input signal are (pseudo-)random. Chirp signals Chirp signals are sinusoidal sequences whose frequency f varies along time. Apart from identification, they are commonly used in sonar and radar applications. This work only focuses in chirp signals whose frequency increase linearly with time, such that u(t) = sin 2πff−f0 Tt+ f0t,(A.1) where f0≤ffare the initial and final frequencies, respectively, and Tis the length of the signal. The initial frequency can be taken approximately equal to the crossover frequency of the system, obtained closing the loop with a relay at a given set-point. The final frequency can be taken two or three orders of magnitude higher than the initial frequency. An example of both signals is shown in Figure A.1, as well as their application to the SISO tank presented in Section 1.3.2. The PRS’s amplitude is chosen from a random uniform distribution from 0 to 1, and the same for the length of each step, between 0.2 and 30 s. The initial frequency of the chirp signal is calculated with a relay test over (ys, us) = (0.25 m,0.5 m3s−1), yielding f0= 4.3 mHz. The final frequency is set to ff= 0.43 Hz. Both signals last T= 100 s. Note that, in general, pseudorandom sequences (Fig. A.1c) generate data points which tend to the set of equilibrium points of the system (cf. Figure 1.2b), while the chirp signals are more likely to generate rounded trajectories in the input-output space (see Figure A.1d). 134 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
A.2. DIMENSIONALITY REDUCTION 0 20 40 60 80 0 0.5 0 20 40 60 80 0 0.5 1 (a) PRS 0 20 40 60 80 100 0 0.5 0 20 40 60 80 100 0 0.5 1 (b) Chirp 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 (c) Data set generated with PRS 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 (d) Data set generated with chirps Figure A.1: PRS and chirp signals applied to a water tank, and the data set that they generate. A.2 Dimensionality reduction High dimensional input spaces may seriously harm the applicability of databased algorithms when learning the ground truth function. Dimensionality reduction is a sub-field of machine learning that aims to reduce the number of independent variables in consideration, that is, to decrease the number of input dimensions nw, while worsening the performance at least as possible. In this thesis, we consider a method that applies a linear transformation to the input space, named principal component analysis; but there exists many others, like discriminant analysis [14] or non-negative matrix factorization [75]. 135 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
A.2. DIMENSIONALITY REDUCTION A.2.1 Principal component analysis Principal component analysis (PCA) was invented by K. Pearson in 1901 [123]. In linear algebra it is called single value decomposition [45]. The principal component of an n-dimensional data set is the (n−1)-dimensional line that minimizes the average squared distance from the points to the line (similar in spirit to the least squares). Then, the second principal component is the line, orthogonal to the first, that minimizes the same value; and continuing until the n-th line. The resulting set conforms an orthogonal basis called the principal components [69]. Given a data set WD={wi∈Rnw, i ∈ I ND 1} ∈ W, the coefficients that transform Winto the space of principal components of WDare gathered in the matrix T. The transformed space is denoted V ⊆ Rnw. Analogously, the transformed input data set is denoted VD∈Rnw×ND. They are obtained applying the linear transformation V(D)=TW(D)={v=Tw, w ∈ W(D)}.(A.2) Example A.1. Figure A.2 shows an input data set WDof ND= 30000 data points with nw= 4, from the continuously stirred tank reactor (CSTR) plant described in Chapter 2, obtained applying a sequence of chirp signals. It is a SISO system with na= 2 and nb= 0, so w(k)=(y(k), y(k−1), y(k−2), u(k)). Both inputs and outputs are scaled to range from 0to 1. Figure A.3 shows the transformed input data set VD. Notice that, so far, the principal component analysis does not reduce the dimensionality of the input space. It only applies a linear transformation. In order to decrease dimensionality, an appropriate number nh≤nwof components (or dimensions) is selected. To this end, the standard deviation of VD is obtained as Si(VD) = v u u t1 ND−1 ND X j=1 |vi,j −mi|2,(A.3) where Siis the standard deviation of the i-th component of VD,vi,j its i, j-th element and miits mean, for i= 1, . . . , nw: mi(VD) = 1 ND ND X j=1 vi,j.(A.4) 136 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
A.2. DIMENSIONALITY REDUCTION Figure A.2: Input data set WDobtained applying chirp signals to the CSTR. Note the high correlation among y(k), y(k−1) and y(k−2). Figure A.3: Transformed data set VDapplying PCA to the data set represented in Figure A.2. Note the smaller ranges in vias iincreases. Then, the nhcomponents whose standard deviation is lower than certain threshold are discarded, yielding a lower-dimensional space. The threshold may be selected given by some percentage σof the maximum standard deviation (the one of the first component of VD), that is, nh= max{i:Si(VD)≥σS1(VD)},(A.5) for some 0 ≤σ≤1 and for all i∈ I nw 1. Table A.1 represents the standard deviations of the example shown in Figure A.3. The ranges in VDare represented in Figure A.4. In this case, it is sensible to choose σ= 10%, and thus, a lower dimensional input data set is obtained with the first two components of the transformed variables. 137 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
A.3. DATA PROCESSING i1 2 3 4 Si0.3966 0.2829 0.0219 0.0012 100Si S1100 71.339 5.5298 0.3132 Table A.1: Standard deviation of the case example shown in Figure A.3. (a) WD(b) VD Figure A.4: Data set ranges by input dimension, before and after the PCA transformation. Note that this process of dimensionality reduction is fully accessible and reversible by the algorithm, provided that the transformation matrix Tis invertible. Hence, the reduction of input dimensions may be applied only in some parts of the control process. For example, in the learning method proposed in Chapter 2, queries are labelled as close to the data set according to the lower-dimensional space, but then the prediction is computed using the whole original input space. A.3 Data processing In general, the process of obtaining data and learning with it involves a large amount of suboperations, constituting a whole field of machine learning, from sampling [46] to data mining [57], which is sometimes referred to as the knowledge discovery in databases [49]. This thesis is not focused in this field, and only two further operations are performed with the raw observed data set, apart from constructing the regressor presented in (1.11), given past input-output trajectories and applying PCA. These operations are scaling and data pruning. 138 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6
A.3. DATA PROCESSING When working with databases, it is desirable (if not compulsory) to work with scaled data, such that every value ranges along the same order of magnitude, typically between 0 and 1. This is how the algorithms of this thesis are going to work. Recall that the regressor (1.11) groups output and input terms. The ranges and units of the inputs usually differ from the outputs, or even among them, so the process of computing metrics among data points does not make sense without a scaling. This is done such that for a set of N measurements of the vector v∈Rn, vs,i =vi−minivi maxivi−minivi ,∀i∈ I N 1.(A.6) Note also that after a principal component analysis, another scaling is required (see Figure A.4). A.3.1 Data pruning The size of the data set of observations, i.e., ND×nw, plays a crucial role in the applicability of any machine learning method, specifically in those considered in this thesis. Section A.2 addressed the problem of dimensionality reduction, giving methods to decrease nw. This section considers reducing ND, that is, the number of points in the data set. To this end, if the data set is too large, it may be decimated, pruning points that are of less interest to the algorithm. A trade-off between prediction accuracy and computational complexity comes up, and it is the engineer’s duty to set the threshold on how large the cardinality of the data set is allowed to be. There exists many techniques for data pruning. T. Alpcan [8] and H. Kingravi [72] studied the dependence of a data point w.r.t. the data stored in D, in the context of Gaussian processes for learning-based control. A procedure to prune uninformative sample points given a Lipschitz constant estimation is given in [34, Appendix C3]. A simple procedure that is applied in this thesis consists in decimating points that are close to each other. The term close is defined with a threshold τp>0 of the distance between data points, such that the j-th sample is considered close to the data set (and therefore uninformative, so it can be removed) if min(kwj−wik)≤τp,∀wi6=wj∈ WD.(A.7) 139 Código seguro de Verificación : GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 | Puede verificar la integridad de este documento en la siguiente dirección : https://sede.administracionespublicas.gob.es/valida ÁMBITOPREFIJO CSV FECHA Y HORA DEL DOCUMENTO GEISER GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6 22/06/2020 12:57:34 Horario peninsular Nº registro DIRECCIÓN DE VALIDACIÓN O00008745e2000025200 https://sede.administracionespublicas.gob.es/valida GEISER-eb19-3fbe-5649-472c-860f-92a2-5b60-e3b6