# Annular Directed Distributed Algorithm for Energy Internet

Authors: Jiaming Tan; Hongyou Li
Journal: International Transactions on Electrical Energy Systems
 (2022)

---

## Abstract

This article investigates an annular directed distributed double optimal algorithm to manage many we-energy frameworks in energy management of energy Internet (EI). The we-energy (WE) is an integrated energy hub containing varied energy devices of different functions including multi-energy production, consumption, and conversion. On this basis, all WE models cooperate to search for a minimum value of an objective function. Energy management in EI has two main goals. On the one hand, it needs to attain the optimality of economy with influence about the fluctuation of distributed renewable energy and randomness of terminal users. On the other hand, the EI should protect the privacy of terminal users well. Besides, discovering optimality value in the oscillation near convergence point, EI also needs a decrease in communication frequency and refraining of Zeno behavior. Zeno behavior means some operation is triggered infinite times in finite times of iteration. For realizing these proposes, this literature establishes an EI system that transfers cyber information in an annular directed path. The algorithm in this EI system adopts a novel annular distributed double-control price guiding strategy. In addition, this algorithm employs other two methods including the alternating direction method of multipliers method and the Newton-downhill method to optimize economy and reach convergence, respectively. Meanwhile, that algorithm adopts a small positive constantw to avoid Zeno behavior. The performance of that algorithm is demonstrated through simulation results. Moreover, the optimality, convergence analysis, and avoiding Zeno behaviors are strictly proved by convex optimization and the monotone-bounded convergence theorem.


---

## Body

## 1. Introduction

Due to increasing concern about cosmopolitan environmental problems, the strategy of carbon emission reduction, and the utilization of environmentally friendly energy resources, the concept of energy internet (EI) [1] is developing at a marvelous speed in recent years. EI is a multi-energy system with a large number of advanced technology containing the theory of multiagent systems; the cyber communication and the physical energy transmission method [2]; the intelligent strategy of energy management; game and synergy theory about the multi-energy resource in a certain region; demand response between every energy units; optimal operation with various energy devices with different functions, among others. EI, a new pattern of energy system, devotes itself to satisfying terminal users’ various kinds of energy demand and enhancing the renewable energy’s utilization efficiency [3]. However, due to the fluctuation of renewable energy resources [4–6], the randomness of terminal users, and complex but strong coupling in various of energy, the optimal dispatch in EI is exceedingly difficult and the method about that is exceedingly unadvanced for complex multi-energy distributed systems and eagerly needed to be reinforced nowadays.We can divide the recent research concerning optimal dispatch of EI into two categories approximately: one is the traditionally centralized optimal strategy, and the other is the newly developing distributed optimal strategy. There are several documents about the centralized optimal strategy. The work in [7] solved the energy management with a non-convex object, and [8] introduced a two-stage multi-objective optimal scheduling in EI. The centralized optimal strategy has high quality in many fields including optimal performance, speed of convergence, and conquering disturbance caused by energy generators and terminal users. So traditional single energy systems including traditional power systems generated by terminal users always adopt the centralized optimal strategy. However, the centralized optimal strategy needs a large-scale energy generator that naturally disperses renewable energy resources that could not build. Centralized EI also brings too much compute and control pressure to its central processing unit. What’s more, the communication frequency in centralized EI is too high, and local malfunctions in the centralized optimal strategy influence global systems strikingly. Besides, it is impossible for the centralized optimal strategy to protect users’ privacy because of its concentrated physical framework. Therefore, the centralized optimal strategy will be replaced by the distributed optimal strategy sooner or later [9]. The literature [10, 11] proposed non-iterative algorithms for the distributed solution of multiagent optimal dispatch problems. The distributed optimal strategy can greatly decrease the frequency of communication and reduce the impact of local malfunction on global systems. The distributed optimal strategy is plug and play, so that it is very easy for distributed EI to expand its scale. And the distributed optimal strategy can protect users’ privacy to a certain extent in the present researches. Most important of all, it is impossible for natural disperse renewable energy resources to build large-scale centralized energy generators, so the only way to utilize renewable energy is the distributed optimal strategy. All advantages above are greatly significant to EI and can only be satisfied by the distributed optimal strategy and not by the centralized optimal strategy. For the foregoing reasons, it is no wonder that concern in the study of EI transfers from the centralized optimal strategy to the distributed optimal strategy in recent years.However, the difficulty of energy management in the distributed optimal strategy is much higher than that in the centralized optimal strategy on account of renewable energy resources’ fluctuation and terminal users’ randomness. Due to the long distance between energy producers and terminal users in the centralized optimal strategy, although it will strongly increase the cost of energy transmission, the fluctuation of renewable energy resources and the randomness of terminal users scarcely affect the EI system. But in the distributed optimal strategy, the distance between generators and terminal users may be very short, so if the fluctuation of renewable energy resources and the randomness of terminal users couple together, it will impact the stability of EI strongly. For solving the problem, various literatures use various strategies. The literature [6] revealed the disturbance caused by renewable energy resources’ fluctuation in power systems. The literature [12] presented an unsupervised algorithm to extract the EV charging loads non-intrusively from the smart meter data. The literature [13] optimized power trading by Stackelberg game. But methods below aimed at power systems, and the type of energy was only power. These methods cannot solve energy management of EI with strong coupling between different types of energy.Generally, we resolve a big complex system to some small and single systems and study small systems respectively to study the big system which is hard to study directly. So, we used agents to study multiagents [12], and we used microgrids to study power systems [14]. But what can we use to study EI? Different people had different ideas. Swiss scholars proposed the energy hub [15], which resolves EI into some small unit for the first time, and interpreted that energy can transform into another kind of energy in EI firstly. Subsequently, Li.L proposed prosumer in [16] to move forward a new step about establishing a small unit, which explained basic energy unit in EI is not only an energy producer and an energy consumer. The prosumer makes much progress not only in coupling between generators and terminal users but also in cooperating among different energies [17]. Sun et al. summed up all models above and referenced the theory of multiagents [18] and proposed we-energy (WE) [19]. We-energy is a full-duplex, hole distributed without a center, intelligent and peer-to-peer energy unit in EI. Compared with other models, WE is most suitable for the current EI for the several reasons below. Firstly, the WE is a full-duplex model, while other models including energy hubs are half-duplex models. Secondly, the WE is equal to the hole network. In other words, the WE is completely a selfish and rational model. So, the decision of the WE is very hard to be influenced by other we-energies, which is very suitable for the distributed algorithms. Thirdly, the WE is a point-to-point model, but other models are point-to-plane models. The WE communicates to other Wes, while others communicate to the hole network. So, the WE is more suitable for the algorithm in this article. Because of these reasons, this study adopts the WE as a basic energy unit.Since the distributed framework and the large-scale system of EI, the suitable method of energy management between WE is the distributed algorithm realized by multiagent. The literature [20] used the distributed algorithm to optimize residential WE and designed an operation method that maps the infeasible solutions into the feasible region. The literature [21] mixed the alternating direction method of multipliers into the distributed algorithm for the first time. The literature [22] converted synchronous communication into out of synchronous communication by event-triggered in the distributed algorithm. But the distributed algorithm in [21, 23] considered energy conversions as must-run energy load and flexible energy load, which made energy conversion could only supply all load in one type of energy demand. For terminal users, the energy that the producer generated and conversion devices transformed made no difference, so EI by no means differentiated load into must-run load and flexible load. In addition, the literature [21, 23] ignored Zeno behavior in the algorithm, which could increase iteration times closed to infinity value under some initial values.Although the distributed optimal strategy and the distributed algorithm reduce communication frequency and protect users’ privacy in EI to a certain degree, the communication frequency’s reduction and privacy’s protection still need to be reinforced. Nowadays, the energy that terminal users need, generators produce, and energy transformers convert is becoming more and more random, varied, and unstable. Furthermore, time delay [22], controllability deficiency [24], and Zeno behavior that always consists in the Newton–Raphson algorithm [23] do serious harm to EI. What is most terrible is Zeno behavior, because it will increase times of iteration to a large value, even to a measureless value. Moreover, traditional researches about privacy protecting only protected privacy from unneighborly agents, but they did nothing about neighbor agents [25–27]. That method had disadvantages. On the one hand, network sparse caused by only-neighbor communication destroyed the controllability of EI due to the theory in [24]. On the other hand, some users were not willing to trust their neighbor agents. To solve these problems, this literature based on the alternating direction method of multipliers in EI [28], Newton–Raphson algorithm [23], and auction-based algorithm [29] proposes an annular directed distributed algorithm, which leads WEs to communicate in an annular route. Compared with traditional distributed optimization methods, the annular directed distributed algorithm can not only observably reduce communicating intensity but also protect all users’ data that cannot be acquired by other agents including neighbor agents. However, in conventional distributed optimization methods, some ratios, power output, or estimated prices have to be shared among neighbored agents so they cannot protect all data. In addition, this study establishes a small number to avoid Zeno behaviors about Newton-downhill factors which traditional researches do not consider. The contributions of this study are summarized as follows:(1)
Reinforce the protection of users’ privacy. In the annular directed distributed algorithm, devices’ data can only be extracted by themselves. By the annular directed distributed algorithm, all we-energies only transfer a power-heat-gas energy flow to another we-energy. Besides, all devices inside the we-energy only transfer that power-heat-gas flow, too. The energy flow is the summation of previous energy manufacture, energy conversion, and encryption factors. Firstly, each energy device does not know where the energy flow begins (the start place is randomness and that of each time of iteration is different), so the energy device does not know the operating condition of each device. Secondly, the privacy of the started device is protected by encryption factors. When the energy flow starts, it will add several values to each type of energy. With this method, the second device cannot know the operating condition of the first device. However, traditional research cannot protect several types of data among neighbor agents [25–27]. In traditional researches, some ratios, power output, or estimated prices have to be shared among neighbored agents so they cannot protect all data.(2)
Reduce the intensity of communication and change the information in communication from last time data to updated data. Firstly, the communication in the annular directed distributed algorithm is single-directed communicating and the communication in traditional researches [27] is non-directed. The non-directed communication means the EI needs to communicate in both ways. So, the intensity of single-directed communication is half of the non-directed communication. Secondly, communication in traditional research is one-to-many communication, while the communication in this article is one-to-one communication. There are several contributions to one-to-one communication. Firstly, the encryption of energy flow described hereinbefore needs that. Secondly, the data one-to-many communication transfers are the data in last time all above, but the data one-to-one communication transfers are the data that have already been updated by the new energy price in this time of iteration, so the information in one-to-one communication is newer than that in one-to-many communication.(3)
Avoid abnormal energy conversion successfully. It is worth noting that the price guiding the alternating direction method of multipliers is unsuitable for the energy conversion because the price guiding always leads to the overshoot in energy conversion. Firstly, the cost functions of energy conversion devices are always linear but the cost functions of energy manufacture devices are always convex so too much energy will increase the energy-producing price but will not change the energy conversion price. Secondly, because the energy sold price need to be changed in the control center of the hole EI, a distributed we-energy cannot change the energy price even though too much energy changes to another type of energy. For the reasons, above one type of cheap energy will change all of themselves to another expensive energy even though the energy conversion may consume all this type of energy and produce too many other energies. For avoiding that, this study enlightened by the idea of mistest in reinforcement learning. The literature [30] proposed a mistesting method to solve the problem of energy conversion and use downhill factors to ensure the astringency of the algorithm.(4)
Avoid Zeno behavior in the Newton-downhill method successfully. Zeno behavior means some operation is triggered infinite times in finite times of iteration. In this literature, Zero behavior means Newton-downhill factors are adjusted too many times in one time of iteration. Theoretically speaking, the Newton-downhill method will go converge sooner or later, so the Zeno behavior will never appear. However, in the realistic project, if the downhill factor changes too many times in one time of iteration, we will regard it as Zeno behavior. This article solved that issue.(5)
Prove that the annular directed distributed algorithm holds asymptotic convergence when the Zeno coefficient reaches maximum value by the monotone-bounded convergence theorem. Meanwhile, optimal performance of the equilibrium point is certified by difference theorem and convex optimization. And by the theory of finite and infinite, avoiding Zeno behavior in certain accuracy is proved, too.The rest of this study is organized as follows: Section2 presents the WE model, energy device models, constraints of main features, and requirements about energy devices. Section 3 introduces some basic knowledge about graph theory; presents the annular directed distributed algorithm for the first time; and proves astringency, optimality, and avoiding Zeno behaviors of it. In Section 4, several simulation results are presented to prove the effectiveness of that. Section 5 concludes the study.


## 2. The We and Device Models in EI

### 2.1. System Modeling and Test Systems

The distributed energy producer includes distributed renewable generators, distributed renewable heat devices, distributed coal or oil combined heat and power devices, and equivalent distributed gas producers. There are two things worth noting: one is that there are no distributed fuel generators and distributed fuel heat devices because they are both more inefficiently than distributed combined heat and power devices. The other is that there are two kinds of combined heat and power devices in EI: one is the distributed coal or oil combined heat and power devices belonging to the distributed energy producer, which produces power and heat by coal or oil, and the other is the distributed gas combined heat and power devices belonging to the distributed energy conversion devices, which produces power and heat by gas. The distributed energy conversion devices include distributed power to gas devices, distributed electric boilers, and distributed gas combined heat and power devices. The distributed energy storage devices include distributed power storage devices, distributed heat storage devices, and distributed gas storage devices. The distributed terminal users and the distributed energy transform devices cannot be divided into smaller devices.The distributed energy producer satisfies the constraints below:(1)PiDP=PiDRG+PiDCOC,HiDP=HiDRHD+HiDCOC,GiDP=GiDGP,where i is the serial number of WE. If there are n WEs in EI, i belongs to 1 to n. PiDP, HiDP, and GiDP are the total power, heat rate, and gas rate produced by the distributed energy producer devices in ith WE, respectively. PiDRG and PiDCOC are the power rate of the distributed renewable generators and the distributed coal or oil combined heat and power devices in ith WE, respectively. HiDRHD andHiDCOC are the heat rate of DRHD and the distributed coal or oil combined heat and power devices in ith WE. GiDGP is the gas rate of the distributed gas producers in ith WE.The distributed energy conversion devices satisfy the constraints below:(2)PiDCD=PiDCD−PiDP2G−PiDEB,HiDCD=HiDEB+HiDGC,GiDCD=GiDP2G−GiDGC,where PiDCD, HiDCD, and GiDCD are the total power, heat rate, and gas rate exchanged by the distributed energy conversion devices in ith WE, respectively—plus or minus of them represents the output or input. PiDGC, PiDP2G, andPiDEB are the exchanging power rate in the distributed gas combined heat and power devices, the distributed power to gas devices, and the distributed electric boiler in ith WE, respectively. HiDEB and HiDGC are the exchanging heat rate in the distributed electric boiler and the distributed gas combined heat and power devices in i WE, respectively. GiDCD, GiDP2G, and GiDGC are the exchanging gas rate in the distributed energy conversion devices, the distributed power to gas devices, and the distributed gas combined heat and power devices in ith WE, respectively.


### 2.2. WE Model

WE is a basic energy unit that can have devices all above. The model of WE in power-gas-heat EI is as follows:(3)PiHiGi=λi,P,Pλi,P,Hλi,P,Gλi,H,Pλi,H,Hλi,H,Gλi,G,Pλi,G,Hλi,G,GPi′Hi′Gi′,where P, H, andG are power, heat, and gas, respectively. Vector PiHiGiT and Pi′Hi′Gi′T are energy flow in terminal side and network side—plus or minus and an absolute value of the element in vector represent the direction and rate of energy flow, respectively. Besides, the following constraints are limited to micro-power-system in WE:(4)Pi′=PiDP+PiDTD,Pi−Pi′=PiDSD,where PiDTD is the power flow between ith WE and others—plus or minus and the absolute value of it are the direction and rate of power flow, respectively. The following constraints in micro-heat-system and micro-gas-system in WE are similar to constraints above. All P are changed into H or G correspondingly. Matrix λ represents energy flow’s proportion of allocation and efficiency of energy conversion. How to design matrix λ? Different research has different points of view in it. Ref. [8] models matrix λ as the dot product of Hadamard matrices and efficiency matrices. But Ref. [8] prohibits circumflex in WE. The circumflex is harmful to smart grid, but may not be harmful to the distributed energy conversion devices in EI because the distributed energy conversion devices could never shorten out. Although circumflex will waste certain energy due to loss in energy conversion, because of a lot of start-stop constraints and ramping rate limits in the distributed energy conversion devices, circumflex is needed in energy management. So, this study models matrix λ as follows:(5)λ=K×A+K′×A′,where A′ and A are Hadamard matrices. An element in A is 1 or 0 represent that there is or is not corresponding energy flowing from network side to terminal side. Oppositely, that in A′ is 1 or 0 represents that there is or is not corresponding energy flowing from terminal side to network side. Elements in K and K′ are corresponding energy flow’s proportion of allocation and efficiency of corresponding energy conversion. For example, if an element in K is κa,b,(6)κa,b=va,bηa,b.There will beva,b proportion of energy a from network side converting into energy b in terminal side, and conversion efficiency is ηa,b. So, κa,b is a positive number smaller than 1. Oppositely, if an element in K′ is κa,b′,(7)κa,b′=1va,b′ηa,b′.There will beva,b′ proportion of energy a from terminal side converting into energy b in network side, and conversion efficiency is ηa,b′. So, κa,b′ is a number larger than 1. It is worth nothing that, if va,b′ equals to zero, κa,b′ will be meaningless. On this occasion, we stipulate va,b′ηa,b′ as a small positive number ζ in order to run the distributed algorithm. Constraints in A′ ensure that ζ could not influence precision of results. So, in power-gas-heat EI, the WE matrix is as follows:(8)PHG=A′∘1vPP′SSTηPP′SST1vPH′DEBηPH′DEB1vPG′DP2GηPG′DP2G1ξ11ξ1vGP′DGCηGP′DGC1vGH′DGCηGH′DGC1vGG′EWC+A∘vPPSSTηPPSST0vGPDGCηGPDGCvPHDEBηPHDEB1vGHDGCηGHDGCvPGDP2GηPGDP2G0vGGEWCP′H′G′,where SST is a solid-state transformer and EWC is corresponding energy without conversion. Energy flow in the network side satisfies the constraints below:(9)∑i=1nPi′=0,∑i=1nHi′=0,∑i=1nGi′=0,where n is the total of WEs.


### 2.3. The Distributed Energy Producer Device Models and Cost Functions

#### 2.3.1. The Distributed Renewable Generators and the Distributed Renewable Heat Device Models

The distributed renewable generators include wind-driven generators and solar-driven generators. But DRHD only includes solar-driven heat devices. Wind and solar are free, and the operating cost of taking the advantage of renewable energy resources is such little that it can be ignored. Models and cost functions of the distributed renewable generators and DRHD are as follows:(10)Pi,tDRG=Pi,twind+Pi,tsolar,Hi,tDRHD=Hi,tsolar,Ci,tDRG=0,Ci,tDRHD=0,where wind and solar are kinds of renewable energy resources. C is the operating cost function of corresponding device. t is hours from 0 to24. P and H are the rate of power and heat. It is worth noting that energy management in this study is hourly dispatch, so all-time in this study is in hour units.Constraints of the distributed renewable generators and DRHD are as below:(11)Piwind−min≤Pi,twind≤Piwind−max,Pisolar−min≤Pi,tsolar≤Pisolar−max,Hisolar−min≤Hi,tsolar≤Hisolar−max,where min and max are corresponding minimum rate and corresponding maximum rate of corresponding energy in corresponding WE.


#### 2.3.2. The Distributed Coal or Oil Combined Heat and Power Device Models

The distributed coal or oil combined heat and power devices adopt coal fuel to produce power and heat. Models and cost functions of the distributed renewable generators and DRHD are as follows:(12)Pi,tDCOC=ηiDCOCFi,tDCOC1φi,tDCOC+1,Hi,tDCOC=ηiDCOCFi,tDCOCφi,tDCOCφi,tDCOC+1,where Pi,tDCOC and Hi,tDCOC are power and heat rate of the distributed coal or oil combined heat and power devices generating, respectively, and Fi,tDCOC is the thermal rate of coal the distributed coal or oil combined heat and power devices consuming. ηiDCOC is the energy conversion of the distributed coal or oil combined heat and power devices. φi,tDCOC is the ratio of heat and power. It is worth noting that, coal does not belong to power-heat-gas network of EI. So, energy management cannot influence coal’s price and the production of coal. People cannot produce coal or oil, after all.(13)φi,tDCOC=PrhPrp,where Prh and Prp are the price of heat and power obviously. The operating cost function of the distributed coal or oil combined heat and power devices is(14)Ci,tDCOC=aiDCOCPiDCOC2+biDCOCPiDCOC+αiDCOCHiDCOC2+βiDCOCHiDCOC+ciDCOCPiDCOCHiDCOC+χiDCOC,where aiDCOCbiDCOCαiDCOC, βiDCOC, ciDCOC, and χiDCOC are constants, and aiDCOCαiDCOC, and χiDCOC are positive constants. The distributed coal or oil combined heat and power devices satisfy the constraints below:(15)−PiDCOC−ramp≤Pi,tDCOC−Pi,t−1DCOC≤PiDCOC−ramp,diDCOCPi,tDCOC+eiDCOCHi,tDCOC+fiDCOC≥0,diDCOCPi,tDCOC+eiDCOCHi,tDCOC≤giDCOC,hiDCOC≤φi,tDCOC≤jiDCOC,where PiDCOC−ramp is the power ramping constrain of the distributed coal or oil combined heat and power devices. diDCOC, eiDCOC, fiDCOC, giDCOC, hiDCOC, and jiDCOC are constants.


#### 2.3.3. The Distributed Gas Producer Models

The operating cost function of the distributed gas producers is(16)Ci,tDGP=ai,tDGPGi,tDGP2+bi,tDGPGi,tDGP+ci,tDGP,where bi,tDGP, and ci,tDGP are constants and ai,tDGP is a positive constant.G represents the rate of gas. The limit of the distributed gas producers is(17)0≤Gi,tDGP−min≤Gi,tDGP≤Gi,tDGP−max,where Gi,tDGP−min and Gi,tDGP−max are lower and higher limits of the distributed gas producers, respectively.


### 2.4. The Distributed Energy Conversion Device Models and Cost Functions

#### 2.4.1. The Distributed Power to Gas Devices and the Distributed Electric Boiler Models

The models and cost functions of the distributed power to gas devices are as follows:(18)Gi,tDP2G=ηi,tDP2GPi,tDP2G,Ci,tDP2G=θi,tDP2GPi,tDP2G,where ηi,tDP2G is the energy conversion efficiency of the distributed power to gas devices. θi,tDP2G is a positive constant. The limit of the distributed power to gas devices is(19)Pi,tDP2G−min≤Pi,tDP2G≤1n×Pnetwork,where Pi,tDP2G−min is the start-stop constraint of the distributed power to gas devices, Pnetwork is the surplus power rate in the network of EI. The models, operating cost functions, and limits of the distributed electric boiler are similar to that of the distributed power to gas devices. We only need to replace gas with heat.


