## Aggregate Power Control of Heterogeneous TCL Populations Governed by Fokker-Planck Equations

## Jun Zheng , Gabriel Laparra , Guchuan Zhu , Meng Lia 1 2 2 3

1

School of Mathematics, Southwest Jiaotong University, Chengdu, Sichuan, P. R. of China 611756 2 Department of Electrical Engineering, Polytechnique Montr´ eal, P . O. Box 6079, Station Centre-Ville, Montreal, QC, Canada H3T 1J4 3 School of Information Science and Technology, Zhejiang Sci-Tech University, Zhejiang, China

## Abstract

This paper addresses the modeling and control of heterogeneous populations of thermostatically controlled loads (TCLs) operated by model predictive control (MPC) schemes at level of each TCL. It is shown that the dynamics of such TCLs populations can be described by a pair of Fokker-Planck equations coupled via the actions in the domain. The technique of input-output feedback linearization is used in the design of aggregate power control, which leads to a nonlinear system in closed loop. Well-posedness analysis is carried out to validate the developed control scheme. The closed-loop stability of the system is assessed by a rigorous analysis. A simulation study is contacted, and the obtained results confirm the validity and the effectiveness of the proposed approach.

Key words: Fokker-Planck equation; thermostatically controlled loads; aggregate load control; well-posedness, stability

## 1 Introduction

Aggregates of large populations of thermostatically controlled loads (TCLs) can be managed to offer auxiliary services, such as frequency control, load following, and energy balancing, which can contribute to maintaining the overall stability of power networks [1,2,3,4]. TCLs can also provide a means for absorbing the fluctuations of renewable energy generated by, e.g., wind turbines and solar photovoltaic plants [1,5]. Moreover, due to the fact that most of the TCLs, including space heaters, air conditioners, hot water tanks, and refrigerators, exhibit flexibilities in power demand for their operation and elasticities in terms of performance restrictions, they are considered to be one of the most important Demand Response (DR) resources that can provide such features as power peak shaving and valley filling and enable dynamic pricing schemes in the context of the Smart Grid [6,7,8,9,10,11,12,13,14,15,16,17]. Indeed, control of aggregated TCL populations is a long-time standing problem, which continues to attract much attention in the recent literature [3,18,19,20].

The aim of this paper is to develop strategies for aggregate power control of heterogeneous TCL populations based on continuum models described by partial differential equations (PDEs). The PDE aggregate model of TCL populations was originally introduced in [21]. Under the assumption of homogenous TCL populations with all the TCLs modeled by thermostat-controlled scalar stochastic differential equations, the dynamics of a population are expressed by two coupled Fokker-Planck equations describing the evolution of the distribution of the TCLs in ON and OFF states, respectively, over a rang of temperature. In about a quarter century later, this PDE aggregate model has been elaborated in the work [1] to incorporate control actions to manipulate the aggregate of TCL populations, such as the total power consumption. This approach is adopted later by, e.g., [6,14,15,22] for different DR applications. A more generic stochastic hybrid system model applicable to a wider class of responsive loads is developed in [23]. PDE aggregate models can also be derived in the framework of deterministic systems. Inspired by the concept of fluid-flow dynamics, a semi-linear transport PDE has been derived in [24] for homogeneous TCL populations. It is

Email addresses: zhengjun@swjtu.edu.cn (Jun Zheng 1 ), gabriel.laparra@polymtl.ca (Gabriel Laparra 2 ), guchuan.zhu@polymtl.ca (Guchuan Zhu ), 2 lmbuaa@gmail.com (Meng Lia ). 3

proposed later in [25] that this PDE model can be extended to describe heterogeneous TCL populations by adding a diffusive term, which indeed results in the same type of PDEs as the one developed in the framework of stochastic systems, i.e., a pair of coupled Fokker-Planck equations. Moreover, a heterogenous population with high diversity can be divided into a finite number of subsets, each of which represents a population with limited variation in its heterogeneity [26].

Another widely adopted model for describing the dynamics of TCL populations is the one originally proposed in [27]. In this model, a temperature deadband is divided into a finite number of segments, called state bins, and the evolution of the distribution of TCLs over the temperature under deadband control is presented by a linear system for which the state transactions can be derived based on Markov chains. This model has been further extended to TCLs of more generic dynamics, e.g., higher order systems, with different control mechanisms for different applications [2,3,4,17,28,29].

In this work, we consider the TCL populations for which the individual TCLs are operated under model predictive control (MPC), which is one of the most frequently adopted techniques in practice due to its ability to handle constraints, time varying processes, delays, and uncertainties, as well as disturbances. It should be noted that as with MPC schemes, TCL switching may occur arbitrarily over the time at any point inside a range of temperature. Consequently, the state-bin model based on the transaction of Markov chains is not applicable. For this reason and being inspired by the work reported in the recent literature, we concentrate on PDE-based aggregate modeling and control of large TCL populations in this paper.

The approach adopted in most of the work on PDE-based aggregate power control of TCL populations is based on the technique of early-lumping with which the control design is carried out on first discretizing the PDE models over the space [1,15,24,26,30]. There are few exceptions, such as the control scheme for output regulation of TCL populations proposed in [31], which are based on the approach of late-lumping and hence, the original PDE model is used in control design without approximations. One of the main advantages of this approach is that it can preserve the basic property of the original system, in particular the closed-loop stability, when the designed control is applied [32].

The control design presented in this work adopts the approach of late-lumping , in which the original PDE model is used. Specifically, the technique of input-output linearization (see, e.g., [33,34,35]), which is an extension of the same technique for finite dimensional nonlinear system control [36,37,38,39], is used in the design of aggregate power (the output) regulation of a TCL population governed by a pair of coupled Fokker-Planck equations. It should be pointed out that this method is greatly inspired by the work presented in [35] while being different in that the control design is carried out directly with the original Fokker-Planck equations. Consequently, a guarantee on closed-loop stability is of primary importance from both theoretical and practical perspectives. To assess the validity of the developed control scheme, an analysis on the well-posedness and the stability of the closed-loop system with nonlinear nonlocal terms is performed. It is worth noting that to the best of our knowledge, the result obtained in this paper is novel for PDE-based aggregate power control of large TCL populations. Particularly, the results on the well-posedness and the stability of Fokker-Planck equations having nonlinear nonlocal terms distributed both in the domain and on the boundary and involving integration in the denominator of nonlocal terms are novel.

The main contributions of this paper include:

- · incorporating MPC-based control schemes at the TCL level into an aggregate model of large TCL populations described by a pair of coupled Fokker-Planck equations;
- · developing an aggregated power control scheme based on the technique of input-output linearization for PDE systems;
- · assessing the well-posedness of the solutions of the closed-loop system with nonlinear nonlocal terms by the technique of iteration and the generalization of the approach developed in [40];
- · conducting a rigorous stability analysis of the closed-loop system composed of finite dimensional input-output dynamics and infinite dimensional internal dynamics.

In the rest of this paper, Section 2 presents the dynamical model of individual TCLs under MPC and the derivation of aggregate model of homogeneous and heterogeneous TCLs populations. Control design is given in Section 3. Wellposedness assessment by the technique of iteration and stability analysis are detailed in Section 4 and Section 5, respectively. A simulation study to evaluate the behavior of the closed-loop system with the developed control scheme is presented in Section 6. Some concluding remarks are provided in Section 7. Finally, details on some technical development are presented in the appendices.

Notations: If there is no ambiguity, we use the notations ∂ m u ∂x m , D u m x , or ux · · · x ︸ ︷︷ ︸ , or u ( m ) to denote the m th-order derivative

m

times

(or weak derivative) of a function u w.r.t. its argument x , where m is a nonnegative integer.

For T ∈ (0 + , ∞ ) , let Q T = ( x, x ) × (0 , T ) with Q T = [ x, x ] × [0 , T ] . Denote Q ∞ = ( x, x ) × (0 , ∞ ) and Q ∞ = [ x, x ] × [0 , ∞ ) . Let l be a nonintegral positive number. The spaces H l ([ x, x ]) with the norm | · | ( ) l ( x,x ) and H l, l 2 ( Q T ) with the norm | · | ( ) l Q T are defined in Appendix 8.1.

For T ∈ (0 + , ∞ ] , the space L ∞ ((0 , T ); L ∞ ( x, x )) consists of all strongly measurable functions u : (0 , T ) → L ∞ ( x, x ) with the norm

$$\| u \| _ { \infty } \vdots = \| u \| _ { L ^ { \infty } ( ( 0, T ) ; L ^ { \infty } ( \underline { x }, \overline { x } ) ) } \vdots = \text{vrai} \sup _ { 0 < s < T } \| u ( \cdot, s ) \| _ { L ^ { \infty } ( \underline { x }, \overline { x } ) } \vdots = \text{vrai} \sup _ { 0 < s < T } \| \text{vrai} \sup _ { \underline { x } < \overline { x } } | u ( x, s ) | < + \infty.$$

Finally, for simplicity, we denote by ‖ u ‖ p the norm ‖ u ‖ L p ( x,x ) for a function u ∈ L p ( x, x ) , where 1 ≤ p ≤ + ∞ . Particularly, we use ‖ u ‖ to denote ‖ u ‖ 2 .

## 2 Modeling Aggregate of TCL Populations

## 2.1 Dynamics and Control of Individual TCLs

We consider a large population of thermostatically controlled loads (TCLs) operated in ON/OFF mode. Denote by x i the temperature of the i -th load in a population of size N . The dynamics of x i are described by a first-order ordinary differential equation (ODE):

$$\frac { d x _ { i } } { d t } = \frac { 1 } { R _ { i } C _ { i } } \left ( x _ { i e } - x _ { i } + s _ { i } u _ { i } R _ { i } P _ { i } \right ), \, i = 1, \dots, N,$$

where x ie is the external temperature, R i is thermal resistance, C i is thermal capacitance, P i is thermal power, s i is defined as

$$s _ { i } = \begin{cases} 1, & \text{heating system;} \\ - 1, & \text{cooling system}, \end{cases}$$

and the switching signal u i :

$$u _ { i } = \begin{cases} 1, & O N ; \\ 0, & \text{OF}, \end{cases}$$

is generated by an MPC scheme specified below. Note that this hybrid dynamic model is widely adopted in the literature, while deadband-based switching is the most used control scheme in different works (see, e.g., [1,5,6,14,15,24,25,26]).

The model used in MPC design is a discrete-time system that can be derived from (1) by using the standard scheme of zeroorder hold in a sampling period T s , which takes the form:

$$x _ { i } ( k + 1 ) = A _ { i } x _ { i } ( k ) + B _ { i } u _ { i } ( k ) + E _ { i } x _ { i e } ( k ),$$

where the system parameters are given by A i = e -1 Ri C i T s , B i = s i C P i i ∫ Ts 0 e -1 Ri C i τ d τ , and E i = 1 R C i i ∫ Ts 0 e -1 Ri C i τ d τ .

Let x ref be the reference signal. Denote by e i ( k ) = x i ( k ) -x ref ( k ) the tracking error. Then for an output regulation problem corresponding to the i -th subsystem the MPC design is an optimization problem expressed by [41]:

$$\min _ { u _ { i } ( 1 ), \dots, u _ { i } ( M ) } \sum _ { k = 1 } ^ { M } \left ( \| e _ { i } ( k ) \| _ { Q _ { m p x } } ^ { 2 } + \| u _ { i } ( k ) \| _ { R _ { m p x } } ^ { 2 } \right ), \\. \quad \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \dots \nu \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \ delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta \delta 0$$

$$\text{s.l.\ } x _ { i } ( k + 1 ) = A _ { i } x _ { i } ( k ) + B _ { i } u _ { i } ( k ) + E _ { i } x _ { i e } ( k ),$$

x

i

( )

t

→

x

i

(0)

,

u

i

( )

t

→

u

i

(0)

,

t &gt;

0

,

(3c)

$$u _ { i } ( k ) \in \{ 0, 1 \},$$

$$x _ { r e f } ( k ) - \underline { \Delta x } \leq x _ { i } ( k ) \leq x _ { r e f } ( k ) + \overline { \Delta x },$$

