CC BY
82
УДК 519.8 ББК Ж 30
ПАРАМЕТРИЧЕСКИЙ СИНТЕЗ ОПТИМАЛЬНОГО В СРЕДНЕМ ДРОБНОГО ПИД-РЕГУЛЯТОРА В ЗАДАЧЕ УПРАВЛЕНИЯ ПОЛЕТОМ
Пантелеев А. В.1, Летова Т. А.2, Помазуева Е. А.3
(Московский авиационный институт (национальный исследовательский университет), Москва)
Поставлена задача поиска параметров дробного ПИД-регулятора, оптимального на множестве начальных состояний замкнутой системы «объект-регулятор» и множестве типовых входных воздействий, как задача поиска минимума функции многих переменных. Сформированы порядок решения задачи и соответствующее программное обеспечение. Приведен пример решения поставленной задачи поиска параметров оптимального в среднем дробного ПИД-регулятора для управления движением самолета по углу тангажа.
Ключевые слова: дробный ПИД-регулятор, критерий оптимизации, множество начальных состояний, множество входных воздействий, метод имитации отжига.
1. Введение
Использованию дробного интегро-дифференциального исчисления для описания динамики процессов управления посвящено большое число работ [3], в которых обсуждаются вопросы устойчивости, управляемости и качества дробных динамических систем. Эти работы условно можно разделить на два класса.
1 Андрей Владимирович Пантелеев, доктор физико-математических наук, профессор (avpanteleev@inbox.ru, тел. (499)158-48-11).
2 Татьяна Александровна Летова, кандидат технических наук, доцент (тел. (495)735-85-08).
3 Екатерина Александровна Помазуева, студент (kate-420@rambler.ru, тел. (499)492-88-88).
1. Работы, посвященные математическому моделированию и качественному анализу динамических систем, описываемых дробными дифференциальными уравнениями. Так, в [10, 11, 12] сформулированы и доказаны необходимые и достаточные условия устойчивости систем с дробно-степенной передаточной функцией (теорема Матиньона). В [18] проведен анализ поведения интегральной оценки (интеграл от модуля и квадрата ошибки) дробных систем управления, на вход которых подается воздействие в виде ступенчатой функции, и сформулированы необходимые и достаточные условия, налагаемые на параметры дробной передаточной функции, при которых может быть обеспечено экспоненциальное затухание ошибки. В [15] рассмотрены вопросы управляемости дробной линейной стационарной системы.
Поиску оптимального управления дробными динамиче -скими системами посвящены работы [ 3, 7, 14, 16]. В этих работах критерий оптимальности представляет собой функционал, зависящий от переменных состояния и управления, а объект управления описывается системой дробных дифференциальных уравнений. В работе [7] содержатся результаты моделирования дробных управляемых динамических систем и сформированы рекомендации по построению высокоэффективных численных схем.
2. Работы, посвященные формированию дробных ПИД-регуляторов различного типа [3, 8, 10, 13], используемых для управления динамическими системами и описываемыми дифференциальными уравнениями дробного и целого порядка. Работы [3, 10] содержат исследования по применению ТГО-контролле-ров для управления динамическими системами. Возможности применения СЯОКЕ-контроллеров для обеспечения высокой работоспособности по отношению к вариации коэффициентов передачи управления описаны в [3, 13].
Решение задачи о выборе параметров наиболее сложного Р1аВД-контроллера, имеющего передаточную функцию
ж =к+ту,
в которой предметом поиска являются параметры К, Д, а, Т-, Тй, описано в [9, 13, 17]. Использование методов многомерной оп-
тимизации для определения параметров Р/°©^-контроллера описано в [3]. В работах [3, 8, 13] содержится информация о применении дробных контроллеров для решения практических задач, связанных с управлением роботами, системами водоснабжения и электропитания.
В настоящей работе рассмотрены задачи поиска параметров дробного ПИД-регулятора с конечной памятью, оптимального в среднем на заданном множестве О начальных состояний системы и множестве О типовых входных воздействий.
2. Общая постановка задачи
Рассматривается замкнутая система «объект - дробный ПИД-регулятор», представленная на рис. 1, в которой:
Рис. 1.
• математическая модель объекта управления представляется уравнениями возмущенного движения
(1) x = f(x(t),u(t),t),x( 0) = хоеЦ
где х = (xi, ..., x„)T - вектор состояния; u = (ui, ..., Um)T, m < n - вектор управления; t - время, t e [0, T], момент T соответствует окончанию процесса; fx, u, t) - заданная непрерывная функция; xo - начальный вектор состояния; Q - множество начальных состояний;
• g(t) - входное воздействие, принадлежащее множеству G типовых входных воздействий;
• £(0 = g(t) - х(0 - ошибка.
Дробный ПИД-регулятор формирует управляющее воздействие н(0 на основании ошибки в виде
(2) и (*) = К£) + К а£) + КвВ(£(*), где
• Кпр, К/, Ки - коэффициенты регулятора;
• Л-д(г)£(?) - интеграл дробного порядка;
17, 0 < * <А,
А(*) = ^ - «память» регулятора; * — А(*) - нижний
предел интегрирования; t - верхний предел интегрирования; а > 0 - дробный порядок интеграла;
• Эо£(0 - дробная производная; ( > 0 - дробный порядок производной.
Интеграл дробного порядка может вычисляться:
а) по формуле Грюнвальда-Летникова [1], когда
(3) /^е=*а 1 (—1У Г(а+1) 4*—(/—А«и
'^ Ц( ) ГС/ + 1)Г(а — / +1) ^ ^у И ) /
где [х] - целая часть; И - шаг; Г(а) - гамма-функция;
б) по формуле Римана-Лиувилля [1]:
(4) /Гдт£(0 = [ £(Т) сТ, 0 <а< 1.
() *—А(* »() Г (а) *—А(* )(* — т)1—а ,
При численной реализации вычисление интеграла в (4) производится по формулам прямоугольников и трапеций. Дробная производная может быть вычислена: а) по формуле Грюнвальда-Летникова:
Г * 1
(5) D$e(t) = hP Y (-1)j-Г(Р +1-s(t - j • s(t - jh)),
0 r(j + 1)Г(Д-j +1)
где P > 0 - порядок дробной производной; б) по формуле Гёльдера [14]:
(ел пР^гл £(t) — s(t - h)
(6) Dps(t) =-—-.
Структура дробного ПИД-регулятора, в отличие от рассмотренной в [1], учитывает «память» системы, так как вычисление интеграла осуществляется не с начального момента времени до текущего значения t, а с момента t - Д(0 до t.
Будем в дальнейшем обозначать дробный ПИД-регулятор, реализующий закон (2), ПИяДд -регулятором [1].
Предположим, что в задаче стабилизации замкнутая система «объект - ПИяДд-регулятор» подвергается следующим воздействиям:
- возможные начальные состояния х(0) = xo, xo е О, где О -множество возможных начальных состояний, в котором выбирается p характерных начальных состояний xok, k = 1, ..., p;
- пробные типовые входные воздействия gl(t), ..., gN(t), образующие множество входных воздействий G, где N - заданное целое число.
Будем искать значения параметров ПИяДд-регулятора Kпр, К/, Ко, а, Д оптимизируя критерий качества управления J, который формируем по следующему правилу, содержащему три этапа.
1. Интегральную ошибку на интервале [0, 7] при заданном входном воздействии gj(t), 1 = 1, ..., N, и фиксированном начальном состоянии х(0) = хо будем оценивать величиной
о
2. Средняя интегральная ошибка J1 на множестве О возможных начальных состояний при заданном входном воздействии &{{), 1 = 1, ., N, равна
где mes Q - мера множества Q.
3. Среднюю интегральную ошибку J на множестве G всех N пробных типовых входных воздействий и множестве возможных начальных состояний найдем по формуле
T
I(x0) = \s2(t)dt, j = 1,...,N.
Jj = Q
, j = 1,...,N,
mes Q
j=i
(7) J(Knp, K, KD, a, ß) = N
Замечание.
1. В частном случае можно считать, что множество Q возможных начальных состояний задается параллелепипедом со сторонами, параллельными координатным осям. По каждой из координат задается равномерная сетка с некоторым шагом. В результате множество Q представляется объединением p непересекающихся элементарных подмножеств Qj, j = 1, ..., p. Координатами центра каждого из подмножеств Qj определяется начальное состояние xok, k = 1, ...,p. Тогда
. p
\Г ( x0)dx0 = ^Г ( xk )mes Q k,
Q k =1
где mes Qk - мера элементарного подмножества (объем).
В результате имеем
f^P'X )mes Qk £Г(х*)
ji ~ _0=1__ _0=1_
mes Q p
2. Значение критерия (7) характеризует поведение пучка траекторий, исходящего из множества начальных состояний Q, при входных воздействиях из множества G.
Требуется при заданном множестве начальных состояний Q и множестве G входных воздействий g(t) определить параметры Кпр, Ki, Kd, a, ß дробного ПИaДß -регулятора из условия минимума критерия J, т.е. среди всех возможных значений Кпр, Ki, Kd, a, ß (a > 0, ß > 0) требуется найти такие К*пр, K*i, K*d, a*, ß*, при которых критерий J принимает наименьшее значение:
J * = min J (Knp, K, KD ,a,ß).
Knp 'KI ,KD ,a,ß
Искомый дробный ПИдДд -регулятор с параметрами К*пр, КI, К d, a, ß будем называть оптимальным в среднем, поскольку при вычислении значения критерия J производится осреднение значений функционала I j(xo) по множествам Q и G.
Поставленная задача представляет собой задачу минимизации функции J(Kпр, К/, Кэ, а, ¡5), которая может быть решена с использованием одного из численных методов. Структура целевой функции 3(Кпр, Кг, Кэ, а, ¡5) и отсутствие информации о виде ее поверхностей уровня не позволяют использовать для эффективного решения задачи методы первого и второго порядка. Поэтому в работе используются два метода оптимизации нулевого порядка: метод адаптивного случайного поиска [4], позволяющий найти локальный минимум функции J с заданной точностью, и метод имитации отжига, относящийся к метаэвристиче-ским методам поиска глобального экстремума [6] и применяющийся для уточнения полученных первым методом результатов (при необходимости), а также для обеспечения возможности нахождения глобального минимума функции ^
3. Порядок решения задачи
1. Сформировать на множестве О множество характерных начальных состояний хок, к = 1, ..., р.
2. Сформировать N пробных типовых входных воздействий gj(t), ] = 1, ..., N (множество О).
3. Выбрать формулы для вычисления дробного интеграла и производной (таблица 1).
Таблица 1._
№ Формулы, используемые для вычисления дробного интеграла и дробной производной
1. - Интеграл вычисляется по формуле Грюнвальда-Летникова. - Производная вычисляется по формуле Грюнвальда-Летникова.
2. - Интеграл вычисляется по формуле Грюнвальда-Летникова. - Производная вычисляется по формуле Гёльдера.
№ Формулы, используемые для вычисления дробного интеграла и дробной производной
3. - Интеграл Римана-Лиувилля по формуле прямоугольников. - Производная вычисляется по формуле Грюнвальда-Летникова.
4. - Интеграл Римана-Лиувилля по формуле прямоугольников. - Производная вычисляется по формуле Гёльдера.
5. - Интеграл Римана-Лиувилля по формуле трапеций. - Производная вычисляется по формуле Грюнвальда-Летникова.
6. - Интеграл Римана-Лиувилля по формуле трапеций. - Производная вычисляется по формуле Гёльдера.
4. Задать пробные начальные значения К°пр, К0/, К°э, а0, 50 параметров ПИЯД( -регулятора, при которых замкнутая система «объект - ПИЯД(-регулятор» будет устойчива (рекомендуется), и:
а) проинтегрировать с использованием численных методов уравнения замкнутой системы «объект-регулятор», где X = /ШМОЛФ) = х0\ и( о = Кпрг(Г)+ВДД((Х(8(0 +вдр8(0
на интервале [0, Т] для Х0к, к = 1, р, у = 1, N;
б) вычислить значения интегральной ошибки / ](хко), у = 1, ..., N к = 1, ..., р;
в) вычислить значения средней интегральной ошибки
р
Т1' (4)
■Я = —-, ] = 1,-, N;
р
г) вычислить значение критерия
N
■0 = ^ = ■ К,К°,К°в,а°, () .
5. Задать начальные значения параметров метода адаптивного случайного поиска: коэффициент сжатия 0 < в < 1; ко-
эффициент растяжения ар > 1; М - максимальное число выполненных испытаний на текущей итерации; Но - начальную величину шага; Я - минимальную величину шага; V - максимальное количество итераций. Решить задачу поиска минимума ЛКПр, К.1, Ко, а, ¡3) методом адаптивного случайного поиска [4]. Для подсчета критерия./использовать операции а)-г) на шаге 4.
Результатом решения задачи будут параметры Кпр,К[,Кв,а, [3,
соответствующие (с заданной точностью) точке минимума критерия
6. Задать начальные значения глобального параметра То (температура), параметра закона распределения Больцмана с > 0, параметра у е [0,8; 0,99] и максимальное число итераций V. Решить задачу поиска глобального минимума критерия ^ методом имитации отжига [6], положив в качестве начальных значений искомых параметров значения Кпр,К1,КЕ>,а, [5 .
4. Пример
Рассмотрим в качестве примера задачу управления самолетом по углу тангажа в режиме горизонтального полета. Уравнения возмущенного движения самолета в режиме горизонтального полета имеют следующий вид [2]: Г( р + )а - р3 = 0,
(8) Г 2
К« Р + «32 )а + (р + пъър)& = -пъ8ъ, где а - отклонение угла атаки; 3 - отклонение угла тангажа; 5ь -отклонение руля высоты; р = й / - оператор дифференцирования.
Будем в дальнейшем записывать систему (8) в нормальной форме Коши, обозначив ^ — а, %2
= 3,хъ=3. Тогда имеем
— ^22*^1 Х3'
(9) I х2 =х3,
.Х^ — ^32 (п0 ^ ^33 )хз ^^Ъ^Ъ '
Управление самолетом по углу тангажа (3 = Х2) осуществляется ПИяДр-регулятором, формирующим отклонение руля высоты Зь(0 по закону
(10) ^ (0 = ^Р^д " + Кр2 Х3(?) + К, СА(, х, (Х^д - Х2(?)) +
+К^х^),
где х]2зад = (у - 1) / N - 1) - желаемое значение угла тангажа, у = 1, ..., N здесь учитывается, что Хззад = 0 - желаемое значение скорости изменения угла тангажа; вектор начальных состояний х0 = (Х10^, Х20к, Х30^)Г, к = 1, р.
Значения х]2зад = (у - 1) / ^ - 1) и х]ззад = 0 определяют постоянные входные воздействия на замкнутую систему:
g2(t) = х2(0 = Х2зад , gз(t) = х3(г) = хЗзад = 0, ] = 1,..., N. Критерии I/(х0к), J/ J имеют вид
V (х^,хко,х3о) = |(х21зад - Х2(0)2Л,] = 1,...,N;
^1 (х10, х20, х30)
3 = ^-, ] = 1,...,N;
Р
N
3 =
N
При формировании критерия J1 предполагалось, что множество О, возможных начальных состояний по координатам х1,х2,х3 (а,3, 3) задается параллелепипедом с заданными сторонами, который содержит в себе р элементарных параллелепипедов. За точку (хюк, х20к, хз0к)Т берется центральная внутренняя точка элементарного параллелепипеда.
Для решения задачи поиска минимума функции J(Kпр\, Кпр2, К Кэ, а, Р) на множестве О начальных состояний х10, х20, хз0 и множестве входных воздействий хУ2 зад(0, У = 1, ..., N разработана программа на языке С#, имеющая модульную структуру и позволяющая осуществлять решение задачи для различных способов вычисления дробного интеграла и дробной производной (таблица 1).
0
Работа с программой начинается с задания коэффициентов уравнений возмущенного движения (9), которые могут быть либо введены вручную, либо перенесены автоматически из предварительно заполненной таблицы (таблица 2) [2] в соответствии с типом самолета (легкий, тяжелый) в ходе ответа на вопрос «Выберите тип самолета» (рис. 2).
Таблица 2.
Коэффициенты Легкий самолет, H = 15 км Тяжелый самолет, H = 8 км
П0 0,7 1,17
n "22 2,5 3
n 32 16 42
n "33 2,2 2,5
nb 100 28
Задаются границы параллелепипеда, определяющие множество Q возможных начальных состояний, начальные значения Х1, Х2, Х3 (угла атаки, угла тангажа и скорости изменения угла тангажа), начальное входное воздействие g(t) е G, G = {g(t): g(t) = Х2 3ad(t) = р- l(t), 0 < р< 1}, значение А.
Рис. 2.
Затем следует выбрать формулы (рис. 3) для вычисления дробного интеграла и дробной производной в соответствии с таблицей 1 и перейти к заданию начальных значений параметров ПИдДд-регулятора (рис. 2).
Рис. 3.
Начальные значения искомых параметров ПИаДд-регулятора задаются на основании предварительных испытаний замкнутой системы «объект - ПИаДд-регулятор» на устойчивость, которые можно провести, интегрируя уравнения (9), (10) (рис. 2). Интегрирование уравнений системы стабилизации (9), (10) осуществляется методом Рунге-Кутты 4-го порядка с шагом к = 0,01. Результат интегрирования системы приведен на рис. 4.
Интегрирование системы (9), (10) осуществляется многократно с учетом множества О начальных состояний системы и множества О типовых входных воздействий g(t) при р = 8, N = 6. В процессе интегрирования вычисляются критерии Р(хок), З, З. Начальное значение критерия З1 = 0,1634, а также переходные процессы для заданного начального состояния из множества О и выбранного входного воздействия Х2 зад = 1 из множества О, приведены на рис. 4.
00 vo
Тип самолета - Легкий самолет. Н = 15 Начальное состояние Угол атаки = 0.1 Угол тангажа — 0.05
Скорость изменения утла тангажа = 0.05 Входное воздействие = 1
Начальные параметры управления
К ,= 1
к = л kd= о
а =0.9
J = 0.1634340650527
Для решения задачи оптимизации:
J -> min
задайте параметры метода адаптивного случайного поиска
Посмотреть описание параметров адаптивного случайного поиска
к= 2 е = 0.8
м= 10 ар 12
R = 0.01 : - 50
Решить задачу оптимизации методом адагтшБного случайного поиска
Результат интегрирования системы уравнений при заданном начальном состоянии
0.1 0.2 о.з
0.5 t
- Угол атаки
- Угол тангажа
■ Нижняя граница трубки
■ Верхняя граница трубки
- Входной сигнал
Рис. 4.
Тип самолета - Легкий. Н = 15
Начальное состояние
Угол атаки = 0.1
Угол тангажа — 0.05
Скорость изменения угла тангажа = 0.05
Входное воздействие = 1
Назальные параметры управления
Кщ,2 = 1
К= -1
к = о
а = 0.9 0= 1
I = 0.1634340650527
Параметры управления после решения методом случайного
Кщ>1= -3.40853169436983
К^2= 0.636372253841674
Кг= -2.79606565744435
Кп= -0.432487272290073
а = 3.66911531645308 0= 0.0382510161546668 I = 0.0189668965976935 Задайте параметры метода случайного поиска
200 0.8
Перед методом случайного поиска
■ После метода случайного поиска
■ После метода имитации отжига
■ Входной сигнал
■ Нижняя граница трубки
■ Верхняя граница трубки
Решить задачу оптимизации методом имитации отжига
Параметры управления после решения, методом имитации отжига
Кф! = -19.1045898262898 - 0.127898046387717 к = 3.94894267906112
I
Кв= 0.198409485304756
а - 9.20617285349323 Р = 0.13027928168875 J = 0.00514651465586379
Строка состояния
Рис. 5.
Решение задачи оптимизации осуществляется в два этапа.
На первом этапе используется метод адаптивного случайного поиска. Для его реализации следует задать параметры метода (рис. 4). Решение задачи оптимизации позволяет определить значения параметров ПИаДр -регулятора, соответствующие с заданной точностью точке минимума критерия J, и построить график переходного процесса по регулируемой координате -углу тангажа (рис. 5).
Такое построение можно осуществлять для любого (из числа заданных) возможного начального состояния и входного воздействия g(t) из заданного множества О. Если качество переходного процесса будет удовлетворять предъявленным техническим требованиям, то дальнейший поиск можно прекратить.
На втором этапе осуществляется поиск глобального минимума методом имитации отжига. Для его реализации следует задать параметры метода (рис. 5). Результатом решения будут значения параметров ПИаДр-регулятора Кр1, Кр2, К Ки, а, р и минимальное значение У* функции J. На рис. 5 представлен переходный процесс по углу тангажа, соответствующий этим параметрам.
На базе разработанной программы были проведены расчеты оптимальных параметров ПИаДр -регулятора с использованием всех приведенных (таблица 1) формул для вычисления дробного интеграла и дробной производной для двух типов самолетов (легкого и тяжелого) при различных значениях А. Настройки ПИаДр-регулятора и средние по пучку траекторий показатели качества переходного процесса по углу тангажа для легкого и тяжелого самолета (ст- перерегулирование и tp - время переходного процесса) приведены в таблицах 3-10 для всех шести наборов формул из таблицы 1. Настройки ПИаДр-регулятора получены в результате решения задачи оптимизации методом имитации отжига, а формулы для вычисления дробного интеграла Римана-Лиувилля методом прямоугольников и трапеций соответственно имеют вид
г
Таблица 3. Настройки ПИа Дв -регулятора для легкого самолета при Д = 0,1
№ Кпр1 Кпр 2 КI а Р
1 -19,104 0,127 0,198 3,948 9,206 0,130
2 -11,370 1,172 -0,940 -19,084 0,181 0,086
3 -22,074 -0,117 1 -28,878 0,623 0,020
4 -61,687 1,513 -0,075 -39,706 0,134 0,253
5 -11,368 0,583 -0,074 -5,655 0,292 0,328
6 -29,987 0,557 -0,017 14,450 0,401 0,766
Таблица 4. Настройки ПИа Дв -регулятора для легкого самолета при Д = 0,2
№ К , пр1 Кпр 2 К а Р
1 -9,02 0,838 -0,301 -32,107 3,835 0,221
2 -12,41 0,454 -0,246 -22,941 13,983 0,163
3 -7,22 0,081 0,047 -7,838 0,814 0,351
4 -19,18 0,509 0,047 -4,759 0,937 0,131
5 -6,69 0,048 0,314 -0,087 0,401 0,108
6 -17,36 2,760 -22,320 -14,079 0,404 0,298
I а5(г) =
£ 1
Г(а) Гг-дот ](п - ] )1
0 < а < 1,
Ге(,) = т±1Ш?кЛ < а < 1, п =
2Г(а) , Г 1 (п - 7 )
а
Таблица 5. Показатели качества переходного процесса по углу тангажа для легкого самолета при А = 0,1
Таблица 6. Показатели качества переходного процесса по углу тангажа для легкого самолета при А = 0,2
№ о, % *Р, с 3
1 10,6 0,09 0,0051
2 1,38 0,14 0,0083
3 5,1 0,08 0,0056
4 4,1 0,04 0,0026
5 6,7 0,08 0,0079
6 10,9 0,77 0,0034
№ о, % *р, с 3
1 6 0,16 0,0099
2 0 0,06 0,0066
3 21,8 0,28 0,0143
4 7,4 0,08 0,0057
5 0,8 0,12 0,0144
6 0 0,33 0,0132
Таблица 7. Настройки ПИа Дв -регулятора для тяжелого самолета при А = 0,1
№ Кпр1 Кпр 2 К а Р
1 -60,49 2,075 -0,679 9,956 64,360 0,220
2 -25,07 1,477 1,314 55,257 16,573 0,189
3 -21,64 3,850 -3,596 41,665 0,106 0,013
4 -15,39 0,755 -0,486 -35,079 0,584 0,390
5 -162,28 8,773 -4,122 -43,444 0,758 0,016
6 -18,66 0,626 -0,969 -2,044 0,444 0,028
Таблица 8. Настройки ПИа Дв -регулятора для тяжелого самолета при А = 0,2
№ Кпр1 Кпр 2 К: а Р
1 -53,68 -1,013 1,960 21,003 30,790 0,108
2 -26,66 1,429 0,029 -75,245 37,492 0,982
3 -10,36 -8,926 8,929 1,055 0,217 0,008
4 -8,77 0,266 -1,179 -0,700 0,227 0,151
5 -4,82 0,663 -0,102 -3,214 0,493 0,614
6 -68,10 5,509 -0,963 -15,304 0,461 0,339
Таблица 9. Показатели качества переходного процесса по углу тангажа для тяжелого самолета при А = 0,1
Таблица 10. Показатели качества переходного процесса по углу тангажа для тяжелого самолета при А = 0,2
№ о, % ¿р, с 3
1 21,2 0,21 0,0081
2 3 0,11 0,0165
3 10,7 0,24 0,0153
4 11,7 0,13 0,0103
5 2,6 0,06 0,0053
6 8 0,15 0,0114
№ о, % ¿р, с 3
1 6,8 0,09 0,0062
2 3,7 0,11 0,0151
3 5,6 0,19 0,0092
4 11,8 0,29 0,0156
5 0 0,35 0,0258
6 0 0,15 0,0120
На рис. 6, 7 приведены переходные процессы по углу тангажа для легкого и тяжелого самолетов при А = 0,1, когда производная и интеграл вычисляются по формуле Грюнвальда-Летникова.
На рис. 8, 9 приведены переходные процессы по углу тангажа для легкого и тяжелого самолетов при А = 0,1, когда для вычисления интеграла Римана-Лиувилля используется формула трапеций, а дробная производная вычисляется по формуле Грюнвальда-Летникова.
Легко видеть, что использование интеграла Римана-Лиувилля при одной и той же формуле вычисления дробной производной позволяет значительно уменьшить перерегулирование. В [5] рассмотрено решение задачи параметрического синтеза обобщенного ПИД-регулятора, формирующего сигнал управления на основании ошибки £(У) в виде
г
иЦ) = КпреЦ) + КвёЦ) + Кс,ёЦ) + К[ | е(т)^т.
г-Д(г)
Структура обобщенного регулятора включает дополнительное слагаемое, учитывающее изменение второй производной ошибки, а в последнем слагаемом предлагается находить значение интеграла от ошибки не на всем промежутке времени от момента ¿0 начала процесса регулирования до текущего момента времени, а лишь на части этого промежутка, характеризующей
«память» системы. Описанные изменения в структуре регулятора по сравнению с классической могут расширить его возможности по обеспечению требуемого качества переходных процессов.
При решении задачи параметрического синтеза обобщенного ПИД-регулятора для легкого и тяжелого самолетов [2] были получены наименьшие значения критерия 3. для легкого самолета при А = 0,1 3 = 0,0102; для тяжелого самолета 3 = 0,0159 при А = 0,1.
Сравнение результатов моделирования для обоих типов самолетов при использовании обобщенного ПИД-регулятора [5] и ПИдДр -регулятора при одном и том же значении А (таблицы 5 и 9) говорит о том, что использование ПИаДр-регулятора позволяет снизить значение критерия 3 для легкого самолета в 3,92 раза, а для тяжелого - в 3 раза.
Время решения задачи оптимизации на этапе использования адаптивного случайного поиска составляет 30-45 с при числе итераций V = 50. Время решения задачи оптимизации на этапе использования метода имитации отжига составляет от 3 до 5 минут при числе итераций не менее 800.
1.2-1 3
Рис. 6. Переходные процессы по углу тангажа для легкого самолета при вычислении дробной производной и дробного интеграла по формулам из п. 1 таблицы 1, А = 0,1: 1 - до начала оптимизации параметров; 2 - после оптимизации параметров методом адаптивного случайного поиска; 3 - после применения метода имитации отжига
Рис. 7. Переходные процессы по углу тангажа для тяжелого самолета при вычислении дробной производной и дробного интеграла по формулам из п. 1 таблицы 1, А = 0,1: 1 - до начала оптимизации параметров; 2 - после оптимизации параметров методом адаптивного случайного поиска; 3 - после применения метода имитации отжига
1.2—| __3
Рис. 8. Переходные процессы по углу тангажа для легкого самолета при вычислении дробной производной и дробного интеграла по формулам из п. 5 таблицы 1, А = 0,1: 1 - до начала оптимизации параметров; 2 - после оптимизации параметров методом адаптивного случайного поиска; 3 - после применения метода имитации отжига
Рис. 9. Переходные процессы по углу тангажа для тяжелого самолета при вычислении дробной производной и дробного интеграла по формулам из п. 5 таблицы 1, А = 0,1: 1 - до начала оптимизации параметров; 2 - после оптимизации параметров методом адаптивного случайного поиска; 3 - после применения метода имитации отжига
5. Выводы
1. Поставлена задача о параметрическом синтезе оптимального в среднем дробного ПИаДр-регулятора на заданном множестве О начальных состояний системы и множестве О возможных входных воздействий g(t) с учетом памяти регулятора
А(0.
2. Поставленная задача сформулирована как задача минимизации функции многих переменных Т(Кпр, К, Кэ, а, Р).
3. Предложен порядок решения задачи, позволяющий найти минимум критерия ^ с последовательным использованием методов адаптивного случайного поиска и метода имитации отжига.
4. Разработано программное обеспечение, поддерживающее предложенный порядок решения задачи.
5. Приведен пример решения задачи поиска параметров оптимального в среднем дробного ПИаДр -регулятора для управления самолетом по углу тангажа в режиме горизонтального полета.
6. Настройки дробного ПИаДр -регулятора в задаче
управления самолетом по углу тангажа позволяют существенно
снизить оптимальное значение критерия J по сравнению с его
значениями в случае синтеза обобщенного ПИД-регулятора [5].
Литература
1. АВСИЕВИЧ А.В., АВСИЕВИЧ В В. Моделирование систем автоматического управления с дробным ПИД-регулятором // Вестник Самарского государственного технического университета. Серия «Техническое науки». -2010. - №1(26). - С. 6-59.
2. БОДНЕР В.А. Теория автоматического управления полетом. - М.: Наука, 1961. - 698 c.
3. БУТКОВСКИЙ А.Г., ПОСТНОВ С.С., ПОСТНОВА Е.А. Дробное интегро-дифференциальное исчисление и его приложения в теории управления. II. Дробные динамические системы: моделирование и аппаратная реализация // Автоматика и телемеханика. - 2013. - №5. - С. 3-34.
4. ПАНТЕЛЕЕВ А.В., ЛЕТОВА Т.А. Методы оптимизации в примерах и задачах. - М.: Высшая школа, 2GG8. - 544 с.
5. ПАНТЕЛЕЕВ А.В., ЛЕТОВА ТА., ПОМАЗУЕВА Е.А. Применение методов глобальной оптимизации для параметрического синтеза обобщенного пропорционально-интегрально-дифференциального регулятора в задаче управления полетом // Труды МАИ. - 2015. - №79. - 23 с.
6. ПАНТЕЛЕЕВ А.В., МЕТЛИЦКАЯ Д.В., АЛЕШИНА Е.А. Методы глобальной оптимизации. Метаэвристические стратегии и алгоритмы. - М.: Вузовская книга, 2013. -244 c.
7. CAI X., LIU F. Numerical Simulation of the Fractional-Order Control System // J. Appl. Math. Comp. - 2007. - Vol. 23, No. 1-2. - P. 229-241.
8. CAPONETTO R., DONGOLA G., FORTUNA L.. et al. Fractional Order Systems. Modelling and Control Applications. -Singapore: World Scintific, 2G1G. - 2GG p.
9. CASTILLO F.J., FELIU V., RIVAS R., et al. Desing of a Class of Fractional Controlles from Fraquency Spetifications with
Guaranteed Time Domane Behavior // Comput. Math. Appr. -2010. - Vol. 50. - P. 1656-1666.
10. CHEN Y.Q., PETRAS I., XUE D. Fractional Order Control - A Tutorial // Proc. 1009 Amer. Control Conf. St. Louis, 2009. -P. 1397-1411.
11. HARLEY T.T., LORENZO C.F. Dynamics and Control of Initialized Fractional-Order Systems // Nonlin. Dyn. - 2002. -Vol. 29, No. 1-4. - P. 201-233.
12. MATINGNON D. Stability Properties for Generalized Fractional Differential Systems // Proc. Colloq. Fractional Differential Systems: Models, Methods and Applications. - Vol. 5. -Paris, 1998. - P. 145-158.
13. MONJE C.A., CHEN Y.Q., VINAGRE B.M., et al. Fractional-Order Systems and Controls: Fundamentals and Applications. -London: Springer-Verlag, 2010. - 409 p.
14. MOPHOU G.M. Optimal Control of Fractional Diffusion Equation // Comput. Math. Appl. - 2011. - Vol. 61. - P. 68-78.
15. MOZYRSKA D., TORRES D.F.M. Modified Optimal Energy and Initial Memory of Fractional Continuous-Time Linear Systems // Signal Proc. - 2011. - Vol. 91. - P. 379-385.
16. ÖZDEMIR N., KARADENIZ D., ISKENDER B.B. Fractional Optimal Control Problem of a Distributed System in Cylindrical Coordinates // Phis. Lett. A. - 2009. - Vol. 373. - P. 221-226.
17. PADULA F., VISIOLI A. Tuning Rules for Optimal PID and Fractional-Oder PID Controlles // J. Proc. Contr. - 2011. -Vol. 21. - P. 69-81.
18. TAVAZOEI M.S. Notes on Integral Performance Indices in Fractional-Order Control Systems // J. Proc. Contr. - 2010. -Vol. 20. - P. 285-291.
PARAMETRIC SYNTHESIS OF AVERAGE OPTIMAL FRACTIONAL PID-CONTROLLER FOR FLIGHT CONTROL PROBLEM
Andrey Panteleev, Moscow Aviation Institute, Moscow, Doctor of Physical and Mathematical Sciences, professor (avpantele-ev@inbox.ru, (499)158-48-11).
Tatiana Letova, Moscow Aviation Institute, Moscow, Candidate of Engineering Sciences, associate professor ((495)735-85-08). Ekaterina Pomazueva, Moscow Aviation Institute, Moscow, student (kate-420@rambler.ru, (499)492-88-88).
Abstract: We consider a problem of fractional PID-controller tuning to optimize it in average over a set of initial states of a closed system object-controller and on the set of typical input actions. The problem is reduced to multivariate optimization problem. We suggest an approach to the solution and implement it algorithmically. The approach is illustrated with finding the average optimal parameters of a fractional PID controller to pitch attitude guidance of an aircraft.
Keywords: fractional PID-controller, optimization criterion, set of initial states, set of input signals, simulated annealing.
Статья представлена к публикации членом редакционной коллегии А.Г. Кушнером.
Поступила в редакцию 24.02.2015.
Опубликована 31.07.2015.
