电力系统程序

1、%本程序的功能是用牛顿拉夫逊法进行潮流计算 B1矩阵：1、支路首端号；2、末端号；3、支路阻抗；4、支路对地电纳5支路的变比；6K1，10% B2矩阵:1、该节点发电机功率；2、该节点负荷功率;3、节点电压初始值%clear；4PV 节点电压 V 的给定值；5、节点所接的无功补偿设备的容量6节点分类标号n=11;%input(请输入节点数：n=); nl=11；%input（请输入支路数：nl=); isb=1；%input（请输入平衡母线节点号：isb=)； pr=0。00001；input(请输入误差精度：pr=); B1=1217+40.2i0.000278i10；1413.6+32.1

2、6i。0002224i10;1615.3+36。18i0.0002502i10;1817+40.2i。000278i10;231.78+53.885i00.9091；2411.9+28。14i0.0001946i11;45。49+40.334i00。9091;67。781+53。885i00.9091;6813.6+32.16i。0002224i10；89。494+40.339i00.9091；8B2=000。i810；t（请输入由支路参数形成的矩阵：023123101;B1=)；0022000060+37.184i220002;0022000040+24.789i220002；0022000

3、060+37。184i220002;00220002;040+24.789i220002;00110请输入各节点参数形成的矩阵： ； n；1，)；s（1，n； S1=zeros(nl)；%修改部分 ym=1；SB=100；UB=220;ym=input（您输入的参数是标么值?(若不是则输入一个不为零的数值）; ifym=0(请输入功率基准值：；t（请输入电压基准值：=；YB=SB。/UB。/UB；BB1=B1；BB2=B2;fori=1：nl13）i，3）44end1矩阵=； disp（B1)fori=1:ni1）2（iB； 22）2。2i，3）B； 2i，4）B； 2i，)。B；end（2矩

4、阵=）; disp(B2） % -fori=1：nl%支路数ifB1（i，6)=0%左节点处于低压侧elsei，1；2;i，2）1（1；侧endY(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，1（i，4）2；%对角元1侧%求导纳矩阵disp（Y=)；disp（Y)-Y）g（Y；分解出导纳阵的实部和虚部fori=1:n%给定各节点初始电压的实部和虚部（i，3）;f（i）=imag(B2(i，3）);V(i）=B2(i，4);PV节点电压给定模值endfori=1:n给定各节点注入功率S(i)=B2(i,1

5、）-B2（i,2);%i节点注入功率 SG-SLB（i，i）=B（i,i）+B2(i，5);%i节点无功补偿量end%=P=real(S）;Q=imag(S); whileIT2=0；1； fori=1：nifi=isbD（forj1=1:n（i)1%（11（1（;jj）endi（i（i节点功率P计算i（j(Gijfj+Bijej）i（ii（)节点功率Q计算i（jei(Gij*fj+Bij*ej)%P,Qf（2；电压模平方= 以下针对非PV节点来求取功率差及 Jacobi矩阵元素 =ifB2(i,6）=3%非PV节点1；节点有功功率 ;%节点无功功率差以上为除平衡节点外其它节点的功率计算 =

6、Jacobi矩阵 =forj1=1：nifj1=isb&j1=i非平衡节点&非对角元)i）B（i，);%B)e（i）G1）f（)；%2%f ;%e i21pN）1；（mJ1；J（N）P；1；4J（2；elseifj1=i&j1=isb非平衡节点对角元Ci，i）e)f（)； i）i）f（i；%DB（e)i）f（i）; =i，i）f（; 1q）3；J（N）;%扩展列Q m=p+1;else=X1；q=q+1；J（p，q）=X4;J(m,N)=DP;扩展列P=X2；%= PVJacobi矩阵的元素 =;%V节点有功误差2；%V节点电压误差forj1=1:nifj1=isbj1=i%非平衡节点非对角元

7、G1）e（i）B（1）f（;%1）e（i）1）)；%X5=0;X6=0;J（pN）mJq）1；P1；J（p，6；2；elseifj1=i&j1=isb%非平衡节点对角元=i）i，ii，i）i）;%i，i）e（i）G（i）f（i）;%=e（i；2)；21；p，J（V；1；J（m，N）J（6；endend=X2；以上为求雅可比矩阵的各个元素 =fork=3：N0N0=2n从第三行开头,第一、二行是平衡节点）N 1即 1扩展列P、Qfork2=k1:N1扩展列P、Q2）k，2k；%非对角元规格化endJk，k）1；%对角元规格化ifk=3=fork3=3:k4% k3k4行消去fork2=k1：N1

8、% k3行后各行下三角元素（J（2）3，)J（k，2）;消去运算endendifk=N0break；end=for k3=k1：N0 fork2=k1:N1end3elsefor k3=k1:N0 fork2=k1:N132）3，2）J（k）2）;%消去运算endend（3，0；%=上面是用线性变换方式将 Jacobi 矩阵化成单位矩阵=for k=3:2：N01=1）L，N;修改节点电压实部k1=k+1;Lf（)J（N）;%修改节点电压虚部end-修改节点电压-fork=3:N0sJkN）ifDET=pr电压偏差量是否满足要求1不满足要求的节点数加1endendend（w=-JV”disp（

9、迭代次数:）; disp（ICT1）;p（disp（ICT2);fork=1：n（k）k）f（k）2； E（k）=e(k）+f(k）j；end计算各输出量 =p（各节点的实际电压标幺值E为(节点号从小到大排列)disp(E)；EE=E*UB；disp（EE);disp(-）;p（各节点的电压大小V为(节点号从小到大排列)）; disp(V)；VV=VUB;disp(VV);disp（-）;disp（各节点的电压相角 sida 为（节点号从小到大排列)：）;disp(sida);forp=1:nC(p）=0;forq=1:nCp）C（p）j（Y（pq）q）end p（各节点的功率S为（节点号从小

10、到大排列):；S；p（）; SS=SSB;disp（SS)；disp（- );p（各条支路的首端功率i为（挨次同您输入1时全都）; for i=1：nli1）ifB1（i,6）=0q（)（)（i4（)（i，)）q）j（1。i，3）)）；p，q；else（p)（4（)（5）q）j（i5）)）; p，q；end disp（Si（p，q） ; S（，r（)r（q）,）=i（disp(ZF);disp(SSi(p，q）;（-)；enddisp（各条支路的末端功率 Sj为(挨次同您输入 B1 时全都)：）;for i=1：nl（1；i，2; ifB1(i，6）=0jp）)jE（q）j（i，4。2）j（E

11、（q。（i,5)）conj（E(p）)*conj(1。/(B1（i,3)B1(i，5)；yj（)；elseSj（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)）; yi）q，)；end SSj(q，p)=Sj（q，p)SB;S（r（q）,r（p,=，r（)）； （F；%disp（SSj(q，p);disp(-)；endp（各条支路的功率损耗S为（挨次同您输入1时全都：）; for i=1:nl11（i，)； ）i（p，q，p； disp(DS（i） ; DDS(i）=DS(i)*SB；S（r（,，r（q,）=r（； disp（ZF）;%disp（DDS（i）)；disp（-）;end subplot（3，3，1）;plot(V);xlabel(节点号);ylabel(电压标幺值)； grid on;plot（sida）;(）l（电压角度grid on；subplot（3，3，3)； bar(real(S）); ylabel（节点注入有功)； grid on;subplot(3,3,4); r（)； ylabel(节点注入无功)；grid on; subplot（3，3，5)； r（)； ylabel（支路首端有功);gr