where M is the prediction horizon, Q mpc and R mpc are, respectively, the weighting factors applied to penalize the tracking errors and the control efforts, and ∆ x and ∆ x are, respectively, the lower and upper temperature bounds. Switching sequences are determined by the solution of the above problem. At each time step of running, the controller compute the solutions until P steps, while only the first one is applied. This process is repeated at the next sampling time, constituting a moving horizon scheme. Note that in the implementation of MPC, the control sequence is applied to the continuous-time system (1) through a mechanism of zero-order hold. Consequently, the hybrid nature of TCLs with their dynamics governed by ODEs while switching between ON and OFF states will be preserved.

## 2.2 PDE Model for Heterogeneous Populations of TCLs

Let v x, t ( ) and w x, t ( ) be the distribution of loads [number of loads/ ◦ C] at temperature x and time t , over the OFF and ON states, respectively. For notational simplicity, we consider hereafter only the setting for cooling systems. The models for the setting of heating systems can be established in the same way.

First, consider the distribution of loads over the ON state. Let F x, t ( ) denote the flow of loads at the point ( x, t ) . Then, based on the model describing the microscopic dynamics of individual TCLs given in Section 2.1, we have for a population of homogeneous loads:

$$F ( x, t ) = \alpha _ { 1 } ( x ) w ( x, t ),$$

where

$$\alpha _ { 1 } ( x ) = \frac { 1 } { R C } \left ( x _ { e } - x - R P \times 1 \right ),$$

with x e , R C , , and P representing the external temperature, thermal resistance, capacitance, and power, respectively, of the population.

Let δ 0 → 1 ( x, t ) denote the added flow due to the switch of loads from OFF state to ON state at ( x, t ) , which is induced by the operation of MPCs at the TCL level. Note that as pointed out in [25], the aggregate power consumption of a heterogenous TCL population exhibits a damped response to step set-point changes. This property can be represented by the diffusivity added to the macroscopic fluid-flow TCL population model (see, e..g. [24]), which leads to the following Fokker-Planck equation:

$$\frac { \partial w } { \partial t } ( x, t ) = & \frac { \partial } { \partial x } \left ( \beta \frac { \partial w } { \partial x } ( x, t ) - \alpha _ { 1 } ( x ) w ( x, t ) \right ) + \delta _ { 0 \rightarrow 1 } ( x, t ), \, ( x, t ) \in Q _ { \infty },$$

where β is the diffusion coefficient.

Similarly, we can derive the dynamic model of the distribution of loads over the OFF state, which is given by

$$\frac { \partial v } { \partial t } ( x, t ) = & \frac { \partial } { \partial x } \left ( \beta \frac { \partial v } { \partial x } ( x, t ) - \alpha _ { 0 } ( x ) v ( x, t ) \right ) + \delta _ { 1 \rightarrow 0 } ( x, t ), \ ( x, t ) \in Q _ { \infty },$$

where

$$\alpha _ { 0 } ( x ) = \frac { 1 } { R C } \left ( x _ { e } - x \right ),$$

and δ 1 → 0 ( x, t ) denotes the added flow due to the switch of loads from ON-state to OFF-state at ( x, t ) .

Suppose that in the considered problem the operation is mass-conservative, that is the size of the population is constant. Indeed, it is reasonable to assume that no TCLs will be added to or removed from the population during a sufficiently long period of operation. Then, the flow due to the switch of loads at ( x, t ) must verify the compatibility condition δ 0 → 1 ( x, t ) =

-δ 1 → 0 ( x, t ) := δ x, t ( ) . Therefore, the dynamics of the evolution of the distribution of loads can be expressed by the following coupled Fokker-Planck equations:

$$\frac { \partial w } { \partial t } ( x, t ) = & \frac { \partial } { \partial x } \left ( \beta \frac { \partial w } { \partial x } ( x, t ) - \alpha _ { 1 } ( x ) w ( x, t ) \right ) + \delta ( x, t ), \ ( x, t ) \in Q _ { \infty },$$

$$\frac { \partial v } { \partial t } ( x, t ) = & \frac { \partial } { \partial x } \left ( \beta \frac { \partial v } { \partial x } ( x, t ) - \alpha _ { 0 } ( x ) v ( x, t ) \right ) ^ { - } \delta ( x, t ), \, ( x, t ) \in Q _ { \infty }.$$

To assure the validity of the PDE aggregate model given in (9), we impose the following boundary and initial conditions

$$\beta \frac { \partial w } { \partial x } ( x, t ) - \alpha _ { 1 } ( x ) w ( x, t ) = 0,$$

$$\beta \frac { \overset { x } { \partial w } } { \underset { \Delta.. } { \partial x } } ( \overline { x }, t ) - \alpha _ { 1 } ( \overline { x } ) w ( \overline { x }, t ) = 0,$$

$$\beta \frac { \partial v } { \partial x } ( \underline { x }, t ) - \alpha _ { 0 } ( x ) v ( \underline { x }, t ) = 0,$$

β

∂v

∂x

(

x, t

)

-

α

0

(

x v x, t

)

(

) = 0

,

(10d)

$$w ( x, 0 ) = w ^ { 0 }, \, v ( x, 0 ) = v ^ { 0 },$$

which guarantees the conservativity of mass of the problem under consideration, as stated below. Proposition 1 Let the size of the population modeled by (9) and (10) be defined as:

$$N _ { a g g } ( t ) = \int _ { \underline { x } } ^ { \overline { x } } \left ( w ( x, t ) + v ( x, t ) \right ) d x, \, t \geq 0.$$

We have then

$$N _ { a g g } ( t ) = N _ { a g g } ( 0 ) = \int _ { \underline { x } } ^ { \overline { x } } \left ( w ^ { 0 } ( x ) + v ^ { 0 } ( x ) \right ) d x, \, \forall t \geq 0.$$

$$\begin{smallmatrix} \\ \end{smallmatrix}$$

Proof The result can be obtained directly by taking integration of (9) and using the boundary conditions (10) .

Remark 1 In the formulation proposed in this work, the reference temperature is the same for all the TCLs in the population. However, it is important to note that the range around the reference which each TCL attempts to restrict to (also called the comfort zone) can be different from other TCLs depending on its performance requirement and power consumption consideration. This feature is significantly different from that in thermostat-based deadband control schemes where the turn-on and turn-off temperatures are identical for all the TCLs in the same population.

For notational simplicity, we may ignore hereafter the domain on which the system is defined and the arguments of different functions if there is no ambiguity.

## 2.3 Controlled PDE Model for Power Consumption Control

As mentioned earlier, power consumption control of a TCL population can be achieved by moving the control volume with the set-point temperature [1,24,26]. By taking the time-derivative of the set-point as the control input: u t ( ) = ˙ x ref , the dynamics of the evolution of the distribution of a population of heterogenous loads can be expressed by the following forced Fokker-Planck equations:

$$\frac { \partial w } { \partial t } ( x, t ) = \frac { \partial } { \partial x } \left ( \beta \frac { \partial w } { \partial x } ( x, t ) - ( \alpha _ { 1 } ( x ) - u ( t ) ) w ( x, t ) \right ) + \delta ( x, t ), \ ( x, t ) \in Q _ { \infty },$$

$$\frac { \partial v } { \partial t } ( x, t ) = & \frac { \partial } { \partial x } \left ( \beta \frac { \partial v } { \partial x } ( x, t ) - ( \alpha _ { 0 } ( x ) - u ( t ) ) v ( x, t ) \right ) ^ { \cdot } - \delta ( x, t ), \ ( x, t ) \in Q _ { \infty },$$

with the following boundary and initial conditions:

$$\beta \frac { \partial w } { \partial x } ( \underline { x }, t ) - ( \alpha _ { 1 } ( \underline { x } ) - u ( t ) ) w ( \underline { x }, t ) = 0,$$

$$\beta \frac { \partial w } { \partial x } ( \overline { x }, t ) - ( \alpha _ { 1 } ( \overline { x } ) - u ( t ) ) w ( \overline { x }, t ) = 0,$$

$$\beta \frac { \tilde { \partial v } } { \partial x } ( \underline { x }, t ) - ( \alpha _ { 0 } ( \underline { x } ) - u ( t ) ) v ( \underline { x }, t ) = 0,$$

$$\beta \frac { \overset { x } { \partial v } } { \partial x } ( \overline { x }, t ) - ( \alpha _ { 0 } ( \overline { x } ) - u ( t ) ) v ( \overline { x }, t ) = 0, \\ \overset { x } { \sim } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda } \hat { \lambda }$$

$$w ( x, 0 ) = w ^ { 0 }, \, v ( x, 0 ) = v ^ { 0 }.$$

For the validity of the PDE aggregate model given in (11) and (12), it is expected that the closed-loop system preserves the property of mass conservativity, which is guaranteed by Proposition 2 given later in Section 4.

## 3 Control Design

For the purpose of controlling the power consumption of the whole population of TCLs, we take the weighted total power consumption as the system output:

$$y ( t ) = \frac { P } { \eta } \int _ { x _ { 2 } } ^ { \overline { x } } ( a x + b ) w ( x, t ) d x,$$

glyph[negationslash]

where a and b are constants with a = 0 . The purpose of introducing the weighting function ax + b in system output defined above is to guarantee that the input-output dynamics of the system are well-defined in terms of characteristic index that is an analogue of the concept of relative degree for finite-dimensional systems [33,34]. The reason to use a specific weighting function in (13) is mainly due to the fact that in general it is very difficult, if not impossible, to conduct well-posedness and stability analysis for nonlinear PDEs with generic nonlinear terms. Note that the weighting function is not unique, and it can be chosen and tuned in control design and implementation depending on control objectives and performance requirements.

We consider in the present work a set-point control problem where the objective is to drive the power consumption of the whole population to a desired level y d representing the demand from, e.g., energy suppliers.

Denoting by e t ( ) = y t ( ) -y d the regulation error, we can derive from (11a) and (13) that

$$\frac { \mathrm d e } { \mathrm d t } = & \frac { P } { \eta } \int _ { \underline { x } } ^ { \overline { x } } ( a x + b ) \frac { \partial w } { \partial t } \mathrm d x \\ = & \frac { P } { \eta } \int _ { \underline { x } } ^ { \overline { x } } ( a x + b ) \frac { \partial } { \partial x } \left ( \beta \frac { \partial w } { \partial x } - ( \alpha _ { 1 } ( x ) - u ( t ) ) w \right ) \mathrm d x + \frac { P } { \eta } \int _ { \underline { x } } ^ { \overline { x } } ( a x + b ) \delta ( x, t ) \mathrm d x.$$

Performing integration by parts and applying the boundary conditions (12), we get

$$\frac { d e } { d t } = - \frac { P } { \eta } \int _ { \underline { x } } ^ { \overline { x } } a \left ( \beta \frac { \partial w } { \partial x } - ( \alpha _ { 1 } ( x ) - u ( t ) ) w \right ) \mathrm d x + \frac { P } { \eta } \int _ { \underline { x } } ^ { \overline { x } } ( a x + b ) \delta ( x, t ) \mathrm d x.$$

Let

$$u ( t ) = - \, \frac { \int _ { \underline { x } } ^ { \overline { x } } \left ( \beta \frac { \partial w } { \partial x } - \alpha _ { 1 } ( x ) w \right ) \text{d} x + \frac { \eta } { a P } \phi ( t ) } { \int _ { \underline { x } } ^ { \overline { x } } w \text{d} x } _ { \, e ^ { \overline { x } } } \, _ { \, e ^ { \overline { x } } } \, _ { \, \dots }$$

$$= & - \frac { \beta ( w ( \overline { x }, t ) - w ( \underline { x }, t ) ) - \int _ { \underline { x } } ^ { \overline { x } } \alpha _ { 1 } ( x ) w \mathrm d x + \frac { \eta } { \ a P } \phi ( t ) } { \int _ { \underline { x } } ^ { \overline { x } } w \mathrm d x },$$

$$\Gamma ( t ) = & \frac { P } { \eta } \int _ { \underline { x } } ^ { \overline { x } } ( a x + b ) \delta ( x, t ) d x.$$

Then the regulation error dynamics become:

$$\frac { d e } { d t } = \phi ( t ) + \Gamma ( t ), \ e ( 0 ) = e _ { 0 },$$

where φ t ( ) is an auxiliary control, and e 0 is the initial regulation error.

Remark 2 It can be seen from (14) that the characteristic index of the input-output dynamics of the system is 1 if at any time the power consumption of the whole population is not null. This constraint is a nature property for a TCL population of sufficiently large size and can be guaranteed by a simple operation in practice.

Remark 3 Note that the term Γ( ) t represents the effect of switching induced by the MPCs at the TCL level. An advantage of treating Γ( ) t as a disturbance in the proposed aggregate power control scheme is that the TCLs do not need to signal the instantaneous switching operations, which will allow greatly simplifying the implementation and considerably relaxing the performance requirement on communication systems. Moreover, the developed control algorithm is a partial state-feedback scheme depending only on w x, t ( ) , which is exactly the same information needed for computing the system output [1,15,24,26, ? ? , ]. Finally, communicating the current operational state and the output (the temperature) of the TCLs can be enabled by, e.g., advanced metering infrastructures in the context of the Smart Grid [1].

The control synthesis will then be completed by finding a φ t ( ) that can robustly stabilize the regulation error dynamics in the presence of disturbances. For this aim, we apply a standard and very effective method, called nonlinear damping (see, e.g., [36,37]). Specifically, it is reasonable to assume that Γ( ) t is uniformly bounded, i.e., | Γ( ) t | ≤ Γ ∞ for all t &gt; 0 with Γ ∞ a positive constant, because at any moment only a limited number of TCLs change their operational state. Then, it is obvious that the control

$$\phi ( t ) = - k _ { 0 } e ( t )$$

with any constant gain k 0 &gt; 0 will globally exponentially stabilize (16) when Γ( ) t ≡ 0 for all t &gt; 0 . Furthermore, V = e 2 is a control Lyapunov function for the disturbance-free counterpart of the system (16). Therefore, due to [36, Lemma 14.1, p. 589], the control given in (17) guarantees that the trajectory of the system (16) is globally uniformly bounded. Moreover, the amplitude of regulation error e t ( ) can be rendered arbitrarily small if the control gain k 0 is sufficiently high. This property is explicitly described in Theorem 3 in Section 5.

The schematic diagram of the considered aggregate power control system for TCL populations is shown in Fig. 1, which has the same structure as most of the systems developed in the literature (see, e.g., [1,3,4,5,14,15,24,26,31], etc.), except that MPC is used at the level of each TCL in a fully distributed manner. The computation of PDE control law requires only the value of the temperature of the TCLs and their current operational state. Whilst, the reference for MPCs is generated from the PDE control signal as x ref ( ) = t x ref (0) + ∫ t 0 ˙ x ref ( τ ) d τ = x ref (0) + ∫ t 0 u τ ( ) d τ .

Fig. 1. Schematic diagram of aggregate power control of a TCL population operated by MPC at the TCL level.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00302v2-figure-cc601dc0-ad89-4227-acda-584330235abd.jpg)

## 4 Well-posedness of the Closed-loop System

In this section, we establish the well-posedness of the closed-loop system composed of (11), (12), (13), (14), (15), (16), and (17). To this aim, we note that α , α 1 0 ∈ C 1 ([ x, x ]) and P, β, η are positive constants, and assume that

$$w ^ { 0 }, v ^ { 0 } \in \mathcal { H } ^ { 2 + \theta } ( [ x, \overline { x } ] ), \delta \in L ^ { \infty } ( ( 0, \infty ) ; L ^ { \infty } ( x, \overline { x } ) ),$$

with

$$& \left | \int _ { \underline { s } } ^ { \overline { x } } w ^ { 0 } ( x ) d x \right | > 0, & & ( 1 9 a ) \\ & \left | \int _ { \underline { s } } ^ { \overline { x } } w ^ { 0 } ( x ) d x + \int _ { 0 } ^ { t } \int _ { \underline { s } } ^ { \overline { x } } \delta ( x, s ) d x d s \right | \geq \delta _ { 0 } > 0, \forall t \in [ 0, + \infty ), & & ( 1 9 b )$$

where θ ∈ (0 1) , , δ 0 &gt; 0 are constants. Note that the conditions given in (19a) and (19b) mean that at any time the TCLs are not all in OFF-state.

Theorem 1 Consider the closed-loop system. For any T &gt; 0 , there exists ( w,v ) ∈ H 2+ θ, 1+ θ 2 ( Q T ) × H 2+ θ, 1+ θ 2 ( Q T ) that satisfies (11) a.e. in Q T .

It should be pointed out that with the proposed control scheme, the closed-loop system is governed by parabolic PDEs with nonlinear nonlocal terms, which brings the following two obstacles for well-posedness analysis:

- (i) Due to the appearance of a nonlinear nonlocal term in the form of w x ( x, t ) ∫ α 1 ( x w x, t ) ( ) d x/ ∫ w x, t ( ) d x , one cannot deal with the original PDEs directly and apply the Lyapunov method to obtain the a prior estimates of the solutions. Consequently, the Galerkin method cannot be directly used for well-posedness analysis.
- (ii) The corresponding PDEs may have a singularity due to the appearance of the term ∫ w x, t ( ) d x in the denominator of the nonlocal term. Therefore, the classical Leray-Schauder fixed-point theorem cannot be applied directly in well-posedness analysis.

To overcome the difficulties brought by the nonlinear nonlocal character, we use in the present work the technique of iteration for well-posedness analysis.

Specifically, for any T &gt; 0 , let { δ n } be a sequence of functions on Q T satisfying:

- · δ n ( x, t ) is H¨lder continuous in o x with exponent θ ;
- · δ n ( x, t ) is C 1 -continuous in t with ∣ ∣ ∣ ∣ ∂δ n ∂t ( x, t ) ∣ ∣ ∣ ∣ ≤ ̂ δ in Q T , where ̂ δ is a positive constant;
- · δ n ( x, t ) → δ x, t ( ) a.e. in Q T , as n →∞ ;
- · ‖ δ n ‖ ∞ ≤ 1 + ‖ ‖ δ ∞ , ∀ n ≥ 1 ;
- · for any n ≥ 1 and t ∈ [0 , T ] , it holds:

$$\left | \int _ { \underline { x } } ^ { \overline { x } } w ^ { 0 } ( x ) d x + \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, s ) d x d s \right | \geq \frac { \delta _ { 0 } } { 2 } > 0.$$

Consider the following two sets of iterating equations in Q T :

$$\frac { \partial w _ { n } } { \partial t } ( x, t ) = - \frac { \partial } { \partial x } \left ( \alpha _ { 1 } ( x ) w _ { n } ( x, t ) \right ) + \beta \frac { \partial ^ { 2 } w _ { n } } { \partial x ^ { 2 } } ( x, t ) - \frac { \partial w _ { n } ( x, t ) } { \partial x } G _ { w _ { n - 1 } } ( t ) + \delta _ { n } ( x, t ),$$

$$\partial t ^ { \stackrel { \sim } { \nu \nu } } \quad \partial x ^ { \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { < } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \stackrel { \sim } { \nu \nu } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{\ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text { } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text
\text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ }
\text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{\ }\text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{  \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{   \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{\} \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{} \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{\{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \end{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \text{ } \begin{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{array} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{array} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrices} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end{matrix} { c c } \end$$

$$\frac { \partial v _ { n } } { \partial t } ( x, t ) & = - \frac { \partial } { \partial x } \left ( \alpha _ { 1 } ( x ) v _ { n } ( x, t ) \right ) + \beta \frac { \partial ^ { 2 } v _ { n } } { \partial x ^ { 2 } } ( x, t ) - \frac { \partial v _ { n } ( x, t ) } { \partial x } G _ { w _ { n } } ( t ) + \delta _ { n } ( x, t ), \\ \beta \frac { \partial v _ { n } } { \partial x } ( \overline { x }, t ) & - \alpha _ { 0 } ( \overline { x } ) v _ { n } ( \overline { x }, t ) - v _ { n } ( \overline { x }, t ) G _ { w _ { n } } ( t ) = 0,$$

and

$$& \beta \frac { \partial v _ { n } } { \partial x } ( \underline { x }, t ) - \alpha _ { 0 } ( \underline { x } ) v _ { n } ( \underline { x }, t ) - v _ { n } ( \underline { x }, t ) G _ { w _ { n } } ( t ) = 0, \\ & v _ { n } ( x, 0 ) = v ^ { 0 } ( x ),$$

where

$$G _ { w _ { n - 1 } } ( t ) = & \frac { \int _ { \underline { x } } ^ { \overline { x } } \left ( \beta \frac { \partial w _ { n - 1 } } { \partial x } - \alpha _ { 1 } w _ { n - 1 } \right ) \mathrm d x } { \int _ { \underline { x } } ^ { \overline { x } } w _ { n - 1 } \mathrm d x } - \frac { k _ { 0 } } { a } \frac { \int _ { \underline { x } } ^ { \overline { x } } ( a x + b ) w _ { n - 1 } \mathrm d x - \frac { \eta } { P } y _ { d } } { \int _ { \underline { x } } ^ { \overline { x } } w _ { n - 1 } \mathrm d x }, \, n \geq 1,$$

with w 0 ( x, t ) := w 0 ( x ) .

We claim that for any n , there exists a unique solution w ,v n n ∈ H 2+ θ, 1+ θ 2 ( Q T ) to the iterating equations (21) and (22) respectively.

Lemma 1 For any n , (21) and (22) admit a unique solution w n ∈ H 2+ θ, 1+ θ 2 ( Q T ) and v n ∈ H 2+ θ, 1+ θ 2 ( Q T ) respectively, satisfying for any t ∈ [0 , T ] :

$$\int _ { \underline { x } } ^ { \overline { x } } w _ { n } ( x, t ) d x & = \int _ { \underline { x } } ^ { \overline { x } } w ^ { 0 } ( x ) d x + \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, s ) d x d s \geq \frac { \delta _ { 0 } } { 2 }, \\ \int _ { \underline { x } } ^ { \overline { x } } v _ { n } ( x, t ) d x & = \int _ { \underline { x } } ^ { \overline { x } } v ^ { 0 } ( x ) d x - \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, s ) d x d s \geq \frac { \delta _ { 0 } } { 2 },$$

and

$$\| w _ { n } ( \cdot, t ) \| _ { 1 } \leq \| w ^ { 0 } \| _ { 1 } + \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, s ) s g n ( w _ { n } ) d x d s,$$

$$r ^ { t } \ r ^ { \overline { x } }$$

$$\| v _ { n } ( \cdot \, t ) \| _ { 1 } \leq \| v ^ { 0 } \| _ { 1 } - \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, s ) s g n ( v _ { n } ) d x d s.$$

Moreover, ( w ,v n n ) has the following uniform estimates:

$$\max _ { \overline { Q } _ { T } } | w _ { n } | + \max _ { \overline { Q } _ { T } } | w _ { n z } | + | w _ { n } | _ { Q _ { T } } ^ { ( 2 + \theta ) } \leq C _ { 1 }, \\ \dots \.$$

$$\max _ { \overline { Q } _ { T } } | v _ { n } | + \max _ { \overline { Q } _ { T } } | v _ { n z } | + | v _ { n } | _ { Q _ { T } } ^ { ( 2 + \theta ) } \leq C _ { 2 },$$

where C 1 and C 2 are positive constants independent of n .

Proof The proof is proceeded in four steps as detailed in Appendix 8.2.

Remark 4 The main mathematical framework used in the proof of Lemma 1 is as follows:

- (i) The existence of a unique solution ( w ,v n n ) ∈ H 2+ θ, 1+ θ 2 ( Q T ) × H 2+ θ, 1+ θ 2 ( Q T ) to the iterating equations (21) and (22) is proved by induction combining with the well-posedness and regularity theory for quasilinear parabolic equations in [40, Chap. V].

(ii) The L 1 -estimates of w n and v n in (23) are established by utilizing appropriate test functions and some basic integrating techniques. Note that since ‖ δ n sgn ( w n ) ‖ ∞ ≤ 1+ ‖ ‖ δ ∞ and ‖ δ n sgn ( v n ) ‖ ∞ ≤ 1+ ‖ ‖ δ ∞ , (23) implies the uniform boundedness of ‖ w n ‖ 1 and ‖ v n ‖ 1 , respectively.

(iii) The uniform boundedness of ‖ w n ‖ 1 and ‖ v n ‖ 1 and the continuity of w n and v n guarantee the uniform boundedness of w n and v n , which, along with the well-posedness and regularity theory for quasilinear parabolic equations in [40, Chap. V], guarantees the uniform estimates of w n and v n in (24) .

Based on Lemma 1, we can establish the existence and the regularity of a solution ( w,v ) to the coupled equations (11) in closed loop formed by (11), (12), (13), (14), (15), (16), and (17).

Proof (Proof of Theorem 1) By Lemma 1, ( w ,v n n ) is bounded in H 2+ θ, 1+ θ 2 ( Q T ) × H 2+ θ, 1+ θ 2 ( Q T ) . Then, by choosing a subsequence, there is a pair of functions ( w,v ) ∈ H 2+ θ, 1+ θ 2 ( Q T ) × H 2+ θ, 1+ θ 2 ( Q T ) satisfying ( w ,v n n ) → ( w,v ) in H 2+ θ, 1+ θ 2 ( Q T ) × H 2+ θ, 1+ θ 2 ( Q T ) as n → ∞ . Noting that δ n → δ a.e. in Q T , we conclude that ( w,v ) satisfies (11) a.e. in Q T .

It should be noticed that Theorem 1 can still not guarantee the uniqueness of a solution of the closed-loop system, because the solution is obtained by the technique of iteration and approximation. Moreover, in general it is not trivial to establish the uniqueness of a solution to an equation with strong nonlinearities. Nevertheless, the uniqueness of the solution of the closedloop system can be guaranteed in a weak sense as stated below.

Theorem 2 Consider the closed-loop system. Let ( w , v 1 1 ) , ( w , v 2 2 ) ∈ H 2+ θ, 1+ θ 2 ( Q T ) ×H 2+ θ, 1+ θ 2 ( Q T ) be two solutions of (11) in Q T . If w 1 ( x, t ) -w 1 ( x, t ) = w 2 ( x, t ) -w 2 ( x, t ) in [0 , T ] , then ( w , v 1 1 ) = ( w , v 2 2 ) in Q T . Proof See Appendix 8.3.

## 5 Stability Analysis

We assess first the stability of the error dynamics (16) with the control given in (17).

Theorem 3 The regulation error e t ( ) is determined by

$$e ( t ) = e ( 0 ) e ^ { - k _ { 0 } t } + \int _ { 0 } ^ { t } \Gamma ( t ) e ^ { - k _ { 0 } ( t - s ) } d s.$$

Furthermore, if | Γ( ) t | ≤ Γ ∞ with a positive constant Γ ∞ for all t &gt; 0 , then

$$| e ( t ) | \leq | e ( 0 ) | e ^ { - k _ { 0 } t } + \frac { \Gamma _ { \infty } } { k _ { 0 } } \left ( 1 - e ^ { - k _ { 0 } t } \right ), \forall t > 0.$$

Proof By (16) and (17) , it follows

$$\frac { d e } { d t } = - k _ { 0 } e ( t ) + \Gamma ( t ), \ e ( 0 ) = e _ { 0 }.$$

Solving this first-order linear ordinary differential equation, we can obtain the result given in (25) .

To assure the closed-loop stability, we need to prove that the internal dynamics composed of (11), (12), (13), and (14) are stable provided Theorem 3 holds true. Towards this aim, we need the properties stated in the following 2 propositions. The first one is the conservativity of mass of the considered problem given below.

Proposition 2 Consider the closed-loop system. For any T &gt; 0 , let ( w,v ) ∈ H 2+ θ, 1+ θ 2 ( Q T ) ×H 2+ θ, 1+ θ 2 ( Q T ) satisfy (11) a.e. in Q T . Then for any t ∈ [0 , T ] , it holds:

$$\begin{array} { l l } ( i ) & \int _ { \underline { x } } ^ { \overline { x } } w ( x, t ) d x = \int _ { \underline { x } } ^ { \overline { x } } w ^ { 0 } ( x ) d x + \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta ( x, s ) d x d s ; \\ ( i i ) & \int _ { \underline { x } } ^ { \overline { x } } v ( x, t ) d x = \int _ { \underline { x } } ^ { \overline { x } } v ^ { 0 } ( x ) d x - \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta ( x, s ) d x d s. \end{array}$$

Moreover, it holds:

$$N _ { a g g } ( t ) = \int _ { \underline { x } } ^ { \overline { x } } ( w ( x, t ) + v ( x, t ) ) d x = \int _ { \underline { x } } ^ { \overline { x } } ( w ^ { 0 } ( x ) + v ^ { 0 } ( x ) ) d x.$$

Proof The result can be obtained directly by taking integration of (11) and using the boundary conditions (12) .

In order to establish the stability of the closed-loop system, we need the following property on the solutions of (11), (12), (13), (14), (15), (16), and (17).

Proposition 3 Consider the closed-loop system. For any T &gt; 0 , let ( w,v ) ∈ H 2+ θ, 1+ θ 2 ( Q T ) ×H 2+ θ, 1+ θ 2 ( Q T ) be a solution of (11) in Q T . Then for any t ∈ [0 , T ] , the following estimates hold

$$( i ) \ \| w ( \cdot, t ) \| _ { 1 } & \leq \| w ^ { 0 } \| _ { 1 } + \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta ( x, s ) s g n ( w ) d x d s ; \\ ( i i ) \ \| v ( \cdot, t ) \| _ { 1 } & \leq \| v ^ { 0 } \| _ { 1 } - \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta ( x, s ) s g n ( v ) d x d s,$$

, which leads to the desired result.

where sgn ( η ) = 1 if η &gt; 0 , sgn ( η ) = -1 if η &lt; 0 , and sgn ( η ) = 0 if η = 0 . Proof By Lemma 1 and the regularity of w n and v n , it suffices to let n →∞ in (23)

As it is guaranteed by Proposition 2 that the closed-loop system is also conservative in terms of the total number of TCLs in the population, there must be a positive constant M such that ∣ ∣ ∫∫ ∣ Q δ x, t ( ) d x t d ∣ ∣ ∣ &lt; M for any Lebesgue measurable set Q ⊂ Q ∞ . We have then the following result regarding the stability of the internal dynamics in closed loop.

Theorem 4 Consider the internal dynamics of the system. Let ( w,v ) ∈ H 2+ θ, 1+ θ 2 ( Q ∞ ) × H 2+ θ, 1+ θ 2 ( Q ∞ ) satisfy (11) a.e. in Q ∞ . Assume that there exists a positive constant M such that ∣ ∣ ∫∫ ∣ Q δ x, t ( ) d x t d ∣ ∣ ∣ &lt; M for any Lebesgue measurable set Q ⊂ Q ∞ . Then for any t ∈ [0 , + ∞ ) , it holds:

(i) ‖ w ,t ( · ) ‖ 1 ≤ ‖ w 0 ‖ 1 +2 M &lt; + ∞ ; (ii) ‖ v ( · , t ) ‖ 1 ≤ ‖ v 0 ‖ 1 +2 M &lt; + ∞ .

Proof We asses only the claim (i). The proof of the claim (ii) can follow the same line. Indeed, for t = 0 , we have yet ‖ w , ( · 0) ‖ 1 = ‖ w 0 ‖ 1 ≤ ‖ w 0 ‖ 1 + 2 M &lt; + ∞ . For any t &gt; 0 , let Q + = { ( x, s ) ∈ Q t ; w x, s ( ) &gt; 0 } and Q -= { ( x, s ) ∈ Q t ; w x, s ( ) &lt; 0 } . Then by assumptions and Proposition 3, we have

$$\chi _ { t }, w ( x, s ) & \prec \mathfrak { J }. \, \i n e n { v y a s u m p o n s u n a r r o p o s u n o } { J, w e n a v e } \\ \| w ( \cdot, t ) \| _ { 1 } \leq & \| w ^ { 0 } \| _ { 1 } + \int _ { 0 } ^ { t } \int _ { x } ^ { \overline { x } } \delta ( x, s ) s g n ( w ) d x d s \\ = & \| w ^ { 0 } \| _ { 1 } + \iint _ { Q ^ { + } } \delta ( x, s ) d x d s - \iint _ { Q ^ { - } } \delta ( x, s ) d x d s \\ \leq & \| w ^ { 0 } \| _ { 1 } + \left | \iint _ { Q ^ { + } } \delta ( x, s ) d x d s \right | + \left | \iint _ { Q ^ { - } } \delta ( x, s ) d x d s \right | \\ \leq & \| w ^ { 0 } \| _ { 1 } + 2 M < + \infty.$$

## 6 Simulation Study

To validate the developed control scheme and illustrate the performance of the system, we conduct a simulation study on a benchmark problem, corresponding to a population of air-conditioners, which has been used in the validation of different solutions in the literature, e.g., [1,25,31,42]. The parameters of the system and the controllers are listed in Table 1. The variation of thermal capacitance creates a heterogenous TCL population.

The setting used in the simulation is a population of 1000 appliances controlled by MPC schemes. The simulation is performed in a normalized time scale with respect to the nominal time constant of the appliances represented by RC (= 20h). The results related to power consumptions are also presented in quantities normalized by the maximal total power consumption of the population. Note that taking b = -20 in the weighting function used in the output function defined in (13) has the effect of moving the data around the center of the operation range. Our experiment has shown that this choice would facilitate the controller tuning.

First, it is known that due to the phenomenon of cold load pickup [43], a step variation of temperature reference in a large population of TCLs may induce high power demand peaks. To avoid this problem, we have adopted a decentralized desynchronization mechanism that combines adding independently a weighting coefficient and a delay to the reference signal at the level of MPC for individual TCLs [44]. To illustrate the effect of this desynchronized MPC scheme, we tested the response of the population to a step change of the reference from 20.5 ◦ C to 19 ◦ C occurred at t = 15 . The MPCs are with a weighting coefficient R mpc from a normally distributed random variable restricted in the interval [6, 14] and an additional delay uniformly distributed in [0, 5] to the reference signal. It can be seen from Fig. 2 that the desynchronized MPC scheme can considerably reduce the peak power demand.

The second part of simulation study is to determine the diffusion coefficient, β , of the Fokker-Planck equations based on the experiments. The diffusivity is estimated by applying the algorithm developed in [42] to the data generated by numerical

Table 1 System parameters

| Parameters   | Description [Unit]                            | Value      |
|--------------|-----------------------------------------------|------------|
| R            | Thermal resistance [ ◦ C/kW]                  | 2          |
| C            | Thermal capacitance [kWh/ ◦ C]                | N (10 , 3) |
| P            | Thermal power [kW]                            | 14         |
| x ie         | Ambient temperature [ ◦ C]                    | 32         |
| η            | Coefficient of performance                    | 2.5        |
| β            | Diffusivity                                   | 0.1        |
| Q mpc        | Weighting coefficient on tracking error e mpc | 100        |
| R mpc        | Weighting coefficient on control input u mpc  | N (10 , 2) |
| M            | Prediction horizon                            | 5          |
| ∆ x , ∆ x    | Temperature deadband width [ ◦ C]             | 0.5        |
| a            | Weighting function coefficient                | -1         |
| b            | Weighting function coefficient                | -20        |
| k 0          | Control gain                                  | 5600       |

Fig. 2. Comparison of aggregate power of the population controlled by standard MPC and desynchronized MPC schemes at the TCL level.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00302v2-figure-fd14891f-e527-44f4-97b3-685eb49fb199.jpg)

simulations, which gives a value β = 0 1 . . To validate this result, we proceed by first generating the data δ x, t ( ) for the population of TCLs specified in Table 1 with MPC at the TCL level. The numerical solutions of the Fokker-Planck equations with the generated data δ x, t ( ) as the inputs corresponding to different values of β are obtained by the forward Euler method [45]. The solution surfaces of the Fokker-Planck equations, w x, t ( ) and v x, t ( ) , with β = 0 1 . are shown in Fig. 3(a). We compute then the normalized mean error in its absolute value of the number of the TCLs in ON and OFF states between the generated data and that given by the solution of the Fokker-Planck equations with different diffusion coefficient varying from 0 to 0.2. As can be seen from Fig. 3(b) that β = 0 1 . is an adequate estimate. It can also be observed that the errors are not sensitive to diffusion coefficient if it is smaller than 0.1. This is in accordance with the results reported in the literature (see, e.g., [42]).

In the last part of simulation study, we test the developed aggregate power control algorithm with a case where the desired aggregate power changes from 0.41 to 0.47 at t = 15 (all in the normalized scale), which corresponds to a variation of the reference temperature from 20.5 ◦ C to 19 ◦ C. [1,24,25]. The desired aggregate power and the actual power consumption of the population are shown in Fig. 4(a). The control signal and the regulation errors are depicted in Fig. 4(b) and Fig. 4(c), respectively. The results show that the developed control scheme performances well with regulation errors lower than 2.5% of the maximal total power consumption at the steady regime. In order to avoid the chartering induced by fast reference variations, the reference sent to the TCLs is the average of the values over the last 10 iterations. The obtained temperature reference is shown in Fig. 5(a). Figure 5(b) illustrates the temperature evolution of 100 randomly selected TCLs. It can be seen that all the TCLs follow well the reference while respecting the temperature deadband. Finally, it is worth noting that in the simulation we used a population containing 1000 TCLs. As the bigger is the size of the population, the more accurate is the continuum model, along with the fact that the switching of the TCLs is considered unknown as specified in (16), the simulation results conformed the robustness of the developed control scheme.

Fig. 3. Experimental validation of the estimate of diffusion coefficient: (a) solution surfaces of the distribution of TCLs in ON and OFF states: w x, t ( ) and v x, t ( ) ; (b) normalized mean errors of the number of TCLs in ON and OFF states between the measured value and the one obtained by the solution of Fokker-Planck equations.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00302v2-figure-dd8c1715-c897-4b7a-b207-34a7a5fb3a30.jpg)

0.6

0.2

0.2

10

0.6

(a)

(b)

Fig. 4. Dynamics of the TCL population: (a) aggregate power consumptions; (b) control signal; (c) regulation error.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00302v2-figure-33ae9522-0df0-430d-a810-74c2b744afb1.jpg)

## 7 Conclusion

In this work, we have developed an input-output linearization-based scheme for aggregate power control of heterogenous thermal thermostatically controlled load (TCL) populations governed by a pair of Fokker-Planck equations coupled via indomain actions representing the operation of model predictive control (MPC) at the level of individual TCLs. Well-posedness and stability analysis on the closed-loop system is conducted. A simulation study is carried out with a benchmark setting, and the results illustrated the validity and the performance of the proposed method.

It should be noted that the present work considered only a set-point regulation problem. While as the proposed linearization control leads to a system with finite-dimensional linear input-output dynamics, it is straightforward to extend this scheme to aggregate power tracking by using the well-established techniques for tracking control of finite-dimensional systems [36,37,38,39]. Finally, the control algorithm developed in this work can be applied to TCL populations operated by other control schemes at the level of individual TCLs, in particular the variations of deadband control. The main difference will be the way with

Fig. 5. Temperature evolution: (a) reference temperature; (b) temperature variations of 100 samples of TCLs.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00302v2-figure-d01ae793-7ca1-4281-9712-bc203860e8e5.jpg)

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00302v2-figure-bb154cf2-99b8-4977-865e-e7624f3a6fa5.jpg)

which the PDEs are coupled, which will lead to settings with different boundary conditions and in-domain inputs. Hence, the well-posedness and the stability of the closed-loop systems have to be assessed, which can also be conducted following the approach used in this work.

## Acknowledgements

The authors thank the anonymous reviewers for their constructive comments, which helped to improve the quality of this work and the presentation of the paper.

## 8 Appendix

$$8. 1 \ N o t a t i o n s \, o n \, \mathcal { H } ^ { l } ( [ \underline { x }, \overline { x } ] ) \, a n d \, \mathcal { H } ^ { l, l / 2 } ( \overline { Q } _ { T } )$$

We introduce first the following definitions of function spaces given in [40, Chapter I], which are used in the subsequent development. Let T &gt; 0 and l be a nonintegral positive number. H l ([ x, x ]) denotes the Banach space whose elements are continuous functions u x ( ) in [ x, x ] , having continuous derivatives up to order [ l ] inclusively and a finite value for the quantity

$$\begin{matrix} 1 \\ < \\ < \end{matrix}$$

$$& \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \\ & | u | _ { ( x, \overline { \lambda } ) } ^ { ( l ) } = \langle u \rangle _ { ( x, \overline { \lambda } ) } ^ { ( l ) } + \sum _ { j = 0 } ^ { [ l ] } \langle u \rangle _ { ( x, \overline { \lambda } ) } ^ { ( j ) }, \\ & \langle u \rangle _ { ( x, \overline { \lambda } ) } ^ { ( 0 ) } = | u | _ { ( x, \overline { \lambda } ) } ^ { ( 0 ) } = \sup _ { ( x, \overline { \lambda } ) } | u |, \\ & \langle u \rangle _ { ( x, \overline { \lambda } ) } ^ { ( j ) } = \sum _ { ( g ) } | D _ { x } ^ { j } u | _ { ( x, \overline { \lambda } ) } ^ { ( 0 ) }, \langle u \rangle _ { ( x, \overline { \lambda } ) } ^ { ( l ) } = \sum _ { ( [ l ] ) } | D _ { x } ^ { [ l ] } u | _ { ( x, \overline { \lambda } ) } ^ { ( l - [ l ] ) }, \\ & \quad \quad \quad \..$$

where (26) defines the norm | u | ( ) l ( x,x ) in H l ([ x, x ]) .

H l,l/ 2 ( Q T ) denotes the Banach space of functions u x, t ( ) that are continuous in Q T , together with all derivatives of the form D D r t s x for 2 r + s &lt; l , and have a finite norm

$$| u | _ { Q _ { r } } ^ { ( l ) } = \langle u \rangle _ { Q _ { r } } ^ { ( l ) } + \sum _ { j = 0 } ^ { [ l ] } ( u ) _ { Q _ { r } } ^ { ( j ) },$$

with

$$\langle u \rangle _ { Q _ { T } } ^ { ( 0 ) } = & | u | _ { Q _ { T } } ^ { ( 0 ) } = \max _ { \overline { Q } _ { T } } | u |, \\ \langle u \rangle _ { Q _ { T } } ^ { ( j ) } = & \sum _ { ( 2 r + s = j ) } | D _ { t } ^ { r } D _ { x } ^ { s } u | _ { Q _ { T } } ^ { ( 0 ) },$$

$$\langle u \rangle _ { Q _ { T } } ^ { ( l ) } = & \langle u \rangle _ { x, Q _ { T } } ^ { ( l ) } + \langle u \rangle _ { t, Q _ { T } } ^ { ( l / 2 ) }, \\ \langle u \rangle _ { x, Q _ { T } } ^ { ( l ) } = & \sum _ { ( 2 r + s = [ l ] ) } \langle D _ { t } ^ { r } D _ { x } ^ { s } u \rangle _ { x, Q _ { T } } ^ { ( l - [ l ] ) }, \\ \langle u \rangle _ { t, Q _ { T } } ^ { ( l / 2 ) } = & \sum _ { 0 < l - 2 r - s < 2 } \langle D _ { t } ^ { r } D _ { x } ^ { s } u \rangle _ { t, Q _ { T } } ^ { ( \frac { l - 2 r - s } { 2 r } ) }, \\, & \quad ( \alpha ) } | u ( x, t ) - u ( 3 )$$

$$\cdot$$

$$\langle u \rangle _ { t, Q _ { T } } ^ { ( \alpha ) } = \sup _ { \substack { ( x, t ), ( x, t ^ { \prime } ) \in \overline { Q } _ { T }, | t - t ^ { \prime } | \leq \rho _ { 0 } \\ \wedge \cdot \cdot \cdot \cdot } } \frac { | u ( x, t ) - u ( x, t ^ { \prime } ) | } { | t - t ^ { \prime } | ^ { \tau } },$$

$$\langle u \rangle _ { x, Q _ { T } } ^ { ( \tau ) } = \sup _ { t \in ( 0, T ) } \sup \rho ^ { - \tau } \text{osc} _ { x } \left \{ u ( x, t ) ; \Omega _ { \rho } ^ { i } \right \}.$$

$$\langle u \rangle _ { t, Q _ { T } } ^ { ( \alpha ) } & = \sup _ { ( x, t ), ( x, t ^ { \prime } ) \in \overline { Q } _ { T }, | t - t ^ { \prime } | \leq \rho _ { 0 } } \frac { | u ( x, t ) - u ( x, t ^ { \prime } ) | } { | t - t ^ { \prime } | ^ { \tau } }, \\ & 0 < \tau < 1, \\ \langle u \rangle _ { x, Q _ { T } } ^ { ( \tau ) } & = \sup _ { t \in ( 0, T ) } \sup \rho ^ { - \tau } \text{osc} _ { x } \left \{ u ( x, t ) ; \Omega _ { \rho } ^ { i } \right \}.$$

In the last line, the second supremum is taken over all connected components Ω i ρ of all Ω ρ with ρ ≤ ρ 0 , while osc x { u x, t ( ); Ω i ρ } is the oscillation of u x, t ( ) on Ω i ρ , i.e., osc x { u x, t ( ); Ω i ρ } = { vrai sup x ∈ Ω i ρ u x, t ( ) -vrai inf x ∈ Ω i ρ u x, t ( ) } , ρ 0 = x -x , Ω ρ = K ρ ∩ ( x, x ) , and K ρ is an arbitrary open interval in ( -∞ , + ∞ ) of length 2 ρ .

## 8.2 Proof of Lemma 1

The proof of Lemma 1 can be proceeded in four steps.

Step 1: We prove the existence of solutions of (21) and (22) by induction combining with the generalization of the approach developed in [40].

Lemma 2 For any n , (21) and (22) have a unique solution w n ∈ H 2+ θ, 1+ θ 2 ( Q T ) and v n ∈ H 2+ θ, 1+ θ 2 ( Q T ) respectively. Proof First we prove the existence and the uniqueness of the solution w n ∈ H 2+ θ, 1+ θ 2 ( Q T ) of (21) . We rewrite (21) as below:

$$\frac { \partial w _ { n } } { \partial t } ( x, t ) - \beta w _ { n x x } + b _ { n } ( x, t, w _ { n }, w _ { n x } ) = 0, & & ( 2 8 a ) \\ \wedge \quad \, \nu = \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \colon & & \wedge \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \, \wedge & & \wedge \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cd.$$

$$\beta w _ { n x } ( \overline { x }, t ) + \psi _ { n } ( \overline { x }, t, w _ { n } ) = 0,$$

$$\beta w _ { n x } ( x, t ) - \psi _ { n } ( x, t, w _ { n } ) = 0,$$

$$w _ { n } ( x, 0 ) = w ^ { 0 } ( x ),$$

where b n ( x, t, w, p ) = α 1 x w +( α 1 + G w n -1 ( )) t p -δ n ( x, t ) and ψ n ( x, t, w ) = ( x + x -2 x )( α 1 + G w n -1 ( )) t w .

We proceed by induction to show that (28) has a unique solution w n ∈ H 2+ θ, 1+ θ 2 ( Q T ) which satisfies:

$$\max _ { \overline { Q } _ { T } } | w _ { n } | \leq m _ { n }, \max _ { \overline { Q } _ { T } } | w _ { n z } | \leq M _ { n }, | w _ { n } | _ { Q _ { T } } ^ { ( 2 + \theta ) } \leq \overline { M } _ { n },$$

where m n depends only on a, b, k 0 , y d , P, η, β , ‖ α 1 ‖ ∞ , ‖ α 1 x ‖ ∞ , ‖ ‖ δ ∞ , max Q T | w n -1 | , ∣ ∫ ∣ x x w 0 d x ∣ ∣ and n M , n depends only on β x x m , , , n , ‖ α 1 ‖ ∞ , ‖ α 1 x ‖ ∞ , max Q T | w n -1 | , | w 0 (2) | ( x,x ) , and M n depends only on M n , m n , θ , β δ , ̂ , | w 0 (2) | ( x,x ) .

Indeed, for n = 1 , recalling that α , w 1 0 ∈ C 1 ([ x, x ]) and ∣ ∫ ∣ x x w 0 d x ∣ ∣ &gt; 0 , there exists a constant g 0 &gt; 0 depending only on a, b, k 0 , y d , P, η, β , ‖ α 1 ‖ ∞ , max [ x,x ] | w 0 | and ∣ ∫ ∣ x x w 0 d x ∣ ∣ such that | G w 0 ( ) t | ≤ g 0 for any t ≥ 0 . For any w and p , we have

$$- w b _ { 1 } ( x, t, w, p ) \leq & \left ( \| \alpha _ { 1 x } \| _ { \infty } + \frac { \| \alpha _ { 1 } \| _ { \infty } + g _ { 0 } + 1 } { 2 } \right ) w ^ { 2 } + \frac { \| \alpha _ { 1 } \| _ { \infty } + g _ { 0 } } { 2 } p ^ { 2 } + \frac { \| \delta _ { 1 } \| _ { \infty } ^ { 2 } } { 2 } \\ \leq & \left ( \| \alpha _ { 1 x } \| _ { \infty } + \frac { \| \alpha _ { 1 } \| _ { \infty } + g _ { 0 } + 1 } { 2 } \right ) w ^ { 2 } + \frac { \| \alpha _ { 1 } \| _ { \infty } + g _ { 0 } } { 2 } p ^ { 2 } + \frac { 1 + \| \delta \| _ { \infty } ^ { 2 } } { 2 }$$

$$\coloneqq & c _ { 1 0 } p ^ { 2 } + c _ { 1 1 } w ^ { 2 } + c _ { 1 2 }, \\ - w \psi _ { 1 } ( x, t, w, p ) \leq & ( \| \alpha _ { 1 } \| _ { \infty } + g _ { 0 } ) \coloneqq c _ { 1 3 } w ^ { 2 }.$$

Then all the structural conditions of [40, Theorem 7.4, Chap. V] are satisfied. Therefore, (28) has a unique solution w 1 ∈ H 2+ θ, 1+ θ 2 ( Q T ) . Moreover, by [40, Theorem 7.3, 7.2, Chap. V], we have the following estimates:

$$\max _ { \overline { Q } _ { T } } | w _ { 1 } | \leq & \lambda _ { 1 1 } e ^ { \lambda _ { 1 } T } \max \{ \sqrt { c _ { 1 2 } }, w ^ { 0 } ( \overline { x } ), w ^ { 0 } ( x ) \} \coloneqq m _ { 1 },$$

$$\max _ { \overline { Q } _ { T } } | w _ { 1 x } | \leq & M _ { 1 },$$

$$\max _ { \varepsilon } | w _ { 1 z } | \leq M _ { 1 },$$

where λ 11 , λ 1 depend only on β, c 10 , c 11 , and c 13 , M 1 depends only on β, x, x, m 1 , g 0 , ‖ α 1 ‖ ∞ , ‖ α 1 x ‖ ∞ , and | w 0 (2) | ( x,x ) . Furthermore, by [40, Theorem 5.4, Chap. V], we have

$$| w _ { 1 } | _ { Q _ { T } } ^ { ( 2 + \emptyset ) } \leq \overline { M } _ { 1 },$$

where M 1 depends only on M ,m ,θ,β,δ 1 1 ̂ , and | w 0 (2) | ( x,x ) .

Assuming that for n = k (28) has a unique solution w n ∈ H 2+ θ, 1+ θ 2 ( Q T ) satisfying (29) , we prove the claim for n = k +1 , where k &gt; 1 is an integer. Noting that by (21) , for all t ∈ [0 , T ] , it follows

$$\int _ { \underline { x } } ^ { \overline { x } } w _ { k } ( x, t ) d x = \int _ { \underline { x } } ^ { \overline { x } } w ^ { 0 } ( x ) d x + \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta _ { k } ( x, t ) d x d t.$$

By the definition of G w k ( ) t , (19b) and (29) , it follows

$$|$$

$$| G _ { w _ { k } } ( t ) | \leq g _ { k }, \forall t \in [ 0, T ],$$

where g k is a positive constant depending only on T, δ 0 , a, b, x, x, k η 0 P , y d , β, ‖ α ‖ ∞ , max Q T | w k | , ∣ ∫ ∣ x x w 0 d x ∣ ∣ . Then for any w and p , we have

$$- w b _ { k + 1 } ( x, t, w, p ) \leq & c _ { ( k + 1 ) 0 } p ^ { 2 } + c _ { ( k + 1 ) 1 } w ^ { 2 } + c _ { ( k + 1 ) 2 }, \\ - w \psi _ { k + 1 } ( x, t, w, p ) \leq & ( \| \alpha _ { 1 } \| _ { \infty } + g _ { k } ) \coloneqq c _ { ( k + 1 ) 3 } w ^ { 2 },$$

$$\begin{matrix} \end{matrix}$$

$$- ^ { \prime }$$

where c ( k +1)0 , c ( k +1)1 , c ( k +1)2 , and c ( k +1)3 are positive constants depending only on g k , ‖ α 1 ‖ ∞ , ‖ α 1 x ‖ ∞ , and ‖ ‖ δ ∞ .

By [40, Theorem 7.4, 7.3, 7.2, Theorem 5.4, Chap. V], we conclude that (28) has a unique solution w n ∈ H 2+ θ, 1+ θ 2 ( Q T ) for n = k +1 , having the estimates in (29) . Therefore, (21) has a unique solution w n ∈ H 2+ θ, 1+ θ 2 ( Q T ) for every n ≥ 1 .

We have then that for any n , there exists w n satisfying (21) . Recalling the definition of G w n ( ) t , G w n ( ) t is fixed. Therefore (22) is a linear equation. Finally, the existence of a unique solution v n ∈ H 2+ θ, 1+ θ 2 ( Q T ) to (22) is guaranteed by [40, Theorem 7.4, Chap. V]. Moreover, v n has similar estimates as (29) , i.e.,

$$\max _ { \overline { Q } _ { T } } | v _ { n } | + \max _ { \overline { Q } _ { T } } | v _ { n x } | + | v _ { n } | _ { Q _ { T } } ^ { ( 2 + \theta ) } \leq C _ { n },$$

with a certain positive constant C n .

Step 2: We prove a property of ( w ,v n n ) stated below.

Lemma 3 Let w ,v n n ∈ H 2+ θ, 1+ θ 2 ( Q T ) be the solutions to (21) and (22) respectively. For any n and any t ∈ [0 , T ] , the following equalities hold:

$$\int _ { \underline { x } } ^ { \overline { x } } w _ { n } ( x, t ) d x = \int _ { \underline { x } } ^ { \overline { x } } w ^ { 0 } ( x ) d x + \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, s ) d x d s \geq \frac { \delta _ { 0 } } { 2 },$$

$$\int _ { \underline { x } } ^ { \overline { x } } v _ { n } ( x, t ) d x = \int _ { \underline { x } } ^ { \overline { x } } v ^ { 0 } ( x ) d x - \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, s ) d x d s \geq \frac { \delta _ { 0 } } { 2 }.$$

Proof It suffices to integrating over ( x, x ) and using the boundary conditions and (20) .

Step 3: We prove the following property of ( w ,v n n ) which implies that ( w ,v n n ) is uniformly bounded in L 1 -norm. Lemma 4 Let w ,v n n ∈ H 2+ θ, 1+ θ 2 ( Q T ) be the solutions to (21) and (22) respectively. For any n and any t ∈ [0 , T ] , the following estimates hold:

$$( i ) \ \| w _ { n } ( \cdot, t ) \| _ { 1 } & \leq \| w ^ { 0 } \| _ { 1 } + \int _ { 0 } ^ { t } \int _ { x } ^ { \overline { x } } \delta _ { n } ( x, s ) s g n ( w _ { n } ) d x d s, \\ ( i i ) \ \| v _ { n } ( \cdot, t ) \| _ { 1 } & \leq \| v ^ { 0 } \| _ { 1 } - \int _ { 0 } ^ { t } \int _ { x } ^ { \overline { x } } \delta _ { n } ( x, s ) s g n ( v _ { n } ) d x d s. \\ \text{Proof } \ \text{We only prove the first inequality in the lemma. The second ineq.}$$

Proof We only prove the first inequality in the lemma. The second inequality can be assessed similarly. For ε &gt; 0 , define

$$\rho _ { \varepsilon } ( r ) = \begin{cases} | r |, \ | r | > \varepsilon, \\ - \frac { r ^ { 4 } } { 8 \varepsilon ^ { 3 } } + \frac { r ^ { 2 } } { 4 \varepsilon } + \frac { 3 \varepsilon } { 8 }, \ | r | \leq \varepsilon, \end{cases}$$

which is C 2 -continuous in r satisfying ρ ε ( r ) ≥ | r | , | ρ ′ ε ( r ) | ≤ 1 , and ρ ′′ ε ( r ) ≥ 0 . Let ˜ G n ( x, t ) = α 1 ( x ) + G w n ( ) t . By (33) , (32) , (20) , we have that

$$|$$

$$| G _ { n - 1 } ( x, t ) | \leq & \| \alpha _ { 1 } \| _ { \infty } + g _ { n - 1 }$$

$$\cdot \quad \cdot \quad \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \\ \leq \frac { 2 ( \| \alpha _ { 1 } \| _ { \infty } + g _ { n - 1 } ) } { \delta _ { 0 } } \int _ { \underline { x } } ^ { \overline { x } } | w _ { n } | d x,$$

for any ( x, t ) ∈ Q T . Multiplying ρ ′ ε ( w n ) to (21) , and integrating by parts, we have

$$\frac { d } { d t } \int _ { \underline { s } } ^ { \overline { x } } \rho _ { \varepsilon } ( w _ { n } ) d x = - \beta \| w _ { n x } \sqrt { \rho _ { \varepsilon } ^ { \prime \prime } ( w _ { n } ) } \| ^ { 2 } + \int _ { \underline { s } } ^ { \overline { x } } w _ { n } \widetilde { G } _ { n - 1 } ( x, t ) \rho _ { \varepsilon } ^ { \prime \prime } ( w _ { n } ) w _ { n x } d x + \int _ { \underline { s } } ^ { \overline { x } } \delta _ { n } ( x, t ) \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) d x.$$

We estimate the second term on the right-hand-side of (35) . By Young's inequality, it follows

$$\int _ { \underline { x } } ^ { \overline { x } } w _ { n } \widetilde { G } _ { n - 1 } ( x, t ) \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) w _ { n x } d x & \leq \frac { 1 } { 4 \beta } \| w _ { n } \widetilde { G } _ { n - 1 } \sqrt { \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) } \| ^ { 2 } + \beta \| w _ { n x } \sqrt { \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) } \| ^ { 2 } \\ & \leq \frac { 1 } { 4 \beta } \| \widetilde { G } _ { n - 1 } \| _ { \infty } ^ { 2 } \| w _ { n } \sqrt { \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) } \| ^ { 2 } + \beta \| w _ { n x } \sqrt { \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) } \| ^ { 2 }.$$

By (35) and (36) , we have

$$\frac { d } { d t } \int _ { \underline { x } } ^ { \overline { x } } \rho _ { \varepsilon } ( w _ { n } ) d x \leq \frac { 1 } { 4 \beta } \| \widetilde { G } _ { n - 1 } \| _ { \infty } ^ { 2 } \| w _ { n } \sqrt { \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) } \| ^ { 2 } + \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, t ) \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) d x.$$

