УЧЕНЫЕ 3 А П И С К И ЦА Г И
Т о м XX
1989
№ 4
УДК 533.6.071.4.011.55 532.525.011.55
ИССЛЕДОВАНИЕ ТЕЧЕНИЯ ГАЗА В ГИПЕРЗВУКОВЫХ СОПЛАХ ПРИ БОЛЬШИХ ЧИСЛАХ РЕЙНОЛЬДСА НА ОСНОВЕ УПРОЩЕННЫХ УРАВНЕНИЙ НАВЬЕ—СТОКСА
В. Я. Безменов, А. П. Быркин, П. И. Горенбух,
В. А. Сабельников, Т. А. Тимофеева, А. И. Толстых
Приведены результаты численных расчетов турбулентных течений гелия и воздуха в конических и профилированных гиперзвуковых соплах, проведенных на основе упрощенных уравнений Навье—Стокса и однопараметрической модели турбулентности, а также экспериментальное их исследование при числах ¡йе= (0,3-ь 1,4) • 107. Получено удовлетворительное согласование расчетных и экспериментальных полей чисел М в невязком ядре потока.
Предложенная методика может быть использована для диагностики полей чисел М в рабочей части гиперзвуковых аэродинамических труб.
При численном решении задач аэромеханики, отвечающих большим числам 1?е, в общем случае необходимо учитывать эффекты турбулентного переноса. Для описания течения газа при этом часто используются полные или упрощенные усредненные уравнения Навье— Стокса, замкнутые с помощью той или иной модели турбулентности (см., например, [1-3]).
В работах [4, 5] на основе упрощенных уравнений Навье—Стокса численно исследовалось течение вязкого газа в гиперзвуковых соплах (М = 6н-18) при ламинарном режиме течения в пограничном слое, когда характерное число Ие, определенное по параметрам газа в ядре потока на выходе и длине сопла, было Не<10°. При этом получено удовлетворительное соответствие результатов расчета и эксперимента.
Целью данной работы является численное моделирование, с использованием модели турбулентной вязкости [7], течений гелия в конических гиперзвуковых соплах при Яе~107, воздуха в конических и профилированном соплах при 1?е~ (3-н5) • 10е и сопоставление результатов расчетов с экспериментальными данными.
1. Исходная система полных уравнений Навье—Стокса для случая стационарного течения газа имеет вид
О)
где х, у — расстояния, измеряемые вдоль оси и по нормали к оси канала; Б, в, Н, Т — векторы-столбцы:
F = у'
pU pV
ри2-т-р prou.
puv , G=y pv2 +p
Il . + v2 \ р И |А + 2 J P v[h± 2 j
0 Tx
0 T 2
, t =
p Ts
0 T4
Здесь 7\ = 0, выражения для Т2, Т3, Г4 можно найти, например, в работе [5].
В уравнениях (1) все величины безразмерные; и, v — компоненты скорости в системе координат х, у, р — плотность; р — давление; h — энтальпия; значение / = 0 соответствует плоскому, /= 1—осесимметричному случаю.
Газ предполагается совершенным, т. е. уравнение состояния имеет
вид
1
рЛ,
(к — показатель адиабаты), коэффициенты вязкости р, = ^т+ц( и теплопроводности %=Хт+Ь, входящие в выражения для Т2, Т3, Т4, включают составляющие, обусловленные соответственно молекулярным и пульсационным переносом. Предполагается, что число Рг( = Ргт=Рг; число Ке! = р\щуул/\1т1 (индекс 1 отвечает значениям размерных величин на оси сопла в критическом сечении, индекс т — стенке). Зависимость молекулярного коэффициента вязкости цт от температуры предполагается заданной по степенному закону или формуле Сезерленда.
Для решения рассматриваемой задачи в области с криволинейной границей ут(х), где Ую(х)—контур сверхзвуковой части сопла, введем новые независимые переменные
Ут (-*) '
При этом по аналогии с работами [4, 5] будем решать не полную систему уравнений Навье—Стокса (1), а упрощенную (см. [4]):
ÔF . ÔG ___________ ¡г . -I гр
аа дт\ ~Re, -
(2)
получаемую из системы (1) (записанной в переменных т]) при исключении из нее вязких членов со вторыми производными по | и смешанными производными ПО I И Т).
Численное решение системы (2) находится по маршевой схеме [4, 5], использующей компактные аппроксимации третьего порядка [6], при следующих начальных и граничных условиях:
Р = РСЇІ). и = и(ц), v=^v(ti), h = h(i\) при S = 0,
до да dh п п п
--------------------п "■ — О ПОИ 7) = О,
ib)
да dh ___________________„
дт| dtj ’
V ■■
и ■■
V ■
:О, h = hw(l)
при
при
7)=1.
Выписанные граничные условия соответствуют заданию профилей газодинамических величин в начальном сечении | = 0 (являющемся в расчетах критическим сечением сопла), условию симметрии течения относительно оси |, прилипанию газа на стенке и заданному закону изменения по | температуры стенки.
2. Для определения коэффициента турбулентной вязкости ^ используется полуэмпирическое уравнение переноса [7]
дч, , дч,
ри —- + pv-------------------
ґ дх ду
+ apv,
ди
ду
+
ajL + vA.V
дх ду )
Р Гі
+
(fa + vm)
(4)
где vm = р^/р» V/ = !X</P> У — расстояние, отсчитываемое от стенки; «. Р, 7> *> эмпирические постоянные.
В монографии [7] для постоянных рекомендованы следующие значения: а = 0,2; р = 0,06; и = 2; g = 2/3; у = 50. Как показали пробные расчеты, при таком наборе значений постоянных наблюдается значительное расхождение расчетных и экспериментальных данных. Причина расхождения заключается в том, что значение диссипативной постоянной у в уравнении (4) сильно завышено. Это приводило к существенному занижению турбулентной вязкости fit, вследствие чего поля газодинамических параметров в сопле практически не отличались от полей для ламинарного режима течения.
Поэтому была сделана корректировка эмпирических постоянных. Значение постоянной у можно оценить, если рассмотреть течение в вязком подслое турбулентного пограничного слоя несжимаемой жидкости, где, как известно, vi<cvm (см., например, [8]), справедлива асимптотика vt^yn, п = 3 или 4.
Асимптотика уравнения (4) в этой области имеет вид
v ----- — -r z == 0,
m dy* 1 y*
откуда находим, что y = n{n—1), т. e. y = 6 при n = 3 и y=12 при n = 4.
Еще одна связь между эмпирическими постоянными следует из рассмотрения логарифмического профиля пограничного слоя,
в котором vi vm и справедливы соотношения ^ = Ли*у, ^у- = ~~ .
постоянная Карма-
»*=у
( ди \ Vm V ду
скорость трения, & = 0,41 на [8]. В логарифмическом слое уравнение (4) имеет вид
d /- d dy
СО
du
dy
-ГР4?
= 0,
откуда находим, что
* +
(5)
Следует подчеркнуть, что выражения для постоянных у и р получены для идеализированных случаев и, тем самым, носят оценочный характер.
Окончательные значения постоянных .выбирались из условия соответствия результатов расчета с экспериментальными данными.
Эти значения приведены в следующем разделе. Необходимо отметить, что экспериментальные данные получены для ядра потока, на развитие которого турбулентная вязкость оказывает хотя и заметное, но не прямое влияние. Поэтому хорошее соответствие результатов расчета и экспериментальных данных в ядре, вообще говоря, не означает, что рассматриваемая модель турбулентности будет с достаточночной точностью описывать профили газодинамических параметров в пограничных слоях*на стенках сопла.
Проделанные расчеты позволяют заключить, что, по крайней мере, суммарное влияние турбулентных пограничных слоев на формирование течения в невязком ядре сопл учитывается с достаточной точностью. Чтобы сделать вывод о применимости модели турбулентности в пристеночной области, необходимо получение надежных экспериментальных данных в этой области и проведение дополнительных расчетов.
Разностное представление уравнения (4) отвечает первому порядку аппроксимации в продольном и поперечном направлениях. Последнее не должно сильно сказаться на точности выходных результатов, поскольку сам коэффициент турбулентной вязкости определялся в рамках приближенной модели турбулентности.
3. Проведен расчет течения гелия (х = 5/3) в коническом гиперзву-ковом сопле (сопло 1) аэродинамической трубы [9].
Радиус критического сечения сопла г/ю 1 = 6,7 мм, полуугол раствора стенок 0и> = 6°, диаметр выходного сечения 300 мм. В окрестности критического сечения контур Уи?(х) (рис. 1) выполнен в виде дуги окружности радиуса ^=10 уш4.
Проведенное экспериментальное исследование данного сопла при параметрах торможения р = 37-105 Па и 70=ЗООК показало, что на его выходе реализуется течение с числом М в ядре потока, равным примерно 17. При этом число Яе, определенное по параметрам потока газа на оси в выходном Сечении и длине сопла, составляет примерно 1,4-107. При таких числах Ие имеет место турбулентный режим течения в пограничном слое у стенки (см., например, [10]).
Для указанных параметров торможения в рассматриваемом сопле сначала был проведен расчет течения гелия с учетом влияния только молекулярной вязкости цт во всей расчетной области (ц,( = 0). При этом, как и в работе [4], температура стенки сопла предполагалась близкой к температуре теплоизолированной (Тш = 7\. = 0,83 Т0). Расчет проводился в предположении, что в критическом сечении сопла пограничный слой отсутствует и имеет место однородный поток газа с числом М! = 1,01. Зависимость коэффициента молекулярной вязко-
о / 'т0,647\
сти цт от температуры предполагалась степенной ((Ат~ * ), число
Ргт=0,68. При дискретизации уравнений (2) и (4) осуществлялось сгущение узлов сетки в пристеночной области. С этой целью делалось преобразование переменной г] по формуле
С= 1 —
1п [!+»(!-4)1
1п(1 +8)
где 8 = 500.
По переменной £ разностная сетка имеет постоянный шаг. Число узлов в поперечном направлении Л^ = 81. При указанных значениях параметров сетки б и N пристеночный шаг Ат] = 0,00016.
Результаты расчета, проведенного на основе метода [4], получены при одной глобальной итерации по всей длине сопла; в этом случае при Мж<1 + е (Мж — число М в расчетной точке, определенное по продольной составляющей скорости, е = 0,1) член д/д1(уюур) в уравнении импульса в продольном направлении системы (2) полагался равным нулю. Однако, как показали расчеты [5], отличие в результатах, полученных при одной и трех глобальных итерациях, незначительно.
В качестве граничного условия на стенке при 0<х/г/№1<2 задавалось условие непротекания, при 2<.х!ут 1<20 — условие скольжения газа при изменении ыте(х) и кц,(х) по плавным законам, при лг/г/иц> >20 — условие прилипания и равенства ТЮ = ТГ.
По результатам данного расчета для ламинарного режима течения в пограничном слое (цг = 0) на рис. 1 штриховой линией приведено распределение по х числа М на оси — Ж0(х). Для сравнения на рисунке нанесены экспериментальные значения чисел М на оси в трех сечениях сопла 1. (Экспериментальные значения чисел М находились на основании показаний приемника полного давления). Видно, что имеется заметное отличие между расчетными (^г = 0) и экспериментальными значениями чисел М на оси сопла.
Сплошной линией на рис. 1 для сопла 1 представлены результаты расчета течения газа при упомянутых выше предположениях с учетом вклада в коэффициенты вязкости и теплопроводности турбулентной составляющей (Ргг = 0,68, Ту, = Тг = 0,88 Т0). Этот расчет выполнен при
нахождения V; из уравнения (4) со следующими значениями констант: и = 2 и | = 2/3. Значение |3 рассчитывалось по формуле (5), значение у варьировалось и принято окончательно равным 12 из условия согласования расчетных и экспериментальных данных по числу М на оси сопла.
Заметим, что при проведении этого расчета начальные условия для определения величины V« из уравнения (4) задавались в сечении сопла х/у-и, 1= 11,2 в виде постоянного распределения V« по т|^*(т)) =0,1).. Кроме того, величина а на начальном участке интегрирования уравнения (4) от сечения х/уи, 1=11,2 до сечения х/у101 = 2Ъ задавалась изменяющейся по линейному закону от а = 0,05 до а=0,2. Такая процедура «мягкого» подключения турбулентной вязкости обеспечивала устойчивый счет.
На рис. 1 приведены рассчитанные профили и (у) в нескольких сечениях сопла 1, а также сопоставлены расчетные и экспериментальные профили М (у) В окрестности ВЫХОДНОГО сечения (х1ую 1 = 205,6). Видно удовлетворительное их согласование, в том числе по размеру невязкого ядра течения.
На рис. 2 представлены рассчитанные профили газодинамических величин и, V, к и р по т), отвечающие окрестности выходного сечения, на основании которых можно судить о структуре течения газа в сопле. Для более детального рассмотрения структуры потока в пограничном слое на рисунке приведены профили и(г), где 2=1—г]— расстояние, отсчитываемое от стенки сопла. Здесь же показан профиль величины vt/vm(z), дающий представление о соотношении молекулярной и турбулентной вязкости.
Аналогичные расчеты выполнены для течения гелия в коническом сопле с тем же диаметром выходного сечения 300 мм, но с радиусом критического сечения 4,8 мм (сопло 2). Эксперименты показали, что при параметрах торможения р0=39-Ю5 Па и 7’о=ЗООК на срезе этого сопла реализуется течение с числом М в ядре потока примерно 20,4; число Ие»1,1-107.
Приведенные на рис. 1 сопоставления результатов расчетов с экспериментом для сопла 2 также свидетельствуют об удовлетворительном их согласовании. Однако в случае сопла 2 в окрестности границы пограничного слоя наблюдается большее расхождение в профилях Щу).
Представляло интерес провести опробование использованных выше констант применяемой однопараметрической модели турбулентности [7], рассчитав гиперзвуковые течения в соплах с видоизмененными условиями. Например, выполнить расчеты течения воздуха в конических соплах, экспериментальному исследованию потока в которых посвящена работа [10].
Исследованные в работе [10] сопла имели полуугол раствора стенок Эю = 6°, диаметр выходного сечения 150 мм; одно сопло имело радиус критического сечения 2,07 мм (сопло У), другое — 3,24 мм (сопло 2). При их испытании параметры воздуха в форкамере были: р0= = 100-105 Па, Т0= 1040 К. Для обеспечения температуры стенки сопла, близкой к постоянной, сопла помещались в ванну с водой, которая в случае сопла 1 доводилась нагревателем до кипения (Тю^373 К), а в случае сопла 2 пропускалась через ванну в сеть (Тю^290 К). При указанных условиях реализуемые в опытах числа Яе составляли ~3-10® (сопло 1) и -5-106 (сопло 2).
Расчеты течения воздуха в указанных соплах проведены по аналогии с расчетами течений в гелиевых соплах (некоторое отличие состояло лишь в том, что процедура «мягкого» подключения трубулент-ной вязкости начиналась с сечения х!ую1^2,А\ кроме того, зависимость \лт(Т) задавалась по формуле Сезерленда, Ргг = Ргт = 0,71).
На рис. 3 представлены полученные распределения Мо(*), которые сопоставлены с экспериментальными. В нижней части рисунка для нескольких сечений (сопло 1) показаны расчетные профили и (у) и профили М(у), отвечающие выходному сечению. В целом также получено удовлетворительное совпадение расчетных и экспериментальных данных.
На рис. 4 приведены результаты расчетного и экспериментального исследования течения воздуха в профилированном сопле, рассчитанном на число М= 10 (испытание данного сопла проведено А. П. Быркиным, В. П. Верховским, А. Ю. Сосуновым). Радиус критического сечения сопла 2,5 мм, диаметр выходного сечения 150 мм. На рисунке показаны невязкий контур сопла (у(х)), контур сопла, определенный с учетом вязкости (уго(х)), характеристический ромб. Штрихпунктирной линией показано распределение М0(л:), которое задавалось при расчете контура сопла при отсутствии вязкости.
Рис. 4
Эксперимент и расчеты прямой задачи о течении в данном сопле с учетом вязкости по предлагаемой в работе методике отвечают параметрам газа в форкамере р0=80-Ю5 Па, 70= 1000 К, при этом число Re~5-10®. В отличие от предыдущего случая стенка сопла была не-охлаждаемой. При расчете распределение Tw(x) задавалось на основе статистических данных измерения температуры внутренней поверхности конических сопл, работающих в аналогичных условиях. Отметим, что на большей части длины температура стенки неохлаждаемых сопел составляла 300—350 К.
Приведенное на рис. 4 сопоставление результатов расчета с экспериментом в целом свидетельствует о близком их соответствии. Некоторое отличие имеет место лишь в размерах невязкого ядра потока на срезе сопла.
В качестве примера на рис. 5 приведены профили газодинамических величин и, V, h, р и величины ví/vm по г), отвечающие выходному сечению сопла.
В заключение отметим, что изложенная в работе методика может быть использована для диагностики полей чисел М на выходе из сопл гиперзвуковых аэродинамических труб.
ЛИТЕРАТУРА
1. Саху Д., Дэнберг Д. Э. Расчеты трансзвуковых течений на основе уравнений Навье—Стокса и дифференциальной двухпараметрической модели турбулентности. — Аэрокосмическая техника, 1987, № 7.
2. Ш е н г Д., Ш е р р С. Д. Расчет обтекания полной компоновки космического самолета на основе уравнений Навье—Стокса. — Аэрокосмическая техника, 1987, № 7.
3. Schiff L. В., Steger J. L. Numerical simulation of steady supersonic viscons flow.—AIAA Paper, 1979, N 79—0130.
4. Быркин А. П., Толстых А. И. Компактные схемы третьего и четвертого порядков в задачах о внутренних течениях вязкого и невязкого газа. — Ж. вычисл. матем. и матем. физ., 1988, т. 28, № 8.
5. Быркин А. П., Тимофеева Т. А., Толстых А. И. Применение компактных схем третьего-четвертого порядка для расчета течения газа в соплах в большими сверхзвуковыми числами М на основе упрощенных уравнений Навье—Стокса. — Ученые записки ЦАГИ, 1988, т. 19, № 6.
6. Толстых А. И. О неявных схемах повышенной точности для систем уравнений. — Журнал вычисл. матем. и матем. физ., 1981, т. 21, № 2.
7. Теория турбулентных струй. Изд. 2-е. /Под редакцией Г. Н. Абрамовича. — М.: Наука, 1984.
8. Schetz J. A. Foundations of boundari layer theory.-—Prentice —
Hall, Inc. New, Jersey, 1984.
9. Безменов В. Я,- Козырев А. П., Маневич О. М. Гипер-звуковая гелиевая аэродинамическая труба. — Ш-я всесоюзная школа по методам агрофизических исследований. — Сб. докладов, т. 2, Новосибирск,
1982.
10. Тимофеева Т. А., Чистов Ю. И. Определение толщины вытеснения турбулентного пограничного слоя в осесимметричных гипёр-звуковых соплах при постоянной температуре стенки. — Труды ЦАГИ,
1972, вып. 1403.
Рукопись поступила 25/IV 1988 г.