Load area aggregation considering integration of electric vehicles to the system

Current electric power systems have an increasing penetration of electric vehicles, and its effect has to be considered in different studies, such as optimal dispatch or voltage stability, among others. Additionally, considering that power system analysis becomes complex when the number of buses increase, this paper presents a methodology for aggregation of load areas that use a measurementbased load modeling approach based on an evolutionary computational technique and a classical reduction method. This aggregate load area model is proposed to reduce areas that consider electric vehicle (EV) load models. The proposed method provides a static equivalent load model and an equivalent network that can be used to reduce the computational effort required by power system studies. In order to validate the application of the proposed methodology, a 30-bus power system considering several disturbances and levels of penetration of the electric vehicles was used. The results show that the equivalent network model allows the reproduction of different events with an acceptable accuracy when it is compared to the original system behavior.


Introduction
Nowadays, power system studies may become complex as new elements are increasingly connected to the network.One of these special elements are the electric vehicles (EV), which pose new challenges and have become an important area of research.These load elements presented a growing interest in the last years, as shown in Figure 1, due to their environmental and economic benefits, compared to conventional vehicles (Global EV Outlook 2015Update, 2015).
Several studies on grid impacts associated with EV charging were proposed; steady-state voltage stability has been analyzed in (Dharmakeerthi, Mithulananthan & Saha, 2014) (Mitra & Venayagamoortghy, 2010), small signal stability in (Das & Aliprantis, 2008), influence in grid load capacity in (Hadley, 2006) and power quality in (Richardson, Flynn & Keane, 2010) (Garcia & Pecas, 2013), among others.However, as mentioned above, the presence of these loads in studies may imply exigent computational analysis efforts, and therefore the incorporation of EV load model in network reduction methodologies is required, since it was not previously analyzed.Currently, power grids are huge and highly interconnected systems, with thousands of elements that present complex interactions among them.Generally, when a specific study is performed in these systems, analysts are interested in the behavior of a specific part of the system (internal), and therefore it is required to employ tools that allow to reduce the rest of the system (external) that does not require a detailed description.The objective of the network reduction is to obtain an equivalent system (static or dynamic) at the interface to the non-reduced part of the grid (internal system), which presents a similar behavior compared to the original system.This allows an improvement of the computational time and reduces the resources used in several analysis tasks (Sturk, Vanfretti, Milano & Sandberg, 2012).
During the last decades, different techniques on power system reduction were proposed.These methods are suitable for stability analysis, steady state power flow for security assessment and control applications.Generator aggregation based in coherency approach (Chow, Galarza, Accari & Price, 1995) (Milano & Srivastava, 2009), modal analysis (Chaniotis & Pai, 2005) (Ghosh & Senroy, 2012) and identification approach (Ramirez, 1999) were widely used, in contrast to the equivalents of loads that was less analyzed.This has been sustained by the uncertainty in the composition and diverse characteristics of loads that do not allow to define an adequate and accurate model.Additionally, the integration of new load models in the power systems, such as electric vehicles (EV) increases such uncertainties.
Related to load area reduction, the following approaches were proposed: In (Wen, Jiang, Wu & Cheng, 2003) (Wen, Wu, Nutall, Shimmin & Cheng 2003) the authors proposed an aggregated load area model (ALAM).A genetic algorithm was used to find the parameters of an equivalent static load model and ficticious transmission lines for power transient stability analysis.A ZIP load model with frequency dependent terms was used to replace the external area in a 5-bus power system.It was demonstrated that the aggregated load area model can accurately capture the behavior of the load area in both steady-state and transient analysis; however, only static loads were considered in these two references.Later, the same authors, Wen, Wu et al, in (Wen et al., 2004) proposed the incorporation of an aggregate exponential recovery model in order to achieve the dynamic characteristics of the ALAM for long-term voltage stability analysis.Another evolutionary computation technique such as Particle Swarm Optimization was used to obtain an equivalent load area model.For instance, in (Wei, Wang, Wu & Lu, 2005) the parameters of fictitious line (R, X), line shunt admittance together with six coefficients of the ZIP load model were identified using one step identification scheme compared with the two-step scheme used in the previous references (Wen, Jiang, Wu & Cheng, 2003) (Wen, Wu, Nutall, Shimmin & Cheng 2003).
The methodologies above exposed also require parameter estimation of equivalent load and parameters of lines that interconnect it to the boundary nodes.This can represent a major computational effort, due to the several parameters that must be estimated.In this paper, the authors propose the reduction of the effort of the estimation algorithm, using, first, a classical reduction methods such as static Ward reduction for fictitious line parameter determination, and, next, only an optimization technique for the estimation of the equivalent load parameter of the load area that includes EVs.
The structure of this paper is organized as follows.Section II presents the theoretical background related to EV applied load models and the network reduction formulation.In Section III the proposed methodology is exposed.Based on that, the simulation results are shown when the methodology is tested on a modified IEEE 30 bus system in Section IV.Finally, in Section V the main conclusions obtained from the result analysis and the future work tasks are presented.

