Keywords: Centrales solaires thermodynamiques, moteurs Stirling, machines à cycle de Rankine Organique Solar thermodynamic power plants, Stirling engines, Organic Rankine Cycle (ORC) machines, control of energy conversion systems, model predictive control (MPC), sliding mode control, nonlinear model identication

Je tiens avant toute chose à remercier mon

List of Figures

Typical disposition of the vectors eld of the dynamic system (6.1c)-(6.1d) when x 2 = x st 2 (x * 3 ) and x 5 = x r 5 are assumed to be constant. . . . . . . . . . . . . . . 87 6.2 Typical disposition of the vectors eld of the dynamic system (6.1e)-(6.1f) for a given x 4 and for x 7 > x r 5 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90 6.3 Control of the (x 5 , x 6 ) system under constant x 4 . Inuence of the allowable values of x 6 on the settling time. Comparison between the cases x max 6 = 10 and x max 6 = 20. Note that the time axis is in milliseconds. The closed-loop system shows a minimum-time like behavior since it saturate the constraints on the current x 6 = I Lbb in order to accelerate the convergence. . . . . . . . . . . . . . 95 6.4 Control of the (x 3 , x 4 ) system under constant x 5 = x d 5 after a step change of +20% on the desired value x d 4 . Note that the bounds on x 3 dened by (D.2)-(D.3) slow down the response in order to meet the constraints on x 2 . The same scenario is depicted in Figure 6.8 using logarithmic scale in order to show the evolution of the variables during rst milliseconds of the scenario. . . . . . . . . 96 6.5 Control of the (x 3 , x 4 ) system under constant x 5 = x d 5 after a step change of -20% on the desired value x d 4 . Note that the bounds on x 3 dened by (D.2)-(D.3) slow down the response in order to meet the constraints on x 2 . The same scenario is depicted in Figure 6.9 using logarithmic scale in order to show the evolution of the variables during rst milliseconds of the scenario. . . . . . . . . 97 6.6 Control of the (x 3 , x 4 ) system under constant x 5 = x d 5 after a step change of -20% on the desired value x d 4 . Inuence of the constraints on x 2 on the response time of the system. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98 6.7 Control of the (x 3 , x 4 ) system under constant x 5 = x d 5 after a step change of +20% on the desired value x d 4 . Inuence of the constraints on x 2 on the response time of the system. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 6.8 Same scenario as the one depicted on Figure 6.4 using logarithmic scale on the time axis in order to show the beginning of the scenario. This scale enables the saturation on x 4 to be clearly shown. . . . . . . . . . . . . . . . . . . . . . . . . 99 6.9 Same scenario as the one depicted on Figure 6.5 using logarithmic scale on the time axis in order to show the beginning of the scenario. This scale enables the saturation on x 4 to be clearly shown. . . . . . . . . . . . . . . . . . . . . . . . . 100 6.10 Performance of the overall closed-loop system under varying power load demand. Case where x max 2 = 5. This scenario is worth comparing with the scenario depicted in Figure 6.11 where the upper bound x max 2 = 4.8 is used. . . 101 6.11 Performance of the overall closed-loop system under varying power load demand. Case where x max 2 = 4.8. This scenario is worth comparing with the scenario depicted in Figure 6.10 where the upper bound x max Today, more than 1.4 billion people, that corresponds to 300 millions households, have no access to electrical energy. They earn less than 2 dollars per day and spend more than 15 dollars per month and per family for energy access 1 :

