Математическая модель САР скорости линейного асинхронного двигателя (Z1=12) при векторном управлении



Математическая модель САР скорости линейного асинхронного двигателя (Z₁=12) при векторном управлении

Емельянов Александр Александрович, доцент;

Бесклеткин Виктор Викторович, ассистент;

Прокопьев Константин Васильевич, студент;

Бурхацкий Владимир Владимирович, студент;

Ситенков Александр Александрович, студент;

Мальцев Никита Васильевич, студент;

Авдеев Александр Сергеевич, студент;

Габзалилов Эльвир Фиргатович, студент

Российский государственный профессионально-педагогический университет

Главной целью данной работы является использование идеи векторного управления для линейного асинхронного двигателя в доступной для понимания студентами форме. Поэтому, вначале рассмотрим главный канал системы без обратных связей (рис. 1). Целесообразно выделить 9 уровней, так как в некоторых из них имеются повторяющиеся элементы с одинаковой структурой программ.

Программирование каждого уровня с одним или множеством элементов производим в Script таким образом, чтобы каждые последующие продукты программ включали предыдущие. В этом случае мы увидим процесс влияния элементов предыдущих уровней на результаты последующих. Только имея представление о роли каждого элемента в главном канале можно переходить к завершающей стадии – введению обратных связей. В основу данной работы положены математические модели линейного асинхронного двигателя, рассмотренные в работах [1] и [2].

Из многочисленных публикаций по векторному управлению, на наш взгляд, необходимо отметить работы [3], [4], [5]. Кроме того, необходимо рассмотреть все наши статьи за последний год, относящиеся к данной теме.

Элементы системы управления (ЗИ, фильтры, регуляторы) «привязаны» к вращающейся системе координат, а математическая модель линейного асинхронного двигателя к неподвижной трехфазной системе индуктора, поэтому здесь необходимо обеспечить двухэтапное преобразование напряжений:

 u_x, u_y → u_α, u_β – из вращающейся системы, связанной с циклической частотой питающего напряжения к неподвижной декартовой системе α, β;

 u_α, u_β → U_a, U_b, U_c – из неподвижной системы α, β переход к трехфазной a, b, c.

Рис. 1. Главный канал разомкнутой системы

Далее токи i_a,b,c линейного асинхронного двигателя проходят также двухэтапное преобразование: i_a,b,c → i_α,β и далее i_α,β → i_xoc,yoc. Сигналы i_xoc, i_yoc далее в замкнутой системе уже можно подать на сравнивающие устройства с сигналами задания.

Циклическая частота вращения системы координат равна циклической частоте подаваемого напряжения на двигатель, т.е. ω_к = ω_s. Связь циклической частоты ω_к = ω_s с линейной скоростью поля .

Поэтому наиболее перспективным представляется рассмотрение следующей схемы: все процессы в отношении к циклической частоте ω_к = ω_s, в том числе задатчик интенсивности, задавать не в системе v_s = f(t), а в системе ω_к = ω_s = f(t), но тогда линейную скорость подвижной части необходимо преобразовать в ω = f(v₀). На выходе математической модели ЛАД необходимо установить преобразователь линейной скорости в циклическую ω = ω_s∙(1-s), где .

Тогда все регуляторы рассчитываются по классическим формулам для векторных систем с АД [3], [4], [5]. Причем за базовую принимаются паспортные данные асинхронного двигателя с близкими по характеристикам к ЛАД (r_s, r_r, L_s, L_r) и по ним определяются постоянные времени Т_э и Т_r, лежащие в основе определения коэффициентов в регуляторах скорости и тока, и далее производится их корректировка с учетом тормозных моментов от продольного краевого эффекта. В соответствии с исследованиями, например [6], можно предположить, что в многополюсных ЛАД такой корректировки не понадобится. По-видимому, в дальнейшем встанет задача компенсации краевого эффекта на уровне систем управления. Причем компенсация возможна как на уровне вращающейся системы координат (x, y), так и в неподвижной системе a, b, c.

Условимся, что если преобразование из неподвижной системы статора (индуктора) с переменными a, b, c идет в направлении системы координат с переменными α, β и далее во вращающуюся систему с переменными x, y, то будем называть это прямым преобразованием, т. е.

a, b, c → α, β → x, y – прямое преобразование,

тогда

x, y → α, β → a, b, c – обратное преобразование.

По данной схеме преобразования все сигналы с задатчика интенсивности (ЗИ), фильтров (ФЗТ, ФОТ), регуляторов (П, И) идут во вращающейся системе координат, а математическая модель линейного асинхронного двигателя (ЛАД) сделана в неподвижной системе координат (a, b, c), поэтому необходимо произвести двухэтапное преобразование сигналов.

Уровень 1. На этом уровне расположен задатчик интенсивности ЗИ. Он задает линейную синхронную скорость поля:

,

гдеτ – полюсное деление линейного асинхронного двигателя (рассмотрим двухполюсную машину); f_c – частота питающего напряжения нарастающего по линейному закону до времени t_k, затем остающегося постоянным до конца времени переходного процесса.

% Программирование задатчика интенсивности Matlab-Script

dt=0.001; tz=9.769*10^-3;

q=2; tau=3*tz*q;

tk=2;

K=input('Длительность цикла k=');

for k=1:(K+1)

if ((k*dt >= 0) && (k*dt <= tk))

fc=k*dt*30/tk;

end;

if (k*dt > tk)

fc=30;

end;

vs=2*tau*fc;

mass_t(k)=k*dt;

mass_vs(k)=vs;

end;

% Построение графиков

