配电系统三相潮流计算 matlab仿真 结果不正确 求解答

咻咻永远

各位~~我根据下图的一个系统用MATLAB编写了一段程序，是电流型前推回代法的三相配电系统潮流计算，但是现在结果不是我想要的，不收敛。请各位帮我看一下，提提建议。图片见附件了。还有个问题：已知的线路都是三相对称线路的阻抗值，那么线路的阻抗矩阵中的互阻抗应该如何计算，是什么样的形式，应用怎样的公式？这里我不是很清楚。

1. %程序名：qiantuihuitui_I_3.m
   9 i9 C( G- D+ Q' P) j2 O5 w
2. %功能：支路电流前推回推法求解潮流
   
3. clc
   * E' Z# H! c3 Z# L$ n
4. clear all;
   & ~5 o# z" ^& J5 W& M9 \- r
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
   2 ~$ H) N1 o1 `
8.     2 2 3 0.000151+j*0.000188 6.15 1 0
   
9.     3 2 4 0.000066+j*0.000082 0.38+j*0.29 1 1
   8 ?/ a! V4 e( q5 }1 M- i
10.     4 2 5 0.000249+j*0.000310 0 1 0
    & i0 d: {6 T# A: Z+ U
11.     5 2 6 0.000172+j*0.000215 0 1 0
    
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
    2 J6 D3 x8 |$ `" \6 _8 a
14.     8 4 9 0.000345+j*0.000430 0 1 0
    2 ^' f" c9 w9 N5 k6 Q5 D
15.     9 4 10 0.000287+j*0.000358 0 1 0
    8 R7 a. d! g+ R
16.     10  5 11 0.020563+j*0.321594 0 2 0
    ( l4 H9 L9 M& `1 X/ c0 B6 {
17.     11  6 12 0.020563+j*0.321594 0 2 0
    0 `3 u- V+ _2 g3 d
18.     12  9 13 0.020563+j*0.321594 0 2 0
    
19.     13  10 14 0.020563+j*0.321594 0 2 0  
    * X1 Q1 Y- G6 ?0 ?$ ^1 X$ w
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
    + ^) Z: x, v* j1 _; F. I+ _
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);  
    0 A1 z; W6 ~3 Q% P/ g9 z+ c
25. B=[1 sin(2*pi/3)+j*cos(2*pi/3) sin(4*pi/3)+j*cos(4*pi/3)];
    ' b& A% X- e) k4 F- w; t& 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];
    " S: W, G# ^5 y" J+ j  u
27. U(:,1)=B(1,1)*C';
    : V3 M2 l. L; k/ G
28. U(:,2)=B(1,2)*C';
    3 `) i9 Q0 J5 K2 f
29. U(:,3)=B(1,3)*C';
    
30. %-------------------求解潮流-----------------
    " E9 }' y: z, U0 r, J
31. for k=1:15
    
32. %   I(:,k)=((DB(:,5).')*(diag(1./U(:,k))))';
    ! u( u! \& U4 P% P* m
33. for i=n:-1:1
    # e! k& N# G& e7 M8 J- ?
34.   %如果尾节点带恒功率负荷，需计算节点上负荷注入的电流
    
35. if DB(i,7)==1
    
36. c=DB(i,5)/3;
    & y( Y1 ^" P2 L
37. d=c/U(i+1,(3*k-2));
    $ _) N- E% l% m* E9 A' x5 [
38. IL(i,3*k-2)=conj(d);
    
39. d=c/U(i+1,3*k-1);
    7 p" l1 V3 Y) H# t! U1 f+ k
40. IL(i,3*k-1)=conj(d);
    4 ]  x* ?+ L( @
41. d=c/U(i+1,3*k);
    7 u+ R' N0 R8 a0 [' S+ d- P" T% k
42. IL(i,3*k)=conj(d);
    6 E& b7 `) j( H
43. else
    
44. IL(i,3*k-2)=0;
    5 Q/ e! N! ?. y# n8 A
45. IL(i,3*k-1)=0;
    
46. IL(i,3*k)=0;
    % Z: T. c+ G! m5 l: m; b' I  u
47. end
    + |/ {! L- A( ?) z9 X8 w
48. %找出所有与尾节点相连的支路,计算进支电流（末端电流），存IKj
    : h  _3 b; d+ R: I1 m4 E3 M" Q: n
49. A=(find(DB(:,2)==(i+1)));
    : _$ ~$ {, }+ g+ L' f1 @! z
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));
    
53. IKj(i,3*k)=IL(i,3*k)-sum(IKj(A,3*k));
    
54. else
    9 l' B& u1 I' j9 f0 {
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);
    % p& x: N) D' V0 c% x/ ~8 S
58. end
    
59. %计算出支电流（始端电流），存IKi
    
60. a=DB(i,2);
    
61. b=DB(i,3);
    , a+ A% _4 I3 X# {" ^
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)]';
    
65. IKi(i,3*k-2)=F(1,1);
    8 I6 |4 i* o, n+ H
66. IKi(i,3*k-1)=F(2,1);
    / E+ ^9 P! ?; y& ?) C
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);
    
71. YII=YT*eye(3);
    
72. YIII=-YT*eye(3);
    5 f4 I: Z" c; c6 m
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)]'));
    4 j" k1 _9 A) M" Z) z
74. U(a,3*k+1)=D(1,1);
    
75. U(a,3*k+2)=D(2,1);
    
76. U(a,3*k+3)=D(3,1);
    " l8 J8 A9 T( H" \! c$ R
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)]';
    / t7 `4 ^2 d% h
78. IKi(i,3*k-2)=E(1,1);
    
79. IKi(i,3*k-1)=E(2,1);
    
80. IKi(i,3*k)=E(3,1);
    
81. end
    : o% o: Y5 P  [2 h( N! z
82. end
    1 C) `  [" [, |+ W" U
83. %前推电压
    
84. for j=2:n
    
85. U(1,3*k-2)=1.02;
    . q! U! G. j! d$ N* ]! S( w! h
86. U(1,3*k-1)=1.02*(sin(2*pi/3)+j*cos(2*pi/3));
    4 \% g0 C. P  I) e$ |- ^
87. U(1,3*k)=1.02*(sin(4*pi/3)+j*cos(4*pi/3));
    6 P! g$ q' E7 n, O! }- T
88. a=DB(j-1,2);
    
89. if DB(j-1,6)==1
    7 q" A6 e; d1 I# [; V0 a# W$ S- E
90. Yi=1/DB(j-1,4)*eye(3);
    8 x7 N- K' v/ @; U% v
91. G1=[IKi(a,3*k-2) IKi(a,3*k-1) IKi(a,3*k)]';
    ) ~  ^8 ^9 C! n4 l/ R
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);
    # o1 {- w, ^; D/ ~
94. U(j,3*k+2)=G(2,1);
    3 k( e" Y/ P5 |/ Y9 ?
95. U(j,3*k+3)=G(3,1);
    0 J" H% w7 G4 I. s4 Y
96. else
    ! T. R! @. r4 x/ p
97. YT=1/real(DB(i,4))+j*(1/imag(DB(i,4)));
    
98. YI=YT*eye(3);
    
99. YII=YT*eye(3);
    0 l! P) Z' z% Q
100. YIII=-YT*eye(3);
     ' A9 f1 e( F# u9 d
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)]');
     * s5 f6 i) f5 `* O7 ^) z! I4 |0 F2 Y
102. U(j,3*k+1)=H(1,1);
     
103. U(j,3*k+2)=H(2,1);
     + _2 Y8 w1 O/ U# k
104. U(j,3*k+3)=H(3,1);
     # }/ a) R% z8 `2 N1 z& b1 n* |
105. end
     6 @2 {5 U4 h# D: k) d2 U
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;
     ) ^: _5 x$ h/ v8 C  F
112. else            
     2 l* U/ v1 r/ B- \3 h' _
113. k=k+1;
     
114. end
     ; q; J' e9 G8 A4 @
115. end
     
116. U
     
117. IL
     
118. IKj
     
119. IKi
     ) e+ s/ i3 m/ P
120. k
     # w5 a! O" A2 ^8 E( k& I

复制代码

redplum

我运行了楼主的这段代码，发现最后结果矩阵中很多元素都是Nan或Inf，是不是潮流计算阻抗矩阵中多数情况是这样的，还是不应该很多元素是NaN或Inf

咻咻永远

回复 2# redplum
3 _2 A6 _2 K) [0 K! [. a

( n0 v3 b' A4 H" [    不应该是这样的，这是一个迭代计算，结果不收敛，所以还是哪个环节有问题，但是我还没有找到。

redplum

这段代码的原始数据有17行7列，K=1:15计算了15次循环。建议将数据弄少一些，这样才能更好判断到底是哪里不收敛

咻咻永远

回复 4# redplum
; ~4 t5 A- J1 x* m5 ~' R9 F  m* u就是即使我弄成迭代次数为2，结果也不啊~~

魏丹

这个问题解决了吗 求答案~