CN111478333A - 一种提升配电网灾后恢复用并行静态安全分析方法 - Google Patents

Abstract

本发明公开了一种提升配电网灾后恢复用并行静态安全分析方法，通过故障情景生成，得到极端灾害情况下的预想故障集，采用GPU进行并行静态安全分析计算，对并行导纳矩阵形成及修改、雅可比矩阵形成、注入功率计算进行并行化设计，对各极端灾害故障情况下的修正方程组求解采用ILU预处理的GMRES算法作为内部迭代的非精确牛顿法并行计算，找到各极端灾害故障情境下的电力系统主要防护元件。本发明能够给出极端灾害故障情境下系统中的重要薄弱环节，根据计算结果，发现能够有效减少重要负荷停电时间的设备位置，发挥快速实时在线分析的优势，尽可能减少维修过程中电力系统的重要负荷的断电时间。

Description

一种提升配电网灾后恢复用并行静态安全分析方法

技术领域

本发明属于电力系统安全规划运行技术领域，具体涉及一种提升配电网灾后恢复用并行静态安全分析方法。

背景技术

电力系统中，为了避免过负荷、电压越限等引起设备损坏，或故障情况下保护设备动作导致的大面积停电，需要快速的在线实时静态安全分析。随着电网规模不断扩大，需要校验的状态数目巨大，传统的串行计算方法已经无法满足大电网实时计算的要求，同时，近年来，大规模电力系统故障往往涉及多个器件故障，此时需要“N-x”故障分析，预想故障集大大增加，导致计算量更大。因此高效的并行化静态安全分析方法对于现代电力系统安全稳定运行非常重要。

此外，近年来发生的极端情况下大电网的事故往往导致配电网瘫痪，失去重要负荷，以及造成大量符合丧失，研究配电网灾害预防和灾后恢复变得越来越有必要。例如2019年在委内瑞拉连续发生了几次大规模停电事故，其原因是网络系统被黑客攻击(网络系统和电力系统耦合)。而恢复力研究涉及灾前、灾中以及灾后过程中采取多种有效措施来减小切负荷量、灾害结束后能够快速恢复到正常运行状态。应对极端事件的能力逐渐成为电力系统必不可少的能力。

通过故障情景生成，得到极端灾害情况下的预想故障集，采用GPU进行并行静态安全分析计算，找到极端灾害故障集下的电力系统薄弱环节，在灾害后作为重点抢修对象，可以使得在灾害发生时能够发挥快速实时在线分析的优势，尽可能减少维修过程中电力系统的重要负荷的断电时间，以及整体上降低负荷的停电时间。

发明内容

本发明所要解决的技术问题在于针对上述现有技术中的不足，提供一种提升配电网灾后恢复用并行静态安全分析方法，能够给出极端灾害预想事故下，系统中的主要防护元件，根据计算结果，揭示有效减少重要负荷停电时间的关键设备，发挥快速实时在线分析的优势，尽可能减少维修过程中电力系统的重要负荷的断电时间。

本发明采用以下技术方案：

一种提升配电网灾后恢复用并行静态安全分析方法，包括以下步骤：

S1、基于致灾因子建模刻画灾害的性质、强度与概率，依照灾变强度的概率分布生成灾害场景；基于脆弱性建模建立电力系统元件故障率与灾变强度之间的条件依赖关系，结合不可修复元件的可靠性理论将灾害场景转化为电力系统元件的故障场景；综合所生成的故障场景，构建极端灾害下的电力系统故障场景集；

S2、根据电力系统结构信息，采用图形处理器并行化构建电力系统导纳矩阵；根据极端灾害故障集，以核函数Kernel在GPU上进行导纳矩阵并行修正；

S3、根据已修正的导纳矩阵，以核函数Kernel在GPU上并行形成雅可比矩阵、并行计算节点的有功无功注入功率，形成功率偏差，得到极端灾害故障集的修正方程组；

S4、采用非精确牛顿法并行求解步骤S3得到的修正方程组，采用预处理的广义极小化残差法对非精确牛顿法的内部迭代，采用不完全LU分解进行预处理；

