phases2.m

同步发电机空载三相两项短路情况下的各相波形以及励磁电流波形

x_d = 1.75581; x_ad = 1.56711; x_F = 1.73883; x_D = 1.63555; r_F = 0.00113; r_D = 0.00929; x_FD = 1.61805;
x_q = 1.75581; x_aq = 1.56711; x_H = 2.22911; x_Q = 1.64064; r_H = 0.01862; r_Q = 0.01788; r = 0.004;
DELT_t = 0.0005; t_max = 0.5 ;U = 1;sita_0 = 1.57;
T=0:DELT_t:0.5;
I=zeros(6,length(T));                      %%%I(:,j)表示第j点的值的存储
u_d=zeros(1,length(T));
u_q=zeros(1,length(T));
u_q(1)=1;
omg=1.0;
T_B=T*100*pi;
delt_t=DELT_t*100*pi;
I(:,1)=[0;0;u_q(1)/x_ad;0;0;0];

C=[-(x_d+delt_t*r/2)    delt_t*omg*x_q/2               x_ad                 x_ad               -delt_t*omg*x_aq/2         -delt_t*omg*x_aq/2
    -delt_t*omg*x_d/2   -(x_q+delt_t*r/2)     delt_t*omg*x_ad/2       delt_t*omg*x_ad/2          x_aq                     x_aq
    -x_ad               0                       x_F+delt_t*r_F/2         x_FD                    0                         0
    -x_ad               0                           x_FD               x_D+delt_t*r_D/2          0                         0
    0                  -x_aq                        0                    0                      x_H+delt_t*r_H/2          x_aq
    0                  -x_aq                        0                    0                      x_aq                      x_Q+delt_t*r_Q/2];   %%%常矩阵C
C=inv(C);

j=length(T);
D=[delt_t/2 0 0 0 0 0
    0 delt_t/2 0 0 0 0]';    %D_1=D(:,1);D_2=D(:,2);  A_G11=D(1,1);A_G12=D(1,2);A_G21=D(2,1);A_G22=D(2,2);B_G1=E(1);B_G2=E(2)(E=C\E计算后)
% D_1=D(:,1);D_2=D(:,2);
u_F=r_F*I(3,1);

for j=2:length(T_B)
    %B_G1=E(1);B_G2=E(2)(E=C\E计算后)
    e1 = -(x_d-delt_t*r/2)*I(1,j-1) - delt_t*omg*x_q/2*I(2,j-1) + x_ad*I(3,j-1) + x_ad*I(4,j-1) + delt_t*omg*x_aq/2*I(5,j-1) + delt_t*omg*x_aq/2*I(6,j-1) + delt_t*u_d(j-1)/2;
    e2 = delt_t*omg*x_d/2*I(1,j-1) - (x_q-delt_t*r/2)*I(2,j-1) - delt_t*omg*x_ad/2*I(3,j-1) - delt_t*omg*x_ad/2*I(4,j-1) + x_aq*I(5,j-1) + x_aq*I(6,j-1) + delt_t*u_q(j-1)/2;
    e3 = -x_ad*I(1,j-1) + (x_F-delt_t*r_F/2)*I(3,j-1) + x_FD*I(4,j-1) + delt_t*u_F;
    e4 = -x_ad*I(1,j-1) + x_FD*I(3,j-1) + (x_D-delt_t*r_D/2)*I(4,j-1);
    e5 = -x_aq*I(2,j-1) + (x_H-delt_t*r_H/2)*I(5,j-1) +x_aq*I(6,j-1);
    e6 = -x_aq*I(2,j-1) + x_aq*I(5,j-1) + (x_Q-delt_t*r_Q/2)*I(6,j-1);
    %     function [E]=EYZFC(XD,XAD,RXF,RXD,RF,RD,RXFD,XQ,XAQ,RXH,RXQ,RH,RQ,RR,DLTT2,OMG,X,UD,UQ,UF)
    %     E=zeros(6,1);
    %     E(1)=-(XD-DLTT2*RR)*X(1)-DLTT2*OMG*XQ*X(2)+XAD*(X(3)+X(4))+DLTT2*OMG*XAQ*(X(5)+X(6))+DLTT2*UD;
    %     E(2)=DLTT2*OMG*XD*X(1)-(XQ-DLTT2*RR)*X(2)-DLTT2*OMG*XAD*(X(3)+X(4))+XAQ*(X(5)+X(6))+DLTT2*UQ;
    %     E(3)=-XAD*X(1)+(RXF-DLTT2*RF)*X(3)+RXFD*X(4)+DLTT2*2.0*UF;
    %     E(4)=-XAD*X(1)+RXFD*X(3)+(RXD-DLTT2*RD)*X(4);
    %     E(5)=-XAQ*X(2)+(RXH-DLTT2*RH)*X(5)+XAQ*X(6);
    %     E(6)=-XAQ*X(2)+XAQ*X(5)+(RXQ-DLTT2*RQ)*X(6);
    %  end
    E=[e1;e2;e3;e4;e5;e6];
    %求u_d(t)&u_q(t)
    E_p=C*E;D_p=C*D;
    A_G11=D_p(1,1);A_G12=D_p(1,2);A_G21=D_p(2,1);A_G22=D_p(2,2);
    B_G1=E_p(1);B_G2=E_p(2);

    G_11=cos(sita_0+T_B(j)-2*pi/3) - cos(sita_0+T_B(j)+2*pi/3);
    G_12=-sin(sita_0+T_B(j)-2*pi/3) + sin(sita_0+T_B(j)+2*pi/3);
    G_21=A_G11*cos(sita_0+T_B(j)) - A_G21*sin(sita_0+T_B(j));
    G_22=A_G12*cos(sita_0+T_B(j)) - A_G22*sin(sita_0+T_B(j));
    F=-B_G1*cos(sita_0+T_B(j)) + B_G2*sin(sita_0+T_B(j));
    u_d(j)=-G_12*F/(G_11*G_22 - G_12*G_21);
    u_q(j)=G_11*F/(G_11*G_22 - G_12*G_21);
    I(:,j)=D_p(:,1)*u_d(j)+D_p(:,2)*u_q(j)+E_p;
    
end
I_a=I(1,:).*cos(T_B+sita_0)-I(2,:).*sin(T_B+sita_0);
I_b=I(1,:).*cos(T_B+sita_0-2*pi/3*ones(1,length(T_B)))-I(2,:).*sin(T_B+sita_0-2*pi/3*ones(1,length(T_B)));
I_c=I(1,:).*cos(T_B+sita_0+2*pi/3*ones(1,length(T_B)))-I(2,:).*sin(T_B+sita_0+2*pi/3*ones(1,length(T_B)));
figure(1)
plot(T,I_b);
% grid on
set(gca,'XTick',0:0.1:0.5)
set(gca,'YTick',-8:2:4)
xlabel('t(s)');
ylabel('i_b (p.u.)');
title('定子b相电流')
figure(2);
plot(T,I(3,:));
set(gca,'XTick',0:0.1:0.5)
set(gca,'YTick',0:2:2)
xlabel('t(s)');
ylabel('i_F (p.u.)');
title('励磁绕组电流')