Placement and optimal size of DG in the distribution network based on nodal pricing reduction with nonlinear load model using the IABC algorithm

The growing use of distributed generation (DG) at the distribution level has led to a change in the status of distribution networks from a passive network to an active network such as transmission systems. Therefore, transmission network pricing method such as nodal pricing could be used in the distribution network. DG connection to the distribution network affects bus nodal pricing. If the DG presence reduces losses and congestion in the distribution network, nodal pricing will also decrease. This paper presents a method for calculating the optimal size and place of DG in the distribution network based on nodal pricing. This planning is done to maximize the profits of distribution companies that have used DG in their network to meet several advantages. The simulation was performed using the improved artificial bee colony algorithm (IABC). In the IABC algorithm, by exchanging the received information between bees according to Newton and gravity laws, it uses all this algorithm capacity to find the ideal answer by considering the constraints applied to the system. In most DG placement articles, network loads are assumed to be constant. Because loads are often sensitive to voltage and frequency, constant load analysis leads to inaccurate results. Therefore, in this paper, the proposed method is implemented on a 38-bus radial distribution system with a model of real loads sensitive to the voltage and frequency of the system, including residential, commercial, and industrial loads.


Introduction
Today, the structure of the electricity industry in most countries is changing and moving towards creating competition in the purchase and sale of electric power. This issue, along with technological advances, environmental and economic issues, and the construction of small power generation units, has led to the expansion of the use of distributed generation [1]. Distributed generation sources can be defined as sources of electric power generation that are connecting to distribution networks and local consumers. The production capacity of these sources is smaller than other conventional sources of electric power production and the technology used to generate power in them is also very different and very diverse [2,3]. Due to the growing demand for energy in the power system, the use of DGs has many advantages, which include reducing power losses and system costs, improving voltage and power quality profiles, preventing system upgrades, delaying investment in distribution and transmission networks, greenhouse gases reduction, improve the integration, reliability, and efficiency of the system [4].
On one hand, the distribution networks location between production and transmission, and on the other hand, the load centers existence, have made the distribution network relatively passive, but the growing use of distributed generation resources in the distribution network in recent years has to change the status of these networks from a passive network to an active network such as a transmission network. Therefore, transmission pricing mechanisms such as Nodal Pricing could be used in the distribution network. Nodal pricing or locational marginal pricing (LMP) is an economically effective pricing method that shows the marginal price of electricity in network nodes or buses. DG connection to the distribution network effects at the nodal pricing of the bus. If the DG presence reduces losses and *For correspondence congestion in the distribution network, nodal pricing will also decrease [5].
Distribution networks are traditional passive networks that are not designed to install generators. Therefore, when connecting a production program in the distribution network, technical issues such as the phenomenon of steadystate voltage increase have to be considered, which results from the connection of generators at low voltage levels. The allowable voltage limit for these systems is assumed to be between 1000 V and 132 kV (±6%) of the rated voltage [6,7]. Installation of DGs in non-optimal locations leads to increasing losses, the phenomenon of increasing voltage, and increasing distribution costs for network users [8]. Therefore, distributed generation sources should be optimally installed in the distribution network to maximize their benefits.
There are many optimal DG placement methods, most of which consider the system power loss as an objective function [9]. In [10], the optimal placement of DG has been performed to reduce losses in the 33-bus IEEE system with the PSO algorithm, assuming linear load changes. The place and optimal size of DG are using loss sensitivity factor (LSF) and simulated annealing algorithm (SA) in [11]. In Reference [12] Optimal DG allocation has been implemented to reduce losses by the Monte Carlo method. In Reference [13] the analytical method for multiple DG placement is proposed to reduce losses by genetic algorithm (GA). Minimizing power losses is an acceptable goal. However, this goal is not always economic, as it does not include the total cost of operation. References [14,15] suggest optimal DG placement based on LMP. In which single DG allocation has been performing, but it is obvious that multiple placements of small capacity DGs are more beneficial than the single installation of large capacity DG. This paper considers only the allowable voltage limit of the buses as a constraint on the formulation of the problem, while the phenomenon of voltage rise, which is an important obstacle when installing DG in the distribution network, is not taken into account. Optimal DG allocation based on nodal pricing using the genetic algorithm is suggested in [16].
Most research articles have considered the constant loads in the problem of optimal DG placement [17][18][19] while the system loads are uncontrollable and depend on the voltage and frequency of the system [20]. Therefore, optimal DG allocation, assuming a constant load, will lead to contradictory and misleading results. Therefore, the assumption of a constant load must be reconsidered. A combined method based on the imperialist competitive algorithm (ICA) and the GA algorithm for simultaneous placement of DG and the capacitor bank is proposed in [21]. This paper presents an effective method based on nodal pricing using the IABC algorithm with terms of load model for optimal location and size of multiple DG. This planning has been implemented to maximize the profits of distribution companies that have used DG in their network to achieve several advantages. IABC is an intelligent method that has a high convergence rate and is less to get stuck at local minimum points and find the best solution. Here the profit is calculated based on the reduced cost of power losses in the presence of DG, the cost of electricity without DG installation, and by combining DG including the cost of electricity provided by DG. The proposed method is testing on the 38-bus radial distribution network according to the load model. To confirm the proposed method, the simulation results are compared with other references and at the end, by introducing the system operation indicators, the result of implementing the proposed method in single DG and multiple DG placement modes is evaluated. The structure of the article is as follows: in section 2 modeling and problem formulation is presented, in section 3 the proposed algorithm and problem-solving method are introduced, in section 4 the simulation results are shown, and finally in section 5 the conclusions are presented.

