Вычислительные технологии Том 19, № 2, 2014
Критические числа Рейнольдса в сверхзвуковом течении Куэтта колебательно-возбуждённого двухатомного газа*
Ю.Н. Григорьев1, И. В. Ершов2 1 Институт вычислительных технологий СО РАН, Новосибирск, Россия 2Новосибирский государственный архитектурно-строительный университет,
Россия
e-mail: [email protected], [email protected]
Григорьев Ю.Н., Ершов И.В. Критические числа Рейнольдса в сверхзвуковом течении Куэтта колебательно-возбуждённого двухатомного газа // Вычисл. технологии. 2014. Т. 19, № 2. С. 20-32.
В рамках нелинейной энергетической теории гидродинамической устойчивости исследовано сверхзвуковое плоское течение Куэтта колебательно-возбуждённого двухатомного газа. Соответствующая спектральная задача для критических значений числа Рейнольдса Recr, определяющих возможное начало ламинарно-тур-булентного перехода, решалась численно с помощью метода коллокаций и QZ-ал-горитма. Показано, что в сверхзвуковом диапазоне рассчитанные значения Recr могут в пределах двух порядков превышать аналогичные значения для дозвуковых чисел Маха.
Ключевые слова: энергетическая теория, гидродинамическая устойчивость, колебательная релаксация, уравнения двухтемпературной газовой динамики, критическое число Рейнольдса, метод коллокаций, QZ-алгоритм.
Grigoryev Yu.N., Ershov I.V. Critical Reynolds numbers in supersonic Couette flow of vibration excited diatomic gas // Comput. Technologies. 2014. Vol. 19, No. 2. P. 20-32.
Supersonic 2D-Couette flow of vibration excited diatomic gas was investigated in the framework of nonlinear energy theory for hydrodynamic stability. The corresponding spectral problem for critical Reynolds numbers Recr defining possible origin of laminarturbulent transit was calculated using collocation method and QZ-algorithm. It was shown that in supersonic range of the calculated values of Recr may be within two orders of magnitude higher than similar values for subsonic Mach numbers.
Key words: energy theory, hydrodynamic stability, two-temperature gas dynamics equations, critical Reynolds number, collocation method, QZ-algorithm.
Введение
В работе [1] на основе нелинейной энергетической теории исследовалась устойчивость дозвукового плоского течения Куэтта колебательно неравновесного молекулярного газа. Выполненное обобщение теории позволило в реальном для двухатомных газов диа-
* Работа выполнена при финансовой поддержке РФФИ (грант № 11-01-00064).
пазоне параметров режима получить значения критических чисел Рейнольдса Recr. При этом найденные значения Recr по порядку величины совпадали с результатами, полученными в аналогичной постановке для несжимаемого течения [2]. Это подтвердило известное представление о том, что дозвуковое течение Куэтта можно считать практически несжимаемым. Вместе с тем немногочисленные исследования [3 - 5], выполненные в постановке классической линейной теории, показывают, что в рамках данного подхода вопрос гидродинамической устойчивости течения Куэтта сжимаемого совершенного газа до последнего времени не имеет однозначного решения. В частности, в работе [3] констатированы абсолютное стабилизирующее влияние вязкости и отсутствие неустойчивости вплоть до чисел Re = 5 • 106 при сверхзвуковых числах Маха M = 2 ^ 5. В то же время в более поздних публикациях [4, 5] были рассчитаны критические числа Рейнольдса в пределах Recr ~ (2 ^ 5) • 104 при числах M = 3 ^ 12.
Таким образом, исходя из результатов численных расчётов [3-5] имеет место следующая ситуация. В рамках линейной теории устойчивости течение Куэтта совершенного газа устойчиво в ближней сверхзвуковой области M < 3 и может проявлять неустойчивость при дальнейшем возрастании числа Маха. Это определило интерес к продолжению начатых в [1] исследований при сверхзвуковом диапазоне M = 2 ^ 5, которым посвящена настоящая работа. Полученные данные для зависимостей критических чисел Рейнольдса от параметров течения представляют самостоятельный интерес. Кроме того, использованная в работе модель двухтемпературной газодинамики при отсутствии возбуждения колебательной моды переходит в модель совершенного газа, что позволяет провести сравнение полученных результатов с данными линейной теории устойчивости [4, 5].
1. Постановка задачи 1.1. Исходные уравнения
Задача устойчивости течения Куэтта колебательно-возбужденного двухатомного газа рассматривается в расчётной области П, представляющей собой прямоугольный параллелепипед, грани которого параллельны координатным плоскостям декартовой системы (x\, x2, x3), а центр совпадает с началом координат. Непроницаемые пластины, вдоль которых направлено основное течение, перпендикулярны оси x2.
Исходной математической моделью служит система уравнений двухтемпературной газовой динамики. В соответствии с физическими представлениями [6-8] эта система описывает течения колебательно-возбуждённого молекулярного газа, когда диссоциацией и возбуждением верхних колебательных уровней молекул, а также поправками на ангармонизм колебаний можно пренебречь.
Как характерные величины для обезразмеривания выбраны полуширина канала L по оси x2, модуль скорости потока U0 на непроницаемых стенках канала, постоянные плотность р0 и температура T0 основного потока, время т0 = L/U0 и давление p0 = р0 Uq .
В качестве невозмущённого потока рассматривается плоское течение Куэтта с линейным профилем скорости и однородным распределением плотности и температур
UsЫ = (x2, 0, 0), Ts(x2) = Tv,s(x2) = Ps(x2) = 1, Ps(x2) = 1/ (7M2) , где ps — статическое давление, постоянное поперёк канала.
Представление мгновенных значений гидродинамических величин возмущённого течения в виде
Р =1+ р', Щ = + п'г, Т =1+ Т', т = 1+ Т'и, р = 1/ (7М2) + р'
позволяет получить [1] из системы уравнений двухтемпературной газовой динамики уравнения для возмущений р', Щ, р', Т', Т' основного течения без ограничения на их амплитуды:
др' др'
дп'г
-¡Т- + пг^~ + Р^
(дп' . дп' тт дп' . I —- + п'—' + П , —' + Щ
Лит+п дх
дг дПя
дхг
дхг
0,
з,3
дх,
дх,
др' 1 д 2п'г 1
+
дхг Ке дх2 Ке
+ тг- «1 + тг
1\ д2п3
3 / дхгдхj
/дТ' . дТ' дТ'\ ~ дп'
+пд" + П-дх^; + 7(7 - 1)М рдх; = КеРг дх2
д2Т' + Тур(Т - Т') +
7(7 - 1)М2
2ке
дп' дх3
дп"
' + "
дх,-
+ 2
дП,, дх,-
+
дП,
, з
дх'
дп' дх3
дп, "
' + 1 +
дх'
+
(
дПя.
\ дх,
+
дП,
, з
дх'
+ 21 а1 - 21
3 / \ дх'
ър
(дТ'
I "-'и
дг
■
+ п,
дТ'
3 дхз 2
+ П
«,3
дТ'
дх,
20
1Ъ
д 2Т'
ЪРСП - Т')
33 Ке Рг дх2
7М2р' = рТ' + р', г = 1, 2, 3, ] = 1, 2, 3.
(5)
Здесь р', п', р', Т', Т' — возмущения плотности и компонент вектора скорости, давления, статической и колебательной температуры газа соответственно. В (1)-(4) и далее по повторяющимся индексам подразумевается суммирование.
Параметры, входящие в уравнения (1)-(5), определяются следующим образом: а1 — относительный коэффициент объёмной вязкости, 7у — параметр, характеризующий степень неравновесности колебательной моды, 7 — показатель адиабаты, И,е, М, Рг — числа Рейнольдса, Маха и Прандтля несущего потока соответственно. Предполагается, что коэффициент а1 и время т в системе (1)-(5) не зависят от статической и колебательной температуры потока и постоянны.
В качестве краевых условий в задаче устойчивости принималось, что при х1 = ±хо/2 и х3 = ±^0/2 возмущения гидродинамических переменных удовлетворяют периодическим граничным условиям, а на непроницаемых границах х2 = ±1 принимают нулевые значения. В расчётах размеры области П по периодическим (однородным) координатам х1, х3 полагались равными длине волны возмущения по соответствующей координате: х0 = 2п/а, г0 = 2п/8 (а, 8 — модули проекций волнового вектора возмущения к на оси координат х1, х3).
1.2. Уравнение энергетического баланса возмущений
В работе [1] из системы (1)-(5) было выведено уравнение энергетического баланса, которое записывается следующим образом:
(Ш* ¿г
= ф = -
(1 + р')п'п3
дП^
3 дхз
+ * а + р')(Т' _ г+ -ММ2 ^
2
2
2
, , „ А IV /ч /дм- 7(7 - 1)М2 Т'
дм- дм —- +---
дж, дж-
2
+
+2
дм- + дМЛ /дЦя, дж, дж- / \ дж,
+
дЦ,
дж-
+
дЦ,, дж,
+
дЦ,
дж-
+ 2< - 2 £
1
+йё
дм-дж,
+1 +§) (Ш +Р
Квадратичная форма
'дТ'\ V 207, /дТ/х 2
дж,-
+
33 \ дж,-
+
£<(() = 2
р^-2 + Т'2 + 7, Т,'2) +
2
7 М2
определяет полную пульсационную энергию возмущений. В последней строке уравнения (6) содержится группа положительно определённых слагаемых, поэтому с уменьшением числа Рейнольдса, начиная с некоторого критического значения Кесг, правая часть этого уравнения становится отрицательной. При этом < 0 и любые началь-
ные возмущения затухают. Критическое число Рейнольдса Кесг соответствует условию dEt/^ = 0 и вычисляется как минимум функционала в правой части энергетического уравнения.
Дальнейшее упрощение уравнения (6) состояло в пренебрежении корреляциями возмущений четвёртого порядка и в частичном разделении переменных под знаком интеграла [1, 8]. При этом зависимости возмущений скорости, плотности и температуры от периодической координаты ж3 представлялись в виде
м/1 = м'/(ж1, ж2) еов(5ж3), м2 = м'2'(жь ж2) еов(5ж3), м'3 = м'3'(ж1, ж2) вт(5ж3), р' = р''(ж1, ж2) еов(5ж3), Т' = Т''(ж1, ж2) еов(5ж3), Т' = Г''(ж1, ж2) еов(5ж3).
(7)
Подстановка (7) в уравнение (6) и интегрирование по переменной ж3 в интервале [-п/5; п/5] приводят к энергетическому уравнению для функций от переменных ж1, ж2
^ = ф'' = _ dt
дм?\2 / дм'/
дж1
+
дж2
+
дм
2, ^дмЬ'Л2, (дм3'х 2
дж2
дж1
+
+
дм3
дж1
+
+
( дм'3'х 2
/5Т
\ дж
дж2 '' \ 2
2
+ ¿2( м1'2 + м2'2 + м"^ + «1 + ^
1\ /дм!' дм'2'
+ 52т ''2 + 207,
33
дм" дж1
+
7, (Г' - Т'')
дТ ''
у V
дж1 2
+
+
3 / \дж1 дж2 дТ \
7
+ <К
дТ'' 2
Рг I \ дж1
+
+ 7 -
дж2
1 + 7 М2 7 М2
+ 52 Т
22
дм''
- 27(7 - 1)М2Т''( д^1+
дж2
Т'М ^ + ^ + 5м3Ч > dS.
дж1 дж2
В результате варьируемый функционал ф'' в правой части становится квадратичным по функциям возмущений м-', Т'', Т'' и их производным.
2
2
2
2
+
1.3. Спектральная задача
Из условия экстремума функционала Ф'' на множестве допустимых функций следуют уравнения Эйлера - Лагранжа, определяющие обобщённую дифференциальную задачу на собственные значения со спектральным параметром И,е:
А .. . 1\дВ 1ЧЛ1Г2 дТ'' Яе
Ам1 +( а + э) д^ - 7(7 - 1)М2 дж" - "2е
м''
7
1 + 7М2\ дТ
7М2 у дж1
,, 1 1\ дВ ~дТ'' Яе
Ам +'а + э) Щ-7(7- 1)М2- "
,, 1 1+ 7М2\ дТ''
м1- 7-
7М2 у дж2
Д„ - ^+0 Д - " (? -^ ''=»■
^Дт.+7 (7 - м+дм!
7
Яе ~2
т
2
2 Т'' - ТЦ ( 1+ 7М'
2 7,--+ 7--тт^-
7М2
, . Т'' Яе 33 (Т' - Т'') =0
РГДТ" - Т 10Т =°.
Здесь
В
дм" дм2'
дж1 ' дж2 +^^ ^ дж2 ' дж2
д2 д2
+
-52.
1 дж2
Амплитудные функции м-', р'', Т'', Т" при ж1 = ±п/а удовлетворяют периодическим граничным условиям, а на непроницаемых границах ж2 = ±1 принимают нулевые значения.
Для периодических по координате ж1 возмущений
(ж1, ж2) = (м1, м2, м3, Т , Т„ ) = ^(ж2) ехр(гаж1),
где ^ = (м, V, ад, 0, ), а — проекция вещественного волнового вектора на координатную ось ж1, г — мнимая единица, спектральная задача (8) сводится к системе обыкновенных дифференциальных уравнений для амплитуд пульсаций м, V, ад, 0, 0,:
1
м'' + га ( а1 + - ) V' -
а2 ( а1 + 3 ) + 52
м + га5 ( а1 + 1 ) ад - 7(7 - 1)М20'
Яе ~2
V +1 - 7 +
1
7 М2
га0
а1 + 4 ) V'' + га ( а1 + 11 м' + 5 Га1 + 1 ) ад' - ^ - 7(7 - 1)М2 га0 =
Яе
т
м + 1 - 7 +
1
7 М2
0'
ад'' - 5 ( а1 + V' - га5 ^а1 + ^ м -
Яе Л 1
= - 7 +
4
52( а1 + |) + в2
ад
0
0
0
Y {0" - k20) + Y(Y - 1)M2 (u' + iav) =
Здесь и далее штрихи у неизвестных функций означают их прозводные соответствующего порядка по переменной х2. Задача (9) замыкается граничными условиями
Как было показано в [1], спектральная задача (9), (10) имеет вещественный спектр собственных значений И,е(а, 8), симметричный относительно осей а = 0, 8 = 0 на плоскости волновых чисел (а, 8). Кроме того, путём сведения системы (9) к спектральному пучку третьего порядка обыкновенных дифференциальных операторов для функции $ = в — в в длинноволновом приближении была получена асимптотическая формула для зависимости критического числа Рейнольдса от параметров задачи
где ао, Ь0 — некоторые положительные постоянные порядка 0(1). 2. Численный метод
Для произвольных значений волновых чисел а, 8 (длин волн возмущений) спектральная задача (9), (10) решалась численно с помощью метода коллокаций [9], в основе которого лежит алгебраическое интерполирование искомого решения по некоторой чебышевской системе функций. В данном случае использовался инструментарий математического пакета МаИаЬ, где в случае непериодических функций применяется интерполяция полиномами Чебышева
Такой выбор обеспечивает экспоненциально быструю сходимость аппроксимирующего ряда для произвольных граничных условий, по крайней мере, в классе бесконечно дифференцируемых функций. На практике этим достигается высокая точность вычислений даже на грубых сетках.
В качестве узлов коллокации (интерполяции) выбирались точки Гаусса — Лобатто x2,n = cos(nn/N), n = 0, 1, ... , N, в которых полином Чебышева N-й степени имеет экстремумы на отрезке [-1, 1]. Дифференциальные операторы первого порядка, входящие в спектральную задачу (9), (10), аппроксимируются на данной сетке матрицей коллокационных производных ülN размером (N + 1) х (N + 1) [9]:
U| x2=±1 = V| Х2 = ±1 = W| X2=±1 = в\ X2=±1 = 1 X2=±1 = 0.
(10)
(11)
Tk(x2) = cos(k arccos x2), k = 0, 1, 2, ...
(-1)i+jsi/[Sj(z - Zj)], l = 3,
<
-Zj/[2(1 - z|)], (2N2 + 1)/6,
-(2N2 + 1)/6, l = 3 = N,
При этом элементы 1-й строки матрицы ДМ являются коэффициентами разностной аппроксимации первой производной в 1-м узле коллокации на шаблоне {х2,п}. Дифференциальные операторы второго порядка аппроксимируются суперпозицией ДМ = ДМ ДМ.
В терминах введённых аппроксимаций задача (9), (10) сводится к обобщённой задаче на собственные значения (линейному матричному пучку) относительно спектрального параметра Л = И,е:
5М+4
¿(Су — Л^.= 0,
3=0
0, 1, ..., 5М + 4.
12)
Здесь вектор неизвестных q размером 5(Ж + 1) состоит из значений собственных функций в узлах коллокации:
q = (щ, Щ, ..., пм, Уо, VI, ..., Ум, 'о, '1, ..., , во, вь ..., вм, в V,о, в„, 1, ..., в V,м).
Матрицы С, Еразмером 5(Ы +1) х 5(Ы +1) вычисляются с использованием специальной процедуры МаШЬ по формулам
С = А 0 ДМ + в 0 ДМ + С 0 1м, Е = В 2 0 ДМ + С 2 0 1м, где матрицы А, В1, С1, В2, С2 размером 5 х 5 имеют вид
(
А
1 0 0 0 0
0 а1 + 4 3 0 0 0
0 0 1 0 0
0 0 0 а 0
0 0 0 0 а
В1
0
га | а1 + 3
С1
V
а2 | а1 + 3 ) + к2 0
—га8 ( а1 + —
0 0
га ( а1 + 3
Ь 0
0
—к2 0
гаЬ 0
В2
/0 0 0 0 0 \
0 0 0 с/2 0
0 0 0 0 0
0 —с/2 0 0 0
\0 0 0 0 0/
1
С2
—8 (а + 3
0 0
га8 ( а1 + 3
0
1
82 ( а1 + 3 | + к2
0 0
0
8 I а1 + 3
0 0
00
гаЬ 0
0
—ак2 0
( 0 1/2 0
1/2 0 0
0 0 0
—гас/2 0 —8с/2
\ 0 0 0
гас/2 0
—8с/2
1>»/т
—Ь 0 00
00
00
00 \
0 0
—ак2
0 0 0
—Yv/т
а = , Ь = 7 (7 — 1)М2, с =1 — 7 + Рг
— 33/(20 т) 33/(20 т) 1
7 М2 ,
1
0
0
0
In — единичная матрица размером (N +1) х (N + 1), а знак ® означает прямое (тензорное) произведение матриц [10].
Однородные граничные условия (10) учитываются неявно через оператор DlN [9] и на дискретном уровне реализуются заменой матриц Dw (l = 1, 2) на окаймлённые матрицы размером (N — 1) х (N — 1). Последние получаются при выполнении условий
= D^j =0, D, 0 = DIn = 0, i = 0, ...,N, j = 0, ...,N, l = 1, 2.
Для нахождения всех собственных значений и функций обобщённой спектральной задачи (12) использовалась процедура Matlab, реализующая QZ-алгоритм [11], который позволяет одновременным ортогональным преобразованием привести пару матриц G, F к обобщённой верхней треугольной форме. В результате применения данной процедуры для фиксированных значений числа Маха M, коэффициента ai, степени неравновесности колебательной энергии yv, времени колебательной релаксации т и каждой пары волновых чисел (а, 5) получался набор из (N + 1) собственных значений, среди которых находилось минимальное по модулю число Рейнольдса Re(a, 5) = | Amin(a, 5)|. Значение критического числа Рейнольдса Recr для данных M, а1, yv, т принималось равным минимальному значению Re во всем диапазоне волновых чисел Re(a, 5): Recr = min [Re(a, 5)]. Затем вычислялись соответствующие Recr собственные функции u, v, w,
(а, ¿) 0, .
Расчёты спектров собственных значений A(a, 5, M, а1, yv, т) проводились в диапазоне волновых чисел а = 0 ^ 10, 5 = 0 ^ 10 при следующих значениях параметров: Yv = 0.111 ^ 0.667, т = 1 ^ 4, а1 = 0 ^ 2, M = 2 ^ 5, Рг = 3/4, y = 7/5. Шаги изменения волновых чисел были выбраны равными ha = h = 0.1. В большинстве расчётов число узлов коллокации в интервале [—1, 1] принималось равным N +1 = 50. Для проверки точности расчётов на основе симметрии спектра выполнялись расчёты в диапазоне волновых чисел а = —10 ^ 10, 5 = —10 ^ 10 и варьировалось число узлов коллокации: N + 1 = 32^ 100.
3. Результаты расчётов и их обсуждение
Выполненные расчёты показали, что при всех рассмотренных значениях параметров задачи минимальные по модулю собственные значения Re(a, ó) = | Am¡n(a, ó)| достигаются на оси а = 0 при ó = 0 в плоскости волновых чисел (а, ó). Изолинии Re(a,ó) приведены на рис. 1.
Как и в случае дозвуковых чисел Маха [1], наиболее "опасными" являются возмущения продольной моды. С учётом периодичности полученного решения по продольной координате £ эти возмущения представляют собой пары двумерных вихрей, вращающихся в противоположных направлениях, с осями, перпендикулярными несущему потоку. Распределение завихренности вычисляется по формуле
ш(ж1,ж2) 1 \( dur\ ( duA .
-=--avi + -— cos ax1 — avr —-— sin ax1 .
Wo V "£2 / V «Ж2)
Здесь ur (x2), ui(x2), vr (x2), vi(x2) — вещественные и мнимые части собственных функций u, v, а нормировочный множитель ш0 представляет собой безразмерный поток за-
Рис. 1. Изолинии поверхностей Ие(а:,5) для М = 3 (а, б, в) и М = 5 (г, д, е) при а = 0, т = 2 (а, г — Yv =0; б, д — 7V = 0.250; в, е — Yv = 0.667; точки на линии 5 = 0 фиксируют критические значения числа Рейнольдса для данного режима)
вихренности через расчётную область (циркуляцию вектора скорости по границе), вычисляемый по формуле
п/а 1
ш0 = J dx1 J dx2 ш(х1 ,х2).
—п/а —1
На рис. 2 представлены примеры изолиний завихрённостей ш(х1,х2) при различных критических числах Рейнольдса Recr(а, а1, Yv, т, M) и значениях амплитуд возмущений скорости, составляющих 10 % значения модуля скорости несущего потока на непроницаемых границах.
Зависимость числа Рейнольдса для продольных мод возмущений от волнового числа а показана на рис. 3. Здесь штрихпунктирные линии соединяют значения абсолютных минимумов на параметризованных по Yv и т кривых Re(a), что позволяет проследить эволюцию Recr. На рис. 4 приведена зависимость Recr от степени неравновесности Yv. Рисунки 1, 3, 4 позволяют констатировать, что с возрастанием параметров а1, Yv, т, M устойчивость течения (критическое число Recr) также возрастает, а соответствующие значения волновых чисел а сдвигаются в сторону более коротких волн. Сопоставление этих результатов с формулой (11) показывает, что длинноволновая асимптотика на качественном уровне правильно воспроизводит зависимость Recr от параметров течения в области волновых чисел а ~ 0(1).
Критические значения числа Рейнольдса Recr (а1, Yv, т, M) приведены в табл. 1, а соответствующие им значения волнового числа а — в табл. 2. Из таблицы 1 следует, что максимальный диапазон изменения Recr при рассмотренных вариациях параметров задачи приближается к полутора порядкам, что существенно больше, чем было получено в [1] для дозвукового течения. Действительно, если при M < 1 и тех же значениях других параметров режима, как и в настоящей работе, критические числа Рейнольдса находились в пределах Recr ~ (1.5 ^ 3.5) • 101, то в данном случае соответствующие пределы составляют Recr ~ 0.5 • 102 ^ 1.5 • 103. Рассматривая степень влияния каждого параметра на Recr при фиксированных значениях остальных параметров, можно
а
'¿г
\ ^Й^Л \ \4\\\5ll 1 i^O.22\ И
\W\\\0i \ \ -0.2? \п ,\0.32 \ Ж mv \\yyvv0 \ 1 wik \
г
У-2
Рис. 2. Изолинии завихренности критических возмущений ш(х1,х2) для М = 3, а\ =0, т = 2 (а — Yv = 0, Ивсг = 106.1; б — Yv = 0.667, Квсг = 201.1; точки на линии х2 = 0 фиксируют максимальное и минимальное значения ш для данного режима)
4 а
4 а
Re 1100
900
700
/ >
/ I i -- ✓
) ( ^^ "2 / 1
б а
6 а
Рис. 3. Зависимость Re(a) для продольных мод возмущений при числах Маха M = 3 (а, б) и M = 5 (в, г) (а, в - ai =0; б, г - ai = 2; 1, 1' - yu = 0.250; 2, 2' - yu = 0.429; 3, 3' -Yv = 0.667; сплошные линии — т = 1, штриховые — т = 3, штрихпунктирные — зависимость критического числа Рейнольдса Recr от волнового числа а)
Recr 1200
600
_______ /
—- —' 2'
" --
1
Recr 1400
б
700
—'—
■
i
0.3
0.6
0.3
0.6 Уъ
а
в
г
а
Рис. 4. Зависимость критического числа Рейнольдса Кесг от степени неравновесности колебательной моды у, (а — «1 =0; б — «1 =2; 1, 1' - М = 2; 2, 2' — М = 3; 3, 3' — М = 4; 4, 4' — М = 5; сплошные линии — т = 1, штриховые — т = 3)
Таблица 1. Критические значения числа Рейнольдса Кесг (а1, 7,, т, М)
М т = 1 т = 4
= 0.111 7, = 0.250 7, = 0.667 = 0.111 7, = 0.250 7, = 0.667
а1 = 0
2 48.22 53.02 59.82 59.82 69.42 82.98
3 139.70 153.62 173.28 173.28 201.08 240.41
4 328.80 361.54 407.81 407.81 473.25 565.80
5 651.37 716.19 807.86 807.86 937.51 1120.85
а1 = 2
2 62.07 68.25 76.99 76.99 89.34 106.81
3 179.83 197.72 223.03 223.03 258.82 309.44
4 423.22 465.34 524.90 524.90 609.14 728.26
5 838.40 921.84 1039.83 1039.83 1206.70 1442.69
Таблица 2. Значения волнового числа а, соответствующие критическим значениям числа Рейнольдса Кесг (а1, 7,, т, М)
т = 1 т = 4
М = 0.111 7, = 0.250 7, = 0.667 = 0.111 7, = 0.250 7, = 0.667
а1 = 0
2 0.672 0.733 0.855 0.855 1.099 1.587
3 0.718 0.784 0.914 0.914 1.175 1.698
4 0.868 0.947 1.105 1.105 1.421 2.052
5 1, 063 1,160 1.353 1.35319 1.739 2.513
а1 = 2
2 1.548 1.689 1.970 1.970 2.533 3.658
3 1.656 1.806 2.107 2.107 2.709 3.913
4 2.001 2.183 2.547 2.547 3.275 4.731
5 2.450 2.673 3.118 3.119 4.009 5.791
заметить, что наибольшее воздействие на возрастание Кесг оказывает рост числа Маха (сжимаемость). При этом в диапазоне М = 2 ^ 5 критическое число Рейнольдса увеличиваются более чем на порядок. В то же время при изменении числа Маха в дозвуковом диапазоне М = 0.2 ^ 0.8 возрастание Кесг лежит в пределах 10 %. Вместе с тем степень влияния коэффициента возбуждения 7, и времени релаксации т, определявших основное воздействие при М < 1, при переходе к сверхзвуковому режиму остаются на прежнем уровне. Тем не менее сделанный в [1] вывод о возможности управления потоком с помощью лазерного возбуждения колебательной моды остается в силе, так как в расчётных пределах изменения 7, число Кесг возрастает в этой работе приблизительно на 30 %.
В заключение отметим, что полученные значения Кесг остаются более чем на порядок ниже критических чисел Рейнольдса, рассчитанных в рамках линейной теории устойчивости для совершенного газа [5]. Кроме того, имеется качественное различие в зависимостях Кесг (М). Если в рассмотренном случае Кесг с ростом числа Маха в диапазоне М = 2 ^ 5 монотонно возрастает, то в рамках линейной теории Кесг в этом диапазоне, наоборот, убывает.
Список литературы
[1] Григорьев Ю.Н., Ершов И.В. Критические числа Рейнольдса в течении Куэтта колебательно возбужденного двухатомного газа. Энергетический подход // ПМТФ. 2012. Т. 53, № 4. С. 57-73.
[2] Гольдштик М.А., Штерн В.Н. Гидродинамическая устойчивость и турбулентность. Новосибирск: Наука. Сиб. отд-ние, 1977.
[3] Duck P.W., Erlebacher G., Hussaini M.Y. On the linear stability of compressible plane Couette flow // J. Fluid Mech. 1994. Vol. 258. P 131-165.
[4] Hu S., Zhong X. Linear stability of viscous supersonic plane Couette flow // Phys. Fluids. 1998. Vol. 10, No. 3. P. 709-729.
[5] Malik M., Dey J., Alam M. Linear stability, transient energy growth, and the role of viscosity stratification in compressible plane Couette flow // Phys. Rev. E. 2008. Vol. 77, iss. 3. P. 036322(15).
[6] Жданов В.М., АлиЕвский М.Е. Процессы переноса и релаксации в молекулярных газах. М.: Наука, 1989.
[7] Нагнибеда Е.А., КустовА Е.В. Кинетическая теория процессов переноса и релаксации в потоках неравновесных реагирующих газов. СПб.: Изд-во С.-Петербургского гос. ун-та, 2003.
[8] Григорьев Ю.Н., Ершов И.В. Устойчивость течений релаксирующих молекулярных газов. Новосибирск: Изд-во СО РАН, 2012.
[9] Canuto C., Hussaini M.Y., Quarteroni A., Zang T.A. Spectral Methods in Fluid Dynamics. Springer Ser. in Comput. Phys. Berlin: Springer-Verlag, 1988.
[10] Корн Г., Корн Т. Справочник по математике. М.: Наука, 1973.
[11] Moler C.B., Stewart G.W. An algorithm for generalized matrix eigenvalue problems // SIAM J. Numer. Anal. 1973. Vol. 10, No. 2. P. 241-256.
Поступила в 'редакцию 30 октября 2013 г.