# Modelling of Hydrothermal Unit Commitment Coordination Using Efficient Metaheuristic Algorithm: A Hybridized Approach

Authors: Suman Sutradhar; Nalin B. Dev Choudhury; Nidul Sinha
Journal: Journal of Optimization
 (2016)

---

## Abstract

In this paper, a novel approach of hybridization of two efficient metaheuristic algorithms is proposed for energy system analysis and modelling based on a hydro and thermal based power system in both single and multiobjective environment. The scheduling of hydro and thermal power is modelled descriptively including the handling method of various practical nonlinear constraints. The main goal for the proposed modelling is to minimize the total production cost (which is highly nonlinear and nonconvex problem) and emission while satisfying involved hydro and thermal unit commitment limitations. The cascaded hydro reservoirs of hydro subsystem and intertemporal constraints regarding thermal units along with nonlinear nonconvex, mixed-integer mixed-binary objective function make the search space highly complex. To solve such a complicated system, a hybridization of Gray Wolf Optimization and Artificial Bee Colony algorithm, that is,h-ABC/GWO, is used for better exploration and exploitation in the multidimensional search space. Two different test systems are used for modelling and analysis. Experimental results demonstrate the superior performance of the proposed algorithm as compared to other recently reported ones in terms of convergence and better quality of solutions.


---

## Body

## 1. Introduction

Hydrothermal unit commitment coordination (HTUC) problem concerns utilizing hydro potential satisfying hydro constraints in such a way that the cost of produced electricity from thermal resources during a scheduling period of time is lowermost and with lowest environmental impact. Many researchers’ have forwarded various solution techniques to solve such problem but due to extreme nonlinearity and compounded structure of hydro and thermal units, the process of finding feasible and optimal solution is quite tough and time-consuming. In addition, nonlinearity of search space increases greatly with increasing complexity of the compounded structure of the system. Most conventional optimization techniques do not work well for problems like HTUC with such nonlinearity and multimodality. The concerned problem involves not only nonlinearity but also nonconvexity. The thermal power generation (due to quadratic cost function), hydro power generation, and emission are nonlinear in nature. In case of thermal power generation the nonconvexity is due to consideration of valve point loading effect for thermal units. In Figure1, the difference between without valve point loading effect and with valve point loading effect on thermal units on production cost is shown. The effect of valve point loading effect is quite considerable on optimization aspect as the problem becomes highly nonconvex in nature. Also due to reliability and security requirements, hydro and thermal constraints such as load balance, generation limits, and water discharge, starting and ending storage volumes of water, intertemporal limitations, and ramp rate limitations have to be included. The limited storage capability of water reservoirs, the stochastic nature of available water, and different limitations regarding commitment status of thermal units make the problem more difficult for hydrothermal systems. In addition, the environmental impacts are also a great issue considering thermal power plant which bound the authors to include the emission as a second objective function, increasing system complication even more.Figure 1
The waveform of production cost of two thermal units (a) with valve point loading effect and (b) without valve point loading effect wherex- and y-axis represent the thermal power generation of respective thermal units and z-axis represents the total production cost.
(a) (b)Several researchers have attempted solving HTUC problem using different solution procedures. Initially most researchers proposed conventional mathematical techniques, such as dynamic programming [1, 2], decomposition techniques [3], and concept of nonlinear network flow [4, 5] to solve such a complex problem. The development and improvement of different soft computing techniques and their reported impressive performance made the researchers incline towards the use of soft computing techniques to handle such highly nonlinear complex problems. Those techniques, such as Genetic Algorithm [6, 7], Simulated Annealing [8], Evolutionary Programming (EP) [9–12], Bender Decomposition [13–15], Harmony Search Algorithm [16], Particle Swarm Optimization (PSO) [17, 18], Evolution Strategy [19], Differential Evolution (DE) algorithm [20, 21], and mixed-integer programming (MIP) [22], and their modifications are also used to solve HTUC problems as these methods are reported to be efficient in finding near global optimal solutions. HTUC problem can be divided into two subproblems, hydro subproblem and unit commitment coordination subproblem. In literature [23] a detailed survey on unit commitment problem is described.Catalão et al. solved hydro subsystem using a nonlinear optimization method [24] and thermal subsystem using ε-constraining method [25]. Sinha et al. [11] solved a test system, initially proposed in [6] with valve point loading effect using fast EP and demonstrated better results. Basu [12] proposed interactive fuzzy satisfying method based on EP technique to solve HTUC with four hydro and three thermal units with valve point loading effect and emission function coefficients. Modified DE is used to solve the HTUC problems [20]. But this modification lacks enough diversity to provide better solution quality in search space. Rebennack [14] solves the problem by a decomposition algorithm, based on Benders decomposition. Although the problems considered here are nonlinear and various nonlinear constraints are accounted for in the formulation, linear approximation and convexification approaches are used to simplify the problem which in turn is bounded to be erroneous up to some extent. Afkousi-Paqaleh et al. [16] solved unit commitment using Harmony Search Algorithm. Mandal and Chakraborty [21] proposed a new test system for HTUC with hydro subsystem and ten thermal plants with quadratic cost function including valve point loading effects. They have solved the problem using DE technique. However, the proposed method failed to provide superior solution quality for larger system. Hota et al. [17] solved the HTUC problem using an improved PSO. Farhat and El-Hawary [26] present a survey of literature on various optimization techniques applied on HTUC. Recently Zhang et al. [27] proposed a new solution technique known as Small Population based PSO (SPPSO) technique. A new algorithm has been found known as CS-MODE [28] where chaotic sequence is used to increase diversity but the proposed method helps increase local optima entrapment by forcefully directing the swarms to another position without exploring properly. Most of these techniques consume significant amount of time with vast multidimensional search space and many researchers have not considered many reliability limitations, necessary to design a practical system. Recently Predator Prey algorithm [29] is evolved which increases the diversity using accelerating operation but with a cost of depth search. Recently Norouzi et al. [30] have solved HTUC problem using Lexicographic Optimization technique with ε-constraint. Clonal Selection Algorithm (CSA) [31] is recently proposed to demonstrate lesser computational efforts with HTUC. Here, the concept of higher population is used for better convergence but there is a strong possibility of local optima entrapment. It is noticed from the recent publications that these improved heuristic algorithms are able to search better solution from a nonlinear, complex search space. But with greater variable and complicated structure, the behavior of search space becomes more nonlinear [31] with higher chances of premature convergence.Therefore, the objective of this research focuses on developing a suitable algorithm for HTUC problem based on hybridized Artificial Bee Colony (ABC) and Gray Wolf Optimization (GWO) algorithm (h-ABC/GWO), initially developed by Mirjalili et al. [32, 33]. The standard GWO is based on the hunting behavior of gray wolves [32, 33]. In general the gray wolves are considered at the top of the food chain. The gray wolves mostly hunt in groups. Their hunting instinct is mostly divided into three steps. Firstly they create a hierarchy among themselves; secondly they encircle the prey; and thirdly they hunt. In mathematical point of view, the exploitation (which mainly tries to improve its position in terms of optimization or local search) is very much familiar with this hunting instinct of gray wolves. The social hierarchy provides the most efficient wolf in the group which can be termed as the closest individual of the population to the global best. The encircling provides the best individual to understand the prey’s position after evaluating its neighborhood. This phenomena provides the algorithm great exploitation or local search capability. Inh-ABC/GWO algorithm, the onlooker bee phase is replaced by this GWO phase as the onlooker bee phase of ABC algorithm is not very efficient against the problems like HTUC. As for the global search, the initial population and the employed bee phase provide sufficient diversity in the population to provide required global search. Also the scout bee phase is used for global search. In this phase all the solutions that are not improved with increase in iteration are used to explore the search space globally until they find another better position in the search space. This exploration capability is very effective for a complicated and vast problem like HTUC. Moreover, the amount of global and local search in the search space can easily be controlled by a parameter used in the algorithm according to the nature of search space of the problem. Overall these processes provide a great balance between global and local search or in other words enhanced exploration and exploitation capability.On the other hand, in GWO algorithm alone the global-or-local search operation on individual position reduces the exploration capability of GWO, which makes GWO fail to provide sufficient diversity in the search space for a complex problem like HTUC with a further deficiency regarding local optima entrapment or premature convergence. For this reason in this research theh-ABC/GWO is used where two algorithms (i.e., ABC and GWO) are merged together to design a more efficient algorithm to solve HTUC problem.The performance of the proposedh-ABC/GWO approach is demonstrated on two test systems and compared with GWO [32, 33], ABC [34], CSA [31], and SPPSO [27] and recent state-of-the-art mixed-integer nonlinear programming (MINLP) [32]. From the computational results, it is observed that the proposed approach yields better convergence than other mentioned algorithms.


