标准的matlab程序
附录1
使用牛顿拉夫逊法进行潮流计算的Matlab程序代码
% 牛拉法计算潮流程序
%-----------------------------------------------------------------------
% B1矩阵:1、支路首端号;2、末端号;3、支路阻抗;4、支路对地电纳 % 5、支路的变比;6、支路首端处于K侧为1,1侧为0
% B2矩阵:1、该节点发电机功率;2、该节点负荷功率;3、节点电压初始值 % 4、PV节点电压V的给定值;5、节点所接的无功补偿设备的容量 % 6、节点分类标号:1为平衡节点(应为1号节点);2为PQ节点;3为PV节点;
%------------------------------------------------------------------------
clear all;
format long;
n=input('请输入节点数:nodes=');
nl=input('请输入支路数:lines=');
isb=input('请输入平衡母线节点号:balance=');
pr=input('请输入误差精度:precision=');
B1=input('请输入由各支路参数形成的矩阵:B1=');
B2=input('请输入各节点参数形成的矩阵:B2=');
Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1=zeros(nl); %------------------------------------------------------------------
for i=1:nl %支路数
if B1(i,6)==0 %左节点处于1侧
p=B1(i,1);q=B1(i,2);
else %左节点处于K侧
p=B1(i,2);q=B1(i,1);
end
Y(p,q)=Y(p,q)-1./(B1(i,3)*B1(i,5)); %非对角元
Y(q,p)=Y(p,q); %非对角元
Y(q,q)=Y(q,q)+1./(B1(i,3)*B1(i,5)^2)+B1(i,4); %对角元K侧
Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4); %对角元1侧
end
%求导纳矩阵
disp('导纳矩阵 Y=');
disp(Y)
%-------------------------------------------------------------------
G=real(Y);B=imag(Y); %分解出导纳阵的实部和虚部
for i=1:n %给定各节点初始电压的实部和虚部
e(i)=real(B2(i,3));
f(i)=imag(B2(i,3));