Integrating Electric Vehicles into Power System Operation Production Cost Models

:


Introduction
Currently, electric vehicles represent a small fraction of the transportation sector. Nevertheless, different market analysts predict a high penetration of electric cars, which will increase up to 33% by 2040 [1] and 50% by 2050 [2]. The advancement propelling the increase in the penetration of electric vehicles is the battery, via both an increase in the energy density (kWh/kg) and a reduction in the cost per energy (USD/kWh). In 2018, costs were in the range of 209 USD/kWh [3]. In 2020, the battery pack achieved a price of 137 USD/kWh [4]. The International Renewable Energy Association (IRENA) anticipates a cost of 150 USD/kWh in the 2020s, which will make EVs a competitive transportation option [5].

•
The paper focuses on the bulk power system rather than the distribution side, which usually involves AC power flow simulations. • This paper is the first to describe the incorporation of EV modeling of coordinated charging directly into a production cost model (PCM). Uncoordinated and coordinated EV charging models have been devised to bookend the positive and negative impacts of charging strategies. • This work successfully demonstrates the implementation of both uncoordinated and coordinated charging models in a PCM for a simple power system.

Batteries and Their Effect on EV Costs
A reduction in the cost of batteries for electric vehicles decreases the cost of new electric cars. Battery packs account for about half the cost of modern EVs. Recently, the cost of batteries for EVs has been decreasing year by year, reducing 73% from 2010 to 2016 [5]. The price reduction from 2010 to 2020 is estimated at 89% [4]. The price drivers for EV lithium-ion (Li-ion) batteries are the chemistry, manufacturing capacity, battery size, and charging speed [13]. Battery cells for fast charging are more complex and thus more expensive to produce. As there has been an increase in the production and size of Li-ion batteries, the manufacturing cost has decreased. The chemistry of the battery is one key to their price reduction. Reducing the percentage of cobalt in the battery cathode reduces their cost.
The reduction in the price of EVs is also due to the cost of batteries. In this regard, the Vehicle and Technology Office at the U.S. Department of Energy (DOE) set a goal to reach 100 USD/kWh by 2033, and 80 USD/kWh by 2038 [30]. Regarding the capital cost of new vehicles, the U.S. National Renewable Energy Laboratory (NREL) projects substantial, yet achievable, reductions in battery costs that would lead to the capital cost of EVs matching the capital cost of equivalent conventional vehicles before the year 2030 [13].
It is essential to consider the operating cost of new vehicles, as this compensates for a higher capital cost. The Electric Power Research Institute (EPRI) showed that the total cost of ownership is lower for EVs than for conventional vehicles and hybrids [31]. In this comparison, EPRI used 2014 cost data (manufacturer, taxes, credits, destination charges, and EV chargers for Palo Alto, CA) and assumed 12 years or 150,000 miles of ownership. Concerning the lower cost of driving EVs, NREL anticipates that the levelized cost of driving (LCOD) [30] (without tax credits or other incentives) will be in parity with conventional vehicles as soon as 2021. For example, using the vehicle cost calculator from the DOE's alternative fuel data center, the Tesla Model 3 Standard Range has an equivalent cost of driving to the Toyota Camry XLE after eight years of ownership [32]. With cost parity between EVs and standard vehicles, there is a high chance that EVs will become a substantial portion of the car fleet in the years to come; as much as 75% of the sales of light-duty vehicles will be electric in 2050 [33].

Production Cost Model
A power system PCM is a software tool that mimics the operation of the power system by minimizing operation costs while adhering to operational constraints, ultimately yielding operational costs and locational marginal prices (LMPs), as well as other useful information. PCM software solves the unit commitment (UC) and economic dispatch (ED) optimization problems to minimize power system operational costs. ED is an optimization problem in which the objective function is to minimize the total operation cost of the system from a predetermined set of online resources to meet the load in the system. ED considers min/max power dispatch, transmission, ramping, and fuel constraints. Linear programming optimization is a common method to solve this type of problem. UC is also an optimization problem, where the goal is to determine the electrical generators that must be available for use, i.e., be committed, to meet the forecasted power demand at the minimum operation cost in a system. Mathematically, the state of a particular generation unit is either off (not committed) or on (committed) in a time step (i.e., the decision variable is binary). A solution technique appropriate for solving both UC and ED is mixed-integer linear programming (MILP), which adds integer multipliers in the objective function to account for the state of each generator as on or off. The MILP formulation of UC includes more integer variables (constraints) than ED, also called "security-constrained unit commitment" (SCUC) [34]. For example, it includes start-up costs, start-up time, minimum up/downtime, reserve requirements, transmission, ramping, and fuel constraints. In this research, statistical models of EV charging appropriate for incorporation into a PCM are presented. The specific PCM employed to demonstrate the method is Polaris Systems Optimization's power system optimizer (PSO), which uses MILP to solve the optimization problems [35].

Methods
This section presents statistical models for EV charging that are adaptable to any regional mobility pattern, based on data from the EV Project. The EV Project data came from over 10,000 public and residential charging stations over three years (2011 to 2013), covering 4 million charging events in 18 different metropolitan areas across the U.S., with more than 4200 Nissan Leaf (24 kWh) and 1800 Chevrolet Volt (16 kWh) vehicles. In addition to describing the charging characteristics, a model of a simple 5-bus grid for use in testing the statistical models, the "NAU 5-bus system", is also presented. Further details on the methods described below are reported in two M.S. theses [36,37].

Battery and Charging Level
This study assumed a battery capacity of 24 kWh (Nissan Leaf 2013) as an average representative of all EVs. For statistical modeling, the EV and its charging characteristics are arbitrary. Regarding the type of charger, a survey study performed in the Salt River Project service area (in Arizona, U.S.) by the EPRI showed that 66.7% of charging events occur using Level 2 chargers [38]. Level 1 and DC fast-charging event occurrences were 32.5% and 0.9%, respectively. Therefore, we assumed that the EV fleet will employ the percentage of each charger portrayed in Figure 1. For the statistical models being created, the battery capacity and distribution of charger types were definable variables (i.e., they are not constrained to the numbers mentioned above).
World Electr. Veh. J. 2021, 12, x FOR PEER REVIEW 4 of 23 more than 4200 Nissan Leaf (24 kWh) and 1800 Chevrolet Volt (16 kWh) vehicles. In addition to describing the charging characteristics, a model of a simple 5-bus grid for use in testing the statistical models, the "NAU 5-bus system", is also presented. Further details on the methods described below are reported in two M.S. theses [36,37].

Battery and Charging Level
This study assumed a battery capacity of 24 kWh (Nissan Leaf 2013) as an average representative of all EVs. For statistical modeling, the EV and its charging characteristics are arbitrary. Regarding the type of charger, a survey study performed in the Salt River Project service area (in Arizona, U.S.) by the EPRI showed that 66.7% of charging events occur using Level 2 chargers [38]. Level 1 and DC fast-charging event occurrences were 32.5% and 0.9%, respectively. Therefore, we assumed that the EV fleet will employ the percentage of each charger portrayed in Figure 1. For the statistical models being created, the battery capacity and distribution of charger types were definable variables (i.e., they are not constrained to the numbers mentioned above).

Battery State of Charge
This research assumed that the charging behavior of the EV owners in the EV Project data can be modeled statistically and scaled to be representative of any number of EV owners. Four statistical distributions of the state-of-charge (SoC) data were tested: normal (Gaussian) distribution, binomial distribution, Poisson distribution, and beta distribution. The chi-square goodness-of-fit test was applied to each distribution for the fitted data. For both the initial and final SoC data, the best fit was found using the beta distribution. At a 95% confidence in the differences between the observed value and the predicted value from a probability distribution, the best fit is the distribution with the smallest chi-square test value. For the initial SoC, the beta distribution with a chi-square test value of 357.76 ranked the lowest.
The beta distribution is a general type of continuous probability distribution. This distribution is defined on the interval (0, 1) parametrized by two positive shape parameters and . The probability density function of a random variable "x" that follows the beta distribution is determined as follows [39]: where B is the beta function, a normalizing constant to keep the total probability equal to 1.

Battery State of Charge
This research assumed that the charging behavior of the EV owners in the EV Project data can be modeled statistically and scaled to be representative of any number of EV owners. Four statistical distributions of the state-of-charge (SoC) data were tested: normal (Gaussian) distribution, binomial distribution, Poisson distribution, and beta distribution. The chi-square goodness-of-fit test was applied to each distribution for the fitted data. For both the initial and final SoC data, the best fit was found using the beta distribution. At a 95% confidence in the differences between the observed value and the predicted value from a probability distribution, the best fit is the distribution with the smallest chi-square test value. For the initial SoC, the beta distribution with a chi-square test value of 357.76 ranked the lowest.
The beta distribution is a general type of continuous probability distribution. This distribution is defined on the interval (0, 1) parametrized by two positive shape parameters α and β. The probability density function of a random variable "x" that follows the beta distribution is determined as follows [39]: where B is the beta function, a normalizing constant to keep the total probability equal to 1. A histogram of the initial SoC data from the EV Project and beta distribution is shown in Figure 2. The least-squares method for the beta distribution took values of α = 3.28 and β = 3.27. A histogram of the initial SoC data from the EV Project and beta distribution is shown in Figure 2. The least-squares method for the beta distribution took values of = 3.28 and = 3.27. A histogram showing EV battery state of charge at the end of charging events (i.e., target state of charge) taken from the EV Project [40], along with a fitted beta distribution, is shown in Figure 3. Note the two peaks in the data, one at 70-80% SoC and the other at 90-100% SoC. The double peak occurs because 80% or lower was the target SoC. A faster charging rate happens below 80% for lithium-ion batteries. However, more than twothirds of the charging ended with an SoC above 80%, and more than half above 90%. Beta distributions were once again found to fit the data best, and in this case, two beta distributions were used, one for the range 0-80% of SoC ( = 5.51 and = 0.16) and another for the range 80-100% of SoC ( = 14.59 and = 0.67). Again, at 95% confidence, the lowest difference from the observed and predicted values for the final SoC came from using two beta distributions with a chi-square test value of 371.   A histogram showing EV battery state of charge at the end of charging events (i.e., target state of charge) taken from the EV Project [40], along with a fitted beta distribution, is shown in Figure 3. Note the two peaks in the data, one at 70-80% SoC and the other at 90-100% SoC. The double peak occurs because 80% or lower was the target SoC. A faster charging rate happens below 80% for lithium-ion batteries. However, more than two-thirds of the charging ended with an SoC above 80%, and more than half above 90%. Beta distributions were once again found to fit the data best, and in this case, two beta distributions were used, one for the range 0-80% of SoC (α = 5.51 and β = 0.16) and another for the range 80-100% of SoC (α = 14.59 and β = 0.67). Again, at 95% confidence, the lowest difference from the observed and predicted values for the final SoC came from using two beta distributions with a chi-square test value of 371.
A histogram of the initial SoC data from the EV Project and beta distribution is shown in Figure 2. The least-squares method for the beta distribution took values of = 3.28 and = 3.27. A histogram showing EV battery state of charge at the end of charging events (i.e., target state of charge) taken from the EV Project [40], along with a fitted beta distribution, is shown in Figure 3. Note the two peaks in the data, one at 70-80% SoC and the other at 90-100% SoC. The double peak occurs because 80% or lower was the target SoC. A faster charging rate happens below 80% for lithium-ion batteries. However, more than twothirds of the charging ended with an SoC above 80%, and more than half above 90%. Beta distributions were once again found to fit the data best, and in this case, two beta distributions were used, one for the range 0-80% of SoC ( = 5.51 and = 0.16) and another for the range 80-100% of SoC ( = 14.59 and = 0.67). Again, at 95% confidence, the lowest difference from the observed and predicted values for the final SoC came from using two beta distributions with a chi-square test value of 371.