Modeling and formulation of Problem
This paper defines the objective function of the problem based on nodal pricing reduction. This optimization is performed by assuming nonlinear loads that are sensitive to the voltage and frequency of the network so that the practical frequency of the system is assuming to be 98% of the nominal frequency of the network. The PQ model is used for DG modeling; In other words, DGs are considered as negative loads and the production cost of DGs is assuming to be 45 $/MWh while the cost of electricity received from the network (price in a slack bus) is assuming to be 44.5 $/MWh.

Nodal pricing
The node price represents the final or marginal price of electricity in the network nodes [22] or, in other words, the node price reflects the marginal cost of delivering one more MWh to each bus in the system. This marginal cost takes into account the effect of power injection and delivery on marginal losses (incremental changes in losses on power). In line loading, this cost increases due to consumption and increased distance from Power Supply Point (PSP) to load. Thus, nodal prices are lower for buses close to the PSP and higher for buses far from the PSP, until they show an increase in marginal losses due to distance from the PSP [23]. Thus, in the distribution network, nodal pricing for active and reactive power is obtained from Equations (4) and (5), respectively: Where q Pi and q Qi marginal loss components (MLC) are active and reactive at each node i in the network. These components are defined as changes in the total active power loss L on the marginal change in the consumption/production of active power Pi and reactive power Qi at each node i in the network. q Qi for i is a PV node is zero and q Pi and q Qi for i is the slack node will be zero. k also shows the price of active power in the PSP.

Modeling of distributed generation resources
DG sources could be modeled by considering the type of DG, operating, and network connection method to PV or PQ buses [24]. In inverter-based DGs, the converter control method determines the type of modeling, while in machinebased DGs, the machine operating mode determines the type of modeling. A DG unit in the PQ node mode is modeling into three different types: In the first type of modeling, the DG is simply modeled as a source of constant active and reactive power, which could be considered as a negative load; In the second type of modeling, DG units have certain values of active power and power factor and are modeling as a machine with constant power factor. In the third type of modeling, DG sources are modeled into variable reactive power generators (such as an induction machine-based wind farm). When DG is modeling like a PV machine, DG units have a certain amount of active power output and bus voltage size [25,26]. To keep the voltage constant at a certain value, the required bus voltage changes must be zeroed by injecting the required reactive power. Since DGs are usually smaller compared to traditional power sources, in this paper a fixed PQ model is used to analyze the power flow of the distribution system, in which case the DG will be considered as a negative load and the production cost of DGs is assuming be 45 $/MWh.

