CC BY
57
УЧЕНЫЕ ЗАПИСКИ Ц А Г И Том XVI 1985
№ 1
УДК 532.525.2: 532.517
ЧИСЛЕННЫЙ РАСЧЕТ НЕИЗОБАРИЧЕСКИХ СВЕРХЗВУКОВЫХ ВЯЗКИХ СТРУЙ, ИСТЕКАЮЩИХ В СПУТНЫЙ СВЕРХЗВУКОВОЙ поток
Н. Ф. Борисов
С помощью укороченных уравнений Навье—Стокса решается задача об истечении сверхзвуковых недорасширенных ламинарных и турбулентных двумерных струй в спутный сверхзвуковой поток. Для численного расчета поля течения используется разработанная явная разностная схема 2-го порядка аппроксимации. Анализируется влияние вязкости на основные характерные участки поля течения неизобарических спутных струй.
На основе модельной задачи о расчете слоя смешения между двумя однородными полубесконечными потоками определяется величина погрешности используемой системы уравнений для расчета сильно недорасширенных струй.
В работе рассматривается задача об истечении сверхзвуковой недорасширенной ламинарной или турбулентной осесимметричной струи в спутный сверхзвуковой поток.
Ранее, в работе [1] было проведено численное исследование двумерной ламинарной спутной струи на основе решения полной системы уравнений Навье—Стокса. Данный подход имеет тот недостаток, что требует для расчета достаточно протяженного участка сильно недорасширенной струи с требуемой для практики точностью на существующих ЭВМ весьма значительного времени счета. В случае, когда число Рейнольдса Ие^>1, достаточную точность обеспечивают укороченные уравнения Навье—Стокса, в которых сохранены все невязкие члены, а в вязких слагаемых не учитываются члены, содержащие производные вдоль оси струи. Указанный подход позволяет перейти от краевой задачи к задаче Коши, что сокращает время счета на ЭВМ на два-три порядка по сравнению с временем, необходимым для расчета полных уравнений Навье—Стокса. Данный метод был использован в [2, 3] для расчета двумерных ламинарных струй, а в [4, 5] — для расчета осесимметричных турбулентных струйных течений.
Следует отметить, что в [5] для расчета струи использовалась разностная схема, имеющая порядок аппроксимации 0(т+/г) (т, к — интервалы расчетной сетки вдоль продольной и поперечной осей), а в [4] расчет проводился на основе численного метода, имеющего точность
О (а + /г2), и уравнения для осредненных параметров замыкались с помощью приближенной интерполяционной формулы для турбулентной вязкости. Кроме того, в |[4, 5] отсутствует исследование влияния начальных данных на структуру поля течения струи.
В настоящей работе для расчета вязких двухкомпонентных спут-ных струй с помощью укороченных уравнений Навье—Стокса используется разработанная явная разностная схема, имеющая порядок аппроксимации 0(т2+/г2). Исходная система уравнений замыкается с помощью дифференциальной модели турбулентной вязкости, предложенной в [6]. На основе разработанного численного алгоритма были проведены расчеты для следующих диапазонов начальных параметров струи и спутного потока: число Маха струи 2<МС<4, число Маха спутного потока 2<МП<10, нерасчетность 1<«<103, температурный фактор 0,1 < Гп/Гс < 1.
1. Для исследования структуры осесимметричной сверхзвуковой струи, истекающей в спутный сверхзвуковой поток, введем в меридиональной плоскости течения ортогональную систему координат г, г так, чтобы ось г совпадала с осью струи. Расчет вязкого течения будем проводить на основе укороченных уравнений Навье—Стокса, которые для двухкомпонентной смеси газов струи и спутного потока в безразмерной форме имеют вид:
д(ра) | д(ру) ду
дх
■1^
У
д (рас) , д(рг>с) дРх
дх ' ду ду
дх
д (риу) , д(р1/24-р)
дх 1 ду
д (риН) д(руН)
ду
_4_
3
ду
дРг | I I р [_ ~ду~ + ~ 1 3
Р, -
(Р1 — рис)-,
I 0)
дх
ду
Р.-Р®Я);
Н-
г 1/1 \ 1 I гг2 -4-1»2
■[ст(1 — с)д] Л------------±—;
(тг-1)М2с 1“ 1 '* ‘'''п ' 2 ’ Р
Р Т
1с Мс
где
р _|_ _|^э_ _дс_ Ш р_________________________________(*»_ ди_ . р _ &и_ _
1 1?е 5сэ ду ’ 2 Re ду ’ 3 Ие ду ’
1
Ие V У
■и_ _^э\ . р _
2 оу ) ’ 6 1?е
а V»
1 дН , 1 Л 1
[_ Ргэ ду + 2 I1 Ргэ
3 ду
_3_______1_
4 Ргэ
1
Бет
+
Ргэ РГТ ^ (хэ \ Рг
1
Ргт > Т„*п
7 — 1 1 /? ’
*п *с с
мс
1-<7
(Тс— 1)М? ' ЗСэ ыс
Ргэ
где и, и — компоненты вектора скорости вдоль осей г, г соответственно, р — плотность, Т — температура, р — давление, Я—полная удельная энтальпия, с—массовая концентрация газа, вытекающего из сопла,
е — коэффициент турбулентной вязкости, |д, — коэффициент молекулярной вязкости, [хт = ре, (хэ = |я + 1?е(хт; Рг и Ргт — ламинарное и турбулентное числа Прандтля, Эс и 5ст — ламинарное и турбулентное числа Шмидта, Я — газовая постоянная, М — число Маха, / = 0; 1—в плоском и осесимметричном течениях соответственно, у — показатель Р «с га
адиабаты, Ке — —число Рейнольдса, рс, ис, Тс — значения га-
Iх (•' с)
зодинамических параметров струи в начальной плоскости на оси симметрии и га — начальный радиус струи, х = г/г&, у = г/г&. Нижние индексы «с» и «п» обозначают величины, относящиеся к струе и спутному потоку соответственно.
В системе (1) при введении безразмерных величин компоненты вектора скорости отнесены к ис, плотность и температура — к рс и Тс, а давление и полная энтальпия — к рс и\ и и\. Область течения, на которой рассчитывалась система (1), была выбрана в виде полуполосы 0<х, 0<у<г/м- На оси у = 0 задавались условия осевой симметрии, т. е.
^ ди др дТ д
ду ду ду '
На верхней границе расчетной области в качестве граничных условий задавались параметры однородного спутного потока. Начальные данные на оси х = 0 выбирались следующим образом. В интервале 0<у<1 задавались начальные значения параметров струи в предположении, что течение в этой области можно рассматривать как течение от сферического источника. В области Ку<ум задавались начальные данные для спутного потока, который предполагался однородным и текущим ВДО'ЛЬ ОСИ X.
В исходную систему (1) из всех коэффициентов переноса входит только коэффициент эффективной вязкости цэ. Процессы теплопроводности и диффузии учитываются с помощью чисел Прандтля и Шмидта. В настоящем исследовании полагалось, что ламинарное и турбулентное числа Прандтля и Шмидта постоянны во всем поле течения. При этом выбор коэффициентов переноса сводится к определению величин |1 и е.
В расчетах функция ц(Т) задавалась или в виде степенной зависимости или определялась по формуле Сатерленда. Для определения е, как уже отмечалось, использовалась модель турбулентной вязкости, предложенная в [6], определяемая из уравнения:
д (рме) д (рие)
дх
ду
д
ду
ду
+ «2 Ре
ди
~ду~
+ аз £
<Х1 ре +
др
де
~ду
] +
и + V
дх ду
(2)
где а1 = 2, а2 = 0,2, а3 = 0,5.
Отметим также, что в последнее время опубликованы работы [7, 8], в которых используется модификация модели турбулентности 1[6], учитывающая влияние числа Маха.
Анализ системы (1), (2) показывает, что постановка для нее задачи Коши корректна только в случае, когда во всем течении струи М*>1, где Мх = и/а (а — скорость звука).
2. Использование приближения пограничного слоя для расчета неизобарических вязких струй на два-три порядка уменьшает время
2 «Ученые записки» № 1
17
счета на ЭВМ по сравнению с расчетом полных уравнении Навье— Стокса. Однако для правомерности использования приближения пограничного слоя система координат должна быть выбрана так, чтобы в области вязкого слоя смешения струи и спутного потока градиенты всех величин вдоль «поперечной» координатной линии у были бы много больше градиентов вдоль «продольной» координатной линии х. Только в этом случае в вязких членах можно пренебречь производными по оси л:. Поскольку при выводе исходной системы (1), (2) предполагалось, что направление координатных осей жестко фиксировано, то из вышесказанного следует, что используемая система уравнений пригодна, строго говоря, только для расчета таких струй, в которых слой смешения имеет достаточно малый угол с осью х. Практически это означает, что степень нерасчетности п = рс/ра не должна быть слишком велика.
Определим диапазон изменения характерных значений угла наклона слоя смешения г*)1 к оси симметрии в сильно недорасширенных сверхзвуковых спутных струях. Экспериментальные исследования и расчетные данные показывают, что в струях с п^> 1 максимальное значение угла наклона Фм слоя смешения реализуется в окрестности кромки сопла. При смещении вниз по потоку угол $ уменьшается и в области изобарического течения струи 0<с1. Оценим значение величины 0М в начальной области слоя смешения. В данной работе рассматривается модель струйного течения, в которой предполагается, что головной скачок уплотнения присоединен к кромке сопла. Поэтому величину дм в этой области можно приближенно определить из задачи о сверхзвуковом невязком обтекании тупого угла, равного я—#м- При этом величина определяется из условия существования присоединенного к основанию угла косого скачка уплотнения. Решение данной задачи известно [9]. Из него следует, что '&м^45°. Проанализируем некоторые особенности, которые могут возникнуть при использовании системы (1), (2) для расчета слоя смешения, который образует заметный угол
с осью х. Рассмотрим для простоты модельную задачу о расчете в декартовой системе координат х, у изобарического слоя смешения между двумя плоскими однородными несжимаемыми полубесконечными потоками (см. рис. 1). Будем считать, что векторы скоростей потоков и У2 параллельны между собой и образуют угол 0 с осью х. Укороченные уравнения Навье—Стокса в системе коордиат х, у имеют вид
„ди. ди , 1 др д2 и ди , dv п
и —ч—- -j- v — -4- ■----------------------v— — v ■■ . ; —к—■ -4— — = 0;
дх ' ду 1 р дх дуг дх 1 ду ’
dv dv . 1 dp 4 д2 к ( ^ ^
дх ду р ду 3 ду2 ’
где V — кинематический коэффициент вязкости.
Преобразуем уравнения (3) к системе координат х', у' (рис. 1), которая повернута на угол 0 относительно системы х, у и «естественно» связана с геометрией слоя смешения. Перейдем в уравнениях (3) к переменным х', у' и величинам и', и', которые являются проекциями вектора скорости на оси х', у'. Тогда, используя в системе х', у' приближение пограничного слоя, из (3) получим
ди’ , ду' _П- ди' , „'ди' , 1.в!п2й\„.,яО,в'
d?- + ar = 0; u'-^r+v'w^V+ Tsin29Jcos26-^r-
Второе уравнение движения в системе координат х', у' в приближении пограничного слоя эквивалентно уравнению р(у')= const.
В то же время, если записать уравнения Навье—Стокса в системе х', у' и упростить их в приближении пограничного слоя, получатся уравнения в виде
ди’ . dv’ „ , ди’ , ди' д2 и'
dJP- + W = 0; 11 + V W = (5)
Из рассмотрения системы (4) видно, что эффективная вязкость \’э, равная
v. = v/(6); /(б) = (l + -Slyi) COS2 6, (6)
при больших углах 0 может быть существенно меньше физической вязкости V.
Определим величину погрешности вычисления толщины слоя смешения при расчете течения с помощью системы (4), которая эквивалентна системе (5) с точностью до замены v3 на v. Обозначим через b и Ь' толщину слоя смешения в системах координат х, у и х', у' соответственно (см. рис. 1). Рассмотрим сначала ламинарный слой смешения. Известно [10], что в этом случае система (5) при условии v = =const имеет автомодельное решение, из которого следует, в частности, что
b' = k \/r^x'jV^ ; k = const.
Аналогично для системы (4) имеем
&'(e) = 6v4(0)*7K- (7)
Подставляя (6) в (7), получим для относительной ошибки вычисления Ь' из системы (4) следующее выражение:
Ь' (х', 0 ) — Ь’(х’ 6) , л l/1 г sin29 /о\
Ь’\х'. 0)--L=1-cos6]/ 1 + —з- • (8)
Из соотношения (8) видно, что погрешность а возрастает с увеличением угла 0. Рассмотрим далее случай турбулентного слоя смешения. Будем определять турбулентную вязкость при расчете по уравнениям (3) в системе х, у по формуле
v# т = %Ь («! — и2); * — const, (9)
а в случае расчета системы (5) — в виде
V, = *Ь'0(У1- V,).
Система уравнений (5) [а следовательно и уравнения (4)] с учетом (10) также имеет автомодельное решение [10], из которого следует, что
Полагаем, что в системе координат х', у' в области слоя смешения выполняются соотношения вида:
и' V'-, х' ^>уг.
Тогда имеем очевидные геометрические соотношения (см. рис. 1)
Подставляя (9), (10) в (11), получим следующее выражение для относительной ошибки вычисления величины Ь':
расчета в области начального слоя смешения струи при использовании уравнений (1), (2). Из проведенных ранее оценок следует, что при 0 = 45° максимальные значения а и ат равны а = 23,6%, ат = 41,7%.
Необходимо отметить, что ошибка при расчете спутных струй с яЗ> 1 будет меньше, чем в данной модельной задаче, поскольку на начальном участке струи заметное влияние на изменение толщины слоя смешения оказывает также продольный градиент давления. Кроме того, при уменьшении 0 значения величин а и ат уменьшаются весьма быстро; так, например, при 0 = 20° а = 4,2%, ат = 8,3|%.
3. В настоящей работе численный расчет уравнений (1), (2) проводится на неоднородной в плоскости х, у разностной сетке, которая выбирается следующим образом. Все поле течения условно разбивается на две области: область струи 0сг/<г/о и область спутного потока Уо<У<Уы■ В качестве границы у0{х) между этими областями выбрана верхняя граница слоя смешения. Сетка вдоль оси у задается таким образом, чтобы она имела достаточно малый интервал в области струи и в то же время была достаточно крупной в области спутного потока. Это позволяет при общем числе расчетных узлов Л^о — 100 по оси у обеспечить большую протяженность исследуемой области течения по оси у, уменьшить схемную вязкость в области слоя смешения, а также увеличить точность расчета системы скачков уплотнения в центральном ядре струи. Для удобства расчета введена новая переменная £ = £(*/) к уравнения рассчитываются на плоскости х, £ в области 0<х; 0<£<£м(г/м). Использование переменной £ позволяет видоизменять в процессе численного счета сетку на оси у при сохранении в плоскости х, 2; однородной по оси £ разностной сетки. Это упрощает численную схему и сокращает время счета на ЭВМ.
На начальной прямой л;=0 функциональная связь у(£) определяется следующим образом. В области струи 0<^<^о функция = у(£,)
имеет вид
х,?т- ; &! = соп51. (11)
у 1
м = и' сое 6; Ь' = Ь соБб.
(12)
(13)
Выражения (8), (13) позволяют оценить величину погрешности
где Ну, Н0 — постоянные интервалы между узлами сеток вдоль осей у
и £■
Из (14) видно, что в области 0<г/<г/0 сетка является однородной и содержит заданное количество И1 расчетных узлов. В области ьо<£<£м функция г/'(£) выбирается в виде
у = Ж2 + 5С + С, ум = у(Си). (15)
Стыковка зависимостей (14) и (15) осуществляется из условия непрерывности величин у(£) и йу/й^ при $ = ?о-
В процессе счета через каждые АЫ шагов по оси л: (ЛЛ^—10) определяется новое значение г/о и производится пересчет сетки на основе функциональных зависимостей (14), (15). Одновременно с этим параметры потока с помощью интерполяционных формул пересчитываются на новую сетку. В данной работе для расчета исходной системы уравнений (1), (2) выбрана явная разностная схема типа предиктор—корректор, имеющая второй порядок аппроксимации. Для модельного уравнения вида
ди , ду (и)
+ "
д
К
+ /
дх . ОС,
используемая разностная схема на расчетной сетке
2 <= {-*' =ут, у = 0, 1, ...; С, = г'А0, I = 0, 1,
может быть представлена следующим образом:
(16)
^о}
,/+Т
а) предиктор
+ и[)
0,5т
+
+ 0
Т(+Г
+ (!-<>)
?{-!
ЗЛ0
4+1 +л .
(17)
Ы1+1-
б) корректор
■ и1 + 4- [“г+2 + 4-2 - 4 (“/+1 + “/-!> + 6“(:1
I /+т /+-
3 - - ф
* 1—
/+т /+т
/. I2- f. 1
,+т
/+•
— а
+
./+т ./+т
“ г-т
+
(18)
где Ф(~Ф(х<, С,), т — интервал расчетной сетки по оси
Разностная схема (17), (18) аппроксимирует уравнение (16) с точностью до 0(х2+Ио) при любых значениях параметров б, а, оа. Варьирование значений этих величин позволяет менять характер схемной вязкости в зависимости от поведения решения в окрестности рассматриваемой точки. Отметим, что при 6 = 0, о=1, со = 0 схема (17), (18) совпадает с разностной схемой, предложенной в [11]. Анализ устой-
чивости системы (17), (18), проведенный при 6 = 0, сг=1 с помощью метода малых возмущений, показал, что область допустимых значений параметра со
0 <;—•(!—х2) при *-<1,
4
(19)
где X = -=—
Л0 <1и
Проведенные методические расчеты модельной задачи о взаимодействии двух плоских полубесконечных потоков идеального газа показали, что наименьшие осцилляции в области тангенциального разрыва и скачка уплотнения получаются при о = 9/8, 6=1/16, со = ОД.
4. На основе описанной выше численной методики составлена программа для ЭВМ и проведен ряд расчетов сверхзвуковых осесимметричных ламинарных и турбулентных спутных струй. Начальные данные течения изменяются в следующих диапазонах:
На рис. 2 и 3 в качестве примера приведены результаты расчета поля течения турбулентной струи для случая, когда начальные данные имеют следующий вид:
где <р0 — полуугол раствора сопла.
Параметры потока в области начального слоя смешения изменяются от значений, соответствующих течению в струе, к параметрам спутного потока по линейному закону. Толщина начального слоя смешения задается равной Аг/ = 0,1. Турбулентная вязкость на начальной прямой у = 0 определяется в соответствии с формулой Прандтля.
На рис. 2 приведены результаты, соответствующие неизобарическому участку струйного течения. Здесь кривая 1 соответствует головному скачку перед струей, кривая 7 — висячему скачку. Кривые 2—6
2<МС<4, 2<МП<10, 1</г<103, 0,1 <7’п/7’с< 1.
Тп/Тс = 0,1; « = 100; Мс = 4; Мп = 3;
Тс = Тп= 1.4; А/д = 250;. N, = 100; Ргт = 0,75; (20)
8ст = 0,625; ¥0 = 0.
0
11 16 ?
а)
б)
Рис. 3
являются линиями постоянной безразмерной избыточной температуры торможения
Тр(у) — Т0 (0)
Т'оп-Т'о(О)
Т*(у)
для значений Г* = 0,99; 0,8; 0,5; 0,2; 0,01 (кривые 2—6 соответственно). Здесь Та — температура торможения. Точками 8 на рис. 2 показана граничная линия тока, выходящая из кромки сопла. Кроме того, для сравнения в виде точек 1' и 8' показано положение скачков уплотнения и аналогичной линии тока для невязкого течения. Величины | и г], отложенные по обеим осям на рис. 2, выбраны в виде:
1 = х/|/п\ у] — у/Уп.Течение в струе в области первой «бочки», которая для данного варианта расчета расположена в интервале 0<|^10, вдоль оси у можно условно разбить на несколько участков, а именно: ядро струи — расположенное между осью симметрии и висячим скачком, сжатый слой — расположенный между висячим скачком и нижней границей слоя смешения (кривая 6), вязкий слой смешения, ограниченный кривыми 2 и 6, и наконец, область возмущенного невязкого спутного потока, которая находится между верхней границей слоя смешения (кривая 2) и головным скачком уплотнения.
Рассмотрим основные особенности поля течения и проанализируем характер влияния вязкости на неизобарический участок струи. Видно, что уже на длине первой «бочки» область вязкого слоя смешения может занимать значительную часть всего поперечного сечения струи. Как следует из рис. 2, наличие толстого слоя смешения может приводить к сильному вязко-невязкому взаимодействию в периферийной области струи и оттеснению висячего скачка к оси симметрии. Аналогичный эффект вытеснения головного скачка выражен значительно слабее. После отражения от оси симметрии висячий скачок пересекает струю и уходит в сторону спутного потока. При прохождении через отраженный скачок линии тока начинают удаляться от оси струи. Вследствие этого в области течения \ 55 8 давление в центральной области струи падает и формируется вторая «бочка», которая, однако, выражена значительно слабее, чем первая «бочка». В области второй «бочки» наличие вязкости наряду с инерционными силами оказывает существенное влияние практически на все поле течения струи. Взаимодействие отраженного скачка и линий тока приводит к тому, что их положение становится существенно отличным от случая невязкой струи. Например, удаление от оси струи граничной линии тока в вязком течении при £=16 в 1,4 раза больше, чем в случае невязкой струи.
На рис. 3,а приведены радиальные распределения величин = = Р/Рп, и, Т при | = 5. Здесь же, для сравнения, пунктирными линиями показаны радиальные профили аналогичных величин для невязкой струи (точки 1, 2, 3 служат для обозначения величин и, р%, Т соответственно) .
Из рис. 3, а видно, что используемая численная методика позволяет достаточно аккуратно рассчитывать скачки уплотнения и тангенциальные разрывы, т. е. области с большими градиентами газодинамических параметров.
На рис. 3,6 показаны радиальные распределения величин р.л., и, е, Т, р*=р/рп при |= 10. Данное сечение соответствует той области течения струи, где отраженный от оси висячий скачок проходит через слой смешения. Особый интерес представляет здесь приосевая область 0<г}<0,15, характеризуемая наличием максимума для Т и минимума для величин и и р*. Анализ результатов расчета показывает, что за точкой отражения висячего скачка от оси струи в приосевой зоне образуется узкая область нагретого газа, которая может сохраняться достаточно далеко вниз по потоку. Наличие указанной области согласуется с данными работы [12], в которой показано существование в зоне регулярного отражения висячего скачка приосевой области с большими градиентами параметров течения.
На рис. 3, б также приведены данные о распределении Т при £=10 для других расчетных сеток. Кривые 1 и 2 соответствуют значениям-*^, равным 50 и 150. Видно, что данные расчета, полученные на разных сетках, достаточно хорошо согласуются между собой за исключением небольшой области около оси симметрии. Это показывает, что данные счета достаточно быстро сходятся к точному решению при уменьшении интервала расчетной сетки.
Из рис. 3, а и б видно, что на неизобарическом участке течения струи радиальные профили практически всех параметров носят сложный немонотонный характер. По мере смещения вниз по потоку течение в струе выравнивается и радиальные распределения параметров в струе приобретают плавный монотонный характер.
На рис. 3, в показаны радиальные профили для величин р%, и, Т, е для сечения | = 35, которое для данного режима истечения струи
соответствует вышеуказанному участку течения. Здесь же для сравнения в виде кривой 1 приведена зависимость и (у) для расчета, в котором турбулентная вязкость задавалась по формуле Прандтля.
Интересно отметить, что различные параметры потока с разной скоростью релаксируют к «автомодельному» профилю. Для данного варианта расчета быстрее всего релаксирует величина и (при |«*15). Примерно одинаково (при |«=:22) происходит релаксация величин р и Т. Изменение радиального распределения Т с увеличением | происходит следующим образом. Сначала уменьшается приосевой максимум температуры, который окончательно исчезает при |=15. Затем основной максимум Т, расположенный при |=10 в точке т] = 0,5 (см. рис. 3, б), начинает смещаться к оси струи вплоть до слияния с ней. В последнюю очередь завершается релаксация турбулентной вязкости, которая происходит при £~30.
На рис. 4 показано влияние показателя адиабаты ус на параметры течения на начальном участке струи. Кривые 1, 2, 3 на рис. 4 показывают положение граничной линии тока и скачков уплотнения для
Рис. 4
значений ус=1,25; 1,4; 1,667 соответственно. Видно, что с увеличением ус уменьшается диаметр первой «бочки» и несколько увеличивается ее длина.
Время расчета данного варианта на ЭВМ БЭСМ-6 при Д^=100 до £ = 50 равно 4 часам.
ЛИТЕРАТУРА
1. Ковалев Б. Д., Мы шенков В. И. Расчет вязкой сверхзвуковой струи, истекающей в спутный поток. — Ученые записки ЦАГЙ,
1978, т. IX, № 3.
2. Бондарев Е. Н., Горина А. Н. Решение задачи о сверхзвуковой ламинарной нерасчетной струе в спутном потоке разностным методом. — Изв. АН СССР, МЖГ, 1968, № 4.
3. М ы ш е н к о в В. И. Расчет течения вязкой ламинарной сверхзвуковой струи в спутном потоке. — Журнал вычислительной математики и математической физики, 1969, вып. 19, № 2.
4. Б о н д а р е в Е. Н., Л и с и ч к о И. Д. Распространение недорас-ширенной турбулентной струи в спутном сверхзвуковом потоке. — Изв.
АН СССР, МЖГ, 1974, № 4.
5. К о п ч е н о в В. И. Метод численного решения задачи о распространении сверхзвуковой недорасширенной турбулентной струи в спутном сверхзвуковом потоке.—’Ученые записки ЦАГИ, 1980, т. XI, № 4.
6. С е к у н д о в А. Н. Применение дифференциального уравнения для турбулентной вязкости к анализу неавтомодельных течений. — Изв. АН СССР, МЖГ, 1971, № 5.
7. Мещеряков Е. А., Сабельников В. А. Горение водорода в сверхзвуковом турбулентном потоке в канале при спутной подаче горючего и окислителя.ФГВ, т. 17, 1981, № 2.
8. К о з л о в В. Е. Метод расчета слабонеизобарической сверхзвуковой турбулентной струи в дозвуковом спутном потоке. — В сб.: Сверхзвуковые газовые струи. — Новосибирск, Наука, 1983.
9. Абрамович Г. Н. Прикладная газовая динамика.—М.: Наука, 1976.
10. Вулис Л. А., Кашка ров В. П. Теория струй вязкой жидкости.— М., Наука, 1965.
11. Lax P., Wendroff'B. Difference schemes for hyperbolic equation with high order of accuracy. — Comm. Pure and Applied Math, 1964, vol. 17, N 3.
12. Мельников Д. А. Отражение скачков уплотнения от оси симметрии. — Изв. АН СССР, Механика и машиностроение, 1962, № 3.
Рукопись поступила 3)111 1983 Переработанный вариант поступил 15/VI 1984