## 2. Problem Formulation

### 2.1. Objective Function

HTUC involves optimization of nonlinear objective function, with a mixture of linear/nonlinear and dynamic network flow constraints. The problem difficulty is compounded by a number of practical considerations and unless several simplifying assumptions are made, the problem is difficult to solve for practical power systems. The basic optimal HTUC involves minimizing overall thermal cost,F, and total emission,E over a given scheduling period of time,T,which can be described as(1a)Minimize,F=∑t=1T∑i=1NSfiPSi,t+SD+ST,(1b)Minimize,EΩ=∑t=1T∑i=1NSαi+βi×PSi,t+γi×PSi,t2+ηi×exp⁡δi×PSti,t,(2)fiPSi,t=ai·PSi,t2+bi·PSi,t+ci+di·sin⁡ei·PSi,tmin-PSi,t,where ai, bi, ci, di, and ei are coefficients for thermal power generation and ST and SD represent the startup and shutdown cost. αi,βi,γi,ηi,δi represent the emission coefficient of ith unit.In Figure2, the nature of production cost and emission based on thermal power generation for a particular unit is demonstrated. From the figure it is clear that both the waveforms are highly nonlinear in nature and the thermal production cost is not only nonlinear but nonconvex as well.Figure 2
The waveform of production cost (see (1b)) and generated emission (see (2)) for a particular thermal unit.In this paper, one day is considered as the scheduling period (T), consisting of 24 discretized intervals of 1 hour. Concerning the thermal units, shutdown cost, SD is quite low in comparison with startup cost (ST) and fuel cost given in (1a), so SD is neglected in this work. Along with quadric cost function the valve point loading effect [12] is also considered as the second part of (2) suggests, which incorporates nonconvexity into the search space of HTUC along with nonlinearity.(a) The Continuity Equation for Reservoir Network. The storage volume of ith reservoir at the end of interval t can be calculated using continuity equation [27](3)Vhi,t-1=Vhi,t+Ihi,t-Qhi,t-Shi,t+∑m=1RuiQhm,t-τm+Shm,t-τm,i=1,2,…,NS,t=1,2,…,T.(b) Physical Limitations on Reservoir Storage Volume and Discharge Rate. The water content of a reservoir has to be within its maximum and minimum reservoir limit, which is represented by (4). The water discharges of hydro units also have to be within its maximum and minimum limit [35](4)Vhmini,t≤Vhi,t≤Vhmaxi,t,i=1,2,…,Nh,t=1,2,…,T,(5)Qhmini,t≤Qhi,t≤Qhmaxi,t,i=1,2,…,Nh,t=1,2,…,T.(c) Initial and Final Reservoir Storage Volume. Water content of each reservoir at initial and end interval of the scheduled period is prespecified [5](6)Vhi,tt=0=Vhbegini,Vhi,tt=T=Vhendi,i=1,2,…,Nh,where Vhbegini and Vhendi are initial and final reservoir volume.(d) Generator Capacity. The power generation of thermal and hydro plants has to be within its maximum and minimum generation limit [6](7)Phmini,t≤Phi,t≤Phmaxi,t,i=1,2,…,Nh,t=1,2,…,T,Ii,tPSmini,t≤Ii,tPSi,t≤Ii,tPSmaxi,ti=1,2,…,NS,t=1,2,…,T.(e) Hydro Power Generation Equation. The hydro power generation is generally a concave nonlinear function in nature [36]. In Figure 3 the waveform of hydro discharge to power generation of a hydro unit for different water volume or water content is shown. From Figure 3 it is evident that the power generation from hydro units is a function of both hydro discharge and water volume. For this reason in this research this quadratic concave function is represented by (8) [11, 27, 37](8)Phi,t=c1·Vhi,t2+c2·Qhi,t2+c3·Vhi,t·Qhi,t+c4·Vhi,t+c5·Qhi,t+c6,i=1,2,…,Nh,t=1,2,…,T,where, c1,c2,c3,c4,c5, and c6 are power generation coefficients of respective hydro units.Figure 3
The waveform of hydro power generation cost (see (8)) with respect to hydro discharges for a particular thermal unit. Hydro discharge Qhi,t.(f) Ramp Rate Limitation. The power generation on a particular thermal unit cannot be increased or decreased beyond ramp up or ramp down limit [23](9)PSi,t-PSi,t-1≤URi,PSi,t-1-PSi,t≤DRi,i=1,2,…,NS,t=1,2,…,T.(g) Minimum Up/Down Time Limitation [23]. If a unit is forced to change its state, it cannot change its state again for a certain amount of time, called minimum up/down time limit(10)
X
O
N
i
,
t
-
1
-
T
O
N
i
·
I
i
,
t
-
1
-
I
i
,
t
≥
0
,
X
O
F
F
i
,
t
-
1
-
T
O
F
F
i
·
I
i
,
t
-
I
i
,
t
-
1
≥
0
,
i
=
1,2
,
…
,
N
S
,
t
=
1,2
,
…
,
T
.(h) System Load Balance. The overall power demand at every hour has to be provided by available thermal and hydro power plants [6](11)∑i=1NsPSi,t+∑i=1NhPhi,t=PDi,tt=1,2,…,T.(i) Startup Cost. Depending on the duration of being off before start, the execution of cold start cost CSi or hot start cost HSi is decided for each start of an individual unit. If a unit starts after a certain period of time, TCOLDi [23], cold start cost CSi is considered or otherwise hot start cost HSi is taken into account for calculation. The summation of all startup costs of all units (see (13)) is incorporated in the cost function given by (1a)(12)SCi,t=CSi×Ii,t-Ii,t-1·Ii,tifXOFFi,t-TCOLDi≥0,SCi,t=HSi×Ii,t-Ii,t-1·Ii,tifXOFFi,t-TCOLDi<0,(13)ST=∑t=1T∑i=1NSSCi,t.


## 2.1. Objective Function

