Reliability Analysis of Electric Power Systems
Using an Objectoriented Hybrid Modeling Approach
Markus Schläpfer, Tom Kessler, Wolfgang Kröger\endnames\namesSwiss Federal Institute of Technology\endnames\namesZurich, Switzerland\endnames\namesschlaepfer@mavt.ethz.ch\endnames
The ongoing evolution of the electric power systems brings about the need to cope with increasingly complex interactions of technical components and relevant actors. In order to integrate a more comprehensive spectrum of different aspects into a probabilistic reliability assessment and to include timedependent effects, this paper proposes an objectoriented hybrid approach combining agentbased modeling techniques with classical methods such as Monte Carlo simulation. Objects represent both technical components such as generators and transmission lines and nontechnical components such as grid operators. The approach allows the calculation of conventional reliability indices and the estimation of blackout frequencies. Furthermore, the influence of the time needed to remove line overloads on the overall system reliability can be assessed. The applicability of the approach is demonstrated by performing simulations on the IEEE Reliability Test System 1996 and on a model of the Swiss highvoltage grid.
Reliability analysis, Monte Carlo simulation, blackout frequency distribution, operator response time
1 Introduction
\PARstartThe ongoing evolution of the electric power systems due to market liberalization and the integration of distributed generation is leading to increasingly complex and hardtopredict interactions of technical components, relevant actors and the operating environments. Furthermore, recent largearea blackouts in North America and Europe demonstrated the potential consequences of inadequate operator response times to contingencies (e.g. [1]). In recent years several advanced methods have been developed to assess the reliability of electric power systems in general and to model and analyze cascading blackouts (e.g. [2, 3]). However, these approaches do not explicitly simulate the evolution of the events in time and represent the operator intervention to contingencies by using highly simplified models not taking into account the time needed for the corrective action. While Anghel et al. [4] introduce a timedependent probabilistic approach incorporating a model for the utility response to line overloads, the influence of the response time on the occurrence of cascading line outages remains neglected.
The contribution of this paper is to present a basic modeling framework which allows the explicit integration of highly nonlinear, timedependent effects and nontechnical factors into a probabilistic reliability assessment. Therefore, we apply an objectoriented hybrid approach combining agentbased modeling techniques [5] with classical methods such as Monte Carlo simulation [6]. Objects represent both technical components such as generators and transmission lines and nontechnical components such as grid operators. They interact with each other directly (e.g. via the generator dispatch) or via the physical power flows on the network. By means of longterm simulations the statistical data is gathered for the calculation of system reliability indices and for the estimation of blackout frequencies.
The paper is organized as follows. Section 2 introduces the conceptual basics of the modeling framework and the derivation of the different component models. In section 3 we present the results of applying the model to the IEEE Reliability Test System 1996 and to the Swiss highvoltage system. Section 4 concludes.
2 Modeling Framework
2.1 Conceptual Basics
The conceptual modeling framework consists in the abstraction of the relevant technical and nontechnical components of the electric power system as individual interacting objects. Each object is modeled by attributes and rules of behavior. An example for an attribute is a technical component constraint such as the rating of a transmission line. The rules of behavior are represented by using finite state machines (FSM) and include both deterministic and stochastic timedependent, discrete events. A deterministic event is, for instance, the outage of a component when reaching a failure threshold, while stochastic processes are probabilistic component failure models using Monte Carlo techniques. The integration of nontechnical components is demonstrated by modeling the behavior of the grid operators in case of line overloads. For the corresponding interactions between the operators and the technical components we make use of agentbased modeling techniques. Furthermore, we account for the possible division of the power system into several control areas. To each control area a distinct grid operator and a distinct control object are assigned. The control object is not an abstraction of a technical component as such but rather represents an implementation construct which controls the balance between generation and load within the corresponding control area. The model captures the system behavior over an operational period of one year.
2.2 Component Models
The components of the power system as being modeled as objects are loads, generators, transmission lines, busbars and grid operators.
Loads
The power demand trajectory of load is described by:
(1) 
The demand factor is the actual timedependent percentage of the peak demand and follows a chronologically changing load profile over the predefined time period of one year. The percent deviation represents stochastic demand fluctuations and is sampled hourly assuming a normal distribution with and standard deviation according to [7]. The value of is assumed to be equal for all loads within the same control area. The term represents the actual amount of partially shed load.
Figure 2.2.1 shows the FSM as implemented for the load objects. With the exception of the restoration process all the transitions of the fourstate model are externally governed by the control object of the corresponding control area.
Partial load shedding occurs only when the control object sheds load due to an operator action for removing a line overload. As soon as the transmission system can be operated within its security margins again, the load object receives the signal to cancel the partial load shedding. The load gets totally disconnected if there is not enough generation capacity available within the entire system to cover its demand or as a consequence of system splitting (see section 2.3). If several loads have to be disconnected within one control area all loads are given the same priority to be shed and are therefore selected randomly.
The restoration process is started once enough generation capacity is available again to cover the disconnected demand, and is modeled by a queue technique. The load which has been disconnected first is also restored first, the subsequent one waits until the previous is reconnected. Based on [8] we assume an incremental overall restoration rate for four different restoration stages according to table 2.2.1, where is the elapsed time measured from the start of the overall restoration process . Hence, the time needed to reconnect a specific load, , is dependent on the actual overall restoration stage.
[min]  [MW/min] 

