Multi-objective day-ahead scheduling of microgrids using modified grey wolf optimizer algorithm

Abstract

Investigation of the environmental/economic optimal operation management of a microgrid (MG) as a case study for applying a novel modified multi-objective grey wolf optimizer (MMOGWO) algorithm is presented in this paper. MGs can be considered as a fundamental solution in order for distributed generators’ (DGs) management in future smart grids. In the multi-objective problems, since the objective functions are conflict, the best compromised solution should be extracted through an efficient approach. Accordingly, a proper method is applied for exploring the best compromised solution. Additionally, a novel distance-based method is proposed to control the size of the repository within an aimed limit which leads to a fast and precise convergence along with a well-distributed Pareto optimal front. The proposed method is implemented in a typical grid-connected MG with non-dispatchable units including renewable energy sources (RESs), along with a hybrid power source (micro-turbine, fuel-cell and battery) as dispatchable units, to accumulate excess energy or to equalize power mismatch, by optimal scheduling of DGs and the power exchange between the utility grid and storage system. The efficiency of the suggested algorithm in satisfying the load and optimizing the objective functions is validated through comparison with different methods, including PSO and the original GWO.

Keywords

Multi objective optimal operation management pareto optimal solution modified grey wolf optimizer micro-grid renewable energy sources

1 Introduction

Recent increases in energy prices and environmental concerns have led to the penetration of renewable distribution generators (DGs) in distribution systems [1]. Microgrid (MG) is a concept resulted from the need of reliable power systems with clean energy sources which can make it easy to have a satisfactory communication and optimal energy management of the power system [2]. A review of modelling, planning and energy management of an MG is presented in [3]. Since control and operation management of MGs are midway through improvement, diverse techniques are proposed in order for optimization of these networks. Consequently, more precise energy source scheduling in MGs considering different objectives seems to be required. Various researches have been conducted dealing with optimal operation scheduling considering different constraints and objectives [4 –12].

The choice of optimization technique depends on several different factors; therefore, different methods such as mathematical programming based optimization approaches and meta-heuristic algorithms can be proposed to solve problems [13, 14]. The lambda iterative method, gradient projection method, linear/non-linear programming interior point methods, dynamic programming, etc. can be named as examples of mathematical programming based methods. On the other hand, meta-heuristic algorithms are proper alternatives to mathematical programming. In the case of solving multi-objective optimization problems, the major motivation of employing meta-heuristic algorithms instead of mathematical programming approaches is the ability of meta-heuristics in finding different solutions in the Pareto optimal front in just one execution of the algorithm, while the mathematical programming methods apply a sequence of independent executions. Furthermore, meta-heuristic algorithms are not sensitive to the continuity and formation of the Pareto front [13].

In the MGs’ environmental∖economic management a certain number of DG units are supposed to supply the load while minimum levels of cost and emission are satisfied under considered constraints [12]. When minimization of both environmental pollutants’ emission and energy cost comprise the objectives, the problem will be multi-objective, since these objective functions are conflicting. Multi-objective optimization approaches are developed for the sake of achieving Pareto optimal solutions of the conflicting objectives such that the operator is capable of making a trade-off among the set of optimal solutions. Some papers are devoted to the multi-objective economic∖emission optimal operation of MGs, applying different meta-heuristic algorithms [5 , 15–19]. A krill herd (KH) algorithm is suggested in [5] for stochastic optimal operation management of a grid-connected MG. In [15] optimization of the environmental economic problem in MG is considered and multi-objective mesh adaptive direct search is presented to minimize the total cost. A multi-objective bi-level optimal operation model for distribution network and MGs is suggested in [16] and the problem is solved using a self-adaptive genetic algorithm and nonlinear programming. In [17] a θ-PSO algorithm is applied to deal with the MG’s energy management problem. Authors in [18] proposed a θ-krill herd algorithm. A multi-objective PSO is used in [19] in order to optimize MG’s short-term performance. The superiority of an optimization algorithm in solving multi-objective problems is revealed from its robustness and fast convergence to a well-distributed Pareto-optimal front over the course of time. In comparison with other meta-heuristic algorithms, the grey wolf optimizer (GWO) algorithm, introduced by Mirjalili et al, manifests approximately proper search speed and convergence in solving some optimization problems [20]. In this paper, in order to investigate the multi-objective optimal operation management (MOOM) problem, a modified multi-objective GWO (MMOGWO) algorithm is proposed. Most of algorithms are needy to onerous tuning process of control parameters which is not required in MMOGWO that makes it superior among all other algorithms. Accordingly, the following objectives are met in the paper:

A novel algorithm is proposed comprising three modifications. The first modification is in the size of population which leads to a variable population. Consequently, trapping in local optima is avoided and the algorithm’s convergence speed is increased. Two other modifications are imposed in the mutations as is described in Section 4.3, which lead to the increase of the accuracy and convergence capability of the algorithm.

A novel method for controlling the size of repository is applied such that the algorithm’s speed improved while the search space becomes immense which leads to finding the optimum global (best compromised) solution faster and more precisely.

Since the ON/OFF states of DGs are taken in to account, a mixed-integer problem is solved in the article. Two different objective functions are assumed. An exquisite Pareto front of optimal solutions is achieved while the computational time is very low. The effectiveness of the proposed approach is affirmed by applying it on the typical MG of [12]. Additionally, the robustness of the algorithm, as its other outstanding feature, is highlighted through the simulation results.

2 Problem formulation

In the considered MOOM problem, optimal allocation of power generation set points besides proper On∖Off states of DG units are defined such that the objective functions, namely the operating cost and emission of the MG, are minimized while several constraints are satisfied [12].

2.1 Cost and emission minimization

In order to consider total operation cost as the first objective function the following should be satisfied:

$\begin{matrix} F_{1} (\vec{X}) = \sum_{t = 1}^{T} {Cost}^{t} = \sum_{t = 1}^{T} {\sum_{i = 1}^{N_{DG}} [u_{i}^{t} . (P_{DGi}^{t}) . B_{DGi}^{t} + \\ {SUC}_{DGi} . u_{i}^{t} . (1 - u_{i}^{t - 1}) + {SDC}_{DGi} . u_{i}^{t - 1} . (1 - u_{i}^{t})] \\ + \sum_{r = 1}^{N_{RES}} [P_{RE S_{r}}^{t} . B_{RE S_{r}}^{t}] + \sum_{s = 1}^{N_{Batt}} [P_{Bat t_{s}}^{t} . B_{{Batt}_{s}}^{t}] \\ + P_{Grid}^{t} . B_{Grid}^{t}} \end{matrix}$ (1)

$\begin{matrix} \vec{X} = {[{\vec{U}}_{DG}, {\vec{P}}_{DG}, {\vec{P}}_{Batt}, {\vec{P}}_{Grid}]}_{1 \times n} \\ n = [((2 \times N_{DG}) + N_{Batt}) + 1] \times T \end{matrix}$ (2)

${\vec{U}}_{DG} = [{\vec{u}}_{1}, {\vec{u}}_{2}, \dots, {\vec{u}}_{N_{DG}}] = {{\vec{u}}_{i}}_{1 \times N_{DG}} \in {0, 1}$ (3)

${\vec{u}}_{i} = [u_{i}^{1}, u_{i}^{2}, \dots, u_{i}^{t}, \dots, u_{i}^{T}]; \begin{matrix} i = 1, 2, \dots, N_{DG} \end{matrix}$ (4)

${\vec{P}}_{DG} = [{\vec{P}}_{D G_{1}}, {\vec{P}}_{D G_{2}}, \dots, {\vec{P}}_{D G_{N_{DG}}}]$ (5)

$\begin{matrix} {\vec{P}}_{DGi} = [P_{DGi}^{1}, P_{DGi}^{2}, \dots, P_{DGi}^{t}, \dots, {PD}_{Gi}^{T}]; \\ i = 1, 2, \dots, N_{DG} \end{matrix}$ (6)

${\vec{P}}_{Batt} = [{\vec{P}}_{Bat t_{1}}, {\vec{P}}_{Bat t_{2}}, \dots, {\vec{P}}_{Bat t_{N_{Batt}}}]$ (7)

$\begin{matrix} {\vec{P}}_{Batt, s} = [P_{Batt, s}^{1}, P_{Batt, s}^{2}, \dots, P_{Batt, s}^{t}, \dots, P_{Batt, s}^{T}]; \\ s = 1, 2, \dots, N_{Batt} \end{matrix}$ (8)

${\vec{P}}_{RES} = [{\vec{P}}_{RE S_{1}}, {\vec{P}}_{RE S_{2}}, \dots, {\vec{P}}_{RE S_{N_{RES}}}]$ (9)

$\begin{matrix} {\vec{P}}_{RE S_{r}} = [P_{RE S_{r}}^{1}, P_{RE S_{r}}^{2}, \dots, P_{RE S_{r}}^{t}, \dots, P_{RE S_{r}}^{T}]; \\ r = 1, 2, \dots, N_{RES} \end{matrix}$ (10)

${\vec{P}}_{Grid} = [P_{Grid}^{1}, P_{Grid}^{2}, \dots, P_{Grid}^{t}, \dots, P_{Grid}^{T}]$ (11)

in which $\vec{X}$ is the vector of design variables, and n is the number of design variables. T is the total number of hours. The total number of dispatchable generations and storage units (battery) along with the number of RESs are N_DG, N_Batt and N_RES, respectively. $u_{i}^{t}$ is applied to imply the ON/OFF states of ith dispatchable DG during each hour of the day. $P_{DGi}^{t}$ , $P_{RE S_{r}}^{t}$ and $P_{Batt, s}^{t}$ represent the real output powers (kWh) of the ith DG, rth RES and sth storage at time t, respectively. The active power which is bought (sold) from (to) the utility at time t is demonstrated by $P_{Grid}^{t}$ , $B_{DGi}^{t}$ , $B_{{RES}_{r}}^{t}$ , $B_{Batt, s}^{t}$ and $B_{Grid}^{t}$ are respectively bids of dispatchable DGs, RESs, storage devices and the utility grid at hour t (€ /kWh). SUC_DGi and SDC_DGi are the start-up and shut-down cost for ith dispatchable DG. $(P_{DGi}^{t} . B_{DGi}^{t})$ , $(P_{Bat t_{s}}^{t} . B_{Bat t_{s}}^{t})$ and $(P_{Grid}^{t} . B_{Grid}^{t})$ represent operational cost of dispatchable DGs, battery and cost of power exchange between the MG and utility (€), respectively. It should be mentioned that DG demonstrates the dispatchable units including fuel cell (FC) and micro-turbine (MT), while Grid and Batt are abbreviated forms of the utility grid and the battery, respectively. The WT and PV are shown with RES.

Note that the utility has to buy all electrical power produced by RES units, consequently RESs’ output powers ( $P_{RES, r}^{t}$ ) are not included in the design variables’ vector.

As the second objective, the environmental pollutants should be minimized as the following [12]:

$\begin{array}{l} F_{2} (\vec{X}) = \sum_{t = 1}^{T} E m i s s i o n^{t} = \sum_{t = 1}^{T} { \sum_{i = 1}^{N_{D G}} [u_{i}^{t} . (P_{D G i}^{t}) . E_{G i}^{t}] + \\ \sum_{s = 1}^{N_{B a t t}} [P_{B a t t_{s}}^{t} . E_{B a t t_{s}}^{t}] + P_{G r i d}^{t} . E_{G r i d}^{t}} \end{array}$ (12)

where $E_{DGi}^{t}, E_{Bat t_{s}}^{t}$ and $E_{Grid}^{t}$ are the amount of pollutants emission (kg/kWh) for each generator, storage device and utility at hour t, respectively. These variables are described as follows [12 –19]: $\begin{matrix} E_{DGi}^{t} = {CO}_{2_{D G_{i}}}^{t} + {SO}_{2_{D G_{i}}}^{t} + {NO}_{x_{D G_{i}}}^{t} \\ E_{Bat t_{s}}^{t} = {CO}_{2_{Batt, s}}^{t} + {SO}_{2_{Batt, s}}^{t} + {NO}_{x_{Batt, s}}^{t} \\ E_{Grid}^{t} = {CO}_{2_{Grid}}^{t} + {SO}_{2_{Grid}}^{t} + {NO}_{x_{Grid}}^{t} \end{matrix}$ (13)

where ${CO}_{2_{D G_{i}}}^{t}$ , ${SO}_{2_{D G_{i}}}^{t}$ and ${NO}_{x_{D G_{i}}}^{t}$ are the amounts of CO₂, SO₂ and NO_x emission from ith DG source at hour t; ${CO}_{2_{Batt, s}}^{t}$ , ${SO}_{2_{Batt, s}}^{t}$ and ${NO}_{x_{Batt, s}}^{t}$ are the amounts of CO₂, SO₂ and NO_x emission from the sth storage unit at hour t of the day, and ${CO}_{2_{Grid}}^{t}$ , ${SO}_{2_{Grid}}^{t}$ and ${NO}_{x_{Grid}}^{t}$ are the amounts of CO₂, SO₂ and NO_x emission from the utility at hour t, respectively.

2.2 Constraints

- Power balance

One of the most important requirements in MG management is the balance of electricity demand and supply, hence [12]: $\sum_{i = 1}^{N_{DG}} [u_{i}^{t} . P_{DGi}^{t}] + \sum_{r = 1}^{N_{RES}} [P_{RE S_{r}}^{t}] + \sum_{s = 1}^{N_{Batt}} [P_{Bat t_{s}}^{t}] + P_{Grid}^{t} = P_{LD}^{t}$ (14)

where $P_{LD}^{t}$ is the total MG load at hour t.

- Battery limits

In order to consider the limitation on charge and discharge rates of the storage devices during each time intervals, along with limits on the state of charge (SOC) of each storage device the following equation and constraints are mentioned for a typical battery [21]: $W_{{Batt}_{s, min}} \leq W_{{Batt}_{s}}^{t} \leq W_{{Batt}_{s, max}}$ (15)

where is the SOC of the sth storage at the end of one-hour interval which is associated to the time t as the following: $\begin{matrix} W_{Batt}^{t} = W_{Batt}^{t - 1} + η_{c} P_{ch} . Δ t - \frac{1}{η_{d}} P_{disch} . Δ t \\ P_{ch}^{t} \leq P_{ch, max}; \\ P_{disch}^{t} \leq P_{disch, max} \end{matrix}$ (16)

in which $W_{Batt}^{t}$ and $W_{Batt}^{t - 1}$ are the amounts of energy stored inside the battery at hours t and t - 1, respectively. During a definite period of time (Δt), P_ch (P_disch) is the permitted rate of charge (discharge), while η_c (η_d) is the efficiency of the battery during charge (discharge) process. The lower and upper limits of amounts of energy storage inside the battery are W_{Batt
_s,min} and W_{Batt
_s,max}, and P_ch,max (P_disch,max) is the maximum rate of battery charge (discharge) during each time interval Δt.

- Real power constraint

Power generations for each dispatchable DG are limited as: $u_{i}^{t} . P_{D G_{i, min}}^{t} \leq P_{D G_{i}}^{t} \leq u_{i}^{t} . P_{D G_{i, max}}^{t}$ (17)

The power exchange with utility grid is constrained as follows: $P_{Grid, min}^{t} \leq P_{Grid}^{t} \leq P_{Grid, max}^{t}$ (18)

Constraints on the rate of charge and discharge of the battery during an hour are considered as the following: $P_{Bat t_{s}, min}^{t} \leq P_{{Batt}_{s}}^{t} \leq P_{{Batt}_{s}, max}^{t}$ (19)

3 Fundamentals of multi-objective optimization

3.1 Characterization of multi-objective optimization

In a typical multi-objective problem, a number of objective functions are simultaneously optimized. In most of the situations, these objective functions are conflict which is a barrier to select an optimum solution for all the objective functions. In most cases, multi-objective problems have more than one optimal solution, which are called the non-dominated solutions. Within the whole search space, the non-dominated solutions are expressed as Pareto-optimal which establish the Pareto-optimal set or Pareto-optimal front. Considering a multi-objective minimization problem, while meeting a number of equality and inequality constraints, a solution X₁ dominates X₂ if the following conditions are satisfied [13, 14]: $\begin{matrix} \forall j \in {1, 2, \dots, n}, f_{j} (X_{1}) \leq f_{j} (X_{2}) \\ \exists k \in {1, 2, \dots, n}, f_{k} (X_{1}) < f_{k} (X_{2}) \end{matrix}$ (20)

where f_i (X) is the ith objective function, X is the vector of the optimization variables, and n is the number of objective functions [13, 14]. In order to assess the suitable particles to be stored in the repository of non-dominated solutions, the concept of Pareto-optimal is employed [13, 14].

3.2 Controlling the size of repository

Since the size of repository is limited, a finite number of solutions can be accumulated. In this paper, repository size is controlled through a novel distance-based method.

When the number of non-dominated solutions in the repository exceeds a predefined value, namely N_L, they are sorted ascending according to one of the objective functions, and the first and the last of the sorted non-dominated solutions are assumed as A and B, respectively. The pseudo-code of the proposed method for controlling the size of repository is as shown in Table 1. In this table N_non-dom is the number of non-dominated solutions, f_i is the objective function and M is the number of objective functions.

Table 1
The pseudo-code of the proposed method for controlling the size of repositor

$epsilon = \frac{\sqrt{\sum_{m = 1}^{M} {(f_{A_{m}} - f_{B_{m}})}^{2}}}{N_{L}}$

Fori = 2 toN_non-dom

$distance = \sqrt{\sum_{m = 1}^{M} {(f_{1_{m}} - f_{i_{m}})}^{2}}$

End

l = 1

Fori = 1 toN_L

κ= find((I-1)epsilon < distance and distance <iepsilon)

If length κ∼= 0

Y =κ (end)

l = l+1

End

End

3.3 Selecting the best compromised solution

As was mentioned, the Pareto front should be scattered uniformly. Additionally, considering each objective function to be minimized independently, the best solution of each should be obtained. In order to achieve these goals, particles should be refrained to accumulate in populated domains. Hence, the proposed method of [21] is applied to improve the properties of Pareto optimal front and to select the global best compromised solution. In this technique, objective functions are normalized such that relatively equal significance is provided to the both objectives as follows:

$M i n φ^{M} = ω_{1} [\frac{\sum_{i = 1}^{N_{_{n o n - d o m}}} f 1_{i} - f_{1}^{^{\min}}}{f_{1}^{^{\max}} - f_{1}^{^{\min}}}] + ω_{2} [\frac{\sum_{i = 1}^{N_{_{n o n - d o m}}} f 2_{i} - f_{2}^{^{\min}}}{f_{2}^{^{\max}} - f_{2}^{^{\min}}}]$ (21)

where f₁ and f₂ are objective functions cost and emission, respectively. Initial guess for ω₁ and ω₂ is equal to 0.5.

4 Modified GWO algorithm

4.1 A brief overview on original GWO algorithm

A swarm-intelligence algorithm impressed by the hierarchical hunting manner of grey wolves, namely GWO, is presented by Mirjalili et al. [20]. The social hierarchy of grey wolves is classified in four groups including the group leader, alpha (α), as the first level, beta (β) as the second level who cooperates alpha in decision making, while the third best solution is called delta (δ) which comes after alpha and beta but leads the fourth hierarchy, omega (ω), which should defer to other three dominant levels. The hunting mechanism of grey wolves is as follows [20]:

Tracking, chasing, and approaching the prey.

Pursuing, encircling, and harassing the prey until it stops moving.

Attacking the prey

In the optimization process ωs randomly update their positions around the prey according to the estimated position of the prey by α, β and δ. The encircling is then performed as follows:

$D = | C . X_{p, k} - X_{k} |$ (22)

$X_{k + 1} = X_{p, k} - A . D$ (23)

k indicates the current iteration, X_p,k is the position vector of the prey, X_k is the position vector of a grey wolf, D is the distance between X_p,k and X_k, while A and C represent the discrimination weight coefficient of search agent, and random mutation coefficients, respectively and are calculated as following:

$A = (2 r_{1} - 1) a$ (24)

$C = 2 r_{2}$ (25)

where r₁ and r₂ are random vectors in [0, 1], a is linearly decreased from 2 to 0 over the course of iteration.

Hunting the prey as the last step of the procedure is guided by α, while β and δ participate in this step. α, β and δ have better knowledge about the potential position of the prey. The position of each grey wolf is updated as follows:

$D_{α} = | C_{1} . X_{α, k} - X_{k} |$ (26)

$D_{β} = | C_{2} . X_{β, k} - X_{k} |$ (27)

$D_{δ} = | C_{3} . X_{δ, k} - X_{k} |$ (28)

$\begin{matrix} X_{1} = X_{α} - A_{1} . (D_{α}), X_{2} = X_{β} - A_{2} . (D_{β}), \\ X_{3} = X_{δ} - A_{3} . (D_{δ}) \end{matrix}$ (29)

$X_{k + 1} = \frac{1}{3} (X_{1} + X_{2} + X_{3})$ (30)

where A₁, A₂ and A₃ are respectively decisive weight coefficients of α, β and δ, values of which decrease progressively.

4.2 Multi-objective GWO (MGWO)

In this paper in order to investigate a bi-objective problem to optimize both the total cost and emission simultaneously, a multi-objective GWO (MOGWO) should be applied. Accordingly, X_α, X_β and X_δ are not fixed for each member and are selected randomly from the repository as is described. If the number of non-dominated solutions in the repository is more than three, after being sorted according to one of the objective functions the repository is divided into three equal sections. X_α is selected randomly from the first section which is from the first non-dominated solution in the repository up to the ${\frac{N_{L}}{3}}^{t} h$ , X_β is selected from ${(\frac{N_{L}}{3} + 1)}^{t} h$ non-dominated solution up to ${\frac{2 N_{L}}{3}}^{t} h$ non-dominated solution, and X_δ is chosen from the last section which is from the ${(\frac{2 N_{L}}{3} + 1)}^{t} h$ non-dominated solution to the last one.

This procedure leads to an intelligent selection of X_α, X_β and X_δ where all three sections of the repository take part in the output Pareto front and consequently the actual Pareto optimal set is obtained.

4.3 Modified multi-objective GWO algorthim (MMOGWO)

In order to achieve a real Pareto optimal front in the MOOM problem, the original GWO should be enhanced. Consequently, three modifications, one in the size of population and the two others in the mutations, are augmented to improve the convergence ability and the accuracy of the approach.

- Modification I

To increase the convergence speed of the algorithm, the size of population is considered variable and changes as follows:

$N = round (\frac{(N_{max} - N_{min}) \times iteration}{ite r_{max}} + N_{min})$ (31)

where N_min and N_max are the minimum and maximum populations respectively, and iter_max is the maximum number of iteration. A prominent feature of the proposed modification is that the selection of a variable population size leads to the increase of the population in each iteration, therefore, the population size is not fixed which helps to avoid trapping in local optima. As a result, the accuracy and convergence capability of the algorithm will improve.

- Modification II

In addition to the described modification in the previous section, the second modification is applied to improve the accuracy of the proposed approach. Five constants k₁ ≠ k₂ ≠ k₃ ≠ k₄ ≠ k₅, unequal to i, are chosen randomly from the population, while three constants $k_{1}^{'} \neq k_{2}^{'} \neq k_{3}^{'}$ are selected from the repository as follows: If the number of non-dominated solutions in the repository is more than three, after being sorted according to one of the objective functions the repository is divided into three sections. $k_{1}^{'}$ is selected randomly from the first section which is from the first non-dominated solution in the repository to the ${\frac{N_{L}}{3}}^{t} h$ , $k_{2}^{'}$ is selected from ${(\frac{N_{L}}{3} + 1)}^{t} h$ non-dominated solution up to ${\frac{2 N_{L}}{3}}^{t} h$ non-dominated solution, and $k_{3}^{'}$ is chosen from the last section which is from the ${(\frac{2 N_{L}}{3} + 1)}^{t} h$ non-dominated solution to the last one. Four mutations $({\vec{X}}_{mu t_{l}}, l = 1, 2, 3, 4)$ are defined as:

${\vec{X}}_{mut 1} = {\vec{X}}_{k_{1}} + rand 1 (.) \times ({\vec{X}}_{k_{2}} - {\vec{X}}_{k_{3}})$ (32)

${\vec{X}}_{mut 2} = {\vec{X}}_{mut 1} + rand 2 (.) \times ({\vec{X}}_{best} - {\vec{X}}_{worst})$ (33)

${\vec{X}}_{mut 3} = {\vec{X}}_{k_{4}} + rand 3 (.) \times ({\vec{X}}_{best} - {\vec{X}}_{k_{5}})$ (34)

${\vec{X}}_{mut 4} = {\vec{X}}_{{k^{'}}_{1}} + rand 4 (.) \times ({\vec{X}}_{{k^{'}}_{2}} - {\vec{X}}_{{k^{'}}_{3}})$ (35)

where after sorting the repository, the first non-dominated solution is selected as X_best and in each iteration the best non-dominated solutions are selected in turn as X_best such that all repository members take part in the population generation, while X_worst is the most dominated solution in the population.

- Modification III

In order to increase the convergence speed of the algorithm the third modification is augmented according to the pseudo code of Table 2.

Table 2

The pseudo-code of modification III

iter = 1

Whileiter ≤ iter_max

For i = 1 to N

1 \leq iter \leq \frac{1}{3} {iter}_{\max}

Select one of the mutations of Equation (21) randomly and apply it to the ith member

IfX_new1 dominates X_i,iter, or non dominates each other

SaveX_new1 in the repository

SetX_i,iter = X_new1

End

Else if

\frac{1}{3} {iter}_{\max} + 1 \leq iter \leq \frac{2}{3} {iter}_{\max}

Select two of the mutations of Equation (21) randomly and apply them to the ith member

Determine the non-dominated solution among X_new1 and X_new2 and X_i,iter

Save the non-dominated solution in the repository

IfX_new1 and X_new2 are the non-dominated solutions

Choose one of them based on Equation (12) and update X_i,iter

Else ifX_new1 or X_new2 is the non-dominated solution

SetX_i,iter = X_new,j (j = 1 or 2)

End

Else if

\frac{2}{3} {iter}_{\max} + 1 \leq iter \leq {iter}_{\max}

Apply all four mutations

Determine the non-dominated solution among X_new1, X_new2, X_new3, X_new4, and X_i,iter

Save the non-dominate solution in the repository

If the number of non-dominated solutions among X_new1, X_new2, X_new3, X_new4 is more than one

Choose one of them based on Equation (12) and update X_i,iter

Else if one of X_new,j (j = 1 to 4), is the non-dominated solution

SetX_i,iter = X_new,j (j = 1 or 2 or 3 or 4)

End

iter = iter+1

End

5 Application of the proposed method

In order to apply the proposed algorithm on the MOOM problem in an MG, the following steps should be taken.

Step 1. Initialize population size, number of design variables and termination criterion. Problem information including MG properties, beside bids and power information of DGs, storages and utility, hourly wind-turbine (WT) and photovoltaic (PV) power forecasts, emission coefficients are specified. The initial charge of the battery is also defined in this step.

Step 2. Since a mixed integer problem is considered in this paper, two types of variables, binary and continuous, are assumed. For states of generators as binary variables, U_i s according to Equation (2) are generated as follows:

$u_{i}^{t} = round (rand (.) \times (u_{i, max}^{t} - u_{i, min}^{t}) + u_{i, min}^{t})$ (36)

However, in order to consider the states of all units, $U_{i}^{t}$ should satisfy the following condition for all hours:

$\begin{matrix} \sum_{i = 1}^{N_{DG}} [u_{i}^{t} . P_{D G_{i, min}}^{t}] + \sum_{r = 1}^{N_{RES}} [P_{RE S_{r, min}}^{t}] + \\ \sum_{s = 1}^{N_{Batt}} [P_{Bat t_{s, min}}^{t}] + P_{Gri d_{min}}^{t} \leq P_{Load}^{t} \\ \leq \sum_{i = 1}^{N_{DG}} [u_{i}^{t} . P_{D G_{i, max}}^{t}] + \sum_{r = 1}^{N_{RES}} [P_{RE S_{r, max}}^{t}] + \\ \sum_{s = 1}^{N_{Batt}} [P_{Bat t_{s, max}}^{t}] + P_{Gri d_{max}}^{t} \end{matrix}$ (37)

If Equation (23) is satisfied, a random population for continuous variables based on the achieved U_i s and according to Equation (2) must be generated as Equation (24).

$P_{D G_{i}}^{t} = rand 5 (.) \times (P_{D G_{i, max}}^{t} - P_{D G_{i, min}}^{t}) + P_{D G_{i, min}}^{t}$ (38) $P_{Bat t_{s}}^{t} = rand 6 (.) \times (P_{Bat t_{s, max}}^{t} - P_{Bat t_{s, min}}^{t}) + P_{Bat t_{s, min}}^{t}$ (39)

$P_{Grid} = rand 7 (.) \times (P_{Grid, max} - P_{Grid, min}) + P_{Grid, min}$ (40)

while the power constraints in Equations (8–10) should be satisfied.

Since in the considered problem some limitations, such as battery constraints, depend on previous and future hours, constraints change in different hours of the day as follows:

$\begin{matrix} W_{disch, Batt}^{t + 1} = max {(\frac{W_{Batt}^{t} - Δ {tP}_{Batt}^{t}}{η_{d}}), W_{Batt, min}}, \\ t = 1, \dots, T \end{matrix}$ (41)

$\begin{matrix} W_{ch, Batt}^{t + 1} = min {(W_{Batt}^{t} - Δ {tP}_{Batt}^{t}, η_{c}), W_{Batt, max}}, \\ t = 1, \dots, T \end{matrix}$ (42)

$\begin{matrix} P_{ch, max} \leq P_{Batt}^{t} \leq P_{disch, max}, \\ t = 1, \dots, T \end{matrix}$ (43)

$\begin{matrix} P_{disch, max}^{t} = \min {(P_{Batt, max}, \frac{(W_{Batt}^{t} - W_{Batt, min}) η_{d}}{Δ t})}, \\ t = 1, \dots, T \end{matrix}$ (44)

$\begin{matrix} P_{ch, max}^{t} = max {(P_{Batt, min}, \frac{(W_{Batt}^{t} - W_{Batt, max})}{(η_{c} Δ t)})}, \\ t = 1, \dots, T \end{matrix}$ (45)

The limitation of the released energy in the discharging mode and the stored energy in charging mode are mentioned in Equations (25 and 26), respectively, while Equation (28) demonstrates the power discharged by the battery in the discharge mode, and the power charged by the grid in the charging mode. Using Equations (28 and 29), the maximum and minimum rates of charging and discharging mode can be calculated [22].

Step 3. Check the power balance violation for each particle. In order to meet the equality constraint, following steps should be carried out:

Step 3.1. t = 1

Step 3.2. Calculate power balance violation:

$\begin{matrix} P_{VLTD}^{t} = \\ (\sum_{i = 1}^{N_{DG}} [u_{i}^{t} . P_{DGi}^{t}] + \sum_{s = 1}^{N_{Batt}} [P_{Batt}^{t}] + P_{Grid}^{t}) - \sum_{d = 1}^{N_{d}} P_{LD}^{t} \end{matrix}$ (46)

Step 3.3. If $P_{VLTD}^{t} ≅ 0$ and t ≤ T then set t = t + 1 and go to Step 3.2; otherwise

If $P_{VLTD}^{t} ≅ 0$ and t = T, go to Step 4; otherwise

If $P_{VLTD}^{t} \neq 0$ select a generated unit randomly (battery, DGs or grid).

Step 3.4. Add $P_{VLTD}^{t}$ to the power of the selected generated unit.

Step 3.5. Check the result with upper and lower limits of the units’ powers, if it violates the upper limit, fix it to the upper limit and if it violates the lower limit, fix it to the lower limit. Go to Step 3.2.

Step 4. Calculate the objective functions for the initial population.

Step 5. Determine the non-dominated solutions in the population and store them in the repository.

Step 6. Choose ‘X_best’ and ‘X_worst’ as described in Section 4.3.

Step 7. Initialize the size of population based on Equation (20).

Step 8. Select ith individual from population.

Step 9. Apply GWO according to Equations (13–19).

Step 10. Apply mutations as described in Section 4.3.

Step 11. If all members are selected, go to Step 13, otherwise set i = i + 1 and got to Step 9.

Step 12. Update the number of population according to Equation (20).

Step 13. Determine the non-dominated solutions in the new population.

Step 14. Update the repository based on new and old non-dominated solutions.

Step 15. If the number of non-dominated solutions is more than a predefined value, N_L, control the size of repository as mentioned in Section 3.2.

Step 16. Control the termination criterion, if it is satisfied, terminate the algorithm, otherwise, set iteration = iteration+1 and turn to Step 6.

Step 17. Report the best compromised solution which is achieved using the technique proposed in Section 3.3.

The flowchart of the proposed algorithm is illustrated in Fig. 1.

6 Simulation results

Fig.1

Flowchart of the proposed algorithm.

The effectiveness of the proposed MMOGWO algorithm is verified in this section where two cases are investigated. The MG of Fig. 2 is considered as the test system. A 24-hour scheduling scheme is assumed for the analysis of the simulated system in order to clarify the performance of each power unit. Besides, the unity power factor is considered for all DGs, thus they just produce active power. The decision about power exchange between the MG and the utility, which is allowed at any hour in a day in order to more profitably exploit the market, is taken by MG central controller (MGCC). The data for the hourly active power of PV and WT, forecasted load demand and the utility power production bid, besides the entire bid data for all DGs along with the power market are available in [12]. PV and WT units do not consume any fuel at the times they produce electrical power during the day, consequently, the utility has to buy all electrical power produced by these units [12]. A Pareto-optimal set is attained for the two incompatible objectives (cost and emission) in each case. It is worth mentioning that applying the proposed approach for controlling the size of repository leads to the extraordinary fast convergence of the proposed MMOGWO that makes it superior among all other existing algorithms. The proposed method was implemented in MATLAB 8.1 and solved in a laptop with Core i5 CPU and 4GB RAM. The number of population and maximum iteration are both considered 100.

Fig.2

A typical MG test system [12].

6.1 First case

In this section, it is assumed that PV and WT (as non-dispatchable units) are in service at all hours during the day and are accepted to exploit at their maximum available output powers. This is while the ON∖OFF states of dispatchable DGs (i.e. FC and MT) are considered. Therefore, in the algorithm process, the solutions for dispatchable units are compared with their minimum limit powers and if lower, the power will be put to zero. In this case no limits on the battery’s initial charge is considered. The Pareto optimal front of this case is revealed in Fig. 3, and results for three algorithms, including original GWO, PSO, and the proposed MMOGWO, are compared. The best compromised solution along with the values of the points where cost and emission are minimum are depicted. It is observed that in solving the MOOM problem, in addition to fast convergence, the proposed MMOGWO algorithm is able to properly find the points where the objective functions are in their minimum values, and the Pareto optimal front maintains between these points, while one deficiency of the two other algorithms in dealing with the considered MG energy management problem is the incapability to detect these points. Evidently, when the cost function is at its minimum value, 281.3 €, the emission is 862.8 kg. Besides, when emission decreases to 455.9 kg, cost equals 857.2 €. The best compromised solution is calculated according to the procedure described in section 3.3 where the cost and emission objective functions are 373.4 € and 566 kg, respectively.

In order to justify the robustness and effectiveness of the proposed MMOGWO algorithm, the program is executed four times and results of these four different executions in the first case are revealed in Fig. 4. Obviously, the achieved Pareto optimal sets in all four runs are approximately similar. It can be concluded that since the numbers of non-dominated solutions saved in the repository in different runs of the program, which are 150, 150, 147 and 149 respectively for the first, second, third and fourth runs, are very close, the proposed algorithm is robust. Thereupon, not only can the proposed algorithm increase the convergence speed but also it decreases quiescence which results in getting away from local optimums. Consequently, the proposed MMOGWO proposes more robust and qualified solutions.

Fig.3

Comparison of the objective functions Pareto optimal fronts of MMOGWO, original GWO and PSO algorithms for the first case

Fig.4

Four different executions of MMOGWO algorithm (first case).

The power dispatch in the best compromised solution obtained using MMOGWO in this case is presented in Tables 3 in details. From these results, it is concluded that all equality and inequality constraints are satisfied. In the dispatch of the battery, when the battery is charging the values of power are negative, while during the discharging hours the values of power are positive. However, for the utility the negative values are representative of delivering energy to the upstream network, while the positive values are related to the times when energy is purchased from the upstream network. According to Table 3, since there is no limit on the battery charge, and the bid of the battery is less than other units, it sells power up to its maximum value in all hours of the day. Additionally, since the price of FC is less than MT, the purchased power from FC in different hours of the day is more. The MT power is limited on its permissible minimum value in most of hours because of the high price of power in comparison with other units.

Table 3

MMOGWO best compromised solution of the first and second cases

Power (kWh)	First case				Second case
	MT	FC	Battery	Utility	MT	FC	Battery	Utility
hour
1	6.027	13.825	30	0.364	30	30	−10.2	0.5
2	6	12.232	30	−0.017	29.8	30	−12.3	0.7
3	6.184	12.8	30	−0.769	29.8	30	−16.1	4.6
4	6	13.22	30	−0.005	29.8	30	−21.6	11
5	6	18.205	29.997	0.013	29.8	29.7	−20.8	15.6
6	6.258	26.024	30	−0.197	30	30	−15.6	17.7
7	9.398	29.184	30	−0.366	29.8	29.9	−0.4	8.9
8	13.513	30	30	−0.018	29.5	30	13.1	0.9
9	10.708	29.995	29.931	−0.169	30	30	10.4	0.1
10	30	30	30	−20.62	30	30	30	−20.6
11	28.776	30	29.999	−30	30	30	−0.5	−0.7
12	21.64	30	30	−30	30	29.9	20.4	−28.7
13	6	29.951	29.954	−21.72	29.8	30	−15.5	−0.1
14	18.58	30	30	−30	30	30	18.6	−30
15	29.974	30	30	−23.63	29.9	30	6.3	0.1
16	29.996	29.991	29.999	−15.52	30	30	14.1	0.4
17	22.678	29.999	29.987	0.0004	29.9	30	0	22.7
18	26.628	30	30	−0.413	30	30	0	26.2
19	28.648	29.999	29.99	0.0601	30	30	0	28.7
20	26.174	29.882	29.928	−0.768	30	30	0	25.2
21	16.699	30	30	0	30	30	0	16.7
22	9.706	30	30	−0.007	30	30	0	9.7
23	6	28.094	29.999	−0.008	29.9	30	0	4.2
24	6.046	19.381	29.974	−0.016	25.6	29.8	0	−0.1

As is expected, it is evidently obvious from Table 4 that when cost is the only considered objective function, it is significantly lower (which is equal to 269.85 €;) in comparison with the case that the ON/OFF states of dispatchable units is not considered where it is 278.25 €. According to Table 4, in the first 8 hours of the day, it is more economical when MT is in OFF mode.

Table 4

MMOGWO optimal solution achieved in the first case where the cost is reduced to 269.85 €

Power (kWh)
	MT	FC	Battery	Utility
hour
1	0	20.171	0.0445	30
2	0	17.924	0.2911	30
3	0	18.145	0.0699	30
4	0	19.129	0.1158	29.9706
5	0	24.159	0.0881	29.9679
6	0	29.906	2.1815	29.9976
7	0	29.992	8.3458	29.8773
8	0	30	18.6006	24.8944
9	29.991	29.999	30	−19.5251
10	30	30	30	−20.615
11	28.817	30	29.9576	−30
12	21.787	29.856	29.9964	−29.9997
13	14.186	29.999	30	−29.9992
14	18.729	29.855	29.9959	−29.9997
15	29.993	30	30	−23.6533
16	30	30	29.9995	−15.5295
17	30	30	30	−7.335
18	0	30	30	26.215
19	0	30	28.7922	29.9058
20	26.429	30	30	−1.2142
21	30	30	30	−13.3005
22	29.905	29.942	29.9868	−20.1343
23	0	30	4.085	30
24	0	25.188	0.2147	29.982

6.2 Second case

In the second scenario the battery initial charge is expected to be zero. The Pareto optimal solutions for this case are illustrated in Fig. 5, where the best compromised solution along with the points where the two objective functions are minimum for the proposed MMOGWO are revealed, besides a comparison between three algorithms, original GWO, PSO, and MMOGWO, is performed. The minimum cost is equal to 484 €, while emission is 1114 kg; in Emission-minimum, the emission is reduced to 950.9 kg, while cost equals 944.4 €. The conflict of the two objectives can be concluded from these obtained values. Consequently, the detailed dispatch of the best compromised solution for the second case is tabulated in Table 3. In this case, the best compromised solution is when the cost and emission objective functions are 643.7 € and 1026 kg, respectively. It is obvious from Table 3 that in the first 7 hours of the day, the local load is supplied by purchasing from the utility grid and even from MT, while the energy surplus is stored in the battery. At hour 8, the battery starts to discharge for the load being supplied, where in this case, the energy surplus is sold to the upstream network. At hour 13 energy is purchased from the upstream network since the market price is low and as the market price is at maximum value at hour 14 this stored energy is sold to the utility at the next hour. Then, at hour 16 the battery is completely discharged. From hour 17 to 24 during which the market price is approximately low, the power is purchased from the upstream network and the battery is not charged nor discharged. Since the price of FC is less than other units, it sells power up to its maximum value in all hours of the day.

Fig.5

Comparison of the objective functions Pareto optimal fronts of MMOGWO, original GWO and PSO algorithms for the second case.

In Fig. 6a comparison between four executions of the algorithm in the second case is shown. Obviously, the achieved Pareto optimal sets in all four runs are approximately similar. As is obvious from Fig. 6, the numbers of non-dominated solutions saved in the repository in different runs of the program, which are 180, 178, 179 and 177 respectively for the first, second, third and fourth runs, are very close. Consequently, the robustness of the proposed MMOGWO algorithm in the second case can be justified through the results of Fig. 6.

Fig.6

Four different executions of MMOGWO algorithm (second case).

When the initial charge of the battery is zero, a strict restriction on the battery charge and discharge is put on. Consequently, the most adjustable scenario is when the battery has initial charge. Because of the more realistic performance of the battery in the second case, more power is purchased from MT, comparatively. In both cases, in most hours a large amount of power is bought from FC since it is less expensive. In first hours of the day where the market price is lower, the system operator purchases power from the utility. This power can be utilized in order to supply local loads or can be stored in the storage device. However, in peak hours, the stored power can be sold to the utility in a much higher price.

The superiority of an optimization algorithm in solving multi-objective problems is concluded from the appropriate and fast convergence, while attaining exact Pareto front. Moreover, the robustness of an optimization algorithm can be another criterion in order to prove the effectiveness of the method. Consequently, as is concluded from the obtained results and by comparing the Pareto optimal sets of three algorithms in solving the MOOM problem, the priority of the proposed MMOGWO algorithm along with its robustness and accuracy is justified. Accordingly, it is proved that the proposed MMOGWO has successfully fulfilled these assumed criteria, and it can significantly achieve the exceptional solution in comparison with other methods in solving MOOM problem of an MG.0

7 Conclusions

The new MMOGWO algorithm was proposed in order to deal with the multi-objective optimal operation management in a typical MG. Three modifications were added to the original GWO which result in more accurate and faster performance of the suggested approach. The variable population size resulted from the first modification avoids trapping in local optima, while the two other modifications were augmented to the mutation procedure in order to increase the convergence speed and robustness of the algorithm. Additionally, a novel method was applied for controlling the size of repository which generates a well-distributed Pareto front in a very low computational time. As a result, the accuracy and speed of the algorithm will improve. It is obvious from the results that an exceptional Pareto optimal set was achieved while applying the presented method comparing with original GWO and PSO algorithms. Two different scenarios were considered in order to justify the effectiveness of MMOGWO. Simulation results manifest that the proposed method is able to deal with mixed-integer problems. Future works can include the following:

Investigating the stochastic MG optimal operation management while considering the uncertainties of renewable resources using the proposed algorithm.

Studying the effects of elements of the future smart grid, such as electric vehicles, in the considered MG energy management.

Footnotes

Acknowledgments

This work was supported in part by Royal Academy of Engineering Distinguished Visiting Fellowship under Grant DVF1617∖6∖45.

References

Vosoogh

, Kamyar

, Akbari

, A novel modification approach based on MTLBO algorithm for optimal management of renewable micro-grids in power systems, J Intell Fuzzy Syst 27(1) (2014), 465–473.

Baghaee

H.R.

, Mirsalim

, Gharehpetian

G.B.

, Multi-objective optimal power management and sizing of a reliable wind/PV microgrid with hydrogen energy storage using MOPSO, J Intell Fuzzy Syst 32(3) (2017), 1753–1773.

, Wu

, Bo

, Liu

, Zhou

, Chen

et al., Modeling, planning and optimal energy management of combined cooling, heating and power microgrid: A review, , Int J Electr Power Energy Syst 54 (2014), 26–37.

Javidtash

, Jabbari

, Niknam

, Nafar

, A novel mixture of non-dominated sorting genetic algorithm and fuzzy method to multi-objective placement of distributed generations in Microgrids, J Intell Fuzzy Syst (2017), 1–8.

Kamankesh

H.R.

, Vassilios

G.A.

, A sufficient stochastic framework for optimal operation of micro-grids considering high penetration of renewable energy sources and electric vehicles, J Intell Fuzzy Syst 32(1) (2017), 373–387.

Abbasi

A.R.

, Abbasi

, Ansari

, Rahmani

, Effect of plug-in electric vehicles demand on the renewable micro-grids, J Intell Fuzzy Syst 29(5) (2015), 1957–1966.

Nguyen

T.A.

, Crow

M.L.

, Stochastic optimization of renewable-based microgrid operation incorporating battery operating cost, IEEE Trans Power Syst 31 (2016), 2289–2296.

K.W.

, Liaw

C.M.

, Incorporated operation control of DC microgrid and electric vehicle, IEEE Trans Indust Electron 63 (2016), 202–215.

Marzband

, Yousefnejad

, Sumper

, Domínguez-García

J.L.

, Real time experimental implementation of optimum energy management system in standalone microgrid by using multi-layer ant colony optimization, Int J Electr Power Energy Syst 75 (2016), 265–274.

10.

Jabbari-Sabet

, Moghaddas-Tafreshi

S.-M.

, Mirhoseini

S.-S.

, Microgrid operation and management using probabilistic reconfiguration and unit commitment, Int J Electr Power Energy Syst 75 (2016), 328–336.

11.

Moradi

M.H.

, Abedini

, Hosseinian

S.M.

, Optimal operation of autonomous microgrid using HS–GA, Int J Electr Power Energy Syst 77 (2016), 210–220.

12.

Parhoudeh

, Baziar

, Mazareie

, Kavousi-Fard

, A novel stochastic framework based on fuzzy cloud theory for modeling uncertainty in the micro-grids, Int J Electr Power Energy Syst 80 (2016), 73–80.

13.

Cheng

, Rodemann

, Fischer

, Olhofer

, Jin

, Evolutionary many-objective optimization of hybrid electric vehicle control: From general optimization to preference articulation, IEEE Trans Emerging Topics Comput Intell 1(2) (2017), 97–111.

14.

Zaman

, Elsayed

S.M.

, Ray

, Sarker

R.A.

, Evolutionary algorithms for dynamic economic dispatch problems, IEEE Trans Power Syst 31 (2016), 1486–1495.

15.

Mohamed

F.A.

, Koivo

H.N.

, Multiobjective optimization using mesh adaptive direct search for power dispatch problem of microgrid, Int J Electr Power Energy Syst 42 (2012), 728–735.

16.

, Ai

, Zhao

, A bi-level multi-objective optimal operation of grid-connected microgrids, Electric Power Syst Res 131 (2016), 60–70.

17.

Baziar

, Kavousi-Fard

, Considering uncertainty in the optimal energy management of renewable micro-grids including storage devices, Renew Energy 59 (2013), 158–166.

18.

Kavousi-Fard

, Abunasri

, Zare

, Hoseinzadeh

, Impact of plug-in hybrid electric vehicles charging demand on the optimal energy management of renewable micro-grids, Energy 78 (2014), 904–915.

19.

Aghajani

, Shayanfar

, Shayeghi

, Presenting a multi-objective generation scheduling model for pricing demand response rate in micro-grid energy management, Energy Convers Manage 106 (2015), 308–321.

20.

Mirjalili

, Mirjalili

S.M.

, Lewis

, Grey wolf optimizer, Adv Eng Software 69 (2014), 46–61.

21.

Rajasomashekar

, Aravindhababu

, Biogeography based optimization technique for best compromise solution of economic emission dispatch, Swarm Evol Comput 7 (2012), 47–57.

22.

Sharma

, Bhattacharjee

, Bhattacharya

, Grey wolf optimisation for optimal sizing of battery energy storage device to minimise operation cost of microgrid, IET Gen Transmission Distrib 10 (2016), 625–637.