CN103248049A - 含dfig风电场电力***电压稳定域切平面的获取方法 - Google Patents

Abstract

本发明公开了一种计算含DFIG风电场电力***电压稳定域切平面的获取方法，涉及电力***领域，包括对DFIG风电场进行数学建模后，对风电场的状态参数进行约束；采用了预测-校正的连续潮流算法确定一条静态电压稳定边界轨迹，在运用连续潮流法时，对预测步以及校正步的雅可比矩阵均进行了修改，将DFIG风电场的约束函数式加入了雅可比矩阵；最后通过对与该二维参数呈线性关系参变量的代换以及切平面点法方程的计算，构建出该含DFIG电力***安全域切平面的表达式，从而使DFIG风电场全注入后***电压稳定性分析直观、可靠，为用安全域分析电网静态电压稳定性提供数据支持。

Description

含DFIG风电场电力***电压稳定域切平面的获取方法

技术领域

本发明属于电力***分析及规划领域，涉及DFIG风电场、电力***静态电压稳定域切平面的获取方法。

背景技术

风力发电是风能利用的重要形式，风能是可再生、无污染、能量大、前景广的能源，大力发展清洁能源是世界各国的战略选择。风具有随机性和不确定性，这就使电网大规模消纳风电具有一定难度。而应用安全域的方法应用到含风电场的电力***中，可通过追踪***运行点与安全域边界来提供安全裕度和最优控制信息，从而使电力***在线实时安全监视、防御与控制更科学和更有效。文献《电力***安全域方法研究述评》概述了安全域分类与分层关系，以及安全域方法的相关应用。文献《Calculating Steady-State Operating Conditions forDoubly-Fed Induction Generator Wind Turbines》提出了一种双馈感应发电机（Doubly Fed Induction Generator，以下简称DFIG）模型的建立。目前，基于DFIG的风电场已经开始并网发电，但基于安全域切平面方法研究DFIG风电场注入后电力***电压稳定性的甚少，而电网对注入后电力***电压稳定性有很高要求，所以研究DFIG风电场电力***静态电压稳定域切平面对预警电网安全有着现实意义。

发明内容

（一）要解决的技术问题

本发明的目的在于解决大规模DFIG风电机组注入后，电力***电压稳定性研究的问题。

（二）技术方案

现实中，绝大部分风电场所使用的发电机机种一致，如该风电场使用DFIG风电机组，则该风电场将使用同一类型风电机组，而不会与其它类型发电装置混用。本发明提出了一种含DFIG电力***静态电压稳定安全域切平面求解的方法，该方法首先对DFIG进行数学建模，利用七个DFIG状态变量对模型进行约束；然后利用预测与校正的连续潮流法对电力***静态电压稳定的二维边界进行追踪，得到一条静态电压稳定边界轨迹，为确定含DFIG电力***静态电压稳定安全域切平面打下基础；最后通过对与该二维参数呈线性关系参变量的代换以及切平面点法式方程的计算，构建出该含DFIG电力***安全域切平面的表达式。

（三）有益效果

本发明解决了大规模DFIG风电机组注入后，电力***电压稳定性研究的问题。

附图说明

图1是本发明的稳定域切平面获取方法计算流程示意图。

具体实施方式

下面结合附图和实施例，对本发明的具体实施方式做进一步描述。以下实施例仅用于说明本发明，但不用来限制本发明的范围。

如图1所示的计算流程。

本发明包括以下步骤：

1）确定DFIG风电场的电力***数学模型

1.1）首先认为所有的双馈风电机组为相同型号的并且风电场周围的风速均匀分布，所以可以用一个双馈风机模型表示风电场。其中，双馈风电机的机械能输出表达式为

$P_m=\frac{1}{2}\mathrm{\ρ}{\mathrm{AC}}_p{U}^3---\left(1\right)$

其中，ρ为空气密度；A为风力机的扫掠面积；U为风速；C_p为风力机的风能利用系数，表明风轮从风中获得的有用风能比例，与尖速比T_tsr有关。

1.2）以下为含DFIG风电场电力***含参潮流方程式：

$f_1=(1+\mathrm{\λ})P_i+P_{\mathrm{im}}-V_i\underset{j\∈i}{\mathrm{\Σ}}{V}_j(G_{\mathrm{ij}}\mathrm{ocs}{\mathrm{\δ}}_{\mathrm{ij}}+B_{\mathrm{ij}}\sin {\mathrm{\δ}}_{\mathrm{ij}})=0$

$f_2=(1+\mathrm{\λ})Q_i+Q_{\mathrm{im}}-V_i\underset{j\∈i}{\mathrm{\Σ}}{V}_j(G_{\mathrm{ij}}\sin {\mathrm{\δ}}_{\mathrm{ij}}+B_{\mathrm{ij}}\cos {\mathrm{\δ}}_{\mathrm{ij}})=0$

f₃=V_scos(θ_s)+R_sI_scos(φ_s)-(X_s+X_m)I_ssin(φ_s)+X_mI_rsin(φ_r)=0

f₄=V_ssin(θ_s)-R_sI_ssin(φ_s)+(X_s+X_m)I_scos(φ_s)-X_mI_rcos(φ_r)=0

f₅=V_rcos(θ_r)-sX_mI_ssin(φ_s)-R_rI_rcos(φ_r)+s(X_s+X_m)I_rsin(φ_r)=0     （2）

f₆=V_rsin(θ_r)+sX_mI_scos(φ_s)-R_rI_rsin(φ_r)-s(X_s+X_m)I_rcos(φ_r)=0

f₇=P_m-V_sI_scos(θ_s-φ_s)+V_rI_rcos(θ_r-φ_r)=0

f₈=Q_m-V_sI_ssin(θ_s-φ_s)=0

$f_9=P_m+I_s^2{R}_s+I_r^2{R}_r+P_{\mathrm{GB},\mathrm{Loss}}-k{(1-s)}^3=0$

式中，V_i为节点电压大小，δ_i为节点电压相角，V_s为定子电压大小，θ_s为定子电压相角，V_r为转子电压大小，θ_r为转子电压相角，I_r为转子电流大小，φ_r为转子电流相角，I_s为定子电流大小，φ_s为定子电流相角，s为滑差。

将含有V_r,θ_r,I_S,φ_S,I_r,φ_r,s七个未知数的七个式子，即式（2）中f₃-f₉加到***的含参潮流方程中后，对其雅可比矩阵修改，再进行预测-校正的计算，以得到含DFIG电力***的静态临界稳定点。

2）用连续潮流法计算含DFIG电力***静态电压稳定边界初始点

在引入双馈风机后，***的含参潮流模型可以记为：

f(V_i,δ_i,V_r,θ_r,I_s,φ_s,I_r,φ_r,s,λ,P_m)=0   (3)