A ⊗ B :=      A(1, 1) • B A(1, 2) • B • • • A(1, m) • B A(2, 1) • B A(2, 2) • B • • • A(2, m
• 500 millions live in India.

• 100 millions live in Indonesia

• 100 millions live in Bangladesh

• 70 millions live in Nigeria.

• more than 220 millions live in Sub-Saharian Africa.

These populations generally live in areas where the electrical grid is not available, generally due to its high cost, or is unstable. However, the need for energy of these populations is increasingly growing since energy is the basis of a sustainable development2 :

• Economical development via local entrepreneurship.

• Mechanization of agricultural production and rural depopulation limitation.

• Access to healthcare. • Reduction of health risks associated with kerosene and biomass use.

The fuel based generators (mainly Diesel engines) are typically used in this context but are not well adapted because of their low reactivity to changes in the power demand level, their high pollution, noise, and their low eciency [START_REF] Kalantar | Dynamic behaviour of stand-alone hybrid power generation system of wind turbine, microturbine, solar array and battery storage[END_REF]. Moreover, the cost for fuel routing to the concerned areas is extremely expensive compared to the daily mean salary of the consumers.

In order to bring solutions to the energy access for those populations, Schneider Electric launched a program called BipBop 3 (Business, Innovation and People at the Base of the Pyramid) that aims at bringing safe, clean and aordable electricity to the "Base Of the Pyramid" (largest, but poorest socio-economic group) focusing on the following three axes 4 depicted in Figure 1.1:

• Development of a dedicated oer for new energy distribution solutions.

• Providing nancial resources to promote the creation of innovative companies in the eld of energy.

• Training young people from the "Base Of the Pyramid" in electrical elds.

Context of the thesis 3

Within the framework of the BipBop program, the "MICROSOL" project was launched within the "Advanced Mechatronics and Control" team (Innovation center) of Schneider Electric. This project aims at developing O-grid solar thermodynamic micro power plants producing a minimum of 150kW h elec /day and functioning 24/24h thanks to an adequate environmentally friendly heat storage. The innovative solution proposed by MICROSOL is competitive and well adapted to the populations of the "Base Of the Pyramid" that do not have access to the electrical grid but located in areas with high sunlight potential.

Preliminary studies prior to MICROSOL have permitted to determine an adequate architecture for the solar thermodynamic micro power plant with a pressurized water tank, harmless for the environment, as a heat storage device. The choice of water as heat transfer and storage uid is due to its low environmental impact in addition of being low cost. However, this choice limits the hot operational temperature and therefore the global eciency of the plant. The innovations proposed by Schneider Electric and its partners are then crucial to overcome these limitations. In particular, ecient control strategies are essential to bring the plant at its maximum eciency operating point regardless of the operating conditions (ambient temperature, weather conditions...).

In the beginning of the MICROSOL project, and taking into account the technological and industrial risks, two variants of the solar micro power plant were considered by Schneider Electric:

• A rst variant uses solar collectors with xed vacuum tubes, provided by SAED, heating water between 110 • C and 150 • C pressurized at a maximum pressure of 10 bars. A Stirling engine provided by STIRAL (http://www.stiral.com) is used in this conguration to transform heat into mechanical work driving an electrical generator while an energy conversion system provided by Schneider Electric adapts the produced power to the O-grid loads.

• The second variant uses parabolic trough collectors with a system of sun tracking provided by EXOSUN (http://www.exosun.fr/). The collectors heat a pressurized water (at a maximum pressure of 15 bars) between 140 • C and 180 • C. An Organic Rankine Cycle (ORC) system converts the water heat into mechanical work thanks to a Scroll turbine (see [Quoilin, 2011]) driving an electrical generator. An adequate power conversion system provided by Schneider Electric is interfaced between the turbine and the O-grid loads.

Concerning the rst variant, in the beginning of the project, the Stirling engine was in the design phase within STIRAL laboratories. It has been then decided to buy a Stirling engine prototype from the company "Cool Energy" (http://coolenergy.com/), that was not a project partner, in order to start the preliminary tests on the energy conversion system associated to this engine and developed by Schneider Electric, this to be ready when the Stirling engine prototype provided by STIRAL will be available. Because of some technical diculties encountered by STIRAL and nancial issues encountered by SAED, this rst variant of the solar micro power plant was abandoned by Schneider Electric and the second variant was selected as a candidate for the future solar micro power plants commercialized by Schneider Electric and its project partners. However, the tests on the Stirling engine prototype, provided by "Cool Energy", with its associated energy conversion system, developed by Schneider Electric, continued in collaboration with the CEA5 .

Aims of the thesis

This CIFRE 6 Ph.D. thesis was prepared within the Gipsa-lab laboratory in partnership with Schneider-Electric in the scope of the MICROSOL project. The general problem under consideration in this thesis is the development of advanced control strategies for the energy management of the two kinds of solar thermodynamic micro power plants considered in the MICROSOL project namely the Stirling based and the ORC based solar micro power systems. This thesis is then divided in two parts, each one focusing on a particular kind of solar micro power plant as follows:

Part I. This part is concerned with the control of an O-grid Stirling based solar micro power plant focusing on the control of the energy conversion system associated to the Stirling engine hybridized to a supercapacitor that serves as an energy buer during load power demand transients. For the considered system, the main issue to solve is the size reduction of the energy buer consisting in the supercapacitor. This is done by controlling the Stirling engine through its electrical part such that to move from an operating point to another (during load power demand changes) as fast as possible while respecting the system constraints. In the literature, one can nd many works that deal with the control of a similar power conversion system topology as the one used in the plant under consideration. However, the resulting controllers found in the literature do not induce optimal performances in term of settling time and do not handle explicitly the physical system constraints. Moreover, one of their major disadvantages is the fact that they are built upon cascaded controllers that necessitate a lot of parameters that could be tricky to tune especially when the system constraints (like actuator saturation) are not explicitly considered.

In this part, two dierent control architectures that tackle the aforementioned issues are proposed. After deriving the nonlinear model (validated experimentally) of the Stirling engine under consideration in Chapter 3 and the nonlinear dynamic model of its associated energy conversion system in Chapter 4, a rst novel control strategy for the Stirling based power plant is proposed in Chapter 5. The closed loop system under the rst proposed control strategy is endowed with near-to-optimal performances in term of settling time in addition of necessitating few design parameters. The control strategy was validated experimentally, rst on a dedicated electrical motor testbed and nally on a real Stirling based power plant. A second control strategy, which is actually a control framework that holds for a family of power circuits, is proposed in Chapter 6. The novelty of the proposed solution is the explicit handling of the physical constraints of the energy conversion system while ensuring near-to-1.2. Aims of the thesis 5 optimal performances of this last one.

Part II. This part focuses on the control of an ORC based solar micro power plant hybridized with a battery bank as energy buer to accommodate for the load power demand fast transients. Here again, the main issue to solve is the size reduction of the battery bank by optimally controlling the ORC system in order to move from an operating point to another as fast as possible while respecting the ORC system constraints. Contrary to the Stirling based power plant where the control is performed acting on its associated electrical part, in the ORC system under consideration, the control is performed through the thermo-hydraulic part acting on the pumps and fan associated to the ORC system in order to dynamically control the produced torque of the turbine to meet some desired electrical power. In this manuscript, we propose an innovative solution that responds to this need by rst deriving in Chapter 7 nonlinear dynamic models of the ORC power plant using experimental data and then synthesizing in Chapter 7 the model based real time constrained optimal controller that regulates in an optimal manner the produced electrical power of the ORC system while respecting the physical constraints.

The works presented in this thesis have been the subject of the following publications:

• Journals:

1. M. Alamir, M.A. Rahmani, and D. Gualino. Constrained control framework for a stand-alone hybrid (stirling engine)/supercapacitor power generation system. Applied Energy, Vol. 118, pp. 192-206, 2014. 2. M. Alamir, M.A. Rahmani, and D. Gualino. Generic model identication structure for thermodynamic engines for use in hybrid power stations control and simulation.

To appear in Journal of Process Control, 2014. Two patents have also been submitted within the framework of the this thesis:

1. M. A. Rahmani, D. Gualino and M. Alamir. Method for producing electricity via an installation comprising a motor and an energy stocking device allowing the latency of motor load changes to be palliated, FR1353642 (patent pending).

2. M. A. Rahmani, D. Gualino and M. Alamir. Optimal Real time control of an o-grid Organic Rankine Cycle (ORC) power plant, patent pending approval from Schneider Electric.

Introduction

Solar energy constitutes the greatest energetic potential available on earth and the idea of using this energy dates back to antiquity, however, the rst experimental solar thermal power plant for practical use (pumping and irrigation) was built in Méadi (Egypt) only in 1913 by Franck Schumann [Mathieu, 2012]. Afterwards, the researches aiming at exploiting this energy experienced ups and downs especially because of the emergence of ecient and low cost fossil fuels. Indeed, the evolution of oil prices have often conditioned the evolution of the researches on solar thermal energy. Hence, because of the attractive cost of oil after the second world war, the researches on solar energy have almost stagnated until 1974 where a sharp rise on oil prices have revived the researches on renewable energies especially solar energy where many experimental solar facilities emerged in many industrialized countries. However, after the fall of oil's price in 1986 and its stabilization in the nineties, the researches on solar energy have also fallen.

Since the last fteen years, global warming has aroused the attention and concern of the international community and public opinion. This, and because of the dramatic decrease of the accessible fossil energy stock, have revived the researches on alternative energies with a particular attention on solar energy since it is the greatest energy resource technically accessible on earth, especially in desert regions around equator, that can respond to the worldwide energetic demand in a sustainable manner. Indeed, it is estimated that if 1% of the surface in 8 Chapter 2. Solar thermodynamic micro power plants arid regions rich in sunlight (like in Sahara) is covered by solar power plants, the worldwide demand in electricity can be fullled [Mathieu, 2012] as it is illustrated in Figure 2.1.

The standard technological solution based on solar energy only consists in using photovoltaic panels associated to a battery bank as an energy storage. This solution leads to a price for a kW h elec within 20 years lower than that of a solution consisting in using fuel based generators even if its initial investment cost is much higher. However, this solution suers from many drawbacks:

• This solution needs a large amount of batteries to supply the loads in O-grid areas during the night and the non sunny days.

• The price for such battery bank is high.

• The batteries lifetime is limited, especially in hot conditions, in addition of being signicantly aected by the irregular charge/discharge due to the load power demand changes.

• The most commonly used batteries are lead-acid ones that are harmful for the environment and for the populations especially in isolated areas lacking of recycling networks.

Currently, the most commonly adopted solution that takes advantage of solar energy while necessitating less batteries consists in hybridizing fuel based generators with photovoltaic panels (associated to a smaller battery bank compared to the previous solution). This solution is currently adopted within the framework of the BipBop program since it is a good trade-o between the fuel saving and the number of electrochemical batteries.

The objective of the BipBop program with the MICROSOL project is to propose a "greener" solution for the electricity production exploiting the thermal energy of the sun through solar thermodynamic micro power plants with a thermal storage as a principal storage device instead of batteries that are used in small amounts during load power demand transients. The general architecture of a such typical solar thermodynamic micro power plants is described in the next section.

General description of a solar thermodynamic micro power plant

Figure 2.2 illustrates the general architecture of a typical solar thermodynamic micro power plant. It consists essentially of a solar collector eld that collects the sunlight to heat some heat transfer uid owing inside solar absorbers, using a volumetric pump, between the solar collectors and a heat storage tank (containing the heat transfer uid). The power plant contains also a cooling loop consisting of a circulation pump and a cooling fan that cools down an appropriate coolant (water-glycol mixture for the thermodynamic power plant under study). The heat transfer uid and the coolant enter a thermodynamic machine that performs a thermodynamic cycle during which it absorbs heat from the heat transfer uid through a heat exchanger and rejects heat to the coolant through another heat exchanger resulting in a positive motor torque that drives an electrical generator producing electricity. The electrical power at the output of the generator is then conditioned through an appropriate energy conversion system to supply the electrical loads (and the system auxiliaries). In O-grid applications, the solar thermodynamic micro power plants are hybridized to an energy buer (batteries and/or supercapacitor) linked to the energy conversion system associated to the thermodynamic engine. This energy buer source responds in a fast manner to the load power demand transients while the thermodynamic machine adapts (through an appropriate control strategy) its produced power in a rather slower manner due to its mechanical and thermal inertia.

In the next sections some of the important components of a solar thermodynamic micro power plant, that greatly inuence its price and eciency, are briey described.

The solar collector eld and solar absorbers

There are many technological solutions for transforming sunlight into useful heat. One can nd a good literature review of these techniques in [Mathieu, 2012]. In a general manner, the transfer uid to be heated ows inside solar absorbers (receivers) whose casing has to be designed such that to maximize the absorption of radiation while minimizing the thermal losses. Many manufacturers provide such absorbers like the vacuum tube based SCHOTT PTR 701 solar absorber depicted in Figure 2.3. The absorbers used in the ORC based solar The solar receivers can absorb the solar radiation using:

• a direct reception: when the absorbers are exposed directly to the solar radiation using simple colectors or advanced ones (see [Mathieu, 2012]).

• a reception by concentration: occurs when mirrors with paraboloid proles are used to concentrate the sunlight into the mirror's focal.

For the ORC based power system under consideration, parabolic trough mirrors, endowed with one directional solar tracking system, are used to concentrate the solar radiations into the absorbers located in the focal of the mirrors as depicted in Figure 2.4.

Chapter 2. Solar thermodynamic micro power plants

The thermodynamic machines

Many thermodynamic cycles can be used to convert the heat of the transfer uid, stored in the heat storage tank, into mechanical work. In this thesis, two thermodynamic cycles have been considered:

• A Stirling cycle engine presented in Chapter 3. A prototype of this engine was provided by "Cool Energy" and depicted in Figure 2.5.

• Organic Rankine Cycle (ORC) presented in Chapter 7. A prototype of this engine was provided by EXOES and depicted in Figure 2.6.

The heat transfer uid

The choice of an adequate heat transfer uid greatly inuences the thermal eciency of the solar collector eld and the thermodynamic engines.

Among the most used heat transfer uids for solar thermal applications [Mathieu, 2012]: water, glycogen polymers, mineral oils, silicone based oils and synthetic oils. For the solar thermodynamic micro power system under study, the selected heat transfer uid is water because it respects the following two conditions that are crucial when the rural electrication market is targeted:

• it has a low environmental impact and it is not dangerous (toxicity, inammability...).

• it is low cost and easily transportable and can be safely manipulated by local sta. Although liquid water is an excellent heat transfer uid, its drawback is its low evaporating temperature at atmospheric conditions. In order to keep it in the liquid thermodynamic state, its pressure should be increased which increases also the price of the hydraulic components especially the pipes and the heat storage tank. For the considered application, the water pressure was limited to 10 bars at 180 • C.

This limitation of the hot water temperature induces technical challenges on the design of the thermodynamic engines (especially their heat exchangers) and on the control strategies used to dynamically optimize the thermoelectric eciency of the power plants.

Conclusion

In this chapter a brief description of solar thermodynamic micro power plants has been presented focusing on their advantages. The main components of the solar thermodynamic micro power plants have also been discussed together with the technological choices of the MI-CROSOL project. In the following chapters, the two kinds of solar thermodynamic micro power plants considered within the framework of this project will be addressed in details focusing on the development of innovative and ecient control strategies for the electrical energy production in O-grid mode.

Introduction

In this part, the control of the o-grid hybrid Stirling engine/Supercapacitor power plant under study is addressed and solved. The proposed controllers drive the system through its electrical part to satisfy the varying load power demand. Two novel controllers have been proposed to reach the control objectives. A rst controller requiring a few tuning parameters is proposed in Chapter 5 based on a backstepping controller and a one step predictive controller. The closed loop system under this rst proposed control is endowed with quasi optimal behavior. The control methodology was then experimentally validated rst on a dedicated electrical motor testbed and then on the real Stirling based power system.

The second proposed control methodology presented in Chapter 6 tackles explicitly the physical constraints of the energy conversion system (associated to the Stirling engine) while ensuring quasi optimal behavior of this last one. This proposed methodology is more than a simple controller, it is a control framework that holds for a family of energy conversion systems.

A quick state of the art on controlling similar hybrid systems Many hybrid power generation systems have been studied in the literature either for grid or o-grid applications. They combine dierent renewable energy sources essentially wind, solar photovoltaic, hydraulic and fuel cells. Dierent power electronics architectures associated to those hybrid power systems have been proposed and their corresponding control strategies developed. In [START_REF] Haque | A novel control strategy for a variable speed wind turbine with a permanent magnet synchronous generator[END_REF], a variable speed wind turbine associated to a simple power electronics have been used to supply AC loads, an MPPT (Maximum Power Point Tracking) algorithm was used to extract maximum power from the wind using a DC/DC boost converter while an inverter was used to create the AC grid and to regulate the DC bus, however because of the random nature of the wind, the AC voltage was perturbed. In [START_REF] Mendis | A novel control strategy for stand-alone operation of a wind dominated raps system[END_REF], a similar system combined to a Diesel Engine and a battery bank was studied. The batteries act as a buer source (through a buck-boost converter) to regulate the DC bus absorbing power in case of high wind and discharging in case of scarce wind. Linearized models were used to control the system at some operating points. However for many applications where the loads and power sources can vary signicantly, such controllers are not suitable. In [START_REF] Haruni | Dynamic operation and control of a hybrid wind-diesel stand alone power systems[END_REF]], a similar system provided by dump loads and its associated switch was studied and PI based controllers were proposed. In [START_REF] Ahmed | Development of an ecient utility interactive combined wind/photovoltaic/fuel cell power system with mppt and dc bus voltage regulation[END_REF] a hybrid wind/photovoltaic/fuel cell power system for grid-connected applications was proposed. An MPPT techniques were used to extract maximum power from the wind and PV panels while the fuel cells system suppress the grid power uctuations. The objectives of such hybrid systems is to act as a battery-less uninteruptible power source (UPS). [START_REF] Valenciaga | Passivity/sliding mode control of a stand-alone hybrid generation system[END_REF], [START_REF] Valenciaga | Variable structure system control design method based on a dierential geometric approach: application to a wind energy conversion subsystem[END_REF], [START_REF] Valenciaga | High-order sliding control for a wind energy conversion system based on a permanent magnet synchronous generator[END_REF]) proposed an interesting technique to control a hybrid wind/photovoltaic (and even extra power sources) with wind as a primary source and a battery bank on the DC bus as an energy storage device. The control design proposed by Valenciaga et al are based on passivity/sliding mode techniques (and a second order sliding mode technique in [START_REF] Valenciaga | High-order sliding control for a wind energy conversion system based on a permanent magnet synchronous generator[END_REF]) for dierent wind regimes.

Supervisory controllers for such systems were proposed in [START_REF] Valenciaga | Supervisor control for a stand-alone hybrid generation system using wind and photovoltaic energy[END_REF], [Qi et al., 2011a] and [START_REF] Qi | A two-time-scale framework to supervisory predictive control of an integrated wind/solar energy generation and water desalination system[END_REF]. All the techniques previously cited don't handle explicitly the constraints of the system's state variables and controllers which can lead to some undesirable behaviors. Moreover, they do not necessarily induce optimal performance of the system in term of settling time which can lead to over-dimensioning of the energy buer source especially when it consists of supercapacitors which is the case for our system because of their many advantages (see [Langlois, 2006] for more details). The aim of this chapter is to study the Stirling thermodynamic engine associated to the Stirling based o-grid power system under consideration. We rst give some generalities about Stirling engines focusing on the alpha type motor since we have in our disposal a 3 kW alpha type low temperature Stirling motor (from "cool energy" company ). Then, a static model of this device taking into account dierent heat losses will be developed and validated using real data gathered on the prototype. The static model will serve to simulate the behavior of the Stirling engine when connected to an adequate power electronics developed within Schneider Electric (a static model is sucient since the dynamics of the thermodynamic part is much slower than that of the electrical part) and on which we designed a real time control strategy presented in Chapters 5 and 6. The static model of the Stirling engine will also serve to construct a cartography that, according to the operating conditions and load power demand, computes the hot and cold mass ow rates (of the external uids owing in the hot and cold heat exchangers) references that optimize the net thermoelectric eciency of the plant.

Generalities about Stirling engines

The rst prototype of Stirling engines was created by Robert Stirling who patented it in 1816 before his brother James industrialized it in 1843. This engine belongs to the family of external combustion engines in the sense that the needed heat is generated externally to the engine. The Stirling engine has many advantages, compared to other kinds of thermodynamic engines, among them [Bonnet, 2005]:

• The Stirling engine accepts any kind of heat sources particulary: the combustion of any kind of fuels (solid, liquid, gaseous, biomass...), solar energy, nuclear source and geothermal sources.

• Stirling engines have comparable performances to Diesel engines having the same nominal power.

• If a combustion of a fuel is used to generate heat, the process is continuous and external to the engine which implies:

the combustion is perfectly controlled and the rejected pollutants can be considerably reduced.

no explosion happens inside the engine.

the produced motor torque is regular and the engine mechanical parts are less solicited which enhance the engine lifetime.

extremely soundless engine operation.

One of the key elements of Stirling engines is the "regenerator" matrix (metal mesh or foam for example) which is an internal heat exchanger playing the role of temporary heat buer placed between the hot and cold spaces of the engine such that the working uid passes through it rst in one direction, taking heat from the uid, then the other, releasing heat to the uid. The regenerator benets from high surface area, high heat capacity, low conductivity and low ow friction. Its function is to retain within the system that heat which would otherwise be exchanged with the environment at temperatures intermediate to the maximum and minimum cycle temperatures, thus enabling the thermal eciency of the cycle to approach the limiting Carnot eciency (http://en.wikipedia.org/wiki/Stirling_engine) leading to higher power output for a given absorbed heat power.

Stirling engines congurations

There are three main kinds of Stirling engines, depending on their geometrical congurations, these are: the alpha, beta and gamma types. The alpha type engines are bi-cylinder motors (see Figure 3.1) where each cylinder contains a piston linked mechanically to a rotating shaft and where the two pistons are out of phase by π 2 in such a way that their kinematic coupling enables the internal working uid to circulate from one cylinder to another through a "regenerator" (which is one of the major element in the engine that greatly inuences the engine eciency) when the shaft is rotating. One cylinder is linked through an adequate heat exchanger to a heat source whereas the other is linked by another heat exchanger to a cold sink to evacuate heat.

Beta type Stirling engine

The beta type engine uses a displacer-piston arrangement (see Figure 3.2) where both the displacer and the piston are situated in an in-line cylinder system. The displacer enables the working uid to be in contact with the hot source and the cold sink (through heat exchangers) when the piston, which is linked to the rotating shaft, is moving. The gamma type engine uses also a displacer-piston arrangement but both the displacer and the piston are situated in separate cylinders (see Figure 3.3).

3.1.2 Stirling engine theoretical Pressure-Volume (P-V) diagram Starting from point 1 (according to Figure 3.4), the working uid undergoes an isothermal expansion at temperature T h until point 2, the two pistons move down and the total volume increases. In order to ensure this isothermal behavior, the working uid absorbs heat Q exp from the external heat source. During this step, the engine produces a positive mechanical work W exp . Between points 2 and 3, the working uid undergoes an isochoric cooling where the uid is discharged from the hot cylinder to the cold one at constant volume (isochoric behavior). During this step, the uid releases a certain amount of heat Q reg c to the regenerator.

Between points 3 and 4, the working uid undergoes an isothermal compression at temperature T c where the two pistons move up simultaneously and the total volume decreases at constant temperature. During this step, the uid consumes a mechanical work W comp and releases heat Q comp to the cold sink to ensure the isothermality of this stage.

The last step occurs between points 4 and 1, to close the cycle, where the working uid undergoes an isochoric heating. During this last step, the uid is discharged back to the hot cylinder at constant volume and absorbs a certain quantity of heat Q reg h from the regenerator which theoretically corresponds exactly to Q reg c .

The net produced mechanical work corresponds to the area delimited by the closed curve 1 -2 -3 -4.

Model with losses of an alpha type Stirling engine

To derive a basic model for the alpha type Stirling engine, we will base our study on the isothermal analysis presented in [START_REF] Urieli | Stirling Cycle Engine Analysis[END_REF] which was perfectly summarized in http://www.ohio.edu/mechanical/stirling/. Based on this analysis, a modied static model for the Stirling engine that takes into account dierent losses will be derived.

Ideal Isothermal Analysis

Consider the simplied diagram of an alpha type Stirling engine depicted on the top of Figure 3.5. As it can be seen from this gure, the Stirling engine consists of a varying volume compression space followed by a xed volume cold side heat exchanger, a regenerator and a hot side heat exchanger linked to a varying volume expansion space. In the ideal isothermal analysis, it is a assumed that the working gas in the compression space and the cold side heat exchanger space is maintained at the same lower bound temperature T c while the gas in the expansion space and in the hot side heat exchanger space is at the same upper bound temperature T h . The temperature in the regenerator is assumed to be varying linearly between T c and T h (see the bottom of Figure 3.5).

The total mass of gas M g inside the whole working volume of the Stirling engine is the sum of gas masses inside each working space as follows:

M g = m comp + m c + m r + m h + m exp (3.1)
where m comp , m c , m r , m h , m exp denote respectively the mass in the compression space, the cooler space, the regenerator space, the heater space and the expansion space. Using the fundamental equation of perfect gases:

M = P.V R.T (3.2)
where M , P , V , R, T denote respectively the mass, pressure, volume, perfect gas constant and temperature, one obtains together with (3.1) the following equation for the total mass of the working gas M g :

M g = P. Vcomp Tc + Vc Tc + Vr Tr + V h T h + Vexp T h R (3.3)
where V comp , V c , V r , V h , V exp denote respectively the compression space, the cooler space, the regenerator space, the heater space and the expansion space volumes. In (http://www.ohio. edu/mechanical/stirling/), a mean eective regenerator temperature was dened by:

T r = T h -T c ln( T h Tc ) (3.4)
This temperature T r takes into account the fact that the the working gas mass distribution along the regenerator length is not uniform but depends on the temperature prole that is assumed to be linear as depicted in Figure 3.5.

Combining (3.3) and (3.4), one obtains:

P (θ, T h , T c ) = M g .R T h -Tc + V h T h + Vexp T h (3.5)
where θ is the shaft angle. Note that the gas pressure given by (3.5) depends on θ since the volumes V comp and V exp depend also on this variable as it will be explained in Section 3.2.2.

The net work produced by the Stirling engine over a complete cycle is then given by the cyclic integration of the product P.dV as follow:

W tot = W exp + W comp = (P (θ, T h , T c ).dV exp ) + (P (θ, T h , T c ).dV comp ) = P (θ, T h , T c ).( dV exp dθ + dV comp dθ ) dθ (3.6)
By performing an energy analysis over a complete cycle for each cell of Figure 3.5 (see (http://www.ohio.edu/mechanical/stirling/) for more details) one obtains the following results:

Q exp = W exp (3.7a) Q comp = W comp (3.7b) Q c = 0 (3.7c) Q r = 0 (3.7d) Q h = 0 (3.7e)
Equations (3.7a)-(3.7b) give the relation between the absorbed (resp. the rejected) heat and the work produced (resp. consumed) which will serve further to compute the produced electrical power. According to (3.7c) and (3.7e), both the hot and cold heat exchangers seem to be redundant since the cyclic heat rejected from the cooler space Q c and transferred to the heater space Q h are set to zero. This is a direct result of the isothermal assumption for the compression and expansion spaces that implies maintaining these chambers at the respective cooler and heater temperatures.

The main issue in the presented model consists in computing the cyclic integral (3.6) which depends on the engine conguration. This is presented in the following section based on the Schmidt analysis.

Schmidt Analysis

In order to compute (3.6), Gustav Schmidt proposed in 1871 an analysis of the cycle assuming a sinusoidal volume variation for the compression and expansion chambers as illustrated in Figure 3.6. This leads to the following equations for the compression and expansion spaces volumes:

V comp = V clc + V swc .(1 + cos θ)/2
(3.8a)

V exp = V cle + V swe .(1 + cos(θ + α))/2 (3.8b)
where V clc , V swc refer respectively to the clearance and swept volumes for the compression chambers (see Figure 3.6) while V cle , V swe refer respectively to the clearance and swept volumes for the expansion chambers.

The analysis proposed by Schmidt yields the closed form solution for the cycle given by (3.9) and (3.10) where the superscript sch refers to "Schmidt" (refer to (http://www.ohio. edu/mechanical/stirling/) for the detailed proof).

W sch comp = Q sch comp = π.V swc .P mean . sin β. ( √ 1 -b 2 -1) b (3.9a) W sch exp = Q sch exp = π.V swh .P mean . sin(β -α). ( √ 1 -b 2 -1) b (3.9b) W sch tot = Q sch tot = W sch comp + W sch exp (3.9c)
where

β = arctan( V swh . sin α T h V swh . cos α T h + Vswc Tc ) (3.10a) c = 1 2 . ( V swh T h ) 2 + 2. V swc T h . V swc T c + ( V swc T c ) 2 (3.10b) s = V swc 2.T c + V clc T c + V c T c + V r .ln( T h Tc ) T h -T c + V h T h + V swh 2.T h + V clh T h (3.10c) b = c s (3.10d) P mean = M g .R s. √ 1 -b 2 (3.10e)
where P mean is the mean working gas pressure inside the Stirling engine cycles.

To complete the previous ideal isothermal analysis, one needs to take into account the heat losses that occur inside the engine in order to be as close as possible from the reality. This is performed in the next section.

Heat losses analysis

The dierent heat losses that occur inside the engine are essentially due to the regenerator imperfection, the heat leaks, the friction losses and the joule eect inside the synchronous generator connected to the Stirling engine. In this section the mathematical model of each of the aforementioned heat losses is derived.

Regenerator losses

The regenerator is a cyclic device where during the rst part of the cycle the hot gas ows from the heater to the cooler through the regenerator releasing some quantity of heat Q reg c to the regenerator matrix. During the second part of the cycle, the cold gas comes back to the cooler through the regenerator and absorbs a quantity of heat Q reg h from the regenerator matrix. In an ideal regenerator |Q reg c | = |Q reg h | that leads to a zero cyclic heat transferred externally to the regenerator: Q r = 0. However, since the regeneration process is not perfect, when the working gas ows from the cooler to the heater (during the rst part of the cycle) through the regenerator, it will have at the regenerator output a temperature lower than that of the heater. This induces an extra heat being externally supplied to increase the working gas temperature and bring it to the heater temperature T h to comply with the isothermal analysis requirement. Analytically, this implies:

Q exp = Q sch exp + Q reg h .(1 -ξ reg ) (3.11)
where ξ reg is the regenerator eciency (that takes a value of 1 for ideal regenerator and 0 for a non regenerative process).

Similarly, during the second part of the cycle, when the working gas comes back from the heater to the cooler through the regenerator, an extra heat has to be rejected from the cooler to decrease the temperature at the output of the regenerator to set it at T c . This implies:

Q comp = Q sch comp -Q reg c .(1 -ξ reg ) (3.12)
The ideal heat power transfered to the regenerator during each part of the cycle depends on the temperature dierence between the heater and the cooler according to:

Qreg h = Qreg c = N rot 60 .M g .c gas p .(T h -T c ) (3.13)
where N rot is the shaft rotational speed in RPM (Rotation Per Minute) and c gas p is the working gas thermal capacity.

In (http://www.ohio.edu/mechanical/stirling/), it has been suggested that the regenerator eciency can be expressed in term of the Number of Transfer Units (NTU) as follow:

ξ reg = 1 1 + 1 N T U (3.14)
This N T U variable, in addition of being dependent on the regenerator dimensions, depends on the working gas mass ow rate inside the regenerator which is proportional to N rot . Thus, we suggested to set the following form for the regenerator eciency:

ξ reg = 1 1 + N T U 0 .( Nrot N 0 ) αreg (3.15)
where N T U 0 , N 0 and α reg are some parameters to be identied using experimental data.

Heat leaks

Since the Stirling engine enclosure is not totally adiabatic some heat leak exists from the hot source to the cold sink whose power Qleak is set according to the following equation:

Qleak = K leak .(T h -T c ) (3.16)
where K leak is an unknown parameter to be identied. This leak power has to be added in the expression of the absorbed power Qexp and substracted from the expression of the rejected power Qcomp .

Friction losses

It is well known that the power of the mechanical friction losses are proportional to the rotational speed and its square value according to:

Qmech = (A + B.N rot ).N rot (3.17)
where A and B are some parameters to be identied. Qmech has also to be respectively added and substracted from Qexp and Qcomp .

Heat loss by Joule eect

Since the Synchronous generator windings have a resistive eect (through the line winding resistance R s ), an extra heat has to be rejected from the cooler that is given by:

Qres = 3.R s .I 2 mot (3.18)
where I mot is the line current of the electrical generator connected to the Stirling engine.

Heat transfer through the heat exchangers

The Stirling engine is connected to the heat source and cold sink through heat exchangers to respectively absorb heat from the hot uid (heated by mean of the solar collector eld and stored inside the heat tank) and reject heat to the cooling loop. The hot uid enters the Stirling engine at temperature T h in with a mass ow rate ṁh and exits the engine at temperature T h out transferring a certain quantity of heat Q exp to the working gas (inside the Stirling cycles) whose heat power Qexp is given by:

Qexp = ṁh .c h p .(T h in -T h out ) (3.19)
where c h p is the hot uid thermal capacity.

Suppose a mean temperature for the hot uid inside the heat exchanger (hot uid side):

T H = T h in + T h out 2
The heat transfer between the hot uid and the working gas that occurs across the heat exchanger walls can be written in term of heat power as follows:

Qexp = K h (N rot , ṁh ).( T h in + T h out 2 -T h ) (3.20)
where K h (N rot , ṁh ) is the thermal conductance of the heat exchanger which is due to the convective heat transfer between the hot uid and the hot side heat exchanger wall (that depends on ṁh ), the thermal conductance through the heat exchanger walls and the convec- tive heat between the working gas inside the cycle and the other heat exchanger wall (that depends on N rot ). Neglecting the conductive heat transfer through the heat exchanger walls, K h (N rot , ṁh ) can be expressed as:

K h (N rot , ṁh ) = 1 1 K f luid h + 1 Kg (3.21)
where K f luid h is the heat conductance by convection between the hot uid and the hot side heat exchanger wall which is proportional to the hot uid mass ow rate and therfore can be set to:

K f luid h = K f luid,nom h .( ṁh ṁnom h ) α h (3.22)
where K f luid,nom h , ṁnom h and α h are some parameters to be identied.

K g is the heat conductance by convection between the working gas and the other heat exchanger wall that is inversely proportional to the rotational speed (when the rotational speed increases, the working uid spends less time in contact with the heat exchanger wall which decreases the heat transfer) and can be set as:

K g = K nom g .( N nom rot N rot ) αg (3.23)
where K nom g , N nom rot and α g are some parameters to be identied.

A similar analysis can be performed for the heat transfer between the cold sink and the compression chamber. Thus, the coolant enters the Stirling engine at temperature T c in with a mass ow rate ṁc and exits the engine at a higher temperature T c out absorbing a certain quantity of heat Q comp whose power Qcomp is given by:

Qcomp = ṁc .c c p .(T c in -T c out ) (3.24)
where c c p is the thermal capacity of the cooling uid.

Also, supposing a mean temperature of the cooling uid inside the cold uid side heat exchanger:

T C = T c in + T c out 2
The heat power rejected to the cold sink can be expressed as follows:

Qcomp = K c (N rot , ṁc ).( T c in + T c out 2 -T c ) (3.25)
where K c (N rot , ṁc ) is the heat conductance for the cold side heat exchanger given by: is the heat conductance by convection between the cold uid and the cold side heat exchanger wall which is proportional the cold uid mass ow rate ṁc and therefore can be set to:

K c (N rot , ṁh ) = 1 1 K f luid
K f luid c = K f luid,nom c .( ṁc ṁnom c ) αc (3.27)
where K f luid,nom c , ṁnom c and α c are some parameters to be identied.

Produced power and cycle eciency

The Stirling engine cyclic produced power P cyclic elec denoted for convenience P elec is given by the sum of the cyclic absorbed power Ẇcomp and the cyclic produced power Ẇexp which corresponds to the area of the closed curve 1 -2 -3 -4 of the Stirling engine P-V diagram (see Figure 3.4) multiplied by the engine frequency Nrot 60 as follow:

P elec = N rot 60 .(W comp + W exp ) = Ẇcomp + Ẇexp (3.28)
The cycle eciency is by denition given by the ratio of the produced power P elec and the absorbed power Qexp :

η cycle = P elec Qexp (3.29)

Summary of the Stirling engine equations

According to the analysis performed in previous sections, the Stirling engine equations that take into account the heat losses are summarized below:

Qexp = n cyl . N rot 60 .Q sch exp + n cyl .M g .c g p .(1 -ξ reg ).(T h -T c ) + K leak .(T h -T c ) (3.30a) Qcomp = n cyl . N rot 60 .Q sch comp -n cyl .M g .c g p .(1 -ξ reg ).(T h -T c ) -K leak .(T h -T c )(3.30b) -(A + B.N rot ).N rot (3.30c) Qexp = ṁh .c h p .(T h in -T h out ) (3.30d) Qexp = K h (N rot , ṁh ).( T h in + T h out 2 -T h ) (3.30e) Qcomp = ṁc .c c p .(T c in -T c out ) (3.30f) Qcomp = K c (N rot , ṁc ).( T c in + T c out 2 -T c ) (3.30g) P elec = Qexp + Qcomp (3.30h) η cycle = P elec Qexp (3.30i)
where:

ξ reg = 1 1 + N T U 0 .( Nrot N 0 ) αreg (3.31a) K h (N rot , ṁh ) = 1 1 K f luid h + 1 Kg (3.31b) K c (N rot , ṁc ) = 1 1 K f luid c + 1 Kg (3.31c)
and:

K g = K nom g .( N nom rot N rot ) αg (3.32a) K f luid h = K f luid,nom h .( ṁh ṁnom h ) α h (3.32b) K f luid c = K f luid,nom c .( ṁc ṁnom c ) αc (3.32c)
Q sch exp and Q sch comp are the indicated exchanged heat coming from the Schmidt analysis and given by (3.9a)-(3.9b) and (3.10).

Note from (3.30a)-(3.30c) that a multiplicative term n cyl has been introduced to account for the number of cylinder pairs that behave identically and linked to the rotating shaft.

In the previous equations the following set of parameters have to be identied experimentally (as described in section 3.3):

K leak , A, B, N T U 0 , N 0 , α reg , K nom g , N nom rot , α g , K f luid,nom h , ṁnom h , α h , K f luid,nom c
, ṁnom c and α c .

Stirling engine parameter identication

This section validates experimentally the nonlinear static model of the Stirling engine developed in Section 3.2 using experimental data gathered on a real testbed that is described below.

Description of the experimental testbed

Figure 3.7 gives an overview of the experimental Stirling based power system under consideration. The cooling loop depicted in Figure 3.8 that consists essentially of a circulating pump and a cooling fan serves as a cold sink to evacuate heat from the Stirling engine. The electrical heater (see Figure 3.7) serves to heat the oil (therminol 66) that enters the hot side Stirling engine heat exchanger. The power conversion system depicted in Figure 3.7 that is used to drive dynamically the Stirling engine acting on the electrical part will be presented in Chapter 4. To generate the data used to identify the unknown parameters of the Stirling engine

Identication of the unknown parameters of the Stirling engine model

Each static operating point of the 3kW elec (maximum electrical output at 300 • C for the hot uid temperature) alpha type Stirling engine is determined by xing a value for: T h in , T c in , ṁh , ṁc , P mean and N rot . To identify the unknown parameters of the Stirling engine model, we have collected measurements on the Stirling testbed corresponding to dierent operating conditions (covering all the operating range of the machine) by varying T h in , T c in , ṁh , ṁc , P mean and N rot over their admissible range with some discretization such that to have a good trade-o between the number of obtained operating points and the time spent to obtain these points. After eliminating the outlier points, it resulted 187 operating points used during the identication process.

Initially, we set the problem of nding the unknown parameters in the form of a nonlinear Figure 3.8: Zoom on the Stirling engine and the cooling loop curve-tting (data-tting) problem in least-squares sense. We tried to use the nonlinear curve tting routines of MATLAB R like "lsqcurvet" but we could not nd the optimal parameters because of the presence of multiple local minima for the objective function to be minimized. Afterwards, we used the functions of the "Global Optimization toolbox" of MATLAB R especially the "MultiStart" functions that enable to nd a global minimum among multiple local minima leading to the determination of the optimal parameters and hence solving the identication problem. The unknown Stirling engine model parameters given below:

K leak =2.6×10 -14 W/K A=2.22 × 10 -14 J B= 2.22 × 10 -5 J.s N T U 0 =0.1 N 0 =200 RPM α reg =-0.653 K nom g =563.5 W/K N nom rot =545.2 RPM α g =-1.699 K f luid,nom h =69.9 W/K ṁnom h =0.176 Kg/s α h =0.78 K f luid,nom c =611.45 W/K ṁnom c =0.971 Kg/s α c =0.493
These parameters make the Stirling engine model ts in a least squares sense the experimental measurements as illustrated in Figure 3.9.

To assess the parameter identication procedure quality and the model validity, a comparison between some of the Stirling engine model outputs and their corresponding measured signals is performed. Hence Figure 3.10 compares the experimental output electrical power, the absorbed and rejected heat powers with their corresponding model outputs for all the operating points while

Stirling engine torque expression

In order to derive a simpler expression of the Stirling engine motor torque T mot that will be used during the control design of the energy conversion system (that will be discussed in Chapters 5 and 6), the following expression for T mot (that has been initially suggested in [Meyer, 2011]) has been tted to the experimental data:

T mot = α( ṁh , ṁc , T h in , T c in ).Ω(t) + β( ṁh , ṁc , T h in , T c in ) (3.33)
where α(.) and β(.) are the functions identied from the experimental measurements, using the Eureqa-Formulize R equation tter software1 ( [START_REF] Schmidt | Distilling free-form natural laws from experimental data[END_REF]), and given by:

α( ṁh , T h in , T c in ) = 0.2498 + 0.004109.T h in + 3.13 × 10 -5 .T h in .T c in . ṁh -1.812. ṁh -1.655 × 10 -5 .T h in .T c in -0.007562.T c in ṁh β(T h in ) = 0.1.T h in -30.03
Figure 3.12 compares the tted motor torque (3.33) with the experimental measurements. One can see from this gure that the derived motor torque expression nicely ts to experimental measurements even if for some operating conditions the model errors can reach 20%. These identications errors will be dynamically estimated using a Luenberger observer as it will be explained in Chapter 5 (see Section 5.3.2) Figure 3.12: Comparison between the tted Stirling engine torque equation and the torque derived from the experimental data. each point in the gure corresponds to an operating condition dened by N rot , ṁh , ṁc , T h in and T c in

Optimal Stirling engine motor torque controller

The main control objective of the Stirling based power plant is to track some desired electrical power reference value. Since the power produced by the Stirling engine is the product of its motor torque T mot by its rotational speed Ω(t), one can control these two quantities to achieve the control objectives. Since these two quantities evolve in two separate time scales, it is possible to design the two controllers in a separate manner. The control of the Stirling engine rotational speed will be discussed in Chapters 5 and 6 when designing the control strategy for the energy conversion system. Concerning the control of the motor torque, this is done by acting on the hot and cold mass ow rates namely ṁh and ṁc (since the cooling fan speed is set to a constant value in our application and considering T h in and T c in as exogenous signals) through controlling the hot side and the cold side circulation pumps. By changing the setpoints of the hot and cold mass ow rates, the pumps will consume more or less electrical power. Since our main control objective is to supply the net electrical power requested by the loads (excluding the pumps), the mass ow rates setpoints are chosen such that to maximize the net thermoelectric eciency of the Stirling based power plant leading to the following formula:

( ṁopt h , ṁopt c ) = arg max ṁh , ṁc η net thermoelectric = P load elec -P pumps elec Qexp (3.34)
Since (3.34) is rather complex to solve online (due to the Stirling engine model complexity), Figure 3.13: Stirling engine optimal torque controller the optimization is performed o-line to construct two cartographies that give, according to the operating conditions (T h in and T c in ) and the load power demand P L , the optimal references for the hot and cold mass ow rates namely ṁopt h and ṁopt c as illustrated in Figure 3.13. The inputs of the cartographies, namely T h in , T c in and P L , are discretized such that to minimize the number of points in these cartographies (to minimize the time needed to construct these cartographies and the time needed during the interpolations) while having a good enough resolution. We have chosen for T h in and T c in a step size of 10 • C and for P L a step size of 100W .

Conclusion

In this chapter, the Stirling engine based power system under study was presented. Then, a static model of the engine taking into account the dierent losses was derived, identied and validated using experimental data gathered on the prototype. Finally, a cartography based controller was developed to optimally control the Stirling engine motor torque. In this chapter, the energy conversion system associated to the Stirling engine (see Figure 3.13) is presented and its dynamic model derived leading to a state space representation of this system that will be used in Chapters 5 and 6 to design advanced control strategies that optimally drive the power plant. Generator (PMSG) connected to the o-grid load to be supplied through an appropriate power conversion system (voltage conditioning stage). The Stirling engine is hybridized with a supercapacitor which serves as fast energy buer to dynamically supply the energy requested by the o-grid loads during fast transients while the Stirling engine adapts its produced power (by an appropriate control Strategy developed in Chapters 5 and 6) in a rather slower manner due to its mechanical and thermal inertia.

Energy conversion system architecture

The general principle of functioning of such system is as follows: the Stirling engine transforms the absorbed heat from the heat source into mechanical work (refer to Chapter 3) by performing a Stirling thermodynamic cycle at each shaft rotation. The produced Stirling torque is used to drive the PMSG that reacts electrically by producing an electromagnetic torque which together with the motor torque (and the friction torque) impose the rotational speed of the shaft. Because of the variable engine speed, the voltage at the output of the PMSG is of variable amplitude and frequency which is not directly usable by the loads. Hence, this voltage is rst rectied through a diode bridge before being conditioned through a DC/DC full bridge converter whose output voltage is maintained (through an appropriate control strategy) to a constant bus voltage V bus = 50V which is compatible with a good functionning of the inverter (linked to the DC bus) that supplies the o-grid load with 230V and 50Hz AC voltage.

The supercapacitor is connected to the DC bus through a bidirectional DC/DC converter that is controlled though its duty ratio u 2 = α bb to tightly regulate the DC bus voltage V bus . The control u 1 = α f b (duty ratio of the DC/DC full bridge converter) has then to be used to control the output current of the DC/DC full bridge converter namely I Lf b at some desired reference value computed so as to meet the load power demand and to enable the charge/discharge of the supercapacitor in order to recover its desired voltage V ref sc .

Dynamic modeling of the energy conversion system

In this section the dynamic equations for each component of the energy conversion system depicted in Figure 4.1 are presented and when needed proved in the appendices.

Dynamic model of the Permanent Magnet Synchronous Generator (PMSG)

In the literature, one can nd many papers and books that treat the modeling and simulation of dierent kinds of synchronous machines (see for instance [Chatelain, 1983], [START_REF] Mouni | Synchronous generator modelling and parameters estimation using least squares method[END_REF] and [START_REF] Barakat | Analysis of synchronous machine modeling for simulation and industrial applications[END_REF]). These models range from the most complex ones (detailed behavioral models) to simpler models. To simulate the dynamic behavior of the synchronous machine under consideration, the simplied Behn-Eshenbourg (see [START_REF] Bracikowski | Modélisation multi-physique par modèles à constantes localisées ; Application à une machine synchrone à aimants permanents en vue de son dimensionnement[END_REF]) model is considered whose per-phase schematic is depicted in Figure 4.2. According to this model, the 3-phase output voltage in the (abc) frame is given by:

[V ] 3 = [E v ] 3 -[R s ] 3 .[I] 3 -[L s ] 3 . d[I] 3 dt (4.1)
where: The 3-phase electromotive force (EMF) of the PSMG denoted by [E v ] 3 is proportional to the shaft rotational speed according to the following equation:

[V ] 3 = V a V b V c T is
[E v ] 3 =   E va E vb E vc   = p.φ f .Ω.   sin(p.Ω.t) sin(p.Ω.t -2.π 3 ) sin(p.Ω.t -4.π 3 )   (4.2)
where p is the number of pair poles of the PMSG, φ f is the EMF constant of the motor and Ω(t) is the shaft rotational speed in rad/s. The inductance and resistance matrices respectively [L s ] 3 and [R s ] 3 are given by:

[L s ] 3 =   L s M s M s M s L s M s M s M s L s   and [R s ] 3 =   R s 0 0 0 R s 0 0 0 R s  
where L s and R s are respectively the line inductance and line resistance of the PMSG and M s is the mutual inductance.

Since the electromagnetic torque of the PMSG corresponds to the ratio of its produced electrical power by its rotational speed, its equation is given by:

T em = E va .I a + E vb .I b + E vc .I c Ω(t)

Model of the full wave diode bridge

The diode bridge depicted in Figure 4.3 is an uncontrolled DC rectier whose mean average output voltage V red depends on the input voltage amplitude

V max = V max a = V max b = V max c
according to the following equation ( [Skvarenina, 2002]):

V red = 3. √ 3 π .V max

Modeling the association of the PMSG and the diode bridge

The models of the Synchronous generator and the diode bridge presented respectively in Sections 4.2.1 and 4.2.2 are used for the simulation of these components. In order to simplify further the model of the system consisting in the synchronous generator connected to the diode bridge, its equivalent DC machine model, which is a control oriented model, is considered which supposes an equivalent serially connected inductance (neglecting the mutual inductance eect), resistance and EMF voltage source of the association PMSG-diode bridge forming an equivalent DC machine electrical circuit whose equations are given by (see [Abdelli, 2007] for more details):

dI red dt = - R s L s .I red - 3.p 2.π .Ω(t).I red + 3. √ 3 2.π.L s .E m (t) - 1 2.L s .V red (4.3a) E m (t) = p • Ω(t) • φ f (4.3b) T em (t) = p. 3. √ 3 π .φ f .I red (t) (4.3c)
where V red is the rectied voltage, I red is the output current of the diode bridge, T em designates the electromagnetic torque of the PMSG, E m (t) its EMF (electromotive force) and φ f is the EMF constant of the motor.

Rotor dynamics

The dynamic equation of the shaft rotational speed results from the motor torque applied by the Stirling engine T mot whose equation is given by (3.33), the electromagnetic torque of the PMSG given by (4.3c) and the friction torque proportional to the speed according to T f r = D f r .Ω(t) where D f r is the friction coecient. This leads to the following well known equation for the speed dynamics:

J rot . dΩ(t) dt = T mot (t) -T em (t) -D f r .Ω(t) (4.4)
where J rot is the total inertia of the rotating parts. the Stirling Engine

Modeling the DC/DC full bridge converter

The DC/DC full bridge converter illustrated in Figure 4.4 is an isolated (due to the isolation transformer) buck transformer consisting in two IGBT (Insulated Gate Bipolar Transistors) legs at the primary of the transformer and a diode full-wave rectier at the secondary of the transformer (see [START_REF] Foch | Convertisseurs continu-continu isolées[END_REF] for more details). In order to model such switching device the mean model formalism is used [Skvarenina, 2002]. In this modeling approach, the averaged model (averaging of the signals and control variable) of the converter over a switching period T sw is considered in which the control variable becomes the duty ratio of the converter which is the ratio of the period on which some appropriate IGBT are enabled over the switching period. By applying this methodology on the considered DC/DC full bridge converter, its mean average model is derived (as explained in detail in Appendix B) and given by the following equations:

L f b . dI Lf b dt = k.V red .α f b -V bus (4.5a) C f . dV red dt = I red -I in f b (4.5b) I in f b = k.I Lf b .α f b (4.5c)
where α f b is the duty ratio of the Full bridge converter, L f b is the output inductance of the DC/DC full bridge, C f is the value of the input capacitor, k is the transformer winding ratio, I Lf b is the output current of the DC/DC full bridge and I in f b is the current at the input of the converter.

Modeling the bidirectional DC/DC converter

The bidirectional DC/DC converter illustrated in Figure 4.5 is a current reversible DC/DC converter that interfaces the supercapacitor to the DC bus in order to adequately control the supercapacitor delivered power. This converter behaves as a buck converter when its output current (current of the output inductance L bb ) denoted by I Lbb ows from the supercapacitor to the DC bus (positive current value) and behaves as a boost converter when the output current ows on the opposite direction (negative current value). Here also, we consider the average modeling approach of this device with a preference for the complementary control (see [Skvarenina, 2002]) that avoids the discontinuous conduction mode of the converter by activating both IGBT (but not simultanously) during a switching period T sw . Assuming that such precaution is taken in the hardware, the mean model used in the control design can be given by (see Appendix B for the detailed proof):

L bb . dI Lbb dt = V sc .α bb -V bus (4.6a) C sc . dV sc dt = -I Lbb .α bb (4.6b)
where α bb is the duty ratio of the bidirectional DC/DC converter, C sc is the supercapacitor capacitance, V bus is the DC bus voltage and V sc is the supercapacitor voltage.

DC bus equation

The equation of the DC bus is given by (see Figure 4.1):

C tot . dV bus dt = I Lf b + I Lbb -I inv (4.7a) C tot = C bb + C f b + C inv (4.

Modeling of the single phase inverter

The single phase inverter depicted in Figure 4.1 used in the testbed under consideration (as explained in Chapter 5) is a commercial product with its own integrated controller. Since we do not have any degree of freedom on this device, we will simply model it by an equation reecting the power conversion between its input and its output as follows:

P L = η inv .V bus .I inv (4.8)
where η inv is the electrical eciency of the inverter and P L is the load power demand.

Summary of the energy conversion system equations

The dynamic equations governing the behavior of the energy conversion system are summarized as follows: The set of physical equations (4.9) for the energy conversion system can be rewritten in a simpler state space form, to facilitate the control design, as follows:

J rot . dΩ dt = α( ṁh , T h in , T c in ).Ω + β(T h in ) -p. 3. √ 3 π .φ f .I red -D f r .Ω (4.9a) dI red dt = - Rs Ls .I red - 3.p 2.π .Ω.I red + 3. √ 3.p.φ f 2.π.L s .Ω - 1 2.Ls .V red (4.9b) C f . dV red dt = I red -k.I Lf b .α f b (4.9c) L f b . dI Lf b dt = k.V red .α f b -V bus (4.9d) C tot . dV bus dt = I Lf b + I Lbb - P L η inv .V bus (4.9e) L bb . dI Lbb dt = V sc .α bb -V bus (4.9f) C sc . dV sc dt = -I Lbb .
ẋ1 = a 1 .x 1 -a 3 .x 2 + a 2 (4.10a) ẋ2 = -a 4 .x 2 -a 5 .x 1 .x 2 + a 6 .x 1 -a 7 .x 3 (4.10b) ẋ3 = a 8 .x 2 -a 8 .k.x 4 .u 1 (4.10c) ẋ4 = -a 9 .x 5 + k.a 9 .x 3 .u 1 (4.10d) ẋ5 = a 10 .(x 4 + x 6 ) - a 10 η inv . P L x 5 (4.10e) ẋ6 = -a 11 .x 5 + a 11 .x 7 .u 2 (4.10f) ẋ7 = -a 12 .x 6 .u 2 (4.10g)
where

x 1 = Ω (rd/s), x 2 = I red (A), x 3 = V red (V), x 4 = I Lf b (A), x 5 = V bus (V), x 6 = I Lbb (A), x 7 = V sc (V).
The control variables are: u 1 = α f b and u 2 = α bb corresponding to both duty ratios of the DC/DC Full Bridge and bidirectional DC/DC converters respectively.

The coecients a i used in the state equations are given by:

a 1 = α( ṁh ,T h in ,T c in )-D f r Jrot , a 2 = β(T h in ) Jrot , a 3 = p.3. √ 3.Φ f π.Jrot , a 4 = Rs Ls , a 5 = 3.p 2.π , a 6 = p.3. √ 3.Φ f 2.π.Ls , a 7 = 1 2.Ls , a 8 = 1 Cf , a 9 = 1 L f b , a 10 = 1 Ctot , a 11 = 1 Lbb , a 12 = 1 Csc .

Control objectives and constraints

Based on the system description explained in Section 4.1 for the energy conversion system, we can translate the control objectives in terms of the notation given in Section 4.3.1 as follows:

• regulate in a sti manner x 5 =V bus around x st 5 =V ref bus =50V which is compatible with the good functioning of the commercial inverter.

• maintain x 7 =V sc around x st 7 =V ref sc , in a rather slower manner, to always have an adequate amount of energy in the supercapacitor. The energy conversion system is subject to the following system constraints: positivity constraints: x i ≥ 0 except for x 6 : Indeed, the bidirectional DC/DC converter is a current reversible converter, therefore the current x 6 = I Lbb can be either positive or negative to enable the charge/discharge of the supercapacitor. However, since the DC/DC Full Bridge converter is not current-reversible, its output current I Lf b is always positive. Moreover, its average current value must be greater than half of the ripple current in order to avoid discontinuous conduction mode, in which case the dierential equations for this converter become more complicated (see for example [START_REF] Sun | Averaged modeling of pwm converters operating in discontinuous conduction mode[END_REF]) in addition of having to switch between the equations for the continuous and the discontinuous conduction mode according to the value of the current.

strong saturation on the control variables since:

• u 1 ∈ [0, 1]: duty ratio of the DC/DC Full Bridge.

• u 2 ∈ [0, 1]: duty ratio of the bidirectional DC/DC converter.

Moreover, the energy conversion system is also subject to the following input disturbances:

• P L : load power demand that can vary in a fast manner.

• T h in and T c in : the input temperatures respectively at the hot side and cold side of the Stirling engine. These temperatures, that enter in the denition of the parameters a 1 and a 2 in (4.10a), have a slow variations compared to the dynamics of the electric stage.

Conclusion

In this chapter, the energy conversion system associated to the Stirling engine has been presented and its dynamic equations derived leading to a state space representation of the system that will be used in Chapters 5 and 6 to design ecient control strategies that will be validated in a real Stirling engine based power system prototype. Now that we have set in Section 4.3.1 the state space equations for the energy conversion system associated to the Stirling based power plant (see Figure 4.1), we present in this chapter a rst model based control strategy that enables the system to reach its control objectives as stated in Section 4.3.2 by controlling the electrical part of the power plant. Among the advantages of the proposed control strategy, the limited number of controller parameters to be tuned and the interesting way to generate the references for the key state variables that will be regulated in near to optimal manner. We present also in this chapter the experimental results that validate the performances of the developed method rst on a dedicated electrical motor testbed (Section 5.2) that reproduces the functioning of the Stirling engine thanks to a torque-controlled DC motor and then on the real Stirling based power system (Section 5.3).

Chapter 5. Development of a rst control strategy for the Energy Conversion

System associated to the O-Grid Stirling based power Plant.

Analysis and control of the energy conversion system

In this section the state space equations (4.10) of the energy conversion system will be analyzed and used to synthesize the control strategy that drives the Stirling based power system through its electrical part to track some desired reference load power demand P L of the O-grid load. These equations are recalled below for convenience:

ẋ1 = a 1 .x 1 -a 3 .x 2 + a 2 (5.1a) ẋ2 = -a 4 .x 2 -a 5 .x 1 .x 2 + a 6 .x 1 -a 7 .x 3 (5.1b) ẋ3 = a 8 .x 2 -a 8 .k.x 4 .u 1 (5.1c) ẋ4 = -a 9 .x 5 + k.a 9 .x 3 .u 1 (5.1d) ẋ5 = a 10 .(x 4 + x 6 ) - a 10 η inv . P L x 5
(5.1e) ẋ6 = -a 11 .x 5 + a 11 .x 7 .u 2 (5.1f) ẋ7 = -a 12 .x 6 .u 2 (5.1g)

where

x 1 = Ω, x 2 = I red , x 3 = V red , x 4 = I Lf b , x 5 = V bus , x 6 = I Lbb , x 7 = V sc , u 1 = α f b and u 2 = α bb .
The coecients a i used in the state equations are given by:

a 1 = α( ṁh ,T h in ,T c in )-D f r Jrot , a 2 = β(T h in ) Jrot , a 3 = p.3. √ 3.Φ f π.Jrot , a 4 = Rs Ls , a 5 = 3.p 2.π , a 6 = p.3. √ 3.Φ f 2.π.Ls , a 7 = 1 2.Ls , a 8 = 1 Cf , a 9 = 1 L f b , a 10 = 1 Ctot , a 11 = 1 Lbb , a 12 = 1 Csc .
First, a global system analysis will be performed leading to some system decomposition and then a controller will be designed for each subsystem while ensuring the stability of the global closed loop system.

Note that the parameters a 1 and a 2 are time varying with a very slow dynamics compared to the dynamics of the electrical part. Moreover these parameters are known since they depend on measured thermodynamic variables.

Global analysis of the system

The rst step in the system global analysis consists in performing some parametric study by replacing the system parameters by the values of the real testbed under consideration. By doing so, one can show that the term corresponding to the diode bridge voltage drop "a 5 .x 1 .x 2 " can be neglected in (5.1b) compared to the other terms of the equation.

By observing the structure of the state space equations (5.1), one can also see that it is possible to split the global system into two subsystems as follow:

1) The rst subsystem (subsystem 1) given by (5.1a)-(5.1d) involves only u 1 as a control variable and is driven by the exogenous signal a 2 (corresponding to the Stirling engine 5.1. Analysis and control of the energy conversion system 55 torque). This subsystem depends also on x 5 , but assuming that this variable is perfectly regulated around its reference (this will be proved later), x 5 can be replaced by its stationary (desired) value x st 5 . Therefore, the equations of subsystem 1 (5.1a)-(5.1d) can be rewritten in the following condensed form:

ż = A(u 1 ).z + B.
a 2 x st

5

(5.2)

where z = x 1 x 2 x 3 x 4 T is a reduced state vector.

2) The second subsystem (subsystem 2) corresponding to equations (5.1e)-(5.1g)

where u 2 is the only control variable being involved.

Note that the two subsystems evolve in a separate time scales. Indeed the equations of subsystem 1 involve the mechanical dynamics associated with the Stirling engine and the motor shaft which makes the evolution of subsystem 1 slower than that of subsystem 2 which involves only the electrical dynamics. Then, in the next subsections, equations of subsystem 2 will be used to tightly regulate x 5 =V bus at x st 5 =V ref bus =50V in a fast manner using the control variable u 2 . These equations will also serve to dene an appropriate reference value for x 4 = i Lf b namely x st 4 that will adapt the state of charge of the supercapacitor by maintaining x 7 =V sc around its setpoint value x st 7 =V ref sc . Finally, the equations of subsystem 1 will be used to derive an appropriate stationary reference value for u 1 and z namely u st 1 and z st (u st 1 ) respectively that will be tracked acting on the control variable u 1 and using the equations of subsystem 1.

Analysis and control of subsystem 2

As mentioned in Section 5.1.1 equations of subsystem 2 (5.1e)-(5.1g) are used to reach one of the control objectives consisting in maintaining tightly x 5 around its setpoint value x st 5 =V ref bus =50V. Since (5.1e)-(5.1f) are in the strict-feedback form, a simple backstepping controller (see [Khalil, 1996]) will be designed for this purpose as follows:

Step 1: In (5.1e), we use the variable x 6 as a virtual control to regulate x 5 =V bus around

x st 5 .
The regulation error being 5 = x 5x st 5 , we use the following positive denite function:

V 5 = 1 2 .
2 5 as a Lyapunov function where its time derivative is given by:

V5 = 5 . ˙ 5 = 5 .(a 10 .(x 4 + x 6 ) - a 10 η inv . P L x 5 )
By setting V5 = -ρ 5 . 2 5 < 0, with ρ 5 > 0 is a design parameter, we get :

x ref 6 = P L η inv .x 5 -x 4 + ρ 5 a 10 .(x st 5 -x 5 ) (5.3)
Step 2: By setting V6 = -ρ 5 . 2 5ρ 6 . 2 6 < 0, with ρ 6 another design parameter, we get:

u 2 (x) := 1 a 11 .x 7
.(a 11 .x 5a 10 .(x 5 -

x st 5 ) + ẋref 6 -ρ 6 .(x 6 -x ref 6 )) (5.4)
Hence, (5.3) and (5.4) locally stabilizes x 5 around its setpoint x st 5 =V ref bus in a fast manner leading rapidly to the following stationary equation ( ẋ5 ≈ 0):

x 6 = -x 4 + P L η inv .x st 5 (5.5)
In order to meet the second control objective consisting in maintaining the supercapacitor voltage x 7 around its reference value x st 7 in a rather slower manner, one can choose the following stationary value for x 6 which together with equation (5.1g) obviously stabilizes x 7 around x st 7 :

x st 6 = k 6 . tanh(β.(x 7x st 7 ))

(5.6)

where k 6 and β are some design parameters.

Replacing (5.6) in (5.5), the following stationary desired value for x 4 namely x st 4 is derived:

x st 4 = -k 6 . tanh(β.(x 7 -x st 7 )) + P L η inv .x st 5 (5.7)
In the next subsection, an appropriate control strategy for subsystem 1 will be designed to track x st 4 which enables to regulate the supercapacitor voltage and to respond to the load power demand.

Analysis and control of subsystem 1

Now that we have dened according to (5.7) the desired stationary value for x 4 namely x st 4 , one can determine online the remaining components of the desired stationary state vector z st and the stationary value for u 1 namely u st 1 by solving analytically the following equation corresponding to ż = 0 (see (5.2)): This leads to the following equations:

z st (u st 1 ) := -[A(u st 1 )] -1 B.
x st 2 = -B+ √ ∆ 2.A , x st 1 = a 3 .x st 2 a 1 , x st 3 = x st 4 .x st 5 x st 2 , u st 1 = x st 2 k.x st 4
where: A = a 3 .a 6 a 1a 4 , B = -a 2 .a 6 a 1 , C = -a 7 .x st 4 .x st 5 and ∆ = B 2 -4.A.C.

Figure 5.1 shows the components of the stationary vector z st with respect to u st 1 .

By computing the eigenvalues of the dynamic matrix A(u st 1 ) for some values of u st 1 ∈ [0, 1], it appears that subsystem 1 is open-loop stable (see Figure 5.2 where the stability indicator of A(u st 1 ), that corresponds to its maximum eigenvalues real part for dierent u st 1 ∈ [0, 1], is plotted) but contains highly oscillatory modes as it can be observed from Figure 5.3 where some open-loop trajectories for the components of the vector z are plotted for three dierent values of u st 1 . Figure 5.4 illustrates the same behavior in the phase plane (x 3 , x 4 ). One can conclude from Figures 5.3 and 5.4 that starting from some initial state, in order to reach some desired state vector z st compatible with some desired x st 4 , applying directly its corresponding u st 1 would almost systematically lead to positiveness constraints violation of the state variables (see Section 4.3.2) that are additionally oscillatory with uncontrolled transients.

Then, it is necessary to design an appropriate control strategy other than u 1 = u st 1 to reach some desired stationary vector z st compatible with some desired x st 4 . For this purpose, we opted for a one step (because of the limitation on the controller computation time that is xed to 100 µs) predictive control [START_REF] Mayne | Constrained model predictive control:stability and optimality[END_REF]] that leads to near to optimal behavior of the System associated to the O-Grid Stirling based power Plant. First, the subsystem 1 equations given by (5.2) are discretized based on the second order approximation of the exponential term leading to the following discrete version of subsystem 1 equations:

z + (u 1 ) = A d (u 1 )z + B d (u 1 ).
a 2 x st

5

(5.9)

Then, to track some desired stationary reference value z st compatible with the desired x st 4 and its corresponding u st 1 , the following one step predictive controller is used:

u opt 1 = arg min u 1 ∈[0,1] J(u 1 ) := z + (u 1 ) -z st 2 P d (u st 1 )
(5.10)

where P d (u st 1 ) is the Lyapunov stability matrix of subsystem 1 (depending on u st 1 ) used as a weighting matrix in the previous cost function.

From preliminary simulations of the closed-loop subsystem 1 under the control u 1 given by (5.10), it appears that x 4 undergoes an undershoot (or an overshoot) during positive (or negative) step changes in x st 4 , which can lead to positiveness constraints violation. To tackle this problem, the reference value x st 4 given by (5.7) is ltered (using an ad hoc rst order lter) to smooth the convergence of x 4 . In addition to this, a corrective term (slow integrator) is added to x st 4 (before the ltering) to correct the static error on x 4 (which is a crucial state variable that adapts the state of charge of the supercapacitor) due to some model uncertainties. This leads to the following modied reference value for x 4 namely x4 st in the discrete form:

x4 st (k + 1) = α f .(x st 4 (k) + (k)) + (1 -α f ). x4 st (k) (5.11a) (k + 1) = k aw . (k) + η corr .(x st 4 (k) -x 4 (k)) (5.11b)
where α f is a ltering parameter, η corr is the integrator gain and k aw is an anti windup gain.

Global control of the system

Figure 5.5 summarizes the global control architecture of the energy conversion system derived from the combination of the subsystem 1 controller (Section 5.1.3) and subsystem 2 controller (Section 5.1.2). The resulting closed loop system is stable and its performances assessed in the next sections.

Closed loop simulations

In order to assess by simulation the performances of the proposed control strategy for the energy conversion system, a load power demand prole is applied to the system as depicted in Figure 5.6. One can observe from Figures 5.7a and 5.7b that the controller tightly regulates the DC bus voltage to its setpoint enabling the good functioning of the inverter that supplies the AC load with a 230V 50Hz AC voltage. According to Figure 5.7c, the proposed controller also maintains the supercapacitor voltage to its reference value enabling having sucient amount of energy to palliate to the load power demand transients. and its instantaneous value I Lf b (t). One can see from this gure that this current is always positive enabling a continuous conduction mode of the device. Figure 5.8a plots the shaft rotational speed Ω(t) that is always stabilized to its setpoint derived from x st 4 .

In Figure 5.8c, one can see the output current of the bidirectional DC/DC converter I Lbb . At t = 2s, this current becomes positive to compensate for an increase in the load power demand (Figure 5.6) that appears between time t = 2s and t = 15s resulting in a discharge of the supercapacitor (see Figure 5.7c). Between t = 2s and t = 15s, I Lbb takes negative and positive values to respectively charge and discharge the supercapacitor in order to maintain its voltage at the desired value. The same remarks hold for the time interval [15s, 30s] where we rst have a decrease and then an increase in the load power demand.

Figures 5.9a and 5.9b give the time proles of the control variables u 1 and u 2 . One can observe from Figure 5.9b that some peaks appear on u 2 , however by zooming on these peaks one observe the same behavior as in Figure 5.7b which is admissible for the system. System associated to the O-Grid Stirling based power Plant.

Experimental Validation on a dedicated electrical motor testbed

Before testing the proposed control strategy on the real Stirling based power system, an experimental electrical motor testbed was built to serve as a rst experimental validation platform in which the behavior of the Stirling engine is reproduced by a DC machine torquecontrolled by mean of a commercial drive and linked mechanically to the synchronous generator which is connected electrically to the energy conversion system we want to control.

The need for this intermediate platform is due to two main reasons. The rst one being that we did not have in our disposal the real Stirling engine at the time we developed the control strategy but we wanted to experimentally validate it in another manner. The second reason is that we wanted to rst validate the control strategy on a more exible testbed especially for motor torque variations and the electrical motor testbed seemed to be the best choice in order to ensure this exibility.

Description of the electrical motor testbed

Figure 5.10 gives the schematic of the electrical motor testbed under consideration. The energy conversion system which is the same as in Figure 4.1 is connected electrically to a synchronous generator linked to a DC machine (reproducing the functioning of the Stirling engine) torque-controlled by mean of a commercial drive. A programmable AC load connected at the output of the inverter is used to generate a variable load prole. Figure 5.11 gives an overview of the experimental testbed for which the following system parameters hold:

PMSG nominal power: 1.8kW, J rot =0.0672Kg.m 2 , α( ṁh , ṁc )=-0.0111, β(T H in , T C in )=45, φ f =0.6014, R s =5, L s =0.28 mH, p=8, C f =2400 µF , L f b =150µH, L bb =223µH, C f b =50µF , C bb =136mF , C sc =63 F .
This leads to the following numerical values for a i parameters: a 1 = -0.183, a 2 = 558.11, a 3 = 118.4453, a 4 = 9615.4, a 5 = 1.3712,a 6 = 5101.1, a 7 = 641.02, a 8 = 425.53, a 9 = 6666.7, a 10 = 7.34, a 11 = 4484.3, a 12 = 0.0159, η inv = 0.95 and k = 0.5.

Because of some physical limitations on the DC/DC full bridge converter used in the testbed, the duty ratio of this device belongs to the interval α f b ∈ [0, 0.9]. This duty ratio is rst conditioned before it is applied to drive the DC/DC Full Bridge converter as it is illustrated in Figure 5.12. In this gure, α f b is rst divided by two before generating its corresponding electrical signal α elec f b in order to drive the two pairs of Insulated Gate Bipolar Transistors (IGBT) (T r 1 ,T r 3 ) and (T r 2 ,T r 4 ), during a cutting frequency T sw = 20 kHz, with the same period but 180 • out of phase (see Section B.1 for more details).

Note also that the bidirectional DC/DC converter used in the testbed does not work in a complementary mode but during a switching period only one IGBT is activated at a time. However, the mean model (5.1f)-(5.1g) is still valid for control design with a slight System associated to the O-Grid Stirling based power Plant. modication before driving the converter. Indeed, the real duty ratio applied to the converter α real bb ∈ [-0.95, 0.95] is derived from α bb as follow:

α real bb = α bb - V bus V sc (5.12)
Finally, the duty ratio is transformed to a PWM signal using standard PWM techniques similar to the one of Figure 5.12 where T r 1 (see Figure 4.5) is activated during a fraction of the switching period corresponding to | α real bb | .T sw if α real bb is positive otherwise T r 2 will be activated during the same fraction of the switching period. Concerning the real time controller on which the developed control strategy is implemented, we opted for a National Instrument controller based on NI PXI-8109 (http://www.ni.com/pdf/manuals/373172a.pdf) real time target (2.66 GHz Dual-Core embedded controller). The control strategies are transcribed in MATLAB/Simulink R environment and ".dll" le embedding a C code is generated thanks to Simulink coder R . The generated le is then deployed on the National Instrument real time target using the National Instrument Veristand R platform (installed on a host PC) that is also used to monitor the system and make data acquisitions.

Experimental results on the electrical motor testbed

In order to assess the performances of the developed control strategy on the experimental electrical motor testbed, a load power demand prole is generated acting on the electronic programmable load (see Figure 5.10), the resulting power prole measured at the input of the inverter is shown in Figure 5.13. Figure 5.14 shows that despite the sharp changes in the load power demand, the DC bus voltage is well regulated around its reference value of 50V, acting on u 2 = α bb (depicted in the bottom of Figure 5.18), which guarantees the good functioning of the inverter that supplies the AC load with uninterruptible 230V 50Hz voltage.

The supercapacitor voltage is also maintained to its setpoint according to Figure 5.15. This enables the system to always have enough energy (and free energy stock) to supply (or absorb) during future increase (or decrease) in the load power demand.

The reference stationary value for the output current of the DC/DC full bridge converter namely x st 4 needed to adapt the SOC (State Of Charge) of the supercapacitor (see Section 5.1.2) is also tracked with a good settling time according to that are tracked as illustrated in Figure 5.17 acting the duty ratio of the DC/DC full bridge converter u 1 = α f b depicted in the top of Figure 5.18. From Figure 5.17, one can observe some small static errors in the tracking of

x ref 1 = Ω ref ,x ref 2 = I ref red and x ref 3 = V ref red
which is due to the modication of the reference signal for x 4 according to (5.11). However this static error will not aect the global performances of the system since the crucial objective being the regulation of x 4 = I Lf b .

Experimental Validation on the Stirling engine based power system

Now that we have validated on the dedicated electrical motor testbed (see Section 5.2) the proposed control strategy for the energy conversion system, we will discuss in this section the experimental results obtained on the real Stirling engine based power plant prototype. First, a brief description of this prototype will be given focussing on the energy conversion system, then a slight modication of the developed control strategy will be addressed taking into account a corrective term for the Stirling engine torque obtained through an additional state observer. Finally, the experimental results on the Stirling engine power plant will be given and analyzed.

Energy conversion system prototype description

The Stirling based power system prototype under consideration is the one depicted in Figure 3.7 whose thermodynamic part was described in Section 3.3.1. The power conversion system whose schematic is depicted in Figure 4.1 and that was used in the experimental electrical System associated to the O-Grid Stirling based power Plant. The second electrical cabinet depicted in Figure 5.21 contains among others the diode bridge and its associated over-voltage protection system (windy boy protection box1 ). This cabinet contains also a commercial variable speed drive (ATV322 variable speed drive) used during the starting procedure (see Appendix C) to initially bring the Stirling engine to some desired speed in order to enable motor torque production.

Corrective term observer for the Stirling engine torque

In Section 3.3.3 an expression of the Stirling engine torque T mot was identied using experimental data resulting in equation (3.33) for the motor torque that was then used in equation (4.9a) leading to the state equation (5.1a) in which a 1 and a 2 contain the terms of the motor torque. Because of some model errors in the torque identication process, it has been added a simple Luenberger observer [Besançon, 2007] that estimates the model error in the Stirling torque identication ∆T mot that has to be added to the calculated torque T calc mot given by (3.33). Preliminary results without adding this observer showed that the closed loop system under the Chapter 5. Development of a rst control strategy for the Energy Conversion System associated to the O-Grid Stirling based power Plant. proposed controller was subjected to some tracking errors on the state variables (x 1 ,x 2 ,x 3 ,x 4 ). This justies the design of the proposed observer that starts from rewriting (5.1a) which is recalled below: ẋ1 = a 1 .x 1a 3 .x 2 + a 2

where:

a 1 = α( ṁh , ṁc)-Dfr Jrot , a 2 = β(T H in ,T C in ) Jrot , a 3 = p.3. √ 3.Φ f π * Jrot , in the following form: ẋ1 = a 1 .x 1 + a calc 2 + γ (5.13)
where:

a 1 = -D f r Jrot , a calc 2 = T calc mot
Jrot , γ = -a 3 .x 2 + ∆a 2 and ∆a 2 = ∆Tmot

Jrot

Considering the slow dynamics on x 2 (and hence on γ), the following augmented system is considered:

ẋ1 = a 1 .x 1 + a calc 2 + γ
(5.14a) γ = 0 (5.14b)

y obs = x 1 (5.14c)
that can be rewritten in the following condensed form considering X obs = x 1 γ T :

X obs = A obs .X obs + B obs .a calc 2
(5.15a)

y obs = C obs .X obs (5.15b)
where A obs = a 1 1 0 0 , B obs = 1 0 and C obs = 1 0

Then the following Luenberger observer is applied on the discrete version of (5.15) to estimate X obs : Finally, the corrective term on the motor torque is obtained as follow:

Xobs (k + 1) = A d
∆a 2 = γ + a 3 .x f 2 (5.17)
where x f 2 is obtained by ltering x 2 using a rst order lter in order to get just the slow dynamics of this variable which is compatible with the dynamics of the Stirling engine torque.

Experimental results on the Stirling based power system prototype

In this section, we present the experimental results on the real Stirling engine prototype depicted in Figure 3.7. Similarly to Section 5.2.2 that validates the developed control strategy on the dedicated electrical motor testbed, a load power prole depicted in Figure 5.22 is generated acting on the electronic load. One can appreciate from Figure 5.23 that the bus voltage is tightly regulated to its reference value despite the step changes in the load power demand. According to Figure 5.24 the supercapacitor voltage is also maintained around its reference value in a fast manner. The measurements disturbances on the DC bus voltage (Figure 5.23) are due to electromagnetic compatibility problems that occurred after setting the DC/DC converters in the electrical cabinets. A low pass lter (with high cutting frequency) has been introduced to limit these measurement disturbances even if these ones have not been totally attenuated. However, the DC bus voltage remains within its admissible range (44V to 56V). Note from Figure 5.22 that between t = 24s and t = 27s the load has been disconnected which explains that during this time interval the supercapacitor voltage keeps increasing (see Figure 5.24). During this time interval and according to Figures 5.25 and 5.26 the state variables x 1 = Ω, x 2 = I red , x 3 = V red and x 4 = I Lf b do not track anymore their respective references but are limited to some lower bound due to the duty ratio α f b saturation as it can be seen from the top of Figure 5.27. If the load is kept disconnected, the supercapacitor voltage will keep increasing until it reaches its maximum value leading to the system shut down (or activating dump loads). One can also see from Figure 5.26 that the measurement of the rotational speed x 1 = Ω is also disturbed in high speed (high frequency) due to electromagnetic compatibility problems. These disturbances have been greatly attenuated using a low pass lter.

The obtained time responses of the closed loop system are very satisfying for the considered application. The reason is that the developed control strategy is based upon an optimal criterion (even if the prediction horizon is one step) and the time constant of the rst order lter given by (5.11a), that limits the undershoots and overshoots eects, has been chosen as low as possible.

Concerning the oscillations observed in Figures 5.25 and 5.26, they are due to the rapid changes in the load power demand (Figure 5.22) which is an extreme behavior for the system. Indeed, at some time instants, before the system reaches its steady state it sees another step change in the load power demand. Despite this extreme behavior, the closed loop system remains stable with a very satisfying time constants. System associated to the O-Grid Stirling based power Plant.

Conclusion

In this chapter, a novel control strategy for an O-Grid hybrid Stirling engine/supercapacitor power plant was developed relying on a specic way to generate the references for the key state variables that were tracked combining a backstepping controller and a one step predictive controller. The performances of the closed loop system were assessed rst by simulations before their implementation on a dedicated electrical motor testbed. Once the control methodology was validated on this testbed, the nal step consisting in implementing the control strategy on the real Stirling based power system was successfully achieved. In Chapter 5, a novel control strategy that drives the Stirling based power system through its electrical part was presented and experimentally validated. Although the proposed methodology endows the closed-loop system with near to optimal performances in term of settling time, it does not handle explicitly the system constraints. Instead, saturations and ad hoc lters has been added to palliate these constraints and enlarge the domain of applicability of the proposed solution.

In this chapter, a constrained control methodology for the energy conversion system is presented. The resulting feedback laws are real time compatible and they handle explicitly the constraints on the state variables. Moreover, a nice feature in the proposed solution is that it suggests a general control design methodology for a family of power circuits despite the fact that its success may depend on the quantitative design of the circuit components. The proposed design underlines the conditions under which the framework is successful and can therefore be also used to guide the design step. 82 Chapter 6. Constrained control framework for the Energy Conversion System associated to the O-Grid Stirling based thermodynamic power Plant

Deeper in the system analysis

This section provides additional results to complement the ones provided in Section 5.1. These results will serve in further sections to construct the constrained feedback laws for the energy conversion system whose state space equations given by (4.10a-4.10g) are recalled thereafter for convenience: where

ẋ1 = a 1 .x 1 -a 3 .x 2 + a 2 (6.1a) ẋ2 = -a 4 .x 2 -a 5 .x 1 .x 2 + a 6 .
x 1 = Ω, x 2 = I red , x 3 = V red , x 4 = I Lf b , x 5 = V bus , x 6 = I Lbb , x 7 = V sc , u 1 = α f b and u 2 = α bb .
The coecients a i used in the state equations are given by:

a 1 = α( ṁh ,T h in ,T c in )-D f r Jrot , a 2 = β(T h in ) Jrot , a 3 = p.3. √ 3.Φ f π.Jrot , a 4 = Rs Ls , a 5 = 3.p 2.π , a 6 = p.3. √ 3.Φ f 2.π.Ls , a 7 = 1 2.Ls , a 8 = 1 Cf , a 9 = 1 L f b , a 10 = 1 Ctot , a 11 = 1 Lbb , a 12 = 1 Csc .
First of all, note that since the Stirling torque variations are assumed to be very slow compared to the electric stage characteristic times and since we are interested in controlling optimally the energy conversion system, the Stirling engine torque T mot is supposed hereafter to be constant resulting in a constant parameters a 1 and a 2 in (6.1a).

Note also that since the supercapacitor used in the experimental testbed presented in Section 5.3 was oversized, we consider here for the simulations a value of 5F for the supercapacitor in order to assess more the performances of the proposed methodology especially for the supercapacitor voltage recovery.

For recall, the control objectives consist in regulating x 5 =V bus and x 7 =V sc around their respective values x r 5 = x st 5 =V ref bus and x st 7 =V ref sc during load power demand P L transients while respecting the following set of constraints:

x i ∈ [x min i , x max i ] ; i ∈ {2, 3, 4, 6} (6.2)
The constraints on the regulated variables x 5 and x 7 will be naturally respected by the success of the control loop while the controlled excursion on x 1 is implicitly imposed through the constraints on both x 2 and x 3 as it can be inferred from equations (6.1a)-(6.1b).

The same system decomposition presented in Section 5.1.1, for the energy conversion system, holds for this chapter: "subsystem 1" which is given by (6.1a)-(6.1d) where u 1 is the control variable and x 5 is assumed to be regulated around its setpoint (by mean of the control u 2 using "subsystem 2" equations) is rewritten in the condensed form (5.2) recalled thereafter:

ż = A(u 1 ).z + B. a 2 x st 5 (6.3)
where z = x 1 x 2 x 3 x 4 T being the reduced state vector for this subsystem and the "subsystem 2" equations are given by (6.1e)-(6.1g). As mentioned before, the "subsystem 2" equations will serve to regulate x 5 to its setpoint using u 2 and to dene an appropriate reference value for x 4 namely x st 4 that will be tracked using "subsystem 1" equations acting on u 1 .

The following assumptions are formulated following [START_REF] Alamir | Constrained control framework for a stand-alone hybrid (stirling engine)/supercapacitor power generation system[END_REF] and apply throughout this chapter. Assumption 6.1

The linear system in the state (x 1 , x 2 ) dened by (6.1a)-(6.1b) (steered by the exogenous signal x 3 ) shows two separated negative real eigenvalues λ f and λ s such that |λ f | |λ s |.

For the studied system, this assumption is satised with |λ f | ≈ 9552 Hz |λ s | ≈ 63.4 Hz. This feature is likely to be satised on any motor testbed as it reects the co-existence of mechanical (x 1 ) and electrical (x 2 ) variables which generally evolve among separated time scales. Assumption 6.2

The lower bound for x 3 has to satisfy x min 3 > x st 5 k Indeed, this condition has to be satised to ensure the controllablity of x 4 using u 1 ∈ [0, 1] as it can be inferred from (6.1d). Assumption 6.3

The bounds x min This assumption can be easily checked through Figure 5.1 where the steady states of (6.1a)-(6.1d) (in which x 5 = x st 5 ) namely x st i (u 1 ), i ∈ {1, . . . , 4} are plotted. Moreover, one can observe from this gure that a wide variation range of x 1 , x 3 and x 4 corresponds to a small variation range of x 2 denoted hereafter by 2 dened as follow:

ε 2 := max u 1 ∈[0,1] x st 2 (u 1 ) -min u 1 ∈[0,1]
x st 2 (u 1 ) (6.4)

The central value of x st 2 (•) is denoted hereafter by x st,c 2 . Consequently,

x st 2 (u 1 ) ∈ I 2 := x st,c 2 + [- ε 2 2 , + ε 2 2 ].
Note that the parameters ε 2 and x st,c associated to the O-Grid Stirling based thermodynamic power Plant and does not depend on any control design. For the system under consideration, x st,c 2 ≈ 4.82 while ε 2 ≈ 0.12 which corresponds to a relative variation of less than 2.5% (see Figure 5.1).

Another feature that can be observed in Figure 5.1 is that all the maps are monotonic so that the following notations can be used:

x st i (x st j ), u st 1 (x st j ) (i, j) ∈ {1, . . . , 4}
to denote the stationary values of x i and u 1 that correspond to a given stationary value x j of the j-th component of the state vector z. This is extensively used in the sequel.

According to Assumption 6.3, it is possible to ensure the fulllment of the constraints on x 2 = I red by monitoring the excursion of the exogenous signal x 3 (exogenous with respect to the subsystem dened by (6.1a)-(6.1b)). This is enunciated in the following lemma: Lemma 6.1

There are easily computable maps x 3 (x 1 , x 2 ) and x3 (x 1 , x 2 ) such that all trajectories of (6.1a)-(6.1b) starting from an admissible initial conditions (regarding the constraints on x 2 ) and steered by a prole x 3 (•) meeting the following condition:

x 3 (t) ∈ x 3 x 1 (t), x 2 (t) , x3 x 1 (t), x 2 (t) (6.5)
satisfy the saturation constraints on x 2 .

Proof: See Appendix D.1 where explicit expressions (D.2)-(D.3) are given for x 3 (x 1 , x 2 )

and x3 (x 1 , x 2 ).

The next two lemmas constitute a formalization of the results obtained in Section 5.1.2. Lemma 6.2

Any appropriate control that asymptotically steers x 5 to its reference value x r 5 leads to the following property:

lim t→∞ [x 6 (t) -x r 6 (t)] = 0 (6.6)
where x r 6 (t) is given by:

x r 6 (t) = P L η inv • x r 5 -x 4 (t) (6.7)
where P L is the load power demand.

Proof: This directly follows from equation (6.1e) in which x 5 is replaced by x r 5 while ẋ5 = 0 is used.

It comes from Lemma 6.2 together with equation (6.1g) that it is possible to regulate the supercapacitor voltage x 7 by choosing the following reference value for x 4 :

x st 4 = -k 6 • tanh(β • (x 7 -x st 7 )) + P L η inv • x st 5 ; k 6 > 0 (6.8)
since according to (6.6)-(6.7) this asymptotically infers the following dynamics to x 7 :

ẋ7 = a 12 • u 2 • k 6 • tanh(β • (x 7 -x st 7
)) (6.9) which obviously stabilizes x 7 around x st 7 .

Note however that for (6.9) to be compatible with the constraints on x 4 , the following conditions are needed:

x min 4 + k 6 ≤ P L η inv • x r 5 ≤ x max 4 -k 6 (6.10)
This obviously introduces constraints on the possible load power levels that are compatible with a sustainable operation. This is rigorously stated in the following assumption: Assumption 6.4

There is a positive real k 6 such that all possible values P L of the load power demand satisfy the inequalities (6.10).

Lemma 6.3

The state feedback law given by:

u * 2 (x) := 1 a 11 .x 7 .(a 11 .x 5 -a 10 .(x 5 -x r 5 ) + ẋref 6 -ρ 6 .(x 6 -x ref 6 )) (6.11)
where:

x ref 6 = P L η inv • x 5 -x 4 + ρ 5 a 10 (x r 5 -x 5 ) (6.12)
in which ρ 5 and ρ 6 are positive design parameters, locally asymptotically steers the pair (x 5 (t), x 6 (t)) to the reference (x r 5 , x r 6 (t)), namely: (6.13) Proof: This result was proved in Section 5.1.2.

lim t→∞ x 5 (t) x 6 (t) - x r 5 x r 6 (t) = 0
In the following sections the previous results will serve to design constrained sliding mode controllers based on vector elds analysis for each subsystem.

6.2 Constrained Sliding Mode control of the energy conversion system 6.2.1 Analysis and constrained sliding mode control for "subsystem 1"

Now that we have dened the reference value for the output current of the DC/DC full bridge converter namely x st 4 = I ref red given by (6.8) leading to the asymptotic stabilization of the

u 1 = 1 u 1 = 1 (x ⇤ 3 , x ⇤ 4 ) u ⇤ 1 u 1 = 1 u 1 = 1 ū1 A B C D E x 3 x 4 x 4 = x max 4 x 4 = x min 4 ū1 u 1 = 0 u 1 = 0 u 1 = 0 u 1 = 0 R 34 L 34
Figure 6.1: Typical disposition of the vectors eld of the dynamic system (6.1c)-(6.1d) when x 2 = x st 2 (x * 3 ) and x 5 = x r 5 are assumed to be constant.

3. DC is the trajectory under the control u 1 = 0 that starts at D and ends at C (this denes D).

4. BC is the trajectory under the control u 1 = 1 that starts at B and ends at C (this denes B).

5. The disposition of the vector elds on the line AB and DE are such that the control ū1 (x) given by: (6.15) steers the state towards D (resp. B) when x 4 = x max 4 (resp. x min 4 ) while keeping x 4 constant. 6. The line ABCDE is attractive under the control u 1 = 1 for any starting point that lies to the right of the curve BCD (this region is denoted by R 34 on Figure 6.1).

ū1 (x) :=          1 kx 3 x 5 + ρ 4 (x max 4 -x 4 ) on DE 1 kx 3 x 5 + ρ 4 (x min 4 -x 4 ) on AB
7. The line ABCDE is attractive under the control u 1 = 0 for any starting point that lies to the left of the curve BCD (this region is denoted by L 34 on Figure 6.1) associated to the O-Grid Stirling based thermodynamic power Plant (x r 5 , x r 6 )

x 5

x 6 = x min

6 x 6 = x max 6 x 6 u 2 = 1 u 2 = 0 u 2 = 0 u 2 = 1 ū2 ū2 u ⇤ 2 A B C D E L 56 R 56
Figure 6.2: Typical disposition of the vectors eld of the dynamic system (6.1e)-(6.1f) for a given x 4 and for x 7 > x r 5 .

6.2.2 Analysis and constrained sliding mode control for "subsystem 2"

By observing the structure of (6.1e)-(6.1f) it appears that the same graphical analysis as in Section 6.2.1 can be done in the 2D plan (x 5 , x 6 ) for a given value of x 4 . This analysis is relevant under (6.14) since as it is explained earlier, equation (6.14) leads to rather slow variations of (x 3 , x 4 ) inducing hence two separated time scales, a slow one for x i , i ∈ {1, 2, 3, 4, 7} on one hand and a fast one for x 5 and x 6 on the other hand.

When considering the vector elds of (6.1e)-(6.1f) for a given value of x 4 , it can be shown that the following assumption holds for the studied system: Assumption 6.6 For a given admissible x 4 and a supercapacitor voltage x 7 > x r 5 , the disposition of the vector elds around the desired pair (x r 5 , x r 6 ) where x r 6 is dened by (6.7) is the one given in Figure 6.2, more precisely:

1. The horizontal lines AB and DE correspond to the lines x 6 = x min 6 and x 6 = x max 6 respectively.

2. C is the desired steady point (x r 5 , x r 6 ) while u * 2 is the locally stabilizing feedback invoked in Lemma 6.3: 3. DC is the trajectory under the control u 2 = 0 that starts at D and ends at C (this denes 6.2. Constrained Sliding Mode control of the energy conversion system 91 D). 4. BC is the trajectory under the control u 2 = 1 that starts at B and ends at C (this denes B).

5. The disposition of the vector elds on the line AB and DE are such that the controls ū2 given by: (6.20) steers the state towards D (resp. B) when x 6 = x max 6 (resp. when x 6 = x min 6 ) while keeping x 6 constant. 6. The line ABCDE is attractive under the control u 2 = 1 for any starting point that lies to the left of the curve BCD (this region is denoted by L 56 on Figure 6.2).

ū2 (x) :=          1 x 7 x 5 + ρ 6 a 11 (x max 6 -x 6 ) on DE 1 x 7 x 5 + ρ 6 a 11 (x min 6 -x 6 ) on AB
7. The line ABCDE is attractive under the control u 2 = 0 for any starting point that lies to the right of the curve BCD (this region is denoted by R 56 on Figure 6.2).

Note that the same comments that followed assumption 6.5 can be reproduced here regarding the items invoked in Assumption 6.6. In particular, there is here also the possibility to check them through explicit conditions similar to (D.4)-(D.7).

Based on the observation of Figure 5, an interpolation-based state feedback K 2 (x|P L ) can be dened in a similar way to the one already used in Denition 6.1 to dene K 1 (x|P L ).

Denition 6.2

Let be given a load power demand P L . Consider the corresponding desired pair (x r 5 , x r 6 ) using

x r 6 dened by (6.7). Compute the following quantities based on the current point p = (x 5 , x 6 ) T :

• Compute u DE 2 (x) and u AB 1 (x) according to:

u DE 2 (x) := ϕ( x 5 (D) -x 5 x 5 (D) , 0, ū2 (x)) (6.21) u AB 2 (x) := ϕ( x 5 -x 5 (B) x 5 (B) , 1, ū2 (x)) (6.22)
where ū2 (x) is given by (6.20).

• Compute u BCD 2 (x) according to:

u BCD 2 (x) :=        ϕ( x max 6 -x 6 x max 6 , 0, u * 2 (x)) if x 6 > x r 6 ϕ( x 6 -x min 6 |x min 6 | , 1, u * 2 (x)) if x 6 ≤ x r 6
where u * 2 (x) is given by (6.11).

u v 2 (x) :=        ϕ( x max 6 -x 6 x max 6 , u DE 2 (x), 1) if x 5 < x5 ϕ( x 6 -x min 6 |x min 6 | , u AB 2 (x), 0) if x 5 ≥ x5 • The feedback u 2 = K 2 (x|P L ) is obtained by horizontal interpolation between u v 2 (x) and u BCD 2 (x)
according to:

u 2 (x|P L ) = K 2 (x|P L ) := ϕ( |x 5 -x5 | |x 5 | , u BCD 2 (x), u v 2 (x)) (6.23)

The State Feedback law

Note that by denition of the interpolated feedback laws K 1 (x, P L ) and K 2 (x, P L ), the resulting closed-loop trajectories passes through neighborhoods V 34 (x * 3 , x * 4 ) and V 56 (x r 5 , x r 6 ) of the targeted pairs (x * 3 , x * 4 ) and (x r 5 , x r 6 ) of size O(1/λ) after a nite time. In other words by choosing λ > 0 (used in the interpolation map (6.16)) suciently high, the trajectories cross neighborhoods that are as small as required. In particular, λ > 0 can be taken suciently high so as to make these neighborhoods entirely contained in the basins of attraction of the locally asymptotically stable dynamics under the local control laws u st 1 (x * 3 ) and u * 2 (x). If then, the feedbacks u st 1 (x * 3 ) and u * 2 (x) are maintained, then asymptotic stability follows.

The preceding discussion suggests a dual-mode controller. Namely, when the state (x 3 , x 4 ) [resp. (x 5 , x 6 )] is far from the desired pair, the interpolated law K 1 [resp. K 2 ] is used to steer the state to small neighborhood of the desired pair. When the state is inside this neighborhood, a switch to the locally asymptotically stabilizing law

u st 1 (x * 3 ) [resp. u * 2 (x)] is operated.
It is worth underlying also that if the neighborhoods are taken suciently small (that is λ is taken suciently high), then the transient excursions still meet all the saturation constraints on the states. The overall control law is therefore admissible.

Note however that such a dual-mode control design has to be done with care in order to avoid undesirable chattering-like behavior (the state enters and leaves the small neighborhood with a high frequency). This is generally avoided by dening the size of the nal neighborhoods using the Lyapunov functions level of the local feedback. This is more rigorously dened as follows: Denition 6.3

For each admissible x * 3 , the Lyapunov matrix reecting the stability of A(u st 1 (x * 3 )) is denoted ). Moreover, for any given pair (x 3 , x 4 ), the following function is dened:

d 34 (x, P L ) :=     0 0 x 3 -x * 3 x 4 -x st 4 (x * 3 )     S 34 (x * 3 ) (6.24)
where x * 3 is given by (6.14).

Denition 6.4

Dene S 56 ∈ R 2 to be the Lyapunov matrix corresponding to the following Hurwitz matrix:   -

a 10 P L η inv (x r 5 ) 2 a 10 0 -ρ 6   (6.25)
Moreover, for any pair (x 5 , x 6 ), the following function can be dened:

d 56 (x, P L ) := x 5 -x r 5 x 6 -x r 6 (x, P L ) S 56 (6.26)
where x r 6 is dened by (6.7).

By now, we have everything we need to dene the complete state feedback: Denition 6.5

Let λ > 0 be the positive scalar invoked in the denition of the interpolation map (6.16). Let γ 34 > 0 and γ 56 > 0 be some positive reals. Let P L be a load power demand. The state feedback law is dened as follows:

u 1 := K 1 (x, P L ) if d 34 (x, P L ) > γ 34 u st 1 (x * 3 (x, P L )) otherwise (6.27) u 2 := K 2 (x, P L ) if d 56 (x, P L ) > γ 56 u * 2 (x, P L ) otherwise (6.28)
Based on the discussions above, the following result can be proved: Proposition 6.1

If Assumptions 6.1-6.6 are satised, then there exist suciently small γ 34 > 0, γ 56 > 0 and suciently high λ × max(γ 34 , γ 56 ) > 0 such that the closed-system obtained by applying the state feedback laws (6.27)-(6.28) to the dynamic system (6.1a)-(6.1g) asymptotically stabilizes the system's state at a steady value that is compatible with the desired values x r 5 and x r 7 of x 5 and x 7 respectively while meeting the problem's constraints (6.2) and using admissible control u ∈ [0, 1] 2 .

Proof: See Appendix D.3. Chapter 6. Constrained control framework for the Energy Conversion System associated to the O-Grid Stirling based thermodynamic power Plant

Simulation-based validation

In this section some simulations are proposed to illustrate several facts regarding the proposed state feedback design. More precisely, the validation is split into two stages:

1. In the rst stage, referred hereafter as the decoupled simulations, each feedback law is validated separately by decoupling the two subsystems. This is done by simulating:

the system with state (x 1 , x 2 , x 3 , x 4 ) T in which x 5 = x d 5 is used. the system with state (x 5 , x 6 , x 7 ) with a given and constant value of x 4 .

These decoupled simulations enable us to show the impact of the saturation on the dierent state variables on the settling time of the corresponding controlled system to be shown. This stage obviously shows that there are two separate time scales which explains the success of the whole scheme.

2. In the second stage, the whole closed-loop system is precisely simulated under varying power demand in order to assess the ability of the control law to tightly regulate the bus voltage x 5 under varying power demand P L while regulating in a slower time scale the super-capacitor voltage x 7 . ∈ {50, 60}.These scenarios clearly shows that closed-loop behavior shows a sort of minimum-time behavior since it saturates the current in order to accelerate the convergence towards the desired value. Figure 6.4 shows how the state x 4 (and therefore x 3 ) are steered to the desired steady pair (x d 3 , x d 4 ) which corresponds to +20% increase in the value of x d 4 while meeting the constraints on all the state variables. In particular, a minimum value of x min 2 = 4 A is considered. This gures clearly shows how the bounding values dened by (D.2)-(D.3) which guarantee the fulllment of the constraints on x 2 slow down the possible evolution of the pair (x 3 , x 4 ). Figure 6.5 shows the same kind of results when a decrease of the desired value x d 4 of -20% is considered. Note that in this case, the upper bound x max 2 = 5 A is saturated in order to accelerate the convergence.

Decoupled Simulations

It is worth noting that during the scenarios of Figures 6.4 and 6.5, the state variable x 4 itself does hit the constraints during the rst milliseconds. This is dicult to see because of the time scale. This is the reason why the same evolutions are plotted in Figures 6.8 and 6.9 using a logarithmic time scale. 3) slow down the response in order to meet the constraints on x 2 . The same scenario is depicted in Figure 6.8 using logarithmic scale in order to show the evolution of the variables during rst milliseconds of the scenario. 3) slow down the response in order to meet the constraints on x 2 . The same scenario is depicted in Figure 6.9 using logarithmic scale in order to show the evolution of the variables during rst milliseconds of the scenario.

Global Simulation

The decoupled set of simulations proposed in section 6.4.1 clearly shows that the time scale of the two subsystems are quite separated. Indeed, the response time of the bus voltage related subsystem is in the range of milliseconds while the response time of the current x 4 = I Lf b that is used to charge the supercapacitor voltage is in the range of seconds. This justies the decoupled design. Nevertheless, the global control design is validated in the present section through closed-loop simulation of the whole system under varying power demand scenario. Figures 6.10, 6.11 and 6.12 show three dierent scenarios that investigate the overall closedloop performance under time varying power load demand proles. In particular, Figures 6.10 and 6.11 show two scenarios that dier only by the maximal allowed value x max 2 on the current x 2 = I red which is set to x max 2 = 5 A and x max 2 = 4.8 A respectively. Note how this slight dierence has a great impact on the amplitude of the excursion of the supercapacitor voltage x 7 = V sc . This suggests that the maximum allowable current x 2 which heavily depends on the choice of the component has to be jointly designed with the choice of the supercapacitor.

All the simulations clearly show that the control objectives are met in the sense that the .9: Same scenario as the one depicted on Figure 6.5 using logarithmic scale on the time axis in order to show the beginning of the scenario. This scale enables the saturation on x 4 to be clearly shown.

bus voltage is tightly regulated while the supercapacitor voltage is restored as soon as the constraints on the dierent control and state variables enable it.

Conclusion

In this chapter, a global constrained control framework is proposed for the stand-alone hybrid (stirling engine)/supercapacitor power generation system. The control architecture explicitly handles saturations on the control variables as well as on the components of the state vector.

The resulting closed-loop behavior can be used in the design step where the choice of the components and the class of power demand proles have to be investigated. This is because such choices highly impact the size of the supercapacitor that can aord the dierent power demand proles.

Beyond the specic system that has been considered in the present paper, it is our claim that the control design methodology can be extended to a wide variety of power conditioning stages that can be encountered even when other thermodynamic cycles and/or dierent storing devices are used.

Introduction

As discussed in Section 1.1, the nal prototype of the solar thermodynamic micro power plant selected by Schneider Electric consists in having an Organic Rankine Cycle (ORC) power plant as thermodynamic machine transforming the water heat into electricity. This ORC turbine was regulated by EXOES (which is one of the MICROSOL project partners) using PI based controllers associated to static maps (for feedforward terms) determined experimentally oline. Two main drawbacks were identied on the control methodology used by EXOES. The rst drawback is the slow response time of the turbine while moving from an operating point to another during load power demand changes. Indeed, this impacts signicantly the sizing of the energy buer consisting in the battery banks that adapts the power consumed by the loads during transients. The other drawback is the way the maps that determine the feedforward terms are computed. Indeed, the methodology consists in browsing experimentally all the possible static operating points (obviously by mean of an adequate step values) and selecting the optimal ones in term of turbine eciency and minimum auxiliaries consumption. although this methodology is the most commonly used by the thermal engineers since it is experience-based, it is very time consuming and it has to be repeated each time a serious modication is performed on the ORC cycle.

In the next chapters, we present a control methodology that tackles the aforementioned drawbacks. The methodology consists in rst identifying a control oriented nonlinear dynamic model of the ORC power plant using real data gathered on the prototype, then an optimal controller in terms of settling time and minimum consumption respecting the system constraints will be designed based on the identied models. Among the advantages of the proposed methodology, the response time of the system is signicantly decreased which has a very good impact on the battery sizing. Indeed, the time response obtained by EXOES using their PI-based control methodology was about 30 s (to respond to a 1kW electrical power demand) while this response time was reduced by half according to simulation results obtained using the proposed constrained optimal control strategy. The proposed methodology can also be viewed as a quasi systematic method for regulating this kind of thermodynamic engines while reducing the necessary time for developing the control laws compared to the traditional methods.

Chapter 7

Nonlinear Dynamic Model Identication of an ORC power plant In this chapter, an identication framework recently developed in [Alamir et al., 2014c], that was used to derive a nonlinear black-box model of a beta-type Stirling engine with rhomicdrive mechanism, is presented. The simulation results presented in [Alamir et al., 2014c] have shown good agreement between the knowledge-based model (used as a data generator) and the identied model all the more that this last one was very fast to simulate which is a very interesting feature for the identied nonlinear models. This identication framework will be applied in this chapter to derive nonlinear maps for each variable of interest for the ORC power plant under study. Then, a state space representation is derived from these nonlinear maps and a linearized model that approximates the nonlinear one is obtained and used in the controller design presented in Chapter 8.

In the literature, one can nd many identication frameworks that deal with nonlinear models like Wiener-Hammerstein, Volterra Serries based formulations (see [START_REF] Doyle | Identication and Control Using Volterra Models[END_REF]) and many others. However, the adopted identication framework stands out from the others by its simplicity in term of model complexity and by the fact that it results in solving a 109 110 Chapter 7. Nonlinear Dynamic Model Identication of an ORC power plant Figure 7.1: Picture of the ORC power plant prototype Quadratic Programming (QP) problem which is quite easy to solve online, then enabling the models to be updated during the system lifetime which is an interesting feature while performing online system diagnosis for example. This chapter is organized as follows: rst, the ORC system under consideration, its principle of functioning, the control objectives and system constraints are presented, then a brief review on modeling such systems is given. Thereafter, the adopted Wiener-based identication framework is discussed and the identication results on the ORC power plant prototype are given and analyzed to assess the performance of the identication framework. Finally, a state space model of the ORC system under consideration is derived.

Presentation of the ORC power plant prototype

The solar thermodynamic micro power plant under consideration that is depicted in Figure 7.1 is a prototype that aims at producing both electricity (thanks to an ORC system) and drinking water (thanks to an appropriate water desalination system provided by a project partners). A set of photovoltaic (PV) panels are also present as a second energy source (which is crucial during the starting procedure). In this work, we focus on the study of the electrical energy production thanks to the ORC process depicted in Figure 2.6 and whose schematic As for the Stirling engine based power plant, the ORC system consists in a hot water loop (hot source) and a cooling loop (cold sink) interfaced through heat exchangers to an internal loop wherein ows the working uid (R245Fa organic uid).

The hot water loop (top of Figure 7.2) whose elements are not represented in Figure 7.2 consists of a pump that makes the hot water circulating between the evaporator and a storage tank. Another loop (which is out of the scope of this work) enables hot water to circulate between the storage tank and the solar collector eld. The cooling loop, on the other hand, enables a water-glycogen mixture to circulate in order to absorb heat from the internal loop (cold sink). The main elements of this loop are the circulation pump and the cooling fan.

The principle of operation of this ORC power plant can be described using both the T-S diagram of Figure 7.3 and the schematic of Figure 7.2. Staring at the volumetric pump inlet (point 1 of Figure 7.2 and Figure 7.3), the working uid is in subcooled liquid state. This uid undergoes an isentropic (meaning some pump eciency) increase in pressure through the volumetric pump to end up at point 2. The uid then absorbs through the heat recuperator (that enhances the cycle eciency) between points 2 and 2 a certain quantity of heat released 112 Chapter 7. Nonlinear Dynamic Model Identication of an ORC power plant by the same uid at the output of the turbine and then enters the evaporator where it absorbs heat from the hot water (Hot water loop) between points 2 and 3. During this step, the working uid is heated and brought to its liquid saturation line where it continues absorbing heat in two-phase thermodynamic state (mixture of liquid and vapor in thermodynamic equilibrium). In the last stage of this step, the uid is superheated to avoid the presence of droplets in the turbine (which can be damageable). At point 3 the working uid is at the inlet of the "scroll" turbine (see [Quoilin, 2011] for more details) and undergoes an isentropic expansion between points 3 and 4 resulting in an output mechanical work W turbine used to drive an asynchronous generator producing an electrical power P elec that is conditioned through an appropriate power electronics (not represented in Figure 7.2). Finally the working uid begins its condensation by rst releasing a quantity of heat between points 4 and 4 through the heat recuperator (the same quantity of heat absorbed between points 2 an 2 ) and enters the condenser where it releases heat to the external cooling loop to end up with a subcooled liquid at the input of the volumetric pump.

As it can be seen in Figure 7.2, the ORC system prototype consists in two scroll turbines: a big turbine (15 kW nominal power) and a small turbine (4 kW nominal power and 5kW maximum power). The two turbines can work either alone or in parallel. In the next sections, we focus our study (model identication and control) on the small turbine (the big turbine being disconnected from the hydraulic circuit). Since the two turbines are from a similar nature, the same analysis can be performed with the big turbine. 7.2. Control objectives and system constraints 113

Control objectives and system constraints

The ORC power plant is designed to operate either in grid-connected mode or in o-grid mode. A supervisor (which is out of the scope of this study) that manages the load consumption and the transition between the two modes of operation generates an electrical power reference P ref elec , based on the batteries SOC (State Of Charge), that the ORC power plant has to track acting on the pumps (volumetric pump, hot water pump and circulation pump) and the cooling fan. The response time of the ORC system in tracking this power setpoint is crucial for the battery sizing especially when the system performs in o-grid mode since the battery banks dynamically meet the load power demand while the ORC turbine adapts its produced power P elec in a rather medium time scale.

The ORC system is subject to a set of constraints in term of minimum and maximum temperatures, pressures and actuator limitations. These constraints are as follows:

• The uid at the inlet of the turbine has to be in superheated vapor state with a superheating 3 • C < ∆T vap < 15 • C such that:

∆T vap = T vap -T sat (P vap ) (7.1)
where [8bar, 25bar] are respectively the temperature and the pressure of the working uid at the outlet of the evaporator (inlet of the turbine) and T sat (P vap ) is the saturation temperature of the working uid (which is a function of the evaporation pressure) whose values will be given by thermodynamic tables for the R245Fa uid provided by the REFPROP R software (http://www.nist.gov/srd/ nist23.cfm). This constraint is crucial since it avoids droplets to be present at the input of the turbine.

T vap ∈ [85 • C, 135 • C] and P vap ∈
Note that T sat (P vap ) can be calculated online using the following fourth order polynomial tted by REFPROP R data (see Figure 7.4):

T sat (P vap ) = -0.0003.P 4 vap + 0.023.P 3 vap -0.79.P 2 vap + 14.32.P vap + 279.49 (7.2)

• The temperature and pressure at the outlet of the turbine respectively T eM ot and P eM ot has to be such that T eM ot < 110 • C and 1.7bar < P eM ot < 7bar. The upper bound on P eM ot is due to the mechanical restrictions on the volumetric pump while the lower bound insures that the the working uid at the outlet of the condenser is in subcooled state to avoid pump cavitation.

• The temperature of the cooling uid (water glycol) at the output of the cooling fan T f r must be such that 30

• C < T f r < 80 • C.
Moreover, the ORC power plant is subject to the following input measured disturbances (exogenous signals): • T aCh: temperature in Kelvin of the hot water at the inlet of the evaporator (hot water side).

• T amb : ambient temperature in Kelvin.

The control variables (manipulated variables) used to achieve the control objectives while respecting the previous constraints are dened below together with their domain of variation:

• f pump ∈ [1, 50] (Hz): rotational frequency of the motor pump (volumetric pump).

• Vhot ∈ [0, 6000] (L/h): volume ow rate setpoint for the hot water side pump.

• V aero rot ∈ [10, 100]: rotational speed of the ventilator in (%) relative to its maximal value.

• V circ rot ∈ [5, 100]: rotational speed of the circulation pump (cold side) in (%) relative to its maximal value.

These control variables have also to be penalized in the expression of the cost function in order to minimize the electrical consumption of the auxiliaries.

To summarize, the control objective for the ORC power plant is to track some dened electrical reference power P ref elec while minimizing the consumed power in the auxiliaries (pumps 7.3. Brief state of the art on physical and semi physical ORC power plants dynamic models 115

and fans) and respecting the following set of constraints:

3 • C ≤ ∆T vap ≤ 15 • C (7.3a) 85 • C ≤ T vap ≤ 135 • C (7.3b) 8bar ≤ P vap ≤ 25bar (7.3c) T eM ot ≤ 110 • C (7.3d)
1.7bar ≤ P eM ot ≤ 7bar (7.3e)

30 • C ≤ T f r ≤ 80 • C (7.3f) 1 ≤ f pump ≤ 50Hz (7.3g) 0 ≤ Vhot ≤ 6000l/h (7.3h) 10% ≤ V aero rot ≤ 100% (7.3i) 10% ≤ V circ rot ≤ 100% (7.3j)
Before designing a controller that achieves the previous control objectives while respecting the constraints, a model for the ORC power plant under consideration has to be derived. In the literature, one can nd many papers that deal with such problem as it is discussed in Section 7.3.

Brief state of the art on physical and semi physical ORC power plants dynamic models

Because of their simplicity and the availability of their components, ORC systems have been widely used in many applications essentially in solar thermal and waste heat recovery (WHR) applications. One can nd in the literature many papers that deal with the optimal sizing of such systems (see [Quoilin, 2011[START_REF] Quoilin | Performance and design optimization of a low-cost solar organic rankine cycle for remote power generation[END_REF] for example) but few documents exist that deal with dynamic modeling and control of ORC power plants. In [Quoilin et al., 2011a], a physical modeling of a small scale ORC power plant for waste heat recovery was proposed focusing on the dynamic modeling of a simple evaporator which is one of the main element (and most dicult to model) in the ORC system. Since the model developed in [Quoilin et al., 2011a] is too computationally heavy, it was not dedicated to real time control but was implemented in Modelica R to serve as a simulator. A static model was also derived to compute feedforward terms that optimize the static ORC eciency while simple tunable PID based controllers were used to reach the control objectives. Since the heat exchangers, especially the evaporator, is the most challenging part to model in the ORC power plant, many papers can be found that deal with this problem. Two main paradigms exist to model the evaporator: the nite volume (FV) paradigm and the Moving Boundary (MB) paradigm.

In [Quoilin et al., 2011a] the nite volume paradigm was used to model a simple evaporator consisting of a hot uid and a cold uid (working uid) channels separated by a wall. The evaporator was decomposed in a set of nite volumes where the mass, energy and heat balances were applied to each nite volume resulting in a large number of nonlinear ordinary dierential equations. In [START_REF] Feru | Heat exchanger modeling and identication for control of waste heat recovery systems in diesel engines[END_REF], the same paradigm was applied to model a simple evaporator Figure 7.5: Moving Boundary paradigm applied to a condenser [START_REF] Shah | Dynamic Modeling And Control of Single and Multi-Evaporator Subcritical Vapor Compression Systems[END_REF].

(for diesel engines WHR applications), but instead of using external thermodynamic tables (as in [Quoilin et al., 2011a]) data tted polynomials were used.

Concerning the Moving Boundary (MB) paradigm depicted in Figure 7.5, it consists in decomposing the evaporator (or the condenser) into three regions (liquid, two phase and superheated) whose length is made time dependent and for each region a mass and energy balance is applied to end up with the following nonlinear state space representation (see [START_REF] Shah | Dynamic Modeling And Control of Single and Multi-Evaporator Subcritical Vapor Compression Systems[END_REF], Rasmussen and Alleyne, 2006, Peralez et al., 2012] for more details):

Z(x, u). ẋ = f (x, u) (7.4)
In [START_REF] Shah | Dynamic Modeling And Control of Single and Multi-Evaporator Subcritical Vapor Compression Systems[END_REF]Alleyne, 2006] for example, the authors decided to take for the condenser model a state vector x(t) that depends on the pressure, the regions lengths, the output enthalpy and the wall temperatures for each region. In [START_REF] Mckinley | A switched system model for heat exchangers using a moving boundary method[END_REF], the author proposed a switching model methodology that tackles the problem of model instability that occurs when one of the thermodynamic state for the working uid becomes much smaller compared to the other thermodynamic phases (which makes Z(x, u) of 7.4 non invertible).

Concerning the models of the other components, [Quoilin, 2011, Twomey et al., 2013] proposed a static models for the scroll expander and the pump that resulted in a large number of parameters to be identied oine using experimental data.

In our approach, we do not consider such physical or semi physical modeling approaches that typically lead to complex and parameter sensitive models. We prefer using a model identication approach that leads to a parsimonious model which means a model that incorporates only the necessary elements needed to design the feedback law. These models can also serve to build fast simulators for the ORC power plants.

Hence, Section 7.4 presents the proposed nonlinear identication framework used to derive nonlinear relationships that link each of the variables of interest: regulated variable P elec and constrained variables (∆T vap , T vap , P vap , T eM ot, P eM ot and T f r ) to some or all of the system control variables (f pump , Vhot , V aero rot and V circ rot ) and the exogenous signals (uncontrolled measured input disturbances) T aCh and T amb . In this section the identication framework introduced in ( [Alamir, 2013], [Alamir et al., 2014a] and [Alamir et al., 2014c]) will be used to derive nonlinear dynamic models for the ORC power system under study. One of the major advantages of this identication methodology is that it results in solving a QP problem and ending up with an easy to use nonlinear model as it can be seen in later sections.

Presentation of the methodology

Consider the identication problem consisting in nding a nonlinear map E(.) that links some desired identied variable q(i) to some: state vector x(i), control vector u(i) and a vector of measured uncontrolled signals w(i). Then, one has to nd E(.), such that:

q(i.τ s ) = q(i) ≈ E(x(i), u(i), w(i)) (7.5)
where τ s is the sampling period of the signals.

As for many identication methods, one has to dene some regressor Z(i) that contains the following extended vector: U (i) = x(i) u(i) w(i) ∈ R n U and its (N -1) past values (N being the regressor buer length) as follows:

Z(i) =      U (i) U (i -m) . . . U (i -(N -1).m)      ∈ R nz ; n z = N.n U
then (7.5) can be rewritten in term of a new nonlinear map F (.), that has to be found, as follows:

q(i) ≈ F (Z(i)) (7.6)
where n U is the length of the vector U (i) and m is some under sampling integer. Suppose now, that one has N data learning data samples used during the identication process.

The identication problem can be dened as follows:

Denition 7.1

Given the learning data:

(q(i), Z(i))

i=N data i=1 ; K = {1, . . . , N data }
the identication problem consists in nding a map F : R nz → R that approximately maps Z to q in a least squares sense, namely In the literature, one can nd many nonlinear structures for the map F (.). In the present work, we constrain F (.) to take the following form:

F (Z) = Γ -1 (Z T .L); Γ(.) streactly increasing (7.8)
Where L ∈ R nz is a parameter vector to be identied together with the nonlinear function Γ(.). The restriction regarding the monotonicity of Γ(.) comes from the need for this function to be invertible. Thus, applying the function, Γ(.) to each side of (7.8) results in:

Γ(q) ≈ Z T .L (7.9)

In [Alamir et al., 2014a], it has been chosen to parametrize Γ(.) using the following structure:

Γ(q) = n b j=1
B j (ξ(q)) .µ j ; ξ(q) := qq min q maxq min ∈ [0, 1] (7.10)

where q min and q max are the minimum and maximum values of q over the learning data and µ j the set of parameters to be identied .The set of function basis B j (.) dened on [0, 1] are given by:

{B j } j=n b j=1 = {1} ∪ {B l 1 } nm-1 l=1 ∪ {B l 2 } nm l=1 (7.11)
in which n b = 2n m is the number of functions in the basis and:

B l 1 (η) = (1 + α l ). η 1 + α l .η ; B l 2 (η) = η 1 + α l .(1 -η) (7.12)
where the coecients α l are given by: α l = exp β.(l-1) -1 for some constant β > 0. Figure 7.6 shows the allure of these maps for β = 0.5.

By putting (7.10) in the following matrix form:

Γ(q) = B(ξ(q)) .µ (7.13)

Where B(ξ(q)) ∈ R 1×n b is the vector whose elements are all the {B j (ξ(q))} j=n b j=1 and µ ∈ R n b is the vector of parameters to be identied, (7.9) can be transformed into:

B(ξ(q)), -Z T .   µ L   ≈ 0 (7.14)
which has to be solved in a least squares sense under the following constraints on µ:

• Constraints ensuring that Γ(.) is strictly increasing: (7.15) for some > 0

∀ξ ∈ [0, 1], dB(ξ) dξ .µ ≥
• Normalization constraint to avoid the trivial solution: µ = 0 and L = 0:

1 0 B(ξ)dξ .µ = q min + q max 2 (7.16)
By dening a suciently dense grid of ξ values on [0, 1], (7.15)-(7.16) can be transformed into a nite set of linear constraints which together with (7.14) dene the QP problem that has to be solved to derive the optimal parameters namely: L ∈ R nz and μ ∈ R n b

Note that the total number of model parameters for each nonlinear map is:

n p = n b + n z = 2n m + n z (7.17)
Note also that a good indicator for the degree of compactness of the learning data into the identied nonlinear map is given by the compression factor which is equal to:

compressionf actor = N data n p (7.18)
7.4.2 Application of the identication framework to the ORC power system

Now that we have dened in Section 7.4.1 the identication framework, we will apply this methodology to the case of the ORC power plant described in Section 7.1. But before going further into details, it is important to emphasize an interesting feature regarding our ORC system that holds also for any thermodynamic machine in general. This feature concerns the inherent open-loop stability of such systems. Indeed, for a given value of the control vector, 120 Chapter 7. Nonlinear Dynamic Model Identication of an ORC power plant the ORC system will always converge to some steady state, not necessarily the desired one, but a stable state anyway. Then, one can exploit an interesting feature of open-loop stable systems which is the asymptotic disappearance of the initial state eect as time progresses which implies that the ORC system state vector x(i) can be approximated as a function of the input proles (controlled and uncontrolled signals) during the last N samples, provided that N is chosen suciently high and the sampling period suciently small, as follows:

x(i) = G u(i), u(i -1)...u(i -N + 1) (7.19)
This implies that for each variable of interest y to be identied which could be the regulated variable P elec or one of the constrained variables (∆T vap , T vap , P vap , T eM ot, P eM ot and T f r ), one can identify a nonlinear relationship that links y(i) to some regressor Z y (i) that depends only on the vector u y (i) (n y u dimentional vector containing some controlled and uncontrolled measured signals) and its (N y -1) past values such that:

Z y (i) =      u y (i) u y (i -m) . . . u y (i -(N y -1).m)      ∈ R n y z ; n y z = N y .n y u (7.20)
To perform the identication task, we have in our disposal experimental data sampled at τ s = 1s coming from two experiments namely "Experiment 1" and "Experiment 2" performed on the ORC prototype. For each variable y to be identied, we combine a data subset of "Experiment 1" and a data subset of "Experiment 2" to generate the set of learning data (y(i), Z y (i))

i=N data i=1 while the whole data set will serve as validation data. These learning data have to contain the maximum and minimum values of the whole data set in order to make the identied models valid on this range. Figure 7.7 and Figure 7.8 show the input variables for Experiment 1 and Experiment 2 respectively. These variables will be used in the denition of the regressor for each identied variable.

Identication Results

For each variable of interest to be identied, namely P elec , T vap , P vap , T f r , P eM ot, T eM ot and ∆T vap , one has to build a regressor as dened in (A.8) used in the identication process. The following input vectors were chosen to construct the regressors for the identied variables:

u Tvap (i) = f pump Vhot T aCh T amb T (7.21a) u Pvap (i) = u Tvap (i) (7.21b) u P elec (i) = f pump Vhot V aero rot V circ rot T aCh T amb T (7.21c) u T f r (i) = u P eM ot (i) = u T eM ot (i) = u P elec (i) (7.21d)
By performing several identication tests, the parameters of Table 7.1 gave the best iden- and7.12-7.15 that a good matching is obtained between the experimental data and the identied nonlinear models for the variables of interest. Note that the learning data used to derive the nonlinear maps during the identication process are plotted in red color.

The models resulting from the identication of the maps for T f r , P eM ot, T eM ot and ∆T vap were linear. However, concerning the variables T vap , P vap and P elec , nonlinear models were needed to capture the nonlinear character of the relationship for these variables as illustrated in where the shapes of the nonlinear maps and their corresponding gradients (to assess the degree of nonlinearity) for the nonlinear models of T vap ,P vap and P elec respectively were plotted.

In the next section, the identied nonlinear maps for each variable of interest will serve to derive a nonlinear state space model for the ORC power plant under consideration. Then a linearized model of the system will be obtained to be used in the design of the constrained optimal controller as explained in Chapter 8. Consider the generic variable y(i) that could be one of the following variables of interest: regulated variable P elec or constrained variables (∆T vap , T vap , P vap , T eM ot, P eM ot and T f r ) that have been previously identied by the following nonlinear map: 128 Chapter 7. Nonlinear Dynamic Model Identication of an ORC power plant where:

y(i.τ s ) = y(i) = F y (Z y (i)) = Γ -1 y (Z T y (i).L y ) (
Z y (i) =      u y (i) u y (i -m) . . . u y (i -(N y -1).m)      ∈ R n y u .N y
The vector u y (i) used in the denition of the previous regressor can be one of the vectors dened in (7.21).

Since all the identied variables are based upon regressors of dierent buer lengths and dierent control vectors, one needs to dene a common regressor for all the variables of interest in order to construct a common state vector and state space equation (as shown in next paragraphs). Then, using the identied nonlinear maps, the relation with the common state vector will be established for each of the identied variables.

To construct the common regressor, consider the following regressor length: N = max{N ∆Tvap , N Tvap , N Pvap , N P elec , N T f r , N T eM ot , N P eM ot } (7.23) and the following input vector that contains all actuation and exogenous variables:

u(i) = f pump (i) Vhot (i) V aero rot (i) V circ rot (i) T aCh(i) T amb (i) T ∈ R nu (7.24)
Now, one can dene a common regressor for all the variables of interest as follow:

Z(i) =      u(i) u(i -m) . . . u(i -(N -1).m)      ∈ R nu.N (7.25)
If we consider the following common state vector (a part of the common regressor Z(i)) :

x(i) =      u(i -m) u(i -2.m) . . . u(i -(N -1).m)      ∈ R (N -1).nu (7.26)
then the state space equation, updated at a period τ m = mτ s , that is common for all the variables of interest will be given by:

x(i + m) =      u(i) u(i -m) . . . u(i -(N -2).m)      = A.x(i) + B.u(i) (7.27)
7.5. State Space representation of the ORC power plant identied dynamic models 129

where:

A =        O nu O nu • • • • • • O nu I nu O nu • • • • • • O nu O nu I nu O nu • • • O nu . . . . . . . . . . . . . . . O nu • • • O nu I nu O nu        ∈ R (N -1
).nu×(N -1).nu and

B =      I nu O nu . . . O nu      ∈ R (N -1).nu×nu
O nu ∈ R nu×nu and I nu ∈ R nu×nu denote respectively the zero and identity matrices.

By making the index change i = km, the new state space equation becomes:

x((k + 1)τ m ) = x(k + 1) = A.x(k) + B.u(k) (7.28)
where the state vector x(k) is given by:

x(k) =      u(k -1) u(k -2) . . . u(k -(N -1))      ∈ R (N -1).nu (7.29)
Similarly, one can rewrite the expressions of the common regressor Z(k) and the generic variable of interest y(k) in term of the new index k as follows:

Z(k) =      u(k) u(k -1) . . . u(k -(N -1))      ∈ R nu.N (7.30) y(k.τ m ) = y(k) = F y (Z y (k)) = Γ -1 y (Z T y (k).L y ) (7.31)
Note that since u(k) can be split into:

u(k) = u T con (k) u T unc (k) (7.32)
where

u con (k) = f pump (k) Vhot (k) V aero rot V circ rot T ∈ R n con u (7.33)
is the vector of controllable actuators and

u unc (k) = T aCh(k) T amb (k) T ∈ R n unc u (7.34)
130 Chapter 7. Nonlinear Dynamic Model Identication of an ORC power plant is the vector of measured uncontrolled (operational) variables, the common state space equation (7.28) can be rewritten as follows:

x(k + 1) = A.x(k) + B con .u con (k) + B unc .u unc (k) (7.35) where B = B con B unc such that B con ∈ R (N -1).nu×n con u and B unc ∈ R (N -1).nu×n unc u .

Now since for each variable of interest y(k) (that could be either T vap (k), P vap (k), T f r (k), P elec (k), P eM ot(k), T eM ot(k) or ∆T vap ) the corresponding input vector u y (k) (used in the regressor denition) can be related to the common control input using the following relation:

u y (k) = R y ed .u(k) (7.36)
where R y ed ∈ R n y u ×nu (the notation R ed comes from the fact that u y (k) is a "reduced" vector deduced from u(k)) is given for each variable of interest by the following set of matrices:

R Pvap ed = R Tvap ed = R ∆Tvap ed =     1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 1     , R P elec ed = R T f r ed = R T eM ot ed = R P eM ot ed = I 6
Then it is possible to relate each regressor Z y (k) for each variable of interest y(k) to the common regressor Z(k) by rst expressing Z y (k) in terms of the common input vector u(k) as follows:

Z y (k) =      u y (k) u y (k -1)
. . . (7.37) where RED y Z = I N y ⊗ R y ed ∈ R N y .n y u ×N y .nu and the symbol ⊗ denotes the Kronecker product (see Denition 0.1) Then making Z(k) appear as follow:

u y (k -(N y -1))      =      R y ed .u(k) R y ed .u(k -1) . . . R y ed .u(k -(N y -1))      = RED y Z .      u(k) u(k -1) . . . u(k -(N y -1))     
     u(k) u(k -1) . . . u(k -(N y -1))      = Π (nu,N ) 1→N y .              u(k) u(k -1) . . . u(k -(N y -1)) u(k -N y ) . . . u(k -(N -1))              = Π (nu,N ) 1→N y .Z(k) (
= I N y O N y ×(N -N y ) ⊗ I nu ∈ R N y .nu×N.nu
By combining (7.37) and (7.38), the relation between Z y (k) and Z(k) becomes:

Z y (k) = RED y Z .Π (nu,N ) 1→N y .Z(k) (7.39)
Knowing that: (see (7.30), (7.29) and (7.32))

Z(k) =   u con (k) u unc (k) x(x)  
and making the following matrix decomposition:

RED y Z .Π (nu,N ) 1→N y = M y ucon M y uunc M y x (7.40) where : M y ucon ∈ R N y .n y u ×n con u , M y uunc ∈ R N y .n y u ×n unc u and M y x ∈ R N y .n y u ×(N -1
).nu , one obtains the relation between each variable related regressor Z y (k) and the common state vector x(k), the common control vector u con (k) and the vector of measured uncontrolled signals u unc (k):

Z y (k) = M y x .x(k) + M y ucon .u con (k) + +M y uunc .u unc (k) (7.41)
Now by combining (7.31) and (7.41), one obtains the output equation of each variable of interest y(k) in terms of x(k), u con (k) and u unc (k) :

y(k) = Γ -1 y (L T y .M y x ).x(k) + (L T y .M y ucon ).u con (k) + (L T y .M y uunc ).u unc (k) (7.42)
Then, the state space model for each variable of interest is given by:

x(k + 1) = A.x(k) + B con .u con (k) + B unc .u unc (k) y(k) = Γ -1 y (L T y .M y x ).x(k) + (L T y .M y ucon ).u con (k) + (L T y .M y uunc ).u unc (k)
This model is obviously observable since the state vector is based upon the past values of the input vector according to (7.29).

ORC Model Linearization

According to Section 7.4.3, the identication of the maps for T f r , P eM ot, T eM ot and ∆T vap resulted in linear models unlike the identication of the maps for T vap , P vap and P elec where nonlinear models were needed to capture the nonlinear character of the relationship for these variables. However, according to the shape of the nonlinear maps for T vap , P vap and P elec (illustrated in Figures 7.9b-7.11b), it is possible to nd a straight line that approximates in a least squares sense these nonlinear map resulting in the end in a linear model for the ORC system that optimally approximates the nonlinear model. The linear model will then be used in the controller design and the model errors will be corrected through slow integrators.

Let σ be dened as follows (according to (7.42)):

σ = (L T y .M y x ).x(k) + (L T y .M y ucon ).u con (k) + (L T y .M y uunc ).u unc (k) (7.44)
Hence, by approximating in a least squares sense the nonlinear map Γ -1 y (σ), given by (7.42) for each variable of interest y(k), as follows:

y(k) = Γ -1 y (σ) = a y .σ + b y (7.45)
and using (7.44), one obtains the nal output equation for each variable of interest y(k):

y(k) = C y .x(k) + D y con .u con (k) + D y unc .u unc (k) + E y (7.46)
Where: C y = a y .L T y .M y x , D y con = a y .L T y .M y ucon , D y unc = a y .L T y .M y uunc and E y = b y .

Figure 7.16 gives the linearization results (equation (7.45)) for the nonlinear maps Γ -1 Tvap (.), Γ -1 Pvap (.) and Γ -1 P elec (.) to assess how close the linear model approximates the nonlinear one for these variables. However, according to the discussion of Section 7.4.3, Γ -1 T f r (.), Γ -1 P eM ot (.), Γ -1 T eM ot (.) and Γ -1 ∆Tvap (.) were found to be linear, so it is not worthwhile to show their corresponding linearization results since there is an exact matching between Γ -1 (.) for these variables and their corresponding linear approximations.

Filtering control variables

Considering the approximated linear state space model for each variable of interest y(k) that was dened by equations (7.35) and (7.46), it is interesting to add a rst order lter on the vector of controlled actuator namely u con (k). Indeed adding this lter will induce the following benets:

• Eliminates the direct feedthrough matrix D y con in the output equation for each variable of interest y(k).

• Avoids the need for introducing rate constraints on control signals (since the rate of variation for the control signals will be dened by the time constants of the lter). This reduces the complexity of the controller.

Hence, the state space model after ltering the control vector becomes as follow:

x(k + 1) = A.x(k) + B con .u con (k) + B unc .u unc (k) (7.47a) u con (k + 1) = A f .u con (k) + B f .u f (k) (7.47b) y(k) = C y .x(k) + D y con .u con (k) + D y unc .u unc (k) + E y (7.47c)
Figure 7.16: Linearization of the nonlinear maps Γ -1 Tvap (.), Γ -1 Pvap (.) and Γ -1 P elec (.).

With: B f = I n con u -A f , where A f is a diagonal ltering matrix whose diagonal elements are set according to the desired lter time constants corresponding to each element of the control vector u con (k). The new control vector becomes u f (k) which has the same dimension as the ltered control vector u con (k).

By extending the state vector with the vector of controlled variables u con (k) as follow:

X(k) = x(k) u con (k) ∈ R nx where n x = n x + n con u (7.48)
the new state space equations for each variable of interest become:

X(k + 1) = Ā. X(k) + Bcon .u f (k) + Bunc .u unc (k) (7.49a) y(k) = Cy . X(k) + D y unc .u unc (k) + E y (7.49b) Where: Ā = A B con O A f , Bcon = O B f
and Cy = C y D y con

Conclusion

In this chapter a recently developed identication framework was applied to derive nonlinear models for the ORC power plant. The nonlinear models so obtained were validated using experimental data gathered on a real prototype. Then a global state space model of the ORC 134 Chapter 7. Nonlinear Dynamic Model Identication of an ORC power plant power plant was derived from the identied models and will be used in the next chapter to derive an optimal control strategy to reach the system objectives while respecting the physical constraints.

Chapter This chapter addresses the problem of designing a constrained optimal control strategy for the ORC power plant based on the state space model (7.49) given for each variable of interest. The presented controller optimally makes the ORC system track some reference electrical power namely P ref elec while respecting the system constraints in term of temperatures and pressures for some key thermodynamic variables. Prior to the control design step, an analysis of the admissible desired steady states is performed to derive feedforward terms that enhance the controller performances.

Study of the admissible steady states

Now that we have dened the state space equations (7.49) and before developing the constrained optimal control strategy for the ORC power plant, it is necessary to study the admissible desired steady states i.e. the steady states that are compatible with some desired reference electrical power P ref elec leading to the respect of the system constraints set in Section 7.2 while penalizing the control vector in the expression of the cost function.

At equilibrium, (7.49a) becomes:

Xd (k + 1) = Ā. Xd (k) + Bcon .u f d (k) + Bunc .u unc (k) = Xd (k) (8.1)
where the superscript "d" stands for "desired".

136

Chapter 8. Constrained Optimal Control of the ORC based power plant Now, consider some vector of regulated variables y r (k) (that obviously contains at least P elec (k) as one of its components) that will be dened later in this section. At equilibrium and according to (7.49b), the resulting desired vector of regulated variables will be given by:

y d (k) = Cyr . Xd (k) + D yr unc .u unc (k) + E yr (8.2)
where Cyr , D yr unc and E yr will be formed by concatenating line by line the matrices Cy , D y unc and E y respectively corresponding to each output variable y(k) which is part of the vector y r (k).

Combining (8.1) and (8.2) in matrix form results in:

( Ā -I nx ) Bcon Cyr O . Xd (k) u f d (k) + Bunc D unc .u unc (k) + O E yr = O I .y d (8.3)
Note from (8.3) that the steady state equation is time dependent since the vector of the measured uncontrolled signals u unc (k), which is time dependent with a slow dynamics compared to the other variables of interest (because it depends on T aCh and T amb that vary slowly), is involved in this equation.

By denoting:

A c = ( Ā -I nx ) Bcon Cyr O , B unc c = Bunc D unc , E c = O E yr , B c = O I and Z d = Xd (k) u f d (k)
the steady state equation (8.3) becomes:

A c .Z d + B unc c .u unc (k) + E c = B c .y d (8.4)
Furthermore, as mentioned in the beginning of this section, the desired steady states must lead to the respect of the ORC system constraints. It is then necessary to express the vector of constrained variables at equilibrium in term of Z d . This begins by expressing the output equation of each variable of interest y(k) as follow:

y(k) = Cy . X(k) + D y unc .u unc (k) + E y = C y c . X(k) u f (k) + D y unc .u unc (k) + E y (
Y min c,Z d -E Z d -D Z d .u unc (k) ≤ C Z d c .Z d ≤ Y max c,Z d -E Z d -D Z d .u unc (k) (8.8)
where

Y min c,Z d = y min c Z min d ; Y max c,Z d =      C Pvap c C Tvap c C ∆Tvap c C T f r c C T eM ot c C P eM ot c I (nx+n con u )                 
Note that equation (8.8) can be rewritten as:

C Z d T .Z d ≤ b Z d T 1 + b Z d T 2 .u unc (k) (8.9)
Where:

C Z d T = C Z d c -C Z d c , b Z d T 1 = Y max c,Z d -E Z d -Y min c,Z d + E Z d and b Z d T 2 = -D Z d D Z d
Finally, the extended steady state vector Xd and the steady state control u f d are given by the steady desired vector Z d = Xd u f d which is the solution of the following Constrained Linear Programming (LP) problem:

min Z d f T .Z d (k) (8.10)
under the following equality and inequality constraints:

C Z d T .Z d (k) ≤ b Z d T 1 + b Z d T 2 .u unc (k) A c .Z d (k) = B c .y d -B unc c .u unc (k) -E c
where f = 1 (nx+n con u ) is a weighting vector. The LP problem (8.10) that enables to derive the admissible extended steady vector Xd and its corresponding steady control vector u f d obviously needs to be fed with the desired steady vector y d whose only a priori known component is the desired electrical power P ref elec . Since this LP problem might be quite heavy to solve online due to the important size of its associated matrices (for the considered application Z d ∈ R 62 and C Z d T ∈ R 130×62 ) and since we will have additionally to solve online a QP problem corresponding to our constrained optimal controller (as it will be explaines in Section 8.2), we proposed another alternative based on an oine computed cartography to determine the desired admissible steady state vectors as it is explained hereafter.

First of all, note that if the A c matrix is invertible (full ranked square matrix) and one can determine an admissible (leading to the respect of the ORC system constraints) desired steady vector y d , the desired admissible steady vector Xd and its corresponding steady control vector u f d will be derived using the following equation: (see (8.4)):

Z d (k) = Xd u f d = A -1 c .(B c .y d -B unc c .u unc (k) -E c ) (8.11)
In order to make the A c matrix invertible, the size of y d (and therefore the size of y r (k)) must be equal to n con u = 4. Then, we have to choose in addition of P elec three more critical variables to add to the vector y r (k). In this case, it would be judicious according to the degree of importance of the constrained variables to choose ∆T vap , T vap and P eM ot as additional variables in y r (k). This leads to the denition of the following matrices:

y r (k) =        P elec (k) ∆T vap (k) T vap (k) P eM ot(k)        , y d =         P ref elec ∆T d vap T d vap P eM ot d         and Cyr =        CP elec C∆Tvap CTvap CP eM ot       
Then, a cartography carto ORC (T aCh vect , T vect amb , P vect elec , ∆T vect vap , T vect vap , P eM ot vect ) is built ofine. This cartography is associated to a set of vectors of interpolation points corresponding to the operational uncontrolled variables namely T aCh and T amb and the set of regulated variables namely P vect elec , ∆T vect vap , T vect vap and P eM ot vect . The bounds of the vectors of interpolation points have to be chosen such that they entirely dene all the admissible domain for the operational (exogenous) and regulated variables.

By scanning all the possible values of the vectors of interpolation points, the LP problem (8.10) is solved oine for each of the values T aCh(i), T amb (j), P vect elec (k), ∆T vect vap (l), T vect vap (m) and P eM ot vect (n) (where i, j, k, l, m and n are the indices of the vectors of interpolation points corresponding to some scanned value). If the LP problem has a solution, then the value 1 is put in the cartography at the indice (i, j, k, l, m, n) otherwise a 0 is put instead.

Then, by measuring the operational variables T aCh and T amb and knowing the reference electrical power P ref elec (power consumed by the o-grid loads), one has to extract in real time from carto ORC (.) the admissible reference vector y d according to the diagram depicted in Figure 8.1. Note in this diagram that a minimum admissible value of ∆T vect vap (l), T vect vap (m) and P eM ot vect (n) is chosen in order to maximize the motor eciency.

Finally, knowing y d , the desired admissible steady vector Xd and its corresponding steady control vector u f d are derived using (8.11). These matrices are determined by tuning giving more importance (and hence more weight) to the regulation respectively of the superheating ∆T vap , T vap , P eM ot and nally P elec in order to let the constrained variables, that intervene in the denition of y r , converge rapidly towards their calculated references. The fact there is a great dispersion in the values of the matrices Q and R is due to a normalization eect.

Note that since some of the constraints in (8.15) are hard, constraints softening methods are used to soften the remaining constraints to guarantee the feasibility of the optimization problem.

The control architecture of the ORC power plant is depicted in Figure 8.2.

Simulation Results

In order to assess the performance of the proposed control architecture (Figure 8.2), a simulation scenario is proposed in which we vary in time the operating conditions (with time constants close to the reality) in term of the hot source temperature prole T aCh and the ambient temperature T amb as depicted in Figure 8.3a. Figure 8.3b gives the tracking results for a reference electrical power P ref elec prole consumed by the o-grid loads. One can appreciate from Figure 8.3b the performances in term of settling time and disturbance (measured uncontrolled signals) rejection of the proposed controller since the latter one makes the ORC system converge rapidly toward the desired reference power with time responses reduced by more than half compared to the time responses obtained by EXOES (one of the MICROSOL project partners) using their PI-based control methodology (indeed, EXOES announced a time response of 30s to respond to a 1kW electrical power demand). The time responses obtained with the proposed control methodology suggest that a signicant reduction in the battery bank size can be performed. Figure 8.4 gives the time prole of the controlled variables (output of the controller) and their corresponding maximum and minimum admissible values that are perfectly respected. The time proles of the constrained thermodynamic variables of the ORC system are depicted in Figure 8.5 together with their maximum and minimum admissible values (according to Section 7.2). One can observe that these constraints are respected (and sometimes active) even though at some time instants some of the constraints slightly overshoot their limitations due to the application of the constraint softening techniques that enable the MPC solver convergence.

Conclusion

In this chapter, the global state space model of the ORC power plant was used to design a constrained optimal control strategy to reach the system objectives while respecting its constraints. The predictive controller is real time implementable and shows very good performances in simulation. Future works will focus on implementing this controller on a real time target and assess experimentally the global performance of the proposed solution.

General Conclusion and Perspectives

This thesis, prepared within the Gipsa-lab laboratory in partnership with Schneider Electric, is concerned with the modeling, simulation and control of two kinds of solar thermodynamic micro power plants that were considered as candidate power generation systems within the framework of the MICROSOL project led by Schneider Electric.

In the case of the rst solar thermodynamic micro power plant, which is based on a Stirling thermodynamic engine linked to a dedicated energy conversion system developed within Schneider Electric laboratories, a nonlinear static model of the Stirling engine has been derived and its unknown parameters identied using experimental data gathered on the prototype. This model has been used to simulate the behavior of the Stirling engine with its power conversion system and to develop a cartography-based torque controller for the engine that maximizes the net thermoelectric eciency of the system. Then, a state space model of the energy conversion system associated to the Stirling engine has been derived and used to develop a rst control strategy, for the hybrid Stirling engine/supercapacitor power plant, that was validated experimentally rst on a dedicated electrical motor testbed and then on the real Stirling based power system. This control strategy in addition of necessitating few design parameters endows the closed-loop system with near-to-optimal behavior in terms of settling time thus reducing the size of the supercapacitor needed as energy buer to palliate the load power demand fast transients. Even if this developed control strategy was ecient and oers the expected performances, it does not handle explicitly the system constraints. For this reason, a second control strategy for the energy conversion system, that holds for a family of power circuits, has been developed. In addition of handling explicitly the saturation on control variables as well as on the component of the state vector, it also endows the closed-loop system with near to optimal performances in term of settling time.

In a second part, the problem of controlling optimally the energy of the second solar thermodynamic micro power plant, which is based on Organic Rankine Cycle (ORC) machines, has been addressed and solved. Prior to the control design step, one of the hurdles consisting in deriving nonlinear models for this engine has been solved using a recently developed nonlinear identication framework. Then based on linearized version of the identied model for the ORC machine, a model predictive controller has been developed that drives the system in a fast manner during load power demand transients, while respecting the physical system constraints, which enables the size reduction (at least by half) of the energy buer consisting in the battery banks associated to the power plant. Preliminary tests on implementing the developed control strategy on a real time controller have proved the real time compatibility of the control algorithm.

The experimental data used to identify nonlinear models of the ORC power plant under consideration represent some specic operating range corresponding to tuning and calibration tests performed on the prototype. When exploring other operating ranges of the plant, one may obtain nonlinear models with a higher degree of nonlinearity. In this situation, either we 148 Conclusion decide to apply the same control methodology, developed during this thesis, using a linearized model representing all the operating ranges and accepting some degradation of the closed loop system performances and less margins on the system constraints (however this approach has the advantage of resulting in an MPC controller whose implementablity on a real time target has been validated), or we decide to apply constrained nonlinear controllers (like NMPC controllers) that has to be validated and that will necessitate more computation resources to be implementable. This constitutes a very interesting research perspective to be investigated.

As other perspectives following the work done during this thesis, the experimental validation, on the real prototype, of the developed control strategy for the ORC based power generation system rst using a fast real time target to prove the concept of the developed methodology and then implementing the control algorithms on the PLC (Programmable Logic Controller) that currently controls the ORC power plant prototype. Another interesting perspective for the present work is the investigation of ORC based power systems functioning with a variable speed turbine. Indeed, by using the speed reference of the turbine as an additional control variable, we will be able to design control strategies that make the ORC based power plant tracks some electrical reference power in a faster manner,compared to the constant speed case, thus reducing further the size of the battery bank associated to the power plant.

Since the proposed methodology for controlling the ORC based power system is quite generic and can be applied to several kinds of thermodynamic engines, it would be also interesting to investigate the applicability of this methodology to other thermal power generation systems like WHR (Waste Heat Recovery) systems which is currently the focus of many researches. L'objectif de cette thèse étant l'élaboration de lois de commande performantes pour la gestion de l'énergie électrique de ces deux types de micro centrales solaire thermodynamiques. Etant donné que la puissance fournie par la machine Stirling est égale à son couple moteur fois la vitesse de rotation de l'arbre moteur et que ces deux grandeurs évoluent dans deux échelles de temps diérentes, la synthèse des deux contrôleurs pour ces 2 grandeurs d'une manière indépendante est justiée. Ainsi un premier contrôleur du couple du moteur Stirling a été développé en se basant sur un modèle statique de la machine validé expérimentalement.

Pour réguler la vitesse de rotation de l'arbre moteur, on agit sur le système de conversion de l'énergie. Ainsi deux lois de commande ont été développées à cet eet. La première loi de commande qui a été développée et qui a été validée expérimentalement gère explicitement les saturations sur les variables de contrôle du système de conversion d'énergie et dote le système en boucle fermé de performances quasi optimales en terme de temps de réponse. Cependant, cette première loi de commande ne gère pas explicitement les contraintes sur tous les états du système. Pour pallier à cette limitation, une deuxième loi de commande a été développée en se basant sur le principe des modes glissants. Pour modéliser en régime statique la machine Stirling, nous avons utilisé un modèle issu de la littérature qui se base sur l'analyse de Schmidt [START_REF] Urieli | Stirling Cycle Engine Analysis[END_REF] qui établit un bilan de puissance parfait dans la machine et qui dépend des paramètres dimensionnant de Ainsi un contrôleur par backstepping (voir [Khalil, 1996]) sert à réguler x 5 autour de x st 5 d'une manière assez raide en utilisant les équations du sous-système 2. Pour réguler la tension du supercondensateur x 7 autour de sa valeur stationnaire x st 7 , on impose à x 6 en régime stationnaire la valeur suivante:

x st 6 = k 6 . tanh(β.(x 7x st 7 )) (A.4) où k 6 and β sont des paramètres de réglage.

On aboutit ainsi à la dénition de la valeur de référence x ref 4 ci-dessous qui devra être suivie en utilisant les équations du sous-système 1:

x ref st issue du ltrage de x 4 (pour lisser certains eets d'undershoot et d'overshoot) et de l'ajout d'un intégrateur lent (pour corriger certaines erreurs du modèle) et on résout les équations correspondantes à ż = 0 pour déterminer le vecteur z st qui sera régulé en utilisant un contrôleur prédictif à un pas appliqué aux équations du sous-système 1 (Figure A.8). Il est à noter, d'après la Figure A.8, l'ajout d'un observateur de Luenberger (voir [Besançon, 2007]) qui sert à estimer les erreurs d'identication du couple de la machine Stirling intervenant dans le paramètre a 2 de l'équation d'état (A.2a).

Avant de valider cette première loi de commande sur le prototype à base de moteur Stirling (Figure A.4), nous l'avons d'abord validé sur un banc moteur électrique illustré à la Figure Dans le cas de la machine ORC considérée dans cette thèse, la vitesse de rotation des turbines est xe. Ainsi, pour modier la puissance fournie par la machine ORC, on ne peut agir que sur son couple moteur produit à travers la partie thermo-hydraulique. Là encore, l'objectif de contrôle pour ce système est de faire passer la machine ORC d'un point de fonctionnement à un autre le plus rapidement possible tout en respectant ses contraintes physiques an de diminuer la taille des batteries (source d'énergie tampon) associées à la machine ORC et pour cela des lois de commande performantes doivent être développées. Pour ce faire, nous avons au préalable identié des modèles dynamiques de la machine ORC en utilisant des données expérimentales recueillis sur un prototype fonctionnel illustré à la In order to get the average model of the DC/DC full bridge converter that will be used during the control design, we apply the following moving average formula over a cutting period T sw to equations B.4:

x(t) =
dV bus dt = K bus1 .(V ref bus -V bus )
where K bus1 is some design parameter. This leads to:

I ref Lf b = -I Lbb + I inv + C tot .K bus1 .(V ref bus -V bus ) (C.2)
To overcome some model errors, measurements errors and undesirable sharp transients during starting procedure, I ref Lf b is ltered using a rst order lter and a corrective term is added to

1. 1

 1 BipBop program axes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 2.1 Representation of the needed Saharian surfaces to produce the currently consumed energy in MENA (Middle East and North Africa), EU (European Union) and the World. (provided by the German Aerospace Center (DLR), 2005.) . . . 8 2.2 Schematic of a solar thermodynamic micro power plant . . . . . . . . . . . . . . 10 2.3 SCHOTT PTR 70 solar absorber . . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.4 Parabolic trough mirrors constituting the solar collector eld of the ORC based power generation system . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.5 Stirling engine provided by "Cool Energy" . . . . . . . . . . . . . . . . . . . . . 12 2.6 ORC prototype provided by "EXOES" . . . . . . . . . . . . . . . . . . . . . . . 13 3.1 Diagram of an alpha type Stirling engine (http://www.ohio.edu/mechanical/ stirling/) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 3.2 Diagram of a beta type Stirling engine (http://www.ohio.edu/mechanical/ stirling/) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 3.3 Diagram of a gamma type Stirling engine (http://www.ohio.edu/mechanical/ stirling/) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 3.4 Pressure-Volume (P-V) diagram of the Stirling Cycle and pistons position for an alpha-type Stirling engine during a complete cycle . . . . . . . . . . . . . . . 23 3.5 temperature distribution for an alpha type Stirling engine in the ideal isothermal model (http://www.ohio.edu/mechanical/stirling/). . . . . . . . . . . 25 3.6 Sinusoidal volume variation (Schmidt analysis) for an alpha type Stirling engine (http://www.ohio.edu/mechanical/stirling/). . . . . . . . . . . . . . . . . . 28 3.7 Overview of the Stirling based power system . . . . . . . . . . . . . . . . . . . . 34 3.8 Zoom on the Stirling engine and the cooling loop . . . . . . . . . . . . . . . . . 35 3.9 Identication procedure for the Stirling engine unknown parameters . . . . . . 36 3.10 Comparaison between the measured and the theoretical output electrical power, absorbed and rejected heat powers . . . . . . . . . . . . . . . . . . . . . . . . between the measured and theoretical output hot and cold temperatures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 3.12 Comparison between the tted Stirling engine torque equation and the torque derived from the experimental data. each point in the gure corresponds to an operating condition dened by N rot , ṁh , ṁc , T h in and T c in . . . . . . . . . . . 38 3.13 Stirling engine optimal torque controller . . . . . . . . . . . . . . . . . . . . . . 39 3.14 Isolines of the optimal hot side mass ow rate ṁref h for P ref elec =100 W and P ref elec =200. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 3.15 Isolines of the optimal cold side mass ow rate ṁref c for P ref elec =300 W and P ref elec =400. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 3.16 Optimal hot and cold mass ow rates references for the following operating conditions : T h in = 170• C and T c in = 20 • C. . . . . . . . . . . . . . . . . . . . 41 4.1 Overview of the Stirling based thermodynamic power plant and its associated power conversion system . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 4.2 Behn-Eshenbourg simplied model for the PMSG . . . . . . . . . . . . . . . . . 45 4.3 Full wave diode bridge. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 4.4 DC/DC Full Bridge converter. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 4.5 bidirectional DC/DC converter. . . . . . . . . . . . . . . . . . . . . . . . . . . . 49 5.1 Evolution of the stationary states of subsystem 1 as function of u st 1 . . . . . . . . 57 5.2 Stability indicator of A(u st 1 ) for u st 1 ∈ [0, 1]. . . . . . . . . . . . . . . . . . . . . 58 5.3 Open loop trajectories for the components of the vector z(t) corresponding to u st 1 = 0.2, u st 1 = 0.4 and u st 1 = 0.9. . . . . . . . . . . . . . . . . . . . . . . . . . . 58 5.4 Open loop trajectories in the phase plane (x 3 ,x 4 ) corresponding to u st 1 = 0.2, u st 1 = 0.4 and u st 1 = 0.9. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 5.5 Global control architecture for the complete system. . . . . . . . . . . . . . . . 60 5.6 Load power demand. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61 5.7 regulated variables. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 5.8 rotor mechanical speed and currents injected in the DC bus by the bidirectional DC/DC converter and the Full bridge converter. . . . . . . . . . . . . . . . . . 63 List of Figures xi 5.9 control variables. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.10 Schematic of the electrical motor testbed together with its Power electronics. . 5.11 Electrical motor testbed used for experimental validation and performance assessment of the proposed control strategy. . . . . . . . . . . . . . . . . . . . . . 5.12 DC/DC Full Bridge converter. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.13 Load power demand P L at the input of the inverter for the electrical motor testbed. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.14 Bus voltage V bus and its reference value V ref bus for the electrical motor testbed. . 5.15 Supercapacitor voltage V sc and its reference value V ref sc for the electrical motor testbed. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.16 Output current of the DC/DC Full Bridge converter I Lf b and its reference value I ref Lf b for the electrical motor testbed. . . . . . . . . . . . . . . . . . . . . . . . . 5.17 Rotational speed of the engine Ω, the rectied voltage V red and the rectied current I red with their respective reference values Ω ref , V ref red and I ref red for the electrical motor testbed. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.18 Duty ratios of the DC/DC Full Bridge converter (upper gure) and bidirectional DC/DC converter (lower gure) for the electrical motor testbed. . . . . . . . . 5.19 Electrical cabinets containing the energy conversion system . . . . . . . . . . . 5.20 Overview of the rst electrical cabinet . . . . . . . . . . . . . . . . . . . . . . . 5.21 Overview of the second electrical cabinet . . . . . . . . . . . . . . . . . . . . . . 5.22 Load power demand P L at the input of the inverter for the Stirling based power plant. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.23 Bus voltage V bus and its reference value V ref bus for the Stirling based power plant. 5.24 Supercapacitor voltage V sc and its reference value V ref sc for the Stirling based power plant. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.25 Output current of the DC/DC Full Bridge converter I Lf b and its reference value I ref Lf b for the Stirling based power plant. . . . . . . . . . . . . . . . . . . . . . . . 5.26 Rotational speed of the engine Ω, the rectied voltage V red and the rectied current I red with their respective reference values Ω ref , V ref red and I ref red for the Stirling based power plant. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.27 Duty ratios of the DC/DC Full Bridge converter (upper gure) and bidirectional DC/DC converter (lower gure) for the Stirling based power plant. . . . . . . .

2= 5

 5 is used. . . . 102 List of Figures xiii 6.12 Performance of the overall closed-loop system under sinusoidal power load demand. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7.1 Picture of the ORC power plant prototype . . . . . . . . . . . . . . . . . . . . . 7.2 Schematic of the Organic Rankine Cycle Power Plant. . . . . . . . . . . . . . . 7.3 Temperature-Entropy (T-S) diagram of the Organic Rankine Cycle Power Plant.112 7.4 comparaison between the tted polynomial for T sat (P vap ) with the data provided by REFPROP R . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7.5 Moving Boundary paradigm applied to a condenser [Shah et al., 2003]. . . . . . 7.6 Allure of the maps B j (.) for β = 0.5 and n m = 10. . . . . . . . . . . . . . . . . 7.7 Control and exogenous variables proles during Experiment 1. . . . . . . . . . . 7.8 Control and exogenous variables proles during Experiment 2. . . . . . . . . . . 7.9 Identication results for T vap . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7.10 Identication results for P vap . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7.11 Identication results for P elec . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7.12 Identication results for T f r . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7.13 Identication results for P eM ot. . . . . . . . . . . . . . . . . . . . . . . . . . . 7.14 Identication results for T eM ot. . . . . . . . . . . . . . . . . . . . . . . . . . . 7.15 Identication results for ∆T vap . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7.16 Linearization of the nonlinear maps Γ -1 Tvap (.), Γ -1 Pvap (.) and Γ -1 P elec (.). . . . . . . . 8.1 Principle diagram for extracting y d from carto ORC (T aCH vect , T vect amb , P vect elec , ∆T vect vap , T vect vap , P eM ot vect ) based on the knowlege of P ref elec . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.2 Diagram of the ORC control System. . . . . . . . . . . . . . . . . . . . . . . . . 8.3 Tracking results and disturbance rejection for the ORC turbine output electrical power P elec . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.4 ORC system actuators: f pump , Vhot , V aero rot and V circ rot . The red dashed lines indicate the minimum and maximum admissible values. . . . . . . . . . . . . . 8.5 ORC system constrained variables. The red dashed lines indicate the minimum and maximum admissible values. . . . . . . . . . . . . . . . . . . . . . . . . . . xiv List of Figures A.1 Schéma global d'une micro centrale solaire thermodynamique . . . . . . . . . . 150 A.2 Schéma global d'une micro centrale solaire thermodynamique à base d'un moteur Stirling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 151 A.3 diagramme Pression-Volume (P-V) du cycle de Stirling et position des pistons pour un moteur Stirling de type alpha pendant un cycle complet. . . . . . . . . 152 A.4 Vue d'ensemble du système de production d'électricité à base de machine Stirling.153 A.5 Procédure d'identication des paramètres inconnus du modèle de la machine Stirling. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 154 A.6 Comparaison entre les valeurs mesurées et théoriques de: la puissance électrique produite, de la puissance thermiques absorbée et de la puissance thermique rejetée.154 A.7 Comparaison entre les valeurs mesurées et théoriques des températures chaude et froide en sortie de la machine Stirling. . . . . . . . . . . . . . . . . . . . . . . 154 A.8 Architecture de la première loi de commande pour le système de conversion d'énergie associé à la machine Stirling. . . . . . . . . . . . . . . . . . . . . . . . 157 A.9 Banc moteur électrique utilisé pour la validation expérimentale de la première loi de commande pour le système de conversion d'énergie associé au moteur Stirling. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 158 A.10 Résultats expérimentaux validant la validation de la première loi de commande pour le système de conversion d'énergie associé au moteur Stirling. . . . . . . . 159 A.11 Disposition typique des champs de vecteurs pour le système dynamique (A.2c)-(A.2d) quand x 2 = x st 2 (x * 3 ) et x 5 = x r 5 sont supposés constants. . . . . . . . . . 160 A.12 Contrôle du système (x 3 , x 4 ) sous x 5 = x d 5 après un changement de +20% pour la valeur désirée x d 4 . Il est à remarquer que les limitations sur x 3 ralentissent le temps de réponse du système an de garantir le respect des contraintes sur x 2 . 161A.13 Schéma global d'une micro centrale solaire thermodynamique à base d'une machine ORC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 162 A.14 Schéma du système de production d'électricité à base d'une machine à cycle de Rankine organique (ORC). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 163 A.15 diagramme Température-Entropie (T-S) du système de production d'électricité à base d'une machine à cycle de Rankine organique (ORC). . . . . . . . . . . . 163 A.16 Prototype du système de production d'électricité à base de machine ORC . . . 164 A.17 Résultats d'identication pour P elec . . . . . . . . . . . . . . . . . . . . . . . . . 167 List of Figures xv A.18 Architecture de contrôle du système ORC. . . . . . . . . . . . . . . . . . . . . . A.19 Résultats de la régulation et du rejet de perturbations pour la puissance électrique P elec en sortie de la turbine du système ORC. . . . . . . . . . . . . . . . A.20 Actionneurs du système ORC: f pump , Vhot , V aero rot and V circ rot . Les lignes en pointillés rouge représentent les valeurs minimum et maximum admissibles. . . A.21 Variables sous contrainte du système ORC. Les lignes en pointillés rouge représentent les valeurs minimum et maximum admissibles. . . . . . . . . . . . . . . . . B.1 temporal diagram of the DC/DC Full Bridge. . . . . . . . . . . . . . . . . . . . B.2 Operational modes of the bidirectional DC/DC converter according to the current ow . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.3 Output current I Lbb of the bidirectional DC/DC converter in the complementary control mode. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . C.1 Overview of the Stirling based power plant and its associated power conversion system including the ATV32 variable speed drive . . . . . . . . . . . . . . . . . C.2 Starting procedure owchart . . . . . . . . . . . . . . . . . . . . . . . . . . . . . D.1 Illustration of the proof of Lemma 6.1: Under the conditions of the Lemma, the position of the pic of y = x 2 is almost independent of x 3 . . . . . . . . . . . List of Tables 7.1 The identied parameters used to dene the nonlinear maps for ∆T vap , T vap , P vap , P elec , T f r , P eM ot and T eM ot . . . . . . . . . . . . . . . . . . . . . . The transpose of a vector or matrix v v any norm of the vector v v P weighted norm of the vector v by the matrix P |v| element wise absolute value of v A ⊗ B Kronecker product of matrices A and B (see Denition 0x→y f (x) the limit of f of x, as x approaches y min(.) minimum value max(.) maximum value arg min x f (x) argument x that minimizes the function f (x) tanh(x) Hyperbolic Tangent of x sin(.),cos(.) sine and cosine functions arctan(.) arc tangent function I n Identity matrix of dimension n xix xx

  Figure 1.1: BipBop program axes

Figure 2

 2 Figure 2.1: Representation of the needed Saharian surfaces to produce the currently consumed energy in MENA (Middle East and North Africa), EU (European Union) and the World. (provided by the German Aerospace Center (DLR), 2005.)

Figure 2 . 2 :

 22 Figure 2.2: Schematic of a solar thermodynamic micro power plant

Figure 2

 2 Figure 2.5: Stirling engine provided by "Cool Energy"

Figure 2

 2 Figure 2.6: ORC prototype provided by "EXOES"

20

  Chapter 3. Nonlinear Static Model of the Stirling Engine

Figure 3 . 1 :

 31 Figure 3.1: Diagram of an alpha type Stirling engine (http://www.ohio.edu/mechanical/ stirling/)

Figure 3

 3 Figure 3.2: Diagram of a beta type Stirling engine (http://www.ohio.edu/mechanical/ stirling/)

Figure 3 .

 3 Figure 3.4 depicts the theoretical (indicated) Pressure-Volume diagram of the Stirling cycle and the disposition of the cylinders during each step of the cycle(here we considered an alpha type conguration that corresponds to the prototype we have in our disposal).

Figure 3

 3 Figure 3.3: Diagram of a gamma type Stirling engine (http://www.ohio.edu/mechanical/ stirling/)

Figure 3

 3 Figure 3.5: temperature distribution for an alpha type Stirling engine in the ideal isothermal model (http://www.ohio.edu/mechanical/stirling/).

Figure 3

 3 Figure 3.6: Sinusoidal volume variation (Schmidt analysis) for an alpha type Stirling engine (http://www.ohio.edu/mechanical/stirling/).

Figure 3 . 7 :

 37 Figure 3.7: Overview of the Stirling based power system

  Figure 3.11 compares the measured and theoretical output hot and cold temperatures. From Figure 3.10 and Figure 3.11 one can appreciate the good model quality since most of the points in the plane (measured variables-output model variables) lie near the identity line.

Figure 3 . 9 :

 39 Figure 3.9: Identication procedure for the Stirling engine unknown parameters

Figures 3 .

 3 Figures 3.14 and 3.15 give the picture of the o-line computed cartographies by depicting the isolines of ṁref h and ṁref c for a xed value of P ref elec . Figure 3.16 gives an example of the resulting optimal mass ow rates for T h in = 170 • C and T c in = 20 • C.

Figure 3 .

 3 Figure 3.14: Isolines of the optimal hot side mass ow rate ṁref h for P ref elec =100 W and P ref elec =200.

Figure 3 .

 3 Figure 3.15: Isolines of the optimal cold side mass ow rate ṁref c for P ref elec =300 W and P ref elec =400.

Figure 3 .

 3 Figure 3.16: Optimal hot and cold mass ow rates references for the following operating conditions : T h in = 170 • C and T c in = 20 • C.

Figure 4

 4 Figure 4.1 shows an overview of the solar thermodynamic micro power plant that consists of an alpha type Stirling engine linked mechanically to a Permanent Magnet Synchronous

Figure 4

 4 Figure 4.1: Overview of the Stirling based thermodynamic power plant and its associated power conversion system

  the 3-phase output voltage of the PMSG and [I] 3 = I a I b I c T is the 3-phase line current.

Figure 4

 4 Figure 4.3: Full wave diode bridge.

Figure 4

 4 Figure 4.4: DC/DC Full Bridge converter.

Figure 4

 4 Figure 4.5: bidirectional DC/DC converter.

  7b) where C bb and C f b are the values of the output capacitors of the bidirectional DC/DC converter and DC/DC full bridge converter respectively. C inv is the value of the input capacitor of the inverter and C tot is the equivalent capacitance of the three capacitors in parallel. I inv denotes the input current of the inverter. 50 Chapter 4. Dynamic Model of the Energy Conversion System associated with the Stirling Engine

  Figure 5.1: Evolution of the stationary states of subsystem 1 as function of u st 1 .

Figure 5 . 2 :

 52 Figure 5.2: Stability indicator of A(u st 1 ) for u st 1 ∈ [0, 1].

Figure 5 . 3 :

 53 Figure 5.3: Open loop trajectories for the components of the vector z(t) corresponding to u st 1 = 0.2, u st 1 = 0.4 and u st 1 = 0.9.

Figure 5 . 4 :

 54 Figure 5.4: Open loop trajectories in the phase plane (x 3 ,x 4 ) corresponding to u st 1 = 0.2, u st 1 = 0.4 and u st 1 = 0.9.

60

  Chapter 5. Development of a rst control strategy for the Energy ConversionSystem associated to the O-Grid Stirling based power Plant.

Figure 5 . 5 :

 55 Figure 5.5: Global control architecture for the complete system.

Figure 5

 5 Figure 5.6: Load power demand.

  Figure 5.8b gives the online computed desired stationary value of the DC/DC full bridge converter output current namely x st 4 = I ref Lf b

  Figure 5.7: regulated variables.

Figure 5

 5 Figure 5.8: rotor mechanical speed and currents injected in the DC bus by the bidirectional DC/DC converter and the Full bridge converter.

Figure 5 .

 5 Figure 5.10: Schematic of the electrical motor testbed together with its Power electronics.

Figure 5 .

 5 Figure 5.11: Electrical motor testbed used for experimental validation and performance assessment of the proposed control strategy.

Figure 5 .

 5 Figure 5.12: DC/DC Full Bridge converter.

  Figure 5.16. As discussed in Section 5.1.3, x st 4 is used to generate the reference values x ref 1 = Ω ref ,x ref 2 = I ref red and x ref 3 = V ref red System associated to the O-Grid Stirling based power Plant.

Figure 5 .

 5 Figure 5.13: Load power demand P L at the input of the inverter for the electrical motor testbed.

Figure 5 .

 5 Figure 5.14: Bus voltage V bus and its reference value V ref bus for the electrical motor testbed.

Figure 5 .

 5 Figure 5.15: Supercapacitor voltage V sc and its reference value V ref sc for the electrical motor testbed.

Figure 5 .

 5 Figure 5.16: Output current of the DC/DC Full Bridge converter I Lf b and its reference value I ref Lf b for the electrical motor testbed.

Figure 5 .

 5 Figure 5.17: Rotational speed of the engine Ω, the rectied voltage V red and the rectied current I red with their respective reference values Ω ref , V ref red and I ref red for the electrical motor testbed.

Figure 5 .

 5 Figure 5.18: Duty ratios of the DC/DC Full Bridge converter (upper gure) and bidirectional DC/DC converter (lower gure) for the electrical motor testbed.

Figure 5 .

 5 Figure 5.19: Electrical cabinets containing the energy conversion system

Figure 5 .

 5 Figure 5.21: Overview of the second electrical cabinet

Figure 5 .

 5 Figure 5.22: Load power demand P L at the input of the inverter for the Stirling based power plant.

Figure 5 .

 5 Figure 5.23: Bus voltage V bus and its reference value V ref bus for the Stirling based power plant.

Figure 5 .

 5 Figure 5.25 shows good results concerning the tracking of the output current of the DC/DC full bridge converter I Lf b to its reference value x ref 4 = I ref Lf b calculated according to (5.11) from

=

  Figure 5.27 plots the control signals u 1 = α f b , used to track x ref 1 Lf b , and u 2 = α bb used to regulate x 5 = V bus to its setpoint x st 5 = V ref bus .

Figure 5 .

 5 Figure 5.25: Output current of the DC/DC Full Bridge converter I Lf b and its reference value I ref Lf b for the Stirling based power plant.

Figure 5 .

 5 Figure 5.26: Rotational speed of the engine Ω, the rectied voltage V red and the rectied current I red with their respective reference values Ω ref , V ref red and I ref red for the Stirling based power plant.

Figure 5 .

 5 Figure 5.27: Duty ratios of the DC/DC Full Bridge converter (upper gure) and bidirectional DC/DC converter (lower gure) for the Stirling based power plant.

3

  and x max 3 invoked in (6.2) are compatible with the constraints on x 2 , namely the steady values of x 2 belong to [x min 2 , x max 2 ] for all constant value of x 3 ∈ [x min 3 , x max 3 ].

Figure 6 .

 6 Figure 6.3 shows how the voltage x 5 is regulated around the desired value x d 5 = 50 V while the desired value x d 4 = 14.7 is considered. Two scenarios are plotted in this Figure for two dierent values of the maximum current x max 6

Figures

  Figures 6.6 and 6.7 shows the impact of the saturation level x min 2

Figure 6

 6 Figure 6.5: Control of the (x 3 , x 4 ) system under constant x 5 = x d 5 after a step change of -20% on the desired value x d 4 . Note that the bounds on x 3 dened by (D.2)-(D.3) slow down the response in order to meet the constraints on x 2 . The same scenario is depicted in Figure6.9 using logarithmic scale in order to show the evolution of the variables during rst milliseconds of the scenario.

98

  Figure 6.6: Control of the (x 3 , x 4 ) system under constant x 5 = x d 5 after a step change of -20% on the desired value x d 4 . Inuence of the constraints on x 2 on the response time of the system.

Figure 6 . 7 :

 67 Figure 6.7: Control of the (x 3 , x 4 ) system under constant x 5 = x d 5 after a step change of +20% on the desired value x d 4 . Inuence of the constraints on x 2 on the response time of the system.

Figure 6 . 8 :

 68 Figure 6.8: Same scenario as the one depicted on Figure 6.4 using logarithmic scale on the time axis in order to show the beginning of the scenario. This scale enables the saturation on x 4 to be clearly shown.

=

  Figure 6.10: Performance of the overall closed-loop system under varying power load demand. Case where x max 2 = 5. This scenario is worth comparing with the scenario depicted in Figure 6.11 where the upper bound x max 2

  Contents 7.1 Presentation of the ORC power plant prototype . . . . . . . . . . . . . 110 7.2 Control objectives and system constraints . . . . . . . . . . . . . . . . 113 7.3 Brief state of the art on physical and semi physical ORC power plants dynamic models . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115 7.4 Wiener Models based Identication Structure . . . . . . . . . . . . . . 117 7.4.1 Presentation of the methodology . . . . . . . . . . . . . . . . . . . . . . . 117 7.4.2 Application of the identication framework to the ORC power system . . 119 7.4.3 Identication Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120 7.5 State Space representation of the ORC power plant identied dynamic models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122 7.5.1 Nonlinear State Space Model of the ORC power plant . . . . . . . . . . . 122 7.5.2 ORC Model Linearization . . . . . . . . . . . . . . . . . . . . . . . . . . . 131 7.5.3 Filtering control variables . . . . . . . . . . . . . . . . . . . . . . . . . . . 132 7.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133

Figure 7 . 2 :

 72 Figure 7.2: Schematic of the Organic Rankine Cycle Power Plant.

Figure 7

 7 Figure 7.3: Temperature-Entropy (T-S) diagram of the Organic Rankine Cycle Power Plant.

Figure 7

 7 Figure 7.4: comparaison between the tted polynomial for T sat (P vap ) with the data provided by REFPROP R .

  Figure 7.6: Allure of the maps B j (.) for β = 0.5 and n m = 10.

Figure 7 . 7 :

 77 Figure 7.7: Control and exogenous variables proles during Experiment 1.

Figure 7

 7 Figure 7.8: Control and exogenous variables proles during Experiment 2.

7. 5

 5 State Space representation of the ORC power plant identied dynamic models 7.5.1 Nonlinear State Space Model of the ORC power plant

ξ

  7.22) 7.5. State Space representation of the ORC power plant identied dynamic models 123 (a) Time proles of the identied and the experimental T vap on the learning and validation data.Gradient of Γ Tvap (ξ)(b) The nonlinear map (given by (7.10)) and its gradient for the T vap -related identication.

Figure 7 . 9 :ξ

 79 Figure 7.9: Identication results for T vap .

Figure 7 .ξ

 7 Figure 7.10: Identication results for P vap .

Figure 7 .

 7 Figure 7.11: Identication results for P elec .

Figure 7 .

 7 Figure 7.13: Identication results for P eM ot.

  Figure 7.14: Identication results for T eM ot.

Figure 7 .

 7 Figure 7.15: Identication results for ∆T vap .

  the following vector of constrained variables: y c = P vap T vap ∆T vap T f r T eM ot P eM ot T ∈ R 6 that has to satisfy according the steady vector Z d is based upon the vector of controlled variables and its past values, it has to satisfy the following inequality constraints:Z min d ≤ Z d ≤ Z max d (8.7)8.1. Study of the admissible steady states 137Thus combining (8.6) and (8.7) and taking into account (8.5), the steady values have to satisfy the following inequality constraints:

Figure

  Figure 8.1: Principle diagram for extracting y d from carto ORC (T aCH vect , T vect amb , P vect elec , ∆T vect vap , T vect vap , P eM ot vect ) based on the knowlege of P ref elec .

  (a) Temperature prole of the measured uncontrolled signals T aCH and T amb . (b) ORC turbine output electrical power P elec and its reference signal P ref elec .

Figure 8 . 3 :

 83 Figure 8.3: Tracking results and disturbance rejection for the ORC turbine output electrical power P elec .

Figure 8

 8 Figure 8.4: ORC system actuators: f pump , Vhot , V aero rot and V circ rot . The red dashed lines indicate the minimum and maximum admissible values.

Figure 8

 8 Figure 8.5: ORC system constrained variables. The red dashed lines indicate the minimum and maximum admissible values.

Figure

  Figure A.1: Schéma global d'une micro centrale solaire thermodynamique

A. 2

 2 Micro Centrale solaire thermodynamique à base d'une machine Stirling Cette section résume les principaux résultats de la Partie I du manuscrit de thèse qui porte sur la gestion de l'énergie électrique d'une centrale solaire thermodynamique utilisant une machine à cycle de Stirling associée à un système de conversion d'énergie approprié. Le schéma d'une telle centrale est illustré à la Figure A.2 Lorsque la charge électrique varie, le moteur Stirling doit adapter sa puissance produite an de compenser celle consommée par la charge électrique mais étant donné que la machine Stirling adapte sa puissance lentement vu son inertie mécanique et thermique, l'hybridation

Figure A. 2 :

 2 Figure A.2: Schéma global d'une micro centrale solaire thermodynamique à base d'un moteur Stirling

A

  .2.1 Modèle statique non linéaire d'une machine Stirling et contrôle de son couple moteur La machine Stirling que nous avons à disposition et qui est étudiée dans le cadre de cette thèse est une machine de type alpha 1 qui est constituée d'un ensemble de 3 cycles liés mécaniquement à un arbre moteur. Chaque cycle est constitué d'un ensemble de deux cylindres-pistons connectés au travers d'un régénérateur. Un des cylindres est connecté grâce à un premier échangeur de chaleur à la source chaude tandis que l'autre cylindre est connecté à la source froide à travers un autre échangeur de chaleur. La machine Stirling fournit un travail mécanique, qui sert à entrainer une Génératrice Synchrone à Aimants Permanents (MSAP), en eectuant un cycle thermodynamique fermé (qui est subit par le gaz de travail dans les cylindres et le régénérateur) à chaque rotation complète de l'arbre moteur comme illustré à la Figure A.3 dans le plan Pression-Volume (P-V).

Figure A. 3 :

 3 Figure A.3: diagramme Pression-Volume (P-V) du cycle de Stirling et position des pistons pour un moteur Stirling de type alpha pendant un cycle complet.

Figure A. 4 :où

 4 Figure A.4: Vue d'ensemble du système de production d'électricité à base de machine Stirling.

Figure A. 8 :

 8 Figure A.8: Architecture de la première loi de commande pour le système de conversion d'énergie associé à la machine Stirling.

4=

  -k 6 . tanh(β.(x 7x st 7 nouvelle valeur de référence x4

Figure A. 13 :

 13 Figure A.13: Schéma global d'une micro centrale solaire thermodynamique à base d'une machine ORC

  Figure A.14: Schéma du système de production d'électricité à base d'une machine à cycle de Rankine organique (ORC).

Figure

  Figure A.15: diagramme Température-Entropie (T-S) du système de production d'électricité à base d'une machine à cycle de Rankine organique (ORC).

ξ

  (a) Proles temporels des puissance électrique P elec mesurée et identiée sur l'ensemble des données d'apprentissage et de validation. Gradient of Γ P elec (ξ) (b) Fonction non linéaire et son gradient pour l'identication de P elec .

Figure A. 17 :

 17 Figure A.17: Résultats d'identication pour P elec .

Figure

  Figure A.18: Architecture de contrôle du système ORC.

  or T sw 2 + α f b . Tsw 2 < t < T swThe DC/DC full bridge converter equations become:L f b . dI Lf b dt = k.V red .u f b -V bus (B.4a) k.I Lf b .u f b (B.4c)

  Figure C.1: Overview of the Stirling based power plant and its associated power conversion system including the ATV32 variable speed drive

I

  Figure C.2: Starting procedure owchart

  Table of symbols and acronyms1 n Ones matrix of dimension n, 1 n ∈ R n O n i ×n j Zeros matrix of dimension n i × n j , O n i ×n j ∈ R n i ×n j O n Zeros matrix of dimension n × n, O n×n ∈ R n×n

	Acronyms
	Control Theory
	PID	Proportional Integral Derivative (controller)
	MPC	Model Predictive Control
	QP	Quadratic Programming (Optimization problem)
	SQP	Sequential Quadratic Programming (Optimization problem)
	Thermodynamics
	ORC	Organic Rankine Cycle
	NTU	Number of Transfer Units
	P-V	Pressure-Volume (diagram)
	T-S	Temperature-Entropy (diagram)
	WHR	Waste Heat Recovery
	FV	Finite Volume
	MB	Moving Boundary
	Electrical engineering
	MPPT	Maximum Power Point Tracking

AC

Alternating Current DC Direct Current

Table of

 of

		symbols and acronyms	xxi
	PV	photovoltaic
	UPS	Uninteruptible Power Source
	SOC	State Of Charge
	RPM	Rotation Per Minute
	PMSG	Permanent Magnet Synchronous Generator
	EMF	Electromotive Force
	IGBT	Insulated Gate Bipolar Transistor
	PWM	Pulse Width Modulation
	Other

BipBop Business, Innovation and People at the Base of the Pyramid PLC Programmable Logic Controller NI National Instrument R DLL Dynamic Link Library Denition 0.1 Let A ∈ R n×m and B ∈ R p×q . The Kronecker product A ⊗ B is dened by:

  Context of the thesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Aims of the thesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4

				Chapter 1
		General Introduction
	Contents 1.1
	. . .	. . .	. . .	. . .

) • B 1.1 Context of the thesis

  In (5.1f), to regulate x 6 around x ref 56 Chapter 5. Development of a rst control strategy for the Energy Conversion System associated to the O-Grid Stirling based power Plant. The regulation error being 6 = x 6x ref 6 , we use the following positive denite function: V 6 = V 5 + 1 2 . 2 6 as a Lyapunov function for subsystem 2 where its time derivative is given by: V6 = 5 . ˙ 5 + 6 . ˙ 6 = 5 .(a 10 .x 4 + a 10 .(x ref 6 + 6 ) -.(-a 11 .x 5 + a 11 .x 7 .u 2 -ẋ6 ref ) = -ρ 5 . 2 5 + + 6 .(-a 11 .x 5 + a 11 .x 7 .u 2 -ẋ6 ref + a 10 . 5 )

	a 10 η inv	.	P L x 5	)
	+ 6

92

  Chapter 6. Constrained control framework for the Energy Conversion System associated to the O-Grid Stirling based thermodynamic power Plant • Compute the abscissa x5 of the horizontal projection of p on BCD • Compute the vertical interpolation u v 2 (x) according to

Table 7 .

 7 1: The identied parameters used to dene the nonlinear maps for ∆T vap , T vap , P vap , P elec , T f r , P eM ot and T eM ot n m β N m n u n z n p CompressionF actor = N data Section 7.4.1 for the parameters denition) as it can be seen from Figures 7.9a-7.11a and 7.12-7.15 that compare the time proles of the experimental and identied variables of interest (namely: T vap , P vap , P elec , T f r , P eM ot, T eM ot and ∆T vap ) on the whole data set for both experiments. One can observe from

	np

  Study of the admissible steady states . . . . . . . . . . . . . . . . . . . 135 8.2 Constrained Model Predictive Control of the ORC power plant . . . 140 8.2.1 Measured uncontrolled signals handling . . . . . . . . . . . . . . . . . . . 140 8.2.2 MPC controller design . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 141

	8
	Constrained Optimal Control of the
	ORC based power plant
	Contents 8.1

8.2.3 Simulation Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 142 8.3 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 145

BipBop program description: http://www2.schneider-electric.com/sites/corporate/en/group/ sustainable-development-and-foundation/access-to-energy/presentation.page

Commissariat à l'énergie atomique et aux énergies alternatives http://www.cea.fr/

SCHOTT website:http://www.schott.com/csp/french/download/schott_csp-datenblatt_ ptr70-classic_fra.pdf

http://www.eureqa.com/

, the control variable u 2 is used.

http://ust.su/upload/iblock/892/WBP-Box-IEN103320.pdf

http://www.schneider-electric.com/download/hk/en/details/2619521-Altivar-32-variable-speed-drives-catalogu ?showAsIframe=false&reference=DIA2ED2100401EN

are determined by the dimensioning of the system

http://www.ohio.edu/mechanical/stirling/

http://www.schneider-electric.com/download/hk/en/details/2619521-Altivar-32-variable-speed-drives-catalogu ?showAsIframe=false&reference=DIA2ED2100401EN 181

Remerciements A.4 Conclusion générale

Cette thèse a été dédiée au développement de lois de commande performantes pour deux types de micros centrales solaires thermodynamiques considérées dans le cadre du projet MI-CROSOL.

La première centrale solaire étudiée utilise une machine thermodynamique à cycle de Stirling associée à un système de conversion d'énergie adéquat. Un premier contrôleur optimal du couple de la machine Stirling a été développé en se basant sur un modèle validé expérimentalement. Pour le contrôle du système de conversion d'énergie, deux lois de commandes ont été développées à cet eet. La première loi de commande développée gère explicitement les contraintes sur les saturations des variables de contrôle tout en étant simple à régler. Elle a été validée expérimentalement sur un prototype fonctionnel à base de machine Stirling et ore les performances attendues. La seconde loi de commande développée gère explicitement les contraintes sur tous les états du système de conversion d'énergie tout en dotant le système en boucle fermé de performances optimales en termes de temps de réponse.

Par la suite, le problème du contrôle optimal sous contraintes d'une micro centrale solaire utilisant une machine ORC a été abordé et résolu en ayant une approche qui part des données expérimentales recueillis sur un prototype fonctionnel an d'identier des modèles dynamiques de la centrale. Ces modèles nous ont par la suite servis à élaborer un contrôleur MPC qui fait passer la machine ORC d'un point de fonctionnement à un autre en puissance le plus rapidement possible tout en respectant les contraintes physiques du système Les temps de réponse obtenus en simulation avec notre approche suggèrent qu'une réduction signicative du nombre de batteries, nécessaires comme source d'énergie tampon, pourrait être envisagée.

Chapter 6. Constrained control framework for the Energy Conversion System associated to the O-Grid Stirling based thermodynamic power Plant supercapacitor voltage x 7 = V sc to its reference value, in this section a feedback law for the control variable u 1 that tacks x st 4 in near-to-optimal manner while respecting the system constraints is derived using "subsystem 1" equations given by (6.3) in which x 5 = x r 5 = x st 5 .

Because of the limitations on x 3 given by (6.5) regarding the fulllment of the saturation constraints on x 2 (refer to Lemma 6.1), a new reference value for x 4 namely x * 4 = x st 4 (x * 3 ) has to be derived and tracked acting on u 1 , where x * 3 is dened to be the closest value to x st 3 (x st 4 ) [see (6.8)] that is compatible with the requirement on the excursion of x 2 . More precisely:

x * 3 (x, P L ) := arg min

|z 3x st 3 (x st 4 (x, P L ))| (6.14) where x st 4 (x, P L ) is given by (6.8) depending on the current value of the supercapacitor voltage x 7 and the present load demand P L . Note that by doing so, the rate with which the supercapacitor voltage is regulated around its reference value is moderated by the dynamic of the mechanical part of the system which is directly linked to the Stirling engine performance. This is done through the denition of the admissible interval x 3 x 1 (t), x 2 (t) , x3 x 1 (t), x 2 (t) . Remark 6.1 Note that in Chapter 5 (especially in Section 5.1.3) where a rst control strategy for the energy conversion system has been developed, it has been observed that in order to avoid negative excursions on some of the system's variables, the variation of the set-point x st 4 needs to be ltered or slowed down. However, contrary to Chapter 5 where this feature has been observed, the discussion above gives a deeper interpretation and a more systematic, ecient and dynamic way to address this limitation. Indeed, while the lter introduced in Chapter 5 operates continuously slowing down the system unconditionally, the values that result from (6.14) slow down the system only when necessary. This is a general fact making explicit handling of the constraint preferable to worst-case-based unconstrained design.

Based on the observation made in Section 6.1 regarding the excursion of the stationary values of x 2 , it can be argued that if the conclusion of Lemma 6.1 is used to dene bounds on the variations of x 3 such that x 2 remains in a tight interval around x st,c 2 , then it becomes relevant to analyze the behavior of the dynamic system (6.1c)-(6.1d) in which the state is (x 3 , x 4 ) while x 2 and x 5 are considered to be constant. Such an analysis can then be performed by visualizing the vector elds in the plan (x 3 , x 4 ). When doing this analysis using the numerical values of the system under consideration, it comes that the following assumption is satised: Assumption 6.5

For any admissible desired pair (x * 3 , x st 4 (x *

3 )), the picture of the vector elds of the dynamic system (6.1c)-(6.1d) in which x 2 = x st 2 (x * 3 ) ∈ I 2 and x 5 = x r 5 is the one depicted in Figure 6.1, more precisely:

1. The horizontal lines AB and DE correspond to the lines x 4 = x min 4 and x 4 = x max 4 respectively.

2. C is the desired steady point (x * 3 , x st 4 (x *

3 )) while u *

3 ) is the corresponding steady control.

Chapter 6. Constrained control framework for the Energy Conversion System associated to the O-Grid Stirling based thermodynamic power Plant Note that the rst 4 items of Assumption 6.5 are only denitions. The item 5) can be easily checked by investigating the disposition of the vector elds on the horizontal lines AB and DE. It is also possible to formulate algebraic conditions [see (D.4)-(D.7) ] to check items 5), 6) and 7) above as it is shown in Appendix D.2.

Although the conditions (D.4)-(D.7) might seem tricky to check, it is likely that the topology of the vector elds depicted in Figure 6.1 holds for a wide class of such voltage conditioning systems. The above conditions are only given for completeness.

Note also that the disposition of the vector elds obviously suggests that a sliding modes controller can be dened in which the curves ABCDE is the sliding surface while the equivalent control on the sliding surface steers the system towards the desired point C. Here we prefer to use continuous control because the equations (6.1a)-(6.1d) represent the mean dynamic and in order to derive an easy to use closed-loop systems avoiding the numerical diculties associated to the sliding modes controllers simulation when long tuning and assessing scenarios are required as it is the case here (see later).

Rather than using the sliding modes formulation, a 2D interpolation techniques is used in the (x 3 , x 4 ) plan. Namely, an interpolated expression of the feedback is derived in which, the expressions u 1 = 1 or u 1 = 0 are obtained when the current state (x 3 , x 4 ) is far from the curve ABCDE while u 1 takes the appropriate values when (x 3 , x 4 ) lies on the curve ABCDE. The appropriate value is itself obtained by appropriate interpolation between the values depicted on Figure 6.1 according to the position of the current point (x 3 , x 4 ) on the curve ABCDE.

In order to properly dene the 2D interpolated expression, the following auxiliary function is needed:

which simply takes the value v 1 for r = 0 and v 2 when r 0 with a transition rate dened by the positive constant λ.

Based on the above discussion, the following explicit feedback can be dened for a given target point (x * 3 , x st 4 (x * 3 )):

)) be dened. Given a current position p := (x 3 , x 4 ) such that

], compute the following quantities:

• Compute u DE 1 (x) and u AB 1 (x) according to:

, 0, ū1 (x) (6.17)

where ū1 (x) is given by (6.15).

Constrained Sliding Mode control of the energy conversion system 89

• Compute u BCD 1 (x) according to:

• Compute the abscissa x3 of the horizontal projection of p on BCD.

• Compute the vertical interpolation u v 1 (x) according to:

3 ) by horizontal interpolation between u v 1 (x) and u BCD 1 (x) according to:

Indeed, for a given value of x 4 , continuity w.r.t x 3 is guaranteed by the continuity w.r.t x 3 of the expressions (6.17), (6.18) and (6.19) while the continuity w.r.t x 4 for a given x 3 is guaranteed by the continuity w.r.t x 4 of the expressions of u BCD 1 (x) and u v 1 (x).

Note that the success of the state feedback u 1 (x|x * 3 ) in steering the state (x 3 , x 4 ) to the desired (x * 3 , x * 4 ) is conditioned by the fact that the excursion of x 2 is small. Indeed, under this assumption, the picture depicted in Figure 6.1 remains relevant while x 2 is slightly moving. Note however that for this excursion to be small, the following two requirements are needed:

• The range of variation ε 2 given by (6.4) of the steady values of x 2 is small.

• The variations of x 3 remains inside the moving interval

The rst requirement mentioned above depends on the choice of the system's components and has to be satised by construction. The second requirement can be satised by an appropriate management of the desired value x * 3 according to (6.14).

The preceding discussion completely denes the feedback law

The next section is dedicated to the denition of the feedback law u 2 = K 2 (x, P L ) through the analysis of "subsystem 2" equations.

Chapter 8. Constrained Optimal Control of the ORC based power plant

Constrained Model Predictive Control of the ORC power plant

Now that we have dened the control-oriented mathematical model for the ORC power plant in Section 7.5 and translated mathematically in Section 8.1 the control objectives and system constraints in terms of the admissible desired steady states, the last task is to design a control strategy that optimally tracks some desired reference electrical power P ref elec while respecting the system constraints. In this work and according to the system model structure, an MPC (Model Predictive Control) controller is the best suited option to achieve the control objectives since it handles explicitly the system constraints while insuring an optimal behavior. The principle of an MPC controller is to solve online a nite horizon open-loop constrained optimization problem that depends on the current system state (rather than determining an o-line feedback law that holds for all the states) to derive at each sampling period an optimal control sequence from which the rst control is applied to the system and the procedure is repeated at each sampling time to cope with the system disturbances.

MPC controllers rely on the use of a model for the system to be controlled. When the model of the plant is linear, the optimization problem consists in solving online a classical Quadratic Programming (QP) problem at each sampling time (refer to [START_REF] Mayne | Constrained model predictive control:stability and optimality[END_REF] where a large literature review on MPC is presented).

Among the numerous advantages of Model Predictive Control:

• The possibility to handle explicitly the system constraints regardless their nature: control saturation, states and output limitations...etc.

• It is well adapted to control multi-variable coupled dynamical systems.

• The fact that MPC design relies on solving an optimization problem endows the methodology by some intuitiveness .

All these advantages made the MPC design widely used in the industry with thousands of succesful implementations (see [START_REF] Qin | A survey of industrial model predictive control technology[END_REF] for some examples of industrial applications of MPC controllers).

In Section 8.2.1, the ORC power plant model is reformulated to handle the system disturbances, then the mathematical formulation of the MPC control problem is presented in Section 8.2.2. Simulation results of the resulting closed loop system will be discussed in Section 8.2.3.

Measured uncontrolled signals handling

If we suppose that the vector of operational variables remains constant during the whole prediction horizon (which is a good assumption since the dynamic of T aCh and T amb is very 8.2. Constrained Model Predictive Control of the ORC power plant 141 slow compared to the dynamics of the other variables), the state space equations for the ORC power plant (7.49) become:

By extending the state vector with the vector of the uncontrolled variables u unc (k) as follow: Based on the extended model (8.13) for the ORC power plant, an MPC controller is designed to track the reference electrical power P ref elec (that will serve to dene y d according to Figure 8.1) while respecting the system constraints.

The objective function to be minimized is given by:

under the constraints dened by:

where y c = P vap T vap ∆T vap T f r T eM ot P eM ot T ∈ R 6 is the vector of constrained variables. Q and R are the weighting matrices of the MPC controller, N p is the prediction horizon and

T is the decision variable from which we extract u f (k) the reference of the rst order lter (7.47b) whose output denes the control vector u con (k) applied to the system.

Appendix A

Résumé en français

Cette partie est un résumé en français du manuscrit de thèse. Etant donné que ce résumé n'est pas exhaustif, le lecteur est cordialement prié de se reporter à la version anglaise pour plus de détails techniques.

A.1 Introduction

Aujourd'hui, plus de 1.4 milliards de personnes n'ont pas accès à l'énergie. Ces populations, qui gagnent moins de 2 dollars par jours et dépensent plus de 15 dollars par mois et par famille pour avoir accès à l'électricité, vivent essentiellement en Inde, Indonésie, Bangladesh, Nigéria et dans certaines régions de l'Afrique sub-saharienne.

An d'apporter une solution pour ce problème d'accès à l'énergie pour ces populations, Schneider Electric a lancé un projet collaboratif appelé MICROSOL qui ÷uvre pour le développement de micros centrales solaires thermodynamiques utilisant l'énergie thermique du soleil pour produire de l'électricité. Le marché visé par le projet MICROSOL étant celui de l'électrication rurale pour des sites "isolés" qui ne sont pas connectés au réseau électrique ou bien qui disposent d'un réseau électrique très peu able. Ces centrales solaires doivent fournir au minimum 150kW h elec par jour et fonctionner 24h/24 grâce à un stockage d'énergie thermique approprié. De plus, ces centrales doivent être ables, respectueuses de l'environnement et accessibles à des populations ne disposant pas de beaucoup de ressources nancières mais qui sont situés dans des zones ayant un fort potentiel d'ensoleillement.

La Figure A.1 donne le schéma d'une micro centrale solaire thermodynamique typique. L'un de ses principaux éléments est le champ de collecteurs solaires qui utilise l'énergie thermique du soleil pour chauer un uide caloporteur qui, dans le cadre du projet MICROSOL, est simplement de l'eau pressurisée an d'être respectueux de l'environnement. Cette eau ainsi chauée sera par la suite stockée dans une cuve qui doit assurer un minimum de pertes thermiques. La centrale solaire se compose également d'une boucle de refroidissement qui permet de refroidir un autre uide caloporteur. Les deux uides, chaud et froid, entrent à l'intérieur d'une machine thermodynamique qui utilise le gradient de température (entre la source chaude et la source froide) pour fournir un travail mécanique en eectuant un cycle thermodynamique adéquat. Ce travail mécanique généré sert à entrainer une génératrice électrique produisant de l'électricité qui sera conditionnée grâce à un système de conversion d'énergie adéquat qui s'interface entre la génératrice électrique et la charge électrique isolée que l'on veut alimenter. La première étape pour contrôler ce système a consisté à établir son modèle dynamique qui est donné par les équations suivantes mises sous forme d'état:

Les paramètres a i dépendent de certaines grandeurs physiques du système (k et η inv sont quand à eux des grandeurs physiques). En particulier a 1 et a 2 dépendent entre autre du couple moteur de la machine Stirling dont l'équation a été identiée expérimentalement.

Les objectifs de contrôle pour ce système sont:

• La régulation d'une manière assez raide de la tension du bus continu x 5 =V bus autour de x st 5 =V ref bus =50V compatible avec un bon fonctionnement de l'onduleur commercial autorégulé.

• La régulation de la tension du supercondensateur x 7 =V sc autour de x st 7 =V ref sc , d'une manière plus lente, an d'avoir une quantité d'énergie tampon adéquate.

Ce système est soumis aux contraintes suivantes: Contraintes de positivité: x i ≥ 0 à l'exception de x 6 .

Fortes saturations sur les variables de contrôle:

• u 1 ∈ [0, 1]: rapport cyclique du convertisseur DC/DC à pont complet.

• u 2 ∈ [0, 1]: rapport cyclique du convertisseur DC/DC bidirectionnel en courant.

De plus, le système de conversion d'énergie est sujet aux perturbations suivantes:

• P L : Puissance électrique consommée par la charge isolée qui peut varier d'une manière très rapide. Vu la structure du modèle du système de conversion d'énergie donné par (A.2a)-(A.2g), il est possible, en négligent le terme "a 5 .x 1 .x 2 " correspondant à l'empiétement du pont de diodes, de décomposer le système en deux sous-systèmes qui évoluent dans deux échelles de temps diérentes:

1) Le premier sous-système (sous-système 1) donné par (A.2a)-(A.2d) ne dépend que de u 1 comme variable de contrôle et de a 2 comme signal exogène qui dépend du couple moteur du Stirling. Ce premier sous système dépend également de x 5 qui est supposé parfaitement régulé autour de sa valeur désirée x st 5 en utilisant les équations du deuxième sous-système. Ainsi le sous-système 1 peut se récrire comme suit:

2) Le deuxième sous-système (sous-système 2) correspond aux équations (A.2e)-(A.2g) qui ne dépendent que de u 2 comme variable de commande.

Les deux sous-systèmes sont liés grâce à la variable d'état x 4 . La Figure A.8 donne l'architecture de la première loi de commande développée pour le système de conversion d'énergie associé à la machine Stirling. • Le fait qu'elle nécessite peu de paramètres de réglage.

• Le fait qu'elle est compatible temps réel.

Cependant cette loi de commande ne gère pas explicitement les contraintes sur les variables d'état du système, à la place, un ltre ad-hoc a été introduit pour élargir le domaine de validité de la loi de commande ce qui conduit à ralentir inconditionnellement le système. Pour pallier à cette limitation une deuxième loi de commande a été développée en se basant sur le principe des modes glissants. Préalablement à cela, des limitations dynamique sur x 3 , qui tiennent compte des valeurs minimum et maximum admissibles pour x 2 et x 3 , ont été établi d'une manière explicite (Se référer au Chapitre 6 pour plus de détails concernant la dénition de ces limitations dy- namiques): • Si un point dans le plan de phase (x 3 , x 4 ) se situe sur la surface de glissement, en appliquant la commande appropriée u 1 , on atteint en un minimum de temps la paire désirée (x * 3 , x * 4 ).

• La surface de glissement ABCDE est attractive.

Pour atteindre une paire désirée (x st 5 , x r 6 ), la même analyse graphique dans le plan de phase (x 5 , x 6 ) peut être eectuée.

An d'illustrer les performances de cette deuxième loi de commande développée pour le système de conversion d'énergie associé à la machine Stirling, des simulations ont été eectuées. La

1.7bar ≤ P eM ot ≤ 7bar (A.6e)

où ∆T vap est la valeur de la surchaue en entrée de la turbine donnée par:

avec T sat la température de saturation du uide de travail.

f pump , Vhot , V aero rot , et V circ rot sont les 4 variables de contrôle du système et représentent respectivement la fréquence de rotation du moteur de la pompe volumétrique de la boucle interne, le débit volumique de la pompe de la boucle d'eau chaude, la vitesse de rotation relative de l'aérotherme et la vitesse de rotation relative du moteur de la pompe de circulation de la boucle froide.

De plus la machine ORC est soumise aux perturbations (variables exogènes non controlables) mesurables suivantes:

• T aCh: température en Kelvin de l'eau chaude en entrée de l'évaporateur.

• T amb : température ambiante en Kelvin.

A.3.3 Identication d'un modèle dynamique pour la machine ORC

Le problème d'identication consiste à trouver une fonction non linéaire F y qui lie chaque variable d'intérêt y(k.τ s ) = y(k) (où τ s est la période d'échantillonnage des signaux), qui peut être la puissance électrique régulée P elec ou bien l'une des variables sous contraintes (∆T vap , T vap , P vap , T eM ot, P eM ot et T f r ), à un régresseur associé Z y (k) déni par:

. . .

où n y u est la dimension du vecteur d'entré u y (k) (qui peut contenir aussi bien des variables de contrôle que des variables non contrôlables mais mesurables) et N y est la profondeur du régresseur. La fonction non linéaire F y (Z) a été restreinte à avoir la forme suivante:

y (Z T y .L y ); Γ y (.) strictement croissante (A.9) où L y est un vecteur à identier. Le lecteur pourra trouver l'allure des fonctions non linéaires Γ y utilisées en se référant au manuscrit de thèse.

La procédure d'identication s'éectue sur un ensemble de données d'apprentissage Une que fois les modèles dynamique pour chaque variable d'intérêt y ont été identié, il est possible après un certain nombre d'opérations (se reporter au manuscrit de thèse pour plus de détails) d'aboutir à un modèle linéarisé au sens des moindre carré pour chacune des variables d'intérêt sous la forme suivante:

où X est un vecteur d'état augmenté qui contient le vecteur d'état d'origine x ainsi que le vecteur de contrôle ltré u con . Le paramètre N est le maximum des profondeurs des régresseurs pour toutes les variables thermodynamiques d'intérêt.

A.3.4 Etude des états stationnaires désirés admissibles

Avant de développer le contrôleur optimal sous contrainte pour la machine ORC, il est intéressant d'étudier tous les états stationnaires admissibles (c'est-à-dire qui respectent les contraintes du système) et qui sont compatibles avec une certaine puissance électrique désirée P ref elec . Pour cela, partant des équations d'état linéarisées (A.11) du système, si on prend comme seule sortie régulée la puissance électrique, et vu qu'on a 4 actionneurs, alors pour chaque puissance électrique désirée il peu y avoir plusieurs états stationnaires compatibles mais qui ne respectent pas forcément tous les contraintes du système. Pour ne converger vers qu'un seul régime stationnaire que l'on maîtrise, on augmente le vecteur de sortie régulé y r avec 3 autres variables thermodynamiques qui sont critiques pour le système: ∆T vap , T vap et P eM ot. On aboutit aux équations d'état suivantes:

Ainsi, il est possible de reformuler le problème qui consiste à trouver les états stationnaires admissible Xd et le vecteur de commande associé en régime stationnaire u f d sous la forme problème d'optimisation LP 5 où l'on cherche à trouver le vecteur d'état augmenté

sous les contraintes suivantes: .15) où f = 1 (nx+n con u ) est un vecteur de pondération. L'équation (A.14) correspond à l'équation d'état stationnaire et à l'équation de sortie stationnaire ou y d est la vecteur désiré correspondant à y r . L'inégalité (A.15) reète quand à elle le fait que tous les états stationnaires doivent respecter les contraintes du système.

Etant donné que ce problème d'optimisation LP est assez lourd pour être résolu en ligne vu la taille des matrices associés et vu qu'on ne connait pas à priori les valeurs des autres composantes du vecteur désiré y d compatibles avec P ref elec , on a préféré utiliser une cartographie calculée hors ligne (dénotée Carto ORC ) qui balaye toutes les valeurs admissibles de y d ainsi que toutes les valeurs admissibles des variables exogènes et qui contient des 0 et des 1 suivant la faisabilité du problème LP. On extrait par la suite en ligne à partir de la cartographie et de ses vecteurs d'appui associés le vecteur y d correspondant à une puissance électrique désirée et à une certaine valeur des conditions opératoires (variables exogènes mesurées) en Here we present the methodology for deriving mean models for the controlled DC/DC full bridge and the bidirectional DC/DC converters which constitute a part of the energy conversion system associated to the Stirling based solar power plant.

B.1 Average modeling of the DC/DC full bridge converter

The DC/DC full bridge converter depicted in Figure 4.4 is controlled acting on its associated IGBTs namely T r 1 , T r 2 , T r 3 and T r 4 in such a way that during a cutting period T sw of this converter, the two IGBT pairs (T r 1 , T r 3 ) and (T r 2 , T r 4 ) are activated during the same controllable period t f but 180 • out of phase according to the the diagram of

one can set its corresponding electrical equations during a cutting frequency as follow:

According to Figure 4.4 and the state of the IGBTs, the equations of the converter are:

One obtains the same equations as (B.2)

One obtains the same equations as (B.3) By dening the following control variable u f b (t) over the cutting frequency: In order to get a continuous conduction mode of the bidirectional DC/DC converter when the output current I Lbb is close to zero, it is interesting to control this converter in the complementary mode in which T r 1 is activated (while T r 2 is disabled) during a period α bb .T sw , where α bb ∈ [0, 1] is the new duty ratio of the converter, and T r 2 is activated (while T r 1 is disabled) during the rest of the cutting period T sw . The equations governing this mode are as follow:

When where:

In order to get the mean model of the bidirectional DC/DC converter used during the control design, we apply the same moving average equation (B.5) to the set of equations (B.9) leading to the following average model of the converter:

Starting procedure for the Stirling based power system prototype Hereafter, we present the starting procedure that brings the Stirling based power plant prototype depicted in Figure 3.7 to its nominal operating mode in which the control strategy developed in Chapter 5 holds. The presented starting procedure holds also for the dedicated motor testbed depicted in Figure 5.11 with some slight modications.

As mentioned in Section 5.3.1, the Stirling based power plant illustrated in Figure 3.7 contains an ATV32 1 speed drive (not represented in Figure 4.1 because this gure concerns the nominal operating mode of the prototype in which the speed drive is disconnected) in addition of the energy conversion system as illustrated in Figure C.1 This speed drive enables to bring the Stirling shaft to some preset speed value Ω start to initially enable the movement of the Stirling engine pistons hence creating some motor torque (when the temperature of the heat source is high enough and the valves inside the Stirling engine between the cylinders and the pressure vessel are closed).

Figure C.2 gives the starting procedure owchart of the Stirling based power system prototype. In this gure the variable "start_system" takes the value 1 when the user pushes a "start system" button on the supervisor (supervision PC of Figure 3.7). In this owchart, the "Switch" corresponds to the one in Figure C.1 which is initially closed to enable the speed drive to control the shaft speed of the PMSG. The variable "Mode_Regule_FB" indicates the regulation mode of the DC/DC full bridge converter as explained hereafter.

C.1 Case where (Mode_Regule_FB)=1

This case happens when the supercapacitor voltage V sc is greater than the bus voltage V bus . In this case, the DC/DC full bridge converter regulates the DC bus voltage V bus to its reference value V ref bus =50V compatible with the good functioning of the inverter.

The regulation is performed using standard PI based cascaded current-voltage loops using where K f b1 is some design parameter. This leads to the following equation for the control variable α f b :

C.2 Case where (Mode_Regule_FB)=2

This case happens when the supercapacitor voltage V sc is less than or equal to the bus voltage V bus . In this situation the diode D2 of the bidirectional DC/DC converter (see Figure 4.5) becomes conductive and then, the supercapacitor C sc becomes in parallel with the DC bus capacitors C tot resulting in greater DC bus total capacitor C tot = C sc + C tot . In this mode the DC/DC full bridge converter regulates also the DC bus voltage V bus to its reference value V ref bus =50V using the same equations developed in Section C.1 (for the reference value of I Lf b and the control α f b ) but replacing C tot by C tot .

C.3 Case where (Mode_Regule_FB)=3

This mode corresponds to the nominal operational mode in which both converters are enabled and for which the control equations were developed in Chapter 5.

Appendix D Proofs related to chapter 6

In this appendix, we present the proofs of: Lemma 6.1, some items of Assumption 6.5 and Proposition 6.1 all of them used during the design step of the constrained sliding mode controller for the energy conversion system associated to the Stirling engine.

D.1 Proof of Lemma 6.1

According to Assumption 6.1, the linear subsystem in (x 1 , x 2 ) dened by (6.1a)-(6.1b) steered by the exogenous signal x 3 shows two eigen-modes: a fast one denoted by z f and a slow one denoted by z s . The output we are interested in for this subsystem is y = x 2 and is given, according to the denition of the eigen-modes, by:

The objective of this proof is to characterize all constant proles on x 3 leading to the respect of the saturation constraints (6.2) on x 2 starting from an admissible initial condition.

To do this, consider the temporal generic structure of the trajectory y(•) depicted in Figure D.1 where the rst order behaviors of the eigen-modes are approximated by straight lines. According to this gure, the position in time of the singular point (zero time-derivative) is close to the settling time (t * ) of the fastest modes and this regardless of the value of the constant steering signal x 3 . This is because after this value t * , the derivative of y becomes monotonic and cannot change its sign anymore. Now for this once-for-all xed value t * , the expression of the output y(t * ) = x 2 (t * ) can be given in terms of the current value of (x 1 , x 2 ) and the input x 3 using constant and o-line computable matrices M * , N * and K * according to:

This obviously leads to the following expressions for x 3 and x3 :

x 3 (x 1 , x 2 ) := max x min 3 , Now if x 3 is not constant but satises at each instant the constraint (6.5), the result can be proved by obvious receding-horizon way by looking at the rst order two eigen-systems. This obviously ends the proof of Lemma 6.1.

D.2 Mathematical conditions to check items 5), 6) and 7) of Assumption 6.5

The Item 5) of Assumption 6.5 can be checked through the following conditions that can be checked by scalar optimization in the decision variable x 3 : min (x 3 ,x min

4

)∈AB

x st 2 (x * 3 ) -

x min 4 x r 5

x 3 > 0 (D.4) max (x 3 ,x max

4

)∈DE

x st 2 (x * 3 ) -

x max 4

x r 5

x 3 < 0 (D.5) since the analyzed expression is nothing but the resulting derivative of x 3 when x 2 , x 4 and x 5 take the appropriate values.

As for the items 6) and 7) regarding the atractivity of the set BCD, it can be checked through the following conditions (in which the notation p = (x 3 , x 4 ) T is used for compactness): where the subsets R ∈ R 2 and L ∈ R 2 are the sets that denote the values of (x 3 , x 4 ) that lie to the right and to the left of BCD respectively. p(t, p 0 , u 1 = 1) denotes the trajectories of D.3. Proof of proposition 6.1 187 (6.1c)-(6.1d) starting from p 0 under the control u 1 = 1 and when x 5 = x r 5 and x 2 = x st 2 (x * 3 ) are used. Similar denitions hold for p(t, p 0 , u 1 = 0) except that the control u 1 = 0 is used. The notation d(p, ABCDE) denotes the distance between p and the curves ABCDE. Finally, the nal time t f is some suciently long computation time in order for the trajectory to cross the curve ABCDE. D.3 Proof of proposition 6.1

Note that by denition of the interpolated feedback laws K 1 (x, P L ) and K 2 (x, P L ), the resulting closed-loop trajectories passes through neighborhoods V 34 (x * 3 , x * 4 ) and V 56 (x 5 r, x r 6 ) of the targeted pairs (x * 3 , x * 4 ) and (x r 5 , x r 6 ) of size O(1/λ) after a nite time.

Let rst γ 34 > 0 and γ 56 > 0 be taken suciently small to guarantee that the corresponding level sets belong to the regions of attraction of the local controls u st 1 (x * 3 (x, P L )) and u * 2 (x, P L ). Once this is done, take λ > 0 suciently high to guarantee that each of the conditions d 34 (x, P L ) < γ 34 and d 56 (x, P L ) < γ 56 will be reached in nite time so that the local controllers can be red. Now since the switching conditions in (6.27)-(6.28) are dened using the Lyapunov level sets associated to the local stabilizing control laws, no chattering behavior occurs and the state is asymptotically steered to the stationary state x st dened by:

x st i (x d 4 ) i=1,. . . ,4 (D.8)

x st 5 = x r 5 (D.9)

x st 6 = 0 (D.10)

x st 7 = x r 7 (D.11) where x d 4 is given by:

x d 4 := P L η inv x r 5 (D.12) Indeed, this last equation comes from (6.8) because of the asymptotic convergence of x 7 towards x r 7 that results from (6.9). This with (6.7) shows that x 6 converges to 0. This obviously ends the proof.