UUUMエンジニアブログ

UUUMのエンジニアによる技術ブログです

カオスの鑑賞

はじめ

こんにちは。UUUMのシステムユニットの柿木です。
今日はカオスを鑑賞していきたいと思います。
用意するのは、python3の実行環境です。
turtleライブラリを利用してカオスを鑑賞します。
* 今回用意したコードは結構PC唸るので、実行は自己責任でお願いします。

そもそもカオスとは?

wikipediaによると

次の3つの条件を満たす写像f: V → Vは、Vの上でカオス的であるといえる[35]。
初期条件に鋭敏に依存する。
位相的に推移的である。
周期点はVにおいて稠密である。

という定義があります。
つまりカオスということですね。

前提

開発用ログファイル 座標など見たい時に使います。 以下ファイルは全て同じ階層です。

IS_LOG_SHOW = False

def log(s):
    if IS_LOG_SHOW:
        print(s)

ロジスティック写像

参考URL
定義がよく分からないので、さっそく実例を見ていきましょう。
ロジスティック写像は

 f(x) = \alpha X(1-X)

で、定義される関数があるとき

 a_1=f(x), a_2=f(f(x)), a_3=f(f(f(x))) \cdots

という風に定義される数列のことをいいます。
これを実装したのが以下のファイルです。 \alpha=4です。
*3000回以上やると描画の関係か重くなってしまうので実行の際は注意してください。

from log import log
from turtle import *

SPEED = 10
SHAPE = 'circle'
LINE_COLOR = 'red'
TIRTLE_COLOR = 'yellow'
SCALE_X = 150
SCALE_Y = 500
N = 1000

def do():
    shape(SHAPE)
    speed(SPEED)
    color(LINE_COLOR, TIRTLE_COLOR)

    y = 0.6

    for x in range(N):
        if x == 0:
            continue
        y = f(y)
        x = symmetric_x(x, 5)
        log(f'(x, y) = {x},{y}')
        setpos(x * SCALE_X, symmetric_y(y) * SCALE_Y)
    
    done()


def f(x):
    return 4 * x * (1 - x)

def symmetric_x(x, l):
    if int(x /l)%4 == 0:
        return x%l
    elif int(x /l)%4 == 1:
        return l - x%l
    elif int(x /l)%4 == 2:
        return -(x%l)
    else:
        return - (l - x%l)

def symmetric_y(y):
    return y - 0.5

これをmain関数から呼び出すと
(以下同じような呼び出しなので次からは省略)

import logistic

def main():
    logistic.do()

main()

以下のような軌跡を描きます。

logistic

はい、まさにカオスですね。
 f(x) = 4X(1-X)という単純な規則を持ちつつも、動きがランダムで0と1の間でしか動いていません。
とても面白いですね。

ラングドンのアリ

参考URL
これはカオスなのか?って思いましたが、見た目がカオスっぽいので入れました。
ラングドンのアリの考案者であるクリストファー・ラングドンさんはカオス研究していたみたいなので、きっとこれはカオスなのでしょう。
さて本題に入ります。
ラングドンのアリはシンプルな2つの規則により構成されます。

白いマスにアリがいた場合、90°右に方向転換し、そのマスの色を反転させ、1マス前進する。
黒いマスにアリがいた場合、90°左に方向転換し、そのマスの色を反転させ、1マス前進する。  

これを12000回動かしてみましょう。
コードは以下になります。

from log import log
from turtle import *

SPEED = 10
PEN_SIZE = 10
SHAPE = 'turtle'
LANGTONS_LIMIT_X = 1001
LANGTONS_LIMIT_Y = 1001
N = 12000

def do():
    pen(fillcolor='black', pencolor='black')
    pensize(PEN_SIZE)
    speed(SPEED)
    shape(SHAPE)

    langtons_map = init_langtons_map(LANGTONS_LIMIT_X, LANGTONS_LIMIT_Y)
    x, y = 501, 501
    direction = 1

    for i in range(N):
        if langtons_map[x][y]:
            langtons_map[x][y] = False
            left(90)
            pen(pencolor='white')
            forward(10)
            direction = (direction - 1 + 4)%4
            x, y = set_x_y(x, y, direction)   
        else:
            langtons_map[x][y] = True
            right(90)
            pen(pencolor='black')
            forward(10)
            direction = (direction + 1)%4
            x, y = set_x_y(x, y, direction)
            
        log(f'x is {x}, y is {y}, direction is {direction}')
    done()

