CN105391057B - A kind of GPU threaded design methods that direction of energy Jacobi battle array calculates - Google Patents
A kind of GPU threaded design methods that direction of energy Jacobi battle array calculates Download PDFInfo
- Publication number
- CN105391057B CN105391057B CN201510809809.3A CN201510809809A CN105391057B CN 105391057 B CN105391057 B CN 105391057B CN 201510809809 A CN201510809809 A CN 201510809809A CN 105391057 B CN105391057 B CN 105391057B
- Authority
- CN
- China
- Prior art keywords
- msub
- mrow
- matrix
- nodes
- node
- 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.)
- Expired - Fee Related
Links
- 239000011159 matrix material Substances 0.000 claims abstract description 83
- 238000004364 calculation method Methods 0.000 claims abstract description 32
- 230000001174 ascending Effects 0.000 claims description 15
- 230000000875 corresponding Effects 0.000 claims description 9
- 238000003860 storage Methods 0.000 claims description 9
- 238000002347 injection Methods 0.000 claims description 6
- 239000007924 injection Substances 0.000 claims description 6
- 210000001503 Joints Anatomy 0.000 claims description 2
- 238000006073 displacement reaction Methods 0.000 claims description 2
- 238000000034 method Methods 0.000 abstract description 5
- 238000009826 distribution Methods 0.000 description 4
- 230000001133 acceleration Effects 0.000 description 3
- 230000000694 effects Effects 0.000 description 2
- 238000004519 manufacturing process Methods 0.000 description 2
- 238000011030 bottleneck Methods 0.000 description 1
- 238000004422 calculation algorithm Methods 0.000 description 1
- 238000005094 computer simulation Methods 0.000 description 1
- 230000003247 decreasing Effects 0.000 description 1
- 230000005611 electricity Effects 0.000 description 1
- 230000005520 electrodynamics Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 241000894007 species Species 0.000 description 1
- 230000003068 static Effects 0.000 description 1
- 230000002123 temporal effect Effects 0.000 description 1
- 230000001550 time effect Effects 0.000 description 1
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
-
- 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]
Abstract
The invention discloses a kind of GPU threaded design methods that direction of energy Jacobi battle array calculates, it is characterised in that：Methods described includes：Input electric network data, pretreatment node admittance battle array Y；CPU distinguishes calculated sub-matrix nonzero element and node admittance battle array Y nonzero element position mapping tables；CPU distinguishes calculated sub-matrix nonzero element and Jacobian matrix nonzero element position mapping table；GPU interior joint injecting power kernel functions S calculates all node injecting powers；Jacobi submatrix calculates kernel function and distinguishes nonzero element in calculated sub-matrix and be stored in Jacobian matrix in GPU.The present invention is calculated in Jacobian matrix non-zero entry calculating process by submatrix to judge the branched structure of the affiliated submatrix region process of element when avoiding and directly being calculated using single kernel function, lifts execution efficiency；And by centralized calculation there is nonzero element in the submatrix of identical calculations formula to solve the unbalanced problem of threads load, lift parallel efficiency calculation.
Description
Technical field
The invention belongs to High performance computing in power system application field, more particularly to a kind of direction of energy Jacobi battle array to calculate
GPU threaded design methods.
Background technology
Load flow calculation is most widely used, most basic and most important a kind of electric computing in power system.In power train
In the research of the method for operation of uniting and programme, it is required for carrying out Load flow calculation to compare the method for operation or plan power supply plan
Feasibility, reliability and economy.Meanwhile for the running status of real-time electric power monitoring system, it is also desirable to carry out a large amount of and fast
The Load flow calculation of speed.Therefore, in the method for operation of programming and planning and schedule system, using offline Load flow calculation；In electricity
In the real-time monitoring of Force system running status, then calculated using online power flow.
And in actual production process, offline trend and online power flow calculating all have this to compare the calculating speed of trend
High requirement.It is being related to planning and designing and complexity situations such as in the offline trend for arranging the method for operation, landing scheme because of equipment, is needing
Want the species of simulation run more, Load flow calculation amount is big, and single Load flow calculation time effects integrally emulate duration；And in power system
The online power flow carried out in operation is calculated to calculating temporal sensitivity height, it is necessary to provide calculation of tidal current in real time, is such as being envisioned
In accident, the equipment Load flow calculation out of service to the influence of static security, system needs to calculate trend under a large amount of forecast accidents
Distribution, and the method for operation Adjusted Option of anticipation is made in real time.
In traditional Newton-Laphson method Load flow calculation, the calculating time of Jacobian matrix accounts for trend and amounts to evaluation time
30%, the calculating speed of Jacobian matrix influences the overall performance of program.And slowing down with the lifting of CPU calculating speeds, existing rank
The single Load flow calculation calculating time of section has reached a bottleneck.The accelerated method of Load flow calculation, which is concentrated mainly on, at present makes
More trends are carried out with coarseness acceleration with cluster and multiple-core server, it is actual into being ground to single trend internal arithmetic acceleration in production
Study carefully less.
GPU is a kind of many-core parallel processor, will be considerably beyond CPU in the quantity of processing unit.GPU traditionally is only
Responsible figure renders, and CPU has all been given in most processing.Method battle array is a kind of multinuclear to present GPU, multithreading, is had
There are powerful calculating ability and high bandwidth of memory, programmable processor.Under universal computer model, associations of the GPU as CPU
Processor works, and is decomposed by task reasonable distribution and completes high-performance calculation.
Jacobian matrix, which calculates, has concurrency, and the calculating of each of which nonzero element is separate, does not rely on pass
System, it can naturally be handled by parallel calculating, be adapted to GPU to accelerate.Therefore can be fast by rational scheduling between CPU and GPU
Speed completes the calculating of Jacobian matrix, and the method that domestic and foreign scholars have begun to carry out GPU Jacobi acceleration is studied,
But without deep optimization threaded design, computational threads design is studied from the distribution of amount of calculation merely, to thread calculating side
Formula, data directory mode are not furtherd investigate, and program can not be made to give full play to GPU advantage.
It would therefore be highly desirable to solve the above problems.
The content of the invention
Goal of the invention：In view of the shortcomings of the prior art, when Jacobian matrix calculating can be greatly decreased in present invention offer one kind
Between and can be lifted Load flow calculation speed a kind of direction of energy Jacobi battle array calculate GPU threaded design methods.
Technical scheme：The present invention proposes a kind of a kind of tide of electric power based on GPU that computing resource is distributed using mapping relations
Flow the GPU threaded design methods that Jacobian matrix calculates.
Load flow calculation：Electrodynamic noun, refer in given power system network topology, component parameters and generating, load parameter
Under the conditions of, calculate the distribution of active power, reactive power and voltage in power network.
Parallel computation：It is a kind of algorithm that once can perform multiple instruction, it is therefore an objective to improve and calculate relative to serial arithmetic
Speed, and by expanding problem solving scale, solve large-scale and complicated computational problem.
GPU：Graphics processor (English：GraphicsProcessingUnit, abbreviation：GPU).
Admittance matrix：Established based on the Equivalent admittance of system element, each node voltage of description electric power networks and
The matrix of relation between Injection Current.
The invention discloses a kind of GPU threaded design methods that direction of energy Jacobi battle array calculates, methods described includes：
(1) electric network data, pretreatment node admittance battle array Y are inputted；
(2) CPU calculates the sub- matrix non-zero elements of H, N, J, L tetra- respectively and node admittance battle array Y nonzero element positions map
Relation table MH2Y、MN2Y、MJ2Y、ML2Y, circular comprises the following steps：
(2.1) the node admittance battle array Y of COO forms is scanned, chooses and wherein corresponds to the element addition set that H-matrix calculates
HY, generation inquiry table MH2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to active reactive PQ nodes or active voltage PV node
When, remember that i=pqpv [x], wherein pqpv are PQ, the index of PV node, x is call number；Row number j is scanned again, when j belongs to PQ or PV
During node, remember j=pqpv [y], choose the element to be added to set HY；
B) element is often added to set HYWhen, set of records ends HYCurrentElement quantity k, and the element is matrix Y's
Position h in COO formsk；Record i position x in PQ, PV nodek, j position y in PQ, PV nodek；Set HYMiddle element representation
For { k, hk, xk, yk}；
C) k → h is generatedkMapping table MH2Y；
(2.2) the node admittance battle array Y of COO forms is scanned, chooses and wherein corresponds to the element addition set that N matrix calculates
NY, generation mapping table MN2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to PQ or PV node, remembers i=pqpv [x], then sweep
Row number j is retouched, when j belongs to PQ nodes, remembers that j=pq [y], wherein pq are the index of PQ nodes, y is call number, chooses the element
It is added to set NY；
B) element is often added to set NYWhen, set of records ends NYCurrentElement quantity k, and the element is matrix Y's
Position n in COO formsk；Record i position x in PQ, PV nodek, j position y in PQ nodesk；Set NYMiddle element representation is
{ k, nk, xk, yk}；
C) k → n is generatedkMapping table MN2Y；
(2.3) the node admittance battle array Y of COO forms is scanned, chooses and wherein adds set corresponding to the element of J matrix computations
MY, generation mapping table MJ2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to PQ nodes, remembers i=pq [x], then scan row number
J, when j belongs to PQ or PV node, remember j=pqpv [y], choose the element to be added to set MY；
B) element is often added to set JYWhen, set of records ends NYCurrentElement quantity k, and the element is matrix Y's
Position m in COO formsk；Record i position x in PQ nodesk, j position y in PQ, PV nodek；Set JYMiddle element representation is
{ k, jk, xk, yk}；
C) k → j is generatedkMapping table MJ2Y；
(2.4) the node admittance battle array Y of COO forms is scanned, chooses and wherein adds set corresponding to the element of L matrix computations
LY, generation inquiry table ML2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to PQ nodes, remembers i=pq [x], then scan row number
J, when j belongs to PQ nodes, remember j=pq [y], choose the element to be added to set LY；
B) element is often added to set LYWhen, set of records ends LYCurrentElement quantity k, and the element is matrix Y's
Position l in COO formsk；Record i position x in PQ nodesk, j position y in PQ nodesk；Set LYMiddle element representation for k,
lk, xk, yk}；
C) k → l is generatedkMapping table ML2Y；
(3) CPU calculates the sub- matrix non-zero elements of H, N, J, L tetra- respectively and the mapping of Jacobian matrix nonzero element position is closed
It is table MH2Ja、MN2Ja、MJ2Ja、ML2Ja, circular comprises the following steps：
(3.1) according to the H obtained in step 2Y、NY、JY、LYIn xk, ykObtain matrix H, N, J, L sparse format
(3.2) definition set Cmp { x, y, NO, zone }, wherein x are Jacobian matrix line number, y is Jacobian matrix row
Number, NO be element position in H, N, J, L Matrix C OO forms, equivalent to the k in step 2, in zone H, N, J, L subregion
The identifier of element, PQ, PV node number are npqpv；
(3.3) element in matrix H, N, J, L is inserted in Cmp, wherein element zone is set to 1 in H battle arrays；Element y values in N battle arrays
Npqpv, zone is added to be set to 2；Element x value adds npqpv, zone to be set to 3 in J battle arrays；Element x value adds npqpv in L battle arrays, and y values add
Npqpv, zone are set to 4；
(3.4) set Cmp is reordered by x ascending orders, y ascending orders, obtains the Jacobi Matrix C OO forms of row ascending order；
In set Cmp add variable num, represent Jacobian matrix in position of the element in standard COO forms, and in order assign 1,
2、3…；
(3.5) to Cmp { x, y, NO, zone, num } by zone ascending orders, NO ascending sorts, the num in each zone regions
Sequence is inquiry table MH2Ja、MN2Ja、MJ2Ja、ML2Ja；
(4) GPU interior joints injecting power kernel function S calculates all node injecting powers；
(5) Jacobi submatrix calculates kernel function H, M, N, L and calculates non-zero in tetra- submatrixs of H, N, J, L respectively in GPU
Element is simultaneously stored in Jacobian matrix.
Wherein, the nodes of power network are n in the step (1), the node admittance battle array Y { i, j, G, B } of power network with line number i,
Row number j, conductance G, susceptance B Coordinate sparse formats storage, by relative skew of the node admittance battle array each row of data in Y
Position is stored in line displacement vector R, and definition node voltage phase angle vector is θ, and voltage magnitude vector is V.
Wherein, the k threads that the calculation of step (4) the interior joint injecting power is kernel function S calculate k nodes
The injection of active and reactive power, calculation formula is：
Wherein,
K is parallel thread number in GPU Programmings；
PkInjected for the active power of k nodes；
QkInjected for the reactive power of k nodes；
kbIndexed for the starting of admittance battle array Y row k nonzero elements, kb=R [k]；
keIndexed for the termination of admittance battle array Y row k nonzero elements, ke=R [k+1] -1；
J, G, B are the data of l-th of position storage in admittance battle array Y；
VkFor the voltage magnitude of k nodes；
VjFor the voltage magnitude of j nodes；
θkFor the voltage phase angle of k nodes；
θjFor the voltage phase angle of j nodes.
Wherein, in the step (5), the computational methods of nonzero element are in H submatrixs：
According to the H submatrixs nonzero element and node admittance battle array Y mapping table M obtained in step (2)H2Y, H core letters
When several t threads calculate, inquiry table MH2YT positions take out Y battle arrays call number kH2Y, take out Y [kH2Y] corresponding in element
Data { i, j, G, B } are calculated, and the calculation formula of t-th of nonzero element of Jacobian matrix neutron matrix H is：
Wherein,
T is parallel thread number in GPU Programmings；HtFor the nonzero element of t-th of position in H submatrix COO forms；
I, j, G, B are kth in admittance battle array YH2YThe data of individual position storage；
QiFor the injection reactive power of node i；
According still further to the H submatrixs nonzero element and matrix J a mapping relations M obtained in step (3)H2Ja, inquire about MH2JaT
Number position take out HtStorage call number k of the result of calculation in Ja matrixesH2Ja, then by HtResult of calculation deposit Ja
[kH2Ja]。
Beneficial effect：Compared with the prior art, beneficial effects of the present invention are：First, the present invention uses GPU processors to refined
Gram accelerated than matrix computations, compared to the speed-up ratio that CPU obtains 16 times；Secondly, threaded design method of the invention is to GPU
The performance of program is optimized, and is calculated in Jacobian matrix non-zero entry calculating process by submatrix to avoid using single
The branched structure of the affiliated submatrix region process of element, lifting GPU branches execution efficiency are judged when kernel function directly calculates；And lead to
Crossing centralized calculation, there is nonzero element in the submatrix of identical calculations formula to solve the unbalanced problem of threads load, and lifting is parallel
Computational efficiency.
Brief description of the drawings
Fig. 1 is the flow chart for the GPU threaded design methods that a kind of direction of energy Jacobi battle array of the present invention calculates
Embodiment
As shown in figure 1, the present invention proposes a kind of a kind of tide of electric power based on GPU that computing resource is distributed using mapping relations
The GPU threaded design methods that Jacobian matrix calculates are flowed, the method includes the steps of：
If the nodes of power network are n, the node admittance battle array Y { i, j, G, B } of power network is with line number, row number, conductance, susceptance
COO sparse formats store, and index in COO forms are converted into CompressedSparseRow (CSR) row compressed format, row is inclined
Shifting array vector is R, and node voltage phase angle vector is θ, and voltage magnitude vector is V；
Step 1：Generate inquiry table MH2Y、MN2Y、MJ2Y、ML2Y
1) the node admittance battle array Y of COO forms is scanned, chooses and wherein corresponds to the element addition set H that H-matrix calculatesY, it is raw
Into inquiry table MH2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to active reactive PQ nodes or active voltage PV node
When, remember that i=pqpv [x], wherein pqpv are PQ, the index of PV node, x is call number；Row number j is scanned again, when j belongs to PQ or PV
During node, remember j=pqpv [y], choose the element to be added to set HY；
B) element is often added to set HYWhen, set of records ends HYCurrentElement quantity k, and the element is matrix Y's
Position h in COO formsk；Record i position x in PQ, PV nodek, j position y in PQ, PV nodek；Set HYMiddle element representation
For { k, hk, xk, yk}；
C) k → h is generatedkMapping table MH2Y；
2) the node admittance battle array Y of COO forms is scanned, chooses and wherein corresponds to the element addition set N that N matrix calculatesY, it is raw
Into mapping table MN2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to PQ or PV node, remembers i=pqpv [x], then sweep
Row number j is retouched, when j belongs to PQ nodes, remembers that j=pq [y], wherein pq are the index of PQ nodes, y is call number, chooses the element
It is added to set NY；
B) element is often added to set NYWhen, set of records ends NYCurrentElement quantity k, and the element is matrix Y's
Position n in COO formsk；Record i position x in PQ, PV nodek, j position y in PQ nodesk；Set NYMiddle element representation is
{ k, nk, xk, yk}；
C) k → n is generatedkMapping table MN2Y；
3) the node admittance battle array Y of COO forms is scanned, chooses and wherein adds set M corresponding to the element of J matrix computationsY, it is raw
Into mapping table MJ2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to PQ nodes, remembers i=pq [x], then scan row number
J, when j belongs to PQ or PV node, remember j=pqpv [y], choose the element to be added to set MY；
B) element is often added to set JYWhen, set of records ends NYCurrentElement quantity k, and the element is matrix Y's
Position m in COO formsk；Record i position x in PQ nodesk, j position y in PQ, PV nodek；Set JYMiddle element representation is
{ k, jk, xk, yk}；
C) k → j is generatedkMapping table MJ2Y；
4) the node admittance battle array Y of COO forms is scanned, chooses and wherein adds set L corresponding to the element of L matrix computationsY, it is raw
Into inquiry table ML2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to PQ nodes, remembers i=pq [x], then scan row number
J, when j belongs to PQ nodes, remember j=pq [y], choose the element to be added to set LY；
B) element is often added to set LYWhen, set of records ends LYCurrentElement quantity k, and the element is matrix Y's
Position l in COO formsk；Record i position x in PQ nodesk, j position y in PQ nodesk；Set LYMiddle element representation for k,
lk, xk, yk}；
C) k → l is generatedkMapping table ML2Y；
Step 2：Generate inquiry table MH2Ja、MN2Ja、MJ2Ja、ML2Ja
1) according to the H obtained in step 1Y、NY、JY、LYIn xk, ykObtain matrix H, N, J, L sparse format
2) definition set Cmp { x, y, NO, zone }, wherein x are Jacobian matrix line number, y is Jacobian matrix row number, NO
For element in H, N, J, L Matrix C OO forms position, equivalent to the k in step 1, zone H, N, J, L subregion interior element
Identifier, PQ, PV node number are npqpv；
3) element in matrix H, N, J, L is inserted in Cmp, wherein element zone is set to 1 in H battle arrays；Element y values add in N battle arrays
Npqpv, zone are set to 2；Element x value adds npqpv, zone to be set to 3 in J battle arrays；Element x value adds npqpv, y values plus npqpv in L battle arrays,
Zone is set to 4；
4) set Cmp is reordered by x ascending orders, y ascending orders, obtains the Jacobi Matrix C OO forms of row ascending order；Collecting
Close Cmp in add variable num, represent Jacobian matrix in position of the element in standard COO forms, and in order assign 1,2,
3…；
5) to Cmp { x, y, NO, zone, num } by zone ascending orders, NO ascending sorts, the num sequences in each zone regions
As inquiry table MH2Ja、MN2Ja、MJ2Ja、ML2Ja。
Step 3：Calculate node injecting power
S kernel functions calculate the active reactive injection of all nodes, and setting kernel function start-up parameter is<<<n/256+1,256>
>>, that is, start one-dimensional piece of block quantity [n/256]+1, open 256 thread thread in each block, amount to n
threads.Required data are passed in GPU in advance, the k threads of S kernel functions calculate the active and reactive power note of k nodes
Enter, thread calculation formula is：
Wherein,
kbIndexed for the starting of admittance battle array Y row k nonzero elements, kb=R [k]；
keIndexed for the termination of admittance battle array Y row k nonzero elements, ke=R [k+1] -1；
J, G, B are the data of l-th of position storage in admittance battle array Y.
Step 4：Calculate nonzero element value in Jacobian matrix
Jacobian matrix is calculated and is divided into four parts, respectively calculated sub-matrix H, M, J, L, each subfunction in calculating process
Thread by mapping directly read calculate data calculated.
By taking the calculating of submatrix H non-zero entries value as an example,
H1Kernel function is responsible for calculated sub-matrix H nonzero element, sets its start-up parameter to be<<<Hnum/256+1,256>>
>, that is, start one-dimensional block quantity [Hnum/ 256]+1, each block is interior to open 256 threads, amounts to HnumIt is individual
Threads, wherein HnumFor submatrix H non-zero entry number.Required data are passed in GPU in advance, H1The t lines of kernel function
When journey calculates, inquiry table MH2YT positions take out Y battle arrays call number kH2Y, take out Y [kH2Y] corresponding data in element i, j,
G, B } calculated, t-th of nonzero element that content is Jacobian matrix neutron matrix H is calculated, calculation formula is：
Wherein,
I, j, G, B are kth in admittance battle array YH2YThe data of individual position storage；
QiFor the injection reactive power of node i.
Establish H2Global the thread number t and matrix J acobi of kernel function mapping relations MH2Ja, H2The t threads of kernel function
After calculating result, M is inquired aboutH2JaT positions take out storage call number k of the result of calculation in Jacobi matrixesH2Ja, then
By result of calculation deposit Jacobi [kH2Ja]。
The calculating of other submatrixs is similar therewith.
Claims (4)
1. a kind of GPU threaded design methods that direction of energy Jacobi battle array calculates, it is characterised in that：Methods described includes：
(1) electric network data, pretreatment node admittance battle array Y are inputted；
(2) CPU calculates the sub- matrix non-zero elements of H, N, J, L tetra- and node admittance battle array Y nonzero element position mapping relations respectively
Table MH2Y、MN2Y、MJ2Y、ML2Y, circular comprises the following steps：
(2.1) the node admittance battle array Y of COO forms is scanned, chooses and wherein corresponds to the element addition set H that H-matrix calculatesY, generation
Inquiry table MH2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to active reactive PQ nodes or active voltage PV node,
Remember that i=pqpv [x], wherein pqpv are PQ, the index of PV node, x is call number；Row number j is scanned again, when j belongs to PQ or PV sections
During point, remember j=pqpv [y], choose the element to be added to set HY；
B) element is often added to set HYWhen, set of records ends HYCurrentElement quantity k, and the element is in matrix Y COO
Position h in formk；Record i position x in PQ, PV nodek, j position y in PQ, PV nodek；Set HYMiddle element representation is
{ k, hk, xk, yk}；
C) k → h is generatedkMapping table MH2Y；
(2.2) the node admittance battle array Y of COO forms is scanned, chooses and wherein corresponds to the element addition set N that N matrix calculatesY, generation
Mapping table MN2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to PQ or PV node, remembers i=pqpv [x], then scan columns
Number j, when j belongs to PQ nodes, remember that j=pq [y], wherein pq are the index of PQ nodes, y is call number, chooses the element to add
To set NY；
B) element is often added to set NYWhen, set of records ends NYCurrentElement quantity k, and the element is in matrix Y COO
Position n in formk；Record i position x in PQ, PV nodek, j position y in PQ nodesk；Set NYMiddle element representation for k,
nk, xk, yk}；
C) k → n is generatedkMapping table MN2Y；
(2.3) the node admittance battle array Y of COO forms is scanned, chooses and wherein adds set M corresponding to the element of J matrix computationsY, generation
Mapping table MJ2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to PQ nodes, remembers i=pq [x], then scans row number j, when
When j belongs to PQ or PV node, remember j=pqpv [y], choose the element to be added to set MY；
B) element is often added to set JYWhen, set of records ends NYCurrentElement quantity k, and the element is in matrix Y COO
Position m in formk；Record i position x in PQ nodesk, j position y in PQ, PV nodek；Set JYMiddle element representation for k,
jk, xk, yk}；
C) k → j is generatedkMapping table MJ2Y；
(2.4) the node admittance battle array Y of COO forms is scanned, chooses and wherein adds set L corresponding to the element of L matrix computationsY, generation
Inquiry table ML2Y；
A) the line number i of each element in COO forms is scanned, when it belongs to PQ nodes, remembers i=pq [x], then scans row number j, when
When j belongs to PQ nodes, remember j=pq [y], choose the element to be added to set LY；
B) element is often added to set LYWhen, set of records ends LYCurrentElement quantity k, and the element is in matrix Y COO
Position l in formk；Record i position x in PQ nodesk, j position y in PQ nodesk；Set LYMiddle element representation is { k, lk,
xk, yk}；
C) k → l is generatedkMapping table ML2Y；
(3) CPU calculates the sub- matrix non-zero elements of H, N, J, L tetra- and Jacobian matrix nonzero element position mapping table respectively
MH2Ja、MN2Ja、MJ2Ja、ML2Ja, circular comprises the following steps：
(3.1) according to the H obtained in step 2Y、NY、JY、LYIn xk, ykObtain matrix H, N, J, L sparse format
(3.2) definition set Cmp { x, y, NO, zone }, wherein x are Jacobian matrix line number, y is Jacobian matrix row number, NO
For element in H, N, J, L Matrix C OO forms position, equivalent to the k in step 2, zone H, N, J, L subregion interior element
Identifier, PQ, PV node number are npqpv；
(3.3) element in matrix H, N, J, L is inserted in Cmp, wherein element zone is set to 1 in H battle arrays；Element y values add in N battle arrays
Npqpv, zone are set to 2；Element x value adds npqpv, zone to be set to 3 in J battle arrays；Element x value adds npqpv, y values plus npqpv in L battle arrays,
Zone is set to 4；
(3.4) set Cmp is reordered by x ascending orders, y ascending orders, obtains the Jacobi Matrix C OO forms of row ascending order；Collecting
Close Cmp in add variable num, represent Jacobian matrix in position of the element in standard COO forms, and in order assign 1,2,
3…；
(3.5) to Cmp { x, y, NO, zone, num } by zone ascending orders, NO ascending sorts, the num sequences in each zone regions
As inquiry table MH2Ja、MN2Ja、MJ2Ja、ML2Ja；
(4) GPU interior joints injecting power kernel function S calculates all node injecting powers；
(5) Jacobi submatrix calculates kernel function H, M, N, L and calculates nonzero element in tetra- submatrixs of H, N, J, L respectively in GPU
And it is stored in Jacobian matrix.
2. the GPU threaded design methods that direction of energy Jacobi battle array according to claim 1 calculates, it is characterised in that：Institute
The nodes for stating power network in step (1) are n, and the node admittance battle array Y { i, j, G, B } of power network is with line number i, row number j, conductance G, susceptance
B Coordinate sparse formats storage, line displacement is stored in by relative deviation post of the node admittance battle array each row of data in Y
In vectorial R, definition node voltage phase angle vector is θ, and voltage magnitude vector is V.
3. the GPU threaded design methods that direction of energy Jacobi battle array according to claim 1 calculates, it is characterised in that：
The k threads that the calculation of step (4) the interior joint injecting power is kernel function S calculate k nodes active and
Reactive power is injected, and calculation formula is：
<mrow>
<msub>
<mi>P</mi>
<mi>k</mi>
</msub>
<mo>=</mo>
<msub>
<mi>V</mi>
<mi>k</mi>
</msub>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>l</mi>
<mo>=</mo>
<msub>
<mi>k</mi>
<mi>b</mi>
</msub>
</mrow>
<msub>
<mi>k</mi>
<mi>e</mi>
</msub>
</munderover>
<msub>
<mi>V</mi>
<mi>j</mi>
</msub>
<mo>&lsqb;</mo>
<mi>G</mi>
<mi> </mi>
<mi>c</mi>
<mi>o</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>&theta;</mi>
<mi>k</mi>
</msub>
<mo>-</mo>
<msub>
<mi>&theta;</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>+</mo>
<mi>B</mi>
<mi> </mi>
<mi>s</mi>
<mi>i</mi>
<mi>n</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>&theta;</mi>
<mi>k</mi>
</msub>
<mo>-</mo>
<msub>
<mi>&theta;</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
</mrow>
<mrow>
<msub>
<mi>Q</mi>
<mi>k</mi>
</msub>
<mo>=</mo>
<msub>
<mi>V</mi>
<mi>k</mi>
</msub>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>l</mi>
<mo>=</mo>
<msub>
<mi>k</mi>
<mi>b</mi>
</msub>
</mrow>
<msub>
<mi>k</mi>
<mi>e</mi>
</msub>
</munderover>
<msub>
<mi>V</mi>
<mi>j</mi>
</msub>
<mo>&lsqb;</mo>
<mi>G</mi>
<mi> </mi>
<mi>s</mi>
<mi>i</mi>
<mi>n</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>&theta;</mi>
<mi>k</mi>
</msub>
<mo>-</mo>
<msub>
<mi>&theta;</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mi>B</mi>
<mi> </mi>
<mi>c</mi>
<mi>o</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>&theta;</mi>
<mi>k</mi>
</msub>
<mo>-</mo>
<msub>
<mi>&theta;</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
</mrow>
Wherein,
K is S kernel functions parallel thread number in GPU Programmings；
PkInjected for the active power of k nodes；
QkInjected for the reactive power of k nodes；
kbIndexed for the starting of admittance battle array Y row k nonzero elements, kb=R [k]；
keIndexed for the termination of admittance battle array Y row k nonzero elements, ke=R [k+1] -1；
J, G, B are the data of l-th of position storage in admittance battle array Y；
VkFor the voltage magnitude of k nodes；
VjFor the voltage magnitude of j nodes；
θkFor the voltage phase angle of k nodes；
θjFor the voltage phase angle of j nodes.
4. the GPU threaded design methods that direction of energy Jacobi battle array according to claim 1 calculates, it is characterised in that：
In the step (5), the computational methods of nonzero element are in H submatrixs：
According to the H submatrixs nonzero element and node admittance battle array Y mapping table M obtained in step (2)H2Y, the t of H kernel functions
When number thread calculates, inquiry table MH2YT positions take out Y battle arrays call number kH2Y, take out Y [kH2Y] corresponding data in element
{ i, j, G, B } is calculated, and the calculation formula of t-th of nonzero element of Jacobian matrix neutron matrix H is：
<mrow>
<msub>
<mi>H</mi>
<mi>t</mi>
</msub>
<mo>=</mo>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<mo>-</mo>
<msub>
<mi>V</mi>
<mi>i</mi>
</msub>
<msub>
<mi>V</mi>
<mi>j</mi>
</msub>
<mrow>
<mo>&lsqb;</mo>
<mrow>
<mi>G</mi>
<mi> </mi>
<mi>sin</mi>
<mrow>
<mo>(</mo>
<mrow>
<msub>
<mi>&theta;</mi>
<mi>i</mi>
</msub>
<mo>-</mo>
<msub>
<mi>&theta;</mi>
<mi>j</mi>
</msub>
</mrow>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mi>B</mi>
<mi> </mi>
<mi>cos</mi>
<mrow>
<mo>(</mo>
<mrow>
<msub>
<mi>&theta;</mi>
<mi>i</mi>
</msub>
<mo>-</mo>
<msub>
<mi>&theta;</mi>
<mi>j</mi>
</msub>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>&rsqb;</mo>
</mrow>
<mo>,</mo>
</mrow>
</mtd>
<mtd>
<mrow>
<mi>i</mi>
<mo>&NotEqual;</mo>
<mi>j</mi>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msubsup>
<mi>V</mi>
<mi>i</mi>
<mn>2</mn>
</msubsup>
<mi>B</mi>
<mo>+</mo>
<msub>
<mi>Q</mi>
<mi>i</mi>
</msub>
<mo>,</mo>
</mrow>
</mtd>
<mtd>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mi>j</mi>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
</mrow>
Wherein,
T is H kernel functions parallel thread number in GPU Programmings；
HtFor the nonzero element of t-th of position in H submatrix COO forms；
I, j, G, B are kth in admittance battle array YH2YThe data of individual position storage；
QiFor the injection reactive power of node i；
ViFor the voltage magnitude of i nodes；
VjFor the voltage magnitude of j nodes；
θiFor the voltage phase angle of i nodes；
θjFor the voltage phase angle of j nodes；
According still further to the H submatrixs nonzero element and matrix J a mapping relations M obtained in step (3)H2Ja, inquire about MH2JaT positions
Put and take out HtStorage call number k of the result of calculation in Ja matrixesH2Ja, then by HtResult of calculation deposit Ja [kH2Ja]。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201510809809.3A CN105391057B (en) | 2015-11-20 | 2015-11-20 | A kind of GPU threaded design methods that direction of energy Jacobi battle array calculates |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201510809809.3A CN105391057B (en) | 2015-11-20 | 2015-11-20 | A kind of GPU threaded design methods that direction of energy Jacobi battle array calculates |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN105391057A CN105391057A (en) | 2016-03-09 |
| CN105391057B true CN105391057B (en) | 2017-11-14 |
Family
ID=55423022
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201510809809.3A Expired - Fee Related CN105391057B (en) | 2015-11-20 | 2015-11-20 | A kind of GPU threaded design methods that direction of energy Jacobi battle array calculates |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN105391057B (en) |
Families Citing this family (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106157176B (en) * | 2016-07-26 | 2019-07-12 | 东南大学 | A kind of LU decomposition method for the direction of energy Jacobian matrix that GPU accelerates |
| CN106684863B (en) * | 2016-12-29 | 2018-02-27 | 华中科技大学 | A kind of discrimination method of power distribution network bus admittance matrix |
| CN107368455A (en) * | 2017-06-22 | 2017-11-21 | 东南大学 | Trigonometric equation group back substitution method on the direction of energy that a kind of GPU accelerates |
| CN108879691B (en) * | 2018-06-21 | 2020-09-04 | 清华大学 | Large-scale continuous power flow calculation method and device |
| CN109344361B (en) * | 2018-08-27 | 2022-05-20 | 南昌大学 | Method for quickly forming Jacobian matrix in power system load flow calculation |
| CN111478333B (en) * | 2020-04-14 | 2021-11-30 | 广东电网有限责任公司广州供电局 | Parallel static security analysis method for improving power distribution network recovery after disaster |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101976835A (en) * | 2010-10-11 | 2011-02-16 | 重庆大学 | Parallel computation method for Newton power flow of large-scale electric power system |
| WO2012061674A3 (en) * | 2010-11-04 | 2013-07-04 | Siemens Corporation | Stochastic state estimation for smart grids |
| CN103617150A (en) * | 2013-11-19 | 2014-03-05 | 国家电网公司 | GPU (graphic processing unit) based parallel power flow calculation system and method for large-scale power system |
| CN103793590A (en) * | 2012-11-01 | 2014-05-14 | 同济大学 | GPU-based computation method for quickly solving power flow in distribution networks |
| CN104967121A (en) * | 2015-07-13 | 2015-10-07 | 中国电力科学研究院 | Large-scale electric power system node load flow computing method |
-
2015
- 2015-11-20 CN CN201510809809.3A patent/CN105391057B/en not_active Expired - Fee Related
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101976835A (en) * | 2010-10-11 | 2011-02-16 | 重庆大学 | Parallel computation method for Newton power flow of large-scale electric power system |
| WO2012061674A3 (en) * | 2010-11-04 | 2013-07-04 | Siemens Corporation | Stochastic state estimation for smart grids |
| CN103793590A (en) * | 2012-11-01 | 2014-05-14 | 同济大学 | GPU-based computation method for quickly solving power flow in distribution networks |
| CN103617150A (en) * | 2013-11-19 | 2014-03-05 | 国家电网公司 | GPU (graphic processing unit) based parallel power flow calculation system and method for large-scale power system |
| CN104967121A (en) * | 2015-07-13 | 2015-10-07 | 中国电力科学研究院 | Large-scale electric power system node load flow computing method |
Non-Patent Citations (3)
| Title |
|---|
| GPU-based power flow analysis with Chebyshev preconditioner andconjugate gradient method;Xue Li，Fangxing Li;《Electric Power Systems Research》;20140614;全文 * |
| 基于GPU的电力系统并行潮流计算的实现;夏俊峰 等;《电力系统保护与控制》;20100916;第38卷(第18期);全文 * |
| 基于道路树分层的大电网潮流并行算法及其GPU优化实现;陈德扬;《电力系统自动化》;20141125;第38卷(第22期);全文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN105391057A (en) | 2016-03-09 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN105391057B (en) | A kind of GPU threaded design methods that direction of energy Jacobi battle array calculates | |
| CN106874113A (en) | A kind of many GPU heterogeneous schemas static security analysis computational methods of CPU+ | |
| CN106157176B (en) | A kind of LU decomposition method for the direction of energy Jacobian matrix that GPU accelerates | |
| CN105576648A (en) | Static security analysis double-layer parallel method based on GPU-CUP heterogeneous computing platform | |
| CN102842917B (en) | A kind of general parallel net type photovoltaic generating system machine-electricity transient model | |
| CN104967121B (en) | A kind of tidal current computing method of large-scale electrical power system node | |
| CN103455948B (en) | A kind of distribution system multi-dimensional multi-resolution Modeling and the method for analysis | |
| CN103337861A (en) | Power distribution network reactive power optimization method based on gold chaotic ecological niche particle swarm algorithm | |
| CN106026107B (en) | A kind of QR decomposition method for the direction of energy Jacobian matrix that GPU accelerates | |
| CN106407158A (en) | GPU accelerated method for performing batch processing of isomorphic sparse matrixes multiplied by full vectors | |
| CN106058863A (en) | Random optimal trend calculation method based on random response surface method | |
| CN101976838B (en) | Newton-process power flow calculation method for study purpose | |
| CN105680473A (en) | Physical fusion modeling method for general electromechanical transient information of photovoltaic power generation system | |
| CN103279824A (en) | Modeling method for relay protection setting calculation system | |
| CN108258725B (en) | Doubly-fed wind turbine dynamic equivalence method based on equivalent power angle coherence | |
| CN106354479B (en) | A kind of GPU acceleration QR decomposition method of a large amount of isomorphism sparse matrixes | |
| CN104734148A (en) | Three-phrase power-distributing network continuation power flow analysis of distributed power supply | |
| CN105896558B (en) | VSC-based UPFC electromechanical transient modular modeling method | |
| CN106294022B (en) | A kind of Jacobian matrix redundancy storage method for static security analysis | |
| CN107368454A (en) | A kind of GPU of the sparse lower trigonometric equation group of a large amount of isomorphisms pushes away method before accelerating | |
| CN106410811B (en) | Iteration small impedance branches endpoint changes the tidal current computing method of Jacobian matrix for the first time | |
| CN109698516A (en) | The maximum capacity computing system and method for renewable energy access power distribution network | |
| Mao et al. | Improved fast short-term wind power prediction model based on superposition of predicted error | |
| Bai et al. | Heuristic optimization for wind energy integrated optimal power flow | |
| CN104182909A (en) | Multi-core parallel successive approximation method of hydropower system optimal scheduling |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20171114 Termination date: 20181120 |
|
| CF01 | Termination of patent right due to non-payment of annual fee |