HTUC involves optimization of nonlinear objective function, with a mixture of linear/nonlinear and dynamic network flow constraints. The problem difficulty is compounded by a number of practical considerations and unless several simplifying assumptions are made, the problem is difficult to solve for practical power systems. The basic optimal HTUC involves minimizing overall thermal cost,F, and total emission,E over a given scheduling period of time,T,which can be described as(1a)Minimize,F=∑t=1T∑i=1NSfiPSi,t+SD+ST,(1b)Minimize,EΩ=∑t=1T∑i=1NSαi+βi×PSi,t+γi×PSi,t2+ηi×exp⁡δi×PSti,t,(2)fiPSi,t=ai·PSi,t2+bi·PSi,t+ci+di·sin⁡ei·PSi,tmin-PSi,t,where ai, bi, ci, di, and ei are coefficients for thermal power generation and ST and SD represent the startup and shutdown cost. αi,βi,γi,ηi,δi represent the emission coefficient of ith unit.In Figure2, the nature of production cost and emission based on thermal power generation for a particular unit is demonstrated. From the figure it is clear that both the waveforms are highly nonlinear in nature and the thermal production cost is not only nonlinear but nonconvex as well.Figure 2
The waveform of production cost (see (1b)) and generated emission (see (2)) for a particular thermal unit.In this paper, one day is considered as the scheduling period (T), consisting of 24 discretized intervals of 1 hour. Concerning the thermal units, shutdown cost, SD is quite low in comparison with startup cost (ST) and fuel cost given in (1a), so SD is neglected in this work. Along with quadric cost function the valve point loading effect [12] is also considered as the second part of (2) suggests, which incorporates nonconvexity into the search space of HTUC along with nonlinearity.(a) The Continuity Equation for Reservoir Network. The storage volume of ith reservoir at the end of interval t can be calculated using continuity equation [27](3)Vhi,t-1=Vhi,t+Ihi,t-Qhi,t-Shi,t+∑m=1RuiQhm,t-τm+Shm,t-τm,i=1,2,…,NS,t=1,2,…,T.(b) Physical Limitations on Reservoir Storage Volume and Discharge Rate. The water content of a reservoir has to be within its maximum and minimum reservoir limit, which is represented by (4). The water discharges of hydro units also have to be within its maximum and minimum limit [35](4)Vhmini,t≤Vhi,t≤Vhmaxi,t,i=1,2,…,Nh,t=1,2,…,T,(5)Qhmini,t≤Qhi,t≤Qhmaxi,t,i=1,2,…,Nh,t=1,2,…,T.(c) Initial and Final Reservoir Storage Volume. Water content of each reservoir at initial and end interval of the scheduled period is prespecified [5](6)Vhi,tt=0=Vhbegini,Vhi,tt=T=Vhendi,i=1,2,…,Nh,where Vhbegini and Vhendi are initial and final reservoir volume.(d) Generator Capacity. The power generation of thermal and hydro plants has to be within its maximum and minimum generation limit [6](7)Phmini,t≤Phi,t≤Phmaxi,t,i=1,2,…,Nh,t=1,2,…,T,Ii,tPSmini,t≤Ii,tPSi,t≤Ii,tPSmaxi,ti=1,2,…,NS,t=1,2,…,T.(e) Hydro Power Generation Equation. The hydro power generation is generally a concave nonlinear function in nature [36]. In Figure 3 the waveform of hydro discharge to power generation of a hydro unit for different water volume or water content is shown. From Figure 3 it is evident that the power generation from hydro units is a function of both hydro discharge and water volume. For this reason in this research this quadratic concave function is represented by (8) [11, 27, 37](8)Phi,t=c1·Vhi,t2+c2·Qhi,t2+c3·Vhi,t·Qhi,t+c4·Vhi,t+c5·Qhi,t+c6,i=1,2,…,Nh,t=1,2,…,T,where, c1,c2,c3,c4,c5, and c6 are power generation coefficients of respective hydro units.Figure 3
The waveform of hydro power generation cost (see (8)) with respect to hydro discharges for a particular thermal unit. Hydro discharge Qhi,t.(f) Ramp Rate Limitation. The power generation on a particular thermal unit cannot be increased or decreased beyond ramp up or ramp down limit [23](9)PSi,t-PSi,t-1≤URi,PSi,t-1-PSi,t≤DRi,i=1,2,…,NS,t=1,2,…,T.(g) Minimum Up/Down Time Limitation [23]. If a unit is forced to change its state, it cannot change its state again for a certain amount of time, called minimum up/down time limit(10)
X
O
N
i
,
t
-
1
-
T
O
N
i
·
I
i
,
t
-
1
-
I
i
,
t
≥
0
,
X
O
F
F
i
,
t
-
1
-
T
O
F
F
i
·
I
i
,
t
-
I
i
,
t
-
1
≥
0
,
i
=
1,2
,
…
,
N
S
,
t
=
1,2
,
…
,
T
.(h) System Load Balance. The overall power demand at every hour has to be provided by available thermal and hydro power plants [6](11)∑i=1NsPSi,t+∑i=1NhPhi,t=PDi,tt=1,2,…,T.(i) Startup Cost. Depending on the duration of being off before start, the execution of cold start cost CSi or hot start cost HSi is decided for each start of an individual unit. If a unit starts after a certain period of time, TCOLDi [23], cold start cost CSi is considered or otherwise hot start cost HSi is taken into account for calculation. The summation of all startup costs of all units (see (13)) is incorporated in the cost function given by (1a)(12)SCi,t=CSi×Ii,t-Ii,t-1·Ii,tifXOFFi,t-TCOLDi≥0,SCi,t=HSi×Ii,t-Ii,t-1·Ii,tifXOFFi,t-TCOLDi<0,(13)ST=∑t=1T∑i=1NSSCi,t.


## 3. Hybrid-Artificial Bee Colony (ABC)/Grey Wolf Optimization (GWO) Algorithm

The main drawback regarding metaheuristic algorithms is extensive random search where there is a high possibility of searching the same position repeatedly while few places in the vast multidimensional search space did not get any attention at all. Superior nonlinearity and nonconvexity increase the effect of this phenomena even further. In order to search such a complex search space it is necessary to balance between the local and global search capability of the acting algorithm or in other words the exploitation and exploration capability according to the search space behavior. For this purpose a hybridized algorithm, that is, hybrid-Artificial Bee Colony/Grey Wolf Optimization (h-ABC/GWO) algorithm, is adopted. In this algorithm the superior exploitation capability of GWO algorithm is merged with diversified exploration capability of ABC algorithm.The challenge of an algorithm in constrained environment is much more complicated than in an unconstrained environment. Algorithm performance also greatly depends on the technique used to handle the acting constraints. Before explaining the function evaluation method and constraint handling process of HTUC problem, the key arrangement ofh-ABC/GWO algorithm is discussed in this section.


### 3.1. Initialization

The initial population will be created using (14), where i and j signify the dimension of the problem and population size(14)xij=li+rand0,1×ui-li.


### 3.2. Employed Bee Phase

The employed bees will search the neighborhoods of the initial positions, in hopes of better position using (15) and modify the initial population with better positions(15)vij=xij+ϕij×xij-xil.


### 3.3. Grey Wolf Phase

Social Hierarchy. In GWO [33] algorithm gray wolves are categorized among α-alpha wolf (fittest solution of the population), β-Beta wolf (second best solution of the population), δ-delta wolf (third best solution of the population), and ω-omega wolves (all other solutions). The hunting for prey (the optimum solution) is guided by alpha, beta, and delta wolves with the help of omega wolves [24, 33].Encircling the Prey. In the next stage, according to gray wolves the prey will be encircled and these specific behaviours can mathematically be replicated by(16)Dij=Cit×Xijpt-Xijt,Xijt+1=Xijt-Ait×Dit,where(17)Ait=2×ait×ri1-ait,Cit=2×ri2.The components of Ait are reduced from 2 to 0 with respect to iterations and ri1, ri2 are vectors of random number within [0,1].Hunting. In optimization perspective, the prey actually refers to the global optima. But in practical scenario it is almost impossible to find global best solution. So in order to imitate the hunting mechanism of grey wolf the best, second, and third best solution are considered as alpha, beta, and delta wolf and used for hunting. The first step is to create distance vectors using α, β, and δ wolves by (18) followed by (19) which gives three different positions for a particular individual of the population. The modified position will be evaluated using (20)(18)Dijα=Ci1×Xijα-Xij;Dijβ=Ci2×Xijβ-Xij;Dijδ=Ci3×Xijδ-Xij;(19)Xij1=Xijα-Ai1×Dijα;Xij2=Xijβ-Ai2×Dijβ;Xij3=Xijδ-Ai3×Dijδ;(20)Xijt+1=Xij1+Xij2+Xij33.


### 3.4. Scout Bee Phase

Consistent individuals that do not improve their position will be converted to scout based on a predefined parameter, limited and forced to search the search space without any guidance using(21)wim=li+rand0,1×ui-li.


### 3.5. Selection

The best solution so far is taken as the global best for the next iteration. This process continues until the iteration number reaches the maximum cycle number (MI). The flow chart forh-ABC/GWO is shown in Figure 4.Figure 4
The detailed flow chart of hybrid-Artificial Bee Colony (ABC)/Grey Wolf Optimization (GWO) (h-ABC/GWO) algorithm.Multiobjective Optimization. As the problem concerned is multiobjective in nature, the algorithm is made suitable for handling multiple objectives simultaneously. A typical multiobjective problem can be stated as(22)Minimize,fnXn=1,2,…,obj;Subject toApX=0p=1,…,XBqX≤q=01,…,YConstraints,where fn(x) is the nth objective function and obj is the total number of objectives. In order to find best compromise solution the first step is to create Pareto optimal front consisting all nondominated individual, which can be found by(23)fiX≤fiY∀i∈1,2,…,obj,fjX≤fjX∃j∈1,2,…,obj.The individualsX and Y satisfying (21) state that X dominates Y. So, in other words solution x is nondominated over solution y. Among all the solutions, the solutions that are not dominated by any other solution are termed as nondominated solution. The Pareto optimal front is created by nondominated solutions only.Best Compromise Solution. Among the individuals of Pareto optimal front there is one individual that has to be termed as best individual. For this purpose, a membership function is used to fit all the objectives within unity bound so that obj, number of objectives, values can be compared accurately. The membership function, gnj forjth individual with respect to nth objective function can be calculated as(24)gnj=1forfnj≤fnminfnmax-fnjfnmax-fnminforfnmin<fnj<fnmax0forfnj≤fnmax,where M is the size of Pareto optimal front and fnmin and fnmax are the lowest and highest function value of nth objective function. In order to find the best compromise solution another normalized function for each individual is found using(25)hj=∑n=1objgnj∑j=1M∑n=1objgnj.The individual with minimumhj is taken as the best compromise solution.


## 3.1. Initialization

The initial population will be created using (14), where i and j signify the dimension of the problem and population size(14)xij=li+rand0,1×ui-li.


## 3.2. Employed Bee Phase

The employed bees will search the neighborhoods of the initial positions, in hopes of better position using (15) and modify the initial population with better positions(15)vij=xij+ϕij×xij-xil.


## 3.3. Grey Wolf Phase

Social Hierarchy. In GWO [33] algorithm gray wolves are categorized among α-alpha wolf (fittest solution of the population), β-Beta wolf (second best solution of the population), δ-delta wolf (third best solution of the population), and ω-omega wolves (all other solutions). The hunting for prey (the optimum solution) is guided by alpha, beta, and delta wolves with the help of omega wolves [24, 33].Encircling the Prey. In the next stage, according to gray wolves the prey will be encircled and these specific behaviours can mathematically be replicated by(16)Dij=Cit×Xijpt-Xijt,Xijt+1=Xijt-Ait×Dit,where(17)Ait=2×ait×ri1-ait,Cit=2×ri2.The components of Ait are reduced from 2 to 0 with respect to iterations and ri1, ri2 are vectors of random number within [0,1].Hunting. In optimization perspective, the prey actually refers to the global optima. But in practical scenario it is almost impossible to find global best solution. So in order to imitate the hunting mechanism of grey wolf the best, second, and third best solution are considered as alpha, beta, and delta wolf and used for hunting. The first step is to create distance vectors using α, β, and δ wolves by (18) followed by (19) which gives three different positions for a particular individual of the population. The modified position will be evaluated using (20)(18)Dijα=Ci1×Xijα-Xij;Dijβ=Ci2×Xijβ-Xij;Dijδ=Ci3×Xijδ-Xij;(19)Xij1=Xijα-Ai1×Dijα;Xij2=Xijβ-Ai2×Dijβ;Xij3=Xijδ-Ai3×Dijδ;(20)Xijt+1=Xij1+Xij2+Xij33.


## 3.4. Scout Bee Phase

Consistent individuals that do not improve their position will be converted to scout based on a predefined parameter, limited and forced to search the search space without any guidance using(21)wim=li+rand0,1×ui-li.


## 3.5. Selection

The best solution so far is taken as the global best for the next iteration. This process continues until the iteration number reaches the maximum cycle number (MI). The flow chart forh-ABC/GWO is shown in Figure 4.Figure 4
The detailed flow chart of hybrid-Artificial Bee Colony (ABC)/Grey Wolf Optimization (GWO) (h-ABC/GWO) algorithm.Multiobjective Optimization. As the problem concerned is multiobjective in nature, the algorithm is made suitable for handling multiple objectives simultaneously. A typical multiobjective problem can be stated as(22)Minimize,fnXn=1,2,…,obj;Subject toApX=0p=1,…,XBqX≤q=01,…,YConstraints,where fn(x) is the nth objective function and obj is the total number of objectives. In order to find best compromise solution the first step is to create Pareto optimal front consisting all nondominated individual, which can be found by(23)fiX≤fiY∀i∈1,2,…,obj,fjX≤fjX∃j∈1,2,…,obj.The individualsX and Y satisfying (21) state that X dominates Y. So, in other words solution x is nondominated over solution y. Among all the solutions, the solutions that are not dominated by any other solution are termed as nondominated solution. The Pareto optimal front is created by nondominated solutions only.Best Compromise Solution. Among the individuals of Pareto optimal front there is one individual that has to be termed as best individual. For this purpose, a membership function is used to fit all the objectives within unity bound so that obj, number of objectives, values can be compared accurately. The membership function, gnj forjth individual with respect to nth objective function can be calculated as(24)gnj=1forfnj≤fnminfnmax-fnjfnmax-fnminforfnmin<fnj<fnmax0forfnj≤fnmax,where M is the size of Pareto optimal front and fnmin and fnmax are the lowest and highest function value of nth objective function. In order to find the best compromise solution another normalized function for each individual is found using(25)hj=∑n=1objgnj∑j=1M∑n=1objgnj.The individual with minimumhj is taken as the best compromise solution.


## 4. Modelling HTUC inh-ABC/GWO Framework

In this section, an algorithm based onh-ABC/GWO method for solving hydrothermal unit commitment coordination problem is described. According to Section 3(i), the initial subpopulations are created where each string of subpopulations (X) in the context of HTUC can be represented by(26)X=QPT,where(27)Q=Qh1,1Qh1,2⋯Qh1,TQh2,1⋯QhNh,T-1QhNh,TT,P=PS1,1PS1,2⋯PS1,TPS2,1⋯PSNS,T-1PSNS,TT.According to the algorithm strategy the hydro discharges (Q) and thermal generations (P) are generated randomly; that is why it is natural for these generated values which may not satisfy most of the limitations, specially equality limitations (end volume limitation and power balance limitation) stated in Section 2. In order to counter the abovementioned problem, constraint handling techniques based on rectification method are used, which are shown in Section 5.Constraint handling technique used to satisfy hydro limitations is as follows:(1)
Considering the initial volume as reference, the end volume ofith unit can be calculated using(28a)Vh0i,T=Vhi,1+∑t=1TIhi,t-∑t=1TQhi,t+∑t=1T∑m=1RuiQhm,t-τm.So, the end volume mismatch is(28b)misi=Vh0i,T-Vhi,T.(2)
A dependent time variable,z, is randomly generated provided z is a natural number and 1≤z≤24.(3)
In order to minimize the mismatch, the discharge atzth hour of the concerned unit is set to a preferable value within its upper and lower bound. If discharge exceeds boundary limit, it is set to its limit. And misi will be calculated for a second time using (28a) and (28b).(4)
Ifmisi is still unacceptable, z will be generated again and processes 2 and 3 continue until misi decreases to an acceptable value.This process is applied on every hydro unit present in the system from least depended to most depended units. Normally a piecewise linear function is used along with Lagrangian multiplier to handle the end volume constraint but this process is time-consuming and erroneous as linearization compromises the accuracy of the system. For this reason in this proposed heuristic approach a rectification based method is used to solve the end volume limitation.Once all the mentioned limitations concerning hydro network are satisfied, the hydro power generation at each hour has to be calculated (see (8)) in order to find the rest of the demanded power that has to be supplied through available thermal units.The constraint handling technique for power balance is described below:(1)
As the minimum up/down time limitations are considered, the units that are available on a particular time interval have to be found based on previous hours and the value (generation) of unavailable units is set to zero. The required power is distributed among the available units.(2)
A dependent variable,zs, is randomly generated, provided zs is a natural number and zS∈NSo, where NSo is the set of available units.(3)
Power balance mismatch attth hour is(29)mist=Preq-∑i=1NSoPSi,t,where(30)Preq=PDi,t-∑i=1NhPhi,t.(4)
The power generation ofzsth unit is set to preferable value within its boundary in order to reduce the mismatch, but due to ramp rate limitation the boundary limits are reduced according to previous hours and upper/lower limit and mist will be calculated.(5)
Ifmist is not within acceptable limit, once again zs will be generated and Steps 2 to 4 continue until mist reduces to acceptable limit.This process will be implemented for every hour to balance generated power at each hour.In any metaheuristic algorithm the constraint handling techniques are very crucial factor. But the problem adopted in this research is not only nonlinear but also highly nonconvex in nature as shown in Figures1, 2, and 3. The basic techniques like Lagrangian multiplier can be used to handle such complex constraints but with the cost of accuracy as this technique requires linearization which additionally increases the time consumption. On the other hand penalty factor method can also be used which forcibly brings the infeasible solutions into feasible region. The selection of penalty factor is very complex considering the complexity and vastness of the problem. Also due to highly restricted feasible region of the search space which is actually scattered around the search space, the penalty factors are not that effective as the proposed heuristic technique shown above. The proposed rectification based technique efficiently brings the infeasible solutions into the feasible region compared to that by penalty factor method.Algorithm Implementation (a)
After creating the subpopulations according to Section3 (i), the positions are modified according to constraint handling techniques given in previous two sections and evaluated for fitness value, using(31)fitXj=11+fXj, whereXj is thejth individual of the population. f(Xj) is the function value evaluated by (1a) and (1b) of the concerned position. In case of multiple objectives, fitness value will be found using (31) for both objectives followed by membership function evaluation using (24) and best compromised solution will be found by (25). In every step this procedure needs to be followed to accurately find the best compromised solution. (b)
For each individual according to the grey wolf optimization phase the better individuals will be selected modifying the population. (c)
The new positions will be evaluated using (31) after constraint handling techniques and based on better fitness value another subpopulation (modified subpopulation) will be created. (d)
Each individual of population will create another position according to (21). The new positions are evaluated for fitness value. (e)
Based on the fitness values of new found positions and their concerned individual of population found in Step(c), a final population of same size will be created. (f)
The best solution of final population will be considered as global best. (g)
Finally, an iteration counter increases.Steps (b)–(g) will continue until iteration counter reachesMI.For multiple objectives, instead of best solution best compromised solution needed to be found. In that case the global best compromised solution can be calculated using pareto optimal front given in (23). After that based on the membership function global best solution will be evaluated.


## 5. Numerical Result

In this section, two test systems are considered for experimentation of the performance of proposed algorithm on HTUC problem. The scheduled time period for both systems is taken as 24 hours. The objective of these studies is to showcase the effectiveness ofh-ABC/GWO to handle HTUC under various scenarios.


### 5.1. Parameter Selection

The parameters of theh-ABC/GWO are tuned manually and reported in Table 1 along with ABC and GWO algorithm. The developed algorithms are run independently with different maximum cycle number (MI) keeping other parameters constant. It is observed from the convergence characteristics that the algorithms converge before 3500 iterations for Test System I and the same is for Test System II as well. So, in order to compare the convergences of different algorithms MI is set to the value as in Table 1. The rest of the other parameters of individual algorithms are found by trial and error method.Table 1
Parameters used for test systems in algorithms.
Test System I Test System II Maximum cycle number (MI) = 3500 Maximum cycle number (MI) = 3500 Limit = 25 Limit = 35 Population size = 40 Population size = 80


### 5.2. Results

The mentioned algorithms are coded in MATLAB (Version: 8.1.0.604 (R2013a)) environment and simulated on Dell XPS15 (2760QM®, Intel configured), 3rd Generation Quad-Core i7 Processor ~3.2 GHz processor speed and 12 GB RAM.


#### 5.2.1. Test System I

Test System I [27] consists of four hydro and ten thermal units. The hydro units are connected in a multichain cascaded formation. The detailed data for hydro units are taken from literature [6] and essential data regarding thermal units and required load are taken from [21]. Valve point loading effect for thermal units is taken into consideration in the fuel cost function.This system is solved in two different scenarios. In first scenario (Case I) all the hydro limitations are considered along with thermal generation limitation (see (9)–(13)) while in second case (Case II) all the limitations of Section 2 are considered.The number of dimensions in the multidimensional search space or decision variables in Case I is24∗(4+10)=336, with 4+24=28 equality and 96 inequality constraints while in second case due to accumulation of binary variables, the total decision variable increases to 576 with added complexity due to ramping restriction and intertemporal limitations.Single Objective EnvironmentCase I. In this case the thermal power is distributed by conventional load dispatch [34] technique without any influence of intertemporal constraints. The mentioned algorithms are simulated 10 times independently and successively and the best solution obtained is$165921.53 byh-ABC/GWO. The solutions of finest individuals found by mentioned algorithms are given in Table 2 Case I column with respect to best, worst, and average solution (of 35 simulations) along with the result of mixed-integer nonlinear programming (MINLP), small population based PSO, and improved DE published in literature [27, 35] and the convergence characteristics forh-ABC/GWO, ABC, GWO, SPPSO, IDE, and MINLP are shown in Figure 5(a). The MINLP modelling of this test case is performed on GAMS software using BONMIN solver.Table 2
Performance comparison for Test System I (Case I and Case II).
Algorithm Case I Case II Best Worst Time Best Worst Time h-ABC/GWO 165921.53 170201.83 202.59 Sec 271736.07 475124.26 394.61 Sec ABC [34] 168386.23 183044.84 213.89 Sec 281645.19 489542.39 1023.64 Sec GWO [32] 180356.68 195635.44 242.65 Sec 325689.26 524586.69 536.98 Sec SPPSO [27] 167710.56 — — — — — IDE [35] 170576.50 — — — — — MINLP [37] 208706.20 — 105.23 Sec — — —Figure 5
(a) Convergence characteristic of different algorithms for Test System I, Case I. (b) Convergence characteristic of different algorithms for Test System I, Case II.
(a) (b)From Table2 and Figure 5, it is observed that there is a significant improvement in the performance of the proposed algorithm as compared to the normal GWO due to added diversity through hybridization in terms of convergence and quality of solutions. It is also shown that the convergence rate for MINLP is far better thanh-ABC/GWO. But the optimal solution of MINLP is very poor with respect to other solutions, which prove higher nonconvexity in the search space. The position (X) in optimal condition, comprising hydro discharges and thermal generations, is shown in Table 3 in hour basis. The water content of every reservoir at each hour can be found using (3). It is evident thath-ABC/GWO algorithm has much faster convergence rate than GWO and ABC. It also avoids premature convergence.Table 3
Performance comparison for Test System II for both Cases I and II in single objective environment.
Algorithm Case I Case II Best Worst Time Best Worst Time h-ABC/GWO 3411361.47 3814805.07 424.36 Sec 4604300.83 6653645.64 394.61 Sec ABC [34] 3923430.55 4811261.91 443.59 Sec 4828932.22 6453793.64 1506.64 Sec GWO [32] 4259836.27 5037176.67 797.53 Sec 5315693.34 7037176.67 624.64 Sec CSA [31] 3503527.75 — — 4823256.68 — — MINLP [37] 4023565.36 — 212.56 Sec — — —Case II. In this scenario, after hydro power production, the rest of required power is distributed among available thermal units by unit commitment coordination [16], where the thermal power distribution is performed using unit commitment, which calls for binary variables in order to define commitment status of thermal units. Due to the mentioned intricacy and additional intertemporal limitations, the search space becomes more nonlinear and vast. The algorithms mentioned in Section 5.1 are implemented on this scenario too and the best solution comes out to be$271736.07 found byh-ABC/GWO. GWO and ABC algorithm are also used to solve the same and the comparison with respect to best, worst, and average solution for 35 consecutive simulations is given in Table 2 in Case II column. Figure 5(b) shows the convergence characteristics of different algorithms proving better solution quality forh-ABC/GWO.


#### 5.2.2. Test System II

In Test System II the number of dimensions is even higher than Test System I and contains almost all the features of a practical hydrothermal system. It incorporates 44 hydro units and 54 thermal units. The hydro system [31] comprises a large multichain hydro network connected together in a cascaded formation. The divisions of water from 20th, 24th, 28th, and 30th hydro units to their downstream units are taken as decision variable.The thermal system [38] comprises two subsystems of 40 and 13 units. The feature of 54th unit is the same as 53rd unit.This system is solved in two different scenarios like Test System I. In Case I the number of variables is24∗44+24∗54+4 (for the diversion of water) = 2356 as thermal power distribution is performed by economic load dispatch while in Case II the number of decision variables is 24∗44+24∗54+24∗54+4=3652 due to inclusion of unit commitment as subproblem. The complexity of this problem is higher than Test System I not only due to additional limitations but also for more complicated structure of the hydro network.Case I. In this case all hydro limitations are considered with thermal power generation limitation for thermal subproblem. Table 1 shows different parameters for the algorithms (h-ABC/GWO, ABC, and MINLP) used. The optimized cost for this scenario is found out to be$3,411,361.47.In order to verify, the results are compared with [31] in Table 3 of Case I column. The convergence characteristics of mentioned algorithms are shown in Figure 6(a). It is observed that there is a significant difference between the performances ofh-ABC/GWO and GWO algorithm.h-ABC/GWO has outperformed GWO and ABC in terms of convergence rate and quality of solutions obtained. It also gives better result over the newly proposed Clonal Selection Algorithm (CSA) reported in [31].Figure 6
(a) Convergence characteristic of different algorithms for Test System II, Case I. (b) Convergence characteristic of different algorithms for Test System II, Case II.
(a) (b)The convergence characteristics show that theh-ABC/GWO is far better than the other mentioned algorithms. The reason behind this vast difference in terms of solution quality for these two algorithm is firstly the lack of diversity in the multidimensional search space and secondly the attraction towards local optimum resulting premature convergence, which is basically a reverse effect of high intensity depth search. Due to lower enormity and inferior nonlinearity in the search space of first problem (with comparison to the second), GWO could reach desired diversity to provide acceptable solution where the effect of premature convergence can be countered with higher population number. But to find an acceptable solution for second problem where the search space is massive with higher nonlinearity requires an algorithm having the capability to provide massive diversity in the multidimensional search space along with controlled depth search.Case II. In this scenario Test System II is solved considering all the limitations of Section 2. The required demand in this case is taken 1 GW less than in Case I for providing sufficient feasible region in nonlinear, nonconvex search space. The problem is also solved byh-ABC/GWO, GWO, and ABC algorithm and optimal solutions are shown in Table 3, Case II column. However data used to solve the problem is not provided in the paper for limited space but will be provided if required. The optimum solution in this case is found to be$4,604,300.83. From the convergence characteristics of Figure 6(b) the superiority ofh-ABC/GWO over GWO and ABC in terms of solution quality is proved.Multiobjective Environment. In multiobjective environment the production cost and emission are considered as the two objectives as mentioned earlier. Both the objectives are needed to be reduced as much as possible.


#### 5.2.3. Test System I

This problem is solved using the same parameter shown in Section5.2.1 for both the cases. The solutions with respect to best and worst solutions are given in Table 4 including time consumption with various other established algorithms like, ABC, GWO, IDE, and recent state of the art MINLP algorithm. And the convergence characteristics for both cases are shown in Figure 7.Table 4
Performance comparison for Test System I for both Case I and Case II in multiobjective environment.
Algorithm Best Worst Cost Emission Time Cost Emission Time h-ABC/GWO 113295.63 19536.65 452.69 123569.65 23326.65 463.59 ABC [34] 167386.23 21831.56 745.36 135643.96 23549.95 745.36 GWO [32] 178546.35 22469.49 721.49 154937.19 24987.48 456.19 IDE [35] 134649.96 234697.49 748.59 194853.15 31497.19 489.46 MINLP [37] 123494.46 200349.45 451.74 — — —Figure 7
Convergence characteristic of (a) production cost and (b) emission for Test System I.
(a) (b)From the convergence characteristics, it is clear that theh-ABC/GWO algorithm is better than other mentioned algorithms in terms of convergence, solution quality, and time consumption.


#### 5.2.4. Test System II

Like the previous section in this case also Test System II is solved under the influence of the parameter used in Section5.2.2. The comparative performance of proposed algorithm and convergence characteristics of the same with respect to other established algorithms is shown in Table 5 and Figure 8 along with time consumption. From these results also the superiority ofh-ABC.GWO can further be consolidated.Table 5
Performance comparison for Test System II for both Case I and Case II in multiobjective environment.
Algorithm Best Worst Cost Emission Time Cost Emission Time h-ABC/GWO 2259863.56 438592.35 1523.32 2456965.58 440015.20 1589.26 ABC [34] 2667386.23 438312.56 1689.36 2785952.63 445963.56 1702.25 GWO [32] 3498830.99 671839.44 1559.45 3925696.54 675895.68 1517.26 MINLP [37] 5939758.65 542569.65 936.36 — — —Figure 8
Convergence characteristic of (a) production cost and (b) emission for Test System II.
(a) (b)In addition to this it is seen that the convergence and time consumption forh-ABC/GWO algorithm are considerably large in case of Test System II compared to those in Test System I for both single and multiobjective environment, which suggest that Test System II is much more complicated than Test System I. But even so the proposed algorithm efficiently solves the problem within acceptable time duration.The reason for such efficiency of the proposed algorithm is due to using the more sophisticated and efficient exploitation technique of GWO algorithm instead of onlooker bee phase of ABC algorithm. These phenomena give the algorithm necessary attribute to efficiently search the globally best solution from such complicated and vast search space. Also the effect of proposed constraint handling technique mentioned in Section4 is considerable. Due to the use of this technique more than 90% of the solutions of the population are feasible for Test System I and more than 80% of solution gets feasible for Test System II. This observation points to the fact that Test System II is much more complicated than first test system.


## 5.1. Parameter Selection

The parameters of theh-ABC/GWO are tuned manually and reported in Table 1 along with ABC and GWO algorithm. The developed algorithms are run independently with different maximum cycle number (MI) keeping other parameters constant. It is observed from the convergence characteristics that the algorithms converge before 3500 iterations for Test System I and the same is for Test System II as well. So, in order to compare the convergences of different algorithms MI is set to the value as in Table 1. The rest of the other parameters of individual algorithms are found by trial and error method.Table 1
Parameters used for test systems in algorithms.
Test System I Test System II Maximum cycle number (MI) = 3500 Maximum cycle number (MI) = 3500 Limit = 25 Limit = 35 Population size = 40 Population size = 80


## 5.2. Results

The mentioned algorithms are coded in MATLAB (Version: 8.1.0.604 (R2013a)) environment and simulated on Dell XPS15 (2760QM®, Intel configured), 3rd Generation Quad-Core i7 Processor ~3.2 GHz processor speed and 12 GB RAM.


### 5.2.1. Test System I

Test System I [27] consists of four hydro and ten thermal units. The hydro units are connected in a multichain cascaded formation. The detailed data for hydro units are taken from literature [6] and essential data regarding thermal units and required load are taken from [21]. Valve point loading effect for thermal units is taken into consideration in the fuel cost function.This system is solved in two different scenarios. In first scenario (Case I) all the hydro limitations are considered along with thermal generation limitation (see (9)–(13)) while in second case (Case II) all the limitations of Section 2 are considered.The number of dimensions in the multidimensional search space or decision variables in Case I is24∗(4+10)=336, with 4+24=28 equality and 96 inequality constraints while in second case due to accumulation of binary variables, the total decision variable increases to 576 with added complexity due to ramping restriction and intertemporal limitations.Single Objective EnvironmentCase I. In this case the thermal power is distributed by conventional load dispatch [34] technique without any influence of intertemporal constraints. The mentioned algorithms are simulated 10 times independently and successively and the best solution obtained is$165921.53 byh-ABC/GWO. The solutions of finest individuals found by mentioned algorithms are given in Table 2 Case I column with respect to best, worst, and average solution (of 35 simulations) along with the result of mixed-integer nonlinear programming (MINLP), small population based PSO, and improved DE published in literature [27, 35] and the convergence characteristics forh-ABC/GWO, ABC, GWO, SPPSO, IDE, and MINLP are shown in Figure 5(a). The MINLP modelling of this test case is performed on GAMS software using BONMIN solver.Table 2
Performance comparison for Test System I (Case I and Case II).
Algorithm Case I Case II Best Worst Time Best Worst Time h-ABC/GWO 165921.53 170201.83 202.59 Sec 271736.07 475124.26 394.61 Sec ABC [34] 168386.23 183044.84 213.89 Sec 281645.19 489542.39 1023.64 Sec GWO [32] 180356.68 195635.44 242.65 Sec 325689.26 524586.69 536.98 Sec SPPSO [27] 167710.56 — — — — — IDE [35] 170576.50 — — — — — MINLP [37] 208706.20 — 105.23 Sec — — —Figure 5
(a) Convergence characteristic of different algorithms for Test System I, Case I. (b) Convergence characteristic of different algorithms for Test System I, Case II.
(a) (b)From Table2 and Figure 5, it is observed that there is a significant improvement in the performance of the proposed algorithm as compared to the normal GWO due to added diversity through hybridization in terms of convergence and quality of solutions. It is also shown that the convergence rate for MINLP is far better thanh-ABC/GWO. But the optimal solution of MINLP is very poor with respect to other solutions, which prove higher nonconvexity in the search space. The position (X) in optimal condition, comprising hydro discharges and thermal generations, is shown in Table 3 in hour basis. The water content of every reservoir at each hour can be found using (3). It is evident thath-ABC/GWO algorithm has much faster convergence rate than GWO and ABC. It also avoids premature convergence.Table 3
Performance comparison for Test System II for both Cases I and II in single objective environment.
Algorithm Case I Case II Best Worst Time Best Worst Time h-ABC/GWO 3411361.47 3814805.07 424.36 Sec 4604300.83 6653645.64 394.61 Sec ABC [34] 3923430.55 4811261.91 443.59 Sec 4828932.22 6453793.64 1506.64 Sec GWO [32] 4259836.27 5037176.67 797.53 Sec 5315693.34 7037176.67 624.64 Sec CSA [31] 3503527.75 — — 4823256.68 — — MINLP [37] 4023565.36 — 212.56 Sec — — —Case II. In this scenario, after hydro power production, the rest of required power is distributed among available thermal units by unit commitment coordination [16], where the thermal power distribution is performed using unit commitment, which calls for binary variables in order to define commitment status of thermal units. Due to the mentioned intricacy and additional intertemporal limitations, the search space becomes more nonlinear and vast. The algorithms mentioned in Section 5.1 are implemented on this scenario too and the best solution comes out to be$271736.07 found byh-ABC/GWO. GWO and ABC algorithm are also used to solve the same and the comparison with respect to best, worst, and average solution for 35 consecutive simulations is given in Table 2 in Case II column. Figure 5(b) shows the convergence characteristics of different algorithms proving better solution quality forh-ABC/GWO.


### 5.2.2. Test System II

In Test System II the number of dimensions is even higher than Test System I and contains almost all the features of a practical hydrothermal system. It incorporates 44 hydro units and 54 thermal units. The hydro system [31] comprises a large multichain hydro network connected together in a cascaded formation. The divisions of water from 20th, 24th, 28th, and 30th hydro units to their downstream units are taken as decision variable.The thermal system [38] comprises two subsystems of 40 and 13 units. The feature of 54th unit is the same as 53rd unit.This system is solved in two different scenarios like Test System I. In Case I the number of variables is24∗44+24∗54+4 (for the diversion of water) = 2356 as thermal power distribution is performed by economic load dispatch while in Case II the number of decision variables is 24∗44+24∗54+24∗54+4=3652 due to inclusion of unit commitment as subproblem. The complexity of this problem is higher than Test System I not only due to additional limitations but also for more complicated structure of the hydro network.Case I. In this case all hydro limitations are considered with thermal power generation limitation for thermal subproblem. Table 1 shows different parameters for the algorithms (h-ABC/GWO, ABC, and MINLP) used. The optimized cost for this scenario is found out to be$3,411,361.47.In order to verify, the results are compared with [31] in Table 3 of Case I column. The convergence characteristics of mentioned algorithms are shown in Figure 6(a). It is observed that there is a significant difference between the performances ofh-ABC/GWO and GWO algorithm.h-ABC/GWO has outperformed GWO and ABC in terms of convergence rate and quality of solutions obtained. It also gives better result over the newly proposed Clonal Selection Algorithm (CSA) reported in [31].Figure 6
(a) Convergence characteristic of different algorithms for Test System II, Case I. (b) Convergence characteristic of different algorithms for Test System II, Case II.
(a) (b)The convergence characteristics show that theh-ABC/GWO is far better than the other mentioned algorithms. The reason behind this vast difference in terms of solution quality for these two algorithm is firstly the lack of diversity in the multidimensional search space and secondly the attraction towards local optimum resulting premature convergence, which is basically a reverse effect of high intensity depth search. Due to lower enormity and inferior nonlinearity in the search space of first problem (with comparison to the second), GWO could reach desired diversity to provide acceptable solution where the effect of premature convergence can be countered with higher population number. But to find an acceptable solution for second problem where the search space is massive with higher nonlinearity requires an algorithm having the capability to provide massive diversity in the multidimensional search space along with controlled depth search.Case II. In this scenario Test System II is solved considering all the limitations of Section 2. The required demand in this case is taken 1 GW less than in Case I for providing sufficient feasible region in nonlinear, nonconvex search space. The problem is also solved byh-ABC/GWO, GWO, and ABC algorithm and optimal solutions are shown in Table 3, Case II column. However data used to solve the problem is not provided in the paper for limited space but will be provided if required. The optimum solution in this case is found to be$4,604,300.83. From the convergence characteristics of Figure 6(b) the superiority ofh-ABC/GWO over GWO and ABC in terms of solution quality is proved.Multiobjective Environment. In multiobjective environment the production cost and emission are considered as the two objectives as mentioned earlier. Both the objectives are needed to be reduced as much as possible.


### 5.2.3. Test System I

This problem is solved using the same parameter shown in Section5.2.1 for both the cases. The solutions with respect to best and worst solutions are given in Table 4 including time consumption with various other established algorithms like, ABC, GWO, IDE, and recent state of the art MINLP algorithm. And the convergence characteristics for both cases are shown in Figure 7.Table 4
Performance comparison for Test System I for both Case I and Case II in multiobjective environment.
Algorithm Best Worst Cost Emission Time Cost Emission Time h-ABC/GWO 113295.63 19536.65 452.69 123569.65 23326.65 463.59 ABC [34] 167386.23 21831.56 745.36 135643.96 23549.95 745.36 GWO [32] 178546.35 22469.49 721.49 154937.19 24987.48 456.19 IDE [35] 134649.96 234697.49 748.59 194853.15 31497.19 489.46 MINLP [37] 123494.46 200349.45 451.74 — — —Figure 7
Convergence characteristic of (a) production cost and (b) emission for Test System I.
(a) (b)From the convergence characteristics, it is clear that theh-ABC/GWO algorithm is better than other mentioned algorithms in terms of convergence, solution quality, and time consumption.


### 5.2.4. Test System II

Like the previous section in this case also Test System II is solved under the influence of the parameter used in Section5.2.2. The comparative performance of proposed algorithm and convergence characteristics of the same with respect to other established algorithms is shown in Table 5 and Figure 8 along with time consumption. From these results also the superiority ofh-ABC.GWO can further be consolidated.Table 5
Performance comparison for Test System II for both Case I and Case II in multiobjective environment.
Algorithm Best Worst Cost Emission Time Cost Emission Time h-ABC/GWO 2259863.56 438592.35 1523.32 2456965.58 440015.20 1589.26 ABC [34] 2667386.23 438312.56 1689.36 2785952.63 445963.56 1702.25 GWO [32] 3498830.99 671839.44 1559.45 3925696.54 675895.68 1517.26 MINLP [37] 5939758.65 542569.65 936.36 — — —Figure 8
Convergence characteristic of (a) production cost and (b) emission for Test System II.
(a) (b)In addition to this it is seen that the convergence and time consumption forh-ABC/GWO algorithm are considerably large in case of Test System II compared to those in Test System I for both single and multiobjective environment, which suggest that Test System II is much more complicated than Test System I. But even so the proposed algorithm efficiently solves the problem within acceptable time duration.The reason for such efficiency of the proposed algorithm is due to using the more sophisticated and efficient exploitation technique of GWO algorithm instead of onlooker bee phase of ABC algorithm. These phenomena give the algorithm necessary attribute to efficiently search the globally best solution from such complicated and vast search space. Also the effect of proposed constraint handling technique mentioned in Section4 is considerable. Due to the use of this technique more than 90% of the solutions of the population are feasible for Test System I and more than 80% of solution gets feasible for Test System II. This observation points to the fact that Test System II is much more complicated than first test system.


## 5.2.1. Test System I

Test System I [27] consists of four hydro and ten thermal units. The hydro units are connected in a multichain cascaded formation. The detailed data for hydro units are taken from literature [6] and essential data regarding thermal units and required load are taken from [21]. Valve point loading effect for thermal units is taken into consideration in the fuel cost function.This system is solved in two different scenarios. In first scenario (Case I) all the hydro limitations are considered along with thermal generation limitation (see (9)–(13)) while in second case (Case II) all the limitations of Section 2 are considered.The number of dimensions in the multidimensional search space or decision variables in Case I is24∗(4+10)=336, with 4+24=28 equality and 96 inequality constraints while in second case due to accumulation of binary variables, the total decision variable increases to 576 with added complexity due to ramping restriction and intertemporal limitations.Single Objective EnvironmentCase I. In this case the thermal power is distributed by conventional load dispatch [34] technique without any influence of intertemporal constraints. The mentioned algorithms are simulated 10 times independently and successively and the best solution obtained is$165921.53 byh-ABC/GWO. The solutions of finest individuals found by mentioned algorithms are given in Table 2 Case I column with respect to best, worst, and average solution (of 35 simulations) along with the result of mixed-integer nonlinear programming (MINLP), small population based PSO, and improved DE published in literature [27, 35] and the convergence characteristics forh-ABC/GWO, ABC, GWO, SPPSO, IDE, and MINLP are shown in Figure 5(a). The MINLP modelling of this test case is performed on GAMS software using BONMIN solver.Table 2
Performance comparison for Test System I (Case I and Case II).
Algorithm Case I Case II Best Worst Time Best Worst Time h-ABC/GWO 165921.53 170201.83 202.59 Sec 271736.07 475124.26 394.61 Sec ABC [34] 168386.23 183044.84 213.89 Sec 281645.19 489542.39 1023.64 Sec GWO [32] 180356.68 195635.44 242.65 Sec 325689.26 524586.69 536.98 Sec SPPSO [27] 167710.56 — — — — — IDE [35] 170576.50 — — — — — MINLP [37] 208706.20 — 105.23 Sec — — —Figure 5
(a) Convergence characteristic of different algorithms for Test System I, Case I. (b) Convergence characteristic of different algorithms for Test System I, Case II.
(a) (b)From Table2 and Figure 5, it is observed that there is a significant improvement in the performance of the proposed algorithm as compared to the normal GWO due to added diversity through hybridization in terms of convergence and quality of solutions. It is also shown that the convergence rate for MINLP is far better thanh-ABC/GWO. But the optimal solution of MINLP is very poor with respect to other solutions, which prove higher nonconvexity in the search space. The position (X) in optimal condition, comprising hydro discharges and thermal generations, is shown in Table 3 in hour basis. The water content of every reservoir at each hour can be found using (3). It is evident thath-ABC/GWO algorithm has much faster convergence rate than GWO and ABC. It also avoids premature convergence.Table 3
Performance comparison for Test System II for both Cases I and II in single objective environment.
Algorithm Case I Case II Best Worst Time Best Worst Time h-ABC/GWO 3411361.47 3814805.07 424.36 Sec 4604300.83 6653645.64 394.61 Sec ABC [34] 3923430.55 4811261.91 443.59 Sec 4828932.22 6453793.64 1506.64 Sec GWO [32] 4259836.27 5037176.67 797.53 Sec 5315693.34 7037176.67 624.64 Sec CSA [31] 3503527.75 — — 4823256.68 — — MINLP [37] 4023565.36 — 212.56 Sec — — —Case II. In this scenario, after hydro power production, the rest of required power is distributed among available thermal units by unit commitment coordination [16], where the thermal power distribution is performed using unit commitment, which calls for binary variables in order to define commitment status of thermal units. Due to the mentioned intricacy and additional intertemporal limitations, the search space becomes more nonlinear and vast. The algorithms mentioned in Section 5.1 are implemented on this scenario too and the best solution comes out to be$271736.07 found byh-ABC/GWO. GWO and ABC algorithm are also used to solve the same and the comparison with respect to best, worst, and average solution for 35 consecutive simulations is given in Table 2 in Case II column. Figure 5(b) shows the convergence characteristics of different algorithms proving better solution quality forh-ABC/GWO.


