一种简单实用的潮流计算程序设计

姚玉斌， 吴志良

(大连海事大学 轮机工程学院， 辽宁 大连 116026)

一种简单实用的潮流计算程序设计

姚玉斌， 吴志良

(大连海事大学 轮机工程学院， 辽宁 大连 116026)

潮流计算是电力系统分析的一项基本计算，也是其它电力系统分析的基础。编制并熟练掌握潮流计算程序不仅有利于深入理解潮流计算的本质，熟悉电力系统的运行状态，还可为研究其它电力系统分析课题打下基础。潮流计算算法难度适中，涉及算法设计和编程技巧却比较全面，是电力系统专业学生学习编制电力系统程序的理想案例。

潮流计算；牛顿法；导纳矩阵；雅可比矩阵；高斯消元法

0 引言

电力系统潮流计算是研究电力系统稳态运行的一项基本计算，根据计算结果，可以评价电力系统现行运行方式或计划运行方式的合理性和经济性；潮流计算也是其它电力系统分析的基础，如安全分析、暂态稳定分析等都要用到潮流计算。所以，潮流计算是电力系统分析中使用得最广泛、最基本和最重要的电力系统分析工具。

潮流计算程序也是一种较简单的电力系统计算程序。通过编制潮流程序可以更好地掌握潮流计算的基本原理，也为编制其它电力系统软件打下基础。

众多不同的潮流计算算法中，牛顿-拉夫逊法(也称牛顿法)及其改进方法(快速分解法)由于具有较满意的收敛性、占用内存适中、计算速度较快的优点，得到了广泛的应用[1-3]。电力公司使用的商用潮流计算软件都使用稀疏矩阵和节点优化编号等技术，编程难度较大，不适合初学者学习；而有些学生设计的潮流计算又过于简单，使用不灵活，如要求PV节点的编号一定编在后面等。为了降低编程难度、突出潮流计算的物理概念，又不失实用性，文献[4]设计了一种不使用稀疏矩阵技术和节点优化编号技术的快速潮流计算方法。本文在文献[4]介绍的方法的基础上进一步完善,提出一种实用潮流计算算法。该算法采用C语言编程，流程清晰，编程简单，是电力系统专业学生学习编制电力系统程序的理想方案。测试结果表明本程序功能齐备、速度较快，完全能够满足研究和在线应用。

1 极坐标牛顿法潮流计算原理

采用计算机算法分析各种实际问题，首先要理清要分析问题的流程。潮流计算实质上是一组非线性方程的求解问题，求解非线性方程组最常用最有效的方法就是牛顿法。

电力系统各节点的功率方程为

(1)

式中，n为节点数，G和B是导纳矩阵元素的实部和虚部。

通过式(1)求出的各节点功率与注入功率应该相等。如果计算电压值与真实电压有偏差，则节点计算功率与注入功率就不等，称为不平衡量，为

(2)

式中，Pis为节点i的有功注入，即节点i的有功出力和有功负荷之差。

令式(2)等于0就得到了潮流计算的一组非线性方程。由于不是所有节点的功率都是给定的，需要求解的方程个数要小于2n，为n-1+m个，m为PQ节点数，即n-1个有功不平衡量方程和m个无功不平衡量方程。

对式(2)线性化，得到潮流计算的修正方程为

(3)

式中，J为雅可比矩阵。

求出电压幅值和相角修正量后，用下式对节点电压进行修正:

(4)

上述潮流计算中，所涉及的电压用幅值和相角表示，则对应的潮流计算方法被称为极坐标牛顿法潮流计算。

2 潮流计算流程

根据极坐标牛顿法潮流计算的原理可以写出潮流计算流程图如图1所示。

图1 极坐标牛顿法潮流计算的流程图

从图1可见，按模块化设计极坐标牛顿法潮流计算算法应该主要包括：输入原始数据模块、设置电压初值模块、形成节点导纳矩阵模块、计算节点功率模块、计算功率不平衡量模块、形成雅可比矩阵模块、方程求解模块、计算支路功率模块。在这些模块中，形成导纳矩阵模块、形成雅可比矩阵模块、方程求解模块、计算支路功率模块是程序设计难点。另外对不熟悉程序设计的学生来说，利用数据文件输入数据也是比较陌生和困难的。

