CC BY
56
Том XL
УЧЕНЫЕ ЗАПИСКИ ЦАГИ
2009
№ 5
УДК 532.526.5.011.7
О ВТОРИЧНОМ ОТРЫВЕ ПРИ СВЕРХЗВУКОВОМ ОБТЕКАНИИ
УГЛА СЖАТИЯ
В. В. ШВЕДЧЕНКО
Проведено численное исследование отрывного течения при сверхзвуковом (M, = 5) обтекании угла сжатия. Детально рассмотрено явление зарождения и развития вторичного отрыва. Классифицированы стадии отрыва по параметру подобия.
Ключевые слова: угол сжатия, отрыв, triple deck, сверхзвуковой поток.
Исследованию проблем отрывного течения в сверхзвуковом потоке вязкого газа около плоской пластины, когда отрыв вызван отклонением на угол 0 задней части пластины, посвящено большое количество экспериментальных и теоретических работ. Обзор результатов можно найти, например, в [1 — 4].
Для развитых отрывных течений, в которых существует область с почти постоянным давлением (область «плато»), обычно применяется подход с использованием критерия Чепмена — Корста [5 — 7]. Для зарождающихся зон отрыва (или малых зон) используется подход, связанный с использованием интегральных уравнений пограничного слоя. После создания асимптотической теории «свободного взаимодействия» [8 — 11] (в литературе часто используется термин triple deck) в рамках этой теории были получены решения с многослойной структурой, например [7, 12]. Расчеты обтекания «угла сжатия» с углом поворота 0 ~Re-14 в рамках асимптотической теории проводились во многих работах [11 — 14]. Наиболее полное изложение результатов теоретических исследований для малых зон отрыва содержится в монографиях [15, 16].
Современные вычислительные методы решения уравнений Навье — Стокса позволяют провести численное моделирование при параметрах течения, когда отрыв не описывается асимптотической теорией свободного взаимодействия. В работе [17] изучено влияние температурного фактора Tw (отношения температуры тела к температуре торможения набегающего потока) на структуру отрывного течения, которое возникает при обтекании сверхзвуковым потоком вогнутого угла. Выявлено сильное влияние температурного фактора на длину зоны отрыва и на создаваемые потоком аэродинамические характеристики. При достаточно больших углах (при фиксированном значении числа Рейнольдса Re), когда течение не может описываться теорией свободного взаимодействия, наблюдались вторичный отрыв (отрыв возвратного течения внутри первичного отрыва), значительные поперечные градиенты давления и вихревые структуры. В этой работе также было показано значительное влияние разрешения сетки на качество получаемых результатов, особенно при больших размерах отрывной зоны. Данная работа является продолжением подробного численного исследования внутренней структуры отрывной зоны, начатого в работе [17].
1. Методика. Численное исследование проводилось на основе пакета программ решения уравнений Навье — Стокса методом установления по времени, разработанного в ЦАГИ [18 — 20].
Уравнения Навье — Стокса в случае двумерной задачи (плоскопараллельное течение), решаемой в произвольной криволинейной системе координат п, где х = х(^, п), У = y(^, п) — декартовы координаты, записываются в дивергентной форме:
д О д Е д С п
—- +-+-= 0 .
дг д^ дп
Здесь О — вектор консервативных зависимых переменных задачи, Е и С — векторы потоков в криволинейной системе координат. Векторы О, Е, С связаны с соответствующими векторами О,,, Ес, Сс в декартовой системе координат формулами:
о=Тос, Е=л Ее £+с , £ 1, С=I (Ее да+Сс д*
дх ду) I дх ду
где 1= д( X, у) — якобиан преобразования. д(^, П)
Декартовы компоненты векторов Ос, Ес, Сс для двумерных уравнений Навье — Стокса имеют вид:
Ос =
р ри Pv
ри , Ее = ри 2 + р — Т хх , С с = р^ — т ху 2
Pv р^ — т ху PV + р — Туу
е риН — ит хх — vт ху + Чх р^ — ит ху — ^ уу + Чу
где р — плотность; и, V — декартовы компоненты вектора скорости V; р — давление; е = р(^Т + (и2 + V 2)/2) — полная энергия на единицу объема; Н = СрТ + (и2 + V 2)/2 — полная энтальпия; Ср и с^, — удельные теплоемкости при постоянном давлении и объеме; X — коэффициент теплопроводности; ц — коэффициент динамической вязкости; т — тензор вязких напряжений с компонентами
, ди 2
,---(
дх 3
(ди дv
т XX = Ц|2 —~ ^ I, Т ху = т ух = Ц ^ + дХ)' Туу = Ц [ 2 ду " ;
, дv 2
дТ
дТ
а — вектор теплового потока с компонентами ах = —X—, ау = —X . 4 Р Чх дх Чу ду
Система уравнений замыкается уравнением состояния для совершенного газа:
р = рЯТ / т .
Здесь Я — универсальная газовая постоянная, т — молярный вес газа. Коэффициенты переноса определяются следующим образом: динамический коэффициент вязкости в зависимости от температуры изменяется по степенному закону ц/цте = (Т/ Тте)и, а число Прандтля принимается постоянным Рг = ЦСр/Х.
Безразмерное напряжение трения после введения безразмерных переменных (используемых _ _ _ _ _ _ _ 2 _ 2
далее) х = х/Ь , у = у /Ь , и = и /, V = V/, ^ = /Ь, р = р / , р = р / , = 1/уМте,
Т = Т/, ц = ц/ц^ =(Т/Тте)И=Т , где Ь — расстояние от начала пластины до угловой точки
записывается в виде с^ = тху =
"ху
ц (ди дv ) Ь
- =—I — + — I, а число Яе =-2 Яе I ду дх
р~и~ Яе ^ ду дх)
Начально-краевая задача решалась численно методом установления по времени на основе интегроинтерполяционного метода (метода конечного объема). При аппроксимации конвективной составляющей векторов потоков в полуцелых узлах использована неявная монотонная схема
Рис. 1. Поле плотности р (Яе = 106, Ту = 1) для угла 6 = 5° (а), 10° (б), 20° (в, г):
области: 1 — невозмущенного течения; 2 — скачка от передней кромки; 3 — скачка от точки отрыва; 4 — зоны отрыва; 5 — скачка от точки присоединения; 6 — область 5 после взаимодействия с областью 3; 7 — пограничный слой; г — внутренняя область отрыва при нестационарном вторичном отрыве (Ту = 0.1,
6-я стадия)
типа Годунова [21] и приближенный метод Роу [22] решения задачи Римана о распаде произвольного разрыва. Для повышения порядка аппроксимации до второго при интерполяции зависимых переменных на грань элементарной ячейки использован принцип минимальных производных [23]. При аппроксимации диффузионной составляющей векторов потоков на грани элементарной ячейки применена разностная схема типа центральных разностей второго порядка точности. Для решения нелинейных сеточных уравнений использовался модифицированный метод Ньютона — Рафсона. Решение системы линейных алгебраических уравнений осуществлялось при помощи метода минимальных невязок ОМЯЕ8(к) [24].
Решалась задача двумерного обтекания угла отклоненного навстречу сверхзвуковому пото-
4 7
ку газа = 5) в широком диапазоне значений числа Рейнольдса (Яе = 10 — 10) и других определяющих параметров. В расчетной области (рис. 1) начало координат (х = 0) совпадает с началом неотклоненной части пластины, точка отклонения угла расположена при х = 1, правая граница расчетной области х = 5 (х, у — декартовы координаты). Максимальный размер расчетной области в поперечном направлении был ~1 (при х = 5). На левой границе ставились граничные условия невозмущенного набегающего потока. Верхняя граница расчетной области выбиралась так, чтобы граничные условия также являлись условиями набегающего потока. Правая граница бралась на достаточно большом расстоянии, чтобы ошибка в мягких граничных условиях на этой границе (Эи/Э£, = 0, Эу/Э£, = 0, Эр/Э£, = 0, ЭТ/Э2, = 0, где — нормаль к правой границе) практически не влияла на решение задачи в окрестности зоны отрыва. На твердой поверхности выполнялись условия прилипания для скорости и = 0, V = 0, а для давления и температуры граничные условия Эр/Эп = 0, Т = Ту/, где п — нормаль к поверхности, а Ту — значение температурного фактора.
В начале пластины производилось специальное сгущение сетки, чтобы корректно отследить резкий градиент давления в носике. Заметим, что погрешности в начале пластины (если начало зоны отрыва не находится в этой области) не имеют большого влияния на решение ниже по течению и быстро затухают при удалении от носика пластины.
Расчеты проводились как на ортогональных сетках с программным сгущением вблизи поверхности (аналитических) [20], так и на подстраивающихся под особенности решения сетках (адаптивных) [25], построенных методом одномерного эквираспределения [26]. Был проведен тщательный анализ получаемых решений на сетках с различным количеством узлов (размерности). Решения, полученные на аналитических сетках малой размерности (низкого разрешения) 200 X 200, 400 X 200, 800 X 200, значительно между собой различались и заметно занижали размеры зоны отрыва (вплоть до его отсутствия на очень грубых сетках). На адаптивных сетках [17] длина зоны отрыва и положение слоя смешения практически не менялись при увеличении разрешения сетки, а уточнялась только внутренняя структура области отрыва. Решения, полученные на аналитических и адаптивных сетках большой размерности (1600 X 200, 1600 X 400), между собой практически не различались, но использование адаптивных сеток позволяло получать относительно хорошие результаты на сетках малой размерности, что существенно ускоряло процесс получения (уточнения) окончательного решения. Применение методики получения «относительно хорошего решения» на сетках малой размерности особенно актуально в трехмерном случае, когда практически нет возможности применять сетки большой размерности.
При решении сеточным методом уравнений Навье — Стокса, кроме вязкой (навье-стоксов-ской) диссипации, всегда присутствует дополнительная сеточная диссипация (вязкость), которую в расчетах всегда стремятся минимизировать. Численные исследования показали наличие заметной дополнительной сеточной диссипации, связанной с грубостью сетки, особенно при больших значениях числа Яе.
При малых значениях числа Яе, когда вязкая диссипация велика, влияние сеточной диссипации практически незначительно. По мере уменьшения навье-стоксовской диссипации при увеличении числа Яе вклад в решение сеточной диссипации на грубых сетках становится значительным. Единственным способом ее снижения является увеличение общего количества узлов или сгущение сетки (повышение разрешения) в принципиально важных областях течения (пограничный слой, слой смешения, область отрыва и присоединения, область вторичного отрыва).
В расчетах использовались квазиравномерные сетки, т. е. сетки с медленно меняющейся функцией распределения узлов (Ai/Ai + 1 = 1 + а, где |а| << 1). Построение квазиравномерных сеток на заданной расчетной области осуществлялось следующим образом. Сначала каким-либо образом строилась ортогональная квазиравномерная сетка. Путем программных сгущений [20] строилась использовавшаяся аналитическая сетка. При построении адаптивной сетки задавался коэффициент к, регулирующий степень адаптации (при к = 0 получалась первоначальная квазиравномерная сетка). К адаптивной сетке применялись такие же программные сгущения, как и для аналитической. Полученная адаптивная сетка отличалась от аналитической только дополнительным адаптивным сгущением в области слоя смешения и некоторых других областях.
В данной задаче применение адаптивных сеток изначально направлено на разрешение слоя смешения, положение которого практически не меняется на адаптивных сетках различной размерности. При этом в ходе решения появилась необходимость изучения внутренней структуры области отрыва, где возникают определенные затруднения и для адаптивных сеток. Внутри зоны отрыва возникали области достаточно интенсивных возвратных течений (М ~ 1). После применения процедуры адаптации сетка автоматически сгущалась и в этих областях. При этом новое решение (возвратная струя) начинало «выпадать» из этих сгущений. На построенной адаптивной сетке получалось решение, отличное от предыдущего, и все различия были сосредоточены во внутренней части первичного отрыва (во вторичном отрыве). Процедура динамического построения адаптивной сетки (построение новой адаптивной сетки после некоторого временного шага) позволила преодолеть эту неоднозначность. Таким образом, после малого изменения в поле строится новая сетка, мало отличающаяся от предыдущей. Использование постоянно меняющихся (динамических) адаптивных сеток позволило «держать» сгущение в этих областях, повысить разрешение сетки в области интенсивных течений в зоне отрыва и заметно снизить сеточную диссипацию возвратной струи. Это значительно повысило качество как стационарного, так и нестационарного решений.
Повышенное разрешение сетки применялось в области слоя смешения, точки присоединения и области вторичного отрыва. В области присоединения достаточно было программного сгущения в продольном и адаптивного сгущения в поперечном направлении. В области вторичного отрыва повышение разрешения достигалось в основном путем увеличения количества узлов сетки в продольном направлении и перераспределения их в эту область. В области слоя смешения использовалось как программное, так и адаптивное сгущение в поперечном направлении (для аналитических и адаптивных сеток соответственно).
В целом, в области отрыва расположено примерно 75% от общего количества узлов сетки как в продольном, так и в поперечном направлении. Из них на слой смешения приходится около 20% (для адаптивной сетки), а на программное сгущение в пристеночном слое порядка 10% от общего количества узлов. Так, например, для области отрыва самого большого размера (0 = 20°, Ту/ = 1) продольный размер области отрыва составляет ~2, а поперечный ~0.17, что для сетки 1600 X 200 соответствует примерно 1200 узлам в продольном направлении и 150 узлам в поперечном. Средний продольный и поперечный размер ячеек сетки внутри отрыва (с учетом неравномерности) составляет А ~ 1 — 2 • 10 . Для области отрыва меньших размеров соответственно размеры ячеек меньше. Поперечный размер ячеек в пристеночном слое составлял А~10 4 — 10 5.
Описанный способ адаптации приводил к искажениям первоначально ортогональных ячеек вдали от поверхности. Специально проведенное тестирование на модельных задачах с различными искажениями типа «параллелограмм» на сетках с разным количеством узлов показало, что это не ухудшает полученное решение, а в отдельных случаях даже улучшает. Об этом также свидетельствует сходимость решений на различных сетках в отсутствие вторичного отрыва.
Отметим, что результаты с использованием граничных условий для твердой поверхности на градиент давления др/дп =0 и др/дп = ± дт^/д^, где т^ — поверхностное напряжение трения; П — нормаль к поверхности; В — касательная к поверхности, практически между собой не различались.
Для нестационарных вариантов была проверена независимость получаемого решения от величины расчетного шага по времени А^ = А^/к, где к = 1, 2, 3, 4, 6.
При получении решения на более подробной сетке в качестве начального приближения использовалось решение, полученное на более грубой сетке. Разницу между начальным приближением и решением, полученным методом установления, можно рассматривать как отклонение от решения (возмущение). Процедуру получения решения на грубых сетках можно рассматривать как способ подавления этих возмущений с использованием дополнительной сеточной вязкости. Переход к более подробным сеткам означает подавление более слабых возмущений навье-сток-совской диссипацией и меньшей сеточной вязкостью при условии, что окончательное решение является стационарным. Следовательно, если на грубых сетках решение стационарно, а на более подробных сетках нестационарно, то «косвенным» источником нестационарности являются эти возмущения. Настоящей причиной данной нестационарности, возможно, является некоторый незатухающий (или слабозатухающий) автоколебательный процесс. При этом не исключается причина, связанная с особенностями расчетного метода.
Для всех нестационарных вариантов проводилась достаточно длительная процедура счета с целью получения стационарного решения. В нестационарном случае в области отрыва наблюдалось квазипериодическое течение при практически неизменном внешнем течении. Небольшие возмущения от взаимодействия вихрей и слоя смешения распространялись вниз по потоку. Детальный анализ обнаруженной нестационарности выходит за рамки данной статьи и является темой отдельного дополнительного исследования.
Были проведены тестовые расчеты (на нестационарных вариантах) для проверки отсутствия влияния продольных размеров расчетной области на полученное решение. Размеры расчетной сетки изменялись путем удаления некоторого количества последних сеточных линий. Таким образом, в неудаленной части поля сетки полностью совпадали. Для всех тестовых вариантов был проведен расчет с некоторым одинаковым интервалом времени, когда полученное решение заметно отличалось от первоначального из-за нестационарности. Все решения между собой полностью совпадали (за исключением небольшой области вблизи выходной границы) при условии, что граница расчетной области располагалась правее резкого градиента давления в области присоединения.
2. Результаты. Основные расчеты проведены в широком диапазоне изменения параметров газа: M = 5, у = 5/3, Re = 104 — 107, Tw = 10-3 — 1, ю = 0.5 (| ~ Iю), Pr = 2/3, 0 = 0 — 20° на сетках c последовательным удвоением числа узлов по продольной координате: 200X 200, 400X 200, 800 X 200, 1600 X 200. Расчеты проводились с использованием как аналитических, так и адаптивных сеток. Для методических целей также выборочно проводились расчеты на сетках 800 X 400, 1600 X 400.
На рис. 2 приведены зависимость х-координаты точек отрыва и присоединения от числа Re. При возрастании угла 0, температурного фактора Tw и числа Re размеры области отрыва увеличиваются. Положение точки отрыва при увеличении числа Re монотонно смещается к началу пластины, и для больших углов 0 = 10 — 20° отрыв начинается с носика пластины. Отметим наличие максимума на зависимостях координаты точки присоединения от числа Re, при котором происходит зарождение вторичного отрыва. Значение х-координаты максимума уменьшается при уменьшении угла 0, а максимум достигается при больших значениях числа Re.
Рис. 2. Зависимость х-координаты точек отрыва (нижняя ветвь) и присоединения (верхняя ветвь) от числа Яе для углов 8 = 2.5, 5, 7.5, 10, 15, 20° (кривые 1, 2, 3, 4, 5, 6 соответственно) при:
Tw = 0.1; б — Tw
1
а
При всех различиях в размерах и формах отрыва для разных углов 0 и чисел Яе результаты удобно разбить по характерным признакам на стадии, каждую из которых можно разделить на начальную и развитую:
а) безотрывное течение;
б) отрыв без вторичного отрыва: 1-я стадия — возникновение и 2-я стадия — развитие отрывной зоны;
в) стационарный вторичный отрыв: 3-я стадия — появление и 4-я стадия — развитие вторичного отрыва;
г) нестационарный вторичный отрыв: 5-я стадия — появление и 6-я стадия — развитие нестационарности.
Эти шесть стадий хорошо коррелируют по параметру подобия = 0 Яе1/4 [8, 9], что отражено в табл. 1, 2. То есть для различных углов 0 значение параметра однозначно определяет момент (число Яе) возникновения отрыва (^о ~ 2), зарождения вторичного отрыва (^о ~ 4) и появления нестационарности внутри отрыва (^о ~ 7).
Таблица 1
Классификация стадий отрыва
№ ^0 Характеристика Особенности
0 <1.5 Безотрывное течение Гладкое распределение давления
1 ~2 Отрыв с гладким трением Практически отсутствует полка давления
2 ~3 Отрыв с искажением по трению Начало формирования полки давления
3 4 — 5 Зарождение вторичного отрыва Сходимость по сеткам, смещение полки давления от угловой точки, появление провала в распределении давления по поверхности
4 5 — 6 Стационарный вторичный отрыв Формирование первых вихрей, появление различий для аналитических сеток с разным количеством узлов, применение адаптивных сеток
5 7 — 8 Зарождение нестационарности вторичного отрыва Рост количества вихрей, значительные различия решений для аналитических сеток с разным количеством узлов, грубые сетки — стационарные, подробные — нестационарные, динамические адаптивные сетки
6 ~10 Нестационарный отрыв На всех сетках нестационарность, усложнение структуры вихрей
Таблица 2
Значения номера стадии отрыва и параметра (в скобках) при различных значениях угла 0 и числа Ке для температурного фактора: слева Т„ = 0.1, справа Т„ = 1
Яе / 0° 2.5 3.75 5 7.5 10 15 20
1 • 103 0 (1.47) 0 0 (1.96) 1
3 • 103 — (1.29) — 0 (1.94) 1 1 (2.58) 1
1 • 104 — (131) — 0 (1.75) 1 1 (2.62) 1 2 (3.49) 2
3 • 104 — (1.15) — 0 (1.72) 1 1 (2.30) 1 2 (3.45) 2 3 (4.59)3
1 • 105 — (1.16) — 0 (1.55) 0 1 (2.33) 1 2 (3.10) 2 3 (4.66) 3 4 (6.21) 4*
3 • 105 — (1.02) — — (1.53) — 0 (2.04) 1 2 (3.06) 2 3 (4.08) 3 4 (6.13) 4* 5 (8.17) 5*
1 • 106 0 (1.38) 0 — (2.07) — 1 (2.76) 1 3 (4.14) 3 4 (5.51) 4* 5 (8.28) 5* 6 (11.0) 6
3 • 106 0 (1.82) 1 — (2.72) — 2 (3.63) 2 4 (5.45) 4* 5 (7.26) 5* —(10.9)— — (14.5) —
1 • 107 1 (2.45) 1 — (3.68) — 3 (4.91) 3 5 (7.36) 5* 6 (9.81) 6 —(14.7)— — (19.6) —
3 • 107 2 (3.23) 2 — (4.84) —
1 • 108 3 (4.26) 3 — (6.55) —
« — » расчеты не проводились.
Для первых стадий (1 — 3) решения мало зависят от размерности сеток. Чем ниже стадия, тем меньше эта зависимость, и для начальных стадий возможно использование сеток с малым количеством узлов. Для более высоких стадий (4 — 6) использовались сетки с большим количеством узлов (800 X 200, 1600 X 200), а узлы перераспределялись так, чтобы максимально повысить разрешение сетки в области вторичного отрыва (в несколько раз выше, чем на равномерной сетке). Размерность сетки выбиралась такой, чтобы решение не зависело от количества узлов по продольному (Nx) и поперечному (Ny) направлениям. Если при больших (4 — 6) стадиях отрыва для адаптивных сеток достаточно Ny = 200, то для аналитических необходимо Ny = 400 и выше. Необходимое количество узлов по продольной координате (Nx = 200 — 1600) в большей степени зависит от номера стадии отрыва, чем от типа сетки.
В табл. 2 приведены номера этих стадий для различных значений угла 0, числа Re и Tw. Отметим, что эволюцию можно просматривать как по числам Re, так и по углам 0, а также наличие
априорной информации — по параметру даже без расчетов. Заметим, что для угла 0 = 2.5° за-
7 8
рождение первичного (Re ~ 10 ) и вторичного (Re ~ 10 ) отрыва происходит при очень больших
значениях числа Рейнольдса.
На рис. 3 приведены распределения давления и напряжения трения по поверхности угла
4 7
0 = 10°, для которого (табл. 2) в диапазоне чисел Re = 10 ■ 10 присутствуют все перечисленные стадии отрыва. Отметим значительную «размазанность» полки в распределении давления вплоть до стадии образования вторичного отрыва, смещение полки давления от угловой точки и характерный провал на кривой распределения давления при образовании вторичного отрыва. При увеличении числа Re «размазанность» переходной зоны распределения давления в области присоединения уменьшается и трансформируется в конкретное значение координаты области присоединения для каждого температурного фактора (это отмечено в [17]). Графики распределения поверхностного напряжения трения хорошо дополняют картину эволюции вторичного отрыва, при котором трение повторно обращается в ноль. Для 5-й и 6-й стадий отрыва наблюдаются характерные осцилляции, связанные с вихревыми структурами.
На рис. 4 для углов 0 = 2.5, 5, 7.5, 10, 15, 20° приведены распределения давления и поверхностного напряжения трения для чисел Re, при которых значение параметра ^0 = 4 — 5 соответствует 3-й стадии развития отрыва (зарождение вторичного отрыва). Для малых углов 0 образование вторичного отрыва происходит при малых размерах области отрыва.
Для больших номеров стадий (4 — 6) характерно формирование вихрей внутри зоны отрыва (аналогично «bubble» в работе [27]), что сопровождается осцилляциями в распределении давления и напряжения трения по поверхности (рис. 3). На рис. 5 приведена эволюция течения внутри области отрыва при различных числах Re и Tw для угла 0 = 20°, для которого отрыв начинается практически с носика (угол 20° выбран из-за большого поперечного размера отрыва и достаточно малых чисел Re). Подобные картины течения реализуются и для меньших углов 0 при больших числах Re и меньших размерах области отрыва.
При зарождении (3-я стадия) вторичный отрыв носит слоистый характер (рис. 5). Для малых температурных факторов (Tw ~ 10-3 — 10-1) зарождающийся вторичный отрыв расположен заметно левее угловой точки (x ~ 0.95, 0 = 20°), а для больших температурных факторов (Tw > 0.1) в угловой точке. Для 4-й стадии (Tw < 0.3) характерно появление первых стационарных вихревых структур. С увеличением числа Re (5-я стадия) количество вихрей возрастает, и они становятся нестационарными. Следует отметить, что если сетки недостаточно подробные, то картина течения существенно зависит от сетки, и могут получаться разные решения: стационарные и нестационарные. Эта неоднозначность разрешается применением динамических адаптивных сеток, для которых все результаты становятся нестационарными. Дополнительные тестовые расчеты с различными временными шагами приводили к одинаковым результатам с нестационарной картиной течения. Такая чувствительность решения по сеткам хорошо объясняется концепцией дополнительной сеточной диссипации. Когда навье-стоксовская диссипация достаточно велика (малые стадии), то малая дополнительная сеточная диссипация незначительно влияет на картину течения. Однако когда навье-стоксовская диссипация становится малой, то становится важно, насколько мала сеточная диссипация. И если она недостаточно мала (грубая сетка), то это может стабилизировать решение и сделать его стационарным.
0.02 -1--
0.1 —.-->1
Рис. 3. Распределение давления р (слева) и поверхностного напряжения трения суЯе1/2 (справа) по х-координате поверхности угла 0 = 10° при числах Яе, соответствующих 1 — 6-й стадиям отрыва (сверху вниз), для температурного фактора Tw = 10-3, 10-1, 1 (соответственно кривые 1, 2, 3)
Рис. 4. Распределение давления р (слева) и поверхностного напряжения трения суКе1/2 (справа) по х-координате поверхности угла 8 = 2.5, 5, 7.5, 10, 15, 20° при числах Яе, соответствующих 3-й стадии отрыва (зарождение вторичного отрыва, ^0 = 4 — 5), для температурного фактора Т=10-3,10-1,1 (соответственно кривые 1, 2, 3)
Рис. 5. Поле чисел Маха (в диапазоне значений от нуля до 0.5) внутри области отрыва угла 0 = 20° для 3 — 6-й стадий (числа Яе: 3-104,1 • 105, 3-105, 1 • 106 соответственно) при температурном факторе Т„ = 0.1, 0.3, 1
- \ 6
- ^-т--1
!_
1
0 0.2 Ч 0.4
Рис. 6. Профили температуры Т по нормали п в угловой точке (0 = 20°) при числе Яе = 3 • 104 (3-я стадия — зарождение вторичного отрыва) для различных значений температурного фактора Т„ = 10-3, 10-2, 0.1, 0.2, 0.3, 1 ( кривые 1, 2, 3, 4, 5, 6 соответственно)
Для 6-й стадии характерны значительное усложнение вихревых структур и их нестационарность (и для ^ = 1) практически для всех не очень грубых сеток, как адаптивных, так и аналитических. При этом необходимо использовать сетки с максимальным разрешением как в продольном, так и в поперечном направлении (1600 X 200, 1600 X 400). При уменьшении температурного фактора интенсивность осцилляций распределения давления и напряжения трения в области отрыва заметно увеличивается (рис. 3, 5).
Если при малом температурном факторе (^ < 0.3) для 4-й и 5-й стадий характерно образование вихрей и появление нестационарности вторичного отрыва, то при большом температурном факторе (^ = 1) сохраняется стационарность протяженного вторичного отрыва с небольшими
осцилляциями возвратной струи в области присоединения (отмечено «*» в табл. 2). Расчеты на более подробной сетке 1600 х 400 не изменили картину течения. Было предположение, что это связано с большим значением вязкости при Tw = 1 (ц ~ T™). Специально проведенные расчеты (табл. 3) для ю = 0 (ц = const) и ю = 1 (Pr = 2/3) показали заметное влияние закона вязкости на процессы первичного, вторичного отрывов и вихреобразования (смещение стадий по числу Re), но качественно не изменили ситуацию. Вероятно, это связано с тем, что для образующегося вторичного отрыва внешним является течение внутри первичного отрыва, которое существенно зависит от температурного фактора. Для Tw = 1 вся зона отрыва равномерно прогрета до слоя смешения высокой температурой стенки, тогда как для Tw < 0.3 температура в области отрыва значительно меньше, чем в слое смешения (рис. 6).
Таблица 3
Значение номера стадии отрыва для угла 0 = 10° в зависимости от числа Ке для различных значений ш, у и
температурного фактора Т„
Re T 1 w 0.1 ^ = 0Re% T = w 1
Y = 5/3 Y = 1.4 Y = 5/3 Y = 1.4
ю = 0 ю = 0.5 ю = 1 ю = 1 ю = 0 ю = 0.5 ю = 1 ю = 1
1 • 104 0 0 0 0 1.75 1 1 0 0
3 • 104 1 1 0 0 2.30 1 1 1 1
1 • 105 2 2 1 1 3.10 2 2 1 1
3 • 105 3 3 2 2 4.08 3 3 2 2
1 • 106 4 4 3 3 5.51 4* 4* 3 3
3 • 106 5 4 7.26 5* 4*
1 • 107 6 9.81 6
Результаты расчетов для у = 1.4 (табл. 3, ю = 1, Pr = 1) также приводят к смещению стадий отрыва к большим значениям числа Re, но качественно картина развития отрыва сохраняется.
На рис. 7 приведены поля давления внутри отрывной области, а на рис. 8 — профили давления в угловой точке для различных стадий отрыва. Для Tw = 0.1 видно, что если для 2-й (и отчасти 3-й) стадии отрыва (рис. 7, а, б) давление поперек отрыва остается постоянным, то для 4 — 6-й стадий (рис. 7, в, г, д) перестает быть постоянным (рис. 8). Для Tw = 1 давление поперек отрыва в угловой точке перестает быть постоянным при заметно большем числе Re (рис. 7, к), при котором и происходит образование вихрей. Аналогичная ситуация реализуется при изменении угла 0 при фиксированном числе Re [17, 28].
Определенный интерес представляет сравнение полученных результатов с решениями [14], полученными в рамках уравнений Прандтля (р = const, др/ду = 0). Так же, как и в [14], зарождение вторичного отрыва происходит при значении параметра подобия а ~ ^0 ~ 4 — 5 (в работе [14] а ~ 4.5). Для диапазона значений параметра а = 4 ■ 8 величина поверхностного напряжения трения (нормированного, как в [14]) не опускается ниже значения -1, тогда как в [14] она гораздо ниже. Значительно меньше, чем в работе [14], аналогичные провалы давления. Эти расхождения растут по мере роста поперечных градиентов давления в решении уравнений Навье — Стокса (рис. 8).
Отметим существенные отличия решений уравнений Навье — Стокса для различных температурных факторов. Для параметра подобия а = 7.5 в работе [14] были получены решения с большим размером вторичного отрыва (~1/2 размера первичного отрыва), а в угловой точке также наблюдался небольшой третичный отрыв. Похожая ситуация с вторичным отрывом большого размера реализуется для 4-й и 5-й стадий отрыва при высоком температурном факторе Tw = 1, где практически во всей области протяженного вторичного отрыва отсутствуют поперечные градиенты давления, и только в области присоединения наблюдается развитие поперечных градиентов давления и вихрей (рис. 7, з, и). Для малого температурного фактора Tw < 0.3 зарождение и развитие вторичного отрыва (3 — 5 стадии) сопровождается появлением поперечных градиентов давления и вихрей во всей области вторичного отрыва (рис. 7, б, в).
Рис. 7. Поле давления внутри области отрыва для угла 0 = 20°, температурного фактора Т„ = 0.1 (а, б, в, г, д) и Т„ = 1
05, 106 (соответственн нулевая линия тока
(е, ж, з, и, к) для чисел Яе = 104, 3 ■ 104, 105, 3 ■ 105, 106 (соответственно 2, 3, 4, 5, 6 стадии отрыва). Сплошная линия —
Рис. 8. Профили давленияр по нормали п в угловой точке для угла 0 = 20°, температурного фактора Тк = 0.1 и 1 при Яе = 104, 3 • 104, 105, 3 • 105, 106 (кривые 2, 3, 4, 5, 6 для соответствующих стадий отрыва)
Для понимания сложной картины результатов с вторичным отрывом была проведена тщательная верификация решений на различных сетках (рис. 9). Для высоких (4 — 6) стадий отрыва при недостаточном разрешении аналитические сетки значительно занижали размеры области отрыва. Чем выше температурный фактор и номер стадии отрыва, тем значительнее эти расхождения. Адаптивные сетки лишены этого недостатка. Однако правильная картина течения внутри отрыва (вихри) получается только на сетках большой размерности (1600 X 200, 1600 X 400). Поэтому сначала на сетках малой размерности устанавливаются внешнее течение и размеры зоны отрыва, а на сетках большой размерности уточняется структура отрывного течения (рис. 10). Результаты на аналитических сетках и адаптивных сетках большой размерности (1600 X 400) близки между собой. Однако для аналитических сеток это требует тщательного подбора сгущений (которые заранее неизвестны) и очень длительной процедуры счета (сетки максимально большой размерности). На грубых сетках получаются заведомо неверные размеры области отрыва.
На более подробных сетках не уточняется решение, а каждый раз получается новое, отличное от предыдущего. Методика применения адаптивных сеток лишена этих недостатков, а результаты получаются достаточно быстро и естественно. Применение динамических адаптивных сеток значительно упрощает и ускоряет процедуру получения решения, особенно для высоких (нестационарных) стадий развития отрыва.
Значительное расхождение решений, полученных на аналитических сетках аналогично результатам, приведенным на рис. 9, наблюдается для всех углов 0 при числах Яе, соответствующих 5-й и 6-й стадиям отрыва. При этом для решений на адаптивных сетках значительных расхождений не наблюдается. Для меньших номеров стадий развития отрыва эти расхождения соответственно меньше и наблюдаются для более грубых аналитических сеток.
Отметим особый характер сеточной диссипации и ее влияние на получаемое решение. При наличии в области отрыва крупномасштабных вихревых структур (для 5-й и особенно
Рис. 10. Поля плотности р внутри области отрыва (5-я стадия) для угла 0 = 20°, Яе = 3 • 105, = 0.1, полученные на динамических адаптивных сетках: 200 X 200 (а), 400 X 200 (б), 800 X 200 (в), 1600 X 200 (г)
для 6-й стадии) и при недостаточном разрешении получаемое решение сильно зависит от размерности сетки (рис. 10). Поэтому и для адаптивных сеток при получении решения в 5-й (рис. 10, г) и 6-й (см. рис. 1, г) стадиях из-за сложной картины вихревых структур внутри области отрыва необходимо применять сетки c большим количеством расчетных узлов.
Заключение. На основе решения уравнений Навье — Стокса проведено детальное численное исследование отрывного сверхзвукового (М„ = 5) течения в угле сжатия. Показано, что параметр ^0 = 0ЯеШ является определяющим при образовании как первичного, так и вторичного отрыва. Классифицированы стадии отрыва. При больших значениях параметра ^0 происходит
Рис. 9. Распределение давления по х-координате угла 0 = 10° для Яе = 3 • 106 и Т„ = 1 (5-я стадия):
сплошные линии 1 — для адаптивных сеток 801 x 201, 1601 x 201, 801 x 401, 1601 x 401; пунктир — для аналитических сеток 1601 x 401, 801 x 401, 1601 x 201, 801 x 201 (кривые 2, 3, 4, 5 соответственно)
значительное усложнение течения в области отрыва: образование вторичного отрыва, появление поперечных градиентов давления, стационарных и нестационарных вихревых структур. Выявлено существенное влияние температурного фактора на структуру вторичного отрыва.
Автор выражает признательность за творческое обсуждение и ценные замечания Нейланду В. Я. и Соколову Л. А., а также Егорову И. В. и Иванову Д. В. за предоставленные коды вычислительных программ.
Работа выполнена при поддержке грантов РФФИ (№ 07-08-000124 и № 06-08-01558а).
ЛИТЕРАТУРА
1. Нейланд В. Я., Куканова Н. И. Исследование течений со срывными зонами // Обзор БНИ ЦАГИ. — 1965. № 129.
2. Лапин Ю. В., Лойцянский Л. Г., Лунькин Ю. П., Нейланд В. Я., Сычев В. В., Тирский Г. А. Динамика вязких жидкостей и газов, теория ламинарных и турбулентных пограничных слоев: В сб. «Механика в СССР за 50 лет». — М.: Наука, 1970. Т. 2.
3. Charwat A. E. Supersonic flows with imbedded separated regions // Adv. Heat Transfer. N.Y.L. Acad. Press. 1970. V. 6.
4. Чжен П. Отрывные течения. — М.: Мир, 1972. Т. 1; 1973. Т. 2. 1973. Т. 3.
5. Chapman D. R. An analysis of base pressure at supersonic velocities and comparison with experiment // NASA Rep. 1951. № 1051.
6. Korst H. H. A theory for base pressure in transonic and supersonic flow // J. Appl. Mech. 1956. V. 23, N 4.
7. Нейланд В. Я. К асимптотической теории плоских стационарных сверхзвуковых течений со срывными зонами // Изв. АН СССР. МЖГ. 1970. № 3.
8. Нейланд В. Я. Сверхзвуковое течение вязкого газа вблизи точки отрыва: Сб. аннотаций докладов 3-го Всесоюзного съезда по теорет. и прикл. механике. — М.: Наука, 1968.
9. Нейланд В. Я. К теории отрыва ламинарного пограничного слоя в сверхзвуковом потоке // Изв. АН СССР. МЖГ. 1969. № 4.
10. Stewartson K., Williams P. G. Self-induced separation // Proc. Roy. Soc. — London. Ser.A. 1969. V. 312, N 1509.
11. Stewartson K. On laminar boundary layers near corners // Quart J. Mech. Appl. Math. 1970.V. 23, N 2.
12. Нейланд В. Я. Течение за точкой отрыва пограничного слоя в сверхзвуковом потоке // Изв. АН СССР. МЖГ. 1971. № 3.
13. Smith F. T., Khorrami A. F. The interactive breakdown in supersonic ramp flow // J. Fluid. Mech. 1991. V. 224.
14. K o ro l e v G. L., Gaj j ar J. S. B., Rub an A. I. Once again on the supersonic flow separation near a corner // J. Fluid Mech. 2002. V. 463.
15. Нейланд В. Я., Боголепов В. В., Дудин Г. Н., Липатов И. И. Асимптотическая теория сверхзвуковых течений вязкого газа. — М.: Физматлит, 2004.
16. N e il and V. Y a., B o g o lep o v V. V., D udin G. N., L ip ato v 1.1. Asymptotic theory of supersonic viscous gas flows // The Netherland: Elsevier; Oxford, 2007.
17. Нейланд В. Я., Соколов Л. А., Шведченко В. В. Влияние температурного фактора на структуру отрывного течения в сверхзвуковом потоке газа // Изв. РАН. МЖГ. 2008. № 5.
18. Егоров И. В., Зайцев О. Л. Об одном подходе к численному решению двумерных уравнений Навье — Стокса методом сквозного счета // ЖВММФ. 1991. Т. 31, № 3.
19. Бабаев И. Ю., Башкин В. А., Егоров И. В. Численное решение уравнений Навье — Стокса с использованием итерационных методов вариационного типа // ЖВММФ. 1994. Т. 34, № 11.
20. Башкин В. А., Егоров И. В., Иванов Д. В. Применение метода Ньютона к расчету внутренних сверхзвуковых отрывных течений // ПМТФ. 1997. Т. 38, № 1.
21. Годунов С. К. Разностный метод численного расчета разрывных течений уравнений гидродинамики // Мат. c6. 1959. Т. 47, № 3.
22. Roe P. L. Approximate Rieman solvers, parameter vectors, and difference scheme // J. Comp. Phys. 1981. V. 43.
23. Колган В. П. Применение принципа минимальных значений производной к построению конечно-разностных схем для расчета разрывных решений газовой динамики // Ученые записки ЦАГИ. 1972. Т. 3, № 6.
24. Saad Y. and Shultz M. H. GMRES: a generalized minimal residual algorithm for solving nonsymmetric linear systems // SIAM J. Sci. and Stat. Comp. 1986. N 6.
25. Гильманов А. Н. Методы адаптивных сеток в задачах газовой динамики. — М.: Физматлит, 2000.
26. Андерсон Д., Таннехилл Дж., Плетчер Р. Вычислительная гидродинамика и теплообмен. — М.: Мир, 1990. T. 1, 2.
27. Stemmer C. and Adams N. A. Investigation of supersonic boundary layers by DNS // ECCOMAS 2004 Proceedings. 2004. V. II.
28. Neyland V. Ya., Sokolov L. A.,Shvedchenko V. V. Temperature factor effect on separated flow features in supersonic gas flow // Boundary & Interior Layers: Comp. & Asymptotic Methods. — Ireland. Jul 2008.
Рукопись поступила 29/XII2008 г.
