電力系統(tǒng)計算編程

牛頓-拉夫遜法進行潮流計算

如圖所示的一個電力網(wǎng)絡，求系統(tǒng)中的潮流分布。計算精度要求各節(jié)點電壓修正量

不大于10-5o

G

G

134

1.250-j3.7510.000-j30.000

5.00-j15.00

1.250-j3.75050

1.667-j5.000

1.667-j5.000

2.500-j7.500

25

■

已知該系統(tǒng)中，節(jié)點1為平衡節(jié)點，保持S'LM+A為定值，其余四個節(jié)點都

是PQ節(jié)點，給定的注入功率分別為：

3t-0.20+J0.2Q3--0.45-JO.15--QL40-J0.05M--0.60-J0.10

計算步驟：

1.形成節(jié)點導納矩陣YB

2.計算各節(jié)點功率的不平衡量

運用如下公式：

稽嘰筑GRq鏟-4泗?泗qT+4號，]

Q?！龉磴魆曖-4加-f”+用e*]

計算各節(jié)點功率.③、Q?

于是△甲6/；

AQ<\QY

3.計算雅可比矩陣中的各元素

先計算各節(jié)點注入電流：

月⑻一雙⑼上⑺一雙⑼

凡⑹工力嚴

然后計算雅可比矩陣各元素：

(1)當iWj時雅可比矩陣的各個元素分別為：

耳=-馬0+q6；

=-qq-4fi；

G=-耳?+q力；

(2)當i=j時雅可比矩陣的各個元素分別為：

國="弓

J.+4；

1?

列出雅可比矩陣J

4.解修正方程式求各節(jié)點電壓

代，以此循環(huán)往復，直到滿足各節(jié)點電壓修正量不大于10-5的要求。

5.計算平衡節(jié)點功率5*1和線路功率與工

平衡節(jié)點功率務務：

?J.2?J-J5

氏Z（%Z（..

工工二……nUj）

線路功率曷羽：

與禹=（Z應一"）Nu

，禹=（區(qū)一%庇-z）%L

6.將各節(jié)點的電壓以極坐標的形式表示

用C語言編程的程序如下：

#includc

#include

doubledivRe(doublebl,doubleb2,doubleb3,doubleb4

(

doubleair;

alr=(b1*b3+b2*b4/(b3*b3+b4*b4;

)

doubleMax(doubleaf],intn

{inti;

doublemax;

for(i=0;i

if(a[i]>a[i+l]

{max=a[i];a[i]=a[i4-l];a[i+l]=max;}

rcturn(max;

}〃找尋△所國和^fi[i]中的最大值

voidmain(

(

inti,j,k,h,km;

doublet;

double

cps,sumpi1,sumpi2,sumqi1,sumqi2,max,sumir,sumii,sumid,sumiq,11r,11i,l1d,l1

doublepi0[5],qi0[5],detpi[5],detqi[5],aii0[5],bii0[5],J0[8][8],detsi[8],detui[8],J0ni[8|[8]

,H[4][4|,N[4][4]J[4][4],L[4][4],euij[5][5],fuij[5][5];

staticdoubleybr[5][5]={{6.250,-5.000,-1.250,0,0},{-5.000,10.834,-1.667,-1.667,-2.500},

{-1.250,-1.667,12.917,-10.000,0),{0,-1.667,-10.000,12.917,-1.250},

{0,-2.500,0,-1.250,3.750)};

staticdoubleybi[5][5]=[(-18.750,15.000,3.750,0,0},(15.000,-32.500,5.000,5.000,7,500),

{3.750,5.000,-38.750,30.000,0),(0,5.000.30.000.-38.750,3.750),

{0,7.500,0,3.750,-11.2501};

doubleei0[5]={1.06,1.0,1.0,1.0,1.0);

doublefi0[5]={0,0A0,0);

doublepi|5|={0,0.2,-0.45,-0.4,-0.6};

doubleqi[5]={0,0.2,-0.15,-0.05,-0.1};

k=0;

km=6;

eps=0.00001;

for(i=l:i<5:i++

{printf("eiO[%d]=%f",i+LeiO[i];

printf("fi1[%d]=%f\n'\i+1,fi0[i];

)

do

{k+=1;

printf("Nowstart..An";

printf("The%dtimes\n",k-l;

sumpi2=0;

sumqi2=0;

for(i=l;i<5;i++〃第一個節(jié)點為平衡節(jié)點，從第二個節(jié)點開始求pi0[i]

{for(j=0;j<5;j++

(

sumpil=(eiO[i]*(ybr[i]U]*eiO|j]-ybi[i皿*血5+砌單挑巾皿*fiO[j]+ybi[i][j]*eiOU];

sumpi2+=sumpil;

1

pi0[i]=sumpi2;

printf("piO[%d]=%fn,i+l,pi0[i];

sumpi2=0;

}〃完成piOE的求解

for(i=l;i〈5;i++〃第一個節(jié)點為平衡節(jié)點，從第二個節(jié)點開始求qiO[i]

{for(j=0;j<5;j++

(

sumqil=(fiO[i]*(ybr[i][j]*eiO[j]-ybi[i][j]*fiO[j]-eiO[i]*(ybr[i][j]*fiO[j]+ybi[i][j]*eiO|j|;

sumqi2+=sumqil;

)

qi0|i]=sumqi2;

printf(,,qiO[%d]=%fu,i+l,qiO[i];

sumqi2=0;

}〃完成qiO[i]的求解

printf("\n";

for(i=l;i<5;i++

(

detpi[i]=pi[i]-piO[i];

detqi[i]=qi[i]-qiO[i];

printf(**detpi[%d]=%f",i+l,detpi[i];

printf("dctqi(%d]=%f\n",iII,dctqi(i];

}〃計算△piO[i],△qiO|i|

fbr(i=l;i<5;i++

{aiiO[i]=divRc(piO[i],qiO[i],ciO[i],fiO[i];

biiO[i]=divIm(piO[i],qiO[i],ciO[i],fiO[i];

printf("aiiO[%d]=%fH,i,aiiO[i];

printf(HbiiO[%d]=%f\nH,i,biiO[i];

}〃求解aiiO[i],biiO|i|

for(i=0;i<4;i++

{for(j=0;j<4;j++

if(i==j{

H[i][j]=-ybi[i+l][j+l]*eiO[i+l]+ybr[i+l][j+l]*fiO[i+l]+biiO[i+l];

N[i][j]=ybr[i+l][j+l]*eiO[i+l]+ybi[i+l]|j+l]*fiO[i+l]+aiiO[i+1];

J[i][j]=-ybr[i+l][j+l]*eiO[i+l]-ybi[i+l][i+l]*fiO[i+l]+aiiO[i+l]；

L[i][j]=-ybi[i+l][j+l]*eiO[i+l]+ybr[i+l][j+l]*fiO[i+l]-biiO[i+1]:

1

else(

=ybr[i+l][j+l]*fiO：i+l|-ybi[i+l][j+l]*eiO[i+l];

N[i][j]=ybr[i+l][j+l]*eiO[i+l]+ybi[i+l][j+l]*fiO[i+l];

J[i][j]=-ybi[i+l]|j+l]*fiO[i+l]-ybr[i+l]|j+ll*eiO|i+l|;

L[i][j]=ybr[i+1][j+1]*fiO[i+1]-ybi[i+1][j+1]*eiO[i+1];

)

}〃求取雅可比矩陣中的各元素

for(i=0;i<8;i++

for(j=0;j<8;j++

(

if(i%2==0&&j%2==0J0j]|j]=H[i/2]fj/2|;

elseif(i%2==0&&j%2!=0J0[i][j]=N[i/2][(j-l/2];

elseif(i%2!=0&&j%2==0J0[i][j]=J[(i-l/2]|j/2];

elseJ0[i]U]=L[i/2][G-l/2];

}〃形成雅可比矩陣

for(i=0;i<8;i++

for(j=0;j<8;j++

printf("JO[%d]|%d]=%3.3f,,,i,j,J0[i]|j|;

printf(”\n”；

}〃顯示雅可比矩陣JOMJ]

for(i=0;i<8;i++

{if(i%2==0detsi[i]=dctpi[(i+2/2];

elsedetsi[i]=detqi[(i+l/2];

printf("detsi[%d]=%f\n",i,detsi|i];

}〃形成的組合矩陣

for(i=0;i<8;i++〃雅可比矩陣JO口皿求逆矩陣

for(j=0;j<8;j++

(

if(i!=j

JOni[i]|j]=O;

else

JOni|i]|j|=l;

}〃逆矩陣JOni[iMj]初始化，單位矩陣

for(i=0;i<8;i++〃將雅可比矩陣JO|i||j]化簡為對角陣逆矩陣跟隨變化

1

for(j=0;j<8;j++

(

if(i!=j

t=JO[j][i]/JO[i][i];

for(h=0;h<8;h++

JO[j][h]-=JO[i][h]*t;

JOni[j][h]-=JOni[i][h]*t;

)

)

)

}

for(i=0;i<8;i++〃雅可比矩陣對角元素化為1

if(JO|i||i|!=l

(

t=JO[i][i]；

for(j=0;j<8;j++

JOni[i][j]=JOni[i]U]/t;

}〃求得逆矩陣JOni[i][j]

for(i=0;i<8;i++

{for(j=0y<8;j++

{printf(',J0ni[%d][%d]=%1.5fH,i,j,J0ni[i][j];}

printf(n\nH;

}〃顯示逆矩陣JOni[i][j]

for(i=0;i<8;i++

(

dctui[i]=O;

for(j=0;j<8;j++

detui[i]-=JOni[i][j]*detsi[i];

detui[i]*=(-l;

for(i=0;i<8;i++

{printf("detui[%d]=%f*i,detui[i];}〃求得△ui|i)

printf("\n";

for(i=l;i<5;i++

(

fiO[i]+=detui[2*i-2];

ei0|i]+=detui[2*i-l];

}〃求得每一次迭代后各節(jié)點的電壓

for(i=l;i<5;i++

{printf("eiO[%d]=%f",i,eiO(i];

printf("fi1[%d]=%fui",i,fiO[i];

)

for(i=l;i<5;i++

(pi[i]=detpi[i]+piO[i];

qi[i]=dctqi[i]+qiO[i];

)

max=Max(dctui,8;

printf("max=%l\n'\max;

)while(max>eps&&k判斷電壓修正量是否不大于10-5

printf("Alldo%dtimes\nu,k;

sumir=0;

sumii=0;

fbr(i=O;i<5;i++

11r=mu!Re(ybr[O][i],-ybi[0][i],ei0[i];

Hi=muHm(ybr[O][i],-ybi[O][i],eiO[i],.fiO[i];〃復數(shù)運算第一步,兩個復數(shù)相乘

sumir+=Ilr;

sumii+=Ili;

)

pi[O]=mulRe(eiO[O],fiO[O],sumir,sumii;

qi[O]=mul!m(eiO[O],fiO[O],sumir,sumii;〃復數(shù)運算第二步

printf("Sl=%f+j%f\n”,pi[O],q