matlab潮流计算程序
MATpower潮流计算软件的使用

MATpower软件的使用一、MATpower软件的使用方法在MATLAB软件中的命令窗口输入runpf(‘程序名’)就可以通过MATpower已经编好的程序计算潮流,而函数runpf的参数是相应需计算潮流的数据文件。
数据文件主要用来定义和返回一下四个变量。
1、baseMVA baseMVA是一个标量,用来设置基准容量。
对于计算中采用有名值,可以根据实际情况设置,如设置100MV A;对于计算中采用标幺值,一般设置为12、bus bus变量是一个矩阵,用来设置电网中各节点参数,该矩阵内的参数如下:%% bus data%bus_i type Pd Qd Gs Bs area Vm Va baseKV zone Vmax Vmin其中,第1列参数即bus_i用来设置母线编号,范围为1~29997;第2列参数type用来设置母线类型,1为PQ节点,2为PV节点,3为平衡节点;第3列参数Pd用来设置母线注入负荷的有功功率;第4列参数Qd用来设置母线注入负荷的无功功率;第5列参数Gs用来设置与母线并联的电导;第6列参数Bs用来设置与母线并联的电纳;第7列参数area用来设置电网断面号,可设置范围为1~100,一般设置为1;第8列参数Vm用来设置母线电压的幅值初值;第9列参数Va用来设置母线电压的相角初值;第10列参数baseKV用来设置该母线的基准电压;第11列参数zone用来设置省耗分区号,可设置范围为1~999,一般设置为1;第12列参数Vmax用来设置工作时母线电压最高幅值;第13列参数Vmin用来设置工作时母线电压最低幅值。
3、gen gen变量是一个矩阵,用来设置接入电网的发电机参数,该矩阵的参数如下:%% generator data%bus Pg Qg Qmax Qmin Vg mBase status Pmax Pmin 其中,第1列参数bus用来设置接入发电机的母线编号;第2列参数Pg用来设置接入发电机的有功功率,注意功率输入的是有名值;第3列参数Qg用来设置接入发电机的无功功率;第4列参数Qmax用来设置接入发电机的无功功率的最大允许值;第5列参数Qmin用来设置接入发电机的无功功率的最小允许值;第6列Vg用来设置接入发电机的工作电压,注意输入的是标幺值;第7列mBase用来设置接入发电机的功率基准;第8列status用来设置发电机的工作状态,1表示投入运行,2表示投出运行;第9列Pmax用来设置接入发电机的无功功率的最大允许值;第10列参数Pmin用来设置接入发电机的无功功率的最小允许值。
基于MATLAB的潮流计算源程序代码(优.选)

%*************************电力系统直角坐标系下的牛顿拉夫逊法潮流计算**********clearclcload E:\data\IEEE014_Node.txtNode=IEEE014_Node;weishu=size(Node);nnum=weishu(1,1); %节点总数load E:\data\IEEE014_Branch.txtbranch=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);endif 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); %对地导纳endY(s,s)=Y(s,s)-j*B/2;Y(e,e)=Y(e,e)-j*B/2; %自导纳的计算情形endfor t=1:nnum;Y(t,t)=-sum(Y(t,:))+Node(t,12)+j*Node(t,13);%求支路自导纳endG=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)==3blancenode=m; %平衡节点编号else if Node(m,2)==0pq=pq+1;pqnode(1,pq)=m; %PQ 节点编号else if Node(m,2)==2pv=pv+1;pvnode(1,pv)=m; %PV 节点编号endendendend%*****************************设置电压初值********************************** Uoriginal=zeros(1,nnum); %对各节点电压矩阵初始化for n=1:nnumUoriginal(1,n)=Node(n,9); %对各点电压赋初值if Node(n,9)==0;Uoriginal(1,n)=1; %该节点为非PV节点时,将电压值赋为1endendPresion=input('请输入误差精度要求:Presion=');disp('该电力系统节点数:');disp(nnum);xiumax=0.1;counter=0;while xiumax>Presion%****************************计算不平衡量***********************************e=real(Uoriginal); %取初始电压的实部f=imag(Uoriginal); %取初始电压的虚部deta=zeros(2*pq+2*pv,1); %构造储存功率变化量的列矩阵n=1;for m=1:pq;Pi=0;Qi=0;for t=1:nnum;Pi=Pi+e(1,pqnode(1,m))*(G(pqnode(1,m),t)*e (1,t)-B(pqnode(1,m),t)*f(1,t))+f(1,pqnode(1,m ))*(G(pqnode(1,m),t)*f(1,t)+B(pqnode(1,m),t) *e(1,t));%计算该PQ节点的负荷有功Qi=Qi+f(1,pqnode(1,m))*(G(pqnode(1,m),t)*e (1,t)-B(pqnode(1,m),t)*f(1,t))-e(1,pqnode(1,m ))*(G(pqnode(1,m),t)*f(1,t)+B(pqnode(1,m),t) *e(1,t));%计算该PQ节点的负荷无功endS1(1,pqnode(1,m))=Pi+j*Qi;P=(Node(pqnode(1,m),7)-Node(pqnode(1,m), 5))/Sj-Pi;%计算该PQ节点的实际有功功率deta(n,1)=P;%在该列向量中储存有功功率n=n+1;Q=(Node(pqnode(1,m),8)-Node(pqnode(1,m),6))/Sj-Qi;%计算该PQ节点的实际无功功率deta(n,1)=Q;%在该列向量中储存无功功率n=n+1;endfor m=1:pv;Pv=0; Qv=0;for t=1:nnum;Pv=Pv+e(1,pvnode(1,m))*(G(pvnode(1,m),t)* e(1,t)-B(pvnode(1,m),t)*f(1,t))+f(1,pvnode(1, m))*(G(pvnode(1,m),t)*f(1,t)+B(pvnode(1,m), t)*e(1,t));%计算该PV节点的负荷有功Ui=e(1,pvnode(1,m))^2+f(1,pvnode(1,m))^2; %计算该节点的负荷电压值Qv=Qv+f(1,pqnode(1,m))*(G(pqnode(1,m),t)* e(1,t)-B(pqnode(1,m),t)*f(1,t))-e(1,pqnode(1, m))*(G(pqnode(1,m),t)*f(1,t)+B(pqnode(1,m), t)*e(1,t));endS1(1,pvnode(1,m))=Pv+j*Qv;P=(Node(pvnode(1,m),7)-Node(pvnode(1,m), 5))/Sj-Pv; %计算该节点的实际有功功率deta(n,1)=P; %储存该有功功率n=n+1;U=Node(pvnode(1,m),3)^2-Ui; %计算电压变化量deta(n,1)=U; %储存该电压变化量n=n+1;enddeta;cerate=zeros(pq+pv,1);for k=1:pqcerate(k,1)=pqnode(1,k);endfor v=1:pvcerate(pq+v,1)=pvnode(1,v);end%******************************雅克比矩阵****************************** Jacob=ones(2*nnum-2);L=0;J=0;H=0;N=0; R=0;S=0;n=1;k=1;form=1:pq; %m表示雅克比矩阵中pq节点的行数for u=1:pq+pv; %u 表示雅克比矩阵中pq节点的列数t=cerate(u,1); %t为中间变量,用来标记雅克比矩阵中指定元素的个数if pqnode(1,m)~=t %非对角元素的情况H=G(pqnode(1,m),t)*f(1,pqnode(1,m))-B(pqn ode(1,m),t)*e(1,pqnode(1,m));N=G(pqnode(1,m),t)*e(1,pqnode(1,m))+B(pq node(1,m),t)*f(1,pqnode(1,m));L=H;J=-N;elseifpqnode(1,m)==t %对角线元素时的情况I=0;for g=1:nnumI=Y(t,g)*Uoriginal(1,g)+I; %计算节点的注入电流endaii=real(I);bii=imag(I);H=-B(t,t)*e(1,pqnode(1,m))+G(t,t)*f(1,pqnode (1,m))+bii;N=G(t,t)*e(1,pqnode(1,m))+B(t,t)*f(1,pqnode (1,m))+aii;L=-B(t,t)*e(1,pqnode(1,m))+G(t,t)*f(1,pqnode (1,m))-bii; J=-G(t,t)*e(1,pqnode(1,m))-B(t,t)*f(1,pqnode( 1,m))+aii;endendJacob(n,k)=H;k=k+1;Jacob(n,k)=N;k=k-1;n=n+1;Jacob(n,k)=J;k=k+1;Jacob(n,k)=L;n=n-1; k=k+1; %按照雅克比矩阵的排列规则排列pq节点的雅克比元素endk=1; n=2*m+1; %将光标定位于下一个待排列PQ节点元素的第一个位置endn=2*pq+1; k=1; %定位于PV节点的第一个位置处for m=1:pv;for u=1:pq+pv;t=cerate(u,1); %t为中间变量,用来标记雅克比矩阵中指定元素的位置if pvnode(1,m)~=t %非对角线元素情况H=G(pvnode(1,m),t)*f(1,pvnode(1,m))-B(pvno de(1,m),t)*e(1,pvnode(1,m));N=G(pvnode(1,m),t)*e(1,pvnode(1,m))+B(pvn ode(1,m),t)*f(1,pvnode(1,m));R=0; S=0;endifpvnode(1,m)==t %对角线元素情况I=0;for g=1:nnumI=Y(t,g)*Uoriginal(1,g)+I; %计算PV节点的注入电流endaii=real(I);bii=imag(I);H=-B(t,t)*e(1,pvnode(1,m))+G(t,t)*f(1,pvnode (1,m))+bii;N=G(t,t)*e(1,pvnode(1,m))+B(t,t)*f(1,pvnode( 1,m))+aii;R=2*f(1,pvnode(1,m));S=2*e(1,pvnode(1,m));endJacob(n,k)=H; k=k+1;Jacob(n,k)=N; k=k-1;n=n+1;Jacob(n,k)=R; k=k+1;Jacob(n,k)=S;n=n-1;k=k+1; %按照雅克比矩阵的排列规则排列PV节点的雅克比元素endk=1;n=n+2; %定位于下一个待排列PV节点的雅克比元素第一个位置end%*************************电压变化量化的计算与存储************************************ Detau=inv(Jacob)*deta; %构建电压的变化量的列向量f=zeros(1,nnum); %给电压实部赋初值0e=zeros(1,nnum); %给电压虚部赋初值0for p=1:(pq+pv);f(1,cerate(p,1))=j*Detau(2*p-1,1);%将电压变量的奇数行赋值给fe(1,cerate(p,1))=Detau(2*p,1); %将电压变量的偶数行赋值给eendt=e+f;xiumax=abs(Detau(1,1)); %将电压变化量的第一个元素赋值给最大允许误差for n=2:2*nnum-2;if abs(Detau(n,1))>xiumaxxiumax=abs(Detau(n,1)); %找出最大的电压误差endendUoriginal=Uoriginal+t; %迭代修正后的电压值counter=1+counter; %统计迭代次数enddisp('迭代次数counter:');disp(counter);%**************************平衡节点功率及显示**********************************m=blancenode;t=0;for n=1:nnum;t=t+(G(m,n)-j*B(m,n))*(real(Uoriginal(1,n))-j*i mag(Uoriginal(1,n)));endS1(1,m)=Uoriginal(1,m)*t;%**************************直角坐标下各节点电压及显示****************************U=zeros(1,nnum);for n=1:nnumUi(n,1)=Node(n,1);U(1,n)=real(Uoriginal(1,n))+i*imag(Uoriginal(1 ,n)); %将电压值由极坐标转化为直角坐标形式Ui(n,2)=U(1,n);Ui(n,3)=S1(1,n);enddisp('各节点电压直角坐标形式及节点注入功率:');disp(' 节点号节点电压值节点注入功率');disp(Ui);disp('修正电压的最大误差:')disp(xiumax);%**************************功率损耗************************************* **for m=1:bnum %支路功率及损耗startnode=branch(m,1);endnode=branch(m,2); %终止节点y=sum(Y,2);Sij=Uoriginal(1,startnode)*((conj(Uoriginal(1,s tartnode))*conj(y(startnode,1)))+(conj(Uorigi nal(1,startnode))-conj(Uoriginal(1,endnode))) *conj(-Y(startnode,endnode)));Sji=Uoriginal(1,endnode)*((conj(Uoriginal(1,e ndnode))*conj(y(endnode,1)))+(conj(Uoriginal (1,endnode))-conj(Uoriginal(1,startnode)))*co nj(-Y(endnode,startnode)));S(m,1)=startnode; %始节点S(m,2)=endnode; %末节点S(m,3)=Sij; %始节点流入的功率S(m,4)=Sji; %末节点流入的功率S(m,5)=Sij+Sji; %线路损耗enddisp('支路功率及损耗:');disp('节点号(I) 节点号(J) 支路功率(I-J)支路功率(J-I)线路损耗(delta_S)')disp(S); 最新文件---------------- 仅供参考--------------------已改成word文本--------------------- 方便更改。
电力系统潮流计算matlab程序

