whole.m

无刷双馈电机建模仿真的matlab、simulink代码。控制绕组的频率一定时

R=[0.81,0,0,0,0,0;0,0.81,0,0,0,0;0,0,0.81,0,0,0;0,0,0,0.81,0,0;0,0,0,0,0.00157,0;0,0,0,0,0,0.00157];
G=[0,0.24,0,0,0,0.00267;-0.24,0,0,0,-0.00267,0;0,0,0,0.63,0,0.0043;0,0,-0.63,0,0.0043,0;0,0,0,0,0,0;0,0,0,0,0,0];
L1=[177.62,0,-101.31,0,-14842,0;0,177.62,0,101.31,0,-14842;-101.31,0,63.74,0,9106.1,0;0,101.31,0,63.74,0,-9106.1;-14842,0,9106.1,0,1.3341e+006,0;0,-14842,0,-9106.1,0,1.3341e+006];
X=[0;0;0;0;0;0];
X1=[0;0;0;0;0;0];
U=[380;0;0;0;0;0];
H=0.0001;%步长


Tem=0;
ala=0;
sun=0;
TE=0;
 k=0;
 v=8;
 c=0;  % present beid angle
 wr0=0;

 ang0=0.1;
 beid0=0;
 c=0;
 kls=1000; Bls=100;Jr=1; Jg=1; ng=4.488; Tb=0.2;
  wg0=wr0*ng;
 a11=-(Bls)/Jr;
 a12=Bls/(Jr*ng);
 a14=-kls/Jr;
 b1=1/Jr;
 a21=Bls/(ng*Jg);
 a22=-Bls/(ng*ng*Jg);
 a24=kls/(ng*Jg);
 b2=1/(Jg);
 a33=-1/Tb;
 b3=1/Tb;
 a41=1;
 a42=-1;
 R0=4;
 Ta=0;
 Tem=0;
  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   
  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  fc=10
for i=0:H:40
   
    if (i>20), c=2, end 
    if (i>40), fc=-5, end
    if (i>40), fc=-10, end
    k=k+1;
   
   wr=H*(a11*wr0+a12*wg0+a14*ang0+b1*Ta)+wr0;
 
   wg=H*(a21*wr0+a22*wg0+a24*ang0+b2*Tem)+wg0;
   beid=H*(a33*beid0+b3*c)+beid0;
   ang=(wr0-wg0/ng)*H+ang0;
      
    sun=sun+(wg+wg0)/2*H;
  %SH=314*i-3*sun;
    U(1,1)=380*cos(314*i-3*sun);
    U(2,1)=-380*sin(314*i-3*sun);
   % U(3,1)=7*cos(sun+2*pi*fc*i+pi/3);
   % U(4,1)=7*sin(sun+2*pi*fc*i+pi/3); this is the direct case
   % U(3,1)=0*cos(sun+2*pi*fc*i+pi/3);% this is the short circuit case
   % U(4,1)=0*sin(sun+2*pi*fc*i+pi/3);  above case function
  % if (fc>0 and fc<5 ) 
   %    U(3,1)=(7+20*abs(fc))*cos(sun+2*pi*fc*i+pi/3);
   % U(4,1)=(7+20*abs(fc))*sin(sun+2*pi*fc*i+pi/3);
   %end 
if (abs(fc)>=5 )
     
    %  U(3,1)=450*cos(sun+2*pi*fc*i+pi/3);
       % U(4,1)=450*sin(sun+2*pi*fc*i+pi/3);
  % if (i<=20)
     U(3,1)=35*abs(fc)*cos(sun+2*pi*fc*i+pi/3);
     U(4,1)=35*abs(fc)*sin(sun+2*pi*fc*i+pi/3);
 %  else 
  % %     U(3,1)=345*cos(sun+2*pi*fc*i+pi/3);
  %   U(4,1)=345*sin(sun+2*pi*fc*i+pi/3);
 %  end 
elseif (fc<0)
     U(3,1)=(7+14*abs(fc))*cos(sun+2*pi*fc*i+pi/3);
     U(4,1)=(7+14*abs(fc))*sin(sun+2*pi*fc*i+pi/3);
 else
     U(3,1)=(7+20*abs(fc))*cos(sun+2*pi*fc*i+pi/3);
     U(4,1)=(7+20*abs(fc))*sin(sun+2*pi*fc*i+pi/3);
end

    X1=L1*(U-wg*G*X-R*X);
    X=X+H*X1;
    Tem=3*0.00089*(X(1,1)*X(6,1)-X(2,1)*X(5,1))+0.0043*(X(3,1)*X(6,1)+X(4,1)*X(5,1));
    Tls(k)=kls*ang+Bls*(wr-wg/ng);
    s1(k)=X(1,1);
     
    
         s3(k)=X(3,1);
        
        
            s6(k)=X(6,1);
     M=wr*30/3.14;
     N=wg*30/3.14;
    b(k)=N;
    d(k)=M;
    a(k)=i;
    
    lam=wr*R0/v;
    lam1=1/(1/(lam-0.02*c)-0.003/(1+c^3));
    m= 0.73*(151/lam1-0.58*c-0.002*c^2.14-13.2)*exp(-18.4/lam1) ;
     if (m<0)
     cp(k)=0;
      else
     cp(k)=m;
       end
     cpp=cp(k);
    Ta=cpp*0.5*1.25*3.14*R0*R0*v^3/wr;
        T(k)=Ta; 
        Pmec(k)=Ta*wr;
    Pp(k)=3*1.5*(U(1,1)*X(1,1)+U(2,1)*X(2,1));
    Pc(k)=1.5*(U(3,1)*X(3,1)+U(4,1)*X(4,1));
    P1(k)=Tem*wg;
    Q(k)=3*1.5*(U(1,1)*X(2,1)-U(2,1)*X(1,1));
    an(k)=ang;
    
    % Q(k)=
    
   wr0=wr;
   wg0=wg;
   beid0=beid;
   ang0=ang;
  
           
end
%i=0:0.0001:3
plot(a,b,'r')