Научная статья на тему 'Численная схема монотонного типа для параболизованных уравнений Навье — Стокса'

Численная схема монотонного типа для параболизованных уравнений Навье — Стокса Текст научной статьи по специальности «Физика»

CC BY
335
79
i Надоели баннеры? Вы всегда можете отключить рекламу.
Область наук
Ключевые слова
TVD-СХЕМА / ЗАДАЧА РИМАНА / СОБСТВЕННЫЕ ЗНАЧЕНИЯ / КОНЕЧНО-РАЗНОСТНАЯ СХЕМА / ПАРАБОЛИЗОВАННЫЕ УРАВНЕНИЯ НАВЬЕ — СТОКСА / TVD SCHEME / RIEMANN PROBLEM / EIGENVALUES / FINITE-DIFFERENCE SCHEME / PARABOLIZED OF THE NAVIER-STOKES EQUATIONS

Аннотация научной статьи по физике, автор научной работы — Макашева Алтын Прмагамбетовна, Найманова Алтыншаш Жамаевна

Построена разностная схема с TVD-аппроксимацией для решения трёхмерных параболизованных уравнений Навье — Стокса. В качестве базовой используется схема с разностями против потока первого порядка точности (Лакса — Фридрихса). Приведены матрицы левых и правых собственных векторов, позволяющие проводить расчёты единым образом в сверх- и дозвуковых областях течений

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Numerical scheme of a monotonous type for the parabolized form of the Navier-Stokes equations

A scheme with TVD approximation for solving the three-dimensional parabolized Navier-Stokes equations is constructed. The scheme relies on the upstream differencing scheme of the first order of accuracy (Laxe-Fridrix). The matrixes of the left and right eigenvalues are given to enable calculations in both superand subsonic region flows in a uniform procedure.

Текст научной работы на тему «Численная схема монотонного типа для параболизованных уравнений Навье — Стокса»

Вычислительные технологии

Том 17, № 5, 2012

Численная схема монотонного типа для параболизованных уравнений Навье — Стокса

А. П. Макашева, А. Ж. Найманова Институт математики Министерства образования и науки, Алматы, Казахстан

e-mail: ked@math.kz

Построена разностная схема с ТУБ-аппроксимацией для решения трёхмерных параболизованных уравнений Навье — Стокса. В качестве базовой используется схема с разностями против потока первого порядка точности (Лакса — Фридрихса). Приведены матрицы левых и правых собственных векторов, позволяющие проводить расчёты единым образом в сверх- и дозвуковых областях течений.

Ключевые слова: ТУБ-схема, задача Римана, собственные значения, конечно-разностная схема, параболизованные уравнения Навье — Стокса.

Введение

Истечение сверхзвуковых струй в спутный поток представляет существенный интерес для областей, связанных с проектированием сопел прямоточного воздушно-реактивного двигателя (ПВРД), газодинамического регулирования вектора тяги маршевых двигателей, а также для изучения теплового воздействия выхлопных струй на хвостовое оперение или фюзеляж самолета. Для таких струйных течений характерны взаимодействия сверхзвуковой зоны течения с дозвуковыми областями, наличие скачков уплотнений и слоев смешения. Известно, что при истечении системы сверхзвуковых струй в спутный сверхзвуковой (дозвуковой) поток в осевом направлении диффузионный перенос потоков массы, количества движения и тепла пренебрежимо мал по сравнению с конвективным переносом, а поле давления вниз по потоку оказывает незначительное влияние на распределения параметров вверх по течению. Наличие этого преимущественного направления развития течения позволяет для исследования струйных течений использовать параболизованные, или упрощенные, уравнения Навье — Стокса, которые по сложности занимают промежуточное положение между полными уравнениями Навье — Стокса и уравнениями пограничного слоя. Многочисленные расчёты [1-3] показали, что приближённая модель на основе параболизованных уравнений Навье — Стокса с достаточной степенью точности воспроизводит параметры потока и требует меньших вычислительных ресурсов. Несмотря на то что в последние десятилетия активно развиваются математические модели на основе полных уравнений Навье — Стокса, модели на основе параболизованных уравнений могут быть использованы для качественного исследования течений с преобладающим направлением. Основными трудностями численного решения задач истечения системы сверхзвуковых струй в спутный сверхзвуковой (дозвуковой) поток с привлечением как параболизованных, так и полных уравнений Навье — Стокса являются наличие осцилляций в решениях и разрывы в рассчитываемых функциях. Эти особенности требуют принятия определённых мер,

z

О

о е о о

о

о е о о

ООО

р. о о

к_

Рис. 1. Схема течения

компенсирующих нежелательные эффекты. При решении в разностные уравнения вводятся дополнительные диссипативные слагаемые, которые видоизменяют разностные операторы в окрестности разрыва. Так, например, в работе [4] для подавления осцил-ляций вблизи сильных ударных волн вводилась искусственная вязкость, представленная членом четвёртого порядка. Следовательно, применение схем повышенного порядка приводит к введению искусственной вязкости, которая может существенно влиять на точность получаемого решения. Квазимонотонные схемы, свободные от указанного выше недостатка, были предложены в [5, 6] и получили название TVD (Total Variation БтшвЫ^)-схем. Основная идея построения TVD-схем с высокой разрешающей способностью на базе TVD-схемы первого порядка точности заключается в модификации векторов потоков, осуществляемой так, что получаемые схемы имеют второй порядок точности в областях плавных изменений решения и первый порядок точности вблизи разрывов решения. На основе TVD-схемы исследованы взаимодействия скачка уплотнения с пограничным слоем при ламинарном течении вязкого сжимаемого газа и обтекание профиля трансзвуковыми и сверхзвуковыми потоками. Численные эксперименты показывают, что расчёты по TVD-схеме хорошо согласуются с экспериментальными данными и демонстрируют наибольшую точность в окрестности ударной волны. Построение TVD-схемы для двумерных и трёхмерных параболизованных уравнений Навье — Стокса приведены в работах [7, 8], в которых рассматриваются задачи обтекания пластины в случае совершенного газа. Однако в них в пограничном слое, где появляются дозвуковые области, применяется схема с направленными разностями первого порядка точности. В связи с этим основная цель настоящей работы — применение TVD-схемы на основе параболизованных уравнений Навье — Стокса, позволяющее эффективно проводить расчёты течений со сверх- и дозвуковыми областями. На основе разработанной методики численно моделируется продольный вдув системы пространственных сверхзвуковых турбулентных струй, истекающих из круглых сопел в спутный сверхзвуковой (дозвуковой) поток (рис. 1). Во всей области течения газ считается совершенным, вязким, а режим течения — турбулентным.

1. Постановка задачи

Исходной является система трёхмерных осреднённых по Фавру параболизованных уравнений Навье — Стокса в декартовой системе координат, полученных из полных уравнений Навье — Стокса путём отбрасывания несущественных вязких членов — вторых производных по продольной координате [9]:

5Е + д (Е - Е?) + д (О - О?)

где

Е

дх

ри ри2 + р рии риш (Е + р) и )

дУ

Е

дг

0,

ри рии ри2 + р риш V (Е + р) и У

О

рш

риш

риш 2

рш2 + р

V (Е + р) ш у

Г,

Е? [0, Тжу, Туу, Ту2, итжу + итуу + шту2 °у ]

+ иТу^ + шТ^ — ]Т

УУ

'У2

2 Со \ 2 /о N

3 ке(2— ш), т** = 3 ке(2ш—%)

Т

ху

Ие

и,

^ / , ч к* — (и2 + ) , Оу = —----к-

Ие 1 2 уП Чу (7 — 1) М^РгИе

Т

-I- 41

= ^ Ке

(7 — 1) М^РгИе

Т

1

рс?

р = (7 — 1)

1

1

Е — ^ (ри2 + ри2 + рш2)

Е — ^ (ри2 + ри2 + рш2)

С?> -

7 (7 — 1) М;

Система уравнений (1) записана в безразмерной консервативной форме. В качестве безразмерных параметров приняты характеристики на срезе сопла ро, ио, То, при этом полная энергия и р соответствуют р0и2. Граничные условия данной системы следующие:

х

в струе: и — Т — р — 1, и — ш — 0,

Т 1 Ма ^ 1 в потоке: Т —1, и — —— л/Т, р —-о,

и

ш

ди дш др др

х>^ ду

ди ди др др

— — — — — — — — ш — 0 при г — 0,Ь.

дг дг дг дг

Здесь р — плотность; и, и, ш — продольная и поперечные составляющие скорости соответственно; Е — полная энергия; п — р0/рте — степень нерасчётности; р0, — давление в струе и в потоке; Т — температура; 7 — ср/с? — показатель адиабаты; ср, с? — удельная теплоёмкость при постоянном давлении и объеме; Ма, — число Маха струи и потока; ^ — коэффициент турбулентной вязкости; к — коэффициент теплопроводности; Ие — число Рейнольдса; Рг — число Прандтля.

В качестве замыкания исходной системы уравнений (1) используется алгебраическая модель турбулентности Болдуина — Ломакса [10].

О

2

1

0

0

2. Методика решения

Система уравнений (1) имеет гиперболическо-параболический тип в сверхзвуковых областях. В этом случае задача Коши корректна и её можно проинтегрировать с помощью маршевых методов. В дозвуковых областях система (1) эллиптическая и для организации устойчивого маршевого счёта ее необходимо регуляризовать, обеспечив гиперболическо-параболический тип во всей области течения. Это можно получить, умножив продольную компоненту градиента давления на некоторый регуляризирую-щий параметр и, т. е. вектор потока Е представляется следующим образом:

Е = [ри, ри2 + ир, рии, рит, (Е + р) и]Т.

Как следует из вышеприведенного, продольный градиент давления в уравнении отсутствует, если и = 0, и полностью сохранен, если и = 1. При и = 1 все собственные значения матрицы С = д¥/дЕ* вещественны и положительны, когда и2 + и2 > а2 или М > 1 (в сверхзвуковых областях). В работе [9] показано, что в дозвуковых областях течений собственные значения матрицы будут оставаться вещественными и положительными при сохранении в уравнении хотя бы части градиента давления (0 < и < 1), если

7 М?

и

1 + (7 - 1) М

(2)

В [9] также приводится зависимость и, обозначенная через / (М?), от местного числа Маха М? (рис. 2, линия 1), которая равна единице для М? = 1 и больше единицы для М? > 1. Кривая 1 является граничной линией устойчивости решения, т.е. выбор значений и ниже этой линий обеспечивает устойчивое, а выше — расходящееся решение. Отсюда следует, что при М? > 1 продольный градиент давления может быть полностью включён в уравнение движения. Однако при М? < 1, чтобы собственные значения оставались вещественными и положительными, следует учитывать только часть этого градиента.

Таким образом, при полном сохранении продольного градиента давления в дозвуковых областях решения параболизованных уравнении Навье —Стокса маршевым мето-

Рис. 2. Кривая, ограничивающая градиент давления в дозвуковых областях потока

дом приводят к неустойчивости (расходящиеся решения). На данный момент имеются различные способы преодоления этой трудности. Отметим наиболее известные из них.

Самый простой, предложенный рядом авторов, — вообще не учитывать градиент давления в дозвуковых зонах. Это позволяет интегрировать исходную систему уравнений устойчивым маршевым методом. Однако в дозвуковых зонах для полей течений с большими продольными градиентами давления маршевая схема приводит к ошибкам. Другой метод учёта продольного градиента давления предложен в работе [9]. Здесь в дозвуковой вязкой зоне одна часть продольного градиента давления ш (др/дх) в уравнении сохраняется, а остальная — (1 — ш) (др/дх) — либо опускается, либо рассчитывается на явном слое при помощи разностей назад. Согласно этому методу вектор потока Е расщепляется на две части

Е = Е* + Ер,

где

(

ри

\

р —

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

ри2 + шр Е* = рии , Е

ри'Ш \ (Е + р) и )

Параметр ш вычисляется по уравнению (2) с некоторым коэффициентом запаса о:

( 0 \

(1 — ш) р 0 0 0

ш=

о! МХ

1 + (7 — 1) М;

(4)

Авторы работы [9], исследуя устойчивость схемы с использованием метода Фурье, обнаружили: если рассматриваемый на явном слое градиент давления опущен, то маршевый по пространственной координате метод будет всегда давать устойчивое решение, так как уравнения остаются гиперболическо-параболическими. Если же этот член остается, то возникает неустойчивость при А х < (А х)т1п.

Влияние коэффициента о на устойчивость решения представлено на рис. 2. Здесь приведены кривые функции (4), рассчитанные при о = 0.71 (линия 2), о = 0.9 (3) и о = 1.2 (4). Видно, что кривые 2 и 3 расположены ниже граничной линии устойчивости 1, тогда как кривая 4 находится выше её. Следовательно, выбор параметра о < 1 обеспечивает устойчивость расчёта, а влияние его на точность результатов необходимо анализировать в процессе численных экспериментов. Использование для вычисления ш более упрощённой формулы

ш = f (Мх) = о7МХ (5)

приводит к тому, что в данном случае линия 5 (о = 0.71) находится ниже предельной, а линия 6 (о = 0.9) — выше её.

Таким образом, для применения устойчивого маршевого метода в качестве функций ш достаточно принять любую функцию, линии которой расположены ниже граничной линии 1. В соответствии с этим в настоящей работе используется метод Виньерона и в качестве параметра ш — формула (5). Здесь выбор коэффициента о осуществляется в ходе численных расчётов. Тогда система (1) с учётом (5) запишется в виде

дЕ* дЕР д (Е — Е^) д (О — О^)

дх дх

дУ

дг

0.

(6)

Система (6) решается численно, аппроксимация конвективных слагаемых производится на основе ТУВ-схемы второго порядка точности. Метод построения этой схемы показан ниже на примере одномерного невязкого случая. Далее схема обобщается на двумерный случай. Для этого рассматривается следующая задача:

дЕ* д Е дЕ* дЕ* п

7- + ^ —0, или — + С— — 0, (7)

дх ду дх ду

дЕ

где С — — ЯЛЯ 1 — матрица Якоби, Я, Я 1 — матрицы, составленные из правых и

левых собственных векторов С, Л — диагональная матрица, состоящая из собственных значений матрицы С.

Конечно-разностная аппроксимация (7) запишется в виде

— --:— Л ^¿+1/2- 1/2) , W — Я-1Е*.

Ах

дуЛ ^¿+1/2— УУ;_1/2

Здесь поперечные потоки могут рассматриваться как на п-м, так и на п + 1-м слое; целые индексы относятся к граням ячеек, полуцелые — к центрам ячеек. Схема (8) переписывается следующим образом:

W^1 — Wn — АхЛ+ ^1/2 — ^^ — АхЛ_ (^^¿+1/2 — 1/2) . (9)

Здесь ЛW¿±1/2 — Л+WL±1/2 + Л_WR±1/2, W¿, Wд — значения W соответственно слева и справа от границы с номером г ± 1/2.

Далее внутри каждой ячейки принимается кусочно-параболическое распределение W вида [1]

/д^\п ,( д2^п W (у) — Wn + (У — У;) ] ^ + ^ (У — У;)^] ^ , (10)

где ^ < 1 — коэффициент схемы. Используя для производных в (10) центральные разности, получим

^1/2 — Wn Д + W; + ^ ,

^^¿+1/2 — wn+l — Д+Wi — А+Wi+l, (11)

где

д+wг — Wг+l — Wi, д_wг — Wг — Wi_l, А+Wг+l — ^^¿+2 — Wi+l.

Аналогично находятся предельные значения ^¿_1//2.

Чтобы при построении схемы удовлетворить условиям монотонности, вводится нелинейная функция ограничения — на приращения в виде

А^ — - (г±) А^, г± — А^Т. (12)

После подстановки (11), (12) и перехода от переменных W — инвариантов Римана — к переменным Е* схема (9) запишется как

Ах

Е*п+1 _ Е*п_

; — ; Ау

Н;+1/2 + Я Л+тгпто^ (Я_1А±Е*п, Ь Я_1АТ Е*п) —

1 — К Л^тгптой (Д-1Д±Е*±!1, Ь К-1Д±Е*га) +

4

- Н"-1/2 - К Л±тгпто<1 (К-1Д^Е*га, Ь К-1 Д^Е*^) + К Л-тгптоА (К-1 Д±Е*га, Ь К-1 ДтЕ*га)

