# Decentralized Optimization Algorithms for Variable Speed Pumps Operation Based on Local Interaction Game

Authors: Shi-Qiang Wang; Jian-Chun Xing; Zi-Yan Jiang; Yun-Chuang Dai
Journal: Journal of Control Science and Engineering
 (2018)

---

## Abstract

A fully distributed optimal control strategy for the parallel variable speed pumps in heating, ventilation, and air-conditioning (HVAC) systems is proposed. Compared with the traditional centralized method, the efficient control and coordination are scattered to every updated smart pump without the need for a monitoring host. Similar to the structure, mechanism, and characteristics of biological communities, a smart pump can communicate with adjacent nodes and operate collaboratively to complete pumps group operation with the least total power consumption under load demand and system constraints. And a decentralized optimization method that is decentralized estimation of distribution algorithm (DEDA) under local interaction games framework has been transplanted to the proposed structure to optimize the pumps working in parallel mode. Besides, convergence property of the two novel mechanisms is analyzed theoretically. Finally, simulation studies have been conducted based on the models of a physical pumps system, and the performance of the proposed algorithm is compared with centralized algorithms in terms of both effectiveness and stability. The results provide support for the validity of the proposed algorithms and control structure.


---

## Body

## 1. Introduction

A pumping system, which generally consists of several circulation pumps, is one of the main energy consuming parts in a chilled water based HVAC system. According to the International Energy Agency, electrical motors have consumed 46% of global electricity demand by end-use [1]. Furthermore, this rate is somewhat higher in the industrial area: up to almost 65% of generated electrical power [2]. And for the HVAC system, 60–70% of the power consumption is attributed to the water pumps and fan-based energy transmission in the water distribution system. These data has reflected the fact that there exists huge potential for energy conserving. And energy efficiency of a pumping system depends not only on the design but also, even more, on its operation. Thus an effective optimization should be proposed to improve the operating efficiency.Currently, frequency converters are widely used with induction motor-driven pumps, enabling the pumps speed variation to control the head and flow rate of pumps under different working conditions in an energy-efficient way, that is, variable speed pumps (VSPs). However, it makes the moderating process more complicated in that not only the combination of pumps but also the frequency of each pump needs to be determined compared to the nonvariable-speed pumps.On the one hand, the physical models of pumps operation are difficult to be designed accurately and the pumps unit commitment is always a mixed integer programming problem, which leads to conventional optimization methods difficult to be used. In this case, the majority of existing researches focus on adopting different kinds of heuristic algorithms to solve such problems. Kecebas et al. utilized the artificial neural network to optimize the geothermal district heating systems [3]. Edson Bortoni had used the dynamic programming algorithm for optimizing a complex pumping system with a set of parallel centrifugal pumps and the efficiency of the method is verified through experiment[4]. Wang et al. introduced an enhanced evolutionary algorithm for biobjective pump scheduling in water supply [5]. Arun Shankar has conducted a comprehensive review of different methods for energy efficiency enhancement initiatives in centrifugal pumping system in recent years, and several practical evolutionary methods were also introduced [6]. However, the theory analysis of heuristic algorithms is not easy so that its convergence characteristic has not been studied well. Therefore, a more valid solving method needs to be applied to the pumps control in some critical application.On the other hand, in most recent studies, the keynote has been put on updating the existing optimization methods or developing more effective ones[7, 8]. However, most of the current algorithms are developed based on a centralized structure, in which all the sensors and actuators are connected to one master controller and the field bus technology is usually used to handle the communication. Figure 1 has shown a typical physical structure of an HVAC system (Figure 1(a)) and a corresponding centralized control system (Figure 1(b))[9]. The sensors and actuators of each separate equipment are gathered to the master controller level by level through the linkers (RS485/Ethernet <-> analog quantities) on the field bus.Figure 1
Structure of the centralized control system of a typical HVAC system.
(a)
Structure of the physical system (b)
Structure of the centralized control systemFrom the engineering perspective, the centralized control structure is not consistent with the physical system, so a mapping mechanism between the two systems is needed, such as configuring and commissioning. Therefore, the development incurs high maintenance and labor costs and somehow diminishes the flexibility and universality of a construction.After that the key procedure is to write the control model and algorithm into the monitoring center. The mathematical model for the total system is always complex and difficult to be obtained. In the algorithm writing process, a lot of work is needed such as system recognizing, software developing, and operating debugging. Moreover, traditional centralized methods require information transmission to the supervising computer during operations, which can cause severe link congestion and operational lag. Besides, since all the devices are managed by the master controller, malfunctions in one device will also bring about influence on other devices which will cause the entire control system running abnormally or even paralysis.Recently, distributed multiagent control systems which are being focused on by many learners have provided a possible solution to this problem [10]. And Andrew Windham has analyzed the application prospects of the multiagent in HVAC control system[11]. However, the structures still require several local central nodes which are not fully distributed in our sense. And, in fact, the challenge is how to achieve the same performance as in a centralized manner where only local state information is available via each smart unit and is not considered in the above studies[12].In this paper, a novel distributed control structure where each pump and chiller are updated to a smart unit by an embedded controller is proposed. These independent units can work independently and cooperate together to fulfill the global optimizing task.To achieve the global optimization in a fully distributed way, a local interaction game is introduced into cognitive ratio networks in [13]. And using the log-linear algorithm[14–18], the network throughput is maximized. However, the log-linear algorithm is not applicable to the pumps optimal load distribution. Therefore, learning from the local state interaction games concept, a decentralized estimation of distribution algorithm (DEDA) is proposed in this paper. Simulation studies have been conducted based on physical models of the chilled water pumps system in a refrigeration station. The performance of DEDA is validated through comparison with centralized manner.The rest of this paper is outlined as follows: In Section2 the pumps system model is presented. Besides, potential game is introduced. Then a decentralized optimization model adaption is formulated in Section 3. And a decentralized algorithm is derived to solve the adapted mathematical model in this section. Section 4 focuses on performance assessment through simulations and experiments. Conclusions are drawn in Section 5.


