An Effective Multi-Solution Approach for Power System Islanding

Controlled islanding is known as a last resort to prevent power system blackouts. The most important challenge in this context is the appropriate selection of separation points in a very short time considering the practical constraints. In this regard, several islanding methods have been proposed until now. Among them, the methods that provide the appropriate solutions are usually very complicated and difficult to implement in the real-time application of power system separation. Also, all methods result in one islanding solution, which may not be optimal due to the limitation of the commonly used objective function. To overcome these limitations and reach a comprehensive solution, this paper proposes a straightforward multi-solution approach through a suggested hierarchical spectral clustering algorithm. In this concept, the most desirable islanding scenario could be selected based on secondary criteria to reach more sustainable islands. In the proposed method, the hierarchical clustering algorithm, which has good records in other applications, is improved such that the generator coherency constraint can be considered in the clustering process. Meanwhile, the transmission lines without remote-controllable circuit breakers could be easily excluded from the islanding solutions. The proposed method is tested using the model of the IEEE 39-bus test system. Furthermore, to evaluate the computational effectiveness and accuracy of the method in a large-scale grid, the model of Khorasan Regional Electric Company (KREC) power system (which is the biggest part of Iran power system) is used. The comparative analysis with the state-of-the-art methods verifies the superiority of the proposed approach.

The associate editor coordinating the review of this manuscript and approving it for publication was Chi-Seng Lam .

φ(S)
Quality of an island. µ(S) Quality index of separation to n islands.
The general objective function for power system islanding. Vertices of the graph.

Vol(S)
The volume of the subgraph. w ij Weight of the edge between the two vertices v i and v j . W Weights of the edges in the graph. At present, electrical power systems are more prone to instability. Recent blackout reports around the world indicate that the occurrence of a short circuit at a transmission line could sometimes result in cascading outages, which in turn may lead to power system blackout. Due to time constraints, it is not typically easy to carry out corrective actions manually in case of successive failure occurrence. Hence, to inhibit the blackout and make easier system restoration, controlled islanding of the power grids through a specific protective system is highly recommended [1], [2]. To ensure sustainable islands, several restrictions must be met, including load-generation balance and coherency of the generators in each island, voltage stability and angle stability of the islands, remote controllability of the circuit breakers, and thermal limitation of the electrical equipment. Therefore, it is hard and cumbersome to find an appropriate approach to be applied in the real-time applications of the power system islanding for meeting all the aforementioned restrictions. Accordingly, a subgroup of these limitations, such as the load-generation balance and coherency of the generators can be utilized as a base for separation [3], [4], [5]. Other constraints can also be met using secondary corrective actions, e.g., load/generation shedding in the subsequent step [5]- [7].
The graph partitioning theory is principally the basis of recently introduced approaches for power system islanding.
These methods are categorized into two general groups based on their objective functions. The two objective functions in this regard are: 1) minimal power imbalance in the islands and 2) minimal power flow disruption. In [8], the power system was separated based on ordered binary decision diagram (OBDD) to satisfy the first aforementioned objective function. In [9], the researchers developed the previous algorithm to decrease the computational complexity. In [10], the grid was separated into three stages: 1) defining a domain for each generator, 2) determining the primary separation points assuming the grouping of coherent generators, and 3) adjusting the splitting points to achieve the minimum power imbalance in the islands. In [11] a new method is proposed based on PMU data to minimize the overall load shedding. Meanwhile, the generator coherency constraint is also considered in this method. In summary, the methods based on the first objective function involved some limitations: 1) the concept of optimization is based on searching techniques that lead to a high computational burden. 2) problem simplification which is sometimes proposed to reduce the running time prevents obtaining a globally optimal solution [5].
Recent researches have mostly focused on the strategy with minimum power flow disruption. For example, a two-stage islanding technique was introduced according to spectral clustering in [5]. In the aforementioned research, the first stage includes determining the coherent groups of generators through normalized spectral clustering, which is based on the dynamic model of the generators. The second stage identifies the separation lines through constrained spectral clustering to minimize power flow disruption and meet the restrictions of generator coherency. This technique is applicable just for two islands. This restriction was overcome by a recursive twophase approach in [5]. The disadvantage of this method is its intensive computational complexity [12], [13]. A novel approach was, therefore, offered based on constrained spectral clustering in [12], to reduce the computational complexity and time. It determines separation points to yield the minimum disruption of the power flow. The main disadvantage of the method is the unchangeable number of islands equaling to those of coherent groups of generators. In [14], the application of the hierarchical spectral clustering algorithm in the power system islanding context is introduced. The aforementioned paper focused on the mathematical basis and the power system islanding constraints, e.g., generator coherency was not considered in it. Accordingly, [14] could be considered as a starting point for future works in this concept. The k-medoids method in spectral clustering applied in [15] offered just a separation scenario similar to [12]. In [16] a new controlled islanding model is presented for the hybrid AC/VSC-HVDC grid. To obtain more sustainable islands, [17] employed a multi-layer clustering approach to minimize the reactive power disruption as well as the active power. The aforementioned algorithm uses the k-mean method for clustering which leads to only one solution. Since the reactive power control is usually a local issue, complicating the algorithm in [17] does not result in VOLUME 8, 2020 significant improvement in the separation outcomes. This method is improved in [18] by adding the correction coefficients between bus frequency components in the algorithm.
The approaches that apply the minimum power imbalance seek those separation lines which result in the least power imbalance in the islands. Such approaches typically utilize the searching strategies and consequently require a long time, which makes them unsuitable to be practically applied in power system islanding. The second group of approaches tries to discover separation lines to minimize the variation of power flowing through the transmission lines [5]. The above techniques enhance transient stability in the islands, decline the possibility of overload in the transmission lines, and facilitate grid restoration [7]. Besides, it is easier to implement this objective function through graph partitioning theory [5].
On the other hand, as the optimization procedure does not directly consider the minimum power imbalance concern, it may be necessary to do some complementary operations (for example, load/generation shedding), following separation to obtain stable islands.
Considering the advantages of the second objective function, this paper introduces an innovative algorithm based on the minimum active power flow disruption approach. Meanwhile, implementing the proposed algorithm results in several candidate solutions. So, to satisfy the operational requirements, e.g., the minimum power imbalance on each island, the most desirable option can be selected depending on a preferred secondary objective function. The basis of the proposed method is the hierarchical spectral clustering algorithm. Indeed, the conventional hierarchical spectral clustering method is improved such that the generator coherency constraint also can be added to the clustering process. Also, the transmission lines without remote-controllable circuit breakers could be easily excluded from the islanding solutions.
The main contribution of this paper could be summarized as follows: • Proposing a straightforward islanding method that is easy to implement: In contrast to the state-of-the-art methods ( [17], [18]), which use active and reactive power and bus frequency, the proposed algorithm is based only on the active power flow in transmission lines and consequently, it is much simpler.
• Introducing the multi-solution approach in the power system islanding context: In the proposed method, several islanding scenarios are obtained simultaneously through the proposed hierarchical clustering algorithm. The most desirable option could be chosen by defining secondary criteria to reach more sustainable islands. Any desired concern related to the post-islanding operation of the power grid could be included in the secondary objective function.
• Improving the hierarchical clustering algorithm: Hierarchical clustering is known as an effective method in data science [19], which is easy to implement. However, it is not possible to include the restriction on the aggregation of vertices in the optimization process.
To make the hierarchical clustering method suitable for power system islanding, new similarity criteria are proposed in the aggregation tree such that the generator coherency constraints could be satisfied in the clustering to reach stable islands. The effectiveness of the proposed method is tested using the model of the IEEE 39-bus test system. Moreover, to evaluate the computational efficiency and the accuracy of the method in a large-scale grid, the model of Khorasan Regional Electric Company (KREC) power system as the most significant part of the Iran power grid is used. The comparative analysis with the state-of-the-art methods confirms the superiority of the proposed approach.
The rest of this paper is organized as follows. Section II describes the graph partitioning theory. The graph spectral clustering technique is discussed in Section III. In Section IV, the objective function of the power system islanding is introduced. The suggested algorithm is presented in Section V. Then its performance is evaluated in Section VI using the test power systems. In Section VII the performance of the proposed algorithm is compared with the existing methods. Finally, the conclusion is provided in Section VIII.

II. GRAPH PARTITIONING THEORYINES
Consider a graph G = (V , E, W ), in which V , E, and W are vertices, edges, and edge weights, respectively. The goal of graph partitioning is to discover a set of vertices (i.e., a cluster) strongly connected and weakly connected with vertices of other sets. To evaluate the quality of an island, two quantities are defined in graph theory, i.e., the volume and the boundary of the subgraph. The boundary is the sum of edge weights between vertices inside S (i.e., the cluster) and vertices outside S [20]: where w ij is the weight of an edge between i-th and j-th vertices. Meanwhile, subgraph S is a group of vertices of the original graph, and i ∈ S shows a single vertex inside the subgraph. The volume of a subgraph is determined to be the sum of weighted degrees of its vertices [20]: where d i is the weighted degree of the i-th vertex (v i ) that needs to be calculated as follows [20]: where N is the number of vertices.
To measure the quality of an island, the following quantity is defined [12].
The smaller the value of φ(S) for a partition, the better the island. That is, an appropriate clustering method yields an island containing vertices strongly connected and weakly with other vertices outside the island. To partition a graph into n independent islands, the total islanding quality is determined by the highest quality of every single island: To discover an optimal solution, µ(S) should be calculated for all potential islanding scenarios, each one contains n independent islands. The scenario with the smallest µ(S) should be selected. Hence, the general objective function for the islanding is defined as As an NP-hard problem, the optimal solution cannot be computed for a big graph [21]. Therefore, it is recommended to employ Cheeger's inequality and spectral clustering to obtain a proper approximate solution [20].

III. GRAPH SPECTRAL CLUSTERING
Spectral clustering is among the efficient graph partitioning techniques, which rely on eigenvectors and eigenvalues of the Laplacian matrix of the graph [20].

A. GRAPH LAPLACIAN MATRIX
The Laplacian matrix L is widely used in graph analysis. For where N is the number of vertices [20], determined as follows: Also, the normalized Laplacian matrix is calculated as follows [22]: where D is a diagonal matrix with non-zero diagonal elements (d i ). The normalized Laplacian matrix can be computed as: which is scale-independent and preferred for clustering applications [13].

B. NORMALIZED LAPLACIAN MATRIX EIGENVALUES
Laplacian matrix eigenvalues are characterized by nonnegative real numbers, and the number of zero eigenvalues is indicative of the number of islands in the graph [23]. For normalized Laplacian matrices, the following inequality holds for eigenvalues: By employing k smallest eigenvalues together with their respective eigenvectors of the normalized Laplacian matrix, spectral clustering gives an approximate solution to solve the optimization problem in (6). This technique is preferred over direct analytical methods owing to its considerably lower computational complexity [24]. Using Cheeger's inequality, one can find out the extent to which the approximate solution is near-optimal [13]: Thus, choosing smaller λ k values result in a smaller ρ G (k) value, yielding better islanding quality.

C. SPECTRAL CLUSTERING THEORY
To determine the vertices in a k-dimensional Euclidean space of R k , spectral clustering employs k eigenvectors of the Laplacian matrix, referred to as spectral k-embedding. According to previous studies, more appropriate solutions are obtained by utilizing the normalized Laplacian in spectral clustering [13], [23]. Following spectral k-embedding, the vertices should be clustered using an appropriate algorithm in Euclidean space. The k-mean and the k-medoids methods, which are often used for this purpose, have some advantages, however, they have some limitations: 1) The number of clusters should be known, and 2) connections of vertices in the graph are not considered [20]. To prevail these limitations, hierarchical spectral clustering is utilized [24], [25], which is outlined in the following subsection.

D. HIERARCHICAL SPECTRAL CLUSTERING THEORY
Hierarchical clustering is well-known in data science and its outperformance is verified in several applications [19]. This approach involves creating a hierarchy of clusters [25]. The agglomerative hierarchical clustering technique (a bottom-up approach) could be employed for power grid islanding which is described as follows: • In a graph with N vertices, two vertices of the most similar (according to the similarity matrix) are chosen to build a cluster, thereby forming a novel graph with N −1 clusters. That is to say, a novel cluster is formed by merging two vertices, and the other clusters each contain one vertex.
• In the freshly-formed graph, a novel cluster is formed by merging the two most similar clusters. Under such circumstances, the newly formed graph consists of N −2 clusters.
• This procedure continues until novel higher-order clusters are formed. The hierarchical spectral clustering results are visualized in the form of a dendrogram. According to the dendrogram (Fig. 1), the similarity of the vertices is shown in a hierarchical form and the clusters could be easily obtained by choosing the number of clusters (n). The hierarchical spectral clustering approach has several advantages than the other VOLUME 8, 2020 clustering methods (e.g., k-mean and k-medoids): 1) changing the number of clusters without additional calculations, 2) providing an overall overview of the similarity between the clusters, and 3) easy implementation.

IV. POWER SYSTEM SEPARATION AND GRAPH THEORY A. GRAPH REPRESENTATION OF POWER SYSTEM
Power grids with N buses can be represented with an undirected weighted graph in the form of G = (V , E, W ). Here, V and E are vertices and edges, respectively, indicating buses and transmission lines in the electrical grid. Hence, This graph is simple without multiple loops and edges due to the nature of the electric power system. W represents the edge weights, i.e., the values of branch power flow. Let's assume that there are no network losses, then w ij = P ij = P ji , in which P ij represents active power flow between buses i and j.

B. THE OBJECTIVE FUNCTION OF POWER SYSTEM ISLANDING
In this paper, minimum power flow disruption is selected as the objective function due to the simple implementation and low computational complexity. In so doing, the largest islands with highly interconnected nodes (the largest power flow in corresponding branches) and loosely connected with other islands (with the smallest power flow in tie lines) should be determined. According to the graph partitioning theory, which is illustrated in Section III, in the graph representation of the power system, the boundary (∂(S)) represents the sum of power flowing through tie-lines connected to the island. Also, the volume (vol(S)) represents the internal power flow of the island plus boundaries. Therefore, the objective function for the power system islanding could be regarded as ρ G (k) which is defined in (6). It is worth noting that there are several practical constraints in power grid islanding that need to be satisfied, e.g., each island must have coherent generators. Furthermore, some of the circuit breakers are not remote-controllable, and thus the relevant transmission lines should be excluded from the final solution of islanding.

V. THE PROPOSED HIERARCHICAL CLUSTERING ALGORITHM
In this paper, to improve the post-islanding operating condition of the power grid, a multi-solution approach through a proposed hierarchical spectral clustering algorithm is introduced. In this concept, the most favorable islanding scenario could be selected based on a desired secondary objective function to reach more sustainable islands. In the proposed method, the hierarchical spectral clustering algorithm is improved such that the generator coherency constraint can be considered in the clustering process. The generator coherency grouping is assumed to be available through the existing PMU based methods. Meanwhile, the transmission lines without remote-controllable circuit breakers could be easily excluded from the islanding solutions.
The proposed approach is described as follows: 1) Identifying the necessity for intentional islanding based on power system vulnerability assessment (according to the existing methods [26]). 2) Employing online power system data: 2-1) Updating the adjacency matrix of the power system weighted graph (G) based on the new configuration of power grids and the results of power flow before the occurrence of disturbance. 2-2) Identifying coherent generators using techniques based on PMU data [27]- [29] and classifying them into k independent group. 3) Determining M lines to be excluded from separation and merging buses at both ends of such lines to build a novel N ×N adjacency matrix, in which N = N −M . 4) Calculating the normalized Laplacian matrix of the power grid graph (L N ) by the new adjacency matrix. 5) Determining eigenvectors corresponding to the smallest k eigenvalues, in which k is the number of coherent groups: 6) Calculating new coordinates of the vertices in the novel k-dimensional Euclidean space R k : 7) Normalizing coordinate vectors of the vertices: 8) Calculating the similarity matrix:   The elements of the similarity matrix for two neighboring vertices in the graph need to be calculated as follows: where u i , u j ∈ R k are the normalized coordinates of i-th and j-th vertices. Moreover, . implies the Euclidean distance. For two vertices that are not adjacent, the similarity is determined as the shortest path between the two in the respective graph. It is calculated according to the Floyd-Warshall algorithm [30]: 9) Determining the hierarchical structure and providing the dendrogram visualization as illustrated in Fig. 2  10) Choosing the desired numbers of clusters (n ≥ k) and identifying the lines that need to be separated.

VI. SIMULATION RESULTS
In this section, the proposed method is tested on the IEEE 39-bus test grid. Also, to evaluate the computational efficiency and the accuracy of the method in a large-scale grid, the KREC power system, which is the largest part of Iran power system, is used. All numerical calculations are performed by MATLAB software [31].

A. THE SIMULATION ON THE IEEE 39-BUS TEST GRID
To evaluate the proposed algorithm in a small grid, the model of the IEEE 39-bus test system is used [32]. This grid contains 10 generators and 46 transmission lines, as depicted in Fig. 3. It is supposed that a three-phase to ground short circuit occurred on the line between 16 and 17 buses and it was cleared after 150 ms by protective relays. Based on the simulation results in DIgSILENT software [33], oscillations happen severely in the grid, moving it towards the instability when an appropriate controlling mechanism is absent. Hence, it is necessary to separate the power system as a final solution for the prevention of wide-area instability.
Simulation results show that the generators oscillate in two coherent groups as follows: Following the number of coherent groups, embedding space dimension (k) is chosen as two. According to the proposed algorithm, the power grid graph and the associated Laplacian matrix (that is normalized) are updated. Then, the eigenvectors of this matrix are calculated, and the new geometrical coordinates of the vertices are determined. By calculating the similarity matrix and running the clustering algorithm, the dendrogram is obtained, as shown in Fig. 4. It is obvious that by choosing the number of required islands (n) on the dendrogram, the vertices of every island could be determined. Also, to achieve another islanding scenario, it is easy to change the number of islands, and no additional calculations are required. Table 1 shows the results of the suggested algorithm. According to the table, for two islands, two transmission lines (3-4 and 8-9) should be disconnected. Meanwhile, the values of φ(S i ); indicate the quality of the partitions; are 1.47% and 0.63% for islands one and two respectively. Hence, the total islanding quality (which is the maximum of φ(S))φ(S) is 1.47%. By comparing this value with the lowest possible value corresponding to the Cheeger inequality (0.6%), one can see the high quality of the simulation results. The power shortage on the first island (97 MW) is insignificant compared to total grid power (6100 MW) and is easily manageable using load shedding. It is worth noting that the numerical calculations are done within 1.2 milliseconds.

B. SIMULATION ON A REAL POWER GRID
To verify and show the effectiveness of the proposed method on a real-scale power system, the static model of Khorasan  Regional Electric Company (KREC) power system is used. This grid is located in Northeast of the Iran electric grid and includes 13 generators, 73 buses, and 109 transmission lines at 400 kV and 132 kV voltage levels and about 2500 MVA load. The Khorasan grid accounts for about 10 percent of the peak load of the Iran grid. However, this grid has a high VOLUME 8, 2020   potential for insatiability in some areas due to the following reasons: 1) The Khorasan grid connects to the Iran grid, only using a long 400 kV tie line (270 km long). 2) This grid is the most extensive grid of Iran. The high density of loads in Northeast of this grid has resulted in locating more than 90% of generation units of Khorasan grid in that area. Hence, long-length transmission lines are used for power transmission to other areas. According to our knowledge about KREC, following a short circuit in the 400 kV transmission line between the 5 th and 86 th buses connecting the south and north of the grid, and the consequent line outage by the protective relay, the generators begin to oscillate in two following groups: Based on the vulnerability analysis, it is necessary to separate the grid to avoid a blackout. To find the splitting lines, the proposed method is carried out and the results are tabulated in Table 2. According to the results, 246 MW (less than 10% of total load) load shedding is required for the bisection case. Also, the simulations take about 5.2 ms. Therefore, the proposed method can be used for real-time applications of power system islanding.

VII. COMPARATIVE ANALYSIS WITH THE STATE-OF-THE-ART METHODS
In this section, the performance of the proposed algorithm is compared with two state-of-the-art methods on the IEEE 39-bus model ( [17], [18]). The islanding methods in both papers are based on the spectral clustering approach (similar to the proposed method). In [17], the main objectives are to minimize the active and reactive powers (P, Q) flow disruptions. In [18], correlation coefficients between bus frequency (f) components are also considered in the objective function beside active and reactive powers (P, Q).
For a fair comparison, we used the same assumptions (initial fault lines and generator coherency grouping) for our algorithm and each competing method. As the computational burden of all spectral clustering methods is in the same order and they easily can be used in real-time applications, the performance of approaches is much important to be compared. According to the simulation results in Table 3, the methods of [17] lead to 853 MW and 619 MW total power shortage in the islands. On the other hand, this value is much less (199 MW) in the proposed method for three islands.
The comparison analysis with [18] is also summarized in Table 4. As shown, the total power shortage is the same for both methods (199 MW). However, it worth mentioning that the proposed method only uses the real power (P) while the method in [18] uses both real (P) and reactive power (Q) and also the frequency (f) as its input.
It should be noted that the two state-of-the-art methods lead to only one islanding solution (the number of islands is equal to the number of coherent groups) while in the proposed algorithm, the number of islands could be chosen equal or more than the number of coherent groups.

VIII. CONCLUSION
This paper introduces a novel power system islanding algorithm based on the minimum active power flow disruption approach. To overcome the limitation of the objective function and reach more sustainable islands, a multi-solution approach is introduced. In this concept, the most convenient islanding scenario could be selected based on a desired secondary objective function to satisfy the operational requirements. The proposed algorithm is based on the hierarchical clustering theory. The conventional hierarchical spectral clustering theory is improved to have the ability to accompany the generator coherency constraint in the clustering process. Also, the transmission lines without remote-controllable circuit breakers could be easily excluded from the islanding solutions. The effectiveness of the proposed method is tested using the model of the IEEE 39-bus test system. Moreover, to evaluate the computational efficiency and the accuracy of the method in a large-scale grid, the model of Khorasan Regional Electric Company (KREC) power system as the most significant part of the Iran power grid is used. The comparative analysis with the state-of-the-art methods shows the outperformance of the proposed approach.
MAHDI AMINI received the B.Sc. degree in elec-