Modeling of load
In traditional power flow analysis, the active and reactive loads in the nodes are usually assumed constant, regardless of the voltage and frequency of the system, while in the practical operation of the power system the real load models (i.e. residential, industrial and commercial loads) depend on voltage and frequency of the systems. The characteristics of different types of loads through exponential load models are shown in Equations (6) and (7): Where P Di and Q Di are active and reactive power in nodes i, P 0i and Q 0i are the amounts of active and reactive power at nominal voltage in each node i, V i is the amount of voltage in bus i, c and s Voltage indicators for active and reactive powers, f and f 0 real and nominal frequency of the network, a and b are frequency coefficients for active and reactive powers. If it is s = c = a = b = 0 in Equation (6) and Equation (7), the proposed load model will be a constant load model used in traditional power flow. The values of s, c, a, and b for different types of loads are showing in table 1 [27]. These load coefficients can be extracted by studying the network over a specified period schedule. Of course, it should be noted that these coefficients may vary in different networks.

Problem formulation
Assume that C i a (no-DG) and C i r (no-DG) are the prices of the active and reactive power nodes in terms of $/MWh and $/MVAr in each node i without DG. P L (no-DG) is the active power loss in MW without a DG connection. P Di and Q Di demand active and reactive power in each node i in terms of MW and MVAr. k The price of grid-supplied electricity or power supplied to the PSP is in $/MWh. The p (no-DG) electricity bill price without DG for each period Dt is obtained according to Equation (8) [28]: Where n is the number of buses and Dt is the period in hours. Since the losses are not dependent on any particular node, the cost of losses is calculated based on the power price in the PSP (ie k) from Equation (3). Assume that C i a (DG) and C i r (DG) are the prices of the active and reactive power nodes in terms of $/MWh and $/ MVAr in each node i with DG. P L (DG) is the active power loss in MW with a DG connection. P Di and Q Di active and reactive power provided by DG in each node i in terms of MW and MVAr. C(DG) the electricity prices supplied by DG in $/MWh. The p (DG) electricity bill price with DG for each period Dt is obtained according to Equation (9): DG is considered a negative load. In other words, the PQ model is used. P DGi and Q DGi are zero in all buses except DG buses. So : Profit is calculated based on savings on electricity bills in the presence of DG, including the price of electricity supplied by DG. Therefore, the problem formulation could be expressed as Equation (11): The following constraints are also included in the placement issue: V 0;i V max 8 secondary side of LTC transformers ð12Þ V end;ijmax load;noDG ! V min 8 nodes at feeders ends ð13Þ S ði;jÞ S ði;jÞmax ð17Þ P DGmin P DGi P DGmax ð18Þ Where V 0,i is the voltage per PSP and V max is the maximum allowable voltage range. V end, i|max load, no DG the voltage at the end of the feeders at maximum load and without DG and V min indicates the minimum allowable voltage range. V DG|min load, max DG indicates the voltage at the DG place with the minimum load and maximum DG penetration, and the voltage of each bus i is denoted by V i . Equations (12) to (16) are about voltage constraints that state the problem of voltage rise in each feeder set with Load tap-changer (LTC). S (i,j) is the apparent power passing, and S (i,j)max is the power flow capacity between nodes i and j in terms of MVA. Equations (18) and (19) also express the limit of DG production capacity.

Methodology
The first step in solving the problem is awareness of the status of the distribution network and find the weaknesses of the network in terms of voltage stability. By solving the power flow problem, the required parameters of the problem could be included: voltage amplitude and angle in all buses, the passing current of each line, losses of each line, total feeder input power, total feeder losses, active and reactive power of each load, the apparent power of the feeders obtained. Therefore, according to previous mentions, for power flow analysis, the following items are important. (1) Determining a bus as a reference bus along with specifying its voltage.

Optimal power flow
In terms of structure and topology of transmission and distribution networks, some differences have caused the method of solving the problem of optimal power flow (OPF) in these two networks to be different. The structure of the transmission network rings, but the structure of the distribution network is radial. In this paper, to find the weakest bus of the network in terms of voltage stability and network losses, a voltage stability index (VSI) is introduced. By equating the n-bus distribution network of figure 1a in the form equal to two buses of figure 1b, the voltage stability index of the bus i is extracted by Equation (20).
Where |V 1 | is the magnitude of the busbar voltage and R eq and X eq are equivalent line resistances and reactances, respectively, calculated by Equation (21). Where R is the sum of the path from a substation to bus i, I i-1,i represents the current from bus i-1 to bus i, |V i | voltage magnitude bus i, and u i the angle difference between the voltage and the load demand current in the bus i. It is assumed that bus 1 is the substation as the reference bus and d 1 = 0. The voltage vector relations of figure 1c will be Equations (23) and (24).
Where d i is the angle of the bus i voltage. Equation (23) and Equation (24) could be rewritten by multiplying their sides in |V i | by Equation (26) and Equation (27), respectively.
The VSI index in Equation (28) depends only on the voltage of substation and bus i also on the bus i voltage angle. In this method, it is no longer necessary to calculate the bus admittance matrix, bus impedance matrix, and equivalent impedance to node i to substation (i.e., Equation (21)). The above advantages will make this index suitable for real-time calculations and the calculation of this index requires only one power flow. Also, the stability limit of the bus i voltage is when VSI i = 0. Therefore, by solving the problem of power flow without DG, the weaknesses of the network could be identified with the help of the VSI index and considered as candidate points for installing DG in the network, and then the problem of power flow with DG could be solved. Figure 2 shows the flowchart of the problem-solving method. As you can see, the optimal solution to the problem is obtained in two steps. In the first stage, the candidate places for installing DG in the network have been extracted and in the second stage, the problem has been solved with the presence of DG in the network. Finally, the output of the problem will include the maximum profit, the optimal place and size of the DG, the node price of the active and reactive power of the network, the reduction of losses, and the improvement of the network voltage profile.

Improved artificial bee colony algorithm
In this paper, the IABC algorithm is used for the optimal placement of DG in the distribution network. The ABC algorithm is a technique for solving optimization problems based on the behavior of bees. In this method, each of the bees, by direct cooperation and sharing information, tries to get the best answer according to the rules of probability. Bees of a complex system to find out the place and quality of food sources outside their hive. First, a set of food sources is selected randomly (initial answers), then the worker bees refer to the sources and check the amount of honey and their quality and return to the hive and give their information to the spectator bees. Then each bee moves to the place and, based on the information, selects a source in its neighborhood, that is, the bee decides whether to stay in the new place or to go to the previous place based on the type of flower and the amount of honey. When the resource is depleting, they move to a new source found by the bees, and this process is repeated until the needs are meting. The algorithm encoding is as follows.
Initialization: The value of X ij is selected as the initial answer in the search space of the algorithm and the value of their competence is examined based on the objective function studied, each of which is introduced as a source of honey used. The choice of these answers was random in the search space and represents worker bees.
Movement of explorer bees: The probability for the selected sources is based on Equation (29) and the selection of a food source using the choice of the roulette wheel for each explorer bee and determining the honey amount for each of them with the model developed based on the reciprocal gravity force between the explorer bees is obtaining in Equation (30) to Equation (34).
The mutual force between two masses in m2, m1 is in the form of Equation (30). Which is shown in figure 3 of the interaction of these forces.
b r 21 ¼ r 2 À r 1 r 2 À r 1 j j ð31Þ Where F 12 , r 21 , and G are mutual forces, unit vector, and constant gravity. In the same way for bees, we express the following equations based on their suitability.
In Equation (32) F(h i ) and F(h j ) are the amount of competency introduced for worker bees, respectively, which results in its effect on the new feed source in the form of Equation (33). Also, considering the interaction of all bees, Equation (33) is expressed by expanding as Equation (34) [15].
In fact, in the IABC algorithm, the criterion used is instead of the random value, which increases the efficiency of this algorithm.
The motion of Pioneer Bees: If the Fitness value of the function is not corrected in the next iterations of the algorithm, it is called Limit and those resources are called abandoned, which is used by the movement of pioneer bees to retrieve and replace new resources. The movement of these bees is in the form of Equation (35).
Placement: If the food source found in later stages is better than before, this amount will enter the bees' memory. And finally, this program continues to meet the optimal answer. Figure 4 shows the flowchart of the IABC algorithm. Finally, the method presented in figure 2 is implemented by the IABC algorithm as follows: Step 1: Create the initial population randomly. Each member of the initial population is the network bus voltage.
Step 2: Check the voltage and other constraints for each member of the population at all network nodes.
Step 3: Extract the network node voltages and DG installation candidate places.
Step 4: Create new members. Each member of the new population includes the size of the DG and the candidate places for DG installation on the network.
Step 5: Apply each member of the new population to the distribution network, i.e. the size of DG produced in each member at the relevant bus place to the network.
Step 6: Calculate the value of the objective function for each member of the population. The member that yields the value of the optimal objective function will be the optimal member. Step 7: If the convergence criterion is fulfilled, stop the operation and show the results, otherwise repeat with the new members until the optimal answer is reached.