030  10.0 
3060  33.3 
6090  66.6 
90  83.3 
Generators
The commitment of the generating units is continuously governed by the control object in order to cover the actual demand within the respective control area . Being constrained by the maximum power outputs , the commitment and economic dispatch follows a heuristic priority list method according to [9] and is implemented in the control object. By using a recursive algorithm and starting with the highest priority, is equally distributed among the units with the same priority. As their maximum capacity is reached the algorithm proceeds to the next lower priority and so forth. In case is larger than the available generation capacity within a control area, the control object analogously commits available generating units from the other control areas of the system. As a simplification, ramp rates and maintenance are not considered at this stage of our work. The FSM for the generator object is made up of a twostate model for repairable forced failures being treated as random events as shown in figure 2.2.2.
The repairable forced failures are modeled by an independent stochastic updownup cycle assuming stationarity and constant failure and repair rates and respectively. Hence, this alternating renewable process is characterized by the cumulative distribution functions of the failurefree times and repair times respectively, and by the probability that the generating unit is in upstate at [10]:
(2) 
(3) 
(4) 
where and are the time spans measured from the moment of entering the upstate and forced down state respectively. All state transitions are reported to the control object of the corresponding control area.
Transmission lines
The time variant line flows are calculated by the DC power flow approximation with , where is the active power flow on line connecting busbar with busbar , having reactance and phase angles and . The approximate solution of the power flow problem does not allow to analyze voltage disturbances. Nevertheless, we assume the DC model to be appropriate for analyzing cascading events due to line overloads and for showing the feasibility of the proposed modeling concept.
A fivestate model for the basic behavior of the transmission line is used considering outages triggered by its protection device and by independent random failures, see figure 2.2.3. Thereby, the protection device is modeled by a separate FSM.
In our model the protection device has two functions. Firstly, it continuously measures the power flow and sends an alarm message (”overload”) to the operator of the control area if becomes equal to or larger than the line rating . Secondly, if reaches it disconnects the line. However, as a consequence of the stochastic timedependent system behavior or the intervention of the operator, may again fall back to less than before reaching and the protection device returns to the idle state. By following the assumptions made by Zima and Andersson [11] the probability for the line outage increases linearly with the power flow, being zero below . Thus, we assume being uniformly distributed in the interval with =1.4. The line is either reconnected if the phase angle difference becomes smaller than the preset value or after a time delay of which models the time until a manual attempt to reclose the breakers. The parameter is used to avoid an immediate recurrence of the overload, potentially resulting in a persistently repeating state change cycle, and is set to =0.9. Analogous to the probabilistic failure model of the generating units (equations (24)), the time to permanent outage and the time to repair follow an exponential distribution with failure rate and respectively.
Busbars
Every busbar object continuously calculates its phase angle relative to its neighboring busbars:
(5) 
where is the net power injection at busbar , to which several loads and generating units might be connected. This distributed approach allows avoiding time consuming matrix calculations in case of network decompositions and restorations due to line outages and reconnections. Potential random outages of busbars are not considered.
Grid operator
The grid operator becomes active in case of transmission line overload contingencies, trying to remove the overload by redispatching the generators or by shedding load if necessary. The basic model for the operator behavior is illustrated for the overload of a tieline between two control areas, see figure 2.2.5. If a tieline becomes overloaded the protection device sends an alarm message to the two operators of both control areas (compare figure 2.2.3). Having the alarm received the neighboring operators try to contact each other with a time delay . The operator which has been assigned responsible for the line then needs some time to find a solution to the overload problem, which is modeled by a time delay .
The corrective action to remove the overload is subsequently formulated as a conventional optimal power flow (OPF) problem [9] and implemented within the control object by using the linear programming (LP) method minimizing potential load shedding, , and the change in generation, , subject to the transmission line constraints and the power balance:
(6) 
subject to
(7) 
(8) 
(9) 
(10) 
where is the busbar specific distance weighting factor and set to =1 for the two busbars at each end of the overloaded line, =10 for the busbars being one line further away and =100 for all other busbars within the overall system. The weighting factor =10000 lets partial load shedding be more expensive relative to the generator redispatch. The linear line sensitivity factor with respect to busbar is dependent on the network connectivity at the model time and is calculated using the conventional matrix method as described in [9]. Equation (10) holds for all lines within the two neighboring control areas. Similar to the model for the reconnection of a failed line, the parameter is used to delay the potential recurrence of an overload and set to =0.8.
The procedure for line overloads within a single control area is basically the same, but without the interaction of the operators and by restricting the load control variables to the busbars and equation (10) to the transmission lines within the control area. In order to prioritize the generator redispatch within the control area concerned, the distance weighting factor is set to =1 for busbars inside and =100 for busbars outside the control area.
2.3 System Splitting
The splitting of the network due to transmission line outages usually leads to an imbalance between load and generation within the separated subsystems. Further, depending on the total inertia within the separated parts, on the frequency control performance and the protection device behavior of the generators, and on implemented automatic load shedding procedures, this imbalance comes along with a frequency deviation potentially leading to stability problems [12]. The consequences range from small load losses to a total collapse of the subsystem (e.g., [1]). In order to include load outages as a consequence of a network splitting while avoiding a complicated model with a high amount of parameters to be estimated we make use of a highly simplified approach. Thereby, an excess of demand within a separated subsystem leads to the immediate disconnection of a minimum number of randomly selected loads so that the excess is at least reduced to zero. An excess of generation leads to the immediate reduction of the generator outputs in order to reestablish the balance and implies no load outages. This strong simplification might be inadequate to represent the real system behavior and the amount of disconnected load thus has to be viewed as a rather indicative value for the system vulnerability regarding the splitting of the network.
2.4 Blackout Frequency Distributions
By means of longterm simulations (i.e. repeatedly over the operation period of one year) the necessary statistical data is gathered for the calculation of conventional reliability indices such as the Expected Energy Not Supplied (EENS). Moreover, frequency distributions of expected blackouts per year are derived. Therefore, let be a random variable counting the number of blackouts per year greater than a specified size . The size is thereby classified by the unserved energy or the maximum amount of demand not being supplied in the course of an event. The expectation is approximated by generating realizations of and calculating their empirical mean, which represents the observed complementary cumulative frequency of events related to one year, denoted by :
(11) 
where denotes the number of simulated years.
Assuming that follows a Poisson distribution, the confidence interval for can be constructed by using the central Chisquare distribution [13]:
(12) 
where is the confidence level, is the probability of error, and and are the degrees of freedom. The blackout events can further be classified into the three outage causes as implemented in the model:

Generation inadequacy:
Loads are disconnected as not enough generation capacity is available to cover the actual demand within the overall system or within a previously separated subsystem. 
System splitting:
Loads are disconnected as a consequence of the separation of the system. 
Operator intervention:
Load is partially shed in order to remove transmission line overloads.
3 Case Studies
3.1 Application to the IEEE Reliability Test System 1996
System layout and model parameters
The threearea IEEE Reliability Test System 1996 (RTS96) has 73 busbars, 120 transmission lines and 96 generating units [14]. We use the yearlong load data with an hourly resolution provided in [14] for modeling the demand trajectories . The three areas have a base case peak load of 2850 MW each and are treated as three single control areas with three corresponding control objects and grid operators. The priorities given to the different generator types are shown in table 3.1.1.
Unit Type  [MW]  Priority 

Hydro  50  1 
Nuclear  400  2 
Coal/Steam  350  3 
Coal/Steam  155  4 
Coal/Steam  76  5 
Oil/Steam  197  6 
Oil/Steam  100  7 
Oil/Steam  12  8 
Oil/CT  20  9 
The failure and repair rates for the generators and the transmission lines are taken from [14]. The parameter value for the time until the manual attempt to reclose the breaker of a disconnected line is assumed to be =1h. Regarding the operator model is set to 2min.
Computational results
The results of two parameter variation studies are presented and discussed: 1) the sensitivity of the blackout frequency to an increase of the system loading without any operator intervention and 2) the influence of the operator response time on the Expected Energy Not Supplied (EENS). Concerning the first experiment we increment the system loading level . The maximum generator outputs, , are augmented by the same factor. In order to gain statistically significant results (i.e. 1000) about 50 hours of simulation time are needed on a single desktop computer (Dell Optiplex GX260 with a Pentium 4 CPU of 2.66GHz and 512MB of RAM). This time was considerably reduced by running several simulations in parallel. Figure 3.1.2 shows the resulting complementary cumulative blackout frequencies with respect to the unserved energy per event, , for four different values of . Regarding the two lower system loading levels (=1.0 and =1.1) the observed complementary cumulative frequencies follow approximately an exponential curve. However, increasing to 1.2 already leads to a remarkable increase of large events, while the shape of the curve in the range of the smaller events (up to about MWh) stays qualitatively the same. The value of =1.37 represents the maximum system loading level where no line overloads would occur without any stochastic component outages. This loading level can be characterized by a high frequency of large blackouts predominantly in the range between MWh and MWh.
In order to further analyze the differences between the overall frequency curves the distributions of the underlying power outage causes have to be considered. The logarithmic histograms of figure 3.1.2 report the impact of increasing the system loading from =1.0 to =1.37 on the absolute frequency of blackouts caused by generation inadequacy (left hand side) and system splitting (right hand side). System splitting is the predominant cause of the observed blackouts for both loading levels. In comparison to generation inadequacy the absolute frequencies for this outage mode show a stronger increase and a stronger shift towards larger events when it comes to an increase of the system loading. Hence, the substantial increase of large blackouts as shown in figure 3.1.2 is mainly the result of an increased frequency of line overloads and subsequent system splitting.
The results of our second parameter variation study are presented in figure 3.1.2, showing the influence of the operator response time on the EENS broken down into the different outage causes for the system loading level =1.37. For the interpretation of the results it should be reminded that thermal aspects of the line overloads are not taken into consideration. Under our model assumptions an operator intervention with a delay of 5 hours after the occurence of the overload still reduces the EENS due to system splitting by about 30%. On the other hand, an increase of the response time from 15min to 30min leads to a significant increase of the EENS due to system splitting of about 26%. The EENS due to generation inadequacy is increasing with the response time as the system is more often separated which, in turn, reduces the redundancy of the generators within the splitted subsystems. The values for the EENS due to the operator intervention are negligible.
3.2 Application to the Swiss HighVoltage Grid
System layout and model parameters
The Swiss electric power system consists of a single control area with an annual energy consumption of GWh and a peak load of about 10 GW. The energy production and installed capacity total to GWh and 12 GW respectively, consisting of 42.2% nuclear, 52.4% hydro and 5.4% conventional thermal generation [15]. The number of components as used in our model for the 380/220 kV transmission grid are shown in table 3.2.1.
99  34  229  161  1 

