GPNBI-inspired MOSFA for Pareto operation optimization of integrated energy system

GPNBI-inspired MOSFA for Pareto operation optimization of integrated energy system

Energy Conversion and Management 151 (2017) 524–537 Contents lists available at ScienceDirect Energy Conversion and Management journal homepage: www...

1MB Sizes 0 Downloads 50 Views

Energy Conversion and Management 151 (2017) 524–537

Contents lists available at ScienceDirect

Energy Conversion and Management journal homepage: www.elsevier.com/locate/enconman

GPNBI-inspired MOSFA for Pareto operation optimization of integrated energy system

MARK



Huaizhi Wanga, Rongquan Zhanga, Jianchun Penga, , Guibin Wanga, Yitao Liua, Hui Jiangb, Wenxin Liua a b

College of Mechatronics and Control Engineering, Shenzhen University, Shenzhen 518060, China College of Optoelectronic Engineering, Shenzhen University, Shenzhen 518060, China

A R T I C L E I N F O

A B S T R A C T

Keywords: Integrated energy system Combined cooling heat and power Firefly algorithm Normal boundary intersection Multi-objective optimization

With the rapid growth of penetration of distributed generations and the strong interdependencies of multiple energy sources, the optimal operation of integrated energy systems (IESs) is becoming more important than ever before. However, current optimal operation strategies for IESs may be impractical because their potential negative impacts on utility grids are generally not considered. Therefore, in this study, an original multiobjective optimization model for IES operation is formulated to minimize the operational cost, primary energy consumption, and carbon dioxide emission of IESs and to optimally reduce the power loss and voltage magnitude deviation of the utility grid. Then, a novel Pareto optimization algorithm, called multiobjective strength firefly algorithm (MOSFA), is proposed to solve the multiobjective operation problem of IES. The strength firefly algorithm (SFA) uses a Boltzmann distribution and a chaotic-sequence-based population selection process to facilitate local exploitation and global exploration. In addition, a generalized piecewise normal boundary intersection (GPNBI) method is developed to transform the multiobjective operation problem into a series of highly constrained single-objective optimization sub-problems that can be effectively solved by the SFA. Finally, a hyper-plane-based decision making strategy is introduced to identify the best compromise solution for the obtained Pareto frontiers. The GPNBI-inspired MOSFA was comprehensively evaluated on a novel IES powered by natural gas, wind and photovoltaic generators. The standard IEEE 39-bus system is considered as the utility grid. The numerical results demonstrated that the proposed MOSFA exhibits competitive performance when compared to the algorithms of the state, which means that the optimized operation strategy provides a better tradeoff between all objectives considered in this study.

1. Introduction In the last few decades, global population growth and worldwide industrial development have led to a exponentially surge in the demand of energy [1]. Thereupon, energy crisis and global warming effect are becoming an serious issue required to be addressed urgently. Current energy consumption structure is significantly dependent on fossil fuels [2], which makes the contradiction between energy consumption and environmental protection more and more severe. Therefore, the microgrid, a decentralized combination of interconnected distributed generators and loads, has attracted significant attention in recent years for improving the utilization efficiency of multiple energy sources and alleviating their environmental impact [3]. Traditionally, the microgrid consists of local generation, energy storage, consumption and the point of common coupling. As an independent controllable entity, the

microgrid is capable of operating in both the grid-connected and island modes to provide good solutions for the power supply in the case of power interruptions in the main grid. In recent years, a typical microgrid, termed as integrated energy system (IES) that simultaneously generates electricity, heating, and cooling from the combustion of fuel or solar heat collector, is gradually becoming a hot research topic worldwide [4]. IES is generally applied into large apartment buildings, hotels and stores because of their high requirements of heating and cooling demands. The main advantage of IES is its efficiency, which is higher than that of conventional microgrid. Another obvious advantage of the IES is that renewable energy resources, such as wind and photovoltaic power generation, can be easily embedded into the IES. In general, the IES contains an electricity subsystem and a combined cooling, heating and power (CCHP) subsystem [5]. To effectively dispatch the energy between multiple energy-generation sources, IES



Corresponding author. E-mail addresses: [email protected] (H. Wang), [email protected] (R. Zhang), [email protected] (J. Peng), [email protected] (G. Wang), [email protected] (Y. Liu), [email protected] (H. Jiang). http://dx.doi.org/10.1016/j.enconman.2017.09.005 Received 2 June 2017; Received in revised form 30 August 2017; Accepted 1 September 2017 0196-8904/ © 2017 Elsevier Ltd. All rights reserved.

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

urban CCHP system was designed and multiobjective evolutionary algorithm was applied to create the PFs. In addition, non-dominated sorting genetic algorithm (NSGA) was used to simultaneously optimize multiple primary objectives, such as the life cycle carbon footprint and the life cycle cost of a building refurbishment system [26]. NSGA was also applied for optimization of operational cost and carbon dioxide emission for microgrid energy management system with high uncertainty [27]. Meanwhile, particle swarm optimization (PSO) has also been implemented to optimize various energy management systems, such as hybrid vehicle systems [28] and grid-connected renewable energy systems [29]. Furthermore, multiple group search optimizer (MGSO) has also been proposed recently to deal with multiobjective operational problems for energy management in CCHP [30] and classical electric power systems [31]. The most significant disadvantage of Pareto evolutionary algorithms is that they cannot produce evenly-distributed PFs with high-quality because they may get clustered stochastically. Therefore, these algorithms are far from successful. Comparatively, multilevel programming methods, such as the normal constraint [32] and normal boundary intersection (NBI) methods [33], have the ability to generate a PF with better uniformity by transforming the multiobjective problem into a series of single-objective sub-problems. In [34], a Pareto adaptive penalty-based NBI method was proposed to address multiobjective operation problem in hybrid renewable energy systems. A reduced wide normalized normal constraint method was applied to optimize both the location and capacity of dynamic reactive power sources in energy systems [35]. In [36], a multiobjective operation and maintenance strategy of a microgrid was determined by applying the mesh adaptive direct search method. In [37], a hybrid CCHP system was established and an analytic hierarchy process was introduced to determine the best compromise solution. However, the uniformity of PFs obtained from multilevel programming methods varies a lot with the Pareto surface shape [38]. Thus, the ability of these methods to generate an evenly-distributed PF is largely limited. Extensive overviews regarding operation strategy optimization analysis for IESs and energy systems can be found in [39,40] and [41,42], respectively. From the above analysis, it is clear that operational optimization for IESs and energy management systems has been extensively reported in the literature. However, there are still some problems requiring immediate solutions. Firstly, the microgrid with sustainable energy integration remains an active research topic [43]. In [44], two classical soft-computing techniques, i.e., PSO and ant colony optimization, were introduced to optimally manage the power flow in a microgrid. In [45], the intelligent energy management of a typical microgrid was innovatively formulated. A thorough review on sustainable energy integrated microgrid could be found in [46]. Nevertheless, to date, IESs and CCHP systems are generally powered by natural gas [47], while sustainable energy integrated IESs and CCHP systems attract little attentions. In practice, wind and solar power can be used directly to satisfy the electric demand or transformed to thermal energy by a regenerative electric boiler indirectly to meet the heating and cooling demand. Therefore, with the rapid growth of the penetration of sustainable energy in IESs, there is a pressing need to study how to make full use of sustainable energy to achieve multiobjective energy dispatch tasks among multiple energy sources. Secondly, current strategies for IES operation may be impractical because they may, to some extent, increase the power loss and voltage magnitude deviation of the utility grid. Hence, these negative impacts of IES operation on the utility grid should be considered. Thirdly, so far, multiobjective optimization algorithms applied in IESs have failed to produce evenly-distributed PFs. It has been demonstrated in [38] that evenly-distributed PFs are of great importance for system operators to meet various particular requirements in real-world. Therefore, an advanced multiobjective optimization algorithm that can produce evenly-distributed PFs is highly required. In summary, to overcome the above challenges, this study focuses on investigating multiobjective strategies for IES operation. The main contributions of this paper are as follows:

operation strategies have been extensively researched in recent years, and can be divided into two basic categories [6]: simple-rule based operation strategies and optimal/suboptimal operation strategies. Simple-rule based operation strategies consist of the following electric load (FEL) model and following thermal load (FTL) model [7]. The FEL model requires that the prime mover in electricity subsystem should be loaded according to the electrical demand [8]. The waste heat generated in the loading process is recycled to meet thermal requirements. Auxiliary heating subsystem will not be operated unless the recovered heat could not satisfy the thermal demand. In FTL model, sufficient recovered heat must be provided in priority to meet the heating and cooling requirements [9]. While, the generated electricity in this case may be not adequate for satisfying electric requirement. The overall performance and effectiveness of the two operation strategies have been extensively studied in [10,11]. However, simple-rule based operation strategies fail to guarantee that the IES is always operated in an optimal even suboptimal operation state [12]. Consequently, many optimal and suboptimal operation strategies based on various optimization techniques were proposed innovatively. Optimal/suboptimal operation strategies comprise single-objective and multiobjective optimization models [13]. A single-objective optimization model minimizes the operational cost, primary energy consumption or carbon dioxide emission. In [14], a linear programming model was proposed to optimally reduce the overall operational cost for CCHP systems. It was found that optimal operation is closely related to the fulfillment of the load conditions. In [15], a minimax regret nonlinear model based on cross-entropy was developed to alleviate the impact of disturbing uncertainties on CCHP economics. In [16], a genetic algorithm was applied to investigate the effects of key parameters on the thermodynamic performance of CCHP systems. In [17], two novel economic models for grid-connected and standalone CCHP systems were formulated and solved by using a colonial competitive algorithm. In [18], a two-stage operation strategy based on a genetic algorithm and mixed-integer linear programming was developed to comprehensively evaluate the CCHP economic or environmental performance. However, single-objective optimization models cannot simultaneously address the operational cost, energy consumption and greenhouse gas emission of the IES. Therefore, multiobjective operation strategies have been extensively explored to optimize the tradeoff between the conflicting objectives in CCHP systems and integrated energy management systems. In general, multiobjective optimization methods are divided into three classes: no-preference, priori and posteriori methods [19]. In nopreference methods, such as the method of global criterion, the decision maker does not explicitly articulate any preference information. However, no-preference methods focus on obtaining a neutral compromise solution instead of forming Pareto frontiers (PFs) and are, therefore, rarely applied in multiobjective energy management optimization. In addition, priori methods, such as the weighted sum method and ∊-constraint method, require that adequate preference information should be obtained before the solution process [20]. In [21], the weighted sum method and genetic algorithm were combined to simultaneously optimize multiple CCHP operation objectives. In [22], the multiobjective energy consumption scheduling problem was modeled as a mixed-integer linear programming model with the ∊-constraint method. In [23], a weighted-sum technique based fuzzy satisfying approach was applied to solve the cost-emission problem of a hybrid battery/Photovoltaic (PV)/fuel cell energy system. At last, posteriori methods try to generate all the Pareto solutions in the PFs. In general, posteriori methods can be categorized into Pareto evolutionary methods, in which one run of the algorithm generates a representative set of Pareto solutions, and multilevel programming methods, in which each run of the algorithm only generates one Pareto solution [24]. Pareto evolutionary methods use mechanisms inspired from biological evolution, such as crossover, mutation and selection, to guide the potential candidates towards their optimal states. In [25], a small-scale distributed 525

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