Electric vehicle load model
According to the previous scientific related documents review, the constant power load model has been widely used to represent the EV load in distribution systems, but it was found that this assumption might not be accurate (Haidar & Muttaqi, 2013).Recently in (Dharmakeerthi et al., 2014), the authors propose a combination of constant power and an exponential load model with negative a for static voltage stability purposes.This negative characteristic is due to the influence of the parasitic resistance of the input filter inductor.Equation (1) represents this formulation.
Where parameters a, b and α depend on resistance values of charger configuration and P o represents the power consumption of the charger at the reference voltage (V o ).This load model was used in this paper.

Network reduction
For power system reduction purposes, it is necessary to divide the power system into three parts: an external area, which is susceptible to being replaced by an equivalent representation (for this paper, it corresponds to the loading area); an internal area, which is the rest of the system under study and remains unaltered during the reduction process; and a set of boundary buses which links these areas.
The existing equivalent methods including the Thevenin equivalent method, the traditional Ward equivalent method (Ward, 1974) and the extended Ward equivalent method (Monticelli et al.,1979) were widely used.For this paper, the classical Ward equivalent was used.The Ward equivalent methodology replaces the external system for a set of transmission lines (that interconnect boundary nodes Y mn ) and a set of power injections and shunt admittances Y 0i connected to boundary nodes.Figure 2 shows the configuration of the system.Initially, consider the full system nodal admittance matrix Y bus (2), that can be rewritten in terms of subscript e and is related to a set of nodes of external subsystem; and subscript k, which corresponds to set of boundary and internal nodes subsystem.
Where I bus is the nodal injection current vector, V bus is the nodal voltage vector.In order to obtain the network reduction with the Ward approach (Deckmann et al, 1980), the Gaussian elimination method is applied to Y bus , as shown in Equation (3).
Where I e×e is an upper triangular matrix of size e×e and Y net is the reduced system´s nodal admittance matrix according to Equation (4).
Subscript f is related to the boundary nodes set and subscript i is related to the internal nodes set.Matrix Y net includes equivalent shunt admittances Y 0j and equivalent branches connecting boundary-boundary buses Y mn of Figure 2, according to Equations ( 5) and ( 6).
Shunt admittances connected to each boundary bus Y 0j are calculated by subtracting the shunt capacitance of lines of the internal system to the value obtained from Equation ( 6).
Finally, the Ward equivalent method considers equivalent power injections at boundary buses (Monticelli et al., 1979), as shown in Equation ( 7).
Where P j LOAD , Q j LOAD are the real and reactive power of load elements connected in boundary node j respectively; P j GEN , Q j GEN are the real and reactive power injections in boundary j, respectively; G jp and B jp are the real and imaginary parts of the j-p element pertaining to the reduced admittance matrix Y net ; V j , θ j correspond to the voltage magnitude and the voltage angle value of node j, respectively.

