Efficient Creation of Datasets for Data-Driven Power System Applications

1 Ef ficient Creation of Datasets for Data-Dri v en Po wer System Applicati ons Andreas V enzke, Daniel K. Molzahn, and Spyros Chatziv a s ileiadis Abstract —Advances in data-driven methods have sparked renewed interest fo r applications in power systems. Creating datasets f or successful application of these methods has prov en to be very challenging, especially wh en considering power system security . This paper p roposes a computationally efficient method to create datasets of secure and insecure operating points. W e propose an infeasibility certificate based on separ ating hyper - planes that can a-priori characterize lar ge parts of the input space as in secure, thus significantly reducing both computation time and p roblem size. Our method can handle an order of magnitude more control variables and creates balanced datasets of secure and insecure operating p oints, which is essential for data-driven appli cations. W h ile we focus on N- 1 security and uncertainty , our method can extend to dynamic security . For PGLib-OPF networks up to 500 buses and up to 125 control variables, we demo nstrate drastic reductions in unclassifi ed input space volumes and computation time, crea te balanced d atasets, and e valuate an illu strativ e data-driv en application. Index T erms —Conv ex relax ation, data-driv en, machine learn- ing, optimal power flow , power system operation. I . I N T R O D U C T I O N Recent advances in data-d riv en metho ds have sh own sub- stantial po tential for power system app lications includ ing secu- rity ass essment under u ncertainty [ 1]–[6], e. g., b y rapidly esti- mating line flo ws [2], training accurate security classifiers [3], and applying these classifiers in the context o f d ata-driven security-co n strained optimal p ower flow [4 ] and deep learning toolboxes [5] . T h e performance of these methods, h owe ver, relies on the quality of the underlying dataset. As h istor ical data is often limited and do es not co ntain many abn ormal situations, the datasets h av e to be enrich e d throu gh simulation. This, however , is a highly co mputation a lly dem anding task. The resulting datasets shou ld be ba la n ced b e tween secur e an d insecure samples to improve classifier per formanc e, take into consideratio n all degrees o f fr eedom of the system, and be able to accurately represent the security bou ndary . In this work , we prop ose an e fficient method to create d atasets with the se proper ties for data- driven ap plications in power systems. The steady-state op erational constraints are de scribed by the AC optimal power flow (A C-OPF) pro blem. The degrees of f reedom of the system, i.e., the inpu ts characterizing each operating poin t, are defined b y the con trol variables, wh ich in A. V enzke and S. Chatzi v asilei adis are with the Department of Electrica l Engineeri ng, T echnica l Univ ersity of Denmark, 2800 Kgs. L yngby , Denmark e-mail: { andven, spchatz } @el ektro.dtu.dk. D. K. Molzahn is with the School of Electrica l and Computer Engineering, Georgi a Institute of T echnology , Atlan ta, GA 30313, USA as well as the Energy Systems Division, Ar gonne Nat ional Labora tory , Lemont, IL 60439 USA, e-mail: molzahn@ga tech.edu. This work is supported by the m ultiDC project , funded by Innov ation Fund Denmark, Grant Agreement No. 6154-00020B. the AC-OPF problem are generator active power and voltage set-points. By definin g these, th e r e maining state variables are determined by solving the AC power flow equ ations [7 ]. Even for medium- sized systems, the resu lting n umber of control variables render s the task of creating datasets covering a wide range of operating points very c o mputation ally challen g ing. T o address this c hallenge, we can directly classify operating points th at are infeasible w ith respect to th e AC-OPF p roblem as insecure a n d av oid any furth e r stability or static secur ity assessment. Ref. [8] f ormulated infeasibility certificates with respect to th e AC-OPF prob lem that are based on h ypersph eres which ce r tify a wide range of op erating points a- priori as insecure. Inspired by [8], our pr e vious work in [ 9] used such certificates to ge n erate large d atasets, r educing the inp ut space and decrea sin g computation time, while considerin g both N-1 security an d small-signal stability . Both work s [8] and [9] consider systems with up to 1 1 contro l variables. Instead of hy perspher es, this p aper pro poses the use of sep arating hyperp lanes, which, amo ng other important benefits, allows us to consider n u mbers of contr o l variables that are a t least an order of magnitude greater than pre vious methods (up to 12 5 in our test cases). Another popular approa c h to create such datasets is throu g h importan ce sampling, e.g., [1 0], [1 1]. I n power systems, how- ev er , the initial sam p ling space is largely un balanced, i.e., the volume of insecur e space is se veral o rders o f ma gnitude larger than the secure sp ace, an d, as we obser ved in [ 9], it can b e challengin g to obtain an a dequate number of secure samp les. In this work, we sho w how o ur proposed method can lead to a balanced d ataset, as it enables us to samp le from inside the secure space. A related strand of research uses histor ica l data that is enrich ed throug h sam pling m ethods such as co mposite modelling ap proaches an d v ine-copu las [1 2], [1 3]. Howe ver , this ca n neglect parts of th e secure space or might not capture abnorm al opera tin g r egions. T o create representative datasets for data- driven power system applications, we pro pose a c o mputation ally efficient method wh ich a) can d eal with high inpu t dimensionality (our test cases have up to 125 control variables), b ) provides a detailed description of the security bo undary , a n d c) creates balance d classes. W e apply this method to AC-OPF pro blems including N-1 security and uncertainty in power injections. The main contributions of our work are: 1) W e prop ose an infe asibility c ertificate based on sepa- rating hyp erplanes. This certificate is co mputed using conv ex relaxations of A C-OPF prob lems and c onsiders both N-1 security and uncertain ty . Compared to the hypersp h ere-based method pro posed in [ 8], our algo- rithm shows two key improvements: First, separ a ting hy- 2 perplane s allow the classification of substantially larger parts of the input space as insecure. Second, as these hyperp lanes form a co n vex polytop e, efficient methods to sample unif ormly fr om in side the remainin g unclassi- fied sp ace ar e available. Based on these, we propose an efficient alg orithm to m aximize the v olume of the input space classified a-prio ri as insecu re. 2) W e ev alu ate this algorithm on PGLib- OPF network s with up to 500 buses an d num b ers of control variables up to 1 25. Compar ed to initial normalized input space vol- umes of 1 (i.e., 10 0 ) based o n sp ecified control variable bound s, the infeasibility certificates reduce the unclas- sified inp ut space volumes significantly , with reduc e d volumes ranging fr om 10 − 2 up to 10 − 40 . 3) W e pr opose a c o mputation ally efficient method to cre- ate datasets for d ata-driven power system applications which can handle systems where the number o f co ntrol variables is at least one order of magnitude greater th an state-of-the- art meth ods (e.g., [9]). Com puting infeasi- bility certificates allows us to efficiently chara c terize the security boundar y in detail and sample f rom inside the secure space. W e create balanced d a tasets for PGLib- OPF networks u p to 50 0 buses an d train ne ural ne twork classifiers as an illustrativ e data-driven application. This pap e r is structured as fo llows: In Section II, we describe the A C-OPF pro blem including N-1 security and un certainty , and its conve x relaxation. I n Section III, we ou tline our pro- posed method o logy to create d atasets, in cluding the infeasibil- ity certificate, b oundar y descriptio n, an d sampling from inside the secure space. Section IV pre sen ts simulation results on PGLib-OPF networks up to 500 b uses. Section V conclu des. I I . O P T I M A L P O W E R F L OW F O R M U L A T I O N This section presen ts the A C-OPF formulations n ecessary for deri ving the data set creation methodolog y . In pa r ticular , we fo rmulate the N-1 security - constrained pr ev entive A C-OPF problem considering uncertainty in p ower injections, and its quadra tic conve x (QC) relaxa tion. For a detailed su r vey on A C-OPF and co n vex relaxations of the AC-OPF , th e rea d er is referred to [7], [14]. Here, fo r brevity , we b uild ou r fo r mu- lation upon the AC-OPF for mulation of [1 5] to facilitate the deriv ation of the QC r elaxation. W e use the QC r e laxation as it represents a good trade -off between comp utational co mplexity and tightness of the relaxation [1 5]. Note that the f ollowing deriv ations cou ld be readily extended via the m any o ther conv ex relaxation s of the power flow equations [ 7]. A. Security-Con strained A C-OPF under Uncertainty A power system is d efined by its set N of b uses. A subset of those buses, which are denoted by G , hav e a c o ntrollable generato r conn ected. A second subset de n oted by U , wh ich can be either genera tion or lo a d buses, are subjec t to uncertain power injectio ns. It is a ssum ed that all b uses of the power system are connected by a set ( i, j ) ∈ L of power lines from bus i to bus j . T o ensure the N-1 secur ity c r iterion d uring operation , we consid e r the potential outag e of a list o f c ritical candidate lines defined b y th e set C ⊂ L . Note that we define the first en try of C to correspo nd to th e intact system state { 0 } , i.e., no transmission line is outaged . The term L c denotes the set of intact p ower lines for ou tage c ∈ C . For th e intact system state, the set L 0 correspo n ds to the set L . The o ptimization variables in the security-co nstrained A C- OPF a r e the co mplex bus voltages V c k for each bus k ∈ N an d contingen cy c ∈ C , the comp lex power dispatch of generator S c G k for each bus k ∈ G and con tingency c ∈ C , and the uncertain com p lex power injection s S U k for each b us k ∈ U . The uncertain power injections do n ot change upon ou tage of system compon ents, i.e., S U = S c U , ∀ c ∈ C . W e assume that the uncertain reacti ve power injection Q U = ℑ{ S U } is determined th rough a fixed power factor c os φ in relation to the un certain activ e power in jection P U = ℜ{ S U } , i.e. , Q U = q 1 − cos 2 φ cos 2 φ P U . If the p ower facto r is not constant, then th e reactive power in jection Q U can be mod e lled as a separate variable, i.e., as a sep arate d egree of freedom in th e data set creation m ethod. The f ollowing co nstraints must be satisfied for the intact system and for each contingen cy c ∈ C : ( V min k ) 2 ≤ V c k ( V c k ) ∗ ≤ ( V max k ) 2 ∀ k ∈ N ( 1a) S min G k ≤ S c G k ≤ S max G k ∀ k ∈ G (1b) | S c ij | ≤ S max ij ∀ ( i, j ) ∈ L c (1c) S c G k − S D k + S U k = X ( k,j ) ∈L c S c kj ∀ k ∈ N (1d) S c ij = ( Y c ij ) ∗ V c i ( V c i ) ∗ − ( Y c ij ) ∗ V c i ( V c j ) ∗ ∀ ( i, j ) ∈ L c (1e) S min U k ≤ S U k ≤ S max U k ∀ k ∈ U (1f) θ min ij ≤ 6 ( V c i ( V c j ) ∗ ) ≤ θ max ij ∀ ( i, j ) ∈ L c (1g) The bus voltage magnitudes are con stra in ed in (1a) by upper and lower limits V min k and V max k . The superscript ∗ deno tes the complex conjugate. Similarly , the genera tors’ complex power outputs ar e limited in (1b) by upper and lower bo unds S min G k and S max G k . The inequ a lity con straints for complex variables are defined a s bounds on the rea l and imag inary parts. The apparen t power flow S ij on the line fro m i to j is u pper bound ed in (1c) by S max ij . The nodal complex p ower b alance (1d) including the load S D , generatio n S G and uncertain injections S U has to hold for eac h b us. The app arent power flow S ij on the line fro m i to j is defined in (1 e). The term Y denotes the admittance matrix of the p ower gr id . Con stra int (1f) mo dels minim u m and max imum bou nds S min U k , S max U k on the u ncertain in je c tio ns. The flo w o n the line from i to j is limited in (1g) b y a lower and upp er limit on angle d ifference s θ min ij and θ max ij , respectively . Please note that for most instan ces the following ho lds: θ min ij = − θ max ij . W e consider preventi ve action s in the security- constrained A C-OPF fo r mulation, i.e., the gen e rator set-poin ts remain fixed during an outage. As a result, we include the following linking constraints between the intact system state and the outaged system states: | V 0 k | = | V c k | ∀ k ∈ G , ∀ c ∈ C \{ 0 } (2a) P 0 G k = P c G k ∀ k ∈ G \{ slack } , ∀ c ∈ C \{ 0 } (2b) The active power dispatch is d enoted as P G , i. e ., P G = ℜ{ S G } . The first co nstraint sets the generato r voltage set- points | V k | o f the outaged system states to the values from 3 the intact system state. T he second constraint do es the same for the active power generation d ispatch, excluding the slack bus which com pensates the difference in active power losses. B. Quadratic Con vex (QC) Relaxation The QC relax ation prop o sed in [15] u ses conve x env elopes of the polar represen tation of the AC -OPF prob lem to relax the depend encies among voltage v ariables. As prop osed in [15], [16], an ad ditional auxiliar y matrix variable W c is intro duced for th e intact system state and each continge n cy c ∈ C , wh ich denotes the produ ct of the co mplex bus voltages: W c ij = V c i ( V c j ) ∗ ∀ c ∈ C (3) This allows refo rmulation of (1a), (1e), (1g), and (2a) as: ( V min k ) 2 ≤ W c kk ≤ ( V max k ) 2 ∀ k ∈ N (4a) S c ij = ( Y c ij ) ∗ W c ii − ( Y c ij ) ∗ W c ij ∀ ( i, j ) ∈ L c (4b) S c j i = ( Y c ij ) ∗ W c j j − ( Y c ij ) ∗ ( W c ij ) ∗ ∀ ( i, j ) ∈ L c (4c) tan( θ min ij ) ≤ ℜ{ W c ij } ℑ{ W c ij } ≤ tan( θ max ij ) ∀ ( i, j ) ∈ L c (4d) W 0 kk = W c kk ∀ k ∈ G , ∀ c ∈ C \{ 0 } (4e) The non-conve xity is encapsulate d in the voltage product (3). T o obtain a conve x r elaxation, the n on-conve x constraint (3) is removed f rom the o ptimization p roblem a n d variables for voltages, v c i 6 θ c i ∀ i ∈ N ∀ c ∈ C , and squared current flows , l c ij ∀ ( i, j ) ∈ L c ∀ c ∈ C , are added . The following conv ex constraints an d en velopes are intr oduced for the intact sy stem state and each contingen cy c ∈ C [15]: W c kk =  v 2 k  T ∀ k ∈ N (5a) ℜ{ W c ij } = D  v c i v c j  M  cos( θ c i − θ c j )  C E M ∀ ( i, j ) ∈ L c (5b) ℑ{ W c ij } = D  v c i v c j  M  sin( θ c i − θ c j )  S E M ∀ ( i, j ) ∈ L c (5c) S c ij + S c j i = Z c ij l c ij ∀ ( i, j ) ∈ L c (5d) | S c ij | 2 ≤ W c ii l c ij ∀ ( i, j ) ∈ L c (5e) The su p erscripts T , M , C, S de n ote conve x en velopes for the square, bilinear produ c t, cosine, and sine functio ns, respec- ti vely . The term Z ij denotes the line imp edance. Refer to [15] for the co mplete QC fo rmulation . The resultin g relaxation of the preventi ve secu rity-con strained A C-OPF under un c ertainty is a seco n d-ord er con e prog r am (SOCP) that minimizes an objective fu n ction, e.g., gener ation cost, subject to (1 b)–(1d), (1f), (2b), (4), and (5). I I I . M E T H O D O L O G Y T O C R E A T E D A TA S E T S The go a l of th e following m ethodolo g y is to create a dataset which maps op erating p oints describe d by th e inp ut vector x to a power sy stem security classification, e.g., secure or insecure. The dataset should be balanced b e tween secure an d insecure samples, take into consider ation th e degrees of f reedom of the system, and have a d etailed d e scr iption o f the security b ound- ary . The power system secur ity classification we consid e r is feasibility with respect to the N-1 security-con strained A C- OPF p r oblem under unc ertainty defined in (1) and (2). The resulting d ataset can be complemen ted with f urther assessment of dyn amic security criteria, e.g . small-signal stability [9]. The input vector x , i. e . , the contro l variables that d efine the r elev an t degrees o f freedom, is defined as follows: x =   P 0 G i | V 0 j | P U k   ∀ i ∈ G \ { slack } , ∀ j ∈ G , ∀ k ∈ U (6) Using the input x , all other states in the AC-OPF p roblem can be determ ined b y solving the non-linear AC power flow equations. The minimu m and m aximum bound s on inpu t vector x max and x min are de fin ed in (1a), ( 1b), and (1f). Please note that in the formulation of the QC relaxation, we use the variable v 0 j instead of the eliminated variable | V 0 j | . The main ch allenge in crea tin g a represen tati ve an d bal- anced dataset is the large number of contro l variables. The dimensiona lity of the in put vector x grows sub stantially with increasing sy stem size. For in stan ce, the IEEE 118 -bus system has 72 control variables, i.e., the dimensionality | x | is 7 2. A na¨ ıve app roach to create a d ataset would be to sample with a prespecified discretizatio n in ter val, e.g., b y spec if ying 10 steps in each dim ension o f the con trol variables, x 1 , x 2 , x 3 , . . . . For the 118 -bus system, this w ould requ ire power flow solutions for 10 72 operating poin ts, which is comp u tationally intractable. Further, as we will empirically show in Section IV -C, large parts of the in put space x ∈ [ x min , x max ] are infe a sible. As a result, identify ing secure samples by na¨ ıvely sampling fro m the entire input space is not possible for larger test cases. T o ad dress these challen ges, we pr esent an efficient method for c reating such da tasets. First, to a-priori classify large parts of the inp ut space as insecure, we prop ose an in feasibility certificate based on separating hyperp lanes in Section III- A. Focusing on th e u nclassified regions, we the n characterize the security bo undary in detail in Section III -B. Finally , we sample inside the secure space in Section III-C. A. Constructing Infeasibility Certificates W e prop ose an infea sib ility certificate which can a-priori certify region s in which th e non -conve x security-co nstrained A C-OPF prob lem und er u ncertainty is infeasible. This exploits the fo llowing p roperty of a convex relaxation : if a re laxation is inf e asible for a giv en operating point, th e origin al non- conv ex problem is a lso guara nteed to be infeasible for that operating point. The pro posed in feasibility c ertificate has th ree compon ents: First, we employ bound tightening to tighten both the QC relaxation and th e in put bou nds; th is better approx imates the secure region , while also redu cing the sample space. Second, we propo se an in feasibility certificate based on separ ating hyper p lanes. Third, we pre sent an efficient algorithm to maxim ize the input region classified as in feasible. 1) Boun d T ig htening Algorithms: Th e tightness of the QC relaxation relies on the tightness of the en velop e s used in (5) including the en velopes on co sine and sine term s. These in turn d epend on the tightness of the bou n ds on the voltage magnitud es a n d ang le differences. The go a l of bound tight- ening is to iterati vely tighten voltage magnitu des and a ngle differences, and, as a result, obtain a tigh ter r elaxation. In the co ntext of our work, the benefits of boun d tigh tening are twofold: First, it tightens the QC relax a tion, i. e., shrinks its 4 ˆ x + + x ∗ | − → n | = R ∗ Hypersphere Hyperpla ne Input s pace x ∈ [ x min x max ] Feasible space (QC relaxat ion) Feasible space (A C-OPF) Fig. 1. Illustra ti ve example of the dif ferences between the infea sibility certific ates using hypersphere s and hyperplanes. For a giv en infeasible point ˆ x , the closest point x ∗ is computed which is feasi ble to th e QC relaxat ion. The normal vector − → n is perpendic ular to the feasible space of the QC relaxation. The feasible space of the non-con vex AC-OPF problem is contained w ithin that of the con vex QC relaxat ion. All points inside the hypersphere or all points that are on the left side of the hype rplane are guarantee d to be infeasible with respect to both the QC relaxation and the non-con ve x AC-OPF problem, respect iv ely . Note that the sets are not drawn to scale. feasible space, makin g the in f easibility certificate based o n separating hyperp la n es m ore e ffecti ve, and seco nd, it allows us to directly tighten the bound s on the in p ut vector x . W e use two bound tigh tening algorithms from the litera- ture: First, we rely on a comp utationally lightweight bound tightening tec hnique for the bran ch ang le differences θ min ij and θ max ij in (1g) from [ 17]. Secon d, we u se an optimization- based boun d tighte n ing alg orithm from [18] which tightens the voltage ma g nitude bound s at each bus V max , V max in (1a), and fur ther tightens the an gle differences for each line θ min ij and θ max ij in (1 g). For this pur pose, we iteratively solve c o n vex optimization problems to calcu late th e ma ximum and mini- mum v alues that the optimization variable under study , i.e., a voltage magnitude or a voltage angle d ifference, can obtain in the relaxed pro b lem. Note that b y tightenin g on e variable bound , it may be possible to f u rther tighten a previously tightened bound . This proce dure can be executed fo r a defined number of iterations or until a fixed point is reached. As a final step in the bou nd tighten in g, we com pute th e tig h tened bound s for the in p ut vector x , i.e., the bounds on active power o f gen erators and uncertain injections. All inpu ts x which are outside the tighten ed min imum an d max imum inp ut bound s x BT , min , x BT , max are gu aranteed to be infeasible with respect to the non-conve x AC-OPF p r oblem. W e calculate the volume of the rem a ining un classified inp ut space volume V B T , normalized by the originally specified bound s on x : V B T = Y k ∈X x BT ,max k − x BT ,min k x max k − x min k (7) The input set X is defined as X : {G \{ slack } , G , U } . 2) Sepa rating Hyperplanes: W e next propose a n infeasi- bility certificate based on separating hyperp la n es. Consider a particular oper ating p oint ˆ x that is infeasib le with respec t to the no n-conve x security-c o nstrained A C-OPF . W e so lve the fo llowing optimization pro blem to compute the clo sest dispatch x ∗ which is f easible to the c on vex QC r elaxation: min x,S,S U ,S G ,v , θ ,l,W,R R (8a) s.t. (1b)–(1d) , (1f) , ( 2b) , ( 4) , (5) , ( 6) (8b) s X k ∈X ( x k − ˆ x k ) 2 ≤ R (8c) If the obtained radius R ∗ is grea te r th an zero, i.e. , th e operatin g point ˆ x is in feasible with r espect to the relaxation , no operating point x exists wh ich is closer to ˆ x than the ob tained p oint x ∗ . This prop e r ty has been used in [8] to con struct in feasibility certificates in the form of hyp erspheres and ellip ses by as- signing different weights to th e compo nents in (8c). Here, we propo se to use hyper planes as infeasibility ce rtificates in order to significantly enlarge the volume classified as infeasible: Pr opo sition 1: For a given infeasible po int ˆ x , if th e solu tion to (8) yields a non-zero radius R ∗ and optimal solution x ∗ , all vectors x which fulfill the following criterio n are infeasible with respect to the A C-OPF constraints (1) and (2): − → n T ( x − x ∗ ) < 0 (9) The normal of the hyperp la n e is defined as − → n := x ∗ − ˆ x and the operator T d enotes the transpose. Pr oof of P r oposition 1: Proof by contrad iction: Assume there exists a feasible point ˜ x that is inside the region classified as infeasible by the hyp erplane: − → n T ( ˜ x − x ∗ ) < 0 . As th e feasible space of o ptimization prob le m (8) is con vex, it must hold that any linear com bination between ˜ x and x ∗ is also feasible: λ ˜ x + (1 − λ ) x ∗ , λ ∈ [0 , 1] . Then, there exists a point ˜ x ∗ = λ ∗ ˜ x + (1 − λ ∗ ) x ∗ which has a radiu s ˜ R ∗ to the initial in feasible point ˆ x that is smaller than R ∗ . Since the optimization p roblem (8) is co n vex, we obtaine d the globally optimal solution x ∗ with the smallest rad ius R ∗ . As a re su lt, there cann ot exist an input ˜ x that has a smaller radius tha n R ∗ . W e have sh own b y contradictio n that ther e cannot exist a feasible point ˜ x that is inside the r egion classified as inf easible by th e hype rplane. The infeasibility cer tificate with respect to the non - conv ex A C-OPF problem (1) and (2) follows fro m the pr operty that in feasibility with respect to the QC relaxation co nstraints (8b) is sufficient to en sure infea sib ility with respe c t to (1) and ( 2). Altern a tively , we can show that Proposition 1 is tru e by taking the first-order T aylor expansion of con stra in t (8 c) at the o ptimal solution x ∗ . For co n vex sets, first-ord er T aylor exp a nsions o f nonlinear constraints are always sep arating hyperplan es [19]. An illustrative comparison of both infea sibility certificates is shown in Fig. 1. By solvin g the same o ptimization prob lem, it is evident that the infeasibility certificate based on hype r planes is able to classify significantly larger spaces as infeasible . This is quantitativ ely analy sed throu g h simulation studies in Section IV -B. 3) An Efficien t Algorithm to Minimize th e Unc lassified Input Spac e: Using the in feasibility certificate, we propo se a n efficient algor ithm to max im ize th e p o rtion o f the input space that can b e classified a-priori a s inf easible. Our algor ithm relies on an insig ht related to the hyper planes: togethe r with the initial input space restriction, subsequ ent h yperplan es form a c on vex poly to pe which can be described as Ax ≤ b . W e 5 Algorithm 1 Computing Infeasibility Certificates 1: Run bou n d tighten ing and obtain x BT , min and x BT , max 2: Set iteration count: k ← 0 3: Initialize unclassified region A (0) x ≤ b (0) : A (0) := [ I | x |×| x | − I | x |×| x | ] T b (0) : = [( x BT , max ) T ( x BT , min ) T ] T 4: while k ≤ N 1 do 5: draw random x ( k ) from inside A ( k ) x ≤ b ( k ) 6: solve (8) with ˆ x := x ( k ) and obtain x ∗ 7: if R > 0 then 8: reduce unclassified region by adding hyperpla n e: A ( k +1) = [( A ( k ) ) T − → n ] T b ( k +1) = [( b ( k ) ) T − → n T x ∗ ] T 9: end if 10: k ← k + 1 11: end while . can write the r ow of A and entr y in b correspond ing to the hyperp lane in (9) as A k := − → n T and b k := − → n T x ∗ . Efficient methods to sample u niformly from inside a conve x p o lytope are available, e.g., “Hit-an d-Run” sampling [20]. This allows us to iterativ ely construct hy perplanes while sampling only inside the currently u nclassified region. Thus, the hyper plane certificates facilitate a sig n ificant improvement on the “rejec- tion” sampling a p proach used with hy perspher e certificates in [8], [9]. The steps of the algorithm to com pute in feasibility cer tifi- cates are detailed in Algor ith m 1. W e start with a description of the conve x polytope restricted to the tighten ed inpu t boun ds. W e iteratively samp le unifo rmly f rom in side the c o n vex poly- tope and add identified hyperplane s until w e reach an upper iteration limit of N 1 samples. This ensures that only samples which have n ot yet been classified as infe asible b y previously added hyp erplanes ar e considered in optimization pr oblem (8). In Section IV -C, we will demonstrate the perfor mance of this algorithm on a range of PGLib-OPF networks u p to 500 b uses b y calcu la tin g the remaining un classified volume as th e volume of th e co n vex poly to pe A ( N 1 ) x ≤ b ( N 1 ) . This shows substantial reductio ns of unc la ssified input space volumes. An alter n ativ e ap p roach to the proposed algorithm using separating hy p erplanes could b e to directly construc t a linear ou ter ap p roximatio n of th e conv ex feasibility set defined by ( 8b). This cou ld be ach ie ved by applying e.g. “Hit-and - Run” sampling to (8b). This is subject of our futur e work. B. Security Bound ary I dentificatio n After com puting the infeasibility certificates, we perform sampling a n d directed walks, similar to [9], to obtain a d etailed description o f th e security boun dary . For this p urpose, we first u niformly dr aw a large numb er N 2 of sam p les f r om the conv ex polytope describing the remaining unclassified input region: A ( N 1 ) x ≤ b ( N 1 ) . F or ea c h sample, we first run A C power flo ws for the intact and the o utaged system states an d check if any o f the con straints in (1) ar e violated. If n ot, we add the cu rrent p oint to the dataset as a fe a sible p oint, otherwise as an inf easible p oint. I f constraints ar e v iolated, we run additiona l AC power flo ws for which we enforce the reactive power limits o f gener ators, i.e. , if any generato r violates its reactive power limit it is co n verted f rom a P V to a P Q b us in th e power flow . This is based on th e observation that reactive power limits are o f ten the only constra in ts violated . If the o btained power flow solution s satisfy all constraints in (1), the v oltage-ad justed point is added to the dataset as feasible sample, i.e., the voltage set-points of generators in x are up dated accor dingly . If bo th stages are not feasible, we solve th e following non- con vex optimizatio n p r oblem which computes the closest feasible dispa tch to the non-c o n vex A C- OPF problem in (1) and (2): min x,V , S U ,S G ,S,R R (10a) s.t. (1) , (2) , (6) , (8 c) (10b) W e add the o btained loc a lly optimal po int x ∗ to the dataset as feasible poin t. W e repeat this pr ocedure f or all N 2 samples and obtain as a result a detailed security boun dary d escription. C. Sampling fr om Insid e the Secur e Space T o obtain a more detailed descr iption of the en tire secur e space, we fit a multiv ar iate no rmal d istribution N to th e feasible points obtaine d. For this purpose, we estimate both the m ean µ and th e covariance matrix Σ from the fe a sib le data points. T o bias the sampling to wards in side the boundary , we reduce th e mag nitude of all entries of th e covariance matrix, i.e., Σ red = s red · Σ , by a co nstant scaling factor s red < 1 . W e draw a large numb er , denoted N 3 , of sam ples fro m N ( µ, Σ red ) . For each of these samples, we first ru n A C po wer flows for the intact and th e outaged system states, check f easibility with respect to all A C-OPF constrain ts, and add the sample with the corresp o nding classification to the dataset. If the sample is infeasible, we run a second ro und of A C p ower flows in which we enforce the gen erators’ rea c ti ve limits and again evaluate the f easibility with respect to all A C-OPF constraints. If the sample is feasible, we add it to the dataset with the gener- ator voltage set- points adjusted accordin gly . Our simulations indicate th a t sampling fr o m a multi variate no rmal distribution N ( µ, Σ red ) results in identification of feasib le samples inside the secure space. W e did not ob serve im p rovements by fitting a Gaussian mixture model. I V . S I M U L A T I O N A N D R E S U LT S W e analy se the perfor mance o f our pr oposed methodolo gy for a r a nge o f test cases from the PGLib-OPF networks. Fir st, we compare th e proposed in feasibility certificate b ased on separating hyperplan es with the certificate based on hyper- spheres from [ 8]. Second, we compute th e volume of the unclassified inp ut space using the infeasibility certificates a n d show substantial reductio ns. Third, we cre a te balanced da tasets and demon strate their applicability using an illustrative data- driven application . A. Simulation Setup In the fo llowing, we first evaluate our proposed methods on 13 PGLib-OPF networks (v 19.05) [21] up to 500 buses 6 for which w e do not con sider N-1 security and uncertainty , i.e., we use the test cases as specified in [2 1]. Seco nd, we use two test cases for which we includ e both N-1 security and uncertainty . W e use ca se39 epri and case16 2 ieee dtc with the following line conting e ncies C = { 0 , 7 , 22 , 24 , 36 , 43 } and C = { 0 , 6 , 8 , 2 4 , 50 , 128 } , resp e cti vely . These lin es co rrespond to the f ollowing b us pairs { – , 3-18, 12-1 3, 14-15, 22-23 , 26- 28 } and { –, 2 - 7, 3- 14, 8 -13, 16-17 , 50-12 5 } , respectively . W e assume th e same parameter s for the outaged system state as f o r the intact sy stem state. Furthe r more, we place three wind farms with rated power of 500 MW and con sider three uncertain loads with ± 50 % variability , i.e. , a to tal of six uncer ta in power injections, at buses U = { 3 , 2 1 , 2 7 , 4 , 2 5 , 2 8 } for case39 epri and U = { 60 , 90 , 145 , 3 , 8 , 52 } fo r ca se162 ieee dtc . For all uncertain injections, we assume a power factor cos φ = 1 . Note that all inputs x ar e n ormalized with re sp ect to th e ir maximum and minimum lim its as specified in [2 1], i.e., if x has d imension | x | , then x ∈ [0 , 1] | x | . Th is nor malization step is standard practice f or many data-driven applicatio ns including neural networks and improves perf ormance [2 2]. For bo th A C power flow an d A C optimal power flow co m- putations, we rely on M A T P O W E R [23] with the I P O P T solver for AC-OPF pro blems [24]. For th e bo und tighten ing, we use the implementation s in [17], [18]. Note tha t we adap ted the implementatio ns in bo th [17] an d [ 1 8] to include u ncertainty in power injections by modeling them as gen erators with active power limits correspon ding to the defined uncertainty set and no r eactiv e power ca p ability , i.e., the lower and upper reactive power limits are set to ze ro (as the po wer factor is assumed to be 1, cos φ = 1 ). W e only tighten the bound s of the intact system state, i.e., the bound s of the o utaged system states are not tig h tened, and we ru n the optimization-based bound tighten ing fo r up to thr ee iter ations. Extension of these toolboxes to the full N-1 case is a dire ction for future work. W e use MOSEK [25] to solve the QC relaxation. T o ap proxim ate the volumes of the conve x polytopes de- scribing the remaining u nclassified input space, we use a volume approx im ation too lbox in C++ [26] which handles floating poin t pr ecision issues. Note that an exact volume computatio n is considered intracta b le fo r dimensiona lity 10 o r higher [ 26]. The relativ e appro x imation error thresh o ld is set to be less than one order o f magnitude, which is sufficiently accurate for our purp o ses since we compute volumes of spaces se veral orders of magnitud es smaller than the initial volume. B. Comparison of Infeasiblity Certificates W e compare the infeasibility ce rtificate b ased on hyp er- spheres p r oposed in [8] with the inf e a sibility certificate b ased on h yperplan es proposed in Section III-A2. The main m etric for comparison is the volume of the rem aining u nclassified space after applying the infeasibility certificates. W e co nsider case39 epri a nd case162 ieee dtc . W e use the QC relaxation for both certificates with all the bo unds tightened as described in Sectio n I II-A1. W e only consider the activ e power g ener- ation in the inp u t variables x and we d o not consider N-1 security o r uncertainty . For the certificate b ased on hyper- planes, we follow the algorithm outlined in Algorith m 1. For 0 10 20 30 40 50 10 − 5 10 − 4 10 − 3 10 − 2 10 − 1 10 0 Iterat ions Unclass. V olume Hypersphere s (39 bus) Hyperpla nes (39 bus) Hypersphere s (162 bus) Hyperpla nes (162 bus) Fig. 2. For case 39 epri and case162 iee e dtc , we compare the remaining unclassifie d volume betwe en an infeasi bility certific ate based on hyperspheres from [8] and the proposed certificat e in Section III-A2 based on hyperplanes. Note tha t these result s cannot be directly compared to those of T able I, as here only the acti ve generator set-points are assumed as degrees of freedom. In T able I, gene rator voltage set -points and uncer tain inje ctions are also considere d. the certificate based on hyperspheres, we assume th at in each iteration we draw a random samp le from the entir e input space, and if it is infeasible we compute the closest feasible input x by so lving (8). If the distance is non- zero, we hav e ob tained an inf easibility certificate. W e use Monte Carlo sampling with 10 6 samples to estimate the volume of the unclassified space not covered by the hy perspheres. Fig. 2 sho ws the variation of the unc lassified inpu t v olume with up to 50 iterations. W e make the following o bservations from the results. First, the hyperpla n e certificates shrink the unclassified r egion by f our ord ers of magnitu de more tha n the hypersp h eres, i.e., th e unclassified volumes evaluate to 10 − 5 versus 10 − 1 compare d to the initial unit hypercube’ s no rmal- ized v olume of 1 . Sec o nd, the algorithm u sing hyperplanes requires significantly fewer iteratio n s. After the first iteratio n, the hype r planes classify a substan tially larger space as infea- sible than the hypersp heres after 5 0 iterations. The reasons for this are twofold: First, as evident in Fig. 1, ce r tificates based on separating h yperplan es cover a larger volume than hypersp h eres for the same sample and second, th e hyper planes enable the use of efficient me th ods for sampling uniform ly from inside the associated con vex polytop e [ 20], [26]. C. Estimating Uncla ssified V o lumes for PGLib-OPF Networks In the following, we compu te infeasibility ce r tificates an d the volume of the re m aining unclassified inp u t space for a range of test cases. For this purp ose, we run Algorith m 1 with the numbe r of iterations N 1 set to 1000 . W e ev aluate the remaining estimated volume for the b ound tightenin g V B T accordin g to (7) and for the separating hyp erplanes describ ed as a co n vex polytope by ru nning the volume approx imation algorithm in [26]. In T able I, the dim ensionality | x | , th e reduced unclassified v olume V B T after bound tightening, the number o f h y perplane s | H P | , and the r e duced u nclassified volume V H P enclosed by the sepa r ating hyp erplanes is listed. Note that both volumes are defined with respect to th e unit hyperc u be x ∈ [0 1 ] | x | normalized by the origin al p ower system limits with volume 1 , i.e . , 10 0 . W e make several observations. First, the bound tightenin g results in a moderate reduction in in put dimensionality of se veral orders of mag n itude ( 10 − 1 to 10 − 4 ) for m ost test cases. Second, th e inf easibility ce r tificates based on hyper- planes en able further substantial reductio ns in the unclassified 7 T ABLE I U N C L A S S I FI E D I N P U T V O L U M E S F O R P G L I B - O P F N E T W O R K S Po wer system case | x | V BT | H P | V H P − log 10 ( V ) | x | A C-OPF without N-1 securit y and without uncertai nty case3 lmbd 4 6.3e-02 28 3.3e-02 37.0% case5 pjm 7 1.0e +00 99 6.9e-03 30.9% case14 ieee 6 2.4e-01 54 6.9e-04 52.7% case24 ieee rts 20 9.2e-01 184 2.3e-06 28.2% case30 ieee 7 6.2e-03 61 8.8e-06 72.2% case39 epri 19 9.9e-02 203 7.0e-08 37.7% case57 ieee 10 3.8e-02 231 4.9e-06 53.1% case73 ieee rts 62 1.0e +00 608 6.1e-16 24.5% case118 ieee 72 1.7e-02 1000 1.6e-17 23.3% case162 ieee dtc 2 3 6.1e-04 371 1.5e-11 47.1% case200 tamu 69 9.3e-01 1000 6.0e-11 14.8% case300 ieee 125 1.0e-12 1000 3.4e-40 31.6% case500 tamu 111 8.6e-02 1000 5.4e-26 22.8% A C-OPF consideri ng N-1 security and uncertai nty case39 epri 25 2.6e-01 271 2.0e-05 18.8% case162 ieee dtc 2 9 2.2e-04 394 6.0e-10 31.8% Median all cases 23 8.6e-02 271 7.0e-08 31.6% volume. As a result, the total unclassified v olume compared to the unit h ypercu be is redu ced between 2 and 40 orders o f magnitud e ( 1 0 − 2 to 10 − 40 ). The med ian of the un classified volume is 1 0 − 8 . T his means that in order to identify o ne sample inside th e unclassified volume, we would h av e to unifor m ly draw 10 8 samples from the origin al bounds on the input x . This highlights the n ecessity of fir st co mputing the infeasibility certificates to b e a ble to identify the secur e space. The n umber of hyperp lanes is below 100 0 fo r most test cases, indicating that Algo rithm 1 has obtain e d a g ood estimation of the u nclassified volume. For th e four test cases f or which 1000 hyperp lanes are added, th e unclassified input volume cou ld be further reduced by increasing N 1 . T o allow for comp arability between test cases with different number of degrees of fr eedom, we p ropose to use a metric defined as − log 10 ( V ) | x | . The metr ic is m otiv ated as follows: If one wants to sample 1 0 steps in each dime n sion, i.e., 10 | x | , then this metric quan tifies by h ow much the exponent is reduced . Note that the value obtained in percent is not the dimensiona lity reduction itself but relates to the redu ction in the or ders of magnitu des of the dim ensionality . This value is between 14.8 % and 72 . 2% for all test cases, showcasing the general applicability of the p roposed infeasibility certificate for A C-OPF problems. D. Dataset Cr eation for PGLib-OPF Networks W e create d atasets o f o p erating p oints classified based on their feasibility with respect to A C-OPF p r oblems including N-1 security an d un certainty . T o this end , we first draw a number of samp les N 2 = 10 4 from the inside of the re- T ABLE II C R E A T E D D AT A S E T S F O R AC - O P F P RO B L E M S Po wer system case Boundary Inside Overa ll Overa ll (MVND) secure point s N 2 = 10 4 N 3 = 10 5 A C-OPF without N-1 securit y and without uncertai nty case3 lmbd 69.5% 36.5% 40.6% 114’389 case5 pjm 68.6% 69.4% 69.3% 125’432 case14 ieee 73.3% 59.0% 61.0% 147’047 case24 ieee rts 66.8% 44.3% 48.7% 131’158 case30 ieee 75.0% 50.2% 54.0% 124’944 case39 epri 57.2% 29.9% 33.9% 154’635 case57 ieee 58.9% 35.2% 38.9% 150’865 case73 ieee rts 63.9% 51.1% 52.7% 222’730 case118 ieee 53.2% 47.0% 47.6% 209’996 case162 ieee dtc 50.0% 40.1% 41.7% 129’165 case200 tamu 50.2% 36.6% 38.1% 177’023 case300 ieee 50.0% 32.6% 34.7% 163’087 case500 tamu 50.0% 35.4% 37.1% 174’774 A C-OPF consideri ng N-1 security and uncertai nty case39 epri 58.2% 78.2% 75.2% 139’756 case162 ieee dtc 50.0% 17.9% 23.2% 121’358 A vera ge all cases 59.7% 44.2% 46.5% 152’424 maining u nclassified volume describ e d in T ab le I and obtain a detailed security bo undary description following th e approach in Section III - B. W e fit a multivariate normal distribution w ith s red = 0 . 25 and classify N 3 = 10 5 samples as secu re or insecure following the ap proach in Section III-C. In T able II, we list the characteristics of the obtained datasets. First, note that in the bound ary identification stage, if th e perc e ntage of secure points is ab ove 50%, then sampling directly from the remaining unclassified volume results in ide n tifying secure operating p o ints. This is th e case for the m ajority of test c a ses, demonstra tin g tha t the infeasibility certificate is able to provide a tig ht app roximation of th e secu re spaces of n on-convex A C- OPF pr oblems. For the test cases wh ere the sampling did not fin d any secur e samples, the numb er of iteration s for the feasibility certificate could b e enlarged or other relaxations such as m oment-ba sed r elaxations described in [7] could be used to furth e r reduce th e unclassified space in Algo rithm 1. Second, the results show that sampling from a multiv ariate normal d istribution fitted to th e b ounda r y samples results in id entification o f a large nu mber of secure samples. The resulting datasets are well b alanced with on a verage 4 6.5% secure samples. Note that this is an impo rtant m etric for the successful app lication o f data-driven methods such as neu ral networks [1]. The num ber of overall points is depend ent on the number of ad ditional f easible samp les identified b y en forcing the g enerators’ reactive power limits in the AC power flo ws and differs between the test cases. Regarding the computatio nal tractability , all simulations were carried o u t on a lapto p and the dataset creation for the largest test c a ses too k a fe w hour s, with the m o st c o m- 8 T ABLE III T E S T S E T AC C U R AC Y O F N E U R A L N E T W O R K C L A S S I F I E R S Po wer system case Full datase t Only boundary A C-OPF without N-1 securit y and without uncertai nty case14 ieee 78.2% 60.5% case39 epri 74.6% 38.5% case162 ieee dtc 84.4% 49.8% A C-OPF includi ng N-1 security and uncertai nty case39 epri 81.0% 80.4% case162 ieee dtc 93.4% 31.9% putationally intense task being the A C-OPF e valuations in the bounda r y identificatio n an d the optimization-b ased bou nd tightening [18]. By usin g hig h-perfo rmance computin g and parallelizing b oth the bou n dary identification a nd the A C power flow co mputation s, we expect that our approac h can scale to systems with thousan d s of buses. Th e nu mber of sam- ples ch osen for each stage of the dataset creatio n metho d n eeds to be adjusted for the da ta- driven applicatio n at hand , and depend s amo ng o ther factors on th e proble m dimen sio nality , the chosen classifier , an d the desired prediction accuracy . A common a p proach is to train an d ev aluate the performan ce of a data-driven metho d on datasets of different sizes. E. T raining Neural Network Classifiers As an illustrati ve data-dr i ven ap plication, we ev aluate the perfor mance of a neural network classifier trained o n several of th e cr eated datasets. The neural network predicts a binar y classification, i.e., whether th e input x be longs to the class “secure” or “insecure ”. W e choose neural network structures with fi ve h idden layers wh ere the numbe rs o f neuron s o f each hidden layer selected to b e 10 times the inpu t dimension | x | . W e split the d ataset into a training set consisting of 85% of all samples and a test set of th e r emaining 15 %. N o te th at the classifier has no informatio n of th e test set during training, and its perf ormance is evaluated on the test set only . This giv es a metric for h ow well the classifier generalizes to unseen data. W e train the neural networks u sing T ensorFlow [ 22] with standard training parame te r s and 25 0 epochs. T able I II shows the test set accu racy , i.e., the share of correctly predicted labels f o r th e test set. First, we use 8 5% of the fu ll dataset for trainin g and 15% of the f ull dataset for testing. Second , we on ly use the bo u ndary samples f r om Section III-B as trainin g data and then test on 15% of the full d ataset. T his gives us an estimation of th e b enefit of the additional sampling from the fitted multivariate d istribution in Section III-C. W e observe that the neu r al n etwork classifier is able to generalize from the training to the test set an d achieve high accur acy whe n using the fu ll dataset. T o further increase the classification accuracy , deep er neu ral network s or a deep autoenco der to identify lower-dimensional features cou ld be used. W e observe that only r elying on the boun dary samples for prediction is n ot suf ficient for most test cases, higlig h ting the impor tance of obtain ing a representative dataset. V . C O N C L U S I O N Successful app lication of d ata-driven methods in power sy s- tems requires d atasets o f sufficient size, covering a wide range of operating po ints. Creatin g a dataset that character izes the security boun dary and sufficiently covers b o th secure and inse- cure o perating regions is a h ighly computatio nally d emanding task, even for medium-sized systems, as we showed in [9 ]. In this paper , we prop o se an effi cient meth od to cr eate such datasets. W e focus on AC-OPF f e asibility and N-1 security , as any o perating point sho uld first satis fy static security criteria. Future work will extend this to in c lu de dynamic security criteria, similar to [9]. W e de velop an infeasibility certificate based on separating hyp erplanes wh ich is able to classify large portion s of the input space as insecu re. W e show th at the infeasibility certificates reduce the unclassified input space volumes significantly , by up to 10 − 40 compare d to an initial normalized input space volume of 1 (i.e., 1 0 0 ) based on defined control variable bou nds. Althoug h the secure ope rating region is a very small po rtion of th e origin al input space, our m e th od is ab le to pro duce balanc e d datasets of secure and insecure operating p o ints, a prop erty desired fo r successful applications of data-d riv en meth ods. As an illu strati ve application, we used th e g e n erated datasets to assess the perf ormance of neural network classifiers. Future work is directed to wards (i) utilizing conve x r e strictions from [27], [2 8] to characterize secure spaces and (ii) exploiting hig h-perfo rmance computin g. A C K N O W L E D G E M E N T S The author s would like to thank the two reviewers for their insightful comments. R E F E R E N C E S [1] L. A. W ehenkel , Automatic learning tech niques in power systems . Springer Science & Business Media, 2012. [2] B. Donnot, I. Guyon, M. Schoenaue r , A. Marot, and P . Panciatic i, “Fast power s ystem security analysis with guided dropout, ” Europea n Symposium on Artifici al Neural Networks, Computatio nal Intell igen ce and Machi ne Learning (ESA NN) , April 2018. [3] J. L. Cremer , I. Konstant elos, S. H. Tinde mans, and G. Strbac, “Data- dri ven power system operatio n: Exploring the balance between cost and risk, ” IEEE T ransactions on P ower Systems , vol . 34, no. 1, pp. 791–801, Jan. 2019. [4] L. Halilbaˇ si ´ c , F . Thams, A. V enzke, S. Chatzi vasi leiadi s, and P . Pinson, “Data-d ri ven security-constra ined A C-OPF for operations and markets, ” in 20th P ower Systems Computation Confer ence (P SCC) , Dublin, Ire- land, June 2018. [5] J.-M. H. Arteaga, F . Hanc harou, F . Thams, and S. Chatzi v asilei adis, “Deep learning for power system security assessment, ” in 13th IEEE P owerT ech , Milan, Italy , June 2019. [6] L. Duchesne , E. Karangel os, and L. W ehenk el, “Using machine learn ing to enable probabilisti c reliabil ity assessment in operation planning, ” in 2018 P ower Systems Computation Confere nce (PSCC) , Dublin, Ireland, June 2018. [7] D. K. Molz ahn and I. A. Hisk ens, “A Surve y of Relaxations and Approximati ons of the Powe r Flo w E quatio ns, ” F oundations and T re nds in Electric Ener gy Systems , vol. 4, no. 1-2, pp. 1–221, F eb. 2019. [8] D. K. Molzahn , “Computing the feasible spaces of optimal power flow problems, ” IEEE T ransactio ns on P ower Systems , vol. 32, no. 6, pp. 4752–4763, Nov . 2017. [9] F . Thams, A. V enzk e, R. E riksson, and S. Chatzi v asileia dis, “Effic ient databa se generatio n for data -dri ven securit y assessment of po wer sys- tems, ” IEE E T ransactions on P ower Systems , vol. 35, no. 1, pp. 30–41, Jan. 2020. 9 [10] V . Krishnan, J. D. McCalley , S. Henry, and S. Issad, “Ef ficient database generat ion for dec ision tree based power system securit y assessment, ” IEEE T ransacti ons on P ower Systems , vol. 26, no. 4, pp. 2319–2327, Nov . 2011. [11] C. Hamon, M. Perninge , and L. S ¨ oder , “ An importance sampling techni que for probabil istic security assessment in power systems with larg e amounts of wind power , ” Electri c P ower Systems Resear ch , vol. 131, pp. 11–18, 2016. [12] Mingyang Sun, I. Konstan telos, S. T indemans, and G. Strbac, “Eval u- ating composit e appr oaches to modelli ng high -dimensional stochastic v ariable s in powe r systems, ” in 19th P ower Systems Computation Confer ence (PSCC) , Genoa, Italy , June 2016. [13] I. Konsta ntelos, M. Sun, S. H. T indemans, S. Issad, P . Panc iatici , and G. Strbac, “Using vine copulas to generate representa ti ve system states for mac hine learnin g, ” IEEE T ransactions on P ower Syste ms , vol. 34, no. 1, pp. 225–235, Jan. 2019. [14] M. B. Ca in, R. P . Oneill, and A. Castillo, “History of optimal po wer flo w and formulations (OPF Pa per 1), ” F eder al Energy Re gulatory Commission , vol. 1, pp. 1–36, 2012. [15] C. Coffri n, H. L. Hijaz i, and P . V an Hentenryc k, “The QC relaxation: A theoret ical and computat ional study on optimal po wer flow , ” IEEE T ransactions on P ower Systems , vol. 31, no. 4, pp. 3008–301 8, July 2016. [16] J. L av aei and S. H. Lo w , “Zero duality gap in optimal power flo w problem, ” IEEE T ransactions on P ower Systems , vol . 27, no. 1, pp. 92–107, Feb . 2012. [17] D. Shchetinin , T . T . De Rubira, and G. Hug, “Effici ent bound tighteni ng techni ques for con vex relaxations of A C opti mal power flow , ” IEEE T ransactions on P ower Systems , vol . 34, no. 5, pp. 3848– 3857, Sept. 2019. [18] K. Sundar , H. Nagaraj an, S. Misra, M. Lu, C. Coffrin, and R. Bent, “Optimiz ation-ba sed bound tightening using a strengthened QC- relaxa tion of the optimal power flow problem, ” , Sept. 2018. [19] S. Boyd, S. P . Boyd, and L. V andenber ghe, Con vex optimization . Cambridge U nive rsity Press, 2004. [20] D. P . Kroese, T . T aimre, and Z. I. Botev , Handbook of Monte Carlo methods . John W iley & Sons, 2013, vol. 706. [21] IEEE PES T ask Force on Benchmarks for V alidation of Emergin g Po wer System Algorithms, “The power gri d libra ry for benchmarking A C opti mal po wer flo w algorith ms, ” arXiv:1908.027 88 , Aug. 2019. [Online]. A va ilable : https://git hub . com/po we r- grid- lib/pglib- opf [22] M. Abadi et al. , “T ensorFlo w: Large-scal e machine learning on heterog eneous systems, ” 2015, softw are av ailable from tensorflo w .org. [Online]. A va ilable : https://www .tensorflo w .org/ [23] R. D. Zimmerman, C. E. Murillo-S ´ anchez, and R. J. T homas, “MA TPOWER: Ste ady-state operations, planning, and analysis tools for po wer systems research and educ ation, ” IEEE T ransactions on P ower Systems , vol. 26, no. 1, pp. 12–19, Feb. 2010. [24] A. W ¨ achter and L. T . Biegler , “On the implementat ion of an interi or- point filter li ne-search algorithm for lar ge-scale nonlinear programming, ” Mathemat ical Pr ogrammi ng , vol. 106, no. 1, pp. 25–57, 2006. [25] MOSEK ApS, MOSEK 8.0.0.37 , 2018. [26] I. Z. Emiris and V . Fisik opoulos, “Prac tical polytop e volume approxi- mation, ” ACM T ransactions on Mathematic al Softwar e (TOMS) , vol. 44, no. 4, p. 38, 2018. [27] D. L ee, H. D. Nguyen, K. Dvijotha m, and K. Tur itsyn, “Con vex restrict ion of power flow feasib ility sets, ” IEE E T ransactions on Contr ol of Network Systems , vol. 6, no. 3, pp. 1235–1245, Sept. 2019. [28] B. Cui and X. A. Sun, “Solva bility of po wer flow equatio ns through exi stence and uniquen ess of complex fixed point, ” , April 2019.