Математическое моделирование элементов технологии гиперзвукового полета тема автореферата и диссертации по механике, 01.02.05 ВАК РФ
Латыпов, Альберт Фатхиевич
АВТОР
|
||||
доктора физико-математических наук
УЧЕНАЯ СТЕПЕНЬ
|
||||
Новосибирск
МЕСТО ЗАЩИТЫ
|
||||
2009
ГОД ЗАЩИТЫ
|
|
01.02.05
КОД ВАК РФ
|
||
|
На правах рукописи
2 7 АЫ 2В0а
Латыпов Альберт Фатхиевич
Математическое моделирование элементов технологии гиперзвукового полета
01.02.05 - Механика жидкости, газа и плазмы
Автореферат диссертации на соискание ученой степени доктора физико-математических наук
Новосибирск - 2009
003475853
Работа выполнена в Институте теоретической и прикладной механики им. С.А. Христиановича СО РАН.
Официальные оппоненты:
академик РАН, доктор физико-математических наук
профессор Левин В.А.
доктор физико-математических наук
профессор Иванов М.С.
доктор технических наук
Серманов В.Н.
Ведущая организация: Экспериментальный машиностроительный завод им. В.М. Мясищева, г. Жуковский.
Защита состоится 13 ноября 2009 года в 10 часов на заседании диссертационного совета Д 003.035.02 по защите диссертаций на соискание ученой степени доктора наук при Институте теоретической и прикладной механики им. С.А. Христиановича по адресу: 630090, Новосибирск, ул. Институтская, д. 4/1.
Отзыв на автореферат в двух экземплярах, заверенный печатью учреждения, просьба направлять на имя ученого секретаря диссертационного совета Д 003.035.02.
С диссертацией можно ознакомиться в библиотеке Института теоретической и прикладной механики им. С.А. Христиановича СО РАН.
Автореферат разослан ^ ^ 2209 2009г.
Учёный секретарь диссертационного совета
доктор технических наук Засыпкин И.М.
Общая характеристика работы
Актуальность. Для создания воздушно-космического самолета (ВКС) как средства транспортировки грузов на трассе «земля - околоземная орбита - земля», несмотря на значительные достижения, необходимо решить множество проблем. Первыми в этом ряду являются взаимосвязанные проблемы прямоточного воздушно-реактивного двигателя (ПВРД) и аэродинамического качества конфигурации летательного аппарата (ЛА), т.к. характеристики силовой установки и аэродинамические характеристики компоновки J1A определяют главным образом необходимые затраты топлива для выведения на орбиту. В настоящее время значительное внимание уделяется решению задачи активного управления обтеканием тел посредством энергетического и/или силового воздействия на набегающий поток, в частности, посредством подвода тепла перед телом в сверхзвуковом потоке. Изучению этой проблемы посвящено значительное число работ. Необходима оценка эффективности такого способа управления обтеканием тел. Его целесообразность может быть установлена методом функционального моделирования. В соответствии с методом JIA, как сложная иерархическая система, должна быть представлена в виде взаимосвязанной совокупности подсистем, определяемых по функциональным признакам. Необходимо построить функциональные математические модели подсистем JIA (аэродинамика, силовая установка, траектория полета) и аппарата в целом. Сравнение расходов топлива на разгон с нагревом воздуха перед JIA и без нагрева должно производиться на оптимальных траекториях полета, которые описываются системой обыкновенных дифференциальных уравнений (ОДУ) с функциями управления. При численном решении задач оптимизации необходимо многократное вычисление функционала. Поэтому разработка алгоритмов с возможно малым числом операций всегда будет важной. Наиболее ресурсоёмкими по количеству операций являются задачи решения систем ОДУ. И чем более адекватна математическая модель физическому процессу, тем более жесткими являются входящие в её состав дифференциальные уравнения, что указывает на важность совершенствования методов их решения. Правые части системы ОДУ, описывающей траекторию полета ВКС, являются разрывными функциями. Существует также проблема расчета глобальной ошибки при численном интегрировании систем ОДУ. Поскольку основой большинства математических моделей являются дифференциальные уравнения и их системы, совершенствование имеющихся численных методов интегрирования систем ОДУ и разработка новых с учетом роста возможностей вычислительных средств была, есть и, вероятно, долго будет актуальной задачей.
Для увеличения эффективного удельного импульса комбинированной силовой установки необходимо увеличить число Маха полета, при котором возможна работа ПВРД. Отсутствуют достоверные экспериментальные результаты, свидетельствующие о сохранении сверхзвукового течения в канале при подводе энергии с эквивалентным коэффициентом избытка воздуха а = 1, тем более, при ограничении статической температуры продуктов сгорания. Это ограничение важно при гиперзвуковых числах Маха полета и связано с ограничением степени диссоциации продуктов сгорания, т.к. диссоциация уменьшает эксергию потока газа. Необходимо определить условия, при которых было бы возможным организовать подвод тепла с учетом названных факторов.
Экспериментальные исследования аэрогазодинамических характеристик моделей гиперзвуковых летательных аппаратов часто проводятся в аэродинамических трубах кратковременного действия. При этом для определения силовых характеристик используются тензометрические весы. Требуется разработка метода для восстановления переменных во времени действующих на модель сил и моментов.
/
Цели работы:
- исследование нестационарных течений в модельном канале ПВРД при импульсно-периодическом подводе энергии, определение условий формирования структуры течений;
- разработка эксергетического метода анализа и оценки характеристик ПВРД;
- разработка новых эффективных численных методов интегрирования систем обыкновенных дифференциальных уравнений на основе многозвенных интерполяционных полиномов Эрмита, их программная реализация;
- разработка функциональных математических моделей элементов гиперзвукового ЛА и оценка эффективности подвода энергии в набегающий поток;
- разработка метода восстановления действующих сил и моментов при испытаниях аэродинамических моделей и моделей ПВРД в импульсных трубах кратковременного действия.
Теоретическое значение и научная новизна работы определяются следующим:
1. Выполнен анализ квазиодномерного и двухмерного квазистационарных течений в канале переменного сечения, описываемых уравнениями Эйлера и формирующихся при импульсно-периодическом энергетическом воздействии при больших числах Струхаля. Получено, что при этом устанавливается периодический режим течения с малыми амплитудами колебаний параметров. Течения устойчивы в среднем на периоде. Получены условия, определяющие структуру течений.
Предложена конфигурация канала, в котором подвод тепла к сверхзвуковому потоку осуществляется с учетом ограничения статической температуры газа. Импульсно-периодический подвод энергии в таком канале позволяет увеличить число Маха полета до значений, при которых возможно использование прямоточного канала в составе комбинированного двигателя для увеличения эффективного удельного импульса.
В результате численного моделирования нестационарного двухмерного течения в канале переменного сечения при подводе тепловой энергии в локальных зонах в импульсно-периодическом режиме получена экспериментально наблюдаемая перестройка начального сверхзвукового течения, определяемая условием подвода заданного количества энергии.
2. Разработаны новые численные методы для решения следующих задач:
2.1. Построено семейство А-, Ь- и £(8)-устойчивых методов решения задачи Коши для жестких систем обыкновенных дифференциальных уравнений, основанные на представлении правых частей системы на шаге к в виде трёх точечных интерполяционных полиномов Эрмита: ЬМЯ{Ь,М, Л,я) - алгоритмы. Определена погрешность методов. Дано определение /.((^-устойчивости одношаговых методов с малым параметром 8. Разработан алгоритм расчета глобальной ошибки и алгоритм решения задачи Коши для систем ОДУ с разрывными правыми частями.
2.2. Разработан метод решения системы линейных интегральных уравнений Вольтерра 1-го рода с разностным аргументом для варианта, когда исходная информация (значения измеряемых функций и ядра) заданы в дискретных точках с известной ошибкой. Решение определяется в классе кусочно-постоянных и кусочно-линейных функций с использованием условия равенства нулю средних значений невязок уравнений на интервалах, на которые разбивается область определения решения. Число интервалов и распределение их длин определяются посредством минимизации среднеквадратичной невязки уравнений.
2.3. Разработан метод восстановления действующих нагрузок в классе кусочно-постоянных функций при испытаниях моделей в аэродинамических трубах кратковременного действия. Приведены примеры решения задач по определению
аэродинамических характеристик эталонной модели НВ-2, демонстратора Ares и модели Expert по результатам испытаний в аэродинамической трубе АТ-303 ИТПМ СО РАН.
2.4. Для решения задач математического программирования методом штрафных функций из условия локального минимума вспомогательного функционала получена оценка для коэффициента штрафа к ~ £.
2.5. Созданы комплексы программ для решения названных классов задач.
3. Разработан эксергетический метод оценки характеристик и анализа ПВРД. Для графического отображения возможных схем подвода тепла в канале ПВРД предложена диаграмма в координатах "полная температура - эксергия". Получено выражение для изменения эксергии в термодинамической системе при подводе тепла и наличии необратимых процессов.
4. Разработана методика оценки эффективности подвода тепла перед JIA при полете со сверхзвуковой скоростью. Полет происходит на границе раздела сред различной плотности (режим глиссирующего полета). Показана значительная эффективность такого способа управления обтеканием JIA как при крейсерском полете, так и при полете с ускорением.
Методика исследований. Проведенные исследования опираются на численные методы механики сплошной среды, методы вычислительной математики, методы условной численной минимизации функционалов.
Обоснованность и достоверность полученных результатов обеспечена использованием известных моделей механики жидкости и газа, методов вычислительной математики, приведением достоверных оценок точности разработанных методов, доказательствами необходимых теорем. Выполнено также тестирование разработанных методов, сравнение полученных результатов с результатами известных методов, с точными решениями и экспериментальными данными.
Практическая значимость работы. Разработанные методы и результаты могут найти применение при: исследованиях динамики газовых потоков при импульсно-периодических воздействиях при больших числах Струхаля; разработке ПВРД для гиперзвуковых чисел Маха полета; разработке активных способов управления аэродинамическими характеристиками летательных аппаратов; интерпретации экспериментальных результатов аэрофизических исследований; оценке эффективности термодинамических систем; решении задач динамики систем и оптимального управления; обратных задач и идентификации параметров; разработке новых численных методов решения уравнений механики жидкости и газа.
Апробация работы. Результаты работы по мере их получения докладывались на следующих российских и международных конференциях:
2-я Международная школа по моделям механики сплошной среды, Владивосток, 1991; Международная конференция "Задачи со свободными границами в механике сплошной среды", Новосибирск, 1991; Всероссийская школа-семинар по комплексам программ математической физики, Новосибирск, 1992; XI Международная конференция по автоматически пилотируемым летательным аппаратам, Англия, Бристоль, 1994; International Aerospace Congress, Moscow, 1994; Международная конференция "Исследование гиперзвуковых течений и гиперзвуковые технологии", Жуковский, 1994; International Congress on Instrumentation in Aerospace Simulation Facilities, Ohio, USA, 1995; AIAA Sixth International Aerospace Planes and Hypersonic Technologies Conference, Chattanooga, 1995; AIAA 8th Intern. Space Planes and Hypersonic Systems and Technologies Conference. Norfolk, USA, 1998; Международная конференция "Фундаментальные исследования для гиперзвуковых технологий", Жуковский, 1998; Всероссийская научная конференция "Краевые задачи и их
приложения", Казань, 1999; Конференция "Юбилейные Чаплыгинские чтения", Новосибирск, 1999; Вторая Международная конференция "Устойчивость и управление для нелинейных трансформирующихся систем", Москва, 2000; The 3rd, 4th, 5th Workshop on Magneto-Plasma-Aerodynamics in Aerospace Applications, Moscow, 2001, 2002, 2003; Международная конференция "Математические модели и методы их исследования", Красноярск, Институт вычислительного моделирования, 2001; IV Международная конференция по неравновесным процессам в соплах и струях, Санкт-Петербург, 2002; International Federation of Automatic Control Workshop - Modeling and Analysis of Logic Controlled Dynamic Systems, Irkutsk, 2003; III Всероссийская конференция "Математика, информация, управление", Иркутск, 2004; European Conference for Aerospace Sciences (EUCASS). Moscow, 2005; II-XIV International Conference on the Methods of Aerophysical Research, Novosibirsk, 1990-2008; West East High Speed Flow Field Conference (WEHSFF), Moscow, 2007; Международная конференция "Обратные и некорректные задачи математической физики", посвященная 75-летию академика М.М. Лаврентьева, Новосибирск, 2007; Седьмой, восьмой и девятый всероссийские съезды по теоретической и прикладной механике; Семинары ИТПМ СО РАН им. С.А. Христиановича.
Публикации. Основное содержание диссертации опубликовано в 43 научных работах.
Личный вклад автора. Автор являлся ведущим разработчиком всех представленных направлений исследований. Из совместных публикаций в диссертацию включены результаты, полученные автором или при его непосредственном участии. Содержание диссертации и автореферата обсуждено и согласовано с соавторами.
Структура и объем работы. Диссертация состоит из введения, четырёх глав, выводов, списка цитируемой литературы из 151 наименований и 61 рисунков. Общий объем работы 199 страниц, включая рисунки.
На защиту выносятся следующие научные результаты:
1. Выполнен анализ квазиодномерного и двухмерного квазистационарных течений в канале переменного сечения, описываемых уравнениями Эйлера и формирующихся при импульсно-периодическом энергетическом воздействии при больших числах Струхаля. Течения устойчивы в среднем на периоде. Получены условия, определяющие структуру течений: максимально допустимое значение энтропии для каждого сечения канала и условие перехода через скорость звука в квазиодномерном случае при подводе энергии и наличии диссипации кинетической энергии.
Получено, что при больших значениях числа Струхаля устанавливается периодический режим течения с малыми амплитудами колебаний параметров. Так как при этом энергия подводится при постоянном объеме, то этот режим обеспечивает максимальное значение эксергии потока и, следовательно, тяги двигателя.
Предложена конфигурация канала, в котором подвод тепла к сверхзвуковому потоку осуществляется с учетом ограничения статической температуры газа. Импульсно-периодический подвод энергии в таком канале позволяет увеличить число Маха полета до значений, при которых возможно использование прямоточного канала в составе комбинированного двигателя для увеличения эффективного удельного импульса.
В результате численного моделирования нестационарного двухмерного течения в канале переменного сечения при подводе тепловой энергии в локальных зонах в импульсно-периодическом режиме получена экспериментально наблюдаемая перестройка начального сверхзвукового течения, определяемая условием подвода заданного количества энергии.
2. Разработаны новые численные методы для решения следующих задач:
2.1. Построено семейство A-, L- и £(5)-устойчивых методов решения задачи Коши для жестких систем обыкновенных дифференциальных уравнений (ОДУ), основанные на представлении правых частей системы на шаге h в виде трёх точечных интерполяционных полиномов Эрмита (LMR(L,M,R,s) —алгоритмы, se [0.5,1.0)—координата внутренней на интервале h точки коллокации). Погрешность методов — ^L+M+R+4 да1Ю определение ¿(«^-устойчивости одношаговых методов с малым параметром S. Для методов LMR(L,—],R,]) определены условия А- и ¿-устойчивости, для методов LMR(0,0,0, i),
0,1, j), £Л/Л(1,1,1, s)—условия А- и ¿(б)-устойчивости. Разработан
алгоритм расчета глобальной ошибки и алгоритм решения задачи Коши для систем ОДУ с разрывными правыми частями.
2.2. Разработан метод решения системы линейных интегральных уравнений Вольтерра 1-го рода с разностным аргументом для варианта, когда исходная информация (значения измеряемых функций и ядра) заданы в дискретных точках с известной ошибкой. Решение определяется в классе кусочно-постоянных и кусочно-линейных функций с использованием условия равенства нулю средних значений невязок уравнений на интервалах, на которые разбивается область определения решения. Число интервалов и распределение их длин определяются посредством минимизации среднеквадратичной невязки уравнений.
2.3. Разработан метод восстановления действующих нагрузок в классе кусочно-постоянных функций при испытаниях моделей в аэродинамических трубах кратковременного действия. Приведены примеры решения задач по определению аэродинамических характеристик эталонной модели НВ-2, демонстратора Ares и модели Expert по результатам испытаний в аэродинамической трубе АТ-ЗОЗ ИТПМ СО РАН.
2.4. Для решения задач математического программирования методом штрафных функций из условия локального минимума вспомогательного функционала получена оценка для коэффициента штрафа к — £.
2.5. Созданы комплексы программ для решения названных классов задач.
3. Разработан эксергетический метод оценки характеристик и анализа прямоточного воздушно-реактивного двигателя. Для графического отображения возможных схем подвода тепла в канале ПВРД предложена диаграмма в координатах "полная температура-эксергия". Получено выражение для изменения эксергии в термодинамической системе при подводе тепла и наличии необратимых процессов.
4. Разработана методика оценки эффективности подвода тепла перед ЛА при полете со сверхзвуковой скоростью. Полет происходит на границе раздела сред различной плотности (режим глиссирующего полета). Показана значительная эффективность такого способа управления обтеканием JLA как при крейсерском полете, так и при полете с ускорением.
Содержание работы
Во введении сформулированы цели и задачи работы, обоснованы актуальность и практическая значимость развиваемого направления исследований, даны описания применяемых методов. Отмечается важная роль функционального математического моделирования на начальном этапе исследований. Приводится краткая характеристика состояния исследований с соответствующими ссылками на работы в части прямоточного воздушно-реактивного двигателя, энергетического метода управления обтеканием летательных аппаратов, методов решения жестких систем обыкновенных
дифференциальных уравнений, методов восстановления действующих сил и моментов на аэродинамические модели и модели ПВРД при их испытаниях в аэродинамических трубах кратковременного действия.
Глава 1. Численные методы. Состоит из четырех разделов. 1. Численные методы решения задачи Коши для жестких систем обыкновенных дифференциальных уравнении на основе многозвенных интерполяционных полиномов Эрмита.
1.1. Постановка задачи. Рассматривается задача Коши для системы обыкновенных дифференциальных уравнений
уЬ)=/(у), хе[0,Х], У о = д40)> (11)
п- число уравнений. Предполагается, что выполнены условия существования и единственности решения задачи. Пусть для хе [0,х(] решение получено. На отрезке [х,,х, + И < X] представим систему (1.1) в следующем виде
= у{0) = у(х,), $ = (*-*()/й,4е[0,1]. (1.2)
1.2. Интерполяционный полином Эрмита по трём точкам. Пусть функция /(х) задана на отрезке [х,,хк ] в К точках х, < хг < ... < Хк вместе со своими производными до порядка Q — Ь + М 4- Л + 3, Ь,М,Я- заданные целые числа. На отрезке Х1 = [х:,х^2] значения функции и её производных
точках -+1»Х(+21 будем обозначать нижним индексом соответственно "О", "5", "1". Функция /(х) на рассматриваемом отрезке представлена интерполяционным полиномом Эрмита следующего вида, где сделана замена переменной
£*=(х-х1)/КИ=хм-х08 = (хм-х1)/К!;еШ, яе(0.5Д), 7(4) -приближенное значение,
Я£>=1 Л1)хо,,(%)+1 Пт>I Лг)х^У- (!-3)
/=0 т-0 г=0
Для ошибки аппроксимации получена оценка
5, =тах /(х)-/(х) <----ъ-
(1.4)
¿.=тах/ш,(*)
' хеХ. * I
1.3. Семейство ЬМК-методов. Зададимся целыми положительными величинами Ь, М и Н. Функции правых частей Ф(£,) представим на шаге Н интерполирующими полиномами Эрмита вида (1.3), подставим в (1.2) и проинтегрируем на отрезках [0,5], и [0,1]. Получим следующую систему алгебраических уравнений
I м я
У(*)=У, + £/0.,(У < + £
;=0 т-0 г=0 ..
I М Л ^ ' '
у(у)=у,+£/0/и< +
1=0 т=0 г=0
ООО Методы решения задачи Коши для систем ОДУ, основанные на пошаговом представлении функций правых частей в форме (1.3) и решении систем
алгебраических уравнений (1.5), назовём LMR(L, А/, R, .sj-методами. Погрешность методов удовлетворяет оценке
Решение линеаризованной по у системы ОДУ (1.2) можно использовать для получения первого приближения в итерационных методах решения алгебраической системы (1.5). Приведены схемы трех алгоритмов из данного семейства: LMR(0,0,Q,s), LMR(\,0,\,s), LMR(\,\,\,s), и доказаны теоремы об устойчивости. Метод LMR(0,0,0,S) пригоден для интегрирован™ неавтономных систем ОДУ.
1.4. Устойчивость методов. Для исследования устойчивости метода используется линейное уравнение
z'(x) = Xz(x), z(0) = zB, х > 0, Re(X)<0,n = Xh (1.6)
Решение на шаге для одношагового метода записывается так z(|l) = z„q(|x).
Определение. Одношаговый метод назовем £ (б)-устойчивым с параметром
5 е (ОД), если метод А -устойчив и lim |q(n| = 5.
М-**
Теорема 1. LMR(0,0,0,s)-метод Л-устойчив при ie[0.5,1.0) и £(б)-устойчив с параметром 8 = (l — s)/s при 5 6 (0.5,1.0).
Теорема 2. ¿М/?(1,0,1,5)-метод А -устойчив при J6 [0.5,1.0) и ¿(б)-устойчив с параметром б = (l - s)/s при s е (0.5,1.0).
Теорема 3. LMR(l,\,l,s)-метод Л -устойчив при £б[0.5,1.0) и ¿(б)-устойчив с параметром ö = (l — s)2/s2 при S6 (0.5,1.0).
Величина S может быть выбрана достаточно малой. В этом случае ¿(5)-устойчивость, как и ¿-устойчивость, является полезным свойством метода для эффективного численного решения сильно жестких систем.
Теорема 4.
LMR{L,-\,R,1) -методы А -устойчивы при L< R< L + 2 и L -устойчивы при L < R < L + 2.
1.5. Расчет локальной и глобальной ошибки. Ошибка 5у= у- у на шаге h определяется из решения квазилинейной системы ОДУ
Для приближенного значения используется форма остаточного члена
аппроксимации функции интерполяционным полиномом Эрмита. Для LMR(0,0,0,s)-метода g(^)« ß(^) = - — l). Для других методов выражение строится аналогично. Коэффициент С, определяется из условия )=ß(^.)> гдс некоторое выбранное значение, например, i;. = 0.5.
Якобиан аппроксимируется квадратичной функцией от ^ + + — <Ü)D. Матрицы A,B,D определяются из условий J(0) = J0,J(1) = J1,J(s) = Js . Интегрирование системы уравнений (1.7) производится методом LMR{ 1,—l,l) при начальном значении ошибки бу0, равном глобальной ошибке интегрирования на предыдущем шаге. Локальная ошибка равна fyloc $>о > используется для регулирования шага h.
1.7. Интегрирование системы ОДУ с разрывными функциями правых частей. Рассматривается следующая задача Коши для системы ОДУ
/(*)=/,(Л*), хф,Х), у0 =^(0), 1.8)
Здесь вектор-функция £ (д>, х) является г -м столбцом матрицы функций определяемой ниже. Изменение номера столбца происходит в точках обращения в нуль задаваемых задачей индикаторных функций Н(у)\[Хт\ (точки разрыва функций правых частей уравнений (1.8)). Предполагается, что Н(у) — непрерывные функции с константой Липшица Ь. В начальной точке х = 0 каждой функции Н1 (_)>0), j = 1, т, ставится в соответствие целочисленный индикатор I) такой, что /^=1, если Н1 (_у0)^0, иначе =2. Составляется целочисленная матрица М1[тхМ\М = Т
', столбцы которой /,,г = \,М, образуют все возможные комбинации значений индикаторов. Исходя из содержания задачи, каждому столбцу 1,,1 = 1,М ставится в соответствие вектор-функция /¡{у,х) матрицы .Р. Если при каком-либо значении х = х происходит смена знака функции Н] (у), то меняется значение соответствующего индикатора по правилу = 3 — IJ. Это определяет
выбор нового столбца матрицы Р для продолжения интегрирования. При этом, используя термины геометрической оптики, возможно отражение траектории от границы, преломление и полное внутреннее отражение (движение вдоль границы). Минимальный интервал разрешения точки разрыва по х задаётся. Для регулирования шага интегрирования задается также какая-либо монотонно убывающая функция К(Г1) такая, что К( 1) = 1. Приводится схема алгоритма интегрирования системы ОДУ
(1.8) с разрывными правыми частями методом ЬМК(0,0,0,5).
Тестирование методов произведено на пяти жестких системах ОДУ и двух системах с разрывными функциями правых частей. Сравнение выполнено с тремя методами из числа доступных программных пакетов: метод Гира, использующий автоматически регулируемый переменный порядок локальной точности выше 7-го; неявный метод Рунге - Кутгы 5-го порядка точности, использующий вторую производную от решения (программная версия 1995 г.); стандартный метод Рунге -Кутты 4-го порядка точности. Сравнительные результаты решений демонстрируют высокую эффективность рассмотренных методов - число обращений к процедуре вычисления функций правых частей существенно меньше - и позволяют получать решения с высокой точностью.
Замечания. Аппроксимация функций многозвенными полиномами Эрмита в форме (1.3) удобна для использования во многих приложениях. Так как коэффициенты полиномов являются значениями функции и ее производных в узловых точках, то не требуется записывать условия сопряжения аппроксимирующих функций смежных интервалов. Примеры применений: 1) для вычисления определенных интегралов, а также для их приближенного аналитического представления; 2) для вычисления функции с заданной точностью в любой точке отрезка, если известны точные значения функции и её производных нужного порядка на некотором множестве реперных точек. Например, для тригонометрических функций можно использовать точки {я/2,0; л/3, л/6; я/4; я/8, Зя/8;...}; 3) для построения кривых и поверхностей; 4) для задания управляющих функций в задачах оптимального управления.
2. Метод интервального осреднения для определения квазнрешсння линейных интегральных уравнений Вольтерра 1-го рода с разностным ядром.
Математическая модель многих измерительных приборов описывается системой линейных интегральных уравнений Вольтерра 1 -го рода с разностным ядром вида
= = /е [О,Г]. (2.1)
о
Уравнение (2.1) описывает линейную динамическую систему, имеющую т входов и п выходов (т<п). Матрица X т] является откликом на ступенчатое
единичное воздействие (нормальная реакция). На матрицу нормальных реакций наложены условия, соответствующие физической системе (скорость распространения сигнала конечна)
1У(0И,^в = С/(0) = 0. (2.2)
Л
Вектор-функции К(/)[лХ1] - реакция системы на воздействие вектор-функции Как правило, измерения функции У(?) производятся в дискретных точках с известной ошибкой. В этих же точках известны значения элементов матрицы {/(/) и также с ошибкой. Функции С(?) подлежит определению.
В задаче (2,1) - (2.2) условия существования и единственности решения не выполнены, задача некорректна. Для решения некорректных задач используются методы регуляризации. Однако любая регуляризация портит исходные уравнения.
В данной работе для определения квазирешения задачи предложен метод на основе выполнения уравнения (2.1) в среднем на каждом интервале. Квазирешения строятся на множествах ступенчатых и кусочно-линейных функций. Отрезок [О, Т] разбивается на N интервалов. Условие равенства нулю средней на интервале невязки исключает проблему выбора точек коллокации. При этом, по крайней мере в одной точке внутри интервала, уравнения выполняются точно; усредняются также ошибки измерений. И возможен выбор такого разбиения отрезка [О,?1], что осредненные на
интервалах матрицы не будут вырождены. Значения функций на интервалах в
случае ступенчатых функций и значения функций в узлах в случае кусочно-линейных функций определяются из рекуррентных соотношений. Во втором случае значение С0 =(/(0) предполагается известным (для измерительных приборов (70 = 0) либо может быть определено из условия минимума среднеквадратичной невязки уравнений на первом интервале. При необходимости порядок аппроксимации может быть увеличен посредством многозвенных двухточечных интерполяционных полиномов Эрмита, коэффициенты которых определяются из условия минимума среднеквадратичной невязки уравнений на интервалах.
2.1. Аппроксимация кусочно-постоянными функциями. Выполняя
осреднение (оператор (•)) в уравнении (2.1) на г-ом интервале, получим решение (штрих сверху означает транспонирование)
с, = (< и, >' ■ < и, >)~1 • < и, >' •[< у, > + < а >}
12(14)= I\и({-1м)-Щ1-11)]-ё,:1=Ш1л ' (2'3)
1=0
Если матрица <Ut> вырождена на любом разбиении, то это означает, что имеются идентичные каналы реакции или воздействия. Квазирешение Gni = 0,N — 1 зависит от количества интервалов N и распределения их длин /?, = — , которые определяются из решения задачи минимизации функционала среднеквадратичной невязки
, . 1 Г
«К^Л^Ь «и/ ■= \z'{t)-z(t)dt. (2.4)
ЛГ,А, ^ О
Для получения точного решения задачи (2.4) требуется выполнить громадный объем вычислений. Приближенное решение в данной работе определяется посредством алгоритма циклического покоординатного сканирования функционала в заданном числе точек с выбором лучшей точки для последовательности значений N.
2.2. Аппроксимация G(t) кусочно-линейными функциями. Принято
G0 = 0. Производя операции, как в предыдущем разделе, получим значение Gi+i
AGi = (<У/>'-< Vt >)"< V/ >' •[< % > + < Qt■ > - < Uj > -G, ]
вМ= I j[v{t~tM)GM-Û{t-ti)G,}-
1=0
1
- j U{t-f)dT
-04(2.5)
i t _
Vi (t, ?,) = — ! U{t - r) dt; Gi+l =Gj+ AG, ;i = 0,N-1;
воздействии в модельной
"II,
Метод тестирован на задаче восстановления динамической системе при т — п = Ъ.
2.3. Динамический метод определения аэродинамических характеристик моделей по результатам экспериментов в аэродинамических трубах кратковременного действия. Разработана методика, учитывающая динамику движения модели и непостоянство параметров потока. Предполагается, что регистрируемые во времени реакции тензометрических весов описываются обыкновенными дифференциальными уравнениями с постоянными коэффициентами, т.е. система «модель + державка + тензометрические весы + система крепления» является линейным динамическим объектом. Основания: действующие нагрузки таковы, что связь между деформациями и напряжениями описывается законом Гука; слабая зависимость показаний тензометрических весов от функции распределения
нагрузок, действующих на модель, при фиксированных их интегральных значениях; это подтверждается, например, результатами испытаний модели прямоточного воздушно-реактивного двигателя и его элементов в стационарной аэродинамической трубе при числах Маха потока М = 2;4;6 при одновременном измерении сил и моментов механическими и тензометричес-кими весами; размеры датчиков малы относительно характерного размера Рис. 1.1. Схема динамической тарировки изменения температур тела
модели НВ-2.
модели и датчиков за время эксперимента незначительно. В этом случае, получив экспериментально методом разгрузки реакции системы на единичные нагрузки (динамическая тарировка, рис. 1.1), определяется система интегральных уравнений вида (2.1) для восстановления действующих нагрузок С(() по регистрируемым во времени реакциям тензометрических весов К(/). Динамическая тарировка моделей должна производиться в аэродинамической трубе в рабочей конфигурации.
Эталонная модель НВ-2. Аэродинамические эксперименты проводились в аэродинамической трубе АТ-303 ИТПМ СО РАН, в которой были установлены конические сопла с углом полураствора /7 = 8°. В таблице 1.1 приведены параметры потока, при которых проводились испытания моделей: Мехр - экспериментальное значение числа Маха, единичное число Рейнольдса, £>с - диаметр сопла в месте расположения носика модели. В таблице 1.2 приведена статическая тарировочная матрица IV тензометрических весов в системе координат, связанных с моделью. Наблюдается достаточно сильное влияние на показания по каналу силы X - силы У, продольного момента М2 и момента Му.
Таблица 1.1
Ос 300 544 300 590
мехр 10.0 9.9 11.7 11.8 13.5 13.9 17.1 18.0
Ле тт тах тт тах тт тах тт тах
Ю-6 Кех 5.8-7.4 31-36 2.1-2.3 11-13 5.1-5.3 26-32 4.1-4.7 11-12
Таблица 1.2
X У 1 Мх му Мг
X 1 -0.0776 0290 0.0022 -0.0579 -0.0687
У 0.0056 1 -0.0137 0.0066 0.0184 0.0194
г -0.0053 -0.0008 1 0.0013 -0.0396 0.0027
мх 0.0224 0.0138 0.0276 1 -0.0840 0.0023
му 0.0006 0.0007 -0.0555 0.0019 1 -0.0027
мг 0.0030 -0.0041 0.0049 0.0003 0.0050 1
На рис. 1.2 приведены нормальные реакции по продольной компоненте Х-11у, нормальной компоненте и продольному моменту М2-1на единичные
нагрузки по 5 каналам для эталонной модели НВ-2, аэродинамические характеристики которой стандартизованы. Канал Му в процессе испытаний вышел из строя. Видно,
что в пределах временного окна ?е [0,60- 70] [мс] рабочего режима аэродинамической трубы влияния нагрузок по всем каналам на показания по каналу сопротивления весьма значительны. На рис. 1.3 приведен пример зависимости от времени значений измеренной продольной силы Сх и вычисленного коэффициента сопротивления сг
Рис. 1.2. Нормальные реакции тестовой модели НВ-2
Рис. 1.3. Зависимость от времени экспериментальных значений продольной силы Сх (а) и соответствующего коэффициента сопротивления С, (б) при числе Маха набегающего потока М^ =12 и угле атаки й = 12° модели НВ-2.
Для сопоставления этих результатов с данными, полученными в стационарных трубах, производится осреднение по времени. На рис. 1.46 приведены значения с, с исключением влияния коничности потока по методу Ньютона для различных значений числа Маха, Рейнольдса и угла атаки. Различие результатов с данными ОЫЕЯА при близких значениях параметров потока составляет около 3% (в пределах погрешности эксперимента, рис. 1.4а).
Рис. 1.4. Осредненные значения коэффициента продольной силы ct.
Получены также аэродинамические характеристики аэрокосмического демонстратора Ares (рис. 1.5) и модели Expert (рис. 1.6).
Рис. 1.5. Модель Ares Рис. 1.6. Модель Expert
Замечания. 1. Представленная методика восстановления действующих на модель переменных во времени сил не накладывает ограничений на массу модели. 2. Могут быть ослаблены требования на качество проектирования и изготовления тензометрических весов в части предельной минимизации взаимовлияния, т.к. методика учитывает взаимовлияние регистрируемых сигналов по всем каналам. 3. Так как деформации чувствительных элементов тензометрических весов происходят при любых нагрузках, то даже при необходимости получения информации о числе компонент обобщенной силы меньше шести требуется использовать шестикомпонентные тензометрические весы.
4. Задача математического программирования. Для решения задачи методом штрафных функций получена оценка коэффициента штрафа из условия локального минимума составного функционала: к ~ £. В этом случае изменение минимизируемого функционала и составляющей штрафа при нарушении ограничений на величину £ будут одного порядка, что позволяет производить эффективный численный поиск решения. Необходимое значение к для конкретной задачи устанавливается в процессе поиска: если значения ограничений в конечной точке поиска превышают заданные значения, то необходимо увеличить значение коэффициента штрафа и продолжить поиск.
Глава 2. Эксергетическнй метод оценки характеристик прямоточного воздушно-реактивного двигателя. Эксергия - термодинамическая функция, определяющая максимум удельной работы, могущей быть произведенной газом. Для потока газа эксергия определяет максимальную скорость истечения во внешнюю среду. Поэтому для анализа термодинамического цикла авиационных двигателей и
расчета их характеристик целесообразно применение этой функции. Методы, используемые для расчета характеристик ПВРД, требуют задания некоторого множества определяющих величин, зависящих от газодинамических и геометрических параметров. Эти методы малопригодны при функциональном моделировании, когда отсутствует конструктивная схема моделируемого объекта. Представляется, что эксергетический метод является наиболее подходящим инструментом, т.к. в его основе лежит оценка потерь работоспособности газа в элементах двигателя вследствие необратимости процессов независимо от их природы. Приращение энтропии в каком-либо элементе двигателя может быть задано на основе оценки верхнего и эталонного значений. Для получения верхней оценки необходимо знание условий, при выполнении которых возможно подвести к потоку данное количество энергии. Требуется также знание максимально допустимого значения диссипации энергии. Если построить какой-либо базовый процесс, то на его основе можно вычислить эталонную оценку. Рабочее значение приращения энтропии в элементе задается как взвешенная сумма этих оценок (может быть задано и при наличии только верхней оценки). Их функциональные зависимости от входных и внешних условий дают основания для предположения о слабой зависимости весовых коэффициентов от режима работы двигателя. Эксергетический метод позволяет минимизировать количество определяющих параметров и строить иерархическую структуру математических моделей для вычисления их значений.
2.1. Эксергия термодинамической системы. Пусть в термодинамическую систему поступает в единицу времени количество тj какого-либо вещества с параметрами Р\,Т\,V\ и выходит такое же количество с параметрами ,?2> (р -давление, Т - температура, V - скорость). Параметры внешней среды р^.Т^. К системе подводится количество тепла Q, и системой совершается работа А. Для изменения эксергии Ле в термодинамической системе из определения эксергии е и закона сохранения энергии получим выражение
Ae=Q-A-T„AS-Se. (2.1)
Приращение энтропии Л? обусловлено необратимыми процессами в системе и может быть вычислено при задании их физических моделей. Дефект эксергии де связан с различием составов входящих и выходящих веществ. Для рабочих газов в ПВРД эта величина мала. Относительная степень необратимости процесса в термодинамической системе характеризуется коэффициентом потери эксергии S— T„AS/в2тах, либо
S = T„AS/Aemax.
Замечание. В книге «В.А. Кирилин, В.В. Сычев, А.Е.Шейдлин. Техническая термодинамика. М., Наука, 1979» эксергия потока тепла Q, отдаваемого телом с температурой Т, определяется следующим образом eg =Q(l — T^/T). Применение
этого соотношения для расчета эксергии продуктов сгорания в ПВРД затруднено, т.к. температура в процессе сгорания топлива - переменная величина и зависит от процесса.
2.2. Условие подвода заданного количества тепла. Для получения верхней оценки приращения энтропии рассмотрим одномерное стационарное течение в канале с переменной площадью сечения F(x) при подводе тепла, а также наличии процессов диссипации энергии. Значения параметров потока (число Маха , давление , температура Т^) во входном сечении F„ заданы. Определим условия, при выполнении которых возможно подвести к потоку заданное количество тепла Q{x)
между входным сечением Fx и сечением F(x). Обозначено также: скорость звука ; температура торможения Tq^ ; критическое сечение . Относительные величины: Т0х = Т0/Т„ ,AS = AS/R,Q(x)= Q(x)/al,
e(x)=(y-l)Q(x)/T0„, F(x) = F(x)/F„,f(x) = F(x)/Ft„. Из соотношений одномерного стационарного движения совершенного газа следует, что максимально допустимое приращение энтропии достигается в сечении при числе Маха М — 1
= + в{х))+1п fix). (2.2)
2 у-\
Второй член в (2.2) определяет максимально допустимое приращение энтропии вследствие необратимых процессов при отсутствии подвода энергии.
2.3. Математическая модель ПВРД. Основные элементы ПВРД -воздухозаборник, камера горения и сопло. Обозначения: коэффициент избытка воздуха в камере сгорания Ctj — Щ при (Х\ > 1 и #2 =' ПРИ 0 < СС] < 1; L0 -стехиометрический коэффициент; Ни - калорийность водорода; gT = 1 /(c?i L0) -удельный расход топлива; у/ - полнота сгорания; а„ - скорость звука невозмущенного потока; R„,Rg - газовые постоянные воздуха и продуктов сгорания; Vc - скорость истечения продуктов сгорания из сопла; Vce = (l + gT )(FC / V„) -эффективная относительная скорость истечения; Fq - площадь входного сечения двигателя; AS/(i = l,2 - приращения энтропии в воздухозаборнике, камере горения; <рс- коэффициент скорости сопла. Из (2.1) следует соотношение для эффективной скорости истечения из сопла
V2
г се
= 0 +
(i-^+ig-ö-^)
Mi
Параметры определяются следующим образом:
<% 2AS, г чЖ2 — w Ни ~ ,
}Mi }Q аг L^ai
Ol
(2.3)
(2.4)
Полнота сгорания топлива существенно зависит от числа Маха полета и определяется многими факторами. Здесь используется термодинамическая оценка y/T(M„) = Q[ZHu -, Q[Zj) тепловой эффект реакций при нормальных условиях
для равновесного состава Zj, включающего перечисленные компоненты: (О2, Н2,
Н20, ОН, Н02, N2> NO, N02, N20, NH, NH2, NH3, HNO, О, H, N, Ar). Эта оценка верхняя, т.к. неравномерность параметров потока, неравновесность процессов уменьшают полноту сгорания. В таблице 2.1 приведены значения Ц/j- (Л/м) при (Х\= 1 и параметрах потока на входе в камеру горения, соответствующих числам Маха полета . Эти значения используются в последующих оценках характеристик ПВРД.
Таблица 2.1
Термодинамическая оценка полноты сгорания
3 4 5 6 8 10 12 14
¥т 0.94 0.92 0.90 0.89 0.85 0.79 0.72 0.64
По данным различных авторов наблюдается значительный разброс экспериментальных оценок полноты сгорания. Расчетные оценки в работе «Rogers R.C., Capriotti D.P., Guy R.W. Experimental supersonic combustion research at NASA Langley // AIAA-98-2506, p. 1-23» дали следующие результаты. При параметрах потока
на входе в камеру сгорания Мк =2.0-3.0,Щ =1940°К;(МХ =6 + 8) полнота сгорания составляет {// = 0.80 + 0.05. С увеличением температуры величина у/ уменьшается. При коэффициенте избытка воздуха СС\<\ величина у/ ~ const и равна значению при Ctj = 1, а при Щ=2 достигает предельного значения у = 1. Такое же качественное поведение величины у/т.
Коэффициенты потери эксергии в элементах ПВРД. Величины оцениваются по экспериментальным данным. Величина S\ определяется через коэффициент восстановления полного давления V в воздухозаборнике <5) = — 2 In и слабо зависит от числа Маха полета - ее оценка для
относительно малых чисел Маха полета лежит в пределах =0.040-*-0.055.
Значение ôx ограничено величиной ôXmax = -2ln fw/{yM^),fw = Fw(F*„)~1.
На нерасчетных режимах воздухозаборника с использованием соотношения (2.2) получена оценка для коэффициента расхода X ■ Для сопротивления жидкого контура
используется оценка Cxj = 2(l — х)^ > ® ~~ предельный угол отклонения потока при сжатии (¿У2 ~ 0. l). Оценка величины 8} вычисляется через коэффициент скорости
сопла ЖРД: <5"3 = 1 - ç2c, <рс = 0.945 + 0.975.
Эффективность подвода тепла зависит от многих факторов: смешение струй, трение, наличие отрывов потока и ударных волн, неравномерность параметров, нестационарность, состав и скорости химических реакций. Моделирование этих процессов при разумных затратах вычислительных ресурсов весьма проблематично. Однако по соотношению (2.2) вычисляется максимально допустимое приращение энтропии ÄS2max в камере сгорания при заданной степени расширения F2 = F2/Fw по значению f2 = Fjfw и известному относительному подводу тепла 9. Нижняя грань AS2mi„ определяется при подводе тепла при постоянном давлении (базовый процесс). Температура воздуха на входе в камеру сгорания вычисляется по процессу сжатия воздуха в воздухозаборнике при заданных значениях Fw,ôj. Рабочее значение приращения энтропии в камере горения задается как взвешенная сумма AS2 = (I ~ £2)^2 от/л — 1 ~ задаваемый параметр. Значение
параметра 82 - по соответствующему соотношению из (2.4).
Тяга, удельный импульс. Далее определяются коэффициент тяги Cfcjf, тяга R и удельный импульс I ПВРД
С,, ff = 2Z(Vce -1)-cX];R = cReff q„F0;I = C^M~a".
На рис. 2.1 представлены значения удельного импульса различных двигателей в зависимости от числа Маха набегающего потока. Оценки удельного импульса ПВРД
Рис. 2.2
на водороде по представленной модели согласуются с расчетными оценками других авторов.
2.4. Тй,ё-диаграмма. Для графического отображения процесса преобразования энергии в двигателе предложена Т0,ё-диаграмма в координатах «полная температура, эксергия» (нормированы соответственно на температуру торможения Тйоо и эксергию набегающего потока). Получены дифференциальные соотношения для изобар полного давления р0 = р0/р0оа и изотерм статической температуры Т = Т/Т0ао. На рис. 2.2 представлены возможные схемы процесса подвода тепла в камере сгорания ПВРД.
Процесс si соответствует сверхзвуковому течению в камере горения, d -дозвуковому для относительно малых чисел Маха полета Мж =Зн-7. Вариант d иллюстрирует псевдоскачковый режим подвода энергии, в котором средняя температура подвода тепла больше, чем в варианте si. Поэтому возможно большее значение эксергии при одинаковом количестве подводимой энергии. При достижении максимально допустимой статической температуры подвод тепла происходит при Т = const (вариант s2, М_ >8). При импульсно-периодическом подводе энергии при числах Струхаля Sh »1 может быть реализован режим, близкий к подводу тепла при постоянном объеме. Этот режим наиболее эффективен, т.к. при этом максимальна средняя температура подвода тепла (условно показан линией v = const.).
Глава 3. Численное моделирование нестационарного течения в канале переменного сечения при распределенном импульсно-периодическом подводе энергии. С целью определения условий для возможности использования прямоточного канала при гиперзвуковых числах Маха полета Л/м >6 — 7 рассмотрены нестационарные течения в модельном канале ПВРД при импульсно-периодическом подводе энергии.
3.1. Квазиодномерное нестационарное течение. В данном разделе приводятся результаты численного моделирования квазиодномерного нестационарного течения в канале, моделирующем элемент ПВРД и состоящем из участков с постоянным и расширяющимся сечениями. Конфигурация канала показана на рис. 3.1 (х1,Х2,Хт1,0 -варьируемые параметры). В классической схеме подвод энергии в камере сгорания осуществляется за счет сгорания топлива в некотором политропном процессе. В данном случае энергия подводится к потоку газа в импульсно-периодическом режиме
Рис. 2.1. Значения удельного импульса различных двигателей в зависимости от числа Маха набегающего потока. Топливо водород: 1 - ПВРД изменяемой геометрии; 2 - ГПВРД; 8 - ГПВРД по данной математической модели (знак «о»); 6 - ЖРД; углеводородное горючее: 3 - ТРД; 4 - ПВРД; 5-ГПВРД; 7-ЖРД.
Рис. 3.1.
равномерно в заданном диапазоне [х^л^]' ^ - период подвода энергии. Исключение из рассмотрения процессов смешения позволило определить непосредственное влияние параметров подводимой энергии (мощности, частоты импульсов, распределения источников по длине канала) на характеристики течения. Число Маха потока во входном сечении канала варьировалось в диапазоне значений = 2.4 -ь 4.0. Соответствующие числа Маха полета =6-¡-12. Для расчета параметров течений использовались
нестационарные квазиодномерные уравнения Эйлера. В качестве меры мощности подводимой энергии принята величина максимальной мощности энергии, выделяющейся при сгорании водорода в воздухе. Предполагается, что длительность импульса настолько мала, что изменением плотности газа и его скорости за соответствующий промежуток времени можно пренебречь. На входе в канал задаются параметры невозмущенного потока, на выходе используется линейная экстраполяция. Для решения задачи в промежутках между моментами подвода энергии применяется метод Р. Маккормака с искусственной вязкостью четвертого порядка малости.
При импульсно-периодическом подводе энергии в начальный сверхзвуковой поток в канале после завершения переходного процесса при больших значениях числа Струхаля устанавливается периодический режим течения с малыми амплитудами колебаний параметров потока. Такой сценарий развития течений реализуется во всех рассмотренных вариантах. При периодическом режиме течения параметры газа, осредненные за период, стационарны. Почти стационарны также параметры, при которых осуществляется импульсный подвод энергии, и для анализа характеристик течений применимы стационарные соотношения. Относительное увеличение температуры газа при однократном подводе энергии определяется выражением
1]{х) =
ЛТ(х)=Г(у-1рй(х) 0 =
Т*{х) БИТ,(х) ' и0Л1
м1
Ц>ао 1У1 о
где й(х) = и(х)/щ,Т,(х) = Т,(х)/Т0 - распределения скорости и температуры газа в зоне [х1, ] > ЗЬ - число Струхаля, к < 1. При 5А » 1 нагрев газа осуществляется малыми порциями {?]{х)«1), соответственно малы и импульсы давления. В противном случае, при возникновении скачков уплотнения, генерируемых большими импульсами давления, в них будут происходить дополнительные потери эксергии. Здесь расчеты выполнены при числе Струхаля 5Л ~ 30 + 200.
Изучено влияние параметров задачи М0, 57), , > , на структуру формирующихся периодических течений. Получена частота "насыщения" импульсов энергии, при превышении которой структура и осредненное за период распределение параметров периодического течения не меняются. Эта частота слабо зависит от параметров задачи, и соответствующий период равен Д1 = 10~3. Стабилизируется также осредненная за период удельная сила (рис. 3.2).
0,4 0.2 0,0 -0,2
о
10
Рис. 3.2. Зависимость осредненной за период удельной силы / от времени
/={{р+р«2)у]А(р+ри2)уШриу\
М0=ЗД к = 0.75, р = 2\ х} =2, =1.0, =1.6
Рис. 3.3. Распределение числа Маха по длине канала(М0 =2.4, х}=2):
1 - к = 0.3, XI = 0.50, х2 = 0.52, /7 = 2°;
2- ¿ = 0.3,*, =1.0,=1.4,>3 = 2°;
3 - А = 0.4,*, = 1.0,=1.4,= 2°;
4- ¿ = 0.5,*, =1.0,х, =1.4,/? = 45°.
Приведем некоторые результаты. При начальном значении числа Маха Мд = 2.4 и значении энергии, подводимой в поток в цилиндрической части канала, соответствующем значению к = 0,3, всюду сохраняется сверхзвуковое течение (кривая 1 на рис. 3.3). В зоне подвода энергии увеличиваются давление и температура. Аналогичный характер течения наблюдается при подводе такой же энергии в широкой области расширяющейся части канала (кривая 2). Структура течения сохраняется до значения к = 0.39. При значении к = 0,4 нарушается условие АБ<АБтах), и в цилиндрической части канала формируется скачок уплотнения, за которым поток является дозвуковым (кривая 3). Скачок уплотнения совершает колебательные движения малой амплитуды с частотой, равной частоте подвода энергии. Кривая 4 соответствует течению в канале с большим углом раскрытия [Р — 45°), особенностью которого является наличие скачка уплотнения внутри зоны подвода энергии. Скорость потока за этим скачком дозвуковая, в конце зоны число Маха М = 1. Этот вариант приведен для того, чтобы показать, что такие течения возможны. Результаты расчета при подводе максимальной в принятой мере энергии к = 1 приведены на рис. 3.4.
В этом случае также формируется скачок уплотнения, но в расширяющейся части
2,5-. 2.01,51.00,50,0-
Рис. 3.4. Распределение числа Маха по длине канала при максимальном количестве подводимой энергии М0 = 2.4, к = 1, /? = 6°, хъ =10,
*1=5.0, *2=5.4
м
а
б
3,0 2,5 2,0 1,5 1,0 0,5 0,0-
2-
8
6
4
X
о
X
0,0 0,5 1,0 1,5 2,0
0,0 0,5 1,0 1,5 2,0
4,0
3.53,02,5 2,01,51,0-
0.5-
т
в
Рис. 3.5. Распределение числа Маха (М), давления (Р) и температуры (Т) по длине канала(М0 =3.0, /) = 2°, х3=2, х, =1.0, х2 =1.6): 1 - А = 0,75, 5Й = 200;2- ¿ = 0,80; 1 - 1 - ¿ = 0,75, 57г = 4
х
' 0,0 0,5 1.0 1,5 2,0
канала. Энергия подводится в дозвуковой поток. Число Маха в малой окрестности правой границы зоны подвода энергии также равно единице. Далее поток имеет сверхзвуковую скорость.
На рис. 3.5 представлены распределения параметров по длине канала при значениях к = 0.15 (кривые 1), к = 0.80 (кривые 2) для начального значения числа Маха Mq =3.0. При к = 0.15 скорость потока всюду сверхзвуковая, а при к = 0.80 энергия подводится в дозвуковой поток. Скачок уплотнения расположен в цилиндрической части канала. На рис. 3.5в показано влияние уменьшения числа Струхаля до значения Sh = 4 на распределение параметров потока (представлена только зависимость Т(х), кривая 1 — 1), зависимости М{х), р{х),р(х) аналогичны. Так как подводимая энергия одна и та же, то при однократном воздействии приращение температуры в этом случае существенно больше. Поэтому и наблюдаются значительные колебания параметров потока. При значениях числа Маха Mq =3.6/4.0 при подводе максимальной энергии в поток всюду сохраняется сверхзвуковая скорость течения.
Для пояснения структуры формирующихся течений продолжим начатое в главе
2 рассмотрение одномерных стационарных течений.
Одномерное стационарное течение (продолжение). В главе 2 получено соотношение (2.2) для максимального значения энтропии при подводе тепла и наличии диссипации кинетической энергии. Если AS {х,) - значение энтропии в некотором сечении F(x*), то возможны следующие случаи. 1. AS(х») = ASmax{x*). Стационарное решение единственно, число Маха в сечении равно единице М = 1. 2. Ж(х*)< AS тах (х»). Возможны два стационарных решения. В зависимости от предыстории потока и граничных условий реализуется течение с дозвуковой или сверхзвуковой скоростью. 3. ASmax(x*). Стационарное решение
отсутствует. Но возможна такая перестройка потока, что реализуется первый случай. Уменьшение суммарного приращения энтропии при заданном значении 9 может быть достигнуто увеличением средней температуры подвода тепла, т.е. подводом тепла
полностью или частично в дозвуковой поток. Возникает прямой скачок уплотнения. Положение скачка уплотнения определяется из условия единственности решения, т.е. равенства суммарного приращения энтропии в скачке уплотнения и при подводе тепла максимально допустимому значению в рассматриваемом сечении А5(х5,) = /Ктах (х,), х5- координата скачка уплотнения.
Условия перехода через скорость звука. Из законов сохранения для элементарного объема с!х следует условие для возможности продолжения непрерывного стационарного решения в точке х = х«, в которой М (х* ) = 1,
СрТ ЯТ Г(ХФ) Н(х) — функция диссипации энергии, штрих - дифференцирование по х.
Соотношения (2.2) и (3.1) определяют структуру формирующихся периодических течений при больших числах Струхаля. В данной задаче за скачком уплотнения Н'(х) = 0, и, поскольку ^'(х) > 0, сечение Р{х,) всегда расположено внутри зоны подвода энергии (как правило X* = х2 — £, £ ~ 0). Течение за ним может быть как сверхзвуковым, так и дозвуковым, и зависит от предыстории течения и граничных условий.
При числах Маха полета М_ > 8 необходимо контролировать максимальное значение температуры газа Ттах. Представленные результаты показывают, что для гиперзвуковых скоростей полета даже при обеспечении сверхзвуковой скорости потока в зоне подвода энергии наблюдаются большие значения температуры газа. Для обеспечения
Т < Ттах необходимо изменить конфигурацию канала, так
Рис. 3.6. выполнения условия
чтобы энергия подводилась в нескольких обособленных зонах (секциях), разделенных участками расширения, и сохранялась при этом всюду сверхзвуковая скорость течения, в том числе при подводе максимальной энергии. Возможная конфигурация канала показана на рис. 3.6.
3.2. Двухмерное нестационарное течение. В двухмерном случае появляется дополнительная степень свободы для размещения источников энергии. Разнообразие связных и несвязных зон подвода энергии в сочетании с различными распределениями мощности источников энергии порождают разнообразие в структурах течений. Однако необходимость подвода энергии в количестве, соответствующем тепловыделению при сгорании водорода с коэффициентом избытка воздуха СС ~ 1, является существенным ограничивающим условием.
Постановка задачи. Моделируется нестационарное течение в плоском канале переменного сечения с распределенным подводом энергии. Решаются двухмерные нестационарные уравнения Эйлера в консервативной форме для газа с постоянным показателем адиабаты у. О форме и размерах канала можно судить по рисункам, приведенным ниже. Минимальные зоны подвода энергии имеют приблизительно прямоугольные формы, из которых формируются зоны произвольных конфигураций. Подвод энергии осуществляется каждый раз настолько быстро, что изменением плотности газа и его скорости за соответствующий очень малый промежуток времени пренебрегается. При этом плотность энергии газа е в г-й зоне увеличивается на
величину Ае(х,у) = Ае,х^х < х < х^у^ < у < у, (номер зоны не связан с номерами координатных узлов сетки). Для определения меры предполагается, что мощность подводимой энергии эквивалентна мощности тепловыделения при сгорании водорода с локальным коэффициентом избытка воздуха <Х,
, риО . _ , Ям
Ах
2
Ь0а0
, где ¿4х = л,
- - протяженность г -й зоны подвода
энергии, к = ОС , - период. Для решения уравнений на входе в канал задаются параметры невозмущенного течения, на выходе при сверхзвуковых скоростях применяется процедура экстраполяции решения, на стенке канала ставится условие непротекания, в плоскости симметрии задаются условия симметрии; начальными условиями являются параметры стационарного течения газа при отсутствии подвода энергии. Для нахождения начального стационарного течения и решения в промежутках между моментами подвода энергии используется конечно-объемная схема (ТУО-реконструкция). Интегрирование по времени проводится по методу Рунге — Кутты третьего порядка точности. Расчетная сетка с числом узлов 200x140 геометрически адаптивна к контуру канала.
Результаты расчета. Основные результаты получены для следующего варианта.
На входе в канал задан равномерный поток
11!
ю
Рис. 3.7.
Маха М0=и(0,у^) = 2, давлением р(0,у,() = 1, плотностью и показателем адиабаты газа
p(0,y,t)=y = 1.33;At = 0Л. На рис. 3.7 представлено распределение числа Маха при периодическом подводе
энергии в одной зоне, перекрывающей все сечение канала (границы зоны 9.48 < х < 10.0). Для лучшего визуального разрешения течения масштаб по оси у на рисунках увеличен. Сформировано течение со скачком уплотнения, расположенным в начальной слабо расширяющейся части канала. Скорость газа за ударной волной дозвуковая. Подвод энергии происходит в дозвуковой области течения. В зоне подвода энергии поток разгоняется до звуковой скорости. Изолиния М = 1 примыкает к правой границе зоны. Сформировавшееся течение - почти периодическое: в моменты времени непосредственно перед подводом энергии положение прямого скачка и распределение параметров течения в канале меняются незначительно, расход газа через выходное сечение канала колеблется в пределах нескольких процентов.
На рис. 3.8 показано поле давления после п = 104 периодов при подводе энергии в трех изолированных зонах почти прямоугольной формы, сумма сечений которых составляет примерно половину сечения канала. Количество подводимой энергии эквивалентно суммарному коэффициенту избытка воздуха
Рис. 3.8.
« = 2,5'/; = 10. Прямой скачок уплотнения вышел из «камеры сгорания» и стабилизировался в начальном слабо расширяющемся участке канала. Качественная структура течений сохраняется и при увеличении длин зон по координате X при постоянной суммарной мощности подводимой энергии. Индикатором реализуемости стационарного течения служит двухмерный аналог условия (2.2)
АЗ(Р)<АБтах(Р), (3.2)
где приращения энтропии вычисляются интегрированием по сечению . Если определить осредненные при условиях сохранения расхода, энергии и энтропии параметры течения, то нарушение условия (3.2) в каком-либо сечении означает, что этот поток не может быть пропущен через данное сечение. Это приводит к перестройке течения такой, чтобы по возможности выполнить условие (3.2). Проведенные расчеты иллюстрируют сказанное. Подвод энергии в дозвуковом потоке, т.е. при большей температуре, уменьшает суммарное приращение энтропии. В экспериментах [15, 16] наблюдается такая перестройка: начальный сверхзвуковой поток при горении водорода перестраивается так, что через 25-30 мс горение происходит в дозвуковом потоке.
Приведенные результаты моделирования показывают, что импульсно-периодический подвод энергии позволяет увеличить число Маха полета до значений, при которых возможно использование прямоточного канала в составе комбинированного двигателя для увеличения эффективного удельного импульса.
Глава 4. Оценка энергетической эффективности подвода тепла перед летательным аппаратом при сверхзвуковой скорости полета. Получена оценка увеличения коэффициента дальности Бреге для крейсерского полета гиперзвукового летательного аппарата и оценка экономии топлива на разгонном участке траектории полета с использованием комбинированного двигателя, включающего прямоточный воздушно-реактивный двигатель (ПВРД) и жидкостный ракетный двигатель (ЖРД). Принято, что подводимая энергия вырабатывается за счет эквивалентного уменьшения эксергии продуктов сгорания топлива.
4.1. Математическая модель процесса. Из работ по расчету структуры течения за пульсирующим источником тепла, следует, что статическое давление в следе достаточно быстро становится равным давлению в окружающей среде. Поэтому предполагается, что подвод тепла в набегающий поток осуществляется при постоянных значениях давления и скорости, так что реализуется бесконечный тепловой след с параметрами р = р„,У = = р/р^ = Р^/Р = £ .
Принятые обозначения: р - давление, V - скорость, Т — температура, Р - сечение следа, р - плотность; индекс " соответствует параметрам газа на бесконечности; £- задаваемый параметр. Тепловая мощность следа перед телом равна б», = р„¥„РтСрТ„(1-е), где Рт - мидель ЛА, ср - теплоемкость воздуха при постоянном давлении. Традиционно эффективность подвода тепла в установившемся
Л\- А 7^-1 %,г .
полете оценивается величинои 77 = —-=-суа7и„, где Ап - исходная
2
мощность двигателя, А - мощность при тепловом воздействии. В этой оценке не учитываются полный энергетический баланс, функциональное назначение летательного аппарата и существенные параметры процесса: степень нагрева газа, несущие свойства ЛА, характеристики двигателя, эффективность преобразования энергии топлива в энергию излучения. В данной работе эффективность подвода тепла в набегающий поток оценивается с учетом перечисленных факторов.
4.2. Аэродинамическая поляра. Предполагается, что полет ЛА происходит на границе раздела сред высокой и низкой плотности в режиме глиссирования: подъемная сила Y в основном создается нижней поверхностью ЛА, обтекаемой невозмущенным потоком воздуха. Для аэродинамической поляры получено
соотношение cx=cy(cyj~ (второй член - сопротивление части ЛА,
находящейся в тепловом следе). При нагреве набегающего потока значение коэффициента подъемной силы, при котором достигается максимальное значение аэродинамического качества, меньше аналогичного значения без нагрева. Поэтому для реализации максимальной эффективности крейсерский полет должен происходить на меньших высотах.
4.3. Силовая установка (СУ). Предполагается, что СУ ГЛА является комбинированной и включает ПВРД и ЖРД. Воздух, поступающий в ПВРД, не подогревается. Сила тяги Щ и удельный импульс 1\ ПВРД зависят от числа Маха полета М„ и коэффициента избытка воздуха СС\ в камере горения, Щ - функция управления. Сила тяги ЖРД - R2 также функция управления, а его удельный импульс /2 = const - заданная величина. Для определения характеристик ПВРД используется эксергетический метод, изложенный в главе 2.
4.4. Удельный импульс при выработке энергии Q„. Выработка энергии, подводимой в набегающий поток, осуществляется за счет эксергии продуктов сгорания топлива при дополнительном увеличении энтропии ASq . Получены
выражения для эффективной скорости истечения при выработке энергии Q^ в тракте ПВРД и ЖРД. Параметр Zj определяет долю энергии, расходуемой на нагрев набегающего потока, от полной энергии, выделяющейся в камере горения ПВРД.
4.5. Крейсерский режим полета. Эффективность нагрева набегающего потока крейсерском режиме полета определялась по увеличению коэффициента дальности
2 he
Бреге Вг£ = —- = Вг
1 + е 1Х
Результаты расчетов. На рис. 4.1 представлена зависимость входного сечения двигателя Р0 (М „). При увеличении крейсерского числа Маха полета потребное значение значительно увеличивается, и при числе Маха = 12 достигает значения ~ 0.1. При решении задачи разгона принято это значение, т.к. при меньших значениях ПВРД не будет доставлять требуемой тяги при больших числах Маха полета.
Рис. 4.1
На рис. 4.2 представлена зависимость относительного коэффициента дальности от степени подогрева набегающего потока для различных значений крейсерского числа Маха. Такое управление обтеканием ЛА является весьма эффективным. При этом даже для значительной степени подогрева доля потребной энергии не превышает 8% (рис. 4.3).
Сравнение расходов топлива на разгон с нагревом воздуха перед ЛА и без нагрева должно производиться на оптимальных траекториях.
1,00-1---- м~
4 6 8 10 12
Рис. 4.2. Относительный коэффициент дальности в зависимости от числа Маха крейсерского полета и относительной плотности набегающего потока.
Значения в: 1 - 0.4, 2 - 0.5, 3 - 0.6, 4-0.7,5-0.8,6-0.9.
4 в 8 10 12
Рис. 4.3. Доля энергии, расходуемой на подогрев набегающего потока, от выделяемой в камере сгорания энергии.
Значения Б : 1 - 0.4, 2 - 0.5,3 - 0.6, 4 -0.7,5-0.8,6-0.9.
4.6. Режим разгона. Построение оптимальной траектории. Расход топлива на разгон ЛА определяется из решения задачи
Етх = 1 - _ та* т{щ),
dm
-= -7И-
dw L
RbRi,cyeD
, (l -w2)cosd 1 + Л-L
+ Л (w) nvK
= -m-f, (4.1)
sind =A(w)-rtv, A(w) = ----5-,we \w0,wx\mü =1
(1 + c)y„Mi
Значение с = 0 соответствует полету при qx = const; режим с —> °° соответствует полету при = const. в - угол наклона траектории ЛА, nv =Voe/g -относительное продольное ускорение. D - допустимая область управления.
Замечание. Уравнение (4.1) является обобщением формулы К.Э. Циолковского для полета в атмосфере. При nv —> с —» I = const получим Vтах = -11пт)с.
Для решения задач оптимального управления широко используется метод Ритца, состоящий в замене исходной задачи некоторой вложенной последовательностью конечномерных задач минимизации, решения которых образуют минимизирующую последовательность. Построение минимизирующей последовательности - весьма трудоемкий процесс. Однако в данном случае оказалось возможным упростить решение задачи (4.1). Показано, что при разгоне с использованием только ПВРД оптимальной является траектория q^ =qmax, при включении ЖРД - nv —nvmax-Получены условия включения ЖРД, выключения ПВРД, полета при максимально допустимом значении угла атаки (Xд = ОСд тах. 4.7. Конечная формулировка задачи
= min f,we = 1 ,nv б (О,пУтах]. (4.2)
aw %
Правая часть в дифференциальном уравнении (4.14.1) - разрывная функция. Точки разрыва определяются условиями:
\.л=
(1 + Я)пу
Л
+
т
с Л 1-^1
/2 {\ + А(и'))пуК
К
;3./,=/2
А
Условие 1 соответствует при > 0 полету при <7„ = Чертах> иначе ал = ОСАтах. Условие 2 - моменту включения ЖРД. Условие 3 - моменту выключения ПВРД. Интегрирование производится по методу, изложенному в главе 1. Для вычисления величин, входящих в выражение /, используются приведенные соответствующие математические модели. При отсутствии нагрева набегающего потока в начальной точке полета заданы значения Кт0 и ^.
4.8. Результаты расчета. Исходные данные:
М0 =А,М\=ПХтй= 4Д = 0. \,Е = 0.3,^ = 1;
с 0 = 0.7,/2 = 4700[м/сек],щ
Утах
в 8 Рис. 4.4.
10 12 1« 18
Рис. 4.5.
0,4 0,3 0,2 0,1 0,0
/
4 6 8 10 12 14 16
Рис. 4.6.
8 10 12 14 16
Рис. 4.7.
8 10 12
Рис. 4.8
На рис. 4.4 и рис. 4.5 представлены значения качества и удельного импульса ПВРД по траектории полета без подвода (кривые 1) и с подводом энергии (кривые 2) в набегающий поток. Аэродинамическое качество увеличивается значительно при малом уменьшении удельного импульса. На рис. 4.6 представлено значение продольной перегрузки с указанием моментов включения ЖРД и выключения
ПВРД. На рис. 4.7 - расход топлива. Видно, что нагрев набегающего потока позволяет значительно уменьшить расход топлива на разгон: Д§т = —0.06 при умеренных значениях затрачиваемой на нагрев энергии (рис. 4.8).
Выводы
1. Выполнен анализ квазиодномерного и двухмерного квазистационарных течений в канале переменного сечения, описываемых уравнениями Эйлера и формирующихся при импульсно-периодическом подводе энергии при больших числах Струхаля. Получено, что при этом устанавливается периодический режим течения с малыми амплитудами колебаний параметров. Течения устойчивы в среднем на периоде. Так как энергия подводится при постоянном объеме, то этот режим обеспечивает максимальное значение эксергии потока и, следовательно, тяги двигателя.
Получены условия, определяющие структуру течений: максимально допустимое значение энтропии для каждого сечения канала и условие перехода через скорость звука в квазиодномерном случае при подводе энергии и наличии диссипации кинетической энергии.
Предложена конфигурация канала, в котором подвод тепла к сверхзвуковому потоку осуществляется с учетом ограничения статической температуры газа. Импульсно-периодический подвод энергии в таком канале позволяет увеличить число Маха полета до значений, при которых возможно использование прямоточного канала в составе комбинированного двигателя для увеличения эффективного удельного импульса силовой установки.
В результате численного моделирования нестационарного двухмерного течения в канале переменного сечения при подводе тепловой энергии в локальных зонах в импульсно-периодическом режиме получена экспериментально наблюдаемая перестройка начального сверхзвукового течения, определяемая условием подвода заданного количества энергии.
2. Разработаны новые численные методы для решения следующих задач:
2.1. Построено семейство А-, Ь- и Д6)-устойчивых методов решения задачи Коши для жестких систем обыкновенных дифференциальных уравнений, основанное на представлении правых частей системы на шаге к в виде трёх точечных интерполяционных полиномов Эрмита (LMR.iL, МЛ я)-методы). Погрешность
/„¿+А/+Л+4 „ »
методов ~ п . Дано определение ¿(о)-устоичивости одношаговых методов с
малым параметром 5. Для методов LMR(L,-l,R,i) определены условия А- и
¿-устойчивости; для методов .у), условия
А- и ¿(б)-устойчивости. Разработаны алгоритм расчета глобальной ошибки и алгоритм решения задачи Коши для систем ОДУ с разрывными правыми частями.
2.2. Разработан метод решения системы линейных интегральных уравнений Вольтерра 1-го рода с разностным аргументом для варианта, когда исходная информация (значения измеряемых функций и ядра) задана в дискретных точках с известной ошибкой. Решение определяется в классе кусочно-постоянных и кусочно-линейных функций с использованием условия равенства нулю средних значений невязок уравнений на интервалах, на которые разбивается область определения решения. Число интервалов и распределение их длин определяются посредством минимизации среднеквадратичной невязки уравнений.
2.3. Разработан метод восстановления действующих нагрузок в классе кусочно-постоянных функций при испытаниях моделей в аэродинамических трубах кратковременного действия. Приведены примеры решения задач по определению
аэродинамических характеристик эталонной модели НВ-2, демонстратора Ares и модели Expert по результатам испытаний в аэродинамической трубе АТ-303 ИТПМ СО РАН.
2.4. Для решения задач математического программирования методом штрафных функций из условия локального минимума вспомогательного функционала получена оценка для коэффициента штрафа к~ S.
2.5. Созданы комплексы программ для решения названных классов задач.
3. Разработан эксергетический метод оценки характеристик и анализа ПВРД. Для графического отображения возможных схем подвода тепла в канале ПВРД предложена диаграмма в координатах "полная температура - эксергия". Получено выражение для изменения эксергии в термодинамической системе при подводе тепла и наличии необратимых процессов.
4. Разработана методика оценки эффективности подвода тепла перед ЛА при полете со сверхзвуковой скоростью. Полет происходит на границе раздела сред различной плотности (режим глиссирующего полета). Показана значительная эффективность такого способа управления обтеканием ЛА как при крейсерском полете, так и при полете с ускорением.
Основные результаты диссертации опубликованы в следующих работах:
1. Латыпов А.Ф., Никуличев Ю.В. Численный метод решения задач оптимального управления // Аэрофизические исследования. Изд. ИТПМ СО АН СССР. 1972. с. 102-103.
2. Латыпов А.Ф., Никуличев Ю.В. Об одном методе поиска минимума функции многих переменных // Изв. СО АН СССР, серия технических наук. 1972. вып. 3, № 13. с. 93-95.
3. Латыпов А.Ф. О решении экстремальных задач с ограничениями // Изв. СО АН СССР, серия технических наук. 1974, вып. 3, №13. с.49-50.
4. Латыпов А.Ф. Об одной модификации метода наискорейшего спуска // Изв. СО АН СССР, серия технических наук. 1974. т. 2, № 8. с. 87-89.
5. Латыпов А.Ф., Хенкин П.В. Параметрический анализ прямоточного воздушно-реактивного двигателя на водороде со сверхзвуковым горением // Вопросы газодинамики. Новосибирск: Изд. ИТПМ СО АН СССР, 1975. с.197-201.
6. Дулов В.Г., Латыпов А.Ф., Пупышев С.Б., и др. Эффективность крейсерского полета гиперзвуковых летательных аппаратов // Исследования по гиперзвуковой аэродинамике. Новосибирск: Изд. ИТПМ СО АН СССР, 1978. с.151-172.
7. Латыпов А.Ф. О математическом моделировании летательных аппаратов на этапе выработки концепции // Численные методы механики сплошной среды. 1979. т. 10, №3. с. 105-110.
8. Латыпов А.Ф., Тенетов В.П. Функциональная математическая модель силовой установки гиперзвукового летательного аппарата: Препринт ИТПМ СО АН СССР №4-83, Новосибирск, 1983. с. 1-31.
9. Латыпов А.Ф., Никуличев Ю.В. Специализированный комплекс программ оптимизации: Препринт ИТПМ СО АН СССР № 15-85. Новосибирск, 1985. с. 1-41.
10. Латыпов А.Ф. Функциональная математическая модель прямоточного и ракетно-прямоточного двигателей // Труды V Школы по методам аэрофизических исследований. Новосибирск: Изд. ИТПМ СО АН СССР, 1990. с. 97-103.
11. Latypov A.F., Niculichev Yu.V. New Methods Based on Hermit Approximation for Solving Problems of Guidance, Forecasts and Optimization // Xl-th International Conference on Remotely Piloted Vehicles. Bristol, 1994. Conference papers, 10 p.
12. Perrier P., Stofflet В., Rostand P., Baev V.K., Latypov A.F., Shumsky V.V., Yaroslavtsev M.L. Integration of an hypersonic airbreathing vehicle: assesment of overall aerodinamic perfomances and of uncertaities. AIAA, 6-th International Aerospace Planes and Hypersonic Technologies Conference, Chattanooga, 1995. AIAA-95-6100, 15 p.
13. Latypov A.F., Yaroslavtsev M.L, Zudov V.N. Application of pulse tube for the test of engines hypersonic aircraft // International Congress on Instrumentation in Aerospace Simulation Facilities: Proceedings. Dayton, Ohio, 1995. p. 185-195.
14. Aulchenko S.M., Latypov A.F. Constructing plane curves by means of parametric polynomials of the fourth degree // Сотр. Math. Phys. 1995. Vol. 35, No. 7. p. 913-915.
15. Chalot F., Rostand P., Perrier P., Goonko Y., Kharitinov A., Latypov A., Mazhul I., Yaroslavtsev M. Validation of global aero propulsive characteristics of integrated configurations. AIAA, 8-th International Space Planes and Hypersonic Systems and Technologies Conference. Norfolk, 1998. AIAA Paper, 1998, No. 98-1624,8 p.
16. Adamov N., Goonko Y„ Kliaritinov A., Latypov A., Mazhul I., Yaroslavtsev M. Chalot F., Rostand P., Perrier P. Study on drag-thmst forces of a scramjet model in blow-down and hot-shot wind tunnels // International Conference on the Methods of Aerophysical Research: Proceedings, Pt. 3. Novosibirsk, 1998. p. 3-9.
17. Аульченко C.M., Латыпов А.Ф., Никуличев Ю.В. Метод численного интегрирования систем обыкновенных дифференциальных уравнений с использованием интерполяционных полиномов Эрмита // Журнал вычислительной математики и математической физики. 1998. т. 38, № 10. с. 1665-1670.
18. Латыпов А.Ф., Фомин В.М. Опыт функционального моделирования воздушно-космических систем выведения грузов на околоземную орбиту // Устойчивость и управление для нелинейных трансформирующихся систем: Тезисы докл. Второй Международной конференции. М., 2000. с.14.
19. Аульченко С.М., Латыпов А.Ф., Никуличев Ю.В. Построение кривых и поверхностей с помощью параметрических полиномов //Автометрия. 2000. К» 4. с. 60-76.
20. Аульченко С.М., Латыпов А.Ф., Никуличев Ю.В. Семейство одношаговых А- и ¿-устойчивых алгоритмов решения задачи Коши для систем обыкновенных дифференциальных уравнений // Математические модели и методы их исследования. Изд. Института вычислительного моделирования СО РАН, Красноярск, 2001. Т. 1. с. 54-58.
21. Латыпов А.Ф. Оценка энергетической эффективности подвода тепла перед телом в сверхзвуковом потоке газа // Восьмой Всероссийский съезд по теоретической и прикладной механике: Тезисы докл. Пермь, 2001. с. 392.
22. Латыпов А.Ф., Фомин В.М. Оценка энергетической эффективности подвода тепла перед телом в сверхзвуковом потоке// Прикладная механика и техническая физика. 2002. т. 43, № 1. с.71 -75.
23. Замураев В.П., Калинина А.П., Латыпов А.Ф. Оценка тяги ПВРД при импульсном подводе энергии // Теплофизика и аэромеханика. 2002. т. 9, №3. с. 405-410.
24. Goonko Y.P., Latypov A.F., Mazhul I.I., Kharitonov A.M., Yaroslavtsev M.I., Rostand P. Structure of flow over a hypersonic inlet with side compression wedges II AIAA Journal. 2003. V. 41, N 3. p. 436448.
25. Kalinina A.P., Latypov A.F., Zamuraev V.P. Modeling of nonstationary flow in variable cross section flat duct at a distributed pulse periodical energy supply // European Conference for Aerospace Sciences (EUCASS): Proceedings, Moscow, 2005. CD-ROM. 6 p.
26. Latypov A.F. Dynamic method used for determining the aerodynamic forces effecting the models tested in short duration wind tunnel // European Conference for Aerospace Sciences (EUCASS): Proceedings, Moscow, 2005. CD-ROM. 6 p.
27. Латыпов А.Ф. Динамический метод определения аэродинамических характеристик моделей по результатам экспериментов в аэродинамических трубах кратковременного действия // Прикладная механика и техническая физика. 2006. т. 47, №5. с.47-55.
28. Замураев В.П., Калинина А.П., Латыпов А.Ф. Моделирование нестационарного течения в канале переменного сечения при распределенном импульсно-периодическом подводе энергии // Механика жидкости и газа. 2006. №2. с. 149-156.
29. Латыпов А.Ф., Никуличев Ю.В. Численные методы решения задачи Коши для жестких систем обыкновенных дифференциальных уравнений на основе многозвенных интерполяционных полиномов Эрмита // Журнал вычислительной математики и математической физики. 2007. т. 47, №2. с. 234-244.
30. Latypov A.F. Numerical methods for solving the linear integral Fredholm and Volterra equations of the first kind // Международная конференция «Обратные и некорректные задачи математической физики», посвященная 75-летию академика М.М. Лаврентьева. Новосибирск, 2007. CD-ROM, 5 р.
31. Латыпов А.Ф. Оценка энергетической эффективности подвода тепла перед телом при полегге с ускорением. Часть 1. Математическая модель // Теплофизика и аэромеханика. 2008. т. 15, №4. с. 573-584.
32. Латыпов А.Ф. Оценка энергетической эффективности подвода тепла перед телом при полете с ускорением. Часть 2. Математическая модель разгонного участка траектории. Результаты расчетов // Теплофизика и аэромеханика. 2009. т. 16, №1. с. 1-12.
33. Латыпов А.Ф. Численное моделирование течения в канале переменной площади сечения при импульсно-периодическом подводе энергии // Прикладная механика и техническая физика. 2009. т. 50, №1. с. 3-11.
34. Латыпов А.Ф. Эксергетическнй анализ прямоточного воздушно-реактивного двигателя // Теплофизика и аэромеханика. 2009. т. 16, №2. с. 319-330.
35. Патент РФ №2347098. Способ работы сверхзвукового пульсирующего прямоточного воздушно-реактивного двигателя и сверхзвуковой пульсирующий прямоточный воздушно-реактивный двигатель / Латыпов А.Ф., Фомин В.М. // 20.02.2009. Заявка №2007122144, приоритет от 13.06.2007.
Ответственный за выпуск А.Ф. Латыпов
Подписано в печать 20.07.2009 Формат бумаги 60 x84/16, Усл. печ. л. 2.0, Уч.-изд. л. 2.0, Тираж 150 экз., Заказ №9
Отпечатано на ризографе ЗАО "ДОКСЕРВИС" 630090, Новосибирск, Институтская, 4/1
Введение
Глава 1. Численные методы.
1. Численные методы решения задачи Коши для жестких систем обыкновенных дифференциальных уравнений на основе многозвенных интерполяционных полиномов Эрмита.
1.1. Постановка задачи.
1.2. Интерполяционный полином Эрмита по трём точкам.
1.3. Семейство ZMR-методов.
1.3.1. LMR(000s) метод. Квазилинейная система ОДУ.
1.3.2. ZA4R(101j) метод. метод.
1.4. Устойчивость методов.
1.5. О решении систем алгебраических уравнений.
1.6. Расчет локальной и глобальной ошибки.
1.7. Интегрирование системы ОДУ с разрывными функциями правых частей.
1.8. Тестовые расчеты. Замечания.
2. Метод интервального осреднения для определения квазирешения линейных интегральных уравнений Вольтерра 1-го рода с разностным ядром.
2.1. Аппроксимация G(t) кусочно-постоянными функциями. Алгоритм минимизации функционала невязки уравнений.
2.2. Аппроксимация G(t) кусочно-линейными функциями.
2.3. Модельная задача.
3. Динамический метод определения аэродинамических характеристик моделей по результатам экспериментов в аэродинамических трубах кратковременного действия.
3.1. Динамический метод.
3.2. Статическая тарировка.
3.3. Исходные уравнения.
3.4. Определение нормальных реакций (динамическая тарировка).
3.5. Определение сил и моментов.
3.6. Эталонная модель НВ-2.
3.7. Модель аэрокосмического демонстратора ARES.
3.8. Модель EXPERT.
4. Оценка коэффициента штрафа в методе штрафных функций для решения задачи математического программирования.
Выводы к главе 1.
Глава 2. Эксергетический метод оценки характеристик прямоточного воздушно- реактивного двигателя.
1. Эксергия термодинамической системы.
2. Одномерное стационарное течение.
Условие подвода заданного количества энергии.
3. Математическая модель ПВРД.
4. Т0,е —диаграмма. Выводы к главе 2.
Глава 3. Численное моделирование нестационарного течения в канале переменного сечения при распределенном импульсно-периодическом подводе энергии.
3.1. Квазиодномерное нестационарное течение. Одномерное стационарное течение. Условия перехода через скорость звука. Результаты расчетов квазиодномерных нестационарных течений. Конфигурация сверхзвуковой камеры сгорания.
3.2. Двухмерное нестационарное течение. Выводы к главе 3.
Глава 4. Оценка энергетической эффективности подвода тепла перед летательным аппаратом при сверхзвуковом полете.
4.1. Математическая модель процесса.
4.2. Математическая модель аэродинамической поляры.
4.3. Силовая установка (СУ).
4.4. Удельный импульс при выработке энергии
4.5. Крейсерский режим полета.
4.6. Уравнения движения на разгонном участке траектории полета.
4.7. Угол наклона траектории.
4.8. Вариация m.
4.9. Построение оптимальной траектории (начальная формулировка).
4.10. Вариация qx.
4.11. Угол атаки аА.
4.12. Угол наклона траектории при аА = const.
4.13. Вариация nv.
4.14. Конечная формулировка задачи. Результаты расчета. Выводы к главе 4.
Выводы.
Предстоящее интенсивное освоение Космоса потребует доставки на околоземную орбиту большого количества грузов. При последующем создании на орбите разнообразных производств будет большим также обратный поток продукции. Существующие ракетные системы не в состоянии обеспечить перемещение грузов в больших объемах. Причины: высокая стоимость, длительное время стартовой подготовки, малое количество стартовых комплексов. Необходимо создание принципиально новых систем выведения грузов на околоземную орбиту. Такие системы могут быть созданы на основе воздушно-космических самолетов (ВКС) с комбинированной силовой установкой, включающей прямоточный воздушно-реактивный двигатель (ПВРД) на водороде и жидкостный ракетный двигатель (ЖРД). Использование воздуха для создания подъемной силы и атмосферного кислорода для окисления топлива на части траектории разгона позволит значительно уменьшить затраты топлива и стартовую массу ВКС.
Преимущества ВКС перед ракетными системами: горизонтальный старт с любого аэродрома; существенно меньшая стоимость доставки на орбиту 1 кг груза; эффективное достижение любых околоземных орбит; малое время подготовки к старту; возможность проведения быстрых спасательных операций; малая зона отчуждения.
Исследования по разработке технологии гиперзвукового полета с ПВРД на водороде ведутся в ряде зарубежных стран в течение последних -50 лет: США («DCR», «GDE», «HXFE»), Франция, Германия («JAPHAR»), Франция («PROMETHEE»), Франция, Россия («WRR»), Япония («Е1», «ATREX»), Англия («Hotol»), Китай, Австралия. Подобные исследования проводились также в СССР («Спираль», «Буран»). Обзор рассматриваемых в мире проектов приведен, например, в [0.1].
Для создания ВКС, несмотря на значительные достижения, необходимо решить множество проблем. Первыми в этом ряду являются взаимосвязанные проблемы ПВРД и аэродинамического качества конфигурации летательного аппарата (ЛА), т.к. характеристики силовой установки и аэродинамические характеристики компоновки JIA определяют главным образом необходимые затраты топлива для выведения на орбиту. По оценкам, используя достигнутые в экспериментах на моделях значения аэродинамического качества конфигураций JIA и удельного импульса ПВРД, масса горючего, необходимого для разгона ВКС до 1-ой космической скорости, составляют около 70% от его стартовой массы [4.11]. Такая же оценка следует из уравнения dm mV° Г (l-w2)cos& dw I, e 2 пЖ
0.1)
A(w) - -——-,sin6 =2,(w)-nv,mQ =1 l + c)r» где m = m/m0 - относительная масса ВКС, w=V/V° - относительная скорость полета, V0 - орбитальная скорость, 1е - эффективный удельный импульс силовой установки, К — аэродинамическое качество, nv — продольная перегрузка, 0 - угол наклона траектории, — показатель адиабаты. Вывод уравнения приводится в главе 4. Траектория разгона разделяется на два участка. На 1-ом участке аппарат разгоняется до скорости Wj = 0.7 с использованием аэродинамической подъемной силы и силовой установки, включающей последовательно ВРД, ПВРД, ПВРД+ЖРД. На 2-ом участке ВКС разгоняется до конечной скорости w2 = 1 с использованием только ЖРД. Расход массы на этом участке определяется из формулы Э.К. Циолковского, которая, в частности, следует из уравнения (0.1) в предположении отсутствия внешних сил.
Из соотношения массового баланса получаем выражение для стартовой массы т0 аппарата ш0=--(0.2)
1 - mT - mK
Здесь mpn - масса полезной нагрузки, тт, Шк - относительные массы горючего и конструкции. Из (0.2) следует, что на массу конструкции накладываются весьма жесткие требования, и значение стартовой массы очень чувствительно к вариации относительной массы горючего Sm0 8тг то тРп
Уменьшение rriT позволяет уменьшить ш0 и ослабить требования к конструкции.
Большая относительная масса конструкции допускается для многоступенчатых систем, в частности, при сбросе по траектории полета отработавших элементов конструкции. Однако при этом усложняются условия эксплуатации, и увеличивается стоимость.
Уменьшение расхода горючего может быть осуществлено увеличением аэродинамического качества ВКС и удельного импульса силовой установки. Многочисленные экспериментальные исследования аэродинамических характеристик гиперзвуковых летательных аппаратов (ГЛА) свидетельствуют о том, что их максимальное аэродинамическое качество в гиперзвуковом диапазоне скоростей составляет около Кт « 4 (рис.0.1) [0.2]. При реальных числах Рейнольдса Кт «6. Это значение не удаётся увеличить посредством аэродинамического конструирования конфигураций ГЛА. Поэтому в настоящее время значительное внимание уделяется решению задачи активного управления обтеканием тел посредством энергетического и/или силового воздействия на набегающий поток, в частности, посредством подвода тепла перед телом в сверхзвуковом потоке. Изучению этой проблемы посвящено значительное число работ (например, [0.3-0.8]). Для технической реализации предполагается использование лазерного и СВЧ-излучения, электрического разряда. В большинстве теоретических и экспериментальных исследованиях изучается задача уменьшения аэродинамического сопротивления. Эффект уменьшения аэродинамического сопротивления связывается, главным образом, с уменьшением плотности газа в набегающем потоке, что подтверждается расчетами [0.09-0.11] и непосредственными измерениями [0.12-0.14]. Дополнительные эффекты возможны из-за изменения режима обтекания вследствие уменьшения числа Маха, изменения числа Рейнольдса, ионизации потока. В [0.15] на примере обтекания гиперзвуковым потоком газа трапециевидного профиля показано значительное влияние ступенчатого распределения температуры в набегающем потоке на аэродинамическую подъемную силу. Установлено, что при условии максимального аэродинамического качества оптимальным является режим глиссирования.
Традиционно эффективность подвода тепла в установившемся полете оценивается величиной
77 = AzA (0.3) где А0 - исходная мощность двигателя, А - мощность при тепловом воздействии, Qw — мощность подвода тепла [0.6-0.8, 0.16-0.18]. Эффективность, оцениваемую выражением (0.3), в приближении бесконечного теплового следа можно представить в виде т7 = ^^cx0Mi. (0.4)
Коэффициент аэродинамического сопротивления cxq определен по миделю летательного аппарата (ДА). В (0.4) не учитываются полный энергетический баланс, функциональное назначение летательного аппарата и существенные параметры процесса: степень нагрева газа, несущие свойства JIA, характеристики двигателя, эффективность преобразования энергии топлива в энергию излучения.
Целесообразность использования такого способа управления обтеканием может быть установлена методом функционального моделирования [0.19,4.8-4.10]. В соответствии с методом JIA, как сложная иерархическая система, представляется в виде взаимосвязанной совокупности подсистем, определяемых по функциональным признакам. Математическая модель (MM) JIA строится так, что отражаются функциональные характеристики и связи элементов и объекта в целом без жесткой привязки к конкретным реализующим устройствам. Последующее изучение такого типа ММ позволяет установить принципиальную возможность достижения цели и определить характеристики эффективности, критические режимы работы, влияние вариаций характеристик элементов на функциональные свойства системы. Решение задач оптимальной компенсации возмущений, накладываемых на базовые значения характеристик элементов, позволяет выработать требования к точности определения характеристик.
Важная особенность функционального моделирования заключается в том, что синтез и анализ объекта производятся при незначительном объеме начальной информации. Это определяет, во-первых, итерационный характер построения ММ, проявляющийся как в количестве включаемых подсистем, так и в точности их описания; и, во-вторых, необходимость построения ММ с минимальным числом задаваемых входных параметров и функциональных характеристик элементов, что уменьшает степень неопределенности определения характеристик ДА. Второе обстоятельство стимулирует поиск обобщенных представлений функциональных свойств элементов. Эти представления, естественно, должны коррелировать с предполагаемым множеством реализующих устройств, однако необходимо различать содержание работы и ее реализацию. Выбор конкретных устройств, их разработка - следующий этап.
Весьма существенна роль учитываемых функциональных ограничений, совокупность которых во многом способствует синтезу технически реализуемых объектов.
В данной работе MM JIA состоит из следующих функциональных блоков: аэродинамических характеристик, тяги и удельного импульса двигателя, траектории полета, функциональных ограничений, решения задачи оптимального управления.
Аэродинамические характеристики определяются в предположении, что полет JIA происходит на границе раздела сред высокой и низкой плотности.
Для определения характеристик ПВРД существующими методами требуется задание некоторого множества определяющих величин, зависящих от газодинамических и геометрических параметров и определяемых, как правило, экспериментально (напр., [0.20, 0.21]). Поэтому эти методы малопригодны при функциональном моделировании. При функциональном моделировании необходимо установить минимальную совокупность определяющих параметров, относительно слабо меняющихся в процессе функционирования системы и с определяемыми границами возможных значений. Ее успешное решение во многом определяет качество получаемых результатов. Здесь для расчета характеристик ПВРД применяется эксергетический метод [0.22 -0.24, 2.3], в основе которого лежит оценка потерь работоспособности газа вследствие необратимости процессов в элементах двигателей. С этой целью для установления приращения энтропии в каком-либо элементе двигателя определяются верхняя и эталонная оценки. Для максимального значения приращения энтропии получено аналитическое выражение; для определения эталонного значения строится базовый процесс. Рабочее значение приращения энтропии в элементе определяется как взвешенная сумма этих оценок. Функциональные зависимости этих оценок от входных и внешних условий дают основания для предположения о слабой зависимости весовых коэффициентов от режима работы двигателя.
В рамках такого описания построена функциональная ММ силовой установки, которая позволяет получать оценки коэффициента тяги и удельного импульса ПВРД и комбинации ракетного и прямоточного двигателей с возможным отбором энергии для целей управления внешним обтеканием. В состав ММ ПВРД входят следующие элементы: воздухозаборник, камера сгорания, сопло.
В [4.13, 4.14] впервые предпринята попытка оценки эффективности подвода тепла перед телом в сверхзвуковом потоке с учетом выше перечисленных факторов. Для получения подводимой мощности увеличивался расход топлива, и вводился коэффициент преобразования энергии топлива в энергию излучения. Однако оценка этого коэффициента весьма затруднена. Эквивалентное подводимой энергии уменьшение эксергии продуктов сгорания топлива устраняет эту неопределенность. В [4.15, 4.16] выполнены оценки увеличения коэффициента дальности Бреге крейсерского полета гиперзвукового JIA и уменьшения затрат топлива при разгоне ВКС с использованием комбинированного двигателя, включающего ПВРД и жидкостный ракетный двигатель (ЖРД).
Сравнение расходов топлива на разгон с нагревом воздуха перед JIA и без нагрева должно производиться на оптимальных траекториях полета, которые описываются системой обыкновенных дифференциальных уравнений (ОДУ) с функциями управления. Для решения задач оптимального управления широко используется метод Ритца, состоящий в замене исходной задачи некоторой вложенной последовательностью конечномерных задач минимизации, решения которых образуют минимизирующую последовательность [0.25].
При численном решении задач оптимизации необходимо многократное вычисление функционала. Поэтому разработка алгоритмов с возможно малым числом операций всегда будет важной. Наиболее ресурсоёмкими по количеству операций являются задачи решения дифференциальных уравнений и их систем. По численным методам решения систем ОДУ написано множество книг. Укажем лишь некоторые, суммирующие достижения в этой области. Обстоятельное изложение теоретических и практических аспектов численного интегрирования систем ОДУ содержит коллективная монография под редакцией Дж. Холла и Дж. Уатта [0.26]. Вопросы устойчивости, точности, сходимости достаточно большого количества методов проанализированы в монографии X. Штеттера [0.27].
Автор согласен с утверждением в [0.28], что чем более адекватна математическая модель физическому процессу, тем более жесткими являются входящие в её состав дифференциальные уравнения, что указывает на важность совершенствования методов их решения.
Понятие жесткости в применении к системам ОДУ, по видимому, впервые было введено в работе [0.29] в 1952 году. Качественно сущность этого явления можно объяснить тем, что в уравнениях требуется учитывать в каждой точке отрезка наблюдения одновременно несколько величин, скорость изменения которых существенно различаются по величине, о чем свидетельствует плохая обусловленность матрицы Якоби. Решение таких систем затруднено накоплением ошибок в процессе вычислений. Ошибки возникают из-за представления чисел в компьютере конечным числом значащих цифр, наличия, как правило, итерационных процедур, дискретизации задачи конечным числом элементов и т.п. Ограниченное быстродействие ЭВМ также требует разработку экономичных методов.
Упомянем некоторые популярные методы интегрирования жестких систем ОДУ. Д. Дальквист [0.30, 0.26] ввел понятие А - устойчивости методов и предложил способ для исследования метода на А -устойчивость. Б.Л.Эль [0.31] ввел понятие L — устойчивости. L — устойчивый метод предоставляет более широкие возможности для увеличения шага интегрирования при сохранении заданной точности интегрирования. Г. Розенброк [0.32] описывает класс новых неявных А -устойчивых методов. Б.Л. Эль [0.33], С.С. Артемьев и Г.В. Демидов [0.34] развивают и совершенствуют методы типа Г. Розенброка. С.О. Фатунла [0.35] приводит описание нового метода с использованием экспоненциальных аппроксимаций, С.В.Гир [0.36] предлагает, так называемый, BDF - метод или метод обратного дифференцирования, использующий аппроксимацию производной искомой функции на шаге интегрирования в виде линейной комбинации значений функции в к точках, где к - порядок метода. Е. Хайрер и Д. Ваннер [0.37] публикуют Фортран-программу для интегрирования жестких систем ОДУ, основанную на неявном методе Рунге-Кутты 5-го порядка. Этот краткий перечень методов показывает, что исследователи продолжают работу по совершенствованию методов численного интегрирования жестких систем ОДУ.
Существует также проблема расчета глобальной ошибки при численном интегрировании систем ОДУ. Теоретические вопросы, связанные с расчетом глобальной ошибки, обсуждаются в [0.27, §1.5], однако алгоритмы не приводятся. Алгоритмический аспект не рассматривается и в упомянутой монографии [0.26]. Автору не удалось найти описания практически используемого метода.
Правые части системы ОДУ, описывающей траекторию полета ВКС, являются разрывными функциями. Разрывы функций обусловлены включением в некоторой точке траектории ЖРД, выключением в другой точке ПВРД и возможным достижением предельного угла атаки. В [0.38] дается качественная теория системы ОДУ с разрывной правой частью. В [0.39] исследуется поведение явных методов Рунге-Кутты при решении систем обыкновенных дифференциальных уравнений с разрывной правой частью. Однако в указанных работах алгоритмы решения не приводятся.
Сравнительная оценка эффективности методов является не простой задачей. Достаточно объективно можно оперировать только с опубликованной программной реализацией, т.к. даже для одного численного метода могут быть различные по своей эффективности программные реализации из-за различий, например, в методах решения вспомогательных задач, в стратегии выбора шагов интегрирования. Обычно относительная эффективность оценивается сравнительными расчетами для нескольких тестовых задач, моделирующих требуемые свойства решений. Количество вычислений правой части системы и якобиана, количество и величины шагов интегрирования достаточно объективно характеризуют эффективность метода.
Математическая модель - численные методы - алгоритмы и программы — компьютерное обеспечение - циклически взаимосвязанная система. Увеличение мощности компьютера позволяет усложнять, увеличивать степень адекватности математических моделей, что, в свою очередь, заставляет совершенствовать численные методы и алгоритмы, а новые математические модели и достижения в численных методах выдвигают новые требования к мощностям вычислительных систем и т.д.
Поскольку основой большинства математических моделей являются дифференциальные уравнения и системы, разработка новых и совершенствование имеющихся численных методов интегрирования систем ОДУ с учетом роста возможностей вычислительных средств была, есть и, вероятно, долго будет чрезвычайно актуальной задачей.
В диссертации приведено новое семейство методов численного решения систем ОДУ, которые основаны на аппроксимации функций правых частей на шаге h трехточечными интерполяционными полиномами Эрмита. Эти методы просты, легко программируемы и гибки по выбору параметров, что позволяет «подгонять» их к особенностям конкретной задачи. И особенно важно, эти методы пригодны для решения жёстких систем ОДУ. В рамках метода оказалось возможным построение алгоритма расчета глобальной ошибки при интегрировании систем ОДУ. Разработан также алгоритм интегрирования систем ОДУ с разрывной правой частью.
Идея прямоточного воздушно-реактивного двигателя для полета с гиперзвуковыми скоростями (числа Маха полета М > 6 - 7 ) предполагает сгорание топлива в сверхзвуковом потоке воздуха в канале. При этом количество сгорающего топлива должно быть достаточным для получения требуемой тяги (коэффициент избытка воздуха а ~ 1). В [0.40, 0.41] предложены схемы возможных способов впрыска топлива в поток и возникающих при этом течений. Однако данные, свидетельствующие об их реализации, отсутствуют. Отмечается, что процессы смешения и горения должны рассматриваться совместно.
Теоретическим и экспериментальным исследованиям ПВРД в течение более 40 лет посвящено значительное число работ. ЦАГИ и ЦИАМ регулярно выпускаются обзоры выполненных работ [0.42]. В [0.43] дан обзор исследований, выполненных в России, в [0.44] - в NASA.
В многочисленных экспериментах в аэродинамических трубах моделей ПВРД при сгорании топлива наблюдается значительная перестройка начального сверхзвукового потока. В переходном процессе увеличивается давление не только в камере сгорания, но и в слабо расширяющемся канале, примыкающем к выходному сечению воздухозаборника. В экспериментах [3.4-3.7] на модели двигателя (рис.0.2) с размерами (102 х 154x860) при числах Маха набегающего потока М00=6-8 время переходного процесса составляло 50-60 мс.
Дозвуковой поток в камере сгорания формируется раньше. На рис.0.3 приведены графики типичных зависимостей от времени относительных давлений в характерных точках воздухозаборника при горении водорода с коэффициентом избытка воздуха а ~ 1 при числе Маха набегающего потока Мод = 6.
На рис.0.4 показано распределение давления по тракту модели для различных моментов времени для этих же условий, характерное для псевдоскачкового режима горения [0.45, 0.46]. Давление воздуха в изоляторе постепенно повышается в системе косых скачков уплотнения, генерируемых в процессе сжатия воздуха в воздухозаборнике и отрывом пограничного слоя. Система косых скачков замыкается прямым скачком уплотнения. Горение топлива происходит в дозвуковом потоке, давление при этом уменьшается. Такой режим течения в канале с подводом тепла наблюдается во многих экспериментальных работах (например, [0.47 — 0.53]). В [0.48] эксперименты проводились при числе Маха во входном сечении канала М=2.75. Псевдоскачковый режим горения наблюдался даже при коэффициенте избытка воздуха а«2, полнота сгорания колеблется в пределах \j/ = 0.6 ч- 0.8.
В [0.54] приведены схема течения (рис. 0.5), распределение числа Маха в канале (рис. 0.6) и распределение давления (рис. 0.7). Некоторые результаты летных испытаний Х-43А при числе Маха полета 10 приведены в [0.55, 0.56]. На рис. 0.8 показана функция распределения давления по каналу двигателя без указания масштабов на осях координат. Идентифицировать структуру течения затруднительно.
Диагностика потоков со сгоранием топлива чрезвычайно трудна вследствие неравномерности распределения параметров течений и неравновестности процессов. Распределения по длине канала статического давления и теплового потока, получаемые в эксперименте, недостаточны для идентификации течений.
Таким образом, подводя итоги краткого рассмотрения результатов в данной области, отметим следующее. Отсутствуют достоверные экспериментальные результаты, свидетельствующие о сохранении сверхзвукового течения в канале при подводе энергии с эквивалентным коэффициентом избытка воздуха а = 1, тем более, при ограничении статической температуры продуктов сгорания. Это ограничение важно при гиперзвуковых числах Маха полета и связано с ограничением степени диссоциации продуктов сгорания, т.к. диссоциация уменьшает эксергию потока газа. Необходимо определить условия, при которых было бы возможным организовать подвод тепла с учетом названных факторов. Необходимо также определить условия и факторы, влияющие на формирование течения определенной структуры.
В классической схеме прямоточного воздушно-реактивного двигателя энергия сгорания топлива подводится к потоку газа в некотором политропном процессе. В данной работе приводятся результаты расчетов одномерного и двумерного нестационарного течений совершенного газа в канале переменного сечения при импульсно-периодическом режиме подвода энергии в локальных зонах. В этой модели процессы смешения отсутствуют. Задача смешения сверхзвуковых реагирующих потоков газа требует отдельного изучения. Исключение из рассмотрения процессов смешения позволило определить непосредственное влияние параметров подводимой энергии (мощности, частоты импульсов, распределения источников по длине канала) на характеристики течения. Математическое моделирование таких течений может восполнить недостающую информацию и установить возможные причины описанной выше наблюдаемой экспериментально перестройки начального сверхзвукового потока в канале при подводе энергии.
Экспериментальные исследования аэрогазодинамических характеристик моделей гиперзвуковых летательных аппаратов часто проводятся в аэродинамических трубах кратковременного действия. При этом для определения силовых характеристик используются тензометрические весы. Для определения переменных во времени действующих на модель сил и моментов применяются следующие методы обработки результатов измерений: 1. метод осреднения; 2. аналитический метод; 3. статистический метод; 4. упрощенный статистический метод [0.58]. При использовании метода осреднения требуется достаточно длительное время работы трубы, чтобы успели затухнуть колебания, вызванные начальными ударными нагрузками, а также достаточно большой промежуток времени, в течение которого параметры потока не меняются. Однако даже при выполнении этих условий сохраняются неконтролируемые систематические ошибки. В трех других методах реализуется принцип компенсации инерции. Для этого в некоторых точках дополнительно измеряются ускорения акселерометрами. Методы основаны на следующих предположениях: модель представляет собой жесткое тело; колеблющаяся система состоит из модели и некоторой части весов, которую нужно определить; произведения скоростей вращения вокруг осей системы достаточно малы, и ими можно пренебречь по сравнению с угловыми ускорениями.
В аналитическом методе основная трудность заключается в установлении части системы, участвующей в движении, с последующим определением массы, центра массы и матрицы моментов инерции. Точность восстановления действующих аэродинамических нагрузок зависит также от координат мест расположения акселерометров.
Статистический метод основан на линеаризации соотношений между силами инерции и ускорениями, измеряемыми акселерометрами, с использованием следующего дополнительного предположения: на выбранном интервале времени, определяющем рабочий режим аэродинамической трубы, аэродинамические коэффициенты постоянны. При этом количество уравнений (6) меньше количества неизвестных (42). Но поскольку при испытаниях система уравнений справедлива в любой момент времени она решается на некотором интервале методом наименьших квадратов.
В упрощенном статистическом методе предполагается линейная зависимость между ускорением и соответствующей силой, что приводит к уменьшению количества неизвестных. Используются три акселерометра, по одному в каждом направлении. Оси измерений совпадают с направлениями связанной системы координат. Неизвестные коэффициенты пропорциональности зависят от геометрических и массовых характеристик модели, которые необходимо определять в каждом эксперименте для модели, установленной в рабочей части аэродинамической трубы. Система возбуждается импульсным ударником, что вызывает ее колебания вокруг нулевой линии без воздействия аэродинамических сил. По результатам измерений ускорений державки и весовых измерений сил и моментов находятся коэффициенты, которые используются в дальнейшем при определении аэродинамических характеристик на некотором интервале.
Решение, получаемое статистическими методами, существенно зависит от положения и протяженности временного интервала. Предположение о постоянстве аэродинамических коэффициентов может не выполняться на необходимой для решения задачи длине интервала. Требуется также согласование периодов колебаний системы с промежутком времени, в течение которого поток стационарен, что накладывает дополнительное ограничение на геометрические и массовые характеристики модели.
В данной работе предложена методика, учитывающая динамику модели и непостоянство параметров потока. В методике отсутствуют указанные выше недостатки. Предполагается, что регистрируемые во времени реакции тензометрических весов описываются обыкновенными дифференциальными уравнениями с постоянными коэффициентами, т.е. система «модель-державка-тензометрические весы — система крепления» является линейным динамическим объектом. В этом случае, получив экспериментально реакции системы на единичные нагрузки, определяется система интегральных уравнений Вольтерра 1-го рода с разностным аргументом для восстановления действующих нагрузок по регистрируемым во времени реакциям тензометрических весов. Поскольку получаемая информация дискретна во времени и задана с ошибкой, то задача восстановления является некорректной. Ее решение существенно осложняется тем, что в начальный момент матрица нормальных реакций и ее производная по времени равна нулю. Для решения некорректных задач используются методы регуляризации [0.59]. Достаточно полный перечень методов и их описания изложены в [0.60]. Однако эффективного метода для решения названной задачи не представлено. В данной работе предлагается метод, основанный на условии равенства нулю средних значений невязок уравнений на интервалах, на которые разбивается область определения функции правой части. Это условие позволяет получать решение в выбранном классе функций без использования стабилизирующего функционала.
Диссертация состоит из введения, четырех глав, выводов, списка литературы и рисунков. Краткое содержание глав.
Выводы.
1. Выполнен анализ квазиодномерного и двухмерного квазистационарных течений в канале переменного сечения, описываемых уравнениями Эйлера и формирующихся при импульсно-периодическом энергетическом воздействии при больших числах Струхаля. Течения устойчивы в среднем на периоде. Получены условия, определяющие структуру течений: максимально допустимое значение энтропии для каждого сечения канала и условие перехода через скорость звука в квазиодномерном случае при подводе энергии и наличии диссипации кинетической энергии.
Получено, что при больших значениях числа Струхаля устанавливается периодический режим течения с малыми амплитудами колебаний параметров. Так как при этом энергия подводится при постоянном объеме, то этот режим обеспечивает максимальное значение эксергии потока и, следовательно, тяги двигателя.
Предложена конфигурация канала, в котором подвод тепла к сверхзвуковому потоку осуществляется с учетом ограничения статической температуры газа. Импульсно-периодический подвод энергии в таком канале позволяет увеличить число Маха полета до значений, при которых возможно использование прямоточного канала в составе комбинированного двигателя для увеличения эффективного удельного импульса.
В результате численного моделирования нестационарного двухмерного течения в канале переменного сечения при подводе тепловой энергии в локальных зонах в импульсно-периодическом режиме получена экспериментально наблюдаемая перестройка начального сверхзвукового течения, определяемая условием подвода заданного количества энергии.
2. Разработаны новые численные методы для решения следующих задач:
2.1. Построено семейство А-, L- и Ь{б)~ устойчивых методов решения задачи Коши для жестких систем обыкновенных дифференциальных уравнений (ОДУ), основанные на представлении правых частей системы на шаге h в виде трёх точечных интерполяционных полиномов Эрмита (LMR(L,M,R,s) - алгоритмы, s е [0.5,1.0)- координата внутренней на интервале h точки коллокации). Погрешность методов ~ hL+M+R+4, Дано определение Ь{д) - устойчивости одношаговых методов с малым параметром д. Для методов LMRiL-l^X) определены условия А— и L- устойчивости, для методов LMi?(0,0,0, s),LMR(\$,\, s), LMR(l,l,l35)-условия А- и L(S)~ устойчивости. Разработан алгоритм расчета глобальной ошибки и алгоритм решения задачи Коши для систем ОДУ с разрывными правыми частями.
2.2. Разработан метод решения системы линейных интегральных уравнений Вольтерра 1-го рода с разностным аргументом для варианта, когда исходная информация (значения измеряемых функций и ядра) заданы в дискретных точках с известной ошибкой. Решение определяется в классе кусочно-постоянных и кусочно-линейных функций с использованием условия равенства нулю средних значений невязок уравнений на интервалах, на которые разбивается область определения решения. Число интервалов и распределение их длин определяются посредством минимизации среднеквадратичной невязки уравнений.
2.3. Разработана методика восстановления действующих нагрузок в классе кусочно-постоянных функций при испытаниях моделей в аэродинамических трубах кратковременного действия. Приведены примеры решения задач по определению аэродинамических характеристик эталонной модели НВ-2, демонстратора Ares и модели Expert по результатам испытаний в аэродинамической трубе АТ-303 ИТПМ СО РАН.
2.4. Для решения задач математического программирования методом штрафных функций из условия локального минимума вспомогательного функционала получена оценка для коэффициента штрафа к~е.
2.4. Созданы комплексы программ для решения названных классов задач.
3. Разработан эксергетический метод оценки характеристик и анализа ПВРД. Для графического отображения возможных схем подвода тепла в канале ПВРД предложена диаграмма в координатах "полная температура— эксергия". Получено выражение для изменения эксергии в термодинамической системе при подводе тепла и наличии необратимых процессов.
4. Разработана методика оценки эффективности подвода тепла перед JIA при полете со сверхзвуковой скоростью. Полет происходит на границе раздела сред различной плотности (режим глиссирующего полета). Показана значительная эффективность такого способа управления обтеканием JIA как при крейсерском полете, так и при полете с ускорением.