Centralized Model Predictive Control DMHE Distributed Moving Horizon Estimation DMPC Distributed Model Predictive Control Energy IN TIME (EiT) Simulation Based Control for Energy Efficiency Operations and maintenance EU European Union FAC Fault Adaptive Control FCU Fan Coil Unit FDD Fault Detection and Diagnosis FP7

Acknowledgment

I would like to thank my supervisors Frederic Hamelin, CRAN-CNRS Nancy France and Joseph Yamé, CRAN-CNRS Nancy France for giving me the opportunity and their support in overcoming the obstacles. I am also thankful to Didier Theilliol, CRAN-CNRS Nancy France who guided me since the start of my research Journey by redirecting my application to my supervisors till the finalization step by his valuable advice ever so often. I take this opportunity to thank Joseba Quevedo, UPC Barcelona, and Didier Dumur, CentralSupelec Paris for their kind patience and time to rigorously evaluate the research work. I deeply appreciate the effort by Estelle Courtial, Orlean France for taking interest in my research work and the kind consideration to examine my research. My special thanks to Christoph Simon, CRAN-CNRS Nancy France for his always cheerful conversations and I am fortunate to have a mentor like him. Million thanks to Sabine Huraux for making all the administrative processes so smooth and her kindness that adds joys in the lab environment.

I especially thank Prof. Vicenç Puig UPC Barcelona, without his moral as well as technical support, I would not have imagined the day of completion of my thesis. I will not hesitate to define his working methodology as an illustration of the ideal and inspiring supervision. His absolute availability, mutual respect, and positivity are the worth lessons I could learn from him with his excellent technical guidance. With all the gratitude, I dedicate this research work to him as a token of acknowledgment to his kindness and exceptional guidance. I am grateful to Stefano Reverso, Senior Scientist UTRC, Ireland for countless discussions and his valuable advice on numerous hurdles without which this journey would have been onerous. I am overwhelmed by realizing his confidence on my capabilities and expertise, that he encouraged me to perceive the opportunity in UTRC and indeed due to his hand, I could secure the position. I certainly owe him tremendous gratitude for his unconditional trust in me.

After a long day at work, it is always heartwarming to share meals and laughs with loved ones. I consider the loved ones are essential fuels to keep the show on. Meera Magnani Ottowa Canda, Savitha Vijay Nancy France, Saikat Muzumdar Vancouver Canada, Hector Sanchez Bonne Germany and Joseph Panach and Helem Corrales, UPC Barcelona are the ones who were always there to make sure I stay sane in all the process. Receiving calls from my (outside country) family Joelle Mallet and Nathalie Mauguin Arles France and hearing their concern to find out if I am doing well, always gave me strength and positive energy. My fellow lab mates and friends, Krishnan Shrinivasrangan and Marc-label Bisch made the working experience joyful and on whom I rely for a hand for any issue at a given time.

I am grateful to my family for freeing me from the conservative social boundaries by letting me stir my ship to embark my journey and encouraging me to think on my feet and make decisions. I realize the friction they would have faced back home for setting their bird free. I consider my family as certainly a good grace.

Finally, I will never miss mentioning Louis Thabault Nancy France for his infinite patience, especially for last few months. I thank him a million for his unconditional care and love.

Motivation

The energy crisis of the early seventies has been certainly one of the strong impetus in the changes of the building industry towards more energy-efficient buildings without sacrificing comfort. Furthermore, in the last several years the awareness of global warming has been also an incentive to radically transform the world economy so that growth may be still possible while halving the emission of greenhouse gas by 2050, including in energy which accounts for about 70 % of emissions. As energy requirement and fuel consumption of heating, ventilation and air-conditioning (HVAC) systems have a direct impact on the cost of operating a building as well as an impact on the environment, building energy performance has become an important issue in many countries, particularly in Europe where the building sector accounts for approximately 40% utilization of overall primary sources (more details on www.ec.europa.eu/energy).

Taking initiative towards energy efficiency for the HVAC building systems, multiple projects are instigated under the FP7 framework by the European Union. Energy in TIME is one of these projects focused on reducing the energy consumption in the operational stage of existing non-residential buildings, resulting in energy savings of up to 20%. The aim of the project is to develop Smart Energy Simulation Based Control methods that are energy efficient without compromising the thermal comfort of the occupants in the buildings. This project brings together a total of 13 partners from 8 different European countries. The project is well organized and divided into nine workpackages. Every workpackage consists of several modules assigned to the group of partners based on their expertise in the respective fields. These workpackages are organized to incorporate different frameworks as Energy Equipment, Control Systems (algorithm and equipment), ICT (communications, sensor and data repositories) and OMM (Operations, Management and Maintenance). The designed control tool is validated on the existing non-residential buildings in the different Europe locations with different climates. These demonstration sites consist of four buildings with different topologies including an airport, offices and test labs, a commercial and office building, and a hotel. Finally, the building energy management systems are controlled automatically and remotely for the given demonstration sites. This serves as proof of concept of the Energy IN TIME solution.

Furthermore, University of Lorraine (UL) mainly deals with the OMM framework involving the tasks like Fault Diagnosis and Adaptive Control applied to HVAC building systems. This thesis work is the part of Energy IN TIME project fulfilling the objectives assigned to UL. The principle objectives delivered by the UL are as: i) Control and Monitoring for the economic building operability within the user defined performance requirements, ii) Maintenance by enabling the early detection of equipment malfunctions and defective system behavior followed by the appropriate corrective action to continue the normal building operability, iii) Scalability of the proposed solution for different classes of HVAC systems This work is motivated to attain these objectives by proposing innovative solutions. To address the energy efficient building control, we propose an economic model predictive control formulation for the cost-effective operation without compromising the thermal comfort of the occupants inside the buildings. We analyze the common fault events occurring in these demo sites and suggest a fault detection, isolation and diagnosis method. Furthermore, we present a method to adapt the faulty condition to maintain the normal building operation while ensuring the energy efficiency. Finally, the concern about the scalability of the solution to the large-scale buildings is answered by proposing the new scheme of distributed control. Lastly, every stage of the proposed solutions are validated on the simulation platforms and the simulation results are provided.

Energy IN TIME Project

The research work carried out in this thesis is supported by the Energy IN TIME (EiT) project funded by the European Union (EU) within the 7th Framework Program FP7-NMP, Sub-program EeB.NMP.2013-4: Integrated control systems and methodologies to monitor and improve building energy performance.

The aim of the project is to develop a Smart Energy Simulation Based Control method which will reduce the energy consumption in the operational stage of existing non-residential buildings, resulting in energy savings of up to 20%. The Table 1 The detailed description regarding the project objectives, workplan and the motivation for the research work is presented in the next sections.

Project Objectives

Buildings Operational stage represents 80% of building's life-cycle cost of which 50% is consequence of the energy use. Up to 90% of the buildings' life cycle carbon emissions occur during their operational phase, mainly as consequence of the HVAC, lighting and appliances' energy use. Therefore, energy and cost saving strategies addressing this building operation phase will have a major impact in the building life cycle cost.

Energy IN TIME project goes beyond existing building control techniques, developing an integrated control & operation approach, that will combine state of the art modeling techniques with the development of an innovative simulation-based control technique with the overarching objective of automating the generation of optimal operational plans tailored to the actual building and users requirements. This approach will allow reducing system inefficiencies and contributing to improve building energy efficiency and comfort.

The target for Energy IN TIME solution will be existing non-residential buildings, which present the building typologies that guarantees higher impact and room for improvement due to the variety and quantity of facilities and equipment covered and the operational management model used in them. A control tool will be implemented in the building energy management systems to be automatically and remotely operated. The methodology for the enhancing solution implementation will be defined for existing buildings and for its implementation in new buildings since its initial commissioning.

Work Packages and Planning

The Energy IN TIME consortium is made up of 13 partners from 8 different EU countries listed below: The project is well planned organizing all the objectives into the nine work packages assigned for the different groups of the partners based on their respective expertise.

Work Packages and Their Description

The work packages and their summary is as follows:

1. WP1-Requirements and System Architecture :

In this WP, target end-users and buildings requirements and needs will be analyzed and identified. The system architecture and Energy IN TIME specifications will also be defined. A data acquisition methodology and the guidelines for the composition of a communication platform will be created.

WP2-Simulation Reference Model:

This Work Package is focused on the development of advanced simulation tools for the analysis and characterization of buildings: identifying their systems, use, occupancy as a prior phase in their management and optimized control. A forecast tool to determine user behaviour evolution and occupancy schedules will also be designed that will be used to design the system models.

WP3-Whole building Intelligent Control System

This Work Package aims to develop an intelligent control system for whole-building energy optimization algorithms that have the capability to adjust set-points for energy equipment to be adapted to occupancy and weather loads through their predictions, minimizing energy consumption while maintaining indoor comfort conditions. Furthermore, the control system will take into account the variations in building operating conditions to provide corrective actions to counteract the effect of component and system-level fault.

WP4-Diagnosis and Continuous Commissioning

This Work Package will develop procedures and algorithms for the detection and diagnosis of faults at system and equipment level. A methodology for Continuous Commissioning will be defined. The procedures for the establishment of a predictive maintenance system will be developed.

WP5-Energy Decision Support Tool

The Energy IN TIME project will be a powerful tool for energy decision making because of the large amount of data gathered from different buildings, systems, devices, facilities, etc. All this information will be integrated in an information repository and a decision support tool designed through analysis and data exploitation.

WP6-System Integration and Pilot-Scale Validation

The objective of this Work Package will be the integration of the developments obtained from the previous WPs to generate the Energy IN TIME solution, including software tools and equipment. A system validation at laboratory scale will be developed in this WP.

WP7-Demonstration

The objective of this WP is to implement the designed and developed systems in various different types of buildings located in regions with different climate conditions. A period of two years has been estimated for solution implementation and monitoring data processing and performance evaluation of the systems.

WP8-Dissemination & Exploitation

The aim of this WP is to ensure that the project meets its goal and overall objectives of exploitation and dissemination activities. Dissemination activities include media, conferences and publications. An Exploitation Plan will be developed to allow the sustainability of the project beyond its life. The focus for the dissemination of the results of the project will be European countries.

WP9-Coordination

The aim of this WP is to ensure an effective coordination allowing the project to meet its goals and overall objectives

The overall Energy IN TIME architecture comprises the various modules where every module is led by the group of different partners. Figure 1 shows the schematic modules distribution that has 15 modules divided into different frameworks. The affinity among the modules and its relationships are established in detail for the smooth execution of the project.

Demonstration Sites

Energy IN TIME will be validated at existing non-residential buildings in different European locations with different climates. The demonstration sites consist of four buildings with different typologies and building uses, including an airport, offices and test labs, a commercial and office building, and a hotel. The variety of scenarios will serve as a "proof of concept" of the Energy IN TIME solution, in terms of validity and adaptability to different typologies, climate conditions and users behavior. The stated demonstration sites are as follows:

1. Airport -FARO (Portugal) 2. Offices and Test Labs -Bucharest (Romania) 3. Commercial and Office Building -Helsinki (Finland) 4. Hotel -Levi-Lapland (Finland)

Description of demonstration sites

We discuss each demonstration site giving the keynotes about the building type, its characteristics and their HVAC configurations. Please not that for better understanding, we classify the HVAC system into two units as primary and secondary. Primary unit contains the unit operations like boilers, heat exchangers, heat pumps, cooling towers etc. Secondary units are mainly directly related to thermal comfort of occupants and they consist of equipments like Variable Air Volume (VAV) boxes, Fan Coil Units (FCU), Constant Air Volume (CAV) units etc. It is important to note that we focus on the energy efficient operations of secondary units only and maintaining the occupants thermal comfort.

Airport, FARO Portugal

Faro Airport is located 7km away from Faro, the capital city of the Algarve, Portugal. Though Faro accommodates around 30 buildings, the airport terminal is the main building of Faro Airport, where both passengers departures and arrivals are processed. Only this building area is accessible under the project purview.

(a) Building Description

The airport terminal has five different floors assigned to different purposes such as check-in and boarding gets, warehouses, electrical distribution rooms, food drinks supply, emergency operations and thermal plants etc.. Following are the peripheral characteristics of the airport:

-Ground build area: 41.000m 2 ; Built 1989 (2001 last refurbishment) -Open spaces with big flows of people at certain times of the day (e.g. flight arrivals or departures). -The main energy source is electrical.

-High level operational plans must be implemented, aimed at maintaining pre-defined environmental conditions (in terms of comfort, etc.). There are two thermal plants responsible for climate control of the open spaces and offices. These thermal plants includes ice bank, boilers, heat pumps like primary HVAC units. These unit operations are precisely controlled in terms of energy efficiency and control performance. Hence, we focus on the building level secondary HVAC units. On specific applications such as server rooms, and to make sure there is redundancy of climate control. Faro Airport has independent systems like split and close control units. Since 1989, Building Management System (BMS) software operates its centralized HVAC system.

Offices and Test Labs, Bucharest, Romania

The ICPE office Building in Romania accepted the Energy IN TIME challenge and therefore devotes all the available efforts for the implementation of this project and for achieving a useful tool for the end-user.

(a) Building Description The entire ICPE area contains four buildings with four floors each. Although, for the project assignments, only two buildings are made available as test bench. This available building part comprise mainly offices with the characteristics are as follows:

-Total built area: 17.384m 2 Built: 1982 -Closed and distributed spaces with constant flows of people and scheduled occupancy. -Strong presence of solar energy (thermal and photovoltaic). District heating system is the base system for heating, covering -Fixed-schedule operational plans for indoor conditions.

-PLC for Monitoring and Control. This ICPE building has a heating substation which provided hot/cold water supply to the building-wide secondary HVAC units. This heating substation consists of heat exchangers and distribution network. The source hot water supply to this substation is through the district heating system. Then, the temperature of this water is manipulated according to the requirements and circulated to the building secondary units. FCUs a form the secondary unit circuit to maintain the thermal comfort inside the buildings.

Commercial building, Helsinki, Finland

Sanomatalo is the headquarters building for the largest newspaper group in Finland.

The building is located in Kluuvi, the commercial center of Helsinki. Thermal comfort maybe not be trivial due to the extreme climate and the high interactions due to the glass walls. The building is connected to the energy district heating network of Helsinki.

In the substation there are three heat exchangers. Rooms are heated mainly by VAV and CAV units.

Hotel, Levi-Lapland, Finland

The building type is of a hotel with the capacity of 388 beds accessed by around 1000 users. The building is connected to the energy district heating network of Adven Oy.

There is possibility to control the temperature in each room ±2 • C with the temperature regulator through the FCU units.

UL contribution

As mentioned before, the research work for this thesis is aligned with the UL contribution to the project. UL mainly deals with the module 8, 10 and 11 titles as Fault adaptive control, Fault Detection & Diagnostics and Predictive Maintenance respectively. The main objectives under this modules are listed as:

-Maintain building operability within the user-specific performance requirements which includes the thermal comfort of occupants under economic building operation.

-Enable efficient detection, localization and diagnostics of faults in the operation of Building system -Reconfigurable control layer to adapt the control system parameters and objective despite of the presence of faults or performance deviation within its specified energy and comfort performance requirements -Scalability of proposed solution for different building environments

The main focus by UL is therefore to develop methods to address the above objectives. These proposed approaches are then demonstrated on the simulator platform based on the architectures of the HVAC systems provided by the demonstration sites in the project.

Thesis Relevance

This research work aims at meeting the project requirements through the following subobjectives:

-Receding Horizon Control to ensure the economic building operation -Fault detection and diagnosis methodology for the common fault of damper stuck in the VAV boxes -Fault adaptive control based on the receding horizon control for faulty event of damper stuck in the VAV boxes -Distributed model predictive control methods to address the economic building operation of large-scale buildings

To illustrate the developed novel solutions, a brief discussion on the simulation results is provided. Please note that the benchmark building simulator developed in this work is motivated by the common HVAC configurations provided on the demonstration sites. These developed simulators have following properties as:

-The HVAC configurations are based on the Variable Air Volume (VAV) and the Fan Coil Unit (FCU) systems.

-The real weather data for winter for Nancy, France is used to validate the simulation results.

-The occupancy schedule is considered according to the type of building e.g. for the office type building, we consider that the occupants are present from 8.00AM to 6.00PM with a lunch break of 2 hours starting from 12.00PM.

We have developed the mathematical models for these benchmark buildings that are extensively used in validating all the novel approaches presented in this work.

Chapter 1

Introduction

Literature review

In this section, we provide a brief overview of the important literature related to energy optimization for HVAC building systems. We streamline the discussion as per the main areas dealt in this thesis work as follows:

Economic MPC for HVAC

Even though the HVAC building systems are complex in nature, there are numerous efforts to address issues related to energy optimization under the thermal comfort constraints of the occupants. The available literature can be roughly classified depending on the three elements as HVAC configurations, mathematical model and focused control attributes as shown in the Figure 1.1. The HVAC units based on convection principle as VAVs, CAVs, FCUs are most favored due to their simplicity in the installation, control, and energy efficient operations. [START_REF] Kreider | Heating and cooling of Buildings Design for Effiiciency[END_REF] promotes the energy-efficient design of these HVAC equipments. Nevertheless, our focus is concerned with the control aspects of the HVAC systems, especially the advancements in the use of MPC. This is due to the inherent advantages of MPC to consider the disturbances, the system constraints and its commercial availability for implementation [START_REF] Afram | Theory and applications of HVAC control systems e A review of model predictive control (MPC)[END_REF]. Needless to say that the mathematical model and the definition of the cost function are essential requisites for the MPC. There are several approaches found in the literature that investigate the ways of mathematical modeling of different HVAC units [START_REF] Afram | Review of modeling methods for hvac systems[END_REF] e.g. based on RC networks [START_REF] Deng | Building thermal model reduction via aggregation of states[END_REF], [START_REF] Goyal | A method for model-reduction of non-linear thermal dynamics of multi-zone buildings[END_REF] or the based on the thermodynamic behavior of a multizone building [START_REF] Bourhan Tashtousha | Dynamic model of an hvac system for control analysis[END_REF].

Work like [START_REF] Ma | Application of economic mpc to the energy and demand minimization of a commercial building[END_REF] develops an autoregressive exogenous models and interesting formulation of the cost function is proposed to represent the total daily electricity expense, which is a combination of energy and demand costs. Furthermore, some researchers link the mathematical modeling problem with the controller design objectives as in [START_REF] Freire | Predictive controllers for thermal comfort optimization and energy savings[END_REF]. This shows a standard MPC with a quadratic cost function that minimizes the energy and the results are shown on the building with conventional heat radiators. Focusing on the formulation of economic MPC aspects, various ways are discovered persistently by the researchers. These techniques are approximately organized according to the mathematical formulation of cost functions, the evaluation of disturbances and To summarize the attempts in the present literature, the energy optimization and the maintenance of the thermal comfort can be handled on the hierarchical level or in a single control layer with multiobjective MPC formulation. But unfortunately, the issue of maintenance is rarely addressed in the context of economic control. The possibility of accounting the maintenance-aware control alongside the energy efficient operation certainly possesses the strong potential for the better economic operation of HVAC systems. This implies that the economic bundling operation can be further improved if maintenance aspects regarding the HVAC units are also included in the control objectives. This encourages the research to emphasize on the MPC formulation allowing to gather various objectives under the single roof of MPC strategy such as minimizing energy consumption with an additional maintenance awareness and controlling the thermal environment.

Fault diagnosis and fault adaptive strategies for HVAC systems

Inadequate control performances and malfunction of the sensors or actuators may result in a poor quality of thermal comfort inside the building. Sometimes the effects can be severe as e.g excessive energy consumption or the complete shutdown of the HVAC units etc... if the faults are not diagnosed. In spite of the great developments in the field of fault diagnosis in the past few decades, there is still lack of research in its application on the It is straightforward that to take some corrective actions after the fault occurrence, some necessary information from fault diagnosis module should be provided to the controller. There are not sufficient advancements in the literature on fault diagnosis and fault adaptive techniques applied to HVAC building systems. However, we acknowledge the efforts addressing the fault adaptive techniques for VAV type HVAC system in [BLS + 15]. In this article, the data based graphical models are presented to detect the potential faults and fault adaptive control is realized by update the model predictive control algorithm constraints.

Although, in the present literature, system level faults are paid more attention despite the equipment level faults may result in serious consequences. Equipment level faults need to be tackled at the local level instead of the hierarchical layer structure. In this work, we consider a most common of damper stuck fault in the VAV type HVAC systems to provide the fault diagnosis and fault adaptive strategy for the same.

MPC for large-scale HVAC systems

All the solutions discussed above assume a centralized control architecture. Their application in the case of large-scale systems may face various challenges due to the various factors like the size of the system, the internal couplings, the communication network difficulties, etc. Generally, this is addressed by means of the decomposition of largescale systems into the subsystems and the coordination between controllers applied to each subsystem (see Figure 1.2). Several architectures and decomposition methods are proposed in various articles in the literature to address the above issues. Siljack [START_REF] Siljak | Decentralized Control of Complex Systems[END_REF] has summarized various decomposition methods. He has introduced the overlapping and non-overlapping decomposition methods based on the sharing of the variables between the subsystems. Some widely used methods of partitioning large-scale systems into the subsystems are based on bipartite graph theory, decomposition and the inclusion principle [Sca09a] [Lun92]. Other interesting approaches such as relative gain array matrix [START_REF] Haggblom | Partial relative gain: a new tool for control structure selection[END_REF] and Grammians [START_REF] Kariwala | Block relative gain: properties and pairing rules[END_REF] use input-output mapping. Sometimes, due to the structural properties of large-scale systems, it is straightforward to derive subsystem models by a system identification procedure, instead of partitioning the centralized large-scale model.

The coordination between the controllers of subsystems is closely related to the degree of interaction between the subsystems [START_REF] Venkat | Distributed Model Predictive Control: Theory and Applications[END_REF]. If the interaction between subsystems is negligible, each subsystem is controlled independently without any coordination between subsystem controllers. This is named as decentralized model predictive control [START_REF] Lunze | Feedback Control of Large Scale Systems[END_REF]. If the interaction between the subsystems is strong, then the coordination between the subsystem controllers improves overall performance of the system. These coordination strategies may differ depending on the several ways of information exchange between controllers mainly classified as noncooperative and cooperative type [START_REF] Rawlings | Model Predictive Control: Theory and Design[END_REF]. In noncooperative architecture, the subsystem controller optimizes locally the MPC problem using information of other subsystems and it reaches to Nash equilibrium. On the other hand, in the cooperative architecture, the subsystem controller optimizes the global objective and it achieves a Pareto optimal solution [START_REF] Venkat | Distributed Model Predictive Control: Theory and Applications[END_REF]. The coordination can be in a hierarchical architecture [START_REF] Scattolini | Architectures for distributed and hierarchical model predictive control a review[END_REF], where master level optimization problem comprises the shared variables or constraints. The optimal solutions are sent to lower level i.e. the subsystem level controllers as coordination parameters. In some articles, DMPC problem is viewed as the partitioning of the Centralized Model Predictive (CMPC) Problem. This is motivated by some decomposition methods of large-scale convex optimization problems e.g. Dantzig Wolfe decomposition, Benders decomposition [START_REF] Morosan | Distributed model predictive control based on benders decomposition applied to multisource multizone building temperature regulation[END_REF], primal and dual decomposition techniques [START_REF] José | Distributed Model Predictive Con-trol Made Easy[END_REF] [START_REF] Pflaum | Comparison of a primal and a dual decomposition for distributed mpc in smart districts[END_REF].

When dealing with the DMPC problem, the decomposition of large-scale systems into subsystems and the coordination between the subsystem controllers are addressed independently. In this work, we propose a novel approach of addressing the system decomposition and controller coordination issues in two distinct stages. In first stage of this proposed method, the optimality conditions of a large-scale optimization problem are formulated and then decomposed to obtain the subsystems. Further, in second stage, the idea of coordination among the controllers is presented using an optimality condition decomposition approach [START_REF] Antonio | Decomposition Techniques in Mathematical Programming[END_REF]. We implement this distributed control scheme on a given multizone building without compromising the main objectives of energy efficient operations.

Distributed moving horizon estimation

In the final part of this thesis, we concentrate on the estimation techniques in the context of large-scale buildings. This focus on the estimation is motivated due to some obvious advantages like the fault detection and isolation techniques, the possibilities of replacing the measurements with an estimation in case of faulty sensors and the minimization of

Contributions

the number of sensors to save the capital cost of the installation, etc.

The Kalman filter is still viewed as the best available strategy for state space estimation. Despite the current developments in the Kalman filter estimation methodologies, other state estimation techniques are also investigated by the researchers. We can find different type of observer designs [START_REF] Patton | Robust Model-Based Fault Diagnosis for Dynamic Systems[END_REF] and moving horizon estimations (MHE) [START_REF] Bagterp | Moving Horizon Estimation and Control[END_REF]. MHE is becoming popular as it is essentially formulated as an optimization problem that facilitates the inclusion of some physical constraints. This optimization problem over horizon N is solved in a receding horizon manner (similarly to MPC) allowing to estimate the states minimizing the errors introduced by the disturbances and noises. To minimize the uncertainty in the initial states, an extra term is introduced in the objective function. This procedure to solve this problem is repeated at each time instant using a sliding window of N values, hence it is termed as moving horizon estimation. The idea of MHE is extended for nonlinear systems in [START_REF] Victor | A fast moving horizon estimation algorithm based on nonlinear programming sensitivity[END_REF] [START_REF] Alessandri | Advances in moving horizon estimation for nonlinear systems[END_REF]. Also, the extension for MHE considering systems with bounded disturbances is presented in [JRH + 16]. Nevertheless, the size of the optimization problem in MHE for large-scale systems may increase exponentially as the size of the system increases so the application of MHE for large-scale systems is an emerging topic. Different schemes of partitioning centralized MHE problem for large-scale systems are proposed in [START_REF] Farina | Moving horizon partition-based state estimation of large-scale systems[END_REF] while [START_REF] Farina | Distributed moving horizon estimation for sensor network[END_REF] suggests the implementation of MHE for each sensor fulfilling convergence properties. [START_REF] Schneider | Convergence and stability of a constrained partition-based moving horizon estimator[END_REF] proposes a sensitivity based partition technique with a detailed discussion on convergence and stability. However, much work still has to be done in the regard of Distributed MHE (DMHE).

Contributions

This work highlights the important issue of economic building operation while maintaining the healthy environments for occupants inside multizone buildings. As discussed earlier, the main contributions of this work are the responses for the objectives assigned under the Energy IN TIME project by providing novel and efficient solutions. These proposed solutions involve energy efficient building performance on the centralized and distributed architectures along with investigation on the fault diagnosis and fault adaptive control methods. The contributions of this thesis work are ordered in the following steps:

-The demonstration sites provided in the project are of various characteristics such as the type of building, the use of the building, the available HVAC configurations and their functional details, etc. An analysis of demonstration sites shows that the VAV and the FCU type HVAC systems are most favored due to their energy efficient operations and the design flexibility. This further motivates the design of a benchmark building based on these types of HVAC units. We derive the mathematical models for both types of HVAC systems using the thermodynamic relations. These evaluated mathematical models are well-suited for the further validation of proposed solutions.

-Economic building operation of HVAC building systems is the foremost requisite of this work, which is achieved by enforcing the controller structure to minimize the energy consumption and the maintenance cost. This is done by introducing a novel 1.3. Outline approach to MPC formulation which considers: i) the dynamics of most energy consuming equipments in the HVAC system to attain the energy efficient approach, ii) a regularization term to suppress the fluctuations in the setpoints provided to the actuators that lowers the actuators fatigue minimizing the maintenance cost. This is finally validated on the benchmark HVAC building system.

-Sometimes a malfunction of the actuators or the sensors may cause problems on various levels e.g. loosing thermal comfort of the occupants, excessive energy consumption and permanent damage of the HVAC equipment if these faults are not detected early, etc. In this work, we propose an Unknown Input Observer based fault diagnosis approach for some given common faults and that is supported by the simulation results.

-Furthermore, the fault diagnosis information plays a key role in the design of fault adaptive methodologies in order to take the corrective action in the faulty events ensuring the thermal comfort of the occupants. Considering the importance of fault adaptive methods, we propose an MPC-based fault adaptive control strategy, which uses the information of the fault diagnosis module developed in the earlier section.

-Final contribution of this work resides in the innovative control approaches applied to large-scale buildings to overcome the challenges of implementing the centralized control architecture. The proposed methods allow us to decompose large-scale systems into subsystems and to design the distributed control architecture based on the optimality condition decomposition method and sensitivity analysis techniques.

The proposed distributed MPC methods are validated on the benchmark building systems.

-Along the lines of the economic building performance, we address the issue of estimating the system outputs (zone temperatures) in the case of limited measurements. This proposed method can be seen as a novel contribution in the area of moving horizon estimation applied to the large-scale buildings. This method holds a great potential for the fault detection and isolation assignments. We implement the proposed distributed estimation technique on the given benchmark building and discuss its effectiveness with some simulation results.

Outline

The remainder of the thesis is structured as follows:

From the analysis of the demonstration sites provided in the Energy IN TIME project, we describe the most commonly used HVAC configurations and their mathematical models in the Chapter 2. These HVAC configurations are essential and used as test benches to validate the solutions developed in this thesis work.

A reminder of the thesis mainly covers the centralized and distributed control aspects. Chapter 3 present the maintenance-aware economic control and the fault diagnosis and fault adaptive control aspects based on a centralized control system architecture. We 1.3. Outline discuss the detailed formulation of the centralized model predictive control ensuring the economic operation of the building that maintains the thermal comfort inside the zones. To address the fault diagnosis issue, an approach using dedicated unknown input observer bank is proposed and verified through the simulation results. Lastly, we present a fault adaptive control strategy that accommodates the fault and maintains the normal operations of the building HVAC system.

Scalability issue regarding the energy efficient control and the fault detection is addressed for larger systems by extending the proposed centralized control solutions. Chapter 4 explores the techniques of distributed model predictive control for the specified benchmark buildings. Two distinct distributed control methods based on optimality condition decomposition and a sensitivity analysis are discussed in depth. Rigorous simulation results provide the validation of the proposed distributed control approaches. The final part presents a great development for moving horizon estimators. Also, its scope is broadened to large-scale systems. The simulation results prove the applicability of the proposed estimation method in the case where fewer measurements are available in the large-scale building.

Finally, in Chapter 4 we present conclusions and future perspectives.

Chapter 2

Building Frameworks and Mathematical Modeling

The solutions developed under Energy IN TIME project are validated on the existing commercial buildings. These buildings are of various topologies that include a hotel, an airport, and offices with different locations throughout the Europe. These developed solutions are deployed on these demonstration buildings to conduct a proof of concept and in this thesis, we analyze the available demonstration sites and their HVAC configurations. Further, based on these configurations, we design a benchmark building prototype and the simulator for a realistic imitation. This simulator serves as a test bench to investigate the proposed solutions offered in this thesis work to fulfill the project module objectives assigned to UL. The building thermal behavior changes according to the occupancy schedule, the type of HVAC system and the coupling among the zones. From the detailed synthesis of the demonstration sites, we summarize the characterization of the different building topologies with respect to the mentioned attributes. The Table 2.1 shows summary of the demonstration sites.

In the existing HVAC systems, the VAV and FCU type configurations are more popular due to their simplicity and energy efficient performances. Hence we design the benchmark building prototypes based on these HVAC configurations. Further, the proposed solutions under the project e.g. energy efficient control methods, fault diagnosis, and fault adaptive techniques etc. are built to suit these benchmark prototypes. Nevertheless, adequate flexibility in the proposed solutions is preserved so that it will apply to a broad range of building types. Further, we present the VAV and FCU type HVAC configurations and their mathematical models.

Thermal Zone Model

It is essential to understand the thermal behavior of the general air conditioning system that uses the air circulation to maintain the thermal comfort of the occupants. This will allow us to derive the mathematical model for the VAV and FCU type HVAC systems respectively. Please note that the thermal comfort of occupants comprises various indi- The thermal behavior of the space heating is well explained based on the first law of thermodynamics. Let us consider n is the total number of zones in the building. For each zone i, (i = 1, ..., n), we denote the temperature of the zone by T i . The mass flow rate of the supply air entering in the i-th zone is represented by ṁi and the supply air temperature by T si . For the winter scenarios where the space heating is required, we write the energy balance equation using first law of thermodynamics as the i-th zone:

ρV i c p dT i dt = ṁi c p (T si -T i ) -Qh,load i (2.1)
with V i the volume of zone i, ρ the air density and c p the air specific heat coefficient. The rate of heat Qh,load i is the sensible heating load of zone i, that is, the net amount of energy that needs to be added to the zone to maintain a specified zone condition. Clearly, the sensible heating load is the sum of all heat losses and the internal heat gains. The heat loss for a zone i is due to the heat transfer from zone i to adjacent zones j and to the outside environment. Using thermal resistances, and denoting by q i all internal heat gains due to the occupancy, electronic devices and solar gain to zone i, T oa the outside air temperature, the sensible heating load for zone i is

Qh,load i = 1 R ext i (T i -T oa ) + n j=1,j =i 1 R ij (T i -T j ) -q i (2.2)
where R ij = R ji is the thermal resistance between zone i and zone j, and R ext i is the thermal resistance between zone i and the exterior of the building. Now, substitute the heating load into equation (2.1), we get

C i dT i dt = ṁi c p (T si -T i ) - 1 R ext i (T i -T oa ) - n j=1,j =i 1 R ij (T i -T j ) + q i (2.3)
where we set C i = ρV i c p as the thermal capacitance of zone i.

It is obvious that to maintain the space heating requirements, the user should manipulate the heat gain through the supply air balancing the sensible heat loads. This implies that to control zone temperature T i in the thermal comfort range, the user may manipulate the supply air temperature T si or the supply air mass flow rate ṁi . This is a key that differentiates the thermal behavior of the HVAC systems based on VAV and FCU units respectively which will be clear in the following sections.

VAV type HVAC system

In large scale non-residential and commercial buildings, the HVAC system must meet the varying needs of different spaces since different zones of the building may have different heating and cooling needs. In that respect, VAV systems were developed to be more energy-efficient and to meet the varying heating and cooling needs of different building zones. All the VAVs receive the supply airflow from a central air handling unit (AHU). Then, VAVs control the supply airflow into the zones by adjusting the damper position to maintain the thermal comfort. Examples of equipments are shown in Figure 2.2. We explain the working principle of a VAV type HVAC system with its mathematical model based on the general thermal behavior given in previous section.

Working Principle

A VAV type HVAC system works as follows [START_REF] Howell | Principles of Heating Ventilation and Air-Conditioning[END_REF]: each VAV terminal box in each zone receives primary air from a central AHU at the same constant temperature, called the supply air-temperature. The VAV terminal box has a primary-air damper which regulates the volume of hot or cold primary air delivered to the box according to the needs of the Figure 2.3: VAV type HVAC system layout spaces. Note that since each VAV box regulates its primary air volume independently, the total volume of primary air delivered by the central AHU varies according to the demands of the VAV boxes in the building. Therefore the central supply fan must vary its output in order to meet the needs of all VAV units. The speed of the central supply fan is consequently controlled to meet the changing demands of the building.

As shown in Figure 2.3, the AHU contains a mixer, a heating coil and a supply fan. The mixer mixes the fresh air and the return air from all the zones. The heating coil is a water to air heat exchanger, which controls the temperature of supply airflow by varying hot water flow of constant temperature, supplied by a boiler. A VAV box at each zone changes the supply airflow varying a damper position with a local PI type controller. We derive a mathematical model of the thermal behavior of the zone and AHU, that is effectively used in the control design.

Mathematical Modeling

It is clear that AHU is a central unit providing the supply air at a constant temperature and the VAV units are present at every zone controlling the supply flow rate to maintain the zone temperature in a required range. Let us refer to the general thermal behavior of the building described in (2.3). The supply air temperature T si are kept constant and supply air flow ṁi is manipulated to maintain zone temperature T i where i = 1, 2 . . . , n. Now, we write explicitly (2.3) for i = 1, ...n as follows,

C 1 dT 1 dt = -α 11 T 1 + 1 R 12 T 2 + 1 R 13 T 3 + ... + 1 R 1n T n -c p T 1 ṁ1 + c p T s ṁ1 + 1 Rext 1 T oa + q 1 C 2 dT 2 dt = -α 22 T 2 + 1 R 21 T 1 + 1 R 23 T 3 + ... + 1 R 2n T n -c p T 2 ṁ2 + c p T s ṁ2 + 1 Rext 2 T oa + q 2 . . . . . . . . . C n dTn dt = -α nn T n + 1 R n1 T 1 + 1 R n2 T 2 + ... + 1 R n,n-1 T n-1 -c p T n ṁn + c p T s ṁn + 1 Rext n T oa + q n (2.4) with α ii = 1 R ext i + n j=1,j =i 1 R ij , i = 1, ..., n
Note that, the supply air temperature T si for i = 1, 2, . . . , n is denoted by a single variable as T s . The n first-order differential equations in (2.4) represents the the model for the overall n-zones building.

The equation system in (2.4) in the matrix form reads as,

      C 1 0 • • • 0 0 C 2 • • • 0 . . . . . . . . . . . . 0 • • • 0 C n             Ṫ1 Ṫ2 . . . Ṫn       =       -α 11 1 R 12 • • • 1 R 1n 1 R 21 -α 22 • • • 1 R 2n . . . . . . . . . . . . 1 R n1 1 R n2 • • • -α nn             T 1 T 2 . . . T n       +       (T s -T 1 ) c p ṁ1 (T s -T 2 ) c p ṁ2 . . . (T s -T n ) c p ṁn       +        1 Rext 1 1 Rext 2 . . . 1 Rext n        T oa +       q 1 q 2 . . . q n      
(2.5)

As the diagonal matrix in the left-hand side of (2.5) is not singular, the above equation can be written in the state space form

ẋ = A v x + u 0 v I n -diag(x) T B v u v + G w + q (2.6)
where (.) T denotes transposition,

A v =       C 1 0 • • • 0 0 C 2 • • • 0 . . . . . . . . . . . . 0 • • • 0 C n       -1        -α 11 1 R 12 • • • 1 R 1n 1 R 21 -α 22 • • • 1 R 2n . . . . . . . . . . . . 1 R n1 1 R n2 • • • -α nn        B v = c p       C 1 0 • • • 0 0 C 2 • • • 0 . . . . . . . . . . . . 0 • • • 0 C n       -1 G =       C 1 0 • • • 0 0 C 2 • • • 0 . . . . . . . . . . . . 0 • • • 0 C n       -1        1 Rext 1 1 Rext 2 . . . 1 Rext n        x = T 1 T 2 • • • T n T u v = ṁ1 ṁ2 • • • ṁn T u 0 v = T s w = T oa (2.7)
and I n is the unity matrix of dimension n. The subscript v denotes the building system dynamics related to a VAV type HVAC system. For constant supply air temperature, which is typical to VAV based air-conditioning systems, the state equation is a bilinear controlled system,

ẋ = A x + diag(x v )B 1 u v + B 2 u v + G w + q (2.8)
where the expressions of matrices B 1 and B 2 are trivial. For our purpose, equation (2.8) is linearized around an operating point (x (0) , u (0) v ) and discretized with a sampling period t s to yield:

x(k + 1) = A v x(k) + B v u v (k) + Gd(k) y(k) = x(k) (2.9)
where A v , B v and G are the resulting discrete-time system matrices of appropriate dimensions. Vector d = [w, q T ] T is the disturbance which accounts for the outside temperature T oa and internal gains q. y is the vector of system outputs implying all the zone temperature measurements are available. Please note that we use the same variables to avoid cumbersome notations that is, x, u v and d in equation (2.9) denote now small variations around their operating point values.

Air Handling Unit

With reference to This implies that the return temperature in the duct is completely determined by:

T r = n i=1 ṁi T i ṁa (2.12)
Next, let T m denote the temperature at the output of the mixer. The mixer mixes the return air at temperature T r with fresh outdoor air of flowrate ṁoa at temperature T oa . Then, again, writing the energy conservation law for the mixer, we have:

ṁr T r + ṁoa T oa = ṁa T m (2.13)
The conservation of mass at the inputs and output of the mixer implies ṁr + ṁoa = ṁa . The return mass airflow rate ṁr is a fraction δ (0 ≤ δ ≤ 1) of the total mass airflow rate , i.e., ṁr = δ ṁa , which implies that ṁoa = (1 -δ) ṁa . Then, equation (2.13) reads as:

T m = δT r + (1 -δ) T oa = δ n i=1 ṁi T i ṁa + (1 -δ) T oa (2.14)
Note that T m depends solely on all zone temperatures T i and on the outside temperature T oa . Above relations are significant in the controller design assignments in coming chapters.

FCU type HVAC system

FCUs are widely used in various types of building topologies due to their simplicity. Most favored traits of the FCUs units are ease of installation and operation, lower noise levels and versatility in type of mounting (floor or ceiling) etc... Numerous forms of FCU units are available in the market such as the example shown in Figure 2.4. A typical FCU unit comprises components as a heating coil, supply fan, filters, mixer and noise attenuation and their forms can vary depending on their internal arrangement inside the FCU unit.

We consider a blow type FCU unit where, as the name suggests, the supply fan is placed before the heating coil (See Figure 2.5). To derive a mathematical model, it is necessary to understand the working principle of the FCU system which is explained next.

Figure 2.4: FCU unit

Working Principle

Let us refer to the general building topology based on FCU system is shown in Figure 2.5. Each zone is equipped with a temperature sensor, a FCU to supply airflow and a return air plenum. The return air plenum recirculates the fraction of the return air to the FCU. The mixer combines the return air flow from the plenum with the outside air. The supply fan maintains constant supply air flow through the heating coil. The heating coil is a water-to-air heat exchanger which maintains the supply air at required temperature by manipulating hot water flow from a boiler or heat pump. This temperature control is achieved with embedded PID controllers. Thus, the supply air of constant flow is circulated into the zone at a suitable temperature. Zone temperature is controlled by modulating supply air temperature directly with the heating coil. Supply air flow is viewed as the control signal for achieving the zone temperature. The supply air temperature is controlled through the heating coil by controlling hot water flow, depending on the temperature of space to which FCU serves.

Mathematical Modeling

To derive a mathematical model for a FCU equipped building system with n zones, a similar approach is extended as that VAV type HVAC system. In FCU systems, unlike VAV systems, supply air flow rate is kept constant and supply air temperature is varied to accommodate the thermal requirements of the zone. The supply air flows ṁi (i = Figure 2.5: FCU based HVAC system -four zones building 1, 2, . . . , n) are viewed as system parameters. The thermal energy balance equations from (2.3) are written explicitly for n zones as,

C 1 dT 1 dt = -α 11 T 1 + 1 R 12 T 2 + 1 R 13 T 3 + ... + 1 R 1n T n -c p T 1 ṁ + c p T s1 ṁ + 1 Rext 1 T oa + q 1 C 2 dT 2 dt = -α 22 T 2 + 1 R 21 T 1 + 1 R 23 T 3 + ... + 1 R 2n T n -c p T 2 ṁ + c p T s2 ṁ + 1 Rext 2 T oa + q 2 . . . . . . . . . C n dTn dt = -α nn T n + 1 R n1 T 1 + 1 R n2 T 2 + ... + 1 R n,n-1 T n-1 -c p T n ṁ + c p T sn ṁ+ 1 Rext n T oa + q n (2.15)
which, in matrix form, reads as,

      C 1 0 • • • 0 0 C 2 • • • 0 . . . . . . . . . . . . 0 • • • 0 C n             Ṫ1 Ṫ2 . . . Ṫn       =       -β 11 1 R 12 • • • 1 R 1n 1 R 21 -β 22 • • • 1 R 2n . . . . . . . . . . . . 1 R n1 1 R n2 • • • -β nn             T 1 T 2 . . . T n       +       T s1 T s2 . . . T sn       γ +        1 Rext 1 1 Rext 2 . . . 1 Rext n        T oa +       q 1 q 2 . . . q n       (2.16) with γ = c p ṁ, β ii = α ii + γ, i = 1, ..., n
The above system (2.16) is written in state space form as follows,

ẋ = A f x + B f u f + G w + q (2.17)
where,

A f =       C 1 0 • • • 0 0 C 2 • • • 0 . . . . . . . . . . . . 0 • • • 0 C n       -1        -β 11 1 R 12 • • • 1 R 1n 1 R 21 -β 22 • • • 1 R 2n . . . . . . . . . . . . 1 R n1 1 R n2 • • • -β nn        B f = γ       C 1 0 • • • 0 0 C 2 • • • 0 . . . . . . . . . . . . 0 • • • 0 C n       -1 u f = T s1 T s2 • • • T sn T w = T oa (2.18)
The subscript f denotes the building system dynamics related to the FCU type HVAC system. It is worth noting that the system (2.17) is linear in the original variables. This linear system equations are discretized with the sampling period t s and represented in the following discrete-time state space form,

x(k + 1) = A f x(k) + B f u f (k) + Gd(k) y(k) = x(k) (2.19)
where, A f , B f and G are discretized zero-order-hold (ZOH) matrices from A f , B f and G .

Concluding Remarks

As seen from the discussion of the demonstration sites, the VAV and FCU type HVAC configurations are popular due to their simplicity in operation and their energy efficient performances. We have presented the mathematical models for the VAV and FCU types of HVAC building systems. These mathematical models will be used to propose novel approaches to centralized and distributed model predictive control. The energy consuming components of these HVAC systems are mainly supply fan and the heating coil. Hence, the control objectives in the next chapters will be designed around these components. These derived models are generic enough to accommodate different types of buildings and weather conditions.

Chapter 3

Centralized control

3.1 Centralized Control -Fault Free Case

Introduction

Model Predictive Control (MPC) has been very popular in recent years and is very favored control technique in various fields e.g. chemicals, food processing, automotive, and aerospace applications. There are great recent developments in the MPC techniques as its basic framework provides simplicity in the understanding and is flexible according to the nature of the application. At a given time instant, MPC generates the numerical sequence for future manipulated variables by predicting the system behavior. The ability to handle the soft and hard system constraints and the robustness towards the disturbance make MPC more beneficial. Nevertheless, in the literature, considerable progress can be seen in the implementation of MPC for HVAC building systems. The recent publications provide various aspects of MPC applied to HVAC building systems e.g. for the energy efficient building operations, for the improvement in the reliability of HVAC equipments, maintenance of thermal comfort inside the buildings etc. As already discussed, sophisticated technological schemes concerning economic control methodologies are now being developed for large scale buildings, based on various control theories like predictive control, fuzzy control [START_REF] Afram | Modeling and Control Design of Residentila HVAC system for Operating cost reduction[END_REF] etc, which maintain thermal comfort while minimizing the operational energy consumption. On the other hand, there has been a noticeable progress in maintenance approaches of HVAC systems categorized as reactive, preventive and predictive maintenance that facilitates efficient working of HVAC elements [ash15]. In the available literature of a HVAC building control strategies relating to the economic operation and maintenance are dealt at different control layers. The noticeable works lack to address an economic HVAC building operation considering the awareness about the equipment maintenance as an objective.

In this chapter, a MPC based strategy is proposed to achieve this cumulative objective of minimizing the operational energy consumption as well as the maintenance cost. The first part of the cumulative objective is quite straightforward. This is obtained through formulating the cost function in the MPC by introducing the term related to the energy consumption by the HVAC equipments.

The related part of cumulative objective as minimizing the maintenance cost of HVAC equipments is achieved through ensuring the smooth functioning of equipment parts. Generally the setpoint trajectories are calculated by the controller and are implemented by the actuators. If these setpoint trajectories contain significant fluctuations, they certainly deteriorate the actuator performance leading increased wear and tear of the actuator parts. This may cause equipment failure, a loss in the control of acceptable indoor environment as well as an increase of the the system maintenance cost. Hence, we propose a methodology to generate the smooth setpoint trajectories by reducing the fluctuations while maintaining the economic operational interests. This proposed methodology is motivated from the one-norm regularization analysis extensively used in the statistics and in the machine learning field. It was originally introduced as the least squares and has been extended to the various statistical models. Moreover, it is also interpreted in other fields of geometry and convex analysis [START_REF] Boyd | Convex Optimization[END_REF].

Brief details about the proposed formulation of MPC describing the cumulative objectives are presented in the next sections. This approach is illustrated on a six-zone HVAC building framework using the weather data from Nancy, France for the winter season.

Maintenance-aware Economic Model Predictive Control

Various interpretations of the economic MPC are available as discussed in the literature survey. In this work, we formulate the desired objectives considering the benchmark HVAC systems provided in Chapter 1. The proposed MPC details are explained considering the VAV type HVAC building system. Note that, referring to the mechanism of FCU type HVAC configuration given in Chapter 1, the cost function formulation can be effortlessly extended for this type of HVAC system.

Mathematical Framework

There has been a great evolution in the formulation of MPC in the literature. We consider the discrete-time state space model (2.9) which reads as:

x(k + 1) = Ax(k) + Bu(k) + Gd(k) (3.1) y(k) = Cx(k)
where x(k) ∈ R nx are the states representing the zone temperatures, u(k) ∈ R nu are the inputs of the system denoting the supply airflow, n x and n u are number of states and number of inputs respectively. d(k) ∈ R n d are the disturbances consisting the weather temperature and heat flux due to occupants and k is the discrete time. A, B and G are system dynamic matrices with appropriate dimensions.

Further, we present the optimization problem for the standard MPC as below, minimize

U k J (U k , x(k)) subject to x (k + j + 1| k) = Ax ( k + j| k) + Bu (k + j|k) +Gd (k + j|k) j = 0, . . . , N -1 x min ≤ x ( k + j| k) ≤ x max j = 0, . . . , N -1 u min ≤ u (k + j|k) ≤ u max j = 0, . . . , N -1 x (k| k) = x(k) (3.2)
where J is the functional representing the overall cost function, N is the prediction horizon, U k = {u(k|k), . . . , u(k + N -1|k)} is the sequence of predicted control inputs at time k. The bounds u min , u max on the input vector u, i.e. on the supply airflow rate represent the damper limits in the VAV box. The bounds on states x min , x max represent the soft bounds on the zone temperature to maintain a thermal comfort. Problem (3.2) is solved repetitively at each time k using the receding horizon principle for the current measured state x(k) along with the predicted state variables {x ( k + j| k)} N j=1 . Assume that the forecast for the disturbances {d(k|k), . . . , d(k + N -1|k)} is available a prior. then the corresponding optimal sequence U k = {u (k) , ..., u (k + N -1)} is obtained and the first element u (k) of this sequence is applied to the system. The procedure is repeated at time k + 1, based on the new measured state x (k + 1).

Please note that an extensive discussion about the mathematical model of benchmark HVAC building system is presented in Chapter 2, which enables us to derive a discrete-time state space model and to formulate the constraints in the MPC problem. Furthermore, the cost function J (U k , x(k)) is the key element and so the novelty of the proposed MPC formulation.

Formulation of Cost Function

Recalling from the previous discussion, the control objectives for the HVAC building operation involve energy efficient performances without compromising the thermal comfort and a reduction in the maintenance cost of actuators. We propose the formulation of the cost function that (i) minimizes the economic operational cost, (ii) maintains the thermal comfort in the zones and (iii) generates smoother control actions by eliminating the fluctuations that improves the actuator life-time.

Energy efficient Operation

The energy efficient building operation is guaranteed by reducing the energy consumed by the building components. In a VAV type HVAC building system, supply fan and heating coil in the AHU are the main energy consuming components. Let us represent the mathematical formulation of the energy consumed by these equipment in order to formulate the cost function. Let J e be the total energy cost for a time interval

[t 0 , t f ] J e = (J h + J f an ) (3.3)
where J h and J f an are the energy associated to the heating coil and to the supply fan in the AHU.

(a) Energy consumed by the heating coil The power or heat transfer rate ( Qcoil ) in the AHU required at the heating coil to deliver an airflow at temperature T s is directly obtained from writing the energy conservation law

Qcoil = c p (T s -T m ) ṁa = c p (T s -T m ) n i=1 ṁi (3.4)
where T m is the temperature of air at the output of the mixer and ṁi is the supply airflow for i-th zone. c p is the air heat capacities. Then, the energy consumption due to heating is simply given by

J h = c 1 t f t 0 Qcoil dt (3.5)
where c 1 represents the related cost per kWh power consumption.

(b) Energy consumed by the supply fan The VAV requires a certain total mass airflow depending on each local (zone) heating load. This mass airflow is discharged by the power fan which is driven by a variable speed drive. The power fan characteristics for a AHU is given by a cubic law [START_REF] Howell | Principles of Heating Ventilation and Air-Conditioning[END_REF], that is:

Ẇfan = α n i=1 ṁi 3 (3.6)
where α is a proportionality constant. With the above power characteristics, the energy consumption for a supply fan is as follows,

J f an = c 2 t f t 0 Ẇfan dt (3.7)
where c 2 corresponds to cost per kWh kWh energy consumption.

Thus, the total energy demand from the AHU can be summarized from (3.3), (3.5) and (3.7). Recalling the system dynamics of VAV system from Chapter 2, the supply airflow rates ṁi for i = 1, 2, . . . , n are represented as the system input vector u. Considering a discrete-time setting, we discretize the integral J e with Euler method using sampling time t s . The total energy cost on an interval [t 0 , t N ] of N sampling intervals is given by,

J e e = N -1 k=0 k (3.8)
where

k = t k+1 t k c 1 c p (T s -T m ) n i=1 ṁi + c 2 α( n i=1 ṁi ) 3 dt (3.9)

Thermal Comfort

With the above economic cost function, it is crucial to maintain the thermal comfort for the occupants. In order to ensure the thermal comfort in the zones, the zone temperatures should be controlled in the comfort range of [x min , x max ].

x min ≤ x(k) ≤ x max (3.10)
These thermal bounds are enforced as the boundary constraints in the MPC formulation.

Actuator fluctuations Reduction [DYH16c]

The above cost function (3.8) and the constraint (3.10) address the energy optimization and thermal comfort aspects. In addition, we introduce a term which indirectly addresses the maintenance cost.This is achieved by reducing the oscillations in the setpoint trajectories calculated by the controller. Indeed, smooth setpoint trajectories reduce the fatigue in the actuators, lowering the system maintenance cost. This term is a regularization term that is formulated as 1-norm ( 1 ) over a total variation (TV) of control signal shown below,

re (U k ) = λ N j=1 u(k + j|k) -u(k + j -1|k) 1 (3.11)
where λ is a regularization parameter with λ > 0.

Here, we present a brief analysis of the above total variation (TV) term. The origin of 1 regularization term can be found in the field of machine learning and signal processing. This is achieved by adding 1-norm regularization term to the existing objective function to be minimized. This results in penalizing non zero elements in the solution. This has been used to develop several motivating examples in data analysis, compressed sensing, audio processing and biomedical expanding beyond the original signal recovery formulation [START_REF] Hastiea | The Elements of Statistical Learning[END_REF]. The introduction of the regularization in the control field is quite recent [START_REF] Gallieri | Lasso-MPC Predictive Control with l1 -Regularised Least Squares[END_REF].

regularization and Sparsity Property

Let us consider the standard 1 regularized optimization problem given as below,

z * = argmin z f (z) criterion + λ z 1 penalty (3.12)
where z ∈ R n is the vector of optimization variables.

For better understanding of 1 regularization, we consider an illustrative example where the criterion is

f (z) = z -z op 2 Q , z = R 2 , Q = [1, 0.1; 1, 1] and z op = [2; 0.4].
As shown in Figure 3.1, the least squares criterion has elliptical contours, centered at z op . Solutions z * were computed for increasing values of λ from 0 to ∞ and solution path is derived. This test is carried out for 1 and 2 regularization as shown in Figure 3.1a and 3.1b respectively. Let λ * be the value of λ, where in ( 1 ) regularization 2 in the solution equal to zero. Above example shows that 1 regularization yields the solution with some zero elements for λ > λ * , this characteristic is described as the sparsity of the solution.

The idea of 1 regularization has been extended in several ways among which an useful extension is obtained by replacing z 1 with Dz 1 where D ∈ R (N -1)×N n is the bidiagonal matrix

D =          -1 T n -1 T n 0 . . . 0 0 0 0 -1 T n -1 T n . . . 0 0 0 . . . . . . . . . . . . . . . . . . . . . 0 0 0 . . . -1 T n -1 T n 0 0 0 0 . . . 0 -1 T n -1 T n         
where 1 T n is the unity row vector of dimension n. This represents the variation in a variable z as ( (z k+j|k -z k+j-1|k ) 1 ) for j = 1, 2...N . This construction attenuates the variation in variable z hence it is also termed as Total Variation regularization [START_REF] Hastiea | The Elements of Statistical Learning[END_REF]. Total variation reconstruction assigns less penalty to large variations compared to small variations ensuring the smoothness of signal [START_REF] Boyd | Convex Optimization[END_REF].

In order to suppress the fluctuation in the control signal provided to the VAV box, an additional term in the objective is considered as follows,

re (U k ) = λ DU k 1 = λ N -1 i=0 u(k + i + 1) -u(k + i) 1 (3.13)
According to the previous argument, the total variation term in the above equation allows some elements r if r = [u(k + j|k) -u(k + j -1|k)] (j = 1, . . . , N ) equal to zero herewith reducing variations in setpoint trajectories provided to the damper. Moreover, the reduction of the fluctuations in the actuators boosts service time cycle of the damper and expands overall life-span with minimizing maintenance cost of a building. Now, the total cost from (3.8) and (3.13) reads as

J (U k , x(k)) = α e e (U k ) + α re re (U k ) (3.14)
where α e and α re are appropriate weights defined by the user. To implement the centralized MPC on a given HVAC building, we solve the optimization problem (3.2) with the cost function defined in (3.14) subject to the system dynamics as equality constraints with usual thermal comfort constraints and actuator limitations.

Simulation Results

We consider an illustrative building example as shown in Figure 3.2 with six zones. Every zone is equipped with a temperature sensor, a VAV box to supply airflow from AHU.

We derive a mathematical model of the thermal behavior of the zones and AHU, that is effectively used in control design. Finally, the proposed centralized MPC has been applied to this simulated six-zones building benchmark. The simulator uses the numerical data and operating point values shown in Table 3.1 and 3.2.

C i 4.5 kJ/s R ext 6 W/ • C R ij 18 W/ • C c p 1.005 kJ/kg • C T min 22 • C T max 24 • C ṁmin 0.192 m 3 /s ṁmax 0.42 m 3 /s Table 3.1: Numerical values T 0 oa 5 • C T 0 s 28 • C ṁ0 i 0.192 m 3 /s q 0 i 0.65 kW T 0 i 23 • C N p 24 h Table 3.2: Operating Point
The state space matrices are given in (3.15) for the sampling time t s = 15min.

A =          0.
         B =          3.
         G =          0.
                                                                                    (3.15)
The thermal comfort band is 23 • C ± 0.5. The supply airflow to each zone varies between 0.0192kg/s and 0.41kg/s with a constant supply air-temperature fixed at 26 • C. A test is performed during 5 days while two scenarios are simulated for comparative analysis of economic MPC, with and without TV regularization. The internal heat gain due to the occupants in each zone and the outside temperature variation during these 5 days are shown in Figure 3.3.

To avoid the possible conflict of infeasible solution in the available numerical solvers (YALMIP in this case) due to the hard thermal comfort as constraints mentioned in (3.10), we relax these constraints by a relaxation parameter ζ as shown below,

-ζ + x min ≤ x(k + j) ≤ ζ + x max (j = 0, . . . , N -1) (3.16)
where 0 ≤ ζ ≤ 0.5 and it is considered as an additional optimization variable by adding a penalty term tc k (ζ),

tc (ζ) = ζ 2 (3.17)
and the this relaxation parameter is treated as the additional optimization variable. This term (3.17) is added in the cost function (3.14) and this gives the total cost as, In first scenario, the classical economic MPC solves the problem considering the cost function concerned with the energy consumption over a prediction horizon of 24 hours. The optimal setpoint values for the supply airflows are computed to maintain the zone temperatures in the comfort range. Figures 3.4 and 3.5 show the temperature responses and the computed optimal setpoints for the supply airflows of the six zones from the AHU. As it is clearly seen that the fluctuations in the airflow setpoints may affect the damper life-span in the VAV box as it causes excessive movement resulting in the fatigue of the equipment. Eventually, it may speed up the wear and tear of the damper mechanism with louder noises and loss of effectiveness in the control and that influences the comfort of occupants inside the zone. Likewise, this decreases the period of systematic maintenance time cycle of the VAV box that increases the economic cost of the building operation.

J (U k , ζ) = α e e (U k ) + α tc tc (ζ) + α re re (U k ) (3.18)
Further, the second scenario is simulated that comprises the economic model predictive control with the TV regularization term. The results are shown in Figures 3.6 and 3.7. The airflow setpoint trajectories are smoothed by eliminating the small variations in the setpoint trajectories of the supply airflows. It is clear that small fluctuations are suppressed while large changes are preserved in the response, motivating stable functioning of the damper to avoid above consequences. As a result, it improves the life expectancy of the VAV components decreasing the service requirements and so the maintenance cost in the long run. Also, this ensures the thermal comfort of occupants, shown in Figure 3.6 being the zone temperature in the comfort range of 23 • C ± 0.5 with no compromise in the indoor environment. These simulation results show the efficiency of the proposed method of the building economic operation with an improvement in the actuator performance.

Concluding Remarks

This proposed MPC formulation is an effective strategy which accomplishes the objective of decreasing the overall building life-cycle cost, by minimizing operational energy consumption while increasing the overall actuator life-span. This strategy is based on a model predictive control scheme where the objective of minimizing the operational energy consumption is incorporated in cost function while constraints are managed to define thermal comfort of the zones. The key feature is to minimize the variation in the setpoints of actuators generated by the model predictive controller by adding one-norm regularization as a penalty term. Due to sparsity property of this regularization term, it assists the actuator setpoint trajectories to retain the large variations and to eliminate the small variations. This follows to decrease the wear and tear of dampers in VAV box increasing the service time cycle, leading reduced maintenance cost in long run.

Centralized Control -Faulty case

Fault Diagnosis Module [DYH16a]

Despite their great flexibility in space conditioning, HVAC building systems are prone to failures due to their complexity and embedded electronics. Such possible failures or Referring to the available literature, there has been some progress on the fault diagnosis (FD) aspects of HVAC systems. Basically, a fault diagnosis module which provides information about the occurred faults, and a controller redesign system to adapt the control law to the faulty mode. Except for few advancements in the fault diagnosis for HVAC building systems discussed in Literature review section, many works are still needed.

