基于慢模式特征值对于线路参数灵敏度的弱连接识别方法

摘要

本发明公开了属于电力系统分析和控制技术领域的一种基于慢模式相关特征值对于线路参数灵敏度的弱连接识别方法，其特征在于：首先计算系统的慢模式相关特征值组对各条线路参数变化灵敏度，再选取灵敏度数值比较大的线路作为系统的弱连接。主要按照以下几个步骤进行：计算机初始化，输入包括系统网络参数、节点电压和功率信息的初始化数据；运用牛顿拉夫逊法计算系统潮流；形成系统线性化的动态方程；运用QR方法求取系统特征值，并根据大小进行排序，选取系统慢模式相关特征值组；计算慢模式相关特征值对于系统各条线路参数变化的灵敏度，并按大小进行排序；选取系统弱连接。本发明的有益效果是相比于以前首先处理节点的方法，能够直接处理线路。

主权项

1.一种慢模式相关特征值对于线路参数灵敏度的弱连接识别方法，其特征在于，具体算法步骤如下：步骤(1)计算机初始化；1.1)输入以下参数：电力系统的网络参数，包括：输电线路的串联电阻、串联电抗、并联电导和并联电纳；变压器的变比和阻抗；并联在输电线路上的电容器和电抗器的阻抗；发电机的动态参数包括发电机的惯性时间常数；发电机的暂态电抗；1.2)设定以下参数：各母线节点的电压和有功功率和无功功率的实时量测值；步骤(2)根据步骤(1)的参数运用牛顿拉夫逊方法计算系统潮流，并记录至少包括系统各个节点电压的幅值和相位当前潮流断面信息；步骤(3)根据步骤(1)、步骤(2)的结果形成多机系统线性化的动态模型方程其中关键是计算矩阵A，其中表示发电机功角偏差向量的二阶导数；形成按以下思路进行：为了计算和推导的方便，本发明采用了简化的电力系统模型，具体模型如下：研究电力系统同调时有一个基本的假设：发电机模型的复杂程度不改变同调机群之间基本的振荡模态；所以本发明中选取很简单的电力系统模型来设计算法，其中发电机取经典二阶模型表同时忽略阻尼；而负荷模型则采用恒阻抗负荷模型；发电机取经典二阶模型并忽略阻尼的情况下，电力系统的机电振荡模型如下：$\left\{\begin{array}{c}2H\frac{d^2\delta }{{\mathrm{dt}}^2}=P_m-P_e\\ P_e=f(E_d,U_G,x_d^{\prime },\delta ,{\theta }_{\mathrm{UG}})\end{array}\right.$第一个等式中：δ为所有发电机功角组成的列相量，$\delta ={\left[\begin{array}{cccc}{\delta }_1& {\delta }_2& ···& {\delta }_{n_G}\end{array}\right]}^T,$n_G为系统中发电机的个数；H为所有发电机惯性时间常数组成的对角矩阵，$H=\mathrm{diag}\left\{\begin{array}{cccc}{H}_1& H_2& ···& H_{n_G}\end{array}\right\};$P_m为所有发电机的机械输入功率组成的列向量，$P_m=\left[\begin{array}{cccc}{P}_{m1}& P_{m2}& ···& P_{{\mathrm{mn}}_G}\end{array}\right],$为常数相量；P_e为所有发电机的电磁输出功率组成的列向量，$P_e=\left[\begin{array}{cccc}{P}_{e1}& P_{e2}& ···& P_{e{n}_G}\end{array}\right];$第二个等式中：$P_{\mathrm{ei}}=\frac{E_{\mathrm{di}}{U}_{\mathrm{Gi}}}{x_{\mathrm{di}}^{\prime }}\sin ({\delta }_i-{\theta }_{\mathrm{UGi}}),$E_di为发电机i的暂态电抗后电势，x′_di为发电机i的暂态电抗，U_gi为发电机i的机端电压；U_gi为发电机i的机端电压，θ_Ugi为发电机i的机端电压的相角；所以f(E_d，U_G，x′_d，δ，θ_UG)是一个矩阵表达式，其形式如下：$f(E_d,U_G,x_d^{\prime },\delta ,{\delta }_{\mathrm{UG}})={\left[\begin{array}{cccc}\frac{E_{d1}{U}_{g1}}{x_{d1}^{\prime }}\sin ({\delta }_1-{\delta }_{\mathrm{Ug}1})& \frac{E_{d2}{U}_{g2}}{x_{d2}^{\prime }}\sin ({\delta }_2-{\delta }_{\mathrm{Ug}2})& ···& \frac{E_{\mathrm{dnG}}{U}_{\mathrm{gnG}}}{x_{\mathrm{dnG}}^{\prime }}\sin ({\delta }_{\mathrm{nG}}-{\delta }_{\mathrm{Ug}{n}_G})\end{array}\right]}^T$其由$E_d={\left[\begin{array}{cccc}{E}_{d1}& E_{d2}& ···& E_{{\mathrm{dn}}_G}\end{array}\right]}^T,$$x_d^{\prime }={\left[\begin{array}{cccc}{x}_{d1}^{\prime }& x_{d2}^{\prime }& ···& x_{d{n}_G}^{\prime }\end{array}\right]}^T,$$U_G={\left[\begin{array}{cccc}{U}_{g1}& U_{g2}& ···& U_{{\mathrm{gn}}_G}\end{array}\right]}^T,$${\theta }_{\mathrm{UG}}={\left[\begin{array}{cccc}{\theta }_{\mathrm{Ug}1}& {\theta }_{\mathrm{Ug}2}& ···& {\theta }_{{\mathrm{Ugn}}_G}\end{array}\right]}^T;$每台发电机的暂态电抗后电势及当前功角按照如下思路确定：假设发电机节点i的功率注入为${\stackrel{*}{S}}_{\mathrm{gi}}=P_{\mathrm{gi}}+j{Q}_{\mathrm{gi}},$同时机端电压模制和相角为U_gi，δ_Ugi，那么首先求出此节点的等效注入电流${\stackrel{*}{I}}_{\mathrm{gi}}=\left(\frac{S_{\mathrm{gi}}}{U_{\mathrm{gi}}\angle {\theta }_{\mathrm{gi}}}\right),$那么：$E_{\mathrm{di}}\angle {\delta }_i=U_{\mathrm{gi}}\angle {\theta }_{\mathrm{Ugi}}+{\stackrel{*}{I}}_{\mathrm{gi}}*\left({\mathrm{jx}}_{\mathrm{di}}\right),$下文中某一个变量有上标星号(*)表示此变量是一个复数矢量；对于系统网络的处理，对于发电机和负荷分别进行简化处理；对于发电机，经典二阶模型下的发电机相当于有内阻的电压源(内阻是x′_di，而电压是E_di)，那么经过诺顿等值将这个电压源转换成一个电流源：阻抗仍为x′_di，$I_{\mathrm{gi}}=\frac{E_{\mathrm{di}}}{x_{\mathrm{di}}^{\prime }},$对于负荷节点k进行如下处理：假设当前潮流下，节点的电压模值为V_{Load_k}，而负荷有功和无功分别为P_k和Q_k，那么等效负荷阻抗为$Z_{\mathrm{Load}\_k}={\left(\frac{V_{\mathrm{Load}\_k}^2}{P_k+{\mathrm{jQ}}_k}\right)}^{*}=R_{\mathrm{Load}\_k}+j{X}_{\mathrm{Load}\_k},$相应得负荷等效导纳为$Y_{\mathrm{Load}\_k}=\frac{1}{Z_{\mathrm{Load}\_k}}=G_{\mathrm{Load}\_k}+j{B}_{\mathrm{Load}\_k};$把负荷节点的等效阻抗以及发电机的暂态电抗加以考虑之后，系统的节点导纳矩阵变为：Y＝G_Y+jB_Y；其中：Y、G_Y、B_Y都是n×n维的矩阵；令$\stackrel{*}{U}=U_x+j{U}_y$表示个节点电压列向量，每个元素为各个节点的电压矢量；$\stackrel{*}{I}=I_x+j{I}_y$代表各个节点的注入电流向量，下文中，某个列向量上标是一个星号(*)都表示此列向量每个元素都是复数矢量；其中：U_x为各个节点电压实部组成的向量，U_y为各个节点电压虚部组成的向量，I_x为各个节点注入电流实部组成的向量，I_y为个节点注入电流虚部组成的向量；经过上述简化之后，只有发电机节点具有注入电流，其他的节点没有电流注入；令$Y^{\prime }=\left[\begin{array}{cc}{G}_Y& -B_Y\\ B_Y& G_Y\end{array}\right],$$U^{\prime }=\left[\begin{array}{c}{U}_x\\ U_y\end{array}\right],$$I^{\prime }=\left[\begin{array}{c}{I}_x\\ I_y\end{array}\right],$那么有：$I^{\prime }=\left[\begin{array}{c}{I}_x\\ I_y\end{array}\right]=Y^{\prime }U=\left[\begin{array}{cc}{G}_Y& -B_Y\\ B_Y& G_Y\end{array}\right]\left[\begin{array}{c}{U}_x\\ U_y\end{array}\right]$令系统所有节点的电压模值向量为U＝[U₁ U₂ ... U_N]，相角向量为θ_U＝{θ_u1 θ_u2 ... θ_uN}^T，N为系统所有节点的个数，那么：U_x＝[U₁cosθ_u1 U₂cosθ_u2 ... U_ncosθ_uN]^T，U_y＝[U₁sinθ_u1 U₂sinθ_u2 ... U_nsinθ_uN]^T令：$Z^{\prime }=Y^{\prime -1}={\left[\begin{array}{cc}{G}_Y& -B_Y\\ B_Y& G_Y\end{array}\right]}^{-1}$那么有：$\left[\begin{array}{c}{U}_x\\ U_y\end{array}\right]=Z\left[\begin{array}{c}{I}_x\\ I_y\end{array}\right]---\left(2\right)$I_x，I_y都只有发电机节点对应元素非零，其他元素均为零；令：$I_{\mathrm{Gx}}={\left[\begin{array}{cccc}{I}_{G1}\sin {\delta }_1& I_{G2}\sin {\delta }_2& ···& I_{{\mathrm{Gn}}_G}\sin {\delta }_{n_G}\end{array}\right]}^T$(1，...，n_G表示发电机序号)$I_{\mathrm{Gy}}={\left[-,\begin{array}{cccc}{I}_{G1}\cos {\delta }_1& I_{G2}\cos {\delta }_2& ···& I_{{\mathrm{Gn}}_G}\cos {\delta }_{n_G}\end{array}\right]}^T$(1，...，n_G表示发电机序号)那么I_x，I_y可以用以下的形式表示：I_x＝T×I_Gx，I_y＝T×I_Gy其中：T是一个N×n_G维的矩阵；T的每一行对应系统的一个节点(编号为k)，如果本节点没有发电机接入，即非发电机节点，那么本行元素全为0，即T(k，：)＝[0...0]；如果本节点有发电机接入，即为发电机节点，假设此发电机序号为i，那么T(k，i)＝1，本行其他元素均为0，即T(k，：)＝[0...1...0]；对(2)式两边线性化得：$\left[\begin{array}{cc}\cos {\theta }_U& -U\sin {\theta }_U\\ \sin {\theta }_U& U\cos {\theta }_U\end{array}\right]\left[\begin{array}{c}\mathrm{\Delta U}\\ \Delta {\theta }_U\end{array}\right]=Z\left[\begin{array}{c}{I}_{\mathrm{x\delta }}\\ I_{\mathrm{y\delta }}\end{array}\right]\mathrm{\Delta \delta }---\left(3\right)$其中：ΔU，Δθ_U，Δδ分别为U，θ_U，δ的偏列差量，N和n_G分别为系统所有节点的个数以及发电机节点的个数，表达式中的各个矩阵的形式如下：cosθ_U＝diag{cosθ_U1 cosθ_U2...cosθ_UN}；sinθ_U＝diag{sinθ_U1 sinθ_U2...sinθ_UN}；Ucosθ_U＝diag{U₁cosθ_U1 U₂cosθ_U2...U_Ncosθ_UN}；Usinθ_U＝diag{U₁sinθ_U1 U₂sinθ_U2...U_Nsinθ_UN}；$I_{\mathrm{x\delta }}=\frac{{\partial I}_x}{\partial \delta },$$I_{\mathrm{y\delta }}=\frac{\partial I_y}{\partial \delta }$均为N×n_G维的矩阵，可知$I_{\mathrm{x\delta }}=\frac{\partial I_x}{\partial \delta }=T\times \frac{{\partial I}_{\mathrm{Gx}}}{\partial \delta },$$I_{\mathrm{y\delta }}=\frac{\partial I_y}{\partial \delta }=T\times \frac{{\partial I}_{\mathrm{Gy}}}{\partial \delta },$而：$\frac{{\partial I}_{\mathrm{Gx}}}{\partial \delta }=\mathrm{diag}\left\{\begin{array}{cccc}{I}_{G1}\cos {\delta }_1& I_{G2}\cos {\delta }_2& ...& I_{\mathrm{GnG}}\cos {\delta }_{n_G}\end{array}\right\}$$\frac{{\partial I}_{\mathrm{Gy}}}{\partial \delta }=\mathrm{diag}\left\{\begin{array}{cccc}{I}_{G1}\sin {\delta }_1& I_{G2}\sin {\delta }_2& ...& I_{\mathrm{GnG}}\sin {\delta }_{n_G}\end{array}\right\}$令：$D=\left[\begin{array}{cc}\cos {\theta }_U& -U\sin {\theta }_U\\ \sin {\theta }_U& U\cos {\theta }_U\end{array}\right],$$E=\left[\begin{array}{c}{I}_{\mathrm{x\delta }}\\ I_{\mathrm{y\delta }}\end{array}\right],$E是2N×n_G维的，那么：$\left[\begin{array}{c}\mathrm{\Delta U}\\ \Delta {\theta }_U\end{array}\right]=D^{-1}\mathrm{ZE\Delta \delta }---\left(4\right)$对于式(1)所示的电力系统机电振荡模型在平衡点处进行线性化，可得：$\frac{d^2\mathrm{\Delta \delta }}{{\mathrm{dt}}^2}=B\left[\begin{array}{c}\Delta U_G\\ \Delta {\theta }_{\mathrm{UG}}\end{array}\right]+\mathrm{C\Delta \delta }---\left(5\right)$其中：ΔU_G，Δθ_UG，Δδ为U_G，θ_UG，δ的偏差量，$B={\left(2H\right)}^{-1}\left[\begin{array}{cc}-\frac{\partial g}{\partial U_G}& 0\\ 0& -\frac{\partial f}{\partial {\theta }_{\mathrm{UG}}}\end{array}\right],$$C=-{\left(2H\right)}^{-1}\frac{\partial f}{\partial \delta }$有：ΔU_G＝T₂×ΔU，Δθ_UG＝T₂×Δθ_U，代入(5)式有：$\frac{d^2\mathrm{\Delta \delta }}{{\mathrm{dt}}^2}=B\left[\begin{array}{c}\Delta U_G\\ \Delta {\theta }_{\mathrm{uG}}\end{array}\right]+\mathrm{C\Delta \delta }=B\left[\begin{array}{cc}{T}_2& 0\\ 0& T_2\end{array}\right]\left[\begin{array}{c}\mathrm{\Delta U}\\ \Delta {\theta }_u\end{array}\right]+\mathrm{C\Delta \delta }={\mathrm{BT}}_2^{\prime }\left[\begin{array}{c}\mathrm{\Delta U}\\ \Delta {\theta }_U\end{array}\right]+\mathrm{C\Delta \delta }---\left(6\right)$其中：$T_2^{\prime }=\left[\begin{array}{cc}{T}_2& 0\\ 0& T_2\end{array}\right],$T₂＝T^T带入(4)式就得到了多机系统线性化模型为：$\frac{d^2\mathrm{\Delta \delta }}{{\mathrm{dt}}^2}=({\mathrm{BT}}_2^{\prime }{D}^{-1}\mathrm{ZE}+C)\mathrm{\Delta \delta }=\mathrm{A\Delta \delta }---\left(7\right)$步骤(4)求取矩阵A的特征值并选择慢模式相关特征值组；4.1)采用QR算法求取A的全部特征值；4.2)把A的特征值按照绝对值由小到大的顺序进行排序；假设排序后的特征值顺序如下：{λ₁ λ₂ ...}，其中n_G为发电机的台数，也正是A的维数；值得注意的是λ₁＝0；4.3)选取慢模式相关特征值组，关键是确定选取慢模式相关特征值的个数r；按照如下的方法选择r：比较所有的λ_i/λ_i+1，i＝2，3，...，n_G-1，令：${\lambda }_r/{\lambda }_{r+1}=\underset{i=\mathrm{2,3},..n_G-1}{\min }{\lambda }_i/{\lambda }_{i+1},$这样就确定了慢模式相关特征值的个数r；这种选择的依据是电力系统多时间尺度特性原理，λ_r/λ_r+1的值越小表示系统的多时间尺度特性越明显，系统的区域性也越明显，系统也越适合解列；经过这一步骤，就确定了慢模式相关特征值组σ_r＝{λ₁，λ₂，..λ_r}；步骤(5)求取慢模式相关特征值组对于每条线路参数变化灵敏度；线路就是指电力系统中实际的输电线路。首先对于所有节点进行编号：1，2，...，N，再对线路进行编号：1，2，...，Nl。一条线路主要有如下元素特征：起始节点编号、终了节点编号、线路电阻r、线路电抗x，线路电导g和线路电纳b，值得注意的是在选取线路起止节点时，哪一端点作为起始节点时无所谓的，所以在确定起始节点时随便指定其中一个端点作为起始节点即可。某条线路的编号是L，其起始节点编号为I，终了节点编号为J，线路的电阻和电抗分别为r和x，电导和电纳分别为g和b，求取慢模式相关特征值组对于线路L参数变化的灵敏度，可按以下两步进行：5.1)求取慢模式相关特征值组σ_r＝{λ₁，λ₂，..λ_r}中每个特征值对线路L参数变化的灵敏度，表示为{S_λ1^L S_λ2^L...S_λr^L}，S_λi^L表示特征值λ_i对于线路L参数变化的灵敏度；特征值λ_i对于线路L参数变化的灵敏度可按照以下思路进行：λ_i为A的一个特征值，λ_i对应的左右特征向量分别为φ_i(n_G×1)和ψ_i(1×n_G)；那么λ_i对于线路L的导纳b变化的灵敏度表达式如下：其中：$N^{,}=\left[\begin{array}{cc}0& M^T\\ M^T& 0\end{array}\right],$$M^{\prime }=\left[\begin{array}{cc}M& 0\\ 0& -M\end{array}\right],$M＝[...1...-1...]^T，其中1的位置是I，而-1的位置是J；λ_i对于线路L的导纳g变化的灵敏度表达式如下：其中：$N^{\prime \prime }=\left[\begin{array}{cc}0& M^T\\ M^T& 0\end{array}\right],$$M^{\prime \prime }=\left[\begin{array}{cc}M& 0\\ 0& M\end{array}\right];$假设线路参数变化时，r和x同比例变化，那么由于$g+\mathrm{jb}=\frac{1}{r+\mathrm{jx}}=\frac{r-\mathrm{jx}}{r^2+x^2},$线路的g和b是按照-r/x的比例变化的，以b的变化为基准(b＞＞g)，总结出λ_i对于线路L参数变化的灵敏度变化表达式如下：$S_{\mathrm{\lambda i}}^L=-S_{\mathrm{\lambda i}}^{\mathrm{Lb}}+\frac{r}{x}{S}_{\mathrm{\lambda i}}^{\mathrm{Lg}}---\left(10\right)$5.2)确定慢模式相关特征值组σ_r＝{λ₁，λ₂，..λ_r}对于线路L参数变化的灵敏度，按下式确定：$S_{\mathrm{\sigma r}}^L=\underset{i=1}{\overset{r}{\Sigma }}{S}_{\mathrm{\lambda i}}^L---\left(11\right)$S_σr^L即为慢模式相关特征值组σ_r＝{λ₁，λ₂，..λ_r}对于线路L参数变化的灵敏度，用这个量来衡量线路参数变化对于慢模式组的影响；经过这一步骤就能够得到所有线路的慢模式相关特征值组对于系统所有线路的参数的灵敏度信息；步骤(6)根据步骤5得到的慢模式相关特征值组对于线路参数变化的灵敏度信息来确定弱连接；6.1)把慢模式相关特征值组对于系统所有的Nl条线路的灵敏度按照绝对值从大到小的顺序进行排序，排序后的顺序为{S_σr¹S_σr²...S_σr^Nl}；6.2)选取弱连接线路，基本思路就是选取灵敏度比较大的一组线路作为系统弱连接线路；至于弱连接线路的条数wl可以由用户自己按照需要设定，也可以自动生成，本发明推荐一种自动生成wl的方法；对于一个多时间尺度特性比较明显的系统，弱连接线路与非弱连接线路的灵敏度数值差异会比较大；这种方法首先对于系统所有线路的灵敏度进行从大到小的排序，然后按照上述原理来选择弱连接线路的条数wl，选取wl的公式表达形式如下：$\mathrm{abs}(S_{\mathrm{\sigma r}}^{\mathrm{wl}}/S_{\mathrm{\sigma r}}^{\mathrm{wl}+1})=\underset{i=\mathrm{1,2},..\mathrm{Nl}}{\max }\mathrm{abs}(S_{\mathrm{\sigma r}}^i/S_{\mathrm{\sigma r}}^{i+1})---\left(12\right)$其中，符号abs表示取一个实数的绝对值。