processingの備忘録　-マンデルブロ集合- - プログラミングの備忘録

こんにちは。

今回は「マンデルブロ集合を描く」という題です。

マンデルブロ集合とは
プログラム化する
まとめ
参考

マンデルブロ集合とは

マンデルブロ集合は、

$$ z_{n+1} = z_n^{2} + c $$

$$ z_0 = 0 $$

で表される複素数の数列で、$n \to \infty$にしたときに発散しないような複素数$c$の集合のことだそうです。

具体的には、ある複素数$c = a+bi$について先に挙げた式にあてはめて計算していき、

$$ \begin{align} z_0 & = 0 \\ z_1 & = c \\ & = a+bi \\ z_2 & = c^{2} + c \\ & = a^{2} - b^{2} + a + \left( 2ab + b \right)i \\ & \vdots \end{align} $$

$n \to \infty$、つまり$z_\infty$のときの値が無限大にならないようなとき、そのときの複素数$c$がマンデルブロ集合の要素ということになります。

そして、複素数$c$を$a+bi$としたときに、座標上の(a, b)の位置に点を打っていくとマンデルブロ集合の図になります。

プログラム化する

マンデルブロ集合

やっていることは$z_n$の値を2乗して$c$を足すことで次の$z_{n+1}$を計算する、ということを繰り返しているだけなので、プログラム上は、

z = z^2 + c

のような形になります。

さて、マンデルブロ集合を描く上では複素数の扱いが必須になりますが、わざわざ複素数のライブラリやなんかを使わなくてもできます。

実部と虚部に分けて計算すれば良いわけです。
つまり先のマンデルブロ集合の式は、複素数$z_n$を$c_n + d_n i$とすると、

$$ \begin{align} c_n {}_+ {}_1 & = c_n^{2} - d_n^{2} + a \\ d_n {}_+ {}_1 & = 2 c_n d_n + b \end{align} $$

$$ \begin{align} c_0 & = 0 \\ d_0 & = 0 \end{align} $$

となります。

ゴチャゴチャ書いても分かりづらいので、プログラムを書いてみます。

とりあえずは、発散しない部分を黒く塗るような描き方でプログラム化します。

size(500, 500);
background(255);

loadPixels(); //画面の読み込み
for(int i = 0; i < width; i ++){ //画面の全範囲に対して
for(int j = 0; j < height; j ++){

    //描画領域（x:-2~2、y:-2~2）に合わせる
    float x = map(i, 0, width, -2, 2);
    float y = map(j, 0, height, -2, 2);
    
    //複素数cの実部･虚部は座標上の位置になる
    float a = x; 
    float b = y;
  
    int n = 1;
    int n_max = 100; //n -> ∞の「∞」を設定

    int inf = 100; //zの値がどこまでいったら発散したと判断するかを設定
    
    while(n < n_max){ //n_maxまでの間

        //z^2を計算
        float xx = x*x-y*y;
        float yy = 2*x*y;
      
        //z^2にcを足してzを更新
        x = xx + a;
        y = yy + b;
  
        if(sqrt(x*x+y*y) > inf){ //もしzの絶対値が無限大になった（発散した）ら
            break; //ループから抜ける
        }
        
        n++; //nを増やす
        
        if(n == n_max){ //計算回数が最大になった（発散しなかった）ら
            int pos = i + j * width; //画面上の位置を計算
            pixels[pos] = color(0); //色を黒くする
        }
    }
} 
}
updatePixels(); //画面を更新

実行すると、

よく見るマンデルブロ集合の図が描けました。

先程のプログラムでは、nは1~100まで計算するという設定にしていましたが、発散したら全て白く塗っていました。
ですが、もしかしたらn=2で発散しているかもしれないし、n=99まで発散しなかったかもしれません。

ということで、発散する速度によって色を変えてみます。

if(n == n_max){
    int pos = i + j * width;
    pixels[pos] = color(0);
}

の部分を、

int pos = i + j * width;
for(int k = 1; k < n_max; k ++){
    if(n == k){ //発散したときのnによって場合分け
        colorMode(HSB, 360, 100, 100); //カラーモードの変更
        float c = map(k, 1, n_max, 250, 0); //kの値を1~n_maxの範囲から250~0の範囲へ変換
        pixels[pos] = color(c, 80, 80); //nに対応した色に塗る
    }
    if(n == n_max){ //発散しなかったら
        pixels[pos] = color(0, 0, 0); //黒に塗る
    }
}

に変えます。

実行すると、

カラフルになりました。

このコードではカラーモードを変えて、1つの変数の違いによってグラデーション的に色を変えられるようにしています。
カラーモードはデフォルトではRGBですが、HSBというモードに変更しています。
HSBはHue (色相)、Saturation (彩度)、Brightness (明度) によって色を変えるモードで、「色相」は「色相環」における位置を指し、これによって赤から青まで色を変えられます。

よって、nの1~100を色相の250~0（青から赤）に対応させて色を変化させているということになります。