CC BY
26
ЭКСПЕРИМЕНТАЛЬНАЯ РЕАЛИЗАЦИЯ МОДЕЛИ ЛОРЕНЦА КОНВЕКТИВНОЙ НЕУСТОЙЧИВОСТИ ЖИДКОСТИ В ВЕРТИКАЛЬНОЙ ТОРОИДАЛЬНОЙ ЯЧЕЙКЕ
Р.А. Браже, О.Н. Куделин
Экспериментально исследованы режимы устойчивой и неустойчивой конвекции глицерина в вертикальной тороидальной ячейке. Представлены результаты фурье-анализа, БРА, вейвлет- и корреляционного анализа особенностей движения жидкости в тороиде. Построен хаотический аттрактор с признаками аттрактора Лоренца.
1. Требования к эксперименту
С тех пор как Э. Лоренц (1963, [1]) при построении математической модели турбулентности в атмосфере получил систему из трех сравнительно простых обыкновенных дифференциальных уравнений, допускающих детерминированные непериодические решения (детерминированный хаос), не прекращаются попытки экспериментальной реализации такой модели. Вот эти уравнения:
X = о (У - X),
у = -хг + гх - у, (1)
г = ху - ьг.
Здесь о = v/х - число Прандтля, V - кинематическая вязкость, х - температуропроводность; г = И,а/Кас, И,а - число Рэлея, И,ас - критическое число Рэлея; ь - параметр, зависящий от формы конвективной ячейки (в случае тороидальной ячейки [2] Ь = 1); X - скорость движения жидкости; У - температура в средней по высоте части тора; г - температура в нижней (нагреваемой) части тора.
Устойчивое конвективное течение (по часовой стрелке или против часовой стрелки) в ячейке наблюдается в диапазоне следующих значений г [1]:
1 <г <г*, (2)
где
= о (о + Ь + 3) = о (о + 4) ' о - Ь -1 о - 2 • (3)
При г > г* конвекция становится неустойчивой.
Результаты многочисленных ранних попыток экспериментального получения неустойчивых, в особенности стохастических, режимов конвекции жидкости в подогреваемых снизу вертикальных ячейках [2-6] не соответствовали модели Лоренца. Как показано в [7], это связано с нарушением предположения об одномерном квази-пуазейлевском характере течения жидкости в канале и изменением температуры по сечению канала. Таким образом, в этих экспериментах не реализовался случай трех степеней свободы, заложенный в основу модели Лоренца (1).
С.М. Дроздов [7] вывел два критерия, обеспечивающих физическую реализацию модели Лоренца для неустойчивых течений, а именно, ограничение по величине:
1) числа Тейлора
/ \ 1/2
Те = Ие( < 84, (4)
где Ие - число Рейнольдса, т0 и Кс - соответственно радиус сечения канала и радиус самого тора;
2) числа Рэлея
Иа < 544. (5)
Чтобы удовлетворить этим условиям, ему пришлось изготовить сложную конструкцию [7] в виде многослойного тора с то = 2.0 мм и Кс = 100 мм, заполненного ртутью. Рабочий канал был окружен последовательными слоями пенопласта и алюминия, и вся система помещалась в воду, служащую для выравнивания температуры на внешней поверхности тора. Нагрев ртути в нижней части канала осуществлялся пропусканием через нее электрического тока.
На наш взгляд, предложенная конструкция вызывает серьезные замечания. Во-первых, протекание электрического тока в проводящей среде провоцирует ее дополнительное течение по каналу и, как минимум, мешает конвекции. Во-вторых, ртуть обладает очень большой температуропроводностью, так что для нее число Прандт-ля о & 0.008 и, согласно (3), т* = —0.016. Таким образом, для ртути выполнение условия (2) для исследования устойчивой конвекции невозможно.
В нашей работе тороидальный канал радиусом Кс = 30 мм и радиусом сечения то = 2.0 мм был заполнен глицерином. Выбор жидкости обусловлен тем, что для наблюдения конвекции на границе порядка и хаоса, когда
адЕ3ЛТ *
т =-> т ,
радиус тороидального канала
/ 8я%у т* \1/3
ЛТ,
Здесь а - температурный коэффициент объемного расширения, д - ускорение свободного падения, лТ - разность температур в подогреваемой и охлаждаемой частях тора. Для глицерина а = 5.0 • 10_4К_1; х = 9.36 • 10"8м2/с; V = 1.10 • 10"3 м2/с; о = 1.18 • 104 и, согласно (3), т* & о = 1.18 • 104. Тогда из (6) для ЛТ = 30 К получаем Ес = 3.1 • 10_2м. Использование критерия (4) дает значение Тес = 84 при максимальной скорости вращения жидкости г>тах = 0.2 м/с. Следовательно, использование такой вязкой жидкости, как глицерин, вполне обеспечивает выполнение требования (4). Условие (5), накладываемое на число Рэлея, также выполняется: в нашем случае Иа = 1.28.
Ъ ^ ^ ЛТ) . (6)
2. Экспериментальная установка и методика измерений
Схема экспериментальной установки показана на рис. 1. В качестве тепловых резервуаров использованы дюралюминиевые бруски шириной 25 мм. В них просверлены сквозные отверстия диаметром 4 мм, равным внутреннему диаметру полипропиленовой трубки. Для возможности крепления трубок с каждой стороны бруска были сделаны дополнительные углубления, соосные каналам, диаметром,
равным внешнему диаметру трубки, что позволило точно совместить внутренний канал в трубке и в бруске. Верхний брусок охлаждался проточной водой. Нижний брусок был прикреплен к дюралюминиевой подставке для уменьшения колебания температуры при нагреве. Обе секции были разнесены с помощью крепежных пластин с низкой теплопроводностью. Радиус тороида в установке равнялся Кс = 30 мм. Измерение температуры осуществлялось с помощью термопар типа хромель-алюмель. В тепловых резервуарах они помещались в специально высверленные отверстия малого диаметра. Сигнал с термопар снимался 4-канальным измерителем температуры Сеи1ег-309, который в реальном масштабе времени отправлял значения измеренной температуры на СОМ-порт компьютера. Разрешение прибора составляет 0.1 К. Абсолютная погрешность измерения температуры складывается из погрешности измерения прибора и погрешности, вносимой термопреобразователем. Общая абсолютная погрешность измерения температуры составляла ±1 К. Отсчеты температуры производились через каждые 5 секунд. На нижней секции поддерживалась постоянная температура при помощи нагревателя, управляемого компьютером. Верхняя секция охлаждалась проточной водой из городского водопровода.
После включения установки верхняя секция охлаждалась до температуры холодной воды в городской магистрали, нижняя же нагревалась до заданной температуры, которая поддерживалась с точностью до 1 градуса включением и выключением нагревателя. Температура воды зависела от времени суток и интенсивности расхода воды в водопроводе. На рис. 2 представлены зависимости ДТ12 = Т2 - Т1 и ДТав = Та - Тв от времени наблюдения, которое в различных экспериментах составляло от 4 часов до 1.5 суток. Тем не менее точность поддержания разности температур ДТ12 в течение эксперимента оказалась вполне удовлетворительной.
В серии из более чем 20 опытов установлено, что однонаправленная (вправо или влево по случайному принципу) конвекция в описанной установке наблюдается до значений ДТ12 ~ 47 К. При больших значениях ДТ12 наблюдалась неустойчивая конвекция с неожиданными, не связанными с колебаниями ДТ12, инверсиями на-
Рис. 1. Схема эксперимента; 1, 2 и А, В - точки подключения термопар
Рис. 2. Наиболее характерные виды зависимостей ЛТ12 и ЛТАВ от времени наблюдения, демонстрирующих неустойчивую конвекцию; средние значения ЛТ12 =78 К (а); ЛТ12 = 97, 5 К (б); ЛТ12 = 107 К (в)
правления течения жидкости в кольце. О скорости течения и его направлении можно было судить по значениям ЛТав . Такого рода скачки становились более частыми и масштабными по достижении ЛТ12 & 90 К.
3. Анализ экспериментальных результатов
На рис. 3 представлены результаты фурье-анализа спектров мощности АТ12 и АТаб для АТ12 = 107 К (см. рис. 2, в). В двойном логарифмическом масштабе оба спектра обнаруживают линейную зависимость от частоты. Из этого можно сделать вывод, что соответствующие временные ряды характеризуются распределением с «тяжелыми хвостами» в виде гиперболической зависимости мощности спектральных гармоник Р от частоты / (распределение Леви - Парето - Ципфа)
Р (/) - /р • (7)
Спектральный показатель в в (7), найденный по тангенсу угла наклона графиков для АТ12 и АТаб, соответственно равен ^12 = 1-4 и $аб = 1-7. Повышение удельного веса низкочастотных гармоник во временной зависимости АТ12 по сравнению с белым шумом (в = 0) связано с включением и выключением время от времени нагревателя. Это производится с целью остановки дрейфа АТ12, вызванного, в основном, колебаниями температуры воды в водопроводе, и поддержания этой разности температур в заданных пределах: АТ12 = АТ12 ± 8. Рис. 2 показывает, что наиболее заметные понижения температуры воды в магистрали, и соответственно рост АТ12, происходит с периодичностью, близкой к 12 часам. Это связано с ростом потребления воды в утренние и вечерние часы. Добавление этих низкочастотных (по сравнению с прочими, случайными факторами) процессов во временные ряды трансформирует их спектр вышеуказанным образом. Полученные значения в показывают, что оба временных ряда обладают признаками розового шума, причем интересующий нас процесс конвекции по своим шумовым характеристикам приближается к коричневому (броуновскому) шуму, для которого в = 2.
Рис. 3. Фурье-спектры мощностей временных зависимостей АТ12 (вверху) и АТав (внизу) для случая, изображенного на рис. 2, в
Степень отклонения значений временного ряда от случайных изменений и наличие долговременной памяти в системе позволяет оценить херст-анализ [8]. Однако, как показано в [9, 10], большей устойчивостью к шуму, более эффективным вычитанием трендов и большей статистической достоверностью при сигналах одинаковой длины обладает модификация херст-анализа - DFA (detrended fluctuation analysis - флуктуационный анализ с удалением трендов). В этом методе вычисляется среднеквадратичная ошибка линейной аппроксимации случайных блужданий k
У (k) = \z (i) — z \ исследуемого временного ряда z(i) = 1, 2,... ,N
i=l
F (n) =
\
1 N
^E \y (k) - yn (k)\5 k=l
и считается, что при наличии масштабной инвариантности (скейлинга) ^ (п) ~ па, где а - так называемая скейлинговая экспонента. Она связана со спектральным показателем в соотношением
в = 2а - 1. (8)
На рис. 4, а показаны результаты применения этого метода к тем же временным рядам, которые использовались выше для фурье-анализа. По тангенсу угла наклона соответствующих прямых найдены а12 = 1.2 и алв = 1-5. Тогда, по формуле (8), в12 = 1-4 и влв = 2.0. Это, в целом, подтверждает выводы, сделанные из фурье-анализа. Для сравнения на рис. 4, б приведены результаты ЭРА хаотизированных (с помощью генератора случайных чисел) временных рядов, полученных из исходных рядов. Так как в хаотизированных рядах значения а оказались равными 0.5, то можно сделать вывод, что исходные ряды имели долговременную память, которая разрушилась в результате рандомизации.
Ни фурье-анализ, ни ЭРА не позволяют ответить на вопрос о степени влияния колебаний разности температур ДТ12 в ячейке на процесс конвекции. Для этого следует сравнить характерные частоты в рядах ДТ12 и ДТлв и их временную динамику, что позволяет сделать вейвлет-анализ [11]. На рис. 5 представлены результаты этого анализа при использовании в качестве базисных вейвлетов гауссианов (рис. 5, а) и вейвлетов Морле (рис. 5, б). Первые лучше локализованы во временном пространстве и поэтому точнее отображают низкочастотные особенности временного ряда. Вторые, напротив, лучше локализованы в частотном пространстве и поэтому точнее отображают его высокочастотные особенности. По оси абцисс отложены временные отсчеты, по оси ординат - масштаб скейлинга а (в логарифмическом виде). Распределение значений вейвлет-преобразования функции f (¿), в качестве которой выступают ДТ12 и ДТлв,
1 ( — А — X
wf <^=а/ Ч t—)>(t) dt
a
—00
условно показано разными оттенками серого цвета.
Анализ скейлограмм приводит к заключению, что на высоких частотах, соответствующих временным интервалам в несколько минут, во временной зависимости
АТ12 разности температур нагреваемой и охлаждаемой частей ячейки присутствует больше шумовых составляющих, чем в зависимости АТаб , характеризующей скорость конвективного движения жидкости. На более низких частотах, с характерными временами в десятки минут, возмущения обеих зависимостей имеют близкие пери-
Рис. 4. а - Результаты БРА для временных зависимостей АТ12 (вверху) и АТав (внизу) из рис. 2, в; б - то же для рандомизированных рядов из этих же измерений
оды, что свидетельствует о возможном влиянии условий нагрева на характер нарушения стационарности конвективного движения. Но наиболее низкочастотные изменения, с характерными временами в несколько часов, не согласованы друг с другом. Отсюда следует, что скачки в направлении и скорости вращения жидкости в то-
О 3 6 9 12 15 18 21 24 27 30 33
Время, час
Рис. 5. Гауссиан- (а) и Морле- (б) скейлограммы зависимостей ДТ12 (вверху) и ДТав (внизу) из рис. 2, в
роидальной ячейке обусловлены не флуктуациями нагрева, а являются следствием неустойчивостей в характере самой конвекции.
В дополнение к вейвлет-анализу изучаемых временных рядов проведен их корреляционный анализ. Связь между значениями одного и того же случайного процесса в разные моменты времени (наличие и продолжительность памяти системы о своих прошлых состояниях) позволяет выявить автокорреляционная функция (АКФ) [12, 13], вычисляемая здесь по формуле
Е (xi - x) (Xi+k - x)
rx (k) =
E (Xi - x)
где X - выборочное среднее временного ряда, а k = 0,1, 2... - сдвиг (lag) между отсчетами. Степень связи различных случайных процессов (в нашем случае AT 12 и ATab) позволяет установить кросскорреляционная функция (ККФ)
' xy
(k) =
Е (Xi - x) (yi+k - y)
E(xi- x)2./E (yi+k- y)
20.0
На рис. 6, а представлены АКФ для рассмотренных выше временных рядов АТ12, АТаб . АКФ разности температур между нижней и верхней точками тороида г 12 (к) имеет вид, близкий к гиперболе, но с тремя локальными максимумами, соответствующими всплескам на графике зависимости АТ12 (¿) (см. рис. 2, в). Наибольшая корреляция между текущими и предшествующими значениями этой зависимости наблюдается перед первым всплеском. Для доверительного интервала по уровню 5% совпадения величина сдвига составляет около 103 отсчетов или 1.4 часа. АКФ разности температур между боковыми точками тороида г2 (к) имеет провал в диапазоне значений к = 0-8+2-0 в область отрицательных значений, что соответствует изменению направления движения жидкости в тороиде.
Кросскорреляционная функция г 12ав зависимостей АТ12 и АТаб показана на рис. 6, б, из которого видно, что наибольшая корреляция обеих зависимостей наблюдается до и после скачка АТаб (соответственно г12ав (к) достигает значений 0.4 и 0.2). Наоборот, в области скачка обе функции максимально нескоррелированы.
^ 0.0 10.0 20.0 О к л О-3
Рис. 6. а - АКФ зависимостей АТ12 и АТав из рис. 2, в; б - ККФ обеих зависимостей; к - сдвиг между отсчетами
Полученные результаты кросскорреляционного анализа говорят, что, хотя и имеется слабая взаимная корреляция (в пределах 20-40%) неустойчивостей в конвективном движении жидкости в тороиде и флуктуаций разности температур в его нижней и верхней точках, тем не менее крупномасштабные нарушения стационарности процесса связаны с внутренними свойствами системы. Это вполне согласуется с выводами вейвлет-анализа.
По результатам эксперимента построен хаотический аттрактор в системе координат, осями которой являются: температура нижней части тора , температура в одной из боковых точек Тл и разность температур в боковых точках тороида ДТлв, пропорциональная угловой скорости вращения жидкости (рис. 7). В данном эксперименте за 36 часов измерений удалось наблюдать лишь одну инверсию потока, причем с несоразмерным с начальным значением периодом вращения. Тем не менее на рисунке отчетливо видно различие в плоскостях изменения фазовых траекторий: (Тл, ДТлв) для вращений в одну сторону и (Тл, Т^) для вращений в другую сторону. Это свойство восстановленной по экспериментальным данным проекции аттрактора на трехмерное пространство сближает его с аттрактором Лоренца. Об истинной размерности вложения данного аттрактора судить сложно из-за сильной зашумленности исходного временного ряда.
Заключение
Осуществлена экспериментальная реализация модели Лоренца на основе тороидальной ячейки радиусом 30 мм с глицерином. С помощью различных методов математической обработки временных рядов измерений (фурье-анализ, ЭРА, вейвлет-и корреляционный анализы) были исследованы возможные режимы конвекции жидкости в ячейке. Выяснилось, что при разности температур в нижней и верхней точках ячейки ДТ12 < 47 К конвективное течение имеет преимущественно регулярный характер: по или против часовой стрелки. Причем, в выборе направления вращения нет никакой закономерности, оно носит случайный характер. Начиная с ДТ12 > 47 К, и в особенности с 80-90 К, конвекция становится неустойчивой, с неожиданными инверсиями направления вращения жидкости. Было доказано наличие долговременной памяти в рядах значений ДТ12, постоянство которых поддерживалось с точностью ±1 К, и в рядах самой измеренной величины ДТлв, пропорциональной скорости вращения жидкости. Фурье-анализ и ЭРА показали, что в условиях нагревания, близкого к розовому шуму, хаотические свойства конвективного движения жидкости близки к броуновскому шуму. Вейвлет-анализ и корреляционный анализ позволили установить, что конвективный процесс слабо коррелирует с особенностями условий нагрева ячейки: все долгопериодические закономерности, характеризующие
Рис. 7. Хаотический аттрактор, построенный по экспериментальным данным из рис. 2, в
неустойчивые режимы конвекции, развиваются самостоятельно. По результатам измерений можно построить хаотический аттрактор, похожий на аттрактор Лоренца. Правда, за 36 часов наблюдений удалось наблюдать лишь однократный перескок из одной плоскости фазового пространства в другую.
Таким образом, в лабораторном макете модели Лоренца на основе вертикального тороида с глицерином можно наблюдать ряд теоретически предсказываемых данной моделью случаев регулярной и хаотической конвекции, в том числе явление детерминированного хаоса.
Библиографический список
1. Lorenz ^.^.Deterministic nonperiodic flow // J. Atmos. Sci. 1963. Vol. 20, № 2. P. 130.
2. Welander P. On the oscillatory instability of a differentially heated loop // J. Fluid Mech. 1967. Vol. 29, Pt 1. P. 17.
3. Creveling H.F. et al. Stability characteristics of a single-phase free convection loop // J. Fluid Mech. 1975. Vol. 67, Pt 1. P. 65.
4. Gorman M., Widman P.J., Robins K.A. Chaotic flow regimes in a convective loop // Phys. Rev. Lett. 1984. Vol. 52, № 25. P. 2241.
5. Wang Y., Singer I., Bau H. Controlling chaos in thermal convecting loop // J. Fluid Mech. 1992. Vol. 237. P. 479.
6. Дроздов С.М.Экспериментальное исследование конвекции жидкости в замкнутом тороидальном канале // Изв. РАН. МЖГ. 1995, № 4. С. 20.
7. Дроздов С.М. Моделирование возникновения нестационарности и хаоса в гидродинамической системе, управляемой небольшим числом степеней свободы // Изв. РАН. МЖГ. 2001, № 1. С. 31.
8. Шредер М.Фракталы, хаос, степенные шумы. Миниатюры из бесконечного рая. Ижевск: НИЦ «Регулярная и хаотическая динамика», 2001.
9. Peng C.-K., Buldyrev S.V., Havlin S., Simons M., Stanley H.E., Goldberger A.L. Mosaic organization of DNA nucleotides // Phys. Rev. E. 1994. Vol. 49. P. 1685.
10. Peng C.-K., Havlin S., Stanley H.E., Goldberger A.L. Quantification of scaling exponents and crossover phenomena in nonstationary heartbeat time series // Chaos. 1995. Vol. 5. P. 82.
11. Астафьева Н.М. Вейвлет-анализ: основы теории и примеры применения // УФН. 1996. Т. 166, № 11. С. 1145.
12. Hamilton J.D. Time series analysis. Princeton University Press, 1994.
13. Box G.E.P., Jenkins G.M., Reinsel G.C. Time series analysis: Forecasting and control. Third edition. Prentice Hall, 1994.
Ульяновский государственный Поступила в редакцию 27.04.2006
технический университет После доработки 4.07.2006
ФНЦП ОАО «Марс», Ульяновск
EXPERIMENTAL REALIZATION OF LORENZ MODEL OF LIQUID'S CONVECTIVE INSTABILITY IN VERTICAL TOROIDAL LOOP
R.A. Brazhe, O.N. Kudelin
Stable and unstable regimes of glycerine convection in vertical toroidal loop are investigated experimentally. The results of Fourier-analysis, DFA, wavelet-, and correlation analysis of liquid's motion peculiarities are presented. Chaotic attractor with Lorenz-attractor signs is constructed.
Браже Рудольф Александрович - родился в селе Ковалево Омской области (1948), окончил Саратовский государственный университет (1971) и там же аспирантуру по специальности «радиофизика» (1974). С 1974 года работает в Ульяновском государственном техническом университете; с 1993 года - заведующий кафедрой физики. Защитил диссертацию на соискание ученой степени кандидата физико-математических наук в СГУ (1975) и доктора физико-математических наук в УлГТУ (2002) в области радиофизики и математического моделирования нелинейных процессов. Автор более 250 научных трудов, в том числе одной монографии и нескольких учебных пособий. Член-корреспондент РАЕН.
Куделин Олег Николаевич - родился в Ульяновске (1979), окончил Ульяновский государственный университет по специальности «физика твердого тела» (2001). Работает на ФНЦП ОАО НПО «Марс» инженером-программистом, автор 15 научных работ.