figure(1);

plot(mass_t,mass_vs);

grid on;

Результат первого уровня приведен на рис. 2.

Рис. 2. Сигнал задатчика интенсивности

Уровень 2. Сигнал задания v_s поступает на регулятор скорости, на выходе которого получаем сигнал задания по проекции yi_{sy задания}. Регулятор скорости принимаем пропорциональным.

Уровень 3. Фильтры в цепи заданий и обратных связей по проекциям имеют одинаковую структуру, поэтому дадим общий вывод.

Пусть задана передаточная функция фильтра с постоянной T_μ в Simulink (рис. 3):

Рис. 3. Передаточная функция фильтра

Переведем это выражение в Script:

.

Перейдем от изображения к оригиналу:

;

.

Решаем дифференциальные уравнения методом конечных разностей Эйлера:

.

Конкретно к проекциям x и y:

;

.

Уровень 4. Пропорциональные и интегральные части регуляторов тока по проекциям одинаковы.

Пропорциональная часть K_Ix = K_Iy.

Интегральная часть:

 Проведем переход из Simulink в Script (рис. 4):

Рис. 4. Передаточная функция интегральной части регулятора тока

 Переходим к оригиналу:

 Выражаем через конечные разности:

;

;

;

;

;

;

;

;

.

Математическое моделирование u_x и u_y производится в Script. Результаты расчета даны на рис. 5.

% Программирование ux, uy в Matlab-Script

dt=0.001; tz=9.769*10^-3; tk=2;

kc=20; q=2; tau=3*tz*q;

Tm=0.005; ixzad(1)=0; iyzad(1)=0;

ixb=6; Ki=8; Ti=24;

K=input('Длительность цикла k=');

for k=1:(K+1)

if ((k*dt >= 0) && (k*dt <= tk))

fc=k*dt*30/tk;

end;

if (k*dt > tk)

fc=30;

end;

vs=2*tau*fc;

vsum=vs;

iyb=vsum*kc;

ixzad(k+1)=ixzad(k)+(ixb-ixzad(k))*dt/Tm; %Фильтр ФЗТ x

iyzad(k+1)=iyzad(k)+(iyb-iyzad(k))*dt/Tm; %Фильтр ФЗТ y

ixsum(k+1)=ixzad(k+1);

iysum(k+1)=iyzad(k+1);

u1(k+1)=ixsum(k+1)*Ki; %Пропорциональная часть x

u2(k+1)=u2(k)+ixsum(k+1)*dt/Ti; %Интегральная часть x

ux=u1(k+1)+u2(k+1);

u3(k+1)=iysum(k+1)*Ki; %Пропорциональная часть y

u4(k+1)=u4(k)+iysum(k+1)*dt/Ti; %Интегральная часть y

uy=u3(k+1)+u4(k+1);

mass_t(k)=k*dt;

mass_ux(k)=ux;

mass_uy(k)=uy;

end;

% Построение графиков

figure(2);

plot(mass_t,mass_ux,'b',mass_t,mass_uy,'r');

grid on;

Рис. 5. Напряжения u_x, u_y на выходе регуляторов тока

Уровень 5. На этом уровне необходимо произвести математическое моделирование по обратному преобразованию из вращающейся системы x, y в неподвижную α, β. В соответствии с [3]:

где определяется по следующим зависимостям:

– угловая скорость вращения системы координат;

– угол поворота.

В конечных разностях Эйлера:

или .

% Программирование ualfa, ubeta в Matlab-Script

dt=0.001; tz=9.769*10^-3; q=2; tau=3*tz*q;

tk=2; kc=20; Tm=0.005; ixzad(1)=0;

iyzad(1)=0; Ki=8; Ti=24; u2(1)=0;

u4(1)=0; ixb=6; teta(1)=0;

K=input('Длительность цикла k=');

for k=1:(K+1)

if ((k*dt >= 0) && (k*dt <= tk))

fc=k*dt*30/tk;

end;

if (k*dt > tk)

fc=30;

end;

vs=2*tau*fc;

vsum=vs;

iyb=vsum*kc;

wk=vs*pi/tau;

teta(k+1)=teta(k)+wk*dt;

ixzad(k+1)=ixzad(k)+(ixb-ixzad(k))*dt/Tm; %Фильтр ФЗТ x

iyzad(k+1)=iyzad(k)+(iyb-iyzad(k))*dt/Tm; %Фильтр ФЗТ y

ixsum(k+1)=ixzad(k+1);

iysum(k+1)=iyzad(k+1);

u1(k+1)=ixsum(k+1)*Ki; %Пропорциональная часть x

u2(k+1)=u2(k)+ixsum(k+1)*dt/Ti; %Интегральная часть x

ux=u1(k+1)+u2(k+1);

u3(k+1)=iysum(k+1)*Ki; %Пропорциональная часть y

u4(k+1)=u4(k)+iysum(k+1)*dt/Ti; %Интегральная часть y

uy=u3(k+1)+u4(k+1);

rox=cos(teta(k+1));

roy=sin(teta(k+1));

% 1 ступень обратного преобразования x,y -> alfa,beta

usalfa=rox*ux-roy*uy;

usbeta=roy*ux+rox*uy;

mass_t(k)=k*dt;

mass_usalfa(k)=usalfa;

mass_usbeta(k)=usbeta;

end;

% Построение графиков

figure(3);

plot(mass_t,mass_usalfa,'b',mass_t,mass_usbeta,'r');

grid on;

Результаты пятого уровня представлены на рис. 6.

Рис. 6. Напряжения u_α, u_β на выходе преобразователя координат (x, y → α, β)

Уровень 6. На этом уровне необходимо преобразовать сигналы u_α, u_β в трехфазную систему U_a, U_b и U_c, непосредственно подаваемых в математическую модель ЛАД.

Уравнения преобразования, в соответствии с [3], имеют следующий вид:

% Программирование Ua, Ub, Uc в Matlab-Script

dt=0.001; tz=9.769*10^-3; q=2; tau=3*tz*q;

tk=2; kc=20; Tm=0.005; ixzad(1)=0;

iyzad(1)=0; Ki=8; Ti=24; u2(1)=0;

u4(1)=0; ixb=6; teta(1)=0;

K=input('Длительность цикла k=');

for k=1:(K+1)

if ((k*dt >= 0) && (k*dt <= tk))

fc=k*dt*30/tk;

end;

if (k*dt > tk)

fc=30;

end;

vs=2*tau*fc;

vsum=vs;

iyb=vsum*kc;

wk=vs*pi/tau;

teta(k+1)=teta(k)+wk*dt;

ixzad(k+1)=ixzad(k)+(ixb-ixzad(k))*dt/Tm; %Фильтр ФЗТ x

iyzad(k+1)=iyzad(k)+(iyb-iyzad(k))*dt/Tm; %Фильтр ФЗТ y

ixsum(k+1)=ixzad(k+1);

iysum(k+1)=iyzad(k+1);

u1(k+1)=ixsum(k+1)*Ki; %Пропорциональная часть x

u2(k+1)=u2(k)+ixsum(k+1)*dt/Ti; %Интегральная часть x

ux=u1(k+1)+u2(k+1);

u3(k+1)=iysum(k+1)*Ki; %Пропорциональная часть y

u4(k+1)=u4(k)+iysum(k+1)*dt/Ti; %Интегральная часть y

uy=u3(k+1)+u4(k+1);

rox=cos(teta(k+1));

roy=sin(teta(k+1));

% 1 ступень обратного преобразования x,y -> alfa,beta

usalfa=rox*ux-roy*uy;

usbeta=roy*ux+rox*uy;

% 2 ступень обратного преобразования alfa,beta -> a,b,c

Ua=usalfa;

Ub=-(1/2)*usalfa+(sqrt(3)/2)*usbeta;

Uc=-(1/2)*usalfa-(sqrt(3)/2)*usbeta;

mass_t(k)=k*dt;

mass_Ua(k)=Ua;

mass_Ub(k)=Ub;

mass_Uc(k)=Uc;

end;

% Построение графиков

figure(4);

plot(mass_t,mass_Ua,'b',mass_t,mass_Ub,'r',mass_t,mass_Uc,'g');

grid on;

Результаты шестого уровня представлены на рис. 7.

Рис. 7. Напряжения U_a, U_b, U_c

Уровень 7. Исследование реакции линейного асинхронного двигателя (i_a, i_b, i_c) при воздействии напряжений U_a, U_b, U_c.

%Программирование ia, ib, ic на выходе ЛАД в разомкнутой системе

Rb=0.1003*10^7; rs=9.5; LsA=0.074; LsB=0.076; LsC=0.07; rr=4.6345*10^-5; Lr=0.0372*10^-5; dt=0.001;

As=rs+LsA/dt; Bs=rs+LsB/dt; Cs=rs+LsC/dt; tz=9.769*10^-3; q=2; tau=3*tz*q; m=121.6; v0=0;

wn=200; UA=wn/dt; F=0; tk=2; kc=20; Tm=0.005; Ki=8; Ti=24;

ixb=6; ixzad(1)=0; iyzad(1)=0; u2(1)=0; u4(1)=0; teta(1)=0; Fc=10;

X=zeros(24,1);

K=input('Длительность цикла k=');

for k=1:(K+1)

if ((k*dt >= 0) && (k*dt <= tk))

fc=k*dt*30/tk;

end;

if (k*dt > tk)

fc=30;

end;

v(1,k)=v0; % Создание вектор-строки для графика скорости

f(1,k)=sum(F); % Создание вектор-строки для графика усилия

i0(1,k)=X(24);

i_a(1,k)=X(21);

i_b(1,k)=X(23);

i_c(1,k)=X(22);

vs=2*tau*fc;

vsum=vs;

iyb=vsum*kc;

wk=vs*pi/tau;

teta(k+1)=teta(k)+wk*dt;

ixzad(k+1)=ixzad(k)+(ixb-ixzad(k))*dt/Tm; %Фильтр ФЗТ x

iyzad(k+1)=iyzad(k)+(iyb-iyzad(k))*dt/Tm; %Фильтр ФЗТ y

ixsum(k+1)=ixzad(k+1);

iysum(k+1)=iyzad(k+1);

u1(k+1)=ixsum(k+1)*Ki; %Пропорциональная часть x

u2(k+1)=u2(k)+ixsum(k+1)*dt/Ti; %Интегральная часть x

ux=u1(k+1)+u2(k+1);

u3(k+1)=iysum(k+1)*Ki; %Пропорциональная часть y

u4(k+1)=u4(k)+iysum(k+1)*dt/Ti; %Интегральная часть y

uy=u3(k+1)+u4(k+1);

rox=cos(teta(k+1));

roy=sin(teta(k+1));

%1 ступень обратного преобразования x,y -> alfa,beta

usalfa=rox*ux-roy*uy;

usbeta=roy*ux+rox*uy;

%2 ступень обратного преобразования alfa,beta -> a,b,c

Ua=usalfa;

