M Alexandre Dolgui

M Nicanor Quijano

M Alain Gauthier

M Naly Rakoto

M Carlos Ocampo-Martinez

M Hervé Gueguen

Members of Jury

Keywords:

Since the complexity and scale of systems have been growing in the last years, distributed approaches for control and decision making are becoming more prevalent. This dissertation focuses on an important problem involving distributed control and decision making, the dynamic resource allocation in a network. To address this problem, we explore a consensus-based algorithm that does not require any centralized computation, and that is capable to deal with applications modeled either by dynamical systems or by memoryless functions. The main contribution of our research is to prove, by means of graph theoretical tools and passivity analysis, that the proposed controller asymptotically reaches an optimal solution without the need of full information.

In order to illustrate the relevance of our main result, we address several engineering applications including: distributed control for energy saving in smart buildings, management of the customers of an aggregating entity in a smart grid environment, and development of an exact distributed optimization method that deals with resource allocation problems subject to lower-bound constraints.

Finally, we explore resource allocation techniques based on classic population dynamics models. In order to make them distributed, we introduce the concept of non-well-mixed population dynamics. We show that these dynamics are capable to deal with constrained information structures that are characterized by non-complete graphs. Although the proposed non-well-mixed population dynamics use partial information, they preserve similar properties of their classic counterpart, which uses full information. Specifically, we prove mass conservation and convergence to Nash equilibrium.

LIST OF FIGURES

5.2 Graph representation of: a. The well-mixed population in Figure 5.1a. b. The non-well-mixed population in Figure 5.1b. Nodes 1,2, and 3 correspond, respectively, to strategies "circle," "triangle," and "square" of the proposed example. . . . . . . . . . . . . . . . . . . .

List of Figures

Introduction

Control and decision making are required in many engineering problems, especially in those that seek to optimize certain criteria, such as increasing economic benefits, improving the quality of service, or reducing energy consumption. In several applications, traditional centralized approaches for control and decision making are no longer appropriate since the complexity and scale of systems have increased in the last years.

Therefore, distributed schemes are becoming more prevalent [START_REF] Schneeweiß | Distributed decision making[END_REF]. The growing interest in these techniques has led to the development of a large number of distributed methods (e.g., see [5,6,[START_REF] Bemporad | Networked control systems[END_REF][START_REF] Cortes | Distributed line search via dynamic convex combinations[END_REF]9, 10] and the references therein), among which we find the algorithms based on multi-agent systems, where the agents 1 make decisions based on local information and coordinate with each other to obtain a desirable global behavior; for instance, the optimization of an overall objective function. An important issue within this field is the dynamic resource allocation over networks of agents. This problem arises when we have a limited amount of a certain resource (electric power, computing capacity, execution time), and it is necessary to establish an optimal distribution policy between some entities (e.g., loads, processors, controllers) that are managed by a set of agents connected by a communication network. Resource allocation is widely used (a comprehensive survey on this topic can be found in [11]), and it has been applied in several fields, such as economics [START_REF] Conrad | Resource economics[END_REF], smart energy systems [START_REF] Pantoja | A population dynamics approach for the dispatch of distributed generators[END_REF], distributed computation [14], communication systems [15], and virus spread mitigation [16].

An extensive literature describes distributed methods for resource allocation based on multi-agent systems: more traditional approaches use decomposition techniques [START_REF] Daniel | A tutorial on decomposition methods for network utility maximization[END_REF], 1 Entities or elements that are capable to make decisions. 1

INTRODUCTION

where the basic idea is to exploit the structure of the overall problem in order to divide it into subproblems that are then assigned to each agent in the network. Other classical techniques are the so-called incremental gradient schemes [START_REF] Sra | Optimization for machine learning[END_REF], where the solution moves following the direction of one component of the objective function's gradient by iteration, and this search direction is alternate by using a token ring scheme. We can also find some methodologies based on heuristics, such as those inspired by market mechanisms [19], [START_REF] An | Characterizing contract-based multiagent resource allocation in networks[END_REF], which use concepts from economic theory to model the interactions between the participants of the resource allocation process. More recently, there has been a widespread interest in game theoretical approaches to solve this kind of problems. Current developments on this field can be found in [21], [START_REF] Pantoja | Distributed optimization using population dynamics with a local replicator equation[END_REF] for the non-cooperative case, and in [START_REF] Saad | Coalitional game theory for cooperative micro-grid distribution networks[END_REF] for the cooperative one. Among all these methods, the appropriate coordination of agents is a crucial issue because it avoids converging to suboptimal solutions.

To ensure this coordination, some methods require the inclusion of a centralized agent; for instance, in classical decomposition techniques, the Lagrange multiplier related to the "price" of the resource is centrally adjusted to get the optimum. Nonetheless, in some real situations, especially when we have a large number of agents and the amount of communications allowed between them is strongly limited, the implementation of a centralized coordinator is impractical, too expensive, or even infeasible [START_REF] Bemporad | Networked control systems[END_REF]. By contrast, other methods are fully decentralized (e.g., [24], [START_REF] Zhu | On distributed convex optimization under inequality and equality constraints[END_REF]). These decentralized methods exploit the communication capabilities of each agent to coordinate its decisions based on the information received from its neighbors. Moreover, decentralized methodologies have important advantages, among which we highlight the increased autonomy and resilience of the whole system in the face of possible failures. This research focuses on distributed resource allocation methods among networks of agents that do not require a central coordinator. Specifically, we are interested in consensus-based approaches since their low computational requirements and appealing results on robustness, convergence rate, and flexibility (see, e.g., [START_REF] Olfati-Saber | Consensus and cooperation in networked multi-agent systems[END_REF][START_REF] Cortés | Distributed algorithms for reaching consensus on general functions[END_REF][START_REF] Moreau | Stability of multiagent systems with time-dependent communication links[END_REF][START_REF] Cao | An overview of recent progress in the study of distributed multi-agent coordination[END_REF]30]).

Regarding consensus-based algorithms for solving resource allocation problems, several works have been conducted. Most of them deals with static optimization, where the objective is to minimize a global cost function that depends directly on the resource allocated to the entities that comprise the whole system. Results on optimality and convergence of these resource allocation methods are not new. Indeed, the first reported algorithm, due to Ho et al., dates back to 1980 [24]. Since then, a large number of 1.1 Resource Allocation among Networks of Dynamical Systems variations of the cited algorithm and contributions for analyzing the proposed methods have been introduced in the literature. Most of these results are largely influenced by the progress of consensus theory [START_REF] Moreau | Stability of multiagent systems with time-dependent communication links[END_REF]31,[START_REF] Olfati-Saber | Consensus problems in networks of agents with switching topology and time-delays[END_REF][START_REF] Ren | Consensus seeking in multiagent systems under dynamically changing interaction topologies[END_REF][START_REF] Solmaz S Kia | Dynamic average consensus under limited control authority and privacy requirements[END_REF]. We highlight the work in [START_REF] Xiao | Optimal scaling of a gradient method for distributed resource allocation[END_REF],

where the authors use linear matrix inequalities in order to derive sufficient conditions that guarantee the convergence of the algorithm proposed in [24] to the optimal solution of a given problem. These conditions have been recently relaxed in [START_REF] Cherukuri | Distributed generator coordination for initialization and anytime optimization in economic dispatch[END_REF] by means of convex optimization theory. Other related results can be found in [START_REF] Johansson | Distributed non-smooth resource allocation over a network[END_REF] for non-smooth objective functions, and in [START_REF] Binetti | A distributed auction-based algorithm for the nonconvex economic dispatch problem[END_REF] for non-convex objective functions. Although resource allocation for static problems is a well studied area, there still exists a research gap when dealing with dynamic problems, especially if the dynamics are modeled by non-linear differential equations.

Resource Allocation among Networks of Dynamical Systems

For the static case, a resource allocation process is equivalent to an equalization mechanism among the marginal costs associated with each entity that conforms the entire system (this equivalence follows from the application of the Karush-Kuhn-Tucker first order necessary conditions). Therefore, if we model each entity as a memoryless system whose input is the resource allocated in the entity, and whose output is equal to the associated marginal cost, a resource allocation algorithm can be viewed as an outputconsensus protocol. Following the same ideas, if we replace the memoryless model by a dynamical model, it is possible to extend the resource allocation concept to dynamic problems, i.e., problems in which the desired objective depends on the outputs of a collection of dynamical systems that are influenced by the resource allocated to the entities.

We can find a large number of applications that fits into this framework. For instance, building temperature regulation [39], control of water distribution systems [START_REF] Ramírez-Llanos | A population dynamics approach for the water distribution problem[END_REF], bandwidth allocation [START_REF] Poveda | Dynamic bandwidth allocation in wireless networks using a Shahshahani gradient based extremum seeking control[END_REF] and access control [START_REF] Tembine | Evolutionary games in wireless networks[END_REF][START_REF] Tembine | Bioinspired delayed evolutionary game dynamics with networking applications[END_REF] in communication networks, dispatch of distributed generators for frequency regulation [START_REF] Mojica-Nava | Dynamic population games for optimal dispatch on hierarchical microgrid control[END_REF], energy saving in lighting systems [START_REF] Pantoja | Distributed optimization using population dynamics with a local replicator equation[END_REF], control for urban drainage systems [START_REF] Barreiro-Gómez | Decentralized Control for Urban Drainage Systems Via Population Dynamics: Bogotá Case Study[END_REF], and so forth.

From the above description, we formulate the resource allocation problem addressed in this dissertation as follows: given a set of dynamical systems, design a distributed

INTRODUCTION

resource allocation policy that drives the systems' outputs to consensus. Thus, our problem is twofold. On the one hand, we require an output-consensus algorithm. On the other hand, this algorithm must satisfy a given resource constraint (i.e., the constraint associated with the limited resource that has to be allocated). Although the literature on output-consensus algorithms is rich (see, e.g., [START_REF] Li | Consensus of multiagent systems and synchronization of complex networks: a unified viewpoint[END_REF][START_REF] Ren | Distributed coordination of multi-agent networks: emergent problems, models, and issues[END_REF][START_REF] Xi | Stable-protocol output consensus for high-order linear swarm systems with time-varying delays[END_REF][START_REF] Chopra | Output synchronization on strongly connected graphs[END_REF][START_REF] Xi | Output consensus analysis and design for high-order linear swarm systems: partial stability method[END_REF]51,52,53], and the references therein), few of these algorithms address the required resource constraint. Among those that take into account the considered constraint, we highlight the methods based on population dynamics [START_REF] Obando | Building temperature control based on population dynamics[END_REF][START_REF] Barreiro-Gómez | Constrained distributed optimization based on population dynamics[END_REF][START_REF] Poveda | Shahshahani gradient-like extremum seeking[END_REF][START_REF] Obando | Replicator Dynamics under Perturbations and Time-Delays[END_REF], where the resource is related to an invariant population mass. Nevertheless, these methods require a specific communication architecture (either the inclusion of a centralized coordinator, or the existence of non-trivial cliques (i.e., complete sub-graphs)). Other algorithms are based on distributed receding horizon techniques [START_REF] Wakasa | Decentralized model predictive control via dual decomposition[END_REF][START_REF] Baotić | Fast Coordinated Model Predictive Control of Large-Scale Distributed Systems with Single Coupling Constraint[END_REF][START_REF] Lefort | Hierarchical control method applied to energy management of a residential house[END_REF], where the treatment of the problem dynamics is made by taking a finite sequence of static problems. Even though this formulation is very flexible (e.g., it allows us to include transient performance criteria and other constraints different from the resource limitation), receding horizon algorithms impose a high computational burden since an optimization problem is solved at each algorithm's step. Finally, we find heuristic methods. For instance, market-based control [19,[START_REF] Clearwater | Saving energy using market-based control. Market-Based Control: A Paradigm for Distributed Resource Allocation[END_REF][START_REF] Ygge | Power load management as a computational market[END_REF][START_REF] Ettore | Market-based control in emerging distribution system operation[END_REF]. These heuristic methods generally lack a theoretical support that guarantees the proper achievement of the required goals.

To tackle the resource allocation problem among networks of dynamical systems, we use a continuous-time version of the center-free algorithm first introduced in [24] and extended in [START_REF] Xiao | Optimal scaling of a gradient method for distributed resource allocation[END_REF]. This technique belongs to the class of fully-distributed algorithms studied in the seminal work of Tsitsiklis [64]. In the framework proposed by Tsitsiklis, each system of the network is associated with an agent that manages the allocated resource taking into account a given objective. Starting with a feasible allocation, each agent updates the corresponding resource by following a given function that depends on the decisions of neighboring agents. The simple idea behind the method makes it appropriate for applications where the agents have elementary computational capabilities.

Moreover, the use of a continuous-time approach leads to a natural interaction of the considered method with dynamical systems modeled by differential equations.

The main contribution of this dissertation is to obtain sufficient conditions that guarantee asymptotical convergence of the distributed resource allocation algorithm to a desired solution, even if the decisions of each agent influence the behavior of a nonlinear 1.2 Resource Allocation Problems with Lower-Bound Constraints dynamical system. Our results are inspired by the ideas developed in [START_REF] Olfati-Saber | Consensus problems in networks of agents with switching topology and time-delays[END_REF], due to the close relationship between the center-free algorithm and the continuous-time consensus protocols analyzed in [START_REF] Olfati-Saber | Consensus problems in networks of agents with switching topology and time-delays[END_REF]. We notice that convergence conditions can be directly derived from [START_REF] Li | Consensus of multiagent systems and synchronization of complex networks: a unified viewpoint[END_REF][START_REF] Xi | Stable-protocol output consensus for high-order linear swarm systems with time-varying delays[END_REF] for systems with linear dynamics. Nonetheless, the nonlinear case is more challenging. In order to address the nonlinear case, we use passivity theory [START_REF] Wyatt | Foundations of Nonlinear Network Theory, Part I: Passivity[END_REF][START_REF] Willems | Dissipative dynamical systems part I: General theory[END_REF].

Specifically, we prove that it is possible to drive a passive system (conceptually, a passive system is a multi-port system that can store but cannot generate energy) to a desired state by using a resource-allocation-based controller. Similar results are derived in [START_REF] Chopra | Output synchronization on strongly connected graphs[END_REF][START_REF] Chopra | Passivity-based control of multi-agent systems[END_REF][START_REF] Arcak | Passivity as a design tool for group coordination[END_REF][START_REF] Bai | Cooperative control design: a systematic, passivity-based approach[END_REF]. However, in these works, resource limitation is not taken into account (a more comprehensive discussion regarding the differences of our contributions and the results given in the cited works is developed in Section 3.2.5).

To illustrate the applicability of our results, we present two engineering problems related to the management of smart energy systems. The first application deals with the control of heating, cooling, and air conditioning (HVAC) systems of large buildings that are not over-designed. In this case, there exists a trade-off between the comfort of the building's occupants and the available heating/cooling power at critical load hours. We propose a distributed approach that maximizes the comfort of the building's occupants under several power constraints. The second application is related to load aggregation for demand response [START_REF] Domınguez | Distributed algorithms for control of demand response and distributed energy resources[END_REF]. This problem is motivated by the fact that active participation of demand is a central issue in the smart grid concept. In the proposed approach, participants are encouraged to stay involved in demand response programs, by reducing the impact that they perceive when it is necessary to curtail their electric load. In all the cases, we provide simulation results that illustrate the efficiency of the proposed method.

Resource Allocation Problems with Lower-Bound Constraints

Once the key contribution on the convergence of the distributed resource allocation algorithm is set, we propose a novel technique based on a passive dynamical system, that optimally allocates a certain resource among a subset of nodes in a connected graph.

Then, this technique is applied to solve static resource allocation problems with lower bound-constraints. The inclusion of lower-bound constraints is crucial in a large number of practical applications, e.g., in [START_REF] Conrad | Resource economics[END_REF]15,[START_REF] Pantoja | Distributed optimization using population dynamics with a local replicator equation[END_REF], where it is required to capture the

INTRODUCTION

non-negativity of the resource allocated to the different entities that are involved in the problem formulation. We use a passivity-based analysis that relies on the tools previously developed in order to prove that the proposed algorithm asymptotically converges to the optimal solution under some mild assumptions related to the convexity of the cost function, and the connectivity of the graph that represents the communication topology.

In order to illustrate our theoretical results, we perform some simulations and compare our method with other techniques reported in the literature.

Although there exists an extensive literature regarding distributed methods for solving static resource allocation problems, this field still attracts considerable research attention [16,[START_REF] Cherukuri | Distributed generator coordination for initialization and anytime optimization in economic dispatch[END_REF][START_REF] Obando | Building temperature control based on population dynamics[END_REF][START_REF] Pantoja | Dispatch of distributed generators under local-information constraints[END_REF][START_REF] Tan | Consensus based approach for economic dispatch problem in a smart grid[END_REF]. A large number of solution methods belong to the general class of NUM algorithms (a survey can be found in [START_REF] Daniel | A tutorial on decomposition methods for network utility maximization[END_REF]). However, most of these methods require either the inclusion of a centralized agent or the use of restrictive information structures (as it is pointed out in [START_REF] Mosk-Aoyama | Fully distributed algorithms for convex optimization problems[END_REF]). Other approaches that solve the static resource allocation problem are inspired on the algorithm presented in [START_REF] Xiao | Optimal scaling of a gradient method for distributed resource allocation[END_REF]. The key difference of the approaches based on [START_REF] Xiao | Optimal scaling of a gradient method for distributed resource allocation[END_REF] compared to our technique is that they do not allow the explicit inclusion of lower bounds on the decision variables, unless these lower bounds are added by means of barrier functions (either logarithmic or exact [START_REF] Cherukuri | Distributed generator coordination for initialization and anytime optimization in economic dispatch[END_REF]). The problem of using barrier functions is that they can adversely affect the convergence time (in the case of using exact barrier functions) and the accuracy of the solution (in the case of using classic logarithmic barrier functions), especially for large-scale problems [74]. There are other methods that consider lower bound constraints in the problem formulation.

For instance, the authors in [START_REF] Tan | Consensus based approach for economic dispatch problem in a smart grid[END_REF][START_REF] Dominguez-Garcia | Decentralized optimal dispatch of distributed energy resources[END_REF], have developed a decentralized technique based on broadcasting and consensus to optimally distribute a resource considering capacity constraints on each entity in the network. Nonetheless, compared to our algorithm, the approach in [START_REF] Tan | Consensus based approach for economic dispatch problem in a smart grid[END_REF][START_REF] Dominguez-Garcia | Decentralized optimal dispatch of distributed energy resources[END_REF] is only applicable to quadratic cost functions. On the other hand, authors in [START_REF] Pantoja | Distributed optimization using population dynamics with a local replicator equation[END_REF] propose a novel methodology based on population dynamics. The main drawback of this technique is that its performance is seriously degraded when the number of communication links decreases. We point out the fact that other distributed optimization algorithms can be applied to solve static resource allocation problems, as those presented in [START_REF] Nedic | Constrained consensus and optimization in multi-agent networks[END_REF], [77], and [START_REF] Johansson | Distributed non-smooth resource allocation over a network[END_REF]. Nevertheless, the underlying idea in these methods is different from the one used in our work, i.e., the authors in [START_REF] Nedic | Constrained consensus and optimization in multi-agent networks[END_REF], [77], and [START_REF] Johansson | Distributed non-smooth resource allocation over a network[END_REF] use consensus steps to refine an estimation of the system state, while in our approach, consensus is used to equalize a quantity that depends on both the marginal cost perceived

Distributed Population Dynamics

by each agent in the network, and the Karush-Kuhn-Tucker (KKT) multiplier related to the corresponding resource's lower bound. In this regard, it is worth noting that the method studied in this paper requires less computational effort than the methods mentioned above. Finally, there are other techniques based on game theory and mechanism design [78, 79] that decompose and solve the original problem. Nonetheless, they need that each agent broadcasts a variable to all the other agents, i.e., a communication topology given by a complete graph is required. In contrast, the method developed in our dissertation only uses a communication topology given by a connected graph, which generally requires lower infrastructure.

Newly, in order to illustrate the performance of the proposed technique, we present an application that focuses on energy efficiency in smart buildings. More precisely, we present a distributed solution to the optimal chiller loading problem in multiple chiller systems [START_REF] Lee | Optimal chiller loading by differential evolution algorithm for reducing energy consumption[END_REF], which are widely used in large air-conditioning systems. The goal is to distribute the cooling load among the chillers that comprise the plant to minimize the total amount of power used by them.

Distributed Population Dynamics

Finally, this dissertation explores resource allocation methods based on population dynamics. Population dynamics have been widely used in the design of learning and control systems for networked engineering applications, where the information dependency among elements of the network has become a relevant issue. Classic population dynamics (e.g., replicator, logit choice, Smith, and projection) require full information to evolve to the solution (Nash equilibrium). The main reason is that classic population dynamics are deduced by assuming well-mixed populations, which limits the applications where this theory can be implemented. In this research, we extend the concept of population dynamics for non-well-mixed populations in order to deal with distributed information structures that are characterized by non-complete graphs. Although the distributed population dynamics proposed in this document use partial information, they preserve similar characteristics and properties of their classic counterpart. Specifically, we prove mass conservation and convergence to Nash equilibrium. To illustrate the performance of the proposed dynamics, we show some applications in the solution of optimization problems, classic games, and the design of distributed controllers.

INTRODUCTION

Research Questions and Outline of the Contributions

The purpose of this doctoral thesis can be summarized by the following research questions:

Q1. Can consensus-based resource allocation algorithms be extended to networks of nonlinear dynamical systems? Under what conditions, classic results on distributed resource allocation hold for this kind of networks?

Q2. How to address lower-bound constraints in static resource allocation problems by using distributed exact algorithms?

Q2. Is it possible to modify classic resource allocation methods based on population dynamics in order to avoid the need of full information? What conditions must be satisfied to guarantee that modified population dynamics that use local information preserve the optimality characteristics exhibited by their classic counterpart?

The above questions are tackled in this document as follows.

Chapter 2

In this chapter, we present some basic definitions and preliminary results on algebraic graph theory and passivity theory. We also introduce an original lemma, related to some properties of Laplacian matrices of connected graphs (Lemma 2.1.1), that is a key tool in subsequent analyses.

Chapter 6

Conclusions and future directions are discussed in Chapter 6.

Other Contributions

The following publications had a significant influence on some of the thesis contributions, but they are not covered in this document.

• G. Obando, J. I. Poveda, and N. Quijano, "Replicator Dynamics Under Perturbations and Time Delays," Submitted to Mathematics of Control, Signal, and Systems, 2015.

INTRODUCTION

Preliminaries

First, we describe the notation used throughout the document and presents some preliminary results on graph theory and passivity that are used in the proofs of our main contributions.

Graph Theory

In the multi-agent framework considered in this dissertation, we use a graph to model the communication network that allows the agents to coordinate their decisions. A graph is mathematically represented by the triplet G = (V, E, A), where V = {1, . . . , n} is the set of nodes, E ⊆ V × V is the set of edges connecting the nodes, and the adjacency matrix A is an n × n nonnegative matrix that satisfies: a ij = 1 if and only if (i, j) ∈ E, and a ij = 0 if and only if (i, j) / ∈ E. Each node of the graph corresponds to an agent of the multi-agent system, and the edges represent the available communication channels (i.e., (i, j) ∈ E if and only if agents i and j can share information). We assume that there is no edges connecting a node with itself, i.e., a ii = 0, for all i ∈ V; and that the communication channels are bidirectional, i.e., a ij = a ji . The last assumption implies that G is undirected. Additionally, we denote by N i = {j ∈ V : (i, j) ∈ E}, the set of neighbors of node i, i.e., the set of nodes that are able to receive/send information from/to node i.

Let us define the n × n matrix L(G) = [l ij ], known as the graph Laplacian of G, as

PRELIMINARIES

follows:

l ij =    j∈V a ij if i = j -a ij if i = j.
(2.1)

Properties of L(G) are related to connectivity characteristics of G as shown in the following theorem. We remark that a graph G is said to be connected if there exists a path connecting any pair of nodes.

Theorem 2.1.1. (Adapted from [START_REF] Godsil | Algebraic graph theory[END_REF]) An undirected graph G of order n is connected if and only if rank(L(G)) = n -1.

From Equation (2.1), it can be verified that

L(G)1 = 0, (2.2)
where 1 = [1, . . . , 1] , 0 = [0, . . . , 0] . A consequence of this fact is that L(G) is a singular matrix. However, we can modify L(G) to obtain a nonsingular matrix as shown in the following lemma, which is the first result of our dissertation.

Lemma 2.1.1. Let L kr (G) ∈ R (n-1
)×n be the submatrix obtained by removing the kth row of the graph Laplacian L(G), and let L k (G) ∈ R (n-1)×(n-1) be the submatrix obtained by removing the kth column of L kr (G). If G is connected, then L k (G) is positive definite. Furthermore, the inverse matrix of L k (G) satisfies L k (G) -1 l kr k = -1, where l kr k is the kth column of the matrix L kr (G).

Proof. First, notice that L(G) is a symmetric matrix because G is an undirected graph. Moreover, notice that according to Equation (2.1), L(G) is diagonally dominant with non-negative diagonal entries. The same holds for L k (G) since this is a sub-matrix obtained by removing the kth row and column of L(G). Thus, to show that L k (G) is positive definite, it is sufficient to prove that L k (G) is nonsingular.

According to Theorem 2.1.1, since G is connected, L(G) has exactly n -1 linearly independent columns (resp. rows). Let us show that the kth column (resp. row) of L(G) can be obtained by a linear combination of the other columns (resp. rows), i.e., the kth column (resp. row) is not linearly independent of the rest of the columns (resp. rows).

Since L(G)1 = 0, notice that l ik = -j∈V,j =k l ij , for all i ∈ V, i.e., the kth column can be obtained by a linear combination of the rest of the columns. Furthermore, since L(G) is a symmetric matrix, the same occurs with the kth row. Therefore, the submatrix L k (G) is nonsingular since its n -1 columns (resp. rows) are linearly independent.

Graph Theory

Now, let us prove that L k (G)

-1 l kr k = -1. In order to do so, we use the fact that L k (G)

-1 L k (G) = I, where I is the identity matrix. Hence, by the definition of matrix multiplication, we have that

n-1 m=1 lk im l k mj = 1 if i = j 0 if i = j , (2.3)
where l k ij and lk ij are the elements located in the ith row and jth column of the matrices L(G) and L k (G)

-1 , respectively. Thus,

n-1 m=1 lk im l k mi = 1, for all i = 1, . . . n -1. (2.4)
Let l kr km be the mth entry of the vector l kr k . Notice that, according to the definition of L k (G) and since L(G)1 = 0, l k mi = -n-1 j=1,j =i l k mjl kr km Replacing this value in Equation (2.4), we obtain

- n-1 j=1,j =i n-1 m=1 lk im l k mj - n-1 m=1 lk im l kr km = 1, for all i = 1, . . . n -1. According to Equation (2.3), n-1 j=1,j =i n-1 m=1 lk im l k mj = 0. This implies that n-1 m=1 lk im l kr km = -1, for all i = 1, . . . , n -1. Therefore, L k (G) -1 l kr k = -1.
Theorem 2.1.1 and Lemma 2.1.1 will be used in the analysis of the methods proposed in this dissertation, which must satisfy the information constraint imposed by the graph G. This constraint can be formally stated by means of the concept of distributed map.

A distributed map over the graph G is defined as follows [START_REF] Cortés | Distributed algorithms for reaching consensus on general functions[END_REF].

Definition 2.1.1. Given two sets R, Z, a function O : R n → Z n is distributed over G if there exist functions O 1 , . . . , O n : R × B(R) → Z (B(R) denotes the collection of all subsets of R), such that O i (r 1 , . . . , r n ) = O i (r i , {r j |j ∈ N i }) ,
for all (r 1 , . . . , r n ) ∈ R n , and all i = 1, . . . , n.

Roughly speaking, we say that a map is distributed over the graph G if it can be computed by each node using only local information.

PRELIMINARIES

Passivity

Passivity theory [START_REF] Willems | Dissipative dynamical systems part I: General theory[END_REF] is a useful tool in the analysis of nonlinear systems. In this section we provide some definitions regarding this field.

Memoryless Systems

Let us consider the following memoryless system

y = ϑ(x), (2.5)
where ϑ : R n → R n is continuous, x ∈ R n is the system input, and y ∈ R n is the system output.

Definition 2.2.1. (Adapted from [START_REF] Hassan | Nonlinear systems[END_REF]) We say that the memoryless system given in Equation (2.5) is passive if x y ≥ 0, for all x ∈ R n .

Passivity concept can be extended to dynamical systems as follows.

Dynamical Systems

Let us consider a dynamical system of the form

ξ = ζ(ξ, x) y = ϑ(ξ, x), (2.6)
where ξ ∈ R m is the system state; x ∈ R n is the system input; y ∈ R n is the system output; ζ : R m × R n → R m is locally Lipschitz; and ϑ : R m × R n → R n is continuous.

Furthermore, ζ(0, 0) = 0 and ϑ(0, 0) = 0. Notice that we assume that the number of inputs and outputs of the system is the same.

Definition 2.2.2. (Adapted from [START_REF] Hassan | Nonlinear systems[END_REF]). We say that the dynamical system given in Equation (2.6) is passive if there exists a continuously differentiable storage function

V (ξ) ≥ 0 such that x y ≥ V (ξ), for all (ξ, x) ∈ R m × R n .
Additionally, we say that this system is:

• Lossless if x y = V (ξ), for all (ξ, x) ∈ R m × R n .
• Strictly passive if x y ≥ V (ξ) + ψ(ξ), for some positive definite function ψ, and

for all (ξ, x) ∈ R m × R n .
Passivity theory and Lyapunov stability are closely related. Indeed, there are some results that link passivity characteristics of systems connected in feedback configurations to stability of the closed loop response. In this regard, we notice that the definitions of storage function and Lyapunov function are similar. A key difference is that Lyapunov functions have to be positive definite whereas storage functions only need to be positive semidefinite. Nonetheless, there exists a result due to Krasovskii and LaSalle that relaxes the positive definiteness requirement on Lyapunov functions. This result is described as follows.

Consider the system ξ = (ξ) (2.7) where : R m → R m is a locally Lipschitz map.

Theorem 2.2.1. Invariance Principle (Adapted from [START_REF] Hassan | Nonlinear systems[END_REF]). Let Ω ⊂ R m be a compact set that is positively invariant with respect to (2.7). Let V : R m → R be a continuously differentiable function such that V (ξ) ≤ 0 in Ω. Let S be the set of all points in Ω where V (ξ) = 0. Let M be the largest invariant set in S. Then every solution of (2.7) starting in Ω approaches M as t → ∞.

In this document, we use passivity theory and the invariance principle of Theorem 3

Distributed Resource Allocation Among Dynamic Networks

Network Dynamics

In several engineering problems, we have to deal with large-scale systems, which are composed by dynamical entities that interact with each other. We can find a large number of examples of this kind of systems, especially in the emergent field of smart grids and smart energy systems. In many of those problems, the main goal is associated with certain control criteria that depend on the state of the overall system. This overall-state dependency is challenging since, in some real situations, the geographical distribution of the entities, the limited power computation, or some privacy issues make the implementation of traditional centralized controllers (i.e., controllers that know the information of the whole system) not viable. Under these scenarios, an information constraint emerges.

Hence, it is necessary to seek for distributed solutions as the one represented by the multi-agent approach depicted in Figure 3.1.

The bottom layer of Figure 3.1 shows a plant (e.g., a smart building) comprised by n dynamical entities or sub-systems, which can be modeled by the following set of differential equations,

Σ p i : Ṫi = g i (T, x i ) y i = h i (T, x i ) , for i = 1, . . . , n, (3.1)
where Σ p i denotes the ith sub-system; T i ∈ R m i is the sub-system state, which can represent a physical quantity (e.g., the power generated by a micro-grid, or the temper- ature of a room in a building); T = [T 1 , . . . , T n ] ; x i ∈ R is the sub-system input; and y ∈ R n is the sub-system output. Furthermore, g i : R m i × R → R p i is locally Lipschitz, and h : R m i × R → R is continuous. Notice that the overall system dynamics are not necessarily decoupled since g i depends on the global state T . The set of sub-systems (3.1) can be written in compact form as follows,

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

G Plant Multi-agent based controller x i y i Σ p i x i {y i , y j : j ∈ N i } Σ c i ith entity ith agent
Σ p : Ṫ = g(T, x) y = h(T, x) , (3.2)
where x = [x 1 , . . . , x n ] , y = [y 1 , . . . , y n ] , g = [g 1 , . . . , g n ] , and h = [h 1 , . . . , h n ] .

On the other hand, the upper layer of Figure 3.1 shows the distributed control scheme based on a multi-agent architecture. In this class of scheme, we have a set of agents (i.e., entities that are capable to make decisions [START_REF] Tsitsiklis | Problems in Decentralized Decision making and Computation[END_REF]. Physically, they can be proces-

Control Objective and Resource Allocation Dynamics

The goal of the agents is to drive the plant to a desired global state. To do this, each agent is responsible for managing one sub-system using only partial information.

Specifically, the ith agent has information on the output of its sub-system and the outputs of its neighboring sub-systems, i.e., the ith agent knows the values Y N i := y i , {y j : j ∈ N i } , where N i = {j ∈ V : (i, j) ∈ E}. These values can be obtained by the ith agent by means of measurements, or by using the information received from its neighbors via the communication network G. Available information is employed by the agents to take decisions, which are characterized by the control outputs x 1 , . . . , x n . Thus, the multi-agent network is mathematically described by the following set of differential equations

Σ c i : ẋi = u i (Y N i ) , for i = 1, . . . , n, (3.3)
where u i : R |N i |+1 is locally Lipschitz. In compact form, (3.3) is denoted by the following system,

Σ c : ẋ = u(y), (3.4)
where u = [u 1 , . . . , u n ] is a distributed map over the graph G (cf., Definition 2.1.1), i.e., u i can be computed by the ith agent employing only local information.

Control Objective and Resource Allocation Dynamics

Control Objective

Several applications require that the outputs of the sub-systems that comprise the whole system reach the same value. These applications include, but are not limited to,... Problems requiring that the sub-systems outputs reach the same value are called output consensus problems. Formally, output consensus is defined as follows (adapted from [START_REF] Chopra | Passivity-based control of multi-agent systems[END_REF]), Definition 3.2.1. Consider the set of sub-systems given in Equation (3.2). We say that output consensus is reached if lim t→∞ |y i (t)y j (t)| = 0, for all i, j = 1, . . . , n, where y i (t) is the output of the sub-system Σ p i at time t.

In this document, we study the case when output consensus has to be achieved under a multi-variable constraint on the control inputs x 1 , . . . , x n . This constraint is associated with the preservation of a certain resource, and it is formulated as follows, n i=1

x i = X, (3.5)
where X ∈ R is the available resource. Notice that, if x i is the energy provided by the ith controller, then (3.5) is related to a general limitation of the total energy that can be applied to the plant. The inclusion of resource constraints is frequent in many control and decision-making problems. For instance, building temperature regulation [39], control of water distribution systems [START_REF] Ramírez-Llanos | A population dynamics approach for the water distribution problem[END_REF], bandwidth allocation [START_REF] Poveda | Dynamic bandwidth allocation in wireless networks using a Shahshahani gradient based extremum seeking control[END_REF] and access control [START_REF] Tembine | Evolutionary games in wireless networks[END_REF][START_REF] Tembine | Bioinspired delayed evolutionary game dynamics with networking applications[END_REF] in communication networks, dispatch of distributed generators for frequency regulation [START_REF] Mojica-Nava | Dynamic population games for optimal dispatch on hierarchical microgrid control[END_REF], energy saving in lighting systems [START_REF] Pantoja | Distributed optimization using population dynamics with a local replicator equation[END_REF], control for urban drainage systems [START_REF] Barreiro-Gómez | Decentralized Control for Urban Drainage Systems Via Population Dynamics: Bogotá Case Study[END_REF], and so forth.

Summarizing, the control goal of the multi-agent network (3.4) is twofold:

O1. Satisfying the constraint (3.5).

O2. Driving (3.2) to output consensus.

Resource Allocation Dynamics

In order to reach the objectives described above, it is necessary to design a set of local control laws, u 1 , . . . , u n , that can be applied by the multi-agent network (3.4). In this regard, we propose to use the following resource allocation algorithm

u i (Y N i ) = j∈N i a ij (y j -y i ), for all i = 1, . . . , n. (3.6)
We recall that N i = {j ∈ V : (i, j) ∈ E} is the set of neighbors of node i, and a ij > 0 is the weight of the link that connects nodes i and j.

Before starting the analysis of this algorithm, it is important to notice the fact that the right hand side of Equation (3.6) is a distributed map over the graph G since the ith agent only requires local information to compute the quantity j∈N i a ij (y jy i ).

We also notice the fact that the resource allocation algorithm corresponds to a classic consensus protocol, whose effectiveness has been proven in the literature. For instance, there exists a large number of results regarding the use of (3.6) to solve consensus problems when y i is the output of a memoryless system that directly depends on x i (cf.,

Control Objective and Resource Allocation Dynamics

Equation(2.5)) [START_REF] Cortés | Distributed algorithms for reaching consensus on general functions[END_REF][START_REF] Moreau | Stability of multiagent systems with time-dependent communication links[END_REF][START_REF] Olfati-Saber | Consensus problems in networks of agents with switching topology and time-delays[END_REF][START_REF] Xiao | Optimal scaling of a gradient method for distributed resource allocation[END_REF][START_REF] Cherukuri | Distributed generator coordination for initialization and anytime optimization in economic dispatch[END_REF]. However, our research deals with the dynamic case (3.2).

We are interested in proving that our control goals O1 and O2 can be achieved by using (3.6). First, let us focus on O1, i.e., satisfaction of the resource constraint (3.5).

Lemma 3.2.1. Assume that the multi-agent network (3.4) uses the resource allocation algorithm (3.6). If the initial conditions of (3.4) satisfy the resource constraint (3.5), i.e., n i=1 x i (0) = X, then x(t) satisfies (3.5) for all t ≥ 0.

Proof. It is sufficient to prove that the quantity ρ := n i=1 x i is positively invariant. To do that, let us calculate the time derivative of this quantity, which is given by ρ = n i=1 ẋi . Replacing ẋi by the dynamics of the multi-agent network (3.4) when the resource allocation algorithm (3.6) is used, we have that ρ = n i=1 j∈N i a ij (y jy i ). In matrix form, this expression is given by ρ = -1 L(G)y, where 1 is the column vector of ones. Notice that, since the graph Laplacian L(G) is symmetric, 1 L(G) = L(G)1 = 0 (cf., Equation (2.2)). Therefore, ρ = 0. This implies that the quantity ρ = n i=1 x i is positively invariant.

Therefore, the resource allocation algorithm attains O1 under a proper initialization.

Notice that this result holds regardless of the topology of the communication network G and the class of sub-systems dynamics (3.2). In the next section, we provide conditions on the graph G and on the dynamics of the sub-systems to guarantee that O2 is also attained by applying (3.6).

Convergence to Output Consensus

The set of sub-systems (3.2) controlled via the multi-agent network (3.4) can be viewed as the feedback interconnection shown in Figure 3.2. If the agents apply the resource allocation algorithm described in Equation (3.6), any equilibrium point (T * , x * ) of this feedback interconnection must satisfy the following property. Proposition 3.2.1. Let (T * , x * ) be an equilibrium point of the feedback connection shown in Figure 3.2, and let y * = h(T * , x * ) be the steady state output of Σ p . If u(y) is given by the resource allocation algorithm (3.6) and the communication graph G is connected, then y * i = y * j , for all i, j = 1, . . . , n, where y * i is the ith element of the vector y * .

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

Σ p :    Ṫ = g(T, x) y = h(T, x) Σ c : ẋ = u(y)
x y Proof. If u(y) is given by (3.6), then the dynamics of Σ c can be written in matrix form as ẋ = -L(G)y. Besides, since (T * , x * ) is an equilibrium point of the closed loop system, we have that -L(G)y * = 0. This implies that y * belongs to the null space of the graph Laplacian L(G). On the other hand, combining the result stated in Theorem 2.1.1, the property given in Equation (2.2), and the fact that G is connected by assumption, we can conclude that the null space of L(G) is equal to the linear span of 1, which is given by span(1) = {ξ ∈ R n : ξ = κ1}, where κ ∈ R is a constant. Therefore, y * = κ1. Hence, y * i = y * j , for all i, j = 1, . . . , n.

According to Definition 3.2.1, the statement in Proposition 3.2.1 implies that if the equilibrium point (T * , x * ) is asymptotically stable, then output consensus is reached. In order to analyze the stability of (T * , x * ), let us express the dynamics of Σ p and Σ c in error coordinates.

The set of sub-systems Σ p given in Equation (3.2) can be written in error coordinates as follows, ėT = g e (e T , e x ) e y = h e (e T , e x ),

where e T = T -T * , e x = xx * , and e y = yy * . Moreover, g e (e T , e x ) = g(T, x) and h e (e T , e x ) = h(T, x)h(T * , x * ), for all (T, x) ∈ R m × R n . Notice that h e (0, 0) = 0, and g e (0, 0) = 0 since (T * , x * ) is an equilibrium point of (3.2). Let us make the following assumption on (3.7).

Assumption 3.2.1. Consider the dynamical system (3.7). If g e (0, e x ) = 0, then e x = 0.

This assumption prevents the original system (3.2) to have other rest points of the form (T * , x * ), where x * = x * . It is worth noting that systems satisfying Assumption 3.2.1

Control Objective and Resource Allocation Dynamics

are common in several applications, including mobile-sensors networks, multi-machine power systems, and social swarms (see e.g., [START_REF] Arcak | Passivity as a design tool for group coordination[END_REF]).

On the other hand, the dynamics of the multi-agent network Σ c (3.4) that applies the resource allocation algorithm (3.6) is expressed in error coordinates as follows (using the compact matrix form), ėx = -L(G)e y .

This is true given the fact that L(G)(y -y * ) = L(G)y. The next original result establishes an important property of this system.

Proposition 3.2.2. Assume that x * satisfies the resource constraint (3.5), i.e., n i=1 x * i = X. If x(0) satisfies (3.5) and G is connected, then the multi-agent systems in error coordinates given in Equation (3.8) is passive and lossless from the input e y to the output -e x .

Proof. Since G is connected, we can use the following storage function,

V x (e x ) = 1 2 e k x L k (G) -1 e k x , (3.9)
where e k x is the vector obtained by removing the kth element from vector e x , and L k (G), is defined in Lemma 2.1.1. Although the inverse matrix of L k (G) is positive definite according to Lemma 2.1.1, the storage function V x (e x ) is positive semidefinite since it does not depend on the kth element of e x . The derivative of V x (e x ) along the trajectories of (3. Let us express ėk x in terms of the matrix L k (G) as follows,

ėk x = -L k (G)e k y -l kr k e y k , (3.11)
where e y k is the kth element of the vector e y , e k y is the vector obtained by removing the kth element from vector e y , and l kr k is defined in Lemma 2.1.1. Replacing (3.11) in (3.10), we have that

Vx (e x ) = -e k x e k y -e k x (L(G)) -1 l kr k e y k .
Taking into account that L k (G) -1 l kr k = -1 (cf., Lemma 2.1.1), we obtain

Vx (e x ) = -e k x e k y + 1 e k x e y k .

(3.12)

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

Notice that 1 e k x = i∈ν,i =k e x i , where e x i is the ith element of the vector e x . On the other hand, according to Lemma 3.2.1, Let us take V (e x , e T ) = V x (e x )+V T (e T ), where V x (e x ) is the storage function of (3.8) that is given in Equation (3.9), and V T (e T ) is the storage function of (3.7). Notice that, since conditions of Proposition 3.2.2 hold and (3.7) is strictly passive, we have that V (e x , e T ) = Vx (e x ) + VT (e T )

e x k (t) = - i∈ν,i =k e x i (
≤ -e x e y + e x e yψ(e T ) = -ψ(e T ), where ψ(e T ) is a positive definite function. Furthermore, let us consider the compact set Ω c = {(e x , e T ) : e x k = -i∈ν,i =k e x i , V x (e x ) + V T (e T ) ≤ c}, where c is any positive real number. Notice that Ω c is bounded since V T (e T ) is positive definite (this is true since (3.7) is strictly passive) and radially unbounded, and V x (e x ) is positive definite with respect to {e x i : i ∈ V, i = k} and radially unbounded. The set Ω c is also positively invariant with respect to (3.7)-(3.8) since Vx (e x ) + VT (e T ) ≤ 0 and i∈V e x i (t) = 0, for all t ≥ 0 (cf., Lemma 3.2.1). Additionally, we define the set S as the set where V (e x , e T ) = 0, i.e., S = {(e x , e T ) : e T = 0}. Let M be an invariant set of S. Notice that if (e T , e x ) ∈ M, then ėT = 0. This implies that g e (0, e x ) = 0. Hence, e x = 0 (given the fact that Assumption 3.2.1 holds). Thus, the only invariant set of S is M = {(0, 0)}. Since Theorem 2.2.1 states that every solution of (3.7)-(3.8) starting in Ω c approaches M as t → ∞, we can conclude that the origin of (3.7)-(3.8) is asymptotically stable.

The result in Theorem 3.2.1 establishes sufficient conditions on the communication graph (A1), initialization of the local control laws (A2), and characteristics of the controlled plant (A3) that guarantee that the control objective O2 is reached by using the resource allocation algorithm (3.6). We notice the fact that, if the controlled plant does not satisfy the requirement given in A3, then it is possible to implement a first-level controller that makes the plant passive.

Summarizing, we have designed a set of local control laws (3.6) that allows passive systems to attain the control goals O1 and O2.

Extension to Memoryless Systems

We can extend the result in Theorem 3.2.1 to deal with passive memoryless systems. To this end, we consider that, in the resource allocation algorithm (3.6), y = [y 1 , . . . , y n ] is the output of the following memoryless system y = ϑ(x).

(3.14)

where ϑ : R n → R n is a continuous function. In order to simplify the analysis, let us write the above system in error coordinates as follows, e y = ϑ e (e x ). (3.15) where e x = xx * , ϑ e (e x ) = ϑ(x)ϑ(x * ), and x * is an equilibrium point of the closed loop system given by Equations (3.4) and (3.14). Similarly to the dynamic case, we make the following assumption on (3.15).

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

Assumption 3.2.2. Consider the memoryless system (3.15). If e x ϑ e (e x ) = 0, then e x = 0.

This requirement is satisfied if, for instance, the vector ϑ is comprised of monotonic functions.

The following theorem states that the results obtained for the dynamic case also hold for the memoryless case.

Theorem 3.2.2. Let x * be an equilibrium point of the system given by (3.4) and (3.14).

Assume that the following conditions hold:

A1. The communication graph G of the multi-agent system (3.4) is connected.

A2. x * and x(0) satisfy the resource constraint (3.5).

A3. The memoryless system (3.14) expressed in error coordinates with respect to x * satisfies Assumption 3.2.2. Moreover it is passive from the input e x to the output e y .

Then O1 is attained, and the output y of the memoryless system (3.14) reaches consensus.

Proof. We proceed as follows,

• O1: The proof of the resource constraint satisfaction i.e., n i=1 x i (t) = X, for all t ≥ 0 follows from Lemma 3.2.1.

• Output consensus: To prove that the system output y reaches consensus, it is sufficient to show that x * is asymptotically stable under the system given by (3.4) and (3.14) (this equivalence follows from the fact that any equilibrium point of (3.4),(3.14) must satisfy output consensus (cf., Proposition 3.2.1)), or equivalently that the origin of (3.8),(3.15) is asymptotically stable.

Let us take V x (e x ) equal to the storage function given in Equation (3.9). Since (3.8) is lossless from e y to -e x (cf., Proposition 3.2.2), we have that Vx (e x ) = -e x e y . Replacing e y by the output of the memoryless system (3.15), we obtain Vx (e x ) = -e x ϑ e (e x ). Since (3.15) is passive and satisfies Assumption 3.2.2, we conclude that Vx (e x ) is negative definite. Finally, we define the positively invariant compact set Ω c = {e x : e x k = -i∈ν,i =k e x i , V x (e x ) ≤ c}, where c is any positive real number, and employ the invariance principle given in Theorem 2.2.1 to prove that any solution of (3.8), (3.15) starting in Ω c asymptotically converges to the origin.

Control Objective and Resource Allocation Dynamics

The statement in Theorem 3.2.2 is not surprising since the steady state behavior of passive dynamical systems is similar to the behavior of passive memoryless functions.

In fact, we can take advantage of this similarity to derive a general case that merges the results obtained in Theorems 3.2.1 and 3.2.2. This general case is formulated as follows.

Consider the function

f : R m × R n → R n defined by f (T, x) = h(T, x) + ϑ(x), (3.16)
where h(T, x) is the output of the dynamical system (3.2), and ϑ(x) is the output of the memoryless system (3.14). Assume that our problem is to attain the objective O1 and reach consensus on f (T, x). To address this problem, we propose to use the multi-agent based controller (3.4) with the following resource allocation algorithm.

u i = j∈N i f j -f i , for all i = 1, . . . , n. (3.17)
To analyze the performance of the proposed controller, we use a similar procedure to the one developed in Section 3.2.3. Thus, we begin by expressing the closed loop system

* -x * ) (ϑ(x * ) -ϑ(x * )) = 0.
If (3.18) satisfies some conditions, then we can guarantee an appropriate performance of the proposed control scheme. This set of conditions is specified in the following result.

Theorem 3.2.3. Let (T * , x * ) be an equilibrium point of the feedback interconnection of systems (3.2) and (3.4), where u is given by the resource allocation algorithm (3.17).

Assume that the following conditions hold:

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

A1. The communication graph G of the multi-agent system (3.4) is connected.

A2. x * and x(0) satisfy the resource constraint (3.5).

A3.

The system (3.2) expressed in error coordinates with respect to (T * , x * ) is strictly passive from the input e x to the output h e (e T , e x ) with radially unbounded storage function.

A4. The memoryless system (3.14) expressed in error coordinates with respect to x * is passive from the input e x to the output ϑ e (e x ).

A5. Assumption 3.2.3 holds.

Then O1 is attained, and the variable f (T, x) given in Equation (3.16) reaches consensus.

Proof. We proceed as follows:

• O1: The proof of the resource constraint satisfaction i.e., n i=1 x i (t) = X, for all t ≥ 0 is straightforward by using the same steps described in the proof of Lemma 3.2.1.

• Consensus on f (T, x): In order to prove that the variable f (T, x) reaches consensus, it is sufficient to show that (T * , x * ) is asymptotically stable under the system (3.2),(3.4) , with u defined in (3.17) (this equivalence follows from the fact that any equilibrium point of this system must satisfy consensus on f (T, x)), or equivalently that the origin of (3.18) is asymptotically stable. To do this, we take the function V (e x , e T ) = V x (e x ) + V T (e T ), where V x (e x ) is the storage function of (3.8) that is given in Equation (3.9), and V T (e T ) is the storage function of (3.7). Then, we use the same steps as in the proof of Theorem 3.2.1 to obtain that the derivative of V (e x , e T ) along the trajectories of (3.18) is given by

V (e x , e T ) = -ψ(e T ) -e x ϑ e (e x ),
where ψ(e T ) > 0, and e x ϑ e (e x ) ≥ 0 since A4 holds.

Furthermore, let us consider the compact set

Ω c = {(e x , e T ) : e x k = - i∈ν,i =k e x i , V x (e x ) + V T (e T ) ≤ c},
where c is any positive real number. Notice that Ω c is bounded since V T (e T ) is positive definite (this is true since (3.7) is strictly passive) and radially unbounded,

Control Objective and Resource Allocation Dynamics

and V x (e x ) is positive definite with respect to {e x i : i ∈ V, i = k} and radially unbounded. The set Ω c is also positively invariant with respect to (3.18) since Vx (e x ) + VT (e T ) ≤ 0 and i∈V e x i (t) = 0, for all t ≥ 0. Additionally, we define the set S as the set where V (e x , e T ) = 0, i.e., S = {(e x , e T ) : e T = 0, e x ϑ e (e x ) = 0}. Let M be an invariant set of S. Notice that if (e T , e x ) ∈ M, then ėT = 0. This implies that g e (0, e x ) = 0. Therefore, if (e T , e x ) ∈ M, then g e (0, e x ) = 0 and e x ϑ e (e x ) = 0. Hence, e x = 0 (given the fact that A5 holds). Thus, the only invariant set of S is M = {(0, 0)}. Since Theorem 2.2.1 states that every solution of (3.18) starting in Ω c approaches M as t → ∞, we can conclude that the origin of (3.18) is asymptotically stable.

In conclusion, we have shown that the resource allocation algorithm (3.17) allows the multi-agent based controller (3.4) to satisfy the resource constraint (3.5) and drive a specified quantity to consensus, as long as this quantity depends on the output of a passive system and/or the output of a passive memoryless function. The inclusion of memoryless terms is important; for instance, to capture actuator constraints as we will

show in Section 3.3.4.

Discussion

In the literature, significant research attention has been paid to the control of dynamical systems through consensus-based protocols (e.g., see [START_REF] Ren | Distributed coordination of multi-agent networks: emergent problems, models, and issues[END_REF][START_REF] Chopra | Passivity-based control of multi-agent systems[END_REF][START_REF] Arcak | Passivity as a design tool for group coordination[END_REF][START_REF] Bai | Cooperative control design: a systematic, passivity-based approach[END_REF]. Most of this research effort is directed towards achieving O2 whereas O1 is somewhat disregarded.

However, in a large number of applications O1 plays a fundamental role since it can be directly linked to the satisfaction of technical constraints on the actuators. For instance, a limitation on the available energy used to control the plant. Some works address O1 by assuming that the dynamics of the controlled plant are linear. In this case, conditions that guarantee a proper performance of the closed-loop system can be derived by employing linear matrix inequalities [START_REF] Xi | Stable-protocol output consensus for high-order linear swarm systems with time-varying delays[END_REF][START_REF] Xi | Output consensus analysis and design for high-order linear swarm systems: partial stability method[END_REF].

For the nonlinear case, passivity based theory is used in [START_REF] Chopra | Passivity-based control of multi-agent systems[END_REF][START_REF] Arcak | Passivity as a design tool for group coordination[END_REF][START_REF] Bai | Cooperative control design: a systematic, passivity-based approach[END_REF]. Nonetheless, these papers do not deal with the same problem studied in this dissertation. Once again, the key difference is that they do not take into account the control goal O1, i.e., the resource constraint (3.5). We highlight the work in [START_REF] Chopra | Passivity-based control of multi-agent systems[END_REF], where the authors propose a

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

control scheme similar to the one discussed in our research. They analyze the closed-loop response of systems of the following form Plant:

   Ṫ = g(T, x) y = h(T )
, Controller:

   x = j∈N i (y j -y i ) . (3.19)
We notice that the main difference with the distributed resource allocation algorithm

(3.3),(3.6) is that in (3.19
), the control law does not have internal dynamics, but it is the output of a memoryless system. The authors in [START_REF] Chopra | Passivity-based control of multi-agent systems[END_REF] require similar conditions to those provided in Theorem 3.2.1 to guarantee that output consensus is reached.

Comparing the proposed control scheme (3.2),(3.4) with the one shown in (3.19), it can be noticed that the latter can be adapted to take the form (3.2),(3.4) . For instance, let us consider the following transformation Plant:

   Ṫ = g(T, x) ẋ = x y = h(T ) , (3.20)
where we have augmented the state of the original plant by adding the vector x. Hence, one might thing that our result in Theorem 3.2.1 is covered by the result in [START_REF] Chopra | Passivity-based control of multi-agent systems[END_REF]. This is not true since the transformation shown in (3.20) does not preserve passivity, i.e., we cannot guarantee that if the plant in (3.19) is passive, then the augmented plant in (3.20) is also passive. To clarify this statement, we propose the following example.

Example 3.2.1. Let us consider a linear system and its corresponding transformation given by Plant:

     Ṫ = -T + x y = T , Augmented Plant :      Ṫ = -T + x ẋ = x y = T , (3.21)
where T, x, x ∈ R. The original plant is strictly passive from the input x to the output y since the time derivative of the storage function

V (T ) = 1 2 T 2 is equal to V (T ) = -T 2 + xy.
To verify passivity of the augmented plant, let us employ the following lemma and the fact that a minimal time invariant realization (as the one that describes the augmented plant in (3.21)) is passive if and only if it is positive real (see [START_REF] Wassim | Nonlinear dynamical systems and control: a Lyapunov-based approach[END_REF])

3.3 Distributed Building Temperature Control Under Power Constraints Lemma 3.2.2. Positive Real (adapted from [START_REF] Hassan | Nonlinear systems[END_REF]). Let G(s)C(sI -A) -1 B be a transfer function matrix, where (A, B.C) is a minimal realization. Then G(s) is positive real if and only if there exist matrices P = P > 0, and L such that

P A + A P = -L L (3.22a) P B = C (3.22b)
If we apply this lemma to the augmented plant in (3.21), we have that any matrix P must have the following form (using (3.22b))

P = p 1 1 0 ,
where p ∈ R. However, det(P ) = -1, i.e., P is not positive definite, which contradicts the requirements in Lemma 3.22. Thus, the augmented plant is not positive real. In this way, we have shown that, although the original plant is passive, the augmented plant does not preserve this characteristic. In conclusion, this example shows that the result given in Theorem 3.2.1, cannot be obtained directly from the contribution developed in [START_REF] Chopra | Passivity-based control of multi-agent systems[END_REF].

The first part of this chapter has presented the general framework studied in our dissertation, and the design of local control laws (3.6) that allow passive systems (dynamical or memoryless) to attain the objectives O1 and O2. To illustrate the relevance of our results, in the next sections we present two application examples. In the first example, we address the distributed temperature control in large-scale buildings whose heating systems are not over-designed (i.e., energy efficient buildings). The second example deals with a widely studied problem involving distributed decision making, the optimal resource allocation in a network of agents. Specifically, we present an application related to the distributed management of the customers of an aggregating entity in a smart grid environment.

Distributed Building Temperature Control Under Power Constraints

Motivation

Due to a widespread interest on environmental and economic issues, energy efficiency is becoming a research trending topic in the last time; an example of this fact is the

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

increasing number of scientific publications on that subject (e.g., see [START_REF] Fereidoon | Energy Efficiency: Towards the End of Demand Growth[END_REF] and the references therein). In this regard, several studies have shown that most of the total energy consumption is done by commercial and residential buildings, e.g., according to [START_REF] Alahmad | Information Technology and the Smart Grid-A Pathway to Conserve Energy in Buildings[END_REF],

buildings in Europe account for around 40% of the total energy, and over 70% for the electricity usage. Additionally, buildings are the main responsible for greenhouse gas emissions (almost 40% of the total [START_REF] Nof | Springer Handbook of Automation[END_REF]), causing harmful effects on the environment.

Therefore, it is not surprising to find in the literature new strategies, guidelines and standards that tend to maximize the buildings energy efficiency. Most of them focus on the optimization of heating, ventilating, and air conditioning (HVAC) systems, given the fact that these systems consume the largest amount of energy in a building. For instance, in [START_REF] Ma | Predictive control for energy efficient buildings with thermal storage: Modeling, stimulation, and experiments[END_REF] a model predictive control (MPC) that uses the rooms' thermal storage capacity is presented, where the main idea is to pre-heat (cool) the rooms depending on the occupation pattern to maintain the temperature in a desirable range, while the cost of the energy is minimized. Similarly, the authors in [START_REF] Lefort | Hierarchical control method applied to energy management of a residential house[END_REF] use MPC to manage the energy in residential buildings. Several constraints and models of the buildings' equipment and energy market are considered. In addition, a hierarchical structure is proposed in order to reduce the underlying computational complexity. Other applications can be found in [START_REF] Goyal | Occupancy-based zone-climate control for energy-efficient buildings: Complexity vs. performance[END_REF] and [START_REF] Zhao | An energy management system for building structures using a multi-agent decisionmaking control methodology[END_REF].

The applications mentioned above consider that the available power in the corresponding HVAC system is always enough to reach the control goals. This is a reasonable assumption, since HVAC systems are, in general, designed to handle worst case loads.

However, this over-design of HVAC equipment is one of the main reasons for energy inefficiency [START_REF] Nof | Springer Handbook of Automation[END_REF]. Recent standards (e.g., [START_REF] Long | Evaluation of ANSI/ASHRAE/USGBC/IES Standard[END_REF]) are tending to solve the problem by means of the implementation of new guidelines, such as the net-zero energy concept (where buildings do not demand energy from the electrical grid, but they are capable to supply their load by means of distributed generators). Additionally, demand response policies are becoming more common in building management systems [START_REF] Alahmad | Information Technology and the Smart Grid-A Pathway to Conserve Energy in Buildings[END_REF]. One of those strategies, which leads to the reduction of energy costs, is the limitation of the power consumption during hours of peak electricity demand. Hence, it is possible that at certain hours of the day, some building's services or comfort levels are affected.

For these reasons, it is important to consider scenarios in which the available power of the HVAC system is constrained. In this regard, we can cite the works in [START_REF] Obando | Building Temperature Control based on Population Dynamics[END_REF] and [START_REF] Clearwater | Saving energy using market-based control. Market-Based Control: A Paradigm for Distributed Resource Allocation[END_REF],

where the authors address the problem using a centralized dynamic resource allocation 3.3 Distributed Building Temperature Control Under Power Constraints approach based on population dynamics and market mechanisms, respectively. Centralized strategies have shown an appropriate performance to control buildings' HVAC systems under power constraints. However, since a building is a large scale system comprised by a large number of rooms (e.g., typical commercial buildings have more than one hundred rooms), the implementation of a centralized controller is not efficient, or even it could be unfeasible. One possibility is to use distributed methods with multiple controllers that make decisions in a coordinated way, where the coordination is made by means of a communication network, instead of a centralized agent. Thus, we propose to use a distributed methodology, which is based on the multi-agent scheme described in (3.1), to regulate the temperature inside the rooms of a building that has a central heating system (which is a common configuration in several residential buildings). We prove that the proposed control method optimally allocates the available heating power among the rooms, even when that power is strongly limited and the setpoints of the rooms are unreachable. In order to do that, we use the passivity-based analysis given in Theorem 3.2.1. For analysis purposes, before describing the controller design, let us introduce the building thermal model.

Thermal Model of a Building

In general terms, the mathematical description of the thermal behavior of building systems is complex since the temperature in the zones of a building is not homogeneous.

However, the problem is usually simplified by dividing the system into sections in which a uniform temperature is assumed, and heat transfer occurs only from section to section.

In [START_REF] Maasoumy | Model-based hierarchical optimal control design for HVAC systems[END_REF][START_REF] Obando | Population dynamics applied to building energy efficiency[END_REF], the authors outline a general model that considers the two most important components that constitute a building: rooms and walls. We use the ideas in [START_REF] Maasoumy | Model-based hierarchical optimal control design for HVAC systems[END_REF][START_REF] Obando | Population dynamics applied to building energy efficiency[END_REF] in order to model the thermal performance of a building consisting of n rooms, where each of them is enclosed by a certain number of walls that can be arranged according to different topologies (i.e., the spatial location of rooms is, in general, arbitrary). Figure 3.3 shows an example of the configuration of a building in which the rooms are arranged in a row. This type of configuration is common in several types of structures, as in livestock buildings.

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

Thermal Model for a Wall

A wall can be divided into layers with uniform temperature (as shown in Figure 3.4).

The temperature T w j,k of the layer k within the wall j is given by (

T w j,k-1 T w j,k+1 T w j,k T r i T w j,k
ρ w j,k c w j,k V w j,k ) Ṫ w j,k = K j k,k+1 (T w j,k+1 -T w j,k )+ K j k,k-1 (T w j,k-1 -T w j,k ) (3.23)
where ρ w j,k , c w j,k , V w j,k are, respectively, the density, specific heat, and volume of the layer k (the term ρ w j,k c w j,k V w j,k is associated with the thermal capacitance of the layer). The thermal conductance K j k,k±1 associated with the link between layers k and k ± 1, can be estimated by harmonic mean as:

K j k,k±1 = A j (L j,k±1 /2)/λ j,k±1 + (L j,k /2)/λ j,k
where A j is the area of the jth wall, L j,k denotes the thickness of the layer k, and λ j,k its thermal conductivity. Subscripts k + 1 and k -1 are related to the adjacent layers.

For the layers which are in direct contact with a room or the external environment, there are some conditions that can be described as follows:

• For a layer k within the wall j that is in contact with the room i and does not have an adjacent layer k ± 1, the term K j k,k±1 (T w j,k±1 -T w j,k ) in Equation (3.23) must 3.3 Distributed Building Temperature Control Under Power Constraints be replaced with h r i A j (T r i -T w j,k ). T r i represents the temperature of the ith room, and h r i its coefficient of convective heat transfer.

• For a layer k within the wall j that is in contact with the external environment and does not have an adjacent layer k ± 1, the term K j k,k±1 (T w j,k±1 -T w j,k ) must be replaced with h a A j (T a j,k -T w j,k ). T a j,k represents the ambient temperature perceived by the layer k of the wall j, and h a the external coefficient of convective heat transfer.

Thermal Model for a Room

A schematic representation of a room surrounded by walls is shown in Figure 3.4. In our model, we assume a uniform temperature in each room. Moreover, in order to formulate an appropriate thermal control strategy, we consider that each of the rooms has a sensor and an actuator (heater). With these assumptions, the temperature T r i of the ith room can be modeled as

(ρ a c a V r i ) Ṫ r i = j∈Ω r i h r i A j (T w j,k -T r i ) + x i + d i , (3.24)
where Ω r i is the set of walls adjacent to room i, T w j,k is the temperature of the layer within the wall j that is in direct contact with the room i; A j is the area of wall j; ρ a , c a are, respectively, the density and the specific heat of the air; V r i is the volume of the ith room, and h r i its coefficient of convective heat transfer; x i is the heating power supplied by the actuator; and d i is a thermal disturbance (e.g., due to opening a window or because of the presence of people who are generating heat). Similarly to Equation (3.23), the terms ρ a c a V r i and h r i A j in Equation (3.24) are associated with the thermal capacitance of the room i and with the thermal conductance of the junction between the room i and the wall j, respectively.

Unified Model

If we define the state vector

T = [T r 1 , ..., T r N , T w 1,1 , ..., T w 1,m 1 , ..., T w M,1 , ..., T w M,m M ] (3.25)

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

where N is the number of rooms of the building, M is the number of walls, and m j is the number of layers of the jth wall, Equations (3.23) and (3.24) can be unified into a single expression as follows [START_REF] Obando | Distributed Building Temperature Control with Power Constraints[END_REF].

θ i Ṫi = W j=1 α i,j (T j -T i ) + α i,a (T a i -T i ) + v i (x i + d i ), for i = 1, . . . , W, (3.26)
where T i is the ith element of the vector T (T i corresponds to the temperature of a component of the building, that can be either a room or a wall's layer); θ i > 0 is its thermal capacitance, α i,j ≥ 0 is the thermal conductance of the junction between elements i and j (we notice that α i,j = 0 if the junction between elements i and j does not exist. Moreover, α i,j = α j,i ); α i,a ≥ 0 is the thermal conductance of the junction between the ith element and the outside environment (again, α i,a = 0 if this junction does not exist); the variable v i takes the values of 1 if i = 1, . . . , n, and 0 otherwise; and T a i is the ambient temperature perceived by the ith element when it is in direct contact with the outside environment, if the ith element is not in direct contact with the outside environment, then T a i = 0. Finally, W denotes the dimension of the vector T, i.e., W is the total number of elements (rooms/layers) that comprise the building.

Since in the next sections we plan to use an analysis based on the result given in

1 θ W      L(Gb)+    α 1,a . . . α W,a      eT +    1 θ 1 . . . 1 θ W    I n 0 (W -n)×n e x (3.27)
where

e T = [T 1 -T * 1 , . . . , T W -T * W ] , e x = [x 1 -x * 1 , . . . , x n -x * n ] , L(G b ) is the graph Laplacian of G b , I n is the n × n identity matrix, and 0 (W -n)×n is a (W -n) × n matrix
with all its entries being zero. We notice that L(G) is symmetric given the fact that α i,j = α j,i . Before stating an important property of this system, let us introduce the following assumption. Assumption 3.3.1. For all i ∈ V b , there exists some j ∈ V b with α j,a > 0 (i.e., the element that corresponds to node j is in direct contact with the external environment), such that there is a path between i and j.

This implies that no component of the building is completely isolated from the outside environment. Without this assumption, one building's component could store energy indefinitely, which is not realistic. Taking into account this assumption, the system given in Equation (3.27) has the following property.

x = [x 1 -x * 1 , . . . , x n -x * n ] to the output e y = [T 1 -T * 1 , . . . , T n - T * n ] .
Proof. We use the following positive definite storage function

V 1 (e T ) = 1 2 W i=1 θ i e 2 T i (3.28)
where e T i is the ith element of the vector e T , i.e., e T i = T i -T * i , for i = 1, . . . , W . The derivative of V 1 along the trajectories of the system stated in Equation (3.27), is given by

V1 (e T ) = -e T L(G b )e T + W i=1 α i,a e 2 T i + n i=1 e T i e x i (3.29)

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

where e x i is the ith element of the vector e x , i.e., e x i = x ix * i , for i = 1, . . . , n. Let us prove that the expression in square brackets is positive definite.

Since L(G b ) is symmetric, we know that e T L(G b )e T = 1 2 i∈V b
j∈V b α i,j (e T je T i ) 2 [START_REF] Olfati-Saber | Consensus problems in networks of agents with switching topology and time-delays[END_REF], which means that both e T L(G b )e T and W i=1 α i,a e 2 T i are positive semidefinite functions. Furthermore, let us note that:

1. e T L(G b )e T = 0 if and only if e T i = e T j , for all α i,j = 0. Hence, e T L(G b )e T = 0

implies that e T i = e T k if there exists a path between i and k.

2.

W i=1 α i,a e 2 T i = 0 if and only if e T i = 0, for all α i,a = 0.

Mixing conditions 1 and 2, and given the fact that Assumption 3.3.1 holds, we conclude that the expression in brackets is positive definite, i.e., e T L(G b )e T + W i=1 α i,a e 2 T i = ψ 1 (e T ), where ψ 1 (e T ) > 0. Moreover, since n i=1 e T i e x i = e x e y , we have that V1 (e T ) = -ψ 1 (e T ) + e x e y , and the system is strictly passive (cf., Definition 2.2.2).

Let us note that we have taken the error of the control signals as the input of our system, and the error of the rooms' temperatures as the output (this fact is since, according to Equation (3.25), the ith output of the system corresponds to the error between the temperature of the ith room and its respective setpoint, for i = 1, . . . , n).

We do not include the temperatures of the walls' layers in the output because they are not relevant for any control strategy. The property pointed out in Proposition 3.3.1 is appealing because passivity guarantees the stability of the closed loop response when we add the multi-agent based controller described in (3.4), (3.6) . Indeed, we notice that the building thermal system (3.27) satisfies A3, i.e, it satisfies all the assumptions imposed on the controlled plant in Theorem 3.2.1.

Problem Statement

The main goal of a building temperature control strategy is to keep the temperature of the building's rooms close to a reference, we denote the reference of the ith room by T s i . The setpoints are not necessarily the same for all the rooms, given the fact that they depend on the preferences of the occupants, the kind of room and its usage. In order to reach the aforementioned goal, we assume that each room of the building has 3.3 Distributed Building Temperature Control Under Power Constraints an actuator (heater). These actuators have some constraints depending on the kind of heating system installed in the building.

A typical configuration, which is common in several buildings, is the central heating system. In this configuration, the total heating power Q is generated in one place (e.g., a furnace room) and then it is distributed to the N rooms that comprise the building.

Hence, we have the following constraint, n i=1

x i + x n+1 = Q, (3.30)
where x i is the heating power consumed in the ith room, for i = 1, . . . , n; and x n+1 is the remaining power, i.e., the amount of heating power that is not consumed in any room. We notice that the value of each x i is nonnegative. Moreover, it must not exceed the capacity of the corresponding actuator. Hence, we have that

0 ≤ x i ≤ x i , for i = 1, . . . , n x n+1 ≥ 0 (3.31)
where x i denotes the maximum capacity of the actuator of the ith room. If the amount of generated heating power is large enough, it is possible to reach all the temperature setpoints in the building's rooms by using individual controllers, for instance, we can use a proportional-integral (PI) controller per room. However, this strategy does not work properly if the available power is scarce, because the power allocation is not performed in an optimal way (e.g., see [19,[START_REF] Obando | Building Temperature Control based on Population Dynamics[END_REF]). This condition is rare but could occur in buildings with non over-designed heating systems, i.e., when the heating system is not designed to handle worst case loads.

In the latter situation, a good alternative is to use a dynamic resource allocation strategy, where the objective is to distribute the generated heating power among the rooms, in such a way that the building's occupants perceive the same welfare. A measure of the success of the allocation process is given in [19,[START_REF] Clearwater | Saving energy using market-based control. Market-Based Control: A Paradigm for Distributed Resource Allocation[END_REF][START_REF] Obando | Building Temperature Control based on Population Dynamics[END_REF], where the authors propose to use the variance of the differences between the setpoints and the steady state temperatures across the rooms, as follows,

Var = 1 n n i=1 ((T s i -T * i ) -T s i -T * i ) 2 , (3.32)
where • denotes mean value. The smaller the variance, the more successful the allocation done by the control strategy.

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

We recall that buildings are generally large scale systems, which are comprised by a large number of rooms. Hence, the application of centralized control strategies is not suitable. In this regard, the distributed method proposed in Section 3.2 is an appealing alternative to address the proposed problem. Under this method, each room has a controller that makes decision in a coordinated way, and the coordination is made by means of a communication network (i.e., we do not require the inclusion of a centralized agent). In the following section, we provide a detailed description of the design of the distributed controller.

Distributed Control Strategy

The performance criterion given in Equation (3.32), suggests that a good strategy to solve the stated problem is to lead an agreement among the temperature errors in the rooms. In fact, if we assume that all the heating power Q is supplied to the rooms (i.e., n i=1 x i = Q), and if we omit the constraints (3.31), then our building temperature control problem is solved by applying the multi-agent approach (3.4),(3.6) , with y i = T i -T s i , and initial conditions satisfying n i=1 x i (0) = Q. In this case, notice that Theorem 3.2.1 guarantees that at equilibrium all the temperature errors T i -T s i are equal; therefore, the variance (3.32) goes to zero. Nonetheless, the assumptions described above are not realistic. First, in situations where there is enough power to reach the temperature references in all the rooms, the power supplied to the rooms is generally lower than Q, i.e., n i=1 x i < Q. Second, constraints (3.31) cannot be removed since they are closely related to the operation of the actuators. Hence, we need to devise a control strategy that addresses all the requirements stated in the problem formulation.

In order to do this, we propose to preserve the multi-agent based controller (3.4), but with some modifications on the resource allocation protocol (3.6).

The first modification is associated with the inclusion of constraints (3.31) by means of memoryless functions. To this end, we change the consensus variable from

y i = T i -T s i to f i (T i , x i ) = (T i -T s i ) + ϕ i (x i ), for all i = 1, . . . , n, (3.33)
where > 0 is a small positive constant, and ϕ i (x i ) has the following properties:

• i) ϕ i (x i
) is a monotonically increasing continuous function defined in (0, x i ),

• ii)

ϕ i (x i ) → -∞, when x i → 0, 3.3 Distributed Building Temperature Control Under Power Constraints • iii) ϕ i (x i ) → +∞, when x i → x i .
We can think of ϕ i (x i ) as the derivative of a convex barrier function that penalizes the control signal when it tries to deviate outside its feasible domain (we refer the reader to [START_REF] Nocedal | Numerical optimization[END_REF] for a detailed description on barrier functions). In this regard, the constant should be designed in such a way that the effect of the function ϕ i (x i ) is minimized when the control signal is far from the boundary of its allowable limits.

The second modification allows the rooms to consume less heating power than Q, i.e., n i=1 x i ≤ Q. This is done by including an additional agent (n + 1) that manages the remaining power x n+1 . Thus, the modified multi-agent system is characterized by

the graph G c = {V c , E c , A c }, where V c = {1, . . . , n + 1}, E c ⊂ V c × V c , and A c = [a ij ] is the adjacency matrix.
After these changes, we can rewrite the dynamics of the multi-agent based controller as follows, ẋi =

j∈N i a ij (f j -f i ), for all i = 1, . . . , n + 1, (3.34)
where

N i = {j ∈ V c : (i, j) ∈ E c } is the neighborhood of node i in G c .
We notice the fact that N i depends on the communication topology (given by G c ), and not on the building's structure (given by G b ). Furthermore, since we have included the additional node (n + 1), we define the corresponding function f n+1 as follows,

f n+1 (x n+1 ) = ϕ n+1 (x n+1 ) (3.35)
where ϕ n+1 (x n+1 ) is a strictly increasing continuous function defined in (0, +∞), which tends to -∞ when x n+1 tends to 0. Notice that, since the consensus variables f i , . . . , f n+1 depends on the output of the building thermal system, and on the memoryless term If the initial conditions of the controller satisfy the power constraint (3.30), i.e., n+1 i=1 x i = Q, and the communication graph G c is connected, then the variable f asymptotically reaches consensus. Moreover, the control signals satisfy (3.30) for all t ≥ 0.

ϕ i (x i ) that

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

Proof. It is sufficient to show that the closed-loop system (3.26), (3.34) , with f = [f 1 , . . . , f n+1 ] given by Equations (3.33)

= -diag 1 θ 1 , . . . , 1 θn     L(G b )+     α 1,a . . . α W,a         A b T * +diag 1 θ 1 , . . . , 1 θn I n 0 (W -n)×n x * + d, (3.36)
where d is a vector of constants that describes the influence of the ambient temperature and the external disturbances on the building thermal system. Hence,

T * = A -1 b I n 0 (W -n)×n x * + A -1 b diag -1 1 θ 1 , . . . , 1 θ n d d . (3.37)
The matrix A b is invertible since it is positive definite (see the proof of Proposition 3.3.1).

On the other hand, the multi-agent network (3.34) must satisfy at equilibrium

-L(G c )f (T * , [x * , x * n+1 ] ) = 0, (3.38)
where

f (T * , [x * , x * n+1 ] ) = I n 0 n×(W -n) 0 1×n 0 1×(W -n) T * -[T s 1 , . . . , T s n , 0] + ϕ 1 (x * 1 ), . . . , ϕ n+1 (x * n+1
) . Since G c is connected by assumption, the null space of L(G c ) is equal to the linear span of 1 (cf., Theorem 2.1.1). Therefore, Equation (3.38) is satisfied if and only if

ξ1 = I n 0 n×(W -n) 0 1×n 0 1×(W -n) T * -[T s 1 , . . . , T s n , 0] + ϕ 1 (x * 1 ), . . . , ϕ n+1 (x * n+1 ) ,

Distributed Building Temperature Control Under Power Constraints

where ξ ∈ R. Replacing the value of T * given in Equation (3.37), we obtain

ξ1 = A -1 b n 0 1×n x * + ϕ 1 (x * 1 ), . . . , ϕ n+1 (x * n+1 ) + d , (3.39)
where (A -1 b ) n is the submatrix formed by the n first rows and the n first columns of matrix A b , and d is a vector of constants, which is defined by

d = I n 0 n×(W -n) 0 1×n 0 1×(W -n) d -[T s 1 , . . . , T s n , 0] .
According to the Karush-Kuhn-Tucker first order necessary conditions, it follows from (3.39), that finding an equilibrium point of (3.26),(3.34) that satisfies the power constraint (3.30) is equivalent to solve the following optimization problem min

(x * ,x * n+1 ) x * A -1 b n x * + d [x * , 0] + Φ(x * , x * n+1 ) subject to n+1 i=1 x i = Q (3.40)
where Φ(x * , x * n+1 ) satisfies ∇Φ = [ϕ 1 (x * 1 ), . . . , ϕ n+1 (x * n+1 )] , and (A -1 b ) n is positive definite since it is a strictly diagonally dominant matrix (see Equation (3.36)). Notice that, the objective function of (3.40) is strictly convex since (A -1 b ) n is positive definite and the gradient of Φ is comprised of monotonically increasing functions. Therefore, the optimization problem (3.40) has a unique solution. This implies that (3.26), (3.34) has an equilibrium point that satisfies the power constraint (3.30). Moreover, this equilibrium point is unique.

• A1 holds by assumption.

• A2 holds since the closed loop system (3.26), (3.34) has an equilibrium point satisfying n+1 i=1 x * i = Q, and the initial conditions satisfy n+1 i=1 x i (0) = Q by assumption.

• A3: The fact that (3.26) expressed in error coordinates is strictly passive, from the

input [x 1 -x * 1 , . . . , x n+1 -x * n+1 ] to the output [T 1 -T * 1 , . . . , T n -T * n , 0
] , follows from Proposition 3.3.1. Moreover, the storage function is radially unbounded.

• A4: Notice that the memoryless term ϕ i (x i ) satisfies the following property inside the feasible domain of x i :

(ϕ i (x i ) -ϕ i (x * i ))(x i -x * i ) ≥ 0, for all i = 1, . . . , n + 1.
(3.41)

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

This is true since ϕ i (x i ) is monotonically increasing. Therefore, the memoryless function [(ϕ 1 (x 1 )ϕ 1 (x * 1 )), . . . , (ϕ n+1 (x n+1 )ϕ n+1 (x * n+1 ))] is passive with respect to the error vector [x 1x * 1 , . . . , x n+1x * n+1 ] .

• A5: Since the ϕ i (x i ) is monotonically increasing, for all i = 1, . . . , n + 1, then

ϕ i (x) -ϕ i (x * i ) = 0 if and only if x i = x * i . Thus, Assumption 3.2.3 holds.
In conclusion, all assumptions of Theorem 3.2.3 hold. This proves the statement in Proposition 3.3.2.

Let us analyze the importance of Proposition 3.3.2.

Power Constraints

A direct consequence of Proposition 3.3.2 is that the proposed resource allocation algorithm prevents the calculated control signal to violate the constraint given in Equation (3.30), i.e, the algorithm guarantees that the sum of the powers allocated in the rooms of the building plus the remaining power (which is not used to heat the rooms) is always equal to the generated heating power Q.

On the other hand, ϕ i (•) acts as a barrier function that forces the ith control signal x i to satisfy the constraint (3.31). In order to explain the effect of ϕ i (•), let us consider the following situation: if x i is close enough to its upper bound x i , then the value of the corresponding function f i is higher than the value of the other functions, this happens given the fact that ϕ i (x i ) (which is a component of f i ) is a monotonically increasing continuous function that tends to +∞ when x i approximates to its upper bound. When the above occurs, ẋi is negative (according to Equation (3.34)), and therefore x i decreases and does not exceed the upper bound x i . The opposite happens when x i is close to its lower bound.

Summarizing, satisfaction of constraints (3.30) and (3.31) is guaranteed under the proposed control strategy.

Equilibrium Points

In this subsection, we use the result given in Proposition 3.3.2 to characterize the equilibrium points of the building thermal system controlled via the algorithm proposed in Proposition 3.3.2 entails that, at equilibrium, all the functions f 1 , . . . , f n+1 , are equal. For a better understanding of the implications of this feature, we define fi (T i ) =

f i (T i , x i ) -ϕ i (x i ) = T i -T s
i , as a function that only depends on the temperature of the ith room and its corresponding setpoint. Furthermore, following the same definition, we have that fn+1 = f n+1 (x n+1 )ϕ n+1 (x n+1 ) = 0. Using the auxiliary functions fi , . . . , fn+1 , let us analyze the properties of the equilibrium point (T * , [x * , x * n+1 ] ).

• Normal conditions: if for all i = 1 . . . , n + 1, x * i is inside its feasible domain and far from its allowable limits (which are given in Equation (3.31)), then the value of ϕ i (x * i ) can be neglected (given the fact that is small). Therefore, we have that f i ≈ fi , for all i = 1 . . . , n + 1; and according to Proposition 3.3.2,

f i (T * i , x * i ) = f n+1 (x * n+1
), for all i = 1, . . . , n. Hence, fi (T * i ) = T * i -T s i ≈ 0, for all i = 1, . . . , n. This means that the temperature setpoints can be reached by using the proposed control strategy, as long as the power required in each room satisfies its corresponding constraints.

• Scarce heating power: another scenario emerges when, at steady state, all the control signals are inside their feasible domains and far from their allowable limits, except the remaining power, i.e., x * n+1 → 0. Let us note that, in this case, ϕ i (x * i ) can be neglected for i = 1, . . . , n, but not for i = n + 1, where we have that ϕ n+1 (x * n+1 ) < 0. Therefore, applying Proposition 3.3.2, we can conclude that

f i = T * i -T s i ≈ ϕ n+1 (x * n+1
) < 0, for all i = 1, . . . , n. This implies that, despite the temperatures of the rooms do not reach the corresponding setpoints due to lack of heating power, all the steady state errors (T * i -T s i , for each i = 1, . . . n) are almost equal, minimizing the variance given in Equation (3.32).

• General case: Finally, let us assume a general case, where in steady state, all the control signals are inside their feasible domains, but only a set of them, denoted by X, is far from the allowable limits, i.e., X = {x * i : x * i is not close to 0, and x * i is not close to x i };

while other set of signals, denoted by X, is close to the lower bound, i.e.,

X = {x * i : x * i → 0};

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

and other set of signals, denoted by X, is close to the upper bound, i.e.,

X = {x * i : x * i → x i }.
As in the previous scenarios,

ϕ i (x * i ) ≈ 0 if x * i ∈ X.
However, this is not the case if x * i belongs to the other sets. Using the properties of the function ϕ i (•), it can be verified that

ϕ i (x * i ) < 0 if x * i ∈ X, while ϕ i (x * i ) > 0 if x * i ∈ X.
Hence, by using the result in Proposition 3.3.2, we can conclude that the following relationships hold at equilibrium:

-fi ≈ fj if x * i ∈ X and x * j ∈ X.
fi < fj if x * i ∈ X and x * j ∈ X.

fi > fj if x * i ∈ X and x * j ∈ X.

Moreover, since fi (T * i ) = T * i -T s i , for all i = 1, . . . , n, the above relationships imply that no heating power is assigned to the rooms that are warmer to the average, while the colder rooms receive the maximum allowable power. In this regard, the heating power is optimally allocated.

Summarizing, the equilibrium point of the closed loop system (3.26),(3.34) exhibits optimal properties under different scenarios of the building temperature control problem.

Besides, since this equilibrium point is asymptotically stable according to Proposition 3.3.2, we can conclude that the multi-agent approach (3.34) properly solves the building temperature control problem described in Section 3.3.3.

Simulation Results

In order to illustrate the performance of the proposed control technique, we simulate the thermal behavior of a building comprised of 50 rooms that are arranged as in Figure 3.3.

The building has a central heating system of 130 KW, with one actuator per room. The maximum power that each actuator can provide is 3.25 KW. Moreover, the temperature setpoints are: 19 • C for rooms 1-17, 20 • C for rooms 18-34, and 21 • C for rooms 34-50.

In the simulation, we assume that the outdoor temperature changes according to Figure of the temperatures in each room from 7h to 20h, Figure 3.5c shows the total heating power used by the rooms, and Figure 3.5d shows the performance index related to the variance of the error across the rooms (cf., Equation (3.32)). It can be noticed that, from 7h-13h it is not possible to reach the setpoints despite the fact that all the available power is used to heat the rooms. However, the controller minimizes the variance of the error across the rooms to a value of zero, which implies that the building's occupants are perceiving the same welfare (i.e., the scarce power situation described in Section 3.3.4.2).

When the weather is warmer, i.e., around 13h-17h, the temperature in each room rises and reaches the corresponding setpoint, while the power consumption is reduced. In the last part of the simulation, we have a similar behavior as the one described in the first part, due to the decrease of the outdoor temperature. Finally, we compare the performance of the proposed technique with the most widely used method in building temperature control, i.e., PI control [START_REF] Nof | Springer Handbook of Automation[END_REF]. In order to do that, we implement one PI controller per room in the same scenario described above. Each

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

PI controller seeks to regulate the temperature of its corresponding room without any information on the state of the other rooms, i.e., their decisions are not coordinated.

Moreover, we notice the fact that the maximum power that each PI controller can use depends on both the actuator capacity and the consumption of the other rooms. The variance of the errors across the rooms with this technique is depicted in dashed line in Figure 3.5d. It can be noticed that individual PI controllers do not work properly when the outdoor temperature is lower and the available heating power is not enough to reach all the setpoints, i.e., from 7h-11h, and from 19h-20h; in that time intervals, the variance has a large value. This adverse effect is due to an inefficient power allocation.

Discussion

When the generated heating power is enough to reach the temperature setpoints of the building's rooms, the behavior of the closed loop system controlled via the proposed algorithm is the same as the behavior achieved with individual PI controllers, i.e., all the temperature references are reached without steady state error. However, under a scarce power scenario, PI controllers do not allocate the available resource in an optimal way, as shown in Figure 3.5d, where the variance of the temperature errors across the rooms has a high value at critical load hours. This fact is due to a lack of coordination among the controllers. The proposed algorithm does not have that problem because the communication network allows the agents to share information on the state of the neighboring rooms' temperatures. In this way, each agent is able to make decisions that improve the performance of the overall system, rather than the performance of a single room.

Optimal results can also be obtained by using model predictive control (centralized or distributed). However, although that technique offers wide formulation flexibility and its predictive features can enhance the energy saving, it suffers from high computational requirements and model dependence, which reduces its scalability, especially for large buildings with not negligible thermal coupling among rooms [START_REF] Petru-Daniel Moros ¸an | A distributed MPC strategy based on Benders' decomposition applied to multisource multi-zone temperature regulation[END_REF]. In this regard, the distributed methodology that we propose has the advantage of being formulated as a set of "quasi-linear" differential equations, whose implementation does not demand high computation resources. Hence, its scalability is substantially higher than that provided by receding horizon techniques.

3.4 Optimal Resource Allocation in a Network of Agents

Optimal Resource Allocation in a Network of Agents

The classical resource allocation problem arises when we have a scarce resource, and it is necessary to distribute it among a set of entities in order to optimize a given objective function. We talk about resource allocation over networks of agents when the resource allocated to each entity is managed by an agent that only has partial information of the overall optimization problem. In order to reach the optimal allocation, the agents coordinate their decisions by using a communication network with a given topology.

If we have a set of n entities and we assign one agent per entity, the aforementioned problem can be expressed mathematically as follows [11],

min φ(x) subject to n i=1 x i = X (3.42)
where φ : R n → R is a strictly convex and differentiable cost function, X ∈ R is a given constant related to the total amount of resource, x i is the amount of resource allocated to the ith entity, and x = [x 1 , . . . , x n ] . In the literature on distributed optimization [START_REF] Bemporad | Networked control systems[END_REF], it is common to assume that the ith agent knows the marginal cost associated with its corresponding entity, i.e., ∂φ ∂x i . Moreover, this information can be shared by the agents using the available communication network, which is modeled by the graph G = (V, E, A).

The local information structure imposed by the graph G must be considered as part of the problem statement. Therefore, the main objective is to formulate a distributed algorithm that can be implemented by each agent to solve (3.42) using only its available information. In this case, notice that the available information of the ith agent is given by ∂φ ∂x i , ∂φ ∂x j : j ∈ N i . A well-known method to address (3.42) that satisfies the local information constraint is formulated as follows [24,[START_REF] Xiao | Optimal scaling of a gradient method for distributed resource allocation[END_REF].

ẋi = j∈N i a i,j ∂φ ∂x j - ∂φ ∂x i , for all i = 1, . . . , n. (3.43)
Although the optimality of (3.43) has been proved in [START_REF] Cherukuri | Distributed generator coordination for initialization and anytime optimization in economic dispatch[END_REF], we are interested in analyzing this method form the perspective of passive memoryless systems given in Section x i = X. To this end, we use the following facts:

• Since φ is strictly convex, then this satisfies the inequality φ(x * ) > φ(x)+∇φ(x)(x *x), for all x = x * (see, e.g., [START_REF] Boyd | Convex optimization[END_REF]). Given the fact that x * is the global optimum of (3.42), then ∇φ(x)(xx * ) > 0, for all x = x * that satisfies n i=1 x i = X.

• According to the Karush-Kuhn-Tucker first order necessary conditions, ∇φ(x * ) belongs to the linear span of 1, i.e., ∇φ(x * ) = λ1, where λ ∈ R. Thus, ∇φ(x * )(x -

x * ) = λ(1 x -1 x * ). Since both x and x * satisfy the resource constraint, we can conclude that ∇φ(x * )(xx * ) = 0.

Combining these two facts, we obtain (∇φ(x) -∇φ(x * )) (xx * ) > 0, for all x = x * .

Thus, the considered memoryless system is passive. Furthermore, notice that this system satisfies Assumption 3.2.2 since (∇φ

(x) -∇φ(x * )) (x -x * ) = 0 only if x = x * .
In conclusion, the result in Theorem 3.2.2 can also be used to guarantee the convergence of (3.43) to the optimum of (3.42) as long as n i=1 x i (0) = X.

Inclusion of Inequality Constraints by means of Barrier Functions

We can take advantage of the convexity of φ to include in (3.42) inequality constraints of the form

x i ≤ x i ≤ x i , for all i = 1, . . . , n, (3.44)
where x i and x i are, respectively, the lower and upper bound of x i . This is done by using convex barrier functions that keep the decision variables strictly inside the feasible set. Barrier functions are continuous functions defined inside the feasible domain whose value increases to infinity when the variable approaches to the boundary of its feasible domain. In our specific case, we use a strictly convex barrier function for each inequality constraint (3.44), as follows [START_REF] Nocedal | Numerical optimization[END_REF]:

b i (x i ) = -ε ln (x i -x i ) + ln (x i -x i ) , ∀i = 1, . . . , n, (3.45)
where ε > 0 is a constant that minimizes the effect of the barrier function when the solution is far from the boundary of the feasible set. These barrier functions are added 52 3.4 Optimal Resource Allocation in a Network of Agents to the objective function φ(x) in Equation (3.42), i.e., our new objective function becomes

φ b (x) = φ(x) -ε n i=1 (ln (x i -x i ) + ln (x i -x i )).
Notice that, if φ is strictly convex, φ b remains strictly convex. Moreover, although the solution of the optimization problem that includes barrier functions is suboptimal for the original problem, this solution can be sufficiently close to the optimum depending on the value of ε. Furthermore, since b i only depends on x i , the algorithm (3.43) remains distributed if we replace ∂φ ∂x i by ∂φ b ∂x i , for all i = 1, . . . , n.

Initial Feasible Solution

The algorithm given in (3.43) needs an initial feasible solution of the corresponding optimization problem. Although finding this starting point is straightforward if the optimization problem does not include the inequality constraints (3.44) (in this case, we can take (x k (0) = X, x i (0) = 0, for all i = 1 . . . , n, i = k) as the initial feasible solution), it is not trivial if they are considered. If this is the case, a preliminary phase that provides us the required feasible starting point is required. To implement the preliminary phase, we use a continuous-time version of the the distributed constraint satisfaction algorithm described in [START_REF] Domınguez | Distributed algorithms for control of demand response and distributed energy resources[END_REF].

In order to maintain the distributed nature of our resource allocation problem, we assume that only one agent knows the total amount of available resource X. The main idea of the preliminary phase is to solve two auxiliary problems that use local information, and then combining the solutions to get the starting point for the main problem (3.42), (3.44) . Let β 1 i and β 2 i , two auxiliary variables associated with the ith node.

Using these variables, we establish the following auxiliary problems,

Problem 1 : min φ 1 (β 1 ) = β 1 β 1 , s. t. 1 β 1 = d - n i=1 x i Problem 2 : min φ 2 (β 2 ) = β 2 β 2 s. t. 1 β 2 = n i=1 (x i -x i ) (3.46)
where

β 1 = [β 1 1 , . . . , β 1n ] and β 2 = [β 2 1 , . . . , β 2n ] .
It can be shown that the optimal solutions of Problems 1 and 2 are,

β * 1 = d-n i=1 x i n 1 and β * 2 = n i=1 (x i -x i ) n
1, respectively.

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

Moreover, since φ 1 and φ 2 are strictly convex, these optimal solutions can be obtained by using the distributed algorithm described in Equation (3.43) with the following initial feasible solutions:

β 1 i (0) = d -x i if i = k -x i otherwise , β 2 i (0) = x i -x i , for all i = 1, . . . , n.
where we assume, without loss of generality, that the kth agent knows the value X. We notice the fact that the agents only require local information to obtain β 1 (0) and β 2 (0).

Using the optimal solutions of Problems 1 and 2 in Equation (3.46), we can obtain a feasible solution of the main problem (3.42), (3.44) , as follows,

x i (0) = x i + (x i -x i ) β * 1,i β * 2,i
, for all i = 1, . . . , n.

Let us note that if 0 ≤

β * 1,i β * 2,i
≤ 1, then this initial solution satisfies all the constraints of our resource allocation problem (3.42), (3.44) . Moreover if

β * 1,i β * 2,i < 0 or β * 1,i β * 2,i > 1,
the problem is infeasible.

Application to Electric Load Aggregation for Demand Response

In this section, we present an application of the optimal resource allocation in networks of agents. The considered application deals with the problem of electric load aggregation for demand response [START_REF] Obando | Distributed Resource Allocation Over Stochastic Networks: An Application in Smart Grids[END_REF]. This problem is motivated by the fact that the active participation of demand is a central issue in the smart grid concept, in this regard, aggregators play a fundamental role. Aggregators are mediators between consumers, who provide their demand flexibilities, and the other participants of electricity markets (e.g., distribution/transmission system operators, retailers, electricity brokers). So, when a demand deviation is required (e.g., due to imbalances in the distribution network), aggregators can provide it by acting on the loads of their clients. In return to these possible load curtailments, aggregators' clients pay lower electricity prices. Moreover, since each aggregator generally has a large number of clients, distributed management strategies are becoming more relevant.

The operation of an aggregator can be viewed as a resource allocation process, where the total amount of load deviation that the aggregator has to provide at certain time,

Optimal Resource Allocation in a Network of Agents

needs to be split between their clients (consumers). One possibility is to do an allocation of the necessary load curtailment considering only technical constraints as it is described in [START_REF] Domınguez | Distributed algorithms for control of demand response and distributed energy resources[END_REF]. However, a better opportunity emerges given the fact that the aggregator can take advantage of the aforementioned allocation process to generate other benefits. For instance, it is important to encourage the clients to stay involved in demand response programs by minimizing the impact that they perceive due to the load curtailments.

To formally describe this situation, let us assume that the aggregator manages a set of n consumers. The authors in [START_REF] Ygge | Power load management as a computational market[END_REF], propose to model the benefit perceived by the ith consumer as follows,

u i (p i ) = k i,1 -k i,2 e -k i,3 p i -k i,4 p i , (3.47)
where p i is the electric power being used by the ith consumer, k i,1 , k i,2 , k i,3 , and k i,4

are positive parameters that depend on the consumer's preferences. Moreover, if the consumer suffers a load curtailment x i , the new benefit that the consumer perceives is

u i (p i -x i ).
Therefore, the cost of the load curtailment for the ith consumer is given by

c i (x i ) = u i (p i ) -u i (p i -x i ). (3.48)
It can be shown that c i is a strictly convex function of x i .

In this regard, the aggregator must satisfy the required load deviation D L , minimizing at the same time, the total cost perceived by its clients, i.e., φ(x) = n i=1 c i (x i ). Mathematically, this problem is formulated as follows

min φ(x) = n i=1 c i (x i ) subject to n i=1 x i = D L 0 ≤ x i ≤ p i , for all i = 1, . . . , n, (3.49)
where x i is the load curtailed from the ith consumer. The inequality constraints model the fact that the aggregator cannot curtail more load than that used by the consumer.

As previously mentioned, an aggregator can have a large number of clients, e.g., small consumers. Hence, a centralized operation is not the best option. Moreover, in a smart grid environment, it is assumed that each consumer has computation and communication capabilities. For these reasons, the optimization problem given in Equation (3.49) can be appropriately approached by means of a multi-agent system that uses the distributed protocol given in Equation (3.43). In order to do that, we assign one agent per client.

Each agent knows (or estimates) the function describing the benefit of the corresponding

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

client, and its current consumption. The coordination of the agents' decisions is done through a communication network with a topology given by a connected graph G.

To illustrate the performance of the distributed resource allocation algorithm, we simulate the operation of an aggregator that manages 8 clients, and uses a communication network with single-path topology (see Figure 3 x 5 (0) = 0.30 MW, x 6 (0) = 0.41 MW, x 7 (0) = 0.23 MW, x 8 (0) = 0.12 MW. It can be noticed that the steady state solution satisfies the technical constraints of the problem, i.e., n i=1 x * i = D L = 2.7 MW, and 0 ≤ x * i ≤ p i , for i = 1, . . . , n. Furthermore, Figure 3.9 shows that the total cost function decreases and finally reaches a reduction of 13.7%.

In this regard, it is worth noting that the impact of load curtailment on the aggregator's clients is reduced by minimizing the cost function stated in Equation (3.48). This can be seen as an incentive for clients to keep involved in demand response programs.

Discussion

The classic resource allocation problem (3.42) have been widely studied in the literature on optimization (see e.g., [11] and the references therein). In fact, the properties of the algorithm stated in Equation (3.43) are well-known [24, [START_REF] Xiao | Optimal scaling of a gradient method for distributed resource allocation[END_REF][START_REF] Cherukuri | Distributed generator coordination for initialization and anytime optimization in economic dispatch[END_REF]. Therefore, the optimality analysis of (3.43) cannot be considered as an original contribution of Section 3.4. The purpose of this section has been to show a novel perspective of the classic resource allocation problem from passivity theory. We have proven that there exists an equivalence between (3.42) and a passive memoryless system that is given by the gradient of the cost function as long as this cost function is convex. Thus, we have shown that the result given in Theorem 3.2.2 is applicable to the solution of the optimal resource allocation problem via the algorithm (3.43). On the other hand, Section 3.4.2 have described a constraint satisfaction algorithm that can be employed to initialize (3.43) when inequality constraints of the form (3.44) are added to the problem formulation.

Notice that, the same initialization phase is required in the multi-agent scheme proposed to address the building temperature control problem (see Section 3.3). Regarding the inequality constraints (3.44), it has been shown that they can be satisfied by introducing convex barrier functions. However, this approach decreases the accuracy of the solution.

DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS

Indeed, in the next chapter, we will show that the adverse effect caused by the influence of barrier function is more pronounced in large-scale applications. A solution to this problem will also be proposed in the next chapter. Finally, we have described an application to the distributed management of the customers of an aggregating entity in a demand response environment. We notice the fact that demand response programs, where consumers allow an aggregator to curtail part of their load, are not exclusively applied by commercial and industrial users nowadays.

They have been extended to residential consumers. This fact has greatly increased the number of clients that aggregators need to manage. In this regard, the development and analysis of distributed techniques, as the one shown in Section 3.4.3, is important to properly deal with large number of users. Hence, these techniques facilitate the inclusion of new participants in demand response programs.

Optimal Resource Allocation with Lower-Bound Constraints

In this chapter, we go deeper into the study of the optimal resource allocation problem stated in Section 3.4. Specifically, we focus our attention on solution methods for problems that include lower-bound constraints in the resource allocated to each entity.

The capability of handling this kind of constraints is crucial in a large number of practical applications, e.g., in [START_REF] Conrad | Resource economics[END_REF]15,16,[START_REF] Cherukuri | Distributed generator coordination for initialization and anytime optimization in economic dispatch[END_REF][START_REF] Obando | Building temperature control based on population dynamics[END_REF][START_REF] Pantoja | Dispatch of distributed generators under local-information constraints[END_REF][START_REF] Tan | Consensus based approach for economic dispatch problem in a smart grid[END_REF], where it is required to capture the non-negativity of the resource allocated to each entity.

The algorithm proposed in this chapter is motivated by the fact that, although lowerbound constraints can be addressed by using barrier functions (cf., Section 3.4.1), the inclusion of these functions affect the accuracy of the solution, especially for large-scale problems [74]. Therefore, it is necessary to develop other methods whose accuracy does not depend on the problem scale.

Problem Statement

We recall that the resource allocation problem arises when there is a limited amount of a certain resource (e.g., electric power, computing capacity, or execution time), and it is necessary to establish an optimal distribution policy between some entities (e.g., loads, processors, or controllers). If we impose lower-bound constraints to the resource allocated to each entity, the problem is formulated as follows [11],

min x φ(x) := n i=1 φ i (x i ) (4.1a) subject to n i=1
x i = X (4.1b)

x i ≥ x i , for all i = 1, . . . , n, (4.1c)
where x i ∈ R is the resource allocated to the ith entity; x = [x 1 , . . . , x n ] ; φ i : R → R is a strictly convex and differentiable cost function; X is the available resource; and x i , is the lower bound of x i , i.e., the minimum amount of resource that has to be allocated to the ith entity.

Given the fact that we are interested in distributed algorithms to solve the problem stated in Equation (4.1), we consider a multi-agent network, where the ith agent is responsible for managing the resource allocated to the ith entity. Moreover, we assume that the agents have limited communication capabilities, so they can only share information with their neighbors. This constraint can be represented by a graph G = {V, E, A} as it was explained in the previous chapter.

Avoiding the lower-bound constraints (4.1c), Karush-Kuhn-Tucker (KKT) conditions establish that at the optimal solution x * = [x * 1 , . . . , x * n ] of the problem given in Equation (4.1a-4.1b), the marginal costs φ i (x i ) = dφ i dx i must be equal, i.e., φ i (x * i ) = λ, for all i = 1, . . . , n, where λ ∈ R. Hence, a valid alternative to solve (4.1a-4.1b) is the use of consensus methods. For instance, we can employ the algorithm presented in Section 3.4, which is described as follows:

ẋi = j∈N i φ j (x j ) -φ i (x i ) , for all i ∈ V. (4.2)
This algorithm has two main properties: i) at equilibrium, φ i (x * i ) = φ j (x * j ) if the nodes i and j are connected by a path; ii)

n i=1 x * i = n i=1 x i (0)
, where x i (0) is the initial condition of x i . Therefore, if the graph G is connected and the initial condition is feasible (i.e., n i=1 x i (0) = X), x asymptotically reaches the optimal solution of (4.1a-4.1b) under (4.2). However, the same method cannot be applied to solve (4.1) (the problem that considers lower bounds in the resource allocated to each entity) since some feasibility issues related with the constraints (4.1c) arise.

Center-Free Resource Allocation Algorithm

In the following section, we propose a novel method that extends the algorithm in Equation (4.2) to deal with the individual inequality constraints given in Equation (4.1c).

Moreover, we show that this method satisfies some passivity properties. Therefore, its analysis can be done by using the result given in Theorem 3.2.3.

Center-Free Resource Allocation Algorithm

Resource Allocation Among a Subset of Nodes in a Graph

First, we consider the following subproblem: let G = {V, E, A} be a graph comprised by a subset of active nodes V a and a subset of passive nodes V p , such that V a V p = V.

A certain amount of resource X has to be split among those nodes to minimize the cost function φ(x) subject to each passive node is allocated with its corresponding lower bound x i . Mathematically, we formulate this subproblem as follows.

min x φ(x) (4.3a) subject to n i=1
x i = X (4.3b)

x i = x i , for all i ∈ V p . (4.3c)
Feasibility of (4.3) is guaranteed by making the following assumption.

Assumption 4.2.1. At least one node is active, i.e., V a = ∅.

According to KKT conditions, the active nodes have to equalize their marginal costs at the optimal solution. Therefore, a consensus among the active nodes is required to solve (4.3). Nonetheless, classic consensus algorithms, as the one given in Equation (4.2), cannot be used directly. For instance, if all the nodes of G apply (4.2) and G is connected, the marginal costs of both passive and active nodes are driven to be equal in steady state. This implies that the resource allocated to passive nodes can violate the constraint (4.3c). Besides, if the resource allocated to passive nodes is forced to satisfy (4.3c) by setting x * i = x i , for all i ∈ V p , there is no guarantee that the new solution satisfies (4.3b). Another alternative, is to apply (4.2) to only active nodes in this case, the neighborhood of node i ∈ V a in Equation (4.2) has to be taken as {j ∈ V a : (i, j) ∈ E}, and the initial condition must satisfy i∈Va x i (0) = X -i∈Vp x i . However, the subgraph formed by the active nodes is not necessarily connected although G is connected.

OPTIMAL RESOURCE ALLOCATION WITH LOWER-BOUND CONSTRAINTS

Hence, marginal cost of active nodes are not necessarily equalized at equilibrium, which implies that the obtained solution is sub-optimal. In conclusion, modification of (4.2) to address (4.3) is not trivial. In order to deal with this problem, we propose the following algorithm:

ẋi = j∈N i (y j -y i ), for all i ∈ V (4.4a) ẋi = (x i -x i ) + j∈N i (y j -y i ), for all i ∈ V p (4.4b) y i = φ i (x i ) if i ∈ V a φ i (x i ) + xi if i ∈ V p .
(4.4c)

In the same way as in (4.2), the variables {x i , i ∈ V} in Equation (4.4) correspond to the resource allocated to both active and passive nodes. Notice that we have added auxiliary variables {x i , i ∈ V p } that allow the passive nodes to interact with their neighbors taking into account the constraint (4.3c). On the other hand, the term j∈N i (y jy i ), in Equations (4.4a)-(4.4b), leads to a consensus among the elements of the vector y = [y 1 , . . . , y n ] , which are given in Equation (4.4c). For active nodes, y i only depends on the marginal cost φ i (x i ), while for passive nodes, y i depends on both the marginal cost and the state of the auxiliary variable xi . Therefore, if the ith node is passive, it has to compute both variables x i and xi . Furthermore, it can be seen that, if all the nodes are active, i.e., (V a = V), then the proposed algorithm becomes the one stated in Equation (4.2).

Notice that the ith node only needs to know y i and the values {y j : j ∈ N i } to compute j∈N i (y jy i ) in (4.4a)-(4.4b). In other words,

L(G)y =   j∈N 1 (y j -y 1 ), . . . , j∈Nn (y j -y n )  
is a distributed map over the graph G (cf., Definition 2.1.1). This implies that the dynamics given in Equation (4.4) can be computed by each node using only local information.

In fact, the message that the ith node must send to its neighbors is solely composed by the variable y i .

Center-Free Resource Allocation Algorithm

Feasibility

Let us prove that, under the multi-agent system proposed in Equation (4.4), x(t) satisfies the first constraint of the problem given by Equation (4.3), for all t ≥ 0, provided that n i=1 x i (0) = X.

1 ẋ = -1 L(G)y. Since G is undirected, 1 L(G) = L(G)1 = 0. Therefore, ∆ = 0.
The above lemma does not guarantee that x(t) is always feasible because of the second constraint in Equation (4.3), i.e., x i = x i , for all i ∈ V p . However, it is possible to prove that, at equilibrium, this constraint is properly satisfied.

Equilibrium Point

The next proposition characterizes the equilibrium point of the multi-agent system given in Equation (4.4). Proposition 4.2.1. Let Assumption 4.2.1 holds. If G is connected, the system in Equation (4.4) has an equilibrium point (x * , {x * i , i ∈ V p }) that satisfies the following conditions:

φ i (x * i ) = λ, for all i ∈ V a , where λ ∈ R is a constant; x * i = x i , for all i ∈ V p ; and n i=1 x * i = X. Moreover, x * i = λ -φ i (x * i ), for all i ∈ V p .
Proof. Since (4.3) is feasible (due to the fact that Assumption 4.2.1 holds), there exists x * that minimizes φ(x) and satisfies (4.3b)-(4.3c), i.e., n i=1 x * i = X, and x * i = x i , for all i ∈ V p . Notice that, applying the KKT first order necessary conditions to (4.3), x * must also satisfy φ i (x * i ) = λ, for all i ∈ V a , and some λ ∈ R. Let x * i = λφ i (x * i ), for all i ∈ V p . To finish the proof, we need to show that the point (x * , {x * i , i ∈ V p }) defined above is an equilibrium of (4.4).

It follows from Equation (4.4a) that y * i = λ, for all i ∈ V, where y * i is the value of y i at (x * , {x * i , i ∈ V p }). Thus, replacing this value in (4.4a), we have that, at (x * , {x * i , i ∈ V p }), ẋi = 0, for all i ∈ V. Besides, replacing the value of y * i in Equation (4.4b), and using the fact that x * i = x i , for all i ∈ V p , we obtain that ẋi = 0, for all i ∈ V p .

OPTIMAL RESOURCE ALLOCATION WITH LOWER-BOUND CONSTRAINTS

Remark 4.2.1. Proposition 4.2.1 states that the system (4.4) has an equilibrium point (x * , {x * i , i ∈ V p }) for which the active nodes equalize their marginal costs, while each passive node is allocated with an amount of resource equal to its corresponding lower bound. Moreover, x * satisfies the constraints (4.3b)-(4.3c). In conclusion, this proposition guarantees that x * minimizes the optimization subproblem given in Equation (4.3). Additionally, notice that the values {x * i , i ∈ V p } are equal to the KKT multipliers associated with the constraint (4.3c).

Convergence Analysis

According to Remark 4.2.1, to show that the algorithm proposed in Equation (4.4) converges to the optimal solution of the subproblem (4.3), it is sufficient to prove that the equilibrium point given in Proposition 4.2.1 is asymptotically stable (AS). Before doing so, let us notice that the algorithm (4.4) is a dynamical system that can be modeled by the feedback interconnection of (4.4a) and (4.4b)-(4.4c). Furthermore, since Equation (4.4a) is equivalent to the model of the multi-agent network (3.4) that applies the resource allocation algorithm (3.6), the convergence analysis can be performed by using the result on passivity stated in Proposition 3.2.2 of the previous chapter.

After these preliminaries, let us prove that the dynamics in Equation (4.4) converge to the optimal solution of (4.3), provided that each φ i (x i ) is strictly convex. Proposition 4.2.2. Assume that φ i (x i ), is a strictly convex cost function, for all i ∈ V. If G is connected, n i=1 x i (0) = X, and Assumption 4.2.1 holds, then x(t) converges to x * under Equation (4.4), where x * is the solution of the optimization problem stated in Equation (4.3), i.e, x * is the same given in Proposition 4.2.1. Furthermore, xi converges to x * i , for all i ∈ V p .

Proof. Let us express our multi-agent system in error coordinates with respect to the equilibrium point

(x * , {x * i , i ∈ V p }), as follows ė = -L(G)e y ėi = e i -(L(G)e y ) i , for all i ∈ V p e y i = φ i (x i ) -φ i (x * i ) if i ∈ V a φ i (x i ) -φ i (x * i ) + êi if i ∈ V p , (4.5)
where L(G) is the graph Laplacian of G; e i = x ix * i , and e y i = y iy * i , for all i ∈ V; êi = xix * i , for all i ∈ V p ; e = [e 1 , . . . , e n ] ; e y = [e y 1 , . . . , e yn ] ; and (L(G)e y ) i represents the ith element of the vector L(G)e y . Since Assumption 4.2.1 holds, V a = ∅. Let k be an active node, i.e., k ∈ V a , and let e k , e k y be the vectors obtained by removing the kth element from vectors e and e y , respectively. We notice that, according to Lemma 4.2.1, e k (t) = -i∈ν,i =k e i (t), for all t ≥ 0. Therefore, Equation (4.5) can be expressed as

ėk = -L k (G)e k y -l kr k e y k e k = -i∈ν,i =k e i ėi = e i -L k (G)e k y + l kr k e y k i , for all i ∈ V p e y i = φ i (x i ) -φ i (x * i ) if i ∈ V a φ i (x i ) -φ i (x * i ) + êi if i ∈ V p , (4.6)
where L k (G), and l kr k are defined in Lemma 2.1.1. In order to prove that the origin of the above system is AS, let us define the following Lyapunov function:

V = 1 2 e k L k (G) -1 e k + 1 2 i∈Vp (e i -êi ) 2 . ( 4.7)
The function V is positive definite since G is connected (the reason of this fact is that, according to Lemma 2.1.1, L k (G) and its inverse are positive definite matrices if G is connected). Notice that the first term of (4.7) is equivalent to the storage function (3.9) defined in the proof of Proposition 3.2.2. Therefore, using the result stated in that proposition, we obtain that the derivative of V along the trajectories of the system stated in Equation (4.6) is given by, V = -e e y -i∈Vp e i (e iêi )

Replacing e y by the expression given in (4.6), we have that

V = - n i=1 e i (φ i (x i ) -φ i (x * i )) - i∈Vp e i êi + i∈Vp e i (ê i -e i ) = - n i=1 (x i -x * i )(φ i (x i ) -φ i (x * i )) - i∈Vp e 2 i ,
where φ i is strictly increasing given the fact that φ i is strictly convex, for all i ∈ V. Therefore, (

x i -x * i )(φ i (x i ) -φ i (x * i )) ≥ 0, for all i ∈ V, and thus V ≤ 0. Since V does not depend on {ê i , i ∈ V p }, it is negative semidefinite. Let S = {{e i , i ∈ V}, {ê i , i ∈ V p } : V = 0}, i.e., S = {{e i , i ∈ V}, {ê i , i ∈ V p } : e i = 0
, for all i ∈ V}. Given the fact that G is connected and V = V p (by Assumption 4.2.1), then ė = 0 iff e y = 0 (see Equation (4.5)). Therefore, the only solution that stays identically in S is the trivial solution, i.e., e i (t) = 0, for all i ∈ V, êi (t) = 0, for all i ∈ V p . Hence, we can conclude that the origin is AS by applying the invariance principle given in Theorem 2.2.1.

In summary, we have shown that the algorithm described in Equation (4.4) asymptotically solves the subproblem in Equation (4.3), i.e., (4.4) guarantees that the resource allocated to each passive node is equal to its corresponding lower bound, while the remaining resource X -i∈Vp x i is optimally allocated to active nodes.

Optimal Resource Allocation with Lower Bounds

Now, let us consider our original problem stated in Equation (4.1), i.e., the resource allocation problem that includes lower bound constraints. Let x * = [x * 1 , . . . , x * n ] be the optimal solution of this problem. Notice that, if we know in advance which nodes will satisfy the constraint (4.1c) with strict equality after making the optimal resource allocation process, i.e., I := {i ∈ V : x * i = x i }, we can mark these nodes as passive and reformulate (4.1) as a subproblem of the form (4.3). Based on this idea, we propose a solution method for (4.1), which is divided in two stages: in the first one, the nodes that belong to I are identified and marked as passive; in the second one, the resulting subproblem of the form (4.3) is solved by using (4.4). Protocol (4.4) can also be used in the first stage of the method as follows: in order to identify the nodes that will satisfy (4.1c) with strict equality at the optimal allocation, we start marking all nodes as active and apply the resource allocation process given by (4.4). The nodes that are allocated with an amount of resource below their lower bounds at equilibrium are marked as passive, and then (4.4) is newly applied (in this way, passive nodes are forced to meet (4.1c)). This iterative process is performed until all nodes satisfy their lower bound constraints. Notice that the last iteration of this procedure corresponds to solving a subproblem of the form (4.3) where the set of passive nodes is equal to the set I. Therefore, this last iteration is equivalent to the second stage of the proposed method.

Summarizing, our method relies on an iterative process that uses the continuoustime protocol (4.4) as a subroutine. The main idea of this methodology is to identify in each step the nodes that have an allocated resource out of their lower bounds. These nodes are marked as passive, so they are forced to satisfy their constraints in subsequent iterations, while active nodes seek to equalize their marginal costs using the remaining resource. In the worst case scenario, the classification between active and passive nodes requires |V| iterations, where |V| is the number of nodes in the network. This fact arises when only one active node becomes passive at each iteration.

Center-Free Resource Allocation Algorithm

The proposed method is formally described in Algorithm 1. Notice that this algorithm is fully decentralized since steps 4-6 can be computed by each agent using only local information. Step 4 corresponds to solve Equation (4.4), while Steps 5 and 6 describe the conditions for converting an active node into passive. Let us note that steps 4-6 have to be performed |V| times since we are considering the worst case scenario.

Therefore, each agent needs to know the total number of nodes in the network. This requirement can be computed in a distributed way by using the method proposed in ( [START_REF] Garin | A survey on distributed estimation and control applications using linear consensus algorithms[END_REF], p.90). We also notice the fact that the agents have to be synchronized (as usual in several distributed algorithms [START_REF] Cortés | Distributed algorithms for reaching consensus on general functions[END_REF][START_REF] Xiao | Optimal scaling of a gradient method for distributed resource allocation[END_REF][START_REF] Garin | A survey on distributed estimation and control applications using linear consensus algorithms[END_REF]) in order to apply the step 4 of Algorithm 1, i.e, all agents must start solving Equation (4.4) at the same time.

Algorithm 1: Resource allocation with lower bounds

Input:

-Parameters of the problem in Equation (4.1).

-An initial value x (0) , such that n i=1 x

i =X. Output: Optimal allocation x * 1 Mark all nodes as active, i.e., Ṽa,0 ← V, Ṽp,0 ← ∅.

2 xi,0 ← x (0) i , for all i ∈ V. 3 for l ← 1 to |V| do 4 xi,l ← x i (t l )
, for all i ∈ V, where x i (t l ) is the solution of Equation (4.4a) at time t l , with initial conditions x(0) = [x 1,l-1 , . . . , xn,l-1 ] , V a = Ṽa,l-1 , V p = Ṽp,l-1 , and {x i (0) = 0, ∀i ∈ V p }.

5

Ṽp,l ← Ṽp,l-1 {i ∈ Ṽa,l-1 : xi < x i }, and Ṽa,l ← Ṽa,l-1 \{i ∈ Ṽa,l-1 : xi < x i }. 7 return x * According to the reasoning described at the beginning of this subsection, we ideally require to know the steady state solution of Equation (4.4) at each iteration of Algorithm 1 (since we need to identify which nodes are allocated with an amount of resource below their lower bounds in steady state). This implies that the time t l in step 4 of Algorithm 1 goes to infinity. Under this requirement, each iteration would demand infinite time and the algorithm would not be implementable. Hence, to relax the infinite time condition, we state the following assumption on the time t l . Assumption 4.2.2. Let x * i,l be the steady state of x i (t) under Equation (4.4), with initial conditions x(0) = xi,l-1 , V a = Ṽa,l-1 , V p = Ṽp,l-1 , and {x i (0) = 0, ∀i ∈ V p }1 .

For each l = 1, . . . , |V| -1, the time t l satisfies the following condition: x i (t l ) < x i if and only if x * i,l < x i , for all i ∈ V.

According to assumption 4.2.2, for the first |V| -1 iterations, we only need a solution of (4.4) that is close enough to the steady state solution. We point out the fact that, if the conditions of Proposition 4.2.2 are met in the lth iteration of Algorithm 1, then

x i (t) asymptotically converges to x * i,l , for all i ∈ V, under Equation (4.4). Therefore, Assumption 4.2.2 is satisfied for large values of t 1 , . . . , t |V|-1 .

Taking into account all the previous considerations, the next theorem states our main result regarding the optimality of the output of Algorithm 1.

Theorem 4.2.1. Assume that G is a connected graph. Moreover, assume that φ i is a strictly convex function for all i = 1, . . . , n. If t 1 , . . . , t |V|-1 satisfy Assumption 4.2.2, and the problem stated in Equation (4.1) is feasible, then the output of Algorithm 4.1 tends to the optimal solution of the problem given in Equation (4.1) as t |V| → ∞.

Proof. The ith component of the output of Algorithm 1 is equal to xi,|V| = x i (t |V| ), where x i (t |V| ) is the solution of Equation (4.4a) at time t |V| , with initial conditions [x 1,|V|-1 , . . . , xn,|V|-1 ] , V a = Ṽa,|V| , and V p = Ṽp,|V| . Hence, it is sufficient to prove that {x * 1,|V| , . . . , x * n,|V| } solves the problem in Equation (4.1). In order to do that, let us consider the following premises (the proof of each premise is written in square brackets).

• P1: {x 1,l , . . . , xn,l } satisfies (4.1b), for all l = 1, . . . , |V| [this follows from Lemma 4.2.1, and form the fact that n i=1 xi,0 = X].

• P2: x * i,l = x i , for all i ∈ Ṽp,l-1 , and for all l = 1, . . . , |V| [this follows directly from Proposition 4.2.2].

• P3: Ṽp,l = Ṽp,l-1 {i ∈ Ṽa,l-1 : x * i,l < x i }, and Ṽa,l = Ṽa,l-1 \{i ∈ Ṽa,l-1 : x * i,l < x i }, for all l = 1, . . . , |V| [this follows from step 5 of Algorithm 1, and from Assumption 4.2.2].

• P4: If for some l, Ṽp,l = Ṽp,l-1 , then Ṽp,l+j = Ṽp,l-1 , for all j = 0, . . . , |V|l [this can be seen from the fact that if the set of passive nodes does not change from one iteration to the next, the steady state of Equation (4.4a) is the same for both iterations]. since P1 and P9 hold, x * i,l has the characteristics given in Proposition 4.2.1, for all i ∈ V, and for all l = 1, . . . , |V|. Hence, φ i (x * i,l ) has the same value for all i ∈ Ṽa,l-1 , and φ i (x * i,l+1 ) has the same value for all i ∈ Ṽa,l . Moreover, since Ṽa,l ⊂ Ṽa,l-1 (according to P3), we have that

φ i (x * i,l ) < φ i (x * i,l+1
), for all i ∈ Ṽa,l . Thus, x * i,l < x * i,l+1 , for all i ∈ V a , l, because φ i is strictly increasing (this follows from the fact that φ i is strictly convex by assumption). Therefore, i∈ Ṽa,l x * i,l < i∈ Ṽa,l x * i,l+1

: k ∈ Ṽa,l , for all l = 1, . . . , |V|. Let φ k (x * k,|V| ) = λ, where λ ∈ R. Moreover, let us define V 0 = {j ∈ V : x * i,|V| > x i }, and V 1 = {j ∈ V : x * i,|V| ≤ x i }. If i ∈ V 0 , then i ∈ Ṽa,|V|-1 (given the fact that, if i / ∈ Ṽa,|V-1| ⇒ i ∈ Ṽp,|V-1| ⇒ x * i,|V| = x i ⇒ i / ∈ V 0 ). Hence, φ i (x * i,|V| ) = φ k (x * k,|V| ) = λ (
this follows from the fact that φ j (x * j,l ) has the same value for all j ∈ Ṽa,l-1 , which in turn follows directly from step 4 of Algorithm 1, and Proposition 4.2.2).

If

i ∈ V 1 , then either i ∈ Ṽa,|V|-1 or i ∈ Ṽp,|V|-1 . In the first case, φ i (x * i,|V| ) = φ k (x * k,|V| ) = λ (following the reasoning used when i ∈ V 0 ). In the second case, ∃l : i ∈ ( Ṽp,l \ Ṽp,l-1 ); hence, φ i (x * i,l ) = φ k (x * k,l
) (this follows from the fact that, if i ∈ ( Ṽp,l \ Ṽp,l-1 ), then i ∈ Ṽa,l-1 ). Furthermore, since i ∈ ( Ṽp,l \ Ṽp,l-1 ), x * i,l < x i (see P3), and given the fact that φ i is strictly increasing, we have that φ i (x * i,l ) < φ i (x i ). Moreover, according to P10,

φ k (x * k,l ) ≥ φ k (x * k,|V| ). Hence, φ i (x i ) > φ k (x * k,|V| ) = λ. In conclusion, if i ∈ V 1 , then φ i (x * i,|V| ) ≥ λ. Thus, we can choose µ i ≥ 0, for all i ∈ V, such that φ i (x * i,|V| ) -µ i = λ, where µ i = 0 if i ∈ V 0 . Hence, let us note that, ∂φ ∂x i x i =x * i,|V| -µ i -λ = 0, for all i ∈ V, where ∂φ ∂x i x i =x * i,|V| = φ i (x * i,|V|
). Therefore, {x * 1,|V| , . . . , x * n,|V| , µ 1 , . . . , µ n , -λ} satisfies the KKT conditions for the problem given in Equation (4.1). Furthermore, since φ(x) is a strictly convex function by assumption, then {x * 1,|V| , . . . , x * n,|V| } is the optimal solution to that problem.

Early Stopping Criterion

Notice that, if the set of passive nodes does not change in the kth iteration of Algorithm 1 because all active nodes satisfy the lower bound constraints (see step 5), then the steady state solutions x * i,k and x * i,k+1 are the same, for all i ∈ V, which implies that the

OPTIMAL RESOURCE ALLOCATION WITH LOWER-BOUND CONSTRAINTS

where, φ b (x) is the new cost function, and > 0 is a constant that minimizes the effect of the barrier function when the solution is far from the boundary of the feasible set.

With this modification, the distributed algorithm is described by the following equation, ẋi =

j∈N i φ b j (x j ) -φ b i (x i ) , for all i ∈ V, (4.8)
where,

φ b i (x i ) = dφ i dx i -db i dx i , i.e., φ b i (x i
) is equal to the marginal cost plus a penalty term induced by the derivative of the corresponding barrier function.

Distributed Replicator Dynamics (DRD)

This methodology is based on the classical replicator dynamics from evolutionary game theory. In the DRD, the growth rate of a population that plays a certain strategy only depends on its own fitness function and on the fitness of its neighbors. Mathematically, the DRD is given by ẋi =

j∈N i (x i -x i )(x j -x j )(v i (x i ) -v j (x j )), v i = -φ i (x i ), for all i ∈ V, (4.9)
where v i is the fitness perceived by the individuals that play the ith strategy. In this case, the strategies correspond to the nodes of the network, and the fitness functions to the negative marginal costs (the minus appears because replicator dynamics are used to maximize utilities instead of minimize costs). On the other hand, it can be shown that, if the initial condition x(0) is feasible for the problem given in Equation (4.1), then x(t) remains feasible for all t ≥ 0, under the DRD.

Distributed Interior Point Method with Exact Barrier Functions (DIPe)

This technique follows the same reasoning of the DIP algorithm. The difference is that DIPe uses exact barrier functions [START_REF] Bertsekas | Necessary and sufficient conditions for a penalty method to be exact[END_REF] to guarantee satisfaction of the lower bound constraints. The exact barrier function for the ith node is given by:

b e i (x i ) = 1 ε [x i -x i ] + , where [•] + = max(•, 0), 0 < ε < 1 2 max x∈F ∇φ(x) ∞
, and F = {x ∈ R n : n i=1 x i = 1, x i ≥ x i } is the feasible region of x for the problem (4.1). Using these exact barrier functions, the augmented cost function can be expressed as,

φ e b (x) = φ(x) + n i=1 b e i (x i ).
The DIPe algorithm is given in terms of the augmented cost function and its generalized

gradient ∂φ e b (x) = [∂ 1 φ e b (x), . . . , ∂ n φ e b (x)] as follows: ẋi ∈ j∈N i (∂ j φ e b (x) -∂ i φ e b (x)), for all i ∈ V, (4.10)
where

∂ i φ e b (x) =    {φ i (x i ) -1 ε } if x i < x i φ i (x i ) -1 ε , φ i (x i ) if x i = x i {φ i (x i )} if x i > x i
In [START_REF] Cherukuri | Distributed generator coordination for initialization and anytime optimization in economic dispatch[END_REF], the authors show that the differential inclusion (4.10) converges to the optimal solution of the problem (4.1) provided that x(0) is feasible.

Comparison

In order to compare the performance of our algorithm with the three methods described above, we use the following simulation scenario: a set of n nodes connected as in Figure 4.1 (we use this topology to verify the behavior of the different algorithms in the face of few communication channels since previous studies have shown that algorithms' performance decreases with the number of available communication links); a nonlinear cost function φ(x) = n i=1 e a i (x i -b i ) + e -a i (x i -b i ) , where a i and b i are random numbers that belong to the intervals (1, 2) and -1 2 , 1 2 , respectively; a resource constraint X = 1; and a set of lower bounds {x i = 0 : i ∈ V}. For each n, we generate 50 problems with We notice that the algorithm proposed in this chapter always reaches the maximum reduction, regardless of the number of nodes that comprise the network. The same happens with the DIPe algorithm. This is an important advantage of our method compared to other techniques. In contrast, the algorithm based on the DRD performs far from the optimal solution. This unsatisfactory behavior is due to the small number of links of the considered communication network. In [START_REF] Pantoja | Distributed optimization using population dynamics with a local replicator equation[END_REF], the authors prove the optimality of the DRD in problems involving well connected networks; however, they also argue that this technique can converge to suboptimal solutions in other cases. On the other hand, the DIP method provides solutions close to the optimum. Nonetheless, its performance decreases when the number of nodes increases. This tendency is due to the influence of barrier functions on the original problem. Notice that, the larger the number of nodes, the bigger the effect of the barrier functions in Equation (4.8).

Regarding the computation time, although convergence of the proposed method is

Application to Optimal Chiller Loading

In this section, we present a distributed solution to the optimal chiller loading problem in multiple chiller systems [START_REF] Chang | Optimal chilled water temperature calculation of multiple chiller systems using Hopfield neural network for saving energy[END_REF], by using the approach developed in this chapter. This problem arises in decoupled chilled-water plants, which are widely used in large airconditioning systems. The goal is to distribute the cooling load among the chillers that comprise the plant for minimizing the total amount of power used by them. For a better understanding of the problem, below we present a brief description of the system. A decoupled chilled-water plant comprised by n chillers is depicted in Figure 4.2.

The purpose of this plant is to provide a water flow f T at a certain temperature T s to the rest of the air-conditioning system. In order to do that, the plant needs to meet a

OPTIMAL RESOURCE ALLOCATION WITH LOWER-BOUND CONSTRAINTS

cooling load C L that is given by the following expression

C L = mf T (T r -T s ), (4.11)
where m > 0 is the specific heat of the water, and T r is the temperature of the water returning to the chillers. Since there are multiple chillers, the total cooling load C L is split among them, i.e., C L = n i=1 Q i , where Q i is the cooling power provided by the ith chiller, which, in turn, is given by

Q i = mf i (T r -T i ), (4.12)
where f i > 0 and T i are, respectively, the flow rate of chilled water and the water supply temperature of the ith chiller. As it is shown in Figure 4.2, we have that f T = n i=1 f i . In order to meet the corresponding cooling load, the ith chiller consumes a power P i that can be calculated using the following expression [START_REF] Chang | Optimal chilled water temperature calculation of multiple chiller systems using Hopfield neural network for saving energy[END_REF]

P i = k 0,i + k 1,i mf i T r + k 2,i (mf i T r ) 2 + + k 3,i -k 1,i mf i -k 4,i mf i T r -2k 2,i (mf i ) 2 T r T i + k 5,i + k 6,i mf i + k 2,i (mf i ) 2 T 2 i , (4.13)
where k j,i , for j = 0, . . . , 6, are constants related to the ith chiller. If we assume that the flow rate f i of each chiller is constant, then P i is a quadratic function of the temperature T i . The optimal chiller loading problem involves the calculation of the chillers' water supply temperatures that meet the total cooling load given in Equation (4.11), and minimize the total amount of power consumed by the chillers, i.e., n i=1 P i . Moreover, given the fact that each chiller has a maximum cooling capacity, we have to consider the following additional constraints

mf i (T r -T i ) ≤ Q i for i = 1, . . . , n, (4.14)
where Q i is the maximum capacity (rated value) of the ith chiller.

Summarizing, the optimal chiller loading problem can be expressed as follows min T 1 ,...,Tn

J = n i=1 P i (T i ) s.t. n i=1 mf i (T r -T i ) = C L T i ≥ T r -Q i mf i , for all i = 1, . . . , n. (4.15)
Now, let us consider that we want to solve the aforementioned problem in a distributed way by using a multi-agent system, in which each chiller is managed by an 4.4 Application to Optimal Chiller Loading agent that decides the value of the water supply temperature. We assume that the ith agent knows (e.g., by measurements) the temperature of the water returning to the chillers, i.e., T r , and the flow rate of chilled water, i.e., f i . Moreover, agents can share their own information with their neighbors through a communication network with a topology given by the graph G. If each P i (T i ) is a convex function, then the problem can be solved by using the method proposed in Algorithm 1 (we take, in this case,

x i = f i T i ).
The main advantage of this approach is to increase the resilience of the whole system in the face of possible failures, due to the fact that the plant operation does not rely on a single control center but on multiple individual controllers without the need for a centralized coordinator.

Illustrative Example

We simulate a chilled-water plant comprised by 7 chillers1 , the cooling capacity and the water flow rate of each chiller are, respectively, Q i = 1406.8 kW, and f i = 65 kg.s -1 , for i = 1, . . . , 7; the specific heat of the water is m = 4.19 kW.s.kg -1 .C -1 ; the supply temperature of the system is T s = 11 C; and the coefficients k j,i of Equation (4.13) are given in Table 4.2. We operate the system at two different cooling loads, the first one is 90% of the total capacity, i.e., C L = 0.9 n i=1 Q i , and the second one is 60% of the total capacity, i.e., C L = 0.6 n i=1 Q i . The P i -T i curves are shown in Figure 4.3a for both cases, it can be noticed that all functions are convex. In order to apply the distributed solution presented in Algorithm 1, we use an agent per chiller (i.e., the ith agent controls the supply temperature T i of the ith chiller) and the communication network shown in Figure 4.1. In all cases the initial conditions of the chillers' supply temperatures are T i (0) = T s , for i = 1, . . . , 7. The results for the first cooling load, i.e., C L = 8862.8 kW, are depicted in Figure 4.3b, where it is shown that the cooling load is properly allocated among the chillers by adjusting the supply temperatures. More efficient chillers (i.e., chiller 3, chiller 6, and chiller 7 in Figure 4.3a) are more loaded than the less efficient ones (i.e., chiller 2 and chiller 5). This can be noticed from the fact that their supply temperatures, in steady state, reach the minimum value. Furthermore, the energy consumption is minimized and power saving reaches to 2.6%. The results for the second cooling load, i.e., C L = 5908.6 kW, are shown in Figure 4.3c, where it can be noticed a similar performance to that obtained with the first cooling load. However, in this case, it is not necessary that the supply temperatures reach the minimum value to meet the required load. Newly, energy consumption is minimized and power saving reaches to 2.8%. As it is stated in Section 4.2, convergence and optimality of the method is guaranteed under the conditions given in Theorem 4.2.1. In both cases we use the early stopping criteria given in Section 4.2.

Although other techniques have been applied to solve the optimal chiller loading problem, e.g., the ones in [START_REF] Chang | Optimal chilled water temperature calculation of multiple chiller systems using Hopfield neural network for saving energy[END_REF], they require centralized information. In this regard, it is worth noting that the same objective is properly accomplished by using our approach, which is fully distributed.

Discussion

The method developed in this chapter solves the problem of resource allocation with lower bounds given in Equation (4.1). The main advantage of the proposed technique is its distributed nature; indeed, our approach does not need the implementation of a centralized coordinator. This characteristic is appealing, especially in applications where communications are strongly limited. Moreover, fully distributed methodologies increase the autonomy and resilience of the system in the face of possible failures.

In Section 4. allows the nodes to share information has few channels. In these cases, the behavior of our approach is better than that of other techniques found in the literature, such as the DIP method, or the DRD. Moreover, it is worth noting that our technique addresses the constraints as hard. This fact has two important consequences: i) in all cases, the solution satisfies the imposed constraints, and ii) the objective function (and therefore the optimum) is not modified (contrary to the DIP method that includes the constraints in the objective function decreasing the quality of the solution as shown in Section 4.3.4).

Other advantage of the method proposed in this chapter is that it does not require an initial feasible solution of the resource allocation problem (4.1). Similarly to the DIPe technique, our method only requires that the starting point satisfies the resource constraint (4.1b), i.e., we need that n i=1 x i (0) = X. Notice that an initial solution x(0) that satisfies (4.1b) is not hard to obtain in a distributed manner. For instance, if we assume that only the kth node has the information of the available resource X, we can use x k (0) = X, {x i (0) = 0 : i ∈ V, i = k} as our starting point. Thus, an initialization phase is not required. In contrast, other distributed methods, such as DIP and DRD needs an initial feasible solution of the problem (4.1), i.e., a solution that satisfies (4.1b) and (4.1c). Finding this starting point is not a trivial problem for systems involving a large number of variables. Therefore, for these methods, it is necessary to employ distributed constraint satisfaction algorithm (as the one described in Section 3.4.2) as a first step.

On the other hand, we notice that to implement the early stopping criterion presented at the end of Section 4.2, it is required to perform an additional min-consensus step in each iteration. Despite this fact, if the number of nodes is large, this criterion saves computational time, because in most of the cases, all passive nodes are identified during the first iterations of Algorithm 1.

Distributed Population Dynamics 1

In the previous chapter, we briefly discussed the distributed replicator dynamics. We showed that, although this model does not reach the optimal solution of resource allocation problems that are subject to sparse graphs, its implementation is appealing to address applications that require fast computation and tolerate suboptimal solutions.

The distributed replicator dynamics were proposed in [81] based on the classic replicator equation [START_REF] Peter | Evolutionary Stable Strategies and Game Dynamics[END_REF], which is a model that belongs to the field of population dynamics [START_REF] William | Population games and evolutionary dynamics[END_REF]. The purpose of this chapter is to extend the results in [81] for dealing with other models, different from the replicator equation. In general terms, we propose a procedure to obtain distributed population dynamics.

Motivation

Population dynamics [START_REF] William | Population games and evolutionary dynamics[END_REF][START_REF] Weibull | Evolutionary Game Theory[END_REF][START_REF] Hofbauer | Evolutionary games and population dynamics[END_REF] describe the dynamical process that a population experiences when there is a strategic interaction among the agents that comprise the population. The agents involved in this dynamical process evolve to an equilibrium according to a revision protocol, which establishes the individual decision rules that agents apply to choose the best strategies (i.e., those strategies earning higher payoffs). Population dynamics properties (e.g., passivity [START_REF] Fox | Population games, stable games, and passivity[END_REF]) can be exploited to design solutions for a variety of engineering problems. For instance, designing of control and learning systems 5.2 Preliminaries by using alternative revision protocols. Besides, we show that a well-mixed population obeys a structure given by a complete graph, whereas a non-well-mixed population has many different possible structures that are generally given by non-complete graphs. In this sense, the proposed approach is versatile to be implemented in a large variety of problems with different information structures. Moreover, we show that the distributed population dynamics exhibit similar stability and invariance properties as their classic counterpart.

Preliminaries

Graph Representation of Population Interactions

We use graph-theoretical tools to describe the constraints on the interactions among agents according to the population type, i.e., well-mixed or non-well-mixed. In this regard, let G = (V, E, A) be a weighted graph, where V = {1, . . . , n} is the set of nodes representing the strategies in the population; E ⊂ {(i, j) : i, j ∈ V} is the set of links representing possible interaction among agents playing different strategies; and A = [a ij ] is an n × n weighted adjacency matrix whose elements satisfy the following property:

a ij > 0 if (i, j) ∈ E, a ij = 0 otherwise.

Population Games

Consider a population of mass m > 0 comprised of a finite and large number of agents in a strategic interaction. Throughout this chapter, we assume without loss of generality, that the mass of the population is equal to one, i.e., m = 1. As was stated before, the set of available strategies for the agents is given by V = {1, ..., n}. The scalar x i ≥ 0 represents the fraction of the population mass that corresponds to the agents choosing the strategy i ∈ V. The vector x ∈ R n + is the state of the population with dimension n whose entries are nonnegative real numbers. Similarly, x ∈ R n ++ denotes a vector of dimension n with strictly positive entries. The set of possible states of the population, which corresponds to all possible distributions of agents among the strategies, is given by the following simplex:

∆ = x ∈ R n + : i∈V x i = 1 .
(5.1)

Mean Dynamics

Taking into account the previously discussed considerations, this section formally describes the evolutionary process of a non-well-mixed population involved in a strategic game. Suppose that the population is composed by M agents, and each of them receives a revision opportunity that is given by an exponential distribution with rate R. Hence, during a time dt, the revision opportunity received by each agent is given by Rdt. Since we assume that the mass of the population is equal to one, the scalar x i is equal to the portion of agents playing the ith strategy, and M x i is the total amount of agents playing strategy i ∈ V. Consequently, the expected number of revision opportunities received by agents playing the ith strategy is approximately M x i Rdt during dt (notice that x i may vary during dt; however, this variation is negligible if dt is small). Agents playing i ∈ V switch to strategy j ∈ V with a probability that depends on the revision protocol, the probability distribution of receiving a revision opportunity, and the encounter probability between strategies i and j (given by the population structure, which is represented by the graph G), i.e., a ij ρ ij /R. Finally, the expected number of agents switching from

strategy i ∈ V to strategy j ∈ V during time dt is M x i a ij ρ ij dt.
Now, if we consider all possible strategies in the population, the expected number of agents switching to strategy i ∈ V is given by M j∈V x j a ji ρ ji dt, and the expected number of agents playing strategy i ∈ V changing to other strategies is given by

M x i j∈V a ij ρ ij dt.
Therefore, the variation of the proportion of agents playing the ith strategy is deduced by a mass balance as follows,

ẋi = j∈V x j a ji ρ ji -x i j∈V a ij ρ ij .
This equation corresponds to the Distributed Mean Dynamics, or mean dynamics for non-well-mixed populations. Since G is undirected, notice that the distributed mean dynamics can be rewritten as follows, ẋi =

j∈N i x j ρ ji -x i j∈N i ρ ij .
(5.2)

DISTRIBUTED POPULATION DYNAMICS 1

For complete graphs, i.e., for well-mixed-populations, we have that N i = V, obtaining the classic Mean Dynamics [START_REF] William | Population games and evolutionary dynamics[END_REF].

Distributed Population Dynamics

Distributed mean dynamics allow the inference of population dynamics involving nonwell-mixed populations comprised of agents that are programmed with a specific revision protocol. This section shows the deduction of different distributed population dynamics using (5.2). The deduced dynamics are named after the classic population dynamics that are generated with the corresponding revision protocol.

Distributed Replicator Dynamics (DRD)

The distributed replicator dynamics are obtained from the distributed mean dynamics using the pairwise proportional imitation protocol (Table 5.1), as follows, ẋi =

j∈N i x j x i [f i -f j ] + -x i j∈N i x j [f j -f i ] + = j∈N i x j x i (f i -f j ).
Finally, the distributed replicator dynamics are given by,

ẋi = x i   f i j∈N i x j - j∈N i x j f j   .
(5.3)

Distributed Smith Dynamics (DSD)

In this case, we use the pairwise comparison protocol (Table 5.1). Substituting this revision protocol in (5.2), we get ẋi =

j∈N i x j [f i -f j ] + -x i j∈N i [f j -f i ] + .
(5.4)

Notice that (5.4) can be written as, ẋi =

j∈N i 1 2 (1 -φ ij )x i + (1 + φ ij )x j [f i -f j ],
where φ ij = sgn(f if j ).

Invariant Set and Stability Analysis

Distributed Logit Dynamics (DLD)

The deduction of the distributed logit dynamics is based on the logit choice protocol (Table 5.1). However, notice that this protocol requires full information since the sum at the denominator is taken over all the strategies. In order to satisfy the information constraint given by the graph G, we modify the protocol as follows,

ρ ij (F, x) = e η -1 f j , η > 0.
Distributed logit dynamics are obtained by replacing the above protocol in the distributed mean dynamics, i.e., ẋi = j∈N i

x j e η -1 f ix i j∈N i e η -1 f j .

(5.5)

Distributed Projection Dynamics (DPD)

The projection dynamics use the modified pairwise comparison protocol, i.e., ẋi =

j∈N i x j [f i -f j ] + x j -x i j∈N i [f j -f i ] + x i , = j∈N i (f i -f j ) .
Thus, the distributed projection dynamics are given by ẋi = |N i |f i -

j∈N i f j , (5.6)
where |N i | denotes the cardinality of the set N i , i.e., the number of neighbors of the ith node.

Invariant Set and Stability Analysis

Invariant Set Analysis

As was mentioned in Section 5.2, the population mass does not vary over time. Hence, all possible states generated during the evolution of the population should belong to the simplex ∆ given in (5.1). This section shows that the simplex ∆ is an invariant set under the distributed population dynamics deduced in the previous section.

DISTRIBUTED POPULATION DYNAMICS 1

Theorem 5.4.1. The simplex ∆ is an invariant set under: the distributed replicator dynamics (5.3), the distributed Smith dynamics (5.4), and the distributed logit dynamics (5.5).

Proof. According to Equation (5.1), ∆ has to conditions: i) i∈V x i = 1 (mass conservation); ii) x i ≥ 0, for all i ∈ V (non-negativeness).

First, we prove that DRD, DSD, and DLD satisfy condition i). Notice that this is equivalent to show that i∈V ẋi = 0 under the distributed mean dynamics (5.2). These dynamics can be written by using the adjacency matrix of the graph G as follows ẋi

