Pythonでカオス・フラクタル
昔,JAVAやらMATLABやらで色々と描いた気がする.ただ,それをそのまま移植すると,ループを回すだけになってしまうので,できる限りPythonicなコードを考えたい.でも,MATLABを使っている時からそうだったけど,ベクトル演算が本当に苦手.色々と参考にしながら,パッと出てくる様になりたい.
参考:
[1] python Mandelbrot マンデルブロ を python で描画する-(仮)タイトル 備忘録
[2] Pythonでカオス・フラクタルを見よう!-Qiita
[3] Python + gnuplot でマンデルブロ集合を描く-剥いだ布団とメモランダム.old
画像(2次元データ)として表したり,等高線で表現したり,散布図として頑張ったり,様々な方法があって面白い.Pythonの場合,複素数はそのまま扱えるし,ブロードキャストを巧く扱えるかどうかが肝になってくるんだろう.x, yをnp.newaxisでブロードキャストする方法を最初考えたけど,そういう時にnp.ogird[]というインスタンス(関数じゃない)が便利らしい.「x = y = np.linspace()」みたいになっている時は,これに置き換えられないか考えれば良さそう.
import numpy as np import matplotlib.pyplot as plt %matplotlib inline x, y = np.ogrid[-2.0:2.0:1000j,-2.0:2.0:1000j] c = x + 1j*y z = 0 for i in range(50): z = z**2 + c mask = np.abs(z) < 4 _, ax = plt.subplots(figsize=(8, 6)) ax.imshow(mask.T,extent=[-2, 2, -2, 2])
閃かなかったのでそのまま書き下した.
import matplotlib.pyplot as plt import seaborn as sns %matplotlib inline sns.set_style('white') fxy = lambda a, b, x, y : [1.0 - a * x**2 + y, b * x] a = 1.4 b = 0.3 xtemp = ytemp = 0 x = [xtemp] y = [ytemp] for i in range(10000): xtemp, ytemp = fxy(a, b, xtemp, ytemp) x.append(xtemp) y.append(ytemp) _, ax = plt.subplots(figsize=(8, 6)) ax.scatter(x, y, color='c')
参考:
[1] Pythonでカオス・フラクタルを見よう!-Qiita
import matplotlib.pyplot as plt import seaborn as sns %matplotlib inline sns.set_style('white') def logistic(a): x = [0.8] for i in range(400): x.append(a * x[-1] * (1 - x[-1])) return x[-100:] al = [2.0 + i * 0.002 for i in range(1000)] _, ax = plt.subplots(figsize=(8, 6)) for a in al: ax.plot([a]*len(logistic(a)), logistic(a), "c.", markersize=1.7)
・樹木曲線(フラクタル・ツリー)
Turtleを使ってもすごく簡単そう.
参考:
[1] Fractal Tree-Python3 codes
import math import matplotlib.pyplot as plt import seaborn as sns %matplotlib inline sns.set_style('white') _, ax = plt.subplots(figsize=(8, 6)) xs = 250 ys = 0 x = [xs] y = [ys] def drawTree(x1, y1, angle, depth): if depth: x2 = x1 + math.cos(math.radians(angle)) * depth * 10.0 y2 = y1 + math.sin(math.radians(angle)) * depth * 10.0 x.append(x2) y.append(y2) drawTree(x2, y2, angle - 15, depth - 1) drawTree(x2, y2, angle + 15, depth - 1) else: ax.plot(x, y) drawTree(xs, ys, -90, 12)
失敗.
import math import matplotlib.pyplot as plt import seaborn as sns %matplotlib inline sns.set_style('white') _, ax = plt.subplots(figsize=(8, 6)) xs = 250 ys = 0 def drawTree(x1, y1, angle, depth): if depth: x2 = x1 + math.cos(math.radians(angle)) * depth * 10.0 y2 = y1 + math.sin(math.radians(angle)) * depth * 10.0 ax.plot([x1, x2], [y1, y2], 'c-') drawTree(x2, y2, angle - 15, depth - 1) drawTree(x2, y2, angle + 15, depth - 1) drawTree(xs, ys, -90, 12)
参考[1]が凄い綺麗.コードそのままで,実行結果の画像が以下.深さで色を変えるだけでここまで綺麗になるんだね.後,matplotlibは遅いのか,やっている事は一緒の筈だけど,これの方が遥かに速い.(Pillowが高速というのは知っていたけど,ここまで違うとは思わなかった.)
取り敢えず,深さで色を変えてみた.
import math import matplotlib.pyplot as plt import matplotlib.cm as cm import seaborn as sns %matplotlib inline sns.set_style('white') _, ax = plt.subplots(figsize=(8, 6)) xs = 250 ys = 0 def drawTree(x1, y1, angle, depth): if depth: x2 = x1 + math.cos(math.radians(angle)) * depth * 10.0 y2 = y1 + math.sin(math.radians(angle)) * depth * 10.0 ax.plot([x1, x2], [y1, y2], color=cm.hsv(depth/12.0)) drawTree(x2, y2, angle - 15, depth - 1) drawTree(x2, y2, angle + 15, depth - 1) drawTree(xs, ys, -90, 12)
今ハマっている春色とか.
PAGE_BREAK: PageBreak
import math import matplotlib.pyplot as plt import matplotlib.cm as cm import seaborn as sns %matplotlib inline sns.set_style('dark') _, ax = plt.subplots(figsize=(8, 6)) ax.invert_yaxis() ax.set_xticklabels([]) ax.set_yticklabels([]) xs = 250 ys = 0 def drawTree(x1, y1, angle, depth): if depth: x2 = x1 + math.cos(math.radians(angle)) * depth * 10.0 y2 = y1 + math.sin(math.radians(angle)) * depth * 10.0 ax.plot([x1, x2], [y1, y2], color=cm.spring(depth/12.0), linewidth=depth) drawTree(x2, y2, angle - 15, depth - 1) drawTree(x2, y2, angle + 15, depth - 1) drawTree(xs, ys, -90, 12)
まとめ
import math import numpy as np import matplotlib.pyplot as plt import matplotlib.cm as cm import seaborn as sns %matplotlib inline def mandelbrot(): sns.set_style('white') x, y = np.ogrid[-2.0:2.0:1000j,-2.0:2.0:1000j] c = x + 1j*y z = 0 for i in range(50): z = z**2 + c mask = np.abs(z) < 4 fig = plt.figure(1, figsize=(8, 6)) ax = fig.add_subplot(111) ax.imshow(mask.T,extent=[-2, 2, -2, 2]) plt.show() def henon_map(): sns.set_style('white') def fxy(a, b, x, y): return 1.0 - a * x**2 + y, b * x #pep8的に駄目な書き方 #fxy = lambda a, b, x, y : [1.0 - a * x**2 + y, b * x] a = 1.4 b = 0.3 xtemp = ytemp = 0 x = [xtemp] y = [ytemp] for i in range(10000): xtemp, ytemp = fxy(a, b, xtemp, ytemp) x.append(xtemp) y.append(ytemp) fig = plt.figure(2, figsize=(8, 6)) ax = fig.add_subplot(111) ax.scatter(x, y, color='c') plt.show() def logistic_map(): sns.set_style('white') fig = plt.figure(3, figsize=(8, 6)) ax = fig.add_subplot(111) def logistic(a): x = [0.8] for i in range(400): x.append(a * x[-1] * (1 - x[-1])) return x[-100:] for a in np.linspace(2.0, 4.0, 1000): x = logistic(a) ax.plot([a]*len(x), x, "c.", markersize=1.7) plt.show() def fractal_tree(): sns.set_style('dark') fig = plt.figure(4, figsize=(8, 6)) ax = fig.add_subplot(111) ax.invert_yaxis() ax.set_xticklabels([]) ax.set_yticklabels([]) xs = 250 ys = 0 def drawTree(x1, y1, angle, depth): if depth: x2 = x1 + math.cos(math.radians(angle)) * depth * 10.0 y2 = y1 + math.sin(math.radians(angle)) * depth * 10.0 ax.plot([x1, x2], [y1, y2], color=cm.spring(depth/12.0), linewidth=depth) drawTree(x2, y2, angle - 15, depth - 1) drawTree(x2, y2, angle + 15, depth - 1) drawTree(xs, ys, -90, 12) plt.show() def main(): mandelbrot() henon_map() logistic_map() fractal_tree() if __name__ == '__main__': main()
ピンバック: フラクタル・ツリー | 粉末@それは風のように (日記)
ピンバック: numpy + PIL | 粉末@それは風のように (日記)