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 PDF

Info

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
Application number
CN201510809809.3A
Other languages
Chinese (zh)
Other versions
CN105391057A (en
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.)
State Grid Corp of China SGCC
Southeast University
China Electric Power Research Institute Co Ltd CEPRI
State Grid Fujian Electric Power Co Ltd
Original Assignee
State Grid Corp of China SGCC
Southeast University
China Electric Power Research Institute Co Ltd CEPRI
State Grid Fujian Electric Power Co Ltd
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 State Grid Corp of China SGCC, Southeast University, China Electric Power Research Institute Co Ltd CEPRI, State Grid Fujian Electric Power Co Ltd filed Critical State Grid Corp of China SGCC
Priority to CN201510809809.3A priority Critical patent/CN105391057B/en
Publication of CN105391057A publication Critical patent/CN105391057A/en
Application granted granted Critical
Publication of CN105391057B publication Critical patent/CN105391057B/en
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • HELECTRICITY
    • H02GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02JCIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J3/00Circuit arrangements for ac mains or ac distribution networks
    • HELECTRICITY
    • H02GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02JCIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J2203/00Indexing scheme relating to details of circuit arrangements for AC mains or AC distribution networks
    • H02J2203/20Simulating, 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

A kind of GPU threaded design methods that direction of energy Jacobi battle array calculates
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>&amp;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>&amp;lsqb;</mo> <mi>G</mi> <mi> </mi> <mi>c</mi> <mi>o</mi> <mi>s</mi> <mrow> <mo>(</mo> <msub> <mi>&amp;theta;</mi> <mi>k</mi> </msub> <mo>-</mo> <msub> <mi>&amp;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>&amp;theta;</mi> <mi>k</mi> </msub> <mo>-</mo> <msub> <mi>&amp;theta;</mi> <mi>j</mi> </msub> <mo>)</mo> </mrow> <mo>&amp;rsqb;</mo> </mrow>
<mrow> <msub> <mi>Q</mi> <mi>k</mi> </msub> <mo>=</mo> <msub> <mi>V</mi> <mi>k</mi> </msub> <munderover> <mo>&amp;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>&amp;lsqb;</mo> <mi>G</mi> <mi> </mi> <mi>s</mi> <mi>i</mi> <mi>n</mi> <mrow> <mo>(</mo> <msub> <mi>&amp;theta;</mi> <mi>k</mi> </msub> <mo>-</mo> <msub> <mi>&amp;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>&amp;theta;</mi> <mi>k</mi> </msub> <mo>-</mo> <msub> <mi>&amp;theta;</mi> <mi>j</mi> </msub> <mo>)</mo> </mrow> <mo>&amp;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>&amp;lsqb;</mo> <mrow> <mi>G</mi> <mi> </mi> <mi>sin</mi> <mrow> <mo>(</mo> <mrow> <msub> <mi>&amp;theta;</mi> <mi>i</mi> </msub> <mo>-</mo> <msub> <mi>&amp;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>&amp;theta;</mi> <mi>i</mi> </msub> <mo>-</mo> <msub> <mi>&amp;theta;</mi> <mi>j</mi> </msub> </mrow> <mo>)</mo> </mrow> </mrow> <mo>&amp;rsqb;</mo> </mrow> <mo>,</mo> </mrow> </mtd> <mtd> <mrow> <mi>i</mi> <mo>&amp;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]。
CN201510809809.3A 2015-11-20 2015-11-20 A kind of GPU threaded design methods that direction of energy Jacobi battle array calculates Expired - Fee Related CN105391057B (en)

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)

* Cited by examiner, † Cited by third party
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)

* Cited by examiner, † Cited by third party
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

Patent Citations (5)

* Cited by examiner, † Cited by third party
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)

* Cited by examiner, † Cited by third party
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