3节电力系统状态估计(WLS算法)

1、3 節電系統狀態估計報告:打开的开关:打开的刀闸G : 发电机: 母线:线路:负荷:连接线（没有阻抗）3 節點系統主接線圖任務：1、 采用最小二乘狀態估計算法,所有量測權重都取1.0,編寫狀態估計程序（C/Matlab） 。2、 按量測類型，列出量測方程（每一類寫出一個方程）3、 畫出程序流程4、 提交源程序，程序中每個函數作用5、 提交計算輸出結果（屏幕拷貝）系統參數：功率基值：100MW電壓基值：230 kV線路阻抗參數（標麼值）：線路首節點末節點R（標麼值）X（ 標麼值）半充電電納Line1-2120.0170.0920Line1-3130.00850.0720Line2-3230. 0

2、1190.10080線路量測（流出母線為正）：線路首端量測首端真值末端量測末端真值Line1-2-61.3+j1.2-63.32+j1.2460+j2.463.99+j2.35Line1-3-46.7-j14.8-47.69-j14.7945.9+j16.547.9+j16.54Line2-324-j6.624.93-j6.6-24+j7.2-24.86+j7.24母線電壓量測：母線電壓量測真值（幅值/角度）1232233/02234.56235.56/3.233236.46236.46/1.82負荷量測（流出母線為正）：負荷量測真值Load1111+j13.50111.02+j13.55發電

3、量測（流入母線為正）：發電機量測真值unit288-j4.2488.92-j4.24unit323+j2423.04+j23.78注：量測存在誤差【數據預處理】首先根據基值將已知量測值均轉換為標么值，並將功率值轉換為流入量，得到如下數據：線路導納參數（標麼值）：線路首節點末節點G（標麼值）B（標麼值）半充電電納Line1-2121.9422-10.51070Line1-3131.6171-13.6980Line2-3231.1551-9.78430線路注入功率量測（標么值）：線路首端量測首端真值末端量測末端真值Line1-20.613-j0.0120.6332-j0.0124-0.60-j0.

4、024-0.6399-j0.0235Line1-30.467+j0.1480.4769+j0.1479-0.459-j0.165-0.479-j0.1654Line2-3-0.24+j0.066-0.2493+j0.0660.24-j0.0720.2486-j0.0724負荷點注入功率量測（標么值）：負荷量測真值Load1-1.11-j0.1350-1.1102-j0.1355發電機節點注入量測（流入母線為正）：發電機量測真值unit20.88-j0.04240.8892-j0.0424unit30.23+j0.240.2304+j0.2378母線電壓量測（標么值）母線電壓量測真值（幅值/角度

5、）11.00871.0130/021.01981.0242/3.2331.02811.0281/1.82【量測方程】0 度，以 vi 表示誤差值。選擇節點1 電壓相角為參考，為1） 節點 1 電壓量測方程：Vi=Vi+v1即 1.0087=V1+v12） 1-3 支路 1 號節點處注入有功功率功率：2Pij=Vi gij-ViVj(gijcos +bijsin)+v20.613=V12g13-V 1V3(g13cos+b13sin)+v2即0.613=-1.6171V12-V1V3(-1.6171cos+13.698sin)+v23） 1 號節點注入功率：2Pi=Vi Gii+Gijcos +

6、Bijsin +v32P1=V 1 G11+G1jcos +B1jsin +v3即 -1.11=3.5613V12+V1V2(-1.9442cos -10.5107sin )+V1V3(-1.6171 cos -13.698 sin )+v3其中 iterations 為迭代次數，可見本例迭代次數為4， 收斂較快，狀態估計得到節點1、 2、 3 電壓分別為：234.01 44444444444444444444444444444444444444444444【程序說明】1、計算h 矩陣函數cal_hfunction h=cal_h(V,th0,B,G) % 其中， V為節點電壓估計值，th0