## 2. System Model and Problem Formulation

### 2.1. Pumps Model and Its Decentralized Control System

The pump head and flow rate of a single VSPi (i=1,2,…,N) running at nonrated speed can be scaled from those of the pump running at the rated speed. According to the affinity laws of centrifugal pumps, the similarity working conditions of a VSP can be determined as(1)QinQin0=nn0,HinHin0=nn02,ηinηin0=1wheren0 is the rated pump speed andn is the pump speed under nonrated conditions.Qi is the flow rate,Hi is the pump head of each branch, andηi is the efficiency.The speed ratiowi for each pumpi is described in (2) and w0 is the lower limitation of the speed ratio.(2)wi=nn0w0<wi<1And, the fitting formulas of pumps running at nonrated speed can be expressed as(3)Hset=Hi=aiQi2+biwiQi+ciwi2,i=1,2,…,N,(4)ηi=jiQiwi2+kiQiwi+li,i=1,2,…,N,(5)W=∑i=1NWi=∑i=1NρgHiQiηi,whereHset is the demanded pressure difference on the main pipe.W is the total energy consumption of the pumps andWi=ρgHiQi/ηi.ρ is the density of water, g is the gravitational acceleration,N is the quantity of pumps.ai,bi,ci andji,ki,li are static physical characteristics of different branches for pumpi(i=1,2,…,N).If the detected pressure differenceHo is not equal toHset, the control system will calculate the demanded total flow rateQp as shown in (6). After that the moderating process will be conducted to ensure the pressure difference by moderatingwi and make the calibration of eachQi satisfyQp.(6)Ho=SQo,(7)Hset=S∑Qi2=SQp2whereSis the resistance of the current pipe net andHo andQo are the current pressure difference and flow rate detected by the sensors.In the largest pump systems, the pump groups are composed by several identical big pumps and one or two small pumps for energy saving. As shown in Figure2, the decentralized control system consists ofN smart updated pump by embedded microprocessors. A flow meter and a pressure meter will be connected to one of the smart units to detect the current state of the pipe net.Figure 2
Decentralized control structure of a parallel pump system.Each smart pump is connected hand to hand with its neighbor building up a chain-like form which corresponds to its actual engineering installation layout as shown in Figure2. The model of each pump is stored in its smart unit as local variables and there is no need to model for the entire system as that in the centralized structure. And the control algorithm and instructions will be written into each smart pump and are identical. The load difference between supply and demand will be transmitted back and forth from one smart unit to another on the chain. And the operating state of each smart pump will be updated with the load difference between supply and demand transmitted from its neighbors. In such way, each smart device can be removed or substituted with little configuration making the system convenient to update and extend.


### 2.2. Problem Formulation

The operation object can be described as finding a combination of the running pumps and different pump can adopt a number of speed ratio strategies that can satisfyQp with the least energy consumption. Thus, it can be built up as(8)minWs.t.Qidown<Qi<Qiup,∑i=1NQi=Qpi=1,2,…,Nwhere Qidown and Qiup denote the power outputs limit of pumpi.Implementation of the centralized model (8) requires a central processor with the total system information. The fully distributed approaches remove the requirement of monitoring host and enable each updated smart DG to complete optimal load distribution task only based on local and neighboring information. In this case, the centralized model (8) is split on each smart pump to be executed in a fully distributed and parallel manner across every processor core in each pump. Each parti(i=1, 2, …,N) is a suboptimization problem of (8) as decentralized optimization.Definition 1 (decentralized optimization for pumps load distribution).
The decentralized optimization problem for each nodei (i=1, 2, …,N) is defined as(9)minWi(10)s.t.Qidown<Qi<Qiup,∑i=1NQi=Qpi=1,2,…,N.Indeed, (9) corresponds to part of the objective function (5). Note that each energy consumptionWi of the pumpi (i=1, 2,…,N) is locally known by smart pumpionly. The penalty term is constructed by load difference and added to the utility function as (11) for pumpi (i=1, 2,…,N).(11)minUi=Wi+Pi=Wi+βi∑i=1NQi-Qp2=Wi+βiΔQ2where β is the penalty coefficient. The parameter β balances the tradeoff between exploration and exploitation. Smaller β implies that the smart pumps are more willing to choose a suboptimal action to explore, whereas higher β implies that they are prone to choose the best response action. Therefore, it is advisable that, at the beginning phase, the value of β is set to small number and it keeps increasing as the optimization iterates[19, 20]. ΔQis the supply-demand difference. Then, (11) means that each pumpi should select a certain speed ratio strategy from its speed space to satisfyQp with the least energy consumption.Thus, the key problem lies in designing the mechanism to achieve global optimization via local interaction among nodes (the terms of agent, node, and smart pump refer to the same concept in this paper) only with direct communication links. As we all know, it is local interaction game. In this manner, optimal load distribution can be completed on the field level. In this paper, any method that solves (9)-(10) without using a central node aggregating data is called decentralized algorithm.


### 2.3. Potential Game

Coordination and combinatorial optimization, when viewed from a game theory perspective, can form a subset of the class of potential games.Definition 2 (potential games, [21]).
Player action sets{Si}i=1N, and ⋃i=1NSi=S, with player objective functions;{Ui:S→R}i=1N is a potential game if, for some potential function Φ:S→R,(12)Uisi,s-i-Uisi′,s-i=ϕsi,s-i-ϕsi′,s-i∀si′,s-i∈Si∀i∈1,2,…,Nfor every playerPi∈P where P ={P1,…,PN} is a group of players (or agents) ands-i is the complimentary set ofsi.And a number of learning algorithms can be found for the potential game. Most of the learning algorithms for potential games guarantee convergence to a (pure) Nash equilibrium.Definition 3 (Nash equilibrium).
A joint strategy profile,s∗, such that no individual agent has an incentive to change to a different strategy, that is,(13)Uisi∗,s-i∗-Uisi,s-i∗≥0∀si∈Si∀i∈1,2,…,N,is a Nash equilibrium.It is easy to see that, in potential games, any action profile maximizing the potential function is a pure Nash equilibrium. Hence, every potential game possesses at least one such equilibrium.Consequently, each smart pump can be treated as a player and the energy consumption (11) is the utility function in potential game. And a potential game, as previously defined, requires perfect alignment between the global objective and the players’ local objective functions in the following sense: if a player unilaterally changed its action, the change in its objective function would be equal to the change in the potential function. The pumps total cost function (7) is taken as the potential function in the optimization. Obviously, (9) subjected to (10) is a potential game.However, there may also exist suboptimal pure Nash equilibrium that do not maximize the potential function[22]. Therefore, the key problem lies in how to design the mechanism to achieve the global optimization in a fully distributed way for the VSPs.