#### 2.4.2. The Distributed Gas Combined Heat and Power Device Models

The model and operating cost functions of the distributed gas combined heat and power devices are as follows:(20)Pi,tDGC=ηiDGCGi,tDGC1φi,tDGC+1,Hi,tDGC=ηiDGCGi,tDGCφi,tDGCφi,tDGC+1,φi,tDGC=PrhPrp,Ci,tDGC=θi,tDGCGi,tDGC,where Pi,tDGC and Hi,tDGC are the power and heat rates of the distributed gas combined heat and power devices generation, Gi,tDGC is the gas rate of coal of the distributed gas combined heat and power devices consumption, ηiDGC is the energy conversion of the distributed gas combined heat and power devices, and φi,tDGC is the ratio of heat and power. The constraints of the distributed gas combined heat and power devices is: Gi,tDGC−min≤Gi,tDGC≤1/n×Gnetwork, which is similar to that in the distributed power to gas devices. The constraints of φi,tDGC is similar to that in the distributed coal or oil combined heat and power devices.


### 2.5. The Distributed Energy Storage Devices Models and Cost Functions

#### 2.5.1. The Distributed Energy, the Distributed Power Store Device, the Distributed Heat Storage Devices, and the Distributed Gas Storage Device Models

The distributed power storage devices have an optimal condition about reserve of stored power. If the stored power is less than optimal condition too much, it will harm batteries and other devices. If stored power is more than optimal condition too much, the stored power will loss too much. So, the optimal performance of the distributed power store device is(21)Oi,tDPSD=aiDPSDPi,tS−DPSDPi,tS−DPSD−2μiDPSD+biDPSD,where O is the optimal function of the reserved power in the distributed power store device, Pi,tS−DPSD is the stored power in time t, μiDPSD is the optimal reserve power in the distributed power store device, aiDPSD and biDPSD are constants, and aiDPSD is a negative constant. The cost function of the distributed power store device is(22)Ci,tDPSD=Oi,t−1DPSD−Oi,tDPSD+θi,tDPSDPi,tDPSD2,where θi,tDPSD is a positive constant and Pi,tDPSD is the power rate of the distributed power store device. If Pi,tDPSD is positive, it expresses power output, vice versa. So, we can know(23)Pi,tDPSD=Pi,t−1S−DPSD−Pi,tS−DPSD.The distributed power store device needs to satisfy the limits below:(24)−Piin−DPSD≤Pi,tDPSD≤Piout−DPSD,Pimin−S−DPSD≤Pi,tS−DPSD≤Pimax−S−DPSD,where Piin−DPSD and Piout−DPSD are the maximum limits of power input and output rate of the distributed power store device, and Pimin−S−DPSD andPimax−S−DPSD are capacity minimum and maximum values, respectively.The models, operating cost functions, and limits of the distributed heat storage devices and the distributed gas storage devices are similar to that of the distributed power store device. We only need to replace power with heat or gas.The energy loss in the distributed energy transform devices is very loss, for convenience, so we ignored that.


### 2.6. Energy Load Models

The must-run load and controllable load which traditional research including Ref. [20, 21] adopts may be unsuitable for EI because of several reasons that I proposed in the preceding part of this study. So, we only adopt Pi,tELHi,tEL, and Gi,tEL to represent power heat and gas loads. If EI supplies less energy to WEs than which they need, the terminal users’ energy demand will not be satisfied. But if EI supplies more energy to WEs than which they need, it will waste a lot of energy. So, the operating cost of EL depends on the unbalance between energy supply and demand. The cost function of power EL is as follows:(25)CtP−EL=θi,tP−ELPi,tsupply−Pi,tEL22,θi,tP−EL=θ1,tP−EL,Pi,tsupply≤Pi,tEL,θ2,tP−EL,Pi,tsupply≥Pi,tEL,where θ1,tP−EL and θ2,tP−EL are two positive constants. The reason for dividing θtP−EL into θ1,tP−EL and θ2,tP−EL is that the harm of surplus and lack of power is different. The cost functions of heat and gas are similar to that of power. The only difference is the type of energy.


### 2.7. Object Functions

This study focuses on the economical optimization by energy management in WEs and EI. So, the object function is used to co-planning all WEs to realize the maximize of the economy in (22):(26)maxF=∑i=1nWi,t,(27)Wi,t=∑j=1mPrp×ΔPij+Prh×ΔHij+Prg×ΔGij−Prf×Fij−Cij,where Wi,t is the economy function of each WE; j is the serial number of each devices in WE; Prg and Prf are price of gas and coal; ΔPij, ΔHij, andΔGij are energy flow rate of power, heat, and gas, respectively. If they are positive, it means that the devices output the corresponding energy, vice versa. Of course, they can be zero. Fij is consumption of coal, which is non-negative. Cij is the operating cost function of devices.


## 2.1. System Modeling and Test Systems

The distributed energy producer includes distributed renewable generators, distributed renewable heat devices, distributed coal or oil combined heat and power devices, and equivalent distributed gas producers. There are two things worth noting: one is that there are no distributed fuel generators and distributed fuel heat devices because they are both more inefficiently than distributed combined heat and power devices. The other is that there are two kinds of combined heat and power devices in EI: one is the distributed coal or oil combined heat and power devices belonging to the distributed energy producer, which produces power and heat by coal or oil, and the other is the distributed gas combined heat and power devices belonging to the distributed energy conversion devices, which produces power and heat by gas. The distributed energy conversion devices include distributed power to gas devices, distributed electric boilers, and distributed gas combined heat and power devices. The distributed energy storage devices include distributed power storage devices, distributed heat storage devices, and distributed gas storage devices. The distributed terminal users and the distributed energy transform devices cannot be divided into smaller devices.The distributed energy producer satisfies the constraints below:(1)PiDP=PiDRG+PiDCOC,HiDP=HiDRHD+HiDCOC,GiDP=GiDGP,where i is the serial number of WE. If there are n WEs in EI, i belongs to 1 to n. PiDP, HiDP, and GiDP are the total power, heat rate, and gas rate produced by the distributed energy producer devices in ith WE, respectively. PiDRG and PiDCOC are the power rate of the distributed renewable generators and the distributed coal or oil combined heat and power devices in ith WE, respectively. HiDRHD andHiDCOC are the heat rate of DRHD and the distributed coal or oil combined heat and power devices in ith WE. GiDGP is the gas rate of the distributed gas producers in ith WE.The distributed energy conversion devices satisfy the constraints below:(2)PiDCD=PiDCD−PiDP2G−PiDEB,HiDCD=HiDEB+HiDGC,GiDCD=GiDP2G−GiDGC,where PiDCD, HiDCD, and GiDCD are the total power, heat rate, and gas rate exchanged by the distributed energy conversion devices in ith WE, respectively—plus or minus of them represents the output or input. PiDGC, PiDP2G, andPiDEB are the exchanging power rate in the distributed gas combined heat and power devices, the distributed power to gas devices, and the distributed electric boiler in ith WE, respectively. HiDEB and HiDGC are the exchanging heat rate in the distributed electric boiler and the distributed gas combined heat and power devices in i WE, respectively. GiDCD, GiDP2G, and GiDGC are the exchanging gas rate in the distributed energy conversion devices, the distributed power to gas devices, and the distributed gas combined heat and power devices in ith WE, respectively.


## 2.2. WE Model

WE is a basic energy unit that can have devices all above. The model of WE in power-gas-heat EI is as follows:(3)PiHiGi=λi,P,Pλi,P,Hλi,P,Gλi,H,Pλi,H,Hλi,H,Gλi,G,Pλi,G,Hλi,G,GPi′Hi′Gi′,where P, H, andG are power, heat, and gas, respectively. Vector PiHiGiT and Pi′Hi′Gi′T are energy flow in terminal side and network side—plus or minus and an absolute value of the element in vector represent the direction and rate of energy flow, respectively. Besides, the following constraints are limited to micro-power-system in WE:(4)Pi′=PiDP+PiDTD,Pi−Pi′=PiDSD,where PiDTD is the power flow between ith WE and others—plus or minus and the absolute value of it are the direction and rate of power flow, respectively. The following constraints in micro-heat-system and micro-gas-system in WE are similar to constraints above. All P are changed into H or G correspondingly. Matrix λ represents energy flow’s proportion of allocation and efficiency of energy conversion. How to design matrix λ? Different research has different points of view in it. Ref. [8] models matrix λ as the dot product of Hadamard matrices and efficiency matrices. But Ref. [8] prohibits circumflex in WE. The circumflex is harmful to smart grid, but may not be harmful to the distributed energy conversion devices in EI because the distributed energy conversion devices could never shorten out. Although circumflex will waste certain energy due to loss in energy conversion, because of a lot of start-stop constraints and ramping rate limits in the distributed energy conversion devices, circumflex is needed in energy management. So, this study models matrix λ as follows:(5)λ=K×A+K′×A′,where A′ and A are Hadamard matrices. An element in A is 1 or 0 represent that there is or is not corresponding energy flowing from network side to terminal side. Oppositely, that in A′ is 1 or 0 represents that there is or is not corresponding energy flowing from terminal side to network side. Elements in K and K′ are corresponding energy flow’s proportion of allocation and efficiency of corresponding energy conversion. For example, if an element in K is κa,b,(6)κa,b=va,bηa,b.There will beva,b proportion of energy a from network side converting into energy b in terminal side, and conversion efficiency is ηa,b. So, κa,b is a positive number smaller than 1. Oppositely, if an element in K′ is κa,b′,(7)κa,b′=1va,b′ηa,b′.There will beva,b′ proportion of energy a from terminal side converting into energy b in network side, and conversion efficiency is ηa,b′. So, κa,b′ is a number larger than 1. It is worth nothing that, if va,b′ equals to zero, κa,b′ will be meaningless. On this occasion, we stipulate va,b′ηa,b′ as a small positive number ζ in order to run the distributed algorithm. Constraints in A′ ensure that ζ could not influence precision of results. So, in power-gas-heat EI, the WE matrix is as follows:(8)PHG=A′∘1vPP′SSTηPP′SST1vPH′DEBηPH′DEB1vPG′DP2GηPG′DP2G1ξ11ξ1vGP′DGCηGP′DGC1vGH′DGCηGH′DGC1vGG′EWC+A∘vPPSSTηPPSST0vGPDGCηGPDGCvPHDEBηPHDEB1vGHDGCηGHDGCvPGDP2GηPGDP2G0vGGEWCP′H′G′,where SST is a solid-state transformer and EWC is corresponding energy without conversion. Energy flow in the network side satisfies the constraints below:(9)∑i=1nPi′=0,∑i=1nHi′=0,∑i=1nGi′=0,where n is the total of WEs.


## 2.3. The Distributed Energy Producer Device Models and Cost Functions

### 2.3.1. The Distributed Renewable Generators and the Distributed Renewable Heat Device Models

The distributed renewable generators include wind-driven generators and solar-driven generators. But DRHD only includes solar-driven heat devices. Wind and solar are free, and the operating cost of taking the advantage of renewable energy resources is such little that it can be ignored. Models and cost functions of the distributed renewable generators and DRHD are as follows:(10)Pi,tDRG=Pi,twind+Pi,tsolar,Hi,tDRHD=Hi,tsolar,Ci,tDRG=0,Ci,tDRHD=0,where wind and solar are kinds of renewable energy resources. C is the operating cost function of corresponding device. t is hours from 0 to24. P and H are the rate of power and heat. It is worth noting that energy management in this study is hourly dispatch, so all-time in this study is in hour units.Constraints of the distributed renewable generators and DRHD are as below:(11)Piwind−min≤Pi,twind≤Piwind−max,Pisolar−min≤Pi,tsolar≤Pisolar−max,Hisolar−min≤Hi,tsolar≤Hisolar−max,where min and max are corresponding minimum rate and corresponding maximum rate of corresponding energy in corresponding WE.