Simulation results
In this paper, the proposed method is implemented on a 38-bus radial distribution network. The single-line diagram of the network is shown in figure 5. Network information is available at [29]. The price of electricity on the PSP is assumed to be 44.5 $/MWh. The proposed method is performed in two modes of simulation of optimal single DG allocation and multiple DG in terms of load model. Similar to reference [30][31][32][33] is assuming that the distribution company (DISCO) owns the DG and the DISCOs profit is calculated based on the savings in the electricity bill, which includes the cost of electricity supplied by the DG. The cost of electricity supplied by DG is equal to 45 $/MWh.

Single DG placement
The proposed method is implemented by a single DG placement for different types of loads in the 38-bus network. The simulation is first implemented with the same type of total loads (such as, the total loads are industrial) and then in the mixed load mode. The results are also compared with the constant load's assumption, which is shown in table 2 of the single DG placement results. As could be seen, the simulation results differ in the presence of voltage and frequency sensitive loads with a constant load mode. The results should be reviewed assuming a constant load. In figures 6 and 7 the price of each unit of active and reactive power in the reference bus (power supplied by the network) in each node without the DG presence and with the optimal placement of DG is calculated. These results are shown in the case of mixed loads. The nodal pricing of the active and reactive power of each node has decreased significantly with the presence of DG.

Multiple DG placement
The proposed method could be extended to more DG placements. Implemented the simulation to placements 2 and 3 DG units on the test system. Tables 3 and 4 show the optimal place and size of the respective DG and total power losses by installing 2 and 3 DG units with the nonlinear load's presence (mixed loads). It is observed that multiple DG placement is more beneficial because the profit increases by 25.3% in two-DG placements and 33.2% in three-DG placements compared to single-DG placements.
On the other hand, active and reactive losses, as seen in tables 3 and 4, are further reduced compared to single-DG placement (losses compared to single-DG placement, in the placement of two-DG units 16.1% and the placement of the three-DG units decreased by 28.2%). The amount of reduction in the active power node price in the buses with the connection of one DG, two DG, and three DG compared to without DG in figure 8. It could be seen that the reduction in the nodal pricing in the buses is better in the case of multiple DG placement than in the single DG placement. Figure 9 shows the total network line losses for the mixed load mode and figure 10 shows the network voltage profile. As you can see, with the help of the stability index introduced in section 3 and the network weaknesses identification and with the DG installation, the results of voltage profile improvement from method [29] are also better.     According to figure 10, the lowest VSI of the system belongs to bus 18-33-37, which is equal to 0.67912, and the voltage amplitude of this bus is 0.9110 p.u. Now, to improve the voltage profile and reduce losses, and then reduce nodal pricing, the tangible effect of DG installation in the network on the voltage stabilization index of bus 18 is obtained according to figure 11. Figure 11 shows that the installation of three-DG in the network was able to increase the VSI of the network buses, especially buses 18-33-37. Usually, in optimization problems, the goal is to minimize the objective function, so in this paper, Equation (11) is optimizing as a (1/max(p no-DG -p DG )) function by the IABC algorithm. Figure 12 shows the convergence rate of the objective function for DG installation in the network.