Availability-Arrival and Departure Time
Charging "availability" is defined here as the period when the charger for an EV is connected to a vehicle and is able to charge. The EV Project data include charging availability at home, at work, and at public locations. A chart showing the charging availability data for each hour of the day on an average weekday is provided in Figure 4 [41]. Of the charging events, 85% occurred at home, with just above 40% availability during the night-time hours, and between 10% and 15% during daytime working hours. As the distribution of charging availability data can vary significantly between different regions of a country, or between different countries, it was not fitted to a statistical distribution but was directly used. In the EV charging models implemented below, a distribution such as this is necessary to determine the charging load profile.

Availability-Arrival and Departure Time
Charging "availability" is defined here as the period when the charger for an EV is connected to a vehicle and is able to charge. The EV Project data include charging availability at home, at work, and at public locations. A chart showing the charging availability data for each hour of the day on an average weekday is provided in Figure 4 [41]. Of the charging events, 85% occurred at home, with just above 40% availability during the nighttime hours, and between 10% and 15% during daytime working hours. As the distribution of charging availability data can vary significantly between different regions of a country, or between different countries, it was not fitted to a statistical distribution but was directly used. In the EV charging models implemented below, a distribution such as this is necessary to determine the charging load profile. In addition to charging availability, the EV Project data included both the arrival and departure times of the EVs. Departure times, which represent an end to an EV's availability for charging, are shown in the histogram in part (a) of Figure 5. Arrival times, the times at which an EV is available for charging, are presented in part (b) of Figure 5. Note the preponderance of departure times in the morning hours and arrival times in the evening hours, with a low and fairly constant number of arrivals and departures during the other hours of the day. The arrival data were used as the start time of charging if there was no TOU price incentive to start charging later in the day. The given distributions for the arrival and departure times were used with a Monte Carlo sampling technique to reproduce the arrival and departure data for a population of vehicles, thus allowing for easy adaptation to other regional or country profiles.  In addition to charging availability, the EV Project data included both the arrival and departure times of the EVs. Departure times, which represent an end to an EV's availability for charging, are shown in the histogram in part (a) of Figure 5. Arrival times, the times at which an EV is available for charging, are presented in part (b) of Figure 5. Note the preponderance of departure times in the morning hours and arrival times in the evening hours, with a low and fairly constant number of arrivals and departures during the other hours of the day. The arrival data were used as the start time of charging if there was no TOU price incentive to start charging later in the day. The given distributions for the arrival and departure times were used with a Monte Carlo sampling technique to reproduce the arrival and departure data for a population of vehicles, thus allowing for easy adaptation to other regional or country profiles.

Modeling the Charging of Electric Vehicles
A survey of modeling techniques used to estimate EV load profiles showed that most researchers implemented bottom-up and agent-based simulations to generate a load profile, and that these were primarily for distribution system load studies [19][20][21][22][23][24][25][26][27]. Bottom-up simulation requires a transition probability at every step, which is unnecessary when creating a load profile for transmission-level studies. Consequently, this research used a Monte Carlo sampling technique, which can generate the load profile by accounting for the charging state of an EV.

Uncoordinated Charging of Electric Vehicles
The uncoordinated charging of EVs occurs when EV owners charge their vehicles according to convenience and without any information about, or communication with, other EV owners or the electric utility supplying the electricity. This section presents a load profile model for reproducing the uncoordinated charging behavior of EV users.
The previous section demonstrated that the various deciding factors of electric vehicle charging may follow some kind of probability distribution. This application of statistical models to address our problem is helpful, as they can be scaled to create load profiles for different numbers of EVs. By definition, Monte Carlo is a method for randomly sampling from a probability distribution to create a random process for a problem [42]. For a chosen number of EVs, a random sampling of the values of the influencing factors (SoC, availability, etc.) for each EV from these distributions can be achieved. In this way, it is possible to produce a load profile from a population of vehicles that adheres to the overall shape of the distribution, scales to the number of vehicles, and gives daily load profiles that differ from day to day.
The sampled values of the influencing factors, as explained above, are processed to obtain the total charging profile of EVs. The flowchart in Figure 6 describes the step-bystep procedure for determining an uncoordinated load profile. The charging profiles for each EV in the fleet are consecutively created in the loop structure. Once all individual EV load profiles are created, they are summed to generate the overall load profile.

Modeling the Charging of Electric Vehicles
A survey of modeling techniques used to estimate EV load profiles showed that most researchers implemented bottom-up and agent-based simulations to generate a load profile, and that these were primarily for distribution system load studies [19][20][21][22][23][24][25][26][27]. Bottom-up simulation requires a transition probability at every step, which is unnecessary when creating a load profile for transmission-level studies. Consequently, this research used a Monte Carlo sampling technique, which can generate the load profile by accounting for the charging state of an EV.

