（电力系统及其自动化专业论文）电力系统暂态稳定约束关于控制量梯度的新算法.pdf

第一章绪论 1 1 引言 第一章绪论 随着社会的发展文明的进步，电力工业在人类社会成为越来越不可或缺的 一部分。“西电东送、南北互供”以及全国联网“】的实施使得我国的电力系统正 逐步形成以大机组、大电网、超高压、长距离、大区联网和交直流混合输电为 特征的复杂大系统。这在带来一系列经济效益的同时，也使得电网结构和运行 方式也变得日益复杂，使电力系统安全稳定问题变得更加突出，对电力系统的 可靠及运行管理水平提出了更高的要求。而一旦系统发生稳定破坏事故，若未 能及时采取有效的控制措施，则可能波及到其他互联网络中，甚至引发大面积 停电事故，造成严重的经济损失和恶劣的社会影响 2 。7 1 。 同时由于我国电力市场的启动和实施，虽然在整体上提高了电力系统的经 济效益和社会效益，但是这也使得电力系统的动态行为变得越来越复杂，区域 间的功率交换更加频繁，线路的承载度大大提高，系统的运行越来越接近稳定 极限。这些因素都使得电力系统的安全稳定问题面临新的挑战孙1 4 】。 电力系统的互联使得在广阔的地域内进行资源的优化配置，互通有无，相 互支援成为可能。但是，在紧密相连的互联电力系统中，一个局部故障能迅速 向整个系统传播，如果说小型电力系统发生事故造成的停电只是局部的，而这 样的大规模电力系统一旦发生稳定性破坏事故，造成的连锁反应，就会引起大 面积、长时间的停电，后果将是非常严重的。近年来国内外发生了多起停电事 故，如：1 9 9 6 年7 8 月美国西部电网( w s c c ) 接连2 次大停电事故，美国总 统甚至认为停电事故己“危及国家安全”：1 9 9 8 年，台湾省的大停电事故，其 损失电量甚至超过了w s c c 的大停电事故；2 0 0 3 年下半年在北美和加拿大、英 国伦敦、瑞典一丹麦、意大利都先后发生过大面积停电事故，震惊世界。特别 是，2 0 0 3 年8 月1 4 曰美加大停电波及5 0 0 0 万人口的供电范围，造成重大经济 损失，是美国历史上最严重的停电事故。 值得注意的是，这几次大停电事故均发生在发达国家或地区，它们电网的 网架结构相当强壮。这些事实再一次表明了，不论电网发展到什么阶段，都应 第一一章绪论 当将确保电力系统的安全稳定运行作为电力生产的第一要务【4 。4 1 。同时这也表 明了，要使电网大的扰动冲击下，不至于发生一系列的连锁故障而崩溃，仅依 靠一一次系统的网络结构是不够的，还必须在实际运行中加强对系统的动态安全 分析，检验系统在合理预想事故扰动下保持暂态稳定性的能力，并针对事故采 取各种紧急控制措施来提高系统的稳定运行水平【2 _ 3 ，1 5 ，2 3 1 。 当前在我国，由于电力二 业发展相对落后于国家经济的发展，以至于出现 了严重缺电和限电现象，这些问题都是难以在短时间内通过电网建设、备用容 量的增加以及新型装置的投入来解决。因而在电网实际运行中，可以通过合理 地运用相对经济的紧急控制和预防控制等手段来提高系统的稳定运行水平 【2 0 ，2 ”。工程实践表明，预防控制措施在提高系统暂态稳定水平上还是比较有效， 相应的算法也在过去几十年问得到广泛的研究【2 2 _ 5 ”。 预防控制是指在网络参数、故障变化以及负荷水平给定的情况，通过调节 系统的控制变量来合理安排运行方式，以提高系统的暂态稳定水平并维持系统 在预想事故扰动下的暂态稳定，防止局部事故扩大引发大规模的停电事故。在 有些文献中，这种控制措施称为动态安全经济调度或考虑暂态稳定性约束的发 电计划再调整m 】等。从这层意义来考虑，预防控制可以作为一个典型的寻优问 题来处理，而系统运行规划人员面临的问题就是如何快速准确地从运行条件所 确定的可行解集合中寻找一种最优意义上的预防控制措施，在有些文献中称其 为考虑暂态稳定性约束的最优潮流( o p t i m a lp o w e rn o wc o n s t r a i t l e db yt r a n s i e n t s t a b i l i t y ，缩写为0 t s l 【4 5 - 5 ”。 作为寻优问题的预防控制所具有的技术经济意义是传统的最优潮流 ( o p 廿m a lp o w e rf l o w ) 所无法实现的。对于实际系统，由常规o p f 所确定的运行 方式往往无法满足系统在暂态过程中的稳定性要求，因而有必要在进行暂态稳 定性分析的基础上作相应的调整；在阻前由于缺乏有效的分析方法，工程人员 往往是在运行经验的基础上采用试探的方法来完成。但是这种解决方法无法将 系统的经济性与安全性统一起来加以考虑，在电力市场环境中，这种方法是无 法保证市场参与者之间的公平竞争的，甚至还有可能破坏市场的调节机制。在 这种情况下，如何统一协调地来考虑系统运行的安全性与经济性成为日益迫切 的问题【4 2 】。预防控制由于能将系统的经济性以及安全性纳入到同一框架中分析， 第一章绪论 因而成为解决这一问题的有效工具。而现在阶段，有关预防控制的算法还不完 善，至少在实用化方面还存在很多的问题，比如暂态约束形式的处理，相关梯 度的计算，预想故障的处理以及优化问题的规模等，这些相关问题的解决对于 预防控制的实用化都是非常关键的。 同时，电力系统安全控制领域中的许多关键性问题，如紧急控制、最大传 输容量、动态阻塞管理等都同时涉及到电力系统的经济性和安全性两个方面， 都是可以用同预防控制相似的数学模型来表示的，即都可以描述为考虑暂态稳 定性约束的优化问题，对这类问题的数学模型及求解方法的研究有着非常重要 工程意义。求解这类问题的关键在于快速准确的求解暂态稳定约束关于控制量 的梯度，这也正是本文研究的出发点。 1 2 国内外的研究现状 由于电力系统是一个典型的大规模非线性动态系统，其稳定性分析和优化 控制不仅在电力系统方面，即使是在更为广泛的数学规划以及非线性动力学领 域也属于难点问题。为对这种复杂系统进行控制管理，研究者们展开了大量的 工作。对于暂态过程这样一组微分代数方程，怎样准确地描述动态过程的约束， 如何给出这种约束关于控制量的梯度表示，怎样合理地处理预想事故以减小问 题的规模和计算都是很值得探讨的问题，而有些学者甚至开始尝试着将预防控 制同紧急控制结合起来考虑1 6 。8 1 。本文将从以下几个方面综述预防控制的研究 现状。 1 2 1 约束的描述 电力系统预防控制是一个典型的大规模非线性动态系统优化问题，但是在 通常的处理方法中往往是将其转化为一个静态优化问题来处理，在这个优化问 题中，最主要的部分就是如何来处理暂态过程中出现的约束，包括功角稳定约 束，电压约束，甚至发电机出力的上下限约束口。，3 2 l 。这样对暂态问题的不同研 究手段也就相应地产生了对暂态约束的不同处理方式。 在现阶段，考虑暂态问题主要存在着数值积分法，直接法以及这两者结合 第一章绪论 的混合方法，而每种方法又可以根据不同的处理方式加以细分。这样在预防控 制问题中，相应地存在很多种约束的表示方法。对于直接法来说，暂态稳定约 束表示为要求系统在加速过程中产生的能量低于稳定的能量阂值，这个包括 e e a c 、b c u 以及p e b s 等具体的处理方式；对于数值积分法来说，由于对暂态稳 定的描述不存在一个很准确的表示形式，往往是按照工程经验设定在某个时段 内功角之间的差值满足某种关系来表示。 下面将只对基于数值积分法的常用暂态稳定约束处理方法作一个说明。 约束的离散化方法 这种方法主要以文献【3 5 ，4 6 4 8 】为代表。在这种方法中，由于是将整个 暂态过程离散化来处理的，为了表示暂态过程中的满足功角约束，学者们是将 暂态过程中每个时步的各个功角与惯量中心角之差不超过某个限值来表示的， 一般表示如下： 增 h i 赋 群- 气_ 一 一 = 1 ( 1 1 ) ( i = 1 ，2 ，n g ；k _ 1 ，2 ，n e n d ) 其中点。为限值，一般取为1 0 0 度。在这种差分化的处理中，当取的步长 比较小时，产生的不等式约束数目将是相当大的。这样处理暂态稳定约束所形 成的大量不等式约束使得优化计算难以收敛。为了减少不等式的数目，在文献 【4 6 】中约束被简化成下面的形式： hj 或 影- ( 1 2 ) 一 j = 1 同时它还将考虑约束的时间段进一步地加以减小，具体做法是：先计算出 上面的这种差值，然后找到最大差值所在的时间点，在考虑不等式约束时只计 及这个时间点以前约束，而将以后的不等式约束加以忽略。 应该说，这种处理也是合理的，既在一定程度上减小了计算量又有效地控 制了不等式约束的数目。 4 第一章绪论 约束的转化表示方法 这种方法以文献【4 9 5 l 】为代表，其实质是将暂态过程中约束的越限情 况用积分的形式加以记录，以此来表示整个暂态过程中约束的越限程度。这种 处理方法的好处是对分段越限及越限的总体效果有一个比较好的表示。通常用 如下的形式来表示： 五；( i ，歹，手) ；c 【m a x ( o ，带( 戈o ) ，y o ) ，手) ) 2 击 ( 1 3 ) 其中掰为各个时步中约束表达式的值与参考值之间的差值，为正时表示约 束越限，此时在由m a ) 【函数来计及越限的影响；否则彬为负时，表示约束在限 值范围类，则不记录相应的差值影响。这样将整个过程用积分的形式表示出来， 就相应地记录下整个过程中的越限情况。 这种处理的方法能够比较好地满足工程中对越限情况的要求，因为在处理 暂态过程中的电压约束时，要求的往往不是简单的越不越限，而是越限的持续 时间不超过一给定值，这样用积分的形式将约束表示出来就可以将约束越限的 持续时间所产生的效应计算进来。不过这种处理方法也只能是一种近似的处理， 首先是暂态过程的离散表示使得上面的积分只能转化为各个时步越限情况的叠 加；其次工程上所要求的是暂态过程中电压幅值越限的持续时间长度不超过某 个值，而要将这种时间量转化为关于控制量的准确的数学表达是比较困难的， 上面的方法最终也只是将其转化为约束在一定时间段内的越限总和不大于某个 闽值来处理，其实质还是作为约束越限的大小情况来处理的。 约束的极值表示方法 对于经典模型，往往只是用来考虑暂态的第一摆稳定问题，因而相应地， 文献【4 5 】在考虑约束时只是从第一摆是否稳定的角度出发，将暂态稳定约束 作为是在某个指定的终端时刻砟功角是否越限来处理。具体的情况是在终端时 刻功角相互之间的最大差值不超过某个阈值，可以用下面的形式表示： 目( “) = m a x 【点( i “) 一c i ( 瑶i 甜) 】2 ) _ 髭。o( 1 4 ) ( i ，j = 1 ，2 ，n g ) 第璋绪论 其中氏。一般是根据实际经验得到，不过在相关的文献中一般取值为1 8 0 度。 这种形式跟上面的两种约束表示相比，无论是形式上还是规模上都要简单很多。 不过还可以将约束相应地表示为整个暂态过程中的最大功角之差不超过某个闽 值。 应该说，这是由于时域类方法本身对暂态稳定的描述不够准确才导致了约 束在表示上的多样性。要改变这种状况，必须要对暂态过程有一个更深入的认 识和研究。不过，时域类方法同直接类方法相比能够更直接地表示出暂态过程 中除功角稳定约束外的其他约束，更符合工程上的需要。 1 2 2 梯度算法的研究 一般来讲，梯度算法往往是同具体的约束表示联系在一起的。但由于暂态 过程是由一组微分代数方程来描述的，而且各个状态变量之间的存在着耦合关 系，为了能够得到比较准确的暂态过程约束梯度，在求取时往往都要同状态变 量关联起来。对应上面的各种暂态约束的不同表示形式，主要存在着以下的几 种梯度表示方法。 离散化的梯度表示 应该说这种方法对于暂态约束的梯度表示是最为直观和简单的。由于这种 方法是将暂态过程中的状态量进行离散化处理的，故作为等式约束的第k + 1 步 的微分代数方程中只包含有第k 和k + 1 步的状态量和代数量，相应的梯度阵中 只有第k 和k + 1 步的梯度子块不为o ，其余的梯度子块为o ，其中对于相应的梯 度子块具有如下的形式： 爵嚣= 一i 一生iqo p 一ip e p e 00oo o000 磷= 一皇， o o i p e 譬p e 瓮 ogb obg 其中瑶嚣和瑶分别为第k + l 步等式约束对第k 步和第k + l 步状态变量的梯度 矩阵，具体的含义参见式( 2 3 3 b ) 和( 2 3 3 c ) 。 而对于暂态功角稳定不等式约束，由于其中只含有发电机的功角占，这些 6 第一章绪论 梯度关系上并未深入研究。在其文章只是针对一个很简单系统的经典模型，用 龙格库塔方法求相应的梯度。而在采用相同约束处理形式的其他文献中，更多 的是将精力放在了约束的表示和用不同的仿真方法来处理暂态约束，而对具体 的如何对状态量的梯度求解涉及得很少。 由于这种方法要同时涉及到状态变量的初值和控制量，故在考虑梯度关系 时，最好是从时域仿真的角度来考虑这个问题。值得注意的是，在进行仿真计 算后整个暂态过程中的状态量已经确定，相应的各个时刻的梯度关系实际上是 一种确定系数的等式关系，从这个角度来看或许会有些简便的形式。不过由于 状态变量耦合情况的存在，使得在求解约束梯度的过程中需要将所有的状态变 量的梯度求解出来，即使约束只涉及到比较少的状态量。 按照时域仿真中对代数方程的不同处理形式，相应的梯度求解也存在两种 情况： 一种是对应微分方程和代数方程交替迭代的求解过程，同样也可以将梯度 关系按照这种迭代求解的形式来处理。 另外一种是对应微分方程和代数方程联立的隐式梯形积分形式，在时域仿 真的第k + 1 步完成时，如果不考虑收敛误差，将会存在一组差分化的等式关系， 利用这种等式关系同样也是可以推导出控制量的梯度。 对这样的梯度关系的探讨将是本文的一个部分，将在文章的正文部分具体 地介绍。 最优控制原理的梯度表示 应用最优参数选取的原理可求出式( 1 4 ) 所表示的暂态稳定约束的梯度。 这种梯度方法最先是针对切机和切负荷提出的一种紧急控制算法4 4 5 2 】，应用最 优参数选取推导出如下的梯度表示： v 目( 。) ：掣：。f 塑殴地业世监型趔出( 1 6 ) 伽 。 砌 式中，h ( x ，弘“，a ，p ) 为约束函数口( “) 所对应的h 口m i ，o n j 册函数。 日( 而弘甜，五，) = 兄7 9 ( x ，弘“) + 7 中( 知 x 第一章绪论 稳定性的量化指标，对各种事故按这种指标值进行排序，计算系统在最严重 故障下所需要的预防性控制措施，确定出系统相应的新运行点，然后重复稳 定性排序、确定新运行点这样一个过程直到系统在所有的预想事故下能够保 持暂态稳定。这种方法作为多预想事故的处理方法是最先被提出来的，通常 与基于规则的启发式预防控制算法结合，并且它的效果已经在多个实际规模 的电力系统上得到验证。 2 在计算过程中同时考虑多个严重故障的方法口1 3 4 ，4 9 _ 5 0 】。这种方法在处理过程 上与上面仅考虑最严重故障的方法大致是一致的，所不同的是对事故的稳定 性进行量化排序后，需要同时考虑最严重的几个故障，将其加入到预防控制 问题中形成新的松弛问题。这种方法在近期的研究中比较受重视，并与并行 计算技术相结合，取得了很好的效果。而在文献【5 0 】中，作者把这种方法 进一步地加以发挥，将上面所要考虑的严重故障加以区分，即在进行稳定性 排序的过程中将其划分为“起作用故障”和“不起作用故障”，在预防控制 问题中仅考虑“起作用故障”，这种试探取得了比较好的效果。 应该说，这两种方法很难在计算量和收敛性能上进行定量的分析，只是第 一种方法由于每次在稳定性排序 x 第一章绪论 在着协调问题的蓁蠢萋蚕冀羹= 薹鬓荔鍪季一崖弼谭毒哩呵咭 隧曜浇灌蔗膨啄懈舞涛瑰河熨二彳j 蚕惯擞熟鞑朔爿。醴弛憩。瓢琵鸶型w 咎荔掣 态。鼎型铂邕蜀鲞甾g j 雒型牲磐美型驯帮蕉到犁寸镧显毙咎；露弛醵竖 划蛆盟型鲢鞋鞋蝴码鹳篓熊积翻墅醪墅蓦曳! 黔酏沁娶萋楚垂跫f 剜蹦踅妻阿 静j 粕藿整蕃爨毯必醛g 如珀孙羁；酏趔业鳞i | l 和歉飘科i 鬻发。 由于暂态过程中的状态变量x 是由( 2 3 ) 式的微分代数方程所决定，即不 存在x关于控制量z ( 或者说u ) 的显示表达式，给梯度的计算带来了非常大的 困难。在约束转化方法的文献中提到重新用积分的方法求状态变量的梯度4 q 5 们， 不过对于具体的形式没有表示；而差分化方法在梯度的表示上相对来说要简单 很多，但是这种方法，即使是不考虑由约束和变量过多所引起的收敛问题，在 优化计算的处理上也不占优势。 既然考虑问题是从时域仿真的角度出发，不妨利用时域仿真的特点来处理 状态变量x 关于控制变量u 约束的梯度。由于在时域仿真中两个连续时步的状 态量之间存在着差分化的等式关系，如在采用隐式梯形积分方法进行仿真时， 可以利用这种关系由前一个时步状态变量关于控制变量的梯度推导出下一个时 步状态变量的梯度。同时，由于状态变量的初值对整个暂态过程同样也是存 在影响的，故在考虑梯度关系时也应计及状态变量初值对梯度的影响。这样对 应终端时刻珏的状态变量x ，关于控制变量u 的梯度可以表示如下： 挈：挈罢+ 挈誓 ( 2 - 8 ) d“ozd ”d x 。d “ 要得到上面的梯度，可以按照下面两个步骤来进行：首先利用仿真过程中隐式 积分的等式关系来得到状态变量x 关于z 和而的梯度善和善，然后利用初 以 蕊。 值方程( 2 4 ) 推导出z 和关于控制变量u 的梯度娶和粤。这样，进一步就 可以很容易地得到式( 2 5 ) 所示的功角稳定约束的梯度。下面将按照这种思路 x 第一章绪论 制需要计及故障发生前的系统状态情况，而紧急控制则只需要针对故障清除后 的系统状态。 这样看来，将预防控制和紧急控制协调在一起来考虑将会是一个混合型的 非线性规划问题，需要从全局优化的角度来处理这个稳定控制决策问题。其实 文献【1 6 】在比较早的时候就提出了一种解决的思路，并给出了协调两者的数 学模型和算法，相应的稳定控制决策空问包含了连续的预防控制子空间和离散 的紧急控制子空间，故两种控制措施的协调优化必然会在连续空间和离散空间 中进行迭代。在其后续的文献【1 8 】中对相应的算法也做了探讨。 但是上面的文献，对预防控制与紧急控制的协调处理只是分阶段来考虑的， 两者之间的协调只是两个阶段各自的最优情况，而不一定是整体情况的最优。 将两者协调作为一个准确的全局最优化问题来解决，还存在着两个问题：需要 对预防控制代价和紧急控制代价的给出一种统一的衡量标准，需要将概率事件 合理地引入目标函数。但是还要注意，这个问题的解决必将涉及到比较多的可 供选择的措施，如果没有可靠的量化指标及相关的灵敏度分析技术，以及一个有 效的搜索策略，那么对这种优化控制问题进行求解将是非常困难的。 1 3 本文的内容及章节安排 本文从时域仿真的角度对预防控制的算法进行了研究，工作主要集中在两 个方面：对暂态约束形式的选取和状态变量梯度的求解。为了保证系统的安全 运行，工程上对包括暂态功角稳定性在内的系统暂态性能提出了各种要求，怎 样用一个合理而又便于优化方法进行求解的表达式来描述暂态约束一直就存在 着多种说法；其次在状态量的梯度求解方面，众多的有关预防控制的文献往往 将重点放在了对暂态约束的描述以及如何来处理预想故障上，对梯度的算法涉 及得相对较少。而预防控制算法要真正地得到改善，暂态过程约束梯度的处理 是关键的一个环= 爷。本文围绕着这两个主要问题，开展了大量的研究工作，主 要包含以下几个方面： 1 在对比几种常用的功角稳定约束表示后，参考文献【4 5 】的约束表 第二章预防控制的约束表示与梯度求解 2 1 数学模型 对于常规的最优潮流其数学模型一般表示为如下的形式： m i n 厂( ，“)( 2 一l a ) s t g ( ，“) = o( 2 - 1 b ) ( ，甜) 0( 2 - l c ) 其中，u 为系统的控制变量，一般包括系统中发电机节点的有功功率( 参考节 点除外) 和电压幅值等；y 。为过程变量，一般包括发电机节点的无功出力、非 发电机节点的电压幅值和所有发电机节点的相角等：g 为潮流等式约束；h 为 不等式约束，包括各个节点的电压幅值的上下限约束，发电机有功出力和无功 出力的上下限约束，线路的传输容量约束等。 上面的模型将变量划分成控制变量和过程变量是与实际的工程情况相符 的，不过由于这样的处理在表示相关梯度时不太直观，故本文为了分析方便， 将上面出现的变量都作为u 来处理，这个虽然与实际的工程情况不相符，但是 这对整个优化问题的分析是没有影响的。同时，这样做也便于后面进行结论分 析。其表示形式如下： m m 厂 )( 2 - 2 a ) s t g ) = 0( 2 2 b ) h ) 兰0( 2 - 2 c ) 在最优潮流中考虑暂态稳定问题时，一般要对如下的微分代数方程进行处 理，或者将其作为等式约束引入，或者从这种等式关系中推导出暂态约束同控 制变量u 的关联： 戈5 f ，y ，z )( 2 - 3 a ) g ( x ，弘z ) = o( 2 3 b ) 其中，( 2 - 3 a ) 式为转予 x 第二章预防控制的约束表示与梯度求解 ( 其中i ，j = l ，g ；最表示第k 步的功角) 其中，吒。的选择根据工程经验来得到，同时为了考察算法的收敛情况也可以 人为地进行设定。相似地，这里的功角约束日还可以按照相对于惯量中心角的 差值的形式来表示。 应该说式( 2 5 a ) 的约束表示同式( 1 1 ) 和式( 1 3 ) 两种表示形式从意义 上来说是一样的，即在整个暂态过程中功角约束都不越限。图2 1 和图2 2 表 明了这种形式的暂态稳定约束的有效性，两个系统的数据参见第四章的数据说 明：( 注，以下图形是各个功角相对惯量中心角的差值表示) 图表2 - 1l 旺e - 3 9 算例6 m x 取1 5 0 度时得到的0 t s 仿真曲线 在这两个图示中，功角约束的限制效果得到了比较好的反映，其中图2 一l 在1 2 6 秒达到功角差的最大值1 4 9 9 2 度，此时的6 m a ) ( 取为1 5 0 度；而图2 2 则在o 4 6 秒达到了最大功角差1 1 9 9 4 ，相应的6 m a x 为1 2 0 度。这两个例子结 果都能很好地符合约束限制的要求。 另外这样处理约束还有一个好处，就是在求解状态变量梯度的过程中不需 要将所有时步的状态变量x 关于控制变量u 的梯度计算出来，而只需要进行到 最大功角差所在的时刻，后面时刻的状态变量由于没有在约束中出现也就不需 要计算相应的梯度。从这两个图示的最大功角差的位置来看，这对减小梯度求 兰三! 堡堕丝型塑丝墨塞重量塑堕查堡 一一一 来推导相关的梯度关系。 2 3 1 状态变量x 关于z 和x 0 的梯度 在采用隐式梯形积分法进行时域仿真时，由( 2 3 ) 式司知，赡+ 1 是由坼拨 照下面的等式关系得到的： n h + 1 = 嘶+ 三( e + 1 + e ) ( 2 9 a ) g ( x ，儿十1 ，z ) = o ( 2 - 9 b ) 其中h 为积分步长。即在时域仿真得到第k + 1 步的状态变量坼+ 。后，可以认为 上面的等式关系是成立的。虽然实际计算时，上面的等式关系还存在着一个容 许误差，但是这对最后得到的梯度的影响是非常小的，这一点可以从后面的数 据分析中得到说明。 在暂态过程中，由于x 是由z 和所决定，不妨将x 看作z 和函数。这 样，用z 和对( 2 - 9 a ) 式的两边求导： 警= 誓+ 害c 譬+ 争 陋， 0 za z2 、d z d z 。 錾：磬+ 鲁( 孥+ ( 2 _ 1 0 b ) 2 、氐 氐。 上面的表示在第k + 1 步和第k 步状态变量关于z 和的梯度之间建立起等式关 系，如果上面表达式f 的偏导关系能够得到，则由第k 步的梯度拿和鍪就能 0 2 够推出第k + l 步的梯度冬丛和安生。 0 2 由( 2 3 a ) 式可知，f 是关于x 、y 和z 的函数，而y 又是x 和z 的函数， 则f 对z 和的偏导分别为： 堡：f 堡十堡鬈l 鬟摹秦鏊毒薹萋 x 第二章预防控制的约束表示与梯度求解 其中，婴和粤的偏导形式可以由潮流方程( 2 3 b ) 加以推出；由隐函数 d x0 2 存在定理可以得到下面的关系，这种方法其实在有关常规最优潮流的简化梯度 算法中就有体现，形式如下： 要+ 票罢：o ( 2 1 2 。) 僦 钾d 譬 誓+ 要害：o ( 2 - 1 2 b ) 洲0 2 这样，将式( 2 1 1 ) 代入( 2 - 1 0 ) 式中，就可以得到梯度第k + 1 步状态变 量的梯度专 和豢 同第k 步状态变量梯度冬和晏的明确关系： o z o x n 0 z 硪、 e “誓：誓+ 鲁( 譬+ e 十1 ) ( 2 - 1 3 。) 1 出色2 、出 。7、一7 砖“譬：导+ 鲁磐 ( 2 _ 1 3 b ) 4 2 氐 、一7 其中，e “和“的具体形式为： 小要争一冬警磐1 4 a ) 2 + - 2 饥+ ，+ 。 、一 ( 2 1 4 b ) 这样黾+ 。和稚关于过程控制量z 和初值的梯度关系就在( 2 1 3 ) 式中确定了， 又由于求梯度的第k + 1 步与第k 步是与时域仿真过程相对应的，则( 2 一1 3 ) 式 中的系数是可以按照( 2 1 4 ) 式来得到的。这样看来，现在最主要的是要得到 梯度计算的初始值，即得到誓和鍪，这样就可以推出暂态过程中与时域仿真 【z 0 相对应的任意时刻的梯度。 由暂态计算的初值方程( 2 - 4 ) 可知，和z 都是由u 所决定的，是两个互 相独立的变量，相互之间不存在偏导关系，即： 堕：n 瑟 ”g 。 ( 2 - 1 5 a ) 盟瑟 鳖瑟 盟 = 呼 第二章预防控制的约束表示与梯度求解 堕：i ( 2 1 5 b ) 其中o 。g n g 为n g 阶的o 矩阵，i 为n g 阶的单位矩阵。这样在迭代计算的初始值 已知的情况下，由( 2 1 3 ) 式就可以通过迭代计算得到终端时刻的状态变量x ， 关于z 和的梯度竽和善。 可以看出，在上面梯度关系的推导过程中没有涉及到具体的网络变化，这 个不难解决，只需要在用式( 2 1 2 ) 进行计算时按照相应的导纳矩阵来处理即 可。同时，与时域仿真的状态量的求解相对应，在得到第k 步的状态量的梯度 后，如果在第k + 1 步发生网络变化，只需要在( 2 1 3 ) 式中按照变化后的导纳 矩阵来处理相应的系数。 2 3 2 状态变量梯度的另一说明 其实式( 2 1 0 ) 的这种梯度关系可以通过对( 2 3 a ) 式两边直接求导，然后 按照隐式梯形积分的形式展开来得到。在这里暂时只考虑用z 求偏导的情况， 对的情况相似处理，简要推导如下。 用z 对( 2 3 a ) 两边求导，可有： 婺：掣 ( i i - l a ) 由于z 在暂态过程中不是时间变量，则有 嬖：喜 ( 1 1 _ 1 b ) 恐 a z 铁 可以将( i i 一1 a ) 表示如下： 竺：掣 ( i i _ 1 c ) 8 z a娩 为了表示方便，不妨先将( i i - 1 c ) 式中的表达式用下面的符号来代替： x 2 = 鍪 ( i i 一2 a ) ，：娑 ( i i _ 2 b ) 进而将( i i - l c ) 式按隐式梯形积分的形式展开，对第k + l 步有： 第二章预防控制的约束表示与梯度求解 + l = + 兰( 巧+ 磺1 ) ( i i - 3 a ) 即： ( 鼽+ 1 _ ( 鼽+ ； ( 瓢+ ( 瓢+ 1 ( i i _ 3 b ) 而( 2 - l o a ) 式的形式为： 纽：盟+ 生f 盟+ 鸟 a z出2 出沈。 可以看出，两者的主要差别在于下标所处的位置不一样，当然表示的意义 也不一样。但是，如果对两个计算表达式来说，有初值和冬相同，巧和孕 的展开形式和系数相同( 最，和竺霉的情况相似) ，应该就可以说明这两个表 达式等价。 首先来看磊和誓的关系：在计算时候，作为计算的初始值和拿都为o ， 出 “ 出 这可以肯定。 而对于巧，按照f 对z 的导数关系来看： 巧：f 罢+ 罢罢+ ( 罢+ 罢罢) 妻1 ( i i 4 a ) 出 砂出 m 砂出出j t 由于在计算巧时，如果不涉及到具体的某个时刻，状态量z 是未知的，相应的 式( i i - 4 a ) 中的系数也无法得到。因而在上面的梯度计算中，下标必须与时域 仿真的第k + 1 步和第k 步相对应，否则计算将无法进行。同时，这也表示得到 的结果是第k + 1 步状态变量在当前数值情况下的梯度，才是符合计算要求的； 而在处理约束过程中，如果约束与对应的梯度不匹配，优化计算的结果是会有 很大的差别，甚至会产生收敛问题，多次的试验计算表明了这一点。 在与时域仿真对应时，( i i 4 a ) 式变成如下的形式： 野：誓+ 婺誓+ ( 誓+ 要挈) 冬( i i 舶) o z k 0 2 出k k 。x k 而( 一4 b ) 式的右边正好是孕的展开形式，表明巧和孕在形式和系数上是 船出 等价的。这样，按照迭代计算得到的和善就是相同的。也就表明( i i 3 b ) 和( 2 1 0 a ) 两式在计算状态变量的梯度上是等价的，进一步表明了，按照( 2 1 3 ) 第二章预防控制的约束表示与梯度求解 式来计算功角的梯度是具有比较高的精度的，在论文的后面部分将通过具体的 数据来予以说明。 2 3 3 z 和关于控制变量u 的梯度 在上面推导出状态变量x 关于过程控制量z 和初值梯度后，为了得到关 于u 的梯度，如式( 2 - 8 ) 所示，还需推导出过程控制量z 和初值关于控制变 量u 的梯度昙和拿。 d 材d “ 考虑初值方程( 2 4 ) ，用所有变量对等式两边求微分，可以得到： 譬如+ 罢氐+ 譬如：o ( 2 - 1 6 ) 船c 瑚 在当前的数值情况下，出和出。是一组关于幽的线性方程，求解( 2 1 6 ) 式就 可以得到当前数值情况下z 和关于u 的梯度昙和拿。 d “a w 这只是一个求解这种梯度转换关系的思路，对具体的发电机模型而言， ( 2 4 ) 式的形式和变量都不相同，因而在梯度转换的计算上也会有不同的处理。 下面将按经典模型来加以叙述，其他的模型在按照相似的方式来处理。 经典模型下的梯度表示 在经典模型下，z 中包含的变量为五+ 和己，x 中包含变量占和国。相应地， 初值方程( 2 - 4 ) 的具体形式为： 耻耻玻+ 竽x ： ( 2 - 1 7 a ) ，= o ( 2 1 7 b ) 只。= 弓 ( 2 一1 7 c ) ( 其中的下标i 表示第i 台发电机) 在( 2 1 7 a ) 式中只存在着瓯，和e 两个变量，则其关于u 的梯度关系在这 个式子中可以确定，先将式子的实部和虚部分开，可得： 巧2 一g c o s ( 瓯，一q ) + 矗纬= o ( 2 1 8 a ) 第二章预防控制的约束表示与梯度求解 e ks i n ( 磊，一只) - 氍弓= o ( 2 1 8 b ) 对上面的两个式子按照全微分的形式展开，整理可得： f c 0 蛾一曰) 谢国一目) 1 f 趔1 一f o2 一e 峨。曰) 叫巧s i i 蛾q ) l 巧姐国曰) 驾s ( 国一日) 八弼i 五 s i ，( 瓦一目)曩巧c 0 毛一句) 镌 ( 2 。1 9 a ) 对于上面的式子，由于倒和d 磊，前面的系数矩阵的逆存在，并且等式的右边为 幽的展开形式，故由( 2 1 9 a ) 可以得到矗和皖。关于u 的梯度。 另外，由( 2 1 7 b ) 和( 2 1 7 c ) 可以得到如下的梯度转换关系： 挚= o( 2 _ 1 9 b ) 堡：堡 抛劫 ( 2 1 9 c ) 其中( 2 1 9 b ) 为0 向量，( 2 1 9 c ) 只有一个对应元素为l 其他元素为0 。这样 在分别得到z 和各个分量关于u 的梯度后，最后通过整理就可以得到经典模 型下整个的梯度矩阵形式。 2 3 4 功角稳定约束的梯度 在上面两节中推导出计算需要的梯度! 、竽、昙和粤，这样就可以 d ”h 根据( 2 8 ) 式来得到终端时刻珏状态变量工，的梯度。 考虑形如( 2 5 a ) 式的功角稳定约束，由于约束中只出现了两台发电机的 功角，则其梯度在计算时只需要! 翌中的两个对应的量，则相应的梯度表示为： 筹- 2 ( j ；一咧竽一筝 ( 2 - 2 0 a ) d “ 。 d “d “ 而对于形如( 2 5 b ) 所示的功角稳定约束，不妨假设在第k 步达到最大功角差， 则相应的状态变量的下标为k ，并且由于约束中只包含这个时步的状态变量， 故相关的梯度计算只需要计算到这一步。其梯度表示为： 第二章预防控制的约束表示与梯度求解 娑：2 ( 一倒) ( 娑一掣) ( 2 - 2 0 b ) d “d “ d “ 至此，得到了功角稳定约束关于控制变量u 的梯度。 这样在2 1 节至2 3 节中就将优化计算所需的数学模型、暂态稳定约束以及 暂态稳定约束的梯度全部表示出来，可以根据需要选择相应的优化方法来考虑 对上面的o t s 问题进行求解。 下面这个小节只是对基于隐式积分形式的相关梯度求解方法一个探讨。 2 4 暂态约束梯度的其他形式 由式( 2 3 ) 可以看出各发电机和负荷只通过网络相互影响，它们之间无直 接联系，因此微分方程式( 2 3 a ) 在各个发电机和各个负荷感应电动机之间没 有直接耦合关系。利用这个特点，工程应用的时域仿真计算往往是采用交替求 解的方法来处理式( 2 3 ) 所示的微分代数方程。与此对应，也可以将这种方法 用来求状态变量的梯度，不妨将其称为电压的附加处理方法。而在上面2 3 1 节 中介绍的梯度方法由于是将暂态过程中的电压作为控制量z 和状态变量x 的函 数来处理不妨称之为电压的函数处理方法。而与此类似，也可以将电压作为变 量来处理，即将电压变量和状态变量联立起来，统一看作x ，这种形式与将动 态方程和网络方程联立的隐式积分是相对应的；这样虽然使得梯度的维数增加 了很多，但是可以看到相关的矩阵具有很高的稀疏度，另外在梯度的表示上也 会简单很多，不妨将其称为电压的变量处理方式。 下面将分别对这两种求梯度的方法进行说明。为了简化表达形式，说明将 只针对经典模型进行。 2 4 1 电压的附加处理方法 这种方法在处理梯度时，所用的方法是与相应的仿真方法是一样的，即通 过微分方程来得到状态变量梯度的预测值，而通过网络方程来对电压部分的梯 度加以修正，进而代入到微分方程中得到修正的状态变量的梯度，这样交替迭 第二章预防控制的约束表示与梯度求解 代直到满足收敛条件。以下将按照这种过程给出求解形式。 在经典模型下，进行暂态计算的所用的动态方程和网络方程分别如下。 动态方程：( 注以下的符号表示在不作说明的情况下按照向量或矩阵处理) a j 一= ，p 卅p 。 ( 2 - 2 1 ) a p m p e 、1j a 弓 网络方程 旧 其中 撒心 觑：型 朋 秒一等 g 和b 分别为导纳阵的实部和虚部。 ( 2 - 2 2 ) ( 2 2 3 a ) ( 2 2 3 b ) ( 2 2 3 c ) 在时域仿真时，由下面的关系来得到第k + l 步状态变量的修正量，其中计 算是的差值和j a c o b i a n 矩阵分别如下： 匕鼽习= 一瓯一言( + q ) 婢+ 。一婢一害( 2 砌+ 1 - ) 婢+ i 一婢一瓦( 2 砌+ 1 - ) ( 2 2 4 ) 其中，i 为n g 阶单位矩阵；船( d ) 为一对角矩阵，其对角元由如下向量组成 。= 隽( 啪一哪) 进而得到第k + 1 步状态变量的预测值： ( 咎+ ( 翱 ( 2 2 4 a ) ( 2 2 5 ) 然后通过网络方程( 2 2 9 ) 来修正发电机的机端电压。这样进行交替迭代计算 以达到收敛要求。 一 = 第二章预防控制的约束表示与梯度求解 2 5 本章小结 这一章给出了进行预防控制计算的0 t s 模型，并对暂态稳定约束的形式进 行了探讨，结合相应的计算结果对其合理性和有效性加以说明；接下来的部分 针对这种形式的约束给出了一种新的考虑状态变量梯度的方法，并进一步对相 关的梯度求解方法进行说明。关于这种约束的处理方法和梯度求解方法的有关 结论将会在第四中结合具体的算例进一步地加以说明。 第三章梯度算法的分析及计算量的比较 第三章梯度算法的分析及计算量的比较 3 o 引言 上一章的梯度算法，在负荷处理上还存在需要修正的地方：例如将负荷在暂 态过程中作为恒定阻抗处理时，每次负荷的等效计算是在得到新的控制变量“川 后进行，这样由于控制变量“。中电压幅值v 的变化必然导致每一一步循环中负荷 等效阻抗不一样，从这种意义上来说导纳矩阵也是控制变量u 的函数；因而在考 虑暂态过程中状态变量对控制变量u 的梯度时需要计及这种影响的。只是在第二 章中为了将整个求梯度的思路表达清楚而没有对这个问题进行讨论，同时对负荷 影响的考虑应该根据不同的负荷模型来具体处理。 另外在本章中还将对这种梯度求解方法的计算量做一个分析，并与差分化的 方法的计算量做一个分析比较。从某种意义上来说，第二章的梯度求解方法是差 分化方法的退化形式，两者在计算量上差别是不大的，只是这种算法在表示形式 和内存的占用上具有优势。 3 1 负荷等效对梯度的影响 正如在引言中所说的那样，负荷等效阻抗的计算是在每次得到新的控制变量 u 后进行的，这样必然导致每一步循环中的导纳矩阵是变化的，可以将导纳矩阵 看成控制变量u 的函数。相应地只有在梯度的计算中考虑到这种影响才会使梯度 的计算更准确。否则，负荷节点电压幅值的梯度将与实际情况有较大的差异；而 当发电机节点存在负荷时，由于相关梯度准确性的影响，将会对整个优化算法的 收敛有着比较大的影响。 这种差异性可以从下面的数据中看出。对3 机9 节点系统【5 3 】用2 _ 3 节的方法 得到的梯度如下： 表格3 - 13 机9 节点系维不考虑负荷影响时发电机节点部分的梯度 第三章梯度算法的分析及计算量的比较 p g ( b u s 2 )一8 7 9 8 6 8 4 5 7 64 6 8 3 5 8 4 9 7 1 2 55 8 2 9 9 5 5 4 p g ( b u s 3 ) 0 5 4 9 7 0 9 3 82 2 9 9 1 6 6 7 1 21 4 2 5 3 2 7 7 5 4 q g ( b u s l ) 0 8 3 8 鹤8 3 50 5 3 1 1 2 0 5 2 50 5 0 0 9 0 5 2 1 2 0 9 ( b u s 2 )0 8 8 5 0 8 7 0 9 74 0 0 5 5 4 7 8 8 82 1 3 4 2 4 0 0 1 9 0 9 ( b u s 3 ) 15 7 1 0 0 8 4 9 5 6 1 0 5 4 5 8 4 1 l2 8 2 8 8 9 3 7 5 5 v ( b u s l )一1 4 ( ) 8 8 5 1 7 3 28 8 1 2 3 2 9 3 1 28 2 7 8 8 3 0 3 8 v ( b u s 2 j7 7 5 2 3 2 4 1 2 83 4 5 9 7 3 3 7 8 71 8 5 3 7 1 1 9 3 5 v ( b u s 3 )8 4 2 8 0 4 0 9 l一3 3 4 8 4 9 1 6 1 9一1 5 6 9 4 6 5 9 3 8 e ( b u s l )0 8 0 1 4 1 6 2 8 32 3 5 6 0 7 3 3 4 33 0 6 2 6 6 3 8 8 7 0 ( b u s 2 )2 8 2 7 6 5 7 1 5 21 0 0 9 7 6 8 5 7 45 9 5 1 2 6 1 7 2 6 e ( b u s 3 )3 0 2 6 2 4 0 8 6 96 7 4 1 6 1 2 3 9 9一1 8 8 8 5 9 7 8 3 9 相同的数据，由差商代微商的方法得到的梯度如下表： 表格3 - 23 机9 节点系统由差商代微商得到的发电机节点的梯度 6 对u 梯度6 ( b u s l )6 ( b u s 2 )6 ( b u s 3 ) p g ( b u s l ) 1 1 4 8 8 4 6 1 2 26 9 4 6 1 4 0 6 4 o 0 1 9 2 8 3 5 5 l p g ( b u s 2 )一8 7 9 8 6 7 3 3 5 74 6 8 3 5 7 7 5 4 32 5 5 8 3 0 1 1 2 8 p g ( b u s 3 ) 一05 4 9 7 1 4 9 6 6 2 2 9 9 1 6 7 4 5 21 4 2 5 3 2 9 8 9 6 q g ( b u s l ) 一0 8 3 8 6 8 8 3 9 805 3 1 1 2 0 5 0 6一o 5 0 0 9 0 5 1 9 8 q g ( b u s 2 ) o 8 8 5 0 8 7 0 0 3 4 0 0 5 5 4 8 6 6 82 1 3 4 2 3 9 6 0 1 0 e ( b u s 3 ) 1 5 7 1 0 0 8 5 861 0 5 4 6 0 6 8 5 2 8 2 8 8 9 0 7 2 2 v ( b u s l )一1 4 0 8 8 5 2 9 6 78 8 1 2 3 2 2 0 3 38 2 7 8 8 3 6 5 9 6 v ( b u s 2 )7 7 5 2 3 1 8 4 2 93 4 5 9 7 3 9 5 5 81 85 3 7 0 9 2 5 9 v ( b u s 3 ) 8 4 2 8 0 4 4 1 l 一3 3 4 8 4 9 8 6 0 2一1 5 6 9 4 5 6 9 0 8 0 ( b u s l )0 8 0 1 4 1 5 2 0 12 3 5 6 0 6 9 3 4 5 3 0 6 2 6 5 9 8 8 6 o ( b u s 2 )一2 8 2 7 6