Methodology
The purpose behind the network reduction is to develop an aggregate load model to facilitate static analysis of the power system when penetration of EV loads is considered.
The proposed approach for load area reduction is based on strategies that divide the power system under study in three subsystems (area to be reduced or external subsystem, the boundary buses and the internal subsystem) as depicted in Figure 3.For this study, load area is considered as the external subsystem (including EV load models and ZIP load models connected to external nodes), which is finally represented by means of equivalent ZIP load models connected to the boundary buses, whose parameters are calculated using a measurement-based load modeling process.
Only the impedances that interconnect boundary buses are obtained from a Ward equivalent method.This allows the reduction of the effort of the estimation algorithm and improves the aggregation load process.
The general procedure is depicted in Figure 3, where initially a Ward equivalent of the external system is calculated; this area is supposed to contain only loads and EV models.Subsequently, measurements of active and reactive power flowing in the external area are obtained, in order to determine an equivalent representation of power consumed by the load area; this procedure is achieved using Particle Swarm Optimization (PSO) for the minimization of the differences between the measured data and model simulated data.

Load model parameterization using measurements
A measurement-based load modeling procedure is applied to obtain a representation by means of a load model of the power supplied to the load area.This procedure requires of the initial selection of a model structure, which is expected to accurately represent the response of elements aggregated by the load model; this model has a defined structure of unknown parameters.Once a model structure is chosen, a fitting process is used to determine a set of these parameters that best fit the model response to the measurements acquired from the power system.

