91夜夜人人揉人人捏人人添-91一区二区三区四区五区-91伊人久久大香线蕉-91在线电影-免费a网址-免费v片网站

當前位置:網站首頁 >> 作文 >> 最新電力系統潮流分析報告(五篇)

最新電力系統潮流分析報告(五篇)

格式:DOC 上傳日期:2023-01-11 15:25:32
最新電力系統潮流分析報告(五篇)
時間:2023-01-11 15:25:32     小編:zdfb

在經濟發展迅速的今天,報告不再是罕見的東西,報告中提到的所有信息應該是準確無誤的。那么報告應該怎么制定才合適呢?下面我給大家整理了一些優秀的報告范文,希望能夠幫助到大家,我們一起來看一看吧。

電力系統潮流分析報告篇一

#include “stdafx.h” #include

#include

#include

#include“complex.h” #include“wanjing.h” #include“gauss.h” using namespace std;

int _tmain(int argc, _tchar* argv[]){

int i;

//i作為整個程序的循環變量

int n=bus::scanfbusno();//輸入節點個數

int l=line::scanflineno();//輸入支路個數

if((l&&n)==0){return 0;} //如果找不到兩個文件中的任意一個,退出

line *line=new line[l];//動態分配支路結構體

line::scanflinedata(line);//輸入支路參數

line::printflinedata(line,l);//輸出支路參數

bus *bus=new bus[n];//動態分配結點結構體

for(int i=0;i

bus[i].=0;

bus[i].=0;}

bus::scanfbusdata(bus);//輸入節點參數

bus::printfbusdata(bus,n);//輸出結點參數

complex **x;x=new complex *[n];for(i=0;i

bus::jisuannodednz(x,line,bus,l,n);//計算節點導納矩陣

bus::printfnodednz(x,n);//輸出節點導納矩陣

int nn=(n-1)*2;double **jacaug;jacaug=new double *[nn];for(i=0;i

double *x;x=new double[nn];int count=1;

loop:

bus::jisuannodei(x,bus,n);//計算節點注入電流

bus::jisuannodescal(x,bus,n);//計算節點功率

bus::jisuannodescal(x,bus,n);//計算節點功率

bus::jisuannodesdelta(bus,n);//計算節點功率差值

bus::printfnodescal(x,bus,n);//輸出節點功率差值

int icon=wehcon1(bus,n);//whether converbence看迭代是否結束

if(icon==1){

cout<<“icon=”<

bus::jisuanjacaug(jacaug,x,bus,n);//計算雅可比增廣矩陣 // bus::printfjacaug(jacaug,n);

gauss::gauss_slove(jacaug,x,nn);//解方程組求出電壓差值

bus::revisenodev(bus,x,n);//修正節點電壓

// bus::printfnodev(bus,n);

count++;

goto loop;}

else

{

for(i=0;i

{

int statemp,endtemp;

complex aa,bb,cc,dd,b;

=0;

=-line[i].b;

statemp=line[i].start;

endtemp=line[i].end;

aa=complex::productcomplex(complex::getconj(bus[statemp-1].v), b);

bb=complex::subcomplex

(complex::getconj(bus[statemp-1].v), complex::getconj(bus[endtemp-1].v));

cc=complex::productcomplex(bb , complex::getconj(line[i].y));

dd=complex::caddc(aa,cc);

line[i].stoe=complex::productcomplex(bus[statemp-1].v,dd);

aa=complex::productcomplex(complex::getconj(bus[endtemp-1].v), b);

bb=complex::subcomplex

(complex::getconj(bus[endtemp-1].v), complex::getconj(bus[statemp-1].v));

cc=complex::productcomplex(bb , complex::getconj(line[i].y));

dd=complex::caddc(aa,cc);

line[i].etos=complex::productcomplex(bus[endtemp-1].v,dd);

}

cout<<“icon=”<

bus::jisuannodescal(x,bus,n);//計算節點功率

for(i=0;i

{

bus[i]. = bus[i]. + bus[i].;//發電機功率=注入功率+負荷功率

bus[i].= bus[i].+ bus[i].;

bus[i].v=complex::rec2polar(bus[i].v);

}

cout<<“====節點電壓===============發電機發出功率======”<

for(i=0;i

{

cout<<“節點”<<(i+1)<<'t';

complex::printfcomplex(bus[i].v);

coutt(bus[i].);

coutt(bus[i].);

cout<

}

cout<<“======線路傳輸功率==========”<

for(i=0;i

{

int statemp,endtemp;

statemp=line[i].start;

endtemp=line[i].end;

cout<

complex::printfcomplex(complex::comspanrea(line[i].stoe,0.01));

complex::printfcomplex(complex::comspanrea(line[i].etos,0.01));

cout<

} }

return 0;}

#include “stdafx.h” #include

#include

#include

#include“complex.h” #include“wanjing.h” #include“gauss.h” using namespace std;

int _tmain(int argc, _tchar* argv[]){

int i;

//i作為整個程序的循環變量

int n=bus::scanfbusno();//輸入節點個數

int l=line::scanflineno();//輸入支路個數

if((l&&n)==0){return 0;} //如果找不到兩個文件中的任意一個,退出

line *line=new line[l];//動態分配支路結構體

line::scanflinedata(line);//輸入支路參數

line::printflinedata(line,l);//輸出支路參數

bus *bus=new bus[n];//動態分配結點結構體

for(int i=0;i

bus[i].=0;

bus[i].=0;}

bus::scanfbusdata(bus);//輸入節點參數

bus::printfbusdata(bus,n);//輸出結點參數

complex **x;x=new complex *[n];for(i=0;i

bus::jisuannodednz(x,line,bus,l,n);//計算節點導納矩陣

bus::printfnodednz(x,n);//輸出節點導納矩陣

int nn=(n-1)*2;double **jacaug;jacaug=new double *[nn];for(i=0;i

double *x;x=new double[nn];int count=1;

loop:

bus::jisuannodei(x,bus,n);//計算節點注入電流

bus::jisuannodescal(x,bus,n);//計算節點功率

bus::jisuannodescal(x,bus,n);//計算節點功率

bus::jisuannodesdelta(bus,n);//計算節點功率差值 bus::printfnodescal(x,bus,n);//輸出節點功率差值

int icon=wehcon1(bus,n);//whether converbence看迭代是否結束

if(icon==1){

cout<<“icon=”<

bus::jisuanjacaug(jacaug,x,bus,n);//計算雅可比增廣矩陣

// bus::printfjacaug(jacaug,n);

gauss::gauss_slove(jacaug,x,nn);//解方程組求出電壓差值

bus::revisenodev(bus,x,n);//修正節點電壓

// bus::printfnodev(bus,n);

count++;

goto loop;}

else

{

for(i=0;i

{

int statemp,endtemp;

complex aa,bb,cc,dd,b;

=0;

=-line[i].b;

statemp=line[i].start;

endtemp=line[i].end;

aa=complex::productcomplex(complex::getconj(bus[statemp-1].v), b);

bb=complex::subcomplex

(complex::getconj(bus[statemp-1].v), complex::getconj(bus[endtemp-1].v));

cc=complex::productcomplex(bb , complex::getconj(line[i].y));

dd=complex::caddc(aa,cc);

line[i].stoe=complex::productcomplex(bus[statemp-1].v,dd);

aa=complex::productcomplex(complex::getconj(bus[endtemp-1].v), b);

bb=complex::subcomplex

(complex::getconj(bus[endtemp-1].v), complex::getconj(bus[statemp-1].v));

cc=complex::productcomplex(bb , complex::getconj(line[i].y));

dd=complex::caddc(aa,cc);

line[i].etos=complex::productcomplex(bus[endtemp-1].v,dd);

}

cout<<“icon=”<

bus::jisuannodescal(x,bus,n);//計算節點功率

for(i=0;i

{

bus[i]. = bus[i]. + bus[i].;//發電機功率=注入功率+負荷功率

bus[i].= bus[i].+ bus[i].;

bus[i].v=complex::rec2polar(bus[i].v);

}

cout<<“====節點電壓===============發電機發出功率======”<

for(i=0;i

{

cout<<“節點”<<(i+1)<<'t';

complex::printfcomplex(bus[i].v);

coutt(bus[i].);

coutt(bus[i].);

cout<

}

cout<<“======線路傳輸功率==========”<

for(i=0;i

{

int statemp,endtemp;

statemp=line[i].start;

endtemp=line[i].end;

cout<

complex::printfcomplex(complex::comspanrea(line[i].stoe,0.01));

complex::printfcomplex(complex::comspanrea(line[i].etos,0.01));

cout<

} }

return 0;}

#include

using namespace std;

class complex//定義復數類 { public: double real;double image;int

recpolar;//0表示直角坐標,1表示極坐標

static complex caddc(complex c1,complex c2);//求兩個復數和

static complex subcomplex(complex c1,complex c2);//求兩個復數差

static complex productcomplex(complex c1,complex c2);//求兩個復數積

static complex spanidecomplex(complex c1,complex c2);//求兩個復數商

static complex comspanrea(complex c1,double r2);//除數

static complex getconj(complex c1);//求一個復數共軛

static complex getinverse(complex c1);//取倒數

static double getcomplexreal(complex c1);//求一個復數實部

static double getcompleximage(complex c1);//求一個復數虛部

static void printfcomplex(complex c1);//顯示一個復數

static void printfmulticomplex(complex c,int n);//顯示多個復數

static void zerocomplex(complex c1);//將復數復零

static complex rec2polar(complex c1);//取極坐標

complex(){

recpolar=0;} };

complex complex::rec2polar(complex c1)//極坐標表示 { complex node;=sqrt(*+*);=atan2(,)*180/3.1415926;ar=1;return node;}

complex complex::caddc(complex c1,complex c2)//復數加法 {

complex node;

=+;

=+;

return node;}

complex complex::subcomplex(complex c1,complex c2)//復數減法 {

complex node;

=;

=-;

return node;}

complex complex::productcomplex(complex c1,complex c2)//復數乘法 {

complex node;

=**;

=*+*;

return node;} complex complex::spanidecomplex(complex c1,complex c2)//復數除法 {

complex node;

=(*+*)/(pow(,2)+pow(,2));=(*-*)/(pow(,2)+pow(,2));return node;} complex complex::comspanrea(complex c1,double r1)//復數除數 { complex node;=/(r1);=/(r1);return node;} complex complex::getconj(complex c1)//取共軛 {

complex node;

=;=-;

return node;}

complex complex::getinverse(complex c1)//取倒數 { complex node;=1;=0;node=(complex::spanidecomplex(node,c1));return node;}

double complex::getcomplexreal(complex c1)//取實部 {

return ;}

double

complex::getcompleximage(complex c1)//取虛部 {

return ;}

void

complex::printfcomplex(complex c1)//按直角坐標輸出 { if(ar==0){ ion(6);

(8);

(ios::right);

cout<

”;

ion(6);

(8);

(ios::left);

cout<

”;} else {

cout<

complex::zerocomplex(complex c1)//清零 { =0;=0;}

class gauss { public: static void gauss_slove(double **a,double *x,int nn);static void gauss_output();};

void gauss::gauss_slove(double **a,double *x,int nn){

int n,i,j,k,*pivrow,**pivcol;double eps,pivot,sum,aik,al;

n=nn;pivrow=new int[n];pivcol=new int *[n];

for(i=0;i

pivot= fabs(a[k][k]);

pivrow[k]=k;//行

pivcol[k][0]=k;pivcol[k][1]=k;//列n*2矩陣

for(i=k;i

{

for(j=k;j

{

if(pivot

{

pivot=fabs(a[i][j]);

pivrow[k]=i;//行

pivcol[k][1]=j;//列

}

}

}

if(pivot

{

cout<<“error”<

getchar();

exit(0);

}

if(pivrow[k]!=k)//行變換

{

for(j=k;j<(n+1);j++)

{

al=a[pivrow[k]][j];

a[pivrow[k]][j]=a[k][j];

a[k][j]=al;

}

}

if(pivcol[k][1]!=k)//列變換

{

for(i=0;i

{

al=a[i][pivcol[k][1]];

a[i][pivcol[k][1]]=a[i][k];

a[i][k]=al;

}

}

if(k!=(n-1))//將矩陣化為上三角形

{

for(i=(k+1);i

{

aik=a[i][k];

for(j=k;j<(n+1);j++)

{

a[i][j]-=aik*a[k][j]/a[k][k];

}

}

} } x[n-1]=a[n-1][n]/a[n-1][n-1];//解方程

for(i=(n-2);i>=0;i--){

sum=0;

for(j=(i+1);j

{

sum +=a[i][j]*x[j];0.182709

0.016894-0.0310701

-0.0402051 0.156702

-0.0355909-0.0668055

-0.00703229-0.0886481

-0.0129814-0.0390805

-0.0135062-0.1023

-0.0460568

-0.0342827

-0.00382402-0.102896

-0.0184062

}

x[i]=(a[i][n]-sum)/a[i][i];} for(k=(n-2);k>=0;k--){

al=x[pivcol[k][1]];

x[pivcol[k][1]]=x[pivcol[k][0]];

x[pivcol[k][0]]=al;}

cout<<“節點電壓修正量”<

cout<

}

====節點功率計算值==== 0.935261

-0.159048 0.573909

0.0789973-0.00289889

-0.00796623-0.0791247

-0.0168362-0.436255

-0.0580392 0.0359139

-0.0106592-0.229118

-0.0885419-0.136179

-0.148207 0.0446243

0.0111298-0.0223764

-0.00695775-0.0237482

-0.198318

-5.24266e-015

-0.0354071

-0.0925078

-1.05629e-015

-0.0391348

0.014529

0.00158644

-0.0258771

-0.109514

icon=1進行第2次迭代 節點電壓修正量

=================-0.00164889-0.000540034-0.00261067-0.00532027-0.00235315-0.00600971-0.00189677-0.00643874-0.0023631-0.00650659-0.00170949-0.0074907-0.00164545-0.00485415-0.00493977-0.0119042-0.00331285-0.0175611-0.00207908

-0.00347744-0.0869347-9.48909e-015-0.0110778-0.0538236-7.53784e-016-0.0168097 7.049e-005-0.00146487-0.00458276 0.00251645

-0.00336375-0.00530645-0.0147816-0.000326161-0.00640487-0.00251701-0.0169829-0.00175286-0.0174333-0.0239063

-0.0119192-0.076014

-0.0160104-0.441997

-0.0750285 0.000250012

3.72542e-005-0.228052

-0.108844-0.100078

-0.105634 0.000410707

0.000378067-0.057497

-0.0195879 0.200039

0.0582563-0.00307326-0.0163809-0.00232773-0.0175806 8.74293e-005-0.0192018 0.000558996-0.0197776-0.000247851-0.0193784-0.00115346-0.0185848-0.00127275-0.0186244-0.00010108-0.0188966 0.000553585-0.0200901-3.76315e-005-0.0208303 0.00308341-0.0219386-0.00195916-0.0205356-0.00184757-0.0076401 0.00197593-0.0245534 0.00434657-0.027534

====節點功率計算值==== 0.98623

-0.134163 0.583136

0.166278-0.111173

0.199792

-0.0621041

-0.0821379

-0.0350785

-0.0902383

-0.0320461

-0.0951562

-0.0220362

-0.175458

4.72557e-015

-0.0320661

-0.0871134

-7.03489e-017

-0.0350769

0.000273455

1.51804e-005

-0.0240417

-0.10604

icon=1進行第3次迭代 節點電壓修正量

=================-2.67079e-005-2.30128e-006-2.20543e-005-6.00686e-005-2.33043e-005-6.85601e-005-3.22294e-005-2.61107e-005-2.80198e-005-6.6167e-005-2.34528e-005

-0.0739846 0.0227868-0.0158709-0.0248173-0.0179447-0.0578368-0.00890719-0.0337091-0.00693706-0.111601 1.21429e-014-0.0159145-0.0667319 9.24355e-016-0.0228592 7.10354e-005-6.6188e-006-0.00889343-0.0184098

-5.66132e-005-4.4646e-005-1.74668e-005-4.50947e-005-0.000181763-3.81763e-006-0.000286581-6.68993e-005-1.28441e-005-5.17172e-005-0.000223284-4.54717e-005-2.47586e-005 4.32335e-007-0.000258494 1.82635e-005-0.000272051-6.95195e-006-0.000251969 1.11318e-005-0.000279418 5.74737e-005-0.000307368 6.86998e-005-0.000320274 5.38112e-005-0.00031447 3.59531e-005-0.00030494 3.37607e-005-0.000307449 5.26532e-005-0.000310721 6.92761e-005-0.000350373 5.60942e-005-0.00040977 0.000123641-0.000440259 1.36149e-005-0.000426973-1.70227e-005-9.37794e-005 0.000113675-0.000544011 0.000176034-0.000636202

====節點功率計算值====

0.986878

-0.133979 0.583

0.167193-0.024

-0.012-0.076

-0.016-0.442

-0.0748606

1.43501e-008

1.07366e-008-0.228

-0.109

-0.0999999

-0.104049 4.51318e-008

8.98835e-008-0.0579999

-0.0199999 0.2

0.0591018-0.112

-0.0749997 0.2

0.0242519-0.062

-0.016-0.082

-0.025-0.035

-0.018

-0.0900001

-0.058-0.032

-0.00899997-0.095

-0.0339999-0.022

-0.00699998-0.175

-0.112

-6.07156e-015

-1.19217e-014-0.032

-0.016-0.087

-0.0669999

7.03078e-017

-9.23979e-016-0.035

-0.0229999

1.09492e-007

4.45699e-008 1.54958e-009

-2.01531e-010-0.024

-0.00899994-0.106

-0.0189996

icon=0,迭代結束。

====節點電壓===============發電機發出功率======

節點1

1.05

0。

98.6878-13.3979

節點2

1.045

-1.846。

29.4193

節點3

1.02384-3.83352。

0

點25 1.01216-9.68486。

0

0 0 節點4

1.01637-4.55698。

0

點26 0.994393

-10.1089。

0 0

0 節點5

1.01

-6.48617。

節 點27 1.02012-9.42025。

0

11.5139 0 節點6

1.01332-5.38073。

0

點28 1.00992-5.86244。

0

0 0 節點7

1.00489-6.38368。

0

點29 1.00022-10.6579。

0

0 節點8 19.5951 節點9 0 節點10 0 節點11 5.91018 節點12 0 節點13 2.42519 節點14 0 節點15 0 節點16 0 節點17 0 節點18 0 節點19 0 節點20 0 節點21 0 節點22 0 節點23 0 節點24 0 1.01

-5.62974。

1.03905-6.78143。

1.03595-8.69362。

-4.5962。

1.04711-7.80323。

1.05

-6.34392。

1.03242-8.7401。

1.02788-8.86784。

1.03458-8.45044。

1.03051-8.83678。

1.01845-9.5141。

1.01604-9.70326。

1.02022-9.50938。

1.0237-9.17478。

1.02432-9.17024。

1.01802-9.36719。

1.01339-9.68362。

0 20

節 點30 0.988705

-11.5464。

0

0 0

======

線路傳輸功率========== 2to1

-57.7373

5.41674i

58.3454

0

-15.1827i

3to1

-39.659

-7.75964i

40.3424

1.78481i

4to2

-30.87

-9.74186i

31.4153

0

3.58352i

4to3

-37.0772

-7.78596i

37.259

6.55964i

5to2

-44.3717

-9.78456i

45.2968

0

4.84242i

6to2

-38.4766

-8.22625i

39.3252

0

2.87667i

6to4

-34.946

1.92384i

35.0885

0

-3.28202i

7to5

-0.16304

-6.41767i

0.171702

0

2.2985i

7to6

-22.637

-4.48233i

22.7745

0

1.44238i

8to6

-11.8939

-5.48098i

11.913

0

3.70557i

6to9

12.3737

-12.3826i

-12.3737

0

13.0033i

6to10

10.9107

-3.80907i

-10.9107

0

4.53223i

11to9

5.91018i

0

-5.08963i

10to9

-32.652

-2.3712i

32.652

0

3.46974i

4to12

23.5411

-11.5375i

-23.5411

0

13.2407i

13to12

2.42519i

1.05

-1.90978i 1.66484i 14to12

-7.9019

-2.06732i

7.97894

30to29

-3.6702

-0.542564i

3.70398

2.22749i 0.606393i 15to12

-18.254

-5.74885i

18.4835

28to8

-1.89152

-3.79982i

1.89395

6.20089i-4.9239i 16to12-7.53872

-2.90237i

7.59633

28to6

-14.7868

-2.82565i

14.8234

3.02352i 0.294601i 15to14-1.69544

-0.461488i

1.70189

請按任意鍵繼續...0.467323i 17to16-4.03014 1.10238i 18to15-6.08074 1.46028i 19to18-2.87549 0.478389i 20to19

6.6418-2.93222i 20to10

-8.8418 3.85077i 17to10-4.96987 4.76656i 21to10-16.1562 9.42843i 22to10-7.87782 4.21401i 22to21

1.34443-2.01837i 23to15-5.59369 2.25006i 24to22-6.48186 2.08163i 24to23-2.38596 0.579814i 25to24-0.167617 0.281364i 26to25

-3.5 2.3674i 27to25

3.39433-2.08638i 28to27

16.1446 3.13006i 29to27-6.10398 1.67047i 30to27-6.92979-1.07089i-1.37839i-0.467767i

2.96679i-3.66679i-4.72911i-9.18162i-4.10132i

2.01969i-2.17981i-2.00141i-0.56401i

-0.28102i-2.29999i

2.11848i-2.10093i-1.50639i

-1.3574i

4.03872

6.12096

2.88074

-6.62452

8.9242

4.98423

16.2709

7.93248

-1.34378

5.62846

6.53339

2.39369

0.167814

3.54513

-3.37751

-16.1446

6.19083

7.09313

高等電力系統分析 ieee30節點潮流程序

班級:電研114班

姓名:王大偉

學號:2201100151

電力系統潮流分析報告篇二

前 言

電力系統的潮流分布是描述電力系統運行狀態的技術術語,它表明電力系統在某一確定的運行方式和接線方式下,系統中從電源經網絡到負荷各處的電壓、電流、功率的大小和方向的分布情況。電力系統的潮流分布,主要取決于負荷的分布、電力網絡參數以及和供電電源間的關系。對電力系統在各種運行方式下進行潮流分布計算,以便確定合理的供電方案,合理地調整負荷。

迄今,電子計算機的運用已經十分普遍,而運用電子計算機計算、分析、研究電力系統時,往往離不開計算其中的潮流分布。本文將以一兩機五節點模型(如圖1-1,參數給定)為基礎,結合matlab軟件,通過牛頓拉夫遜算法和pq分解法分析潮流分布情況。

關鍵詞:潮流計算、matlab、牛頓拉夫遜法、pq分解法

i

目錄

前 言.........................i

第一章 電力系統潮流計算.........................1第一節 電力系統潮流計算簡介.......................1

第二節 電力網絡方程.........................1

第二章 復雜系統的計算機算法........................1

第一節 牛頓拉夫遜算法.....................1

第二節 牛頓拉夫遜法潮流計算的基本步驟.................3第三節 計算機算法程序(見附件)...............3

第三章pq分解法.........................3

第一節 pq分解法潮流計算時的修正方程式...............3

第二節 pq分解法潮流計算基本步驟....................4第三節 pq分解法的matlab實現(見附件)...............4

第四章 兩種算法的比較......................4

第一章 電力系統潮流計算

第一節 電力系統潮流計算簡介

電力系統潮流計算是電力系統穩態運行情況的一種計算,它是基于給定的運行條件及系

統接線方式,確定整個電力系統個部分的運行狀態。在電力系統規劃的設計和現有電力系統運行方式的研究中,通過潮流計算,可以定量的分析比較供電方案和運行方式的合理性。通過潮流計算,還可以發現系統中的薄弱環節,檢查設備、元件是否過負荷,各節點是否符合要求,以便提出必要的改進措施,實施相應的調壓措施,保證電力系統的電能質量,并使整個電力系統獲得最大的經濟性。

第二節 電力網絡方程

在電路的理論學習中,已導出了運用節點導納矩陣的節點電壓方程

電力系統潮流分析報告篇三

南 京 理 工 大 學

《電力系統穩態分析》

課程報告

姓名

xx

學 號: 5*** 自動化學院 電氣工程

基于牛頓-拉夫遜法的潮流計算例題編程報學院(系): 專

業: 題

目: 任課教師 碩士導師 告

楊偉 xx

2015年6月10號

基于牛頓-拉夫遜法的潮流計算例題編程報告

摘要:電力系統潮流計算的目的在于:確定電力系統的運行方式、檢查系統中各元件是否過壓或者過載、為電力系統繼電保護的整定提供依據、為電力系統的穩定計算提供初值、為電力系統規劃和經濟運行提供分析的基礎。潮流計算的計算機算法包含高斯—賽德爾迭代法、牛頓-拉夫遜法和p—q分解法等,其中牛拉法計算原理較簡單、計算過程也不復雜,而且由于人們引入泰勒級數和非線性代數方程等在算法里從而進一步提高了算法的收斂性和計算速度。同時基于matlab的計算機算法以雙精度類型進行數據的存儲和運算, 數據精確度高,能進行潮流計算中的各種矩陣運算,使得傳統潮流計算方法更加優化。

一 研究內容

通過一道例題來認真分析牛頓-拉夫遜法的原理和方法(采用極坐標形式的牛拉法),同時掌握潮流計算計算機算法的相關知識,能看懂并初步使用matlab軟件進行編程,培養自己電力系統潮流計算機算法編程能力。

例題如下:用牛頓-拉夫遜法計算下圖所示系統的潮流分布,其中系統中5為平衡節點,節點5電壓保持u=1.05為定值,其他四個節點分別為pq節點,給定的注入功率如圖所示。計算精度要求各節點電壓修正量不大于10-6。

二 牛頓-拉夫遜法潮流計算 1 基本原理

牛頓法是取近似解x(k)之后,在這個基礎上,找到比x(k)更接近的方程的根,一步步地迭代,找到盡可能接近方程根的近似根。牛頓迭代法其最大優點是在方程f(x)=0的單根附近時誤差將呈平方減少,而且該法還可以用來求方程的重根、復根。電力系統潮流計算,一般來說,各個母線所供負荷的功率是已知的,各個節點的電壓是未知的(平衡節點外)可以根據網絡結構形成節點導納矩陣,然后由節點導納矩陣列寫功率方程,由于功率方程里功率是已知的,電壓的幅值和相角是未知的,這樣潮流計算的問題就轉化為求解非線性方程組的問題了。為了便于用迭代法解方程組,需要將上述功率方程改寫成功率平衡方程,并對功率平衡方程求偏導,得出對應的雅可比矩陣,給未知節點賦電壓初值,將初值帶入功率平衡方程,得到功率不平衡量,這樣由功率不平衡量、雅可比矩陣、節點電壓不平衡量(未知的)構成了誤差方程,解誤差方程,得到節點電壓不平衡量,節點電壓加上節點電壓不平衡量構成節點電壓新的初值,將新的初值帶入原來的功率平衡方程,并重新形成雅可比矩陣,然后計算新的電壓不平衡量,這樣不斷迭代,不斷修正,一般迭代三到五次就能收斂。2 基本步驟和設計流程圖

形成了雅克比矩陣并建立了修正方程式,運用牛頓-拉夫遜法計算潮流的核心問題已經解決,已有可能列出基本計算步驟并編制流程圖。由課本總結基本步驟如下:

1)形成節點導納矩陣y;

2)設各節點電壓的初值,如果是直角坐標的話設電壓的實部e和虛部f;如果是極坐標的話則設電壓的幅值u和相角a;

3)將各個節點電壓的初值代入公式求修正方程中的不平衡量以及修正方程的系數矩陣的雅克比矩陣;

4)解修正方程式,求各節點電壓的變化量,即修正量; 5)計算各個節點電壓的新值,即修正后的值;

6)利用新值從第(3)步開始進入下一次迭代,直至達到精度退出循環; 7)計算平衡節點的功率和線路功率,輸出最后計算結果; ① 公式推導

② 流程圖

matlab編程代碼

clear;

% 如圖所示1,2,3,4為pq節點,5為平衡節點

y=0;

% 輸入原始數據,求節點導納矩陣

y(1,2)=1/(0.07+0.21j);

y(4,5)=0;y(1,3)=1/(0.06+0.18j);

y(1,4)=1/(0.05+0.10j);

y(1,5)=1/(0.04+0.12j);

y(2,3)=1/(0.05+0.10j);

y(2,5)=1/(0.08+0.24j);

y(3,4)=1/(0.06+0.18j);

for i=1:5

for j=i:5

y(j,i)=y(i,j);

end

end

y=0;

% 求節點導納矩陣中互導納

for i=1:5

for j=1:5

if i~=j

y(i,j)=-y(i,j);

end

end

end

% 求節點導納矩陣中自導納

for i=1:5

y(i,i)=sum(y(i,:));

end

y

% y為導納矩陣

g=real(y);

b=imag(y);% 輸入原始節點的給定注入功率

s(1)=0.3+0.3j;

s(2)=-0.5-0.15j;

s(3)=-0.6-0.25j;

s(4)=-0.7-0.2j;

s(5)=0;

p=real(s);

q=imag(s);

% 賦初值,u為節點電壓的幅值,a為節點電壓的相位角

u=ones(1,5);

u(5)=1.05;

a=zeros(1,5);

x1=ones(8,1);

x2=ones(8,1);

k=0;

while max(x2)>1e-6

for i=1:4

for j=1:4

h(i,j)=0;

n(i,j)=0;

m(i,j)=0;

l(i,j)=0;

op(i)=0;

oq(i)=0;

end

end

% 求有功、無功功率不平衡量

for i=1:4

for j=1:5

op(i)=op(i)-u(i)*u(j)*(g(i,j)*cos(a(i)-a(j))+b(i,j)*sin(a(i)-a(j)));

oq(i)=oq(i)-u(i)*u(j)*(g(i,j)*sin(a(i)-a(j))-b(i,j)*cos(a(i)-a(j)));

end

op(i)=op(i)+p(i);

oq(i)=oq(i)+q(i);

end

x2=[op,oq]';

% x2為不平衡量列向量

% 求雅克比矩陣

% 當i~=j時,求h,n,m,l

for i=1:4

for j=1:4

if i~=j

h(i,j)=-u(i)*u(j)*(g(i,j)*sin(a(i)-a(j))-b(i,j)*cos(a(i)-a(j)));

n(i,j)=-u(i)*u(j)*(g(i,j)*cos(a(i)-a(j))+b(i,j)*sin(a(i)-a(j)));

l(i,j)=h(i,j);

m(i,j)=-n(i,j);

end

end

end

% 當i=j時,求h,n,m,l

for i=1:4

for j=1:5

if i~=j h(i,i)=h(i,i)+u(i)*u(j)*(g(i,j)*sin(a(i)-a(j))-b(i,j)*cos(a(i)-a(j)));n(i,i)=n(i,i)-u(i)*u(j)*(g(i,j)*cos(a(i)-a(j))+b(i,j)*sin(a(i)-a(j)));

m(i,i)=m(i,i)-u(i)*u(j)*(g(i,j)*cos(a(i)-a(j))+b(i,j)*sin(a(i)-a(j)));

l(i,i)=l(i,i)-u(i)*u(j)*(g(i,j)*sin(a(i)-a(j))-b(i,j)*cos(a(i)-a(j)))

end

end

n(i,i)=n(i,i)-2*(u(i))^2*g(i,i);

l(i,i)=l(i,i)+2*(u(i))^2*b(i,i);

end

j=[h,n;m,l]

% j為雅克比矩陣

x1=-((inv(j))*x2);

% x1為所求△x的列向量

% 求節點電壓新值,準備下一次迭代

for i=1:4

oa(i)=x1(i);

ou(i)=x1(i+4)*u(i);

end

for i=1:4

a(i)=a(i)+oa(i);

u(i)=u(i)+ou(i);

end

k=k+1;

end

k,u,a

% 求節點注入功率

i=5;

for j=1:5

p(i)=u(i)*u(j)*(g(i,j)*cos(a(i)-a(j))+b(i,j)*sin(a(i)-a(j)))+p(i);

q(i)=u(i)*u(j)*(g(i,j)*sin(a(i)-a(j))-b(i,j)*cos(a(i)-a(j)))+q(i);

end

s(5)=p(5)+q(5)*sqrt(-1);

s

% 求節點注入電流

i=y*u'

運行結果

節點導納矩陣

經過五次迭代后的雅克比矩陣

迭代次數以及節點電壓的幅值和相角(弧度數)

節點注入功率和電流

五 結果分析

在這次學習和實際操作過程里:首先,對電力系統分析中潮流計算的部分特別是潮流計算的計算機算法中的牛頓-拉夫遜法進行深入的研讀,弄明白了其原理、計算過程、公式推導以及設計流程。牛頓-拉夫遜法是求解非線性方程的迭代過程,其計算公式為?f?j?x,式中j為所求函數的雅可比矩陣;?x為需要求的修正值;?f為不平衡的列向量。利用x(*)=x(k+1)+?x(k+1)進行多次迭代,通過迭代判據得到所需要的精度值即準確值x(*)。六 結論

通過這個任務,自己在matlab編程,潮流計算,word文檔的編輯功能等方面均有提高,但也暴漏出一些問題:理論知識儲備不足,對matlab的性能和特點還不能有一個全面的把握,對word軟件也不是很熟練,相信通過以后的學習能彌補這些不足,達到一個新的層次。

電力系統潮流分析報告篇四

電力系統潮流計算發展史

對潮流計算的要求可以歸納為下面幾點:

(1)算法的可靠性或收斂性(2)計算速度和內存占用量(3)計算的方便性和靈活性

電力系統潮流計算屬于穩態分析范疇,不涉及系統元件的動態特性和過渡過程。因此其數學模型不包含微分方程,是一組高階非線性方程。非線性代數方程組的解法離不開迭代,因此,潮流計算方法首先要求它是能可靠的收斂,并給出正確答案。隨著電力系統規模的不斷擴大,潮流問題的方程式階數越來越高,目前已達到幾千階甚至上萬階,對這樣規模的方程式并不是采用任何數學方法都能保證給出正確答案的。這種情況促使電力系統的研究人員不斷尋求新的更可靠的計算方法。

在用數字計算機求解電力系統潮流問題的開始階段,人們普遍采用以節點導納矩陣為基礎的高斯-賽德爾迭代法(一下簡稱導納法)。這個方法的原理比較簡單,要求的數字計算機的內存量也比較小,適應當時的電子數字計算機制作水平和電力系統理論水平,于是電力系統計算人員轉向以阻抗矩陣為主的逐次代入法(以下簡稱阻抗法)。

20世紀60年代初,數字計算機已經發展到第二代,計算機的內存和計算速度發生了很大的飛躍,從而為阻抗法的采用創造了條件。阻抗矩陣是滿矩陣,阻抗法要求計算機儲存表征系統接線和參數的阻抗矩陣。這就需要較大的內存量。而且阻抗法每迭代一次都要求順次取阻抗矩陣中的每一個元素進行計算,因此,每次迭代的計算量很大。

阻抗法改善了電力系統潮流計算問題的收斂性,解決了導納法無法解決的一些系統的潮流計算,在當時獲得了廣泛的應用,曾為我國電力系統設計、運行和研究作出了很大的貢獻。但是,阻抗法的主要缺點就是占用計算機的內存很大,每次迭代的計算量很大。當系統不斷擴大時,這些缺點就更加突出。為了克服阻抗法在內存和速度方面的缺點,后來發展了以阻抗矩陣為基礎的分塊阻抗法。這個方法把一個大系統分割為幾個小的地區系統,在計算機內只需存儲各個地區系統的阻抗矩陣及它們之間的聯絡線的阻抗,這樣不僅大幅度的節省了內存容量,同時也提高了節省速度。

克服阻抗法缺點的另一途徑是采用牛頓-拉夫遜法(以下簡稱牛頓法)。牛頓法是數學中求解非線性方程式的典型方法,有較好的收斂性。解決電力系統潮流計算問題是以導納矩陣為基礎的,因此,只要在迭代過程中盡可能保持方程式系數矩陣的稀疏性,就可以大大提高牛頓潮流程序的計算效率。自從20世紀60年代中期采用了最佳順序消去法以后,牛頓法在收斂性、內存要求、計算速度方面都超過了阻抗法,成為直到目前仍被廣泛采用的方法。

在牛頓法的基礎上,根據電力系統的特點,抓住主要矛盾,對純數學的牛頓法進行了改造,得到了p-q分解法。p-q分解法在計算速度方面有顯著的提高,迅速得到了推廣。

牛頓法的特點是將非線性方程線性化。20世紀70年代后期,有人提出采用更精確的模型,即將泰勒級數的高階項也包括進來,希望以此提高算法的性能,這便產生了保留非線性的潮流算法。另外,為了解決病態潮流計算,出現了將潮流計算表示為一個無約束非線性規劃問題的模型,即非線性規劃潮流算法。

近20多年來,潮流算法的研究仍然非?;钴S,但是大多數研究都是圍繞改進牛頓法和p-q分解法進行的。此外,隨著人工智能理論的發展,遺傳算法、人工神經網絡、模糊算法也逐漸被引入潮流計算。但是,到目前為止這些新的模型和算法還不能取代牛頓法和p-q分解法的地位。由于電力系統規模的不斷擴大,對計算速度的要求不斷提高,計算機的并行計算技術也將在潮流計算中得到廣泛的應用,成為重要的研究領域。

電力系統潮流分析報告篇五

電力系統潮流計算程序設計

姓名:韋應順

學號:2011021052 電力工程學院

牛頓—拉夫遜潮流計算方法具有能夠將非線性方程線性化的特點,而使用matlab語言是由于matlab語言的數學邏輯強,易編譯。

【】【】程序12

function tisco %這是一個電力系統潮流計算的程序 n=input(‘n請輸入節點數:n=’); m=input(‘請輸入支路數:m=’);ph=input(‘n請輸入平衡母線的節點號:ph=’); b1=input(‘n請輸入支路信號:b1=’);%它以矩陣形式存貯支路的情況,每行存貯一條支路 %第一列存貯支路的一個端點 %第二列存貯支路的另一個端點 %第三列存貯支路阻抗

%第四列存貯支路的對地導納

%第五列存貯變壓器的變比,注意支路為1 %第六列存貯支路的序號

b2=input(‘n請輸入節點信息:b2=’); %第一列為電源側的功率 %第二列為負荷側的功率 %第三列為該點的電壓值

%第四列為該點的類型:1為pq,2為pv節點,3為平衡節點 a=input(‘n請輸入節點號及對地阻抗:a=’); ip=input(‘n請輸入修正值:ip=’); %ip為修正值);y=zeros(n);

y(p,q)=y(p,q)-1./(b1(i3)*b1(i5);e=zeros(1,n);

y(p,q)=y(p,q);f=zeros(1,n);

no=2*ph=1; y(q,q)=y(q,q)+1./b1(i3)+b1(i4)/2;

end for i=1:n

g=real(y);if a(i2)=0

b=imag(y);p=a(i1);

y(p p)=1./a(i2);for i=1:n end e(i)=real(b2(i3));end f(i)=imag(b2(i3));for i=1:m s(i)=b2(i1)-b2(i2);p=b1(i1);v(i)=b2(i3);p=b1(i2);end y(p,p)=y(p,p)+1./(b1(i3)*b1(i5)^2+b1(i4)./2p=real(s);q=imag(s);[c,d,df]=xxf(g,b,e,f,p,q,n,b2,ph,v,no);j=jacci(y,g,b,p,q,e,f,v,c,d,b2,n,ph,no);[de,di]=hxf(j,d,f,ph,n,no);t=0;while

max(abs(de))>ip&max(abs(dfi)>ip

t=t+1;

e=e+de;

f=f+df;

[c,d,df]=xxf(g,b,e,f,p,q,n,b2,ph,v,no);

j=jacci(y,g,b,p,q,e,f,v,c,d,b2,n,ph,no);

[de,df]=hxf(j,df,ph,n,no);end v=e+f*j;for i=1:n hh(i)=conj(y(ph,i)*v(i));end s(ph)=sum(hh)*v(ph);b2(ph,1)=s(ph);v=abs(v);

jd=angle(v)*180/p;resulte1=[a(:,1),real(v),imag(v),v,jd,real(s’),imag(s’),real(b2(:1)),imag(b2(:1)),real(b2(:2)),imag(b2(:,2))];for i=1:m

a(i)=conj((v(b1(i1))/b1(i5)-v(b1(i2))/b1(i3));

b(i)=v(b1(i1))*a(i)-j*b1(i4)*v(b1(i))^2/2;

c(i)=-v(b1(i2))*a(i)-j*b1(i4)*v(b1(i2))^2/2;end result2=[b1(:,6),b1(:,1),b1(:,2),real(b’),imag(b’),real(c’),imag(c’), real(b’+c’),imag(b’+c’)];printcut(result1,s,b,c,result2);type resultm function [c,d,df]=xxf(g,b,e,f,p,q,n,b2,ph,v,no)%該子程序是用來求取df for i=1:n

if

i=ph

c(i)=0;

d(i)=0;

for j=i:n

c(i)=c(i)+g(i,j)*e(j)-b(i,j)*f(j);d(i)=d(i)+g(i,j)*f(j)+b(i,j)*e(j);end

p1=c(i)*e(i)+d(i)*f(i);q1=c(i)*f(i)-d(i)*e(i);v1=e(i)^2+f(i)^2;if

b2(i4)=2 p=2*i-1;

df(p)=p(i)-p1;p=p+1;else p=2*i-1;

df(p)=p(i)-p1;p=p+1;

df(p)=q(i)-q1;end end end df=df’;if ph=n df(no?=[];end

function [de,df]=hxf(j,df,ph,n,no)%該子函數是為求取de df dx=jdf;dx1=dx;

x1=length(dx1);if ph=n dx(no)=0;dx(no+1)=0;

for i=(no+2):(x1+2)dx(i)=dx1(i-2);end else

dx=[dx1,0,0];end k=0;

[x,y]=size(dx);for i=1:2:x k=k+1;

df(k)=dx(i);de(k)=dx(i+1);end end case 2 function for j=1:n j=jacci(y,g,b,pq,e,f,v,c,d,b2,n,ph,no)x1=g(i,j)*f(i)-b(i,j)*e(i);

x2=g(i,j)*e(i)+b(i,j)*f(i);%該子程序是用來求取jacci矩陣

for i=1:n x3=0;switch b2(i4)x4=0;case 3 p=2*i-1;continue q=2*j-1;case 1 j(p,q)=x1;for j=1:n m=p+1;if

j=&j=ph j(m,q)=x3;x1=g(i)*f(i)-b(i,j)*e(i);q=q+1;x2=g(i,j)*e(i)+b(i,j)*f(i);j(p,q)=x2;x3=-x2;j(m,q)=x4;x4=x1;x1=d(i)+g(i,j)*f(i)-b(i,j)*e(i);p=2*i-1;x2=c(i)+g(i,j)*e(i)+b(i,j)*f(i);q=2*j-1;x3=0;j(p,q)=x1;x4=0;m=p+1;p=2*i-1;j(p,q)=x2;q=2*j-1;j(m,q)=x4;j(p,q)=x1;else if j=&j=jph m=p+1;x1=d(i)+g(i,j)*f(i)-b(i,j)*e(i);j(m,q)=x3;x2=c(i)+g(i,j)*e(i)+b(i,j)*f(i);q=q+1;x3= c(i)+g(i,j)*e(i)-b(i,j)*f(i);j(p,q)=x2;x4= c(i)+g(i,j)*f(i)-b(i,j)*e(i);j(m,q)=x4;p=2*i-1;end q=2*j-1;end j(p,q)=x1;end m=p+1;end j(m,q)=x3;if ph=n q=q+1;j(no:)=[];j(p,q)=x2;j(no:)=[];j(m,q)=x4;j(:,no)=[];end j(:,no)=[];end

2實例驗證 【例題】設有一系統網絡結線見圖1,各支路阻抗和各節點功率均已以標幺值標示于圖1中,其中節點2連接的是發電廠,設節點1電壓保持u1=1.06定值,試計算其中的潮流分布,請輸入節點數:n=5 請輸入支路數:m=7 請輸入平衡母線的節點號:ph=l 請輸入支路信息:

bi=[ l 2 0.02+0.06i o l 1;1 3 0.08+0.24i 0 1 2;2 3 0.06+0.18i 0 l 3: 2 4 0.06+0.18i o l 4: 2 5 0.04+0.12i 0 l 5: 3 4 0.01+0.03i 0 l 6: 4 5 0.08+0.24i o 1 7] 請輸入節點信息:

b2=[ 0 0 1.06 3;0.2+0.20i 0 1 1;一o.45一o.15i 0 l l;一0.4-0.05i 0 l 1;一0.6—0.1i 0 1 l] 請輸入節點號及對地阻抗: a=[l 0;2 0;3 0;4 0;5 o ] 請輸入修正值:ip=0.000 0l

參考文獻

[1]陳珩.電力系統穩定分析[m].北京:中國電力出版社,2002:139—187.

[2]鄭阿奇.matlab實用教程[m].北京:電子工業出版社,2005:1-243.

[3] 束洪春,孫士云,等.云電送粵交商流混聯系統全過 程動態電壓研究[j】.中國電力,2008,4l(10):l-4. shu hong—ch吼,sun shi-yun,et a1.research on fun prc'cess dyn鋤ic voltage stabil時of hybrid ac/dc power tmnsmission system舶m yu衄an province to g啪gdong province【j】.electric power,2008,4l(10): l-4.

[4] 朱新立,湯涌,等.大電網安全分析的全過程動態仿 真技術[j】.電網技術,2008,32(22):23—28. song xin—ii,tang yof唱,et a1. full dyn鋤ic simulation for the stabilhy a眥lysis of large power system【j】.power system融ilriolo影,2008,32(22): 23.28.

[5]roytelm鋤i,shallidehpour s m.a comprehcnsivc long teml dynaiilic simulation for powcr system recovery【j】. ieee transactions 0n power systems,1994,9(3). [6] 石雩梅,汪志宏,等.發電機勵磁系統數學模型及參 數對電網動態穩定性分析結果影響的研究[j】.繼電 器,2007,35(21):22-27.

shi xue.mei,wang zlli-hon舀et a1.iksearch on the innuence of g鋤e翰to璐bascd ∞de詛iled excitation system models柚d parameters t0 power鏟id dyn鋤ic stabil時【j】.relay,2007,35(2 1):22-27.

[7] 方思立,朱方.快速勵磁系統對系統穩定的影響[j】.中 國電機工程學報,1986,6(1):20.28.

fang si.1i,zhu fang.the effbct of f弧t.respon∞

excitation system on the stability of power netwofk【j】. proceedings ofthe csee,1986,6(1):20-28.

[8] 劉?。娏ο到y穩定性及發電機勵磁控制[m】.北京: 中國電力出版社,2007.

liu qu.power system s詛bility鋤d generator excitation control【m】.beuing:chiila electric powef press,2007. [9] dallachy j l,anderson t.experience with rcplacing ro詛ting exciters wim static exciters【j】.1k institution of electrical engineers,1 996.

[10] 陳利芳,陳天祿.淺談自并勵勵磁系統在大容量機組 中的應用【j】.繼電器,2007,35(1):8l培4. chen li-f抽島chen tian—lil application of 辯l仁exci組tion mode in large capacity髫memtor unit【j】. reiay'2007,35(1):81-84.

[11] 方思立,劉增煌,孟慶和.大型汽輪發電機自并勵勵 磁系統的應用條件【j].中國電力,1994,27(12):61.63. fang si.ii,liu zeng-hu鋤g,meng qin爭hc.m application conditions of large turbine generator self-excitation system【j】.electric powef,1994,27(12): 61.63.

[12]梁小冰,黃方能.利用emtdc進行長持續時間過程 的仿真研究【j】.電網技術,2002,26(9):55.57. liang xiao-bing,huang fan爭眥ng.how to cany out simulalion of long dul‘ation processes by use of emtdc【j】.power system 11echnology,2002,26(9): 55-57.

[13]王卉,陳楷,彭哲,等.數字仿真技術在電力系統中 的應用及常用的幾種數字仿真工具【j】.繼電器,2004,32(21):7l一75.

wang hui,chen kai,peng zhe,et a1.application of digital simulation眥hniques棚d severaj simulation tools in power system[j】.relay,2004,32(21):71·75.

[14]ieee power engmeering socie哆.ieee std 421.5.2005 ieee玎ccommended practice for excitation system models for power system stabii時studies【s】.

全文閱讀已結束,如果需要下載本文請點擊

下載此文檔
a.付費復制
付費獲得該文章復制權限
特價:5.99元 10元
微信掃碼支付
已付款請點這里
b.包月復制
付費后30天內不限量復制
特價:9.99元 10元
微信掃碼支付
已付款請點這里 聯系客服
主站蜘蛛池模板: 69国产成人综合久久精 | 深爱婷婷网| 日韩中文字幕一在线 | 一级在线观看视频 | 人人爽天天爽夜夜爽曰 | 两个人看的www免费 两个人看的www高清免费视频 | 成人午夜免费剧场 | 国产成人精品s8p视频 | 91国内精品久久久久免费影院 | 免费在线视频成人 | 99久久亚洲| 午夜免费看片 | 欧美激情综合色综合啪啪五月 | aa一级黄色片 | 天天干天天爽 | 国产日韩亚洲不卡高清在线观看 | 中国二级毛片 | 欧美一级高清片欧美国产欧美 | 亚洲成a人片在线网站 | www.福利视频 | 日韩精品亚洲专区在线观看 | 亚洲 欧美 中文 日韩专区 | 又黄又爽视频好爽视频 | 国产成人刺激视频在线观看 | 在线a人片免费观看国产 | jizzjizz老师太多水日本 | a黄色一级片 | 亚洲日本va中文字幕在线不卡 | 成年人视频免费 | 亚洲精品美女视频 | 日韩国产欧美成人一区二区影院 | 韩国一级黄色大片 | 毛片免费观看网址 | 日韩黄色毛片 | 一级一片| 国产成人精品一区二三区在线观看 | 欧美日本高清一本二本三本 | 噜噜色综合 | 一级毛片免费视频网站 | 久久yy| 成人免费高清 |