CC BY
65
Т о м IX
УЧЕНЫЕ ЗАПИСКИ Ц А Г И
197 8
М 2
УДК 532.525.2.011.5.032
РАСЧЕТ ВЯЗКОЙ СВЕРХЗВУКОВОЙ СТРУИ, ИСТЕКАЮЩЕЙ В ЗАТОПЛЕННОЕ ПРОСТРАНСТВО
Б. Д, Ковалев, В. И. Мышенков
Рассматривается истечение сверхзвуковой осесимметричной струи вязкого совершенного газа в затопленное пространство. Предлагается метод расчета, достаточно эффективный в широких диапазонах изменения характерных параметров, основанный на решении краевой задачи типа Дирихле для полных уравнений Навье — Стокса методом установления с помощью неявной разностной схемы расщепления, что позволяет рассматривать как сверхзвуковые, так и дозвуковые струи. Исследуется изменение газодинамических параметров в начальном участке струи и проводится сравнение результатов расчета с имеющимися расчетными данными для идеальной струи и данными экспериментальных исследований.
Задаче истечения недорасширенной струи в затопленное пространство посвящено много работ как экспериментальных, так и теоретических (см., например, [1—5] и приведенные в них библиографии). В большинстве опубликованных работ исследуется течение на начальном участке сверхзвуковой струи, имеющем характерную бочкообразную форму, ограниченную висячей и замыкающей центральной (диск Маха) ударными волнами. В результате проведенных исследований найдены определяющие безразмерные параметры задачи и получены зависимости характерных размеров струи от числа М на срезе сопла Ма, показателя адиабаты х, степени нерасчетности я, числа Рейнольдса Ие. Обнаружена автомо-дельность по п струйного течения при больших п.
Однако, как показывают эксперименты [1, 2, 5], вязкое перемешивание у границы струи приводит (особенно с уменьшением чисел Ие) к существенной перестройке профилей плотности, давления и температуры в сжатом слое газа между висячей ударной волной и границей струи, изменяет положение ударных волн и границы струи. Образование различных вихревых течений над струей может привести к изменению характерных геометрических размеров струи до 20% [2]. В этой связи для определения газоди-
намических параметров в струе возникает необходимость построения метода расчета истечения вязкой струи.
Предлагаемый в настоящей работе метод расчета истечения струи вязкого газа в затопленное пространство в сравнении с методами расчета струй идеального газа [3, 4] характеризуется тем, что здесь (помимо добавления характерных параметров Ке, Рг, Таоо) решается задача типа Дирихле, в которой учитываются все возможные влияния одних областей течения на другие, а поэтому такая постановка является более полной и строгой. Она, очевидно, тем предпочтительней, чем меньше число Ие. (Аналогичная постановка задачи использовалась в работе [6], где решение было получено в непосредственной окрестности сопла).
1. Рассмотрим двумерную (осесимметричную или плоскую) задачу истечения струи вязкого совершенного газа из сопла в затопленное пространство при условии, что параметры потока заданы на срезе сопла (плотность ра, продольный и поперечный компоненты скорости иа и внутренняя энергия еа) и в окружающем затопленном пространстве на бесконечности (рм, и^, v<x, есо).
Задачу будем решать в рамках уравнений Навье — Стокса ко-нечно-разностным методом установления, предполагая существование и единственность решения. В качестве характерных параметров примем параметры на срезе сопла, а в качестве характерного размера — радиус среза сопла га.
Систему уравнений, описывающую течение вязкого сжимаемого теплопроводного совершенного газа, запишем в безразмерных переменных, аналогично работе [7], в виде
где 2 = 2 9г — дифференциальный матричный оператор,
и
г=1 дх
О и
О
О
О
__4__д_
дх ЗИе р дх ^дх
О О
О О
_д___1__д_
дх Ие р дх ^ дх
О
О О
д
и -Г-
дх Р*еРгр дх г дх
а, =
V— ООО
ду
0 + о о
д_ ду
О 0 г»
ду
4 / д _д_ д 1 п
' 3 Ие р ( ду ^ ду у ^ ду I и
ООО ч*-—^— + » ^
ду ИеРгр \дугду'у^ду}
24 =
0 д р^зг 0 0
е д -- Р дх 0 0 (х- -1) д дх
0 0 0 0
0 -\)е- д дх 0 0
0 0 Р д ду 0
0 0 0 0
!) * а Р ду 0 0 -1) д ду
0 0 (*- \)е д Ту 0
/=(р, и, V, е);
"0 0 0 0
1 0 0 у, о
Ке р 0 <р2 <Рз 0
_0 <Р4 ?5 0
?1=
<Э 2 д д <о Г <Э 2 д , . дх ду
д д
<р4 = -
д■ V дх
Г-
1 д«
а
дх
2 (?>- д* ] + Р [(-|у)
2 <о дц
~ 3~ У ~ду '
Л Й 1 (о 4
дх ду I ус
4
Ш'
1 ди д
+
2 Л* ду 4
1 ду дх
д
+
У 3
(т ~ ;
^(^-¿77- («-и*)-
Здесь приняты следующие обозначения; t — время, х, у—^прямоугольная система координат, ц, — коэффициент вязкости, Рг — число
Прандтля, Ие =
Ра "а Га Ра
число Рейнольдса, о> — коэффициент,
равный нулю в случае плоского и 1 в случае осесимметричного течения. Давление газа р исключено из системы уравнений посредством уравнения состояния /> = (*—1)ре. Зависимость коэффициентов вязкости ц и теплопроводности X от температуры Т принята в виде [8]
_£_ _ А. _ /_М3/2 Та+Тя
Ра ^а ^ Та )
Т+ Т3
где 7^= 110 К при Та = 300К — постоянная Сазерленда. Будем искать решение стационарной задачи
(е+Л-_!!-*)/= о
при граничных условиях
//= 0,
(1)
(2)
где / —граничный оператор в области, лежащей в верхней полуплоскости ху и ограниченной снизу осью (плоскостью) симметрии, а сверху поверхностью ¿Л, расположенной от оси на достаточном удалении, чтобы там можно было поставить условия, эквивалентные условиям на бесконечности (фиг. 1). С левой стороны область интегрирования ограничивается срезом сопла и выше него до поверхности Ьн твердой стенкой, на которой задаются условия
прилипания потока и температура стенки Тт. Правая граница области счета находится от среза сопла на расстоянии нескольких (более 15) радиусов га, и на ней задаются гладкие условия сопряжения типа д//дх — 0, либо д2//дх2 — 0. На срезе сопла почти во всех приведенных ниже расчетах параметры течения принимались постоянными, равными /а\ ра = иа=1; г>а = 0; еа = 1/х(х— 1)Ма.
Начальные условия можно задавать довольно произвольно. В настоящих расчетах примыкающая к оси часть рассчитываемой области заполнялась параметрами, равными параметрам на срезе сопла, а оставшаяся — параметрами на бесконечности.
2. Для численного решения поставленной задачи дискретизи-руем ее, сопоставив системе дифференциальных уравнений (1) некоторую разностную систему. Для этого в прямоугольном цилиндре £>н = £>Х Н (прямоугольная область, где ищется решение, Я=[0, Т], Г = (? X Н — граница) введем разностную сетку
М N
■*т = 2 Л = ¿ = тХ£; 6 = 0, ..., К. Здесь т — ите-
т—0 п=0
рационный параметр, кх{т) и й2(га) — шаги сетки по х, у.
Введем сеточную функцию /ь(?ктп, иктп, юктп, ектп) и разностные
операторы 2/, Ан, Вн, 1Н, аппроксимирующие исходные дифференциальные матричные и граничный операторы с первым или вторым порядком точности. В настоящей работе использовалась схема, где конвективные члены уравнений аппроксимировались односторонними разностями с первым порядком точности, а все остальные члены — симметричными разностями со вторым порядком точности. Расчеты по схеме второго порядка дали совпадающие результаты. Наконец, поставим в соответствие дифференциальной задаче (1) —(2) разностную задачу: в цилиндре Он найти решение уравнений
Фиг. 1
с граничными условиями
/»д=о.
(4)
Полученную разностную краевую задачу (3) — (4) будем решать методом установления с помощью неявной схемы расщепления, предложенной в работе [7]:
^ + -у В* + 6*+1/4 = х + Ан + у В*)/*;
(Я +тдйз)|А+3/4=^+2/4; + х2л4) = |*+3/4; /Г1 = +
где £*+л4 (/=1, 2, 3, 4) — вспомогательный вектор.
Исключение дробных шагов дает схему универсального алгоритма, обладающую свойством полной аппроксимации:
7 = 1 У У=1
Устойчивость приведенной схемы исследовалась в работе [7? на модельных уравнениях. Показана ее абсолютная устойчивость в линейном случае, а для нелинейных уравнений ее устойчивость экспериментально показана для чисел Куранта К 10, где
К== М + ^М + 'тахС!^!, |Л,|)),
здесь с — скорость звука.
Таким образом, система разностных уравнений (3) расщепляется на четыре системы более простых уравнений (5), каждая из которых решается посредством скалярной прогонки. Один полный шаг по времени Ь состоит из четырех дробных шагов.
Решение считаем установившимся, если вектор продольной составляющей скорости и в некотором выбранном поперечном сечении удовлетворяет условию |д£//с^|<8, где 8 < Ю-3. Расчеты с различными значениями 8<10~3 и Ю-4 показали практическое совпадение конечных распределений р, и, V, е.
Для проверки правильности работы программы с ее помощью была решена задача о течении в следе за пластиной при М=0,3 и Не = 100 в постановке работы [9]. Полученные распределения параметров течения хорошо согласуются с данными работы [9].
3. Задача истечения осесимметричной струи вязкого совершенного газа в затопленное пространство была решена при следующих значениях характерных параметров: ^ = 1,4, Рг = 0,71, = = Та/Т<х=1 и 10. Числа М, Ие и степень нерасчетности варьировались в диапазонах 1<Ма<4, 20< Ие < 104(105), 0,5<я<103.
Область счета содержала более 1700 расчетных точек (42X41). В окрестности среза сопла шаги расчетной сетки и Л2 примерно на расстоянии одного-двух радиусов среза сопла были постоянны и равны Л,-=Дл:, А3= Ду(Дл; = га/12, Ду = га/'16). Далее с увеличением х и у шаги кх и Л2 увеличивались по арифметической
прогрессии с шагами Ах, Ду соответственно. Такое сжатие области счета оказалось вполне удовлетворительным для большинства расчетов. При больших числах М и степенях нерасчетности применялись и другие способы сжатия координат, обеспечивающие достаточное увеличение области счета.
Исследование сходимости результатов расчета в зависимости от шага сетки проводилось для одних и тех же размеров области решения на следующих сетках: 10X 10, 25X25, 30X30, 35X35,
40 X 40. Результаты расчетов (для Ма=1, Не=103, п = 10) показывают хорошую сходимость распределений давления и плотности
0,08
ОМ
- Ма=7 п Чэ ?и 712;п =78,5
#
о 2 Ч- 6 8 70 х Фиг. 3
почти во всей рассчитываемой области с увеличением числа точек сетки, начиная с сетки 25 X 25.
Влияние числа М и степени нерасчетности на параметры истекающей струи в настоящей работе исследовалось при Не = 103, л =10 для Ма = уаг и Ма = 4 для п = уаг. Поскольку основные закономерности .качественного и количественного влияния этих характеристик течения на параметры струи уже выявлены (см. [1—4]), мы не будем обсуждать подробно этот вопрос. Отметим лишь, что полученные в настоящей работе данные достаточно хорошо согласуются с имеющимися результатами. На фиг. 2 приведены, в частности, распределения давления в различных поперечных сечениях струи, полученные для вязкого потока при Ке = 103, Ма=1, п =10 (сплошные линии) в сравнении с данными расчета для течения идеального газа при Мо=1,05, »=10 из работы [3] (пунктирные линии). Некоторое расхождение в распределении параметров у среза сопла и границ струи естественно и объясняется влиянием вязкого взаимодействия.
Для сравнения с экспериментальными данными работы [5] были проведены расчеты течения струи с Мв=1, п ==18,5 для числа Ке = 112. Результаты расчетов (сплошная линия) и экспериментальные данные (звездочки) по распределению плотности, отнесенной к плотности газа в форкамере (плотности торможения р0), приведены на фиг. 3. Как видно из графика, соответствие данных вполне удовлетворительное. Небольшое расхождение результатов вблизи сопла, очевидно, объясняется некоторой неравномерностью экспериментального профиля скоростей на срезе сопла, тогда как в расчете он был равномерным. Расхождение за центральным скачком определяется влиянием правой границы рассчитываемой области течения.
Наиболее интересными представляются результаты расчетов с различными числами Ие, поскольку в большинстве имеющихся
экспериментальных работ исследуются лишь внешние проявления воздействия числа Йе (изменение внешних характерных геометрических размеров) и описывается качественная картина его влияния. Лишь в отдельных работах [2, 5] приведены некоторые количественные данные по распределению плотности и давления.
Влияние числа Ие исследовалось при «=10, Ма=1 и 4 в диапазоне чисел 20<Ие<105. Расчеты показали, что распределения параметров течения- вдоль оси струи при числах Ие>103 почти не отличаются друг от друга, что объясняется, по-видимому, действием только сеточной вязкости. С уменьшением числа 1?е внешняя граница струи расширяется. Висячий скачок уплотнения располагается ближе к оси струи, чем в идеальном случае, а изменение параметров в нем и в сжатом слое более плавное. О границах этих областей можно судить по распределению статического давления, причем довольно условно, поскольку при сквозном счете (без выделения поверхностей разрыва) градиенты давления в сжатом слое и висячем скачке оказываются соизмеримыми и трудно различимыми на используемых сетках. Определение внешней границы струи для вязкого течения также затруднительно в связи с образованием слоя смешения и вторичных течений в затопленном пространстве. Приблизительное положение ее, по-видимому, можно найти по распределению статической температуры. На фиг. 1 схематически приведены характерные геометрические размеры вязкой струи для Ма=1, Ие^Ю8, «=10 (сплошные линии) в сравнении с данными идеальной струи (пунктирные линии). Здесь цифрами обозначены: /—внутренняя граница области влияния вязкости, 2— висячая ударная волна, 3—отраженная ударная волна, 4—р — граница струи (линия, выше которой р = Роо), 5 —Г —граница струи (линия внутри слоя смешения, где Г=Ггаах), 6 — Т — граница слоя смешения (над которой Т = 7^), 7—форма ядра струи, полученная экспериментально [I], 8— форма первой бочки идеальной струи [3]'
Для Ма = 1, как видно из фиг. 1, используемая сетка позволила получить картину течения, включающую область первой бочки. Интересно отметить, что в этом случае даже при Ие = 103 не было получено нерегулярного отражения висячего скачка (диск Маха) с характерной для него дозвуковой зоной. Отражение получалось регулярным также и для малых степеней нерасчетности п<С1 при Мй = 4, Не = 103. Аналогичная структура течения наблюдалась и экспериментально, например, в работах [5, 10]. Далее за центральным скачком прослеживается четко выраженный протяженный кольцевой слой сжатого газа, на который уже обращалось внимание в работе [2]. С уменьшением числа Ие вся эта сложная картина течения сглаживается, скачок приближается к срезу сопла и при Ие=20 бочкообразная структура течения исчезает вообще.
Эволюция течения с уменьшением числа Ие для М0 = 4 существенно отличается от эволюции для Ма = 1. Так, например, для Ие = 50 при Ма = 4 на выходе из сопла повышается температура, давление и плотность (скорость уменьшается), чего не происходит приМа=1. Это изменение термодинамических функций обусловлено вязким взаимодействием задаваемого сильно завихренного слоя с основным потоком у среза сопла. Интенсивность и характер этого взаимодействия определяются степенью завихренности 2, толщиной вихревого слоя 8 в. основном потоке перед срезом сопла и величиной числа Ие.
С уменьшением числа Яе интенсивность взаимодействия возрастает. Так, задание равномерных профилей р, и, V, е с прилипа-
/ ^ ди
нием на кромке среза при этом 2 =
ду
16, 8 =
приводит
при Ие — 50 к сильному возрастанию (до 50%) указанных функций р, р, Т в небольшой окрестности кромки среза сопла, в то время
Фиг. 5
как у оси струи это возрастание не превышает 8%. При уменьшении завихренности, задаваемой на входе, это различие нивелируется, распределение параметров в поперечном сечении струи становится более равномерным. Для параболического профиля скоростей, задаваемого на входе (2=1,8, 8 = 1), указанное повышение параметров на начальном участке струи наблюдалось вплоть до Ие = 103. При Мв = 1 подобные эффекты не проявляются явно, поскольку возникающее в этом случае повышение р, р, Т компенсируется их снижением в волне разрежения, начинающейся сразу от среза сопла. При Мв = 4 и Ие = 50 воздействие волны разрежения начинает проявляться лишь на некотором расстоянии от среза сопла. Далее по потоку вдоль оси струи давление монотонно падает, в то время как температура после небольшого минимума, соответствующего начальному участку волны разрежения, сильно возрастает, достигая максимального для слоя смешения значения, равного 1,5 Та и остающегося практически неизменным по всей длине зоны смешения (фиг. 4).
Высокая чувствительность температуры к изменениям числа Яе была характерной для всех расчетов. Это позволило по ее отклонениям сравнительно легко следить за перестройкой картины течения в струе. При Ма = 4 для больших чисел Ие статическая температура имеет перед висячей ударной волной значительный минимум и лишь за скачком уплотнения повышается до уровня Т^. При снижении числа Ие этот минимум уменьшается по абсолютной величине и при Ие-<100 исчезает совсем. Бочка заметно сужается, вязкий слой быстро растет, что выражается, в первую очередь, в значительном росте температуры по обе стороны от границы струи. При Ма=1 утолщение слоя смещения по температуре носит иной характер. Из фиг. 5 видно, что этот слой имеет тенденцию разрастаться с уменьшением числа Ие в направлении к оси струи, а не в сторону затопленного пространства. В целом при Ма=1 изменение температуры в поперечных: сечениях струи имеет более плавный характер.
Расчеты показали, что над струей возникает инициируемое ею слабое движение окружающей среды, меняющееся с изменением числа Ие и имеющее вихреобразную (циркуляционную) структуру.
Исследование влияния температуры струи при п =10, Ма=1, Ие = 103 показало, что изменение температуры струи может привести к существенным изменениям газодинамических параметров вязкой струи и перестройке ее геометрии. Так, например, при увеличении относительной температуры Та х с 1 до 10 при п= 10 длина бочки примерно на 25% увеличивается. В поперечном сечении струи, наиболее сильным деформациям подвергается сжатый слой, в то время как в волне разрежения характер течения в обоих случаях одинаков. Такое воздействие температуры струи обусловлено уменьшением сил вязкости в слое смешения у внешней границы струи.
Задание неравномерного профиля скоростей на срезе сопла в виде параболы (Ма=1, Ке=103, я =10), наоборот, сокращает длину бочки примерно на 20% и приводит к такому распределению газодинамических параметров, которое эквивалентно распределению при истечении струи с равномерным профилем скоростей, но меньшего (примерно на 20%) радиуса.
Отсутствие в большинстве расчетов сколько-нибудь значительных поперечных возмущений параметров затопленного пространства в окрестности границы струи свидетельствует, что верхняя граница области счета /,Л, имитирующая бесконечность, выбрана правильно, т. е. достаточно удалена от оси течения и не оказывает заметного влияния на положение границы струи. Лишь при степенях нерасчетности, больших 50, вертикальный размер рассчитываемой области оказался недостаточно большим, что вызвало деформацию границы струи и висячего скачка уплотнения, однако на распределении параметров потока у оси струи это не отразилось.
Распределение давления на стенке, расположенной в сечении среза сопла („донной" стенке), как показывают расчеты, зависит от значений чисел М, Ие и нерасчетности. С уменьшением числа Ие у кромки сопла образуется большая область пониженного давления, достигающая, например, при Ие = 50 величины 0,2 га и га при Ма = 1 и 4 соответственно. В качественном отношении это явление аналогично смещению точки отрыва потока вниз по „донной" стенке при обтекании ступеньки [9]. Далее вверх по стенке распределение давления полностью определяется вторичным течением, инициируемым струей. С увеличением степени нерасчетности: при Ма = 4 эти эффекты несколько усиливаются, а с уменьшением при «< 1 исчезают.
ЛИТЕРАТУРА
1. АвдуевскийВ. С., Иванов А. В., Карпман И. М., Трасковский В. Д., ЮделовичМ. Я. Течение в сверхзвуковой вязкой недорасширенной струе. „Изв. АН СССР, МЖГ", 1970. № 3.
2. Волчков В. В., И в а н о в А. В., К и с л я к о в Н. И., Р е б-ров А. К., Сухнев В. А., Шарафутдинов Р. Г. Струи низкой плотности за звуковым Соплом при больших перепадах давления. ПМТФ, 1973, № 2.
2—Ученые записки № 2
17
3. А в е р е н к о в а Г. И., А ш р а т о в Э. А., В о л к о н с к а я Т. Г., Дьяконов Ю. Н., Егорова М. И., Мельников Д. А., Росляков Г. С., Усков В. И. Сверхзвуковые струи идеального газа, ч. 1, 2, М., Изд-во МГУ, 1970.
4. Жохов В. А., Хомутский А. А. Атлас сверхзвуковых течений свободнорасширяющегося идеального газа, истекающего из осесимметричного сопла. Труды ЦАГИ, вып. 1224, 1970.
5. Шарафутдинов Р, Г. Исследование структуры недорас-ширенных струй низкой плотности с помощью пучка электронов. Дисс. канд. физ.-мат. наук. Новосибирск, Ин-т теплофиз. СО АН СССР, 1969.
6. Кузнецова Л. В., Павлов Б. М. О расчете струйных течений вязкого газа. В сб. „Вычислительные методы и программирование", вып. 23. М., Изд-во МГУ, 1974.
7. Бе резин Ю. А., Ковеня В. М., Яненко Н. Н. Об одной неявной схеме течения вязкого теплопроводного газа. В сб. .Численные методы механики сплошной среды", т. 3, № 4, 1972.
8. Голубев И. Ф. Вязкость газов и газовых смесей. М., Физ-матгиз, 1959.
9. Мышенков В. И. Численное исследование течения вязкого газа в следе плоского тела. Журн. вычисл. матем. и матем. физ., т. 12, № 3, 1972.
10. Driftmyer R. Т. A correlation of freejet data. .AIAA J"., vol. 10, N 8, 1972.
Рукопись поступила 3fV 1977
