基于广义特征值的时滞稳定上限计算系统及其计算方法.pdf

上传人:罗明 文档编号:6148271 上传时间:2019-04-19 格式:PDF 页数:22 大小:1.95MB
返回 下载 相关 举报
摘要
申请专利号:

CN201410067208.5

申请日:

2014.02.26

公开号:

CN103838965A

公开日:

2014.06.04

当前法律状态:

授权

有效性:

有权

法律详情:

授权|||实质审查的生效IPC(主分类):G06F 19/00申请日:20140226|||公开

IPC分类号:

G06F19/00(2011.01)I

主分类号:

G06F19/00

申请人:

华北电力大学

发明人:

马静; 李俊臣; 高翔; 丁秀香; 王增平

地址:

102206 北京市昌平区回龙观朱辛庄2号

优先权:

专利代理机构:

北京麟保德和知识产权代理事务所(普通合伙) 11428

代理人:

周恺丰

PDF下载: PDF下载
内容摘要

本发明公开了电力系统稳定分析技术领域中的一种基于广义特征值的时滞稳定上限计算系统及其计算方法。系统包括顺序相连的数据采集模块、时滞系统处理模块、时滞上限求解模块和结果输出模块;方法包括:采集建立时滞系统状态方程所需的网络结构参数、发电机频率和发电机功角;建立时滞系统状态方程,并对时滞系统状态方程中的参数矩阵进行降阶处理,得到降阶后的时滞系统状态方程;生成基于改进自由权矩阵的时滞稳定判据;利用时滞稳定判据求解时滞稳定上限。本发明能够有效降低时滞稳定上限计算过程中的保守性,具有很好的正确性和有效性。

权利要求书

权利要求书
1.  一种基于广义特征值的时滞稳定上限计算系统,其特征是所述系统包括顺序相连的数据采集模块、时滞系统处理模块、时滞上限求解模块和结果输出模块;
所述数据采集模块用于采集建立时滞系统状态方程所需的网络结构参数、发电机频率和发电机功角,并将采集的数据发送至时滞系统处理模块;
所述时滞系统处理模块用于建立时滞系统状态方程,并对时滞系统状态方程中的参数矩阵进行降阶处理;
所述时滞上限求解模块用于生成基于改进自由权矩阵的时滞稳定判据,并利用时滞稳定判据求解时滞稳定上限;
所述结果输出模块用于输出时滞稳定上限结果。

2.  一种基于广义特征值的时滞稳定上限计算方法,其特征是所述方法包括:
步骤1:采集建立时滞系统状态方程所需的网络结构参数、发电机频率和发电机功角;
步骤2:建立时滞系统状态方程,并对时滞系统状态方程中的参数矩阵进行降阶处理,得到降阶后的时滞系统状态方程;
步骤3:生成基于改进自由权矩阵的时滞稳定判据;
步骤4:利用时滞稳定判据求解时滞稳定上限。

3.  根据权利要求2所述的计算方法,其特征是所述降阶后的时滞系统状态方程为
其中,x(t)为降阶后的电力系统状态向量;
Α为降阶后的电力系统状态矩阵;
Αd为降阶后的电力系统时滞矩阵;
为降阶后的电力系统状态量对应的状态值;
d(t)为时滞,0≤d(t)≤h且
h为时滞稳定上限且h>0;
μ为时滞最大变化率。

4.  根据权利要求3所述的计算方法,其特征是所述基于改进自由权矩阵的时滞稳定判据为:
ΦhNhShMhAcT(Z1+Z2)hNT-hZ1000hST0-hZ100hMT00-hZ20h(Z1T+Z2T)Ac000-h(Z1+Z2)<0;]]>
其中,Φ=Φ1+Φ2+Φ2T;]]>
Φ1=PA+ATP+Q+RPAd0AdTP-(1-μ)Q000-R;]]>
Φ2=[N+M-N+S-M-S];
Ac=[AAd0];
N、M和S为改进自由权矩阵;
改进自由权矩阵N满足2ζ1T(t)N[x(t)-x(t-d(t))-&Integral;t-d(t)tx&CenterDot;(s)ds]=0;]]>
改进自由权矩阵S满足2ζ1T(t)S[x(t-d(t))-x(t-h)-&Integral;t-ht-d(t)x&CenterDot;(s)ds]=0;]]>
改进自由权矩阵M满足2ζ1T(t)M[x(t)-x(t-h)-&Integral;t-htx&CenterDot;(s)ds]=0;]]>
ζ1(t)=[xT(t)xT(t-d(t))xT(t-h)]T;
P、Q、R、Z1和Z2为待定矩阵。

5.  根据权利要求4所述的计算方法,其特征是所述利用时滞稳定判据求解时滞稳定上限包括:
子步骤A1:将时滞稳定判据等价变换为
ΦNSMAcT(Z1+Z2)NT-Y1000ST0-Y100MT00-Y20(Z1T+Z2T)Ac000-(Y1+Y2)<0;]]>
其中,Y1和Y2为附加矩阵,并且Y1=Y1T≥0,Y2=Y2T≥0,
Y100Y2<vZ100Z2,v=1/h;]]>
子步骤A2:以v最小为目标,以时滞稳定判据等价变换后的矩阵不等式ΦNSMAcT(Z1+Z2)NT-Y1000ST0-Y100MT00-Y20(Z1T+Z2T)Ac000-(Y1+Y2)<0]]>以及Y100Y2<Z100Z2]]>为约束条件,计算待定矩阵P、Q、R、Z1、Z2、Y1和Y2,进而得到时滞稳定上限h。

说明书

