Том XXXVI
УЧЕНЫЕ ЗАПИСКИ ЦАГИ
2 00 5
№ 3 — 4
УДК 629.735.33.016
ОСОБЕННОСТИ РАСЧЕТА ТРАЕКТОРИЙ БОЛЬШОЙ ДАЛЬНОСТИ В ЗАДАЧЕ О ПЕРЕЛЕТЕ САМОЛЕТА С МИНИМАЛЬНЫМ РАСХОДОМ ТОПЛИВА
А. В. ТЕГИН
Рассматриваются особенности численного решения задачи оптимизации траекторий самолета по критерию минимума расхода топлива применительно к большим значениям конечной дальности полета. Описана разработанная на основе принципа максимума и энергетического подхода итерационная процедура для расчета оптимального программного движения самолета при полете на большую дальность в вертикальной плоскости.
Одной из важных задач современных разработок в области гражданской авиации является повышение экономичности воздушных перевозок, которая определяется в основном их топливной эффективностью. Методы оптимального траекторного управления позволяют решать эту задачу, однако их практическое применение требует разработки эффективных вычислительных схем, которые мало освещены в известных публикациях [1] — [4].
Для разработки эффективных алгоритмов планирования полетов и их реализации в бортовых системах управления оказалось недостаточным использование теории оптимизации динамических систем на основе принципа максимума Л. С. Понтрягина. Развитие получило относительно новое направление в математической теории оптимального управления применительно к задачам управления полетами самолетов на транспортных режимах — «метод сингулярных возмущений» [4]. При этом открываются возможности для практической реализации близкого
к оптимальному управления полетом, обеспечивающего при использовании имеющихся аэродинамических и расходно-тяговых возможностей самолетов наибольшую экономию топлива.
Система дифференциальных уравнений движения центра масс самолета в вертикальной плоскости имеет 5-й порядок. Если задача оптимального управления решается с использованием принципа максимума, возникает двухточечная краевая задача 10-го порядка (кроме исходной появляется сопряженная система уравнений 5-го порядка). Она представляет трудности для решения, в том числе по причине плохой численной обусловленности, что отмечалось и ранее (см. [4]). Поэтому необходимо по возможности упростить методику расчета оптимальной траектории.
Для этого, в частности, используется энергетический подход (см. [1], [4]), который позволяет уменьшить число фазовых уравнений до трех. При использовании энергетического подхода принимается, что угол наклона траектории и ее кривизна малы. Это обычно допустимо для транспортных траекторий. Для решения рассматриваемой проблемы данный подход в сочетании с упоминавшимся выше методом сингулярных возмущений является конструктивным. В этом методе переменные динамической системы разбиваются на «быстрые» и «медленные» и в соответствии с таким разбиением производится декомпозиция задачи. Метод служит не только для упрощения расчета траекторий полета, его применение является условием успешного решения задачи в случае траекторий достаточно большой дальности.
В статье рассматривается процедура решения задачи оптимизации траектории самолета в вертикальной плоскости по критерию минимума расхода топлива при больших значениях дальности полета и свободном времени. Вопросы реализации программных траекторий, получаемых при решении подобного рода задач, рассмотрены в работе [5] — [7].
1. Постановка задачи. Задача перелета с минимальным расходом топлива на заданную дальность формулируется в рамках энергетического подхода. Предполагается, что маневрирование по курсу незначительно и слабо влияет на оптимальный профиль полета. Углом наклона траектории и его производной можно пренебречь, при этом нормальная скоростная перегрузка полагается равной 1. Движение самолета в вертикальной плоскости описывается следующими уравнениями:
Е = пхУ, х = V, т = -/ (И, V, Т). (1)
V 2
Здесь Е = И +----— удельная механическая энергия; V — скорость; И — высота полета;
2 Я
Я — ускорение свободного падения; пх — тангенциальная перегрузка; х — продольная координата; т — масса самолета; / — расход топлива в единицу времени; Т — тяга силовой установки.
Требуется найти такую траекторию, чтобы при заданных начальных
Е ( 0 ) = Е0, х ( 0 ) = х0, т ( 0 ) = т0
и конечных условиях
Е (/ ) = Е/, хЫ =1
обеспечивался минимальный расход топлива
1/
3 =| /йх. (2)
0
Нижний индекс / обозначает конечный (незаданный) момент времени.
В качестве управляющих переменных будем рассматривать высоту полета И и тягу силовой установки Т, меняющиеся в пределах ограничений:
Итт ( Е )< И < Иmax (Е),
Ттт (И Е )< Т < Tmax (К Е ).
Полученная в результате преобразования
И = И - Иmin (Е)
Иmax (Е) Итт ( Е )
переменная И уже не ограничивается зависящими от фазовой переменной Е функциями, а изменяется на фиксированном интервале:
0 < И < 1.
2. Условия оптимальности. Рассматриваемая оптимизационная задача решается с использованием процедуры принципа максимума. В этом случае гамильтониан выглядит как
H = / + ХхУ + ХЕпхУ + V! (Т - ^ (И, Е)) + У2 (Тпш (И, Е) - Т), (3)
где V!, V2 — множители Лагранжа, соответствующие ограничениям на максимальную Т^ах и минимальную ТшП величины тяги силовой установки; ХЕ, Xх — сопряженные переменные, соответствующие переменным состояния Е, х и удовлетворяющие условиям
д H
д H
к E =-------------, к x = —
Лі ї—т ■'Л
дE ox
(4)
При решении рассматриваемой здесь задачи обычно допускается, что изменение массы самолета не влияет на функционал (2) (см., например, [8]).
В связи с тем, что фазовая координата x не входит в гамильтониан, переменная кх сохраняет на траектории постоянное значение:
к x = const.
Согласно принципу максимума функция Гамильтона (3) на оптимальной траектории достигает минимума по управляющим переменным и при свободном (незаданном) времени полета удовлетворяет условию
minH = 0.
h, T
(5)
В расчетах по оптимизации траекторий, проводившихся автором, в аппроксимированном виде использовались аэродинамические и расходно-тяговые характеристики
среднемагистрального пассажирского самолета типа Ту-204. Нахождение оптимальных значений управляющих переменных осуществлялось с помощью численной процедуры. Для поиска минимума гамильтониана использовалась двумерная (по двум переменным управления) итерационная процедура, основанная на методе «золотого сечения» [9]. Точность поиска подбиралась такой, чтобы ее повышение не приводило к заметному изменению траектории.
Величина сопряженной переменной ХЕ выбирается из условия равенства нулю гамильтониана в начальной или конечной точке (в зависимости от направления интегрирования). В дальнейшем эта величина используется как краевое условие при интегрировании первого из уравнений (4).
3. Исследование чувствительности решения задачи к конечным условиям. На рис. 1 показано изменение конечной дальности полета в зависимости от переменной Xх. Начальная и конечная величины удельной энергии выбирались равными 910 м, что соответствует высоте полета 400 м и скорости 100 м/с. Начальная масса самолета составляла 72 т.
Кривая, приведенная на рис. 1, имеет существенно нелинейный характер, и необходимая точность задания переменной Хх с ростом конечной дальности полета также нелинейно
Рис. 1. Зависимость конечной дальности от параметра Xх при «сквозной»
оптимизации
увеличивается. Эту точность можно оценить, исходя из потребной точности расчета конечной дальности Ь и производной дЬ/дХх для этой дальности. Оценка производной дЬ/дХх показывает, что при увеличении дальности полета в диапазоне от 0 до 500 км ее значение возрастает приблизительно с 104 км2/т до 105 км2/т, в диапазоне дальностей от 500 до 650 км — до 109 км2/т, в диапазоне от 650 до 950 км — до 3 • 1012 км2/т и в диапазоне от 950 км до 1150 км — до
13 2
5 • 10 км /т. Численные исследования при выполнении настоящей работы проводились с использованием написанной на алгоритмическом языке Си++ вычислительной программы, в которой применялся формат двойной точности для всех переменных с плавающей запятой (точность мантиссы переменных данного типа составляет 15 значащих цифр). Таким образом, для реализации в расчетах траекторий с конечной дальностью полета 1100—1200 км, уже необходимо применять точность задания сопряженной переменной Хх, близкую к максимальной. Налицо высокая чувствительность решения задачи к конечным условиям.
Вместе с тем видно, что даже использование самой высокой точности при вычислении величины Хх не решает проблемы расчета траектории для больших значений конечной дальности, так как оптимизация без декомпозиции профиля полета на отдельные участки («сквозная» оптимизация) обеспечивает удовлетворительные результаты решения краевой задачи только для траекторий с конечной дальностью полета не более 1200 км. При этом возможности самолета, математическая модель которого применялась в расчетах, позволяют совершать перелеты с дальностью 5000—6000 км. Поэтому, чтобы получать в результате численной оптимизации траектории с большей, чем это возможно при «сквозном» расчете, дальностью, необходимо применение метода сингулярных возмущений.
4. Выбор методики расчета траекторий самолета с большой конечной дальностью. В процедуре оптимизации траекторий большой дальности используется их декомпозиция. Стандартный профиль полета самолета содержит три участка: набор высоты, участок крейсерского режима и снижение.
Особенностью крейсерского режима полета является его вырожденность, которая выражается в отсутствии уравнения (1) для удельной энергии в исходной и первого из уравнений (4) для ХЕ в сопряженной системе уравнений. Таким образом, уменьшается порядок как исходной, так и сопряженной систем уравнений. В терминах метода сингулярных возмущений крейсерский, или стационарный, режим полета является внешним решением задачи, а примыкающие к нему слева и справа нестационарные режимы — внутренним решением.
В работе [3] показано, что для оптимизации крейсерского участка полета можно пользоваться следующей формулой:
. • /
X х =-mm —
х к, V V
0. (6)
пх =0,
Пу =1
Условие (6) отражает известный факт, что решением задачи оптимизации крейсерского режима полета самолета с незаданным временем является полет с минимальным километровым расходом топлива. Сопряженная переменная Хх при переходе от внутреннего к внешнему решению и обратно не претерпевает разрывов, т. е. она имеет постоянную на всех участках величину [10]. Для траекторий с декомпозицией, содержащих крейсерский участок, она может быть вычислена согласно (6) и применяться для расчета нестационарных режимов полета.
При декомпозиции задачи расчет параметров каждого участка полета производится отдельно, но так, чтобы вместе они составляли единый профиль полета, удовлетворяющий граничным условиям. Выполнить такой расчет можно различными способами, приведенными ниже.
Один из способов, использовавшийся ранее, в том числе и автором, заключается в подборе такого значения Xх, которое обеспечило бы плавный переход от участка набора высоты к крейсерскому режиму полета. Переход происходит на дальности, где параметры участков набора
высоты и крейсерского режима совпадают с заданной точностью. Отметим, что при данном способе расчета следует производить одновременный расчет двух участков полета.
Предлагаемый альтернативный способ расчета заключается в следующем. Оптимизация профиля полета начинается в точке перехода от участка набора высоты к крейсерскому участку. Расчет участка набора высоты осуществляется в обратном времени. Значение Хх вычисляется
с использованием формулы (6). Условие на массу самолета в начальной точке удовлетворяется методом последовательных приближений. Итерационный процесс сходится, как правило, за 2 —4 итерации.
Последовательность этапов расчета оптимальной траектории, при которой сначала определяется внешнее решение задачи, а затем с ним сопрягаются полученные отдельно участки внутреннего решения, является фундаментальной особенностью метода сингулярных возмущений и объясняет второе его название «метод сращиваемых асимптотических разложений».
На взгляд автора, такой подход к разработке процедуры расчета профилей полета является по ряду причин предпочтительным. Во-первых, при определении в начале расчета траектории параметров внешнего решения с помощью формулы (6) мы имеем возможность достаточно легко вычислить величину Х х. При другом подходе, упоминавшемся выше, для выхода в крейсерские условия необходимо подобрать значение параметра Хх с высокой точностью. Во-вторых, нет необходимости вычислять параметры крейсерского участка параллельно с расчетом набора высоты.
5. Процедура расчета траекторий с декомпозицией профиля полета. Здесь описана разработанная автором итерационная процедура оптимизации траекторий с большой заданной дальностью полета.
В данной процедуре расчет участка набора высоты выполняется в обратном времени, а участков стационарного полета и снижения — в прямом. Точки перехода между участками подбираются итерационно.
Расчет начинается в точке перехода из набора высоты на крейсерский участок. В соответствии с условиями пх = 0, Пу = 1 ив предположении малости угла атаки тяга на
стационарном режиме равна аэродинамическому сопротивлению, а подъемная сила — весу самолета. Таким образом, крейсерская тяга определяется весом самолета, и так как секундный расход топлива / является функцией тяги, высоты и скорости полета, то для нахождения оптимальных параметров крейсерского режима достаточно при известной массе самолета применить формулу (6).
Перед началом процедуры выбирается нулевое приближение для значения массы самолета, при которой происходит переход с участка набора высоты на крейсерский участок. Далее выполняется оптимизация крейсерских скорости и высоты полета согласно (6). В результате нахождения оптимальных параметров вычисляется величина Хх и производится переход к
расчету участка набора высоты. Следует упомянуть о следующей особенности данного перехода.
Оптимизация нестационарного и стационарного режимов полета при значении удельной энергии стационарного режима дает одни и те же значения управляющих переменных. При нулевом значении тангенциальной перегрузки, полученном для нестационарного режима, в соответствии с уравнением (1) самолет остается на крейсерском уровне энергии. Выход в условия нестационарного режима полета обеспечивается за счет отрицательного приращения к значению удельной энергии, полученному в результате оптимизации скорости и высоты согласно (6).
При численном решении задачи для рассмотренного типа самолета абсолютная величина приращения принималась равной 100 м. Выбор указанного приращения определяется двумя факторами. С одной стороны, чем меньше его величина, тем продолжительнее участок набора высоты. В то же время резкий переход от расчета одного участка программы к другому увеличивает промежуточный участок сглаживания при отслеживании программной траектории с использованием полных уравнений движения самолета в вертикальной плоскости. Вопросы сопряжения участков программы полета, в числе других аспектов отслеживания траекторий, рассмотрены в статье [7].
Интегрирование уравнений движения для участка набора высоты производится в обратном времени до выполнения начального условия на величину удельной энергии. Начальное значение массы удовлетворяется итерационно. Дальность участка наборы высоты при этом остается свободной.
По результатам каждой итерации вычисляется Ьт0) — разность между заданным значением начальной массы самолета и полученным в результате расчета набора высоты; назовем эту величину невязкой начальной массы самолета. Значение массы самолета в момент перехода от набора высоты к крейсерскому участку (в обратном времени) тнв — кр на /-й итерации
корректируется на величину невязки:
т^ —т(/ 1) + 8т(/ 1)
н.в ——кр н.в — кр 0
Т аблица 1
№ итерации Дальность участка, км Время набора высоты, с Начальная масса самолета, т Конечная масса на наборе высоты, т
0 272.51 1260 87.1301 85
1 259.46 1205 89.9739 87.870
2 259.40 1205 90.0007 87.896
В табл. 1 приведены результаты процедуры расчета участка набора высоты. Начальная масса самолета была задана равной 90 т, удельная энергия соответствовала начальной высоте полета 400 м, начальной скорости — 100 м/с. На нулевой итерации значение массы самолета в момент перехода от наборы высоты к крейсерскому участку назначалось равным 85 т. По условиям процедуры итерационный цикл завершался, если отличие начальной массы от заданной не превышало 1 кг. Как видно из таблицы, в этом случае заданная точность достигнута за 2 итерации.
Крейсерский участок полета заканчивается в тот момент, когда разность между заданной конечной и текущей дальностью равна величине, необходимой для последующего снижения и выхода в конечные условия.
Момент перехода от крейсерского участка к снижению может быть определен итерационным путем. На начальной итерации дальность полета, на которой начинается снижение, выбирается достаточно произвольно. Для самолета рассматриваемого типа в качестве начального приближения дальность участка снижения можно принять равной 300 км. Интегрирование траектории происходит до выполнения конечного условия на величину удельной энергии. По окончании каждой итерации полученная дальность полета сравнивается с заданной и вычисляется невязка конечной дальности:
Ы1 — Ь - х/г.
С помощью итерационной процедуры абсолютная величина невязки уменьшается, пока заданная дальность полета не будет достигнута с необходимой точностью. При этом значение дальности Ьсх, при которой нужно сходить с крейсерского участка, корректируется на (/ + 1)-й итерации на величину невязки, вычисленной на /-й итерации:
Ь х+1} —4Х +м.
На рис. 2 показана траектория с дальностью полета 1500 км, для которой момент перехода на снижение определяется по описанной методике. В данном случае условие на конечную дальность было удовлетворено за 2 итерации. Результаты расчета данного участка снижения приводятся в табл. 2.
1 ©
\V 1 J
\\ !
1 ®Vfr ]
I k\ 1
I . I \— ; ! \.\]
1000 1100 1200 1300 1400 1500 L. км
Рис. 2. Участки снижения в итерационном процессе удовлетворения условия
на дальность полета
Т аблица 2
№ итерации Дальность участка, км Время снижения, с Расход топлива на снижение, т Средний расход топлива на единицу пути, кг/км Невязка конечной дальности на предыдущей итерации, км
1 466.586 2279.98 1.2092 2.592 113.891
2 296.660 1537.69 0.6143 2.071 57.452
В статье описана итерационная процедура для расчета оптимальных, в смысле минимума расхода топлива, траекторий с большой заданной дальностью полета. Процедура используется для расчета траекторий, содержащих три участка полета: набор высоты, крейсерский участок и снижение. При этом на участке набора высоты корректируется невязка начального значения массы самолета, на крейсерском участке и снижении — невязка конечной дальности путем выбора точки начала снижения. Опыт расчетов с использованием процедуры показывает ее эффективность по затратам машинного времени.
ЛИТЕРАТУРА
1. Zagalsky N. R., Irons R. P., Schultz R. L. Energy state approximation and min-imum-fuel fixed range trajectories // J. Aircraft. — 1971. Vol. 8, N 6.
2. Илларионов В. Ф., Пашинцев В. Т. Качественный анализ семейства оптимальных траекторий в задаче полета самолета на максимальную дальность //Труды ЦАГИ. — 1974. Вып. 1591.
3. Гревцов Н. М., Ефимов О. Е., Мельц И. О., Трубецкой А. Б. Соотношение условий оптимальности стационарного и нестационарного режимов полета в методе сингулярных возмущений // Ученые записки ЦАГИ // 1995. Т. XXVI, № 1 — 2.
4. Chacravarty A., Vagners J. Application of singular perturbation theory to onboard aircraft trajectory optimization // AIAA 19th Aerospace Sciences Meeting, St., Louis, Missouri. — January 12 — 15, 1981.
5. T e g i n A. V. 4D navigation problem — features of organization of numerical solution procedure / Proceedings of IV Seminar on Recent Research and Design Progress in Aeronautical Engineering. — Warsaw. — 2000.
6. Тегин А. В. Вычислительная процедура решения оптимизационной задачи о перелете самолета на заданные дальность и время / Сб. трудов 6-го Международного симпозиума «Авиационные технологии 21 века: новые рубежи авиационной науки». —
Жуковский. — 14 — 17 августа 2001.
7. Гревцов Н. М., Тегин А. В. Формирование управления самолетом для отслеживания траектории в задаче 4-мерной навигации // Ученые записки ЦАГИ. — 2000.
Т. XXXI, № 1 — 2.
8. Erzberger Н., Lee H. Constrained optimum trajectories with specified range //
J. Guidance and Control. — 1980. Vol. 3, N 1.
9. Полак Э. Численные методы оптимизации. — М.: Мир. — 1974.
10. Брайсон А., Хо Ю - Ши. Прикладная теория оптимального управления. — М.: Мир. — 1972.
Рукопись поступила 5/12004 г.