CC BY
65
_________УЧЕНЫЕ ЗАПИСКИ ЦАГИ
Том XXVIII 199 7 "
№3-4
УДК 533.6.011.8
ЧИСЛЕННОЕ ИССЛЕДОВАНИЕ ОБТЕКАНИЯ ТЕЛ С ПРОТОКОМ ГИПЕРЗВУКОВЫМ ПОТОКОМ РАЗРЕЖЕННОГО ГАЗА
А. И. Ерофеев, В. Д. Перминов
Методом Монте-Карло рассчитывается обтекание модели треугольного крыла с гондолой, имитирующей канал прямоточного двигателя. Исследуется зависимость суммарных аэродинамических характеристик, коэффициента расхода газа через канал и коэффициента сопротивления канала от числа Рейнольдса и угла атаки. Расчеты проводились при числе Маха М* = 20,2 и температурном факторе /*, = 0,1 для случаев открытого и закрытого каналов. Для исследования детальной картины течения газа в канале проведено также решение двумерной задачи об обтекании двух плоских пластин различной длины, образующих канал.
1. Постановка задачи. Метод расчета. В настоящее время имеется достаточно большое число исследований аэродинамических характеристик тел при их обтекании гиперзвуковым потоком разреженного газа (см., например, [1—3]). Спектр этих исследований чрезвычайно широк: от изучения обтекания простых тел (пластина, конус, сфера и т.д.) до расчета течений около летательных аппаратов типа «Space Shuttle» или «Буран»; состав потока при этом варьируется от простого одноатомного газа до смеси газов с учетом физико-химических превращений. Однако теоретические исследования обтекания летательных аппаратов с каналами, имитирующими, например, канал прямоточного двигателя, пока весьма немногочисленны. В работах [4], [5] рассмотрена двумерная задача о движении газа в канале, образованном двумя плоскими параллельными пластинами, и основное внимание уделено вопросам интерференции. В работе [6] рассмотрена более сложная задача об обтекании двух пластин разной длины; особое внимание уделено случаю, когда скачок уплотнения, сформировавшийся при обтекании более длинной пластины, падает на переднюю кромку короткой затупленной пластины. В отмеченных работах расчеты проводились методом Монте-Карло при некоторых фиксированных значениях параметров подобия: чисел Маха Мж, Кнудсена Кпда, температурного фактора
tw =TW / Tq и не рассматривался вопрос о расходе газа через канал. В работе [7] на примере задачи обтекания двух плоских параллельных пластин одинаковой длины проведено исследование течения газа в канале в зависимости от величины параметра разреженности — числа Рейнольдса Re0 = pJU^L / m-(7q) (здесь р*,, U„ — плотность и скорость газа в набегающем потоке, ц.(Г0) — коэффициент вязкости, рассчитанный при температуре торможения Tq, L — характерный линейный размер) для нескольких значений числа Маха и температурного фактора. Было показано, что расход газа через канал является немонотонной функцией числа Reg, а также отмечено сильное влияние температурного фактора на формирование течения газа в канале, что необходимо принимать во внимание, например, при моделировании натурных условий обтекания аппаратов с воздухозаборниками в лабораторных условиях.
В данной работе рассматривается трехмерная задача обтекания гиперзвуко-вым потоком разреженного газа простой модели летательного аппарата с каналом прямоточного двигателя — плоского треугольного крыла с гондолой (рис. 1). Рассматривается обтекание модели как с открытым, так и с закрытым каналом, причем в последнем случае канал с входной стороны закрывался пластиной, расположенной под углом 26,5° к плоскости крыла (в соответствии с профилем боковых стенок гондолы), а со стороны выхода — поперечной пластиной. Модель имеет следующие геометрические размеры: отношение размаха крыла b к длине L равно 2/3, длина нижней поверхности гондолы / составляет 1/3 длины крыла, отношение высоты гондолы А к / равно 1/4. ,
Решение задачи проводилось методом прямого статистического моделирования, подробное описание которого дано в [8]. Здесь отметим только, что в этом методе реальное течение газа моделируется движением ансамбля частиц в некоторой расчетной области. Расчетная область разбивается на ячейки, размер которых должен быть меньше местной длины свободного пробега частиц.В начальный момент времени область течения заполняется частицами, поступательные скорости и внутренняя энергия которых определяются по начальной функции распределения, соответствующей состоянию газа в невозмущенном потоке. Затем последовательно на каждом шаге по времени At проводятся свободное перемещение частиц и столкновение между ними, причем сталкиваться могут лишь частицы, находящиеся в одной геометрической ячейке. При движении частиц в расчетной области они могут
а) Треугольное кршо с гондолой
б) Канал из dSyx плоских пластин ____________________It_____________
_____________Л
Рис. 1
сталкиваться с поверхностью летательного аппарата, а также вылетать за пределы области. Вылетевшие частицы исключаются из дальнейшего рассмотрения, а на каждом временном шаге с границ области проводится вбрасывание частиц в соответствии с граничной функцией распределения. После некоторого количества шагов в системе устанавливается квазистационарное состояние, и с этого момента начинается сбор необходимой информации о поле течения, о потоках импульса и энергии на поверхности и других выходных параметрах.
В данной работе решение проводилось одним из вариантов метода Монте-Карло [9], в котором вероятность столкновения частиц в ячейке за время А/ определяется соотношением р = a(g)gAt/V, где ст(#) — полное сечение столкновения, g — относительная скорость пары частиц, V — объем ячейки. Сечение столкновения определялось по модели сфер переменного диаметра [10], основывающейся на степенном потенциале взаимодействия частиц 1/(г) ~г~у. Для этой модели зависимость коэффициента вязкости от температуры описывается соотношением ц(Г) ~ , где со = 1/2 ч- 2/у; в дальнейшем будет принято
V = 10; при этом значении V теоретическая зависимость коэффициента вязкости ц(Г) хорошо аппроксимирует экспериментальную для азота при температурах от 300 до 20000 К: погрешность аппроксимации справочных данных [11] не превышает 6%. Одним из основных требований при применении метода Монте-Карло является условие ё <Х, где с1 — размер ячейки, а X — местная средняя длина свободного пробега частиц. Практика применения метода показала, что результаты расчетов не зависят от размеров ячейки при й < Х/3; кроме того, выяснилось [12, 13], что такое жесткое требование необходимо выполнять в направлении наибольших градиентов функции распределения или ее моментов. Поскольку в различных зонах течения величины X могут сильно различаться, желательно иметь методику, позволяющую в процессе счета «подстраивать» расчетную сетку под местную длину свободного пробега частиц. Такая «подстройка» в данной работе проводится на основе данных о поле плотности (см. также [13], [14]). Действительно, средняя длина свободного пробега частиц для модели сфер переменного диаметра в равновесных условиях X / Хх ={Т / 7,00)2//'’(л00 / п). Если температура в поле течения Т >ТХ, то для оценок можно принять, что X / Х„ = п„ / и, и использовать это соотношение для определения размеров ячеек. Ниже адаптация размеров ячеек под поле плотности проводилась следующим образом: в начальный момент времени устанавливался некоторый базовый размер ячеек (расчеты проводились на декартовой сетке) и на текущей итерации определялось поле плотности; на следующей итерации, в случае если плотность газа в ячейке превышает значение п* (величина п* принималась равной 1,5 пх), то проводится дробление базовой ячейки в направлении той координатной оси, вдоль которой градиент плотности максимален. Если же ячейка примыкает к твердой поверхности, то дробление проводится по той координатной оси, вдоль которой максимальна проекция нормали к поверхности. Такое дробление
в определенном смысле обеспечивает и адаптацию сетки к обтекаемым поверхностям. Как правило, двух-трех итераций обычно бывает достаточно для получения сходящихся результатов.
, Расчеты проводились для газа, имеющего вращательные степени свободы. Обмен энергией между поступательной и вращательной модами описывался моделью Ларсена — Боргнакке [15] с постоянным параметром ф = 0,2 и 0,35. В этой модели параметр ф связан с параметром релаксации zr = tr / tc, где tr — время вращательной релаксации^ /е — среднее время свободного пробега молекул, соотношением ф(гг) = 2,06 / zr- Так что величина zr, характеризующая число столкновений, необходимых для установления равновесия между поступательными и вращательными степенями свободы, принималась равной 10 или 6. Контрольные расчеты при одних и тех же значениях числа Мм и температурного фактора показали, что влияние параметра zr (или ф) на аэродинамические характеристики и на расход газа через канал несущественно.
Поскольку исследование детальной картины трехмерных течений представляет собой достаточно трудоемкую задачу, для прояснения некоторых особенностей поля течения в трехмерном случае была рассмотрена также двумерная задача об обтекании двух плоских пластин различной длины (см. рис. 1, б). Рассматривались два случая: отношение длины верхней пластины lt к длине нижней пластины ld равно 2:1 и 2,5:1; расстояние между пластинами hc в обоих случаях составляло ld / 4.
2. Результаты расчетов. Расчеты двумерной и трехмерной задач проводились при числе = 20,2 и температурном факторе tw = 0,1. В качестве параметра разреженности использовалось число Reg. В качестве характерной длины Z в трехмерном случае принималась длина треугольной пластины, а при сопоставлении трехмерной и двумерной задач — длина нижней пластины поскольку длина верхней пластины в двумерном случае варьировалась. Аналбгично при определении аэродинамических характеристик за характерную площадь принималась площадь треугольной пластины в трехмерном Случае и площадь нижней пластины — при сопоставлении результатов. Рассмотрим сначала зависимость интегральных характеристик от параметров задачи. На рис. 2, 3 представлены результаты расчетов аэродинамических характеристик модели треугольного крыла с гондолой для случаев открытого (сплошные линии) и закрытого каналов (штриховые линии) для углов атаки а = 0, 5 и 30°. Из приведенных данных следует, что при малых углах атаки различие в коэффициенте сопротивления сХа, весьма значительное в свободномолекулярном пределе, указанном на графиках чертой слева, с увеличением числа Re0 практически исчезает. При а = 30°,напротив, с увеличением числа Re0 наблюдается некоторое расслоение зависимостей сХа (Reg) для закрытого и открытого каналов. Видно, что коэффициент сопротивления модели с открытым каналом
1,2
0,8
ОЛ
-о-
а—О
-0,5
0,5 1
Рис. 2
1,5 ІіЯе0
при малых углах атаки является немонотонной функцией числа 11е0, что характерно для обтекания аэродинамических острых тел [16]. В случае закрытого канала вклад в сХа пластины, закрывающей вход в
канал, оказывается значительным при малых числах К.е0 и немонотонность сХа (Яво) лишь слабо проявляется при а = 0. Что же касается коэффициента подъемной силы сУа (рис. 3), то при а = 5° во всем диапазоне расчетных чисел Ке0 для модели с закрытым каналом он примерно на 50% выше, чем для модели с открытым каналом. Для а = 30° ве-очень близкие в свободномолекулярном пределе, с
личины
ьУа-
увеличением числа Ке0 расслаиваются, причем сУа для модели с закрытым каналом становится примерно на 10% выше.
Известно, что при газодинамическом режиме течение в канале существенным образом зависит от того, попадет ли скачок уплотнения, образующийся при обтекании верхней пластины, внутрь канала или проходит ниже входного сечения канала. С другой стороны, из-за отсутствия скачка при свободномолекулярном режиме имеется только один тип течения. Ниже делается попытка исследовать особенности течения в канале при промежуточном режиме. При этом ввиду трудоемкости решения трехмерной задачи (особенно при сравнительно больших числах Рейнольдса) исследование проводится с привлечением результатов решения двумерной задачи обтекания двух пластин разной длины. Поскольку результаты решения трехмерной задачи показывают, что линия максимальных значений плотности при 11е0 = 28 (линия 1 на рис. 4) в плоскости симметрии проходит ниже входного сечения канала, отношение I, / 1а для двумерного случая было выбрано таким, чтобы при близком значении числа Ке0 в одном случае (линия 2, Яе0 =40, -/,//</= 2:1) эта линия пересекала нижнюю пластину канала, а в другом (ли- ' "—3
ния 3, Яе0 = 31,6, 1( / 1а = 2,5:1) — прохо- 7
дила ниже входного сечения канала. На Рис. 4
а) а=0
1?Яс0 1
Рис. 6
1рЯев 2
Рис. 5
рис. 5, 6 приведены расчетные зависимости коэффициента сопротивления сх канала и коэффициента расхода газа через канал Cj=J/ от числа Кео для трехмерной и двумерной задач (/«, — расход газа через трубку тока в невозмущенном потоке с поперечным размером, равным проекции площади входного сечения канала на плоскость, перпендикулярную вектору скорости). На этих рисунках линия 1 соответствует решению трехмерной задачи, 2 — двумерной для /, / 1а =2:1, 3 — /г / 1й =2,5:1, черточками слева с соответствующими маркерами обозначены результаты для свободномолекулярного случая (для а = 5° свободномолекулярные пределы оказались практически одинаковыми для трехмерной и двумерной задач). Напомним, что при проведении такого сравнения в качестве характерных величин длины и площади принимались длина и площадь нижней поверхности канала. Из приведенных данных следует, что при а = 0 и 5° данные для двумерной задачи при отношении длин 7,//^ =2,5:1 во многих случаях хорошо согласуются с результатами трехмерной задачи не только качественно, но и количественно, в то время как при /,//^ = 2:1 различия весьма существенны, особенно для коэффициента расхода газа через канал. При а = 30° хорошее согласование расчетных данных дву- и трехмерной задач имеет место при /, / 1а = 2:1.
На рис. 7 представлены поля плотности в нескольких сечениях внутри канала и вблизи входа в канал (двумерный случай, а = О, Ке0 = ЮО). Видна несимметричность поля профиля плотности внутри канала, в случае когда скачок уплотнения проходит внутрь канала. В этом случае, как показывает анализ полей скорости и плотности, поток газа дважды «отражается» от стенок канала — сначала от нижней поверхности, а затем от верхней. Этот вывод подтверждается распределением локальных силовых характеристик на внутренних поверхностях канала (рис. 8). Напротив, когда в канал попадает более разреженный поток газа, происходит большая симметризация течения, так что на выходе поток газа становится почти симметричным относительно сре-
-/ xfhe'0 7 Z 3 4
lt:la=Z,S:I
0 Z Чп/п.
■ ■ ■ ' ■ '
0 2 Чп/Яс
» | « » * '
1
Ж x/h=l 7 1 1 J Г 4 1
~1,S
Рис. 7
lt:ld=Z,5.1
Рис. 8
динной поверхности канала. На рис. 8 приведены эпюры трения Су и нормальной к поверхности силы с„, отнесенные к скоростному напору, при обтекании плоского канала при а = 0 (показана только область самого канала, линии 1 соответствуют нижней поверхности канала; 2 — верхней; штриховка нанесена только для сп; масштаб величин отмечен в нижней части рисунка вертикальным отрезком и относится к обеим величинам.
Аналогичная ситуация имеет место и при а = 5°.
При а = 30° картина обтекания существенно изменяется по сравнению с малыми углами атаки (рис. 9). В двумерной задаче (lt / ld = 2:1) хорошо видна структура ударной волны, на которой достигается отношение плотностей л/=5,3 (рис. 9, а — Re0 = 40, рис. 9,6 — Re0 = 100). Представленные на рис. 9, а профили плотности для трехмерного случая вблизи плоскости симметрии (крестики, Reg = 28) и двумерного случая хорошо согласуются между собой. Из этих данных следует, что в обоих случаях в канал попадает плотный газ, прошедший
О 5 10 п/п.
1,1,1 I I '
а) Не 0=М
О $ 10л/п.
Ц I I | '
Рис. 9
через ударную волну, что и определяет согласование результатов расчетов коэффициентов расхода газа через канал и аэродинамических коэффициентов. Внутри канала, гак и при малых углах атаки, дважды происходит «отражение» ударной волны от стенок канала, что проявляется в увеличении нормальной силы сп (рис. 9, в) на стенках канала: сначала на нижней, а затем на верхней поверхности. Как видно из рис. 9, поток газа на выходе практически симметричен относительно срединной поверхности канала.
В заключение отметим основные качественные результаты, полученные в работе. Расчетные данные доказывают, что аэродинамические характеристики треугольного крыла с открытым и закрытым каналами воздухозаборника в переходном режиме обтекания могут заметно отличаться. В большей степени это относится к коэффициенту подъемной силы; что же касается коэффициента сопротивления, то наибольшие отличия имеют место при малых углах атаки в свободномолекулярном режиме обтекания. Показано, что на режимах обтекания, когда сформировался головной скачок уплотнения, течение в канале существенно зависит от положения фронта ударной волны относительно входного сечения канала. Практически во всем рассмотренном диапазоне чисел Рейнольдса при правильном моделировании этого положения результаты решения двумерной задачи обтекания двух пластин различной длины не только качественно, но и количественно согласуются с решением трехмерной задачи, описывающей течение в канале.
1. Bird G. A. Monte Cario simulation in an engineering context//Progrcss in Aeronautics and Astronautics.— 1991. Vol. 74.
2. Иванов М. С. Статистическое моделирование гиперзвуковых течений разреженного газа//Докт. дисс., ИТФ СО РАН.— 1993.
3. Erofeev A. I., Nikolaev К. V., Perminov V. D. Some feature of hypersonic rarefied gas 3D-flows. Proc. 19-th Intern. Symp. Rarefied Gas Dynamics. Vol. 2. Oxford.— 1995.
4. Wilmoth R. G. Interference effects on the hypersonic rarefied flow about a flat plate. Proc. 16-th Intern. Symp. Rarefied Gas Dynamics. Pasadena.— 1988.
5. Yasuhara М., Nakamura Y., Tanaka J. Monte Cario simulation of flow into channel with sharp leading edge. Proc. 16-th Intern. Symp. Rarefied Gas Dynamics. Pasadena.— 1988.
6. Wilmoth R. G. Adaptive domain decomposition for Monte Carlo simulations on parallel processors. Proc. 17-th Intern. Symp. Rarefied Gas Dynamics. Aachen.—1990.
7. Ерофеев А. И. Обтекание двух плоских пластин гиперзвуковым потоком разреженного газа//Ученые записки ЦАГИ.—1996. Т. XXVII, № 1-2.
8. Берд Г. Молекулярная газовая динамика,— М.: Мир.— 1981.
9. БелоцерковскийО. М., ЯницкийВ. Е. Статистический метод частиц в ячейках для решения задач динамики разреженного газа// Ж. вычисл. матем. и матем. физ,— 1975. Т. 15, № 5, 6.
10. Е р о ф е е в А. И. О моделировании межмолекулярного взаимодействия при решении уравнения Больцмана методом Монте-Карло//Изв. АН СССР, МЖГ.-1977, № 6.
11. Гордеев О. А., Калинин А. П., Комов А. Л., Люстер-ник В. Е., Самуйлов Е. В., Соколова А. А., Фокин Л. Р. Потенциалы взаимодействия, упругие течения, интегралы столкновений компонентов воздуха для температур до 20 ООО К (Методы определения, рекомендуемые данные)//Обзоры по теплофизическим свойствам веществ, № 5 (55).— М. ИВТ АН СССР.- 1985.
12. Николаев К. В. Аэродинамические и тепловые характеристики обтекания затупленных тел разреженным газом//Канд. дисс.— М.: МФТИ.— 1990.
13. Тау 1 о г J. С., М о s s J. N., Н ass a n Н. A. Study of hypersonic flow past sharp cones//AIAA-89-1713.
14. Olynic D. P., Hass an H. A., Moss J. N. Grid generation and adaptation for direct simulation Monte Cario method//AIAA Thermophysics, Plasmadynamics and Lasers Conference, June 27—29.— 1988, San Antonio, Texas.
15. Larsen P. S., Borgnakke C. Statistical collision model for simulating polyatomic gas with restricted energy exchange.— Proc. 9-th Intern. Symp. Rarefied Gas Dynamics, DFVLR-Press, Potz-Wahn. Vol. 1.— 1974.
16. Гусев В. H., Ерофеев А. И., Климова Т. В., Перепу-хов В. А., Рябов В. В., Толстых А. И. Теоретические и экспериментальные исследования обтекания тел простой формы гиперзвуковым потоком разреженного газа//Труды ЦАГИ.— 1977. Вып. 1855.
Рукопись поступила 17/11996 г.