In this work, we address a common fault of damper stuck encountered in the VAV type HVAC building system. Recalling the working principle of VAV box from Chapter 2, it essentially controls the supply airflow entering into the zones by manipulating the damper position. Due to the stiction phenomenon which is the combining effect of stick and friction, the damper may slower the movements causing a permanent stick at particular positions in extreme conditions. This may adversely affect the thermal comfort of the occupants inside the affected zone. Also, it can cause the excessive energy consumption. Hence we focus on the detection and diagnosis of the VAV damper stuck fault for the provided benchmark HVAC building system. In the next sections, the thorough discussion about these two sequential steps and its illustration for the six zones VAV type HVAC building system is provided.

FDI based on Unknown Input Observer Design

Consider u is the optimal values of supply airflow calculated by the controller and to be implemented by the VAV box while u A is the actual supply airflow discharge, (see Figure 3.10(a)). Referring to the Figure 3.10(b), a dedicated UIO is designed for given thermal zone that generates the residual based on the controller setpoints u, disturbance information and the actual zone temperature. Let us recall the discrete time state space model for a VAV type HVAC building system from the previous section,

x(k + 1) = Ax(k) + Bu(k) + Gd(k) (3.19)
where x(k) ∈ R nx are the states representing the n x zone temperatures, u(k) ∈ R nu are the n u inputs of the system denoting the supply airflowS. d(k) ∈ R n d is the disturbance vector consisting of the weather temperature and the heat flux due to occupants and k is the discrete time. Please note that the benchmark building under consideration satisfies n x = n u = n, hence for the simplification, we will use n as the dimension. A, B and G are system dynamic matrices with appropriate dimensions. With the proposed scheme, simultaneous VAV damper stucks detection and isolation is made possible. A key point for the design of the unknown input observer design for this building system is that the outside temperature w = T oa and the occupancy q are known disturbance signals provided respectively by an external temperature sensor and the occupancy profile of the zones.Then, under damper stucks of the VAVs, the actual input to the plant can be modeled as

u A = u + f where signal f = f 1 f 2 • • • f n T
is the vector of faults on the n VAV dampers. Under these faults, the discrete time state space model from (3.19) of the building reads as

x (k + 1) = Ax (k) + Bu (k) + B f f (k) + Gd (k) (3.20) with B f = B u 0 1×n I n
where 0 1×n is the (1 × n)-dimensional zero matrix and I n the ndimensional identity matrix. Thanks to the structure of the control input u, the nominal free-fault model (3.19) can be written as

x (k + 1) = Ax (k) + B υ i υ i (k) + B (-i) u (-i) (k) (3.21)
for all i = 1, ..., n, where the vectors υ i (k) and u (-i) (k) are defined by:

υ i = u i , d T T u (-i) = u 1 , ..., u i-1 u i+1, ..., u n T (3.22)
and the matrices B υ i and B (-i) u by:

B υ i = b i | b 0 | B d B (-i) = b 1 | ... b i-1 | b i+1 | ... b n (3.23)
Now, viewing the vector u (-i) in model (3.21) as the unknown input to the i-th zone dynamics, the structure of the i-th residual generator for detection of the i-th damper stuck is therefore described by the dynamical system [JR99]

z i (k + 1) = F i z i (k) + T i B υ i υ i (k) + K i y (k) r i (k) = (I n -H i ) y(k) -z i (k) (3.24)
where matrices F i , T i , K i , H i are the parameters of the i-th generator to be determined.

The measured outputs being all the states (i.e., the zone temperatures) it follows that the existence of these parameters for the building model is always ensured as the two conditions below are satisfied 1. rank

I n B (-i) =rank B (-i) 2. (T i A, I n ) is detectable
In fact, it turns out that (T i A, I n ) has the stronger property of being observable. The parameters of the UIRGs for all the actuators i = 1, . . . , n are computed through the followings steps:

Algorithm 1 Computation of UIOs parameters Input Data: A; B (-i) u , i = 1, ..., n Results : F i , T i , K i , H i for i = 1 : n H i ← B (-i) B (-i) t B (-i) -1 B (-i) t T i ← I n -H i choose K i,1 such that F i = T i A -K i,1 is a stable matrix K i ← K i,1 + F i H i end for
The decision logic under a fault-free damper or a damper stuck in the i-th VAV is then simply given by e T n,i r i (k) < th i for fault-free mode e T n,i r i (k) ≥ th i for faulty case where th i is a threshold which, when exceeded, indicates the occurrence of a damper stuck and e n,i is the i-th unit vector (i = 1, ..., n) of the standard basis of the n-dimensional Euclidian space.

Estimation of the damper stuck fault

The above UIRGs scheme performs detection and isolation of VAV damper faults, however in order to achieve a complete diagnosis of a VAV damper stuck, an estimation of the stuck is needed. Consider the i-th filter (3.24) where the i-th control variable u i is explicitly exhibited:

z i (k + 1) = F i z i (k) + T i Gd (k) + T i b i u i (k) + K i y (k) r i (k) = (I n -H i ) y(k) -z i (k) (3.25)
Next, instead of the control signal u i , consider the filter (3.24) with the actual i-th input to the plant, i.e., the actual actuator output u A,i .

ζ i (k + 1) = F i ζ i (k) + T i Gd (k) + T i b i u A,i (k) + K i y (k) ρ i (k) = (I n -H i ) y(k) -ζ i (k) (3.
(k) = (I n -H i ) y(k).
The algebraic computation leads to the estimation u A,i of the i-th VAV damper stuck

u A,i (k) = (T i b i ) + × {ζ i (k + 1) -F i ζ i (k) -T i Gd(k)K i y(k)} (3.27)
where (T i b i ) + is the pseudo inverse given by (b

T i T T i T i b i ) -1 b T i T T i .
Practically, the estimation of the damper stuck is triggered at least one sampling period after the fault detection stage. If k denotes the time at which the estimation (3.27) is computed, then the fault occurred at or before time k -1. At time k, the damper stuck being a constant signal, we have u A,i (k) = u A,i (k -1) = u A,i . Then the estimation value at time k by (3.27) is equivalently given by

u A,i (k) = (T i b i ) + × { I n -H i |-F i (I n -H i ) -K i y(k) y(k -1) -T i Gd(k -1)} (3.28)
Further, a demonstration of suitably designed fault estimators is presented on the six zones benchmark building system in the next section.

Simulation Results

We consider the benchmark building shown in Figure 3.2. As discussed earlier, in this benchmark building, every zone is equipped with the temperature sensor and the VAV box to maintain the thermal comfort. PID Controllers are tuned for the closed loop control of the VAV dampers to adjust the airflow rates to accommodate for the disturbances and occupants schedule. The zone temperatures are maintained at the setpoints of 23 • C ± 0.5 controlling the suppply airflow rates. For demonstration purpose, damper stuck fault is considered in zone-1. Figure 3.11 explains this damper stuck scenario. At the beginning of the third day, in zone-1, dampers in the VAV box is considered to be stuck allowing constant supply air flow rate of 0.33kg/s, while dampers in the VAV boxes from another zones are considered in healthy condition. As shown in Figure 3.12, zone temperature behavior is presented. Due to damper stuck, temperature in the first zone changes abruptly and leaves the comfort zone, while in temperature in another zones are maintained in the thermal comfort range. As per the formulation, Unknown Input Observer receives the information about the current zone temperature, disturbances and the control input calculated by the MPC controller (see

F 1 =          
0.1393 -0.0325 -0.0275 -0.0058 -0.0272 -0.0028 -0.0915 -0.3894 -0.0089 -0.0599 -0.0612 -0.0036 -0.0858 -0.0090 0.1397 -0.0664 -0.0052 -0.0596 -0.0065 -0.0598 -0.0663 -0.0875 -0.0035 -0.0596 -0.1451 -0.1209 -0.0087 -0.0072 0.1273 -0.0006 -0.0067 -0.0071 -0.1191 -0.1191 -0.0005 0.3868

          T 1 =          
0.9925 -0.0489 -0.0504 0.0015 -0.0504 0.0008 -0.0489 0.0024 0.0025 -0.0001 0.0025 -0.0000 -0.0504 0.0025 0.0026 -0.0001 0.0026 -0.0000 0.0015 -0.0001 -0.0001 0.0000 -0.0001 0.0000 -0.0504 0.0025 0.0026 -0.0001 0.0026 -0.0000 0.0008 -0.0000 -0.0000 0.0000 -0.0000 0.0000

                                                           (3.
H 1 =          
0.0075 0.0489 0.0504 -0.0015 0.0504 -0.0008 0.0489 0.9976 -0.0025 0.0001 -0.0025 0.0000 0.0504 -0.0025 0.9974 0.0001 -0.0026 0.0000 -0.0015 0.0001 0.0001 1.0000 0.0001 -0.0000 0.0504 -0.0025 -0.0026 0.0001 0.9974 0.0000 -0.0008 0.0000 0.0000 -0.0000 0.0000 1.0000

          K 1 1 =          
0.3710 0.0664 0.0598 0.0073 0.1191 0.0071 0.0664 0.3877 0.0073 0.0598 0.1191 0.0071 0.0598 0.0073 -0.1414 0.0664 0.0071 0.1191 0.0073 0.0598 0.0664 0.0875 0.0071 0.1191 0.0596 0.0596 0.0035 0.0035 -0.1289 0.0005 0.0035 0.0035 0.0596 0.0596 0.0005 -0.3868

                                                           (3.30)
where

K 1 = K 1 1 + F 1 H 1 .
Note that, subscript 1 denoted the design matrices for the UIO concerning the VAV box in the zone-1. The residuals calculated using (3.24) and presented in the Figure 3.13.

Value of threshold th 1 is 0.5 and the fault occurrence is detected if the residual crosses the boundaries of the threshold. As indicated in Figure 3.13, residual leaves the band of threshold and fault can be detected easily. Also it is seen, residual for another zones remain unaffected and stays within the threshold limit. This plays key role in isolating the fault. This is an advantage of designing dedicated UIOs for every actuator, hence the isolation of fault becomes simpler. This also enables the to detect the multiple faults at a given time. Moreover, in order to estimate the stuck value, the detection of damper stuck plays a key role by the formulation in (3.28). Using this simple algebraic relation, the supply flow rate in the affected zone is calculated. This is verified for the simulation example where u A,1 = 0.3209kg/s. This approach is very convenient as it effectively detects and estimates the damper stuck fault without changing the actual configuration of the closed loop system.

Fault Adaptive Control Module [DYH16a]

Early fault diagnosis and the adaptive action against the fault become necessary as they certainly prevent excessive energy consumption. Sometimes the fault adaptive actions avoid the immediate impact of the faults retaining the occupants functional under the fault occurrence state. The immediate corrective action against the fault prolongs the normal building operation providing extra time for the maintenance team to engage the rearing and servicing tasks of the faulty equipments. Fault adaptive strategies work in twofold as: an FDD system which provides information about the occurred faults, and a controller redesign system to adapt the control law to the faulty mode. We have described FDD module in depth in the previous section for the fault under consideration of the VAV type HVAC building system. Nevertheless, we emphasize the fault adaptive method in this section for considered damper stuck fault.

Furthermore, our objective is to go beyond the mere recovery of control performance by considering explicitly the energy consumed by the overall building as a key performance indicator to be minimized under comfort constraints during faulty modes. We focus on the damper stuck fault in the VAV box for which we have investigated the diagnosis method in detail in the previous chapter. In particular, we show that by cascading an FDD module with a suitably designed model predictive controller, tightly connected to the FDD module, fault adaptive control to VAV damper stuck can be achieved and makes possible the reduction of energy consumption of the building at levels that are achievable for the faulty conditions. The model predictive controller computes in real-time the setpoint settings of the local controllers of the AHU and VAV boxes under new constraints adapted to the indicated faulty mode to minimize an overall predicted energy consumption cost-functional within the thermal comfort bounds. The fault-adaptive control scheme is illustrated on the same six zones benchmark building.

As we know, in the healthy condition, controller varies the supply air flow rate in order to maintain the zone temperature around setpoint. In case of damper stuck, it is no longer possible for controller to change the supply air flow rate. Hence, a corrective action is taken by varying supply air temperature. This is illustrated in the schematic shown in the Figure 3.14. As FDD diagnoses fault in damper and the position at which it is stuck, the controller reconfigures to calculate required supply air temperature setpoint for providing constant supply air flow rate through VAV.

MPC based FAC

A model predictive control problem which aims at minimizing the energy consumption in the nominal as well as in the faulty modes of the HVAC benchmark system is now

Nominal mode

Faulty mode : i-th VAV damper stuck

C 1 T min i ≤ T i ≤ T max i T min i ≤ T i ≤ T max i C 2 ṁmin i ≤ ṁi ≤ ṁmax i 1.
for the i:thVAV : ṁi = ṁstuck i 2. for all j = i : ṁmin

j ≤ ṁj ≤ ṁmax j C 3 T s = T (0) s T (0) s -εT η ≤ T s ≤ T (0)
s + εT η Table 3.3: MPC constraints Reconfiguration in the faulty conditions formulated. As the fault-adaptive control is achieved through online modification of the constraints on the decision variables under damper stuck failures occurrence, details are presented of these constraints and the energy cost functional.

Under the nominal and faulty modes, the system states and control variables are subject to constraints imposed by the control performance requirements and the VAVs damper limits.

Table 3.3 summarizes the MPC configuration when damper stuck fault occurs. Constraints C 1 relates to the comfort bounds, C 2 to the VAV dampers limits and constraints C 3 to the supply air temperature. Note that in the nominal mode, the supply air-temperature T s is constant as required for energy-efficiency of the HVAC system but under a faulty mode, this variable is used as a redundant control variable which is allowed to vary within an interval centered around its nominal value, T (0) s -εT η ≤ T s ≤ T (0) s + εT η , where T η is a constant temperature and ε is a scalar variable which is the magnitude of the relaxation of the constraint. The constraints C 2 with regards to the i-th VAV under the faulty mode is automatically set to the stuck value ṁstuck i by the FDD algorithm. The above constraints are easily translated in terms of the variation signals x, u of the linearized discrete-time model (3.2). The constrained finite-time optimal control problem at time k is formulated as

Minimize U k ,ζ,ε J (U k , x k ) + ρε 2 subject to x (k + j + 1| k) = Ax ( k + j| k) + Bu (k + j) +Gd (k + j) j = 0, . . . , N -1 x min ≤ x ( k + j| k) ≤ x max j = 0, . . . , N -1 u min ≤ u (k + j|k) ≤ u max j = 0, . . . , N -1 x (k| k) = x(k) (3.31)
where U k = {u (k|k) , ...u (k + N -1|k)} is the set of predicted control inputs at time k and d ( k + j| k) is the predicted disturbance at time k +j . The scalar ε, the magnitude of the relaxation of the supply air-temperature constraint, is treated as an optimization variable. Under the nominal mode, ρ is set to zero in (3.31) but whenever a fault is detected and isolated the optimization problem is adapted to the faulty mode by modifying accordingly the constraints in (3.31) with ρ = 0. Problem (3.31) is solved repetitively at each time k for the current measurement x k and for the corresponding mode at that time as indicated by the FDD algorithm along with the predicted state variables{x ( k + j| k)} N j=1 and the predicted disturbances d ( k

+ j| k) N j=1
. Then, the corresponding optimal se-quence U k = {u (k) , ...u (k + N -1)} is obtained and the first element u (k) of this sequence is applied to the system. The procedure is repeated at time k + 1, based on the new state x (k + 1) and the new information provided by the FDD. It follows naturally from this procedure a fault-adaptive control algorithm to able to cope with damper stuck failures.

Simulation Results

The FAC algorithm has been applied for the six zones VAV type benchmark HVAC building system from the Figure 3.2. We refer to the previous sections for detailed description of benchmark building, its mathematical model and the numerical data used.

We use the same internal gain profiles due to the occupants in a zone and weather temperature profile over a period of 5 days are shown in Figure 3.3. As discussed earlier, for nominal case, the thermal comfort band is 23.5 • C ± 0.5 and the supply airflow vary between 0.0192kg/s and 0.41kg/s with a constant supply air temperature of 26 • C. MPC solves the problem stated in equation (3.31) over a prediction horizon of 6 hours, where optimal setpoint values for supply airflows are computed to maintain zone temperatures in the comfort range. As shown in Figure 3.11, the damper in the VAV box of zone 1 is considered to be stuck at the beginning of the third day and the fault diagnosis module helps us to detect and to estimate the stuck value which is 0.33kg/s. This triggers the constraints modification for the supply airflow in zone 1 and the supply air-temperature from the AHU.

It is worth observing that the supply air-temperature is modified in a way such that that the temperature in zone 1 is effectively maintained within the comfort band whereas the supply airflow setpoints in other zones are automatically adjusted accordingly to this new supply air-temperature without affecting the thermal comfort in the respective zones, presented in the Figure 3.15, 3.16 and 3.17. This test shows that the fault-adaptive control algorithm successfully compensates for VAV damper stuck failures with the benefits of minimizing the energy consumption of the HVAC system under such failure modes.

Concluding Remarks

This work presents a potential approach to detect and diagnose a damper stucks in VAV boxes in the HVAC building systems. The novelty in the proposed approach is the design of the dedicated bank of unknown input residual generators for fault detection/isolation. This detection step is followed by the suitably designed fault estimators that evaluates the damper stuck value. Furthermore, we also present a novel fault-adaptive control method for the VAV damper stuck fault. The adaptive control module receives the information from the proposed fault diagnosis module. The key feature of this fault-adaptive method is to modify the constraints online in the nominal MPC controller. Though, in the general MPC-based fault-adaptive control, it is not a trivial problem to know how to change the problem formulation under the occurrence of faults, the constraint modification performed here under VAV damper stuck was made possible thanks to the designed fault estimation filters from the fault diagnosis module. Finally, the simulation results are provided for the six zone building benchmark building.

Chapter 4 Distributed Control and Estimation

A building can be a complicated system depending on the building topology, its use and the interconnections among the zones. Normally, a centralized controller for the whole building system results in an optimal energy performance maintaining the thermal comfort requirements. To obtain the optimal performance, it is necessary for the centralized controller to receive the building-wide sensor information and communicate the control signals back to the actuators installed in the building-wide network. Specially, if the controller is based on the MPC, the mathematical model with respect to the entire building structure must be provided a priori. However, in the case of large scale buildings these requirements to implement the centralized control may appear as a great challenge due to communication network failures, difficulty in deriving the centralized mathematical model, isolation during fault occurrences. Hence, in the recent years the theory of decentralized control and distributed control for large-scale HVAC building systems have attracted considerable attention of many researchers. The basic notion behind the decentralized and distributed control is to manage the subsystem performance instead of considering the overall system. In the decentralized control framework, the individual controllers maintain the respective subsystems performance neglecting the interaction information among the subsystems. If these subsystems interactions are significant and largely contribute in the overall building behavior then the decentralized control may result in suboptimal performance. But unlike decentralized control, in the distributed control, the communication between the subsystem controllers may ensure the optimal performance considering the interactions between the subsystems. This distributed control framework offers various advantages e.g. possibility of isolation of the subsystem in case of faulty events, modular maintenance scheme or plugging new subsystem in the existing distributed control framework etc.. Due to the obvious advantages of the distributed control framework, the Energy IN TIME project considers the importance of the scalability of the proposed energy efficient solutions given under the centralized schemes. In this chapter, we propose methodologies of distributed model predictive control (DMPC) for the benchmark HVAC building structures. A detailed discussion about the mathematical formulation of the proposed methods with their validation on the simulator platforms are presented in the next sections.

Optimality Conditions Decomposition based Distributed Model Predictive Control

Introduction

We consider the CMPC problem for the targeted large scale system is available. This CMPC comprises the knowledge of control objectives as cost function, system dynamics as equality constraints and physical limitations on actuators as bounds on the variables.

We propose an approach of distributed model predictive control that addresses the system decomposition and the coordination strategy between subsystem controller in two distinct steps. In the first step of system decomposition into subsystems, we formulate the Karush-Kuhn-Tucker (KKT) system of the CMPC. This KKT system is carefully analyzed and modified in an efficient manner to obtain separable subsystems. The resulting separable KKT subsystems allow to obtain the structure of the decomposed subproblems. The degree of modification in the original KKT system will define the convergence of the decomposed solution to the centralized solution. In second step, the coordination between the subsystem controllers is proposed based on the Optimality Condition Decomposition method. Finally, the proposed approach is illustrated using a benchmark HVAC building case study.

Motivation

We consider a discrete-time state space model representation for the building system as,

x(k + 1) = Ax(k) + Bu(k) + Gd(k) (4.1)
where x(k) ∈ R nx , u(k) ∈ R nu and d(k) ∈ R n d are the states, inputs and disturbances of the system (recall that there are as many input variables as state variables for the zones, as n x = n u = n). Let the formulation of the optimization problem for model predictive control for the above large scale system be, minimize

U k J (U k , x(k)) subject to x (k + j + 1| k) = Ax ( k + j| k) + Bu (k + j|k) +Gd (k + j|k) j = 0, . . . , N -1 x min ≤ x ( k + j| k) ≤ x max j = 0, . . . , N -1 u min ≤ u (k + j|k) ≤ u max j = 0, . . . , N -1 x (k| k) = x(k) (4.2)
recalling, U k = {u(k|k), . . . ,u(k + N -1|k)} is the predicted control input sequence while x(k) is the current measured state, {d(k|k), . . . , d(k + N -1|k)} is the forecast of the disturbances and N is the prediction horizon. Although this centralized systemwide MPC achieves the best attainable solution, it may suffer from some limitations. Computational burden and inadequate reliability are major drawbacks of centralized MPC [START_REF] Reverso | Distributed and Plug and -Play Control for Constrained Systems[END_REF]. Hence, as mentioned in the introduction, these drawbacks motivate different schemes for distributed control architectures. The class of proposed distributed architectures assume two distinct stages: i) decomposition of a large scale system into the subsystems ii) the coordination scheme between subsystem controllers. We propose novel approaches to address both the stages and obtain a innovative distributed model predictive control scheme.

The DMPC scheme proposed in this work is motivated by the decomposition technique for large scale nonlinear programming introduced by Conejo [START_REF] Antonio | Decomposition Techniques in Mathematical Programming[END_REF]. The basic idea behind this technique is based on the decomposition of optimality conditions of the original CMPC problem. Let us rewrite problem (4.2) into an equivalent general mathematical form to avoid cumbersome notations and simplify our systemwide CMPC presentation,

F : Minimize z f (z) subject to h(z) = 0 z min ≤ z ≤ z max (4.3)
where z ∈ R 2N n is the vector of the optimization variables, f is the objective function and the h = [h 1 , ..., h n ] T are equality constraints representing the dynamics of the system. To be more specific in translating (4.2) into (4.3), the optimization variables are defined as:

z =              u(0|0) x(1|0) u(1|0) x(2|0) . . . u(N -1|0) x(N |0)              (4.4)
for more concise notation, k = 0 is the present time and u(k + j|k) denotes the predicition of the signal u for k + j at time k. The bounds on the control variables and states translates as:

z min = 1 N ⊗ u min x min z max = 1 N ⊗ u max x max (4.5)
The system dynamics is described by the equality constraint as:

h(z) = Az -B = 0 (4.6) A =       B -I n 0 0 0 0 0 A B -I n 0 0 . . . 0 0 0 0A B -I n       B =       Ax(0|0) 0 . . . 0       +       G 0 0 0 0 G 0 0 . . . 0 0 0 G             d(0|0) d(1|0) . . . d(N -1|0)       (4.7)
where 1 N is vector of dimension N whose all components are equal to 1, the symbol ⊗ represent the Kronecker product. The vector B is known considering the forecast of disturbances is available. We assume the following for the problem (4.3): Assumption 1 The cost function f (z) is convex and twice differentiable Assumption 2 The constraints functions h(z) (h 1 , . . . , h n ) are linear. Now, we aim at deriving the distributed control problems from the given centralized problem in (4.3). To achieve this objective, we identify the following challenges: i) Decomposition of the vector z into p-subvectors as zi (i = 1, 2, .., p) ii) Defining the coordination between subproblems based on the coupling information.