= j∈V a ij ρ ji x j - j∈V a ij ρ ij x i . Hence, i∈V ẋi = i∈V j∈V a ij ρ ji x j - i∈V j∈V a ij ρ ij x i . Since G is undirected (i.e., a ij = a ji ), we have i∈V ẋi = i∈V j∈V a ji ρ ji x j - j∈V i∈V a ji ρ ji x j = 0.
Second, we prove that each dynamic satisfies condition ii):

• DRD: Non-negativeness of each x i is satisfied given the fact that ẋi = 0 if x i = 0 under distributed replicator dynamics. Thus, if x i (0) ≥ 0, then x i (t) ≥ 0 for all t ≥ 0.

• DSD: According to (5.4), notice that when x i = 0 for any i ∈ V, then ẋi ≥ 0.

Hence, the non-negativeness of x i is satisfied under distributed Smith dynamics.

• DLD: Notice that ẋi ≥ 0 when x i = 0 under distributed logit dynamics (5.5).

Therefore, if x(0) ∈ ∆, then x i (t) ≥ 0 for all t ≥ 0. Proposition 5.4.1. The set ∆ = {x ∈ R n :

i∈V x i = 1} is invariant under the distributed projection dynamics (5.6).

Invariant Set and Stability Analysis

Proof. The distributed projection dynamics can be written by using the adjacency matrix of the graph G as follows ẋi