Note that

$$\| w _ { n } \sqrt { \rho _ { \varepsilon } ^ { \prime \prime } ( w _ { n } ) } \| ^ { 2 } & = \int _ { \underline { x } } ^ { \underline { x } } w _ { n } ^ { 2 } \rho _ { \varepsilon } ^ { \prime \prime } ( w _ { n } ) \chi _ { \{ | w _ { n } | > \varepsilon \} } ( x ) d x + \int _ { \underline { x } } ^ { \underline { x } } w _ { n } ^ { 2 } \rho _ { \varepsilon } ^ { \prime \prime } ( w _ { n } ) \chi _ { \{ | w _ { n } | \leq \varepsilon \} } ( x ) d x \\ & = \int _ { \underline { x } } ^ { \underline { x } } w _ { n } ^ { 2 } \rho _ { \varepsilon } ^ { \prime \prime } ( w _ { n } ) \chi _ { \{ | w _ { n } | \leq \varepsilon \} } ( x ) d x \\ & = \int _ { \underline { x } } ^ { \underline { x } } w _ { n } ^ { 2 } \frac { 3 ( \varepsilon ^ { 2 } - w _ { n } ^ { 2 } ) } { 2 \varepsilon ^ { 3 } } \chi _ { \{ | w _ { n } | \leq \varepsilon \} } ( x ) d x$$

$$\leq _ { \overline { 2 } } ^ { 3 } \varepsilon ( \overline { x } - \overline { x } ).$$

By (34) , (37) , and (38) , it follows

$$\Sigma \times \cdot \sum \nolimits \nolimits \dots \cdot \sum \nolimits \nu \nolimits \dots \cdot \sum \\ \frac { d } { d t } \int _ { \underline { x } } ^ { \underline { x } } & \rho _ { e } ( w _ { n } ) d x \leq \frac { 3 \sigma ( \overline { x } - \underline { x } ) } { 4 \beta } \frac { \| \alpha _ { 1 } \| _ { \infty } + g _ { n - 1 } } { \delta _ { 0 } } \left ( \int _ { \underline { x } } ^ { \overline { x } } | w _ { n } | d x \right ) ^ { 2 } + \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, t ) \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) d x \\ & \leq \frac { 3 \varepsilon ( \overline { x } - \underline { x } ) } { 4 \beta } \frac { \| \alpha _ { 1 } \| _ { \infty } + g _ { n - 1 } } { \delta _ { 0 } } \left ( \int _ { \underline { x } } ^ { \overline { x } } \rho _ { e } ( w _ { n } ) d x \right ) ^ { 2 } + \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, t ) \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) d x \\ & \coloneqq c _ { n } \varepsilon \left ( \int _ { \underline { x } } ^ { \overline { x } } \rho _ { e } ( w _ { n } ) d x \right ) ^ { 2 } + \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, t ) \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ) d x. \\ \nu \cdot \cdot \cdot \cdot \cdot \sum \nolimits \nolimits \dots \cdot \cdot \cdot \cdot \cdot \cdot \cdot$$

We infer from Gronwall's inequality:

$$\int _ { \underline { \varepsilon } } ^ { \overline { \varepsilon } } \rho _ { \varepsilon } ( w _ { n } ) d x & \leq \int _ { \underline { \varepsilon } } ^ { \overline { \varepsilon } } \rho _ { \varepsilon } ( w ^ { 0 } ) d x \times e ^ { \epsilon _ { n } \varepsilon \int _ { 0 } ^ { t } \int _ { \underline { \varepsilon } } ^ { \overline { \varepsilon } } \rho _ { \varepsilon } ( w _ { n } ( x, s ) ) d x d s } \\ & \quad + \int _ { 0 } ^ { t } \int _ { \underline { \varepsilon } } ^ { \overline { \varepsilon } } \delta _ { n } ( x, s ) \rho _ { \varepsilon } ^ { \prime } ( w _ { n } ( x, s ) ) d x \times e ^ { \epsilon _ { n } \varepsilon \int _ { 0 } ^ { t } \int _ { \underline { \varepsilon } } ^ { \overline { \varepsilon } } \rho _ { \varepsilon } ( w _ { n } ( x, t ) ) d x d r } d s.$$

Note that w n is smooth and bounded on [ x, x ] × [0 , T ] . Letting ε → 0 , we get

$$\int _ { \underline { x } } ^ { \overline { x } } | w _ { n } ( x, t ) | d x \leq \int _ { \underline { x } } ^ { \overline { x } } | w ^ { 0 } ( x ) | d x + \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta _ { n } ( x, s ) s g n ( w _ { n } ) d x d s.$$

Step 4: We prove that every constant of the estimates of w n and v n in the proof of Lemma 4 are independent of n , i.e., we have the result stated below.

Lemma 5 Let w ,v n n ∈ H 2+ θ, 1+ θ 2 ( Q T ) be the solutions to (21) and (22) respectively. For any n , it holds:

$$\max _ { \overline { Q } _ { T } } | w _ { n } | \leq m, \max _ { \overline { Q } _ { T } } | w _ { n x } | \leq M, | w _ { n } | _ { Q _ { T } } ^ { ( 2 + \theta ) } \leq \overline { M }, \\ \dots \.$$

$$\max _ { \overline { Q } _ { T } } | v _ { n } | \leq l, \max _ { \overline { Q } _ { T } } | v _ { n z } | \leq L, | v _ { n } | _ { Q _ { T } } ^ { ( 2 + \theta ) } \leq \overline { L },$$

where m,M,M,M,l,L ˜ , and L are positive constants independent of n . Proof By the proof of Lemma 4 and ‖ δ n ‖ ∞ ≤ 1 + ‖ ‖ δ ∞ , ‖ w n ( · , t ) ‖ 1 is uniformly bounded in n . Then by the continuity of w n , we deduce that w n is uniformly bounded in n on Q T . Note that ∫ x x βw nx d x = β w ( n ( x, t ) -w n ( x, t )) . By (32) , (19b) , Lemma 4, and the a.e. convergence of δ n to δ , we have that G w n ( ) t is uniformly bounded in n on Q T . Then m ,M ,M n n n in (29) are all independent of n . Therefore (41a) holds. (41b) is a consequence of the uniform boundedness of G w n ( ) t and [40,

Theorem 5.4, 7.2-7.4, Chap. V].

## 8.3 Proof of Theorem 2

Note that as (11b) is linear in v , it suffices to prove that w 1 = w 2 in Q T . Letting ˜ w = w 1 -w 2 , it follows

$$\frac { \partial \widetilde { w } } { \partial t } = \frac { \partial } { \partial x } \left ( \beta \frac { \partial \widetilde { w } } { \partial x } - \alpha _ { 1 } ( x ) \widetilde { w } - \widetilde { w } G _ { w _ { 2 } } ( t ) - w _ { 1 } ( G _ { w _ { 1 } } ( t ) - G _ { w _ { 2 } } ( t ) ) \right ), \\ \beta \frac { \partial \widetilde { w } } { \partial x } ( x, t ) - \alpha _ { 1 } ( x ) \widetilde { w } ( x, t ) - \widetilde { w } ( x, t ) G _ { w _ { 2 } } ( t ) - w _ { 1 } ( x, t ) ( G _ { w _ { 1 } } ( t ) - G _ { w _ { 2 } } ( t ) ) = 0, \\ \beta \frac { \partial \widetilde { w } } { \partial x } ( \widetilde { x }, t ) - \alpha _ { 1 } ( \widetilde { x } ) \widetilde { w } ( \widetilde { x }, t ) - \widetilde { w } ( \widetilde { x }, t ) G _ { w _ { 2 } } ( t ) - w _ { 1 } ( \widetilde { x }, t ) ( G _ { w _ { 1 } } ( t ) - G _ { w _ { 2 } } ( t ) ) = 0, \\ \widetilde { w } ( x, 0 ) = 0.$$

Then we have

$$\frac { 1 } { 2 } \frac { d } { d t } \| \widetilde { w } \| ^ { 2 } = - \int _ { \Xi } ^ { \overline { x } } \frac { \partial \widetilde { w } } { \partial x } \left ( \beta \frac { \partial \widetilde { w } } { \partial x } - \alpha _ { 1 } ( x ) \widetilde { w } - \widetilde { w } G _ { w _ { 2 } } ( t ) \right ) \mathrm d x + \int _ { \Xi } ^ { \overline { x } } \frac { \partial \widetilde { w } } { \partial x } w _ { 1 } ( G _ { w _ { 1 } } ( t ) - G _ { w _ { 2 } } ( t ) ) \mathrm d x.$$

Note that there exist positive constants C 0 and C 1 such that

$$| \alpha _ { 1 } ( x ) + G _ { w _ { 2 } } ( t ) | \leq C _ { 0 }, \forall ( x, t ) \in \overline { Q } _ { T },$$

and

$$\left | G _ { w _ { 1 } } ( t ) - G _ { w _ { 2 } } ( t ) \right | & = \frac { \left | \int _ { \underline { x } } ^ { \overline { x } } \left ( \alpha _ { 1 } ( x ) + \frac { k _ { 0 } } { a } \left ( \alpha x + b \right ) \widetilde { w } d x \right | } { \left | \int _ { \underline { x } } ^ { \overline { x } } w _ { 2 } d x \right | } \leq C _ { 1 } \| \widetilde { w } \| _ { 1 },$$

where we used ˜ w x, t ( ) = ˜ w x, t ( ) , H¨lder's inequality, and the fact that by Proposition 2 o

$$\left | \int _ { \underline { x } } ^ { \overline { x } } w _ { 1 } \mathrm d x \right | = \left | \int _ { \underline { x } } ^ { \overline { x } } w _ { 2 } \mathrm d x \right | = \left | \int _ { \underline { x } } ^ { \overline { x } } w ^ { 0 } ( x ) \mathrm d x + \int _ { 0 } ^ { t } \int _ { \underline { x } } ^ { \overline { x } } \delta ( x, s ) \mathrm d x \mathrm d s \right | \geq \delta _ { 0 }.$$

By (42), (43), (44), H¨lder's inequality, and Young's inequality, we have o

$$\frac { 1 } { 2 } \frac { \mathrm d } { \mathrm d t } \| \widetilde { w } \| ^ { 2 } & \leq - \beta \| \widetilde { w } _ { x } \| ^ { 2 } + \frac { \varepsilon } { 2 } \| \widetilde { w } _ { x } \| ^ { 2 } + \frac { \varepsilon } { 2 } \| \widetilde { w } _ { x } \| ^ { 2 } + \frac { C _ { 0 } ^ { 2 } } { 2 \varepsilon } \| \widetilde { w } \| ^ { 2 } + \frac { \varepsilon } { 2 } \| \widetilde { w } _ { x } \| ^ { 2 } + \frac { 1 } { 2 \varepsilon } \| w _ { 1 } \| ^ { 2 } \cdot C _ { 1 } ^ { 2 } \| \widetilde { w } \| _ { 1 } ^ { 2 } \\ & \leq - \left ( \beta - \frac { 3 \varepsilon } { 2 } \right ) \| \widetilde { w } _ { x } \| ^ { 2 } + C _ { 2 } \| \widetilde { w } \| ^ { 2 },$$

where ε &gt; 0 , C 2 is a positive constant depending only on x, x, ‖ w 1 ‖ , C 0 , C 1 , and ε . Choosing ε = 2 3 β , it follows

$$\frac { \mathrm d } { \mathrm d t } \| \widetilde { w } \| ^ { 2 } \leq C _ { 3 } \| \widetilde { w } \| ^ { 2 },$$

where C 3 is a positive constant. By Gronwall's inequality, we have

$$\| \widetilde { w } ( \cdot, t ) \| ^ { 2 } \leq \| \widetilde { w } ( \cdot, 0 ) \| ^ { 2 } \cdot { \mathrm e } ^ { C _ { 3 } t } = 0,$$

based on which and by the continuity of ˜ w , it yields ˜ w ≡ 0 on Q T .

## References

- [1] D. S. Callaway, 'Tapping the energy storage potential in electric loads to deliver load following and regulation, with application to wind energy,' Energy Conv. Manag. , vol. 50, no. 5, pp. 1389-1400, 2009.
- [2] G. S. Ledva, E. Vrettos, S. Mastellone, G. Andersson, and J. L. Mathieu, 'Managing communication delays and model error in demand response for frequency regulation,' IEEE Trans. Power Syst. , vol. 33, no. 2, pp. 1299-1308, March 2018.
- [3] M. Liu, Y. Shi, and X. Liu, 'Distributed MPC of aggregated heterogeneous thermostatically controlled loads in smart grid,' IEEE Trans. Ind. Electron. , vol. 63, no. 2, pp. 1120-1129, Feb. 2016.
- [4] M. Liu and Y. Shi, 'Model predictive control of aggregated heterogeneous second-order thermostatically controlled loads for ancillary services,' IEEE Trans. Power Syst. , vol. 31, no. 3, pp. 1963-1971, May 2016.
- [5] N. Mahdavi, J. H. Braslavsky, M. M. Seron, and S. R. West, 'Model predictive control of distributed air-conditioning loads to compensate fluctuations in solar power,' IEEE Trans. Control Syst. Technol. , vol. 8, no. 6, pp. 3055-3065, Nov. 2017.
- [6] D. Angeli and P.-A. Kountouriotis, 'A stochastic approach to 'dynamic-demand' refrigerator control,' IEEE Trans. Control Syst. Technol. , vol. 20, no. 3, pp. 581-592, May 2012.
- [7] M. Brandstetter, A. Schirrer, M. Mileti´ c, S. Henein, M. Kozek, and F. Kupzog, 'Hierarchical predictive load control in smart grids,' IEEE Trans. Smart Grid , vol. 8, no. 1, pp. 190-199, Jan. 2017.

| [8]                 | L. D. Collins and R. H. Middleton, 'Distributed demand peak reduction with non-cooperative players and minimal communication,' IEEE Trans. Smart Grid , vol. 10, no. 1, pp. 153-162, Aug. 2019.                                                                                                                                  |
|---------------------|----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------|
| [9]                 | W. Gu, Z. Wang, Z. Wu, Z. Luo, Y. Tang, and J. Wang, 'An online optimal dispatch schedule for CCHP microgrids based on model predictive control,' IEEE Trans. Smart Grid , vol. 8, no. 5, pp. 2332-2342, Sept. 2017.                                                                                                             |
| [10] J.             | Ma, H. Chen, and L. Song, 'Residential load scheduling in smart grid: A cost efficiency perspective,' IEEE Trans. Smart Grid , vol. 7, pp. 771-784, March 2016.                                                                                                                                                                  |
| [11]                | K. X. Perez, M. Baldea, and T. F. Edgar, 'Integrated HVACmanagement and optimal scheduling of smart appliances for community peak load reduction,' Energy and Buildings , vol. 123, pp. 34-40, 2016.                                                                                                                             |
| [12]                | J. Qi, Y. Kim, C. Chen, X. Lu, and J. Wang, 'Demand response and smart buildings: A survey of control, communication, and cyber-physical security,' ACM Transactions on Cyber-Physical Systems , vol. 1, no. 4, p. 18, Oct. 2017.                                                                                                |
| [13]                | W. H. Sadid, S. A. Abobakr, and G. Zhu, 'Discrete-event systems-based power admission control of thermal appliances in smart buildings,' IEEE Trans. Smart Grid , vol. 8, no. 6, pp. 2665-2674, Nov. 2017.                                                                                                                       |
| [14]                | S. H. Tindemans, V. Trovato, and G. Strbac, 'Decentralized control of thermostatic loads for flexible demand response,' IEEE Trans. Control Syst. Technol. , vol. 23, no. 5, pp. 1685-1700, Sept. 2015.                                                                                                                          |
| [15]                | L. C. Totu, R. Wisniewski, and J. Leth, 'Demand response of a TCL population using switching-rate actuation,' IEEE Trans. Control Syst. Technol. , vol. 25, no. 5, pp. 1537-1551, Sept. 2017.                                                                                                                                    |
| [16] J.             | Yao, G. T. Costanzo, G. Zhu, and B. Wen, 'Power admission control with predictive thermal management in smart buildings,' IEEE Trans. Ind. Electron. , vol. 62, no. 4, pp. 2642-2650, April 2015.                                                                                                                                |
| [17] W.             | Zhang, J. Lian, C.-Y. Chang, and K. Kalsi, 'Aggregated modeling and control of air conditioning loads for demand response,' IEEE Trans. Power Syst. , vol. 28, no. 4, pp. 4655-4664, Nov. 2013.                                                                                                                                  |
| [18] A.             | C. Kizilkale, R. Salhab, and R. Malham´, e 'An integral control formulation of mean field game based large scale coordination of loads in smart grids,' Automatica , vol. 100, pp. 312-322, 2019.                                                                                                                                |
| [19]                | T. Voice, 'Stochastic thermal load management,' IEEE Trans. Autom. Control , vol. 63, no. 4, pp. 931-946, April 2018.                                                                                                                                                                                                            |
| [20] L.             | Zhao, W. Zhang, H. Hao, and K. Kalsi, 'A geometric approach to aggregate flexibility modeling of thermostatically controlled loads,' IEEE Trans. Power Syst. , vol. 32, no. 6, pp. 4721-4731, Nov. 2017.                                                                                                                         |
| [21] R.             | Malham´ e and C.-Y. Chong, 'Electric load model synthesis by diffusion approximation of a high-order hybrid-state stochastic system,' IEEE Trans. Autom. Control , vol. 30, no. 9, pp. 854-860, Sept. 1985.                                                                                                                      |
| [22]                | N. Beeker, P. Malisani, and N. Petit, 'Modeling populations of electric hot water tanks with Fokker-Planck equations,' in 2nd IFAC Workshop on Control of Systems Governed by Partial Differential Equations , Bertinoro, Italy, June 13-15 2016, pp. 66-73.                                                                     |
| [23] L.             | Zhao and W. Zhang, 'A unified stochastic hybrid system approach to aggregate modeling of responsive loads,' IEEE Trans. Autom. Control , vol. 63, no. 12, pp. 4250-4263, Dec. 2018.                                                                                                                                              |
| [24] S.             | Bashash and H. K. Fathy, 'Modeling and control of aggregate air conditioning loads for robust renewable power management,' IEEE Trans. Control Syst. Technol. , vol. 21, no. 4, pp. 1318-1327, July 2013.                                                                                                                        |
| [25] S.             | Moura, V. Ruiz, and J. Bendtsen, 'Modeling heterogeneous populations of thermostatically controlled loads using diffusion-advection PDEs,' in ASME Dynamic Systems and Control Conference , Palo Alto, California, USA, Oct. 2013, p. V002T23A001.                                                                               |
| [26] A.             | Ghaffari, S. Moura, and M. Krsti´, c 'Modeling, control, and stability analysis of heterogeneous thermostatically controlled load populations using partial differential equations,' J. Dyn. Syst. Meas. Control , vol. 137, 2015, p. 101009.                                                                                    |
| [27] S. in          | Koch, J. L. Mathieu, and D. S. Callaway, 'Modeling and control of aggregated heterogeneous thermostatically controlled loads for ancillary services,' Proc. Power Syst. Comput. Conf. , Stockholm, Sweden, Aug. 2011, pp. 1-8.                                                                                                   |
| [28] S.             | E. Z. Soudjani and A. Abate, 'Aggregation and control of populations of thermostatically controlled loads by formal abstractions,' IEEE Trans. Control Syst. Technol. , vol. 23, no. 3, pp. 975-990, May 2015.                                                                                                                   |
| [29] H. no.         | Hao, B. M. Sanandaji, K. Poolla, and T. L. Vincent, 'Aggregate flexibility of thermostatically controlled loads,' IEEE Trans. Power Syst. , vol. 20, 1, pp. 189-198, Jan. 2015.                                                                                                                                                  |
| [30] C.             | Le Floch, E. Kara, and S. Moura, 'PDE modeling and control of electric vehicle fleets for ancillary services: A discrete charging case,' IEEE Trans. Smart Grid , vol. 9, no. 2, pp. 573-581, March 2018.                                                                                                                        |
| [31] M.             | Ghanavati and A. Chakravarthy, 'Demand-side energy management by use of a design-then-approximate controller for aggregated thermostatic loads,' IEEE Trans. Control Syst. Technol. , vol. 26, no. 4, pp. 1439-1448, July 2018.                                                                                                  |
| [32] M.             | J. Balas, 'Active control of flexible dynamic systems,' Journal of Optimization Theory and Applications , vol. 25, no. 3, pp. 415-436, 1978.                                                                                                                                                                                     |
| [33] P.             | D. Christofides and P. Daoutidis, 'Feedback control of hyperbolic PDE systems,' AIChE Journal , vol. 42, no. 11, pp. 3063-3086, 1996.                                                                                                                                                                                            |
| [34] P. [35] A. no. | D. Christofides, Nonlinear and Robust Control of PDE Systems: Methods and Applications to Transport-Reaction Processes . Boston: Birkhuser, 2001. Maidi and J.-P. Corriou, 'Distributed control of nonlinear diffusion systems by input-output linearization,' Int. J. Robust Nonlinear Control , vol. 24, 3, pp. 389-405, 2014. |
| [36] [37]           | H. K. Khalil, Nonlinear Systems , 3rd ed. Englewood Cliffs, NJ: Prentice-Hall, 2002. M. Krstic, I. Kanellakopoulos, and P. Kokotovic, Nonlinear and Adaptive Control Design . New York: John Wiley &Sons, Inc., 1995.                                                                                                            |
| [38] A.             | Isidori, Nonlinear Control Systems , London: Springer-Verlage, 1995.                                                                                                                                                                                                                                                             |
| [39] J.             | L´vine, e Analysis and Control of Nonlinear Systems: A Flatness-based Approach . Berlin: Springer-Verlag, 2009.                                                                                                                                                                                                                  |
| [40] O.             | A. Ladyzenskaja, V. A. Solonnikov, and N. N. Uralceva, Linear and Quasi-linear Equations of Parabolic Type . Providence, RI: American Mathematical Society, 1968.                                                                                                                                                                |

- [41] S. A. Abobakr, H. Sadid, Waselul, and G. Zhu, 'A game-theoretic decentralized model predictive control of thermal appliances in discrete-event systems frameworks,' IEEE Trans. Ind. Electron. , vol. 65, no. 8, pp. 6446-6456, Aug. 2018.
- [42] S. Moura, J. Bendtsen, and V. Ruiz, 'Parameter identification of aggregated thermostatically controlled loads for smart grids using PDE techniques,' Int. J. Control , vol. 87, no. 7, pp. 1373-1386, 2014.
- [43] S. Ihara and F. C. Schweppe, 'Physically based modelling of cold load pickup,' IEEE Trans. Power Appart. Syst. , vol. PAS-100, no. 9, pp. 4142-4150, September 1981.
- [44] G. Laparra, M. Li, G. Zhu, and Y. Savaria, 'Desynchronized model predictive control for large populations of thermostatically controlled loads,' IEEE Trans. Smart Grid , 2019, Early Access, doi:10.1109/TSG.2019.2923668.

[45] F. B. Hildebrand, Introduction to Numerical Analysis . Courier Corporation, 1987.
Jun Zheng Gabriel Laparra Guchuan Zhu Meng Li
Aggregate Power Control of Heterogeneous TCL Populations Governed by Fokker-Planck Equations
This paper addresses the modeling and control of heterogeneous populations of
thermostatically controlled loads (TCLs) operated by model predictive control
(MPC) schemes at level of each TCL. It is shown that the dynamics of such TCLs
populations can be described by a pair of Fokker-Planck equations coupled via
the actions in the domain. The technique of input-output feedback linearization
is used in the design of aggregate power control, which leads to a nonlinear
system in closed loop. Well-posedness analysis is carried out to validate the
developed control scheme. The closed-loop stability of the system is assessed
by a rigorous analysis. A simulation study is contacted, and the obtained
results confirm the validity and the effectiveness of the proposed approach.