电力系统稳定计算软件的开发

摘要

随着世界各国的电力系统发展的越来越庞大，电力系统的运行也随之越来越复杂，系统中发生的事故越来越难以用传统的分析方法预测，对电力系统仿真技术提出了更高的要求。目前，国际上有多种电力系统分析软件，在各国的电力系统的实践和研究中发挥了很大的作用。

本文研究了基本的电力系统暂态分析方法。开发的程序可以计算计及励磁系统和调试系统的包含10台发电机的小电力系统的机电暂态过程。

电力系统的暂态过程分析和电力系统各元件的动态特性有密切的关系。因此本文的主要思路是在分析电力系统暂态稳定计算基本原理的基础上，建立电力系统各元件的数学模型，并根据元件之间的拓扑关系将这些模型组合为全系统模型，然后采用适当的数学方法来进行求解。在本文中，电力系统模型的分析求解采用的是隐式积分法。最后采用matlab 绘图得到各发电机端电压，电磁功率，机械功率，绝对功角，相对功角曲线，用来判断系统在扰动后是否能够恢复稳定状态。

本程序通过算例的验证，可以计算小系统的稳定计算。

关键词：隐式积分法电力系统模型网络化简

The development of Power system stability

calculation software

Abstract

Along with the development of power system, and also the operation is more and more complicated, the accident happened in the system has become more and more difficult to predict with traditional analysis method, which puts forward a higher request to the power system simulation technology. At present, there are many international power system analysis softwares, that play an important roal in the power system of practice and study.

In this paper we research the basic power system transient analysis method. This program can calculate the electrical transient process of small system. This system contains excitation system and debugging systems, also contains 10 generators.

Power system transient process analysis and power system dynamic characteristics of each component has close relationship. The main ideas of the paper are based on analysising the basic principle of the power system transient stability calculation. Firstly established the mathematical model of components. Then put them to a system model according to the network topology. Again, using appropriate mathematical method to solve the problem. In this paper,we can use implicit integral method . Finally using matlab drawing the electromagnetic power, voltage generator, mechanical power, absolute Angle, relative Angle curve to verify the stability of the system.

This program is verified by an example. It can solve the stability calculation of small systems.

Keywords: Implicit integration scheme;Power system model;Network simplification

摘要 . ..................................................................... Ⅰ Abtrsact ................................................................... Ⅱ

1 绪论 . .................................................................... 1

1.1 课题背景 . .............................................................. 1

1.2 电力系统仿真研究现状与发展 . ............................................ 2

1.2.1 仿真算法的选择 . ...................................................... 2

1.3 电力系统数字仿真的发展与展望 . .......................................... 3

1.4 本文任务 . .............................................................. 4

2 电力系统暂态稳定计算的基本原理与元件数学模型 . ............................ 5

2.1 基本原理 . .............................................................. 5

2.2 基本假设 . .............................................................. 6

2.3 电力系统中各元件的数学模型 . ............................................ 7

2.3.1 同步发电机组的数学模型 . .............................................. 7

2.3.1.1 同步发电机组转子运动方程 . .......................................... 7

2.3.1.2 发电机电动势变化方程 . .............................................. 7

2.3.1.3 同步发电机定子方程 . ................................................ 8

2.3.2 发电机励磁系统的数学模型 . ............................................ 8

2.3.3 负荷的数学模型 . ..................................................... 10

2.3.4 原动机的数学模型 . ................................................... 10

2.3.5 电力网络的数学模型 . ................................................. 12

2.3.6 机网衔接 . ........................................................... 12

3 解题过程中采用的数学方法 . ............................................... 14

3.1 牛顿—拉夫逊法潮流计算 . ............................................... 14

3.2 隐式积分法解微分方程组 . ............................................... 16

3.3 xy dq 坐标变换 ....................................................... 19

3.4 化简网络中的无源节点 . ................................................. 20

4 机电暂态仿真计算流程 . ................................................... 21 . ...........................................................................

5 算例验证及结果分析 . ..................................................... 25

6 结论 . ................................................................... 29

7 参考文献 . ............................................................... 30

致谢 . ..................................................................... 31

1 绪论

1.1 课题背景

1831年，法拉第发现电磁感应原理，奠定了发电机的理论基础。科学的发现，引起了技术的发明。1866年，维·西门子发明了励磁电机。接着，1876年，贝尔发明了电话；1879年，爱迪生发明了电灯。这三大发明开启了人类电气化的道路，引起电力技术革命。电力技术和电力工业的出现和发展改变了人们的生产和生活面貌，促使经济以前所未有的速度向前发展。

随着电子技术、电子计算机技术和自动化技术的发展，电力工业自动化迅速向前发展。以大机组、大电厂、高电压、大电网、高度自动化为特点的现代化电力工业在不同的国家已经形成或正在形成。

时至今日，世界各国的电力系统发展得越来越庞大，电力系统的运行也随之越来越复杂。由原来的简单的输电网络变成了先进的可以说是人类历史上构造最复杂的工业系统之一，它是一个强非线性、高维数、分层分布的动态大系统。一旦发生事故，对于国民经济造成的损失是十分巨大的。

由于现代电力系统具有这样的特点，系统中发生的事故越来越难以用传统的分析方法预测。在我国，随着三峡电站的建设、西电东送工程的实施和全国联网工程的推进，我国大陆部分在本世纪初叶将形成规模巨大的全国性的互联电网，而且，全国联网系统中既有交流线路，又有直流线路，还包括诸如无功静止补偿器之类的电力电子设备和FACTS 装置，从而使电网的安全稳定运行控制变得十分复杂。

为了保证电力系统可靠、安全、经济的运行，在规划、分析和研究电力系统时必须确切完整的考察实际电力系统的静态特性和动态特性。然而电力系统的特点决定了难以采用试验的方法来实现，必须采用仿真的手段。

电力系统数字仿真是对电力系统进行分析和研究的重要手段，不受原型系统规模和结构复杂性的限制。能够保证被研究系统的安全性，并具有良好的经济性、方便性等诸多优点，已成为电力系统规划、运行、设计以及安全稳定分析等不可或缺的工具。按照研究对象、模型和算法的不同，电力系统数字仿真又可以分为机电暂态仿真和电磁暂态仿真。机电暂态仿真主要研究电力系统受到大扰动后的暂态稳定性和受到小扰动后的静态稳定性，而电磁暂态仿真侧重研究电力系统元件中电场和磁场以及相应的电压和电流的变化情况。本文要研究的主要是电力系统机电暂态仿真[1]。机电暂态仿真主要研究电力系统受到大扰动后的暂态稳定和受到小扰动后的静态稳定性能，需要联立求解电力系统微分方程组和代数方程组，以获得物理量的时域解。现在的机电暂态仿真软件对仿真规模一般没有限制，对称分量法的引入，使其也可以计算系统的不对称故障。因此在实际工程中特别是对于大型电力系统的稳定研究中，得到了广泛的应用。

1.2 电力系统仿真研究现状与发展

电力系统数字仿真的基本要求之一是所仿真系统与实际系统的结果尽可能一致，而两种情况的吻合程度在很大程度上取决于仿真所采用模型的准确性[2]。在满足数字仿真精度的同时，要实现数字仿真的快速性，尽可能缩短仿真所花费的时间，这就要求对仿真的模型作适当地简化。因此，在不同的仿真要求和条件下，同时兼顾仿真的快速性和准确性，就必须对同一装置或设备采用多种不同复杂程度的电力系统数学模型。

