|
|
楼主 |
发表于 2015-6-16 15:02:40
|
显示全部楼层
回复 3# 玉门关山
: j' e8 x; T9 `% H( W" m* b
$ h* w U' ^" T9 E5 }6 p S1 x7 l J+ h% C H8 g( ~4 f
因为我不知道问题出在哪里了。。2 B) v, L0 ?7 G+ j) `
do
6 x' c" Q$ t" o8 A; d {
- K8 W2 u B: R //求解不平衡量
$ k! ~% y- {! R, h1 a/ q for(i=0;i<nB;i++)8 E- @* f* r8 n. r: ]2 j" x/ M
{
% V9 h- F! q. u, \" d if(sB.Type!=2)//假如不是平衡节点1 v2 B# |; i6 a. }5 Q* Z" i, G+ ~3 b
{ $ |* }1 [! e8 z% @: }8 R
DP=sB.GenP-sB.LoadP;: V% G: I# J. Z) u
DQ=sB.GenQ-sB.LoadQ;
( O! o! \4 l) |- b& A( N9 |
+ O- m0 l* I( j4 p- v for(j=0;j<nB;j++)* g. n/ ]3 t, Y4 \" e
{
. t8 F$ Q( r& E1 K) k A=sB.Phase-sB[j].Phase;
1 _! ?& L1 T% b5 ]) ^ DP-=sB.Volt*(sB[j].Volt*(g[j]*cos(A)+b[j]*sin(A)));
+ J3 R5 J) {6 l* I 1 y0 `4 _% O% j+ E/ X$ E f8 {: i
if(sB.Type==0)//PQ节点+ m5 B r' f i1 K2 k1 [
DQ-=sB.Volt*(sB[j].Volt*(g[j]*sin(A)-b[j]*cos(A)));
" v! s$ M* C4 N3 r: ~# {
- |' q6 b0 x) B8 x9 {7 F else if(sB.Type==1)//PV节点! f, n9 R: I# Y3 F( c, X+ T
DQ=0;
# b+ r% n y4 l( J# M# T! D& g) h }/ ^: e. Q$ P( y, ^2 h' l" v! e
}
* Q# W- w! [' F' O else if(sB.Type==2)//平衡节点- P- m1 s, N, e% h: [+ `4 g8 A3 F* }" ]- |
DP=DQ=0;* W. c+ M2 p2 ]( X
}
5 {/ m: L# U, k) ~7 h9 O8 N //for(i=0;i<nB;i++). g8 {- g- n$ I0 B* a( X: Z- F
// printf("DP[%d]===%f,DQ[%d]===%f\n",i,DP,i,DQ);" M7 l* g$ h) _# Y, B
# g4 p* r- u) z: z. `
//求解修正方程
# f5 p6 D( W' r& E; Q for(i=0;i<nB-1;i++)" T. K+ A7 n" V* s4 \
AA1=DP[i+1]/sB[i+1].Volt;
* x, J: P% [6 v for(i=0;i<nB-1-count_PVnode;i++)( @( w5 m" r( g8 a% K& c
AA2=DQ[i+1+count_PVnode]/sB[i+1+count_PVnode].Volt;$ n/ p/ @+ x9 Y! q4 ]3 o
calculate_gaosi((double **)b1,BB1,AA1,NBUS-1);//AA是不平衡量,BB是解向量( a7 N6 D6 U) }$ r+ U
calculate_gaosi((double **)b2,BB2,AA2,NBUS-1);
5 }& Z: Z- b9 m4 H) ]6 ~
; s) B- C1 ^- \" j" `7 K. T max1=fabs(AA1[0]);( ^/ L/ |/ K, R* L2 E8 i
for(i=1;i<nB-1;i++)
1 c7 J9 A4 r, y* ~ if(max1<fabs(AA1)) % z* T2 D. M: l# M+ J. ]
max1=fabs(AA1);
$ K9 t/ X$ G" e max2=fabs(AA2[0]);; q5 B5 o* q! C8 h! N. z
for(i=1;i<nB-1-count_PVnode;i++)6 ?5 Z; |7 J- b; g
if(max2<fabs(AA2))
& @0 l5 v# f# h7 x- u/ g0 D8 B* J# f) x max2=fabs(AA2);
3 u8 p% y+ ^' Y/ e6 ~* ~& Z8 a for(i=0;i<nB-1;i++)
8 Q8 T/ d% D9 w: }5 S9 L/ P0 T1 p sB[i+1].Phase+=BB1/sB.Volt;
7 o! A; {9 A6 u1 u H5 e% T for(i=0;i<nB-1-count_PVnode;i++)
9 z# w# K+ O; h! A2 K0 T sB[i+1+count_PVnode].Volt+=BB2;( o' C- H3 w7 `+ |6 }
for(i=0;i<nB;i++)' p& l) L8 L! x( N: M7 V0 V9 L( O0 X
{
; w1 ?" z# B% C9 ^ @ printf("sB[%d].Volt=%f,sB[%d].Phase=%f\n",i,sB.Volt,i,sB.Phase*180/PI);5 P# I+ P) f3 R- W) {7 o5 Y' k
! J, {, ] Y% M$ t0 z4 w" c }
) E/ y2 X4 Z6 H0 N/ L/ P printf("\n");
9 Y4 e! T. Q9 `/ o1 h: h8 z; ? ci++;& ?- ~# N: ~; R, ?2 d: r2 j
}
' D" j, I3 s p' Q2 e while(fabs(max1)>0.00001&&fabs(max2)>0.00001&&ci<40);8 O1 b; ]& W2 @+ O9 k' x
这是我求潮流的程序,用的PQ分解法,最后得到的结果是只能精确到小数点后第二位,第三位就不对了。 |
|