matlab潮流计算Word文件下载.docx

1、%-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=B

2、1（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侧%求导纳矩阵disp（导纳矩阵 Y=disp（Y）%-G=real（Y）;B=imag（Y）; %分解出导纳阵的实部和虚部n %给定各节点初始电压的实部和虚部e（i）=real（B2（i,3）

3、;f（i）=imag（B2（i,3）; V（i）=B2（i,4）; %PV节点电压给定模值n %给定各节点注入功率 S（i）=B2（i,1）-B2（i,2）; %i节点注入功率SG-SL B（i,i）=B（i,i）+B2（i,5）; %i节点无功补偿量%-P=real（S）;Q=imag（S）; %分解出各节点注入的有功和无功功率ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0; %迭代次数ICT1、a；不满足收敛要求的节点数IT2while IT2=0 % N0=2*n 雅可比矩阵的阶数；N=N0+1扩展列 IT2=0;a=a+1;n if i=isb %非平衡节点C（i）=0;

4、D（i）=0;for j1=1: 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） 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）%求i节点有功和无功功率P,Q的计算值 V2=e（i）2+f（i）2; %电压模

5、平方%以下针对非PV节点来求取功率差及Jacobi矩阵元素- if B2（i,6）=3 %非PV节点 DP=P（i）-P1; %节点有功功率差 DQ=Q（i）-Q1; %节点无功功率差%以上为除平衡节点外其它节点的功率计算-%求取Jacobi矩阵- 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

6、*i-1;q=2*j1-1; J（p,q）=X3;J（p,N）=DQ;m=p+1; % X3=dQ/de J（p,N）=DQ节点无功功率差 J（m,q）=X1;J（m,N）=DP;q=q+1; % X1=dP/de J（m,N）=DP节点有功功率差J（p,q）=X4;J（m,q）=X2; % X4=dQ/df X2=dp/dfelseif 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,

7、i）*f（i）; % dQ/de X4=-C（i）+G（i,i）*e（i）+B（i,i）*f（i）;% dQ/dfJ（p,q）=X3;%扩展列Q m=p+1;%扩展列Pelse %下面是针对PV节点来求取Jacobi矩阵的元素- % PV节点有功误差 DV=V（i）2-V2; % PV节点电压误差 % dP/de % dP/df X5=0;X6=0;J（p,q）=X5;J（p,N）=DV; % PV节点电压误差J（p,q）=X6; % PV节点有功误差 X5=-2*e（i）; X6=-2*f（i）;%以上为求雅可比矩阵的各个元素及扩展列的功率差或电压差- for k=3:N0 % N0=2*n

8、 （从第三行开始，第一、二行是平衡节点） k1=k+1;N1=N; % N=N0+1 即 N=2*n+1扩展列P、Q 或U for k2=k1:N1 % 从k+1列的Jacobi元素到扩展列的P、Q 或U J（k,k2）=J（k,k2）./J（k,k）;% 用K行K列对角元素去除K行K列后的非对角元素进行规格化 J（k,k）=1; % 对角元规格化K行K列对角元素赋1 %回代运算- if k=3 % 不是第三行 k 3 k4=k-1; for k3=3:k4 % 用k3行从第三行开始到当前行的前一行k4行消去N1 % k3行后各行上三角元素 J（k3,k2）=J（k3,k2）-J（k3,k）*

9、J（k,k2）;%消去运算（当前行k列元素消为0） end %用当前行K2列元素减去当前行k列元素乘以第k行K2列元素 J（k3,k）=0; %当前行第k列元素已消为0 if k=N0 %若已到最后一行break; %前代运算- for k3=k1:N0 % 从k+1行到2*n最后一行N1 % 从k+1列到扩展列消去k+1行后各行下三角元素%消去运算 else %是第三行k=3 %第三行k=3的前代运算-N0 %从第四行到2n行（最后一行）N1 %从第四列到2n+1列（即扩展列）%消去运算（当前行3列元素消为0） end %用当前行K2列元素减去当前行3列元素乘以第三行K2列元素 %当前行第3列元素已消为0%上面是用线性变换方式高斯消去法将Jacobi矩阵化成单位矩阵%-for k=3:2:N0-1 L=（k+1）./2; e（L）=e（L）-J（k,N）;