def init_langtons_map(n, m):
    langtons_map = []
    for i in range(n):
        temp_list = []
        for j in range(m):
            temp_list.append(False)
        langtons_map.append(temp_list)
    return langtons_map

def set_x_y(x, y, direction):
    if direction == 0:
        y += 1
    elif direction == 1:
        x += 1
    elif direction == 2:
        y -= 1
    else:
        x -= 1
    return x, y

実行結果が以下になります。

langtons

本来はマス目の塗りつぶしなので綺麗になるのですが、turtleライブラリはあくまで軌跡を描くものになるのでちょっと塗り残しが発生してしまします。
しかし、これはこれで趣がありますね。
また、上の方でまっすぐ道を作ってるのも興味深いです。

ローレンツ方程式

参考URL
ローレンツ方程式は、ローレンツさんが大気現象のモデルを研究していた時に発見した方程式です。 方程式は以下になります。

 {\displaystyle {\frac {dx}{dt}}=-px+py}
 {\displaystyle {\frac {dy}{dt}}=-xz+rx-y}
 {\displaystyle {\frac {dz}{dt}}=xy-bz}

左辺は時間微分です。 なので、コードは以下のようになります。
 p=10, r=28, b=\frac{8}{3}, h=0.01です。

from log import log
from turtle import *

SPEED = 10
SHAPE = 'circle'
LINE_COLOR = 'red'
TIRTLE_COLOR = 'yellow'
SCALE_X = 10
SCALE_Y = 10
N = 10000

def do():
    speed(SPEED)
    shape(SHAPE)
    color(LINE_COLOR, TIRTLE_COLOR)

    p = 10
    r = 28
    b = 8/3
    h = 0.01
    x = 1
    y = 1
    z = 1

    setpos(x * SCALE_X, y * SCALE_Y)
    for i in range(N):
        x, y, z = f(x, y, z, p, r, b, h)
        log(f'x is {x}, y is {y}')
        setpos(x * SCALE_X, y * SCALE_Y)
    
    done()

def f(x, y, z, p, r, b, h):
    x_next = x - (p * h * x) + (p * h * y)
    y_next = y - (h * x * z) + (r * h * x) - (h * y)
    z_next = z + (h * x * y) - (b * h * z)
    return x_next, y_next, z_next

実行結果は以下のようになります。 xy平面の断面図です。

lorents

今まで動かしてきたのはカクカクとしたものだったので、カーブが描かれるとまた新鮮に感じますね。

レスラー方程式

参考URL
レスラー方程式は、上記のローレンツ方程式を参考に単純化されたものになります。
方程式は以下になります。

 {\displaystyle {\frac {dx}{dt}}=-y-z}
 {\displaystyle {\frac {dy}{dt}}=x+ay}
 {\displaystyle {\frac {dz}{dt}}=b+xz+cz}

ローレンツ方程式と比べると少し簡単になっていますね。
左辺は時間微分となっているので、ローレンツ方程式の時と同様に求めます。
なので、コードは以下のようになります。
 a=0.2, b=0.2, c=5.7, h=0.1です。

from log import log
from turtle import *

SPEED = 10
SHAPE = 'circle'
LINE_COLOR = 'red'
TIRTLE_COLOR = 'yellow'
SCALE_X = 20
SCALE_Y = 20
N = 10000

def do():
    shape(SHAPE)
    speed(10)
    color(LINE_COLOR, TIRTLE_COLOR)

    a = 0.2
    b = 0.2
    c = 5.7
    h = 0.1
    x = 10
    y = 1
    z = 1
    
    for i in range(N):
        x, y, z = f(x, y, z, a, b, c, h)
        setpos(x * SCALE_X, y * SCALE_Y)
    
    done()

def f(x, y, z, a, b, c, h):
    x_next = x - (y * h) - (z * h)
    y_next = y + (x * h) + (a * y * h)
    z_next = z + b + (x * z * h) - (c * z * h)
    return x_next, y_next, z_next 

実行結果は以下のようになります。 xy平面の断面図です。

rossler
ローレンツ方程式よりは単純な動きをしていますね。

最後に

いかがでしたでしょうか?
私の拙い言葉でどれだけカオスの魅力をどれだけ伝えられたか分かりませんが、1ミリでもカオスの魅力が伝わったのなら幸いです。
今回紹介したコードでは回数、描くスピード、色、縮尺など変えれるようにしたのでぜひ遊んでみてください。