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

咻咻永远

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

1. %程序名：qiantuihuitui_I_3.m
   
2. %功能：支路电流前推回推法求解潮流
   * V" S& ^( H, X, d3 _4 \
3. clc
   
4. clear all;
   , H% P) E: i+ Z
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.     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
    
11.     5 2 6 0.000172+j*0.000215 0 1 0
    $ Y# Q6 o9 s7 ]4 ?9 Y) Q! `
12.     6 4 7 0.000156+j*0.000195 6.06 1 0
    ; L$ V2 s* \. `
13.     7 4 8 0.000162+j*0.000202 6.04 1 0
    - p4 B" S4 h4 l# H, V
14.     8 4 9 0.000345+j*0.000430 0 1 0
    
15.     9 4 10 0.000287+j*0.000358 0 1 0
    , @5 w% Y" W: I9 a
16.     10  5 11 0.020563+j*0.321594 0 2 0
    
17.     11  6 12 0.020563+j*0.321594 0 2 0
    ' S% w' [, }3 ], Q) v
18.     12  9 13 0.020563+j*0.321594 0 2 0
    
19.     13  10 14 0.020563+j*0.321594 0 2 0  
    . ^' d- q0 D1 _$ N4 P3 {& m3 j+ m
20.     14  11 15 0.000237+j*0.000408 5.72+j*0.12 1 1
    " L% a2 C+ _4 P4 Y9 B# @1 w# m) r
21.     15  12 16 0.000237+j*0.000408 5.76+j*0.09 1 1
      z8 z' ]3 N' u. \
22.     16  13 17 0.000292+j*0.000502 5.86+j*0.11 1 1
    1 F/ W, Q6 G2 m. `
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];
    2 }8 e9 }% Q1 x' Z
27. U(:,1)=B(1,1)*C';
    ; X* W; C, C* O; H9 ]$ D7 f/ B2 n
28. U(:,2)=B(1,2)*C';
    
29. U(:,3)=B(1,3)*C';
    , f" n* w5 |4 c& z0 @
30. %-------------------求解潮流-----------------
    
31. for k=1:15
    4 E) W; t9 u2 z6 a1 ?. s
32. %   I(:,k)=((DB(:,5).')*(diag(1./U(:,k))))';
    8 j& Q" u$ ~( Q/ |4 @/ w1 _
33. for i=n:-1:1
    2 B0 ^9 ^) K7 a8 U9 N2 }! l
34.   %如果尾节点带恒功率负荷，需计算节点上负荷注入的电流
      @# O. P& G) a" p; d5 a/ w+ S
35. if DB(i,7)==1
    
36. c=DB(i,5)/3;
    , J2 \( _. v$ i6 g# w) g
37. d=c/U(i+1,(3*k-2));
    
38. IL(i,3*k-2)=conj(d);
    
39. d=c/U(i+1,3*k-1);
    6 O9 q) D6 i# L+ }8 ]  @' M8 f
40. IL(i,3*k-1)=conj(d);
    2 N2 N/ Z5 e6 N7 z$ u9 c3 m
41. d=c/U(i+1,3*k);
    
42. IL(i,3*k)=conj(d);
    1 Y0 E( M% {5 r6 E3 x# \5 h
43. else
    
44. IL(i,3*k-2)=0;
    3 ^9 D1 ?) U! X/ G2 T8 M4 u
45. IL(i,3*k-1)=0;
      I. v5 e4 Y' I; ?  ~: x, k+ O3 c+ p
46. IL(i,3*k)=0;
    + g/ R- t) b: |" k
47. end
    
48. %找出所有与尾节点相连的支路,计算进支电流（末端电流），存IKj
    
49. A=(find(DB(:,2)==(i+1)));
    
50. if isempty(A)~=1
    8 }$ M, b- {. b; g0 k4 b$ q
51. IKj(i,3*k-2)=IL(i,3*k-2)-sum(IKj(A,3*k-2));
    7 T5 k1 |6 ?% c& ^4 {, _
52. IKj(i,3*k-1)=IL(i,3*k-1)-sum(IKj(A,3*k-1));
    6 y% f! ?8 e, ^+ H) b2 c% y) s
53. IKj(i,3*k)=IL(i,3*k)-sum(IKj(A,3*k));
    
54. else
    
55. IKj(i,3*k-2)=IL(i,3*k-2);
    0 Q2 [- ?! q0 N5 h+ t
56. IKj(i,3*k-1)=IL(i,3*k-1);
    ) {0 X2 Z2 p# t9 _
57. IKj(i,3*k)=IL(i,3*k);
    
58. end
    0 P! l) T  P* d: Y( H9 o7 ~3 ?; ]1 s
59. %计算出支电流（始端电流），存IKi
    6 M+ \# s, i3 P% t& y; X/ L
60. a=DB(i,2);
    % a2 H' M4 F' X( \3 j: f
61. b=DB(i,3);
    4 L! @0 P( q' ^: U
62. Yi=1/DB(i,4)*eye(3);
    8 X; {! e2 Y; q$ g, e  A8 R5 g
63. if DB(A,6)==1
    6 L: C. ^8 W& y) }
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)]';
    6 m8 Z. t; g+ S6 b& ]( f
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);
    ! d2 s7 Y  {9 Y0 G, k
68. else
    
69. YT=1/real(DB(i,4))+j*(1/imag(DB(i,4)));
    6 {, f: G% p6 k% y, ?- B; S$ C; r0 O
70. YI=YT*eye(3);
    
71. YII=YT*eye(3);
    4 L: h4 w; o6 n; a/ m
72. YIII=-YT*eye(3);
    3 d5 ~) b- k- A- a
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)]'));
    6 e! ]0 Q! X. k. D# w1 p8 T
74. U(a,3*k+1)=D(1,1);
      Y8 C6 g4 v! i* R+ H/ u
75. U(a,3*k+2)=D(2,1);
    
76. U(a,3*k+3)=D(3,1);
    
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);
    ' b3 z# V6 w- L3 B: l
79. IKi(i,3*k-1)=E(2,1);
    
80. IKi(i,3*k)=E(3,1);
    
81. end
    
82. end
    
83. %前推电压
    
84. for j=2:n
    
85. U(1,3*k-2)=1.02;
    0 d- m6 S$ H7 V. @/ O, X
86. U(1,3*k-1)=1.02*(sin(2*pi/3)+j*cos(2*pi/3));
    ! D* T  w. v4 Q$ M5 Z$ Z
87. U(1,3*k)=1.02*(sin(4*pi/3)+j*cos(4*pi/3));
    . o: ^0 p4 D" I! o
88. a=DB(j-1,2);
      ~* x) q. Q% z; ^
89. if DB(j-1,6)==1
    ! t. ]- z/ U3 E6 M% U" L
90. Yi=1/DB(j-1,4)*eye(3);
      F2 o5 H9 _) }
91. G1=[IKi(a,3*k-2) IKi(a,3*k-1) IKi(a,3*k)]';
    8 Z5 |. w# N$ _: p- _' o7 K) {& 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)]');
    
93. U(j,3*k+1)=G(1,1);
    
94. U(j,3*k+2)=G(2,1);
    
95. U(j,3*k+3)=G(3,1);
    / r( B; d: h! G3 w& o
96. else
    
97. YT=1/real(DB(i,4))+j*(1/imag(DB(i,4)));
    
98. YI=YT*eye(3);
    5 M- F7 N5 z/ {& g* X) J
99. YII=YT*eye(3);
    
100. YIII=-YT*eye(3);
     4 t6 e# H6 X( F, Q
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)]');
     & B5 C# R8 K! k$ y7 ]% |
102. U(j,3*k+1)=H(1,1);
     6 e2 E* A% m$ o) S$ @& w4 k
103. U(j,3*k+2)=H(2,1);
     
104. U(j,3*k+3)=H(3,1);
     
105. end
     
106. end
     ! T5 \# S: b3 S( C  s
107. dU1=U(:,3*k+3)-U(:,3*k);
     ' Q- a. \4 r7 Q
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)
     - ?* @: v( J6 F+ m6 A! W5 [4 c
111. break;
     7 @7 d% x8 m) r
112. else            
     8 x6 j2 y; Z, e% j3 i4 H- X
113. k=k+1;
     
114. end
     
115. end
     
116. U
     
117. IL
     1 Y7 V  S, k& J( A
118. IKj
     
119. IKi
     & m4 V  i' |8 _$ m
120. k
     

复制代码

redplum

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

咻咻永远

回复 2# redplum


    不应该是这样的，这是一个迭代计算，结果不收敛，所以还是哪个环节有问题，但是我还没有找到。

redplum

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

咻咻永远

回复 4# redplum
/ ^3 |: x. M% T0 T8 l就是即使我弄成迭代次数为2，结果也不啊~~

魏丹

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