Optimal Siting and Sizing of Distributed Generation Based on Improved Nondominated Sorting Genetic Algorithm II

With the development of distributed generation technology, the problem of distributed generation (DG) planning become one of the important subjects. This paper proposes an Improved non-dominated sorting genetic algorithm-II (INSGA-II) for solving the optimal siting and sizing of DG units. Firstly, the multi-objective optimization model is established by considering the energy-saving benefit, line loss, and voltage deviation values. In addition, relay protection constraints are introduced on the basis of node voltage, branch current, and capacity constraints. Secondly, the violation constrained index and improved mutation operator are proposed to increase the population diversity of non-dominated sorting genetic algorithm-II (NSGA-II), and the uniformity of the solution set of the potential crowding distance improvement algorithm is introduced. In order to verify the performance of the proposed INSGA-II algorithm, NSGA-II and multiple objective particle swarm optimization algorithms are used to perform various examples in IEEE 33-, 69-, and 118-bus systems. The convergence metric and spacing metric are used as the performance evaluation criteria. Finally, static and dynamic distribution network planning with the integrated DG are performed separately. The results of the various experiments show the proposed algorithm is effective for the siting and sizing of DG units in a distribution network. Most importantly, it also can achieve desirable economic efficiency and safer voltage level.


Introduction
In recent years, distributed generation systems using renewable energy technologies such as wind power generation and photovoltaic power generation have become one of the hotspots of research at home and abroad. After the distributed generation (DG) units are connected to the distribution network (DN), the structure, operation mode, and control strategy of the DN will undergo tremendous changes [1,2]. The research shows that DGs provide more flexibility and expansibility for distribution network. When DGs are connected to the distribution network (DN), the performances of DGs are most important for the power quality, reliability, and security of DN. In addition, DGs with renewable energy technology can effectively reduce system line loss and transmission congestion [3][4][5]. Despite the above advantages, if the placement and sizing of the DGs are improperly selected, it may cause a series of power system safety hazards such as power flow reverse, insufficient voltage stability, and malfunction of the protection device [6,7]. Therefore, the placement and sizing of the DGs has become one of the important subject. subjectivity of weight selection makes the calculation result unsatisfactory. Aiming at the above problems, multi-objective optimization algorithm is introduced to solve the DG planning problem [21]. J. Neale solved the reconfiguration problem of the DN with the integrated DG by using non-dominated sorting genetic algorithm-II (NSGA-II) [22]. This method can effectively improve the operation performance of DN. A. Ameli used the multi-objective particle swarm optimization (MOPSO) [23] with adaptive grid method to optimize the distribution of energy supply system electric vehicle (EV) in DG technology. The results show the effectiveness of MOPSO programming; R. S. Maciel proposed an evolutionary particle swarm optimization algorithm (MEPSO) [24] as a multi-objective optimization algorithm for DG, of which the planning scheme can ultimately improve the objectives. However, in the multi-objective optimization algorithm, there is a common problem; that is, the evaluation information of solution set is insufficient due to crowding distance, which leads to the elimination of potential high-quality solutions in the truncation process. The distribution of the final algorithm planning results is uneven.
In addition, the mutation operator of NSGA-II will gradually lose the ability of local search as the dimension of solution variable increases, so the population updating strategy in the fireworks algorithm [25] is introduced to simultaneously improve the mutation operator and the diversity of solution set. Based on the effective Pareto-optimal set of multi-objective optimization problems, decision makers cannot get a unique solution. Therefore, this paper introduces an unbiased compromise strategy [26] to choose a best compromise from the Pareto-optimal set.
According to the above analysis, a lack of the operation state research planning of DN with the integrated DG considers the relay protection. Actually, when the DG is connected to the DN, the DN that is originally powered by a single system power becomes a network with multiple power structures, causing a change in the magnitude and direction of the short-circuit current during the failure. The action affects the safe and stable operation of the DN. Therefore, this paper increases the short-circuit current constraint condition and compares the result with the optimization result without considering the short-circuit constraint, so as to prove the necessity of considering the relay protection condition. At present, with the development of technology, distributed energy resources are bound to penetrate into the fields of industry, commerce, and urban and rural residents. The problem of DG planning with timing characteristics needs to be solved urgently. Therefore, the paper studies the planning of various examples and plans the static state and dynamic running state of the IEEE 33-bus system separately to provide a reasonable DG planning for decision makers.
The main contributions of this work are as follows: (1) An improved non-dominated sorting genetic algorithm (INSGA-II) for placement and sizing of DGs is proposed. The violation constrained index and improved mutation operator are proposed to increase the population diversity of NSGA-II, and the uniformity of the solution set of the potential crowding distance improvement algorithm under the 3D objective function is introduced. (2) The energy-saving benefit of DGs is considered, and a multi-objective optimization model is established. (3) The restraint condition of relay protection is added; that is, the size of short-circuit current is limited, so that the protection device will not malfunction. (4) Convergence metric and spacing metric are employed to evaluate the performance of INSGA-II, NSGA-II, and MOPSO. (5) The effectiveness of the algorithm in different cases is verified.
The proposed algorithm is applied to solve the static and dynamic planning of the DN with the integrated DG units. The rest of this article is organized as follows. Section 2 presents the mathematical programming model. Section 3 briefly introduces the algorithms needed in the planning scheme and proposes INSGA-II. Section 4 provides the comparison of the simulation results of the proposed algorithm in various cases and the performance of the algorithm. Section 5 concludes our work.