Model structure definition
Equivalent loads connected to each boundary bus are represented using a ZIP load model as shown in Figure 3.This structure is defined according to Equations ( 8) and ( 9): where P, Q are the active and reactive power consumed by the load, P 0 , Q 0 are the initial active and power the load (obtained from a power flow a 0 , a 1 , a 2 are coefficients related to constant power, constant current and constant impedance components of active load, and b 0 , b 1 , b 2 are coefficients related to constant power, constant current and constant impedance components of the reactive load, respectively.

Estimation of load parameters
The estimation process can be formulated as an optimization process, where a set of load model parameters is determined in order to minimize the differences between measurements and load response.The optimization process applied for this study is Particle Swarm Optimization (Del Valle et al., 2008) (Hong, 2009) (Yoshida et al., 2000).PSO is a population-based method where solutions, called particles, fly around the search space according to its own experience and swarm experience.Each particle p is composed of four elements, as indicated in Equation ( 10).
With the remaining parameters, a 1 , b 1 are calculated according to Equations ( 8) and (9).A modified particle swarm optimization method that uses a mechanism of speed restart is applied for estimation (Rodríguez, Pérez and Mora, 2013).

Test results
In this section, the proposed aggregation load model approach is validated in a modified IEEE 30-bus test system (PSTCA, 2015) that considers EV load model incorporation.
The power system is divided into the following external, boundary and internal networks in order to reduce the loading area.The proposed division is shown in Table 1 and is depicted in Figure 4.This division is proposed to define an area of a high concentration of loads, which is separated from the main sources of power generation.Loads in the internal system are represented using ZIP load models.EV load models exposed in Equation ( 1) are also considered for loads connected to buses 18, 21 and 24, equivalent to 40 %, 20 % and 30 % of total active load connected to the bus.The parameters of EV load models are selected according to (Dharmakeerthi, Mithulananthan and Saha, 2014).The selected parameters consider the worst-case scenario for voltage stability, which corresponds to the most negative value for the exponent of EV load model.

Validation of the proposed methodology
In order to test the capabilities of the aggregation strategy, a cross-validation-based method is proposed.A set of three disturbances is simulated and measurements of voltage magnitude at boundary buses, active power and reactive power flowing into the external area are also stored.The considered disturbances include outage of line connected between buses 1 and 2 (dataset 1), three-phase fault at bus 7, fault resistance of 0.1 p.u. and fault time of 200 ms (dataset 2) and three-phase fault at bus 30, fault resistance of 0.1 p.u. and fault time of 250 ms (dataset 3).
These disturbances are simulated to consider different levels of voltage variation, which have a direct effect on the load model ability to represent the aggregated load response.To illustrate the severity of the disturbances in the power system, Figure 6 shows the voltage variation at bus 10 during these disturbances.In this case, the second has the greatest voltage variation of the disturbances, which may have implications on load area equivalent generalization.
The equivalent system obtained for dataset 1 is shown in Figure 6, where three equivalent transmission lines are included, L10-12, L12-24 and L10-24, and three equivalent ZIP load models at the boundary buses.The results obtained with the network reduction methodology are shown in Table 2.    Figures 7 and 8 show the results for active and reactive power of the power injection at bus 10; these results include the validation of the equivalent response to the same disturbance used for estimation, indicated as Disturbance 1.However, it is important to analyze if this model is also adequate to represent the power injection under different contingencies; in this case, a cross validation method is useful to analyze the generalization capability of the model.In general terms, if a set of N measurements are stored, a load equivalent is obtained using one of these datasets and validated with the remaining datasets.An error indicator is determined for each validation case, as for active and reactive power.An average relative error is calculated according to (11): In relation to the generalization capabilities of the equivalent load model presented in Table 2, the cross validation process was applied.With exception of load area active power response to disturbance 2, measured data from power system and simulated data from load area equivalent are overlapped in Figures 7 and 8.The differences observed for the validation process considering disturbance 2 are due to the fact that load equi-valent was obtained with data of a voltage disturbance from 1.03 p.u. to 0.99 p.u.; however, it was tested to represent a disturbance with a voltage variation from 1.03 to 0.68.
Load equivalents are more accurate when these are tested with disturbances that have a similar voltage variation compared to the disturbance used for estimation.Tables 4 and 5 summarize cross-validation results for active and reactive power, respectively.
For all models and buses, estimation errors are less than 1 %.The load area equivalent obtained using the Dataset 3 (Model 3) has considerable higher errors when it is validated using Dataset 1 and Dataset 2. This is due to the fact that Model 3 has a better response when disturbances have a similar response to the measurements used for obtaining the load model.However, it is evident that estimation errors are still low and load equivalents are valid under these circumstances.

Equivalent validation considering changes in penetration levels of EV load
In order to verify the equivalent response after changes in load composition, loads of EV load models are increased to 1.5 times its rated value.Equivalent load area is obtained using Dataset 2 (Model 2).Initial power of ZIP loads needs to be recalculated using load flow information; however, ZIP load model parameters are not modified.Fault event of dataset 1 is simulated for the modified system.Active and reactive power results for each bus are depicted in Figures 9-11.Despite variation in load composition due to the increase of EV active load, the load area equivalent response is adequate for the representation of the load area.

Conclusions
In this paper, a methodology based on measurement based load modeling approach and the classical reduction method for aggregation of load area model is proposed.The equivalent ZIP load model and the fictitious lines and shunt admittances obtained, can substitute the external system that considers EV load models, preserving the internal system behavior.The equivalent model system parameters can be obtained using an evolutionary computational technique based on system measurements and the Ward equivalent method.Satisfactory results were obtained when the equivalent load model had been tested under several perturbation and different penetration levels of EV load.
This methodology can be easily implemented when measurements at boundary nodes are available.This will allow a continuous updating process of the equivalent load model.
Due to the fact that large scale deployment of EVs could bring major challenges to power system analysis, the methodology proposed in this paper can provide an alternative tool to reduce the associated efforts.Future work will involve other equivalent load models when a large number of EVs in the power grid are considered.

Figure 5 .
Figure 5. Voltage variation at bus 10 for the simulated

Figure 6 .
Figure 6.Equivalent power system for modified IEEE 30-bus system.
Figures 7 and 8 also show the results for active and reactive power response of load area equivalent obtained from Dataset 1 at bus 10, when it is validated using information of Datasets 1 to 3. Dotted lines indicate Measured data obtained from datasets (M), while solid lines indicate Simulated data obtained from load area equivalent response (S).

Figure 7 .
Figure 7. Active power response of load area equivalent under three different scenarios.

Figure 8 .
Figure 8. Active power response of load area equivalent under three different scenarios.

Figure 9 .
Figure 9. Active and reactive power response of power flowing through bus 10.

Figure 10 .
Figure 10.Active and reactive power response of power flowing through bus 12.

Figure 11 .
Figure 11.Active and reactive power response of power flowing through bus 24.

Table 1 .
Power system division for network reduction.

Table 2 .
Equivalent line parameters and shunt admittances for dataset 1.

Table 3
summarizes the estimated parameters for equivalent ZIP loads connected to boundary buses.

Table 3 .
Estimated ZIP parameters for equivalent loads.