## 2.1. Pumps Model and Its Decentralized Control System

The pump head and flow rate of a single VSPi (i=1,2,…,N) running at nonrated speed can be scaled from those of the pump running at the rated speed. According to the affinity laws of centrifugal pumps, the similarity working conditions of a VSP can be determined as(1)QinQin0=nn0,HinHin0=nn02,ηinηin0=1wheren0 is the rated pump speed andn is the pump speed under nonrated conditions.Qi is the flow rate,Hi is the pump head of each branch, andηi is the efficiency.The speed ratiowi for each pumpi is described in (2) and w0 is the lower limitation of the speed ratio.(2)wi=nn0w0<wi<1And, the fitting formulas of pumps running at nonrated speed can be expressed as(3)Hset=Hi=aiQi2+biwiQi+ciwi2,i=1,2,…,N,(4)ηi=jiQiwi2+kiQiwi+li,i=1,2,…,N,(5)W=∑i=1NWi=∑i=1NρgHiQiηi,whereHset is the demanded pressure difference on the main pipe.W is the total energy consumption of the pumps andWi=ρgHiQi/ηi.ρ is the density of water, g is the gravitational acceleration,N is the quantity of pumps.ai,bi,ci andji,ki,li are static physical characteristics of different branches for pumpi(i=1,2,…,N).If the detected pressure differenceHo is not equal toHset, the control system will calculate the demanded total flow rateQp as shown in (6). After that the moderating process will be conducted to ensure the pressure difference by moderatingwi and make the calibration of eachQi satisfyQp.(6)Ho=SQo,(7)Hset=S∑Qi2=SQp2whereSis the resistance of the current pipe net andHo andQo are the current pressure difference and flow rate detected by the sensors.In the largest pump systems, the pump groups are composed by several identical big pumps and one or two small pumps for energy saving. As shown in Figure2, the decentralized control system consists ofN smart updated pump by embedded microprocessors. A flow meter and a pressure meter will be connected to one of the smart units to detect the current state of the pipe net.Figure 2
Decentralized control structure of a parallel pump system.Each smart pump is connected hand to hand with its neighbor building up a chain-like form which corresponds to its actual engineering installation layout as shown in Figure2. The model of each pump is stored in its smart unit as local variables and there is no need to model for the entire system as that in the centralized structure. And the control algorithm and instructions will be written into each smart pump and are identical. The load difference between supply and demand will be transmitted back and forth from one smart unit to another on the chain. And the operating state of each smart pump will be updated with the load difference between supply and demand transmitted from its neighbors. In such way, each smart device can be removed or substituted with little configuration making the system convenient to update and extend.


## 2.2. Problem Formulation

The operation object can be described as finding a combination of the running pumps and different pump can adopt a number of speed ratio strategies that can satisfyQp with the least energy consumption. Thus, it can be built up as(8)minWs.t.Qidown<Qi<Qiup,∑i=1NQi=Qpi=1,2,…,Nwhere Qidown and Qiup denote the power outputs limit of pumpi.Implementation of the centralized model (8) requires a central processor with the total system information. The fully distributed approaches remove the requirement of monitoring host and enable each updated smart DG to complete optimal load distribution task only based on local and neighboring information. In this case, the centralized model (8) is split on each smart pump to be executed in a fully distributed and parallel manner across every processor core in each pump. Each parti(i=1, 2, …,N) is a suboptimization problem of (8) as decentralized optimization.Definition 1 (decentralized optimization for pumps load distribution).
The decentralized optimization problem for each nodei (i=1, 2, …,N) is defined as(9)minWi(10)s.t.Qidown<Qi<Qiup,∑i=1NQi=Qpi=1,2,…,N.Indeed, (9) corresponds to part of the objective function (5). Note that each energy consumptionWi of the pumpi (i=1, 2,…,N) is locally known by smart pumpionly. The penalty term is constructed by load difference and added to the utility function as (11) for pumpi (i=1, 2,…,N).(11)minUi=Wi+Pi=Wi+βi∑i=1NQi-Qp2=Wi+βiΔQ2where β is the penalty coefficient. The parameter β balances the tradeoff between exploration and exploitation. Smaller β implies that the smart pumps are more willing to choose a suboptimal action to explore, whereas higher β implies that they are prone to choose the best response action. Therefore, it is advisable that, at the beginning phase, the value of β is set to small number and it keeps increasing as the optimization iterates[19, 20]. ΔQis the supply-demand difference. Then, (11) means that each pumpi should select a certain speed ratio strategy from its speed space to satisfyQp with the least energy consumption.Thus, the key problem lies in designing the mechanism to achieve global optimization via local interaction among nodes (the terms of agent, node, and smart pump refer to the same concept in this paper) only with direct communication links. As we all know, it is local interaction game. In this manner, optimal load distribution can be completed on the field level. In this paper, any method that solves (9)-(10) without using a central node aggregating data is called decentralized algorithm.


## 2.3. Potential Game

