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

咻咻永远

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

1. %程序名：qiantuihuitui_I_3.m
   
2. %功能：支路电流前推回推法求解潮流
   
3. clc
     W" c4 C0 n7 n3 s
4. clear all;
   ' Z) X$ h9 H$ z4 ~+ {4 N' y6 a
5. %--------------输入网络参数--------------
   
6. %1-支路编号，2-首节点，3-尾节点，4-自阻抗，5-尾节点复功率，6-支路性质（1-馈线段支路，2-变压器支路）,7-尾节点是否带负荷
   + g5 P. J) D8 k( k7 {
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
   1 O  J) L! z8 N" u& X  x0 J
9.     3 2 4 0.000066+j*0.000082 0.38+j*0.29 1 1
   5 ^3 t3 W3 d2 Y8 D& d; G& ^) {
10.     4 2 5 0.000249+j*0.000310 0 1 0
    # p4 ?8 b3 R$ z7 {% n- o5 S
11.     5 2 6 0.000172+j*0.000215 0 1 0
    ( {: L. [6 ?2 G- v$ B& f
12.     6 4 7 0.000156+j*0.000195 6.06 1 0
    - |0 [4 B# O1 F% v; @. e  p+ |. y
13.     7 4 8 0.000162+j*0.000202 6.04 1 0
    8 S& \& l& Y. ~2 q: N( S2 Y' j
14.     8 4 9 0.000345+j*0.000430 0 1 0
    * o$ G* ~3 T  s) G
15.     9 4 10 0.000287+j*0.000358 0 1 0
    # ^2 r2 F% _0 n9 Q2 T) a
16.     10  5 11 0.020563+j*0.321594 0 2 0
    - [5 w1 m" U' x" ~
17.     11  6 12 0.020563+j*0.321594 0 2 0
    ) s: T& l5 m( N! m. Q7 p/ `* v2 ?
18.     12  9 13 0.020563+j*0.321594 0 2 0
    
19.     13  10 14 0.020563+j*0.321594 0 2 0  
      |! p; j, w* Z4 E3 A9 `
20.     14  11 15 0.000237+j*0.000408 5.72+j*0.12 1 1
    ' c  N( R* X  F  Q% N2 x
21.     15  12 16 0.000237+j*0.000408 5.76+j*0.09 1 1
    2 ?7 `; E9 k2 g  d
22.     16  13 17 0.000292+j*0.000502 5.86+j*0.11 1 1
    
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)];
    7 D8 D8 U' e7 l1 T
26. C=[1.02 1 1.02 1 1 1 1.02 1.02 1 1 1 1 1 1 1 1 1 1];
    5 C7 G' ]) A5 d5 s
27. U(:,1)=B(1,1)*C';
    5 L' R% f) s5 t. |
28. U(:,2)=B(1,2)*C';
    $ Z: P3 G0 x) w& y  L6 |- }
29. U(:,3)=B(1,3)*C';
    
30. %-------------------求解潮流-----------------
    " O* z# w4 P' x
31. for k=1:15
    , J  V" q7 Y& q3 f: K9 {
32. %   I(:,k)=((DB(:,5).')*(diag(1./U(:,k))))';
    
33. for i=n:-1:1
    " [, K2 K/ s* `
34.   %如果尾节点带恒功率负荷，需计算节点上负荷注入的电流
    % V- H" Z0 E3 X# Q4 }& \$ ~
35. if DB(i,7)==1
    
36. c=DB(i,5)/3;
    ! F) ~: M) D. f- s( m. |
37. d=c/U(i+1,(3*k-2));
    1 [4 Y3 }5 q; p7 p* M+ B
38. IL(i,3*k-2)=conj(d);
    1 [  x8 D2 b% q/ @6 g4 N5 }
39. d=c/U(i+1,3*k-1);
    0 j8 d" G' H* c, g  @# T0 }4 ?* O5 P
40. IL(i,3*k-1)=conj(d);
      z. u! N2 c7 A7 s
41. d=c/U(i+1,3*k);
    - j  j$ B% g9 ?3 C
42. IL(i,3*k)=conj(d);
    3 d* U" H  i0 ^3 i: K( z
43. else
      U7 X: O" f" P6 B5 n7 m, U5 F# \) N
44. IL(i,3*k-2)=0;
    ) y$ H; Y# N5 W
45. IL(i,3*k-1)=0;
    ; }: t' R1 B3 ?/ k" D
46. IL(i,3*k)=0;
    ' e$ j3 D2 O" M5 Q! p; e
47. end
    
48. %找出所有与尾节点相连的支路,计算进支电流（末端电流），存IKj
    
49. A=(find(DB(:,2)==(i+1)));
    
50. if isempty(A)~=1
    
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));
    ; o' H' r% a  V8 _+ E
53. IKj(i,3*k)=IL(i,3*k)-sum(IKj(A,3*k));
      Q  f( y( y/ V$ R* u# w
54. else
      Z3 f' V3 G6 ~: v% _
55. IKj(i,3*k-2)=IL(i,3*k-2);
    
56. IKj(i,3*k-1)=IL(i,3*k-1);
    
57. IKj(i,3*k)=IL(i,3*k);
    " S& m0 T' d4 o1 }4 [3 ], v
58. end
    
59. %计算出支电流（始端电流），存IKi
      O3 [5 W1 J& Q# q+ l/ t
60. a=DB(i,2);
    
61. b=DB(i,3);
    
62. Yi=1/DB(i,4)*eye(3);
    
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)]';
    ) b6 b/ {& g4 B! }% Y' [; S
65. IKi(i,3*k-2)=F(1,1);
    
66. IKi(i,3*k-1)=F(2,1);
    
67. IKi(i,3*k)=F(3,1);
    
68. else
    
69. YT=1/real(DB(i,4))+j*(1/imag(DB(i,4)));
    
70. YI=YT*eye(3);
    ; M3 r2 |  R8 H4 O/ D! J
71. YII=YT*eye(3);
    
72. YIII=-YT*eye(3);
    2 A. U! o, x1 j
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)]'));
    # j4 b0 Y. X% j! o( W
