[PYTHON] Numerical calculation with Phython while learning nonlinear dynamics and chaos [1]

Contents of this series

This article is by Steven H. Strogatz ["Strogatz Nonlinear Dynamics and Chaos"](https://www.amazon.co.jp/%E3%82%B9%E3%83%88%E3%83%AD% E3% 82% AC% E3% 83% 83% E3% 83% 84-% E9% 9D% 9E% E7% B7% 9A% E5% BD% A2% E3% 83% 80% E3% 82% A4% E3 % 83% 8A% E3% 83% 9F% E3% 82% AF% E3% 82% B9% E3% 81% A8% E3% 82% AB% E3% 82% AA% E3% 82% B9-Steven-H -A summary of nonlinear dynamics and chaos created based on Strogatz / dp / 4621085808) and implementation code using Python. To give a brief self-introduction, the author of this article is a mechanical student at a local national university, and the mathematics that appears in this article is basically self-taught. (So I'm worried about the accuracy.) I think it's better to think about the equations that appear in this book as a memorandum implemented by an amateur while thinking alone. We plan to add the contents of other chapters one by one. If you find any mistakes, I would be very grateful if you could point them out!

0. Execution environment

All the code in the article is the one executed by Python3.8.2.

1. One-dimensional flow

In this article, in accordance with Strogats' book, We will deal with one-dimensional equations → two-dimensional equations → three-dimensional equations (chaos appearance) in this order. Let's first consider a one-dimensional equation, that is, a function that can be expressed in the following form.

\frac{dx}{dt} = f(x)

1.1 Fixed points and their stability

To find the solution of *** $ x'= f (x) $ *** starting from any initial condition *** $ x_0 $ ***, use a virtual particle *** $ x_0 Place it in $ *** and observe how it is carried by the flow. Over time, this point moves along the *** $ x $ *** axis according to some function *** $ x (t) $ ***. This function is called the ** orbit of the solution ** starting from *** $ x_0 $ *** and represents the solution of the differential equation starting from the initial condition *** $ x_0 $ ***. The one that represents all the qualitatively different activations of the system is called a ** phase diagram **. In this phase diagram, the point where *** $ x'= 0 $ *** is called ** fixed point **.

[Example 1.1] Find all the fixed points of
*** $ x'= x ^ 2-1 $ *** and classify them.
[Sample code 1.1]

example1_1.py


import sympy as sym
from sympy.plotting import plot
x = sym.Symbol('x')
dx = x**2 - 1
plot(dx, (x, -10, 10), aspect_ratio='auto', title='example1.1', xlim = (-5, 5), ylim = (-3, 25),ylabel="x'")
[Execution result 1.1] example1_1.png
[Answer 1.1]
As shown in the figure, there are two fixed points *** $ x = -1, x = 1 $ ***. The former is stable because the restoring force works even for small perturbations near the fixed point. On the other hand, the latter is not stable (called ** sink **) because the point moves in the direction of amplifying a small perturbation. (Called ** source **)

1.2 Linear stability analysis

Now, let us consider a method for analyzing stability near a fixed point as shown in the previous example. *** $ x ^ \ ast $ *** as a fixed point, *** $ \ eta (t) = x (t)-x ^ \ ast $ *** as *** $ x ^ \ ast $ Let it be a small perturbation from ***.

\begin{eqnarray}
\dot{\eta} &=& \frac{d}{dt}(x-x^\ast)=\dot{x}\\
&=& f(x)\\
&=& f(x^\ast+\eta)\\
&=& f(x^\ast)+\eta \cdot f'(x^\ast) + o(\eta^2)\\
&\simeq& \eta \cdot f'(x^\ast)
\end{eqnarray}

Therefore, the linearization around a fixed point is

\dot{\eta} \simeq \eta \cdot f'(x^\ast)

It can be seen that it can be expressed as. (Taylor expansion was used in the transformation from the 3rd line to the 4th line in the expression expansion) Based on this result

* $ f'(x ^ \ ast)> 0 $ * will grow exponentially, If
* $ f'(x ^ \ ast) <0 $ *, it will be attenuated.
* However, these are valid in the case of * $ f'(x ^ \ ast) \ neq0 $ * You can see that.

1.3 Existence and uniqueness of the solution

From here, a very important theorem for considering non-linear dynamical systems.

*** $ \ dot {x} = f (x) , x (0) = x_0 $ ***. Assume that *** $ f (x) $ *** and *** $ f'(x) $ *** are continuous on the open interval R on the *** $ x $ *** axis. It is assumed that *** $ x_0 $ *** is a point in R. Then, this initial value problem is solved by *** $ (-\ gamma, \ gamma) $ *** while walking around *** $ t = 0 $ ***. *** $ x (t) $ With ***, this solution is unique. In other words, if *** $ f (x) $ *** is smooth enough, then the solution exists and is unique </ font>.