Based on a particular system snapshot taken on a January morning the fluctuating power injections at the different nodes are derived by using publicly available statistical data [15]. For each hydro power generator a different production capacity is assigned for the winter halfyear and the summer halfyear respectively. The failure and repair rates for all hydro generators are equally set to and , and for all nuclear units to and . Regarding the transmission lines the failure model parameters are chosen as and . As the phase shifting transformers have considerable influence on the power flows corrective injections were calculated for the nodes adjacent to a phase shifting transformer. The energy exchange with the neighboring countries is modeled by independent positive or negative power injections at the surrounding boundary nodes. The parameter values for the time until the manual attempt to reclose the breaker of a disconnected line and for the operator response time are assumed to be =1h and =15min, respectively.
Computational results
It should be noted that the intention of the analysis was primarily to investigate the applicability of the proposed modeling method to a real system. The computational results thus make no claim to quantify the reliability of the Swiss highvoltage grid in absolute terms.
The estimated blackout frequencies and the histogram of the different outage causes both with respect to the unserved energy per event are depicted in figure 3.2.2. The model potentially overestimates the duration of the events and thus the unserved energy as switching operations on lower voltage levels for the reconnection of deenergized loads and the possibility to import extra power from neighboring countries to overcome generation shortages are not taken into consideration.
The complementary cumulative blackout frequency follows an exponential curve. Generation inadequacy is the dominant factor regarding the larger events while load shedding for line overload relief becomes important in the range of the smaller events. The influence of load disconnections due to system splitting is significant but the frequency of this outage cause never exceeds the frequency of load disconnections due to generation inadequacy or load shedding due to the operator action. Hence, under our model assumptions, it can be concluded that the system reliability is somewhat more sensitive to generation outages than to transmission line failures.
The benefit of the operator response to line overloads is shown in figure 3.2.2 where the frequencies of the events with and without operator action are compared. The event size is thereby measured by the maximum unserved demand.
The impact of the operator intervention becomes significant in the range of the larger events where a high fraction of blackouts with a size greater than 200 MW is prevented. A large number of disconnected loads due to system splitting and thus a high value for the unserved demand generally needs a high number of subsequently disconnected lines due to overload. Such a sequence of events potentially gives the operator a higher chance to intervene in comparison to a disconnection of a single load due to the outage of a few lines without further cascading failures.
The relative overload frequencies for each transmission line, , are reported in figure 3.2.2. About 15% of all overload contingencies are occurring on only two lines. Furthermore, several groups of adjacent lines can be identified as being prone to overloads, helping to highlight the most critical system regions.
4 Conclusions
We presented an objectoriented hybrid modeling framework for a comprehensive reliability analysis of electric power systems. The main advantages are the explicit integration of highly nonlinear, timedependent effects and the possibility to include nontechnical factors. The chosen level of modeling detail allows analyzing a multitude of different (timedependent) reliability aspects such as the identification of weak points and the assessment of system upgrades. Although several model refinements need to be further developed, the results of the case studies performed on the IEEE RTS96 and on a model of the Swiss highvoltage grid confirm the applicability of the approach with respect to midperiod power system planning purposes. Optimizing the technical implementation of the models [16] together with the evolution of both hard and software will fasten up the simulation speed. Gaining experience in applying the proposed approach will give insight in the parameters to be used, thus lessen the problem of the high number of parameters to be estimated.
Acknowledgements
The authors would like to thank ”swisselectric research” for cofinancing the work, Swissgrid AG for the fruitful collaboration and for providing the operational data of the Swiss electric power system, and Walter Sattinger (Swissgrid AG) for his helpful feedback on the manuscript.
References
 U.S.Canada Power System Outage Task Force, ”Final Report on the August 14, 2003 Blackout in the United States and Canada: Causes and Recommendations”, 2004.
 Chen, J., Thorp, J. S. and Dobson, I., ”Cascading Dynamics and Mitigation Assessment in Power System Disturbances via a Hidden Failure Model”, Electrical Power and Energy Systems, Vol. 27, No. 4, pp. 318326, 2005.
 Rios, M. A., Kirschen, D. S., Jayaweera, D., Nedic, D. P. and Allan, R. N., ”Value of Security: Modeling TimeDependent Phenomena and Weather Conditions”, IEEE Trans. Power Syst., Vol. 17, No. 3, pp. 543548, 2002.
 Anghel, M., Werley, K. A. and Motter, A. E., ”Stochastic Model for Power Grid Dynamics”, 40th Hawaii International Conference on System Sciences, 2007.
 D’Inverno, M., Luck, M., ”Understanding Agent Systems”, Springer, Berlin, 2004.
 Billinton, R., Li, W., ”Reliability Assessment of Electric Power Systems Using Monte Carlo Methods”, Plenum Press, New York, 1994.
 Dai, Y., McCalley, J. D., AbiSamra, N., Vittal, V., ”Annual Risk Assessment for Overload Security”, IEEE Trans. Power Syst., Vol. 16, No. 4, pp. 616623, 2001.
 Kirschen, D. S., Bell, K. R. W., Nedic, D. P., Jayaweera, D. and Allan, R. N., ”Computing the value of security”, IEE Proc.Gener. Transm. Distrib., Vol. 150, No. 6, pp. 673678, 2003.
 Wood, A. J. and Wollenberg, B. F., ”Power Generation, Operation, and Control”, ed., John Wiley & Sons, New York, 1996.
 Birolini, A., ”Reliability Engineering, Theory and Practice”, ed., Springer, Berlin, 2007.
 Zima, M. and Andersson, G., ”On security criteria in power systems operation”, Proc. IEEE Power Engineering Society General Meeting, San Francisco, June 2005.
 Anderson, P. M. and Mirheydar, M., ”An adaptive method for setting underfrequency load shedding relays”, IEEE Trans. Power Syst., Vol. 7, No. 2, pp. 647655, 1992.
 Cowan, G., ”Statistical Data Analysis”, Clarendon Press, Oxford, 1998.
 IEEE RTS Task Force of APM Subcommittee, ”The IEEE Reliability Test System  1996”, IEEE Trans. Power Syst., Vol. 14, No. 3, pp. 10101020, 1999.
 Swiss Federal Office of Energy, ”Swiss Electricity Statistics 2006” (german), Bern 2007.
 Schläpfer, M. and Mancarella, P., ”Probabilistic modeling and simulation of transmission line temperatures under fluctuating power flows”, IEEE Trans. Power Deliver., Vol. 26, No. 4, pp. 22352243, 2011.