### 2.3.2. The Distributed Coal or Oil Combined Heat and Power Device Models

The distributed coal or oil combined heat and power devices adopt coal fuel to produce power and heat. Models and cost functions of the distributed renewable generators and DRHD are as follows:(12)Pi,tDCOC=ηiDCOCFi,tDCOC1φi,tDCOC+1,Hi,tDCOC=ηiDCOCFi,tDCOCφi,tDCOCφi,tDCOC+1,where Pi,tDCOC and Hi,tDCOC are power and heat rate of the distributed coal or oil combined heat and power devices generating, respectively, and Fi,tDCOC is the thermal rate of coal the distributed coal or oil combined heat and power devices consuming. ηiDCOC is the energy conversion of the distributed coal or oil combined heat and power devices. φi,tDCOC is the ratio of heat and power. It is worth noting that, coal does not belong to power-heat-gas network of EI. So, energy management cannot influence coal’s price and the production of coal. People cannot produce coal or oil, after all.(13)φi,tDCOC=PrhPrp,where Prh and Prp are the price of heat and power obviously. The operating cost function of the distributed coal or oil combined heat and power devices is(14)Ci,tDCOC=aiDCOCPiDCOC2+biDCOCPiDCOC+αiDCOCHiDCOC2+βiDCOCHiDCOC+ciDCOCPiDCOCHiDCOC+χiDCOC,where aiDCOCbiDCOCαiDCOC, βiDCOC, ciDCOC, and χiDCOC are constants, and aiDCOCαiDCOC, and χiDCOC are positive constants. The distributed coal or oil combined heat and power devices satisfy the constraints below:(15)−PiDCOC−ramp≤Pi,tDCOC−Pi,t−1DCOC≤PiDCOC−ramp,diDCOCPi,tDCOC+eiDCOCHi,tDCOC+fiDCOC≥0,diDCOCPi,tDCOC+eiDCOCHi,tDCOC≤giDCOC,hiDCOC≤φi,tDCOC≤jiDCOC,where PiDCOC−ramp is the power ramping constrain of the distributed coal or oil combined heat and power devices. diDCOC, eiDCOC, fiDCOC, giDCOC, hiDCOC, and jiDCOC are constants.


### 2.3.3. The Distributed Gas Producer Models

The operating cost function of the distributed gas producers is(16)Ci,tDGP=ai,tDGPGi,tDGP2+bi,tDGPGi,tDGP+ci,tDGP,where bi,tDGP, and ci,tDGP are constants and ai,tDGP is a positive constant.G represents the rate of gas. The limit of the distributed gas producers is(17)0≤Gi,tDGP−min≤Gi,tDGP≤Gi,tDGP−max,where Gi,tDGP−min and Gi,tDGP−max are lower and higher limits of the distributed gas producers, respectively.


## 2.3.1. The Distributed Renewable Generators and the Distributed Renewable Heat Device Models

The distributed renewable generators include wind-driven generators and solar-driven generators. But DRHD only includes solar-driven heat devices. Wind and solar are free, and the operating cost of taking the advantage of renewable energy resources is such little that it can be ignored. Models and cost functions of the distributed renewable generators and DRHD are as follows:(10)Pi,tDRG=Pi,twind+Pi,tsolar,Hi,tDRHD=Hi,tsolar,Ci,tDRG=0,Ci,tDRHD=0,where wind and solar are kinds of renewable energy resources. C is the operating cost function of corresponding device. t is hours from 0 to24. P and H are the rate of power and heat. It is worth noting that energy management in this study is hourly dispatch, so all-time in this study is in hour units.Constraints of the distributed renewable generators and DRHD are as below:(11)Piwind−min≤Pi,twind≤Piwind−max,Pisolar−min≤Pi,tsolar≤Pisolar−max,Hisolar−min≤Hi,tsolar≤Hisolar−max,where min and max are corresponding minimum rate and corresponding maximum rate of corresponding energy in corresponding WE.


## 2.3.2. The Distributed Coal or Oil Combined Heat and Power Device Models

The distributed coal or oil combined heat and power devices adopt coal fuel to produce power and heat. Models and cost functions of the distributed renewable generators and DRHD are as follows:(12)Pi,tDCOC=ηiDCOCFi,tDCOC1φi,tDCOC+1,Hi,tDCOC=ηiDCOCFi,tDCOCφi,tDCOCφi,tDCOC+1,where Pi,tDCOC and Hi,tDCOC are power and heat rate of the distributed coal or oil combined heat and power devices generating, respectively, and Fi,tDCOC is the thermal rate of coal the distributed coal or oil combined heat and power devices consuming. ηiDCOC is the energy conversion of the distributed coal or oil combined heat and power devices. φi,tDCOC is the ratio of heat and power. It is worth noting that, coal does not belong to power-heat-gas network of EI. So, energy management cannot influence coal’s price and the production of coal. People cannot produce coal or oil, after all.(13)φi,tDCOC=PrhPrp,where Prh and Prp are the price of heat and power obviously. The operating cost function of the distributed coal or oil combined heat and power devices is(14)Ci,tDCOC=aiDCOCPiDCOC2+biDCOCPiDCOC+αiDCOCHiDCOC2+βiDCOCHiDCOC+ciDCOCPiDCOCHiDCOC+χiDCOC,where aiDCOCbiDCOCαiDCOC, βiDCOC, ciDCOC, and χiDCOC are constants, and aiDCOCαiDCOC, and χiDCOC are positive constants. The distributed coal or oil combined heat and power devices satisfy the constraints below:(15)−PiDCOC−ramp≤Pi,tDCOC−Pi,t−1DCOC≤PiDCOC−ramp,diDCOCPi,tDCOC+eiDCOCHi,tDCOC+fiDCOC≥0,diDCOCPi,tDCOC+eiDCOCHi,tDCOC≤giDCOC,hiDCOC≤φi,tDCOC≤jiDCOC,where PiDCOC−ramp is the power ramping constrain of the distributed coal or oil combined heat and power devices. diDCOC, eiDCOC, fiDCOC, giDCOC, hiDCOC, and jiDCOC are constants.


## 2.3.3. The Distributed Gas Producer Models

The operating cost function of the distributed gas producers is(16)Ci,tDGP=ai,tDGPGi,tDGP2+bi,tDGPGi,tDGP+ci,tDGP,where bi,tDGP, and ci,tDGP are constants and ai,tDGP is a positive constant.G represents the rate of gas. The limit of the distributed gas producers is(17)0≤Gi,tDGP−min≤Gi,tDGP≤Gi,tDGP−max,where Gi,tDGP−min and Gi,tDGP−max are lower and higher limits of the distributed gas producers, respectively.


## 2.4. The Distributed Energy Conversion Device Models and Cost Functions

### 2.4.1. The Distributed Power to Gas Devices and the Distributed Electric Boiler Models

The models and cost functions of the distributed power to gas devices are as follows:(18)Gi,tDP2G=ηi,tDP2GPi,tDP2G,Ci,tDP2G=θi,tDP2GPi,tDP2G,where ηi,tDP2G is the energy conversion efficiency of the distributed power to gas devices. θi,tDP2G is a positive constant. The limit of the distributed power to gas devices is(19)Pi,tDP2G−min≤Pi,tDP2G≤1n×Pnetwork,where Pi,tDP2G−min is the start-stop constraint of the distributed power to gas devices, Pnetwork is the surplus power rate in the network of EI. The models, operating cost functions, and limits of the distributed electric boiler are similar to that of the distributed power to gas devices. We only need to replace gas with heat.


### 2.4.2. The Distributed Gas Combined Heat and Power Device Models

The model and operating cost functions of the distributed gas combined heat and power devices are as follows:(20)Pi,tDGC=ηiDGCGi,tDGC1φi,tDGC+1,Hi,tDGC=ηiDGCGi,tDGCφi,tDGCφi,tDGC+1,φi,tDGC=PrhPrp,Ci,tDGC=θi,tDGCGi,tDGC,where Pi,tDGC and Hi,tDGC are the power and heat rates of the distributed gas combined heat and power devices generation, Gi,tDGC is the gas rate of coal of the distributed gas combined heat and power devices consumption, ηiDGC is the energy conversion of the distributed gas combined heat and power devices, and φi,tDGC is the ratio of heat and power. The constraints of the distributed gas combined heat and power devices is: Gi,tDGC−min≤Gi,tDGC≤1/n×Gnetwork, which is similar to that in the distributed power to gas devices. The constraints of φi,tDGC is similar to that in the distributed coal or oil combined heat and power devices.


## 2.4.1. The Distributed Power to Gas Devices and the Distributed Electric Boiler Models

The models and cost functions of the distributed power to gas devices are as follows:(18)Gi,tDP2G=ηi,tDP2GPi,tDP2G,Ci,tDP2G=θi,tDP2GPi,tDP2G,where ηi,tDP2G is the energy conversion efficiency of the distributed power to gas devices. θi,tDP2G is a positive constant. The limit of the distributed power to gas devices is(19)Pi,tDP2G−min≤Pi,tDP2G≤1n×Pnetwork,where Pi,tDP2G−min is the start-stop constraint of the distributed power to gas devices, Pnetwork is the surplus power rate in the network of EI. The models, operating cost functions, and limits of the distributed electric boiler are similar to that of the distributed power to gas devices. We only need to replace gas with heat.


## 2.4.2. The Distributed Gas Combined Heat and Power Device Models

The model and operating cost functions of the distributed gas combined heat and power devices are as follows:(20)Pi,tDGC=ηiDGCGi,tDGC1φi,tDGC+1,Hi,tDGC=ηiDGCGi,tDGCφi,tDGCφi,tDGC+1,φi,tDGC=PrhPrp,Ci,tDGC=θi,tDGCGi,tDGC,where Pi,tDGC and Hi,tDGC are the power and heat rates of the distributed gas combined heat and power devices generation, Gi,tDGC is the gas rate of coal of the distributed gas combined heat and power devices consumption, ηiDGC is the energy conversion of the distributed gas combined heat and power devices, and φi,tDGC is the ratio of heat and power. The constraints of the distributed gas combined heat and power devices is: Gi,tDGC−min≤Gi,tDGC≤1/n×Gnetwork, which is similar to that in the distributed power to gas devices. The constraints of φi,tDGC is similar to that in the distributed coal or oil combined heat and power devices.


## 2.5. The Distributed Energy Storage Devices Models and Cost Functions

### 2.5.1. The Distributed Energy, the Distributed Power Store Device, the Distributed Heat Storage Devices, and the Distributed Gas Storage Device Models