S5、设置重要元件判据，根据计算结果，在GPU中并行判断各支路功率越限情况和排序、机组状态的统计和排序，重要元件排序后，输出作为极端灾害下的主要防护元件。

具体的，步骤S2中，创建核函数adMatrixKernel以GPU线程级别形成对应电力系统导纳矩阵的互导纳，一个线程形成一个互导纳元素，当电力系统支路数超过1024*240时，开启多个网格；在已形成导纳矩阵的基础上，根据步骤S1中形成的极端灾害下的电力系统故障场景集，开启新的核函数Y_XZ_Kernel进行导纳矩阵的并行修正；每个线程上的导纳矩阵修正对应每个故障情况下元件故障导致的对应连接节点的导纳矩阵修正，当一个元件故障时，修正的导纳元素包括两端节点的自导纳和他们的互导纳共4个元素，当N_element个元件同时故障时，修正的是4*N_element个导纳矩阵元素。

具体的，步骤S3中，k号线程计算k号节点的有功注入功率P_k和无功注入功率Q_k具体为：

其中，θ_k为k号节点的电压相角，θ_j为j号节点的电压相角，G为导纳矩阵实部，V_j为j号节点的电压幅值，V_k为k号节点的电压幅值，B为导纳矩阵虚部，k为核函数Make_S_Kernel在GPU中的线程号；NYseq[k]为导纳矩阵第k行非零元素起始索引；NYseq[k+1]-1为导纳矩阵第k行非零元素终止索引。

具体的，步骤S3中，雅可比子矩阵为：H_n×n、N_n×m、M_m×n、L_m×m，其中，m为PQ节点个数、n为系统总结点数-1。

进一步的，首先在CPU上形成雅可比子矩阵H、N、M、L与导纳矩阵的对应关系H2Y、N2Y、M2Y、L2Y，作为并行形成雅可比子矩阵的下标索引；然后在CPU上计算雅可比子矩阵H、N、M、L与雅可比矩阵的元素对应关系H2Ja、N2Ja、M2Ja、L2Ja；根据与导纳矩阵不同节点类型的对应关系，在GPU上创建核函数Make_H_Kernel、Make_N_Kernel、Make_M_Kernel、Make_L_Kernel并行计算雅可比子矩阵H、N、M、L，H子矩阵的核函数线程数为n×n，N、H子矩阵的核函数线程数为n×m，L子矩阵的核函数线程数为m×m。

具体的，步骤S4中，通过求解极小化J(y)得到y_m向量，确定线性方程组的近似解，通过预处理矩阵M，将线性方程组Ax＝b转化为同解的易于迭代求解的新的方程组M^-1Ax＝M^-1b，选取不完全LU分解(ILU)进行预处理；通过相对残差rr_k+1限制内部GMRES迭代求解算法求解当前步骤解s_k的精度值，通过内部迭代求解得到的相对残差与预设的阈值比较，决定内部迭代是否收敛；当相对残差小于预设的强迫因子η时，内部迭代停止。

进一步的，相对残差rr_k+1具体为：

其中，e_k为误差范数，F(x_k)为，F′(x_k)为雅可比矩阵，s_k为修正量。

具体的，步骤S5中，在CPU上进行重要元件排序，重要元件判据根据相对越限程度进行判断，定义重要元件判据为：

其中，F_ij为故障i情况下支路j的支路功率或故障i情况下节点j有功、无功；

为系统中支路j的支路功率或节点j标准有功、标准无功。

与现有技术相比，本发明至少具有以下有益效果：

本发明一种提升配电网灾后恢复用并行静态安全分析方法，采用并行算法，提高实时在线分析速度非常重要；首先，本发明能够根据输入条件由灾害建模得到发生概率较大的极端灾害故障集，根据GPU的结构特征，在粗粒度并行中进行极端灾害故障集下的静态安全分析，根据CUDA建立每个极端故障情况下的线程级潮流计算，使静态安全分析的计算速度提高倍数。然后，对节点注入功率和每个潮流计算过程中雅可比矩阵的形成进行了细粒度的并行设计，以提高单次潮流的计算速度，达到较好的加速度效果。对于修正方程的求解，由于目前电网和配电网规模较大，电网和配电网较为复杂。此外，修正方程多为稀疏方程，因此采用迭代法求解。由于每次牛顿迭代过程修正方程的求解不需要很精确，这里引入非精确的牛顿方法，以提高内部迭代计算速度。针对网络中存在相移器而导致的Jacobian矩阵不对称的情况，内部迭代方法采用广义最小化残差法(GMRES)，通过不完全的LU预处理进一步加速内部迭代收敛，进一步提高细粒并行性的效率。

进一步的，在进行了极端灾害故障集潮流的粗粒度并行后，需要提高单次潮流计算速度，单次潮流计算中，修正方程的求解消耗时间占比约60％-70％，雅可比矩阵的形成消耗时间占比约30％左右，此外剩余的时间消耗主要在于导纳矩阵的形成等准备步骤，在这些步骤中，导纳矩阵的形成具有很好的并行性，因此导纳矩阵的并行形成步骤对于提升单次潮流计算速度有一定的积极作用；在此基础上，对于每个极端灾害故障情况，如果都进行一次导纳矩阵形成将耗费巨大的时间，但是采用并行修正导纳矩阵的形式将大大节省这部分时间。

进一步的，较常规潮流计算的顺序求解雅可比矩阵各元素，耗费的时间占潮流计算的30％左右，因此将雅可比矩阵的形成过程并行化设计，将使得这部分时间大大节省，每次计算消耗的时间是计算单个雅可比矩阵元素的四倍左右，注入功率的并行计算同样可以进一步降低单次潮流求解的时间。之后形成的基于线程计算的修正方程则是静态安全分析粗粒度并行的结果，可以大大节省计算时间。

进一步的，占据潮流求解时间70％的修正方程组求解需要采用更快速的计算方法，尽管LU分解已经可以达到很快的速度，但是对于系统规模更大的情况求解耗时仍然很大，并且难以实现并行化，还需要寻找更高效的求解算法。对于大规模稀疏线性方程组，非精确牛顿法是一种改进的牛顿法，在求解大规模稀疏线性方程组时十分有效，且由于其自身适合与迭代法配合，因此可以应用于GPU并行化计算，大大降低计算时间；这里采用的广义极小化残差法(GMRES)适用于系统中存在移相器等情况导致的雅可比矩阵非对称正定，不完全LU分解法作为预处理进一步加快内部迭代收敛速度，常数强迫因子的选取降低了对修正方程求解精度的限制，可以避免直接法出现的过解现象。

进一步的，设置的重要元件判据原则是相对与元件自身的越限百分比，这样的选取原则可以避免元件自身标准值过大导致的绝对偏差过大的情况，实现元件排序的合理性；在GPU中并行判断各故障下的越限情况，可以有效降低大量数据在GPU和CPU之间传输过程的时间消耗，同时也降低了在CPU中串行判断和排序消耗的大量时间。同时，这种判据和排序的方式保证了重要元件作为薄弱环节的合理性。

综上所述，本发明能够给出极端灾害故障情境下系统中的重要薄弱环节，根据计算结果，发现能够有效减少重要负荷停电时间的设备位置，发挥快速实时在线分析的优势，尽可能减少维修过程中电力系统的重要负荷的断电时间。

下面通过附图和实施例，对本发明的技术方案做进一步的详细描述。

附图说明

图1为整体算法流程图；

图2为GPU在CUDA框架下的内部结构图；

图3为非精确牛顿法迭代求解流程图；

图4为预处理广义极小化残差法并行化计算流程设计图。

具体实施方式

本发明提供了一种提升配电网灾后恢复用并行静态安全分析方法，通过故障情景生成，得到极端灾害情况下的预想故障集，采用GPU进行并行静态安全分析计算，对并行导纳矩阵形成及修改、雅可比矩阵形成、注入功率计算进行并行化设计，对各极端灾害故障情况下的修正方程组求解采用ILU预处理的GMRES算法作为内部迭代的非精确牛顿法并行计算，找到各极端灾害故障情境下的电力系统主要防护元件。

请参阅图1，本发明一种提升配电网灾后恢复用并行静态安全分析方法，包括以下步骤：

S1、基于致灾因子建模刻画灾害的性质、强度与概率，依照灾变强度的概率分布生成灾害场景；基于脆弱性建模建立电力系统元件故障率与灾变强度之间的条件依赖关系，结合不可修复元件的可靠性理论将灾害场景转化为电力系统元件的故障场景；综合所生成的故障场景，构建极端灾害下的电力系统故障场景集；

S2、根据电力系统结构信息，利用NVIDIA推出的通用并行计算架构(CUDA，ComputeUnified Device Architecture)，采用图形处理器(GPU)并行化构建电力系统导纳矩阵；根据极端灾害故障集，以核函数Kernel在GPU上进行导纳矩阵并行修正；

请参阅图2，为GPU的内部结构示意图，根据NVIDIA推出的通用并行计算架构CUDA对GPU内部结构的定义，GPU设备内部软件结构包括众多的网格(grid)，网格(grid)的下一级结构是线程块(Thread Block)，线程块(Thread Block)的下一级结构为线程(Thread)。以往的串行计算过程，都是按照顺序计算分别计算每种故障情况下潮流，但是这种静态安全分析存在着很大的粗粒度并行性，按照每个线程(Thread)求解一种故障情况下的潮流，这样整体的计算速度可以达到N_fault倍的提升，这里的N_fault为极端灾害故障集的故障个数。这里可以通过CPU调用GPU并行执行的方式，开启并行计算，采用核函数Kernel实现这一功能。

以C语言联合CUDA实现本发明的功能。

首先创建核函数“adMatrixKernel”以GPU线程级别形成对应电力系统导纳矩阵的互导纳，一个线程形成一个互导纳元素，由于通常GPU一个网格(grid)包含240个线程块(Thread Block)，一个线程块(Thread Block)最大为1024个线程，因此当电力系统支路数超过1024个时，需要开启多个线程块，当电力系统支路数超过1024*240时，需要开启多个网格。互导纳的线程并行实现以支路编号为索引对应的线程编号threadIdx.x，编程实现第k个线程计算第k条支路的互导纳；

在已经形成导纳矩阵的基础上，根据步骤S1中形成的极端灾害预想故障集，开启新的核函数“Y_XZ_Kernel”进行导纳矩阵的并行修正。这里新的核函数个数与极端灾害预想故障集中的故障个数N_fault相同，同时也是接下来要创建的非精确牛顿法并行求解修正方程组的线程个数，也就是实现极端灾害故障的粗粒度并行；

每个线程上的导纳矩阵修正对应着每个故障情况下的元件故障导致的对应连接节点的导纳矩阵修正，当一个元件故障时，需要修正的导纳元素包括两端节点的自导纳和他们的互导纳共4个元素，当N_element个元件同时故障时，需要修正的是4*N_element个导纳矩阵元素。同时，在每种故障情况下，需要由CPU建立故障元件连接节点与线程编号的对应关系，方便并行修正导纳矩阵元素。

S3、根据已修正的导纳矩阵，以核函数Kernel在GPU上并行形成雅可比矩阵、并行计算节点注入功率，形成功率偏差，得到极端灾害故障集的修正方程组；

节点注入功率的计算，采用CUDA创建一个核函数“Make_S_Kernel”，根据下式进行计算，k号线程计算k号节点的有功无功注入功率：

其中，k为核函数“Make_S_Kernel”在GPU中的线程号；NYseq[k]为导纳矩阵第k行非零元素起始索引；NYseq[k+1]-1为导纳矩阵第k行非零元素终止索引。

雅可比矩阵的并行设计时，需要考虑节点类型对于导纳矩阵与雅可比矩阵中的元素索引的影响，雅可比子矩阵：

H_n×n、N_n×m、M_m×n、L_m×m

其中，m为PQ节点个数、n为系统总结点数-1。

首先在CPU上形成雅可比子矩阵H、N、M、L与导纳矩阵的对应关系“H2Y、N2Y、M2Y、L2Y”，将其作为并行形成雅可比子矩阵的下标索引；