Assessment of the proposed method
The following performance indicators are introduced to verify the proposed method and to investigate the effect of the load model on distributed generation planning.
Active power loss indices are defined according to Equation (36) [31]: P LDG is active power loss with DG, P L is an active loss without DG. A decrease in the value of this index indicates a reduction in system losses. The penetration level of the DG affects the network voltage profile. Optimal DG installation leads to an improved grid voltage profile. VPI is calculating according to Equation (37) [31]: If V DDG and V D are the mean system voltage deviations from the reference with and without DG installation, the voltage deviation index (VDI) is obtained from Equation (40) [32]: Reducing the VDI value indicates improved voltage regulation in the bus. A value of 1 indicates a state without DG. To present the effect of DG on the nodal price of network buses, the NPI index according to Equation (41) is introduced: The simulation results for all modes (without DG, one, two, and three DG) for the mixed load model are shown in figures 13 and 14. It could be seen that all evaluation indicators have improved with DG. This change is more noticeable in the multiple placements of DG so that compared to the case without DG in the network, in the placement positions of one, two, and three DG units with LPI index, 43.8%, 52.7% and 59.6%, respectively, VDI index 38.6% and 50%, 57.1%, VPI index, 31.5% and 46.4%, 48.4% and NPI index 45.7% and 68.2%, 70.9% decreased. The percentage change of indices between two and three DG units is insignificant. Therefore, the placement of more than three DG units in the sample distribution network is not recommended. Table 5 shows a comparison between some of the methods presented so far in the power network considering several criteria and the percentage of LMP reduction in each method. As you can see from the results, the percentage of LMP reduction of the proposed method is higher than other methods. However, it should be noted that the type of network studied, network size, load model, and type   In [3,6,[37][38][39], a higher reduction was obtained, and this is in the case that the largest network of choice in these methods is the IEEE 33-bus network, and sometimes for the transmission network, this percentage is reduced. This percentage will certainly decrease as the network gets bigger. Therefore, the proposed method has a good performance in reducing LMP. Table 5 shows the dependence of each method on quantities, indices, parameters in reducing LMP. The proposed method uses the network voltage stability index to the DG placement and combines loads (including residential, industrial, and commercial) to make the results real.
Finally, the features of the proposed method can be listed as follows: • Nodal Pricing Reduction in the Distribution Network using the IABC Algorithm. • Consider different load models.
• Calculate the price of active and reactive power in the whole network node with optimal DG placement. • Provide voltage stability index (VSI) to solve the OPF problem. • Extraction of optimal DG installation place in the network, based on VSI.

Conclusion
In this paper, an effective method based on nodal pricing for optimal place and size of multiple DGs in the distribution network by considering nonlinear loads using the IABC algorithm is proposed. The efficiency of the proposed method has been tested on a 38-bus radial distribution network. As observed, the simulation results are different from the assumed constant load by considering the nonlinear load model. It has been observed that an optimized small-capacity DG produces more profit than a non-optimized large-capacity DG. Therefore, placement DG is essential to meet maximum profit and minimum loss. It is also clear that multiple DG placement with small sizes is both technically and economically more profitable than single placement DG with large capacity, as nodal prices, losses, and voltage profiles have been much improved. On the other hand, according to the proposed method evaluation, the indices of losses, voltage profiles, nodal pricing, and voltage deviation in the case of multiple DG placement have changed much compared to the single placement of DG and without DG.

Declaration
Conflict of interest: The authors declare that they have no conflict of interest.