电力系统分析潮流计算最终完整版文档格式.docx
《电力系统分析潮流计算最终完整版文档格式.docx》由会员分享，可在线阅读，更多相关《电力系统分析潮流计算最终完整版文档格式.docx（21页珍藏版）》请在冰点文库上搜索。
%B2为节点参数矩阵
%第一列为节点注入发电功率参数
%第二列为节点负荷功率参数
%第三列为节点电压参数
%第四列
%第五列
%第六列为节点类型参数，“T为平衡节点，“2”为PQ节点，“3”为PV节点参数
%X为节点号和对地参数矩阵
%第一列为节点编号
%第二列为节点对地参数
%默认算例
%n=4;
%n1=4;
%isb=4;
%pr=0.00001;
%B1=[120.1667i00.88641;
130.1302+0.2479i0.0258i10;
140.1736+0.3306i0.0344i10;
34
0.2603+0.4959i0.0518i10];
%B2=[001002;
0-0.5-0.3i1002;
0.201.05003;
0-0.15-0.1i1.05001];
%X=[10;
20.05i;
30;
40];
clear;
clc;
num=input（'
是否采用默认数据？
（1-默认数据;
2-手动输入）'
）；
ifnum==1
n=4;
n1=4;
isb=4;
pr=0.00001;
B1=[120.1667i00.88641;
340.2603+0.4959i0.0518i10];
B2=[001002;
X=[10;
else
n=input（'
请输入节点数:
n='
）;
n1=input（'
请输入支路数:
n1='
isb=input（'
请输入平衡节点号：
isb='
pr=input（'
请输入误差精度:
pr='
B1=input（'
请输入支路参数：
B1='
B2=input（'
请输入节点参数:
B2='
X=input（'
节点号和对地参数:
X='
end
Times=1;
%迭代次数
%创建节点导纳矩阵
Y=zeros（n）;
fori=1:
n1
ifB1（i,6）==0%不含变压器的支路
p=B1（i,1）;
q=B1（i,2）;
Y（p,q）=Y（p,q）"
B1（i,3）;
Y（q,p）=Y（p,q）;
Y（p,p）=Y（p,p）+侣1（i,3）+0.5*B1（i,4）;
Y（q,q）=Y（q,q）+侣1（i,3）+0.5*B1（i,4）;
else%li■有变压器的支路
Y（p,q）=Y（p,q）-B1（i,5）/B1（i,3）;
Y（p,p）=Y（p,p）+B1（i,5）/B1（i,3）+（1-B1（i,5））/B1（i,3）;
Y（q,q）=Y（q,q）+B1（i,5）/B1（i,3）+（B1（i,5）*（B1（i,5）-1））/B1（i,3）;
Y（i,i）=Y（i,i）+X（i,2）;
%计及补偿电容电纳
disp（'
导纳矩阵为：
’）；
disp（Y）;
%显示导纳矩阵
咖始化OrgSDetaS
OrgS=zeros（2*n-2,1）;
DetaS=zeros（2*n-2,1）;
%创建OrgS,用于存储初始功率参数
h=0;
j=0;
n%对PQ节点的处理
ifi~=isb&
B2（i,6）==2M是平衡点&
是PQ点
h=h+1;
forj=1:
n
昭式8-74
%Pi=ei*（Gij*ej-Bij*fj）+fi*（Gij*fj+Bij*ej）
%Qi=fi*（Gij*ej-Bij*fj）-ei*（Gij*fj+Bij*ej）
OrgS（2*h-1,1）=OrgS（2*h-1,1）+real（B2（i,3））*（real（Y（i,j））*real（B2（j,3））-imag（Y（i,j））*imag（B2（j,3）））+imag（B2（i,3））*（real
（Y（i,j））*imag（B2（j,3））+imag（Y（i,j））*real（B2（j,3）））;
OrgS（2*h,1）
=OrgS（2*h,1）+imag（B2（i,3））*（real（Y（i,j））*real（B2（j,3））-imag（Y（i,j））*imag（B2（j,3）））-real（B2（i,3））*（real（Y（i,j））*imag（B2（j
3））+imag（Y（i,j））*real（B2（j,3）））;
n%对PV节点的处理，注意这时不可再将h初始化为0
B2（i,6）==3%不是平衡点&
是PV点
forj=1:
昭式8-75-a
OrgS（2*h,1）=OrgS（2*h,1）+imag（B2（i,3））*（real（Y（i,j））*real（B2（j,3））-imag（Y（i,j））*imag（B2（j,3）））-real（B2（i,3））*（real（Y（i,j
））*imag（B2（j,3））+imag（Y（i,j））*real（B2（j,3）））;
%创建PVU用于存储PV节点的初始电压
PVU=zeros（n-h-1,1）;
t=0;
ifB2（i,6）==3
t=t+1;
PVU（t,1）=B2（i,3）;
endend
%创建DetaS,用于存储有功功率、无功功率和电压幅值的不平衡量h=0;
B2（i,6）==2h=h+1;
DetaS（2*h-1,1）=real（B2（i,2））-OrgS（2*h-1,1）;
%delPiDetaS（2*h,1）=imag（B2（i,2））-OrgS（2*h,1）;
%delQiendendt=0;
B2（i,6）==3h=h+1;
t=t+1;
DetaS（2*h-1,1）=real（B2（i,1））-OrgS（2*h-1,1）;
%delPi
DetaS（2*h,1）=real（PVU（t,1））A2+imag（PVU（t,1））A2-real（B2（i,3））A2-imag（B2（i,3））A2;
%delUi
endend%DetaS
%创建I,用于存储节点电流参数i=zeros（n-1,1）;
h=0;
fori=1:
nifi~=isbh=h+1;
I（h,1）=（OrgS（2*h-1,1）-OrgS（2*h,1）*sqrt（-1））/conj（B2（i,3））;
%conj求共轴endend
%创建Jacbi碓可比矩阵）Jacbi=zeros（2*n-2）;
k=0;
n%对PQ节点的处理ifB2（i,6）==2h=h+1;
forj=1:
nifj~=isbk=k+1;
ifi==j%对角元素的处理Jacbi（2*h-1,2*k-1）=-imag（Y（i,j））*real（B2（i,3））+real（Y（i,j））*imag（B2（i,3））+imag（I（h,1））;
Jacbi（2*h-1,2*k）=real（Y（i,j））*real（B2（i,3））+imag（Y（i,j））*imag（B2（i,3））+real（I（h,1））;
Jacbi（2*h,2*k-1）=-Jacbi（2*h-1,2*k）+2*real（I（h,1））;
Jacbi（2*h,2*k）=Jacbi（2*h-1,2*k-1）-2*imag（I（h,1））;
else%^对角元素的处理
Jacbi（2*h-1,2*k-1）=-imag（Y（i,j））*real（B2（i,3））+real（Y（i,j））*imag（B2（i,3））;
Jacbi（2*h-1,2*k）=real（Y（i,j））*real（B2（i,3））+imag（Y（i,j））*imag（B2（i,3））;
Jacbi（2*h,2*k-1）=-Jacbi（2*h-1,2*k）;
Jacbi（2*h,2*k）=Jacbi（2*h-1,2*k-1）;
ifk==（n-1）%务用于内循环的指针置于初始值，以确保雅可比矩阵换行
k=0;
n%对PV节点的处理
ifj~=isb
k=k+1;
ifi==j%对角元素的处理
Jacbi（2*h-1,2*k-1）=-imag（Y（i,j））*real（B2（i,3））+real（Y（i,j））*imag（B2（i,3））+imag（I（h,1））;
Jacbi（2*h-1,2*k）=real（Y（i,j））*real（B2（i,3））+imag（Y（i,j））*imag（B2（i,3））+real（I（h,1））;
Jacbi（2*h,2*k-1）=2*imag（B2（i,3））;
Jacbi（2*h,2*k）=2*real（B2（i,3））;
else%^对角元素的处理
Jacbi（2*h,2*k-1）=0;
Jacbi（2*h,2*k）=0;
初始雅可比矩阵为：
disp（Jacbi）;
%求解修正方程，获取节点电压的不平衡量
DetaU=zeros（2*n-2,1）;
DetaU=inv（Jacbi）*DetaS;
%inv矩阵求逆
%DetaU
%修正节点电压j=0;
n%对PQ节点处理
ifB2（i,6）==2
j=j+1；
B2（i,3）=B2（i,3）+DetaU（2*j,1）+DetaU（2*j-1,1）*sqrt（-1）;
j=j+1;
%B2
"
始循环**********************************************************************
whileabs（max（DetaU））>
pr
B2（i,6）==2
B2（i,6）==3
OrgS（2*h-1,1）=OrgS（2*h-1,1）+real（B2（i,3））*（real（Y（i,j））*real（B2（j,3））-imag（Y（i,j））*imag（B2（j,3）））+imag（B2（i,3））*（real（Y（i,j））*imag（B2（j,3））+imag（Y（i,j））*real（B2（j,3）））;
end%OrgS
%创建DetaS
DetaS（2*h-1,1）=real（B2（i,2））-OrgS（2*h-1,1）;
DetaS（2*h,1）=imag（B2（i,2））-OrgS（2*h,1）;
%DetaS（2*h-1,1）=real（B2（i,2））-OrgS（2*h-1,1）;
%DetaS
%创建I
i=zeros（n-1,1）;
ifi~=isb
I（h,1）=（OrgS（2*h-1,1）-OrgS（2*h,1）*sqrt（-1））/conj（B2（i,3））;
%I
%创建Jacbi
Jacbi=zeros（2*n-2）;
ifi==j
Jacbi（2*h,2*k-1）=-Jacbi（2*h-1,2*k）+2*real（I（h,1））;
Jacbi（2*h,2*k）=Jacbi（2*h-1,2*k-1）-2*imag（I（h,1））;
ifk==（n-1）
end%Jacbi
%修正节点电压
Times=Times+1;
%迭代次数加1
迭代次数为：
'
disp（Times）;
收敛时电压修正量为：
：
disp（DetaU）;
fork=1:
E（k）=B2（k,3）;
e（k）=real（E（k））;
f（k）=imag（E（k））;
V（k）=sqrt（e（k）A2+f（k）A2）;
sida（k）=atan（f（k）./e（k））*180./pi;
%===============计算各输出量===========================
各节点的实际电压标幺值E为（节点号从小到大排列）：
disp（E）;
%显示各节点的实际电压标幺值E用复数表示
'
）
各节点的电压大小V为（节点号从小到大排列）：
disp（V）;
%显示各节点的电压大小V的模值
各节点的电压相角sida为（节点号从小到大排列）：
disp（sida）;
%显示各节点的电压相
forp=1:
C（p）=0；
forq=1:
C（p）=C（p）+conj（Y（p,q））*conj（E（q））;
%计算各节点的注入电流的共轴值
S（p）=E（p）*C（p）;
%计算各节点的功率S=电压X注入电流的共轴值
各节点的功率S为（节点号从小到大排列）：
disp（S）;
%显示各节点的注入功率
Sline=zeros（n1,5）;
各条支路的首端功率Si为（顺序同您输入B1时一致）：
Sline（i,1）=B1（i,1）;
Sline（i,2）=B1（i,2）;
ifB1（i,6）==0
Si（p,q）=E（p）*（conj（E（p））*conj（B1（i,4）./2）+（conj（E（p）*B1（i,5））-conj（E（q）））*conj（1./（B1（i,3）*B1（i,5））））;
Siz（i）=Si（p,q）;
Si（p,q）=E（p）*（conj（E（p））*（（1-B1（i,5））/B1（i,3））+（conj（E（p））-conj（E（q）））*（B1（i,5）/B1（i,3）））;
SSi（p,q）=Si（p,q）;
Sline（i,3）=Siz（i）;
ZF=['
S（'
num2str（p）,'
'
num2str（q）,'
）='
num2str（SSi（p,q））];
disp（ZF）;
各条支路的末端功率Sj为（顺序同您输入B1时一致）：
Sj（q,p）=E（q）*（conj（E（q））*conj（B1（i,4）./2）+（conj（E（q）./B1（i,5））-conj（E（p）））*conj（1./（B1（i,3）*B1（i,5））））;
Sjy（i）=Sj（
升级会员 