然后在CPU上计算雅可比子矩阵H、N、M、L与雅可比矩阵的元素对应关系“H2Ja、N2Ja、M2Ja、L2Ja”；

接下来，根据与导纳矩阵不同节点类型的对应关系，开始在GPU上创建核函数“Make_H_Kernel”、“Make_N_Kernel”、“Make_M_Kernel”、“Make_L_Kernel”并行计算雅可比子矩阵H、N、M、L，其中需要注意形成H子矩阵的核函数线程数为n×n，N、H子矩阵的核函数线程数为n×m，L子矩阵的核函数线程数为m×m。

具体程序中，每个线程计算雅可比子矩阵一个子元素，并根据与雅可比矩阵的对应关系存入稀疏雅可比矩阵对应的位置。

根据雅可比矩阵和节点注入功率，并行形成极端故障情况下的修正方程组。

S4、采用非精确牛顿法并行求解修正方程组，其中非精确牛顿法的内部迭代采用预处理的广义极小化残差法(GMRES)，预处理的方法采用不完全LU分解(ILU)；

在步骤S3中得到基于线程的各故障情况的修正方程，且均为大规模稀疏线性方程组，非精确牛顿法是一种改进的牛顿法，在求解大规模稀疏线性方程组时十分有效，且由于其自身适合与迭代法配合，因此可以应用于GPU并行化计算，大大降低计算时间；同时也降低了对修正方程求解精度的限制，可以避免直接法出现的过解现象。

下式是非精确牛顿法的基本原理：

F′(x_k)s_k＝-F(x_k)+e_k,||e_k||/||F(x_k)||＜η (3)

x_k+1＝x_k+s_k (4)

其中，F′(x_k)为雅可比矩阵；s_k为修正量；e_k为误差范数；η为强迫因子。

相比较于牛顿法，非精确牛顿法在求解线性方程组的时候加入了一定的偏差量e_k，因此在求解该方程组时，得到的线性解s_k并非完全精确的。

通过近似求解修正方程，可以大大降低计算量，提高计算速度。非精确牛顿法是通过在每次迭代中引入偏差以减少每次迭代的计算量，达到增加计算效率的目的，但是带来的负面效果是增加了外部迭代的次数。

因此，非精确牛顿法的关键在于选取合适的内部迭代算法进行非精确求解，以及选取合适的收敛准则，使得正面效果大于负面效果。非精确牛顿法的求解流程如图3所示；外部循环是牛顿循环，用来求解非线性潮流方程；内部迭代算法求解线性方程组。

由于实际电力系统中移相器的存在等，导致雅可比矩阵非对称正定，因此采用广义极小化残差算法(GMRES)作为内部迭代算法批处理求解修正方程组，本方法仅要求系数矩阵非奇异。GMRES算法是Krylov子空间迭代算法的典型代表，是一种斜投影方法。

Krylov-GMRES子空间法就是基于Krylov子空间上，选择约束子空间

使得

的斜投影方法，即搜索子空间κ为Krylov子空间：

κ_m(A,v)＝span{v,Av,A²v,…,A^m-1v} (5)

其中，v初始残差v＝r₀＝b-Ax₀；A线性方程组系数矩阵，即修正方程式的系数矩阵F′(x_k)，也就是雅克比矩阵。

Krylov-GMRES基本原理如下：

对于空间M中的任意向量x可以以基向量的形式表示出来，如下式

x＝x₀+V_my (6)

其中，y是长度和系数矩阵规模相同的向量。

定义J(y)如下

J(y)＝||b-Ax||₂＝||b-A(x₀+V_my)||₂ (7)

由正交基的Arnoldi的构造过程，可以做如下推导:

由于V_m+1列向量是正交的，因此

J(y)＝||b-A(x₀+V_my)||₂＝||βe₁-H_my||₂ (9)

通过求解极小化J(y)，可以得到y_m向量，这样线性方程组的近似解可以表示为：

x_m＝x₀+V_my_m (10)

为加快内部迭代收敛速度，提高非精确求解的效率，需要采取合适的预处理手段对系数矩阵进行预处理，改善系数矩阵的条件数。

系数矩阵预处理的核心思想是，通过预处理矩阵M，将线性方程组Ax＝b转化为同解的易于迭代求解的新的方程组M^-1Ax＝M^-1b，这里选取不完全LU分解(ILU)进行预处理，以提高内部迭代的求解和收敛速度。

ILU预处理GMRES算法流程如下：

在选取内部迭代算法的基础上，选取合适的内部迭代收敛准则对于计算效果有很大影响。下式为相对残差，在程序设计时，相对残差rr_k+1限制内部GMRES迭代求解算法求解当前步骤解s_k的精度值；通过内部迭代求解得到的相对残差与预设的阈值比较，以此决定内部迭代是否收敛。

这个阈值即为强迫因子，当相对残差小于预设的强迫因子η时，内部迭代停止。

强迫因子η本质上刻画了内部牛顿方程求解的精确程度，它影响着非精确牛顿法的收敛性以及计算速度。图4为内部迭代算法并行化设计示意图，其中预处理矩阵在CPU上得到，雅可比矩阵和偏差量在GPU中并行计算得到，修正方程整体求解过程在GPU上并行计算。

S5、根据潮流计算结果，在GPU中并行判断各支路功率越限情况和排序、机组状态的统计和排序。

设置重要元件判据，根据越限情况将元件排序，并将得到的重要元件和越限数据传输到CPU显示，根据输出结果作为极端灾害下的主要防护元件；

由于从故障情景生成开始，主体计算保持在GPU中进行，因此计算的结果需要输出到CPU中，由于故障情景很多，同时不仅需要输出该种故障模式下运行方式的越限情况，还要输出该种故障模式下，具体线路发电机的越限情况，导致需要传输的数据量过大、GPU存储资源占用过多。因此此处的越限判断设置在GPU内部执行，并行越限判据包括支路功率越限、和发电机无功越限。在图1潮流计算收敛或终止后，直接进行越限判据，在判断过程中将未越限的故障情况和未越限的节点直接舍弃，仅存储越限的故障情况形成越限故障集

以及某故障情况下的越限节点集τ_node。将得到的越限故障集

和越限节点集τ_node及相关数据传输到CPU上。

在CPU上进行重要元件排序，重要元件判据根据相对越限程度进行判断，定义重要元件判据为：

其中，F_ij为故障i情况下支路j的支路功率或故障i情况下节点j有功、无功；

为系统中支路j的支路功率或节点j标准有功、标准无功。重要元件排序后，输出作为极端灾害下的主要防护元件。

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。通常在此处附图中的描述和所示的本发明实施例的组件可以通过各种不同的配置来布置和设计。因此，以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围，而是仅仅表示本发明的选定实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

计算实验设置：

计算实验在一台配备NVIDIA GeForce GTX 1660Ti GPU上进行，该GPU有24个流多处理器，每个流处理器有64个CUDA核、共计1536个CUDA核心，计算能力7.5，线程束(wrap)大小为32，每个流处理器中的最大线程数1024个，每个线程块(block)上面的最大线程数1024。

服务器拥有6核Intel(R)Core(TM)i7-9750H CPU@2.60GHz和16GB内存。

CUDA驱动程序版本是10.2，Visual Studio的版本是2017版。

操作系统是64位Windows 10操作系统。

为了保持与Matlab相同的精度，CUDA的所有计算都是基于双精度浮点运算。强迫因子选取常数强迫因子0.001，迭代停止准则是相对残差(rr)小于0.001。

算例系统选择MATPOWER中的case 2383wp是波兰系统示例之一，共计2383个节点、2896条支路，展示几组本地连接的子系统。首先通过转换程序，将MATPOWER中的case2383wp.m文件转换为C语言输入格式数据，其次将下表中的断线支路设置为极端灾害故障，共设置4种极端故障情况。

表1极端灾害故障集

根据四种极端灾害的情况，创建4条线程并行执行潮流计算，导纳矩阵的形成中互导纳需要创建2896个线程；每次潮流计算中雅可比矩阵子阵和注入功率的并行计算最多需要创建2382个线程；每次潮流计算修正方程的求解最多需要创建的线程不超过2*2383个。而本算例的计算实验平台的最大可创建线程数为24*1024个，因此硬件配置可以完成本算例。