System decomposition

The Lagrange function for the problem (4.3) with Lagrange multipliers λ (λ 1 , ..., λ n ) reads as,

L(z) = f (z) + n j=1 λ j h j (z) (4.8)
The bounds on the variables will not affect the decomposition so we ignore the associated inequalities. Thank to this Lagrangian, we define the KKT matrix by [START_REF] Boyd | Convex Optimization[END_REF],

KKT = ∇ 2 z L ∇ zλ L ∇ T zλ L 0 (4.9)
where its block-matrices are described in the next subsection. The proposed system decomposition method exploits the special structures in the KKT (Karush-Kuhn-Tucker) matrix. Essentially, the KKT matrix system is transformed into block-diagonal structure. These decomposed block structures help to identify the partition of the given vector z. This KKT matrix can be computed by various methods, it essentially represents the optimality conditions. Let us denote the KKT matrix as KKT cent for the centralized optimization problem (4.3):

KKT cent = KKT = ∇ 2 z L ∇ zλ L ∇ T zλ L 0 (4.10)

Significance of KKT matrix in the decomposition

A closed look at the KKT matrix from equation (4.10) reveals crucial information about the system structure. Due to the symmetric nature of the matrix, two distinct significant blocks are identified as ∇ 2 z L and ∇ zλ L. Their significance is given below:

1. Hessian of the Lagrangian Block w.r.t. z (∇ 2 z L) The upper triangular block from the KKT matrix (4.10) is the Hessian of Lagrange function with respect to variable z, i.e. (∇ 2 z L). This block contains information about the separability of the cost function with respect to the variable z. If ∇ 2 z L is in block-diagonal form, then the cost function is separable. In that case, every sub-block represents the structure of cost function for the decomposed subproblems.

Sensitivity matrix for the dynamics of the system (∇ zλ L)

The notation (∇ zλ L) denotes the hessian of L with respect to the vector z λ . This matrix represents the sensitivity of the system dynamics with respect to vector z

    ∇ z 1 h 1 (z) . . . ∇ z 1 h n (z) . . . . . . . . . ∇ zn h 1 (z) . . . ∇ zn h n (z)     (4.11)
The off-diagonal coefficients in this matrix represent the degree of coupling between the original components of vector z. This block holds the structural information about the system dynamics and plays a key role in achieving the distributed structure. For building applications, this block is sparse that makes it convenient to identify the group of variables that share strong interactions.

Formulating KKT matrix

The KKT matrix reveals the significant system structure which is exploited to derive the distributed structure. The standard ways to formulate this KKT matrix are: i) using developments in the primal dual interior point method or ii) gradient based approach [START_REF] Boyd | Convex Optimization[END_REF]. Both procedures are explained briefly in the following.

Primal-Dual Approach

Let us reconsider a general optimization problem (4.3) and its Lagrange function (4.8), where z and λ are primal and dual variables, respectively. We assume the problem is strictly feasible, then the KKT optimality conditions are written as,

∇L(z * , λ * ) = ∇f (z * ) + λ * T ∇h(z * ) = 0 (4.12) h(z * ) = 0 (4.13)
The above set of KKT conditions can be solved for (2N n) variables with to (2N n) nonlinear equations. Interior point method solves the equation (4.3) or KKT conditions (4.12) and (4.13) by applying the iterative Newton's method [START_REF] Boyd | Convex Optimization[END_REF]. This procedure seeks the search direction for both primal and dual variables. Hence, it is termed as Primal-Dual Interior Point Method. Nevertheless, the hierarchical procedure is described here to obtain KKT system as follows. Let us define the residual for the given iteration as,

r(z, λ) = r dual r prim = ∇f (z) + λ T ∇h(z) h(z) (4.14)
where r dual is dual residual, r prim is primal residual. Consider residual r(y) where y = (z, λ) with search direction ∇r(y), then a Newton's step ∆y = (∆z, ∆λ) is characterized by the linear equation r(y + ∆y) = r(y) + ∇r(y)∆y (4.15)

Now from (4.14) and (4.15) we obtain

∇ 2 f (z) + λ T ∇ 2 h(z) ∇h(z) ∇h(z) T 0 ∆z ∆λ = - r dual r prim (4.16)
The values of primal and dual search directions depend on feasible initial values of z and λ. Note that Newton's method can be extended in the case where the initial point is not feasible. Observe that the matrix in the left hand side of (4.16) is actually the KKT matrix.

Gradient Descent based Approach

An alternative approach to obtain the KKT matrix is based on the evaluation of gradient descent of Lagrange function (4.8). Let us reconsider the original problem (4.3) and its Lagrange function (4.8). We rewrite the first order optimality condition for Lagrange function at optimal point (y * ) at the given iteration of the gradient descent approach using the same notation y = (z, λ), ∇L(y * ) = 0 (4.17)

If y -1 is an initial feasible value, a step of gradient descent in the direction

∇ 2 L(y) is y -y -1 = ∆y = -(∇ 2 L(y)) -1 ∇L(y) (4.18)
Substituting y = (z, λ) and ∆y = (∆z, ∆λ) into (4.18) leads to,

∇ 2 z L ∇ zλ L ∇ T zλ L 0 ∆z ∆λ = - ∇ z L ∇ λ L (4.19)
The expansion of the above equation system (4.19) with the original Lagrange function from equation (4.8) gives the same result as the KKT system from equation (4.16). Therefore, the system given in (4.16) and (4.19) are equivalent.

Starting from the centralized KKT matrix (KKT cent ) in (4.10), it is always possible to transform it into the following equivalent form

KKT cent ∼             KKT 1 * . . . *
            (4.20)
where the symbol ∼ stands for equivalent. By the very transformation of KKT cent into the form (4.20), it naturally appears, in the transformed matrix, p KKT block-matrices on the diagonal. These KKT block matrices KKT i (i = 1, . . . , p) may be viewed as those of some p subsystems (S 1 , S 2 , . . . , S p ) composing the overall system. The (*) in (4.20) denotes sparse block matrices which bear the interaction between the subsystems. If these sparse block matrices are identically zero matrices then the KKT cent reduces to the following form

KKT dec =               KKT 1 0 . . . 0 0 KKT 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . KKT p-1 0 0 0 . . . KKT p               (4.21)
This KKT system (4.21) is regarded as a decentralized structure. Looking closely at (4.20), this equivalent centralized form may be viewed also as a distributed structure in which p subsystems are clearly exhibited with the interactions between these subsystems through the *-blocks. For this reason, we denote this matrix in (4.20) as KKT dist . Now we aim at deriving practically the KKT dist from KKT cent .

Transformation of KKT cent into KKT dist

It is worth to note that transforming the centralized KKT system into a distributed KKT system achieves the goal of decomposing the problem (4.2) on the control level as well as on the dynamics/constraints level. In the literature, various methods describe the transformation of a symmetric matrix into block-diagonal matrix form. For example, Dulmage-Mendelsohn [START_REF] Pothen | Computing the block triangular form of a sparse matrix[END_REF] technique or Reachability matrix techniques [START_REF] Siljak | Decentralized Control of Complex Systems[END_REF] are efficient methods in bipartite graph theory. But due to the large size of KKT matrix resulting from MPC problem (4.2), the complexity of bipartite graph methods increases significantly. Moreover, for computational purposes, these methods replace non zeros matrix coefficients by ones and transform the KKT matrix into an equivalent binary matrix. This may undermine the coupling information between the variables in the KKT system. The nested decomposition method could also be used [START_REF] Siljak | Decentralized Control of Complex Systems[END_REF]. This method is based on the graph theory and is very popular in the matrix decomposition literature. In this method, matrix elements that are less than are replaced by zeros. Then, the modified matrix is reordered to obtain a block diagonal form. Often, this procedure is carried out iteratively by augmenting for every iteration till the block diagonal form is achieved. Though this method seems promising, it compromises the coupling information linking system variables due to elimination of some matrix coefficients.

To avoid the loss of information, it is advised to reorder the original KKT matrix into the block-diagonal structure. In this paper, we use the Sparse Reverse Cuthill-McKee (CM) ordering algorithm [START_REF] Alan | Computer Solution of Large Sparse Positive Definite Systems[END_REF]. This algorithm permutes the rows/columns of sparse symmetric matrix to result in a narrow bandwidth band matrix. With this algorithm, it is possible to extract overlapping block diagonal matrices from the KKT matrix. This algorithm is illustrated with an example as shown in Figure 4.1, where M is the sparse symmetric matrix and p is the permutation of rows/columns by sparse reverse Cuthill-McKee method.

Algorithm 2 summarizes the approach of decomposition of the large-scale system into subsystems.

Algorithm 2 System decomposition based on KKT matrix Input Data: { h (h 1 , . . . , h n ); z 0 (z 0 1 , . . . , z 0 2n );f (z) } Result :z i (i = 1, . . . , p) 1. Formulate KKT cent matrix for F (z) 2. Transform KKT cent system into KKT dist using sparse reverse CM methods 3. Identify separable blocks from modified KKT dist matrix 4. Determine the sub-vectors zi (i = 1, . . . , p) based on the identified separable blocks

Optimality Condition Decomposition

With the provided system decomposition method, one can successfully partition a largescale system into p-subsystems. Now, raises the issue of synthesizing the controllers for the identified p-subsystems and their coordination. The main idea and tool used to address this issue is an optimization procedure suggested by Conejo. This approach is based on the so-called optimality condition decomposition which can be viewed as a modified Lagrangian relaxation.

Thanks to the decomposition of z into p-subvectors z = [z T 1 , zT 2 , . . . , zT p ] T , the centralized problem (4.3) can be written explicitly

F : Minimize z1 ,z 2 ,...,zp f (z 1 , z2 , . . . , zp ) subject to h(z 1 , z2 , . . . , zp ) = 0 zmin i ≤ zi ≤ zmax i (i = 1, 2, . . . , p) (4.22)
Since the constraints h(z 1 , z2 , . . . , zp ) = 0 is actually is the overall system dynamics, it naturally decomposes into p constraints describing the dynamics of the p subsystems S i (i = 1, 2, . . . , p) and how they interact with the other subsystems S j (j = i) (j = 1, 2, . . . , p). These p constraints h i (i = 1, 2, . . . , p) are actually coupling constraints that if they were removed, the resulting optimization problem would be easier to solve. The iterative procedure proposed by Conejo towards solving problem (4.24) known as the Optimality Condition Decomposition (OCD), is to enforce "some separable approximation" for both f and h i (i = 1, 2, . . . , p) by fixing some of the variables in these functions to their last computed values (indicated by overlining) to obtain the following p subproblems solved in parallel:

h(z 1 , z2 , . . . , zp ) = 0 ⇔          h 1 (z 1 ,
F i : Minimize zi f ( z1 , . . . , zi-1 , zi , zi+1 , . . . , zp ) + p j=1 j =i λj h j ( z1 , . . . , zi-1 , zi , zi+1 , . . . , zp ) (4.25) subject to h i ( z1 , . . . , zi-1 , zi , zi+1 , . . . , zp ) = 0 z min i ≤ zi ≤ z max i (i = 1, 2, . . . , p)
In (4.25), λj (j = i, j = 1, 2, . . . , p) denotes the dual variables (Lagrange multipliers) associated to the constraint h j and fixed at its last computed value. Note that the procedure looks like Lagrangian relaxation, but it has the computational advantage to perform a single iteration for each subproblem before updating the variables ( z1 , . . . , zp ). The benefits of this formulation rely in removing the coupling in the constraints from the original large-scale problem. We present the detailed steps of the Algorithm 3.

H.2

The Jacobian of the constraints ∇h has full row rank at z * H.3 The matrix (4.21), KKT dec , is non-singular at the optimal point z * The convergence of the distributed control solution to the centralized optimal solution is given by following proposition,

Proposition 1 Under the above assumptions (H.1, H.2 and H.3), if at the optimal point at z * it holds that

ρ * = ρ(I -(KKT dec ) -1 KKT dist ) < 1 (4.26)
then the procedure used in the Algorithm 3 converge locally to z * with linear rate at least equal to ρ * .

Proof 1 see Conejo [Con06]

Illustrative Example

An example to illustrate the proposed approach for better understanding is introduced.

Let us consider the following optimization problem

Minimize

z 1 ,z 2 -(20z 1 + 16z 2 -2z 2 1 -z 2 2 ) subject to z 2 1 + z 2 2 = 12.5572 z 1 + z 2 = 5 (4.27)
The centralized solution of the problem is,

z * = z * 1 z * 2 = 2.33 2.66 (4.28)
By comparing (4.28) with (4.3),

f (z) = -(20z 1 + 16z 2 -2z 2 1 -z 2 2 ) (4.29) h(z) = z 2 1 + z 2 2 -12.5572 z 1 + z 2 -5 = 0 (4.30)
and by applying Algorithm (2), the KKT system partition is graphically represented in the Figure 4.2. The following subproblems are formulated, minimize

z 1 -(20z 1 + 16z 2 -2z 2 1 -z2 2 ) + λ2 (z 1 + z2 -5) subject to z 2 1 + z2 2 = 12.5572 (4.31) minimize z 2 -(20z 1 + 16z 2 -2z 2 1 -z 2 2 ) + λ1 (z 2 1 + z 2 2 -12.5572) subject to z1 + z 2 = 5 (4.32)
where zi , λi (i = 1, 2) are initial feasible values. The subproblems (4.31) and (4.32) are solved using available numerical solvers e.g. interior point method. The solution of the subproblems is as follows,

z * 1 λ * 1 = 2.33 10.66 z * 2 λ * 2 = 2.66 1.7 (4.33)
The second order derivative of the constraint h (h 1 , h 2 ) is null matrix which leads to the optimal solution of subproblems is identical as optimal solution of centralized problem.

Simulation Results

In this section, we demonstrate the proposed approach on the complex building system of 6-zones as shown in the Figure 3.2. Each zone is equipped with a temperature sensor, a Fan Coil Unit (FCU) to supply airflow and a return air plenum. The return air plenum recirculates the fraction of the return air to the FCU. The extensive details about working principles of FCU unit are explained in Chapter 2. FCU controls the supply airflow into the zones to maintain the zone temperatures at the setpoint. The supply airflow temperatures are manipulated in the heating coil that is air-to-water heat exchanger (see Figure 2.5). A supply fan in the FCU units maintain the constant supply airflow into the zones. The heating coil is the energy consuming equipment in the FCU unit. We focus to formulate the optimization problem in the MPC to minimize the energy consumed by heating coil while maintaining the thermal comfort in the zones. The heat exchange between the zones can be described as the coupling between the subsystems. We initially define the CMPC problem, followed by the derivation of its KKT matrix. Then, with appropriate matrix transformation we modify this KKT matrix into block diagonal form. Finally, we extract the groups of variables from the modified KKT matrix to formulate the DMPC problems. The efficiency of this approach is validated with the implementation results of the extracted DMPC problems and its comparison with the centralized control performance.

CMPC Formulation

Mathematical model and cost function formulation are the requisites to implement the model predictive control. In this section we recall the centralized model predictive control problem described in detail in Chapter 2 referring to the mathematical model details presented in Chapter 1.

Mathematical Model

The system dynamics for the above benchmark building is based on the thermal balance equation of the zones. It is very straightforward to represent a discretetime state space model for the six zones FCU type HVAC building system from the thorough analysis given in Section 2.3. The A, B and G matrices of the model (2.19) of this benchmark building are:

A =          
0.5266 0.0664 0.0598 0.0073 0.0596 0.0035 0.0664 0.5266 0.0073 0.0598 0.0596 0.0035 0.0598 0.0073 0.5266 0.0664 0.0035 0.0596 0.0073 0.0598 0.0664 0.5266 0.0035 0.0596 0.1191 0.1191 0.0071 0.0071 0.4702 0.0005 0.0071 0.0071 0.1191 0.1191 0.0005 0.4702

          B =           3.
          G =           0.
          (4.34)

Formulation of the optimization problem

The control objectives is to minimize the consumed energy by the building with the constraints of keeping thermal comfort in the zones. Thanks to the particular control structure of the FCU units, the energy minimizing cost J e (3.9) has stage costs given by

l k (u(k)) = α1 T n u(k) (k = 0, 1, . . . , N -1) (4.35)
where α is a constant and 1 T n is the row vector of dimension n whose components are all equal to 1. The control sequence

u 0 =       u(0|0) u(1|0) . . . u(N -1|0)       (4.36)
can be expressed as u 0 = Ez using the padding matrix

E =       I 0 0 . . . 0 0 0 0 I . . . 0 0 . . . 0 0 . . . 0 I 0       (4.37)
so that the cost function reads simply as

J e = f (z) = c T z (4.38) where c T = 1 T n(N -1) E and 1 T n(N -1) = [1 T n , 1 T n , . . . , 1 T n ] N -1 . F : Minimize z f (z) subject to h(z) = Az -B = 0 z min ≤ z ≤ z max (4.39)
where the matrix A and vector B are given by (4.7). From this formulation, the KKT matrix is given by,

KKT = H A T A 0 (4.40)
where H is the hessian matrix of f (z) which is here the zero matrix of appropriate dimension. With regards to the coupling dynamics of the system, it is sufficient to look at this KKT matrix for an horizon N = 2 i.e.

KKT =           H 0 B T 0 -I T n A T 0 H 0 B T 0 -I T n B -I n 0 0 0 0 0 A B -I n 0 0 0           (4.41)

Implementation of DMPC based on optimality Conditions Decomposition

We can write vector z explicitly as:

z =      u(0|0) x(1|0) u(1|0) x(2|0)      (4.42)
where,

u(0|0) =           u 1 (0|0) u 2 (0|0) u 3 (0|0) u 4 (0|0) u 5 (0|0) u 6 (0|0)           u(1|0) =           u 1 (1|0) u 2 (1|0) u 3 (1|0) u 4 (1|0) u 5 (1|0) u 6 (1|0)           x(1|0) =           x 1 (1|0) x 2 (1|0) x 3 (1|0) x 4 (1|0) x 5 (1|0) u 6 (1|0)           x(2|0) =           x 1 (2|0) x 2 (2|0) x 3 (2|0) x 4 (2|0) x 5 (2|0) x 6 (2|0)           (4.43)
Now, we transform the KKT matrix from (4.41) using Sparse Reverse Cuthill-McKee (CM) ordering algorithm giving two separable blocks. These blocks provides the following partitions as:

z1 =                          u 1 (0|0) x 1 (1|0) u 2 (0|0) x 2 (1|0) u 5 (0|0) x 5 (1|0) u 1 (1|0) x 1 (2|0) u 2 (1|0) x 2 (2|0) u 5 (1|0) x 5 (2|0)                          z2 =                          u 3 (0|0) x 3 (1|0) u 4 (0|0) x 4 (1|0) u 6 (0|0) x 6 (1|0) u 3 (1|0) x 3 (2|0) u 4 (1|0) x 4 (2|0) u 6 (1|0) x 6 (2|0)                          (4.44)
To implement the distributed model predictive controllers for these two identified subsystems, the prediction horizon is N =24 hours. We have implemented the algorithm of distributed control structure given in Algorithm 3. The simulation results for the zone temperatures and supply airflow temperatures are shown in Figure 4.3 and 4.4 respectively.

Figure 4.3 shows the temperature in the three zones over five working days period. The zone temperatures and the control inputs for zone 1, 2 and 5 are similar to that of zone 3, 4 and 6 respectively. Hence, to avoid repetition, we have shown the zone temperatures and supply air temperatures for zone 1, 2 and 5. Thermal comfort is maintained between 23 • C ± 0.5 while occupants are present in the building. The corresponding control actions are presented in Figure 4.4 showing the supply airflow temperature setpoint trajectories that are implemented by the heating coils in the FCUs.

Comparative analysis

We compare the performance of the proposed DMPC algorithm with the decentralized and centralized approaches. The CMPC problem for proposed benchmark building can be implemented by solving the centralized problem (4.2). It can be clearly seen that, in CMPC, if the number of zones increases, the problem complexity and computational load increases exponentially. On the other hand, decentralized model predictive control can be applied by considering no communication between subsystem controllers. Also, the couplings between subsystems are ignored. This obviously deteriorates the performance with respect to DMPC and CMPC and results comparatively consuming more energy. The solution for centralized architecture and the proposed distributed architecture are very close as stated in proposition 1. Energy consumed by the system over 5 days is compared in Figure 4.5 and at the end of the fifth working day, total energy consumption is states in Table 4.1.

Concluding Remarks

In this section, a new approach of DMPC is proposed that realized in two steps. In first step, we decompose the system into subsystems by partitioning the KKT matrix obtained from the CMPC problem. In second step, the design of subsystem controllers and their coordination strategy are established using the Optimality Conditions Decomposition introduced by Conejo. Finally, this proposed control architecture is applied to the sixzones benchmark building. The comparison between the proposed strategy, decentralized control and centralized control is analyzed.

For the illustration purpose, we consider a benchmark HVAC building system with variable-air-volume (VAV) systems. We consider every zone is provided with a VAV box in which a damper manipulates the airflow of supply air with the constant temperature into the zones to maintain the thermal comfort. This supply air with constant temperature is provided by an Air Handling Unit (AHU). We propose the DMPC approach to achieve the same performance as centralized control architecture without compromising the thermal comfort of the occupants.

System Decomposition

Here, we address the issue of a system decomposition into subsystems that principally relies on the system dynamic behavior. Consider a building with n number of zones where every zone has its VAV box to provide the thermal comfort for the occupants (see Chapter 2 Section 2.2 for more details of the VAV type HVAC system). For each zone i, (i = 1, ..., n), we denote the temperature of the zone by T i , the mass flow rate at the output of the i-th VAV box by ṁi and the supply air temperature by T s . Then, the first law of thermodynamics applied to each zone is (as presented in (2.3))

C i dT i dt = ṁi c p (T s -T i ) - 1 R ext i (T i -T oa ) - n j=1,j =i 1 R ij (T i -T j ) + q i (4.45)
where C i is the thermal capacitance of zone i, R ij = R ji is the thermal resistance between zone i and zone j and R ext i is the thermal resistance between zone i and the exterior of the building. T oa is the outside temperature and q i is the heat flux due to occupancy and electronic devices. The system dynamics (4.45) is linearized around the operating point (x 0 , u 0 ) and discretized with a sampling period t s . This results in:

s i (x i , x j , u i ) = ẋi = a i x i + n j=1,j =i a ij x j + b i u i + g i w + q i (4.46)
where state x i is the zone temperature T i , input u i is the supply air flowrates ṁi and w is the outside temperature.

We propose an approach based on the global sensitivity of the system motivated by Sobieski [START_REF] Sobieszczanski-Sobieski | On the sensitivity of complex and internally coupled systems[END_REF]. He suggested to obtain the system sensitivity equations to evaluate the internal couplings and system behavior related to variable changes. This approach has been used for distributing the computing task of mathematical model design into various engineering disciplines in the 90s, especially for the aircraft wing design problems. Here, we extend this notion in decomposing the large-scale system into the subsystems.

Sensitivity Matrix

It is worth to note that the thermal balance equation for i-th zone as shown in (4.45) represents the change in the zone temperature with respect to the inputs i.e. supply flow rates (u i ) (i = 1, . . . , n) and temperatures of the neighboring zones x j (j = i, j = 1, . . . , n). For example, the coefficients a ij in (4.46) represents the sensitivity of the i-th zone temperature with respect to j-th zone temperature. The values of the coefficients b i represents the sensitivity of i-th temperature zone with respect to the i-th input (u i ). Note that inputs from the neighboring zones (u j ) will affect the the i-th zone temperature through j-th zone temperature (x j ). This will be accounted in the coefficient a ij . We write the sensitivity equations (4.46) for all n zones and represent them in matrix form as follows,

S gs =              ∂s 1 ∂x 1 ∂s 2 ∂x 1 . . . ∂s n-
             (4.47)
where the i-th column represents s i (i = 1, . . . , n) denoting the thermal balance for i-th zone. The rows represent the variables with respect to which the sensitivity is calculated. This sensitivity matrix contains the information about the system couplings with the states and the inputs. The following section explains the methodology to exploit this information in the system decomposition. The off-diagonal coefficients in this block matrix represent the degree of the sensitivity of the state variables x (x 1 , . . . , x n ) with respect to other state variables and inputs u (u 1 , . . . , u n ). The basic idea behind the decomposition is to partition the matrix (4.47) into p separable blocks. Every block will represent the group of zones representing the corresponding subsystem. The methodology of the matrix partition ensuring the minimal loss of information is explained in detail next.

Partitioning based on sensitivity

Note that, the sensitivity matrix obtained in (4.47) is a large-scale sparse matrix. There are various methods proposed in the literature to transform a sparse matrix into the block diagonal form [GL96] [START_REF] Pothen | Computing the block triangular form of a sparse matrix[END_REF]. In this work, we use the nested decomposition method [START_REF] Siljak | Decentralized Control of Complex Systems[END_REF]. This method is based on the graph theory and is very popular in the matrix decomposition literature. In this method, matrix coefficients that are less than are replaced by zeros. Then, the modified matrix is reordered to obtain a block diagonal form. Often, this procedure is carried out iteratively by augmenting such that ( k < k+1 )

where k represents the iteration till the block diagonal form is achieved. Let S k gs be the matrix after eliminating matrix elements less than k at k-th interval. This matrix S k gs is permuted to obtain a diagonal form S k gs using existing algorithms as e.g. reverse Cuthill-McKee algorithms [START_REF] Golub | Matrix Computation[END_REF]. To ensure minimal loss of the information in the modified sensitivity matrix S k gs the following condition should be verified [START_REF] Golub | Matrix Computation[END_REF],

(I -(S gs ) -1 S k gs ) < 1 (4.48)
where S gs is the sensitivity matrix S gs after applying the same permutation applied to the S k gs . denote spectral radius of the matrix. The condition (4.48) should be verified for each iteration. The detailed procedure of the decomposition of sensitivity matrix into block diagonal form is stated in the Algorithm 4.

Algorithm 4 Decomposition of global sensitivity matrix

Input Data: S gs , k Result : p-subsystems 1. Replace S gs (ij) by zero if S gs (ij) < k where i, j denote the number of row and column respectively. We obtain S k gs .

2. Permute S k gs system using sparse reverse CM methods into the matrix S k gs .

3. Verify the condition (I -(S gs ) -1 S k gs ) < 1 is satisfied.

4. Identify separable blocks from S k gs . If still is not trivial then augment k to k+1 and repeat step 1. 5. Otherwise determine the subsystems based on identified separable blocks from S k gs matrix

Sensitivity based cooperation in DMPC

Using sensitivity analysis, we decompose the system into p-subsystems. Next, we will investigate the issue of coordination between the controllers for corresponding subsystems. To address this issue, we introduce a sensitivity-based coordination mechanism [START_REF] Scheu | Sensitivity-based coordination in distributed model predictive control[END_REF]. This is achieved by adding the linear approximation of the objectives of the other subsystems into the objective of the local controller. We represent the problem (4.2) in simplified general mathematical form as mentioned in (4.3) recalling the definitions of variables from (4.6).

J : minimize z f (z) subject to h(z) = Az -B = 0 (4.49) z min ≤ z ≤ z max
From the subsystems derived using from the Algorithm 4, we have the p-subvectors from the vector z as z1 , . . . , zp .

Let us assume the cost function is separable with respect to the obtained partition as given below:

f (z) = f (z 1 , z2 , . . . , zp ) = f 1 (z 1 ) + f 2 (z 1 ) + . . . + f p (z p ) (4.50)
The equivalent representation for the constraints is:

h(z 1 , z2 , . . . , zp ) = 0 ⇔          h 1 (z 1 ,
z min ≤ z ≤ z max
The Lagrangian function for the above optimization problem is written as,

L = f (z) + λ T h(z) (4.53) = p i=1 f i (z i ) + λ T 1 h 1 (z 1 , z2 , . . . , zp ) + . . . + λ T p h p (z 1 , z2 , . . . , zp ) (4.54)
From the obtained partition of vector z, we group the Lagrangian functions as:

L = {f 1 (z 1 ) + λ T 1 h 1 (z 1 , z2 , . . . , zp )} + . . . + f 1 (z 1 ) + λ T p h p (z 1 , z2 , . . . , zp ) = p i=1 f i (z i ) + λ T i h i (z 1 , z2 , . . . , zp ) (4.55) If L i = f i (z i ) + λ T i h i (z 1 , z2 , .
. . , zp ) for i = 1, . . . , p then (4.55) becomes:

L = p i=1 L i (4.56)
Let Li be defined as,

Li = L i + p j=1 i =j ∂L j ∂ zi (z i -zi ) (4.57)
This function may be viewed as the perturbation of the Lagrangian term of the i-th subsystem L i (z i , λ i ) by the linear approximation of the Lagrangian terms of the other subsystems. Then the first order optimality condition applied to Li with respect to decision variables zi is,

∂ Li ∂ zi = ∂L i ∂ zi + p j=1 i =j ∂L j ∂ zi = 0 (4.58)
We aggregate the optimality conditions for i = 1, . . . , p,

p i=1 ∂ Li ∂ zi = p i=1 ∂L i ∂ zi + p j=1 i =j ∂L j ∂ zi = 0 (4.59)
Looking closely at (4.59), it appears that the optimality condition for (4.56) is equivalent to that of (4.59). Writing Li explicitly as,

Li = f (z i ) + λ T i h i (z 1 , z2 , . . . , zp ) + p j=1 i =j ∂f ∂ zi + p j=1 j =i λT j ∂h j ∂ zi (z i -zi ) (4.60)
This suggests the formulation of corresponding optimization problem of the i-th subsystem as: Note that the terms f ( z1 , z2 , . . . , zp ) is a constant and will not affect the optimal solution. Thus, f (z) can be approximated into separable functions as f 1 (z 1 ), . . . , f p (z p ).

J i : Minimize zi f (z i ) + ∂f ∂ zi + p j=1 j =i λT j ∂h j ∂ zi (z i -zi ) subject to h i ( z1 , . . . , zi-1 , zi , zi+1 , zp ) = 0 z min i ≤ zi ≤ z max i (4.

Simulation Results

We describe a benchmark school building as shown in Figure 4.6 used to demonstrate the proposed sensitivity based DMPC approach. The building has two floors with 16 zones having a total area 648m 2 . The cross sectional layout for the benchmark building is shown VAV terminal provides supply air flow to each zone in order to maintain the thermal comfort which is recirculated to the AHU. Please refer to the Chapter 2 where the VAV type of HVAC system are explained briefly. Thermal balance equations (4.46) are evaluated for every zone in the benchmark building. The numerical data used in simulating the case study school building is obtained from the Table 3.1. We assume that the number of subsystems or partitions p = 2 is decided by the user. The sensitivity matrix (4.47) is calculated and partitioned to obtain the subsystems as groups of zones. From the building layout presented in Figure 4.7, we obtain two groups p = 2 as {1, 2, 3, 4} and {5, 6, 7, 8}.

It is worth to note that, the decomposition obtained by this method is identical to the decomposition based in the physical properties (such as partitioning based on building floor or architecture, etc). So, the sensitivity based system decomposition can be viewed as a mathematical explanation of the manner of partitioning the building system dynamics. Applying Algorithm for sensitivity based DMPC controller is simulated considering the obtained subsystems. The occupants are present in the building from 8.00AM to 6.00PM, the corresponding heat flux due to occupants is shown in Figure 4.8. Also, the winter weather profile is plotted over the period of five working days. To represent the per- formance of the implemented architecture, the temperature response and corresponding supply airflow for all the are shown in Figure 4.9 and Figure 4.10, respectively. Note that the thermal range and actuators limits can be different for each zone and it is convenient to define in the proposed distributed architecture. The performance of the sensitivity based DMPC is equivalent to the CMPC framework. This is clear from the temperature and supply airflow behavior of the respective control architectures. We also compare the control performances with decentralized MPC architecture [START_REF] Siljak | Decentralized Control of Complex Systems[END_REF]. In decentralized MPC, local subsystem controller operates independent of the other subsystems i.e. without any coordination or data exchange between the controllers. Also, the dynamics of the local subsystems in decentralized control architecture completely ignores the coupling with the other subsystems. To support the argument, we also compare the energy consumed by the benchmark HVAC building over the five working days in the Figure 4.11. As, in the building system, the coupling between the zones are effective and if ignored, it results in consuming more energy and poor control performance.

Concluding Remarks

In this section, we propose an approach that addresses two stages of DMPC for the VAV based HVAC building system as i) the decomposition of the building system into subsystems and ii) the coordination between obtained subsystem controllers. The proposed method of system decomposition uses the sensitivity matrix derived from the system dynamics. This sensitivity matrix contains the coupling information between the variables. Based on this information, we have partitioned the sensitivity matrix into separable blocks that represent the corresponding subsystems. Moreover, the second stage of defining the coordination between the subsystems controllers is based on the notion of sensitivity of other subsystem controllers with respect to the local controller. This is achieved by adding a linear approximation of the objectives of other subsystems to the objective of the local controller. The decomposition and coordination strategies are demonstrated on the VAV type HVAC building system.

Distributed Estimation

Introduction

Under the Energy in TIME project, the overall building life cycle cost is minimized with the energy efficient operation by means of economic control and fault adaptive approaches as discussed in earlier chapters. We have demonstrated the energy efficient control method coupled with early fault diagnosis and adaptive control to achieve the cost-effective building operation. In large-scale systems, the complexity of fault diagnosis increases exponentially [START_REF] Farina | Distributed moving horizon estimation for sensor network[END_REF]. The number of required instrumentation and communication network grows with the size of the system. A great amount of attention is needed to address these issues in case of large scale HVAC building systems. In recent years, Moving Horizon Estimation (MHE) has been an important research topic in optimal state estimation of dynamical systems. MHE aims at estimating the states and outputs of the system with the knowledge of the system dynamics and the past available data. There are many variations and developments in this technique presented in the literature. In this section, we tackle MHE problem for the large-scale building applications. We propose a novel method for distributed MHE that decomposes the centralized MHE problem into the several local MHE problems. A coordination strategy between these local MHE problems is proposed which is based on the coupling information in the system dynamics. This method essentially relies on the partitioning of the optimality conditions of the centralized MHE problem. We demonstrate this proposed method on the six-zones benchmark building.

Moving Horizon Estimation Formulation

Consider the following discrete-time state space system,

x(k + 1) = Ax(k) + Bu(k) + Gd(k) + ω(k) (4.65) y(k) = Cx(k) + η(k)
where x ∈ R nx , u ∈ R nu and y ∈ R ny are states, inputs and outputs, respectively. The disturbances d ∈ R nd are known. Sensor noise is denoted by η ∈ R ny while ω represents the uncertainty in the system dynamics. A, B, C and G are the system matrices with the appropriate dimensions and k is discrete time instant.

The detailed MHE formulation for the above generalized state space system is presented in the following according to [START_REF] Christopher | Constrained linear state estimation -a moving horizon approach[END_REF]. The general form of the optimization problem for MHE reads as follows, minimize

X J(X, U m , Y m ) subject to x(j + 1) = Ax(j) + Bu m (j) + Gd m (j) j = 2, . . . , N y m (j) = C x(j) j = 1, . . . , N ω(1) = x m (1) -x(1) (4.66)
where the uncertainty ω and the sensor noise η are considered bounded The past control actions and output measurements over the horizon N is denoted by U m and Y m . The estimated state sequence is denoted by X.

U m =       u m (1) u m (2) . . . u m (N )       Y m =       y m (1) y m (2) . . . y m (N )       X =       x(1) x(2) . . . x(N )       (4.67)
where subscript m corresponds to the available data. The initial state value x m (1) is approximately known while d m (j) (j = 1, . . . , N ) denotes the available disturbance measurements. The solution of optimization problem (4.66) is the sequence X from which we consider only the last element x(N ) discarding other elements in the solution sequence.

Formulation of cost function

The cost function is multiobjective and aims i) to minimize the residual between measured and estimated outputs ii) to minimize the uncertainty in the system dynamics and in the measurements and iii) to minimize the uncertainty in the initial condition x(1). The formulation for each objective is as follows, where x 2 P denotes x T P x and P is a positive definite matrix. The cost (4.68) can be represented as,

η = P 1 X -(-Y m ) 2 P N (4.69)
where P 1 and P N are give as follows:

P 1 =          -C 0 0 0 0 0 -C 0 0 0 . . . 0 0 0 -C 0 0 0 0 0 -C          P N =         
P 0 0 0 0 0 P 0 0 0 . . . 0 0 0 P 0 0 0 0 0

P          (4.70)
ii) Error associated to input disturbances As mentioned in (4.65), disturbance ω are unknown but bounded. To minimize the error associated to input disturbances the following objective is considered:

ω = N j=1 x(j + 1) -Ax(j) -Bu m (j) -Gd(j) 2 Q (4.71)
where the weight matrix Q associated to the disturbance uncertainty level. To represent the above cost function in quadratic form as:

ω = Q 1 X -Q 2 2 Q N (4.72)
where Q 1 , Q 2 and Q N are give as follows:

Q 1 =          -A I n 0 0 0 0 0 -A I n 0 0 0 . . . 0 0 0 -A I n 0 0 0 0 0 -A I n          Q N =          Q 0 0 0 0 0 Q 0 0 0 . . . 0 0 0 Q 0 0 0 0 0 Q          Q 2 =          B 0 0 0 0 0 B 0 0 0 . . . 0 0 0 B 0 0 0 0 0 B          U m +          G 0 0 0 0 0 G 0 0 0 . . . 0 0 0 G 0 0 0 0 0 G          D m (4.73)
where I n is the identity matrix of dimension n. Note D m can be defined with the same notion as in (4.67)

D m = [d T (1), . . . , d T (N )] T .
iii) Error associated to the unknown initial condition To penalize the error in the initial condition,

icx = x(1) -x(1) 2 R (4.74)
Transform (4.74) into the explicit form:

icx = R 1 X -R 2 2 R N (4.75) R 1 =       -I n 0 0 0 0 0 . . . . . . . . . 0 0 0       R 2 =       -I n 0 . . . 0      
x m (1) (4.76)

R N =       R 0 0 0 0 0 . . . . . . . . . 0 0 0      
Collecting all previous objectives, the cumulative cost function J from (4.68), (4.71) and (4.74) can be written as,

J = η + ω + icx = P 1 X -P 2 2 P N + Q 1 X -Q 2 2 Q N + R 1 X -R 2 2 R N (4.77)
Now we write the optimization problem (4.66) associated to MHE for the system (4.65) with the cost function (4.77) as, Φ : minimize

X J(X) subject to g(X) = AX -B = 0 (4.78)
where,

A = Q 1 P 1 B = Q 2 -Y m (4.79)
The optimization problem (4.78) is a quadratic programming problem that can be efficiently solved with available numerical solvers. But in case of large scale systems, the size of the centralized MHE increases exponentially with increasing number of variables. To deal with the curse of dimentionality, we propose a decompsoition of the centralized MHE into different small sized MHE problems. A coordination strategy is introduced between the decomposed MHE problems so that the distributed MHE solution is comparable with the centralized MHE solution.

Distributed Moving Horizon Estimation