在本文中，发电机一般采用dq0坐标下的基本方程，即派克方程；负荷及其他系统元件采用以额定频率同步旋转的直角坐标，即xy0坐标。发电机的数学模型经坐标变换与网络方程衔接即可进行全系统方程组的联立求解。

在确定所要仿真系统的模型后，就要着手仿真算法的研究。仿真算法类型繁多，不同算法的适应性不同，需根据具体的仿真系统与仿真要求地选用。一般应依据以下原则：(1)算法的数值稳定性要高；(2)算法的健壮性要好（所谓健壮性，指一个算法除了能对合法的输入数据得到正确的结果外，还应对非法的输入数据做出正确合理的处理，健壮性体现了思维的缜密）；(3)算法要简单且收敛性好。

1.2.1 仿真算法的选择

早期的暂态稳定计算的程序大都采用显式数值积分法，包括欧拉法、改进欧拉法和龙格—库塔法[3]。欧拉法和改进欧拉法由于计算精度较低且数值稳定性差，通常只适用于当所采用的元件数学模型比较简单而且所要求计算的暂态过程持续时间不长时的暂态稳定计算中。龙格—库塔法的计算精度优于改进欧拉法，所采用的数学模型可以较为详细，目前仍应用于一些暂态稳定计算中。

目前仍在使用的一些仿真算法有：(1)单步法如欧拉法与改进欧拉法、龙格—库塔法；

(2)线性多步法如哈明算法，阿达姆斯算法；(3)针对病态方程的Gear 法、向后差分法、隐式及半隐式龙格——库塔法等。单步法的优点是可以自起步，即由稳态值出发，用通用公式一步步求取数值解，而多步法则不能自起步，因而需要多于一步的状态信息。电力系统仿真计算一般采用单步法。

目前在电力系统数字仿真软件中广泛采用的是隐式梯形积分法，它具有数值稳定性好，对刚性微分方程的适应性好，可采用大步长，并可通过联立求解同步发电机的微分方程和网络的代数方程来消除接口误差等优点。

而且由于隐式积分法具有较好的数据计算特性，因而在用在电力系统稳定计算时，可任意采用较详细的元件数学模型。对于发电机，采用考虑阻尼绕组时的数学模型，并考虑励磁系统、原动机及其调速系统的影响；对于负荷，考虑其动态和静态特性模型等。

1.3 电力系统数字仿真的发展与展望

现代化的电网迫切需要现代化的分析手段，随着全国电网及电力市场化的发展，越来越多新技术的投入使用，计算机技术的日新月异，新的数值计算方法的出现，电力系统数字仿真也将会有新的发展和变化[2]。

1 计算机技术的发展

计算机硬件技术的发展使得计算机对运行环境的要求更低，使得硬件的维护、升级更方便，使得软件的可移植性更好。

计算机软件技术的发展推动电力系统仿真的前进。以往电力系统仿真软件的开发基本上是结构化的理论体系，即采用面向功能的设计思路，该方法使系统组织方式、用户界面和计算功能的组织是专用的、固定的，不同系统之间缺乏共同性和相互综合的可能性。因此系统的通用性和可维护性较差。当系统要进行扩展或修改时，修改的工作量较大，容易导致系统冗余增加，结构变得复杂而且利用率下降。面向对象技术是解决这一问题的有效途径。程序中类、对象、继承、成员、结构、封装等机制的引入为进一步改善模块的结构和管理，以及模块之间的信息传递提供了便利，其封装性的特点可以保证系统某一功能模块不断扩展和维护，但仍然有许多问题需要进一步研究和探讨。另外，多媒体技术的发展和应用，为开发友好的人机界面提供了可能。

2 新元件的建模及原有模型的改进

随着各种新型电力元件的出现，需要增添新的元件模型。如柔性交流输电系统装置的应用在高压互联电网中已成为必然的趋势，因此需建立相应的动态元件模型，以分析其对系统的作用。

3 模型定义方式的改进

采用模块化的建模思路，按照标准模块搭建系统模型。另外，允许用户自定义模块和采用动态模型库也是一种有效的方法。这些技术使模型的修改很方便，而且对新出现的电力元件也可容易地建立模型。

4 数字仿真数值方法的新变化

除了现阶段提出的经典数值积分方法外，正在发展一些新算法。如现代数值仿真算法，主要是离散近似法、可调整数值积分法、增广矩阵法、采样控制系统仿真法和面向结构图仿真法等。无论使用哪种仿真算法，都是按环节离散化或按环节使用数值积分公式，无法自动推出整个系统的状态方程。这样，对于整个系统的仿真误差，已不能用截断误差的方法进行分析。

5 统一的电磁暂态与机电暂态混合仿真

机电暂态仿真与电磁暂态仿真在仿真步长、模型变量的表示、仿真规模和应用领域等方面都存在差异。电力系统发生故障或操作后，将产生复杂的电磁暂态过程和机电暂态过程，两者同时发生并相互影响。由于这两个暂态过程的变化速度相差很大，通常近似地对它们分别进行仿真。随着电力电子设备和高压直流输电技术的广泛使用，机电暂态过程仿

真中使用这些设备的准稳态模型或简化模型，会导致仿真结果存在较大误差；另一方面，虽然电磁暂态仿真能够准确表达电力电子设备模型，由于受模型与算法的限制，其仿真规模不大，难以适应于现代大电力系统。电磁暂态过程与机电暂态过程的混合仿真可兼得各自的优点。混合仿真的关键在于两仿真程序之间的接口问题。

1.4 本文任务

限于能力与时间的限制，本文对电力系统分析的研究将仅限于小系统机电暂态的过程，建立一个数学模型来解决小系统机电暂态的问题，最终形成的程序可以计算计及励磁系统和调试系统的包含10台发电机的小电力系统的机电暂态过程。。

2 电力系统暂态稳定计算的基本原理

与元件数学模型

2.1 基本原理

电力系统运行稳定性的问题就是当系统在某一正常运行状态下受到某种干扰后，能否经过一段时间后回到原来的运行状态或者过渡到一个新的稳定运行状态的问题[4]。如果能够，则认为系统在该正常运行状态下是稳定的。反之，若系统不能回到原来的运行状态或者不能建立一个新的稳定运行状态，则说明系统的状态变量没有一个稳态值，而是随着时间不断增大或振荡，系统是不稳定的。通常所考虑的扰动包括各种短路故障、切除大容量发电机或输电设备以及某些负荷的突然变化等。

在正常的稳定运行情况下，电力系统中各发电机组输出的电磁功率和原动机的机械功率保持一致，因此所有发电机的转子速度保持恒定。但是当系统受到较大的扰动后，系统中将经历暂态过程。即由于系统的结构或参数发生了较大的变化，使得系统的潮流和各发电机的输出功率发生变化，破坏电磁功率与机械功率的平衡，在发电机转轴上产生不平衡转矩，导致转子加速或者减速。一般来说，在多发电机的复杂电力系统中发生大扰动时，各发电机输出的电磁功率将按扰动后网络特性重新分配。这样有些发电机因电磁功率小于原动机机械功率而加速，有的则因为发电机电磁功率大于原动机机械功率而减速，从而使各发电机之间的相对角度发生较大的变化。这种角度的变化又反过来影响各发电机的输出功率，从而使各个发电机的功率，转速和转子之间的相对角度继续发生变化。

同时，由于发电机端电压和定子电流的变化，也将引起励磁系统的调节作用；由于机组转速的变化，引起发电机调速系统的调节过程。这些调节过程也将影响发电机的功率平衡状况。

