Total Optimization of Energy Networks in a Smart City by Multi-Population Global-Best Modified Brain Storm Optimization with Migration

This paper proposes total optimization of energy networks in a smart city by multi-population global-best modified brain storm optimization (MP-GMBSO). Efficient utilization of energy is necessary for reduction of CO2 emission, and smart city demonstration projects have been conducted around the world in order to reduce total energies and the amount of CO2 emission. The problem can be formulated as a mixed integer nonlinear programming (MINLP) problem and various evolutionary computation techniques such as particle swarm optimization (PSO), differential evolution (DE), Differential Evolutionary Particle Swarm Optimization (DEEPSO), Brain Storm Optimization (BSO), Modified BSO (MBSO), Global-best BSO (BSO), and Global-best Modified Brain Storm Optimization (GMBSO) have been applied to the problem. However, there is still room for improving solution quality. Multi-population based evolutionary computation methods have been verified to improve solution quality and the approach has a possibility for improving solution quality. The proposed MS-GMBSO utilizes only migration for multi-population models instead of abest, which is the best individual among all of sub-populations so far, and both migration and abest. Various multi-population models, migration topologies, migration policies, and the number of sub-populations are also investigated. It is verified that the proposed MP-GMBSO based method with ring topology, the W-B policy, and 320 individuals is the most effective among all of multi-population parameters.


Introduction
In various environmental problems, global warming is a serious problem all over the world. It is considered that an increase of greenhouse gas causes global warming. Hence, many countries have conducted smart city demonstration projects for reduction of total consumption energies and CO 2 emission [1,2]. Smart city is an eco-city which can realize low carbon emission and energy consumption by Internet of Things (IoT), renewable energies, storage batteries, and so on. For example, in Japan, after the Great East Japan Earthquake, introduction of the smart city was investigated, especially in Tohoku area [3].
Basically, it is difficult to evaluate the actual CO 2 emission reduction in the actual smart city. Therefore, a smart city model including various sector models should be utilized for the evaluation. Static models which can treat various energy balances and dynamic models which can treat dynamic behaviors have been developed separately in each sector [4][5][6][7]. However, a smart city model, which can calculate various environmental and energy loads among all sectors considering environmental and energy flow among various sectors, had not been developed. Considering these backgrounds, a smart city model was developed in order to evaluate energy costs and CO 2 emission of the whole smart city considering environmental and energy flow among various sectors in Japan [8][9][10].
Using the smart city model, the authors have proposed total optimization of energy network in a whole smart city for minimization of energy costs, electric power loads at peak load hours, and CO 2 emission by Particle Swarm Optimization (PSO) [11], Differential Evolution (DE) [12], Differential Evolutionary PSO (DEEPSO) [13], BSO [14], MBSO [15], GBSO [16], and GMBSO [17]. In addition, considering energy facility characteristics, energy load characteristics, cost characteristics, and continuity of weekday operation of various energy facilities, the authors have proposed reduction methods of search space in order to solve the problem effectively [11,12]. Therefore, GMBSO considering the reduction of search space can obtain the highest quality solution so far [17]. However, there is still room for improvement of quality of solution.
Multi-population based evolutionary computation methods have been verified to improve quality of solution [18][19][20][21][22][23][24][25][26]. Using the methods, one population with numbers of individuals is divided into several sub-populations with small numbers of individuals. Individuals at each sub-population perform solution search separately. Interaction models of multi-population can be divided into three groups: the migration model [18,19,21,[23][24][25], the abest model [20,22], in which abest is the best searching point so far among all of sub-populations, and a model using both migration and abest. To the best of the authors' knowledge, since a model using both migration and abest has not been proposed, the authors propose the model [26]. However, the most suitable model has to be investigated for each evolutionary computation method.
Considering these backgrounds, this paper proposes total optimization of energy networks in a smart city by MP-GMBSO with only migration. MP-GMBSO with migration is newly proposed in order to realize improvement of quality of solution. The proposed MP-GMBSO with migration-based method is applied to a model of mid-sized city such as Toyama city. The results of the proposed method are compared with those of the GMBSO based method with a single population, and the MP-GMBSO based methods with only abest, and with both migration and abest model. Various numbers of sub-populations, various topologies (trigonal pyramid with four sub-populations, cube with eight sub-populations, hyper-cube with 16 sub-populations, and ring with 2, 4, 8, 16 sub-populations), various policies, various interaction models (only using migration, only using abest, and using both migration and abest), and various numbers of individuals have been investigated in simulation results.
The following gives a summary of the paper contributions: -A proposal of a new evolutionary computation method, namely MP-GMBSO with migration, in order to realize improvement of quality of solution, -An application of MP-GMBSO with migration-based method to total optimization of energy networks in a smart city, -Verification of efficacy of the conventional GMBSO based method for total optimization of smart city by comparing with the conventional DEEPSO, BSO, MBSO, and GBSO based methods, -Verification of efficacy of the proposed MP-GMBSO based method with migration for total optimization of smart city by comparing with the original GMBSO (GMBSO with one population) based method, and the MP-GMBSO based methods with various interaction model (only using migration, only using abest, and using both migration and abest), various policies, various topologies, various numbers of individuals, and various numbers of sub-populations, -It is verified that quality of solution is the most improved by the proposed MP-GMBSO with migration-based method using the ring topology with 16 sub-populations and 320 individuals, and the W-B policy (the worst individual of a sub-population is substituted by the best individual of other sub-populations) among all of multi-population parameters.
The rest of the paper is organized as follows. Section 2 explains a concept of the smart city model. Section 3 introduces problem formulation. Section 4 explains the details of the conventional

Concept of The Whole Smart City Model
The smart city model has been developed in order to quantitatively evaluate the smart city and the model can calculate energy costs and the amount of CO2 emission of whole smart city [8]. There are various sectors in the model such as electric power utility, natural gas utility, drinking water treatment plant and waste water treatment plant, industry, building, residence, and railroad (see Figure 1). There are energy supply-side and demand-side groups. The supply-side group includes natural gas utility, electric power utility, and drinking water treatment plant and waste water treatment plant sectors. The demand-side group includes industry, building, residence, and railroad sectors. The model quantitatively calculates energy cost, the amount of energy supply and demand, and CO2 emission in the whole smart city using energy flows among various sectors.

Sector Models in The Supply-Side Group
In the supply-side group, various energies such as electric power, natural gas, and drinking water are supplied interactively to the demand-side sector models [9]. The amount of supplied energies by supply-side sectors should be equal to the amount of energy loads of demand-side sectors.
The natural gas supply sector model can be modelled as a natural gas energy source. Namely, required amount of natural gas is supplied to various demand-side sectors by the model.
There are various electric power generation plants such as nuclear power, thermal power, hydroelectric power, and renewable energy generation in the electric power utility sector model. Costs and amount of CO2 emission of the power generation plants can be evaluated by the model. Power generation fuels can be purchased from outside of the city and fuel costs can be calculated. The amount of renewable energy and hydroelectric power generation output can be input to the model per an hour as fixed numerical values. Summation of the ratios of these power generation plants are set to 1.
In the drinking water treatment plant sector model, the amounts of drinking water loads, renewable energy generation, and required demand response are fixed. In the waste water treatment plant sector model, the amounts of inflow to waste water pipes, renewable energy generation, and required demand response are fixed.

Sector Models in the Supply-Side Group
The group interactively supply energies (natural gas, electric power, and drinking water). Namely, supplied energies by supply-side sectors and energy loads of demand-side sectors should be the same. The details are shown in [9,27].

Sector Models in the Demand-Side Group
Various energy supply facilities inside each sector are modeled. Various energy loads such as electric power, hot and cold heat, and steam loads are also dealt with. The energy supply facilities supply tertiary energies to the various energy loads when operations of facilities are decided. Hence, required secondary energies are calculated, and the amount is supplied from the energy supply group when various hourly energy loads of one day (24 points) are obtained and operations of facilities (decision variables) are determined. Figure 2 shows an industrial model as an example. The details are shown in [10,27]. boilers, and turbo and steam refrigerators are treated. The model also deals with various energy loads such as electric power load, heat load, and steam load (see Figure 2) [10]. When decision variables are determined, the energy supply facilities supply tertiary energies to various energy loads. Therefore, required secondary energy values are made clear and the values are supplied from the energy supply group when various hourly energy load values of one day (24 points) are given and decision variables are determined. Large buildings and large shopping malls are considered as the building sector model. The model deals with energy supply facilities such as gas turbine generators, turbo refrigerators, and steam refrigerator, and so on. The sector considers offices and shopping malls as energy loads and treats hourly loads of various energies and the amount of required demand response.
Condominiums and detached houses are treated as the residential sector model. The model treats both condominiums and detached houses with the same model using different parameters. The model deals with various facilities such as a fuel cell, a storage battery, a natural gas stove, a heat storage tank, a tank-less water heater, a heat pump water heater, and so on. These facilities supply electric power, hot water, and heat energies.

Problem Formulation of Total Optimization of Energy Networks in a Smart City
The smart city optimization considers the smart city which is a closed city such as an industrial park or a local government of a city, and minimizes energy cost, actual electric power at high load hours, and CO2 emission. Therefore, when the electric power companies or gas companies buy fuel from resource companies, the cost is not considered in the problem.

Decision Variables
Decision variables of various sectors are shown below: (a) Drinking water treatment plant model: The amounts of river inflow, water inflow to a service reservoir, electric power output by co-generators, and charged/discharged electric power of storage batteries. (b) Waste water treatment plant model: The amounts of waste water pumping to the waste water treatment plant, electric power output by co-generators, and charged/discharged electric power of storage batteries. (c) Industrial model: The amounts of electric power output by gas turbine generators, heat output by turbo refrigerators, heat output by steam refrigerators, and charged/discharged electric power of SBs. (d) Building model: The amounts of electric power output by gas turbine generators, heat output by TRs, and heat output by steam refrigerators. (e) Residential model: The amounts of electric power output by fuel cells, heat output by heat pump

Problem Formulation of Total Optimization of Energy Networks in a Smart City
This paper considers smart cities such as industrial parks and local governments of cities aiming at minimization of energy cost, actual electric power at high load hours, and CO 2 emission. Purchase costs of various fuels by electric power and gas companies from resource companies are not included in the problem.

Decision Variables
Decision variables in considered sectors are shown in Table 1. Each decision variable has to consider 24 hourly decision variables per day. Therefore, when the number of each facility at each sector is set to 1, the total number of decision variables can be 816. Hence, the problem can be considered as one of the large-scale optimization problems.

Objective Function
The objective function of the smart city problem considers three terms. The function considers minimization of energy costs of whole smart city as the first term. Actual electric power at high load hours are minimized by the second term of the function. The function also considers minimization of CO 2 emission in the whole city as the third term. The functions are shown in the following Equation (1): NumSec n=1 T j=1 pg n j × GU nj + pe n j × EU n j + w 2 where NumSec is the number of sectors except electric power and natural gas utilities, T is the number of hours per day (=24), pg nj is natural gas purchase at hour j of sector n, GU nj is a natural gas purchase unit cost at hour j of sector n, pe nj is electric power purchase at hour j of sector n, EU nj is an electric power purchase unit cost at hour j of of sector n, el nm is an actual electric power load at hour m of sector n, PS is the start hour of peak load hours of actual electric power loads, PE the end hour of peak load hours of actual electric power loads, GC is a coefficient between the purchased natural gas and CO 2 emission, EC is a coefficient between purchased electric power and CO 2 emission, and w 1 , w 2 , and w 3 are weighting coefficients (w 1 + w 2 + w 3 = 1).
Authors include one of a member of a smart city model committee in IEE of Japan. In the committee, a scope of a smart city has been discussed in order to consider a closed community such as industrial parks and local governments. In such a case, purchase costs of various fuels by electric power and gas companies from resource companies are not out of the scope of the city and we follow the scope in (1). After decision variables are calculated, the amounts of natural gas and electric power purchase are calculated as dependent variables. Hence, when the dependent variables violate limits of the constraints, penalty function values are calculated and they are added to the objective function values.

Constraints
(1) Energy balances: Electric power, hot and cold heat, and steam energy balances are considered.
These energy balances are expressed using the following equation: where g nr (y i , z i ) is a balance equation of energy r in sector n, y i is a startup or shutdown status of a facility for decision variable i (∈ {0, 1}), and z i is an input or output real value of a facility for decision variable i (∈ R), NumE n is the number of energies in sector n, and NumDim is the number of decision variables.
where h nq (y i , z i ) is an efficiency function of facility, or upper and lower bounds of facility q of sector n, and NumF n is the number of facilities in sector n. Efficiency of facility should be sometimes expressed with nonlinear functions. Hence, the problem is considered as one of mixed-integer nonlinear optimization problems (MINLPs) and evolutionary computation methods should be utilized in order to treat the problem.

The Proposed Multi-Population Global-Best Modified Brain Storm Optimization with Migration
MP-GMBSO is based on GMBSO which has been proposed using MBSO and GBSO [17]. Since MP-GMBSO is based on BSO variations, this section summarizes BSO, MBSO, GBSO and GMBSO.

Overview of BSO
BSO has developed by Shi in 2011 [28]. The BSO's main procedure is shown below.
Step 1 Initialization: Randomly generate NumInd individuals and calculate the objective function values of NI individuals.
Step 2 Clustering: The k-means method is applied to divide NumInd individuals into K clusters.
Step 3 Generation of New individual: One or two clusters are randomly selected, and new individuals are generated. Step 4 Selection of individuals: Individuals which are newly generated are compared with the current individuals which have the same individual indices of the newly generated individuals. Keep the better ones and the individuals are stored as the current individuals.
Step 5 Evaluation of individuals: The newly stored NumInd individuals are evaluated.
Step 6 The procedure can be stopped and go to Step 7 when the iteration number reaches the maximum iteration number which is pre-determined. Otherwise, go to Step 2 and repeat the procedures.
Step 7 The objective function value and the finally obtained variables are output as a set of the final solution.
Details of Step 2 and 3 are explained below.

Clustering of BSO
BSO utilizes the k-means method as a clustering method for dividing search space into small regions. A probability p clustering is utilized in order to change a cluster center to a randomly generated searching point. The method can avoid fast convergence, namely premature convergence, and individuals can escape from the local optima. The clustering algorithm is shown below: Step 2-1 Clustering: The k-means algorithm is applied to divide NumInd individuals into K clusters.
Step 2-3 Individuals are ranked in each cluster.
If r clustering ≥ p clustering (a pre-determined probability), the best individual is set as the cluster center in each cluster, Otherwise, One individual is randomly selected in each cluster, and the selected individual is set as the cluster center.

Generation of New Individual of BSO
The following equations are utilized in order to generate new individuals based on one or two current individuals: where, x new ij is decision variable j of the ith new individual, x old ij is decision variable j of the ith current individual, ξ(t) is a step size function, ITER is the number of maximum iteration, iter is the current number of iteration, c is a coefficient to change slope of the log-sigmoid transfer function, x old aj and x old bj are decision variable j of two individuals which are randomly selected two cluster centers or randomly selected two individuals of the selected two clusters. Several conditions utilized to determine x old ij is shown below: If p generation > rand(1, 0), one cluster is randomly selected.
If p OneCluster > rand(1, 0), the cluster center in the selected one cluster is set as x old ij and one new individual is generated using (4) and (5). Otherwise, one individual is randomly set as x old ij in the selected one cluster and one new individual is generated using (4) and (5).
If p generation ≤ rand(1, 0), two clusters are randomly selected. If p TwoCluster < rand(1, 0), Set the two cluster centers as x old aj and x old bj and the two centers are combined using (6). one new individual is generated using the combined x old ij with (4) and (5). Otherwise, Randomly select two individuals from the two selected clusters and the selected individuals are set as x old aj and x old bj and combined using (6). One new individual is generated using the combined x old ij with (4) and (5).

Overview of MBSO
MBSO has been proposed in order to improve two procedures of BSO [29]: Clustering and generation of new individuals. The k-means method is utilized for clustering in BSO algorithm and it takes a long time. Hence, simple grouping method (SGM) is applied in MBSO for reduction of calculation time. In addition, the method utilizes a different method for generation of new individuals for increase of individual diversification.

Clustering of MBSO
In MBSO, SGM is utilized as a clustering method. The algorithm of SGM is shown below: Step 1 K different individuals are randomly selected from the current generation as group centers of K groups.
Step 2 Calculate distances between the individuals and each group center. Distances to all group centers are compared. The individuals are assigned to the closest group.
Step 3 Individuals are ranked in each cluster.
If r clustering ≥ p clustering (a pre-determined probability), the best individual is set as the cluster center in each cluster, Otherwise, One individual is randomly selected in each cluster, and the selected individual is set as the cluster center.

Generation of New Individuals in MBSO
MBSO utilizes the following different equation in order to diversify individuals: where, pr is a probability utilized to determine which equation is utilized, L j is the lower bound of decision variable i, H j is the upper bound of decision variable i, and x s1 ij and x s2 ij are the randomly selected individuals for decision variable j of individual i.

Overview of GBSO
GBSO has been proposed in order to improve BSO [30]. The fitness-based grouping method is utilized as a clustering method in GBSO. In addition, the global-best idea information is added to individuals when a condition is satisfied. Details of the fitness-based grouping and the way to update individuals considering the global-best idea are explained in the next subsections. GBSO utilizes fitness-based grouping as a clustering method. An algorithm of the fitness-based grouping method is explained as follows: Step 1 Individuals are ranked using calculated values of the objective function.
Step 2 NumInd individuals are divided into K groups using (8).
where, g(i) is the selected group number of individual i, r(i) is a ranking of an individual i.

Generation of New Individuals of GBSO
GBSO utilizes the same BSO algorithm for generation of new individuals. However, information of "gbest", which is the best individual so far among all individuals, is applied to individuals when a certain condition is satisfied. In such a case, using (11), gbest information is added to x old ij of a selected one individual or x old ij which is combined with (6) using two selected individuals before new individual is generate in GBSO. Generally, when evolutionary computation is applied, diversification should be focused at the early searching iterations and intensification should be focused at the final searching iterations. Therefore, the authors change the condition (9) [16].
where, C is a probability and utilized to determine whether information of "gbest" is applied or not, Cmax is the maximum bound of C, Cmin is the minimum bound of C.
When the condition (9) is satisfied, information of "gbest" is applied to x old ij using the following equation:

Overview of GMBSO
The GMBSO has been proposed by the authors as a combined method of MBSO [29] and GBSO [30]. For clustering in GBSO, the Fitness-grouping method is utilized, which is proposed in [30]. For new individual generation, an update equation developed for MBSO is utilized which is proposed in [29] considering diversification. In addition, information of "gbest", is applied to individuals when a certain condition is satisfied. In such a case, using (11), information of "gbest" is applied to x old ij of a selected one individual or x old ij which is combined with (6) using two selected individuals before new individual is generate in GMBSO [17]. The procedure of GMBSO is shown below: Step 1 Initialization: NumInd individuals are randomly generated and evaluated.
Step 2 Clustering: NumInd individuals are divided into K clusters by Fitness-based grouping explained in 4.3.
Step 3 Generation of new individuals: Randomly select one cluster or two clusters. When the condition (9) is satisfied, information of "gbest" is applied to x old ij using (11). Then, new individuals are generated using Equation (7) explained in 4.2.
Step 4 Selection: The individuals which are newly generated are compared with the current individuals with the same individual indices. The better one is kept and stored as the current individual.
Step 5 Evaluation: The NumInd individuals are evaluated.
Step 6 The procedure can be stopped and go to Step 7 if the number of current iteration reaches the maximum number of iteration which is pre-determined. Otherwise, go to Step 2 and repeat the procedures.
Step 7 The objective function value and the finally obtained variables are output as a final solution.

Overview of the Proposed MP-GMBSO
Conventionally, the migration model [18,19,21,[23][24][25] (see Figure 3), the abest model [20,22] (see Figure 4), and a model using both migration and abest [26,31] (see Figure 5) has been considered in multi-population based evolutionary computation methods. The migration and abest model has been proposed by the authors because the smart city problem can be considered as one of the large-scale optimization problems [26]. Evolutionary computation methods with multiple searching points utilize only one population. On the contrary, the multi-population models divide searching individuals of one population into several sub-populations, and different searching process is performed at each sub-population. In addition, multi-population models exchange and replace searching individuals among sub-populations at certain intervals. The model shares the information of the best searching point of all sub-populations (abest) among all sub-populations in the abest model, and the combined migration and abest model. The various parameters of the models are shown below.

Update Equations of The Proposed MP-GMBSO with Migration
This paper proposes MP-GMBSO with migration. Equations for Fitness-based grouping, updating x old ij in order to apply information of "gbest" to x old ij , and equations for new individual generation are shown below.
Following Equation (12) is utilized for fitness-based grouping of MP-GMBSO.
where g s (i) is a selected group number of individual i in a sub-population s, R s (i) is the rank of individual i in a sub-population s, K s is the number of cluster in a sub-population s, and NumSubPop is the number of sub-populations.

Update Equations of The Proposed MP-GMBSO with Migration
This paper proposes MP-GMBSO with migration. Equations for Fitness-based grouping, updating in order to apply information of "gbest" to , and equations for new individual generation are shown below.
Following Equation (12) is utilized for fitness-based grouping of MP-GMBSO.
where is a selected group number of individual in a sub-population , is the rank of individual in a sub-population , is the number of cluster in a sub-population , and is the number of sub-populations.

Update Equations of The Proposed MP-GMBSO with Migration
This paper proposes MP-GMBSO with migration. Equations for Fitness-based grouping, updating in order to apply information of "gbest" to , and equations for new individual generation are shown below.
Following Equation (12) is utilized for fitness-based grouping of MP-GMBSO.
where is a selected group number of individual in a sub-population , is the rank of individual in a sub-population , is the number of cluster in a sub-population , and is the number of sub-populations.

Update Equations of The Proposed MP-GMBSO
This paper proposes MP-GMBSO. Equations for Fitness-based grouping, updating in order to add gbest information to , and new individual generation are shown below. Following Equation (12) is utilized for fitness-based grouping of MP-GMBSO.
where ( ) is group number of individual in a sub-population , ( ) is the rank of individual in a sub-population , is the number of cluster in a sub-population , and is the number of sub-populations.
An equation for updating in order to add gbest information to is expanded for multipopulation. There are two equations. One is for only migration model. Another one is for only abest model, and both migration and abest model.
where is the th individual of decision variable in sub-population , and is the best individual of decision variable in sub-population . -only abest model (see Figure 4), and a model using both migration and abest (see Figure 5): where is the best individual of decision variable among all individuals. (14) is utilized for comparison in Simulations.
An equation for new individual generation is expanded for multi-population as shown below.
where is the th newly generated individual of decision variable in sub-population , and and are the randomly selected two individuals for decision variable of individual in sub-population s.

Cutout Transformation Function
The variables are composed of discrete and continuous variables. This paper utilizes a cutout transformation function because decision variables have to deal with both startup/shutdown status and output values of various facilities in operation. The function can transform both discrete and continuous variables to only continuous variables. Using the function, the problem can be solved using only continuous variables. In addition, various evolutionary computation techniques can be applied easily by the function. Discrete variables ( = 0, … , ) express startup/shutdown status An equation for updating x old ij in order to apply information of "gbest" to x old ij is expanded for multi-population. There are two equations. One is for the only migration model. Another one is for the only abest model, and the combined migration and abest model. - the proposed MP-GMBSO with only migration model (see Figure 3): is decision variable j of the best individual in sub-population s. -only abest model (see Figure 4), and a model using both migration and abest (see Figure 5): where x abest j is decision variable j of the best individual among all individuals. (14) is utilized for comparison in Simulations.
An equation for new individual generation is expanded for multi-population as shown below.

Cutout Transformation Function
Both discrete and continuous variables are dealt as the decision variables of the problem. In order to easily utilize evolutionary computation method, decision variables should be only continuous variables. Therefore, the authors utilize a cutout transformation function which changes both discrete and continuous variables to only continuous variables. In addition, the function can reduce the number of decision variables. Startup and shutdown status of facilities can be expressed as discrete variables y i (i = 0, . . . , NumDim) (0 and 1). Operational values of facilities are expressed as continuous variables z i . Discrete variable y i is expressed as follows: Shutdown status (y i = 0) of facilities cannot be searched effectively considering (16), because y i is equal to 0 if and only if z i is equal to 0. Hence, new continuous variables x i (i = 1, . . . , NumDim) are utilized and the relationship between z i and x i are expressed as follows: where l i is a minimum bound of z i (the minimum input/output value of a facility), u i is an maximum bound of z i (the upper input/output bound of a facility), β i is a parameter of x i for shutdown status of a facility, α i is a lower limit of x i , and γ i is an upper limit of x i . Both discrete and continuous variables are utilized in the problem. However, the function can change both the discrete and continuous variables to only continuous variables and the evolutionary computation can be easily utilized.

Reduction of Search Space
Conventionally and usually, from the current operating conditions of various facilities, search space is calculated using the maximum decrease/increase of energy output/input values of facilities, upper/lower bounds of energy output/input values of facilities. The authors have proposed a reduction method of search space considering not only characteristics of various facilities but also considering electric and natural gas hourly cost and hourly load characteristics, and smooth transition to the next day facility operation in weekdays. As an example, the reduction of search space of the heat storage tank is explained in this section. The upper and lower bounds of the stored heat energies are calculated from the initial condition (HST initial ) [11,12]. The bounds are calculated by the equations shown below: ForMin i = HST i−1 − HinMin i + Hout i + Loss where HST i is the stored heat energies at hour i, ForMax i is the upper bound of stored heat energies in the heat storage tank at hour i, ForMin i is the lower bound of stored heat energies in the heat storage tank at hour i, HinMax i is the upper bound of heat input value at hour i, HinMin i is the lower bound of heat input value at hour i, Hout i is heat output at hour i, and Loss is hourly heat loss. Solid lines in Figure 9 show ForMax i and ForMin i and dash-dotted lines in Figure 9 show BackMax i and BackMin i . Using the same calculation concept of (20) and (21), BackMax i and BackMin i can be calculated using the opposite direction from T hour considering smooth transition to the next day facility operation in weekdays. Finally, considering the upper and lower bounds explained above, the following (22) and (23) shows the reduced search space at each hour.
where TMax i is the upper bound of the stored heat energies in the heat storage tank in the reduced search space at hour i, TMinx i is the lower bound of the stored heat energies in the heat storage tank in the reduced search space at hour i, HSTMax is an upper bound of heat energies in the heat storage tank, and HSTMin a lower bound of heat energies in the heat storage tank.
The reduced search space is shown with a shaded area in Figure 9 considering not only characteristics of various facilities but also smooth transition to the next day facility operation in weekdays and heat hourly load characteristics. Consequently, the reduction of search space can realize effective search to the target problem.
Algorithms 2018, 11, x FOR PEER REVIEW 13 of 25 energies in the heat storage tank in the reduced search space at hour , is an upper limit of heat energies in the heat storage tank.  Figure 9). However, the search space can be reduced to a shaded area in Figure 9 considering not only facility characteristics but also continuity of weekday facility operation and heat load characteristics. The large-scale optimization problem of the whole smart city can be effectively solved using the reduction of search space.

The Proposed Total Optimization Algorithm of Smart City by MP-GMBSO
The proposed total optimization algorithm by MP-GMBSO is shown below: Step 1 Initialization: All individuals are divided into sub-populations. Initial searching points of individuals are generated at each sub-population inside the reduced search space for a smart city.
Step 2 Energy cost, actual electric power load at peak load hour, and CO2 emission for the objective function are calculated at individuals in sub-populations. If the operational variables are not inside the reduced search space, penalty values are added to the objective function values.
Step 3 Clustering: Step 3-1 Generate clusters using in all sub-populations.
Step 3-2 Calculate objective functions of all individuals in each cluster in all subpopulations.
Step 3-4 The highest rank individuals at each cluster of all sub-populations are set as cluster centers. If > (1,0), randomly generate a new individual and replace a cluster center with the newly generated individual.
Step 4 New individual generation: select one or two individuals. When condition (9) is satisfied, gbest information is added to using (13) or (14). Generate new individuals considering several conditions explained in 4.1.2 using Equation (15).
Step 5 Selection: Step

The Proposed Total Optimization Algorithm of Smart City by MP-GMBSO with Migration
The proposed total optimization algorithm by MP-GMBSO with migration is shown below: Step 1 Initialization: Divide all individuals into NumSubPop sub-populations. Generate initial individuals at each sub-population considering the reduced search space for a smart city.
Step 2 Calculate the objective function at all individuals in sub-populations.
Step 3-2 Calculate objective functions of all individuals in each cluster in all sub-populations.
Step 3-4 The highest rank individuals at each cluster of all sub-populations are set as cluster centers. If p clustering > rand(1, 0), randomly generate a new individual and replace a cluster center with the newly generated individual.
Step 4 Generation of new individual: When condition (9) is satisfied, information of "gbest" is applied to x old ijs using (13) or (14). New individuals are generated considering several conditions explained in 4.1.2 using Equation (15).
Step 5 Selection: Step 5-1 The objective function values are calculated for all individuals.
Step 5-2 The new individual is compared with the current individual with the same individual index. Keep the better one and the individual is stored as the current individual in all sub-populations.
Step 6 Evaluation: Calculate objective function values of individuals in all sub-populations. The best individual is updated when the objective function value of the individual is better than the current best individual. Step 7 Individuals are migrated when the current iteration number reaches the migration interval which is pre-determined.
Step 8 the whole procedure is stopped and go to Step 9 when the current number of iterations reaches the maximum number of iterations which is pre-determined. Otherwise, go to Step 3 and repeat the procedures.
Step 9 The finally obtained objective function value and the best operational values are output as a final solution.

Simulation Conditions
The proposed MP-GMBSO with migration is applied to a typical mid-sized city model referring to data of Toyama city in Japan. The numbers of models in each sector in the model are shown below. They are determined so that sector load ratios of the model are approximately the same as the sector load ratios of Toyama city [32]: Industry model: 15, Building model: 50, Residential model: 45,000, Railroad: 1, Drinking water treatment plant: 1, Wastewater treatment plant: 1.
Usually, each smart city has each goal. Therefore, this paper considers the following three cases in order to consider different goals of smart cities: The maximum iteration number for DEEPSO based method is set to 1000 The conventional DEEPSO based method evaluates the objective function twice at each iteration for both new and clone individuals. Therefore, the maximum iteration number for DEEPSO is set to be half in order to keep the same number of objective function evaluation among all methods. Initial searching points are randomly generated for all methods. The simulation software has been developed using C language (gcc version 4.4.7) on a PC (Intel Xeon E5-2660 (2.20 GHz)).

Simulation Results
Firstly, the comparison of the mean, the minimum, the maximum, and the standard deviation of the objective function values for Case 1 among DEEPSO, BSO, MBSO, GBSO, and GMBSO based methods are shown in Table 2. Eighty individuals are utilized in Table 2. The number of individuals is set to 80 because the condition is the most severe for all methods. The results in the table are percentages of the objective function value when the mean of the objective function value by DEEPSO is set to 100%. It is verified that the mean, the minimum, the maximum, and the standard deviation of GMBSO based method can be the most reduced among all methods. Table 3 shows results of the mean ranks for Friedman test and p-value among the conventional DEEPSO, BSO, GBSO, MBSO, and GMBSO based methods. It is verified that there are significant differences at 0.05 significance level among all methods and the mean ranks by the GMBSO based method for each case are the best in Table 2.  Secondly, the best multi-population model is investigated among three models, namely, only migration model, only abest model, or a model using both migration and abest. Table 4 shows comparison of the mean, the minimum, the maximum, the standard deviation, and the mean rank of the objective function value among various numbers of sub-populations, and topologies using both migration with the W-B policy and abest, only migration with the W-B policy, and only abest with 640 individuals for Case 1 through 50 trials and a p-value by Friedman test. Six hundred and forty individuals are utilized in the table because there is a possibility to obtain better results by the MP-GMBSO based method through previous studies. The results of Table 4 are calculated when the mean of the objective function value using one sub-population is set to 100%. It is verified that the mean, the minimum, the maximum, and the standard deviation can be the most reduced by the only migration model with 16 sub-populations (bold numbers). It can be considered that the GMBSO based method can focus on intensification more than the other methods including the multi-swarm DEEPSO based method with a model using both migration and abest [26]. Therefore, the only migration model is the best model especially for the proposed MP-GMBSO based model in order to balance diversification and intensification. In addition, the mean rank by the only migration model with 16 sub-populations is the best among all parameters in the table. Table 5 shows the mean, the minimum, the maximum, the standard deviation, and the average rank of the objective function value of the optimal objective function values of each migration model with 640 individuals through 50 trials for Case 1 and a p-value by Friedman test. It is verified that the mean, the minimum, the maximum value, and the average rank using only migration model are the best among all multi-population models. It is also verified that there are significant differences at 0.05 significance level among all parameters in Table 4. Table 4. Comparison of the mean, the minimum, the maximum, the standard deviation, and average rank of the objective function value of the optimal objective function values among various numbers of sub-populations, and topologies using both migration with the W-B policy and abest, only migration with the W-B policy, and only abest with 640 individuals through 50 trials for Case 1 and a p-value by Friedman test. 1

Model
Mig. Policy # of sub-pop.

Ring
Trigonal  Next, the best migration policy and topology are investigated in Table 6 when the only migration model is utilized. Table 6 shows comparison of the mean, the minimum, the maximum, the standard deviation values, and the average rank of the optimal objective function values using only migration model among various migration policies, various topologies, and various numbers of sub-populations with 640 individuals through 50 trials for Case 1 and a p-value by Friedman test. The results of the table are calculated when the average of the objective function value using one sub-population is set to 100%. It is studied that, using only migration, the results can be the most reduced by the "W-B" policy with 16 sub-populations and the ring topology (bold numbers), and the average rank by the "W-B" policy with 16 sub-populations and the ring topology is the best among all parameters in the table.  Table 7 shows the best values of the mean, the minimum, and the maximum, the standard deviation values, and the average rank of the results calculated using optimal objective function values using only migration model of each migration policy with 640 individuals through 50 trials for Case 1 and a p-value by Friedman test. It is studied that the mean, the minimum, and the maximum values, and the average rank using W-B policy are the best among all migration policies (bold numbers). It is also studied that there are significant differences at 0.05 significance level among all parameters in Table 6. Table 7. The best values of the mean, the minimum, the maximum, the standard deviation values, and the average rank of the optimal objective function values using only migration model of each migration policy with 640 individuals through 50 trials for Case 1 from Table 6  The best topology, the best number of individuals, and the best number of sub-populations are investigated in Table 8. Table 8 explains the comparison of the mean, the minimum, the maximum the standard deviation values, and average rank of the results calculated using objective function values using only migration model with the W-B policy among various topologies, various numbers of individuals, and various number of sub-populations and a p-value by Friedman test. When the number of individuals is 40 with 4 or 8 sub-population, the results include penalty value because dependent variables are out of allowable ranges. On the contrary, when the number of individuals is more than 80 or equal to 80, the results do not include penalty because dependent variables are inside allowable ranges. Especially, when the number of individuals is set to 320 with 16 sub-populations using ring topology, the mean, the minimum, and the maximum values can be the most reduced (bold numbers). Namely, the proposed method with such conditions can balance diversification and intensification the most effectively for the problem. In addition, the average rank is the best when the number of individuals is set to 320 with 16 sub-populations using ring topology among all parameters in Table 8.  Table 9 shows the best values of the mean, the minimum, the maximum, the standard deviation, and average rank of the objective function value of the optimal objective function values of each number of individuals with W-B policy using only migration and ring topology through 50 trials for Case 1 and a p-value by Friedman test. It is verified that the mean, the minimum, the maximum value, and the average rank with 320 individuals and 16 sub-populations are the best among all multi-population models (bold numbers). It is also verified that there are significant differences at 0.05 significance level among all methods.
Tables 10 and 11 investigate the above results (Tables 2-9) from a facility operation point of view. In order to clarify these results, Table 10 shows comparison of the optimal operation of the facilities in an industrial model of Case 2 by MP-GMBSO using migration model, and W-B migration policy with ring topology when the number of individuals is set to 320 among various numbers of sub-populations. Column A indicates electric power output of a GTG and column B indicates electric power which is purchased from electric power utility at each hour. From 8 to 22 h, in the model, electric power which is purchased from electric power utility is expensive and electric power which is output of a GTG is inexpensive and purchased electricity is expensive. Hence, electric power which is purchased from electric power utility should be reduced and electric power which is output of a GTG should be increased from 8 to 22 h for reduction of energy cost. It is investigated that the electric power which is purchased from electric power utility can be the most reduced and electric power which is output of a GTG can be the most increased from 8 to 22 h by the proposed method (MP-GMBSO using migration model) with 16 sub-populations (bold numbers). Table 11 shows comparison of the best facility operation of an industrial model for Case 3 using only migration model, and W-B migration policy with ring topology when the number of individuals is set to 320 among various numbers of sub-populations. Electric power which is purchased from electric power utility should be reduced and electric power which is output of a GTG should be increased a whole day for reduction of CO 2 emission in the model. It was verified that electric power which is purchased from electric power utility can be reduced and electric power output of a GTG can be increased at most a whole day by the proposed method with 16 (bold numbers).
Consequently, the proposed MP-GMBSO with only migration-based method with the W-B policy, the ring topology, 320 individuals, and 16 sub-populations is the most effective to the total optimization of energy networks in the smart city problem. Table 9. The best values of the mean, the minimum, the maximum, the standard deviation, and average rank of the objective function value of the optimal objective function values of each number of individuals with WB policy using only migration and ring topology through 50 trials for Case 1 from

Conclusions
This paper proposes total optimization of energy networks in a smart city by multi-population global-best modified brain storm optimization (MP-GMBSO) with migration. A new evolutionary computation technique, called MP-GMBSO is proposed in order to obtain better results and the MP-GMBSO based method is proposed for total optimization of energy networks in a smart city. Effectiveness of the conventional GMBSO based method for total optimization of smart city by comparing with the conventional DEEPSO, BSO, MBSO, and GBSO based methods is verified. Moreover, effectiveness of the proposed MP-GMBSO based method for total optimization of smart city with migration is verified by comparing with the original GMBSO (MP-GMBSO with a sub-population) based method, and the MP-GMBSO based methods with various migration model (only using migration, only using abest, and using both migration and abest), various policies, various topologies, various numbers of individuals, and various numbers of sub-populations. It is verified that the quality of solution is the most improved by the proposed MP-GMBSO with migration-based method using the ring topology with 16 sub-populations and 320 individuals, and the W-B policy (the worst individual of a sub-population is replaced by the other sub-populations) among all multi-population parameters. Solutions by the proposed method can be changed by tuning weighting coefficients of the objective function considering the purpose of the target smart city.
As future works, several challenges can be considered. The smart city problem can be considered as one of the large-scale MINLP problems and it is hard to obtain solutions effectively. Therefore, firstly, state-of-the-art techniques of evolutionary computation will be applied in order to improve quality of solution. Secondly, various uncertainties of the problem will be considered. For example, uncertainty of renewable energy output depends on weather conditions. Uncertainty of various energy loads depends on various sector conditions. Hence, stochastic approaches such as a Monte Carlo simulation