Coordination and combinatorial optimization, when viewed from a game theory perspective, can form a subset of the class of potential games.Definition 2 (potential games, [21]).
Player action sets{Si}i=1N, and ⋃i=1NSi=S, with player objective functions;{Ui:S→R}i=1N is a potential game if, for some potential function Φ:S→R,(12)Uisi,s-i-Uisi′,s-i=ϕsi,s-i-ϕsi′,s-i∀si′,s-i∈Si∀i∈1,2,…,Nfor every playerPi∈P where P ={P1,…,PN} is a group of players (or agents) ands-i is the complimentary set ofsi.And a number of learning algorithms can be found for the potential game. Most of the learning algorithms for potential games guarantee convergence to a (pure) Nash equilibrium.Definition 3 (Nash equilibrium).
A joint strategy profile,s∗, such that no individual agent has an incentive to change to a different strategy, that is,(13)Uisi∗,s-i∗-Uisi,s-i∗≥0∀si∈Si∀i∈1,2,…,N,is a Nash equilibrium.It is easy to see that, in potential games, any action profile maximizing the potential function is a pure Nash equilibrium. Hence, every potential game possesses at least one such equilibrium.Consequently, each smart pump can be treated as a player and the energy consumption (11) is the utility function in potential game. And a potential game, as previously defined, requires perfect alignment between the global objective and the players’ local objective functions in the following sense: if a player unilaterally changed its action, the change in its objective function would be equal to the change in the potential function. The pumps total cost function (7) is taken as the potential function in the optimization. Obviously, (9) subjected to (10) is a potential game.However, there may also exist suboptimal pure Nash equilibrium that do not maximize the potential function[22]. Therefore, the key problem lies in how to design the mechanism to achieve the global optimization in a fully distributed way for the VSPs.


## 3. Decentralized Optimization Based on Probability Parallel Variable-Frequency Speed Pumps

In fact, there exist many evolutionary and swarm intelligent algorithms for coordination and combinatorial optimization in engineering[23, 24], but they might not be applicable to decentralized pumps control system due to a number of reasons. Firstly, the existing evolutionary or swarm intelligent algorithm refers to quite identical individualities. That is, every population in swarm intelligent algorithm or chromosomes in evolutionary algorithm involves a global variable, but only local information is available in the decentralized control structure of a parallel pump system. As an exploration, an example of multiple autonomous robots following simple local rules and collectively provided by [25] is quite similar to the principles of the decentralized control system in this paper. However, the challenge is how to achieve the same performance as in a centralized manner where only local information is available.


### 3.1. Decentralized Estimation of Distribution Algorithm

Estimation of distribution algorithm (EDA) is an evolutionary optimization method that guides the search for the optimum by building and sampling explicit probabilistic models of promising candidate solutions[26]. Provided that a solutionx hasN factorsxi,i= 1, 2, …,N, the probability distribution of observingx is estimated by the product of density for each factor as follows:(14)Px=∏i=1NPθixi,whereθi is the vector of parameters for the probability density of factorxi.The EDA is adapted into decentralized fashion in this study. DEDA is different from traditional evolutionary or swarm intelligent algorithms, where the populations of each smart pump only cover its own variable as shown in Figure3.Figure 3
Comparison between EDA and decentralized optimization algorithm.And the detailed description of DEDA is illustrated in Figure4. Further, from the potential game theory perspective, an adjustment mechanism can typically produce sequential or parallel actions by agents (smart pump). Observed by the characteristic of decentralized parallel pumps control structure, sequential random schedule with local interaction is adopted in the present study. And an arbitrary agent can be picked as the indicator in the pumps optimal load distribution.Figure 4
Flow chart of the DEDA in each smart pump and executed in serial mode.


### 3.2. Convergence Proof

The convergence proof of the DEDA is presented in Proposition4 and illustrated in the following paragraph.Proposition 4.
The solution xde= [x1de, …, xide, …,xNde] for the decentralized optimization (11) generated by DEDA can converge to the optimal Nash equilibrium with a probability of 1.Proof.
Assume the global solution space is defined asΩ=Ω1 ∩ Ω2 ∩ …∩ΩN, whereΩi,i=1, 2, …,N, is the space of theith agent and its corresponding variablexi. An elite-preserving operator is used in the presented algorithm. This operator favors theg elites of a population by giving them an opportunity to be directly carried over to the next generation. Thus, the new population consists of (L−g) new chromosomes, withg<L.
Firstly, the optimality convergence of each agenti(i=1, 2, …,N) is discussed. The present study assumes that the search spaceΩi of problem (14) is a closed and bounded region. The objective function is taken as the fitness functionJi(xi), and the optimal solution set is expressed as(15)D0=xi∈Ωi∣Jixi-Ji∗<εwhere Ji∗ = minJi(xi): xi∈Ωi.
For agenti (i=1, 2,…,N), theL chromosomes of EDA can be divided into the following properties:
(I) At least one chromosome belongs toD0, which is written asS0.
(II) AllL chromosomes belong toD1, which is written asS1, where(16)D1=Ω∖D0.
And ifqij (i,j= 0,1) represents the probability that the (k+1)th generationxi(k+1) transfers toSj when thek-th generationxi(k) keeps the state ofSi, then the following conditions hold:
(I)xi (k) atS0 andq00 =1.
Condition (I) is obviously valid because the elite-preserving operator is used in this paper.
(II)xi (k) atS1 andq11 ≤c, wherec∈(0,1).
Conditions (II) is analyzed as follows: For anyε> 0, ifxi0∈{xi∈Ωi∣Jixi-Ji∗<ε/2}, there existsr > 0 satisfying(17)Jixi-Jixi0<ε2when(18)xi-xi0∞=maxxi-xi0≤r.
Define(19)Qxi0,r=xi∈Ωi∣xi-xi0∞≤r,so that(20)Qxi0,r⊂D0due to(21)Jixi-Ji∗≤Jixi-Jixi0+Jixi0-Ji∗<ε2+ε2=ε.
The probability thatxi transfers toS0 can be defined by(22)Pxi+∈Qxi0,r=Pxi+-xi0≤r,where xi+ is the chromosome of a new generation andi= 1, 2, …,N. Given that the evolutionary computation of EDA is based on probability, each element of the chromosome generates a series of the probability function g(xi), withi=1, 2, …,N.
Forxi+, (23)xi+~gxiand thus(24)Pxi+∈Qxi0,r=∫xi0-rxi0+rgxidxi.
For simplicity,Pxi+∈Qxi0,r can be denoted by(25)P1xi+=Pxi+∈Qxi0,r,and then(26)0<P1xi+<1.
Define(27)P1yi0=min⁡P1xi+,so that(28)P1yi0≤P1xi+≤q10.
As(29)q11+q10=1,the following result can be achieved:(30)q11=1-q10≤c0<c<1,where the constantc is(31)c=1-P1yi0.
Further, for anyε>0, set(32)pk=PJixi∗k-Ji∗≥εand(33)pk=0∃T∈0,1,2,…,k,xi∗T∈D0pk¯,xi∗m∉D0m=0,1,2,…,k,and, thus, the following equation is obtained:(34)pk¯=Pxi∗m∉D0,m=0,1,2,…,k=q11k≤ck,so that(35)∑k=1∞pk¯≤∑k=1∞ck=c1-c<∞.
Then the following theorem is introduced to complete the proof process.
T h e o r e m   5 (Borel-Cantelli theorem).The present study assumes that A1, A2,… is the event sequence based on probability and definespk= P{Ak}. If(36)∑k=1∞pk<∞,
then (37) P ⋂ m = 1 ∞ ⋃ k ≥ m A k = 0 .
If (38) ∑ k = 1 ∞ p k = ∞
and A k are independent of one another, then (39) P ⋂ m = 1 ∞ ⋃ k ≥ m A k = 1 .According to Theorem 5,(40)P⋂m=1∞⋃k≥mJixi∗k-Ji∗≥ε=0is true. Therefore,xit(m), wheret= 1, 2, …,L, is the population of themth generation of EDA for agenti (i=1,2,…,N), which consists ofL chromosomes. If the fitness functionJi(xi) is continuous on the bounded regionΩi, the optimal solution xi∗(m) converges to the optimum for agentiwith a probability of 1, and(41)xi∗m=argmin1≤t≤LJixitm.Then, the validity of the global optimization is deduced in the subsequent paragraphs.According to the above analysis, the solutionxide can converge to the optimum xi∗de for agentiwith a probability of 1. Define the probability aspi, that is,(42)xide⟶pixi∗de,and(43)pi=P⋂m=1∞⋃k≥mJxidek-Ji∗<ε=1.According to (24), the global combined probability density can be described as(44)P=∏i=1Npi=1.For the coupled variables, a similar result is also obtained:(45)P=∏i=1Npipij=1.Then the present study directly obtains(46)xde⟶Px∗dewhere x∗de = [x1∗de, …, xi∗de, …,xN∗de] is a Nash equilibrium of problem (11). In conclusion,xde converges to an optimal solution of the centralized optimization problem with a probability of 1 as the iteration stepm→∞. Then Proposition 4 can be obtained as follows.


## 3.1. Decentralized Estimation of Distribution Algorithm

Estimation of distribution algorithm (EDA) is an evolutionary optimization method that guides the search for the optimum by building and sampling explicit probabilistic models of promising candidate solutions[26]. Provided that a solutionx hasN factorsxi,i= 1, 2, …,N, the probability distribution of observingx is estimated by the product of density for each factor as follows:(14)Px=∏i=1NPθixi,whereθi is the vector of parameters for the probability density of factorxi.The EDA is adapted into decentralized fashion in this study. DEDA is different from traditional evolutionary or swarm intelligent algorithms, where the populations of each smart pump only cover its own variable as shown in Figure3.Figure 3
Comparison between EDA and decentralized optimization algorithm.And the detailed description of DEDA is illustrated in Figure4. Further, from the potential game theory perspective, an adjustment mechanism can typically produce sequential or parallel actions by agents (smart pump). Observed by the characteristic of decentralized parallel pumps control structure, sequential random schedule with local interaction is adopted in the present study. And an arbitrary agent can be picked as the indicator in the pumps optimal load distribution.Figure 4
Flow chart of the DEDA in each smart pump and executed in serial mode.


## 3.2. Convergence Proof

The convergence proof of the DEDA is presented in Proposition4 and illustrated in the following paragraph.Proposition 4.
The solution xde= [x1de, …, xide, …,xNde] for the decentralized optimization (11) generated by DEDA can converge to the optimal Nash equilibrium with a probability of 1.Proof.
Assume the global solution space is defined asΩ=Ω1 ∩ Ω2 ∩ …∩ΩN, whereΩi,i=1, 2, …,N, is the space of theith agent and its corresponding variablexi. An elite-preserving operator is used in the presented algorithm. This operator favors theg elites of a population by giving them an opportunity to be directly carried over to the next generation. Thus, the new population consists of (L−g) new chromosomes, withg<L.
Firstly, the optimality convergence of each agenti(i=1, 2, …,N) is discussed. The present study assumes that the search spaceΩi of problem (14) is a closed and bounded region. The objective function is taken as the fitness functionJi(xi), and the optimal solution set is expressed as(15)D0=xi∈Ωi∣Jixi-Ji∗<εwhere Ji∗ = minJi(xi): xi∈Ωi.
For agenti (i=1, 2,…,N), theL chromosomes of EDA can be divided into the following properties:
(I) At least one chromosome belongs toD0, which is written asS0.
(II) AllL chromosomes belong toD1, which is written asS1, where(16)D1=Ω∖D0.
And ifqij (i,j= 0,1) represents the probability that the (k+1)th generationxi(k+1) transfers toSj when thek-th generationxi(k) keeps the state ofSi, then the following conditions hold:
(I)xi (k) atS0 andq00 =1.
Condition (I) is obviously valid because the elite-preserving operator is used in this paper.
(II)xi (k) atS1 andq11 ≤c, wherec∈(0,1).
Conditions (II) is analyzed as follows: For anyε> 0, ifxi0∈{xi∈Ωi∣Jixi-Ji∗<ε/2}, there existsr > 0 satisfying(17)Jixi-Jixi0<ε2when(18)xi-xi0∞=maxxi-xi0≤r.
Define(19)Qxi0,r=xi∈Ωi∣xi-xi0∞≤r,so that(20)Qxi0,r⊂D0due to(21)Jixi-Ji∗≤Jixi-Jixi0+Jixi0-Ji∗<ε2+ε2=ε.
The probability thatxi transfers toS0 can be defined by(22)Pxi+∈Qxi0,r=Pxi+-xi0≤r,where xi+ is the chromosome of a new generation andi= 1, 2, …,N. Given that the evolutionary computation of EDA is based on probability, each element of the chromosome generates a series of the probability function g(xi), withi=1, 2, …,N.
Forxi+, (23)xi+~gxiand thus(24)Pxi+∈Qxi0,r=∫xi0-rxi0+rgxidxi.
For simplicity,Pxi+∈Qxi0,r can be denoted by(25)P1xi+=Pxi+∈Qxi0,r,and then(26)0<P1xi+<1.
Define(27)P1yi0=min⁡P1xi+,so that(28)P1yi0≤P1xi+≤q10.
As(29)q11+q10=1,the following result can be achieved:(30)q11=1-q10≤c0<c<1,where the constantc is(31)c=1-P1yi0.
Further, for anyε>0, set(32)pk=PJixi∗k-Ji∗≥εand(33)pk=0∃T∈0,1,2,…,k,xi∗T∈D0pk¯,xi∗m∉D0m=0,1,2,…,k,and, thus, the following equation is obtained:(34)pk¯=Pxi∗m∉D0,m=0,1,2,…,k=q11k≤ck,so that(35)∑k=1∞pk¯≤∑k=1∞ck=c1-c<∞.
Then the following theorem is introduced to complete the proof process.
T h e o r e m   5 (Borel-Cantelli theorem).The present study assumes that A1, A2,… is the event sequence based on probability and definespk= P{Ak}. If(36)∑k=1∞pk<∞,
then (37) P ⋂ m = 1 ∞ ⋃ k ≥ m A k = 0 .
If (38) ∑ k = 1 ∞ p k = ∞
and A k are independent of one another, then (39) P ⋂ m = 1 ∞ ⋃ k ≥ m A k = 1 .According to Theorem 5,(40)P⋂m=1∞⋃k≥mJixi∗k-Ji∗≥ε=0is true. Therefore,xit(m), wheret= 1, 2, …,L, is the population of themth generation of EDA for agenti (i=1,2,…,N), which consists ofL chromosomes. If the fitness functionJi(xi) is continuous on the bounded regionΩi, the optimal solution xi∗(m) converges to the optimum for agentiwith a probability of 1, and(41)xi∗m=argmin1≤t≤LJixitm.Then, the validity of the global optimization is deduced in the subsequent paragraphs.According to the above analysis, the solutionxide can converge to the optimum xi∗de for agentiwith a probability of 1. Define the probability aspi, that is,(42)xide⟶pixi∗de,and(43)pi=P⋂m=1∞⋃k≥mJxidek-Ji∗<ε=1.According to (24), the global combined probability density can be described as(44)P=∏i=1Npi=1.For the coupled variables, a similar result is also obtained:(45)P=∏i=1Npipij=1.Then the present study directly obtains(46)xde⟶Px∗dewhere x∗de = [x1∗de, …, xi∗de, …,xN∗de] is a Nash equilibrium of problem (11). In conclusion,xde converges to an optimal solution of the centralized optimization problem with a probability of 1 as the iteration stepm→∞. Then Proposition 4 can be obtained as follows.


## 4. Case Study

### 4.1. Experimental Setup

This case study is carried out to validate the performance of DOABP by applying it to the model of a physical pump system under different load demands. The physical parameters (PAR) corresponding to performance curves of the pumps are shown in Table1, and these parameters are recognized from the practical operating data of a district freezing station. As seen from Table 1, the performances of pumps of the same type are close to each other. The pump system is composed of four pumps working in parallel mode including two big pumps (KP1020-3 457mm in impeller diameter) and two smaller pumps (KP1020-3 406mm in impeller diameter) manufactured by Grundfos. A determined pressure difference control is applied here where the designed pressure difference and flow rate demand for the entire system are 45m and 3234m3/h, respectively.Table 1
Pumps parameters.
PAR PUMP 1.KP1020-3 406 2.KP1020-3 406 3.KP1020-3 457 4.KP1020-3 457 a -6.07e-6 -5.36e-6 -8.54e-6 -8.41e-6 b -0.00822 -0.00949 -0.00048 -0.00046 c 60.41767 60.61789 74.54932 73.55836 j -4.08e-7 -3.96e-7 -3.87e-7 -3.82e-7 k 0.000940 0.000911 0.000999 0.000985 l 0.316357 0.328377 0.249898 0.259137


### 4.2. Simulation Result