Uncoordinated Charging of Electric Vehicles
The uncoordinated charging of EVs occurs when EV owners charge their vehicles according to convenience and without any information about, or communication with, other EV owners or the electric utility supplying the electricity. This section presents a load profile model for reproducing the uncoordinated charging behavior of EV users.
The previous section demonstrated that the various deciding factors of electric vehicle charging may follow some kind of probability distribution. This application of statistical models to address our problem is helpful, as they can be scaled to create load profiles for different numbers of EVs. By definition, Monte Carlo is a method for randomly sampling from a probability distribution to create a random process for a problem [42]. For a chosen number of EVs, a random sampling of the values of the influencing factors (SoC, availability, etc.) for each EV from these distributions can be achieved. In this way, it is possible to produce a load profile from a population of vehicles that adheres to the overall shape of the distribution, scales to the number of vehicles, and gives daily load profiles that differ from day to day.
The sampled values of the influencing factors, as explained above, are processed to obtain the total charging profile of EVs. The flowchart in Figure 6 describes the step-by-step procedure for determining an uncoordinated load profile. The charging profiles for each EV in the fleet are consecutively created in the loop structure. Once all individual EV load profiles are created, they are summed to generate the overall load profile.

Coordinated Charging of Electric Vehicles
The coordinated charging model of EVs in this research employed an aggregator scheme. An overview of the algorithm used to model the aggregator is described in the flowchart in Figure 7. The first step in coordinated charging is to determine the daily energy requirement for charging the fleet of EVs. The projected maximum power draw at each time of the day comes from the charging availability distribution, the number of EVs (e.g., 1000 or 10,000, etc.), and the charging capacity, as shown in Figure 8. This profile defines a dispatchable EV load during each period of the day. The time in Figure 8 is in hours, but any time increment desired can be created, depending on the time resolution of the source data, which in this case was 15 min. As the size of the EV fleet increases to higher penetrations, EV charging impacts both the total load profile and the correspondent local marginal pricing (LMP) of electricity. Consequently, trying to allocate a coordinated charging profile by using only the day-ahead forecast of LMP with no EVs neglects the influence of EV charging in the LMP profile. Therefore, the coordinated charging profile is determined in the PCM optimization as a dispatchable load. This is necessary so that the charging time series can respond to LMP changes as the load increases or decreases.

Start
Enter the EV fleet total daily energy requirement in MWh.
Enter the time series for the projected maximum power in MW calculated using the charging availability times the number of cars times the average charging power Setup the EVs as a load injector (dispatchable load) in the model of the power system. Solve for UC and ED in a PCM Create EV load profile as part of PCM optimization, subject to constraints identified above The first step in coordinated charging is to determine the daily energy requirement for charging the fleet of EVs. The projected maximum power draw at each time of the day comes from the charging availability distribution, the number of EVs (e.g., 1000 or 10,000, etc.), and the charging capacity, as shown in Figure 8. This profile defines a dispatchable EV load during each period of the day. The time in Figure 8 is in hours, but any time increment desired can be created, depending on the time resolution of the source data, which in this case was 15 min. As the size of the EV fleet increases to higher penetrations, EV charging impacts both the total load profile and the correspondent local marginal pricing (LMP) of electricity. Consequently, trying to allocate a coordinated charging profile by using only the day-ahead forecast of LMP with no EVs neglects the influence of EV charging in the LMP profile. Therefore, the coordinated charging profile is determined in the PCM optimization as a dispatchable load. This is necessary so that the charging time series can respond to LMP changes as the load increases or decreases. The goal of creating the EV load models described above was to use them in studying the effects of EVs on the operation of an actual power system (e.g., the Western Interconnect in the U.S., a regional utility power system model, etc.). However, to show the effec- Hour of a weekday

NAU 5-Bus System
The goal of creating the EV load models described above was to use them in studying the effects of EVs on the operation of an actual power system (e.g., the Western Interconnect in the U.S., a regional utility power system model, etc.). However, to show the effectiveness of the models, it is only necessary to use a simple transmission system model that contains loads, generators, and transmission linkages that are representative of an actual power system. Therefore, a simple 5-bus system, namely the "NAU 5-bus system" (NAU is the acronym for Northern Arizona University), was created, and is depicted in Figure 9. The generator types and capabilities, the renewable energy generation time series, and the load time series are representative of those found in the southwestern U.S.; it includes a subset of the variable renewable energies and generator characteristics from the Reliability Test System of the Grid Modernization Laboratory Consortium (RTS-GMLC) [43,44].

NAU 5-Bus System
The goal of creating the EV load models described above was to use them in studying the effects of EVs on the operation of an actual power system (e.g., the Western Interconnect in the U.S., a regional utility power system model, etc.). However, to show the effectiveness of the models, it is only necessary to use a simple transmission system model that contains loads, generators, and transmission linkages that are representative of an actual power system. Therefore, a simple 5-bus system, namely the "NAU 5-bus system" (NAU is the acronym for Northern Arizona University), was created, and is depicted in Figure  9. The generator types and capabilities, the renewable energy generation time series, and the load time series are representative of those found in the southwestern U.S.; it includes a subset of the variable renewable energies and generator characteristics from the Reliability Test System of the Grid Modernization Laboratory Consortium (RTS-GMLC) [43,44]. The NAU 5-bus system has an installed generation capacity that is representative of a utility in the southwestern U.S. (see Figure 10). The researchers scaled the load from region 1 of the RTS-GMLC down to 1580 MW and applied it to the NAU 5-bus system. This peak load of 1580 MW matches the peak load in the "Pennsylvania-New Jersey-Maryland (PJM) 5-bus system", from which the NAU 5-bus system drew its inspiration [45]. Similarly, VRE sources were also scaled down from region 1 in the RTS-GMLC to the NAU  4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 Figure 9. Northern Arizona University (NAU) 5-bus system layout that portrays buses (bold vertical lines), transmission lines (connecting the buses), loads (inverted triangles), and generators (circles and rectangles).
The NAU 5-bus system has an installed generation capacity that is representative of a utility in the southwestern U.S. (see Figure 10). The researchers scaled the load from region 1 of the RTS-GMLC down to 1580 MW and applied it to the NAU 5-bus system. This peak load of 1580 MW matches the peak load in the "Pennsylvania-New Jersey-Maryland (PJM) 5-bus system", from which the NAU 5-bus system drew its inspiration [45]. Similarly, VRE sources were also scaled down from region 1 in the RTS-GMLC to the NAU 5-bus system. The capacity of the generators used in the NAU 5-bus system is shown in Table 1. The time-series data for load and VRE are the actual real-time (RT) of a 1-h resolution in the RTS-GMLC. 5-bus system. The capacity of the generators used in the NAU 5-bus system is shown in Table 1. The time-series data for load and VRE are the actual real-time (RT) of a 1-hour resolution in the RTS-GMLC. Refer to Table 2 for the load distribution per node in the NAU 5-bus system, which is similar to that of the PJM 5-bus system, except for the additional EV load in node C. A summary of the dispatchable generation fleet of the NAU 5-bus system is presented in Table 3. The prices of various fuels are considered for the year 2025. The transmission parameters for the NAU 5-bus system are the same as those in the PJM 5-bus system [45].  Refer to Table 2 for the load distribution per node in the NAU 5-bus system, which is similar to that of the PJM 5-bus system, except for the additional EV load in node C. A summary of the dispatchable generation fleet of the NAU 5-bus system is presented in Table 3. The prices of various fuels are considered for the year 2025. The transmission parameters for the NAU 5-bus system are the same as those in the PJM 5-bus system [45]. Table 2. Load distribution presented in the NAU 5-bus system.

Production Cost Model
As mentioned, this research aimed to create models of EV load that are suitable for studying the integration of EVs into power system operations via a power system production cost model. Thus, the time resolution in the models is 1 h to replicate the day-ahead electricity market integration of EV charging. The objective of the PCM is to minimize the electricity generation cost, subject to a host of power system constraints. Production costs depend on a variety of factors, including the characteristics of the load and the generation on each node, the capabilities of the interconnecting transmission elements, the reserve requirements, the accuracy of load and VRE forecasts, etc. Consequently, the electricity prices at each node, i.e., the LMPs, often differ from one another. LMPs are an important output from a PCM and are typically higher during the peak load hours and lower during the off-peak hours. Due to the variation in the LMP with the time of day, a price-based charging strategy is a good choice for minimizing operating costs and reducing charging costs.

Power System Optimizer (PSO)
An appropriate PCM tool for this research requires capabilities similar to those typical in power system planning and operations. Polaris Systems Optimization's power system optimizer (PSO) is a tool with the capabilities that are necessary for an accurate simulation of power system operations that mimics the actions and timing of power system planners and operators. PSO uses MILP, has multilevel nested time intervals (overlapping time frames), and simulates uncertainties for both load and VRE resource forecasting [35]. Using multilevel nested time intervals means that the solver can commit units at multiple points in the time ahead of the hour of operation (e.g., commit units with long start-up time days ahead, units with 12-h start-up times a day ahead, units with a 1-h start-up time 1 h ahead, etc.). PSO also employs a "look-ahead" feature that permits commitment decisions to be made in full consideration of forecasted load and VRE beyond the current decision cycle. This multilevel nested time interval capability is important because it emulates the forecast uncertainties that system operators face in committing generators and ancillary services.
For this simulation, the planning cycle for committing units was a day ahead. It had a collection of horizons consisting of 24 1-h periods with a look-ahead of 24 h, and the time increment for dispatch (i.e., the period resolution) was hourly. The look-ahead ensures that the commitment decisions in every horizon are related to the expected conditions in the next horizon. The PCM outputs analyzed to show the effect of the EVs include total production cost (USD), generation (energy produced in megawatt hours; MWh), revenue (USD), and electricity price (LMP, USD/MWh). The revenue (USD) is calculated as the LMP (USD/MWh) times the generation (MWh). Therefore, the EV charging cost is equivalent to the negative of the revenue for the EV load. The model simulated seven days, using RTS-GMLC data from Wednesday 31 January to Tuesday 4 February. Net loads were low at times when the VRE was relatively high. The model run started on Wednesday and finished on Tuesday to ensure that the commitment on Sunday (the day for which analyzed the results) was in accordance with the previous and future power system conditions (i.e., there was sufficient time in advance of Sunday to start the slowest starting coal units). In the NAU 5-bus system on Sunday 2 February, the VRE reached its highest hourly penetration of 18% of the load between 11 a.m. and 1 p.m.
For each time increment of the model data, the "net load" that the conventional, dispatchable generation must serve is found by subtracting the VRE from the load on each bus. During system operation, dispatchable generators are committed and dispatched to serve the net load. The heat curves (heat rates) of the dispatchable generators listed in Table 3 (taken from RTS-GMLC) are multiplied by their fuel costs to calculate the dispatch curves (cost curves), which have units of USD/MWh (refer to Figure A1).