Ub=-(1/2)*usalfa+(sqrt(3)/2)*usbeta;

Uc=-(1/2)*usalfa-(sqrt(3)/2)*usbeta;

ТелопрограммыЛАД [2];

% Электромагнитное усилие

F(1)=X(2)*Ir(1)/(2*tz);

for n=1:18

F(n+1)=(X(n+2)-X(n))*Ir(n+1)/(2*tz);

end;

F(20)=-X(19)*Ir(20)/(2*tz);

% Скорость

v0=v0+((sum(F)-Fc)/m)*dt;

mass_t(k)=k*dt;

mass_i_a(k)=i_a(1,k);

mass_i_b(k)=i_b(1,k);

mass_i_c(k)=i_c(1,k);

end;

% Построение графиков

figure(5);

plot(mass_t,mass_i_a,'b',mass_t,mass_i_b,'r',mass_t,mass_i_c,'g');

grid on;

Результат седьмого уровня представлен на рис. 8.

Рис. 8. Токи i_a, i_b, i_c на выходе ЛАД в разомкнутой системе

Уровень 8. Для того чтобы управлять этими токами i_a, i_b, i_c, необходимо над ними произвести двухэтапное прямое преобразование во вращающуюся систему координат.

На восьмом уровне выходные токи i_a, i_b и i_c ЛАД поступают на первый преобразователь координат, реализующий следующие математические зависимости [3]:

% Программирование ialfa, ibeta в Matlab-Script

Rb=0.1003*10^7; rs=9.5; LsA=0.074; LsB=0.076; LsC=0.07; rr=4.6345*10^-5; Lr=0.0372*10^-5; dt=0.001;

As=rs+LsA/dt; Bs=rs+LsB/dt; Cs=rs+LsC/dt; tz=9.769*10^-3; q=2; tau=3*tz*q; m=121.6; v0=0;

wn=200; UA=wn/dt; F=0; tk=2; kc=20; Tm=0.005; Ki=8; Ti=24;

ixb=6; ixzad(1)=0; iyzad(1)=0; u2(1)=0; u4(1)=0; teta(1)=0; Fc=10;

X=zeros(24,1);

K=input('Длительность цикла k=');

for k=1:(K+1)

if ((k*dt >= 0) && (k*dt <= tk))

fc=k*dt*30/tk;

end;

if (k*dt > tk)

fc=30;

end;

v(1,k)=v0; % Создание вектор-строки для графика скорости

f(1,k)=sum(F); % Создание вектор-строки для графика усилия

i0(1,k)=X(24);

i_a(1,k)=X(21);

i_b(1,k)=X(23);

i_c(1,k)=X(22);

vs=2*tau*fc;

vsum=vs;

iyb=vsum*kc;

wk=vs*pi/tau;

teta(k+1)=teta(k)+wk*dt;

ixzad(k+1)=ixzad(k)+(ixb-ixzad(k))*dt/Tm; %Фильтр ФЗТ x

iyzad(k+1)=iyzad(k)+(iyb-iyzad(k))*dt/Tm; %Фильтр ФЗТ y

ixsum(k+1)=ixzad(k+1);

iysum(k+1)=iyzad(k+1);

u1(k+1)=ixsum(k+1)*Ki; %Пропорциональная часть x

u2(k+1)=u2(k)+ixsum(k+1)*dt/Ti; %Интегральная часть x

ux=u1(k+1)+u2(k+1);

u3(k+1)=iysum(k+1)*Ki; %Пропорциональная часть y

u4(k+1)=u4(k)+iysum(k+1)*dt/Ti; %Интегральная часть y

uy=u3(k+1)+u4(k+1);

rox=cos(teta(k+1));

roy=sin(teta(k+1));

%1 ступень прямого преобразования a,b,c -> alfa,beta

i_alfa=(1/3)*(2*i_a(1,k)-i_b(1,k)-i_c(1,k));

i_beta=(1/sqrt(3))*(i_b(1,k)-i_c(1,k));

%1 ступень обратного преобразования x,y -> alfa,beta

usalfa=rox*ux-roy*uy;

usbeta=roy*ux+rox*uy;

%2 ступень обратного преобразования alfa,beta -> a,b,c

Ua=usalfa;

Ub=-(1/2)*usalfa+(sqrt(3)/2)*usbeta;

Uc=-(1/2)*usalfa-(sqrt(3)/2)*usbeta;

ТелопрограммыЛАД [2];

% Электромагнитное усилие

F(1)=X(2)*Ir(1)/(2*tz);

for n=1:18

F(n+1)=(X(n+2)-X(n))*Ir(n+1)/(2*tz);

end;

F(20)=-X(19)*Ir(20)/(2*tz);

% Скорость

v0=v0+((sum(F)-Fc)/m)*dt;

mass_t(k)=k*dt;

mass_i_alfa(k)=i_alfa;

mass_i_beta(k)=i_beta;

end;

% Построение графиков

figure(6);

plot(mass_t,mass_i_alfa,'b',mass_t,mass_i_beta,'r');

grid on;

Результат восьмого уровня представлен на рис. 9.

Рис. 9. Токи i_α, i_β в разомкнутой системе

Уровень 9. На этом уровне сигналы i_α, i_β с помощью нового прямого преобразователя координат переводятся во вращающуюся систему координат ω_k = ω_s или в синхронно бегущую систему координат .

Математическая модель прямого преобразователя имеет следующую форму [3]:

где и совпадают с ранее рассмотренным при обратном преобразовании координат.

% Программирование ixoc, iyoc в Matlab-Script