As in the distributed model predictive control proposed in Section 4.1, the distributed estimation is carried out in two distinct stages as: i) System decomposition based on the KKT matrix decomposition ii) Coordination between subsystem estimators. As described earlier, the system decomposition into subsystems is achieved by partitioning the KKT matrix of centralized MHE problem into p blocks. These blocks represent the corresponding subsystems. To coordinate the estimators of these subsystems we use the algorithm [START_REF] Antonio | Decomposition Techniques in Mathematical Programming[END_REF] based on Optimality Condition Decomposition.

System Decomposition

The Lagrange function for the problem (4.78) with Lagrange multipliers λ (λ 1 , ..., λ n ) can be written as, L (X) = J(X) + λ T g(X) (4.80)

The bounds on the variables will not affect the decomposition so we will ignore the associated inequalities. The KKT matrix for the problem reads as [START_REF] Boyd | Convex Optimization[END_REF],

KKT M HE = ∇ 2 z L ∇ zλ L ∇ T zλ L 0 (4.81)
With regard to system decomposition, it is sufficient to consider N =2, and the KKT matrix becomes,

KKT M HE =    H Q T 1 P T 1 Q 1 0 0 P 1 0 0    =      H 1 -A T C T H 2 I T n 0 -A I n 0 0 C 0 0 0      (4.82)
It is worth noting that, the cost function is separable and the interactions are present in the constraints. Using the Algorithm 2, we decompose this KKT matrix (4.82) into p block diagonal blocks. (For more details please refer to Section 4.1.3). These diagonal-blocks represent the partition of vector X into p-subvectors as X 1 , X 2 , . . . , X p .

Distributed Estimation scheme using Optimality Condition Decomposition

We obtain the distributed estimator structure from the available centralized moving horizon estimation problem (4.78) based on the Optimality Condition Decomposition method. Thanks to the system decomposition, The cost function and the constraints can be written as:

J(X) = J(X 1 , X 2 , . . . , X p ) (4.83) g(X) = g(X 1 , X 2 , . . . , X p ) = 0 ⇔          g 1 (X 1 , X 2 , . . . , X p ) = 0 . . . g p (X 1 , X 2 , . . . , X p ) = 0 (4.84)
Let us write the decomposed i-th estimation as:

Φ i : minimize X i J( X1 , . . . , Xi-1 , X i , Xi+1 , . . . , Xp ) + p j=1 j =i λT j g j ( X1 , . . . , Xi-1 , X i , Xi+1 , . . . , Xp ) subject to g i ( X1 , . . . , Xi-1 , X i , Xi+1 , . . . , Xp ) = 0 (4.85)
where ( X1 , . . . , Xk+1 , . . . , Xp ) and λj (j = 1, ..., p; j = k) represents feasible initial values.

Concluding Remarks

In this section we presented a method for DMHE for linear systems. We formulated the CMHE by delineating the associated centralized cost function. With this formulation of CMHE, the optimality conditions are written in order to formulate the KKT matrix. This KKT matrix is partitioned into separable blocks where these blocks provide the description of the subsystems. The distributed scheme of estimation is implemented using the optimization procedure suggested by Conejo. To illustrate the proposed method, simulation results are presented for the six-zones HVAC building system considering the weather data for Nancy France.

Chapter 5

Conclusions and Future Perspectives

This thesis work was performed within the Energy IN TIME project that is focused to reduce the energy consumption in the operational stage of existing non-residential buildings. Four questions were raised and addressed in this work: i) the issue of dynamic modeling of building with different types of air conditioning terminal units ii) the issue of centralized thermal comfort control of the buildings iii) concerns related to fault diagnosis and adaptive control to the equipment failures and iv) the issue of distributed control of large-scale buildings. This chapter describes the main original contributions in solving the above raised issues and it gives suggestions for future work, mainly oriented towards the improvements which might be required for real-buildings implementations. Firstly, we have established a benchmark building prototype based on the common VAV or FCU type HVAC units to maintain the thermal comfort inside the building. Mathematical models for the given building frameworks based on the thermodynamic behavior of the zones have been devised. The key feature of these models is that they are generic enough to accommodate different building thermal dynamics. They are well suited to test the control solutions provided in this work. For the simulation purposes, we assume that the forecast of the weather and the occupancy profile are available a priori. The winter weather data for the location Nancy, France is used for the simulation purposes. The typical office type occupancy profile is used where occupants are considered to be working from 8 AM to 6 PM with the two hours lunch break starting from 12 PM.

Next, we have addressed the issue of the economic operation of the buildings by proposing a novel approach to model predictive control. In the literature, usually, model predictive control for buildings considers standard quadratic objectives providing thermal comfort to the occupants. No much work has been reported on minimization of energy consumed by HVAC equipments. In this work, we have tackled the energy minimization problem through an economic perspective under the constraints of thermal comfort of the occupants. An innovative formulation of the overall control problem has been stated with an additional objective related to the awareness of the control strategy to the HVAC equipment maintenance. To the best of our knowledge, this formulation is new. This objective is achieved by adding the 1 regularization term that eliminates the fluctuations in the setpoint trajectories provided to the actuators. This reduces the wear and tear of HVAC equipments and improves their overall life span. The proposed centralized model predictive control strategy has been validated on the six zones benchmark building system. An important issue in the building operation is the occurrences of faults that may affect the economic performance of the building operation by consuming excessive energy and by compromising thermal comfort. One class of such faults is the actuator failures within HVAC terminal units. In this work, we dealt with the common fault in VAV boxes that is the damper stuck due to friction and stiction. If this damper stuck fault is not diagnosed and corrective actions are not taken on time, then such faults may impact severely the energy consumed and the thermal comfort in the building. We have proposed a fault detection and isolation module based on a dedicated unknown input observers bank. Moreover, to estimate the damper stuck value we have devised a new methodology for such estimation. To take corrective actions, a novel fault adaptive control approach based on model predictive control reconfiguration is designed. This fault adaptive module receives the fault diagnosis information and the controller achieves the fault tolerance through the constraint adaptation without compromising the thermal comfort of the occupants. The efficiency of the proposed fault diagnosis and adaptive control modules is shown with the simulation results on the six zones benchmark building.

The scalability of the proposed solutions under centralized control is addressed through the design of distributed architecture. A novel approach to distributed model predictive control is proposed based on the Optimality Condition Decomposition methodology. This approach of distributed model predictive control is realized in two steps: i) decomposition of the system into subsystems ii)the design of subsystem controllers and their coordination. A deep insight into the structure of KKT matrix allows us to decompose the system into subsystems. The decomposed blocks of this KKT system represent the subproblems in the distributed control architecture. To implement the MPC for the subsystems and establish their coordination, we used the Optimality Condition Decomposition algorithm introduced by Conejo. This approach is based on the modified Lagrangian relaxation that ensures the distributed solution is equivalent to the centralized solution. Under some conditions, simulation results on a FCU type HVAC building system show the potential of this distributed control scheme at approximating closely the centralized control performance.

A second approach to distributed model predictive control based on the sensitivity analysis is presented. As in the first DMPC proposed approach, the sensitivity-based distributed model predictive control proceeds in two steps: i)system decomposition into subsystems and ii) design and coordination of subsystem controllers and coordination between the controllers. The decomposition stage uses the sensitivity information from the system dynamics. This sensitivity information presented in the matrix form is then partitioned into several blocks where each block denotes the structure of the corresponding subsystems. In the second stage, the controllers for the subsystems are designed based on the sensitivity with respect to other controllers. Simulations on the benchmark building systems based on VAV type of HVAC systems are given to illustrate the performance of this scheme.

Finally, an estimation technique aiming at addressing the fault detection and isolation issues for the distributed frameworks is proposed. This distributed estimation scheme relies on moving horizon techniques and uses the optimality condition decomposition tool already used in distributed model predictive control. The centralized moving horizon estimation problem formulates a KKT matrix. We use the same KKT system decomposition technique to identify the subsystems. Once subsystems are identified, the dis-tributed estimation is realized by using Optimality Condition Decomposition Algorithm. To demonstrate this estimation approach, we have considered the case in which few measurements are available in the benchmark building system and shown that all the outputs are estimated using the proposed approach.

Despite the proposed solutions to the economic building operations in this work, there are still several issues which need to be addressed. The simulation validation has been limited to the terminal units with certain HVAC configurations, the proposed solutions can be applied to the central equipments as well as other types of HVAC configurations. In the case of VAV type terminal units, we have used the linearized models for the synthesis of control problems. An interesting approach can be considering the bilinearity in the VAV type terminal units. In this work, we assume the knowledge of weather forecast and occupancy schedule is available. It may be interesting to consider the unavailability of the forecast and occupancy schedule as future research topics.

The novel approach proposed for distributed model predictive control performs the decomposition of the large-scale system into subsystems. It also designs the subsystem controllers and their coordination. However, the issues related to the stability and controllability properties still need to be addressed. This distributed model predictive control approach has been extended to the distributed moving horizon estimation. Nevertheless, an interesting future topic would be to adapt the distributed moving horizon estimation techniques in the field of fault detection and isolation.

The proposed solutions under the project e.g. energy efficient control methods, fault diagnosis, and fault adaptive techniques etc. are built to suit these benchmark prototypes. Nevertheless, adequate flexibility in the proposed solutions is preserved so that it will apply to a broad range of building types. Further, we present the VAV and FCU type HVAC configurations and their mathematical models. These mathematical models will be used to propose novel approaches to centralized and distributed model predictive control. The energy consuming components of these HVAC systems are mainly supply fan and the heating coil. Hence, the control objectives in the next chapters will be designed around these components. These derived models are generic enough to accommodate different types of buildings and weather conditions.

Centralized Control :

A novel maintenance-aware model predictive control (MPC) method is proposed in this work. In this proposed control scheme, the MPC is designed considering multiple objectives as minimization of consumption of energy, maintenance of thermal comfort and reducing the maintenance cost of HVAC equipments. The first part of the cumulative objective is quite straightforward. This is obtained through formulating the cost function in the MPC by introducing the term related to the energy consumption by the HVAC equipments. The related part of cumulative objective as minimizing the maintenance cost of HVAC equipments is achieved through ensuring the smooth functioning of equipment parts. Generally the setpoint trajectories are calculated by the controller and are implemented by the actuators. If these setpoint trajectories contain significant fluctuations, they certainly deteriorate the actuator performance leading increased wear and tear of the actuator parts. This may cause equipment failure, a loss in the control of acceptable indoor environment as well as an increase of the the system maintenance cost. Hence, we propose a methodology to generate the smooth setpoint trajectories by reducing the fluctuation while maintaining the economic operational interests. This proposed methodology is motivated from the one-norm regularization analysis extensively used in the statistics and in the machine learning field. It was originally introduced as the least squares and has been extended to the various statistical models. Moreover, it is also interpreted in other fields of geometry and convex analysis. Brief details about the proposed formulation of MPC describing the cumulative objectives are presented in the next sections. This approach is illustrated on a six-zone HVAC building framework using the weather data from Nancy, France for the winter season.

Furthermore, as in buildings, fault occurrences may result in deteriorating the energy efficiency as well as the thermal comfort for the occupants inside the buildings. To address this issue, we design a fault diagnosis and fault adaptive control techniques based on the model predictive control and demonstrate the simulation results on the benchmark building. To achieve this, we present a potential approach to detect and diagnose a damper stucks in VAV boxes in the HVAC building systems. The novelty in the proposed approach is the design of the dedicated bank of unknown input residual generators for fault detection/isolation. This detection step is followed by the suitably designed fault estimators that evaluates the damper stuck value. Furthermore, we also present a novel fault-adaptive control method for the VAV damper stuck fault. The adaptive control module receives the information from the proposed fault diagnosis module. The key feature of this fault-adaptive method is to modify the constraints online in the nominal MPC controller. Though, in the general MPC-based fault-adaptive control, it is not a trivial problem to know how to change the problem formulation under the occurrence of faults, the constraint modification performed here under VAV damper stuck was made possible thanks to the designed fault estimation filters from the fault diagnosis module. Finally, the simulation results are provided for the six zone building benchmark building.

Distributed Control and Estimation :

In the case of large scale buildings these requirements to implement the centralized control may appear as a great challenge due to communication network failures, difficulty in deriving the centralized mathematical model, isolation during fault occurrences. Hence, in the recent years the theory of decentralized control and distributed control for large-scale HVAC building systems have attracted considerable attention of many researchers. The basic notion behind the decentralized and distributed control is to manage the subsystem performance instead of considering the overall system. In the decentralized control framework, the individual controllers maintain the respective subsystems performance neglecting the interaction information among the subsystems. If these subsystems interactions are significant and largely contribute in the overall building behavior then the decentralized control may result in suboptimal performance. But unlike decentralized control, in the distributed control, the communication between the subsystem controllers may ensure the optimal performance considering the interactions between the subsystems. This distributed control framework offers various advantages e.g. possibility of isolation of the subsystem in case of faulty events, modular maintenance scheme or plugging new subsystem in the existing distributed control framework etc.. Due to the obvious advantages of the distributed control framework, the Energy IN TIME project considers the importance of the scalability of the proposed energy efficient solutions given under the centralized schemes. In this chapter, we propose methodologies of distributed model predictive control (DMPC) for the benchmark HVAC building structures.

(a) Optimality Condition Based Distributed Model Predictive Control : A detailed discussion about the mathematical formulation of the proposed methods with their validation on the simulator platforms. a new approach of DMPC is proposed that realized in two steps. In first step, we decompose the system into subsystems by partitioning the Karush-Kuhn-Tucker (KKT) matrix obtained from the CMPC problem. This KKT matrix contains very significant information about the system and the cost function. We arrange this KKT matrix into block-diagonal form where blocks at the diagonal represent the subsystems and the non-diagonal elements provide the interaction information between these subsystems. This transformation of KKT matrix into block-diagonal form is performed using Cuthil-McKee algorithm. We have shown the system decomposition exercise based on KKT matrix transformation applied to a six-zone VAV based building layout. In second step, the design of subsystem controllers and their coordination strategy are established using the Optimality Conditions Decomposition introduced by Conejo. We discuss the convergence properties of the proposed distributed control scheme i.e. the equivalence between the centralized solution and distributed solution. Finally, this proposed control architecture is applied to the six-zones benchmark building. The comparison between the proposed strategy, decentralized control and centralized control is analyzed.

(b) Optimality Condition Based Distributed Moving Horizon Estimation :

The proposed distributed scheme is successfully extended to the field of state estimation based on moving horizon method. We present this method for Distributed Moving Horizon Estimation (DMHE) for linear systems. Firstly, we formulated the Centralized Moving Horizon Estimation (CMHE) problem by delineating the associated centralized cost function. Please note, this cost function is multiobjective as it minimizes the error in estimation, uncertainty in the initial state value and the output noise. Inspired from the optimality condition decomposition based distributed control, we write the the optimality conditions are written in order to formulate the KKT matrix from the given formulation of CMHE. This KKT matrix is partitioned into separable blocks where these blocks provide the description of the subsystems. The distributed scheme of estimation is implemented using the optimization procedure using Optimality Condition Decomposition. The obtained solutions are breaodcasted to the network and are used as initial values for solving next optimization problems, like the distributed control scheme. To illustrate the proposed method, simulation results are presented for the six-zones HVAC building system considering the weather data for Nancy France. In this example we demonstrate that the temperatures for the zones where sensors not available can be estimated using this proposed scheme. This notion of distributed estimation has great potential for the extension in the field of fault detection and diagnosis.

(c) Sensitivity based Distributed Model Predictive Control : Moreover, we also propose another approach of distributed control based on sensitivity analysis. This method is also addressed in two stages of DMPC for the VAV based HVAC building system as i) the decomposition of the building system into subsystems and ii) the design and coordination between obtained subsystem controllers. The proposed method of system decomposition uses the sensitivity matrix derived from the system dynamics. This sensitivity matrix contains the coupling information between the variables. Based on this information, we have partitioned the sensitivity matrix into separable blocks that represent the corresponding subsystems. Moreover, the second stage of defining the coordination between the subsystems controllers is based on the notion of sensitivity of other subsystem controllers with respect to the local controller. This is achieved by adding a linear approximation of the objectives of other subsystems to the objective of the local controller. The convergence of this proposed method is stated to demonstrate the equivalence of the distributed solution to the centralized solution. Finally, the decomposition and coordination strategies are demonstrated on the VAV type HVAC building system.

Keywords

Contrôle centralisé:

Une nouvelle méthode de contrôle prédictif basée sur la maintenance (MPC) est proposée dans ce travail. Dans ce schéma de contrôle proposé, le MPC est conçu avec de multiples objectifs à l'esprit, tels que la réduction de la consommation d'énergie, le maintien du confort thermique et la réduction des coûts de maintenance pour les équipements CVC. La première partie de l'objectif cumulatif est assez simple. Ceci est réalisé en formulant la fonction de coût dans le MPC en introduisant le terme lié à la consommation d'énergie par les équipements HVAC. La partie connexe de l'objectif cumulatif tel que la minimisation du coût d'entretien de l'équipement CVC est obtenue en assurant le bon fonctionnement des pièces d'équipement. Généralement, les chemins de consigne sont calculés par le contrôleur et sont mis en oeuvre par les actionneurs. Si ces chemins de consigne contiennent de grandes fluctuations, ils détériorent certainement les performances de l'actionneur, entraînant une usure accrue des parties de l'actionneur. Cela peut entraîner une défaillance de l'équipement, une perte de contrôle de l'environnement intérieur et une augmentation des coûts de maintenance du système. Par conséquent, nous proposons une méthodologie pour générer les trajectoires de consigne en réduisant la fluctuation tout en maintenant les intérêts économiques opérationnels. Cette méthodologie proposée est motivée par l'analyse de régularisation à une norme largement utilisée dans les statistiques et dans le domaine de l'apprentissage automatique. Il a été présenté à l'origine comme les moindres carrés et a été étendu à différents modèles statistiques. En outre, il est également interprété dans d'autres domaines de la géométrie et de l'analyse convexe. De brefs détails sur la formulation de MOC proposée décrivant les objectifs cumulatifs sont présentés dans les sections suivantes. Cette approche est illustrée sur un cadre de construction de CVC à six zones utilisant des données météorologiques de Nancy, en France, pour la saison d'hiver.

De plus, comme dans les bâtiments, les pannes peuvent entraîner une dégradation de l'efficacité énergétique et du confort thermique des occupants à l'intérieur des bâtiments. Pour résoudre ce problème, nous concevons un diagnostic de défaut et des techniques de contrôle adaptatif des défaillances basées sur le contrôle prédictif du modèle et démontrons les résultats de la simulation sur le bâtiment de référence. Pour ce faire, nous présentons une approche potentielle pour détecter et diagnostiquer un verrou d'amortisseur dans VAV boîtes dans les systèmes de construction HVAC. La nouveauté dans l'approche proposée est la conception à partir de la banque dédiée des générateurs résiduels d'entrée inconnue pour la détection / l'isolement des défauts. Cette étape de détection est suivie par des estimateurs de défauts conçus de manière appropriée la valeur de l'amortisseur de l'amortisseur. En outre, nous introduisons également une nouvelle méthode adaptative de contrôle de défaut pour la défaillance de l'amortisseur VAV. Le module de contrôle adaptatif reçoit les informations du module de diagnostic de panne proposé. La principale caractéristique de cette méthode de défaillance adaptative est de modifier les contraintes en ligne dans le contrôleur MPC nominal. Bien que, dans le contrôle général adaptatif basé sur MPC, ce n'est pas un problème trivial de savoir comment changer la formulation du problème sous l'occurrence de fautes, la modification de contrainte effectuée ici sous l'amortisseur VAV a été rendue possible. estimation du module de diagnostic de panne. Enfin, des résultats de simulation sont fournis pour le bâtiment de référence du bâtiment à six zones. (c) Sensibilité basée sur le model prédictif du mode distribué: De plus, nous proposons également une autre approche du contrôle distribué basée sur l'analyse de sensibilité. Cette méthode est également abordée en deux étapes du DMPC pour le système de bâtiment HVAC basé sur VAV: i) la décomposition du système de construction en sous-systèmes et ii) la conception et la coordination entre les contrôleurs de sous-systèmes obtenus. La méthode de décomposition du système proposée utilise la matrice de sensibilité dérivée de la dynamique du système. Cette matrice de sensibilité contient les informations de couplage entre les variables. Sur la base de ces informations, nous avons partitionné la matrice de sensibilité en blocs séparables qui représentent les sous-systèmes correspondants. De plus, la deuxième étape de la définition la coordination entre les sous-systèmes contrôleurs est basée sur la notion de sensibilité de d'autres contrôleurs de sous-système par rapport au contrôleur local. Ceci est réalisé en ajoutant une approximation linéaire des objectifs d'autres sous-systèmes par rapport à l'objectif de manette. La convergence de cette méthode proposée est indiquée pour démontrer l'équivalence de la solution distribuée à la solution centralisée. Enfin, les stratégies de décomposition et de coordination sont démontrées sur le VAV système de construction de type HVAC.