3 数据结构

著名计算机科学家沃思提出“数据结构+算法=程序”的见解，可见数据结构的设计对程序设计是很重要的[5]。潮流计算涉及的原始数据有支路(包括输电线路和变压器支路)参数、节点类型、功率注入数据、无功补偿参数等。为了清晰起见，采用结构体设计这些数据的数据结构，如表1～表3所示。三者分别为支路参数、节点数据、无功补偿设备的数据结构。

表1 支路参数数据结构

表1对应的变压器模型如图2所示，输入数据时，节点号前应该添一个负号，以区分输电线路和变压器支路，以及判断哪一侧是非标准变比侧。

图2 变压器模型

表3 无功补偿参数数据结构

电压的幅值和相角、导纳矩阵元素的实部和虚部都可设计为结构体。当然表1～表3表示的结构体数据结构也可以用普通变量的数组表示。

潮流计算需要定义一些变量来表示节点数、支路数、无功补偿数、收敛精度和最大迭代次数等。

4 主要模块设计

潮流计算采用模块化设计，涉及很多模块。主要有：

4.1 输入原始数据模块

首先用FILE定义文件指针，指向数据文件。

FILE *fp1, *fp2;

再用fopen函数打开数据文件，指定读写属性。

fp1 = fopen(“data.txt”, “r”);

fp2 = fopen(“flowout.txt”, “w”);

以上2个语句把相关数据文件赋予文件指针，并规定文件的读写属性，读写函数使用文件指针对数据文件进行操作。

用fscanf函数读数据，读整型数(int)用%d格式，读单精度浮点数(float)用%f，读双精度浮点数(double)用%lf，注意这里不是数字1而是字母l。变量和格式要一一对应，变量前需要取地址符(&)。

fscanf(fp1, “%d %lf”, &n, &eps);

用fprintf函数写数据，写整型数用%d格式，可以指定数据宽度，如下例中“5”；写浮点数用%f，可以指定数据宽度和小数位数，如下例中“7.4”表示数据宽度为7，其中小数位为4。注意变量前没有取地址符(&)。

fprintf(fp2, “%5d %7.4f”, n, eps);

4.2 形成导纳矩阵模块

牛顿法潮流计算是基于导纳矩阵的方法，因此需要先根据支路参数和无功补偿参数形成导纳矩阵。

手工计算导纳矩阵元素一般都是按节点进行，一次计算一个节点的自导纳和相关互导纳。计算机计算导纳矩阵的过程也可以采用手工计算的过程，但编程较复杂，计算效率低。一般采用按支路计算，每扫描一条支路，计算出支路两端节点间的互导纳，以及这两个节点自导纳的一部分再与上次值累加起来，待所有支路都扫描完，完整导纳矩阵就形成了。

另外，还要根据无功补偿参数追加相应节点的导纳矩阵对角元的虚部。

对多回支路，可以在输入数据后处理，进行合并，但需要增加模块，增加计算量，也无法计算各并联支路的支路功率。其实在本模块稍微改动即可解决多回支路问题，只需将互导纳计算也改为累加方式即可。

4.3 形成雅可比矩阵模块

雅可比矩阵是非线性潮流计算方程线性化得到的线性方程组的系数矩阵，其元素是功率不平衡量对电压幅值和相角增量的偏导。

雅可比矩阵是一个(n-1+m) ×(n-1+m)阶矩阵，很多程序设计规定，PQ节点的编号为前m个，平衡节点放在最后。然而很多情况下电力系统数据的节点编号不满足这一要求，这样就需要增加一个重新编号模块，各节点在计算过程需要换号，大大增加编程难度，也容易出错。文献[4]设计了一种新的处理方法，把雅可比矩阵设计为 (2n) ×(2n)阶矩阵，每个节点对应2行2列，PV节点的ΔQ对应的行及ΔU对应的列空着不用，置为0。这样对节点编号就没有特殊要求，也不用重新编号；还能处理多平衡节点系统，当然形成雅可比矩阵和解方程时要特殊设计。

