qh1.m

一个三相异步电动机（感应电机）设计程序

function imeval
clear;
% imdesign.m contains dimension, winding, & rating data
krc=2; % 1=copper; 2=aluminum
% Winding data
VL=540; f=60; p=6; m=3; HP=50;
S1=54; S2=51; Cs=8; pitch=8; a=1;
% Conductor data
sa1=0.01563; sb=0.250*0.700; sr=0.700*0.531;
d2c=0.673; b2c=0.250;
% General dimensions
D=12.600; Dr=12.530; Do=17.62;
Dsh=4.00; dc=1.00; la=7.00; SF=0.96;
% Stator slot dimensions
bs1=0.375; ds0=0.185; ds1=0.900;
% Rotor slot dimensions
br2=0.171; br1=0.3125;
dr2=0.116; dr1=1.016;
% Squirrel cage data
Der=11.398; der=0.700; wer=0.531;
lb=7.00; br=0.00; Kring=1.00;
V1=VL/sqrt(3)+0*j; ws=2/p*2*pi*f;
Ns1=S1/p/m; Cphi=Ns1*Cs*p/a; del=(D-Dr)/2;
lam1=pi*D/S1; lam2=pi*Dr/S2;
% Fringing coefficients
qty=lam1*(5*del+bs1); ks=qty/(qty-bs1^2);
qty=lam2*(4.4*del+0.75*br2); kr=qty/(qty-br2^2);
% Stator winding factor
gam=p*pi/S1; kd=sin(Ns1*gam/2)/Ns1/sin(gam/2);
rho=pitch*p/S1*pi; kp=cos(pi/2-rho/2); kw1=kd*kp;
% Flux density checks
phim=0.97e08*a*V1/2.22/kd/kp/p/Ns1/f/Cs;
ts1=lam1-bs1; Bstm=pi*p*phim/2/S1/ts1/la/SF;
Bscm=phim/(Do-D-2*ds1)/la/SF;
tr2=pi*(Dr-2*dr1)/S2-br1; Brtm=pi*p*phim/2/S2/tr2/la/SF;
Brcm=phim/(Dr-Dsh-2*dc-2*dr1)/la/SF;
clc; disp(' '); disp(' ');
disp('FLUX DENSITY CHECKS - kilolines/sq. in'); disp(' ');
disp(['Bstm = ', num2str(Bstm/1000)]);
disp(['Bscm = ', num2str(Bscm/1000)]);
disp(['Brtm = ', num2str(Brtm/1000)]);
disp(['Brcm = ', num2str(Brcm/1000)]);
% R1 calculation
alf=asin((bs1+0.050)/lam1); tau=pi*D/p;
lalf=rho*tau/(2*pi*cos(alf)); lc=la+1+2*(lalf+ds1);
R1=1.056e-06*lc*Cphi/a/sa1; R1=1.05*R1; % ac value
% Xm calculation
dele=ks*kr*del;
Xm=1.6713e-07*f*D*la/p/dele*(kw1*S1*Cs/2/a/m)^2;
% Rc calculation
Bt13=pi*p*phim/2/(pi*(D+2*ds1/3)-S1*bs1)/la/SF;
Bscc=0.866*phim/(Do-D-2*ds1)/la/SF;
Wt=(pi/4*((D+2*ds1)^2-D^2)-S1*bs1*ds1)*la*SF*0.283;
Pt=Wt*Pc(Bt13)*4;
Wy=pi/4*(Do^2-(D+2*ds1)^2)*la*SF*0.283;
Py=Wy*Pc(Bscc)*4;
Rc=3*(0.97*VL/sqrt(3))^2/(Pt+Py);
% X1 components calculation
thetask=p*pi/S1;
Xsk=Xm*thetask^2/12;
pupitch=p*pitch/S1;       % Phase factor
if pupitch <= 0.5; Ks=pupitch;
elseif pupitch <= 2/3; Ks=-0.25+1.5*pupitch;
else; Ks=0.25+0.75*pupitch; end
Fs=(ds1-ds0)/3/bs1+ds0/bs1;
Xss=2.0055e-07*m*f*Cphi^2*la*Ks*Fs/S1;
if S2<S1; ks1=1; else; ks1=S1/S2; end
a1=abs(pi*D/S1-bs1-br2)/2;
Fz=a1/4/del+dr2*(1-ks1/2)/br2;
Xz1=2.0055e-07*m*f*Cphi^2*la*Ks*Fz/S1;
Fes=1+4*log(3.94*D/S1/sqrt(bs1*ds1));
Xes=0.3192e-07*f*(kw1*Cphi)^2*(0.5+lalf)*Ks*Fes/p;

% X2 components calculation
if S1<S2; kr1=1; else; kr1=S2/S1; end
a2=abs(pi*Dr/S2-bs1-br2)/2;
q2=S2/m/p; gamr=p*pi/S2; kdr=sin(q2*gamr/2)/q2/sin(gamr/2);
rhor=fix(S2/p)*p*pi/S2; kpr=cos((pi-rhor)/2); kwr=kdr*kpr;
Fr=(dr1-dr2)/3/br1+br2/br2;
Xrs=2.0055e-07*m*f*(kw1*Cphi)^2*la*Fr/kwr^2/S2;
Fz2=a2/4/del+ds0*(1-kr1/2)/bs1;
Xz2=2.0055e-07*m*f*(kw1*Cphi)^2*la*Fz2/kwr^2/S2;
fer1=S2*Der/m/p*(1+4*log(3.94*Dr/S2/sqrt(wer*der)));
fer2=4*br*p*(1+4*log(3.94*Dr/S2/sqrt(br1*dr1)));
Xer=7.98e-09*f*(kw1*Cphi)^2/p^2*(fer1+fer2);

