マンデルブロ集合

ここでは、書籍「フラクタルの美　複素力学系のイメージ」のp.191-192に解説されている「実験その２：マンデルブロー集合」を参考にして、初めてのマンデルブロ集合による作品を作ってみました。

Contents [非表示]

1 マンデルブロ集合

今回描いたマンデルブロ集合による作品は以下の図形となります。

よく見かける、マンデルブロ集合のかっこいい図形を描くことができました。

マンデルブロ集合の考え方

マンデルブロ集合とは複素平面上の点の集まりで、 $z_0=0$ として漸化式

$z_{n+1} \mapsto z_n^2 + c \ \ ( n = 0,1,2,\cdots )$

で定義される複素数列を考えたときに $n \to \infty$ の極限で無限大に発散しないという条件を満たす複素数 $c$ の全体が作る集合を指します。

このマンデルブロ集合をまともに計算することはできません。そのため、マンデルブロ集合を用いた図形を描く場合、

$|z_n|$ はある値 $M$ を越えたときに複素数列は発散したとみなす
$n$ が適当な上限値 $K$ に達したとき、 $z_K < M$ であれば複素数列は発散しなかったとみなす

として複素数 $c$ を分類し、発散したとみなされた場合はそのときの $n$ の値、発散しなかったとみなされた場合は上限値 $K$ に対応させて、複素数 $c$ の位置にある画素に色を与えるようにして描いています。

プログラムの流れ

マンデルブロ集合による作品を作る際のプログラムの流れを記載しておきます。

Step 0

キャンバスサイズ $w, h$ を決めます。
複素数 $c$ の範囲を $p_\mathrm{min} =-2.25, \ \ q_\mathrm{min} =-1.5, \ \ p_\mathrm{max} =0.75, \ \ q_\mathrm{max} =1.5$ 、判定条件 $M=200, \ \ K=256$ を決めます。
$\Delta p = (p_\mathrm{max} – p_\mathrm{min})/w, \ \ \Delta q = (q_\mathrm{max} – q_\mathrm{min})/h$ とおきます。

以下の処理はキャンバスのすべての画素 $(n_p, n_q)$ に対して行います。

Step 1

$p=p_\mathrm{min} + n_p \cdot \Delta p, \ \ q=q_\mathrm{min} + n_q \cdot \Delta q, \ \ k=0, \ \ x_0=y_0=0$ とおきます。

Step 2（反復ステップ）

$(x_k,y_k)$ をもとにして $(x_{k+1},y_{k+1})$ を漸化式

$x_{k+1} = x_k^2-y_k^2+p, \ \ y_{k+1}=2x_k y_k + q$

にしたがって計算します。このとき、カウンター $k$ の値に $1$ を加えます。

Step 3（判定ステップ）

$r = \sqrt{x_k^2+y_k^2}$ を計算します。

(i) $r>M$ ならば、 $k$ の値に応じてその画素の色を割り当てます。

(ii) $k=K$ ならば、その画素の色に $0$ （黒色）を割り当てます。

(iii) $r \leq M$ かつ $k<K$ ならば、Step 2に行きます。

プログラムコード

最後に、このマンデルブロ集合による作品を描いた際のプログラムコード（Processing）を示しておきます。

void setup(){
  size(1000,1000);
  
  PImage img = createImage(width, height, RGB);
    
  int Kmax = 256;
  float M = 200;
  float pmax = 0.75;
  float pmin = -2.25;
  float qmax = 1.5;
  float qmin = -1.5;

  float delta_p = (pmax - pmin) / width;
  float delta_q = (qmax - qmin) / height;
  
  float x, y, next_x, next_y, r;
  float p, q;
  int k;
  for(int i=0; i<width; i++){
    for(int j=0; j<height; j++){
      p = pmin + i * delta_p;
      q = qmax - j * delta_q;
      x = 0.0;
      y = 0.0;
      r = sqrt(x*x+y*y);
      k=0;
      img.set(i,j,color(0,0,0));
      while(k < Kmax){
        if(r > M){
          img.set(i,j,color(10*k,0,0));
          break;
        }
        next_x = x*x - y*y + p;
        next_y = 2.0 * x * y + q;
        x = next_x;
        y = next_y;
        r = sqrt(x*x+y*y);
        k++;
      }
    }
  }
  
  image(img,0,0,width,height);
}