CN114566967A - Fast decomposition method load flow calculation method suitable for research purpose - Google Patents

Fast decomposition method load flow calculation method suitable for research purpose Download PDF

Info

Publication number
CN114566967A
CN114566967A CN202210109034.9A CN202210109034A CN114566967A CN 114566967 A CN114566967 A CN 114566967A CN 202210109034 A CN202210109034 A CN 202210109034A CN 114566967 A CN114566967 A CN 114566967A
Authority
CN
China
Prior art keywords
node
delta
nodes
reactive power
unbalance
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202210109034.9A
Other languages
Chinese (zh)
Inventor
姚玉斌
李震
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Dalian Maritime University
Original Assignee
Dalian Maritime University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Dalian Maritime University filed Critical Dalian Maritime University
Priority to CN202210109034.9A priority Critical patent/CN114566967A/en
Publication of CN114566967A publication Critical patent/CN114566967A/en
Pending legal-status Critical Current

Links

Images

Classifications

    • H—ELECTRICITY
    • H02—GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02J—CIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J3/00—Circuit arrangements for ac mains or ac distribution networks
    • H02J3/04—Circuit arrangements for ac mains or ac distribution networks for connecting networks of the same frequency but supplied from different sources
    • H02J3/06—Controlling transfer of power between connected networks; Controlling sharing of load between connected networks
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00—Details relating to the application field
    • G06F2113/04—Power grid distribution networks
    • H—ELECTRICITY
    • H02—GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02J—CIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J2203/00—Indexing scheme relating to details of circuit arrangements for AC mains or AC distribution networks
    • H02J2203/10—Power transmission or distribution systems management focussing at grid-level, e.g. load flow analysis, node profile computation, meshed network optimisation, active network management or spinning reserve management
    • H—ELECTRICITY
    • H02—GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02J—CIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J2203/00—Indexing scheme relating to details of circuit arrangements for AC mains or AC distribution networks
    • H02J2203/20—Simulating, e g planning, reliability check, modelling or computer assisted design [CAD]
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02E—REDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
    • Y02E60/00—Enabling technologies; Technologies with a potential or indirect contribution to GHG emissions mitigation

Abstract

The invention discloses a quick decomposition method load flow calculation method suitable for research purposes, and provides a Newton method load flow calculation method which is easy to modify and maintain for scientific research personnel who further research on the basis of quick decomposition method load flow calculation. For an electric power system with n nodes, the node numbers are independent of the node types and can be randomly numbered, the coefficient matrixes B 'and B' of load flow calculation are stored according to n multiplied by n, the elements of the row corresponding to the delta P and the column corresponding to the delta theta in the B 'are all 0, the elements of the row corresponding to the PV node or the delta Q and the column corresponding to the delta U in the B' are all 0, but the diagonal elements are not cleared, and the original calculated values are reserved. Therefore, the linear equation trigonometric decomposition method with excellent performance provided by the programming language can be used for solving the correction equation of the load flow calculation, the programming difficulty is reduced, the calculation speed is increased, and the reliability of the algorithm is ensured.

Description

Fast decomposition method load flow calculation method suitable for research purpose
Technical Field
The invention relates to a fast decomposition method load flow calculation method for an electric power system, in particular to a fast decomposition method load flow calculation method suitable for research purposes.
Background
The power system load flow calculation is a basic calculation for studying the steady-state operation of a power system, and determines the operation state of the whole power system according to the given operation condition and network structure of the power system. The power flow calculation is also the basis of other analyses of the power system, such as safety analysis, transient stability analysis and the like. Because of the advantages of reliable convergence, fast calculation speed and less memory requirement, the fast decomposition method becomes one of the mainstream methods of current flow calculation, and researchers often make further research on the basis of fast decomposition method flow calculation. Practical commercial software adopts high-level technologies such as sparse matrix technology, node optimized numbering and the like. Although the technologies can greatly improve the speed of tidal current calculation and reduce the memory occupation amount, programming is very troublesome, modification and maintenance are difficult, and new functions are not easy to increase, so that the technology is not suitable for being used for research purposes by scientific researchers.
The power equation of the load flow calculation node by the rapid decomposition method is as follows:
Figure BDA0003494688160000011
in the formula, Pi、QiRespectively the active power and the reactive power of the node i; u shapei、UkThe voltage amplitudes of the node i and the node k are respectively; thetaik＝θi-θk，θiAnd thetakVoltage phase angles of the node i and the node k respectively; g ik、BikRespectively, node admittance matrix elements YikReal and imaginary parts of (c); n is the number of nodes.
Setting No. 1-m nodes as PQ nodes, No. m + 1-n-1 nodes as PV nodes, No. n nodes as balance nodes, and the equation of the power unbalance is as follows:
Figure BDA0003494688160000021
in the formula,. DELTA.Pi、ΔQiRespectively the active power unbalance amount and the reactive power unbalance amount of the node i; pis、QisRespectively giving an injection active power and an injection reactive power for the node i; m is the number of PQ nodes.
The fast decomposition method load flow calculation is a process of solving equation (2) to obtain the amplitude and phase angle of each node voltage and further calculating the power flowing through each branch.
The basic equations of the load flow calculation are nonlinear equations, and a successive linearization method is generally adopted for iterative solution. The linearized equation is called a correction equation and is used for solving the correction quantity of the voltage amplitude and the phase angle. The quick decomposition method correction equation is obtained by decoupling and improving on the basis of a polar coordinate Newton method load flow calculation correction equation.
The quick decomposition method correction equation is as follows:
-B′Δθ＝ΔP/U (3)
-B″ΔU＝ΔQ/U (4)
in the formula, delta P/U is a column vector obtained by dividing the active power unbalance by the voltage amplitude, and the dimension is n-1; delta Q/U is a column vector obtained by dividing the reactive power unbalance by the voltage amplitude, and the dimension is m; delta theta is a column vector of the n-1 dimensional voltage phase angle correction quantity; delta U is a m-dimensional voltage amplitude correction quantity column vector; b 'is the imaginary part of the simplified admittance matrix, branch resistance, admittance to the ground and nonstandard transformation ratio are not considered when the simplified admittance matrix is calculated, and the B' comprises the related rows and columns of the PQ node and the PV node and is a (n-1) x (n-1) order matrix; b' is the imaginary part of the admittance matrix, including only the rows and columns associated with the PQ nodes, which is an m by m order matrix.
As shown in fig. 1, the conventional fast decomposition method flow calculation method mainly includes the following steps:
A. raw data and an initialization voltage are input.
According to the characteristics of the nodes of the power system, the nodes of the power system are divided into 3 types by load flow calculation: the node with known active power and reactive power and unknown node voltage amplitude and voltage phase angle is injected into the node and is called as a PQ node; the node with known node injection active power and voltage amplitude and unknown node injection reactive power and voltage phase angle is called PV node; the node with known node voltage amplitude and voltage phase angle and unknown node injection active power and reactive power is called a balance node.
The voltage initialization adopts flat start, namely the voltage amplitudes of the PV node and the balance node are set values, and the voltage amplitude of the PQ node is 1.0; the voltage phase angles of all nodes take 0.0. The voltage phase angle is here in units of radians, other quantities take per unit values.
B. A node admittance matrix is formed.
C. Coefficient matrices B' and B "of the correction equation are formed and triangulated.
D. Set iteration count t to 0 and set Δ Qmax＝10εQ。
E. Calculating the active power unbalance amount delta P of the node(t)Calculating the maximum unbalance amount delta P of the active powermax。
The balance nodes do not participate in P-theta iterative calculation, and the active power unbalance of the nodes does not need to be calculated.
And solving the maximum value of the absolute value in the active power unbalance amount of each node, which is called as the active power maximum unbalance amount.
F. Judging absolute value | delta P of maximum active power unbalancemaxWhether | is less than convergence precision ∈P(ii) a If less than the convergence accuracy εPGo to step H; otherwise, go to step G.
G. Solving equation (3) to obtain Delta theta(t)Correcting the voltage phase angle according to the formula (5), and turning to the step I;
Figure BDA0003494688160000031
in the formula, superscript (t) represents the t iteration; delta thetaiIs the voltage phase angle correction column vector for node i.
H. Judging absolute value | delta Q of maximum unbalance amount of reactive powermaxWhether | is less than convergence precision εQ(ii) a Such asResult less than convergence precision epsilonQTurning to the step N; otherwise, go to step I.
I. Calculating the node reactive power unbalance amount delta Q(t)Calculating the maximum unbalance amount Delta Q of the reactive powermax。
The balance node and the PV node do not participate in Q-U iterative computation, and the reactive power unbalance of the nodes does not need to be computed.
And the maximum absolute value of the reactive power unbalance of each node is called as the maximum reactive power unbalance.
J. Judging absolute value | delta Q of maximum unbalance amount of reactive powermaxWhether | is less than convergence precision εQ(ii) a If less than the convergence accuracy εQGo to step L; otherwise, go to step K.
K. Obtaining Delta U by solving the formula (4)(t)Correcting the voltage amplitude according to the formula (6), and turning to the step M;
Figure BDA0003494688160000041
in the formula, superscript (t) represents the t-th iteration; delta UiThe column vector is modified for the voltage magnitude at node i.
L, judging the absolute value | delta P of the maximum unbalance of the active powermaxWhether | is less than convergence precision ∈P(ii) a If less than the convergence accuracy εPTurning to the step N; otherwise, go to step M.
And M, making t equal to t +1, and returning to the step E for the next iteration.
And N, calculating the power of the balance node, the reactive power of the PV node and the branch power.
And O, outputting a calculation result.
In the linear equation shown in the formula (3), the PQ node and the PV node have Δ P equations, and the balance node does not have Δ P equations. Therefore, there is not a Δ P equation for every node, nor is it necessary to solve for Δ θ for every node. In the linear equation shown in equation (4), the PQ node has a Δ Q equation, and the balance node and the PV node do not have a Δ Q equation. Therefore, there is not a Δ Q equation for every node, nor is it necessary to solve for Δ U for every node. Therefore, the corresponding relation between the node and the equation needs to be determined according to the type of the node, the type is changed, the corresponding relation is also changed, and great troubles are brought to programming and debugging. Although the first m nodes can be defined as PQ nodes, nodes m +1 to n-1 are PV nodes, and node n is a balance node. However, this provides that when data is input, it may be necessary to renumber the nodes according to this specification; in addition, when the node type changes, the corresponding relationship between the node and the equation needs to be readjusted.
Therefore, the Chinese patent CN201010585176.X provides a rapid decomposition method flow calculation method for an electric power system, and provides a rapid decomposition method flow calculation algorithm which is easy to modify and maintain for scientific research personnel who perform further research on the basis of rapid decomposition method flow calculation. The B ' matrix and the B ' matrix of the patent trend calculation are stored according to n multiplied by n, elements of a row corresponding to a balance node delta P in the B ' and a column corresponding to a delta theta are all 0, elements of a row corresponding to a balance node delta Q in the B ' and elements of a column corresponding to a row corresponding to a delta U in the B ' are all 0, so that although the memory requirement is increased, the corresponding relation between the nodes and the lines and columns of the equation coefficient matrix can be simplified, the programming difficulty is greatly reduced, and the calculation amount is not increased. When the B 'matrix and the B' matrix are subjected to triangular decomposition, a row with a main diagonal element of 0 in the coefficient matrix is skipped through judgment (the main diagonal element is 0 indicates that the row elements are all 0, and no corresponding equation exists), and no processing is performed.
Although the method of the chinese patent cn201010585176.x simplifies the corresponding relationship between nodes and rows and columns of equation coefficient matrixes, and greatly reduces the programming difficulty, when solving the load flow calculation correction equation, special processing is required, and the solution cannot be performed by using a trigonometric decomposition method with excellent performance provided by a programming language.
Disclosure of Invention
In order to solve the problems in the prior art, the invention provides a fast decomposition method flow calculation method suitable for research purposes, and provides a fast decomposition method flow calculation algorithm which is easy to modify and maintain for scientific researchers who carry out further research on the basis of fast decomposition method flow calculation.
To achieve the above object, the present invention modifies the coefficient matrices B' and B ″ as follows: coefficient matrixes B 'and B' are still stored according to n multiplied by n, in B ', elements of a row corresponding to a balance node delta P and a column corresponding to delta theta are all 0, in B', elements of a row corresponding to a PV node or a balance node delta Q and a column corresponding to delta U are all 0, but diagonal elements are not cleared, and an original calculated value is reserved. Therefore, the correction equation of the power flow calculation can be solved by using the excellent performance trigonometric decomposition method provided by the programming language.
The technical scheme of the invention is as follows: a fast decomposition method flow calculation method suitable for research purposes comprises the following steps:
A. raw data and an initialization voltage are input.
B. A node admittance matrix is formed.
C. Coefficient matrices B' and B "of the correction equation are formed and triangulated.
B 'is the imaginary part of the simplified admittance matrix, branch resistance, admittance to the ground and nonstandard transformation ratio are not considered when the simplified admittance matrix is calculated, and the B' comprises the related rows and columns of the PQ node and the PV node and is a (n-1) x (n-1) order matrix; b' is the imaginary part of the admittance matrix, including only the rows and columns associated with the PQ nodes, which is an m by m order matrix.
The steps for forming the coefficient matrices B' and B "are as follows:
c1, calculating B 'elements by all nodes according to the PQ node types to form an n multiplied by n order B' matrix;
c2, setting the node count i to 1;
c3, judging whether the node i is a balanced node, if not, turning to the step C10;
c4, setting count k to 1;
c5, judging whether k is true or not, if so, going to step C8;
c6, let coefficient matrix B 'ith row and kth column element B'i,k＝0；
C7, line B'k,i＝0；
C8, let k be k + 1;
c9, judging whether k is larger than n, if k is larger than n, turning to the step C10; otherwise, returning to step C5;
c10, let i ═ i + 1;
c11, judging whether i is larger than n, if i is larger than n, turning to the step C12; otherwise, return to step C3.
C12, all nodes are used as PQ nodes to take the imaginary part of the admittance matrix to form a B 'matrix, and an n multiplied by n order B' matrix is formed;
c13, setting the node count i to 1;
c14, judging whether the node i is a PQ node, if so, turning to the step C21;
c15, setting count k to 1;
c16, judging whether k is true or not, if so, going to step C19;
c17, let coefficient matrix B "row i, column k element Bi,k＝0；
C18, order B ″k,i＝0；
C19, let k be k + 1;
c20, judging whether k is larger than n, if k is larger than n, turning to the step C21; otherwise, returning to step C16;
C21, let i ═ i + 1;
c22, judging whether i is larger than n, if i is larger than n, turning to step D; otherwise, return to step C14.
D. Set iteration count t to 0 and Δ Qmax＝10εQ
E. Calculating the active power unbalance amount delta P of the node(t)Calculating the maximum unbalance amount delta P of the active powermax。
And (4) calculating the active power unbalance of the node according to the formula (7).
Figure BDA0003494688160000071
In the formula,. DELTA.PiThe active power unbalance amount of the node i is obtained; pisThe injected active power is given for node i.
The balance nodes do not participate in P-theta iterative calculation, and the active power unbalance of the nodes does not need to be calculated.
Delta P is n-dimensional vector, and the active power unbalance quantity delta P corresponding to the balance node iiAnd setting 0.
The maximum absolute value of the active power unbalance of each node is calculated and used as the maximum active power unbalance delta Pmax。
F. Judging absolute value | delta P of maximum unbalance amount of active powermaxWhether | is less than convergence precision εP(ii) a If less than the convergence accuracy εPGo to step H; otherwise, go to step G.
G. Solving P-theta iterative correction equation to obtain delta theta(t)Correcting the voltage phase angle according to the formula (9), and turning to the step I;
the quick decomposition method P-theta iteration correction equation is as follows:
-B′Δθ＝ΔP/U (8)
in the formula, delta P/U is a column vector obtained by dividing the active power unbalance by the voltage amplitude, and the dimension is n; and delta theta is a n-dimensional voltage phase angle correction quantity column vector.
The voltage phase angle correction formula is as follows:
Figure BDA0003494688160000072
in the formula, superscript (t) represents the t-th iteration; delta thetaiIs the voltage phase angle correction column vector for node i.
H. Judging absolute value | delta Q of maximum unbalance amount of reactive powermaxWhether | is less than convergence precision εQ(ii) a If less than the convergence accuracy εQTurning to the step N; otherwise, go to step I.
I. Calculating the node reactive power unbalance amount delta Q(t)Calculating the maximum unbalance amount Delta Q of the reactive powermax。
And (5) calculating the node reactive power unbalance according to the formula (10).
Figure BDA0003494688160000081
In the formula,. DELTA.QiIs a section ofThe reactive power unbalance amount of the point i; qisThe node i is given the injected reactive power.
The balance node and the PV node do not participate in Q-U iterative computation, and the reactive power unbalance of the nodes does not need to be computed.
Delta Q is n-dimensional vector, and the reactive power unbalance quantity Delta Q corresponding to the balance node i or PV node iiAnd setting 0.
Calculating the maximum absolute value of the reactive power unbalance of each node as the maximum reactive power unbalance delta Qmax。
J. Judging absolute value | delta Q of maximum unbalance amount of reactive powermaxWhether | is less than convergence precision εQ(ii) a If less than the convergence accuracy εQGo to step L; otherwise, go to step K.
K. Solving the Q-U iterative correction equation to obtain delta U(t)Modifying the voltage amplitude according to the formula (12), and turning to the step M; the quick decomposition method Q-U iteration correction equation is as follows:
-B″ΔU＝ΔQ/U (11)
In the formula, delta Q/U is a column vector obtained by dividing the reactive power unbalance by the voltage amplitude, and the dimension is n; Δ U is an n-dimensional voltage amplitude correction quantity column vector.
The voltage amplitude correction formula is as follows:
Figure BDA0003494688160000082
in the formula, superscript (t) represents the t iteration; delta UiThe column vector is modified for the voltage magnitude at node i.
L, judging the absolute value | delta P of the maximum unbalance amount of the active powermaxWhether | is less than convergence precision εP(ii) a If less than the convergence accuracy εPTurning to the step N; otherwise, go to step M.
And M, making t equal to t +1, and returning to the step E for the next iteration.
And N, calculating the power of the balance node, the reactive power of the PV node and the branch power.
And O, outputting a calculation result.
Compared with the prior art, the invention has the following beneficial effects:
1. the coefficient matrixes B 'and B' of load flow calculation are stored according to n multiplied by n, so that the corresponding relation between the nodes and the rows and columns of the equation coefficient matrixes can be simplified, the programming difficulty is greatly reduced, and the calculation amount is not increased.
2. For an equation without a balance node, other elements except diagonal elements in a corresponding row of the coefficient matrix B' are all 0 and are used for representing that no corresponding equation exists; for an equation where no PV node or balance node exists, the other elements in the row corresponding to the coefficient matrix B ″ except for the diagonal elements are all 0, which is used to indicate that no corresponding equation exists. Because the diagonal element is not 0, a trigonometric decomposition method provided by a programming language can be used for solving a correction equation of load flow calculation, and the programming difficulty is greatly reduced. The trigonometric decomposition method provided by the programming language is optimized, the calculation speed is high, the stability is high, and the calculation speed and the program stability of the load flow calculation can be greatly improved.
3. The nodes of the load flow calculation are flexibly numbered, the numbering is irrelevant to the node type, the balancing node is the last node without requiring the PQ node to be numbered in the front, and the original numbering of the system nodes is not required to be changed according to the program design requirement.
Drawings
The invention is shown in figure 2. Wherein:
fig. 1 is a flowchart of a conventional fast decomposition flow calculation.
FIG. 2 is a block flow diagram of the formation of coefficient matrices B 'and B' in accordance with the present invention.
Detailed Description
The invention is further described below with reference to the accompanying drawings. According to the flow chart of the fast decomposition method power flow calculation shown in the figure 1 and the flow chart formed by the coefficient matrixes B 'and B' of the invention shown in the figure 2, the power flow calculation is carried out on a practical large-scale power grid. The actual large-scale power grid is provided with 445 nodes which contain a large number of small impedance branches, and in order to enable the conventional load flow calculation method to calculate, the small impedance branches are changed into normal branches. Convergence accuracy epsilon of load flow calculationPAnd epsilonQAre all 0.00001.
For comparison, the following 3 methods are adopted to perform load flow calculation on the actual large-scale power grid:
the conventional method comprises the following steps: a conventional fast decomposition method load flow calculation method;
the patented method comprises the following steps: the method of patent cn201010585176. x;
The comparison method comprises the following steps: the method of patent cn201010585176.x is adopted, but the solution of the correction equation of the load flow calculation adopts a triangular decomposition method of Matlab.
All the load flow calculation methods are programmed by adopting M files of Matlab, and the calculation time of different load flow calculation methods is shown in a table 1.
TABLE 1 calculation time(s) for different load flow calculation methods
Method Conventional methods Patented process Comparison method The invention
Calculating time 0.424438 0.399523 The coefficient matrix is singular and cannot be solved 0.151596
As can be seen from Table 1, for the modified 445-node actual power system calculation example, the calculation time of the conventional fast decomposition method load flow calculation method and the patent method is close, and the method of the invention adopts the trigonometric decomposition method provided by Matlab programming language for calculation when solving the equation, so that the calculation speed is high. Some diagonal elements in the coefficient matrix of the comparison method are 0, so that the coefficient matrix of the correction equation is singular and cannot be solved by a trigonometric decomposition method provided by a programming language.
The present invention can be implemented using any programming language and programming environment, such as the M-file programming languages of C, C + +, FORTRAN, Delphi, MATLAB, and the like. The development environment may employ Visual C + +, Borland C + + Builder, Visual FORTRAN, MATLAB, and the like.
The present invention is not limited to the embodiment, and any equivalent idea or change within the technical scope of the present invention is to be regarded as the protection scope of the present invention.

Claims (1)

1. A fast decomposition method flow calculation method suitable for research purposes is characterized in that: the method comprises the following steps:
A. inputting original data and an initialization voltage;
according to the characteristics of the nodes of the power system, the nodes of the power system are divided into 3 types by load flow calculation: the node with known active power and reactive power and unknown node voltage amplitude and voltage phase angle is injected into the node and is called as a PQ node; the node with known active power and voltage amplitude injected by the node and unknown reactive power and voltage phase angle injected by the node is called a PV node; the node with known node voltage amplitude and voltage phase angle and unknown node injection active power and reactive power is called a balance node;
the voltage initialization adopts flat start, namely the voltage amplitudes of the PV node and the balance node are set values, and the voltage amplitude of the PQ node is 1.0; the voltage phase angles of all the nodes are 0.0; the unit of the voltage phase angle is radian, and other quantities adopt per unit values;
B. forming a node admittance matrix;
C. forming coefficient matrixes B 'and B' of a correction equation and carrying out triangular decomposition;
b 'is the imaginary part of the simplified admittance matrix, branch resistance, admittance to the ground and nonstandard transformation ratio are not considered when the simplified admittance matrix is calculated, and the B' comprises the related rows and columns of the PQ node and the PV node and is a (n-1) x (n-1) order matrix; b' is the imaginary part of the admittance matrix, including only the rows and columns associated with the PQ nodes, which is an m × m order matrix; n is the number of nodes; m is the number of PQ nodes;
The steps for forming the coefficient matrices B' and B "are as follows:
c1, calculating B 'elements by all nodes according to the PQ node types to form an n multiplied by n order B' matrix;
c2, setting the node count i to 1;
c3, judging whether the node i is a balanced node, if not, turning to the step C10;
c4, setting count k to 1;
c5, judging whether k is true or not, if so, going to step C8;
c6, let coefficient matrix B 'ith row and kth column element B'i,k＝0；
C7, line B'k,i＝0；
C8, let k be k + 1;
c9, judging whether k is larger than n, if k is larger than n, turning to the step C10; otherwise, returning to step C5;
c10, let i ═ i + 1;
c11, judging whether i is larger than n, if i is larger than n, turning to the step C12; otherwise, returning to step C3;
c12, all nodes are used as PQ nodes to take the imaginary part of the admittance matrix to form a B 'matrix, and an n multiplied by n order B' matrix is formed;
c13, setting the node count i to 1;
c14, judging whether the node i is a PQ node, if so, turning to the step C21;
c15, setting count k to 1;
c16, judging whether k is true or not, if so, going to step C19;
c17, let coefficient matrix B "row i, column k element Bi,k＝0；
C18, order B ″k,i＝0；
C19, let k be k + 1;
c20, judging whether k is larger than n, if k is larger than n, turning to the step C21; otherwise, returning to step C16;
C21, let i ═ i + 1;
c22, judging whether i is larger than n, if i is larger than n, turning to step D; otherwise return to step C14;
D. set iteration count t to 0 and set Δ Qmax＝10εQ
E. Calculating the active power unbalance amount delta P of the node(t)Calculating the maximum unbalance amount Delta P of the active powermax；
Calculating the active power unbalance of the node according to the formula (1);
Figure FDA0003494688150000021
in the formula,. DELTA.PiThe active power unbalance amount of the node i is obtained; pisInjecting active power given for the node i; piIs the active power of node i; u shapei、UkThe voltage amplitudes of the node i and the node k are respectively; thetaik＝θi-θk，θiAnd thetakVoltage phase angles of the node i and the node k respectively; gik、BikAre respectively node admittance matrix elements YikThe real and imaginary parts of (c);
the balance nodes do not participate in P-theta iterative calculation, and the active power unbalance of the nodes does not need to be calculated;
delta P is n-dimensional vector, and the active power unbalance quantity delta P corresponding to the balance node iiSetting 0;
calculating the maximum value of the absolute value in the active power unbalance amount of each node as the maximum active power unbalance amount delta Pmax；
F. Judging absolute value | delta P of maximum unbalance amount of active powermaxWhether | is less than convergence precision εP(ii) a If less than the convergence accuracy εPGo to step H; otherwise, go to step G;
G. solving P-theta iterative correction equation to obtain delta theta (t)Correcting the voltage phase according to equation (3)Turning to the step I;
the quick decomposition method P-theta iteration correction equation is as follows:
-B′Δθ＝ΔP/U (2)
in the formula, delta P/U is a column vector obtained by dividing the active power unbalance by the voltage amplitude, and the dimension is n; delta theta is a n-dimensional voltage phase angle correction sequence vector;
the voltage phase angle correction formula is as follows:
Figure FDA0003494688150000031
in the formula, superscript (t) represents the t iteration; delta thetaiA voltage phase angle correction quantity column vector of a node i;
H. judging absolute value | delta Q of maximum unbalance amount of reactive powermaxWhether | is less than convergence precision εQ(ii) a If less than the convergence accuracy εQTurning to the step N; otherwise, go to step I;
I. calculating the node reactive power unbalance amount delta Q(t)Calculating the maximum unbalance amount Delta Q of the reactive powermax；
Calculating the reactive power unbalance of the nodes according to the formula (4);
Figure FDA0003494688150000032
in the formula,. DELTA.QiIs the reactive power unbalance of the node i; qisInjecting reactive power given for node i; qiRespectively the reactive power of the node i;
the balance node and the PV node do not participate in Q-U iterative computation, and the reactive power unbalance of the nodes does not need to be computed;
delta Q is n-dimensional vector, and the reactive power unbalance quantity Delta Q corresponding to the balance node i or PV node iiSetting 0;
calculating the maximum absolute value of the reactive power unbalance of each node as the maximum reactive power unbalance delta Q max；
J. Judging absolute value | delta Q of maximum unbalance of reactive powermaxWhether | is less than convergence precision ∈Q(ii) a If less than the convergence accuracy εQGo to step L; otherwise, go to step K;
K. solving the Q-U iterative correction equation to obtain delta U(t)Correcting the voltage amplitude according to the formula (6), and turning to the step M;
the quick decomposition method Q-U iteration correction equation is as follows:
-B″ΔU＝ΔQ/U (5)
in the formula, delta Q/U is a column vector obtained by dividing the reactive power unbalance by the voltage amplitude, and the dimension is n; delta U is an n-dimensional voltage amplitude correction quantity column vector;
the voltage amplitude correction formula is as follows:
Figure FDA0003494688150000041
in the formula, superscript (t) represents the t iteration; delta UiA voltage amplitude correction quantity column vector of a node i;
l, judging the absolute value | delta P of the maximum unbalance amount of the active powermaxWhether | is less than convergence precision εP(ii) a If less than the convergence accuracy εPTurning to the step N; otherwise, go to step M;
m, making t equal to t +1, and returning to the step E for the next iteration;
n, calculating the power of a balance node, the reactive power of a PV node and the branch power;
and O, outputting a calculation result.
CN202210109034.9A 2022-01-29 2022-01-29 Fast decomposition method load flow calculation method suitable for research purpose Pending CN114566967A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210109034.9A CN114566967A (en) 2022-01-29 2022-01-29 Fast decomposition method load flow calculation method suitable for research purpose

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210109034.9A CN114566967A (en) 2022-01-29 2022-01-29 Fast decomposition method load flow calculation method suitable for research purpose

Publications (1)

Publication Number Publication Date
CN114566967A true CN114566967A (en) 2022-05-31

Family

ID=81713523

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210109034.9A Pending CN114566967A (en) 2022-01-29 2022-01-29 Fast decomposition method load flow calculation method suitable for research purpose

Country Status (1)

Country Link
CN (1) CN114566967A (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116760039A (en) * 2023-08-16 2023-09-15 山东大学 Node type expanded AC/DC power flow calculation method, system, equipment and medium

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116760039A (en) * 2023-08-16 2023-09-15 山东大学 Node type expanded AC/DC power flow calculation method, system, equipment and medium
CN116760039B (en) * 2023-08-16 2023-11-28 山东大学 Node type expanded AC/DC power flow calculation method, system, equipment and medium

Similar Documents

Publication Publication Date Title
CN106356859B (en) A kind of rectangular coordinate Newton load flow calculation method based on Matlab
Ghatak et al. A fast and efficient load flow technique for unbalanced distribution system
CN106532711A (en) Newton method power flow calculation method which changes Jacobi matrix with iteration and node types
CN109617080B (en) Rectangular coordinate Newton method load flow calculation method based on improved Jacobian matrix
CN103632046A (en) Power grid load flow calculation method
CN114566967A (en) Fast decomposition method load flow calculation method suitable for research purpose
CN106856327B (en) A kind of compensation of line series containing small impedance branches algorithm quicksort tidal current computing method
CN106410811B (en) Iteration small impedance branches endpoint changes the tidal current computing method of Jacobian matrix for the first time
CN116166915A (en) Improved GMRES-based power system power flow calculation method and system
CN106532712A (en) Rectangular coordinate Newton method load flow calculation method for small-impedance-branch-containing power grid based on compensation method
CN104967121A (en) Large-scale electric power system node load flow computing method
CN107846022B (en) Large-scale power distribution network power flow analysis method based on ILUTP preprocessing parallel iteration method
CN106712029B (en) The Newton load flow calculation method of small impedance branches PQ endpoint change Jacobian matrix
CN113890042B (en) Load flow calculation method based on active power distribution network
CN106529089B (en) Penalty method algorithm quicksort tidal current computing method for the power grid containing small impedance branches
CN111049144A (en) Active power distribution network load flow calculation PV node processing method based on constant change rate matrix
CN108123434B (en) Method for calculating slope of PV curve to obtain operating point of PV curve
CN110649624B (en) Power system load flow parallel computing method
CN110336287B (en) Three-phase load flow calculation method of power distribution system based on Jacobian element extraction
CN114566968A (en) Polar coordinate Newton method load flow calculation method suitable for research purpose
CN110336288B (en) Three-phase power flow calculation method for power distribution network based on extraction of Jacobian elements through matrix operation
CN107658880B (en) The algorithm quicksort coefficient matrix calculation method of operation based on correlation matrix
CN114566969A (en) Rectangular coordinate Newton method load flow calculation method suitable for research purpose
CN108183487A (en) A kind of distribution power flow rapid analysis method based on linear mathematical model
CN107944682B (en) Load flow calculation admittance matrix calculation method based on Matlab matrix operation

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination