Engineering, Technology & Applied Science Research Vol. 8, No. 3, 2018, 2901-2906 2901 www.etasr.com Le and Nguyen: Stochastic Sizing of Energy Storage Systems for Wind Integration Stochastic Sizing of Energy Storage Systems for Wind Integration Dinh Duong Le Faculty of Electrical Engineering The University of Danang—University of Science and Technology Danang, Vietnam ldduong@dut.udn.vn Nhi Thi Ai Nguyen Faculty of Electrical Engineering The University of Danang—University of Science and Technology Danang, Vietnam ntanhi@dut.udn.vn Abstract—In this paper, we present an optimal capacity decision model for energy storage systems (ESSs) in combined operation with wind energy in power systems. We use a two-stage stochastic programming approach to take into account both wind and load uncertainties. The planning problem is formulated as an AC optimal power flow (OPF) model with the objective of minimizing ESS installation cost and system operation cost. Stochastic wind and load inputs for the model are generated from historical data using clustering technique. The model is tested on the IEEE 39- bus system. Keywords-energy storage system; ESS; OPF; sizing; stochastic; wind I. INTRODUCTION Renewable energy has been increasingly integrated into power systems as a result of the effort to reduce CO2 emissions and build a future power grid economically feasible and environmentally sustainable. Particularly, according to the Blue Map scenario for power supply, electricity generation from renewable energy will provide a share of 22% of global electricity generation in 2050, which grows almost threefold compared to the baseline scenario [1]. Along with this growing share of renewable technologies, greater interest has been attracted to the use of ESSs due to the variable nature of most renewable energy sources. ESSs can accommodate renewable generation in time-shifting its energy to match demand and avoid power curtailment. They can also be used to mitigate transmission congestion and hedge forecast errors, etc. In this context, appropriate sizing of storage systems is of importance not only for power system operation but also for economic consideration. In recent years, there has been extensive study on optimal sizing of ESSs for different applications with wind generation [2-10]. Analytical techniques are developed in [2, 3] to determine optimal ESS capacity with wind integration. Authors in [4] propose a methodology to size ESS for power balancing of variable energy generation. They use discrete Fourier transform to decompose the required balancing power into different time-varying periodic components, which can be used to quantify storage capacity for different types of energy storages. To deal with the variability and uncertainty of wind power, author in [5] decides optimal ESS capacity by using chance-constrained programming. He uses GA combined with Monte-Carlo simulation to solve the optimization problem in order to minimize energy cost while ensuring the difference between wind/ESS output and a predefined profile within a certain limit. In [6], a dynamic sizing based on statistical scenario forecasts and real market situation is presented to assess necessary ESS capacity in each period for compensating deviations from proposed bids of wind power. In [7], a methodology is presented to optimize storage capacity for mitigating wind prediction errors. Typical seven-scenario approximation of wind power forecast errors are modeled and used in the system tests. Authors in [8] focus on storage sizing for voltage support in low voltage distribution networks. The model is formulated as a stochastic optimization problem and a scenario reduction procedure is proposed to reduce the size of the problem. A probabilistic sizing method is proposed in [9] to reduce wind forecast uncertainty. The sizing method is based on persistence-based forecast data and verified with real-world forecast data. Authors in [10] apply two-stage stochastic optimization to optimally size ESSs while minimizing generation cost and maximizing the use of renewable generation. However, they do not incorporate transmission constraints and only focus on an intra-hourly economic dispatch. Stochastic programming has been widely employed to handle uncertainties in some areas of power system planning such as power capacity expansion planning [11-13] and transmission planning [14-16]. Thus, in this paper, we propose a two-stage stochastic model to determine the optimal capacity of ESSs in time-shifting wind generation considering wind and load variability. The planning problem is formulated as a non- linear programming AC OPF model, aiming to minimize ESS installation cost and system operation cost. Our contribution is providing a tool for long-term sizing of ESSs under the stochastic behavior of both wind generation and load. This planning tool allows to optimally size ESSs for time-shifting application considering investment costs, generation costs, taking into account network constraints. Extensive tests are performed on IEEE 39-bus system. The remaining of the paper is organized as follows: In section II, the mathematical formulation of the model is introduced. In section III, we Engineering, Technology & Applied Science Research Vol. 8, No. 3, 2018, 2901-2906 2902 www.etasr.com Le and Nguyen: Stochastic Sizing of Energy Storage Systems for Wind Integration describe the stochastic wind and load inputs. In section IV, we discuss a case study and simulation results and section V concludes the paper. II. MATHEMATICAL FORMULATION Generally, a two-stage stochastic programming problem has the form [17]: min ( ) + ( , )   s.t. g(x) ≤ 0 ( , ) is the optimal value of the second-stage problem: min ( , )  s.t. G(y, ) ≤ 0 where, x are first-stage decision variables, y are second- stage decision variables, ξ is a vector of random parameters, and the expectation operator at the first-stage problem (1) is taken with respect to probability distribution of the random parameters ξ. For random parameters with discrete distribution, i.e., ξ has a finite number of realizations (scenarios) ξω with respective probabilities pω, ω = 1,...,Nscen. Then: ( , ) ∑ = 1 ( , )  where, Nscen is total number of scenarios of the random parameters. In this paper, the two-stage stochastic programming model has the goal of minimizing ESS installation cost in the first-stage and total generation cost in the second-stage; therefore, the objective function is formulated as: min, ∑ + + ∑ ∑ ∑ + , + , + =1 , + ℎ ℎ ,   In (4), and , which are respectively energy rating and power rating of the ESS, are control variables of the first-stage. , , , , , are variables of the second-stage. They are respectively real generation power of generating unit at bus i in hour t and scenario ω, including power from both conventional and wind generators, and discharging and charging power of ESS at bus i in hour t and scenario ω. , , are cost coefficients of the generating unit at bus i. and are discharging and charging cost of the ESS at bus i. Ns and Ng are total number of ESSs to be installed and total number of generating units, respectively. T is the optimization period considered for the second-stage problem. In this model, the second-stage problem is based on a daily basis, and hence T=24 hours. C [$/kWh/day] and C [$/kW/day] are respectively ESS daily energy-related and power-related capital cost. Assuming a life time of N years, ESS total capital cost is converted into daily capital cost by multiplying the energy-related cost and power-related cost by the daily capital recovery factor [18]: = ( )( ) (5)  = ( )( )  (6) where, CB[$/kWh] and CR[$/kW] are energy-related and power-related capital costs of ESSs, r is annual interest rate and Nday is the number of days in a year, which is 365 days. The 3rd component in the objective function (4) is the implementation of complementary constraint to make sure an ESS is not charged and discharged simultaneously. This constraint is incorporated by applying suitable values for charging and discharging costs (cch and cd) of the ESS. Since the operational cost of charging cch is the locational marginal price (LMP) at the ESS bus, it is set to zero. To prevent simultaneous charging and discharging, the operational cost of discharging cd is set to a very small quantity, i.e., cd =10-2 [19]. This objective function has to fulfill both network constraints (7)-(13) and constraints for ESSs (14)-(20):  Power balance equations Include equations for real and reactive power at each node i in each time period t for each scenario ω: , − , + , − , = ∑ , , cos , −, + sin , − , (7) , − , + , − , = ∑ , , sin , −, − cos , − , (8) where, Nb is the total number of buses in the system. , and , denotes real and reactive power of load at bus i in hour t and scenario ω. , defines reactive power of generating unit at bus i in hour t and scenario ω. , and , are reactive charging and discharging power of ESS at bus i in hour t and scenario ω. , and , are voltage magnitudes of bus i and k in hour t and scenario ω. , and , are voltage angles of bus i and k in hour t and scenario ω. Gik and Bik are line conductance and line susceptance of branch ik.  Upper and lower limits for voltage magnitudes: ≤ , ≤ (9) where, and are lower and upper limits of voltage magnitude at bus i.  Bounds on real and reactive generation powers: Engineering, Technology & Applied Science Research Vol. 8, No. 3, 2018, 2901-2906 2903 www.etasr.com Le and Nguyen: Stochastic Sizing of Energy Storage Systems for Wind Integration ≤ , ≤ (10)  ≤ , ≤   where, and are lower and upper limits of real power generation of generating unit at bus i. and are lower and upper limits of reactive power generation of generating unit at bus i.  Branch current limits: 0 ≤ , ≤ (12) or 0 ≤ , ≤ (13) where, , and , are magnitude of currents flowing from bus i to bus j and from bus j to bus i in hour t and scenario , respectively. and are upper limits of currents flowing from bus i to bus j and from bus j to bus i, respectively.  ESS energy balance equations Include energy balance equation for each ESS in each period, considering charging and discharging efficiencies: , = , + ( , + , ⁄ )∆ (14) where, , and , are energy levels of ESS at bus i in hour t and t-1 in scenario . and are ESS charging and discharging efficiencies, respectively.  ESS energy continuity This constraint is added to make sure energy of ESS i at the end of the day is equal to energy at the beginning of the day B0. , = , (15) where, , and , are respectively energy level at the beginning of the day and at hour 24 of ESS at bus i in scenario .  ESS charging/discharging power bounds: ≤ , ≤ (16)  ≤ , ≤    ≤ , ≤    ≤ , ≤   where, is lower limit of real charging/discharging power of ESS at bus i. and are lower and upper limits of reactive charging/discharging power of ESS at bus i.  ESS energy limits: ≤ , ≤ (20) where, , is energy level of ESS at bus i in hour t and scenario . is lower limit of energy of ESS at bus i. The proposed two-stage stochastic model is a non-linear and non-convex problem. It is implemented in GAMS on a PC with Intel Core i7 – 3.4GHz CPU and 8.0GB of memory, using IPOPT solver with optimality gap of 0.5%. III. STOCHASTIC INPUT DATA Stochastic programming is a programming approach involving uncertain parameters which can be characterized by discrete distributions with a set of scenarios or realizations. In stochastic programming, deterministic mathematical formulation of any model is solved with this set of scenarios of uncertain parameters. Therefore, a necessary step in applying stochastic programming approach is generating a set of scenarios that realistically represents parameter uncertainties. For planning problems, the long-term planning horizon may result in a very large number of scenarios, which leads to intractability of the problems. To attain tractability with reasonable computation time for such planning problems, there is a need to approximate the original scenario set into a smaller subset, which still preserves essential features of the original one by using different scenario reduction techniques. In this paper, we consider the stochastic nature of both wind and load. Scenarios for these parameters are generated from historical hourly data using PCA-guided K-means clustering technique [20]. Daily wind and load data are clustered into distinct groups and the resulting representative wind and load profiles of each cluster along with their probabilities are obtained. K-means is known as one of the most popular and efficient clustering techniques in many applications because of its simple and efficient implementation. However, for high dimensional data set such as wind power and load, it is extremely difficult to identify their coherent patterns. Therefore, applying directly K-means technique to these data may result in clusters with less accuracy. A possible solution to this problem is the use of a dimension reduction technique such as principle component analysis (PCA) [21], which projects data to a lower dimensional subspace, where their patterns are more easily identified. In PCA-guided K-means clustering technique, K-means is performed in the PCA subspace to initialize centroids for clustering in the full data space. This PCA-guided algorithm is proved to be effective in obtaining better solutions and computation time, especially for such high dimensional data as wind and load. To prepare stochastic inputs for the model, we use historical 5-minute data of wind power and load from 2007 to 2013 [22].These wind and load data are averaged to create hourly data and then scaled down to a suitable level for the test Engineering, Technology & Applied Science Research Vol. 8, No. 3, 2018, 2901-2906 2904 www.etasr.com Le and Nguyen: Stochastic Sizing of Energy Storage Systems for Wind Integration system (Figure 1). The PCA-guided K-means clustering technique is applied to these data by first performing PCA analysis on the set of wind and load data to create another set of the same size of linearly independent variables called principle components (PCs). The first few components, which contain the highest amount of information (up to more than 90%) of the original data, are selected. The original data set is thus reduced to another with smaller dimensionality called PCA subspace. The selected PCs of the above wind and load data can be seen in Figure 2. K-means clustering is then performed in the PCA subspace to obtain initial cluster centroids and performed in the full space afterwards using the above initial cluster centroids. Finally, the representative scenario of each cluster is determined by averaging the members of the cluster. a) Wind data b) Load data Fig. 1. Input data for IEEE 39-bus case study. a) Wind PCs b) Load PCs Fig. 2. PCs of the input data. After performing clustering on wind and load data, we obtain a set of wind scenarios and load scenarios as shown in Figure 3, each of them associated with a probability. These wind and load scenarios are used to create sets of wind-load scenarios for the problem. Each set of scenarios consists of one scenario of wind and another scenario of load. These sets are created by listing all possible combinations of wind scenarios and load scenarios. The probability of any combined wind-load scenario is obtained by convolving the two probabilities of wind and load scenarios. IV. CASE STUDY AND DISCUSSION We performed tests on a modified IEEE 39-bus system (Figure 4), which has 9 conventional generators with a total capacity of 2200MW and system load with peak value of 2007MW. Wind farms were connected to buses 32 and 35. In order to test the model with different wind penetration levels, which is the ratio between wind farm installed capacity and the peak load, we assume the total installed capacity of both wind farms as 400MW, 600MW, and 800MW, accounting for 20%, 30%, and 40% of wind penetration respectively. An ESS system, with capacity to be determined, is installed at a wind bus, for example, bus 32. Three different ESS technologies are considered, i.e., battery energy xtorage (BES), compressed air energy storage (CAES) and pumped hydro storage (PHS). Parameters of these storage technologies are shown in Table I, assuming ηch=ηd=η [23–25]. For the sake of simplicity, we assume a correlation coefficient of 1 for the wind farms. If the correlation coefficient is different from 1, wind power of each wind farm will be represented by each set of wind scenarios, which results in an increase in the total number of scenarios. This will only increase the computation time and does not affect the accuracy of the proposed model. The daily capital cost of each ESS technology is calculated according to (5) and (6) is shown in Table II. a) Wind scenarios b) Load scenarios Fig. 3. Wind and load scenaios for IEEE 39-bus case study. Fig. 4. IEEE 39-bus system. TABLE I. PARAMETERS OF THE THREE ESS TECHNOLOGIES CB ($/kWh) CR ($/kW) N (years) η r (%) BES 330 400 15 0.85 5 CAES 5 700 30 0.79 PHS 14 1000 30 0.87 Engineering, Technology & Applied Science Research Vol. 8, No. 3, 2018, 2901-2906 2905 www.etasr.com Le and Nguyen: Stochastic Sizing of Energy Storage Systems for Wind Integration TABLE II. DAILY CAPITAL COST OF THE THREE ESS TECHNOLOGIES CBd ($/kWh/day) CRd ($/kW/day) BES 87.1 105.6 CAES 0.89 124.7 PHS 2.5 178.2 From the 6 wind scenarios and 5 load scenarios (each associated with a probability) shown in Figure 3, we combines each wind scenario with the 5 load scenarios to create a set of 30 wind-load scenarios as system input. The probability of each scenario in this set is obtained by convolving probabilities of the corresponding wind and load scenarios. We run the model with this set of 30 wind-load scenarios at different wind penetration levels. The test is performed first with BES system, and then repeated for the remaining ESS technologies considered, i.e., CAES and PHS. It turns out that due to the high investment cost of BES, the model does not choose to install this technology for any wind penetration level. In other words, installing BES in this case does not yield any benefit for the system since its capital cost overweighs the profit from the combined operation of wind and BES. The model, instead, chooses to install CAES and PHS. Optimal capacities of the CAES and PHS are shown in Table III. ESS of both technologies is installed at wind penetration level from 30% to 40%. At 20% wind penetration level, which is 400MW for both wind farms, since the amount of excess wind power to be time- shifted is also not enough to cover the investment cost of either CAES or PHS, the model does not install the storage device. At higher wind penetration levels (30% and 40%), the model decides to install either CAES or PHS system and capacity of the storage device is gradually increased with the wind penetration. Capacity of CAES is larger than that of PHS. This capacity difference is higher at higher wind penetration levels. As noted from Table I, even PHS technology has higher efficiency than CAES, its higher daily capital cost still results in smaller capacities of the storage device. TABLE III. OPTIMAL CAPACITY OF ESS IN 39-BUS SYSTEM Wind 20% 30% 40% PHS CAES PHS CAES PHS CAES Bmax (MWh) 0.0 0.0 74.4 83.1 109.3 162.8 Rmax (MW) 0.0 0.0 12.2 15.0 18.3 29.5 We also performed tests for the case when there is no ESS installed to compare the total yearly cost of the system in both cases, with and without the storage device (Table IV). This cost includes the yearly capital cost of the ESS and the expected value of the yearly generation cost. TABLE IV. TOTAL COST OF 39-BUS SYSTEM WITH AND WITHOUT ESS Wind 30% 40% No ESS 4.6355e+8 4.4588e+8 CAES 4.6329e+8 4.4526e+8 PHS 4.6343e+8 4.4547e+8 As can be seen from Table IV, the total yearly cost of the system is lower when ESS is installed. Moreover, the case of PHS results in a slightly higher yearly cost compared to that of CAES. The cost difference between these three cases is higher at higher wind penetration, when higher ESS capacity is installed. That is, the CAES gains a cost reduction of 0.06% at 30% wind penetration and 0.14% at 40% wind penetration while the PHS yields 0.03% cost reduction at 30% wind penetration and 0.09% at 40% wind penetration. Daily operation of the CAES in all scenarios at 30% wind penetration level, can be seen in Figure 5 as an example. Basically, in many scenarios (corresponding to load scenarios with higher peaks in Figure 3b), the CAES is charged at off-peak periods (hours 1 to 5 and hours 13 to 16) and then discharged at peak periods (hours 7 to 12 and hours 18 to 20). In other scenarios (corresponding to load scenarios with lower peaks in Figure 3b), it is charged at off-peak periods (hours 1 up to 7) and gradually discharged at high load periods (hours 9 to 11, 12 to 16, and 17 to 21). In other words, the storage has been fully employed to time-shift wind energy for matching the demand in all scenarios of wind and load. Fig. 5. Energy level of the CAES in daily operation. We also performed tests on a reduced model by removing the constraints on branch current limit, i.e., constraints (12) - (13), for branches with current flow less than 70% of the branch limit. The size of the model is thus reduced considerably and the computation time is also significantly less than that of the full model. In particular, computation time of the IEEE 39-bus system is reduced from 17 minutes to only 2 minutes. V. CONCLUSIONS In this paper, we propose a two-stage stochastic model to optimally size ESSs in power systems with high wind integration aiming to minimize ESS capital cost and system generation cost. The model is incorporated into an AC OPF problem, taking into account network constraints and ESS constraints. Stochastic input data of wind and load are generated from historical data using PCA-guided K-means clustering technique. The model is extensively tested on IEEE 39-bus with three ESS technologies and different wind penetration levels. Optimal sizes of the ESSs, energy rating and power rating, are explicitly determined. Simulation results show that the ESSs are efficiently employed in time-shifting wind energy to meet the demand. The model is solved using direct method and tests on the case study give reasonable computation time. Engineering, Technology & Applied Science Research Vol. 8, No. 3, 2018, 2901-2906 2906 www.etasr.com Le and Nguyen: Stochastic Sizing of Energy Storage Systems for Wind Integration ACKNOWLEDGMENT This work was supported by the University of Danang, University of Science and Technology, code number of project: T2017-02-95. REFERENCES [1] P. Taylor, Energy Technology Perspectives 2010 -Scenarios and Strategies to 2050, International Energy Agency, Paris, France, 2010 [2] H. T. Le, T. Q. Nguyen, “Sizing energy storage systems for wind power firming: An analytical approach and a cost-benefit analysis”, IEEE Power and Energy Society General Meeting – Conversionand Delivery of Electrical Energy in the 21st Century Pittsburgh, USA, July 20-24, 2008 [3] Z. Y. Gao, P. Wang, L. Bertling, J. H. Wang, “Sizing of energy storage for power systems with wind farms based on reliability cost and wroth analysis”, IEEE Power and Energy Society General Meeting, San Diego, USA, July 24-29, 2011 [4] Y. V. Makarov, P. Du, M. C. W. Kintner-Meyer, C. Jin, H. F. Illian, “Sizing Energy Storage to Accommodate High Penetrationof Variable Energy Resources”, IEEE Transaction on Sustainable Energy, Vol. 3, No. 1, pp. 34-40, 2012 [5] L. Li; L. Yang, “A chance-constrained programming based energystorage system sizing model considering uncertainty of wind power”, International Conference on Sustainable Power Generation and Supply (SUPERGEN 2012), Hangzhou, China, September 8-9, 2012 [6] P. Pinson, G. Papaefthymiou, B. Klockl, J. Verboomen, “Dynamic Sizing of Energy Storage for Hedging Wind Power Forecast Uncertainty”, IEEE Power & Energy Society General Meeting, PES 2009, Calgary, AB, Canada, July 26-30, 2009 [7] T. Hou, Y. Cui, Y. Li, W. Zhang, X. Zhou, “Sizing battery energy storage for wind farms based on wind power forecast uncertainty in the bulk power system”, IEEE PES Asia-Pacific Power and Energy Engineering Conference (APPEEC), Xi'an, China, pp. 1548-1552, 2016 [8] M. Bucciarelli, S. Paoletti, A. Vicino, “A scenario reduction approach for optimal sizing of energy storage systems in power distribution networks”, IEEE 56th Annual Conference on Decision and Control (CDC), Melbourne, Australia, pp. 4497-4502, December 12-15, 2017 [9] H. Bludszuweit, J. A. Dominguez-Navarro, “A Probabilistic Method for Energy Storage Sizing Based on Wind Power Forecast Uncertainty”, IEEE Transaction on Power Systems, Vol. 26, No. 3, pp. 1651-1658, 2011 [10] K. Baker, G. Hug, X. Li, “Optimal storage sizing using two-stage stochastic optimization for intra-hourly dispatch”, North American Power Symposium (NAPS), Pullman, USA, September 7-9, 2014 [11] A. P. Sanghvi, I. H. Shavel, “Investment planning for hydrothermal power system expansion: Stochastic programming employing the Dantzig-Wolfe decomposition principle”, IEEE Transaction on Power Systems, Vol. 1, No. 2, pp. 115-121, 1986 [12] B. G. Gorenstin, N. M. Campodonico, J. P. Costa, M. V. F. Pereira, “Power system expansion planning under uncertainty”, IEEE Transaction on Power Systems, Vol. 8, Vo. 1, pp. 129-136, 1993 [13] A. Marn, J. Salmern, “Electric capacity expansion under uncertain demand: Decomposition approaches”, IEEE Transaction on Power Systems, Vol. 13, No. 2, pp. 333-339, 1998 [14] M. Carrin, J. M. Arroyo, N. Alguacil, “Vulnerability- constrainedtransmission expansion planning: A stochastic programming approach”, IEEE Transaction on Power Systems, Vol. 22, No. 4, pp. 1436-1445, 2007 [15] H. Yu, C. Y. Chung, K. P. Wong, J. H. Zhang, “A chance constrained transmission network expansion planning method with consideration of load and wind farm uncertainties”, IEEE Transaction on Power Systems, Vol. 24, No. 3, pp. 1568-1576, 2009 [16] T. Akbari, S. Zolfaghari, A. Kazemi, “Multi-stage stochastic transmission expansion planning under load uncertainty using Benders decomposition”, International Review of Electrical Engineering, Vol. 4, No. 5, pp. 976-984, 2009 [17] A. Shapiro, D. Dentcheva, A. Ruszczynski, Lectures on Stochastic Programming: Modeling and Theory, Society for Industrial and Applied Mathematics and Mathematical Programming Society, Philadelphia, USA, 2009 [18] C. S. Park, Fundamentals of Engineering Economics, Pearson Education, Inc., New Jersey, 2004 [19] A. Castillo, D. F. Gayme, “Profit maximizing storage allocation in power grids”, IEEE 52nd Annual Conference on Decision and Control (CDC), pp. 429-435, 2013 [20] Q. Xu, C. Ding, J. Liu, B. Luo, “PCA-guided search for K-means”, Pattern Recognition Letters, Vol. 54, pp. 50-55, 2015 [21] D. D. Le, G. Gross, A. Berizzi, “Probabilistic Modeling of Multisite Wind Farm Production for Scenario-based Applications”, IEEE Transactions on Sustainable Energy, Vol. 6, No. 3, pp. 748-756, 2015 [22] Bonneville Power Administration, Wind Generation &Total Load in The BPA Balancing Authority, available at: http://transmission.bpa.gov/ Business/Operations/Wind/ (date of last access: February 21, 2018) [23] S. M. Schoenungand, W. V. Hassenzahl, Long- vs. Short-Term Energy Storage Technologies Analysis, Sandia National Laboratories, U.S. Department of Commerce, 2003 [24] S. M. Schoenungand, Characteristics and Technologies for Long- vs. Short-Term Energy Storage, Sandia National Laboratories, U.S. Department of Commerce, 2003 [25] S. M. Schoenungand, Energy Storage Systems Cost Update, Sandia National Laboratories, 2011