______ученые записки а а г и
Том XXIV 1993
№ 1
УДК 532.525.011.5 533.6.071.4
К РАСЧЕТУ ПРОФИЛИРОВАННЫХ ГИПЕРЗВУКОВЫХ СОПЛ АЭРОДИНАМИЧЕСКИХ ТРУБ С УЧЕТОМ ВЛИЯНИЯ ВЯЗКОСТИ
А. П. Быркин, В. П. Верховский
Предложена методика решения обратной задачи о построении с учетом эффектов вязкости контура профилированного гиперзвукового сопла, обеспечивающего в невязком ядре на выходе однородное течение газа с заданным числом М. Методика основана на использовании итерационного процесса, состоящего в последовательном решении прямой задачи для вязкого и обратной задачи для невязкого течений. Приведен пример конкретного расчета в случае ламинарного режима течения. В рассмотренном примере для окончания итерационной процедуры потребовалось всего две итерации.
1. Метод учета влияния вязкости. Рассматривается обратная задача расчета с учетом влияния вязкости контура уш(х) профилированного сверхзвукового сопла (см. рис. 1), реализующего однородное течение газа с заданным числом М в невязком изоэнтропическом ядре потока.
При условии, что толщина пограничного слоя б мала по сравнению с радиусом (/м(6<С(/и>), задача о профилировании сопла сводится к нахождению поправки к исходному невязкому контуру у(х), равной толщине вытеснения ламинарного или турбулентного пограничного слоя б*, (в зависимости От режима течения). Контур у(х) рассчитывается методом характеристик и обеспечивает заданное распределение числа М вдоль оси сопла—Мо(х) и однородный поток газа на его выходе с потребным числом М. В случае гиперзвукового сопла, ко1*да толщина пограничного слоя становится соизмеримой с радиусом (б~//ш), необходимо учитывать эффекты поперечной кривизны, скольжения и скачка температуры газа на стенке. Решению данной задачи для ламинарного течения посвящены работы [1, 2].
В настоящей работе решение обратной задачи с учетом влияния вязкости, так же как и в [1, 2], строится исходя из требования реализации заданного распределения числа М вдоль оси невязкого ядра в сопле с искомым контуром ут (х), которое совпадает с распределением Мо(х), отвечающим течению невязкого газа в сопле с контуром у(х). Указанное условие означает совпадение характеристик обоих названных течений в пределах невязкого ядра потока.
Отметим, что предложенный метод применим и в случае турбулентного течения в пограничном слое. При этом предполагается использование той или иной модели турбулентности с соответствующими эмпирическими константами, адекватно описывающей рассматриваемый класс течений.
Искомый контур уш(х) определяется в процессе последовательных приближений следующим образом. Так, приняв в первом приближении для функции А(х)=уи>(х) — —У(х) распределение по длине толщины вытеснения пограничного слоя на стенке сопла 6*(х) (на основе численных расчетов в случае ламинарного и экспериментальных корреляционных зависимостей в случае турбулентного режимов течения), получим приближенное решение рассматриваемой задачи:
(*) = у (*)+ Д0> (х) = у(х) + Ь* (х).
Решая далее с использованием упрощенных уравнений Навье—Стокса прямую задачу о течении вязкого газа в сопле с найденным контуром у^ {х)при заданных исходных предположениях (однородное течение газа и отсутствие пограничного слоя в критическом сечении сопла, заданная температура стенки и др.), определим распределение числа М вдоль оси (х) вплоть до правой границы характеристического
ромба. Для полученного распределения (х), решая затем обратную задачу для
течения невязкого газа, найдем соответствующий ему невязкий контур сопла уО )(.£).
При этом примем во втором приближении
где
„(2)
цу
Д<2>
Уш(х)=* У (•*).+ д(2) М.
У(1,(Х).
По мере необходимости процесс уточнения распределенной поправки ^п\х) повторяется. Критерием его окончания должно быть совпадение с приемлемой точностью заданного М0(х) и полученного У\^п> (х) распределений.
Предложенная, методика апробирована на примере расчета контура сопла на число М=18 для ламинарного режима течения совершенного газа (х=1,4), результаты которого приводятся ниже.
2. Алгоритмы расчетов течений вязкого и невязкого газа. При реализации изложенного выше метода последовательных приближений численное решение прямой задачи о течении вязкого газа в сопле с заданным контуром у^\х) проводилось по методу [3].
Использовалась упрощенная система уравнений Навье—Стокса для случая течения в осесимметричном сопле. Переходя от цилиндрической системы координат х, у к системе координат %,—х, т] = у/уы(х) и сохраняя составляющие скорости и и и вдоль осей х я у, упрощенные уравнения Навье—Стокса получим из полных уравнений при отбрасывании в них вязких членов со вторыми производными по | и смешанными производными по £ и г). Полученные таким образом уравнения содержат полностью все члены уравнений Эйлера и вязкие члены с производными по 11 и имеют вид
дР дй
+ -Щ=Н+ 1ІГ7''
(1>
где /\ б, Н, Т — векторы-столбцы
УУ<с
ра
РИ2 -1- ГУ р« (Л +
И3 + Vі
а = у
Р (« - У^а)
Р (и - Ут г1и) и~У т Р (и — Ут Ї1М) V + Р
Р
0 Ті
0 — Т2
р т = Тз
о
н =У«
В уравнении (1) все величины беразмерные; р — плотность, р — давление, Л — энтальпия, — динамическая вязкость, Кеі = р1мі(/г0 і/(Хі (индекс 1 отвечает значениям размерных величин на оси сопла в критическом сечении, индекс ни — на стенке); Уно—производная контура ущ(х); Ті = 0, составляющие Т2, Т3, Г4 учитывают вязкие члены с производными по г|. Газ предполагается совершенным, т. е. уравнение состояния имеет вид
Р Л.
и счцтэется заданной зависимости ^(Г),
Граничные условия для системы (1) следующие:
при £ — О,
р = р(т;), u = u(fi), v — v^), h = h (ї]) д[і ди dh
d-ц д'ч df\
и = v = 0, h = h (і)
= 0, t» = О при ■>1=0, при Т) = 1.
Для численного решения системы (1) с граничными условиями (2) применяется метод глобальных итераций, когда на каждой итерации используется маршевый алгоритм, а аппроксимация продольной составляющей градиента давления строится с учетом найденных на предыдущей итерации величин давления вниз по потоку от рассматриваемого сечения | = const. Применение такого подхода сводится к тому, что при М*<1+е (М*— число М в расчетной точке, определенное по продольной составляющей скорости, е>0) член д/д%{ууюр) в уравнении импульсов в направлении ответственный за передачу возмущений вверх по потоку, представляется в разностном виде следующим образом:
д. . 1я (УУ«Р)і+1]~(УУ^Р)1.і
til {УУ™Р) I/. / =--------Д£..........•
где 1, / — номера узловых точек расчетной области в направлении | и т], п — номер текущей глобальной итерации.
В схеме глобальных итераций необходимо, кроме того, задание граничного условия для давления на границе расчетной области.
При дискретизации уравнения (1) предварительно осуществлялось сгущение узлов в пристеночной области, достигавшееся преобразованием переменной т) и введением разностной сетки с постоянным шагом.
Отметим, что при численном решении задачи о течении в соплах при числах Яе~(5—10)-10е необходимо учитывать эффекты турбулентного переноса (число Ие определено по параметрам газа в невязком ядре потока на срезе сопла и его длине). Опыт расчета таких течений в соплах (прямая задача, см. [4]) при использовании однопараметрической модели турбулентности [5] (с учетом коррекции эмпирических констант для рассматриваемого случая) дает основание надеяться, что, по крайней мере, суммарное влияние турбулентного, пограничного слоя на формирование течения в невязком ядре сопл учитывается с достаточной точностью. При этом на основе модели [5] определяется турбулентная вязкость эффективная вязкость ц принимается равной сумме молекулярной и турбулентной составляющих (|х = цт + [х(), Ргт = Рг, = Рг.
При реализации метода последовательных приближений для нахождения искомого контура сопла </»(*) на каждом шаге требуется также решать обратную задачу для течения невязкого газа. Решение указанной обратной задачи находилось методом характеристик [6].
3. Пример расчета. По предложенной в п. 1 методике проведен расчет контура профилированного сопла с учетом влияния вязкости ут(х), обеспечивающий на выходе однородный поток газа с числом М=18. Невязкий контур сопла у(х) в конце области разгона имеет конический участок, на выходе которого реализуется радиальное течение газа. Полуугол раствора конического участка 0=15°. Рабочий газ — воздух, параметры торможения р0=Ы07 Па, Го = 2000 К, радиус критического сечения у к 1 = 1 мм, температура стенки Тт = 377 К. Для указанных условий значение Ие1 = = 1,56-105, при этом значение Не, определенное по параметрам газа на выходе и длине сопла, равно — 1,9 -106. Течение предполагалось ламинарным, причем расчеты невязкого и вязкого течений проводились как для совершенного газа (и = 1,4), с тем чтобы исключить влияние эффектов неравновесности и основное внимание уделить влиянию эффектов вязкости. Коэффициент вязкости определялся по формуле Сазерленда (постоянная С = 110 К), число Рг = 0,71. Влияние скольжения и температурного скачка на стенке для рассмотренных условий оказалось пренебрежимо малым.
На рис. 1 показаны невязкий контур сопла у(х), распределения чисел М на оси — Мо(л-) и невязком контуре—М(х), характеристический ромб.
В первом приближении в качестве поправки Д(х) к невязкому контуру использовалось распределение толщины вытеснения пограничного слоя б* по х, определенной расчетным путем по методу [7] с применением асимптотических уравнении Прандтля для заданных условий. Найденная поправка к невязкому контуру у(х) в в виде где Щв (л) = с1у№х, представлена на рчс 1,
J74
Рис. 1
Рис. 2
После решения прямой задачи о течении вязкого газа в сопле с контуром У®* (■*) = У (■*) + 8* (*) получено распределение М^дс), показанное на рисунке.
Видно, что найденная уже в первом приближении зависимость M^IJ (лг) оказалась весьма близкой к потребной зависимости М0(*), хотя значение б* на выходе сопла примерно равно радиусу невязкого ядра. Данный факт является примечательным.
По данным расчетов вязкого течения в первом приближении на рис, 1, 2 приводятся поперечные профили чисел М в различных сечениях сопла и профили газодинамических величин ц, v, h, р в окрестности вмаддиргд сечения,
Уточненный контур Ую\х) 00 втором приближении обеспечивает в пределах точности расчетов практически однородное течение в области характеристического ромба.
Представленные на рис. 1 данные, в частности, показывают, что радиус невязкого ядра потока в сопле уя в сечении .*=1500 и *=1980 имеет примерно одно и то же значение ~75, при этом отношение уя/ут соответственно равно —0,45 и 0,38. Полученный результат означает существование оптимальной длины сопла, отвечающей максимальному значению у„/ую. Последнее находится в соответствии с идеей использования оптимальных в указанном смысле (или укороченных) сопл (см., например, [8]).
ЛИТЕРАТУРА
1. Михайлов В. В. Метод расчета сверхзвуковых сопл с учетом влияния вязкости //Изв. АН СССР. МЖГ. — 1969, № 1.
2. Денисенко О. В. Метод расчета сверхзвуковых сопл при сильном влиянии вязкости//Ученые записки ЦАГИ.— 1992. Т. 13, № 4.
3. Б ы р к и н А. П., Тимофеева Т. А., Толстых А. И. Применение компактных схем третьего—четвертого порядка для расчета течения газа в соплах с большими сверхзвуковыми числами М на основе упрощенных уравнений Навье—Стокса//Ученые записки ЦАГИ.— 1988. Т. 19, № 6.
4. Безменов В Я., Б ы р к и и А. П., Г о р е н б у х Г1. И., Сабельников В. А., Тимофеева Т. А., Толстых А. И. Исследование течения газа в гиперзвуковых соплах при больших числах Рейнольдса на основе упрощенных уравнений Навье—Стокса//Ученые записки ЦАГИ. — 1989. Т. 20, № 4.
5. Теория турбулентных струй/Под ред. Г. Н. Абрамовича. — М.: Наука, 1984.
6. Кацкова О. Н., Наумова И. Н., Шмыглевский 10. Д., Ш у л и ш н и н а М. П. Опыт расчета плоских и осесимметричных течений методом характеристик. — М.: ВЦ АН СССР, 1961.
7. Б ы р к и н А. П., Щ е н н и к о в В. В. Об одном численном методе расчета ламинарного пограничного слоя//Ж. вычисл. матем. и матем. физ. — 1970. Т. 10, № 1.
8. М е ж и р о в И. И. Исследование течений в гиперзвуковых соплах аэродинамических труб//Труды ЦАГИ.-—1981. Вып. 2119.
Рукопись поступила 15/1У 1991 г.