CN107204617A - 基于线性规划的直角坐标形式的区间潮流计算方法 - Google Patents
基于线性规划的直角坐标形式的区间潮流计算方法 Download PDFInfo
- Publication number
- CN107204617A CN107204617A CN201710059967.0A CN201710059967A CN107204617A CN 107204617 A CN107204617 A CN 107204617A CN 201710059967 A CN201710059967 A CN 201710059967A CN 107204617 A CN107204617 A CN 107204617A
- Authority
- CN
- China
- Prior art keywords
- node
- voltage
- power
- interval
- active
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000205 computational method Methods 0.000 title claims abstract description 13
- PXFBZOLANLWPMH-UHFFFAOYSA-N 16-Epiaffinine Natural products C1C(C2=CC=CC=C2N2)=C2C(=O)CC2C(=CC)CN(C)C1C2CO PXFBZOLANLWPMH-UHFFFAOYSA-N 0.000 claims abstract description 92
- 238000004364 calculation method Methods 0.000 claims abstract description 14
- 239000000243 solution Substances 0.000 claims description 59
- 238000002347 injection Methods 0.000 claims description 52
- 239000007924 injection Substances 0.000 claims description 52
- 239000011159 matrix material Substances 0.000 claims description 47
- 238000004422 calculation algorithm Methods 0.000 claims description 28
- 238000000034 method Methods 0.000 claims description 25
- 239000013598 vector Substances 0.000 claims description 25
- 230000005611 electricity Effects 0.000 claims description 19
- 230000008859 change Effects 0.000 claims description 17
- 239000000203 mixture Substances 0.000 claims description 15
- 230000003321 amplification Effects 0.000 claims description 14
- 238000003199 nucleic acid amplification method Methods 0.000 claims description 14
- 230000015572 biosynthetic process Effects 0.000 claims description 7
- 230000006835 compression Effects 0.000 claims description 7
- 238000007906 compression Methods 0.000 claims description 7
- 230000010349 pulsation Effects 0.000 claims description 6
- 238000003825 pressing Methods 0.000 claims description 4
- 230000035945 sensitivity Effects 0.000 claims description 3
- 238000003786 synthesis reaction Methods 0.000 claims description 3
- 239000011800 void material Substances 0.000 claims description 3
- 238000000342 Monte Carlo simulation Methods 0.000 description 15
- 238000004458 analytical method Methods 0.000 description 6
- 238000009826 distribution Methods 0.000 description 4
- 241000208340 Araliaceae Species 0.000 description 2
- 240000002853 Nelumbo nucifera Species 0.000 description 2
- 235000006508 Nelumbo nucifera Nutrition 0.000 description 2
- 235000006510 Nelumbo pentapetala Nutrition 0.000 description 2
- 235000005035 Panax pseudoginseng ssp. pseudoginseng Nutrition 0.000 description 2
- 235000003140 Panax quinquefolius Nutrition 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000005206 flow analysis Methods 0.000 description 2
- 235000008434 ginseng Nutrition 0.000 description 2
- 238000004519 manufacturing process Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 241001269238 Data Species 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 239000003990 capacitor Substances 0.000 description 1
- 238000011217 control strategy Methods 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000004870 electrical engineering Methods 0.000 description 1
- 239000000686 essence Substances 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 230000008707 rearrangement Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
Classifications
-
- H—ELECTRICITY
- H02—GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
- H02J—CIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
- H02J3/00—Circuit arrangements for ac mains or ac distribution networks
- H02J3/04—Circuit arrangements for ac mains or ac distribution networks for connecting networks of the same frequency but supplied from different sources
- H02J3/06—Controlling transfer of power between connected networks; Controlling sharing of load between connected networks
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16Z—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS, NOT OTHERWISE PROVIDED FOR
- G16Z99/00—Subject matter not provided for in other main groups of this subclass
-
- H02J3/382—
-
- H02J3/383—
-
- H02J3/386—
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02E—REDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
- Y02E10/00—Energy generation through renewable energy sources
- Y02E10/50—Photovoltaic [PV] energy
- Y02E10/56—Power conversion systems, e.g. maximum power point trackers
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02E—REDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
- Y02E10/00—Energy generation through renewable energy sources
- Y02E10/70—Wind energy
- Y02E10/76—Power conversion electric or electronic aspects
Landscapes
- Engineering & Computer Science (AREA)
- Power Engineering (AREA)
- Supply And Distribution Of Alternating Current (AREA)
Abstract
本发明公开了一种基于线性规划的直角坐标形式的区间潮流计算方法，包括以下步骤：步骤1、建立直角坐标形式的潮流方程，计算区间中点潮流解。假设风电机组的有功和无功出力分别可以表示成区间。建立直角坐标形式下的潮流方程，对区间中点值进行潮流计算，得到潮流解。步骤2、估计潮流所在的波动区间。步骤3、建立电压的直角坐标仿射坐标形式。步骤4、计算节点功率的仿射坐标形式。将电压仿射坐标形式代入直角坐标潮流方程，得到节点功率的仿射坐标形式。步骤5、采用线性规划压缩潮流解区间。构造一个线性规划，分别最小化和最大化噪声元，得到区间长度最小的电压实部和虚部。
Description
技术领域
本发明涉及求解电力系统区间潮流计算的技术，尤其涉及一种基于线性规 划的直角坐标形式的区间潮流计算方法，该方法在涉及区间潮流方程求解时， 将电压相量表示成直角坐标形式，并采用仿射算术将直角坐标下的区间潮流方 程问题转化为线性规划问题的求解，通过求解线性规划来获取区间潮流解的上 下边界的领域。
背景技术
新能源机组并入电网后，其出力的不确定性使电压频繁地波动，给电网安 全运行带来威胁。为分析新能源出力对电网潮流(主要是电网电压)的影响， 传统的做法是采用蒙特卡洛(Monte Carlo)方法对新能源机组出力进行模拟，在每 一个随机生成的场景下进行一次潮流计算，得到电网各节点的电压，最后将所 有场景下潮流计算得到的电压进行统计，得到每个节点电压的最大和最小值， 即电压的分布区间，若所有电压的分布区间都在电压上下限内，则在当前的新 能源机组出力下电网是安全的，否则，应采取相应的无功电压控制策略来控制 越限的节点电压，使其处于正常运行范围。随着现代电力系统的不断发展，电 网规模的不断扩大，以及越来越多新能源的接入，传统的蒙特卡洛方法在抽样 数目上会大大增加，潮流计算的规模也增大，计算时间急剧增长，不再适合用 来分析新能源出力对电网电压的影响。因此，需要探索新的分析方法，区间潮 流分析为解决这一问题提供了新的思路，但其估计的区间范围过于保守且收敛 性得不保证，基于线性规划和仿射算术的区间潮流方法很好地解决这一问题， 但极坐标形式下的潮流方法存在三角函数等非仿射形式的区间计算，对区间的 估计精度造成了影响。因此，本发明提出了基于线性规划的直角坐标形式的区 间潮流方法，避免了三角函数的仿射逼近运算，缩小了区间潮流解的范围。
区间分析计算最早由Moore于1966年提出，后来在工程领域得到了广泛应 用。近年来，区间分析逐渐被用于区间潮流计算，例如，文[1]将Krawkzyk-Moore 算子用于配电网三相潮流计算，该方法既可以计算负荷确定的潮流，也可以计 算含不确定性负荷的区间潮流。文[2]将Krawkzyk算法用于直流区间潮流的计 算，用于处理高压输电系统潮流计算中的不确定性问题，包括网络参数、机组 出力和负荷的不确定性。文[3]将基于仿射算术的Krawkzyk-Moore算法用于电力 系统潮流计算,有效缩小了区间潮流的解的范围。文[4]提出了基于线性规划的极 角坐标形式的区间潮流计算方法，此方法在计算时间、收敛性能和精度都取得 了突破。
参考文献：
文[1]王成山,王守相.基于区间算法的配电网三相潮流计算及算例分析[J]. 中国电机工程学报,2002,22(3):58-62.
文[2]郑志杰,王守相,赵龙,等.基于Krawczyk算法的直流潮流区间算法 [J].电力系统自动化,2012,36(20):50-53.
文[3]丁涛,崔翰韬,顾伟,等.基于区间和仿射运算的不确定潮流算法[J]. 电力系统自动化,2012,36(13):51-55.
文[4]Vaccaro A,Canizares C,Villacci D.An affine arithmetic-basedmethodology for reliable power flow analysis in the presence of datauncertainty[J]. Power Systems,IEEE Transactions on,2010,25(2):624-632.
发明内容
本发明的目的在于克服现有技术的缺点与不足，提供一种基于线性规划的 直角坐标形式的区间潮流计算方法，该方法在求解区间潮流方程时，先将电压 相量表示成直角坐标形式，进而将潮流方程改写成直角坐标形式。然后，将电 压的实部和虚部表示成仿射形式，此仿射形式的噪声包含了各新能机组出力对 电压实部或虚部产生的变化量(采用偏导数乘以新能源机组的变化量表示)。为 了让初始的仿射形式下的电压实部和虚部包含实际潮流所在的区间，在电压实 部和虚部的噪声元素系数上乘以一个放大系数，再将电压的实部和虚部代入直 角坐标形式下的潮流方程中，构造线性规划，分别最小化和最大化电压的实部 和虚部的范围，得到压缩后区间，作为最终区间算法估计的区间。由于潮流方 程是直角坐标形式，避免了极坐标形式中三角函数的仿射计算，简化了仿射计 算，提高了计算的精度。本发明的目的通过下述技术方案实现：
本发明的目的通过以下技术方案实现：一种基于线性规划的直角坐标形式 的区间潮流计算方法，具体可以包括以下步骤：
步骤1、建立直角坐标形式的潮流方程，计算区间中点潮流解；假设风电机 组的有功和无功出力分别可以表示成区间和(i∈nW)，nW代表所有 风电机组组成的集合，Pi 为有功出力下限，为有功出力上限，Q i为无功出力下 限，为无功出力上限，方便本方法的阐述，不考虑风电机组运行的控制模式； 建立直角坐标形式下的潮流方程，对区间中点值进行潮流计算，得到潮流解 X(0)＝[e(0),f(0)]，其中，e(0)为节点电压的实部，f(0)为节点电压的虚部；
步骤2、估计潮流所在的波动区间；根据潮流方程的雅可比矩阵的逆矩阵可 以计算电压实部e和虚部f在X(0)＝[e(0),f(0)]点对各个风电机组出力的灵敏度矩 阵；即可以将电压实部e和虚部f的变化量表示成：其中 X＝[e,f]不包含平衡节点的电压实部和虚部，W为风电机组的有功、无功出力以 及节点电压平方构成的向量，ΔW为节点功率和电压平方的变化量，为潮流 解对节点功率和电压平方的偏导数，W(0)为有功出力、无功出力以及节点电压平 方的区间中点值，同时为了使估计的电压区间包含实际潮流的电压波动的范围， 需要对变化量乘以一个放大系数α(大约1-10较合适)，这样，可以得到初始估 计的电压波动区间为：X(0)为初始潮流解，ΔX为潮流解的变化量；
步骤3、建立电压的仿射坐标形式；由于电压主要变化主要受到风电机组的 影响，因此，我们可以将风电出力看成噪声的来源，进而将电压的实部和虚部 可以表示成含噪声的仿射坐标形式；
步骤4、计算节点功率的仿射坐标形式；将步骤3中得到的电压仿射坐标形 式代入直角坐标潮流方程中，得到了节点功率的仿射坐标形式；
步骤5、采用线性规划压缩潮流解区间；将得到的节点功率的仿射形式转化 成区间形式，构造一个线性规划，分别最小化和最大化噪声元，进而得到区间 长度最小的电压实部和虚部。
在上述步骤1中，建立直角坐标形式的潮流方程，计算区间中点潮流解的 步骤具体为：
1)将节点电压写成直角坐标形式，即n为系统节点个 数，表示第i个节点的节点电压，ei和fi分别为实部和虚部；
对于PQ节点，列写潮流方程如下：
式中，PQ为PQ节点集合，ΔPi为节点i的有功不平衡量，PLi为节点的有功负荷， Pi为节点j的有功出力，Gij为导纳矩阵的第i行第j列元素的虚部，Bij为导纳矩 阵的第i行第j列元素的虚部，ΔQi节点i的无功不平衡量，Qi为节点i的无功出力， QLi为节点i的无功负荷，ei和ej为节点i和节点j的电压的实部，fi和fj分别节点 i和节点j的电压的虚部；
对于PV节点，列写潮流方程如下：
式中，PV表示PV节点集合，ΔPi为节点i的有功不平衡量，PLi为节点的有功负 荷，Pi为节点j的有功出力，Gij为导纳矩阵的第i行第j列元素的虚部，Bij为导 纳矩阵的第i行第j列元素的虚部，ei和ej为节点i和节点j的电压的实部，fi和fj分别节点i和节点j的电压的虚部；
2)将风电机组的出力用区间的中点值代替，即和 i∈nW，nW代表所有风电机组组成的集合，Pi 为有功出力下限， 为有功出力上限，Q i为无功出力下限，为无功出力上限；采用牛顿拉夫逊法 便可以得到区间中点的潮流解X(0)＝[e(0),f(0)]，e(0)为节点电压的实部，f(0)为节点 电压的虚部。
所述步骤2中的估计潮流所在的波动区间，需要计算区间中点潮流解处的 雅可比矩阵的逆矩阵；具有步骤如下：
1)已知区间中点潮流解为X(0)＝[e(0),f(0)]，可以得到功率变化量的表达式：
式中，ΔS＝[ΔP,ΔQ,ΔV2]T为节点注入功率和节点电压平方的变化量，ΔP为有功功 率的变化量，ΔQ为无功功率的变化量，ΔV2为节点电压幅值平方的变化量，为雅可比矩阵，为节点有功注入功率对电压实部的偏导数， 为有功注入功率对电压实部的偏导数，为节点电压的平方对电压实部的 偏导数，为节点有功注入功率对电压虚部的偏导数，为有功注入功率对电 压虚部的偏导数，为节点电压的平方对电压虚部的偏导数，ΔX＝[Δe,Δf]T为 潮流变化量，Δe为节点电压实部的变化量，Δf为节点电压虚部的变化量；
2)根据1)中的式子，可以推出相应的电压实部和虚部变化量的表达式如下：
式中，H＝J-1为雅可比矩阵的逆矩阵，ΔS＝[ΔP,ΔQ,ΔV2]T为节点注入功率和节点 电压平方的变化量，ΔP为有功功率的变化量，ΔQ为无功功率的变化量，ΔV2为 节点电压幅值平方的变化量，ΔX＝[Δe,Δf]T为潮流变化量，Δe为节点电压实部的 变化量，Δf为节点电压虚部的变化量，S(0)为节点注入功率和节点电压平方所在 区间的中点值，同时我们可以将H写成：
式中，为电压实部对节点注入有功功率的偏导数，为电压虚部对节点有功注入功率的偏导数，为电压实部对节点注入无功功率的偏导数，为电压虚 部对节点注入无功功率的偏导数，为电压实部对节点电压的平方的偏导数， 为电压虚部对节点电压的平方的偏导数；若只考虑风电机组的出力发生波 动，即忽略节点电压的变化，则可以将电压实部和虚部的变化量写成 W＝[P,Q]，P和Q分别为风电机组的有功和无功出力为风电机 组的有功和无功出力，ΔW为节点功率的变化量，为潮流解对节点功率的偏 导数，W(0)为有功和无功出力区间中点值；
3)为了使估计的潮流区间包含实际潮流解的波动范围，在变化量上加放大 系数α，放大系数的选取比较宽松，但α太大会导致估计的区间太大，导致最终 获取的潮流区间的放大，α太小会导致估计的区间不能完全包含实际潮流解，一 般需要根据风电波动区间的幅度选取(一般可选1-10之间)，进而可以得到估计 的潮流区间X(0)为初始潮流解，ΔX为潮流解的变化量。
所述步骤3中建立电压的仿射坐标形式的步骤具体如下：
1)确定噪声元的个数；由于只有风电机组的出力具有波动性，它们是引起 电压波动的根源，因此噪声元的个数包含了2k个，k为风电机组的个数；
2)形成仿射坐标。步骤2中已经计算出了电压实部与虚部变化量的表达式， 因此，可以将电压实部和虚部表示成如下的仿射形式：
式中，ΔPj为节点j的注入有功功率变化量，ΔQj为节点j的注入无功功率变化量，为 节点j注入有功功率区间的中点值，为节点j注入无功功率区间的中点值， 为节点i的电压实部对节点j的注入有功功率的偏导数，为节点i的电压 实部对节点j的注入无功功率的偏导数，为节点i的电压虚部对节点j的注入 有功功率的偏导数，为节点i的电压虚部对节点j的注入无功功率的偏导数； 为节点的节点注入有功功率对应的噪声元，为节点的节点注入 有功功率对应的噪声元。
所述步骤4中计算节点功率的仿射坐标形式的步骤具体如下：
1)将步骤3中电压的仿射坐标形式代入到潮流方程中；根据仿射算术的则 运算法则，即：
xy＝x0y0+(x1y0+x0y1)ε1+(x2y0+x0y2)ε2+…+(xny0+x0yn)εn+zkεk，
式中，εi为第i个噪声元，xi和yi为区间和对应的第i个位置的噪声系数， zk为新增加的噪声元εk对应的系数，x0和y0分别为区间和对应的区间中点值；
便可以计算出潮流方程的仿射形式：
和
式中，分别为有功功率、无功功率和电压平方的仿射形式， 均为仿射计算中有功功率、无功功率和电压平方表达式出现的常数项，为仿射计算后噪声元的系数，分别为新 产生噪声元的对节点注入有功功率、无功功率和电压平方对应的系数，nN为新 产生噪声元的集合，nW代表所有风电机组组成的集合；
2)将常数项替换；将功率表示成仿射形式是为了计算出现功率的波动范围， 而产生的常数项实际上应该是原来波动区间的中点值，因此，我们需要将相应 的功率和电压仿射式中产生的常数进行替换，公式为：
式中，分别为有功功率、注入无功功率和电压区间中点值， 若为常数，则区间中点为本身；分别为节点的注入有功功率、注入 无功功率和电压平方。PLi和QLi分别为节点i的有功负荷和无功负荷。
所述步骤5中采用线性规划压缩潮流解区间的具体步骤如下：
1)将新增加的噪声元部分，按照仿射运算的逆运算，将其转化成区间形式； 例如：其中，即为所有系数的绝对值的乘积；
2)将潮流方程进行排序，依次为功潮流方程、无功潮流方程和电压潮流方 程，将相应的表达式表示成矩阵向量形式：
S＝AY+B，
式中，A矩 阵中分别代表由组成元素形成的矩阵，为仿射计算后噪声元的系数； B向量为新增噪声元合成的区间向量，[-BP,BP]为有功功率合成区间向量， [-BQ,BQ]为无功功率合成区间向量，[-BV,BV]为节点电压合成区间向量；为有 功功率对应的噪声元，为无功功率对应的噪声元；分 别表示有功功率，无功功率和节点电压平方对应的区间；
3)压缩噪声元向量；对噪声元向量初始时我 们将它们的限定在区间[-1,1]中，并乘以了相应的放大系数，目的是包含实际的 潮流区间其中，Pmin、Qmin和分别为注入有功功率、无功功 率以及电压平方的下限，Pmax、Qmax和分别为注入有功功率、无功功率以及电 压平方的上限；但实际上Y可以在包含潮流区间的基础上进一步被压缩，直到使 AY+B＝fSP，因此我们建立如下线性规划：
以及
式中，C＝fSP-B，inf(·)表示求下限，sup(·)表示求上限，nW代表所有风电 机组组成的集合，n为系统的节点总数，为节点的节点注入有功功率对 应的噪声元，为节点的节点注入有功功率对应的噪声元，为实数， 取值范围为-1到1之间，Aij和Aik为系数矩阵A的元素；
一共需要求解4k个单目标线性规划，才能得到所有压缩的噪声元，即：进一步，我们将其代入了电压实部和虚部的仿 射形式中，得到：
式中，和分别为电压实部和虚部的仿射形式，和分别为电压实部和虚部的区间中点值，和分别为电压实部和虚部噪声元对节点注入有功功率的 对应的系数，和分别为电压实部和虚部噪声元对节点注入无功功率的对应 的系数，和分别为经过仿射后的新的噪声元，nW代表所有风电机组组成 的集合；
采用仿射逆运算，便可以得到电压实部和虚部所在区间，即和emin为实部的下限，emax为实部的上限，fmin为虚部的下限，fmax为 虚部的上限；
4)求电压幅值和相角范围；对PV节点，我们不需要求电压幅值的范围， 只需要求其相角的范围。对于PQ节点，电压幅值和相角的范围都需要重新确定； 对电压幅值，有对于电压相角，arctan为正切三角 函数的反函数，为电压幅值的区间，和分别为电压实部和虚部所在区间。
本发明的目的也可以通过以下技术方案实现：一种基于线性规划的直角坐 标形式的区间潮流计算方法，考虑风电出力都在相应的区间内变化，可以包括 以下步骤：
步骤1、建立直角坐标形式的潮流方程，计算区间中点潮流解。假设风电机 组的有功和无功出力分别可以表示成区间和(i∈nW)，nW代表所有 风电机组组成的集合，Pi 为有功出力下限，为有功出力上限，Q i为无功出力下 限，为无功出力上限，方便本方法的阐述，不考虑风电机组运行的控制模式。 建立直角坐标形式下的潮流方程，对区间中点值进行潮流计算，得到潮流解 X(0)＝[e(0),f(0)]，其中，e(0)为节点电压的实部，f(0)为节点电压的虚部。
步骤2、估计潮流所在的波动区间。根据潮流方程的雅可比矩阵的逆矩阵可 以计算电压实部e和虚部f在X(0)＝[e(0),f(0)]点对各个风电机组出力的灵敏度矩 阵。即可以将电压实部e和虚部f的变化量表示成：其中 X＝[e,f]不包含平衡节点的电压实部和虚部，W为风电机组的有功、无功出力以 及节点电压平方构成的向量，ΔW为节点功率和电压平方的变化量，为潮流 解对节点功率和电压平方的偏导数，W(0)为有功出力、无功出力以及节点电压平 方的区间中点值，同时为了使估计的电压区间包含实际潮流的电压波动的范围， 需要对变化量乘以一个放大系数α(大约1-10较合适)，这样，可以得到初始估 计的电压波动区间为：X(0)为初始潮流解，ΔX为潮流解的变化量。
步骤3、建立电压的仿射坐标形式。由于电压主要受到风电机组的影响，因 此，我们可以将风电出力看成噪声的来源，进而将电压的实部和虚部可以表示 成含噪声的仿射坐标形式。
步骤4、计算节点功率的仿射坐标形式。将步骤3中得到的电压仿射坐标形 式代入直角坐标潮流方程中，得到了节点功率的仿射坐标形式。
步骤5、采用线性规划压缩潮流解区间。将得到的节点功率的仿射形式转化 成区间形式，构造一个线性规划，分别最小化和最大化噪声元，进而得到区间 长度最小的电压实部和虚部。
在上述步骤1中，建立直角坐标形式的潮流方程，计算区间中点潮流解的 步骤具体为：
1)将节点电压写成直角坐标形式，即n为系统节点个 数，表示第i个节点的节点电压，ei和fi分别为实部和虚部。对于PQ节点，列 写潮流方程如下：
式中，PQ为PQ节点集合，ΔPi为节点i的有功不平衡量，PLi为节点的有功负荷， Pi为节点j的有功出力，Gij为导纳矩阵的第i行第j列元素的虚部，Bij为导纳矩 阵的第i行第j列元素的虚部，ΔQi节点i的无功不平衡量，Qi为节点i的无功出力， QLi为节点i的无功负荷，ei和ej为节点i和节点j的电压的实部，fi和fj分别节点 i和节点j的电压的虚部。
对于PV节点，列写潮流方程如下：
式中，PV表示PV节点集合，ΔPi为节点i的有功不平衡量，PLi为节点的有功负 荷，Pi为节点j的有功出力，Gij为导纳矩阵的第i行第j列元素的虚部，Bij为导 纳矩阵的第i行第j列元素的虚部，ei和ej为节点i和节点j的电压的实部，fi和fj分别节点i和节点j的电压的虚部。
2)将风电机组的出力用区间的中点值代替，即和 i∈nW，nW代表所有风电机组组成的集合，Pi 为有功出力下限， 为有功出力上限，Q i为无功出力下限，为无功出力上限。采用牛顿拉夫逊法 便可以得到区间中点的潮流解X(0)＝[e(0),f(0)]，e(0)为节点电压的实部，f(0)为节点 电压的虚部。
所述步骤2中的估计潮流所在的波动区间，需要计算区间中点潮流解处的 雅可比矩阵的逆矩阵。具有步骤如下：
1)已知区间中点潮流解为X(0)＝[e(0),f(0)]，可以得到功率变化量的表达式：
式中，ΔS＝[ΔP,ΔQ,ΔV2]T为节点注入功率和节点电压平方的变化量，ΔP为有功功 率的变化量，ΔQ为无功功率的变化量，ΔV2为节点电压幅值平方的变化量，为雅可比矩阵，为节点有功注入功率对电压实部的偏导数， 为有功注入功率对电压实部的偏导数，为节点电压的平方对电压实部的 偏导数，为节点有功注入功率对电压虚部的偏导数，为有功注入功率对电 压虚部的偏导数，为节点电压的平方对电压虚部的偏导数，ΔX＝[Δe,Δf]T为 潮流变化量，Δe为节点电压实部的变化量，Δf为节点电压虚部的变化量。
2)根据1)中的式子，可以推出相应的电压实部和虚部变化量的表达式如下：
式中，H＝J-1为雅可比矩阵的逆矩阵，ΔS＝[ΔP,ΔQ,ΔV2]T为节点注入功率和节点 电压平方的变化量，ΔP为有功功率的变化量，ΔQ为无功功率的变化量，ΔV2为 节点电压幅值平方的变化量，ΔX＝[Δe,Δf]T为潮流变化量，Δe为节点电压实部的 变化量，Δf为节点电压虚部的变化量，S(0)为节点注入功率和节点电压平方所在 区间的中点值，同时我们可以将H写成：
式中，为电压实部对节点注入有功功率的偏导数，为电压虚部对节点 有功注入功率的偏导数，为电压实部对节点注入无功功率的偏导数，为电 压虚部对节点注入无功功率的偏导数，为电压实部对节点电压的平方的偏导 数，为电压虚部对节点电压的平方的偏导数。若只考虑风电机组的出力发生 波动，即忽略节点电压的变化，则可以将电压实部和虚部的变化量写成 W＝[P,Q]，P和Q分别为风电机组的有功和无功出力，ΔW为 节点功率的变化量，为潮流解对节点功率的偏导数，W(0)为有功和无功出力 区间中点值。
3)为了使估计的潮流区间包含实际潮流解的波动范围，在变化量上加放大 系数α，放大系数的选取比较宽松，但α太大会导致估计的区间太大，导致最终 获取的潮流区间的放大，α太小会导致估计的区间不能完全包含实际潮流解，一 般需要根据风电波动区间的幅度选取(一般可选1-10之间)，进而可以得到估计 的潮流区间X(0)为初始潮流解，ΔX为潮流解的变化量。
所述步骤3中建立电压的仿射坐标形式的步骤具体如下：
1)确定噪声元的个数。由于只有风电机组的出力具有波动性，它们是引起 电压波动的根源，因此噪声元的个数包含了2k个，k为风电机组的个数。
2)形成仿射坐标。步骤2中已经计算出了电压实部与虚部变化量的表达式， 因此，可以将电压实部和虚部表示成如下的仿射形式：
式中，ΔPj为节点j的注入有功功率变化量，ΔQj为节点j的注入无功功率变化量，为 节点j注入有功功率区间的中点值，为节点j注入无功功率区间的中点值， 为节点i的电压实部对节点j的注入有功功率的偏导数，为节点i的电压 实部对节点j的注入无功功率的偏导数，为节点i的电压虚部对节点j的注入 有功功率的偏导数，为节点i的电压虚部对节点j的注入无功功率的偏导数。 为节点的节点注入有功功率对应的噪声元，为节点的节点注入 有功功率对应的噪声元。
所述步骤4中计算节点功率的仿射坐标形式的步骤具体如下：
1)将步骤3中电压的仿射坐标形式代入到潮流方程中。根据仿射算术的则 运算法则，即：
式中，εi为第i个噪声元，xi和yi为区间和对应的第i个位置的噪声系数，zk为 新增加的噪声元εk对应的系数，x0和y0分别为区间和对应的区间中点值。
便可以计算出潮流方程的仿射形式：
和
式中，分别为有功功率、无功功率和电压平方的仿射形式， 均为仿射计算中有功功率、无功功率和电压平方表达式出现的常数项，为仿射计算后噪声元的系数，分别为新 产生噪声元的对节点注入有功功率、无功功率和电压平方对应的系数，nN为新 产生噪声元的集合，nW代表所有风电机组组成的集合。
2)将常数项替换。将功率表示成仿射形式是为了计算出现功率的波动范围， 而产生的常数项实际上应该是原来波动区间的中点值，因此，我们需要将相应 的功率和电压仿射式中产生的常数进行替换，公式为：
式中，分别为有功功率、注入无功功率和电压区间中点值， 若为常数，则区间中点为本身；分别为节点的注入有功功率、注入 无功功率和电压平方。PLi和QLi分别为节点i的有功负荷和无功负荷。
所述步骤5中采用线性规划压缩潮流解区间的具体步骤如下：
1)将新增加的噪声元部分，按照仿射运算的逆运算，将其转化成区间形式。 例如：其中，即为所有系数的绝对值的乘积。
2)将潮流方程进行排序，依次为功潮流方程、无功潮流方程和电压潮流方 程，将相应的表达式表示成矩阵向量形式：
S＝AY+B，
式中，A矩阵中分别代表由组成元素形成的矩阵，为仿射计算后噪声元的系数；B向量为新增噪声元合成的区间向量，[-BP,BP]为有功功率合成区间向量， [-BQ,BQ]为无功功率合成区间向量，[-BV,BV]为节点电压合成区间向量；为有 功功率对应的噪声元，为无功功率对应的噪声元；分 别表示有功功率，无功功率和节点电压平方对应的区间；
3)压缩噪声元向量。对噪声元向量初始时我 们将它们的限定在区间[-1,1]中，并乘以了相应的放大系数，目的是包含实际的 潮流区间其中，Pmin、Qmin和分别为注入有功功率、无功功 率以及电压平方的下限，Pmax、Qmax和分别为注入有功功率、无功功率以及电 压平方的上限。但实际上Y可以在包含潮流区间的基础上进一步被压缩，直到使 AY+B＝fSP，因此我们建立如下线性规划：
以及
式中，C＝fSP-B，inf(·)表示求下限，sup(·)表示求上限，nW代表所有风电 机组组成的集合，n为系统的节点总数，为节点的节点注入有功功率对 应的噪声元，为节点的节点注入有功功率对应的噪声元，为实数， 取值范围为-1到1之间，Aij和Aik为系数矩阵A的元素。
一共需要求解4k个单目标线性规划，才能得到所有压缩的噪声元，即：进一步，我们将其代入了电压实部和虚部的仿 射形式中，得到：
式中，和分别为电压实部和虚部的仿射形式，和分别为电压实部 和虚部的区间中点値，和分别为电压实部和虚部噪声元对节点注入有功功 率的对应的系数，和分别为电压实部和虚部噪声元对节点注入无功功率的 对应的系数，和分别为经过仿射后的新的噪声元，nW代表所有风电机组 组成的集合。
采用仿射逆运算，便可以得到电压实部和虚部所在区间，即和emin为实部的下限，emax为实部的上限，fmin为虚部的下限，fmax为 虚部的上限。
4)求电压幅值和相角范围。对PV节点，我们不需要求电压幅值的范围， 只需要求其相角的范围。对于PQ节点，电压幅值和相角的范围都需要重新确定。 对电压幅值，有对于电压相角，arctan为正切三角 函数的反函数，为电压幅值的区间，和分别为电压实部和虚部所在区间。
本发明相对于现有技术具有如下的优点及效果：
(1)本发明可用于分析风电、光伏等新能源机组出力和负荷为不确定区间 时的潮流，确定电压和平衡节点出力的波动区间，为调度运行工作人员提供电 网运行的安全信息。
(2)本发明采用的仿射算术考虑了区间运算过程中区间相关性的处理，可 以大幅缩减区间宽度。
(3)本发明在区间潮流的仿射计算中采用了直角坐标形式，避免了三角函 数的计算，无需采用切比雪夫近似，简化了计算过程，可以进一步提高区间潮 流的精度。
(4)本发明利用线性规划来获取区间潮流的上下边界，简化了区间计算。
附图说明
图1是修改后IEEE30节点系统的接线图，图中增加了3台风电机组。表 示发电机节点，表示负荷节点，表示平衡机节点，表示风电机组， 表示变压器，表示电容器。
图2是直角坐标区间算法和Monte Carlo模拟方法得到的各节点电压实部的 区间情况，纵坐标单位为标幺值(p.u.)，其中e1l和e1u为区间算法得到的电压实 部的下限和上限，e2l和e2u为Monte Carlo得到的电压实部的下限和上限。可 见，区间算法所得的电压实部的区间包含了Monte Carlo方法得到的电压实部的 区间。
图3是直角坐标区间算法和Monte Carlo模拟方法得到的各节点电压虚部的 区间情况，纵坐标单位为标幺值(p.u.)，其中f1l和f1u为区间算法得到的电压虚 部的下限和上限，f2l和f2u为Monte Carlo得到的电压虚部的下限和上限。可见， 区间算法所得的电压虚部的区间包含了Monte Carlo方法得到的电压虚部的区 间。
图4是直角坐标区间算法和Monte Carlo模拟方法得到的各节点电压幅值的 区间情况，纵坐标单位为标幺值(p.u.)，其中V1l和V1u为区间算法得到的电压 幅值的下限和上限，V2l和V2u为Monte Carlo得到的电压幅值的下限和上限。 可见，区间算法所得的电压幅值的区间包含了Monte Carlo方法得到的电压幅值 的区间。
图5是直角坐标区间算法和Monte Carlo模拟方法得到的各节点电压相角的 区间情况，纵坐标单位为度(deg.)，其中T1l和T1u为区间算法得到的电压相角 的下限和上限，T2l和T2u为Monte Carlo得到的电压相角的下限和上限。可见， 区间算法所得的电压相角的区间包含了Monte Carlo方法得到的电压相角的区 间。
具体实施方式
下面结合实施例及附图对本发明作进一步详细的描述，但本发明的实施方 式不限于此。
实施例
为便于理解本发明，下面结合附图进行阐述。
采用修改后的IEEE30节点系统进行测试，系统的接线图见附图中的图1， 该系统有24条传输线路，4台变压器，2个无功补偿点，9台发电机组，参数见 表2，其中3台为风电机组，1台为平衡机组，24个负荷节点。三台风电机组的 出力的区间如表1(表1为风电机组出力参数表)：
表1
所有参数的计算都采用标幺制，基准功率取100M V·A。
下面具体说明直角坐标区间潮流计算的算法步骤：
第一步，读取IEEE30节点数据，包括了发电机、负荷、线路、变压器和接 地电容参数。并设置风电有功出力和无功出力波动的区间范围。
第二步，利用支路追加法形成导纳矩阵。
第三步，取风电有功功率出力和无功功率出力的中点值进行潮流计算，采 用牛顿拉夫逊法计算，得到中点处的潮流解为X(0)＝[e(0),f(0)]。
第四步，计算雅可比矩阵的逆矩阵，并采用步骤3中的电压虚部和实部的 变化量计算公式，计算得到电压实部和虚部的估计区间和仿射形式的系数 并取放大系数α＝1。
第五步，采用区间的仿射算术将潮流方程转化为仿射形式，得到相关的仿 射系数即可得到步骤5中的矩阵A和向量B等线 性规划的相关参数，从而构造求解区间潮流的两类线性规划。
第六步，采用MATLAB的linprog函数求解线性规划，得到压缩后的噪声 元，其代入了电压实部和虚部的仿射形式中，便可得到它们的波动区间。
第七步，根据电压的实部和虚部的变化区间，计算出电压幅值和相角的波 动的区间。
为了进一步验证算法的有效性，我们采用了Monte Carlo法(MC)对风电 出力区间进行模拟潮流计算，假设其出力在区间内服从均匀分布，并且各风电 机组出力的样本相互独立，抽取了5000个样本，统计电压实部、电压虚部、电 压幅值和电压相角的最大和最小值。为方便作图，对节点编号进行重排，1-5号 为常规发电机节点，30号为平衡节点，其余节点按原来编号从小到大排列，并 将得到的结果和直角坐标的区间潮流计算结果进行比较，得到各节点电压实部 和虚部的区间分布情况如图2和图3所示，可见区间算法所得的电压实部和虚 部所在的区间包含了Monte Carlo方法得到的电压实部和虚部所在的区间。电压幅值和相角所在的区间分布情况如图4和图5所示，区间算法所得的电压幅值 和电压相角的区间包含了Monte Carlo方法得到的电压幅值和相角的区间。从以 上分析可知，区间潮流算法所得到的潮流区间包含了Monte Carlo模拟方法得到 的潮流区间，两者相差很小，符合区间算法的计算要求和有效性。但计算时间 方面，直角坐标区间潮流算法只需要3s,而Monte Carlo模拟需要大约2min，在 计算效率上远不及直角坐标区间潮流算法，进一步验证了它的有效性。
上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实 施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、 替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。
Claims (6)
1.基于线性规划的直角坐标形式的区间潮流计算方法，其特征在于，包括以下步骤：
步骤1、建立直角坐标形式的潮流方程，计算区间中点潮流解；假设风电机组的有功和无功出力分别可以表示成区间和nW代表所有风电机组组成的集合，P i为有功出力下限，为有功出力上限，Q i为无功出力下限，为无功出力上限，方便本方法的阐述，不考虑风电机组运行的控制模式；建立直角坐标形式下的潮流方程，对区间中点值进行潮流计算，得到潮流解X(0)＝[e(0),f(0)]，其中，e(0)为节点电压的实部，f(0)为节点电压的虚部；
步骤2、估计潮流所在的波动区间；根据潮流方程的雅可比矩阵的逆矩阵可以计算电压实部e和虚部f在X(0)＝[e(0),f(0)]点对各个风电机组出力的灵敏度矩阵；即可以将电压实部e和虚部f的变化量表示成：其中X＝[e,f]不包含平衡节点的电压实部和虚部，W为风电机组的有功、无功出力以及节点电压平方构成的向量，ΔW为节点功率和电压平方的变化量，为潮流解对节点功率和电压平方的偏导数，W(0)为有功出力、无功出力以及节点电压平方的区间中点值，同时为了使估计的电压区间包含实际潮流的电压波动的范围，需要对变化量乘以一个放大系数α(大约1-10较合适)，这样，可以得到初始估计的电压波动区间为：X(0)为初始潮流解，ΔX为潮流解的变化量；
步骤3、建立电压的仿射坐标形式；由于电压主要变化主要受到风电机组的影响，因此，我们可以将风电出力看成噪声的来源，进而将电压的实部和虚部可以表示成含噪声的仿射坐标形式；
步骤4、计算节点功率的仿射坐标形式；将步骤3中得到的电压仿射坐标形式代入直角坐标潮流方程中，得到了节点功率的仿射坐标形式；
步骤5、采用线性规划压缩潮流解区间；将得到的节点功率的仿射形式转化成区间形式，构造一个线性规划，分别最小化和最大化噪声元，进而得到区间长度最小的电压实部和虚部。
2.根据权利要求1所述的基于线性规划的直角坐标形式的区间潮流计算方法，其特征在于：
在上述步骤1中，建立直角坐标形式的潮流方程，计算区间中点潮流解的步骤具体为：
1)将节点电压写成直角坐标形式，即n为系统节点个数，表示第i个节点的节点电压，ei和fi分别为实部和虚部；
对于PQ节点，列写潮流方程如下：
式中，i∈PQ，PQ为PQ节点集合，ΔPi为节点i的有功不平衡量，PLi为节点的有功负荷，Pi为节点j的有功出力，Gij为导纳矩阵的第i行第j列元素的虚部，Bij为导纳矩阵的第i行第j列元素的虚部，ΔQi节点i的无功不平衡量，Qi为节点i的无功出力，QLi为节点i的无功负荷，ei和ej为节点i和节点j的电压的实部，fi和fj分别节点i和节点j的电压的虚部；
对于PV节点，列写潮流方程如下：
式中，i∈PV，PV表示PV节点集合，ΔPi为节点i的有功不平衡量，PLi为节点的有功负荷，Pi为节点j的有功出力，Gij为导纳矩阵的第i行第j列元素的虚部，Bij为导纳矩阵的第i行第j列元素的虚部，ei和ej为节点i和节点j的电压的实部，fi和fj分别节点i和节点j的电压的虚部；
2)将风电机组的出力用区间的中点值代替，即和nW代表所有风电机组组成的集合，P i为有功出力下限，为有功出力上限，Q i为无功出力下限，为无功出力上限；采用牛顿拉夫逊法便可以得到区间中点的潮流解X(0)＝[e(0),f(0)]，e(0)为节点电压的实部，f(0)为节点电压的虚部。
3.根据权利要求1所述的基于线性规划的直角坐标形式的区间潮流计算方法，其特征在于：
所述步骤2中的估计潮流所在的波动区间，需要计算区间中点潮流解处的雅可比矩阵的逆矩阵；具有步骤如下：
1)已知区间中点潮流解为X(0)＝[e(0),f(0)]，可以得到功率变化量的表达式：
式中，ΔS＝[ΔP,ΔQ,ΔV2]T为节点注入功率和节点电压平方的变化量，ΔP为有功功率的变化量，ΔQ为无功功率的变化量，ΔV2为节点电压幅值平方的变化量，为雅可比矩阵，为节点有功注入功率对电压实部的偏导数，为有功注入功率对电压实部的偏导数，为节点电压的平方对电压实部的偏导数，为节点有功注入功率对电压虚部的偏导数，为有功注入功率对电压虚部的偏导数，为节点电压的平方对电压虚部的偏导数，ΔX＝[Δe,Δf]T为潮流变化量，Δe为节点电压实部的变化量，Δf为节点电压虚部的变化量；
2)根据1)中的式子，可以推出相应的电压实部和虚部变化量的表达式如下：
式中，H＝J-1为雅可比矩阵的逆矩阵，ΔS＝[ΔP,ΔQ,ΔV2]T为节点注入功率和节点电压平方的变化量，ΔP为有功功率的变化量，ΔQ为无功功率的变化量，ΔV2为节点电压幅值平方的变化量，ΔX＝[Δe,Δf]T为潮流变化量，Δe为节点电压实部的变化量，Δf为节点电压虚部的变化量，S(0)为节点注入功率和节点电压平方所在区间的中点值，同时我们可以将H写成：
式中，为电压实部对节点注入有功功率的偏导数，为电压虚部对节点有功注入功率的偏导数，为电压实部对节点注入无功功率的偏导数，为电压虚部对节点注入无功功率的偏导数，为电压实部对节点电压的平方的偏导数，为电压虚部对节点电压的平方的偏导数；若只考虑风电机组的出力发生波动，即忽略节点电压的变化，则可以将电压实部和虚部的变化量写成W＝[P,Q]，P和Q分别为风电机组的有功和无功出力为风电机组的有功和无功出力，ΔW为节点功率的变化量，为潮流解对节点功率的偏导数，W(0)为有功和无功出力区间中点值；
3)为了使估计的潮流区间包含实际潮流解的波动范围，在变化量上加放大系数α，放大系数的选取比较宽松，但α太大会导致估计的区间太大，导致最终获取的潮流区间的放大，α太小会导致估计的区间不能完全包含实际潮流解，一般需要根据风电波动区间的幅度选取(一般可选1-10之间)，进而可以得到估计的潮流区间X(0)为初始潮流解，ΔX为潮流解的变化量。
4.根据权利要求1所述的基于线性规划的直角坐标形式的区间潮流计算方法，其特征在于：
所述步骤3中建立电压的仿射坐标形式的步骤具体如下：
1)确定噪声元的个数；由于只有风电机组的出力具有波动性，它们是引起电压波动的根源，因此噪声元的个数包含了2k个，k为风电机组的个数；
2)形成仿射坐标。步骤2中已经计算出了电压实部与虚部变化量的表达式，因此，可以将电压实部和虚部表示成如下的仿射形式：
式中，i＝1,2,…,n-1， ΔPj为节点j的注入有功功率变化量，ΔQj为节点j的注入无功功率变化量，为节点j注入有功功率区间的中点值，为节点j注入无功功率区间的中点值，为节点i的电压实部对节点j的注入有功功率的偏导数，为节点i的电压实部对节点j的注入无功功率的偏导数，为节点i的电压虚部对节点j的注入有功功率的偏导数，为节点i的电压虚部对节点j的注入无功功率的偏导数；为节点的节点注入有功功率对应的噪声元，为节点的节点注入有功功率对应的噪声元。
5.根据权利要求1所述的基于线性规划的直角坐标形式的区间潮流计算方法，其特征在于：
所述步骤4中计算节点功率的仿射坐标形式的步骤具体如下：
1)将步骤3中电压的仿射坐标形式代入到潮流方程中；根据仿射算术的则运算法则，即：
式中，εi为第i个噪声元，xi和yi为区间和对应的第i个位置的噪声系数，zk为新增加的噪声元εk对应的系数，x0和y0分别为区间和对应的区间中点值；
计算出潮流方程的仿射形式：
和
式中，i∈PQ，i∈PV，分别为有功功率、无功功率和电压平方的仿射形式，均为仿射计算中有功功率、无功功率和电压平方表达式出现的常数项，为仿射计算后噪声元的系数， 分别为新产生噪声元的对节点注入有功功率、无功功率和电压平方对应的系数，nN为新产生噪声元的集合，nW代表所有风电机组组成的集合；
2)将常数项替换；将功率表示成仿射形式是为了计算出现功率的波动范围，而产生的常数项实际上应该是原来波动区间的中点值，因此，我们需要将相应的功率和电压仿射式中产生的常数进行替换，公式为：
式中，分别为有功功率、注入无功功率和电压区间中点值，若为常数，则区间中点为本身；分别为节点的注入有功功率、注入无功功率和电压平方。PLi和QLi分别为节点i的有功负荷和无功负荷。
6.根据权利要求1所述的基于线性规划的直角坐标形式的区间潮流计算方法，其特征在于：
所述步骤5中采用线性规划压缩潮流解区间的具体步骤如下：
1)将新增加的噪声元部分，按照仿射运算的逆运算，将其转化成区间形式；例如：其中，即为所有系数的绝对值的乘积；
2)将潮流方程进行排序，依次为功潮流方程、无功潮流方程和电压潮流方程，将相应的表达式表示成矩阵向量形式：
S＝AY+B，
式中，A矩阵中分别代表由组成元素形成的矩阵，为仿射计算后噪声元的系数；B向量为新增噪声元合成的区间向量，[-BP,BP]为有功功率合成区间向量，[-BQ,BQ]为无功功率合成区间向量，[-BV,BV]为节点电压合成区间向量；为有功功率对应的噪声元，为无功功率对应的噪声元； 分别表示有功功率，无功功率和节点电压平方对应的区间；
3)压缩噪声元向量；对噪声元向量初始时我们将它们的限定在区间[-1,1]中，并乘以了相应的放大系数，目的是包含实际的潮流区间其中，Pmin、Qmin和分别为注入有功功率、无功功率以及电压平方的下限，Pmax、Qmax和分别为注入有功功率、无功功率以及电压平方的上限；但实际上Y可以在包含潮流区间的基础上进一步被压缩，直到使AY+B＝fSP，因此我们建立如下线性规划：
以及
式中，C＝fSP-B，inf(·)表示求下限，sup(·)表示求上限，nW代表所有风电机组组成的集合，n为系统的节点总数，为节点的节点注入有功功率对应的噪声元，为节点的节点注入有功功率对应的噪声元，为实数，取值范围为-1到1之间，Aij和Aik为系数矩阵A的元素；
一共需要求解4k个单目标线性规划，才能得到所有压缩的噪声元，即：进一步，我们将其代入了电压实部和虚部的仿射形式中，得到：
式中，i＝1,2,…,n-1，和分别为电压实部和虚部的仿射形式，和分别为电压实部和虚部的区间中点值，和分别为电压实部和虚部噪声元对节点注入有功功率的对应的系数，和分别为电压实部和虚部噪声元对节点注入无功功率的对应的系数，和分别为经过仿射后的新的噪声元，nW代表所有风电机组组成的集合；
采用仿射逆运算，便可以得到电压实部和虚部所在区间，即和emin为实部的下限，emax为实部的上限，fmin为虚部的下限，fmax为虚部的上限；
4)求电压幅值和相角范围；对PV节点，我们不需要求电压幅值的范围，只需要求其相角的范围。对于PQ节点，电压幅值和相角的范围都需要重新确定；对电压幅值，有对于电压相角，arctan为正切三角函数的反函数，为电压幅值的区间，和分别为电压实部和虚部所在区间。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710059967.0A CN107204617B (zh) | 2017-01-24 | 2017-01-24 | 基于线性规划的直角坐标形式的区间潮流计算方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710059967.0A CN107204617B (zh) | 2017-01-24 | 2017-01-24 | 基于线性规划的直角坐标形式的区间潮流计算方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN107204617A true CN107204617A (zh) | 2017-09-26 |
| CN107204617B CN107204617B (zh) | 2020-06-19 |
Family
ID=59904867
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710059967.0A Active CN107204617B (zh) | 2017-01-24 | 2017-01-24 | 基于线性规划的直角坐标形式的区间潮流计算方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107204617B (zh) |
Cited By (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108306344A (zh) * | 2018-02-12 | 2018-07-20 | 北京省鑫知识产权运营服务有限公司 | 一种混合供电系统电力平衡控制系统 |
| CN108549985A (zh) * | 2018-04-13 | 2018-09-18 | 深圳供电局有限公司 | 一种求解区间直流潮流模型的改进蒙特卡洛方法 |
| CN109586300A (zh) * | 2018-12-21 | 2019-04-05 | 深圳供电局有限公司 | 一种获取风电潮流模型中潮流变量变化区间的方法及系统 |
| CN109861231A (zh) * | 2019-02-20 | 2019-06-07 | 武汉大学 | 一种基于凸多边形的电力系统区间潮流方法 |
| CN111799799A (zh) * | 2020-07-13 | 2020-10-20 | 福州大学 | 一种基于区间泰勒展开法的交直流混合配电网区间潮流计算方法 |
| CN112510714A (zh) * | 2020-07-07 | 2021-03-16 | 广西电网有限责任公司南宁供电局 | 一种考虑风电场相关性的区间潮流求解方法及系统 |
| CN112886597A (zh) * | 2021-01-22 | 2021-06-01 | 河海大学 | 一种考虑新能源出力不确定性的辐射状配电网仿射潮流解析计算方法 |
| CN117767314A (zh) * | 2023-12-21 | 2024-03-26 | 国网山东省电力公司潍坊供电公司 | 含分布式电源配电网的改进复仿射潮流计算方法及装置 |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102709955A (zh) * | 2012-05-30 | 2012-10-03 | 中国电力科学研究院 | 一种基于多断面潮流控制的方法 |
| CN102946098A (zh) * | 2012-10-23 | 2013-02-27 | 四川大学 | 基于网络拓扑聚类的电力系统主动解列方法 |
-
2017
- 2017-01-24 CN CN201710059967.0A patent/CN107204617B/zh active Active
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102709955A (zh) * | 2012-05-30 | 2012-10-03 | 中国电力科学研究院 | 一种基于多断面潮流控制的方法 |
| CN102946098A (zh) * | 2012-10-23 | 2013-02-27 | 四川大学 | 基于网络拓扑聚类的电力系统主动解列方法 |
Cited By (12)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108306344A (zh) * | 2018-02-12 | 2018-07-20 | 北京省鑫知识产权运营服务有限公司 | 一种混合供电系统电力平衡控制系统 |
| CN108306344B (zh) * | 2018-02-12 | 2020-10-30 | 科沃智能工业(赣州)有限公司 | 一种混合供电系统电力平衡控制系统 |
| CN108549985A (zh) * | 2018-04-13 | 2018-09-18 | 深圳供电局有限公司 | 一种求解区间直流潮流模型的改进蒙特卡洛方法 |
| CN108549985B (zh) * | 2018-04-13 | 2022-04-19 | 深圳供电局有限公司 | 一种求解区间直流潮流模型的改进蒙特卡洛方法 |
| CN109586300A (zh) * | 2018-12-21 | 2019-04-05 | 深圳供电局有限公司 | 一种获取风电潮流模型中潮流变量变化区间的方法及系统 |
| CN109861231A (zh) * | 2019-02-20 | 2019-06-07 | 武汉大学 | 一种基于凸多边形的电力系统区间潮流方法 |
| CN109861231B (zh) * | 2019-02-20 | 2022-08-05 | 武汉大学 | 一种基于凸多边形的电力系统区间潮流方法 |
| CN112510714A (zh) * | 2020-07-07 | 2021-03-16 | 广西电网有限责任公司南宁供电局 | 一种考虑风电场相关性的区间潮流求解方法及系统 |
| CN111799799A (zh) * | 2020-07-13 | 2020-10-20 | 福州大学 | 一种基于区间泰勒展开法的交直流混合配电网区间潮流计算方法 |
| CN111799799B (zh) * | 2020-07-13 | 2022-03-08 | 福州大学 | 一种基于区间泰勒展开法的交直流混合配电网区间潮流计算方法 |
| CN112886597A (zh) * | 2021-01-22 | 2021-06-01 | 河海大学 | 一种考虑新能源出力不确定性的辐射状配电网仿射潮流解析计算方法 |
| CN117767314A (zh) * | 2023-12-21 | 2024-03-26 | 国网山东省电力公司潍坊供电公司 | 含分布式电源配电网的改进复仿射潮流计算方法及装置 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN107204617B (zh) | 2020-06-19 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107204617A (zh) | 基于线性规划的直角坐标形式的区间潮流计算方法 | |
| CN103150606B (zh) | 一种分布式电源最优潮流优化方法 | |
| CN105633948B (zh) | 一种分布式能源接入电力系统随机模糊潮流算法 | |
| Eltamaly et al. | Load flow analysis by gauss-seidel method; a survey | |
| CN111082427A (zh) | 一种基于全纯函数的微电网潮流计算方法 | |
| CN108683191B (zh) | 一种下垂控制型孤岛微电网的三相潮流分析方法 | |
| CN107093901A (zh) | 一种分布式潮流控制器的机电暂态模型与仿真方法 | |
| CN104993491B (zh) | 一种计及电压和无功的线性化潮流计算方法 | |
| CN108448631A (zh) | 基于模型降阶的含分布式电源接入配电网的动态相量建模方法 | |
| CN101882785A (zh) | 一种含分布式电源的智能配电网络潮流分析装置及方法 | |
| CN104113061A (zh) | 一种含分布式电源的配电网三相潮流计算方法 | |
| Zhao et al. | Distributed continuation power flow method for integrated transmission and active distribution network | |
| CN107039981A (zh) | 一种拟直流线性化概率最优潮流计算方法 | |
| CN106786629A (zh) | 一种风场内部无功电压协调控制方法 | |
| CN109861231B (zh) | 一种基于凸多边形的电力系统区间潮流方法 | |
| CN105119269B (zh) | 一种考虑多端统一潮流控制器的随机潮流计算方法 | |
| CN105226685B (zh) | 主动配电网三相有功无功协调的电压相量校正控制方法 | |
| CN110120673A (zh) | 基于戴维南等值参数辨识的分布式输配协同无功优化方法及系统 | |
| CN106294022B (zh) | 一种用于静态安全分析的雅可比矩阵冗余存储方法 | |
| CN108718091A (zh) | 一种应用于主动配电网的三相极坐标系线性潮流计算方法 | |
| CN110336299B (zh) | 一种考虑综合能源系统小干扰稳定的配网重构方法 | |
| CN109586300A (zh) | 一种获取风电潮流模型中潮流变量变化区间的方法及系统 | |
| Zhu et al. | Integrated power flow solution for interconnected transmission and distribution network | |
| CN109309385A (zh) | 一种有源配电网中混合有源滤波器优化配置方法 | |
| Bian et al. | Identification and improvement of probabilistic voltage instability modes of power system with wind power integration |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |