02 matlab搜題買粉絲(matlab用牛頓法計算潮流需要在命令窗口輸入什么)

改進牛頓法通過對傳統法進行一定的近似，將J陣寫成UDUT 的形式。U僅由網絡拓撲決定，是一個上三角矩陣；D是一個對角矩陣。在牛拉法中，需要對J陣因子分解與前代回代，改進法則只有前推回代的計算過程。它很好地改善了傳統法以及前推回推法。經過算例計算結果證明，改進法可以避免J陣病態，且擁有前推回代法的收斂速度、精度，又由于它屬于牛頓型算法，所以該算法已經得到了廣泛的運用[18]。

下面附帶電力系統分析牛頓法算例及matlab程序：

網絡結構如下：系統結構圖

系統參數如下：

在上圖所示的簡單電力系統中，系統中節點1、2為PQ節點，節點3為PV節點，節點4為平衡節點，已給定P1s+jQ1s=-0.30-j0.18 P2s+jQ2s=-0.55-j0.13 P3s=0.5 V3s=1.10 V4s=1.05∠0°

容許誤差ε=10-5

節點導納矩陣：

導納矩陣

各節點電壓：

節點 e f v ζ

1.0.984637 -0.008596 0.984675 -0.500172

2.0.958690 -0.108387 0.964798 -6.450306

3.1.092415 0.128955 1.100000 6.732347

4.1.050000 0.000000 1.050000 0.000000

各節點功率：

節點 P Q

1-0.300000 -0.180000

2–0.550000 -0.130000

3 0.500000 -0.551305

4 0.367883 0.264698

matlab程序如下：

// 牛頓法潮流計算matlab程序

clc;

Y=[1.042093-8.242876i -0.588235+2.352941i 3.666667i -0.453858+1.891074i;

-0.588235+2.352941i 1.069005-4.727377i 0 -0.480769+2.403846i;

3.666667i 0 -3.333333i 0;

-0.453858+1.891074i -0.480769+2.403846i 0 0.934627-4.261590i];

%導納矩陣

e=[1 1 1.1 1.05];%初始電壓

f=zeros(4,1);

V=zeros(4,1);%節點電壓

Ws=[-0.3 ; -0.18 ; -0.55 ; -0.13 ; 0.5 ; 1.1];%初始功率

W=zeros(6,1);

n=length(Y);%節點數

J=zeros(2*(n-1));%雅可比矩陣

delta_v=zeros(1,6);

delta_w=Ws;

G=real(Y);

B=imag(Y);

S=zeros(4,2);

c=0;%循環次數

m=input('請輸入PQ節點數：');

while max(abs(delta_w))>10^-5

for i=1:(n-1)%以下為求取雅可比矩陣

for j=1:(n-1)

if (i~=j)

J(2*i-1,2*j-1)=-(G(i,j)*e(i)+B(i,j)*f(i));

J(2*i,2*j)=-J(2*i-1,2*j-1);

J(2*i-1,2*j)=B(i,j)*e(i)-G(i,j)*f(i);

J(2*i,2*j-1)=J(2*i-1,2*j);

end

end

end

for j=1:(n-2)

J(6,2*j-1)=0;

J(6,2*j)=0;

end%以上為非對角線元素

s1=0;

s2=0;

for i=1:(n-1)

for j=1:n

s1=s1+(G(i,j).*e(j)-B(i,j).*f(j));

s2=s2+(G(i,j).*f(j)+B(i,j).*e(j));

end

J(2*i-1,2*i-1)=-s1-G(i,i) *e(i)-B(i,i)*f(i);

J(2*i-1,2*i)=-s2+B(i,i) *e(i)-G(i,i)*f(i);

s1=0;

s2=0;

end

for i=1:m

for j=1:n

s1=s1+G(i,j).*f(j)+B(i,j).*e(j);

s2=s2+(G(i,j).*e(j)-B(i,j).*f(j));

end

J(2*i,2*i-1)=s1+B(i,i) *e(i)-G(i,i)*f(i);

J(2*i,2*i)=-s2+G(i,i) *e(i)+B(i,i)*f(i);

s1=0;

s2=0;

end

J(6,5)=-2*e(3);

J(6,6)=-2*f(3);%對角線元素求解

for i=1:m

for j=1:n

s1=s1+e(i)*(G(i,j).*e(j)-B(i,j).*f(j))+f(i)*(G(i,j).*f(j)+B(i,j).*e(j));

s2=s2+f(i)*(G(i,j).*e(j)-B(i,j).*f(j))-e(i)*(G(i,j).*f(j)+B(i,j).*e(j));

end

delta_w(2*i-1)=Ws(2*i-1)-s1;

delta_w(2*i)=Ws(2*i)-s2;

W(2*i-1)=s1;

W(2*i)=s2;

s1=0;

s2=0;

end

for j=1:n

s1=s1+e(3)*(G(3,j).*e(j)-B(3,j).*f(j))+f(3)*(G(3,j).*f(j)+B(3,j).*e(j));

end

delta_w(5)=Ws(5)-s1;

delta_w(6)=(Ws(6)^2-(e(3)^2+f(3)^2));

W(5)=s1;

W(6)=sqrt(e(3)^2+f(3)^2);%以上求功率差值

delta_v=-inv(J)*delta_w;

for i=1:(n-1)

e(i)=e(i)+delta_v(2*i-1);

f(i)=f(i)+delta_v(2*i);

end%求電壓差值

c=c+1;

end

for x=1:4

V(x)=e(x)+f(x)*1i;

end%節點電壓

s1=0;

for x=3:4

for j=1:4

s1=s1+買粉絲nj(Y(x,j))*買粉絲nj(V(j));

end

S(x,1)=real(V(x)*s1);

S(x,2)=imag(V(x)*s1);

s1=0;

end%PV與平衡節點功率

for x=1:2

S(x,1)=W(2*x-1);

S(x,2)=W(2*x);

end%節點功率

c

J

V

S

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

100

101

102

103

104

105

106

1

2

3

4

5

6

7

8

9

10

11

12

1

很赞哦!（23125）