7、為節點電壓相角估計% 值， B 為節點電導矩陣，G 為節點電納矩陣b=-B;% 線路電導矩陣g=-G;% 線路電納矩陣P=zeros(3,1);% 初始化，節點注入功率Q=zeros(3,1);PP=zeros(3,3); % 線路注入功率QQ=PP;th=0;th0;% 節點 1 電壓相角為0for i=1:3P_P=0;Q_Q=0; for j=1:3 if (j=i)P_P=P_P+V(i)*V(j)*(G(i,j)*cos(th(i)-th(j)+B(i,j)*sin(th(i)-th(j);Q_Q=Q_Q+V(i)*V(j)*(G(i,j)*sin(th(i)-th(j)-B(i,j

8、)*cos(th(i)-th(j); PP(i,j)=(V(i)2)*g(i,j)-V(i)*V(j)*(g(i,j)*cos(th(i)-th(j)+b(i,j)*sin(th(i)-th(j);QQ(i,j)=-(V(i)2)*b(i,j)-V(i)*V(j)*(g(i,j)*sin(th(i)-th(j) -b(i,j)*cos(th(i)-th(j);endendP(i)=(V(i)2)*G(i,i)+P_P;Q(i)=-(V(i)2)*B(i,i)+Q_Q;endVV=V(1);V(2);V(3);h=P;Q;PP(1,2);PP(2,3);PP(3,1);QQ(1,2);QQ(2,

9、3);QQ(3,1);PP(1,3);PP(2,1);P P(3,2);QQ(1,3);QQ(2,1);QQ(3,2);VV;2、計算H 矩陣函數cal_HHfunction H=cal_HH(V,th0,G,B,P,Q)% 其中，P， Q為根據電壓估計值計算得到節點% 注入電壓b=-B;g=-G;PV=zeros(3,3); % 節點注入功率對電壓幅值偏導數QV=zeros(3,3);Pth=zeros(3,3); % 節點注入功率對電壓相角偏導數Qth=zeros(3,3);PPV=zeros(3,3); %P ij對 Vj偏導數QQV=zeros(3,3); %Q ij對 Vj 偏導數P

10、Pth=zeros(3,3); %P ij對 th j偏導數QQth=zeros(3,3); %Q ij對 th j偏導數PPV1=zeros(3,3); %P ij對 Vi偏導數QQV1=zeros(3,3); %Q ij對 Vi偏導數PPth1=zeros(3,3); %P ij對 th i偏導數QQth1=zeros(3,3); %Q ij對 th i偏導數VV=eye(3);Vth=zeros(3,2);th=0;th0;for i=1:3for j=1:3 if (i=j)PV(i,j)=V(i)*(G(i,j)*cos(th(i)-th(j)+B(i,j)*sin(th(i)-th

11、(j); QV(i,j)=V(i)*(G(i,j)*sin(th(i)-th(j)-B(i,j)*cos(th(i)-th(j); PPV(i,j)=-V(i)*(g(i,j)*cos(th(i)-th(j)+b(i,j)*sin(th(i)-th(j); QQV(i,j)=-V(i)*(g(i,j)*sin(th(i)-th(j)-b(i,j)*cos(th(i)-th(j); PPV1(i,j)=2*V(i)*g(i,j)-V(j)*(g(i,j)*cos(th(i)-th(j)+b(i,j)*sin(th(i)-th(j);QQV1(i,j)=-2*V(i)*b(i,j)-V(j)*(g

12、(i,j)*sin(th(i)-th(j)-b(i,j)*cos(th(i)-th(j);if (j=1)Pth(i,j)=V(i)*V(j)*(G(i,j)*sin(th(i)-th(j)-B(i,j)*cos(th(i)-th(j);Qth(i,j)=-V(i)*V(j)*(G(i,j)*cos(th(i)-th(j)+B(i,j)*sin(th(i)-th(j);PPth(i,j)=-V(i)*V(j)*(g(i,j)*sin(th(i)-th(j)-b(i,j)*cos(th(i)-th(j);QQth(i,j)=-V(i)*V(j)*(-g(i,j)*cos(th(i)-th(j)-

13、b(i,j)*sin(th(i)-th(j);endif(i=1)PPth1(i,j)=V(i)*V(j)*(g(i,j)*sin(th(i)-th(j)-b(i,j)*cos(th(i)-th(j);QQth1(i,j)=-V(i)*V(j)*(g(i,j)*cos(th(i)-th(j)+b(i,j)*sin(th(i)-th(j);endelsePV(i,j)=(G(i,i)*(V(i)2)+P(i)/V(i);QV(i,j)=(Q(i)-(V(i)2)*B(i,i)/V(i); if (j=1)Pth(i,j)=-B(i,i)*(V(i)2)-Q(i);Qth(i,j)=P(i)-(V

14、(i)2)*G(i,i);endendendendH=PV,Pth(:,2:3);QV,Qth(:,2:3);.PPV1(1,2),PPV(1,2),0,PPth(1,2),0;.0,PPV1(2,3),PPV(2,3),PPth1(2,3),PPth(2,3);.PPV(3,1),0,PPV1(3,1),0,PPth1(3,1);.QQV1(1,2),QQV(1,2),0,QQth(1,2),0;.0,QQV1(2,3),QQV(2,3),QQth1(2,3),QQth(2,3);.QQV(3,1),0,QQV1(3,1),0,QQth1(3,1);.PPV1(1,3),0,PPV(1,3)

15、,0,PPth(1,3);.PPV(2,1),PPV1(2,1),0,PPth1(2,1),0;.0,PPV(3,2),PPV1(3,2),PPth(3,2),PPth1(3,2);.QQV1(1,3),0,QQV(1,3),0,QQth(1,3);.QQV(2,1),QQV1(2,1),0,QQth1(2,1),0;.0,QQV(3,2),QQV1(3,2),QQth(3,2),QQth1(3,2);. VV,Vth;3、主程序calculate_all.m 文件format longG=3.5613,-1.9442,-1.6171;.-1.9442,3.0993,-1.1551;.-1.6

16、171,-1.1551,2.7722;B=-24.2087,10.5107,13.698;.10.5107,-20.295,9.7843;.13.698,9.7843,-23.4832;%B 為節點電導矩陣%G 為節點電納矩陣P=-1.11;0.88;0.23;Q=-0.135;-0.0424;0.24;節點注入功率量測值PP=0.613;-0.24;-0.459;QQ=-0.012;0.066;-0.165;% 線路 1-2,2-3,3-1 注入功率在首端量測值% 線路 1-3,2-1,3-2 注入功率在首端量測值PP1=0.467;-0.6;0.24;QQ1=0.148;-0.024;-0

17、.072;% 節點電壓幅值量測值V=1.0087;1.0198;1.0281;R=diag(ones(21,1);Z=P;Q;PP;QQ;PP1;QQ1;V;權重都取為1% 量測值矩陣V0=1;1;1;% 初值th0=0;0;delta=100; iterations=0; % 迭代次數 while delta>0.000001iterations=iterations+1;h=cal_h(V0,th0,B,G); % 計算 h矩陣H=cal_HH(V0,th0,G,B,h(1:3,1),h(4:6,1);% 計算H矩陣A=H'*inv(R)*H; b=H'*inv(R)*(Z-h); d=Ab; % 求解修正值 delta=max(abs(d);V0=V0+d(1:3,1);% 修正估計值th0=th0+d(4:5,1); end iterations V0=V0*230;% 轉換為有名值th0=th0*180/pi;% 轉換為度for i=1:3j=num