EV Penetration and Modeling
Regarding the quantity of EVs, we tested 1000 EVs. Given the assumed EV characteristics and their availability profile, this represented a maximum possible EV load of 2.8 MW during the morning hours, and about 0.8 MW during the midday hours (refer to Figure 8). This was a small but realistic load penetration compared to the peak load of~550 MW in the results to be presented in the next section. During periods of low load (e.g., night-time), fewer dispatchable generators will be online. Consequently, it is expected that any VRE or dispatchable load will be a more decisive factor influencing the system operation cost.
Within the PSO, the uncoordinated charging model was implemented as an "injector" with a discrete time series of EV load data points. These data resulted from the uncoordinated model previously described. To apply the coordinated charging model, we employed the built-in library "Energy Limited" (EL) in the PSO. The setup for this library requires an estimate of the amount of energy that needs to be served in a decision horizon (i.e., over the course of a day); this energy can be determined from the total daily energy supplied in the uncoordinated charging model, or by using the fleet's daily energy consumption or total charging load in MWh. The dispatch constraint for the EVs is their charging availability, which represents the number of vehicles plugged into the grid; this number caps the power that can be drawn from the power network to charge the EVs in any time period (see Figure 8). The EL library ensures that the total amount of the charging load is served, while honoring the EV availability constraints. For a detailed examination of the input data in PSO format (comma-separated values format), refer to the Supplementary Materials.

Results
This section presents the results of running the PSO production cost model for the NAU 5-bus system for the week of interest. As a representative sample, the results of the simulation are presented for one day, in this case, 22 February, as it had the largest VRE penetration. This section starts by showing the results for the NAU 5-bus system with no EVs, followed by the results when charging 1000 EVs under the uncoordinated and coordinated charging scenarios.

System Operation without EVs
The generation stack and system LMPs for each hour of the day are presented in Figure 11 for a typical day of operation in the period analyzed. The system LMP is represented by the dark dashed line and refers to the scale on the right side of the graph. The LMP varied between 24.6 USD/MWh at low net load, up to 35.4 USD/MWh during the peak load hours. The bars show the generation stack and refer to the scale on the left-hand side of the graph. The nuclear generation is represented by the green color on the bottom of each column and is the base load generation, which does not vary from hour to hour. The next generator stacked on each bar is the combined cycle (CC) gas generation (light blue color). This represents the marginal generator, and its output varies from hour to hour. Between hours 5 and 17, there was substantial VRE generation, represented by the small bars (dark blue, yellow, and red, respectively) stacked on the top of each bar. Note that the LMP increased and decreased as the CC generation changed, according to its cost curve (see Figure A1 frame (e)). The LMP followed the expected behavior and had two peaks, one in the early morning and a higher second peak in the early evening. The LMP valley in the middle of the day was a consequence of high solar and wind power that reduced the net load. During this day, the total operation cost was calculated as 359,060 USD, the daily average system LMP was 29.65 USD, and the peak load was 551 MW at 18 h (6 p.m.).
small bars (dark blue, yellow, and red, respectively) stacked on the top of each bar. Note that the LMP increased and decreased as the CC generation changed, according to its cost curve (see Figure A1 frame (e)). The LMP followed the expected behavior and had two peaks, one in the early morning and a higher second peak in the early evening. The LMP valley in the middle of the day was a consequence of high solar and wind power that reduced the net load. During this day, the total operation cost was calculated as 359,060 USD, the daily average system LMP was 29.65 USD, and the peak load was 551 MW at 18 h (6 p.m.).

System Operation with Uncoordinated EV Charging
Applying the uncoordinated charging model scaled to 1000 EVs resulted in the hourly EV load profile shown in Figure 12

System Operation with Uncoordinated EV Charging
Applying the uncoordinated charging model scaled to 1000 EVs resulted in the hourly EV load profile shown in Figure 12 Running the PCM PSO for the NAU 5-bus system with the addition of the EV load resulted in an operating cost of 359,390 USD (an increase of only 330 USD). In addition, the system LMP did not change; the system dispatch stack and LMP are similar to those shown in Figure 11 because the EV charging profile was small compared with the rest of the injectors.

System Operation with Coordinated EV Charging
Implementing the coordinated charging modeling in the PSO resulted in the daily charging profile shown in Figure 13. Recall that the coordinated charging model used the LMP in combination with EV availability and battery SoC in determining the time to charge as part of the goal to minimize operating costs over the course of a day. Compared with the uncoordinated charging profile (Figure 12), most of the EV charging was moved to the low-LMP hours of the day (hours 1, 2, and 11 to 14). The total daily energy used in charging was the same as in the uncoordinated case (10.75 MWh), but the charging timing was shifted in a way that complied with EV charging availability, battery SoC, and daily Hour of a weekday Running the PCM PSO for the NAU 5-bus system with the addition of the EV load resulted in an operating cost of 359,390 USD (an increase of only 330 USD). In addition, the system LMP did not change; the system dispatch stack and LMP are similar to those shown in Figure 11 because the EV charging profile was small compared with the rest of the injectors.

System Operation with Coordinated EV Charging
Implementing the coordinated charging modeling in the PSO resulted in the daily charging profile shown in Figure 13. Recall that the coordinated charging model used the LMP in combination with EV availability and battery SoC in determining the time to charge as part of the goal to minimize operating costs over the course of a day. Compared with the uncoordinated charging profile (Figure 12), most of the EV charging was moved to the low-LMP hours of the day (hours 1, 2, and 11 to 14). The total daily energy used in charging was the same as in the uncoordinated case (10.75 MWh), but the charging timing was shifted in a way that complied with EV charging availability, battery SoC, and daily charging requirement.
Running the PCM PSO for the NAU 5-bus system with the addition of the EV load resulted in an operating cost of 359,390 USD (an increase of only 330 USD). In addition, the system LMP did not change; the system dispatch stack and LMP are similar to those shown in Figure 11 because the EV charging profile was small compared with the rest of the injectors.

System Operation with Coordinated EV Charging
Implementing the coordinated charging modeling in the PSO resulted in the daily charging profile shown in Figure 13. Recall that the coordinated charging model used the LMP in combination with EV availability and battery SoC in determining the time to charge as part of the goal to minimize operating costs over the course of a day. Compared with the uncoordinated charging profile (Figure 12), most of the EV charging was moved to the low-LMP hours of the day (hours 1, 2, and 11 to 14). The total daily energy used in charging was the same as in the uncoordinated case (10.75 MWh), but the charging timing was shifted in a way that complied with EV charging availability, battery SoC, and daily charging requirement. Hour of a weekday In the case of coordinated charging, Figure 14 shows the generation stack and LMPs. The system operation cost for the day was 359,330 USD. Relative to the plot shown in Figure 11, the system LMP increased between 1 a.m. and 2 a.m., which means that as the marginal generator increased its power dispatch due to the EV load, it increased a step in its cost curve. Figure 15 directly compares this point, plotting the EV charging profiles and the system LMPs. Regardless of this increase in LMP, the overall production cost was less than for the case of uncoordinated charging, as the charging occurred during the lower-cost (i.e., lower LMP) hours during the day. In the case of coordinated charging, Figure 14 shows the generation stack and LMPs. The system operation cost for the day was 359,330 USD. Relative to the plot shown in Figure 11, the system LMP increased between 1 a.m. and 2 a.m., which means that as the marginal generator increased its power dispatch due to the EV load, it increased a step in its cost curve. Figure 15 directly compares this point, plotting the EV charging profiles and the system LMPs. Regardless of this increase in LMP, the overall production cost was less than for the case of uncoordinated charging, as the charging occurred during the lower-cost (i.e., lower LMP) hours during the day.

Comparing the EV Cost of Charging
The cost to charge the EVs was determined by multiplying the magnitude of energy used to charge the EVs each hour by the LMP that hour. The total cost to charge was then

Comparing the EV Cost of Charging
The cost to charge the EVs was determined by multiplying the magnitude of energy used to charge the EVs each hour by the LMP that hour. The total cost to charge was then computed by summing the hourly charging costs over all hours of the day. The uncoordinated charging approach led to a cost to charge of 334 USD, while the coordinated charging cost was 294 USD (see Figure 16). Using the coordinated charging approach reduced the cost of charging by 40 USD, equal to 12%. Interestingly, the system production cost decreased by 60 USD between the uncoordinated and coordinated cases, which was slightly different from the reduction in cost to a charge of 40 USD. This phenomenon is due to the coordinated EV charging profile changing the LMP in the early morning hours. Putting the charging numbers into a per vehicle context, the average EV will charge 10.75 kWh. The average uncoordinated charging cost was close to 0.33 USD per day per vehicle, and the coordinated charging cost was close to 0.29 USD, representing roughly 12% in savings. These numbers are representative of the low energy cost experienced by EV owners. computed by summing the hourly charging costs over all hours of the day. The uncoordinated charging approach led to a cost to charge of 334 USD, while the coordinated charging cost was 294 USD (see Figure 16). Using the coordinated charging approach reduced the cost of charging by 40 USD, equal to 12%. Interestingly, the system production cost decreased by 60 USD between the uncoordinated and coordinated cases, which was slightly different from the reduction in cost to a charge of 40 USD. This phenomenon is due to the coordinated EV charging profile changing the LMP in the early morning hours. Putting the charging numbers into a per vehicle context, the average EV will charge 10.75 kWh. The average uncoordinated charging cost was close to 0.33 USD per day per vehicle, and the coordinated charging cost was close to 0.29 USD, representing roughly 12% in savings. These numbers are representative of the low energy cost experienced by EV owners.  Figure 17 shows a graph of the system and the EV loads for the day of 2 February, the same day as is presented in the previous figures with the dispatch stacks. The system load is shown by the light-blue bars (refers to the left-hand scale), and the EV load in the red (coordinated) and yellow (uncoordinated) bars, both referring to the right-hand scale and having a much lower magnitude than the system load. Compared to the system load alone, the uncoordinated EV charging increased the system's peak load. In our particular simulation, 1000 EVs led to an increase in peak demand of 1 MW. Coordinated EV charging occurred during low-LMP hours and therefore did not increase the peak load.  Figure 17 shows a graph of the system and the EV loads for the day of 2 February, the same day as is presented in the previous figures with the dispatch stacks. The system load is shown by the light-blue bars (refers to the left-hand scale), and the EV load in the red (coordinated) and yellow (uncoordinated) bars, both referring to the right-hand scale and having a much lower magnitude than the system load. Compared to the system load alone, the uncoordinated EV charging increased the system's peak load. In our particular simulation, 1000 EVs led to an increase in peak demand of 1 MW. Coordinated EV charging occurred during low-LMP hours and therefore did not increase the peak load. It is interesting to compare operating costs and LMPs for each of the cases. In Figure  18, the operation cost is shown in the blue columns associated with the primary y axis on the left. In addition, in Figure 18, the average LMP is shown by the orange line related to the secondary y axis on the right. In this figure, the left column is for the system operation with no EV with a cost of 359,060 USD, the central column is for the system operation with uncoordinated EV charging with a cost of 359,390 USD, and the right column is for the system operation with coordinated EV charging with a cost of 359,330 USD. Thus, as shown by comparing the central and right columns with the left one, the addition of the EV load increased the system operation cost. However, the EV coordinated charging profile (right blue column) led to lower operating costs than the uncoordinated charging profile (central column). At the same time, the daily average LMP increased slightly (in 0.47 USD) with coordinated charging (right line). This is not uncommon in power system operation; average LMPs can increase while system operating costs decrease.  4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22   It is interesting to compare operating costs and LMPs for each of the cases. In Figure 18, the operation cost is shown in the blue columns associated with the primary y axis on the left. In addition, in Figure 18, the average LMP is shown by the orange line related to the secondary y axis on the right. In this figure, the left column is for the system operation with no EV with a cost of 359,060 USD, the central column is for the system operation with uncoordinated EV charging with a cost of 359,390 USD, and the right column is for the system operation with coordinated EV charging with a cost of 359,330 USD. Thus, as shown by comparing the central and right columns with the left one, the addition of the EV load increased the system operation cost. However, the EV coordinated charging profile (right blue column) led to lower operating costs than the uncoordinated charging profile (central column). At the same time, the daily average LMP increased slightly (in 0.47 USD) with coordinated charging (right line). This is not uncommon in power system operation; average LMPs can increase while system operating costs decrease. It is interesting to compare operating costs and LMPs for each of the cases. In Figure  18, the operation cost is shown in the blue columns associated with the primary y axis on the left. In addition, in Figure 18, the average LMP is shown by the orange line related to the secondary y axis on the right. In this figure, the left column is for the system operation with no EV with a cost of 359,060 USD, the central column is for the system operation with uncoordinated EV charging with a cost of 359,390 USD, and the right column is for the system operation with coordinated EV charging with a cost of 359,330 USD. Thus, as shown by comparing the central and right columns with the left one, the addition of the EV load increased the system operation cost. However, the EV coordinated charging profile (right blue column) led to lower operating costs than the uncoordinated charging profile (central column). At the same time, the daily average LMP increased slightly (in 0.47 USD) with coordinated charging (right line). This is not uncommon in power system operation; average LMPs can increase while system operating costs decrease.  4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 Figure 18. Summary of the system operation cost and LMP change with the introduction of EV charging, gross comparison against the system with no EVs.