电力系统遭受大干扰后发生的机电暂态过程可能有两种不同的结局。一种是各发电机转子之间的相对角随时间的变化呈振荡的状态，且振幅值逐渐衰减，各发电机之间的相对运动将逐渐消失，从而系统过渡到一个新的稳定运行状态，各发电机仍然保持同步运行。这时，我们称电力系统是暂态稳定的。另一种结局是在暂态过程中某些发电机转子之间始终存在着相对运动，使得转子间的相对角度随时间不断增大，最终失步。严重情况下甚至可能导致系统的解列或崩溃。

电力系统正常运行的必要条件是所有发电机保持同步。因此，电力系统在大干扰下的稳定性分析，就是分析遭到大干扰后系统中各发电机能否维持同步运行的能力，这就是电力系统稳定性分析。

判断电力系统是否具有暂态稳定性的具体方法就是看各发电机转子之间相对角的变化特性。在相对角中，只要有一个相对角随时间的变化趋势是不断增大（或不断减小）时，系统就是不稳定的；如果所有的相对角经过振荡之后都能稳定在某一值，则系统是稳定的。

因为绝对角是发电机相对于同步旋转轴的角度，因此，若绝对角δi 随时间不断增大，则意

味着第i 台发电机的转速高于同步转速；若若δi 随时间不断减小，则意味着第i 台发电机的转速低于同步转速。当所有发电机的绝对角都随时间不断增大时，系统仍然可能是稳定的，它只意味着在新的稳定运行状态下，系统频率高于额定值。

判断电力系统在预想事故下能否稳定运行，需要进行暂态稳定分析。当系统不稳定时，还需要研究提高系统稳定性的有效措施；当系统发生重大稳定破坏事故时，需要进行事故分析，找出系统的薄弱环节，并提出相应的对策。

2.2 基本假设

由大扰动引起的电力系统暂态过程，是一个电磁暂态过程和发电机转子间机电暂态过程交织在一起的复杂过程。当计及原动机调速器、发电机励磁系统等调节设备的动态过程时，暂态过程将更加的复杂[3]。

精确地确定所有的电磁参数和机械运动参数在暂态过程中的变化时困难的，对于解决一般的工程实际问题也是不必要的。通常，暂态分析计算的目的在于确定系统在给定的大干扰下发电机能否保持同步运行。因此，只需要考虑暂态过程中对转子机械运动起主要影响的因素，对于影响不大的因素，则予以忽略或近似考虑。

1、忽略发电机定子电流的非周期分量和与它相对应的转子电流的周期分量

在大扰动中，定子的非周期分量电流将在定子回路电阻中产生有功损耗，增加发电机转轴上的电磁功率，在某些情况下，可能使发电机的电磁功率大于原动机的机械功率，从而使发电机减速运动。然而，由于定子非周期分量电流衰减时间常数很小，而且转子的机械惯性较大，因而对转子整体相对运动的影响很小。

采用这一假设即忽略电力网络各元件的电磁暂态过程。此时，定子绕组中只包含基频分量，定子电压方程就变成了代数方程。

2、不计零序和负序电流对转子运动的影响

对于零序电流来说，一方面，由于发电机接在升压变压器的三角形侧，故障发生在高压侧时，零序电流并不通过发电机；另一方面，即使发电机流通零序电流，由于定子三相绕组在空间对称分布，零序电流对转子运动也没有影响。

对于负序电流，其与转子绕组直流电流相互作用所产生的转矩是以近二倍同步频率交变的转矩，其平均值接近于零，对转子运动的总趋势影响很小。

采用这个假设时，发电机输出的电磁功率仅由正序分量确定。故障时确定正序分量的等效电路与正常运行时的等效电路的不同之处在于故障处接入由故障类型确定的附加电抗不同。

3、忽略暂态过程中发电机的附加损耗

这些附加损耗对转子的加速运动有一定的制动作用，但其数值不大。

4、不考虑频率变化对系统参数的影响

在一般的暂态过程中，发电机的转速偏离同步转速不多，可以不考虑频率变化对系统参数的影响，各元件参数值都按额定频率计算。

2.3 电力系统中各元件的数学模型

在电力系统分析中，各元件所采用的数学模型，不但与稳定分析结果的正确性直接相关，而且对稳定分析的复杂性的很大的影响。因此，选用适当的数学模型描述各元件的特性，使得稳定分析的结果满足合理的精度要求并且计算简单，是电力系统稳定分析中一个至关重要的问题。

2.3.1 同步发电机组的数学模型

为了便于对电力系统的稳定性问题进行准确的分析和计算，必须首先建立描述发电机转子运动的方程[3] [4] [5]。

2.3.1.1 同步发电机组转子运动方程

根据资料[3]与资料[4]可得发电机转子运动方程的标幺值形式为：

⎧d δ

=(ω-1) ω0

⎪dt ⎪

⎨d ω1

⎪=(P T -P e -D (ω-1)) ⎪⎩dt T J

(2-1)

其中： T J ：发电机组惯性时间常数 P m ：发电机机械功率标幺值 P e ：发电机电磁功率标幺值

δ：功角（发电机电动势响亮与同步转轴的夹角） ω0：同步转轴的角速度，其标幺值为1 ω：发电机的角速度

此方程组反应了功角δ在发电机同步转速发生变化时随时间变化的规律：当发电机机械功率大于电磁功率时，发电机加速，功角变大。

2.3.1.2 发电机电动势变化方程

虽然发电机中装有自动励磁系统，但是由于励磁调节系统和发电机转子回路的暂态过程，暂态电动势E q ' 和空载电动势E q 都不会一直保持一个常数。若在分析系统的稳定性时要求计及电动势的实际变化，则必须分析电动势的变化过程[3] [4] [6]。

同步电机的暂态电势E q ' 和E d ' 与转子绕组的暂态电磁过程有密切关系，描述E q ' 和

E d ' 电势变化规律的微分方程如下[4] [5]：

dE q ' ⎧

E =E +T ' ⎪qe q d 0

dt ⎨

⎪E q =E q ' +I d (x d -x d ') ⎩其中： E qe ：强制空载电动势

(2-2)

T d 0' ：励磁绕组时间常数，也可表示为T f

上式表明电动势E q 和E q ' 是受控于强制空载电动势E qe 即励磁电压u f 的，而励磁电压

u f 则是由自动调节励磁系统控制的。 2.3.1.3 同步发电机定子方程

电力系统暂态分析中发电机的模型一般用微分方程的阶次来说明，常用的有三阶模型和五阶模型[5]。三阶模型包括转子运动的两阶微分方程和励磁绕组电磁暂态的一阶微分方程。三阶模型再加上纵轴和横轴阻尼绕组各一阶微分方程，便构成五阶模型。

本文中采用的是三阶模型。即发电机用E q ' 表示时，按照固定在本机转子上的d,q 坐标系建立定子回路电压平衡方程，即：

⎧⎪E q ' =U Gq +I Gq r G +I Gd X d '

⎨

0=U -I X +I r ⎪Gd Gq q Gd G ⎩

(2-3)

由于发电机中的电压电流量是基于dq 坐标轴，而电网计算中的电压电流量是基于xy 坐标轴，因此需要把dq 坐标轴的量转换为xy 坐标轴的量，其中xy 坐标轴中常取稳态潮流计算中的平衡节点的电压相位为零作为参考，将其取做x 轴，具体的转换方法见下文。

2.3.2 发电机励磁系统的数学模型

发电机的自动调节励磁系统调节发电机励磁绕组两端的励磁电压u f ，因而影响发电机的电动势，对发电机的电磁功率和系统的稳定性均有较大的影响。因而在暂态稳定计算中有必要考虑励磁系统的作用[4] [5] [6]。

自动励磁调节器的作用是接收、处理和放大输入的控制信号，以形成合适的励磁调节。控制信号有发电机电压、电流等。自动励磁调节器中一般包括的测量、放大、稳定作用以及限幅等环节。自动励磁调节器的种类很多，目前多采用晶闸管和微机控制系统。

晶闸管励磁调节器传递函数的框图[5]如下所示：

图 2-1 励磁系统传递函数框图

其主要参数的取值范围如表2—1[5]所示：

表2—1 晶闸管励磁调节系统主要参数的数值范围

在近似的简化分析中，可以用一个等值的一阶惯性环节来表示励磁系统，如下所示：

图 2-2 励磁系统简化传递函数框图

其表达式为：

-∆U G 1+K T e p =∆u f

在上式中，负号表示端电压下降时励磁电压增加，K e 和T e 为等效的放大倍数和时间常数。对于闭路调节回路，包括励磁机的常规励磁系统，由于自动励磁调节器的滤波、综合

放大环节的时间常数较小，T e ≈T ff 。

由于励磁电压u f 和强制空载电动势E qe 间为线性关系，在标幺制中若取它们的基准值满足比例关系，则u f 和E qe 的标幺值相等，于是，励磁系统的表达式可以改写为：

d ∆E qe T e =-∆E qe -K e ∆U G

dt (2-4)

其中

∆E qe =E qe -E qe 0；∆U G =U G -U G 0

上述为正常励磁时描述强制空载电动势E qe 变化规律的微分方程。

当发电机机端电压由于系统发生短路等故障而大幅下降时，启用强行励磁，励磁电压

立刻升值其极限值u f max ，则应有：

du f T e =u f max -u f

相应的强行励磁时描述强制空载电动势E qe 变化规律的微分方程为：

dE qe T e =E qe max -E qe

dt 用。

(2-5)

以上两式即可描述发电机强制空载电动势的变化规律，表明了发电机励磁系统的作

2.3.3 负荷的数学模型

负荷是电力系统的重要组成部分，在电力系统暂态分析中，需要知道反映某一节点的全部负荷性能的数学模型[3]。

负荷特性可以分为静态特性和动态特性，前者指电压或频率进入稳态时负荷功率与电压或频率的关系，后者指电压或频率变化过程中负荷功率与电压或频率的关系。系统电压和频率的变化直接影响负荷从同吸收的有功与无功功率，后者又反过来影响系统的状态，由此可见，负荷特性对系统的稳定性有相当大的影响。

较为常用的负荷模型有三种:恒阻抗模型、综合负荷静态模型和综合负荷动态模型。由于综合负荷由各种不同种类的负荷组成，其组成情况不尽相同且随时变化，因此要准确的获得负荷的数学模型是比较困难的。本文只考虑负荷的恒阻抗模型来简化计算过程。

负荷的等效阻抗数值由潮流计算的结果得到。若负荷用电阻，电抗的并联支路来模拟，则：

其中： P ：节点有功功率 Q ：节点无功功率 R L ：负荷的等效电阻 X L ：负荷的等效电抗 U ：负荷节点电压

这是最简单的负荷处理办法，实际上就是近似认为负荷从系统吸收的功率总是正比于负荷节点电压的平方。采用这种方法在仿真时处理简单，计算速度快，精度较差。

P =U 2/R L Q =U /X L

2.3.4 原动机的数学模型

电力系统中提供机械功率或机械能的装置一般统称为原动机。目前电力系统中的原动机主要是水轮机和汽轮机[5] [7]。为了控制原动机想发电机输出的机械功率，并保持电网的正常运行频率，以及在并列运行的发电机之间更加合理的分配负荷，基本上所有的原动机都配备了调速器系统。通过改变汽轮机的汽门开度或水轮机的导水叶开度来实现原动机输出功率的调节。

通过调速系统对原动机输出机械功率的调节来实现对整个电力系统稳定性的调节作用。

在我国的电力系统中，水电和火电机组的原动机配置了不同类型的调速器。大型汽轮机普遍采用了高速弹簧片和旋转阻尼液压调速器。大型水轮机除了其中一部分仍然采用离心飞摆机械调速器外，近年来有较多的机组配置了电气液压调速器。为了适应大型中间再热汽轮机调节的要求，功率-频率电气液压调速器（简称功频电液调速器）的应用正在发展中。

从电力系统计算的角度出发，按照它们的作用原理，把各种调速器分为水轮机调速器

和汽轮机调速器两种基本类型。其中汽轮机调速器又可以分为液压调速器和功频电液调速器两种。本文重点介绍的是汽轮机液压调速器。

汽轮机液压调速器有旋转阻尼液压调速器和高速弹簧片液压调速器两种类型。旋转阻尼液压调速器的结构图如（2-3）所示：

其调节过程简述如下：当外界负荷增加时，汽轮机的转速降低。此时测速部件——旋转阻尼油压室的一次脉动油压下降，通过油压放大器的作用，是二次脉动油压上升，继动器活塞下移，关小碟阀下面的排油口，使得错油门克服下端弹簧的作用而下移，于是压力油进入油动机的下腔室，油动机活塞上移，开大调速气门，增加汽轮机的进气量，是汽轮机的出力调高并与外界负荷的增加相平衡，转速得以回升。在油动机活塞上移的同时，通过杠杆的作用，反馈弹簧将继动器活塞上拉，与转速回升的作用一起，使错油门回到中间位置，于是调节过程结束。外界负荷减少时，其调节过程与上相反。

图2-3 旋转阻尼液压调速器结构图

汽轮机液压调速系统的传递函数框图如图(2-4)所示：

图2-4 汽轮机液压调速系统传递函数框图

调速系统中各环节参数的范围见资料[17]。

本文中对于原动机的模型不进行重点研究，直接采用现成的原动机模型。其输入量为

∆P e 和∆ω，输出量为∆P T 。

2.3.5 电力网络的数学模型

在分析电力系统暂态过程中，忽略电力网络中的电磁暂态过程，则计算动态过程的网络模型采用节点导纳矩阵的形式。

其中： Y ：节点导纳矩阵 U ：节点电压列向量 I ：节点注入电流列向量

负荷采用恒阻抗模型时，可将负荷阻抗并入节点导纳矩阵中，消去负荷节点，只保留发电机节点。

当系统发生故障或切出故障时，即网络的节点导纳矩阵发生变化；发生故障的瞬间，节点注入电流保持不变。

YU =I

(2-6)

2.3.6 机网衔接

将式2-3表示的电压电流方程写成矩阵的形式为[3]：

⎡E q ' ⎤⎡U Gq ⎤⎡r G X d ' ⎤⎡I Gq ⎤

⎥+⎢-X ⎢⎥ ⎢⎥=⎢⎥

⎣0⎦⎣U Gd ⎦⎣q r G ⎦⎣I Gd ⎦

对式2-7进行坐标变换，由dq 轴到xy 轴，具体方法见下文，得结果：

⎡I Gx ⎤⎡G x B x ⎤⎡U Gx ⎤⎡C x ⎤

⎢I ⎥=-⎢B G ⎥⎢U ⎥+⎢C ⎥E q '

y ⎦⎣Gy ⎦⎣Gy ⎦⎣y ⎣y ⎦其中：

