电力系统及其自动化--潮流计算--咻咻永远

咻咻永远

各位，我想请教一个问题，我现在做的船舶交流电力系统，想要仿真一个18节点破冰船舶电力系统，主要参考的文献为《船舶电力系统拓扑分析和潮流计算研究-周容华》，需要的是3相潮流计算，但是现在编写的程序不收敛，想请大家帮我看看问题主要出在哪里？源代码在下边。

1. %程序名：qiantuihuitui_I_3.m
   + ]8 ^  n( {# b8 J& N
2. %功能：支路电流前推回推法求解潮流
   
3. clc
   , }, }) K8 @: w; _$ A2 ~: `
4. clear all;
   5 O( P! k9 W8 G: i0 s# p) T
5. %--------------输入网络参数--------------
   
6. %1-支路编号，2-首节点，3-尾节点，4-自阻抗，5-尾节点复功率，6-支路性质（1-馈线段支路，2-变压器支路）,7-尾节点是否带负荷
   
7. DB=[1 1 2 0.000167+j*0.000208 0.42+j*0.31 1 1
   
8.     2 2 3 0.000151+j*0.000188 6.15 1 0
   9 A  Z: `7 ~3 B2 R0 k$ ~
9.     3 2 4 0.000066+j*0.000082 0.38+j*0.29 1 1
   
10.     4 2 5 0.000249+j*0.000310 0 1 0
    $ r0 ~0 {6 o/ S8 S6 E. I
11.     5 2 6 0.000172+j*0.000215 0 1 0
    0 Y# \6 P6 f( K! C! Z" h8 L, V8 w
12.     6 4 7 0.000156+j*0.000195 6.06 1 0
    
13.     7 4 8 0.000162+j*0.000202 6.04 1 0
    ' q4 S) S# s3 d" o  `) a
14.     8 4 9 0.000345+j*0.000430 0 1 0
    
15.     9 4 10 0.000287+j*0.000358 0 1 0
    8 }  N& _5 `$ j6 E, K
16.     10  5 11 0.020563+j*0.321594 0 2 0
    , l, f& c( U" X/ w
17.     11  6 12 0.020563+j*0.321594 0 2 0
    . |0 J: c  L# H( S' n2 Y
18.     12  9 13 0.020563+j*0.321594 0 2 0
    $ L: a; k9 m+ [" I' V
19.     13  10 14 0.020563+j*0.321594 0 2 0  
    
20.     14  11 15 0.000237+j*0.000408 5.72+j*0.12 1 1
    
21.     15  12 16 0.000237+j*0.000408 5.76+j*0.09 1 1
    
22.     16  13 17 0.000292+j*0.000502 5.86+j*0.11 1 1
    , X  r/ b/ a' `- i
23.     17  14 18 0.000274+j*0.000470 5.81+j*0.14 1 1];
    
24. [n,m]=size(DB);  
    
25. B=[1 sin(2*pi/3)+j*cos(2*pi/3) sin(4*pi/3)+j*cos(4*pi/3)];
    
26. C=[1.02 1 1.02 1 1 1 1.02 1.02 1 1 1 1 1 1 1 1 1 1];
    4 X+ j6 k" X3 A: X" `% h
27. U(:,1)=B(1,1)*C';
    % e  N4 Q7 @  S$ N4 o/ l- s0 G1 g4 j3 p, J
28. U(:,2)=B(1,2)*C';
    6 E5 _" F9 o6 X
29. U(:,3)=B(1,3)*C';
    
30. %-------------------求解潮流-----------------
    
31. for k=1:15
    
32. %   I(:,k)=((DB(:,5).')*(diag(1./U(:,k))))';
    5 @2 X% _" x- ~. ~" y' a  V6 A
33. for i=n:-1:1
    & |6 M6 M% E# g" y, o
34.   %如果尾节点带恒功率负荷，需计算节点上负荷注入的电流
    
35. if DB(i,7)==1
    
36. c=DB(i,5)/3;
    
37. d=c/U(i+1,(3*k-2));
    
38. IL(i,3*k-2)=conj(d);
    7 x: J+ G+ L# q7 R
39. d=c/U(i+1,3*k-1);
    , C3 }6 G9 O# x9 F8 i
40. IL(i,3*k-1)=conj(d);
    1 M1 u! [0 c6 O2 j. Y8 a% j- y4 R
41. d=c/U(i+1,3*k);
    - n. G' m$ j% t0 y  y1 L; \1 @
42. IL(i,3*k)=conj(d);
    
43. else
    
44. IL(i,3*k-2)=0;
    
45. IL(i,3*k-1)=0;
    
46. IL(i,3*k)=0;
    - j# r" L- g# h4 H$ Y
47. end
    8 M; Q% _) h7 e$ g" F8 z
48. %找出所有与尾节点相连的支路,计算进支电流（末端电流），存IKj
    
49. A=(find(DB(:,2)==(i+1)));
    # t6 G- c5 d2 b4 M7 U
50. if isempty(A)~=1
    : ]3 t& w, D: [  o( I0 y* }+ Q2 Y
51. IKj(i,3*k-2)=IL(i,3*k-2)-sum(IKj(A,3*k-2));
    
52. IKj(i,3*k-1)=IL(i,3*k-1)-sum(IKj(A,3*k-1));
    0 G" ^) \2 |: Z6 S& F) u
53. IKj(i,3*k)=IL(i,3*k)-sum(IKj(A,3*k));
    / C# v! {; r( h! q# j2 @8 {
54. else
    
55. IKj(i,3*k-2)=IL(i,3*k-2);
    ) \. f- l- ~9 n' B% h
56. IKj(i,3*k-1)=IL(i,3*k-1);
    
57. IKj(i,3*k)=IL(i,3*k);
    % ~5 N  M2 f+ _: w3 E
58. end
    
59. %计算出支电流（始端电流），存IKi
    . T. O: `; b5 b/ [
60. a=DB(i,2);
    
61. b=DB(i,3);
    9 _, Q1 d: Z% y' P0 _, R. B- j
62. Yi=1/DB(i,4)*eye(3);
    8 l: e3 _6 h- X) G1 q7 Q! m
63. if DB(A,6)==1
    
64. F=0.5*Yi*[(U(a,3*k-2)+U(b,3*k-2)) (U(a,3*k-1)+U(b,3*k-1)) (U(a,3*k)+U(b,3*k))]'+[IKj(i,3*k-2) IKj(i,3*k-1) IKj(i,3*k)]';
    8 H; G4 E& N# _. m
65. IKi(i,3*k-2)=F(1,1);
    & x% W6 h  h$ Z2 r: M% R! Z
66. IKi(i,3*k-1)=F(2,1);
    6 }# c+ w& m4 a3 n! I
67. IKi(i,3*k)=F(3,1);
    
68. else
    # @! j' P. s& r' |  E
69. YT=1/real(DB(i,4))+j*(1/imag(DB(i,4)));
    9 `  h% U& S8 B/ r0 _5 ]1 p! ?
70. YI=YT*eye(3);
    9 J! ~8 U$ \8 \% V) q6 b! z
71. YII=YT*eye(3);
    
72. YIII=-YT*eye(3);
    
73. D=inv(YIII)*(-[IKj(i,3*k-2) IKj(i,3*k-1) IKj(i,3*k)]'-YII*([U(b,3*k-2) U(b,3*k-1) U(b,3*k)]'));
    
74. U(a,3*k+1)=D(1,1);
    6 _9 P1 _( ]( ], P+ k6 X
75. U(a,3*k+2)=D(2,1);
      {8 }+ ]( M  {) A6 A
76. U(a,3*k+3)=D(3,1);
    2 N( K4 a& v  @) r& E8 p, q' T
77. E=YI*[U(a,3*k+1) U(a,3*k+2) U(a,3*k+3)]'-YIII*[U(b,3*k-2) U(b,3*k-1) U(b,3*k)]';
    
78. IKi(i,3*k-2)=E(1,1);
    / Z) D9 H; |! Y7 O5 H
79. IKi(i,3*k-1)=E(2,1);
    , `7 x' h0 G' i/ p: i( T8 H
80. IKi(i,3*k)=E(3,1);
    6 I/ @( u. o, q9 ~. d$ ^6 \+ Y
81. end
    + z5 ~2 K) ~2 u+ e
82. end
    $ N* m, Z; h5 |" Y/ c
83. %前推电压
    ) B( F/ ]" r$ T
84. for j=2:n
    
85. U(1,3*k-2)=1.02;
    ! `# n7 S0 d( O! I
86. U(1,3*k-1)=1.02*(sin(2*pi/3)+j*cos(2*pi/3));
    " C' T1 i# h) x7 y" l
87. U(1,3*k)=1.02*(sin(4*pi/3)+j*cos(4*pi/3));
    : L- J- {: E7 L: |
88. a=DB(j-1,2);
    
89. if DB(j-1,6)==1
    : \  q3 a( ~# a- r6 M
90. Yi=1/DB(j-1,4)*eye(3);
    
91. G1=[IKi(a,3*k-2) IKi(a,3*k-1) IKi(a,3*k)]';
    
92. G=[U(a,3*k-2) U(a,3*k-1) U(a,3*k)]'-DB(j-1,4)*eye(3)*(G1-0.5*Yi*[U(a,3*k-2) U(a,3*k-1) U(a,3*k)]');
    ' c- a) N5 B1 I5 b/ ]: F
93. U(j,3*k+1)=G(1,1);
    3 x& v# |- D2 s" q2 Y' `- ~
94. U(j,3*k+2)=G(2,1);
    6 @8 f6 a) j4 n( d3 R
95. U(j,3*k+3)=G(3,1);
    
96. else
    
97. YT=1/real(DB(i,4))+j*(1/imag(DB(i,4)));
    
98. YI=YT*eye(3);
    ( f1 C/ R: W% z$ y
99. YII=YT*eye(3);
    
100. YIII=-YT*eye(3);
     / n5 i: e' L- E7 N7 P, O
101. H=inv(YIII)*([IKi(a,3*k-2) IKi(a,3*k-1) IKi(a,3*k)]'-YI*[U(a,3*k+1) U(a,3*k+2) U(a,3*k+3)]');
     
102. U(j,3*k+1)=H(1,1);
     
103. U(j,3*k+2)=H(2,1);
     
104. U(j,3*k+3)=H(3,1);
     3 e; W, W6 R( p6 U! {& A
105. end
     
106. end
     
107. dU1=U(:,3*k+3)-U(:,3*k);
     
108. dU2=U(:,3*k+2)-U(:,3*k-1);
     
109. dU3=U(:,3*k+1)-U(:,3*k-2);
     
110. if (max(abs(dU1))<0.001&&max(abs(dU2))<0.001&&max(abs(dU2))<0.001)
     
111. break;
     6 u# k' P0 V+ ?
112. else            
     + S% q) y7 S7 f( J  l
113. k=k+1;
     ) v, B) C& B% D- U" Q2 F( b- a+ S
114. end
     ! U+ E2 I- r- h5 a% j  Y# ~6 W/ o
115. end
     7 o  \! x4 T  R  r. D7 F/ V4 M
116. U
     8 R" g/ @( h4 K3 V
117. IL
     
118. IKj
     9 D( R- v' i- ~+ F5 w% L: C& c
119. IKi
     6 V& }& O3 F' T3 ?
120. k
     7 G+ I8 _) Y! A

复制代码

sun710504

学习学习