Rb=0.1003*10^7; rs=9.5; LsA=0.074; LsB=0.076; LsC=0.07; rr=4.6345*10^-5; Lr=0.0372*10^-5; dt=0.001;

As=rs+LsA/dt; Bs=rs+LsB/dt; Cs=rs+LsC/dt; tz=9.769*10^-3; q=2; tau=3*tz*q; m=121.6; v0=0;

wn=200; UA=wn/dt; F=0; tk=2; kc=20; Tm=0.005; Ki=8; Ti=24;

ixb=6; ixzad(1)=0; iyzad(1)=0; u2(1)=0; u4(1)=0; teta(1)=0; Fc=10;

X=zeros(24,1);

K=input('Длительность цикла k=');

for k=1:(K+1)

if ((k*dt >= 0) && (k*dt <= tk))

fc=k*dt*30/tk;

end;

if (k*dt > tk)

fc=30;

end;

v(1,k)=v0; % Создание вектор-строки для графика скорости

f(1,k)=sum(F); % Создание вектор-строки для графика усилия

i0(1,k)=X(24);

i_a(1,k)=X(21);

i_b(1,k)=X(23);

i_c(1,k)=X(22);

vs=2*tau*fc;

vsum=vs;

iyb=vsum*kc;

wk=vs*pi/tau;

teta(k+1)=teta(k)+wk*dt;

ixzad(k+1)=ixzad(k)+(ixb-ixzad(k))*dt/Tm; %Фильтр ФЗТ x

iyzad(k+1)=iyzad(k)+(iyb-iyzad(k))*dt/Tm; %Фильтр ФЗТ y

ixsum(k+1)=ixzad(k+1);

iysum(k+1)=iyzad(k+1);

u1(k+1)=ixsum(k+1)*Ki; %Пропорциональная часть x

u2(k+1)=u2(k)+ixsum(k+1)*dt/Ti; %Интегральная часть x

ux=u1(k+1)+u2(k+1);

u3(k+1)=iysum(k+1)*Ki; %Пропорциональная часть y

u4(k+1)=u4(k)+iysum(k+1)*dt/Ti; %Интегральная часть y

uy=u3(k+1)+u4(k+1);

rox=cos(teta(k+1));

roy=sin(teta(k+1));

% 1 ступень прямого преобразования a,b,c -> alfa,beta

i_alfa=(1/3)*(2*i_a(1,k)-i_b(1,k)-i_c(1,k));

i_beta=(1/sqrt(3))*(i_b(1,k)-i_c(1,k));

% 2 ступень прямого преобразования alfa,beta -> x,y

ixoc(k)=rox*i_alfa+roy*i_beta;

iyoc(k)=-roy*i_alfa+rox*i_beta;

% 1 ступень обратного преобразования x,y -> alfa,beta

usalfa=rox*ux-roy*uy;

usbeta=roy*ux+rox*uy;

% 2 ступень обратного преобразования alfa,beta -> a,b,c

Ua=usalfa;

Ub=-(1/2)*usalfa+(sqrt(3)/2)*usbeta;

Uc=-(1/2)*usalfa-(sqrt(3)/2)*usbeta;

ТелопрограммыЛАД[2];

% Электромагнитное усилие

F(1)=X(2)*Ir(1)/(2*tz);

for n=1:18

F(n+1)=(X(n+2)-X(n))*Ir(n+1)/(2*tz);

end;

F(20)=-X(19)*Ir(20)/(2*tz);

% Скорость

v0=v0+((sum(F)-Fc)/m)*dt;

mass_t(k)=k*dt;

mass_ixoc(k)=ixoc(k);

mass_iyoc(k)=iyoc(k);

end;

% Построение графиков

figure(7);

plot(mass_t,mass_ixoc,'b',mass_t,mass_iyoc,'r');

grid on;

Результат девятого уровня представлен на рис. 10.

Рис. 10. Токи i_{x oc}, i_{y oc} в разомкнутой системе

Наконец, приступим к главной части – замыканию всех обратных связей.

Математическое моделирование САР скорости линейного асинхронного двигателя.

Замкнутая САР скорости линейного асинхронного двигателя показана на рис. 11.

% Математическая модель САР скорости ЛАД с укладкой статорной обмотки классическим способом (z=12) с нулевым проводом

% Исходные данные асинхронного двигателя

Rb=0.1003*10^7; rs=9.5; LsA=0.074; LsB=0.076; LsC=0.07; rr=4.6345*10^-5; Lr=0.0372*10^-5; dt=0.001; v0=0;

As=rs+LsA/dt; Bs=rs+LsB/dt; Cs=rs+LsC/dt; tz=9.769*10^-3; q=2; tau=3*tz*q; m=121.6; wn=200; UA=wn/dt;

F=0; tk=2; kc=20; Tm=0.005; ixzad(1)=0; iyzad(1)=0; ix(1)=0; iy(1)=0; Ki=8;

Ti=24; u2(1)=0; u4(1)=0; ixb=6; teta(1)=0; Fc=10;

X=zeros(24,1);

K=input('Длительность цикла k=');

for k=1:(K+1)

if ((k*dt >= 0) && (k*dt <= tk))

fc=k*dt*30/tk;

end;

if (k*dt > tk)

fc=30;

end;

v(1,k)=v0; % Создание вектор-строки для графика скорости

f(1,k)=sum(F); % Создание вектор-строки для графика усилия

Рис. 11. Замкнутая САР скорости линейного асинхронного двигателя

i0(1,k)=X(24);

i_a(1,k)=X(21);

i_b(1,k)=X(23);

i_c(1,k)=X(22);

vs=2*tau*fc;

vsum=vs-v0;

iyb=vsum*kc;

wk=vs*pi/tau;

teta(k+1)=teta(k)+wk*dt;

rox=cos(teta(k+1));

roy=sin(teta(k+1));

% 1 ступень прямого преобразования a,b,c -> alfa,beta

i_alfa=(1/3)*(2*i_a(1,k)-i_b(1,k)-i_c(1,k));

i_beta=(1/sqrt(3))*(i_b(1,k)-i_c(1,k));

% 2 ступень прямого преобразования alfa,beta -> x,y

ixoc(k)=rox*i_alfa+roy*i_beta;

iyoc(k)=-roy*i_alfa+rox*i_beta;

ixoc1(k+1)=ixoc(k)+(ix-ixoc(k))*dt/Tm; %Фильтр ФОТ x

ixzad(k+1)=ixzad(k)+(ixb-ixzad(k))*dt/Tm; %Фильтр ФЗТ x

iyoc1(k+1)=iyoc(k)+(iy-iyoc(k))*dt/Tm; %Фильтр ФОТ y

iyzad(k+1)=iyzad(k)+(iyb-iyzad(k))*dt/Tm; %Фильтр ФЗТ y

ixsum(k+1)=ixzad(k+1)-ixoc1(k+1);

iysum(k+1)=iyzad(k+1)-iyoc1(k+1);

u1(k+1)=ixsum(k+1)*Ki; %Пропорциональная часть x

u2(k+1)=u2(k)+ixsum(k+1)*dt/Ti; %Интегральная часть x

ux=u1(k+1)+u2(k+1);

u3(k+1)=iysum(k+1)*Ki; %Пропорциональная часть y

u4(k+1)=u4(k)+iysum(k+1)*dt/Ti; %Интегральная часть y

uy=u3(k+1)+u4(k+1);

% 1 ступень обратного преобразования x,y -> alfa,beta

usalfa=rox*ux-roy*uy;

usbeta=roy*ux+rox*uy;

% 2 ступень обратного преобразования alfa,beta -> a,b,c

Ua=usalfa;

Ub=-(1/2)*usalfa+(sqrt(3)/2)*usbeta;

Uc=-(1/2)*usalfa-(sqrt(3)/2)*usbeta;

% Формирование матрицы А

A=zeros(24);

B=2*Rb*(rr+Lr/dt)+1/dt;

B1=6*Rb*(rr+Lr/dt)+(-4*Rb)*Lr*v0/(2*tz)+1/dt;

B2=55*Rb*(rr+Lr/dt)+(-45*Rb)*Lr*v0/(2*tz)+1/dt;

B3=550*Rb*(rr+Lr/dt)+(-450*Rb)*Lr*v0/(2*tz)+1/dt;

B4=1000*Rb*(rr+Lr/dt)+1/dt;

B5=550*Rb*(rr+Lr/dt)+450*Rb*Lr*v0/(2*tz)+1/dt;

B6=55*Rb*(rr+Lr/dt)+(45*Rb)*Lr*v0/(2*tz)+1/dt;

B7=6*Rb*(rr+Lr/dt)+(4*Rb)*Lr*v0/(2*tz)+1/dt;

C=-Rb*(rr+Lr/dt)+(2*Rb*Lr+1)*v0/(2*tz);

C1=-Rb*(rr+Lr/dt)+(6*Rb*Lr+1)*v0/(2*tz);

C2=-5*Rb*(rr+Lr/dt)+(55*Rb*Lr+1)*v0/(2*tz);

C3=-50*Rb*(rr+Lr/dt)+(550*Rb*Lr+1)*v0/(2*tz);

C4=-500*Rb*(rr+Lr/dt)+(1000*Rb*Lr+1)*v0/(2*tz);

C5=-500*Rb*(rr+Lr/dt)+(550*Rb*Lr+1)*v0/(2*tz);

C6=-50*Rb*(rr+Lr/dt)+(55*Rb*Lr+1)*v0/(2*tz);

C7=-5*Rb*(rr+Lr/dt)+(6*Rb*Lr+1)*v0/(2*tz);

D=-Rb*Lr*v0/(2*tz);

D1=5*D;D2=50*D;D3=500*D;

E=-Rb*(rr+Lr/dt)-(2*Rb*Lr+1)*v0/(2*tz);

E1=-5*Rb*(rr+Lr/dt)-(6*Rb*Lr+1)*v0/(2*tz);

E2=-50*Rb*(rr+Lr/dt)-(55*Rb*Lr+1)*v0/(2*tz);

E3=-500*Rb*(rr+Lr/dt)-(550*Rb*Lr+1)*v0/(2*tz);

E4=-500*Rb*(rr+Lr/dt)-(1000*Rb*Lr+1)*v0/(2*tz);

E5=-50*Rb*(rr+Lr/dt)-(550*Rb*Lr+1)*v0/(2*tz);

E6=-5*Rb*(rr+Lr/dt)-(55*Rb*Lr+1)*v0/(2*tz);

E7=-Rb*(rr+Lr/dt)-(6*Rb*Lr+1)*v0/(2*tz);

T=-wn*Lr*v0/(2*tz);

Y=-wn*(rr+Lr/dt);

M=Y+T;

N=Y-T;

W1=-wn*Lr/dt;

P=-Rb*Lr/dt;

Q=(2*Rb*Lr+1)/dt;

Q1=(6*Rb*Lr+1)/dt;

Q2=(55*Rb*Lr+1)/dt;

Q3=(550*Rb*Lr+1)/dt;

Q4=(1000*Rb*Lr+1)/dt;

for n=1:3

A(2*n+2,n+20)=(-1)^(n+1)*T;

A(2*n+3,n+20)=(-1)^(n+1)*M;

A(2*n+4,n+20)=(-1)^(n+1)*N;

A(2*n+5,n+20)=(-1)^n*T;

A(2*n+8,n+20)=(-1)^n*T;

A(2*n+9,n+20)=(-1)^n*M;

A(2*n+10,n+20)=(-1)^n*N;

A(2*n+11,n+20)=(-1)^(n+1)*T;

end;

for n=1:3

A(24,n+20)=1;

end;

A(24,24)=-1;

for n=1:12

A(n+4,n+4)=B;

A(n+5,n+4)=E;

A(n+3,n+4)=C;

end;

for n=1:13

A(n+2,n+4)=D;

A(n+5,n+3)=-D;

end;

A(1,1)=B4; A(1,2)=C5; A(1,3)=D2; A(2,1)=E4; A(2,2)=B5; A(2,3)=C6; A(2,4)=D1; A(3,1)=-D3;

A(3,2)=E5; A(3,3)=B6; A(3,4)=C7; A(4,2)=-D2; A(4,3)=E6; A(4,4)=B7; A(5,3)=-D1; A(5,4)=E7;

A(16,17)=C1; A(16,18)=D1; A(17,17)=B1; A(17,18)=C2; A(17,19)=D2; A(18,17)=E1; A(18,18)=B2; A(18,19)=C3;

A(18,20)=D3; A(19,17)=-D1; A(19,18)=E2; A(19,19)=B3; A(19,20)=C4; A(20,18)=-D2; A(20,19)=E3; A(20,20)=B4;

for n=1:2

A(21,n+4)=UA;

A(21,n+10)=-UA;

A(22,n+8)=UA;

A(22,n+14)=-UA;

A(23,n+6)=-UA;

A(23,n+12)=UA;

end;

A(21,21)=As;

A(22,23)=Bs;

A(23,22)=Cs;

% Матрица свободных членов

S=[ Q4*X(1)+P*( 500*X(2)); %1

Q3*X(2)+P*(500*X(1)+50*X(3)); %2

Q2*X(3)+P*(50*X(2)+5*X(4)); %3

Q1*X(4)+P*(5*X(3)+X(5)); %4

W1*X(21)+Q*X(5)+P*(X(4)+X(6)); %5

W1*X(21)+Q*X(6)+P*(X(5)+X(7)); %6

(-1)*W1*X(22)+Q*X(7)+P*(X(6)+X(8)); %7

(-1)*W1*X(22)+Q*X(8)+P*(X(7)+X(9)); %8

W1*X(23)+Q*X(9)+P*(X(8)+X(10)); %9

W1*X(23)+Q*X(10)+P*(X(9)+X(11)); %10

(-1)*W1*X(21)+Q*X(11)+P*(X(10)+X(12)); %11

(-1)*W1*X(21)+Q*X(12)+P*(X(11)+X(13)); %12

W1*X(22)+Q*X(13)+P*(X(12)+X(14)); %13

W1*X(22)+Q*X(14)+P*(X(13)+X(15)); %14

(-1)*W1*X(23)+Q*X(15)+P*(X(14)+X(16)); %15

(-1)*W1*X(23)+Q*X(16)+P*(X(15)+X(17)); %16

Q1*X(17)+P*(X(16)+5*X(18)); %17

Q2*X(18)+P*(5*X(17)+50*X(19)); %18

Q3*X(19)+P*(50*X(18)+500*X(20)); %19

Q4*X(20)+P*500*X(19); %20

UA*(X(5)+X(6)-X(11)-X(12))+(LsA/dt)*X(21)+Ua; %21

UA*(X(9)+X(10)-X(15)-X(16))+(LsB/dt)*X(23)+Ub; %22

UA*(X(13)+X(14)-X(7)-X(8))+(LsC/dt)*X(22)+Uc; %23

0]; %24

% Решение методом Гаусса-Жордана

Z=rref([A S]); %Приведение расширенной матрицы к треугольному виду

X=Z(1:24,25:25); %Выделение последнего столбца из матрицы

% Ток в роторе

Ir=[ 1000*Rb*X(1)-Rb*(500*X(2)); %1

550*Rb*X(2)-Rb*(500*X(1)+50*X(3)); %2

55*Rb*X(3)-Rb*(50*X(2)+5*X(4)); %3

6*Rb*X(4)-Rb*(5*X(3)+X(5)); %4

-wn*X(21)+2*Rb*X(5)-Rb*(X(4)+X(6)); %5

-wn*X(21)+2*Rb*X(6)-Rb*(X(5)+X(7)); %6

(-1)*(-wn)*X(22)+2*Rb*X(7)-Rb*(X(6)+X(8)); %7

(-1)*(-wn)*X(22)+2*Rb*X(8)-Rb*(X(7)+X(9)); %8

-wn*X(23)+2*Rb*X(9)-Rb*(X(8)+X(10)); %9

-wn*X(23)+2*Rb*X(10)-Rb*(X(9)+X(11)); %10

(-1)*(-wn)*X(21)+2*Rb*X(11)-Rb*(X(10)+X(12)); %11

(-1)*(-wn)*X(21)+2*Rb*X(12)-Rb*(X(11)+X(13)); %12

-wn*X(22)+2*Rb*X(13)-Rb*(X(12)+X(14)); %13

-wn*X(22)+2*Rb*X(14)-Rb*(X(13)+X(15)); %14

(-1)*(-wn)*X(23)+2*Rb*X(15)-Rb*(X(14)+X(16)); %15

(-1)*(-wn)*X(23)+2*Rb*X(16)-Rb*(X(15)+X(17)); %16

6*Rb*X(17)-Rb*(X(16)+5*X(18)); %17

55*Rb*X(18)-Rb*(5*X(17)+50*X(19)); %18

550*Rb*X(19)-Rb*(50*X(18)+500*X(20)); %19

1000*Rb*X(20)-Rb*(500*X(19))]; %20