(2-7)

(2-8)

r G +(X q -X d ')sin δcos δ⎧G =⎪x 2

r +X q *X d ' G ⎪

⎪r G -(X q -X d ')sin δcos δ⎪G y =2

r G +X q *X d ' ⎪

⎪2⎪B =X d ' +(X q -X d ')sin δ

2⎪x

r G +X q *X d ' ⎪

⎨2

X d ' +(X q -X d ') cos δ⎪B =-2⎪y

r G +X q *X d '

⎪⎪r G cos δ+X q sin δ⎪C x =2

r +X q *X d ' G ⎪

⎪r sin δ-X q cos δ⎪C y =G

⎪r G +X q *X d ' ⎩

这些系数都是绝对角δ的函数，它们必须根据各个时间段的功角值，不断地加以修改。而网络方程（即为正常潮流计算的网络方程），以发电机端点、负荷端点和联络节点为系统节点，用网络元件的参数形成节点导纳矩阵。节点i 的方程为

或

I =∑Y ij U *j

j =1

I ix +jI iy =∑(G ij +jB ij )(U jx +jU jy )

j =1

写成矩阵的形式为：

⎡I ix ⎤n ⎡G ij -B ij ⎤

⎥⎢I ⎥=∑⎢

B G ⎢ij ⎥ij ⎦⎣iy ⎦j =1⎣式2-8与式2-9联立可得

⎡I ix ' ⎤n ⎡G ij -B ij ⎤

⎥⎢I ' ⎥=∑⎢

B G j =1iy ⎢ij ⎥⎣⎦j ≠i ⎣ij ⎦其中

⎡U jx ⎤

⎢⎥(i =1, 2,......, n ) ⎢U jy ⎦⎥⎣

(2-9)

⎡U jx ⎤⎡G ii +G ix -B ii +B ix ⎤⎡U ix ⎤⎢⎥+⎢⎥⎢U ⎥

B +B G +G U iy ii iy ⎦⎣iy ⎦⎢⎣jy ⎥⎦⎣ii

(2-10)

⎡I ix ' ⎤⎡C ix ⎤

⎢I ' ⎥=⎢C ⎥E qi ' ⎣iy ⎦⎣iy ⎦

由此可见，只要用(I ix ' I iy ' )代替(I ix I iy )，并且用式2-8中的二阶导纳矩阵修改网络节点导纳矩阵中对应的二阶对角子阵，就可以把发电机模型纳入网络方程。

暂态稳定计算中，解微分方程算出该时刻各发电机的δ后，有式2-8中的公式计算出各发电机的系数，修正各发电机节点的自导纳二阶子矩阵元素，由式2-10得出各发电机的计算用注入电流I ix ' 和I iy ' ，求解网络方程，得到发电机的端电压U x 和U y ，再由式2-9或式2-6求出发电机的定子电流I x 和I y ，即可得到发电机的电磁功率：

P e =U x I x +U y I y +(I x +I y ) r G

通过以上过程即实现了发电机与网络之间衔接的过程。

(2-11)

3 解题过程中采用的数学方法

3.1 牛顿—拉夫逊法潮流计算

牛顿—拉夫逊法是常用的解非线性方程组的方法，也是当前广泛采用的计算潮流的方法[7] [9]。

牛顿型潮流计算的核心问题是修正方程式的建立和求解。在电力系统中存在多种节点，根据给定变量的不同可以分为三类：第一类是PQ 节点，这些节点的注入功率P 和Q 是给定值；第二类是PV 节点，这些节点的电压和有功功率是已知的；第三类是平衡节点，对这种节点，它的电压保持恒定，通常是担负调整系统频率任务的发电厂母线。在修正方程式中对这些节点分别进行了研究。

为了更好的说明修正方程式，先做如下约定：

1、网络中共有n 个节点，编号为1，2，3，…，n ，其中包含一个平衡节点，编号为s 。 2、网络中有(m -1) 个PQ 节点，编号为1，2，3，…，m ，其中编号为s 的是平衡节点。 3、网络中有(n -m ) 个PV 节点，编号为m+1，m+2，…，n 。所以可以建立修正方程式如下所示：

⎡∆P ⎤⎡H N 11H 12N 12⎢1⎥⎢11

L 12⎢∆Q 1⎥⎢J 11L 11J 12

⎢⎥⎢∆P ⎢2⎥⎢H 21N 21H 22N 22⎢∆Q 2⎥⎢J 21L 21J 22L 22⎢⎥⎢

⎢⎥=⎢

⎢∆P ⎥⎢H N p 1H p 2N p 2⎢p ⎥⎢p 1

2⎥⎢∆U p ⎢R p 1S p 1R p 2S p 2

⎢⎥⎢∆P ⎢n ⎥⎢H n 1N n 1H n 2N n 2⎢2⎥⎢R S n 1R n 2S n 2∆U ⎣n ⎦⎣n 1

H 1p J 1p H 2p J 2p

H pp R pp H np R np

N 1p L 1p N 2p L 2p N pp S pp N np S np

H 1n J 1n H 2n J 2n H pn R pn H nn R nn

N 1n ⎤

⎥L 1n ⎥N 2n ⎥

⎥L 2n ⎥

⎥ ⎥N pn ⎥

⎥S pn ⎥

⎥N nn ⎥S nn ⎥⎦⎡∆f 1⎤

⎢⎥∆e ⎢1⎥⎢∆f ⎥⎢2⎥⎢∆e 2⎥⎢⎥⎢⎥ ⎢∆f p ⎥⎢⎥∆e ⎢p ⎥⎢⎥∆f n ⎢⎥⎢⎣∆e n ⎥⎦

(3-1)

其中：∆P i ，∆Q i 和∆U i 2是注入功率和节点电压平方的不平衡量

j =n

⎧

⎪∆P i =P i -∑[e i (G ij e j -B ij f j ) +f i (G ij f j +B ij e j )]

j =1

⎪

j =n ⎪⎪

⎨∆Q i =Q i -∑[f i (G ij e j -B ij f j ) -e i (G ij f j +B ij e j )]

j =1⎪

⎪∆U 2=U 2-(e 2+f 2)

i i i

⎪i ⎪⎩雅克比矩阵的元素为：

当i ≠j 时

(3-2)

∂P i ⎧H =⎪ij ∂f =-B ij e i +G ij f i

⎪⎪∂P i N =⎪ij ∂e =G ij e i +B ij f i

i ⎪

⎪∂Q i J ==-G ij e i -B ij f i =-N ij ⎪ij

∂f i

⎪

⎨∂Q ⎪L ij =i =-B ij e i +G ij f i =H ij

∂e i ⎪⎪2

⎪R =∂U i =0⎪ij ∂f i ⎪

∂U i 2⎪

⎪S ij =∂e =0

i ⎩

(3-3)

当i =j 时

∂P i ⎧H =⎪ii ∂f =-B ij e i +G ij f i +b ii

⎪⎪∂P i N =⎪ii ∂e =G ij e i +B ij f i +a ii

i ⎪j =n

⎧⎪∂Q i

⎪a ii =G ii e i -B ii f i +∑(G ij e j -B ij f j ) =-G ij e i -B ij f i +a ii ⎪J ii =j =1⎪∂f i j ≠i ⎪⎪

其中 ⎨⎨j =n ∂Q ⎪L ii =i =-B ij e i +G ij f i -b ii ⎪b =G f +B e +(G f +B e )

∑ii i ii i ij j ij j

