アルゴリズム

[Python] マンデルブロ集合のフラクタル図形を描画する方法

2025-01-18更新日: 2025-01-18

マンデルブロ集合のフラクタル図形をPythonで描画するには、複素数の計算と反復処理を用います。

各点がマンデルブロ集合に属するかを判定するために、複素数 $c$ に対して反復式 $z_{n + 1} = z_{n}^{2} + c$ を使用し、一定回数の反復後に発散するかどうかを確認します。

発散しない点は集合に属し、発散する点は属しません。

描画には、matplotlibやPILなどのライブラリを使い、各点の発散速度に応じて色を付けることでフラクタル図形を生成します。

目次から探す

マンデルブロ集合とは
マンデルブロ集合の数学的背景
Pythonでマンデルブロ集合を描画する準備
マンデルブロ集合の描画手順
完全なサンプルコード
高速化のための工夫
応用例
まとめ

マンデルブロ集合とは

マンデルブロ集合は、複素数の特定の集合であり、フラクタル図形の一例です。

この集合は、反復関数 $z_{n + 1} = z_{n}^{2} + c$ において、初期値 $z_{0} = 0$ と複素数 $c$ を用いて生成されます。

マンデルブロ集合に属する点は、反復を行った際に発散しない複素数 $c$ の集合です。

視覚的には、非常に複雑で美しいパターンを形成し、数学や芸術の分野で広く注目されています。

特に、マンデルブロ集合の境界部分は無限に細かい構造を持ち、ズームインすることで新たなフラクタルパターンが現れることが特徴です。

マンデルブロ集合の数学的背景

反復関数 $z_{n + 1} = z_{n}^{2} + c$ の説明

マンデルブロ集合は、反復関数 $z_{n + 1} = z_{n}^{2} + c$ に基づいて定義されます。

ここで、 $z_{n}$ は複素数で、初期値は $z_{0} = 0$ です。

複素数 $c$ は、描画したい点を表します。

この関数を繰り返し適用することで、各点がマンデルブロ集合に属するかどうかを判断します。

具体的には、反復を続けた結果が無限大に発散するかどうかが重要です。

発散と収束の判定方法

マンデルブロ集合の点 $c$ が集合に属するかどうかは、反復の結果が発散するか収束するかによって決まります。

具体的には、次の条件を用います：

反復の過程で $| z_{n} |$ が 2 を超えた場合、その点 $c$ は発散すると判断します。
反復を一定回数(例えば、最大反復回数)行っても発散しなかった場合、その点は収束すると見なされ、マンデルブロ集合に属します。

発散の閾値と最大反復回数

発散の閾値は、一般的に $| z_{n} | > 2$ と設定されます。

この値を超えると、反復は無限大に向かうため、その点はマンデルブロ集合には含まれません。

また、最大反復回数は、計算の効率を考慮して設定されます。

例えば、1000回の反復を行い、その間に発散しなかった場合、その点はマンデルブロ集合に属すると判断します。

カラーリングの基準(発散速度に基づく色付け)

マンデルブロ集合の描画では、発散の速度に基づいて色付けを行います。

具体的には、各点の反復回数に応じて異なる色を割り当てることで、視覚的に美しいフラクタルを生成します。

例えば、早く発散した点には明るい色を、遅く発散した点には暗い色を使用することが一般的です。

この方法により、マンデルブロ集合の複雑な構造を視覚的に表現することができます。

Pythonでマンデルブロ集合を描画する準備

必要なライブラリのインストール

マンデルブロ集合を描画するためには、Pythonのいくつかのライブラリが必要です。

主に使用するのは、数値計算を行うための numpy と、グラフを描画するための matplotlib です。

これらのライブラリをインストールすることで、効率的にマンデルブロ集合を描画できます。

matplotlibのインストール

matplotlib は、Pythonでグラフや図を描画するためのライブラリです。

以下のコマンドを使用してインストールできます。

pip install matplotlib

numpyのインストール

numpy は、数値計算を効率的に行うためのライブラリです。

以下のコマンドを使用してインストールできます。

pip install numpy

描画領域の設定

マンデルブロ集合を描画する際には、描画領域を設定する必要があります。

描画領域は、複素数平面の特定の範囲を指定することで決まります。

通常、横軸(実数部)と縦軸(虚数部)の範囲を設定します。

例えば、以下のように設定することが一般的です：

実数部：-2.0 から 1.0
虚数部：-1.5 から 1.5

この範囲を設定することで、マンデルブロ集合の全体像を描画することができます。

複素数平面の範囲を決める

複素数平面の範囲を決める際には、描画する解像度も考慮する必要があります。

解像度が高いほど、より詳細なフラクタルを描画できますが、計算時間も増加します。

例えば、1000×1000ピクセルの解像度で描画する場合、以下のように範囲を設定します：

実数部の範囲： $[- 2.0, 1.0]$
虚数部の範囲： $[- 1.5, 1.5]$

この設定により、マンデルブロ集合の詳細な構造を視覚的に捉えることが可能になります。

マンデルブロ集合の描画手順

複素数平面の各点を計算する

マンデルブロ集合を描画するためには、複素数平面の各点に対して反復計算を行います。

まず、指定した範囲内の実数部と虚数部の値を用いて、複素数 $c$ を生成します。

以下のサンプルコードでは、1000×1000ピクセルの解像度で複素数平面の各点を計算します。

import numpy as np
# 描画領域の設定
x_min, x_max = -2.0, 1.0
y_min, y_max = -1.5, 1.5
width, height = 1000, 1000
# 複素数平面の各点を生成
x = np.linspace(x_min, x_max, width)
y = np.linspace(y_min, y_max, height)
C = np.array(np.meshgrid(x, y)).T.reshape(-1, 2)

発散判定の実装

次に、各点に対して反復計算を行い、発散するかどうかを判定します。

発散の閾値は $| z_{n} | > 2$ とし、最大反復回数を設定します。

以下のコードでは、発散判定を実装しています。

# 発散判定の設定
max_iter = 1000
Z = np.zeros(C.shape[0], dtype=np.complex128)
divergence = np.zeros(C.shape[0], dtype=int)
for i in range(max_iter):
    mask = np.abs(Z) <= 2  # 発散しない点をマスク
    Z[mask] = Z[mask]**2 + (C[mask, 0] + 1j * C[mask, 1])  # 反復計算
    divergence[mask] = i  # 発散しない点の反復回数を記録

色付けの実装

発散判定の結果を基に、各点に色を付けます。

発散の速度に応じて色を変えることで、視覚的に美しいフラクタルを生成します。

以下のコードでは、色付けを実装しています。

import matplotlib.pyplot as plt
# 色付けの実装
colors = plt.cm.viridis(divergence / max_iter)  # カラーマップを使用

matplotlibを使った描画

色付けが完了したら、matplotlib を使って描画します。

以下のコードでは、マンデルブロ集合を描画しています。

# 描画
plt.figure(figsize=(10, 10))
plt.imshow(colors.reshape(height, width, 4), extent=(x_min, x_max, y_min, y_max))
plt.title("Mandelbrot Set")
plt.xlabel("Re(c)")
plt.ylabel("Im(c)")
plt.axis("off")  # 軸を非表示にする
plt.show()

描画結果の保存

描画した結果を画像ファイルとして保存することもできます。

以下のコードでは、描画結果をPNG形式で保存しています。

# 描画結果の保存
plt.imsave("mandelbrot_set.png", colors.reshape(height, width, 4), extent=(x_min, x_max, y_min, y_max))

この手順を通じて、マンデルブロ集合を描画することができます。

各ステップを実行することで、複雑で美しいフラクタル図形を生成することが可能です。

完全なサンプルコード

以下に、マンデルブロ集合を描画するための完全なサンプルコードを示します。

このコードを実行することで、マンデルブロ集合のフラクタル図形を生成し、表示および保存することができます。

import numpy as np
import matplotlib.pyplot as plt
# 描画領域の設定
x_min, x_max = -2.0, 1.0
y_min, y_max = -1.5, 1.5
width, height = 1000, 1000
# 複素数平面の各点を生成
x = np.linspace(x_min, x_max, width)
y = np.linspace(y_min, y_max, height)
C = np.array(np.meshgrid(x, y)).T.reshape(-1, 2)
# 発散判定の設定
max_iter = 1000
Z = np.zeros(C.shape[0], dtype=np.complex128)
divergence = np.zeros(C.shape[0], dtype=int)
# 発散判定の実装
for i in range(max_iter):
    mask = np.abs(Z) <= 2  # 発散しない点をマスク
    Z[mask] = Z[mask]**2 + (C[mask, 0] + 1j * C[mask, 1])  # 反復計算
    divergence[mask] = i  # 発散しない点の反復回数を記録
# 色付けの実装
colors = plt.cm.viridis(divergence / max_iter)  # カラーマップを使用
# 描画
plt.figure(figsize=(10, 10))
plt.imshow(colors.reshape(height, width, 4), extent=(x_min, x_max, y_min, y_max))
plt.title("Mandelbrot Set")
plt.xlabel("Re(c)")
plt.ylabel("Im(c)")
plt.axis("off")  # 軸を非表示にする
plt.show()
# 描画結果の保存
plt.imsave("mandelbrot_set.png", colors.reshape(height, width, 4), extent=(x_min, x_max, y_min, y_max))

このコードを実行すると、マンデルブロ集合のフラクタル図形が表示され、同時に mandelbrot_set.png という名前で画像ファイルとして保存されます。

必要に応じて、描画領域や解像度、最大反復回数を変更することで、異なる結果を得ることができます。

高速化のための工夫

マンデルブロ集合の描画は計算量が多いため、処理を高速化するための工夫が重要です。

以下に、いくつかの方法を紹介します。

NumPyを使ったベクトル化処理

NumPyを使用することで、ループを避けてベクトル化処理を行うことができます。

これにより、計算速度が大幅に向上します。

例えば、複素数の反復計算を一度に行うことができるため、以下のように実装できます。

# 発散判定のベクトル化処理
Z = np.zeros(C.shape[0], dtype=np.complex128)
divergence = np.zeros(C.shape[0], dtype=int)
for i in range(max_iter):
    mask = np.abs(Z) <= 2  # 発散しない点をマスク
    Z[mask] = Z[mask]**2 + (C[mask, 0] + 1j * C[mask, 1])  # 反復計算
    divergence[mask] = i  # 発散しない点の反復回数を記録

このように、NumPyの配列操作を活用することで、Pythonの標準ループよりもはるかに高速に計算を行うことができます。

マルチスレッドやマルチプロセスの活用

計算をさらに高速化するために、マルチスレッドやマルチプロセスを利用することも有効です。

Pythonの concurrent.futures モジュールを使用することで、複数のスレッドやプロセスを使って計算を並行して行うことができます。

以下は、マルチプロセスを使用した例です。

from concurrent.futures import ProcessPoolExecutor
def mandelbrot_chunk(start, end):
    Z = np.zeros((end - start, ), dtype=np.complex128)
    divergence = np.zeros((end - start, ), dtype=int)
    # 反復計算のロジックをここに実装
    return divergence
# プロセスプールを作成
with ProcessPoolExecutor() as executor:
    results = list(executor.map(mandelbrot_chunk, range(0, height, chunk_size)))

この方法により、計算を複数のCPUコアに分散させることができ、全体の処理時間を短縮できます。

GPUを使った高速化(numbaやcupyの利用)

GPUを利用することで、さらに大規模な計算を高速化することが可能です。

numba や cupy といったライブラリを使用することで、GPU上での計算を簡単に実行できます。

numbaを使用した例

