close

Вход

Забыли?

вход по аккаунту

?

Главный инвариант для подсеточных моделей турбулентности.

код для вставкиСкачать
МАТЕМАТИКА
N. Rozov
BUFFERNESS PHENOMENON IN THE MATHEMATICAL MODELS
OF SCIENSE AND TECHNOLOGY
A mathematical model of the distributed van der Pol generator is under consideration. The model is described by the linear system of telegraph equations
with non-linearity in one of the boundary conditions; the parameters of the model
represent physical properties of the generator. The bufferness phenomenon appears in the model under consideration, i.e. simultaneous existence of any a priori defined finite number of stable time-periodic solutions of the considered
boundary value problem is discovered under the appropriate choice of parameter
values.
А. М. Балонишников, А. В. Копыльцов
ГЛАВНЫЙ ИНВАРИАНТ
ДЛЯ ПОДСЕТОЧНЫХ МОДЕЛЕЙ ТУРБУЛЕНТНОСТИ
На основе упрощенного анализа уравнений для мелкомасштабных компонент скоростей турбулентного потока несжимаемой жидкости и соображений размерности получено новое выражение для подсеточного коэффициента турбулентной вязкости. Это выражение в отличие от
известной модели Смагоринского определяется некоторым инвариантом
тензора градиента крупномасштабной скорости, который явно содержит
завихренность крупномасштабной скорости движения, что приводит к
уменьшению коэффициента турбулентной вязкости.
Введение
Подсеточное моделирование было введено как попытка преодолеть трудности прямого численного моделирования развитой турбулентности на основе
исходных уравнений Навье—Стокса, поскольку современным компьютерам
не хватает быстродействия для расчета динамики вихрей [1]. В подсеточных
моделях выбирают некоторый масштаб h — шаг разностной сетки по пространству — и рассматривают вихри, диаметр которых превосходит h. Из-за
нелинейности уравнений Навье—Стокса, описывающих движение несжимаемой жидкости, возникают существенные трудности. В уравнения, описывающие крупные вихри, входят подсеточные напряжения Рейнольдса, которые
определяются скоростями мелкомасштабных вихрей, чей диаметр меньше,
чем h. Подсеточные напряжения Рейнольдса зависят от масштаба h и градиентов (и/или более высоких производных) скоростей больших вихрей и описывают переход энергии от крупных вихрей к более мелким и с последующим
переходом в тепло. Масштаб h не должен превышать диаметра энергосодержащих вихрей [2]. В идеальном случае выражение для подсеточных напряжений Рейнольдса желательно выводить аналитически из уравнений для скоростей мелкомасштабных вихрей с последующим осреднением произведения
компонент мелкомасштабной скорости, которое определяет напряжения Рейнольдса, по объему кубической ячейки с длиной ребра h. В этом случае подсеточное моделирование было бы свободно от эмпирических констант, и такую
модель можно было бы применять для всех типов течений, включая внутренние течения, свободную турбулентность, ламинарно-турбулентный переход. К
96
Аналитическое решение граничных задач для семейства БГК-уравнений…
настоящему моменту исследователи далеки от решения этой задачи. Различные подходы к подсеточному моделированию отражены, например, в работе
[3]. Простейшей и наиболее используемой параметризацией является следующее выражение Смагоринского для подсеточного коэффициента турбулентной вязкости ν t [4], ν t = ch 2 S , величина которого определяет напряжения
Рейнольдса, где c — некоторая константа, h — шаг разностной сетки по пространству при дискретизации уравнений гидродинамики, S — инвариант тензора деформации крупномасштабной скорости:
S = 2 S ij S ji , S ij = 12 (∂ iU j + ∂ j U i ) ,
где U i (i = 1, 2, 3) — компонента крупномасштабной скорости. В действительности коэффициент турбулентной вязкости может зависеть и от антисимметричной части Ωij = 12 (∂ iU j − ∂ jU i ) тензора градиента скорости. Имеется простая связь этой антисимметричной части с завихренностью
крупномасштабного движения Ω i (Ω = rotU ): Ωi = 2ε ijk Ω jk , где ijk — тензор
ε
Леви—Чевитта. В работе [5] показано, что всего для несжимаемой жидкости
имеется пять независимых инвариантов симметричной и антисимметричной
части тензора градиента деформации крупномасштабной скорости:
Tr{S 2 }, Tr{Ω 2 }, Tr{S 3 }, Tr{Ω 2 S}, Tr{Ω 2 S 2 } ,
где Tr — след матрицы. Получить зависимость коэффициента турбулентной
вязкости от этих пяти инвариантов эмпирическим путем весьма затруднительно [6].
Зная коэффициент турбулентной вязкости как функцию производных
крупномасштабной скорости, можно определить подсеточные напряжения
Рейнольдса и тем самым получить замкнутое описание динамики крупномасштабных вихрей.
Вывод выражения для коэффициента подсеточной турбулентной вязкости из упрощенного уравнения для мелкомасштабной скорости и соображений размерности
Покажем, что исходя из максимально упрощенной теории устойчивости
имеется единый основной инвариант, который и должен определять подсеточный коэффициент турбулентной вязкости. Если сделать предположение,
что крупномасштабную скорость и ее тензор градиента можно считать постоянными при рассмотрении динамики мелкомасштабной компоненты скорости
[7], т. е. отказаться от широко используемого предположения Праудмена—
Бэтчелора (Proudman—Batchelor) [8] о линейности крупномасштабной или
средней скорости (что дает вклад в спектр скорости мелкомасштабных вихрей), то уравнение для Фурье-компонент мелкомасштабной скорости будет
иметь следующий вид [9]:
(∂ t + νk 2 )uα (k , t ) = − 2i Pαβγ (k )∑ u β (p, t )uγ (k − p, t ) + Pαβ ∂ mU β u m (k , t ) ,
(1)
p
где Pij = δ ij − k i k j / k 2 , Pαβγ = k β Pαγ + kγ Pαβ , индексы в этих формулах пробегают
значения 1, 2, 3 (здесь и всюду по тексту по повторяющимся индексам подра97
МАТЕМАТИКА
зумевается суммирование). Если жидкость или газ для простоты считать несжимаемыми ki ui = 0 , то можно трансформировать уравнение (1) к двум независимым компонентам мелкомасштабной скорости, направленным вдоль единичных векторов ε 1j и ε 2j (j = 1, 2, 3), ортогональным вектору k, как это было
сделано для изотропной турбулентности в отсутствие градиентов крупномасштабной скорости dU [10]:
∑Φ
(∂ t + νk 2 )u γ (k , t ) = a γµ u µ +
γαβ
u α (p, t )u β (q, t ),
(2)
p ,q , p + q = k
где a γµ = ε γj ε mµ ∂ mU j , Φ γαβ (k , p, q) = −ik m ε γj (k )ε αj (p)ε mβ (q), верхние индексы пробегают значения 1, 2, а нижние — значения 1, 2, 3.
Входящие в эти формулы вектора имеют следующие декартовы компоненты, выраженные через углы Эйлера [10]:
k = (k cosθ cosη , k sin θ cosη , k sin η ) ,
где cosθ = k1 / k " , sin θ = k 2 / k " , cosη = k " / k , sin η = k 3 / k , k " = k12 + k 22 .
Если выбрать единичный вектор e = (cosθ , sin θ ,0), тогда единичный вектор ε 1 (k ) определяется соотношением
ε 1 = e × k / | e × k |,
или
ε 1 = (sin θ ,− cosθ ,0) ,
а единичный вектор ε 2 (k ) определяется соотношением
ε 2 = k × ε 1 (k ) / | k × ε 1 |,
или
ε = (cosθ sin η , sin θ sin η ,− cosη ).
2
При инверсии вектора k справедливы соотношения
ε 1 (k ) = −ε 1 (−k ), ε 2 (k ) = ε 2 (−k ).
Верно соотношение ε 1 × ε 2 = k / k . Связь обычных Фурье-компонент мелкомасштабной скорости с ее поляризационными компонентами выражается
соотношением
ui (k , t ) = ε iµ (k )u µ (k , t ) ,
(3)
причем
k ⋅ ε µ = 0, ε µ (k ) ⋅ ε λ (k ) = δ µλ , ε iµ (k )ε µj (k ) = Pij (k ) .
(4)
С помощью замены переменных u = Bv, где матрица B составлена из собственных векторов матрицы A = {aαβ } , можно привести матрицу A к диагональной форме (в общем случае к жордановой матрице) [11], а само уравнение для поляризационных компонент мелкомасштабной скорости примет вид
(∂ t , + νk 2 )v λ = J λµ v µ + ( B −1 ) λγ ∑ Φ γαβ ( Bv)α ( Bv) β ,
98
(5)
Аналитическое решение граничных задач для семейства БГК-уравнений…
где J λµ — элемент диагональной матрицы, на диагонали которой стоят собственные числа матрицы A : λ1 и λ 2 , которые определяются как
λ1, 2 =
P
±
2
P2
− Q,
4
(6)
где P = ε 1jε m1 ∂ mU j + ε 2j ε m2 ∂ mU j ,
Q = (ε 1jε m1 ∂ mU j )(ε l2ε p2 ∂ pU l ) − (ε 1jε m2 ∂ mU j )(ε l2ε 1j ∂ jU l ) .
Если разложить dU на симметричную и антисимметричную части
∂ mU j = S jm + Ω jm и ввести единичный вектор n c компонентами ni = ki / k , то,
используя свойства векторов ε λ и условие несжимаемости, можно выражения
для P и Q упростить:
P = − ni n j S ij ,
Q = (ε 1j ε m1 S mj )(ε l2ε p2 S lp ) − (ε 1j ε m2 S mj ) 2 + 14 (n ⋅ Ω) 2 .
(7)
(8)
Покажем, что в действительности величина Q определяется непосредственно ориентацией вектора k (или единичным вектором n), исключив вспомогательные вектора ε λ :
Q = (ε m1 ε l2 − ε l1ε m2 )ε 1j ε 2p S mj S lp + 14 (n ⋅ Ω) 2 .
Введем антисимметричный тензор τ ml = ε m1 ε l2 − ε l1ε m2 . Выписав члены в Q ,
отличные от 0 с m ≠ l , получим после очевидных преобразований
Q = (τ 12 S1 j S 2 p + τ 13 S1 j S3 p + τ 23 S 2 j S3 p )τ jp + 14 (n ⋅ Ω) 2 .
Нетрудно убедиться, что τ 12 = n3 , τ 23 = n1 , τ 13 = −n2 . Поэтому
Q = n32 ( S11S 22 − S12 S 21 ) + n22 ( S11S33 − S13 S31 ) + n12 ( S 22 S33 − S 23 S32 ) + 2n1n2 ( S 23 S13 −
(9)
S12 S 33 ) + 2n1 n3 ( S12 S 23 − S 22 S13 ) + 2n2 n3 ( S12 S13 − S11 S 23 ) + 14 (n ⋅ Ω) 2 .
В инвариантном виде
Q = − 12 Tr[(n × S) 2 ] + 14 (n ⋅ Ω) 2 ,
где Tr[(n × S) 2 ] = ε ijk n j S kpε pab na Sbi .
Окончательно имеем следующее выражение для собственных чисел матрицы линеаризированных уравнений для мелкомасштабной скорости (без вязких членов):
λ1, 2 = − 12 Tr (CS) ±
1
4
[Tr (CS)]2 + 12 Tr[(n × S) 2 ] − 14 (n ⋅ Ω) 2 ,
(10)
где C ij = n i n j .
Согласно принципу подчинения Хакена [12], неустойчивые моды подчиняют себе устойчивые. Если имеется набор неустойчивых мод, подчиняющихся одному и тому же уравнению, то следует ожидать, что наиболее существенный вклад в подсеточные напряжения Рейнольдса дадут наиболее
неустойчивые моды. Если дискриминант квадратного уравнения, определя99
МАТЕМАТИКА
ющий собственные числа λ1 и λ2 , положителен, то эти собственные числа
действительны. Если оба собственных числа положительны, или одно положительно, а другое отрицательно, то, по-видимому, основной вклад дадут моды с наибольшим числом λ1 .
Среди ориентаций вектора k следует выбрать то направление, которое дает
максимальное значение λ1 .
Для нахождения этого максимума
по угловым переменным будем счи-
тать, что предварительно был осуществлен переход с помощью ортогонального преобразования в систему координат, в которой тензор S диагонален. Тогда
выражения (7) и (10), определяющие величины P и Q, существенно упрощаются:
P = − n12 S1 − n22 S 2 − n32 S3 ,
Q = n32 S1S 2 + n22 S1 S3 + n12 S 2 S3 +
(11)
1
(n ⋅ Ω) 2 ,
4
(12)
где S1 = S11 , S 2 = S 22 , S 3 = S 33 — собственные числа тензора S, S1 ≥ S 2 ≥ S3 .
Условие несжимаемости будет иметь вид:
S1 + S 2 + S 3 = 0 .
(13)
Компоненты единичного вектора n удовлетворяют соотношению
n12 + n22 + n32 = 1 .
(14)
Используя соотношения (11)–(14), выражение (10) для величины λ1 принимает вид:
λ1 = − 12 [(2n12 + n22 − 1) S1 + (2n22 + n12 − 1) S 2 ] + Z,
(15)
где
Z=
1
4
[(2n12 + n22 − 1)S1 + (2n22 + n12 )S2 ]2 + [n12 S 2 (S2 + 2S1 ) + n22 S1 (S1 + 2S2 ) − S1S2 ] − R ,
R = (n1Ω1 + n2 Ω 2 + n3Ω 3 ) 2 , n3 = ± 1 − n12 − n22 .
(16)
Необходимыми условиями максимума функции λ1 (n) являются
∂λ1
= 0;
∂n1
∂λ1
=0.
∂n2
Используя соотношение (15), получим систему двух алгебраических уравнений:
(2 S1 + S 2 )n1 Z − (2 S1 + S 2 ) n1[(2n12 + n22 − 1) S1 + ( 2n22 + n12 − 1) S 2 ] − n1S 2 (2 S1 + S 2 ) +
+ (Ω1 +
∂n3
Ω 3 )(n1Ω1 + n2 Ω 2 + n3 Ω 3 ) = 0,
∂n1
(17)
(2S 2 + S1 )n2 Z − (2S 2 + S1 )n2 [(2n1 + n22 − 1) S1 + (2n22 + n12 − 1) S 2 ] − n2 S1 (2S 2 + S1 ) +
+ (Ω 2 +
∂n3
Ω 3 )(n1Ω1 + n2 Ω 2 + n3Ω 3 ) = 0 ,
∂n2
100
(18)
Аналитическое решение граничных задач для семейства БГК-уравнений…
где
∂n3
∂n3
= ∓ n1 / 1 − n12 − n22 ,
= ∓ n2 / 1 − n12 − n23 .
∂n2
∂n1
Система уравнений (17) и (18) не становится проще при переходе к угловым переменным η , θ . Однако в случае Ω = 0 можно получить решение
λ1 max = 2S 1 . В других случаях необходимо использовать приближенные численные методы, например метод Ньютона.
Таким образом, соответствующий основной инвариант будет I = λ1 max , если подкоренное выражение в (15) положительно, где максимум берется по
различным направлениям вектора k . Если подкоренное выражение отрицаI = sup{Re(λ1 (n))} при sup(Re(λ1 )) > 0
или
I = 0 при
тельно, то
sup{Re(λ1 )} ≤ 0 , где Re — действительная часть числа. В этом случае возникновение неустойчивостей возможно лишь за счет нелинейных взаимодействий, а диссипацию будут обеспечивать вязкие члены. Физически первый случай (подкоренное выражение в (15) положительно) соответствует
превалированию деформации крупномасштабным полем скорости мелкомасштабных вихрей над эффектом вращения. Второй случай (подкоренное выражение в (15) равно 0 или отрицательно) соответствует обратной ситуации. Во
втором случае ячейка со стороной h находится в ядре когерентного вихря [13],
где использовался более грубый метод индентификации когерентных структур: инвариант q тензора dU положителен, где q = 12 (Ωij Ω ji − Sij S ji ) .
В итоге из соображений размерностей можно считать, что коэффициент
подсеточной турбулентной вязкости равен
ν t = ch 2 I ,
(19)
где с — константа, которую следует определить путем сравнения численных
экспериментов с натурными экспериментами, h — шаг разностной сетки.
Обсуждение результатов
Отметим, что эффект уменьшения коэффициента турбулентной вязкости
при увеличении крупномасштабной завихренности — явление достаточно известное. Так, в [14] была предложена крупномасштабная феноменологическая
модель с учетом наличия крупномасштабной завихренности с коэффициентом
турбулентной вязкости:
ν t = ν 0 (k 2 / ε )[1 /(1 + 0,36kΩ / ε )],
где ν 0 — константа, k — кинетическая энергия турбулентности, ε — удельная скорость диссипации турбулентной энергии, Ω = WijW ji — модуль крупномасштабной завихренности, Wij — антисимметричная часть тензора градиента скорости dU. Похожее выражение использовалось в [15] для
коэффициента турбулентной вязкости ν t :
ν t = cµ f µ kτ 1 ,
101
МАТЕМАТИКА
где f µ — коэффициент, зависящий от числа Рейнольдса, τ 1 — характерное
время оборота энергосодержащих вихрей: τ 1 = k / ε + ν / ε , ν — кинематический коэффициент молекулярной вязкости,
cµ =
1
,
A1 + Asτ 1 max(S , Ω)
где A1 = 8, As = 3 cos Φ, Φ = 13 arccos( 6 x), x = 21.5
Sij S jk S ki
.
S3
В этих двух работах ν t → 0 как 1 / Ω при Ω → ∞. В нашем подходе, как
следует из (10), ν t пропорционален aS + b 2 S 2 − Ω 2 , если Ω ≤ bS , где a, b —
константы; при Ω > bS ν t пропорционален S, т. е. коэффициент турбулентной
вязкости в первом приближении перестает зависеть от завихренности.
Какая ситуация более адекватна действительности, в настоящий момент
сказать трудно. Величины k , ε могут, кроме того, непосредственно зависеть
от крупномасштабной завихренности Ω через свои уравнения переноса.
Таким образом, получено явное выражение основного инварианта I для
турбулентного течения несжимаемой жидкости, которое в первую очередь
определяет подсеточные эффекты произвольного неоднородного турбулентного потока. Впервые учтено явным образом влияние крупномасштабной завихренности течения на подсеточный коэффициент турбулентной вязкости.
Определение эмпирической константы c требует проведения дополнительных численных расчетов на основе предложенной модели для различных типов турбулентных течений. Теоретическое определение этой эмпирической
константы и, возможно, непосредственное определение компонент тензора
Рейнольдса через градиенты крупномасштабной скорости требуют более глубокого анализа нелинейных уравнений для мелкомасштабной скорости (2).
Предложенный подход может быть распространен на магнитогидродинамическую турбулентность, сжимаемые течения, течения с химическими реакциями и т. п.
БИБЛИОГРАФИЧЕСКИЕ ССЫЛКИ
1. Piomelli U. Large-eddy simulation: Achievements and challenges // Prog. Aerosp. Sci.
1999. Vol. 35. Р. 335.
2. Balonishnikov A. M. Extended local balance model of turbulence and Couette-Taylor flow
// Physical Review E. 2000. Vol. 61. № 2. Р. 1390–1394.
3. Voelkl T., Pullin D. I., Chang D. C. A physical-space version of the stretched-vortex subgrid-stress model for large-eddy simulation // Physics of Fluids. 2000. Vol. 12. № 7. Р. 1825.
4. Smagorinsky J. General circulation experiments with the primitive equations // Mon.
Weather Rev. 1963. Vol. 91. Р. 99–105.
5. Pope S. B. A more general effective viscosity hypothesis // J. Fluid Mech. 1975. Vol. 72.
P. 331–340.
6. Jorgen T., Gatski T. B. General explicit algebraic stress relations and best approximation
for three-dimensional flows // Int. J. Eng. Science. 1998. Vol. 36. P. 739–763.
7. Скворцов Г. Е., Тимохов Л. А. К теории турбулентности // Вестник ЛГУ. 1980.
Вып. 3. № 13. С. 106–110.
8. Batchelor G. K., Proudman J. The effect of rapid distorsion on a fluid in turbulent motion
// Q. J. Mech. Appl. Math. 1954. Vol. 7. P. 83–90.
102
Аналитическое решение граничных задач для семейства БГК-уравнений…
9. Балонишников А. М. Однопараметрическая модель неоднородной гидродинамической турбулентности // Системы автоматизации в науке и производстве. — СПб., 1984.
С. 50–54.
10. Lee J. The triad-interaction representation of homogeneous turbulence // J. Math. Phys.
1975. № 7. Р. 1359–1366.
11. Брюно А. Д. Локальный метод нелинейного анализа дифференциальных уравнений.
— М., 1979.
12. Хакен Г. Синергетика. — М., 1980.
13. Dubief Y., Delcayre F. On coherent-vortex identification in turbulence // Journal of Turbulence. 2000. Vol. 1. P. 2–22.
14. Zhou Y. A phenomenological treatment of rotating turbulence // Physics of Fluids. 1995.
Vol. 7. Р. 2094.
15. Merci B., De Langhe C., Vierendeels J., Dick E. A quasi-Realizable Cubic Low-Reynolds
Eddy-Viscosity Turbulence Model with a New Dissipation Rate Equation // Flow, Turbulence
and Combustion. 2001. Vol. 64. P. 133–160.
А. Balonishnikov, A. Kopyltsov
THE MAIN INVARIANT FOR SUB-GRID MODELS OF TURBULENCE
A new expression for the sub-grid coefficient of turbulent viscosity based on a
simplified analysis of equations for small-scale components of the incompressible
turbulent flow and the dimensional analysis has been obtained. This expression in
contrast with the Smagorinsky model is determined by an invariant of the tensorgradient of the large-scale velocity, which contains explicitly the large-scale vortices, which diminishes the coefficient of turbulent viscosity.
103
Документ
Категория
Без категории
Просмотров
6
Размер файла
370 Кб
Теги
главные, инвариантов, моделей, турбулентность, подсеточной
1/--страниц
Пожаловаться на содержимое документа