∂e i ⎪⎪ii

j =1

⎪j ≠i ⎪⎩2

∂U i ⎪R ==2f i ii

⎪∂f i ⎪

∂U i 2⎪

⎪S ii =∂e =2e i

i ⎩

牛顿—拉夫逊法潮流计算的基本步骤为： ① 形成节点导纳矩阵；

② 设各节点电压的初值e i (0)、f i (0)或U i (0)、δi (0)；

(3-4)

③ 将各节点电压初值代入式3-2，求修正方程式中的不平衡量∆P i ，∆Q i 和∆U i 2； ④将各电压初值代入式3-3与式3-4，得到雅克比矩阵中各个元素值；

⑤解修正方程式，求各节点电压的变化量，即修正量∆e i (0)，∆f i (0)或∆U i (0)，∆δi (0)； ⑥计算各节点电压的新值， e i (1)=e i (0)+∆e i (0); f i (1)=f i (0)+∆f i (0)

或

U i (1)=U i (0)+∆U i (0); δi (1)=δi (0)+∆δi (0)

⑦运用各节点电压的新值进入下一步迭代。

3.2 隐式积分法解微分方程组

为了得到故障时发电机功角随时间变化的特性δ(t ) ，必须求解一系列的方程组。由于在这个方程组中存在非线性的常微分方程式。在一般的情况下不能求得解析解，只能用数值方法近似求解。非线性微分方程式的数值解法，是计算数学的研究课题之一，目前已有各种不同精度的算法。例如，欧拉法、改进欧拉法和龙格库塔法等显式积分法[3] [10] [11] [12]。

欧拉法和改进欧拉法由于计算精度较低且数值稳定性差，通常只适用于当所采用的元件数学模型比较简单而且所要求计算的暂态过程持续时间不长时的暂态稳定计算。龙格库塔法德计算精度优于改进欧拉法，所采用的元件数学模型可以较为详细，目前仍用于一些暂态稳定计算程序中。

与显式积分法相比，隐式积分法具有良好的数值稳定性和对刚性微分方程组的适应性，从而可以采用较大的积分步长，并可以模拟时间常数较小的环节。目前认为隐式梯形法是计算暂态稳定性比较理想的方法。

由于隐式积分法具有良好的数据计算性，因此，用在电力系统暂态稳定性计算时，一般使用较详细的数学模型：对于发电机，采用考虑阻尼绕组时的数学模型，并考虑励磁系统，原动机及其调速系统的作用。

用隐式积分法进行暂态稳定计算时，要将微分方程差分化，形成差分方程。由发电机转子运动方程, 即式2-1：

⎧d δ

=(ω-1) ω0

⎪dt ⎪

⎨d ω1

⎪=(P T -P e -D (ω-1)) ⎪⎩dt T J

在时刻t 的角度δ(t ) ，相对转速ω(t ) ，电磁功率P e (t ) ，机械功率P T (t ) 已知时，由上式两侧同时积分可得：

即为：同理可得到

δ(t +∆t ) =δ(t ) +⎰

δ(t +∆t ) =δ(t ) +

t +∆t

