Численная оптимизация ПИД-регуляторов с использованием детектора правильности движения в целевой функции
Жмудь В. А.1' 2'3, Ядрышников О. Д.:
1ФГБОУ ВПО НГТУ,
2НИУ НГУ, 3ОАО «НИПС», Россия
Аннотация. Оптимизация регулятора в системах, склонных к колебаниям вследствие специфичности модели объекта трудно осуществить численными методами даже при наличии хорошо апробированной методики. В статье предложены способы модификации целевой функции, позволяющие эффективно решить эту задачу. Результативность предложенного подхода продемонстрирована на примерах.1
Ключевые слова: численная оптимизация, регуляторы, автоматика, моделирование, динамические системы, точность управления
ВВЕДЕНИЕ
Высокотехнологическое производство
научные исследования и другие отрасли технической деятельности требуют
совершенствования методов управления в контуре с отрицательной обратной связью. В таких системах для достижения требуемой точности и качества управления применяются ПИД-регуляторы или их модификации [1-5]. Вычисление параметров настройки регуляторов может быть осуществлено, например, методами численной оптимизации [2-9]. Ряд объектов обладают математической моделью, которая затрудняет применение этого способа. К таким особенностям можно отнести высокую склонность к колебаниям, особенно, при замыкании отрицательной обратной связи.
В данной статье рассматривается несколько приемов, развивающих методику численной оптимизации регуляторов для таких систем.
I. ПОСТАНОВКА ЗАДАЧИ
Общий вид ПИД-регулятора задается уравнением:
Работа выполнена по заданию Министерства образования и науки РФ, проект №7.599.2011, Темплан, НИР № 01201255056.
Г^) = КР + К / s + Кс s. (1)
Здесь КР, Кь Кс - коэффициенты усиления различных каналов регулятора, 5 - аргумент преобразования Лапласа или символический оператор дифференцирования в случае использования дифференциальных уравнений. Многие статьи предлагают определение коэффициентов методом Зигеля-Никольса [4]. Этот выбор очевидно не оптимален, поскольку этот метод не самый эффективный для достижения наилучшего качества управления.
Нами ранее разработана методика расчета этих коэффициентов путем моделирования в программе VisSim и численной оптимизации по критерию, содержащему интеграл от модуля ошибки управления, умноженному на время с момента начала переходного процесса.
У (Т) = 11 е(1) 1Ш.
(2)
Применение такого метода не дает результата по следующим причинам.
1. Процедура оптимизации прерывается вследствие достижения какой-либо вычисленной величиной максимально допустимого для программы значения 10303. Данный вид прерывания не дает возможности получения хотя бы промежуточного результата оптимизации.
2. Если удается осуществить оптимизацию без проблемы, указанной в предыдущем пункте, результат оптимизации дает систему с неудовлетворительным переходным процессом, содержащим слишком большое перерегулирование и (или) множество колебаний.
II. МЕТОДЫ МОДИФИКАЦИИ ЦЕЛЕВОЙ ФУНКЦИИ
В целевую функцию (2) предлагается внести ограничение ее величины по модулю большим значением, например, 108. Это предотвратит прерывание по пункту 1.
0
^(Г, У) =
108, 108; У,-108 < 108; -108, У<-108.
(3)
В целевую функцию (2) также предлагается внести под интеграл произведение ошибки на ее производную, умноженную на большой коэффициент усиления. Логика этой модификации заключается в том, что если ошибка и ее производная совпадают по знаку, это означает неправильное развитие переходного процесса: ошибка возрастает, когда должна убывать или убывает, когда должна возрастать. Наличие таких участков переходного процесса нежелательно, такие участки возникают при перерегулировании и при колебаниях. т
У 2(Г) = |{Я[в(г)йе(г)/ Ж]+1 е(г) | г}йг. (4)
ГI, I > 0;
Я[/] = \ (5)
и [ 0,1 < 0.
Уз(Г) = ^(Г, (6)
Здесь Я - детектор колебаний.
Пример 1. Рассмотрим объект с передаточной функцией
^(Б) = (- 6 Б + 3) / (5 Б2 + 3 Б + 1). (7)
На рис. 1 показана схема моделирования и оптимизации системы с объектом (3) и регулятором (1). Применена целевая функция (4) с ограничением.
Результатом оптимизации является следующий набор коэффициентов регулятора: КР = 0,67; К = 0,1259; Кс = 0,787. Целевая функция при этом равна 1,757 на интервале до 100 с.
Этот результат не является робастным, то есть при небольших изменениях параметра объекта или регулятора устойчивость нарушается. В переходном процессе имеется небольшое перерегулирование, процесс происходит с несколькими колебаниями.
г-Н!]-ИЗ)-1 +>)
—1|-*0>-КЖ1——^—►
'—ИЗ)-М ^епуаНуе |—^
-бэ+З
5Б'+ЗВ+1
0.02
I 0 035 Ь
ГоТЦ-
-►I рагап^егЦпкшмт |-►П>Т-
-►I рагате1егЦпкгижп |-И~П—
-►I рагате1егЦпкгимт |-НЗ-
1.75751
.670041
.125972
.787021
-►[сёЦ]
-+Ш-
20 1.5
1 0
^ -5
О
-.5
-1 0
\\
20 40 60 "Пте (э)
80
100
Рис. 1. Структура для моделирования и оптимизации задачи Примера 1 с ограничением целевой функции и результат оптимизации
Рис. 2. Структура с модификацией путем введения детектора колебаний и результат оптимизации
0
На рис. 2 показана структура для решения той же задачи с применением детектора колебаний на основе произведения ошибки и ее производной (6). В результате оптимизации получены другие коэффициенты регулятора: КР = 0,0203; К = 0,0334; Кс = -0,199. Длительность переходного процесса увеличивается, но полученный результат является робастным, то есть в меньшей степени устойчивость и качество переходного процесса зависит от точности реализации коэффициентов регулятора.
г+0-Ке>-1 |—
—i—hil>—кт——►
L-КЗ)-H derivative I-^ -—
-6343 5s*+3s+1
г _-= -
I 0.02 I-И [jaran"ieterUnknown |-H P h
I 0 035 I-И EarameterUnknow"n~l-КП-
I -0.2 I-[jaran"ieterUnknown |-H d h
■СИШ] -+CZEMMI2]
-+C
- 133354 I
кз-
-W abs
derivative
^ 20000>—1
_ r-H cost I
—КЖН
331 02 I
■ Р1о I 111®
1.2 1.0 .8 0 г = ® CL 1 4 .2 0 п
/
.......
.......
J.....
"0 20 40 80 80 100 Time [s)
Рис. 3. Структура с модификацией по Примеру 2
г-КЗ-Кв>-1
—UKD—км}——►
*-И derivative |—^
-6s+3
5s2+3s+1
□
^3—,
I 0 0212 I-Ц para meterUnkn own I-ЙП-
-►I paraineterllnknc.-Yn I-КП—
2 12e-2 I
I 0.0248 h ГоП-
2 43e-2 I
-H parameterUnknown I-KdT-
m
-КЗ"
~H abs
derivative
>1 20000>—1
_ r-Ц cost I
—КЖ1—1
1313 33 |
■ Plot НЕ
1.4 12 1 0 . 8 "О f 6 Ё < .4 .2 0 -.2 [
__________ А
.
I
If
J
) 20 40 60 80 100 Time (з)
Рис. 4. Структура с модификацией по Примеру 3 с запаздыванием от 1 до 5 с
функция которого дана соотношением:
Пример 2. Введем в объект из Примера 1 звено запаздывания. С регулятором, рассчитанным вторым способом в Примере 1, система остается устойчивой. Процедура оптимизации дает уточненный результат, показанный на рис. 3.
Робастность этого результата
подтверждается серией переходных процессов, полученных при различных значениях запаздывания от 1 до 5 с.
Пример 3. Введем в тот же объект запаздывание 8 с. Получаем результат, показанный на рис. 4. Робастность подтверждается для значений запаздывания 12 с, система остается устойчивой и при значениях до 20 с, однако перерегулирование в этом случае возрастает до величины 50 %.
Для значений запаздывания от 16 до 24 с также может быть рассчитан робастный регулятор, результаты этого расчета показаны на рис. 6.
Пример 4. Рассмотрим объект, передаточная
W^s) = 1 / (10 s3 + 4 s2 + 2 s + 1). (3)
Объект характеризуется сильной склонностью к колебаниям. Переходный процесс самого объекта без системы при воздействии ступенчатого скачка показан на рис. 7, это - колебания с возрастающей амплитудой.
В результате оптимизации получены следующие коэффициенты регулятора: КР = -0,59; К = 0,0429; Кс = -0,1537. Мы видим, что переходный процесс вначале идет в противоположную сторону, то есть ошибка вначале возрастает. Эта ситуация не желательна, однако, система устойчива и успешно завершает переходный процесс, сводя ошибку к нулю.
Рис. 5. Структура с модификацией по Примеру 3 с запаздыванием от 8 до 20 с
Рис. 6. Структура с модификацией по Примеру 3 с запаздыванием от 16 до 24 с
Рис. 7. Переходный процесс в объекте из Примера 4 без регулятора
Для оценки этого результата рассмотрим один из альтернативных методов решения этой задачи из Примера 4. Альтернативный метод состоит в введении в регулятор звена, которое компенсирует нежелательный полином знаменателя передаточной функции объекта введением такого же полинома в числитель звена, включаемого на входе объекта. При этом в знаменателе такого звена помещается полином, обеспечивающий более
благоприятный переходный процесс.
Для иллюстрации действия этого компенсирующего звена на рис. 9 показаны два переходных процесса: переходный процесс с
колебаниями соответствует передаточной функции объекта, а переходный процесс без колебаний - передаточной функции модели с благоприятным полиномом в знаменателе. На рис. 10 показан переходный процесс компенсирующего звена.
Полученный объект совместно с компенсатором можно рассматривать как новый объект, для которого можно осуществить численную оптимизацию регулятора. Структура для такой оптимизации и ее результат показаны на рис. 11.
Если бы моделирование осуществлялось только на интервале 90 с, можно было бы сделать ошибочный вывод, что результат оптимизации соответствует поставленной задаче. Моделирование на дальнейшем интервале показывает, что система неустойчива. После достижения нулевого значения ошибки выходной сигнал системы начинает колебаться около равновесного состояния с возрастающей амплитудой. При отличии только одного коэффициента компенсатора на 5% в любую сторону (коэффициент при первой степени б) колебания начинаются раньше и возрастают быстрее, как показано на рис. 12.
г*Н—Н1>—
-*П>—►TvsT-
-H derivative"
1
10s +4s'+2s+1
>
-+Ш—^
I -0 530242 I-И parameterlln known I-ИП-
I 0 0453732 I-И parameterlln known I-ИЗ—
I 0 15231= I-И parameterUn known I-HjLb
-+C
- 530063 I
-+ЕПМИШ
-+C
153742 I
+Ш-
~H abs
I derivative
>1 20000V-1
Г-Н cost I
32S.11 I
7*©—ИТ/Тр
g Plot ш
2.0 1.5 1.0 "О f 0 Е < ,5 10 15 -2.0 (
К. .Л
/
\ I
\J"
) 20 40 60 80 100 120| Time [s)
Рис. 8. Структура с модификацией по Примеру 4 и результат оптимизации регулятора
Ш-
Рис. 9. Структура и переходные процессы звеньев объекта и улучшенной модели
Юз +4з +29+1
Юз' -На 47з+1
т plo R
Amplitude
/
I
0 20 40 60 80 100 120l Time (s)
Рис. 10. Структура и переходный процесс компенсационного звена
Рис. 11. Структура с модификацией по Примеру 4 с компенсационным звеном и результат оптимизации регулятора
Таким образом, следует признать, что метод с компенсирующим элементом ненадежен, получаемый результат неудовлетворителен.
При моделировании этой же задачи с таким же решением в программе МЛТЬЛБ результат может оказаться положительным, то есть в случае полного совпадения полиномов знаменателя объекта и числителя компенсирующего элемента неустойчивое движение не возникнет. Однако практической ценности такой положительный результат не имеет, поскольку достижение точного равенства двух полиномов, один из которых известен лишь из экспериментальных оценок, следовательно, с ограниченной точностью, невозможно. И даже при полном совпадении этих полиномов их действие осуществляется в разных точках, поэтому результат будет соответствовать результату, получаемому в программе У^Бт, то есть даже точное совпадение моделей не даст устойчивого решения. Неустойчивость может проявиться в конце переходного процесса и при моделировании за малое время никак не проявиться.
Рис. 12. Переходные процессы в системе при отличии одного коэффициента компенсирующего элемента на 0,5 %
ЗАЛЮЧЕНИЕ
Таким образом, предложенная модификация целевых функций позволяет успешно вычислять коэффициенты регулятора методом численной оптимизации в программе У^Бт. Данная работа продолжает и пополняет результаты исследований авторского коллектива в данном направлении [3-5 ].
Работа выполнена по заданию Министерства образования и науки по проекту «Исследование предельных точностей оптических методов измерения параметров движения и мехатронных методов управления движением и разработка новых робототехнических и электромеханических систем», Темплан, проект № 7.559.2011, НИР № 01201255056.
ЛИТЕРАТУРА
[1] Ziegler, J.G and Nichols, N.B. Optimum Settings for Automatic Controllers. Transactions of the ASME 64, 1942, pp. 759-768.
[2] Zavorin A.N. et al. The modification of the quality characteristics of system of control with feedback by the using of PI2D2-regulators. Collection of science works of NSTU. 2010. 4(62). P.41 - 50. (In Russian: Заворин А.Н., Ядрышников О. Д., Жмудь В. А. Усовершенствование качественных характеристик систем управления с обратной связью при использовании ПИ2Д2 -регулятора. Сборник научных трудов НГТУ. Новосибирск. 2010. 4 (62). С.41-50).
[3] Voevoda A.A. et al. Comparative analysis of the optimization methods with the use of programs MATLAB and VisSim. Mechatronics, automation and control. 2012. 9. 37-43. (In Russian: ., А.Н. Заворин, О.Д. Ядрышников. Сравнительный анализ методов оптимизации регуляторов с использованием программных средств VisSim и MATLAB // Мехатроника, автоматизация и управление. № 9, 2012. с. 37-43)
[4] Poller B. V. et al. The design of robust regulator with the method of double iterative numerical optimization. Science Bulletin of NSTU. 2012. 2. P. 196200. (In Russian: Синтез робастного регулятора методом двойной итеративной параллельной численной оптимизации / Б. В. Поллер, В. А. Жмудь, С. П. Новицкий, А. Н. Заворин, О. Д. Ядрышников // Научный вестник НГТУ. - 2012. - № 2 . С 196 -200).
[5] Zhmud V.A. et al. The method of the designing of adaptive control systems for the controlling of non-stationary object with delay. Science Bulletin of NSTU. 2012. 3. (In Russian: Метод проектирования адаптивных систем для управления нестационарными объектами с запаздыванием. / В. А. Жмудь, А.Н. Заворин, Полищук А.В., О. Д. Ядрышников // Научный вестник НГТУ. 2012. №3).
Вадим Жмудь - заведующий кафедрой Автоматики в НГТУ, профессор, доктор технических наук, автор 200 научных статей., главный научный сотрудник Института лазерной физики СО РАН. Область научных интересов и компетенций -теория автоматического управления, электроника, измерительная техника.
E-mail: [email protected] Олег Ядрышников, аспирант кафедры Автоматики НГТУ, автор более 10 научных статей. Область научных интересов и компетенций - теория автоматического управления, оптимальные и адаптивные системы, оптимизация, многоканальные системы. E-mail: [email protected]