表2导纳矩阵计算测试结果

算例	CPU(ms)	GPU(ms)	加速比
				case 2383wp	5.30	0.95	5.58

表3雅可比矩阵计算测试结果

算例	CPU(ms)	GPU(ms)	加速比
				case 2383wp	201.74	24.84	8.12

表4静态安全分析计算测试结果

算例	CPU(ms)	GPU(ms)	加速比
				case 2383wp	706.30	220.03	3.21

根据表2的测试结果，可以看出导纳矩阵并行计算后可以达到5.58倍的加速比，说明并行化导纳矩阵可以达到很好的效果。但是对比表4整体的计算时间消耗，导纳矩阵计算消耗的时间很少，几乎可以忽略不计。根据表3的测试结果，可以看出并行雅可比矩阵的计算可以实现很好的加速效果，加速比达到8.12倍，同时可以看到雅可比矩阵的计算耗时占总计算耗时很大一部分比例，因此很有必要进行并行化设计。

根据表4的测试结果，可以看出，针对case 2383wp算例系统，本方法相对与常规的静态安全分析能够实现明显的加速效果，加速比达到3.21。理论上的加速比可以达到4左右，这里的原因可能是由于CPU与GPU之间频繁的数据传输导致了额外的时间损耗，因此难以达到理想的加速效果，但是目前的加速比已经很大了，足以说明本发明的效果。

本发明一种提升配电网灾后恢复用并行静态安全分析方法，具有以下特点：

1)依照灾变强度的概率分布生成灾害场景，基于脆弱性建模建立电力系统元件故障率与灾变强度之间的条件依赖关系，结合不可修复元件的可靠性理论将灾害场景转化为电力系统元件的故障场景，构建极端灾害下的电力系统故障场景集；

2)采用图形处理器(GPU)并行化构建电力系统导纳矩阵；根据极端灾害故障集，以核函数Kernel在GPU上进行导纳矩阵并行修正；

3)根据已修正的导纳矩阵，以核函数Kernel在GPU上并行形成雅可比矩阵、并行计算节点注入功率，形成功率偏差，得到极端灾害故障集的修正方程组；

4)采用ILU预处理的GMRES算法非精确牛顿法并行求解修正方程组；

5)采取重要元件判据判断，在GPU中并行判断各支路功率越限情况和排序、机组状态的统计和排序，输出到CPU中，根据输出结果作为极端灾害下的主要防护元件。

综上所述，本方法可以给出极端灾害故障情境下，系统中的重要薄弱环节，根据计算结果，发现能够有效减少重要负荷停电时间的设备位置，发挥快速实时在线分析的优势，尽可能减少维修过程中电力系统的重要负荷的断电时间。

以上内容仅为说明本发明的技术思想，不能以此限定本发明的保护范围，凡是按照本发明提出的技术思想，在技术方案基础上所做的任何改动，均落入本发明权利要求书的保护范围之内。

Claims (8)

1.一种提升配电网灾后恢复用并行静态安全分析方法，其特征在于，包括以下步骤：

S1、基于致灾因子建模刻画灾害的性质、强度与概率，依照灾变强度的概率分布生成灾害场景；基于脆弱性建模建立电力系统元件故障率与灾变强度之间的条件依赖关系，结合不可修复元件的可靠性理论将灾害场景转化为电力系统元件的故障场景；综合所生成的故障场景，构建极端灾害下的电力系统故障场景集；

S2、根据电力系统结构信息，采用图形处理器并行化构建电力系统导纳矩阵；根据极端灾害故障集，以核函数Kernel在GPU上进行导纳矩阵并行修正；

S3、根据已修正的导纳矩阵，以核函数Kernel在GPU上并行形成雅可比矩阵、并行计算节点的有功无功注入功率，形成功率偏差，得到极端灾害故障集的修正方程组；

S4、采用非精确牛顿法并行求解步骤S3得到的修正方程组，采用预处理的广义极小化残差法对非精确牛顿法的内部迭代，采用不完全LU分解进行预处理；

