实验三 电力系统暂态稳定分析
电力系统暂态稳定计算实际上就是求解发电机转子运动方程的初值问题，从而得出δ-t和ω-t的关系曲线。每台发电机的转子运动方程是两个一阶非线性的常微分方程。因此，首先介绍常微分方程的初值问题的数值解法。
常微分方程的初值问题
(一)问题及求解公式的构造方法
我们讨论形如式(3-1)的一阶微分方程的初值问题
(3-1)
设初值问题(3-1)的解为，为了求其数值解而采取离散化方法，在求解区间[]上取一组节点
称()为步长。在等步长的情况下，步长为
用表示在节点处解的准确值的近似值。
设法构造序列所满足的一个方程(称为差分方程)
(3-2)
作为求解公式，这是一个递推公式，从(，)出发，采用步进方式，自左相右逐步算出在所有节点上的近似值()。
在公式(3-2)中，为求只用到前面一步的值，这种方法称为单步法。在公式(3-2)中的由明显表示出，称为显式公式。而形如(3-3)
(3-3)
的公式称为隐式公式，因为其右端中还包括。
如果由公式求时，不止用到前一个节点的值，则称为多步法。
由式(3-1)可得
= (3-4)
两边在[，]上积分，得
(3-5)
由此可以看出，如果想构造求解公式，就要对右端的积分项作某种数值处理。这种求解公式的构造方法叫做数值积分法。
(二)一般的初值问题的解法
欧拉法和改进欧拉法
对于初值问题(3-1)，采用数值积分法，从而得到(3-5)。对于(3-5)右端的积分用矩形公式(取左端点)，则得到
进而得到(3-1)的求解公式(3-2)
(=0，1，2，n-1) (3-6)
此公式称为欧拉(Euler)格式。
如果对式(3-5)右端的积分用梯形公式
则可以得到初值问题(3-1)的梯形求解公式如式(3-7)
(=0，1，2，n-1) (3-7)
式(3-7)是个隐式公式。可以采取先用欧拉格式求一个的初步近似值，记作，称之为预报值，然后用预报值替代式(3-7)右端的，再计算得到，称之为校正值，这样建立起来的预报－校正方法称为改进欧拉格式
(3-8)
龙格—库塔方法
在单步法中，应用最广泛的是龙格－库塔(Runge-kutta)法，简称R－K法。下面直接给出一种四阶的龙格－库塔法的计算公式(3-9)
(3-9)
它也称为标准(古典)龙格－库塔法。
例3-1 研究下列微分方程的初值问题
解：
这是一个特殊的微分方程，其解的解析式可以给出，为
应用龙格－库塔法，取=0.25，根据式(3-9)编写一段程序，由零开始自左相右逐步算出在所有节点上的近似值。计算结果见表3-1。计算结果表明，四阶龙格－库塔方法的精度是较高的。
表3-1
2.004.3e-54.002.5e-66.003.7e-78.009.2e-8实际上，MATLAB为常微分方程提供了很好的解题指令，使得求解常微分方程变得很容易，并且能将问题及解答表现在图形上。因此，我们可以不用根据式(3-9)编写较复杂的程序，而只需应用MATLAB提供的常微分方程解题器来解决问题。下面给出用MATLAB编写的解题程序。
首先编写描述常微分方程的ODE文件，文件名为ˊmyfunˊ，便于解题器调用它。
function dy = myfun(x,y)
dy = zeros(1,1);
dy=1/(1+x^2)-2*y^2;
再编写利用解题器指令求解y的程序。
clear
x0=0;
for i=1:4
xm=2*i;
y0=0;
[x,y] = ode45('myfun',[x0 xm],[y0]);
format long
y(length(y))
end
plot(x,y,'-')
运行上述程序，在得到几个点的函数值的同时，也得到函数y的曲线，如图3-1所示。
图3-1 根据运算结果画出y的曲线
简单电力系统的暂态稳定性
(一)物理过程分析
某简单电力系统如图3-2(a)所示，正常运行时发电机经过变压器和双回线路向无限大系统供电。发电机用电势作为其等值电势，则电势与无限大系统间的电抗为
(3-10)
这时发电机发出的电磁功率可表示为
(3-11)
如果突然在一回输电线路始端发生不对称短路，如图3-2(b)所示。故障期间发电机电势与无限大系统之间的联系电抗为
(3-12)
在故障情况下发电机输出的电磁功率为
(3-13)
在短路
