5. Динамічні системи: використання SciPy¶

Зміст глави

  • Презентація модуля SciPy
  • Виконання розширеної обробки даних за допомогою SciPy
  • Генерація вдосконалених фігур за допомогою Matplotlib
  • Напишіть документацію функції

У рамках наукових модулів Python є вражаюча колекція функцій, які виконують різні операції над цифровими масивами NumPy. Перш за все, їх слід використовувати, а не винаходити колесо, перекодуючи (можливо, менш добре) функціональність, яка вже існує в інших місцях !

В NumPy, ми вже бачили, що існує багато операцій, що дозволяють, наприклад,

  • генерувати конкретні масиви: np.arange, np.ones, np.linspace, .
  • виконувати операції зі значень масиву: np.sum, np.sin, np.histogram тощо.
  • виражати (булеві) умови на елементах масиву або витягувати частини масиву: np.logical_not, np.in1d
  • змінити розташування значень масиву, його форму або навіть створити новий масив з кількох інших: np.reshape, np.concatenate .
  • тощо. .

Модуль SciPy тут набір інструментів числові для масивів NumPy. Ми знаходимо в SciPy класичні операції з обробки/обробки цифрових даних, але специфічні для певного типу додатків (лінійна алгебра, статистика тощо). Отже, це функції більш “високого рівня”, ніж функції NumPy.

SciPy - це стабільний, добре перевірений та відносно добре задокументований модуль.

Модуль SciPy виконує різні операції з числовими масивами (ndarray) numpy. Тому ми можемо безпосередньо використовувати ці масиви як аргументи для різних функцій

Щоб продемонструвати використання SciPy, ми зупинимось на інтегруванні диференціальних рівнянь, розглядаючи динамічні системи на основі механічних маятників.

5.1. Простий маятник¶

системи

Написавши фундаментальне відношення динаміки (збереження імпульсу), ми отримуємо просте рівняння маятника, задане формулою

де - кут маятника відносно вертикалі, і ми позначимо точкою похідну від часу.

Для малих коливань ми можемо зробити апроксимацію. Коли наближення не є дійсним, це диференціальне рівняння потрібно інтегрувати чисельно, щоб отримати еволюцію положення та кутову швидкість маятника з часом.

Отже, нам потрібно мати інтегратор диференціальних рівнянь, який можна було б очікувати у SciPy. Але яка відповідна функція? Давайте спробуємо в Google знайти “scipy інтегрувати диференціальне рівняння”, або перегляньте Зміст довідки http://docs.scipy.org/doc/scipy/reference/index.html. Існує інтегрований підмодуль, який сам містить функцію odeint

Для початку ми повинні поставити диференціальне рівняння 2-го порядку у вигляді системи рівнянь першого порядку

Тепер ми можемо інтегрувати траєкторію з початкової умови

Як і функція odeint, існує велика кількість функцій, що дозволяють виконувати наукові обчислювальні роботи високого рівня.

Наприклад, ми можемо перевірити збереження механічної енергії з часом

Щоб сформувати малюнок вище, ми використали ряд команд з модуля matplotlib .

Для побудови різних елементів фігури (кривої, міток, легенди тощо) існують спеціальні функції matplotlib, які можна "виявити" за допомогою документації matplotlib http://matplotlib.org/. Зокрема, приклад галереї http://matplotlib.org/gallery.html дуже корисний для того, щоб побачити, як генерувати різні типи фігур.

Напишіть сценарій python для побудови простої діаграми фаз маятника, показаної нижче. Для цього вам потрібно

  • інтегрувати диференціальне рівняння для різних початкових умов від 0 до .
  • представляють як функцію для різних рішень.

5.2. Примусовий маятник¶

Можна додати ступінь свободи динамічній системі маятника, наклавши силу, що залежить від часу, наприклад, періодичну в часі. Ми також розглянемо демпфування маятника через тертя повітря. Потім ми розглянемо наступне рівняння для вимушеного маятника

Як і раніше, ми нормалізуємо час за допомогою фіксації. Тоді є три параметри, за якими ми можемо грати: A, q та. Тому ми визначаємо нову функцію, що визначає динамічну систему першого порядку

Для зручності ми також визначимо функцію, що дозволяє інтегрувати траєкторію з початкової умови за певний інтервал часу:

Тепер ми можемо інтегрувати траєкторію:

Для цієї траєкторії еволюція періодичного залишку з тим самим періодом, що і форсуючий період

Коли амплітуда форсування збільшена, рух залишається періодичним після досить тривалого перехідного режиму. З іншого боку, ми спостерігаємо подвоєння періоду, амплітуда коливань коливається між двома періодами форсування

Нарешті, коли ми додатково збільшуємо амплітуду форсування, рух перестає бути періодичним і здається "хаотичним":

Ми підтвердимо спостереження за періодом подвоєння періоду завдяки дискретному перетворення Фур'є сигналів. Для дискретного перетворення Фур'є існує також підмодуль SciPy:

Спочатку потрібно вибрати ту частину сигналу, яка здається періодичною, тобто видалити перехідний стан. Для цього ми використовуємо логічну маску

Оскільки нульова частота нас не цікавить, ми віднімаємо середнє значення сигналу

Зверніть увагу, що ми можемо відняти скаляр із масиву ndarray: тоді віднімання виконується елемент за елементом.

Обчислюємо дискретне перетворення Фур'є:

Для представлення TF нам також потрібна таблиця частот, що відповідає таблиці fft_theta_0. Для цього є спеціальна функція в модулі fftpack, fftfreq:

Тоді ми можемо представити модуль перетворення Фур'є як функцію від частоти

Помічено, що основна частота сигналу відповідає частоті примусу. Ми також спостерігаємо потрійну гармоніку в .

Для проміжного форсування ми повторюємо ту ж серію операцій

Потім ми спостерігаємо новий пік на частоті половин, що відповідає подвоєнню періоду.

Щоб піти далі, ми можемо візуалізувати еволюцію спектру для великої кількості параметрів за допомогою кольорового кодування амплітуди модуля перетворення Фур'є:

Потім ми спостерігаємо каскад переходу до хаосу, з хаотичними вікнами та іншими періодичними виданнями.

Однією з властивостей хаотичних динамічних систем є їх велика чутливість до початкових умов: дві спочатку дуже близькі траєкторії з часом швидко віддаляться.

У хаотичному режимі ми перевіримо, що дві спочатку дуже близькі траєкторії дуже швидко віддаляються. Для параметрів

ініціалізувати дві траєкторії, розділені дистанцією eps (початкові умови та, наприклад), і представляти їх відстань у фазовому просторі з часом

Щоб перевірити, чи відстань зростає з часом в геометричній прогресії, зробіть графічне представлення в напівлогарифмічних координатах (plt.semilogy).

Для подальшого

Ми можемо представити в 3-D граничний цикл, що описує динаміку маятника.

5.3. Вхід/вихід: як зберегти масиви та завантажити файли¶

Наразі ми ще не зберегли жодного результату! Ви можете легко зберегти масив NumPy у двох форматах:

  • у тексті ascii з np.savetxt (якщо вам потрібно переглянути таблицю в текстовому редакторі) для таблиць розмірів np.save (для менших файлів у пам'яті)

Для масивів розміром три або більше використовуйте np.save

Так само є дві функції для завантаження масиву numpy з файлу

Спеціалізовані модулі можуть зберігати таблиці в інших форматах: csv, hdf5 тощо.