% Электромагнитное усилие

F(1)=X(2)*Ir(1)/(2*tz);

for n=1:18

F(n+1)=(X(n+2)-X(n))*Ir(n+1)/(2*tz);

end;

F(20)=-X(19)*Ir(20)/(2*tz);

% Скорость

v0=v0+((sum(F)-Fc)/m)*dt;

mass_t(k)=k*dt;

mass_ixb(k)=ixb;

mass_ixoc(k)=ixoc(k);

mass_iyoc(k)=iyoc(k);

mass_ux(k)=ux;

mass_uy(k)=uy;

mass_usalfa(k)=usalfa;

mass_usbeta(k)=usbeta;

mass_Ua(k)=Ua;

mass_Ub(k)=Ub;

mass_Uc(k)=Uc;

mass_i_a(k)=i_a(1,k);

mass_i_b(k)=i_b(1,k);

mass_i_c(k)=i_c(1,k);

mass_i_alfa(k)=i_alfa;

mass_i_beta(k)=i_beta;

end;

% Построение графиков

figure(1);

plot(mass_t,mass_ixb,'b',mass_t,mass_ixoc,'r');

grid on;

figure(2);

plot(mass_t,mass_iyb,'b',mass_t,mass_iyoc,'r');

grid on;

figure(3);

plot(mass_t,mass_ux,'b',mass_t,mass_uy,'r');

grid on;

figure(4);

plot(mass_t,mass_usalfa,'b',mass_t,mass_usbeta,'r');

grid on;

figure(5);

plot(mass_t,mass_Ua,'b',mass_t,mass_Ub,'r',mass_t,mass_Uc,'g');

grid on;

figure(6);

plot(mass_t,mass_i_a,'b',mass_t,mass_i_b,'r',mass_t,mass_i_c,'g');

grid on;

figure(7);

plot(mass_t,mass_i_alfa,'b',mass_t,mass_i_beta,'r');

grid on;

figure(8);

plot(mass_t,mass_ixoc,'b',mass_t,mass_iyoc,'r');

grid on;

figure(9);

k=0:K;

subplot(2,1,1);

plot(k*dt,v);

title('Линейная скорость');

xlabel('t,c');

ylabel('v,m/c');

grid on;

subplot(2,1,2);

plot(k*dt,f);

title('Электромагнитное усилие');

xlabel('t,c');

ylabel('F,H');

grid on;

Результаты расчетов в замкнутой системе приведены на рис. 12, …, 20.

Рис. 12. Токи i_{x зад.} и i_{x oc} на входе фильтров ФЗТ и ФОТ

Рис. 13. Токи i_{y зад.} и i_{y oc} на входе фильтров ФЗТ и ФОТ

Рис. 14. Напряжения u_x, u_y на выходе регуляторов тока

Рис. 15. Напряжения u_α, u_β на выходе преобразователя координат (x, y → α, β)

Рис. 16. Напряжения U_a, U_b, U_c на выходе преобразователя координат (α, β → a, b, c)

Рис. 17. Токи i_а, i_b, i_c на выходе ЛАД

Рис. 18. Токи i_α, i_β на выходе преобразователя (a, b, c → α, β)

Рис. 19. Токи i_{x ос}, i_{y ос} на выходе прямого преобразователя координат

Рис. 20. График скорости и момента в замкнутой системе

Литература:

1. Сарапулов Ф. Н., Емельянов А. А., Иваницкий С. В., Резин М. Г. Исследование электромеханических переходных процессов линейного асинхронного короткозамкнутого двигателя // Электричество. — 1982. — № 10. — С. 54–57.
2. Емельянов А. А., Кобзев А. В., Козлов А. М., Бесклеткин В. В., Бочкарев Ю. П., Авдеев А. С., Киряков Г. А., Чернов М. В., Габзалилов Э. Ф., Иванин А. Ю. Программирование линейного асинхронного двигателя (Z1=12) с трехфазной обмоткой индуктора с нулевым проводом // Молодой ученый. — 2014. — № 3. — С. 28–47.
3. Шрейнер Р. Т. Математическое моделирование электроприводов переменного тока с полупроводниковыми преобразователями частоты. — Екатеринбург: УРО РАН, 2000. — 654 с.
4. Шрейнер Р. Т., Емельянов А. А., Медведев А. В. Оптимизация перемежающихся режимов работы частотно-регулируемого асинхронного электропривода // Известия ТулГУ. Технические науки. Вып. 3: в 5 ч. Туда: Изд-во ТулГУ, 2010. Ч 4. — 262 с.
5. Шрейнер Р. Т., Емельянов А. А., Медведев А. В. Ресурсы энергосбережения в повторно-кратковременных режимах работы асинхронного привода // Промышленная энергетика. — 2011. — № 11. — С. 22–27.
6. Веселовский О. Н., Коняев А. Ю., Сарапулов Ф. Н. Линейные асинхронные двигатели. — М.: — Энергоатомиздат, 1991. — 256 с.