S5、设置重要元件判据，根据计算结果，在GPU中并行判断各支路功率越限情况和排序、机组状态的统计和排序，重要元件排序后，输出作为极端灾害下的主要防护元件。

2.根据权利要求1所述的提升配电网灾后恢复用并行静态安全分析方法，其特征在于，步骤S2中，创建核函数adMatrixKernel以GPU线程级别形成对应电力系统导纳矩阵的互导纳，一个线程形成一个互导纳元素，当电力系统支路数超过1024*240时，开启多个网格；在已形成导纳矩阵的基础上，根据步骤S1中形成的极端灾害下的电力系统故障场景集，开启新的核函数Y_XZ_Kernel进行导纳矩阵的并行修正；每个线程上的导纳矩阵修正对应每个故障情况下元件故障导致的对应连接节点的导纳矩阵修正，当一个元件故障时，修正的导纳元素包括两端节点的自导纳和他们的互导纳共4个元素，当N_element个元件同时故障时，修正的是4*N_element个导纳矩阵元素。

3.根据权利要求1所述的提升配电网灾后恢复用并行静态安全分析方法，其特征在于，步骤S3中，k号线程计算k号节点的有功注入功率P_k和无功注入功率Q_k具体为：

其中，θ_k为k号节点的电压相角，θ_j为j号节点的电压相角，G为导纳矩阵实部，V_j为j号节点的电压幅值，V_k为k号节点的电压幅值，B为导纳矩阵虚部，k为核函数Make_S_Kernel在GPU中的线程号；NYseq[k]为导纳矩阵第k行非零元素起始索引；NYseq[k+1]-1为导纳矩阵第k行非零元素终止索引。

4.根据权利要求1所述的提升配电网灾后恢复用并行静态安全分析方法，其特征在于，步骤S3中，雅可比子矩阵为：H_n×n、N_n×m、M_m×n、L_m×m，其中，m为PQ节点个数、n为系统总结点数-1。

5.根据权利要求4所述的提升配电网灾后恢复用并行静态安全分析方法，其特征在于，首先在CPU上形成雅可比子矩阵H、N、M、L与导纳矩阵的对应关系H2Y、N2Y、M2Y、L2Y，作为并行形成雅可比子矩阵的下标索引；然后在CPU上计算雅可比子矩阵H、N、M、L与雅可比矩阵的元素对应关系H2Ja、N2Ja、M2Ja、L2Ja；根据与导纳矩阵不同节点类型的对应关系，在GPU上创建核函数Make_H_Kernel、Make_N_Kernel、Make_M_Kernel、Make_L_Kernel并行计算雅可比子矩阵H、N、M、L，H子矩阵的核函数线程数为n×n，N、H子矩阵的核函数线程数为n×m，L子矩阵的核函数线程数为m×m。

6.根据权利要求1所述的提升配电网灾后恢复用并行静态安全分析方法，其特征在于，步骤S4中，通过求解极小化J(y)得到y_m向量，确定线性方程组的近似解，通过预处理矩阵M，将线性方程组Ax＝b转化为同解的易于迭代求解的新的方程组M^-1Ax＝M^-1b，选取不完全LU分解(ILU)进行预处理；通过相对残差rr_k+1限制内部GMRES迭代求解算法求解当前步骤解s_k的精度值，通过内部迭代求解得到的相对残差与预设的阈值比较，决定内部迭代是否收敛；当相对残差小于预设的强迫因子η时，内部迭代停止。

7.根据权利要求6所述的提升配电网灾后恢复用并行静态安全分析方法，其特征在于，相对残差rr_k+1具体为：

其中，e_k为误差范数，F(x_k)为，F′(x_k)为雅可比矩阵，s_k为修正量。

8.根据权利要求1所述的提升配电网灾后恢复用并行静态安全分析方法，其特征在于，步骤S5中，在CPU上进行重要元件排序，重要元件判据根据相对越限程度进行判断，定义重要元件判据为：

其中，F_ij为故障i情况下支路j的支路功率或故障i情况下节点j有功、无功；

为系统中支路j的支路功率或节点j标准有功、标准无功。

Images