基于GPU计算平台的大规模电力系统暂态稳定计算.pdf

第４１卷第４期２０１３年２月１６日电力系统保护与控制ＰｏｗｅｒＳｙｓｔｅｍＰｒｏｔｅｃｔｉｏｎａｎｄＣｏｎｔｒｏｌＶｌ０１．４１Ｎｏ．４Ｆｅｂ．１６．２０１３基于ＧＰＵ计算平台的大规模电力系统暂态稳定计算江涵，江全元（浙江大学电气工程学院，浙江杭州３１００２７）摘要：为满足对大规模互联电网进行快速暂态稳定仿真的需求，提出了一种基于图形处理器（ＧｒａｐｈｉｃｓＰｒｏｃｅｓＳｉｎｇＯｎｉｔ，ＧＰＵ）计算平台的暂态稳定并行算法。算法依据暂态稳定联立矩阵的双层对角加边结构，将整体计算分解为三部分：１）动态元件相关计算；２）分区系统计算；３）边界系统计算。第１）、２）部分被分配到计算平台中的多核ＣＰＵ上进行处理。第３）部分则采用可完全并行化的稳定双共轭梯度法在ＧＰＵ上计算，并且为了减少迭代次数使用了稀疏近似逆预处理技术。针对一个含１２８２３节点、１４３１台发电机组的算例进行了仿真测试，验证了所提算法的有效性和实用性。结果显示，算法并行加速比可达到７．０１倍，仿真速度快于实际暂态过程，为暂态稳定并行计算提供了新的解决思路。关键词：暂态稳定；并行计算；图形处理器；双层对角加边结构；稳定双共轭梯度；稀疏近似逆预处理ＡｐａｒａＨｅｌｔｒａｎｓｉｅｎｔｓｔａｂｉｌｉｔｙａｌｇｏｒｉｔｈｍｆｏｒｌａｒｇｅ－ｓｃａｌｅｐｏｗｅｒｓｙｓｔｅｍｂａｓｅｄｏｎＧＰＵｐｌａｔｆｏｒｍ—ＪＩＡＮＧＨａｒｔ，ＪＩＡＮＧＱｕａｎｙｕａｎ（ＣｏｌｌｅｇｅｏｆＥｌｅｃｔｒｉｃａｌＥｎｇｉｎｅｅｒｉｎｇ，ＺｈｅｊｉａｎｇＵｎｉｖｅｒｓｉｔｙ，Ｈａｎｇｚｈｏｕ３１００２７，Ｃｈｉｎａ）Ａｂｓｔｒａｃｔ：Ｉｎｏｒｄｅｒｔｏｓａｔｉｓｆｙｔｈｅｄｅｓｉｒｅｆｏｒｆａｓｔｔｒａｎｓｉｅｎｔｓｔａｂｉｌｉｔｙｓｉｍｕｌａｔｉｏｎ，ａｎｏｖｅｌｍｉｘｅｄａｌｇｏｒｉｔｈｍｂａｓｅｄｏｎＧｒａｐｈｉｃｓＰｒｏｃｅｓｓｉｎｇＵｎｉｔ（ＧＰＵ）ｐｌａｔｆｏｒｍｉｓｉｎｔｒｏｄｕｃｅｄ．Ａｃｃｏｒｄｉｎｇｔｏｔｗｏ－ｌｅｖｅｌｂｌｏｃｋｂｏｒｄｅｒｅｄｄｉａｇｏｎａｌｆｏｒｍ（ＢＢＤＦ）ｏｆｔｈｅｍａｔｒｉｘ，ｔｈｅａｌｇｏｒｉｔｈｍ—ｄｅｃｏｍｐｏｓｅｓｔｈｅｗｈｏｌｅｃｏｍｐｕｔａｔｉｏｎｉｎｔｏｔｈｒｅｅｐａｒｔｓ，ａｎｄｔｈｅｙａｒｅｔｈｅｃｏｍｐｕｔａｔｉｏｎｏｆｄｙｎａｍｉｃｕｎｉｔｓａｎｄｉｎｊｅｃｔｉｏｎｃｕｒｒｅｎｔｓ，ｓｕｂｎｅｔｗｏｒｋｃｏｍｐｕｔａｔｉｏｎａｎｄｔｈｅｂｏｕｎｄａｒｙｎｅｔｗｏｒｋｃｏｍｐｕｔａｔｉｏｎ．Ｔｈｅｆｏｒｍｅｒｔｗｏｐａｒｔｓａｒｅｃｏｍｐｕｔｅｄｏｎｍｕｌｔｉ．ｃｏｒｅＣＰＵｓ．ＴｈｅｌａｓｔｐａｒｔｉＳｓｏｌｖｅｄｏｎＧＰＵｂｙｃｏｍｐｌｅｔｅｌｙｐａｒａｌｌｅｌｂｉｃｏｎｊｕｇａｔｅｄｇｒａｄｉｅｎｔｓｔａｂｉｌｉｚｅｄｍｅｔｈｏｄ．Ｔｏｄｅｃｒｅａｓｅｔｈｅｔｉｍｅｓｏｆｉｔｅｒａｔｉｏｎ，ｔｈｅｓｐａｒｓｅａｐｐｒｏｘｉｍａｔｅｉｎｖｅｒｓｅｐｒｅｔｒｅａｔｍｅｎｔｔｅｃｈｎｉｑｕｅｉＳａｄｏｐｔｅｄ．Ｔｈｅｓｉｍｕｌａｔｉｏｎｒｅｓｕｌｔｓｓｈｏｗｔｈａｔｆｏｒａｌａｒｇｅｎｅｔｗｏｒｋｗｉｔｈ１２８２３ｎｏｄｅｓａｎｄ１４３１ｇｅｎｅｒａｔｉｏｎｓ．ｔｈｅｐｒｏｐｏｓｅｄｆｌｇｏｒｉｔｈｍｒｕｎｓ７．Ｏ１ｔｉｍｅｓｆａｓｔｅｒｔｈａｎｏｎａｓｉｎｇｌｅＣＰＵｃｏｒｅ．ａｎｄｔｈｅｓｉｍｕｌａｔｉｏｎｔｉｍｅｉＳｌｅｓｓｔｈａｎｔｈｅａｃｔｕａｌｔｒａｎｓｉｅｎｔｐｒｏｃｅｓｓｔｉｍｅ．Ｉｔｐｒｏｖｉｄｅｓａｎｅｗｓｏｌｕｔｉｏｎｓｔｒａｔｅｇｙｆｏｒｌａｒｇｅ．ｓｃａｌｅｐｏｗｅｒｓｙｓｔｅｍｓｉｍｕｌｍｉｏｎ．ＴｈｉｓｗｏｒｋｉｓｓｕｐｐｏｓｅｄｂｙＮａｔｉｏｎａｌＨｉ【ｇｈＴｅｃｈｎｏｌｏｇｙＲｅｓｅａｒｃｈａｎｄＤｅｖｅｌｏｐｍｅｎｔＰｒｏｇｒａｍｏｆＣｈｉｎａ（８６３Ｐｒｏｇｒａｍ）（Ｎｏ．２０１１ＡＡ０５Ａ１１８），ＮａｔｉｏｎａｌＰｒｏｇｒａｍｏｎＫｅｙＢａｓｉｃＲｅｓｅａｒｃｈＰｒｏｊｅｃｔ（９７３Ｐｒｏｇｒａｍ）（Ｎｏ．２０１２ＣＢ２１５１０６），ａｎｄＮａｔｉｏｎａｌＮａｔｕｒａｌＳｃｉｅｎｃｅＦｏｕｎｄａｔｉｏｎｏｆＣｈｉｎａ（Ｎｏ．５０９７７０８２）．—Ｋｅｙｗｏｒｄｓ：ｔｒａｎｓｉｅｎｔｓｔａｂｉｌｉｔｙ；ｐａｒａｌｌｅｌｃｏｍｐｕｔａｔｉｏｎ；ＧＰＵ；ｔｗｏｌｅｖｅｌＢＢＤＦ；ｂｉｃｏｎｊｕｇａｔｅｇｒａｄｉｅｎｔｓｔａｂｉｌｉｚｅｄｍｅｔｈｏｄ；ｓｐａｒｓｅａｐｐｒｏｘｉｍａｔｅｉｎｖｅｒｓｅｐｒｅｔｒｅａｔｍｅｎｔｔｅｃｈｎｉｑｕｅ中图分类号：ＴＭ７４４文献标识码：Ａ文章编号：１６７４．３４１５（２０１３）０４－００１３．０８Ｏ引言电力系统暂态稳定仿真是进行电力系统规划、安全稳定分析和电力科学研究的重要工具。电网规模的扩大和仿真模型的日益精细使其计算量不断增长。与此同时，为进一步提高工作效率以及满足对电网进行更快的安全控制，对仿真速度的要求也不基金项目：国家８６３高技术基金项目（２０１１ＡＡ０５Ａ１１８）；国家９７３重点基础研究发展规划项目（２０１２ＣＢ２１５１０６）；国家自然科学基金项目（５０９７７０８２）断提升。在此背景下，并行计算技术被引入暂态稳定仿真领域，并不断发展【ＪＪ。近几十年来，随着并行计算技术的革新和进步，众多科研工作者提出了多种暂态稳定并行算法，一般可分为空间并行算法，时间并行算法和引自ＶＬＳＩ领域的波形松弛算法三类［３】。算法的数学本质都是使用并行方法解决～个大规模微分．代数方程组初值问题，对于空间并行算法，一个主要思路是基于串行联立算法，首先将整个系统形成一个统一的非线性方程组，进而通过牛顿法将方程组线性化，再进行并行处理，代表性的算法有并行ＶＤＨＮ法，．１４．电力系统保护与控制—ＳＯＲ－Ｎｅｗｔｏｎ法，ＭａｃｌａｕｒｉｎＮｅｗｔｏｎ法ｌ等。另一条思路是基于交替算法，直接将微分与代数方程分开处理，两部分交替计算，其并行研究主要集中在任务的平衡分配和网络矩阵的处理【６Ｊ。这些算法的共同之处是都需要处理一个线性代数方程组。对线性方程组的求解，直接法（如ＬＵ分解）占据了主流地位。但是，直接法由于其序列性，不能完全地并行化。在很多科研领域，迭代法已经显示出其优秀的应用性质。在众多迭代算法中，共轭梯度（ＣｏｎｊｕｇａｔｅＧｒａｄｉｅｎｔＭｅｔｈｏｄ，ＣＧ）及其改进算法【８】具有快速收敛、存储量小、不必预先估计参数等特点，其计算过程主要是矩阵向量相乘。此类方法近些年在电力系统潮流、仿真等计算中得到了很多研究应用。文献【９］在超立方计算机上实现了预处理ＣＧ算法，应用于暂态稳定计算。文献［１Ｏ１针对分块式电力系统暂态稳定并行算法，采用乔列斯基分解的预处理共轭梯度法用于边界系统方程组的求解。文献【ｌｌ】将ＩＬＵ共轭梯度法与矩阵方程直接求解法相结合，求解暂态稳定计算中的整个雅可比矩阵，获得了较好的加速效果。文献［１２１在Ｂｅｏｗｕｌｆ集群平台上，提出一种采用多项式预处理的共轭梯度迭代法求解电力系统分析计算中的稀疏线性方程组，收到了一定的效果。但是，共轭梯度算法针对的是实对称矩阵，很多电力系统矩阵并不具备这个特点（如包含移相调节变压器的网络导纳矩阵），会影响算法的收敛性，导致迭代次数增加。针对求解不对称矩阵问题，文献［１３］提出的稳定双共轭梯度法（ＢｉＣｏｎｊｕｇａｔｅＧｒａｄｉｅｎｔＳｔａｂｉｌｉｚｅｄＭｅｔｈｏｄ，ＢｉＣＧＳＴＡＢ）是一种有效的解决方案。基本的并行计算硬件平台主要基于ＰＣ集群和共享内存的多核计算机等多ＣＰＵ体系结构。近年来，图形处理器（ＧｒａｐｈｉｃｓＰｒｏｃｅｓｓｉｎｇＵｎｉｔ，ＧＰＵ）发展迅猛。自１９９９年诞生后，其发展的速度是ＣＰＵ更新速度的３倍多【１训。广泛应用于图像与视频处理等领域。２０１０年，文献［１５】首次将ＧＰＵ引入电力系统计算，在包含ＧＰＵ的多核计算平台上（本文简称为ＧＰＵ计算平台）使用直接法计算电力系统暂态稳—定性。我国学者同时也展开了相应研究：文献［１６１７１分析了ＧＰＵ的特性，并使用其计算电力系统潮流和进行谐波分析，提升了计算速度，提高了分析的实用性。本文提出了一种应用于ＧＰＵ计算平台的电力系统暂态稳定仿真并行算法。首先，在联立算法框架下，依据非诚实牛顿法（ＶｅｒｙＤｉｓｈｏｎｅｓｔＮｅｗｔｏｎＭｅｔｈｏｄ，ＷＨＹ）原理，将联立矩阵解耦，使动态元件与网络方程可以分别计算。进一步，将网络矩阵重排为对角加边（ｂｌｏｃｋｂｏｒｄｅｒｅｄｄｉａｇｏｎａｌｆｏｒｍ，ＢＢＤＦ）形式，使此部分计算分解为子分区系统和边界系统两部分。子分区系统与相连的动态元件交由多核ＣＰＵ并行处理。对于边界系统，本文在ＧＰＵ上使用基于稀疏近似逆预处理的ＢｉＣＧＳＴＡＢ法求解。针对一个１２８２３节点的大型电力系统进行测试，验证了算法的准确性和有效性。结果显示算法提高了边界系统计算速度，可获得较高加速比。１暂态稳定并行算法１．１联立算法基本原理电力系统暂态稳定模型的数学本质为大规模微分代数方程初值问题＝ｆ（ｘ，Ｖ）（１）—０＝ＩＶ（２）式中：为状态向量，描述发电机组等元件的动态行为；为电压向量，描述节点电压变化；Ｊ为电流向量，描述节点注入电流变化。ｙＮ为输电网络导纳矩阵。公式（１）描述电网各种动态元件的行为，如发电机组。各动态元件仅由网络连接，因此，式（１）由多组相互独立的子集组成。采用联立算法时，求解过程如下：１）应用数值积分算法将微分方程组差分化，获得联立代数方程组。本文采用隐式梯形积分法：当ｆ时步状态向量，电压向量已知时，令ｈ为积分步长，由式（１）得１时步的状态向量Ｘ表达式为】＋】＝＋芸，）＋十ｌ，）］（３）由式（２）得，１时步网络方程为’ｙＮ＋１＝，＋ｌｖｔ＋１）（４）定义Ｆ（ｘｔ＋ｌ，）＝＋ｌ一一－￣Ｌｆ（ｘ，，Ｖｔ）＋ｆ（ｘｔ＋】，Ｖｔ＋１）ｌ（５）日（＋１，＋１）ｙＮ＋】一，（＋１，）（６）２）联立计算ＪＦ（ｘｔ＋ｌ，＋１）＝０…ＩＨ（ｘｔ＋１，＋１）＝０式（７）为非线性代数方程组，＋。与＋。为未知量，可采用Ｎｅｗｔｏｎ法求解。第ｋ次迭代计算公式为『Ｉ１Ｉ畿：寓ｌ：－１至０『Ｉ一『＿书ｌ麓ｌ（８）ｌ毗，嘛）广ｌ一【ｃＩ峨【＝－Ｊｔｅ；＇－ｅ￣ｌＹＤ为输电网络导纳ｙＮ（固定不变）与动态元件自导纳ｙＧ之和。（《１，。）与日（．ｌ，。）为第ｋ江涵，等基于ＧＰＵ计算平台的大规模电力系统暂态稳定计算．１５．△次迭代时的残差向量。域，为迭代修正量。雅可比矩阵仅在迭代次数ｋ达到设定阈值时才更新。式（８）迭代收敛后，进入２时步计算。当计及发电机凸极效应时，如图１所示。ｆ（《．）＝＝＝图１发电机网络接口图Ｆｉｇ．１Ｉｎｔｅｒｆａｃｅｏｆｔｈｅｇｅｎｅｒａｔｏｒｗｉｔｈｔｈｅｎｅｔｗｏｒｋ图中，ｙＧＮ、ｓ分别为ｙＧ中的固定部分和可变部分，即ＹＧ＝ＹＧ￣－ｙＧｓ。ｆ（）为发电机注入电流Ｊ，）中与机端电压无关部分。由图及ＶＤＨＮ原理引，式（８）的计算可简化为—馘１＝－Ａ～ｎａｙ／＋ｌＡＦ（《１，１）（９）Ｊ＝，（）一ｌ＋ＣＡＦ（《ｌ，１）（１０）略＝ｊ（１１）式中，Ｙ＝（ｙＮ＋ｙＧＮ），其在计算过程中为固定值，仅需在网络拓扑变化时重新形成，十分利于算法并行化。１．２双层ＢＢＤＦ并行算法由电力系统特点，可对上述联立算法进行双层ＢＢＤＦ并行分解。第ｌ层：电力系统中，发电机组通过电网连接，这种特点使式（８）的具体结构如式（１２）所示。ｌＡ。…ｌｌ…００．．Ｉ△《１．１……一一…。。…一…一…‘…。一……‘。｛ｃｌ△《１．Ｇ△ｌ（１２）…Ｉ…Ｆ。Ｊ其中：Ｇ代表发电机组总数。ｇ代表第ｇ台发电机；，，Ｃｇ，域１，，分别为，Ｂ，Ｃ，嗽１，的子矩阵。定义虞为Ｂ中非零向量，为中非零向量。因此，动态元件（发电机组）相关计算，即式（９）、式（１０）可分解为Ａｇｚ￣Ｘｔｋ＋ｌ，ｇ＝一—ＢｇＡＥ￣Ｉ，ｇＦ３．．１，ｇ（１３）＝（《）一，咝，＋４（１４）∈其中：ｇ［１，Ｇ］，，为中元素，代表一个节点的电压值；为ｊ中元素，代表一个节点的注入电流值。此部分的计算具有内在并行性，可分组后并行计算。第２层：由于整个电网为大型稀疏网络，可将电网分为印个子分区系统和一个边界系统，使式（１１）成如下结构【９］．譬嗡：－喘嘴Ｉｑ。一Ｉｂ（１５）∈其中：ｂ代表边界系统；／代表第子分区，Ｊ［１，ｑ】；，，Ｖ／＋分别表示第ｕ，分区注入电流向量与电压向量。因此，式（１１）可分解为边界系统计算∑［一Ｉ），Ｔ－∑＝一（１６）ｊ＝ｌｊ＝ｌ子分区系统计算，＝，ｆ一ｆ（１６）本文中子分区系统式（１６）的计算与相连的发电机组分配到相同ＣＰＵ核心并行计算。同时需要说明的是，为便于处理式（１６），本文使用实数导纳。１。３基于稳定双共轭梯度的边界网络方程计算应用１．２节并行算法后，式（１６）成为唯一需要串行计算的部分，传统方法使用ＬＵ分解法进行计算。文献［１０，１９１中引入ＣＧ算法进行并行迭代求解，简化式（１６）成式（１８）形式。＝（１７）其中：ｙ代表左边项【一ＪｙＪ一］；代表ｊ＝ｌ［ｊｆ６一ｊ，】；代表，方程阶数为。ｊ＝ｌ可将其构造成无约束优化问题＝０．５Ｖ一７，ｆ（ｖ）的极小值点便是式（１７）的解。其局限性是ＣＧ法本质上适合于处理对称实数矩阵，而式ｆ１６）中，ｌ，是经各系统导纳矩阵运算后的合成矩阵，并且当网络中包含移相调节变压器时为非对称实矩阵ｆ２。】。其在某些情况下也不再具有对角矩阵占优这个通常电力系统导纳矩阵具备的特点［１们，条件数会增大，甚至使方程出现病态，使求解时ＣＧ方法的迭代次数■蔓．．．．Ｌ．．．．Ｌ．．．一¨～一～～～一芎¨一＝：一一¨～一．ｈ一．“＿ｌ＿¨～～一ｒ～～～～ｈ一．一～～．＝：～一一一¨＿＿＿一～一～一：＝～一～～～～．１６．电力系统保护与控制增加。针对不对称实矩阵的迭代求解，ＢｉＣＧＳＴＡＢ法被证明是一种十分有效的改进方案，其基于双边Ｌａｎｃｚｏｓ算法和残差正交子空间的理论引，具有双正交和平方收敛的特点，运算量小，收敛性优于ＣＧ算法［２１１。同时，从算法流程上，其继承了ＣＧ法涉及大量的行（列）的乘积与求和的特点，并行性很高。因此，非常适合用于边界系统的并行求解。对式（１７）进行求解时，矩阵Ｉ，的条件数也将影响算法的收敛性。而预处理技术可以降低矩阵的条件数。因此，ＢｉＣＧＳＴＡＢ算法一般与预处理算法相结合，算法一般流程如下。步骤１）给定系数矩阵】，，向量Ｊ，及预处理矩阵Ｍ（的近似矩阵）。设初始值，最大迭代次数ｆｍ，最大容许误差８ｍ，计算ｒ０＝一。，并令ｉ＝１，＝ｒｏ。步骤２）若ｆ／ｍ且ｓｓ，则转到３）；否则，得到的值，结束迭代；步骤３）Ｐｆ＿１＝一１步骤４）若ｆ：１，那么Ｐｉ＝ｒｏ；否则＝（Ｐｊ＿ｌ／ｐｆ＿２）（ａ／ｃｏｉ一１）—ｐｆ＝一１＋（一１０９ｉ一１Ｈｉ一１）步骤５）计算／ｇｉ＝，０［＝Ｐｆ一１／（ｆ），＝一１一。步骤６）￣＝ｌｌｓｌｌ，若ｓ＜ｓ，则Ａｖ＝，１，＝Ｖｆ＿１＋Ａｖｆ，结束迭代；否则进入下一步。∞步骤７）计算＝Ｍｓ，ｔ＝，＝ｔＴｓ／ｔｔ，ｌ，ｆ＝ｌ，ｆ一１—＋ｆ＋ｑ，ｒｉ＝ＣＯｆｔ，ｓ＝ｌＩ，令ｆ＝ｆ＋１，转到２。此计算过程都为矩阵向量间乘加操作，可完全并行化。其中，的构造目前较常用的是ＩＬＵ预处理方法，但并行性有限。近年来，稀疏近似逆预处理技术以其天然的并行性而受到重视，得到迅速的发展。而且，其具有较少的自定义的参数，不确定性小。其原理为：若能构造两个相互关于Ｉ，共…轭的向量集合｛ｚ＿１．，…和｛ｆ｝ｆ＿１＿组成两个矩阵……Ｚ＝【７，１，ｚ２，，ｚ］，Ｗ＝［１，２，，¨＇】，则……ｙｚ＝Ｄ＝ｄｉａｇ（ｄ１，，，，）（１８）≠其中，，：死，０。那么～＝∑ＺＤＷ：（）／（１９）／＝１ｚ，，Ｄ的完整构造过程如文献【２４】所述。我们为了得到一个尽量稀疏的预条件子，ｚ，，可被不完全地计算，即扔掉那些绝对值小于给定值的元≈≈≈素，即２Ｚ，Ｗ，Ｄ，则稀疏近似逆预条件子Ｍ２ｂ。２ＧＰＵ计算平台简介及算法实现２．１计算平台架构图２显示了ＧＰＵ计算平台的架构。计算平台包含两颗主频为２．６６Ｈｚ的ＩｎｔｅｌＸｅｏｎＥ５６５０处理器，每颗ＣＰＵ包含６个核４￣，（ｃｏｒｅ），共享２４Ｇ内存。ＧＰＵ与ＣＰＵ通过ＰＣＩｅ２．０总线相连。ＧＰＵ型号为ＧｅＦｏｒｃｅＧＴＸ４６０，其包含多个流处理器组（ｓｔｒｅａｍｍｕｌｔｉ．ｐｒｏｃｅｓｓｏｒ，ＳＭ），每个ＳＭ中包含３２（或４８）个流处理器（ｓｔｒｅａｍＰｒｏｃｅｓｓｏｒ，ＳＰ）。ＳＰ是ＧＰＵ上进行数据计算的基本单元。——在ＮＶＩＤＩＡ公司开发编程环境统一计算设备架构（ＣｏｍｐｕｔｅＵｎｉｆｉｅｄＤｅｖｉｃｅＡｒｃｈｉｔｅｃｔｕｒｅ，ＣＵＤＡ）下，如图３所示，一个程序分为两部分：ｈｏｓｔ端和ｄｅｖｉｃｅ端。Ｈｏｓｔ端是指在ＣＰＵ上执行的部分，Ｈｏｓｔ端可使用ＭＰＩ等技术并行。图２ＧＰＵ计算平台硬件架构Ｆｉｇ．２ＨａｒｄｗａｒｅａｒｃｈｉｔｅｃｔｕｒｅｏｆＧＰＵｐｌｕｇｇｅｄｉｎｔｏａＰＣｍｏ￣ｅｒｂｏａｒｄ图３ＣＵＤＡ计算模式Ｆｉｇ．３ＨｉｅｒａｒｃｈｙｏｆｃｏｍｐｕｔｉｎｇｓｔｒｕｃｔｕｒｅｏｆＣＵＤＡｄｅｖｉｃｅ端是在ＧＰＵ芯片上执行的部分，称为内核（ｋｅｒｎｅ１）函数。在ｋｅｒｎｅｌ中，当一个指令被触发时，将被大量ＧＰＵ线程并行执行，每个线程被分配给一个ｓＰ，而且线程的创建、撤销和切换３种操作Ｉ江涵，等基于ＧＰＵ计算平台的大规模电力系统暂态稳定计算．１７．仅需要极少的系统开销［２引。每个线程都具有唯一的编号（ｉｎｄｅｘ），供编程时使用。为便于组织和使用这些ＧＰＵ线程，ＣＵＤＡ中还定义了两层结构：网格（ｇｒｉｄ）和块（ｂｌｏｃｋ）。一个ｇｒｉｄ包含多个ｂｌｏｃｋ，每个ｂｌｏｃｋ管理多个线程（Ｔｈｒｅａｄ）。２．２算法实现２．２．１基于ＧＰＵ的边界系统计算本文中边界系统由１－３节所述ＢｉＣＧＳＴＡＢ算法计算，此方法计算过程包含大量的矩阵．矢量，矢量．矢量计算。以两个矢量相加为例，＝＋Ａｖ。在ＣＰＵ上，此类计算的实现模式是ｆｏｒ（ｉｎｄｅｘ：０；ｉｎｄｅｘ＜，ｚ；ｉｎｄｅｘ＋＋）ｖｉ［ｉｎｄｅｘ】÷＿－Ｖｉ一１［ｉｎｄｅｘ】＋Ａｖｉ［ｉｎｄｅｘ］循环计数器ｉｎｄｅｘ从０增加到ｎ，向量的相应元素顺序进行加运算。但在此计算过程中，每个元素的计算实际是相互独立的，没有先后依存关系。因此在ＧＰＵ上，整个向量一次性使用ｎ个ＧＰＵ线程同时操作，如下所示ｖｉ［ｉｎｄｅｘ】÷…＿－Ｖｉ一１［ｉｎｄｅｘ】＋Ａｖｉ［ｉｎｄｅｘ］，ｉｎｄｅｘ＝１，，ｎ其中，ｉｎｄｅｘ既表示向量元素的编号，同时定位了ＧＰＵ线程，所有的数据是同步并行操作的。ＧＰＵ可对成百乃至上万的元素并行操作，这是ＣＰＵ所不能比拟的。２．２．２基于ＣＰＵ的分区系统及动态元件计算的分配根据１．２节算法，使用ＭＰＩ将式（１３）、式（１４）、式（１６）的计算分配到多个ＣＰＵ核心中计算。在多核平台上，为了获得较高的并行效率，需根据计算量，进行任务分配。本文依据文献［７】所述方法进行任务分配，分配指标为ｍｉｎ｛Ｍａｘｊ（Ｃｏ．￣Ｃｏｓｔｊ）＋ＣｏｍｐＣｏｓ６｝（２０）其中：为任务数（＝ｃＰｕ核心数＝分区总数），ＣｏｍｐＣｏｓｂ代表了任务／的总计算量，必须说明的是此计算量由网络子分区和相连的发电机组计算量共同确定，即式（１３）、式（１４）、式（１７）；ＣｏｍｐＣｏｓｔＢ代表了边界系统的计算量。由于本文中使用共享内存的多核计算机，其通讯损耗小于ＰＣ集群，因此本文算法相当于放宽了对边界系统的约束，而更加重视子分区计算量的平均分配。从并行计算的角度，算法一次迭代流程如图４所示。３算例分析３．１测试方案本文以联立算法为框架，在ＧＰＵ计算平台上实现了并行暂态稳定算法。ＣＰＵ函面———————协调系统电压ｌ计算子分区电压———————ｒ＿一更新子分区内发电机组状态变量Ｖ／＼更新阈值／更新Ｊａｃｏｂｉａｎ矩阵Ｉ击ｌ飘ｌ注入电流计算并向协调系统映射ｌ子分区ｑ计算子分区电压二二］［＝二更新子分区内发电机组状态变量／是否达到＼翌星土更新Ｊａｃｏｂｉａｎ矩阵＝＝＝＝］［二＝＝＝子分区注入电流计算并向协调系统映射否ＧＰＵ由近似逆预处理ＢｉｃＧｓＴＡＢ算法计算边界系统图４并行算法流程Ｆｉｇ．４Ｆｌｏｗｃｈａｒｔｏｆａｎｉｔｅｒａｔｉｏｎｏｆｔｈｅｐａｒａｌｌｅｌａｌｇｏｒｉｔｈｍ本文算法本质上是一种组合算法，为验证其有效性，使用２种算法组合进行效率比较测试，两种算法外层框架都使用双层ＢＢＤＦ并行算法，边界系统分别采用基于稀疏近似逆预条件子的ＣＧ法（简称为ＰＣＧ）和基于稀疏近似逆预条件子的ＢｉＣＧＳＴＡＢ法（简称为ＰＢｉＣＧＳＴＡＢ）计算，如表１所示。表１算法列表Ｔａｂｌｅ１Ｌｉｓｔｓｏｆａｌｇｏｒｉｔｈｍｓ每时步迭代收敛条件均为１Ｏ～。算法步长为０．０１Ｓ，ＣＧ与ＢｉＣＧＳＴＡＢ法收敛精度为１０一。测试方案为仿真１０Ｓ暂态过程，０．０Ｓ系统中某条线路近母线处发生三相接地短路，Ｏ．０８３Ｓ故障切除。加速比公式为ｓｐ＝，其中表示单核』ｐａｒａｌｌｅｌ串行计算时间，。。表示并行计算时间。３．２测试结果及分析本文使用某１２８２３节点系统，系统由１４３１台发电机组（发电机采用６阶模型，并考虑励磁控制及原动机调速系统），１２３６５条线路，５６１２台变压器（包含部分移相调节变压器）组成，串行算法仿真时间１为６８．４１Ｓ。首先，选取系统中３个节点电压曲线，对本文电力系统保护与控制算法（即算法２）精度进行对比分析，从图５中可以看出，串行算法与并行算法电压曲线基本相同，可以保证计算精度。ｊ型１口匝（一一线为串行算法电压曲线，一线为算法２电压曲线）图５算法电压曲线比较Ｆｉｇ．５Ｃｏｍｐａｒｉｓｏｎｏｆｖｏｌｔａｇｅｃｕｒｖｅｓｂｅｔｗｅｅｎｄｉｆｆｅｒｅｎｔａｌｇｏｒｉｔｈｍｓ算法１、２并行框架同为双层ＢＢＤＦ法，因此边界系统规模是相同的。边界系统的划分结果及计算时间如表２所示，其中表示每时步内ＰＣＧ法求解的平均迭代次数，Ｉ＿厕表示每时步内ⅣＰＢｉＣＧＳＴＡＢ法求解的平均迭代次数。Ｈ为渤与户ＣＧ之比，为ＰＢｉＣＧＳＴＡＢ法与ＰＣＧ法对边界系统计算总时间之比。表２边界系统计算结果对比Ｔａｂｌｅ２Ｃｏｍｐａｒｉｓｏｎｏｆｔｈｅｓｉｍｕｌ￣ｉｏｎｒｅｓｕｌｔｓｏｆｔｈｅｂｏｕｎｄａｒｙｓｙｓｔｅｍ由表２可知，分区数从２增加到１２的过程中，边界系统由２节点增长到１３１节点。从迭代次数的角度，ＰＣＧ法的增加较快，平均次数由２．５次增加到７．１次，而ＰＢｉＣＧＳＴＡＢ法迭代次数仅由２．１次增加到３．１次。从实际计算时间的角度，当边界系统增大到１３１节点时，ＰＢｉＣＧＳＴＡＢ的时间仅为ＰＣＧ法的０．６２。仿真整体加速比和总时间如图６、图７所示。从原理上，两种算法的动态元件与网络子分区部分的计算策略是相同的，而边界系统算法不同。因此，图中两种算法在并行分区较少时加速比几乎相同，但随着分区的增多，算法１加速比接近饱和，算法２的加速比逐渐高于算法１。这是由于算法２较好的处理了分区增多导致的边界系统增大的情况。算法１目算法２６．Ｏ８５．３４７Ｏ１２４８１２ＣＰＵ核心数图６算法加速比对比Ｆｉｇ．６Ｃｏｍｐａｒｉｓｏｎｏｆｓｐｅｅｄ－ｕｐｂｅｔｗｅｅｎｔｗｏａｌｇｏｒｉｔｈｍｓ图７算法总计算时间对比Ｆｉｇ．７Ｃｏｍｐａｒｉｓｏｎｏｆｃｏｍｐｕｔａｔｉｏｎｔｉｍｅｂｅｔｗｅｅｎｔｗｏａｌｇｏｒｉｔｈｍｓ在串行环境下，仿真时间大于实际暂态过程。在并行环境下，当使用１２个核时，算法２｝ＪＨ速比达到了７．０１倍。仿真时间是实际暂态过程的９７．６％，小于实际暂态过程，说明了了算法２的有效性和一定的实用性。４结论本文在联立算法的基础上，依据ＶＤＨＮ原理和双层ＢＢＤＦ算法将整体暂态稳定计算分为３部分：１）动态元件相关计算。２）子分区系统计算。３）边界系统计算。在高性能ＧＰＵ计算平台上，首先将１），１１２）交由多核ＣＰＵ负责的计算。然后，重点研究了边界系统计算，针对其方程组左边项为不对称实矩阵的特点，引入了ＢｉＣＧＳＴＡＢ算法，并进一步通过稀疏近似逆预处理技术加快收敛速度。最后，通过算例测试表明该算法性能优于稀疏近似逆预处理ＣＧ并行算法，有效的降低了边界系统的计算时间，获得了较高的并行加速比，对于１２８２３节点算例达到了实时仿真。参考文献［１］石恒初．基于ＰＣ机群的电力系统暂态稳定评估【Ｊ］．电８７６５４３２ｌ０景江涵，等基于ＧＰＵ计算平台的大规模电力系统暂态稳定计算．１９一力系统保护与控制，２００９，３７（１０）：５－１４．—ＳＨＩＨｅｎｇｃｈｕ．Ｐｅｒｓｏｎａｌｃｏｍｐｕｔｅｒｃｌｕｓｔｅｒｂａｓｅｄｐｏｗｅｒｓｙｓｔｅｍｔｒａｎｓｉｅｎｔｓｔａｂｉｌｉｔｙａｓｓｅｓｓｍｅｎｔ［Ｊ］．ＰｏｗｅｒＳｙｓｔｅｍＰｒｏｔｅｃｔｉｏｎａｎｄＣｏｎｔｒｏｌ，２００９，３７（１０１：５．１４．［２］汪芳宗，张磊．基于辛Ｒａｄａｕ方法的暂态稳定性并行计算方法［Ｊ］．电力系统保护与控制，２０１１，３９（２３）：…３５．３８．—ＷＡＮＧＦａｎｇｚｏｎｇ，ＺＨＡＮＧＬｅｉ．ＰａｒａｌｌｅｌａｌｇｏｒｉｔｈｍｂａｓｅｄｏｎｓｙｍｐｌｅｃｔｉｃＲａｄａｕｍｅｔｈｏｄｆｏｒｔｒａｎｓｉｅｎｔｓｔａｂｉｌｉｔｙｓｉｍｕｌａｔｉｏｎ［Ｊ］．ＰｏｗｅｒＳｙｓｔｅｍＰｒｏｔｅｃｔｉｏｎａｎｄＣｏｎｔｒｏｌ，２０１１，３９（２３）：３５．３８．［３］薛巍，舒继武，王心丰，等．电力系统暂态稳定仿真并行算法的研究进展［Ｊ］．系统仿真学报，２００２，１４（２）：１７７．１８２．—ＸＵＥＷｅｉ，ＳＨＵＪｉｗｕ，ＷＡＮＧＸｉｎ．ｆｅｎｇ，ｅｔａ１．Ａｄｖａｎｃｅｏｆｐａｒａｌｌｅｌａｌｇｏｒｉｔｈｍｆｏｒｐｏｗｅｒｓｙｓｔｅｍｔｒａｎｓｉｅｎｔｓｔａｂｉｌｉｔｙｓｉｍｕｌａｔｉｏｎ［Ｊ］．ＪｏｕｒｎａｌｏｆＳｙｓｔｅｍＳｉｍｕｌａｔｉｏｎ，２００２，１４（２）：１７７．１８２．１４ＪＣｈａｉＪＳ，ＺｈｕＮ，ＢｏｓｅＡ，ｅｔａ１．ＰａｒａｌｌｅｌＮｅｗｔｏｎｔｙｐｅｍｅｔｈｏｄｓｆｏｒｐｏｗｅｒｓｙｓｔｅｍｓｔａｂｉｌｉｔｙａｎａｌｙｓｉｓｕｓｉｎｇｌｏｃａｌａｎｄｓｈａｒｅｄｍｅｍｏｒｙｍｕｌｔｉｐｒ０ｃｅｓｓｏｒｓ［Ｊ】．ＩＥＥＥＴｒａｎｓｏｎＰｏｗｅｒＳｙｓｔｅｍｓ，１９９１，６（４）：１５３９．１５４５．［５］ＣｈａｉＪＳ，ＢｏｓｅＡ．Ｂｏｔｔｌｅｎｅｃｋｓｉｎｐａｒａｌｌｅｌａｌｇｏｒｉｔｈｍｓｆｏｒｐｏｗｅｒｓｙｓｔｅｍｓｔａｂｉｌｉｔｙａｎａｌｙｓｉｓ［Ｊ］．ＩＥＥＥＴｒａｎｓｏｎＰｏｗｅｒＳｙｓｔｅｍｓ，１９９３，８（１）：９－１５．［６］李亚楼，周孝信，吴中习．一种可用于大型电力系统数字仿真的复杂故障并行计算方法［Ｊ］．中国电机工程学报，２００３，２３（１２）：１－５．——ＬＩＹａ－ｌｏｕ，ＺＨＯＵＸｉａｏｘｉｎ，ＷＵＺｈｏｎｇｘｉ．Ａｐａｒａｌｌｅｌ—ｃｏｍｐｌｅｘｆａｕｌｔｃｏｍｐｕｔｍｉｏｎａｌｇｏｒｉｔｈｍｆｏｒｌａｒｇｅｓｃａｌｅｐｏｗｅｒｓｙｓｔｅｍｄｉｇｉｔａｌｓｉｍｕｌａｔｉｏｎ［Ｊ］．ＰｒｏｃｅｅｄｉｎｇｓｏｆｔｈｅＣＳＥＥ，２００３，２３（１２）：１－５．—［７］ＳＨＵＪｉ－ｖｖ２１，ＷＥＩＸｕｅ，ＺＨＥＮＧＷｅｉｍｉｎ．Ａｐａｒａｌｌｅｌｔｒａｎｓｉｅｎｔｓｔａｂｉｌｉｔｙｓｉｍｕｌａｔｉｏｎｆｏｒｐｏｗｅｒｓｙｓｔｅｍｓ［Ｊ］．ＩＥＥＥＴｒａｎｓｏｎＰｏｗｅｒＳｙｓｔｅｍｓ，２００５，２０（４）：１７０９．１７１７．１８ｊｄｅＰｉｌｌｉｓＬＧＡｃｏｍｐａｒｉｓｏｎｏｆｉｔｅｒａｔｉｖｅｍｅｔｈｏｄｓｆｏｒｓｏｌｖｉｎｇｎｏｎｓｙｍｍｅｔｒｉｃｌｉｎｅａｒｓｙｓｔｅｍｓ［Ｊ］．ＡｃｔａＡｐｐｌｉｃａｎｄａｅＭａｔｈｅｍａｔｉｃａｅ，１９９８，５１（２）：１４１．１５９．１９ＪＤｅｃｋｅｒＩＣ，ＦａｌｃａｏＤＭ，ＫａｓｚｋｕｒｅｗｉｃｚＥ．Ｃｏｎｊｕｇａｔｅｇｒａｄｉｅｎｔｍｅｔｈｏｄｓｆｏｒｐｏｗｅｒｓｙｓｔｅｍｄｙｎａｍｉｃｓｉｍｕｌａｔｉｏｎｏｎｐａｒａｌｌｅｌｃｏｍｐｕｔｅｒｓ［Ｊ］．ＩＥＥＥＴｒａｎｓｏｎＰｏｗｅｒＳｙｓｔｅｍｓ，１９９６，１１（３）：１２１８－１２２７．［１０］韩晓言，韩祯祥．预处理并轭梯度法在电力系统暂态稳定分析并行算法中的应用研究［Ｊ］．电力系统及其自动化学报，１９９６，８（２）：１－６．ＨＡＮＸｉａｏ－ｙａｎ，ＨＡＮＺｈｅｎ－ｘｉａｎｇ．Ｒｅｓｅａｒｃｈｏｎｔｈｅａｐｐｌｉｃａｔｉｏｎｏｆｐａｒａｌｌｅｌａｌｇｏｒｉｔｈｍｏｆｐｒｅｃｏｎｄｉｔｉｏｎｉｎｇｃｏｎｊｕｇａｔｅｄｇｒａｄｉｅｎｔｏｆｐｏｗｅｒｓｙｓｔｅｍｔｒａｎｓｉｅｎｔｓｔａｂｉｌｉｔｙ—ａｎａｌｙｓｉｓ［Ｊ］．ＰｒｏｃｅｅｄｉｎｇｓｏｆｔｈｅＣＳＵＥＰＳＡ，１９９６，８（２）：１．６．林济铿，仝新宇，李杨春，等．基于预处理共轭梯度法的电力系统机电暂态仿真［Ｊ］．电工技术学报，２００８，—２３（５）：９３９９．—ＬＩＮＪｉ－ｋｅｎｇ，ＴＯＮＧＸｉｎｙｕ，ＬＩＹａ—ｎｇｃｈｕｎ，ｅｔａ１．Ｅｌｅｃｔｒｉｃａｌ－ｍａｃｈａｎｉｃａｌｔｒａｎｓｉｅｎｔｓｉｍｕｌａｔｉｏｎｏｆｐｏｗｅｒｓｙｓｔｅｍｂａｓｅｄｏｎｐｒｅｃｏｎｄｉｔｉｏｎｅｄｃｏｎｊｕｇａｔｅｇｒａｄｉｅｎｔｍｅｔｈｏｄ［Ｊ］．ＴｒａｎｓａｃｔｉｏｎｓｏｆＣｈｉｎａＥｌｅｃｔｒｏｔｅｃｈｎｉｃａｌＳｏｃｉｅｔｙ，２００８，２３（５）：９３・９９．［１２］刘洋，周家启，谢开贵，等．基于Ｂｅｏｗｕｌｆ集群的大规模电力系统方程并行ＰＣＧ求解［Ｊ］．电工技术学报，—２００６，２１（３）：１０５１１１．ＬＩＵＹａ—ｎｇ，ＺＨＯＵＪｉａｑｉ，ＸＩＥＫａｉ－ｇｕｉ，ｅｔａ１．ＰａｒａｌｌｅｌＰＣＧｓｏｌｕｔｉｏｎｏｆｌａｒｇｅｓｃａｌｅｐｏｗｅｒｓｙｓｔｅｍｅｑｕａｔｉｏｎｓｂａｓｅｄｏｎＢｅｏｗｕｌｆｃｌｕｓｔｅｒ［Ｊ］．ＴｒａｎｓａｃｔｉｏｎｓｏｆＣｈｉｎａＥｌｅｃｔｒｏｔｅｃｈｉｃａｌＳｏｃｉｅｔｙ，２００６，２１（３）：１０５－１１１．［１３］ｖａｎｄｅｒＶｏｒｓｔＨＡ．ＢＩ．ＣＧＳＴＡＢ：ａｆａｓｔａｎｄｓｍｏｏｔｈｌｙｃｏｎｖｅｒｇｉｎｇｖａ—ｒｉａｎｔｏｆＢＩＣＧｆｏｒｔｈｅｓｏｌｕｔｉｏｎｏｆｎｏｎｓｙｍｍｅｔｒｉｃｌｉｎｅａｒｓｙｓｔｅｍｓ［Ｊ］．ＳＩＡＭＪｏｕｒｎａｌｏｎＳｃｉｅｎｔｉｆｉｃａｎｄＳｔａｔｉｓｔｉｃａｌＣｏｍｐｕｔｉｎｇ，１９９２，１３（２）：６３１．６４４．［１４］张浩，李利军，林岚．ＧＰＵ的通用计算应用研究［Ｊ］．计算机与数字工程，２００５（１２）：６０．６２．ＺＨＡＮＧＨａｏ，ＬＩＬｉ－ｊｕｎ，ＬＩＮＬａｎ．Ｇｅｎｅｒａ１．ｐｕｒｐｏｓｅｃｏｍｐｕｔａｔｉｏｎｏｎｇｒａｐｈｉｃｓｐｒｏｃｅｓｓｏｒｓ［Ｊ］．Ｃｏｍｐｕｔｅｒ＆ＤｉｇｉｔａｌＥｎｇｉｎｅｅｒｉｎｇ，２００５（１２）：６０－６２．——［１５］ＪａｌｉｌｉＭａｒａｎｄｉＶＤｉｎａｖａｈｉＶＳＩＭＤ．ｂａｓｅｄｌａｒｇｅｓｃａｌｅｔｒａｎｓｉｅｎｔｓｔａｂｉｌｉｔｙｓｉｍｕｌａｔｉｏｎｏｎｔｈｅｇｒａｐｈｉｃｓｐｒｏｃｅｓｓｉｎｇｕｎｉｔ［Ｊ］．ＩＥＥＥＴｒａｎｓｏｎＰｏｗｅｒＳｙｓｔｅｍｓ，２０１０，２５（３）：ｌ５８９．１５９９．［１６］夏俊峰，杨帆，李静，等．基于ＧＰＵ的电力系统并行潮流计算的实现［Ｊ］．电力系统保护与控制，２０１０，３８（１８１：１００－１０３．—ＸＩＡＪｕｎｆｅｎｇ，ＹＡＮＧＦａｎ，ＬＩＪｉｎｇ，ｅｔａ１．ＩｍｐｌｅｍｅｎｔａｔｉｏｎｏｆｐａｒａｌｌｅｌｐｏｗｅｒｆｌｏｗｃａｌｃｕｌａｔｉｏｎｂａｓｅｄｏｎＧＰＵ［Ｊ］．ＰｏｗｅｒＳｙｓｔｅｍＰｒｏｔｅｃｔｉｏｎａｎｄＣｏｎｔｒｏｌ，２０１０，３８（１８）：—１００１０３．［１７］韩志伟，刘志刚，鲁晓帆，等．基于ＣＵＤＡ的高速并行小波算法及其在电力系统谐波分析中的应用［Ｊ］．电力自动化设备，２０１０，３８（１）：９８．１０１．—ＨＡＮＺｈｉｗｅｉ，ＬＩＵＺｈｉ－ｇａｎｇ，ＬＵＸｉａｏ．ｆａｎ，ｅｔａ１．Ｈｉｇｈ－ｓｐｅｅｄｐａｒａｌｌｅｌｗａｖｅｌｅｔａｌｇｏｒｉｔｈｍｂａｓｅｄｏｎＣＵＤＡ．２０．电力系统保护与控制ａｎｄｉｔｓａｐｐｌｉｃａｔｉｏｎｉｎｐｏｗｅｒｓｙｓｔｅｍｈａｒｍｏｎｉｃａｎａｌｙｓｉｓ［Ｊ］．ＥｌｅｃｔｒｉｃＰｏｗｅｒＡｕｔｏｍａｔｉｏｎＥｑｕｉｐｍｅｎｔ，２０１０，３８（１）：９８．１０１．江涵，江全元．一种可变步长的暂态稳定自适应修正牛顿组合算法【Ｊ］．中国电机工程学报，２０１１，３１（３４）：１０５．１１２．ＪＩＡＮＧＨａｎ，ＪＩＡＮＧＱｕａｎ・ｙｕａｎ．Ｐｏｗｅｒｓｙｓｔｅｍｔｒａｎｓｉｅｎｔ—ｓｔａｂｉｌｉｔｙｓｉｍｕｌａｔｉｏｎｂａｓｅｄｏｎＳｈａｍａｎｓｋｉｉＶＤＨＮａｌｇｏｒｉｔｈｍｗｉｔｈｖａｒｉａｂｌｅｓｔｅｐｓｉｚｅ［Ｊ］．ＰｒｏｃｅｅｄｉｎｇｓｏｆｔｈｅＣＳＥＥ，２０１１，３１（３４）：１０５一ｌ１２．ＤｅｃｋｅｒＩＣ，ＦａｌｃａｏＤＭ，ＫａｓｚｋｕｒｅｗｉｃｚＥ．Ｐａｒａｌｌｅｌｉｍｐｌｅｍｅｎｔａｔｉｏｎｏｆａｐｏｗｅｒｓｙｓｔｅｍｄｙｎａｍｉｃｓｉｍｕｌａｔｉｏｎｍｅｔｈｏｄｏｌｏｇｙｕｓｉｎｇｔｈｅｃｏｎｊｕｇａｔｅｇｒａｄｉｅｎｔｍｅｔｈｏｄ［Ｊ］．ＩＥＥＥＴｒａｎｓｏｎＰｏｗｅｒＳｙｓｔｅｍｓ，１９９２，７（１）：４５８－４６５￣．昆德．电力系统稳定与控制［Ｍ］．北京：中国电力出版社，２００１．Ｋｕｎｄｕｒ．Ｐｏｗｅｒｓｙｓｔｅｍｓｔａｂｉｌｉｔｙａｎｄｃｏｎｔｒｏｌ［Ｍ］．Ｂｅｉｊｉｎｇ：ＣｈｉｎａＥｌｅｃｔｒｉｃＰｏｗｅｒＰｒｅｓｓ，２００１．ＴｏｐｓａｋａｌＥ，ＫｉｎｄｔＲ，ＳｅｒｔｅｌＫ，ｅｔａ１．ＥｖａｌｕａｔｉｏｎｏｆｔｈｅＢＩＣＧＳＴＡＢ（１）ａｌｇｏｒｉｔｈｍｆｏｒｔｈｅｆｉｎｉｔｅ－－ｅｌｅｍｅｎｔ／ｂｏｕｎｄａｒｙ・－ｉｎｔｅｇｒａｌｍｅｔｈｏｄ［Ｊ］．ＡｎｔｅｎｎａｓａｎｄＰｒｏｐａｇａｔｉｏｎＭａｇａｚｉｎｅ，ＩＥＥＥ，２００１，４３（６）：１２４－１３１．谷同祥，迟学斌，刘兴平．稀疏近似逆与多层块ＩＬＵ预条件技术［Ｊ］．应用数学和力学，２００４，２５（９）：９２７－９３４．ＧＵＴｏｎｇ－ｘｉａｎｇ，ＣＨＩＸｕｅ－ｂｉｎ，ＬＩＵＸｉｎｇ－ｐｉｎｇ．ＡＩＮＶａｎｄＢＩＬＵＭｐｒｅｃｏｎｄｉｔｉｏｎｉｎｇｔｅｃｈｎｉｑｕｅｓ［Ｊ］．Ａｐｐｌｉｅｄ—ＭａｔｈｅｍａｔｉｃｓａｎｄＭｅｃｈａｎｉｃｓ，２００４，２５（９）：９２７９３４．［２３］汪芳宗，何一帆，叶婧．基于稀疏近似逆预处理的牛顿．广义极小残余潮流计算方法［Ｊ］．电网技术，２００８，—３２（１４）：５０５３．—ＷＡＮＧＦａｎｇ－ｚｏｎｇ，ＨＥＹｉｆａｎ，ＹＥＪｉｎｇ．ＬｏａｄｆｌｏｗｃａｌｃｕｌａｔｉｏｎｏｆＮｅｗｔｏｎ－ＧＭＲＥＳｍｅｔｈｏｄｗｉｔｈｓｐａｒｓｅａｐｐｒｏｘｉｍａｔｅｉｎｖｅｒｓｅｐｒｅｃｏｎｄｉｔｉｏｎｅｒｓ［Ｊ】．ＰｏｗｅｒＳｙｓｔｅｍ—Ｔｅｃｈｎｏｌｏｇｙ，２００８，３２（１４）：５０５３．［２４］ＢｅｎｚｉＭｉｃｈｅｌｅ，ＴｕｍａＭｉｒｏｓｌａｖ．Ａｓｐａｒｓｅａｐｐｒｏｘｉｍａｔｅｉｎｖｅｒｓｅｐｒｅｃｏｎｄｉｔｉｏｎｅｒｆｏｒｎｏｎｓｙｍｍｅｔｒｉｃｌｉｎｅａｒｓｙｓｔｅｍｓ［Ｊ］．ＳＩＡＭＪｏｕｒｎａｌｏｎＳｃｉｅｎｔｉｆｉｃＣｏｍｐｕｔｉｎｇ，１９９８，１９（３）：９６８－９９４．［２５］ＯｗｅｎｓＪＤ，ＨｏｕｓｔｏｎＭ，ＬｕｅｂｋｅＤ，ｅｔａ１．ＧＰＵＣｏｍｐｕｔｉｎｇ［Ｊ］．ＰｒｏｃｅｅｄｉｎｇｓｏｆｔｈｅＩＥＥＥ，２００８，９６（５）：８７９。８９９．—收稿日期：２０１２－０４２３；—修回日期：２０１２－０５２９作者简介：江涵（１９８４－），男，博士研究生，主要从事电力系统稳定分析、电力系统高性能计算的研究；Ｅ．ｍａｉｌ：ｊｉａｎｇｈａｎ＠ｚｊｕ．ｅｄｕ．ｃａ江全元（１９７５－），男，教授，主要从事电力系统稳定与控制、电力系统高性能计算，电力系统最优化的研究。∞¨ｎｎ