Objective Functions
Three objectives should be taken into account when establishing multi-criteria optimization model of DN with the integrated DG units: (1) Improving the energy-saving benefit of the system; (2) reducing system line losses; and (3) reducing the node voltage deviation.

Maximizing annual energy-saving benefit
Maximizing annual energy-saving benefit of DN with the integrated DG units can be expressed as follows: where Z NODG cost is the total annual cost of DN without DGs and Z NODG cost is the total annual cost of DN with DGs.
Total annual costs of DN excluding DGs can be expressed as follows: where C loss is the annual loss of DN without DGs and C b is the annual purchase cost of DN without DGs.
where k is the number of branches in the DN; C p is the unit price of electricity consumed per unit ($/kWh); τ max is the annual maximum load loss hour (h) of the DN; and ∆P j is the active power loss (kW) of the jth branch.
where P load is the DN total active load and T max is the DN annual maximum load utilization hours. The total annual cost of DN with DGs can be expressed as follows: where C dgm is the annual investment and maintenance cost of DGs; C ploss is the annual cost of DN with DGs (calculated in the same way as C loss ); C B is the annual purchase cost of DN with DGs; and C sub is the financial subsidy of new energy generation.
where r is the annual interest rate; t is the planning period; C dgi is the ith distributed generation equipment investment cost ($/kWh); C mi is the annual operation and maintenance cost ($/kWh) for the ith distributed generation; and P dgi is the installed capacity (kW) of the ith distributed generation.
In order to reflect the benefits of DGs in environmental protection, the government has put forward policy support to DGs, that is, financial subsidies for distributed generation. It can be expressed as follows: where C sp is the subsidy amount ($) of distributed generation unit.

Minimize line losses
The existence of line losses will lead to line heating, which will accelerate the aging of insulated lines, will reduce the insulation of lines, and will ultimately lead to the risk of leakage. Reducing line losses can improve energy efficiency of electrical equipment or processes.
Lowering the line losses can improve the energy utilization efficiency of the energy-using equipment or process, and it is also one of the important measures of energy saving. The objective function can be expressed as follows [27]: where n is the total number of DN nodes; g ij is the admittance of the branch (i, j); U i and U j are the voltage amplitudes of branches i and j respectively; and θ ij is the voltage phase angle.

Minimize the total voltage deviation
With the increase of load, the system voltage stability will deteriorate gradually, even voltage collapse, resulting in a system in danger. Therefore, the voltage deviation is one of the important indexes to evaluate the operation safety and power quality of the system. In addition, the increase of network node voltage can effectively reduce the reactive power loss of the system. The objective function can be expressed as follows: where U i is the node ith voltage real value of the DN when DG is incorporated and U speci f ied i is the voltage rated value. This paper assumes that the rated voltage is 1 (p.u).

Influence of DG on Relay Protection of DN
DN is generally equipped with three-stage current protection. After the DG is connected to the DN, the system changes from single power supply to multiple power supply. When the transmission line is short-circuited, the magnitude and direction of the short-circuit current will change. Figure 1 shows a typical simple DN. L 1 and L 2 represent feeders; T 1 represents a transformer; PD 1 , PD 2 , PD 3 , and PD 4 are protective devices of the corresponding feeders; and DG is connected to the DN from node C.
where sp C is the subsidy amount ($) of distributed generation unit.

Minimize line losses
The existence of line losses will lead to line heating, which will accelerate the aging of insulated lines, will reduce the insulation of lines, and will ultimately lead to the risk of leakage. Reducing line losses can improve energy efficiency of electrical equipment or processes.
Lowering the line losses can improve the energy utilization efficiency of the energy-using equipment or process, and it is also one of the important measures of energy saving. The objective function can be expressed as follows [27]: where n is the total number of DN nodes; ij g is the admittance of the branch ( i , j ); i U and j U are the voltage amplitudes of branches i and j respectively; and ij θ is the voltage phase angle.

Minimize the total voltage deviation
With the increase of load, the system voltage stability will deteriorate gradually, even voltage collapse, resulting in a system in danger. Therefore, the voltage deviation is one of the important indexes to evaluate the operation safety and power quality of the system. In addition, the increase of network node voltage can effectively reduce the reactive power loss of the system. The objective function can be expressed as follows: where i U is the node ith voltage real value of the DN when DG is incorporated and specified i U is the voltage rated value. This paper assumes that the rated voltage is 1 (p.u).

Influence of DG on Relay Protection of DN
DN is generally equipped with three-stage current protection. After the DG is connected to the DN, the system changes from single power supply to multiple power supply. When the transmission line is short-circuited, the magnitude and direction of the short-circuit current will change. Figure 1 shows a typical simple DN. L1 and L2 represent feeders; T1 represents a transformer; PD1, PD2, PD3, and PD4 are protective devices of the corresponding feeders; and DG is connected to the DN from node C.  Assuming that the DG is connected, the short-circuit fault can be divided into the following three cases: 1. The adjacent feeder connected to the DG is short-circuited When the short-circuit fault F 1 occurs at the end of the adjacent feeder, the power source S and DG provide the overlapping short-circuit current to PD 1 and PD 2 . Compared without DG connection, the short-circuit current increases and its value may be larger than the original I segment setting value of PD 1 , which makes PD 1 malfunction, resulting in the loss of coordination between PD 1 and PD 2 . Simultaneously, PD 3 will also flow through the reverse current provided by DG. When the capacity of DG is larger, the value of reverse current may be larger than the original I segment setting value of PD 3 , resulting in PD 3 malfunction.
2. The upstream of DG access is short-circuited When the DG upstream feeder AC is short-circuit fault F 2 , PD 4 can operate normally. However, when PD 3 is activated, the island downstream of the DG will form an island operation. At the same time, the DG will always provide short-circuit current to the fault point, which will affect the reliability of the system.