电力系统潮流计算matlab程序电力系统潮流计算是电力系统运行和规划中的重要环节,它用于计算电力系统中各节点的电压、功率和电流等参数。
随着电力系统规模的不断扩大和复杂性的增加,传统的手工计算方法已经无法满足需求,因此,利用计算机编程进行潮流计算成为了一种必要的选择。
Matlab是一种功能强大的科学计算软件,它提供了丰富的数学函数和工具箱,可以方便地进行电力系统潮流计算。
下面我将介绍一下如何使用Matlab编写电力系统潮流计算程序。
首先,我们需要建立电力系统的节点模型。
节点模型是电力系统中各节点的电压、功率和电流等参数的数学表示。
在Matlab中,我们可以使用矩阵来表示节点模型。
假设电力系统有n个节点,我们可以定义一个n×n的复数矩阵Y来表示节点之间的导纳关系,其中Y(i,j)表示节点i和节点j之间的导纳。
同时,我们还需要定义一个n×1的复数向量V来表示各节点的电压,其中V(i)表示节点i的电压。
接下来,我们需要编写潮流计算的主程序。
主程序的主要功能是根据节点模型和潮流计算算法,计算出各节点的电压、功率和电流等参数。
在Matlab中,我们可以使用循环语句和矩阵运算来实现潮流计算。
具体的计算过程可以参考电力系统潮流计算的算法。
在编写主程序之前,我们还需要定义一些输入参数,如电力系统的节点数、发电机节点和负荷节点等。
这些参数可以通过用户输入或者读取文件的方式获取。
同时,我们还需要定义一些输出参数,如各节点的电压、功率和电流等。
这些参数可以通过矩阵运算和循环语句计算得到,并输出到文件或者显示在屏幕上。
最后,我们需要进行程序的测试和调试。
可以通过输入一些测试数据,运行程序并检查输出结果是否正确。
如果发现程序有错误或者结果不准确,可以通过调试工具和打印调试信息的方式进行调试。
总之,利用Matlab编写电力系统潮流计算程序可以提高计算效率和准确性,为电力系统的运行和规划提供有力的支持。
当然,编写一个完整的潮流计算程序需要考虑很多细节和特殊情况,这需要有一定的电力系统和编程知识。
潮流计算MATLAB程序

0 0.343+0.21256i 1 0 0 2]
;%input('请输入各节点参数形成的矩阵: B2=');
Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1=zeros(nl);
J(p,q)=X4;J(m,q)=X2; % X4=dQ/df X2=dp/df
elseif j1==i&j1~=isb %非平衡节点&对角元
X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);% dP/de
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)
J(m,q)=X2;
end
end
else
%=============== 下面是针对PV节点来求取Jacobi矩阵的元素 ===========
8 9 0.03512+0.08306i 0.13455i 1 0]
B2=[0 0 1.1 1.1 0 1;
0 0 1 0 0 2;
Байду номын сангаас%================= 求取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
MATLAB电力系统PQ潮流计算程序设计

MATLAB电力系统PQ潮流计算程序设计1 绪论1.1潮流计算1.1.1 潮流计算概述电力系统潮流计算是研究电力系统稳态运行情况的一种计算,它根据给定的运行条件及系统接线情况确定整个电力系统各部分的运行状态:各母线的电压,各元件中流过的功率,系统的功率损耗等等。
在电力系统规划的设计和现有电力系统运行方式的研究中,都需要利用潮流计算来定量地分析比较供电方案或运行方式的合理性、可靠性和经济性。
此外,电力系统潮流计算也是计算系统动态稳定和静态稳定的基础。
所以潮流计算是研究电力系统的一种很重要也很基础的计算。
电力系统潮流计算也分为离线计算和在线计算两种,前者主要用于系统规划设计和安排系统的运行方式,后者则用于正在运行系统的随时监视及实时控制。
利用计算机进行电力系统潮流计算从50年代中期就已经开始。
在这20年内,潮流计算曾采用了各种不同的方法,这些方法的发展主要围绕着对潮流计算的一些基本要求进行的。
对潮流计算的要求可以归纳为下面几点:(1)计算方法的可靠性或收敛性;(2)对计算机内存量的要求;(3)计算速度;(4)计算的方便性和灵活性。
电力系统潮流计算问题在数学上是一组多元非线性方程式求解问题,其解法都离不开迭代。
因此,对潮流计算方法,首先要求它能可靠地收敛,并给出正确答案。
由于电力系统结构及参数的一些特点,并且随着电力系统不断扩大,潮流计算方程式的阶数也越来越高,对这样的方程式并不是任何数学方法都能保证给出正确答案的。
这种情况成为促使电力系统计算人员不断寻求新的更可靠方法的重要因素。
在用数字计算机解电力系统潮流问题的开始阶段,普遍采取以节点导纳矩阵为基础的逐次代入法。
这个方法的原理比较简单,要求的数字计算机内存量比较低,适应50年代电子计算机制造水平和当时电力系统理论水平。
但它的收敛性较差,当系统规模变大时,迭代次数急剧上升,在计算中往往出现迭代不收敛的情况。
这就迫使电力系统计算人员转向以阻抗矩阵为基础的逐次代入法。
基于MATLAB进行潮流计算

基于MATLAB进行潮流计算本文介绍了基于MATLAB软件的潮流计算方法。
电力系统潮流计算方法分为手算潮流和计算机潮流计算两类。
手算潮流主要适用于规模较小的辐射型电力潮流计算,而计算机潮流计算有两种途径:编程实现网络方程的迭代求解和借助电力系统分析仿真软件搭建系统模型完成潮流计算。
MATLAB具有强大的矩阵运算功能和电力系统仿真平台,可以为实现潮流计算提供更便捷的手段。
本文采用极坐标形式牛顿─拉夫逊法进行潮流计算,为其他形式的潮流计算提供借鉴。
Abstract: The power flow n method can be divided into two categories: hand n of tidal current and computer power flow XXX simplified equivalent circuits。
making it XXX: programming XXX ns。
or using power system XXX system model for power flow n。
MATLAB are has strong matrix ns and its power system XXX-Raphson method of power flow n in polar coordinates with MATLAB are。
and can serve as a reference for other forms of power flow n.1.电力系统中的牛顿法潮流计算是一种常用的电力系统分析方法。
该方法基于节点电压的相等条件和潮流方程的等式条件,通过迭代求解电压和相位的不平衡量,最终得到各节点的电压、相位和功率等参数。
2.牛顿法潮流计算的步骤包括输入系统原始数据、形成节点导纳矩阵、给定各节点电压初值、计算功率偏差向量、判断收敛条件、计算雅克比矩阵、解修正方程、计算节点电压和相位的修正值、迭代计算直至满足收敛条件、计算各节点功率等参数并输出计算结果。
matlab电力系统潮流计算程序
matlab电力系统潮流计算程序电力系统潮流计算是电力系统分析的关键步骤之一,用于确定电力系统各节点的电压和相角分布。
以下是一个简单的MATLAB电力系统潮流计算的基本步骤和代码示例:1.定义电力系统参数:-定义系统节点数量、支路数据、发电机数据、负荷数据等电力系统参数。
```matlab%电力系统参数busdata=[1,1.05,0,0,0,0,0,0;2,1.02,0,0,0,0,0,0;%...其他节点数据];linedata=[1,2,0.02,0.06,0.03;%...其他支路数据];gendata=[1,2,100,0,999,1.05,0.95;%...其他发电机数据];loaddata=[1,50,20;%...其他负荷数据];```2.构建潮流计算矩阵:-利用节点支路导纳、节点负荷和发电机功率等信息构建潮流计算的阻抗矩阵。
```matlabYbus=buildYbus(busdata,linedata);```3.迭代求解潮流方程:-利用迭代算法(如牛顿-拉夫森法)求解潮流方程,更新节点电压和相角。
```matlab[V,delta]=powerflow(Ybus,gendata,loaddata,busdata);```4.结果分析和可视化:-分析计算结果,可视化电压和相角分布。
```matlabplotVoltageProfile(busdata,V,delta);```这只是一个简单的潮流计算示例。
具体的程序实现可能涉及更复杂的算法和工程细节,取决于电力系统的复杂性和精确性要求。
您可能需要根据实际情况和数据格式进行调整和改进。
在实际工程中,也可以考虑使用专业的电力系统仿真软件。
电力系统潮流计算的MATLAB辅助程序设计-潮流计算程序
电力系统潮流计算的MATLAB辅助程序设计潮流计算,通常指负荷潮流,是电力系统分析和设计的主要组成部分,对系统规划、安全运行、经济调度和电力公司的功率交换非常重要。
此外,潮流计算还是其它电力系统分析的基础,比如暂态稳定,突发事件处理等。
现代电力系统潮流计算的方法主要:高斯法、牛顿法、快速解耦法和MATLAB的M语言编写的MATPOWER4.1,这里主要介绍高斯法、牛顿法和快速解耦法.高斯法的程序是lfgauss,其与lfybus、busout和lineflow程序联合使用求解潮流功率。
lfybus、busout和lineflow程序也可与牛顿法的lfnewton程序和快速解耦法的decouple程序联合使用。
(读者可以到MATPOWER主页下载MATPOWER4.1,然后将其解压到MATLAB目录下,即可使用该软件进行潮流计算)一、高斯—赛德尔法潮流计算使用的程序:高斯—赛德法的具体使用方法读者可参考后面的实例,这里仅介绍各程序的编写格式:lfgauss:该程序是用高斯法对实际电力系统进行潮流计算,需要用到busdata和linedata两个文件。
程序设计为输入负荷和发电机的有功MW和无功Mvar,以及节点电压标幺值和相角的角度值。
根据所选复功率为基准值将负荷和发电机的功率转换为标幺值。
对于PV节点,如发电机节点,要提供一个无功功率限定值。
当给定电压过高或过低时,无功功率可能超出功率限定值。
在几次迭代之后(高斯—塞德尔迭代为10次),需要检查一次发电机节点的无功出力,如果接近限定值,电压幅值进行上下5%的调整,使得无功保持在限定值内。
lfybus:这个程序需要输入线路参数、变压器参数以及变压器分接头参数。
并将这些参数放在名为linedata的文件中。
这个程序将阻抗转换为导纳,并得到节点导纳矩阵.busout:该程序以表格形式输出结果,节点输出包括电压幅值和相角,发电机和负荷的有功和无功功率,以及并联电容器或电抗器的有功和无功功率。
潮流计算MATLAB 粗略程序
%========================================================================== %========================================================================== %========================================================================== %潮流计算MATLAB 粗略程序 C.zhou 2009.3.27%========================================================================== %========================================================================== %========================================================================== %creat a new_datat=0;s=0;r=0;w=0;number=input('How many node are there=');% Convert Pq to a new arrayfor ii=1:numberif data(ii,4)==1t=t+1;for jj=1:14new_data1(t,jj)=data(ii,jj);end;a(1,t)=ii;s=s+1; %record the number of the PQ node end;end;%Convert pv to a new arrayfor ii=1:numberif data(ii,4)==2t=t+1;for jj=1:14new_data1(t,jj)=data(ii,jj);end;a(1,t)=ii;r=r+1; %record the number of the PV node end;end;%Convert set_v to a new arrayfor ii=1:numberif data(ii,4)==3t=t+1;for jj=1:14new_data1(t,jj)=data(ii,jj);end;a(1,t)=ii;w=w+1;end;end;%creat a new_data2[x,y]=size(data2)for ii=1:xfor jj=1:2for mm=1:numberif data2(ii,jj)==a(1,mm)new_data2(ii,jj)=mm;end;end;end;end;for ii=1:xfor jj=3:14new_data2(ii,jj)=data2(ii,jj);end;end;%creat a YY=zeros(number,number);YY=zeros(number,number);yy=zeros(number,number);for ii=1:x% for jj=1:14iii=new_data2(ii,1);jjj=new_data2(ii,2);if new_data2(ii,5)==2sub=new_data2(ii,6)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6 ))-new_data2(ii,7)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6))*i;Y(iii,jjj)=-sub./new_data2(ii,14);YY(iii,jjj)=sub./new_data2(ii,14);Y(jjj,iii)=-sub/new_data2(ii,14);YY(jjj,iii)=sub./new_data2(ii,14);yy(iii,jjj)=(1.-new_data2(ii,14))./(new_data2(ii,14).*new_data2(ii,14)).*sub;yy(jjj,iii)=(new_data2(ii,14)-1)./(new_data2(ii,14)).*sub;elseY(iii,jjj)=-new_data2(ii,6)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2 (ii,6))+new_data2(ii,7)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6))*i;YY(iii,jjj)=new_data2(ii,6)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data 2(ii,6))-new_data2(ii,7)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6))*i;Y(jjj,iii)=-new_data2(ii,6)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2 (ii,6))+new_data2(ii,7)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6))*i;YY(jjj,iii)=new_data2(ii,6)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6))-new_data2(ii,7)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6))*i;yy(iii,jjj)=new_data2(ii,8)./2.*i;yy(jjj,iii)=new_data2(ii,8)./2.*i;end;%end;end;for iii=1:numberY(iii,iii)=0;end;%for ii=1:x% for jj=1:14for iii=1:numberfor jj=1:number% if iii~=jjY(iii,iii)=Y(iii,iii)+YY(iii,jj)+yy(iii,jj);% end;end;end;%creat B, Gfor ii=1:numberfor jj=1:numberG(ii,jj)= real(Y(ii,jj));B(ii,jj)= imag(Y(ii,jj));end;end;%creat Initial_P Initial_Q Initial_Vfor ii=1:(s+r)set_P(ii,1)=(new_data1(ii,9)-new_data1(ii,7))./100;end;for ii=1:s;set_Q(ii,1)=(new_data1(ii,10)-new_data1(ii,8))./100;end;for ii=1:rset_V(ii,1)=new_data1(ii+s,12).*new_data1(ii+s,12);%try to modify for sike of correcting end;Initial_p_q_v=[set_P;set_Q;set_V];disp(Initial_p_q_v);%creat Initial_e,Initial_ffor ii=1:number-1e(ii,1)=1;f(ii,1)=0.0;%change f to test used to be 1.0end;e(number,1)=new_data1(number,12);f(number,1)=0;% e(64,1)=0.88;%test 118ieee% f(64,1)=0.39395826829394;% f(14,1)=0;% e(10,1)=1.045;%e(11,1)=1.01;%e(12,1)=1.07;%e(13,1)=1.09;%//////////////////////////////////////////////////////////////////////////%/////////////////////////////////////////////////////////////////////////%//////////////////////////////////////////////////////////////////////////%//////////////////////////////////////////////////////////////////////////% Start NEWTOWN CALULATIONfor try_time=1:25%Creat every node consume P Q and Un=s;m=r;for ii=1:(n+m)sum1=0;for jj=1:(n+m+1)sum1=sum1+e(ii,1).*(G(ii,jj).*e(jj,1)-B(ii,jj).*f(jj,1))+f(ii,1).*(G(ii,jj).*f(jj,1)+B(ii,jj).*e(jj,1));end;p(ii,1)=sum1;end;for ii=1:nsum2=0;for jj=1:(n+m+1)sum2=sum2+f(ii,1).*(G(ii,jj).*e(jj,1)-B(ii,jj).*f(jj,1))-e(ii,1).*(G(ii,jj).*f(jj,1)+B(ii,jj).*e(jj,1));end;q(ii,1)=sum2;end;disp('q=');disp(q);u=zeros((n+m),1);for ii=(n+1):(n+m)u(ii,1)=e(ii,1).*e(ii,1)+f(ii,1).*f(ii,1);end;for ii=n+1:(n+m)extra_u((ii-n),1)=u(ii,1);end;disp('extra_u=');disp(extra_u);sum=[p;q;extra_u];disp(sum)disp(s);disp(p);%creat Jacobiandisp(n);disp(m);for ii=1:(n+m)for jj=1:(n+m)if (ii~=jj)PF(ii,jj)=B(ii,jj).*e(ii,1)-G(ii,jj).*f(ii,1);PE(ii,jj)=-G(ii,jj).*e(ii,1)-B(ii,jj).*f(ii,1);elsess=0;qq=0;for num=1:(n+m+1)ss=ss+G(ii,num).*f(num,1)+B(ii,num).*e(num,1);qq=qq+G(ii,num).*e(num,1)-B(ii,num).*f(num,1);end;PF(ii,jj)=-ss+B(ii,jj).*e(ii,1)-G(ii,jj).*f(ii,1);%TEST+1PE(ii,jj)=-qq-G(ii,jj).*e(ii,1)-B(ii,jj).*f(ii,1);%TEST+1end;end;end;copy=3.14159;disp('================copy================')for ii=1:nfor jj=1:m+nif (ii~=jj)QE(ii,jj)=B(ii,jj).*e(ii,1)-G(ii,jj).*f(ii,1);%TEST+1QF(ii,jj)=G(ii,jj).*e(ii,1)+B(ii,jj).*f(ii,1);%TEST+1elsess=0;qq=0;for num=1:(n+m+1)ss=ss+G(ii,num).*f(num,1)+B(ii,num).*e(num,1);qq=qq+G(ii,num).*e(num,1)-B(ii,num).*f(num,1);end;QF(ii,jj)=-qq+G(ii,jj).*e(ii,1)+B(ii,jj).*f(ii,1);%TEST+1QE(ii,jj)=ss+B(ii,jj).*e(ii,1)-G(ii,jj).*f(ii,1);%TEST+1end;end;end;%disp('QF');%disp(QF);%disp('QE');%disp(QE);UE=zeros((n+m),(n+m));UF=zeros((n+m),(n+m));for ii=n+1:n+mfor jj=1:(n+m)if (ii~=jj)UE(ii,jj)=0;UF(ii,jj)=0;elsess=0;qq=0;for num=1:(n+m+1)ss=ss+G(ii,num).*f(num,1)+B(ii,num).*e(num,1);qq=qq+G(ii,num).*e(num,1)-B(ii,num).*f(num,1);end;UF(ii,jj)=-2.*f(ii,1);UE(ii,jj)=-2.*e(ii,1);end;end;end;for ii=(n+1):(n+m)for jj=1:(n+m)extra_UE((ii-n),jj)=UE(ii,jj);extra_UF((ii-n),jj)=UF(ii,jj);end;end;%disp('extra_UE');%disp(extra_UE);%disp('extra_Uf');%disp(extra_UF);Jacobian=[PF,PE;QF,QE;extra_UF,extra_UE];%disp('Jacobian=');%disp(Jacobian);%creat substract resultsubstract_result=Initial_p_q_v-sum;%disp('substract_result');%disp(substract_result);%calculate delta_f_edelta_f_e=-inv(Jacobian)*substract_result; %disp(delta_f_e);for ii=1:number-1;f(ii,1)=f(ii,1)+delta_f_e(ii,1);e(ii,1)=e(ii,1)+delta_f_e(ii+number-1,1); end;if max(substract_result)<1e-4break;end ;end;%disp('substract_result');%disp(substract_result);%disp('e=');%disp(e);%disp('f=');%disp(f);for ii=1:numberuuu(ii,1)= e(ii,1).*e(ii,1)+f(ii,1).*f(ii,1);U_RESULT(ii,1)=sqrt(uuu(ii,1));end;for ii=1:numberfor jj=1:numberif ii==a(1,jj)Old_Uresult(ii,1)=U_RESULT(jj,1)end;end;end;for ii=1:numberOld_Uresult(ii,2)=ii;end;%disp('U_result');%disp(U_RESULT);disp('====================================='); disp('The last result is :')disp('===========U===================BUS-NO.'); disp('U=')disp(Old_Uresult);%calculate the anglePI=3.141592for ii=1:numberAngle(ii,1)=atan(f(ii,1)./e(ii,1))./PI*180;end;for ii=1:numberfor jj=1:numberif ii==a(1,jj)Old_Angle(ii,1)=Angle(jj,1);Old_Angle(ii,2)=ii;end;end;end;disp('=======Angle===================BUS-NO.'); disp('Angle=');disp(Old_Angle);disp('=====Try-times=======================') disp('Try-times=')disp(try_time);%disp('p====================');%disp(p);% for jj=1:number% if a(1,jj)==118% man=jj% end;%end;%disp('man=========');%disp(man)sum4=0;for jj=1:numberY_conj(number,jj)=conj(Y(number,jj));sum4=sum4+Y_conj(number,jj).*(e(jj,1)-f(jj,1)*i);end;%sum4=sum4*e(number,1);disp('===============Balance P Q=========BUS-NO');%disp(sum4);Blance_Q(1,1)=imag(sum4)*100;Blance_Q(1,2)=a(1,number);Blance_P(1,1)=real(sum4)*100;Blance_P(1,2)=a(1,number);disp('Q Of the Balance node= ');disp(Blance_Q);disp('P Of the Balance node= ')disp(Blance_P);disp('=================================BUS-NO');%calculate the Q of the P-V nodeQ_PV_node=zeros(number,2);Y_conj=conj(Y);for ii=(s+1):(s+r)for jj=1:numberQ_PV_node(ii,1)=Q_PV_node(ii,1)+(e(ii,1)+f(ii,1)*i).*(Y_conj(ii,jj).*(e(jj,1)-f(jj,1)*i));end;end;for ii=(s+1):(s+r);Q_PV_node(ii,1)=Q_PV_node(ii,1).*100+new_data1(ii,8)*i;end;disp('This program is from /breadwinner') ;for ii=1:numberold_number=a(1,ii);Q_PV_node_old(old_number,1)=Q_PV_node(ii,1);end;for ii=1:numberQ_PV_node_old(ii,1)=imag(Q_PV_node_old(ii,1));end;for ii=1:numberQ_PV_node_old(ii,2)=ii;end;disp('Q gen=');disp(Q_PV_node_old);。
22节点潮流计算matlab仿真程序
【22节点潮流计算matlab仿真程序详解】引言主题指定的22节点潮流计算matlab仿真程序是电力系统领域中的重要内容。
在本文中,我将以深入、广泛的方式来探讨这一主题,帮助您全面了解和掌握相关知识。
一、22节点潮流计算的概念1. 22节点潮流计算的基本介绍22节点潮流计算是电力系统分析中的重要部分,它是对电力系统节点之间功率、电压、相角等参数进行计算和分析的过程。
通过22节点潮流计算,可以有效评估电力系统的稳定性和可靠性。
2. 22节点潮流计算的基本原理在进行22节点潮流计算时,需要考虑节点之间的功率平衡方程、电压平衡方程等基本原理。
这些原理是潮流计算的理论基础,也是实际仿真程序设计的重要依据。
二、22节点潮流计算matlab仿真程序的设计与实现1. matlab在电力系统仿真中的应用matlab作为一种功能强大的科学计算软件,在电力系统领域有着广泛的应用。
通过matlab,可以方便地进行潮流计算、稳定性分析、拓扑优化等工作。
2. 22节点潮流计算matlab仿真程序的设计要点在设计22节点潮流计算matlab仿真程序时,需要考虑程序的模块化、可扩展性、运行效率等要点。
通过合理的程序设计,可以提高仿真程序的稳定性和可靠性。
三、22节点潮流计算matlab仿真程序的应用实例1. 22节点潮流计算matlab仿真程序的基本框架通过一份完整的22节点潮流计算matlab仿真程序,来实际演示程序的结构和实现过程。
这样可以让读者更直观地理解程序的设计和应用。
2. 22节点潮流计算matlab仿真程序的实际应用案例以一个具体的电力系统实例为例,演示22节点潮流计算matlab仿真程序的实际应用过程。
这样可以让读者更清晰地了解程序在实际工程中的价值和作用。
结束语通过本文的深度探讨和广度展示,相信您已经对22节点潮流计算matlab仿真程序有了更全面的了解。
也希望您能够在实际工作中灵活运用所学知识,不断提升自己在电力系统领域的技术水平。