• A new IES structure consisting of an electricity subsystem and a





generated from the heat recovery system, solar heat collector and two boilers. In addition, the compression chiller and electric boiler are designed as backups for providing heating and cooling services. In this study, the IEEE 39-bus system was employed as the utility grid and the proposed IES was connected to bus 16 to evaluate the feasibility and effectiveness of the whole IES system, as shown in Fig. 1. The IEEE system topology and parameters are available in [49]. The IES system proposed in this subsection is inspired from a conventional CCHP system [50] and the IES presented in [5]. Comparatively, the proposed IES has at least three advantages. The first is that the wind and PV energy are considered as distributed generations that reduce carbon dioxide emission and further improve the utilization efficiency of multiple energy sources. The second advantage is that the IES does not consider a steam turbine that uses exhausted high-pressure steam to generate electricity. This is because, in practice, the generated power accounts for a small percentage of the total electric demand. Therefore, the profit would not recover the investment and operational cost. The last advantage is that the electric boiler is only powered by a visual power generation unit (PGU) rather than various distribution generators because the visual PGU can provide more stable and smooth power. More on the relationships between the subsystems and the operation principles of IESs can be found in [1,2,5], and [51].

CCHP subsystem is proposed. Wind turbines, solar photovoltaic cells and gas turbines are integrated in the electricity subsystem. Compared to classical CCHP systems, the CCHP subsystem considered in this work has an extra solar heat collector and a heat storage tank. Both of them are designed to largely improve the energy utilization efficiency. A novel multiobjective optimization model for IES operation is formulated. In this model, two types of objectives are considered for optimization, i.e., conventional IES operation objectives and the negative impacts of IES operation on the utility grid. The former includes the operational cost, primary energy consumption, carbon dioxide emission and the integrated performance criterion (IPC) [1]. The latter consists of power loss and deviations in voltage magnitudes. A multiobjective strength firefly algorithm (MOSFA) inspired by the piecewise normal boundary intersection (GPNBI) is developed to solve the multiobjective optimization model for IES operation. The MOSFA is capable of producing more uniformly distributed PFs compared to Pareto evolutionary algorithms. In addition, a hyperplane-based decision making strategy is proposed to recognize the best-compromise solution for the obtained PFs.

This paper is organized as follows. In Section 2, the proposed IES structure and its corresponding multiobjective optimization model are presented in detail. In Section 3, we describe the GPNBI-inspired MOSFA. Subsequently, the feasibility and effectiveness of the proposed IES structure, multiobjective optimization model, and MOSFA are comprehensively evaluated in Section 4 by using real energy demand data from an office building in Shenzhen, China. Lastly, the conclusions are drawn in Section 5.

2.2. Multiobjective optimization model In general, the conventional optimization goals of IESs and CCHP systems are to minimize the operational cost, primary energy consumption, and carbon dioxide emission. However, conventional optimization models may not provide operation strategies that are suitable for practical implementation because they may deteriorate the economics of the connected utility grid. For example, the operation strategy may increase the power loss or deviation of the bus voltage magnitudes. Therefore, both the power loss and voltage profile objectives are considered in this study for minimizing the impacts of IESs on the utility grid. (1) Operational Cost Objective: Since the financial year of 1990, many countries, such as Finland, Poland, New Zealand, Mexico, and Japan, have begun to implement carbon tax [52]. From 2012, China also plans to collect carbon tax. Therefore, the carbon tax should be considered in the operational cost of IESs without the loss of generality, although the tax amounts to a small percentage of the total operational cost at present. Coupled with the electric power cost and fuel cost, the operational cost of IES can be described as follows,

2. Problem formulation In this section, a novel IES structure consisting of an electricity subsystem and a CCHP subsystem is proposed and presented. In addition, an original multiobjective optimization model of IES operation is formulated. 2.1. IES structure description The IES structure proposed in this paper is composed of an electricity subsystem, which is powered by natural gas, wind, and photovoltaic energy, and a CCHP subsystem. Note that the IES considered here is different from the microgrid presented in [48], although they look similar. Firstly, a battery-based energy storage ensemble system was considered in the microgrid [48] but not considered in the proposed IES. In addition, two thermal storages were used in the microgrid to store hot and cold water. However, only one heat storage tank exists in our IES. Last but not least, wind energy is considered in our IES but not in [48]. The proposed IES structure is illustrated in Fig. 1. In the electricity subsystem, PV cells and wind turbines are used to transform renewable solar and wind energy into electricity. In addition, natural gas is employed as a backup source to fire gas turbines and produce electricity. Then, the electricity from the PV cells, wind turbines, and gas turbines are combined together to satisfy the electric load in the IES. However, owing to the chaotic nature of the weather, wind and solar power always exhibit strong intermittency and randomness. Therefore, the proposed IES is connected to a utility grid that is used as a controllable unit that absorbs extra power from the IES when solar and wind power are in surplus or sends the electricity into the IES to maintain the electric balance in the IES when a power shortage exists. In the CCHP subsystem, the heating and cooling demand are mainly satisfied by the heating coil and absorption chiller, respectively, using the thermal energy from the heat storage system. The heat storage system is used to store the high-pressure steam

f1 = [Egrid Ce + (FAB + FGT )(Cf + uf Cc )−Eex Cs]/ SP1

(1)

where f1 is the operational cost objective, Egrid represents the electricity provided by utility grid, Eex denotes the electricity flowed back to utility grid, FGT and FAB are the fossil fuels required in gas turbine and auxiliary boiler, respectively. The locations of Eex, Egrid, FGT and FAB are labeled as shown Fig. 1. In addition, Ce, Cf, Cc and Cs are the unit price of electricity from utility grid, natural gas, carbon tax and electricity sold back to utility grid, uf is emission conversion factor for natural gas, SP1 represents the cost of separation production system [1]. (2) Primary Energy Consumption Objective: This objective is aimed at minimizing the primary energy consumption and can be formulated as a linear function [53], as follows,

f2 = [Egrid ke + (FAB + FGT ) kf ]/ SP2

(2)

where f2 is primary energy consumption objective, ke and kf are site-toprimary energy conversion factors for electricity from gird and for natural gas, SP2 represents the energy consumption of separation production system. (3) Carbon Dioxide Emission Objective: Owing to the increasing public awareness of global warming and environmental protection, carbon dioxide emission is generally considered as a criterion for effectively reducing pollutant emission, as follows, 526

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

Fig. 1. The IES model with wind and photovoltaic energy.

f3 = [Egrid ue + (FAB + FGT ) uf ]/ SP3

where EWT, EPV and EGT are electric power generated by wind turbine, PV cell and gas turbine, EEB and ECC are power consumptions in electric boiler and compression chiller, EED is the electric demand. EWT, EPV, EGT, EEB and ECC are shown in Fig. 1 and described as follows [3],

(3)

where f3 is carbon dioxide emission objective and ue represents the emission conversion factor for electricity from grid, SP3 represents the emission of separation production system. Obviously, all the units of f1, f2 and f3 should be per unit (p.u.), because the operational cost, primary energy consumption and carbon dioxide emission in IESs are divided by the corresponding cost, energy consumption and emission, respectively, in the separation production system. In addition, compared with f1, f3 seems not to be necessary because the carbon tax is already considered in f1. Nevertheless, up to now, the unit price of the carbon tax is much smaller than that of natural gas. Therefore, the minimization of f1 would have little effect on the minimization of f3. Thus, in this study, carbon dioxide emission is considered as a primary objective to be optimized independently. (4) Power Loss Objective: This objective is to minimize the transmission power loss of the utility grid, subjected to several operation constraints. Newton-Raphson method is generally used to construct the loss objective, as follows,

EWT



Gi [Vi2 + V j2−2Vi Vj cos(δi−δj )]

i=1

(4)

where f4 denotes power loss, NL and Gi represent number of transmission lines and conductance of the ith line, Vi and Vj are the voltage magnitudes at bus i and bus j, respectively. In addition, δi and δj correspond to the voltage angles, respectively. (5) Voltage Profile Objective: This objective aims to minimize the sum of absolute deviations in voltage magnitudes at all system buses. Generally, this objective can be expressed as:



|Vi −Vi,rate |

i=1

(8)

EGT = FGT εηGT / τ

(9)

EEB = QEB / ηEB

(10)

ECC = QCC / ηCC

(11)

(12)

where Qcool is the cooling load in IES, QCD is the thermal energy from absorption chiller and can be calculated as follows,

QCD = {[QEB + Fsl ηsl + FGT ε (ηGTC −ηGT ) ηHRS + FAB ηAB ] ηHST −QHC } ηAC (13)

(5)

where Fsl and ηsl are the energy from solar collector and its corresponding energy conversion factor, ηAB, ηHST and ηAC represent the working efficiency of auxiliary boiler, heat storage tank and absorption chiller, QHC is the thermal energy flowed to heating coil, ηHRS and ηGTC are the recovery efficiency of heat recovery system and the combined heat and power conversion efficiency of wind generator. Obviously, Eq. (13) describes the thermal energy that flows from energy sources to cooling load. (3) Heating Demand Balance: This constraint demands that the thermal energy generated from heating coil should meet the heating load in IES, described as follows,

where NB represents the number of system buses, Vi,rate is the voltage rating at the ith bus. Generally, Vi,rate is set to 1. 2.3. Constraints The constraints required to be considered in IES contain equality constraints that describe the balance between supply and demand, and inequality constraints that specify the acceptable operating ranges for every module. (1) Electric Demand Balance: This constraint specifies that the electric power from utility grid, wind turbine, PV cell and gas turbine should be equal to the sum of electric demand in IES and power consumptions in electric boiler and compression chiller, mathematically described as follows,

Egird + EWT + EPV + EGT = EEB + ECC + EED

EPV = SPV ηPV θ

QCC + QCD = Qcool

NB

f5 =

(7)

where v, vin, vout, and vrat are real-time wind speed, cut-in, cut-out and rated wind speed, respectively. EWT,rat represents rated wind power. In addition, SPV, ηPV and θ denote the area of PV cells, working efficiency and solar irradiation intensity, respectively. QEB and QCC are thermal energy from electric boiler and energy from compression chiller, ηEB and ηCC correspond to their efficiency. Moreover, ɛ, τ and ηGT are the electrical efficiency parameters with respect to gas turbine. (2) Cooling Demand Balance: This constraint requires that the thermal energy from compression chiller and absorption chiller should be equal to the cooling load, as follows,

NL

f4 =

0, 0 ⩽ v ⩽ vin or v>vout ⎧ ⎪E ( )( ), v v v v vin ⩽ v ⩽ vout − − WT , rat in rat in = ⎨E vrat ⩽ v ⩽ vout ⎪ WT ,rat , ⎩

QHC / ηHC = Qheat

(14)

where ηHC represents the working efficiency of heating coil. (4) Operational Range Inequality Constraints: This type of inequality constraints specify that these module variables in (1)–(14) should be

(6) 527

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

initialized for multiobjective optimization because the locations of the optimums are unknown. According to [57], Halton point set is able to produce a uniform distribution over the entire search space and also exhibits a degree of irregularity. Meanwhile, it has been numerically demonstrated that Halton point set may improve the algorithm performance [57]. Therefore, in this paper, Halton point set based initialization process is introduced, described as follows,

within their lower and upper limits. More detail information with respect to the inequality constraints within the modules in IES could be found in [6]. 2.4. Pareto optimization formulation Taking all the five optimization objectives and three types of constraints into consideration, the optimal operation strategy problem can be formulated as a compact multiobjective Pareto optimization model, as follows,

min Xi0,j = Ximin + Haltoni (i,j ) × (Ximax ,j ,j −Xi,j )

i = 1,2,…,Np; j = 1,2,…,D

Ximax ,j

and represent the upper and lower bound for the jth where decision variable in ith firefly, Xi0,j corresponds to the initialized value for the variable. In addition, D is the size of variables, Haltoni(i,j) represents the pseudo-random value obtained from Halton code [57]. (2) Attraction Selection: According to the above analysis, it is obvious that in standard FA, there are too many attractions between fireflies that may lead them to oscillate in the solution space and would thus deteriorate the algorithm performance. Therefore, this study proposes an attraction selection method based on the Boltzmann distribution to reduce the number of attractions. Considering a given firefly j, we assume that Na different attractions are randomly chosen from the attraction pool based on the Boltzmann distribution, as follows,

min[f1 (X ),f2 (X ),…,fNobj (X )] s. t .

⎧ pj (X ) = 0 j = 1,2,…,Neq ⎨ qk (X ) ⩽ 0 k = 1,2,…,Nineq ⎩

(16)

Ximin ,j

(15)

where X represents the set of decision variables, Nobj, Neq and Nineq are the numbers of objectives, equality and inequality constraints. In addition, pj and qk denote the jth equality constraint and kth inequality constraint. The primary objective of the multiobjective model is to find a set of Pareto-optimal solutions that are distributed as evenly as possible so that none of these solutions dominates any others for all objectives in the solution set.

Na

i probBoltz = exp(fopt / θ ·f i )/ ∑ exp(fopt / θ ·f i )

3. GPNBI inspired MOSFA

(17)

i=1

where probi Boltz represents the Boltzmann probability of the attraction connected to the ith firefly, i ∈ Np & i ≠ j, θ is the temperature factor. In addition, f i is the normalized fitness value of the ith firefly, and fopt corresponds to the optimal fitness value that has been found among all the fireflies in the searching process. fopt is determined after the population initialization process and updated in each iteration. From (17), it is obvious that for the given firefly, a smaller fitness value indicates a higher probability. This means that all the fireflies would be more inclined to retain the attractions that are connected to fireflies with better fitness values and discard other attractions. Therefore, it can be concluded that the Boltzmann-distribution-based attraction selection method exhibits a wide diversity range at the beginning of the searching process because the fitness values of most fireflies are similar. Meanwhile, the iterative searching process would gradually reduce any deviations because all the fireflies are forced to move toward the global optimum. These features give the SFA a capability of convergence acceleration. (3) Firefly Update: At each iteration, a firefly is required to move toward a new destination according to its attractions. In this study, a chaotic-theory-based firefly update rule is proposed. A chaotic sequence is employed to guide the fireflies toward the global optimum. In addition, [58] reported that chaotic theory may enhance the algorithm performance. Accordingly, the firefly update rule can be described as follows:

In order to effectively solve the Pareto multiobjective optimization model (15), a hybrid multiobjective algorithm consisting of a strength firefly algorithm (SFA), GPNBI, and a hyper-plane-based decisionmaking strategy is proposed. 3.1. Strength firefly algorithm In nature, fireflies are a family of insects that emit light in the summer in tropical and wet areas as a signal system to attract each other. Inspired by their flashing phenomenon, a new population-based meta-heuristic algorithm, termed as the firefly algorithm (FA), was proposed by Yang in 2008 [54]. For simplicity, the standard FA was formulated by assuming ideal flashing behaviors, namely: (1) All fireflies are unisexual, which means that any individual firefly will be attracted by all other fireflies; (2) the attractiveness of fireflies is proportional to their brightness; thus, less-bright fireflies are forced to move towards brighter fireflies. In addition, a firefly moves randomly if no individual is brighter than the given firefly; and (3) the brightness of fireflies should be associated with an objective function. Up to date, the FA has been used in a variety of practical fields [55]. From these assumptions, it can be seen that the searching behavior of the FA is determined by the attractions between fireflies and each firefly is attracted by all other brighter fireflies. This searching process is termed as the full attraction model. However, it should be noted that too many attractions exist in the full attraction model. As analyzed in [56], the average number of attractions for each firefly is (Np − 1)/2, where Np is the number of fireflies. In reality, these attractions would indeed help the fireflies to find new promising candidates. However, for any problem, too many attractions would more or less make the fireflies, i.e., candidates, oscillate in the solution space, resulting in a higher computational cost, and more importantly, a relatively slow convergence speed. Even worse, the proposed multiobjective optimization model (15) exhibited nonlinearity and was highly constrained, thus, the size of the fireflies in the problem domain was relatively large. This would again hugely increase the attractions between fireflies. Therefore, to effectively surmount this demerit, a new strength FA is proposed to randomly select part of the attractions so as to avoid these unnecessary oscillations during the searching process and to guide the fireflies toward global optimality. (1) Initialization: In general, the candidates in FA should be evenly-

2

k Xik + 1 = Xik + β0 e−γrij τik (XBolt −Xik ) + ψik ζ

(18)

where rij is the distance between Xi and Xj, β0 represents the attractiveness at the distance r = 0, γ is the light absorption coefficient, ζ is a random number between 0 and 1, XkBolt is the candidate that are stochastically selected from (17). In addition, τik is a time series chaotic number, and ψik is a non-monotonic decreasing number. τik and ψik are updated as follows,

τik + 1 = μ (τik )2sin(πτik ) ψik =

(19)

ϑ 1−ζ 1 − (k / NG)

(20) ψik,

where NG is the number of generation, ϑ is the decay rate of μ is a control factor that is used to determine the behavior of chaotic number. At every iteration, the fireflies are regenerated based on the above attraction selection method and the firefly update rule until the algorithm 528

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

f3

⎡ η1 Φ11+…+ηi Φ1i +…+ηm Φ1m ⎤ ⎥ ⎢ ⋮ ⎥ ⎢ P (η1,…,ηi ,…,ηm) = ⎢ η1 Φj1+…+ηi Φji +…+ηm Φjm ⎥ ⎥ ⎢ ⋮ ⎥ ⎢ +…+ +…+ η Φ η Φ η Φ m 1 mi mm i m ⎥ ⎢ 1 ⎦ ⎣

A3

' M 2,3

M2,3

A2

where η is a set of factors that determine the location of the point and m meanwhile satisfy ∑1 ηi = m-1 and 0 ⩽ ηi ⩽ 1, Φji is the normalized value of the jth row and ith column that can be computed based on (21). (2) Piecewise Decomposition: At first, two anchor points on the created normalized utopia plane are randomly chosen and the midpoint of the two anchor points is obtained. Then, the Pareto solution corresponding to the midpoint can be obtained by solving the following single objective NBI problem:

M1,3 ' M1,3

M1,2 mid Eud1,2

f2

' M1,2

A1

f1

min−Eudimid ,j mid

⎧ Φ·ηi,j ⎪ s. t . : ⎨ ⎪ ⎩

Fig. 2. Illustration of the piecewise decomposition process.

converges. Then, the solution with the best fitness value is considered as the optimum for the problem.

⋯ f1 (Xm∗ ) ⎤ ⋮ ⎥ ⎥ ⋯ fi (Xm∗ ) ⎥ ⋱ ⋮ ⎥ ⎥ ⋯ f m∗ (Xm∗ ) ⎥ ⎦

fi (X )−f i∗ (Xi∗) ∗ max[fi (X1 ),…,f i∗ (Xi∗),…,fi (Xm∗ )]−f i∗ (Xi∗)

k = 1,2,…,Neq

qk (X ) ⩽ 0

k = 1,2,…,Nineq

(24)

3.3. Hyper-plane based decision making strategy After obtaining a set of Pareto solutions, the best compromise solution must be identified and implemented by using a decision-making strategy. In general, fuzzy-based decision-making strategy is used to select the best compromise solution from the obtained PFs [59]. However, considering a realistic scenario, such as that the case that the obtained PFs are incomplete and partly known owing to the limited capability of the optimization algorithm, compromise solutions obtained with the fuzzy-based strategy may exhibit high instability, as demonstrated latterly in Section 4.3. Therefore, in order to overcome the disadvantage of high instability, a novel hyper-plane-based decision making strategy is proposed as follows:

(21)

Afterwards, a normalization transformation is employed to normalize the payoff table because the objective functions have different order of magnitudes. The normalization transformation is presented as follows,

fi (X ) =

i≠j

pk (X ) = 0

is a vector of η in which ηi = ηj = 0.5 and others are equal to where 1, n = [n1, n2,…,nm] represents the normal unit vector F (X ) = [f1 (X ),…,fi (X ),…,fm (X )]. In addition, Eudimid is the distance be,j tween the obtained midpoint and the Pareto surface. Repeat the above two steps until all Pareto points with regard to the midpoints of any two anchors have been found. It is clear that there are Cm2 Pareto solutions with regard to Cm2 midpoints when there has m objectives. The entire Pareto surface can be decomposed into (Cm2 + 1) Pareto sub-surfaces by using the newly-generated Cm2 Pareto points and the m original anchor points. Obviously, each of these Pareto sub-surfaces can be further piecewise decomposed by searching the Pareto points with regard to the midpoints of any two anchors on the utopia sub-plane. The piecewise decomposition process can be repeated until enough number of Pareto points have been created. The proposed GPNBI method is graphically illustrated in Fig. 2, in which m is set to 3. Therefore, there are C32 = 3 midpoints, that are M1,2, M1,3 and M2,3. Afterwards, 3 Pareto points corresponding to the ′ , three midpoints can be found by solving (24) and are denoted by M1,2 ⌢ ′ and M2,3 ′ . Then, the original Pareto surface A1 A2 A3 can be decomM1,3 ⌢ ′ M1,3 ′ , posed into C32+1 = 4 Pareto sub-surfaces, that are A1 M1,2 ⌢ ⌢ ⌢ ′ ′ ′ ′ ′ ′ ′ A2 M1,2 M2,3 , A3 M1,3 M2,3 and M ,12 M1,3 M2,3 . It can be easily concluded that each of the created four sub-surfaces can be further piecewise decomposed. Therefore, from above analysis, the multiobjective optimization problem can be piecewise transformed into a series of single objective NBI sub-problems and then effectively addressed by using the proposed SFA algorithm in Section 3.1.

The primary objective of multi-objective optimization is to generate a family of Pareto solutions in which none of these solutions can outperform others for all objectives. Meanwhile, these optimal solutions are required to be evenly-distributed on the Pareto profile so as to meet particular requirements in energy dispatch [38]. However, the solutions obtained from Pareto evolutionary algorithms may readily get clustered and the solutions from conventional multilevel programming methods are largely dependent on the real Pareto surface to the problem being optimized. Consequently, the uniformity of the PFs obtained from these two types of methods is hugely limited. Therefore, in order to produce a PF as uniform as possible, a GPNBI method inspired from piecewise idea [38] is originally proposed. In this method, normalized payoff table is firstly required to be formed. Then, a Pareto solution producing process based on GPNBI is presented, described as follows. (1) Normalized Payoff Table: Regarding a multiobjective optimization problem with m objectives, the individual optimum fi∗(Xi∗) (1 ≤ i ≤ m) for the ith objective can be computed. Meanwhile, the values of other objectives with Xi∗ can be assessed and denoted by f1(Xi∗),…, fi-1(Xi∗), fi+1(Xi∗),…, fm(Xi∗). Then, we get the payoff table Φ, described as follows,

⋯ f1 (Xi∗) ⋱ ⋮ ⋯ f i∗ (Xi∗) ⋮ ⋯ fm (Xi∗)

+ Eudimid ,j n = F (X )

ηimid ,j

3.2. Generalized piecewise normal boundary intersection

∗ ∗ ⎡ f1 (X1 ) ⎢ ⋮ ⎢ ∗ Φ = ⎢ fi (X1 ) ⎢ ⋮ ⎢ ∗ ⎢ fm (X1 ) ⎣

(23)

(22)

Apparently, with the normalization transformation, the multiobjective problem (15) can be solved in a non-dimensional space. In this space, any normalized point P(η1,…,ηi,…,ηm) can be expressed as follows:

Nobj

max

∑ i=1

529

j

j

j

φ

φ

φ

i )) Eud ((f 1 ,…,f i ,…,f Nobj ),(f 1 i ,…,f i i ,…,f Nobj

(25)

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

Fig. 3. Flowchart of GPNBI inspired MOSFA.

Start Randomly select a midpoint of two anchor points and determine its coordinates ƐƚĂďůŝƐŚĂ'WE/ŵŝŶŝŵŝnjĂƟŽŶƐƵďͲƉƌŽďůĞŵďĂƐĞĚŽŶ;ϮϰͿ Initialization (16) EŽ EŽ

zĞƐ

i
j
Attraction selection (17) Firefly update (18) YzĞƐ EŽ

ůůŽĨƚŚĞWĂƌĞƚŽƐŽůƵƟŽŶƐ ŽŶĐƵƌƌĞŶƚƐƵďͲƐƵƌĨĂĐĞ ůĂLJĞƌŚĂǀĞďĞĞŶĂĐŚŝĞǀĞĚ͍

zĞƐ

ŶŽƵŐŚWĂƌĞƚŽƐŽůƵƟŽŶƐ ŚĂǀĞďĞĞŶŐĞŶĞƌĂƚĞĚ͍

zĞƐ



Piecewise decompose current Pareto-surface into a set of Pareto sub-surfaces Hyper-plane based decision-making strategy End

ŽŶǀĞƌŐĞŶĐĞĞƌƌŽƌŽĨ ŽƉĞƌĂƟŽŶĂůĐŽƐƚ;ΨͬŚŽƵƌͿ

Fig. 4. The plot of convergence error for operational cost.

/ƚĞƌĂƟŽŶƐ

ŽŶǀĞƌŐĞŶĐĞĞƌƌŽƌŽĨƉƌŝŵĂƌLJ ĞŶĞƌŐLJĐŽŶƐƵŵƉƟŽŶ;ŬtŚͬŚŽƵƌͿ

Fig. 5. The plot of convergence error for primary energy consumption.

/ƚĞƌĂƟŽŶƐ

530

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

ŽŶǀĞƌŐĞŶĐĞĞƌƌŽƌŽĨĐĂƌďŽŶ ĚŝŽdžŝĚĞĞŵŝƐƐŝŽŶ;dŽŶƐͬŚŽƵƌͿ

Fig. 6. The plot of convergence error for carbon dioxide emission.

/ƚĞƌĂƟŽŶƐ Convergence error of primary energy conƐƵmpƟon ;ŬtŚͬŚoƵrͿ

Fig. 7. Plot of statistical convergence errors of primary energy consumption.

^imƵůaƟon ƟmeƐ

Convergence error of carďon ĚiodžiĚe emiƐƐion ;donƐͬŚoƵrͿ

Fig. 8. Plot of statistical convergence errors of carbon dioxide emission.

^imƵůaƟon ƟmeƐ

Convergence error of operaƟonaů coƐƚ ;ΨͬŚoƵrͿ

Fig. 9. Plot of statistical convergence errors of operational cost.

^imƵůaƟon ƟmeƐ

j

f i = fi j /(max(fi1 ,…fi j ,…,f iNPF ))

maximize the distance between the solution and the hyper-plane consisting of all the least-wanted-Pareto-solutions found so far for each objective. The benefits of the hyper-plane based decision making strategy can be found in tion 4.3.

(26)

where NPF represents the number of Pareto solutions, Eud means Euclidean distance, φi is the index of a solution with the worst fitness on the ith objective. In addition, fij denotes the fitness of the jth PF point for ith objective. From (25), it can be seen that the hyper-plane based decision making strategy is designed to obtain a compromise solution that can

3.4. Solution framework The flow chart of the GPNBI inspired MOSFA is illustrated in Fig. 3. 531

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

including a dual-objective investigation and a tri-objective investigation.

Table 1 Daily operational cost, primary energy consumption and carbon dioxide emissions. Time

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24

Real operational cost ($)

Real primary energy consumption (kWh)

Real carbon dioxide emission (Ton)

IES

CCHP

IES

CCHP

IES

CCHP

29.84 28.89 28.39 27.94 25.82 14.72 15.36 18.92 20.84 22.52 24.03 26.43 26.81 27.16 29.83 32.51 30.63 30.89 29.96 27.04 23.86 35.05 34.47 32.97

40.55 31.94 40.36 32.96 32.87 19.17 19.62 22.49 24.70 85.21 100.20 68.64 94.96 91.60 69.26 56.79 34.86 35.19 60.31 61.52 61.57 38.91 39.51 36.46

74.18 72.38 71.23 70.18 69.12 67.86 69.12 74.64 77.40 81.70 84.55 87.98 89.44 90.59 94.23 97.87 94.38 95.17 93.87 89.14 84.71 83.31 81.90 79.78

79.31 77.40 76.78 75.26 74.48 72.53 74.66 79.35 82.94 89.94 94.49 99.12 101.9 101.9 102.3 104.0 99.50 100.0 99.46 94.91 90.59 88.43 86.92 84.77

15.54 15.16 14.92 14.70 14.47 14.20 14.47 15.64 16.21 17.11 17.72 18.43 18.74 18.98 19.75 20.52 19.78 19.94 19.67 18.68 17.75 17.45 17.16 16.72

16.63 16.21 16.11 15.77 15.60 15.19 15.63 16.62 17.37 18.95 19.93 20.84 21.49 21.47 21.52 21.85 20.85 20.97 20.89 19.95 19.04 18.53 18.21 17.76

4.1. Experimental settings The test system used in this study is a hybrid of the proposed IES and the IEEE 39-bus system. The topology of the hybrid test system is shown in Fig. 1. The IEEE 39-bus system consisting of 10 generators and 46 branches is used as the utility grid, and the entire IES is connected to bus 16. The model parameters of the IEEE 39-bus system are presented in [49]. The optimization objectives are the operational cost, primary energy consumption, carbon dioxide emission of the IES, the power loss and voltage magnitude deviation of the 39-bus system, as defined in (1)-(5). The constraints of the multiobjective optimization model are detailed in Section 2.3. The IES model parameters are mainly obtained from [1,7] and [53]. Specifically, the unit price of natural gas, i.e., Cf, is set to 3.23$/kWh and the unit price of the carbon tax is set to 0.00002$/kWh. It is clear that the price of the carbon tax is much smaller than that of natural gas. In addition, the unit price of the electricity sold back to the utility grid is set to 0$/kWh [1,53]. This is because the feedback of the electric power from the microgrid and IES to the utility grid is not allowed in many regional distribution systems, such as Sanmenxia in Henan province, China. Furthermore, the unit price of electricity from the unity grid is set to 0.93$ at 6:00am-21:00 pm and 0.55$ at other times. The emission conversion factors ue and uf are set to 968 g/kWh and 220 g/ kWh, respectively. The PEC conversion factor of the utility grid ke is set to 3.336 p.u., and that of natural gas kf is set to 1.047 kWh/m3. The daily electricity demand EED, heating demand Qheat and cooling demand Qcool are estimated by using EnergyPlus software. Moreover, the working efficiency coefficients ηEB, ηCC, ηHC, ηAB, ηAC, ηGTC, ηGT, ηHST, ηHRS and ηsl are set to 3.0, 0.45, 3.2, 3.5, 0.5, 0.9, 0.8, 0.98, 0.97, and 0.13, respectively. The real-time wind speed data and the related parameters, such as the cut-in and cut-out wind speed, can be found in [60,61]. Other parameters related to solar power generation and the gas turbine generator are reported in [7]. The eight parameters, θ, β0, γ, μ, ϑ, NG, Na and Np are crucial in the implementation of the GPNBI inspired SFA for IES operation and should be set according to the following guidelines. The temperature of the Boltzmann distribution, i.e., θ, determines how the probabilities of attraction selection are assigned. A larger θ tends to allocate the same probability to each attraction and a smaller θ, such as one close to 0, would tend to assign a high probability to attractions with better fitness values. The primary purpose for applying the Boltzmann distribution is to select promising attractions among the attraction pool. Therefore, a smaller θ is preferable. Our simulations show that a value of θ ranging from 0.05 to 0.4 is acceptable, and here it is set to 0.15. In addition, the attractiveness at distance = 0, i.e., β0, is the initial time step for candidate updating. Our simulations show that a value in the range of 0.1–2 works well because a larger value would affect the stability of the algorithm and a smaller value would slow down the convergence speed. Here, β0 is set to 1. The light absorption coefficient γ characterizes the variation of the attractiveness and thus plays a crucial role in determining the convergence speed and the behavior of the SFA. Based on the analysis in [62], γ is set to 1 in this paper. In addition, the factor μ is used to maintain the updating sequence chaotic and fixed at 0.7 [58]. The decay rate ϑ is used as a coefficient to control the magnitude of the exploration component in (18). Generally, a value within the range from 0.1 to 1 is acceptable. In this study, it is set to 0.9. Furthermore, Np is the number of fireflies. Generally, Np should be within 20–40 when the problem to be optimized is relatively simple and 60–100 when the problem is complex and highly-constrained. A smaller Np would result in an algorithm lacking diversity and a larger Np would degrade the algorithm convergence. Therefore, Np is set to 60 based on this analysis. The number of candidates Na is employed to determine the scale of the attractions. Theoretically, Na can be any value from 1 to Np. However,

Table 2 Statistical span metric of various evolutionary algorithms.

Min Max Average Standard variance

Proposed algorithm (p.u.)

MGSO (p.u.)

MOEA (p.u.)

NSGAII (p.u.)

1.3340 1.4142 1.3661 0.0245

1.1567 1.4790 1.3179 0.2278

0.6372 0.9815 0.8629 0.1955

1.1924 1.2271 1.2098 0.0439

Table 3 Statistical convergence metric of various evolutionary algorithms.

Min Max Average Standard variance

Proposed algorithm (p.u.)

MGSO (p.u.)

MOEA (p.u.)

NSGAII (p.u.)

0.0203 0.0251 0.0232 4.1558E-4

0.0244 0.0349 0.0297 7.0451E-4

0.0391 0.1280 0.0773 1.3781E-3

0.0482 0.2136 0.1530 4.0695E-3

Table 4 Statistical spacing metric of various evolutionary algorithms.

Min Max Average Standard variance

Proposed algorithm (p.u.)

MGSO (p.u.)

MOEA (p.u.)

NSGAII (p.u.)

0.0051 0.0106 0.0067 7.69E-05

0.0084 0.0257 0.0169 5.39E-04

0.0103 0.0311 0.0179 2.97E-03

0.0083 0.0392 0.0269 2.36E-03

4. Case studies The feasibility of the IES structure, the multiobjective optimization model, and the overall performance of the GPNBI-inspired MOSFA were comprehensively tested and evaluated using the following case studies, 532

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

Power loss (p.u.)

Fig. 10. The Pareto profile for carbon dioxide emission and power loss.

Carbon dioxide emission (p.u.)

Power loss (p.u.)

Fig. 11. The Pareto profile for integrated performance criterion and power loss.

IPC (p.u.)

solƚage magniƚude deviaƟon (p.u.)

Fig. 12. The Pareto profile for primary energy consumption and voltage magnitude deviation.

Primary energy consumpƟon (p.u.)

4.2. Investigations of dual-objective optimization

as discussed in Section 3.1, the number of attractions should be limited. Therefore, Na is set to 4. Lastly, the number of generations NG is set to 70 because the algorithm converged after 70 iterations. Namely, the moving average convergence error is below 10−6. It is required to be noted that some of these parameters, including θ, Na, μ and ϑ, are not sensitive to the optimization problems. In other words, these four parameters do not need to be tuned again for different problems or different constraint sets. The other four parameters, β0, γ, NG and Np, may require to be varied for changing problems or different scenarios. The above analysis can be used as guidelines for the determination of these parameters. The settings of these parameters and coefficients have been heuristically adjusted by a large amount of comparative simulations and studies. Meanwhile, these parameters remain unchanged in all the following case studies, except for a special setting. The proposed algorithm is implemented in MATLAB R2014a and carried out on a personal computer with an Intel(R) Core i5-4590 3.3-GHz CPU and 16.00 GB of RAM.

(1) Case 1: This case is aimed at evaluating the overall performance of the SFA in terms of optimization objectives. This performance is measured by convergence error, and comprehensively compared to classical differential evolution (DE), PSO, and the standard firefly algorithm (FA). The convergence errors of the operational cost, primary energy consumption and carbon dioxide emission for various algorithms are presented in Figs. 4–6. Apparently, after about 50 iterations, the convergence errors of SFA for all three objectives approach to the value of zero and errors of PSO, FA and DE for each objective always exist. This means SFA has the ability to find much better global optimums than PSO, FA and DE. It also can be seen that the convergence errors of SFA falls very quickly, indicating a fast convergence rate. This is due to the adopted Boltzmann distribution and chaotic sequence that provide a general balance between global exploration and local exploitation. Especially note that the standard FA converges much slower

533

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

Table 5 Daily gas fuels of three extreme solutions in resultant Pareto profile. Gas fuel (MWh) time (h)

Proposed IES

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 Average

Table 6 Statistics of the four metrics for case 2.

Conventional CCHP

Index

Method

Best

Worst

Average

Proposed MGSO Proposed MGSO Proposed MGSO Proposed MGSO

0.0755 0.0613 16.8826 25.6920 0.0715 0.0506 0.1199 0.2362

0.0394 0.0220 69.4581 101.1362 0.1419 0.1207 0.2980 0.5594

0.0603 0.0481 39.1630 62.9144 0.0973 0.0810 0.1671 0.3783

Power loss

VMD

IPC

Power loss

VMD

IPC

Area (p.u.2)

6.819 6.579 6.515 6.349 6.380 6.234 6.327 6.816 7.157 7.586 7.736 8.181 8.273 8.385 8.537 8.668 8.693 8.651 8.592 8.203 7.753 7.573 7.533 7.300 7.535

6.775 6.494 6.416 6.369 6.179 6.174 6.323 6.896 7.021 7.426 7.603 7.973 8.134 8.211 8.476 8.845 8.598 8.592 8.596 7.993 7.654 7.458 7.423 7.340 7.457

6.595 6.415 6.468 6.377 6.287 6.158 6.334 6.616 7.086 7.385 7.615 7.892 8.096 8.061 8.464 8.740 8.601 8.566 8.554 8.032 7.519 7.444 7.424 7.080 7.409

7.048 7.006 6.913 6.848 6.664 6.623 6.636 7.158 7.500 7.909 8.369 8.841 9.083 8.984 9.060 9.170 8.831 8.956 8.868 8.471 8.039 7.891 7.876 7.563 7.929

7.034 6.837 6.874 6.850 6.704 6.610 6.721 7.246 7.384 7.905 8.320 8.769 8.948 9.065 9.069 9.218 8.890 8.977 8.796 8.412 8.035 7.823 7.794 7.582 7.911

7.113 6.861 6.876 6.742 6.703 6.514 6.560 7.121 7.402 8.021 8.301 8.753 8.936 8.975 9.135 9.264 8.871 8.996 8.742 8.349 8.079 7.907 7.756 7.587 7.899

lmax/lmin (p.u.) l (p.u.) σ (p.u.)

convergence errors for each objective are considered to be 100% converged and the biggest errors are considered to be 0% converged. Then, with respect to primary energy consumption, the SFA algorithm is averagely 99.99% converged with a variance of 1.31E-9. The average convergence of PSO, FA and DE are 97.89% with a variance of 4.29E-4, 95.44% with a variance of 1.72E-5 and 16.67% with a variance of 1.90E-2. Considering carbon dioxide emission, the convergence of SFA has a minimal of 99.97% and a maximum of 100% with an average of 99.99% and a variance of 4.56E-09. The convergences of PSO, FA and DE are ranged from 94.32% to 99.91%, from 95.04% to 96.17% and from 0% to 48.62%, respectively. Regarding operation cost, the convergences of SFA, PSO, FA and DE are varied from 99.99% to 100%, from 96.52% to 98.55%, from 97.12% to 98.30%, and from 0% to 87.60%, respectively. From these statistical results, it is clear that the proposed SFA algorithm exhibits the lowest convergence errors and lowest variances when compared to PSO, FA and DE. Therefore, these results demonstrate that the proposed SFA can not only significantly improve the algorithm performance, but also exhibit high stability and robustness. (2) Case 2: To illustrate the feasibility and efficiency of the proposed IES system, a series of simulations under the proposed operation strategy are carried out according to a representative days’ electric load, thermal and cooling demands from [53]. The daily operational cost, primary energy consumption and carbon dioxide emission are presented in Table 1 and compared to standard CCHP system [1]. It can be seen that the real operational cost of the IES system ranges from 14.72$ to 35.05$, with an average of 26.87$ and a standard variance of 5.42$. Apparently, the real operation cost of the IES is more economic than the CCHP benchmark. In addition, the real primary energy consumption of the IES ranges from 67.86kWh to 97.87kWh, with an average of 82.28kWh. While, the average of the primary energy consumption for CCHP is 88.79kWh. Furthermore, the averages of carbon dioxide emission for IES and CCHP are 17.24 Ton and 18.64 Ton, respectively. From these numerical results, it is clear that the IES is operated in a higher efficiency with lower operational cost, energy consumption and emission. This can be explained by two factors. One is that wind and PV generators are integrated in the proposed IES system

than SFA. This is because the attractions in SFA are hugely reduced by using Boltzmann distribution, i.e., 4 out of 59 attractions are only retained. The reduction in the number of attractions would avoid the fireflies, or candidates, being oscillated during the searching process. In other words, the Boltzmann distribution is used to retain the promising attractions in attraction pool and discard other attractions. Each more attraction means that the search direction for each firefly would be changed once and cause the firefly to be oscillated in the solution space. And, therefore, the SFA exhibits a faster convergence rate than standard FA. In addition, it can be seen that DE has the worst convergence property with bigger errors and not smoothness. This is because DE is found to be stage in the fine-tuning process [63]. Generally, the performance of heuristic optimization algorithms, such as the FA, is known to change from one simulation to another owing to the variability of the model parameters. Therefore, twenty independent runs for all algorithms are carried out to demonstrate the relatively high-stability of the proposed SFA. The algorithm parameters of FA, PSO and DE are well-tuned based on the guidelines in [62,64] and [65]. These parameters remain unchanged in the following simulation environment. For each optimization objective, Figs. 7–9 present the plots of statistical convergence errors at 20 trials. In order to better explain the meaning expressed in Figs. 7–9, we assume that the lowest

Fig. 13. The 3-dimentional Pareto profile considering power loss, IPC and VMD.

Integrated Performance Criterion

0.91 0.905 0.9 0.895 0.89 0.885 0.88 0.45

0.77 0.78 0.5

0.55

0.6

0.79

Powe

0.65

r Loss

0.7

0.8 0.75

0.8

0.81

Ma ltage

Vo

534

oĮle

de Pr

gnitu

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

Table 7 The compromise solutions from various decision making strategies. PF size

Algorithms

Fuzzy decision making

Nash decision making

Proposed decision making

Power loss (p.u.)

IPC (p.u.)

VMD (p.u.)

Power loss (p.u.)

IPC (p.u.)

VMD (p.u.)

Power loss (p.u.)

IPC (p.u.)

VMD (p.u.)

45

MOSFA MGSO

0.7135 0.6970

0.9021 0.9177

0.7983 0.8006

0.7054 0.6970

0.9014 0.9177

0.7978 0.8006

0.6821 0.6607

0.8995 0.8957

0.7971 0.7980

153

MOSFA MGSO

0.6352 0.6785

0.8952 0.8930

0.7916 0.7928

0.6295 0.6472

0.8906 0.8963

0.7954 0.7951

0.6405 0.6258

0.8959 0.8933

0.7903 0.7892

set to 60 and 160, respectively. Other parameters are set based on the analysis in Section 4.1. From Fig. 10, it is clear that the two endpoints of the whole Pareto profile are solid red circle. This means the length between the two red endpoints obtained from GPNBI inspired MOSFA is the longest among the four benchmark algorithms in Fig. 10. Because the physical meaning of span metric is the distance between two normalized endpoints, it can be concluded that the proposed GPNBI inspired MOSFA algorithm has the best span metric among the four algorithms. In addition, it is obvious that the Pareto points from MOEA are concentrated in the middle of the whole Pareto profile. Meanwhile, most of the Pareto points from NSGAII are scattered on the right of the whole Pareto profile, and most Pareto points from MGSO are distributed on the left of the whole Pareto profile. Comparatively, the solid red circles are evenly-distributed on the whole Pareto profile, which means the proposed algorithm exhibits a better uniformity, i.e., spacing metric. Moreover, convergence metric is computed as the sum of the distance between the Pareto point in the algorithm and corresponding points in Pareto reference. Therefore, we can conclude from the distributions of Pareto points of the four algorithms that the GPNBI inspired MOSFA algorithm has the best convergence metric. The same conclusions could also be obtained from the Pareto profiles in Figs. 11 and 12. Therefore, the proposed algorithm performs the best among the four algorithms in terms of span, spacing and convergence metrics. This confirms the high potential of the proposed algorithm to find Pareto optimality.

and the electrical requirement can thus be partly provided by wind and solar energy. The other factor is the IES structure consisting of a solar heat collector, a heat recovery system and a heat storage tank. These elements collaborate with each other to satisfy the cooling and heating demands, resulting in a higher efficiency. The two factors correspond to the reasons why the proposed IES system offers a higher working efficiency. (3) Case 3: This case is to evaluate the overall performance of the proposed GPNBI inspired MOSFA in terms of widely-used span, spacing and convergence metrics [31]. The span metric is computed as the normalized distance of two terminals in the obtained Pareto profile. Convergence metric is applied to assess the closeness between the obtained PF and PF reference. The spacing metric is used to evaluate the crowding distance of adjacent Pareto solutions. Note that all the span, convergence and spacing metrics are per unit values. In addition, power loss of the utility grid and operational cost of the IES are considered as two optimization objectives. The MGSO [31], non-sorted genetic algorithm II (NSGAII) [66] and multiobjective evolutionary algorithm (MOEA) [67] are used for performance comparison. The parameters in MGSO, NSGAII and MOEA are well-tuned in a trial-and-error manner based on the parameter adjustment guidelines in [31,66,67]. For fair comparison, the population size Np and total number of generations of GPNBI inspired MOSFA, NSGAII, MGSO and MOEA are all set to 60 and 160, respectively. The size of PF is set to 17. Moreover, for each algorithm, ten independent runs are carried out because the performance of evolutionary algorithms is generally known to change from one simulation to another. The span, spacing and convergence metrics are statistically presented in Table 2–4, respectively. From Table 2, it can be seen that the span metric of GPNBI inspired MOSFA ranges from 1.3340 p.u. to 1.4142 p.u., with an average of 1.3661 p.u. and a standard variance of 0.0245 per unit. It is evident that GPNBI inspired MOSFA performs the best when compared to MGSO, NSGAII and MOEA. The results indicate that the proposed algorithm can effectively find the non-dominated solutions in the feasible solution space. From Table 3, it is obvious that the proposed algorithm outperform MGSO, NSGAII and MOEA in term of convergence metric. This means that the proposed algorithm has the best searching ability for Pareto solutions. In addition, from Table 4, it can be seen that the spacing metric from GPNBI inspired MOSFA performs better than the spacing metrics from MGSO, NSGAII and MOEA. This means that the variance of the distance of adjacent PF solutions obtained with the proposed GPNBI inspired MOSFA is minimum and the PF solutions are thus more uniformlydistributed. Furthermore, in order to comprehensively evaluate the overall performance of the proposed GPNBI inspired MOSFA, three other simulations were carried out using different combinations of two objectives. The results of the utility grid power loss and carbon dioxide emission are graphically presented in Fig. 10. The Pareto profile of the power loss and IPC [1], which is an index that simultaneously addresses operational cost, primary energy consumption and carbon dioxide emission, are illustrated in Fig. 11. The Pareto profile in combination with voltage magnitude deviation (VMD) and primary energy consumption are presented in Fig. 12. The population size Np and total number of generations NG of all algorithms used for comparison are also

4.3. Investigations on tri-objective optimization (1) Case 1: The purpose of this case is to test the feasibility and high efficiency of the proposed IES system. The IPC, power loss and VMD are chosen as the three optimization objectives. A standard CCHP system is used for performance comparison. The number of PF points is set to 153. The total number of generations NG is set to 300 that is bigger than the numbers adopted in dual-objective optimization. Table 5 presents the gas fuels of the three extreme solutions for the proposed IES system and the standard CCHP system over a typical day. From the results, it can be seen that the gas fuel of the IES under the power loss objective ranged from 6.234 MWh to 8.693 MWh with an average of 7.535 MWh. Compared to standard CCHP, the gas fuel consumption is reduced by 4.97% on average. Considering the VMD and IPC, the gas fuels are evenly reduced by 5.74% and 6.20%, respectively. Therefore, it can be concluded that the proposed IES system exhibits a high working efficiency compared to the standard CCHP system. (2) Case 2: This case is to comprehensively assess the overall performance of the proposed multiobjective optimization algorithm. MGSO is selected as the only benchmark algorithm because of its comparatively better performance in Section 4.2. In this case, the uniformity is evaluated by four indices [38], i.e., lmax/lmin, l, σ, and area. The index lmax/lmin denotes the ratio of maximum distance of adjacent Pareto solutions to the minimum distance. The indices l and σ represent the mean and standard deviation of these adjacent distances. Area is computed as the normalized area of the three endpoints in the resultant Pareto surface. These four indices are used to systematically measure the overall performance of the GPNBI inspired MOSFA algorithm. The 535

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

higher stability and robustness. The hyper-plane based decision-making takes the least-wanted solutions as the baseline, then measures the distance between the solutions in PF and the hyper-plane, and so the solution with the largest distance is considered as the best compromise solution. This kind of measures of distance endow the proposed decision-making strategy with robustness and high-stability.

three optimization objectives are the same as case 1. The number of generations NG is set to 300. In addition, ten independent runs of GPNBI inspired MOSFA and MGSO are carried out respectively to calculate the four indices. The number of Pareto optimal points is fixed to 153. The Pareto profile obtained from the proposed algorithm is graphically plotted in Fig. 13. The four indices, i.e., lmax/lmin, l, σ, and area, are tabulated in Table 6, respectively. From Table 6, it can be seen that the area of the proposed algorithm has a minimum of 0.0394 p.u.2 to 0.0755 p.u.2 with an average of 0.0603 p.u.2, while, the area of the MGSO is from 0.0220 p.u.2 to 0.0613 p.u.2. Evidently, the proposed algorithm has a better normalized area. Therefore, it can be concluded that the proposed multiobjective algorithm has the ability to find solutions with much better fitness for all three objectives compared to MGSO. In addition, it is clear that the proposed algorithm has a minimum of lmax/lmin and also a minimal of σ, respectively. This means the normalized distance between adjacent PF points obtained with the proposed algorithm shows little variation. Thus, the PF solutions are more uniformly-distributed. The uniformlydistributed PF profile plays a crucial role for dispatchers of energy systems because any particular requirement in multiobjective operation strategy would be met. Moreover, the proposed algorithm exhibits a larger l. This could be explained by the fact that the mean value of normalized adjacent distance increases as the area expands when the number of the PF solutions is fixed. From these numerical results and analysis, it is evident that the proposed algorithm performs better than the MGSO. This confirms a high potential of the proposed algorithm in searching global Pareto optimality. (3) Case 3: This case is aimed at demonstrating the feasibility and stability of the hyper-plane based decision-making strategy. The Pareto frontiers are repeatedly produced with different number of Pareto optimal points based on the presentations in case 2. For ease of analysis, the numbers of Pareto optimal points are set to 45 and 153, respectively. The PFs with 45 Pareto points mean these PFs are incomplete and partly known. On the contrary, the PFs with 153 points indicate that these PFs are complete and thorough. The two kinds of PFs are used to test the feasibility and stability of decision-making approaches. The numerical results from the produced PF samples are detailed in Table 7. The fuzzy decision-making strategy and Nash-based decisionmaking approach are used for comparison. It can be seen from Table 7 that the power loss of the compromise solutions obtained from the fuzzy decision-making strategy varies from 0.6352 per unit to 0.7135 per unit, IPC varies from 0.8930 per unit to 0.9177 per unit, and VMD varies from 0.7983 per unit to 0.8006 per unit, respectively. Considering Nash-based decision making, the power loss, IPC and VMD of the compromise solution vary from 0.6295 per unit to 0.7054 per unit, from 0.8906 per unit to 0.9177 per unit, and from 0.7951 per unit to 0.8006 per unit, respectively. However, in the hyper-plane based decision-making, the power loss, IPC and VMD range from 0.6258 per unit to 0.6821 per unit, from 0.8933 per unit to 0.8995 per unit, and from 0.7892 per unit to 0.7980 per unit, respectively. Evidently, the proposed hyper-plane based decision-making strategy performs best in terms of power loss, IPC and VMD. Therefore, it can be concluded that the proposed decision-making strategy is a better choice for energy dispatchers from an viewpoint of overall performance. In addition, with respect to the proposed MOSFA algorithm, the power loss obtained from the proposed decision-making strategy varies from 0.6821 per unit to 0.6405 per unit as the number of Pareto points is changed from 45 to 153, i.e., from an incomplete PF to a complete PF. Under the same condition, the power loss from fuzzy varies from 0.7135 per unit to 0.6352 per unit and loss from Nash varies from 0.7054 per unit to 0.6295 per unit. Hence, the change of power loss based on hyper-plane decision-making strategy is minimum. From Table 7, it can also be seen that in most cases, the changes of compromise solutions from hyper-plane based decision making are smaller than that from fuzzy decision making and Nash decision making. These results demonstrate that the proposed decision-making strategy exhibits

5. Conclusions In this paper, a new IES structure consisting of electricity subsystem with sustainable energy resources and CCHP is firstly proposed to achieve a high energy utilization efficiency. Then, a multiobjective optimization model for IES operation is formulized. Both conventional IES operation objectives and the impacts of IES on utility grid are considered. At last, A hybrid multiobjective algorithm is originally developed to effectively solve the multiobjective optimization model. The algorithm is a hybrid of SFA, GPNBI and hyper-plane based decision-making strategy. The feasibility and efficiency of the proposed IES system, multiobjective model and hybrid algorithm are comprehensively evaluated in the case studies where IEEE 39-bus system is considered as the utility grid. The numerical results show that the IES system exhibit a higher operation efficiency than conventional CCHP systems. Meanwhile, the proposed MOSFA algorithm inspired from GPNBI is demonstrated to exhibit better PF profile in terms of several metrics when compared to the frequently-used benchmarks in literature, which could benefit the dispatchers of energy systems for a specific use. It also demonstrates that the performance of the proposed algorithm is very stable and robust. The superiority of the proposed optimization algorithm is attributed to the techniques that used, including chaotic sequence, Boltzmann distribution, generalized piecewise decomposition and the hyper-plane based decision-making. It is, therefore, apparent that the proposed algorithm has a high potential for practical implementations in electrical power and energy systems. Acknowledgement The authors would like to thank the anonymous reviewers for their insightful and detailed comments, as these comments led us to an improvement of this paper. This work was jointly supported by Natural Science Foundation of China (51707123, 51507103, 51477104), Natural Science Foundation of Guangdong Province (2017A030310061, 2016A030313041), and the Foundations of Shenzhen Science and Technology Committee (JCYJ20170302153607971, JCYJ20160520164316902). References [1] Fang F, Wang QH, Shi Y. A novel optimal operational strategy for the CCHP system based on two operating modes. IEEE Trans Power Syst 2012;27(2):1032–41. [2] Hashemi R. A developed offline model for optimal operation of combined heating and cooling and power systems. IEEE Trans Energy Convers 2009;24(1):222–9. [3] Jin D, et al. Toward a cyber resilient and secure microgrid using software-defined networking. IEEE Trans Smart Grid 2017;8(5):2494–504. [4] Gu Wei, Wang Jun, Lu Shuai, Luo Zhao, Wu Chenyu. Optimal operation for integrated energy system considering thermal inertia of district heating network and buildings. Appl Energy 2017;199:234–46. [5] Ju L, Tan Z, Li H, Tan Q, Yu X, Song X. Multi-objective operation optimization and evaluation model for CCHP and renewable energy based hybrid energy system driven by distributed energy resources in China. Energy 2016;111(15):322–40. [6] Hu M, Cho H. A probability constrained multi-objective optimization model for CCHP system operation decision support. Appl Energy 2014;116(1):230–42. [7] Zhu Q, Luo X, Zhang B, Chen Y. Mathematical modeling and optimization of a largescale combined cooling, heat, and power system that incorporates unit changeover and time-of-use electricity price. Energy Convers Manage 2017;133(1):385–98. [8] Wang J, Lu Y, Yang Y, Mao T. Thermodynamic performance analysis and optimization of a solar-assisted combined cooling, heating and power system. Energy 2016;115(Part 1):49–59. [9] Li M, Mu H, Li N, Ma B. Optimal design and operation strategy for integrated evaluation of CCHP (combined cooling heating and power) system. Energy 2016;99(15):202–20. [10] Jalalzadeh-Azar A. A comparison of electrical- and thermal-load-following CHP systems. ASHRAE Trans 2004;110(2):85–94.

536

Energy Conversion and Management 151 (2017) 524–537

H. Wang et al.

[39] Cho Heejin, Smith Amanda D, Mago Pedro. Combined cooling, heating and power: a review of performance improvement and optimization. Appl Energy 2014;136(31):168–85. [40] Wei Gu, Wu Zhi, Bo Rui, Liu Wei, Zhou Gan, Chen Wu, Wu Zaijun. Modeling, planning and optimal energy management of combined cooling, heating and power microgrid: a review. Int J Electr Power Energy Syst 2014;54:26–37. [41] Fadaee M, Radzi MAM. Multi-objective optimization of a stand-alone hybrid renewable energy system by using evolutionary algorithms: a review. Renew Sustain Energy Rev 2012;16(5):3364–9. [42] Hina Fathima A, Palanisamy K. Optimization in microgrids with hybrid energy systems – a review. Renew Sustain Energy Rev 2015;45:431–46. [43] Moghaddam Amjad Anvari, Seifi Alireza, Niknam Taher. Multi-operation management of a typical micro-grids using Particle Swarm Optimization: a comparative study. Renew Sustain Energy Rev 2012;16(2):1268–81. [44] Colson CM, Nehrir MH, Pourmousavi SA. Towards real-time microgrid power management using computational intelligence methods. IEEE PES General Meeting: Minneapolis, MN; 2010. p. 1–8. [45] Chaouachi A, Kamel RM, Andoulsi R, Nagasaka K. Multiobjective intelligent energy management for a microgrid. IEEE Trans Industr Electron 2013;60(4):1688–99. [46] Justo Jackson John, Mwasilu Francis, Lee Ju, Jung Jin-Woo. AC-microgrids versus DC-microgrids with distributed energy resources: a review. Renew Sustain Energy Rev 2013;24:387–405. [47] Chang Huawei, Wan Zhongmin, Zheng Yao, Chen Xi, Shu Shuiming, Tu Zhengkai, Chan Siew Hwa. Energy analysis of a hybrid PEMFC–solar energy residential microCCHP system combined with an organic Rankine cycle and vapor compression cycle. Energy Convers Manage 2017;142(15):374–84. [48] Zhang X, Sharma R, He Yanyi. Optimal energy management of a rural microgrid system using multi-objective optimization. 2012 IEEE PES Innovative Smart Grid Technologies (ISGT), Washington, DC 2012;2012:1–8. [49] Zimmerman R, Gan D. MATPOWER: A Matlab Power System Simulation Package. [Online]. Available: < http://www.pserc.cornell.edu.matpower > . [50] Mohammadi Amin, Ahmadi Mohammad H, Bidi Mokhtar, Joda Fatemeh, Valero Antonio, Uson Sergio. Exergy analysis of a combined cooling, heating and power system integrated with wind turbine and compressed air energy storage system. Energy Convers Manage 2017;131(1):69–78. [51] Vögelin Philipp, Georges Gil, Boulouchos Konstatinos. Design analysis of gas engine combined heat and power plants (CHP) for building and industry heat demand under varying price structures. Energy 2017;125(15):356–66. [52] World Bank. State and trends of carbon pricing 2016, Washington, D.C., World Bank Group, available at: < http://documents.worldbank.org/curated/en/ 598811476464765822/State-and-trends-of-carbon-pricing > . [53] Fu Xueqian, Huang Shangyuan, Li Rui, Sun Hongbin. Electric power output optimization for CCHP using PSO theory. Energy Procedia 2016;103:9–14. [54] Ibrahim Ibrahim Anwar, Khatib Tamer. A novel hybrid model for hourly global solar radiation prediction using random forests technique and firefly algorithm. Energy Convers Manage 2017;138(15):413–25. [55] Wang Deyun, Luo Hongyuan, Grunder Olivier, Lin Yanbing, Guo Haixiang. Multistep ahead electricity price forecasting using a hybrid model based on two-layer decomposition technique and BP neural network optimized by firefly algorithm. Appl Energy 2017;190(15):390–407. [56] Wang Hui, Wang Wenjun, Zhou Xinyu, Sun Hui, Zhao Jia, Yu Xiang, et al. Firefly algorithm with neighborhood attraction. Inform Sci 2017;382–383:374–87. [57] Price K, Storn RM, Lampinen JA. Differential evolution: a practical approach to global optimization (Natural Computing Series). New York, NY, USA: SpringerVerlag; 2005. [58] Caponetto R, Fortuna L, Fazzino S, Xibilia MG. Chaotic sequences to improve the performance of evolutionary algorithms. IEEE Trans Evol Comput 2003;7(3):289–304. [59] Abido MA. Multi-objective evolutionary algorithms for electric power dispatch problem. IEEE Trans Evol Comput 2006;10(3):315–29. [60] Wang Huai-zhi, Li Gang-qiang, Wang Gui-bin, Peng Jian-chun, Jiang Hui, Liu Yitao. Deep learning based ensemble approach for probabilistic wind power forecasting. Appl Energy 2017;188(15):56–70. [61] Wang HZ, Wang GB, Li GQ, Peng JC, Liu YT. Deep belief network based deterministic and probabilistic wind speed forecasting approach. Appl Energy 2016;182(15):80–93. [62] Gandomi AH, Yang X-S, Talatahari S, Alavi AH. Firefly algorithm with chaos. Commun Nonlinear Sci Numer Simul 2013;18(1):89–98. [63] Adarsh BR, Raghunathan T, Jayabarathi T, Yang X. Economic dispatch using chaotic bat algorithm. Energy 2016;96:666–75. [64] Shayeghi H, Pirayeshnegab A, Jalili A, Shayanfar HA. Application of PSO technique for GEP in restructured power systems. Energy Convers Manage 2009;50(9):2127–35. [65] Basu M. Quasi-oppositional differential evolution for optimal reactive power dispatch. Int J Electr Power Energy Syst 2016;78:29–40. [66] Lotfan S, Akbarpour Ghiasi R, Fallah M, Sadeghi MH. ANN-based modeling and reducing dual-fuel engine’s challenging emissions by multi-objective evolutionary algorithm NSGA-II. Appl Energy 2016;175:91–9. [67] Muhsen Dhiaa Halboot, Ghazali Abu Bakar, Khatib Tamer, Abed Issa Ahmed, Natsheh Emad M. Sizing of a standalone photovoltaic water pumping system using a multi-objective evolutionary algorithm. Energy 2016;109:961–73.

[11] Cardona E, Piacentino A, Cardona F. Matching economical, energetic and environmental benefits: an analysis for hybrid CHCP-heat pump systems. Energy Convers Manage 2006;47(20):3530–42. [12] Yao E, Wang H, Wang L, Xi G, Maréchal F. Multi-objective optimization and exergoeconomic analysis of a combined cooling, heating and power based compressed air energy storage system. Energy Convers Manage 2017;138(15):199–209. [13] Wei Dajun, Chen Alian, Sun Bo, Zhang Chenghui. Multi-objective optimal operation and energy coupling analysis of combined cooling and heating system. Energy 2016;98(1):296–307. [14] Kong XQ, Wang RZ, Huang XH. Energy optimization model for a CCHP system with available gas turbines. Appl Therm Eng 2005;25(2–3):377–91. [15] Wang L, Li Q, Sun M, Wang G. Robust optimisation scheduling of CCHP systems with multi-energy based on minimax regret criterion. IET Gener Transm Distrib 2016;10(9):2194–201. [16] Xia Jiaxi, Wang Jiangfeng, Lou Juwei, Zhao Pan, Dai Yiping. Thermo-economic analysis and optimization of a combined cooling and power (CCP) system for engine waste heat recovery. Energy Convers Manage 2016;128(15):303–16. [17] Karami H, Sanjari MJ, Gooi HB, Gharehpetian GB, Guerrero JM. Stochastic analysis of residential micro combined heat and power system. Energy Convers Manage 2017;138(15):190–8. [18] Guo Li, Liu Wenjian, Cai Jiejin, Hong Bowen, Wang Chengshan. A two-stage optimal planning and design method for combined cooling, heat and power microgrid system. Energy Convers Manage 2013;74:433–45. [19] Marler RT, Arora JS. Survey of multi-objective optimization methods for engineering. Struct Multidiscip Opt 2004;26(6):369–95. [20] Miettinen KM. Nonlinear multiobjective optimization, international series in operations research & management science. Kluwer Academic Publishers; 1999. p. 29–43. [21] Zeng Rong, Li Hongqiang, Jiang Runhua, Liu Lifang, Zhang Guoqiang. A novel multi-objective optimization method for CCHP–GSHP coupling systems. Energy Build 2016;112(15):149–58. [22] Zhang Di, Evangelisti Sara, Lettieri Paola, Papageorgiou Lazaros G. Economic and environmental scheduling of smart homes with microgrid: DER operation and electrical tasks. Energy Convers Manage 2016;110(15):113–24. [23] Majidi Majid, Nojavan Sayyad, Esfetanaj Naser Nourani, Zare Kazem, Zare Kazem. A multi-objective model for optimal operation of a battery/PV/fuel cell/grid hybrid energy system using weighted sum technique and fuzzy satisfying approach considering responsible load management. Solar Energy 2017;144:79–89. [24] Giagkiozis I, Fleming PJ. Methods for multi-objective optimization: an analysis. Inf Sci 2015;293:338–50. [25] Abdollahi Gh, Meratizaman M. Multi-objective approach in thermoenvironomic optimization of a small-scale distributed CCHP system with risk analysis. Energy Build 2011;43(11):3144–53. [26] Schwartz Yair, Raslan Rokia, Mumovic Dejan. Implementing multi objective genetic algorithm for life cycle carbon footprint and life cycle cost minimisation: a building refurbishment case study. Energy 2016;97:58–68. [27] Kuznetsova Elizaveta, Ruiz Carlos, Li Y-Fu, Zio Enrico. Analysis of robust optimization for decentralized microgrid energy management under uncertainty. Int J Electr Power Energy Syst 2015;64:815–32. [28] João P, Antunes Carlos Henggeler. A comparative analysis of meta-heuristic methods for power management of a dual energy storage system for electric vehicles. Energy Convers Manage 2015;95:281–96. [29] Hemmati Reza, Saboori Hedayat, Jirdehi Mehdi Ahmadi. Stochastic planning and scheduling of energy storage systems for congestion management in electric power systems including renewable energy resources. Energy 2017;133:380–7. [30] Zheng JH, Wu QH, Jing ZX. Coordinated scheduling strategy to optimize conflicting benefits for daily operation of integrated electricity and gas networks. Appl Energy 2017;192(15):370–81. [31] Zhou B, Chan KW, Yu T, Chung CY. Equilibrium-inspired multiple group search optimization with synergistic learning for multiobjective electric power dispatch. IEEE Trans Power Syst 2013;28(4):3534–45. [32] Hennen Maike, Voll Philip, Bardow André. An adaptive normal constraint method for bi-objective optimal synthesis of energy systems. Comput Aid Chem Eng 2014;33:1279–84. [33] Roman C, Rosehart W. Evenly distributed Pareto points in multi-objective optimal power flow. IEEE Trans Power Syst 2006;21(2):1011–2. [34] Ming Mengjun, Wang Rui, Zha Yabing, Zhang Tao. Pareto adaptive penalty-based boundary intersection method for multi-objective optimization. Inf Sci 2017;414:158–74. [35] Deng Z, Liu M, Ouyang Y, Lin S, Xie M. Multi-objective mixed-integer dynamic optimization method applied to optimal allocation of dynamic var sources of power systems. In: IEEE Transactions on Power Systems, vol. PP, no. 99, pp. 1–1. [36] Mohamed Faisal A, Koivo Heikki N. System modelling and online optimal management of MicroGrid using Mesh Adaptive Direct Search. Int J Electr Power Energy Syst 2010;32(5):398–407. [37] Yousefi Hossein, Ghodusinejad Mohammad Hasan, Noorollahi Younes. GA/AHPbased optimal design of a hybrid CCHP system considering economy, energy and emission. Energy Build 2017;138(1):309–17. [38] Li Q, Liu M, Liu H. Piecewise normalized normal constraint method applied to minimization of voltage deviation and active power loss in an AC–DC hybrid power system. IEEE Trans Power Syst 2015;30(3):1243–51.

537