The downstream of DG access is short-circuited
When DG downstream feeder CD short-circuit fault F 3 occurs, protections PD 1 and PD 2 can operate normally. Although the short-circuit current through PD 3 is only supplied by the system power S, the fault current will be smaller than before due to DG connection, which may cause the backup protection of PD 3 to refuse to operate. The short-circuit current through PD 4 is provided by the system power S and DG together. The increase of short-circuit current may increase the protection distance of PD 4 , thus losing the selectivity.
When there are multiple DGs in the system, the analysis method is similar. From the above analysis, it can be seen that the access of DG will have adverse effects on the reliability of the relay protection and system of DN, so it is necessary to consider the restraint of relay protection when optimizing the configuration of DG. In the traditional three-stage current protection, the I-stage instantaneous current quick-break protection is the most important, so this paper focuses on the analysis of the impact of DG on the I-stage of the original current protection of the PD in the DN.

Constraint Condition
DG connected to a distribution network has great influence on power flow, voltage distribution, branch current, and power quality of the system. In order to design a reasonable programming model, the following equality constraints and inequality constraints are included in this paper.

Equality constraint
Below are the power flow equations constraints.
where P DGi and Q DGi are the ith node active and reactive power of the DGs injected; P Li and Q Li are the ith node active and reactive powers of load output; U k is the voltage amplitude of all kth nodes connected to ith nodes; G ik is branch conductance; B ik is the branch susceptance; θ ik is the difference between the voltage angles of ith node and kth node; and N 1 is the number of branches associated with the ith node.

Inequality constraints
In order to make the planning scheme meet the power system operation standards, it is usually required that the node voltage, current, power, and permeability satisfy the constraints after DGs are connected to the grid. Node voltage constraints where U i represents the voltage amplitude of the ith node, and U min i and U max i represent the upper and lower limits of the ith nodes voltage amplitude, respectively. b.
Branch current constraints where I ij is the current of branch (i, j) and I max ij is the maximum current allowed of branch (i, j).
c. Single DG access capacity constraints where P DGi and Q DGi are the active power and reactive power of the DG connected to the ith node, respectively. d.
Distributed generation access total capacity constraints where S DGi is the total capacity of distributed generation and S L is the maximum load value of DN. e.
Short-circuit current constraints where K III sen,i is the sensitivity coefficient of current III segment protection of branch i as backup protection; I I set,i , I II set,i , and I III set,i are the setting values of the current I, II, and III protections of branch i, respectively; I I set,j is the setting value of the current I segment protection of the branch j of the feeder line of the DG; I k,j.max are the maximum short-circuit currents through branch i of the terminal three-phase short circuit of branch i and branch j, respectively; and I (2) k,i+1.max is the minimum short-circuit current through branch i of the two-phase short circuit at the end of the downstream of branch i.

Overview formulation
From the above analysis, the multi-objective function and constraints can be described as follows: where f i is the ith objective function, N obj is the number of objective functions, x s and x c denote the state vector and the control vector respectively, h i is the equation constraint, p is the number of equality constraints, g i is the inequality constraint, and q is the number of inequality constraints.
x c is composed of n variables, and n represents the number of nodes of the optimized network, where each variable represents two components: the DG installation placement and capacity. For example, a variable is assigned a value indicating the placement of the variable to installation DG. The value of the variable indicates the installation capacity. If DG is not installed, the corresponding variable value is 0. x c can be illustrated as follows: where L DGi indicates that the ith node is the DG installation location and P DGi represents the DG installation capacity of the ith node.
x s is the network parameter that x c is calculated from the power flow. Variables such as line losses, node voltage, and branch current are used to calculate multi-objective functions and to determine the satisfaction of constraints. x s can be illustrated as follows: where P , Q , U , θ , and I represent the network feeder active power vector, reactive power vector, node voltage vector, voltage argument vector, and branch current vector respectively after DG is connected.

Establish PQ Mode of DG
The DG injects its power output into the grid through power electronics. Typically, for a PQ model, which shows the active power (P) versus reactive power (Q) called a PQ model. DG is modeled as a constant power factor and negative load model.
In this case, the DG is modeled as a constant power source. P DG is the actual power output specified for the DG model. The load at ith node with the DG device is modified as Equations (22)- (23).
where P loadi is the original network load power of ith node, P DGi and Q DGi are the active and reactive power of DG at the ith node, P loadi is the active power after installing DG, and cos ϕ is the power factor. In general, constrained problems can be solved using either deterministic or stochastic algorithms. However, deterministic approaches such as feasible direction and generalized gradient descent require strong mathematical properties of the objective function such as continuity and differentiability. In cases where these properties are absent, evolutionary computation, such as NSGA-II, offers reliable alternative methods [13].

NSGA-II Algorithm Overview
NSGA-II was proposed by Deb, K. in 2002 [28], which is one of the most popular multi-objective optimization algorithms. It uses simulated binary crossover and polynomial variation to introduce non-dominated sorting and crowding distance operator instead of NSGA. The sharing radius of the algorithm ensures the diversity of the population and the uniformity of the Pareto-optimal set and introduces the elite retention strategy and the elimination strategy, which improves the operation speed and stability of the algorithm.