Conclusions
The purpose of this work was to create statistical models for EV charging that were suitable for power system operation. Statistical models for EV charging were created using data from the EV Project. This modeling with the appropriate data resulted in reasonable and scalable EV charging load profiles. The models recognized coordinated and uncoordinated charging loads. A production cost model could integrate the EV charging model and assess the EV effects on the power system cost and loads.
The modeling of uncoordinated and coordinated charging profiles required the battery size and charging level (power kW). In addition, the arrival time and planned energy consumption were necessary (i.e., the initial and final SoC). However, for coordinated charging, both the arrival and departure, or the equivalent availability, plus the energy or charging requirement (MWh) for the whole fleet were required.
Statistical modeling recreated the uncoordinated charging profile. A Monte Carlo simulation technique modeled the uncoordinated charging profile of each vehicle by sampling the initial SoC, final SoC, charging level, and arrival time as the start time of charging. The fleet charging profile was the sum of the individual charging profiles. The statistical analysis chi-square method was used to find the best fit distribution for the initial and final SoC. The best fit for the initial SoC was a beta distribution, and two beta distributions for the final SoC. The least-mean-square method provided the beta parameters to fit the data to beta distributions.
The authors developed an algorithm for modeling the uncoordinated charging and incorporated the result into the PCM model as a time series, similar to the load forecast or VRE forecast. However, coordinated charging modeling must be part of the PCM and modeled as a controllable load with the constraints imposed by the fleet characteristics, energy requirement (MWh), charging availability, and charging level (kW). Otherwise, the modeling will neglect the influence of the EV fleet on the electricity price.
The NAU 5-bus system, a simple 5-bus power system similar to the PJM 5-bus model but with time-series data from the RTS-GMLC model, was created to test the statistical models of the EVs. Polaris Systems Optimization's PSO was the PCM used for this demonstration, and it showed the potential effects of 1000 EVs (maximum demand of~0.2% of the NAU 5-bus peak load) on power system operating costs, LMPs, and dispatch stacks when implementing the two charging strategies mentioned above. The demonstration showed interesting results, and revealed the potential positive or negative effects of EVs on system operation. Coordinated charging exemplified the positive effects of accommodating EV charging coincident with VRE. However, uncoordinated charging demonstrated the negative effect of increasing the system peak load. With scalable statistical models of EVs, such as those presented here, it is possible to investigate the impact of charging strategies for the large-scale adoption of EVs on power system operations and costs. Acknowledgments: Polaris Systems Optimization for providing PSO and support, AIMMS for providing an academic license of AIMMS and CPLEX, and NREL for providing the RTS-GMLC model.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A. Cost Curves NAU 5-Bus System
Dispatch or cost curves of the committable generators presented in the NAU 5-bus system, showing the variable cost of producing a unit of energy with respect to the power output of these generators: Acknowledgments: Polaris Systems Optimization for providing PSO and support, AIMMS for providing an academic license of AIMMS and CPLEX, and NREL for providing the RTS-GMLC model.

Conflicts of Interest:
The authors declare no conflict of interest.

Cost Curves NAU 5-Bus System
Dispatch or cost curves of the committable generators presented in the NAU 5-bus system, showing the variable cost of producing a unit of energy with respect to the power output of these generators: