基于多机电力系统单一发电机电磁功率简洁表达方式的暂态稳定在线分析方法 【技术领域】
本发明涉及电力系统暂态稳定分析技术，特别是多机系统单一发电机电磁功率简洁表达方式及其在电力系统暂态稳定在线分析方面的应用。
背景技术
电力系统稳定分析在寻找电力系统稳定裕度和失稳机理方面取得了大量的研究成果。
在分析方法方面主要有暂态能量函数法，数值仿真法和二者相结合的混合法。能量函数法又包括基于势能边界法与相关不稳定平衡点法的能量函数法，扩展等面积法(EEAC及SIME)，单一发电机等面积(SGEAC)等方法。时域仿真的优点是能考虑任意复杂的系统模型，通过多次计算，可得到基于临界切除时间的稳定裕度，但不能快速得到反应系统稳定裕度的量化指标。对大系统而言，由于很难确定失稳模式及相关不稳定平衡点，基于相关不稳定平衡点法的能量函数法使用起来不是很方便。基于最大势能搜索的势能边界法(PEBS)很难得到稳定裕度量化指标。为了克服势能边界法的能量函数法的不准确性，在系统惯性中心的基础上采用了基于双机等值的修正动能，该修正动能修正掉了与二群各自内部各机组之间的“布朗运动”相对应的，对系统失稳没有影响的动能。EEAC及其改进方法在同步坐标基础上定义了动态二群机组的等值角度与速度，并进一步定义了基于动态二群相对运动且能反应系统失稳裕度的加速面积与减速面积(动能与势能)。
在理论研究的基础上，扩展等面积法EEAC(Extended equal area criterion)或者单一发电机无限大等值机方法SIME(Single Machine Equivalent)动态地将暂态过程中的同步发电机分为前向机群与后向机群。在二群等值的基础上，进一步将二等值机系统简化为一单一发电机无限大系统，等值机的电磁功率表达为前向等值机群惯性中心角与后向机群惯性中心角之差。但本方法需要动态划分前向机群与后向机群，而前向机群与后向机群的划分是非常难于实现的。
多机系统中单一发电机电磁功率曲线非常复杂，理论推导可知，当负荷采用恒定阻抗模拟时，任意单一发电机的输出电磁功率是系统中所有发电机转子角度的多元非线性函数。当负荷采用感应电动机模型时情况更为复杂，任意单一发电机的输出电磁功率是系统中所有发电机转子角度与感应电动机滑差的多元非线性函数。
多机系统的单一发电机等面积定则暂态稳定判据的特点是没有对电力系统动态模型作任何假设，不需要对失稳模式进行任何假设，稳定判别结果与时域仿真完全一致。与时域仿真法相比较，单一发电机等面积法可方便地求得基于能量的稳定与不稳定裕度。
【发明内容】
本发明的目的是，提出一种基于多机电力系统单一发电机电磁功率简洁表达方式的暂态稳定在线分析方法，即在多机系统单一发电机输出电磁功率特性分析的基础上，提出多机电力系统单一发电机电磁功率简洁表达方式，采用时域仿真与数据拟合得到用惯性中心坐标系下发电机本身转子角度表达的功角曲线，将其外推得到发电机的不稳定平衡点角度，采用事后评估与预测评估得到系统的稳定裕度。本发明能够实现预想事故电力系统暂态稳定快速在线分析，用于大规模电力系统在线暂态稳定评估，为电网运行人员提供当前电网运行状态在预想故障情况下系统稳定裕度有多大。
为达到上述发明目的，本发明采取的技术方案是，一种基于多机电力系统单一发电机电磁功率简洁表达方式的暂态稳定在线分析方法，包括如下基础步骤：
(1)由电力系统的数据采集系统采集得到全网数据；
(2)采用状态估计方法对采集到的全网数据进行处理，得到全网当前运行状态；
(3)采用时域仿真对预想事故下当前电网进行仿真，得到故障后电网暂态功角摇摆曲线与电磁功率曲线；计算得到故障切除以后系统惯性中心坐标系COI下单一发电机的转子角度与输出电磁功率，在正向摇摆与负向摇摆群中选定多台临界机组作为临界机群；根据当前时刻前一段时间内各临界机的功角与电磁功率数据，采用数据拟合得到各发电机功角曲线的二次或高次多项式系数；将当前时刻前拟合得到的功角曲线外延得到当前时刻以后的未来时刻的发电机电磁功角曲线，利用单一发电机等面积准则计算得到预测临界机群中的各单一临界机组的稳定裕度；
(4)不稳定故障的系统不稳定裕度为临界机群中的临界机组在不稳定平衡点的稳定裕度最小者；稳定故障的系统稳定裕度为临界机群中各临界机组在最大转子角度处的稳定裕度最小者；
(5)稳定裕度为负时，警示预想事故，转其它预防控制或紧急控制分析模块；
(6)重复进行上述步骤(3)至(5)，实现预想事故下电力系统暂态稳定的在线分析。
其中，步骤3中当前时刻前一段时间一般选取当前时刻的前0.1-0.2秒内。
本发明的工作原理和过程进一步描述如下：
电力系统的数据采集系统SCADA(Supervisory Control And DataAcquisition)采集得到全网相关数据，包括节点电压、网络拓扑、节点注入功率等。再由状态估计进行数据处理后得到全网当前运行状态；状态估计方法可采用当前成熟的、在电力调度中心使用的任意商用软件算法。
采用时域仿真对预想事故下当前电网进行仿真研究，计算系统惯性中心坐标系COI下单一发电机组的转子角度与输出电磁功率，利用数据拟合得到当前时刻前一段时间内单一发电机输出电磁功率简洁表达形式中的多项式系数，将其外延得到预测的单一发电机电磁功率，预测临界机群中地单一发电机稳定裕度；具体过程是：
如果负荷采用恒定阻抗并入电网，由理论分析可知机组i的等效电磁功率可表示为：
Pei′=f(θ~1,...,θ~g)]]>
上面式中的(其中i为1，2，…g”)为发电机i相对系统惯性中心的转子角度，g为发电机台数。由上式可有如下发电机等效电磁功率对系统惯性中心坐标系COI下发电机本身转子位置角度的一阶导数：
dPei′dθ~i=∂f∂θ~1×dθ~1dt×dtdθ~i+...+∂f∂θ~gdθ~gdtdtdθ~i=∂f∂θ~1×dθ~1dt×1dθ~i/dt+...+∂f∂θ~gdθ~gdt1dθ~i/dt]]>进一步对上式微分，可计算出电磁功率对功角的高阶导数。
如果负荷采用考虑机电暂态的异步机负荷，由理论分析可知发电机i的等效电磁功率可表示为：
Pei′=f(θ~1,...,θ~g,s1,...,sl)]]>
上面式中为发电机i相对系统惯性中心的转子角度，sj(其中j为1，2，…l”)为异步感应电动机j的转子滑差，g为发电机台数，l为感应发电机台数。由上式可有如下发电机等效电磁功率对系统惯性中心COI坐标系下发电机本身转子位置角度的一阶导数：
dPei′dθ~i=∂f∂θ~1×dθ~1dt×dtdθ~i+...+∂f∂θ~gdθ~gdtdtdθ~i+∂f∂s1×ds1dt×dtdθ~i+...+∂f∂sldsldt1dθ~i]]>
=∂f∂θ~1×dθ~1dt×dtdθ~i/dt+...+∂f∂θ~gdθ~gdtdtdθ~i/dt=∂f∂s1×ds1dt×1dθ~i/dt+...+∂f∂sldsldt1dθ~i/dt]]>进一步对上式微分，可计算出电磁功率对功角的高阶导数。
任意单一发电机在达到最大角度之前P′e，i与之间是一条连续可微曲线。由泰勒级数可有
Pe,i′=Pe,i′|θ~i=θi0~+dPe,i′dθ~i|θ~i=θi0~(θ~i-θi0~)+...+1k!×dkPe,i′dθik~|θ~i=θi0~(θ~i-θi0~)k]]>
上式可以简写成
Pe,i′=fi(θ~i)=b0,i+Σk=1lbk,i(θ~i)k]]>
上式即为单一发电机电磁功率简洁表达方式，其中b0，i为发电机i的功角曲线常数项；bk，i为发电机功角曲线中第k次多项式的系数。
由于上面高阶导数计算起来十分复杂，利用仿真得到的(P′e，i，)数据，通过拟合得到的一元高次多项式来预测单一发电机从当前功角点到不稳定平衡点之间的P-曲线。由数学分析可知，在展开点附近的一个邻域内，泰勒级数能较好地吻合原曲线。由当前功角点之前的时刻点上功角曲线数据拟合得到的功角曲线函数应该可以很好地预测当前功角点到不稳定平衡点之间的曲线。
不稳定故障的系统不稳定裕度为临界机群中的临界机组在不稳定平衡点的不稳定裕度最小者；
即不稳定情况下基于单一发电机等面积定则的机组不稳定裕度：
uηi=-ω~icr,2ω~iS,2]]>
系统不稳定裕度为：uηs＝minuηi
稳定故障的系统稳定裕度为临界机群中各临界机组在最大转子角度处的稳定裕度最小者；
即稳定情况下基于单一发电机等面积定则的加速机组的系统稳定裕度：
ηi=SDi12Miω~is,2]]>
SDi=∫θ~isθ~icr[Pmi-Pei′]dθ~i]]>
其中是临界发电机i的最大功角点功角；是临界发电机i的不稳定平衡点功角；是故障切除时刻临界机组i的角速度或稳定平衡点时刻临界发电机i的角速度；是不稳定平衡点时刻临界发电机i的角速度。
系统稳定裕度为：ηs＝minηi
通过解下式，可求得发电机i的不稳定平衡点，
fi(θ~i)=Pm,i]]>
式中为故障切除后拟合确定的曲线。
同时，警示稳定裕度为负的预想事故，转其它预防控制或紧急控制分析模块；重复进行上述步骤，实现预想事故下电力系统暂态稳定的在线分析。
采用本发明所述分析方法分析大系统稳定性时，采用单一发电机等面积方法不必计算所有机组的稳定裕度，只需计算前向机群或后向机群中的少数几台临界机组的稳定裕度就可得到系统稳定性。
综上所述，本发明提供的所述基于多机电力系统单一发电机电磁功率简洁表达方式的暂态稳定在线分析方法的优点包括如下几个方面：
(1)本发明所述分析方法得到了多机系统中单一发电机电磁功率的简洁表达方式。将多机系统中单一发电机多元非线性功角函数表达为一用发电机本身功角表达的一元高次多项式函数，使得单一发电机等面积方法能够方便地应用于多机电力系统的暂态稳定分析；
(2)利用时域仿真得到故障切除以后系统的暂态轨迹，得到临界机群中单一发电机在系统惯性中心坐标系COI下的转子角度与输出电磁功率，方法不受模型的限制；
(3)采用事后评估方法评估系统的稳定裕度时，在确定稳定系统的稳定裕度时，利用转子角度达到最大值之前的一段时间内的数据，拟合得到临界机群中各单一发电机功角曲线多项式系数，将其外延得到时域仿真不能得到的功角曲线，计算得到临界机组的稳定裕度；对于不稳定故障，利用不稳定平衡点时刻机组转子的角速度得到机组的不稳定裕度，结果可靠；
(4)采用预测评估时，利用当前时刻前的一段时间内的数据，拟合得到临界机群中各临界机组的功角曲线，外延预测得到未来时刻机组的电磁功角曲线。利用单一发电机等面积定则可以预测到各临界机组的不稳定裕度。可以提前预测得到系统的不稳定裕度与稳定裕度。
本发明所述基于多机电力系统单一发电机电磁功率简洁表达方式的暂态稳定在线分析方法适用于大规模电力系统预想事故在线暂态稳定评估。
【附图说明】
图1为本发明所述分析方法的流程框图；
图2为实施例1中6机系统失稳摇摆曲线图；
图3为实施例1中6机稳定系统摇摆曲线图；
图4为实施例2中IEEE17机系统75号节点三相短路，0.356秒切除故障线路时发电机转子摇摆曲线图；
图5为实施例2中IEEE17机系统中支路25-26靠近26节点侧发生三相接地，0.196秒切除故障线路时摇摆曲线图；
其中：图2-图5中的横坐标为时间，纵坐标为转子摇摆角度(弧度)。
【具体实施方式】
如图1所示，由电力系统的数据采集系统SCADA采集得到全网相关数据，如，发电功率、负荷功率、网络拓扑等，由任意商用状态估计软件进行数据处理后得到全网当前运行状态。
采用时域仿真对预想事故下当前电网进行仿真研究，计算系统惯性中心坐标系COI下单一发电机组的转子角度与输出电磁功率，利用数据拟合得到当前时刻前单一发电机输出电磁功率理论表达形式中的多项式系数，将其外延得到预测的单一发电机电磁功率，预测临界机群中各单一发电机稳定裕度；不稳定故障的不稳定裕度为临界机群中各临界机组在不稳定平衡点的稳定裕度最小者；稳定故障的系统稳定裕度为临界机群中各临界机组在最大转子角度处的稳定裕度最小者；警示稳定裕度为负的预想事故，转其它预防控制或紧急控制分析模块；重复进行上述过程，实现预想事故下电力系统暂态稳定的在线分析。
实施例1：
本实施例采用事后评估与预测评估的方法来评估系统的稳定裕度。对于失稳情况，事后评估是根据仿真结果得到故障切除时刻或稳定平衡点时刻研究发电机的角速度，以及研究发电机在不稳定平衡点时刻的角速度，计算得到相应的稳定裕度。对于稳定情况，事后评估是根据仿真得到的研究发电机到达最大功角之前的功角曲线数据拟合得到最大功角点到不稳定平衡点之间的功角曲线，利用本发明提出的相关计算公式计算得到相应稳定裕度。所谓预测评估就是利用当前功角点之前的功角曲线数据，利用拟合并外延得到相应功角曲线，能够在机组达到最大功角点或者不稳定平衡点之前并计算出相应机组的稳定裕度。
本实施例介绍的是6机系统失稳故障的事后评估和预测评估，6机系统的相关原始数据参见后述附录。
其中6机系统失稳故障的事后评估为：
6机系统在5号节点发生三相瞬时短路，0.36秒故障，摇摆曲线如图2所示。在故障切除以后的最初阶段，根据单一发电机等面积定则，临界发电机G1、G2、G3、G6的稳定裕度如表1所示。表中为故障切除时刻或平衡点时刻发电机的角加速度；为不稳定平衡点时刻发电机的角加速度；tcr为发电机越过不稳定平衡点的时间。由表中单一发电机的一摆稳定裕度可知，系统的稳定裕度为-0.568。
表1：临界单一发电机不稳定裕度
采用EEAC(扩展等面积)评估系统稳定裕度时，不同的机组分群方式有不同的评估结果。由表2数据可知，不同的分群方式会有不同的稳定裕度结果。
表2：不同分群得到的稳定裕度
6机系统失稳故障的预测评估：
在当前功角点时，可以得到稳定裕度计算式为
ηs=∫θ~θ~cr(Pe′-Pm)dθ~0.5Mω~s,2-ω~2ω~s,2]]>
式中与分别为研究发电机在当前时刻的功角与角速度，为研究单一发电机在稳定平衡点时刻的角速度或故障切除时刻的角速度，P′e为研究单一发电机的功角曲线，Pm为研究单一发电机的机械功率，为研究单一发电机的不稳定平衡点功角。随着仿真的推进，不断向不稳定平衡点功角靠近，拟合的P′e得到不断更新，得到的预测稳定裕度不断向事后评估裕度逼近，因而仿真时间越长计算得到的稳定裕度越准确。上面的稳定裕度计算式中，与都是仿真(测量)得到的，只有在与之间的功角曲线及是由拟合或求解得到。裕度计算公式计算稳定裕度时充分利用了已有仿真或测量结果，对拟合P′e的外延区间减到了最小，采用简单的二次函数拟合即可得到令人满意的结果。利用不同样本时间的数据，发电机2的功角曲线拟合结果、预测的不稳定平衡点角度与稳定裕度如表3所示。由时域仿真得到的不稳定平衡点角度为-1.39，随着仿真时间的加长，拟合预测得到的不稳定平衡点角度越来越接近时域仿真得到的准确的不稳定平衡点角度。发电机3和发电机6的预测稳定裕度如表4、5所示。
表3：发电机2的预测稳定裕度
表4：发电机3的预测稳定裕度
表5：发电机6的预测稳定裕度
6机系统稳定场景的稳定评估：
对稳定故障，即可以采用预测评估的方法，也可以采用事后评估的方法进行稳定裕度的量化计算。系统故障为15号节点发生三相短路，0.16秒故障切除。摇摆曲线如图3所示，在观察时段内系统是稳定的。对正向摇摆的临界机组与负向摇摆的临界机组进行了稳定裕度的评估，其结果如表6与表7所示。由时域仿真与曲线拟合得到各临界机组的功角曲线，求得相应的单一发电机不稳定平衡点角度，得到各机组的稳定裕度。表6与表7稳定裕度结果采用的拟合样本时间为各临界机组达到最大功角点之前0.1秒。应当特别指出的是，稳定情况下，机组的不稳定平衡点通过时域仿真是得不到的，将本发明拟合得到的功角曲线外推可以得到仿真无法得到的不稳定平衡点角度。由摇摆曲线与表6及表7所示的数据可知，表中的稳定裕度大小正确反应了机组摇摆的严重程度。
表6：负向摇摆临界机的预测稳定裕度
表7：正向摇摆临界机的预测稳定裕度
实施例2：
IEEE17机系统测试的第一个故障是支路75-9靠近75号节点处三相短路，0.356秒切除故障线路。摇摆曲线如图4所示。故障切除以后，系统一摆是稳定的，二摆失去稳定。利用故障切除以后0.1秒的数据，由拟合得到单一发电机简洁功角曲线及预测的单一发电机一摆稳定裕度如表8所示。随着采样时间接近最大转子角度时间，预测得到的单一发电机稳定裕度越来越接近表9所示的事后评估得到的稳定裕度。预测的处于130节点的发电机二摆稳定裕度为负，系统失去稳定，与仿真情况一致。
表8：样本时间为故障切除后01秒下IEEE17机系统单一发电机功角曲线及稳定裕度
表9：样本时间为最大功角前0.1秒下IEEE17机系统单一发电机功角曲线及稳定裕度
IEEE17机系统的第二个测试故障是支路25-26靠近26节点侧发生三相接地，0.196秒切除故障线路。摇摆曲线如图5所示，由摇摆曲线可知在节点73与76上的临界发电机组是前向机组，系统以二群失稳的模式失去同步稳定。由单一发电机等面积法则得到76号机上的稳定裕度为-0.732。与实际情况相符。
实施例3：
某电力系统有584个节点，531个负荷，87台发电机以及1284条支路。EGBZ2，EGBZ3节点间463号线路靠近EGBZ2侧发生三相短路，0.2秒故障切除。故障切除以后，采用本发明方法得到4个临界机群的一摆稳定裕度分别如表10所示。表中UEP1是样本时间为0.2-0.3秒预测得到的不稳定平衡点，预测稳定裕度是根据拟合得到的单一发电机功角曲线，采用预测评估方法得到的。UEP2为采用时域仿真得到的、准确的不稳定平衡点，但由于196号节点上的发电机一摆稳定，其UEP2是根据功角曲线与原动机功率平衡求解得到。节点190、191及193节点上的发电机单一发电机稳定裕度是根据时域仿真得到的不稳定平衡点时刻发电机转子速度得到的。由于196节点上发电机一摆稳定，其稳定裕度是由拟合得到的功角曲线得到的。系统稳定裕度为单一发电机稳定裕度的最小者。系统稳定裕度与仿真结果完全一致。
由表10数据可知预测的稳定裕度、不稳定平衡点与仿真得到的稳定裕度、不稳定平衡点非常接近。
表10：某87机系统单一发电机功角曲线及稳定裕度
所在节点 UEP1 (预测) 稳定裕度 (预测) UEP2 (实际测量) 事后评估的 稳定裕度 190 1.1729 -1.1397 1.1862 -1.29063 191 1.1832 -1.1398 1.2028 -1.24097 192 1.4512 -1.2937 1.497 -1.14091 196 2.5942 2.9351 2.5048 2.9156
附录：6机系统参数
本附录给出了6机系统标幺值参数，基准值为100MVA。
表A1：支路阻抗
支路Branch(I→J) 电阻R 电抗X 1→18 0.021 0.064 1→23 0.0 0.0795 2→4 0.0521 0.216 2→12 0.0298 0.178 2→24 0.0 0.0397 3→5 0.009 0.02 3→21 0.0 0.148 4→5 0.0 0.171 6→13 0.0 -0.0141 6→14 0.0 0.07 6→15 0.0 0.213 7→8 0.0565 0.311 7→19 0.0 0.039 8→9 0.0 0.0804 8→12 0.039 0.104 9→10 0.013 0.104 10→11 0.013 0.028 11→20 0.0 0.07 12→13 0.0 0.0852 15→16 0.046 0.144 15→18 0.041 0.114 16→17 0.0316 0.0912 17→22 0.0 0.165
表A2：：发电机参数
节点 Node 定子电 阻Ra 同步电 抗Xd 同步电 抗Xq 暂态电 抗X′d 惯性时间常数 Tj(seconds) 有功Pg 无功 Qg 额定容量SN (MVA) 19(G1) 0.0 0.5 0.1 0.1 11.6 1.8 0.1072 100.0 20(G2) 0.0 1.3 1.0 1.0 12.5 0.0 1.0 100.0 21(G3) 0.0 1.5 1.0 1.0 5.6 0.0 0.055 100.0 22(G4) 0.0 1.6 1.0 1.0 5.6 0.25 0.182 100.0 23(G5) 0.0 1.3 1.0 1.0 15.0 0.735 0.066 100.0 24(G6) 0.0 1.2 1.0 1.0 18.0 1.3675 0.225 100.0
表A3：节点电压
节点 电压(p.u.) 角度.(Radian) 节点 电压(p.u.) 角度(Radian) 1 1.0454620 3.848650E-002 13 9.461755E-001 -2.130444E-001 2 1.0414580 -4.964663E-002 14 9.753382E-001 -2.228975E-001 3 1.0422390 -3.114384E-001 15 9.809795E-001 -7.883766E-002 4 1.0527970 -2.015636E-001 16 9.984322E-001 -9.799489E-002 5 1.0411910 -3.109851E-001 17 1.0214000 -8.097351E-002
6 9.358647E-001 -2.228975E-001 18 1.0258040 -3.533235E-003 7 1.0481540 4.260347E-001 19 1.0500000 4.898637E-001 8 9.443007E-001 -1.401441E-001 20 1.0500000 -3.795634E-001 9 9.503843E-001 -2.387923E-001 21 1.0500000 -3.114384E-001 10 9.569663E-001 -3.666559E-001 22 1.0500000 -4.250140E-002 11 9.834906E-001 -3.795634E-001 23 1.0500000 9.174172E-002 12 9.664204E-001 -1.523288E-001 24 1.0500000 .0000000
表A4：负荷功率
节点 有功 无功 5 0.7 0.32 8 0.492 0.251 10 1.07 0.642 13 1.27 0.9 14 0 -0.55(Capacitive compensation) 16 0.33 0.165
表A5：等值并联容性电纳
节点 等值电纳 节点 等值电纳 1 0.05 4 0.354 5 0.215 7 0.48 8 0.559 12 0.206 16 0.15 18 0.1