74. U(a,3*k+1)=D(1,1);
    
75. U(a,3*k+2)=D(2,1);
    * T6 U* e: ?9 P2 n4 h# B
76. U(a,3*k+3)=D(3,1);
    8 {) p* G& c4 ^$ k% x& f6 k: \, `
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);
    
79. IKi(i,3*k-1)=E(2,1);
    ; m% x4 q! i. f- T  T4 V: T& V% ^5 {
80. IKi(i,3*k)=E(3,1);
    
81. end
    4 c' i- u! m: v; R& @0 ]" e
82. end
    
83. %前推电压
    
84. for j=2:n
    
85. U(1,3*k-2)=1.02;
    
86. U(1,3*k-1)=1.02*(sin(2*pi/3)+j*cos(2*pi/3));
    ( [  ~' W) y" w9 W& [* q
87. U(1,3*k)=1.02*(sin(4*pi/3)+j*cos(4*pi/3));
    
88. a=DB(j-1,2);
    3 y# t$ j  A, h. ~* O
89. if DB(j-1,6)==1
    
90. Yi=1/DB(j-1,4)*eye(3);
    ; o# A' |" ?" i& o
91. G1=[IKi(a,3*k-2) IKi(a,3*k-1) IKi(a,3*k)]';
    5 c# m: |% b& W2 E* J+ O4 s
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)]');
    
93. U(j,3*k+1)=G(1,1);
    
94. U(j,3*k+2)=G(2,1);
    ; n1 ]/ e, l* l3 S; z
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);
    0 O8 E, Y) W, g2 b
99. YII=YT*eye(3);
    
100. YIII=-YT*eye(3);
     
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);
     0 c. v& V4 o$ Y7 [: C' [$ M, e
103. U(j,3*k+2)=H(2,1);
     0 c0 w4 U$ d, A  g
104. U(j,3*k+3)=H(3,1);
       y+ `4 f6 j4 R6 ], q
105. end
     4 |  M+ o6 s5 J
106. end
     4 I( t9 Y5 L" E9 v5 h
107. dU1=U(:,3*k+3)-U(:,3*k);
     
108. dU2=U(:,3*k+2)-U(:,3*k-1);
     % O9 W# H/ Y: K; _! v
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;
     4 n' W/ i" F: P" Y# f, m8 H2 \2 Z
112. else            
     " |+ C3 t9 u5 W# v/ |2 H$ A6 k
113. k=k+1;
     
114. end
     
115. end
     
116. U
     2 D' l# n) f% u7 S
117. IL
     9 {9 {: p. b- f- W
118. IKj
     / j: e2 s+ D$ _/ U! k  l3 H
119. IKi
     
120. k
     

复制代码

sun710504

学习学习