CC BY
13
ЭЛЕКТРОННЫЙ НАУЧНЫЙ ЖУРНАЛ «APRIORI. CЕРИЯ: ЕСТЕСТВЕННЫЕ И ТЕХНИЧЕСКИЕ НАУКИ»
№ 1 2016
УДК 519.71
ОПРЕДЕЛЕНИЕ СООТНОШЕНИЯ ОПТИМАЛЬНОСТИ В ЛИНЕЙНО-КВАДРАТИЧНОЙ ЗАДАЧЕ УПРАВЛЕНИЯ С ОБОБЩЕННЫМ ФУНКЦИОНАЛОМ
Афонин Виктор Васильевич
кандидат технических наук Мордовский государственный университет им. Н.П. Огарёва, Саранск
Аннотация. В данной статье рассматриваются особенности программной реализации в системе MATLAB соотношений оптимальности в линейно-квадратичной задаче управления с обобщенным функционалом. Соотношения оптимальности связывают между собой параметры объекта (матрицы А, В), функционала (матрицы Q, R, N и П-регулятора
Ключевые слова: MATLAB; квадратичный функционал; матрицы объекта управления; весовые матрицы функционала; оптимальность.
DETERMINING THE RATIO OF OPTIMALITY IN LINEAR-QUADRATIC CONTROL PROBLEM WITH GENERALIZED FUNCTIONALITY
Afonin Victor Vasilievich
candidate of technical sciences Ogaryov Mordovia State University, Saransk
Abstarct. This article discusses the features of a software implementation in MATLAB relations optimality in the linear-quadratic control problem with the generalized functionality. Optimal ratio between an associated object parameters (matrix A, B), functional (matrix Q, R, N) and the P controller (K).
Key words: MATLAB; quadratic functional; matrices of the control object; the weight matrix functionality; optimality.
В данной статье рассматривается задача оптимальной стабилизации или задача аналитического конструирования оптимальных регуляторов (АКОР), для которой рассматриваются соотношения оптимальности [1-3] в плане особенностей их программной реализации для случая обобщенного функционала. Вывод соотношений оптимальности в случае традиционного квадратичного функционала дан в ранней работе автора [2]. В работе [3] рассмотрены некоторые теоремы относительно соотношений оптимальности. В [1; 2] рассмотрены различные виды традиционных функционалов, в которых задаются две весовые матрицы. Обобщенный квадратичный функционал содержит три весовые матрицы, традиционно их обозначают Q, R, N. В качестве среды и языка программирования выбран MATLAB R2014b, так как он предоставляет собой удобный и гибкий функционал для решения задач самого различного характера. Для поставленной задачи необходимо дополнение к ядру системы MATLAB, а именно пакет расширения (тулбокс) Control System Toolbox.
Особенности программирования можно понять из постановки задачи и полученных результатов. Опираясь на [1-3], приведем необходимые уравнения и собственно задачу управления, чтобы затем отметить особенности программной реализации полученных соотношений оптимальности.
Рассмотрим вполне управляемую модель объекта в виде
dX(t)
= AX(t) + BU(t), (1)
где X(t) - п-мерный вектор состояния, U(t) - /--мерный вектор управления, А - матрица состояния, матрица действительных чисел размера п*п, В - матрица входа, матрица действительных чисел размера п*г. Предполагается, что на управления uj(t) (j = 1,2,...,r) ограничения не
наложены. Управление должно быть выбрано так, чтобы при произволь-
ном начальном условии Х(0) минимизировать квадратичный функционал (обобщенный функционал) вида
00
(2)
О
где О - положительно-определенная симметрическая матрица действительных чисел размера пхп, [ - положительно-определенная симметрическая матрица действительных чисел размера гхг, Т - символ транспонирования, N - матрица действительных чисел, которая должна удовлетворять условию положительной определенности следующей коагулированной матрицы W:
Задача (1), (2) - это задача оптимальной стабилизации [4]. Рассмотрение задач оптимальной стабилизации с функционалом (2) -обобщенным функционалом - можно найти в [5-7]. Решение задачи (1), (2) имеет вид
где Р - (пхп) симметрическая матрица, положительно-определенное решение нелинейного матричного алгебраического уравнения Риккати
и (г) = + ЩХ(х),
(4)
АТР + РА- (РВ + Ы)Я-1(ВТР + Ыт) + (} = 0.
(5)
В(4) обозначим
К = /Г1^7/5 + ЛГ).
(6)
С учетом введенного обозначения (6) перепишем выражение (4) в виде
и = -КХ(0, (7)
где К - постоянная матрица оптимального регулятора размера гхп.
Предполагается, что все переменные состояния доступны для измерения, а система (1) полностью управляема по Калману, матрица вида (3) положительно-определенная.
Объект (1), замкнутый на оптимальное управление (7), имеет следующее описание:
=(А- ВК)Х{0. (8)
В (7) собственные числа матрицы А - В К будут располагаться в левой полуплоскости комплексной плоскости корней.
Объект управления (7) является устойчивым [7]. Если к объекту (7) приложить внешнее воздействие, которое обозначим П ( £) , то получим описание системы управления в виде
сIX (г)
= (А - В К) X (О + В П (0 . (9 )
Получим аналитическое соотношение оптимальности, подобное тому, которое приводится в работах [1, 2, 3].
Теорема. Для задачи (1), (2) с решением (4), (6), (7) и ограниченным внешним постоянным управляющим воздействием справедливо следующее соотношение:
Я = ВТБТ(25В + (Ег + ВТБТКТ)Я(ЕГ + КБВ) - 2ВТ5ТЫ(ЕГ + КБВ), (10)
где (2, Д, Ы- весовые матрицы обобщенного функционала (2), 5 = (л - вку1,
К- оптимальный П-регулятор, определяемый соотношением (6), £г - единичная матрица размера г хг, А, В- матрицы состояния и входа объекта управления (1), т- символ транспонирования.
Доказательство. В силу того, что оптимальная система (8) устойчивая, то при ограниченном воздействии на нее решение системы (9) будет стремиться к установившемуся значению. Следовательно, на достаточно большом отрезке времени (^ да) изменение вектора состояния Х(£) будет равным нулю, т.е.
ахи)
= 0, t->oo. (11)
(И
С учетом (11) уравнение (9) становится алгебраическим. Разрешим его относительно вектора состояния х(£):
0 = (Л -5 / X ( 0 + 5 и ( 0 ; (Л - 5 / X (О = - 5 и (0 ; (Л - 5 / - 1 (Л - 5 / X (О = - (Л - 5 / - х5 и (0;
X (О = - (л - 5 / - ^ и (0 . В последнем соотношении обратную матрицу обозначим через 5, т.е.
5 = (л - вку1. (12)
С учетом (12) будем иметь следующее выражение относительно вектора состояния:
X (г) = - бв П (0 . (1 з )
Соответственно, транспонируя вектор состояния (13), получим
хт (г) = - Пт (О В тБт. (14)
Уравнение (9) может быть получено, если в (1) управление представить в виде
и (0 = - КХ (0 + П (0. (1 5 )
Соответственно, транспонируя выражение (15), получим
П П
Поскольку интегральный функционал (2) задается в пределах от нуля до бесконечности, то в него подставим значения векторов (13)-(16). При этом выполним над интегрантом функционала следующие эквивалентные преобразования, опуская для краткости зависимости вектора состояния и управления от аргумента - времени:
Хт((Х + ПтЯ П + 2 хты и = = ПтВтБт ((БВ П + (Пт + ПтВтБтКт) я (П + кбв П) - 2 ПтВтБтЫ (П + кбв П).
Вынесем за скобки П П (соответственно, влево и вправо), получим следующее выражение интегранта:
UT[BTSTQSB + (Er + ВTSTKT)R(Er + KSB) - 2 ВTSTN(Er + KSB) ]U. ( 1 7)
Поскольку в заданном функционале (2) между ит( г) и и(г) находится весовая матрица я, то с учетом (17), можно принять
Я = 5Т5Т(255 + (Ег + ВТБТКТ)Я(ЕГ + КБВ) - 2ВТБТЫ(ЕГ + КБВ). (18)
Тем самым, полученное соотношение (18) завершает доказательство теоремы.
Очевидно, выполнение (18) связано с численными расчетами обратных матриц и алгебраическими операциями над матрицами. Все это может привести к определенным погрешностям.
Рассмотрим численные расчеты, выполненные в системе MATLAB. Для этого сначала приведем небольшую программу, в которой легко выполнить изменения размерности задачи аналитического конструирования регуляторов. Текст программы приведен в листинге 1.
% Листинг 1. Script - скрипт-файл или М-сценарий clear,clc
format short %% long %% формат вывода числовых значений n = 8; %% размерность вектора состояния r = 6; %% размерность вектора управления rng(1119); %% для воспроизводимости результатов
%% Задание матриц объекта случайно по нормальному и равномерному законам A = randn(n); %% матрица состояния объекта управления B = rand(length(A),r); %% матрица входа if rank(ctrb(A,B)) < n
еппоп('Объект не полностью управляем по состоянию. Расчет прерывается'); end
fprintf('\n Размерность матриц объекта управления.\n')
fprintf(' A:%gx%g, B:%gx%g\n', size(A,1),size(A,2),size(B,1),size(B,2))
% Определение симметрических положительно-определенных матриц % из галереи библиотечных матриц MATLAB Q = 1.23*galleny('toeppd',n); R = galleny('lehmen',n)*6.78;
nng(1062) %% для воспроизводимости результатов и выбора матрицы N
%% Задание матрицы N по равномерному закону из интервала от 0 до 1
N = nand(length(Q), length(R));
W = [Q, N; N', R]; if eq(W, W')
disp(' Матрица W симметрическая') else
disp(' Матрица W не симметрическая')
end
fprintf('\n Размерность матриц функционала.\n')
fprintf(' Q: %gx%g; R: %gx%g; N: %gx%g \n',n,n,r,r, size(N,l),size(N,2))
fprintf('\n R = \n')
disp(R)
disp('============================================')
K = lqn(A,B,Q,R, N); % расчет оптимального регулятора S = inv(A - B*K); Sn = A - B*K;
En = eye(size(R)); % единичная матрица % Прямая кодировка формулы (18)
% R2=B'*S'*Q*S*B+(Er+B'*S'*K')*R*(Er+K*S*B)-2*B'*S'*N*(Er+K*S*B); % Кодировка с учетом интеллектуальной подсказки компилятора MATLAB R2=B'*S'*Q*(Sn\B)+(En+B'*S'*K')*R*(En+K*(Sn\B))-2*B'*S'*N*(En+K*(Sn\B)); fprintf('\n R2 = \n') disp(R2)
stn='\n Максимальная абсолютная ошибка (R - R2): %g\n\n'; fpnintf(stn, max(max(abs(R-R2))))
Расчетное выражение (18) в программном коде обозначено как R2. Следует отметить, что выбор матрицы N функционала (2) не однозначен. С матрицей N коагулированная матрица W должна быть положительной полуопределенной, т. е. все ее миноры главной диагонали не
должны быть отрицательными. Соответственно, собственные числа матрицы (3) также не должны быть отрицательными. В противном случае при расчете оптимального регулятора с помощью библиотечной функции lqr выводится сообщение:
Warning: The [Q N;N' R] matrix should be positive semi-definite. Сначала приведем результат выполнения рассмотренной программы: A:8x8, B:8x6 Размерность матриц функционала. Q: 8x8; R: 6x6; N: 8x6 Матрица Mn симметрическая R =
6 .7800 3 .3900 2 2600 1 .6950 1 3560 1 .1300
3 .3900 6 .7800 4 .5200 3 .3900 2 7120 2 2600
2 .2600 4 .5200 6 7800 5 .0850 4 0680 3 3900
1 .6950 3 .3900 5 0850 6 .7800 5 4240 4 .5200
1 .3560 2 .7120 4 0680 5 .4240 6 7800 5 6500
1 .1300 2 .2600 3 3900 4 .5200 5 6500 6 7800
Wanning: The [Q N;N' R] matrix should be positive semi-definite. Type "help lqr" for more information. R2 =
6 .7800 3 .5271 2 4101 1 .8163 1 2775 1 0008
3 .2529 6 .7800 4 .4415 3 .4422 2 7308 2 2218
2 .1099 4 .5985 6 7800 4 .9405 4 0868 3 3321
1 .5737 3 .3378 5 .2295 6 .7800 5 3155 4 5171
1 .4345 2 .6932 4 .0492 5 .5325 6 7800 5 5767
1 .2592 2 .2982 3 4479 4 .5229 5 7233 6 7800
Максимальная абсолютная ошибка ^ - К2): 0.150149
Теперь рассмотрим пример, когда параметры в программе следующие:
П = 8; г = 1; ^(111);
A = randn(n); B = rand(length(A),r);
Q = 1.23*gallery(,toeppd,,n); R = diag(1:(size(в,2)))*3.4; N = rand(length(Q), length(R));
Результат, выводимый в командное окно MATLAB:
Размерность матриц объекта управления. A:8x8, B:8x1 Матрица Ы симметрическая
Размерность матриц функционала. Q: 8x8; R: 1x1; ^ 8x1
R =
3.4000
R2 =
3.4000
Максимальная абсолютная ошибка (R - R2): 4.21041e-11
В заключение приведем результаты расчетов со следующими параметрами:
П = 15; г = 10; ^(1119)
A = randn(n); B = rand(length(A),r);
Q = 1.23*gallery(,toeppd,,n); R = gallery(,lehmer,,r)*6.78;
Результат без вывода матриц R и R2:
Размерность матриц объекта управления.
A:15x15, B:15x10
Размерность матриц функционала.
Q: 15x15; R: 10x10; N: 15x10
Матрица W симметрическая
Максимальная абсолютная ошибка ^ - R2): 0.418571
Заключение. В статье было выполнено доказательство одного соотношения оптимальности в линейно-квадратичной задаче управления с обобщенным функционалом, в котором три весовые матрицы. Полученное соотношение - формула (18) аналитически связывает между собой параметры линейного стационарного объекта управления (матрицы А, В), весовые матрицы квадратичного функционала (О, М) и матрицу оптимального регулятора (К). Данный результат, например, может быть использован при расчетах на компилируемом языке, которым не является MATLAB, служить критерием останова программы численного расчета оптимального регулятора, как на это указывалось в [2]. В программном коде численных расчетов оптимального регулятора с последующей проверкой соотношения оптимальности были учтены предупреждения компилятора системы MATLAB для случая прямой кодировки формулы (18). Подобный учет предупреждений компилятора не был сделан, например, в [3]. Изменение кодировки формулы (18) позволяет повысить быстродействие расчетов.
Список использованных источников
1. Афонин В.В. Аналитические соотношения оптимальности в линейно-квадратичной задаче управления // Математическое моделирование. РАН. 2000. Т. 12. № 3. С. 9.
2. Афонин В.В. Об одном выводе аналитической оценки решения задачи оптимальной стабилизации // Дифференциальная алгебра и динамика систем: межвуз. сб. науч. тр. Саранск: Изд-во Мордов. унта, 2008. С. 78-84.
3. Афонин В.В., Мурюмин С.М. Соотношения оптимальности в линейно-квадратичной задаче управления // Журнал Средневолжского математического общества. 2014. Т.16. № 2. С. 118-120.
4. Ройтенберг Я.Н. Автоматическое управление. М.: Наука,1978. 552 с.
5. Веремей Е.И. Линейные системы с обратной связью: Учеб. пособие. СПб.: Лань, 2013. 448 с.
6. Гуэрра М. Обобщенные решения особых задач оптимального управления // Современная математика. Фундаментальные направления. 2008. Т. 27. С. 60-184.
7. Мижидон А.Д., Имыхелова М.Б. Предельные возможности в задаче активной виброзащиты // Вестник Бурятского государственного университета. 2011. № 9. С. 26-30.
