CN110336288B - 基于矩阵运算提取雅可比元素的配电网三相潮流计算方法 - Google Patents

基于矩阵运算提取雅可比元素的配电网三相潮流计算方法 Download PDF

Info

Publication number
CN110336288B
CN110336288B CN201910695420.9A CN201910695420A CN110336288B CN 110336288 B CN110336288 B CN 110336288B CN 201910695420 A CN201910695420 A CN 201910695420A CN 110336288 B CN110336288 B CN 110336288B
Authority
CN
China
Prior art keywords
node
phase
matrix
formula
elements
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.)
Active
Application number
CN201910695420.9A
Other languages
English (en)
Other versions
CN110336288A (zh
Inventor
姚玉斌
王俊峰
吴志良
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Dalian Maritime University
Original Assignee
Dalian Maritime University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Dalian Maritime University filed Critical Dalian Maritime University
Priority to CN201910695420.9A priority Critical patent/CN110336288B/zh
Publication of CN110336288A publication Critical patent/CN110336288A/zh
Application granted granted Critical
Publication of CN110336288B publication Critical patent/CN110336288B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • HELECTRICITY
    • H02GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02JCIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J3/00Circuit arrangements for ac mains or ac distribution networks
    • H02J3/04Circuit arrangements for ac mains or ac distribution networks for connecting networks of the same frequency but supplied from different sources
    • H02J3/06Controlling transfer of power between connected networks; Controlling sharing of load between connected networks
    • HELECTRICITY
    • H02GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02JCIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J2203/00Indexing scheme relating to details of circuit arrangements for AC mains or AC distribution networks
    • H02J2203/20Simulating, e g planning, reliability check, modelling or computer assisted design [CAD]
    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02EREDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
    • Y02E60/00Enabling technologies; Technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Power Engineering (AREA)
  • Supply And Distribution Of Alternating Current (AREA)
  • Complex Calculations (AREA)

Abstract

本发明公开了一种基于矩阵运算提取雅可比元素的配电网三相潮流计算方法,利用Matlab具有丰富的函数并且擅长矩阵运算的特点来形成雅可比矩阵有效行标号数组和有效列标号数组。采用Matlab的矩阵乘法运算计算雅可比矩阵与其转置的乘积矩阵A,再用diag函数提取乘积矩阵A的对角元形成数组R,如果数组R的第i个元素为0,则雅可比矩阵第i行内元素全部为0,利用Matlab的find函数查找数组R非0元素的位置,进而形成雅可比矩阵有效行标号数组;同理形成雅可比矩阵有效列标号数组。由于Matlab的内置函数和矩阵运算速度很快,使本发明大大缩短了形成雅可比矩阵有效行数组和有效列数组的时间,提高了计算速度。

Description

基于矩阵运算提取雅可比元素的配电网三相潮流计算方法
技术领域
本发明涉及一种配电系统三相潮流计算方法,特别是一种基于矩阵运算提取雅可比矩阵元素的配电网三相潮流计算方法。
背景技术
电力系统是由生产、变换、输送、分配和使用电能的设备(元件)组成的十分庞大而复杂的统一体,分为输电系统和配电系统两部分。输电系统主要起变换和输送电能的作用,电压等级较高;配电系统主要起变换和分配电能的作用,电压等级较低。
电力系统潮流计算是研究电力系统稳态运行的一项基本计算,它根据给定的运行条件和网络结构确定整个网络的运行状态。潮流计算也是电力系统其他分析的基础,如安全分析、暂态稳定分析等都要用到潮流计算。
配电系统是由配电线路、配电变压器、配电开关、配电电容器和配电负荷等组成的直接向终端用户分配电能的一个网络系统。对配电系统进行电力系统分析首先需要进行潮流计算,得到它的运行状态。与高压输电系统的三相对称运行方式不同,配电系统的负荷和网络都可能不对称,配电系统进行潮流计算时,应考虑三相不对称的特点,建立各元件的三相模型,进行三相潮流计算。
三相潮流计算的数学模型和算法设计都非常复杂,编程和调试工作量都大,三相潮流计算研究迫切需要易于编程和调试的软件开发环境。Matlab平台可以方便地处理各种矩阵和复数运算,有大量常见且实用的函数,给编程带来很大便利,同时还具有非常强大的可视化图形处理和交互式功能,为科学研究以及工程应用提供了一种高效的编程工具,目前已经成为许多科学领域的基本工具和首选平台,在各种科学和工程计算领域得到了广泛的应用。因此Matlab平台很适合三相潮流计算算法的开发和研究。但Matlab为解释型编程语言,通常对C语言等编译型编程语言非常有效的实现方法,用Matlab实现时计算速度则很慢,因此设计潮流计算方法时要充分考虑Matlab的特点。
根据电力系统节点的特点,潮流计算把电力系统节点分成3类:节点有功功率和无功功率已知、节点电压幅值和电压相角未知的节点称为PQ节点;节点有功功率和电压幅值已知、节点无功功率和电压相角未知的节点称为PV节点;节点电压幅值和电压相角已知、节点有功功率和无功功率未知的节点称为平衡节点。
潮流计算的功率方程或电流方程是非线性方程,必须采用迭代方式求解。迭代方法求解潮流计算问题时,需要设定各节点电压的初值。电压初始化一般采用平启动,单相潮流计算采用平启动设定电压初值时,PV节点和平衡节点的节点电压幅值取给定值,PQ节点的节点电压幅值取1.0;所有节点电压的相角都取0.0。三相潮流计算采用平启动设定电压初值时,PV节点和平衡节点的三相节点电压幅值都取给定值,PQ节点的三相节点电压幅值都取1.0;所有节点A相节点电压的相角都取0°,B相节点电压的相角都取-120°,C相节点电压的相角都取120°。这里单位除相角外都采用标幺值。
设配电系统的节点数为n,其中PV节点数为m。为叙述方便,设PV节点的编号为1~m。
三相潮流计算的节点电压方程为:
I=YU (1)
式中,I为三相注入电流相量列向量,Y为3n×3n的导纳矩阵,U为三相节点电压相量列向量。
导纳矩阵由输电线路、变压器支路和无功补偿元件形成,为3n×3n阶矩阵。配电系统的输电线路一般为三相,也有的为两相或单相,因此,导纳矩阵中有些行或列的元素全部为0。
各相注入复功率为:
Figure BDA0002149208760000021
式中,
Figure BDA0002149208760000022
为节点i的A相复功率,
Figure BDA0002149208760000023
为节点i的B相复功率,
Figure BDA0002149208760000024
为节点i的C相复功率,
Figure BDA0002149208760000025
为节点i的A相电压相量,
Figure BDA0002149208760000026
为节点i的B相电压相量,
Figure BDA0002149208760000027
为节点i的C相电压相量,
Figure BDA0002149208760000028
为节点i的A相注入电流相量,
Figure BDA0002149208760000029
为节点i的B相注入电流相量,
Figure BDA00021492087600000210
为节点i的C相注入电流相量,上标*表示共轭,i=1、2、…、n。
由式(1),得:
Figure BDA0002149208760000031
式中,导纳矩阵元素
Figure BDA0002149208760000032
为节点i的A相和节点k的A相之间的互导纳,当k=i时,
Figure BDA0002149208760000033
为节点i的A相的自导纳;
Figure BDA0002149208760000034
为节点i的A相和节点k的B相之间的互导纳;
Figure BDA0002149208760000035
为节点i的A相和节点k的C相之间的互导纳;
Figure BDA0002149208760000036
为节点i的B相和节点k的B相之间的互导纳,当k=i时,
Figure BDA0002149208760000037
为节点i的B相的自导纳;
Figure BDA0002149208760000038
为节点i的B相和节点k的A相之间的互导纳;
Figure BDA0002149208760000039
为节点i的B相和节点k的C相之间的互导纳;
Figure BDA00021492087600000310
为节点i的C相和节点k的C相之间的互导纳,当k=i时,
Figure BDA00021492087600000311
为节点i的C相的自导纳;
Figure BDA00021492087600000312
为节点i的C相和节点k的A相之间的互导纳;
Figure BDA00021492087600000313
为节点i的C相和节点k的B相之间的互导纳;k=1、2、…、n。
由式(2),得:
Figure BDA00021492087600000314
将式(4)减去式(3)得到:
Figure BDA00021492087600000315
式中,
Figure BDA00021492087600000316
为节点i的A相注入电流相量偏差,
Figure BDA00021492087600000317
为节点i的B相注入电流相量偏差,
Figure BDA0002149208760000041
为节点i的C相注入电流相量偏差。
为了提高配电系统运行的可靠性,中性点一般不接地。当中性点不接地时,中性点注入电流为0,即:
Figure BDA0002149208760000042
式中,
Figure BDA0002149208760000043
为节点i的中性点注入电流。
将式(4)代入式(6),得:
Figure BDA0002149208760000044
式中,
Figure BDA0002149208760000045
为节点i的中性点注入电流相量偏差。
对于PQ节点,已知其每一相的有功功率和无功功率或者三相总有功功率和三相总无功功率,求其电压相量。当中性点不接地时,由中性点注入电流等于0可得,3个单相有功功率和3个单相无功功率并不是完全可控的,只能控制其中4个量,另外两个是状态量。根据中性点注入电流为0的条件,多出了中性点注入电流实部和虚部为0的两个方程,共有8个方程和6个状态量。8个方程分别是A、B、C相和中性点N的电流实部偏差和虚部偏差方程,6个状态量为A、B、C相电压实部和虚部,所以需要增加2个状态量,使方程个数和状态量个数一致。控制变量选取总有功功率Pi∑和3个单相无功功率,其值不变,为给定值;A相有功功率Pi A和B相有功功率Pi B作为状态变量,其值待求,为计算值,C相有功功率由Pi∑、Pi A和Pi B求得,也是计算值。当中性点接地时,去掉中性点N的电流实部偏差和虚部偏差方程,以及将A相有功功率Pi A和B相有功功率Pi B作为控制变量,为已知量。为了列写方程方便,中性点接地的节点也列写中性点N的电流实部偏差和虚部偏差方程,且中性点N的电流实部偏差和虚部偏差都置为0,对应方程各变量的系数也都置为0,相应节点的Pi A和Pi B的修正量也为0。
PV节点有足够的可调节的无功功率容量,使得它可以保证电压幅值的恒定,一般将有无功功率储备的发电厂和具有调节无功功率能力的变电所选作PV节点。PV节点电压幅值是给定的,分别为:
Figure BDA0002149208760000046
这样与PQ节点相比,PV节点需要增加3个关于电压幅值偏差的方程;同时需要将控制变量Qi A、Qi B、Qi C变为状态变量,其值为计算值。
电压幅值偏差公式如下:
Figure BDA0002149208760000051
式中,
Figure BDA0002149208760000052
为节点i的A相给定电压幅值,
Figure BDA0002149208760000053
为节点i的B相给定电压幅值,
Figure BDA0002149208760000054
为节点i的C相给定电压幅值;
Figure BDA0002149208760000055
fi A分别为
Figure BDA0002149208760000056
的实部和虚部,
Figure BDA0002149208760000057
fi B分别为
Figure BDA0002149208760000058
的实部和虚部,
Figure BDA0002149208760000059
fi C分别为
Figure BDA00021492087600000510
的实部和虚部。
可见,PV节点比PQ节点多了3个状态变量修正量和3个偏差量。为了设计方便,可以把PV节点的状态变量修正量和偏差量分成两部分,一部分是与PQ节点相同的部分,另一部分则是多出的3个状态变量修正量和3个偏差量。
对PQ节点和PV节点都适用的状态变量修正量为:
Figure BDA00021492087600000511
式中,
Figure BDA00021492087600000512
和Δfi A分别为
Figure BDA00021492087600000513
的实部修正量和虚部修正量,
Figure BDA00021492087600000514
和Δfi B分别为
Figure BDA00021492087600000515
的实部修正量和虚部修正量,
Figure BDA00021492087600000516
和Δfi C分别为
Figure BDA00021492087600000517
的实部修正量和虚部修正量,ΔPi A和ΔPi B分别为Pi A和Pi B的修正量,上标T表示转置。
对PQ节点和PV节点都适用的偏差量为:
Figure BDA00021492087600000518
式中,
Figure BDA00021492087600000519
Figure BDA00021492087600000520
分别为
Figure BDA00021492087600000521
的实部和虚部,
Figure BDA00021492087600000522
Figure BDA00021492087600000523
分别为
Figure BDA00021492087600000524
的实部和虚部,
Figure BDA00021492087600000525
Figure BDA00021492087600000526
分别为
Figure BDA00021492087600000527
的实部和虚部,
Figure BDA00021492087600000528
Figure BDA00021492087600000529
分别为
Figure BDA00021492087600000530
的实部和虚部。
为了列写方程方便,也列写平衡节点的偏差方程,且置对应方程的ΔWi P和各变量系数为0,对应的
Figure BDA00021492087600000531
也为0。
PV节点附加的状态变量修正量为:
Figure BDA00021492087600000532
式中,
Figure BDA00021492087600000533
为节点i的A相无功功率修正量,
Figure BDA00021492087600000534
为节点i的B相无功功率修正量,
Figure BDA00021492087600000535
为节点i的C相无功功率修正量,上标T表示转置。
PV节点附加的偏差量为:
Figure BDA00021492087600000536
式中,
Figure BDA00021492087600000537
为节点i的A相电压幅值偏差,
Figure BDA00021492087600000538
为节点i的B相电压幅值偏差,
Figure BDA0002149208760000061
为节点i的C相电压幅值偏差。
由式(9)和式(11),得到状态变量修正量ΔX为:
Figure BDA0002149208760000062
由式(10)和式(12),得到偏差量列向量ΔW为:
Figure BDA0002149208760000063
平衡节点不参与迭代计算,不需要计算电流偏差或电压偏差,也不需要计算电压修正量和功率修正量,迭代完成后再计算各相功率。
潮流计算的修正方程为:
ΔW=-JΔX (15)
式中,J为雅可比矩阵,ΔW为偏差量列向量,ΔX为状态变量修正量列向量。
将雅可比矩阵写成如下的分块矩阵:
Figure BDA0002149208760000064
式中,J0为8n×8n阶基本雅可比子矩阵,D为8m×3m阶子矩阵,F为3m×8m阶子矩阵,O1为8(n-m)×3m阶零矩阵,O2为3m×8(n-m)阶零矩阵,O3为3m×3m阶零矩阵。
基本雅可比子矩阵J0为:
Figure BDA0002149208760000071
式中,Jik、JDii为分块子矩阵,diag表示对角矩阵。
子矩阵Jik为:
Figure BDA0002149208760000072
式中,
Figure BDA0002149208760000073
分别为导纳矩阵元素
Figure BDA0002149208760000074
的实部和虚部,
Figure BDA0002149208760000075
分别为导纳矩阵元素
Figure BDA0002149208760000076
的实部和虚部,
Figure BDA0002149208760000077
分别为导纳矩阵元素
Figure BDA0002149208760000078
的实部和虚部;
Figure BDA0002149208760000079
分别为导纳矩阵元素
Figure BDA00021492087600000710
的实部和虚部,
Figure BDA00021492087600000711
分别为导纳矩阵元素
Figure BDA00021492087600000712
的实部和虚部,
Figure BDA00021492087600000713
分别为导纳矩阵元素
Figure BDA00021492087600000714
的实部和虚部;
Figure BDA00021492087600000715
分别为导纳矩阵元素
Figure BDA00021492087600000716
的实部和虚部,
Figure BDA00021492087600000717
分别为导纳矩阵元素
Figure BDA00021492087600000718
的实部和虚部,
Figure BDA00021492087600000719
分别为导纳矩阵元素
Figure BDA00021492087600000720
的实部和虚部。
子矩阵JDii为:
Figure BDA00021492087600000721
式中元素分别表示如下:
Figure BDA0002149208760000081
Figure BDA0002149208760000082
式中,Pi A
Figure BDA0002149208760000083
分别为节点i的A相有功功率和无功功率,Pi B
Figure BDA0002149208760000084
分别为节点i的B相有功功率和无功功率,
Figure BDA0002149208760000085
为节点i的C相无功功率,P为节点i的三相总有功功率。
式(19)中元素分别表示如下:
Figure BDA0002149208760000086
Figure BDA0002149208760000087
如果节点i的中性点接地,则式(19)的第4行和第8行、第4列和第8列元素清零。
D为分块对角阵:
D=diag(D1,D2,…,Dm) (24)
Figure BDA0002149208760000088
式中,
Figure BDA0002149208760000089
为节点i的A相电压相量的幅值,
Figure BDA00021492087600000810
Figure BDA00021492087600000811
为节点i的B相电压相量的幅值,
Figure BDA00021492087600000812
Figure BDA00021492087600000813
为节点i的C相电压相量的幅值,
Figure BDA00021492087600000814
如果节点i的中性点接地,则式(25)的第4行和第8行元素清零。
F为分块对角阵:
F=diag(F1,F2,…,Fm) (26)
Figure BDA0002149208760000091
如图1-3所示,现有基于电流偏差型直角坐标牛顿法三相潮流计算方法,主要包括以下步骤:
A、输入原始数据和初始化电压。
B、形成三相导纳矩阵。
根据输入的输电线路、变压器和无功补偿元件参数形成三相导纳矩阵。
C、令t=0。
D、形成偏差量列向量ΔW并计算最大偏差量ΔWmax
根据式(5)和式(7)计算各节点的注入电流相量偏差,根据式(8)计算PV节点的电压相量偏差,形成偏差量列向量ΔW,并计算最大偏差量ΔWmax,即绝对值最大的偏差量。
E、判断|ΔWmax|是否满足收敛精度ε,如果满足,转至步骤L。
F、形成基本雅可比子矩阵J0
计算式(17)所示的基本雅可比子矩阵J0各分块子矩阵,其中J0各分块子矩阵Jik按式(18)计算;J0各分块子矩阵JDii按式(19)-(23)计算。形成基本雅可比子矩阵J0后,平衡节点对应的基本雅可比子矩阵J0行和列元素清零。
G、形成雅可比矩阵J。
在基本雅可比子矩阵J0的基础上追加与PV节点有关的附加雅可比矩阵元素形成完整的雅可比矩阵J,计算过程如下:
按式(25)计算分块对角矩阵D的各分块矩阵,并形成如式(24)所示的对角矩阵D;按式(27)计算分块对角矩阵F的各分块矩阵,并形成如式(26)所示的对角矩阵F。由基本雅可比子矩阵J0、对角矩阵D、对角矩阵F以及零矩阵O1、O2、O3形成完整的雅可比矩阵J。
H、判断t=0是否成立,如果不成立,转至步骤J。
I、形成雅可比矩阵有效行标号数组JR和有效列标号数组JC。
雅可比矩阵J中有许多整行元素全部为0的行或整列元素全部为0的列,为无效的行或列,进行计算时,需要去掉这些行和列,只取有效的行和列。具体有两种方法,第一种方法是根据Matlab按列存储矩阵数据的特点,按列取数据较快,按行取数据较慢,因此可以通过矩阵转置将判断矩阵一行中元素是否全部为0的运算转化为判断矩阵一列中元素是否全部为0的运算;第二种方法是用Matlab的find函数查找并记录雅可比矩阵一行(列)内不为0元素的列(行)号,再用length函数计算其个数,进而判断该行(列)元素是否全部为0。
选择形成雅可比矩阵有效行标号数组和有效列标号数组的方法,如果选择第一种方法,则转步骤I1;如果选择第二种方法,则转步骤I21;
I1、令V=JT
I2、令l=1,s=0,u=0;
I3、令q=1;
I4、判断Vql=0是否成立,如果成立,转步骤I6;
I5、令s=s+1,JRs=l,转步骤I8;
I6、令q=q+1;
I7、判断q大于雅可比矩阵的列数N是否成立,如果不成立,转步骤I4;
I8、令q=1;
I9、判断Jql=0是否成立,如果成立,转步骤I11;
I10、令u=u+1,JCu=l,转步骤I13;
I11、令q=q+1;
I12、判断q大于雅可比矩阵的行数N是否成立,如果不成立,转步骤I9;
I13、令l=l+1;
I14、判断l大于雅可比矩阵的行数N是否成立,如果成立,转步骤J;否则转步骤I3;
I21、令V=JT
I22、令l=1,s=0,u=0;
I23、令A=find(V(:,l));
I24、令q=length(A);
I25、判断q=0是否成立,如果成立,转步骤I27;
I26、令s=s+1,JRs=l;
I27、令A=find(J(:,l));
I28、令q=length(A);
I29、判断q=0是否成立,如果成立,转步骤I31;
I30、令u=u+1,JCu=l;
I31、令l=l+1;
I32、判断l大于雅可比矩阵的行数N是否成立,如果成立,转步骤J;否则转步骤I23;
J、解修正方程并修正状态变量;
考虑去掉雅可比矩阵中元素全部为0的行和列后的修正方程为
ΔWJR=-JJR,JCΔXJC (28)
式中,ΔWJR为按数组JR记录的标号提取偏差量列向量ΔW的元素形成的新列向量;JJR,JC为按数组JR记录的行标号和JC记录的列标号提取雅可比矩阵J的元素形成的新矩阵;ΔXJC为按数组JC记录的标号提取状态变量修正量ΔX的元素形成的新列向量。
利用Matlab的除法运算求解式(28)所示的修正方程,得到状态变量修正量ΔXJC如下:
ΔXJC=-JJR,JC\ΔWJR (29)
式中,\为Matlab的除法运算符号,其左侧是分母,右侧是分子。
按式(30)修正状态变量如下:
Figure BDA0002149208760000111
式中,上标t表示第t次迭代;
Figure BDA0002149208760000112
分别为按数组JC记录的标号提取第t次、第t+1次迭代的状态变量X的元素形成的新列向量。
K、令t=t+1,转至步骤D。
L、输出节点及支路数据,结束。
由于配电系统的规模较大、雅可比矩阵维数很高,在Matlab编程时采用传统第一种方法或第二种方法形成雅可比矩阵有效行标号数组和有效列标号数组时,计算时间很长,难以满足三相潮流计算对速度的要求。
发明内容
为解决现有技术存在的上述问题,本发明要提出一种能大大缩短计算时间,提高计算速度的基于矩阵运算提取雅可比元素的配电网三相潮流计算方法。
本发明的基本思路是充分利用Matlab具有丰富的函数并且擅长矩阵运算的特点来形成雅可比矩阵有效行标号数组和有效列标号数组。采用Matlab的矩阵乘法运算计算雅可比矩阵与其转置的乘积矩阵A,再用diag函数提取乘积矩阵A的对角元形成数组R,如果数组R的第i个元素为0,则雅可比矩阵第i行内元素全部为0,利用Matlab的find函数查找数组R非0元素的位置,进而形成雅可比矩阵有效行标号数组;采用Matlab的矩阵乘法运算计算雅可比矩阵转置与雅可比矩阵的乘积矩阵E,再用diag函数提取乘积矩阵E的对角元形成数组C,如果数组C的第i个元素为0,则雅可比矩阵第i列内元素全部为0,利用Matlab的find函数查找数组C非0元素的位置,进而形成雅可比矩阵有效列标号数组。Matlab的矩阵运算以及diag函数和find函数都是针对整个矩阵或数组的运算,不是具体某个矩阵或数组元素的运算,运算速度非常快。
为了实现上述目的,本发明的技术方案如下:基于矩阵运算提取雅可比元素的配电网三相潮流计算方法,包括以下步骤:
A、输入原始数据和初始化电压。
B、形成三相导纳矩阵。
根据输入的输电线路、变压器和无功补偿元件参数形成三相导纳矩阵。
C、令t=0。
D、形成偏差量列向量ΔW并计算最大偏差量ΔWmax
根据式(5)和式(7)计算各节点的注入电流相量偏差,根据式(8)计算PV节点的电压相量偏差,形成偏差量列向量ΔW,并计算最大偏差量ΔWmax,即绝对值最大的偏差量。
E、判断|ΔWmax|是否满足收敛精度ε,如果满足,转至步骤L。
F、形成基本雅可比子矩阵J0
计算式(17)所示的基本雅可比子矩阵J0各分块子矩阵,其中J0各分块子矩阵Jik按式(18)计算;J0各分块子矩阵JDii按式(19)-(23)计算。形成基本雅可比子矩阵J0后,平衡节点对应的基本雅可比子矩阵J0行和列元素清零。
G、形成雅可比矩阵J。
在基本雅可比子矩阵J0的基础上追加与PV节点有关的附加雅可比矩阵元素形成完整的雅可比矩阵J,计算过程如下:
按式(25)计算分块对角矩阵D的各分块矩阵,并形成如式(24)所示的对角矩阵D;按式(27)计算分块对角矩阵F的各分块矩阵,并形成如式(26)所示的对角矩阵F。由基本雅可比子矩阵J0、对角矩阵D、对角矩阵F以及零矩阵O1、O2、O3形成完整的雅可比矩阵J。
H、判断t=0是否成立,如果不成立,转至步骤J。
I、形成雅可比矩阵有效行标号数组JR和有效列标号数组JC。
雅可比矩阵J中有许多整行元素全部为0的行或整列元素全部为0的列,为无效的行或列,进行计算时,需要去掉这些行和列,只取有效的行和列。形成雅可比矩阵有效行标号数组和有效列标号数组的具体步骤如下:
I1、令矩阵Z为雅可比矩阵J的转置,即:
Z=JT (31)
I2、计算矩阵的乘积,得:
A=J·Z (32)
E=Z·J (33)
I3、利用Matlab的函数diag提取矩阵的对角元形成数组R和C,得:
R=diag(A) (34)
C=diag(E) (35)
式中,diag为Matlab中提取矩阵对角元的函数。
I4、分别用式(36)和式(37)提取数组R和C中不为0元素的标号,得到数组JR、JC:
JR=find(R) (36)
JC=find(C) (37)
式中,find为Matlab中查找不为0的数组元素的标号。
J、解修正方程并修正状态变量。
考虑去掉雅可比矩阵中元素全部为0的行和列后的修正方程为:
ΔWJR=-JJR,JCΔXJC (38)
式中,ΔWJR为按数组JR记录的标号提取偏差量列向量ΔW的元素形成的新列向量;JJR,JC为按数组JR记录的行标号和JC记录的列标号提取雅可比矩阵J的元素形成的新矩阵;ΔXJC为按数组JC记录的标号提取状态变量修正量ΔX的元素形成的新列向量。
利用Matlab的除法运算求解式(38)所示的修正方程,得到状态变量修正量ΔXJC如下:
ΔXJC=-JJR,JC\ΔWJR (39)
式中,\为Matlab的除法运算符号,其左侧是分母,右侧是分子。
按式(40)修正状态变量如下:
Figure BDA0002149208760000141
式中,上标t表示第t次迭代;
Figure BDA0002149208760000142
分别为按数组JC记录的标号提取第t次、第t+1次迭代的状态变量X的元素形成的新列向量。
K、令t=t+1,转至步骤D。
L、输出节点及支路数据,结束。
与现有技术相比,本发明具有以下有益效果:
本发明采用Matlab的矩阵乘法运算计算雅可比矩阵与其转置的乘积矩阵A以及雅可比矩阵转置与雅可比矩阵的乘积矩阵E,用diag函数分别提取乘积矩阵A和乘积矩阵E的对角元形成数组R和数组C,再利用Matlab的find函数查找数组R和数组C的非0元素的位置,进而形成雅可比矩阵有效行标号数组JR和有效列数组JC。由于Matlab的内置函数和矩阵运算速度很快,与传统采用循环方法或逻辑函数方法相比,本发明方法大大缩短了形成雅可比矩阵有效行数组和有效列数组的时间,提高了计算速度,同时也简化了编程。
附图说明
本发明共有附图4张。其中:
图1是现有配电系统三相潮流计算的流程图。
图2是现有形成雅可比矩阵有效行和有效列数组的方法1的流程图。
图3是现有形成雅可比矩阵有效行和有效列数组的方法2的流程图。
图4是本发明方法配电系统三相潮流计算的流程图。
具体实施方式
下面结合附图对本发明进行进一步地说明,按照图4所示配电系统三相潮流计算流程对4个大型实际配电系统算例进行了计算。
算例1:有269个节点,264条线路,7个变压器,220个负荷节点,2个电源节点,0个补偿电容器节点。
算例2:有665个节点,653条线路,13个变压器,605个负荷节点,2个电源节点,1个补偿电容器节点。
算例3:有1447个节点,1437条线路,12个变压器,548个负荷节点,2个电源节点,5个补偿电容器节点。
算例4:有2661个节点,2647条线路,16个变压器,1498个负荷节点,2个电源节点,5个补偿电容器节点。
采用本发明和几种对比方法对这4个实际系统算例进行了计算,各种潮流计算均采用了Matlab的稀疏技术,计算时收敛精度为0.00001。几种对比三相潮流计算算法分别为:
方法1:采用传统方法1形成雅可比矩阵有效行数组和有效列数组的三相潮流计算算法。
方法2:采用传统方法2形成雅可比矩阵有效行数组和有效列数组的三相潮流计算算法。
方法3:本发明方法。
3种对比方法计算各算例时有效数组形成时间和潮流计算时间见表1~表4,潮流计算的计算时间不包括数据读入和输出及支路功率计算的时间。
表1几种三相潮流算法计算算例1时的时间比较
潮流计算方法 有效数组形成时间(s) 潮流计算时间(s) 时间占比
方法1 13.75322 14.76090 93.17%
方法2 0.10973 1.11541 9.84%
方法3 0.00761 1.01254 0.75%
从表1可见,对于具有269个节点的系统,传统方法1形成雅可比矩阵有效行数组和有效列数组的时间较长,占整个潮流计算的计算时间的93.17%;传统方法2形成雅可比矩阵有效行数组和有效列数组的时间虽然大大减少,但占整个潮流计算的计算时间也高达9.84%,对形成雅可比矩阵有效行数组和有效列数组这一简单功能的计算时间也显得太长了。本发明形成雅可比矩阵有效行数组和有效列数组的时间非常短,仅占整个潮流计算的计算时间的0.75%。
表2几种三相潮流算法计算算例2时的时间比较
潮流计算方法 有效数组形成时间(s) 潮流计算时间(s) 时间占比
方法1 127.22290 130.41379 97.55%
方法2 0.46918 3.66107 12.82%
方法3 0.01746 3.20539 0.54%
从表2可见,对于具有665个节点的系统,传统方法1形成雅可比矩阵有效行数组和有效列数组的时间很长,竟高达127s,占整个潮流计算的计算时间的97.55%,无法接受;传统方法2形成雅可比矩阵有效行数组和有效列数组的时间虽然大大减少,但占整个潮流计算的计算时间也高达12.82%,对形成雅可比矩阵有效行数组和有效列数组这一简单功能的计算时间也显得太长了。本发明形成雅可比矩阵有效行数组和有效列数组的时间非常短,仅占整个潮流计算的计算时间的0.54%。
表3几种三相潮流算法计算算例3时的时间比较
潮流计算方法 有效数组形成时间(s) 潮流计算时间(s) 时间占比
方法2 1.823120 14.86029 12.27%
方法3 0.031793 13.06857 0.24%
对于具有1447个节点的系统,传统方法1形成雅可比矩阵有效行数组和有效列数组的时间非常长,无法接受,不得不中断运行,故没有得到计算时间;从表3可见,传统方法2形成雅可比矩阵有效行数组和有效列数组的时间虽然较少,但占整个潮流计算的计算时间也高达12.27%,对形成雅可比矩阵有效行数组和有效列数组这一简单功能的计算时间也显得太长了。本发明形成雅可比矩阵有效行数组和有效列数组的时间非常短,仅占整个潮流计算的计算时间的0.24%。
表4几种三相潮流算法计算算例4时的时间比较
Figure BDA0002149208760000161
Figure BDA0002149208760000171
对于具有2661个节点的系统,传统方法1形成雅可比矩阵有效行数组和有效列数组的时间非常长,无法接受,不得不中断运行,故没有得到计算时间;从表4可见,传统方法2形成雅可比矩阵有效行数组和有效列数组的时间较长,占整个潮流计算的计算时间高达11.41%,对形成雅可比矩阵有效行数组和有效列数组这一简单功能的计算时间也显得太长了。本发明形成雅可比矩阵有效行数组和有效列数组的时间仍然非常短,仅占整个潮流计算的计算时间的0.14%。本发明不仅计算时间短,程序也进一步简化。
本发明可以在任何版本的MATLAB编程语言实现,但建议使用较新版本的MATLAB语言。
本发明不局限于本实施例,任何在本发明披露的技术范围内的等同构思或者改变,均列为本发明的保护范围。

Claims (1)

1.基于矩阵运算提取雅可比元素的配电网三相潮流计算方法,包括以下步骤:
A、输入原始数据和初始化电压;
根据电力系统节点的特点,潮流计算把电力系统节点分成3类:节点有功功率和无功功率已知,节点电压幅值和电压相角未知的节点称为PQ节点;节点有功功率和电压幅值已知,节点无功功率和电压相角未知的节点称为PV节点;节点电压幅值和电压相角已知,节点有功功率和无功功率未知的节点称为平衡节点;
B、形成三相导纳矩阵;
C、设置迭代计数t=0;
D、形成偏差量列向量ΔW并计算最大偏差量ΔWmax
PQ节点和PV节点的电流相量偏差公式为:
Figure FDA0003650269760000011
Figure FDA0003650269760000012
式中,
Figure FDA0003650269760000013
为节点i的A相注入电流相量偏差,
Figure FDA0003650269760000014
为节点i的B相注入电流相量偏差,
Figure FDA0003650269760000015
为节点i的C相注入电流相量偏差;
Figure FDA0003650269760000016
为节点i的中性点注入电流相量偏差;
Figure FDA0003650269760000017
为节点i的A相复功率,
Figure FDA0003650269760000018
为节点i的B相复功率,
Figure FDA0003650269760000019
为节点i的C相复功率;
Figure FDA00036502697600000110
为节点i的A相电压相量,
Figure FDA00036502697600000111
为节点i的B相电压相量,
Figure FDA00036502697600000112
为节点i的C相电压相量;导纳矩阵元素
Figure FDA00036502697600000113
为节点i的A相和节点k的A相之间的互导纳,当k=i时,
Figure FDA00036502697600000114
为节点i的A相的自导纳;
Figure FDA00036502697600000115
为节点i的A相和节点k的B相之间的互导纳;
Figure FDA00036502697600000116
为节点i的A相和节点k的C相之间的互导纳;
Figure FDA00036502697600000117
为节点i的B相和节点k的B相之间的互导纳,当k=i时,
Figure FDA00036502697600000118
为节点i的B相的自导纳;
Figure FDA00036502697600000119
为节点i的B相和节点k的A相之间的互导纳;
Figure FDA0003650269760000021
为节点i的B相和节点k的C相之间的互导纳;
Figure FDA0003650269760000022
为节点i的C相和节点k的C相之间的互导纳,当k=i时,
Figure FDA0003650269760000023
为节点i的C相的自导纳;
Figure FDA0003650269760000024
为节点i的C相和节点k的A相之间的互导纳;
Figure FDA0003650269760000025
为节点i的C相和节点k的B相之间的互导纳;上标*表示共轭;k=1、2、…、n;
PV节点的电压幅值偏差公式如下:
Figure FDA0003650269760000026
式中,
Figure FDA0003650269760000027
为节点i的A相电压幅值偏差,
Figure FDA0003650269760000028
为节点i的B相电压幅值偏差,
Figure FDA0003650269760000029
为节点i的C相电压幅值偏差;
Figure FDA00036502697600000210
为节点i的A相给定电压幅值,
Figure FDA00036502697600000211
为节点i的B相给定电压幅值,
Figure FDA00036502697600000212
为节点i的C相给定电压幅值;
Figure FDA00036502697600000213
fi A分别为
Figure FDA00036502697600000214
的实部和虚部,
Figure FDA00036502697600000215
fi B分别为
Figure FDA00036502697600000216
的实部和虚部,
Figure FDA00036502697600000217
fi C分别为
Figure FDA00036502697600000218
的实部和虚部;
偏差量列向量ΔW为:
Figure FDA00036502697600000219
式中,n为节点数;m为PV节点数;设PV节点的编号为1~m;
式(4)中ΔWi P为:
Figure FDA00036502697600000220
式中,
Figure FDA00036502697600000221
Figure FDA00036502697600000222
分别为
Figure FDA00036502697600000223
的实部和虚部,
Figure FDA00036502697600000224
和Δhi B分别为
Figure FDA00036502697600000225
的实部和虚部,
Figure FDA00036502697600000226
Figure FDA00036502697600000227
分别为
Figure FDA00036502697600000228
的实部和虚部,
Figure FDA00036502697600000229
Figure FDA00036502697600000230
分别为
Figure FDA00036502697600000231
的实部和虚部;上标T表示转置;
式(4)中ΔWi PV0为:
Figure FDA00036502697600000232
在偏差量列向量ΔW中查找绝对值最大的值,得到最大偏差量ΔWmax
平衡节点不参与迭代计算,不需要计算电流相量偏差或电压幅值偏差;
E、判断|ΔWmax|是否满足收敛精度ε,如果满足,转至步骤L;否则,执行步骤F;
F、形成基本雅可比子矩阵J0
基本雅可比子矩阵J0为:
Figure FDA0003650269760000031
式中,Jik、JDii为分块子矩阵,diag表示对角矩阵;
式(7)中子矩阵Jik为:
Figure FDA0003650269760000032
式中,
Figure FDA0003650269760000033
分别为导纳矩阵元素
Figure FDA0003650269760000034
的实部和虚部,
Figure FDA0003650269760000035
分别为导纳矩阵元素
Figure FDA0003650269760000036
的实部和虚部,
Figure FDA0003650269760000037
分别为导纳矩阵元素
Figure FDA0003650269760000038
的实部和虚部;
Figure FDA0003650269760000039
分别为导纳矩阵元素
Figure FDA00036502697600000310
的实部和虚部,
Figure FDA00036502697600000311
分别为导纳矩阵元素
Figure FDA00036502697600000312
的实部和虚部,
Figure FDA00036502697600000313
分别为导纳矩阵元素
Figure FDA00036502697600000314
的实部和虚部;
Figure FDA00036502697600000315
分别为导纳矩阵元素
Figure FDA00036502697600000316
的实部和虚部,
Figure FDA00036502697600000317
分别为导纳矩阵元素
Figure FDA00036502697600000318
的实部和虚部,
Figure FDA00036502697600000319
分别为导纳矩阵元素
Figure FDA00036502697600000320
的实部和虚部;
式(7)中子矩阵JDii为:
Figure FDA0003650269760000041
式(9)中元素分别表示如下:
Figure FDA0003650269760000042
Figure FDA0003650269760000043
Figure FDA0003650269760000044
Figure FDA0003650269760000045
式中,Pi A
Figure FDA0003650269760000046
分别为节点i的A相有功功率和无功功率,Pi B
Figure FDA0003650269760000047
分别为节点i的B相有功功率和无功功率,
Figure FDA0003650269760000048
为节点i的C相无功功率,P为节点i的三相总有功功率;
如果节点i的中性点接地,则式(9)的第4行和第8行、第4列和第8列元素清零;
形成基本雅可比子矩阵J0后,平衡节点对应的雅可比矩阵行和列元素清零;
G、形成雅可比矩阵J;
在基本雅可比子矩阵J0的基础上追加与PV节点有关的附加雅可比矩阵元素形成完整的雅可比矩阵J如下:
Figure FDA0003650269760000051
式中,J0为(8n)×(8n)阶基本雅可比子矩阵;D为(8m)×(3m)阶子矩阵;F为(3m)×(8m)阶子矩阵;O1为(8n-8m)×(3m)阶零矩阵;O2为(3m)×(8n-8m)阶零矩阵;O3为(3m)×(3m)阶零矩阵;
式(14)中D为分块对角阵:
D=diag(D1,D2,…,Dm) (15)
式(15)中元素分别表示如下:
Figure FDA0003650269760000052
式中,
Figure FDA0003650269760000053
为节点i的A相电压相量的幅值,
Figure FDA0003650269760000054
Figure FDA0003650269760000055
为节点i的B相电压相量的幅值,
Figure FDA0003650269760000056
Figure FDA0003650269760000057
为节点i的C相电压相量的幅值,
Figure FDA0003650269760000058
如果节点i的中性点接地,则式(16)的第4行和第8行元素清零;
式(14)中F为分块对角阵:
F=diag(F1,F2,…,Fm) (17)
式(17)中元素分别表示如下:
Figure FDA0003650269760000059
由基本雅可比子矩阵J0、对角矩阵D、对角矩阵F以及零矩阵O1、O2、O3形成完整的雅可比矩阵J;
其特征在于:还包括以下步骤:
I、形成雅可比矩阵有效行标号数组JR和有效列标号数组JC;
雅可比矩阵J中有许多整行元素全部为0的行或整列元素全部为0的列,为无效的行或列,进行计算时,需要去掉这些行和列,只取有效的行和列;形成雅可比矩阵有效行标号数组和有效列标号数组的具体步骤如下:
I1、令矩阵Z为雅可比矩阵J的转置,即:
Z=JT (19)
I2、计算矩阵的乘积,得:
A=J·Z (20)
E=Z·J (21)
I3、利用Matlab的函数diag提取矩阵的对角元形成数组R和C,得:
R=diag(A) (22)
C=diag(E) (23)
式中,diag为Matlab中提取矩阵对角元的函数;
I4、分别用式(24)和式(25)提取数组R和C中不为0元素的标号,得到数组JR、JC:
JR=find(R) (24)
JC=find(C) (25)
式中,find为Matlab中查找不为0的数组元素的标号;
J、解修正方程并修正状态变量;
潮流计算的修正方程为:
ΔW=-JΔX (26)
式中,J为雅可比矩阵;ΔW为偏差量列向量;ΔX为状态变量修正量列向量;
式(26)中ΔX为:
Figure FDA0003650269760000071
式(27)中
Figure FDA0003650269760000072
为:
Figure FDA0003650269760000073
式中,
Figure FDA0003650269760000074
和Δfi A分别为
Figure FDA0003650269760000075
的实部修正量和虚部修正量,
Figure FDA0003650269760000076
和Δfi B分别为
Figure FDA0003650269760000077
的实部修正量和虚部修正量,
Figure FDA0003650269760000078
和Δfi C分别为
Figure FDA0003650269760000079
的实部修正量和虚部修正量;ΔPi A和ΔPi B分别为Pi A和Pi B的修正量;
式(27)中
Figure FDA00036502697600000710
为:
Figure FDA00036502697600000711
式中,
Figure FDA00036502697600000712
为节点i的A相无功功率修正量,
Figure FDA00036502697600000713
为节点i的B相无功功率修正量,
Figure FDA00036502697600000714
为节点i的C相无功功率修正量;
考虑去掉雅可比矩阵中元素全部为0的行和列后的修正方程为:
ΔWJR=-JJR,JCΔXJC (30)
式中,ΔWJR为按数组JR记录的标号提取偏差量列向量ΔW的元素形成的新列向量;JJR,JC为按数组JR记录的行标号和JC记录的列标号提取雅可比矩阵J的元素形成的新矩阵;ΔXJC为按数组JC记录的标号提取状态变量修正量ΔX的元素形成的新列向量;
利用Matlab的除法运算求解式(30)所示的修正方程,得到状态变量修正量ΔXJC
ΔXJC=-JJR,JC\ΔWJR (31)
式中,\为Matlab的除法运算符号,其左侧是分母,右侧是分子;
按式(32)修正状态变量:
Figure FDA00036502697600000715
式中,上标t表示第t次迭代;
Figure FDA0003650269760000081
分别为按数组JC记录的标号提取第t次、第t+1次迭代的状态变量X的元素形成的新列向量;
Figure FDA0003650269760000082
为按数组JC记录的标号提取第t次迭代的状态变量修正量ΔX的元素形成的新列向量;
K、令t=t+1,转至步骤D;
L、输出节点及支路数据,结束。
CN201910695420.9A 2019-07-30 2019-07-30 基于矩阵运算提取雅可比元素的配电网三相潮流计算方法 Active CN110336288B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910695420.9A CN110336288B (zh) 2019-07-30 2019-07-30 基于矩阵运算提取雅可比元素的配电网三相潮流计算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910695420.9A CN110336288B (zh) 2019-07-30 2019-07-30 基于矩阵运算提取雅可比元素的配电网三相潮流计算方法

Publications (2)

Publication Number Publication Date
CN110336288A CN110336288A (zh) 2019-10-15
CN110336288B true CN110336288B (zh) 2022-08-16

Family

ID=68148068

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910695420.9A Active CN110336288B (zh) 2019-07-30 2019-07-30 基于矩阵运算提取雅可比元素的配电网三相潮流计算方法

Country Status (1)

Country Link
CN (1) CN110336288B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111355236B (zh) * 2020-04-10 2023-05-30 大连海事大学 一种计及中性点电压变量的配电网三相潮流计算方法

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106356859A (zh) * 2016-09-29 2017-01-25 大连海事大学 一种基于Matlab的直角坐标牛顿法潮流计算方法

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106356859A (zh) * 2016-09-29 2017-01-25 大连海事大学 一种基于Matlab的直角坐标牛顿法潮流计算方法

Also Published As

Publication number Publication date
CN110336288A (zh) 2019-10-15

Similar Documents

Publication Publication Date Title
Ghatak et al. A fast and efficient load flow technique for unbalanced distribution system
CN106356859B (zh) 一种基于Matlab的直角坐标牛顿法潮流计算方法
CN106532711B (zh) 随迭代和节点类型改变雅可比矩阵的牛顿法潮流计算方法
Wei et al. An innovative fixed-pole numerical approximation for fractional order systems
CN107204617B (zh) 基于线性规划的直角坐标形式的区间潮流计算方法
Wang et al. Three-phase distribution power flow calculation for loop-based microgrids
US11853384B2 (en) Methods of patel loadflow computation for electrical power system
CN107196306B (zh) 基于Matlab稀疏矩阵的快速分解法潮流计算方法
CN111355236B (zh) 一种计及中性点电压变量的配电网三相潮流计算方法
CN106229988B (zh) 一种基于Matlab的极坐标牛顿法潮流计算方法
CN106709243B (zh) 含小阻抗支路电网的补偿法极坐标牛顿法潮流计算方法
CN109617080A (zh) 基于改进的雅可比矩阵的直角坐标牛顿法潮流计算方法
CN110336288B (zh) 基于矩阵运算提取雅可比元素的配电网三相潮流计算方法
CN107168925A (zh) 一种线性系统的状态与未知输入估计方法
CN110336287B (zh) 一种基于雅可比元素提取的配电系统三相潮流计算方法
CN106532712A (zh) 含小阻抗支路电网的补偿法直角坐标牛顿法潮流计算方法
CN109494748B (zh) 基于节点类型和修正的雅可比矩阵的牛顿法潮流计算方法
CN106410811B (zh) 首次迭代小阻抗支路端点改变雅可比矩阵的潮流计算方法
CN110417022A (zh) 矩阵乘法运算提取雅可比元素的配电网三相潮流计算方法
CN107194131A (zh) 基于Matlab稀疏矩阵的极坐标牛顿法潮流计算方法
He et al. Fast short-circuit current calculation of unbalanced distribution networks with inverter-interfaced distributed generators
CN106712029B (zh) 小阻抗支路pq端点变雅可比矩阵的牛顿法潮流计算方法
CN111049144A (zh) 基于定变化率矩阵的有源配电网潮流计算pv节点处理方法
CN114566969A (zh) 一种适合研究目的使用的直角坐标牛顿法潮流计算方法
CN107181260A (zh) 基于Matlab稀疏矩阵直角坐标牛顿法潮流计算方法

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