Этот топик продолжает серию моих статей на Хабре, посвященных исследованию аттрактора Лоренца.

Часть 1. Критический взгляд на аттрактор Лоренца
Часть 2. Динамическая система Лоренца и вычислительный эксперимент
Часть 3. О существовании периодических решений в системе Лоренца
Часть 4. Три цикла в аттракторе Лоренца

Итак, рассмотрим нелинейную систему дифференциальных уравнений, введенную Эдвардом Лоренцом в 1963 году:

$$

где

$$

классические значения параметров системы.

Изучая иностранную литературу, на днях наткнулся на работы [1, 2] профессора Мичиганского университета Дивакара Вишваната (Divakar Viswanath) об итерационном алгоритме вычисления периодических орбит динамических систем, основанном на методе Линдштедта-Пуанкаре (ЛП) (для ознакомления с ним рекомендую книгу [3, с. 408-411]). Преимуществом данного метода является то, что он не требует численного интегрирования дифференциального уравнения, поэтому может быть применён к построению и неустойчивых циклов.

Это третий мой топик на Хабре (часть 1 и часть 2), посвященный динамической системе Лоренца. Я продолжаю заниматься исследованием вопроса о существовании периодических решений (циклов) в этой системе. Удалось получить интересный результат при определенном соотношении ее параметров.

При построении и анализе поведения решений систем обыкновенных дифференциальных уравнений иногда требуется определять матричную экспоненту [1]. Классический метод связан с тем, что приходится рассчитывать большие степени матриц. В данном топике рассматривается алгоритм приближенного вычисления матричной экспоненты, который за фиксированное число матричных операций дает результат с заданной точностью. Проведен вычислительный эксперимент с целью анализа эффективности алгоритма.

Данный пост является продолжением моей статьи [1] на Хабрахабре об аттракторе Лоренца. Здесь рассмотрим метод построения приближенных решений соответствующей системы, уделив внимание программной реализации.

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



где , r и b являются положительными числами. Для исследования поведения решений системы обычно берут классические значения параметров системы, т.е.



Как было отмечено в статье [1], в этом случае в системе (1) имеет место неустойчивость ее решений на аттракторе. По сути, это делает некорректным применение классических численных методов на больших отрезках времени (а на таких отрезках и строятся предельные множества динамических систем). Одним из вариантов преодоления этой проблемы является переход к высокоточным вычислениям, но такой подход ставит исследователя в жесткие рамки: во-первых, малая степень свободы для уменьшения ошибки (изменение величины шага интегрирования и точности представления вещественного числа для управления вычислительным процессом), во-вторых, большой объем вычислений при очень малых .

Intro

Я публикую этот топик как обучающий. Собственно говоря, существенной новизны в материале нет, тема заезжена. Думаю, что интересным будет подход к решению задачи.

Помню, на первом курсе на занятиях по математическому анализу пришел в голову один интеграл. Преподаватель вызвал к доске, но прозвенел звонок. По дороге домой в автобусе сложился «скелет» решения кубического уравнения. Общая схема, конечно, не самая рациональная. Есть более эффективная — тригонометрическая формула Виета. Там сразу выписывается корень по виду уравнения, а, вообще, по объему вычислений все-таки лучше использовать численный метод Ньютона, поскольку степенные ряды для обратных тригонометрических функций сходятся медленно (по ним строятся вычисления таких функций в некоторых библиотеках). Вот что получилось.

В прикладной математике иногда возникает задача построения периодических решений нормальной системы обыкновенных дифференциальных уравнений вида



где функция представляет собой сумму



многомерного многочлена и тригонометрического полинома , являющегося -периодической векторной функцией.

Многие из теорем существования периодических решений системы (1) используют тот фундаментальный факт, что такие решения полностью определяются неподвижными точками оператора сдвига по траекториям системы. Однако использование данных теорем для непосредственного нахождения нужного периодического решения, скорее всего, не представляется возможным.

1. Об аттракторе Лоренца

Эдвард Нортон Лоренц (1917 – 2008) является основателем теории хаоса, очень популярной в науке на сегодняшний день. Он учился в колледже Дартмут штата Нью-Гемпшир США и Гарвардском университете в Кембридже. Во время Второй мировой войны служил метеорологом в авиационном корпусе армии США, потом до конца своих дней работал профессором в Массачусетском технологическом институте.

В 1963 году в журнале «Journal of the Atmospheric Sciences» вышла его статья «Deterministic Nonperiodic Flow» (русский перевод: Лоренц Э. Детерминированное непериодическое течение // Странные аттракторы. — М.: Мир, 1981, с. 88-117), заложившая не только основы теории хаоса, но и изменившая представления о моделировании погодных явлений. В этой работе из системы уравнений Навье-Стокса впервые была получена нелинейная автономная система обыкновенных дифференциальных уравнений третьего порядка (динамическая система), описывающая движение воздушных потоков в плоском слое жидкости постоянной толщины при разложении скорости течения и температуры в двойные ряды Фурье с последующем усечением до первых-вторых гармоник:

В булевой алгебре, используемой для моделирования работы цифровых схем, не учитывается время перехода (или инертность) логических элементов (например, «И», «ИЛИ») из одного состояния в другое (например, из 0 в 1). В случаях, когда время распространения сигнала внутри элемента достаточно мало, задержку перехода можно не учитывать. Но с повышением частоты изменения входных сигналов в реальных схемах начинает сказываться влияние времени распространения сигнала внутри ее элементов. Такие задержки могут порождать неустойчивую работу устройств (т.е. в сигналах после последовательного прохождения через узлы схемы появляются не учитываемые моделью схемы переходы, называемые сбоями). Многие производители современных процессоров держат в тайне, как они борются со сбоями, порождаемыми задержками, на частотах порядка ГГц. Ведь обычный проводник с множеством близких изгибов на плате в таком режиме работы уже превращается в индуктивность.

Чтобы проанализировать простейшие схемы на риски сбоя, обычно применяют метод временных диаграмм, ставший уже классическим. Сигналы в каждом узле прочерчиваются строго друг под другом; где нужно, производится искусственная задержка перехода из одного состояния в другое, а далее по булевому представлению строятся выходные сигналы. Данный метод не хорош тем, что требует работы с графиками, что может внести погрешность в получаемые сигналы. А нам нужно знать, возникнет ли сбой, и какой он будет иметь вид.

В 2005-м году я касался этой темы и изложил свои мысли в данной заметке на английском языке. Однако мой хороший знакомый посоветовал написать в Хабрахабр для обсуждения с пользователями.

В данном топике процедура анализа схем переведена из графического представления в математическое представление. При этом вводится в булеву алгебру параметр «время» посредством замены логических операций арифметическими. Для простоты анализа рассматриваются асинхронные, т.е. неуправляемые внешним (синхронизирующим или тактирующим) сигналом, цифровые схемы.

Начнем…

Хочу рассказать о двух странностях, с которыми мне пришлось столкнуться, программируя вычислительные алгоритмы на языке C.

Итак, первое неожиданное поведение для некоторых программистов. Вот маленькая прога.

#include <stdio.h>
int main()
{
    unsigned char a = 1, b;
    b = ~a >> 1;
    printf("%u\n", b);
    return 0;
}

Разберем ее. Поразрядная операция ~ инвертирует состояние каждого бита байта a, в который изначально записана единица. В результате должны получить 11111110b, то есть 254. Сдвигая этот байт вправо на один бит, должны получить 127. Однако код, который дает, например, компилятор gcc, выводит в консоль число 255?!

Сначала я подумал о том, что дело в приоритете — вдруг у компилятора приоритет операций «косячит»? То есть будто бы сначала делается сдвиг, а потом — инверсия (а что, логично...). Так в чем же дело?