[ω(t ' ) -1]ω0dt '

(3-5)

∆t

[ω(t ) +ω(t +∆t ) -2]ω0 2

ω(t +∆t ) =ω(t ) +

∆t 1

**{P T (t ) -P e (t ) -D [ω(t ) -1]2T J +P T (t +∆t ) -P e (t +∆t ) -D [ω(t +∆t ) -1]}

(3-6)

式3-5与式3-6联立可得

⎧T (t +∆t ) -P e (t +∆t ) -D ω(t +∆t )]+δt 0⎪δ(t +∆t ) =αδ[P

⎨

ω(t +∆t ) =αδ(t +∆t ) +ωωt 0⎪⎩

其中：

(3-7)

ω0∆t 2⎧α=⎪δT *(2)

⎪⎪⎪δt 0=δ(t ) +ω0*∆t *[ω(t ) -1]+αδ*[P T (t ) -P e (t ) -D ω(t ) +2*D ]

⎨

⎪α=2⎪ω∆t *ω0⎪⎪⎩ωt 0=2-αω*δ(t ) -ω(t )

αδ和αω是在∆t 确定后的一个常数，这种常数称为差分常数。δt 和ωt 是一个已知数，

它在不同的时段有不同的值，这是差分方程式区别于代数方程式的主要方面。

同理，由发电机电动势变化方程，即式3-2可得：

' ' ⎧E q (t +∆t ) =α' *[E qe (t +∆t ) -E q (t +∆t )]+E q E q t 0⎪

⎨'

E (t +∆t ) =E q (t +∆t ) +αE q *I d (t +∆t ) ⎪⎩q 其中：

(3-8)

∆t ⎧

α' =⎪E q 2*T '

d 0

⎪⎪' '

⎨E q t 0=αE q ' *[E qe (t ) -E q (t )]+E q (t ) ⎪' ⎪αE q =x d -x d ⎪⎩

发电机励磁系统可以分为正常励磁和强制励磁两种形式。两种形式均为常微分方程式，分别化为差分方程的形式如下所示：

发电机励磁系统为正常励磁，即为式2-4时，化简所得的差分方程为：

E qe (t +∆t ) =αE qe *Ke *U G (t +∆t ) +E qe t 0 其中：

(3-9)

∆t ⎧α=-⎪E qe

2*Te +∆t

⎪

Te ⎨

E (t ) ⎪E qe t 0=-2*αE qe *(E qe 0+Ke *U G 0) +αE qe *(E qe (t ) +Ke *U G (t )) +qe

Te +⎪

⎩2发电机励磁系统为强制励磁，即为式2-5时，化简所得的差分方程为： E qe (t +∆t ) =αE qe E qe (t ) +E qe t 0 其中：

(3-10)

2*Te -∆t ⎧

α=⎪E qe 2*Te +∆t ⎪

∆t ⎨

E ⎪E qe t 0=∆t qe max

Te +⎪

⎩2

将以上差分方程联立（发电机励磁系统采用正常励磁方式为例），将方程3-7,3-8,3-9重写如下：

⎧δ(t +∆t ) =α[P (t +∆t ) -P (t +∆t ) -D ω(t +∆t )]+δ

δT e t 0

⎪

⎪ω(t +∆t ) =αωδ(t +∆t ) +ωt

⎪⎪' '

⎨E q (t +∆t ) =αE q ' *[E qe (t +∆t ) -E q (t +∆t )]+E q t 0

⎪'

⎪E q (t +∆t ) =E q (t +∆t ) +αE q *I d (t +∆t ) ⎪

E (t +∆t ) =αE qe *Ke *U G (t +∆t ) +E qe t 0⎪⎩qe

(3-11)

在以上方程中，由于等式两端都存在待求变量，因此这些方程是隐式方程，在本文中，采用迭代方式求解。

在t ~∆t 时段的计算中，P P U G (t +∆t ) e (t +∆t ) 和I d (t +∆t ) ，T (t +∆t ) 由原动机系统确定，

由代数方程式2-9,2-10,2-11确定，需要迭代求出的量为δ(t +∆t ) ，E q (t +∆t ) ，ω(t +∆t ) ，

E q (t +∆t ) 和E qe (t +∆t ) 。

迭代计算的思路是：先给定的初值δ(0)(t +∆t ) =δ(t ) ；ω(0)(t +∆t ) =ω(t ) ；

'(0)' (0)(0)E q (t +∆t ) =E q (t ) ；E q (t +∆t ) =E q (t ) ；E qe (t +∆t ) =E qe (t ) ，然后代入式3-11求解

⎧δ(1)(t +∆t ) =α[P (t +∆t ) -P (t +∆t ) -D ω(0)(t +∆t )]+δ

δT e t 0

⎪

⎪ω(1)(t +∆t ) =αδ(t +∆t ) +ω

ωt 0

⎪⎪'(1)(0)(0)'

⎨E q (t +∆t ) =αE q ' *[E qe (t +∆t ) -E q (t +∆t )]+E q t 0

⎪(1)

' (0)

⎪E q (t +∆t ) =E q (t +∆t ) +αE q *I d (t +∆t ) ⎪(1)

(0)

E qe (t +∆t ) =αE qe *Ke *U G (t +∆t ) +E qe t 0⎪⎩

检验是否满足下列条件：

⎧|δ(1)(t +∆t ) -δ(0)(t +∆t ) |≤ε

⎪

⎪|ω(1)(t +∆t ) -ω(0)(t +∆t ) |≤ε⎪'(1)'(0)

⎨|E q (t +∆t ) -E q (t +∆t ) |≤ε ⎪(1)(0)|E (t +∆t ) -E (t +∆t ) |≤εq q ⎪

⎪|E (1)(t +∆t ) -E (0)(t +∆t ) |≤ε

qe ⎩qe

其中ε是事先给定的精度指标，ε>0。如果满足条件，则有：

⎧δ(t +∆t ) =δ(1)(t +∆t ) ⎪

⎪ω(t +∆t ) =ω(1)(t +∆t ) ⎪' '(1)

⎨E q (t +∆t ) =E q (t +∆t )

⎪(1)

⎪E q (t +∆t ) =E q (t +∆t ) ⎪E (t +∆t ) =E (1)(t +∆t )

qe ⎩qe

计算结束。若不满足条件，则将新的到的t +∆t 值代入式3-11继续求解，直到满足下列收敛条件：

⎧|δ(k +1) (t +∆t ) -δ(k ) (t +∆t ) |≤ε⎪

⎪|ω(k +1) (t +∆t ) -ω(k ) (t +∆t ) |≤ε⎪'(k +1) '(k )

⎨|E q (t +∆t ) -E q (t +∆t ) |≤ε ⎪(k +1) (k ) |E (t +∆t ) -E q (t +∆t ) |≤ε⎪q

⎪|E (k +1) (t +∆t ) -E (k ) (t +∆t ) |≤ε

qe ⎩qe

迭代结束。可以进行下一时段的求解。

3.3 xy -dq 坐标变换

在电力系统动态分析时，各个发电机和网络的电压电流向量可能采用不同的坐标系，常见的是发电机采用各自的dq 坐标，网络采用xy 坐标，需要用各个发电机的绝对角将它们的坐标统一在xy 坐标系下[3]。

在暂态稳定计算中，将同步旋转的参考轴选为与x 轴重合，则转子的绝对角为发电机转子的q 轴与x 轴之间的夹角。

两种坐标系的变换为下列关系：由xy 轴到dq 轴：

图3-1 xy -dq 坐标变换图

⎡f d ⎤⎡sin δ⎢f ⎥=⎢⎣q ⎦⎣cos δ-cos δ⎤⎡f X ⎤

⎢⎥⎥sin δ⎦⎣f Y ⎦-cos δ⎤⎡f d ⎤⎢f ⎥

sin δ⎥⎦⎣q ⎦

由dq 轴到xy 轴：

⎡f X ⎤⎡sin δ

⎢f ⎥=⎢cos δ⎣Y ⎦⎣

3.4 化简网络中的无源节点

随着电力事业的发展，电网逐步形成巨大的互联系统，网络中的节点数越来越多，网络的拓扑结构也越来越复杂。而机电暂态计算时，必然要面对网络方程的求解。系统的规模越大，网络方程的阶数越高，求解中消耗的内存和时间也越多。这时候我们就需要寻找一种方法，在尽可能保证所研究部分的动态静态特征的前提下，化简网络中不感兴趣的节点。例如，在电力系统中往往存在许多既不接发电机也不接负荷的节点，这些节点的存在不会影响系统暂态稳定性的判断，因而可以舍去，对这类节点进行化简的时候，只需要对节点导纳矩阵做些修改就可以了[14]。

利用高斯消去法化简网络，可以一次消去若干个节点。其方法如下：

令子集I 表示要保留的节点，子集E 表示要消去的节点（注：消去点均为无源点，则有网络的节点导纳矩阵为：

⎡Y II Y =⎢

⎣Y EI

Y IE ⎤

⎥Y EE ⎦

各节点注入电流列向量为：

各节点电压列向量为：即有推出

⎡U I ⎤

U =⎢⎥

⎣U E ⎦

⎡I ⎤I =⎢I ⎥ ⎣0⎦

⎡Y II ⎢Y ⎣EI Y IE ⎤Y EE ⎥⎦⎡U I ⎤⎡I I ⎤

⎢U ⎥=⎢0⎥ ⎣E ⎦⎣⎦

即

⎧Y II *U I +Y IE *U E =I I ⎨

⎩Y EI *U I +Y EE *U E =0

-1IE EE EI

(Y II -Y Y Y )*U I =I I

所以，化简网络中的无源节点以后所得的只含有源节点的网络的节点导纳矩阵为

' -1Y =Y -Y Y Y EI II IE EE

通过高斯消去法得到的简化网络，方程的阶数从节点数下降为发电机数，极大地减少了运算量和运算误差，使仿真更加精确。

4 机电暂态仿真计算流程

图4-1 机电暂态稳定分析流程框图

图4-2 隐式积分法部分框图

图4-3 牛顿-拉夫逊法潮流计算部分框图

如图4-1为机电暂态仿真流程框图[16] [17]。

第一步，输入原始数据，包括系统元件参数，网络拓扑结构，仿真步长，仿真时间，故障点信息等。

第二步，对输入的各元件参数进行归算整理。

第三步，根据输入的信息进行系统全网的潮流计算，得到各节点的电压，电流，并计算出负荷的等效阻抗，为后面的暂态稳定计算做准备。

第四步，根据各元件参数以及全网的拓扑结构形成节点导纳矩阵。注意发电机暂态电抗和负荷等值电抗的处理。并且要在这一步形成暂态稳定计算的初始值，置计数单元t=0。

第五步，判断是否有故障发生。若有，则需要根据故障情况修改节点导纳矩阵和微分方程。若无扰动，直接进入下一步的计算。

第六步，采用隐式积分法求解微分方程，根据T n 时刻系统的状态量求出T n +1系统的状态量。值得注意的是，网络拓扑结构的变化，将引起某些运行变量的突然变化，例如节点电压，发电机输出的电磁功率等。这不符合客观实际。因此，计算时应对这种情况进行特殊处理。例如，在第i 时段网络结构发生变化，造成发电机输出电磁功率从P e ' (i ) 变化至P e '' (i ) ，则在第i 时段末的电磁功率一般取其平均值

' ''

P e (i ) =(P e (i ) +P e (i ) ) /2

此时即完成一个步长的计算。

第七步，判断仿真总时间是否到达，若已到达，则输出结果，结束仿真；若未到达, 计数单元t=t+1，进行下一轮的计算。

其中第五、六步的运算是整个仿真程序的核心，其框图单独列出，如图4-2所示：其公式及其原理可见前文所述。

牛顿—拉夫逊法潮流计算的程序框图如图4-3所示。

5 算例验证及结果分析

选择一个简单多机系统的电路图如图5-1所示

图5-1 简单模型系统电路图

该系统的其它参数为：

节点④处的发电机1：x d =0.369；x q =0.0318；x q =0.02；x 2=0.025；r =0.005；

发电机惯性时间常数T J =1.57；

励磁绕组时间常数T f =8.18；发电机阻尼常数D =10；自动调节励磁系统时间常数T e =0.5；自动调节励磁系统等效放大倍数K e =40；强制空载电动势最大值E qe max =3

' 节点⑤处的发电机2：x d =0.369；x q =0.0318；x q =0.02；x 2=0.025；r =0.01；

发电机惯性时间常数T J =1.97；

节点②④之间变压器：电抗x t =0.015 节点③⑤之间变压器：电抗x t =0.03 线路②③：负序电抗x 2=0.9 线路①②：负序电抗x 2=0.75 线路①③：负序电抗x 2=1.05 接地阻抗Z f =0.001+j0.001

（1）故障发生在线路①②的靠近节点②处中点处，故障类型为三相短路，3秒发生故障，3.25秒切除故障，最终恢复到稳定状态的情况。

图①

图②

图③

图④

图⑤

图5-2 系统切除故障后恢复稳定算例结果（一）

①发电机1机端电压变化过程

②发电机1电磁功率（红色）与有功功率（蓝色）变化过程 ③发电机2机端电压变化过程

④发电机2电磁功率（红色）与有功功率（蓝色）变化过程 ⑤各发电机绝对角（蓝色）与相对角（红色）变化过程

秒切除故障，最终恢复到稳定状态的情况。

图③ 图④ 图① 图②

图⑤ 图5-3 系统切除故障后恢复稳定算例结果（二） ①发电机1机端电压变化过程 ②发电机1电磁功率（红色）与有功功率（蓝色）变化过程 ③发电机2机端电压变化过程 ④发电机2电磁功率（红色）与有功功率（蓝色）变化过程 ⑤各发电机绝对角（蓝色）与相对角（红色）变化过程

障，5秒切除故障，最终不能恢复到稳定状态的情况。

图③ 图④ 图① 图②

图⑤ 图5-4 系统切除故障后不能恢复稳定算例 ①发电机1机端电压变化过程 ②发电机1电磁功率（红色）与有功功率（蓝色）变化过程 ③发电机2机端电压变化过程 ④发电机2电磁功率（红色）与有功功率（蓝色）变化过程 ⑤各发电机绝对角（蓝色）与相对角（红色）变化过程

6 结论

世界各国的电力系统发展得越来越庞大，电力系统的运行也随之越来越复杂电力系统数字仿真技术的重要性也越来越凸显出来。本文研究了基本的电力系统暂态分析方法，在分析电力系统暂态稳定计算基本原理的基础上，考察了电力系统仿真中现存的各种方法，并最终选定了隐式积分法来计算暂态分析过程中产生的微分方程组。

本文的最终成果是在C 语言环境下形成了能够解决小系统暂态稳定分析与计算的程序。

该程序建立的主要过程为：①建立电力系统各元件的数学模型；②根据元件之间的拓扑关系将各元件的数学模型组合为全系统模型；③然后采用适当的数学方法来进行求解。在暂态稳定计算之前的潮流计算中，采用了牛顿—拉夫逊法；在暂态稳定计算中，采用了隐式积分法来解非线性方程，并且采用了高斯消去法来化简网络中的无源节点，大大降低了迭代次数，提高了程序效率，减少了计算时间。

最后采用matlab 绘图得到各发电机端电压，电磁功率，机械功率，绝对功角，相对功角曲线，用来判断系统在扰动后是否能够恢复稳定状态。

从matlab 所绘图中分析可知：当系统发生短路时，发电机端电压与电磁功率下降，发电机间的相对角剧烈波动，系统不稳定；当系统切除故障时，发电机端电压恢复原来的状态，电磁功率到达一个新的平衡状态；当发生三相短路时，系统的稳定性最差，不容易恢复到稳定状态，这个结果与预计相符，说明本程序的设计基本符合设计要求，能够解决预计问题。

在本次毕业设计中，我学到了很多知识，对电力系统机电暂态的了解更加深刻，我认为本次毕业设计是成功的，达到了基本目的。

7 参考文献

[1] 石访. 电力系统机电暂态和电磁暂态数字混合仿真[J].电网与清洁能源.2009，25(12)

[2] 潘学萍．电力系统数字仿真研究综述[J]．江苏电机工程.2005，24（1）

[3] 常鲜戎，赵书强．电力系统暂态过程[M]．机械工业出版社，2010

[4] 李光琦．电力系统暂态分析[M]．北京：中国电力出版社，2007

[5] 西安交通大学，清华大学，浙江大学，湖南大学，成都工学院，水利电力部电网调度研究所．电力系统计算[M]．水利电力出版社，1985

[6] 刁士亮，柳中莲．同步发电机组模型与电力系统稳定分析[M]．北京：华南理工大学出版社，1989

[7] 吴际舜．电力系统静态安全分析．上海：上海交通大学出版社．1984

[8] 于庆广，胥埴伦．电力系统的电磁暂态与机电暂态混合仿真[J]．中国电力，2007，40(11)

[9] 陈珩．电力系统稳态分析[M]．北京：中国电力出版社，2007

[10] Nikolay F.Djagarov ． A method of transient electromechanical processes modeling in power systems[J]．IEEE Xplore．2009．

[11] 马东升, 雷永军．数值计算方法[M]．机械工业出版社，2008

[12] 同济大学应用数学系．线性代数[M]．北京：高等教育出版社，2003

[13] Hui Yang；Fushuan Wen；Liping Wang．Newton-Raphson on power flow algorithm and Broyden Method in the distribution system．IEEE Xplore．2008.

[14] 常鲜戎．电力系统分析[M]．保定：华北电力大学内部交流讲义，2007

[15] X iaoyin Yao ，Karady, G .G . ，Farmer, R.G . ，Agrawal, B.L. ．Estimation of generator excitation-system parameters by tabu search．IEEE Xplore．2001．

[16] 谭浩强．C 语言程序设计[M]．2005

[17] 成颖．C/C++程序设计语言．南京：东南大学出版社．2003

致谢

不知不觉间，我的大学生活马上就要画上句号了，在大学四年的生活中，我学到了很多，也提高了很多。在这里，向所有曾向我提供帮助的老师和同学致以最诚挚的谢意。

毕业设计是对我们知识运用能力和学习能力的一次全面的考核，培养我们综合运用所学知识来分析和解决问题的能力。本人的本科毕业设计一直是在导师常鲜戎老师的悉心指导下进行的。常老师的治学作风严谨细致，一丝不苟，在整个毕业设计的过程中，常老师一直引导我思考问题，并提出新的问题，使我的毕业设计能够比较顺利的深入下去，在此表示诚挚的感谢和由衷的敬意。

在建模期间，针对所要讨论的内容孙一莹同学与我进行了多次深入讨论；并且在程序的调试过程中，孙一莹同学也给了我很大的帮助。在此，对她表示深深地感谢。

感谢宿舍的姐妹们，没有你们，我就不会有这么快乐的大学生活，希望你们在以后的日子里越过越幸福。

感谢所有在大学期间曾经给予我关心，支持我，帮助我的人们，我会永远记住你们，谢谢！

谨以此文来结束我大学四年的快乐生活，我的师长与伙伴们，祝你们永远幸福！！