Mots-clés: Modèle de contrôle prédictif, systèmes HVAC, diagnostic de pannes et contrôle tolérant, Contrôle Distribué 108

  1. ACCIONA, Spain 2. ANA Aeroportos de Portugal, Portugal 3. CIRCE, Spain 4. Cork Institute for Technology, Ireland 5. Université de Lorraine Nancy, France 6. Centre Scientifique et Technique de Batiment, France 7. FUNIBER, NY, USA 8. Institutul de Cercetari Electrotehnice, Bucharest, Romania 9. Integrated Environmental Solutions, Glasgow ,United Kingdom 10. STAM SRL, Genoa, Italy 11. Universidad de Granada, Granada, Spain 12. United Technologies Research Center, Ireland 13. Caverion Helsinki, Finland

Figure 1 :

 1 Figure 1: EiT Project Architecture

Figure 2 :

 2 Figure 2: Faro Airport, Portugal

Figure 3 :

 3 Figure 3: SCADA screenshot of the different FCU units

Figure 4 :

 4 Figure 4: Office Building, Bucharest, Romania

  Figure 5: Commertial building, Helsinki, Finland

  (a) Building Description -Total built area: 42.500m 2 Built: 2010 -Four main heating distribution centers with different usages: hotel (170 rooms), parking, 2 apartment buildings with seasonal and high variable occupation. -The main heating source is district heating (3.630 kWh). A radiant floor installation with a heating power of 660 kW. -Operational Plans distributed in four independent spaces areas. -Building Automation System with external network based interfaces for monitoring provided.

Figure 6 :

 6 Figure 6: Hotel, Levi, Finland

Figure 1

 1 Figure 1.1: Overview of available literature dealing with the energy efficient HVAC building operation using MPC

Figure 1

 1 Figure 1.2: Centralized and distributed control scenarios (a) Control Level Decomposition (b) Mathematical Model Level Decomposition

Figure 2 . 1 :

 21 Figure 2.1: Schematic of air conditioning of a typical thermal zone

Figure 2 . 2 :

 22 Figure 2.2: Example of a Damper and AHU unit

  Figure 2.3, let T r denote the temperature of the return air flow rate at the input of the mixer. Then, assuming that there is no leakage of mass flow rate in the duct, i.e., in the return duct reads as: ṁa T r = ṁ1 T 1 + ... + ṁn T n (2.11)

Figure 3

 3 Figure 3.1: Example with (a) 1-norm regularization (b) 2-norm regularization

Figure 3

 3 Figure 3.2: Six-zones Building Layout

Figure 3

 3 Figure 3.3: Heat flux due to occupancy and outside temperature

  Figure 3.4: Zone temperatures with EMPC

Supply air flow rate for Zone 6 Figure 3

 63 Figure 3.6: Zone temperatures with 1 regularized EMPC

Figure 3

 3 Figure 3.8: Damper Stuck fault in VAV box

Figure 3 .

 3 Figure 3.10: (a)Block-diagram showing the faulty case scenario (b)Overview of the proposed UIO based fault diagnosis scheme

Figure 3

 3 Figure 3.11: Damper stuck scenario for 6 zones benchmark building

Figure 3

 3 Figure 3.12: Temperature profiles for 6 zones benchmark building

Figure 3

 3 Figure 3.13: Residual estimation for 6 zones benchmark building

Figure 3 .

 3 Figure 3.14: Structure of the fault adaptive control system

Supply air flow rate for Zone 6 Figure 3 Figure 3

 633 Figure 3.15: Zone temperatures after FAC implementation

Figure 4 . 1 :

 41 Figure 4.1: Reordering with Sparse Reverse Cuthill-McKee Algorithm

Figure 4 . 2 :

 42 Figure 4.2: Reordering and partitioning of KKT matrix for Illustrative Example

Days

Figure 4 . 3 :

 43 Figure 4.3: Zone temperature responses for 3 zones

Figure 4 . 4 :

 44 Figure 4.4: Supply Air Temperatures for 3 zones

Figure 4 . 6 :

 46 Figure 4.6: Typical School Building

Figure 4

 4 Figure 4.8: Heat Flux due to occupants and weather temperature

Figure 4 . 9 :

 49 Figure 4.9: Temperatures response for all the zones

i)

  Error associated to measurement noise Considering the past measurements (Y m ) over horizon N are available then, minimizing the error between measured and estimated output over horizon N at discrete time k,

Figure 4

 4 Figure 4.12: State estimation for all the zones

  (a) Condition d'optimisation basée sur le model prédictif du mode distribué:Une discussion détaillée sur la formulation mathématique des méthodes proposées avec leur validation sur les plateformes de simulateurs. une nouvelle approche de DMPC est proposée qui s'est réalisée en deux étapes. Dans un premier temps, nous décomposons le système en sous-systèmes en partitionnant la matrice de Karush-Kuhn-Tucker (KKT) obtenue à partir du problème CMPC. Cette matrice KKT contient des informations très importantes sur le système et la fonction de coût. Nous plaçons cette matrice KKT en forme de bloc-diagonale où les blocs de la diagonale représentent les sous-systèmes et les éléments non-diagonaux fournissent les informations d'interaction entre ces sous-systèmes. Cette transformation de la matrice KKT en forme de bloc-diagonale est réalisée en utilisant l'algorithme de Cuthil-McKee. Nous avons montré l'exercice de décomposition du système basé sur la transformation de la matrice KKT appliquée à une disposition de bâtiment basée sur VAV à six zones. En deuxième étape, la conception des contrôleurs de sous-systèmes et leur la stratégie de coordination est établie en utilisant la décomposition des conditions optimales présenté par Conejo. Nous discutons les propriétés de convergence du schéma de contrôle distribué proposé, c'est-à-dire l'équivalence entre la solution centralisée et la solution distribuée. Enfin, cette architecture de contrôle proposée est appliquée au bâtiment de référence à six zones. La comparaison entre la stratégie proposée, décentralisée le contrôle et le contrôle centralisé sont analysés. (b) Optimisation de l'espacement mobile basé sur les conditions d'optimalité: Le schéma distribué proposé est étendu avec succès au domaine de l'estimation d'état basé sur la méthode de l'horizon mobile. Nous présentons cette méthode pour l'estimation de l'horizon mobile distribué (DMHE) pour les systèmes linéaires. Premièrement, nous avons formulé l'estimation de l'horizon mobile centralisé (CMHE) en délimitant la fonction de coût centralisée associée. Veuillez noter que cette fonction de coût est multiobjective car elle minimise l'erreur d'estimation, l'incertitude dans la valeur d'état initial et le bruit de sortie. Inspiré du contrôle distribué basé sur la décomposition de conditions d'optimalité, nous écrivons les conditions d'optimalité sont écrites afin de formuler la matrice KKT à partir de la formulation donnée de CMHE. Cette matrice KKT est partitionnée en blocs séparables où ces blocs fournissent le description des sous-systèmes. Le schéma d'estimation distribué est mis en oeuvre en utilisant la procédure d'optimisation utilisant Optimality Condition Decomposition. Les solutions obtenues sont breaodcasted sur le réseau et sont utilisées comme valeurs initiales pour résoudre les problèmes d'optimisation suivants, comme le schéma de contrôle distribué. Pour illustrer la méthode proposée, les résultats de la simulation sont présentés pour le système de construction de CVAC à six zones en considérant Données cartographiques sur Nancy France Dans cet exemple, nous démontrons que les températures pour les zones où les capteurs non disponibles peuvent être estimés en utilisant ce schéma proposé. Cette notion d'estimation distribuée présente un fort potentiel d'extension dans le domaine de la détection et du diagnostic des fautes.

Table 1 :

 1 Project Details

	gives the

Type of Building Occupancy Coupling among Type of HVAC the zones systems Office

		almost fixed	negligible	VAV or FCU type
	Hotel	variable	negligible	FCU type
	Airport	highly	high coupling	VAV/CAV type
	check-in/waiting area	variable	(glass walls)

Table 2 .

 2 1: Building types and their characteristics cators as zone temperature, humidity and CO 2 concentration inside the thermal zones. However, in this work, we focus only on the zone temperature as a thermal comfort indicator. The air conditioning schematic of the typical thermal zone is shown in Figure2.1.

  26)Clearly, with the true actuator output u A,i , residual ρ i (k) generated by filter (3.26) is close to zero. Note that filters (3.25), (3.26) are completely equivalent in the fault-free case, i.e., when u A,i = u i . Under a VAV damper stuck, the stuck position is no longer known but all signals in the second equation of filter (3.26) are available since from the nullity of the residual ρ i (k), it follows that ζ i

  The function f and the constraints h have second derivatives in an open set containing z

	Algorithm 3 Algorithm for Distributed control architecture based on Optimality Con-
	dition Decomposition
	Initial Data: f , z max ,z max , ĥi , zi , zi , (i = 1, . . . , p)
	Result : z *
	1. Solve optimization problems given in (4.25) and obtain search directions ∆z i and
	∆λ

i for the given iteration 2. Update zi ← zi + ∆z i ; λi ← λi + ∆λ i (i = 1, . . . , p) 3. These updated information is broadcasted on the network 4. Stop if variables do not change significantly in two consecutive iterations, otherwise continue from Step 1

5. Implement final z * ← zi on the system

In order to state the next result, let z * denotes the optimal solution of centralized problem (4.2). For this optimal point, it should hold that:

H.1 * .

  Sensitivity based Distributed Model Predictive AlgorithmInitial Data: f , z max ,z max ,h i , zi , zi , (i = 1, . . . , p)

	Then,
	f (z) f ( z1 , z2 , . . . , zp )+ f1 (z 1 ) +	∂f ∂ z1	z (z 1 -z1 )	+ . . . . . .+ fp (z p ) +	∂f ∂ zp	z (z p -zp )	(4.64)
		f 1 (z 1 )					fp(zp)
								61)
	recall, zi (i = 1, . . . , p) are values from the previous iteration. Following algorithm sum-
	marizes the proposed method:
	Algorithm 5 Result : z *
	1. Solve optimization problems given in (4.61) and obtain search directions ∆z i and
	∆λ i for the given iteration
	2. Update zi ← zi + ∆z i and λi ← λi + ∆λ i , (i = 1, . . . , p)
	3. These updated information is broadcasted on the network
	4. Stop if variables do not change significantly in two consecutive iterations, otherwise
	continue from Step 1
	5. Implement final zi on the system
	In a case where f (z) is a non-separable cost function, to apply above algorithm, it
	may be possible to represent f (z) as follows,
	f (z) = f1 (z 1 ) + . . . + fp (z p ) + f (z 1 , z2 , . . . , zp )	(4.62)
	where fi (z i ) represents the explicit terms in f (z) that are only function of zi while
	f (z 1 , z2 , . . . , zp ) represent the coupling terms in f (z). Let us write the first order ap-
	proximation for this coupling term for the feasible values z = ( z1 , z2 , . . . , zp ) as:
	f (z 1 , z2 , . . . , zp ) = f ( z1 , z2 , . . . , zp ) +	∂f ∂ z1	z (z 1 -z1 ) + . . . +	∂f ∂ zp	z (z p -zp )	(4.63)

:

  Model Predictive Control, HVAC systems, Fault diagnosis and tolerant control, Distributed Control. sites de démonstration, nous formulons des modèles d'espace-état linéaires pour les systèmes CVC donnés. Ces modèles sont ensuite utilisés pour les exercices de contrôle et de diagnostic des pannes.Les solutions proposées dans le cadre du projet, par ex. des méthodes de contrôle éconergétiques, un diagnostic des pannes et des techniques d'adaptation aux fautes, etc. sont construits pour convenir à ces prototypes de référence. Néanmoins, une flexibilité adéquate dans les solutions proposées est préservée afin qu'elle s'applique à un large éventail de types de bâtiments. En outre, nous présentons les configurations HVAC de type VAV et FCU et leurs modèles mathématiques. Ces modèles mathématiques seront utilisés pour proposer de nouvelles approches de contrôle prédictif centralisé et distribué. Les composants énergivores de ces systèmes HVAC sont principalement le ventilateur d'alimentation et le serpentin de chauffage. Par conséquent, les objectifs de contrôle dans les prochains chapitres seront conçus autour de ces composants. Ces modèles dérivés sont suffisamment génériques pour s'adapter à différents types de bâtiments et aux conditions météorologiques.

Benchmark Building Framework: Les solutions développées dans le cadre de ce projet devraient être validées sur les sites de démonstration de différents sites européens fournis dans le cadre du projet Energy IN TIME. Nous concevons un cadre de construction de benchmark général pour émuler le comportement des sites de démonstration. Ce cadre de construction de benchmarks sert de banc d'essai pour la validation des solutions proposées dans ce travail de thèse. A partir de la synthèse détaillée des sites de démonstration, nous résumons la caractérisation des différentes topologies de construction par rapport aux attributs mentionnés. Dans les systèmes existants de ventilation et de climatisation (CVC), les configurations de type volume d'air variable (VAV) et ventilo-convecteur (FCU) sont plus populaires en raison de leur simplicité et de leurs performances écoénergétiques. Par conséquent, nous concevons les prototypes de bâtiment de référence basés sur ces configurations de CVC. Sur la base de la conception de la disposition du bâtiment de référence, nous présentons une formulation de contrôle économique en utilisant un modèle de contrôle prédictif minimisant la consommation d'énergie. De plus, nous présentons le modèle mathématique basé sur le comportement thermique du bâtiment en fonction du programme d'occupation, du type de système HVAC et du couplage entre les zones. En utilisant les données existantes des

Contrôle distribué et estimation:Enfin, une attention particulière est portée au problème de l'estimation dans le cadre de mesures limitées dans les grands bâtiments. Les techniques d'estimation avancées suggérées sont basées sur les méthodologies de l'horizon mobile et sont démontrées sur les systèmes de construction de référence. Dans le cas de bâtiments à grande échelle, ces exigences pour mettre en oeuvre le contrôle centralisé peuvent apparaître comme un grand défi en raison des défaillances du réseau de communication, de la difficulté à dériver le modèle mathématique centralisé, de l'isolement pendant les défaillances. Ainsi, au cours des dernières années, la théorie du contrôle décentralisé et du contrôle distribué pour les systèmes de construction de CVC à grande échelle a attiré l'attention de nombreux chercheurs. La notion de base du contrôle décentralisé et distribué est de gérer les performances du sous-système au lieu de considérer le système global. Dans le cadre de contrôle décentralisé, les contrôleurs individuels maintiennent les performances des sous-systèmes respectifs en négligeant les informations d'interaction entre les sous-systèmes. Si ces interactions de sous-systèmes sont significatives et contribuent largement au comportement global du bâtiment, le contrôle décentralisé peut entraîner des performances sous-optimales. Mais contrairement au contrôle décentralisé, dans le contrôle distribué, la communication entre les contrôleurs de sous-système peut assurer la performance optimale compte tenu des interactions entre les sous-systèmes. Ce cadre de contrôle distribué offre divers avantages, par ex. possibilité d'isolation du sous-système en cas d'événements fautifs, schéma de maintenance modulaire ou branchement d'un nouveau sous-système dans le cadre de contrôle distribué existant etc. En raison des avantages évidents du cadre de contrôle distribué, le projet Energy IN TIME considère l'importance de l'évolutivité des solutions d'efficacité énergétique proposées dans les schémas centralisés. Dans ce chapitre, nous proposons des méthodologies du contrôle prédictif de modèle distribué (DMPC) pour les structures de bâtiments HVAC de référence.

Resumé 1. Motivation: Depuis les deux dernières décennies, il y a eu une prise de conscience croissante sur le changement climatique et le réchauffement de la planète qui a suscité plusieurs initiatives de la Direction de diverses administrations. Ces initiatives concernent principalement la maîtrise des émissions de gaz à effet de serre, l'utilisation de ressources énergétiques non conventionnelles et l'optimisation de la consommation d'énergie dans les systèmes existants. L'Union européenne a proposé de nombreux projets dans le cadre du 7e PC pour réaliser les économies d'énergie jusqu'à 20% d'ici 2020. En particulier, conformément à la directive sur l'efficacité énergétique, les bâtiments sont principalement responsables de 40% des ressources énergétiques en Europe. 36% d'émission de CO2. Ainsi, une catégorie de projets dans le cadre du 7e PC encourage l'utilisation de technologies intelligentes dans les bâtiments et rationalise les règles existantes. Energy IN TIME est l'un des projets axés sur le développement d'une méthode de contrôle basée sur la simulation intelligente de l'énergie qui permettra de réduire la consommation d'énergie au stade opérationnel des bâtiments non résidentiels existants. Essentiellement, cette thèse propose plusieurs solutions originales pour atteindre les objectifs du projet assignés à l'Université de Lorraine.

Sensitivity Based Distributed Model Predictive Control

Introduction

In the previous section, we have proposed an approach to the distributed model predictive control methodology assuming the centralized MPC problem is available. To formulate the centralized MPC problem, overall system dynamics and systemwide cost function are the necessary requisites. In the literature, the possible presented distributed control methods requires the overall system dynamics a priori. In this section, we propose: i) a method to decompose the building system into subsystems based on the sensitivity analysis and ii) a coordination strategy between the subsystem controllers considering the sensitivities of the other subsystem controllers and the coupling information. The above distributed estimation problems can be solved using Algorithm 3. In this algorithm, values for ( X1 , . . . , Xk+1 , . . . , Xp ) and λj (j = 1, ..., p; j = k) are taken as last iterates.

Simulation Results

The proposed approach is illustrated on a Fan Coil Unit (FCU) based heating, ventilation and air conditioning (HVAC) system for a benchmark building shown in Figure 3.2. FCU maintains the thermal comfort in the given area i.e. in a thermal zone. Essentially, FCUs vary the supply air temperature entering in the zones based on the feedback from temperature sensors installed in the zones. In this case study, we demonstrate that the estimation of the temperatures are possible in a case where some of the temperature measurements are not available. For simulation purpose, we consider the temperature measurement for zone 2 and zone 5 are not available. The motivation of the DMHE implementation is to estimate the unmeasured temperature values with the provided system dynamics (2.19). The measurement data is available for the remaining temperature sensors. The formulation of CMHE for the 6 zones building can be easily written in the form of (4.66) for N = 24 hours. Applying Algorithm 2 to the problem (4.66), we evaluate KKT matrix. Partitioning the KKT matrix, we obtain the subproblems with the following two groups (p = 2) of zones as {1, 2, 5} and {3, 4, 6}. The MHE problems are implemented for both the subsystems and the simulation results for the achieved decomposition is given in next section. The type of occupancy considered corresponds to an office building where the working time is from 08:00 to 18:00. The study is carried out during the winter season at Nancy in France. The plots of heat flux due to occupancy and weather temperature over five workings days are shown in Figure 3.3.

Firstly, we simulate the open loop responses for the all states and outputs for the overall system under the defined input values. This establishes the base for the comparison of the available states estimation techniques. Figure 4.12 and 4.13 represents the estimation plots for the states and the outputs over a period of five working days. A comparative analysis with the Kalman filter state estimation is presented alongside the available actual data. Figure 4.12 shows the graphs for state estimation by the proposed method and centralized MHE. This verifies the claim discussed in the convergence section 4.1.3 that the distributed solution and centralized solutions are essentially equivalent. Finally, in Figure 4.13, the output measurements for four zones are shown with all the available state estimation techniques.

Abstract

Motivation :

Since the last two decades, there has been a growing awareness about the climate change and global warming that has instigated several Directorate initiatives from various administrations. These initiatives mainly deal with controlling greenhouse gas emissions, use of non-conventional energy resources and optimization of energy consumption in the existing systems. The European Union has proposed numerous projects under FP7 framework to achieve the energy savings up to 20% by the year 2020. Especially, stated by the Energy Efficiency Directive, buildings are majorly responsible for 40% of energy resources in Europe and 36% of CO2 emission. Hence a class of projects in the FP7 framework promotes the use of smart technology in the buildings and the streamline existing rules. Energy IN TIME is one of the projects focused on developing a Smart Energy Simulation Based Control method which will reduce the energy consumption in the operational stage of existing nonresidential buildings. Essentially, this thesis proposes several novel solutions to fulfill the project objectives assigned to the University of Lorraine.

Benchmark Building Framework :

The developed solutions under this project should be validated on the demonstration sites from various European locations provided under the Energy IN TIME project. We design a general benchmark building framework to emulate the behavior of demonstration sites. This benchmark building framework serves as a test bench for the validation of proposed solutions given in this thesis work. From the detailed synthesis of the demonstration sites, we summarize the characterization of the different building topologies with respect to the mentioned attributes. In the existing Heating Ventilation and Air conditioning (HVAC) systems, the Variable Air Volume (VAV) and Fan Coil Unit (FCU) type configurations are more popular due to their simplicity and energy efficient performances. Hence we design the benchmark building prototypes based on these HVAC configurations. Based on the design of benchmark building layout, we present an economic control formulation using model predictive control minimizing the energy consumption. Further, we present the mathematical model based on the building thermal behavior according to the occupancy schedule, the type of HVAC system and the coupling among the zones. Using existing data from Demonstrations sites, we formulate linear statespace models for the given HVAC systems. These models are then used for the control and fault diagnosis exercises.