《牛顿拉夫逊潮流程序稀疏矩阵》由会员分享，可在线阅读，更多相关《牛顿拉夫逊潮流程序稀疏矩阵（10页珍藏版）》请在金锄头文库上搜索。
1、牛顿拉夫逊潮流程序稀疏矩阵牛顿拉夫逊潮流程序稀疏矩阵function zll()%=原始数据录入模块=%error=1e-4;max_n=7;balance_power_node=1.05+0i; %平衡节点电压line=1120.01920.05750.02642130.04520.18520.0204;3240.057 0.17370.0184;4250.04720.19830.0209;5260.05810.17630.0187;6340.01320.03790.0042;7460.01190.04140.0045;812 14 0.12310.25590;912 15 0.06620.
2、13040;1012 16 0.09450.19870;11570.046 0.116 0.0102;12670.02670.082 0.0085;13680.012 0.042 0.0045;1410 17 0.03240.08450;1510 20 0.09360.209 0;1610 21 0.03480.07490;1710 22 0.07270.14990;18628 0.01690.05990.0065;19828 0.06360.20.0214;2014 15 0.221 0.19970;2115 18 0.107 0.21850;2215 23 0.10.202 0;2316
3、17 0.08240.19320;2418 19 0.06390.12920;2519 20 0.034 0.068 0;2621 22 0.01160.02360;2722 24 0.115 0.179 0;2823 24 0.132 0.270;2924 25 0.18850.32920;3025 26 0.25440.380;3125 27 0.10930.20870;3227 29 0.21980.41530;3327 30 0.32020.60270;3429 30 0.23990.45330;35911 0 0.2080;36910 00.110;3712 13 00.140;38
4、6900.208 0;39610 00.556 0;40412 00.256 0;4128 27 00.396 0;node=1 1.050.9877-0.064900;21.0207000.024 0.012;31.0135000.076 0.016;41.011 0000;51.0531000.228 0.109;61.03750000;71.0339000.058 0.02;81.0458000.112 0.075;91.031 000.062 0.016;10 1.0263000.082 0.025;11 1.0335000.035 0.018;12 1.0284000.090.058
5、;13 1.0167000.032 0.009;14 1.0142000.095 0.034;15 1.0183000.022 0.007;16 1.0216000.175 0.112;17 1.02220000;18 1.016 000.032 0.016;19 1.0108000.087 0.067;20 1.00930000;21 0.9915000.035 0.023;22 1.01710000;23 1.00690000;24 0.9972000.024 0.009;25 0.9856000.106 0.019;26 1.045 0.80.41710.217 0.127;27 1.0
6、10 0.50.166 0.942 0.19;28 1.010 0.20.29310.30.3;29 1.050.20.067100;30 1.050.20.033700;%=节点导纳矩阵形成模块=%mn=size(line); %获得参数矩阵的列 mn(1)和行 mn(2)sum_node=max(max(line(:,2),max(line(:,3); %获得总的节点数目ybus=zeros(sum_node,sum_node);c1=linspace(0,0,mn(1);for i=1:mn(1); y(i)=0;if line(i,2) line(i,3);temp=line(i,2)
7、;line(i,2)=line(i,3);line(i,3)=temp;end;y(i)=complex(line(i,4),line(i,5);y(i)=1/ y(i);%线路阻抗化为导纳c1(i)=complex(0,line(i,6); %对地 1/2电纳 C1 化为复数end;for j=1:sum_node;y_temp=0;for i=1:mn(1);if line(i,2)=j|line(i,3)=j; %如果是对角上的参数则相加，因为行小列大y_temp=y_temp+y(i)+c1(i);end;end;ybus(j,j)=y_temp; %填充对角线end;ybus(10,
8、10)=ybus(10,10)+0.19i;ybus(24,24)=ybus(24,24)+0.04i; %对地并联电容 2 个for i=1:mn(1);if line(i,2)=0 ybus(line(i,2),line(i,3)=-y(i); %填充对角线上面的元素endend;ybus=conj(ybus)+ybus-tril(ybus,0); %-参数再处理-%duidi=linspace(0,0,sum_node);for i=1:mn(1); if line(i,2)=0; duidi(line(i,3)=y(i);end;end;for i=1:sum_node;power(n
9、ode(i,1)=complex(node(i,5),node(i,6); end;star_node=2;for i=2:sum_nodeif node(i,3)=0end_node=node(i,1); %由于没有涉及 PV 节点，暂时用此式表示break;endendend_node=end_node-1;voltage=linspace(1.05+0i,1.05+0i,sum_node); %PQ 节点voltage(1)=balance_power_node; %规定第一个节点为平衡节点pv_star=end_node+1;for i=pv_star:sum_nodeui(i)=co
10、mplex(node(i,2),0);end%-牛顿拉夫逊算法实现-%for diedai=1:max_n;diedai%注入功率计算e=real(voltage);f=imag(voltage);G=real(ybus);B=imag(ybus);for i=1:sum_node;for j=1:end_node ;temp_P(i,j)=e(i)*(G(i,j)*e(j)-B(i,j)*f(j)+f(i)*(G(i,j)*f(j)+B(i,j)*e(j);temp_Q(i,j)=f(i)*(G(i,j)*e(j)-B(i,j)*f(j)-e(i)*(G(i,j)*f(j)+B(i,j)*e
11、(j);end;P(i)=sum(temp_P(i,:);Q(i)=sum(temp_Q(i,:);end;%注入功率误差计算power_n=P+Q*(1i);d_power=power-power_n; ; %功率不平衡量d_power_real=real(d_power);d_power_imag=imag(d_power);d_power_real;%注入电流计算for i=pv_star:sum_nodeui2_n(i)=sqrt(abs(ui(i)-sqrt(e(i)-sqrt(f(i);end voltage_conj=conj(voltage);power_n_conj=conj
12、(power_n);for i=1:sum_node;current(i)=power_n_conj(i)/voltage_conj(i);end;%雅克比修正矩阵a=real(current);b=imag(current);k=1;for i=star_node:end_node; %PQ 节点部分 d_s_T(k:k+1,1)=d_power_real(i);d_power_imag(i);voltage_T(k:k+1,1)=f(i);e(i); k=k+2;m=(i-star_node)*2+1;for j=star_node:sum_node;if i=jH(i,i)=- B(i,
13、i) * e(i)+ G(i,i) *f(i)+b(i);N(i,i)= G(i,i)*e(i)+B(i,i)*f(i)+a(i);J(i,i)=-G(i,i)*e(i)-B(i,i)*f(i)+a(i);L(i,i)=-B(i,i)*e(i)+G(i,i)*f(i)-b(i);R(i,i)= 2*f(i);S(i,i)= 2*e(i);elseH(i,j)=-B(i,j)*e(i)+G(i,j)*f(i);N(i,j)= G(i,j)*e(i)+B(i,j)*f(i);J(i,j)=-N(i,j);R(i,j)= 0;L(i,j)= H(i,j);S(i,j)= 0; end;n=(j-st
14、ar_node)*2+1;Jacobian(m:(m+1),n:(n+1)=H(i,j),N(i,j);J(i,j),L(i,j);end;end;for i=pv_star:sum_node; %PQ 节点部分d_s_T(k:k+1,1)=d_power_real(i);ui2_n(i);voltage_T(k:k+1,1)=f(i);e(i); k=k+2;m=(i-star_node)*2+1;for j=star_node:sum_node;if i=jH(i,i)=-B(i,i)*e(i)+G(i,i)*f(i)+b(i);N(i,i)= G(i,i)*e(i)+B(i,i)*f(i)+a(i);J(i,i)=-G(i,i)*e(i)-B(i,i)*f(i)+a(i);L(i,i)=-B(i,i)*e(i)+G(i,i)*f(i)-b(i);R(i,i)= 2*f(i);S(i,i)= 2*e(i);elseH(i,j)=-B(i,j)*e(i)+G(i,j)*f(i);N(i,j)= G(i,j)*e(i)+B(i,j)*f(i);J(i,j)=-N(i,j);R(i,j)= 0;L(i,j)= H(i,j);S(i,j)= 0; end;n=(j-star_node)*2+1;Jacobian(m:(m+1),n:
