一种基于矢量化技术的特征值解析灵敏度计算方法

技术领域

本发明涉及一种电力系统小干扰稳定分析技术领域，尤其涉及一种基于矢 量化技术的特征值解析灵敏度计算方法。

背景技术

随着跨区域互联电网规模的不断扩大和新能源渗透率的不断提高，电力系 统出现低频振荡的事故越发频繁。为了抑制低频振荡，提高系统阻尼，目前常 采取三种措施：电网规划、加装阻尼控制器和发电重新调度技术。其中发电重 新调度技术由于周期短、成本低、可操作性强，已成为抑制低频振荡的最重要 手段。

研究发现，电力系统发生低频振荡主要跟弱阻尼振荡模式相关，如果能够 找到振荡模式阻尼对系统中发电机出力的灵敏度，即特征值对发电机出力的灵 敏度，就能够利用它进行发电重新调度，这样便能够在短时间内提高电力系统 的稳定裕度，从而达到抑制低频发生的目的。

目前求解特征值对发电机出力灵敏度的常用方法为基于摄动理论的数值微 分法，该方法物理概念清晰、实现过程简单，缺点包括计算精度、效率不足， 平衡节点灵敏度、PV节点无功灵敏度无法求取等。例如，一种在中国专利文献 上公开的“确定电力系统稳定器抑制低频振荡效果的方法及装置”，其公告号 CN105375497B，包括：在电力系统稳定器所在的发电机的功角保持恒定的情况 下，获取所述发电机的励磁参考电压至电磁转矩的频率响应特性曲线；计算不 稳定特征值关于电力系统稳定器增益的灵敏度；根据所述频率响应特性曲线或 所述灵敏度整定所述电力系统稳定器的参数；利用整定后的电力系统稳定器的 参数进行仿真，得到电力系统的闭环特征值。

随着电网规模的不断扩大，计算特征值灵敏度的计算时间大幅增加，为了 满足在线重新调度要求，必须加快特征值灵敏度的计算效率，同时在优化模型 求解过程中特征值的灵敏度精度也同样会影响到模型的收敛性。因此，开发一 种既能大幅度提高计算效率又能提高计算精度的特征值灵敏度方法具有相当大 的必要性和紧迫性。

发明内容

本发明主要解决现有技术计算精度、效率不足的问题；提供一种基于矢量 化技术的特征值解析灵敏度计算方法，利用复合函数求导法则推导出特征值对 发电机有功出力的一阶表达式；结合电力系统实际运行规律，对一阶表达式作 相应的处理，最终获得特征值对发电机出力的解析灵敏度。既能大幅度提高计 算效率又能提高计算精度。同时解决了平衡节点灵敏度、PV节点无功灵敏度无 法求取的问题。

本发明的上述技术问题主要是通过下述技术方案得以解决的：

本发明包括以下步骤：

S1：采集电网实时运行数据以及各元件参数；

S2：建立矢量化的稳定分析模型；

S3：分别计算特征值对系统变量的导数和发电机出力对系统变量的导数， 利用复合函数求导法则得到特征值对发电机出力的导数；

S4：结合电力系统的实际运行规律获取特征值对发电机出力解析灵敏度；

S5：建立含1-范数目标函数的发电重新调度模型进行有效性和高效性的验 证。

电网实时运行数据以及各元件参数主要包括潮流计算的网架结构参数、发 电机出力、负荷大小以及包括发电机、励磁机在内等动态元件的数学模型；本 方案首先获取电力系统网架结构数据并利用牛顿迭代法算出潮流；然后根据发 电机、励磁机等动态元件的数学模型构建相应的基于矢量化技术的微分代数方 程并计算系统各变量初值；利用李雅普诺夫线性化技术获取系统的状态矩阵， 求出特征值及对应的特征向量。最后，利用矢量化技术和复合函数求导法则求 取特征值对发电机出力的导数并作相应处理后获得特征值解析灵敏度。本方案 利用复合函数求导法则和矢量化技术推导出灵敏度解析解以提高灵敏度计算精 度和计算效率，为抑制低频振荡的发电重新调度方案提供关键的参考信息。

