潮流计算程序
最大负荷时
%本程序的功能是用牛顿——拉夫逊法进行潮流计算
% B1矩阵：1、支路首端号；2、末端号；3、支路阻抗；4、支路对地电纳
% 5、支路的变比；6、支路首端处于K侧为1，1侧为0
% B2矩阵：1、该节点发电机功率；2、该节点负荷功率；3、节点电压初始值% 4、PV节点电压V的给定值；5、节点所接的无功补偿设备的容量
% 6、节点分类标号
clear;
n=12;%input('请输入节点数：n=');
nl=14;%input('请输入支路数：nl=');
isb=1;%input('请输入平衡母线节点号：isb=');
pr=0.00001;%input('请输入误差精度：pr=');
B1=[1 2 3.3+11.48i 0 1 0;
1 4 3.87+15.375i 0 1 0;
1 6 5.8+18.8i 0 1 0;
2 3 2.39+34.95i 0 1.02 1;
4 5 1.59+23.3i 0 1.02 1;
2 4 4.67+18.34i 0 1 0;
6 7 2.39+24.35i 0 1.02 1;
8 9 2.39+24.35i 0 1.02 1;
10 11 2.39+24.35i 0 1.02 1;
8 10 1.7+11.2i 0 1 0;
4 12 4.2+13.9i 0 1 0;
6 12 4.2+14.64i 0 1 0;
8 12 2.2+13.9i 0 1 0;
10 12 2.4+15.2i 0 1 0];
%input('请输入由支路参数形成的矩阵： B1=');
B2=[0 0 231 231 0 1;
0 0 220 0 0 2;
0 100+75i 220 0 0 2;
0 0 220 0 0 2;
0 100+75i 220 0 0 2;
0 0 220 0 0 2;
0 150+112.5i 220 0 0 2;
0 0 220 0 0 2;
0 150+112.5i 220 0 0 2;
0 0 220 0 0 2;
0 150+112.5i 220 0 0 2;
450 0 231 231 0 3];%input('请输入各节点参数形成的矩阵： B2='); Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1= zeros(nl);
%－－－－－－－修改部分－－－－－－－－－－－－
ym=1;
SB=100;UB=220;
%ym=input('您输入的参数是标么值？（若不是则输入一个不为零的数值）'); if ym~=0
%SB=input('请输入功率基准值:SB=');
%UB=input('请输入电压基准值:UB=');
YB=SB./UB./UB;
BB1=B1;
BB2=B2;
for i=1:nl
B1(i,3)=B1(i,3)*YB;
B1(i,4)=B1(i,4)./YB;
end
disp('B1矩阵B1=');
disp(B1)
for i=1:n
B2(i,1)=B2(i,1)./SB;
B2(i,2)=B2(i,2)./SB;
B2(i,3)=B2(i,3)./UB;
B2(i,4)=B2(i,4)./UB;
B2(i,5)=B2(i,5)./SB;
end
disp('B2矩阵B2=');
disp(B2)
end
% % %---------------------------------------------------
for i=1:nl %支路数
if B1(i,6)==0 %左节点处于低压侧
p=B1(i,1);q=B1(i,2);
else
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)./2; %对角元K侧
Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2; %对角元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));
V(i)=B2(i,4); %PV节点电压给定模值
end
for i=1:n %给定各节点注入功率
S(i)=B2(i,1)-B2(i,2); %i节点注入功率SG-SL
B(i,i)=B(i,i)+B2(i,5); %i节点无功补偿量
end
%=================================================================== P=real(S);Q=imag(S);
ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0;
while IT2~=0
IT2=0;a=a+1;
for i=1:n
if i~=isb %非平衡节点
C(i)=0;D(i)=0;
for j1=1:n
C(i)=C(i)+G(i,j1)*e(j1)-B(i,j1)*f(j1);%Σ(Gij*ej-Bij*fj)
D(i)=D(i)+G(i,j1)*f(j1)+B(i,j1)*e(j1);%Σ(Gij*fj+Bij*ej)
end
P1=C(i)*e(i)+f(i)*D(i);%节点功率P计算eiΣ(Gij*ej-Bij*fj)+fiΣ(Gij*fj+Bij*ej)
Q1=C(i)*f(i)-e(i)*D(i);%节点功率Q计算fiΣ(Gij*ej-Bij*fj)-eiΣ(Gij*fj+Bij*ej)
%求P',Q'
V2=e(i)^2+f(i)^2; %电压模平方
%========= 以下针对非PV节点来求取功率差及Jacobi矩阵元素=========
if B2(i,6)~=3 %非PV节点
DP=P(i)-P1; %节点有功功率差
DQ=Q(i)-Q1; %节点无功功率差
%=============== 以上为除平衡节点外其它节点的功率计算=================
%================= 求取Jacobi矩阵 ===================
for j1=1:n
if j1~=isb&j1~=i %非平衡节点&非对角元
X1=-G(i,j1)*e(i)-B(i,j1)*f(i); %
dP/de=-dQ/df
X2=B(i,j1)*e(i)-G(i,j1)*f(i); % dP/df=dQ/de X3=X2; % X2=dp/df X3=dQ/de
X4=-X1; % X1=dP/de X4=dQ/df
p=2*i-1;q=2*j1-1;
J(p,q)=X3;J(p,N)=DQ;m=p+1;
J(m,q)=X1;J(m,N)=DP;q=q+1;
J(p,q)=X4;J(m,q)=X2;
elseif j1==i&j1~=isb %非平衡节点&对角元
X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);% dP/de
X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);% dP/df
X3=D(i)+B(i,i)*e(i)-G(i,i)*f(i); % dQ/de
X4=-C(i)+G(i,i)*e(i)+B(i,i)*f(i);% dQ/df
p=2*i-1;q=2*j1-1;J(p,q)=X3;J(p,N)=DQ;%扩展列△Q
m=p+1;
J(m,q)=X1;q=q+1;J(p,q)=X4;J(m,N)=DP;%扩展列△P J(m,q)=X2;
end
end
else
%=============== 下面是针对PV节点来求取Jacobi矩阵的元素===========
DP=P(i)-P1; % PV节点有功误差
DV=V(i)^2-V2; % PV节点电压误差
for j1=1:n
if j1~=isb&j1~=i %非平衡节点&非对角元
X1=-G(i,j1)*e(i)-B(i,j1)*f(i); % dP/de
X2=B(i,j1)*e(i)-G(i,j1)*f(i); % dP/df
X5=0;X6=0;
p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV;
m=p+1;
J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6;
J(m,q)=X2;
elseif j1==i&j1~=isb %非平衡节点&对角元
X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);% dP/de
X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);% dP/df
X5=-2*e(i);
X6=-2*f(i);
p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV;
m=p+1;
J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6;
J(m,q)=X2;
end
end
end
end
end
%========= 以上为求雅可比矩阵的各个元素 =====================
for k=3:N0 % N0=2*n （从第三行开始，第一、二行是平衡节点）
k1=k+1;N1=N; % N=N0+1 即 N=2*n+1扩展列△P、△Q for k2=k1:N1 % 扩展列△P、△Q
J(k,k2)=J(k,k2)./J(k,k); % 非对角元规格化
end
J(k,k)=1; % 对角元规格化
if k~=3 % 不是第三行
%============================================================
k4=k-1;
for k3=3:k4 % 用k3行从第三行开始到当前行前的k4行消去
for k2=k1:N1 % k3行后各行下三角元素
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算
end
J(k3,k)=0;
end
if k==N0
break;
end
%==========================================
for k3=k1:N0
for k2=k1:N1
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算
end
J(k3,k)=0;
end
else
for k3=k1:N0
for k2=k1:N1
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算
end
J(k3,k)=0;
end
end
end
%====上面是用线性变换方式将Jacobi矩阵化成单位矩阵=====
for k=3:2:N0-1
L=(k+1)./2;
e(L)=e(L)-J(k,N); %修改节点电压实部
k1=k+1;
f(L)=f(L)-J(k1,N); %修改节点电压虚部
end
%------修改节点电压-----------
for k=3:N0
DET=abs(J(k,N));
if DET>=pr %电压偏差量是否满足要求
IT2=IT2+1; %不满足要求的节点数加1
end
end
ICT2(a)=IT2;
ICT1=ICT1+1;
end
%用高斯消去法解"w=-J*V"
disp('迭代次数：');
disp(ICT1);
disp('没有达到精度要求的个数：');
disp(ICT2);
for k=1:n
V(k)=sqrt(e(k)^2+f(k)^2);
sida(k)=atan(f(k)./e(k))*180./pi;
E(k)=e(k)+f(k)*j;
end
%=============== 计算各输出量 ===========================
disp('各节点的实际电压标幺值E为(节点号从小到大排列)：');
disp(E);
EE=E*UB;
disp(EE);
disp('-----------------------------------------------------'); disp('各节点的电压大小V为(节点号从小到大排列)：');
disp(V);
VV=V*UB;
disp(VV);
disp('-----------------------------------------------------'); disp('各节点的电压相角sida为(节点号从小到大排列)：');
disp(sida);
for p=1:n
C(p)=0;
for q=1:n
C(p)=C(p)+conj(Y(p,q))*conj(E(q));
end
S(p)=E(p)*C(p);
end
disp('各节点的功率S为(节点号从小到大排列)：');
disp(S);
disp('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~');
SS=S*SB;
disp(SS);
disp('-----------------------------------------------------');
disp('各条支路的首端功率Si为(顺序同您输入B1时一致)：');
for i=1:nl
p=B1(i,1);q=B1(i,2);
if B1(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);
else
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);
end
disp(Si(p,q));
SSi(p,q)=Si(p,q)*SB;
ZF=['S(',num2str(p),',',num2str(q),')=',num2str(SSi(p,q))];
disp(ZF);
%disp(SSi(p,q));
disp('-----------------------------------------------------'); end
disp('各条支路的末端功率Sj为(顺序同您输入B1时一致)：');
for i=1:nl
p=B1(i,1);q=B1(i,2);
if B1(i,6)==0
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(q,p);
else
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(q,p);
end
disp(Sj(q,p));
SSj(q,p)=Sj(q,p)*SB;
ZF=['S(',num2str(q),',',num2str(p),')=',num2str(SSj(q,p))];
disp(ZF);
%disp(SSj(q,p));
disp('-----------------------------------------------------'); end
disp('各条支路的功率损耗DS为(顺序同您输入B1时一致)：');
for i=1:nl
p=B1(i,1);q=B1(i,2);
DS(i)=Si(p,q)+Sj(q,p);
disp(DS(i));
DDS(i)=DS(i)*SB;
ZF=['DS(',num2str(p),',',num2str(q),')=',num2str(DDS(i))];
disp(ZF);
%disp(DDS(i));
disp('-----------------------------------------------------'); end
figure(1);
subplot(2,2,1);
plot(V);
xlabel('节点号');ylabel('电压标幺值');
grid on;
subplot(2,2,2);
plot(sida);
xlabel('节点号');ylabel('电压角度');
grid on;
subplot(2,2,3);
bar(S);
ylabel('节点注入无功');
grid on;
subplot(2,2,4);
bar(Siz);
ylabel('支路首端无功');
grid on;
A—2最小负荷时的潮流程序
%本程序的功能是用牛顿——拉夫逊法进行潮流计算
% B1矩阵：1、支路首端号；2、末端号；3、支路阻抗；4、支路对地电纳
% 5、支路的变比；6、支路首端处于K侧为1，1侧为0
% B2矩阵：1、该节点发电机功率；2、该节点负荷功率；3、节点电压初始值% 4、PV节点电压V的给定值；5、节点所接的无功补偿设备的容量% 6、节点分类标号
clear;
n=12;%input('请输入节点数：n=');
nl=14;%input('请输入支路数：nl=');
isb=1;%input('请输入平衡母线节点号：isb=');
pr=0.00001;%input('请输入误差精度：pr=');
B1=[1 2 3.3+11.48i 0 1 0;
1 4 3.87+15.375i 0 1 0;
1 6 5.8+18.8i 0 1 0;
2 3 2.39+34.95i 0 1.02 1;
4 5 1.59+23.3i 0 1.02 1;
2 4 4.67+18.34i 0 1 0;
6 7 2.39+24.35i 0 1.02 1;
8 9 2.39+24.35i 0 1.02 1;
10 11 2.39+24.35i 0 1.02 1;
8 10 1.7+11.2i 0 1 0;
4 12 4.2+13.9i 0 1 0;
6 12 4.2+14.64i 0 1 0;
8 12 2.2+13.9i 0 1 0;
10 12 2.4+15.2i 0 1 0];
%input('请输入由支路参数形成的矩阵： B1=');
B2=[0 0 231 231 0 1;
0 0 220 0 0 2;
0 80+60i 220 0 0 2;
0 0 220 0 0 2;
0 80+60i 220 0 0 2;
0 0 220 0 0 2;
0 120+90i 220 0 0 2;
0 0 220 0 0 2;
0 120+90i 220 0 0 2;
0 0 220 0 0 2;
0 120+90i 220 0 0 2;
360 0 231 231 0 3];%input('请输入各节点参数形成的矩阵： B2='); Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1= zeros(nl);
%－－－－－－－修改部分－－－－－－－－－－－－
ym=1;
SB=100;UB=220;
%ym=input('您输入的参数是标么值？（若不是则输入一个不为零的数值）'); if ym~=0
%SB=input('请输入功率基准值:SB=');
%UB=input('请输入电压基准值:UB=');
YB=SB./UB./UB;
BB1=B1;
BB2=B2;
for i=1:nl
B1(i,3)=B1(i,3)*YB;
B1(i,4)=B1(i,4)./YB;
end
disp('B1矩阵B1=');
disp(B1)
for i=1:n
B2(i,1)=B2(i,1)./SB;
B2(i,2)=B2(i,2)./SB;
B2(i,3)=B2(i,3)./UB;
B2(i,4)=B2(i,4)./UB;
B2(i,5)=B2(i,5)./SB;
end
disp('B2矩阵B2=');
disp(B2)
end
% % %---------------------------------------------------
for i=1:nl %支路数
if B1(i,6)==0 %左节点处于低压侧
p=B1(i,1);q=B1(i,2);
else
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)./2; %对角元K侧
Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2; %对角元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));
V(i)=B2(i,4); %PV节点电压给定模值
end
for i=1:n %给定各节点注入功率
S(i)=B2(i,1)-B2(i,2); %i节点注入功率SG-SL
B(i,i)=B(i,i)+B2(i,5); %i节点无功补偿量
end
%=================================================================== P=real(S);Q=imag(S);
ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0;
while IT2~=0
IT2=0;a=a+1;
for i=1:n
if i~=isb %非平衡节点
C(i)=0;D(i)=0;
for j1=1:n
C(i)=C(i)+G(i,j1)*e(j1)-B(i,j1)*f(j1);%Σ(Gij*ej-Bij*fj)
D(i)=D(i)+G(i,j1)*f(j1)+B(i,j1)*e(j1);%Σ(Gij*fj+Bij*ej)
end
P1=C(i)*e(i)+f(i)*D(i);%节点功率P计算eiΣ(Gij*ej-Bij*fj)+fiΣ(Gij*fj+Bij*ej)
Q1=C(i)*f(i)-e(i)*D(i);%节点功率Q计算fiΣ(Gij*ej-Bij*fj)-eiΣ(Gij*fj+Bij*ej)
%求P',Q'
V2=e(i)^2+f(i)^2; %电压模平方
%========= 以下针对非PV节点来求取功率差及Jacobi矩阵元素=========
if B2(i,6)~=3 %非PV节点
DP=P(i)-P1; %节点有功功率差
DQ=Q(i)-Q1; %节点无功功率差
%=============== 以上为除平衡节点外其它节点的功率计算=================
%================= 求取Jacobi矩阵 ===================
for j1=1:n
if j1~=isb&j1~=i %非平衡节点&非对角元
X1=-G(i,j1)*e(i)-B(i,j1)*f(i); %
dP/de=-dQ/df
X2=B(i,j1)*e(i)-G(i,j1)*f(i); % dP/df=dQ/de X3=X2; % X2=dp/df X3=dQ/de
X4=-X1; % X1=dP/de X4=dQ/df
p=2*i-1;q=2*j1-1;
J(p,q)=X3;J(p,N)=DQ;m=p+1;
J(m,q)=X1;J(m,N)=DP;q=q+1;
J(p,q)=X4;J(m,q)=X2;
elseif j1==i&j1~=isb %非平衡节点&对角元
X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);% dP/de
X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);% dP/df
X3=D(i)+B(i,i)*e(i)-G(i,i)*f(i); % dQ/de
X4=-C(i)+G(i,i)*e(i)+B(i,i)*f(i);% dQ/df
p=2*i-1;q=2*j1-1;J(p,q)=X3;J(p,N)=DQ;%扩展列△Q
m=p+1;
J(m,q)=X1;q=q+1;J(p,q)=X4;J(m,N)=DP;%扩展列△P J(m,q)=X2;
end
end
else
%=============== 下面是针对PV节点来求取Jacobi矩阵的元素
===========
DP=P(i)-P1; % PV节点有功误差
DV=V(i)^2-V2; % PV节点电压误差
for j1=1:n
if j1~=isb&j1~=i %非平衡节点&非对角元
X1=-G(i,j1)*e(i)-B(i,j1)*f(i); % dP/de
X2=B(i,j1)*e(i)-G(i,j1)*f(i); % dP/df
X5=0;X6=0;
p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV;
m=p+1;
J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6;
J(m,q)=X2;
elseif j1==i&j1~=isb %非平衡节点&对角元
X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);% dP/de
X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);% dP/df
X5=-2*e(i);
X6=-2*f(i);
p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV;
m=p+1;
J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6;
J(m,q)=X2;
end
end
end
end
end
%========= 以上为求雅可比矩阵的各个元素 =====================
for k=3:N0 % N0=2*n （从第三行开始，第一、二行是平衡节点）
k1=k+1;N1=N; % N=N0+1 即 N=2*n+1扩展列△P、△Q for k2=k1:N1 % 扩展列△P、△Q
J(k,k2)=J(k,k2)./J(k,k); % 非对角元规格化
end
J(k,k)=1; % 对角元规格化
if k~=3 % 不是第三行
%============================================================
k4=k-1;
for k3=3:k4 % 用k3行从第三行开始到当前行前的k4行消去
for k2=k1:N1 % k3行后各行下三角元素
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算
end
J(k3,k)=0;
end
if k==N0
break;
end
%==========================================
for k3=k1:N0
for k2=k1:N1
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算 end
J(k3,k)=0;
end
else
for k3=k1:N0
for k2=k1:N1
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算 end
J(k3,k)=0;
end
end
end
%====上面是用线性变换方式将Jacobi矩阵化成单位矩阵===== for k=3:2:N0-1
L=(k+1)./2;
e(L)=e(L)-J(k,N); %修改节点电压实部
k1=k+1;
f(L)=f(L)-J(k1,N); %修改节点电压虚部
end
%------修改节点电压-----------
for k=3:N0
DET=abs(J(k,N));
if DET>=pr %电压偏差量是否满足要求
IT2=IT2+1; %不满足要求的节点数加1
end
end
ICT2(a)=IT2;
ICT1=ICT1+1;
end
%用高斯消去法解"w=-J*V"
disp('迭代次数：');
disp(ICT1);
disp('没有达到精度要求的个数：');
disp(ICT2);
for k=1:n
V(k)=sqrt(e(k)^2+f(k)^2);
sida(k)=atan(f(k)./e(k))*180./pi;
E(k)=e(k)+f(k)*j;
end
%=============== 计算各输出量 ===========================
disp('各节点的实际电压标幺值E为(节点号从小到大排列)：');
disp(E);
EE=E*UB;
disp(EE);
disp('-----------------------------------------------------');
disp('各节点的电压大小V为(节点号从小到大排列)：');
disp(V);
VV=V*UB;
disp(VV);
disp('-----------------------------------------------------');
disp('各节点的电压相角sida为(节点号从小到大排列)：');
disp(sida);
for p=1:n
C(p)=0;
for q=1:n
C(p)=C(p)+conj(Y(p,q))*conj(E(q));
end
S(p)=E(p)*C(p);
end
disp('各节点的功率S为(节点号从小到大排列)：');
disp(S);
disp('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~');
SS=S*SB;
disp(SS);
disp('-----------------------------------------------------');
disp('各条支路的首端功率Si为(顺序同您输入B1时一致)：');
for i=1:nl
p=B1(i,1);q=B1(i,2);
if B1(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);
else
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);
end
disp(Si(p,q));
SSi(p,q)=Si(p,q)*SB;
ZF=['S(',num2str(p),',',num2str(q),')=',num2str(SSi(p,q))];
disp(ZF);
%disp(SSi(p,q));
disp('-----------------------------------------------------'); end
disp('各条支路的末端功率Sj为(顺序同您输入B1时一致)：');
for i=1:nl
p=B1(i,1);q=B1(i,2);
if B1(i,6)==0
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(q,p);
else
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(q,p);
end
disp(Sj(q,p));
SSj(q,p)=Sj(q,p)*SB;
ZF=['S(',num2str(q),',',num2str(p),')=',num2str(SSj(q,p))];
disp(ZF);
%disp(SSj(q,p));
disp('-----------------------------------------------------'); end
disp('各条支路的功率损耗DS为(顺序同您输入B1时一致)：');
for i=1:nl
p=B1(i,1);q=B1(i,2);
DS(i)=Si(p,q)+Sj(q,p);
disp(DS(i));
DDS(i)=DS(i)*SB;
ZF=['DS(',num2str(p),',',num2str(q),')=',num2str(DDS(i))];
disp(ZF);
%disp(DDS(i));
disp('-----------------------------------------------------'); end
figure(1);
subplot(2,2,1);
plot(V);
xlabel('节点号');ylabel('电压标幺值');
grid on;
subplot(2,2,2);
plot(sida);
xlabel('节点号');ylabel('电压角度');
grid on;
subplot(2,2,3);
bar(S);
ylabel('节点注入无功');
grid on;
subplot(2,2,4);
bar(Siz);
ylabel('支路首端无功');
grid on;
A—3 最大负荷时发电厂B与变电站3断线时的潮流程序
%本程序的功能是用牛顿——拉夫逊法进行潮流计算
% B1矩阵：1、支路首端号；2、末端号；3、支路阻抗；4、支路对地电纳
% 5、支路的变比；6、支路首端处于K侧为1，1侧为0
% B2矩阵：1、该节点发电机功率；2、该节点负荷功率；3、节点电压初始值% 4、PV节点电压V的给定值；5、节点所接的无功补偿设备的容量% 6、节点分类标号
clear;
n=12;%input('请输入节点数：n=');
nl=13;%input('请输入支路数：nl=');
isb=1;%input('请输入平衡母线节点号：isb=');
pr=0.00001;%input('请输入误差精度：pr=');
B1=[1 2 3.3+11.48i 0 1 0;
1 4 3.87+15.375i 0 1 0;
1 6 5.8+18.8i 0 1 0;
2 3 2.39+34.95i 0 1.02 1;
4 5 1.59+23.3i 0 1.02 1;
2 4 4.67+18.34i 0 1 0;
6 7 2.39+24.35i 0 1.02 1;
8 9 2.39+24.35i 0 1.02 1;
10 11 2.39+24.35i 0 1.02 1;
8 10 1.7+11.2i 0 1 0;
6 12 4.2+14.64i 0 1 0;
8 12 2.2+13.9i 0 1 0;
10 12 2.4+15.2i 0 1 0];
%input('请输入由支路参数形成的矩阵： B1=');
B2=[0 0 231 231 0 1;
0 0 220 0 0 2;
0 100+75i 220 0 0 2;
0 0 220 0 0 2;
0 100+75i 220 0 0 2;
0 0 220 0 0 2;
0 150+112.5i 220 0 0 2;
0 0 220 0 0 2;
0 150+112.5i 220 0 0 2;
0 0 220 0 0 2;
0 150+112.5i 220 0 0 2;
450 0 231 231 0 3];%input('请输入各节点参数形成的矩阵： B2='); Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1= zeros(nl);
%－－－－－－－修改部分－－－－－－－－－－－－
ym=1;
SB=100;UB=220;
%ym=input('您输入的参数是标么值？（若不是则输入一个不为零的数值）'); if ym~=0
%SB=input('请输入功率基准值:SB=');
%UB=input('请输入电压基准值:UB=');
YB=SB./UB./UB;
BB1=B1;
BB2=B2;
for i=1:nl
B1(i,3)=B1(i,3)*YB;
B1(i,4)=B1(i,4)./YB;
end
disp('B1矩阵B1=');
disp(B1)
for i=1:n
B2(i,1)=B2(i,1)./SB;
B2(i,2)=B2(i,2)./SB;
B2(i,3)=B2(i,3)./UB;
B2(i,4)=B2(i,4)./UB;
B2(i,5)=B2(i,5)./SB;
end
disp('B2矩阵B2=');
disp(B2)
end
% % %---------------------------------------------------
for i=1:nl %支路数
if B1(i,6)==0 %左节点处于低压侧
p=B1(i,1);q=B1(i,2);
else
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)./2; %对角元K侧
Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2; %对角元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));
V(i)=B2(i,4); %PV节点电压给定模值
end
for i=1:n %给定各节点注入功率
S(i)=B2(i,1)-B2(i,2); %i节点注入功率SG-SL
B(i,i)=B(i,i)+B2(i,5); %i节点无功补偿量
end
%=================================================================== P=real(S);Q=imag(S);
ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0;
while IT2~=0
IT2=0;a=a+1;
for i=1:n
if i~=isb %非平衡节点
C(i)=0;D(i)=0;
for j1=1:n
C(i)=C(i)+G(i,j1)*e(j1)-B(i,j1)*f(j1);%Σ(Gij*ej-Bij*fj)
D(i)=D(i)+G(i,j1)*f(j1)+B(i,j1)*e(j1);%Σ(Gij*fj+Bij*ej)
end
P1=C(i)*e(i)+f(i)*D(i);%节点功率P计算eiΣ(Gij*ej-Bij*fj)+fiΣ(Gij*fj+Bij*ej)
Q1=C(i)*f(i)-e(i)*D(i);%节点功率Q计算fiΣ(Gij*ej-Bij*fj)-eiΣ(Gij*fj+Bij*ej)
%求P',Q'
V2=e(i)^2+f(i)^2; %电压模平方
%========= 以下针对非PV节点来求取功率差及Jacobi矩阵元素=========
if B2(i,6)~=3 %非PV节点
DP=P(i)-P1; %节点有功功率差
DQ=Q(i)-Q1; %节点无功功率差
%=============== 以上为除平衡节点外其它节点的功率计算=================
%================= 求取Jacobi矩阵 ===================
for j1=1:n
if j1~=isb&j1~=i %非平衡节点&非对角元
X1=-G(i,j1)*e(i)-B(i,j1)*f(i); %
dP/de=-dQ/df
X2=B(i,j1)*e(i)-G(i,j1)*f(i); % dP/df=dQ/de X3=X2; % X2=dp/df X3=dQ/de
X4=-X1; % X1=dP/de X4=dQ/df
p=2*i-1;q=2*j1-1;
J(p,q)=X3;J(p,N)=DQ;m=p+1;
J(m,q)=X1;J(m,N)=DP;q=q+1;
J(p,q)=X4;J(m,q)=X2;
elseif j1==i&j1~=isb %非平衡节点&对角元
X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);% dP/de
X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);% dP/df
X3=D(i)+B(i,i)*e(i)-G(i,i)*f(i); % dQ/de
X4=-C(i)+G(i,i)*e(i)+B(i,i)*f(i);% dQ/df
p=2*i-1;q=2*j1-1;J(p,q)=X3;J(p,N)=DQ;%扩展列△Q
m=p+1;
J(m,q)=X1;q=q+1;J(p,q)=X4;J(m,N)=DP;%扩展列△P J(m,q)=X2;
end
end
else
%=============== 下面是针对PV节点来求取Jacobi矩阵的元素===========
DP=P(i)-P1; % PV节点有功误差
DV=V(i)^2-V2; % PV节点电压误差
for j1=1:n
if j1~=isb&j1~=i %非平衡节点&非对角元
X1=-G(i,j1)*e(i)-B(i,j1)*f(i); % dP/de
X2=B(i,j1)*e(i)-G(i,j1)*f(i); % dP/df
X5=0;X6=0;
p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV;
m=p+1;
J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6;
J(m,q)=X2;
elseif j1==i&j1~=isb %非平衡节点&对角元
X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);% dP/de
X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);% dP/df
X5=-2*e(i);
X6=-2*f(i);
p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV;
m=p+1;
J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6;
J(m,q)=X2;
end
end
end
end
end
%========= 以上为求雅可比矩阵的各个元素 =====================
for k=3:N0 % N0=2*n （从第三行开始，第一、二行是平衡节点）
k1=k+1;N1=N; % N=N0+1 即 N=2*n+1扩展列△P、△Q for k2=k1:N1 % 扩展列△P、△Q
J(k,k2)=J(k,k2)./J(k,k); % 非对角元规格化
end
J(k,k)=1; % 对角元规格化
if k~=3 % 不是第三行
%============================================================
k4=k-1;
for k3=3:k4 % 用k3行从第三行开始到当前行前的k4行消去
for k2=k1:N1 % k3行后各行下三角元素
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算
end
J(k3,k)=0;
end
if k==N0
break;
end
%==========================================
for k3=k1:N0
for k2=k1:N1
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算
end
J(k3,k)=0;
end
else
for k3=k1:N0
for k2=k1:N1
J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算
电力系统计算程序设计(包含源程序)
电力系统计算程序设计（包含matlab源程序）
广西大学电气工程学院 2007年1月 第一章原始数据 电力系统原始数据是电力系统计算的基础。电力系统每个计算程序都要求输入一定的原始数据，这些数据可以反映电力网络结构、电力系统正常运行条件、电力系统各元件参数和特性曲线。不同的计算程序需要不用的原始数据。 第一节电力网络的描述 电力网络是由输电线路、电力变压器、电容器和电抗器等元件组成。这些元件一般用集中参数的电阻、电抗和电容表示。为了表示电力网络中各元件是怎样互相连接的，通常要对网络节点进行编号。电力网络的结构和参数由电力网络中各支路的特性来描述。 1.1.1 线路参数 在电力系统程序设计中，线路参数一般采用线路的Π型数学模型，即线路用节点间的阻抗和节点对地容性电纳来表示，由于线路的对地电导很小，一般可忽略不计。其等价回路如下： r+jx -jb/2 对于线路参数的数据文件格式一般可写为： 线路参数（序号，节点i，节点j，r，x，b/2） 1.1.2 变压器参数
在电力系统程序设计中，变压器参数一般采用Π型等值变压器模型，这是一种可等值地体现变压器电压变换功能的模型。在多电压级网络计算中采用这种变压器模型后，就可不必进行参数和变量的归算。双绕组变压器的等值回路如下： k Z T k:1 Z T (a)接入理想变压器后的等值电路(b) 等值电路以导纳表示 (c) 等值电路以导纳表示 三绕组变压器的等值回路如下： 综合所述，三绕组变压器的等值电路可以用两个双绕组变压器的等值电
路来表示。因此，对于变压器参数的数据文件格式一般可写为：变压器参数（序号，节点i，节点j，r，x，k0） 其中，k0表示变压器变比。 1.1.3对地支路参数 对地支路参数一般以导纳形式表示，其等价回路如下： i g-jb 对地支路参数的数据文件格式一般可写为： 接地支路参数（序号，节点i，g i，b i） 第二节电力系统运行条件数据 电力系统运行条件数据包括发电机（含调相机）所连接的节点号、有功与无功功率；负荷所连接的节点号、有功与无功功率；PV节点与给定电压值；平衡节点的节点号与给定电压值。 1.2.1节点功率参数 电力系统中有流入流出功率的称为功率节点，有流入功率的称发电节点，一般为各发电站、枢纽变电站等节点；有流出功率的称负荷节点。对于电力系统稳态计算来说，功率节点都用有功功率P和无功功率Q来简单表示。其等价回路如下： Q G P G P L Q L 节点功率参数的数据文件格式一般可写为：
计算机潮流计算程序代码
G(5,5)=6.250;B(5,5)=-18.750;G(5,2)=-5.00;B(5,2)=15.000; G(5,3)=-1.250;B(5,3)=3.750;G(5,4)=0;B(5,4)=0; G(5,1)=0;B(5,1)=0; G(2,5)=-5.00;B(2,5)=15.000;G(2,2)=10.834;B(2,2)=-32.500;G(2,3)=-1.677;B(2,3)=5.00; G(2,4)=-1.667;B(2,4)=5.00;G(2,1)=-2.500;B(2,1)=7.500; G(3,5)=-1.250;B(3,5)=3.750;G(3,2)=-1.667;B(3,2)=5.00;G(3,3)=12.917;B(3,3)=-38.750; G(3,4)=-10.000;B(3,4)=30.000;G(3,1)=0;B(3,1)=0; G(4,5)=0;B(4,5)=0;G(4,2)=-1.667;B(4,2)=5.000;G(4,3)=-10.000;B(4,3)=30.000; G(4,4)=12.917;B(4,4)=-38.750;G(4,1)=-1.250;B(4,1)=3.750; G(1,5)=0;B(1,5)=0;G(1,2)=-2.500;B(1,2)=7.500;G(1,3)=0;B(1,3)=0; G(1,4)=-1.250;B(1,4)=3.750;G(1,1)=3.750;B(1,1)=-11.250; Y=G+j*B delt(1)=0; delt(2)=0; delt(3)=0; delt(4)=0; u(1)=1.0; u(2)=1.0; u(3)=1.0; u(4)=1.0; ps(2)=0.2;qs(2)=0.20;ps(3)=-0.45;qs(3)=-0.15; ps(4)=-0.4;qs(4)=-0.05; ps(1)=-0.6;qs(1)=-0.1; k=1;precision=1; N1=4; while precision>0.0001 delt(5)=0;u(5)=1.06; for m=1:N1 for n=1:N1+1 pt(n)=u(m)*u(n)*(G(m,n)*cos(delt(m)-delt(n))+B(m,n)*sin(delt(m)-delt(n))); qt(n)=u(m)*u(n)*(G(m,n)*sin(delt(m)-delt(n))-B(m,n)*cos(delt(m)-delt(n))); end pi(m)=sum(pt);qi(m)=sum(qt); dp(m)=ps(m)-pi(m); dq(m)=qs(m)-qi(m); end for m=1:N1 for n=1:N1 if m==n H(m,m)=-qi(m)-u(m)^2*B(m,m); N(m,m)=pi(m)+u(m)^2*G(m,m); J(m,m)=pi(m)-u(m)^2*G(m,m); L(m,m)=qi(m)-u(m)^2*B(m,m); JJ(2*m-1,2*m-1)=H(m,m); JJ(2*m-1,2*m)=N(m,m); JJ(2*m,2*m-1)=J(m,m); JJ(2*m,2*m)=L(m,m); else H(m,n)=u(m)*u(n)*(G(m,n)*sin(delt(m)-delt(n))-B(m,n)*cos(delt(m)-delt(n))); J(m,n)=-u(m)*u(n)*(G(m,n)*cos(delt(m)-delt(n))+B(m,n)*sin(delt(m)-delt(n))); N(m,n)=-J(m,n);L(m,n)=H(m,n); JJ(2*m-1,2*n-1)=H(m,n);JJ(2*m-1,2*n)=N(m,n); JJ(2*m,2*n-1)=J(m,n); JJ(2*m,2*n)=L(m,n); end
matlab电力系统潮流计算
华中科技大学 信息工程学院课程设计报告书题目: 电力系统潮流计算 专业：电气工程及其自动化 班级： 学号： 学生姓名： 指导教师： 2015年 11 月 10 日
2015年11月12日
信息工程学院课程设计成绩评定表
摘要 电力系统稳态分析包括潮流计算和静态安全分析。本文主要运用的事潮流计算，潮流计算是电力网络设计与运行中最基本的运算，对电力网络的各种设计方案及各种运行方式进行潮流计算，可以得到各种电网各节点的电压，并求得网络的潮流及网络中的各元件的电力损耗，进而求得电能损耗。本位就是运用潮流计算具体分析，并有MATLAB仿真。 关键词：电力系统潮流计算 MATLAB仿真
Abstract Electric power system steady flow calculation and analysis of the static safety analysis. This paper, by means of the calculation, flow calculation is the trend of the power network design and operation of the most basic operations of electric power network, various design scheme and the operation ways to tide computation, can get all kinds of each node of the power grid voltage and seek the trend of the network and the network of the components of the power loss, and getting electric power. The standard is to use the power flow calculation and analysis, the specific have MATLAB simulation. Key words: Power system; Flow calculation; MATLAB simulation
c语言程序设计课程计算器设计报告
课程设计说明书 题目计算器程序设计 起讫日期 2006 年 7月 3日至 2006 年 8月 6日 所在院系软件学院 专业机械+软件班级 04-2 学生姓名偶偶哦学号 指导教师 2006年 8 月日
摘要 当今社会是信息社会，科技经济高速发展的社会！为了更方便人们的工作生活和加速人们处理信息的速度，计算器应运而生。由于它体积小巧，携带方便，价格便宜，构造简单等诸多的优点成为人们生活中的必备品！ 随着科技的发展计算器的种类变得更多，功能变得更强大，体积变得更小！电脑的出现改变人们的生活习惯，很多事情都可以电脑来完成！电脑的更大一个优点就是可以通过软件的应用无限的延伸电脑功能的外延！下面我们将用我们学习的c语言编写一个简易的计算器程序！实现简单的初步的计算功能！ 本程序的编写基础是Tubro 汉化版,它在tubro c的原有基础上实现了多汉字的支持方便了我们的使用。生成的程序可移植性强兼容性好稳定！现在只实现了加、减、乘、除、求幂、求模，求平方根，求Sin，求Cos，求Log10，以及一个时钟原代码。这个系统是基于软件发展的生命周期来研制的，它可以直接输入数学表达式，不需要任何转换，就可以直接输出数学四则运算的结果。但是，每次只能运算一个表达式。不能运算多个表达式。在程序里面在添加一组选择函数即可。本论文主要介绍了本课题的开发背景，开发的过程和所要完成的功能。重点的说明了系统设计思想，设计的步骤、难点技术和解决方案。 关键词：C语言 Tubro c 汉化版计算器时钟
目录 第一章综述 (1) 1．1 课题的现实意义 (1) 1．2 软件环境 (1) 1．3 硬件环境 (1) 第二章系统设计流程图 (2) 2．1 系统流程图 (2) 2．2 主要功能表 (2) 第三章系统分析和设计 (3) 3．1 图形的绘制和输出 (3) 3．2 文本的输出显示 (3) 3．3 计算函数的调用 (4) 3．4 程序的运行和退出 (5) 第四章系统测试 (6) 4．1 系统测试 (6) 4．2 调试 (6) 4．3 错误原因分析一 (6) 4．4 错误原因分析二 (6) 第五章用户使用说明书 (8)
潮流计算的计算机算法
第四章潮流计算的计算机算法 第一节概述 潮流计算是电力系统最基本、最常用的计算。根据系统给定的运行条件、网络接线及元件参数，通过潮流计算可以确定各母线的电压（幅值及相角），各元件中流过的功率、整个系统的功率损耗等。潮流计算是实现电力系统安全经济发供电的必要手段和重要工作环节。因此潮流计算在电力系统的规划设计、生产运行、调度管理及科学研究中都有着广泛的应用。 电力系统潮流计算分为离线潮流计算和在线潮流计算。前者主要用于系统规划设计和安排系统的运行方式，后者则用于正在运行系统的经常监视及实时控制。本章主要讨论离线潮流计算问题，它的基本算法同样适用于在线潮流计算。 潮流计算在数学上是多元非线性方程组的求解问题，求解的方法有很多种。自从五十年代计算机应用于电力系统以来，当时求解潮流的方法是以节点导纳矩阵为基础的逐次代入法（导纳法），后来为解决导纳法的收敛性较差的问题，出现了以阻抗矩阵为基础的逐次代入法（阻抗法）。到六十年代，针对阻抗法占用计算机内存大的问题又出现了分块阻抗法及牛顿－拉夫逊（Newton－Raphson）法。Newton —Raphson法是数学上解非线形方程式的有效方法，有较好的收敛性。将N－R法用于潮流计算是以导纳矩阵为基础的，由于利用了导纳矩阵的对称性、稀疏性及节点编号顺序优化等技巧，使N－R法在收敛性、占用内存、计算速度方面的优点都超过了阻抗法，成为六十年代末期以后普遍采用的方法。同时国内外广泛研究了诸如非线形规划法、直流法、交流法等各种不同的潮流计算方法。七十年代以来，又涌现出了更新的潮流计算方法。其中有1974年由B、Stott、O、Alsac 提出的快速分解法以及1978年由岩本伸一等提出的保留非线性的高 129
潮流计算matlab程序
clear; %各节点参数：节点编号，类型，电压幅值，电压相位，注入有功，注入无功%类型：1＝PQ节点，2＝PV节点，3＝平衡节点 %本程序中将最后一个节点设为平衡节点 R_1=[1 1 1.0 0 0.2 0.2j; 2 1 1.0 0 -0.45 -0.15j; 3 1 1.0 0 -0.45 -0.05j; 4 1 1.0 0 -0.6 -0.1j; 5 3 1.0 0 0 0]; %支路号首端节点末端节点支路导纳 R_2=[1 5 2 1.25-3.75j; 2 2 3 10.00-30.00j; 3 3 4 1.25-3.75j; 4 1 4 2.50-7.50j; 5 1 5 5.00-15.00j; 6 1 2 1.667-5.00j]; n=5;L=6;%需要改变的到此为止 i=0;j=0;a=0;precision=1;k=0; Y=zeros(n,n);u=zeros(1,n);delt=zeros(1,n);P=zeros(1,n);Q=zeros(1,n); G=[];B=[];PP=[];uu=[];U=[];dp=[];dq=[]; for a=1:L i=R_2(a,2); j=R_2(a,3); Y(i,j)=-R_2(a,4); Y(j,i)=Y(i,j); end for a=1:n for b=1:n if a~=b Y(a,a)=Y(a,a)+Y(a,b); end end end for i=1:n for j=1:n if i==j Y(i,j)=-Y(i,j); end end end Y %形成导纳矩阵 for i=1:n for j=1:n G(i,j)=real(Y(i,j));
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 %给定各节点初始电压的实部和虚部
C语言课程设计 简单计算器程序
课程设计名称：C语言课程设计课程设计题目：简单计算器程序
目录 第1章需求分析 (1) 1.1设计要求 (1) 1.2任务 (1) 第2章总体设计 (2) 2.1设计简介及设计方案论述 (2) 2.2功能模块层次图 (2) 第3章详细设计 (3) 3.3由（后缀）逆波兰表达式计算中缀表达式原理 (8) 3.3.1算法描述 (8) 第4章调试分析 (10) 4.1程序设计中所遇到的错误及犯错的原因 (10) 4.2错误的解决方法 (10) 第5章用户手册 (11) 总结 (15) 参考文献 (16) 附录（程序清单） (17)
第1章需求分析 1.1 设计要求 （1）用 C 语言数据结构实现程序设计； （2）利用结构体、栈、进行相关信息处理； （2）系统的各个功能模块要求用函数的形式实现； （4）界面简单，可操作性高。 1.2任务 （1）定义一个结构体类型数组，输入0~9 及+、--、*等符号的信息，将其信息存储起来； （2）输入简单的加减乘除算术计算式，并在屏幕上显示逆波兰（后缀式）表达式和计算结果； （3）编写代码； （4）程序分析与调试。 说明： 本课程设计将实现一个简单计算器。在功能上尽量模仿windows 的计算器。系统界面不做牵制要求。该程序能实现标准型中+、-、*、/、（、）、.、的混合运算表达式（一般意义上的中缀表达式），将其转换成逆序波兰表达式（后缀表达式）并计算输出结果。在进行运算后可以选择继续运算或者结束当前运算。即时准确地获得需要的计算的结果，充分降低了数字计算的难度和节约了时间，对人们的生活有一定的帮助。
第2章 总体设计 2.1设计简介及设计方案论述 逆波兰表达式又叫做后缀表达式。在通常的表达式中，二元运算符总是置于与之相 关的两个运算对象之间，所以，这种表示法也称为中缀表达式。波兰逻辑学家 J.Lukasiewicz 于 1929 年提出了另一种表示表达式的方法。按此方法，每一运算符都置 于其运算对象之后，故称为后缀表达式。 后缀表达式的优点是显而易见的， 编译器在处理时候按照从左至右的顺序读取逆波 兰表达式，遇到运算对象直接压入堆栈，遇到运算符就从堆栈提取后进的两个对象进行计算，这个过程正好符合了计算机计算的原理。后缀表达式比前缀表达式更加易于转换，并且它的最左面一定为数字，这一点在实 际编程的时候就会体会到它的好处了。 逆波兰表达式有一个更大的优点，就是拆括号，根据运算符的级别将中缀表达式转 换成逆波兰表达式后，运算顺序就已经替代了运算符的级别，这样也避免了括号提高运 算级别的特殊处理。 2.2功能模块层次图 将算术表达式转化为逆波兰表达式 计算逆波兰表达式的值 简单计算器 表 达 式 格 式 转 换 系统 求 值 计 算 系 统
用matlab电力系统潮流计算
题目：潮流计算与matlab 教学单位电气信息学院姓名 学号 年级 专业电气工程及其自动化指导教师 职称副教授
摘要 电力系统稳态分析包括潮流计算和静态安全分析。本文主要运用的事潮流计算，潮流计算是电力网络设计与运行中最基本的运算，对电力网络的各种设计方案及各种运行方式进行潮流计算，可以得到各种电网各节点的电压，并求得网络的潮流及网络中的各元件的电力损耗，进而求得电能损耗。本位就是运用潮流计算具体分析，并有MATLAB仿真。 关键词：电力系统潮流计算 MATLAB Abstract Electric power system steady flow calculation and analysis of the static safety analysis. This paper, by means of the calculation, flow calculation is the trend of the power network design and operation of the most basic operations of electric power network, various design scheme and the operation ways to tide computation, can get all kinds of each node of the power grid voltage and seek the trend of the network and the network of the components of the power loss, and getting electric power. The standard is to use the power flow calculation and analysis, the specific have MATLAB simulation. Key words: Power system; Flow calculation; MATLAB simulation
潮流计算(matlab)实例计算
潮流例题：根据给定的参数或工程具体要求（如图），收集和查阅资料；学习相关软件（软件自选：本设计选择Matlab进行设计）。 2.在给定的电力网络上画出等值电路图。 3.运用计算机进行潮流计算。 4.编写设计说明书。 一、设计原理 1．牛顿-拉夫逊原理 牛顿迭代法是取x0 之后，在这个基础上，找到比x0 更接近的方程的跟，一步一步迭代，从而找到更接近方程根的近似跟。牛顿迭代法是求方程根的重要方法之一，其最大优点是在方程f(x) = 0 的单根附近具有平方收敛，而且该法还可以用来求方程的重根、复根。电力系统潮流计算，一般来说，各个母线所供负荷的功率是已知的，各个节点电压是未知的（平衡节点外）可以根据网络结构形成节点导纳矩阵，然后由节点导纳矩阵列写功率方程，由于功率方程里功率是已知的，电压的幅值和相角是未知的，这样潮流计算的问题就转化为求解非线性方程组的问题了。为了便于用迭代法解方程组，需要将上述功率方程改写成功率平衡方程，并对功率平衡方程求偏导，得出对应的雅可比矩阵，给未知节点赋电压初值，一般为额定电压，将初值带入功率平衡方程，得到功率不平衡量，这样由功率不平衡量、雅可比矩阵、节点电压不平衡量（未知的）构成了误差方程，解误差方程，得到节点电压不平衡量，节点电压加上节点电压不平衡量构成新的节点电压初值，将新的初值带入原来的功率平衡方程，并重新形成雅可比矩阵，然后计算新
的电压不平衡量，这样不断迭代，不断修正，一般迭代三到五次就能收敛。 牛顿—拉夫逊迭代法的一般步骤： （1）形成各节点导纳矩阵Y。 （2）设个节点电压的初始值U和相角初始值e 还有迭代次数初值为0。 （3）计算各个节点的功率不平衡量。 （4）根据收敛条件判断是否满足，若不满足则向下进行。 （5）计算雅可比矩阵中的各元素。 （6）修正方程式个节点电压 （7）利用新值自第（3）步开始进入下一次迭代，直至达到精度退出循环。 （8）计算平衡节点输出功率和各线路功率 2．网络节点的优化 １）静态地按最少出线支路数编号 这种方法由称为静态优化法。在编号以前。首先统计电力网络个节点的出线支路数，然后，按出线支路数有少到多的节点顺序编号。当由n 个节点的出线支路相同时，则可以按任意次序对这n 个节点进行编号。这种编号方法的根据是导纳矩阵中，出线支路数最少的节点所对应的行中非零元素也２）动态地按增加出线支路数最少编号在上述的方法中，各节点的出线支路数是按原始网络统计出来的，在编号过程中认为固定不变的，事实上，在节点消去过程中，每消去一个节点以后，与该节点相连的各节点的出线支路数将发生变化(增加，减少或保持不变)。因此，如果每消去一个节点后，立即修正尚未编号节点的出线支路数，然后选其中支路数最少的一个节点进行编号，就可以预期得到更好的效果，动态按最少出线支路数编号方法的特点就是按出线最少原则编号时考虑了消去过程中各节点出线支路数目的变动情况。 3．MATLAB编程应用 Matlab 是“Matrix Laboratory”的缩写，主要包括：一般数值分析，矩阵运算、数字信号处理、建模、系统控制、优化和图形显示等应用程序。由于使用Matlab 编程运算与人进行科学计算的思路和表达方式完全一致，所以不像学习高级语言那样难于掌握，而且编程效率和计算效率极高，还可在计算机上直接输出结果和精美的图形拷贝，所以它的确为一高效的科研助手。 二、设计内容 1.设计流程图
C语言简易计算器的实现
目录 一.课程设计目的 (1) 二.设计环境 (1) 三.设计内容 (1) 四.设计说明 (2) 五.设计程序流程图 (2) 六.调试 (4) (1)错误原因分析一 (4) (2)语法错误 (5) (3)逻辑错误 (5) 七. 调试结果图 (6) 八. 结论与心得体会 (7) 九.附录 (8) 具体代码实现 (8) 十．参考文献 (18)
一.课程设计目的 1.通过一个学期的学习，我认为要学号C语言程序这门课程，不仅要认真阅读课本知识，更重要的是要通过上机实践来巩固我 们的知识，特别是学计算机专业的，我们更应该注重这一环节， 只有这样我们才能成为一个合格的计算机人才。通过这一个课程 设计，进一步来巩固所学的语句，如：循环，和分支结构的运用。还要熟悉四则运算和函数的算法。 2.通过这次课程设计扩展自己的知识面，课本上的东西是远 远不够的，可以通过上网或去图书馆查资料等方式得到一些新的 知识， 3.通过课程设计，加深对课程化设计思想的理解，能进行一 个系统功能分析，并设计一个合理的模块化结构，提高程序开发 能力。 二.设计环境 1.硬件：一台完整的电脑，包括键盘、鼠标，最小硬盘空间1GHz 2.软件：安装有Microsoft visual c++6.0 三.设计内容 以简易计算器为例，通过对简单应用软件计算器的设计，编制、调试，实现
简单的加，减，乘，除等运算，以学习应用MFC库类编写对话框的原理，加深对C++类的学习及应用。 （1）定义一个结构体类型数组，输入0~9及+、--、*等符号的信息，将其信息存入文件中； （2）输入简单的加减乘除算术计算式，并在屏幕上显示计算结果； （3）画出部分模块的流程图； （4）编写代码； （5）程序分析与调试。 四.设计说明 1)包含的功能有：加、减、乘、除运算，开方、平方等功能。 (2)计算器上数字0—9为一个控件数组，加、减、乘、除为一个控件数组，其余为单一的控件。 (3)输入的原始数据、运算中间数据和结果都显示在窗口顶部的同一个标签中。 (4)计算功能基本上是用系统内部函数。 (5)程序可以能自动判断输入数据的正确性，保证不出现多于一个小数点、以0开头等不正常现象。 (6)“CE”按钮可以清除所有已输入的数据从头计算 五.设计程序流程图
基于MATLAB的潮流计算源程序代码(优.选)
%*************************电力系统直角坐标系下的牛顿拉夫逊法潮流计算********** clear clc load E:\data\IEEE014_Node.txt Node=IEEE014_Node; weishu=size(Node); nnum=weishu(1,1); %节点总数 load E:\data\IEEE014_Branch.txt branch=IEEE014_Branch; bwei=size(branch); bnum=bwei(1,1); %支路总数 Y=(zeros(nnum)); Sj=100; %********************************节点导纳矩阵******************************* for m=1:bnum; s=branch(m,1); %首节点 e=branch(m,2); %末节点 R=branch(m,3); %支路电阻 X=branch(m,4); %支路电抗 B=branch(m,5); %支路对地电纳 k=branch(m,6); if k==0 %无变压器支路情形 Y(s,e)=-1/(R+j*X); %互导纳 Y(e,s)=Y(s,e); end if k~=0 %有变压器支路情形 Y(s,e)=-(1/((R+j*X)*k)); Y(e,s)=Y(s,e); Y(s,s)=-(1-k)/((R+j*X)*k^2); Y(e,e)=-(k-1)/((R+j*X)*k); %对地导纳 end Y(s,s)=Y(s,s)-j*B/2; Y(e,e)=Y(e,e)-j*B/2; %自导纳的计算情形 end for t=1:nnum; Y(t,t)=-sum(Y(t,:))+Node(t,12)+j*Node(t,13); %求支路自导纳 end G=real(Y); %电导 B=imag(Y); %电纳 %******************节点分类************************************* * pq=0; pv=0; blancenode=0; pqnode=zeros(1,nnum); pvnode=zeros(1,nnum); for m=1:nnum; if Node(m,2)==3 blancenode=m; %平衡节点编号 else if Node(m,2)==0 pq=pq+1; pqnode(1,pq)=m; %PQ 节点编号 else if Node(m,2)==2 pv=pv+1; pvnode(1,pv)=m; %PV 节点编号 end end end end %*****************************设置电压初值********************************** Uoriginal=zeros(1,nnum); %对各节点电压矩阵初始化 for n=1:nnum Uoriginal(1,n)=Node(n,9); %对各点电压赋初值 if Node(n,9)==0;
C语言实现计算器功能
实验一多功能计算器 一、问题描述 设计一个多功能计算器,可以完成基本的计算。 设计要求: 1、具备整型数据、浮点型数据的算术(加、减、乘、除)运算功能。依次输入第一个运算数、运算符(+,-,*,/)、第二个运算数,然后输出结果。结果可以作为下一个运算的第一运算数。按‘C’清屏,按‘R’返回菜单。 例如:输入:2 + 5 输出:7 2、实现单运算符表达式计算的功能。输入的操作数可以包含整数或浮点数。输入表达式如下: 例如:输入:2+5 输出:7 二、算法说明 1.数据结构说明(可以图示说明,也可以文字说明) 本程序主要根据选择菜单编写了六个自定义函数,用于在main()函数中调用,在main()中,用一个字符变量num1来记录下菜单选项的标号,根据num1的值来决定调用哪个函数。 程序要完成的功能及所要用到的函数如下:
下面就是整个程序的流程图:
2.算法说明(即函数说明) void suanshuyunsuan() //做算术运算时调用的函数 void suanshuyunsuan2() //选择继续做算术运算调用的函数,将上次运算的结果做为下次算术运算的第一个操作数//判断算术运算就是否继续 void panduan() //判断算术运算就是否继续 void biaodashiyunsuan() //单运算符表达式实现函数 void qingping() //清除屏幕 void fanhuicaidan() //显示菜单 三、测试结果(这部分需文字与图示结合) 1.第一组测试用例 (1)测试输入: 测试目的:测试算术运算的功能 结果输出:
(2)再一次输入:1测试目的:测试算术运算就是否能继续 结果输出: (3)这时输入:0 测试目的:退出算术运算 结果输出:
电力系统分析潮流计算课程序设计及其MATLAB程序设计
电力系统分析潮流计算程序设计报告题目：13节点配电网潮流计算 学院电气工程学院 专业班级 学生姓名 学号 班内序号 指导教师房大中 提交日期 2015年05月04日
目录 一、程序设计目的 (2) 二、程序设计要求 (4) 三、13节点配网潮流计算 (4) 3.1主要流程................................................................................................ 错误！未定义书签。 3.1.1第一步的前推公式如下（1-1）-（1-5）： .................................. 错误！未定义书签。 3.1.2第二步的回代公式如下（1-6）—（1-9）： ................................ 错误！未定义书签。 3.2配网前推后代潮流计算的原理 (7) 3.3配网前推后代潮流计算迭代过程 (8) 3.3计算原理 (9) 四、计算框图流程 (10) 五、确定前推回代支路次序....................................................................................... 错误！未定义书签。 六、前推回代计算输入文件 (11) 主程序： (11) 输入文件清单： (12) 计算结果： (13) 数据分析： (13) 七、配电网潮流计算的要点 (14) 八、自我总结 (14) 九、参考文献 (15) 附录一 MATLAB的简介 (15)
大学计算机c语言计算器源代码
C++语言编写。。 #include
char *base; }stack_char; stack_num S_num;//定义 stack_char S_char;//定义 char fu[18] = {'\n', ')', '+', '-', '*', '/', '%', '^', 'Q', 'L', 'C', 'S', 'T', 'c', 's', 't', '('}; int compare[1000];//表现出各运算符号的优先级 double shu[1000];//存储"数"的数组 double dai_result;//运算的结果，是为了处理M运算(简介函数里有M的定义) int biao = 0;//和dia_result一样，为了处理M运算 char line[SIZE];//输入的所要计算的表达式 void init()//初始化 { compare[fu[0]] = -2;//用数字的大小表现出符号的优先级 compare[fu[1]] = -1; compare[fu[2]] = 2; compare[fu[3]] = 2; compare[fu[4]] = 4; compare[fu[5]] = 4; compare[fu[6]] = 4;
潮流计算的计算机算法
高等电力系统分析 （潮流计算的计算机算法）PQ分解法潮流计算（IEEE14)
目录 一、MATLAB源程序 二、对支路参数（B1）、节点参数（B2）的说明 三、带入数据，运行结果
一、MATLAB源程序 clear close all n=input('请输入节点数:n='); n1=input('请输入支路数:n1='); isb=input('请输入平衡节点号:isb='); pr=input('请输入误差精度:pr='); B1=input('请输入支路参数:B1='); B2=input('请输入节点参数:B2='); n2=input('请输入PQ节点个数:n2='); Y=zeros(n); for i=1:n1 p=B1(i,1); q=B1(i,2); Y(p,q)=Y(p,q)-1/(B1(i,3)+B1(i,4)*1j); %非对角元 Y(q,p)=Y(p,q); Y(p,p)=Y(p,p)+1/(B1(i,3)+B1(i,4)*1j)+B1(i,6)*1j; %对角元 Y(q,q)=Y(q,q)+1/(B1(i,3)+B1(i,4)*1j)+B1(i,6)*1j; end disp('导纳矩阵Y='); disp(Y) %--------------------------------------------- %---------------下面是求P,Q,V,O矩阵--------------- V=zeros(1,n);O=zeros(1,n);P=zeros(1,n);Q=zeros(1,n); G=real(Y);B=imag(Y); for i=1:n P(i)=B2(i,3); Q(i)=B2(i,4); V(i)=B2(i,5); O(i)=B2(i,6); end B3=B(1:n-1,1:n-1); %不含平衡节点，由节点导纳虚部构成 B4=B(1:n2,1:n2); %所有PQ节点 %---------------------------------------------- %---------------下面是求ΔP,ΔQ矩阵--------------- DX=0;ICT=1;Mp=1;Mq=1; while ICT~=0 m1=1;m2=1; for i=1:n
潮流计算C++程序
程序为计算书3-4的过程 程序可以解决开式单直网络和树状网络的计算。树状网络计算时要自己先设定好支路的起始节点和终止节点标号以及计算顺序源代码： #include
潮流计算程序
最大负荷时 %本程序的功能是用牛顿——拉夫逊法进行潮流计算 % B1矩阵：1、支路首端号；2、末端号；3、支路阻抗；4、支路对地电纳 % 5、支路的变比；6、支路首端处于K侧为1，1侧为0 % B2矩阵：1、该节点发电机功率；2、该节点负荷功率；3、节点电压初始值% 4、PV节点电压V的给定值；5、节点所接的无功补偿设备的容量 % 6、节点分类标号 clear; n=12;%input('请输入节点数：n='); nl=14;%input('请输入支路数：nl='); isb=1;%input('请输入平衡母线节点号：isb='); pr=0.00001;%input('请输入误差精度：pr='); B1=[1 2 3.3+11.48i 0 1 0; 1 4 3.87+15.375i 0 1 0; 1 6 5.8+18.8i 0 1 0; 2 3 2.39+34.95i 0 1.02 1; 4 5 1.59+23.3i 0 1.02 1; 2 4 4.67+18.34i 0 1 0; 6 7 2.39+24.35i 0 1.02 1; 8 9 2.39+24.35i 0 1.02 1; 10 11 2.39+24.35i 0 1.02 1; 8 10 1.7+11.2i 0 1 0; 4 12 4.2+13.9i 0 1 0; 6 12 4.2+14.64i 0 1 0; 8 12 2.2+13.9i 0 1 0; 10 12 2.4+15.2i 0 1 0]; %input('请输入由支路参数形成的矩阵： B1='); B2=[0 0 231 231 0 1; 0 0 220 0 0 2; 0 100+75i 220 0 0 2; 0 0 220 0 0 2; 0 100+75i 220 0 0 2; 0 0 220 0 0 2; 0 150+112.5i 220 0 0 2; 0 0 220 0 0 2; 0 150+112.5i 220 0 0 2; 0 0 220 0 0 2; 0 150+112.5i 220 0 0 2; 450 0 231 231 0 3];%input('请输入各节点参数形成的矩阵： B2='); Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1= zeros(nl); %－－－－－－－修改部分－－－－－－－－－－－－ ym=1;