Dominant, Non-Dominated, Pareto-Optimal Set and Crowding Distance
NSGA-II introduces Pareto-optimal set and crowding distance to replace the fitness of traditional intelligent optimization algorithm. The solution set is divided into dominated solution set and non-dominated solution set according to the dominant relationship among the solutions. The rank of the solution is from 1 to n: the better the quality, the smaller the value. Solutions of the same rank are divided into pros and cons by comparing the size of the crowding distance.
Multi-objective optimization problems can be described as follow: where I is a feasible solution space.

Definition 1.
If solution x 1 dominates x 2 , then the following definitions must be met.

Definition 2.
If solution x 1 does not dominate x 2 , then the following definitions must be met.
Definition 3. Pareto-optimal set is composed of mutually non-dominated solution sets, and ranks are composed of mutually dominated solution sets.
Definition 4. Assume that P = {x 1 , x 2 , · · · , x n } is a Pareto-optimal set and that the crowding distance represents the distribution density of other solution sets around the solution. The larger the solution distance at the same layer, the better the solution distribution, that is, the better the solution set diversity. Its calculation formula is as follows: where cd k i is the crowding distance on the kth objective function of x i , m is the number of objective functions, index (x k i ) is the sort index of x i on the kth objective function, f k i+1 is the value of objective function corresponding to the last solution on the axis of the kth objective function, and f k max and f k min represent the maximum and minimum values of the kth objective function, respectively.

INSGA-II: Improved Mutation Operator
Because mutations in genetics often lead to worse results, the probability of mutation operations in genetic algorithm (GA) is usually set to be very small. In GA, crossover and mutation operations have two functions: global search and local search. When dealing with low-dimensional convex problems, only the crossover operator can solve these problems well. However, when dealing with high-dimensional nonconvex problems, it often falls into the local optimal solution. For this reason, the mutation operator is improved by referring to the population-updating method in the fireworks algorithm [25]; the new mutation individual is created by Equations (29)- (31).
where d i denotes the dimension of the ith solution for mutation, vd denotes the variable dimension, rand denotes the random number between [0, 1], S denotes the index for performing the mutation operation, u mu j denotes the jth gene performing the mutation operation, and u j represents the value of the jth gene after the mutation operation.
Crossover and mutation operations often produce solutions beyond the feasible space. The general method is to assign solutions less than the lower limit of feasible region to the lower limit, and the solution that exceeds the upper limit of the feasible region is given the upper limit value. This operation is simple to perform but reduces the diversity of solution sets. According to Equation (32), solutions beyond the boundary can be effectively mapped back to the feasible space.
where u max i and u min i represent the upper limit and lower limit of the gene, respectively, and % is the remainder operation.

INSGA-II: Violation Constrained Index
In view of the subjectivity of penalty factor selection, the violation constrained index (VCI) is proposed when dealing with infeasible solutions.
The method divides the solution set into feasible solution and infeasible solution by calculating the error degree (ED) of the constraint condition of the solution.
The VCI calculation formula is as follows: where ED max j and ED min j represent the maximum and minimum error degree of infeasibility under the jth constraint, respectively, and C indicates the number of constraints. The classification process of the solution set is performed according to Equation (35).
where S_t denotes the type of the solution set. The default index can be effectively used in combination with the non-dominated sorting strategy. For example, when the S_t is 0 or 1/2, the feasible solution part is given the true rank and the crowding distance, and the rank of the infeasible solution is sorted in descending order according to VCI i and continues to be sorted from the feasible solution. Figure 2 shows that, after solutions a-h pass the truncation strategy, a, d, e, g and h are preserved but the eliminated solution c has a better crowding distance than the solution d, which is due to the traditional crowding distance information being missing, and the truncation strategy eliminates the solution set with small crowding distances at one time, resulting in uneven distribution of the final set. The calculation of the traditional crowding distance only considers the spatial distribution of two adjacent solutions on the axis of the objective function without considering that the effect of deleting the neighborhood solution leads to the elimination of potential high-quality solutions. Finally, a potential crowding distance is proposed. The potential crowding distance in the two-dimensional space is simple to calculate. As shown in Figure 2a, when the solution map is re-target functions 1 and 2, the neighborhoods are the same in the solution set and only the order changes. The truncation process is shown in Figure 2b,c:

INSGA-II: Improved Crowding Distance Operator
where _ S t denotes the type of the solution set. The default index can be effectively used in combination with the non-dominated sorting strategy. For example, when the _ S t is 0 or 1/2, the feasible solution part is given the true rank and the crowding distance, and the rank of the infeasible solution is sorted in descending order according to i VCI and continues to be sorted from the feasible solution. Figure 2 shows that, after solutions a-h pass the truncation strategy, a, d, e, g, and h are preserved but the eliminated solution c has a better crowding distance than the solution d, which is due to the traditional crowding distance information being missing, and the truncation strategy eliminates the solution set with small crowding distances at one time, resulting in uneven distribution of the final set. The calculation of the traditional crowding distance only considers the spatial distribution of two adjacent solutions on the axis of the objective function without considering that the effect of deleting the neighborhood solution leads to the elimination of potential high-quality solutions. Finally, a potential crowding distance is proposed. The potential crowding distance in the two-dimensional space is simple to calculate. As shown in Figure 2a, when the solution map is retarget functions 1 and 2, the neighborhoods are the same in the solution set and only the order changes. The truncation process is shown in Figure 2b,c:  Therefore, the formula for calculating the potential crowding distance in two-dimensional space is as follows:

INSGA-II: Improved Crowding Distance Operator
( )   Therefore, the formula for calculating the potential crowding distance in two-dimensional space is as follows: where ∆d 1 i and ∆d 2 i denote the increments of the crowding distance after the adjacent solution of the kth solution is deleted, f 1 i−1 and f 2 i−1 denote the neighborhood solutions mapped to the kth solution on the objective function 1, and pd i denotes the potential crowding distance of the kth solution.
The objective function above three-dimensions does not have the symmetric relationship of the two-dimensional objective function. Each solution no longer has only the same two adjacent solution sets on each objective function axis but has (2-n) different solutions. Therefore, the potential crowding distance of the solution under the 3D objective function is calculated as follows: where Φ denotes the set of adjacent fields mapped to the ith solution on each objective function axis; m denotes the number of objective functions; k∆d i denotes the increment generated by the crowding distance of the ith solution after the kth solution in the neighborhood set is eliminated; and 1 j , 1 l , and 1 r respectively indicate that the kth solution belongs to the neighborhood of the ith solution or belongs to the left neighborhood solution or the right neighborhood solution on the jth objective function.

INSGA-II: Improved Crowding Distance
NSGA-II uses the tournament selection operator for population reproduction, and the selection criteria considers the rank and crowding distance. The specific selection process is described in Reference [27]. Since the improved NSGA-II algorithm increases the potential crowding distance, it needs to involve the crowding distance for the selection operation. After some modifications, the selection rules are as follows: x i , i f cd i > cd j and pd i > pd j x j , i f cd j > cd i and pd j > pd i x i or x j , otherwise where newx k denotes the selected kth descendant solution, and x i and x j respectively denote two different solutions in the parent. Improved selection strategy comprehensively considers the current crowding distance and potential crowding distance. Therefore, it can select a better solution and can avoid the defects of insufficient diversity of the solution set when the truncation strategy is executed.

Unbiased Compromise Strategy
In order to solve the cumbersome scheme of the Pareto-optimal set, this paper uses an unbiased compromise strategy based on fuzzy set, which uses unbiased member parameter ω to evaluate the quality of each individual in the Pareto-optimal set. The calculation formula is as follows: where ω i j denotes the unbiased parameter of the ith solution in the Pareto-optimal set on the jth objective function; f min j and f max j denote the minimum and maximum values of the jth objective function corresponding to the Pareto-optimal set, respectively; f i j denotes the value of the jth objective function corresponding to the ith solution in the Pareto-optimal set; n denotes the number of objective functions; m denotes the number of solutions in the Pareto-optimal set; and ω * denotes the unbiased member parameters of the comprehensive optimal solution.

Optimized Configuration Process Based on Improved NSGA-II
The proposed algorithm optimizes the placement and sizing process configuration of DG, as follows:

Algorithm 1. The algorithm of the proposed INSGA-II
Input: Set IEEE-33 distribution network-related parameters, the number of target functions, the number of variables, the size of the population, the maximum number of iterations, etc. Outputs: Comprehensive optimal solution x * . 1: Initialization population P(x c ) t , set the number of iterations t = 0; 2: P(x s ) t is calculated by forward and backward substitution method, and then, the objective function values f 1 , f 2 , f 3 , and S_t are calculated; 3: Algorithm iteration start, while t < t max do; 4: The solution set is classified by S_t, feasible solution to perform non-dominated sorting strategy, and the infeasible solution to calculate DCI and sorted; 5: Generation of children Q(x c ) t from the parent P(x c ) t by improving selection operations and improving genetic manipulation; 6: Q(x c ) t performs power flow calculation, mixes with P(x c ) t to form a mixed population R(x c ) t and R(x c ) t for non-dominated sorting strategy, and calculates the crowding distance according to Equations (37) and (38); 7: Selecting a new parent population R(x c ) t of size from P(x c ) t+1 by a truncation strategy; 8: t = t + 1; 9: End while; 10: Using unbiased compromise strategy to select comprehensive optimal solution x * ; 11: Return x * .

Experiment Setup and Description
In order to better verify the effectiveness of the proposed algorithm, IEEE 33

IEEE 33-Bus System Case Simulation Experiment
An IEEE 33-bus system topology diagram is shown in Figure 3, where the system includes 33 nodes and 32 branches [30]. The 0 node is assumed to be a balanced node with a voltage reference of 12.66 kV, total power load of 3.17 MW, and reactive power 2.30 MVAr.

Analysis of Optimization Results with Different Parameters
As shown in Table 1, in order to improve the operation speed of the algorithm and the accuracy of the solution, experimental comparisons of the experimental hyperparameters are made. Compared with experiment 1 and experiment 2, the optimization results are improved with the increase of iteration times. Comparing experiments 1, 3, and 4, the optimization effects are improved with the increase of population size but the time increases. The population selection of 100 and the number of iterations of 100 are taken as the hyperparameters of the algorithm.
Through large numbers of experiments, subsequent experiments decided to set improved NSGA-II parameters: population size is 100, number of iterations is 100, and crossover probability is 0.7.