## 5.2.2. Test System II

In Test System II the number of dimensions is even higher than Test System I and contains almost all the features of a practical hydrothermal system. It incorporates 44 hydro units and 54 thermal units. The hydro system [31] comprises a large multichain hydro network connected together in a cascaded formation. The divisions of water from 20th, 24th, 28th, and 30th hydro units to their downstream units are taken as decision variable.The thermal system [38] comprises two subsystems of 40 and 13 units. The feature of 54th unit is the same as 53rd unit.This system is solved in two different scenarios like Test System I. In Case I the number of variables is24∗44+24∗54+4 (for the diversion of water) = 2356 as thermal power distribution is performed by economic load dispatch while in Case II the number of decision variables is 24∗44+24∗54+24∗54+4=3652 due to inclusion of unit commitment as subproblem. The complexity of this problem is higher than Test System I not only due to additional limitations but also for more complicated structure of the hydro network.Case I. In this case all hydro limitations are considered with thermal power generation limitation for thermal subproblem. Table 1 shows different parameters for the algorithms (h-ABC/GWO, ABC, and MINLP) used. The optimized cost for this scenario is found out to be$3,411,361.47.In order to verify, the results are compared with [31] in Table 3 of Case I column. The convergence characteristics of mentioned algorithms are shown in Figure 6(a). It is observed that there is a significant difference between the performances ofh-ABC/GWO and GWO algorithm.h-ABC/GWO has outperformed GWO and ABC in terms of convergence rate and quality of solutions obtained. It also gives better result over the newly proposed Clonal Selection Algorithm (CSA) reported in [31].Figure 6
(a) Convergence characteristic of different algorithms for Test System II, Case I. (b) Convergence characteristic of different algorithms for Test System II, Case II.
(a) (b)The convergence characteristics show that theh-ABC/GWO is far better than the other mentioned algorithms. The reason behind this vast difference in terms of solution quality for these two algorithm is firstly the lack of diversity in the multidimensional search space and secondly the attraction towards local optimum resulting premature convergence, which is basically a reverse effect of high intensity depth search. Due to lower enormity and inferior nonlinearity in the search space of first problem (with comparison to the second), GWO could reach desired diversity to provide acceptable solution where the effect of premature convergence can be countered with higher population number. But to find an acceptable solution for second problem where the search space is massive with higher nonlinearity requires an algorithm having the capability to provide massive diversity in the multidimensional search space along with controlled depth search.Case II. In this scenario Test System II is solved considering all the limitations of Section 2. The required demand in this case is taken 1 GW less than in Case I for providing sufficient feasible region in nonlinear, nonconvex search space. The problem is also solved byh-ABC/GWO, GWO, and ABC algorithm and optimal solutions are shown in Table 3, Case II column. However data used to solve the problem is not provided in the paper for limited space but will be provided if required. The optimum solution in this case is found to be$4,604,300.83. From the convergence characteristics of Figure 6(b) the superiority ofh-ABC/GWO over GWO and ABC in terms of solution quality is proved.Multiobjective Environment. In multiobjective environment the production cost and emission are considered as the two objectives as mentioned earlier. Both the objectives are needed to be reduced as much as possible.


## 5.2.3. Test System I

This problem is solved using the same parameter shown in Section5.2.1 for both the cases. The solutions with respect to best and worst solutions are given in Table 4 including time consumption with various other established algorithms like, ABC, GWO, IDE, and recent state of the art MINLP algorithm. And the convergence characteristics for both cases are shown in Figure 7.Table 4
Performance comparison for Test System I for both Case I and Case II in multiobjective environment.
Algorithm Best Worst Cost Emission Time Cost Emission Time h-ABC/GWO 113295.63 19536.65 452.69 123569.65 23326.65 463.59 ABC [34] 167386.23 21831.56 745.36 135643.96 23549.95 745.36 GWO [32] 178546.35 22469.49 721.49 154937.19 24987.48 456.19 IDE [35] 134649.96 234697.49 748.59 194853.15 31497.19 489.46 MINLP [37] 123494.46 200349.45 451.74 — — —Figure 7
Convergence characteristic of (a) production cost and (b) emission for Test System I.
(a) (b)From the convergence characteristics, it is clear that theh-ABC/GWO algorithm is better than other mentioned algorithms in terms of convergence, solution quality, and time consumption.


## 5.2.4. Test System II

Like the previous section in this case also Test System II is solved under the influence of the parameter used in Section5.2.2. The comparative performance of proposed algorithm and convergence characteristics of the same with respect to other established algorithms is shown in Table 5 and Figure 8 along with time consumption. From these results also the superiority ofh-ABC.GWO can further be consolidated.Table 5
Performance comparison for Test System II for both Case I and Case II in multiobjective environment.
Algorithm Best Worst Cost Emission Time Cost Emission Time h-ABC/GWO 2259863.56 438592.35 1523.32 2456965.58 440015.20 1589.26 ABC [34] 2667386.23 438312.56 1689.36 2785952.63 445963.56 1702.25 GWO [32] 3498830.99 671839.44 1559.45 3925696.54 675895.68 1517.26 MINLP [37] 5939758.65 542569.65 936.36 — — —Figure 8
Convergence characteristic of (a) production cost and (b) emission for Test System II.
(a) (b)In addition to this it is seen that the convergence and time consumption forh-ABC/GWO algorithm are considerably large in case of Test System II compared to those in Test System I for both single and multiobjective environment, which suggest that Test System II is much more complicated than Test System I. But even so the proposed algorithm efficiently solves the problem within acceptable time duration.The reason for such efficiency of the proposed algorithm is due to using the more sophisticated and efficient exploitation technique of GWO algorithm instead of onlooker bee phase of ABC algorithm. These phenomena give the algorithm necessary attribute to efficiently search the globally best solution from such complicated and vast search space. Also the effect of proposed constraint handling technique mentioned in Section4 is considerable. Due to the use of this technique more than 90% of the solutions of the population are feasible for Test System I and more than 80% of solution gets feasible for Test System II. This observation points to the fact that Test System II is much more complicated than first test system.


## 6. Conclusion

The hydrothermal unit commitment coordination (HTUC) is modelled having included various limitations of hydro and thermal subsystems. The analysis is performed based on an efficient optimization algorithm for better performance. A hybridized approach of ABC and GWO is presented in this work in both single and multiobjective environment. This proposed method is tested with two different test systems in different scenarios. Though the proposed method is performing better than other mentioned algorithms through thorough analysis, it is found that, in case of Test System I, where the number of variables is comparatively low, theh-ABC/GWO helps to improve robustness and convergence. But in Test System II (where the number of variables is considerably high) the proposed approach not only improves the convergence but also successfully avoids the local optimum along with better robustness. The computation time and cost function are improved by the proposed method and are compared with other methods. The computation performance is found to be significantly superior to the proposed method using hybridization. Superiority of the proposedh-ABC/GWO algorithm is further consolidated through experimentation with a larger system with faster convergence and better quality of solutions.



---

*Source: 4529836-2016-12-28.xml*