一种计及一次调频不确定性的概率潮流分析方法

专利类型：发明专利 

语 言：中文 

申 请 号：CN201510590799.9 

申 请 日：20150915 

发 明 人：赵霞罗兰余渌绿颜伟余娟 

申 请 人：重庆大学 

申请人地址：400044 重庆市沙坪坝区沙正街174号 

公 开 日：20171103 

公 开 号：CN105207204B 

代 理 人：王翔 

代理机构：重庆大学专利中心 50201 

摘　　要：本发明公开提出了一种计及一次调频不确定性的概率潮流分析方法，构造计及一次调频特性的动态潮流模型，应用点估计方法处理模型中的不确定问题，采用牛拉法求解确定性潮流问题，从而求出扰动不确定情况下，计及系统一次调频特性及其系数不确定性时的系统频率、节点电压等状态变量的统计特征值。在计及一次调频特性的动态潮流模型设计中，考虑发电机和负荷单位调节功率必然存在的不确定性，因此可以得到更接近实际情况的状态变量，为相关的分析计算提供更准确的数据。该方法可以更为有效准确地描述系统的运行状态，为电力系统进一步的分析计算提供基础。同时，点估计法具有速度快，低阶准确率高等优点，为本方法的实现提供了有效途径。 

主 权 项：一种计及一次调频不确定性的概率潮流分析方法，其特征在于，包括以下步骤：1)输入确定参数及随机变量的分布信息；输入常规潮流计算所需要的基础信息、随机变量的分布信息和随机变量的数量；所述输入常规潮流计算所需要的基础信息，包括：网络拓扑数据，支路阻抗，母线并联电导电纳，节点类型，节点负荷功率，PV节点上发电机出力数据，PQ节点上发电机出力数据，发电机的额定容量及功率因素，发电机的调差系数标幺值R_g和系统负荷扰动功率预测值；所述随机变量的分布信息，随机变量包括：节点负荷扰动预测误差ΔP_{L_pre}、负荷的单位调节功率k_L和发电机的单位调节功率k_G；所述节点负荷扰动预测误差ΔP_{L_pre}服从正态分布，确定其正态分布的期望和方差；所述发电机的单位调节功率k_G服从二项分布，确定发电机单位调节功率标幺值的两点取值K_Gg＝1/R_g和K_Gg＝0，及其对应的概率；所述负荷的单位调节功率k_L服从均匀分布，确定其标幺值k_Lg的分布区间[a',b']；所述随机变量的数量为输入系统中随机变量的数量，总个数m，其中随机变量ΔP_{L_pre}的个数为m_p,随机变量k_L的个数为m_kL，随机变量k_G的个数为m_kG；2)计算基准状态下的潮流分布；结合步骤1获得的所述常规潮流基础信息，建立如下常规潮流模型；$\left\{\begin{array}{c}{P}_{Gi}-P_{Li}^0=U_i\underset{j=1}{\overset{j=n}{\Sigma }}{U}_j(G_{ij}{\mathrm{cos\delta }}_{ij}+B_{ij}{\mathrm{sin\delta }}_{ij})\\ Q_{Gi}-Q_{Li}^0=U_i\underset{j=1}{\overset{j=n}{\Sigma }}{U}_j(G_{ij}{\mathrm{sin\delta }}_{ij}-B_{ij}{\mathrm{cos\delta }}_{ij})\end{array}\right.---\left(1\right)$公式1是基于极坐标的常规潮流方程，其中P_Gi,Q_Gi分别表示节点i所接发电机的有功出力和无功出力；分别表示基准状态下，节点i上负荷吸收的有功功率和无功功率；U_i,U_j分别表示节点i的电压幅值和节点j的电压幅值，δ_ij＝θ_i?θ_j表示节点i、j间电压的相角差，其中θ_i和θ_j分别表示节点i和节点j电压的相角；G_ij,B_ij分别表示节点i、j间的电导和电纳，若节点i、j间无支路则对应的导纳为零；运用牛拉法求解上述方程，得到基准状态下发电机有功出力无功出力节点电压幅值和相角3)构建计及系统一次调频特性的动态潮流模型；结合步骤1获得的所述常规潮流基础信息以及步骤2求解得到的所述基准状态变量，构建负荷扰动时，计及负荷和发电机一次调频作用的动态潮流模型如下：$\left\{\begin{array}{c}{P}_{Gi}-P_{Li}=U_i\underset{j=1}{\overset{j=n}{\Sigma }}{U}_j(G_{ij}{\mathrm{cos\delta }}_{ij}+B_{ij}{\mathrm{sin\delta }}_{ij})\\ Q_{Gi}^0-Q_{Li}^0=U_i\underset{j=1}{\overset{j=n}{\Sigma }}{U}_j(G_{ij}{\mathrm{sin\delta }}_{ij}-B_{ij}{\mathrm{cos\delta }}_{ij})\\ P_{Gi}=P_{Gi}^0-k_{Gi}\Delta f\\ P_{Li}=P_{Li}^{\prime }+k_{Li}\Delta f=(P_{Li}^0+{\mathrm{\Delta P}}_{Li})+k_{Li}\Delta f\end{array}\right.---\left(2\right)$公式2中，系统发生扰动后，计及一次调频作用时节点i上发电机的实际有功出力用P_Gi表示，负荷吸收的实际有功功率用P_Li表示；k_Gi表示节点i上发电机的单位调节功率,k_Li表示负荷的单位调节功率；Δf＝f?f₀为系统的频率偏移，其中f₀和f分别表示基准状态的系统频率(50Hz)和扰动后的系统频率；P′_Li表示扰动发生瞬间节点i上的有功，表示扰动前基准状态节点i上的有功负荷,ΔP_Li表示基准状态节点i上发生的实际负荷扰动功率；4)根据分布信息，计算随机变量所需的统计特征值；对概率密度为f_X(x)的变量X,求取下列统计特征值：期望μ_X：${\mu }_X=\underset{-\infty }{\overset{+\infty }{\int }}{\mathrm{xf}}_X\left(x\right)dx---\left(3\right)$标准差σ_X：${\sigma }_X=\sqrt{\left(\underset{-\infty }{\overset{+\infty }{\int }}{(x-{\mu }_X)}^2{f}_X\left(x\right)dx\right)}---\left(4\right)$X的j阶中心距M_j(X)：$M_j\left(X\right)=\underset{-\infty }{\overset{+\infty }{\int }}{(x-{\mu }_X)}^j{f}_X\left(x\right)dx---\left(5\right)$定义λ_{X_j}表示随机变量X的j阶中心距与标准差σ_X的j次方的比值：${\lambda }_{X\_j}=\frac{M_j\left(X\right)}{{\sigma }_X^j}---\left(6\right)$其中，λ_{X_1}＝0，λ_{X_2}＝1，λ_{X_3}和λ_{X_4}分别表示随机变量X的偏度系数和峰度系数；结合公式3?6，并根据步骤1获得的所述随机变量分布信息，计算所述随机变量的期望、均方差、偏度系数和峰度系数；计算过程如下：I)负荷扰动预测误差；实际的扰动功率ΔP_L由负荷扰动的预测值P_{L_pre}及预测误差ΔP_{L_pre}两部分组成，计算表达式如公式7：ΔP_L＝P_{L_pre}+ΔP_{L_pre}?????(7)用变量X1表示预测误差ΔP_{L_pre}，各统计特征值求取如下：利用公式8计算扰动预测误差ΔP_{L_pre}的偏度系数λ_{X1_3}：${\lambda }_{X1\_3}=\frac{M_3\left(X1\right)}{{\sigma }_{X1}^3}=\frac{0}{{\sigma }_{X1}^3}=0---\left(8\right)$利用公式9计算扰动预测误差ΔP_{L_pre}的峰度系数λ_{X1_4}：${\lambda }_{X1\_4}=\frac{M_4\left(X1\right)}{{\sigma }_{X1}^4}=\frac{3{\sigma }_{X1}^4}{{\sigma }_{X1}^4}=3---\left(9\right)$II)发电机的单位调节功率；机组的单位调节功率k_G服从二项分布，p表示参与一次调频的概率，对应的单位调节功率的值为1/R，其中R是机组的调差系数；q＝1?p表示不参与一次调频的概率，对应的单位调节功率的值为0；将步骤1中输入的调差系数标幺值R_g转换为其有名值R，转换公式如下：$R=R_g\frac{f_N}{P_{GN}}---\left(10\right)$其中f_N＝50Hz，为系统额定频率；P_GN是机组的额定有功功率，通过公式11求取：P_GN＝S_GNcosα?????(11)其中S_GN机组的额定容量，α为功率因素；用变量X2表示发电机的单位调节功率k_G，其各统计特征值求取如下：利用公式12计算发电机的单位调节功率k_G的期望μ_X2：${\mu }_{X2}=p\frac{1}{R}+q\times 0=\frac{p}{R}---\left(12\right)$利用公式13计算发电机的单位调节效率k_G的均方差σ_X2：${\sigma }_{X2}=\sqrt{p{(\frac{1}{R}-\frac{p}{R})}^2+q{(0-\frac{p}{R})}^2}=\frac{\sqrt{pq}}{R}---\left(13\right)$利用公式14计算发电机的单位调节功率k_G的偏度系数λ_{X2_3}：${\lambda }_{X2\_3}=\frac{M_3\left(X2\right)}{{\sigma }_{X2}^3}=\frac{p{(1/R-p/R)}^3+q{(0-p/R)}^3}{{(\sqrt{\mathrm{pq}}/R)}^3}=\frac{q^2-p^2}{\sqrt{pq}}---\left(14\right)$利用公式15计算发电机的单位调节功率k_G的峰度系数λ_{X2_4}：${\lambda }_{X2\_4}=\frac{M_4\left(X2\right)}{{\sigma }_{X2}^4}=\frac{p{(1/R-p/R)}^4+q{(0-p/R)}^4}{{(\sqrt{pq}/R)}^4}=\frac{q^3+p^3}{pq}---\left(15\right)$III)负荷的单位调节功率；负荷的单位调节功率标幺值k_Lg在[a',b']区间内均匀分布；与步骤II一样，需要先将标幺值区间[a',b']转换为有名值区间[a,b]，按照公式16转换：$k_L=k_{Lg}\frac{P_{LN}}{f_N}---\left(16\right)$其中P_LN表示额定有功；f_N＝50Hz，为系统额定频率；用变量X3表示负荷的单位调节功率k_L，其各统计特征值求取如下：利用公式17计算负荷的单位调节功率k_L的期望μ_X3：${\mu }_{X3}=\frac{b-a}{2}---\left(17\right)$利用公式18计算负荷的单位调节功率k_L的均方差σ_X3：${\sigma }_{X3}=\sqrt{\underset{a}{\overset{b}{\int }}{(x-{\mu }_{X3})}^2\frac{1}{b-a}dx}=\frac{b-a}{2\sqrt{3}}---\left(18\right)$利用公式19计算负荷的单位调节功率k_L的偏度系数λ_{X3_3}：${\lambda }_{X3\_3}=\frac{M_3\left(X3\right)}{{\sigma }_{X3}^3}=\frac{\underset{a}{\overset{b}{\int }}{(x-{\mu }_{X3})}^3\frac{1}{b-a}dx}{{\sigma }_{X3}^3}=0---\left(19\right)$利用公式20计算负荷的单位调节功率k_L的峰度系数λ_{X3_4}：${\lambda }_{X3\_4}=\frac{M_4\left(X3\right)}{{\sigma }_{X3}^4}=\frac{\underset{a}{\overset{b}{\int }}{(x-{\mu }_{X3})}^4\frac{1}{b-a}dx}{{\sigma }_{X3}^4}=\frac{9}{5}---\left(20\right)$5)根据2m+1点估计法构造随机变量的样本点，确定样本值及其权重系数；首先求取各样本值对应的权重系数；根据点估计原理，对系统中的所有随机变量z_i(z_i∈{z₁,z₂,…z_m})分别按照公式21构造三个样本值z_i,1,z_i,2,z_i,3；z_i,j＝μ_i+ξ_i,jσ_i?i＝1,…m；j＝1,2,3?????(21)其中ξ_i,j表示第i个随机变量的第j个位置系数，通过公式22求取；$\left\{\begin{array}{c}{\xi }_{i,1}=\frac{{\lambda }_{i,3}}{2}+\sqrt{{\lambda }_{i,4}-\frac{3}{4}{\lambda }_{i,3}^2}\\ {\xi }_{i,2}=\frac{{\lambda }_{i,3}}{2}-\sqrt{{\lambda }_{i,4}-\frac{3}{4}{\lambda }_{i,3}^2}\\ {\xi }_{i,3}=0\end{array}\right.---\left(22\right)$其中λ_i3和λ_i4分别表示步骤4中求取的第i个随机变量的偏度系数和峰度系数；基于公式22得到的参数，利用公式23求取各样本值对应的权重系数w_i,1、w_i,2和w_i,3；$\left\{\begin{array}{c}{w}_{i,1}=\frac{1}{{\xi }_{i1}^2-{\xi }_{i1}{\xi }_{i2}}\\ w_{i,2}=\frac{1}{{\xi }_{i2}^2-{\xi }_{i1}{\xi }_{i2}}\\ w_{i,3}=\frac{1}{m}-\frac{1}{{\lambda }_{i,4}-{\lambda }_{i,3}^2}\end{array}\right.---\left(23\right)$其中m表示系统中随机变量的个数；然后，基于得到的负荷扰动预测误差的三个样本点，按照步骤4中的公式7，计算得到负荷扰动的样本；最后构造系统样本点；构造随机变量样本值z_i,j(i＝1,…m；j＝1,2,3)对应的系统样本点时，随机变量z_i＝z_i,j，其余变量取均值，即样本点Z(i,j)＝(μ₁,…,z_i,j,…μ_m)；其中Z(1,3)＝Z(2,3)＝…＝Z(m,3)＝(μ₁,…,μ_i,…μ_m)，只需进行一次计算，故对有m个随机变量的系统，只需完成2m+1次确定性潮流计算即可；6)参数初始化；完成步骤5后，确定潮流计算所需的初始参数；设置最大迭代次数Tmax，收敛精度ε＝10^?5～10^?3；根据步骤2得到的基准状态数据，确定发电机有功出力初值、无功出力初值分别为和节点i电压幅值和相角初值分别为和系统频率初值f_N＝50Hz；利用公式24计算节点i的注入无功$Q_{i\_input}^0=Q_{Gi}^0-Q_{Li}^0---\left(24\right)$7)采用牛拉法计算各样本对应的系统状态变量值；完成步骤6后，基于步骤5得到的由负荷扰动、发电机单位调节功率和负荷单位调节功率组成的系统样本和步骤3构造的动态潮流模型，分别计算3m个样本点对应的确定性潮流，获得该样本点对应的系统频率f，节点电压幅值U和相角θ，发电机的实际有功出力P_G和负荷吸收的实际有功功率P_L，计算过程包括；i)设置迭代次数初值，计算初始不平衡量；首先设置迭代次数初值t＝0；利用公式25计算基准状态发生扰动ΔP_Li瞬间的初始不平衡量；$\left\{\begin{array}{cc}{\mathrm{\Delta P}}_i^0=P_{i\_input}^0-U_i^0\underset{j=1}{\overset{j=n}{\Sigma }}{U}_j^0(G_{ij}{\mathrm{cos\delta }}_{ij}^0+B_{ij}{\mathrm{sin\delta }}_{ij}^0)& (i\in [pq;pv;ref\left]\right)\\ {\mathrm{\Delta Q}}_i^0=Q_{i\_input}^0-U_i^0\underset{j=1}{\overset{j=n}{\Sigma }}{U}_j^0(G_{ij}{\mathrm{sin\delta }}_{ij}^0-B_{ij}{\mathrm{cos\delta }}_{ij}^0)& (i\in pq)\end{array}\right.---\left(25\right)$公式25中为基准状态节点i和节点j电压的相角差，即；${\delta }_{ij}^0={\theta }_i^0-{\theta }_j^0---\left(26\right)$表示扰动后节点i的注入有功，根据公式27求取；$\{\begin{array}{c}{P}_{i\_input}^0=P_{Gi}^0-{P^{\prime }}_{Li}^0\\ {P^{\prime }}_{Li}^0=P_{Li}^0+{\mathrm{\Delta P}}_{Li}=P_{Li}^0+P_{Li\_pre}+{\mathrm{\Delta P}}_{Li\_pre}\end{array}---\left(27\right)$其中和分别表示基准状态下节点i上发电机的有功出力和有功负荷；表示扰动瞬间节点i的有功负荷；P_{Li_pre}和ΔP_{Li_pre}分别表示扰动预测值和扰动预测误差；此时系统的不平衡量为ΔF；ii)计算雅克比矩阵；基于步骤3提出的动态潮流模型，计算雅克比矩阵；$J=\left[\begin{array}{ccc}H& N& M\\ K& L& T\\ W& O& V\end{array}\right]---\left(29\right)$其中H为n_pq+n_pv阶的方阵；$H_{ij}=\frac{\partial {\mathrm{\Delta P}}_i}{\partial {\theta }_j}=\left\{\begin{array}{cc}{U}_i^t{U}_j^t(G_{ij}{\mathrm{sin\delta }}_{ij}^t-B_{ij}{\mathrm{cos\delta }}_{ij}^t)& (i\ne j)\\ -U_i^t\underset{k\ne i}{\overset{n}{\underset{k=1}{\Sigma }}}{U}_k^t(G_{ik}{\mathrm{sin\delta }}_{ik}^t-B_{ik}{\mathrm{cos\delta }}_{ik}^t)& (i=j)\end{array}\right.,(i\in [pq;pv];j\in [pq;pv\left]\right)---\left(30\right)$N为(n_pq+n_pv)×n_pq阶的矩阵；$H_{ij}=\frac{\partial {\mathrm{\Delta P}}_i}{\partial U_j}{U}_j=\left\{\begin{array}{cc}{U}_i^t{U}_j^t(G_{ij}{\mathrm{cos\delta }}_{ij}^t+B_{ij}{\mathrm{sin\delta }}_{ij}^t)& (i\ne j)\\ U_i^t\underset{k=1}{\overset{n}{\Sigma }}{U}_k^t(G_{ik}{\mathrm{cos\delta }}_{ik}^t+B_{ik}{\mathrm{sin\delta }}_{ik}^t)+{\left(U_i^t\right)}^2{G}_{ii}& (i=j)\end{array}\right.,(i\in [pq;pv];j\in pq)---\left(31\right)$M为(n_pq+n_pv)×1阶的矩阵；$M_i=\frac{\partial {\mathrm{\Delta P}}_i}{\partial f}=(k_{Gi}+k_{Li}),(i\in [pq:pv\left]\right)---\left(32\right)$K为n_pq×(n_pq+n_pv)阶的矩阵；

关 键 词： 

法律状态：公开 

IPC专利分类号：H02J3/06