close

Вход

Забыли?

вход по аккаунту

?

2000-0058-0-01

код для вставкиСкачать
МИНИСТЕРСТВО ОБРАЗОВАНИЯ РОССИЙСКОЙ ФЕДЕРАЦИИ
Санкт-Петербургский
государственный университет аэрокосмического приборостроения
А. А. Умнов
ПРОЕКТИРОВАНИЕ
БОРТОВЫХ КОМПЛЕКСОВ
УПРАВЛЕНИЯ
Текст лекций
Санкт-Петербург
2000
УДК 697.7.05
ББК 32.965.5
У54
Умнов А.А.
У54 Проектирование бортовых комплексов управления: Текст лекций/ СПбГУАП. СПб., 2000. 59 с.
В тексте лекций использованы результаты работ по созданию программно-математического обеспечения систем оптимизации режимов
полета и вертикального маневра различных летательных аппаратов. Иллюстративный материал дан для самолета Ту-154 "М". Рассмотрены вопросы построения математических моделей объекта управления, представляющего собой систему планер–двигатель. Приведены методы и алгоритмы определения оптимальных параметров крейсерского полета, расчета оптимальных траекторий набора высоты и снижения, а также оптимального управления при вертикальном маневрировании.
Рецензенты:
кафедра процессов управления
Балтийского государственного университета;
доктор технических наук профессор Е. Ф. Фурмаков
Утверждено
редакционно-издательским советом университета
в качестве текста лекций
©
©
2
СПбГУАП, 2000
А.А.Умнов, 2000
Предисловие
Современный уровень вычислительной техники позволяет в значительной степени автоматизировать процесс управления полетом. Основу бортового комплекса управления представляет вычислительная система самолетовождения (ВСС). Не касаясь задач аппаратной реализации, в тексте лекций автор рассматривает вопросы создания алгоритмического и программно-математического обеспечения ВСС.
Для гражданской авиации одним из основных требований к ВСС
является обеспечение высокой экономической эффективности эксплуатации летательного аппарата (ЛА). Поэтому одной из важных составляющих ВСС является система оптимизации режимов полета (СОРП),
обладающая развитыми свойствами адаптации к изменяющимся в широких диапазонах условиям полета.
3
1. Общая структура адаптивной оптимальной
системы управления полетом
Задача создания многофункциональной оптимальной системы
управления движением ЛА, обладающей свойствами адаптации
к изменяющимся условиям полета, может быть решена на основе сочетания алгоритмов оптимального оценивания параметров
(идентификации) и состояния (фильтрации) управляемого процесса
с алгоритмами оптимизации управления с использованием, например, прогнозирующих моделей. Адаптивная система управления
движением самолета должна предполагать выполнение следующих процедур:
определение динамических характеристик управляемого объекта в
процессе его функционирования;
оценивание состояния управляемого объекта;
формирование управляющих сигналов с использованием информации, получаемой с помощью первых двух процедур.
В общем виде структура многопараметрической адаптивной системы управления представлена на рис.1.1.
Движение летательного аппарата как объекта управления в общем
случае описывается дифференциальным уравнением
(1.1)
x& = F ( x, a, u, t ) + ? x ,
где x – n-мерный вектор состояния объекта; a – r-мерный вектор параметров, определяемый конструкцией ЛА и свойствами среды; u – mмерный вектор управляющих воздействий; ?x – n-мерный вектор возмущений; F – n-мерная вектор-функция, получаемая на основе теоретических и экспериментальных исследований объекта.
Наблюдение за движением объекта осуществляется с помощью комплекса датчиков, измеряющих компоненты вектора состояния и вектора управления:
z = h(x, a, u, t) + ?z ,
(1.2)
4
Объект управления
Возмущения оx
Управление u
x& = F ( x, a, u , t ) + ? x
x
Шумы оz
Датчики информации
Априорная информация об объекте
z = h ( x, a , u , t ) + ? z
Адаптивная СОРП
z
Программная
настройка
a = ?( z, t )
a
Идентификация
параметров
x?
x? = ? ( z , a? , t0 , t )
a? = H ( z , a , t0 , t )
a?
Цель управления и критерии
Оценивание
состояния
x?
u = ? ( x? , a? , t , tk )
Рис 1.1. Функциональная схема адаптивной СОРП
где z – l-мерный вектор наблюдений; ?z – l-мерный вектор аддитивных
шумов; h – l-мерная вектор-функция, известная на основе теоретических и экспериментальных исследований датчиков информации.
Идеология построения интегрированных комплексов предполагает
комплексирование на уровне чувствительных элементов. В этом случае уравнение (1.2) описывает преобразование информации чувствительными элементами.
Результаты измерений поступают в адаптивную оптимальную систему управления, где используются для определения динамических характеристик объекта и оптимального оценивания его состояния. Рассматриваемая структура предусматривает два типа процессов определения характеристик объекта. Первый из них является в значительной
степени традиционным для авиации и реализует программное восстановление основных характеристик объекта непосредственно по сигна5
лам датчиков. В общем случае программа описывается вектор-функцией:
a = ? ( z, t ) ,
(1.3)
где a – r-мерный вектор из области программных значений параметров
(аэродинамических коэффициентов, инерционно-массовых характеристик, коэффициентов уравнений динамики).
Другой процесс определения динамических характеристик объекта
представляет собой параметрическую идентификацию, использующую
входные и выходные сигналы объекта. В качестве начальных оценок
вектора параметров могут использоваться значения, определяемые (1.3)
В общем виде процесс идентификации описывается оператором
a? (t ) = H ( z , a , t 0 , t ) ,
(1.4)
где a? – r-мерный вектор оценок компонент вектора a.
На основе сигналов датчиков и оценок параметров объекта осуществляется оптимальное оценивание состояния объекта, позволяющее в
значительной степени повысить точность информации о векторе x. Оператор, описывающий оценивание, в общем случае имеет вид
x? (t ) = ? ( z , a? , t 0 , t ) ,
(1.5)
где x? – n-мерный вектор оценок компонент вектора x.
Допускается и отсутствие оценивания состояния в случаях, когда
точность измерения вектора состояния (1.2) достаточна для решения
задач управления.
Итоговой процедурой многопараметрической адаптивной системы
оптимизации режимов полета (СОРП) является оптимизация управляющих сигналов на основе задаваемых цели управления и критериев оптимизации. Оператор, формально описывающий формирование вектора
оптимального управления, имеет вид
u = ?( x? , a? , t , t k )
(1.6)
Путь создания адаптивной СОРП содержит в укрупненном плане
следующие этапы:
1 – разработка математической модели объекта;
2 – формулирование критериев оптимальности и синтез законов оптимального управления;
6
3 – разработка алгоритмов адаптации (настройки) законов управления по режимам функционирования объекта;
4 – реализация полученных алгоритмов с помощью БЦВМ.
Математическая модель, описывающая пространственное движение ЛА и положенная в основу разрабатываемого алгоритмического
обеспечения СОРП, определяет “глубину” управляемых и, следовательно, оптимизируемых процессов. Полнота отражения моделью динамических свойств объекта влияет на эффективность и трудоемкость как
процессов идентификации и оценивания, так и процессов оптимизации.
7
2. Полная нелинейная модель
пространственного движения ЛА
2.1. Системы координат
Для описания движения ЛА используют уравнения в проекциях на
оси выбранных систем координат (СК). В динамике полета получили
распространение следующие правые прямоугольные системы координат.
1. Нормальная земная система координат Fg . Начало находится на
поверхности Земли: в определенной точке взлетно-посадочной полосы
(ВПП), в точке расположения ориентира, в центре наземной цели и т.д.
Оси Oo Xg , Oo Zg расположены в горизонтальной плоскости, а ось OoYg
направлена вверх (вдоль местной вертикали). Ось Oo Xg направлена на
север или совпадает с направлением полета.
2. Нормальная система координат Fн . Начало находится в центре
масс ЛА. В дальнейшем будем полагать, что оси нормальной и нормальной земной СК параллельны.
Относительное положение этих СК определяется вектором r .
3. Связанная система координат Fсв . Основная плоскость OXY является плоскостью симметрии самолета. Оси связанной системы совпадают с продольной ОХ, нормальной ОY и поперечной ОZ осями самолета. Начало координат системы располагается в центре масс ЛА. Относительное положение связанной и нормальной СК определяется в общем случае девятью направляющими косинусами, т.е. косинусами девяти углов между соответствующими осями связанной и нормальной
СК. Часто для определения пространственного положения самолета используют эйлеровы углы тангажа ? , крена ? и рыскания ? .
4. Скоростная (поточная) система координат Fс . Эта система используется, в основном, для определения аэродинамических сил, действующих на самолет. Скоростная ось ОХ? направлена вдоль воздушной скорости самолета V , ось подъемной силы помещается в плоскости симметрии самолета и направлена вверх, боковая ось OZ? образует
с осями OX? и OY? правую систему координат.
8
Положение самолета относительно воздушного потока, определяющее величину аэродинамических сил, задается двумя углами ? и ?, определяющими относительное положение связанной и скоростной систем координат.
По отношению к нормальной системе координат скоростная система повернута на углы ?? , ?? , ?? – скоростные углы рыскания, тангажа
и крена, введенные по аналогии с эйлеровыми углами ?, ?, ? для связанной системы.
5. Полусвязанная система координат Fl . Ее ось OXl совпадает с
проекцией вектора воздушной скорости V на плоскость симметрии самолета, ось OZl – с поперечной осью в связанной. Таким образом, система OXl Yl Zl повернута относительно скоростной на угол скольжения ?
вокруг OYl = OY? и относительно связанной на угол атаки ? вокруг
OZl = OZ.
6. Траекторная система координат Fт . Ось OXk – совпадает с направлением земной скорости Vk самолета. Ось OYk помещается в вертикальную плоскость, проходящую через ось OXk , и направлена вверх
от поверхности Земли. Ось OZk образует правую систему координат.
По отношению к нормальной траекторная система координат повернута на углы ? и ?.
2.2. Уравнения пространственного движения
самолета
При выводе уравнений пространственного движения самолета принимаются следующие допущения:
– конструкция самолета считается жесткой;
– масса самолета в процессе моделирования постоянна, жидкое наполнение отсутствует;
– главные оси инерции ЛА совпадают с осями связанной системы
координат;
– сила тяги двигателя лежит в плоскости симметрии ЛА и составляет с осью ОХ угол ?д;
– вращение Земли и кривизна ее поверхности не учитываются;
– аэродинамика ЛА нелинейная по углам атаки и скольжения;
– атмосфера является стандартной по ГОСТ 4401-81.
9
Уравнения движения самолета относительно инерциальной системы
отсчета могут быть получены из основных теорем динамики твердого
тела, представленных векторными уравнениями
dQ
dK
= F;
=M ,
dt
dt
(2.1)
где Q , K – главный вектор и главный момент количества движения твердого тела относительно центра масс; F , M – главный вектор и главный
момент относительно центра масс внешних сил, действующих на твердое тело.
Рассматривая самолет как твердое тело с постоянной массой m, запишем производную вектора количества движения Q = mV в виде
d
d
Q =m V,
dt
dt
(2.2)
где V – абсолютная скорость центра масс.
Для выполнения дифференцирования в уравнениях (2.1) рассмотрим
несколько вспомогательных задач векторного анализа.
Пусть твердое тело имеет наподвижную точку О. Свяжем жестко с
телом систему координат OXYZ (рис.2.1). Система координат OXYZ
однозначно определяет положение рассматриваемого тела по отношению к неподвижной системе отсчета OX1Y1Z1. Положение произвольной точки M твердого тела определяется радиус-вектором r. Если x, y, z
– координаты точки М в подвижной сиZ
1
Z
стеме координат, а i , j , k – единичные
векторы осей этой СК, то
M
r
k
j
r = xi + yj + zk .
Y
O
Y1
i
X1
X
Рис 2.1. Движение тела с
одной закрепленной точкой
10
(2.3)
Координаты x, y, z точки М в подвижной СК являются постоянными величинами, а единичные вектора i , j , k будут функциями времени, так как система координат движется вместе с телом.
Скорость точки М определяется по формуле
V =
d
r,
dt
(2.4)
поэтому, дифференцируя (2.3) по t , получим
V =x
d
d
d
i + y j+z k .
dt
dt
dt
(2.5)
Умножая обе части равенства (2.5) скалярно на i , j , k получим
di
dj
dk
i +y i +z
i;
dt
dt
dt
di
dj
dk
Vy = V ? j = x
j+y
j+z
j;
dt
dt
dt
di
dj
dk
Vz = V ? k = x k + y k + z
k.
dt
dt
dt
Vx = V ? i = x
(2.6)
Поскольку векторы i , j , k взаимно перпендикулярны, то
i 2 = 1,
i ? j = 0,
j 2 = 1, k 2 = 1;
j ? k = 0, k ? i = 0.
(2.7)
Дифференцируя эти равенства по времени, найдем две группы формул
di
i = 0,
dt
dj
j = 0,
dt
dk
k = 0;
dt
di
dj dj
dk
dk
di
j =? i, k =?
j,
i =? k .
dt
dt dt
dt
dt
dt
Введем в рассмотрение вектор
? = ?x i + ? y j + ?zk ,
(2.8)
(2.9)
(2.10)
проекции которого ? x , ? y , ? z на оси x, y, z определим равенствами:
?x =
dj
dk
k =?
j;
dt
dt
11
dk
di
i = ? k;
dt
dt
di
dj
j =? i.
?z =
dt
dt
?y =
(2.11)
Используя соотношения (2.8), (2.9) и (2.11), перепишем уравнения
(2.6) в виде
Vx = ? y z ? ? z y;
V y = ? z x ? ? x z;
Vz = ? x y ? ? y x.
(2.12)
Рассмотрим векторное произведение
i
? Ч r = ?x
x
(
)
j
?y
y
k
?z =
z
(
)
= ? y z ? ? z y i + (? z x ? ? x z ) j + ? x y ? ? y x k .
Проекции этого векторного произведения на оси x,y, z равны правым
частям соотношения (2.12), следовательно, вектор скорости V может
быть записан в виде
(2.13)
V = ?Чr.
Геометрическое место точек, скорость которых равна нулю, находится из уравнения ? Ч r = 0 , определяющего прямую линию, называемую мгновенной осью вращения. Введенный нами вектор направлен
по мгновенной оси вращения.
Рассмотрим теперь задачу дифференцирования вектора, определенного в системе координат, которая может двигаться произвольным образом. Пусть даны: основная система координат и подвижная система
координат, которая совершает произвольное движение. Пусть вектор
a = a(t) определен в подвижной СК, т.е. проекции этого вектора ax , ay , az
на оси подвижной системы – заданные функции времени. Если i , j , k –
единичные векторы подвижной СК, то вектор a может быть представлен в виде
12
a = axi + a y j + az k .
(2.14)
Установим теперь правило нахождения производной в неподвижной
СК (абсолютной производной) по времени этого вектора. Дифференцируя обе части равенства (2.14) по времени, будем иметь в виду, что
векторы i , j , k вследствие движения подвижной СК, меняют свое направление, т.е. являются функциями времени. Таким образом, абсолютная производная a по времени будет равна
da y
da
d a dax
di
d j
dk
=
+ ay
+ az
i +
j + z k + ax
.
(2.15)
dt
dt
dt
dt
dt
dt
dt
Сумма первых трех слагаемых представляет собой производную от
вектора a в подвижной СК. Эта сумма называется относительной про-
~
da
, т. е.
изводной и обозначается
dt
~
da y
da
d a da x
i+
j+ z k .
=
(2.16)
dt
dt
dt
dt
Заменяя в формуле (2.13) радиус-вектор r последовательно на
i , j , k , получим
di
dj
dk
= ?Ч i ,
= ?Ч j,
= ?Чk .
dt
dt
dt
Поэтому сумма последних трех слагаемых (2.15) может быть представлена в виде
(
)
di
dj
dk
+ ay
+ az
= ? Ч ax i + a y j + az k = ? Ч a ,
dt
dt
dt
где ? – угловая скорость подвижной СК, следовательно,
ax
(2.17)
d a d% a
=
+ ?Чa .
(2.18)
dt
dt
Таким образом, абсолютная производная вектора равна сумме относительной производной этого вектора и векторного произведения угловой скорости подвижной системы на этот вектор.
13
Z1
Z2
Q
Z
VA Y
A
Y2
X
O
X2
X1
Y1
Рассмотрим движение свободного твердого тела . Введем неподвижную СК OX1Y1Z1, подвижную
AX2Y2Z2, перемещающуюся поступательно относительно осей
OX1Y1Z1 и связанную с телом только в одной точке А, и подвижную СК
AXYZ, жестко связанную с телом
(рис.2.2). Скорость точки А ( V A )
определяется выражением (2.4).
Ускорение точки А ( W A ) можно оп-
Рис. 2.2. Движение
свободного твердого тела
ределить как производную вектора
скорости на основания выражения (2.18)
dV A
~
= WA = WA + ? Ч VA ,
dt
~
где W A =
(2.19)
dV y
dV x
dV
i+
j + z k,
dt
dt
dt
V x , V y , V z – проекции вектора скорости на связанные оси; ? – угловая
скорость связанной СК.
Подставляя (2.19) в (2.2), получим уравнение динамики самолета
как твердого тела постоянной массы m в произвольной СК, вращающейся с абсолютной угловой скоростью ?
mV&k + m? Ч Vk = F ,
(2.20)
где Vk – вектор земной скорости ЛА.
Отметим, что выражения (2.2) и (2.20) следуют из основного уравнения динамики
mW = F ,
где m – постоянная масса материальной точки, а W =
14
(2.21)
dV
.
dt
Умножим обе части уравнения (2.21) векторно слева на радиус-вектор r , определяющий положение материальной точки относительно
какой-либо точки O, которую будем называть центром:
(2.22)
r Ч mW = r Ч F .
Преобразуем левую часть этого уравнения следующим образом:
r Ч mW = r Ч m
dV
d
dr
= ( r Ч mV ) ?
Ч mV .
dt
dt
dt
Но dr / dt = V и векторное произведение параллельных векторов
V Ч mV равно нулю, поэтому уравнение (2.22) можно записать в виде
d
(r Ч mV ) = r Ч F .
dt
(2.23)
Вектор K 0 = r Ч mV называется моментом количества движения
материальной точки относительно центра O. Вектор M 0 = r Ч F представляет собой момент силы, приложенной к точке относительно центра. Таким образом
dK 0
= M0.
dt
(2.24)
Момент количества движения K материальной системы относительно центра O (главный момент) определяется суммой моментов количеств движений всех материальных точек, входящих в систему:
n
n
k =1
k =1
K = ? K 0 k = ? rk Ч mk Vk .
(2.25)
Если материальная система представляет непрерывно распределенную материальную среду, заполняющую некоторый объем, то сумма
переходит в соответствующий интеграл.
Рассмотрим твердое тело , вращающееся с угловой скоростью ?
вокруг неподвижной оси, например Z (рис.2.3). Выделим в теле элемент
объема М с массой dm и будем рассматривать его как материальную
точку. При вращении тела вокруг неподвижной оси элемент объема М
будет двигаться по окружности с центром в точке О с радиусом, рав15
ным расстоянию h от точки М до оси
Z
вращения. Проекция скорости V элемента объема М на касательную к окружности равна V? = ? z hz , а проекция
o
количества движения на ту же ось будет V? dm = ? z hz dm .
h
M
V dm
V
Рис. 2.3. Вращение
твердого тела вокруг
неподвижной оси
Поскольку плечо вектора V dm относительно оси вращения равно hz, то
момент количества движения элемента
объема М относительно оси равен
V? dmhz = ? z hz2 dm . Для всего тела бу-
дем иметь K z = ? ? z hz2 dm , где интегрирование распространено на
массу всего тела.
Проекция угловой скорости ?z одинакова для всех точек,и, следовательно, ее можно вынести за знак интеграла
K z = ?z ? hz2 dm .
Получившийся интеграл зависит только от характера распределения
массы в теле и не зависит от его кинематического состояния. Он называется моментом инерции тела относительно оси Z и обозначается символом Iz
I z = ? hz2 dm .
(2.26)
В этих обозначениях будем иметь
K z = I z? z .
(2.27)
Рассмотрим зависимости между параметрами, определяющими положение ЛА как твердого тела (например, углами Эйлера), и силами,
приложенными к телу. Для этого воспользуемся уравнениями (2.1), (2.24),
(2.25).
Пусть неподвижная СК OX1Y1Z1 имеет начало в закрепленной точке
тела. Свяжем жестко с телом подвижную СК с началом в той же точке
О (рис.2.4).
16
При движении тела его моменZ
ты инерции и центробежные моY
менты инерции в неподвижной СК
Z
K
будут переменными величинами,
так как тело при своем движении
Y1
k
j
изменяет свое положение относиO
тельно осей системы.
i
В неподвижной же СК, жестко
X1
связанной с телом, моменты инерX
ции и центробежные моменты инерции являются величинами постоянными. Поэтому будем считать векРис. 2.4. К выводу уравнения
тор момента количества движения
абсолютной производной
K определенным в подвижной СК,
момента количества движения
жестко связанной с телом. Тогда, в соответствии с зависимостью
(2.18), уравнение (2.1) можно записать в виде
~
dK
+ ?ЧK = M,
dt
(2.28)
~
dK y
dK z
d K dK x
i+
j+
k,
=
где
dt
dt
dt
dt
(2.29)
а K x , K y , K z – проекции вектора K на оси связанной СК; ? – угловая
скорость тела.
Имея в виду формулу (2.29), а также известное равенство
i
? Ч K = ?x
j
?y
k
?z = (?y K z ? ?z K y ) i +
Kx
Ky
Kz
(2.30)
+(?z K x ? ?x K z ) j + (?x K y ? ?y K x )k ,
перепишем уравнение (2.28) в проекциях на оси связанной СК
dK x
+ (? y K z ? ? z K y ) = M x ;
dt
17
dK y
+ (? z K x ? ? x K z ) = M y ;
dt
dK z
+ (? x K y ? ? y K x ) = M z ,
dt
(2.31)
где ? x , ? y , ? z – проекции угловой скорости ЛА как твердого тела на
оси связанной СК; M x , M y , M z – главные моменты всех внешних сил
относительно тех же осей.
Если координатные оси X, Y, Z являются главными осями инерции
для точки О, то все центробежные моменты инерции будут равны нулю
и согласно (2.27)
K x = I x?x , K y = I y? y , K z = I z ?z .
(2.32)
Уравнения (2.31) примут следующий вид:
d?x
+ ( I z ? I y )? y ? z = M x ;
dt
d?y
+ ( I x ? I z )? z ? x = M y ;
Iy
dt
d?z
Iz
+ ( I y ? I x )? x ? y = M z .
dt
Ix
(2.33)
Уравнения (2.33) называются динамическими уравнениями Эйлера.
18
3. Уравнения продольного движения ЛА
Представим векторное уравнение движения центра масс ЛА (2.20) в
следующем виде:
R A
V&k = ??Ч Vk + + + g ,
(3.1)
m m
где R – сила тяги; A – аэродинамическая сила; g – приведенная сила
тяжести; Vk – земная скорость.
Далее будем рассматривать продольное движение ЛА (движение
в вертикальной плоскости), поскольку это движение связано с изменением высоты и скорости, т. е. энергии ЛА, и задача оптимизации управления этим движением является наиболее актуальной.
Для связанной СК рассмотрим проекции (3.1) на продольную OX и
нормальную OY оси
R 1
V&x = ?zVy + + (Y? sin ? ? X ? cos ?) ? g sin ? ;
m m
(3.2)
1
(3.3)
V&y = ??zVx + ( X ? sin ? + Y? cos ?) ? g cos ? ,
m
где R – сила тяги, задаваемая в связанной СК; X?– лобовое сопротивление и Y? – подъемная сила – составляющие аэродинамической силы A
задаются в полусвязанной СК; g – приведенная сила тяжести, задаваемая в нормальной СК.
Из (2.33) получаем уравнение вращательного движения относительно центра масс в продольной плоскости
& z = Mz ?
?
1
Iz .
(3.4)
Для определения пространственного положения ЛА к динамическим
уравнениям (3.2) – (3.4) добавляются кинематические соотношения
?& = ? ;
z
19
L& = V x cos ? ? V y sin ? ;
H& = V x sin ? + V y cos ? .
(3.5)
Угол атаки ? может быть определен следующим образом:
? = ?arctg
Vy
.
(3.6)
Vx
При исследовании динамики ЛА его масса m может cчитаться неизменной. Например, для самолета Ту-154”М” часовой расход топлива
составляет около 5 т при массе порядка 90 т.
Для анализа длительных процессов систему (3.2) – (3.5) необходимо
дополнить уравнением
m& = ?GТ ,
где GТ – временной расход топлива.
Аэродинамические силы и моменты определяются следующим образом:
X ? = Cx
?(h)V 2
S,
2
? = Cy
?(h)V 2
S,
2
M z = mz
?V 2
S ? ba .
2
Y
(3.7)
где ?(h) – плотность воздуха; V = Vx2 + V y2 – земная скорость; S –
площадь крыла (для Ту-154 “М” S = 180м2); ba – средняя аэродинамическая хорда (САХ) (для Ту-154“М” ba = 5,28м).
Безразмерные аэродинамические коэффициенты подъемной силы
Cy , лобового сопротивления Cx и момента тангажа mz для Ту-154“М”
рассчитываются по формулам:
C x = C x ( M , ?) ?в = 0
+ ?С xи ( M , ?) ? Kи (?и ) + ?Сxш ( M ) + (Cx (?, ?в) ?
?ст = ?6°
?Cx (?,0)) K?в ( M ) + (Cx (?, ?ст ) ? Cx (?, ?6°)) K? ( M ),
ст
20
C y = C y ( M , ?) ?в = 0
+ ?С yи ( M , ?) ? Kи (?и ) + ?С yш ( M ) + (C y (?, ?в) ?
?ст = ?6°
?C y (?,0)) K?в ( M ) + (C y (?, ?ст ) ? C y (?, ?6°)) K? ( M ),
ст
+ ?mzи ( M , ?) ? Kи (?и ) + ?mzш ( M ) +
mz = mz ( M , ?) ?в = 0
?ст = ?6°
+ (mz (?, ?в ) ? mz (?,0)) K ?в ( M ) + (mz (?, ?ст ) ? mz (?, ?6°)) K ? ( M ) +
ст
?
+C y ( X Т ? X То ) + mz z (?)?z + mz?& (?)?& ,
ba
ba
; ?& = ?&
.
57.3V
57.3V
Функции одной или двух переменных, стоящие в правых частях, задаются таблицами, составляющими банк аэродинамики. Текущие значения находятся линейной интерполяцией. Первым управляющим воздействием при движении в вертикальной плоскости является положение руля высоты ?В , оказывающее наибольшее влияние на mz.
Сила тяги, необходимая для осуществления полета, обеспечивается
силовыми установками. Так, самолет Ту-154”М” оснащен тремя двухконтурными турбореактивными двигателями Д-30КУ, схема которых
приведена на рис.3.1.
где X Т – центровка в % САХ; ?z = ?z
Рис. 3.1. Схема двигателя Д-30КУ
21
Этот двигатель относится к типу двухвальных. Управляющим воздействием является положение сектора газа ?сг , выходными параметрами – тяга R, расход топлива GТ и давление в компрессоре высокого
давления Pk (вторая ступень на рис.3.1). Структурная схема двигателя
Д30КУ приведена на рис.3.2.
Динамическая модель двигателя может быть представлена двумя
дифференциальными уравнениями первого порядка относительно оборотов компрессора высокого давления N2 и компрессора низкого давления N1:
1
N& 2 = ( k 2 ? СГ ? N 2 ),
?2
1
N& 1 = ( k1 N 2 ? N1 ),
?1
R = k R N1 , GТ = kG N 2 ,
(3.8)
где параметры модели двигателя k1, k2, kR, kp, kG , постоянные времени
?1 и ?2 являются, как и при вычислении аэродинамических коэффициентов,. таблично заданными функциями скорости (числа Маха М), высоты, температуры и т.д.
Pk max
Kp
Pk
KR
R
KG
GT
17,8
бсг
Н
k2
?2 p + 1
k1
?1 p + 1
G T min
Н
7 7 7 кг/час
Рис. 3.2. Структурная схема двигателя Д30КУ
22
Таким образом, продольное движение ЛА может быть представлено нелинейным векторным уравнением
(3.9)
x& = F ( x, u ) ,
где x = (Vx ,V y , ?z , ?, L, H , N 2 , N1 )T ; u = (?в , ?СГ )T ; F – вектор-функция, содержащая уравнения (3.2)–(3.5) и (3.8).
Моделирование уравнения (3.9) показало, что динамика самолета
Ту-154 “М” в продольном канале характеризуется двумя видами движений: во-первых, короткопериодическими (с периодом 2–3 с) сильнодемпфированными колебаниями относительно центра масс и, во-вторых, длиннопериодическими (с периодом 100–120 с) слабодемпфированными колебаниями центра масс.
23
4. Учет ветра в уравнениях движения
В процессе полета основным возмущающим воздействием на ЛА
является ветер. Поскольку земная скорость ЛА является суммой воздушной скорости и скорости ветра, уравнение движения центра масс
(3.1) можно представить в виде
R A
V& + ?Ч V = + + g ? (W& + ?Ч W ) ,
m m
(4.1)
где V – воздушная скорость; W – скорость ветра; составляющая, определяемая скоростью ветра, (W& + ?Ч W ) проявляется в виде эффективной силы.
Рассмотрим оценку последнего слагаемого в (4.1). Поскольку ветер
задают в земной СК (данные метеостанций и др.), для связанной СК
? Wx ?
? Wx ?
? h?
? ?
н
W = ? Wy ? = Dсв ? Wyh ?
?
?,
?W ?
?
?
W
? z?
? zh ?
(4.2)
н – матрица перехода из нормальной СК в связанную [1];
где Dсв
WxH ,W yH ,WzH – известные составляющие скорости ветра.
Полные производные компонент скорости ветра в земной СК имеют
вид
? ?W ?
? ?W ?
? ?W ?
? ?W ?
W& xH = ? x ? x&H + ? x ? y& H + ? x ? z&H + ? x ? ;
? ?x ? H
? ?z ? H
? ?t ? H
? ?y ? H
? ?W y ?
? ?W y ?
? ?W y ?
? ?W y ?
W& yH = ?
? x&H + ?
? y& H + ?
? z&H + ?
? ;
? ?x ? H
? ?y ? Н
? ?z ? H
? ?t ? H
24
? ?W ?
? ?W ?
? ?W ?
? ?W ?
W& zH = ? Z ? x&H + ? z ? y& H + ? z ? z&H + ? z ? .
? ?x ? H
? ?z ? H
? ?t ? H
? ?y ? H
Тогда для связанной СК имеем
? W& x
? W& x ?
? H
? ?
н
W& = ? W& y ? = Dсв
? W& yH
?? &
? W& ?
? z?
? WzH
?
? Wx
?
? H
н
&
+
D
?
св ? W yH
??
??
?
? WzH
?
?
?
?? .
?
(4.3)
Таким образом, на основе (4.2) и (4.3) можно оценить составляющую (W& + ?Ч W ) , влияющую на составляющие ускорения центра масс
по связанным осям.
Нами было показано, что безразмерный момент тангажа mz зависит
от ?& , а поскольку tg ? = ?
Vy
Vx
?& =
,то
VyV&x ? VxV&y
Vx2 + Vy2
.
(4.4)
Для связанной СК найдем проекции уравнения (4.1) на продольную и
нормальную оси:
(
)
R 1
V&x = + (Y? sin ? ? X ? cos ? ) ? W& x + ?z Vy + Wy ? g sin ? ,
m m
1
(4.5).
V&y = ( x? sin ? + y? cos ? ) ? W& y ? ?z (Vx + Wx ) ? g cos ? .
m
Подставив (4.5) в (4.4), можно оценить влияние порыва ветра на уравнение моментов, а следовательно, на изменение ориентации ЛА относительно центра масс.
25
5. Задачи минимизаций функций и функционалов
При постановке любой инженерной оптимизационной задачи , в том
числе и задач по определению оптимальных режимов полета, возникает
проблема выбора метода оптимизации. Могут применяться, в зависимости от располагаемых вычислительных мощностей, как сравнительно сложные машинные методы в той или иной вариационной постановке, так и методы оптимизации на основе простых инженерных алгоритмов в рамках вырожденных вариационных задач.
Простейший класс задач оптимизации связан с нахождением значений m координат вектора управления u = (u1 ,K , um )T , минимизирующих скалярный критерий качества L(u). Если на возможные значения u
не наложены какие-либо ограничения (связи) и если функция L(u) имеет
первые и вторые частные производные для любого u, то необходимое
условие минимума функции L по u имеет вид
?L
?L ?L
?L
,
,K ,
) = 0.
=(
?u
?u1 ?u2
?um
(5.1)
Более общий класс задач оптимизации связан с определением значений координат вектора управления u = (u1 ,K , um )T , минимизирующих скалярный критерий качества, зависящий от m + n переменных
L(x, u), причем n координат вектора состояния x связаны с координатами управления с помощью соотношения
f (x, u)=0,
(5.2)
где f (x, u) – n-мерная вектор-функция.
Такая задача носит название задачи на условный экстремум. Рассмотрим геометрическую интерпретацию этой задачи. Пусть x и u
принадлежат пространству R1. Возьмем координатные оси x и u и
построим в них функцию f (x, u) = 0. Дадим функции L(x, u) постоянные значения C1<C2<C3 , тогда на плоскости {x,u} получим линии равного уровня функции L(x, u) (рис.5.1). Очевидно, кривая L(x, u) = C2
имеет единственную точку касания с множеством f (x, u) = 0 . Именно
26
u
P
C0
C1
r
?f
C2
L( x,u)=C 2
C3
r
?L
f ( x,u)=0
x
Рис. 5.1. Геометрическая интерпретация задачи
на условный экстремум
в данной точке функция L(x, u) достигает экстремума на заданном ограничении f (x, u) = 0 . Проведем в точке касания функций касательную Р. Эта касательная будет общей как для f (x, u) = 0 , так и для
L(x, u) = C2 . Градиенты функций L(x, u) = C2 и f (x, u) = 0 в точке
касания будут направлены в разные стороны и перпендикулярны к касательной, поэтому они расположены на одной прямой, откуда следует:
r
r
?L = ???f , или
? ?L ?L ?
? ?f ?f ?
? , ? = ?? ? , ? .
? ?x ?u ?
? ?x ?u ?
(5.3)
Постоянная ? выравнивает, согласовывает градиенты функций
L(x, u) = C2 и f (x, u) = 0 в точке касания.
Из (5.3) следует:
?L
?f
?L
?f
+ ? = 0;
+?
= 0.
(5.4)
?x
?x
?u
?u
Но соотношения (5.4) и условие (5.2) определяют безусловный эк-
стремум новой функции H ( x, u , ? ) = L( x, u ) + ?f ( x, u ) .
(5.5)
Для этой функции условия экстремума (5.1) имеют вид:
?H
?H
?H
= 0;
= 0;
=0,
(5.6)
?x
?u
??
что полностью согласуется с (5.2), (5.4).
Постоянная ? получила название множителя Лагранжа. С помощью
этого множителя задача на условный экстремум сводится к задаче на
безусловный экстремум для функции H(x,u,?).
27
Рассмотрим теперь задачу минимизации функционала при дифференциальных ограничениях, описывающих поведение непрерывной системы. Пусть система описывается нелинейным векторным дифференциальным уравнением
(5.7)
x& = f ( x(t ), u (t )), x(t0 ) – задано, t 0 ? t ? t k .
Здесь x(t) – n-мерный вектор состояния, который определяется выбором m-мерного вектора управления u(t). Введем скалярный критерий
качества
tk
?
I = V ( x(tk ), tk ) + L[ x(t ), u (t ), t ]dt.
(5.8)
t0
Задача состоит в том, чтобы найти вектор-функцию x(t), минимизирующую I. Эта задача оптимального программирования управления u(t)
для непрерывной системы относится к задачам вариационного исчисления. Ее можно рассматривать как предельный случай задач оптимального программирования для дискретной многошаговой системы
(рис.5.2), когда интервал времени между шагами становится малым по
сравнению с общим временем движения.
u (0)
x (0)
x (1)
f
0
u (N - 1)
u (1)
x (2)
f
x (N - 1 )
1
x (N )
f
N -1
Рис. 5.2. Блок-схема дискретной
многошаговой системы
Разобьем интервал наблюдения [t 0, t k] на N малых интервалов. Тогда система (5.7) в дискретном варианте примет вид
(5.9)
x(i + 1) = f i ( x(i ), u (i )); x(0) – задано, i = 0,…, N-1.
Эти уравнения представляют собой последовательность условий в
виде равенств, где x(i) – последовательность значений n-мерного вектора состояния, определяемая в свою очередь выбором последовательности значений m-мерного вектора управления u(i). Критерий качества
(5.3) в рассматриваемом случае примет вид
28
N ?1
I = V ( x( N )) + ? Li [ x(i ), u (i )] = L( x, u ),
i =0
(5.10)
где x = ( x (1),K , x( N ))T – nN-мерный вектор; u = (u (0),K , u ( N ? 1))T –
mN-мерный вектор.
Задача состоит в том, чтобы найти последовательность u(i) (вектор u(t)), которая минимизирует функцию L(x, u) (5.10) при условии (5.9),
которое можно представить в виде
(5.11)
f ( x, u ) = ( f 0 ? x1 , f 1 ? x2 ,K , f N ?1 ? xN )T = 0 ,
где f (x,u) – nN-мерная вектор-функция.
Как было показано ранее, данная задача эквивалентна задаче на безусловный экстремум новой функции (5.5)
H ( x, u , ? ) = L ( x, u ) + ? T ? f ( x, u ) ,
(5.12)
необходимыми условиями экстремума которой являются условия (5.6).
Тогда для системы (5.10), (5.11), (5.12) имеем:
1-е необходимое условие
?H ?L
?f ?V ( xN ) ?L0
=
+ ?T
=
+
?x ?x
?x
?x
?x
0
0
K
? ?1
? 1
? ?f
?1 0
K
? ?x1
?
?f 2
?
+ (?1 ,K, ? N ) ? ? 0
?1
K
?x2
?
K
?K K K
?
?f N ?1
? 0
K
K
?
?x N ?1
?
?L1
?LN ?1
+K+
+
?x
?x
0?
?
0?
?
?
?V ( x N )
?
)+
0 ? = (0, 0K 0,
?x N
?
K?
?
?1 ??
?
+
? ?L1 ?L2
?
?f 2
?LN ?1 ? ? ?f 1
,
,K
, 0 ? + ? ?2
? ? 2 ,K , ?? N ? = 0;
? ?1 , ? 3
?
?
? ?x1 ?x2
?x2
?x N ? 1 ?? ?
?x1
?
?
+?
29
отсюда получаем последовательность значений множителя ?(i)
T
T
? ?f i ?
? ?Li ?
? (i ) = ?
? ? (i + 1) + ?
? , i = 0,1,K, N ? 1
? ?x(i ) ?
? ?x(i ) ?
?
?
?
?
с граничным условием
(5.13)
T
? ?V ( x N ) ?
?? .
? ( N ) = ??
? ?x( N ) ?
(5.14)
2-е необходимое условие
?H ?L
=
+ ?T
?u ?u
? ?f 0
0
?
? ?u0
?
?f 1
? 0
?
?u1
??
0
? 0
? K
K
?
?
0
?? 0
?
?t ?L0 ?L1
?LN ?1
=
+
+K +
+ (?1 , ? 2 ,K , ? N ) ?
?u ?u ?u
?u
?
0 K
0 ?
?
?
0 K
0 ?
? ? ?L0 ?L1
?LN ?1 ?
,
,K ,
?+
? = ??
?u N ?1 ??
0 ? ? ?u0 ?u1
K K
K K
K ?
?
?f N ?1 ?
K K
?
?u N ?1 ??
? ?f 0
?f 1
?f N ?1 ?
+ ? ?1
, ?2
,K , ? N
? = 0.
? ?u0
?
u
u
?
?
?
1
N
1
?
?
Отсюда получаем последовательность соотношений для определения значений оптимального управления u(i)
?Li
?f i
T
+ ? (i + 1)
= 0.
?u (i )
?u (i )
(5.15)
Третье необходимое условие приводит к соотношению (5.11) и, в итоге,
к исходной последовательности (5.9).
30
Переходя к пределу при N, стремящемся к бесконечности в соотношениях (5.13) – (5.15), получим решение исходной непрерывной задачи
для системы (5.7) с функционалом (5.8). Итак, для того чтобы найти
вектор управления u(t) , при котором критерий качества (5.8) достигает
минимума, нужно решить систему дифференциальных уравнений
(5.16)
x& = f ( x, u , t );
? ?f ?
? ?L ?
?& = ?? ? ? ? ? ? ,
? ?x ?
? ?x ?
T
T
5.17)
где u(t) определяется из условия
?L
?f
+ ?T
=0.
?u
?u
(5.18)
Граничные условия для уравнений (5.16), (5.17) разделены: одни из
них заданы при t = t0 , другие – при t = tk :
(5.19)
x(t0) – задано,
T
? ?V (t k ) ?
?? .
? (t k ) = ??
(
)
?
x
t
k ?
?
(5.20)
Отметим теперь, что соотношения (5.16) – (5.20) можно получить,
приравнивая нулю первую вариацию нового вспомогательного критерия
качества (лагранжиана), получаемого из (5.7), (5.8):
tk
I = V ( x k , t k ) + ? {L( x, u , t ) + ?T ? [ x& ? f ( x, u, t )]}dt ,
(5.21)
t0
а также введением для удобства вспомогательной скалярной функции
(гамильтониана)
(5.22)
H ( x , u , ? , t ) = L ( x , u , t ) + ?T ? f ( x , u , t ) .
Этот способ будет нами широко использоваться в дальнейшем.
31
6. Прогнозирующие системы оптимального
управления
Рассмотрим задачу оптимального управления, состоящую в отыскании траектории x(t)?Rn и управления u(t)?Rm, доставляющих минимум функционалу
tk
?
J = V(x(tk), tk)+ L( x, u , t ) dt
(6.1)
to
и удовлетворяющих уравнениям движения
x& = F(x, u, t),
причем различают случаи:
а) u = ? (непосредственный привод);
(6.2)
б) u = ?& (интегрирующий привод).
В соответствии с рецептом Лагранжа введем новый управляемый
лангранжиан
~
L (t, x, x& , ? , u, ? ) = L(x, u, t) + ? T( x& `- F(t, x, ? )),
(6.3)
При этом задача (6.1), (6.2) эквивалентна задаче на безусловный экстремум нового функционала
t
k
~
~
J =V(x(tk),tk) + ? L (t , x, x& , ? , u , ? )dt .
t0
~
Рассмотрим теперь вариацию критерия качества J , соответствующую вариациям вектора управления u(t), вектора состояния x(t), начального t0 и конечного tk моментов времени:
случай u = ? :
tk
?J% = ?[V ( x(tk ), tk ] + ?[ L% (t , x, x& , u , ?)dt ] =
?
t0
32
tk
t
= [Vx ?x + (Vx x& + Vt )?t ]tk + L% ? t tk + ( L% x ? x + L% x& ?x& + L%u ?u + L%? ??)dt.
0
?
t0
Возьмем интеграл от второго слагаемого подынтегрального выражения по частям
tk
?
& = L% x& ? x ttk ?
L% x& ?xdt
0
t0
tk
d
? dt ( L%x& )?xdt.
t0
В результате получим
?J% = [(Vx + L% x& )?x + (Vx x& + Vt + L% )?t ]tk ? ( L% ?t + L% x ?x)t0 +
tk
+ ? [{L% x ?
t0
d %
( Lx& ))?x + L%u ?u + L%? ?? ]dt = 0,
dt
так как для экстремума ?J% = 0 .
Таким образом, гладкие экстремали должны удовлетворять уравнениям Эйлера
d %
( Lx& ) ? L%x = 0;
dt
L%? = 0; L%u = 0
и граничным условиям:
((Vx + L% x& )?x + (Vx x& + Vt + L% )?t )tk ? ( L% ?t + L% x ?x )t0 = 0 ,
(6.4)
(6.5)
где индексы tk и t0 означают, что выражения, стоящие в скобках, необходимо вычислить в конечной и начальной точках траектории; ?t – возможные вариации концов интервала [t0, tk]; ?x – вариация траектории x,
определенная, как и x на интервале, содержащем [t0, tk] с некоторой его
окрестностью.
~
При подстановке явного выражения для L (6.3) в (6.4), получим
?& T + ?T Fx ? Lx = 0;
x& ? F = 0;
Lu ? ?T Fu = 0.
(6.6)
33
Соотношения, аналогичные (6.4), (6.6) в случае интегрирующего привода u = ?& , имеют вид
d %
( Lx& ) ? L% x = 0 ;
?& T + ?T Fx ? Lx = 0 ;
dt
L%? = 0 ;
x& ? F = 0 ;
или
d %
( L & ) ? L%? = 0 ;
dt ?
(6.7)
d
( Lu ) + ?T F? = 0 .
dt
Рассмотрим теперь управляемый гамильтониан, отвечающий лаг-
~
ранжиану L .
а) случай u = ? :
H (t , x, ?, u ) = ?T x& ? L% (t , x, x& , u , ? ) = ?T F ? L ( x, u , t ).
Уравнения Эйлера переписываются следующим образом:
(6.8)
x& = H ?Т ,
x& = F ;
T
T
?& T = ? H x ; или ?& = ?? Fx + Lx .
(6.9)
Оптимальное управление u = u0 находится из соотношения:
H (t , x, ?, u0 ) = max H (t , x, ?, u ) .
u?U
(6.10)
Если U – открытая область, то вместо (6.10) можно записать
(6.11)
Hu= 0 или ?T Fu ? Lu = 0 .
Часто бывает полезно так называемое уравнение Гамильтона-Якоби. Зададимся некоторым управлением u(t). На характеристиках, т. е.
решениях системы Гамильтона (6.9), отвечающих выбранному u(t), комбинация
?T dx ? Hdt
оказывается полным дифференциалом, т. е. существует функция
S (t , x) ? R 1 , такая, что
(6.12)
dS = ?T dx ? Hdt , т.е. S x = ?T , St = ? H .
Последнее равенство в (6.12) дает уравнение Гамильтона-Якоби:
34
S t + H (t , x, S x , u ) = 0 .
(6.13)
б) случай u = ?& (интегрирующий привод):
H (t , x, ?, ?, w) = ?T x& + w?& ? L% (t , x, x& , ?, ?& , ? ) =
= ?T F (t , x, ?) + w?& ? L(t , x, ?& ),
(6.14)
~
здесь w = L?& = Lu .
Уравнения Гамильтона примут вид
x& = H ?T ;
x& = F ;
?& T = ? H x ; или ?& T = ??T Fx + Lx ;
?& = H w ;
?& = u;
w& = ? H ?;
w& = ??T F? .
(6.15)
На характеристиках системы (6.15) комбинация
(6.16)
?T dx + wd ? ? Hdt = dS
является полным дифференциалом скалярной функции S(t, x, ? ) , так
что
S x = ? T ; S ? = w ; St = ? H .
Отсюда получаем уравнение Гамильтона-Якоби
(6.17)
(6.18)
St + H (t , x, S x , ?, S? ) = 0 .
Рассмотрим систему с линейным управлением, квадратичным функционалом и прямым приводом.
Для последнего десятилетия характерно бурное развитие общей теории оптимального управления, базирующейся на рассмотрении пространства состояний и функционалов качества управления. Методы этой
теории все в большей мере конкурируют с классическими методами
анализа и синтеза систем регулирования. К числу таких методов относится аналитическое конструирование регуляторов. Наиболее полное
решение задачи оптимизации в общем виде, характерное для аналитического конструирования, обуславливает интерес к этому методу с точки зрения осуществления совмещенного синтеза оптимальных управлений в адаптивной системе управления.
35
За полноту общего решения аналитическое конструирование расплачивается условиями, накладываемыми на динамические свойства управляемого объекта и на структуру критерия, отражающего предъявляемые к управлению требования. Укажем эти условия.
1. Уравнение (6.2), описывающее движение объекта, линейно относительно вектора управляющих воздействий, т. е. вместо (6.2) можно
пользоваться записью
(6.19)
x& = f ( x, t ) + ?( x, t )u ,
где f и ? – известные векторная и матричная функции указанных аргументов.
2. Область возможных значений управляющих воздействий U ? R m
является незамкнутой.
3. Все возможные переходные функции объекта (6.19) непрерывно
дифференцируемы в R n .
4. Минимизируемый функционал (6.1) является квадратичным относительно вектора управления
tk
t
1 k
J = V ( x(t k ), t k ) + ? Q( x, t )dt + ? u T K ?1 (t )udt ,
2 t0
t0
(6.20)
где K T = K ; ?T K ? > 0 , т. е. K (t ) – положительно определенная невырожденная симметричная матрица.
Перечисленные условия в большинстве случаев приводят к несильным ограничениям в практических задачах.
Управляемый лагранжиан задачи (6.19), (6.20) имеет вид
1
L% (t , x, x& , u , ? ) = Q( x, t ) + uT K ?1 (t )u +
2
+?T (t )( x& ? f ( x, t ) ? ?( x, t )u ).
Управляемый гамильтониан, соответствующий (6.21)
H (t , x, ?, u ) = ?T x& ? L% (t , x, x& , u , ? ) =
1
= ?T ( f ( x, t ) + ?( x, t )u ) ? Q( x, t ) ? uT K ?1 (t )u.
2
36
(6.21)
Гладкие экстремали удовлетворяют уравнениям Эйлера:
x& = H ?Т ;
x& = f + ?u;
?& T = ? H x ; или ?& T = ??T ( f x + ? x u ) + Qx ;
H u = 0;
?T ? ? uT K ?1 = 0.
(6.22)
Последнее из уравнений (6.22), определяющее точку u 0 максимума гамильтониана,
max H (t , x, ? , u ) = H (t , x, ? , u o ) ,
u?U
легко разрешимо; учитывая, что K T = K , получим
(6.23)
u0 = K ?T ?.
Выражение (6.23) определяет оптимальное управление как функцию
t , x, ? . Подставив (6.23) в первые два уравнения (6.22), получим:
x& = f + ?K ?T ?;
?& T = ??T f x ? ?T ? x K ?T ? + Q.
(6.24)
Уравнение Гамильтона-Якоби, связанное с (6.24), имеет вид
1
St + S x f + S x ?K ?T S xT = 0.
2
(6.25)
Чтобы синтезировать оптимальное управление u 0 (t ) , достаточно найти x (t ), ?T (t ) = S x (t , x (t )) и воспользоваться формулой (6.23). Определение траектории x(t ) и импульса ?(t) на промежутке [t 0 , t k ] , т. е. в будущем, является прогнозом, поэтому соответствующая система управления носит название прогнозирующей.
По методу А. А. Красовского прогнозируемая характеристика
x(t ) ,?(t) отвечает неуправляемому движению u (t ) = 0 . Алгоритм вычисления оптимального управления имеет вид
1. Интегрируется исходная система (вычисляется x(t ) на интервале [t 0 , t k ] )
37
x& = f ( x, t ), x(t0 ) = x0 .
2. Интегрируется в обратном времени сопряженная система (вы-
числение ?(t) на интервале [t 0 , t k ] )
?& T = ??T f x (t , x ) + Qx (t , x)
с граничным условием
?T (tk ) = S x (tk , x (tk )) = ?Vx (tk ,( x (tk ));
3. Определяется оптимальное управление
u (t ) = K (t )?T (t , x (t ))? (t ).
Достоинство метода Красовского – значительная экономия вычислительных средств. Недостаток – прогнозирующая характеристика,
отвечающая u = 0 , не совпадает с реальной, а потому u (t ) не совпадает, вообще говоря, с u 0 (t ) .
38
7. Линейная модель движения ЛА
Для решения задач управления часто бывает весьма полезно иметь
линеаризованную модель исследуемого объекта, поскольку линеаризация позволяет значительно упростить исходную нелинейную модель.
Методика линеаризации уравнений движения ЛА известна и включает:
– выбор опорного движения ЛА (часто в качестве такого движения
выбирается горизонтальный крейсерский полет);
– разложение уравнений возмущенного движения в ряд Тейлора относительно опорной траектории;
– удержание в рядах Тейлора только линейных членов (на основании
предположения о малости отклонений возмущенного движения от опорного);
– вычитание из уравнений и соотношений возмущенного движения
соответствующих уравнений и соотношений опорного движения.
В результате применения этой методики модель движения ЛА может быть приведена к модели в малых приращениях относительно произвольного опорного движения ЛА. В форме Коши это модель имеет
вид
(7.1)
?x& = A?x + B?u .
Движение ЛА описывается дифференциальным уравнением
x& = F ( x, a, u , t ),
(7.2)
где x – n-мерный вектор состояния; a – r-мерный вектор параметров, определяемых свойствами среды; u – m-мерный вектор положения
рулевых органов; t – текущее время.
Введем в рассмотрение вектор состояния ЛА xоп (t ) и вектор положения рулевых органов uоп (t ) .
Разложим теперь уравнение (7.2) в окрестности xоп и uоп в ряды
Тейлора, полагая, что эти векторы соответствуют опорному движению.
В результате получим
x&оп + ?x&оп = F ( xоп , a, uоп , t ) +
?F
?x
?x +
оп
?F
?u
?u + O 2 ,
оп
(7.3)
39
?F
где ?x = x ? xоп ; ?u = u ? uоп ;
?x
,
оп
?F
?u
оп
– матрицы частных произ-
водных векторной функции F по компонентам векторов x и u , вычисленные при значениях x = xоп, u = uоп ; O 2 – малые величины более высоко-
го порядка.
При достаточно малых отклонениях x, u от опорных (балансировочных) значений xоп и uоп малыми величинами более высокого порядка в
(7.3) можно пренебречь.
Для автоматизации линеаризации уравнений движения самолета используется процедура линеаризации нелинейных дифференциальных уравнений (7.2). Процедура задает отклонения по составляющим векторов
состояния ? x , управления ?u и вычисляет элементы матриц А и В
линеаризованной системы (7.1) по формулам
aij =
bik =
f i ( x1 ,..., x j + ?x j ,...) ? f i ( x1 ,..., x j ? ?x j ,...)
2?x j
;
f i ( x, u1 ,..., u k + ?u k ,...) ? f i ( x, u1 ,..., u k ? ?u k ,...)
.
2?u k
В практике использования методов аналитического конструирования
оптимальных регуляторов при наличии линейной модели объекта (7.1)
наибольшее применение нашли интегральные квадратичные функционалы качества
t
1 T
1 k T
J (u ) = x (t k )? k x(t k ) + ? ( x ?x + u T ?u )dt ,
2
2 t0
(7.4)
здесь ? k , ? – неотрицательно определенные, а ? – положительно определенная матрицы весовых коэффициентов.
Слагаемое x ?x оценивает отклонение фазовых координат от желаемых. При этом система “штрафуется” за большие ошибки намного
T
больше, чем за малые. Слагаемое u T ?u оценивает “стоимость” управления. Терминальный член x T (t k )? k x(t k ) гарантирует малость
40
ошибки в конечный момент времени t k . Его необходимо учитывать в
том случае, когда влияние величины x(t ) в конечный момент времени
особенно важно.
Выбор элементов весовых матриц ?, ? k и ? определяет динамические свойства оптимальной системы и осуществляется с использованием метода последовательных приближений. Для предварительного
выбора значений весовых коэффициентов, в соответствии с принципом
равных вкладов, пользуются следующими рекомендациями [3].
1. Считаем, что максимальные допустимые отклонения фазовых
координат xi в любой момент времени вносят в функционал качества
одинаковый вклад. Если аналогичные рассуждения распространить и
на отклонения сигналов управления u j , то эти условия можно выразить
с помощью формул
?ii xi2max = ?11 x12max , i = 2,..., n;
? jj u 2jmax = ?11u12max , j = 2,..., r ;
?kii xi2max (tk ) = ?k11 x12max (tk ), i
= 2,..., n;
(7.5)
где ximax – максимальное допустимое отклонение i -й координаты на
интервале T = t k ? t 0 ;
u jmax – максимальное допустимое отклонение j -го сигнала управления
на том же интервале;
ximax (t k ) – максимальное допустимое отклонение i -й фазовой координаты в конечный момент времени t k .
2. Суммарный вклад максимальных допустимых отклонений фазовых координат равняется суммарному вкладу максимальных допустимых сигналов управления
n
?T ?ii xi2
i =1
max
r
n
j =1
i =1
= ? T ? jj u jmax 2 = ? ?kii xi2max (tk ).
(7.6)
41
3. Матрицы ? k , ? и ? выбираются диагональными, какой-либо из
элементов назначается произвольно (например ?11 = 1 ), тогда остальные элементы вычисляются из соотношений (7.5) и (7.6).
Если рассматривать задачу стабилизации относительно начала координат пространства отклонений на достаточно большом интервале
времени, то координаты вектора состояния, являющиеся невязками
между желаемыми и действительными параметрами полета, для устойчивого объекта с течением времени стремятся к нулю. Поэтому
терминальный член в такой задаче можно не рассматривать и критерий
качества примет вид
t
J (u ) =
1 k T
( x ?x + u T ?u )dt.
?
2 t0
(7.7)
Для решения этой задачи запишем расширенный функционал
tk
1
1
J (u ) = [ xT ?x + uT ?u + ?T ( x& ? Ax ? Bu )]dt.
2
2
?
t0
Найдем вариацию этого функционала, соответствующую вариациям траектории ?x и управления ?u ,
tk
?J = ? [ xT ??x + uT ??u + ??T ( x& ? Ax ? Bu ) + ?T (?x& ? A?x ? B?u )]dt.
t0
Берем по частям следующий интеграл
tk
??
T
& = ? ?x
?xdt
T
t0
tk
t0
tk
? ?& T ?xdt.
?
t0
Получаем, приравнивая к нулю вариацию расширенного функционала, краевое условие ?T (tk ) = 0 и уравнения Эйлера
x& = Ax + Bu ,
x(t0 ) = x0 ;
?& = ? AT ? + ?x, ?(tk ) = 0;
u = ? ?1BT ?.
(7.8)
Будем искать ? в виде ? = ? Kx. Тогда система (7.8) примет вид
42
x& = Ax + Bu ,
? Kx& = AT Kx + ?x.
Умножим первое уравнение на матрицу K слева и сложим; получим
( KA ? KB? ?1 B T K + AT K + ? ) x = 0.
Поскольку это уравнение должно выполняться для любого x , то
(7.9)
KA + AT K ? KB? ?1 B T K + ? = 0.
Это матричное алгебраическое уравнение Риккати для определения матрицы K. Решив его, можно получить оптимальное управление в
виде обратной связи
(7.10)
u = ?? ?1BT Kx.
Для решения алгебраического уравнения Риккати можно численно интегрировать дифференциальное уравнение Риккати
K& = KA + AT K ? KB? ?1 B T K + ? до получения стационарного значе-
ния. Эта процедура осуществляется достаточно долго, особенно при
высокой размерности рассматриваемой задачи. Поэтому часто используют другой путь – записывают уравнение Риккати в виде матричной
функции матричного аргумента
F ( K ) = KA + AT K ? KB? ?1 B T K + ?,
(7.11)
корень которой нужно найти;
При использовании итеративной процедуры Ньютона, матричную
n(n + 1)
скалярных нели2
нейных уравнений (в силу симметрии матрицы К)
функцию (7.11) представляют в виде системы
?( y ) = 0,
yT = ( K11 ,..., K1n , K 22 ,..., K 2n , K33 ,..., K3n ,..., K nn ),
(7.12)
n(n + 1)
– мерная вектор-функция.
2
Тогда очередное приближение решения уравнения (7.12) находится
из соотношения
где ? ?
yi = yi ?1 ? (? y )?1 ?( yi ?1 ),
43
где ? y – матрица частных производных вектор-функции ? по вектору y в точке y i ?1 ,
?( y ) = ( f11 ,..., f1n , f 22 ,..., f 2 n , f33 ,..., f3n ,..., f nn )T ,
f ij – элементы матричной функции F (K ) , которые выражаются через
параметры исследуемой системы в следующем виде:
f ji =
n
n
k =1
k =1
n
n
? k jk aki + ? a jk kki + ? ji ? ?? k jk mkl kli ,
k =1 l =1
где kij , aij , ?ij – элементы матриц K , A, ? ; mij – элементы матрицы
M = B? ?1 B T .
44
8. Анализ процессов управления ЛА
на различных этапах полета
8.1. Управление в соответствии с руководством
по лётной эксплуатации
На этапе набора высоты двигатели работают в номинальном режиме. Это означает, что сектором газа выдерживаются N2ном = 92,5%
N2max = const. В канале руля высоты используют два режима: максимальный дальности (МД) и максимальный крейсерской скорости (МКР).
Для самолёта Ту-154 “М” первому соответствует постоянная приборная скорость Vпр = 550 км/ч, второму – Vпр = 575 км/ч. Приборная скорость (с некоторыми поправками, определяемыми по таблицам, учитывающим влияние сжимаемости воздуха) равна такой истинной (характеризующей дальность, проходимую ЛА в воздухе за единицу времени), при которой в стандартных условиях у Земли получался тот же скоростной напор
q = ??
Таким образом, Vпр ? Vи ?
V2
.
2
?н
,
?0
В крейсерском полёте используется три режима: МД, МКР и экономичный режим (ЭКР). Скорости определяются по соответствующим
номограммам в зависимости от массы, температуры и ветра.
На этапе снижения двигатели работают в режиме “малого газа”. В
канале руля высоты на режиме МКР выдерживается постоянная приборная скорость Vпр = 575км/ч, а на режиме МД – Vпр = 550 км/час.
Указанные алгоритмы удобны в режиме ручного пилотирования, но наверняка не являются оптимальными.
45
8.2. Расчёт оптимальных параметров крейсерского полёта
Показатель стоимости полёта представляют в виде суммы двух составляющих, зависящих от затрачиваемого топлива mT и времени t.
J = Cт ? mТ + Сt ? t ,
(8.1)
где Сm – стоимость единицы массы топлива; Сt – стоимость единицы
времени полёта.
Рассматривая элементарный участок dL крейсерского полёта, можно
записать
(V + U)dt = dL,
где V – воздушная скорость; U – попутная (встречная) составляющая
скорости ветра.
Тогда критерий (8.1) можно представить в виде
Lкр
С ? q + Ct
?mТ
J = ? (Cт
dL =
+ Ct ) dt = ? т t
V +U
?t
0
0
C
0
Lкр Ст ? qкм + t
V dL ,
= ?
U
0
1+
V
t кр
0
здесь qкм =
(8.2)
qt
– штилевой километровый расход топлива.
V
Ядро интеграла (8.2) представляет собой монотонную функцию и
отражает стоимость одного километра полёта. Поэтому для крейсерских режимов критерий оптимальности принимает вид
J кр =
где I =
0
( m,V , H ) +
qкм
1+
U
V
Ct
– индекс стоимости топлива.
Cm
I
V ,
(8.3)
Рассмотрим алгоритм вычисления оптимальной скорости крейсерского полёта для самолёта Ту-154 “М”. Критерий (8.3) применительно к
этому объекту имеет вид
46
J кр =
q0 [1 + ( a2 + a3 M)?M 2 ] +
I
3.6 ? aн ? M
U
1+
aн М
,
(8.4)
где aн – скорость звука на высоте H; qкм = q0[1+(a2+a3M)?M2] – аппроксимация крейсерского километрового расхода топлива; ?M = М ? М0;
q0 = fq(G, H) – аппроксимация минимального километрового крейсерского расхода топлива; M0 = fм(G, H) – аппроксимация числа М для минимального километрового расхода топлива.
Необходимое условие минимума (8.4)
?J кр
=0
?M
?M
приводит к следующему уравнению относительно x =
:
M0
x ? (1 + x 2 ) ? (1 + B ? x) = A ,
I
4x
? q0U ? q0UM 02 x{2(a2 + a3M 0 ) + 3 x[ a2 + a3M 0 + a3M 0 (1 + )]}
3
A = 3.6
,
3
2(a2 + a3 M 0 )M 0 q0 aн
B=
3a3M 0
.
2(a2 + a3M 0 )
(8.5)
Уравнение (8.5) решается с помощью итеративного алгоритма
y ( xi +1 ) = A( xi ),
xi +1 = f ( A( xi )),
x0 = 0,
где f (?) – функция, обратная функции y (?).
График функции y = x(1+x)2(1+Bx) показан на рис.8.1.
Эта функция имеет экстремумы при x1 = ?1, а также при
x2,3 = ?
1
3
1
B + ± ( B ? )2 + 2
2
2
.
47
y
x*
-B
x
-1
y*
Рис. 8.1. График функции O(N)
Значению x* соответствует глобальный минимум
x* = ?
1
3
1
B + + ( B ? )2 + 2
2
2
.
Значение x* является (рис.8.1) левой границей интервала монотонного возрастания функции y(x), т.е., начиная от точки с координатами
(x*, A*), где A* = x*(1 + x* )2(1 + Bx*) может быть построена обратная
функция x = f(h), вид которой представлен на рис.8.2.
Для того чтобы исследуемая функция была чисто положительной
переходим к новым переменным x1 = x + x* , A1 = A + A*. Вид новой
функции x1 = f1 (h1) приведён на рис.8.3. Поскольку вблизи минимума
(x*, A*) исходной функции y(x) = A (рис.8.1) разложение в ряд Тейлора
не имеет линейных слагаемых, т.е.
?
?1
y*
h
h1
x*
Рис. 8.2. График
обратной функции
48
Рис. 8.3. График
смещенной функции
y = y ( x*) +
1 ''
y ( x*)( x ? x*)2 + O( x ? x*)3 ? y * + a ( x ? x*)2 ,
2
то при y = y*
y ? y*
,
a
т.е. обратная функция f1 (h1) близка к функции “корень квадратный”.
Для “выпрямления” этой зависимости проводится вторая замена переменной
x ? x*+
A 2=
A 1,
при этом зависимость x1 = f2 (A2) становится “почти линейной”. Для
получения искомого x делается третья замена (обратное смещение)
A3 = A 2 ? x* = A + A * ? A*,
а полученная “почти линейная” зависимость
x = f? ( A + A * ? A *)
для повышения точности аппроксимируется дробно-рациональной функцией fx(A, B). Итеративный алгоритм (8.6) в этом случае примет вид
x0 = 0,
x1 = A( x0 ),
xi +1 = f x ( A( xi ), B ), i = 1, 2.
8.3. Оптимизация управления ЛА
на режиме набора высоты
Задача отыскания оптимального управления при наборе высоты (снижении) может быть решена двумя способами. Первый путь предполагает синтез регулятора с обратной связью, переводящего самолёт из
точки начала набора высоты (снижения) с параметрами H0,V0 в точку с
параметрами Н3, V3 заданного крейсерского режима.
Регулятор формирует управляющие воздействия в каналах скорости
и высоты в виде обратных связей по невязкам между текущим и заданным состояниями ЛА. Такие управляющие устройства могут быть получены методами аналитического конструирования оптимальных регуляторов.
49
Другой путь – формирование оптимальной программы набора высоты и стабилизации ЛА относительно этой программы.
Для расчёта оптимальной по минимуму расхода топлива программы
набора высоты (снижения) рассмотрим уравнение баланса топлива на
всей дальности полёта (от L = 0 до L = Lк)
mTк = mTнв + mTкр + mTсн ,
где mTк , mTнв, mTкр , mTсн – массы израсходованного топлива соответственно на всей дальности полета Lк , дальностях набора высоты, крейсерского участка и снижения.
Крейсерский участок полета Lкр , величина которого определяет расход топлива и времени в крейсерском полёте, зависит от протяжённости
участков набора высоты Lнв и снижения Lсн
Lкр = Lк ? Lнв ? Lсн .
Представим уравнение баланса топлива на всей дальности полёта в
форме
где
mTк = mТ*нв + mТ*кр + mТ*сн ,
mТ*кр =
LК
? qкр dL;
0
Lк
mТ*сн = mТсн ?
?
qкр dL;
Lк ? Lсн
mТ*нв = mТнв ?
Lнв
? qкр dL;
0
Таким образом, если крейсерский режим полёта известен, то величина топлива, расходуемого на всей траектории полёта, распадается на
три не зависимые друг от друга составляющие m*Тнв , m*Ткр , m*Тсн :
m*Тнв – часть топлива, расходуемая при наборе высоты только на
увеличение полной энергии самолёта (потенциальной и кинетической),
являющаяся функцией программы набора высоты;
m*Ткр – для заданной дальности полёта Lк является постоянной величиной и представляет собой условную величину топлива, расходуемого
50
ЛА на преодоление всей заданной дальности Lк на крейсерском режиме;
m*Тсн – часть топлива, расходуемая при снижении только на уменьшение полной энергии самолёта, являющаяся функцией программы снижения.
Поскольку m*Тнв и m*Тсн – топливо, расходуемое только на изменение
полной энергии самолёта при наборе высоты или снижении, то m*Тнв и
m*Тсн называются энергетическими частями топлива или энергетическими параметрами.
Таким образом, для выполнения условия минимального расхода топлива на всей дальности Lк необходимо найти программы набора высоты и снижения соответственно
(m*Тнв)min и (m*Тсн )min.
Поскольку в модели двигателя, используемой при синтезе оптимальной программы набора высоты, вычисляется мгновенное значение часового расхода топлива qt , представим минимизируемый критерий в виде
J нв =
тТ*нв
=
tнв
?
0
?тТнв
dt ? qtкр ? tнв ,
?t
(8.7)
здесь крейсерский секундный расход qtкр в первом приближении считается постоянным.
Рассмотрим процедуру расчёта программы набора высоты. Энергетический метод определения оптимальной программы набора высоты предполагает введение новой независимой переменной – энергетической высоты Hэ, объединяющей такие параметры полёта как высота
и скорость
Hэ = H +
V2
; V 2 = Vx2 + V y2 ,
2g
(8.8)
где Vx и Vy – проекции воздушной скорости на связанные оси.
На участке набора высоты (снижения) Нэ монотонно зависит от времени и, следовательно, может использоваться в качестве новой независимой переменной, при этом
d dHэ d
& d .
=
?
= Hэ
dt
dt dHэ
dHэ
51
Критерий (8.7) относительно энергетической высоты может быть
представлен в виде
J нв =
H ЭK
?
H ЭН
=
?тТнв 1
dН ? qt кр
&
?t Hэ
H ЭK
H ЭK
?
H ЭН
1
dHэ =
&
Нэ
N
qt ? qt кр
dHэ ? ? ?i ?Hэ.
&
Нэ
?
(8.9)
i =0
H ЭН
Интервал интегрирования [ H Эн , H Эк ] разбивают на N подынтервалов,
веса каждого из которых di(Vx, Vy, H, N2) и определяют итоговое значение критерия (8.9).
Поэтому на каждом шаге решается задача минимизации
?i (Vx ,V y , H , N 2 ) ?????
?
? min .
V ,V , H , N
x
y
2
Модель продольного движения ЛА в проекциях на связные оси имеет вид
.
Vx = V y ?Z ? g sin ? +
R ? ?V 2 S
+
(C y sin ? ? C x cos ?);
m
2m
.
Vy = ?Vx ?Z ? g cos ? +
.
?z =
? ?V 2S
(Cx sin ? + C y cos ?);
2m
? ? V 2 Sba
?
? mz ;
2
Iz
.
? = ?z ;
.
L = Vx cos ? ? Vy sin ?;
.
H = Vx sin ? + Vy cos ?;
Двигатель:
52
.
N 2 = f N 2 (V , H , ? сг );
.
N 1 = f N1 (V , H , N 2 );
алгебраические соотношения:
V = Vx2 + V y2 ;
? = ?arctg
Vy
Vx
;
R = f R (V , H , N1 ); qt = f q (V , H , N 2 ).
Производная по времени от энергетической высоты, присутствующая в выражении (8.9), записывается следующим образом:
.
.
. 1
.
.
1 .
V V + H = H + (Vx Vx + V y V y ).
(8.10)
g
g
С учётом (8.10) выражение для энергетического топлива ?i примет
вид
HЭ =
?i =
=
qti (Vx ,V y , H i , N 2i ) ? qt
i
i
H& Э
qti (Vxi ,V y , H i , N 2i ) ? qt
i
R (Vx ,V y , H i , N1i )
i
i
mg
кр
=
кр
?( H i )(Vx2 + V y2 )
+
i
i
[Vx (C y sin ? ? C x cos ?) + V y (C x sin ? + C y cos ?)]
i
(8.11)
2mg
,
i
где все параметры взяты из приведённой модели ЛА.
Необходимо отметить, что хотя ?i формально зависит от 4 аргументов, одна из переменных H,Vx,Vy может быть выражена через остальные в соответствии с (8.8).
Приращения высоты ?Hi и скорости по связным осям ?Vxi и ?Vyi
при выполнении условия ?Hэ<<I связаны соотношением
?Hэ ? ?Hi +
Vxi ?1
g
?Vxi ?1 +
V yi ?1
g
?V yi = const.
(8.12)
53
Если принять в качестве независимых координат ?Vxi и ?Vyi то область изменения этих координат будет, в соответствии с (8.12), следующей:
0 ? ?Vxi ?
0 ? ?Vyi ?
?H Э g
,
Vxi ?1
?H Э g Vxi ?1
?
?Vxi ,
V yi ?1 V yi ?1
(8.13)
т. е. ограничение на вторую переменную является зависимым. Третья
переменная ?Hi определяется с необходимостью из соотношения (8.12)
в виде
?H i = ?H Э ?
Vxi ?1
g
?Vxi ?
V yi ?1
g
?V yi .
(8.14)
Четвёртый аргумент функции ?i (обороты двигателя) не зависят от
энергетической высоты и ограничены максимальным и минимальным
допустимыми значениями.
(8.15)
N 2 min ? N 2i ? N 2 max .
Таким образом, решая на каждом шаге ?Hэ задачу минимизации
энергетического расхода ?i ,i = 1,.., N, как функции трех переменных,
область изменения которых задана соотношениями (8.13) – (8.15), находим оптимальное сочетание Vx, Vy, H и N2, а следовательно, можем получить программы V = V(H) и N2= N2(H), задаваемые таблично или
аппроксимационными выражениями.
8.4. Управление ЛА на режиме смены эшелона
Исследования показали, что модель ЛА для крейсерского режима
может быть успешно аппроксимирована линейным приближением
.
? X = A?X + B?U .
Задача состоит в переводе ЛА из крейсерского режима с параметрами (Н0,V0) в крейсерский режим с параметрами (Н3,V3) .
Анализ общей зависимости расхода топлива
qt = qt (V , H , N 2 )
54
(8.16)
от скорости, высоты и оборотов компрессора высокого давления показал, что для крейсерских режимов эта зависимость с погрешностью не
более 5% может быть представлена в виде
qt = qt з + K T ?X ,
(8.17)
где qtз – расход, отвечающий заданному крейсерскому режиму с параметрами (H,V)з; ?X – невязка между текущими и заданными координатами состояния ЛА.
Критерий оптимальности должен минимизировать затраты на режиме смены эшелона и может быть представлен с учётом (8.17) в виде
J=
tK
2
? (qt ? qtз ) dt =
0
tК
? ?X
T
K ? K T ? Xdt.
(8.18)
0
Для обеспечения требуемого качества процесса регулирования (апериодический характер вписывания в заданную крейсерскую траекторию, ограничения по перегрузкам) в ядро функционала (8.18) добавляется квадратичная форма управлений
J=
tR
? (? X
T
? ? ? X + ?U T ? ? ?U )dt ,
(8.19)
0
где ? = К?КТ, ? – подбирается на основе эксперимента с выполнением
требований, сформулированных в разд. 7;
?U = U-Uз – невязка между текущим и балансировочным заданным
управлениями.
Поскольку объект линеен (с матрицами А и В), функционал (8.19)
квадратичен, оптимальный регулятор формирует управление в виде
?U опт = ?? ?1BT S ? ?X ,
где ?U = (?? сг , ?? В )Т ;
X T = (Vx ,V y , ?z , ?, H , N 2 , N1 ); S – решение мат-
ричного алгебраического уравнения Риккати
SA + AT S ? SB? ?1B T S + ? = 0.
В зависимости от полётной обстановки матрица коэффициентов линейных обратных связей
P = ?? ?1BT S
55
меняется, поэтому для различных точек области параметров полётной
обстановки могут быть рассчитаны значения Р, а затем элементы матрицы P аппроксимированы аналитическими выражениями вида
pij = pij ( H ,V , ?Tн , m),
где ?Тн – отклонение температуры от стандартной; m – текущая масса
ЛА.
В связи с изложенным, бортовой алгоритм оптимального управления ЛА на этапе смены эшелона содержит следующие процедуры.
1. Блок вычисления по Нз и Vз (параметры заданного крейсерского
режима) значений заданных координат вектора состояния
X з = (Vx ,V y , ?z , ?, H , N 2 , N1 )Тз
и вектора управления
U з = (? сг , ?в )Тз .
2. Блок формирования сигналов управления
U = U з + Р ( X ? X з ).
Приведённые алгоритмы оптимизации управления на всех этапах
полёта (наборе высоты, крейсерских участках, вертикальных манёврах
и снижении) очерчивают общий подход к алгоритмическому обеспечению бортовых вычислительных процедур, реализуемых в СОРП.
При практической реализации СОРП могут быть разработаны различные модификации для конкретных ЛА, отличающиеся блоками сменных коэффициентов, однако подход к их построению, основные функциональные связи будут идентичны.
56
Библиографический список
1. Буков В. Н. Адаптивные прогнозирующие системы управления полётом. М.: Наука,1987. 230 с.
2. Черкасов Б. А. Автоматика и регулирование воздушно-реактивных двигателей. М.: Машиностроение, 1988. 359 с.
3. Брайсон А., Хо Ю-ши. Прикладная теория оптимального управления. М.: Мир, 1972. 540 с.
4. Скрипниченко С. Ю. Оптимизация режимов полёта по экономическим критериям. М.: Машиностроение, 1988.151 с.
5. Лигум Т. И., Скрипниченко С.Ю., Шишмарёв А.В.: Аэродинамика самолёта ТУ – 154. М.: Транспорт, 1985. 263 с.
57
Оглавление
Предисловие ............................................................................................
1. Общая структура адаптивной оптимальной системы управления
полетом ................................................................................................
2. Полная нелинейная модель пространственного движения ЛА .....
2.1. Системы координат ................................................................
2.2. Уравнения пространственного движения самолета ..........
3. Уравнения продольного движения ЛА .............................................
4. Учет ветра в уравнениях движения ..................................................
5. Задачи минимизаций функций и функционалов................................
6. Прогнозирующие системы оптимального управления ...................
7. Линейная модель движения ЛА ........................................................
8. Анализ процессов управления ЛА на различных этапах полета ...
8.1. Управление в соответствии с руководством по лётной
эксплуатации ............................................................................
8.2. Расчёт оптимальных параметров крейсерского полёта ....
8.3. Оптимизация управления ЛА на режиме набора высоты
8.4. Управление ЛА на режиме смены эшелона ........................
Библиографический список ....................................................................
58
3
4
8
8
9
19
24
26
32
39
45
45
46
49
54
57
Учебное издание
Умнов Александр Алексеевич
ПРОЕКТИРОВАНИЕ
БОРТОВЫХ КОМПЛЕКСОВ
УПРАВЛЕНИЯ
Текст лекций
Редактор В.П. Зуева
Компьютерная верстка Ю. С. Бардуковой
Лицензия ЛР № 020341 от 07.05.97. Сдано в набор 26.09.00. Подписано к печати 20.12.00.
Формат 60Ч84 1/16. Бумага тип. №3. Печать офсетная. Усл. печ. л. 2,9. Усл. кр.-отт. 3,0.
Уч. -изд. л. 3,15. Тираж 100 экз. Заказ №
Редакционно-издательский отдел
Лаборатория компьютерно-издательских технологий
Отдел оперативной полиграфии
СПбГУАП
190000, Санкт-Петербург, ул. Б. Морская, 67
59
щую (W& + ?Ч W ) , влияющую на составляющие ускорения центра масс
по связанным осям.
Нами было показано, что безразмерный момент тангажа mz зависит
от ?& , а поскольку tg ? = ?
Vy
Vx
?& =
,то
VyV&x ? VxV&y
Vx2 + Vy2
.
(4.4)
Для связанной СК найдем проекции уравнения (4.1) на продольную и
нормальную оси:
(
)
R 1
V&x = + (Y? sin ? ? X ? cos ? ) ? W& x + ?z Vy + Wy ? g sin ? ,
m m
1
(4.5).
V&y = ( x? sin ? + y? cos ? ) ? W& y ? ?z (Vx + Wx ) ? g cos ? .
m
Подставив (4.5) в (4.4), можно оценить влияние порыва ветра на уравнение моментов, а следовательно, на изменение ориентации ЛА относительно центра масс.
25
5. Задачи минимизаций функций и функционалов
При постановке любой инженерной оптимизационной задачи , в том
числе и задач по определению оптимальных режимов полета, возникает
проблема выбора метода оптимизации. Могут применяться, в зависимости от располагаемых вычислительных мощностей, как сравнительно сложные машинные методы в той или иной вариационной постановке, так и методы оптимизации на основе простых инженерных алгоритмов в рамках вырожденных вариационных задач.
Простейший класс задач оптимизации связан с нахождением значений m координат вектора управления u = (u1 ,K , um )T , минимизирующих скалярный критерий качества L(u). Если на возможные значения u
не наложены какие-либо ограничения (связи) и если функция L(u) имеет
первые и вторые частные производные для любого u, то необходимое
условие минимума функции L по u имеет вид
?L
?L ?L
?L
,
,K ,
) = 0.
=(
?u
?u1 ?u2
?um
(5.1)
Более общий класс задач оптимизации связан с определением значений координат вектора управления u = (u1 ,K , um )T , минимизирующих скалярный критерий качества, зависящий от m + n переменных
L(x, u), причем n координат вектора состояния x связаны с координатами управления с помощью соотношения
f (x, u)=0,
(5.2)
где f (x, u) – n-мерная вектор-функция.
Такая задача носит название задачи на условный экстремум. Рассмотрим геометрическую интерпретацию этой задачи. Пусть x и u
принадлежат пространству R1. Возьмем координатные оси x и u и
построим в них функцию f (x, u) = 0. Дадим функции L(x, u) постоянные значения C1<C2<C3 , тогда на плоскости {x,u} получим линии равного уровня функции L(x, u) (рис.5.1). Очевидно, кривая L(x, u) = C2
имеет единственную точку касания с множеством f (x, u) = 0 . Именно
26
u
P
C0
C1
r
?f
C2
L( x,u)=C 2
C3
r
?L
f ( x,u)=0
x
Рис. 5.1. Геометрическая интерпретация задачи
на условный экстремум
в данной точке функция L(x, u) достигает экстремума на заданном ограничении f (x, u) = 0 . Проведем в точке касания функций касательную Р. Эта касательная будет общей как для f (x, u) = 0 , так и для
L(x, u) = C2 . Градиенты функций L(x, u) = C2 и f (x, u) = 0 в точке
касания будут направлены в разные стороны и перпендикулярны к касательной, поэтому они расположены на одной прямой, откуда следует:
r
r
?L = ???f , или
? ?L ?L ?
? ?f ?f ?
? , ? = ?? ? , ? .
? ?x ?u ?
? ?x ?u ?
(5.3)
Постоянная ? выравнивает, согласовывает градиенты функций
L(x, u) = C2 и f (x, u) = 0 в точке касания.
Из (5.3) следует:
?L
?f
?L
?f
+ ? = 0;
+?
= 0.
(5.4)
?x
?x
?u
?u
Но соотношения (5.4) и условие (5.2) определяют безусловный эк-
стремум новой функции H ( x, u , ? ) = L( x, u ) + ?f ( x, u ) .
(5.5)
Для этой функции условия экстремума (5.1) имеют вид:
?H
?H
?H
= 0;
= 0;
=0,
(5.6)
?x
?u
??
что полностью согласуется с (5.2), (5.4).
Постоянная ? получила название множителя Лагранжа. С помощью
этого множителя задача на условный экстремум сводится к задаче на
безусловный экстремум для функции H(x,u,?).
27
Рассмотрим теперь задачу минимизации функционала при дифференциальных ограничениях, описывающих поведение непрерывной системы. Пусть система описывается нелинейным векторным дифференциальным уравнением
(5.7)
x& = f ( x(t ), u (t )), x(t0 ) – задано, t 0 ? t ? t k .
Здесь x(t) – n-мерный вектор состояния, который определяется выбором m-мерного вектора управления u(t). Введем скалярный критерий
качества
tk
?
I = V ( x(tk ), tk ) + L[ x(t ), u (t ), t ]dt.
(5.8)
t0
Задача состоит в том, чтобы найти вектор-функцию x(t), минимизирующую I. Эта задача оптимального программирования управления u(t)
для непрерывной системы относится к задачам вариационного исчисления. Ее можно рассматривать как предельный случай задач оптимального программирования для дискретной многошаговой системы
(рис.5.2), когда интервал времени между шагами становится малым по
сравнению с общим временем движения.
u (0)
x (0)
x (1)
f
0
u (N - 1)
u (1)
x (2)
f
x (N - 1 )
1
x (N )
f
N -1
Рис. 5.2. Блок-схема дискретной
многошаговой системы
Разобьем интервал наблюдения [t 0, t k] на N малых интервалов. Тогда система (5.7) в дискретном варианте примет вид
(5.9)
x(i + 1) = f i ( x(i ), u (i )); x(0) – задано, i = 0,…, N-1.
Эти уравнения представляют собой последовательность условий в
виде равенств, где x(i) – последовательность значений n-мерного вектора состояния, определяемая в свою очередь выбором последовательности значений m-мерного вектора управления u(i). Критерий качества
(5.3) в рассматриваемом случае примет вид
28
N ?1
I = V ( x( N )) + ? Li [ x(i ), u (i )] = L( x, u ),
i =0
(5.10)
где x = ( x (1),K , x( N ))T – nN-мерный вектор; u = (u (0),K , u ( N ? 1))T –
mN-мерный вектор.
Задача состоит в том, чтобы найти последовательность u(i) (вектор u(t)), которая минимизирует функцию L(x, u) (5.10) при условии (5.9),
которое можно представить в виде
(5.11)
f ( x, u ) = ( f 0 ? x1 , f 1 ? x2 ,K , f N ?1 ? xN )T = 0 ,
где f (x,u) – nN-мерная вектор-функция.
Как было показано ранее, данная задача эквивалентна задаче на безусловный экстремум новой функции (5.5)
H ( x, u , ? ) = L ( x, u ) + ? T ? f ( x, u ) ,
(5.12)
необходимыми условиями экстремума которой являются условия (5.6).
Тогда для системы (5.10), (5.11), (5.12) имеем:
1-е необходимое условие
?H ?L
?f ?V ( xN ) ?L0
=
+ ?T
=
+
?x ?x
?x
?x
?x
0
0
K
? ?1
? 1
? ?f
?1 0
K
? ?x1
?
?f 2
?
+ (?1 ,K, ? N ) ? ? 0
?1
K
?x2
?
K
?K K K
?
?f N ?1
? 0
K
K
?
?x N ?1
?
?L1
?LN ?1
+K+
+
?x
?x
0?
?
0?
?
?
?V ( x N )
?
)+
0 ? = (0, 0K 0,
?x N
?
K?
?
?1 ??
?
+
? ?L1 ?L2
?
?f 2
?LN ?1 ? ? ?f 1
,
,K
, 0 ? + ? ?2
? ? 2 ,K , ?? N ? = 0;
? ?1 , ? 3
?
?
? ?x1 ?x2
?x2
?x N ? 1 ?? ?
?x1
?
?
+?
29
отсюда получаем последовательность значений множителя ?(i)
T
T
? ?f i ?
? ?Li ?
? (i ) = ?
? ? (i + 1) + ?
? , i = 0,1,K, N ? 1
? ?x(i ) ?
? ?x(i ) ?
?
?
?
?
с граничным условием
(5.13)
T
? ?V ( x N ) ?
?? .
? ( N ) = ??
? ?x( N ) ?
(5.14)
2-е необходимое условие
?H ?L
=
+ ?T
?u ?u
? ?f 0
0
?
? ?u0
?
?f 1
? 0
?
?u1
??
0
? 0
? K
K
?
?
0
?? 0
?
?t ?L0 ?L1
?LN ?1
=
+
+K +
+ (?1 , ? 2 ,K , ? N ) ?
?u ?u ?u
?u
?
0 K
0 ?
?
?
0 K
0 ?
? ? ?L0 ?L1
?LN ?1 ?
,
,K ,
?+
? = ??
?u N ?1 ??
0 ? ? ?u0 ?u1
K K
K K
K ?
?
?f N ?1 ?
K K
?
?u N ?1 ??
? ?f 0
?f 1
?f N ?1 ?
+ ? ?1
, ?2
,K , ? N
? = 0.
? ?u0
?
u
u
?
?
?
1
N
1
?
?
Отсюда получаем последовательность соотношений для определения значений оптимального управления u(i)
?Li
?f i
T
+ ? (i + 1)
= 0.
?u (i )
?u (i )
(5.15)
Третье необходимое условие приводит к соотношению (5.11) и, в итоге,
к исходной последовательности (5.9).
30
Переходя к пределу при N, стремящемся к бесконечности в соотношениях (5.13) – (5.15), получим решение исходной непрерывной задачи
для системы (5.7) с функционалом (5.8). Итак, для того чтобы найти
вектор управления u(t) , при котором критерий качества (5.8) достигает
минимума, нужно решить систему дифференциальных уравнений
(5.16)
x& = f ( x, u , t );
? ?f ?
? ?L ?
?& = ?? ? ? ? ? ? ,
? ?x ?
? ?x ?
T
T
5.17)
где u(t) определяется из условия
?L
?f
+ ?T
=0.
?u
?u
(5.18)
Граничные условия для уравнений (5.16), (5.17) разделены: одни из
них заданы при t = t0 , другие – при t = tk :
(5.19)
x(t0) – задано,
T
? ?V (t k ) ?
?? .
? (t k ) = ??
(
)
?
x
t
k ?
?
(5.20)
Отметим теперь, что соотношения (5.16) – (5.20) можно получить,
приравнивая нулю первую вариацию нового вспомогательного критерия
качества (лагранжиана), получаемого из (5.7), (5.8):
tk
I = V ( x k , t k ) + ? {L( x, u , t ) + ?T ? [ x& ? f ( x, u, t )]}dt ,
(5.21)
t0
а также введением для удобства вспомогательной скалярной функции
(гамильтониана)
(5.22)
H ( x , u , ? , t ) = L ( x , u , t ) + ?T ? f ( x , u , t ) .
Этот способ будет нами широко использоваться в дальнейшем.
31
6. Прогнозирующие системы оптимального
управления
Рассмотрим задачу оптимального управления, состоящую в отыскании траектории x(t)?Rn и управления u(t)?Rm, доставляющих минимум функционалу
tk
?
J = V(x(tk), tk)+ L( x, u , t ) dt
(6.1)
to
и удовлетворяющих уравнениям движения
x& = F(x, u, t),
причем различают случаи:
а) u = ? (непосредственный привод);
(6.2)
б) u = ?& (интегрирующий привод).
В соответствии с рецептом Лагранжа введем новый управляемый
лангранжиан
~
L (t, x, x& , ? , u, ? ) = L(x, u, t) + ? T( x& `- F(t, x, ? )),
(6.3)
При этом задача (6.1), (6.2) эквивалентна задаче на безусловный экстремум нового функционала
t
k
~
~
J =V(x(tk),tk) + ? L (t , x, x& , ? , u , ? )dt .
t0
~
Рассмотрим теперь вариацию критерия качества J , соответствующую вариациям вектора управления u(t), вектора состояния x(t), начального t0 и конечного tk моментов времени:
случай u = ? :
tk
?J% = ?[V ( x(tk ), tk ] + ?[ L% (t , x, x& , u , ?)dt ] =
?
t0
32
tk
t
= [Vx ?x + (Vx x& + Vt )?t ]tk + L% ? t tk + ( L% x ? x + L% x& ?x& + L%u ?u + L%? ??)dt.
0
?
t0
Возьмем интеграл от второго слагаемого подынтегрального выражения по частям
tk
?
& = L% x& ? x ttk ?
L% x& ?xdt
0
t0
tk
d
? dt ( L%x& )?xdt.
t0
В результате получим
?J% = [(Vx + L% x& )?x + (Vx x& + Vt + L% )?t ]tk ? ( L% ?t + L% x ?x)t0 +
tk
+ ? [{L% x ?
t0
d %
( Lx& ))?x + L%u ?u + L%? ?? ]dt = 0,
dt
так как для экстремума ?J% = 0 .
Таким образом, гладкие экстремали должны удовлетворять уравнениям Эйлера
d %
( Lx& ) ? L%x = 0;
dt
L%? = 0; L%u = 0
и граничным условиям:
((Vx + L% x& )?x + (Vx x& + Vt + L% )?t )tk ? ( L% ?t + L% x ?x )t0 = 0 ,
(6.4)
(6.5)
где индексы tk и t0 означают, что выражения, стоящие в скобках, необходимо вычислить в конечной и начальной точках траектории; ?t – возможные вариации концов интервала [t0, tk]; ?x – вариация траектории x,
определенная, как и x на интервале, содержащем [t0, tk] с некоторой его
окрестностью.
~
При подстановке явного выражения для L (6.3) в (6.4), получим
?& T + ?T Fx ? Lx = 0;
x& ? F = 0;
Lu ? ?T Fu = 0.
(6.6)
33
Соотношения, аналогичные (6.4), (6.6) в случае интегрирующего привода u = ?& , имеют вид
d %
( Lx& ) ? L% x = 0 ;
?& T + ?T Fx ? Lx = 0 ;
dt
L%? = 0 ;
x& ? F = 0 ;
или
d %
( L & ) ? L%? = 0 ;
dt ?
(6.7)
d
( Lu ) + ?T F? = 0 .
dt
Рассмотрим теперь управляемый гамильтониан, отвечающий лаг-
~
ранжиану L .
а) случай u = ? :
H (t , x, ?, u ) = ?T x& ? L% (t , x, x& , u , ? ) = ?T F ? L ( x, u , t ).
Уравнения Эйлера переписываются следующим образом:
(6.8)
x& = H ?Т ,
x& = F ;
T
T
?& T = ? H x ; или ?& = ?? Fx + Lx .
(6.9)
Оптимальное управление u = u0 находится из соотношения:
H (t , x, ?, u0 ) = max H (t , x, ?, u ) .
u?U
(6.10)
Если U – открытая область, то вместо (6.10) можно записать
(6.11)
Hu= 0 или ?T Fu ? Lu = 0 .
Часто бывает полезно так называемое уравнение Гамильтона-Якоби. Зададимся некоторым управлением u(t). На характеристиках, т. е.
решениях системы Гамильтона (6.9), отвечающих выбранному u(t), комбинация
?T dx ? Hdt
оказывается полным дифференциалом, т. е. существует функция
S (t , x) ? R 1 , такая, что
(6.12)
dS = ?T dx ? Hdt , т.е. S x = ?T , St = ? H .
Последнее равенство в (6.12) дает уравнение Гамильтона-Якоби:
34
S t + H (t , x, S x , u ) = 0 .
(6.13)
б) случай u = ?& (интегрирующий привод):
H (t , x, ?, ?, w) = ?T x& + w?& ? L% (t , x, x& , ?, ?& , ? ) =
= ?T F (t , x, ?) + w?& ? L(t , x, ?& ),
(6.14)
~
здесь w = L?& = Lu .
Уравнения Гамильтона примут вид
x& = H ?T ;
x& = F ;
?& T = ? H x ; или ?& T = ??T Fx + Lx ;
?& = H w ;
?& = u;
w& = ? H ?;
w& = ??T F? .
(6.15)
На характеристиках системы (6.15) комбинация
(6.16)
?T dx + wd ? ? Hdt = dS
является полным дифференциалом скалярной функции S(t, x, ? ) , так
что
S x = ? T ; S ? = w ; St = ? H .
Отсюда получаем уравнение Гамильтона-Якоби
(6.17)
(6.18)
St + H (t , x, S x , ?, S? ) = 0 .
Рассмотрим систему с линейным управлением, квадратичным функционалом и прямым приводом.
Для последнего десятилетия характерно бурное развитие общей теории оптимального управления, базирующейся на рассмотрении пространства состояний и функционалов качес
Документ
Категория
Без категории
Просмотров
0
Размер файла
347 Кб
Теги
0058, 2000
1/--страниц
Пожаловаться на содержимое документа