:13)

11 3

где Н.,1/2 = -Д+Е,--вдп (КЛК-1) Д+Е,, ф (г+) = тгптой (1, Ьг+), 1 < Ь < --

+ '2 2 1 — — параметр схемы.

Порядок аппроксимации (П.А.) схемы (13)

П.А. = — |^^ (Д х)2 д Е д3Е*

4 Г д Е* ду3 '

Параметр — определяет порядок точности схемы. Третьему порядку соответствует зна-1

чение — = —, второму — — = —1, 0,1. При — = —1 аппроксимация сводится к полностью 3

1

односторонним, при — = 0 — к центральным, при — = —, 1 — к смещённым против

3

потока разностям.

Далее обобщение схемы (13) для двумерной задачи (6) будет иметь вид

Ат

тр*га+1 _ рл*га ^^ Гц- тт , тр* тр* , тр тр ] I

Ег,7 = ЕЧ? Ду [Нг+1/2,^ Н,—1/2,7 + Е г+1/2,^ Е г-1/2,7 + Е^г+1/2,7 Е^г-1/2,7] +

Дт

+ Д [Н,,^'+1/2 — Нг,^'-1/2+ О,7+1/2 — °.-1/2 + °^г,7+1/2 — °^г,7-1/2] —

— (Е. — Е^-1) , (14)

где

Н,± 1/2,7 = — 18дп (Кг± 1/2,7Л,±1/2,.К-±1/2,^ ,

е*+1/2,. = ^К,+1/2,.Л++1/2,.тгптой Д±ЕПП.,ЬЩ^/. ДТЕ".

-1—-Кг+1/2,7Л"+1/2,.тгпт°А (К"+11/2,. Д+ЕП+1,. , Ь К:+1/2,7 Д+Еы ) —

к,+1/2,7 Л:+1/2,7тгпто^ (к-+1/2,. Д-Е".,Ь К:+11/2,. Д+Е?+и).

Аналогично определяются Е*.-1/2, О*7-1/2 и 0*.+1/2, а вязкие члены Е^, О^ — аппроксимируются со вторым порядком точности. Здесь производные с Ер аппроксимируются явно, разностью назад.

Полученная система уравнений (14) решается с применением неявной схемы Би-ма — Уорминга. Численное решение осуществляется согласно принципу расщепления методом матричной прогонки. На каждом шаге по времени в (14) конвективные слагаемые линеаризуются с помощью замены их значений с (п + 1)-го слоя разложениями в ряд Тейлора с известными значениями на предыдущем слое с номером п:

Е*п+1_А™ига+1 е"+1_в™ига+1 О"+1_^гаи™+1 Е"+1_в™ига+1

где А — дЕ*/дU, в — -е/-и ^ — дО/дU, в — дU, ^ — дО^/ди — матрицы

Якоби.

Как видно из построенного алгоритма, для получения решения необходимо знание собственных значений и собственных векторов. Левый собственный вектор в направле-

Оп+1 — ^пип+1, и — [р, ри, ри, рш, Е*]т,

— дЕ/ди, S — дО/ди, в^ — дЕ/ди, й — дО/ди

нии У

Я

1

—^ (с2 + Е2) с2 — Е2 + (и2 + ш2 )(Ь — 1) ш(и2 — ш2) + и2

— Л1с2 + д1и с2 + ш —

(7— 1) и ^(7— 1) и ^(7— 1) — (7— 1)

ш(и2 + ш2) + и2

1

1

и

—Мс2 + Е2)

1

и

—ш— и и ш— и

2

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

— Лбс2 + 05 и с2 + <25и

ш ш— и ш ш— и

05 ш

0

0

05

^ 1,5 — и — иЛ 1,5, Е2 — 2 (7 — 1) (и2 + и2 + ш2) , Ь — 7 (1 — ш) + ш,

01 — ^1(7— 1), 05 — ^5(7— ^ и правый собственный вектор, являющийся обратной матрицей Я-1,

Яу

—УЁА ^2

(—У5 и — шс2и)^1 и^2 — (У5и + с2и)^1 Ы2

р1^э иХ^з иХ^з V

р2^4

и(Я + X )^4

+ X )^4

V,

(шр1--)^3 (шр2--)^4

Ы2

ш

ш

(Л, *р2 — 2шV )^4

У 1,5 — ии (7 (1 — ш) + ш) — Л 1,5а1, Е — а2 — 4а1 а3, к —

Е2

(У1и — шс2и)^5 (У1 и + с2и)^5

(Ь 1

7— 1 2

- (и2 + и2 + ш2) ,

1

с2 ^ 1

— —, — — и(—ши2 — и2), —

с

р1 — X — и2, р — я + X + и2,

с2^ 5

л/Ё'

X

и

2 шс2

(2шс2 + и2 (ш (7— 1) —(7— 1))) , V

и2 (ши2 + и2)

и ( ) Я — — (2 с2 — и2 (ш (7— 1) —(7— 1)))

определяются соответственно из собственных значений матрицы Якоби С, которые име-

ют вид

где

_ — а2 + у/ а2 — 4а1«з

Л 1 — -^-, Л2,3,4

2а1

и, Л 5

и

—а2 — -у/а2 — 4а1а3

2а1

а 1 — — шс — ш (7 — 1) и + 7и , а2 — ш (7 — 1) ии — (7 + 1) ии, а3 — и — с .

В направлении г собственные векторы Я _1, Я 2 и собственные значения матрицы яко-би Е вычисляются аналогично.

с

3. Результаты расчётов и их анализ

Численный расчёт проводился при следующих значениях определяющих параметров: 1 < Ма < 5, 2 < Мто < 10, 1 < п < 10, Т0 — — 1, шаг по маршевой координате варьировался в пределах А х — 0.006 ^ 0.15. Предварительно были проведены тестовые расчёты по влиянию числа узлов на результаты вычислений, вследствие чего выбрана сетка в поперечных направлениях 100 х 100 узлов с шагами А у — А г — 0.15. Ниже описаны результаты численных экспериментов, проведённых с вышеуказанными параметрами.

В качестве тестовой рассматривается задача истечения системы плоских сверхзвуковых струй идеального газа в спутный сверхзвуковой поток при Ма — 3, Мто — 5, п — 10, Т0 — — 1 (рис. 3). Из профилей продольной составляющей скорости (рис. 3, а) и давления (рис. 3, б), рассчитанных по ТУБ-схеме (сплошная линия), методом характеристик в работе [11] (штрих) и по схеме Уорминга — Катлера — Ломакса (штрихпунктир), видно, что разработанная ТУБ-схема, как и метод характеристик, хорошо отслеживает волны возмущения, тогда как в силу влияния искусственных диссипативных членов в схеме Уорминга — Катлера — Ломакса вниз по потоку (х — 7.2, 11.2) профили скорости и давления сглаживаются.

Ниже представлены результаты численного моделирования осесимметричного течения вязкой сверхзвуковой струи в спутный сверхзвуковой поток со следующими параметрами: Ма — 2, Мто — 4, п — 10, Т0 — — 1. На рис. 4 приведено сравнение рассчитанного распределения температуры в поперечном сечении струи при х — 5.5 с данными [12]. Видно, что результаты, полученные по ТУБ-схеме (штриховая линия),

Рис. 3. Расчётные профили продольной составляющей скорости (а) и статического давления (б) сверхзвуковой струи в спутный сверхзвуковой поток: сплошная линия — ТУВ-схема, штрих — метод характеристик, штрихпунктир — схема Уорминга — Катлера — Ломакса; Ма = 3, Мто = 5, = 10, То = = 1

лучше согласуются с расчётными данными [13] для идеальной струи (сплошная линия), при этом данные расчёта [12] с искусственным сглаживанием (кружки) значительно занижены.

На рис. 5 приведены результаты численного моделирования вдува сверхзвуковой струи в спутный сверхзвуковой поток с меньшими значениями параметров струи: Ма = 1.5, М^ = 2, п = 4, Т0 = Т^ = 1, рассчитанные по ТУБ-схеме (жирные линии) и неявной схеме Бима — Уорминга (тонкие линии). Из рисунка следует, что представленные данные практически совпадают. В расчёте по схеме Бима — Уорминга для подавления осцилляции в решение вводятся сглаживающие члены второго и четвёртого порядка с малыми коэффициентами е. Следует отметить, что введённый коэффициент искусственной вязкости по давлению оказался недостаточным фактором для подавления осцилляций при решении задачи вдува сверхзвуковой струи в спутный дозвуковой поток, поэтому был введён дополнительный коэффициент вязкости по температуре. Эта

Рис. 4. Расчётное распределение температуры в поперечном сечении струи при х = 5.5: сплошная линия — для идеальной струи [13], штрих — результаты настоящей работы, о — данные [12]; Ма = 1.5, М^ = 2, п = 4, То = Т^ = 1

Рис. 5. Сравнение режимных параметров на оси струи: а — давление, б — температура; жирные линии — ТУБ-схема, тонкие — схема Бима — Уорминга; Ма = 1.5, М^ = 2, п = 4, То = Т^ = 1

проблема была также указана в работе [14]. Коэффициенты варьировались в таком же диапазоне (0 < £ < 0.3), как и в работе [1]. Таким образом, существенным недостатком численного решения с введением искусственной вязкости является то, что в каждом конкретном случае необходимо подбирать значения коэффициентов вязкости с максимумом на стыке сверх- и дозвуковой областей, т. е. требуется значительное количество численных экспериментов, тогда как предлагаемая ТУБ-схема существенно упрощает расчёт, поскольку не требует введения дополнительных искусственных членов.

Как было показано в методике решения, для обеспечения устойчивого расчёта достаточным является выбор коэффициента запаса а < 1. Вместе с тем его влияние на точность результатов расчёта остается не изученным. Поэтому были проведены численные эксперименты с различными значениями а для задачи истечения сверхзвуковой струи в спутный дозвуковой поток с при Ыа = 3, Ы^ = 0.258, п = 1, Т0 = Тж = 1. На рис. 6

0 50 100 х 150

Рис. 6. Распределение осевой продольной составляющей скорости. Результаты расчётов: а = 0.71 (сплошная линия), а = 0.5 (штриховая), а = 0 (штрихпунктир); результаты экспериментов [15]: о - т = 0, А - т = 0.081, • - т = 0.316; Ма = 3, Мто = 0.258, п = 1, То = Тж = 1

а

0 5 10 15 20 25 30 л- 35

Рис. 7. Изобары в плоскости хОу: а — п =10, б — п = 2; Ма = 3, М^ =4, Т0 = Тж = 1

0.0 7.5 0.0 7.5

Рис. 8. Изотермы в сечениях 13.5 (а), 27.0 (б), 43.9 (в), 60.7 (г); А — То = 2, = 1, Б — То = 1, = 2; Ма = 3, = 5, п = 10

(см. с. 90) приведено сравнение результатов расчёта осевой составляющей скорости при значениях коэффициента запаса о = 0.71 (сплошная линия), о = 0.5 (штрих) и о = 0 (штрихпунктир) с экспериментальными данными [15]. Видно, что уменьшение о приводит к снижению осевой составляющей скорости, что уменьшает согласование с экспериментом. Уменьшение коэффициента запаса в 1.4 раза вызывает снижение скорости в 1.6 раза. Максимальное расхождение между проведенным расчётом и экспериментом имеет место при о = 0. Численное моделирование показывает наилучшее согласование полученных данных с экспериментальными при о = 0.71, что и было принято в дальнейших расчётах.

Ниже представлены результаты влияния основных режимных параметров на характеристики потока.

На рис. 7 (см. с. 90) представлены изотермы в плоскости тОу, рассчитанные при режимных параметрах Ма = 2, = 5, Т0 = = 1 и степени нерасчётности п = 10 (а) и п = 2 (б). Видно, что вследствие разности давлений в струе и в потоке перед струей возникает присоединённая ударная волна, которая распространяется в сторону внешнего спутного потока, достигает оси компоновки и отражается при п = 10 раньше (т ~ 11), чем при п = 2 (т ~ 15). Следовательно, с увеличением п степень смешения струи и потока возрастает.

Поскольку в реальных струях температуры струи и потока различны, то практический интерес представляет исследование их влияния на параметры течения. Результаты расчётов задачи с вдувом сверхзвуковой струи в спутный сверхзвуковой поток (Ма = 3, = 5, п = 10) при То = 2, = 1 (рис. 8, А) и То = 1, = 2 (рис. 8, Б) показывают, что для первого случая температуры струи и потока выравниваются в сечении т = 43.9 (рис. 8, А, г), в то время как для второго случая в этом же сечении (рис. 8, Б, г) наблюдается ярко выраженный слой смешения, который сохраняется почти до т = 100, т.е. при То > смешение струи и потока происходит быстрее, чем при То < Тте.

Таким образом, результаты вычислительных экспериментов показали, что построенная численная модель позволяет проводить расчёты для различных режимных параметров (числа Маха струи и потока, степень нерасчётности, температура струи и потока) и с достаточной точностью отслеживает ударно-волновые процессы.

Список литературы

[1] СинхА Н.Д., Дэш С.М. Расчёты сверхзвуковых течений в каналах при наличии горения, выполненные посредством решения параболизованных уравнений Навье — Стокса // Аэрокосм. техника. 1988. № 7. С. 48-60.

[2] Бондарев Е.Н., Гущин Г.А. Пространственное взаимодействие струй, распространяющихся в спутном сверхзвуковом потоке // Изв. АН СССР. Механика жидкости и газа. 1972. № 6. С. 88-93.

[3] Рогов Б.В., Соколова И.А. Упрощённые уравнения Навье —Стокса для внутренних смешанных течений и численный метод их решения // Изв. РАН. Механика жидкости и газа. 2001. № 3. С. 61-67.

[4] Шифф Л.Б., Стегер Дж,Л. Численный расчёт стационарных сверхзвуковых вязких течений // Ракетная техника и космонавтика. 1980. Т. 18, № 12. С. 16-29.

[5] Йи Г.С., Хартен А. Неявные схемы TVD для гиперболических систем уравнений, записанных в консервативной форме относительно системы криволинейных координат // Аэрокосм. техника. 1987. № 11. С. 11-21.

[6] Harten A. High resolution scheme of the computation of weak solutions of hyperbolic conservation laws // J. of Comput. Phys. 1983. Vol. 49. P. 357-393.

[7] Lawrence S.L., Tannehill J.C. An upwind algorithm for the parabolized Navier — Stokes equations // AIAA J. 1989. Vol. 27, No. 9. P. 1-17.

[8] Lawrens S.L., Chaussee D.S., Tannehill J.C. Application of an upwind algorithm to the three-dimensional parabolized Navier —Stokes equations. AIAA Pap., 1987.

[9] Андерсон Д., Таннехил Дж., Плетчер Р. Вычислительная гидромеханика и теплообмен. Пер. с англ. М.: Мир, 1990.

[10] Иванов М.Я., Крупа В.Г., НигмАтуллин Р.З. Неявная схема С.К. Годунова повышенной точности для интегрирования уравнений Навье — Стокса // Журн. вычисл. математики и матем. физики. 1989. Т. 29, № 6. С. 888-901.

[11] Аверенкова Г.И., Ашратов Э.А., ВолконскАя Т.Г. и др. Сверхзвуковые струи идеального газа // Труды ВЦ МГУ. 1970. Ч. 1. С. 4 -35.

[12] Мышенков В.И. Расчёт течения вязкой ламинарной сверхзвуковой струи в спутном потоке // Журн. вычисл. математики и матем. физики. 1979. Т. 19, № 2. С. 474-485.

[13] Ковалёв Б.Д., Мышенков В.И. Расчёт вязкой сверхзвуковой струи, истекающей в спутный поток // Уч. зап. ЦАГИ. 1978. № 3. С. 125-130.

[14] ДрАммонд Дж.Ф. Численный расчёт вдува звуковой струи водорода перпендикулярно потоку воздуха в канале // Ракетная техника и космонавтика. 1979. Т. 15, № 5. С. 95-97.

[15] Абрамович Г.Н. Теория турбулентных струй. М.: Наука, 1984. 715 с.

Поступила в 'редакцию 29 марта 2011 г., с доработки — 29 июня 2012 г.

i Надоели баннеры? Вы всегда можете отключить рекламу.