此模块涉及的公式较多，逻辑判断也较多，需要分析一下。

(1) 雅可比矩阵元素，利用下式关系编程，可以减少计算量，并使得程序清晰，减少出错概率。

(5)

(2) 节点功率与雅可比矩阵元素存在下式关系，利用此关系可以顺便求出节点功率，不需要再单独计算功率，减少计算量。

(6)

(3) 由于节点较多时，导纳矩阵比较稀疏，即零元素很多，可以判断导纳是否为0，以减少计算量。

4.4 方程求解模块

线性方程组一般采用消去法求解，如高斯消去法、约当消去法，其中高斯消去法计算量较小。

1) 高斯消去法原理

对于一般形式的方程组AX=F，阶数为p，消去运算显然只对A和F进行，与未知向量X无关。为了叙述和程序设计的方便，令F的元素fi=ai,p+1作为第p+1列添加在A之后，构成增广矩阵B。

用高斯消去对矩阵B进行消去后，得到如下矩阵，其中系数矩阵是单位上三角阵：

(7)

2) 算法设计

高斯消去法包括消去和回代两个过程，消去又可分为按行消去和按列消去，回代也可分为按行回代和按列回代。本文采用按行消去和按行回代。

对第i行消去过程公式如下：

(8)

式中，上标(k)表示第k次迭代的值。

回代过程公式如下：

(9)

回代后，增广矩阵最后1列元素就是方程的解。

3) 算法分析

根据高斯消去法原理可知，消去后，式(8)的老元素就没有用了，因此可以用新元素更新老元素。

由式(9)可见，回代时没用到系数矩阵下三角和对角元，因此按式(8)编程，消去后系数矩阵的下三角并不为0，对角元也不为1。如果循环时，消去过程控制变量j从k开始，反而会出错。

由于雅可比矩阵有一些行的元素全部为0，模块中应通过判断对角元跳过这些行。

4.5 迭代模块

潮流计算是一个非线性问题，牛顿法采用逐次线性化的方法求解，是一个反复迭代的过程。最大不平衡量小于收敛精度，即为收敛，停止迭代。为了避免死循环，迭代次数大于最大迭代次数时，终止迭代。

设计迭代模块不难，但调试它有难度，因为如果出现反复迭代，计算结果不正确，就很难确定是哪一步有问题。

4.6 计算支路功率模块

计算支路功率模块，要计算支路首节点到末节点的功率，也要计算末节点到首节点的功率。

计算出电压后，就可以计算支路电流，进一步就可以算出支路功率。

(10)

(11)

式中，b0为支路对地电纳，上标“^”表示共轭。

由于理想变压器没有损耗，所以流入变压器的功率就是流经支路阻抗zij的功率。如图2只要计算出i′点的电压，再按没有对地支路的线路计算即可。

(12)

4.7 其它模块

(1) 电压初值设置模块：潮流计算的电压初值一般采用平启动，即PV节点和平衡节点的电压幅值设为给定值，PQ节点的电压幅值设为1，所有节点的电压相角都设为0。如果考虑多平衡节点，有些平衡节点的相角可能不为0，此时应对平衡节点另行设置。

(2) 节点功率计算模块：节点功率计算可以在计算雅可比矩阵时顺便计算，这样就不需要节点功率计算模块了。

(3)功率不平衡量计算模块：计算功率不平衡量要注意平衡节点的有功和无功不平衡量以及PV节点的无功不平衡量要置0。

(4)节点功率计算模块：先由式(1)计算出待求节点的计算功率，再根据计算功率等于出力与负荷之差就能求出平衡节点的有功和无功出力、PV节点的无功出力。

5 程序调试建议

潮流计算涉及迭代过程，数据量又比较大，因此调试比较困难。有以下建议供参考。

5.1 语法错误

(1) 由于编译器会提示错误信息，语法错误比较容易发现和改正。但注意有些错误不一定在编译器指出错误的行上，可能前1行错误引起的，因此查找错误不要只在指定行查找。

(2) C语言使用的逗号、分号都是英文标点，不能写成中文标点。因二者很相近，容易混淆，应特别注意。

5.2 输入数据错误

原始数据读入的错误是导致计算错误的一个主要因素。原因有：① 数据输入模块编写有误；② 数据文件录入的数据有错误；③ 数据文件的文件名不对；④ 数据文件不存在或不在指定文件夹。读入数据一定要再写到输出文件以检查输入的数据是否正确。

5.3 逻辑错误

(1) 避免使用字母“l”、“O”，因它们与数字“1”、“0”很像，容易混淆，出错后不易发现。

(2) 采用输出信息方式帮助调试，按图1所示流程依次检查各模块的输出。信息较少时可以输出到显示器上，大量信息则写到文件中。

(3) 对数据分析，要掌握技巧。如导纳矩阵，分析互导纳时，先判断支路类型，以确定是输电线路还是变压器支路的流程有误；分析自导纳时，先看只连输电线路的节点，再看仅连变压器支路的节点，最后看既有输电线路也有变压器支路的节点。

(4) 雅可比矩阵元素分“i=j”和“i≠j”两种情况检查，并且检查出错的是H、N、M、L中哪一类。要检查头两次迭代过程中的雅可比矩阵，仅第1次迭代的雅可比正确并不能保证以后迭代过程的雅可比矩阵正确，因为首次迭代时，电压相位为0，sinθij= 0,导致有些错误不能暴露。

6 算例及性能分析

6.1 考核算例

根据文献[1]的算例，进行修改得到图3所示的考核算例，其具有以下特点：① 节点编号随意，没有把PQ节点都排在最前面，平衡节点也不是最后一个节点；② 存在双回线；③ 存在无功补偿设备。

图3 考核算例系统图

6.2 性能分析

采用本文方法设计的潮流计算程序对某445节点规模的电网进行了分析。本算例是在主频为2.27 GHz的Intel Core i3的PC机上进行的，迭代次数为5，计算时间为0.2 s，速度很快，完全满足研究需要。

本程序对节点编号没有特殊要求，节点可以任意编号；能够处理多回线；也可以处理无功补偿设备；还能够处理多平衡节点系统。

7 结语

本文设计不采用稀疏矩阵和节点优化编号技术，使得潮流计算程序简单明了，大大降低了编程难度，更利于理解潮流计算的原理，有利于初学者学习电力系统软件编程和加深对电力系统运行状态的理解。

(姚玉斌等文)

本方法能够考虑多回线、无功补偿设备、多平衡节点的情况，功能齐备，计算速度较快，完全能够满足研究使用，也能满足在线应用。

[1] 西安交通大学，等. 电力系统计算[M]. 北京：水利电力出版社，1978

[2] 何仰赞，温增银，汪馥瑛，周勤慧. 电力系统分析(第3版) [M]. 武汉：华中理工大学出版社，1996

[3] 陈珩. 电力系统稳态分析(第3版) [M]. 北京：中国电力出版社，2007

[4] Yubin Yao, Benhua Gu, Hong Zhao, et al. A Newton Power Flow for Study Purpose[C], Proceedings of the 2012FEEM，Energy Procedia 16, 1372-1379

[5] 谭浩强. C语言设计(第3版) [M]. 北京：清华大学出版社，2005

Design of a Simple and Practical Power Flow Algorithm

YAO Yu-bin, WU Zhi-liang

(SchoolofMarineEngineeringofDalianMaritimeUniversity,Dalian116026,China)

Power flow is a basic tool for power system static analysis, and also base for other power flow analysis. By programming and mastering power flow algorithm, students can not only deeply comprehend the essence of power flow and be familiar with operational states of power system, but also lay foundation for other power flow analysis. Power flow algorithm is ideal case for students to study power system programming, which has moderate level of difficulty and needs comprehensive skill for algorithm designing and programming.

power flow; Newton method; admittance matrixes; Jacob matrixes; Gauss elimination.

2016-04-22;

2016-07- 23

姚玉斌(1965-)，男，教授，博士，主要从事电力系统分析和控制的教学和研究工作，E-mail：yaoyubin@dlmu.edu.cn 吴志良(1963-) ，男，教授，博士，主要从事船舶自动化系统的教学与研究工作。

TM744

A

1008-0686(2017)00-0074-06