式中，V_i,δ_i,V_r,θ_r,I_s,φ_s,I_r,φ_r,s为***的状态向量；λ为反映***发电机有功、无功功率注入和负荷变化的向量；p_m为风电场机械功率注入向量（即***可变参数向量）。

先固定p_m不变，λ为自由变化参数，再通过连续潮流法描绘出一条有功功率和电压的P-U曲线，并且得出一个鞍节分叉点（Saddle NodeBifurcation，以下简称SNB点），即为初始边界点。以下为连续潮流计算方法。

2.1）预测步的计算

首先要确定变量的预测方向。以切线方向为例，对式（3）求微分，得

$f_{V_i}{d}_{V_i}+f_{{\mathrm{\δ}}_i}{d}_{{\mathrm{\δ}}_i}+f_{V_r}{d}_{V_r}+f_{{\mathrm{\θ}}_r}{d}_{{\mathrm{\θ}}_r}+f_{I_s}{d}_{I_s}+f_{{\mathrm{\φ}}_s}{d}_{{\mathrm{\φ}}_s}$

                  （4）

$+f_{I_r}{d}_{I_r}+f_{{\mathrm{\φ}}_r}{d}_{{\mathrm{\φ}}_r}+f_s{d}_s+f_s{d}_s+f_{\mathrm{\λ}}{d}_{\mathrm{\λ}}=0$

如果当前解处在解曲线的平凡位置，即雅可比矩阵非奇异，则以λ作为参数化变量，在临界点之前λ的变化方向为+1，在临界点之后λ的变化方向为-1。可得下式

$\left[\begin{array}{cccccccccc}{f}_{{1V}_i}& f_{{1\mathrm{\δ}}_i}& f_{{1V}_r}& f_{{1\mathrm{\θ}}_r}& f_{{1I}_S}& f_{{1\mathrm{\φ}}_S}& f_{{1I}_r}& f_{{1\mathrm{\φ}}_r}& f_{1s}& f_{1\mathrm{\λ}}\\ f_{{2V}_i}& f_{{2\mathrm{\δ}}_i}& f_{{2V}_r}& f_{{2\mathrm{\θ}}_r}& f_{{2I}_S}& f_{{2\mathrm{\φ}}_S}& f_{{2I}_r}& f_{{2\mathrm{\φ}}_r}& f_{2s}& f_{2\mathrm{\λ}}\\ f_{{3V}_i}& f_{{3\mathrm{\δ}}_i}& f_{{3V}_r}& f_{{3\mathrm{\θ}}_r}& f_{{3I}_S}& f_{{3\mathrm{\φ}}_S}& f_{{3I}_r}& f_{{3\mathrm{\φ}}_r}& f_{3s}& f_{3\mathrm{\λ}}\\ f_{{4V}_i}& f_{{4\mathrm{\δ}}_i}& f_{{4V}_r}& f_{{4\mathrm{\θ}}_r}& f_{{4I}_S}& f_{{4\mathrm{\φ}}_S}& f_{{4I}_r}& f_{{4\mathrm{\φ}}_r}& f_{4s}& f_{4\mathrm{\λ}}\\ f_{{5V}_i}& f_{{5\mathrm{\δ}}_i}& f_{{5V}_r}& f_{{5\mathrm{\θ}}_r}& f_{{5I}_S}& f_{{5\mathrm{\φ}}_S}& f_{{5I}_r}& f_{{5\mathrm{\φ}}_r}& f_{5s}& f_{5\mathrm{\λ}}\\ f_{{6V}_i}& f_{{6\mathrm{\δ}}_i}& f_{{6V}_r}& f_{{6\mathrm{\θ}}_r}& f_{{6I}_S}& f_{{6\mathrm{\φ}}_S}& f_{{6I}_r}& f_{{6\mathrm{\φ}}_r}& f_{6s}& f_{6\mathrm{\λ}}\\ f_{{7V}_i}& f_{{7\mathrm{\δ}}_i}& f_{{7V}_r}& f_{{7\mathrm{\θ}}_r}& f_{{7I}_S}& f_{{7\mathrm{\φ}}_S}& f_{{7I}_r}& f_{{7\mathrm{\φ}}_r}& f_{7s}& f_{7\mathrm{\λ}}\\ f_{{8V}_i}& f_{{8\mathrm{\δ}}_i}& f_{{8V}_r}& f_{{8\mathrm{\θ}}_r}& f_{{8I}_S}& f_{{8\mathrm{\φ}}_S}& f_{{8I}_r}& f_{{8\mathrm{\φ}}_r}& f_{8s}& f_{8\mathrm{\λ}}\\ f_{{9V}_i}& f_{{9\mathrm{\δ}}_i}& f_{{9V}_r}& f_{{9\mathrm{\θ}}_r}& f_{{9I}_S}& f_{9{\mathrm{\φ}}_S}& f_{{9I}_r}& f_{{9\mathrm{\φ}}_r}& f_{9s}& f_{9\mathrm{\λ}}\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 1\end{array}\right]\left[\begin{array}{c}{d}_{V_i}\\ d_{{\mathrm{\δ}}_i}\\ d_{V_r}\\ d_{{\mathrm{\θ}}_r}\\ d_{I_s}\\ d_{{\mathrm{\φ}}_s}\\ d_{I_r}\\ d_{{\mathrm{\φ}}_r}\\ d_s\\ d_{\mathrm{\λ}}\end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ \±1\end{array}\right]---\left(5\right)$

由此可确定预测方向。

如果当前解处在解曲线的位置接近临界点，即雅可比矩阵接近奇异，式（5）的系数矩阵也将接近奇异，则应选择状态变量（例如变化率最大的节点k的电压）作为参数化变化量，而将λ作为普通变化量，状态变量的斜率作为切线方向，此时由下式可确定预测方向：

$\left[\begin{array}{cccccccccc}{f}_{{1V}_i}& f_{{1\mathrm{\δ}}_i}& f_{{1V}_r}& f_{{1\mathrm{\θ}}_r}& f_{{1I}_S}& f_{{1\mathrm{\φ}}_S}& f_{{1I}_r}& f_{{1\mathrm{\φ}}_r}& f_{1s}& f_{1\mathrm{\λ}}\\ f_{{2V}_i}& f_{{2\mathrm{\δ}}_i}& f_{{2V}_r}& f_{{2\mathrm{\θ}}_r}& f_{{2I}_S}& f_{{2\mathrm{\φ}}_S}& f_{{2I}_r}& f_{{2\mathrm{\φ}}_r}& f_{2s}& f_{2\mathrm{\λ}}\\ f_{{3V}_i}& f_{{3\mathrm{\δ}}_i}& f_{{3V}_r}& f_{{3\mathrm{\θ}}_r}& f_{{3I}_S}& f_{{3\mathrm{\φ}}_S}& f_{{3I}_r}& f_{{3\mathrm{\φ}}_r}& f_{3s}& f_{3\mathrm{\λ}}\\ f_{{4V}_i}& f_{{4\mathrm{\δ}}_i}& f_{{4V}_r}& f_{{4\mathrm{\θ}}_r}& f_{{4I}_S}& f_{{4\mathrm{\φ}}_S}& f_{{4I}_r}& f_{{4\mathrm{\φ}}_r}& f_{4s}& f_{4\mathrm{\λ}}\\ f_{{5V}_i}& f_{{5\mathrm{\δ}}_i}& f_{{5V}_r}& f_{{5\mathrm{\θ}}_r}& f_{{5I}_S}& f_{{5\mathrm{\φ}}_S}& f_{{5I}_r}& f_{{5\mathrm{\φ}}_r}& f_{5s}& f_{5\mathrm{\λ}}\\ f_{{6V}_i}& f_{{6\mathrm{\δ}}_i}& f_{{6V}_r}& f_{{6\mathrm{\θ}}_r}& f_{{6I}_S}& f_{{6\mathrm{\φ}}_S}& f_{{6I}_r}& f_{{6\mathrm{\φ}}_r}& f_{6s}& f_{6\mathrm{\λ}}\\ f_{{7V}_i}& f_{{7\mathrm{\δ}}_i}& f_{{7V}_r}& f_{{7\mathrm{\θ}}_r}& f_{{7I}_S}& f_{{7\mathrm{\φ}}_S}& f_{{7I}_r}& f_{{7\mathrm{\φ}}_r}& f_{7s}& f_{7\mathrm{\λ}}\\ f_{{8V}_i}& f_{{8\mathrm{\δ}}_i}& f_{{8V}_r}& f_{{8\mathrm{\θ}}_r}& f_{{8I}_S}& f_{{8\mathrm{\φ}}_S}& f_{{8I}_r}& f_{{8\mathrm{\φ}}_r}& f_{8s}& f_{8\mathrm{\λ}}\\ f_{{9V}_i}& f_{{9\mathrm{\δ}}_i}& f_{{9V}_r}& f_{{9\mathrm{\θ}}_r}& f_{{9I}_S}& f_{9{\mathrm{\φ}}_S}& f_{{9I}_r}& f_{{9\mathrm{\φ}}_r}& f_{9s}& f_{9\mathrm{\λ}}\\ & & & & e_k^t& & & & & 0\end{array}\right]\left[\begin{array}{c}{d}_{V_i}\\ d_{{\mathrm{\δ}}_i}\\ d_{V_r}\\ d_{{\mathrm{\θ}}_r}\\ d_{I_s}\\ d_{{\mathrm{\φ}}_s}\\ d_{I_r}\\ d_{{\mathrm{\φ}}_r}\\ d_s\\ d_{\mathrm{\λ}}\end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ \±1\end{array}\right]---\left(6\right)$

其中，

为第k个元素为+1，其余元素为0的行矢量。应注意式（6）与式（5）的系数矩阵不同。由于选择了状态变量作为参数化变量，即使雅可比矩阵接近奇异（假定其秩是n-1），也可以证明是（6）的系数矩阵是非奇异的。

根据式（5）或式（6）确定的预测方向，可以计算预测点如下：

$\left[\begin{array}{c}{\tilde{V}}_i\\ {\widetilde{\mathrm{\δ}}}_i\\ {\tilde{V}}_r\\ {\widetilde{\mathrm{\θ}}}_r\\ {\tilde{I}}_S\\ {\widetilde{\mathrm{\φ}}}_S\\ {\tilde{I}}_r\\ {\widetilde{\mathrm{\φ}}}_r\\ \tilde{s}\\ \widetilde{\mathrm{\λ}}\end{array}\right]=\left[\begin{array}{c}{V}_{i0}\\ {\mathrm{\δ}}_{i0}\\ V_{r0}\\ {\mathrm{\θ}}_{r0}\\ I_{S0}\\ {\mathrm{\φ}}_{S0}\\ I_{r0}\\ {\mathrm{\φ}}_{r0}\\ s_0\\ {\mathrm{\λ}}_0\end{array}\right]+\mathrm{\σ}\left[\begin{array}{c}{d}_{V_i}\\ d_{{\mathrm{\δ}}_i}\\ d_{V_r}\\ d_{{\mathrm{\θ}}_r}\\ d_{I_S}\\ d_{{\mathrm{\φ}}_S}\\ d_{I_r}\\ d_{{\mathrm{\φ}}_r}\\ d_s\\ d_{\mathrm{\λ}}\end{array}\right]---\left(7\right)$

其中，σ为预测步。

2.2）校正步的计算

在校正步中，如果预测方向是由式（5）得到的，应先固定λ，采用垂直校正方法，以为初值求解式（1）的潮流方程。以牛顿-拉夫逊法为例，迭代格式如下：

$\left[\begin{array}{cccccccccc}{f}_{{1V}_i}& f_{{1\mathrm{\δ}}_i}& f_{{1V}_r}& f_{{1\mathrm{\θ}}_r}& f_{{1I}_S}& f_{{1\mathrm{\φ}}_S}& f_{{1I}_r}& f_{{1\mathrm{\φ}}_r}& f_{1s}& f_{1\mathrm{\λ}}\\ f_{{2V}_i}& f_{{2\mathrm{\δ}}_i}& f_{{2V}_r}& f_{{2\mathrm{\θ}}_r}& f_{{2I}_S}& f_{{2\mathrm{\φ}}_S}& f_{{2I}_r}& f_{{2\mathrm{\φ}}_r}& f_{2s}& f_{2\mathrm{\λ}}\\ f_{{3V}_i}& f_{{3\mathrm{\δ}}_i}& f_{{3V}_r}& f_{{3\mathrm{\θ}}_r}& f_{{3I}_S}& f_{{3\mathrm{\φ}}_S}& f_{{3I}_r}& f_{{3\mathrm{\φ}}_r}& f_{3s}& f_{3\mathrm{\λ}}\\ f_{{4V}_i}& f_{{4\mathrm{\δ}}_i}& f_{{4V}_r}& f_{{4\mathrm{\θ}}_r}& f_{{4I}_S}& f_{{4\mathrm{\φ}}_S}& f_{{4I}_r}& f_{{4\mathrm{\φ}}_r}& f_{4s}& f_{4\mathrm{\λ}}\\ f_{{5V}_i}& f_{{5\mathrm{\δ}}_i}& f_{{5V}_r}& f_{{5\mathrm{\θ}}_r}& f_{{5I}_S}& f_{{5\mathrm{\φ}}_S}& f_{{5I}_r}& f_{{5\mathrm{\φ}}_r}& f_{5s}& f_{5\mathrm{\λ}}\\ f_{{6V}_i}& f_{{6\mathrm{\δ}}_i}& f_{{6V}_r}& f_{{6\mathrm{\θ}}_r}& f_{{6I}_S}& f_{{6\mathrm{\φ}}_S}& f_{{6I}_r}& f_{{6\mathrm{\φ}}_r}& f_{6s}& f_{6\mathrm{\λ}}\\ f_{{7V}_i}& f_{{7\mathrm{\δ}}_i}& f_{{7V}_r}& f_{{7\mathrm{\θ}}_r}& f_{{7I}_S}& f_{{7\mathrm{\φ}}_S}& f_{{7I}_r}& f_{{7\mathrm{\φ}}_r}& f_{7s}& f_{7\mathrm{\λ}}\\ f_{{8V}_i}& f_{{8\mathrm{\δ}}_i}& f_{{8V}_r}& f_{{8\mathrm{\θ}}_r}& f_{{8I}_S}& f_{{8\mathrm{\φ}}_S}& f_{{8I}_r}& f_{{8\mathrm{\φ}}_r}& f_{8s}& f_{8\mathrm{\λ}}\\ f_{{9V}_i}& f_{{9\mathrm{\δ}}_i}& f_{{9V}_r}& f_{{9\mathrm{\θ}}_r}& f_{{9I}_S}& f_{9{\mathrm{\φ}}_S}& f_{{9I}_r}& f_{{9\mathrm{\φ}}_r}& f_{9s}& f_{9\mathrm{\λ}}\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 1\end{array}\right]\left[\begin{array}{c}{\mathrm{\Δ}}_{V_i}\\ {\mathrm{\Δ}}_{{\mathrm{\δ}}_i}\\ {\mathrm{\Δ}}_{V_r}\\ {\mathrm{\Δ}}_{{\mathrm{\θ}}_r}\\ {\mathrm{\Δ}}_{I_s}\\ {\mathrm{\Δ}}_{{\mathrm{\φ}}_s}\\ {\mathrm{\Δ}}_{I_r}\\ {\mathrm{\Δ}}_{{\mathrm{\φ}}_r}\\ {\mathrm{\Δ}}_s\\ {\mathrm{\Δ}}_{\mathrm{\λ}}\end{array}\right]=-\left[\begin{array}{c}{f}_1\\ f_2\\ f_3\\ f_4\\ f_5\\ f_6\\ f_7\\ f_8\\ f_9\\ 0\end{array}\right]---\left(8\right)$

如果上述潮流计算收敛，则可以得到P-U曲线上的一点，然后开始新的预测步计算。如果潮流发散，有两种对应的措施：其一是减小步长σ用式（7）预测新的点，并用垂直校正法重新用式（8）迭代；另一种办法是选择状态变量（例如变化率最大的节点k的电压）作为参数化变量，将λ作为普通变量，采用水平校正方法解潮流方程，迭代格式如下：

$\left[\begin{array}{cccccccccc}{f}_{{1V}_i}& f_{{1\mathrm{\δ}}_i}& f_{{1V}_r}& f_{{1\mathrm{\θ}}_r}& f_{{1I}_S}& f_{{1\mathrm{\φ}}_S}& f_{{1I}_r}& f_{{1\mathrm{\φ}}_r}& f_{1s}& f_{1\mathrm{\λ}}\\ f_{{2V}_i}& f_{{2\mathrm{\δ}}_i}& f_{{2V}_r}& f_{{2\mathrm{\θ}}_r}& f_{{2I}_S}& f_{{2\mathrm{\φ}}_S}& f_{{2I}_r}& f_{{2\mathrm{\φ}}_r}& f_{2s}& f_{2\mathrm{\λ}}\\ f_{{3V}_i}& f_{{3\mathrm{\δ}}_i}& f_{{3V}_r}& f_{{3\mathrm{\θ}}_r}& f_{{3I}_S}& f_{{3\mathrm{\φ}}_S}& f_{{3I}_r}& f_{{3\mathrm{\φ}}_r}& f_{3s}& f_{3\mathrm{\λ}}\\ f_{{4V}_i}& f_{{4\mathrm{\δ}}_i}& f_{{4V}_r}& f_{{4\mathrm{\θ}}_r}& f_{{4I}_S}& f_{{4\mathrm{\φ}}_S}& f_{{4I}_r}& f_{{4\mathrm{\φ}}_r}& f_{4s}& f_{4\mathrm{\λ}}\\ f_{{5V}_i}& f_{{5\mathrm{\δ}}_i}& f_{{5V}_r}& f_{{5\mathrm{\θ}}_r}& f_{{5I}_S}& f_{{5\mathrm{\φ}}_S}& f_{{5I}_r}& f_{{5\mathrm{\φ}}_r}& f_{5s}& f_{5\mathrm{\λ}}\\ f_{{6V}_i}& f_{{6\mathrm{\δ}}_i}& f_{{6V}_r}& f_{{6\mathrm{\θ}}_r}& f_{{6I}_S}& f_{{6\mathrm{\φ}}_S}& f_{{6I}_r}& f_{{6\mathrm{\φ}}_r}& f_{6s}& f_{6\mathrm{\λ}}\\ f_{{7V}_i}& f_{{7\mathrm{\δ}}_i}& f_{{7V}_r}& f_{{7\mathrm{\θ}}_r}& f_{{7I}_S}& f_{{7\mathrm{\φ}}_S}& f_{{7I}_r}& f_{{7\mathrm{\φ}}_r}& f_{7s}& f_{7\mathrm{\λ}}\\ f_{{8V}_i}& f_{{8\mathrm{\δ}}_i}& f_{{8V}_r}& f_{{8\mathrm{\θ}}_r}& f_{{8I}_S}& f_{{8\mathrm{\φ}}_S}& f_{{8I}_r}& f_{{8\mathrm{\φ}}_r}& f_{8s}& f_{8\mathrm{\λ}}\\ f_{{9V}_i}& f_{{9\mathrm{\δ}}_i}& f_{{9V}_r}& f_{{9\mathrm{\θ}}_r}& f_{{9I}_S}& f_{9{\mathrm{\φ}}_S}& f_{{9I}_r}& f_{{9\mathrm{\φ}}_r}& f_{9s}& f_{9\mathrm{\λ}}\\ & & & & e_k^t& & & & & 0\end{array}\right]\left[\begin{array}{c}{\mathrm{\Δ}}_{V_i}\\ {\mathrm{\Δ}}_{{\mathrm{\δ}}_i}\\ {\mathrm{\Δ}}_{V_r}\\ {\mathrm{\Δ}}_{{\mathrm{\θ}}_r}\\ {\mathrm{\Δ}}_{I_s}\\ {\mathrm{\Δ}}_{{\mathrm{\φ}}_s}\\ {\mathrm{\Δ}}_{I_r}\\ {\mathrm{\Δ}}_{{\mathrm{\φ}}_r}\\ {\mathrm{\Δ}}_s\\ {\mathrm{\Δ}}_{\mathrm{\λ}}\end{array}\right]=-\left[\begin{array}{c}{f}_1\\ f_2\\ f_3\\ f_4\\ f_5\\ f_6\\ f_7\\ f_8\\ f_9\\ 0\end{array}\right]---\left(9\right)$

在校正步中，如果预测方向本身就是由式（6）得到的，也应按式（9）解潮流方程。

3）计算SNB点精确值

精确求解SNB初始点，SNB点的约束方程式可描述为

由SNB点构成的电压稳定域边界面可描述为

$\left\{\begin{array}{c}f(V_i,{\mathrm{\δ}}_i,V_r,{\mathrm{\θ}}_r,I_s,{\mathrm{\φ}}_s,I_r,{\mathrm{\φ}}_r,s,\mathrm{\λ},P_m)=0\\ f_{V_i}{|}_{*}y+f_{{\mathrm{\δ}}_i}{|}_{*}y+f_{V_r}{|}_{*}y+f_{{\mathrm{\θ}}_r}{|}_{*}y+f_{I_S}{|}_{*}y+f_{{\mathrm{\φ}}_S}{|}_{*}y+f_{I_r}{|}_{*}y+f_{{\mathrm{\φ}}_r}{|}_{*}y+f_s{|}_{*}y=0\\ y^{t}y=1\end{array}\right.---\left(10\right)$

y^t=[dV_i dδ_i dV_r dθ_r dI_s dφ_s dI_r dφ_r ds]

式中，V_i,δ_i,V_r,θ_r,I_S,φ_S,I_r,φ_r,s为***的n维状态向量；λ为反映***发电机有功、无功功率注入和负荷变化的向量；P_m为风电场机械功率注入向量（即***可变参数向量），此时先将P_m值固定；y为潮流方程雅可比矩阵的右特征向量。以连续潮流得到的解为初值，用牛顿—拉夫逊法进行求解方程变量。精确得出初始点。

4)计算含有DFIG电力***稳定域的SNB点轨迹

该稳定域SNB点轨迹的算法同样是应用预测-校正方法，将式（10）可简记为

φ(z)=0     （11）

此时

z^t=[V_i δ_i V_r θ_r I_S φ_S I_r φ_r s λ P_m]∈R²ⁿ⁺²,

$\mathrm{\φ}=\left[\begin{array}{c}f(V_i,{\mathrm{\δ}}_i,V_r,{\mathrm{\θ}}_r,I_s,{\mathrm{\φ}}_s,I_r,{\mathrm{\φ}}_r,s,\mathrm{\λ},P_m)\\ f_{V_i}{|}_{*}y+f_{{\mathrm{\δ}}_i}{|}_{*}y+f_{V_r}{|}_{*}y+f_{{\mathrm{\θ}}_r}{|}_{*}y+f_{I_S}{|}_{*}y+f_{{\mathrm{\φ}}_S}{|}_{*}y+f_{I_r}{|}_{*}y+f_{{\mathrm{\φ}}_r}{|}_{*}y+f_s{|}_{*}y\\ y^{t}y-1\end{array}\right]:R^{2n+2}\→R^{2n+1}$

由上式所确定的雅可比矩阵的维数为(2n+1)×(2n+2)，因此加入一个参数化方程

(z-z₁)^tv=τ     （12）

式中，v为z₁点的切向量，z₀为初始点，τ为步长。v是由下式确定的：

${\mathrm{\φ}}_z{|}_{z=z_0}v=0---\left(13\right)$

4.1）预测步的计算

首先假定由步骤2）已经计算得到一点z₀在该电力***稳定域SNB点轨迹上。从该点开始进行预测，由下式可确定预测方向：

$\left[\begin{array}{cccccccccccc}{f}_{{1V}_I}& f_{{1\mathrm{\δ}}_i}& f_{{1V}_r}& f_{{1\mathrm{\θ}}_r}& f_{{1I}_S}& f_{{1\mathrm{\φ}}_S}& f_{{1I}_r}& f_{{1\mathrm{\φ}}_r}& f_{1s}& 0& f_{1\mathrm{\λ}}& f_{{1P}_m}\\ f_{{2V}_i}& f_{{2\mathrm{\δ}}_i}& f_{{2V}_r}& f_{{2\mathrm{\θ}}_r}& f_{{2I}_S}& f_{{2\mathrm{\φ}}_S}& f_{{2I}_r}& f_{{2\mathrm{\φ}}_r}& f_{2s}& 0& f_{2\mathrm{\λ}}& f_{{2P}_m}\\ f_{{3V}_i}& f_{{3\mathrm{\δ}}_i}& f_{{3V}_r}& f_{{3\mathrm{\θ}}_r}& f_{{3I}_S}& f_{{3\mathrm{\φ}}_S}& f_{{3I}_r}& f_{{3\mathrm{\φ}}_r}& f_{3s}& 0& f_{3\mathrm{\λ}}& f_{{3P}_m}\\ f_{{4V}_i}& f_{{4\mathrm{\δ}}_i}& f_{{4V}_r}& f_{{4\mathrm{\θ}}_r}& f_{{4I}_S}& f_{{4\mathrm{\φ}}_S}& f_{{4I}_r}& f_{{4\mathrm{\φ}}_r}& f_{4s}& 0& f_{4\mathrm{\λ}}& f_{{4P}_m}\\ f_{{5V}_i}& f_{{5\mathrm{\δ}}_i}& f_{{5V}_r}& f_{{5\mathrm{\θ}}_r}& f_{{5I}_S}& f_{{5\mathrm{\φ}}_S}& f_{{5I}_r}& f_{{5\mathrm{\φ}}_r}& f_{5s}& 0& f_{5\mathrm{\λ}}& f_{{5P}_m}\\ f_{{6V}_i}& f_{{6\mathrm{\δ}}_i}& f_{{6V}_r}& f_{6{\mathrm{\θ}}_r}& f_{{6I}_S}& f_{{6\mathrm{\φ}}_S}& f_{{6I}_r}& f_{{6\mathrm{\φ}}_r}& f_{6s}& 0& f_{6\mathrm{\λ}}& f_{{6P}_m}\\ f_{{7V}_i}& f_{{7\mathrm{\δ}}_i}& f_{{7V}_r}& f_{{7\mathrm{\θ}}_r}& f_{{7I}_S}& f_{{7\mathrm{\φ}}_S}& f_{{7I}_r}& f_{{7\mathrm{\φ}}_r}& f_{7s}& 0& f_{7\mathrm{\λ}}& f_{{7P}_m}\\ f_{{8V}_i}& f_{{8\mathrm{\δ}}_i}& f_{{8V}_r}& f_{{8\mathrm{\θ}}_r}& f_{{8I}_S}& f_{{8\mathrm{\φ}}_S}& f_{{8I}_r}& f_{{8\mathrm{\φ}}_r}& f_{8s}& 0& f_{8\mathrm{\λ}}& f_{{8P}_m}\\ f_{{9V}_i}& f_{{9\mathrm{\δ}}_i}& f_{{9V}_r}& f_{{9\mathrm{\θ}}_r}& f_{{9I}_S}& f_{{9\mathrm{\φ}}_S}& f_{{9I}_r}& f_{{9\mathrm{\φ}}_r}& f_{9s}& 0& f_{9\mathrm{\λ}}& f_{{9P}_m}\\ \frac{\∂g}{{\∂V}_i}& \frac{\∂g}{{\∂\mathrm{\δ}}_i}& \frac{\∂g}{{\∂V}_r}& \frac{\∂g}{{\∂\mathrm{\θ}}_r}& \frac{\∂g}{{\∂I}_S}& \frac{\∂g}{{\∂\mathrm{\φ}}_S}& \frac{\∂g}{{\∂I}_r}& \frac{\∂g}{{\∂\mathrm{\φ}}_r}& \frac{\∂g}{\∂s}& \frac{\∂g}{\∂y}& \frac{\∂g}{\∂\mathrm{\λ}}& \frac{\∂g}{{\∂P}_m}\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& {2y}^t& 0& 0\\ {\mathrm{\ΔV}}_i& {\mathrm{\Δ\δ}}_i& {\mathrm{\ΔV}}_r& {\mathrm{\Δ\θ}}_r& {\mathrm{\ΔI}}_S& {\mathrm{\Δ\φ}}_S& {\mathrm{\ΔI}}_r& {\mathrm{\Δ\φ}}_r& \mathrm{\Δs}& {\mathrm{\Δy}}^t& \mathrm{\Δ\λ}& {\mathrm{\ΔP}}_m\end{array}\right]\left[\begin{array}{c}{d}_{V_i}\\ d_{{\mathrm{\δ}}_i}\\ d_{V_r}\\ d_{{\mathrm{\θ}}_r}\\ d_{I_S}\\ d_{{\mathrm{\φ}}_S}\\ d_{I_r}\\ d_{{\mathrm{\φ}}_r}\\ d_s\\ d_y\\ d_{\mathrm{\λ}}\\ d_{P_m}\end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\\ \±1\end{array}\right]---\left(14\right)$

根据（14）所确定的预测方向，可以计算预测点如下：

$\left[\begin{array}{c}{\tilde{V}}_i\\ {\widetilde{\mathrm{\δ}}}_i\\ {\tilde{V}}_r\\ {\widetilde{\mathrm{\θ}}}_r\\ {\tilde{I}}_S\\ {\mathrm{\φ}}_S\\ {\tilde{I}}_r\\ {\widetilde{\mathrm{\φ}}}_r\\ \tilde{s}\\ \tilde{y}\\ \widetilde{\mathrm{\λ}}\\ {\tilde{P}}_m\end{array}\right]=\left[\begin{array}{c}{V}_{i0}\\ {\mathrm{\δ}}_{i0}\\ V_{r0}\\ {\mathrm{\θ}}_{r0}\\ I_{S0}\\ {\mathrm{\φ}}_{S0}\\ I_{r0}\\ {\mathrm{\φ}}_{r0}\\ s_0\\ y_0\\ {\mathrm{\λ}}_0\\ P_{m0}\end{array}\right]+\mathrm{\τ}\left[\begin{array}{c}{d}_{V_i}\\ d_{{\mathrm{\δ}}_i}\\ d_{V_r}\\ d_{{\mathrm{\θ}}_r}\\ d_{I_S}\\ d_{{\mathrm{\φ}}_S}\\ d_{I_r}\\ d_{{\mathrm{\φ}}_r}\\ d_s\\ d_y\\ d_{\mathrm{\λ}}\\ d_{P_m}\end{array}\right]---\left(15\right)$

式中，V_i,δ_i,V_r,θ_r,I_S,φ_S,I_r,φ_r,s为***的n维状态向量；y为潮流方程雅可比矩阵的右特征向量；λ为反映***发电机有功、无功功率注入和负荷变化的向量；p_m为风电场机械功率注入向量（即***可变参数向量）；τ为预测步长。

4.2）校正步的计算

在校正格式中，增加了一维校正方程与式（10）联立得校正方程

$\left[\begin{array}{cccccccccccc}{f}_{{1V}_I}& f_{{1\mathrm{\δ}}_i}& f_{{1V}_r}& f_{{1\mathrm{\θ}}_r}& f_{{1I}_S}& f_{{1\mathrm{\φ}}_S}& f_{{1I}_r}& f_{{1\mathrm{\φ}}_r}& f_{1s}& 0& f_{1\mathrm{\λ}}& f_{{1P}_m}\\ f_{{2V}_i}& f_{{2\mathrm{\δ}}_i}& f_{{2V}_r}& f_{{2\mathrm{\θ}}_r}& f_{{2I}_S}& f_{{2\mathrm{\φ}}_S}& f_{{2I}_r}& f_{{2\mathrm{\φ}}_r}& f_{2s}& 0& f_{2\mathrm{\λ}}& f_{{2P}_m}\\ f_{{3V}_i}& f_{{3\mathrm{\δ}}_i}& f_{{3V}_r}& f_{{3\mathrm{\θ}}_r}& f_{{3I}_S}& f_{{3\mathrm{\φ}}_S}& f_{{3I}_r}& f_{{3\mathrm{\φ}}_r}& f_{3s}& 0& f_{3\mathrm{\λ}}& f_{{3P}_m}\\ f_{{4V}_i}& f_{{4\mathrm{\δ}}_i}& f_{{4V}_r}& f_{{4\mathrm{\θ}}_r}& f_{{4I}_S}& f_{{4\mathrm{\φ}}_S}& f_{{4I}_r}& f_{{4\mathrm{\φ}}_r}& f_{4s}& 0& f_{4\mathrm{\λ}}& f_{{4P}_m}\\ f_{{5V}_i}& f_{{5\mathrm{\δ}}_i}& f_{{5V}_r}& f_{{5\mathrm{\θ}}_r}& f_{{5I}_S}& f_{{5\mathrm{\φ}}_S}& f_{{5I}_r}& f_{{5\mathrm{\φ}}_r}& f_{5s}& 0& f_{5\mathrm{\λ}}& f_{{5P}_m}\\ f_{{6V}_i}& f_{{6\mathrm{\δ}}_i}& f_{{6V}_r}& f_{6{\mathrm{\θ}}_r}& f_{{6I}_S}& f_{{6\mathrm{\φ}}_S}& f_{{6I}_r}& f_{{6\mathrm{\φ}}_r}& f_{6s}& 0& f_{6\mathrm{\λ}}& f_{{6P}_m}\\ f_{{7V}_i}& f_{{7\mathrm{\δ}}_i}& f_{{7V}_r}& f_{{7\mathrm{\θ}}_r}& f_{{7I}_S}& f_{{7\mathrm{\φ}}_S}& f_{{7I}_r}& f_{{7\mathrm{\φ}}_r}& f_{7s}& 0& f_{7\mathrm{\λ}}& f_{{7P}_m}\\ f_{{8V}_i}& f_{{8\mathrm{\δ}}_i}& f_{{8V}_r}& f_{{8\mathrm{\θ}}_r}& f_{{8I}_S}& f_{{8\mathrm{\φ}}_S}& f_{{8I}_r}& f_{{8\mathrm{\φ}}_r}& f_{8s}& 0& f_{8\mathrm{\λ}}& f_{{8P}_m}\\ f_{{9V}_i}& f_{{9\mathrm{\δ}}_i}& f_{{9V}_r}& f_{{9\mathrm{\θ}}_r}& f_{{9I}_S}& f_{{9\mathrm{\φ}}_S}& f_{{9I}_r}& f_{{9\mathrm{\φ}}_r}& f_{9s}& 0& f_{9\mathrm{\λ}}& f_{{9P}_m}\\ \frac{\∂g}{{\∂V}_i}& \frac{\∂g}{{\∂\mathrm{\δ}}_i}& \frac{\∂g}{{\∂V}_r}& \frac{\∂g}{{\∂\mathrm{\θ}}_r}& \frac{\∂g}{{\∂I}_S}& \frac{\∂g}{{\∂\mathrm{\φ}}_S}& \frac{\∂g}{{\∂I}_r}& \frac{\∂g}{{\∂\mathrm{\φ}}_r}& \frac{\∂g}{\∂s}& \frac{\∂g}{\∂y}& \frac{\∂g}{\∂\mathrm{\λ}}& \frac{\∂g}{{\∂P}_m}\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& {2y}^t& 0& 0\\ {\mathrm{\ΔV}}_i& {\mathrm{\Δ\δ}}_i& {\mathrm{\ΔV}}_r& {\mathrm{\Δ\θ}}_r& {\mathrm{\ΔI}}_S& {\mathrm{\Δ\φ}}_S& {\mathrm{\ΔI}}_r& {\mathrm{\Δ\φ}}_r& \mathrm{\Δs}& {\mathrm{\Δy}}^t& \mathrm{\Δ\λ}& {\mathrm{\ΔP}}_m\end{array}\right]\left[\begin{array}{c}{\mathrm{\Δ}}_{V_i}\\ {\mathrm{\Δ}}_{{\mathrm{\δ}}_i}\\ {\mathrm{\Δ}}_{V_r}\\ {\mathrm{\Δ}}_{{\mathrm{\θ}}_r}\\ {\mathrm{\Δ}}_{I_S}\\ {\mathrm{\Δ}}_{{\mathrm{\φ}}_S}\\ {\mathrm{\Δ}}_{I_r}\\ {\mathrm{\Δ}}_{{\mathrm{\φ}}_r}\\ {\mathrm{\Δ}}_s\\ {\mathrm{\Δ}}_y\\ {\mathrm{\Δ}}_{\mathrm{\λ}}\\ {\mathrm{\Δ}}_{P_m}\end{array}\right]=-\left[\begin{array}{c}{f}_1\\ f_2\\ f_3\\ f_4\\ f_5\\ f_6\\ f_7\\ f_8\\ f_9\\ g\\ y^{t}y-4\\ 0\end{array}\right]---\left(16\right)$

修正后的变量为

$\left[\begin{array}{c}{V}_i^{(n+1)}\\ {\mathrm{\δ}}_i^{(n+1)}\\ V_r^{(n+1)}\\ {\mathrm{\θ}}_r^{(n+1)}\\ I_S^{(n+1)}\\ {\mathrm{\φ}}_S^{(n+1)}\\ I_r^{(n+1)}\\ {\mathrm{\φ}}_r^{(n+1)}\\ s^{(n+1)}\\ y^{(n+1)}\\ {\mathrm{\λ}}^{(n+1)}\\ P_m^{(n+1)}\end{array}\right]=\left[\begin{array}{c}{V}_i^{\left(n\right)}\\ {\mathrm{\δ}}_i^{\left(n\right)}\\ V_r^{\left(n\right)}\\ {\mathrm{\θ}}_r^{\left(n\right)}\\ I_S^{\left(n\right)}\\ {\mathrm{\φ}}_s^{\left(n\right)}\\ I_r^{\left(n\right)}\\ {\mathrm{\φ}}_r^n\\ s^{\left(n\right)}\\ y^{\left(n\right)}\\ {\mathrm{\λ}}^{\left(n\right)}\\ P_m^{\left(n\right)}\end{array}\right]+\left[\begin{array}{c}{\mathrm{\ΔV}}_i^{\left(n\right)}\\ \mathrm{\Δ}{\mathrm{\δ}}_i^{\left(n\right)}\\ \mathrm{\Δ}{V}_r^{\left(n\right)}\\ \mathrm{\Δ}{\mathrm{\θ}}_r^{\left(n\right)}\\ \mathrm{\Δ}{I}_S^{\left(n\right)}\\ \mathrm{\Δ}{\mathrm{\φ}}_S^{\left(n\right)}\\ \mathrm{\Δ}{I}_r^{\left(n\right)}\\ \mathrm{\Δ}{\mathrm{\φ}}_r^{\left(n\right)}\\ \mathrm{\Δ}{s}^{\left(n\right)}\\ \mathrm{\Δ}{y}^{\left(n\right)}\\ \mathrm{\Δ}{\mathrm{\λ}}^{\left(n\right)}\\ \mathrm{\Δ}{P}_m^{\left(n\right)}\end{array}\right]---\left(17\right)$

上述预测-校正方法中函数g的定义及相关计算如下：

$g=f_{V_i}{|}_{*}y+f_{{\mathrm{\δ}}_i}{|}_{*}y+f_{V_r}{|}_{*}y+f_{{\mathrm{\θ}}_r}{|}_{*}y$

           （18）

$+f_{I_S}{|}_{*}y+f_{{\mathrm{\φ}}_S}{|}_{*}y+f_{I_r}{|}_{*}y+f_{{\mathrm{\φ}}_r}{|}_{*}y+f_s{|}_{*}y=0$

由于上述雅可比矩阵含零元素较多，所以为稀疏矩阵。

在电压稳定域边界上的已知一个SNB点为初始点的情况下，利用该预测-校正法可得到一条含有DFIG电力***稳定域的SNB点轨迹。

5）计算含有DFIG电力***稳定域的切平面

在SNB点处将式（10）线性化得到

$f_{V_i}{|}_{*}\mathrm{\Δ}{V}_i+f_{{\mathrm{\δ}}_i}{|}_{*}{\mathrm{\Δ\δ}}_i+f_{V_r}{|}_{*}\mathrm{\Δ}{V}_r+f_{{\mathrm{\θ}}_r}{|}_{*}\mathrm{\Δ}{\mathrm{\θ}}_r+f_{I_S}{|}_{*}\mathrm{\Δ}{I}_S+f_{{\mathrm{\φ}}_S}{|}_{*}\mathrm{\Δ}{\mathrm{\φ}}_S$

                （19）

$+f_{I_r}{|}_{*}\mathrm{\Δ}{I}_r+f_{{\mathrm{\φ}}_r}{|}_{*}\mathrm{\Δ}{\mathrm{\φ}}_r+f_s{|}_{*}\mathrm{\Δs}+f_{\mathrm{\λ}}{|}_{*}\mathrm{\Δ\λ}+f_{P_m}{|}_{*}\mathrm{\Δ}{P}_m=0$

式中，为***在相应SNB点处的雅可比矩阵；f_λ|_*和

分别为式（10）在该点处对λ和P_m的导数矩阵。

对上式（19）左乘该SNB点处雅可比矩阵零特征值对应的左特征向量ω，得

$\mathrm{\ω}{f}_{V_i}{|}_{*}\mathrm{\Δ}{V}_i+\mathrm{\ω}{f}_{{\mathrm{\δ}}_i}{|}_{*}\mathrm{\Δ}{\mathrm{\δ}}_i+\mathrm{\ω}{f}_{V_r}{|}_{*}\mathrm{\Δ}{V}_r+\mathrm{\ω}{f}_{{\mathrm{\θ}}_r}{|}_{*}\mathrm{\Δ}{\mathrm{\θ}}_r+\mathrm{\ω}{f}_{I_s}{|}_{*}\mathrm{\Δ}{I}_s+\mathrm{\ω}{f}_{{\mathrm{\φ}}_s}{|}_{*}\mathrm{\Δ}{\mathrm{\φ}}_s$

               （20）

$+{\mathrm{\ωf}}_{I_r}{|}_{*}{\mathrm{\ΔI}}_r+{\mathrm{\ωf}}_{{\mathrm{\φ}}_r}{|}_{*}{\mathrm{\Δ\φ}}_r+{\mathrm{\ωf}}_s{|}_{*}\mathrm{\Δs}+{\mathrm{\ωf}}_{\mathrm{\λ}}{|}_{*}\mathrm{\Δ\λ}+{\mathrm{\ωf}}_{P_m}{|}_{*}{\mathrm{\ΔP}}_m=0$

式中，“|_*”表示对应系数的取值来自SNB点。

设L为负荷裕度，即λ|_*=λ₀+L。又由于式子中SNB处所以上式可化简为

${\mathrm{\ωf}}_{\mathrm{\λ}}{|}_{*}\mathrm{\ΔL}+{\mathrm{\ωf}}_{P_m}{|}_{*}{\mathrm{\ΔP}}_m=0---\left(21\right)$

其均为零，式（21）可视为L-P_m空间中电压稳定域的局部近似条件，若已知一个SNB点z₁，将该点的对应数代入式（21）中，可得到系数向量

该系数向量就是L-P_m空间中的电压稳定域边界在点z₁处的法向量。设z₁在该空间中的坐标为(L₁,P_m1)，则根据点法式方程可得到电压稳定域边界在z₁处的切线表达式为

${\mathrm{\ωf}}_{\mathrm{\λ}}{|}_{*}^1(L-L_1)+{\mathrm{\ωf}}_{P_m}{|}_{*}^1(P_m-P_{m1})=0---\left(22\right)$

将式（22）给出的L-P_m空间中电压稳定域边界切平面方程映射到发电机、负荷和DFIG风电场全注入空间，考虑如下关系：

$(1+k_{i}L){P_i}^0=P_i{|}_{*}---\left(23\right)$

式中k_i为***各个节点发电机和负荷有功注入的增长方向。将式（23）和式（22）联立，可得到在发电机、负荷和风电场全注入空间中静态电压稳定安全域边界在z₁点的切平面方程

$\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\ωf}}_{\mathrm{\λ}}{|}_{*}^1\frac{(P_i{|}_{*}-P_i{|}_{*}^1)}{k_i{P_i}^0}+{\mathrm{\ωf}}_{P_m}{|}_{*}^1(P_m{|}_{*}-P_m{|}_{*}^1)=0---\left(24\right)$

式（24）就是所求解的含DFIG风电场的电力***静态电压稳定安全域局部切平面解析式。

以上实施方式仅用于说明本发明，而并非对本发明的限制，有关技术领域的普通技术人员，在不脱离本发明的精神和范围的情况下，还可以做出各种变化和变型，因此所有等同的技术方案也属于本发明的保护范畴。

Claims (8)

1.一种含DFIG风电场的电力***电压稳定域切平面的获取方法，其特征在于包括如下步骤：

1）对DFIG风电场进行数学建模；

2）用连续潮流法计算含DFIG电力***静态电压稳定边界SNB点的近似值；

3）以连续潮流计算的SNB点近似值为初值用精确算法计算SNB点精确值；

4）以计算得到的SNB点精确值为初始点，计算含有DFIG电力***稳定域的边界轨迹；

5）计算含有DFIG电力***稳定域的切平面。

2.根据权利要求1所述的含DFIG的风电场的电力***电压稳定域切平面的获取方法，其特征在于，所述步骤1）中DFIG风电场的数学模型需要用含有风电场状态参数的方程式进行约束，建立含该参数的潮流方程。

3.根据权利要求1所述的含DFIG的风电场的电力***电压稳定域切平面的获取方法，其特征在于，所述步骤2）中，使用***的DFIG风电场的状态参数潮流模型，对连续潮流法的预测步以及校正步的雅可比矩阵进行修改，将DFIG风电场的状态参数加入该雅可比矩阵对应元素中，求解潮流方程。

4.根据权利要求1所述的含DFIG的风电场的电力***电压稳定域切平面的获取方法，其特征在于：所述步骤4）中由SNB点构成的电压稳定域边界方程需调用已建立的潮流方程。

5.根据权利要求4所述的含DFIG的风电场的电力***电压稳定域切平面的获取方法，其特征在于：SNB点构成的电压稳定域边界面方程加入参数化方程，构成新方程组。

6.根据权利要求5所述的含DFIG风电场的电力***电压稳定域切平面的获取方法，其特征在于，在计算DFIG电力***稳定域的SNB点轨迹的方程组中，再次应用预测-校正的方法对静态电压稳定边界进行追踪，构成新的预测步与校正步中的雅可比矩阵。

7.根据权利要求1所述的含DFIG风电场的电力***电压稳定域切平面的获取方法，其特征在于，所述步骤5）中，根据点法式方程，求得已知SNB点在电压稳定域边界在该点处的切平面表达式。

8.根据权利要求7所述的含DFIG风电场的电力***电压稳定域切平面的获取方法，其特征在于：确定在发电机、负荷和风电场全注入空间中静态电压稳定安全域边界在SNB已知点的切平面方程，该方程为含DFIG风电场的电力***静态电压稳定安全域局部切平面解析式。