Analysis of DG Installation Capacity Optimization Results
The DG optimization configuration results of the improved NSGA-II are shown in Table 2. As shown in Table 2, Case-1 has no DG, Case-2 is installation without short-circuit current constraints (other constraints are considered), and Case-3 is installation with short-circuit current constraints (other constraints are considered). Compared with Case-1, the planning of Case-2 and -3 obtained by INSGA-II can effectively improve the objective function. Case-2 can produce energysaving benefits of about $0.197 million, the voltage deviation is improved by 73.24%, and the line loss is reduced by 69.60%. Case-3 can produce energy-saving benefits of about $0.193 million, the voltage deviation is improved by 74.83%, and the line loss is reduced by 60.49%. Comparison of Case-2 and -3 show that, under the influence of short-circuit current protection constraints, the capacity of DGs Multi-objective function parameter setting: The maximum installation capacity of each node of DG is 1 MW, the annual maximum load utilization time and the annual maximum line loss hour of DGs are 4500 h, the investment cost of micro-turbine (MT) unit equipment is 0.07 million $/kW, annual operation and maintenance cost is 0.009 $/kWh, wind turbine generator (WTG) group unit equipment investment cost is 0.114 million $/kW, annual operation and maintenance costs are 0.004 $/kWh, MT and WTG unified as a PQ model for processing, power factor is 0.9, the unit price of the loss is 0.071 $/kWh, the financial subsidy per unit of electricity is 0.019 $/kWh, and the conversion factor of the annual equipment investment cost of the distributed power source is obtained by 3% of the annual interest rate. We use the algorithm of the proposed INSGA-II for the following series of analysis.

Analysis of Optimization Results with Different Parameters
As shown in Table 1, in order to improve the operation speed of the algorithm and the accuracy of the solution, experimental comparisons of the experimental hyperparameters are made. Compared with experiment 1 and experiment 2, the optimization results are improved with the increase of iteration times. Comparing experiments 1, 3, and 4, the optimization effects are improved with the increase of population size but the time increases. The population selection of 100 and the number of iterations of 100 are taken as the hyperparameters of the algorithm. Through large numbers of experiments, subsequent experiments decided to set improved NSGA-II parameters: population size is 100, number of iterations is 100, and crossover probability is 0.7.

Analysis of DG Installation Capacity Optimization Results
The DG optimization configuration results of the improved NSGA-II are shown in Table 2. As shown in Table 2, Case-1 has no DG, Case-2 is installation without short-circuit current constraints (other constraints are considered), and Case-3 is installation with short-circuit current constraints (other constraints are considered). Compared with Case-1, the planning of Case-2 and -3 obtained by INSGA-II can effectively improve the objective function. Case-2 can produce energy-saving benefits of about $0.197 million, the voltage deviation is improved by 73.24%, and the line loss is reduced by 69.60%. Case-3 can produce energy-saving benefits of about $0.193 million, the voltage deviation is improved by 74.83%, and the line loss is reduced by 60.49%. Comparison of Case-2 and -3 show that, under the influence of short-circuit current protection constraints, the capacity of DGs decrease and voltage deviation and line losses are improved. In addition, when the DGs are close, the short-circuit current of some branches will be very high, so the Case-3 optimization results do not show that the DG access points are close. The above results confirm the rationality of considering short-circuit current constraints.
As shown in Figure 4a,b, Case-2 and -3 can improve the voltage amplitude of each node and line losses by installing DG appropriately. Because most of the installation nodes of Case-3 are terminal nodes of the system, the supporting effect of the node voltage is stronger [31]. Active network loss is usually determined by node voltage and branch resistance. Distribution network voltage level is lower and R/X value is larger, so it will lead to larger network losses. When DG is reasonably connected, the voltage fading will be improved [31]. In the case of the conditions, Case-3 has a greater effect on node voltage support than Case-2, so the line losses are better. As shown in Figure 4a,b, Case-2 and -3 can improve the voltage amplitude of each node and line losses by installing DG appropriately. Because most of the installation nodes of Case-3 are terminal nodes of the system, the supporting effect of the node voltage is stronger [31]. Active network loss is usually determined by node voltage and branch resistance. Distribution network voltage level is lower and R/X value is larger, so it will lead to larger network losses. When DG is reasonably connected, the voltage fading will be improved [31]. In the case of the conditions, Case-3 has a greater effect on node voltage support than Case-2, so the line losses are better.

Performance Analysis of INSGA-II Algorithm
In order to verify the optimization performance of the proposed algorithm, INSGA-II, NSGA-II, and MOPSO are used to analyze and compare the planning results of the IEEE 33-bus system independently. The parameters of NSGA-II are set as INSGA-II, the inertia factor is 0.8, the local speed factor is 0.1, and the global speed factor is 0.1; other parameters are the same as INSGA-II.
As shown in Figure 5, NSGA-II and MOPSO have uneven distribution of Pareto solutions due to the defects of crowding distance operators. MPSO optimizes by choosing leaders and by archiving mechanism. Similar to NSGA-II, its selection process also depends on a dominant relationship, so

Performance Analysis of INSGA-II Algorithm
In order to verify the optimization performance of the proposed algorithm, INSGA-II, NSGA-II, and MOPSO are used to analyze and compare the planning results of the IEEE 33-bus system independently. The parameters of NSGA-II are set as INSGA-II, the inertia factor is 0.8, the local speed factor is 0.1, and the global speed factor is 0.1; other parameters are the same as INSGA-II.
As shown in Figure 5, NSGA-II and MOPSO have uneven distribution of Pareto solutions due to the defects of crowding distance operators. MPSO optimizes by choosing leaders and by archiving mechanism. Similar to NSGA-II, its selection process also depends on a dominant relationship, so crowding distance is also needed to participate in it. When solving the high-dimensional nonconvex problem, the solution set is easy to fall into the local optimal solution and the population diversity is insufficient. Obviously, the uniformity of Pareto-optimal set obtained by INSGA-II is the best of the three algorithms, which verifies the validity of the potential crowding distance and improves the mutation operator to improve the optimization ability of the algorithm.

