潮流计算的matlab程序实现方法
这是一个电气狗熬两个礼拜图书馆的成果,根据华中科技大学《电力系统分析》中原理编写,可用牛顿-拉夫逊和PQ分解法计算给定标幺值条件的潮流。本人水平有限,仅供参考,欢迎一起找Bug。
2019/11/17 添加算例系统图和基础数据、参考文献。
2019/01/05 添加word文档 潮流计算课程设计。
2018/07/06 说明:由于本人变压器建模与PSASP不同,本人使用模型如下图,参数输入时请按该模型计算。
2018/06/18 主程序更新:增加补偿电容参数
2018/06/18 增加下载地址,详见文章底部。
应用算例系统图:
基础数据(标幺值):
(1)母线数据
母线名 |
1 |
2 |
3 |
4 |
5 |
6 |
7 |
8 |
9 |
基准电压 |
230 |
230 |
230 |
230 |
230 |
230 |
18 |
13.8 |
16.5 |
(2)交流线数据
I侧母线 |
J侧母线 |
电阻 |
电抗 |
电纳的1/2 |
6 |
1 |
0.01 |
0.085 |
0.088 |
1 |
4 |
0.032 |
0.161 |
0.153 |
4 |
3 |
0.0085 |
0.072 |
0.0745 |
3 |
5 |
0.0119 |
0.1008 |
0.1045 |
5 |
2 |
0.039 |
0.17 |
0.179 |
2 |
6 |
0.017 |
0.092 |
0.079 |
(3)变压器数据
I侧母线 |
J侧母线 |
电抗 |
变比 |
9 |
6 |
0.0576 |
1 |
7 |
4 |
0.0625 |
1 |
8 |
5 |
0.0586 |
1 |
(4)发电机数据
母线名 |
母线类型 |
有功功率 |
无功功率 |
电压幅值 |
电压相角 |
9 |
slack |
- |
- |
1.04 |
0 |
7 |
PV |
1.63 |
- |
1.025 |
- |
8 |
PV |
0.85 |
- |
1.025 |
- |
(5)负荷数据
母线名 |
母线类型 |
有功功率 |
无功功率 |
1 |
PQ |
1.25 |
0.5 |
2 |
PQ |
0.9 |
0.3 |
3 |
PQ |
1 |
0.35 |
主程序
% file name:chaoliu_lj.m
% auther: 山东科技大学 罗江
% function:使用牛顿-拉夫逊法、PQ分解法计算潮流
% update:2018/6/18 13:22 增加补偿电容参数
%节点类型 标号
%PQ节点 1
%PV节点 2
%slack节点 3
%能计算给定标幺值网络,有且仅有一个平衡节点的潮流
%注意:母线标号顺序要求:PQ节点-PV节点-平衡节点
%若某元件不存在,其导纳为0,阻抗为infclear %清除工作空间变量
clc %清屏
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%数据输入(标幺值)
SB=100; %基准容量,单位MVA%母线基准电压
Bus=[230 230 230 230 230 230 18 13.8 16.5];%交流线参数:I侧母线 J侧母线 阻抗 1/2接地导纳
Line=[6 1 0.01+0.085i 0.088i;1 4 0.032+0.161i 0.153i; 4 3 0.0085+0.072i 0.0745i;3 5 0.0119+0.1008i 0.1045i;5 2 0.039+0.17i 0.179i;2 6 0.017+0.092i 0.079i];%变压器参数:I侧母线 J侧母线 阻抗 变比 %变压器阻抗归算到I侧
Trans=[9 6 0.0576i 1;7 4 0.0625i 1;8 5 0.0586i 1];%加接地电容器补偿: 母线 导纳
Cap=[2 0.0i];%发电机参数:母线 节点类型 P V/U θ
Gen=[9 3 1.04 0;7 2 1.63 1.025;8 2 0.85 1.025]; %负荷参数:母线 节点类型 P Q
%按参考方向,发电机发出功率(正值),负荷消耗功率(负值)
Load=[1 1 -1.25 -0.5;2 1 -0.9 -0.3;3 1 -1 -0.35];mode=1; %1-极坐标下牛拉法, 2-PQ分解法
Tmax=10; %最大迭代次数
limit=1.0e-8; %要求精度
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%变压器π型等效阻抗参数
Zt=zeros(size(Trans,1),3);
Zt(:,1)=Trans(:,3)./Trans(:,4);
Zt(:,2)=Trans(:,3)./(1-Trans(:,4));
Zt(:,3)=Trans(:,3)./(Trans(:,4).^2-Trans(:,4));
Trans_pi=[Trans(:,1:2) Zt(:,1) 1./Zt(:,2) 1./Zt(:,3)];n=numel(Bus); %总节点数
m=n-1; %PQ节点数
for i=1:size(Gen,1)%数组行数if Gen(i,2)==2 %除去PV节点就是PQ节点m=m-1;end
end
for i=1:size(Load,1)if Load(i,2)==2m=m-1;end
end
%PQ节点包含浮游节点,其PQ=0%提取P,Q,U向量
P=zeros(1,n); %P,Q为原始数据,Pi,Qi为计算结果
Q=zeros(1,n);
U=ones(1,n); %电压初始值由此确定
cita=zeros(1,n); %此处未知节点皆设为1.0∠0 %注意:此处角度单位为度,提取后再转换成弧度,后面计算使用弧度
for i=1:size(Gen,1)if Gen(i,2)==1 %PQ节点P(Gen(i,1))=Gen(i,3);Q(Gen(i,1))=Gen(i,4);endif Gen(i,2)==2 %PV节点P(Gen(i,1))=Gen(i,3);U(Gen(i,1))=Gen(i,4);endif Gen(i,2)==3 %slack节点U(Gen(i,1))=Gen(i,3);cita(Gen(i,1))=Gen(i,4);end
end
for i=1:size(Load,1)if Load(i,2)==1 %PQ节点P(Load(i,1))=Load(i,3);Q(Load(i,1))=Load(i,4);endif Load(i,2)==2 %PV节点P(Load(i,1))=Load(i,3);U(Load(i,1))=Load(i,4);endif Load(i,2)==3 %slack节点U(Load(i,1))=Load(i,3);cita(Load(i,1))=Load(i,4);end
end
disp('初始条件:')
disp('各节点有功:')
disp(P);
disp('各节点无功:')
disp(Q);
disp('各节点电压幅值:')
disp(U);
cita=(deg2rad(cita)); %角度转换成弧度
disp('各节点电压相角(度):')
disp(rad2deg(cita)); %显示依然使用角度%节点导纳矩阵的计算
Y=zeros(n); %新建节点导纳矩阵
y=zeros(n); %网络中的真实导纳
%计算y(i,j)
for i=1:size(Line,1) %与交流线联结的真实导纳ii=Line(i,1); jj=Line(i,2);y(ii,jj)=1/Line(i,3);y(jj,ii)=y(ii,jj);
end
for i=1:size(Trans_pi,1) %与变压器联结的真实导纳ii=Trans_pi(i,1); jj=Trans_pi(i,2);y(ii,jj)=1/Trans_pi(i,3);y(jj,ii)=y(ii,jj);
end
%计算y(i,i)
for i=1:size(Line,1) %与交流线联结的对地导纳ii=Line(i,1); jj=Line(i,2);y(ii,ii)=y(ii,ii)+Line(i,4);y(jj,jj)=y(jj,jj)+Line(i,4);
end
for i=1:size(Trans_pi,1) %与变压器联结的对地导纳ii=Trans_pi(i,1); jj=Trans_pi(i,2);y(ii,ii)=y(ii,ii)+Trans_pi(i,4);y(jj,jj)=y(jj,jj)+Trans_pi(i,5);
end
%算上补偿电容
for i=1:size(Cap,1)ii=Cap(i,1);y(ii,ii)=y(ii,ii)+Cap(i,2);
end
%由y计算Y
ysum=sum(y,1); %每一行求和
for i=1:nfor j=1:nif i==jY(i,j)=ysum(i);elseY(i,j)=-y(i,j);endend
end
disp('节点导纳矩阵:');
disp(Y);G=real(Y); %电导矩阵
B=imag(Y); %电纳矩阵%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%以上为基础数据整理
%接下来是牛拉法的大循环
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%计算功率不平衡量
[dP,dQ,Pi,Qi]=Unbalanced( n,m,P,Q,U,G,B,cita );
disp('有功不平衡量:');
disp(dP);
disp('无功不平衡量:');
disp(dQ);for i=1:Tmaxfprintf('第%d次迭代\n',i);%雅可比矩阵的计算if(mode==1)J=Jacobi( n,m,U,cita,B,G,Pi,Qi );disp('雅可比矩阵');disp(J);end%求解节点电压修正量if(mode==1)[dU,dcita]=Correct( n,m,U,dP,dQ,J ); else[dU,dcita]=PQ_LJ( n,m,dP,dQ,U,B );enddisp('电压、相角修正量:');disp(dU);disp(rad2deg(dcita));%修正节点电压U=U+dU;cita=cita+dcita;disp('节点电压幅值:');disp(U);disp('节点电压相角:');disp(rad2deg(cita));%计算功率不平衡量[dP,dQ,Pi,Qi]=Unbalanced( n,m,P,Q,U,G,B,cita );disp('有功不平衡量:');disp(dP);disp('无功不平衡量:');disp(dQ);if (max(abs(dP))<limit && max(abs(dQ))<limit )break;end%if
end%for%迭代结束,判断收敛
if (max(abs(dP))<limit && max(abs(dQ))<limit )disp('计算收敛');
elsedisp('计算不收敛或未达到要求精度');
end
%打印功率
fprintf('迭代总次数:%d\n', i);
disp('节点电压幅值:');
disp(U);
disp('节点电压相角:');
disp(rad2deg(cita));
disp('有功计算结果:');
disp(Pi);
disp('无功计算结果:');
disp(Qi);
子程序一
% filename:Unbalanced.m
% author: 山东科技大学 罗江
% function: 计算功率不平衡量
function [ dP,dQ,Pi,Qi ] = Unbalanced( n,m,P,Q,U,G,B,cita )
%计算ΔPi有功的不平衡量
for i=1:nfor j=1:nPn(j)=U(i)*U(j)*(G(i,j)*cos(cita(i)-cita(j))+B(i,j)*sin(cita(i)-cita(j)));endPi(i)=sum(Pn);
end
dP=P(1:n-1)-Pi(1:n-1); %dP有n-1个%计算ΔQi无功的不平衡量
for i=1:nfor j=1:nQn(j)=U(i)*U(j)*(G(i,j)*sin(cita(i)-cita(j))-B(i,j)*cos(cita(i)-cita(j)));endQi(i)=sum(Qn);
end
dQ=Q(1:m)-Qi(1:m); %dQ有m个end%func
子程序二
% filename:Jacobi.m
% author:山东科技大学 罗江
% function: 计算雅可比矩阵
function [ J ] = Jacobi( n,m,U,cita,B,G,Pi,Qi )
%雅可比矩阵的计算
%分块 H N K L
%i!=j时
for i=1:n-1for j=1:n-1H(i,j)=-U(i)*U(j)*(G(i,j)*sin(cita(i)-cita(j))-B(i,j)*cos(cita(i)-cita(j)));end
end
for i=1:n-1for j=1:mN(i,j)=-U(i)*U(j)*(G(i,j)*cos(cita(i)-cita(j))+B(i,j)*sin(cita(i)-cita(j)));end
end
for i=1:mfor j=1:n-1K(i,j)=U(i)*U(j)*(G(i,j)*cos(cita(i)-cita(j))+B(i,j)*sin(cita(i)-cita(j)));end
end
for i=1:mfor j=1:mL(i,j)=-U(i)*U(j)*(G(i,j)*sin(cita(i)-cita(j))-B(i,j)*cos(cita(i)-cita(j)));end
end
%i==j时
for i=1:n-1H(i,i)=U(i).^2*B(i,i)+Qi(i);
end
for i=1:mN(i,i)=-U(i).^2*G(i,i)-Pi(i);
end
for i=1:mK(i,i)=U(i).^2*G(i,i)-Pi(i);
end
for i=1:mL(i,i)=U(i).^2*B(i,i)-Qi(i);
end%合成雅可比矩阵
J=[H N;K L];end
子程序三
% filename:Correct.m
% author:山东科技大学 罗江
% function:修正节点电压
function [ dU,dcita ] = Correct( n,m,U,dP,dQ,J )
%求解节点电压修正量
for i=1:mUd2(i,i)=U(i);
end
dPQ=[dP dQ]';
dUcita=(-inv(J)*dPQ)';
dcita=dUcita(1:n-1);
dcita=[dcita 0];
dU=(Ud2*dUcita(n:n+m-1)')';
dU=[dU zeros(1,n-m)];end
子程序四
% filename:PQ_LJ.m
% author:山东科技大学 罗江
% function:使用PQ分解法计算电压修正量
function [ dU,dcita ] = PQ_LJ( n,m,dP,dQ,U,B )
dP_U=dP./U(1:n-1);
dQ_U=dQ./U(1:m);
dUdcita=(-inv(B(1:n-1,1:n-1))*dP_U')';
dcita=dUdcita./U(1:n-1);
dU=(-inv(B(1:m,1:m))*dQ_U')';
dU=[dU zeros(1,n-m)];
dcita=[dcita 0];%补零end
程序下载地址:点击下载程序
参考文献
[1]《电力系统分析(第四版)》 何仰赞 温增银,华中科技大学出版社,2017.
[2]《MATLAB/Simulink电力系统建模与仿真(第2版)》 于群 曹娜,机械工业出版社,2017.
[3]《电力系统分析课程设计指导及示例分析》郭丽萍 顾秀芳,中国水利水电出版社,2016.
[4]《电力系统分析的计算机算法》邱晓燕 刘天琪 黄媛,中国电力出版社,2015.
潮流计算的matlab程序实现方法相关推荐
- 交流潮流matlab程序,大神们,求个电力系统潮流计算的matlab程序。
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 n=input('Please input n\n'); %n表示系统的节点数 d0=input('Please input d0\n'); %d0表示系 ...
- matlab电力系统潮流计算,大神们,求个电力系统潮流计算的matlab程序。
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 n=input('Please input n\n'); %n表示系统的节点数 d0=input('Please input d0\n'); %d0表示系 ...
- 电力系统潮流计算matlab程序,大神们,求个电力系统潮流计算的matlab程序。
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 n=input('Please input n\n'); %n表示系统的节点数 d0=input('Please input d0\n'); %d0表示系 ...
- 电力系统潮流计算程序 matlab,大神们,求个电力系统潮流计算的matlab程序。
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 n=input('Please input n\n'); %n表示系统的节点数 d0=input('Please input d0\n'); %d0表示系 ...
- matlab在电力系统潮流计算程序,大神们,求个电力系统潮流计算的matlab程序。
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 n=input('Please input n\n'); %n表示系统的节点数 d0=input('Please input d0\n'); %d0表示系 ...
- matlab潮流计算求节点自导纳,大神们,求个电力系统潮流计算的matlab程序。
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 n=input('Please input n\n'); %n表示系统的节点数 d0=input('Please input d0\n'); %d0表示系 ...
- 电力系统潮流计算及Matlab编程实现
目录 1. 潮流计算: 2. 潮流计算常用算法: 2.1 牛顿-拉夫逊算法 2.1.1 牛顿-拉夫逊法的基本原理 2.1.2 潮流计算的修正方程 2.1.3 节点电压用极坐标表示时的牛顿-拉夫逊法潮 ...
- 二阶锥松弛在配电网最优潮流计算中的应用(IEEE33节点配电网最优潮流算例matlab程序)(yalmip+cplex)
二阶锥规划在配电网最优潮流计算中的应用IEEE33节点配电网最优潮流算例matlab程序(yalmip+cplex) 参考文献:二阶锥规划在配电网最优潮流计算中的应用 最优潮流计算是电网规划.优化运行 ...
- 牛顿拉夫逊基波潮流计算通用型程序,runpf函数的替换
牛顿拉夫逊基波潮流计算通用型程序,runpf函数的替换,可提供matlab版和python版 ID:4480672886448715
- python-pandapower电力系统潮流计算无法收敛情况解决方法
python-pandapower电力系统潮流计算无法收敛情况解决方法 无法收敛原因 解决办法 示例 数据集 无法收敛原因 pandpower网络上的潮流计算可能会因为各种各样的原因而无法收敛,这通常 ...
最新文章
- 机房收费系统【VB版】——前期准备
- 从DCGAN到SELF-MOD:GAN的模型架构发展一览
- 科技部认定的独角兽名单来了!共164家
- android listView嵌套gridview的使用心得
- 两篇关于MCU的嵌入式应用的文章【ZZ】
- 直播|实时音视频抗弱网技术揭秘
- 210221阶段三线程、信号量、互斥锁
- Redis数据库的连接
- 最大流问题 Edmonds-Karp算法
- 腾讯:人们回归工作导致四季度游戏收入减缓
- 从其他电脑拷mysql到自己电脑_mysql 数据库复制到其他电脑
- async 和 await 关键字
- 淘宝客APP如何配置阿里妈妈sdk详细教程(uniapp配置)
- 积分商城系统积分兑换运营开源架构
- 微信公众平台测试号的申请与使用
- ckplayer ajax,谁能帮我做一个脚本啊?能让这个网页视频播放可以拉动进度条 可以快进...
- linux内核 4g拨号,openwrt 基于qmi的 3G|4G拨号
- android 清理系统缓存文件怎么恢复,文件过期或已被清理怎么恢复(微信如何恢复已清理文件)...
- linux 性能分析工具perf
- 力扣 1598. 文件夹操作日志搜集器
热门文章
- 企业邮箱在outlook登录邮件如何撤回?
- 斐波那契数列 Java 实现。
- ArcGIS单波段影像重分类与批处理
- 微电子计算机是信息技术的,信息技术说第三十五说,计算机微电子技术
- 大型门户网站的商业计划书(包括技术解决方案)
- 初次接触面元法对螺旋桨的性能预报,发现之前很多学者都是用fortran进行编程进行性能预报,为什么不用matlab呢,两者的差异在哪里,建议初学者用这哪个软件呢
- Linux和onenote很像的软件,Microsoft OneNote替代方案?
- Java SE基础教程——Eclipse开发工具的安装与使用
- 超简单的方法找出QQ共同好友
- 3月第4周网络安全报告:境内76.2万个主机感染病毒