A fast decomposition power flow calculation method suitable for research purposes

Abstract

The present invention discloses a fast decomposition method for power flow calculation suitable for research purposes, and provides a Newton method for power flow calculation that is easy to modify and maintain for scientific researchers who conduct further research based on the fast decomposition method for power flow calculation. For an n-node power system, the node number is independent of the node type and can be arbitrarily numbered. The coefficient matrices B′ and B″ of the power flow calculation are stored as n×n. The elements of the rows corresponding to the balance node ΔP and the columns corresponding to Δθ in B′ are all 0, and the elements of the rows corresponding to the PV node or the balance node ΔQ and the columns corresponding to ΔU in B″ are also all 0, but the diagonal elements are not cleared to zero, and the original calculated values ​​are retained. In this way, the linear equation triangular decomposition method with excellent performance provided by the programming language can be used to solve the modified equation for power flow calculation, which reduces the programming difficulty, improves the calculation speed, and ensures the reliability of the algorithm.

Description

Technical Field

[0001] The present invention relates to a power system rapid decomposition method for power flow calculation, in particular to a rapid decomposition method for power flow calculation suitable for research purposes. Background Art

[0002] Power system flow calculation is a fundamental calculation used to study the steady-state operation of power systems. It determines the operating state of the entire power system based on given operating conditions and network structure. Flow calculations also serve as the foundation for other power system analyses, such as security analysis and transient stability analysis. Due to its advantages of reliable convergence, fast computational speed, and minimal memory requirements, the rapid decomposition method has become one of the mainstream methods for current flow calculations, and researchers often use this method as a foundation for further research. Practical commercial software utilizes advanced techniques such as sparse matrix technology and optimized node numbering. While these techniques can significantly increase the speed of flow calculations and reduce memory usage, they are cumbersome to program, difficult to modify and maintain, and difficult to add new functionality, making them unsuitable for research purposes.

[0003] The node power equation for power flow calculation using the fast decomposition method is:

[0004]

[0005] Where, P i , Q i are the active power and reactive power of node i respectively; U i 、U k are the voltage amplitudes of nodes i and k respectively; θ ik =θ i -θ k ,θ i and θ k are the voltage phase angles of node i and node k respectively; G ik 、B ik are the node admittance matrix elements Y ik The real and imaginary parts of ; n is the number of nodes.

[0006] Assume that nodes 1 to m are PQ nodes, nodes m+1 to n-1 are PV nodes, and node n is a balance node. The power imbalance equation is:

[0007]

[0008] Where, ΔP i , ΔQ i are the active power imbalance and reactive power imbalance of node i respectively; P is , Q isare the injected active power and injected reactive power given to node i respectively; m is the number of PQ nodes.

[0009] The rapid decomposition method power flow calculation is the process of solving equation (2) to obtain the amplitude and phase angle of each node voltage, and then calculate the power flowing through each branch.

[0010] The basic equations for power flow calculation are nonlinear and are typically solved iteratively using a successive linearization method. The resulting linearized equations are called correction equations, which are used to determine voltage amplitude and phase angle corrections. The correction equations using the rapid decomposition method are derived by decoupling and improving upon the correction equations for power flow calculation using the polar coordinate Newton method.

[0011] The corrected equation of the rapid decomposition method is:

[0012]

[0013] Where ΔP / U is the column vector of active power imbalance divided by voltage amplitude, with dimension n-1; ΔQ / U is the column vector of reactive power imbalance divided by voltage amplitude, with dimension m; Δθ is the n-1-dimensional column vector of voltage phase angle correction; ΔU is the m-dimensional column vector of voltage amplitude correction; B′ is the imaginary part of the simplified admittance matrix. The simplified admittance matrix calculation does not take into account branch resistance, ground admittance, and non-standard transformation ratios. B′ includes rows and columns related to the PQ node and the PV node and is an (n-1)×(n-1)-order matrix; B″ is the imaginary part of the admittance matrix, which includes only rows and columns related to the PQ node and is an m×m-order matrix.

[0014] like Figure 1 As shown in Figure 2, the conventional rapid decomposition method for power flow calculation mainly includes the following steps:

[0015] A. Input original data and initialization voltage.

[0016] According to the characteristics of power system nodes, power flow calculation divides power system nodes into three categories: nodes where the active power and reactive power injected by the node are known, but the node voltage amplitude and voltage phase angle are unknown are called PQ nodes; nodes where the active power and voltage amplitude injected by the node are known, but the node reactive power and voltage phase angle are unknown are called PV nodes; nodes where the voltage amplitude and voltage phase angle are known, but the node active power and reactive power injected by the node are unknown are called balance nodes.

[0017] Voltage initialization uses a flat start, meaning the voltage amplitudes at the PV and balancing nodes are set to given values, the voltage amplitude at the PQ node is set to 1.0, and the voltage phase angles at all nodes are set to 0.0. The voltage phase angles are expressed in radians, and all other quantities are expressed in per-unit values.

[0018] B. Form a node admittance matrix.

[0019] C. Form the coefficient matrices B′ and B″ of the modified equation and perform triangular decomposition.

[0020] D. Set the iteration count t = 0 and set ΔQ max =10ε Q .

[0021] E. Calculate the node active power imbalance ΔP (t) , calculate the maximum unbalanced active power ΔP max .

[0022] The balanced nodes do not participate in the P~θ iterative calculation, and there is no need to calculate the node active power imbalance.

[0023] The maximum absolute value of the active power imbalance of each node is called the maximum active power imbalance.

[0024] F. Determine the absolute value of the maximum unbalanced active power |ΔP max Is it less than the convergence accuracy ε? P ; If it is less than the convergence accuracy ε P , go to step H; otherwise, go to step G.

[0025] G. Solve equation (3) to get Δθ (t) , correct the voltage phase angle according to formula (5) and go to step I;

[0026]

[0027] Where, the superscript (t) indicates the tth iteration; Δθ i is the column vector of voltage phase angle correction at node i.

[0028] H. Determine the absolute value of the maximum unbalanced reactive power |ΔQ max Is it less than the convergence accuracy ε? Q ; If it is less than the convergence accuracy ε Q , go to step N; otherwise, go to step I.

[0029] I. Calculate the node reactive power imbalance ΔQ (t) , calculate the maximum unbalanced reactive power ΔQ max .

[0030] The balancing nodes and PV nodes do not participate in the Q~U iterative calculation, and there is no need to calculate the node reactive power unbalance.

[0031] The value with the largest absolute value among the reactive power unbalance values ​​of each node is called the maximum reactive power unbalance value.

[0032] J. Determine the absolute value of the maximum unbalanced reactive power |ΔQ max Is it less than the convergence accuracy ε?Q ; If it is less than the convergence accuracy ε Q , go to step L; otherwise, go to step K.

[0033] K, solve equation (4) to get ΔU (t) , correct the voltage amplitude according to formula (6) and go to step M;

[0034]

[0035] Where, the superscript (t) indicates the tth iteration; ΔU i is the column vector of voltage amplitude correction at node i.

[0036] L. Determine the absolute value of the maximum unbalanced active power |ΔP max Is it less than the convergence accuracy ε? P ; If it is less than the convergence accuracy ε P , go to step N; otherwise, go to step M.

[0037] M. Set t=t+1 and return to step E for the next iteration.

[0038] N. Calculate the balancing node power, PV node reactive power, and branch power.

[0039] O. Output the calculation results.

[0040] In the linear equation shown in Equation (3), the PQ nodes and PV nodes have ΔP equations, while the equilibrium nodes do not. Therefore, not every node has a ΔP equation, and not every node needs to solve for Δθ. In the linear equation shown in Equation (4), the PQ nodes have ΔQ equations, while the equilibrium nodes and PV nodes do not. Therefore, not every node has a ΔQ equation, and not every node needs to solve for ΔU. Thus, the correspondence between nodes and equations must be determined based on the node type. When the type changes, the correspondence also changes, which brings great trouble to programming and debugging. Although it is possible to specify the first m nodes as PQ nodes, nodes m+1 to n-1 as PV nodes, and node n as the equilibrium node, this specification may require renumbering the nodes according to this specification when inputting data; in addition, when the node type changes, the correspondence between nodes and equations must be readjusted.

[0041] To this end, Chinese patent CN201010585176.X proposes a method for calculating power system flow using a rapid decomposition method, providing researchers who conduct further research based on rapid decomposition method flow calculation with an algorithm for rapid decomposition method flow calculation that is easy to modify and maintain. The B′ and B″ matrices for flow calculation in this patent are stored as n×n, and the elements of the rows corresponding to the balancing node ΔP and the columns corresponding to Δθ in B′ are all 0. The elements of the rows corresponding to the balancing node and PV node ΔQ and the columns corresponding to ΔU in B″ are all 0. Although this increases the memory requirement, it can simplify the correspondence between the nodes and the rows and columns of the equation coefficient matrix, greatly reducing the programming difficulty and not increasing the amount of calculation. When performing triangular decomposition of the B′ and B″ matrices, the rows whose main diagonal elements in the coefficient matrix are 0 are skipped by judgment (the main diagonal elements being 0 means that the elements in this row are all 0 and there is no corresponding equation) and are not processed.

[0042] Although the method in Chinese patent CN201010585176.X simplifies the correspondence between nodes and rows and columns of the equation coefficient matrix, greatly reducing the programming difficulty, it requires special processing when solving the power flow calculation correction equation, and cannot be solved using the high-performance triangular decomposition method provided by the programming language. Summary of the Invention

[0043] In order to solve the above problems existing in the prior art, the present invention proposes a rapid decomposition method for tidal current calculation suitable for research purposes, and provides a rapid decomposition method for tidal current calculation that is easy to modify and maintain for scientific researchers who conduct further research based on rapid decomposition method for tidal current calculation.

[0044] To achieve the aforementioned objectives, the present invention improves the coefficient matrices B′ and B″ as follows: Coefficient matrices B′ and B″ are still stored as n×n. In B′, the elements in the rows corresponding to the balancing node ΔP and the columns corresponding to Δθ are all zero. In B″, the elements in the rows corresponding to the PV node or balancing node ΔQ and the columns corresponding to ΔU are also all zero. However, the diagonal elements are not cleared to zero, retaining their original calculated values. This allows the use of the high-performance triangular decomposition method provided by the programming language to solve the modified equations for power flow calculations.

[0045] The technical solution of the present invention is as follows: A rapid decomposition method for power flow calculation suitable for research purposes comprises the following steps:

[0046] A. Input original data and initialization voltage.

[0047] B. Form a node admittance matrix.

[0048] C. Form the coefficient matrices B′ and B″ of the modified equation and perform triangular decomposition.

[0049] B′ is the imaginary part of the simplified admittance matrix. The simplified admittance matrix does not take into account branch resistance, admittance to ground, and non-standard transformation ratio. B′ contains rows and columns related to the PQ node and the PV node, and is an (n-1)×(n-1)-order matrix. B″ is the imaginary part of the admittance matrix, which only includes rows and columns related to the PQ node and is an m×m-order matrix.

[0050] The steps to form the coefficient matrices B′ and B″ are as follows:

[0051] C1. To simplify the correspondence between node numbers and rows and columns of matrix B′, i.e., the i-th row of matrix B′ corresponds to the active power unbalance of node i, and the i-th column of matrix B′ corresponds to the voltage phase angle correction of node i, all nodes are calculated according to the PQ node type. The B′ elements form an n×n order B′ matrix;

[0052] C2. Set node count i = 1;

[0053] C3. Determine whether node i is a balanced node. If not, go to step C10.

[0054] C4. Set count k = 1;

[0055] C5. Determine whether k=i. If so, go to step C8.

[0056] C6. Let the coefficient matrix B′ be the i-th row and k-th column element B′ i,k =0;

[0057] C7, let B′ k,i =0;

[0058] C8. Let k = k + 1;

[0059] C9. Determine whether k is greater than n. If k is greater than n, go to step C10; otherwise, return to step C5;

[0060] C10. Let i=i+1;

[0061] C11. Determine whether i is greater than n. If i is greater than n, go to step C12; otherwise, return to step C3.

[0062] C12. To simplify the correspondence between node numbers and rows and columns of the matrix B″, i.e., the i-th row of the matrix B″ corresponds to the reactive power unbalance of node i, and the i-th column of the matrix B″ corresponds to the voltage amplitude correction of node i, all nodes are treated as PQ nodes and the imaginary part of the admittance matrix is ​​taken to form the B″ matrix, forming an n×n order B″ matrix;

[0063] C13. Set node count i=1;

[0064] C14. Determine whether node i is a PQ node. If it is a PQ node, go to step C21.

[0065] C15, set count k = 1;

[0066] C16, determine whether k=i. If yes, go to step C19;

[0067] C17, let the coefficient matrix B″ be the i-th row and k-th column element B″ i,k =0;

[0068] C18, let B″ k,i =0;

[0069] C19, let k = k + 1;

[0070] C20, determine whether k is greater than n, if k is greater than n, go to step C21; otherwise return to step C16;

[0071] C21, let i = i + 1;

[0072] C22. Determine whether i is greater than n. If i is greater than n, go to step D; otherwise, return to step C14.

[0073] D. Set the iteration count t = 0 and set ΔQ max =10ε Q

[0074] E. Calculate the node active power imbalance ΔP (t) , calculate the maximum unbalanced active power ΔP max .

[0075] Calculate the node active power imbalance according to formula (7).

[0076]

[0077] Where, ΔP i is the active power imbalance of node i; P is is the injected active power given by node i.

[0078] The balanced nodes do not participate in the P~θ iterative calculation, and there is no need to calculate the node active power imbalance.

[0079] ΔP is an n-dimensional vector, and the active power imbalance ΔP corresponding to the balancing node i i Set to 0.

[0080] The maximum absolute value of the active power imbalance of each node is taken as the maximum active power imbalance ΔP max .

[0081] F. Determine the absolute value of the maximum unbalanced active power |ΔP max Is it less than the convergence accuracy ε? P ; If it is less than the convergence accuracy ε P , go to step H; otherwise, go to step G.

[0082] G. Solve the P~θ iterative correction equation to obtain Δθ (t) , correct the voltage phase angle according to formula (9) and go to step I;

[0083] The P~θ iterative correction equation of the rapid decomposition method is:

[0084] -B′Δθ=ΔP / U (8)

[0085] Where ΔP / U is the column vector of active power imbalance divided by voltage amplitude, with dimension n; Δθ is the n-dimensional column vector of voltage phase angle correction.

[0086] The voltage phase angle correction formula is:

[0087]

[0088] Where, the superscript (t) indicates the tth iteration; Δθ i is the column vector of voltage phase angle correction at node i.

[0089] H. Determine the absolute value of the maximum unbalanced reactive power |ΔQ max Is it less than the convergence accuracy ε? Q ; If it is less than the convergence accuracy ε Q , go to step N; otherwise, go to step I.

[0090] I. Calculate the node reactive power imbalance ΔQ (t) , calculate the maximum unbalanced reactive power ΔQ max .

[0091] The node reactive power unbalance is calculated according to formula (10).

[0092]

[0093] Where ΔQ i is the reactive power imbalance of node i; Q is is the injected reactive power given by node i. Balancing nodes and PV nodes do not participate in the Q-U iterative calculation, and there is no need to calculate the node reactive power unbalance.

[0094] ΔQ is an n-dimensional vector, and the reactive power unbalance ΔQ corresponding to the balancing node i or PV node i is i Set to 0. Find the maximum absolute value of the reactive power imbalance of each node as the maximum reactive power imbalance ΔQmax .

[0095] J. Determine the absolute value of the maximum unbalanced reactive power |ΔQ max Is it less than the convergence accuracy ε? Q ; If it is less than the convergence accuracy ε Q , go to step L; otherwise, go to step K.

[0096] K. Solve the Q~U iterative correction equation to obtain ΔU (t) , correct the voltage amplitude according to formula (12) and go to step M;

[0097] The Q~U iterative correction equation of the fast decomposition method is:

[0098] -B″ΔU=ΔQ / U (11)

[0099] Where ΔQ / U is the column vector of reactive power unbalance divided by voltage amplitude, with dimension n; ΔU is the n-dimensional column vector of voltage amplitude correction.

[0100] The voltage amplitude correction formula is:

[0101]

[0102] Where, the superscript (t) indicates the tth iteration; ΔU i is the column vector of voltage amplitude correction at node i.

[0103] L. Determine the absolute value of the maximum unbalanced active power |ΔP max Is it less than the convergence accuracy ε? P ; If it is less than the convergence accuracy ε P , go to step N; otherwise, go to step M.

[0104] M. Set t=t+1 and return to step E for the next iteration.

[0105] N. Calculate the balancing node power, PV node reactive power, and branch power.

[0106] O. Output the calculation results.

[0107] Compared with the prior art, the present invention has the following beneficial effects:

[0108] 1. The coefficient matrices B′ and B″ of the power flow calculation of the present invention are stored in n×n format, which can simplify the correspondence between nodes and rows and columns of the equation coefficient matrix, greatly reducing the programming difficulty and not increasing the amount of calculation.

[0109] 2. For equations without a balancing node, all elements except the diagonal elements in the corresponding row of the coefficient matrix B′ are zero, indicating that the corresponding equation does not exist. For equations without a PV node or balancing node, all elements except the diagonal elements in the corresponding row of the coefficient matrix B″ are zero, indicating that the corresponding equation does not exist. Because the diagonal elements are not zero, the triangular decomposition method provided by the programming language can be used to solve the corrected equations for power flow calculations, greatly reducing the programming difficulty. The triangular decomposition method provided by the programming language has been optimized for fast calculation speed and high stability, which can greatly improve the calculation speed and program stability of power flow calculations.

[0110] 3. The node numbering of the power flow calculation of the present invention is flexible and has nothing to do with the node type. The PQ node number is not required to be in the front, and the balance node is the last node. There is no need to change the original numbering of the system nodes according to the program design requirements. BRIEF DESCRIPTION OF THE DRAWINGS

[0111] The present invention has the following Figure 2 Zhang. Among them:

[0112] Figure 1 This is the flow chart of the conventional rapid decomposition method for power flow calculation.

[0113] Figure 2 It is a flowchart of the formation of coefficient matrices B′ and B″ of the present invention. DETAILED DESCRIPTION

[0114] The present invention will be further described below with reference to the accompanying drawings. Figure 1 The flowchart of the fast decomposition method power flow calculation is shown in the figure and Figure 2 The flow chart formed by the coefficient matrices B′ and B″ of the present invention is shown in the figure. A flow chart is used to calculate the power flow of an actual large-scale power grid. The actual large-scale power grid has 445 nodes and contains a large number of small impedance branches. In order to make the conventional power flow calculation method applicable, these small impedance branches are converted to normal branches. The convergence accuracy of the power flow calculation is ε P and ε Q Both are 0.00001.

[0115] For comparison, the following three methods were used to calculate the power flow of the actual large-scale power grid:

[0116] Conventional method: conventional rapid decomposition method for power flow calculation;

[0117] Patented method: Patent CN201010585176.X;

[0118] Comparison method: The method of patent CN201010585176.X is adopted, but the modified equation for the power flow calculation is solved using the triangular decomposition method of Matlab.

[0119] Various power flow calculation methods are programmed using Matlab M files. The calculation time of different power flow calculation methods is shown in Table 1.

[0120] Table 1 Calculation time (s) of different power flow calculation methods

[0121]

[0122] As shown in Table 1, for the modified 445-node real-world power system example, the conventional fast decomposition method and the patented method have similar calculation times. The proposed method, which uses the triangular decomposition method provided by the Matlab programming language to solve the equation, has a faster calculation speed. However, some diagonal elements in the coefficient matrix of the comparison method are zero, resulting in a singular coefficient matrix for the modified equation, making it impossible to solve using the triangular decomposition method provided by the programming language.

[0123] The present invention can be implemented using any programming language and programming environment, such as C language, C++, FORTRAN, Delphi, MATLAB M-file programming language, etc. The development environment can be Visual C++, Borland C++Builder, Visual FORTRAN, MATLAB, etc.

[0124] The present invention is not limited to this embodiment, and any equivalent concepts or modifications within the technical scope disclosed by the present invention are included in the protection scope of the present invention.

Claims

1. A rapid decomposition method for power flow calculation suitable for research purposes, characterized by: The following steps are involved: A. Input original data and initialization voltage; According to the characteristics of power system nodes, power flow calculation divides power system nodes into three categories: nodes with known active power and reactive power injection but unknown voltage amplitude and voltage phase angle are called PQ nodes; nodes with known active power and voltage amplitude but unknown reactive power injection and voltage phase angle are called PV nodes; nodes with known voltage amplitude and voltage phase angle but unknown active power and reactive power injection are called balancing nodes; The voltage initialization adopts flat start, that is, the voltage amplitude of the PV node and the balance node takes the given value, and the voltage amplitude of the PQ node takes 1.0; the voltage phase angle of all nodes is 0.0; the voltage phase angle unit here is radian, and other quantities use per unit value; B. Forming a node admittance matrix; C. Form the coefficient matrices B′ and B″ of the modified equation and perform triangular decomposition; B′ is the imaginary part of the simplified admittance matrix. The simplified admittance matrix does not take into account branch resistance, admittance to ground, and non-standard transformation ratio. B′ contains rows and columns related to the PQ nodes and PV nodes, and is an (n-1)×(n-1)-order matrix. B″ is the imaginary part of the admittance matrix, which only includes rows and columns related to the PQ nodes and is an m×m-order matrix. n is the number of nodes; m is the number of PQ nodes. The steps to form the coefficient matrices B′ and B″ are as follows: C1. To simplify the correspondence between node numbers and rows and columns of matrix B′, i.e., the i-th row of matrix B′ corresponds to the active power unbalance of node i, and the i-th column of matrix B′ corresponds to the voltage phase angle correction of node i, all nodes are calculated according to the PQ node type. The B′ elements form an n×n order B′ matrix; C2. Set node count i = 1; C3. Determine whether node i is a balanced node. If not, go to step C10. C4. Set count k = 1; C5. Determine whether k=i. If so, go to step C8. C6. Let the coefficient matrix B′ have the element B in row i and column k. i ′,k=0; C7. Let B k ′, i = 0; C8. Let k = k + 1; C9, determine whether k is greater than n, if k is greater than n, go to step C10; Otherwise, return to step C5; C10. Let i=i+1; C11. Determine whether i is greater than n. If i is greater than n, go to step C12; otherwise, return to step C3; C12. To simplify the correspondence between node numbers and rows and columns of the matrix B″, i.e., the i-th row of the matrix B″ corresponds to the reactive power unbalance of node i, and the i-th column of the matrix B″ corresponds to the voltage amplitude correction of node i, all nodes are treated as PQ nodes and the imaginary part of the admittance matrix is ​​taken to form the B″ matrix, forming an n×n order B″ matrix; C13. Set node count i=1; C14. Determine whether node i is a PQ node. If it is a PQ node, go to step C21. C15, set count k = 1; C16, determine whether k=i. If yes, go to step C19; C17, let the coefficient matrix B″ have the element B in row i and column k i ′,′ k =0; C18, let B k ″,i=0; C19, let k = k + 1; C20, determine whether k is greater than n, if k is greater than n, go to step C21; otherwise return to step C16; C21, let i = i + 1; C22. Determine whether i is greater than n. If i is greater than n, go to step D; otherwise, return to step C14. D. Set the iteration count t = 0 and set ΔQ max =10ε Q E. Calculate the node active power imbalance ΔP (t) , calculate the maximum unbalanced active power ΔP max ; Calculate the node active power imbalance according to formula (1); Where ΔP i is the active power imbalance of node i; P is is the injected active power given by node i; P i is the active power of node i; U i 、U k are the voltage amplitudes of nodes i and k respectively; θ ik =θ i -θ k ,θ i and θ k are the voltage phase angles of node i and node k respectively; G ik 、B ik are the node admittance matrix elements Y ik The real and imaginary parts of The balancing nodes do not participate in the P~θ iterative calculation, and there is no need to calculate the node active power imbalance; ΔP is an n-dimensional vector, and the active power imbalance ΔP corresponding to the balancing node i i Set to 0; Find the maximum absolute value of the active power imbalance of each node as the maximum active power imbalance ΔP max ; F. Determine the absolute value of the maximum unbalanced active power |ΔP max Is it less than the convergence accuracy ε? P ; If it is less than the convergence accuracy ε P , go to step H; otherwise, go to step G; G. Solve the P~θ iterative correction equation to obtain Δθ (t) , correct the voltage phase angle according to formula (3) and go to step I; The P~θ iterative correction equation of the rapid decomposition method is: -B′Δθ=ΔP / U (2) Where ΔP / U is the column vector of active power imbalance divided by voltage amplitude, with dimension n; Δθ is the n-dimensional column vector of voltage phase angle correction; The voltage phase angle correction formula is: Where, the superscript (t) indicates the tth iteration; Δθ i is the column vector of voltage phase angle correction at node i; H. Determine the absolute value of the maximum unbalanced reactive power |ΔQ max Is it less than the convergence accuracy ε? Q ; If it is less than the convergence accuracy ε Q , go to step N; Otherwise, go to step I; I. Calculate the node reactive power imbalance ΔQ (t) , calculate the maximum unbalanced reactive power ΔQ max ; Calculate the node reactive power unbalance according to formula (4); Where ΔQ i is the reactive power imbalance of node i; Q is is the injected reactive power given by node i; Q i are the reactive power of node i respectively; The balancing nodes and PV nodes do not participate in the Q-U iterative calculation, and there is no need to calculate the node reactive power unbalance; ΔQ is an n-dimensional vector, and the reactive power unbalance ΔQ corresponding to the balancing node i or PV node i is i Set to 0; Find the maximum absolute value of the reactive power imbalance of each node as the maximum reactive power imbalance ΔQ max ; J. Determine the absolute value of the maximum unbalanced reactive power |ΔQ max Is it less than the convergence accuracy ε? Q ; If it is less than the convergence accuracy ε Q , go to step L; otherwise, go to step K; K. Solve the Q~U iterative correction equation to obtain ΔU (t) , correct the voltage amplitude according to formula (6) and go to step M; The Q~U iterative correction equation of the fast decomposition method is: -B″ΔU=ΔQ / U (5) Where ΔQ / U is the column vector of reactive power unbalance divided by voltage amplitude, with dimension n; ΔU is the n-dimensional column vector of voltage amplitude correction; The voltage amplitude correction formula is: Where, the superscript (t) indicates the tth iteration; ΔU i is the column vector of voltage amplitude correction at node i; L. Determine the absolute value of the maximum unbalanced active power |ΔP max Is it less than the convergence accuracy ε? P ; If it is less than the convergence accuracy ε P , go to step N; Otherwise, go to step M; M. Set t = t + 1 and return to step E for the next iteration; N. Calculate the balancing node power and PV node reactive power and branch power; O. Output the calculation results.

Images