[发明专利]一种基于双重并行计算的在线Prony分析方法

摘要

本发明公开了一种基于双重并行计算的在线Prony分析方法，涉及电力系统调度自动化领域和计算机高性能计算领域。该方法针对大电网发生低频振荡情况下，振荡参数辨识采用传统基于串行Prony算法资源利用率低、计算速度慢的弊端，提出了多计算节点的分布式并行，有效地进行任务调度和负载均衡，极大地降低系统的响应时间。首次实现了Prony数学模型上的并行化设计，采用多线程并行计算技术实现了Prony的多线程并行计算。本发明能够同步实现对电网多支路、多个电气量的振荡幅值、频率、初相位和衰减因子等参数的在线辨识，有效提高了计算分析速度，能够更好的适应大电网同步在线计算的要求。

主权项

一种基于双重并行计算的在线Prony分析方法，其特征在于，该方法包含下列步骤：(a)在由I条支路和J个节点组成的交流互联电网中，各支路和所有节点均装设同步测量单元PMU；计算平台由一个客户端、一个管理节点和C个计算节点组成的机群构成，其中C_i表示第i个计算节点，P_i为C_i的核数，C_ij表示C_i中的第j个核，i＝1，2，…，C，j＝1，2，3，…，P_i；(b)当电网中发生低频振荡时，将PMU采集的M条支路E个电气量包括：有功功率ΔP、频率Δω，取同一时间段采样数据为一组发送给计算平台进行Prony分析，则需要分析的任务数为N＝EM，M＝1，2，…，I；每一类电气量中包含Len个采样数据，500≤Len≤3000；(c)多机并行计算过程，首先，进行计算节点能力测试，即客户端将接收到的N个任务发送给管理节点，管理节点从N个任务中任意选取一个任务，同时发送给C_i节点上进行计算，i＝1，2，…，C，记录每个节点对该任务的执行时间，节点C_i的执行时间记为T_i，根据第i节点的执行时间与C个节点总执行时间的比例对各计算节点的计算能力进行判别，根据计算能力的强弱对各计算节点进行排序，作为任务分配的依据，然后，将N个计算任务分配到各计算节点开展多机并行计算，节点C_i分配到的任务数为N_i：$T=\underset{i=1}{\overset{C}{\mathrm{\Σ}}}{T}_i$$N_i=\left[\frac{T_i}{T}N\right],i=\mathrm{1,2},...,C$将N_i个任务发送到计算节点C_i上计算，未能整除余下的任务根据计算能力排序依次分配到各计算节点，由C_i进行单节点多线程并行计算；(d)单节点多核并行计算过程：单节点多线程的并行化采用OpenMP的隐式任务并行和显式任务并行技术实现；当节点C_i接收到的任务数为N_i，每个任务的数据长度为Len，第m个采样数据为X_m，m＝1，2，…，Len；Prony计算过程中取迭代次数k＝Len/2；此步骤包含以下计算过程：(d1)最小二乘计算元素求解$R_{u,v}=\underset{m=k}{\overset{m=N-1}{\mathrm{\Σ}}}{X}_{m-v}{X}_{m-u},u,v=\mathrm{1,2},......,k$$r_u=-\underset{m=k}{\overset{m=N-1}{\mathrm{\Σ}}}{X}_m{X}_{m-v}$采用隐式任务并行方式对R_u,v、r_u进行并行计算，采用静态调度方式进行子任务的分配与调度，将R_u,v、r_u平均划分为P_i个子任务，每个子任务负责[k/P_i]×k个数据；由OpenMP将划分得到的P_i个子任务自动分配到核C_ij上执行，i＝1,2,…，C；j＝1，2，3，…，P_i；并行计算后得到R_u,v和r_u，u＝1,2,…，k；v＝1，2,…，k；(d2)采用高斯约旦消元法求解线性方程组：$\left[\begin{array}{ccccc}{R}_{\mathrm{1,1}}& R_{\mathrm{1,2}}& R_{\mathrm{1,3}}& ...& R_{1,k}\\ R_{\mathrm{2,1}}& R_{\mathrm{2,2}}& R_{\mathrm{2,3}}& ...& R_{2,k}\\ R_{\mathrm{3,1}}& R_{\mathrm{3,2}}& R_{\mathrm{3,3}}& ...& R_{3,k}\\ ...& ...& ...& ...& ...\\ R_{k,1}& R_{k,2}& R_{k,3}& ...& R_{k,k}\end{array}\right]\·\left[\begin{array}{c}{\mathrm{\α}}_1\\ {\mathrm{\α}}_2\\ {\mathrm{\α}}_3\\ ...\\ {\mathrm{\α}}_k\end{array}\right]=\left[\begin{array}{c}{r}_1\\ r_2\\ r_3\\ ...\\ r_k\end{array}\right]$得到α＝(α₁、α₂、...、α_k)^T采用隐式并行策略实现高斯约旦消元法求解线性方程组的并行化；高斯约旦消元法要求进行k次消元，每次消元过程需要对k个方程进行消元，将k个方程的消元过程分配到P_i个核上并行执行，每个核分配[k/P_i]个方程的消元；划分时，使任务划分的数据间隔为一个高速缓存行的整数倍大小，防止伪共享现象，并在调度上采用动态调度方式。划分后的子任务将由OpenMP自动调度到核C_i1、C_i2、…、C_iPi上计算；计算后得到α₁、α₂、…、α_k；(d3)将α₁、α₂、…、α_k代入式X(n)＝‑[α₁X(n‑1)+α₂X(n‑2)+...+α_kX(n‑k)]求解高次方程z^k+α₁z^k‑1+α₂z^k‑2+...+α_k＝0即求解如下矩阵的特征值：$\left[\begin{array}{ccccc}-{\mathrm{\α}}_1& -{\mathrm{\α}}_2& ...& -{\mathrm{\α}}_{k-1}& -{\mathrm{\α}}_k\\ 1& 0& ...& 0& 0\\ ...& ...& ...& ...& ...\\ 0& 0& 1& 0& 0\\ 0& 0& 0& 1& 0\end{array}\right]$采用显式任务并行方式实现矩阵特征值的并行求解；求解过程需要k次循环，每次循环过程分为两个部分的计算，一是与其他任务有数据相关的部分，即对矩阵元素更新部分，二是无数据相关的部分，即计算出特征值部分；主线程只负责执行有数据相关的部分，无数据相关的部分利用task指令显式构造并行任务，添加到任务队列中，在OpenMP的调度下由其他线程执行。核Ci1执行主线程，在其对元素更新之后，将其余计算部分显式构造成子任务添加到任务队列中，进而由空闲的核取子任务并执行；最后显式汇总求得：z₁、z₂、……、z_k。(d4)求解线性方程组$\left[\begin{array}{ccccc}1& 1& 1& ...& 1\\ z_1& z_2& z_3& ...& z_k\\ z_1^2& z_2^2& z_3^2& ...& z_k^2\\ ...& ...& ...& ...& ...\\ z_1^{N-1}& z_2^{N-1}& z_3^{N-1}& ...& z_k^{N-1}\end{array}\right]\·\left[\begin{array}{c}{b}_1\\ b_2\\ b_3\\ ...\\ b_k\end{array}\right]=\left[\begin{array}{c}{x}_1\\ x_2\\ x_3\\ ...\\ x_{N-1}\end{array}\right]$采用隐式任务并行方式完成并行计算，将转置矩阵分为P_i个较小的矩阵块，每个矩阵块大小为(k/P_i)×N，这样就产生P_i个矩阵相乘的子任务，指定调度方式为静态调度，P_i个子任务由OpenMP自动调度到核上完成并行计算；并行计算后得到k个方程，用高斯约旦消元法和隐式任务并行方式进行线性方程组的并行求解，指定调度方式为动态调度；并行计算得到b₁、b₂、…、b_k；(d5)将以上求出的z_u和b_u代入下式$\left\{\begin{array}{c}{A}_u=\left|b_u\right|\\ {\mathrm{\θ}}_u=\arctan (\mathrm{lm}\left(b_u\right)/\mathrm{Re}\left(b_u\right))/2\mathrm{\π\Δt}\\ {\mathrm{\α}}_u=\mathrm{lm}\left(z_u\right)/\mathrm{\Δt}\\ f_u=\arctan (\mathrm{lm}\left(z_u\right)/\mathrm{Re}\left(z_u\right))/2\mathrm{\π\Δt}\end{array}\right.,u=\mathrm{1,2},...k$计算出每个任务的幅值A_u、相位θ_u、频率f_u和衰减因子α_u，u＝1，2，…，k；对分配到的N_i个计算任务重复步骤(d1)至(d5)，完成节点C_i上的N_i个计算任务，得到N_i个任务对应电气量的振幅、频率、初相位、衰减因子的值；(e)在每个节点完成其负责的N_i个计算任务后，将各电气量的振幅、频率、初相位和衰减因子的值发送到管理节点；由管理节点汇总C个计算节点的结果，将结果返回给客户端。