说明书基于广义特征值的时滞稳定上限计算系统及其计算方法
技术领域
本发明属于电力系统稳定分析技术领域,尤其涉及一种基于广义特征值的
时滞稳定上限计算系统及其计算方法。
背景技术
基于广域信息的反馈控制器的设计存在时滞问题,而这必然会造成控制器的控制效果降低,甚至出现负阻尼的情况。因此,迫切需要对系统的时滞稳定上限进行研究。
目前国内外学者在系统时滞稳定上限的研究方面取得了大量有益成果,主要可以分为3类:时域法、频域法和直接法。时域法可判定系统在特定场景下是否稳定,但在稳定程度、时滞稳定上限等信息的获取方面还需要进一步研究。频域法通过在实数空间中搜索时滞系统的关键特征根,能够在一定程度上揭示时滞系统变化规律,但计算量较大,求解速度有待提升。直接法借助Lyapunov理论和线性矩阵不等式(Linear Matrix Inequality,LMI)技术,可同时考虑时滞的随机变动、存在切换环节等情况,适用范围更为广泛,但该方法具有一定的保守性。
针对以上问题,本发明提出一种基于广义特征值的时滞稳定上限计算系统及其计算方法。本发明从电力系统实际出发,能够有效解决直接法求解时滞上限过程中保守性较高的问题。首先利用读入的广域信号数据建立时滞系统状态方程,并在有效保留系统中低频振荡成分的基础上降阶系统。其次,该发明形成基于改进自由权矩阵的时滞稳定判据,在此基础上,将时滞稳定判据等价变 换,利用广义特征值法求解系统的时滞稳定上限。基于IEEE4机11节点系统和IEEE16机68节点系统仿真表明,本发明可以较好的降低传统方法的保守性,具有很好的有效性和正确性。
发明内容
本发明的目的在于,提供一种基于广义特征值的时滞稳定上限计算系统及其计算方法,用于解决直接法求解时滞上限过程中保守性较高的问题。
为了实现上述目的,本发明提出的技术方案是,一种基于广义特征值的时滞稳定上限计算系统,其特征是所述系统包括顺序相连的数据采集模块、时滞系统处理模块、时滞上限求解模块和结果输出模块;
所述数据采集模块用于采集建立时滞系统状态方程所需的网络结构参数、发电机频率和发电机功角,并将采集的数据发送至时滞系统处理模块;
所述时滞系统处理模块用于建立时滞系统状态方程,并对时滞系统状态方程中的参数矩阵进行降阶处理;
所述时滞上限求解模块用于生成基于改进自由权矩阵的时滞稳定判据,并利用时滞稳定判据求解时滞稳定上限;
所述结果输出模块用于输出时滞稳定上限结果。
一种基于广义特征值的时滞稳定上限计算方法,其特征是所述方法包括:
步骤1:采集建立时滞系统状态方程所需的网络结构参数、发电机频率和发电机功角;
步骤2:建立时滞系统状态方程,并对时滞系统状态方程中的参数矩阵进行降阶处理,得到降阶后的时滞系统状态方程;
步骤3:生成基于改进自由权矩阵的时滞稳定判据;
步骤4:利用时滞稳定判据求解时滞稳定上限。
所述降阶后的时滞系统状态方程为其中,x(t)为降阶后的电力系统状态向量;
Α为降阶后的电力系统状态矩阵;
Αd为降阶后的电力系统时滞矩阵;
为降阶后的电力系统状态量对应的状态值;
d(t)为时滞,0≤d(t)≤h且
h为时滞稳定上限且h>0;
μ为时滞最大变化率。
所述基于改进自由权矩阵的时滞稳定判据为:
ΦhNhShMhAcT(Z1+Z2)hNT-hZ1000hST0-hZ100hMT00-hZ20h(Z1T+Z2T)Ac000-h(Z1+Z2)<0;]]>
其中,Φ=Φ1Φ2Φ2T;]]>
Φ1=PA+ATP+Q+RPAd0AdTP-(1-μ)Q000-R;]]>
Φ2=[N+M-N+S-M-S];
Ac=[AAd0];
N、M和S为改进自由权矩阵;
改进自由权矩阵N满足2ζ1T(t)N[x(t)-x(t-d(t))-&Integral;t-d(t)tx&CenterDot;(s)ds]=0;]]>
改进自由权矩阵S满足2ζ1T(t)S[x(t-d(t))-x(t-h)-&Integral;t-ht-d(t)x&CenterDot;(s)ds]=0;]]>
改进自由权矩阵M满足2ζ1T(t)M[x(t)-x(t-h)-&Integral;t-htx&CenterDot;(s)ds]=0;]]>
?1(t)=[xT(t)xT(t-d(t))xT(t-h)]T;
P、Q、R、Z1和Z2为待定矩阵。
所述利用时滞稳定判据求解时滞稳定上限包括:
子步骤A1:将时滞稳定判据等价变换为
ΦNSMAcT(Z1+Z2)NT-Y1000ST0-Y100MT00-Y20(Z1T+Z2T)Ac000-(Y1+Y2)<0;]]>
其中,Y1和Y2为附加矩阵,并且Y1=Y1T≥0,Y2=Y2T≥0,
Y100Y2<vZ100Z2,v=1/h;]]>
子步骤A2:以v最小为目标,以时滞稳定判据等价变换后的矩阵不等式ΦNSMAcT(Z1+Z2)NT-Y1000ST0-Y100MT00-Y20(Z1T+Z2T)Ac000-(Y1+Y2)<0]]>以及Y100Y2<Z100Z2]]>为约束条件,计算待定矩阵P、Q、R、Z1、Z2、Y1和Y2,进而得到时滞稳定上限h。
本发明利用改进自由权矩阵建立时滞系统稳定判据,借助广义特征值法对系统时滞上限进行求解,其能够有效降低时滞稳定上限计算过程中的保守性,具有很好的正确性和有效性。
附图说明
图1是基于广义特征值的时滞稳定上限计算系统结构图;
图2是IEEE4机11节点系统结构图;
图3是IEEE4机11节点SMA降阶前后特征根对比图;其中,(a)是全阶开环状态矩阵和降阶开环状态矩阵的特征根对比图,(b)是全阶闭环状态矩阵和降阶闭环状态矩阵的特征根对比图;
图4是IEEE4机11节点不同延时情况下发电机1-4相对功角动态响应曲线图;
图5是IEEE4机11节点不同延时情况下发电机2-3相对功角动态响应曲线图;
图6是IEEE4机11节点系统G1与G4功角差各时滞时间下的阻尼比结果表;
图7是IEEE4机11节点系统G2与G3功角差各时滞时间下的阻尼比结果表;
图8是IEEE16机68节点系统结构图;
图9是IEEE16机68节点Schur降阶前后频率响应对比图;其中,(a)是全阶开环状态矩阵和降阶开环状态矩阵的特征根对比图,(b)是全阶闭环状态矩阵和降阶闭环状态矩阵的特征根对比图;
图10是IEEE16机68节点不同延时情况下发电机1-16相对功角动态响应曲线图;
图11是IEEE16机68节点不同延时情况下发电机3-14相对功角动态响应曲线图;
图12是IEEE16机68节点不同延时情况下发电机10-15相对功角动态响应曲线图;
图13是16机系统G1与G16功角差各时滞时间下的阻尼比表;
图14是16机系统G3与G14功角差各时滞时间下的阻尼比表;
图15是16机系统G10与G15功角差各时滞时间下的阻尼比表。
具体实施方式
下面结合附图,对优选实施例作详细说明。应该强调的是,下述说明仅仅是示例性的,而不是为了限制本发明的范围及其应用。
实施例1
图1是本发明提供的基于广义特征值的时滞稳定上限计算系统结构图。如图1所示,本发明提供的基于广义特征值的时滞稳定上限计算系统结构图包括顺序相连的数据采集模块、时滞系统处理模块、时滞上限求解模块和结果输出模块。
数据采集模块用于采集网络结构参数、发电机频率和发电机功角,并将采集的数据发送至时滞系统处理模块。其中,网络结构参数是用于建立时滞系统状态方程所需的参数,而发电机频率和发电机功角则是时滞系统状态方程中的状态向量。
时滞系统处理模块用于建立时滞系统状态方程,并对时滞系统状态方程中的参数矩阵进行降阶处理,得到降阶后的时滞系统状态方程。
时滞上限求解模块用于生成基于改进自由权矩阵的时滞稳定判据;并对所述时滞稳定判据进行等价变换,利用广义特征值法求解系统的时滞稳定上限。
结果输出模块用于输出时滞稳定上限结果。
本发明提供的基于广义特征值的时滞稳定上限计算方法包括:
步骤1:采集建立时滞系统状态方程所需的网络结构参数、发电机频率和发电机功角。
网络结构参数包括时滞系统中线路的阻抗、导纳、发电机的内阻抗和负荷的等值阻抗。发电机频率用于计算时滞系统转速,发电机功角和时滞系统转速的变化量为时滞系统状态方程的状态向量。
步骤2:建立时滞系统状态方程,并对时滞系统状态方程中的参数矩阵进行降阶处理,得到降阶后的时滞系统状态方程。
多输入多输出电力系统的状态方程可表示为:
x&CenterDot;(t)=Ax(t)+Bu(t)u(t)=K1x(t)---(1)]]>
其中,x′(t)∈Rn为电力系统状态向量,u′(t)∈Rm为电力系统控制输入向量,A′∈Rn×n为电力系统状态矩阵,B′∈Rn×m为电力系统控制矩阵。
通过状态反馈后得到相应的闭环系统为:
x&CenterDot;(t)=Cx(t)---(2)]]>
其中,C′为闭环状态矩阵,当系统经过状态反馈时,闭环状态矩阵为C′=A′+B′K′1,其中,K′1∈Rm×n为各附加控制器的综合状态反馈矩阵。
在实际电力系统中,控制输入向量通过SCADA/WAMS系统向各控制器传达,信号传递过程必然存在一定时滞,则相应闭环系统可描述为:
x&CenterDot;(t)=Ax(t)+BK1x(t-d(t))---(3)]]>
由式(3)可知,电力系统时滞矩阵A′d=B′K′1。对于含时滞环节的系统,其状态方程有如下形式:

其中,矩阵A′和A′d分别为电力系统状态矩阵和电力系统时滞矩阵,h为时滞稳定上限。公式(4)中时滞d(t)满足条件:
0≤d(t)≤h                          (5)
d&CenterDot;(t)μ---(6)]]>
公式(4)-式(6)即为时滞系统状态方程,μ为时滞最大变化率。
实际系统分析中,仅对低频振荡特征根相应的特征向量中与系统状态量Δω(发电机转速变化量)和Δδ(发电机功角变化量)相对应的元素感兴趣,以便了解在Δωi中所含的该振荡模式分量的相对幅值及相位。常用的系统降阶方法包括SMA降阶方法和Schur平衡降阶方法。由于系统降阶方法已经是本领域常用的方法,因此本发明只以选择模式分析法SMA为例,对系统降阶做简单介绍。
对于系统状态方程将其按下式划分
X&CenterDot;1X&CenterDot;2=A11A12A21A22X1X2---(7)]]>
其中,X1=[ΔωT,ΔδT]为保留变量,X2为其他变量,待消去。
由式(7)可消去X2,得:
X&CenterDot;1=[A11+A12(pI-A22)-1A21]X1---(8)]]>
其中,I为单位阵,p为微分算子。
将上式改写为
X&CenterDot;r=Ar(p)Xr---(9)]]>
其中,Xr=X1为保留变量,Ar(p)为运算形式的“降阶”系统系数阵。
由式(7)-(9)可以得到两个重要性质:
(1)如果p=λ1(i=1,2,…,N)为式(7)相应系统特征根,即|λiI-A|=0,则p=λi也为式(8)或(9)形式上降阶的系统特征根,即亦有|λiI-Ar(λi)|=0,特征根不发生变化,系统模式不变。
(2)对于原系统,λi的特征向量ui,有Aui=λiui。设降阶系统λi相应的特征 向量为uri,即Ar(λi)uri=λiuri,则uri和ui中保留变量Xr相对应的元素相等,即特征向量的相应元素不变。因此,在Xr保留变量处去观察同一模式λi的振荡时,相对幅值即相位不变,或者说模态不变。这样,所关心频带输入输出特性被完整的保留下来。
通过降阶处理,可以得到降阶后的时滞系统状态方程为:

公式(10)中,x(t)是降阶后的电力系统状态向量,Α为降阶后的电力系统状态矩阵,Αd为降阶后的电力系统时滞矩阵,为降阶后的电力系统状态量对应的状态值,d(t)、h和μ的含义同公式(4)且满足式(5)和式(6)。
步骤3:生成基于改进自由权矩阵的时滞稳定判据。
构造如下形式的Lyapunov-Krasovskii泛函(李雅普诺夫-克拉索夫斯基泛函):
V(x)=xT(t)Px(t)+&Integral;t-d(t)txT(s)Qx(s)ds+&Integral;t-htxT(s)Rx(s)ds+&Integral;-h0&Integral;t+θtx&CenterDot;T(s)(Z1+Z2)x&CenterDot;(x)dsdθ---(11)]]>
公式(11)中,P=PT>0,Q=QT≥0,R=RT≥0和Zi=ZiT≥0(i=1,2)都是待定矩阵。
由Newton-Leibniz公式(牛顿-莱布尼兹公式)可知,对于改进自由权矩阵N、S和M,式(12)-式(14)成立:
2ζ1T(t)N[x(t)-x(t-d(t))-&Integral;t-d(t)tx&CenterDot;(s)ds]=0---(12)]]>
2ζ1T(t)S[x(t-d(t))-x(t-h)-&Integral;t-ht-d(t)x&CenterDot;(s)ds]=0---(13)]]>
2ζ1T(t)M[x(t)-x(t-h)-&Integral;t-htx&CenterDot;(s)ds]=0---(14)]]>
其中,ζ1(t)=[xT(t)xT(t-d(t))xT(t-h)]T。
计算式(11)中V(x)关于t的导数为:
V&CenterDot;(x)=2xT(t)Px&CenterDot;(t)+xT(s)Qx(s)-(1-d&CenterDot;(t))xT(t-d(t))Qx(t-d(t))+xT(t)Rx(t)-xT(t-h)Rx(t-h)+hx&CenterDot;T(t)(Z1+Z2)x&CenterDot;(t)-&Integral;t-htx&CenterDot;T(s)(Z1+Z2)x&CenterDot;(s)ds---(15)]]>
将式(10)中第一个式子及式(12)-式(14)代入式(15),并加入必要的松散项可得:
V&CenterDot;(x)xT(t)(PA+ATP)x(t)+xT(s)(Q+R)x(s)-(1-μ)xT(t-d(t))Qx(t-d(t))-xT(t-h)Rx(t-h)+hx&CenterDot;T(t)(Z1+Z2)x&CenterDot;(t)-&Integral;t-htx&CenterDot;T(s)Z1x&CenterDot;T(s)ds-&Integral;t-ht-d(t)x&CenterDot;T(s)Z1x&CenterDot;T(s)ds-&Integral;t-htx&CenterDot;T(s)Z2x&CenterDot;T(s)ds+2ζ1T(t)N[x(t)-x(t-d(t))-&Integral;t-d(t)tx&CenterDot;(s)ds]+2ζ1T(t)S[x(t-d(t))-x(t-h)-&Integral;t-d(t)tx&CenterDot;(s)ds]+2ζ1T(t)M[x(t)-x(t-h)-&Integral;t-htx&CenterDot;(s)ds]ζ1T(t)[Φ+Φs]ζ1(t)-&Integral;t-d(t)tθ1ds-&Integral;t-ht-d(t)θ2ds-&Integral;t-htθ3ds---(16)]]>
Φs=hAcT(Z1+Z2)Ac+hNZ1-1NT+hSZ1-1ST+hMZ1-1MT---(17)]]>
θ1=[ζ1T(t)N+x&CenterDot;(t)Z1]Z1-1[NTζ1(t)+Z1x&CenterDot;(t)]---(18)]]>
θ2=[ζ1T(t)S+x&CenterDot;(t)Z1]Z1-1[STζ1(t)+Z1x&CenterDot;(t)]---(19)]]>
θ3=[ζ1T(t)M+x&CenterDot;(t)Z2]Z2-1[MTζ1(t)+Z1x&CenterDot;(t)]---(20)]]>
上述公式中,Φ=Φ1+Φ2+Φ2T.]]>
Φ1=PA+ATP+Q+RPAd0AdTP-(1-μ)Q000-R,]]>Φ2=[N+M-N+S-M-S]。
Ac=[AAd0]。
NT=N1TN2TN3T,]]>ST=S1TS2TS3T,]]>MT=M1TM2TM3T.]]>
N1、N2和N3为矩阵N的分块矩阵且N1的维度与PA+ATP+Q+R的维度相等,N2的维度与的维度相等,N3的维度与R的维度相等。
S1、S2和S3为矩阵S的分块矩阵且S1的维度与PA+ATP+Q+R的维度相等,S2的维度与的维度相等,S3的维度与R的维度相等。
M1、M2和M3为矩阵M的分块矩阵且M1的维度与PA+ATP+Q+R的维度相等,M2的维度与的维度相等,M3的维度与R的维度相等。
考虑到式(11)中Zi=ZiT≥0,(i=1,2),因此公式(18)-式(20)中θi=θiT≥0,(i=1,2,3)。若Φ+Φs≤0,则式(16)中对于Φ+Φs,利用Schur补后可得式(10)表征的时滞系统稳定判据如下:
对于给定标量h>0和μ,若存在P=PT>0,Q=QT≥0,R=RT≥0和Zi=ZiT≥0(i=1,2),NT=N1TN2TN3T,]]>ST=S1TS2TS3T]]>MT=M1TM2TM3T,]]>使得如下线性矩阵不等式成立:
ΦhNhShMhAcT(Z1+Z2)hNT-hZ1000hST0-hZ100hMT00-hZ20h(Z1T+Z2T)Ac000-h(Z1+Z2)<0---(21)]]>
则对于同时满足条件式(5)和式(6)的时滞系统(10)渐进稳定。
步骤4:利用时滞稳定判据求解时滞稳定上限。
式(21)表征的线性矩阵不等式仅能判定系统是否稳定,而无法获取系统时滞稳定上限等信息。考虑到广义特征值法能够求解优化问题的全局最小值,因此,本发明提出利用广义特征值法计算系统的时滞稳定上限。由于式(21)不是标准的广义特征值形式,无法直接利用广义特征值法进行求解,因此,做如下变换,将式(21)同时左乘与右乘式(22),如下:
I00000I/h00000I/h00000I/h00000I/h>0---(22)]]>
可得:
ΦNSMAcT(Z1+Z2)NT-vZ1000ST0-vZ100MT00-vZ20(Z1T+Z2T)Ac000-v(Z1+Z2)<0---(23)]]>
其中,v=1h,I为单位阵。
为了求解最大时滞稳定上限h,即最小的v,本发明引入附加矩阵Yi=YiT≥0(i=1,2),该矩阵需满足式(24)成立:
Y100Y2<vZ100Z2---(24)]]>
再将式(24)代入式(23)可得:
ΦNSMAcT(Z1+Z2)NT-Y1000ST0-Y100MT00-Y20(Z1T+Z2T)Ac000-(Y1+Y2)<0---(25)]]>
由此,时滞稳定上限问题已经转化为以下优化问题:
以v最小为目标,以ΦNSMAcT(Z1+Z2)NT-Y1000ST0-Y100MT00-Y20(Z1T+Z2T)Ac000-(Y1+Y2)<0]]>Y100Y2<vZ100Z2]]>为约束条件,计算待定矩阵P、Q、R、Z1、Z2、Y1和Y2,进而得到时滞稳定上限h。公式表示如下:
min v
P,Q,R,Zi,Yi
s.t.
Y100Y2<vZ100Z2ΦNSMAcT(Z1+Z2)NT-Y1000ST0-Y100MT00-Y20(Z1T+Z2T)Ac000-(Y1+Y2)<0---(26)]]>
通过求解式(26),可以计算出以矩阵P、Q、R、Z1、Z2、Y1和Y2为变量,以式(25)和式(26)为约束的最小v。最终,利用h=1/v可以求出时滞稳定上限。
实施例2
基于MATLAB仿真软件搭建的图2所示的IEEE4机11节点系统,发电机 采用6阶详细模型,励磁系统采用快速励磁,基准模型下的负荷采用50%恒阻抗和50%恒电流模型。首先,通过模态分析法得到四机系统的状态矩阵,并利用SMA方法分别对开、闭环状态矩阵进行降阶,如图3所示。
图3(a)是全阶开环状态矩阵和降阶开环状态矩阵的特征根对比图,(b)是全阶闭环状态矩阵和降阶闭环状态矩阵的特征根对比图。由图3可以看出,开、闭环状态矩阵在利用SMA进行降阶后,均保留了系统中低频振荡成分。因此,利用降阶后系统矩阵能够求解系统允许的最大时滞。将降阶后状态矩阵Α和时滞矩阵Αd代入式(26),求得最大时滞边界h=281.88ms。其中,降阶后状态矩阵Α如下:
A=000376.9000000376.9000000376.9-0.0730.0650.004-0.7300.2720.0760.058-0.0870.0091.160-0.343-0.1340.0080.011-0.085-0.0200.047-0.554.]]>
降阶后时滞矩阵Αd如下:
Ad=000000000000000000000-0.234-0.8390.0100-0.00110.0010-0.348-1.362-0.13800.001000.049-0.290-0.638.]]>
时滞分别设置为h1=50ms,h2=100ms,h3=288.88ms和h4=350ms。所观察的物理量为发电机G1和G4之间的功角差,以及G2和G3之间的功角差,分别如图4和图5所示。
由图4和图5均可看出,在未加广域阻尼控制时,发电机之间的功角发生 了严重的低频振荡。加入广域阻尼控制后,在不考虑时滞的情况下,低频振荡可得到有效抑制;然而随着时滞的增加,阻尼效果随之减弱,在最大时滞288ms的情况下,虽然仍存在一定的阻尼,但其阻尼比已经降到10%以下,说明此时控制器已经不满足控制要求。
利用prony算法对各时滞下的功角差曲线进行阻尼比分析,结果如图6和图7所示。
实施例3
基于MATLAB仿真软件搭建的图8所示的IEEE16机68节点系统进一步考查基于广义特征值法的最大时滞求解方法的有效性和通用性。该系统重要联络线为区域4和区域5的联络线1-2,1-27和8-9。发电机采用6阶详细模型,励磁采用IEEE-DC1型励磁,负荷模型15%恒有功功率,25%的恒有功电流和15%的恒无功功率,25%的恒无功功率和60%恒阻抗。首先利用Schur平衡降阶方法,在保证所关心频带的输入输出特性保持不变的前提下,对系统进行降阶,如图9所示。
图9(a)中实线为开环全阶系统的频率响应,虚线为开环降阶系统的频率响应,可以看出,降阶系统和全阶系统的输入输出特性相同。图9(b)中实线为闭环全阶系统的频率响应,虚线为闭环降阶系统的频率响应,可以看出,闭环降阶系统同样保留了全阶系统的特性。因此,利用降阶后的系统矩阵求解全阶系统的时滞稳定上限具有有效性和可行性。
降阶后的状态矩阵Αre与降阶后的时滞矩阵Αdre代入式(26),可得最大的时滞上限为h=93ms。时滞分别设置为h1=50ms和h2=93ms。以发电机G1和G16之间的功角差,G3和G14之间的功角差与G10和G15之间的功角差,为所观察的物理量,分别如图10,图11和图12所示。其中,降阶后的状态矩阵Αre如下:
Are=-16.1234191.6215-34.7234-93.6034-46.9813-40.7245148.6905-9.40370.5458-0.1871-8.582273.851138.379628.2807-31.1118-3.53283.0915-0.81650.46700.38718.6901-56.6390-119.0312-154.230389.48598.39817.654547.8142-16.057220.4179-16.3933128.3466100.713891.8190-63.40741.0582-30.2333-5.40161.8643-2.0558-21.5589183.071841.8903-25.8071-3.006515.5680-58.595454.2491-17.982124.67538.9164-85.969941.1758105.8473-47.0510-8.151754.598135.2028-2.44248.9442-3.329627.08186.4952-4.40472.32595.3976-8.772833.8530-9.137913.11600.6362-4.3077-5.9066-6.77253.3248-1.2021-0.6273-11.95335.6377-0.85970.3450-3.42042.05374.4665-1.9128-0.94820.6125-6.28614.34356.76910.4673-4.19030.94373.9651-2.1654-0.74902.0926-2.3038-5.0688-4.4607.]]>
降阶后的时滞矩阵Αdre如下:
Adre=-0.01680.13870.04010.0041-0.01390.0056-0.0470-0.003500.00230.1610-1.3994-0.20930.29680.0753-0.01220.32610.1445-0.02640.00950.0414-0.2064-0.5163-0.61590.2242-0.08370.2789-0.27640.0783-0.0886-0.02470.3531-0.7123-1.18810.4669-0.0484-0.1044-0.58400.1524-0.167700.09470.42360.6927-0.5043-0.09540.60540.4081-0.10200.12510.1195-0.61200.42101.2198-0.9955-0.33451.96670.7972-0.18510.2245-0.1383-0.00100.5268-0.01491.26200.8615-4.0937-0.45670.0720-0.1352-0.50983.6470-2.0410-5.41943.24880.7903-5.5637-3.16910.7561-0.86150.5506-3.81021.91745.4325-3.4277-0.92816.24423.2459-0.76730.88331.0523-7.11844.580911.8740-7.5930-1.953913.12297.0687-1.68381.9525.]]>
随着时滞的增加,阻尼效果随之减弱,在最大时滞93ms的情况下,虽然存在一定的阻尼效果,但是其阻尼比已经降到10%以下,此时控制器同样也不能满足控制的要求。利用prony方法对各时滞下的功角差曲线进行阻尼比分析,结果如图13、图14和图15所示。
以上所述,仅为本发明较佳的具体实施方式,但本发明的保护范围并不局限于此,任何熟悉本技术领域的技术人员在本发明揭露的技术范围内,可轻易想到的变化或替换,都应涵盖在本发明的保护范围之内。因此,本发明的保护范围应该以权利要求的保护范围为准。

基于广义特征值的时滞稳定上限计算系统及其计算方法.pdf_第1页
第1页 / 共22页
基于广义特征值的时滞稳定上限计算系统及其计算方法.pdf_第2页
第2页 / 共22页
基于广义特征值的时滞稳定上限计算系统及其计算方法.pdf_第3页
第3页 / 共22页
点击查看更多>>
资源描述

《基于广义特征值的时滞稳定上限计算系统及其计算方法.pdf》由会员分享,可在线阅读,更多相关《基于广义特征值的时滞稳定上限计算系统及其计算方法.pdf(22页珍藏版)》请在专利查询网上搜索。

1、(10)申请公布号 CN 103838965 A (43)申请公布日 2014.06.04 CN 103838965 A (21)申请号 201410067208.5 (22)申请日 2014.02.26 G06F 19/00(2011.01) (71)申请人 华北电力大学 地址 102206 北京市昌平区回龙观朱辛庄 2 号 (72)发明人 马静 李俊臣 高翔 丁秀香 王增平 (74)专利代理机构 北京麟保德和知识产权代理 事务所 ( 普通合伙 ) 11428 代理人 周恺丰 (54) 发明名称 基于广义特征值的时滞稳定上限计算系统及 其计算方法 (57) 摘要 本发明公开了电力系统稳定分析。

2、技术领域中 的一种基于广义特征值的时滞稳定上限计算系统 及其计算方法。系统包括顺序相连的数据采集模 块、 时滞系统处理模块、 时滞上限求解模块和结果 输出模块 ; 方法包括 : 采集建立时滞系统状态方 程所需的网络结构参数、 发电机频率和发电机功 角 ; 建立时滞系统状态方程, 并对时滞系统状态 方程中的参数矩阵进行降阶处理, 得到降阶后的 时滞系统状态方程 ; 生成基于改进自由权矩阵的 时滞稳定判据 ; 利用时滞稳定判据求解时滞稳定 上限。本发明能够有效降低时滞稳定上限计算过 程中的保守性, 具有很好的正确性和有效性。 (51)Int.Cl. 权利要求书 2 页 说明书 11 页 附图 8 。

3、页 (19)中华人民共和国国家知识产权局 (12)发明专利申请 权利要求书2页 说明书11页 附图8页 (10)申请公布号 CN 103838965 A CN 103838965 A 1/2 页 2 1. 一种基于广义特征值的时滞稳定上限计算系统, 其特征是所述系统包括顺序相连的 数据采集模块、 时滞系统处理模块、 时滞上限求解模块和结果输出模块 ; 所述数据采集模块用于采集建立时滞系统状态方程所需的网络结构参数、 发电机频率 和发电机功角, 并将采集的数据发送至时滞系统处理模块 ; 所述时滞系统处理模块用于建立时滞系统状态方程, 并对时滞系统状态方程中的参数 矩阵进行降阶处理 ; 所述时滞上。

4、限求解模块用于生成基于改进自由权矩阵的时滞稳定判据, 并利用时滞稳 定判据求解时滞稳定上限 ; 所述结果输出模块用于输出时滞稳定上限结果。 2. 一种基于广义特征值的时滞稳定上限计算方法, 其特征是所述方法包括 : 步骤 1 : 采集建立时滞系统状态方程所需的网络结构参数、 发电机频率和发电机功角 ; 步骤 2 : 建立时滞系统状态方程, 并对时滞系统状态方程中的参数矩阵进行降阶处理, 得到降阶后的时滞系统状态方程 ; 步骤 3 : 生成基于改进自由权矩阵的时滞稳定判据 ; 步骤 4 : 利用时滞稳定判据求解时滞稳定上限。 3. 根据权利要求 2 所述的计算方法, 其特征是所述降阶后的时滞系统。

5、状态方程为 其中, x(t) 为降阶后的电力系统状态向量 ; 为降阶后的电力系统状态矩阵 ; d为降阶后的电力系统时滞矩阵 ; 为降阶后的电力系统状态量对应的状态值 ; d(t) 为时滞, 0 d(t) h 且 h 为时滞稳定上限且 h0 ; 为时滞最大变化率。 4. 根据权利要求 3 所述的计算方法, 其特征是所述基于改进自由权矩阵的时滞稳定判 据为 : 其中, 权 利 要 求 书 CN 103838965 A 2 2/2 页 3 2=N+M-N+S-M-S ; Ac=AAd0 ; N、 M 和 S 为改进自由权矩阵 ; 改进自由权矩阵 N 满足 改进自由权矩阵 S 满足 改进自由权矩阵 M。

6、 满足 1(t)=xT(t)xT(t-d(t)xT(t-h)T; P、 Q、 R、 Z1和 Z2为待定矩阵。 5. 根据权利要求 4 所述的计算方法, 其特征是所述利用时滞稳定判据求解时滞稳定上 限包括 : 子步骤 A1 : 将时滞稳定判据等价变换为 其中, Y1和 Y2为附加矩阵, 并且 Y1=Y1T 0, Y2=Y2T 0, 子 步 骤 A2 : 以 v 最 小 为 目 标, 以 时 滞 稳 定 判 据 等 价 变 换 后 的 矩 阵 不 等 式 以及 为 约束条件, 计算待定矩阵 P、 Q、 R、 Z1、 Z2、 Y1和 Y2, 进而得到时滞稳定上限 h。 权 利 要 求 书 CN 10。

7、3838965 A 3 1/11 页 4 基于广义特征值的时滞稳定上限计算系统及其计算方法 技术领域 0001 本发明属于电力系统稳定分析技术领域, 尤其涉及一种基于广义特征值的 0002 时滞稳定上限计算系统及其计算方法。 背景技术 0003 基于广域信息的反馈控制器的设计存在时滞问题, 而这必然会造成控制器的控制 效果降低, 甚至出现负阻尼的情况。因此, 迫切需要对系统的时滞稳定上限进行研究。 0004 目前国内外学者在系统时滞稳定上限的研究方面取得了大量有益成果, 主要可 以分为 3 类 : 时域法、 频域法和直接法。时域法可判定系统在特定场景下是否稳定, 但在 稳定程度、 时滞稳定上限。

8、等信息的获取方面还需要进一步研究。频域法通过在实数空间 中搜索时滞系统的关键特征根, 能够在一定程度上揭示时滞系统变化规律, 但计算量较 大, 求解速度有待提升。直接法借助 Lyapunov 理论和线性矩阵不等式 (Linear Matrix Inequality,LMI) 技术, 可同时考虑时滞的随机变动、 存在切换环节等情况, 适用范围更为 广泛, 但该方法具有一定的保守性。 0005 针对以上问题, 本发明提出一种基于广义特征值的时滞稳定上限计算系统及其计 算方法。本发明从电力系统实际出发, 能够有效解决直接法求解时滞上限过程中保守性较 高的问题。首先利用读入的广域信号数据建立时滞系统状。

9、态方程, 并在有效保留系统中低 频振荡成分的基础上降阶系统。 其次, 该发明形成基于改进自由权矩阵的时滞稳定判据, 在 此基础上, 将时滞稳定判据等价变换, 利用广义特征值法求解系统的时滞稳定上限。基于 IEEE4 机 11 节点系统和 IEEE16 机 68 节点系统仿真表明, 本发明可以较好的降低传统方法 的保守性, 具有很好的有效性和正确性。 发明内容 0006 本发明的目的在于, 提供一种基于广义特征值的时滞稳定上限计算系统及其计算 方法, 用于解决直接法求解时滞上限过程中保守性较高的问题。 0007 为了实现上述目的, 本发明提出的技术方案是, 一种基于广义特征值的时滞稳定 上限计算。

10、系统, 其特征是所述系统包括顺序相连的数据采集模块、 时滞系统处理模块、 时滞 上限求解模块和结果输出模块 ; 0008 所述数据采集模块用于采集建立时滞系统状态方程所需的网络结构参数、 发电机 频率和发电机功角, 并将采集的数据发送至时滞系统处理模块 ; 0009 所述时滞系统处理模块用于建立时滞系统状态方程, 并对时滞系统状态方程中的 参数矩阵进行降阶处理 ; 0010 所述时滞上限求解模块用于生成基于改进自由权矩阵的时滞稳定判据, 并利用时 滞稳定判据求解时滞稳定上限 ; 0011 所述结果输出模块用于输出时滞稳定上限结果。 0012 一种基于广义特征值的时滞稳定上限计算方法, 其特征是。

11、所述方法包括 : 说 明 书 CN 103838965 A 4 2/11 页 5 0013 步骤 1 : 采集建立时滞系统状态方程所需的网络结构参数、 发电机频率和发电机 功角 ; 0014 步骤 2 : 建立时滞系统状态方程, 并对时滞系统状态方程中的参数矩阵进行降阶 处理, 得到降阶后的时滞系统状态方程 ; 0015 步骤 3 : 生成基于改进自由权矩阵的时滞稳定判据 ; 0016 步骤 4 : 利用时滞稳定判据求解时滞稳定上限。 0017 所述降阶后的时滞系统状态方程为其中, x(t) 为降阶后的电力系统状态向量 ; 0018 为降阶后的电力系统状态矩阵 ; 0019 d为降阶后的电力系。

12、统时滞矩阵 ; 0020 为降阶后的电力系统状态量对应的状态值 ; 0021 d(t) 为时滞, 0 d(t) h 且 0022 h 为时滞稳定上限且 h0 ; 0023 为时滞最大变化率。 0024 所述基于改进自由权矩阵的时滞稳定判据为 : 0025 0026 其中, 0027 0028 2=N+M-N+S-M-S ; 0029 Ac=AAd0 ; 0030 N、 M 和 S 为改进自由权矩阵 ; 0031 改进自由权矩阵 N 满足 0032 改进自由权矩阵 S 满足 0033 改进自由权矩阵 M 满足 说 明 书 CN 103838965 A 5 3/11 页 6 0034 ?1(t)=。

13、xT(t)xT(t-d(t)xT(t-h)T; 0035 P、 Q、 R、 Z1和 Z2为待定矩阵。 0036 所述利用时滞稳定判据求解时滞稳定上限包括 : 0037 子步骤 A1 : 将时滞稳定判据等价变换为 0038 0039 其中, Y1和 Y2为附加矩阵, 并且 Y1=Y1T 0, Y2=Y2T 0, 0040 0041 子步骤 A2 : 以 v 最小为目标, 以时滞稳定判据等价变换后的矩阵不等式 以及 为 约束条件, 计算待定矩阵 P、 Q、 R、 Z1、 Z2、 Y1和 Y2, 进而得到时滞稳定上限 h。 0042 本发明利用改进自由权矩阵建立时滞系统稳定判据, 借助广义特征值法对。

14、系统时 滞上限进行求解, 其能够有效降低时滞稳定上限计算过程中的保守性, 具有很好的正确性 和有效性。 附图说明 0043 图 1 是基于广义特征值的时滞稳定上限计算系统结构图 ; 0044 图 2 是 IEEE4 机 11 节点系统结构图 ; 0045 图 3 是 IEEE4 机 11 节点 SMA 降阶前后特征根对比图 ; 其中, (a) 是全阶开环状态 矩阵和降阶开环状态矩阵的特征根对比图, (b) 是全阶闭环状态矩阵和降阶闭环状态矩阵 的特征根对比图 ; 0046 图 4 是 IEEE4 机 11 节点不同延时情况下发电机 1-4 相对功角动态响应曲线图 ; 0047 图 5 是 IE。

15、EE4 机 11 节点不同延时情况下发电机 2-3 相对功角动态响应曲线图 ; 0048 图 6 是 IEEE4 机 11 节点系统 G1 与 G4 功角差各时滞时间下的阻尼比结果表 ; 0049 图 7 是 IEEE4 机 11 节点系统 G2 与 G3 功角差各时滞时间下的阻尼比结果表 ; 0050 图 8 是 IEEE16 机 68 节点系统结构图 ; 0051 图 9 是 IEEE16 机 68 节点 Schur 降阶前后频率响应对比图 ; 其中, (a) 是全阶开环 状态矩阵和降阶开环状态矩阵的特征根对比图, (b) 是全阶闭环状态矩阵和降阶闭环状态 说 明 书 CN 1038389。

16、65 A 6 4/11 页 7 矩阵的特征根对比图 ; 0052 图 10 是 IEEE16 机 68 节点不同延时情况下发电机 1-16 相对功角动态响应曲线 图 ; 0053 图 11 是 IEEE16 机 68 节点不同延时情况下发电机 3-14 相对功角动态响应曲线 图 ; 0054 图 12 是 IEEE16 机 68 节点不同延时情况下发电机 10-15 相对功角动态响应曲线 图 ; 0055 图 13 是 16 机系统 G1 与 G16 功角差各时滞时间下的阻尼比表 ; 0056 图 14 是 16 机系统 G3 与 G14 功角差各时滞时间下的阻尼比表 ; 0057 图 15 。

17、是 16 机系统 G10 与 G15 功角差各时滞时间下的阻尼比表。 具体实施方式 0058 下面结合附图, 对优选实施例作详细说明。 应该强调的是, 下述说明仅仅是示例性 的, 而不是为了限制本发明的范围及其应用。 0059 实施例 1 0060 图 1 是本发明提供的基于广义特征值的时滞稳定上限计算系统结构图。如图 1 所 示, 本发明提供的基于广义特征值的时滞稳定上限计算系统结构图包括顺序相连的数据采 集模块、 时滞系统处理模块、 时滞上限求解模块和结果输出模块。 0061 数据采集模块用于采集网络结构参数、 发电机频率和发电机功角, 并将采集的数 据发送至时滞系统处理模块。其中, 网络。

18、结构参数是用于建立时滞系统状态方程所需的参 数, 而发电机频率和发电机功角则是时滞系统状态方程中的状态向量。 0062 时滞系统处理模块用于建立时滞系统状态方程, 并对时滞系统状态方程中的参数 矩阵进行降阶处理, 得到降阶后的时滞系统状态方程。 0063 时滞上限求解模块用于生成基于改进自由权矩阵的时滞稳定判据 ; 并对所述时滞 稳定判据进行等价变换, 利用广义特征值法求解系统的时滞稳定上限。 0064 结果输出模块用于输出时滞稳定上限结果。 0065 本发明提供的基于广义特征值的时滞稳定上限计算方法包括 : 0066 步骤 1 : 采集建立时滞系统状态方程所需的网络结构参数、 发电机频率和发。

19、电机 功角。 0067 网络结构参数包括时滞系统中线路的阻抗、 导纳、 发电机的内阻抗和负荷的等值 阻抗。发电机频率用于计算时滞系统转速, 发电机功角和时滞系统转速的变化量为时滞系 统状态方程的状态向量。 0068 步骤 2 : 建立时滞系统状态方程, 并对时滞系统状态方程中的参数矩阵进行降阶 处理, 得到降阶后的时滞系统状态方程。 0069 多输入多输出电力系统的状态方程可表示为 : 0070 0071 其中, x (t) Rn为电力系统状态向量, u (t) Rm为电力系统控制输入向量, 说 明 书 CN 103838965 A 7 5/11 页 8 A Rnn为电力系统状态矩阵, B R。

20、nm为电力系统控制矩阵。 0072 通过状态反馈后得到相应的闭环系统为 : 0073 0074 其 中, C 为 闭 环 状 态 矩 阵, 当 系 统 经 过 状 态 反 馈 时, 闭 环 状 态 矩 阵 为 C =A +B K 1, 其中, K1 R mn 为各附加控制器的综合状态反馈矩阵。 0075 在实际电力系统中, 控制输入向量通过 SCADA/WAMS 系统向各控制器传达, 信号传 递过程必然存在一定时滞, 则相应闭环系统可描述为 : 0076 0077 由式 (3) 可知, 电力系统时滞矩阵 A d=B K1。对于含时滞环节的系统, 其状 态方程有如下形式 : 0078 0079 。

21、其中, 矩阵A和Ad分别为电力系统状态矩阵和电力系统时滞矩阵, h为时滞稳 定上限。公式 (4) 中时滞 d(t) 满足条件 : 0080 0 d(t) h (5) 0081 0082 公式 (4)- 式 (6) 即为时滞系统状态方程, 为时滞最大变化率。 0083 实际系统分析中, 仅对低频振荡特征根相应的特征向量中与系统状态量 (发 电机转速变化量) 和 (发电机功角变化量) 相对应的元素感兴趣, 以便了解在 i中所 含的该振荡模式分量的相对幅值及相位。常用的系统降阶方法包括 SMA 降阶方法和 Schur 平衡降阶方法。由于系统降阶方法已经是本领域常用的方法, 因此本发明只以选择模式分 。

22、析法 SMA 为例, 对系统降阶做简单介绍。 0084 对于系统状态方程将其按下式划分 0085 0086 其中, X1=T,T 为保留变量, X2为其他变量, 待消去。 0087 由式 (7) 可消去 X2, 得 : 0088 0089 其中, I 为单位阵, p 为微分算子。 0090 将上式改写为 0091 0092 其中, Xr=X1为保留变量, Ar(p) 为运算形式的 “降阶” 系统系数阵。 0093 由式 (7)-(9) 可以得到两个重要性质 : 说 明 书 CN 103838965 A 8 6/11 页 9 0094 (1) 如果 p=1(i=1,2,N) 为式 (7) 相应系。

23、统特征根, 即 |iI-A|=0, 则 p=i 也为式(8)或(9)形式上降阶的系统特征根, 即亦有|iI-Ar(i)|=0, 特征根不发生变化, 系统模式不变。 0095 (2) 对于原系统, i的特征向量 ui, 有 Aui=iui。设降阶系统 i相应的特征向 量为uri, 即Ar(i)uri=iuri, 则uri和ui中保留变量Xr相对应的元素相等, 即特征向量的相 应元素不变。因此, 在 Xr保留变量处去观察同一模式 i的振荡时, 相对幅值即相位不变, 或者说模态不变。这样, 所关心频带输入输出特性被完整的保留下来。 0096 通过降阶处理, 可以得到降阶后的时滞系统状态方程为 : 0。

24、097 0098 公式 (10) 中, x(t) 是降阶后的电力系统状态向量, 为降阶后的电力系统状态 矩阵, d为降阶后的电力系统时滞矩阵,为降阶后的电力系统状态量对应的状态值, d(t)、 h 和 的含义同公式 (4) 且满足式 (5) 和式 (6) 。 0099 步骤 3 : 生成基于改进自由权矩阵的时滞稳定判据。 0100 构造如下形式的 Lyapunov-Krasovskii 泛函 (李雅普诺夫 - 克拉索夫斯基泛函) : 0101 0102 公式 (11) 中, P=PT0, Q=QT 0, R=RT 0 和 Zi=ZiT 0(i=1,2) 都是待定矩阵。 0103 由 Newto。

25、n-Leibniz 公式 (牛顿 - 莱布尼兹公式) 可知, 对于改进自由权矩阵 N、 S 和 M, 式 (12)- 式 (14) 成立 : 0104 0105 0106 0107 其中, 1(t)=xT(t)xT(t-d(t)xT(t-h)T。 0108 计算式 (11) 中 V(x) 关于 t 的导数为 : 0109 0110 将式 (10) 中第一个式子及式 (12)- 式 (14) 代入式 (15), 并加入必要的松散项可 说 明 书 CN 103838965 A 9 7/11 页 10 得 : 0111 0112 0113 0114 0115 0116 上述公式中, 0117 2=N。

26、+M-N+S-M-S。 0118 Ac=AAd0。 0119 0120 N1、 N2和 N3为矩阵 N 的分块矩阵且 N1的维度与 PA+ATP+Q+R 的维度相等, N2的维度 与的维度相等, N3的维度与 R 的维度相等。 0121 S1、 S2和 S3为矩阵 S 的分块矩阵且 S1的维度与 PA+ATP+Q+R 的维度相等, S2的维度 与的维度相等, S3的维度与 R 的维度相等。 0122 M1、 M2和 M3为矩阵 M 的分块矩阵且 M1的维度与 PA+ATP+Q+R 的维度相等, M2的维度 与的维度相等, M3的维度与 R 的维度相等。 0123 考虑到式 (11) 中 Zi=。

27、ZiT 0, (i=1,2), 因此公式 (18)- 式 (20) 中 i=iT 0, (i=1,2,3)。若 +s 0, 则式 (16) 中对于 +s, 利用 Schur 补后可得式 说 明 书 CN 103838965 A 10 8/11 页 11 (10) 表征的时滞系统稳定判据如下 : 0124 对于给定标量 h0 和 , 若存在 P=PT0, Q=QT 0, R=RT 0 和 Zi=ZiT 0(i=1,2), 和使得如下 线性矩阵不等式成立 : 0125 0126 则对于同时满足条件式 (5) 和式 (6) 的时滞系统 (10) 渐进稳定。 0127 步骤 4 : 利用时滞稳定判据求。

28、解时滞稳定上限。 0128 式 (21) 表征的线性矩阵不等式仅能判定系统是否稳定, 而无法获取系统时滞稳 定上限等信息。 考虑到广义特征值法能够求解优化问题的全局最小值, 因此, 本发明提出利 用广义特征值法计算系统的时滞稳定上限。由于式 (21) 不是标准的广义特征值形式, 无法 直接利用广义特征值法进行求解, 因此, 做如下变换, 将式 (21) 同时左乘与右乘式 (22), 如 下 : 0129 0130 可得 : 0131 0132 其中, v=1h, I 为单位阵。 0133 为 了 求 解 最 大 时 滞 稳 定 上 限 h, 即 最 小 的 v, 本 发 明 引 入 附 加 矩。

29、 阵 Yi=YiT 0(i=1,2), 该矩阵需满足式 (24) 成立 : 0134 0135 再将式 (24) 代入式 (23) 可得 : 说 明 书 CN 103838965 A 11 9/11 页 12 0136 0137 由此, 时滞稳定上限问题已经转化为以下优化问题 : 0138 以 v 最小为目标, 以和 为约束条件, 计算待定矩阵 P、 Q、 R、 Z1、 Z2、 Y1和 Y2, 进而得到时滞 稳定上限 h。公式表示如下 : 0139 min v 0140 P, Q, R, Zi, Yi 0141 s.t. 0142 0143 通过求解式 (26), 可以计算出以矩阵 P、 Q、。

30、 R、 Z1、 Z2、 Y1和 Y2为变量, 以式 (25) 和 式 (26) 为约束的最小 v。最终, 利用 h=1/v 可以求出时滞稳定上限。 0144 实施例 2 0145 基于 MATLAB 仿真软件搭建的图 2 所示的 IEEE4 机 11 节点系统, 发电机采用 6 阶 详细模型, 励磁系统采用快速励磁, 基准模型下的负荷采用 50% 恒阻抗和 50% 恒电流模型。 首先, 通过模态分析法得到四机系统的状态矩阵, 并利用 SMA 方法分别对开、 闭环状态矩阵 进行降阶, 如图 3 所示。 0146 图 3(a) 是全阶开环状态矩阵和降阶开环状态矩阵的特征根对比图, (b) 是全阶闭。

31、 环状态矩阵和降阶闭环状态矩阵的特征根对比图。由图 3 可以看出, 开、 闭环状态矩阵在利 用 SMA 进行降阶后, 均保留了系统中低频振荡成分。因此, 利用降阶后系统矩阵能够求解系 说 明 书 CN 103838965 A 12 10/11 页 13 统允许的最大时滞。将降阶后状态矩阵 和时滞矩阵 d代入式 (26), 求得最大时滞边界 h=281.88ms。其中, 降阶后状态矩阵 如下 : 0147 0148 降阶后时滞矩阵 d如下 : 0149 0150 时滞分别设置为 h1=50ms, h2=100ms, h3=288.88ms 和 h4=350ms。所观察的物理量为 发电机 G1 和。

32、 G4 之间的功角差, 以及 G2 和 G3 之间的功角差, 分别如图 4 和图 5 所示。 0151 由图 4 和图 5 均可看出, 在未加广域阻尼控制时, 发电机之间的功角发生了严重 的低频振荡。加入广域阻尼控制后, 在不考虑时滞的情况下, 低频振荡可得到有效抑制 ; 然 而随着时滞的增加, 阻尼效果随之减弱, 在最大时滞 288ms 的情况下, 虽然仍存在一定的阻 尼, 但其阻尼比已经降到 10% 以下, 说明此时控制器已经不满足控制要求。 0152 利用 prony 算法对各时滞下的功角差曲线进行阻尼比分析, 结果如图 6 和图 7 所 示。 0153 实施例 3 0154 基于 MA。

33、TLAB 仿真软件搭建的图 8 所示的 IEEE16 机 68 节点系统进一步考查基于 广义特征值法的最大时滞求解方法的有效性和通用性。该系统重要联络线为区域 4 和区域 5 的联络线 1-2, 1-27 和 8-9。发电机采用 6 阶详细模型, 励磁采用 IEEE-DC1 型励磁, 负荷 模型 15% 恒有功功率, 25% 的恒有功电流和 15% 的恒无功功率, 25% 的恒无功功率和 60% 恒 阻抗。首先利用 Schur 平衡降阶方法, 在保证所关心频带的输入输出特性保持不变的前提 下, 对系统进行降阶, 如图 9 所示。 0155 图 9(a) 中实线为开环全阶系统的频率响应, 虚线为。

34、开环降阶系统的频率响应, 可 以看出, 降阶系统和全阶系统的输入输出特性相同。图 9(b) 中实线为闭环全阶系统的频率 响应, 虚线为闭环降阶系统的频率响应, 可以看出, 闭环降阶系统同样保留了全阶系统的特 性。因此, 利用降阶后的系统矩阵求解全阶系统的时滞稳定上限具有有效性和可行性。 0156 降阶后的状态矩阵re与降阶后的时滞矩阵dre代入式(26), 可得最大的时滞上 限为 h=93ms。时滞分别设置为 h1=50ms 和 h2=93ms。以发电机 G1 和 G16 之间的功角差, G3 说 明 书 CN 103838965 A 13 11/11 页 14 和G14之间的功角差与G10和。

35、G15之间的功角差, 为所观察的物理量, 分别如图10, 图11和 图 12 所示。其中, 降阶后的状态矩阵 re如下 : 0157 0158 降阶后的时滞矩阵 dre如下 : 0159 0160 随着时滞的增加, 阻尼效果随之减弱, 在最大时滞 93ms 的情况下, 虽然存在一定 的阻尼效果, 但是其阻尼比已经降到 10% 以下, 此时控制器同样也不能满足控制的要求。利 用 prony 方法对各时滞下的功角差曲线进行阻尼比分析, 结果如图 13、 图 14 和图 15 所示。 0161 以上所述, 仅为本发明较佳的具体实施方式, 但本发明的保护范围并不局限于此, 任何熟悉本技术领域的技术人员。

36、在本发明揭露的技术范围内, 可轻易想到的变化或替换, 都应涵盖在本发明的保护范围之内。因此, 本发明的保护范围应该以权利要求的保护范围 为准。 说 明 书 CN 103838965 A 14 1/8 页 15 图 1 图 2 图 3 说 明 书 附 图 CN 103838965 A 15 2/8 页 16 图 4 图 5 说 明 书 附 图 CN 103838965 A 16 3/8 页 17 图 6 图 7 说 明 书 附 图 CN 103838965 A 17 4/8 页 18 图 8 说 明 书 附 图 CN 103838965 A 18 5/8 页 19 图 9 说 明 书 附 图 CN 103838965 A 19 6/8 页 20 图 10 图 11 说 明 书 附 图 CN 103838965 A 20 7/8 页 21 图 12 图 13 图 14 说 明 书 附 图 CN 103838965 A 21 8/8 页 22 图 15 说 明 书 附 图 CN 103838965 A 22 。

展开阅读全文
相关资源
猜你喜欢
相关搜索

当前位置:首页 > 物理 > 计算;推算;计数


copyright@ 2017-2020 zhuanlichaxun.net网站版权所有
经营许可证编号:粤ICP备2021068784号-1