Next Article in Journal
Visualization and Measurement of Swirling Flow of Dry Ice Particles in Cyclone Separator-Sublimator
Next Article in Special Issue
A Comprehensive Evaluation Method and Strengthening Measures for AC/DC Hybrid Power Grids
Previous Article in Journal
Microgrid Energy Management System Based on Fuzzy Logic and Monitoring Platform for Data Analysis
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Parametric Transient Stability Constrained Optimal Power Flow Solved by Polynomial Approximation Based on the Stochastic Collocation Method

1
The PowerChina Huadong Engineering Corporation Limited, Hangzhou 311122, China
2
The College of Electrical Engineering, Zhejiang University, Hangzhou 310027, China
3
The Huadong Branch of State Grid, Shanghai 200120, China
4
The Department of Electrical and Computer Engineering, University of Macau, Macau, China
*
Author to whom correspondence should be addressed.
Energies 2022, 15(11), 4127; https://doi.org/10.3390/en15114127
Submission received: 6 May 2022 / Revised: 26 May 2022 / Accepted: 2 June 2022 / Published: 3 June 2022
(This article belongs to the Special Issue Analysis and Control of Complex Power Systems)

Abstract

:
To better respond to the impact of power system-uncertain parameters on transient stability, a novel model named the parametric transient stability constrained optimal power flow (parametric TSCOPF) is proposed. It seeks the optimal control scheme of transient stability constrained optimal power flow (TSCOPF) expressed by the function of uncertain parameters in power systems. The key difficulty to solve this model lies in that the relationship between the parametric TSCOPF solution and uncertain parameters is implicit, which is hard to derive generally. To this end, this paper approximates the optimal solution of parametric TSCOPF by polynomial expressions of uncertain parameters based on the stochastic collocation method. First, the parametric TSCOPF model includes both uncertain parameters and transient stability constraints, in which the transient stability constraint is constructed as a set of polynomial expressions using the SCM. Then, to derive the relationship between the parametric TSCOPF solution and uncertain parameters, the SCM is applied to the parametric Karush–Kuhn–Tucker (KKT) conditions of the parametric TSCOPF model, so that the optimal solution of the parametric TSCOPF is approximated by using polynomial expressions with respect to uncertain parameters. The proposed parametric TSCOPF model has been tested on a 3-machine, 9-bus system, and the IEEE 145-bus system, which verifies the effectiveness of the proposed method.

1. Introduction

The transient stability constrained optimal power flow (TSCOPF) [1,2,3,4,5] is an efficient tool to reconcile the requirements of economic operation and dynamic security in power systems, and the model of TSCOPF problem has been generally addressed as a deterministic optimization problem. However, renewable resources and nodal loads are uncertain factors in power systems, which have significant influence on the solution of deterministic TSCOPF. So it is necessary to consider the effect of these uncertain parameters in the TSCOPF problem, and implement control schemes to respond to the fluctuations of uncertain factors in power systems.
According to the manners of parameter uncertainty quantification and optimization model construction, there are two commonly used methods of solving optimization problems considering uncertain parameters, i.e., probabilistic optimization [6,7,8] and robust optimization [9,10]. Probabilistic optimization assumes that the uncertain parameters obey a certain determined probability distribution. It minimizes the expectation of generation cost, etc. as the objective function to build a stochastic probabilistic model satisfying a certain probability level constraints. This approach has been used in problems such as optimal siting and sizing of distributed generators [11], security-constrained unit commitment [12,13,14,15] and reserve dispatch with uncertain wind power output [16,17,18]. It has been shown that the strategy derived from probabilistic programming can optimize the expectations of the objective function under uncertainty. However, probabilistic optimization modeling requires the knowledge of complete information about the probability distribution of uncertain parameters. Robust optimization is widely used because it does not need to know the information of probability distribution of uncertain parameters. Robust optimization limits the value of all uncertain parameters within their variation range, and transforms the problem into a deterministic optimization problem to cater to the needs of all extreme scenarios. Robust optimization has been widely used to deal with uncertain parameters in economic dispatch [19,20,21,22], security-constrained unit commitment [23], and cooperative energy and reserve scheduling for multi-microgrids [24]. However, although robust optimization can effectively deal with uncertain parameters, the optimization results are often very conservative due to its neglection of the distribution information of uncertain parameters.
In this paper, a novel approach named parametric optimization is proposed. In the model of the parametric optimization, uncertain parameters are treated as independent variables which vary freely in their range, and the system variables and control parameters are treated as dependent variables of uncertain parameters. So the parametric optimization seeks a solution that characterize the quantitive relationship between the parametric optimization solution and uncertain parameters.
The difficulty of solving the well-known TSCOPF is that the transient equality constraint is a set of differential algebraic equations (DAEs) which are hard to be dealt with in the optimization model. For the parameteric TSCOPF (Pa-TSCOPF), this difficulty will be even more serious because the uncertain parameters are free to change in the optimization model. Therefore, to obtain the parametric solution, the relationship between system variables (containing both state and algebraic variables) in the transient process and parameters (containing both uncertain and control parameters) needs to be reformulated, thus the DAEs can be eliminated from the optimization model. There are three kinds of ways to solve this problem: (1) The Monte Carlo method [25] is a straightforward method to analyze the statistical characteristic of the transient variables. However, the dynamic models of power systems have the characteristics of high dimensionality and strong nonlinearity; thus this method is time-consuming. (2) The analytical trajectory sensitivity method [26,27] is a linearization approximation method, which can be used to approximate the effect of parameters on the transient variables. Although this method is time-saving, it is a local linear approximation that is only accurate in the neighborhood of the expansion point. (3) The polynomial chaos expansion (PCE) method [28,29] uses a linear combination of a set of polynomial basis functions with undetermined coefficients to retain the nonlinearity of the transient process, so it is more accurate compared with the trajectory sensitivity method. This method falls into the category of polynomial approximation, which approximates the relationship between system variables in the transient process and parameters explicitly by polynomial approximation. Because this method has high accuracy and good convergence, it has been applied to many static or dynamic problems in engineering fields [30,31,32,33,34,35,36].
In terms of ways to solve undetermined coefficients, the PCE method can be classified into the stochastic Galerkin method (SGM) [37] and the stochastic collocation method (SCM) [38]. These two methods both can be applied to the parametric transient stability analysis and have high accuracy. The difference is that the SGM is an intrusive method and may become computationally intensive or even too hard to solve for systems with many parameters, whereas SCM is a non-intrusive PCE method that treats the system model as a black box and only needs to acquire samples of parameter-output function at the collocation points from the system model. Therefore, the SCM is more portable than the SGM as it can be easily integrated with existing commercial software such as BPA and PSS/E. Because of this, in this paper, the SCM method is applied to analyze the parametric transient stability by constructing the polynomial approximation of transient variables concerning parameters.
With transient variables been approximated by algebraic polynomial formulation, the DAEs transient constraint can be reformulated as a set of algebraic equations, and hence the differential algebraic Pa-TSCOPF can be reformulated in terms of algebraic parametric nonlinear programming (Pa-NLP). In efforts to solve the algebraic Pa-NLP, Qiu et al. in [39] propose an SGM-based method. However, in this paper, the Pa-TSCOPF model would be even more computationally demanding and too hard to solve because of the large size of the power system and the complexity of the model. As an alternative solution, an SCM-based method is developed to effectively solve this algebraic Pa-NLP by using PCE to approximate the relationship between the Pa-TSCOPF solution and uncertain parameters.
The main contributions of this paper are twofold: (1) the proposed SCM-based method is used to reformulate the original complicated set of stochastic DAEs in the optimization model into a simple set of algebraic equations in polynomial form, which can explicitly describe the relationship between system variables in the transient process and parameters. In this way, the original system of stochastic DAEs can be eliminated from the optimization model, which facilitates the solution of the optimization problem. (2) Compared with the robust TSCOPF and the probabilistic TSCOPF, the proposed Pa-TSCOPF can provide effective preventive control schemes according to parameter variations, and can avoid the conservativeness of the solution while improving the transient stability of the power system. In practice, the variation ranges of uncertain parameters are predicted in advance, then once the information about the uncertain parameters has been measured, the values of the uncertain parameters can be substituted into the parametric optimization solution. Because the computation time required to obtain the preventive control scheme through the parametric optimization solution is negligible, it can ensure the timeliness of the preventive control scheme.
The rest of this paper is organized as follows. The model of Pa-TSCOPF is first introduced in Section 2, Section 3 applies the SCM to the parametric transient stability analysis to reformulate the Pa-TSCOPF model as an algebraic Pa-NLP. Section 4 applies the SCM to solve the algebraic Pa-NLP, and Section 5 validates this method on the 3-machine, 9-bus system and IEEE 145-bus system. Section 6 concludes the paper.

2. Formulation of Parametric Transient Stability Constrained Optimal Power Flow

The proposed parametric TSCOPF (Pa-TSCOPF) is formulated as the following parametric nonlinear programming (Pa-NLP) in the form of DAEs, where the control parameters u are variables to be optimized, such as rescheduling the generators’ active powers or adjusting their terminal voltage, while the uncertain parameters p are treated as variables that vary freely in their range. In the parametric problem, the control variables u , the state variables x and the algebraic variables y depend on the value of uncertain parameters p ; thus they are functions in terms of p . In addition to that, the initial conditions x 0 and y 0 are also affected by the uncertain parameters, so initial conditions are also functions in terms of p . The mathematical models of the parametric problem are as follows:
min x 0 , y 0 , u c ( x 0 ( p ) , y 0 ( p ) , u ( p ) ; p )
subject to:
0 = g ( x 0 ( p ) , y 0 ( p ) , u ( p ) ; p )
0 = f ( x 0 ( p ) , y 0 ( p ) , u ( p ) ; p )
H 0 . mim H 0 ( x 0 ( p ) , y 0 ( p ) , u ( p ) ; p ) H 0 . max
0 = g ( x ( t ) , y ( t ) , u ( p ) ; p ) , t t 0 , t end
d x ( t ) d t = f ( x ( t ) , y ( t ) , u ( p ) ; p ) , t t 0 , t end
H t ( x ( t ) , y ( t ) , u ( p ) ; p ) 0 , t t 0 , t end .
The objective function c : R n x × R n y × R n u × R n p R in Equation (1) can be expressed as the fuel cost, network loss, or the renewable energy consumption capability, etc. g : R n x × R n y × R n u × R n p R n y in Equation (2) and (5) includes the models of network and machine–network interface. f : R n x × R n y × R n u × R n p R n x in Equations (3) and (6) includes dynamic models of generators and loads. Equations (2) and (3) describe the steady-state equality constraints including power flow balances. Equation (4) includes limits on bus voltage magnitude, generator active and reactive power outputs etc. Equations (5) and (6) describe the dynamics of the power system during the transient period t t 0 , t end , and Equation (7) represents the corresponding transient stability constraints such as the limit on relative swing angles in transient process.
Obviously, if p is set to a deterministic value, then Pa-TSCOPF Equations (1)–(7) becomes the traditional TSCOPF [1,2,3,4,5].
In the robust problem, the control parameters u will be influenced by the variation ranges of the uncertain parameters p but not by the particular values of p . But the state variables x , the algebraic variables y , the initial conditions x 0 and y 0 will depend on it. As we do not optimize the objective function for each p of the range independently, but jointly for all p , the objective function is reformulated as its supremum considering the variation ranges of p in the robust problem, such as
min x 0 , y 0 , u sup p D c ( x 0 ( p ) , y 0 ( p ) , u ; p ) s . t . H ( x 0 ( p ) , y 0 ( p ) , u ; p ) 0 , p D ,
where sup is the supremum calculation. The solution of the robust optimization is a set of deterministic optimal control schemes, which are independent of uncertain parameters.
For the probabilistic problem, the issue is similar to the one of the robust problem. The difference is that the control parameters will be influenced by the probability distribution of the uncertain parameters. In addition to that, the objective function is reformulated as its expectation considering the probability distribution of p in the probabilistic problem. Therefore, the mathematical models of the parametric problem are as follows:
min x 0 , y 0 , u E p F [ c ( x 0 ( p ) , y 0 ( p ) , u ; p ) ] s . t . P { H ( x 0 ( p ) , y 0 ( p ) , u ; p ) 0 } β , p F ,
where p F is the probability distribution of p . The solution of the probabilistic optimization is a set of deterministic optimal control schemes, which are independent of uncertain parameters either.
The aim of the above kinds of TSCOPF models is to search an optimal steady-state operating point x 0 , y 0 that is transient stable through regulate the control parameters u . The transient equality constraints Equations (5) and (6) is a set of parametric DAEs. Because this set of equations is difficult to be dealt with in the optimization model, the SCM is applied to reformulate Equations (5) and (6) as a set of parametric polynomial equations in the next section.

3. Reformulation of Transient Stability Constraints Based on the Stochastic Collocation Method

3.1. Polynomial Approximation of Transient Dynamics

As can be seen from the dynamic constraints of (5) and (6), the relationship between parameters p and system variables y is implicitly described by a set of DAEs with initial condition y 0 . To eliminate this set of DAEs from the optimization model to make the optimization model easy to solve, the SCM is introduced to approximate the relationship between y ( t ) and p explicitly with polynomial expression.
Because g / y 0 for t t 0 , t end in general, according to the implicit function theorem, the relationship between y and p can be approximated by
y ( t , p ) y ^ ( t , p ) = k = 1 N b c k y ( t ) Φ k ( p ) ,
where Φ k ( p ) φ k 1 ( p 1 ) φ k 2 ( p 2 ) φ k i ( p i ) φ k d ( p d ) is the k-th orthogonal polynomial basis function. As N is the total degree of Φ ( p ) , then N b = ( N + d ) ! N ! d ! . In this paper, the Legendre polynomial series are selected as basis functions; thus the range of p is transformed to [ 1 , 1 ] by linear transformations.
Then, the undetermined coefficients c k y ( t ) in (10) can be solved by the SCM with the following formula:
c k y ( t ) = 1 χ k m = 1 M α m Φ k ( p ( m ) ) ϕ y ( y 0 , t , p ( m ) ) ,
where χ k = D Φ k ( p ) Φ k ( p ) ω ( p ) d p with ω ( p ) = 1 for Legendre polynomial series; p ( m ) , m = 1 , , M are collocation points (also called integration points). The collocation points and corresponding integration coefficient α m are determined through the Smolyak sparse grid quadrature [28]; ϕ y ( y 0 , t , p ( m ) ) is the trajectory of y with fault information, initial condition y ( t 0 ) = y 0 and p = p ( m ) , which can be obtained through simulation by numerical methods. Therefore, the nonlinearity of the initial condition and the accuracy of numerical simulation results will affect the accuracy of the polynomial approximation. The common numerical integration includes Euler’s method, improved Euler’s method, implicit trapezoidal method, and Runge–Kutta method. Considering that the implicit trapezoidal method has better accuracy and numerical stability, in this paper, the implicit trapezoidal method is utilized. It can be observed in Equation (11) that c k y ( t ) is determined by initial condition y 0 when χ k , α m , p ( m ) are predesigned in SCM.
Note that polynomial approximation Equation (10) cannot reflect the statistical information of uncertain parameters. However, when their probability distribution is known, Equation (10) can be used to retrieve the statistical information of y . For instance, if the uncertain parameters are uniformly distributed, according to the orthogonality of Legendre polynomials, the mean and standard variance of y ( t , p ) can be approximated as [29]
μ y ( t ) = E [ y ( t , p ) ] ( k = 1 N b c k y ( t ) Φ k ( p ) ) d F ( p ) = c 1 y ( t ) ,
σ y = E [ ( y ( t , p ) μ y ( t ) ) 2 ] k = 2 N b [ χ k ( c k y ( t ) ) 2 ] .
Figure 1 is plotted to illustrate Equations (12) and (13). The fault condition of Figure 1 is presented in Section 6.1, and the uncertain parameter p L 5 is uniformly distributed in p L 5 [ 112.5 , 137.5 ] MW. By setting the total degree of polynomial bases to 2, the polynomial approximation of δ G 2 in terms of p L 5 is
δ G 2 δ ^ G 2 = k = 1 3 c k δ ( t ) Φ k ( p L 5 ) .
Coefficients c k δ ( t ) , k = 1 , 2 , 3 in (14) are plotted in Figure 1. As can be seen, the trajectories of c 1 y ( t ) and k = 2 3 χ k ( c k δ G 2 ( t ) ) 2 are coincident with the mean trajectory δ ¯ G 2 and standard variance trajectory σ respectively, which validates Equations (12) and (13).
If uncertain parameters are not uniformly distributed, statistical quantities of y ( t , p ) can also be readily approximated by combining the polynomial approximation Equation (10) with sampling methods such as the Monte Carlo method and Latin hypercube sampling [24]. Note that the computation of polynomial approximation does not require the uncertain parameters to be non-independent. However, when the obtained polynomial approximation is further used to analyze the statistical properties of system variables, it is necessary to consider the correlation of the uncertain parameters.

3.2. Reformulatin of Transient Stability Constraints

After the polynomial approximation of y has been obtained, the polynomial expression can be used to substitute variables y in transient stability constraints (7). Then transient stability constraints are discreted into a series of algebraic equations with numerical discretization, such that:
H t ( k = 1 N b c k p ( t ) Φ k ( y ) , p ) 0 ,
where H t = [ H t 1 , H t 2 , , H t s ] T and s is the number of total integration steps. According to (11) and (15) can be simplified as:
H t ( y 0 , p ) 0 .
Because all transient variables in transient stability constraints (7) have been substituted by polynomial expressions as shown in (10), the dynamic constraints (5) and (6), which define the transient outputs, are implicitly included in the reformulated transient stability constraints (16). Therefore, dynamic constraints (5) and (6) and the corresponding transient variable y can be eliminated from the model of Pa-TSCOPF, and the following algebraic Pa-NLP which eliminates differential equations is formulated:
min y 0 , u c ( y 0 , u , p )
subject to:
0 = G ( y 0 , u , p )
H 0 . mim H 0 ( y 0 , u , p ) H 0 . max
H t ( y 0 , u , p ) 0 .
The elimination of the dynamic constraints (5) and (6) from the parametric TSCOPF model gives the following advantages:
(1)
The effect of uncertain parameters p on the transient process can be explicitly evaluated by the transient stability constraint (20), which is a series of polynomials and easy to be evaluated in the Pa-TSCOPF problem.
(2)
The dynamic constraints (5) and (6) in the general Pa-TSCOPF model is a set of DAEs, which cannot be dealt with directly by existing optimization methods. Differently, the solution for the algebraic Pa-NLP in the reformulated Pa-TSCOPF model is easier to be dealt with.

4. Solution of Algebraic Pa-NLP Model Based on the Stochastic Collocation Method

The existence of uncertain parameters p in the proposed algebraic Pa-NLP model implies the optimal steady state operating point y 0 * ( p ) and the optimal control scheme u * ( p ) depend on uncertain parameters p . The challenge for solving the algebraic Pa-NLP is how to effectively evaluate the effect of uncertain parameters p on the optimal solution. The SCM-based optimization method proposed here is a parametric optimization method for algebraic Pa-NLP. It aims to approximate the relationship between the Pa-TSCOPF solution and uncertain parameters explicitly by polynomial expressions. This method starts with the formulation of parametric Karush–Kuhn–Tucker (KKT) conditions by the primal–dual interior point method (PDIPM) and gets the parametric optimal solution by applying the SCM to solve KKT conditions.

4.1. Parametric KKT Conditions of Algebraic Pa-NLP Model

Based on the classical PDIPM [40], by introducing slack variables for inequality constraints (19) and (20) and appending the logarithmic barrier functions to the objective, algebraic Pa-NLP model (17)–(20) can be reformulated as the following subproblem:
min y 0 , u c ( y 0 , u , p ) μ j = 1 r 0 log ( l ( p ) ) μ j = 1 r 0 log ( υ ( p ) ) μ j = 1 r t log ( υ t ( p ) )
subject to:
0 = G ( y 0 , u , p )
H 0 ( y 0 , u , p ) + υ ( p ) = H 0 . max
H 0 ( y 0 , u , p ) l ( p ) = H 0 . min
H t ( y 0 , u , p ) + υ t ( p ) = 0 ,
where l = [ l 1 ( p ) , , l r ( p ) ] , υ = [ υ 1 ( p ) , , υ r ( p ) ] , υ t = [ υ t 1 ( p ) , , υ t s ( p ) ] are the slack variables for inequality constraints respectively, and r, s are the numbers of steady-state inequality constraints (19) and transient stability inequality constraints (20). As μ tends towards zero, the solution of the reformulated subproblem approaches the solution of the primary optimization model (17)–(20).
The Lagrange function of the subproblem (21)–(25) is
L μ = c ( y ¯ , p ) ξ T ( p ) G ( y ¯ , p ) z T ( p ) [ H 0 ( y ¯ , p ) + υ ( p ) H 0 . max ] w T ( p ) [ H 0 ( y ¯ , p ) + l ( p ) H 0 . min ] γ T ( p ) [ H t ( y ¯ , p ) + υ t ( p ) ] μ j = 1 r 0 log ( l ( p ) ) μ j = 1 r 0 log ( υ ( p ) ) μ j = 1 r t log ( υ t ( p ) ) ,
where ξ ( p ) , z ( p ) , w ( p ) , γ ( p ) are Lagrange multipliers (dual variables) for constraints (22)–(25) respectively. The optimal solution of the subproblem (21)–(25) satisfies KKT first-order conditions:
y 0 L μ y 0 c ( y ¯ , p ) y 0 G ( y ¯ , p ) ξ y 0 H 0 ( y ¯ , p ) ( w + z ) y 0 H t ( y ¯ , p ) γ = 0
ξ L μ G ( y ¯ , p ) = 0
z L μ H 0 ( y ¯ , p ) l ( p ) H 0 . min = 0
w L μ H 0 ( y ¯ , p ) υ ( p ) H 0 . max = 0
γ L μ H t ( y ¯ , p ) + υ t ( p ) = 0
l L μ L Z e μ e = 0
υ L μ U W e + μ e = 0
υ t L μ U t R e = 0 ,
where e = [ 1 , , 1 ] T , L = diag ( l 1 , , l r ) , Z = diag ( z 1 , , z r ) , U = diag ( υ 1 , , υ r ) , W = diag ( w 1 , , w r ) , U t = diag ( υ t 1 , , υ t s ) , R = diag ( γ 1 , , γ s ) .
So far, the algebraic Pa-NLP model (17)–(20) is translated into a group of parametric KKT conditions (27)–(34). We will use the SCM introduced in Section 3.2 to find the polynomial approximations of the implicit functions mapping parameters p to the optimal solution.

4.2. Parametric Solution of Parametric TSCOPF Model

In Section 3.2, the way to apply the SCM to reformulate DAEs as algebraic equations has been presented; therefore, to avoid redundancy, here we will only pay attention to how to solve the parametric KKT conditions (27)–(34) using the SCM.
First, all variables in KKT conditions (27)–(34) can be represented as linear combinations of the polynomial basis functions Φ ( p ) as follows:
ϖ ( p ) ϖ ¯ ^ ( p ) = k = 1 N b c k ϖ Φ k ( p ) ,
where ϖ = [ y ¯ T , ξ T , l T , υ T , υ t T , z T , w T , γ T ] is the compact form of variables in KKT conditions.
Then, by the SCM introduced in Section 3.2 the unknown coefficients c k ϖ can be solved:
c k ϖ ( t ) = 1 χ k m = 1 M α m Φ k ( p ( m ) ) ϖ ( p ( m ) ) ,
where ϖ ( p ( m ) ) is the solution of parametric KKT conditions (27)–(34) with p = p ( m ) .
So far, the parametric optimal solution (35) of the initial parametric TSCOPF (1)–(7) has been obtained by the SCM. The parameterized optimal solution contains polynomial approximation for the optimal control parameters, system states, slack variables, and Lagrange multipliers, which explicitly describes the quantitative relationship between the Pa-TSCOPF solution and uncertain parameters.

5. Procedures of Parametric TSCOPF in Multi-Contingency Case

The preceding sections have introduced how to use the SCM to solve the Pa-TSCOPF considering a single contingency. In practical applications, when multiple contingencies are considered, the corresponding transient stability constraints need to be reformulated separately for each contingency according to the proposed method. The main procedures of the proposed method for solving the parametric TSCOPF problem in multi-contingency case are summarized below and depicted in Figure 2.
Step 1:
Input system data and construct the parametric TSCOPF model. Determine the k max contingencies to be considered in the transient stability analysis and set k = 1 .
Step 2:
For the k-th contingency, construct the polynomial approximation of transient variables as (10), and then the k-th transient stability constraint in parametric TSCOPF model can be reformulted as (20).
Step 3:
Add the k-th trainsient stability constraint into the reformulated model of the parametric TSCOPF, set k = k + 1 .
Step 4:
If k < k max , return to step 2. If k = k max , solve the reformulated model of parametric TSCOPF (17)–(20) by the method proposed in Section 4; the polynomial formulation (35) with coefficient solved by (36) is the parametric optimal solution of the reformulated parametric TSCOPF problem.

6. Case Studies

The 3-machine, 9-bus system case is used to verify the theoretical validity of the proposed method. In addition, because it is difficult to obtain the transient parameters of an actual power system, we use the IEEE 145-bus system to verify the feasibility of the proposed method in practical power systems. All computations are performed on Wolfram Mathematica 12.0 on a laptop with a CPU of Intel i5-6200U 2.30 GHz and a RAM of 4.0 GB, and the time step of numerical integration is 0.01 s .

6.1. 3-Machine 9-Bus System Case

6.1.1. Case Settings

The detailed data of the 3-machine, 9-bus system shown in Figure 3 can be found in [31]. The mechanical power of generators is assumed to remain constant during the transient process, and the dynamic models for transient stability analysis are 4th order models equipped with IEEE Type 1 exciter which is cited from [16]. The slack generator is at bus 1. The controllable parameters are power generation P G 2 [ 130 , 155 ] MW and P G 3 [ 45 , 155 ] MW.
Two fault scenarios are considered in this case study: fault 1 involves a three-phase fault on line 5 7 at the bus 7 end at t = 0.1 s , which is cleared after t cl = 0.07 s by tripping the faulted line 5 7 ; fault 2 is a three-phase fault on line 4 5 at the bus 4 end at t = 0.1 s , which is cleared after t cl = 0.2 s by tripping this faulted line.
There are 3 uncertain parameters allowed to vary, i.e., the real power P L 5 , P L 6 , P L 8 of the load at bus 5, bus 6, and bus 8, whose variation ranges are P L 5 [ 112.5 , 137.5 ] MW , P L 6 [ 81 , 99 ] MW , P L 8 [ 90 , 110 ] MW .
The objective function is defined as the traditional operating cost:
C G = i = 1 3 ( a i P G i 2 + b i P G i + c i ) ,
where P G i is active power of the i-th generator. [ a 1 , a 2 , a 3 ] T = [ 2 , 1 , 1 ] T , [ b 1 , b 2 , b 3 ] T = [ 3 , 1 , 1 ] T and [ c 1 , c 2 , c 3 ] T = [ 3 , 1 , 1 ] T are the fuel cost coefficients of generators.
The steady state constraint is that voltages of PQ buses are required to be in range [ 0.95 , 1.05 ] p.u. The transient stability constraint is defined as the following inequality constraints:
δ G i = δ G i k = 1 N G H k δ G k k = 1 N G H k δ G i . max , t ( t 0 , t end ] ,
where H k is inertia constant of the k-th generator, δ G i is the rotor angle of the i-th generator and δ G i is its relative rotor angle with respect to the center of inertia reference frame. δ G i . max is the upper limit of the rotor angle. It should be noted that other constraints, such as voltage dips or rises and line flow limits can also be included in transient stability constraints.

6.1.2. Approximation Error and Computation Time

To show the accuracy of the proposed method for the approximation of transient variables, fault 1 is taken as an example and the trajectory of δ G 2 is plotted in Figure 4 to compare the approximation error of the proposed stochastic collocation method (SCM) and the trajectory sensitivity method (TSM), where the results of the numerical simulation method are used as a benchmark. The nominal trajectory of δ G 2 , shown as a black line in Figure 4, is obtained with the parameter variation P L = [ P L 5 , P L 6 , P L 8 ] T = 0 by the numerical simulation method (NSM). Approximated trajectories were then formed for P L = ± 7.5 % P L 0 , where P L 0 is the expected value of P L . These approximated trajectories are shown as dashed lines for the SCM with N = 2 and dotted lines for the TSM respectively. For comparison, corresponding exact (fully simulated) trajectories are shown as solid lines. As we can see, the approximations obtained by the proposed SCM presented in this example are very accurate with parameter variation ± 7.5 % P L . Both the proposed method and the trajectory sensitivity-based approximation tracked the exact trajectory quite closely before t = 0.5 s . After that time, the approximation error of the trajectory sensitivities increased obviously. To quantitatively evaluate the approximation error, we use the mean absolute percentage error (MAPE) E y % as the error index:
E y % ( t ) 100 n i = 1 n y ( t ) y ^ ( t ) y ( t ) ,
where n is the number of simulation points.
Figure 5 plots the MAPE of rotor angle E δ G 2 % versus time t, and Table 1 shows the average MAPE E ¯ δ G % and computation time of the TSM and the SCM with different order N. It can be seen from Figure 5 that, no matter what N is, E δ G 2 % with the SCM during the whole time domain is smaller than that of the TSM. However, from Table 1 we can see that though the SCM’s average RMSE is much smaller than that of TSM, it needs more computation time in the respect that it requires a small number of simulations at the collocation points. Moreover, for the SCM, when N = 3 , although the average MAPE is reduced compared to that of N = 2 , more computation time is required. This is because the higher N, the more collocation points are needed for the SCM, which leads to more computation time. So it is suggested that N = 2 is sufficient for analysis of transient stability to retain both accuracy and efficiency in the general case.

6.1.3. Transient Stability Analysis of the Initial Operation Point

Because the polynomial approximation of transient variables concerning the uncertain parameters has been formed, their expectation and standard variance can be easily obtained by using Equations (10) and (11). When no preventive control is applied in the fault 1 and fault 2, the expectation trajectories of relative rotor angles δ G 1 , δ G 2 , δ G 3 with respect to a center of inertia reference frame and their standard deviation are plotted in Figure 6 and Figure 7, which are obtained through Equations (12) and (13).
As we can see, the trajectories of rotor angles will be affected by uncertain parameters in both faults. According to Chebyshev’s inequality, for a variable with any distribution, at least 93.75 % of the data were within 4 standard deviations of the expectation, or a little more conservative, at least 96 % were within 5 standard deviations of the expectation. Therefore, in fault 1, the possible maximum rotor angles is δ G 2 = 87.63 + 5 × 5.08 = 113.035 at t = 0.84 s . In fault 2, the possible maximum rotor angles is δ G 2 = 101.17 + 5 × 5.0 = 126.44 at t = 0.64 s . The maximum rotor angles in the critical stable scenario are δ G 2 = 114.3 in fault 1 and δ G 2 = 138.73 in fault 2 respectively. Therefore, affected by uncertain parameters, the initial operation point without preventive control is close to the critical stable scenario, and hence the implementation of parametric TSCOPF is needed. For the consideration of safety, we set the upper limit of the relative swing angle as 80 for fault 1 and 100 for fault 2 in (38).

6.1.4. Solution of the Parametric TSCOPF Problem

To show how the optimal solution varies with uncertain parameters, we make P L 8 be fixed to 100 MW , and hence objective function becomes a function of P L 5 and P L 6 , which forms a surface as shown in Figure 8. As we can see, the surface obtained by the proposed method is composed of 2 parts, part 1 in dark blue is transient stability constrained as δ G 2 encountering its upper limit in fault 2, whereas part 2 in light blue is transient stability constrained as δ G 2 encountering its upper limit in fault 1. When P L 5 or P L 6 decreases to the dashed line where part 1 and part 2 intersect, the encountering of the transient stability constraint is transformed from fault 1’s to fault 2’s. The RMSE in the whole parameter domain is 0.09 (about 0.83 % compared to the mean value of the objective function), showing good accuracy of the proposed approximation.
Furthermore, 500 pairs of the uncertain parameters are generated to compare the performance of the parametric control scheme obtained by the proposed Pa-TSCOPF with the control scheme obtained by the Pr-TSCOPF and R-TSCOPF. Note that for R-TSCOPF and Pa-TSCOPF with a 15 % probability of allowing constraints violation, δ G 2 80 in fault 1 and δ G 2 100 in fault 2 can not be satisfied. Therefore, the transient stability constraints are relaxed to δ G 2 90 and δ G 2 115 in fault 1 and fault 2 respectively.
Figure 9 plots the results in terms of the cost C G and the maximum rotor angle δ G 2 . max by different methods: (1) When the proposed Pa-TSCOPF is implemented, δ G 2 . max is controlled to be lower than 85 in fault 1 and 105 in fault 2. (2) Assume that the uncertain parameters are uniformly distributed over their range, when the Pr-TSCOPF with β = 0.9 . and R-TSCOPF are implemented, because there does not exist an optimal control scheme such that δ G 2 . max 85 in fault 1 and δ G 2 . max 105 in fault 2 are satisfied simultaneously, the transient constraints in Pr-TSCOPF and R-TSCOPF are relaxed to δ G 2 . max 90 in fault 1 and δ G 2 . max 110 in fault 2. As shown in Figure 9, it can be seen that in comparison with the Pa-TSCOPF and Pr-TSCOPF, although the cost of Pr-TSCOPF is generally lower, the transient stability constraint at many sampling points has been violated and the transient stability of the system under the influence of parameter uncertainty is not guaranteed; whiereas in comparison with the Pa-TSCOPF and R-TSCOPF, although the transient stability margin of the sampling results of R-TSCOPF is larger, the cost of the control scheme is generally higher than that of the Pa-TSCOPF.
The corresponding conclusions can also be obtained from the results in Table 2, which show that, compared with the other two methods, the control scheme provided by Pa-TSCOPF not only minimizes the out-of-limit probability of transient constraints but also makes the average cost of the control scheme relatively low. Therefore, the advantage of the proposed Pa-TSCOPF is to provide the most economical solution while maintaining the transient stability of the system. It should be noted that because the polynomial approximation still has a certain approximation error, the out-of-limit probability of transient constraints in fault 1 and fault 2 is 0.2 % and 0 % , respectively. That is to say, the proposed method may not be able to guarantee 100 % transient stability of the system. Besides, once the parametric solution of the Pa-TSCOPF is obtained, the optimal control scheme obtained by this parametric solution is indeed essentially equal to the deterministic TSCOPF. In practice, they both obtain the solution of TSCOPF based on the measurement results of uncertain parameters, and then apply the optimal control scheme to the power system. However, the computation time required to obtain the preventive control scheme through the solution of Pa-TSCOPF (around 0.0 s in this case) is much less than the computation time required to implement a deterministic TSCOPF (2.296 s in this case). Therefore, compared with the deterministic TSCOPF, the proposed Pa-TSCOPF can provide a preventive control scheme according to the real-time operating state of the power system in a more timely manner.

6.2. IEEE 145-Bus System Case

The detailed data of the IEEE 145-bus system can be found in [41]. In this case, three wind farms are added at bus 101, bus 105, and bus 112 to verify the effectiveness of the proposed method, and each wind farm consists of 100 DFIGs with rated power 1.5 MW, whose parameters are cited from [42]. Therefore, the variation range of the uncertain wind power generation P W 101 , P W 105 , P W 112 is [ 0 , 150 ] MW.
The fault scenario considered in this case is a three-phase fault occurring at one end of line 1–6 near bus 6 at t = 0.1 s , and subsequently cleared by tripping line 1–6 after a time period t cl = 0.12 s . Considering some generators have little influence on transient stability, to avoid adverse effects on the convergence and calculation efficiency of optimization due to too many parameters, the trajectory sensitivity method is used to screen the critical generators. Therefore, the control parameters in this case are the power generation by generators at bus 89, bus 98, bus 100, and bus 104, which are assumed to vary in the range of 80∼ 120 % of their nominal values.
By the SCM, the expectation rotor angle trajectories of several generators affected by wind speed and control parameters in the base case are plotted in Figure 10. Note that the expected trajectories and standard variance in Figure 10 are obtained by combining the polynomial approximation with the sampling-based Monte Carlo method. As we can see, the magnitude of δ G 104 ’s expectation trajectory is the highest among all generators, whose approximated possible maximum rotor angle is δ G 104 = 136.152 + 5 × 23.653 = 254.417 at t = 1.24 s . The maximum rotor angle in the critical stable scenario is δ G 104 = 182.935 . Therefore, affected by uncertain parameters and control parameters, the optimal operation point without considering transient stability constraints can be transiently unstable; thus the implementation of preventive control is required. The parametric TSCOPF model is the same as that of the first case. The upper limit of the relative rotor angle is set to 135 in this case.
The parametric TSCOPF model is solved by using the proposed approach, and the results are shown in Table 3. As we can see, the optimal control scheme will vary with the wind power generation P W 101 , P W 105 , P W 112 . In practice, this set of parametric control scheme can be utilized to adjust the control parameters according to the wind speed or schedule power generation in advance based on wind speed forecasts. In Figure 11, the performance of the parametric control scheme obtained by the proposed method and the control scheme obtained by the normal OPF is compared. The lower bound and upper bound trajectories of δ G 104 are formed by their vertices. As we can see, when the parametric control scheme is implemented, the maximum of δ G 104 is 136.2 . Almost all δ G 104 is controlled to be lower than 135 except 9 realizations, and all realizations are transient-stable under the proposed parametric control scheme; when the normal OPF which does not consider transient stability constraints is implemented, all realizations do not satisfy the transient stability constraint because the lower bound of δ G 104 is larger than 135 since t = 0.4 s and the power system is insecure when the fault occurs.

7. Conclusions

In this paper, a Pa-TSCOPF model is proposed for power system control to enhance the transient stability security, and the SCM-based polynomial approximation is presented to solve the Pa-TSCOPF. Firstly, to reformulate the transient stability constraints in the form of DAEs into algebraic equations, the SCM is utilized to approximate the transient variables by a set of polynomial expressions. Consequently the initial differential algebraic Pa-TSCOPF model is reformulated as an algebraic Pa-NLP model. Then, by applying the SCM to the parametric KKT conditions in an interior point method-like manner, the optimal control scheme is approximated by using polynomial expressions, which explicitly depict the relationship between the Pa-TSCOPF solution and uncertain parameters.
The advantage of the proposed Pa-TSCOPF model is that it can avoid the conservativeness of the control scheme and improves its economy compared with the R-TSCOPF, and can effectively deal with the inconsistency between the actual probability distribution and the predicted results of uncertain parameters compared with the Pa-TSCOPF. Note that based on the characteristics of different methods, suitable methods can be selected for different application scenarios to improve the transient stability of the power system. The proposed method can improve the transient stability of the power system in scenarios wherein the probability distribution of uncertain parameters varies greatly and the cost of preventive control schemes needs to be reduced.
Note that the disadvantages of the proposed method are mainly as follows: (1) since the resulting parametric control scheme needs to be computed based on the real-time values of the uncertain parameters, it leads to a more complicated acquisition of its control scheme than the deterministic control scheme provided by the probabilistic and robust methods; (2) the proposed method may suffer from the problem of large amounts of computational effort when the parameter dimension is high or the polynomial approximation order is high. Therefore, the screening of polynomial basis functions for specific computational scenarios to reduce the computational effort is a future research direction.

Author Contributions

Conceptualization, B.X. and H.W.; methodology, B.X.; software, B.X.; validation, W.Y., L.C. and Y.S.; formal analysis, B.X.; investigation, B.X.; resources, H.W.; data curation, H.W.; writing—original draft preparation, B.X.; writing—review and editing, H.W.; visualization, L.C.; supervision, Y.S.; project administration, W.Y.; funding acquisition, H.W. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by China NSF 51777184, Science and Technology Project SGTYHT/20-JS-221 from East Branch of State Grid Corporation of China.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Conflicts of Interest

The authors declare no conflict of interest.

Nomenclature

c { · } The objective function.
P { * } The probability calculation operator.
x State variables.
y Algebriac variables.
u Control parameters.
p Uncertain parameters.
n x , n y The dimension of x and y
n u , n p The dimension of u and p
f { · } The vector of functions corresponding to the differential part of power system model.
g { · } The vector of functions corresponding to the algebraic part of power system model.
x 0 , y 0 The initial values of x and y .
H 0 { · } The vector of functions corresponding to the steady-state inequality constraints.
H 0 . min The lower bound of H 0 .
H 0 . max The upper bound of H 0 .
H t { · } The vector of functions corresponding to the transient inequality constraints.
H The vector of functions corresponding to the inequality constraints.
DThe domain of p .
E [ * ] The expectation calculation operator.
β A fixed risk security level.
F ( · ) The distribution function of uncertain parameters.
p The compact form of u and p .
D The domain of p .
y The compact form of x and y .
y 0 The compact form of x 0 and y 0 .
c k y The vector of time-varied undetermined coefficients.
MThe number of collocation points.
Φ ( p ) The orthogonal multi-parameter polynomial basis function in terms of p .
φ ( p i ) The orthogonal single-parameter polynomial basis function in i-th dimension of p .
k , j The subscript of k-th or j-th multi-parameter polynomial basis function Φ ( p ) in i-th
dimension of p .
k i , j i The subscript of k-th or j-th single-parameter polynomial basis function φ ( p i ) in i-th
dimension of p .
NThe total degree of { Φ ( p ) } .
dThe dimension of p .
N b The number of polynomial basis functions.
α m The integration coefficient of the m-th collocation point.
ω The density of orthogonal Legendre polynomial series.
μ y The mean of y .
σ y The standard variance of y .
p Li The active load at bus i.
δ Gi The relative swing angle of generator i.
n s The number of total integration steps.
G The compact form of functions g and f .
μ The barrier parameter of interior point method.
y ¯ The compact form of y 0 and u .
ϖ All Lagrange variables in KKT conditions.

References

  1. Tang, L.; Sun, W. An Automated Transient Stability Constrained Optimal Power Flow Based on Trajectory Sensitivity Analysis. IEEE Trans. Power Syst. 2017, 32, 590–599. [Google Scholar] [CrossRef]
  2. Gan, D.; Thomas, R.J.; Zimmerman, R.D. Stability-constrained optimal power flow. IEEE Trans. Power Syst. 2000, 15, 535–540. [Google Scholar] [CrossRef] [Green Version]
  3. Xu, Y.; Dong, Z.Y.; Meng, K.; Zhao, J.H.; Wong, K.P. A Hybrid Method for Transient Stability-Constrained Optimal Power Flow Computation. IEEE Trans. Power Syst. 2012, 27, 1769–1777. [Google Scholar] [CrossRef]
  4. Pizano-Martínez, A.; Fuerte-Esquivel, C.R.; Zamora-Cárdenas, E.A.; Lozano-García, J.M. Directional Derivative-Based Transient Stability-Constrained Optimal Power Flow. IEEE Trans. Power Syst. 2017, 32, 3415–3426. [Google Scholar] [CrossRef]
  5. Yang, Y.; Qin, Z.; Liu, B.; Liu, H.; Hou, Y.; Wei, H. Parallel solution of transient stability constrained optimal power flow by exact optimality condition decomposition. IET Gener. Transm. Distrib. 2018, 12, 5858–5866. [Google Scholar] [CrossRef]
  6. Riaz, M.; Ahmad, S.; Hussain, I.; Naeem, M.; Mihet-Popa, L. Probabilistic Optimization Techniques in Smart Power System. Energies 2022, 15, 825. [Google Scholar] [CrossRef]
  7. Chen, S.; Hu, W.; Du, Y.; Wang, S.; Zhang, C.; Chen, Z. Three-stage relaxation-weightsum-correction based probabilistic reactive power optimization in the distribution network with multiple wind generators. Int. J. Electr. Power Energy Syst. 2022, 141, 108146. [Google Scholar] [CrossRef]
  8. Bian, J.; Wang, H.; Wang, L.; Li, G.; Wang, Z. Probabilistic Optimal Power Flow of an AC/DC System with a Multiport Current Flow Controller. CSEE J. Power Energy Syst. 2021, 7, 9. [Google Scholar]
  9. Huang, Z.; Zhang, Y.; Xie, S. Data-adaptive robust coordinated optimization of dynamic active and reactive power flow in active distribution networks. Renew. Energy 2022, 188, 164–183. [Google Scholar] [CrossRef]
  10. Khojasteh, M.; Faria, P.; Vale, Z. A Robust Model for Aggregated Bidding of Energy Storages and Wind Resources in the Joint Energy and Reserve Markets. Energy 2021, 238, 121735. [Google Scholar] [CrossRef]
  11. Liu, Z.; Wen, F.; Ledwich, G. Optimal Siting and Sizing of Distributed Generators in Distribution Systems Considering Uncertainties. IEEE Trans. Power Deliv. 2011, 26, 2541–2551. [Google Scholar] [CrossRef]
  12. Zhang, H.; Li, P. Chance Constrained Programming for Optimal Power Flow under Uncertainty. IEEE Trans. Power Syst. 2011, 26, 2417–2424. [Google Scholar] [CrossRef]
  13. Wang, J.; Shahidehpour, M.; Li, Z. Security-Constrained Unit Commitment with Volatile Wind Power Generation. IEEE Trans. Power Syst. 2008, 23, 1319–1327. [Google Scholar] [CrossRef]
  14. Wu, L.; Shahidehpour, M.; Li, T. Stochastic Security-Constrained Unit Commitment. IEEE Trans. Power Syst. 2007, 22, 800–811. [Google Scholar] [CrossRef]
  15. Wang, Q.; Guan, Y.; Wang, J. A Chance-Constrained Two-Stage Stochastic Program for Unit Commitment with Uncertain Wind Power Output. IEEE Trans. Power Syst. 2011, 27, 206–215. [Google Scholar] [CrossRef]
  16. Xia, S.; Luo, X.; Chan, K.W.; Zhou, M.; Li, G. Probabilistic Transient Stability Constrained Optimal Power Flow for Power Systems with Multiple Correlated Uncertain Wind Generations. IEEE Trans. Power Syst. 2016, 7, 1133–1144. [Google Scholar] [CrossRef]
  17. Papavasiliou, A.; Oren, S.S.; O’Neill, R.P. Reserve Requirements for Wind Power Integration: A Scenario-Based Stochastic Programming Framework. IEEE Trans. Power Syst. 2011, 26, 2197–2206. [Google Scholar] [CrossRef] [Green Version]
  18. Ahmadi-Khatir, A.; Conejo, A.J.; Cherkaoui, R. Multi-Area Energy and Reserve Dispatch under Wind Uncertainty and Equipment Failures. IEEE Trans. Power Syst. 2013, 28, 4373–4383. [Google Scholar] [CrossRef]
  19. Peng, C.; Xie, P.; Pan, L.; Yu, R. Flexible Robust Optimization Dispatch for Hybrid Wind/Photovoltaic/Hydro/Thermal Power System. IEEE Trans. Smart Grid 2016, 7, 751–762. [Google Scholar] [CrossRef]
  20. Lorca, A.; Sun, X.A. Adaptive Robust Optimization with Dynamic Uncertainty Sets for Multi-Period Economic Dispatch under Significant Wind. IEEE Trans. Power Syst. 2015, 30, 1702–1713. [Google Scholar] [CrossRef] [Green Version]
  21. Wu, W.; Chen, J.; Zhang, B.; Sun, H. A Robust Wind Power Optimization Method for Look-Ahead Power Dispatch. IEEE Trans. Sustain. Energy 2014, 5, 507–515. [Google Scholar] [CrossRef]
  22. Jabr, R.A. Adjustable Robust OPF with Renewable Energy Sources. IEEE Trans. Power Syst. 2013, 28, 4742–4751. [Google Scholar] [CrossRef]
  23. Bertsimas, D.; Litvinov, E.; Sun, X.A.; Zhao, J.; Zheng, T. Adaptive Robust Optimization for the Security Constrained Unit Commitment Problem. IEEE Trans. Power Syst. 2013, 28, 52–63. [Google Scholar] [CrossRef]
  24. Li, Y.; Zhao, T.; Wang, P.; Gooi, H.B.; Wu, L.; Liu, Y.; Ye, J. Optimal Operation of Multimicrogrids via Cooperative Energy and Reserve Scheduling. IEEE Trans. Ind. Inf. 2018, 14, 3459–3468. [Google Scholar] [CrossRef]
  25. Hajian, M.; Rosehart, W.D.; Zareipour, H. Probabilistic Power Flow by Monte Carlo Simulation with Latin Supercube Sampling. IEEE Trans. Power Syst. 2013, 28, 1550–1559. [Google Scholar] [CrossRef]
  26. Hiskens, I.A.; Pai, M.A. Trajectory sensitivity analysis of hybrid systems. IEEE Trans. Circuits Syst. I Regul. Pap. 2000, 47, 204–220. [Google Scholar] [CrossRef] [Green Version]
  27. Hiskens, I.A.; Alseddiqui, J. Sensitivity, Approximation, and Uncertainty in Power System Dynamic Simulation. IEEE Trans. Power Syst. 2006, 21, 1808–1820. [Google Scholar] [CrossRef] [Green Version]
  28. Shen, D.; Wu, H.; Xia, B.; Gan, D. Polynomial Chaos Expansion for Parametric Problems in Engineering Systems: A Review. IEEE Syst. J. 2020, 14, 4500–4514. [Google Scholar] [CrossRef]
  29. Xiu, D. Numerical Methods for Stochastic Computations: A Spectral Method Approach; Princeton University Press: Princeton, NJ, USA, 2010. [Google Scholar]
  30. Hockenberry, J.R.; Lesieutre, B.C. Evaluation of uncertainty in dynamic simulations of power system models: The probabilistic collocation method. IEEE Trans. Power Syst. 2004, 19, 1483–1491. [Google Scholar] [CrossRef] [Green Version]
  31. Xia, B.; Wu, H.; Qiu, Y.; Lou, B.; Song, Y. A Galerkin Method-Based Polynomial Approximation for Parametric Problems in Power System Transient Analysis. IEEE Trans. Power Syst. 2019, 34, 1620–1629. [Google Scholar] [CrossRef]
  32. Zhou, Y.; Wu, H.; Gu, C.; Song, Y. A Novel Method of Polynomial Approximation for Parametric Problems in Power Systems. IEEE Trans. Power Syst. 2017, 32, 32983307. [Google Scholar] [CrossRef] [Green Version]
  33. Qiu, Y.; Lin, J.; Liu, F.; Song, Y. Explicit MPC Based on the Galerkin Method for AGC Considering Volatile Generations. IEEE Trans. Power Syst. 2020, 35, 462–473. [Google Scholar] [CrossRef]
  34. Kaintura, A.; Dhaene, T.; Spina, D. Review of polynomial chaos based methods for uncertainty quantification in modern integrated circuits. Electronics 2018, 7, 30. [Google Scholar] [CrossRef] [Green Version]
  35. Najm, H.N. Uncertainty quantification and polynomial chaos techniques in computational fluid dynamics. Annu. Rev. Fluid Mech. 2009, 41, 35–52. [Google Scholar] [CrossRef]
  36. Kim, K.-K.K.; Shen, D.E.; Nagy, Z.K.; Braatz, R.D. Wiener’s polynomial chaos for the analysis and control of nonlinear dynamical systems with probabilistic uncertainties [historical perspectives]. IEEE Control Syst. Mag. 2013, 33, 58–67. [Google Scholar]
  37. Xiu, D.; Karniadakis, G.E. The Wiener–Askey polynomial chaos for stochastic differential equations. SIAM J. Sci. Comput. 2002, 24, 619–644. [Google Scholar] [CrossRef]
  38. Xiu, D.; Hesthaven, J.S. High-order collocation methods for differential equations with random inputs. SIAM J. Sci. Comput. 2005, 27, 1118–1139. [Google Scholar] [CrossRef]
  39. Qiu, Y.; Wu, H.; Song, Y.; Wang, J. Global Approximation of Static Voltage Stability Region Boundaries Considering Generator Reactive Power Limits. IEEE Trans. Power Syst. 2018, 33, 5682–5691. [Google Scholar] [CrossRef]
  40. Wu, Y.C.; Debs, A.S.; Marsten, R.E. A direct nonlinear predictor-corrector primal-dual interior point algorithm for optimal power flows. IEEE Trans. Power Syst. 1994, 9, 876–883. [Google Scholar]
  41. Power System Test Case Archive, University of Washington (2021, January) IEEE 145-Bus System. Available online: http://labs.ece.uw.edu/pstca/dyn50/pg_tcadd50.html (accessed on 5 May 2022).
  42. Ekanayake, J.B.; Holdsworth, L.; Wu, X.; Jenkins, N. Dynamic modeling of doubly fed induction generator wind turbines. IEEE Trans. Power Syst. 2003, 18, 803–809. [Google Scholar] [CrossRef] [Green Version]
Figure 1. Trajectotires of swing angle and its coefficients in 3-generator, 9-bus system case.
Figure 1. Trajectotires of swing angle and its coefficients in 3-generator, 9-bus system case.
Energies 15 04127 g001
Figure 2. Flowchart of SCM-based polynomial approximation for solving parametric TSCOPF model.
Figure 2. Flowchart of SCM-based polynomial approximation for solving parametric TSCOPF model.
Energies 15 04127 g002
Figure 3. Diagram of 9-buses system.
Figure 3. Diagram of 9-buses system.
Energies 15 04127 g003
Figure 4. Approximated trajectories of δ G 2 with P L = ± 7.5 % P L 0 by the SCM and the TSM in fault 1.
Figure 4. Approximated trajectories of δ G 2 with P L = ± 7.5 % P L 0 by the SCM and the TSM in fault 1.
Energies 15 04127 g004
Figure 5. Mean absolute percentage error of the SCM and the TSM in fault 1.
Figure 5. Mean absolute percentage error of the SCM and the TSM in fault 1.
Energies 15 04127 g005
Figure 6. Expectation and standard variance of rotor angles in fault 1.
Figure 6. Expectation and standard variance of rotor angles in fault 1.
Energies 15 04127 g006
Figure 7. Expectation and standard variance of rotor angles in fault 2.
Figure 7. Expectation and standard variance of rotor angles in fault 2.
Energies 15 04127 g007
Figure 8. Surface of objective function with P L 5 [ 112.5 , 137.5 ] MW, P L 6 [ 81 , 99 ] MW. Part 1 in dark blue: δ G 2 encounters its upper limit in fault 2. Part 2 in light blue: δ G 2 encounters its upper limit in fault 1.
Figure 8. Surface of objective function with P L 5 [ 112.5 , 137.5 ] MW, P L 6 [ 81 , 99 ] MW. Part 1 in dark blue: δ G 2 encounters its upper limit in fault 2. Part 2 in light blue: δ G 2 encounters its upper limit in fault 1.
Energies 15 04127 g008
Figure 9. Control objective versus maximum rotor angle δ G 2 . max of the parametric TSCOPF and normal OPF under 500 random realization of the uncertain parameters in the 3-machine, 9-bus system case.
Figure 9. Control objective versus maximum rotor angle δ G 2 . max of the parametric TSCOPF and normal OPF under 500 random realization of the uncertain parameters in the 3-machine, 9-bus system case.
Energies 15 04127 g009
Figure 10. Expectation and standard variance of rotor angles in the base case of IEEE 145-bus system.
Figure 10. Expectation and standard variance of rotor angles in the base case of IEEE 145-bus system.
Energies 15 04127 g010
Figure 11. The lower bound and upper bound trajecories of δ G 104 obtained by the SCM and TSM.
Figure 11. The lower bound and upper bound trajecories of δ G 104 obtained by the SCM and TSM.
Energies 15 04127 g011
Table 1. Computation time and approxiamtion error of the SCM and TSM in fault 1.
Table 1. Computation time and approxiamtion error of the SCM and TSM in fault 1.
Method E ¯ δ G 1 % 1 E ¯ δ G 2 % 1 E ¯ δ G 3 % 1 N p 2Total Time
TSM 38.29 % 6.91 % 13.11 % -8.548 s
SCM, N = 1 2.45 % 2.78 % 1.84 % 712.125 s
SCM, N = 2 0.88 % 0.98 % 0.73 % 2545.75 s
SCM, N = 3 0.27 % 0.30 % 0.24 % 69136.766 s
1 E ¯ δ G 1 % , E ¯ δ G 2 % , E ¯ δ G 3 % are mean values of E δ G 1 % , E δ G 2 % , E δ G 3 % . 2 N p is the number collocation points of the SCM.
Table 2. Comparison of optimization results of different methods in 3-machine, 9-bus system case.
Table 2. Comparison of optimization results of different methods in 3-machine, 9-bus system case.
MethodMean of δ G 2 . max Mean of CostOut-of-Limit
Probability
Case1Case2Case1Case2
Pr-TSCOPF 87.51 102.44 12.3511.8%10.4%
R-TSCOPF 86.37 97.66 13.112.2%0.4%
Pa-TSCOPF 84.00 102.04 12.340.2%0.0%
Table 3. Parametric Control Scheme of the IEEE 145-bus system case.
Table 3. Parametric Control Scheme of the IEEE 145-bus system case.
Control ParamterParametric Control Scheme
P G 89 202.356 0.0011872 P W 101 0.0066261 P W 105
0.001142 P W 112
P G 98 35016 778.754 P W 101 + 124.31 P W 101 2
4074.73 P W 105 + 1.34181 P W 101 P G 105
+ 126.188 P W 105 2 + 570.413 P W 112
2.06418 P W 102 P W 112 1.85842 P W 105 P W 112
96.1712 P W 112 2
P G 100 801.093 0.0074 P W 101 0.02623 P W 105
0.0045 P W 112
P G 104 507.08 0.003 P W 101 + 0.017 P W 105 0.003 P W 112
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Xia, B.; Wu, H.; Yang, W.; Cao, L.; Song, Y. Parametric Transient Stability Constrained Optimal Power Flow Solved by Polynomial Approximation Based on the Stochastic Collocation Method. Energies 2022, 15, 4127. https://doi.org/10.3390/en15114127

AMA Style

Xia B, Wu H, Yang W, Cao L, Song Y. Parametric Transient Stability Constrained Optimal Power Flow Solved by Polynomial Approximation Based on the Stochastic Collocation Method. Energies. 2022; 15(11):4127. https://doi.org/10.3390/en15114127

Chicago/Turabian Style

Xia, Bingqing, Hao Wu, Wenbin Yang, Lu Cao, and Yonghua Song. 2022. "Parametric Transient Stability Constrained Optimal Power Flow Solved by Polynomial Approximation Based on the Stochastic Collocation Method" Energies 15, no. 11: 4127. https://doi.org/10.3390/en15114127

APA Style

Xia, B., Wu, H., Yang, W., Cao, L., & Song, Y. (2022). Parametric Transient Stability Constrained Optimal Power Flow Solved by Polynomial Approximation Based on the Stochastic Collocation Method. Energies, 15(11), 4127. https://doi.org/10.3390/en15114127

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop