Основное

1. Передаточная функция резонансного регулятора в непрерывной форме. Коэффициент Kr имеет тот же смысл что интегральная составляющая в ПИ-регуляторе

   

2. Резонансный регулятор обеспечивает бесконечно большой коэффициент усиления на частоте f₀, обеспечивая астатическое регулирование одновременно по амплитуде и фазе сигнала с этой фиксированной частотой.

Что это такое и для чего необходимо

Иногда возникает необходимость стабилизировать напряжение или ток преобразователя, например, бесперебойного источника питания, солнечного инвертора и иных устройств где на выходе имеется синусоида.

Для трёхфазных систем обычно используется преобразование Парка-Горева. Вначале трёхфазный сигнал ABC преобразуется в ортогональный αβ0 (плюс нулевая последовательность) посредством преобразования Кларке, который, затем, посредством преобразованием Парка-Горева, становится проекциями на синхронно вращающуюся систему координат, представляющие собой сигналы с постоянной составляющей, содержащие информацию о фазе и амплитуде исходной трёхфазной системы. В идеальном варианте симметричной системы, содержащей только прямую последовательность, эти сигналы не содержат вторую гармонику. Управление по симметричным компонентам - отдельная тема, следует отметить, что резонансный регулятор позволяет реализовать независимое управление по отдельным фазам.

В однофазных системах используется восстановление ортогональной составляющей производной/интегратором/фазовым фильтром либо специальные методы по измерению амплитуды и относительной фазы требуемых величин применяя счётчики и компараторы, однако это не всегда надёжно и "плавно". В этом случае на помощь приходит резонансный регулятор, где изначально в качестве сигнала задания используется синусоидальный сигнал с необходимой амплитудой и фазой.

Анализ замкнутой системы с Р-регулятором в непрерывной форме

Общая структура представлена на рис. 1. Имеется входной синусоидальный сигнал, ПР-регулятор, объект управления в общем виде.

С использованием предложенной методики покажем равенство фазы и амплитуды в установившемся режиме. Общий скрипт для вычислений в системе символьных вычислений Maxima:

/*Передаточная функция ПР-регулятора*/
V1:Kr*p/(p^2+(2*%pi*f_0)^2);
/*Уравнение замкнутой системы*/
V2:[e=x-y,y=Fk(p)*V1*e];

Вначале задаётся вид регулятора, составляется система уравнений для замкнутой системы на рис. 1.

Далее формируется решение, подставляется в него значение входного сигнала x с целью получения переходной характеристики:

/*Решение системы относительно выходной величины*/
V3:factor(y),solve(V2,[e,y]);
/*Задание входного сигнала во временной форме*/
V4:A*sin(2*%pi*f_0*t+ψ_0);
/*Преобразование Лапласа от входного сигнала*/
V5:x=laplace(V4,t,p);
/*Подстановка входного сигнала в решение системы уравнений
- переходная характеристика*/
V6:factor(V3),V5;

Производится умножение на коэффициент перехода к комплексно-частотной точке в установившемся режиме

/*Задание коэффициента перехода к комплексно-частотной точке*/
V7:(p^2+4*%pi^2*f_0^2)/(2*%pi*f_0);
/*Умножение на коэффициент и упрощение*/;
V8:factor(V6*V7);
/*Переход к комплексно-частотной характеристике в заданной точке*/
V9:p = %i*2*%pi*f_0;
/*Комлпексное значение в установившемся режиме*/
V10:factor(V8),V9;

Производится подстановка для перехода к комплексно-частотной точке в установившемся режиме

Легко заметить что получается ничто иное как формула Эйлера. Таким образом, амплитуда выходного сигнала в установившемся режиме равна A и фаза равна ψ₀, что и является непосредственно входным сигналом. Обратить внимание что выходной сигнал не зависит от передаточной функции объекта управления в установившемся режиме по интересующей гармонике.

Чтобы это показать наглядно, зададим необходимую функцию для объекта управления. Пусть это будет апериодическое звено первого порядка вида ФНЧ с постоянной времени τ=0.1 секунды, определяется как 1/(τ•p+1).
Остальные параметры соответственно равны Kr = 1000, ψ₀ = -0.2, f₀ = 50, A = 1.3.

Далее приведены результаты численного анализа для схемы, составленной в среде Scilab+xCos.

Параметры моделирования и исходные данные

Результаты моделирования. В разомкнутой системе амплитуда сигнала достаточно мала и соответствует прохождению через ФНЧ, при этом фаза является отстающей и примерно равна 90°.

Получим тот же результат для проверки с использованием Maxima, продолжим скрипт с использованием вспомогательной библиотеки для численных методов, позволяющих осуществить обратное преобразование Лапласа для функций высших порядков.

/*Задать форму объекта управления Fk(p)*/
V11:Fk(p)=1/(τ*p+1);
/*Подставить объект управления в переходную характеристику*/
V12:factor(V6),V11;
/*Задать численные параметры*/
Num:[A=1.3,Kr=2000,ψ_0=-0.2,f_0=50,τ=0.02];
/*переходная характеристика в численном виде*/
V13:V12,Num;
/*Загрузить пакет для работы с обратным преобразованием
Лапласа в численном виде*/
load("coma");
/*получить обратное преобразование Лапласа */
V14:nilt(V13,p,t);
/*Входной сигнал с численными значениями*/
V15:V4,Num;
/*Построение графика переходного процесса*/
wxplot2d([V15,V14],[t,0,0.4],[color,green,black],[legend,false]),wxplot_size=[480,320];

Задание объекта управления и получение переходной характеристики

После чего можно подставить численные параметры, загрузить пакет "coma" для работы с функциями высших порядков (более 4-го). Команда по-умолчанию ilt работает до 4 порядка включительно, для преодоления этого ограничения была разработана nilt (numerical inverse Laplace transform).

Результат представлен далее, видно, что выходной сигнал стремится к необходимой фазе и амплитуде 1.3.

Таким образом, резонансный регулятор позволяет реализовать астатическое (с нулевой ошибкой) управление как по амплитуде, так и по фазе сигнала управляющего воздействия, обеспечивая в установившемся режиме их равенство.

Построим частотную характеристику регулятора, чтобы показать бесконечную добротность на частоте управления.

/*Переход к комплексно-частотной характеристике (КЧХ) для всех частот*/
V16:p = %i*2*%pi*f;
/*КЧХ регулятора*/
V17:V1,V16;
/*АЧХ и ФЧХ регулятора, предположить что параметры неотрицательные для упрощения не содержащего знак модуля*/;
assume(Kr>0,f>0,f_0>0);
/*Действительная и мнимая части КЧХ*/
V17R:realpart(V17);
V17I:factor(imagpart(V17));
/*АЧХ и ФЧХ регулятора, предположить что параметры неотрицательные для упрощения не содержащего знак модуля*/;
assume(Kr>0,f>0,f_0>0);
V17A:factor(cabs(V17));
V17P:carg(V17);
/*Подстановка численных параметров */
V18A:V17A,Num;V18P:V17P,Num;
/*График АЧХ*/
wxplot2d(V18A,[f,0,400],[y,0,10],[grid2d,true],[color,green,black],[legend,false],[ylabel,"|f(j2πf)|"]),wxplot_size=[480,320];
/*Построение графика ФЧХ*/
wxplot2d(V18P*180/%pi,[f,0,400],[y,0,360],[grid2d,true],[ylabel,"∠f(j2πf)"],[color,green,black],[legend,false]),wxplot_size=[480,320];
/*Построение графика мнимой и действительной частей*/
wxplot2d([V18R,V18I],[f,0,400],[y,-10,10],[grid2d,true],[ylabel,"Re,Im f(j2πf)"],[color,green,black],[legend,false]),wxplot_size=[480,320];

Осуществляется подстановка для оператора Лапласа, осуществляющая переход в комплексно-частотную область для всех значений частоты, далее осуществляется взятие модуля, угла и нахождение действительной и мнимой частей КЧХ, обратить внимание, что КЧХ - чисто мнимая.

На рис.14 представлены результирующие характеристики. Обратить внимание на нулевую действительную часть КЧХ и асимптоты, соответствующие частоте 50 Гц.

Цифровая форма реализации

Для получения регулятора в дискретном виде воспользуемся методом согласования нулей и полюсов. Об этом методе будет отдельная публикация равно как и о реализации в Maxima цифровых фильтров.

Для начала используем следующий код в новом скрипте

/*Инициализация*/
kill(all);ratprint:false;fpprintprec:4;
/*Резонансный регулятор*/
V1:Kr*p/(p^2+(2*%pi*f_0)^2);
/*Знаменатель и числитель*/
V1n:num(V1);V1d:denom(V1);

Общая схема (рис. 15) метода нулей и полюсов представлена на рисунке, здесь и далее fd - частота дискретизации

Далее находятся корни числителя и знаменателя, записывается передаточная функция в общем виде исходя из значения корней. Обратить внимание на использование функции demoivre, которая преобразует комплексные экспоненты в синусы-косинусы согласно формулы де Муавра.

/*Корни знаменателя и числителя*/
Kn:solve(V1n,p);Kd:solve(V1d,p);
/*Нахождение произведения разностей вида z-exp(p[i]/fd) для i-го корня*/
/*Для числителя*/
Zn:product(subst(Kn[i],(z-exp(p/fd))),i,1,length(Kn));
/*Для знаменателя*/
Zd:trigsimp(demoivre(product(subst(Kd[i],(z-exp(p/fd))),i,1,length(Kd))));
/*Передаточная функция в общем виде*/
V2:Krd*Zn/Zd;
/*Сравнение частотных характеристик*/;
/*Переход в КЧХ для непрерывной функции*/
V3:p=%i*2*%pi*f;
/*Переход в КЧХ для дискретной функции*/
V4:z=exp(%i*2*%pi*f/fd);
/*Подстановка в непрерывную ПФ*/
V5:V1,V3;
/*Подстановка в дискретную ПФ*/
V6:factor(V2),V4;

Вывод представлен на рисунке 16 в виде последовательности выражений

Следует отметить, что билинейное преобразование (метод Тастина) смещает частоты, поэтому используется именно метод согласования нулей и полюсов, дающий точные частотные точки, однако, при этом, искажается амплитуда. Для этого необходимо найти вспомогательный коэффициент Krd, который приближает непрерывную функцию к дискретной в заданной частотной точке по амплитуде.

/*Подстановка в непрерывную ПФ*/
V5:V1,V3;
/*Подстановка в дискретную ПФ*/
V6:factor(V2),V4;
/*Потребуем, чтобы скорость изменения составляющих КЧХ была одинаковой
на нулевой частоте. Используется производная*/
DV5:subst(f=0,diff(V5,f));
DV6:subst(f=0,diff(V6,f)); 
/*Решение для коэффициента, удовлетворяющего равенству производных на нулевой частоте*/
V7:solve(DV5=DV6,Krd);
/*Зададим численные параметры*/
Num:[Kr=2000,f_0=50,fd=1000];
/*Нахождение вспомогательного коэффициента*/
V7N:float(V7),Num;
/*Общий вектор численных параметров*/
Num:append(V7N,Num);
/*Упрощение выражений для АЧХ резонансного регулятора 
в непрерывном и дискретном виде*/
assume(f>0,f_0>0,fd>0,Kr>0,Krd>0);
V5A:factor(cabs(V5));V6A:factor(trigreduce(trigsimp(factor(cabs(V6)))));
/*Подстановка численных параметров*/
V5N:float(V5A),Num;V6N:float(V6A),Num;
/*Построение АЧХ*/
wxplot2d([V5N,V6N],[f,0,2000],[y,0,5],[grid2d,true],[ylabel,"|Fd(p)|, |Fz(p)|"],[color,green,black],[legend,false]),wxplot_size=[480,320];

При обычном методе где на нулевой частоте отличное от нуля значение, например, для ФНЧ, вспомогательный коэффициент находится подстановкой p=0 и z=1. Если же необходимо обеспечить одинаковый рост частотной характеристики вблизи этой точки, особенно если она равна нулю, используется производная. Порядок нахождения представлен на рис. 17.

После подстановки этого коэффициента в дискретную передаточную функцию, подстановки численных параметров, можно построить следующий график АЧХ (рис. 18)

Далее реализуем цифровой фильтр и сравним результат во временной форме. Получим сперва результат обратного преобразования Лапласа от воздействия на резонансный регулятор ступенчатой функции, он довольно прост: (20•sin(100•π•t))/π. Видно, что это незатухающие синусоидальные колебания (благодаря бесконечной добротности).

/*Найдём результат обратного преобразования Лапласа
от воздействия ступеньки на резонансный регулятор*/
W1:ilt(V1/p,p,t),Num;
/*Подставляем значение для перехода от дискретной передаточной функции
к цифровому фильтру*/
V7:factor(V2),z=1/zm;
/*Раскрыть скобки у полиномов числителя и знаменателя*/
V8n:expand(num(V7));V8d:expand(denom(V7));
/*Выделить для числителя и знаменателя коэффициенты
при степенях zm*/
V9cn:makelist(coeff(V8n,zm,i),i,0,hipow(V8n,zm));
V9cd:makelist(coeff(V8d,zm,i),i,0,hipow(V8d,zm));
/*Разделить числитель и знаменатель на 
свободное слагаемое знаменателя (при нулевой степени zm)
*/
V10cn:makelist(i/V9cd[1],i,V9cn);
V10cd:makelist(i/V9cd[1],i,V9cd);
/*Подставить численные параметры*/
V11N:float(V10cn),Num;
V11D:float(V10cd),Num;

Цифровому фильтру будет уделено внимание в отдельной публикации, поэтому кратко можно описать следующий процесс, как показано на рисунке 12.

Переход:

• подставить zm = z⁻¹ в дискретную передаточную функцию и упростить

• выделить полиномы числителя и знаменателя

• сформировать массив коэффициентов при zm начиная с нулевой степени (свободное слагаемое)

• разделить коэффициенты числителя и знаменателя на свободное слагаемое знаменателя, тем самым, произведя нормировку цифрового фильтра
  

Далее следует скрипт, который позволяет реализовать цифровой фильтр в Maxima. Более подробно в следующей публикации.

/*Коэффициенты числителя и знаменателя*/
CnumN:V11N;CdenN:V11D;
DigFiltReset(Filter):=block
(
    [Nx,Ny,i],
    Cx:Filter[3],Cy:Filter[4],
    Nx:length(Cx),Ny:length(Cy),
    for i:1 thru Nx do (Filter[1][i]:0),
    for i:1 thru Ny do (Filter[2][i]:0),
    Filter
);
DigFiltInit(Cx,Cy):=block
(
    [Fliter,X,Y,Nx,Ny],
    Nx:length(Cx),Ny:length(Cy),Filter:makelist(0,4),   
    X:makelist(0,Nx),Y:makelist(0,Ny),Cy[1]:0,
    Filter[1]:X,Filter[2]:Y,Filter[3]:Cx,Filter[4]:Cy,
    Filter
);
DigFiltStep(Filter,InputX):=block
(
[Acc,i,T1,T2,Cx,Cy,X,Y],
X:Filter[1],Y:Filter[2],
Cx:Filter[3],Cy:Filter[4],
Nx:length(Cx),Ny:length(Cy),
X[1]:InputX,
Cy[1]:0,Acc:0,
(
  for i:1 thru Nx do
    (Acc:float(Acc+X[i]*Cx[i])),
  for i:1 thru Ny do
    (Acc:float(Acc-Y[i]*Cy[i])),
  Y[1]:Acc,
  T1:X[1],
  for i:1 thru Nx do
    (T2:X[i],X[i]:T1,T1:T2),
 T1:Y[1],        
  for i:1 thru Ny do
    (T2:Y[i],Y[i]:T1,T1:T2)        
),       
Filter[1]:X,Filter[2]:Y,
Filter[3]:Cx,Filter[4]:Cy,
Filter
);
DigFiltCycle(Filter,Function,Niter,Fd):=block
(
[Fst,j,Out,OutIn,Fv],Fst:Filter,
DigFiltReset(Fst),
Out:makelist([0,0],i,Niter),OutIn:makelist([0,0],i,Niter),
for j:1 thru Niter do
(
    Fv:float(ev(subst(t=(j-1)/Fd,Function))),    
    Vhod:1,
    DigFiltStep(Fst,Vhod),
    Out[j][2]:Fst[2][1],
    Out[j][1]:j,       
    OutIn[j][2]:Fv,
    OutIn[j][1]:j           
) ,
[Out,OutIn]           
);
/*Инициализация переменных состояния фильтра и коэффициентов*/
Fa:DigFiltInit(CnumN,CdenN);
Out:DigFiltCycle(Fa,W1,50,subst(Num,fd))$
wxplot2d([[discrete, Out[1]],[discrete, Out[2]]],[style,points,lines],[y,-10,10],[grid2d,true],[ylabel,"Fp(t),Fz(t)"],[color,green,black],[legend,false]),wxplot_size=[480,320];

В итоге получается требуемый отклик во временной области:

Коэффициент цифрового фильтра 2 cos((2•π•f₀)/fd) можно рассчитывать непосредственно как константу или используя f₀ от блока ФАПЧ.