Optimization Models for Advanced Cyber-security in the Smart Grid by Seyedamirabbas Mousavian A dissertation submitted to the Graduate Faculty of Auburn University in partial ful llment of the requirements for the Degree of Doctor of Philosophy Auburn, Alabama August 2, 2014 Keywords: Cyber-security, Risk Mitigation, Phasor Measurement Units, Optimal placement, Mixed Integer Linear Programming, Arti cial neural networks Copyright 2014 by Seyedamirabbas Mousavian Approved by Jorge Valenzuela, Chair, Full Professor of Industrial and Systems Engineering Jianhui Wang, A liated Professor of Industrial and Systems Engineering Saeed Maghsoodloo, Professor Emeritus of Industrial and Systems Engineering Chase Murray, Assistant Professor of Industrial and Systems Engineering Acknowledgments I would rst and foremost like to thank my God because without Him, I would not have completed the process of obtaining a Ph.D. Secondly, I am eternally grateful for my advisor, Dr. Jorge Valenzuela, for continuing to guide my dissertation from afar. He continued to encourage, support, criticize and guide me even when it was uncertain where my dissertation was going. In addition, I would like to especially thank Dr. Jianhui Wang and the Argonne National Laboratory for the research and nancial supports throughout my studies. Also, I would like to thank my committee, Dr. Chase Murray and Dr. Allison Jones-Farmer, for their support in the dissertation process. Last but not the least, I would like to thank my lovely parents and sister, Aliasghar Mousavian and Zahra Jaberansari and Samaneh Mousavian, for being so encouraging and supportive. They continued to push me when I wanted to give up and listened to all of the issues I ran into even when it was not interesting. I am so grateful and thankful for them and would like to dedicate this hard-work to my family who loves regardless. ii Table of Contents Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ii List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . v List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1 Real-time Data Reassurance in Electrical Power Systems based on Arti cial Neu- ral Networks . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.1 Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.3 Optimal power ow model . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 1.4 Analysis of cyber-attacks to the network data . . . . . . . . . . . . . . . . . 8 1.5 Detection model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 1.5.1 Arti cial Neural Network Model . . . . . . . . . . . . . . . . . . . . . 12 1.5.2 ANN Reduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 1.5.3 Detection Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 1.6 Experimental Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 1.6.1 DC-OPF Detection Algorithm . . . . . . . . . . . . . . . . . . . . . . 19 1.6.2 AC-OPF Detection Algorithm . . . . . . . . . . . . . . . . . . . . . . 26 1.6.3 Computational Time Analysis . . . . . . . . . . . . . . . . . . . . . . 28 1.7 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 2 A Probabilistic Risk Mitigation Model for Cyber-Attacks to PMU Networks . . 33 2.1 Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 iii 2.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 2.3 Estimating the Threat Levels . . . . . . . . . . . . . . . . . . . . . . . . . . 37 2.4 Response Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 2.5 Experimental Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 2.5.1 6-bus Test System . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 2.5.2 24-bus Test System . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53 2.5.3 Dealing with Large Power Systems . . . . . . . . . . . . . . . . . . . 57 2.5.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58 2.6 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 3 Investment Decisions on Optimal Allocation of Phasor Measurement Units . . . 64 3.1 Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64 3.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 3.3 PMU Placement Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 3.3.1 Network Observability Rules . . . . . . . . . . . . . . . . . . . . . . . 68 3.3.2 Single-phase Optimal PMU Placement Model . . . . . . . . . . . . . 70 3.3.3 Case Study I: Two-phase Investment Approach . . . . . . . . . . . . 72 3.3.4 Case Study II: The Impact of Transmission Switching on Optimal PMU Placement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74 3.4 Two-phase optimal PMU placement model considering transmission switching 75 3.5 Experimental Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78 3.5.1 Optimal PMU Placement without Transmission Switching . . . . . . 78 3.5.2 Optimal PMU Placement with Transmission Switching . . . . . . . . 79 3.5.3 Computational Time Analysis . . . . . . . . . . . . . . . . . . . . . . 84 3.6 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 iv List of Figures 1.1 The 6-bus Test System . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 1.2 The network model under attack . . . . . . . . . . . . . . . . . . . . . . . . . . 11 1.3 Multi-layer ANN with two hidden layers . . . . . . . . . . . . . . . . . . . . . . 12 1.4 Block diagram of the inputs/outputs of the algorithm . . . . . . . . . . . . . . . 13 1.5 Flowchart of the anomaly detection algorithm . . . . . . . . . . . . . . . . . . . 17 1.6 24-bus test system . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 1.7 Standardized load . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 1.8 Power ows of four lines during the training period . . . . . . . . . . . . . . . . 20 1.9 Fitted Cumulative Distribution Function to the sum of the squared errors for ANN-2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 1.10 Plot of the ANN-1 false alarms . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 1.11 Plot of the ANN-2 false alarms . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 1.12 Plot of the ANN-3 false alarms . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 2.1 The 6-bus Test System . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 2.2 Threat levels in case of no response action . . . . . . . . . . . . . . . . . . . . . 49 v 2.3 E ect of on threat level of PMU2 . . . . . . . . . . . . . . . . . . . . . . . . . 50 2.4 Comparison of two potential responses . . . . . . . . . . . . . . . . . . . . . . . 52 2.5 E ect of computational time on threat level of PMU6 . . . . . . . . . . . . . . 52 2.6 IEEE 24-bus Test System . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53 2.7 Threat levels for no response to compromised PMU7 . . . . . . . . . . . . . . . 55 2.8 Threat levels for no response to compromised PMU20 . . . . . . . . . . . . . . . 55 2.9 E ect of single and multiple attacks on threat levels . . . . . . . . . . . . . . . 56 2.10 Maximum threat levels for two potential responses . . . . . . . . . . . . . . . . 57 3.1 Network Observability Rule 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 3.2 Network Observability Rule 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 3.3 Network Observability Rule 3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 3.4 IEEE 57-bus test system . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72 vi List of Tables 1.1 Load in MW . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 1.2 Power generation output with no cyber-attack . . . . . . . . . . . . . . . . . . . 9 1.3 Flows with no cyber-attack . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 1.4 Flows after the cyber-attack . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 1.5 Con gurations of ANN-1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 1.6 Solution of the MILP model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 1.7 Con gurations of ANN-2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 1.8 Detection results of one anomaly injection . . . . . . . . . . . . . . . . . . . . . 25 1.9 Cyber-attack scenarios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 1.10 DC-OPF detection results in di erent cyber-attack scenarios . . . . . . . . . . . 26 1.11 Con gurations of ANN-3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 1.12 AC-OPF detection results in di erent cyber-attack scenarios . . . . . . . . . . . 27 1.13 Training and Detection Computational Times . . . . . . . . . . . . . . . . . . . 28 2.1 Nodal distances between PMUs in the 6-bus test system . . . . . . . . . . . . . 49 2.2 Candidate responses for the 6-bus test system . . . . . . . . . . . . . . . . . . . 49 2.3 Optimal response action in the 6-bus system . . . . . . . . . . . . . . . . . . . . 50 2.4 Initial observabilites of the 24-bus system . . . . . . . . . . . . . . . . . . . . . 54 2.5 Nodal distances between PMUs in the 24-bus test system . . . . . . . . . . . . . 54 2.6 E ect of response time on the optimal solution in case of multiple attacks . . . 56 2.7 Optimal response for single cyber-attacks to the 24-bus test system . . . . . . . 58 vii 3.1 Optimal number of PMUs for full observability . . . . . . . . . . . . . . . . . . 71 3.2 Alternative placements of the IEEE 57-bus system . . . . . . . . . . . . . . . . 73 3.3 Two-phase investment plan for alternative solutions of the IEEE 57-bus system 74 3.4 Transmission switiching scenarios . . . . . . . . . . . . . . . . . . . . . . . . . . 75 3.5 Two-phase investment plan for the IEEE test systems . . . . . . . . . . . . . . 78 3.6 Optimal number of PMUs in large power systems . . . . . . . . . . . . . . . . . 79 3.7 Transmission switiching scenarios for IEEE 118-bus system . . . . . . . . . . . . 80 3.8 Optimal placements for the transmission switching scenarios of IEEE 118-bus system . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81 3.9 Optimal placement for IEEE 118-bus system considering transmission switching scenarios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82 3.10 Optimal number of PMUs for removed transmission switching scenarios . . . . . 83 3.11 Computational Times (s) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84 viii Introduction The power grid is increasingly dependent on information and communication technolo- gies, which puts more emphasis on the importance of power system security as one of the top priorities. Internal and external factors can put the security of the power system at risk. The external factors include cyber-terrorist attacks, sabotage and environmental im- pacts while the internal factors are inherent to the accuracy of power system applications and their associated input data. As the utility industry becomes more automated and relies more on automated devices, the major threat to the grid is shifting from equipment failures to cyber-security attacks. Hence, improvement of the cyber-security of the power grid along with the reliability of the power system operations are primary focus of the United States government, power industry executives, and the research community. Cyber-attacks to power systems could cause huge nancial losses and substantial dam- ages to the power grid. Hence, the main goal of the research summarized in this dissertation is to highlight the risks associated to the cyber-attacks to power systems and develop al- gorithms and models to improve the safe operations of the electric power grid. To achieve this goal, a detection algorithm is developed to avoid cyber-attacks to the optimal power ow (OPF) software module of the power systems. Since there is no guarantee that detec- tion models detect all potential cyber-attacks, it is necessary to equip power systems with response models which avoid propagation of the cyber-atacks. In this dissertation, cyber- attacks to a speci c network of the power systems, known as the phasor measurement unit (PMU) network, have been discussed, and a risk mitigation model for cyber-attacks to PMU networks is proposed. PMUs provide the system operators with the real-time operating sta- tus of the power grid, which can be further utilized for better cyber-security of the power grid. Since PMUs require considerable capital investments, an investment decision model is 1 developed in this dissertation for the optimal allocation of phasor measurement units. There- fore, this dissertation outlines a research that will consider the OPF and state estimation modules of the electrical power systems and deliver the following outcomes. An algorithm for real-time data reassurance in the OPF module A probabilistic risk mitigation model for cyber-attacks to PMU networks An investment decision model for the optimal allocation of phasor measurement units The remainder of this dissertation is organized as follows. Chapter 1 explains the cyber- security problem of the OPF module and introduces the real-time data reassurance model to improve the cyber-security of the OPF module. Chapter 2 provides background information about cyber-security of the PMU networks and explains the probabilistic risk mitigation model for cyber-attacks to the PMU networks. Chapter 3 discusses the challenges of the optimal placement of PMUs (OPP) problem and describes the developed investment decision model to allocate PMUs with minimum capital costs. 2 Chapter 1 Real-time Data Reassurance in Electrical Power Systems based on Arti cial Neural Networks 1.1 Abstract Power system security is vulnerable to cyber-attacks that may cause signi cant damages to the power grid and result in huge nancial losses. In this paper, we show the risks associ- ated with cyber-attacks and propose an arti cial neural network-based protection approach. The proposed algorithm can monitor the output of power ow calculations and detect data anomalies in real-time. The network observability rules are formulated as a mixed integer linear program (MILP) problem. The results of the MILP problem are used to decrease the amount of data input required by the algorithm while the system stays observable. We run our experiments on the IEEE 24-bus reliability test system. The experimental results show that the developed algorithm is a promising enhancement to ensure data integrity in control centers. keywords: Power system security, Arti cial neural networks, Cyber-security, Network ob- servability Nomenclature Sets k: Set of lines to/from the bus k : Set of zero injection buses k: Set of lines to/from the zero-injection bus k2 : Set of generators 3 Constants K: Number of buses L: Number of lines Lk: Load at bus k (l): From bus of line l (l): To bus of line l Hk;l: Incidence matrix coe cient (-1, 0 or 1) at bus k of line l Sl: Susceptance of line l (Siemens) Cg: Energy cost of generator g ($/MWh) Pming : Minimum generation of generator g (MW) Pmaxg : Maximum generation of generator g (MW) Fmaxl : Maximum capacity of line l (MVA) Nk: Number of transmission lines to/from zero-injection bus k Variables pg: Energy dispatched by generator g2 (MW) fl: Power ow at line l = 1;;L (MVA) k Voltage angle at bus k = 1;::;K (Degrees) xl: Binary variable which equals 1 if line l is chosen x0l: Binary variable which equals 1 if the ow on line l can be computed yk: Binary variable which equals 1 if bus k is chosen y0k: Binary variable which equals 1 if the voltage on bus k can be computed i: Threshold value i SSEt: The mean of the sum of squared errors at time t F : Power ow matrix 4 1.2 Introduction Power system security is one of the top priorities for control center operators. Internal and external factors can put the security of the power system at risk. The external factors include cyber-terrorist attacks, sabotage and environmental impacts [1] while the internal factors are inherent to the accuracy of power system applications and their associated input data. As the utility industry becomes more automated and relies more on automated devices, the major threat to the grid is shifting from equipment failures to cyber-security attacks [2]. According to [3], cyber threats happen when unauthorized users exploit cyber system vulnerabilities. A cyber terrorist could wisely design a malicious data-tampering attack to deliberately in ict major damage on the power grid. The intruder can gain access to the supervisory control of a SCADA system and initiate control actions. Several software modules are used by power system operators to support decision making in the control centers. As a case in point, the state estimation software gives system operators an updated picture of the system status by estimating the actual values of system variables using real- time data. The software estimates the voltage magnitudes and voltage angles at all network buses [4]. The e ectiveness of state estimation can be a ected by bad data stemming from equipment installation problems, localized equipment failures, communication errors, etc. It has been also pointed out in [5] that the cyber-attackers may take advantage of the bad data for nancial arbitrage such as virtual bidding at selected pairs of nodes. Since state estimation procedures already consider that data measurements can be bad, existing procedures have been modi ed to detect malicious attacks to the data [6]. The weighted least squares (WLS) method, which solves Gauss normal equations iteratively, was initially used in state estimation. And consequently, the Bad Data Detection, Identi cation and Elimination (BDDIE) method was proposed to detect data attacks [6]. These orthogonal transformation techniques are not widely used in the power system community due to the required high computation e orts in large systems [7]. The bad data suppression (BDS) algorithm, based on a non-quadratic cost function, was proposed to improve the performance of the WLS 5 technique in the presence of bad data [7]. In this technique, the least normalized residuals (LNR) are used in detection and elimination of bad data. This use of residual analysis and non-quadratic estimation criteria laid the groundwork for the concept of interacting vs. non- interacting bad data and the ability to probabilistically predict false alarms [8]. Valenzuela et al. [9] considered another important software module used by the control centers, the Optimal Power Flow (OPF), that can be a target for a data attack. OPF determines the steady-state operation point which ensures the minimum generation cost while maintaining system constraints on real and reactive power, generator outputs, transmission line ows, bus voltages, etc. The authors pointed out that an undetected cyber-attack on the input data to the OPF module could cause power to be dispatched erroneously, overloading transmission lines and possibly resulting in cascading power outages. Traditional state estimation can detect di erences in database parameters and estimated parameters during the bad data detection process; however, cyber terrorists could carefully design an attack after the most recent state estimation program has run that is undetectable and intentionally designed to hamper the grid. In [10], it was shown that false data injection attacks cannot be avoided in todays SCADA systems. The notion of a false data injection attack was rst introduced in 2009 [11]. Besides, it is discussed in [12] that current bad data detection schemes would not detect all types of parameter changes, especially when branch power ows and power injections at both ends of a branch are critical to estimate the conventional state vector. Hence, it can be inferred that data manipulations would remain undetected and the incorrect database parameters would be used for later decisions. It has been pointed out also in [13] that there is a crucial need for SCADA real-time intrusion detection algorithms to mitigate the risk of cyber threats. To our knowledge, the data supplied to the OPF module is not well- protected against these kinds of attacks and could be an attractive target for cyber-attackers. In a similar research to this paper, Valenzuela et al. [9] proposed a bad data detection algorithm that monitors the AC power ow results of the OPF. The algorithm uses Principal Component Analysis (PCA) to determine whether the power system input data has been 6 compromised. In this paper, we address the same problem as in [9], but we use a di erent approach which allows for several enhancements to our algorithm. Our approach is based on a forecasting technique while the other paper uses a variability monitoring technique. The advantage of our technique is that the threshold value can be computed statistically. The threshold value in [9] is obtained by experimentation. Another di erence is that [9] considers just the AC OPF software while we consider both the AC and DC OPF. Lastly, we show by an example how a cyber-attack to the network data can endanger the physical power system, which is not discussed in [9]. We use arti cial neural network (ANN) to verify the trustworthiness of the results from the OPF. ANN has been used extensively in nonlinear systems such as the single-ended fault location of transmission lines [14], short term load forecasting [15] and power transformer fault diagnosis [16]. However, to the best of our knowledge, it has not been used in detecting data anomalies in power system applications. We also model network observability rules as a mixed integer linear programming (MILP) problem to reduce the dimensionality of the problem while still maintaining the critical variables. The rest of this chapter is organized as follows: Section 1.3 describes the optimal power ow model. Section 1.4 analyzes the cyber-attacks to the network data. Section 1.5 discusses the anomaly detection model and the ANN algorithm. Section 1.6 provides experimental results and Section 1.7 reports our conclusions. 1.3 Optimal power ow model Modern power system control centers run a sophisticated collection of computer appli- cations and maintain huge databases that ensure the economical operations of the power system. The input data to various application modules is cleansed and sampled for compu- tation, storage and further analysis. In this paper we model a system where a cyber terrorist has compromised the integrity of the data supplied to the OPF module. We focus on the OPF module because it plays a signi cant role in power generation and transmission, and an 7 undetected attack to the input data to the OPF could be disastrous to the power grid. OPF is an optimization-based module which minimizes the total generation cost of the system. For simplicity, a simpli ed DC OPF is presented in this section. However, in Section 1.5, we provide the results for both the DC and AC models. z = min X g2 Cgpg (1.1) Subject to LX l=1 Hk;lfl +pk = Lk k = 1;:::;K (1.2) fl Sl( (l) (l)) = 0 l = 1;:::;L (1.3) Fmaxl fl Fmaxl l = 1;:::;L (1.4) Pming pg Pmaxg 8g2 (1.5) The objective function in equation (1.1) is subject to power balance constraints at each bus k given in equation (1.2). The power ow on each line l is shown in equation (1.3). Constraints in equation (1.4) represent thermal ow limits for all lines, and constraints in equation (1.5) are generation capacities for each generator. The variables k are unrestricted. 1.4 Analysis of cyber-attacks to the network data The main goal of this research is to provide the system operators with an approach to protect the power system from attacks against the network model used in the OPF calculation. The network model represents the physical parameters of the network and is stored in a database. We use a 6-bus test system given in [17] to show how a cyber-attack to the network data can endanger the physical power system. The 6-bus test system is shown in Figure 1.1. This test system includes six buses, three generators and eleven transmission lines. TABLE 1.1 shows the load at each bus. 8 B u s 4 B u s 5 B u s 1 B u s 6 B u s 3 B u s 2 L i n e 1 L i n e 2 Line 3 Line 4 L i n e 5 L i n e 6 Line 7 L i n e 8 Line 9 L i n e 1 0 L i n e 1 1 Figure 1.1: The 6-bus Test System Table 1.1: Load in MW Bus 1 2 3 4 5 6 Load (MW) 0 0 0 80 80 80 First, we run the DC version of the OPF software to obtain power generation outputs and ows in the transmission lines assuming that there is no contaminated data. We use MatPower [18] in this paper for the OPF calculations. The power generation results are given in TABLE 1.2, and TABLE 1.3 shows the ows on transmission lines. The last column of TABLE 1.3 shows the percentage of transmission line capacity used for dispatching. Table 1.2: Power generation output with no cyber-attack Generator 1 2 3 Total pg (MW) 46.77 103.18 90.05 240 Next, we assume that intruders attack the database by getting access to the physical parameters of the network and carefully making changes to hamper the power grid. In addition, we assume that the system operator is unaware of these changes and operates the 9 Table 1.3: Flows with no cyber-attack Line Fmaxl fl(MVA) % Line Capacity Used 1 40 -1.37 3.4 2 60 26.67 44.5 3 40 21.46 53.6 4 40 -2.71 6.8 5 60 56.08 93.5 6 30 22.37 74.6 7 90 26.07 29.0 8 70 28.42 40.6 9 80 58.92 73.6 10 20 2.76 13.8 11 40 -4.99 12.5 system using the tampered data. The attacks consist of changing the network topology by modifying the origin and destination buses of some transmission lines. Mathematically, the network topology is represented by Hk;l in equation (1.2) of the OPF model. In the 6-bus test system, line 1 distributes the power from bus 1 to bus 2. Line 5 goes from bus 2 to bus 4. We assume that the attack aims at changing the destination of line 1 to bus 3 and the destination of line 5 to bus 3. Notice that these cyber-attacks will change only the network topology and there will be no change in the physical parameters of the network. Figure 1.2 shows the network model under attack. The dashed lines are the actual lines in the physical network which the attacker modi es. As a result, the network topology stored in the database no longer maps to the real system. Because the network model does not match the physical network, the physical ows will be di erent from the ows obtained by the OPF software. To compute the physical ows after the cyber-attack, we used MatPower [18], set the input voltage angles at all buses and obtained physical power ows. TABLE 1.4 shows the physical power ows on each line when the system is under cyber-attack. As shown in TABLE 1.4, line 9 is obviously overloaded. Overloaded transmission lines are highly undesirable to the grid stability since the overloaded transmission lines may lead to cascading outages and blackouts. In addition to possible overloaded lines, an attack can a ect the economics of power systems by using more expensive generating units or taking 10 B u s 4 B u s 5 B u s 1 B u s 6 B u s 3 B u s 2 L i n e 1 L i n e 2 Line 3 Line 4 L i n e 5 L i n e 6 Line 7 L i n e 8 Line 9 L i n e 1 0 L i n e 1 1 Line 1 - Attacked Line 5 - Attacked Figure 1.2: The network model under attack Table 1.4: Flows after the cyber-attack Line Power Flow % Line Capacity Used 1 0.34 0.8 2 25.58 42.6 3 11.78 29.4 4 -37.00 92.5 5 50.47 84.1 6 11.55 38.5 7 -0.26 0.29 8 48.91 69.9 9 91.98 115.0 10 -3.96 19.9 11 -11.72 29.3 advantage of nancial arbitrage in virtual bidding at di erent nodes [19]. Therefore, it is extremely important to detect this type of cyber-attack to the system. 1.5 Detection model In this section, we describe a model that allows the system operators to detect whether the network model has been compromised. The detection model is based on ANN, which 11 has been used extensively in the area of power systems. In particular, ANN has broad applications in fault diagnosis methods because of its ability to deal with noisy inputs, non- linear function approximation, and adaptive learning. 1.5.1 Arti cial Neural Network Model ANNs are usually trained o ine and then used to detect faults online. As shown in Figure 1.3, an ANN is composed of interconnected layers of arti cial neurons. ANN is a powerful computational model which can adjust the values of the connections, namely the weights, to approximate almost any nonlinear function [20]. Input Layer Output Layer Hidden Layer 2 Hidden Layer 1 Figure 1.3: Multi-layer ANN with two hidden layers In our detection model, ANN is used to estimate the power ows at period t given the load and the power ows at the previous period t 1. To detect anomalies the estimated power ows are then compared to the power ow output from the OPF module in the same period to determine whether an anomaly exists. The ows and the load at the previous period are the input data to the ANN model. The load is included to eliminate the variations stemming from load variability. Figure 1.4 shows a block diagram of the inputs/outputs of the algorithm. The MATLAB Neural Network Toolbox described in [20] is used to train the ANN module of the detection algorithm. 12 Trained ANN Loads at Period t -1 Power Flows at Period t -1 Estimate of Power Flows at Period t OPF Loads at Period t Generators Status at Period t Power Flows at Period t Is there any anomaly? Figure 1.4: Block diagram of the inputs/outputs of the algorithm 1.5.2 ANN Reduction The number of input vector entries to the ANN model is determined by the number of transmission lines plus the number of buses where loads are present. Since one of the factors a ecting the computational e ciency of the detection algorithm is the size of the input vector, we use the concept of network observability de ned in [21] to reduce the size of ANN. The network observability is a measure for how well to infer the state of a system by knowledge of a subset of its external outputs. Hence, we can just select a subset of transmission ows and buses such that the power system can still remain observable with less data. The problem of nding these subsets is formulated as a MILP model where the constraints assure the observability of the power network. The MILP model is formulated as follows: min LX l=1 xl + KX k=1 yk (1.6) Subject to y0 (l) = y (l)xl 8l such that (l) =2 (1.7) x0l = y (l)y (l) 8l such that (l) or (l) =2 (1.8) X l2 k xl ! x0l0 Nk 2 8k2 and8l0 2 k (1.9) xl +x0l 1 l = 1;:::;L (1.10) 13 yk +y0k 1 8k =2 (1.11) X l2 k xl 1 k = 1;:::;K (1.12) xl;x0l;yk;y0k2f0;1g l = 1;:::;L and8k =2 (1.13) In the objective function, equation (1.6), the rst term is the total number of selected transmission lines while the second term is the total number of chosen buses. Constraints equation (1.7) assure that if we know the ow of a branch and the bus voltage on one end, then the bus voltage on the other end can be calculated via the power ow equations. Since the term y (l)xl is the product of two binary variables, these constraints are nonlinear. These constraints can be represented by a group of linear constraints as follows: 2y0 (l) y (l) xl 0 (1.14) y0 (l) y (l) xl 1 (1.15) These two constraints set the binary variable y0 (l) to be the product of y (l) and xl. This is if y (l) = 0 or xl = 0, then y0 (l) = 0. The value of y0 (l) is 1 only if y (l) = 1 and xl = 1. Equation (1.8) assures that if we know the voltages at both buses of a transmission line, the line ow can be calculated. Again, the term y (l)y (l) is a nonlinear term. Similarly, these constraints can be represented by two linear constraints as follows: 2x0l y (l) y (l) 0 (1.16) x0l y (l) y (l) 1 (1.17) Equation (1.9) assures that for a zero-injection node, if just one of the incidence trans- mission line ows is unknown, then it can be calculated by Kirchho s law. Equation (1.10) and equation (1.11) assure that all bus voltages and transmission ows are either directly or 14 indirectly computed. Since the ANN predicts power ows, equation (1.12) assures that at least one transmission line from each bus is included in the reduced set of inputs. 1.5.3 Detection Algorithm The detection algorithm is theoretically grounded on forecasting the power ows of the current hour t based on the power ows of the previous hour t 1. Forecasted power ows at hour t are compared to the power ows generated by the OPF module at hour t. The di erence is computed by the sum of squared errors (SSE) of the forecasted ows and the ows from the OPF. If the SSE is greater than a threshold value, the alarm informs the system operator that the power system may have been compromised. We use ANN as the forecasting tool since it has been shown in other applications [14{16] to have excellent multivariate forecasting capabilities. We describe our methodology by the following pseudo code: Step 0. Simulate the generators status Simulate the generators status (on/o ) from a continuous-time Markov chain Step 1. Generate historical power ows tSim = Simulation Time for t = 1 to 2 tSim+ 1 do Sample generator status Solve the OPF using the load at time t Compute the tth row of the power ow matrix F end for Step 2. Create and Train the Neural Network Create the matrix Training-Input by using the entries 1 to tSim 1 of the ow matrix F as the input to the training function of the neural network Create the Matrix Training-Target by using the entries 2 to tSim of the ow matrix F as the target of the training function of the neural network 15 Create the neural network object Train the ANN using the matrices Training-Input and Training-Target Step 3. Statistically compute the threshold value for t = tSim+ 1 to 2 tSim+ 1 do Obtain an estimate for the ows at time t, ^P, by using the load and the ows at time t 1 as input to the neural network Compute the vector ^U by dividing ^P by the maximum emergency capacity of the transmission lines Compute the tth row of the used-capacity matrix U by dividing the tth row of the actual power ow matrix F by the maximum emergency capacity of the transmission lines Compute the sum of squared error SSEt between the tth row of U and ^U end for Compute a threshold value based on tting a Weibull distribution function to SSEt Step 4. Test Detection Algorithm for t = 2 tSim+ 2 to T do Sample generator status Obtain an estimate for the ows at time t, ^P, by using the load and the ows at time t 1 as input to the neural network Compute the vector ^U by dividing ^P by the maximum emergency capacity of the transmission lines Solve the OPF using the load at time t Compute the vector U by dividing the power ows obtained from the OPF by the maximum emergency capacity of the transmission lines Compute the mean of the sum of squared error SSEt between ^U and U If SSEt > then an anomaly exists Else the system is in control end if 16 end for We have scaled down the power ows (matrix U) to values between 0 and 1 to avoid a computer over ow when computing the sum of the squared errors. A summary of the algorithm is shown as a owchart in Figure 1.5. Create the ANN Train the ANN using the first half of the historical data (Time 1 to tSim) Use the ANN to obtain a prediction for the flows at time t Scale down the predicted power flows and the actual power flows Compute and store the SSE at time t Calculate the threshold value based on the distribution of SSE from time tSim+1 to time 2tSim+1 t = t + 1 t = 2?tSim + 1 Y e s No Train the ANN Determine the Threshold Value Anomaly Detection Sample generator status Store the loads and actual power flows Solve the OPF using the load at time t t = tSim + 1 Use the ANN to obtain a prediction for the flows at time t Scale down the predicted power flows and the actual power flows Compute the SSE at time t Sample generator status Store the loads and actual power flows Solve the OPF using the load at time t t = t + 1 SSE > Threshold Anomaly exists The system is under control No Yes Figure 1.5: Flowchart of the anomaly detection algorithm 1.6 Experimental Results To assess the performance of our detection algorithm, we use data from the IEEE Reliability Test System [22]. We run our experiments on the 24-bus system which consists of 38 transmission lines, 24 buses of which 17 have loads, and 33 generators. A picture of the 24-bus test system is provided in Figure 1.6. We refer the readers to [18] and [22] for additional details. The load data are obtained from the PJM [23]. A standardized load pro le of 744 hours by scaling down the load data based on the peak load is created so that the value 1.0 17 L i n e 1 L i n e 2 Line 3 Line 4 L i n e 5 L i n e 6 Line 7 L i n e 8 L i n e 9 L i n e 1 0 L i n e 1 1 Line 12 Line 13 Line 14 L i n e 1 5 L i n e 1 6 L i n e 1 7 L i n e 1 8 L i n e 1 9 L i n e 2 0 Line 21 L i n e 2 2 Line 23 L i n e 2 4 L i n e s 2 5 , 2 6 Line 27 L i n e 2 8 Line 29 Line 30 L i n e 3 1 L i n e s 3 2 , 3 3 Lines 34,35 L i n e s 3 6 , 3 7 L i n e 3 8 B u s 8 B u s 9 B u s 1 0 B u s 7 B u s 2 2 B u s 1 1 B u s 2 3 B u s 6 B u s 5 B u s 1 B u s 2 4 Bus 21 Bus 20Bus 19 B u s 1 8 Bus 15 B u s 3 B u s 4 B u s 2 B u s 1 4 B u s 1 3 B u s 1 6 B u s 1 7 B u s 1 2 Figure 1.6: 24-bus test system corresponds to the peak load. The load at each hour is calculated by applying the load pro le to the IEEE test system. Although the load at each bus follows the same variability pattern, the actual values are distinct at di erent buses. Figure 1.7 is a plot of the standardized load of the rst 168 hours, used to train the neural network. 18 0 20 40 60 80 100 120 140 1600.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Load (MWh) Time (Hours) Figure 1.7: Standardized load 1.6.1 DC-OPF Detection Algorithm In this section, we develop the detection algorithm for DC-OPF. We rst construct the con guration of the ANN. The loads of buses and power ows of transmission lines (at time t 1) are the input variables to the ANN. The output variables are the forecasted power ows for time t. ANN Con guration We use the multi-layer feed-forward back-propagation technique with an adaptive learn- ing rate and momentum. The learning rate, which controls the rate of convergence, is chosen to be 0.05. The training process is stopped when either the minimum performance gradient is reached or the performance goal is met. We simulate the operation of the power system using the DC-OPF for 168 hours (1/4 of the total hours of available data) to generate the historical power ows which are used to train the neural network. To simulate generator fail- ures, we sample the state of the generators according to a continuous-time two-state Markov 19 chain. Figure 1.8 shows the variations of the power ows on four of the lines in the 24-bus test system during the training period. 0 20 40 60 80 100 120 140 160?300 ?200 ?100 0 100 200 300 Time (Hours) Flows (MWA) Figure 1.8: Power ows of four lines during the training period We design two con gurations, ANN-1 and ANN-2. ANN-1 uses all buses and lines (62 inputs) while ANN-2 uses the reduced network obtained by solving the MILP formulation described in section 1.4. The performance of the ANN depends on its con guration and the best con guration depends on the power network. To nd the best con guration for a given power system, we have developed a computer code. The computer code tries di erent number of neurons for the ANN and computes the mean and standard deviation of the sum of squared errors (MSSE, SSSE) of the target output and the ANNs outputs of each con guration. The performance of an ANN con guration is de ned as MSSE+SSSE and computed by using the data for hours 169 to 337. This process needs to be executed every time a new system is studied. In our experiments, a total of 500 di erent ANN con gurations are generated by chang- ing the number of layers and the number of neurons on each layer. TABLE 1.5 shows the statistics of the ten best con gurations for ANN-1. Based on the criterion mentioned above, 20 the con guration in the rst row is selected for ANN-1, which is a four-layer feed-forward ANN. This network is activated by the tan-sigmoid (2 layers), log-sigmoid (1 layer) and linear (output layer) functions. Table 1.5: Con gurations of ANN-1 Number of Neurons The Sum of Squared Errors Con guration Layer 1 Layer 2 Layer 3 Min Mean Max Standard Deviation 1 17 17 15 0.0002 0.0362 0.2156 0.2154 2 17 17 18 0.0002 0.0362 0.2644 0.2642 3 15 20 15 0.0002 0.0362 0.2973 0.2935 4 15 15 15 0.0002 0.0364 0.2965 0.2963 5 14 18 18 0.0003 0.0377 0.2385 0.0547 6 17 15 16 0.0005 0.0377 0.3722 0.0652 7 7 7 0 0.0005 0.0394 0.5256 0.5251 8 8 10 0 0.0003 0.0405 0.2214 0.2211 9 11 8 0 0.0005 0.0462 0.2653 0.2648 10 13 8 0 0.0002 0.0462 0.3027 0.3025 Similarly, we determine the best con gurations for ANN-2. This network uses as input the load at buses and power ows on lines determined by the mentioned MILP formula- tion. We solve the optimization model for the 24-bus test system. The solution is given in TABLE 1.6, which shows that by applying the network observability rules we decrease the size of the input vector by 42% (from 62 to 36 inputs). For such a small system, the reduction of the number of input/output variables may not be strictly necessary. However, for real world power systems where thousands of buses and transmission lines may exist, without reduction the number of input variables could be extremely large and prohibited for 21 computation. The advantage of using network observability rules is that the crucial variables are maintained in the reduced set. Table 1.6: Solution of the MILP model Included buses in the input vector 1, 2, 3, 5, 6, 8, 9, 10, 14, 15, 16, 19 Included lines in the input vector 3, 4, 5, 7, 8, 11, 14-22, 26-28, 30-32, 34, 36, 38 We train 500 di erent ANN con gurations and select the one with the best performance as is done with ANN-1. TABLE 1.7 shows the ten best con gurations for ANN-2. The con guration in the rst row is chosen. Table 1.7: Con gurations of ANN-2 Number of Neurons The Sum of Squared Errors Con guration Layer 1 Layer 2 Layer 3 Min Mean Max Standard Deviation 1 9 13 10 0.0001 0.0145 0.0976 0.0179 2 9 10 8 0 0.0156 0.0908 0.0181 3 7 11 12 0.0002 0.0162 0.1055 0.0203 4 13 9 7 0.0001 0.0174 0.1021 0.0215 5 12 11 9 0.0001 0.0178 0.1003 0.0217 6 12 7 8 0.0003 0.0194 0.0919 0.0223 7 10 12 7 0 0.0194 0.1158 0.0240 8 12 7 10 0.0001 0.0193 0.1217 0.0245 9 8 10 8 0.0001 0.0193 0.1153 0.0247 10 9 10 9 0.0002 0.0202 0.0966 0.0238 Determining the Threshold Value First, we statistically determine the threshold value 1 for the ANN-1 model and cal- culate the sum of the squared error SSE at each hour for the next 168 hours which are 22 hours from 169 to 337. The values of SSE are tted to a Weibull distribution, which is commonly used in reliability engineering and failure analysis. For ANN-1, the scale and shape parameters are estimated to be 0.031 and 0.80 respectively. Assuming that a 2.5% false alarm rate is acceptable, the 97.5 percentile is considered as the threshold value, i.e. 1 = 0:16. Figure 1.9 shows the cumulative distribution function of the sample data and the hypothesized distribution. Similarly, the scale and shape parameters are estimated to be 0.013 and 0.83 for ANN-2, and the threshold value is set to 2 = 0:05. 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.20 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Sum of Squared Errors Cumulative probability Sum of Squared ErrorsFitted Weibull Distribution Figure 1.9: Fitted Cumulative Distribution Function to the sum of the squared errors for ANN-2 Testing the Detection Algorithm We rst study the performance of the detection algorithm under normal operation con- ditions. We run our simulation for the next 407 hours (hour 338 to hour 744). The nor- mal operation includes load variability and generator outages. Figure 1.10 and Figure 1.11 show the time series of the sum of the squared errors SSEt for ANN-1 and ANN-2, respec- tively. Three false alarms are observed in Figure 1.10 and seven false alarms are observed in Figure 1.11. 23 0 50 100 150 200 250 300 350 400 4500 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2 Time (Hours) Sum of Squared Error Threshold Figure 1.10: Plot of the ANN-1 false alarms 0 50 100 150 200 250 300 350 400 4500 0.02 0.04 0.06 0.08 0.1 0.12 Time (Hours) Sum of Squared Error Threshold Figure 1.11: Plot of the ANN-2 false alarms To study the performance of the detection algorithm under a cyber-attack, we change the database of the power system at one particular point in time. At hour 744, we change the phase shifter angle from 0 to =6 radians of each line, one line at a time. We run the two 24 detection algorithms. TABLE 1.8 shows whether the change of a particular line is detected using either detection algorithm. The results show that ANN-1 and ANN-2 detect 37 out of the 38 injected changes. Table 1.8: Detection results of one anomaly injection Line ANN-1Detection ANN-2Detection Line ANN-1Detection ANN-2Detection Line ANN-1Detection ANN-2Detection 1 Yes Yes 14 Yes Yes 27 Yes Yes 2 Yes Yes 15 Yes Yes 28 Yes Yes 3 Yes Yes 16 Yes Yes 29 Yes Yes 4 Yes Yes 17 Yes Yes 30 Yes Yes 5 Yes Yes 18 Yes Yes 31 Yes Yes 6 Yes Yes 19 Yes Yes 32 Yes Yes 7 Yes Yes 20 Yes Yes 33 Yes Yes 8 Yes Yes 21 Yes Yes 34 Yes Yes 9 Yes Yes 22 Yes Yes 35 Yes Yes 10 Yes Yes 23 Yes Yes 36 Yes Yes 11 Yes No 24 Yes Yes 37 Yes Yes 12 Yes Yes 25 Yes Yes 38 Yes Yes 13 No Yes 26 Yes Yes Next, we create four other cyber-attack scenarios in which we change multiple lines simultaneously. Since multiple lines can be randomly selected, we run 200 replications for each attack scenario. It is assumed that the attacker may change the origin/destination of a transmission line, Hk;l, the normal capacity of a transmission line, Fmaxl or the line availability status. We assume that the attacker may decrease the capacity of a transmission line to 10% of its original value. It is important to mention that under certain cyber-attacks the DC-OPF cannot obtain a feasible solution. Since the operator should be alerted to this condition, this cyber-attack scenario is considered detected. The cyber-attack scenarios are summarized in TABLE 1.9. The performances of the detection algorithms ANN-1 and ANN-2 are given in TABLE 1.10. The results show that the reduced ANN-2 which eliminates redundant information performs slightly better than ANN-1 for scenarios A, B and C. In scenario D the performances of both algorithms are comparable. Hence, the reduced ANN is chosen for experiments with the AC-OPF based detection algorithm. 25 Table 1.9: Cyber-attack scenarios Scenarios Line A B C D 1st ModifyOrigin ModifyOrigin ModifyCapacity ModifyOrigin 2nd ModifyDestination ModifyCapacity ModifyCapacity ModifyOrigin 3rd Remove ModifyCapacity ModifyOrigin Table 1.10: DC-OPF detection results in di erent cyber-attack scenarios %Detection Scenario ANN-1DC-OPF Model ANN-2DC-OPF Model A 92.5 97.5 B 95.0 97.0 C 95.5 97.5 D 99.5 98.5 1.6.2 AC-OPF Detection Algorithm The AC-OPF detection algorithm, ANN-3, is similar to the DC-OPF detection algo- rithm except that the AC-OPF algorithm takes the real AC power ows obtained from AC-OPF as inputs to the ANN and returns the estimate of the real AC power ows of the next period. The same process used for DC-OPF is used to determine the best design for ANN-3. TABLE 1.11 shows the best con gurations of ANN-3. The rst con guration is selected. 26 Table 1.11: Con gurations of ANN-3 Number of Neurons The Sum of Squared Errors Con guration Layer 1 Layer 2 Layer 3 Min Mean Max Standard Deviation 1 13 11 8 0.0001 0.0163 0.1035 0.0223 2 12 8 11 0.0001 0.0204 0.0822 0.0236 3 13 11 12 0.0003 0.0182 0.1337 0.0262 4 9 8 12 0.0005 0.0183 0.1269 0.0261 5 9 9 12 0.0003 0.0165 0.1698 0.0288 For ANN-3, the scale and shape parameters are estimated to be 0.013 and 0.74 respec- tively. Assuming that a 2.5% false alarm rate is acceptable, the 97.5 percentile is considered as the threshold value, i.e. 3 = 0:078. Similarly, to study the performance of the AC-OPF detection algorithm under normal operation conditions, we run our simulation for hours 338 to hour 744. Figure 1.12 shows the time series of the sum of the squared errors. Eight false alarms are observed in Figure 1.12. To study the performance of ANN-3 in detecting anomalies, the same attack scenarios mentioned in TABLE 1.12 are used. The detection results are provided in TABLE 1.12. The results show that the algorithm performance is highly satisfactory. The percentage of detection is greater than 95% under all scenarios. Table 1.12: AC-OPF detection results in di erent cyber-attack scenarios Scenario % Detection of ANN-3 (AC-OPF Model) A 95.0 B 95.0 C 98.5 D 99.0 27 0 50 100 150 200 250 300 350 400 4500 0.5 1 1.5 Time (Hours) Sum of Squared Error Threshold Figure 1.12: Plot of the ANN-3 false alarms 1.6.3 Computational Time Analysis All experiments are performed on a 64-bit laptop with an Intel Core i5 2.4GHz processor and 4GB RAM. The training time of the arti cial neural network and execution time of 100 replications of anomaly detections are given in TABLE 1.12. Table 1.13: Training and Detection Computational Times Arti cial Neural Network Training Time (s) Detection Time (s) ANN1 35.19 7.05 ANN2 17.53 10.55 ANN3 10.77 24.15 28 1.7 Conclusions We have illustrated that cyber-attacks could be as dangerous as physical attacks to the power grid since they could cause major physical losses and damages. Current data cleansing methods are not powerful enough to detect all cyber-attacks. We have proposed an algorithm which uses arti cial neural networks to detect cyber-attacks against the transmission network data of an electric power system, which is an additional security measure to the traditional state estimation software. An alarm from our algorithm would indicate that a parameter was changed deliberately. We have used network observability rules to reduce the size of the inputs to the ANN while maintaining the critical variables. We have tested our algorithm in two software mod- ules, DC-OPF and AC-OPF. We have simulated cyber-attacks by changing the parameters of components and transmission lines. The algorithm was able to detect 92 to 99.5% of the introduced anomalies with a small number of false alarms. The detection capability of the algorithm depends on how the altered parameters change transmission power ows. Although the algorithm is e ective on detecting the presence of an anomaly in the system, it cannot identify, locate, or eliminate the anomaly. The main obstacle for using the algorithm on a much larger power system is the computer processing time, which is highly depending on the number of input/output variables of the ANN. Fortunately, the ANN approach accepts parallel computing and can be easily implemented on a computer with multiprocessors. Certainly, more research needs to be done on studying the performance and the level of scalability of the proposed approach. However, the results from the 24-bus power system have shown that the algorithm is a promising tool for adding an extra level of cyber-security to a power system. 29 References [1] S. Halilcevic, F. Gubina, and A. F. Gubina, \Prediction of power system security levels," IEEE Transactions on Power Systems, vol. 24, pp. 368{377, 2009. [2] I. L.G.Pearson, \Smart grid cybersecurity for europe," Energy Policy, vol. 39, pp. 5211{ 5218, 2011. [3] C. W. Ten, C. Liu, and G. Manimaran, \Vulnerability assessment of cybersecurity for scada systems," IEEE Transactions on Power Systems, vol. 23, pp. 1836{46, 2008. [4] F. C. Schweppe and J. Wildes, \Power system static-state estimation, part i-iii," IEEE Transactions on Powe, vol. PAS-89, pp. 120{135, 1970. [5] L. Xie, Y. Mo, and B. Sinopoli, \Integrity data attacks in power market operations," IEEE Transactions on Smart Grid, vol. 2, pp. 659{666, 2011. [6] V. H. Quintana, A. Simoes-Costa, and M. Mier, \Bad data detection and identi ca- tion techniques using estimation orthogonal methods," IEEE Transactions on Power Apparatus and Systems, vol. PAS-101, pp. 3356{3364, 1982. [7] H. M. Merrill and F. C. Schweppe, \Bad data suppression in power system static state estimation," IEEE Transactions on Power Apparatus and Systems, vol. PAS-90, pp. 2718{2725, 1971. [8] E. Handschin, F. C. Schweppe, J. Kohlas, and A. Fiechter, \Bad data analysis for power system state estimation," IEEE Transactions on Power Apparatus and Systems, vol. 94, pp. 329{337, 1975. 30 [9] J. Valenzuela, J. Wang, and N. Bissinger, \Real-time intrusion detection in power system operations," IEEE Transactions on Power Systems, vol. PP, no. 99, p. 1, 2012. [10] L. Xie, Y. Mo, and B. Sinopoli, \False data injection attacks in electricity markets," in First IEEE International Conference on Smart Grid Communications (SmartGrid- Comm), Gaithersburg, MD, USA, 2010, pp. 226{231. [11] Y. Liu, M. K. Reiter, and P. Ning, \False data injection attacks against state estimation in electric power grids," in 16th ACM Conference on Computer and Communication Security, Chicago, IL, USA, 2009. [12] A. Abur and A. G. Exposito, Power System State Estimation: Theory and Implemen- tation. Marcel Dekker Inc., 2004. [13] C. W. Ten, G. Manimaran, and C. C. Liu, \Cybersecurity for critical infrastructures: Attack and defense modeling," IEEE Transactions on Systems, Man and Cybernetics, Part A: Systems and Humans, vol. 40, pp. 853{865, 2010. [14] Z. C. Z and J. C. Maun, \Arti cial neural network approach to single-ended fault locator for transmission lines," IEEE Transactions on Power System, vol. 15, pp. 370{375, 2000. [15] A. G. Bakirtzis, J. B. Theocharis, S. J. Kiartzis, and K. J. Satsios, \Short term load forecasting using fuzzy neural networks," IEEE Transactions on Power System, vol. 10, pp. 1518{1524, 1995. [16] K. Meng, Z. Y. Dong, D. H. Wang, and K. P. Wong, \A self-adaptive RBF neural network classi er for transformer fault analysis," IEEE Transactions on Power System, vol. 25, pp. 1350{1360, 2010. [17] A. J. Wood and B. F. Wollenberg, Power Generation, Operation and Control. John Wiley & Sons, 1996. 31 [18] R. D. Zimmerman, C. E. Murillo-Sanchez, and R. J. Thomas, \MATPOWER: Steady- state operations, planning, and analysis tools for power systems research and education," IEEE Transactions on Power Systems, vol. 26, pp. 12{19, 2011. [19] L. Xie, Y. Mo, and B. Sinopoli, \Integrity data attacks in power market operations," IEEE Transactions on the Smart Grid, vol. 2, pp. 659{666, 2011. [20] M. H. Beale, M. T. Hagan, and H. B. Demuth, Neural Network Toolbox User Guide?s. The MathWorks, Inc, 2011. [21] J. Peng, Y. Sun, and H. F. Wang, \Optimal pmu placement for full network observability using tabu search algorithm," Electrical Power and Energy Systems, vol. 28, pp. 223{ 231, 2006. [22] Subcommittee, \IEEE reliability test system," IEEE Transactions on Power Apparatus and Systems, vol. 98, pp. 2047{2054, 1979. [23] www.pjm.com. 32 Chapter 2 A Probabilistic Risk Mitigation Model for Cyber-Attacks to PMU Networks 2.1 Abstract The power grid is becoming more dependent on information and communication tech- nologies. Complex networks of advanced sensors such as PMUs are used to collect real time data to improve the observability of the power system. Recent studies have shown that the power grid has signi cant cyber vulnerabilities which could increase when PMUs are used extensively. Therefore, recognizing and responding to vulnerabilities are critical to the se- curity of the power grid. This paper proposes a risk mitigation model for optimal response to cyber-attacks to PMU networks. We model the optimal response action as a mixed in- teger linear programming (MILP) problem to prevent propagation of the cyber-attacks and maintain the observability of the power system. keywords: Cyber-attack, Phasor Measurement Units, Cyber-security, Networks, Observ- ability Nomenclature Sets and Indices : Set of buses : Set of buses equipped with PMUs : Set of PMUs detected as compromised : Set of buses with conventional devices : Set of branches with conventional devices i;j;k: Indices of buses 33 Constants M: Number of detected compromised PMUs Hi;j: Connectivity between buses i and j Ti: Threshold threat level (between 0 and 1) of PMUi Dij: Nodal distance between PMUi and PMUj t: Time that a propagation attempt takes m: Number of t that the system operator takes to respond to the cyber-attack Decision Variables xj: : Binary decision variable which equals 1 if PMUj is kept connected to the network, and 0 otherwise ?i: Observability number, number of times that bus i is observed, which is 1 if bus i is observable yi;j: Binary variable which equals 1 if the measurement from the conventional device at bus j is assigned to compute the unknown voltage phasor of bus i, and 0 otherwise. i;j: Binary variable which equals 1 if the measurement from the conventional device at trans- mission line i j is assigned to compute the unknown voltage phasor of bus i, and 0 otherwise. Random Variables Aj(t): random variable which equals 1 if PMUj is attacked, and 0 otherwise Probabilities ij: The probability that PMUj is attacked through PMUi : The probability that an attack propagates through a router : The probability that an attack propagates to a PMU through a router j(t): The threat level of PMUj at time t 34 2.2 Introduction The power grid is increasingly dependent on information and communication technolo- gies due to the integration of intelligent measurement devices such as phasor measurement units (PMU). Smart grid investment grants and demonstration project investments have sig- ni cantly accelerated the pace of phasor technology deployment [1]. According to [2], PMUs will ultimately replace conventional devices. PMUs can measure in real time synchronized phasors of bus voltages and currents for better observability of the power grid [3]. Synchro- nization is achieved by timing signals from the Global Positioning System (GPS) satellite. A PMU takes about 30 to 120 measurements per second and sends its measurements to a phasor data concentrator (PDC) through a wireless communication network based on the NASPInet architecture [4, 5]. In the NASPInet architecture, PMUs are connected to an IP-based communication network like an Intranet. Although the communication network is dedicated Intranet and isolated from public networks, it is not immune to cyber-attacks [6]. Currently, PMUs transmit their measurements to a pre-de ned PDC in a hierarchical man- ner by using the IP Unicast routing protocol. However, hierarchical architectures su er from drawbacks such as delays of messages. A technical report from CISCO proposes that PMUs should send measurements using the IP Multicast routing protocol [7]. Under this proto- col, a PMU is directly connected to a router and sends out data packets to pre-con gured destinations. The list of these predetermined destinations can be further manipulated by a cyber-attacker to propagate the cyber-attack to other PMUs. The propagation of the cyber- attacks in shared communication networks has been also studied in [8{13]. Furthermore, it has been reported that the communication network shows poor network security and in- su cient software security [14]. Moreover, the authors in [15] studied the spoo ng attack as an optimization problem to maximize the PMU?s receiver clock o set before and after the attack. The authors in [16] mentioned that there is no available defense against GPS spoo ng which is a threat to critical infrastructure applications such as PMUs that rely on the publicly known civilian GPS signal. Under these conditions, cyber-attackers could gain 35 access to the PMU communication network, inject false measurement data and propagate their cyber-attack to the other PMUs to endanger the reliability of the power grid. Dealing with erroneous data has been a concern of state estimation programs since their inception in the late nineteen-sixties [17]. It is shown in [18] that sophisticated attacks may not be detected by conventional state estimation algorithms. Thus, the detection problem is considered to be challenging and several algorithms have been proposed. The authors in [19] presented a Principal Component Analysis (PCA) based approach to detect cyber- attacks in the optimal power ow (OPF) module. The authors in [20] showed how data attacks can endanger the physical structure of the power grid and developed a detection algorithm using the arti cial neural networks. The authors in [21] have proposed new points of view to enrich the detection solutions such as modeling the dynamics of attacker versus defender. The authors in [22] developed greedy algorithms to obtain perfect protection and partial protection against stealth attacks given a limited budget for protection. In [23], the authors used the generalized likelihood ratio test to develop a computationally e cient detection algorithm where the cyber-attacker uses a graph theoretic approach to launch stealthy malicious data attacks. After the detection, actions need to be taken to prevent the propagation of the cyber-attack. The authors of [13] have formulated an optimization model to avoid propagation of the cyber-attacks in open-science computer network where collaboration and communication exist among network sites. The optimization model is a mixed integer linear programming problem, which determines the sites to be disconnected from the network to maximize the number of users connected to the network resources. The decision is constrained to keep the threat levels of the sites below a certain threshold value. However, such a model cannot be directly applied to a PMU network as that model is focused more on maximizing the available connections on the network, which is not a priority for PMUs. Additionally, observability is not an issue in open-science computer networks. 36 In this paper, we propose an optimal response model to cyber-attacks to PMU networks where the state estimation principally relies on PMUs. Our model minimizes the threat levels by disabling known compromised PMUs and PMUs that are likely to be compromised due to the propagation of the cyber-attacks, while keeping the power system observable. Here, a threat level represents the probability of a PMU being contaminated at a certain time. We rst use a probabilistic model to estimate the threat levels for PMUs and then formulate the optimal response as a mixed-integer linear programming problem. Here, a response stands for disconnecting the contaminated PMU buses to ensure the resultant PMU network is secure and the grid is observable. The rest of this chapter is organized as follows: Section 2.3 describes threat level es- timation. Section 2.4 explains our proposed optimization model. Experimental results are provided in Section 2.5 and the conclusions are reported in Section 2.6. 2.3 Estimating the Threat Levels In this section, we calculate the threat levels of the uncompromised PMUs if an intruder were to attacked one or more PMUs. Propagation of the cyber-attack has been studied on the other networks. As a case in point, the authors in [24] and [25] considered worm prop- agation in mobile ad-hoc networks and energy meters in a secondary distribution network, respectively. Similarly in a PMU network, the intruder controls the attacked PMUs which can be used to transmit false measurements. Moreover, the intruder can use the communica- tion links between PMUs to disseminate the attack to the other PMUs via the compromised PMUs. If the attack propagates to more PMUs, it could jeopardize the observability of the power system even further. The attack could propagate via all routers between the compromised and uncompromised PMUs to contaminate the uncompromised PMUs. Let us assume at time t = 0, M PMUs are detected to be compromised while the remaining PMUs are uncompromised. It takes time t to disable the compromised PMUs from the communication network. Naturally, their measurements will no longer be used for 37 the state estimation software. Disabling PMUs can be done automatically by the detection software or manually by the system operator. During this time, there is a chance that the attack could have been propagated to uncompromised PMUs and the detection software has not detected them yet. The reason is that the detection software cannot detect at 100% e ciency [21]. The cyber-attack can propagate to uncompromised PMUs through a path of interconnected routers. If the cyber-attack successfully breaks into all routers between the compromised and uncompromised PMUs, it is likely to contaminate the uncompromised PMUs as well. Moreover, these new compromised PMUs can further infect other PMUs, and so forth. We represent by the probability, j(t), the likelihood of a PMU being compromised (also called threat level) at time t. Notice that the threat levels increase over time as long as the network still contains compromised PMUs. It takes a time period of m t for the system operator to run the optimization model to obtain the optimal response and con rm that the alarm is not a false alarm. Thus, certain PMUs, determined by the optimization model, are disabled at time (m+ 1) t. At time t = 0, equations (2.1) and (2.2) hold. Pr(Aj(0) = 1) = 1 8j2 (2.1) Pr(Aj(0) = 1) = 0 8j =2 (2.2) By time t, all compromised PMUs are disabled. However, due to the possible prop- agation of the cyber-attack, there is a chance that the remaining PMUs could have been compromised and undetected. Therefore at time t, equations (2.3) and (2.4) hold. Pr(Aj( t) = 1) = 0 8j2 (2.3) Pr(Aj( t) = 1) = 1 Y i2 (1 ij) 8j =2 (2.4) 38 where ij is the probability that the attack propagates from compromised PMUi to an uncompromised PMUj during the time t, and it is given by equation (2.5). ij = Dij i2 ; j =2 (2.5) In equation (2.5), is the probability that an attack propagates through a router, is the probability that an attack propagates to another PMU and Dij, called nodal distance, is the minimum number of routers that connect PMUi and PMUj on the communication network. It is likely that there are multiple shortest paths between two PMUs in a large communication network. In this case, the value of ij would increase and would be given by equation (2.6). ij = 1 (1 Dij)Nij (2.6) Where Nij is the number of shortest paths between PMUi and PMUj. There may be other paths in addition to the shortest paths. Considering that it may not be simple to nd all paths between every two PMUs, and that the probability ij decreases exponentially when the distance increases, we use equation (2.6) to estimate the probability that the attack propagates from one PMU to another. This function indicates that the propagation becomes less probable when the nodal distance between the two PMUs becomes larger. We have stated that it takes time t to disable the compromised PMUs. If, during m t, the system operator concludes that the alarm is false, the disabled PMUs would be enabled again. Otherwise, the operator would begin disabling the PMUs as determined by the optimization model discussed in the next section at m t and all contaminated PMUs would be disabled by time (m + 1) t. Thus, we need to calculate the threat levels at time (m+1) t. Since the threat levels are calculated in an iterative process, we need to calculate 39 the threat levels at time 2 t, 3 t, :::, n t. Equations (2.7) and (2.8) hold at time 2 t. Pr(Aj(2 t) = 1) = 0 8j2 (2.7) Pr(Aj(2 t) = 1) = Prf(Aj(2 t) = 1jAj( t) = 0g PrfAj( t) = 0g + Prf(Aj(2 t) = 1jAj( t) = 1g PrfAj( t) = 1)g 8j =2 (2.8) In equation (2.8), we have expressions for all terms except for Prf(Aj(2 t) = 1jAj( t) = 0g, which can be obtained from equation (2.9). PrfAj(2 t) = 1jAj( t) = 0g= 1 PrfAj(2 t) = 0jAj( t) = 0g= 1 Y k;j=2 k6=j (1 Pr(Ak( t) = 1) kj) j;k2 (2.9) We denote Pr(Aj(t) = 1) by j(t). Therefore, we can rewrite equation (2.8) as equation (2.10) in terms of threat levels. j(2 t) = 0 BB B@1 Y k;j=2 k6=j (1 k( t) kj) 1 CC CA (1 j( t)) + 1 j( t) j;k2 (2.10) It can be shown that equation (2.11) is true when n 2 and can be used to calculate threat levels for time 2 t and further. j(n t) = 0 BB B@1 Y k;j=2 k6=j (1 k((n 1) t) kj) 1 CC CA 1 j((n 1) t) ! + j((n 1) t) = 1 Y k;j=2 k6=j (1 k((n 1) t) kj) ! (1 j((n 1) t) j;k2 (2.11) 40 If another cyber-attack is detected by time (m+1) t, the fundamental change due to the new attack is the change on set , meaning that the set of the detected compromised PMUs is adjusted to include the new detected PMUs. In this case, the threat levels after the new cyber-attack are re-calculated by setting the initial threat levels to the threat levels of the previous attack. This is, if the second detection occurs at time S = n t and 1 n (m+1): j;new attack(0) = j;old attack(n t) 8j =2 (2.12) In the next step, equation (2.11) will be used to update the threat levels after time S = n t. 2.4 Response Model The response to cyber-attacks to a PMU network is modeled using mixed integer linear programming. The objective function of the model is the minimization of the maximum threat level of all connected PMUs at time (m + 2) t, which is one t after disabling the PMUs determined by the model. The threat levels from time t = 0 to t = (m + 2) t are summarized in equations (2.13)-(2.19). j(0) = 1 8j2 (2.13) j(0) = 0 8j =2 (2.14) j( t) = 0 8j2 (2.15) j( t) = 1 Y i2 (1 ij) ! 8j =2 (2.16) j(n t) = 0 8j2 ; 2 n m+ 2 (2.17) j(n t) = 1 Y k;j=2 k6=j (1 k((n 1) t) kj) ! 41 (1 j((n 1) t) 8j;k =2 ; 2 n m+ 1 (2.18) j (m+ 2) t = 1 Y k;j=2 k6=j (1 k((m+ 1) t) kjxk) ! (1 j((m+ 1) t)) 8j;k =2 (2.19) Notice the presence of the binary decision variable xk in equation (2.19), which equals 1 if the PMUk is kept connected to the network, and 0 otherwise. It should be also noted that threat levels from time t = 0 to t = (m+ 1) t are all assumed to be constants, and the system operators cannot decrease them due to the physical constraints such as control and communication delays in disabling PMUs from the network. However, disabling of suspicious PMUs occurs at time (m + 1) t, which decreases the threat levels of the remaining PMUs at time (m+2) t. In fact, the response optimization model determines which PMUs should be disabled such that the threat levels at time (m+ 2) t are minimized. To be able to solve the response model more e ciently, we reformulate equation (2.19) by an equivalent linear equation (2.20). lnf1 j (m+ 2) t g= X j;k=2 k6=j ln(1 k (m+ 1) t kj)xk + lnf1 j (m+ 1) t g 8j =2 (2.20) In obtaining equation (2.20), we have used the following equality where K is a constant (notice again that xi is a binary variable). ln (1 Kxi) = xi ln (1 K) (2.21) 42 The objective function of the response model is the minimization of the maximum threat of all connected PMUs at time (m+ 2) t. Z = minx maxj ( j (m+ 2) t xj) 8j =2 (2.22) Since j((m+2) t) and consequently the objective function are not linear, we used its equiva- lent function given in equation (2.23), which can be reformulated linearly in equations (2.24)-(2.29). Z = minx maxj ( ln[1 j (m+ 2) t ] xj ) 8j =2 (2.23) Equation (2.23) can be represented by the following set of equations (2.24)-(2.29). Z = minY (2.24) Subject to: wj xj 8j =2 (2.25) wj ln [1 j (m+ 2) t ] 8j =2 (2.26) wj ln [1 j (m+ 2) t ] (1 xj) 8j =2 (2.27) wj 0 8j =2 (2.28) wj Y 8j =2 (2.29) In the response model, we used equation (2.20) in equations (2.26)-(2.27) to obtain equivalent linear equations. Hence, equations (2.26)-(2.27) can be represented as equations (2.30)-(2.31), respectively. wj X j;k=2 k6=j ln(1 k((m+ 1) t) kj)xk lnf1 j((m+ 1) t)g 8j =2 (2.30) 43 wj X j;k=2 k6=j ln(1 k((m+ 1) t) kj)xk lnf1 j((m+ 1) t)g (1 xj) 8j =2 (2.31) We add equation (2.32) to keep PMUs with threat levels less than a threshold value connected to the network. j((m+ 2) t) >Tj xj 8j =2 (2.32) which can be represented by: ln[1 j (m+ 2) t ] < ln[1 Tj +xj] 8j =2 (2.33) We can reformulate the right-hand side to a linear equivalent equation as: ln[1 j (m+ 2) t ] < (1 xj) ln(1 Tj) +xj ln(2 Tj) 8j =2 (2.34) Using equation (2.20), equation (2.34) can be represented linearly as equation (2.35). X j;k=2 k6=j ln(1 k((m+ 1) t) kj)xk + lnf1 j((m+ 1) t)g < (1 xj) ln(1 Tj) +xj ln(2 Tj) 8j =2 (2.35) Disabling PMUs may a ect the observability of the system, so there should be a set of constraints to ensure full observability. To have a full observable power system, the voltage phasors (state variables) of all buses need to be known. The values of the voltage phasors allow the calculation of all other variables such as current, real power, and reactive power. 44 Hence, the observability function is given in equation (3.2) [26]. ?i = X j2 Hi;jxj + X j2 Hi;jyi;j + X fi;jg2 Hi;j i;j 8i2 (2.36) Where Hi;j shows the connectivity between buses i and j. The rst term of the right hand side, P j2 Hi;jxj, provides observability of the buses through the remaining connected PMUs; the second term, P j2 Hi;jyi;j, calculates the observability of the buses through conventional devices installed at buses, and the third term, P fi;jg2 Hi;j i;j, represents the observability through conventional devices installed on branches. To guarantee the full observability, the observability variable value should be greater than or equal to one as given in (2.37). ?i 1 8i2 (2.37) Conventional devices are used to observe a group of buses when they are not observable by PMUs direct measurements. In such cases, a system of equations is solved to obtain the unknown state variables. To guarantee the solvability of the system of equations, we need to ensure that each conventional measurement is tied to one state variable [26]. Hence, equations (2.38)-(2.39) need to be met. Equation (2.38) ensures that a conventional voltage measurement is tied to one state variable and equation (2.39) ensures that a conventional current measurement is tied to one state variable. When the system is under attack, mea- surements of the conventional devices may be more reliable because the attacker needs to have complete knowledge of system topology and considerable resources in order to con- sistently manipulate conventional measurements. Thus, we use equation (2.40) to force the measurements of conventional devices to observe more buses. X i2 Hi;jyi;j = 1 8j2 (2.38) i;j + j;i = Hi;j 8fi;jg2 (2.39) 45 X j2 Hi;jyi;j + X fi;jg2 Hi;j i;j 1 8i2 (2.40) where i;j is a binary variable which equals 1 if the state variable of bus i is computed by a conventional device at branch i;j, and 0 otherwise. Finally, equation (2.41) determines that the decision variables are binary. xi; yi;j; i;j2f0;1g (2.41) To summarize, the objective function given in equation (2.22) determines the PMUs to be disabled from the network such that the maximum threat level of the PMUs still connected to the network is minimized. If the threat level of a PMU exceeds the threshold value, the objective function would disable that PMU from the network only if it does not a ect the observability. If the threat level of a PMU is less than the threshold value, the PMU is kept connected to the network. This is guaranteed by equation (2.32). Finally, equation (3.2) ensures that all buses are observable at least once by either PMUs or conventional devices. The proposed response model is given by the following mixed integer linear programming problem. Z = minY (2.42) Subject to the constraints given in: Equations (2.13)-(2.18) Equation (2.25) Equations (2.28)-(2.29) Equations (2.30)-(2.31) Equation (2.35) Equations (3.2)-(2.41) 46 Finally, it is important to notice that the cyber-attack could involve blocking the com- munication to and from the PMUs, which may prevent the disabling of the compromised PMUs after detection. This situation is less likely to occur because the intruders would want to be undetected and a fault in the communication network would alert the system operators immediately. Nevertheless, if the communication is blocked, the PMUs would be isolated from the communication network, which has the same e ect as disabling the compromised PMUs. The measurements would not be used by the state estimation software, and there would be no propagation of the cyber-attack to other PMUs. 2.5 Experimental Results We test the performance of our response optimization model on the 6-bus and 24-bus test systems introduced in [27] and [28], respectively. We use the small power system to describe the cyber-attack propagation problem and our methodology. We use the middle size power system to test our approach. This power system includes conventional devices. We intentionally consider a few conventional devices because we want to show the use of our optimization response model where PMUs are the main devices for power system state estimation. In our experiments, we assume that an attack to a PMU propagates through a router with probability = 0:05, and it e ectively compromises another PMU with probability = 0:05. For simplicity, we also assume that there is only one shortest path between PMUs. We set Nij = 1. We set t to 0.1 seconds, t = 0:1(s), and the threshold values to 0.005, Ti = 0:005. In a real case, experiments can be conducted to estimate the values of these parameters. Experimental network infrastructures such as the Virtual Network Infrastructure (VINI), X-Bone, and Violin can be used to run these experiments with real tra c and routing software [29]. 47 2.5.1 6-bus Test System The 6-bus test system includes six buses, three generators and eleven transmission lines. Figure 2.1 shows the 6-bus test system. B u s 4 B u s 5 B u s 1 B u s 6 B u s 3 B u s 2 L i n e 1 L i n e 2 Line 3 Line 4 L i n e 5 L i n e 6 Line 7 L i n e 8 Line 9 L i n e 1 0 L i n e 1 1 PMU PMU PMU PMU PMU Shared Communication Network Supervisory Controller and State Estimations PDC PDC Figure 2.1: The 6-bus Test System Buses 1, 2, 3, 4 and 6 are equipped with PMUs which make the system fully observable. Phasor measurements such as voltage magnitudes and voltage angles are transmitted to the PDCs by the ve PMUs through a digital communication network. The communication network is depicted in Figure 2.1 and it consists of interconnected routers. TABLE 2.1 shows the nodal distances between the PMUs installed in the 6-bus test system. Notice that the greater the nodal distance between two PMUs, the less likely that the cyber-attack can propagate from one PMU to another one. For this case study, we assume that at time t = 0 the system operator is informed that PMU1 and PMU3 have been attacked by a cyber-intruder. In Figure 2.2, we show the threat levels of the three uncompromised PMUs over time when the system operator does not disable the compromised PMUs from the network. Notice that the threat levels increase 48 Table 2.1: Nodal distances between PMUs in the 6-bus test system Bus 1 2 3 4 6 1 0 2 3 1 2 2 2 0 3 3 2 3 3 3 0 1 2 4 1 3 1 0 3 6 2 2 2 3 0 Table 2.2: Candidate responses for the 6-bus test system PMU CandidateResponse Status 1 2 3 4 5 6 7 8 x2 1 1 1 1 0 0 0 0 x4 1 0 0 1 1 0 1 0 x6 1 1 0 0 1 1 0 0 ThreatLevels 2 0.00013 0.00013 0.00013 0.00013 0 0 0 0 4 0.00500 0 0 0.00500 0.00500 0 0.00500 0 6 0.00025 0.00025 0 0 0.00025 0.00025 0 0 Max( j) 0.00500 0.00025 0.00013 0.00500 0.00500 0.00025 0.00500 0 Observaiblity Yes Yes Yes Yes Yes No No No non-linearly until all PMUs become compromised with probability 1. PMU4 is more at risk because it is closer to the compromised PMUs. The nodal distance from PMU3 to PMU4 is 1 while to PMU2 is 3. Figure 2.3 illustrates the e ect of the value of on the threat level of PMU2. To avoid propagation, the system operator should disable the compromised PMU1 and PMU3 and other PMUs which may be compromised because of the propagation. There are 0 0 .1 0 .2 0 .3 0 .4 0 .5 0 .6 0 .7 0 .8 0 .9 1 0 200 400 600 800 1 0 0 0 1 2 0 0 1 4 0 0 T h re a t L e v e l T i m e ( Se c o n d s ) Figure 2.2: Threat levels in case of no response action 49 0 0 .1 0 .2 0 .3 0 .4 0 .5 0 .6 0 .7 0 .8 0 .9 1 0 200 400 600 800 1 0 0 0 1 2 0 0 1 4 0 0 T h re a t L e v e l T i m e ( Se c o n d s ) Figure 2.3: E ect of on threat level of PMU2 eight possible choices which are shown in TABLE 2.2. The smallest threat levels can be obtained when all PMUs are disabled. However, this solution is not feasible since the power system would no longer be observable. The second candidate is to disable PMU4 and PMU6 but the threat level of PMU6 is less than the threshold value, T6 = 0:005, and therefore it should remain connected to the network. The third candidate is to disable PMU4. This action keeps the power network observable and minimizes the maximum threat level of all connected PMUs. In TABLE 2.3, we give the optimal solution, observability number of the buses obtained from equation (3.2), and the threat level of each connected PMU right after disabling PMU4. Table 2.3: Optimal response action in the 6-bus system PMU xj Observability number( i) Threat level( j(3 t)) 1 0 1 0 2 1 2 0.00013 3 0 2 0 4 0 1 0 5 NA 2 NA 6 1 2 0.00025 We show in Figure 2.4 the maximum threat level of the connected PMUs when just the compromised PMU1 and PMU3 are disabled and when PMU4 is also disabled from 50 the network. Notice that the threat levels still increase after the response. However, the reduction in threat levels is considerable if the optimal response action is taken. 51 0 .0 0 7 0 .0 0 9 0 .0 1 1 0 .0 1 3 0 .0 1 5 0 .0 1 7 0 .0 1 9 0 .0 2 1 0 .0 2 3 0 .0 2 5 0 0 .2 0 .4 0 .6 0 .8 1 1 .2 1 .4 1 .6 1 .8 2 T h r e a t L e v e l T i m e ( s ) D i s a b l i n g t h e d e t e c t e d c o m p r o m i s e d P M U s O p t im a l r e s p o n s e a c t io n Figure 2.4: Comparison of two potential responses We have assumed that the optimization results are obtained in 0.1 seconds. To see the e ect of greater computational time, we consider that the optimal results are available at times 0.1, 2, and 5 seconds. In all of these cases, the optimal solution is to disable PMU4. In Figure 2.5, we show the threat level for each case. Certainly, a shorter processing time is desired. 2 .4 9 E -0 4 2 .5 0 E -0 4 2 .5 1 E -0 4 2 .5 2 E -0 4 2 .5 3 E -0 4 2 .5 4 E -0 4 2 .5 5 E -0 4 2 .5 6 E -0 4 2 .5 7 E -0 4 2 .5 8 E -0 4 2 .5 9 E -0 4 2 .6 0 E -0 4 0 5 10 15 20 25 30 35 40 T h re a t L e v e l T i m e ( S e c on d s ) Figure 2.5: E ect of computational time on threat level of PMU6 52 2.5.2 24-bus Test System The proposed response model is tested using the IEEE 24-bus test system which consists of 38 transmission lines, 24 buses and 33 generators. The 24-bus test system is shown in Figure 2.6. We refer the readers to [28] and [30] for additional system data. L i n e 1 L i n e 2 Line 3 Line 4 L i n e 5 L i n e 6 Line 7 L i n e 8 L i n e 9 L i n e 1 0 L i n e 1 1 Line 12 Line 13 Line 14 L i n e 1 5 L i n e 1 6 L i n e 1 7 L i n e 1 8 L i n e 1 9 L i n e 2 0 Line 21 L i n e 2 2 Line 23 L i n e 2 4 L i n e s 2 5 , 2 6 Line 27 L i n e 2 8 Line 29 Line 30 L i n e 3 1 L i n e s 3 2 , 3 3 Lines 34,35 L i n e s 3 6 , 3 7 L i n e 3 8 B u s 8 B u s 9 B u s 1 0 B u s 7 B u s 2 2 B u s 1 1 B u s 2 3 B u s 6 B u s 5 B u s 1 B u s 2 4 Bus 21 Bus 20Bus 19 B u s 1 8 Bus 15 B u s 3 B u s 4 B u s 2 B u s 1 4 B u s 1 3 B u s 1 6 B u s 1 7 B u s 1 2 Figure 2.6: IEEE 24-bus Test System We use the optimal placement of seven PMUs given in [31] and add seven more PMUs and ve conventional measuring devices to increase the observability of the power system. The PMUs are located at buses 1, 2, 3, 4, 8, 10, 11, 15, 16, 17, 21, and 23 and the conventional units are located at buses 1, 19, 22 and transmission lines 9 and 17. The initial observability of the buses is given in TABLE 2.4. TABLE 2.5 gives PMUs? nodal distances randomly generated between 1 and 4. 53 Table 2.4: Initial observabilites of the 24-bus system Bus Observability Bus Observability Bus Observability Bus Observability 1 3 7 2 13 2 19 2 2 3 8 3 14 2 20 3 3 3 9 4 15 3 21 2 4 2 10 3 16 3 22 3 5 3 11 2 17 2 23 2 6 2 12 3 18 2 24 2 Table 2.5: Nodal distances between PMUs in the 24-bus test system PMU Locations 1 2 3 4 7 8 10 11 15 16 17 20 21 23 1 0 1 2 1 2 4 4 2 3 4 1 2 3 2 2 1 0 3 3 2 1 3 2 4 3 3 3 1 3 3 2 3 0 2 1 4 2 3 1 3 3 2 1 2 4 1 3 2 0 3 1 2 2 3 4 3 3 1 1 7 2 2 1 3 0 1 3 1 1 4 1 4 2 3 8 4 1 4 1 1 0 2 3 1 1 1 2 3 2 10 4 3 2 2 3 2 0 4 4 2 3 1 2 2 11 2 2 3 2 1 3 4 0 2 2 3 2 4 1 15 3 4 1 3 1 1 4 2 0 4 1 2 2 4 16 4 3 3 4 4 1 2 2 4 0 1 4 1 1 17 1 3 3 3 1 1 3 3 1 1 0 3 4 4 20 2 3 2 3 4 2 1 2 2 4 3 0 3 3 21 3 1 1 1 2 3 2 3 2 1 4 3 0 1 23 2 3 2 1 3 2 2 1 4 1 4 3 1 0 First, we assume that just one PMU is compromised at time zero. We consider two cases to show the propagation e ect of the cyber-attack. We assume: a)PMU7 is compromised and b) PMU20 is compromised at time zero. In both cases, it is assumed that the compromised PMU stays connected to the network (no response action). Figure 2.7 shows the threat levels of a selected group of PMUs, PMUs 4, 10 and 15, under case (a) while Figure 2.8 shows the threat level for the same PMUs under case (b). Notice that PMU10 is less a ected by a cyber-attack to PMU7 but it is considerably a ected if the attack occurs on PMU20. Next, we study the e ect of simultaneous cyber-attacks. We assume that PMUs 1, 3, 7 and 10 are compromised at time zero and they remain connected to the network (no response action). Figure 2.9 illustrates the threat level of PMU4 under a single attack to PMU7 and 54 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 10 20 30 40 50 60 70 80 T h r eat L evel T i m e (S ec on d s ) Figure 2.7: Threat levels for no response to compromised PMU7 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 10 20 30 40 50 60 70 80 T h r e at L e ve l T im e ( Sec ond s) Figure 2.8: Threat levels for no response to compromised PMU20 the multiple attacks. Notice the increase on the threat level of PMU4 under a multiple attack. To show the e ect of the computational time m t (time to obtain the optimal response) on the threat levels, we use the multiple attacks case. We change the value of m from 1 to 700, i.e. the time that the system operator takes to respond to the cyber-attack varies from 55 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 10 20 30 40 50 60 70 80 T h r e at L e ve l T im e ( Sec ond s) Figure 2.9: E ect of single and multiple attacks on threat levels 0.1 to 70 seconds. We obtain the optimal response for each value of m. The results are given in TABLE 2.6. Table 2.6: E ect of response time on the optimal solution in case of multiple attacks Case Response time (s) Disabled PMUs S1 m t 8:9 1, 3, 7, 10, 17 S2 9:0 m t 10:8 1, 3, 4, 7, 10, 17 S3 10:9 m t 12:7 1, 3, 4, 7, 10, 16, 17 S4 12:8 m t 70:0 1, 3, 4, 7, 10, 16, 17, 23 S5 m t 70:0 Not Feasible In Figure 2.10, we show the maximum threat level of the system (objective function of the optimization problem) when only the compromised PMUs are disabled from the network and when the optimal response is implemented. The gure shows the results under di erent response times. Notice that when the proposed model is used the maximum threat level is reduced for all response times considered. The reduction increases as the response time increases. 56 0.0 00 0 0.0 050 0.0 100 0.0 150 0.0 20 0 0.0 250 9 109 128 Maxim u m T h r e at Leve l m O p tim al re s p on s e Disc on n e c tin g on ly c om p r om ised P M Us Figure 2.10: Maximum threat levels for two potential responses 2.5.3 Dealing with Large Power Systems Since the main contribution of the paper is the formulation of the MILP and the mod- eling of the probabilistic threat levels, we demonstrated our approach using a small and a medium size power system. All experiments were performed on a 64-bit laptop with an Intel Core i5 2.4 GHz processor and 4GB RAM. The computation time consists of two compo- nents, the threat level calculation and the optimization time. For the experiment on the IEEE 24-bus test system, the threat level calculations using MATLAB R2012a took 6.95 seconds and the optimization using the optimization software LINGO 11.0 took 1 second. It is known that MILP solvers su er from the curse of dimensionality and therefore waiting to obtain the optimal solution for a larger power system can threaten the security of the power system more severely. However, considering that a single cyber-attack is more likely to occur, the proposed model can be run o ine for all possible single attacks. The system operator would already know the optimal response when a single attack occurs in the network. In TABLE 2.7, we show the optimal responses within 20 seconds for all single attacks to a PMU on the 24-bus test system. This approach, however, would not work for multiple attacks due to the numerous possible combinations of cyber-attacks to PMUs. In 57 Table 2.7: Optimal response for single cyber-attacks to the 24-bus test system Compromised PMU PMUs to be disabled Compromised PMU PMUs to be disabled 1 1, 4, 8, 16, 17 11 7, 11 2 1, 2, 8, 16, 21 15 4, 8, 15, 16, 17 3 3, 7, 16, 17 16 4, 8, 15, 16, 17, 23 4 1, 4, 8, 16, 21, 23 17 1, 2, 8, 15, 16, 17, 23 7 4, 7, 15, 16, 17 20 20 8 1, 2, 15, 16, 17, 23 21 1, 4, 8, 16, 21, 23 10 10 23 4, 8, 16, 21, 23 this case, the optimization software can be stopped at a predetermined time. A trade-o has to be made between the closeness of the obtained solution at the predetermined time to the optimal solution and the increase of the threat levels over that time. 2.5.4 Discussion We have developed a response optimization model to a cyber- attack to power systems that rely heavily on PMUs for the state estimation. Thus, the proposed model becomes more bene cial as the ratio of PMUs to conventional devices increases. It is expected that in the near future, as PMUs become widespread, the state estimation process will be more in uenced by PMU measurements. The proposed model disables PMUs restricted to secure system observability. However, pseudo measurements that are used to handle measurement unavailability in conventional state estimation could be also used to remove additional PMUs at the cost of loss of observ- ability. This feature can be incorporated to our model by adding a set of new constraints. Furthermore, if information from a detection scheme is available it can serve as an input to our model. For example, given the estimation on the likelihood of infection of all uncompromised PMUs at the time of the detection of the cyber-attack, the initial threat level of uncompromised PMUs, j(0), can be set to this value. For the sake of simplicity, in our experiments we assume that no information comes from the detection software and therefore all initial threat levels are set to zero. 58 2.6 Conclusions Observability of the power system is important for grid operation and control. An attacker could design a cyber-attack to PMUs that endangers the observability of the power system. We showed that in addition to disabling the compromised PMUs, other PMUs should also be disabled to reduce the probability of propagation of the cyber-attack. We developed a mixed integer linear programming model to determine the PMUs that should be disabled under the restriction that the remaining PMUs continue to maintain the observability of the power system. The model minimizes the maximum threat level of the PMUs that remain connected to the network. We have shown experimental results for two power systems. The results in both cases demonstrated a signi cant reduction in the propagation of the cyber-attacks when the solution obtained by the optimization model is implemented. References [1] (2010, October) Real-time application of synchrophasors for improving reliability. North American Electric Reliability Corporation (NERC). [Online]. Available: https://www.naspi.org/File.aspx? leID=519 [2] H. Lin, Y. Deng, S. Shukla, J. Thorp, and L. Mili, \Cyber security impacts on all-PMU state estimator - a case study on co-simulation platform GECO," in IEEE SmartGrid- Comm Symposium - Wide Area Protection and Control (WAMPAC), 2012. [3] D. Dua, S. Dambhare, R. K. Gajbhiye, and S. A. Soman, \Optimal multistage schedul- ing of PMU placement: An ILP approach," IEEE Transactions on Power Delivery, vol. 23, pp. 1812{1820, 2008. [4] W. Wang and Z. Lu, \Cyber security in the smart grid: Survey and challenges," Com- puter Networks, 2013. 59 [5] (2009, May) Data bus technical speci cations for north american syncro-phasor initiative network (NASPInet). North American Syncro-Phasor Initiative Network. [Online]. Available: https://www.naspi.org/File.aspx? leID=587 [6] H. Lin, Y. Deng, S. Shukla, J. Thorp, and L. Mili, \Cyber security impacts on all-PMU state estimator a case study on co-simulation platform GECO," in IEEE Third In- ternational Conference on Smart Grid Communications (SmartGridComm), November 2012, pp. 587{592. [7] (2012) PMU networking with IP multicast. CISCO Public. [Online]. Available: http://www.cisco.com/ [8] T. Liu, X. Guan, Q. Zheng, and Y. Qu, \A new worm exploiting IPv6 and IPv4-IPv6 dual-stack networks: experiment, modeling, simulation, and defense," IEEE Networks, vol. 23, pp. 22{29, 2009. [9] B. Sun, G. Yan, Y. Xiao, and T. A. Yang, \Self-propagating mal-packets in wireless sensor networks: Dynamics and defense implications," Ad Hoc Networks, vol. 7, pp. 1489{1500, 2009. [10] W. Xiaoming and L. Yingshu, \An improved SIR model for analyzing the dynamics of worm propagation in wireless sensor networks," Chinese Journal of Electronics, vol. 18, pp. 8{12, 2009. [11] C. C. Zou, D. Towsley, and W. Gong, \Modeling and simulation study of the propagation and defense of internet e-mail worms," IEEE Transactions on Dependable and Secure Computing, vol. 4, pp. 105{118, 2007. [12] B. Stephenson and B. Sikdar, \A quasi-species model for the propagation and contain- ment of polymorphic worms," IEEE Transactions on Computers, vol. 58, pp. 1289{1296, 2009. 60 [13] M. Altunay, S. Ley er, J. T. Linderoth, and Z. Xie, \Optimal response to attacks on the open science grid," Computer Networks: The International Journal of Computer and Telecommunications Networking, vol. 55, pp. 61{73, January 2011. [14] \NSTB assessments summary report: Common industrial control system cyber security weaknesses," Idaho National Laboratory (INL), Tech. Rep., 2010. [15] J. Xichen, Z. Jiangmeng, B. J. Harding, J. J. Makela, and A. D. Dominguez-Garcia, \Spoo ng GPS receiver clock o set of phasor measurement units," IEEE Transactions on Power Systems, vol. 28, pp. 3253{3262, 2013. [16] D. P. Shepard, T. E. Humphreys, and A. A. Fansler, \Evaluation of the vulnerability of phasor measurement units to GPS spoo ng attacks," International Journal of Critical Infrastructure Protection, vol. 5, pp. 146{153, 2012. [17] F. C. Schweppe, J. Wildes, and D. B. Rom, \Power system static state estimation, parts i,ii, and iii," IEEE Transactions on Power Apparatus and Systems, vol. PAS-89, pp. 120{135, 1970. [18] Y. Liu, P. Ning, and M. K. Reiter, \False data injection attacks against state estimation in electric power grids," in Proceedings of the 16th ACM conference on computer and communications security, 2009, pp. 21{32. [19] J. Valenzuela, J. Wang, and N. Bissinger, \Real-time intrusion detection in power system operations," IEEE Transactions on Power Systems, vol. PP, no. 99, p. 1, 2012. [20] S. Mousavian, J. Valenzuela, and J. Wang, \Real-time data reassurance in electrical power systems based on arti cial neural networks," Electric Power Systems Research, vol. 96, p. 285295, March 2013. 61 [21] C. Shuguang, H. Zhu, S. Kar, T. T. Kim, H. V. Poor, and A. Tajer, \Coordinated data- injection attack and detection in the smart grid: A detailed look at enriching detection solutions," IEEE Signal Processing Magazine, vol. 29, pp. 106{115, September 2012. [22] G. Dan and H. Sanberg, \Stealth attacks and protection schemes for state estimators in power systems," in Proceedings of the International IEEE Conference on Smart Grid Communications, Gaithersburg, MD, October 2010. [23] O. Kosut, L. Jia, R. J. Thomas, and L. Tong, \Malicious data attacks on smart grid state estimation: Attack strategies and countermeasures," in Proceedings of IEEE In- ternational Conference on Smart Grid Communications, Gaithersburg, MD, October 2010. [24] R. G. Cole, N. Phamdo, M. A. Rajab, and A. Terzis, \Requirements on worm mitigation technologies in MANETS," in Workshop on Principles of Advanced and Distributed Simulation, 2005, pp. 207{214. [25] Y. H. Chang, P. Jirutitijaroen, and C.-W. Ten, \A simulation model of cyber threats for energy metering devices in a secondary distribution network," in 5th International Conference on Critical Infrastructure (CRIS), Beijing, 2010, pp. 1{7. [26] S. Azizi, B. G. Gharehpetian, G. Hug-Glanzmann, and A. Dobakhshari, \Optimal in- tegration of phasor measurement units in power systems considering conventional mea- surements," IEEE Transactions on Smart Grid, vol. 4, pp. 1113{1121, 2012. [27] A. J. Wood and B. F. Wollenberg, Power Generation, Operation and Control. John Wiley & Sons, 1996. [28] Subcommittee, \IEEE reliability test system," IEEE Transactions on Power Apparatus and Systems, vol. 98, pp. 2047{2054, 1979. 62 [29] A. Bavier, N. Feamster, M. Huang, L. Peterson, and J. Rexford, \In VINI veritas: realistic and controlled network experimentation," in SIGCOMM ?06 Proceedings of the 2006 conference on Applications, technologies, architectures, and protocols for computer communications, vol. 36, October 2006, pp. 3{14. [30] R. D. Zimmerman, C. E. Murillo-Sanchez, and R. J. Thomas, \MATPOWER: Steady- state operations, planning, and analysis tools for power systems research and education," IEEE Transactions on Power Systems, vol. 26, pp. 12{19, 2011. [31] B. K. S. Roy, A. K. Sinha, and A. K. Pradhan, \An optimal PMU placement technique for power system observability," International Journal of Electrical Power & Energy Systems, vol. 42, pp. 71{77, 2012. 63 Chapter 3 Investment Decisions on Optimal Allocation of Phasor Measurement Units 3.1 Abstract Reliability of the electrical power systems necessitates wide-area monitoring and full observability of the power grid. Phasor measurement units (PMUs), as the state-of-the- art measurement devices, collect synchronized phasors of voltages and currents in real time and are utilized for full observability of the power systems. Due to budget restrictions and considerable cost of installing PMUs, it is not possible to equip all buses with PMUs. In this paper, we discuss the necessity of considering transmission switching and single contingencies in the optimal PMU placement problem. We show that considering transmission switching and single contingencies increases the PMU investment costs substantially and propose an integer linear programming model to determine the optimal PMU placement plan in two investment phases. In the rst phase, PMUs are installed to achieve full observability of the power grid whereas additional PMUs will be installed in the second phase to guarantee the N 1 observability of the power grid. In each phase, power grid observability is guaranteed for all resulting topologies of the power grid stem from transmission switching. Simulation results are provided on several IEEE test systems which show that our proposed approach is a promising enhancement to the methods available for the optimal placement of PMUs. keywords: Phasor measurement unit, Optimal placement, Network observability, Trans- mission switching, Integer linear programming 64 Nomenclature Sets and Indices : Set of buses : Set of power grid topologies stem from transmission switching i;j: Indices of buses : Index of topologies Constants i: Observability number of bus i i : Observability number of bus i after the rst phase of PMU placement under topology i : Observability number of bus i after the second phase of PMU placement under topology Cj: Cost of installing a PMU at bus j : The annual percentage of the change in PMU prices Hi;j: Binary parameter that equals to 1 if i = j or there is a transmission line between bus i and bus j, and 0 otherwise. H i;j: Binary parameter that equals to 1 if i = j or there is a transmission line between bus i and bus j under topology , and 0 otherwise. I: In ation-free interest rate K: Number of years between two phases of investment Decision Variables xj: Binary decision variable which is equal to 1 if bus j is equipped with a PMU in the rst phase, and 0 otherwise. yj: Binary decision variable which is equal to 1 if bus j is equipped with a PMU in the second phase, and 0 otherwise. 65 3.2 Introduction Wide-area monitoring and full network observability of the electrical power systems in real time was impractical until the emergence of phasor measurement units (PMUs). PMUs are power system devices that measure synchronized phasors of voltages and currents in real time [1]. Synchronization is achieved by timing signals from the global positioning system (GPS) satellite with the accuracy in the order of 1 microsecond. In the future, it is expected that the smart grid will consist of at least 10,000 PMUs each taking about 30 to 120 measurements per second [2]. To ensure the observability of the power system, voltage phasors of all buses should be either directly measured or computed from other measurements [3]. Two types of observabil- ity have been addressed, numerical and topological observability. A network is numerically observable if the measurement Jacobian is of full rank [4, 5]. These methods are computa- tionally extensive due to the iterative procedure of matrix manipulations [6]. Alternatively, topological observability considers interconnections of the buses and network observability rules to obtain the states vector of the power system. Unlike conventional measurement devices, a PMU can measure the current phasors of multiple lines and provide measure- ments to compute the voltage phasors of adjacent buses. Thus, there is no need, in terms of observability, to install a PMU at all buses. Recently, the problem concerning Optimal Placement of PMUs (OPP) has been studied by researchers. The OPP problem considers the minimum number of PMUs and their instal- lation locations that makes the power system observable. Many researchers have developed heuristic and meta-heuristic methods to solve the OPP problem. Chakrabarti and Kyri- akides in [7] applied an exhaustive binary search methodology to tackle the OPP problem and nd the associated locations of PMUs. An iterative three-stage heuristic method has been introduced in [8] where in the rst two stages less important and strategically impor- tant buses are determined, and the last stage returns the optimal solution using pruning operation. In [9] and [10], simulated annealing is used to solve the OPP problem. Other 66 meta-heuristic methods such as Tabu search in [11] and binary particle swarm optimization in [12] have been applied to nd the minimum number of PMUs required for full observabil- ity of the power system. Integer programming is used in [13] and [14] to nd the optimal placement of PMUs. The authors in [15] applied integer linear programming (ILP) to solve the OPP problem considering conventional measurement devices. Although the OPP prob- lem has been studied by many researchers, there are still certain practical aspects of the problem that need to be considered. In this paper, we consider new investment decisions on the placement of PMUs and deliver the following outcomes as our contributions. First, failure of any PMU or transmission line may a ect full observability of the power grid. Therefore, considering single contingencies in the allocation of PMUs is essential to meet the reliability requirements of the power system [15]. However, considering single contingencies increases the investment costs substantially, more than twice the initial costs in many cases. We refer the readers to Section 3.3.2 that illustrates this observation. Hence, we propose an integer linear programming problem that minimizes the total investment costs and determines the optimal placement of PMUs in two investment phases. In the rst investment phase, PMUs are installed to achieve full observability of the grid. In the second phase of investment, additional PMUs are placed in service to meet the N 1 reliability requirement, which assures the full observability of the power grid even in the case of single contingencies. Secondly, it has been discussed in the literature that transmission switching can provide additional economical advantages for a power system through changing its topology dur- ing operations [16{20]. Hence, transmission switching should be considered in the optimal placement of PMUs. Otherwise, it may put the observability of the power grid at risk. We refer the readers to Section 3.3.4 for a case study on how transmission switching may a ect power grid observability. 67 Finally, we propose an integer linear programming model for optimal placement of PMUs in two investment phases. The model minimizes the total investment costs and considers transmission switching in the optimal PMU allocations. The rest of this chapter is organized as follows: Section 3.3 describes the general optimal PMU placement model and our developed case studies. Section 3.4 discusses our proposed two-phase optimal PMU placement model. Section 3.5 provides some experimental results and Section 3.6 reports our conclusions. 3.3 PMU Placement Model The Optimal PMU placement problem is de ned as nding the installation location of PMUs required for the observability of the power system such that the total cost is minimized. Network observability rules can be used to avoid installing PMUs at all buses and reduce associated costs signi cantly. 3.3.1 Network Observability Rules The network observability rules for topological observability of the power system given in [12] are described here. 1. If a PMU is installed at bus i, voltage phasor of bus i and current phasors of all incident transmission lines to bus i are known (Figure 3.1). PM U I - M easured I - M easured V - M easured Figure 3.1: Network Observability Rule 1 68 2. If voltage phasor of one end of a transmission line and the current phasor of the transmission line are known, the voltage phasor of the other end of the transmission line can be calculated (Figure 3.2). V - M e a s u r e d I - M easured V - C a l c u l a t e d Figure 3.2: Network Observability Rule 2 3. If voltage phasors of both ends of a transmission line are known, the current phasor of the transmission line can be calculated (Figure 3.3). V - M e a s u r e d I - Calculated V - M e a s u r e d Figure 3.3: Network Observability Rule 3 Measurements obtained by Rule 1 are direct measurements. Rules 2 and 3 provide pseudo measurements. Zero-injection buses, which do not inject currents into the system, have the potential to reduce the number of required PMUs for observability of the power system, but we do not consider zero-injection buses in this paper. We refer the reader to [12] for more information on the network observability rules on zero-injection buses. 69 3.3.2 Single-phase Optimal PMU Placement Model The objective function of the optimal PMU placement problem is to minimize the in- vestment cost of PMUs. The objective function is given in equation (3.1). Z = minx X j2 Cjxj (3.1) A PMU at bus j can make all the adjacent buses observable by measuring the current phasors of the incident lines. Hence, the observability number of a bus is obtained by equation (3.2) using the network observability rules. i = X j2 Hi;jxj 8i2 (3.2) We refer to the number of PMUs that make bus i observable, i, as the observability number of bus i. To have a topologically observable power system, the observability number of all buses should be greater than or equal to 1. Therefore, we write equation (3.3). i 1 8i2 (3.3) Equations (3.1)-(3.3) assure the full observability of the power grid with minimum in- vestment costs. Without loss of generality, it is common to assume that Cj = 1. Thus, the investment cost is interpreted in terms of number of PMUs. Moreover, single contingencies such as failure of a PMU or a transmission line should be taken into account in the placement of PMUs. Therefore, each bus in the power grid should be observable by at least two PMUs to make sure that any single contingency does not a ect the full observability of the network. Mathematically speaking, we write equation (3.4) to ensure that the optimal placement is resilient against any single contingencies. i 2 8i2 (3.4) 70 We denote the ILP model given in Equations (3.1)-(3.4) by single-phase optimal PMU placement model since it is assumed that all PMUs are installed together. The single-phase optimal PMU placement problem with and without considering single contingencies have been solved by di erent heuristic, meta-heuristic and ILP methods. Table 3.1 shows the minimum number of PMUs required to make di erent test systems observable with and without considering single contingencies [8,21]. Table 3.1: Optimal number of PMUs for full observability Power System Minimum Observability N 1 Observability IEEE 14-bus 4 9 IEEE 24-bus 7 14 IEEE 30-bus 10 21 IEEE 39-bus 13 28 IEEE 57-bus 17 33 IEEE 118-bus 32 68 Considering the results provided in Table 3.1, it can be inferred that achieving the N 1 observability placement costs almost twice the minimum observability placement. Regarding the substantial capital cost of installing PMUs, a utility company may prefer to install PMUs in two phases. In the rst phase, PMUs are installed to make the power grid fully observable by PMUs and postpone the N 1 observability placement to the second phase. Then, extra PMUs will be installed in the second phase to achieve N 1 observability. However, installing PMUs in the rst phase should be done wisely to avoid any unnecessary additional investment in the second phase. In the next subsection, we discuss how PMU investment of the second phase may be a ected by the PMU placement of the rst phase. 71 3.3.3 Case Study I: Two-phase Investment Approach To demonstrate the two-phase investment approach, we use the IEEE 57-bus system which consists of 57 buses, 80 transmission lines and 7 generators. The IEEE 57-bus test system is depicted in Figure 3.4. We refer the readers to [22] and [23] for additional system data. Figure 3.4: IEEE 57-bus test system B 5 B 4 B 3 B 2 B 1 B 1 6 B 6 B18 B 1 4 B 4 5 B 1 2 B19 B20 B21 B15 B17 B 1 3 B 2 6 B48 B47 B46 B44B23 B24 B22 B38 B37 B 5 0 B 4 9 B 3 9 B 5 7 B 4 0 B 5 6 B 4 1 B 1 1 B 1 0 B 5 1 B 9 B42 B43 B25 B27 B30 B33 B35 B 2 8 B 3 1 B 3 2 B 3 4 B 7 B 8 B29 B52 B53 B54 B55 B36 From Table 3.1, we know that IEEE 57-bus system requires at least 17 PMUs for full observability and 33 PMUs for N 1 observability. Table 3.2 shows alternative optimal placements for full observability of the IEEE 57-bus test system. Each alternative may be selected in the rst phase since they all require the minimum investment cost. However, it should be noticed that installing PMUs in the rst phase creates new constraints for the optimal allocation of PMUs in the second phase. Therefore, it is likely that the total investment cost is not completely minimized if two phases are not considered together. To achieve N 1 observability in phase-2, we proceed and use the eight alternative solutions (P1 P8) given in Table 3.2. For each alternative solution, the number of PMUs and their locations to achieve the N 1 observability at the second-phase are given in Table 3.3. Notice that seven alternatives ultimately required more number of PMUs. The best investment 72 plan is given by U8 which consists of installing 17 PMUs during the rst phase and 16 PMUs during the second phase. The advantage of the two-phase plan is that the utility can decide whether to install all PMUs in one or in two phases while avoiding potential unnecessary investments. Table 3.2: Alternative placements of the IEEE 57-bus system Placement PMU Locations P1 3, 6, 12, 15, 19, 22, 25, 27, 32,36, 39, 41, 45, 47, 50, 52, 55 P2 2, 6, 12, 19, 22, 25, 27, 32, 36,39, 41, 45, 46, 49, 51, 52, 55 P3 2, 6, 12, 19, 22, 25, 27, 32, 36,39, 41, 45, 46, 49, 50, 52, 55 P4 2, 6, 12, 19, 22, 25, 27, 29, 32,36, 41, 45, 46, 49, 51, 52, 54 P5 2, 6, 12, 15, 19, 22, 25, 27, 32,36, 41, 44, 47, 50, 52, 54, 57 P6 2, 6, 12, 13, 19, 22, 25, 27, 32,36, 39, 41, 44, 47, 50, 52, 54 P7 1, 4, 9, 19, 22, 26, 29, 30, 32,36, 41, 45, 46, 49, 51, 54, 57 P8 1, 4, 9, 19, 22, 26, 29, 30, 32,36, 41, 45, 46, 47, 50, 54, 57 73 Table 3.3: Two-phase investment plan for alternative solutions of the IEEE 57-bus system Investment Plan PMULocations - Phase1 PMUs Phase1 PMULocations- Phase2 PMUs Phase2 Total U1 3,6,12,15,19,22,25,27,32, 36,39,41,45,47,50,52,55 17 1,4,9,20,24,26,29,30,33, 35,38,43,46,51,54,56,57 17 34 U2 2,6,12,19,22,25,27,32,36, 39,41,45,46,49,51,52,55 17 1,4,9,15,20,23,26,29,30,33, 34,37,43,44,47,50,53,56 18 35 U3 2,6,12,19,22,25,27,32,36, 39,41,45,46,49,50,52,55 17 1,4,9,10,14,20,23,26,29,30, 33,35,43,44,48,54,56,57 18 35 U4 2,6,12,19,22,25,27,29,32, 36,41,45,46,49,51,52,54 17 1,4,9,11,15,20,24,28,31,33, 35,38,39,47,50,55,56,57 18 35 U5 2,6,12,15,19,22,25,27,32, 36,41,44,47,50,52,54,57 17 1,4,9,20,24,28,29,31,33, 35,38,39,43,46,51,55,56 17 34 U6 2,6,12,13,19,22,25,27,32, 36,39,41,44,47,50,52,54 17 1,4,9,20,24,26,29,31,33,35, 38,43,45,46,51,55,56,57 18 35 U7 1,4,9,19,22,26,29,30,32, 36,41,45,46,49,51,54,57 17 3,6,12,15,20,24,28,31,33, 35,38,39,43,47,50,53,56 17 34 U8 1,4,9,19,22,26,29,30,32, 36,41,45,46,47,50,54,57 17 3,6,12,15,20,24,28,31,33, 35,37,38,43,51,53,56 16 33 3.3.4 Case Study II: The Impact of Transmission Switching on Optimal PMU Placement Transmission switching reduces the electricity generation cost by temporarily removing ine cient transmission lines out of service. To demonstrate how transmission switching 74 in uences the optimal placement of PMUs, we use IEEE 57-bus system, PMU placement U8 given in Table 3.3 and developed scenarios given in Table 3.4. In each scenario, we assume that given transmission lines are removed temporarily from the power grid. Table 3.4: Transmission switiching scenarios Scenario Removed Transmission Lines S1 No Line S2 7-8, 10-51, 11-13, 13-14 S3 9-12, 19-20, S4 6-7, 7-29 In scenario S1 when there is no removed transmission line, all buses are N 1 observable. In scenario S2, all buses are still N 1 observable although four transmission lines are temporarily removed. In scenario S3, removing two transmission lines a ects the N 1 observability of the power grid since buses 9, 12, 19 and 20 would no longer be N 1 observable. Hence, a single contingency such as failure of a PMU may further a ect the full observability of the grid. Removing lines 6-7 and 7-29 in scenario S4 results in not observability of bus 7. It can be concluded that it is essential to consider transmission switching in the optimal placement of PMUs since it changes the topology of the power grid. 3.4 Two-phase optimal PMU placement model considering transmission switch- ing The two-phase optimal placement of PMUs considering transmission switching is mod- eled using integer linear programming. The optimal placement of PMUs to achieve N 1 observability occurs in two phases. Hence, the objective function is to minimize the present value of the total investment costs. To calculate the present value, we use the concept of constant dollar analysis [24] in which price increases due to the in ation are ignored. Hence, the in ation-free interest rate, I, is used to compare the investment costs of the two phases. 75 It is assumed that the second phase of PMU placement occurs K years after the initial phase. Therefore, the objective function is as follows. Z = minx;y X j2 Cjxj + 1(1 +I)K X j2 KCjyj ! (3.5) Where xj equals 1 if a PMU is installed in bus j in phase 1. Similarly, yj equals 1 if a PMU is installed in bus j in phase 2. Moreover, it is likely that PMU prices decrease annually due to the free market competitions and advancements in manufacturing of PMUs. We use parameter , the annual percentage of the change in PMU prices, in the objective function to consider these likely price decreases. Notice that the objective function favors installing PMUs in the second phase. The goal of the rst phase is to achieve the full observability of the power grid. Therefore, the constraints given in equations (3.6)-(3.7) should be satis ed. i = X j2 Hi;jxj 8i2 (3.6) i 1 8i2 (3.7) However, notice that equations (3.6)-(3.7) ensure the full observability of the power grid only when all transmission lines of the power system are in service. As it is mentioned, transmission switching changes the topology of the network. Therefore, it is necessary to make sure that allocated PMUs make the power grid observable for all potential topologies of the power grid stem from transmission switching. Hence, we modify equations (3.6)-(3.7) and obtain equations (3.8)-(3.9) as follows. 76 i = X j2 H i;jxj 8i2 ;8 2 (3.8) i 1 8i2 ;8 2 (3.9) Where i represents the observability number of bus i after the rst phase of PMU place- ments under topology 2 , and H i;j is the connectivity parameter of buses i and j under topology 2 . Furthermore, N 1 observability should be accomplished in the second phase consid- ering di erent topologies stem from the transmissions switching. Therefore, we consider the constraints given in equations (3.10)-(3.11). i = i + X j2 H i;jyj 8i2 ;8 2 (3.10) i 2 8i2 ;8 2 (3.11) Where i represents the observability number of bus i after the second phase of PMU place- ments under topology 2 . Next, it should be ensured that at most one PMU is installed in each bus. Therefore, the constraint given in equation (3.12) should be satis ed. xi +yi 1 8i2 (3.12) Finally, equation (3.13) determines that the decision variables are binary. xi;yi2f0;1g 8i2 (3.13) 77 3.5 Experimental Results To test the performance of our algorithm, we use data from the IEEE Reliability Test Systems [22]. We run our experiments on the IEEE 14-bus, IEEE 24-bus, IEEE 30-bus, IEEE 39-bus, IEEE 57-bus, IEEE 118-bus, IEEE 300-bus, IEEE 2383-bus and IEEE 3120-bus test systems. In our experiments, we set the cost of installing a PMU at bus j to 1, Cj = 1. Also, we assume that the price of manufacturing a PMU is the same in both phases of investment, = 1. We arbitrarily chose the value of 0:5% for the in ation-free interest rate, I = 0:5%. Finally, we assume that the two phases of investment occur in two consecutive years, K = 1. 3.5.1 Optimal PMU Placement without Transmission Switching We use our two-phase ILP model to nd the optimal locations of the PMUs such that full observability of the power grid is achieved in the rst phase, and N 1 observability is accomplished in the second phase. The optimal placements are provided in Table 3.5. Table 3.5: Two-phase investment plan for the IEEE test systems TestSystem PMULocations- Phase1 PMUsPhase1 PMULocations- Phase2 PMUsPhase2 Total IEEE14-busSystem2,6,7,9 4 1,3,8,11,13 5 9 IEEE24-busSystem3,4,8,10,16,21,23 7 1,2,7,11,15,17,20 7 14 IEEE30-busSystem2,3,6,10,11,12,15,18,25,27 10 1,7,8,9,13,16,19,21,23,26,29 11 21 IEEE39-busSystem2,6,9,10,13,14,17,19,22,23,29,34,37 13 1,3,8,11,16,20,25,26,30,31,32,33,35,36,38 15 28 IEEE57-busSystem1,4,9,19,22,26,29,30,32,36,41,45,46,47,50,54,57 17 3,6,12,15,20,24,28,31,33,35,37,38,43,51,53,56 16 33 IEEE118-busSystem 2,5,9,12,15,17,21,24,26,28,34,37,41,45,49,53,56,62,63,68,71, 75,77,80,85,87,90,94,102,105,110,114 32 3,6,10,11,19,22,27,29,30,32,35,40,44,46,51,54,57,59,64, 66,70,73,79,84,86,89,92,96,100,107,108,111,112,116,117,118 36 68 The results of our proposed method are compared with available PMU placement meth- ods given in Table 3.1. Notice that the total number of PMUs required for N 1 observability in a two-phase plan is the same as that of a single-phase plan for the given test systems. However, there is no guarantee that the two-phase plan always returns the same number of 78 PMUs as the single-phase plan. Fortunately, our method is exible to obtain the optimal single-phase placement of PMUs as well by setting the value of the parameter k to zero. It is one of the advantages of our two-phase optimal PMU placement model that it provides the utility companies with more exibility. They obtain the optimal PMU placements for single-phase and two-phase plans and decide whether to install all PMUs in one or in two phases while avoiding any potential unnecessary investments. Next, we use our model to obtain the number of PMUs required for N 1 observability in larger systems. The results are given in Table 3.6 for IEEE 300-bus, IEEE 2383-bus and IEEE 3120-bus test systems. The results show that there are PMU placements for these large test systems, which are optimal solutions for both single-phase and two-phase PMU placement methods. Table 3.6: Optimal number of PMUs in large power systems Test System PMUs Phase 1 PMUs Phase 2 Total Single Phase IEEE 300-bus 87 115 202 202 IEEE 2383-bus 746 935 1681 1681 IEEE 3120-bus 994 1212 2206 2206 3.5.2 Optimal PMU Placement with Transmission Switching In this section, we use IEEE 118-bus system to test the performance of our ILP model for two-phase optimal placement of PMUs considering the transmission switching. This test system consists of 118 buses, 186 transmission lines and 54 generators. We refer the readers to [22] and [23] for additional system data. 79 We create the scenarios given in Table 3.7 randomly in which mentioned transmission lines are temporarily removed. In real world, the transmission lines which might be taken o the grid for transmission switching are known in advance. Thus, di erent topologies of the power grid stem from transmission switching can be obtained by enumeration. Table 3.7: Transmission switiching scenarios for IEEE 118-bus system Scenario Removed Transmission Lines W1 No Line (No Transmission Switching) W2 17-18, 23-25, 37-39, 78-89, 80-81 W3 24-70, 48-49, 103-110 W4 32-113, 65-66, 83-85, 103-105 W5 23-25, 41-42, 49-66, 71-72, 80-96 Next, we obtain the optimal placements for all given scenarios and provide the results in Table 3.8. At the bottom of Table 3.8, we consider all scenarios together to gure out the PMU placement which makes all resulting topologies N 1 observable if our proposed ILP model considering transmission switching is not employed. Therefore, it requires to place 94 PMUs in service, 61 PMUs in phase-1 and 33 PMUs in phase-2. Notice that the investment cost of the rst phase is increased considerably. Afterwards, we use our proposed ILP model to nd the optimal solution for the men- tioned transmission switching scenarios. The optimal placement is provided in Table 3.9. It shows that the optimal placement for N 1 observability of the IEEE 118-bus system with mentioned transmission switching scenarios requires 76 PMUs, 35 PMUs in phase-1 and 41 PMUs in phase-2. Therefore, our proposed model could decrease a great deal of investment costs comparing to the given placement W1 W5 provided in Table 3.8. Moreover, notice that the optimal solution needs much less number of PMUs in the rst phase which is critical 80 Table 3.8: Optimal placements for the transmission switching scenarios of IEEE 118-bus system ScenarioPMULocations- Phase1 PMUsPhase1 PMULocations - Phase2 PMUsPhase2 Total W1 2,5,9,12,15,17,21,24,26,28,34,37,41,45,49,53,56,62,63,68,71, 75,77,80,85,87,90,94,102,105,110,114 32 3,6,10,11,19,22,27,29,30,32,35,40,44,46,51,54,57,59,64, 66,70,73,79,84,86,89,92,96,100,107,108,111,112,116,117,118 36 68 W2 2,5,9,12,15,19,22,27,30,31,32,35,40,43,45,49,52,56,62,64,68, 70,71,75,77,80,85,87,90,94,101,105,110 33 3,6,10,11,17,18,21,24,25,28,34,37,39,41,46,51,53,57,59,61,67, 73,78,79,81,83,86,89,92,96,100,106,108,111,112,114,116,117,118 39 72 W3 2,5,9,12,15,17,21,23,28,30,36,40,44,46,50,52,56,62,63,68,71, 75,77,80,85,87,90,94,102,105,110,115 32 3,6,10,11,19,20,25,29,35,37,42,43,48,49,51,54,61,64,67,70,72, 73,79,84,86,89,92,96,100,107,108,111,112,113,114,116,117,118 38 70 W4 2,5,9,12,15,17,21,25,29,34,37,40,45,49,53,56,62,64,68,70,71, 76,79,84,87,89,92,96,100,105,110,114 32 3,6,10,11,19,22,26,28,31,35,41,44,46,51,54,57,59,65,67,72,73, 75,77,80,82,85,86,90,94,102,107,108,111,112,113,115,116,117 38 70 W5 2,5,9,12,15,17,21,24,25,29,34,37,40,45,49,53,56,62,63,68,73, 75,77,80,85,87,90,94,101,105,110,114 32 3,6,10,11,19,22,27,28,30,32,35,41,44,46,51,54,57,61,64,66,71, 72,74,78,84,86,89,92,96,100,107,108,111,112,116,117,118 37 69 W1 W5 2,5,9,12,15,17,19,21-32,34-37,40,41,43-46,49,50,52,53,56, 62-64,68,70,71,73,75-77,80,84,85,87,89,90,92,94,96,100-102, 105,110,114,115 61 3,6,10,11,18,20,39,42,48,51,54,57,59,61,65,66,67,72,74,78,81, 82,83,86,106,107,108,111,112,113,116,117,118 33 94 from investment point of view. It is likely that PMU prices decrease later which bene ts a utility by decreasing the investment costs even more. 81 Table 3.9: Optimal placement for IEEE 118-bus system considering transmission switching scenarios Phase PMU Locations Number of PMUs 1 2, 5, 9, 12, 15, 17, 19, 21, 24, 25, 28, 35, 40, 43, 46, 49, 52, 56, 59, 62, 65, 68, 73, 75, 77, 80, 84, 87, 89, 92, 94, 102, 105, 110, 114 35 2 3, 6, 10, 11, 18, 22, 26, 27, 29, 36, 37, 39, 41, 44, 48, 53, 57, 58, 63, 66, 70, 71, 72, 76, 78, 79, 81, 82, 85, 86, 90, 96, 100, 107, 108, 111, 112, 113, 115, 116, 117 41 Furthermore, it is likely that a utility would like to know if they could decrease the PMU investments by not using one or more of the transmission switching scenarios. Table 3.10 provides the optimal number of PMUs if a certain scenario would not be used. As a case in point, it can be concluded that including scenario W2 may be a strategic decision since including W2 itself increases the number of PMUs as many as including all three scenarios W3, W4 and W5. Moreover, notice that the utility needs to install eight more PMUs due to the transmission switching scenarios. Otherwise, the IEEE 118-bus system requires only 68 PMUs for N 1 observability. 82 Table 3.10: Optimal number of PMUs for removed transmission switching scenarios Removed Scenario PMUs Phase 1 PMUs Phase 2 Total None 35 41 72 W2 34 38 72 W3 34 41 75 W4 33 41 74 W5 35 40 75 W2;W3 34 37 71 W2;W4 33 37 70 W2;W5 33 38 71 W3;W4 33 40 73 W3;W5 34 40 74 W4;W5 33 41 74 W2;W3;W4 32 37 69 W2;W3;W5 32 38 70 W2;W4;W5 32 38 70 W3;W4;W5 33 39 72 W2;W3;W4;W5 32 36 68 83 3.5.3 Computational Time Analysis All experiments were performed on a 64-bit laptop with an Intel Core i5 2.4 GHz pro- cessor and 4GB RAM. We used IBM ILOG CPLEX Optimization Studio 12.6 as the opti- mization solver in our experiments. Table 3.11 summarizes the computational time for obtaining the optimal PMU place- ments for di erent cases discussed in Section 3.5.1. Computational time analysis shows that our two-phase optimal PMU placement model is solved very quickly even for large power systems. Table 3.11: Computational Times (s) Test System Single-phase Two-phase IEEE 14-bus 0.81 0.81 IEEE 24-bus 0.82 0.82 IEEE 30-bus 0.82 0.83 IEEE 39-bus 0.82 0.84 IEEE 57-bus 0.87 0.89 IEEE 118-bus 0.97 1.00 IEEE 300-bus 1.36 1.39 IEEE 2383-bus 33.15 33.71 IEEE 3120-bus 57.79 58.88 3.6 Conclusions Observability of the power system is important for grid operation and control. Complex networks of PMUs, as the state-of-the-art measurement devices, are used to collect real time data to improve the observability of the power systems. However, due to budget restrictions and considerable cost of installing PMUs, it is not possible to equip all buses with PMUs. Also, it has been discussed that the reliability requirement of N 1 observability and including transmission switching in the optimal placement of PMUs increase the investment costs signi cantly. Therefore, the network observability rules and PMUs? characteristics should be used to the full extent in order to minimize the PMU investment costs. In this 84 paper, we developed an integer linear programming model that minimizes the investment costs of PMUs and determines the optimal locations of PMUs in two investment phases. In the rst phase, PMUs are installed to achieve full observability of the power grid whereas additional PMUs will be installed in the second phase to guarantee the N 1 observability of the power grid. The proposed model is able to provide utilities with single-phase and two-phase optimal placements, which gives investors more exibility on whether to install all PMUs in one or in two phases while avoiding any potential unnecessary costs. Furthermore, it has been shown that it is critical to consider transmission switching in the optimal PMU placement problem. Transmission switching changes the topology of the power grid and may cause not observability of some buses. In our ILP model, we integrated the transmission switching concept into the optimal PMU placement problem such that the obtained optimal placement is N 1 observable for all topologies of the power grid stem from transmission switching. The performance of the developed ILP model is tested on several IEEE test systems. Experimental results show that our developed model is a promising enhancement to the optimal PMU placement problem and also ensures the observability of the power systems when transmission switching is utilized. Computational times have also been reported to show that our model can be solved very quickly by optimization solvers even for large power systems. 85 References [1] B. Singh, N. Sharma, A. Tiwari, K. Verma, and S. Singh, \Applications of phasor measurement units (PMUs) in electric power system networks incorporated with facts controllers," International Journal of Engineering, Science and Technology, vol. 3, pp. 64{82, 2011. [2] K. P. Briman, L. Ganesh, and R. van Renesse, \Running smart grid control software on cloud computing architectures," in Workshop Computational Needs for the Next Generation Electric Grid. Cornell University, 2011, pp. 1{33. [3] A. Monticelli, State Estimation in Electric Power Systems: A Generalized Approach. Norwell, MA: Kluwer Academic Publishers, 1999. [4] A. B. Antonio, J. R. A. Torreao, and M. B. D. C. Filho, \Meter placement for power system state estimation using simulated annealing," in Proceedings of the IEEE Power Tech Conference, 2001. [5] B. Milosevic and M. Begovic, \Nondominated sorting genetic algorithm for optimal phasor measurement placement," IEEE Transactions on Power Systems, vol. 18, pp. 69{75, 2003. [6] R. Sodhi, S. C. Srivastava, and S. N. Singh, \Optimal PMU placement method for com- plete topological and numerical observability of power system," Electric Power System Research, vol. 80, pp. 1154{1159, 2010. [7] S. Chakrabarti and E. Kyriakides, \Optimal placement of phasor measurement units for power system observability," IEEE Transactions on Power Systems, vol. 23, pp. 1433{1440, 2008. 86 [8] B. K. S. Roy, A. K. Sinha, and A. K. Pradhan, \An optimal PMU placement technique for power system observability," Electrical Power and Energy Systems, vol. 42, pp. 71{ 77, 2012. [9] T. L. Baldwin, L. Mili, J. M. B. Boisen, and R. Adapa, \Power system observability with minimal phasor measurement placement," IEEE Transactions on Power Systems, vol. 8, pp. 707{715, 1993. [10] R. F. Nuqui and A. G. Phadke, \Phasor measurement unit placement techniques for complete and incomplete observability," IEEE Transactions on Power Delivery, vol. 20, pp. 2381{2388, 2005. [11] J. Peng, Y. Sun, and H. F. Wang, \Optimal PMU placement for full network observ- ability using tabu search algorithm," Electrical Power and Energy Systems, vol. 28, pp. 223{231, 2006. [12] A. Ahmadi, Y. Alinejad-Beromi, and M. Moradi, \Optimal PMU placement for power system observability using binary particle swarm optimization and considering measure- ment redundancy," Expert Systems with Applications, vol. 38, pp. 7263{7269, 2011. [13] B. Xu and A. Abur, \Optimal placement of phasor measurement units for state estima- tion," PSERC, Tech. Rep., 2005. [14] S. Azizi, A. S. Dobakhshari, S. A. N. Sarmadi, and A. M. Ranjbar, \Optimal PMU placement by an equivalent linear formulation for exhaustive search," IEEE Transac- tions on Smart Grid, vol. 3, pp. 174{182, 2012. [15] S. Azizi, B. G. Gharehpetian, G. Hug-Glanzmann, and A. Dobakhshari, \Optimal in- tegration of phasor measurement units in power systems considering conventional mea- surements," IEEE Transactions on Smart Grid, 2012. 87 [16] C. Liu, J. Wang, and J. Ostrowski, \Static switching security in multi-period transmis- sion switching," IEEE Transactions on Power Systems, vol. 27, no. 4, pp. 1850{1858, November 2012. [17] G. Granelli, M. Montagna, F. Zanellini, P. Bresesti, R. Vailati, and M. Innorta, \Op- timal network reco guration for congestion management by deterministic and genetic algorithms," Electric Power System Research, vol. 76, p. 549556, April 2006. [18] R. P. O?Neill, R. Baldick, U. Helman, M. H. Rothkopf, and W. Stewart, \Dispatchable transmission in RTO markets," IEEE Transactions on Power Systems, vol. 20, no. 1, p. 171179, February 2005. [19] E. B. Fisher, R. P. O?Neill, and M. C. Ferris, \Optimal transmission switching," IEEE Transactions on Power Systems, vol. 23, no. 3, pp. 1355{1364, August 2008. [20] K. W. Hedman, R. P. O?Neill, E. B. Fisher, and S. S. Oren, \Optimal transmission switching - sensitivity analysis and extensions," IEEE Transactions on Power Systems, vol. 23, no. 3, p. 14691479, August 2008. [21] D. Dua, S. Dambhare, R. K. Gajbhiye, and S. A. Soman, \Optimal zero injection considerations in PMU placement: An ILP approach," IEEE Transactions on Power Delivery, vol. 23, pp. 1812{1820, 2008. [22] Subcommittee, \IEEE reliability test system," IEEE Transactions on Power Apparatus and Systems, vol. 98, pp. 2047{2054, 1979. [23] R. D. Zimmerman, C. E. Murillo-Sanchez, and R. J. Thomas, \MATPOWER: Steady- state operations, planning, and analysis tools for power systems research and education," IEEE Transactions on Power Systems, vol. 26, pp. 12{19, 2011. [24] C. S. Park, Fundamentals of Engineering Economics, 3rd ed. Prentice Hall, 2013. 88