Performance Analysis of INSGA-II Algorithm
In order to verify the optimization performance of the proposed algorithm, INSGA-II, NSGA-II, and MOPSO are used to analyze and compare the planning results of the IEEE 33-bus system independently. The parameters of NSGA-II are set as INSGA-II, the inertia factor is 0.8, the local speed factor is 0.1, and the global speed factor is 0.1; other parameters are the same as INSGA-II.
As shown in Figure 5, NSGA-II and MOPSO have uneven distribution of Pareto solutions due to the defects of crowding distance operators. MPSO optimizes by choosing leaders and by archiving mechanism. Similar to NSGA-II, its selection process also depends on a dominant relationship, so crowding distance is also needed to participate in it. When solving the high-dimensional nonconvex problem, the solution set is easy to fall into the local optimal solution and the population diversity is insufficient. Obviously, the uniformity of Pareto-optimal set obtained by INSGA-II is the best of the three algorithms, which verifies the validity of the potential crowding distance and improves the mutation operator to improve the optimization ability of the algorithm.  In order to verify the performance of the improved algorithm more intuitively, the convergence metric [29] and the spacing metric [32] are used to compare the three optimization algorithms.
(1) C_mertic Convergence metric can be expressed as follows: where X and X are respectively two solution sets, a and a are solutions belonging to two solution sets, p is the dominant symbol, and |X | is the number of solutions in the solution set X . If I C (X , X ) is 1, all solutions in the solution set X dominate the solution in the solution set X , and if I C (X , X ) is equal to 0, all solutions in the solution set X are not dominated by the solution set X . As shown in Table 3, about 32.39% of NSGA-II and 33.04% of MOPSO are dominated by INSGA-II. In contrast, only about 6.41% and 7.92% of INSGA-II are dominated by NSGA-II and MOPSO, respectively. In addition, the computational complexity O (mn 3 ) is similar, so the running time of the three algorithms is almost the same, which proves that the improved crowding distance operator can effectively improve the quality of the solution set and can ensure operation efficiency. (2) S_mertic Spacing metric [27] can be expressed as follows: where N represents the number of solutions, d i represents the shortest distance from the ith individual to the rest of the solution, and d represents the mean of all individuals d i . A smaller spacing metric means that the solution distribution in the Pareto solution is more uniform. The zero value of the interval metric means that all solutions in the Pareto-optimal set are equally spaced.
As shown in Figure 6, (3) Analysis of the relationship between objective functions As shown in Figure 7a, the energy-saving benefit has a nonlinear relationship with line loss, and the degree of line loss improvement will have an extreme point, which will not decrease with the increase of DG capacity but will be damaged. In view of Figure 7b, the voltage increases with the increase of energy-saving benefits. The results of DG optimization show that the capacity of DGmounted nodes is too large to reduce the total voltage deviation correction effect. Therefore, the energy-saving benefits have a linear relationship with the voltage deviation when the permeability of DG is satisfied. Line loss and voltage deviation also show a nonlinear relationship in Figure 7c, and the two restrict each other.  (3) Analysis of the relationship between objective functions As shown in Figure 7a, the energy-saving benefit has a nonlinear relationship with line loss, and the degree of line loss improvement will have an extreme point, which will not decrease with the increase of DG capacity but will be damaged. In view of Figure 7b, the voltage increases with the increase of energy-saving benefits. The results of DG optimization show that the capacity of DG-mounted nodes is too large to reduce the total voltage deviation correction effect. Therefore, the energy-saving benefits have a linear relationship with the voltage deviation when the permeability of DG is satisfied. Line loss and voltage deviation also show a nonlinear relationship in Figure 7c, and the two restrict each other.
As shown in Figure 7a, the energy-saving benefit has a nonlinear relationship with line loss, and the degree of line loss improvement will have an extreme point, which will not decrease with the increase of DG capacity but will be damaged. In view of Figure 7b, the voltage increases with the increase of energy-saving benefits. The results of DG optimization show that the capacity of DGmounted nodes is too large to reduce the total voltage deviation correction effect. Therefore, the energy-saving benefits have a linear relationship with the voltage deviation when the permeability of DG is satisfied. Line loss and voltage deviation also show a nonlinear relationship in Figure 7c, and the two restrict each other.

Case Analyses of IEEE 33-, 69-, and 118-Bus Systems
In order to prove the improvement degree of the NSGA-II algorithm on network optimization results of different node systems, the IEEE 33-, 69-and 118-bus systems are introduced. The parameters of the IEEE 69-and 118-bus systems are described in References [33,34]. Table 4 shows the optimization results of the INSGA-II algorithm on three examples. For the objective function, with the increase of network complexity, the improvement degree of the objective function of the three DNs increases. For example, IEEE 33-bus, 69-bus, and 118-bus are improved by 74.38%, 77.39%, and 77.32%, respectively, in terms of voltage deviation. On the hand, it also proves the importance of reasonable installation of DG. On the other hand, most of the installation nodes are the end of the system, which also proves that DGs can improve the voltage deviation better than other nodes.

Case Analyses of IEEE 33-, 69-, and 118-Bus Systems
In order to prove the improvement degree of the NSGA-II algorithm on network optimization results of different node systems, the IEEE 33-, 69-and 118-bus systems are introduced. The parameters of the IEEE 69-and 118-bus systems are described in References [33,34]. Table 4 shows the optimization results of the INSGA-II algorithm on three examples. For the objective function, with the increase of network complexity, the improvement degree of the objective function of the three DNs increases. For example, IEEE 33-bus, 69-bus, and 118-bus are improved by 74.38%, 77.39%, and 77.32%, respectively, in terms of voltage deviation. On the hand, it also proves the importance of reasonable installation of DG. On the other hand, most of the installation nodes are the end of the system, which also proves that DGs can improve the voltage deviation better than other nodes.

Experiments on Load and DG Output Considering Annual Time Series Changes
With the increase of installations and capacity of DG in the future, considering only static load and DG operation status, it may not be applicable to other load levels of permeability constraints and may even lead to reverse power flow, voltage limit, original system protection measure failures, and other phenomena, so that the system is threatened.
To this end, this section optimizes the configuration scheme to consider the DG annual output change and load change trend and the revised IEEE 33-bus system-related data reference [34]. The data was recorded every hour, so it is assumed that the DG output is equal every hour and that the load demand is constant. The four-season output level and load demand of wind turbines are shown in Figure 8a,b.

Experiments on Load and DG Output Considering Annual Time Series Changes
With the increase of installations and capacity of DG in the future, considering only static load and DG operation status, it may not be applicable to other load levels of permeability constraints and may even lead to reverse power flow, voltage limit, original system protection measure failures, and other phenomena, so that the system is threatened.
To this end, this section optimizes the configuration scheme to consider the DG annual output change and load change trend and the revised IEEE 33-bus system-related data reference [34]. The data was recorded every hour, so it is assumed that the DG output is equal every hour and that the load demand is constant. The four-season output level and load demand of wind turbines are shown in Figure 8a   As seen in Figure 8a,b, the wind speed is higher in spring and winter and the output of wind turbine is also improved. On the contrary, the output of wind turbine in summer and autumn is small and the seasonal data of the load has changed but the overall change is not large. On the one hand, if the optimization is based on the summer and autumn season time data, it may lead to waste of DG capacity installation. On the other hand, optimizing the configuration in spring and winter will destroy the network penetration rate and cause economic and security crisis. Here, we choose to optimize the four seasons separately and consider them comprehensively. The four seasonal DG optimization results using INSGA-II are shown in Table 5. As shown in Table 5, due to the relatively stable output of the wind turbines in summer and autumn, it is possible to seek a relatively high-quality configuration result. The spring and winter seasons are relatively unsatisfactory because of the large fluctuations in output. The combined installation capacity of the dynamic optimization scheme should meet the DG permeability constraint of the remaining time period. In each season, it is quite difficult for a DG configuration scheme to improve the optimization function while satisfying the different constraints of 24 h. Therefore, the installation capacity of each season is almost close to the maximum installation capacity, which is a compromise that is comprehensively considered to optimize the time-series data of each objective function. In order to meet the annual operating constraints, the combination of the smallest installation capacity on each quarter node is selected as the planning scheme.

Conclusions
Considering the influence of DG access on the economy, security, and reliability of DN, this paper establishes a multi-objective optimization model to minimize the line losses and voltage deviation and to maximize the annual energy benefit. For DN with integrated DG units, equality constraints and inequality constraints, which include power-flow equality, nodal voltage, branch current, DG capacity, and short-circuit current, are considered in the optimization model. The principle of installation of DG is summarized through experiments and power system theory, and the violation constraint indicators are proposed to solve the infeasible solution.
The mutation operator, crowding distance operator, and selection operator of traditional NSGA II are improved to increase the population diversity and consistency in the optimal allocation. The introduction of an unbiased, compromised strategy can quickly give decision makers a comprehensive plan to weigh each goal.
The different examples of the IEEE 33-, 69-, and 118-bus systems are analyzed and optimized, and different super parameters are compared experimentally. The performance of the proposed INSGA-II is verified by comparing NSGA-II and MOPSO algorithms. The IEEE 69-and 118-bus systems are selected as the verification case. The proposed algorithm is applied to solve the static and dynamic characteristics of the DN with the integrated DGs. Finally, the IEEE 33-bus system is modified and optimized by using the four season data of wind turbine and load. The results indicate that the proposed method can achieve better precision and diversity and can provide a good configuration plan for decision makers under the premise of meeting the annual penetration rate.
In practice, the choice of the best site may not always be feasible due to many reality constraints. However, the optimization and analysis here suggest that considering multi-objectives helps to decide siting and sizing of DG units for the decision maker. The optimization planning of a multi-distributed generation access distribution system combined with geographic information system (GIS) technology will be further investigated.

Acknowledgments:
The authors would like to thank the editors and the reviewers for their constructive comments and suggestions.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

Abbreviations
Z NODG cost total annual cost of DN without DGs Z NODG cost total annual cost of DN with DGs C loss annual loss of DN without DGs C b annual purchase cost of DN without DGs. C p unit price of electricity consumed per unit τ max annual maximum load loss hour (h) of the DN ∆P j active power loss (kW) of the jth branch P load DN total active load P loadi the active power after installing DG T max DN annual maximum load utilization hours C dgm annual investment and maintenance cost of DGs C ploss annual cost of DN with DGs C B annual purchase cost of DN with DGs C sub financial subsidy of new energy generation C dgi distributed generation equipment investment cost ($/kWh) C mi annual operation and maintenance cost ($/kWh) for the ith distributed generation P dgi installed capacity (kW) of the ith distributed generation C sp subsidy amount ($) of distributed generation unit g ij admittance of the branch (i, j) U speci f ied i voltage-rated value P DGi Q DGi ith node active and reactive power of the DGs injected P Li Q Li ith node active and reactive powers of load output B ik the branch susceptance S DGi the total capacity of distributed generation S L maximum load value of DN K III sen,i the sensitivity coefficient of current III segment protection of branch i as backup protection I I