CN104091092A - Feature value analysis system for small-interference stability of large-scale power system - Google Patents

Abstract

The invention provides a feature value analysis system for small-interference stability of a large-scale power system, and belongs to the technical field of power system simulation and analysis. The feature value analysis system comprises a BPA data interface module for achieving seamless read-write operation of a BPA data file, a sparse matrix calculation module for achieving relevant sparse matrix processing in small-interference stability feature value analysis, a power flow calculation module, a system modeling and linearization module, a state matrix calculation module, an all feature value calculation engine module, a partial feature value calculation engine module and an oscillation mode extraction and analysis module. Omni-bearing and missing-free feature value simulation analysis for small-interference stability of an actual large-scale power grid is achieved by the adoption of TTQRE and IRAM and JDM with multiple calculation schemes.

Description

The Eigenvalues analysis system of large-scale electrical power system small signal stability

Technical field

What the present invention relates to is the system of a kind of electric system simulation and analysis technical field, specifically the Eigenvalues analysis system of the large-scale electrical power system small signal stability of three kinds of different characteristic value-based algorithms of a kind of employing.

Background technology

Power system small signal stability refers to when after electrical network experience small sample perturbations, continues to keep the ability of synchronous operation.Power system small signal stability generally adopts Liapunov first method as criterion.Liapunov first method points out, if the state matrix after system linearization does not occur zero or the eigenwert of positive real part, so just can judge that current system is small interference stability.Therefore, in electric system, eigenwert is calculated be all the assembling addressing of the identification that realizes low frequency power oscillation pattern, all kinds of stability controllers and parameter optimization, operational factor all the time to control parameter sensitivity analysis, detect online important prerequisite and the basic guarantee of the function such as modal information extraction of oscillation data.

Due to the continuous expansion of electrical network scale, from eighties of last century eighties, the researchist of power industry has just dropped into a large amount of energy in eigenwert computing method and analysis on Small Disturbance Stability systematic research and exploitation.Up to now, all Eigenvalues analysis systems that is used for large-scale electrical power system small signal stability in world wide, are there is, for example: U.S.'s Pacific Ocean gas and the SSAT of EISEMAN, the SSSP that U.S.'s DianKeYuan is combined Ontario, Canada hydroelectric board joint development of Utilities Electric Co.'s exploitation, Canadian power technology development in laboratory, NEVA, the PacDyn of Brazilian power science research centre exploitation, PSASP and the PSD ?SSAP etc. of China Electric Power Research Institute's exploitation that Siemens is developed.

Quantity according to eigenwert to be asked is sorted out, and the Eigenvalues analysis method of power system small signal stability is divided into All Eigenvalues analytic approach and partial feature value analytic approach.Although nearly all large-scale electrical power system analysis on Small Disturbance Stability system all comprises All Eigenvalues analysis and partial feature value is analyzed two aspects, but still there is following problem and shortage:

1) what the nucleus module of All Eigenvalues analytic approach still adopted is the dual step displacement implicit expression QR algorithm that the eighties of last century Kublanovskay sixties and Francis propose, thereby cause domestic power industry personage to generally believe: for the analysis on Small Disturbance Stability of large-scale power system, there is low memory in QR algorithm, computing time is very long, the eigenwert error calculating is very large, algorithm such as may not restrain at the problem [Wang Kang, Jin Yuqing, Gande is strong, Deng. the stability analysis of electric system small-signal and control summary [J]. Electric Power Automation Equipment, 2009 (5): 10 ?19. Xue Yu victory, Hao Sipeng, Liu Junyong. about the commentary [J] of low-frequency oscillation analysis method. Automation of Electric Systems, 2009, 33 (3): 1 ?8. China Electric Power Research Institutes, the little interference calculation user manual of PSASP7.0 version [R], Beijing: China Electric Power Research Institute, 2010. China Electric Power Research Institute, PSD ?SSAP analysis on Small Disturbance Stability program user handbook (2.5.2 version) [R], Beijing: China Electric Power Research Institute, 2012.].Undeniable, there is above-described variety of problems in early stage dual step displacement implicit expression QR algorithm really.But, along with the continuous progress of numerical computation method and computer hardware technique, QR algorithm has completed already from the multiple step displacement of the multiple step displacement-chain type of dual step displacement-bulk fritter-with positive differentiation [Francis J G F.The QR transformation a unitary analogue to the LR transformation-Part 1[J] the .The Computer Journal of the multiple step displacement of two step fritters of early-age shrinkage strategy, 1961, 4 (3): 265 ?271.Francis J G F.The QR transformation-part 2[J] .The Computer Journal, 1962, 4 (4): 332 ?345.Bai Z, Demmel J.On a Block Implementation of Hessenberg Multishift QR Iteration[J] .International Journal of High Speed Computing, 1989, 1 (1): 97 ?112.Braman K, Byers R, Mathias R.The multishift QR algorithm.part I:Maintaining well ?focused shifts and level 3 performance[J] .SIAM Journal on Matrix Analysis and Applications, 2002, 23 (4): 929 ?947.Braman K, Byers R, Mathias R.The multishift QR algorithm.Part II:Aggressive early deflation[J] .SIAM Journal on Matrix Analysis and Applications, 2002, 23 (4): 948 ?973.], personal computer is also already from only supporting 32 of 4GB addressing, become and can support 16EB (1EB=2 ³⁰gB) 64 of addressing.Therefore,, based on the All Eigenvalues analytic approach of dual step displacement implicit expression QR algorithm and 32 personal computers, obviously can not meet the computation requirement of current large-scale electrical power system analysis on Small Disturbance Stability;

2) partial feature value analytic approach is the main stream approach of present analysis large-scale electrical power system small signal stability.The nucleus module of partial feature value analytic approach is a class iterative projection method (Iterative projection methods) [Bai Z J, et al.Templates for the solution of algebraic eigenvalue problems:a practical guide[M] .Siam, 2000.].In the eigenwert algorithm research of the large-scale electrical power system analysis on Small Disturbance Stability over nearly 10 years, the iterative projection method that the frequency of occurrences is maximum has two, one is IRAM (the Implicitly Restarted Arnoldi Method under Krylov subspace, implicit restarted Arnoldi algorithm) [Kim D J, Moon Y H.Application of the implicit restarted Arnoldi method to the small ?signal stability of power systems[J] .Journal of Electrical Engineering & Technology, 2007, 2 (4): 428 ?433. second month in a season realize it, Song Xinli, Tang Yong, Deng. the electric system frequency domain character value parallel search algorithm [J] based on multi-process. Automation of Electric Systems, 2010 (21): 11 ?16.], another is Jacobi ?Davidson method (the JDM) [Du Zhengchun under non-Krylov subspace, Liu Wei, Fang Wanliang, Deng. in the analysis on Small Disturbance Stability based on Jacobi ?Davidson method, critical eigenvalue calculates [J]. Proceedings of the CSEE, 2005, 25 (14): 19 ?24.Tsai S H, Lee C Y, Wu Y K.Efficient calculation of critical eigenvalues in large power systems using the real variant of the Jacobi – Davidson QR method[J] .IET generation, transmission & distribution, 2010, 4 (4): 467 ?478.].

Although IRAM has extensively been integrated in each power system small signal stability analysis system, but the function that its algorithm possesses itself could not farthest be brought into play, for example: PSD ?SSAP and PSASP only provide based on Ping Yi ?the frequency sweep IRAM of inverse transformation, and be very suitable for calculating the Cayley conversion IRAM of crucial mode of oscillation, in two kinds of methods, all do not embody.Although frequency sweep IRAM has linear speed of convergence for utmost point eigenwert, its search behavior exists randomness, in the irrational situation of parameter configuration, easily occurs " leakage root ", misses some crucial especially mode of oscillation.In the time of Exact Solution update equation, JDM has the speed of convergence of progressive second order, but the more important thing is, even if do not use spectral transformation, JDM still can be according to the rule of certain appointment (real part maximum, damping ratio minimum), and order restrains those crucial mode of oscillation that to affect power system small signal stability.

But, there is the math library ARPACK increasing income unlike IRAM, because JDM only has the example procedure under Matlab, therefore, JDM is common in academic research always, and is not integrated in any large-scale electrical power system small signal stability Eigenvalues analysis system.

" three China " (Hua Zhong ?Hua Dong ?North China) extra-high voltage AC/DC electrical network interconnected, be bound to occur newly-increased low frequency power oscillation pattern, now " tetanic weak friendship " electrical network in first stage of construction, very likely make some pattern wherein will or present unsure state, therefore, adopt three kinds of diverse eigenwert algorithms of Numerical Principle to analyze the small signal stability of large-scale electrical power system, not only there is theory value, have more realistic meaning.

Summary of the invention

The present invention is directed to prior art above shortcomings, a kind of Eigenvalues analysis system of large-scale electrical power system small signal stability is proposed, according to PSD ?the electric network data file (* .dat and * .swi) of BPA form, adopt TTQRE (Two ?tone small ?bulge multishift QR algorithm with aggressive early deflation, with the multiple step displacement QR of the two step fritters algorithm of positive early-age shrinkage strategy), identify fast and accurately the whole electromechanic oscillation modes that affect power system dynamic stability; Employing has the IRAM of various computing schemes, the crucial electromechanic oscillation mode of part of extraction and analyzing influence power system dynamic stability; Employing has the JDM of various computing schemes, the crucial electromechanic oscillation mode of part of extraction and analyzing influence power system dynamic stability, the electromechanic oscillation mode obtaining by mating and compare three kinds of algorithms, final realization carries out to practical large-scale electrical network mode of oscillation identification and the model analysis that comprehensive, exhaustively little interfere with dynamic is stable.

The present invention is achieved by the following technical solutions:

The present invention relates to a kind of Eigenvalues analysis system of large-scale electrical power system small signal stability, comprise: BPA data interface module, be used for the sparse matrix computing module of the sparse matrix relevant treatment that realizes small signal stability Eigenvalues analysis, trend computing module, system modelling and linearization block thereof, state matrix computing module, All Eigenvalues computing engines module, partial feature value computing engines module and mode of oscillation are extracted and analysis module, wherein: BPA data interface module receives electric network data file, and be connected with system modelling and linearization block thereof with trend computing module respectively and API (Application Program Interface is provided, application programming interfaces), sparse matrix computing module is connected and API is provided with trend computing module, system modelling and linearization block thereof, state matrix computing module and partial feature value computing engines module respectively, the API that trend computing module provides according to BPA data interface module and sparse matrix computing module, realizes the trend of large-scale electrical power system and calculates and export calculation of tidal current to system modelling and linearization block thereof, the calculation of tidal current that the API that system modelling and linearization block thereof provide according to BPA data interface module and trend computing module provide, realize the linearization modeling of large-scale electrical power system, generate the system state matrix of augmentation and export respectively state matrix computing module and partial feature value computing engines module to, the system state matrix of the augmentation that the API that state matrix computing module provides according to sparse matrix computing module and system modelling and linearization block thereof provide, obtains system state matrix and exports All Eigenvalues computing engines module to, All Eigenvalues computing engines module obtains All Eigenvalues and the part left/right proper vector of system state matrix and exports mode of oscillation to and extracts and analysis module from system state matrix according to TTQRE, the API that partial feature value computing engines module provides according to sparse matrix computing module, with RCI (Reverse Communication Interface, contrary communication interface) mode carry out the numerical value iterative process of IRAM or JDM, the partial feature value and the part left/right proper vector that obtain system state matrix also export mode of oscillation to and extract and analysis module, mode of oscillation is extracted and analysis module is specified the eigenwert and the left/right proper vector that merge from the system state matrix of All Eigenvalues computing engines module or partial feature value computing engines module according to user, realize identification and the model analysis of electromechanic oscillation mode, and the pattern information of gained and modal analysis result are exported to user in the mode of Excel form.

The system state matrix J of described augmentation _augmeet: $J_{\mathrm{aug}}=\left[\begin{array}{cc}{J}_A& J_B\\ J_C& J_D\end{array}\right],$ Wherein: J _afor the block diagonal matrix that dynamic element inearized model is spliced, J _bfor the relational matrix of dynamic element and static cell, J _cfor the relational matrix of static cell and dynamic element, J _dfor the system admittance matrix of revising.

Described system state matrix S refers to:

In described All Eigenvalues computing engines module, comprise TTQRE unit, this unit is connected with state matrix computing module, and obtain according to system state matrix after the All Eigenvalues and part left/right proper vector of system state matrix, export mode of oscillation to and extract and analysis module.

Described partial feature value computing engines module comprises: IRAM unit and JDM unit, wherein: IRAM unit is connected with system modelling and linearization block thereof, according to user's selection determine adopt Ping Yi ?inverse transformation in conjunction with the maximum Ritz value selection strategy of mould or Cayley conversion the numerical procedure in conjunction with mould maximum Ritz value selection strategy, and obtain according to the system state matrix of augmentation after the partial feature value and part left/right proper vector of system state matrix, export mode of oscillation to and extract and analysis module; JDM unit is connected with system modelling and linearization block thereof, according to user's selection determine adopt Ping Yi ?inverse transformation in conjunction with the maximum Ritz value of mould selection strategy, Cayley conversion in conjunction with the maximum Ritz value of mould selection strategy, original matrix in conjunction with real part maximum Ritz value selection strategy or original matrix the numerical procedure in conjunction with damping ratio minimum Ritz value selection strategy, and obtain according to the system state matrix of augmentation after the partial feature value and part left/right proper vector of system state matrix, export mode of oscillation to and extract and analysis module.

Described mode of oscillation identification and model analysis refer to: the correlation factor that obtains individual features value according to every a pair of left/right proper vector, obtain the electromechanical circuit relevance ratio of this eigenwert according to correlation factor, if electromechanical circuit relevance ratio be greater than 1 and this eigenwert for plural number, judge that this eigenwert is as electromechanic oscillation mode, corresponding right proper vector is Oscillatory mode shape; The size of correlation factor mould value has been reacted the degree of coupling between Electrical Power System Dynamic element and mode of oscillation, on the unit of degree of coupling maximum, installing power system stabilizer, PSS additional, is the customary practice that suppresses the electric system low frequency power oscillation being brought out by this mode of oscillation; The amplitude of Oscillatory mode shape has determined the degree of waving of unit under this mode of oscillation, angle has determined the relative relation of waving between unit and unit under this mode of oscillation, and Oscillatory mode shape is to judge the influence degree of mode of oscillation to unit, judge whether multiple units can carry out the important indicator of dynamic equivalent.

Technique effect

Compared with prior art means, beneficial effect of the present invention comprises:

1. adopt the TTQRE that can represent the present art, realized the All Eigenvalues analysis of large-scale electrical power system small signal stability.This algorithm, in eigenwert iterative link, adopts chain type fritter to catch up with technology, has eliminated the displacement blooming that multiple step displacement QR piece catches up with, and has solved the problem that algorithm is not restrained; Adopt positive early-age shrinkage strategy, reduced the number of times of QR iteration, improved the computational accuracy of algorithm; Adopt Ju Zhen ?the senior computing of matrix, improved the utilization factor of cpu cache, saved the computing time of algorithm.From the angle of Practical, prove the feasibility that QR algorithm is analyzed for large-scale power system small signal stability All Eigenvalues;

2. the spectral transformation unification of IRAM is described as to operational form, realized have Ping Yi ?inverse transformation in conjunction with the maximum Ritz value of mould selection strategy, Cayley conversion the IRAM in conjunction with two kinds of numerical procedures of mould maximum Ritz value selection strategy.The great advantage of Cayley conversion IRAM is, can effectively identify some crucial mode of oscillation that damping ratio in electric system is less than or equal to certain designated value, compare Ping Yi ?inverse transformation IRAM, substantially there is not randomness in the search behavior of Cayley conversion IRAM, therefore, even if in the not exclusively rational situation of parameter configuration, Cayley conversion IRAM still can identify most crucial mode of oscillation, reduce to a great extent while analyzing the small signal stability of large-scale electrical power system with IRAM, occurred the possibility of " leakage root ";

3. realized JDM in the mode of RCI first at home and abroad.RCI is the perfect separation of having annotated between eigenwert algorithm and electric system application not only, and RCI also makes iterative projection method possess consistent calling rule.Compare the iterative projection method under Krylov subspace, JDM under non-Krylov subspace not only has speed of convergence faster, and in the situation that not using spectral transformation, JDM still can, according to the Ritz value selection strategy of certain appointment (damping ratio minimum, real part maximum), restrain the eigenwert that to want successively.Realize the JDM with 4 kinds of numerical procedures, and it has been integrated into together with IRAM in large-scale electrical power system small signal stability Eigenvalues analysis system, the practical application of JDM has been played to the effect of promotion.

Brief description of the drawings

Fig. 1 is technical scheme enforcement figure provided by the invention;

Fig. 2 a and Fig. 2 b are respectively left-half and the right half part of the numerical value iterative process schematic diagram of IRAM and JDM;

Fig. 3 be Ping Yi ?the effect figure of inverse transformation frequency sweep IRAM;

Fig. 4 is the effect figure of Cayley conversion in conjunction with the maximum Ritz value of mould selection strategy IRAM;

Fig. 5 is the effect figure of original matrix in conjunction with the minimum Ritz value of damping ratio selection strategy JDM;

Fig. 6 is the convergence process figure of original matrix in conjunction with the minimum Ritz value of damping ratio selection strategy JDM;

Fig. 7 is the model analysis figure of a certain inter-area oscillation mode in HD8241 system.

Embodiment

Below embodiments of the invention are elaborated, the present embodiment is implemented under taking technical solution of the present invention as prerequisite, provided detailed embodiment and concrete operating process, but protection scope of the present invention is not limited to following embodiment.

Embodiment 1

As shown in Figure 1, the present embodiment comprises: for realize BPA data file seamless read-write operation BPA data interface module, extract and analysis module for sparse matrix computing module, trend computing module, system modelling and linearization block thereof, state matrix computing module, All Eigenvalues computing engines module, partial feature value computing engines module and the mode of oscillation of the sparse matrix relevant treatment that realizes small signal stability Eigenvalues analysis, wherein:

Described BPA data interface module: because domestic power grid enterprises are partial to deposit actual electric network data by the file layout of BPA, for this reason, based on object based programming thought, adopt C++ to develop the BPA data interface module with dynamic link library form.This module not only has the seamless read-write capability of BPA data file, can also provide API for other modules.The embodiment of this module comprises following 3 key steps:

The Analysis of Hierarchy Structure of <1>BPA data file, is used for determining deposit position in BPA data file of static state/dynamic element parameter of electrical network and deposits order, realizes data card and the effective of control card separates;

The design of <2> class formation framework, according to the requirement of object based programming, designs a set of and BPA program user handbook and has the class formation framework of identical hierarchical relationship, the maintenance of Convenient interface module and expansion;

The generation of <3>BPA data-interface, adopts C++ to the realization of programming of designed class formation framework, and the data-interface of generation not only can provide API, and independently application program more can be provided.

Described sparse matrix computing module: the simulation analysis of any large-scale electrical power system and calculating, all be unable to do without sparse technology.The SuiteSparse that this module has adopted Timothy professor A.Davis of Florida State University to provide, realizes all computings relevant to sparse matrix, provides API taking the form of dynamic link library as other modules.This module is mainly made up of 5 parts:

The I/O of <1> data, realizes the bidirectional data transfers of sparse matrix between from hard disk to internal memory, and this part had both been supported the sparse matrix that 2 system forms are preserved, and met again the demand of 32/64 bit addressing simultaneously;

The computing of <2> fundamental matrix, realizes some conventionally calculations of sparse matrix, for example: transposition, rearrangement, Ju Zhen ?vector multiplication, Ju Zhen ?addition of matrices, Ju Zhen ?matrix multiplication etc.;

The sparse LU of <3> decomposes, and symbolic analysis and the LU numerical value of realizing sparse matrix decompose;

Sparse former generation/the back substitution of <4>, realizes the quick former generation/back substitution of sparse matrix;

The computing of <5> extended matrix, realizes the sparse matrix computing of some expansions, for example: the Xi that compute sparse matrix is mended, blocked about the Schur of certain sub-block dredge Ju Zhen ?vector multiplication etc.

In prior art, SuiteSparse is the sparse matrix solver equally celebrated for their achievements with SuperLU, PARDISO, TAUCS etc., the more important thing is, on the basis of SuiteSparse source code, can develop the sparse matrix computing function that is exclusively used in power system small signal stability Eigenvalues analysis, for example, the function shown in <4> and the <5> part of sparse matrix computing module.The practical application effect of these expansion sparse matrix computing functions, will provide in the back.

Described trend computing module: the API that utilizes BPA data interface module and sparse matrix computing module to provide, realizes the trend of large-scale electrical power system and calculates, for system modelling and linearization block thereof provide calculation of tidal current.

The calculation of tidal current that the API that described system modelling and linearization block thereof provide according to BPA data interface module and trend computing module provide, and [encourage just in conjunction with document, Chen Chen. application Design Mode exploitation analysis on Small Disturbance Stability method [J]. Proceedings of the CSEE, 2002,22 (1): 12 ?16.] the unified element connection modeling method that proposes, analysis on Small Disturbance Stability is carried out to system linearization modeling, generate the system state matrix J of augmentation simultaneously _augmeet: $J_{\mathrm{aug}}=\left[\begin{array}{cc}{J}_A& J_B\\ J_C& J_D\end{array}\right]$ Wherein: J _afor the block diagonal matrix that dynamic element inearized model is spliced, J _bfor the relational matrix of dynamic element and static cell, J _cfor the relational matrix of static cell and dynamic element, J _dfor the system admittance matrix of revising.

Described state matrix computing module: the sparse matrix J that utilizes API that sparse matrix computing module provides and system modelling and linearization block thereof to provide _a, J _b, J _cand J _d, realize the sparse of system state matrix S and solve: because S is dense matrix, and J _a, J _b, J _c, J _dbe sparse matrix, solve S according to the computing method of dense matrix, computing time can be very long on the one hand, and memory cost also can be very large on the other hand, so this module adopts following steps to solve S:

<1> is to J _dcarrying out symbolic analysis and sparse LU decomposes: P _dj _dq _d=L _du _d; Wherein: P _dand Q _dbe respectively the rearrangement matrix of realizing column selection pivot and the sequence of approximate minimum degree.

<2> J _acarry out the S:S=full (J of the dense storage of initialization _a)

The extended matrix calculation function that <3> utilizes sparse matrix computing module to provide, compute sparse matrix J _augabout sub-block J _dschur mend: due to the J in power system small signal stability analysis _band J _cthere is very special sparsity structure, and make the right-hand vector of <3> step formula not need to carry out whole numerical evaluation, therefore, sparse matrix computing module just can be on the source code basis of SuiteSparse, has realized the extended matrix calculation function of " compute sparse matrix is mended about the Schur of certain sub-block ".The practical application effect of this function will provide in the back.

