УДК 533.6.011.55, 533.6.011.6, 532.526.3
В.Н. Бражко1,2, А. В. Ваганов1,2, В. Я. Нейланд1,2, М. А. Стародубцев1,2,
В. И. Шалаев1,2
1 Центральный аэрогидродинамический институт 2 Московский физико-технический институт (государственный университет)
Моделирование особенностей обтекания наветренной стороны треугольного крыла с затупленными передними кромками на основе численного решения уравнений Навье^Стокса
В работе, на основе численного решения методом конечного объема стационарных уравнений Навье-Стокса для сжимаемого вязкого и теплопроводного газа, исследовано гиперзвуковое обтекание конфигурации треугольного крыла с конусом. Рассмотрены ламинарный и турбулентный режимы течения, а также режим с ламинарно-турбулентным переходом. Особый акцент сделан на исследование структур пространственного течения около наветренной поверхности крыла. Результаты расчетов позволяют качественно объяснить наличие локальных зон пиковых тепловых потоков на крыле, наблюдаемых в экспериментах. В работе построено несколько сеточных моделей от 4 млн до 100 млн расчетных узлов. Непосредственное численное интегрирование выполнено в рамках пакета программ .\\SYS СГХ на вычислительном кластере ФАЛТ МФТИ.
Ключевые слова: гиперзвуковое течение, треугольное крыло с затупленными кромками, численное моделирование, тепловые потоки.
1. Введение
Особенности теплообмена на наветренной стороне треугольных крыльев являлись предметом многочисленных экспериментальных исследований разных авторов. Достаточно исчерпывающая библиография приведена в работах [1, 5]. Однако экспериментальные исследования, особенно при гиперзвуковых скоростях потока, сильно ограничены по средствам измерений, и имеющиеся данные не позволяют объяснить природу возникновения наблюдаемых особенностей. Современные средства численного моделирования дают возможность проанализировать тонкую структуру течения, связать ее с данными экспериментов и на этой основе построить физическую модель течения. В ряде случаев это позволяет расширить область параметров и осуществить переход от модельных условий к условиям натурного полета.
В настоящей работе численно исследуется интенсивность теплопередачи и пространственная структура течения вблизи наветренной поверхности треугольного крыла для чисел Маха набегающего потока М = 7.5, 8.3 и 10.5 при углах атаки а = 0 10° — режимы, исследованные экспериментально в аэродинамической трубе Т-117 ЦАГИ [1]. При численном интегрировании уравнений Навье-Стокса применяется метод конечного объема, реализованный в пакете программ СРХ. Рассматривается ламинарный, турбулентный
режим течения и режим, учитывающий ламинарно-турбулентный переход. Особое внимание уделяется структуре течения в области носика крыла в местах сопряжения сферического и цилиндрического затупления, а также вблизи плоскости симметрии на наветренной поверхности треугольного крыла.
Сеточная модель создавалась средствами пакета программ А^УБ 1СЕМ СГТ), где была выполнена блочная структура на 17 блоков и построен ряд сеток от 4 до 100 млн узлов. Сразу отметим, что основная часть работы выполнена на сетке в 50 млн узлов.
Эффекты турбулентности исследуются на модели к — шББТ, а при учете ламинарно-турбулептпого перехода — её модификации. Отметим, что в расчетной модели присутствует
демпфирование генерации турбулентности на скачках уплотнения, а в области присоединения отрывов, наоборот, вводится дополнительная генерация. К сожалению, существующие на сегодняшний день модели перехода, хотя и способны учитывать уровень внешней турбулентности и рост амплитуды волн Толлмина Шлихтинга, не позволяют оценить влияние неустойчивости поперечного течения (cross-flow instability) на переход. Последнее обстоятельство может иметь первостепенное значение при расчете тепловых нагрузок на наветренной стороне исследуемого крыла. Тем не менее в настоящей работе с целью оценки имеющегося подхода приведены некоторые результаты, полученные с помощью модели турбулентности к — wSST с ламинарно-турбулентным переходом.
Численное моделирование выполнено на 24-х узловом кластере ФАЛТ МФТИ.
2. Геометрическая модель
Проектирование компоновки затупленного крыла длиной L = 600 мм с конусным телом на верхней поверхности выполнено в НИО-8 ЦАГИ. Угол стреловидности крыла составляет % = 75°, а вершина и передняя кромка имеют радиус затупления 3 мм. Верхняя плоскость треугольного крыла отсекает конусное тело, образуя параболу. Внешний вид компоновки, включая подветренную и наветренную поверхности модели, представлен на рисунке 1.
Рис. 1. Треугольное; крыло с конусным телом: слева подветренная поверхность: справа наветренная поверхность
3. Сеточная модель
При построении сеточной модели использовался пакет программ АКБУБ 1СЕМ СТБ, куда импортировалась математическая модель ЮЕБ компоновки треугольного крыла с конусным телом. Форма и размер расчетной области выбраны таким образом, чтобы скачки уплотнения на всех рассматриваемых режимах не пересекали внешнюю границу. В работе делается акцент на изучение течения вблизи вершины и передней кромки крыла, поэтому область за телом не рассматривается. Полагается, что течение симметрично и для моделирования течения достаточно рассмотреть лишь одну половину крыла. При построении гекеаэдрической сетки применяется блочный подход, реализованный в модуле НЕХА. Созданная в данном модуле блочная структура состоит из 17 блоков и имеет плотное сгущение сеточных линий вблизи твердых стенок модели (рис. 2). Высота пристеночных ячеек фиксировалась и составляла 0.001 мм. Особое внимание уделено области около вершины крыла и местам стыка цилиндрической поверхности крыла со сферическим затуплением, где сосредоточена половина всех расчетных узлов. Всего выполнено 4 сеточных модели: на
4, 20, 50 и 100 млн узлов. Первые три из них относятся к половине, последняя к полной компоновке треугольного крыла.
Все приводимые в работе результаты относятся к расчетной сетке в 50 млн узлов. Применение более грубых сеток не позволило разрешить тонкие вихревые структуры на наветренной поверхности крыла вблизи плоскости симметрии, а сетка на 100 млн узлов создавалась зеркальным отображением сеточной модели на 50 млн узлов для оценки влияния условий симметрии.
Рис. 2. Расчетная область: слева блочная структура: справа фрагмент сетки
4. Физическая модель
Физическая постановка задачи численного моделирования гинерзвукового обтекания треугольного крыла для чисел Маха 7.5-10.5 (И,е = 0.5 — 1.7 ■ 106) обязана своим появлением результатам аэродинамического эксперимента ЦАГИ в трубе Т-117 [1]. В экспериментах на наветренной стороне крыла были обнаружены зоны аномальной теплопередачи (рис. 3), причины формирования которых пока не имеют объяснения. Целью настоящей работы является численное моделирование течений, реализованных в трубном эксперименте, и построение физической картины обтекания крыла, позволяющей объяснить повышение тепловых потоков в указанных зонах.
Рис. 3. Тепловые нагрузки на наветренной стороне треугольного крыла. Эксперимент в аэродинамической трубе Т-117 ЦАГИ. М = 7.5, Ие^ = 1.7 • 106, а = 10*
Для этого в расчетах применяется модель сжимаемого газа, в которой теплопроводность и вязкость зависят от температуры по закону Сазерленда. Поверхность треугольного крыла и конуса полагается изотермической, температура стенки Тт = 18 °С. Газодинамические параметры набегающего потока фиксировались на входной границе, а на выходе применялись «мягкие» условия. Рассматриваются следующие параметры набегающего потока:
Таблица 1
Режимы течений
№ а М» V», м/с Ро, Атм т», к ReL ■ 106
1 0 7.5 1094 11.42 53.0 0.8
2 10 7.5 1110 11.42 54.5 0.8
3 10 7.5 1123 25.50 55.8 1.7
4 10 8.2 1370 18.76 69.5 0.7
5 10 8.2 1280 32.65 59.7 1.2
6 0 10.5 1512 53.40 52.0 0.6
7 0 10.5 1543 94.67 53.8 1.0
8 10 10.5 1465 52.46 50.8 0.6
9 10 10.5 1455 97.85 47.3 1.0
Пояснения. В таблице используются следующие обозначения параметров набегающего потока: а — угол атаки, М» — число Маха, V» — скорость, Р0 _ полное давление, Т» — температvpa, Re^ — число Рейнольдса вычисленное по длине модели.
Исследуется ламинарный, турбулентный и переходный режимы течения. Для моделирования турбулентности используется модель к — wSST. Необходимо отметить, что она несколько отличается от классического варианта [2], а именно, имеет место демпфирование генерации турбулентности на скачках уплотнения, а в области присоединения потока, наоборот, ее увеличение.
Существующие на сегодняшний день модели ламинарно-турбулентного перехода способны учитывать уровень внешней турбулентности и рост амплитуды волн Толмина-Шлихтинга на положение перехода, однако неустойчивость поперечного течения (cross-flow instability) в них не рассматривается, и поэтому использование таких моделей в численном моделировании пространственных течений сильно ограничено. Тем не менее в настоящей работе приведены некоторые результаты, полученные с помощью одной такой модели, чтобы оценить возможный тип ламинарно-турбулентного перехода [3].
5. Описание численного подхода
Анализ рассматриваемых течений строится на решении полных уравнений Навье-Стокса (усредненных по Фавру в случае турбулентного течения) вязкого теплопроводного сжимаемого газа в стационарной постановке. При интегрировании уравнений применяется хорошо зарекомендовавший себя метод конечного объема, гарантирующий полное выполнения законов сохранения массы, импульса и энергии для любой составленной из контрольных объемов подобласти внутри расчетного домена. Указанное свойство консервативности расчетной схемы является необходимым условием для численного моделировании сильных скачков уплотнения. В работе используется сопряженный решатель программы ANSYS CFX, обеспечивающий второй порядок аппроксимации по пространственным переменным.
6. Некоторые результаты численного моделирования
Основным результатом проведенного численного моделирования является пространственная структура течения, возникающая на рассматриваемых режимах. Пример нескольких сечений поля чисел Маха в случае ламинарного обтекания крыла потоком газа М» = 7.5 под нулевым углом атаки представлен на рисунках 4 и 5.
На рисунках хорошо прослеживается положение ударной волны и разрешающая способность сеточной модели как около вершины (рис. 4), так и ниже по потоку в зоне интерференции скачков от крыла и конуса (рис. 5).
Как в численном, так и в трубном экспериментах на передней кромке треугольного крыла на расстоянии около 10г, где г — радиус затупления крыла, присутствует зона локаль-
ного повышения теплопередачи (рис. 6). Представленные на рис. 6 результаты демонстрируют, что рассчитанные размер и положение этой зоны качественно хорошо согласуются с данными измерений. Следует отметить, что распределение безразмерных тепловых потоков д/до (до — тепловой поток в критической точке) является автомодельным при использовании безразмерной продольной координаты х/К, где К — радиус носового затупления, для которого рассчитывался тепловой поток до-
Рис. 4. Поло чисел Маха в плоскости симметрии
Рис. 5. Поле числе Маха в сечении X = 100 мм и X = 500 мм
Что касается механизма возникновения повышенных тепловых потоков, то его можно объяснить, опираясь на следующие положения.
Ударные волны, формирующиеся перед сферически затупленной вершиной крыла и далее вниз по потоку перед передней его кромкой с цилиндрическим затуплением, отличаются по своей интенсивности, форме и удаленности от поверхности тела. Взаимодействие указанных скачков уплотнения приводит к появлению переходной зоны «пологого» фронта волны и сужению области течения в этой зоне. Это происходит вследствие того, что углы наклона ударной волны от вершины выше этой области меньше угла наклона волны, формирующейся в окрестности кромки ниже по течению. Газ, прошедший через пологий фронт волны, имеет пониженное значение энтропии и, следовательно, обладает большей кинетической энергией. Он прижимает высокоэнтроиийный поток от вершины и пограничный слой к поверхности, что приводит к сужению области сверхзвукового течения, повышению давления и уменьшению скорости течения вдоль кромки.
Из результатов расчетов ноля статической энтропии на рис. 7 видно, что зона максимальных тепловых потоков на поверхности крыла расположена прямо под областью торможения и, по-видимому, именно этот эффект является причиной интенсификации теплообмена. Сужение области течения около кромки является также причиной интенсификации поперечного течения во внешней части течения над поверхностью крыла вне кромки. Представленные результаты в целом совпадают с результатами, представленными в работах [4, 6].
Рис. 6. Тепловые нагрузки вдоль передней кромки крыла: слева численное моделирование: справа эксперимент Т-117
Рис. 7. Поле энтропии потока (верхняя часть) и тепловые потоки вдоль наветренной поверхности крыла (пижняя часть)
Другим важным следствием взаимодействия низкоэнтропийнш'о потока от кромки крыла и ударнш'о слоя, раенроетраняющм'оея вниз но потоку за ех'о вершиной, является формирование продольных вихревых структур на расстоянии около 40 мм от плоскости симметрии и в средней части крыла (рис. 8). Взаимодействие высокоэнтропийных медленных струек тока за сферической головной волной с более быстрым окружающим потоком приводит к образованию вихревой структуры вблизи плоскости симметрии один из вихрей этой структуры показан в левой части рисунка 8. Следует отметить, что подобный вихрь был обнаружен при решении уравнений Эйлера [7], поэтому природа ех'о образования, возможно, имеет невязкий характер. Аналогичная ситуация наблюдается и в зоне, где низкоэнтронийные быстрые струйки тока оказываются в окружении более медленных.
Как и в случае с высокоэнтранийными струйками, там, где проходит низкоэнтроний-ный газ, образуются новые вихревые структуры в средней части крыла (два вихря в правой части рисунка 8). Заметим, что такие структуры не были обнаружены при совместном ре-
шснии параболизованных уравнений Навье Стокса и уравнений Эйлера [6] и причины их образования, по-видимому, определяется более сложными эффектами, нежели централь-но!'о вихря.
Рис. 8. Поло энтропии потока (нижняя часть рисунка) и тепловые потоки вдоль наветренной поверхности крыла (верхняя часть рисунка)
Выделим теперь линии растекания нонеречноі'о вихревоі'о движения и нанесем их координаты на поверхности крыла. Полученные таким образом точки точно ложатся на полученные в трубном эксперименте [1] границы зон повышенных тепловых потоков (рис. 9). Отметим лишь, что линии растекания были получены только для ламинарноі'о режима течения. В турбулентном течении вихревое движение в зоне низкоэнтронийной струйки (средний вихрь) в численном эксперименте не наблюдаются. По-видимому, эффективная вязкость турбулентного потока не дает развиться слабым вихрям, быстро ноглощая низко-энтрониную трубку тока.
У, нм
Рис. 9. Тепловые нагрузки на наветренной стороне треугольного крыла. Эксперимент в аэродинамической трубе Т-117 ЦАГИ
Хорошо предсказывая положение зон повышенной теплопередачи, расчет ламинарного течения количественно сильно занижает пиковое значение теплопередачи (рис. 10). Возможной причиной такого несоответствия может оказаться переход течения от ламинарного режима к турбулентному, что не учитывается в расчетной модели.
Действительно, числа Рейнольдса но параметрам набегающего потока близки к переходным, а вихревые структуры в пограничном слое ускоряют процесс ламинарнотурбулентного перехода. Именно на линии растекания вихревые жгуты посредством поперечного перетекания способны ускорить ламинарно-турбулентный переход, инициировать
образование турбулентных клиньев и тем самым интенсифицировать повышенную теплопередачу.
Рис. 10. Тепловые нагрузки вдоль наветренной поверхности крыла: слева — сечение X = 300 мм; справа — сечение X = 500 мм
Косвенным подтверждением этой идеи является совпадение максимальных значений тепловых потоков при моделировании турбулентного течения с данными трубного эксперимента (рис. 10). Тогда, возникающие участки повышенной теплопередачи можно было бы трактовать как следствие образования вытянутых турбулентных зон, между которыми находятся ламинарные участки. К сожалению, существующие модели ламинарнотурбулентного перехода не учитывают влияние неустойчивости поперечного течения (crossflow instability) и предпринятая авторами попытка учесть переход также не увенчалась успехом. Полученная в расчете интенсивность теплопередачи представлена на рисунке 10. Несмотря на то, что экспериментальная кривая тепловых потоков имеет области, совпадающие с расчетами ламинарного и турбулентного течений, в последнем случае пространственное ноле течения не содержит вихревые структуры в средней части крыла, инициализированные низкоэнтронийной струей и, как следствие, вторая зона повышенной теплопередачи отсутствует.
24-узловой кластер ФАЛТ МФТИ
Производитель: Т-Платформы Операционная система: §ц£Е Ю-О Server SP1 Количество процессоров: 50 Количество ядер: 100
Суммарный объем оперативной памяти: 208 Qb Количество узлов: 24 Дисковая подсистема: 2920 QJ} Коммуникационная сеть: ШшБзШ Транспортная и сервисная сеть: Gigabit Ethernet Производительность: 844.9 GElQ.pp-Linp.acK 1021 СГ орр - Пиковая
Рис. 11. Используемый при численном моделировании кластер
Так, на рисунке 10, где представлено распределение тепловых потоков в сечении
X = 300 мм и X = 500 мм, максимум теплопередачи, полученный при расчете турбулентного течения, в обоих сечениях хорошо совпадает с данными трубного эксперимента, но отсутствует второй пик. С другой стороны, расчеты ламинарного режима правильно предсказывают положение центрального и среднего пиков, но величина их значительно меньше тепловых потоков, наблюдаемых в эксперименте.
Выполнение данной работы потребовало использования значительных расчетных ресурсов. Для этого все расчеты были выполнены на вычислительном кластере ФАЛТ МФТИ (рис. 11). Не вдаваясь в детали работы кластера, отметим, что вычисления на сетке в 100 млн узлов использовали 170 Гб оперативной памяти, а в расчетах были задействовано 96 расчетных ядер.
7. Заключение
В работе численно смоделирован ряд режимов гиперзвукового обтекания модели, состоящей из треугольного крыла с конусным телом на подветренной поверхности. Рассмотрены как ламинарный, так и турбулентный режимы течения. Продемонстрировано образование высоко- и низкоэнтропийных струй, наличие которых приводит к образованию вихревых структур за ударной волной. В случае ламинарного обтекания численные расчеты правильно предсказывают положение центрального и вторичного пика теплового потока в средней части крыла, но величина их существенно меньше данных трубного эксперимента. Расчеты турбулентного режима, наоборот, правильно предсказывают максимальную величину теплового потока в области центрального вихря, но не описывают развитие вторичной вихревой структуры в средней части крыла и связанный с ним второй пик теплового потока.
Учет ламинарно-турбулентного перехода позволил бы провести более корректное сравнение расчетных и экспериментальных данных, но имеющийся в программе расчета механизм не включает в себя влияние неустойчивости поперечного течения (cross-flow instability) и потому не дает полного описания пространственного течения.
Представленные результаты показывают, что рассмотренные в работе явления имеют сложную природу и требуют дальнейшего изучения на базе более совершенных подходов.
Работа выполнена при финансовой поддержке РФФИ (проекты № 09-01-00206а и 13-
01-00573а) и программы совместных фундаментальных исследований ЦАГИ и ИТПМ СО
РАН.
Литература
1. Бражко В.Н., Ваганов А.В., Дудин Г.Н., Ковалева Н.А., Липатов И.И., Скуратов А. С. Экспериментальное исследование особенностей аэродинамического нагревания треугольного крыла при больших числах Маха// Труды МФТИ. — 2009. — Т. 1, № 3. — С. 58-67.
2. Menter F.R. Two-Equation Eddy-Viscosity Turbulence Models for Engineering Applications 11 AIAA Journal. - 1994. - V. 32, N 8. - P. 1598-1605.
3. Langtry R.B., Menter F.R. Transition Modeling for General CFD Applications in Aeronautics 11 AIAA Paper. - 2005. - N 522.
4. Губанова О.И., Землянский Б.А., Лесин А.Б., Лунев В.В. и др. Аномальный теплообмен на наветренной стороне треугольного крыла с затупленным носком при геперзвуковом обтекании// Аэродинамика воздушно-космических систем. Ч. 1, Б.М. — 1992. — С. 188— 196.
5. Лесин А.Б., Лунев В.В. О пиковых тепловых потоках на треугольной пластине с притупленным носком в гиперзвуковом потоке // МЖГ. — 1994. — № 2. — С. 131-137.
6. Власов В.И., Горшков А.Б., Ковалев Р.В., Лунев В.В. Тонкая треугольная пластина с притупленным носком в вязком гиперзвуковом потоке // МЖГ. — 2009. — № 4. —
С.134-145.
7. Ваганов А.В., Нейланд В.Я., Чернов С.В., Юмашев В.Л. Численное исследование особенностей невязкого течения около наветренной стороны треугольного крыла с затупленными передними кромками при гиперзвуковых скоростях потока // Материалы десятой Международной школы-семинара «Модели и методы аэродинамики». - М.: МЦНМО. - 2010. - С. 29-30.
Поступим в редакцию 31.12.2010.