close

Вход

Забыли?

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

?

Фрик П.Г. - Турбулентность- модели и подходы Часть 1(1998).pdf

код для вставкиСкачать
1
М инистерство общ его и профессионального образования Российской Ф едерации
П ермский государственный технический университет
Кафедра математического моделирования систем и процессов
П .Г.Ф рик
ТУРБУЛЕН ТН О СТЬ:
М ОДЕЛИ И П ОДХОДЫ
Курс лекций
Часть I
Рекомендовано учебно-методическим отделом по направлению
«Электроника и прикладная математика» в качестве учебного пособия
для студентов специальности «П рикладная математика»
П ермь 1998
2
УДК 532.517.4
Турбулентность: модели и подходы. Курс лекций. Часть I /
П .Г.Ф рик; П ерм. гос. техн. ун-т. П ермь, 1998. 108 с.
П ервая часть курса лекций вклю чает в себя введение и три из семи
разделов курса «Турбулентность: модели и подходы». П ервый раздел содержит базовые сведения из механики жидкости, необходимые для дальнейш его изложения. Второй посвящ ен вопросам, связанным со стохастическим поведением маломодовых систем гидродинамического типа. В третьем разделе выводятся уравнения для статистических моментов пульсаций
скорости и дается краткий обзор моделей, используемых для их замыкания.
Для студентов и аспирантов физико-математических специальностей.
И л.64. Библиогр. 12 назв.
Рецензенты:
кафедра физики П ермского
государственного технического университета,
д-р физ.-мат.наук, профессор Д.В.Любимов
ISBN
©
1998
П ермский государственный
технический
университет,
3
ВВЕДЕН И Е ........................................................................................................................4
1
ОСН О ВЫ ......................................................................................................................7
1.1
Уравнения движения жидкости..........................................................................................................7
1.2
Устойчивость течений.......................................................................................................................21
1.3
Свободная конвекция несжимаемой жидкости................................................................................26
1.4
Конвективная устойчивость.............................................................................................................31
1.5
М аломодовая модель конвекции (система Лоренца) ......................................................................37
2
ХАОС В Д И Н А М И Ч ЕСК И Х СИСТЕМ АХ ............................................................. 42
2.1
Консервативные и диссипативныесистемы .....................................................................................43
2.2
Бифуркации.......................................................................................................................................50
Как описать переход и хаос ? ..........................................................................................................................52
2.4
Спектры Ф урье..................................................................................................................................58
2.5
Странный аттрактор.........................................................................................................................63
2.6
Ф ракталы...........................................................................................................................................67
2.7
Субгармонический каскад................................................................................................................74
2.8
Н екоторые примеры .........................................................................................................................79
3
ПОЛУЭМ П И РИ Ч ЕСКИ Е М ОДЕЛИ ....................................................................... 92
3.1
Развитая турбулентность..................................................................................................................92
3.2
Уравнения для статистических моментов........................................................................................99
3.3
Турбулентная вязкость................................................................................................................... 102
3.4
Длина пути смеш ения...................................................................................................................... 103
3.5
М одели переноса турбулентной вязкости...................................................................................... 105
3.6
Двухпараметрическиемодели........................................................................................................ 105
4
ВВЕДЕН И Е
Турбулентность остается одним из наиболее сложных объектов исследования механики жидкости и газа. За почти столетню ю историю ее
изучения предложены десятки различных подходов, почти всегда отражаю щ ие наиболее активно развиваемые перспективные направления математики и физики соответствую щ его периода времени. Статистическая физика и теория вероятности, теория размерности, фурье анализ и прямые
численные методы, теория динамических систем, теория фракталов и вейвлет-анализ- вот далеко не полный перечень областей науки, которые давали основные идеи исследователям турбулентности.
Теория турбулентности далека от своего заверш ения. П родолжаю т
появлятся и все новые подходы к ее изучению . Растет число моделей, предлагаемых для лучш его понимания отдельных ее свойств. Дать представление об основных идеях, движущ их этот процесс, продемонстрировать возможности различных подходов и показать проблемы, ими не разреш енные,
представить современные модели, не вош едш ие ещ е в учебники и не ставш ие хрестоматийными - вот цель предлагаемого курса лекций.
Курс предназначен для студентов специальности "прикладная математика", ориентирую щ ихся на работу в научно-исследовательских учреждениях и на кафедрах, в особенности тех, что связаны с реш ением задач механики жидкости и газа. В то же время, в курсе рассматриваю тся и общ ие
подходы к моделированию сложных динамических систем, которые могут
быть полезными специалистам, занимаю щ имся моделированием самых
различных (и не только механических) систем и явлений. Курс рассчитан на
студентов, получивш их ш ирокую базовую подготовку по основным математичеким дисциплинам, вклю чая методы математической физики, функциональный анализ и теорию вероятности, а также прослуш авш их спецкурсы по механике (механику сплош ных сред, теорию определяю щ их соотнош ений).
Курс лекций состоит из двух частей. В первую часть вклю чены три
главы, вклю чаю щ ие в основном сведения, которые можно найти в различных учебниках и монографиях, но собранныевоедино и изложенныевсвете
задач, обсуждаемых в этом курсе. Вторая часть содержит результаты, которые, за редким исклю чением, не вош ли ещ е в книги и могут быть найдены только в оригинальных статьях.
П ервая глава содержит базовые сведения по динамике несжимаемых
жидкостей, вклю чая вывод уравнений движения для идеальной и вязкой
жидкости и примеры задач, имею щ их точные реш ения. Даны основы тео-
5
рии устойчивости, имею щ ей важнейш ее значение в понимании проблем перехода от ламинарных течений к турбулентным. П одробно обсуждаются
две задачи : устойчивость плоского течений П уазейля (задача ОрраЗоммерфельда) и задача Релея о конвективной устойчивости подогреваемого снизу горизонтального слоя несжимаемой жидкости. П оследняя задача предворяется выводом уравнений свободной конвекции в приближении
Буссинеска и обсуждением необходимых условий устойчивости неоднородно нагретой жидкости, находящ ейся в поле сил тяжести. Особое внимание уделяется вопросу о безразмерном представлении уравнений движения,
о законах подобия и о безразмерных параметрах и их роли в описании
процессов перехода к хаотическому поведению . Глава заканчивается выводом маломодовой модели конвекции (модель Лоренца). Этот вывод имеет
методическую цель - показать и обсудить проблему проектирования нелинейных уравнений движения на конечномерный базис и переход от уравнений в частных производных к обыкновенным дифференциальным уравнениям. В то же время подробный вывод модели полезен, так как полученная
система уравнений ш ироко используется в следую щ ей главе, где подробно
обсуждаю тсяеесвойства.
Значительный прогресс в понимании природы и свойств турбулентности произош ел в последние десятилетия благодаря успехам теории динамических систем, позволивш им понять как хаотическое поведение возникает в детерминированных системах. Этим результатам посвящ ена вторая
глава, в которой приводятся базовые сведения из теории динамических
систем и обсуждаю тся некоторые приложения. Вводится понятие фазового
пространства и даны примеры фазовых портретов некоторых простых динамических систем. Обсуждаются особенности эволю ции консервативных
и диссипативных систем. Для диссипативных систем вводится понятие аттрактора, обсуждаются свойства аттракторов стохастических систем. И злагаются краткие сведения из теории фракталов, дается понятие обобщ енной размерности и описаны алгоритмы определения размерности аттракторов стохастических систем. Даны основы теории бифуркаций, рассмотрены некоторые методы исследования перехода к хаосу и характреистики
динамических систем при периодическом и хаотическом поведении (сечения П уанкаре, показатели Ляпунова, энтропия Колмогорова, спектры Ф урье). Описаны и обсуждены основные сценарии перехода от порядка к хаосу: сценарий Ландау, сценарий Рю эля и Таккенса, субгармонический каскад. В заклю чение главы рассматриваются примеры гидродинамических
систем, демонстрирую щ их хаотическое поведение. П роведен подробный
анализ поведения модели Лоренца, уравнения которой выведены в первой
главе. Рассмотрена также простейш ая модель генерации магнитного поля
Земли (динамо Рикитаки), воспроизводящ ая эффект случайных перебросов
направления магнитного поля. П оказаны и обсуждены также результаты
6
экспериментального наблю дения хаотизации конвективного течения в
замкнутой полости.
В третьей главе начинается знакомство с методами описания развитой турбулентности, а именно, с исторически первым и наиболее развитым
подходом к описанию турбулентных потоков. Это подход Рейнольдса и
выросш ие из него многочисленные полуэмпирические модели турбулентности. Н ачинается глава с определения статистических моментов случайных полей, характеризую щ их турбулентный поток. Далее дан вывод уравнения Рейнольдса для средних полей и обсуждаются вопросы, связанные с
появлением в уравнениях тензора напряжений Рейнольдса. П оказано, как
получается цепочка уравнений Ф ридмана-Келлера и формулируется проблема замыкания. Разговор о путях реш ения этой проблемы начинается с
описания гипотезы Буссинеска для тензора напряжений, определения понятия турбулентной вязкости, описания и обсуждения модели пути смеш ения
П рандтля. В последую щ их параграфах рассмотрены более сложные модели: модели переноса турбулентной вязкости и двухпараметрическиемодели
типа k - e модели. П олуэмпирическим моделям в предлагаемом курсе лекций уделено сравнительно скромное место по двум причинам. Во-первых,
именно этот подход наиболее полно освещ ен в литературе и может быть
свободно изучен по учебникам. Во-вторых, основной целью данного курса
является знакомство с методами изучения свойств мелкомасш табной турбулентности (однородной изотропной турбулентности), которая как раз и
остается за полем зрения полуэмпирических моделей. П оэтому описание
этих подходов необходимо только для общ его знакомства с идеологией метода, даю щ его возможность ссылаться на него в дальнейш ем и проводить
необходимыесравнения.
7
1 ОСН О ВЫ
1.1 Уравнения движения жидкости
Гидродинамика - это раздел механики сплош ных сред, описываю щ ий
движение жидкостей и газов в рамках модели сплош ной среды. П оследнее
означает, что рассматриваю тся масш табы l > > l , где l - длина свободного
пробега молекул.
Рассматривается физически
бесконечно малый объем, и вводятся хаr
рактеристики среды: скорость v и две термодинамические величины: давление P и плотность r .
1.1.1 Уравнение непрерывности
Законы движения выводятся иззаконов сохранения. Сначала используется закон сохранения вещ ества. В пространстве фиксируется некоторый
объем V, ограниченный поверхностью S , масса которого равна
m = тrdV .
V
И зменение массы этого объема есть
¶
¶
m = тrdV ,
¶t
¶t V
а вытекаю щ ий из объема поток жидкости
тr
v n dS .
S
Если за положительное направление принять направление движения
израссматриваемого объема, то условие сохранения массы можно записать
в виде
¶
rdV = - тrv n dS .
¶t Vт
S
8
П равая часть равенства преобразуется по теореме ОстроградскогоГаусса
r
тrv dS = тdiv( rv )dV .
n
S
V
Тогда
й¶r
ткл ¶t +
rщ
div ( rv )ъdV = 0 ,
ы
а так как равенство должно быть справедливо для любого объема, то
подынтегральноевыражение должно удовлетворять уравнению
¶r
r
+ div( rv )= 0 ,
¶t
(1.1)
которое называют уравнением непрерывности (неразрывности). Для
несжимаемой жидкости плотность есть величина постоянная ( r = const ) и
уравнение (1.1) упрощ ается:
.
r
div(v )= 0
(1.2)
Важно отметить, что уравнение неразрывности справедливо и для
идеальной,
и
для
реальной
жидкости.
1.1.2 И деальная жидкость
Уравнения для скорости выведем сначала для идеальной жидкости.
И деальная жидкость- это жидкость без вязкости и теплопроводности.
Закон сохранения импульса для движущ егося жидкого объема есть
d
dt
(тrvrdv )= е F ,
i
где в правой части стоит сумма всех сил, действую щ их на выделенный
объем. Ограничиваясь рассмотрением силы тяжести и сил давления, запиш ем
d
r
r
rv dV = тrgdV +
т
dt V
V
т(- P )dS .
S
9
Учитывая, что
dr
тdt dV є 0
(интеграл берется по жидкой частице, то
есть по заданному количеству жидкости, а не по заданному объему), можно
переписать уравнение в виде
d r
r
тr dt (v )= т(rg - С P )dV
V
и, снова исходя из произвольного выбора объема частицы, перейти к
дифференциальной форме
r
dv r С P
=g.
dt
r
(1.3)
r
dv
Входящ ая в уравнение производная
dt
это субстанциональная
производная, которая описывает изменение скорости жидкой частицы.
Рассмотрение движения отдельных жидких частиц называется подходом
Лагранжа к описанию движения жидкости. В больш инстве случаев предпочтительным является подход Эйлера, который заклю чается в описании
характеристик жидкости в заданной точке. Чтобы получить уравнение
движения в форме Эйлера, нужно получить связь между субстанциональной и локальной производными. Запиш ем приращ ение скорости
r
r ¶v
¶v
¶v
¶v
dv =
dt +
dx +
dy +
dz
¶t
¶x
¶y
¶z
и получим изнего связь субстанциональной (полной) производной по
времени с частной производную скорости по времени (изменение скорости
в заданной точке)
r
r
r
dv ¶v ¶v dx ¶v dy ¶v dz ¶v
¶v
¶v
¶v
=
+
+
+
=
+ vx
+ vy
+ vz
,
dt
¶t ¶x dt ¶y dt ¶z dt dt
¶x
¶y
¶z
или
r
r
dv ¶v r r
=
+ (v С )v .
dt
¶t
(1.4)
И спользуя полученное соотнош ение, приходим к уравнению Эйлера,
полученному
им ещ е в 1755 г.:
r
r
¶v r r
1
+ (v С )v = - С P + g .
¶t
r
(1.5)
10
Гидростатическое приближение получается при условии отсутствия
движения, то есть равенства нулю скорости и производной по времени:
r
¶
= 0 и v = 0.
¶t
Таким образом,
-
r
1
Сp+ g = 0,
r
(1.6)
r
илиС p = rg . Учитывая, что сила тяжести направлена вертикально вниз
r
r
и считая, что по вертикали направлена координата z , т.е. g = - ge z , получим
¶P
= - rg ,
¶z
а
P = P0 - rgz .
Запиш ем теперь поток импульса в тензорных обозначениях. Отметим,
что в дальнейш ем мы иногда производную по времени будем обозначать
как ¶t .
¶t ( rvi ) = r¶t vi + ¶t rvi
Уравнение непрерывности перепиш ем в виде
¶t r -
¶( rv k )
= 0,
¶x k
а уравнение Эйлера (1.5) в виде
¶t vi = - v k
¶vi 1 ¶P
.
¶x л r ¶xi
П одставим две последние формулы в выражение для изменения импульса:
¶t ( rvi )= - rv k
= - dik
¶vi ¶P
¶( rv k )
¶
¶P
( rvi vk )=
- vi
=¶x k ¶xi
¶x k
¶xi ¶x k
¶P
¶
( rvi v k )= - ¶ (dik P + rvi vk )
¶x k
¶x k ¶x k
и введем тензор плотности потока импульса, описываю щ ий перенос
i-ой компоненты импульса через площ адку, перпендикулярную k-ой оси
11
Х
ik
= dik P + rvi v k .
(1.7)
Тогда уравнение для изменения импульсазапиш ется в виде
¶ Х ik
¶
( rv i )= ¶x k
¶t
,
(1.8)
а для конечного объема
¶
т¶t rv dV = тi
V
V
¶ Х ik
dV = - тХ
¶x k
S
ik
dS k .
Н ами не использован закон сохранения энергии. Н апомним, что рассматривается идеальная жидкость, а это означает, что в жидкости отсутствую т теплообмен и трение. В таком случае движение адиабатично в каждой жидкой частице. Следовательно, закон сохранения энергии выливается
в утверждение, что энтропия каждого жидкого элемента остается постоянной
dS
= 0.
dt
П ереходя от переменных Лагранжа к переменным Эйлера, получаем
¶S
+ (vС )S = 0 .
¶t
(1.9)
1.1.3 Реальная жидкость
Реальная жидкость - это жидкость с вязкостью (внутренним трением)
и теплопроводностью . Н ачнем с рассмотрения уравнений движения для
изотермической жидкости и для начала ещ е раз напомним, что уравнение
непрерывности (1.1) справедливо и для реальной жидкости, так как его вывод основывался только на законе сохранения вещ ества.Далее воспользуемся уравнением Эйлера, записанным в формезакона для переноса импульса (1.7)-(1.8), и попытаемся дописать в него слагаемые, отвечаю щ ие за перенос импульса в результате действия вязких сил
¶
(rv i )= - ¶ (rv iv k + pdik + поток импульса из - за вязкости)=
¶x k
¶t
=-
¶
(rv iv k + pdik - s ikў )
¶x k
12
где величина s ik = pdik - s ikў называется тензором напряжений, а s ikў тензором вязких напряжений.
Тензор вязких напряжений s ikў должен характеризовать неоднородности поля скорости, которые можно описать производными поля скорости
¶vi
¶2 v i
,
, ...... .
¶x k
¶xi ¶x k
Требуется угадать форму зависимости тензора вязких напряжений от
этих производных. Н а этом этапе делается самое важное ограничение на
пути получения уравнения движения. Оно состоит в том, что учитываю тся
только линейные комбинации первых производных поля скорости. Кажется естественным, что однородное поле скорости не приводит к появлению
вязких напряжений. Н ужно, однако, учесть, что есть специальный случай,
когда поле скорости неоднородно, а вязкие напряжения возникать не
должны. Это случай твердотельного вращ ения жидкости.
Сущ ествуют только две линейные комбинации первых производных,
удовлетворяющ ие этому требованию. Это
ж ¶vi ¶v k ц
з
з¶x + ¶x ч
ч
i ш
и k
и
r ¶v
div v = k
¶x k
(здесь и далее подразумевается суммирование по повторяю щ имся индексам).1
Общ ий вид тензора второго ранга, удовлетворяю щ его поставленным
условиям, есть
ж ¶vi ¶v k ц
r
s ikў = aз
з¶x + ¶x ч
ч + bdik div v .
i ш
и k
П ринята несколько иная форма записи
Убедимся, что эти две комбинации равны нулю при твердотельном вращ ении жидкости.
r r r
v = W ґr
vx = W y z - W z y
1
vy = W z y - W x z
vz = W x y - W y x
¶v x ¶v y
¶v x ¶v y ¶v z
+
= - W z + W z = 0 и т.д.
=
=
=0
¶y
¶x
¶x
¶y
¶z
13
ж ¶vi ¶v k 2
rц
r
s ikў = h з
ч+ xdik div v ,
з¶x + ¶x - 3 dik div v ч
i
ш
и k
(1.10)
удобная тем, что сумма диагональных членов в скобке равна нулю. В
выражении присутствуют два коэффициента:
h -сдвиговая вязкость
x -объемная (вторая) вязкость.
Таким образом, уравнение движения приобретает вид
ж¶vi
¶P
¶ жж ¶vi ¶vk ц 2
¶vi ц
rц ¶ ж ¶vk ц
з
ч
rз
ч.
ч - 3 dik div v ч + ¶x з
зx ¶x ч
ч = - ¶x + ¶x h зз
з¶x + ¶x ч
з ¶t + vk ¶x ч
k ш
i
k
k
i
i
k
ш
ш
и
и
и
и
ш
(1.11)
Коэффициенты вязкости зависят от температуры и давления и не являю тся постоянными вдоль жидкости. Однако, во многих случаях можно
считать эту зависимость слабой и вынося коэффициенты вязкости за операторы дифференцирования, прийти к виду
r
r ж hц
r
й¶v r rщ
r к + (v С )v ъ = - С p + hDv + зx + чgrad div v ,
3ш
и
л¶t
ы
(1.12)
который и принято называть уравнением Н авье-Стокса. Важным частным случаем является случай несжимаемой жидкости. Тогда уравнения
движения (1.1),(1.12) записываются в виде
r
r
¶v r r
1
+ (v С )v = - С P + nDv
r
¶t
r
div v = 0
(1.13)
где n = h / r - коэффициент кинематической вязкости.
Для реш ения конкретной задачи уравнения должны быть дополнены
граничными условиями (например, условие прилипания на твердой границе или условие отсутствия касательных напряжений на свободной границе).
Основные проблемы реш ения уравнений Н авье-Стокса связаны с нелинейным членом. И звестно небольш ое число задач, в которых этот член
обращ ается в нуль и задачи приводят к точным реш ениям. П риведем только два хорош о известных примера таких задач.
Течение Куэтта. Рассматривается течение несжимаемой жидкости
в горизонтальном слое толщ иной d ,
нижняя граница которого неподвижна, а верхняя движется с постоянной
Рис. 1.1.
14
горизонтальной скоростью v 0 , направленной вдоль оси x . Ось z направлена вертикально вверх. И щ ется стационарное реш ение, то есть производная
по времени равна нулю . Считается также, что задача плоская, то есть нет
зависимости от координаты y и нет соответствую щ ей компоненты скорости ( v y = 0 ). Более того, течение горизонтально и v z = 0 . Отсутствует также
горизонтальный градиент давления.
И з уравнения непрерывности немедленно следует, что оставш аяся
компонента скорости не можетзависеть от координаты x :
¶v x
= 0.
¶x
Следовательно, v x = f (z ) , и нелинейный член исчезает
vx
¶v x
¶v
¶v
+ v y x + vz x = 0 .
¶x
¶y
¶z
В результате, от уравнения Н авье-Стокса остается
¶2 v x
= 0 или v x = az + b .
¶z 2
П остоянные интегрирования находятся изграничных условий
v x = 0 при z = 0
v x = v0
при z = d
и получаетсярезультат
v x = v0
z
.
d
П ри этом, отлична от нуля только одна компонента тензора вязких
напряжений
s xz = h
¶v x v0
= h,
¶z
d
с которой просто связана сила, действую щ ая на площ адку поверхности площ адью S
F=
v0hS
.
d
15
ТечениеП уазейля.
Второй
хорош о известный пример задачи о
течении вязкой жидкости, имею щ ей
точное реш ение, является задача П уазейля о течении жидкости в слое с
твердыми границами (или трубе) под
действием приложенной к краям разности давления. Рассмотрим сначала
плоский горизонтальный слой толщ иной 2d и длиной L , на концах которого задано давление P1 и P2 , соответстРис. 1.2.
венно. Как и в предыдущ ей задаче,
ищ ем стационарное реш ение ( ¶t = 0 )
только для компоненты скорости v x ( v y = v z = 0 ) и по тем же причинам
¶x = ¶y = 0 . В этом случае снова исчезает нелинейный член, так как возникаю щ ий градиент скорости направлен перпендикулярно самой скорости.
Тогда уравнение Н авье-Стокса принимаетвид
-
¶2 v
1 ¶P
+ n 2x = 0 ,
r ¶x
¶z
аего реш ение, с учетом граничных условий ( v x = 0 при z = ±d ) есть
vx =
P1 - P2 2
( z - d2 ).
2hL
Для цилиндрической трубы радиуса R задача реш ается аналогично. В
этом случае оператор Лапласа нужно записать в цилиндрической системе
координат
1 d ж dv ц P1 - P2
зr ч =
Lh
r dr и dr ш
и его реш ение примет вид
v=
P1 - P2 2
r + C lnr + B .
4hL
П остоянная интегрирования C = 0 , так как при r = 0 значение скорости
должно быть конечно. Определив вторую константу из условия прилипания на стенкетрубы, получим
16
v=
P1 - P2 2
R - r2 .
4hL
(
)
Важной характеристикой является расход жидкости, протекаю щ ей
через трубу. Для него имеем
R
Q = 2p тrrvdr =
0
pr ( P1 - P2 ) 4
R .
8hL
1.1.4 Число Рейнольдса
П олученные уравнения движения несжимаемой жидкости приведем к
безразмерному виду. И меем (1.13),
r
r
¶v r r
СP
+ (v С )v = + nDv ,
¶t
r
r
div v = 0.
В уравнения входят следую щ ие физические величины: время t , расстояние l , скорость v , плотность r , давление P и вязкость n . Если мы
принимаем систему единиц СИ, то каждая изэтих величин будет иметь следую щ ую размерность:
[t ] = c ;
[l ] = м;
[v ] = м/с;
[ r ] = кг/м3 ;
[P ] = кг/м Чс2 ;
[n ] = м 2 /с.
И дея обезразмеривания состоит в том, чтобы измерять все величины
в единицах, являю щ ихся характерными параметрами конкретной задачи.
Так, например, в качестве единицы измерения длины можно выбрать некий
характерный размер L (это может быть толщ ина слоя жидкости, диаметр
трубы, размер обтекаемого тела и т.д.), за единицу измерения скорости характерную скорость V (скорость верхней пластины в течении Куэтта,
скорость на оси трубы в течениеП уазейля, скорость набегаю щ его потока в
задачах об обтекании тела и т.д.). Единица измерения времени выражается
через две введенные величины и есть L / V , а единицей давления может служить величина rV 2 .
Безразмерные величины (обозначим их буквами с тильдами)
будут связаны со старыми, размерными величинами как
~r r
~
~
~
v = v /V , ~
xi = xi / L , ~
t = tV / L , P = P / rV 2 , С = LС , D = L2 D.
17
П одставляя эти соотнош ения в уравнения движения, получим
~r
V 2 ¶v V 2 ~r ~ ~r
V2 ~~
V ~ ~r
+
v
С
v
=
С P + n 2 Dv
~
L ¶t
L
L
L
~r
div v = 0,
( )
а сокращ ая подобные множители и опуская тильды, приходим к
уравнениям
r
¶v r r
1 r
+ (v С )v = - С P + Dv
¶t
R
r
div v = 0,
где безразмерная величина R =
(1.14)
VL
называется числом Рейнольдса.
n
Это число характеризует отнош ение инерционных сил к вязким (нелинейного члена к вязкому) и именно оно является критерием, определяю щ им
этапы перехода от ламинарных течений к турбулентным.
Важно подчеркнуть, что приведенный способ обезразмеривания
уравнений не является единственно возможным. Н апример, в качестве единицы времени можно взять величину L2 / n , характеризую щ ую время вязкой
диссипации, а в качестве единицы скорости - величину n / L . П ереходя к
безразмерным переменным, в этом случае получим уравнение
r
r
¶v r r
+ (v С )v = - С P + Dv ,
¶t
не содержащ ее каких-либо параметров. Н е следует думать, что таким
образом мы получили уравнение, лиш енное параметра. В действительности, роль числа Рейнольдса выполняет теперь безразмерная скорость. Если
при первом способе обезразмеривания безразмерная скорость по определению лежала в интервале [0,1] (или вблизи него), то при втором способе
единичной является скорость вязкого переноса, а безразмерная скорость
можетдостигать величин порядка
v
VL
v~ =
»
,
n/L n
то есть являетсяаналогом числа Рейнольдса.
С числом Рейнольдса тесно связан вопрос о подобии различных течений, то есть вопрос о том, каким критериям должна удовлетворять модель исследуемого течения. П усть рассматривается определенный тип тече-
18
ний жидкости (например, течение по трубам или обтекание тел определенной формы). Очевидно, что для моделирования движения нужно в первую
очередь обеспечить геометрическое подобие. Тогда геометрические свойства задачи определяются одним линейным размером L .
И з параметров, характеризую щ их жидкость, в уравнения входит
r
только кинематическая вязкость n (поля скорости v и давления, отнесенного к плотности, P / r являются неизвестными функциями, которые необходимо найти). Если рассматривается обтекание тела потоком, то характеристикой течения в целом является скорость потока (на бесконечности) V .
М ы видим, что в рамках заданного типа движений реш ение определяется
тремя параметрами: n ,V , L . И з этих трех размерных величин можно составить только одну безразмерную комбинацию , а именно, введенное выше
число Рейнольдса. И скомые поля (опять же, для заданного типа течений)
должны будут выражаться зависимостями вида
r
r
жr ц
v = Vf 1 з , R ч ,
иL ш
r
жr ц
P = rV f 2 з , R ч .
иL ш
2
Суть закона подобия, сформулированного Рейнольдсом в 1883
году, состоит в том, что течения одного типа с равным числом Рейнольдса
подобны. П одобие двух течений состоит в том, что все поля могут быть
получены друг из друга простым масш табным преобразованием координат
и скорости.
Если в задаче появляется дополнительный параметр, то из
имею щ ихся четырех величин можно составить два независимых безразмерных комплекса и для обеспечения подобия задач потребуется обеспечить
равенство обоих безразмерных параметров. Так, если в рассматриваемом
течении сущ ественно влияние сил тяжести, то в качестве дополнительного
размерного параметра в задачу входит ускорение силы тяжести g . Тогда
новым безразмерным параметром можетслужить число Ф руда
F=
V2
,
Lg
являю щ ееся мерой отнош ения кинетической энергии движущ ейся
жидкости к потенциальной.
1.1.5 Течение в диффузоре
М ы рассмотрели выше два простейш их примера точных реш ений уравнений Н авье-Стокса. И звестно ещ е несколько задач, для которых
19
найдены точные реш ения. Это, например, задача о затопленной струе, задача о течении
вблизи вращ аю щ егося диска, течение в диффузоре и некоторые другие. Н е воспроизводя реш ения задачи, остановимся на течение жидкости в плоском диффузоре (задача Гамеля,
1917г.).
П лоский диффузор образован двумя полу-плоскостями, выходящ ими из начала
координат под углом a (рис.1.3). В начале коРис. 1.3.
ординат находится источник жидкости мощ ностью Q . Если Q < 0 , то источник становится
стоком, а устройство называетсяконфузором.
Реш ение ищ ется в цилиндрической системе координат (r ,j , z )
для чисто радиального течения (vj = v z = 0; v r = v(r , j ) ) . Уравнение непрерывности, записанное в цилиндрических координатах
1 ¶(rv r ) 1 ¶vj ¶v z
+
+
= 0,
r ¶r
r ¶j
¶z
показывает, что (rv ) не зависит от радиуса и может быть только
функцией угла j . Реш ение поэтому ищ етсядля автомодельной переменной
u (j ) =
1
rv .
6n
Вид реш ения, получаю щ егося для конфузора при малых и больш их
числах Рейнольдса, иллю стрируетрис.1.4.
И нтересной особенностью задачи Гамеля является то, что для конфузора (втекание жидкости, Q < 0 ) реш ение сущ ествует для лю бых значений
Рис. 1.4.
20
числа Рейнольдса, которое определяется через расход и есть
R=
|Q|
,
rn
а для диффузора ( Q > 0 ) симметричное расходящ ееся течение сущ ествует только при ограниченных значениях числа Рейнольдса R < Rmax и ограниченных значениях угла раствора
a < a max . П редельные параметры связаны простым соотнош ением
Rmax
ц
жp 2
= 6з - a ч,
ч
зa
ш
и
которое определяет область сущ ествования симметричных реш ений
Рис. 1.5.
на плоскости (R,a ) (см. рис.1.5). П ри
R > Rmax сущ ествуют только несимметричные реш ения, в которых имею тся области возвратных течений. П римеры профилей скорости, соответствую щ их таким реш ениям, приведены на
рис.1.6. Важно отметить, что реш ение в конфузоре при R ® Ґ стремится к
реш ению для идеальной жидкости (столбообразное течение с проскальзыванием на границе), а в диффузоре предельного перехода нет : при
R ® Ґ число перегибов в реш ении неограниченно возрастает.
Задача о диффузоре интересна тем, что является примером задачи, в
которой сущ ествует граничное значение числа Рейнольдса, при превышении которого реш ение данного вида не сущ ествует. Н еследует путать этот
случай с ситуацией, когда реш ение в принципесущ ествует, но не реализуетсяв силу возникаю щ ей неустойчивости. Об этом пойдетречь далее.
Рис. 1.6.
21
1.2 Устойчивость течений
Вопрос об устойчивости того или иного состояния (реш ения, режима)
возникает в самых разных задачах. Достаточно вспомнить простейш ий
пример об устойчивости ш арика, лежащ его на различных поверхностях
(рис.1.7). В первом случае положение ш арика абсолютно устойчиво, то есть
при любом конечном воздействии ш арик по окончании действия возмущ аю щ ей силы возвращ ается в исходное состояние. Во втором случае положение ш арика абсолю тно неустойчиво - любое, сколь угодно малое возмущ ение, безвозвратно уводит его из начального положения. Третий случай
иллю стрирует пример состояния, устойчивого по отнош ению к малым возмущ ениям, но наруш аю щ егося, если возмущ ения превышаю т критическую величину.
Н ас интересует вопрос
об устойчивости стационарных течений. Для конкретРис. 1.7.
ности будем говорить о течении П уазейля. Возмущ ения в реальных течениях сущ ествую т всегда. И х источником служат ш ероховатости стенок, входные участки (бесконечных труб нет), просто флуктуации характеристик самой жидкости и т.д. Н ужно ответить на вопрос о
том, какое возмущ ение является самым опасным и где та граница, при превышении которой это возмущ ение приведет к разруш ению сущ ествую щ его
течения.
И так, имеем течение несжимаемой жидкости, для которой запиш ем
уравнения Н авье-Стокса в безразмерной форме(1.14)
r
¶v r r
1 r
+ (v С )v = - С P + Dv ,
¶t
R
r
div v = 0.
Стационарное реш ение задачи (имеем в виду течение П уазейля, хотя
до определенного этапа все рассуждения не зависят от конкретного вида
реш ения) обозначим как v0 , P0 . Это реш ение, в свою очередь, удовлетворяетуравнениям
(vr0 С )vr0 = - С P0 +
r
div v0 = 0.
1 r
Dv 0 ,
R
(1.15)
22
П оля скорости и давления представим в виде сумм стационарных реш ений и возмущ ений
r
r
r
v ( x, y, z, t ) = v 0 ( z ) + v ў( x, y, z, t ),
P ( x, y, z, t ) = P0 ( z ) + P ў( x, y, z, t ).
(1.16)
Отметим, что в отличие от исследуемого стационарного реш ения,
слагаемые со ш трихами описываю т возмущ ения, которые могут зависеть
от времени и от всех координат. Введенные разложения подставляю тся в
исходныеуравнения
r
r r r r
r r
¶v ў r r
1 r 1 r
+ (v0 С )v0 + (v0 С )v ў+ (v ўС )v 0 + (v ўС )v ў = - С P0 - С P ў+ Dv0 + Dv ў
R
R
¶t
r
r
div v0 + div v ў= 0
(1.17)
и, после вычитания из них уравнений для стационарных реш ений
(1.15), получаем
r
r r
¶v ў r r r r
1 r
+ (v0 С )v ў+ (v ўС )v0 + (v ўС )v ў = - С P ў+ Dv ў,
¶t
R
r
div v ў= 0.
(1.18)
Н аибольш ие трудности в реш ении этих уравнений представляет неr r
линейное по искомым возмущ ениям слагаемое (v ўС )v ў. Следую щ ий, принципиальный ш аг состоит в том, что это слагаемое отбрасывается. Тем самым мы ограничиваем себя рамками линейной теории устойчивости, рассматриваю щ ей эволю цию малых возмущ ений. Это значит, что
r
r
| v ў| << | v0 | .
Линейная теория работает только вблизи порога возникновения неустойчивости. П о прохождению порога, возмущ ения нарастаю т и линейные уравнения перестаю т работать. Тем не менее, поставленная задача при
этом может считаться выполненной, так как требовалось указать именно
сам порог и наиболее опасные возмущ ения, которые начинаю т нарастать в
первую очередь.
Отказавш ись от написания ш трихов, мы придем к системе уравнений,
которую необходимо дополнить граничными условиями для возмущ ений.
Н апример, можно предположить, что на границах возмущ ения равны нулю .
23
r
¶v r r r r
1 r
+ (v 0 С )v + (v С )v 0 = - С P + Dv ,
¶t
R
r
=
div v 0,
r
v Г = 0.
(1.19)
Далее делаю т ещ е ряд сущ ественных упрощ ений. П ервое состоит в
том, что рассматриваются только плоские возмущ ения. Этот ш аг оправдывается теоремой Скваера, которая утверждает, что самыми опасными являю тся именно плоские возмущ ения. Такое предположение означает, что
r
v = (v x ,0,v z ) и
¶
= 0.
¶y
С учетом того, что
(vr0 С )= v0
¶
¶x
(vrС )= v x
и
¶
¶
+ vz
,
¶x
¶z
уравнения движения для оставш ихся двух компонент запиш утся в виде
¶v x
+
¶t
¶v z
+
¶t
¶v x
+
¶x
¶v
¶v x
¶P 1
+ Dv x ,
+ vz 0 = ¶z
¶x R
¶x
¶v
¶P 1
v0 z = + Dv z ,
¶x
¶z R
¶v z
= 0.
¶z
v0
Следую щ ий ш аг состоит в том, что вводится функция тока y , связанная с
компонентами вектора скорости: v x = -
¶y
¶z
vz =
¶y
. Введение функции тока
¶x
позволяет уменьш ить число переменных. П латой за это является повышение порядка дифференциальных уравнений, которые принимаю т вид:
¶ ¶y
¶2y
¶y ¶v0
¶P 1 ¶y
,
- v0
+
=- D
¶t ¶z
¶x¶z ¶x ¶z
¶x R ¶z
¶ ¶y
¶2y
¶P 1 ¶y
+ v0 2 = + D
,
¶t ¶x
¶z R ¶x
¶x
¶2y
¶2y
+
є 0.
¶x¶z ¶z¶x
-
П оследнее уравнение (это уравнение непрерывности) выполняется тождественно. Это не удивительно, так как функция тока вводится именно для
24
несжимаемой жидкости. Следую щ ий ш аг также является общ епринятым для того, чтобы избавиться от давления и получить одно уравнение для
функции тока, необходимо второе уравнение продифференцировать по координате x и вычесть из него первое, продифференцированное по координате z . Результирую щ ееуравнение есть
¶ ж¶2y
¶2y
з
+
2
¶t з
¶x 2
и ¶z
ц
¶ ¶2y
¶2y ¶v0 ¶2y ¶v0 ¶y ¶2 v0 ¶2y ¶v 0
¶ ¶2y
ч
+
+
+
+
=
v
v
0
ч 0 ¶x ¶z 2 ¶x¶z ¶z ¶x¶z ¶z
¶x ¶z 2 ¶x¶z ¶x
¶x ¶x 2
ш
¶2 P ¶2 P 1 ж¶2y
¶2y
+
+ Dз
+
2
¶x¶z ¶x¶z R з
¶z 2
и ¶x
ц
ч
ч.
ш
Сокращ ая подобные члены и учитывая, что ¶v0 ¶x = 0 , приходим к уравнению
¶
¶
1
І ¶y
Dy + v0 Dy - v0
= DDy ,
¶t
¶x
¶x R
(1.20)
которое дополняетсяграничными условиями для функции тока:
при z = ±1 :
¶y
¶y
=
=0
¶x
¶z
Н апомним, что функция тока введена для возмущ ений поля скорости, возr
никаю щ их на фоне стационарного течения v0 . Ш трихами обозначено дифференцирование по вертикальной координате z .
П олученное уравнение (1.20) можно реш ать численно, задавая различные начальные возмущ ения и наблюдая за их эволю цией при различных числах Рейнольдса. Этот путь не снимает, однако, вопроса о выборе
вида возмущ ений. Следуя обычному для теории устойчивости способу, будем рассматривать нормальные возмущ ения, то есть возмущ ения вида
y ( x, z, t )= j ( z )e i (w t - kx ).
Рис. 1.8.
(1.21)
П ри это фактически мы провели
разделение переменных, вклю чив
зависимость от вертикальной координаты z только в амплитуду
возмущ ений j . Зависимость от
продольной координаты и времени принята в виде гармонических
волн, распространяю щ ихся вдоль
оси x ( w - частота, k - волновое
25
число). Частота является величиной комплексной: w = a + ib , что позволяет
переписать выражение для нормальных возмущ ений в виде
y ( x, z, t )j ( z )e i (w t - kx ) = j ( z )e - bt e i ( at - kx ) .
Характер эволю ции колебаний во времени определяется мнимой частью частоты: если b > 0 , то возмущ ения убывают со временем, а если b < 0 ,
то возмущ ения нарастаю т (см. рис.1.8). И менно знак величины b и интересен с точки зрения вопроса об устойчивости течения. Требуется узнать, при
каком значении числа Рейнольдса появляется реш ение с отрицательным b
и какое волновое число k соответствует этому реш ению .
Возмущ ения в нормальной форме подставляются теперь в уравнение
для функции тока. Соответствую щ ие производные определяю тся формулами:
¶y
= iwy ,
¶t
¶y
= - iky ,
¶x
¶y
= j ўеi (w t - kx ) ,
¶z
¶2y
¶2y
Dy = 2 +
= j ўў- k 2j e i (w t - kx ),
2
¶x
¶z
(
[
)
)]
(
(
)
(
)
DDy = j IV - k 2j ўў- k 2 j ўў- k 2j e i (w t - kx ) = j IV - 2k 2j ўў+ k 4j e i (w t - kx ).
П осле подстановки получаем
(
)
(
)
1 IV
й
щ i (w t - kx )
І
2
2
2
4
кiw j ўў- k j - ikv0 j ўў- k j + ikv0 j = R j - 2k j ўў+ k j ъe
л
ы
(iw - ikv0 )j ўў- k 2j + ikv0 Іj = 1 j IV - 2k 2j ўў+ k 4j ,
R
(
)
(
)
а после деления на ik и добавления граничных условий приходим к окончательной форме уравнения, называемого уравнением Орра-Зоммерфельда
(1937г.):
i
wц
ж
І
2
(
j IV - 2k 2j ўў+ k 4j ),
зv 0 - ч(j ўў- k j )- v 0 j =
kш
kR
и
j z = ±1 = 0 j ўz = ±1 = 0.
(1.22)
Задача остается чрезвычайно сложной и впервые для плоского слоя
была реш ена только в 1945 г. Линем. П оучительна история реш ения этого
уравнения. П ервые подходы были связаны с попытками реш ать уравнение
Орра-Зоммерфельда с отброш енной правой частью. Соответствую щ ее
26
уравнение называют уравнением Релея. Отметим, что отбрасывая члены с
четвертой производной j IV , мы лиш аемся возможности использовать все
граничные условия и можем требовать обращ ения в нуль только нормальной компоненты скорости (этому соответствует условие ¶y ¶x = 0 и j = 0 ).
Отбрасывание правой части мотивировалось тем, что она описывает действие вязкости, а вязкость, казалось, должна играть стабилизирую щ ую
роль. Результат реш ения уравнения Релея состоял в том, что оно оказывалось абсолю тно устойчивым.
Линь показал, что фазовая скорость возмущ ений vф = w / k меньш е максимальной скорости потока в центре слоя. Точки, в которых фазовая скорость возмущ ений совпадает со скоростью основного течения, являю тся критическими и именно вблизи этих точек начинается
нарастание возмущ ений. Основной результат
исследования уравнения Орра-Зоммерфельда
качественно иллю стрируется рисунком 1.9, на
Рис. 1.9.
котором представлена так называемая нейтральная кривая, нарисованная на плоскости
k - R . Область неустойчивости заш трихована. Критические параметры отмечены на рисунке звездочками. Н аименьш ее значение числа Рейнольдса,
при котором начинается рост возмущ ений R * = 5700 . Соответствую щ ее ему
критическое значение волнового числа k * » 1 . Это означает, что наиболее
опасными возмущ ениями являю тся возмущ ения с длиной волны, превыш аю щ ей толщ ину слоя приблизительно в 2p раз.
И нтересна ещ е одна особенность нейтральной кривой. П ри некоторых значениях волнового числа в область неустойчивости можно попасть
и двигаясь от больш их чисел Рейнольдса к малым. Это означает, что вязкость можетиграть и дестабилизирую щ ую роль.
1.3 Свободная конвекция несжимаемой жидкости
П од свободной конвекцией понимаю т движения жидкости, возникаю щ ие за счет сил Архимеда при наличии неоднородности плотности
жидкости в поле массовых сил. В основном будем рассматривать термогравитационную конвекцию, т.е. случай, когда неоднородности жидкости
связаны с ее неравномерным нагревом и течение возникает в поле силы
тяжести. П ри этом будем иметь в виду жидкости, плотность которых падает с ростом температуры, т.е. ¶r ¶T < 0 (напомним, что аномальное поведе-
27
ние дает вода в интервале от 0 до 4о С). Считаем, что неоднородность температуры являетсяединственным источником движения и что
Dr << r ,
т.е. рассматривается слабая конвекция. В уравнении движения появляетсяслагаемое, описываю щ еедействие силы тяжести
r
r
r
й¶v r rщ
r к + (v С )v ъ = - С P + hDv + rg
л¶t
ы
и нужно учесть изменения плотности. П оследняя в общ ем случае есть
функция температуры и давления r = r (T , P ), а приращ ение плотности есть
ж¶r ц
ж ¶r ц
dr = з ч dT + з ч dP .
и¶P шT
и¶T шP
Далее делается важное ограничение, состоящ ее в том, что рассматривается
несжимаемая жидкость, означаю щ ее что вторым слагаемым в этом равенстве можно пренебречь. Таким образом, полагается, что плотность зависит
только от температуры: r = r (T ), а приращ ение плотности есть
dr =
1 ж ¶r ц
з чr 0 dT = - br 0 dT .
r 0 и¶T ш
Здесь b - коэффициент объемного расш ирения. Температуру жидкости представим в виде
T = T0 + T ў,
(1.23)
где T0 - средняя температура, а T ў- вариации температуры, малые в том
смысле, что вызываемые ими вариации плотности остаются малыми
( Dr << r ). П лотность представляется, соответственно, в виде r = r0 + rў(T ) ,
где r 0 - плотность жидкости при температуре T0 . И з сказанного выше следует, что
r ў= - r 0 bT ў
или
r = r 0 (1 - bT ў).
(1.24)
П ринятое ограничение слабой конвекции предполагает, что bT ў<< 1 .
Вспомним, что для воды b = 2 Ч10 - 4 , и следовательно приближение годится
28
практически для любых возможных разностей температуры. Для газов
b » 1 273 , что сущ ественно больш е, но также позволяет пользоваться принятыми ограничениями при достаточно больш их разностях температуры.
И зотермической жидкости с температурой T = T0 и соответствую щ ей
этой температуре плотностью r = r 0 отвечает гидростатическое давление
P0 , подчиняю щ ееся уравнению
r
С P0 = r0 g .
П оле давления, устанавливаю щ ееся при конвективном движении,
представим в видесуммы
P = P0 + Pў.
П одставляя в уравнения движения всевведенныеразложения, получаем
r
й¶v r rщ
( r0 + r ў)к + (v С )v ъ = - С P0 - С P ў+ hDvr + r0 gr + r ўgr ,
л¶t
ы
ў
r
r
¶r
+ r 0 div v + div r ўv = 0 .
¶t
Теперь нужно вычесть из первого уравнения уравнение гидростатики и
сделать самое важное допущ ение. Оно состоит в том, что добавкой к плотности r ў, возникаю щ ей за счет изменения температуры, пренебрегаю т всю ду, заисклю чением члена, описываю щ его силу Архимеда. Тогда
r
r
r
й¶v r rщ
r0 к + (v С )v ъ = - С P ў+ hDv - bT ўg .
л¶t
ы
Систему необходимо дополнить уравнением для температуры. Если
пренебречь нагревом жидкости за счет вязкой диссипации, то закон переноса удельной энергии записывается в виде
й¶S r
щ
rT к + (v С )S ъ = kDT ,
л ¶t
ы
где k - коэффициент теплопроводности, а энтропия S связана с температурой и давлением
ж ¶S ц
ж ¶S ц
S = S 0 + з ч T ў+ з ч P ў.
и¶P шT
и¶T шP
И спользуясоотнош ение
29
сp
ж ¶S ц
з ч =
и¶T шP T0
и считая третье слагаемое пренебрежимо малым (это логично сделать, так
как зависимостью плотности от давления уже пренебрегли), приходим к
соотнош ению
S = S0 +
cp
T0
T ў.
П одставляя в уравнение для энтропии и ограничиваясь членами, линейными по T ў, получаем
¶T ў r
k
+ (v С )T ў=
DT ў.
¶t
rc p
Далее, откажемся от написания ш трихов (не забывая при этом, что температура отсчитывается от среднего значения, а давление - от гидростатического давления) и запиш ем результат - систему уравнений для термогравитационной конвекции несжимаемой жидкости в приближении Буссинеска
r
r
r
¶v r r
СP
+ nDv + gbTe z ,
+ (v С )v = r0
¶t
r
¶T
+ (v С )T = cDT ,
¶t
r
div v = 0.
r
(1.25)
r
М ы учли, что g = - ge z и ввели коэффициент температуропроводности
c = k / rc p . Систему необходимо дополнить граничными условиями. Для
r
скорости можно принять, например, условия прилипания ( v | Г = 0 ), а для
температуры - либо задать ее распределение на границе ( T | Г = f1 ( Г) ), либо
теплопоток через границу
¶T
¶n
= f 2 ( Г) .
Г
Обсудим возможные способы представления уравнений свободной
конвекции в безразмерной форме. Особенностью конвективных задач является отсутствие заданной характерной скорости - скорость есть результат
приложенной (заданной) разности температуры. Возможный набор единиц
измерения есть: расстояния - характерный размер L , температуры - характерная разность температур q , скорости - величина n L , времени - L2 n и
давления - r 0n 2 L2 . П ереходя к безразмерным величинам, получаем систему
уравнений
30
r
r
r
¶v r r
+ (v С )v = - С P + Dv + GTe z ,
¶t
r
¶T
1
+ (v С )T = DT ,
s
¶t
r
div v = 0.
(1.26)
В уравнения входят два безразмерных комплекса: число Грассхофа
G=
и число П рандтля
gbJL3
n2
s=
n
.
c
Число Грассхофа характеризует отнош ение архимедовых сил к вязким и
свидетельствует о сильной зависимости конвективных механизмов от размера (в число Грассхофа размер входит в кубе). В отличие от числа Грассхофа, число П рандтля есть физический параметр жидкости, не зависящ ий
от конкретной задачи, и характеризую щ ий отнош ение коэффициентов кинематической вязкости и температуропроводности. П риведем несколько
типичных примеров значений числа П рандтля. Для газов число П рандтля
порядка единицы, у воды s » 7 , у ртути s » 10 - 2 , у глицерина - s » 10 3 . В
жидкостях с малым числом П рандтля теплопередача эффективней конвекции и наоборот, при высоких П рандтлях температура «вморожена» в жидкость и перенос тепла за счет конвекции становится более эффективен, чем
теплопередача.
Н аряду с двумя введенными безразмерными параметрами, в конвективных задачах часто используется число Релея, являю щ ееся произведением чисел П рандтля и Грассхофа
Ra = sG =
gbJL3
.
nc
Если за единицу скорости взять величину c L , оставив все остальные единицы измерения прежними, то мы придем к системе уравнений, содержащ ей число Релея
r
r
r
¶v r r
+ (v С )v = - С P + Dv + RaTe z ,
¶t
r
¶T
s
+ (v С )T = DT ,
¶t
r
div v = 0.
(1.27)
31
За единицу скорости можно выбрать и скорость, приобретаемую жидкой
частицей, перегретой на величину J относительно окружаю щ ей ее жидкости и разгоняю щ ейся на расстоянии L . И з условия rV 2 ~ r ўgL получаем
V ~ gbJL . П ринимая заединицу времени величину L / V , получаем
r
¶v r r
1 r r
+ (v С )v = - С P + Dv + Te z ,
¶t
R
r
¶T
1
+ (v С )T =
DT ,
¶t
sR
r
div v = 0.
(1.28)
В уравнениях появилось число Рейнольдса, что обусловлено введением характерной скорости. И спользуя выражение для введенной единицы скорости, просто получить связь появивш егося числа Рейнольдса с числом
Грассхофа
R=
VL
=
n
gbJL ЧL
= G.
n
1.4 Конвективная устойчивость
Рассмотрим вопрос о том, может ли жидкость оставаться неподвижной при наличии неоднородного распределения температуры. Чтобы убедиться, что равновесие неравномерно нагретой жидкости возможно, достаточно вспомнить ш кольный опыт по кипячению воды в наклоненной пробирке, на дне которой находитсялед, а нагреваетсятолько верхняя часть.
Н айдем необходимое условие механического равновесия жидкости
(при наличии неоднородности температуры). М еханическое равновесие
подразумевает отсутствие скоростей и стационарность:
r
v =0,
¶
= 0.
¶t
С учетом этих условий от уравнений Буссинеска остается
-
r
1
С P + gbTe = 0
r0
DT = 0.
32
Н а первое уравнение подействуем оператором rot. Так как
rot С = 0
r
r
r
rot (Te )= T rot e + С T ґ e = 0,
r
а rot e = 0 , то условие равновесия жидкости сводитсяк требованию
r
СT ґe = 0,
то есть градиент параллелен вертикальной оси и температура можетменяться только по вертикали: T = T ( z ). Это означает, что лю бой горизонтальный градиент температуры приводит к возникновению конвективного
движения.
И з второго уравнения
DT =
¶2 T
=0
¶z 2
следует, что температура может
быть только линейной функцией
высоты:
T = Az + B .
М ы не получили никакой
информации даже относительно
знака градиента температуры.
Опыт подсказывает, что устойчивым может быть нагрев сверху. Более точный ответ состоит в
том, что неустойчивость наступает при подогреве снизу после
превышения некоторого (совсем
небольш ого) критического градиента температуры. Н апример,
Рис. 1.10.
в горизонтальном слое с твердыми границами критическое
число Релея, при котором возникает конвекция, равно 1708. Оценим соответствую щ ую критическую разность температуры, имея в виду для определенности слой воды толщ иной h :
DT = R *
nc
1708 Ч10 - 6 Ч1.4 Ч10 - 7 10 - 7 град
»
» 3
.
м
gbh 3
10 Ч2 Ч10 - 4 h 3
h
33
Таким образом, в слое воды глубиной 1 метр при подогреве снизу неустойчивость возникает уже при вертикальной разности температуры величиной всего 10 - 7 градуса, в слое толщ иной 1 сантиметр критическая разность температуры равна 0.1 градуса, а слой воды толщ иной один миллиметр практически абсолю тно
устойчив.
Задача об устойчивости
горизонтального слоя жидкости при наличии вертикального градиента температуры (задача Релея-Бенара) является
классической задачей о конвективной
устойчивости.
И менно в подогреваемом снизу горизонтальном слое жидкости со свободной верхней
границей Бенар в 1900 году
обнаружил возникновение после превышения критического
Рис. 1.11.
градиента температуры гексагональных структур, получивш их название ячеек Бенара (рис.1.10). Ф отография, взятая из работы
[Koschmieder E.L. Adv.Chem.Phys., 1974, V.26. P.177-212.], иллю стрирует
высокую чувствительность гексагональной структуры к возмущ ениям слабая деформация поверхности медной пластины, образую щ ей дно сосуда, приводит к локальному наруш ению вида ячеек. Течение в слое силиконового масла визуализируется с помощ ью алюминиевой пудры.
Отметим, что гексагональные структуры возникают в слое только
при наличии свободной поверхности и направление циркуляции в жидкостях и газах при этом противоположно. В жидкости горячий поток поднимается в центре ячейки, а в газах наоборот - в центре ячейки холодный поток жидкости направлен вниз. Отметим, что возникновение гексагональных структур связано с действием поверхностного натяжения. П ри твердых
горизонтальных границах возникаю т конвективные валы. Этот вид конвективных течений иллю стрирует рис.1.11, где показана валиковая конвекция в слое силиконового масла в круглом сосуде, закрытом сверху стеклом.
Ф орма сосуда навязывает валам осевую симметрию.
Задача Релея. Теоретически задачу о конвективной устойчивости
жидкости впервые реш ил Релей в 1916 году. Он рассмотрел горизонтальный слой жидкости толщ иной h со свободными, но не деформируемыми
границами (такие не совсем реальные граничные условия дают самую про-
34
стую постановку), на которых поддерживается температура T1 и T2 , соответственно. Уравнения Буссинеска записываются в безразмерной форме(на
этот разединицы измерений выбраны следую щ им образом: единица длины
- h , температуры - ( T1 - T2 ), времени - h 2 / n , скорости - c / h ):
r
r
r
¶v 1 r r
1
+ (v С )v = - С P + Dv + RTe z ,
¶t s
s
r
¶T
s
+ (v С )T = DT ,
¶t
r
div v = 0.
Реш ается двумерная задача в плоскости (x,z) , то есть имеются в виду конвективные валы, направленные вдоль оси y. Граничные условия :
z = 0:
z = 1:
¶v x
= 0, v z = 0, T = 1.
¶z
¶v x , y
= 0, v z = 0, T = 0.
¶z
Температура задается в виде T = J - z , так что величина J описывает отклонение температуры от равновесного (линейного) распределения. Для поля
скорости вводитсяфункция тока
vx = -
¶y
,
¶z
vz =
¶y
.
¶x
Рассмотрение ведется в рамках линейной теории устойчивости, то
есть из уравнений выбрасываются все члены, квадратичные по скорости и
возмущ ениям равновесного профиля температуры. В результате получаю тсялинейные уравнения
¶
¶J
Dy = DDy + R
,
¶t
¶x
¶J
¶y
s
= DJ +
.
¶x
¶t
П оследнее слагаемое во втором уравнении - это остаток от нелинейного слагаемого, так как
r
r
(v С )T = (v С )J - v z .
Граничные условия на верхней и нижней границах имею т одинаковый вид:
35
y = Dy = J = 0 .
Следую щ ий ш аг состоит в использовании нормальных возмущ ений,
которыезадаются в форме периодических возмущ ений с экспоненциальной
зависимостью амплитуды от времени:
y = y 0 e - lt sin(pnz ) sin(pax)
J = J0 e - lt sin(pnz ) cos(pax) .
Учитывая, что
(
)
Dy = - p 2 n 2 + a 2 y
DDy = p 4 n 4 + 2a 2 n 2 + a 4 y ,
(
)
получаем уравнения
(
)
lp 2 n 2 + a 2 y
(
0
(
)
= p 4 n 4 + 2a 2 n 2 + a 4 y
)
- lsq 0 = - p n + a q 0 + pay
2
2
2
0
- Rapq 0
0
представляю щ ие собой систему линейных, однородных уравнений для амплитуд y 0 и J0 :
(
)[
(
)]
p 2 a2 + n2 l - p 2 n2 + a2 y
pay
0
[
(
)]
0
+ Rapq 0 = 0
+ ls - p 2 n 2 + a 2 q 0 = 0 .
Система имеет реш ение, если ееопределитель равен нулю
)[
(
(
p a2 + n2 l - p 2 n2 + a2
pa
)]
Ra
=0
ls - p 2 n 2 + a 2
(
)
Раскрывая определитель, получаем уравнение
(
)[
(
p a 2 + n 2 l2s - lg- slg+ p 4 n 2 + a 2
) ]- pRa
2
2
= 0,
реш ение которого дает значения для декремента l :
(
)
(
)
p 2 (1 + s ) n 2 + a 2
p 4 n 2 + a 2 (1 - s )
Ra 2
l=
±
+
.
2s
4s 2
s a2 + n2
2
2
(
)
(1.29)
36
П о виду реш ения (1.29) можно сделать ряд полезных выводов. Вопервых, видно, что при положительных значениях числа Релея (а при принятых обозначениях положительным числам Релея соответствует нагрев
слоя снизу) подкоренное выражение всегда положительно. Это означает,
что оба корня уравнения являются вещ ественными величинами и, следовательно, возмущ ения эволю ционирую т монотонным образом. П ри этом
один корень всегда положителен, а
второй при некотором значении
R = Rc меняет знак.
Во-вторых, при отрицательных
числах Релея (подогрев сверху) вещ ественная часть обоих корней всегда
положительна. Следовательно, все
возмущ ения при подогреве сверху затухают. В то же время с ростом величины подогрева возникает ситуация,
когда выражение под корнем станоРис. 1.12.
вится отрицательным, то есть появляется два комплексно-сопряженных корня, описываю щ их затухаю щ ие, но
колебательные режимы. Это происходит при
p 4 (n 2 + a 2 ) 3 (1 - s ) 2
R =
.
4sa 2
*
Н а рис.1.12 показан график зависимости вещ ественной части декремента
затухания от числа Релея. Н а графике отмечены три области: I - область затухаю щ их колебательных возмущ ений, II - область монотонно затухаю щ их
возмущ ений и III - область монотонно нарастаю щ их возмущ ений.
Н айдем критическое значение числа Релея, при достижении которого
начинается нарастание возмущ ений. И з
условия l = 0 получаем
Rc =
p 4 (a 2 + n 2 ) 3
.
a2
Так как требуется найти самые опасные
возмущ ения, то нужно определить соответствую щ ие значения a и n . Дифференцированиепо a дает
¶R 2p 4 2
= 3 ( a + n 2 ) 2 ( 2a 2 - n 2 ) = 0
¶a
a
и
Рис. 1.13.
37
ac =
n
2
,
Rc =
27p 4 n 4
.
4
Самые малые критические значения появляю тся при n = 1 , что соответствуетодному слою конвективных валов. Следовательно,
ac =
1
2
,
Rc =
27 p 4
= 657,5 .
4
Вид нейтральной кривой показан на рис.1.13.
1.5 М аломодовая модель конвекции (система Лоренца)
В заклю чение вводной части курса остановимся на выводе простой
динамической системы, описываю щ ей конвективные течения в той же самой задаче Релея о конвекции в подогреваемом снизу горизонтальном слое
несжимаемой жидкости. Эта система стала одной из наиболее известных
динамических систем, иллю стрирую щ их переход к хаосу и возникновение
странных аттракторов (см. следую щ ую главу). Н а данном этапе нас интересует сам процесс получения конечномерных проекций уравнений движения жидкости и переход к системе обыкновенных дифференциальных уравнений. Учитывая общ епринятый вид системы Лоренца, мы сохраним единицы размерности и обозначения его работы (Lorenz E., Deterministic Nonperiodic Flow, Journal of Atmospheric Sciences, 1963, V.20, P.130-141.)
Как и в описанной выше задаче Релея рассматриваются только плоские движения жидкости (конвективные валы). Вектор скорости имеет две
r
компоненты v = (v x ,0, v z ) и уравнения Буссинеска, записанные покомпонентно, имеют вид
¶v x
¶v
¶v
1
+ vx x + vz x = ¶t
¶x
¶z
r0
¶v z
¶v
¶v
1
+ vx z + vz z = ¶t
¶x
¶z
r0
¶T
¶T
¶T
+ vx
+ vz
= cDT ,
¶t
¶x
¶z
¶v x ¶v z
+
= 0.
¶x
¶z
¶P
+ nDv x ,
¶x
¶P
+ nDv z + gbT ,
¶z
(1.30)
38
Далееснова вводится функция тока (мы повторяем вывод уравнений, так
как теперь в них сохраняю тся нелинейныеслагаемые)
-
1 ¶P
¶ ¶y
¶y ¶2y
¶y ¶2y
¶y
+
=- nD
r 0 ¶x
¶t ¶z
¶z ¶x¶z ¶x ¶x¶z
¶z
1 ¶P
¶ ¶y
¶y ¶2y
¶y ¶2y
¶y
=+ nD
+ gbT
2
¶t ¶x
¶z ¶x
¶x ¶x¶z
r 0 ¶z
¶z
¶T ¶y ¶T ¶y ¶T
+
= xDT
¶t
¶z ¶x ¶x ¶z
и после обычной процедуры дифференцирования первого и второго уравнений соответственно по z и по x и вычитания первого из второго, получаем
¶
¶T
Dy + {y , Dy }= nDDy + gb
,
¶t
¶x
¶T
+ {y , T }= xDT ,
¶t
(1.31)
где для упрощ ения записи использованы скобки П уассона
{A, B}= ¶A ¶B ¶x ¶z
¶A ¶B
.
¶z ¶x
Учитывая линейную зависимость равновесной температуры по высоте,
представим, как и ранее, температуру в видесуммы
T =q -
DTz
,
h
где q - есть отклонение температуры от линейного профиля. Тогда
¶
¶q
,
Dy + {y , Dy }= nDDy + gb
¶t
¶x
¶q
DT ¶y
+ {y ,q}= cDq .
¶t
h ¶x
Н а границах:
y = Dy = q = 0 .
(1.32)
39
Дальнейш ий путь состоит в том, что функция тока и температура раскладываются в двойные ряды Ф урьесзависящ ими от времени коэффициентами
жpmx ц жpnz ц
(t ) sin з
чsin з
ч,
и h ш и h ш
жpmx ц жpnz ц
q ( x, z , t ) = е q nm (t ) cosз
чsin з
ч.
и h ш и h ш
y ( x, z , t ) = е y
nm
П одставляя эти разложения в уравнения и приравнивая коэффициенты при
одинаковых функциях от x и z , получают систему обыкновенных дифференциальных уравнений для коэффициентов y nm (t ) и q nm (t ) . Отличительной
особенностью модели Лоренца является то, что в разложениях оставлено
минимальное число членов, сохраняю щ их нелинейность системы, а именно, один член из ряда для функции тока и два - для температуры. Этот выбор был обусловлен результатами численных исследований конечномерных систем, проведенных Сольцменом (Saltzman B. Finite amplitude free
convection as an initial value problem, Journal of Atmospheric Sciences, 1962,
V.19, P.329-341.), в которых было показано, что при некоторых значениях
параметров системы действительно возникают режимы, при которых все
остальные переменные стремятся к нулю , а поведение трех оставш ихся характеризуетсянерегулярными непериодическими колебаниями.
М ы, следуя Лоренцу, сразу оставим в разложениях только эти три
члена, обозначив амплитуды соответствую щ их мод как X , Y и Z . Отметим,
что при этом используется не совсем обычный способ обезразмеривания, в
том смысле, что в единицы измерений входят критические параметры. За
единицы измерения приняты величины: длины - h , времени t = h 2 /(p 2 (1 + a 2 ) c ) , функции тока - h 2 / t , температуры - DT . Вводится обозначение b = 4 /(1 + a 2 ) и нормированное число Релея
r=
gb DT h 3 a 2
R
=
.
Rc cn p 4 (1 + a 2 ) 3
Безразмерная система уравнений примет вид
¶
sb
4sr ¶q
Dy + {y , Dy }=
DDy +
,
2
¶t
4p
ba 2 ¶x
¶q
¶y
b
+ {y ,q}=
+
Dq .
¶t
¶x 4p 2
(1.33)
(1.34)
В эти уравнения подставляю тся разложения для функции тока и для температуры в виде
40
y = X (t )
q=
[
2
sin(pax) sin(pz )
p a
2
1
Y (t ) 2 cos(pax) sin(pz ) - Z (t ) sin( 2pz )
pr
]
В уравнении (1.33) скобки П уассона равны нулю и простые преобразования
приводят к уравнению (производные по времени обозначаем точками)
X& = s (Y - X )
Уравнение (1.34) дает
[
]
1 &
Y 2 cos(pax) sin(pz ) - Z&sin( 2pz ) +
pr
2
X
cos(pax) sin(pz ) Y 2 cos(pax) cos(pz ) - 2Z cos(2pz ) +
p r
2
X
sin(pax) cos(pz )Y 2 sin(pax) sin(pz ) =
p r
2
2
b
X
cos(pax) sin(pz ) - Y
cos(pax) sin(pz ) +
Z sin( 2pz ).
pr
pr
p
[
]
Учитывая, что сумма слагаемых, содержащ их произведение XY , дает
[XY (pr ) - 1 sin(2pz)]можно упростить уравнение
[Y& - rX + Y ]cos(pax) sin(pz) - 2 XZ cos(pax) sin(pz) cos(2pz) =
1 &
[Z - XY + bZ ]sin(2pz).
2
Это уравнение разделяется на два путем последовательного умножения на sin(pz ) и на sin(2pz ) и интегрирования по координате z . Таким образом, система уравнений для амплитуд трех выбранных мод выглядит следую щ им образом
X& = s (Y - X )
Y& = - XZ + rX - Y
Z& = XY - bZ
(1.35)
Н апомним, что система (1.35) имеет отнош ение к реальным конвективным
движениям только при небольш их надкритичностях (относительное число
Релея не намного превосходит единицу). Н есмотря на это, поведение этой
системы оказалось интересным само по себе и многочисленные численные
исследования ее свойств проводились в очень ш ироком диапазоне параметра r . В вычислениях обычно использую т число П рандтля s = 10 , а па-
41
раметр b = 8 / 3 , что соответствует результату Релея для критического значения a = 1 / 2 .
Рекомендуемая литература к первой главе:
1. Ландау Л.Д., Лифш иц Е.М . Гидродинамика. М .: Н аука, 1988. 736с.
2. Лойцянский Л.Г. М еханика жидкости и газа. М .: Н аука, 1978.736с.
3. Герш уни Г.З., Ж уховицкий Е.М . Конвективная устойчивость несжимаемой жидкости. М .: Н аука, 1972. 392с.
4. Валандер С.В. Лекции по гидроаэромеханике. Л.: И зд-во Ленингр.
ун-та, 1978. 296с.
42
2 ХАОС В ДИ Н АМ И Ч ЕСКИ Х СИСТЕМ АХ
Турбулентные течения представляю т собой системы, характеризую щ иеся наличием хаотически распределенных и хаотически осцилирую щ их
структур самого различного масш таба. Турбулентность - это воплощ ение
хаоса, а хаос долгое время ассоциировался с системами, имею щ ими огромное число степеней свободы, и развитая турбулентность считалась лиш енной какого-либо порядка. Однако, начиная с конца 60-х годов наш его века
наметился значительный прогресс в понимании природы турбулентности,
связанный с осознанием природы и структуры хаоса.
Во-первых, была установлена возможность хаотического поведения в
нелинейных системах с совсем небольш им числом степеней свободы. И нтересно, что впервые хаотическое поведение в простых гамильтоновых системах обнаружил А.П уанкаре около ста лет назад, но только после работы
Э.Лоренца (1963г.), в которой исследовалось хаотическое поведение диссипативной системы из трех обыкновенных дифференциальных уравнений
первого порядка (1.35), было оценено значение этого факта и началось активное исследования хаотического поведения динамических систем. П равда, произош ло это тоже не сразу, а только после клю чевой работы Д.Рю эля
и Ф .Таккенса 1971г., в которой было сформулировано понятие странного
аттрактора и указана его роль в формировании нерегулярного поведения
системы.
Во-вторых, было понято, что даже в самом развитом турбулентном
потоке сущ ествуют элементы порядка, а число реально возбужденных степеней свободы значительно меньш е ожидаемого. В 70-80-х годах появляю тся многочисленные работы о когерентных структурах в турбулентных потоках и делаются первые попытки описания турбулентности на языке
фракталов.
И менно в это время сформировались такие науки, как теория катастроф и синергетика, появились первые книги о «детерминированном хаосе»
и «порядке в хаосе». Важно подчеркнуть, что обычно рассматриваемые в
этих книгах проблемы динамических систем невысокого порядка не имеют
прямого отнош ения к развитой турбулентности. В них речь идет о хаотическом во времени поведении небольш ого числа заданных в пространстве
мод (такие течения реально сущ ествуют при небольш их надкритичностях,
то есть вблизи порога неустойчивости), в то время, как «истинная» турбулентность хаотична и в пространстве и во времени. Тем не менее, рассматриваемые в качественной теории динамических систем вопросы чрезвычайно полезны как для понимания путей развития турбулентных течений
43
(сценариев перехода к хаосу), так и для отработки методов описания хаотических (в том числе и турбулентных) систем.
Н еобходимо остановиться на самом понятии детерминированный хаос. П од ним понимают нерегулярное поведение нелинейных систем, эволю ция которых однозначно описывается динамическими уравнениями при заданных начальных условиях. П ри этом нелинейность является необходимым, но не достаточным условием возникновения хаотического поведения,
а его возникновение связано не с наличием источников ш ума или бесконечного числа степеней свободы, а со свойством нелинейных систем экспоненциально быстро разводить реш ения в ограниченной области фазового
пространства.
В данной главе мы остановимся на базовых понятиях теории динамических систем, рассмотрим основные виды бифуркаций и основные сценарии перехода от упорядоченного движения к хаосу. М ы подробно разберем
свойства системы Лоренца, не только сыгравш ей важнейш ую роль в становлении науки о хаосе, но и имею щ ей самое прямое отнош ение к теме наш его курса. Далее мы приведем пример ещ е одной динамической системы,
имею щ ей отнош ение к гидродинамическим системам - это простейш ая модель земного динамо Рикитаке. В заверш ение будут приведены некоторые
результаты лабораторного исследования стохастизации конвективного
движения в замкнутой полости.
2.1 Консервативные и диссипативныесистемы
Любые движения можно разделить на монотонные и колебательные,
а колебательные в свою очередь, на регулярные (периодические) и нерегулярные. Среди периодических колебаний наиболее изучены гармонические
колебания. Это вполне естественно, так как гармонические колебания
чрезвычайно ш ироко распространены в самых различных системах, а также потому, что лю бой колебательный процесс с помощ ью преобразования
Ф урье может быть представлен как сумма гармонических колебаний. Н е удивительно, что знакомство с динамическими системами традиционно начинаю т с рассмотрения простого осцилятора.
Рассмотрим хорош о известный со ш кольной
скамьи математический маятник - точечное тело массой m , подвеш енное на стержне длиной l и находящ ееся в поле силы тяжести, характеризуемой ускорением
свободного падения g (Рис.2.1). М аятник имеет одну
степень свободы, описываемую углом отклонения от
Рис. 2.1.
44
вертикали q . Основной закон механики приводит к уравнению, которое в
цилиндрических координатах имеет вид
g
q&&+ sinq = 0 .
(2.1)
l
Для малых угловых отклонений, когда sinq » q , уравнение (2.1) стано-
витсялинейным уравнением
g
q&&+ q = 0 ,
l
(2.2)
реш ением которого являются гармонические колебания q = q0 sin(w t + j 0 ) с
круговой частотой w = g l .
2.1.1 Ф азовое пространство
Состояние маятника в любой момент времени полностью задается
двумя величинами: положением q (t ) и угловой скоростью q&(t ) . Если мы
введем систему координат, осями которой будут служить эти две величины, то точка на плоскости (q ,q&) будет полностью характеризовать состояние системы, а любому реш ению будет соответствовать та, или иная линия
(траектория).
Ф азовое пространство определим как пространство, в котором осями
координат служат переменные, описываю щ ие состояние системы, в случае
осцилятора - положение и скорость. Ф азовой траекторией называется кривая в фазовом пространстве, описываю щ ая эволю цию системы. Совокуп-
Рис. 2.2.
45
ность фазовых траекторий, описываю щ их эволю цию системы при различных начальных условиях, образуетфазовый портрет системы.
Н а рисунке 2.2 приведен фазовый портрет маятника. Картина периодична по оси q с периодом 2p . В области применимости уравнения (2.2)
фазовые траектории представляют собой окружности с центрами в точках
q& = 0, q = ±2pn, n -целое число. Эти кривые соответствуют гармоническим
колебаниям, частота которых не зависит от амплитуды. С ростом амплитуды колебаний траектории принимают эллиптическую форму и период
колебаний растет. Если энергия колебаний превышает величину 2 g / l , то
колебания переходят во вращ ения вокруг оси. Траектории, точно соответствую щ ие этому значению энергии, проходят через верхнее, неустойчивое
положение равновесия и период колебаний стремится к бесконечности. Эта
траектория разделяет области фазового пространства с различным характером поведения (колебания и вращ ение) и являетсясепаратрисой. Стрелки
на рисунке указываю т направление движения.
2.1.2 Консервативныесистемы
М аятник, описываемый уравнением (2.1) сохраняет энергию. Действительно,
E=
йq&2 g
щ
mv 2
+ mgl( 1 - cos q ) = ml 2 к + ( 1 - cos q )ъ
2
l
л2
ы
и
dE
g
й
щ
= ml 2 кq&&+ sin q ъq& є 0 .
dt
l
л
ы
(2.3)
Это означает, что линии на рисунке 2.2 можно интерпретировать как
линии равной энергии. Энергия с точностью до множителя совпадает с
функцией Гамильтона, а уравнение (2.1) приводится к системе уравнений
первого порядка
¶H
q& =
,
¶p
Здесь H ( p,q ) =
p& =
¶H
.
¶q
p2
+ g (1 - cosq ) , и p = lq&.
2l
Таким образом, рассмотренный маятник относится к гамильтоновым
системам, которые, как известно, консервативны.
И з консервативности (сохранения энергии) следует одно очень важное свойство - сохранение площ адей (в общ ем случае - объема) в фазовом
46
пространстве. Элемент объема в фазовом пространстве можно рассматривать как множество начальных условий. В процессе эволюции это множество преобразуется в другой элемент фазового пространства (каждая точка
следует своей фазовой траектории), объем которого должен оставаться постоянным.
Следует подчеркнуть, что сохранение объема не подразумевает при
этом сохранения формы, так как сохранение объема может достигаться
двумя различными способами. В первом случае элемент фазового объема
переносится вдоль траектории практически без деформации. Во втором
случае происходит экспоненциальное удлинение объема в некотором направлении с одновременным сжатием в перпендикулярном направлении
(также экспоненциальным). Хотя фазовый объем сохраняется в обоих случаях, поведение системы отличается принципиально. В первом случае траектории, близкие в начальный момент времени, остаю тся близкими - траектории (а следовательно, и реш ение) устойчивы. Во втором случае малое
начальное возмущ ение приводит к быстрому расхождению траекторий со
временем - они неустойчивы.
Отметим ещ е одно свойство консервативных систем, состоящ еев том,
что они инвариантны к обращ ению времени (замене t на - t ). В случае маятника это означает, что если его движения заснять на видеофильм, то
фильм можно прокручивать в обоих направлениях и отличить правильное
направление от обратного по воспроизводимым движениям маятника будет невозможно.
2.1.3 Диссипативныесистемы
П римером простейш ей диссипативной системы может служить тот же
простой маятник, но подверженный действию сил трения. Реально силы
трения присутствуют всегда (трение на оси, сопротивление воздуха и т.д.) и
ни один свободный осцилятор несоверш ает колебания неограниченно долго. Для учета действия сил сопротивления нужно добавить в уравнение
(2.1.) слагаемое, например, пропорциональное скорости движения
маятника
g
q&&+ mq& + sinq = 0 ,
l
Рис. 2.3.
(2.4)
где m есть коэффициент трения.
П овторяя вычисления для
скорости изменения энергии, вместо (2.3) получим теперь
47
dE
= - mml 2q&2 .
dt
(2.5)
Таким образом, при лю бом положительном значении коэффициента трения энергия убывает со временем, стремясь в конечном итоге к нулю (отрицательной энергия стать не может). Это означает, что семейство траекторий, представлявш ее собой в отсутствие трения множество концентрических окружностей, превращ ается теперь во множество траекторий, сходящ ихся к началу координат. Н а рисунке 2.3 показаны фазовые портреты маятника с трением для малого (а) и больш ого (б) трения. В первом случае характерное время
затухания значительно превышает период колебаний и траектории представляют собой спирали
с малым ш агом. Соответствую щ ий фазовый
портрет называется фокусом. Во втором случае
затухание происходит за время, меньш ее периоРис. 2.4.
да. Колебания становятся апериодическими, а
портрет называется узлом. В обоих случаях все
фазовые траектории заканчиваю тся в одной точке, которая называетсяпритягиваю щ ей точкой или аттрактором.
Н аличие аттрактора является важнейш им свойством диссипативных
систем. Аттрактор является точкой только в простейш их случаях. В общ ем
случае аттрактор - это притягиваю щ ее множество (линия, поверхность и
т.д.). П редставим, что в рассматриваемом нами осциляторе добавлена вынуждаю щ ая сила (для конкретности представим себе гирю в часахходиках). Теперь, независимо от начальных условий фазовые траектории
сходятся к окружности, радиус которой определяется действую щ ей силой
Рис. 2.5.
(рис.2.4). Эта круговая траектория и является аттрактором (предельным
циклом). Важным является тот факт, что в диссипативной системе пропала
зависимость реш ения от начальных условий (на достаточно больш их временах, когда система выходит на аттрактор).
48
Рассмотренный пример иллю стрирует ещ е одно важнейш ее свойство
диссипативных систем - сжатие площ адей (объема) в фазовом пространстве. Объем любого множества начальных условий уменьш ается в среднем во
времени. Однако, как и в консервативных системах, эволюция множества
может происходить различным образом. И ногда (как в простом маятнике с
трением) это множество равномерно стягивается в точку (или стремится к
предельному циклу) и все траектории сближаю тся со временем. Н о не всегда уменьш ение объема подразумевает неизбежное сокращ ение длин. Растяжение объема в одном направлении может компенсироваться более эффективным сжатием в другом направлении. Эти два сценария сжатия фазового объема показаны на рисунке2.5.
П оследнее принципиальное отличие диссипативных систем от консервативных связано с тем, что они не инвариантны к обращ ению времени.
Если фильм о затухаю щ ем маятнике просматривать в обратном направлении, то маятник станет раскачиваю щ имся.
2.1.4 П ример немеханической системы
П риведем простой пример диссипативной системы из живого мира.
Это модель системы жертва - хищ ник. Система бесспорно диссипативна,
так как в отсутствие пищ и любая биологическая популяция вымирает.
П усть в изолированном лесу обитаю т только зайцы и волки, за популяциями которых мы и собираемся следить ( N - количество волков, n - количество зайцев). Ф азовое пространство есть в этом случае один квадрант
на плоскости (n, N ) , так как отрицательные значения для численности животных не возможны. П остараемся нарисовать фазовый портрет системы,
не выписывая уравнений.
Какие параметры определяю т возможные сценарии развития жизни в
лесу ? Это рождаемость обоих видов, естественная смертность, аппетит
волков. Очевидно, что у каждого вида есть наименьш ее критическое число
(соответственно, nc и N c ), необходимое для того, чтобы вид мог воспроизводиться. Отложим на осях эти критические значения и подумаем, как может развиваться система если начальные условия задаю т старт фазовой
траектории вблизи осей координат. Ясно, что реш аю щ им является число
зайцев. Если количество зайцев не достаточно для поддержания вида, то
вымрут зайцы, а следом с неизбежностью вымрут и волки. Если мало волков ( N < N c ) , а зайцев достаточно, то после вымирания волков численность
популяции зайцев (в упрощ енной модели) будет зависеть только от наличия травы в наш ем лесу (обозначим это число как nm ). Таким образом, в
системе выявились две притягиваю щ ие точки, каждая из которых имеет
свою область притяжения.
49
Если число зайцев и волков достаточно, то наиболее вероятное развитие событий - это возникновение колебаний: размножились волки уменьш ается число зайцев, стало мало зайцев - уменьш ается численность
волков, стало меньш е волков - снова размножаю тся зайцы и т.д. Такой
сценарий немедленно следует и из простейш ей модельной системы
n& = a n - bnN ,
N& = - g
N + dnN ,
(2.6)
где a - рождаемость зайцев, g - смертность волков (смертностью зайцев от
старости пренебрегаем), b и d - коэффициенты, описываю щ ие результат
встречи зайцев с волками (как часто такие встречи заканчиваю тся трагически и сколько волков могут насытиться в результате одной удачной охоты).
Система (2.6) имеет стационарное реш ение: n = a / b , N = g/ d, а линеаризациясистемы вблизи точки равновесия приводит к уравнению
n&&= ag
n,
имею щ им своим реш ением гармонические колебания. Таким образом, если
стационарное реш ение является неустойчивым, то можно ожидать появления в системе предельного цикла. Всесказанноесуммирует рисунок 2.6, где
приведен качественный вид фазового портрета системы зайцы - волки.
Видно, что аттрактор системы вклю чает два узла и предельный цикл, и что
каждый изтрех элементов аттрактора имеет свою область притяжения. Области притяжения разделены сепаратрисами, обозначенными пунктиром.
Рис. 2.6.
50
2.2 Бифуркации
2.2.1 Что такое бифуркация ?
В рассмотренных нами примерах
диссипативных систем с подводом энергии
(маятник, энергия которого поддерживается за счет опускаю щ ейся гири, животные в лесу, питаю щ иеся в конечном итоге
за счет травы) мы обош ли молчанием
важный вопрос о том, как устойчивое реш ение (точка в фазовом пространстве)
становится неустойчивым и сменяется
предельным циклом. Ясно, что поведение
системы зависит от некоторых управляю щ их параметров (масса гири в часах, при
недостатке которой маятник остановится,
рождаемость зайцев и т.д.) и при изменение этого параметра возможны не только
количественные, но и качественные перестройки характера эволю ции системы.
Точка в пространстве параметров,
при которой происходят качественные изменения характера реш ений, называется
точкой бифуркации, а соответствую щ ее
значение параметра называется критическим. Вспомним результаты анализа конвективной устойчивости нагретой жидкости в горизонтальном слое, описанные в
первой главе и представим их на плоскости ( R, A) , где R - число Релея, а A - амплитуда (скорость вращ ения) конвективных
валов (см. рис.2.7). П ри R < Rc , единственным реш ением является устойчивая неподвижная точка (конвекция отсутствует).
В точке R = Rc рождается дополнительная
пара реш ений (это также устойчивые точки), каждое из которых соответствует
вращ ению валов в ту или иную сторону.
Рис. 2.7.
Рис. 2.8.
51
П ри этом прежнее реш ение становится неустойчивым. В этой точке имеет
место бифуркация, называемая вилкой (ответвление пары реш ений в виде
притягиваю щ их точек). Таким образом, точкой бифуркации называется
точка, в которой происходит ветвление реш ений.
2.2.2 Бифуркация Хопфа
Бифуркацией Хопфа называется процесс рождения предельного цикла из точки. П оведение системы вблизи точки бифуркации иллю стрирует
рисунок 2.8. Н а рисунке схематически изображены фазовые траектории
при трех значениях управляю щ его параметра e : e < ec , e = ec , e > ec .
Отметим два важных свойства бифуркации Хопфа. Во-первых, вблизи точки бифуркации период колебаний не зависит от величины надкритичности e - eс . Во-вторых, амплитуда колебаний (амплитуда предельного
цикла) зависит от надкритичности по корневому закону, то есть пропорциональна величине | e - ec | .
И менно с бифуркацией Хопфа связан первый предложенный сценарий перехода от ламинарного течения к турбулентности (Ландау, 1944г.).
Согласно сценарию Ландау переход к турбулентности представляет собой
бесконечную цепочку бифуркаций Хопфа, каждая из которых приводит к
появлению новой частоты. В такой схеме аттрактор представляет собой n мерный тор с n , стремящ имся к бесконечности, и хаос рождается в системе
с очень больш им числом степеней свободы.
2.2.3 Н ормальные и обратныебифуркации
П редставленная на рис.2.7 бифуркационная диаграмма соответствует
нормальной (суперкритической) бифуркации вилки. Это означает, что возникаю щ ая в точке бифуркации пара реш ений ответвляется от начального реш ения
мягко, то есть с нулевой начальной амплитудой, которая монотонно растет по
мере роста надкритичности.
Точно также нормальной (суперкритической) называется бифуркация
Хопфа, если предельный цикл рождается
с нулевой амплитудой и в точке бифур-
Рис. 2.9.
52
кации система находится в состоянии нейтральной устойчивости. П о мере
удаления от точки бифуркации происходит плавное увеличение амплитуды
предельного цикла.
Возможна и другая картина, когда в точке бифуркации происходит
жесткий переход к циклу конечной амплитуды (или, в случае вилки, две новые точки появляю тся на конечном расстоянии друг от друга). Это происходит, когда нелинейные члены в уравнениях стремятся усилить возникаю щ ую неустойчивость. П роходя точку бифуркации справа налево
(рис.2.9) можно видеть, что неустойчивая неподвижная точка превращ ается
в устойчивую неподвижную точку и неустойчивый предельный цикл. Такая
бифуркация называетсяобратной или субкритической.
Важной особенностью обратных бифуркаций является наличие интервала управляю щ его параметра eсў < e < ec , в котором сосущ ествуют два
устойчивых реш ения. Какое из этих реш ений реализуется, зависит от предыстории: при движении слева направо неподвижная точка остается устойчивой до значения e = eс , после чего реш ение перепрыгивает на одну из
двух устойчивых ветвей. П ри движении справа налево реш ение следует
вдоль этой ветви до точки e = eсў, где скачком переходит в устойчивую неподвижную точку на оси.
Такое явление называется гистерезисом и хорош о известно в самых
различных областях физики и механики.
2.3 Как описать переход и хаос ?
2.3.1 Сечения П уанкаре
Рис. 2.10.
И дея метода П уанкаре состоит в
снижении объема обрабатываемой
информации при изучении поведения
фазовых траекторий путем рассмотрения лиш ь дискретного ряда точек на
траектории. Реализуется эта идея путем выбора некоторой (вообщ е говоря, произвольной) плоскости в фазовом пространстве и наблю дения за
точками пересечения этой плоскости
фазовыми траекториями. М етод поясняет рисунок 2.10, где для трехмерного фазового пространства показаны
53
точки пересечения плоскости фазовой траекторией (причем фиксируются только точки,
в которых траектории пересекают плоскость
в одном направлении, в данном случае, сверху вниз).
М ножество точек пересечения Pi образую т сечение П уанкаре, а преобразование,
связываю щ ее последую щ ую точку с предыдущ ей
Pi + 1 = T ( Pi )
(2.7)
называетсяотображением П уанкаре.
П ри переходе от фазовых траекторий к
Рис. 2.11.
сечению П уанкаре происходит снижение
размерности исследуемого множества. П ри
этом рассматривается не система дифференциальных уравнений с непрерывным временем, а отображение (2.7) с дискретным временем и дифференциальные уравнения заменяются разностными. В то же время, сечение П уанкаресохраняет топологические свойства
породивш его его потока. Так для консервативной системы сечение сохраняет, а для диссипативной сокращ ает площ ади на плоскости S .
Если реш ение системы периодическое, характеризуемое частотой f1 ,
то фазовая траектория представляет собой замкнутую кривую и сечение
П уанкаре представляет собой в простейш ем случае одну единственную
точку (или несколько точек, если траектория очень извилистая и/или неудачно выбрана плоскость сечения). Если в реш ении появляется вторая
частота f 2 и аттрактор представляет собой двумерный тор, то точки в сечение П уанкаре ложатся на замкнутую кривую, которая может иметь или
не иметь точек самопересечения (рис.2.11). П ри этом точки могут образовывать на этой кривой конечное множество, если отнош ение частот f1 / f 2
рационально и фазовая траектория представляет собой замкнутую линию ,
или покрывать кривую непрерывным образом, если отнош ение частот иррационально.
П осмотрим, как выглядит проблема устойчивости периодического
реш ения с точки зрения отображения П уанкаре. Вопрос состоит в том, является ли замкнутая траектория устойчивой по отнош ению к малым возмущ ениям. И наче говоря, нужно узнать, как изменится положение точки P
на следую щ ем ш аге, если на данном ш аге внести возмущ ение в ее положение. Ограничиваясь линейным анализом устойчивости, для описания отображения П уанкаре T (P ) вводят матрицу
54
й¶T щ
M = к i ъ, i, j = 1,2 ,
к
л¶x j ъ
ы
(2.8)
называемую матрицей Ф локе. Эта матрица характеризует реакцию отображения T вдоль координаты i на возмущ ение вдоль координаты j . Устойчивость цикла определяется собственными
значениями матрицы (2.8). Смещ ение траектории на следую щ ем витке экспоненциально
убывает со временем, если все собственные
значения лежат внутри единичной окружности на комплексной плоскости. Ели же какоелибо собственное значение становится по модулю больш е единицы, то смещ ения растут со
временем и цикл становится неустойчивым.
Рис. 2.12.
И зучение свойств матрицы Ф локе позволяет не только определить устойчив или
нет предельный цикл, но и узнать вид бифуркации, соответствую щ ей потере устойчивости. П отеря устойчивости, как уже отмечалось выше, происходит при пересечении модуля собственного значения через единичную окружность. Это пересечение может происходить тремя различными способами (рис.2.12).
В первом случае, собственное значение действительно и пересекает
окружность в точке+1. Этот переход соответствует бифуркации узел-седло,
означаю щ ей, что появляется одно неустойчивое направление и периодическое движение разруш ается.
Во втором случае, собственное значение также действительно, но пересекает окружность в точке -1. М омент перехода соответствует ситуации,
когда траектория через раз снова попадает в прежню ю точку. Это так называемая бифуркация удвоения периода (субгармоническая бифуркация).
Она может быть нормальной и обратной. П ри нормальной субгармонической бифуркации реш ение заменяется новым устойчивым периодическим
реш ением с удвоенным периодом (см. параграф 1.7), при обратной бифуркации возникает временная перемежаемость, когда долгие интервалы почти периодического движения сменяю тся хаотическими осциляциями.
Третий тип перехода возникает при комплексных собственных значениях. В этом случае пара комплексно-сопряженныхзначений одновременно
пересекает единичную окружность. Этот переход отвечает бифуркации
Хопфа (возникает блуждание траектории вокруг устойчивой прежде точки). Если бифуркация нормальная, то предельный цикл переходит в тор,
если обратная, то вновь возникаетперемежаемость.
55
2.3.2 П оказатели Ляпунова
Теория Ф локе рассматривает устойчивость замкнутой фазовой траектории, интересуясь при этом только
поведением всего цикла в целом. М ожно поставить вопрос и о локальной устойчивости траектории, независимо от
того, является ли она замкнутой или
нет. И наче говоря, речь идет о характеристике скорости расхождения (схождения) начально близких траекторий
в фазовом пространстве. Количественной мерой расходимости траекторий
являю тся показатели Ляпунова.
Чтобы ввести показатели Ляпунова, необходимо рассмотре
ть эволю r
цию малого возмущ ения dX (t ) фазовой
r
траектории X (t ) . И нтегрируя числено
исследуемую систему уравнений, можно построить матрицу M , связываю щ ую вектор возмущ ений в момент
времени t + dt с вектором в момент времени t :
r
r
dX (t + dt ) = M (dt )dX (t ) .
Для n - мерной системы матрица
M будет иметь размерность n 2 и
n собственных значений. Траектория
устойчива, если модули всех собственных чисел меньш е единицы (или показатели степени при экспоненциальном
представлении собственных чисел отрицательны). Н а практике интерес
представляет наиболее опасное направление и определяется только один,
самый больш ой показатель Ляпунова.
И сходя из того, что на конечных временах возмущ енная траектория уходит
Рис. 2.13.
56
в самом неустойчивом направлении, практическое определение первого
показателя Ляпунова
можно реализовать по следую щ ей схеме.
r
r
В точке X (t ) на заданной траектории вносится возмущ ение dX (t ) , отстоящ ее на расстояние d 0 от основной траектории. Реш ая далее исследуемую систему уравнений для невозмущ енного и возмущ енного реш ения, вычисляют расстояниемежду траекториями d (t ) через промежуток времени t .
Далее, возмущ енную точку снова устанавливают на расстоянии d 0 от основной
траектории, но так, что она остается в том направлении от точки
r
X (t + t ) , что было получено в результате вычислений возмущ енного реш ения. Тем самым на каждом ш аге мы вычисляем скорость расхождения траекторий в наиболее опасном направлении. Считая, что расхождение траекторий подчиняется экспоненциальному закону d (t + t ) = d 0 e l t и многократно повторяя эту процедуру, приходим к следую щ ей формуле для вычисления первого показателя Ляпунова:
1
d
1 m
l1 = lim
еi=1 ln d i
m ® Ґ mt
0
.
2.3.3 Энтропия Колмогорова
Другой важной характеристикой хаотического движения в фазовом
пространстве является энтропия Колмогорова (К-энтропия). Н апомним,
что энтропия есть мера беспорядка (в термодинамике) или мера информации, необходимой для определения положения системы в некотором состояния (в теории информации) и определяется формулой
S=-
е P ln P ,
i
i
i
где Pi есть вероятность нахождения системы в состоянии i .
П усть система эволю ционирует в d - мерном фазовом пространстве,
которое разбивается на ячейки размера l (всего l d ячеек). Состояние системы фиксируется через интервалы времени t и на каждом ш агеrрегистрируется номер ячейки, в которой оказалась фазовая траектория X (t ) . Обозначим Pi ....i совместную вероятность того, что система, стартовав при t = t 0 в
ячейке i0 , прош ла через ячейки i1 ,i2,.... и в момент t = t 0 + nt оказалась в ячейке
in . И нформация, необходимая для определения положения системы на заданной траектории, пропорциональна энтропии Ш енона
0
n
57
Kn = -
е
Pi0 ....in ln Pi0 ....in .
i0 ...in
Тогда, если известно, что система прош ла цепочку состояний i0 ...in , то для
предсказания положения системы на следую щ ем ш аге требуется дополнительная информация K n+ 1 - K n . И наче говоря, эта разность описывает потерю информации на ш аге n + 1 .
Энтропия Колмогорова вводится как характеристика скорости потери информации
K = lim lim lim
t ® 0 l ® 0 m® Ґ
1
mt
m- 1
е
n=0
( K n + 1 - K n ) = - lim lim lim
t ® 0 l ® 0 m® Ґ
1
mt
е
Pi0 ....im ln Pi0 ....im .
(2.9)
i0 ...im
П роцедуру вычисления энтропии иллю стрирует рисунок 2.13 на примере одномерной системы с дискретным временем. Ось абсцисс соответствует времени, разбитому на интервалы длиной t . П ри рассмотрении дискретного времени предел по t не берется. ВероятностьPi = l , а число ячеек,
в которые может попасть система на следую щ ем ш аге пусть остается постоянным и равным N . Тогда вероятность Pi i = l / N , Pi i i = l / N 2 , а
Pi ...i = l / N m . Тогда общ еечисло возможных траекторий есть M = N m / l и
0
01
0
01 2
m
K = - lim lim
l ® 0 m® Ґ
1 M
1
l
Pi0 ....im ln Pi0 ....im = - lim lim M m (ln l - m ln N ) = ln N .
е
l ® 0 m® Ґ m
m 1
N
Н а рис.2.13а показан пример регулярного движения, когда из ячейки
i0 система однозначно переходит в данную ячейку i1 и т.д., а первоначально
близкие траектории остаю тся близкими. В этом случае N = 1 и K = 0 . В случае, показанном на рис.2.13б, близкие траектории расходятся экспоненциально и N = e l . Тогда K = l и, как видим, К - энтропия совпадает в этом
случае с показателем Ляпунова. П оследний случай (рис.2.13в) соответствует случайной системе, в которой на каждом ш аге система с равной вероятностью оказывается в любой ячейке. Это приводит к тому, что N ® Ґ и
K® Ґ.
58
2.4 Спектры Ф урье
2.4.1 Н епрерывноеи дискретное преобразование Ф урье
Анализ Ф урье играет особую роль при исследовании не только периодических, но также квазипериодических и стохастических сигналов. В
контексте задач, рассматриваемых в этой главе, он интересует нас как инструмент, позволяю щ ий отличать периодические режимы от стохастических, но значение метода Ф урье в изучении проблемы турбулентности этим
не исчерпывается. В дальнейш ем мы увидим, насколько он полезен при
численном исследовании турбулентных потоков и при обработке результатов измерений. Все это делает необходимым краткое изложения основных свойств непрерывного и дискретного преобразования Ф урье.
Н апомним, что Ф урье предложил разложение функций в ряд по гармоническим функциям как метод реш ения уравнения теплопроводности,
которое в одномерном случае имеет вид
¶t T = h¶ xxT .
(2.10)
Если задача реш ается на отрезке (0,L) и имеет, например, нулевые
граничные условия, то температура представляется рядом
ж 2pnx ц
T ( x, t ) = е bn (t )sin з
ч.
и L ш
n
(2.11)
П одстановка (2.11) в (2.10), дает уравнение
е
n
2
ж 2pnx ц
ж 2pn ц
ж 2pnx ц
b&n (t )sin з
ч,
ч = - h е bn (t )з
ч sin з
и L ш
и L ш
и L ш
n
(2.12)
которое распадается на отдельные уравнения для каждой гармоники (для
этого достаточно умножить уравнение на sin(2pm / L) и проинтегрировать по
рассматриваемому отрезку)
2
ж 2pm ц
b&m (t ) = - з
ч hbm (t ) .
и L ш
(2.13)
Реш ение поставленной задачи становится в результате тривиальным:
после разложения в ряд для каждой гармоники имеется реш ение (2.13),
имея которые, можно восстановить по (2.11) распределение температуры в
лю бой момент времени.
59
В общ ем случае периодическую функцию f (t ) с периодом T , для ко+ T /2
торой сущ ествуетинтеграл т- T / 2 | f (t ) | dt , можно разложить в ряд Ф урье:
f (t ) =
a0
+
2
Ґ
е (a n cos(nw 0 t ) +
n =0
bn sin( nw 0 t ) )=
Ґ
еce
n=- Ґ
inw 0 t
n
,
(2.14)
где w 0 = 2p / T , а коэффициенты Ф урье определяются выражениями:
T /2
T /2
2
an =
f (t ) cos(nw 0 t )dt ,
T - Tт/ 2
2
bn =
f (t ) sin( nw 0 t )dt ,
T - Tт/ 2
T /2
c n = c -* n =
1
f (t )e - iw 0t dt ,
T - Tт/ 2
(2.15)
(2.16)
где звездочкой обозначено комплексноесопряжение.
Действительную функцию f (t ) можно представить интегралом Ф у+Ґ
рье, если для неесущ ествуетинтеграл т- Ґ | f (t ) | dt . Тогда
+Ґ
f (t ) =
тf?(n )e
-Ґ
+Ґ
f?(n ) =
тf (t )e
2pin t
dn ,
- 2pin t
dt .
(2.17)
(2.18)
-Ґ
Здесь f?(n ) есть фурье-образ функции f (t ) , n - частота (будем также
пользоваться круговой частотой w = 2pn ). Отметим, что когда речь идет о
преобразовании Ф урье от функции координат f (x) , то в преобразовании
вместо частот появляются волновые числа k и g ( k = 2pg, в полной аналогии с частотами).
2.4.2 Основныесвойства фурье-преобразования
П риведем формулировки основных теорем, касаю щ ихся свойств непрерывного фурье-преобразования, помня при этом, что все они имею т
прямой аналогвтерминах дискретного преобразования.
И так, пусть f (x) - действительная функция, для которой сущ ествует
+Ґ
интеграл т- Ґ | f ( x) | dx . Тогда
60
f ( x )=
f?(k )=
или, с учетом связи k = 2pg,
1
2p
1
)
тf (k )e
ikx
тf (x )e
2p
dk
- ikx
dx
(2.19)
(2.20)
)
f ( x )= тf (g)e 2pigx dg
)
f (g)= тf ( x )e - 2pigx dx
И спользуя для преобразования Ф урье обозначение
)
~
f (k )= F [f ( x )],
сформулируем его основныесвойства.
1. Единственность: преобразование (2.19)-(2.20) однозначно.
2. Линейность:
)
)
~
F[
a 1 f1 ( x )+ a 2 f 2 ( x )]= a 1 f1 (k )+ a 2 f 2 ( k ) .
(2.21)
3. Теорема о масш табах:
1 )ж k ц
~
F [f (a x )]= f з ч.
a иa ш
(2.22)
)жk ц
~
F [f ( x + a )]= e ika f з ч.
иn ш
(2.23)
)
)
~
F ( f1 * f 2 )= f1 (k )Чf 2 (k ),
)
)
~
F ( f1 Чf 2 )= f1 (k )* f 2 (k ).
(2.24)
.
4. Теорема о сдвиге:
5. Теорема о свертке2:
6. Теорема о дифференцировании:
)
~
n
F ( f (n )( x ))= (ik ) f (k ).
(2.25)
7. Теорема П арсеваля3:
)
)*
тf (x )f (x )dx = тf (k )f (k )dk .
*
1
2
1
2
(2.26)
f1 ( x) * f 2 ( x) = тf1 ( x - x ў)f 2 ( x ў)dx
)
3 Важным следствием теоремы П арсеваля является сохранение энергии: | f ( x )|2 dx = | f (k )|2 dk
т
т
2
Н апомним, что сверткой называется интегральная операция
61
8. Теорема о комплексном сопряжении:
)
~
F ( f * )= f * (- k ).
(2.27)
)
)
)
Если f - вещ ественное число, то F~ ( f * )= F~ ( f )= f * (- k ) т.е. f (k )= f * (- k )
2.4.3 Спектры
П усть имеется временной сигнал f (t ) , для которого сущ ествуют преобразования (2.17)-(2.18). Для этого сигнала можно ввести корреляциионную функцию (автокорреляцию )
T
1
C (t ) = lim тf (t ) f (t + t )dt .
T® Ґ T
0
(2.28)
Корреляционная функция (2.28) есть среднее произведение двух значений
сигнала, сдвинутых на величину t и характеризует степень зависимости текущ его значения сигнала от его предыдущ их значений.
Спектральной плотностью сигнала f (t ) называется функция
)
F ( n ) =| f ( n ) |2 . Связь спектральной плотности с автокорреляционной функцией устанавливает теорема Хинчина-Винера:
+Ґ
F (n ) = тС (t )e
- 2pin t
dt .
(2.29)
-Ґ
Следует отметить, что обрабатываемые сигналы представляю т собой, как правило, последовательность дискретных точек (по
крайней мере, сигнал становится
таковым на этапе ввода в цифровую вычислительную маш ину). В
Рис. 2.14.
этом случае приходится иметь дело с конечной выборкой и важной
становится теорема Котельникова, утверждаю щ ая, что функция f (t ) ,
спектр которой ограничен конечным интервалом частот - n max < n < n max , однозначно определяется выборкой на дискретном множестве точек с ш агом
62
Dt = 1 / 2n max . Точнее говоря, функция f (t ) восстанавливается по конечной
выборке f n = f (nDt ) с помощ ью соотнош ения
f (t ) =
Ґ
е
n=- Ґ
жt
ц
sin p з - n ч
и Dt
ш.
fn
жt
ц
pз - nч
и Dt
ш
(2.30)
Рис. 2.15.
Другими словами, теорема Котельникова устанавливает предельную
частоту, которая может быть определена по сигналу, регистрируемому с
ш агом Dt .
П ри дискретной выборке, состоящ ей из N равноотстоящ их точек, исходному ряду соответствует ряд фурье-коэффициентов (2.16), которые для действительного сигнала равны
- 2pimn
N
)
fn = е fme N
(2.31)
m =1
Спектральной плотности F (n ) при дискретном представлении соот)
ветствует ряд величин Fn =| f n |2 , называемый спектром мощ ности (а также
энергетическим спектром, или
просто спектром Ф урье).
Остановимся на том,
как выглядят спектры различных типов сигналов. Н ачнем со случая, когда функция
f (t ) есть периодический сигнал с периодом T . В про-
Рис. 2.16.
63
стейш ем случае это гармонический сигнал (синус или косинус) и его спектр
состоит из одной ненулевой компоненты с частотой 1 / T . Для периодического сигнала другой формы в спектре появляю тся кратные гармоники (с
частотами 2 / T , 3 / T , 4 / T ,..... ) (рис.2.14).
Более сложно выглядит спектр квазипериодического сигнала. Как
уже указывалось выше, аттрактор квазипериодического движения представляет собой тор размерности d. Это означает, что у функции сущ ествует
d аргументов, по которым функция периодична с соответствую щ ими периодами Тi . В общ ем случае спектр квазипериодического движения может
иметь довольно сложный вид. П росто он выглядит только тогда, когда
сигнал есть суперпозиция периодических функций и спектр в силу линейности преобразования представляет собой сумму спектров отдельных периодических функций. Если квазипериодическая функция есть нелинейная
комбинация периодических функций, то ее спектр содержит комбинационные частоты типа n1n1 + n2n 2 + ... + ndn d . Н а рисунке 2.15 показаны два спектра
квазипериодических сигналов с двумя частотами n1 и n 2 . П ри этом, на
рис.2.15а показан случай, когда отнош ение частот есть величина иррациональная, а на рис.2.15б это отнош ение рационально и равно 2/3. Во втором
случае все пики в спектре соответствуют гармоникам с частотами, кратными разности частот (n 2 - n1 ) . В обоих случаях спектр сигналов остается дискретным.
Н а рисунке 2.16 показан типичный спектр апериодического сигнала.
В отличие от предыдущ их спектров он непрерывен (сплош ной, или заполненный спектр). Н а практике вопрос о принадлежности спектра апериодическому или квазипериодическому сигналу не всегда прост, так как квазипериодический сигнал с больш им числом частот приближается по своему
виду к спектру стохастического сигнала. П редельный вид стохастического
сигнала называю т белым ш умом. Это сигнал с плоским спектром, корреляционная функция которого есть дельта-функция.
2.5 Странный аттрактор
Теперь вернемся к вопросу о том, каким должен быть аттрактор хаотического движения. М ы уже упоминали выше, что первый сценарий перехода к Хаосу был предложен Ландау и представлял собой бесконечную цепочку бифуркаций Хопфа. Такому движению соответствует аттрактор в
виде тора T Ґ . Н о уже система с тремя степенями свободы дает сплош ной
спектр Ф урье, что являетсяпризнаком хаотического движения.
64
Н еобходим аттрактор, который объясняет хаотическое поведение
системы в фазовом пространстве низкой размерности (для определенности
будем иметь в виду трехмерное фазовое пространство, так как известно,
что в трехмерных нелинейных системах возможно сущ ествование хаотических режимов). Соответствую щ ий аттрактор был предложен Рю элем и
Таккенсом в 1971г. и назван странным аттрактором. Эти же авторы предложили и сценарий перехода к турбулентности, состоящ ий в том, что в системе после двух бифуркаций Хопфа (приводящ их к появлению в спектре
двух независимых частот) происходит третья бифуркация, приводящ ая к
возникновению странного аттрактора (и появлению заполненного спектра).
Важнейш им свойством, которым должен обладать аттрактор хаотического движения является чувствительность к заданию начальных условий
(ЧЗНУ). Это означает, что близкие траектории должны расходиться
(должны быть положительные показатели Ляпунова) или, иными словами,
система должна забывать о начальных условиях благодаря наличию малых
возмущ ений.
В то же время нужно помнить, что речь идет о диссипативных системах, в которых объем в фазовом пространстве сокращ ается и объем аттрактора должен быть равен нулю . П отеря памяти о начальных условиях
обеспечивается и сокращ ением объемов, так как независимо от начальных
условий фазовая траектория выходит на аттрактор. Чтобы объем множества точек был равен нулю , его размерность d должна быть меньш е размерности пространства. Следовательно,
d < 3.
И з требования ЧЗНУ следует, что траектории в фазовом пространстве должны расходиться, однако, система является детерминированной, а
это означает, что в каждой точке должно сущ ествовать единственное реш ение и траектории не должны пересекаться (разве что в конечном числе особых точек). С учетом того, что траектория должна занимать конечную область фазового пространства, на плоскости эти два требования совместить
не возможно и мы приходим ко второму ограничению на размерность аттрактора:
65
d > 2.
Таким образом, апериодический (странный) аттрактор должен:
а) притягивать фазовыетраектории изобласти притяжения;
б) удовлетворять требованию ЧЗНУ
Рис. 2.17.
в) иметь дробную размерность (в конкретном случае размерность между двойкой и тройкой).
Отложим вопрос о дробной размерности до следую щ его параграфа и
приведем несколько качественных соображений, касаю щ ихся возможной
структуры аттрактора стакими свойствами.
М оделью возможного построения странного аттрактора является так
называемая подкова Смейла. Эта модель отражает важное свойство странных аттракторов - они всегда содержат в себе элементы растяжения с последую щ им складыванием.
П остроение подковы Смейла иллю стрирует рисунок 2.17. И меется
прямоугольник, который растягивается в 2 раза вдоль оси x и сжимается в
2h раз вдоль оси y . Коэффициент h > 1 и характеризует степень сжатия
площ ади.
Н а втором ш аге вытянутый прямоугольник складывается в подкову и
возвращ ается таким образом в исходную область пространства. П ри этом
66
он занимает не всю исходную область, так как появились пробелы, обусловленныесжатием.
Третий ш аг повторяет первый и так далее. Отметим, что деформацию
можно характеризовать числами (показателями) Ляпунова. Растяжение по
оси x характеризуется положительным показателем l1 = ln 2 , а сжатие по
оси y - отрицательным показателем l2 = - ln 2h .
Рис. 2.18.
Вертикальное сечение полученного объекта в точности воспроизводит так называемое канторово множество, размерность которого будет определена в следую щ ем параграфе. Здесь же отметим только, что в пределе
слабой диссипации (h ® 1 ) размерность подковы стремится к двум (она занимает почти всю плоскость). В пределе сильной диссипации (h ® Ґ ) на
плоскости остаю тся редкие линии и размерность множества стремится к
единице.
Другую попытку представить возможность сущ ествования аттрактора с требуемыми свойствами представляет рисунок 2.18. Н а первом ш аге
происходит разбегание траекторий (обеспечиваю щ ее ЧЗНУ). Н а втором
происходит складывание и на третьем - сворачивание полученной пространственной структуры в «кольцо» таким образом, что сложенная вдвое
растянутая сторона смыкается с начальной недеформированной. Вспоминая, что траектории не должны при этом пересекаться, мы приходим к выводу, что должна образоватьсямноголистная структура.
67
2.6 Ф ракталы
2.6.1 П онятие фрактала
П усть имеется множество точек, расположенных в некотором пространстве размерностью D . Введем сферу радиуса r (гиперсферу, если
D > 3 ) и будем подсчитывать среднее число точек N , попадаю щ их в сферу
при различных ее положениях в пространстве. Естественно рассчитывать
на то, что зависимость числа точек от радиуса сферы будет иметь степенную форму
N (r ) » r d
(2.32)
и размерность множества есть
d=
ln N (r )
.
ln r
(2.33)
Если точки множества расположены на линии, то d = 1 , если они лежат на плоскости, то d = 2 , а если точки занимаю т все трехмерное пространство, то опять же получаетсяобычная (евклидова) размерность d = 3 .
Ф ракталами называю т объекты с нецелой размерностью. П ростейш им примером фрактального множества является канторово множество,
строящ ееся по следую щ ему правилу. Единичный отрезок разбивается на
три равных части и средняя часть
удаляется. Н а втором ш аге каждый
из оставш ихся двух отрезков снова
делится на три части с последую щ им удалением центральных частей. П роцедура повторяется до
бесконечности (рис.2.19). Таким
образом, получается такое множество, что лю бой сколь угодно малый объем области обязательно
содержит точки, этому множеству
Рис. 2.19.
не принадлежащ ие.
Оценим размерность построенного множества по формуле(2.33). И зпроцедуры построения множества
следует, что при каждом увеличении радиуса сферы в три раза, число точек, в нее попадаю щ их, увеличивается вдвое ( r » 3n , N » 2 n ). Следовательно,
68
d=
ln 2
= 0,63 .
ln 3
Это не единственный способ определения фрактальной размерности.
Н аиболее известна так называемая размерность Хаусдорфа-Безиковича.
Она определяется следую щ им образом. П усть N (l ) - наименьш ее число кубов (сфер) с ребром (диаметром) l , которым можно покрыть все точки
множества. Тогда размерность Хаусдорфа - Безиковича есть
D = lim
l® 0
ln N (l )
.
ln(1 / l )
(2.34)
Оценивая размерность введенного выше канторова множества по
(2.34), мы придем к тому же самому результату, что и при вычислениях по
формуле (2.33). Одинаковый результат получается при оценке размерности
однородных фракталов. Н есколько примеров однородных фракталов и по-
Рис. 2.20.
лучаемые для них размерности приведены на рис.2.20. В общ ем случае неоднородных фракталов размерности d и D могут не совпадать, но всегда
d Ј D (см. п.1.6.3).
Объекты с фрактальными свойствами возникаю т в самых различных
приложениях. Одной изпервых практических задач, приведш их к развитию
69
теории фракталов, была задача об определении длины береговой линии.
П роблема состоит в том, что по мере использования карт с более мелким
разреш ением получаемая длина береговой линии все увеличивается и процесс не сходится. Береговая линия является, таким образом, типичным
фрактальным объектом (сравните со структурой снежинки Коха, рис.2.20).
Ф рактальными свойствами обладают облака, кораллы, растущ ие кристаллы, семейства трещ ин при процессах разруш ения и поле диссипации энергии в развитом турбулентном течении.
К фракталам приводят многие математические задачи. П ростейш ий
пример дает задача о границах областей притяжения рациональных отображений комплексной плоскости в себя. Н апример, рассматривается
уравнение
z3 = 1,
Рис. 2.21.
70
имею щ ее три корня (1, - 1 / 2 + i 3 / 2, - 1 / 2 - i 3 / 2) , и используется итерационный метод Н ью тона для его реш ения. Это значит, что для уравнения
f ( z ) = 0 строится последовательностьзначений z n , таких, что
f ( z n ) + ( z n + 1 - z n ) f ў( z n ) = 0 .
В наш ем случае это приводит к выражению
zn - 1
3
z n+ 1 = z n -
3z n
2
.
(2.35)
И терационный процесс(2.35) стартуетсразличных начальныхзначений z 0
на комплексной плоскости и приводит, в конце концов, к одному изтрех
корней уравнения. Задача состоит в том, чтобы построить границу раздела
трех областей притяжения. Такие границы называю тсямножествами Ж ю лиа (задача Ж ю лиа датируется1918 годом !) и обладаю т замечательным
свойством: каждая точка границы разделяет все три области притяжения.
М ножества Ж юлиа строятсяи для логистического уравнения
z n + 1 = z n2 + C ,
для которого показано (М андельброт, 1980г.), что уравнение сущ ествует
только для определенныхзначений C на комплексной плоскости. П риняв за
линию уровня число итераций, необходимых для попадание в e окрестность реш ения и рисуя разные уровни разными цветами, получаю т
живописные картинки, украш аю щ ие многие книги и журнальные статьи.
М ы не приводим их из-забедности черно-белого представления и отсылаем
к соответствую щ им изданиям (см. список рекомендуемой литературы).
Эстетическое наслаждение можно получить и от рассматривания изображений аттракторов динамических систем, примеры которых можно видеть на рисунке 2.21. (Н а рисунке, взятом из книги Г.Ш устера «Детерменированный хаос», показаны примеры странного аттрактора и сечения П уанкаре, полученные при реш ении уравнения для нелинейных осциляторов.)
Вспоминая, что именно размерность аттракторов динамических систем с
хаотическим поведением заставили насобратиться к фракталам, вернемсяк
вопросу о том, как именно можно измерить размерность аттрактора.
2.6.2 Алгоритм вычисления размерности аттрактора
71
Вопрос об измерении размерности аттрактора становится особенно
сложным при попытках обработки экспериментальных данных, когда даже
вопрос о размерности фазового пространства, то есть вопрос о необходимом числе независимых переменных, остается открытым. П одход к реш ению этой задачи дает так называемая теорема Таккенса, суть которой состоит в следующ ем.
П усть имеется динамическая
система (не слиш ком больш ой
размерности N ), описываемая
системой
дифференциальных
уравнений
первого
порядка.
П ринципиально, от системы N
уравнений первого порядка можно
перейти к дифференциальному
уравнению N -ого порядка, содержащ ему N производных , но одной переменной (например, остается переменная X (t ) и ее производные X&(t ), X&&(t ), X&&&(t ), и т.д.). П ри
Рис. 2.22.
представлении дифференциальных
уравнений в конечных разностях
это
соответствует
одновременному
знанию
величин
X (t ), X (t + t ), X (t + 2t ), X (t + 3t ), и т.д., где t - постоянная. Теорема Таккенса утверждает, что каждая переменная системы X (t ) отражает основные свойства этой системы, а аттрактор, построенный в фазовом пространстве переменных X (t ), X (t + t ), X (t + 2t ), X (t + 3t ),...... , сохраняет основные топологические
свойства аттрактора исходной системы.
П рактически, алгоритм вычисления размерности аттрактора строится следую щ им образом. Для измеряемой величины X (t ) выбирается характерное время сдвига t и строится фазовая траектория на p переменных
X (t ), X (t + t ), ............, X (t + ( p - 1)t ) как показано на рисунке 2.22. Эта траектория
состоит из последовательности точек,
каждая из которых определяется в
r
фазовом пространстве вектором X i . В каждую из этих точек помещ ается
гиперсфера радиуса r и вычисляется число точек фазовой траектории, попавш их в пределы этой сферы. Затем вводится функция
C (r ) = lim
m® Ґ
1
m2
m
е H (r -
i , j =1
r
r
| X i - X j |) ,
(2.36)
72
характеризую щ ая среднее число пар точек, попадаю щ их в сферу заданного
радиуса. Здесь H - функция Хевисайда, равная по определению единице при
положительных и нулю при остальныхзначениях аргумента. Ожидая, что
C (r ) » r d ,
строят эту функцию в двойном логарифмическом масш табе и при наличии
в таком представлении прямолинейного участка определяют его наклон,
равный величине d . Отметим, что степенной закон можно ожидать только
на масш табах r , заметно меньш их размеров области, занимаемой аттрактором.
П роцедура вычисления величины d повторяется для все возрастаю щ их значений размерности используемого фазового пространства p . П ри
этом вычисленные значения d равны p до тех пор, пока размерность используемого пространства остается меньш ей размерности аттрактора. Если
вычисленная размерность d перестает зависеть от p , то это означает, что
она равна размерности самого аттрактора. Н аименьш ее целое число,
больш ее полученной (фрактальной) размерности аттрактора, называется
размерностью вложения и определяет реальное число степеней свободы
рассматриваемой системы. П ример поведения функции C (r ) по мере роста
p , построенная по результатам реальных измерений в конвекции РелеяБенара (из работы Malraison B. et al., Comptes Rendus Acad.Sc.Paris, 1983,
C297, p.209.) приведена на рис.2.23. В этом примере наклон прямых линий
перестает возрастать с p = 4 , хотя предельный наклон прямых есть 2,8 (то
есть размерность вложения равна трем).
Рис. 2.23.
73
2.6.3 Обобщ енная размерность
П усть система эволю ционирует в некотором фазовом пространстве.
Разобьем это пространство на ячейки (n-мерные кубики) с ребром l (всего
M ячеек) и вычислим вероятность попадания системы в каждую i -тую
ячейку
pi = lim
N® Ґ
ni
.
N
Где ni - число точек, попавш их в данную ячейку, а N - общ еечисло рассмотренных точек.
Обобщ енная размерность (размерность Рени) определяется как
M
1
Dq = lim
l® 0 q - 1
ln е pi
q
i =1
ln l
.
(2.37)
Таким образом вводится последовательность величин D q , связанныхссоответствую щ ими моментами распределения вероятности. П осмотрим, какой смысл имеет эта величина при конкретныхзначениях q .
1) q = 0 . Тогда
M
D0 = lim
l® 0
ln е pi
0
i =1
ln l
сумма в числителе равна числу ячеек, в которых оказалась хотя бы одна
точка. Следовательно,
D0 = lim
l® 0
ln N (l )
,
ln(1 / l )
(2.38)
где N (l ) есть число ячеек, содержащ их точки, и (2.38) совпадает, таким образом, с определением размерности Хаусдорфа (2.34).
2) q = 1 . В этом случае возникает проблема деления на ноль. Рассматриваетсяпредел q ® 1 и с помощ ью правила Лопиталя
74
M
1
D1 = lim
lim
l ® 0 ln l q ® 1
ln е pi
M
q
i =1
q- 1
1
= lim
lim
l ® 0 ln l q ® 1
е
i =1
M
q
pi ln pi
M
е
i =1
pi
q
= lim
l® 0
е
i =1
pi ln pi
ln l
(2.39)
Числитель под знаком предела есть энтропия Ш енона, а размерность
D1 называю т информационной размерностью.
3) q = 2 . Теперь в числителе под знаком суммы стоит квадрат вероятности попадания точки в ячейку, то есть совместная вероятность одновременного попадания пары точек. Таким образом,
M
D2 = lim
l® 0
ln е pi
i =1
ln l
2
= lim
l® 0
ln С (l )
,
ln l
(2.40)
где С (l ) есть функция (2.36), а размерность (2.40) называется корреляционной размерностью .
Справедливо общ ее правило: Di і D j , если i < j . Это означает,
что наибольш еезначение всегда имеет Хаусдорфова размерность D0 .
2.7 Субгармонический каскад
Рис. 2.24.
В этом параграфе речь пойдет о
переходе к хаотическому движению по
сценарию , называемому субгармоническим каскадом и представляю щ ему собой последовательность бифуркаций
удвоения периода. М ы уже упоминали
бифуркацию этого типа, разбирая возможные типы потери устойчивости траектории при анализе матрицы Ф локе.
Качественно перестройку фазовой траектории, соответствую щ ую бифуркации
удвоения периода, иллю стрирует рисунок 2.24. П редельный цикл после бифуркации замыкается только на втором
витке, удваивая тем самым период движения системы в фазовом пространстве.
75
П ри это в сечении П уанкаре число точек удваивается, а в спектре Ф урье
появляетсяновая частота, вдвое меньш ая той, что была до бифуркации.
П рекрасной иллю страцией свойств субгармонического каскада является работа Ф ейгенбаума «Универсальное поведение квадратичных отображений» (Feigenbaum M.J., The universal properties of nonlinear transformations, J.Stat.Phys., 1979, V.21, P.669.), содержание которой мы в основном
и постараемсяпересказать.
Рассмотрим отображение первого возвращ ения
(2.41)
x k + 1 = f ( x k ) = 4mx k (1 - x k )
где x О [0,1] и 0 Ј m Ј 1 . Отображение ставит в соответствие каждой точке из
интервала [0,1] другую точку из этого же интервала. m - управляю щ ий параметр.
П ри m < 0,25 сущ ествует только одна точка, в которой x k + 1 = x k . Это
точка x = 0 и она устойчива. Действительно,
f ў( x) = 4m (1 - 2 x)
и
f ў(0) = 4m .
Это означает, что при m < 1 / 4 производная в точке пересечения функции
Рис. 2.25.
Рис. 2.26.
76
f (x) с биссектрисой x k + 1 = x k
остается
меньш е единицы, что обеспечивает устойчивость реш ения (см. рис.2.25).
П ри 0,25 < m < 0,75 реш ение x = 0
становится неустойчивым, но появляетсядругое реш ение
x* = 1 -
которое устойчиво,
0,25 < m < 0,75
1
,
4m
так как при
| f ў( x * ) |= 2 | 1 - 2m |< 1 .
П уть, по которому реш ение выходит в
этом случае на устойчивую точку, показан на рис.2.26. В точке m = m1 = 0,75 и
Рис. 2.27.
эта точка становится неустойчивой.
Характер возникаю щ его реш ения иллю стрирует рисунок 2.27, где показано реш ение для m = 0,8 . В реш ении возникают две неподвижные точки. Это так называемый 2-цикл, при котором
реш ение возвращ ается в данную точку через ш аг. И наче говоря, реш ение
определяется условием: x k + 2 = x k . Запиш ем
xk + 2 = f ( xk + 1 ) = f 2 ( xk ) = g ( xk ) ,
где явный вид функции g есть
Рис. 2.28.
Рис. 2.29.
77
g ( x) = 16m 2 ( x - x 2 - 4mx 2 + 8mx 3 - 4mx 4 ) .
График этой функции показан на рисунке 2.28, а два выделенных квадрата
поясняю т тот факт, что в них воспроизводится картинка, представленная
на рис.2.26. В дальнейш ем все повторяется. Ф ункция g (x) теряет устойчивость при m = m 2 = (1 + 6 ) / 4 = 0,86237... Далеерассматриваетсяфункция
h( x ) = g 2 ( x ) = f 4 ( x ) ,
представленная на рис.2.29. Квадрат
на рисунке снова показывает, что
вблизи каждой устойчивой точки
воспроизводится ситуация рисунка
2.26. Ф ункция h(x) становится неустойчивой при m = m 3 = 0,875 , и т.д. Каждый раз имеет месть бифуркация удвоения периода (период цикла удваивается). Ф ейгенбаум обнаружил два
закона подобия, характеризую щ их
субгармонический каскад. Во первых,
он показал, что последовательность
mi быстро сходится
Рис. 2.30.
m Ґ = 0,892486418... ,
и сущ ествует предел
lim
i® Ґ
mi - mi- 1
= d.
mi+ 1 - mi
Важно, что величина d не зависит от конкретного вида функции f (x) (лю бая выпуклая, непрерывная, дифференцируемая функция с одним максимумом) и равна
d= 4,6692016091....
Это первый закон подобия. Второй закон подобия касается положения устойчивых точек. Н а рисунке 2.30 схематически показана структура реш ений уравнения (2.41). Рассматриваю тся точки пересечения с прямой x = 0,5
до ближайш ей к ней точки на устойчивом 2 n -цикле. Для соответствую щ их
расстояний d n справедливо соотнош ение
78
lim
n® Ґ
dn
=- a
d n+ 1
и вторая константа Ф ейгенбаума
a = 2,5029078750... .
Отметим, что в результате бесконечной последовательности бифуркаций удвоения при m = m Ґ возникает бесконечное множество (аттрактор
Ф ейгенбаума), который имеет размерность Хаусдорфа D = 0,548... . Важно,
что при всех m < m Ґ показатель Ляпунова отрицателен, стремясь при m ® m Ґ
к нулю. Следовательно, аттрактор Ф ейгенбаума неявляетсястранным.
Хаос возникает при m > m Ґ , где показатель Ляпунова в основном
положителен. П оведение в этой области достаточно сложное. Хаотические
области чередую тсяс«окнами периодичности» (светлыезоны на рис.2.31).
Рис. 2.31.
79
2.8 Н екоторые примеры
2.8.1 Система Лоренца
Рассмотрим подробно свойства системы Лоренца, полученной ранее
в параграфе 1.5 как пример маломодовой модели конвекции в подогреваемом снизу слое жидкости. И меем систему (1.35)
X& = s (Y - X ),
Y& = - XZ + rX - Y ,
Z& = XY - bZ .
(2.42)
Н апомним, что управляю щ им параметром является относительное
число Рейнольдса r , а число П рандтля и параметр b для определенности во
всех случаях, когда будут обсуждаться численные результаты, будем полагать s = 10, b = 8 / 3 .
Уравнения (2.42) имею т тривиальное реш ение X 0 = Y0 = Z 0 = 0 , отвечаю щ ее отсутствию конвекции. П роверим это реш ение на устойчивость.
Для этого представим всетри переменные в виде
X = X 0 + xe - lt , Y = Y0 + ye - lt , Z = Z 0 + ze - lt ,
(2.43)
считая x, y, z - малыми возмущ ениями. (2.43) подставляем в (2.42) и
Рис. 2.32.
80
отбрасываем нелинейные по малым возмущ ениям члены. В результате, после сокращ ения на экспоненты, получаем линейную алгебраическую систему
( l - s ) x + sy = 0,
rx + ( l - 1 ) y = 0,
( l + b ) z = 0.
Реш ая задачу на собственные значения, приравниваем нулю определитель системы и получаем
(s + 1) 2 - 4s (1 - r )
s+ 1
l=
±
.
2
2
Видно, что при r > 1 один из двух корней становится отрицательным,
то есть в точном соответствии с результатом Релея (иначе и быть не может)
при r = 1 возникаетконвективное движение.
Система (2.42) имеет и нетривиальное реш ение
X = Y = ± b(r - 1) ,
(2.44)
Z = r - 1.
У переменных X и Y действительная часть появляется при r > 1 . Таким образом, в точке r = 1 имеет место нормальная бифуркация вилки и появляется два устойчивых реш ения, соответствую щ их стационарной валиковой конвекции с противоположным направлением вращ ения конвективных
валов.
П овторяя линейный анализ устойчивости для реш ения (2.44), приходим к кубическому уравнению
Рис. 2.33.
81
Рис. 2.34.
l3 - (s + b + 1)l2 + (r + s )bl - 2sb(r - 1) = 0 ,
в одном из корней которого появляется отрицательная действительная часть при
r=
s (s + b + 3)
.
s- b- 1
П ри s = 10, b = 8 / 3 это выражение дает значение r = 24,74 . В этой
точке имеет место субкритическая
бифуркация Хопфа. Особенность
поведения системы Лоренца в том,
что устойчивый предельный цикл
не возникает в ней вовсе (напомним, что согласно сцерарию Рю эляТаккенса, странный аттрактор возникает после двух бифуркаций
Хопфа) и странный аттрактор возникает сразу после первой (обратной) бифуркации Хопфа. Бифуркационная диаграмма представлена
на рисунке 2.32. Следует отметить,
Рис. 2.35..
82
что «чистый» странный аттрактор сущ ествует в небольш ом интервале числа Релея 24,06 < r < 30,1 . Обратим внимание и на то, что на левом краю этого
интервала сущ ествует гистерезис - при понижении числа Релея странный
аттрактор сущ ествует до r = 24,06 , а не до r = 24,74 . Левее этой границы в
интервале чисел Релея r > 13,93 сущ ествует область так называемого метастабильного хаоса. В этой области малые возмущ ения стационарного реш ения монотонно затухают, но больш ие
возмущ ения приводят к хаотическим
режимам, которые в конечном итоге
также затухают, но успеваю т при этом
выписать в фазовом пространстве
многочисленные хаотические петли,
напоминаю щ ие поведение системы на
странном аттракторе. П ри r > 30,1 диаграмма режимов представляет собой
чередование областей с хаотическим и
периодическим движениями, напоминая поведение отображения Ф ейгенбаума в области m Ґ < m < 1 (рис.2.31).
Рис. 2.36.
П оявлению области с периодическим
аттрактором предш ествует обратный
каскад, а само «окно периодичности» включает субгармонический каскад.
Число «окон периодичности», по-видимому, бесконечно и при больш их
числах Релея их ш ирина растет. П оследнее окно неограниченно и занимает
всю область r > 214,364 .
В своей знаменитой работе Лоренц численно исследовал поведение
системы при r = 28 . Н а рисунке 2.33 показан фрагмент поведения во времени переменной X (t ) при этом значении r , а на рис.2.34 - характерный вид
фазовой траектории системы на странном аттракторе. Н а рис.2.35 - проекции фазовой траектории на плоскости ( X , Z ) .
Н аблю дение за эволю цией фазовой траектории показывает, что траектория описывает
витки вокруг точек, соответствую щ их ставш им
неустойчивыми реш ениям (2.44), переходя случайным образом от вращ ения вокруг одного
фокуса к вращ ению вокруг другого.
Н аблю дая за эволюцией фазовой траектории в плоскости ( X , Z ) , Лоренц сделал важный вывод. Траектория раскручивается воРис. 2.37.
круг одного фокуса, увеличивая на каждом
витке радиус орбиты. Этот процесс происхо-
83
дит до тех пор, пока на очередном витке в точке максимума траектория не
выйдет за значение Z = 38,5 . Как только траектория превысит это значение,
она уходит в область притяжения другого фокуса и все повторяется вновь.
П ри этом число витков, которое соверш ит траектория, зависит от величины превышения траектории над этим критическим значением перед перебросом. Лоренц использовал метод точечных отображений, позволяю щ ий
перейти от системы с непрерывным временем к системе с дискретным временем - вариант сечения П уанкаре, называемый отображением первого
возвращ ения. В качестве отображения использовалось значение величины
Z в текущ ем локальном максимуме, как функция от значения в предыдущ ем
максимуме (рис.2.36). Левая, восходящ ая часть функции соответствует
процессу раскручивания, а переход за пик - перебросу к другому фокусу.
Лоренц предложил простейш ую модель наблю даемого процесса - отображение отрезка [0,1] на себя вида (рис.37)
M n+ 1
м
2M
п
п n
=н
п2(1 - M )
n
п
о
1
2
1
Mn >
2
Mn <
(2.45)
Если рассматривается последовательность, начинаю щ аяся со значения M 0 , то она будетразвиватьсяпо следую щ ей цепочке:
м2M 0
M1 = н
о2 - 2M 0
м4M 0
п2 - 4M
п
0
M2 = н
п4 - 4M 0
п
о- 2 + 4M 0
м8M 0
п2 - 8M
0
п
п4 - 8M 0
п
п6 - 8M 0
M3 = н
п8 - 8M 0
п- 2 + 8M 0
п
п- 4 + 8M 0
п- 6 + 8M
0
о
....... M n = m n ± 2 n M 0 .
Здесь m n - четное число, такое, что оно сдвигает величину 2 n M 0 в интервал [0,1]. Все возможные последовательности можно разделить на три
типа:
П оследовательности, заканчиваю щ иеся в нуле. Таких последовательностей счетное множество и они начинаю тся с элемента видаM 0 = u 2 p , где
u - нечетное целое число. Тогда M p - 1 = 1 2 и M p = 0 .
84
П ериодические последовательности. Они возникают, если M 0 = u 2 p v ,
где u,v - простые числа. Тогда M p + 1+ k
2 p + 1+ k u
2 Ч2 k u
= m±
= m±
. П ростеш ие приv
2pv
меры получаю щ ихсяпериодических последовательностей есть
[2 3, ].....
[2 / 5, 4 / 5, ].....
[2 / 7, 4 / 7, 6 / 7, ]....
[2 / 9, 4 / 9, 8 / 9, ]....
3) Апериодические последовательности.
Эта модель иллю стрирует ещ е одно важное свойство системы - неустойчивость к малым возмущ ениям (ЧЗНУ). Действительно, если рассмотреть последовательность с малым возмущ ением начального элемента
M 0ў = M 0 + e , то после n итераций
M nў = m n ± 2 n ( M 0 ± e) = M n ± 2 n e ,
что свидетельствуетоб экспоненциальном ростевозмущ ений.
Отметим, что модельное отображение (2.45) при всей своей простоте
сохраняет важнейш ее свойство, приводящ ее к ЧЗНУ в диссипативных системах - это растяжение в сочетании со складыванием. Растяжение на каждом ш аге приводит к экспоненциальному росту начального смещ ения (расхождению траекторий), а складывание обеспечивает возвращ ение в ограниченную область (в данном случае интервал).
2.8.2 М одель динамо
Рикитаке
Другой пример динамической
системы со стохастическим поведением дает так называемая модель двухдискового динамо Рикитаке, предложенная в связи с задачей об инверсиях
геомагнитного поля. М агнитное поле
Земли в первом приближении представляет собой диполь, который по
палеомагнитным данным многократ-
Рис. 2.38.
85
но и нерегулярно менял свою полярность. Н а сегодняш ний день ш кала полярности геомагнитного поля восстановлена более чем за 1700 миллионов
лет, что составляет порядка половины возраста Земли. За это время зарегистрировано 593 переброса магнитного поля, причем время между двумя перебросами колеблется в интервале от 10 тысяч до сотен миллионов лет, демонстрируя хаотическое поведение, лиш енное каких-либо периодичностей.
Согласно принятой на сегодня точке зрения, магнитное поле Земли
возбуждается в результате конвективного движения в жидком (электропроводящ ем) ядре. П роцесс возбуждения магнитного поля в движущ ейся проводящ ей среде получил название М ГД-динамо. Земное динамо представляет собой сложный нелинейный магнитогидродинамический процесс, исследование которого находится лиш ь на начальной стадии. Больш ой интерес
представляют поэтому любые упрощ енный модели процесса генерации
магнитного поля, способные приводить к случайным сменам полярности
генерируемого магнитного поля.
Самые простые модели оперируют не потоками проводящ ей жидкости, а движущ имися проводниками. П ервая попытка построить такого рода модель принадлежит Булларду (Bullard E.C., Proc.Cambridge Philos.
Soc.,1955, v.51, p.744.), который предложил однодисковое динамо, но такая
модель не дает смены полярности генерируемого поля. Рикитаке (Rikitake
T., Proc.Cambridge Philos. Soc.,1958, v.54, p.89.) рассмотрел систему двух
дисковых динамо, связанных таким образом, что ток от одного диска питает катуш ку возбуждения другого и наоборот. Эта ситуация изображена на
рис.2.38. Оба диска вращ аю тся без трения и находятся под действием одинаковых моментов сил G , компенсирую щ их омические потери в дисках и
обмотках. Уравнения, описываю щ ие эволю цию токов I 1 , I 2 и угловых скоростей W 1 , W 2 можно записать в виде
LI&1 + RI 1 = MW 1 I 2 ,
LI&2 + RI 2 = MW 2 I 1 ,
CW& = G - MI I ,
1
(2.46)
1 2
CW& 2 = G - MI 1 I 2 ,
где L - коэффициент самоиндукции, R - сопротивление каждой цепи, M коэффициент взаимоиндукции, C - момент инерции диска.
Два последних уравнения (2.46) показываю т, что разность угловых
скоростей есть величина постоянная
W1 - W 2 =
GL
A,
CM
где A - константа. Это позволяет перейти к системетрех уравнений.
86
Система записывается в безразмерном виде. П ри этом за единицу тока принимаю т величину G / M , угловой скорости - GL / CM , а за единицу
времени - величину t et m . Единица времени выражена через два характерных масш таба времени, присущ их системе. Это время t m , за которое диск
под действием приложенного момента сил разгоняется до характерной
скорости R / M ,
tm =
CR
GM
и время электромагнитной диффузии
te =
L
,
R
характеризую щ ее время вырождения магнитного поля при остановке диска. И х отнош ение является безразмерным параметром системы
m=
tm
CR 2
=
.
t e GLM
Обозначая безразмерные токи как X i , а безразмерные угловые скорости
как Yi (в уравнениях остается одна переменная Y , так как Y1 - Y2 = A ), прихо-
Рис. 2.39.
87
дим к системе
X&1 + mX 1 = YX 2 ,
X&2 + mX 2 = (Y - A) X 1 ,
Y& = 1 - X X .
1
(2.47)
2
Система (2.47) имеетстационарные реш ения
X 1 = ±K ,
X 2 = ± K - 1 , Y = Y1 = mK 2 , Y2 = mK - 2 ,
где A = m ( K 2 - K - 2 ).
М ы не будем подробно описывать свойства системы Рикитаке, оставляя ее изучение для самостоятельных работ. Н а рис.2.39 показана только
фазовая траектория системы для случая m = 1,5; K = 2. М ожно видеть, что ее
топология близка аттрактору Лоренца.
2.8.3 Реальная конвекция
Н аибольш ее число экспериментальных работ по исследованию перехода от упорядоченных течений к хаотическим выполнено, пожалуй, в исследованиях конвективных течений. М ы приведем некоторые результаты
исследований перехода от ламинарного движения к турбулентности при
конвекции в кубической полости, взятые из работы: Зимин В.Д., Кетов
А.И . Н адкритические конвективные движения в кубической полости.
И зв.АН СССР, М еханика жидкости и газа, 1974, N.5, С.110.
Рис. 2.40.
88
Рис. 2.41.
Рис. 2.42.
И змерения проводились подогреваемой снизу в кубической полости с
ребром 40 мм, образованной медными стенками. Горизонтальные стенки
термостатировались, обеспечивая заданную разность температуры, а вертикальные обеспечивали равновесный однородный градиент температуры.
Н адкритические течения, возникаю щ ие в кубической полости и имею щ ие
наиболее низкие уровни устойчивости, схематически показаны на рис.2.40,
где стрелками показано направление движения жидкости в верхней части
полости, а знаками «плю с» и «минус» обозначены области, в которых температура
оказывается выше или ниже средней. Критические числа Релея для движений типа А
и Б равны 8224, для В - 9184 и для Г 14032. В полости были установлены дифференциальные термопары, расположенные таким образом, что их показания позволяли выделять движения всех четырех
типов.
Н е останавливаясь на сценариях
развития неустойчивости и переходов от
одного режима движения к другому, приведем лиш ь некоторые данные, иллю стрирую щ ие поведение системы в одночастотном режиме, двухчастотном и стохастическом режимах. Для каждого из трех режимов на рисунках представлены изменения
во времени показаний термопар, соответствую щ их каждому из выделяемых течеРис. 2.43.
ний, проекции фазовых траекторий на
89
Рис. 2.44.
Рис. 2.45.
плоскости, образованные всеми парами термопар и спектры мощ ности
пульсаций температуры, регистрируемой каждой из четырех термопар.
Рисунки 2.41-2.43 относятся к одночастотному режиму, регистрируемому при числе Релея R = 2 Ч10 5 . П ервый рисунок показывает характер колебаний показаний всех четырех термопар, второй - соответствую щ ие этим
колебаниям проекции фазовых траекторий, ясно указываю щ ие на сущ ествование предельного цикла. Об этом
же свидетельствую т и спектры Ф урье
(рис.2.43) состоящ их их одного главного пика на частоте 0,054 Гц и пика на
удвоенной частоте, обусловленный негармонической формой колебаний.
Следую щ ая группа рисунков
представляет результаты для числа Релея R = 2,24 Ч10 5 . Н а рисунке 2.44 показаны пульсации показаний термопар,
на рис.2.45 - соответствую щ ие фазовые
траектории (за время соответствую щ ее
периоду низкочастотных колебаний), а
на рисунке 2.46 - спектры, свидетельствую щ ие о сущ ествовании двухчастотного режима (частоты 0,0451 Гц и 0,304
Гц).
Движение становится стохастическим при R = 2,50 Ч10 5 . П оказания
термопар для этого режима представРис. 2.46.
лены на рис.2.47, фазовые траектории -
90
на рис.2.48, а спектры мощ ности - на рис.2.49. Видно, что фазовые траектории имеют чрезвычайно запутанную структуру, а спектры становятся
сплош ными, сохраняя лиш ь слабые локальные максимумы, свидетельствую щ ие о сохранении периодических составляю щ их.
Рис. 2.47.
Рис. 2.48.
Рис. 2.49.
91
Рекомендуемая литература ко второй главе:
1. БержеП ., П омо И ., Видаль Л. П орядок в хаосе. М осква: М ир. 1991.
366с.
2. Ш устер Г. Детерминированный хаос. М осква: М ир. 1988. 240с.
3. Странные аттракторы. Сборник статей. Серия «М атематика. Н овоев зарубежной науке», выпуск 22. М осква: М ир. 1981. 254с.
92
3 ПОЛУЭМ П И РИ Ч ЕСК И Е М ОДЕЛИ
3.1 Развитая турбулентность
3.1.1 Вводныезамечания
В данной главе мы начинаем рассматривать подходы к описанию
развитой турбулентности, то есть течений, возникаю щ их при значительном
превышении критических значений управляю щ их параметров (числа Рейнольдса, если речь идет об изотермическом течении в отсутствии дополнительных силовых полей). Такие течения характеризую тся наполненными
спектрами Ф урье, причем нетолько временными, но и пространственными.
Н апомним ещ е раз, что именно в этом и есть основное отличие турбулентности от хаоса в динамических системах невысокого порядка: в турбулентном потоке хаос и пространственный, и временной, а хаотическое поведение маломодовых систем (соответствую щ их например конвективным течениям при невысокой надкритичности) представляет собой хаотическую во
времени эволю цию мод с относительно простой пространственной структурой.
П риступая к рассмотрению развитых турбулентных течений, следует
сделать ряд важных замечаний. П ервое из них касается уравнений движения жидкости. В первой главе мы получили уравнения Н авье-Стокса, как
основные уравнения, с помощ ью которых мы описываем в дальнейш ем все
течения жидкости. Снова подчеркнем, что мы действительно продолжаем
считать, что эти уравнения описываю т течения жидкости и в турбулентном
режиме, даже при экстремально больш их значениях безразмерных параметров (более того, мы будем рассматривать только случай несжимаемой
жидкости). Уверенность в том, что это возможно, держится на результатах
многочисленных успеш ных попыток использования этих уравнений для
турбулентных течений. Сама возможность приложения уравнений Н авьеСтокса к турбулентности совсем не очевидна (и продолжает подвергаться
критике), так как при их выводе было сделано достаточно сильное предположение о том, что тензор вязких напряжений включает в себя только линейные комбинации первых производных поля скорости. В ламинарных и
слабо надкритических течениях это предположение кажется разумным и
прекрасно работает, но в сильно нелинейных режимах нельзя исклю чить,
что тензор вязких напряжений будет иметь более сложную зависимость от
структуры поля скорости. Оправданием использованию уравнений движения в принятой формеможетслужить только сопоставление результатов их
реш ения с экспериментальными данными.
93
Далее, пусть уравнения движения справедливы и предположим, что
мы располагаем мощ нейш им компью тером, способным реш ать трехмерные
уравнения движения с лю бой желаемой точностью (например, будем считать трехмерный поток на сетке 1000х1000х1000). Это, однако, не снимает
проблемы описания турбулентности, так как в результате такого реш ения
мы будем иметь огромное количество информации, осознание которой
требует ее представления в некотором виде, а это фактически опять же
предполагает введение определенной модели процесса. П о сути, такой суперкомпью тер отличается от реального турбулентного течения, наблю даемого в лаборатории или природе, только несравненно больш ими возможностями съема информации относительно состояния потока в лю бой точке
и в любой момент времени.
П роблема описания турбулентного движения состоит в выделении
характеристик, описываю щ их свойства системы с огромным числом степеней свободы, а лю бой подход к ее описанию - это тот или иной способ ограничения числа степеней свободы.
Турбулентные поля (скорость, давление, температура и т.д.) представляю т собой случайные поля. В лю бой точке потока можно установить
датчик и зарегистрировать реализацию процесса в данной точке. М ногократно повторяя эту процедуру, принципиально возможно получить плотr
ность вероятности P ( f ) для интересую щ ей нас величины f (r , t ) . В общ ем
случае, плотность вероятности также есть функция координат и времени.
Сущ ествуетряд важных частных случаев, которыемы и перечислим.
Турбулентность является однородной, если плотность распределения
вероятности независит от сдвига
r
r
r
P (t , r + Dr ) = P(t , r ) .
Турбулентное течение называется стационарным, если плотность вероятности независит от времени, то есть
r
r
P (t + t , r ) = P (t , r ) .
П роцесс называется эргодическим, если осреднение по времени эквивалентно для него осреднению по ансамблю реализаций
T
r
r
1
f (r ) = lim тf (t , r )dt .
T® Ґ T 0
Угловыми скобками будем обозначать среднее по ансамблю реализаций. Очевидно, что только стационарный процесс может быть эргодическим. Гипотеза эргодичности ш ироко используется при исследовании ста-
94
ционарных течений, так как на практике измеряю тся именно средние по
времени величины.
В реальных измерениях ш ироко используется и гипотеза Тейлора, позволяю щ ая связать пространственные и временные флуктуации исследуеr
мой величины f (r , t ) . Согласно этой гипотезе, если сущ ествует среднее теr
чение, характеризуемой скоростью U , то справедливо соотнош ение
¶f
¶f
= Ui
.
¶t
¶xi
П ользуясь этой гипотезой, по измерениям в заданной точке пространства определяют пространственные флуктуации исследуемого поля и
их статистические характеристики.
3.1.2 Статистическиемоменты случайных полей
r
Ф ункция распределения плотности вероятности P (r , t ) содержит полr
ную информацию о случайном поле f (r , t ) , однако, ееопределение в полном объеме практически невозможно. И звестно, что заданию плотности
вероятности эквивалентно задание последовательности (в принципе- бесконечной) статистических моментов
Mf m = тf m P ( f )df .
П ри этом момент нулевого порядка равен единице в силу условия нормировки
Mf 0 = тP ( f )df = 1 ,
а момент первого порядка, называемый такжематематическим ожиданием,
даетсреднее значение величины
Mf 1 = тfP ( f )df = f .
Для моментов второго и болеевысоких порядков обычно использую т
центральныемоменты, вычисляемые относительно средних значений
M(f -
f ) m = т( f -
f ) m P ( f )df .
Н апомним, что центральный момент второго порядка называется
дисперсией.
95
С точки зрения описания турбулентных полей, необходимы статистиr
ческие характеристики связи между значениями величины f (r , t ) в различных точках пространства. Это требует введения совместной плотности ве-
Рис.3.1.
r
r
роятности P ( f (r1 ), f (r2 )) и (или) соответствую щ их двухточечных моментов.
Важнейш им среди двухточечных моментов является момент второго порядка, называемый корреляционной функцией
r r
r
r
r
r
B (r1 , r2 ) = т( f (r1 ) - f1 )( f (r2 ) - f 2 ) P ( f (r1 ), f (r2 ))df df 2 = ( f1 - f1 )( f 2 - f 2 ) . (3.1)
1
Если речь идет о векторном поле(например, скорости), то появляется
корреляционный тензор
r r
r
r
r
r
Bij (r1 , r2 ) = (vi (r1 ) - vi (r1 )(v j (r2 ) - v j (r j ) ) .
(3.2)
Для однородной турбулентности (3.1) и (3.2) зависят только от взаr r r
имного расположения двух точек, то есть, если r2 = r1 + r , то
r r
r
Bij (r1 , r2 ) = Bij (r ) .
(3.3)
Важным частным случаем является однородная и изотропная турбулентность, в которой совместная плотность вероятности (а, следовательно,
r
и двухточечныемоменты) не зависят и от направления вектора r . Тогда
r r
r
Bij (r1 , r2 ) = Bij (| r |) = Bij (r ) .
(3.4)
96
Чащ е всего использую т корреляционные функции Bll (r ) и Bnn (r ) , характеризую щ ие корреляцию продольных и поперечных составляю щ их
пульсаций скорости. Здесь индексом l обозначена составляю щ ая скорости
r
r
вдоль линии, соединяю щ ей точки r1 и r2 , а индексом n составляю щ ая, нормальная этой линии. Характерный вид этих функций иллю стрирует рисунок 3.1.
Выш е, в параграфе 2.4.3, указывалось на связь корреляционной
функции со спектрами (теорема Хинчина) в случае временного сигнала.
Аналогичное соотнош ение связывает и пространственные спектры с двухточечными корреляционными функциями. П режде чем написать это соотнош ение, остановимся несколько подробнее на вопросе о пространственных спектрах турбулентности.
3.1.3 П ространственныеспектры
П редположим, что рассматриваемое случайное (турбулентное) поле
занимает ограниченный объем и величина f ( x, y, z, t ) может быть представлена интегралом Ф урье
Ґ
rr r
) r
r
1
f ( r ,t ) = 3 тf ( k ,t )e ik r dk ,
8p - Ґ
(3.5)
где
) r
f ( k ,t ) =
Ґ
rr
r
- ik r r
f
(
r
,
t
)
e
dr ,
т
(3.6)
-Ґ
r
r
r = ( x, y, z ) - радиус-вектор, k = (k x , k y , k z ) - волновой вектор.
Считая рассматриваемую турбулентность стационарной, определим
трехмерный энергетический спектр случайного поля :
r
) r
F ( k ) =< | f ( k ) |2 >
(3.7)
Угловыескобки означают в этом случае осреднение по времени. Трехмерr
ный спектр связан с корреляционной функцией B (r ) (теорема Хинчина)
r
r rr r
1
F (k ) = 3 тB (r )e - ik r dr
8p
(3.8)
В теории турбулентности, говоря о ее спектральных свойствах, обычно
имею т в виду энергетический спектр E (k ) , который характеризует энергию
97
всех гармоник с заданным модулем волнового вектора, независимо от его
направления.
r r
E (k ) = тF (k )dk ,
(3.9)
r
|k |
или, в сферической системе координат,
2p p
E (k ) =
r
ттF (k )k
2
sin J dJ dj .
0 0
r
В важном частном случае изотропной турбулентности, когда F (k ) = F (k ) ,
связь становитсяочень простой:
E (k ) = 4pk 2 F (k ) .
(3.10)
Отметим, что все оценки для спектральных законов развитой турбулентности касаются обычно именно энергетического спектра E (k ) .
Если в турбулентном потоке измерения проводятся вдоль одной прямой, то по этим измерениям можно построить одномерное фурьепреобразование. Ограничиваясь однородной и изотропной турбулентностью , в которой все прямые равноправны, рассмотрим прямую y = z = 0 и
запиш ем
+Ґ
)
f1 ( k x ) = тf ( x , y , z )e - ixk x dx .
-Ґ
Квадрат модуля этой величины есть одномерный энергетический спектр
)
F1 ( k x ) =| f1 ( k x ) |2 .
(3.11)
Чтобы получить связь между одномерным и трехмерным спектрами,
выразим исходную величину на прямой y = z = 0 через обратное преобразованиеФ урье. С одной стороны
f ( x ,0,0 ) =
1
2p
)
тf ( k
1
x
)e ixk x dk x ,
а с другой стороны
1 )
i ( xk + 0 k + 0 k )
f ( k x ,k y ,k z )e x y z dk x dk y dk z =
3 т
8p
) r
1
ж 1
ц
=
f
( k )dk y dk z ч e ixk x dk x .
з
2 т
т
2p
и 4p
ш
f ( x ,0,0 ) =
98
Таким образом,
)
1
f1 ( k x ) =
4p 2
а
F1 (k x ) =
) r
f
т ( k )dk y dk z ,
r
1
F
(
k
)dk y dk z .
16p 4 т
В следую щ их главах, рассматривая структуру мелкомасш табной турбулентности, мы постоянно будем обращ аться к спектрам, описываемым
степенными законами. П окажем, как связаны между собой введенные спектры турбулентности при степенной зависимости энергии от масш таба (волнового числа). П усть имеется однородное изотропное поле скалярной величины, энергетический спектр которой следует степенному закону
E (k ) ~ k a .
Тогда трехмерный спектр
F (k ) ~ k
a- 2
2
= (k x + k y + k z )
2
2
a- 2
2
,
а одномерный
F1 (k x ) ~
тт(k
ж
a- 2
= k x ттз1 +
з
и
= kx
a
тт(1 + h
2
x
жk y
з
зk
и x
2
+ ky + kz )
2
2
2
2
ц жk z ц
ч
з ч
ч
ч + з
ш иk x ш
+ x2
)
a- 2
2
a- 2
2
dk y dk z =
a- 2
2
ц
ч
ч
ш
dhdx ~ k x
dk y dk z =
a
(проведена замена переменных h = k y / k x ; x = k z / k x ).
Таким образом, в однородной изотропной турбулентности энергетический спектр E (k ) и одномерный спектр F1 (k ) следую т одному степенному
закону, а степень убывания трехмерного спектра меньш е на двойку (т.е.
трехмерный спектр значительно круче).
99
3.2 Уравнения для статистических моментов
3.2.1 Уравнение Рейнольдса
Рассмотрим уравнения Н авье-Стокса в тензорных обозначениях
¶t vi + v j ¶ j vi = - r - 1¶i p + n¶2jj vi + f i ,
(3.12)
(3.13)
¶k v k = 0 .
Входящ ие в них величины представим в виде сумм средних полей и
пульсаций:
r
r
r
v i (r , t ) = U i ( r , t ) + u i ( r , t ) ,
r
r
r
p(r , t ) = P (r , t ) + p ў(r , t )
(3.14)
П ри этом, согласно принятым определениям, предполагаются следую щ ие правила осреднения (угловые скобки по-прежнему обозначают осреднение по ансамблю реализаций):
vi = U i ,
p = P,
Ui = Ui ,
P = P,
(3.15)
(3.16)
u i = 0;
p ў = 0;
Разложения (3.14) подставим в исходные уравнения (3.12)-(3.13):
¶tU i + ¶t u i + U j ¶ jU i + U j ¶ j u i + u j ¶ jU i + u j ¶ j u i = - r - 1 (¶i P + ¶i p ў) + n (¶2jjU i + ¶2jj u i ) + f i
(3.17)
¶k U k + ¶k u k = 0 ,
(3.18)
и проведем осреднение
¶t U i + ¶t u i + U j ¶ j U i + U j ¶ j u i + u j ¶ j U i + u j ¶ j u i =
- r - 1 (¶i P + ¶i p ў) + n (¶ 2jj U i + ¶2jj u i ) +
fi
¶k U k + ¶k u k = 0 .
Учитывая правила осреднения (3.15)-(3.16), приходим к уравнению
Рейнольдса:
¶tU i + U j ¶ jU i = - r - 1¶i P + n¶2jjU i - ¶ j u j u i +
fi
,
(3.19)
и уравнению неразрывности длясреднего поля скорости
100
¶k U k = 0 .
(3.20)
В уравнение Рейнольдса для средних полей входит одноточечный
корреляционный тензор пульсаций скорости, называемый тензором напряжений Рейнольдса
t ij = u i u j .
(3.21)
Этот тензор нельзя выразить через осредненные характеристики турбулентных полей. Следовательно, число неизвестных превышает число
имею щ ихся уравнений и система (3.19)-(3.20) являетсяне замкнутой.
3.2.2 Ц епочка уравнений Ф ридмана-Келлера и проблема
замыкания
В уравнении Рейнольдса появилась новая неизвестная величина - тензор напряжений Рейнольдса (3.21), для которого также можно получить
эволю ционноеуравнение. Так как
¶tt ij = ¶t u i u j = u i ¶t u j + u j ¶t u i ,
то сначала требуется получить уравнение для пульсаций скорости,
для чего из уравнения (3.17) необходимо вычесть уравнение (3.19). П олучим
(немые индексы j заменены на k )
¶t u i + U k ¶k u i + u k ¶k U i + u k ¶k u i = - r - 1¶i p ў- ¶k u i u k + n¶2kk u i + f iў.
(3.22)
Аналогичное уравнение получается и для компоненты u j :
¶t u j + U k ¶k u j + u k ¶k U j + u k ¶k u j = - r - 1¶ j p ў- ¶k u j u k + n¶2kk u j + f jў.
(3.23)
Уравнение (3.22) умножается на u j и складывается с уравнением
(3.23), умноженным на u i :
u i ¶t u j + u j ¶t u i =
- U k ¶k (u i u j ) - u j u k ¶k U i - u i u k ¶k U j - u j ¶k (u i u k ) - u i ¶k (u j u k ) - u i ¶k u j u k - u j ¶k u i u k
2
- r - 1 (u i ¶ j p ў+ u i ¶ j p ў) - n (u i ¶kk
u j + u j ¶2kk u i ) + u i f jў+ u j f iў
101
П осле осреднения приходим к уравнению:
(
)
¶t u j u j + U k ¶ k u i u j = - u i u k ¶ k U j + u j u k ¶ k U i - ¶ k u i u j u k
- r - 1 ( u i ¶ j p ў + u j ¶i p ў ) - n ( u i ¶2kk u j + u j ¶ 2kk u i ) + u i f jў + u j f ў .
(3.24)
i
В уравнении для корреляционного тензора пульсаций скорости второго порядка (3.24) появился корреляционный тензор (момент) третьего
порядка u i u j u k и новые моменты второго порядка, описываю щ ие корреляции пульсаций компонент скорости с давлением и скорости со вторыми
производными скорости.
Для вновь появивш ихся статистических моментов также можно написать эволю ционные уравнения типа (3.24), но проблемы это не реш ит,
так как в уравнение для момента третьего порядка войдут момент четвертого порядка и новые моменты третьего порядка и так далее. Система
уравнений для моментов все возрастаю щ их порядков называется цепочкой
уравнений Ф ридмана-Келлера и является незамкнутой в принципе. П роблема обрыва этой цепочки и получения замкнутой системы называется
проблемой замыкания и является центральной проблемой на пути построения моделей турбулентности, предназначенных для описания осредненных
полей скорости (температуры, концентрации примеси и т.д.).
Всеполуэмпирическиемодели основаны на различных искусственных
способах обрыва цепочки уравнений Ф ридмана-Келлера. Всякая процедура
замыкания тем или иным способом выражает моменты порядка n через
моменты низш их порядков с помощ ью неких гипотез. М оделями замыкания первого порядка называю т модели, выражаю щ ие моменты второго порядка через моменты первого порядка. М одели замыкания второго порядка оставляют моменты второго порядка, выражая через них моменты
третьего порядка и т.д. Н азвание полуэмпирические модели отражает тот
факт, что всемодели непременно содержат константы, требую щ ие их определения из опыта.
П роблему замыкания можно проиллю стрировать и на примере уравнения для давления. Как известно, уравнение для определения давления получается из уравнения Н авье-Стокса (3.12) путем применения к последнему
операции С . В результате получаетсяуравнение
Dp = - r (¶i v j ¶ j vi - ¶i f i ) .
(3.25)
В уравнение (3.25) подставляем разложения (3.14)
DP + Dp ў= - r¶ij (U iU j + U i v j + U j vi + vi v j ) - r (¶i Fi + ¶i f i )
2
(3.26)
102
и после осреднения получаем
DP = - r¶ij (U iU j + vi v j ) - r¶i Fi .
(3.27)
2
Таким образом, в уравнении для средних величин снова появился
тензор напряжений Рейнольдса. Для того чтобы выразить статистические
моменты, вклю чаю щ ие пульсации давления (см. уравнение (3.24)), потребуется написать уравнение для величины p ў, что можно сделать, вычтя (3.27)
из(3.26),
[
]
Dp ў= - r ¶ij (U i v j + U j vi + vi v j - vi v j ) - ¶i f i .
2
(3.28)
Это уравнение вклю чает и тензор напряжений Рейнольдса и произведение пульсаций, что неминуемо приведет при попытках написания уравнений для моментов, вклю чаю щ их пульсации давления, к появлению новых
моментов старш их порядков.
3.3 Турбулентная вязкость
Самыми простыми являются модели первого порядка, которые тем
или иным образом выражаю т тензор напряжений Рейнольдса через характеристики среднего поля скорости. П ри этом, практически все модели первого порядка оперирую т понятием «турбулентная вязкость». В наиболее
общ ем виде турбулентная вязкость вытекает изформулы Буссинеска, предложенной для тензора напряжений Рейнольдса по аналогии с выражением
для вязких напряжений, принятом для несжимаемой жидкости (1.10)
t ij =
ж¶U
¶U j
1 2
u i dij - n t з i +
з ¶x
¶xi
3
и j
ц
ч
ч
ш
(3.25)
Важно подчеркнуть, что в отличие от молекулярной вязкости, турбулентная вязкость n t не является свойством жидкости, а зависит от самого
течения и даже для заданного течения может меняться от точки к точке.
Другими словами, концепция турбулентной вязкости основана на рассмотрении некой «турбулентной жидкости», отличной по своим свойствам от
вязкой жидкости в турбулентном течении.
Самый простой подход к рассмотрению турбулентных течений состоит в том, чтобы предположить, что турбулентная вязкость и энергия
103
турбулентных пульсаций k = u i 2 / 2 для данного течения есть величины постоянные, не изменяю щ иеся от точки к точке. В этом случае уравнение
Рейнольдса (3.19) принимает простейш ий вид
¶t U i + U j ¶ jU i = - r - 1¶i P + (n + n t )¶2jjU i +
fi
.
(3.26)
Н есмотря на чрезвычайную грубость такого предположения, оно позволяет в некоторых случаях правдоподобно описывать крупномасш табную структуру турбулентного течения. П олученное реш ения представляет в
этом случае «ламинарный аналог» реального течения, так как получаемые
профили скорости соответствую т ламинарным, а не турбулентным режимам течения. Значения турбулентной вязкости часто превышают при этом
молекулярную вязкость на многие порядки. Так, например, для задач описания крупномасш табных течений в атмосфере принимают значения турбулентной вязкости в диапазоне 10 2 ё 10 4 м 2 / с, в то время как молекулярная
кинематическая вязкость воздуха равна 2 Ч10 - 5 м 2 / c (т.е. различие составляет
7-9 порядков !).
3.4 Длина пути смеш ения
М ногие простые схемы замыкания опираю тся на идею П рандтля о
длине пути смеш ения, характеристике потока, под которой понимают расстояние, проходимое жидкой частицей поперек потока, прежде чем происходит ее смеш ение с окружаю щ ей жидкостью . П онятие пути смеш ения исходит из аналогии между турбулентным перемеш иванием и молекулярным
переносом в газах, когда характеристики молекул остаются постоянными в
промежутках между соударениями.
М одель П рандтля применяется обычно к простым потокам, в которых средняя скорость имеет только одну компоненту (пограничные
слои,
r
каналы, трубы). Для определенности будем считать что U = (U x ,0,0) , а сущ ественным является только градиент средней скорости вдоль оси z . Тогда,
следуя П рандтлю (1925г.), можно написать, что
ж¶U ц
=-l з xч .
и ¶z ш
2
uxuz
2
(3.27)
Ф ормула (3.27) получается и из качественных соображений, использую щ их идею турбулентной вязкости. Действительно, если считать, что ве-
104
личина пульсаций скорости в турбулентном потоке пропорциональна градиенту средней скорости, то из размерных соображений появляется коэффициент с размерностью длины : u i » l ¶ zU x . Логично также предположить,
что турбулентная вязкость тем больш е, чем выше уровень турбулентных
пульсаций. Соображения размерности снова требую т наличия множителя с
размерностью длины: n t » lu . Тогда
n t » l 2 ¶ zU x ,
что в принципеэквивалентно формуле(3.27).
П еречислим некоторые задачи, в которых ш ироко используется гипотезаП рандтля о пути смеш ения.
Свободный слой со сдвигом ш ириной d . В этом случае длина пути
смеш ения считается постоянной
l = Cd ,
где С - эмпирическая константа, величина которой имеет порядок
С » 0,1 .
Турбулентный пограничный слой. П редположение о том, что размер
доминирую щ их вихрей пропорционален расстоянию от стенки z , приводит
к выражению
l = Cz .
В этом случае эмпирическая константа С » 0,4 .
Течение в открытом канале. Для ка??ала глубиной d используется
оценка
l = Cz 1 -
z
.
d
Эта формула применима и для закрытого канала. В этом случае глубина d заменяется на полуш ирину d / 2 . Ф ормула работает и в случае круглой трубы (вместо глубины в ней появляется радиус канала). Значение эмпирической константы в каждом случае свое.
Важно отметить, что определение длины пути смеш ения (длины перемеш ивания), предложенное П рандтлем (3.27) не является единственно
возможным. Ш ироко используются и некоторые другие модели, опираю щ иеся на это понятие. Н апример, Тейлор ввел модель, в которой тензор
напряжений Рейнольдса для одномерного турбулентного потока задается
выражением
105
u x u z = - lU x ¶ zU x .
(3.28)
3.5 М одели переноса турбулентной вязкости
В общ ем случае турбулентная вязкость меняется от точки к точке и
r
может изменяться со временем, то есть n t = n t (t , r ) . К моделям переноса турбулентной вязкости относятся модели, в которых для турбулентной вязкости записывается эволю ционноеуравнение.
Ф ормально, для любой переносимой течением скалярной величины
a , для которой выполняется закон сохранения, можно записать уравнение
вида
r
¶t a + (v С )a = ¶ j q j + G + D ,
(3.29)
r
где q - поток величины a засчет диффузии,
G - слагаемое, характеризую щ еегенерацию величины a ,
D - слагаемое, характеризую щ еедиссипацию этой величины.
Если предположить, что полная вязкость (сумма молекулярной и
турбулентной вязкостей) есть переносимая потоком скалярная величина, то
для нееможно записать уравнение вида (3.29).
П риведем в качестве примера такой модели переноса турбулентной
вязкости уравнение, предложенное Н и и Коважным для плоского пограничного слоя (Nee V., Kovasznay L. Simple phenomenological theory of turbulent shear flow, Phys.Fluids, 1969, V.12, P.473-484.)
¶tn t + U j ¶ jn t = ¶ j ((n + n t )¶ jn t )+ An t ¶ zU x - Bn t (n + n t )
(3.30)
Выражение для потока полной вязкости записано в предположении,
что коэффициент диффузии равен этой же полной вязкости (условие самодиффузии). Уравнение вклю чает две эмпирические константы. П араметр A
характеризует интенсивность генерации турбулентной вязкости за счет
сдвига (авторы модели принимали его значение близким к 0,1) и параметр
B , характеризую щ ий «самосжигание» турбулентной вязкости.
3.6 Двухпараметрическиемодели
106
Больш ую группу моделей составляю т модели, основанные на рассмотрении кинетической энергии пульсаций скорости k = u i 2 2 . В моделях
этого типа обычно появляется и вторая важная характеристика - скорость
диссипации энергии e . Турбулентная вязкость выражается через эти две
величины. Соображения размерности приводят к соотнош ению
nt = C
k2
e
Уравнение для энергии пульсаций скорости можно получить из уравнения (3.24), положив в нем j = i (не путаем в уравнении кинетическую
энергию пульсаций и индекс k ):
й жu 2 p ўц
щ
i
ч - n¶k k ъ + u i f iў. ,
¶t k + U k ¶k k = - u i u k ¶k U i - ¶k к u k з
з2
rч
к
ъ
ш
л и
ы
(3.31)
однако, это уравнение по-прежнему вклю чает неизвестныемоменты и
не снимает проблему замыкания.
Замыкание уравнения (3.31) приводит к ш ирокой группе моделей переноса кинетической энергии. Н е претендуя даже на беглый обзор полуэмпирических моделей этого типа, мы только приведем пример k - e модели
для описания течения в плоском пограничном слое на стенке
Рис.3.2.
¶t k + U x ¶ x k + U z ¶ z k = ¶ z (n t ¶ z k )+ n t (¶ zU x ) - e
2
(3.32)
107
e
e2
2
¶t e + U x ¶ x e + U z ¶ z e = ¶z (n t ¶ z e)+ C1 n t (¶ zU x ) - C 2
k
k
(3.33)
Замкнутую
систему
образую т
при
этом
уравнения
(3.19),(3.20),(3.25),(3.32) и (3.33).
Для иллю страции возможностей полуэмпирических моделей на рисунке3.2, взятом из книги [4], показаны результаты вычислений осесимметричного следа за ш аром в несжимаемой жидкости с помощ ью различных
моделей. Точками на рисунке обозначены экспериментальные данные,
пунктирной линией - результаты расчета с помощ ью однопараметрической
модели, ш трих-пунктирной - результаты расчета с помощ ью k - e модели,
сплош ной - результаты расчета с помощ ью другой двухпараметрической
модели, специально разработанной для свободных течений.
Очевидно, что уравнения для статистических моментов, характеризую щ их более сложное течение, например, турбулентную конвекцию,
должны вклю чать соответствую щ ие моменты для температурных пульсаций и смеш анные моменты, характеризую щ ие корреляции поля скорости и
поля температуры.
В заклю чение ещ е раз отметим, что полуэмпирические модели представляю т наиболее разработанное направление в изучении турбулентных
течений, и что по ним сущ ествует подробная литература. Для начального
систематического знакомства с ними можно порекомендовать удачно подобранные сборники статей под редакцией Ф роста и М оулдена [4] и Кольмана [5].
Рекомендуемая литература к третьей главе:
А.С.М онин, А.М .Яглом, Статистическая гидромеханика. Ч.1. М .:
Н аука, 1965. 639с.
А.С.М онин, А.М .Яглом, Статистическая гидромеханика. Ч.2. М .:
Н аука, 1967. 720с.
А.Дж.Рейнольдс, Турбулентные течения в инженерных приложениях.
М .: Энергия, 1979. 408с.
Турбулентность. П ринципы и применения. П од. ред. У.Ф роста,
Т.М оулдена. М .:М ир, 1980. 536с.
М етоды расчета турбулентных течений. П од. ред. В.Кольмана. М .:
М ир, 1984. 464с.
В этой
точке имеет место субкритическая
бифуркация Хопфа. Особенность
поведения системы Лоренца в том,
что устойчивый предельный цикл
не возникает в ней вовсе (напомним, что согласно сцерарию Рю эляТаккенса, странный аттрактор возникает после двух бифуркаций
Хопфа) и странный аттрактор возникает сразу после первой (обратной) бифуркации Хопфа. Бифуркационная диаграмма представлена
на рисунке 2.32. Следует отметить,
Рис. 2.35..
82
что «чистый» странный аттрактор сущ ествует в небольш ом интервале числа Релея 24,06 < r < 30,1 . Обратим внимание и на то, что на левом краю этого
интервала сущ ествует гистерезис - при понижении числа Релея странный
аттрактор сущ ествует до r = 24,06 , а не до r = 24,74 . Левее этой границы в
интервале чисел Релея r > 13,93 сущ ествует область так называемого метастабильного хаоса. В этой области малые возмущ ения стационарного реш ения монотонно затухают, но больш ие
возмущ ения приводят к хаотическим
режимам, которые в конечном итоге
также затухают, но успеваю т при этом
выписать в фазовом пространстве
многочисленные хаотические петли,
напоминаю щ ие поведение системы на
странном аттракторе. П ри r > 30,1 диаграмма режимов представляет собой
чередование областей с хаотическим и
периодическим движениями, напоминая поведение отображения Ф ейгенбаума в области m Ґ < m < 1 (рис.2.31).
Рис. 2.36.
П оявлению области с периодическим
аттрактором предш ествует обратный
каскад, а само «окно периодичности» включает субгармонический каскад.
Число «окон периодичности», по-видимому, бесконечно и при больш их
числах Релея их ш ирина растет. П оследнее окно неограниченно и занимает
всю область r > 214,364 .
В своей знаменитой работе Лоренц численно исследовал поведение
системы при r = 28 . Н а рисунке 2.33 показан фрагмент поведения во времени переменной X (t ) при этом значении r , а на рис.2.34 - характерный вид
фазовой траектории системы на странном аттракторе. Н а рис.2.35 - проекции фазовой траектории на плоскости ( X , Z ) .
Н аблю дение за эволю цией фазовой траектории показывает, что траектория описывает
витки вокруг точек, соответствую щ их ставш им
неустойчивыми реш ениям (2.44), переходя случайным образом от вращ ения вокруг одного
фокуса к вращ ению вокруг другого.
Н аблю дая за эволюцией фазовой траектории в плоскости ( X , Z ) , Лоренц сделал важный вывод. Траектория раскручивается воРис. 2.37.
круг одного фокуса, увеличивая на каждом
витке радиус орбиты. Этот процесс происхо-
83
дит до тех пор, пока на очередном витке в точке максимума траектория не
выйдет за значение Z = 38,5 . Как только траектория превысит это значение,
она уходит в область притяжения другого фокуса и все повторяется вновь.
П ри этом число витков, которое соверш ит траектория, зависит от величины превышения траектории над этим критическим значением перед перебросом. Лоренц использовал метод точечных отображений, позволяю щ ий
перейти от системы с непрерывным временем к системе с дискретным временем - вариант сечения П уанкаре, называемый отображением первого
возвращ ения. В качестве отображения использовалось значение величины
Z в текущ ем локальном максимуме, как функция от значения в предыдущ ем
максимуме (рис.2.36). Левая, восходящ ая часть функции соответствует
процессу раскручивания, а переход за пик - перебросу к другому фокусу.
Лоренц предложил простейш ую модель наблю даемого процесса - отображение отрезка [0,1] на себя вида (рис.37)
M n+ 1
м
2M
п
п n
=н
п2(1 - M )
n
п
о
1
2
1
Mn >
2
Mn <
(2.45)
Если рассматривается последовательность, начинаю щ аяся со значения M 0 , то она будетразвиватьсяпо следую щ ей цепочке:
м2M 0
M1 = н
о2 - 2M 0
м4M 0
п2 - 4M
п
0
M2 = н
п4 - 4M 0
п
о- 2 + 4M 0
м8M 0
п2 - 8M
0
п
п4 - 8M 0
п
п6 - 8M 0
M3 = н
п8 - 8M 0
п- 2 + 8M 0
п
п- 4 + 8M 0
п- 6 + 8M
0
о
....... M n = m n ± 2 n M 0 .
Здесь m n - четное число, такое, что оно сдвигает величину 2 n M 0 в интервал [0,1]. Все возможные последовательности можно разделить на три
типа:
П оследовательности, заканчиваю щ иеся в нуле. Таких последовательностей счетное множество и они начинаю тся с элемента видаM 0 = u 2 p , где
u - нечетное целое число. Тогда M p - 1 = 1 2 и M p = 0 .
84
П ериодические последовательности. Они возникают, если M 0 = u 2 p v ,
где u,v - простые числа. Тогда M p + 1+ k
2 p + 1+ k u
2 Ч2 k u
= m±
= m±
. П ростеш ие приv
2pv
меры получаю щ ихсяпериодических последовательностей есть
[2 3, ].....
[2 / 5, 4 / 5, ].....
[2 / 7, 4 / 7, 6 / 7, ]....
[2 / 9, 4 / 9, 8 / 9, ]....
3) Апериодические последовательности.
Эта модель иллю стрирует ещ е одно важное свойство системы - неустойчивость к малым возмущ ениям (ЧЗНУ). Действительно, если рассмотреть последовательность с малым возмущ ением начального элемента
M 0ў = M 0 + e , то после n итераций
M nў = m n ± 2 n ( M 0 ± e) = M n ± 2 n e ,
что свидетельствуетоб экспоненциальном ростевозмущ ений.
Отметим, что модельное отображение (2.45) при всей своей простоте
сохраняет важнейш ее свойство, приводящ ее к ЧЗНУ в диссипативных системах - это растяжение в сочетании со складыванием. Растяжение на каждом ш аге приводит к экспоненциальному росту начального смещ ения (расхождению траекторий), а складывание обеспечивает возвращ ение в ограниченную область (в данном случае интервал).
2.8.2 М одель динамо
Рикитаке
Другой пример динамической
системы со стохастическим поведением дает так называемая модель двухдискового динамо Рикитаке, предложенная в связи с задачей об инверсиях
геомагнитного поля. М агнитное поле
Земли в первом приближении представляет собой диполь, который по
палеомагнитным данным многократ-
Рис. 2.38.
85
но и нерегулярно менял свою полярность. Н а сегодняш ний день ш кала полярности геомагнитного поля восстановлена более чем за 1700 миллионов
лет, что составляет порядка половины возраста Земли. За это время зарегистрировано 593 переброса магнитного поля, причем время между двумя перебросами колеблется в интервале от 10 тысяч до сотен миллионов лет, демонстрируя хаотическое поведение, лиш енное каких-либо периодичностей.
Согласно принятой на сегодня точке зрения, магнитное поле Земли
возбуждается в результате конвективного движения в жидком (электропроводящ ем) ядре. П роцесс возбуждения магнитного поля в движущ ейся проводящ ей среде получил название М ГД-динамо. Земное динамо представляет собой сложный нелинейный магнитогидродинамический процесс, исследование которого находится лиш ь на начальной стадии. Больш ой интерес
представляют поэтому любые упрощ енный модели процесса генерации
магнитного поля, способные приводить к случайным сменам полярности
генерируемого магнитного поля.
Самые простые модели оперируют не потоками проводящ ей жидкости, а движущ имися проводниками. П ервая попытка построить такого рода модель принадлежит Булларду (Bullard E.C., Proc.Cambridge Philos.
Soc.,1955, v.51, p.744.), который предложил однодисковое динамо, но такая
модель не дает смены полярности генерируемого поля. Рикитаке (Rikitake
T., Proc.Cambridge Philos. Soc.,1958, v.54, p.89.) рассмотрел систему двух
дисковых динамо, связанных таким образом, что ток от одного диска питает катуш ку возбуждения другого и наоборот. Эта ситуация изображена на
рис.2.38. Оба диска вращ аю тся без трения и находятся под действием одинаковых моментов сил G , компенсирую щ их омические потери в дисках и
обмотках. Уравнения, описываю щ ие эволю цию токов I 1 , I 2 и угловых скоростей W 1 , W 2 можно записать в виде
LI&1 + RI 1 = MW 1 I 2 ,
LI&2 + RI 2 = MW 2 I 1 ,
CW& = G - MI I ,
1
(2.46)
1 2
CW& 2 = G - MI 1 I 2 ,
где L - коэффициент самоиндукции, R - сопротивление каждой цепи, M коэффициент взаимоиндукции, C - момент инерции диска.
Два последних уравнения (2.46) показываю т, что разность угловых
скоростей есть величина постоянная
W1 - W 2 =
GL
A,
CM
где A - константа. Это позволяет перейти к системетрех уравнений.
86
Система записывается в безразмерном виде. П ри этом за единицу тока принимаю т величину G / M , угловой скорости - GL / CM , а за единицу
времени - величину t et m . Единица времени выражена через два характерных масш таба времени, присущ их системе. Это время t m , за которое диск
под действием приложенного момента сил разгоняется до характерной
скорости R / M ,
tm =
CR
GM
и время электромагнитной диффузии
te =
L
,
R
характеризую щ ее время вырождения магнитного поля при остановке диска. И х отнош ение является безразмерным параметром системы
m=
tm
CR 2
=
.
t e GLM
Обозначая безразмерные токи как X i , а безразмерные угловые скорости
как Yi (в уравнениях остается одна переменная Y , так как Y1 - Y2 = A ), прихо-
Рис. 2.39.
87
дим к системе
X&1 + mX 1 = YX 2 ,
X&2 + mX 2 = (Y - A) X 1 ,
Y& = 1 - X X .
1
(2.47)
2
Система (2.47) имеетстационарные реш ения
X 1 = ±K ,
X 2 = ± K - 1 , Y = Y1 = mK 2 , Y2 = mK - 2 ,
где A = m ( K 2 - K - 2 ).
М ы не будем подробно описывать свойства системы Рикитаке, оставляя ее изучение для самостоятельных работ. Н а рис.2.39 показана только
фазовая траектория системы для случая m = 1,5; K = 2. М ожно видеть, что ее
топология близка аттрактору Лоренца.
2.8.3 Реальная конвекция
Н аибольш ее число экспериментальных работ по исследованию перехода от упорядоченных течений к хаотическим выполнено, пожалуй, в исследованиях конвективных течений. М ы приведем некоторые результаты
исследований перехода от ламинарного движения к турбулентности при
конвекции в кубической полости, взятые из работы: Зимин В.Д., Кетов
А.И . Н адкритические конвективные движения в кубической полости.
И зв.АН СССР, М еханика жидкости и газа, 1974, N.5, С.110.
Рис. 2.40.
88
Рис. 2.41.
Рис. 2.42.
И змерения проводились подогреваемой снизу в кубической полости с
ребром 40 мм, образованной медными стенками. Горизонтальные стенки
термостатировались, обеспечивая заданную разность температуры, а вертикальные обеспечивали равновесный однородный градиент температуры.
Н адкритические течения, возникаю щ ие в кубической полости и имею щ ие
наиболее низкие уровни устойчивости, схематически показаны на рис.2.40,
где стрелками показано направление движения жидкости в верхней части
полости, а знаками «плю с» и «минус» обозначены области, в которых температура
оказывается выше или ниже средней. Критические числа Релея для движений типа А
и Б равны 8224, для В - 9184 и для Г 14032. В полости были установлены дифференциальные термопары, расположенные таким образом, что их показания позволяли выделять движения всех четырех
типов.
Н е останавливаясь на сценариях
развития неустойчивости и переходов от
одного режима движения к другому, приведем лиш ь некоторые данные, иллю стрирую щ ие поведение системы в одночастотном режиме, двухчастотном и стохастическом режимах. Для каждого из трех режимов на рисунках представлены изменения
во времени показаний термопар, соответствую щ их каждому из выделяемых течеРис. 2.43.
ний, проекции фазовых траекторий на
89
Рис. 2.44.
Рис. 2.45.
плоскости, образованные всеми парами термопар и спектры мощ ности
пульсаций температуры, регистрируемой каждой из четырех термопар.
Рисунки 2.41-2.43 относятся к одночастотному режиму, регистрируемому при числе Релея R = 2 Ч10 5 . П ервый рисунок показывает характер колебаний показаний всех четырех термопар, второй - соответствую щ ие этим
колебаниям проекции фазовых траекторий, ясно указываю щ ие на сущ ествование предельного цикла. Об этом
же свидетельствую т и спектры Ф урье
(рис.2.43) состоящ их их одного главного пика на частоте 0,054 Гц и пика на
удвоенной частоте, обусловленный негармонической формой колебаний.
Следую щ ая группа рисунков
представляет результаты для числа Релея R = 2,24 Ч10 5 . Н а рисунке 2.44 показаны пульсации показаний термопар,
на рис.2.45 - соответствую щ ие фазовые
траектории (за время соответствую щ ее
периоду низкочастотных колебаний), а
на рисунке 2.46 - спектры, свидетельствую щ ие о сущ ествовании двухчастотного режима (частоты 0,0451 Гц и 0,304
Гц).
Движение становится стохастическим при R = 2,50 Ч10 5 . П оказания
термопар для этого режима представРис. 2.46.
лены на рис.2.47, фазовые траектории -
90
на рис.2.48, а спектры мощ ности - на рис.2.49. Видно, что фазовые траектории имеют чрезвычайно запутанную структуру, а спектры становятся
сплош ными, сохраняя лиш ь слабые локальные максимумы, свидетельствую щ ие о сохранении периодических составляю щ их.
Рис. 2.47.
Рис. 2.48.
Рис. 2.49.
91
Рекомендуемая литература ко второй главе:
1. БержеП ., П омо И ., Видаль Л. П орядок в хаосе. М осква: М ир. 1991.
366с.
2. Ш устер Г. Детерминированный хаос. М осква: М ир. 1988. 240с.
3. Странные аттракторы. Сборник статей. Серия «М атематика. Н овоев зарубежной науке», выпуск 22. М осква: М ир. 1981. 254с.
92
3 ПОЛУЭМ П И РИ Ч ЕСК И Е М ОДЕЛИ
3.1 Развитая турбулентность
3.1.1 Вводныезамечания
В данной главе мы начинаем рассматривать подходы к описанию
развитой турбулентности, то есть течений, возникаю щ их при значительном
превышении критических значений управляю щ их параметров (числа Рейнольдса, если речь идет об изотермическом течении в отсутствии дополнительных силовых полей). Такие течения характеризую тся наполненными
спектрами Ф урье, причем нетолько временными, но и пространственными.
Н апомним ещ е раз, что именно в этом и есть основное отличие турбулентности от хаоса в динамических системах невысокого порядка: в турбулентном потоке хаос и пространственный, и временной, а хаотическое поведение маломодовых систем (соответствую щ их например конвективным течениям при невысокой надкритичности) представляет собой хаотическую во
времени эволю цию мод с относительно простой пространственной структурой.
П риступая к рассмотрению развитых турбулентных течений, следует
сделать ряд важных замечаний. П ервое из них касается уравнений движения жидкости. В первой главе мы получили уравнения Н авье-Стокса, как
основные уравнения, с помощ ью которых мы описываем в дальнейш ем все
течения жидкости. Снова подчеркнем, что мы действительно продолжаем
считать, что эти уравнения описываю т течения жидкости и в турбулентном
режиме, даже при экстремально больш их значениях безразмерных параметров (более того, мы будем рассматривать только случай несжимаемой
жидкости). Уверенность в том, что это возможно, держится на результатах
многочисленных успеш ных попыток использования этих уравнений для
турбулентных течений. Сама возможность приложения уравнений Н авьеСтокса к турбулентности совсем не очевидна (и продолжает подвергаться
критике), так как при их выводе было сделано достаточно сильное предположение о том, что тензор вязких напряжений включает в себя только линейные комбинации первых производных поля скорости. В ламинарных и
слабо надкритических течениях это предположение кажется разумным и
прекрасно работает, но в сильно нелинейных режимах нельзя исклю чить,
что тензор вязких напряжений будет иметь более сложную зависимость от
структуры поля скорости. Оправданием использованию уравнений движения в принятой формеможетслужить только сопоставление результатов их
реш ения с экспериментальными данными.
93
Далее, пусть уравнения движения справедливы и предположим, что
мы располагаем мощ нейш им компью тером, способным реш ать трехмерные
уравнения движения с лю бой желаемой точностью (например, будем считать трехмерный поток на сетке 1000х1000х1000). Это, однако, не снимает
проблемы описания турбулентности, так как в результате такого реш ения
мы будем иметь огромное количество информации, осознание которой
требует ее представления в некотором виде, а это фактически опять же
предполагает введение определенной модели процесса. П о сути, такой суперкомпью тер отличается от реального турбулентного течения, наблю даемого в лаборатории или природе, только несравненно больш ими возможностями съема информации относительно состояния потока в лю бой точке
и в любой момент времени.
П роблема описания турбулентного движения состоит в выделении
характеристик, описываю щ их свойства системы с огромным числом степеней свободы, а лю бой подход к ее описанию - это тот или иной способ ограничения числа степеней свободы.
Турбулентные поля (скорость, давление, температура и т.д.) представляю т собой случайные поля. В лю бой точке потока можно установить
датчик и зарегистрировать реализацию процесса в данной точке. М ногократно повторяя эту процедуру, принципиально возможно получить плотr
ность вероятности P ( f ) для интересую щ ей нас величины f (r , t ) . В общ ем
случае, плотность вероятности также есть функция координат и времени.
Сущ ествуетряд важных частных случаев, которыемы и перечислим.
Турбулентность является однородной, если плотность распределения
вероятности независит от сдвига
r
r
r
P (t , r + Dr ) = P(t , r ) .
Турбулентное течение называется стационарным, если плотность вероятности независит от времени, то есть
r
r
P (t + t , r ) = P (t , r ) .
П роцесс называется эргодическим, если осреднение по времени эквивалентно для него осреднению по ансамблю реализаций
T
r
r
1
f (r ) = lim тf (t , r )dt .
T® Ґ T 0
Угловыми скобками будем обозначать среднее по ансамблю реализаций. Очевидно, что только стационарный процесс может быть эргодическим. Гипотеза эргодичности ш ироко используется при исследовании ста-
94
ционарных течений, так как на практике измеряю тся именно средние по
времени величины.
В реальных измерениях ш ироко используется и гипотеза Тейлора, позволяю щ ая связать пространственные и временные флуктуации исследуеr
мой величины f (r , t ) . Согласно этой гипотезе, если сущ ествует среднее теr
чение, характеризуемой скоростью U , то справедливо соотнош ение
¶f
¶f
= Ui
.
¶t
¶xi
П ользуясь этой гипотезой, по измерениям в заданной точке пространства определяют пространственные флуктуации исследуемого поля и
их статистические характеристики.
3.1.2 Статистическиемоменты случайных полей
r
Ф ункция распределения плотности вероятности P (r , t ) содержит полr
ную информацию о случайном поле f (r , t ) , однако, ееопределение в полном объеме практически невозможно. И звестно, что заданию плотности
вероятности эквивалентно задание последовательности (в принципе- бесконечной) статистических моментов
Mf m = тf m P ( f )df .
П ри этом момент нулевого порядка равен единице в силу условия нормировки
Mf 0 = тP ( f )df = 1 ,
а момент первого порядка, называемый такжематематическим ожиданием,
даетсреднее значение величины
Mf 1 = тfP ( f )df = f .
Для моментов второго и болеевысоких порядков обычно использую т
центральныемоменты, вычисляемые относительно средних значений
M(f -
f ) m = т( f -
f ) m P ( f )df .
Н апомним, что центральный момент второго порядка называется
дисперсией.
95
С точки зрения описания турбулентных полей, необходимы статистиr
ческие характеристики связи между значениями величины f (r , t ) в различных точках пространства. Это требует введения совместной плотности ве-
Рис.3.1.
r
r
роятности P ( f (r1 ), f (r2 )) и (или) соответствую щ их двухточечных моментов.
Важнейш им среди двухточечных моментов является момент второго порядка, называемый корреляционной функцией
r r
r
r
r
r
B (r1 , r2 ) = т( f (r1 ) - f1 )( f (r2 ) - f 2 ) P ( f (r1 ), f (r2 ))df df 2 = ( f1 - f1 )( f 2 - f 2 ) . (3.1)
1
Если речь идет о векторном поле(например, скорости), то появляется
корреляционный тензор
r r
r
r
r
r
Bij (r1 , r2 ) = (vi (r1 ) - vi (r1 )(v j (r2 ) - v j (r j ) ) .
(3.2)
Для однородной турбулентности (3.1) и (3.2) зависят только от взаr r r
имного расположения двух точек, то есть, если r2 = r1 + r , то
r r
r
Bij (r1 , r2 ) = Bij (r ) .
(3.3)
Важным частным случаем является однородная и изотропная турбулентность, в которой совместная плотность вероятности (а, следовательно,
r
и двухточечныемоменты) не зависят и от направления вектора r . Тогда
r r
r
Bij (r1 , r2 ) = Bij (| r |) = Bij (r ) .
(3.4)
96
Чащ е всего использую т корреляционные функции Bll (r ) и Bnn (r ) , характеризую щ ие корреляцию продольных и поперечных составляю щ их
пульсаций скорости. Здесь индексом l обозначена составляю щ ая скорости
r
r
вдоль линии, соединяю щ ей точки r1 и r2 , а индексом n составляю щ ая, нормальная этой линии. Характерный вид этих функций иллю стрирует рисунок 3.1.
Выш е, в параграфе 2.4.3, указывалось на связь корреляционной
функции со спектрами (теорема Хинчина) в случае временного сигнала.
Аналогичное соотнош ение связывает и пространственные спектры с двухточечными корреляционными функциями. П режде чем написать это соотнош ение, остановимся несколько подробнее на вопросе о пространственных спектрах турбулентности.
3.1.3 П ространственныеспектры
П редположим, что рассматриваемое случайное (турбулентное) поле
занимает ограниченный объем и величина f ( x, y, z, t ) может быть представлена интегралом Ф урье
Ґ
rr r
) r
r
1
f ( r ,t ) = 3 тf ( k ,t )e ik r dk ,
8p - Ґ
(3.5)
где
) r
f ( k ,t ) =
Ґ
rr
r
- ik r r
f
(
r
,
t
)
e
dr ,
т
(3.6)
-Ґ
r
r
r = ( x, y, z ) - радиус-вектор, k = (k x , k y , k z ) - волновой вектор.
Считая рассматриваемую турбулентность стационарной, определим
трехмерный энергетический спектр случайного поля :
r
) r
F ( k ) =< | f ( k ) |2 >
(3.7)
Угловыескобки означают в этом случае осреднение по времени. Трехмерr
ный спектр связан с корреляционной функцией B (r ) (теорема Хинчина)
r
r rr r
1
F (k ) = 3 тB (r )e - ik r dr
8p
(3.8)
В теории турбулентности, говоря о ее спектральных свойствах, обычно
имею т в виду энергетический спектр E (k ) , который характеризует энергию
97
всех гармоник с заданным модулем волнового вектора, независимо от его
направления.
r r
E (k ) = тF (k )dk ,
(3.9)
r
|k |
или, в сферической системе координат,
2p p
E (k ) =
r
ттF (k )k
2
sin J dJ dj .
0 0
r
В важном частном случае изотропной турбулентности, когда F (k ) = F (k ) ,
связь становитсяочень простой:
E (k ) = 4pk 2 F (k ) .
(3.10)
Отметим, что все оценки для спектральных законов развитой турбулентности касаются обычно именно энергетического спектра E (k ) .
Если в турбулентном потоке измерения проводятся вдоль одной прямой, то по этим измерениям можно построить одномерное фурьепреобразование. Ограничиваясь однородной и изотропной турбулентностью , в которой все прямые равноправны, рассмотрим прямую y = z = 0 и
запиш ем
+Ґ
)
f1 ( k x ) = тf ( x , y , z )e - ixk x dx .
-Ґ
Квадрат модуля этой величины есть одномерный энергетический спектр
)
F1 ( k x ) =| f1 ( k x ) |2 .
(3.11)
Чтобы получить связь между одномерным и трехмерным спектрами,
выразим исходную величину на прямой y = z = 0 через обратное преобразованиеФ урье. С одной стороны
f ( x ,0,0 ) =
1
2p
)
тf ( k
1
x
)e ixk x dk x ,
а с другой стороны
1 )
i ( xk + 0 k + 0 k )
f ( k x ,k y ,k z )e x y z dk x dk y dk z =
3 т
8p
) r
1
ж 1
ц
=
f
( k )dk y dk z ч e ixk x dk x .
з
2 т
т
2p
и 4p
ш
f ( x ,0,0 ) =
98
Таким образом,
)
1
f1 ( k x ) =
4p 2
а
F1 (k x ) =
) r
f
т ( k )dk y dk z ,
r
1
F
(
k
)dk y dk z .
16p 4 т
В следую щ их главах, рассматривая структуру мелкомасш табной турбулентности, мы постоянно будем обращ аться к спектрам, описываемым
степенными законами. П окажем, как связаны между собой введенные спектры турбулентности при степенной зависимости энергии от масш таба (волнового числа). П усть имеется однородное изотропное поле скалярной величины, энергетический спектр которой следует степенному закону
E (k ) ~ k a .
Тогда трехмерный спектр
F (k ) ~ k
a- 2
2
= (k x + k y + k z )
2
2
a- 2
2
,
а одномерный
F1 (k x ) ~
тт(k
ж
a- 2
= k x ттз1 +
з
и
= kx
a
тт(1 + h
2
x
жk y
з
зk
и x
2
+ ky + kz )
2
2
2
2
ц жk z ц
ч
з ч
ч
ч + з
ш иk x ш
+ x2
)
a- 2
2
a- 2
2
dk y dk z =
a- 2
2
ц
ч
ч
ш
dhdx ~ k x
dk y dk z =
a
(проведена замена переменных h = k y / k x ; x = k z / k x ).
Таким образом, в однородной изотропной турбулентности энергетический спектр E (k ) и одномерный спектр F1 (k ) следую т одному степенному
закону, а степень убывания трехмерного спектра меньш е на двойку (т.е.
трехмерный спектр значительно круче).
99
3.2 Уравнения для статистических моментов
3.2.1 Уравнение Рейнольдса
Рассмотрим уравнения Н авье-Стокса в тензорных обозначениях
¶t vi + v j ¶ j vi = - r - 1¶i p + n¶2jj vi + f i ,
(3.12)
(3.13)
¶k v k = 0 .
Входящ ие в них величины представим в виде сумм средних полей и
пульсаций:
r
r
r
v i (r , t ) = U i ( r , t ) + u i ( r , t ) ,
r
r
r
p(r , t ) = P (r , t ) + p ў(r , t )
(3.14)
П ри этом, согласно принятым определениям, предполагаются следую щ ие правила осреднения (угловые скобки по-прежнему обозначают осреднение по ансамблю реализаций):
vi = U i ,
p = P,
Ui = Ui ,
P = P,
(3.15)
(3.16)
u i = 0;
p ў = 0;
Разложения (3.14) подставим в исходные уравнения (3.12)-(3.13):
¶tU i + ¶t u i + U j ¶ jU i + U j ¶ j u i + u j ¶ jU i + u j ¶ j u i = - r - 1 (¶i P + ¶i p ў) + n (¶2jjU i + ¶2jj u i ) + f i
(3.17)
¶k U k + ¶k u k = 0 ,
(3.18)
и проведем осреднение
¶t U i + ¶t u i + U j ¶ j U i + U j ¶ j u i + u j ¶ j U i + u j ¶ j u i =
- r - 1 (¶i P + ¶i p ў) + n (¶ 2jj U i + ¶2jj u i ) +
fi
¶k U k + ¶k u k = 0 .
Учитывая правила осреднения (3.15)-(3.16), приходим к уравнению
Рейнольдса:
¶tU i + U j ¶ jU i = - r - 1¶i P + n¶2jjU i - ¶ j u j u i +
fi
,
(3.19)
и уравнению неразрывности длясреднего поля скорости
100
¶k U k = 0 .
(3.20)
В уравнение Рейнольдса для средних полей входит одноточечный
корреляционный тензор пульсаций скорости, называемый тензором напряжений Рейнольдса
t ij = u i u j .
(3.21)
Этот тензор нельзя выразить через осредненные характеристики турбулентных полей. Следовательно, число неизвестных превышает число
имею щ ихся уравнений и система (3.19)-(3.20) являетсяне замкнутой.
3.2.2 Ц епочка уравнений Ф ридмана-Келлера и проблема
замыкания
В уравнении Рейнольдса появилась новая неизвестная величина - тензор напряжений Рейнольдса (3.21), для которого также можно получить
эволю ционноеуравнение. Так как
¶tt ij = ¶t u i u j = u i ¶t u j + u j ¶t u i ,
то сначала требуется получить уравнение для пульсаций скорости,
для чего из уравнения (3.17) необходимо вычесть уравнение (3.19). П олучим
(немые индексы j заменены на k )
¶t u i + U k ¶k u i + u k ¶k U i + u k ¶k u i = - r - 1¶i p ў- ¶k u i u k + n¶2kk u i + f iў.
(3.22)
Аналогичное уравнение получается и для компоненты u j :
¶t u j + U k ¶k u j + u k ¶k U j + u k ¶k u j = - r - 1¶ j p ў- ¶k u j u k + n¶2kk u j + f jў.
(3.23)
Уравнение (3.22) умножается на u j и складывается с уравнением
(3.23), умноженным на u i :
u i ¶t u j + u j ¶t u i =
- U k ¶k (u i u j ) - u j u k ¶k U i - u i u k ¶k U j - u j ¶k (u i u k ) - u i ¶k (u j u k ) - u i ¶k u j u k - u j ¶k u i u k
2
- r - 1 (u i ¶ j p ў+ u i ¶ j p ў) - n (u i ¶kk
u j + u j ¶2kk u i ) + u i f jў+ u j f iў
101
П осле осреднения приходим к уравнению:
(
)
¶t u j u j + U k ¶ k u i u j = - u i u k ¶ k U j + u j u k ¶ k U i - ¶ k u i u j u k
- r - 1 ( u i ¶ j p ў + u j ¶i p ў ) - n ( u i ¶2kk u j + u j ¶ 2kk u i ) + u i f jў + u j f ў .
(3.24)
i
В уравнении для корреляционного тензора пульсаций скорости второго порядка (3.24) появился корреляционный тензор (момент) третьего
порядка u i u j u k и новые моменты второго порядка, описываю щ ие корреляции пульсаций компонент скорости с давлением и скорости со вторыми
производными скорости.
Для вновь появивш ихся статистических моментов также можно написать эволю ционные уравнения типа (3.24), но проблемы это не реш ит,
так как в уравнение для момента третьего порядка войдут момент четвертого порядка и новые моменты третьего порядка и так далее. Система
уравнений для моментов все возрастаю щ их порядков называется цепочкой
уравнений Ф ридмана-Келлера и является незамкнутой в принципе. П роблема обрыва этой цепочки и получения замкнутой системы называется
проблемой замыкания и является центральной проблемой на пути построения моделей турбулентности, предназначенных для описания осредненных
полей скорости (температуры, концентрации примеси и т.д.).
Всеполуэмпирическиемодели основаны на различных искусственных
способах обрыва цепочки уравнений Ф ридмана-Келлера. Всякая процедура
замыкания тем или иным способом выражает моменты порядка n через
моменты низш их порядков с помощ ью неких гипотез. М оделями замыкания первого порядка называю т модели, выражаю щ ие моменты второго порядка через моменты первого порядка. М одели замыкания второго порядка оставляют моменты второго порядка, выражая через них моменты
третьего порядка и т.д. Н азвание полуэмпирические модели отражает тот
факт, что всемодели непременно содержат константы, требую щ ие их определения из опыта.
П роблему замыкания можно проиллю стрировать и на примере уравнения для давления. Как известно, уравнение для определения давления получается из уравнения Н авье-Стокса (3.12) путем применения к последнему
операции С . В результате получаетсяуравнение
Dp = - r (¶i v j ¶ j vi - ¶i f i ) .
(3.25)
В уравнение (3.25) подставляем разложения (3.14)
DP + Dp ў= - r¶ij (U iU j + U i v j + U j vi + vi v j ) - r (¶i Fi + ¶i f i )
2
(3.26)
102
и после осреднения получаем
DP = - r¶ij (U iU j + vi v j ) - r¶i Fi .
(3.27)
2
Таким образом, в уравнении для средних величин снова появился
тензор напряжений Рейнольдса. Для того чтобы выразить статистические
моменты, вклю чаю щ ие пульсации давления (см. уравнение (3.24)), потребуется написать уравнение для величины p ў, что можно сделать, вычтя (3.27)
из(3.26),
[
]
Dp ў= - r ¶ij (U i v j + U j vi + vi v j - vi v j ) - ¶i f i .
2
(3.28)
Это уравнение вклю чает и тензор напряжений Рейнольдса и произведение пульсаций, что неминуемо приведет при попытках написания уравнений для моментов, вклю чаю щ их пульсации давления, к появлению новых
моментов старш их порядков.
3.3 Турбулентная вязкость
Самыми простыми являются модели первого порядка, которые тем
или иным образом выражаю т тензор напряжений Рейнольдса через характеристики среднего поля скорости. П ри этом, практически все модели первого порядка оперирую т понятием «турбулентная вязкость». В наиболее
общ ем виде турбулентная вязкость вытекает изформулы Буссинеска, предложенной для тензора напряжений Рейнольдса по аналогии с выражением
для вязких напряжений, принятом для несжимаемой жидкости (1.10)
t ij =
ж¶U
¶U j
1 2
u i dij - n t з i +
з ¶x
¶xi
3
и j
ц
ч
ч
ш
(3.25)
Важно подчеркнуть, что в отличие от молекулярной вязкости, турбулентная вязкость n t не является свойством жидкости, а зависит от самого
течения и даже для заданного течения может меняться от точки к точке.
Другими словами, концепция турбулентной вязкости основана на рассмотрении некой «турбулентной жидкости», отличной по своим свойствам от
вязкой жидкости в турбулентном течении.
Самый простой подход к рассмотрению турбулентных течений состоит в том, чтобы предположить, что турбулентная вязкость и энергия
103
турбулентных пульсаций k = u i 2 / 2 для данного течения есть величины постоянные, не изменяю щ иеся от точки к точке. В этом случае уравнение
Рейнольдса (3.19) принимает простейш ий вид
¶t U i + U j ¶ jU i = - r - 1¶i P + (n + n t )¶2jjU i +
fi
.
(3.26)
Н есмотря на чрезвычайную грубость такого предположения, оно позволяет в некоторых случаях правдоподобно описывать крупномасш табную структуру турбулентного течения. П олученное реш ения представляет в
этом случае «ламинарный аналог» реального течения, так как получаемые
профили скорости соответствую т ламинарным, а не турбулентным режимам течения. Значения турбулентной вязкости часто превышают при этом
молекулярную вязкость на многие порядки. Так, например, для задач описания крупномасш табных течений в атмосфере принимают значения турбулентной вязкости в диапазоне 10 2 ё 10 4 м 2 / с, в то время как молекулярная
кинематическая вязкость воздуха равна 2 Ч10 - 5 м 2 / c (т.е. различие составляет
7-9 порядков !).
3.4 Длина пути смеш ения
М ногие простые схемы замыкания опираю тся на идею П рандтля о
длине пути смеш ения, характеристике потока, под которой понимают расстояние, проходимое жидкой частицей поперек потока, прежде чем происходит ее смеш ение с окружаю щ ей жидкостью . П онятие пути смеш ения исходит из аналогии между турбулентным перемеш иванием и молекулярным
переносом в газах, когда характеристики молекул остаются постоянными в
промежутках между соударениями.
М одель П рандтля применяется обычно к простым потокам, в которых средняя скорость имеет только одну компоненту (пограничные
слои,
r
каналы, трубы). Для определенности будем считать что U = (U x ,0,0) , а сущ ественным является только градиент средней скорости вдоль оси z . Тогда,
следуя П рандтлю (1925г.), можно написать, что
ж¶U ц
=-l з xч .
и ¶z ш
2
uxuz
2
(3.27)
Ф ормула (3.27) получается и из качественных соображений, использую щ их идею турбулентной вязкости. Действительно, если считать, что ве-
104
личина пульсаций скорости в турбулентном потоке пропорциональна градиенту средней скорости, то из размерных соображений появляется коэффициент с размерностью длины : u i » l ¶ zU x . Логично также предположить,
что турбулентная вязкость тем больш е, чем выше уровень турбулентных
пульсаций. Соображения размерности снова требую т наличия множителя с
размерностью длины: n t » lu . Тогда
n t » l 2 ¶ zU x ,
что в принципеэквивалентно формуле(3.27).
П еречислим некоторые задачи, в которых ш ироко используется гипотезаП рандтля о пути смеш ения.
Свободный слой со сдвигом ш ириной d . В этом случае длина пути
смеш ения считается постоянной
l = Cd ,
где С - эмпирическая константа, величина которой имеет порядок
С » 0,1 .
Турбулентный пограничный слой. П редположение о том, что размер
доминирую щ их вихрей пропорционален расстоянию от стенки z , приводит
к выражению
l = Cz .
В этом случае эмпирическая константа С » 0,4 .
Течение в открытом канале. Для ка?
Документ
Категория
Без категории
Просмотров
52
Размер файла
1 155 Кб
Теги
1998, подход, pdf, часть, турбулентность, модель, фрик
1/--страниц
Пожаловаться на содержимое документа