CC BY
131
Мензульский Сергей Юрьевич
аспирант кафедры «Аэрокосмические системы» (МГТУ им. Н.Э. Баумана)
УДК 629.7.015.4:629.016.55
Методика расчета динамических аэроупругих характеристик конструкции, обтекаемой гиперзвуковым потоком воздуха
С.Ю. Мензульский
Описывается методика определения аэроупругих характеристик конструкции гиперзвуковых летательных аппаратов, основанная на использовании трехмерного конечно-элементного моделирования и подходов «модифицированной поршневой теории». Проведена верификация данной методики с результатами эксперимента по определению критической скорости флаттера конструкции.
Ключевые слова: динамика конструкции, гиперзвуковой летательный аппарат, аэроупругие колебания, флаттер.
The article describes a method for estimation of aeroelastic characteristics of a hypersonic flying vehicles construction, based on the use of three-dimensional finite-element modeling and «modified piston theory» approaches. The method was verified against the results of the experiment to determine a critical speed of the flutter.
Keywords: structural dynamics, hypersonic flying vehicle, aeroelastic vibrations, flutter.
Введение
середине прошлого века в мире произошло более 200 катастроф летательных аппаратов, связанных с самовозбуждающимися колебаниями конструкции (флаттером). Основной проблемой была неизученность аэроупругих автоколебательных процессов и отсутствие необходимых расчетных методик.
В 60-х гг. эта проблема по большей части была решена. Количество авиационных происшествий, связанных с флаттером, резко снизилось, однако аэроупругие автоколебания продолжают становиться причиной авиакатастроф (например, ТУ-134, ТУ-22, B-747, «Гжель»).
При освоении сверхзвуковых скоростей полета конструкторы и эксплуатанты вновь столкнулись с автоколебаниями, и из-за неприменимости дозвуковых методик этот вопрос потребовал дополнительного изучения. Благодаря наземной отработке и экспериментам на исследовательских моделях аэроупругие колебания не были серьезной проблемой при создании сверхзвуковых летательных аппаратов. Накопились сведения, позволившие создать их таким образом, что при эксплуатации существенных проблем с флаттером не возникло.
Рис. 1. Общие виды американских ГЛА X-33 и X-34
В настоящее время проводятся работы по созданию гиперзвуковых летательных аппаратов (ГЛА). Создаются и испытываются летающие лаборатории, проводится наземная отработка высокоскоростных двигателей [1], исследуется экономическая эффективность. Для предотвращения авиационных происшествий требуется достаточно простая и точная методика расчета аэроупругих характеристик конструкции в гиперзвуковом потоке воздуха.
Постановка задачи
Традиционно при решении задачи аэроупругости принимается гипотеза базовой плоскости. При этом вся упругая динамическая модель аппарата, как правило, сводится в одну расчетную плоскость. Аэродинамические нагрузки перпендикулярны базовой плоскости. Учитываются лишь поперечные колебания, перпендикулярные продольной оси аппарата.
Для ГЛА этот подход неприменим. Они имеют сравнительно небольшое удлинение, их конструкции могут иметь элементы, произвольно ориентированные в пространстве (рис. 1, [2]). При обтекании наклонных панелей поток поворачивается. Различные элементы конструкции могут совершать малые колебания в различных направлениях. Кроме поперечных возникают и продольные колебания соизмеримой амплитуды. Возникает связанная проблема анализа продольно-поперечных форм колебаний аппарата в потоке. Поэтому
при расчете аэроупругих колебаний ГЛА необходимо использование трехмерных аэродинамических и конечно-элементных моделей.
При расчете аэроупругих характеристик с помощью метода конечных элементов уравнение колебаний системы имеет вид [3]
[ М ]д +[ВА]д + «СаЖС КШ = 0,
(1)
где [М] — матрица масс конструкции; [^А] — матрица аэродинамического демпфирования; [СА] — матрица аэродинамической жесткости; [СК] — матрица жесткости конструкции; {#} — столбец обобщенных координат.
Матрицы [^А] и [СА] описывают влияние внешнего потока газа на упругую конструкцию летательного аппарата. Для описания сверхзвукового потока воздуха широко используется поршневая теория [4]. Согласно этой теории, аэроупругие силы в точке описываются выражением
Ар(х ,7 ) = к
ду(х ) ду(х )
ы
-+ V-
дх
Здесь
к =
2pv
л/м2 -1:
(2)
(3)
где р — плотность воздуха; М — число Маха; V — скорость набегающего потока; х — координата вдоль вектора скорости набегающего потока; у — координата поперек вектора скорости
ду( х ^)
набегающего потока;-— местный угол
дх
атаки (а).
Выражение (2) можно записать в виде Ар(х) = к^^ + ку дУ(Х)
д?
дХ
(4)
Из уравнения (4) наглядно следует, что член ду( Х ,1)
к-определяет аэродинамическое демп-
д?
фирование и входит в матрицу [^А], а член ду( Х ,1)
ку-— аэродинамическую жесткость и вхо-
дХ
дит в матрицу [СА] уравнения (1).
Следует отметить, что расчеты аэроупругих характеристик конструкции с помощью поршневой теории проводятся с применением гипотезы стационарности [5], суть которой заключается в допущении, что коэффициент к не зависит от малых деформаций конструкции в результате колебаний.
Использование поршневой теории в аэроупругих расчетах дает корректные результаты при величинах числа Маха 2 < М< 3. «При больших сверхзвуковых скоростях необходимо принимать во внимание эффекты, связанные с вязкостью, диссоциацией, ионизацией, фазовыми превращениями на границе газа и тела и другими эффектами. Однако локальный характер между возмущениями давления и скосом потока должен, разумеется, сохраниться» [4]. Таким образом, при М > 3 формула (2) останется без изменений, однако коэффициент к потребует коррекции с учетом описанных выше особенностей гиперзвукового потока.
В 80-х гг. прошлого века на предприятии НПО «Машиностроение» (г. Реутов Московской области) В.В. Никитенко, Ю.М. Ватрухи-ным иС.В. Куликом в качестве аэродинамического оператора была предложена модифицированная поршневая теория для расчетов аэроупругих характеристик летательного аппарата при 1,2 > М >2 [6]. Суть данного метода заключалась в том, что в уравнении (2) коэффициент к определялся исходя из результатов натурных продувок изделия по формуле
к = Лр( х) /
ду(х)
>-
дХ
(5)
За прошедшее время модифицированная поршневая теория показала свою адекватность
и хорошее совпадение расчетных аэроупругих характеристик с экспериментальными.
В данной методике было решено адаптировать подход модифицированной поршневой теории для определения аэроупругих характеристик гиперзвукового летательного аппарата (М > 5). Так как расчет аэроупругих характеристик ГЛА должен быть проведен до испытаний, то для определения коэффициента к должны быть использованы не результаты продувок, а результаты аэродинамических расчетов. Численное моделирование параметров гиперзвукового потока можно проводить в нескольких программных продуктах, например в среде ANSYS/Fluent. Этот пакет представляет собой программную реализацию численного интегрирования уравнений газовой динамики методом конечных объемов. Методика расчета основана на численном решении уравнений сохранения массы, импульса и энергии, записанных в интегральной форме. Следует отметить, что данная программа показывает хорошую сходимость результатов расчетов стационарного гиперзвукового обтекания с экспериментальными данными [7].
Программная реализация методики определения динамических аэроупругих характеристик ГЛА
В настоящее время ни один массовый программный пакет не может обеспечить расчет аэроупругих характеристик гиперзвукового летательного аппарата. В данной ситуации логичным является создание специальной программной реализации с максимальным использованием возможностей коммерческих пакетов (рис. 2).
Как сказано выше, система ANSYS/Fluent адекватно моделирует параметры гиперзвукового потока воздуха. Таким образом, расчетное определение коэффициента давления к1 корректно проводить в данном программном пакете.
Для задач анализа колебаний конструкции летательного аппарата хорошо себя зарекомен-
Рис. 2. Алгоритм методики по определению динамических аэроупругих характеристик гиперзвукового
летательного аппарата
довал пакет конечно-элементного моделирования MSC.Nastran. С его помощью удобно провести построение и верификацию упруго-массовой модели конструкции гиперзвукового летательного аппарата. Также в системе MSC.Nastran удобно провести задание распределения коэффициента давления по поверхности ГЛА. Однако интерфейс системы MSC.Nastran не позволяет использовать эти данные для составления матриц [СА] и [^А]. Поэтому составление и анализ уравнения (1) следует проводить в третьем программном пакете, обладающем достаточно гибким интерфейсом и мощными вычислительными возможностями.
Таким пакетом является язык математических вычислений МАТЬАВ. В нем после импортирования составленных системой MSC.Nastran матриц массы и жесткости производится формулировка матриц [СА] и [Д^], составление уравнения (1) для данной расчетной модели, нахождение с помощью функции eigs характеристических корней уравнения и визуализация результатов.
Составление матрицы [Х>А] производится следующим образом. Экспортированный из системы MSC.Nastran столбец значений коэффициента к,,
имеющий вид {к1 ,к2,... ,кп } , преобразуется в диагональную матрицу [/)А], имеющую вид
к1 0 0 к2
0 0
0 0
к
При составлении матрицы [СА] коэффициент умножается на значение скорости набегающего потока V. Полученные данные вводятся в матрицу как внедиагональные элементы. Место элемента соответствует степени свободы, по которой действует неконсервативная сила, и степени свободы, от которой она зависит. В общем виде матрица [СА] имеет вид
[Са] =
0 0 к у 0 ■ 0
0 0 0 ку ■ 0
0 0 0 0 ■ • к-У
0 0 0 0 ■ 0
0 0 0 0 ■ 0
Следует отметить, что, так как, во-первых, реальная расчетная модель имеет несколько десятков тысяч степеней свободы, а, во-вторых, значительная часть элементов столбца значе-
ний коэффициента к1 — нулевые, коэффициентами кп-2 V, кп-1 Vи к„¥при построении матрицы [СА] можно пренебречь.
Верификация расчетных
и экспериментальных аэроупругих характеристик конструкции
Для верификации расчетной методики должно быть проведено сравнение экспериментальных и расчетных динамических аэроупругих характеристик конструкции. К сожалению, в доступной литературе не удалось обнаружить достаточно подробного описания экспериментальных исследований флаттера конструкции при числах М> 3. В то же время разработанная методика применима и при скоростях 1,2 < М< 3. Поэтому в качестве эталонных были приняты аэроупругие характеристики пластины в сверхзвуковом потоке газа, полученные при проведении эксперимента [8].
Экспериментальные пластины изготавливались из стали 1Х18Н9 и дуралюмина Д16АТ размером 300 х 300 мм и 250 х 250 мм различной толщины. Приспособление для крепления образцов в аэродинамической трубе показано на рис. 3. Оно представляет собой плиту, которая двумя кромками крепится к стенкам трубы. Другие две кромки клинообразные (для обтекаемости). Плита имеет квадратную полость под пластиной. В днище полости сделаны дренажные отверстия для быстрого выравнивания давлений и уменьшения демпфирования воздуха в полости. Давление в полости практически равно давлению в потоке. Плита с закрепленной на ней экспериментальной пластиной устанавливается по уровню (горизонтально) в рабочей части аэродинамической трубы. Таким образом, пластина обдувается под нулевым углом атаки, с верхней стороны. С нижней стороны в полости приспособления находится неподвижный воздух. Для определения момента наступления флаттера применялись тензометры сопротивления. Тензометры наклеивались на нижнюю сторону пластины. Провода от тензометров выводились через тело плиты за стенку трубы. При продувке пластина доводилась до
Рис. 3. Пластина с экспериментальной оснасткой
момента разрушения конструкции под действием автоколебательных явлений.
По разработанной методике был проведен расчет динамических аэроупругих характеристик дюралюминиевых пластин 300 х 300 мм толщиной 1,1; 1,18; 1,3 мм. Граничные условия — жесткое защемление по контуру. Давление, плотность и температура невозмущенного потока соответствуют параметрам на уровне моря.
Рис. 4. Конечно-элементная модель пластины в системе MSC.Nastran
На рис. 5 показана эволюция собственных значений характеристического показателя для дюралюминиевой пластины размером 300 х 300 мм и толщиной 1,3 мм. Число М изменялось от МН = 2,5 до МК = 3,5 с шагом 0,05. Как видно из рисунка, флаттер наступает при 3 < М< 3,05. Увеличение шага и уменьшение диапазона рассматриваемых чисел М позволяет определить
Рис. 5. Эволюция собственных значений характеристического показателя 1-го, 2-го и 3-го тонов колебаний на комплексной плоскости в зависимости от числа М потока
Рис. 6. Зависимость критической скорости флаттера от параметров пластины
критическую скорость флаттера МКР = 3,02. Для пластины толщиной 1,18 мм расчет показывает МКР = 2,26, толщиной 1,1 мм — МКР = 1,7.
На рис. 6 даны экспериментальные значения и эмпирическая зависимость границы области устойчивости в работе [8], а также расчетные значения (звездочки), полученные по разработанной методике. Буквами обозначены: h — толщина пластины; a — длина кромки пластины; М — скорость Маха. Кривая построена для дюралюминиевых пластин и давлений, соответствующих уровню моря. Эксперименталь-
ные точки также приводились к этим условиям. Каждая экспериментальная точка соответствует пластинам такой толщины, при которых автоколебания еще возникают. У более толстых пластин флаттера не наблюдалось.
Как видно из приведенного сравнения, значения, полученные расчетным путем, хорошо согласуются с экспериментальными данными, что говорит о корректности разработанной расчетной методики для сверхзвуковых скоростей полета. В настоящее время ведутся работы по верификации расчетных аэроупругих харак-
теристик путем их сравнения с экспериментальными данными, полученными при продувках стендовых моделей ГЛА.
Выводы
1. Предложена методика расчета аэроупругих характеристик гиперзвуковых летательных аппаратов, позволяющая учесть особенности гиперзвукового потока.
2. Разработана программная реализация данной методики, основанная на максимальном использовании возможностей коммерческих расчетных систем.
3. Предложена модель аэродинамического нагружения в задаче гиперзвуковой аэроупругости. Как и в поршневой теории, аэродинамическое давление в точке принимается прямо пропорциональным местному эффективному углу атаки, что позволяет проводить исследование аэроупругих характеристик используя отработанные методики.
4. Проведена верификация расчетной методики путем сравнения с результатами экспериментальных исследований флаттера пласти-
ны. Получено хорошее совпадение расчетных и экспериментальных данных.
Список литературы
1. Marshall L.A., Catherine B.M., Corpening G.P., Sherrill R.P. Overview With Results and Lessons Learned of the X-43A Mach 10 Flight. Online at www.nasa.gov.
2. McNamara J.J., Friedman P.P. Aeroelastic and Aerothermoelastic Analysis of Hypersonic Vehicles: Current Status and Future Trends. AIAA Journal, 2007. № 2013. Р. 3.
3. Колесников К..С, Сухов В.Н. Упругий летательный аппарат как объект автоматического управления. М.: Машиностроение, 268 с. 1974.
4. Болотин В.В. Неконсервативные задачи теории упругой устойчивости. М.: Государственное издательство физико-математической литературы, 1961. С. 250—257.
5. Фершинг Г.Ф. Основы аэроупругости. М.: Машиностроение, 598 с. 1984.
6. Ватрухин Ю.И., Никитенко В.И., Кулик С.В. Методика определения аэродинамических нагрузок на летательный аппарат в динамических задачах аэроупругости. Международная конференция «Авиационные технологии 2000». Тезисы докладов. ЦАГИ. 1997. С. II—22.
7. Kurbatskii К., Montanari F. Application of pressure-based coupled solver to the problem of hypersonic missiles with aerospikes. AIAA. 2007. № 462.
8. Микишев Г.Н. Экспериментальное исследование автоколебаний квадратной пластины в потоке. Известия АН СССР. ОТН. Механика и машиностроение. 1959. № 1. С. 154—157.
Статья поступила в редакцию 24.01.2011 г.