作为优选，所述的矢量化的稳定分析模型表示如下：

消去运行参数得到：

得到状态矩阵的表达式：

其中，x为状态变量，y为运行参数变量，μ为输入变量；

Δx为状态变量；

E

Δμ为输入变量。

建立矢量化的小干扰稳定分析模型，并计算潮流初值和系统变量各初值； 建立矢量化的小干扰稳定模型，最重要的是确认电力系统采用的包含励磁机、 发电机在内的动态元件数学模型，分别对发电机模型、励磁器模型和网络模型 进行建模在对电力系统元件线性化，方便求解。

作为优选，所述的矢量化的稳定分析模型定义一个称为谱横坐标的参数η来 描述系统的安全裕度，所述的安全裕度为：

Re(λ

其中，λ(A)表示状态矩阵A的所有特征值；

Re(λ)表示特征值λ的实部；

λ

把谱横坐标参数作为发电重新调度的稳定指标，在简化模型同时还能够降 低调度执行的时间。此外，谱横坐标极限值应用于优化调度中，可以保证优化 结果中的特征值整体向虚轴左半部分偏移，因此和传统使用阻尼比指标一样具 有反应系统小干扰稳定性的功能。

作为优选，所述的步骤S3包括以下步骤：

S31：特征值对系统变量的求导；

其中，

λ

X为系统变量；

S32：计算发电机出力对系统变量的导数

包括构造潮流方程、网络代数方程以及稳态下的微分方程组；并对系统变 量求导，求解得出发电机出力对系统变量的导数

S33：利用复合函数求导法则得到特征值对发电机出力的导数；

其中，P

δ是矢量化转子角度；

ω是矢量化转子速度；

V为矢量化节点电压幅值；

θ为矢量化节点电压相角；

E′

I

E

V

R

特征值对系统变量的求导的左右特征向量经过归一化处理，简化计算，通 过消掉逆矩阵的偏导形式并把三维矩阵运算转化为二维矩阵运算后得到，提高 计算效率和简化算法的开发难度。

作为优选，所述的步骤S32中构建的方程如下：

完全极坐标形式潮流方程：

网络方程极坐标形式：

定子代数方程的极坐标形式：

取自与系统变量求导相关的稳态下的微分方程组如下：

其中，符号

G

G

G

P

P

Y表示为系统节点导纳矩阵；

X

X′

V

R

diag(V)表示把列向量V转化为对角稀疏矩阵；

E

发电机出力对系统变量的求导方法不需要对矩阵整体求逆，只需要对单个 变量求逆即可，求取方便。

作为优选，所述的步骤S4包括以下步骤：

S41：在步骤S32求导过程中对含自动发电控制ACC的发电机组进行修正， 即对连接含AGC的母线k的偏导数设置为零，具体为

S42：在步骤S32求导过程中对安装有自动化电压调节器AVR的发电机组进 行修正，即含AVR的母线j的偏导数也需要设置为零，具体为

V

V

利用复合函数求导法则得到特征值对发电机出力的导数后需要结合电力系 统的实际运行规律作相应处理，才能获得特征值对发电机出力灵敏度的解析表 达式。

作为优选，所述的1-范数目标函数为：

约束条件为：

其中，

σ

η

S

ΔP

上述方案考虑了平衡节点的灵敏度与PV节点的无功灵敏度。

作为优选，对所述的1-范数目标函数调整，调整后的模型为：

约束条件为：

其中，发电机组出力的调整量ΔP

正调节量

调整后的模型维持1-范数目标函数的优势，又不使其呈现非光滑性。

本发明的有益效果是：

1.在灵敏度的计算过程中去掉了矩阵求逆和三维海森矩阵的运算，计算效 率得到大幅度提升，

2.优化模型为简单的线性规划模型，迭代次数少，收敛速度快。

3.采用含1-范数目标函数的发电重新调度模型进行有效性和高效性的验 证，考虑了平衡节点的灵敏度与PV节点的无功灵敏度。

附图说明

图1是本发明的一种特征值解析灵敏度计算方法的流程图。

图2是本发明的一种WSCC3机9节点系统结构拓扑图。

图3是本发明的一种新英格兰10机39节点系统结构拓扑图。

具体实施方式

下面通过实施例，并结合附图，对本发明的技术方案作进一步具体的说明。

实施例：

本实施例的一种基于矢量化技术的特征值解析灵敏度计算方法，如图1所 示，包括以下步骤：

S1：采集电网实时运行数据以及各元件参数。

电网实时运行数据以及各元件参数主要包括潮流计算的网架结构参数、发 电机出力、负荷大小以及包括发电机、励磁机在内等动态元件的数学模型。

S2：建立矢量化的稳定分析模型，计算潮流初始值和系统变量初始值。

建立矢量化的小干扰稳定模型，最重要的是确认电力系统采用的包含励磁 机、发电机在内的动态元件数学模型。

本实施例以励磁机采用IEEE DC-1型、发电机采用双轴模型为实施例进行 分析说明。

小干扰稳定模型表达为如下公式(1)的线性化微分代数方程：

消去运行参数得到如下公式(2)：

得到状态矩阵的表达式如下公式(3)：

其中，x为状态变量，y为运行参数变量，μ为输入变量。

Δx为状态变量；

E

Δμ为输入变量。

如果状态矩阵A的所有特征值实部都是负的，那么根据李雅普诺夫定理， 该系统在受到小干扰时是稳定的。理论上当特征值出现非负实数时系统才会出 现低频振荡，但是在电力系统的实际运行过程中，需要预留一定的安全裕度。

定义一个称为谱横坐标的参数η来描述系统的安全裕度，表达式为如下公式 (4)：

Re(λ

其中，λ(A)表示状态矩阵A的所有特征值，Re(λ)表示特征值λ的实部，λ

把谱横坐标参数作为发电重新调度的稳定指标，在简化模型同时还能够降 低调度执行的时间。

此外，谱横坐标极限值应用于优化调度中，可以保证优化结果中的特征值 整体向虚轴左半部分偏移。

在本实施例中，小干扰稳定分析采用双轴同步发电机模型。

其矢量化的微分方程表示如下，依次为公式(5)～公式(8)：

其中，符号

δ是矢量化转子角度；

ω是矢量化转子速度；

ω

E′

X

X′

I

T′

T

M表示矢量化惯性时间常数；

D表示矢量化发电机阻尼系数；

E

定子代数方程的极坐标形式为：

表达式依次为公式(9)和公式(10)。

其中，V

R

IEEE-DC I型励磁器模型的矢量化微分方程为：

方程依次为公式(11)～公式(13)。

其中，T

K

V

T

K

V

R

S

网络模型的表达式为：

表达式依次为公式(14)和公式(15)。

其中，S

diag(V)表示把列向量V转化为对角稀疏矩阵；

Y表示为系统节点导纳矩阵；

非发电机组节点的矢量化极坐标潮流方程如下：

表达式依次为公式(16)和公式(17)。

其中，S

对电力系统元件线性化，小干扰稳定分析中的重要一步是线性化其微分代 数方程，对上述方程进行线性化得：

o＝C

O＝C

0＝D

表达式依次为公式(18)～公式(21)。其中：

A

B

B

C

C

D

D

D

D

D

D

D

ΔI

ΔV

Δθ

其中，公式(18)的详细展开如下，得到公式(22)：

其中，场饱和函数的表达式为：

其中，A

E

K

T

公式(18)～公式(21)经过矩阵变换后得到如下公式(23)：

则公式(23)和公式(3)对应的关系如下：

上述表达式依次为公式(24)～公式(27)。上式

上述表达式依次为公式(28)～公式(34)。

其中，E表示单位对角矩阵；

V

能够看出D

S3：分别计算特征值对系统变量的导数和发电机出力对系统变量的导数， 利用复合函数求导法则得到特征值对发电机出力的导数。

S31：特征值对状态变量的求导表示为：

上述表达式为公式(35)。

其中，λ

为了简化灵敏度的计算，需要对左右特征向量进行归一化处理，使得

为了提高计算效率和简化算法的开发难度，须对上式进行转换求解。

通过消掉逆矩阵的偏导形式并把三维矩阵运算转化为二维矩阵运算后得 到：

上述表达式为公式(36)。其中，

对于本实施例采用的发电机、励磁机模型，其系统变量为：

上式为公式(37)。

特征值对每个系统变量求导过程如下所示：

计算

有步骤S2的推导过程得知，E′

在根据公式(22)得到如下公式(40)：

其中，

根据步骤S2的推导过程可知，I

上式中τ

上述公式依次为公式(43)和公式(44)。其中，τ

V和θ仅存在于矩阵C

上述公式依次为公式(45)～公式(47)。其中，

对公式(28)求导得如下公式(48)：

上式中sparse(i，j，s，n，m)表示向量s的稀疏表达，即构造n行m列的稀疏 矩阵，其中横纵坐标i、j对应的数据为s。τ

对公式(29)求导有：

上述公式依次为公式(49)～公式(51)。其中，C

上述公式分别为公式(52)～公式(54)。

对公式(31)求导有如下公式(55)：

上式中sparse(i，j，s，n，m)T表示对稀疏矩阵求转置。 对公式(32)求导有：

上式分别为公式(56)～公式(58)。

对公式(33)求导得：

上述公式依次为公式(59)～公式(61)。

对D

其中，τ

推导潮流方程的矢量化海森矩阵，有如下公式(63)：

上式中的

S32：计算发电机出力对系统变量的导数

电力系统稳态下的代数方程：主要包括潮流方程、网络代数方程以及稳态 下的微分方程式。

完全极坐标形式潮流方程如下公式(64)所示：

网络方程极坐标形式如下公式(65)所示：

定子代数方程的极坐标形式如下公式(66)所示：

取自与系统变量求导相关的稳态下的微分方程组如下公式(67)所示：

G

G

G

P

P

对潮流方程式求导得如下公式(68)所示：

因此有如下公式(69)所示：

同理能够求得

公式(65)～公式(67)对发电机出力P

上式中的

至此，潮流方程、网络代数方程以及稳态下的微分方程对系统变量的求导 全部给出，然后求解公式(70)，即可得到发电机有功出力对系统变量导数

S33：利用复合函数求导法则得到特征值对发电机出力的导数。

根据步骤S31和步骤S32的推导过程，分别得到特征值对系统变量的导数 和特征值对发电机出力的导数，然后利用复合函数求导法则即可获得特征值对 发电机出力的解析灵敏度如下公式(72)所示：

S4：结合电力系统的实际运行规律获取特征值对发电机出力解析灵敏度。

对于步骤S33，利用复合函数求导法则得到特征值对发电机出力的导数后需 要结合电力系统的实际运行规律作相应处理，才能获得特征值对发电机出力灵 敏度的解析表达式。其中需要进行的处理包括：

S41：在步骤S32求导过程中对含自动发电控制AGC的发电机组进行修正， 即对连接含AGC的母线k的偏导数设置为零，具体为

S42：在步骤S32求导过程中对安装有自动化电压调节器AVR的发电机组进 行修正，即含AVR的母线j的偏导数也需要设置为零，具体为

通过上述所有步骤的分析计算，最终可获得特征值对发电机有功出力的解 析灵敏度。

S5：建立含1-范数目标函数的发电重新调度模型进行有效性和高效性的验 证。

目前计算特征值对发电机有功出力的灵敏度的最常用方法为基于摄动理论 的数值微分法。数值型特征值灵敏度通常只考虑谱横坐标灵敏度，其表示如下 公式(73)所示：

η(A)为基态下的谱横坐标，ε为发电机的摄动量，η(A

数值特征值灵敏度的最大优点是算法简单，不过也存在许多缺点：例如， 摄动量ε是人为设置的，所获得的计算结果具有很大的不确定性。此外，平衡节 点的灵敏度和PV节点的无功灵敏度都无法计算。

现分别选取WSSC 3机9节点系统、新英格兰10机39节点系统、IEEE 54 机118节点系统进行仿真计算。

利用公式(72)计算特征值对发电机有功出力的解析灵敏度，利用公式(73) 计算得到特征值数值灵敏度。获得的谱横坐标灵敏度分别如下：

表1 WSCC 3机9节点系统谱横坐标灵敏度

表2新英格兰10机39节点系统谱横坐标灵敏度

表3 IEEE 54机118节点系统谱横坐标灵敏度

数值谱横坐标灵敏度与解析谱横坐标灵敏度存在的共同点便是绝对值最大 灵敏度都位于同一台发电机，并且数值十分相近。其它节点的灵敏度数值有所 差距，这种差距的产生主要在于数值灵敏度的计算中多考虑了平衡节点的灵敏 度与PV节点的无功灵敏度，而解析灵敏度只考虑该节点上的有功灵敏度。它们 两者之间的关系可以表示如下：

谱横坐标可以表达成两个独立变量

谱横坐标为发电机的有功输出

上式为公式(74)。公式(74)按泰勒级数展开并且忽略高阶项后有如下公 式(75)：

其中S

在实际电力系统中，如果仅仅改变第i台发电机有功出力，即使摄动量很 小，平衡节点的有功出力、无功出力以及所有PV节点的无功出力也会随着变化。 因此上式可以进一步表达成如下公式(76)：

其中，ΔP

S

两边同时除以ΔP

式中，

通过观察公式(77)发现，即使在仅保留泰勒级数一阶展开式的情况下，谱横坐标的 数值灵敏度与解析灵敏度仍然存在一些差异，即数值灵敏度比解析灵敏度多考虑了后三项。

为验证本实施例所提方法的有效性和高效性，利用现有的技术，即数值灵 敏度进行对比分析。对于步骤S5，建立含1-范数目标函数的发电重新调度模型 如下：

约束条件为：

为了维持1-范数目标函数的优势，又不使其呈现非光滑性，将发电机组出 力的调整量ΔP

约束条件为：

该模型为简单的线性规划模型，利用Matlab中的函数linprog即可求解。 该模型中谱横坐标灵敏度与发电机有功出力的关系假设为线性关系。在求解过 程中由于调节量较大，第一次优化计算得到的谱横坐标η

根据上述公式对灵敏度进行更新：

复上述模型的优化计算，如此反复，直至满足裕度要求的精度：

其中，ξ为精度。

所有仿真计算都是在处理器内存为8G、主频为3.2GHz的一体机上完成。

特征值数值灵敏度以及基于矢量化技术的特征值解析灵敏度的编程工作皆 利用Matlab软件完成。仿真结果如下：

表4两种不同方法下WSCC 3机9节点系统机组有功出力的调整量

3机9节点系统仿真系统，基态下稳定裕度η

表5两种不同方法下新英格兰10机39节点系统机组有功出力的调整量

10机39节点系统仿真系统，基态下稳定裕度η

表6两种不同方法下IEEE 54机118节点系统机组有功出力的调整量

54机118节点系统仿真系统，基态下稳定裕度η

从上述仿真结果可以看出解析灵敏度比数值灵敏度的效果更优，即在达到 同样的稳定裕度要求下，10机39节点仿真结果表明解析灵敏度需要调整的有功 出力值更小，54机118节点系统仿真结果显示数值灵敏度无法收敛而解析灵敏 度仍然能够收敛。此外，两种方法的仿真时间对比分析如下：

表7基于1-范数有功调度方法的计算时间

从表格7我们可以看出，在计算特征值灵敏度的效率上，解析灵敏度的计 算效率更高，尤其随着网络节点的增加，数值特征值灵敏度的计算效率急剧下 降，当网络节点为118时，其计算时间为解析法的12倍多。

从上面的范例仿真结果可以看出，本发明提出的一种基于矢量化技术的小 干扰稳定的特征值灵敏度计算方法在提高系统的稳定裕度上比现有技术更为有 效和高效。