Impact analysis of time‐varying voltage‐dependent load models on hybrid DG planning in a radial distribution system using analytical approach

Correspondence M. L. Merlin Sajini, EEE Department, Coimbatore Institute of Technology, Civil Aerodrome Post, Coimbatore, India. Abstract This paper investigates the impact of time-varying voltage-dependent load models on hybrid distributed generation planning studies. Firstly, a multi-objective index is constructed with indices namely real and reactive power loss index, voltage deviation index, and mega-volt-ampere capacity index with an objective of reduction of real and reactive power loss, voltage profile improvement, and decreased line loading. The multi-objective index derived is minimised to locate the distributed generation. Secondly, the analytical expression is obtained to size the distributed generation units with the time-varying voltage-dependent load models and probabilistic nature of solar and wind generation. Finally, the impact of load models with the generation effects is separately analysed for PV/wind (hybrid) distributed generation planning. The proposed perspective was validated on the IEEE 33-bus, IEEE 69-bus radial distribution network, and a real 16-bus distribution substation for its effectiveness because the distributed generators are placed in the distribution side of an electrical network. The outcomes reveal that the distributed generation allocation for time-varying voltage-dependent load model and probabilistic nature of solar and wind has significant impacts on the distribution system by reducing real and reactive power loss, improving voltage deviation and reduced mega-volt-ampere intake.


INTRODUCTION
Renewable energy based power generation has entered into electric power generation in this era because of increasing power and decreased carbon emission. A small generation unit with a capacity of few kilowatts to 50 MW is defined as a distributed generation (DG) unit. Solar and wind are the naturally available resources in the living planet and therefore solar and wind-based renewable energy technologies have the potential to meet future energy demand [1]. In an electrical network, the DG unit has to be implemented in the distribution side because of its high R/X ratio in distribution network than in the transmission network. As R/X ratio increases, it results in more energy loss throughout the distribution feeders [2]. Therefore DG planning plays a major role in the distribution network. Also, it is indicated that inappropriate location and sizing of DG leads the network system to high power loss and poor voltage profile [3,4]. Therefore by proper analysis, distributed generators have to This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited. © 2020 The Authors. IET Renewable Power Generation published by John Wiley & Sons Ltd on behalf of The Institution of Engineering and Technology be placed in the distribution network for the upgradation of grid, reliable operation, power loss reduction, voltage profile improvement, and to reduce the running cost. Hybrid DG installation in the distribution has its advantage by providing consistent energy generation [5]. Various techniques are discussed to optimally allocate and size the DG units.
A 'golden rule' or '2/3 rule' for optimal location of the distributed generator in radial distribution network was proposed with a concept of installing 2/3 rating of the distributed generator at a distance of 2/3 length of the line. The proposed technique applies only to an evenly scattered load [6]. Despite its simplicity and easy approach, it cannot be used for a feeder with other types of load.
The DG is considered to be located in the primary distribution system, and the objective of DG placement is to reduce the losses. The cost of DG and the other associated benefits have not been considered while solving the location and sizing problem. The sizing and placement of DG are based on single instantaneous demand at a peak, where the losses are more. An analytical approach is discussed [2] to find the appropriate location and size to reduce only the losses by using the exact loss formula. On the other hand, Wang and Nehrir in [7] discussed an analytical approach considering a DG of unity power factor. The optimal bus is found based on a bus admittance matrix, generator data, and load data. In this approach only the appropriate location, to place a single DG to reduce the loss, was found but the optimal sizing of DG was not discussed.
Apart from analytical approaches, optimisation algorithms have been proposed for the DG allocation and sizing problem, in a distribution network. DG allocation and sizing with multi-objective considering real and reactive power losses and deviation in voltage profile have been discussed in [8][9][10][11], only for a constant load and constant generation model. In [12,13], the optimal sizing and placement of DGs have been proposed based on the Hybrid Grey Wolf optimiser (HGWO) and Artificial Bee Colony technique, considering only real power loss reduction, with equality and inequality constraints. In [14], the optimal location and rating of different cases of DG capable of injecting real and reactive power have been proposed based on particle swarm optimisation to minimise real power loss considering the system constraints. An approach based on genetic algorithm [15] is explained for the optimal allocation of DG considering the objective of power loss reduction by increasing the number of DGs at different load levels.
An evolutionary programming based approach with an objective of power loss reduction is proposed in [16]. But in this approach, the probabilistic nature of solar and wind generation was considered for DG allocation and sizing in a distribution network. A new methodology to increase the annual energy loss reduction by optimally allocating the solar and wind-based DGs considering the probabilistic model of generation and load is discussed in [17]. The allocation and sizing problem with an objective only to improve the voltage stability is proposed in [18] considering the probabilistic nature of solar and wind. A new voltage stability factor with objectives of power loss reduction and voltage profile improvement in radial distribution network considering the solar and wind-based DG have been proposed in [19]. But in the above discussions, the impact of the DG installation in the distribution network is not studied. An approach to optimally allocate and size the DG considering uncertainties based on the minimum life cycle cost is discussed in [20]. Minimisation of life cycle cost includes investment and operational cost, emission taxes, and residual equipment cost at the end of the planning period.
A novel technique based on power loss sensitivity and voltage stability index for optimal allocation and sizing of DG units is discussed in [21]. The improvement in voltage profile, reduction in power loss, apparent power intake by the grid, and savings in energy cost are determined for the radial distribution network. A novel analytical method is discussed in [22] for DG placement and sizing with and without P (real power bus) and PQV buses (bus whose voltage magnitude is controlled by real bus). This method does not use the bus impedance matrix and the objective is to reduce the power loss.
The review of the literature confirms that an analytical approach for hybrid DG planning considering system real and reactive power loss index, voltage deviation index, and megavolt-ampere (MVA) capacity index as indices in the multiobjective function, the effect of voltage-dependent and combined load-generation model, probabilistic nature of solar and wind, and various optimal power factor of the hybrid (solar/ wind) DG on DG planning in a radial distribution system simultaneously has not been proposed in any of the research works in recent years. Also, in many of the papers for DG planning, the main objective is loss minimisation and voltage profile improvement. But the objective to reduce line loading with optimal power factor is not considered through an analytical approach.
Thus, the main contribution of this paper is the optimal allocation of hybrid DG units considering the effect of voltagedependent time-varying load models has been studied. The multi-objective index (MOI ) is a combination of four indices: active power loss index (PLI ), reactive power loss index (QLI ), voltage deviation index (VDI ), and MVA capacity index (MVAI ). Based on this an analytical expression is derived for optimal sizing of hybrid DG at different optimal power factors. Finally, the impact analysis for the different loads with the probabilistic generation of solar and wind-based DG unit is done and evaluated.
The rest of the paper is structured as follows. In Section 2, the modelling of load, solar, and wind are discussed. In Section 3, the problem is formulated with the impact indices related to real and reactive power loss, voltage deviation, and MVA capacity, with the inclusion of hybrid DG and without DG. In Section 4, with these impact indices, a multi-objective function based analytical approach is formulated for the optimal sizing and placement of solar and wind and the computation method. Section 5 explains the test system model (i.e. 33 and 69 radial bus distribution networks) and the obtained results. Section 6 explains the contribution of solar and wind in the distribution network and their impacts on the time-varying loads and the conclusion of the work is summarised in Section 7.

SYSTEM MODELLING
The generation and load model considered for the evaluation of performance indices is shown in Figure 1. The historical data of solar irradiance and wind speed modelled with beta probability density function (pdf) and Weibull pdf, respectively, and the load demand data are used with the DG units to frame a combined generation and load model.

Load modelling
In this proposed work as the radial distribution network has to supply time-varying loads, the dynamic load models are considered. As discussed in [23], the three types of dynamic load models such as industrial, residential, and commercial are considered for optimal placement and sizing of DG units. The normalised daily load curve with a peak load of 1 p.u is shown in Figure 2 for the considered load patterns. The mathematical model of the voltage-dependent timevarying load model for a particular time period is given by where P m and Q m are the real and reactive power injection at bus m, respectively, P om and Q om are the real and reactive load at bus m, respectively, V m is the voltage at bus m, e p and e q are the real and reactive load voltage exponents.
If a constant load model is used for load flow studies, the values of e p and e q are equal to 0. But for commercial, industrial, and residential load patterns, the values of e p and e q are tabulated in Table 1.

Distributed generation modelling
In the proposed DG planning approach, apart from IEEE distribution network considered a real 16-bus distribution substation is also considered from the Coimbatore region. The solar and wind generation for the Coimbatore region is predicted by using the historical data. High-resolution hourly GHI (W/m 2 ) and wind speed (m/s) data are used from the NSRDB data viewer provided by NREL, for the Coimbatore region, with the latitude and longitude values 11.0209 • N and 77.0639 • E, respectively. The Coimbatore region is at an elevation of 415 m above the sea level and the wind speed is measured at 10 m above the ground level.

Solar PV modelling
The solar irradiance for every hour of a day can be modelled by the beta pdf for the historical data which has been collected for four years. For four years, each day data has been categorised into 24-time segments and the beta pdf for all time segments is obtained for each hour. Each hour of the day is assumed to have 20 states of solar irradiance using a step of 0.05 kW/m 2 . From the mean and standard deviation of the hourly global horizontal irradiance data, the beta pdf with 20 states for global horizontal irradiance for each hour of the day is determined. Correspondingly, the expected power output for that particular hour is determined. As discussed in [3], the beta pdf can be employed to explain the probabilistic characteristic of the solar irradiance. This method of solar PV model is adopted in several studies [12,[16][17][18]. The probability of solar irradiance for each day can be expressed as where F (s) is the beta distribution function of s, Γ is the gamma function, s is a random variable of solar irradiance in kilowatt per square metre. The shape parameters of beta distribution function which are calculated by using the mean and standard deviation of solar irradiance are as follows: The cumulative distribution function (CDF) can be used to estimate the time probability of solar irradiance during any specific hour. It can be expressed as follows: where s 1 and s 2 are limits of solar irradiance state. The power output of the photovoltaic module at any state s, as discussed in [3], can be expressed as where N is the number of modules, K i and K v are the current and voltage temperature coefficients (A/ • C and V/ • C), respectively. FF is the fill factor of the solar panel.
The average expected power output of a photovoltaic system over a particular period t, P S (t ) for time period t = 1 h, can be obtained from (6) and (7) as discussed in [3]: From the previous four years of hourly data, the mean ( ) and standard deviation ( ) are found which are tabulated in Table A1. The shape parameters and are calculated with the values of and . Therefore, for each hour with an interval of 0.05 kW/m 2 the pdf is generated and shown in Figure 3.
The pdf for each hour is different because the solar irradiance is time-and weather-dependent. With these pdf and from the specification of the PV module, given in Table A2, the expected output of the PV module for the respective 20 states is calculated using Equation (12) and plotted in Figure 4 from 7:00 to 11:00 AM.  Expected PV output 7:00 -11:00 AM Therefore, the average expected power output for the day in an hour is the area under the curve which is in kilowatt-hour, as the period is considered for 1 h. For example, the average expected power output for 9:00 AM from a single solar panel is 287.79 W and the PV module is expected to output power of 287.79 Wh in a day for 9:00 AM.
Thus, by calculating the average expected output power for each hour, the normalised PV output for a day can be calculated accordingly and is plotted in Figure 5.
With the average power (P avg s ) calculated and maximum power output (P max s ) from the installed PV module the capacity factor of PV module can be calculated as discussed in [3]: In this proposed work, the PV power generation is considered to be capable of delivering active power and exporting Normalised PV output for a day or consuming reactive power. As per the standard IEEE1547, PV inverters used with the system are not allowed to generate reactive power to the grid as discussed in [24]. But [25] if reactive power is inadequate, the voltage regulation is disturbed in the existing distribution system. Conventional devices such as switching capacitors, automatic voltage regulators can be used to compensate for the reactive power, but those conventional devices are not fast enough to compensate for the sudden change in events due to the stochastic nature of the PV generation unit. Therefore to enhance a fast response, the inverter-based photovoltaic topology is preferred to inject or absorb reactive power as per the German code [26]. Thus, the active and reactive power of PV unit at a bus j is given by Equation (2) in [2] as Q sj = j P sj (14) where j = ±tan(cos −1 (p f s j )). The value of j is positive for PV unit supplying reactive power and negative for PV unit consuming reactive power. p f s j is the operating power factor of the PV unit at bus j.

Wind modelling
The wind speed is measured at a height of 10 m above the ground level. The landscape type is considered as small towns, villages, with agricultural land and for this landscape the friction coefficient is 0.4 as discussed in [5]. Thus, the wind velocity at any height can be found by Hellman's equation, as discussed in [27]: where v is the speed to height h, v o is the speed to height h o , and is the friction coefficient or Hellman exponent.
The wind speed of every hour can be modelled by the Weibull pdf. Each hour of the day is assumed to have 16 states of wind speed with a step of 0.05 m/s. As discussed in [3], the Weibull distribution for the wind turbine is given by For the wind speed model, the Rayleigh distribution of Weibull distribution is used in which the value of k is considered as 2 and c is considered as 1.128 v m , where v m is the average wind speed specified in Table A2. Therefore, the probability distribution function of wind speed based on Rayleigh distribution is given by The cumulative density function (CDF) of wind speed is given by The CDF can be used to estimate the time probability of wind speed during any specific hour. It can be expressed as follows: Therefore, from the power curve of a wind turbine, it is understood that the power generated from a wind turbine can be represented as The average expected power output of a wind turbine over a particular period t, P w (t ) for time period t = 1 h, can be obtained from Therefore, for each hour with an interval of 0.05 m/s the pdf is generated and shown in Figure 6. With these pdfs and from the specification of the wind turbine, given in Table A4, the expected output of the wind turbine for the respective 16 states is calculated using Equation (21) and plotted in Figure 7. Therefore, the average expected power output from the wind turbine for the day in an hour is the area under the curve which is in watt-hour, as the period is considered for 1 h. For example, the average expected power output for 8:00 AM is 3.54 kW from a single wind turbine and is expected to output power of 3.54 kWh in a day for 8:00 AM. Thus, by calculating the average expected output power for each hour, the normalised wind output for a day can be calculated accordingly and is plotted in Figure 8.
With the average power calculated (P avg w ) and maximum power output (P max w ) from the installed wind turbine, the capacity factor of the wind turbine can be calculated as The wind generation unit considered is capable of delivering active power and export or consuming reactive power. According to [28] IEEE standard 1094-1991, for wind farm generating stations, the generators used requires reactive power (var) for its operation. If the reactive power is not supplied for the wind generator, it draws reactive power from the utility. The common method of reactive power compensation is to provide syn- Normalised distributed generation output for a day chronous generators or power conditioners and shunt capacitors. Thus, the active and reactive power of wind generation unit at bus j is given by where j = ±tan(cos −1 (p f w j )). The value of j is positive for wind generation unit supplying reactive power and negative for wind generation unit consuming reactive power. p f w j is the operating power factor of the wind generation unit at bus j. Thus, the average expected output for a day from combined solar and wind (hybrid) DG unit is shown in Figure 9.

Optimal power factor
In a distribution system to reduce the system power losses, the power factor of the DG unit that has to be installed should be Thus, the optimal power factor of the hybrid DG unit is given as follows:

Combined generation and load model
The generation model of the DG and the model of the load are employed to obtain a combined generation and load model. The wind speed and solar irradiance have a different probability for a particular hour in a day with 20 states for solar irradiance and 16 states for wind speed, respectively. But the load model is constant for each hour in a day and its probability is 1. Thus, the probability of combined generation and load model is the probability of the solar and wind generation itself.

PROBLEM FORMULATION
Consider a radial distribution network shown in Figure 10(a) with n number of buses and (n − 1) branches. For a distribution network, the source is the substation, in which P and Q are the active and reactive power flow through each branch, respectively. P D and Q D are the active and reactive load power at buses, respectively. If a DG unit is located at bus jP DG j and Q DG j are the active and reactive powers injected by DG unit at bus j, respectively. The impact indices considered are as follows.

Active power loss index (PLI)
For a distribution network with n number of branches, the total power loss is given by [29] Now, if a DG unit is placed at bus m as shown in Figure 10(b), the active and reactive powers flow from source to the DG-injected bus is decreased and the power flow in all other branches remains the same. Thus, the real power loss can be rearranged as In general, the active and reactive power of DG unit at a bus j is given by Substituting (25) and (27) in (26) we get Therefore, the active power loss index (PLI) can be given by the ratio of active power loss with and without a solar or windbased DG unit in the bus given by PLI = P LDG P L (29)

Reactive power loss index (QLI)
The total reactive power loss, Q L in the bus system with n number of branches, is given by If a DG unit is inserted in bus j, the reactive power loss can also be written as Substituting (30) and (27) in (31) we get Therefore, the reactive power loss index (QLI) can be given by the ratio of reactive power loss with and without a solar or wind-based DG unit in the bus given by

Voltage deviation index
The change in voltage along the branch from bus i to bus i + 1, with branch resistance and reactance, R i + jX i can be expressed as denoted in [30]: Therefore, for n branches the cumulative voltage change is given by Squaring this equation on both sides: Now, if a DG unit is injected at bus j, the voltage change can be rewritten as Substituting Equations (23) and (36) in Equation (37): Thus, the voltage deviation index for the distribution system is given by

MVA capacity index
The MVA capacity of a distribution network is the maximum capacity that can be handled by the conductor. While supplying power to the nearby loads, the power flow may decrease in some sections of the distribution network, but in some sections, the power flow may increase above the distribution capacity limit. Therefore, the MVA capacity index discusses the branch's current flow through the network regarding the maximum capacity of the conductors. With this information, the upgradation of the distribution line can be made [31]. If the capacity index value is lesser, it is understood that the distribution line can withstand the increase in branch current flow. If the index value is more than unity, the distribution line has to be upgraded to withstand the increase in branch currents. The power intake from the substation is given by the total active and reactive load attached to each bus and the total power losses in the system.
Without DG, With DG, Thus, the apparent power intake from the substation is as follows: Therefore, MVA capacity index for a distribution line is given by where S i is the MVA flow in line i, P i is the real power flow in line I, and Q i is the reactive power flow in line i. C si is the MVA capacity of line i. Similarly, the MVA flow in each line can be determined. Now, if a DG unit is injected in a bus j, the MVA capacity of the distribution network is given by From this Equation (49), the MVA capacity of the distribution line can be found if a DG unit is placed on any of the buses. Finally, the MVA capacity index for the distribution network if a DG unit is injected is given by From the above analysis, the impact of DG units in a distribution network for different load models is determined by finding the real and reactive power loss, apparent power intake for the distribution network for each hour in a day without DG units and with DG units. Impact indices such as real power loss index, reactive power loss index, voltage change index, and MVA capacity index for each hour are also determined. Finally, the reduction in energy loss with the DG units is determined for each load model for the distribution network.

MULTI-OBJECTIVE INDEX
In most of the analytical approach, a single objective is discussed [2,6,7] with optimal placement and sizing. In this approach, the multi-objective is derived as a single-objective function using the MOI, which is a combination of real and reactive power loss index, voltage change index, and MVA capacity index considering the objective of reduction in power loss, improvement in voltage stability, and decreased line loading. The indices are provided with proper weights: where ∑ 4 i=1 i = 1 and i ∈ [0, 1], i = 1, 2, 3, 4. The impact indices are within the values between 0 and 1, as discussed in [32].
The MOI value is equal to 1 when calculated for a base case system. As discussed in [31,32], the impact indices are given proper weights with the relative importance in DG planning and studies. In this work, the active power loss has a higher weight of 0.4, since it is important in the placement of DG to meet the demand. The MVA capacity index has the second major role and is given a weight of 0.25 as it gives important information regarding the upgradation of conductors in the distribution network. Next, the voltage change index is given a weight of 0.15 and the reactive power loss index receives a weightage of 0.2 each considering its power quality impact. Depending on the importance, the weights of the impact indices can be adjusted. The MOI Equation (51) is minimised subject to the following operational constraints.
Voltage constraint: The operating voltage of all the buses in the distribution network should be within the limit as discussed by the DISCO as V min , V max represent the upper and the lower limits of bus voltage, respectively.
Feeder capacity limit: The apparent power flow through each line in the distribution network should not exceed the maximum capacity of the line: where S i represents the line flow of the ith line and S imax represents the maximum line flow of the ith line.
As the MOI depends on the stochastic nature of photovoltaic and wind generation, the MOI is also stochastic in nature, for each hour. Therefore, the MOI is calculated for the expected generation and load model for any particular period, MOI(t) as given by From this the average MOI can be obtained as where dt is the difference in the time stamp. The AMOI is calculated placing DG unit in each bus, for each hour expected power output. The lowest value of AMOI calculated for each bus implies the best location for DG allocation, thereby the real and reactive power loss is reduced, the voltage profile is improved, and line loading is decreased.
In the analytical approach, the objective is framed with a multi-objective function considering the reduction in power loss, voltage profile improvement, and decreased line loading.

Sizing of the DG unit
To find the optimal size of the DG unit, the framed MOI has to be minimised and equated to 0 by substituting the values of real and reactive power loss index, voltage change index, and MVA capacity index.
To find the minimum MOI, differentiate partially Equation (50) with respect to P DG j and equating to 0,  (28), (32), (38), and (50) with respect to P DG j , the equations can be written as Substituting (59)-(62) in Equation (58): For a particular value of power factor and j , from Equation (27) the reactive power DG can also be calculated. Using Equation (63), the size of the DG unit is found by placing the DG unit in each bus.

Location of distributed generation unit
After calculating, the size of the DG unit in each bus, the expected power output for that particular size of DG unit is found. Now from the MOI, the expected MOI(t) for any specific period has to be found by placing the DG unit on a particular bus. The MOI(t) was calculated for each hour in a day with the combined load and generator model. Finally, the average MOI is calculated from the MOI(t) calculated for each hour on a specific day. The above-mentioned steps have to be done by placing DG in other buses and average MOI in all the buses has to be calculated. The lowest value of AMOI implies the best location to allocate the DG unit.

Method of computation
The method of computation is summarised in the following steps: Step 1: Perform the load flow analysis for the system without incorporating the DG unit.
Step 3: Find the optimal size of the DG unit using Equation (63).
Step 4: Locate the DG unit in a bus and calculate MOI(t) for each hour in a day with the expected power output and load model. Step 5: Find the average value of MOI(t) for a day using the values of MOI(t) for each hour.
Step 6: Repeat the steps (4) and (5) by locating the DG unit in all other buses.
Step 7: Find the minimum value of AMOI and identify the best location to place the DG unit.
Step 8: For the, identified location and the corresponding value of P DG j in the particular bus, find the maximum output power using the relation where P avg DG j is the size of DG unit at the average load level at that bus, CF P DG j is the capacity factor of DG unit, and k DG is the adjustment factor.
Step 9: Find the output of the DG unit for a particular time period as This method of computation was performed for the placement of only one DG unit in a radial distribution network.

Test system
The test system considered for the proposed work was the IEEE 33, IEEE 69 radial bus distribution network, and 110/11 kV distribution substation. Figure 11 shows the IEEE 33bus distribution networks. As discussed in [9], the maximum demand of 33-bus system is 3.7 MW and 2.3 MVAr. The 69bus distribution network is shown in Figure 12. The maximum demand of the 69-bus system is 3.8 MW and 2.7 MVAr. The base values considered for the bus system are 100 MVA and 11 kV. Also, a 16-bus distribution substation located at Irugur, Coimbatore district, Tamil Nadu, has been considered as a test case. The maximum demand of the substation is 36.4 MW. The substation has 16 buses with 12 outgoing feeders and three  Figure 13. To study the impact of time-varying load models, loads assumed are constant, residential, industrial, and commercial load pattern for 24 h in a day is presented in Figure 2. While analysing the impact of residential load model, the load in all the buses are assumed to be residential. Similarly, the analysis is done for the other different load models.

Location of DG at different OPF
The radial bus network is made to run for the load flow analysis with different load models and corresponding MOI is calculated for each bus. The bus with the lowest value of MOI is identified as the best location for DG placement. For a 69-bus system with industrial load for OPF > 0 and OPF < 0 the MOI obtained for each bus is plotted in Figure 14. From the figure it is noted that MOI at bus number 61 is minimum. Thus bus number 61 is selected as the candidate bus for optimal placement of the DG unit for a 69-bus distribution system. Also, the optimal location of 33-bus system is identified as bus number 6. The same procedure is repeated for different load models (commercial, industrial, residential, and mixed). The optimal location of real 16-bus distribution substation is identified as bus number 8. But the optimal location for the different load models may differ depending on the various characteristics of the system and load pattern.

Sizing of DG at different OPF
The optimal size of the DG unit to be located in each bus is calculated using Equation (63) and is plotted in Figure 14 for the 69-bus system for the industrial load model at different optimal power factors. As the location of DG is identified with the minimum MOI, the optimal size is chosen as the DG size corresponding to the minimum MOI location. From the plot, it is noted that 1.85 MW is the optimal size of the DG unit at OPF > 0 that can be placed in bus number 6. Similarly, the average size of the DG unit can be calculated for the different load models and is tabulated in Tables 2, 3, and 5 for 33-and

FIGURE 14
Optimal DG size and MOI values at lagging and leading power factor for the industrial load for 69-bus system 69-bus system, and 16-bus distribution substation, respectively. The part of the optimal size of DG unit is occupied by wind turbine based DG unit and is chosen with respect to the available turbine specification provided in Table A3. The other part is occupied by solar based DG unit because the PV array can be arranged in series and parallel to generate almost the power required depending on the specification provided in Table A1.

FIGURE 15
Expected total power output at bus 61 for 69-bus system with industrial load model at OPF DG Now, for the industrial load model at OPF > 0 and OPF < 0, the total expected power output from the DG unit is plotted in Figure 15. Also, with the average size of the DG unit for different load models, the expected output power with a hybrid DG unit is plotted in Figure 16. The expected power output of DG unit for a 69-bus system, for different load models follow the normalised DG pattern for a day. The MOI value for industrial    Figure 17. It is inferred that the MOI values decrease when there is increased generation. The expected MOI values for different load models are plotted in Figure 18 for the 69-bus system locating a DG unit in the identified bus, respectively. It is understood from the figure that during the night PV generation is zero and the MOI approaches towards unity and during the daytime the MOI is depending on the PV and wind generation.

IMPACT ANALYSIS OF VARIOUS LOAD MODELS
The impact of different load models with the allocation of DG units has been studied by calculating the MVA capacity saving (MVA system ), real power loss (P LDG ), and reactive power loss (Q LDG ) for each hour. Illustration of percentage reduction in power intake from grid with commercial load for 33-bus system For the IEEE 33-bus system, with the installation of DG unit the expected power output, the commercial load pattern for a single day, and the percentage reduction in power intake from the grid is plotted in Figure 19. It is seen that for commercial load after installation of the hybrid DG unit the energy loss reduction for a day is reduced by 58.43%. The energy loss reduction is tabulated in Tables 2, 3, and 5 for the 33-bus system, 69-bus system, and 16-bus real distribution substation, respectively.

Time-varying constant load model
A constant load model, for each hour, consumes constant power from the substation. When a DG unit is placed in the distribution system, the real and reactive power intake is reduced during the period from 9:00 AM to 3:00 PM for 33-, 69-, and 16-bus distribution substation as the solar availability in excess during this period. The real power loss and reactive power loss are also reduced considerably depicted in Figures 20 and 21 for 33-, 69-, and 16-bus distribution system, respectively. Also, there is MVA capacity saving of 2.2, 2.85, and 4.4 MVA at noon for the 33-, 69-, and 16-bus distribution systems, respectively, when the total expected output availability is more as shown in Figure 22. Thus, when a PV and windbased DG unit is installed, the reduction in power loss is more compared to the scenarios, only when a solar and wind-based DG unit is installed because the uncertainty of one renewable source is compensated by other renewable sources.

Time-varying commercial load model
For the commercial load, the integration of solar and windbased DG unit reduces the real and reactive power intake from the grid. The real and reactive power loss reduction is depicted in Figures 20 and 21 for all the test systems. But the reactive power intake is considerably more from DG units compared to other loads. The real and reactive power intake from the substation, and the real and reactive power loss reduction is more significant for commercial load because the commercial load pattern resembles the generation pattern. Also, the MVA capacity saving is significantly large compared to the other load models as shown in Figure 22. This is because the system MVA depends on the active and reactive power intake which wholly depends on the capacity of the DG unit installed in the distribution system. Thus, a larger capacity of DG unit for commercial load model reduces the MVA capacity of substation compared to other load models, with a

Time-varying residential load model
Considering the residential load model, the real power intake is less compared to a constant, commercial, industrial, and mixed load model, however it is more than mixed load model, if DG unit is allocated. Also, the reactive power intake is very much lower than all the load models. The real and reactive power loss reduction is depicted in Figures 20 and 21 for all the test systems, as the real and reactive power loss reduction is less than constant and commercial load but more than industrial and mixed loads. The residential load pattern has its peak demand from 9:00 AM to 3:00 PM and from 6 to 10:00 PM. Thus, the MVA capacity saving of residential load is very less compared to constant, commercial load but there is a marginal difference of MVA capacity savings when compared to mixed and industrial loads as shown in Figure 22.

Time-varying industrial load model
For the industrial load model, the real power intake is very less compared to the constant and commercial load model but higher than the residential load model. The reactive power   intake is also less compared to a constant and commercial load model, however it is more than mixed load model, if the DG unit is allocated. The reduction in reactive power intake is higher compared to residential and mixed load models during the period 11:00 AM to 2:00 PM. This hourly variation of reactive power intake for different load models is due to the irregularity of the hourly demand patterns as shown in Figure 1. The real and reactive power loss and reactive power loss reduction are less than constant, commercial, and residential load models. The industrial load pattern has its peak demand during night-time from 6:00 PM to 12:00 AM. Similarly, the MVA capacity saving is more than residential load model but less than constant and commercial load model as shown in Figure 22.

Time-varying mixed load model
Considering the mixed load model, the real power intake is less compared to a constant, commercial, and industrial load model, however it is more than residential load model if the DG unit is allocated. Also, the reactive power intake is very much lower than all the load models. The real and reactive power loss reduction is depicted in Figures 20 and 21 for all the test systems as the real and reactive power loss reduction is less than constant, commercial, and residential load but more than industrial loads. From Figures 20 and 21, it is demonstrated that for different hours of the day, there exists an irregularity in the reduction pattern of these load models. The variation is because of the avail-ability of expected power output for each hour and hourly load demand pattern. Also, the MVA capacity saving of mixed load is very less compared to constant, commercial load but there is a marginal difference of MVA capacity savings when compared to residential and industrial loads as shown in Figure 22. The real and reactive power loss with different DG units and the performance indices values are tabulated in Table 4 for 33and 69-bus systems, respectively. For the 16-bus substation, all those values are tabulated in Table 5. The indices PLI, QLI, VDI, and IC depend on real power loss, reactive power loss, voltage deviation, and line loading. The performance indices demonstrate that if the values are high, the allocation of DG unit has less impact on the load models, whereas low values of indices demonstrate high impact on the load models. From the average indices plot shown in Figures 23-25, it is observed that PLI and QLI have less value for the commercial load models but industrial and residential loads are more if DG unit is installed because the commercial load pattern ultimately matches the generation pattern owing to the reduction in power loss. And the other reason is industrial and residential load demand is more beyond the daytime. The average voltage deviation index is more for the constant load model and less for the residential load model, as shown in Figures 23-25.
The voltage deviations were considered as an operational constraint as V min ≤ V i ≤ V max and only very small allowable deviations were observed. It is understood that the voltage profile is improved with the DG unit for all the load models. Considering the line loading, the MVA capacity index is almost same for all load models as shown in Figures 23-25. Considering the line loading the apparent power flow through each line in the distribution network should not exceed the maximum capacity of the line. As the apparent power intake from the substation is reduced, MVA capacity index, when the DG unit is allocated, has a value lesser than unity. From this index, it is noted that if the value approaches towards unity, the system line has to be upgraded. Tables 2 and 3 show the summary and comparison of energy loss reduction for both 33-and 69-bus system when a DG unit is allocated in the distribution network. The energy loss for 16-bus substation is tabulated in Table 5. The energy loss is found by calculating the hourly power loss over a day. It is understood that the energy loss reduction is more when a hybrid DG unit is installed because of the availability of both sources throughout the day. Comparing the load models, the energy loss reduction is more for the commercial load than the residential loads and industrial as the hybrid DG unit generation matches the commercial load pattern.

DG penetration level
The penetration of the hybrid DG unit is shown in Figure 26 for both 33-and 69-bus systems. It is understood that the penetration level is 52.06% for commercial load and comparatively the penetration of DG unit is less for industrial load (45.19%) and residential load (45.83%) for 33-bus system. Similarly, for the 69-bus system, the penetration level for the commercial load is 41.96% which is more than the penetration level for the industrial load (42.03%) and residential load (42.49%). For the 16-bus substation the penetration level for the commercial load is (22.81%) which is more than the penetration level for the industrial load is (24.12%) and residential load is (22.33%). The reason behind the high penetration of the DG unit for the commercial load is the generation pattern matches the commercial load pattern. It is also observed that behaviour of the system load plays a major role in DG penetration. Considering the two IEEE distribution networks it is observed that the DG penetration is lower for 69-bus system compared to a 33-bus system as shown in Figure 26.

CONCLUSION
The impact analysis of different voltage-dependent time-varying load models on a single hybrid DG planning in an IEEE 33bus, 69-bus, and a real 16-bus distribution substation is investigated in this approach. Here an MOI based analytical expression is derived with the real power loss, reactive power loss, voltage deviation and MVA capacity to reduce the losses, improving voltage profile and decreased line loading. With the analytical expression derived, the size of the DG unit to be allocated was found. The MOI is minimised to locate the optimal place to install the DG unit. The hybrid DG unit is allocated on the identified bus and load flow analysis is carried out, and the reduction in real and reactive power loss, improvement in voltage profile, and line capacity are observed from the performance indices for different load models. Finally, energy loss reduction with DG unit allocation is also studied. Considering each load models, the percentage reduction of energy loss is more for the commercial load model because the commercial load pattern resembles the generation pattern. In future work, the analytical approach can be extended for optimal placement and sizing of multiple DG units. Further, off-the-shelf optimisation algorithms can be applied to study the practical scenarios.

Nomenclature
AMOI average multi-objective index CDF cumulative distribution function DG distributed generation DISCO distribution company MOI multi-objective index MVAI mega-volt-ampere capacity index NREL National Renewable Energy Laboratory NSRDB National Solar Radiation Database OPF optimal power factor PLI active power loss index PV photo-voltaic QLI reactive power loss index VDI voltage deviation index v wind speed (m/s) V ci cut-in speed of wind turbine (m/s) V av average speed of wind turbine (m/s) V r rated speed of wind turbine (m/s) V co cut-off speed of wind turbine (m/s) P r rated output power of wind turbine (kW) P D active load power at each bus Q D reactive load power at each bus P i active power flow in each branch Q i reactive power flow in each branch R i resistance of each branch X i reactance of each branch V i voltage magnitude of each branch P DG j active power of DG unit Q DG j reactive power of DG unit , shape parameters of beta distribution function mean of the solar irradiance standard deviation of the solar irradiance