Described All Eigenvalues computing engines module: adopt TTQRE, realize solving of All Eigenvalues to state matrix S, part left/right proper vector, comprise following 6 key steps:

<1> matrix equilibrium, adopts matrix to reset and isolates the part real character value being exposed on S diagonal line, adopts circulation convergent-divergent to reduce 2 norms of matrix A after balance: A=(PD) ^-1s (PD); Wherein: P is for resetting matrix, and D is diagonal matrix.

<2>Hessenberg reduction, adopts a series of Householder reflections, and A is approximately turned to upper Hessenberg form H:H=Q ^taQ; Wherein: Q is the orthogonal matrix that a series of Householder reflection accumulations form.

The multiple step displacement implicit expression QR iteration that <3> chain type fritter catches up with, adopting a series of Householder reflections and Givens rotation to carry out chain type QR piece catches up with, adopt contraction window to restrain multiple eigenwerts, by H similarity transformation, for intending upper triangular matrix T, the All Eigenvalues result of calculation of S is all diagonal elements of T: T=Q ^taQ; Wherein: Q is the orthogonal matrix that further accumulation forms, comprise a series of Householder reflections and Givens rotation.

Solving of <4> standard Schur type Partial Feature vector, determines after the eigenwert Λ wanting in T, calculates left eigenvector V and the right proper vector U of T about Λ: $\left\{\begin{array}{c}{V}^{H}T={\mathrm{\&Lambda;V}}^H\\ \mathrm{TU}=\mathrm{U\&Lambda;}\end{array}\right.;$

Q, V and the U that step <4> obtains that the P that the solving of <5> state matrix Partial Feature vector obtains according to step <1> and D, step <3> obtain, push back and calculate left eigenvector V and the right proper vector U of S about Λ: $\left\{\begin{array}{c}U=\mathrm{PDQU}\\ V={\mathrm{PD}}^{-1}\mathrm{QV}\end{array}\right.;$

<6>, by the All Eigenvalues calculating and part left/right proper vector, is input in text with the form of 2 systems.

The numerical evaluation storehouse that can complete above-mentioned calculation procedure in prior art is very many, for example: the reference BLAS/LAPACK of the MKL of Intel, the ACML of AMD, Netlib, ATLAS/LAPACK, the Matlab etc. of SourceForge.For making TTQRE algorithm can give play to best calculated performance on disposed computing machine, the OpenBLAS v0.2.8 that the present embodiment adopts Chinese Academy of Sciences software Suo Zhangxianyi team to provide, and the reference LAPACK v3.5.0 that provides of Netlib.The practical application effect of this module will provide in the above-described embodiments.

OpenBLAS is a set of optimization BLAS subroutine library that is directed to different CPU kernel framework.Its Shi Zhangxianyi team on the basis of the GotoBLAS2 of Hou rattan and luxuriant exploitation, exploitation and perfect multi-threaded parallel BLAS subroutine library gradually.For example, for different CPU micro-architectures: Penryn, Athlon, SandyBridge etc., the computing velocity of OpenBLAS is faster than MKL, ACML, ATLAS, Matlab etc., and matrix size is larger, and the acceleration effect of OpenBLAS is more obvious.

LAPACK is on the basis of BLAS subroutine library, adopts a set of linear algebra method that is exclusively used in dense matrix numerical evaluation of Fortran77 exploitation.LAPACK is early stage linear algebra method LINPACK and the aggregate of eigenvalue method EISPACK, contained that the Solving Linear in LINPACK, linear least-squares approach and EISPACK in the function such as eigenwert calculating.In addition, the LAPACK of 3.1.0 and later version, encapsulation be exactly TTQRE.

Described partial feature value computing engines module: in the mode of RCI function readjustment, the API that the numerical value iterative process of execution IRAM or JDM provides according to sparse matrix computing module, realizes all computings relevant to sparse matrix.This module comprises following 11 key steps:

<1> reads in the J of system modelling and linearization block generation thereof _a, J _b, J _cand J _d;

<2> splices J _augand carry out symbolic analysis, the rearrangement matrix Q that stet is analyzed _aug;

<3> is to J _dcarrying out symbolic analysis and LU numerical value decomposes;

<4> reads in the configuration information of partial feature value computing engines from * .cfg file, comprising: dimension, convergence precision, maximum iteration time, translation point σ and the anti-translation point μ of selected algorithm, spectral transformation type, Ritz value selection strategy, the eigenwert number of search, min/max search subspace;

<5> utilizes translation point σ to revise J _augpart diagonal element, and to revised J _aug– σ I carries out LU numerical value decomposition: P _aug(J _aug-σ I) Q _aug=L _augu _aug; Wherein: P _augand Q _augbe respectively the rearrangement matrix of realizing column selection pivot and the sequence of approximate minimum degree.

If the selected algorithm of <6> is IRAM, proceed to step <7>, otherwise selected algorithm is JDM, proceed to step <8>;

<7> calls the RCI function of IRAM, carries out the numerical value iteration of IRAM, proceeds to step <9>;

<8> calls the RCI function of JDM, carries out the numerical value iteration of JDM, proceeds to step <9>;

<9> judges whether numerical value iteration finishes, and if so, proceeds to step <11>, if not, proceeds to step <10>;

<10> according to the type of spectral transformation operator Op according to the intermediate variable of preserving in abovementioned steps, realize two class Ju Zhen in RCI function readjustment process ?sparse the solving of implicit expression of vector multiplication, after solving, get back to step <6>;

Two described class Ju Zhen ?vector multiplication refer to: wherein: the state matrix that S is system, the RCI readjustment vector that x provides for IRAM or JDM, the RCI readjustment scalar that θ provides for JDM, y be offer IRAM or JDM wait to ask vector.

If <11> algorithm convergence has obtained a part of eigenwert and right proper vector thereof, adopt twice inverse power method to solve left eigenvector corresponding to this partial feature value, and all result of calculation is input in text with the form of 2 systems.

As Fig. 2 a and Fig. 2 b (left-half and the right half part of the numerical value iterative process schematic diagram of IRAM and JDM) have provided the specific implementation process of above-mentioned steps <6>～step <10>, can be found out by Fig. 2 a and Fig. 2 b, the implementation of RCI function readjustment is very complicated, but teach at document [Dongarra J as Jack Dongarra, Eijkhout V, Kalhan A.Reverse communication interface for linear algebra templates for iterative methods[J] .UT, CS ?95 ?291, May, 1995.] described in introduction: " contrary communication is a technology, by it we can by various operations in iterative algorithm realize details give stash.On the one hand, algorithm development person can be without considering which type of data structure user adopts preserve the required matrix of iterative device, and on the other hand, the mode that user can want according to oneself, realizes the desired computing of iterative device ".Although the present embodiment only provides the embodiment of IRAM and two kinds of eigenwert algorithms of JDM, but can find out, for other iterative projection method, be deployed as RCI form as long as slightly make an amendment, just can directly add in partial feature value computing engines module, therefore, the present embodiment has very strong extensibility.

In described step <10>, two class Ju Zhen ?sparse the solving of implicit expression of vector multiplication, be that " Xi blocking dredge Ju Zhen ?vector multiplication " this extended matrix calculation function is realized by calling in sparse matrix computing module.

Be transformed to example with Cayley, i.e. Op (S)=(S – μ I) × (S – σ I) ^?1, IRAM and JDM all adopt following step realize the 1st class Ju Zhen ?sparse the solving of implicit expression of vector multiplication.

<a> comes initialization and J with x _augwith the column vector z=[x of dimension; 0];

<b> carries out sparse former generation/back substitution: z=Q to z _aug× (U _aug(L _aug(P _aug× z)))

The exponent number that <c> note n is S, the front n that takes out z is capable, generates column vector y=z (1:n);

<d> calls 2 grades of BLAS subroutine gemv, calculates y=x+ (σ – μ) × y.

With Ping Yi ?be inversely transformed into example, i.e. Op (S)=(S – σ I) ^?1, JDM adopt following step realize the 2nd class Ju Zhen ?sparse the solving of implicit expression of vector multiplication.

<e> calculates interim vector v: v=Q _d× (U _d(L _d(P _d× (J _c× x))))

<f> calculates interim vectorial u: $u=-\frac{(J_A-\mathrm{\&sigma;I})\&times;x+J_B\&times;v}{\mathrm{\&theta;}}$

<g> utilizes interim translation point α=(1+ θ × σ)/θ to revise J _augpart diagonal element, and to revised J _aug– α I carries out the decomposition of LU numerical value: ${\hat{P}}_{\mathrm{aug}}(J_{\mathrm{aug}}-\mathrm{\&alpha;I})Q_{\mathrm{aug}}={\hat{L}}_{\mathrm{aug}}{\hat{U}}_{\mathrm{aug}}$

<h> comes initialization and J with u _augwith the column vector z=[u of dimension; 0];

<i> carries out sparse former generation/back substitution to z: $Q_{\mathrm{aug}}\&times;({\hat{U}}_{\mathrm{aug}}\backslash ({\hat{L}}_{\mathrm{aug}}\backslash ({\hat{P}}_{\mathrm{aug}}\&times;z)))$

The exponent number that <j> note n is S, the front n that takes out z is capable, generates column vector y=z (1:n).

Described mode of oscillation is extracted and analysis module: the result of calculation that reads All Eigenvalues computing engines module or partial feature value computing engines module, realize identification and the model analysis function of electromechanic oscillation mode, calculate the most at last pattern and the modal information of gained, with the formal output of Excel form.

In order to further illustrate correctness and the validity of the present embodiment system, respectively the small signal stability of 3 ieee standard systems and 3 real systems is carried out to simulation analysis and calculating.Wherein, 3 real systems are taken from respectively the East China Power Grid high layout data of summer of 05,09 and 15 year.The details of these 6 test macros are as shown in table 1.

Table 1 test macro is described

Test macro	Node number	Branch road number	Generator number	State variable number	Algebraically variable number
						IEEE118	118	186	34	340	236
IEEE300	300	411	69	690	600
						IEEE600	600	823	139	1390	1200
HD4171	4171	5318	375	3750	8342
						HD5473	5473	7440	472	4720	10946
HD8241	8241	11657	500	3000	16482

Due to the data of test macro be all stored as PSD ?the file layout of BPA, therefore, by PSD ?SSAP v2.5.2 as the canonical reference of precision, canonical reference by Matlab R2012a as speed, correctness and the validity of the small signal stability Eigenvalues analysis system that checking the present embodiment provides.All tests are all carried out under same computing environment: the double-core CPU E5200 of Windows XP operating system, 4GB internal memory, Intel 2.5GHz.Because this CPU is Penryn micro-architecture, therefore, in the makefile.rule of OpenBLAS file, specify TARGET=PENRYN.Due to PSD ?SSAP be serial computing pattern, to even things up, the present embodiment and Matlab are also appointed as the serial computing pattern of core CPU.

Verification of correctness comprises following 2 aspects:

1) adopt respectively the TTQRE of the present embodiment and PSD ?the QR method of SSAP, calculate All Eigenvalues and the left/right proper vector of 6 test macros, after the electromechanic oscillation mode that identification is obtained mates between two, the relatively maximum absolute deviation of result of calculation, result is as shown in table 2;

2) adopt respectively the present embodiment Ping Yi ?inverse transformation IRAM and PSD ?SSAP Ping Yi ?inverse transformation IRAM, calculate partial feature value and the left/right proper vector of 6 test macros, after the electromechanic oscillation mode that identification is obtained mates between two, the relatively absolute deviation of result of calculation, result is as shown in table 3.Frequency sweep IRAM adopt PSD ?the parameter configuration of SSAP acquiescence.

The verification of correctness that table 2 All Eigenvalues is analyzed

Table 3 based on Ping Yi ?inverse transformation frequency sweep IRAM partial feature value analyze verification of correctness

Can be found out by table 2 and table 3, the present embodiment is correct.Produce mode of oscillation number and occur that the reason of deviation is,

The present embodiment and PSD ?SSAP adopted diverse ways to calculate correlation factor.

Validation verification comprises following 3 aspects:

1) adopt respectively the present embodiment state matrix computing module, PSD ?the sparse matrix computing function of SSAP, Matlab, show the state matrix that calculates 6 test macros, three's computing time is as shown in table 4;

2) adopt respectively the present embodiment TTQRE, PSD ?the QR method, the eig function of Matlab of SSAP, calculate All Eigenvalues and the left/right proper vector of 6 test macros, three's computing time is as shown in table 5;

3) adopt respectively the present embodiment Ping Yi ?inverse transformation IRAM and PSD ?SSAP Ping Yi ?inverse transformation IRAM, calculate partial feature value and the left/right proper vector of 6 test macros, both computing times are as shown in table 6.Similarly, frequency sweep IRAM adopt PSD ?the parameter configuration of SSAP acquiescence.

The validation verification that table 4 state matrix calculates

The validation verification that table 5 All Eigenvalues and proper vector are calculated

Test macro	The TTQRE/s of the present embodiment	PSD ?QR method/s of SSAP	Eig function/s of Matlab
				IEEE118	0.203	0.540	0.408
IEEE300	1.705	4.190	2.029
				IEEE600	13.144	31.115	14.371
HD4171	236.007	694.957	253.122
				HD5473	450.542	1812.102	482.496
HD8241	128.400	376.078	139.143

Table 6 based on Ping Yi ?inverse transformation frequency sweep IRAM partial feature value analyze validation verification

Test macro	The Ping Yi of the present embodiment ?inverse transformation frequency sweep IRAM/s	PSD ?SSAP Ping Yi ?inverse transformation frequency sweep IRAM/s
			IEEE118	1.058	1.206
IEEE300	1.523	1.618
			IEEE600	3.224	3.388
HD4171	63.072	76.291

HD5473	72.184	85.682
			HD8241	24.825	27.071

Table 4 has fully proved in the sparse matrix computing module of the present embodiment, the validity of " compute sparse matrix is mended about the Schur of certain sub-block " this extended matrix calculation function.Table 5 has fully proved in the All Eigenvalues computing engines module of the present embodiment, the validity of TTQRE algorithm.Table 6 fully proved in sparse matrix computing module, the validity of " Xi blocking dredge Ju Zhen ?vector multiplication " this extended matrix calculation function.

For investigate the Cayley conversion that provides in the present embodiment in conjunction with the IRAM unit of the maximum Ritz value of mould selection strategy, original matrix the JDM unit in conjunction with damping ratio minimum Ritz value selection strategy, for the practical application effect of the crucial mode of oscillation of identification large-scale electrical power system, the HD8241 system of getting in above-described embodiment is carried out emulation testing.Crucial mode of oscillation refers to damping ratio and is not more than 5% eigenwert.

Figure 3 shows that Ping Yi ?the electromechanic oscillation mode that identifies of inverse transformation frequency sweep IRAM.Frequency sweep IRA identifies 115 patterns, wherein 72 critical modes that are less than 5% for damping ratio, and all mode has 499, and wherein 87 is critical mode, therefore, there is " leakage root " in frequency sweep IRA method, particularly near 1.2Hz.

Figure 4 shows that the electromechanic oscillation mode that Cayley conversion identifies in conjunction with the maximum Ritz value of mould selection strategy IRAM.Cayley conversion IRAM has iteration altogether 70 times, and 60.058 seconds consuming time, identify the crucial electromechanic oscillation mode that whole 87 damping ratios are less than 5%, therefore, do not occur " leakage root ".

Figure 5 shows that the electromechanic oscillation mode that original matrix identifies in conjunction with the minimum Ritz value of damping ratio selection strategy JDM.The JDM of this numerical procedure has iteration altogether 109 times, and 208.007 seconds consuming time, identify the crucial electromechanic oscillation mode that whole 87 damping ratios are less than 5%, therefore, do not occur " leakage root ".

Figure 6 shows that the numerical value iteration convergence process of original matrix in conjunction with the minimum Ritz value of damping ratio selection strategy JDM.Although the counting yield of JDM is not as good as IRAM in the present embodiment, as can be seen from this figure, JDM can, according to the ascending order of damping ratio, restrain crucial mode of oscillation successively.

The electromechanic oscillation mode identifying due to native system and modal information are all with the formal output of Excel form, and therefore, the function that user both can adopt Excel to carry, can adopt again other mode to use result of calculation.

Still, taking HD8241 system as example, provide two kinds of application scenarioss of native system simulation analysis result:

1) electromechanic oscillation mode in Excel, original matrix being obtained in conjunction with the minimum Ritz value selection strategy JDM of damping ratio and translation ?inverse transformation frequency sweep IRAM according to the arrangement of damping ratio ascending order after, find that it is 0.04303542 that IRAM has missed a damping ratio, frequency is the crucial mode of oscillation of 1.19297526Hz, by the correlation factor of this mode of oscillation according to the descending sort of mould value after, find with the maximally related unit of this mode of oscillation to be " ZCN_1 " and " ZCN_2 ", so, on these two units, install respectively power system stabilizer, PSS additional, can effectively improve the dynamic stability of HD8241 system,

2) modal information that the JDM of reading and saving in Excel file calculates in Matlab, finding damping ratio is 0.02234904, frequency is the crucial mode of oscillation of 0.62502012Hz, because the oscillation frequency of this pattern is lower, can judge that it is inter-area oscillation mode, in order to find out the vibration region under this pattern, the modal information of this pattern is mapped to result as shown in Figure 7, as seen from Figure 7 under this pattern, between the unit in " vibration region 1 " and " vibration region 2 ", relatively wave, therefore, analysis result shows that these two all units in vibration region can carry out dynamic equivalent.

Claims (6)

1. the Eigenvalues analysis system of a large-scale electrical power system small signal stability, it is characterized in that, comprise: comprising: BPA data interface module, be used for the sparse matrix computing module of the sparse matrix relevant treatment that realizes small signal stability Eigenvalues analysis, trend computing module, system modelling and linearization block thereof, state matrix computing module, All Eigenvalues computing engines module, partial feature value computing engines module and mode of oscillation are extracted and analysis module, wherein: BPA data interface module receives electric network data file, and be connected and API is provided with system modelling and linearization block thereof with trend computing module respectively, sparse matrix computing module is connected and API is provided with trend computing module, system modelling and linearization block thereof, state matrix computing module and partial feature value computing engines module respectively, the API that trend computing module provides according to BPA data interface module and sparse matrix computing module, realizes the trend of large-scale electrical power system and calculates and export calculation of tidal current to system modelling and linearization block thereof, the calculation of tidal current that the API that system modelling and linearization block thereof provide according to BPA data interface module and trend computing module provide, realize the linearization modeling of large-scale electrical power system, generate the system state matrix of augmentation and export respectively state matrix computing module and partial feature value computing engines module to, the system state matrix of the augmentation that the API that state matrix computing module provides according to sparse matrix computing module and system modelling and linearization block thereof provide, obtains system state matrix and exports All Eigenvalues computing engines module to, All Eigenvalues computing engines module obtains All Eigenvalues and the part left/right proper vector of system state matrix and exports mode of oscillation to and extracts and analysis module from system state matrix according to TTQRE, the API that partial feature value computing engines module provides according to sparse matrix computing module, carry out the numerical value iterative process of IRAM or JDM in the mode of RCI, obtain partial feature value and the part left/right proper vector of system state matrix and export mode of oscillation to and extract and analysis module, mode of oscillation is extracted and analysis module is specified the eigenwert and the left/right proper vector that merge from the system state matrix of All Eigenvalues computing engines module or partial feature value computing engines module according to user, realize dynamo-electric mode of oscillation identification and model analysis, and the pattern information of gained and modal analysis result are exported to user in the mode of Excel form.

2. system according to claim 1, is characterized in that, the system state matrix J of described augmentation _augmeet: $J_{\mathrm{aug}}=\left[\begin{array}{cc}{J}_A& J_B\\ J_C& J_D\end{array}\right],$ Wherein: J _afor the block diagonal matrix that dynamic element inearized model is spliced, J _bfor the relational matrix of dynamic element and static cell, J _cfor the relational matrix of static cell and dynamic element, J _dfor the system admittance matrix of revising.

3. system according to claim 1, is characterized in that, described system state matrix S refers to:

4. system according to claim 1, it is characterized in that, in described All Eigenvalues computing engines module, comprise TTQRE unit, this unit is connected with state matrix computing module, and obtain according to system state matrix after the All Eigenvalues and part left/right proper vector of system state matrix, export mode of oscillation to and extract and analysis module.

5. system according to claim 1, it is characterized in that, described partial feature value computing engines module comprises: IRAM unit and JDM unit, wherein: IRAM unit is connected with system modelling and linearization block thereof, according to user's selection determine adopt Ping Yi ?inverse transformation in conjunction with the maximum Ritz value selection strategy of mould or Cayley conversion the numerical procedure in conjunction with mould maximum Ritz value selection strategy, and obtain according to the system state matrix of augmentation after the partial feature value and part left/right proper vector of system state matrix, export mode of oscillation to and extract and analysis module; JDM unit is connected with system modelling and linearization block thereof, according to user's selection determine adopt Ping Yi ?inverse transformation in conjunction with the maximum Ritz value of mould selection strategy, Cayley conversion in conjunction with the maximum Ritz value of mould selection strategy, original matrix in conjunction with real part maximum Ritz value selection strategy or original matrix the numerical procedure in conjunction with damping ratio minimum Ritz value selection strategy, and obtain according to the system state matrix of augmentation after the partial feature value and part left/right proper vector of system state matrix, export mode of oscillation to and extract and analysis module.

6. system according to claim 1, it is characterized in that, described mode of oscillation identification and model analysis refer to: the correlation factor that obtains individual features value according to every a pair of left/right proper vector, obtain the electromechanical circuit relevance ratio of this eigenwert according to correlation factor, if electromechanical circuit relevance ratio be greater than 1 and this eigenwert for plural number, judge that this eigenwert is as electromechanic oscillation mode, corresponding right proper vector is Oscillatory mode shape; The size of correlation factor mould value has been reacted the degree of coupling between Electrical Power System Dynamic element and mode of oscillation, on the unit of degree of coupling maximum, installing power system stabilizer, PSS additional, is the customary practice that suppresses the electric system low frequency power oscillation being brought out by this mode of oscillation; The amplitude of Oscillatory mode shape has determined the degree of waving of unit under this mode of oscillation, angle has determined the relative relation of waving between unit and unit under this mode of oscillation, and Oscillatory mode shape is to judge the influence degree of mode of oscillation to unit, judge whether multiple units can carry out the important indicator of dynamic equivalent.