|
软件程序
软件/程序名称: |
牛顿拉夫讯 |
软件/程序大小: |
1k |
软件/程序语言: |
简体中文 |
运行平台: |
Windows XP/2003/Vista |
功能简介: |
牛顿拉夫讯 |
授权方式: |
免费版 |
研究/处理: |
汉化 |
马上加入,结交更多好友,共享更多资料,让你轻松玩转电力研学社区!
您需要 登录 才可以下载或查看,没有账号?立即加入
×
谁能帮忙看看这个牛顿拉夫讯怎么改进的,能用于配电网潮流计算?
& d4 Q# a! T7 r3 h牛顿拉夫讯
; o( r* `" E. R3 [, A0 I$ Adouble PowerFlow(int iNode, int iBalc, ROAD *pRoad, NodePQ *pNodePQ, bool vOutputVoltage = false)
) F+ u* e* ?* n{
: h4 Z k. S. I$ N9 G int M=2*(iNode-1); `+ o0 B$ @# @$ r' g t
int flag1[N]; // PV节点标志
) I6 R( i% t/ ^
, c1 N# O2 I- a' f# Q% o( H5 V- l; M- u double g[N][N]; // 电导; D. {) o$ e" P# o9 ~! ]- a
double b[N][N]; // 电纳9 E3 {5 F2 T# s- A2 L0 ?" T6 u
double dvv[N]; // 电压微分; Q7 G+ i6 Z* T2 @6 [+ |. ^
double angle[N]; // 相角
" a2 s, ]9 [" N& b1 D double p[N],dp[N]; // 有功及其微分6 p( K- K6 m( M/ I5 e0 K
double q[N],dq[N]; // 无功及其微分
* |" E/ M3 |5 ?& x0 Y double e[N]; // 节点电压实部/ Y# N+ d f5 t5 l* V/ P* I7 i9 { c
double f[N]; // 节点电压虚部) @' y3 E. y1 d4 y& B
double *jcbb = new double[M];
/ D) L h1 @7 O; r double *jcb = new double [M][M];
8 {( l j+ I K7 |+ N* k [) q% {; ~, L) | J
int i,j,i1,j1,n,m,t,k1,balc,node/*,flag,npv*/;$ C( [9 P/ f" N
double rr,xx,bb,k,com1,com2,sump=0.0,sumq=0.0;
$ A H* q7 P* o# l1 T' T+ c( R double lossp=0.0, lossq=0.0;
4 ^' T o7 u: o3 V) _7 G node = iNode; // 节点数) m9 d0 G$ I( x9 X! w7 P7 V9 Z1 y
balc = iBalc - 1; // 平衡节点3 m% L# P& ?% v
for(i=0;i<node;i++){
6 b. M6 E0 V: V9 E flag1[i]=0;
- q! D+ V( }6 K" i dp[i]=0.0;& t3 E& D% N: x
dq[i]=0.0;- I7 x- Z" H( m+ u/ H
dvv[i]=0.0;( q; t/ H5 A. K- e' S2 r6 d( b3 F
}
+ W8 [9 t: F) f0 i1 v5 O# B 6 y6 d1 Y8 H# Q! x8 B4 A
for(i=0;i<node;i++) {2 ] |$ R! w/ e; U; ?
for(j=0;j<node;j++){# e% i2 l4 Q+ p
g[i][j]=0.0;
% a9 w1 ^) u8 \) G; n/ ?+ r b[i][j]=0.0; ' ~2 u; Y: @5 ]( }& R4 o
}
* T# z$ w/ c0 t; o* O! ]; d3 Q }
8 c8 A3 U* y% t1 I for(n = 0; n < node-1; n++){ 3 R2 b7 P5 J' |+ ~
i = pRoad[n].head;
- Z8 g# m- e4 p) }7 Y4 U. q6 L j = pRoad[n].tail;
4 y) S, e8 t2 Q! n/ I7 ` rr = pRoad[n].r;
; V" Z+ m4 a. D# R; j5 W xx = pRoad[n].x;
5 F8 [7 v" P5 ^; T2 @' Z$ T// printf("%d\t%d\t%.3f\t%.3f\n", pRoad[n].head, pRoad[n].tail, pRoad[n].r, pRoad[n].x);
& Y' _+ o+ `" J bb = 0;2 `* M- Z# Y6 Y/ i6 W' n
k = 1;) t9 o9 _0 S# n, m4 T0 i# w/ Q
9 y% Z/ g6 Y6 ~* ]2 c- X g[i-1][j-1]+=-(rr/(rr*rr+xx*xx)/k);. _8 v- \5 ]+ p% N9 _3 m
b[i-1][j-1]+=-(-xx/(rr*rr+xx*xx)/k);/ x3 q6 n. a$ t4 G8 W
g[j-1][i-1]+=-(rr/(rr*rr+xx*xx)/k);+ P2 c( C* H7 V
b[j-1][i-1]+=-(-xx/(rr*rr+xx*xx)/k);! d/ S4 T2 e+ Z
g[i-1][i-1]+=rr/((rr*rr+xx*xx)*k*k);
+ k. o5 e3 N. {! F8 d: \7 [0 K2 Y7 _ g[j-1][j-1]+=rr/(rr*rr+xx*xx);: Y9 j' P; V: b0 y; z0 q. w7 E) _. F
b[i-1][i-1]+=bb+(-xx/((rr*rr+xx*xx)*k*k));1 J6 t4 K0 d2 x; ]1 p
b[j-1][j-1]+=bb+(-xx/(rr*rr+xx*xx));# l4 A! M/ {1 F5 f" J! ^
}
7 r. g6 E6 r) k% `# F$ E for(i=0;i<node;i++){
7 F4 ~# a$ C! H$ @- g0 Y8 b angle[i] = 0;- U. g9 \/ B+ n5 ]7 l. g
v[i] = 1;
$ J7 w; m$ h5 T; H9 `2 b p[i] = pNodePQ[i].p;
$ z4 v: S8 F i q[i] = pNodePQ[i].q;
5 x8 X3 j7 L3 {& y4 ]// printf("P:%.3f\t Q:%.3f\n", pNodePQ[i].p, pNodePQ[i].q);4 X0 X& \. ~6 E# a, V
e[i]=v[i]*cos(angle[i]*PI/180);
* H9 D. @) k" C4 F f[i]=v[i]*sin(angle[i]*PI/180);
3 {5 _, M4 v4 ^% {) h4 T }; x" d! K3 `8 k
// flag1[1]=1;4 `6 E0 S) H8 ~/ K! P0 m
. w3 |$ d+ X J L// printf("共有n=%d条支路\n",n);+ v0 R: V( q: y4 G0 g8 w# D8 m
// printf("节点导纳矩阵如下:\n");; x% r; Q6 y/ R8 f
( t8 |# O8 _: K' M// for(i=0;i<node;i++){
7 _: x$ k% R/ [$ ^// for(j=0;j<node;j++)
8 ~; @3 ^. |; z) e, z6 H0 w// printf("%lf+j%lf ", g[i][j],b[i][j]);
9 U3 C# M; _" n; L9 Y+ f1 S// printf("\n\n");}
/ }( B$ T8 q9 e//
! N+ `6 H& W# F' f* v6 g// printf("初始值如下:\n");
& C( T) Y. X% s// printf("节点号 电压 相角 有功 无功 e f\n");
% q) B5 h) z8 b3 M// for(i=0;i<node;i++)) o' }0 Q$ s s# O! J/ _
// printf("%d %lf %lf %lf %lf %lf %lf\n", (i+1), v[i], angle[i], p[i], q[i], e[i], f[i]);
3 d) c+ a! A* v6 p& x3 j
$ H! x4 [/ E% E! E# T for(t=0;t<50;t++){
2 {4 g4 |7 U7 j8 o# d// printf("\n第%d次迭代开始:\n", t);& t0 T: ?$ K* i5 f
for(i=0;i<node;i++){
7 \; T2 }8 N; y: @. C8 a5 L- O if(i==balc) /*如果是平衡节点*/
1 n! H& b* g* M6 E" w) Z, ]% m continue;' L! x3 k5 l0 \( R' q; I) ?( Q
else{ 9 c- W. B- a. J3 ?6 k. Z! `# J- p
com1=0.0, com2=0.0;
+ B3 u% B* h0 o& E6 D/ ^3 p% e for(j=0;j<node;j++){& [2 I; d( L" ~1 J& p% _! ]* @% X& t
com1+=g[i][j]*e[j]-b[i][j]*f[j];. [3 ?) _; F% b% j8 S# T5 z; z" d
com2+=g[i][j]*f[j]+b[i][j]*e[j];0 M8 X+ p% ]3 p& g
}' I( X( h: @8 z8 S
if(flag1[i]!=1){ /*如果是PQ节点*/* O5 \3 b$ }2 B9 B; [# `; b
dp[i]=p[i]-e[i]*com1-f[i]*com2;
0 a" `3 T9 K2 X5 V( D0 M- v( H dq[i]=q[i]-f[i]*com1+e[i]*com2;
2 V* \& q1 F- N( [# J }
" G! A2 }1 `/ d0 A; [0 N else{ V6 [6 V' U: l+ E# n4 V
dp[i]=p[i]-e[i]*com1-f[i]*com2;
( {, ~2 a- z& A. y! y dvv[i]=v[i]*v[i]-(e[i]*e[i]+f[i]*f[i]); k8 j4 Z4 P8 }. e c
}6 X0 B% G0 s- F" C) }0 P0 r* X
}
0 C; q, l0 T! ?" \8 g }
* i$ C" Y# E8 q: K1 ~2 @ ; P( ~* M/ B, r5 }/ G1 D# ?8 C/ o
for(i=0,j=0;i<node;i++){
9 t/ k% T) [0 u# n+ E; m0 Z! ?$ ]& K. O4 G if(i==balc) continue;
3 }& F' y+ r H$ C7 g else if(flag1[i]!=1){' h7 Z6 C# K K; C& J% ?) i
jcbb[j++]=dp[i];9 G/ M* V% ]! @* E
jcbb[j++]=dq[i];, H6 m3 `" ~4 w3 o
// printf("\ndp[%d]=%lf\ndq[%d]=%lf",(i+1),dp[i],(i+1),dq[i]);
- M( A6 B; @- [8 p2 ~: c2 c R4 l }3 v5 H' C) {" N0 z3 A
else{
, f* W" D9 Z4 S; R( T) t- y; N5 v jcbb[j++]=dp[i];
' J# n: k; s. H7 e m4 y# w jcbb[j++]=dvv[i];
; ?8 L! F0 d8 y6 u, f8 `// printf("\ndp[%d]=%lf\ndvv[%d]=%lf",(i+1),dp[i],(i+1),dvv[i]);4 Y$ x4 m) U% z, Q% a) x
}9 { K( p# j( s( O
}2 C2 w% S+ O& R+ h! o2 f
// printf("\n");
6 A# o- P3 Q& a1 z3 n 8 o9 N( ~) p7 T
for(i=0;i<node;i++) {: M5 j4 n2 a) I* s5 t9 x
if(fabs(dp[i])>=E||fabs(dq[i])>=E||fabs(dvv[i])>=E)
6 ^- `* X6 ]. @- O8 C& C2 {% v! v break;* k K; c+ ~/ z# n' ^
}
5 i! d+ p9 ~7 g7 q0 ]$ Q if(i>=node)% [0 O3 X. w/ w( z" ]9 K g
break;
8 N7 B( x! W" E ?/ d, b5 F/ }
' U$ P6 v% C/ I. n
! H: y; Q) w6 p M9 \7 |; `% c( \ for(i=0,i1=0;i<node;i++){
8 s, B2 L, F/ [+ \ if(i==balc)4 v' U+ ]: ~& M1 I+ A. H" j6 [0 H, \
continue;
9 y: A0 I T! x7 {* [ i1++;% @. g ]8 \& H- ]1 n% ^
n=2*i1;
) h+ z4 z# ~* r9 A2 S% |: L* N0 o for(j=0,j1=0;j<node;j++){, K) h/ Z; j/ J" I* ]
if(j==balc) continue;
5 N6 n' k {+ s5 O8 }) o j1++;$ R: L4 p( l$ @" r% G( j
m=j1*2;
' j0 o" v$ ]1 X. l- e" P* H if(i!=j){ 9 C& G% q' X# M( s0 l% u/ E
jcb[n-2][m-2]=-(g[i][j]*e[i]+b[i][j]*f[i]);: G' {8 M* n5 V; @# `, z
jcb[n-2][m-1]=b[i][j]*e[i]-g[i][j]*f[i];
9 j C3 w: Y! f6 r- F if(flag1[i]!=1)
& Z1 Z) \, I9 }# q) G1 K* ]- g {
, L# `8 d# k% C" p/ v9 T& h jcb[n-1][m-2]=jcb[n-2][m-1];5 L" c* y+ K; G4 ^
jcb[n-1][m-1]=-jcb[n-2][m-2];6 q' U% x, Y' y- k
}
$ P* K, e1 p+ i else
( r& E# P' Q! C jcb[n-1][m-2]=jcb[n-1][m-1]=0.0;
( ^) }7 o; V- Y c+ D% F9 q0 ~' E4 q' B& m( ] }
7 i f) r5 T% F/ U; f5 V* J) z, K else{# b4 L6 y k) a% t; a: j# f
com1=0.0, com2=0.0;
6 w2 z+ g, T9 Q7 t for(k1=0;k1<node;k1++)+ ^+ p& T0 y* D; T
{5 i/ A/ W i1 Y+ ]+ u2 ~3 O* e8 I
com1+=g[i][k1]*e[k1]-b[i][k1]*f[k1];
4 G) N/ s/ l. J" o* O com2+=g[i][k1]*f[k1]+b[i][k1]*e[k1];0 h3 H, z& [7 l' ]1 T$ ~7 w" ^4 c
}
, z8 O, v; o% r! w/ i/ { d: \; u jcb[n-2][m-2]=-com1-g[i][i]*e[i]-b[i][i]*f[i];
/ e Q U3 K+ k: C jcb[n-2][m-1]=-com2+b[i][i]*e[i]-g[i][i]*f[i];% w" H7 W1 ^# T C
if(flag1[i]!=1)
2 f+ `! \$ ^, {/ W {& i9 K% w; n1 w+ G* x# a
jcb[n-1][m-2]=com2+b[i][i]*e[i]-g[i][i]*f[i];4 ?2 d; S% l* ]
jcb[n-1][m-1]=-com1+g[i][i]*e[i]+b[i][i]*f[i];9 z: ^; P- |. E: M0 }6 Y' k
}4 P8 S7 g; E1 \% k% W
else j9 K0 n4 M" [) m2 I
{
# m6 k: b5 @7 u5 G) T4 D' w0 H jcb[n-1][m-2]=-2*e[i];
6 h- }) N2 r: D3 \5 ? jcb[n-1][m-1]=-2*f[i];
6 M! x. \# s( A) y% s2 s }! e# U- [$ V. v
}
( s" R9 A6 y& G- f) i1 c }
# A1 r; c( U/ G- Y3 @& O5 W0 _ }
% j# B E G/ s & S5 E4 e$ N) B
// printf("\n第%d次输出-J矩阵如下:\n", t);
& i, } {) s# o2 ? for(i=0;i<2*(node-1);i++){
1 g! @% M8 B" q! Q. a! J2 z5 l for(j=0;j<2*(node-1);j++){
- k( H3 a9 M' H4 r; p3 \' y jcb[i][j]*=-1.0;
; J8 I( t0 |, ?. F7 n* J// printf("%lf ", jcb[i*M+j]);. ^5 [! @5 W$ [3 P( p
}5 P: a q1 ]3 I1 g* i- R
// printf("\n\n");
+ {/ B- j2 q+ [; P; H }1 F$ u6 J$ C3 p/ e0 h, G8 d% Z- Y
& {0 l) N& y. F. x* U# Y
/*算线形方程组*/" z8 A" E: v: d t' z* A |
if (agaus(&jcb[0],&jcbb[0],2*(node-1))!=0) {' k A, n; `/ c3 ^
for (i=0,j=0;i<2*(node-1);i++) {
5 r1 h* ^; i/ d. I/ Y) B- |2 O if(j==balc) , f9 H' T3 K1 t% Y+ |' R' R
j++;
* x1 K. M; B. D. S
6 K6 |$ E b' o if(i%2==0){
+ m$ g/ J+ N6 {. P+ G/ h e[j]+=jcbb[i];
- p( v) [' A( U4 f' f9 N- ^// printf("de[%d]=%lf e[%d]=%lf\n", (j+1),jcbb[i],(j+1),e[j]);# s$ k8 U, e0 Y& z/ E. x) ~( \
}1 _$ l7 i; j8 b {
else{
* i5 ]( ?5 X4 J z* r4 V4 s \ f[j]+=jcbb[i];1 P6 B& f9 N6 d+ V1 A, z
// printf("df[%d]=%lf f[%d]=%lf\n",(j+1),jcbb[i],(j+1),f[j]);. r9 r0 k {! U4 K" t! o2 M5 D
j++;
( ^7 i* j2 ]4 y% a/ r }$ p2 K9 b. r- V/ c' a: }$ y
}* j* ~( n( [/ z
}& V- A4 H; v2 Z1 M
}
6 G9 @% Z; c% O7 L% v
$ p+ W5 T1 C! { I2 r8 j// printf("共经过了%d次迭代\n", t);
/ A3 v( k/ D6 P- x8 u
8 o/ B n" F4 V7 u. j4 e /*求平衡节点P,Q*/' `# ~! w, l9 R
i=balc;
9 o+ ?( l! h: C! |& j3 n com1=0.0, com2=0.0;; F' D; S7 w" {; t# _" }6 Q
for(j=0;j<node;j++){7 T6 {5 F7 g. K; R8 T
com1+=g[i][j]*e[j]-b[i][j]*f[j];: }8 `- b$ E* e' z$ {
com2+=g[i][j]*f[j]+b[i][j]*e[j];1 s- `: b8 c% G0 l8 ]% k4 G0 n
}
# n+ r2 H3 d) X' R9 v9 z6 ^ p[i]=e[i]*com1+f[i]*com2;
, P1 g2 `6 f2 l+ E. [ q[i]=f[i]*com1-e[i]*com2;
& g) b+ j8 D' g R% Q
) i7 w) |3 r/ u' h4 o /*求其他各节点P,Q*/
, {( k! ?" J$ H# @* \5 G2 q: a3 A for(i=0;i<node;i++){
, M& r$ u8 E& H0 q5 D if(i==balc) continue;, r- S1 V) u& |- g! H r/ C4 ?; }* N
else' F& i( ~" ^4 v8 P
{
' ^0 d* [7 V+ c' P com1=0.0, com2=0.0;; }' h0 J5 a4 w- K2 w; `3 ~
for(j=0;j<node;j++)7 V7 K% v9 r" T
{
+ f& l- L' x- q+ e com1+=g[i][j]*e[j]-b[i][j]*f[j];
; F3 k& l8 K+ h9 e+ W' R com2+=g[i][j]*f[j]+b[i][j]*e[j];
) Q0 d4 }( Q' I! h2 f }! n3 D6 W' M7 W# o2 u* d% _6 i/ b
p[i]=e[i]*com1+f[i]*com2;
; S7 Z4 r1 A+ N6 k2 \. Q. b- ^ q[i]=f[i]*com1-e[i]*com2;
. ?! k) v0 S3 S0 X6 Q }
8 i/ L0 p) N. d% F8 ]# ^2 Z1 C }3 i7 k9 G* r, u* z. G8 ~5 ]
: A; \" u8 Z9 O* ^8 a( u
// printf("最终结果为:\n");
6 E8 f" B+ _6 l2 M& m6 X. l# u. K// printf("节点号 电压 相角 有功 无功 e f\n");$ B# ^7 s9 \& b+ J: J, K5 K
for(i=0;i<node;i++)
$ x# f$ f% r, h {
* X" Q( C4 W0 }1 }9 M4 d v[i]=sqrt(e[i]*e[i]+f[i]*f[i]);' x7 i) n/ ?; ?, l' ]1 p
angle[i]=atan(f[i]/e[i])*180/PI;" g) W5 {3 {& k
sump+=p[i];- p) ~' B+ N$ w
sumq+=q[i];( _' v) q1 j @) R2 @6 F
// printf("%d %lf %lf %lf %lf %lf %lf\n", (i+1), v[i], angle[i], p[i], q[i], e[i], f[i]); w1 [0 X' ~1 {
} |
|