= j∈V a ij f i - j∈V a ij f j . Therefore, i∈V ẋi = i∈V j∈V a ij f i - i∈V j∈V a ij f j . Since a ij = a ji because G is undirected, we obtain i∈V ẋi = j∈V i∈V a ji f i - i∈V j∈V a ij f j =0.
Remark 5.4.1. It should be noticed that the distributed projection dynamics satisfy one of the conditions of the original simplex ∆, i.e., i∈V x i = 1 (mass conservation). However, the non-negativeness of x i is not guaranteed. This fact also occurs in the classic projection dynamics. ♦ Remark 5.4.2. Notice that Theorem 5.4.1 and Proposition 5.4.1 do not impose any conditions on the interaction graph G. Thus, the studied distributed population dynamics exhibit simplex invariance under any population structure. ♦

Stability Analysis

Classic population dynamics usually converge to Nash equilibria since Nash equilibria correspond to the expected outcome of games played by rational individuals (i.e., individuals that are trying to maximize their profit). Given a population game F , the set of Nash equilibria is defined as follows [START_REF] William | Population games and evolutionary dynamics[END_REF],

N E(F ) = {x * ∈ ∆ : x * i > 0 ⇒ f i (x * ) ≥ f j (x * ), ∀i, j ∈ V}.
Thus, in a Nash equilibrium, all players perceives the same profit.

This section provides sufficient conditions guaranteeing that a Nash equilibrium x * of the population game F is asymptotically stable under the distributed population dynamics derived in Section 5.3.2. These conditions, which are related to the connectivity of the interaction graph and the characteristics of the Nash equilibrium, are summarized in the following assumptions.

Invariant Set and Stability Analysis

Thus, if x(0) ∈ B, the null space of L (x) is equal to span{1} (we use 1 to denote a vector of dimension n whose entries are all 1) since G (x) is connected (we conclude that G (x) is connected since: i) G (x) and G have the same topology in B, i.e., if x ∈ B, a

ij = 0 only if a ij = 0; ii) G is connected by assumption). In this case, ĖV (x) = 0 if and only if f i = f j , for all i, j ∈ V , i.e., ĖV (x) = 0 only in x * . Therefore, x * is asymptotically stable.

Remark 5.4.3. Proposition 5.4.2 requires that, in steady state, all strategies are played by the individuals involved in the game. Indeed, when any proportion of individuals is extinct at equilibrium (i.e., x * i = 0 for some i ∈ V), then convergence of the distributed replicator equation and the distributed Smith dynamics to the Nash equilibrium is not guaranteed. ♦ Theorem 5.4.3. Let F be a continuously differentiable stable game, let x * ∈ N E(F ), and let ẋ be the distributed projection dynamics (5.6). If Assumptions 5.4.1 and 5.4.2 hold, then x * is asymptotically stable.

Proof. Consider the pairwise comparison protocol ρ ij = [f jf i ] + , and define

ρ ij = ϕ(f j -f i ),
where ϕ(•) = [•] + . Then, consider the Lyapunov function candidate:

V (x) = i∈V j∈V a ij f j-f i 0 ϕ(s)ds.
Since ϕ : R → R + is increasing on [0, +∞) and G is connected, the function V (x) > 0, for all x = x * . Additionally, V (x * ) = 0 since f j (x * ) = f i (x * ), for all i, j ∈ V. Moreover, notice that

∂V (x) ∂x l = i∈V j∈V a ij ∂f j ∂x l - ∂f i ∂x l ϕ(f j -f i ) = i∈V j∈V a ij ϕ(f j -f i ) ∂f j ∂x l - j∈V i∈V a ji ϕ(f i -f j ) ∂f j ∂x l .
Taking into account that a ij = a ji , we obtain

∂V (x) ∂x l = i∈V j∈V a ji (f j -f i ) ∂f j ∂x l = j∈V ∂f j ∂x l i∈N j (f j -f i ).

DISTRIBUTED POPULATION DYNAMICS 1

According to (5.6), notice that i∈N j (f jf i ) = ẋj , where ẋj is the jth element of the distributed projection dynamics ẋ. Hence,

