Особенности численного моделирования гиперзвукового обтекания простых тел

Бурное развитие во второй половине XX в. авиационной и ракетно-космической техники вызвало повышенный интерес к проблемам гиперзвукового обтекания тел. Для изучения процессов теплообмена при проектировании летательных аппаратов, получения стационар­ных аэродинамических характеристик в числе традиционных теоретических и эксперимен­тальных исследований широко применяется численное моделирование. Зачастую недости­жимые в лабораторных условиях возможность воспроизведения натурных условий и сравни­тельная дешевизна стимулировали быстрое создание эффективных численных алгоритмов получения данных о структуре течений около обтекаемых поверхностей.

Задача определения тепловых нагрузок на тела при входе их в плотные слои атмосферы с большими скоростями характеризуется ши­роким диапазоном изменения определяющих параметров по траектории движения; в част­ности, числа Рейнольдса набегающего потока изменяются от десяти до сотен миллионов, температура торможения достигает десятков тысяч градусов, а давление торможения - со­тен атмосфер. При этом высокие температуры инициируют протекание в ударном слое около обтекаемого тела сложных физико-химических процессов, т. к. внутренние степени свободы молекул возбуждаются и сохранение исходно­го состояния газа уже невозможно: молекулы кислорода и азота начинают диссоциировать на атомы. С дальнейшим ростом температуры начинается процесс ионизации с образованием свободных электронов. Образовавшиеся ато­мы, ионы и электроны диффундируют в более холодную область - к поверхности тела. Там происходит обратная реакция - рекомбинация, идущая с выделением тепла, что даёт допол­нительный вклад в нагрев вблизи поверхности тела. При этом изменяются физические свой­ства и состав воздуха, существенно влияющие на величины вязкости, теплопроводности и характеристики сжимаемости. Кроме того, из- за высоких температур значительно меняется удельная теплоёмкость компонентов воздуха. Стоит отметить, что диссоциация и ионизация могут поглощать до 75 % энергии потока, что делает неприменимыми многие результаты га­зовой динамики совершенного газа [1, 2].

Результаты расчётов конвективного те­плообмена в спускаемых аппаратах, входящих в атмосферу Земли с первой космической ско­ростью, в значительной степени зависят от вы­бора газодинамической и физико-химической модели высокотемпературного воздуха. Тепло­вой поток лучистого теплообмена принимает большие значения на малых высотах [3]. В сложной задаче одновременного учёта конвек­тивного и лучистого составляющих теплообме­на конечный результат зависит от правильного выбора математических моделей гиперзвукового обтекания тела и определения границ их ис­пользования. Предварительные исследования показали, что применимость этих моделей по высотам и скоростям зависит от размеров кос­мического спускаемого аппарата (КСА) и его траекторных параметров.

Цель работы состояла в определении воз­можностей и границ применимости числен­ных схем аппроксимации потоков Roe FDS и AUSM+, используемых для решения задач гиперзвуковой аэродинамики и исследований с их помощью течений в возмущённой области перед обтекаемым затупленным телом для установления закономерностей тепловых и газодинамических характеристик, связанных с разработкой необходимой тепловой защиты КСА.

При решении задачи гиперзвукового об­текания лучистый теплообмен не учитывался. Для определения конвективного теплового по­тока на поверхности КСА был применён конеч­нообъёмный метод решения полных уравнений Навье - Стокса со структурированными сетка­ми. В качестве рассматриваемых тел использо­вались тела цилиндрической формы.

Для создания геометрической и сеточной модели потока около цилиндра (рис. 1) исполь­зована свободно-распространяемая программа Gmsh [4].

В представленной сеточной модели раз­меры ячеек распределены следующим образом. Вдоль стенки, принадлежащей цилиндру, раз­меры ячеек распределены равномерно. Закон распределения размера ячеек перпендикулярно стенке (рис. 2) имеет следующее аналитиче­ское выражение:

m - минимальный размер элемента;

M - максимальный размер элемента; b - длина отрезка «крутого спуска» (на рис. 2 - участок кривой при l от 0 до 1 м);

K - отношение максимального размера эле­мента к размеру элемента в начале «крутого спуска»;

L - максимальное значение I.

Для оценки влияния параметров сеточ­ной модели на точность расчёта варьировался только параметр m, т. е. высота 1-й пристеноч­ной ячейки.

Основная проблема, возникающая при численном решении газодинамических уравнений - это устойчивость разностной схемы. Для того чтобы схема была устойчивой, она должна содержать разности, ориентированные против потока [5], т. е. на этапе построения схемы не­обходимо определить, в каком направлении каждая группа волн распространяется по рас­чётной сетке и как взаимодействуют расчётные ячейки. Как правило, используют две модели такого взаимодействия [6].

В первой модели ячейки взаимодейству­ют посредством дискретных волн, определяе­мых с помощью точного или приближённого решения задачи о распаде произвольного разрыва, заданного на границе между ячейками (задачи Римана). Дискретизация в этой модели реализуется на основе схем Годунова [7], Ошера [8], Роу [9].

Во второй модели взаимодействие меж­ду ячейками осуществляется через группы частиц, перемещающихся между ячейками и имеющих заданное распределение скоростей. Для разделения групп частиц на движущиеся «вперёд» и «назад» применяют методы рас­щепления потока. К разностным схемам рас­щепления относятся схемы Ван Лира [10], Лио и Стефана [11]. Группа схем Лио и Стефана получила название AUSM(Advection Upstream Splitting Method).

На данный момент разработано множе­ство модификаций схемAUSM. Первоначально метод Лио и Стефана использовался для расчё­та типичных аэродинамических задач, затем, с целью повышения точности, был усовершен­ствован, а позже обобщён на все скоростные режимы и многофазные течения. Также были предложены модификации этих обобщений.

В данной работе рассматривались два ме­тода расщепления потоков - Roe Flux-Difference Splitting Scheme (Roe FDS) и AUSM+. В резуль­тате анализа существующих методов [12], в том числе и перечисленных, выяснилось, что метод Roe FDS относится к схеме, склонной к неустойчивости - излому отошедшей ударной волны характерного для эффекта «карбункула» (carbunclephenomenon), в то время как AUSM+ является более устойчивой к этому явлению схемой (рис. 3). Характерный пример эффекта «карбункула», в качестве «шипа» на ударной волне, изображен на рис. 3а.

В настоящее время существуют две точ­ки зрения относительно эффекта «карбунку­ла» [13].

Первая из них заключается в следующем: проявление этого эффекта только в определён­ном диапазоне изменения значений определяю­щих параметров, главными из которых являют­ся высота H и скорость V полета (числа Маха M_∞ и Рейнольдса Re_D), прямо указывают на физичность его природы. Ещё более эта пози­ция аргументируется фактом одновременного развития других особенностей аэродинамики гиперзвукового обтекания - в частности, из­менением характеристик пограничного слоя и головной ударной волны, а также образова­нием сомкнутого вязкого ударного слоя. При этом не отрицается влияние вычислительного алгоритма на получаемое численное решение.

Вторая точка зрения на этот эффект про­тивоположна первой и заключается в отри­цании (практически безоговорочном) физики природы «карбункула», сводя его появление к свойствам вычислительного алгоритма. Глав­ным аргументом этой точки зрения является то обстоятельство, что эффект «карбункула» для одной и той же задачи, с полностью иден­тичным набором определяющих параметров, одними вычислительными алгоритмами об­наруживается, а другими не обнаруживается.

В качестве основного метода расщепле­ния потоков было решено использовать метод AUSM+.

Тестирование метода AUSM+ проводи­лось на примере распределения коэффициента давления C_p и удельного теплового потока q₁ по поверхности цилиндра с результатами экс­периментов [14] и эмпирических соотношений [3] при следующих параметрах гиперзвукового обтекания: радиус цилиндра R =0,0381 м; U_∞ = 2108,61 м/с (М_∞ = 15,622); T_∞ = 45,17 К; р_∞ = 23,622 Па; температура стенки T_w = 300,33 К; число Рейнольдса Re_D = 2,045·10⁵. Стоит отме­тить, что приведённые условия эксперимента позволяют определить точность численных схем аппроксимации без учёта турбулентности и физико-химических реакций.

Для дальнейших исследований использо­вались следующие параметры сеточной моде­ли: H=0,8D; M=0,02D; b=0,05D; K=30; m=10^-3D; 10^-4D; 10^-5D.

Расположение и форма ударной волны при поперечном гиперзвуковом обтекании ци­линдра, полученные в результате численного моделирования, практически не отличаются от данных по эмпирическим зависимостям (рис. 4), приведённым в [15]: погрешность расположения ударной волны в области точки торможения минимальна (менее 0,5 %). Этот результат позволяет говорить о том, что чис­ленная схема AUSM+ позволяет с достаточ­ной точностью предсказывать расположение и форму ударной волны при гиперзвуковых режимах течения.

Стоит отметить, что изменение высоты 1-й пристеночной ячейки не влияет на рас­пределение коэффициента давления C_p по по­верхности, а изменение величины q_l в окрест­ности точки торможения составляет порядка 5 % (рис. 5), поэтому в целях экономии вычис­лительных ресурсов достаточно использовать параметр m=10^-3D.

Полученные кривые распределения ко­эффициента давления C_p и удельного тепло­вого потока q_l по поверхности цилиндра (рис. 6, 7) большей своей частью подобны экспе­риментальным (погрешность 5-10 %) и эмпи­рическим (погрешность 3 %) зависимостям. Расхождение полученных данных между физическим и численным экспериментом объ­ясняется тем, что численный эксперимент не учитывал такие факторы, как степень турбу­лентности потока в трубе, неравномерность поля воздушного потока, влияние размеров рабочей части аэродинамической трубы, воз­мущения, вносимые датчиками, и т. д.

Вывод

После проведения численного моделирования и решения задачи теплообмена между газом и обтекаемым телом на примере гиперзвукового обтекания цилиндра выяснилось, что рассмо­тренный метод расщепления потоков AUSM+ способен решать задачи определения тепловых и газодинамических характеристик с прием­лемой точностью с определёнными, в рамках экономии вычислительных ресурсов, параме­трами сеточной модели.