The distributed power storage devices have an optimal condition about reserve of stored power. If the stored power is less than optimal condition too much, it will harm batteries and other devices. If stored power is more than optimal condition too much, the stored power will loss too much. So, the optimal performance of the distributed power store device is(21)Oi,tDPSD=aiDPSDPi,tS−DPSDPi,tS−DPSD−2μiDPSD+biDPSD,where O is the optimal function of the reserved power in the distributed power store device, Pi,tS−DPSD is the stored power in time t, μiDPSD is the optimal reserve power in the distributed power store device, aiDPSD and biDPSD are constants, and aiDPSD is a negative constant. The cost function of the distributed power store device is(22)Ci,tDPSD=Oi,t−1DPSD−Oi,tDPSD+θi,tDPSDPi,tDPSD2,where θi,tDPSD is a positive constant and Pi,tDPSD is the power rate of the distributed power store device. If Pi,tDPSD is positive, it expresses power output, vice versa. So, we can know(23)Pi,tDPSD=Pi,t−1S−DPSD−Pi,tS−DPSD.The distributed power store device needs to satisfy the limits below:(24)−Piin−DPSD≤Pi,tDPSD≤Piout−DPSD,Pimin−S−DPSD≤Pi,tS−DPSD≤Pimax−S−DPSD,where Piin−DPSD and Piout−DPSD are the maximum limits of power input and output rate of the distributed power store device, and Pimin−S−DPSD andPimax−S−DPSD are capacity minimum and maximum values, respectively.The models, operating cost functions, and limits of the distributed heat storage devices and the distributed gas storage devices are similar to that of the distributed power store device. We only need to replace power with heat or gas.The energy loss in the distributed energy transform devices is very loss, for convenience, so we ignored that.


## 2.5.1. The Distributed Energy, the Distributed Power Store Device, the Distributed Heat Storage Devices, and the Distributed Gas Storage Device Models

The distributed power storage devices have an optimal condition about reserve of stored power. If the stored power is less than optimal condition too much, it will harm batteries and other devices. If stored power is more than optimal condition too much, the stored power will loss too much. So, the optimal performance of the distributed power store device is(21)Oi,tDPSD=aiDPSDPi,tS−DPSDPi,tS−DPSD−2μiDPSD+biDPSD,where O is the optimal function of the reserved power in the distributed power store device, Pi,tS−DPSD is the stored power in time t, μiDPSD is the optimal reserve power in the distributed power store device, aiDPSD and biDPSD are constants, and aiDPSD is a negative constant. The cost function of the distributed power store device is(22)Ci,tDPSD=Oi,t−1DPSD−Oi,tDPSD+θi,tDPSDPi,tDPSD2,where θi,tDPSD is a positive constant and Pi,tDPSD is the power rate of the distributed power store device. If Pi,tDPSD is positive, it expresses power output, vice versa. So, we can know(23)Pi,tDPSD=Pi,t−1S−DPSD−Pi,tS−DPSD.The distributed power store device needs to satisfy the limits below:(24)−Piin−DPSD≤Pi,tDPSD≤Piout−DPSD,Pimin−S−DPSD≤Pi,tS−DPSD≤Pimax−S−DPSD,where Piin−DPSD and Piout−DPSD are the maximum limits of power input and output rate of the distributed power store device, and Pimin−S−DPSD andPimax−S−DPSD are capacity minimum and maximum values, respectively.The models, operating cost functions, and limits of the distributed heat storage devices and the distributed gas storage devices are similar to that of the distributed power store device. We only need to replace power with heat or gas.The energy loss in the distributed energy transform devices is very loss, for convenience, so we ignored that.


## 2.6. Energy Load Models

The must-run load and controllable load which traditional research including Ref. [20, 21] adopts may be unsuitable for EI because of several reasons that I proposed in the preceding part of this study. So, we only adopt Pi,tELHi,tEL, and Gi,tEL to represent power heat and gas loads. If EI supplies less energy to WEs than which they need, the terminal users’ energy demand will not be satisfied. But if EI supplies more energy to WEs than which they need, it will waste a lot of energy. So, the operating cost of EL depends on the unbalance between energy supply and demand. The cost function of power EL is as follows:(25)CtP−EL=θi,tP−ELPi,tsupply−Pi,tEL22,θi,tP−EL=θ1,tP−EL,Pi,tsupply≤Pi,tEL,θ2,tP−EL,Pi,tsupply≥Pi,tEL,where θ1,tP−EL and θ2,tP−EL are two positive constants. The reason for dividing θtP−EL into θ1,tP−EL and θ2,tP−EL is that the harm of surplus and lack of power is different. The cost functions of heat and gas are similar to that of power. The only difference is the type of energy.


## 2.7. Object Functions

This study focuses on the economical optimization by energy management in WEs and EI. So, the object function is used to co-planning all WEs to realize the maximize of the economy in (22):(26)maxF=∑i=1nWi,t,(27)Wi,t=∑j=1mPrp×ΔPij+Prh×ΔHij+Prg×ΔGij−Prf×Fij−Cij,where Wi,t is the economy function of each WE; j is the serial number of each devices in WE; Prg and Prf are price of gas and coal; ΔPij, ΔHij, andΔGij are energy flow rate of power, heat, and gas, respectively. If they are positive, it means that the devices output the corresponding energy, vice versa. Of course, they can be zero. Fij is consumption of coal, which is non-negative. Cij is the operating cost function of devices.


## 3. The Annular Directed Distributed Algorithm and Its Certifications

### 3.1. Basic Knowledge of Graph Theory

The graph theory adoptsGraph=V,E,B to represent a graph. V=vi|i=1,2,…,n is a finite nonempty set of nodes in the graph. And the E⊆V×V is the set of sides. The side vi,vj means that there is a side between vi and vj. B=bi,j∈Rm×n is the weighting neighbor matrix of the graph. The weighting neighbor matrix B is used to express the relationship between nodes and sides. The diagonal elements in that are all zero constantly. If the non-diagonal element bi,j>0, vi,vj∈E. But if bi,j=0, vi,vj∉E. In the undirected graph, sides between nodes are not directed, so vi,vj∈E equals to vj,vi∈E. The paths between vi and vj consist of a lot of sides like vi,vi1, vi1,vi2⋯vik,vij. If there is a path between vi and vj, we define that vi andvj are connected. If entire pairs of nodes are connected, we define that the graph is connected. If all elements in B except diagonal elements are positive, the graph is complete.


### 3.2. The Relationship between Graph Theory and EI

Traditional researches [19, 20] use an undirected connected graph to model EI, use nodes to model WE or other agents about energy microgrid, and use sides to model channels which can transfer information and energy because they think the undirected connected graph can not only control all WEs in EI due to the path between entire pairs of nodes but also protect users’ privacy data which could only be acquired by neighbor nodes. They may be right to some extent. Agents in the unconnected graph cannot communicate among all nodes and that in the complete graph cannot protect users’ privacy. But that model has some serious disadvantages. Firstly, agents’ privacy data can be acquired by neighbor agents, which cannot protect users’ privacy entirely. Secondly, only-neighbor communicating causes a lot of trouble for the distributed algorithm in comparison with that in the complete graph. As for solving this problem, some researches pursue the consensus of tiny growth rate about voltage, air pressure in gas pipeline, and flow of hot water [29]. Although the consensus makes the distributed algorithms easily, it may not be the most optimal method for energy management. Thirdly, if an agent is strict with privacy and unwilling to communicate with other agents, the only way is to segment it from the whole system of EI. That agent will be operated in island mode and other agents will be operated in grid-connected mode, which is very bad for unifying dispatch of EI. Fourthly, only-neighbor energy transformation leads to a lot of energy loss and operating costs. Information has privacy, but energy flow does not have, after all! Last but not least, the connected but most complete graph is not well distributed in compactness, which is very bad for controllability of EI in accordance with the theory [24] proposed.For solving these problems, this study proposed a new algorithm named the annular directed distributed algorithm. In the annular directed distributed algorithm, the path of information and energy is separate. The information path is a directed annular path, which is shown in Figure1 and the hole energy path constitutes an undirected complete graph. After ensuring all devices’ operating conditions, how to transmit corresponding energy from the supply side to the demand side will be an elementary math problem that does not need algorithms. So, the most significant problem of energy management in the annular directed distributed algorithm is how to ensure devices’ optimal operating condition in the information side.Figure 1
Gas supply-demand mismatch (kw) in 100 times ofw.


### 3.3. Main Algorithm