from numba import jit
@jit(nopython=True)
def mandelbrot_numba(C, max_iter):
    Z = np.zeros(C.shape[0], dtype=np.complex128)
    divergence = np.zeros(C.shape[0], dtype=int)
    for i in range(max_iter):
        mask = np.abs(Z) <= 2
        Z[mask] = Z[mask]**2 + C[mask]
        divergence[mask] = i
    return divergence

cupyを使用した例

import cupy as cp
# GPU上での計算
C_gpu = cp.asarray(C)
Z_gpu = cp.zeros(C_gpu.shape[0], dtype=cp.complex128)
divergence_gpu = cp.zeros(C_gpu.shape[0], dtype=int)
# 反復計算のロジックをGPU上で実行

これらの方法を活用することで、マンデルブロ集合の描画を大幅に高速化することができます。

特に、GPUを利用することで、数百万の点を同時に計算することが可能になり、描画速度が飛躍的に向上します。

応用例

マンデルブロ集合の描画技術を応用することで、さまざまな興味深いフラクタルを生成することができます。

以下にいくつかの応用例を紹介します。

ズームインして詳細なフラクタルを描画する

マンデルブロ集合の特性の一つは、無限にズームインすることで新たなパターンが現れることです。

特定の領域にズームインすることで、より詳細なフラクタルを描画できます。

以下のように、描画範囲を変更することでズームインを実現できます。

# ズームインする範囲を設定
x_min, x_max = -0.75, -0.7
y_min, y_max = 0.186, 0.188

このように範囲を狭めることで、より複雑なフラクタルの構造を観察することができます。

色のグラデーションをカスタマイズする

マンデルブロ集合の描画では、色のグラデーションをカスタマイズすることで、視覚的に魅力的な結果を得ることができます。

matplotlib のカラーマップを変更することで、異なる色合いを使用できます。

例えば、以下のようにカラーマップを変更することができます。

# 色のカスタマイズ
colors = plt.cm.plasma(divergence / max_iter)  # 'plasma'カラーマップを使用

このように、異なるカラーマップを試すことで、フラクタルの印象を大きく変えることができます。

アニメーションでフラクタルの変化を表現する

マンデルブロ集合の描画をアニメーション化することで、フラクタルの変化を視覚的に表現することができます。

matplotlib のアニメーション機能を使用して、ズームインや色の変化をアニメーション化することが可能です。

以下は、アニメーションの基本的な構造の例です。

from matplotlib.animation import FuncAnimation
# アニメーションの設定
fig, ax = plt.subplots()
def update(frame):
    # ズームインや色の変更を行うロジック
    ax.clear()
    ax.imshow(...)  # 新しいフラクタルを描画
    return ax,
ani = FuncAnimation(fig, update, frames=100, repeat=True)
plt.show()

このように、アニメーションを通じてフラクタルの動的な変化を楽しむことができます。

他のフラクタル(ジュリア集合など)を描画する

マンデルブロ集合の描画技術を応用して、他のフラクタルも描画できます。

例えば、ジュリア集合は、マンデルブロ集合と密接に関連しており、異なる複素数 $c$ を用いて生成されます。

ジュリア集合の描画は、以下のように実装できます。

def julia_set(c, max_iter, x_min, x_max, y_min, y_max, width, height):
    x = np.linspace(x_min, x_max, width)
    y = np.linspace(y_min, y_max, height)
    Z = np.array(np.meshgrid(x, y)).T.reshape(-1, 2)
    
    Z = Z[:, 0] + 1j * Z[:, 1]
    divergence = np.zeros(Z.shape[0], dtype=int)
    for i in range(max_iter):
        mask = np.abs(Z) <= 2
        Z[mask] = Z[mask]**2 + c
        divergence[mask] = i
    return divergence
# ジュリア集合の描画
c = complex(-0.7, 0.27015)  # ジュリア集合のパラメータ
divergence = julia_set(c, max_iter, x_min, x_max, y_min, y_max, width, height)