I1m=zeros(1,100); PF=I1m; TTd=I1m; eta=I1m; Ps=I1m;
wm=linspace(0,ws-0.01,100); npts=length(wm);
for i=1:npts;      %  Speed loop - performance calculations
s=(ws-wm(i))/ws; fr=s*f;
if krc==1; kalf=0.3505; kbet=1; else; kalf=0.2700; kbet=1.64; end
alfd=kalf*sqrt(fr*b2c/br1)*d2c;
K1=alfd*(sinh(2*alfd)+sin(2*alfd))/(cosh(2*alfd)-cos(2*alfd));
J1=3*(sinh(2*alfd)-sin(2*alfd))/(2*alfd*(cosh(2*alfd)-cos(2*alfd)));
%  R2 calculation with skin effect adjustment
R2=kbet*1.169e-06*m*(kw1*Cphi/kwr)^2*(lb*K1/S2/sb+2*Der*Kring/pi/p^2/sr);
X1=Xsk+Xes+Xss+Xz1;
X2=Xer+Xz2+Xrs*J1;      % Adjustment for skin effects
Z1=R1+X1*j; Z2=R2/s+X2*j; ZM=(0+Rc*Xm*j)/(Rc+Xm*j);
Z=Z1*Z2+Z1*ZM+Z2*ZM; I1=V1*(Z2+ZM)/Z; I1m(1,i)=abs(I1);  % Stator current
I2=ZM/(Z2+ZM)*I1; E=I2*Z2;
PF(1,i)=cos(atan2(imag(Z/(Z2+ZM)),real(Z/(Z2+ZM)))); %  Power factor
TTd(1,i)=3*abs(I2)^2*R2/s/ws;                 % Developed torque
pfw=0.03*HP*746*(wm(i)/ws)^(2.8);  %  F&W losses
Ps(1,i)=TTd(1,i)*wm(i)-pfw; if Ps(1,i)<0; Ps(1,i)=0; end
Pin=3*real(V1*conj(I1)); eta(1,i)=Ps(1,i)/Pin*100;   % Efficiency
end;
nm=wm*30/pi; Ps=Ps/746;
% Plot performance
figure(1);
subplot(2,1,1), plot(nm,Ps); grid; title('Output power');
xlabel('Speed, rpm'); ylabel('Power, hp');
subplot(2,1,2), plot(nm,TTd);grid; title('Developed torque');
xlabel('Speed, rpm'); ylabel('Developed torque, N-m');
figure(2);
subplot(2,1,1), plot(nm,I1m);grid; title('Input current');
xlabel('Speed, rpm'); ylabel('Stator current, A');
subplot(2,1,2), plot(nm,PF);grid; title('Input power factor');
xlabel('Speed, rpm'); ylabel('Input PF');
figure(3);
plot(nm,eta);grid; title('Efficiency');
xlabel('Speed, rpm'); ylabel('Efficiency, %');
% Linear interpolation for rated conditions
for i=1:npts-1; if Ps(npts-i) >= HP; k=npts-i; break; end; end;
if i==npts-1; disp('DESIGN CANNOT PRODUCE RATED OUTPUT');
else;
nmR=nm(k+1)-(HP-Ps(k+1))*(nm(k+1)-nm(k))/(Ps(k)-Ps(k+1));
delnm=nm(k+1)-nm(k); delnm1=nm(k+1)-nmR;
I1R=I1m(k+1)+(I1m(k)-I1m(k+1))*delnm1/delnm;
TTdR=TTd(k+1)+(TTd(k)-TTd(k+1))*delnm1/delnm-0.03*HP*746*...
(nmR*pi/30)^(1.8)/ws^(2.8);
etaR=eta(k+1)+(eta(k)-eta(k+1))*delnm1/delnm;
PFR=PF(k+1)+(PF(k)-PF(k+1))*delnm1/delnm;
disp(' '); disp(' ');disp('RATED CONDITIONS');disp(' ');
disp(['Speed(rpm) = ',num2str(nmR)]);
disp(['Current(A) = ',num2str(I1R)]);
disp(['Torque(N_m) = ',num2str(TTdR)]);
disp(['Efficiency(%) = ',num2str(etaR)]);
disp(['Power factor = ',num2str(PFR)]);
end


function y = Pc(Bx)
% B-Wlb values that follow are valid for M-27,24 ga,60 Hz
B=[0 10 19 24 30 35 40 45 50 55 60 64 70 75 80 85 90 ...
95 100 105 110 115 120 125 130 135 138 165 210 ]*1000;
Wlb=[0 0.08 0.15 0.195 0.27 0.35 0.43 0.52 0.625 0.76 ...
0.89 1.0 1.2 1.35 1.55 1.75 2.0 2.25 2.55 2.8 3.15 3.5 ... 
3.85 4.2 4.7 5.2 5.6 11.0 30.0 ];
% Activate to plot B-Wlb curve
% m=29; plot(Wlb(1:m),B(1:m)); grid; pause; % Linear plot
% m=29; semilogx(wlb(2:m),B(2:m)); grid; pause; % Semilog plot    
n=length(B); k=0;
if Bx==0; k=-1; y=0; end
if Bx<0; k=-1; y=0; disp('WARNING - Bx < 0, Wlb = 0 returned'); end
if Bx>B(n); y=Wlb(n); k=-1; disp('CAUTION - Beyond B-Wlb curve');end
for i=1:n
if k==0 & (Bx-B(i))<=0; k=i; break; end
end
if k>0;
y=Wlb(k-1)+(Bx-B(k-1))/(B(k)-B(k-1))*(Wlb(k)-Wlb(k-1));  
else;
end