The proposed management in the annular directed distributed algorithm is to ensure devices’ optimal operating condition in the information side. The main algorithm is shown in Table1. The method in step 12 is:Table 1
Annular distributed algorithm.
Initialize and iterate:k=1,i=1,r=1,dh=1;(1) Set initial value ofPrpPrhPrg, w.(2) Input all reserves of stored energy of each of the distributed energy store devices in timet−1.(3) Set initial value of energy input of all the distributed energy conversion devices in a small value in all WEs, includingPi,tDP2G, Pi,tDEB, and Gi,tDGC.(4) Initialize power-gas-heat energy flowPx,t, Hx,t, andGx,t and power-gas-heat energy load flow PL,tall, HL,tall, and GL,tall. Input them into control center of terminal load in first WE. Turn to first WE.(5) StorePx,tHx,tGx,t, PL,tall, HL,tall, and GL,tall by the designations of Px,tsHx,tsGSx,t, PL,tS−all, HL,tS−all, and GL,tS−all in control center of terminal load in first WE.(6)PL,tall, HL,tall, and GL,tall equals to the sum of corresponding type of energy loads of ith WE and corresponding PL,tall, HL,tall, and GL,tall.(7) Ifi<n, i=i+1, transmit PL,tall, HL,tall, and GL,tall to ith WE, and return step 6, else, continue.(8) Ifi=1, transmit PL,tall, HL,tall, and GL,tall to control center of terminal load in first WE. PL,tall, HL,tall, and GL,tall is equal to themselves subtracting corresponding PL,tS−all, HL,tS−all, and GL,tS−all which are stored in first WE.(9) Compute and store ratios of load toPL,tallHL,tall, and GL,tall by the designations of Pi,tr, Hi,tr, and Gi,tr in control center of terminal load in ith WE.(10) Ifi<5, i=i+1, transmit PL,tall, HL,tall, and GL,tall to ith WE, and return step 9, else, i=1, continue.(11) Update entire operating conditions containing energy input or output volumes of all the distributed energy producer and the distributed energy store devices inith WE to the value that (26) and (27) can get peak value. Because all parts of (27) are linear or non-convex, the work is very easy and does not need complex algorithm.(12) Update entire operating conditions containing energy input or output volumes of all the distributed energy conversion devices inith WE by method described below.(13) Add all energy input and output volumes in all the distributed energy producer, the distributed energy conversion devices, and the distributed energy store devices inith WE to Px,t, Hx,t, andGx,t. If i<n, i=i+1, transmit the volumes of Px,t, Hx,t, and Gx,t to ith WE and return to the step 11, else i=1, transmit the volumes of Px,t, Hx,t, and Gx,t to ith WE, continue.(14) If all types of energy supply-demand balance are reached, energy management finishes, else, if one type of energy supply-demand mismatch absolute value is higher than that in last time of iteration,dh=0.5×dh, return step 11, else, if one type of energy supply-demand mismatch changes symbol, r=0.5×r, return step 3, else, r=1, dh=1k=k+1, return step 11.If we want to change the typea of energy to the type b of energy, the energy input will be updated by the method that:(28)Ia,i,tk+1=Ia,i,tkdhηab×PrbPra+θab−1+1,where I is the input of energya (the starting value of I is initialized by step 3), k is the number of iteration times, ηab is the efficiency of energy conversion, and Pra and Prb are the energy prices of energy a and b. What is noteworthy is that if a device can change one type of energy to two or more types of energies, the price of output energies is the weighted average of all output energies, and weighted factors are the rate ratio of that energy to hole output energies. θab is the energy conversion price from energy a to energy b., dh is a downhill factor whose effect is to adjust devices’ operating condition in case of divergence of the annular directed distributed algorithm. dh does not change between different types of the distributed energy conversion devices, the starting value in each iteration of dh is 1, and the updating conditions and methods of that will be introduced below.The method in step 14 is as follows:(29)Pi,tsupply=Pi,tr×Px,t,Hi,tsupply=Hi,tr×Hx,t,Gi,tsupply=Gi,tr×Gx,t.wherePi,tsupply, Hi,tsupply, and Gi,tsupply are the power, heat, and gas which the energy management supplied to energy load in ith WE, respectively. But they may not be equal to Pi,tEL, Hi,tEL, and Gi,tEL. So, we should adjust energy prices. We use cost functions of ELs as subgradient functions to adjust the energy price. We adjust power price like Figure 2. w in Figure 2 is a small positive constant, which may be different with different types of energy. By the way, the precision of Psd andPsm is no possibility infinite in practical engineering, so if they are very less theoretically, we regard them as zero. r is another Newton-downhill factor about avoiding overshoot, if one type of energy supply-demand mismatch changes symbol, r=0.5×r.Figure 2
Power price adjusting.The method to adjust the price of gas and heat is similar to that of power, so we only need to change the energy type. After price adjustment, we should inspect the astringency of the algorithm. Whichever supply-demand mismatching in all types of energy absolute value ink+1 times of iteration is more than that in k times of iteration, or supply-demand mismatching in all types of energy value changes symbol, the result in step 14 must be invalid. The algorithm will return step 12 or step 3, respectively, and energy management will amend downhill factor dh or r to half of quondam value and repeat processes above until all energy supply-demand mismatching absolute value in k+1 times of iteration is no more than that in k times of iteration. Besides, the value of dh or r in k times of iteration should not be retained into k+1 times of iteration. It will initialize into 1 again in that. What is always ignored by traditional research but is very significant to energy management is there is a problem of Zeno behavior. Zeno behaviors mean some operation is triggered infinite times in finite times of iteration. In the annular directed distributed algorithm, Zeno behaviors are that the frequency of that r or dh change is too high in one time of iteration.It is worth noting that, because accuracy in down-to-earth engineering is limited, if the energy supply-demand mismatch is very less, we can regard energy supply-demand balance as achieved and stop algorithm. Besides, if we do not stop the algorithm in that circumstance, the algorithm result may be choppy near the exact value due to the energy mismatch is less than that the small positive constantw adjusted.Four things need to be added. Firstly, after one time of iteration, the sequence of WEs should be changed. Theith WE will be the i+1th WE, and the first WE will be the last WE. We should transmit all information from the old first WE to the new first WE before iteration. The reason for that is too much control and computing load is at the last WE. Changing the sequence can relieve the load of the last WE. Besides, it can also increase the precision and universality of algorithm results. Secondly, all initialized energy flow and energy load flow is needed to be initialized over and over again in each iteration. The annular directed distributed algorithm has a great advantage in which it can protect each user’s privacy much better than traditional EI or traditional multiagent systems [25–27].Although preceding researchers can protect users’ privacy which cannot be searched by all agents, neighbor agents still can access their private data. But in the annular directed distributed algorithm, all agents can only access total before themselves. However, the privacy of the first agent cannot be protected by this method, hence we use fictitious initial energy flows to protect it. For reinforcing the confidentiality of privacy data, initial energy flows should be changed each iteration. Thirdly, privacy protecting is needed not only among different WEs, but also inside WEs. The information inside WEs transmitted in one directed line whose direction is identical to the direction of the information path between WEs. Fourthly, although the annular information path has numerous advantages, it is exceeding sensitive to malfunctions. For handling this issue, if a device in one WE is broken down, the energy load in that WE will operate by island mode and other WEs will operate by a mode combined to the grid as normal. Finally, in the annular directed distributed algorithm, the information path is different from the energy transmission path. The energy transmission path is independent of privacy protecting, so all devices can transmit energy to others. After realizing a balance between energy supply and energy demand, how to transmit energy in the shortest way is an easy problem which does not need complex algorithms, so this study does not discuss energy transmission path. It is worth noting that the optimization in this article is distributed. The principle is to adjust energy price to change the energy mismatch. If the adjust of energy price is overshoot to the optimization valve, the algorithm will be oscillation or divergency and the energy supply-demand balance will never be reached. So, we should adopt the Newton-downhill method to avoid the overshoot of energy price.We list the data exchange and privacy protection that is required by the proposed distributed algorithm and conventional distributed methods in Tables2 and 3 (✔ for entirely exchanging, ✖ for not exchanging, ● for cryptographically exchanging) (✔ for protecting, ✖ for not protecting).Table 2
Privacy protecting between we-energies.
Type of data exchange among neighbor we-energiesThe proposed distributed algorithmConventional distributed methodsPrivacy protection for that data in proposed distributed algorithmPrivacy protection for that data in conventional distributed methodsPower output or input ratio●✔✔✖Heat output or input ratio●✔✔✖Gas output or input ratio●✔✔✖Power price✔✔✖✖Heat price✔✔✖✖Gas price✔✔✖✖Table 3
Privacy protecting inside we-energies.
Type of data exchange inside each we-energyThe proposed distributed algorithmConventional distributed methodsPrivacy protection for that data in proposed distributed algorithmPrivacy protection for that data in conventional distributed methodsPower output or input ratio●✔✔✖Heat output or input ratio●✔✔✖Gas output or input ratio●✔✔✖Power price✔✔✖✖Heat price✔✔✖✖Gas price✔✔✖✖From that two tables you can know that the proposed algorithm can protect the data of power-heat-gas output or input conditions of each device better. It is worth noting that, all manufactures of one type of energy are commercial competitors, so are all we-energies. So, it is unreasonable for conventional distributed methods to regard some energy manufactures and some we-energies as neighbor agents and exchange data among them. The data cryptographically exchanging method in this article is a better method.


## 3.1. Basic Knowledge of Graph Theory

The graph theory adoptsGraph=V,E,B to represent a graph. V=vi|i=1,2,…,n is a finite nonempty set of nodes in the graph. And the E⊆V×V is the set of sides. The side vi,vj means that there is a side between vi and vj. B=bi,j∈Rm×n is the weighting neighbor matrix of the graph. The weighting neighbor matrix B is used to express the relationship between nodes and sides. The diagonal elements in that are all zero constantly. If the non-diagonal element bi,j>0, vi,vj∈E. But if bi,j=0, vi,vj∉E. In the undirected graph, sides between nodes are not directed, so vi,vj∈E equals to vj,vi∈E. The paths between vi and vj consist of a lot of sides like vi,vi1, vi1,vi2⋯vik,vij. If there is a path between vi and vj, we define that vi andvj are connected. If entire pairs of nodes are connected, we define that the graph is connected. If all elements in B except diagonal elements are positive, the graph is complete.


## 3.2. The Relationship between Graph Theory and EI

Traditional researches [19, 20] use an undirected connected graph to model EI, use nodes to model WE or other agents about energy microgrid, and use sides to model channels which can transfer information and energy because they think the undirected connected graph can not only control all WEs in EI due to the path between entire pairs of nodes but also protect users’ privacy data which could only be acquired by neighbor nodes. They may be right to some extent. Agents in the unconnected graph cannot communicate among all nodes and that in the complete graph cannot protect users’ privacy. But that model has some serious disadvantages. Firstly, agents’ privacy data can be acquired by neighbor agents, which cannot protect users’ privacy entirely. Secondly, only-neighbor communicating causes a lot of trouble for the distributed algorithm in comparison with that in the complete graph. As for solving this problem, some researches pursue the consensus of tiny growth rate about voltage, air pressure in gas pipeline, and flow of hot water [29]. Although the consensus makes the distributed algorithms easily, it may not be the most optimal method for energy management. Thirdly, if an agent is strict with privacy and unwilling to communicate with other agents, the only way is to segment it from the whole system of EI. That agent will be operated in island mode and other agents will be operated in grid-connected mode, which is very bad for unifying dispatch of EI. Fourthly, only-neighbor energy transformation leads to a lot of energy loss and operating costs. Information has privacy, but energy flow does not have, after all! Last but not least, the connected but most complete graph is not well distributed in compactness, which is very bad for controllability of EI in accordance with the theory [24] proposed.For solving these problems, this study proposed a new algorithm named the annular directed distributed algorithm. In the annular directed distributed algorithm, the path of information and energy is separate. The information path is a directed annular path, which is shown in Figure1 and the hole energy path constitutes an undirected complete graph. After ensuring all devices’ operating conditions, how to transmit corresponding energy from the supply side to the demand side will be an elementary math problem that does not need algorithms. So, the most significant problem of energy management in the annular directed distributed algorithm is how to ensure devices’ optimal operating condition in the information side.Figure 1
Gas supply-demand mismatch (kw) in 100 times ofw.


## 3.3. Main Algorithm

The proposed management in the annular directed distributed algorithm is to ensure devices’ optimal operating condition in the information side. The main algorithm is shown in Table1. The method in step 12 is:Table 1
Annular distributed algorithm.
Initialize and iterate:k=1,i=1,r=1,dh=1;(1) Set initial value ofPrpPrhPrg, w.(2) Input all reserves of stored energy of each of the distributed energy store devices in timet−1.(3) Set initial value of energy input of all the distributed energy conversion devices in a small value in all WEs, includingPi,tDP2G, Pi,tDEB, and Gi,tDGC.(4) Initialize power-gas-heat energy flowPx,t, Hx,t, andGx,t and power-gas-heat energy load flow PL,tall, HL,tall, and GL,tall. Input them into control center of terminal load in first WE. Turn to first WE.(5) StorePx,tHx,tGx,t, PL,tall, HL,tall, and GL,tall by the designations of Px,tsHx,tsGSx,t, PL,tS−all, HL,tS−all, and GL,tS−all in control center of terminal load in first WE.(6)PL,tall, HL,tall, and GL,tall equals to the sum of corresponding type of energy loads of ith WE and corresponding PL,tall, HL,tall, and GL,tall.(7) Ifi<n, i=i+1, transmit PL,tall, HL,tall, and GL,tall to ith WE, and return step 6, else, continue.(8) Ifi=1, transmit PL,tall, HL,tall, and GL,tall to control center of terminal load in first WE. PL,tall, HL,tall, and GL,tall is equal to themselves subtracting corresponding PL,tS−all, HL,tS−all, and GL,tS−all which are stored in first WE.(9) Compute and store ratios of load toPL,tallHL,tall, and GL,tall by the designations of Pi,tr, Hi,tr, and Gi,tr in control center of terminal load in ith WE.(10) Ifi<5, i=i+1, transmit PL,tall, HL,tall, and GL,tall to ith WE, and return step 9, else, i=1, continue.(11) Update entire operating conditions containing energy input or output volumes of all the distributed energy producer and the distributed energy store devices inith WE to the value that (26) and (27) can get peak value. Because all parts of (27) are linear or non-convex, the work is very easy and does not need complex algorithm.(12) Update entire operating conditions containing energy input or output volumes of all the distributed energy conversion devices inith WE by method described below.(13) Add all energy input and output volumes in all the distributed energy producer, the distributed energy conversion devices, and the distributed energy store devices inith WE to Px,t, Hx,t, andGx,t. If i<n, i=i+1, transmit the volumes of Px,t, Hx,t, and Gx,t to ith WE and return to the step 11, else i=1, transmit the volumes of Px,t, Hx,t, and Gx,t to ith WE, continue.(14) If all types of energy supply-demand balance are reached, energy management finishes, else, if one type of energy supply-demand mismatch absolute value is higher than that in last time of iteration,dh=0.5×dh, return step 11, else, if one type of energy supply-demand mismatch changes symbol, r=0.5×r, return step 3, else, r=1, dh=1k=k+1, return step 11.If we want to change the typea of energy to the type b of energy, the energy input will be updated by the method that:(28)Ia,i,tk+1=Ia,i,tkdhηab×PrbPra+θab−1+1,where I is the input of energya (the starting value of I is initialized by step 3), k is the number of iteration times, ηab is the efficiency of energy conversion, and Pra and Prb are the energy prices of energy a and b. What is noteworthy is that if a device can change one type of energy to two or more types of energies, the price of output energies is the weighted average of all output energies, and weighted factors are the rate ratio of that energy to hole output energies. θab is the energy conversion price from energy a to energy b., dh is a downhill factor whose effect is to adjust devices’ operating condition in case of divergence of the annular directed distributed algorithm. dh does not change between different types of the distributed energy conversion devices, the starting value in each iteration of dh is 1, and the updating conditions and methods of that will be introduced below.The method in step 14 is as follows:(29)Pi,tsupply=Pi,tr×Px,t,Hi,tsupply=Hi,tr×Hx,t,Gi,tsupply=Gi,tr×Gx,t.wherePi,tsupply, Hi,tsupply, and Gi,tsupply are the power, heat, and gas which the energy management supplied to energy load in ith WE, respectively. But they may not be equal to Pi,tEL, Hi,tEL, and Gi,tEL. So, we should adjust energy prices. We use cost functions of ELs as subgradient functions to adjust the energy price. We adjust power price like Figure 2. w in Figure 2 is a small positive constant, which may be different with different types of energy. By the way, the precision of Psd andPsm is no possibility infinite in practical engineering, so if they are very less theoretically, we regard them as zero. r is another Newton-downhill factor about avoiding overshoot, if one type of energy supply-demand mismatch changes symbol, r=0.5×r.Figure 2
Power price adjusting.The method to adjust the price of gas and heat is similar to that of power, so we only need to change the energy type. After price adjustment, we should inspect the astringency of the algorithm. Whichever supply-demand mismatching in all types of energy absolute value ink+1 times of iteration is more than that in k times of iteration, or supply-demand mismatching in all types of energy value changes symbol, the result in step 14 must be invalid. The algorithm will return step 12 or step 3, respectively, and energy management will amend downhill factor dh or r to half of quondam value and repeat processes above until all energy supply-demand mismatching absolute value in k+1 times of iteration is no more than that in k times of iteration. Besides, the value of dh or r in k times of iteration should not be retained into k+1 times of iteration. It will initialize into 1 again in that. What is always ignored by traditional research but is very significant to energy management is there is a problem of Zeno behavior. Zeno behaviors mean some operation is triggered infinite times in finite times of iteration. In the annular directed distributed algorithm, Zeno behaviors are that the frequency of that r or dh change is too high in one time of iteration.It is worth noting that, because accuracy in down-to-earth engineering is limited, if the energy supply-demand mismatch is very less, we can regard energy supply-demand balance as achieved and stop algorithm. Besides, if we do not stop the algorithm in that circumstance, the algorithm result may be choppy near the exact value due to the energy mismatch is less than that the small positive constantw adjusted.Four things need to be added. Firstly, after one time of iteration, the sequence of WEs should be changed. Theith WE will be the i+1th WE, and the first WE will be the last WE. We should transmit all information from the old first WE to the new first WE before iteration. The reason for that is too much control and computing load is at the last WE. Changing the sequence can relieve the load of the last WE. Besides, it can also increase the precision and universality of algorithm results. Secondly, all initialized energy flow and energy load flow is needed to be initialized over and over again in each iteration. The annular directed distributed algorithm has a great advantage in which it can protect each user’s privacy much better than traditional EI or traditional multiagent systems [25–27].Although preceding researchers can protect users’ privacy which cannot be searched by all agents, neighbor agents still can access their private data. But in the annular directed distributed algorithm, all agents can only access total before themselves. However, the privacy of the first agent cannot be protected by this method, hence we use fictitious initial energy flows to protect it. For reinforcing the confidentiality of privacy data, initial energy flows should be changed each iteration. Thirdly, privacy protecting is needed not only among different WEs, but also inside WEs. The information inside WEs transmitted in one directed line whose direction is identical to the direction of the information path between WEs. Fourthly, although the annular information path has numerous advantages, it is exceeding sensitive to malfunctions. For handling this issue, if a device in one WE is broken down, the energy load in that WE will operate by island mode and other WEs will operate by a mode combined to the grid as normal. Finally, in the annular directed distributed algorithm, the information path is different from the energy transmission path. The energy transmission path is independent of privacy protecting, so all devices can transmit energy to others. After realizing a balance between energy supply and energy demand, how to transmit energy in the shortest way is an easy problem which does not need complex algorithms, so this study does not discuss energy transmission path. It is worth noting that the optimization in this article is distributed. The principle is to adjust energy price to change the energy mismatch. If the adjust of energy price is overshoot to the optimization valve, the algorithm will be oscillation or divergency and the energy supply-demand balance will never be reached. So, we should adopt the Newton-downhill method to avoid the overshoot of energy price.We list the data exchange and privacy protection that is required by the proposed distributed algorithm and conventional distributed methods in Tables2 and 3 (✔ for entirely exchanging, ✖ for not exchanging, ● for cryptographically exchanging) (✔ for protecting, ✖ for not protecting).Table 2
Privacy protecting between we-energies.
Type of data exchange among neighbor we-energiesThe proposed distributed algorithmConventional distributed methodsPrivacy protection for that data in proposed distributed algorithmPrivacy protection for that data in conventional distributed methodsPower output or input ratio●✔✔✖Heat output or input ratio●✔✔✖Gas output or input ratio●✔✔✖Power price✔✔✖✖Heat price✔✔✖✖Gas price✔✔✖✖Table 3
Privacy protecting inside we-energies.
Type of data exchange inside each we-energyThe proposed distributed algorithmConventional distributed methodsPrivacy protection for that data in proposed distributed algorithmPrivacy protection for that data in conventional distributed methodsPower output or input ratio●✔✔✖Heat output or input ratio●✔✔✖Gas output or input ratio●✔✔✖Power price✔✔✖✖Heat price✔✔✖✖Gas price✔✔✖✖From that two tables you can know that the proposed algorithm can protect the data of power-heat-gas output or input conditions of each device better. It is worth noting that, all manufactures of one type of energy are commercial competitors, so are all we-energies. So, it is unreasonable for conventional distributed methods to regard some energy manufactures and some we-energies as neighbor agents and exchange data among them. The data cryptographically exchanging method in this article is a better method.


## 4. Testification of Convergence Analysis, Optimality, and Avoiding Zeno Behaviors

The convergence of the algorithm and the balance between energy supply and energy demand are proved below.Firstly, we establish an assumption that the operating condition of the distributed energy conversion devices ink+1 times of the iteration is equal to those in k times of the iteration. Under that circumstance, all energy supply-demand is decided by formula (27). We discuss power supply-demand as an example. We can calculate power input or output in each of the distributed energy producer and the distributed energy store devices by solving the partial differential equation below.(30)∂Wi,t∂ΔPij=0.And we can disassemble Wi,t by:(31)Q=∑j=1mPrh×ΔHij+Prg×ΔGij−Prf×Fij−Cij,Wi,t=∑j=1mPrp×ΔPij+Q.So, we can transform (28) into(32)∂Q∂ΔPij=−Prp.Due to the non-convex ofQ, ∂Q/∂ΔPij is a decreasing function. If power is excess in EI, the price of that will decrease. So, the solution of (32) will be less, vice versa. Moreover, because of the downhill factor r, if overshoot will never appear, r will reduce until overshoot vanishes. And the change of r only alters the adjusting extent of power price but does nothing about adjusting the direction of that. So r does not affect the monotonicity of power price. It is worth noting that the monotonicity of power price aims at power price when adjusting of r finishes. The power price of the overshoot systems needs to be abandoned so it is meaningless to energy management. Above all, we can draw a conclusion in line with mathematical induction that under the assumption that the operating condition of the distributed energy conversion devices in k+1 times of the iteration is equal to that in k times of the iteration, the power price, which is one of the results of the algorithm, is monotone. However, without the assumption, is the power price monotone? The answer is right. A contradiction method is used to prove it.First of all, we set up an assumption that the power price is not monotone. So, step 12 to step 14 will be an endless loop anddh will be halved over and over again. Hence,(33)dh=limm⟶+∞12m=0.So, the operating condition of the distributed energy conversion devices ink+1 times of the iteration is equal to that in k times of the iteration. On the basis of the conclusions mentioned above, the algorithm is monotone, which is in contradiction with the assumption at the beginning of the contradiction method. So, we can prove that the power price is monotone. It is worth noting that the proof process does not contain the circumstance that when the power supply-demand mismatch tends to a very small value. If that mismatch is very small, we will deem that the power supply-demand is balanced and we will regard it and Psm as zero in practical engineering. Due to the finite precision of measurement and operation in practice, it is useless to discuss the small power supply-demand mismatch.Becausew is small, power supply-demand mismatch cannot change direction without going through a small value. And due to the reason above, if that is small, we will regard Psm as 0 so the power price will not update. So, the power price is unilateral bounded and the bound is the value when that mismatch is exceedingly small.In accordance with the monotone-bounded convergence theorem and conclusions all above, we can prove that the power price is convergence. We can figure out the value of power price when it is convergence by the way that:(34)Prpk+1=Prpk×1−0⋅2×w×Psm,Prpk+1=Prpk.So,Psm=0, the price of power will converge to a value that can make the power supply-demand balance. Similarly, the price of heat and gas will converge to a value that can reach the heat or gas supply-demand balance, too. The proving process of them are analogical to that of power, and we only need to change the types of energy. To sum up, after numerous times of iteration, all types of energy prices are converged and all kinds of energy are balanced between energy supply and energy demand. EI is stabilized.The optimality of the convergence point in the annular directed distributed algorithm is proved below.Similar to convergence proving, the attesting of optimality in energy management needs to prove optimality about power systems without energy conversion first.(35)∂Wi,t∂Pij=Prp−dCPdP.BecauseCP is a convex function, dCP/dP is an increasing function. Step 11 in the annular directed distributed algorithm requires all power-producing devices to choose operating conditions when dCP/dP=Prp. At the convergence point, power supply-demand balance is satisfied. We assume another operating condition of EI as device A produce ΔP power less than its power-producing in the annular directed distributed algorithm. For maintaining power supply-demand balance, device B will produce ΔP power more. Cost decrease in A is ΔP×ΔCAΔP. ΔCAΔP is backward difference gradient in the convergence point of device A. Analogically, cost decrease in B is ΔP×ΔCBΔP. ΔCBΔP is forward difference gradient in the convergence point of device B. Because all CP are convex functions, ΔCAΔP<Prp<ΔCBΔP. So, the convergence point is optimal for power without energy conversion. By the same token, the convergence point is optimal for heat and gas without energy conversion. And we should prove optimality among energy conversion.Becausew is a very small constant, energy price in a large number of iteration only changes a little. Hence, the speed of energy price changing is much lower than that of energy conversion changing so that we can ignore the change of energy price when we study energy conversion. If ηab×Prb is more than Pra+θab, due to (28), there is more and more energy that changes from a to b until ηab×Prb is equal to Pra+θab. Thereafter, ηab×Prb is equal to Pra+θab though the hole iteration. Whenever energy price changes a little, energy conversion will change at a far high speed to reach equality between ηab×Prb and Pra+θab and goes along with slowly energy price changing. That equality is still satisfied until convergence point. By the way, because of the difference between that two speeds, frequency about adjusting of dh is much less than that of r. So, Zeno behavior never happens in dh. However, the Zeno behavior never happens in r, either. And the reason for it will be described in this literature, too. We assume another operating condition of EI as energy a transforms Δa energy less than its energy conversion in the annular directed distributed algorithm. For maintaining power supply-demand balance, energy management should produce Δa energy a less and produce ηab×Δa energy more. The cost decrease in producing and transforming in a is Δa×ΔCaΔa+θab. And the cost increase in producing in b is Δa×ηab×ΔCbΔa×ηab. ΔCaΔa and ΔCbΔa×ηab are difference gradients. Because all C are convex functions, ΔCaΔa+θab<Pra+θab=ηab×Prb<ΔCbΔa×ηab×ηab. So Δa×ΔCaΔa+θab<Δa×ηab×ΔCbΔa×ηab, energy conversion in the convergence point is optimal. In addition, equal relations mentioned above do not consider the influence of constraint. So, it is no wonder that in reality EI project,To sum up, the convergence point in the annular directed distributed algorithm has optimality.Avoiding Zeno behavior of the annular directed distributed algorithm is proved below:On account of the reason in the preceding part of this literature, Zeno behavior never happens indh. But Zeno behavior in r still needs analysis.If overshoot happens, energy overshoot and energy mismatch in the last time of iteration are all finite for the reason of if they are infinitely small we will regard them as zero and energy supply-demand are reached. (Both infinitely small and infinitely great are all infinite.) The adjusting of energy price can be calculated through linear or quadratic functions (all functions in this literature are linear or quadratic) by four fundamental rules, gradient rules, and inverse function rules from finite values (energy overshoot and energy mismatch in last time of iteration). Hence, they are finite, and so adjusting times ofr in one or finite times of iteration is finite, thus Zeno behavior does not happen in r or dh.To sum up, the annular directed distributed algorithm is not only optimal but also convergence. In addition, Zeno behavior never occurs in it.


## 5. Simulation Results

The framework of the EI system used to test the annular directed distributed algorithm containing 5 WEs and data of devices all above are exhibited in Appendix. In this article, if all the types of energy mismatches are less than500kw, we think of the balance of energy-producing and energy-consuming as achieved and regard Psm as zero. After 42 times of iteration, that balance is reached. The simulation results are discoursed belowFigures3–5 are prices of power, heat, and gas. Figures 6–8 are supply-demand mismatches of power, heat, and gas. The Zeno coefficient is times of r adjusting. It shows that after 42 times of the iteration, when the Zeno coefficient is 2, power, heat, and gas mismatch are all lower than 500kw. At this moment kw⋅ht, the power mismatch is 462kw less. The heat mismatch is 387kw more. And the gas mismatch is near zero in experimental precision. Comparing with traditional distributed algorithm 28 which needs hundreds of times to reach energy supply-demand balance, the annular directed distributed algorithm is exceedingly quick. EI system should fix the prices of power heat and gas as 11 ￠ per kw⋅h, 9 ￠ per kw⋅h, and 8 ￠ per kw⋅h. (Due to accuracy of USA dollars, retain prices to integral multiples.) Figures 9 and 10 are times of r and dh adjusting in each time of iteration. Times of those are no more than twice as they show. So, Zeno behavior does not appear. 11 and 12 are the power-heat-gas mismatch in 100 times of w, and Figures 13–15 are the power-heat-gas mismatch in 0.01 times of w. They show that if the value of w is too large, the algorithm will vibrate seriously. If the value of w is too small, the algorithm will converge in a very slow speed. So, the algorithm in this article is very sensitive to the value of w. To this end, the plug and play performance of the proposed distributed method is low because if the EI framework and the energy load changes largely, it is very hard for us to adjust w. What's more, the proposed method ignores energy transmitting cost so some areas which are hard to transmit energy are unsuitable for the proposed method. However, some areas require large demand of privacy protecting, so they must adopt the proposed method.Figure 3
Power price.Figure 4
Heat price.Figure 5
Gas price.Figure 6
Power supply-demand mismatch.Figure 7
Heat supply-demand mismatch.Figure 8
Gas supply-demand mismatch.Figure 9
Maximum value of the Zeno coefficient.Figure 10
Maximum value of the vice-Zeno coefficient.Figure 11
Power supply-demand mismatch (kw) in 100 times ofw.Figure 12
Heat supply-demand mismatch (kw) in 100 times ofw.Figure 13
Power supply-demand mismatch (kw) in 0.01 times ofw.Figure 14
Heat supply-demand mismatch (kw) in 0.01 times ofw.Figure 15
Power supply-demand mismatch (kw) in 0.01 times ofw.EI is a new energy system in scientists’ conceive. Although EI has got large concern among researchers, there is not an unabridged EI system with a high permeability of renewable energy resources, a large scale of physical system, and a highly interconnected network among different kinds of energy nowadays because it has just been proposed. A lot of cases in articles which are imagined by authors may be absurd and unreasonable in actuality. The only thing we can know is the energy demand. We summarize various energy demands in the world in Table4 and analyze if they are fit for the proposed algorithm (✔ for fit, ✖ for unfit). Because there are too many kinds of energy demand in the whole world, our work may be incomplete.Table 4
Application scenario analysis.
Serial numberSceneRegional characteristicsTypical caseFit for traditional distributed methodsFit for the proposed distributed methodA1Suburb typeSmall area with low energy demandUrban✖✔A2Seasonal switching typeNeed heat in winter and cold in summer. The demand of power and gas is large.Cities in Yangtze plain, middle and lower, China✔✖A3Cold plateauThe atmospheric pressure is low because of high altitude, and the temperature is cold for the same reason. All energy demands are low.Qinghai-Tibet plateau, China✖✖A4Island (in sea)Wind and petroleum are richNansha six reefs, China✖✔A5Cold areaHeat load is highNortheast, China✔✖A6MountainousEnergy load is low. Wind and solar are rich.Southwest, China✖✔A7High latitude portDC load of ships is high. Heat load is high because the temperature is coldPort cities in Northern Europe and North America✖✔A8Dispersing areaThe natural gas is rich, and the energy load is dispersedNorthwest, China✖✔A1 is suitable for the proposed distributed method because the business competition in urban is sharp and the privacy protection demand is large. A2 is unsuitable for the proposed method because the energy load is changeable which cause great trouble ofw. A3 is unsuitable for the both traditional and proposed method because the heat and gas transmitting is hard in that area. They need another method about heat and gas transmitting. A4 is suitable for the proposed method because the privacy protecting in that area is important due to political factors. A5 is unsuitable for the proposed method because the privacy protecting in that area is not important. A6 is suitable for the proposed method because the change of that area is very slow and the w is easy to choose. A7 is suitable for the proposed method because the business competition in that area is serious so the privacy protecting is important. A8 is suitable for the proposed method because the change of that area is very slow so the w is easy to choose. To sum up, the proposed method is suitable for some cases with high demand of privacy protecting and some cases which change slowly. However, it is unsuitable for some cases with low demand of privacy protecting and some cases which change quickly.


## 6. Conclusions

Based on the traditional distributed and the alternating direction method of multipliers algorithm in energy market, this article proposed a new communicating and optimizing algorithm of energy management EI. EI in this algorithm communicates information in an annular way, which can not only greatly reduce communicating times, but also can protect the privacy data of all agents that can only be searched by itself. Compared with previous privacy protection in multiagent systems in which privacy data are transmitted to neighbor agents, privacy in this article is protected much more strongly. In addition, this article uses the subgradient method to quicken the annular directed distributed algorithm and solve the issue of astringency in energy transforming. Simulation results and theoretical identification containing the monotone-bounded convergence theorem, theory of limit, and contradiction have demonstrated the effectiveness of the proposed algorithm (see Table5).Table 5
Data.
Distributed renewable generators and distributed renewable heat devicesPi,twind, ×103kwPi,tsolar, ×103kwHi,tsolar, ×102kwPiwind−min, ×102kwPiwind−max, ×104kwPisolar−min, ×102kwHisolar−min, ×102kwHisolar−max, ×104kwWE13.98431.42540.45131.21825.24520.14201.10352.7131WE24.53242.28710.25411.15462.62640.31971.32143.6157WE30.58441.12502.25090.91374.25670.09842.31471.9548WE41.65341.00700.01310.89544.32180.12471.21811.7496WE51.24553.15430.30651.31933.21480.22141.31201.3427Distributed coal or oil combined heat and power devicesaiDCOC, ×10−4biDCOC, ×10−3αiDCOC, ×10−4βiDCOC, ×10−3ciDCOC, ×10−4χiDCOC, ×102diDCOCeiDCOCWE11.651.501.801.201.742.2536.1557.15WE21.350.451.951.501.722.8435.4059.55WE31.801.951.201.501.435.9238.7060.45WE41.651.801.352.101.473.7140.0563.15WE51.952.851.651.651.763.6939.3062.55Distributed coal or oil combined heat and power devicesfiDCOCgiDCOChiDCOCjiDCOCηiDCOCPiDCOC−ramp, kwPi,t−1DCOC, kwWE13203554887600.520.8689%10563908WE23237604797600.540.8784%9873847WE33054605188050.510.9191%8533691WE42980954823550.590.9092%9644058WE52775305463150.580.8895%12084061Distributed gas producersai,tDGP, ×10−4bi,tDGPci,tDGPGi,tDGP−minkwGi,tDGP−maxkwWE11.53.03601485201.35941.40Distributed power to gas devicesηi,tDP2Gθi,tDP2GPi,tDP2G−min, kwInitialized operating conditionkwWE177%0.03520520WE282%0.02610610WE381%0.04710710WE479%0.01610610WE584%0.02600600Distributed electric boilerηi,tDEBθi,tDEBPi,tDEB−min, kwInitialized operating conditionkwWE192%0.02980980WE295%0.04720720WE393%0.01840840WE497%0.0310301030WE591%0.05790790Distributed gas combined heat and power devicesθi,tDGCηiDGChiDGCjiDGCInitialized conditionGi,tDGC−min, kwWE10.0383%0.595.3120750750WE20.0279%0.555.2019690690WE30.0582%0.585.7187840840WE40.0184%0.575.2162720720WE50.0486%0.545.9146710710Distributed power storage devicesaiDPSDbiDPSDθi,tDPSDPiin−DPSD, kwPiout−DPSD, kwPimin−S−DPSD, kwPimax−S−DPSD, kwPi,t−1S−DPSD, kwμiDPSDWE11.121540.1543256727115681100625WE21.091360.1345961422617341025629WE30.931270.1750341219117211328735WE40.821980.2144145014014301327577WE51.331210.114044421421331954518Distributed heat storage devicesaiDHSDbiDHSDθi,tDHSDHiin−DHSD, kwHiout−DHSD, kwHimin−S−DHSD, kwHimax−S−DHSD, kwHi,t−1S−DHSD, kwWE11.302530.1456164117117651207534WE21.257680.1954059025415831011672WE30.713540.1357252536918521319601WE40.964170.1750353034217541021608WE51.041390.1450952324415431055512Distributed gas storage devicesaiDGSDbiDGSDθi,tDGSDGiin−DHSD, kwGiout−DHSD, kwGimin−S−DGSD, kwGimax−S−DGSD, kwGi,t−1S−DGSD, kwμiDGSDWE12.107340.23325371461206435542WE22.862610.14317384711671429531WE31.844580.12349363611087456568WE42.104300.2140131282960463591WE51.614190.19366350691427421601Energy loadPi,tEL×104kwHi,tEL, ×103kwGi,tEL, kwθ1,tP−ELθ1,tH−ELθ1,tG−ELθ2,tP−ELθ2,tH−ELθ2,tG−ELWE11.37906.961256402.82.62.11.01.02.7WE22.33337.342962722.82.62.11.01.02.7WE31.49007.005650802.82.62.11.01.02.7WE41.41386.975267362.82.62.11.01.02.7WE51.76187.114355922.82.62.11.01.02.7



---

*Source: 7717605-2022-07-05.xml*