∂V (x) ∂x l = j∈V ẋj ∂f j ∂x l .
(5.7)

Therefore, the time derivative of the Lyapunov function is

V (x) = (∇V (x)) ẋ = ẋ DF (x) ẋ,
where ẋ DF (x) ẋ ≤ 0 since F is stable.

Remark 5.4.4. As was stated in the proof of Theorem 5.4.2, the fact that a Nash equilibrium x * belongs to int∆ implies that all the fitness functions reach the same value, i.e., f i (x * ) = f j (x * ), for all i, j ∈ V. Therefore, the results given in Theorems 5.4.2 and 5.4.3 are related to the contributions reported in the literature on consensus in multiagent networks (e.g., see [START_REF] Olfati-Saber | Consensus and cooperation in networked multi-agent systems[END_REF][START_REF] Moreau | Stability of multiagent systems with time-dependent communication links[END_REF][START_REF] Ren | Consensus seeking in multiagent systems under dynamically changing interaction topologies[END_REF]). An essential difference is that Theorems 5.4.2 and 5.4.3 show a direct relationship between game-theoretic properties and Lyapunov stability of a population game under distributed dynamics. ♦

The connectivity condition of the graph G in Theorems 3 and 4 is sufficient for Nash equilibrium stability. Regarding this fact, it is interesting to study if this condition is also necessary. The following proposition gives us insights on this issue.

Proposition 5.4.2. Assume that the population game F has a unique Nash equilibrium, which is in the interior of the simplex ∆, i.e., x * ∈ int∆. Let ẋ be the distributed mean dynamics (5.2). If x(t) asymptotically converges to x * , for all x(0) ∈ int∆, then the graph G is connected.

Proof. We prove the contrapositive. Assume that G is non-connected. We can express G as the union of r ≥ 2 connected components (maximal connected sub-graphs) denoted by G p = (V p , E p ), where p = 1, . . . , r, i.e., G = r p=1 G p . We use the arguments in the proof of Theorem 5.4.1 to conclude that, under the distributed mean dynamics, i∈V p x i (t) = i∈V p x i (0), for all p = 1, . . . , r, and for all t ≥ 0. Take two connected components G 1 and G 2 of the graph G. Furthermore, take the following initial condition: . Each decision variable is managed by a node in a network. Furthermore, we impose an information constraint given by the graph shown in Figure 5.3. This graph is obtained by following the Erdös-Rényi model (which is the simplest model of several kind of social and biological networks [START_REF] Bornholdt | Handbook of graphs and networks[END_REF]) with edge generation probability equal to 0.01. Besides, we add a path connecting all nodes to guarantee that the generated graph is connected. The information constraint implies that the ith node only has information about the state of its neighbors. To solve the problem in (5.8), we define a full potential game F = ∂V ∂x 1 , . . . , ∂V ∂x 50 (i.e., the fitness functions correspond to the marginal utilities) and apply the distributed population dynamics derived in Section 5.3.2. Notice that all nodes satisfy the information constraint (this fact is not possible by using the classic population dynamics).

x i (0) =        x * i + ε |V 1 | if i ∈ V 1 x * i -ε |V 2 | if i ∈ V 2
Results are shown in Figure 5.4 considering an initial condition x i (0) = The academical example proposed in Equation (5.8) can be adapted to model practical optimization problems. For instance, let us consider the economic dispatch of distributed generators reported in [81]. This problem is stated as follows: consider a set of n generators that are connected to the electric distribution network. These generators have to supply a certain load denoted by L. Therefore, n i=1 p i = L, where p i ≥ 0 is the power supplied by the ith generator. The cost of production of p i is given by a quadratic function c i (p i ) = α i + β i p i + γ i p 2 i , where α i , β i , γ > 0 are parameters associated with the i-th generator. The goal of the economic dispatch problem is to minimize the total cost of energy production, which is given by J(p 1 , . . . , p n ) := n i=1 c i (p i ). Mathematically, this problem is formulated as follows: (5.9)

Notice that the formulation of problems (5.8) and (5.9) are similar. Indeed, if we let the generators be the strategies of a game, and the power supplied by the ith generator be the amount of population playing the strategy i. Then, it is possible to model the economic dispatch of distributed generators as a population game, where the fitness functions are given by the marginal cost of energy production, i.e., f i = ∂J ∂p i , for all i = 1, . . . , n. In this regard, we can address the economic dispatch of distributed generators by using the same procedure described in the solution of the optimization problem (5.8). This approach is appealing, specially in scenarios where the number of distributed generators is large.

Illustrative Examples

Classic Population Games: Bad Rock-Paper-Scissors with a Twin

To evaluate the behavior of a non-well mixed population involved in a strategic interaction, we use a classic game called "bad rock-paper-scissor with a twin" (adapted from [START_REF] William | Population games and evolutionary dynamics[END_REF]). Moreover, we compare the performance of the distributed dynamics proposed in Section 5.3.2 with their classic counterpart (i.e., assuming a well-mixed population).

The bad rock-paper-scissor with a twin preserves the same rules as the rock-paperscissors game. The difference is that the losing strategies have double penalty, i.e., this is not a zero-sum game. Moreover, another strategy called "twin" is added, which earns the same payoff as scissors. Summarizing, this game can be represented by the following payoff matrix

A =     0 -2 1 1 1 0 -2 -2 -2 1 0 0 -2 1 0 0     ,
where the first, second, third, and fourth column/row corresponds to rock, paper, scissors, and twin, respectively. In this case, the fitness functions are given by F (x) = Ax.

Evolution of the population state is shown in Figure 5.5 under distributed (first row) and classic (second row) dynamics. For the distributed case, we use a path graph.

It can be noticed a similar behavior between the trajectories obtained by using the classic population dynamics and its distributed counterpart. DSD, DRD, and DPD exhibit oscillations around the set of Nash equilibria (depicted in red). Indeed, a limit cycle emerges in each case. The main difference is that the limit cycle lies in the plane orthogonal to the set of Nash equilibria for the classic dynamics while it is not orthogonal for the distributed case. Additionally, the angle formed by the plane containing the limit cycle and the set of Nash equilibria is the same for DSD, DRD, and DPD. On the other hand, both DLD and classic logit dynamics reach the same equilibrium point. However, this rest point does not belong to the set of Nash Equilibria. We notice the fact that this characteristic has been reported before in the literature related to classic logit dynamics (e.g., see [START_REF] William | Population games and evolutionary dynamics[END_REF]).

DISTRIBUTED POPULATION DYNAMICS 1

Distributed Control of Dynamical Systems

Finally, we design a distributed controller for the optimal transportation of drinking water2 (relationship between distributed control and games has been pointed out in several papers, e.g., [START_REF] Marden | Cooperative Control and Potential Games[END_REF][START_REF] Wang | A control perspective for centralized and distributed convex optimization[END_REF]). The system is composed by n coupled tanks as shown in Figure 5.6. The arrows in the graphical representation show how flow directions are. Each tank has an outflow given by an unknown demand considered as a disturbance that is denoted by d i , and an inflow u i from a limited water source, i.e., the control outputs are subject to a constraint given by n i=1 u i ≤ K, where K is the total available resource. The mentioned inflows are controlled by x i that determines a percentage of the total resource K, then u i = Kx i , and n i=1 x i = 1. We assume that there are local controllers at the valves to guarantee the desired inflow according to the signal x i . The dynamics for this system are as follows,

dh 1 dt = u 1 -ρgh 1 -d 1 dh i dt = u i + ρgh i-1 -ρgh i -d i , i = 2, ..., n -1 dh n dt = u n + ρgh n-1 -d n , x 1 d n d 2 d 1 Source • • • • • • x 2
x n Figure 5.6: Simple drinking water system with a unique resource and unknown demands.

where h i is the water level of the ith tank, ρ is the density of the fluid, and g is the gravity. The proposed example considers the case of 4 tanks, and the control objective is to maintain the water level of each tank at a safety value of reference, which is given by the company in charge of the management of the network. This safety value is obtained according to the demand that each tank supplies. For this particular example the safety reference is established at 0.5 m, and the unknown demand profile at each node during two days is shown in Figure 5.7 (adapted from [START_REF] Ocampo-Martinez | Hierarchical and decentralised model predictive control of drinking water networks: application to Barcelona case study[END_REF]). To control this plant using population dynamics, we propose the following analogy:

the population is related to the available water, the strategies correspond to the tanks.

Therefore, the population state x represents the distribution of the water (in percentage) among the tanks. Additionally, the fitness functions are selected to be the error at each 5. DISTRIBUTED POPULATION DYNAMICS 1 tank, i.e., f i = 0.5-h i , for all i = 1, . . . , n. Notice that this fitness is appropriate since: i) more proportion of inflow is assigned to those tanks with larger error, and ii) it is known that h i increases as the inflow controlled by x i increases, i.e., the fitness function f i is decreasing with respect to x i satisfying the condition for a stable game [START_REF] William | Population games and evolutionary dynamics[END_REF]. Furthermore, consider the case in which the quantity of water required to meet the safety levels is less than the available resource. In this situation, if all the available resource is allocated in the tanks, the safety levels would be exceeded. Therefore, it is necessary to introduce an additional strategy whose corresponding population x n+1 is used as a slack variable.

The fitness function for this slack variable is f n+1 = 0. Hence, if the level in the ith tank is higher than the safety level (i.e., f i < 0), it is more profitable to allocate resource to the slack variable than to the ith tank.

Figure 5.8 shows the control performance considering full information in the classic Smith dynamics, i.e., that at each point of the network, the information related to all the system is available to make decisions. Figure 5.9 shows the performance of a distributed controller designed based on DSD. The information graph considered for this example is a path graph, i.e., that the ith tank only has information about the (i -1)th and the (i + 1)th tanks. First, notice that in this particular problem with four control actions and one slack variable, the classical Smith dynamics with full information require 10 communication links. In contrast, the distributed population dynamics approach for a problem with four control actions, one slack variable, and a path configuration just Regarding control performance, results show that the control objective is achieved with both the Smith dynamics controller with full information and with the DSD controller by using partial information despite the hard and unknown outflow disturbances at each tank. Furthermore, it can be noticed that the settling time of the distributed controller is not far from the one obtained with full information. Additionally, once the system achieves the set-point, both controllers have similar behavior showing the well performance of the distributed population dynamics approach.

Discussion

This chapter proposes a methodology to generate distributed population dynamics from revision protocols considering different population structures. The core of the proposed method lies in a generalization of the mean dynamics for non-well-mixed populations.

This novel concept provides us a tool to deal with information constraints related to the strategic interactions of the individuals that conform the population. Allowed agents' interactions within the population are described by means of an undirected graph. Following this convention, we have shown that well-mixed populations are characterized by complete graphs, while non-well-mixed populations are represented by non-complete graphs. The main implication of this feature is the possibility to derive population dynamics that only use local information to evolve, i.e., distributed population dynamics.

DISTRIBUTED POPULATION DYNAMICS 1

The appealing features of distributed population dynamics have been exploited in some application examples that illustrate their usefulness. We have shown that distributed population dynamics can be used in the design of constrained optimization algorithms, and the synthesis of control systems for problems requiring non-centralized schemes. Notice that normal form games can also model similar engineering applications. Therefore, distributed learning strategies might be used to solve these problems.

However, the authors in [START_REF] Marden | State based potential games[END_REF] suggest that this approach fails if the application includes a coupled constraint that involves all the decision variables. Under this scenario, distributed population dynamics become a proper alternative since their trajectories evolve to the optimal solution satisfying the referred coupled constraint as shown in Theorems 2-4.

Simulations have illustrated the equivalence between the steady state behavior of the proposed distributed dynamics and classic population dynamics. Moreover, regarding the transient performance, simulations have also shown that there exists a relationship between the convergence rate and the algebraic connectivity of the underlying communication graph. The larger the algebraic connectivity is, the faster the transient response of the corresponding dynamics. Indeed, the fastest response is obtained by using a complete-graph, i.e., by using classic population dynamics. Similarly to other distributed schemes, the dependence of the performance on the graph connectivity evidences a tradeoff between the number of communication links and how fast the optimal solution is reached.

6

Conclusions and Future Directions 6.1 Distributed Resource Allocation Among Dynamic Networks

We have extended the applicability of consensus-based resource allocation algorithms from static problems to dynamic scenarios. We have proved that, if the systems that comprise the network are passive (either dynamical or memoryless), then a multi-agent controller that uses classic consensus protocols can drive the network to a desired global behavior, while a resource constraint is satisfied for any time. This result has been used to synthesize a distributed technique for controlling the temperature in the rooms of a building that uses a central heating system. We have proved by means of graph theoretical tools and invariance analysis, that the proposed controller is able to deal with several power constraints without the need of full information. Moreover, we have used concepts form passivity theory to show that, under connected graphs, our methodology optimally allocates the available heating power among the rooms, considering different scenarios, e.g., when the temperature references cannot be reached due to critical load conditions.

In this situation we have used the variance of the steady state temperature errors as our performance index. Some simulations that show the advantages of the proposed method compared to other widely used techniques have been presented. Furthermore, we have analyzed static resource allocation problems from the perspective of passive memoryless systems. We have shown that convex problems have associated a passive memoryless 6.3 Distributed Population Dynamics some applications of the proposed dynamics in the design of distributed optimization algorithms, solution of classic games, and the synthesis of controllers for multi-variable system. Simulation results have shown that distributed population dynamics exhibit an appropriate performance and provide optimal solutions despite the lack of full information. Thus, distributed population dynamics are promising to address problems requiring non-centralized information structures. This fact is relevant since distributed problems are becoming more common in many fields, especially when large-scale and complex systems are involved.

A future research direction is the inclusion of dynamics in the fitness functions (in a similar way as in Chapter 3). We notice that this modification would give more flexibility to distributed population dynamics for dealing with a wider range of problems.

3. 1 A

 1 large-scale system controlled by a network of agents. . . . . . . . . . . 3.2 Feedback interconnection of systems (3.2) and (3.4). . . . . . . . . . . . . 3.3 n rooms arranged in a row. . . . . . . . . . . . . . . . . . . . . . . . . . . 3.4 Schematic representation of a room and a wall. The leftmost figure illustrates the cross section of a wall divided into layers, while the rightmost one shows the top view of a room surrounded by walls. . . . . . . . . . . . 3.5 a. Ambient temperature during a day (7h-20h). b. Temperatures of the building's rooms. Rooms 1-17 (red), rooms 18-34 (yellow), rooms 35-50 (blue). c. Heating Power allocated in the building's rooms. d. Variance of the error across the rooms. Distributed controller (solid line), PI controllers (dashed line). . . . . . . . . . . . . . . . . . . . . . . . . . . 3.6 Cost function for each aggregator's client. . . . . . . . . . . . . . . . . . . 3.7 Communication network with single-path topology. . . . . . . . . . . . . . 3.8 Evolution of the decentralized resource allocation algorithm. . . . . . . . . 3.9 Evolution of the total cost perceived by the aggregator's clients. . . . . . . 4.1 Single path topology for n nodes. . . . . . . . . . . . . . . . . . . . . . . . 4.2 Decoupled chilled-water plant with n chillers (adapted from [2]). . . . . . 4.3 a. P i -T i curves for each chiller, b. Evolution of supply temperatures and total power consumed by the chillers, C L = 8862.8 kW, c. Evolution of supply temperatures and total power consumed by the chillers, C L = 5908.6 kW. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.1 a. Well-mixed population. b. Non-well-mixed population. . . . . . . . . v

1. 4 . 4

 44 Chapter 5Q3 is answered in Chapter 5, where we extend the results in[81] in order to avoid the need of full information of classic population dynamics. Our main contribution is the design of a general method that allows us to deduce several distributed population dynamics (i.e., population dynamics that evolve using only local information). To illustrate our methodology, we present a distributed version of the fundamental population dynamics (those obtained by applying classic revision protocols), i.e., the distributed replicator dynamics, the distributed Smith dynamics, the distributed logit dynamics, and the distributed projection dynamics. Finally, we prove that, although distributed population dynamics do not use full information, they preserve similar properties as their classic counterpart, i.e., simplex invariance (cf., Theorem 5.4.1 and Proposition 5.4.1), and convergence to Nash equilibria (cf., Theorems 5.4.2 and 5.4.3). Related Publications: • J. Barreiro-Gómez, G. Obando, N. Quijano, "Distributed Population Dynamics: Optimization and Control Applications," Submitted to the IEEE Transactions on Systems, Man, and Cybernetics: Systems, 2015.

2. 2 .

 2 1 to analyze the stability properties of local control laws.

Figure 3 . 1 :

 31 Figure 3.1: A large-scale system controlled by a network of agents.

Figure 3 . 2 :

 32 Figure 3.2: Feedback interconnection of systems (3.2) and (3.4).

( 3 . 2 )

 32 ,(3.4) , with u defined in(3.17), in error coordinates with respect to the equilibrium point (T * , x * ) as follows. ėT = g e (e T , e x ) ėx = -L(G) h e (e T , e x ) + ϑ e (e x ) .(3.18)Let us make the following assumption on this system. Assumption 3.2.3. Consider the error dynamics (3.18). If g e (0, e x ) = 0 and ϑ e (e x ) = 0, then e x = 0. Different from Assumption 3.2.1, Assumption 3.2.3 allows the system (3.2) to has other equilibrium points of the form (T * , x * ), where x * = x * , as long as these points do not satisfy the equation (x

Figure 3 . 3 :

 33 Figure 3.3: n rooms arranged in a row.

Figure 3 . 4 :

 34 Figure 3.4: Schematic representation of a room and a wall. The leftmost figure illustrates the cross section of a wall divided into layers, while the rightmost one shows the top view of a room surrounded by walls.

  Theorem 3.2.1, it is convenient to write the dynamics of the building thermal system given in Equation(3.26) in error coordinates, i.e., expressing(3.26) in the form (3.7). To do this, let T * = [T * 1 , . . . , T * W ] , x * = [x * 1 , . . . , x * N ] be an equilibrium point of the referred system (considering fixed values for T a 1 , . . . , T a W , d 1 , . . . , d n ). Moreover, in order to simplify the notation, let us define the weighted graph G b = {V b , E b , A b } that represents the connections between the elements that comprise the building, where the set of nodes V b = {1, . . . , W } is related to the components of the building, following the order established in Equation (3.25) (e.g., node 1 corresponds to the first room, node n+1 corresponds to the first layer of the wall 1, etc.). The set of edges E b ⊂ V b ×V b is related to the thermal interactions between the building's components, i.e., (i, j) ∈ E b if and only if α i,j > 0. The weighted adjacency matrix A b is formed by the thermal conductances of the junctions between the components of the building, i.e., the component that corresponds to the ith row and jth column of A b is given by α i,j . Using G b , we can express the 3.3 Distributed Building Temperature Control Under Power Constraints building thermal model in error coordinates as follows,

Proposition 3 . 3 . 1 .

 331 Let T * = [T * 1 , . . . , T * W ] , x * = [x * 1 , . . . , x * n ] be an equilibrium point of the system stated in Equation (3.26). If Assumption 3.3.1 holds, then (3.27) is strictly passive from the input e

Proposition 3 . 3 . 2 .

 332 keeps the control signals feasible, the multi-agent controller (3.34) have the form (3.4),(3.17) . Therefore, we can use the result in Theorem 3.2.3 to study the properties of the proposed control strategy. Consider the building thermal system (3.26) controlled via the multi-agent network (3.34) with f = [f 1 , . . . , f n+1 ] given by Equations (3.33) and (3.35).

3. 3

 3 Distributed Building Temperature Control Under Power Constraints Equation (3.34) under different conditions. Moreover, we show that these equilibrium points solve the problem stated in Section 3.3.3.

3. 5a .

 5a The results obtained with the proposed controller and a communication network with ring topology, are shown in Figures 3.5b-3.5d. Figure 3.5b shows the evolution 3.3 Distributed Building Temperature Control Under Power Constraints

Figure 3 .

 3 Figure 3.5: a. Ambient temperature during a day (7h-20h). b. Temperatures of the building's rooms. Rooms 1-17 (red), rooms 18-34 (yellow), rooms 35-50 (blue). c. Heating Power allocated in the building's rooms. d. Variance of the error across the rooms. Distributed controller (solid line), PI controllers (dashed line).

3. 2 . 4 .

 24 Notice that(3.43) is equivalent to the resource allocation algorithm given by Equations (3.4) and(3.6), where y = ∇φ(x). Let us show that the memoryless system 3. DISTRIBUTED RESOURCE ALLOCATION AMONG DYNAMIC NETWORKS y = ∇φ(x) is passive from the input xx * to the output yy * , where x * is the global optimum of (3.42), y * = ∇φ(x * ), and x satisfies the resource constraint, i.e., n i=1

  Figure 3.6. Let us note that the cost functions for clients 1 and 4 (5 and 8) are different despite the fact that their corresponding benefit functions have the same parameters. This happens due to the different power consumption levels. Moreover, in order to address the inequality constraints of the problem, we use logarithmic barrier functions, as it was explained in Section 3.4.1. The evolution of of the distributed resource allocation algorithm (3.43) is depicted in Figure 3.8. In this case, we use the following initial feasible point: x 1 (0) = 0.35 MW, x 2 (0) = 0.53 MW, x 3 (0) = 0.53 MW, x 4 (0) = 0.23 MW,

Figure 3 . 6 : 56 3. 4

 36564 Figure 3.6: Cost function for each aggregator's client.

Figure 3 . 7 :

 37 Figure 3.7: Communication network with single-path topology.

Figure 3 . 8 :Figure 3 . 9 :

 3839 Figure 3.8: Evolution of the decentralized resource allocation algorithm.

Lemma 4 . 2 . 1 .

 421 The quantity

n

  i=1 x i (t) is invariant under Equation (4.4), i.e., if n i=1 x i (0) = X, then n i=1 x i (t) = X, for all t ≥ 0.Proof. It is sufficient to prove that ∆ = 0, where ∆ = n i=1 x i . Notice that ∆ = n i=1 ẋi = 1 ẋ, where ẋ = [ ẋ1 , . . . , ẋn ] . Moreover, according to Equation (4.4),

6x*

  ← [x 1,l , . . . , xn,l ] .

Figure 4 . 1 :

 41 Figure 4.1: Single path topology for n nodes.

  slower than the one shown by DRD and DIP, it is faster than the convergence of the method based on exact barrier functions, i.e., DIPe. Therefore, among the methods that guarantee optimality of the solution, our technique shows the best convergence speed. Computation time taken by DIPe is affected by the use of penalty terms that generate strong changes in the value of the cost function near to the boundaries of the feasible set. The drastic variations of the generalized gradient of exact barrier functions produces oscillations of numerical solvers around the lower bounds (a visual inspection of the results given in Figure3of[START_REF] Cherukuri | Distributed generator coordination for initialization and anytime optimization in economic dispatch[END_REF] confirms this claim). These oscillations are the main responsible for the low convergence speed shown by DIPe. On the other hand, DRD and DIP exhibit the fastest convergence. Hence, DRD and DIP are appealing to be implemented in applications that require fast computation and tolerate suboptimal solutions.

Figure 4 . 2 :

 42 Figure 4.2: Decoupled chilled-water plant with n chillers (adapted from [2]).

Figure 4 .

 4 Figure 4.3: a. P i -T i curves for each chiller, b. Evolution of supply temperatures and total power consumed by the chillers, C L = 8862.8 kW, c. Evolution of supply temperatures and total power consumed by the chillers, C L = 5908.6 kW.

1 50 i=1 x i = 1 x

 1501 Solving a Distributed Optimization Problem First, we propose the following distributed optimization problem, max V (x) := -x x + b x s.t. i > 0 for all i = 1, . . . , 50, (5.8) where x ∈ R 50 is the vector of decision variables, and b ∈ R 50 is a vector of constants, whose entries are given by b i = 2i 1275 , i.e., b = 1 1275 [2 4 ... 100]

Figure 5 . 3 :

 53 Figure 5.3: Non-complete graph for the distributed optimization example.

50 i=1

 50 Figure 5.4, DRD, DSD, and DPD equalize the fitness functions' values in steady state, i.e., these dynamics converge to a Nash equilibrium. This behavior is consistent with the results stated in Theorems 5.4.2 and 5.4.3 since V (x) corresponds to a strictly concave potential function, i.e., F is a full potential and stable game. Moreover, convergence time varies from one dynamic to another. DLD shows the fastest time response while the convergence of DRD is the slowest.

  min J(p 1 , . . . , p n ) := n i=1 α i + β i p i + γ i p 2 i s.t. n i=1 p i = L p i > 0 for all i = 1, . . . , n.

Figure 5 . 4 :

 54 Figure 5.4: Evolution under different distributed population dynamics: states (1st row), fitness functions (2nd row), objective function (3rd row), equality constraint (4th row).

8 )

 8 Figure 5.5: Distributed (up) and classic (bottom) population dynamics applied to bad rock(R)-paper(P)-scissors(S) with a twin(T).

Figure 5 . 7 :

 57 Figure 5.7: Demand profile (d 1 , d 2 , d 3 and d 4 ) during 2 days for the 4 tanks case.

Figure 5 . 8 :

 58 Figure 5.8: System states evolution for a controller with full information.

102 5Figure 5 . 9 :

 10259 Figure 5.9: System states evolution for a distributed population dynamics based controller.

  and (3.35), meets the conditions of Theorem 3.2.3. • First, let us show that (3.26)(3.34) has an equilibrium point (T * , [x * , x *

	n+1 ] )
	that satisfies n+1 i=1 x * i = Q: According to Equation (3.27), if (T * , [x * , x * n+1 ] ), is
	a rest point of (3.26), then
	0

Table 3 . 1 :

 31 Utility functions' parametersk i,1 k i,2 k i,3 k i,4

	i = 1 2.95 2.95 2.00 1.77
	i = 2 1.77 1.77 2.00 0.53
	i = 3 1.31 1.31 2.00 0.13
	i = 4 2.95 2.95 2.00 1.77
	i = 5 1.30 1.30 5.00 0.39
	i = 6 1.11 1.11 5.00 0.11
	i = 7 1.60 1.60 5.00 0.96
	i = 8 1.30 1.30 5.00 0.39
	1	2	3	4
	8	7	6	5

  , which contradicts P8]. Now, let us prove that {x * 1,|V| , . . . , x * n,|V| } solves the Problem in Equation (4.1). First, the solution {x * 1,|V| , . . . , x * n,|V| } is feasible according to P1 and P7. On the other hand, from P9, it is known that ∃k

Table 4 .

 4

		1: Distributed Algorithms' Performance
		Percentage Decrease, Computation Time
	Number of Nodes Proposed Approach	DIP	DRD	DIPe
	n = 5	100%, 0.08 s.	100%, 0.04 s. 91%, 0.04 s.	100%, 7 s.
	n = 20	100%, 0.7 s.	99%, 0.4 s.	69%, 0.2 s. 100%, 153 s.
	n = 50	100%, 3.2 s.	98%, 1.6 s.	57%, 1.3 s. 100%, 841 s.
	n = 100	100%, 17.8 s.	96%, 9.1 s.	50%, 6.5 s.	-
	n = 200	100%, 181.2 s.	94%, 68.7 s. 46%, 41.6 s.	-

Table 4 .

 4

			2: Chillers' parameters
	i = 1 i = 2 i = 3	i = 4	i = 5 i = 6 i = 7
	k 0,i 113.51 71.70 62.75 112.68 74.13 61.98 76.54
	k 1,i 0.21	-0.45	0.49	0.18	-0.44	0.55	0.34
	k 2,i 0.35	0.48	0.30	0.36	0.44	0.30	0.31
	k 3,i -8.19	-5.13 -4.53	-8.13	-5.30 -4.48 -5.527
	k 4,i 0.43	-0.14	0.71	0.40	-0.13	0.76	0.54
	k 5,i 0.14	0.09	0.08	0.14	0.09	0.08	0.10
	k 6,i 0.01	-0.01	0.02	0.01	-0.01	0.02	0.01

  . . . , 50. First and fourth row of Figure5.4 show that x(t) satisfies the problem constraints for all time, i.e., x i (t) remains nonnegative, for all i = 1, . . . , 50; and

	5.5 Illustrative Examples
	i = 1,
	1 50 , for all

As well as in step 4 of Algorithm 1, we have initialized the auxiliary variables xi to zero by convention. If these variables are initialized to other value, convergence of (4.4) is not affected (cf. Proposition 4.2.2).

Algorithms were implemented in a computer with an Intel Core i5 processor.

Simulation parameters are adapted from[START_REF] Chang | Optimal chilled water temperature calculation of multiple chiller systems using Hopfield neural network for saving energy[END_REF].

G. Obando and J. Barreiro-Gómez contributed equally to this chapter. Results are reported in[START_REF] Barreiro-Gómez | Distributed Population Dynamics: Optimization and Control Applications[END_REF].

Parameters of this academical problem have been scaled from real transportation of drinking water problems.

• P5: Ṽa,l Ṽp,l = V, for all l = 1, . . . , |V| [from P3, we know that Ṽa,l Ṽp,l = Ṽa,l-1 Ṽp,l-1 , for all l = 1, . . . , |V|. Moreover, given the fact that Ṽp,0 = ∅, and Ṽa,0 = V, (see step 1 of Algorithm 1) we can conclude P5].

• P6: Since the problem in Equation Hence, from P2 and P3, we know that x * i,l ≤ x i , for all i ∈ V; moreover, {i ∈ Ṽa,l-1 : x * i,l < x i } = ∅. Therefore, n i=1 x * i,l < n i=1 x i . According to P1, we know that n i=1 x * i,l = X; thus, X < n i=1 x i , which contradicts the feasibility assumption].

• P7: {x * 1,|V| , . . . , x * n,|V| } satisfies the constraints (4.1c) [in order to prove P7, we proceed by contradiction: assume that {x * 1,|V| , . . . , x * n,|V| } does not satisfy the constraints (4.1c). Since P2 holds, this assumption implies that {i ∈ Ṽa,|V-1| : x * i,|V| < x i } = ∅. Therefore, Ṽp,|V| = Ṽp,|V|-1 (see P3). Using P4, we can conclude that Ṽp,|V| = Ṽp,|V|-1 = • • • = Ṽp,0 = ∅, i.e., {i ∈ Ṽa,|V|-j : x * i,|V|-j+1 < x i } = ∅, for all j = 1, . . . , |V|. Thus, according to P3,

• P8:

i∈ Ṽa,l x * i,l ≥ i∈ Ṽa,l x * i,l+1 [we prove P8 as follows: using P1 and the result in Lemma 4.2.1, we know that i∈V x * i,l = i∈V x * i,l+1 = X. Moreover, according to P5, V can be expressed as V = Ṽa,l Ṽp,l , where Ṽp,l = Ṽp,l \ Ṽp,l-1 Ṽp,l-1 (see P3). Thus, we have that i∈ Ṽa,l x * i,l + i∈ Ṽp,l ,i / ∈ Ṽp,l-1 x * i,l + i∈ Ṽp,l-1 x * i,l = i∈ Ṽa,l x * i,l+1 + i∈ Ṽp,l ,i / ∈ Ṽp,l-1 x * i,l+1 + i∈ Ṽp,l-1 x * i,l+1 . Furthermore, since Ṽp,l-1 ⊂ Ṽp,l , and given the fact that P2 holds, we have that i∈ Ṽa,l x * i,l

, where x ix * i,l > 0, for all i ∈ Ṽp,l , i / ∈ Ṽp,l-1 (according to P3). Hence, we can conclude P8].

• P9: There exists k, such that k ∈ Ṽa,l , for all l = 1, . . . , |V| [in order to prove P9, we use the fact that, if k ∈ Ṽa,l , then k ∈ Ṽa,l-j , for all j = 1, . . . , l (this follows from P3). Moreover, according to P5 and P6, | Ṽa,|V| | = 0; hence, there exists k, such that k ∈ Ṽa,|V| . Therefore, P9]. P9 guarantees that Assumption 4.2.1 holds at each iteration.

), for all i ∈ Ṽa,l [we prove P10 by contradiction: assume that φ i (x * i,l ) < φ i (x * i,l+1 ), for some i ∈ Ṽa,l . According to Proposition 4.2.2, and

Simulation Results and Comparison

set of passive nodes also does not change in the (k + 1)th iteration. Following the same reasoning, we can conclude that

, for all i ∈ V. Therefore, in this case, {x * 1,k , . . . , x * n,k } is the solution of our resource allocation problem. Practically speaking, this implies that Algorithm 1 does not need to perform more iterations after the kth one. Thus, it is possible to implement a flag z * i (in a distributed way) that alerts the agents if all active nodes satisfy the lower bound constraints after step 4 of Algorithm 1. A way to do that is by applying a min-consensus protocol [START_REF] Cortés | Distributed algorithms for reaching consensus on general functions[END_REF] with initial conditions z i (0) = 0 if the node i is active and does not satisfy its lower bound constraint, and z i (0) = 1 otherwise. Hence, notice that our flag z * i (i.e., the result of the min-consensus protocol) is equal to one, for all i ∈ V, only if all the active nodes satisfy the lower bound constraints, which corresponds to the early stopping criterion described above.

Simulation Results and Comparison

In this section, we compare the performance of our algorithm with other continuoustime distributed techniques found in the literature. We have selected three techniques that are capable to address nonlinear problems and can handle lower bound constraints: i) a distributed interior point method [START_REF] Xiao | Optimal scaling of a gradient method for distributed resource allocation[END_REF], ii) the local replicator equation [START_REF] Pantoja | Distributed optimization using population dynamics with a local replicator equation[END_REF], iii) a distributed interior point method with exact barrier functions [START_REF] Cherukuri | Distributed generator coordination for initialization and anytime optimization in economic dispatch[END_REF]. The first one is a traditional methodology that uses barrier functions; the second one is a novel technique based on population dynamics; and the third one is a recently proposed method that follows the same ideas as the first one, but replaces classic logarithmic barrier functions by exact penalty functions. Below, we briefly describe the aforementioned algorithms.

Distributed Interior Point (DIP) Method

This algorithm is a variation of the one presented in Equation (4.2) that includes strictly convex barrier functions to prevent the solution to flow outside the feasible region (cf., Section 3.4.1). The barrier functions b i (x i ) are added to the original cost function as follows,

for all i ∈ V, 5. DISTRIBUTED POPULATION DYNAMICS 1 [START_REF] Marden | Game Theory and Distributed Control[END_REF]110], optimization [111,[START_REF] Barreiro-Gomez | Constrained distributed optimization based on population dynamics[END_REF], coordination [START_REF] Tembine | Evolutionary games in wireless networks[END_REF], dynamic resource allocation [39],

and so forth.

When using population dynamics for solving learning, control, and optimization problems, some elements of the problem are associated with "strategies" that agents in the population can adopt, and other elements are associated with "masses" of agents playing each strategy. This analogy has a direct implication in the information required to implement a solution based on a population dynamics algorithm, since the existing algorithms assume that the population is well-mixed [START_REF] William | Population games and evolutionary dynamics[END_REF][START_REF] Weibull | Evolutionary Game Theory[END_REF] (i.e., any pair of agents playing any pair of strategies can interact with each other). A consequence of the wellmixed population structure assumption is that the elements of the problem are allowed to interact each other without any constraint (i.e., following a full-information structure).

Therefore, classic population dynamics are restricted to be implemented in problems characterized by a centralized information scheme. However, the number of problems that require distributed solutions has increased dramatically in the last few years. In this regard, some approaches have been proposed to model the interaction constraints in a population of players. For instance, the authors in [START_REF] Berninghaus | Local Interaction on Random Graphs[END_REF]114] deal with normalform games and propose a graph-based interaction model, where each node in the graph represents an individual that repetitively plays a symmetric game with its neighbors.

However, this approach is not suitable to be applied in the population game framework since, in this framework, it is preferable to avoid the individuation of players [START_REF] William | Population games and evolutionary dynamics[END_REF]. On the other hand, other approaches aim to apply learning algorithms that are capable to deal with information constraints [START_REF] Gharesifard | Distributed convergence to Nash equilibria by adversarial networks with directed topologies[END_REF][START_REF] Boussaton | On the distributed learning of Nash equilibria with minimal information[END_REF]. Similarly, the authors in [81] modify the well known replicator dynamics model to relax the full-information dependency. They propose a distributed replicator equation in which the evolution of each strategy is only governed by the "neighboring" strategies (according to a given topology).

This chapter extends the results in [81]. Our main contribution is the design of a general method that allows us to deduce several distributed population dynamics. The core of the proposed method is the use of the mean dynamics [START_REF] William | Population games and evolutionary dynamics[END_REF] in non-well-mixed populations. To illustrate our methodology, we present a distributed version of the fundamental population dynamics (those obtained by applying classic revision protocols), i.e., the distributed replicator dynamics, the distributed Smith dynamics, the distributed logit dynamics, and the distributed projection dynamics. It is worth noting that the deduction presented in this chapter can be used to generate other distributed dynamics Agents playing the ith strategy obtain a reward given by a fitness function f i (x), i.e., f i : ∆ → R is a continuous map that specifies the payoff associated with the strategy i ∈ V. Notice that a population game is completely characterized by the fitness vector

There are several types of games depending on F (x).

Below, we present the definition of two of them, which are found in a large number of applications [START_REF] Hofbauer | Stable games and their dynamics[END_REF].

Definition 5.2.1. Let F : R n + → R n be a population game with payoffs defined on the positive orthant. If there exists a continuously differentiable potential function

Definition 5.2.2. The population game F : ∆ → R n is a stable game if:

The next theorem, adapted from [START_REF] William | Population games and evolutionary dynamics[END_REF], gives an equivalent characterization of a stable game in terms of the Jacobian matrix of F (x), i.e., DF (x); and the tangent space of the simplex ∆.

Theorem 5.2.1. Let the population game F : ∆ → R n be continuously differentiable.

F is a stable game if and only if ξ DF (x)ξ ≤ 0, for all ξ ∈ T ∆, x ∈ ∆, where T ∆ is the tangent space of the simplex ∆, which is defined by

A population game combined with a revision protocol lead to the emergence of population dynamics [START_REF] William | Population games and evolutionary dynamics[END_REF]. The function ρ : R n ×∆ → R n×n + is known as the revision protocol, and describes the timing and the result of the decisions of agents in the strategic interaction. The revision protocol takes the payoff vector F (x) and a determined population state x ∈ ∆, returning a non-negative matrix, whose element of the ith row and jth column ρ ij (F (x), x) represents the conditional switch rate from strategy i to strategy j, where i, j ∈ V. Depending on the revision protocol used by the individuals, we can find several kinds of population dynamics (see Table 5.1), e.g., replicator dynamics, Smith dynamics, logit dynamics, projection dynamics, etc.

In order to simplify the notation for fitness functions, and revision protocols, we remove their corresponding arguments from now on, i.e, f i = f i (x), and ρ ij = ρ ij (F (x), x).

General dynamics on graphs

Table 5.1: Some revision protocols and their corresponding population dynamics [START_REF] William | Population games and evolutionary dynamics[END_REF] Revision protocol Population dynamics

Pairwise proportional imitation:

Logit choice protocol:

Projection dynamics

General dynamics on graphs

The dynamics describing a population behavior depend on the population structure. In this regard, current literature assumes that the population under consideration is wellmixed, i.e., if we take any portion of the entire population, this contains all the strategies with the same probability. Figure 5.1a illustrates this fact by showing a population composed by a large and finite number of agents involved in a game. Each element in the figure represents an agent, and the shape of the element ("circle," "square," or "triangle") denotes the strategy that the agent has adopted. In population games, all agents have the same probability to receive a revision opportunity. The agent receiving the revision opportunity randomly choose another agent from its neighbors and can change its own strategy by the neighbor's strategy depending on the selected revision protocol. Since the population is well-mixed, the probability that the selected opponent is playing any of the available strategies is the same.

On the other hand, there could be a non-well-mixed population as the one shown in Figure 5.1b. For this population, all agents have the same probability to receive an opportunity to make a revision. However, the probability that the opponent is playing a particular strategy is not equal (e.g., if the strategy played by the agent receiving the revision opportunity is "square," then there is the same probability to select an opponent playing strategy "triangle" or "square;" but the probability to select an opponent playing strategy "circle" is zero since the population structure does not allow it).

DISTRIBUTED POPULATION DYNAMICS 1

a. b. Interactions among agents playing different strategies can be represented by a graph G = {V, E, A}. The set of nodes V is associated with the available strategies and the set of links E is related to the encounter probability between strategies, i.e., there exists a link between two strategies if their encounter probability is different from zero. Hence, the elements of the corresponding adjacency matrix A = [a ij ] are as follows:

denotes that strategies i and j can encounter each other, while a ij = 0 denotes that the population structure makes impossible a matching between strategies i and j. According to this convention, the scenarios associated with well-mixed and non-well-mixed populations can be represented by two kinds of graphs. The well-mixed population case is always represented by a complete graph, whereas a non-well-mixed population is represented by a graph with a specific topology depending on the particular population structure (see Figure 5.2).

In this chapter, we assume that the encounter probability between strategies i and j is the same as the one of strategies j and i, i.e., the graph G is undirected. The non-well-mixed population in Figure 5.1b. Nodes 1,2, and 3 correspond, respectively, to strategies "circle," "triangle," and "square" of the proposed example.

DISTRIBUTED POPULATION DYNAMICS 1

Assumption 5.4.1. The graph G that describes the population structure is connected.

Assumption 5.4.2. The Nash equilibrium x * belongs to the interior of the simplex ∆, i.e., x * ∈ int∆, where int∆ = x ∈ R n ++ : i∈V x i = 1 .

Below, we provide our results on convergence of the distributed population dynamics to a Nash equilibrium.

Theorem 5.4.2. Let F be a full potential game with strictly concave potential function V (x), and let x * ∈ N E(F ). If Assumptions 5.4.1 and 5.4.2 hold, then x * is asymptotically stable under the distributed replicator dynamics (5.3) and the distributed Smith dynamics (5.4).

Proof. Since x * ∈ N E(F ) and x * ∈ int∆, we conclude that f i (x * ) = f j (x * ), for all i, j ∈ V. Moreover, notice that x * = arg max x∈∆ V (x) (applying the Karush-Kuhn-Tucker conditions). Additionally, since V (x) is strictly concave, we can take E V (x) = V (x * ) -V (x) as a Lyapunov function candidate. The derivative of E V (x) along the trajectories of DRD (5.3) and DSD (5.4) is given by

is a matrix whose entries l (x) ij are defined as follows, for DRD: l

Notice that L (x) is the Laplacian of the undirected graph given by the tuple

ij ] is the adjacency matrix whose entries are defined as follows: a

These entries are nonnegative since x ∈ ∆. Thus, L (x) ≥ 0 and ĖV (x) ≤ 0. Therefore, x * is stable under DRD and DSD.

Considering that x * ∈ int∆ is stable, a set B around x * can be defined such that if x(0) ∈ B, then x(t) ∈ int∆, for all t ≥ 0 (it is possible to show that B = int∆ for DRD).

where i ∈ V, and ε > 0. Notice that, for small values of ε, x(0) ∈ int∆ since x * ∈ int∆. Under this initial condition, it is not possible that x(t) converges to the unique Nash equilibrium x * since i∈V

Therefore, if we need to guarantee convergence to the Nash equilibrium from any initial condition inside the simplex ∆, we require that the graph G is connected. However, this condition might not be necessary if the initial conditions x(0) are constrained. For instance, suppose that the graph G in the population game is non-connected. Moreover, G is composed of r connected components (maximal connected sub-graphs) denoted by

, where p = 1, . . . , r, i.e., G = r p=1 G p . Then, it can be shown (following the same reasoning as in proof of Proposition 5.4.2) that the equilibrium point x * ∈ int∆ is asymptotically stable if i∈V p x * i = i∈V p x i (0) for all p = 1, . . . , r. Consequently, in this case, the connectivity condition of G is not necessary.

Illustrative Examples

Classic population dynamics have been applied to solve a large number of engineering problems. For instance, access control in communication networks [START_REF] Tembine | Evolutionary games in wireless networks[END_REF], combinatorial optimization [START_REF] Bomze | Approximating the maximum weight clique using replicator dynamics[END_REF], bandwidth allocation [START_REF] Poveda | Dynamic bandwidth allocation in wireless networks using a Shahshahani gradient based extremum seeking control[END_REF], hierarchical frequency control in microgrids [START_REF] Mojica-Nava | Dynamic population games for optimal dispatch on hierarchical microgrid control[END_REF], dispatch of electric generators [START_REF] Pantoja | A Population Dynamics Approach for the Dispatch of Distributed Generators[END_REF], building temperature control [39], constrained extremum seeking [START_REF] Poveda | Shahshahani gradient-like extremum seeking[END_REF], control of drinking water networks [START_REF] Barreiro-Gomez | Distributed control of drinking water networks using population dynamics: Barcelona case study[END_REF], and so forth.

The same applications can be addressed by using distributed population dynamics. This approach has some benefits related to the information privacy, resilience to central failures, and parallelization of the computations. In order to show the versatility of the distributed population dynamics, three application examples are presented corresponding to the following key areas: distributed optimization, classic games, and distributed control design. function, which is given by the gradient of the corresponding cost function.

An interesting future direction is the extension of the results given in this dissertation to deal with directed and time-varying graphs. This extension is important since directed and time-varying graphs emerge in a number of real applications. For instance, in networks of mobile agents [31], where communications between agents vary depending on agents' relative positions. This proposed future direction is challenging since the analysis of distributed resource allocation algorithms highly depends on the graph structure.

Optimal Resource Allocation with Lower-Bound Constraints

We have developed a distributed method that solves a class of resource allocation problems with lower bound constraints. The proposed approach is based on a multi-agent system, where coordination among agents is done by using a consensus protocol. We have proved that convergence and optimality of the method is guaranteed under some mild assumptions, specifically, we require that the cost function is strictly convex and the graph related to the communication network that enables the agents to share information is connected. The main advantage of our technique is that it does not need a centralized coordinator, which makes the method appropriate to be applied in largescale distributed systems, where the inclusion of centralized agents is undesirable or infeasible. As future work, we propose to use a switched approach in order to eliminate the iterations in Algorithm 1. Moreover, we plan to include upper bound constraints in our original formulation.

Distributed Population Dynamics

We have generalized a methodology to generate distributed population dynamics from the distributed mean dynamics and different revision protocols. To illustrate the method, we have derived four distributed population dynamics, i.e., the distributed Smith dynamics, the distributed replicator dynamics, the distributed projection dynamics and the distributed logit dynamics. We have proved that the distributed population dynamics presented in this work exhibit mass conservation and convergence to Nash equilibria in the same way as classic population dynamics do. Finally, we have presented

Abstract

Since the complexity and scale of systems have been growing in the last years, distributed approaches for control and decision making are becoming more prevalent. This dissertation focuses on an important problem involving distributed control and decision making, the dynamic resource allocation in a network.

To address this problem, we explore a consensus-based algorithm that does not require any centralized computation, and that is capable to deal with applications modeled either by dynamical systems or by memoryless functions. The main contribution of our research is to prove, by means of graph theoretical tools and passivity analysis, that the proposed controller asymptotically reaches an optimal solution without the need of full information.

In order to illustrate the relevance of our main result, we address several engineering applications including: distributed control for energy saving in smart buildings, management of the customers of an aggregating entity in a smart grid environment, and development of an exact distributed optimization method that deals with resource allocation problems subject to lower--bound constraints.

Finally, we explore resource allocation techniques based on classic population dynamics models. In order to make them distributed, we introduce the concept of non--well--mixed population dynamics. We show that these dynamics are capable to deal with constrained information structures that are characterized by non-complete graphs. Although the proposed non--well-mixed population dynamics use partial information, they preserve similar properties of their classic counterpart, which uses full information. Specifically, we prove mass conservation and convergence to Nash equilibrium.

Key Words

Resource allocation, output consensus, distributed algorithms, graph theory, passivity, multi-agent systems, population dynamics.

Resumen

Dado que la complejidad y la escala de los sistemas se han ido incrementando en los últimos años, las técnicas centralizadas de control y toma de decisiones están siendo reemplazadas por métodos distribuidos. Esta tesis se centra en un importante problema que involucra control y toma de decisiones distribuidas: la asignación dinámica de recursos en redes. Para abordar este problema, exploramos un algoritmo basado en consenso que no requiere computación centralizada, y que puede ser usado en aplicaciones modeladas ya sea por sistemas dinámicos o funciones sin memoria. La principal contribución de esta tesis es probar, por medio de teoría de grafos y pasividad, que el algoritmo propuesto alcanza asintóticamente una solución óptima sin la necesidad de usar información completa.

Para ilustrar la relevancia del resultado principal de esta disertación, abordamos varias aplicaciones en ingeniería, incluyendo: el control distribuido en edificios inteligentes orientado a la eficiencia energética, la gestión de los clientes de un agregador en una red inteligente en la que se aplican estrategias de respuesta de la demanda, y el desarrollo de un método de optimización exacto que permite incluir restricciones de límite inferior.

Finalmente, se exploran otras técnicas de asignación de recursos inspiradas en modelos de dinámicas poblacionales. Se introduce el concepto de poblaciones no-bien-mezcladas, y se muestra que las dinámicas asociadas a este tipo de poblaciones cuentan con una estructura de información local, caracterizada por grafos que no son completos. A pesar de que las dinámicas propuestas usan información parcial, ellas preservan características similares a las dinámicas poblacionales clásicas que usan información completa.

Palabras Clave

Asignación de recursos, consenso, algoritmos distribuidos, teoría de grafos, pasividad, sistemas multiagente, dinámicas poblacionales.