The strategy space is discretized, and it becomes a limited strategy space to ensure that the potential game is a finite game. And, the speed ratio (SR) of each pump, as selected action set, is shown in Table2.Table 2
Pumps discreted strategy space.
SR PUMP 1.KP1020-3 406 2.KP1020-3 406 3.KP1020-3 457 4.KP1020-3 457 1 0 0 0 0 2 0.93 0.94 0.825 0.825 3 0.94 0.945 0.85 0.85 4 0.945 0.95 0.875 0.875 5 0.95 0.955 0.9 0.9 6 0.96 0.96 0.915 0.915 7 0.97 0.97 0.925 0.925 8 0.98 0.98 0.95 0.95 9 0.99 0.99 0.975 0.975 10 1 1 1 1 Note. The speed ratio of 0 means that the pump is off.In the case study, the proposed algorithms will be triggered as follows: when the pipe net has been changed, the smart unit will detect that the pressure differenceHo is not equal toHset andQo which can be got through the flow meter on the trunk; thusQi can be solved in (3).DEDA can be downloaded to microprocessor chip in each smart pumpi (i=1,2, 3,4), and then each smart pump can be active. And, each smart pump can communicate and cooperate with other smart pump nodes in serial mode. The corresponding results by the proposed DEDA are shown in Figure 5.Figure 5
Iterative process in decentralized optimization.As shown in Figure5, the probability distribution (PD) on each pump is dynamic during the iterative process in DEDA. And the optimal solution of each pump fluctuates and finally reaches a new equilibrium state under mutual negotiation and adjustment. By about 60 iteration steps, the pumps optimal load distribution task can be achieved. It can be seen that the probability of the optimal output flow rate tends to approach 1. Conversely, the probability of the other solution space approaches 0.The iterative process of each pumpi(i=1,2,3,4) cost is illustrated in Figure 6. Further, the iterative process of total energy cost is shown in Figure 7.Figure 6
Iterative process of cost on each pump in decentralized optimization.Figure 7
Iterative process of total cost in decentralized optimization.Comparing Figure6 and Figure 7, we find that although not all the costs of the four pumps decrease, the total cost decreases during the iterative process. This can validate the effectiveness of DEDA. For a comparison, the iterative process of EDA in centralized method is also illustrated in Figure 8. As shown in Figures 5 and 8, the DEDA can achieve the actual optima which are equivalent to that obtained with the centralized algorithm.Figure 8
Iterative process in centralized optimization.Further, a more detailed calculation result is shown in Table3, where EDA and DEDA are adopted, respectively.Table 3
Comparison between centralized and decentralized results.
#Working condition DEDA EDA Demand: 1496.9 m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2 Pump 3 Pump 4 47.000.00 41.25 0.00 700.170.00 796.71 0.00 0.790.00 0.85 0.00 47.000.00 41.25 0.00 700.170.00 796.71 0.00 0.790.00 0.85 0.00  Power/kW 227.4804 227.4804  Difference 0.0169 0.0169  Demand: 2583.4m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2Pump 3Pump 4 48.500.00 0.00 48.75 886.980.00 0.00 1695.44 0.830.00 0.00 0.82 48.500.00 0.00 48.75 886.980.00 0.00 1695.44 0.830.00 0.00 0.82  Power/kW 392.2045 392.2045  Difference 0.9772 0.9772  Demand: 2952.9m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2Pump 3Pump 4 47.2549.00 0.00 45.00 733.01926.41 0.00 1292.65 0.800.84 0 0.89 48.0047.50 40.00 45.75 827.24745.890.001379.76 0.820.8000.88  Power/kW 435.5667 439.1761  Difference 0.8373 0.0181  Demand: 3234.7m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2Pump 3Pump 4 50.0047.75 45.75 40.00 1054.53777.35 1402.09 0.00 0.850.81 0.87 0.00 48.5047.70 0.00 47.50 886.98777.34 0.00 1569.19 0.830.81 0.00 0.84  Power/kW 475.8191 485.6048  Difference 0.7344 1.1869  Demand: 3858.3m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2Pump 3Pump 4 47.0049.00 41.25 46.25 700.17926.41 796.71 1435.66 0.790.84 0.85 0.87 48.5050.00 41.25 43.75 886.981037.53 796.71 1136.63 0.830.85 0.85 0.89  Power/kW 572.8274 561.5738  Difference 0.6517 0.4537  Demand: 4901.2 m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2Pump 3Pump 4 49.0049.00 50.00 46.25 944.62926.41 1831.72 1435.66 0.840.84 0.78 0.87 48.0050.00 46.25 50.00 827.241037.51 1456.97 1815.93 0.820.85 0.86 0.79  Power/kW 778.3465 777.4647  Difference 0.4889 1.2384


## 4.1. Experimental Setup

This case study is carried out to validate the performance of DOABP by applying it to the model of a physical pump system under different load demands. The physical parameters (PAR) corresponding to performance curves of the pumps are shown in Table1, and these parameters are recognized from the practical operating data of a district freezing station. As seen from Table 1, the performances of pumps of the same type are close to each other. The pump system is composed of four pumps working in parallel mode including two big pumps (KP1020-3 457mm in impeller diameter) and two smaller pumps (KP1020-3 406mm in impeller diameter) manufactured by Grundfos. A determined pressure difference control is applied here where the designed pressure difference and flow rate demand for the entire system are 45m and 3234m3/h, respectively.Table 1
Pumps parameters.
PAR PUMP 1.KP1020-3 406 2.KP1020-3 406 3.KP1020-3 457 4.KP1020-3 457 a -6.07e-6 -5.36e-6 -8.54e-6 -8.41e-6 b -0.00822 -0.00949 -0.00048 -0.00046 c 60.41767 60.61789 74.54932 73.55836 j -4.08e-7 -3.96e-7 -3.87e-7 -3.82e-7 k 0.000940 0.000911 0.000999 0.000985 l 0.316357 0.328377 0.249898 0.259137


## 4.2. Simulation Result

The strategy space is discretized, and it becomes a limited strategy space to ensure that the potential game is a finite game. And, the speed ratio (SR) of each pump, as selected action set, is shown in Table2.Table 2
Pumps discreted strategy space.
SR PUMP 1.KP1020-3 406 2.KP1020-3 406 3.KP1020-3 457 4.KP1020-3 457 1 0 0 0 0 2 0.93 0.94 0.825 0.825 3 0.94 0.945 0.85 0.85 4 0.945 0.95 0.875 0.875 5 0.95 0.955 0.9 0.9 6 0.96 0.96 0.915 0.915 7 0.97 0.97 0.925 0.925 8 0.98 0.98 0.95 0.95 9 0.99 0.99 0.975 0.975 10 1 1 1 1 Note. The speed ratio of 0 means that the pump is off.In the case study, the proposed algorithms will be triggered as follows: when the pipe net has been changed, the smart unit will detect that the pressure differenceHo is not equal toHset andQo which can be got through the flow meter on the trunk; thusQi can be solved in (3).DEDA can be downloaded to microprocessor chip in each smart pumpi (i=1,2, 3,4), and then each smart pump can be active. And, each smart pump can communicate and cooperate with other smart pump nodes in serial mode. The corresponding results by the proposed DEDA are shown in Figure 5.Figure 5
Iterative process in decentralized optimization.As shown in Figure5, the probability distribution (PD) on each pump is dynamic during the iterative process in DEDA. And the optimal solution of each pump fluctuates and finally reaches a new equilibrium state under mutual negotiation and adjustment. By about 60 iteration steps, the pumps optimal load distribution task can be achieved. It can be seen that the probability of the optimal output flow rate tends to approach 1. Conversely, the probability of the other solution space approaches 0.The iterative process of each pumpi(i=1,2,3,4) cost is illustrated in Figure 6. Further, the iterative process of total energy cost is shown in Figure 7.Figure 6
Iterative process of cost on each pump in decentralized optimization.Figure 7
Iterative process of total cost in decentralized optimization.Comparing Figure6 and Figure 7, we find that although not all the costs of the four pumps decrease, the total cost decreases during the iterative process. This can validate the effectiveness of DEDA. For a comparison, the iterative process of EDA in centralized method is also illustrated in Figure 8. As shown in Figures 5 and 8, the DEDA can achieve the actual optima which are equivalent to that obtained with the centralized algorithm.Figure 8
Iterative process in centralized optimization.Further, a more detailed calculation result is shown in Table3, where EDA and DEDA are adopted, respectively.Table 3
Comparison between centralized and decentralized results.
#Working condition DEDA EDA Demand: 1496.9 m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2 Pump 3 Pump 4 47.000.00 41.25 0.00 700.170.00 796.71 0.00 0.790.00 0.85 0.00 47.000.00 41.25 0.00 700.170.00 796.71 0.00 0.790.00 0.85 0.00  Power/kW 227.4804 227.4804  Difference 0.0169 0.0169  Demand: 2583.4m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2Pump 3Pump 4 48.500.00 0.00 48.75 886.980.00 0.00 1695.44 0.830.00 0.00 0.82 48.500.00 0.00 48.75 886.980.00 0.00 1695.44 0.830.00 0.00 0.82  Power/kW 392.2045 392.2045  Difference 0.9772 0.9772  Demand: 2952.9m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2Pump 3Pump 4 47.2549.00 0.00 45.00 733.01926.41 0.00 1292.65 0.800.84 0 0.89 48.0047.50 40.00 45.75 827.24745.890.001379.76 0.820.8000.88  Power/kW 435.5667 439.1761  Difference 0.8373 0.0181  Demand: 3234.7m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2Pump 3Pump 4 50.0047.75 45.75 40.00 1054.53777.35 1402.09 0.00 0.850.81 0.87 0.00 48.5047.70 0.00 47.50 886.98777.34 0.00 1569.19 0.830.81 0.00 0.84  Power/kW 475.8191 485.6048  Difference 0.7344 1.1869  Demand: 3858.3m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2Pump 3Pump 4 47.0049.00 41.25 46.25 700.17926.41 796.71 1435.66 0.790.84 0.85 0.87 48.5050.00 41.25 43.75 886.981037.53 796.71 1136.63 0.830.85 0.85 0.89  Power/kW 572.8274 561.5738  Difference 0.6517 0.4537  Demand: 4901.2 m3/h n(Hz) Q (m 3 /h) η n(Hz) Q (m 3 /h) η  Pump 1Pump 2Pump 3Pump 4 49.0049.00 50.00 46.25 944.62926.41 1831.72 1435.66 0.840.84 0.78 0.87 48.0050.00 46.25 50.00 827.241037.51 1456.97 1815.93 0.820.85 0.86 0.79  Power/kW 778.3465 777.4647  Difference 0.4889 1.2384


## 5. Conclusion

Most of the optimizing methods for parallel pump systems in HVAC system are based on a centralized structure under which all the management tasks have to be settled with one master controller. In this structure, the computation efficiency and stability are relatively low in that all the data and inspection have to be processed by one controller and that one malfunctioning device will influence the whole system. Besides, the control strategy is a case-by-case one which cannot be transplanted easily to other systems.In this study, a decentralized control strategy has been applied to the parallel pump systems in which the traditional master controllers have been substituted by a list of independent smart units embedded in each pump.(i) A new fully distributed mathematical model for pumps optimal distribution and decentralized flat architecture are presented. With the rapid development of the electronic industry, smart hardware has been widely utilized in different fields. According to the vision of the decentralized method, traditional units can be upgraded to smart agents by incorporating a microcontroller chip inside. The decentralized algorithm can be written into the chip. The smart unit communicates with neighboring nodes such that they can work collaboratively to achieve the optimal load distribution task. In this case, the complicated onsite modeling, configuration, commissioning, and other secondary developing work can be simplified to the wiring of communication connection among different nodes because the smart unit can self-organize and plug-and-play.(ii) Decentralized optimization is used to solve pumps optimal load distribution under potential game framework. And, to achieve effectiveness of global goals, DEDA using local interaction is introduced. The convergence property of the novel method is also analyzed theoretically. Simulation results presented in the preceding paragraph indicate that the proposed method can achieve optimal load distribution.This scenario is extremely motivating for further investigation in the field of distributed networked control in other large-scale interconnected nonlinear systems. And, there are still several potential issues of the presented method that needed to be emphasized here. Firstly, the proposed mechanism is executed in sequence. Therefore, the proposed distributed algorithms should in principle be applicable to chain-like networks. Each smart pump is connected hand to hand with its neighbor, but the meshed network should also be considered. Secondly, future work will consider the hardware test and implementation on practical engineering.



---

*Source: 5468398-2018-07-09.xml*