マンデルブロ集合の見た目

マンデルブロ集合(Mandelbrot Set)は、このような図形で表されます。この図の黒い部分がマンデルブロ集合です。

Mandel zoom 00 mandelbrot set.jpg
By Created by Wolfgang Beyer with the program Ultra Fractal 3. - Own work, CC BY-SA 3.0, Link

この黒とそうでない領域の境界は複雑な形をしていて、どれだけ拡大しても複雑な図形が現れ、しかもそれは繰り返さず常に異なる図形が現れるという性質を持ちます。

Mandelbrot sequence new.gif
By Simpsons contributor at English Wikipedia - Transferred from en.wikipedia to Commons by Franklin.vp using CommonsHelper., Public Domain, Link

集合?

なぜこの図形が「集合」なのかというと、数の集まりを表現したものだからです。整数の集合、自然数の集合、実数の集合などは、一次元の数直線上に描くことができます。

  • 整数の集合 $ \mathbb{Z} = { …, -3, -2, -1, 0, 1, 2, 3, … } $
  • 自然数の集合 $ \mathbb{N} = { 1, 2, 3, … } $
  • 実数の集合 $ \mathbb{R} = { \sqrt{2}, e, \pi, \log{2}, 0, 1.5, … } $
  • 偶数の集合 $ = { 0, 2, 4, 8, … } $

マンデルブロ集合は、複素数の集合です。複素数は、虚数単位i を用いて、$ x + yi $ の形で表現できる数です。そして複素数は、\(x\)軸を実数、$y$軸を虚数にとった二次元の複素平面上で描画できます。

Illustration of a complex number.svg
By Никита Воробьев - Own work, Public Domain, Link

偶数の集合は、2で割って余りのない整数の集まりとして定義できます。では、マンデルブロ集合はどのような複素数の集まりなのでしょうか。

複素数の数列 $z_n$ を以下のように定義するとき、$ n \rightarrow \infty $ で発散しなければ複素数$c$ はマンデルブロ集合に含まれます。

\[\begin{eqnarray} \left\{ \begin{array}{l} z_{n+1} = z_n^2 + c \\ z_0 = 0 \end{array} \right. \end{eqnarray}\]

具体的に計算してみます。例えば $c = 0$ はマンデルブロ集合に含まれるでしょうか?

\[\begin{eqnarray} \begin{array}{l} z_0 = 0 \\ z_1 = {z_0}^2 + c = 0^2 + 0 = 0 \\ z_2 = {z_1}^2 + c = 0^2 + 0 = 0 \\ z_3 = {z_2}^2 + c = 0^2 + 0 = 0 \\ \end{array} \end{eqnarray}\]

$z_3$ まで計算して、これは $n$ が大きくなっても 0 のままであろうということがわかりました。$n \rightarrow \infty$ で発散しないので、$c = 0$ はマンデルブロ集合に含まれます。

次に、$c = 1 + i$ はマンデルブロ集合に含まれるでしょうか?

\[\begin{eqnarray} \begin{array}{l} z_0 = 0 \\ z_1 = {z_0}^2 + c = 0 + 1 + i = 1 + i \\ z_2 = {z_1}^2 + c = (1 + i)^2 + 1 + i = 1 + 2i - 1 + 1 + i = 1 + 3i \\ z_3 = {z_2}^2 + c = (1 + 3i)^2 + 1 + i = 1 + 6i - 9 + 1 + i = -7 + 7i \\ z_4 = {z_3}^2 + c = (-7 + 7i)^2 + 1 + i = 49 -98i - 49 + 1 + i = 1 + 99i \\ z_5 = {z_4}^2 + c = (1 + 99i)^2 + 1 + i = 1 + 198i - 9801 + 1 + i = -9799 + 199i \\ \end{array} \end{eqnarray}\]

これは発散しそうです。従って $c = 1 + i$ はマンデルブロ集合に含まれません。図を見ることでもわかります。

Mandelset hires.png
By Connelly (talk· contribs), Public Domain, Link

この図は、x軸が実数軸(Re: Real)、y軸が虚数軸(Im: Imaginary)になっています。$z = 1 + i$ は、$Re = 1$、$Im = 1$ が対応しますが、この点は真っ白であることから、黒で塗りつぶされたマンデルブロ集合に含まれないことがわかります。

JavaScript で描画する

まず複素数 $c = cRe + cIm \cdot i$ を生成する処理です。縦横pixelピクセルの画像を、x = 0 から pixel まで、$y = 0$ から pixel まで二重の for文で回して、各$(x, y)$ をスケールや描画範囲を加味して複素数 $c = cRe + cIm \cdot i$ に変換しています。

for (let x = 0; x < this.pixel; x++) {
    for (let y = 0; y < this.pixel; y++) {
        let cRe = x * this.scale / this.pixel - this.scale / 2 + this.centerX
        let cIm = y * this.scale / this.pixel - this.scale / 2 + this.centerY
        this.drawPoint(x, y, this.color(cRe, cIm))
    }
}

次に、各点の色を計算する処理です。与えられた複素数 $(cRe, cIm)$ が発散する場合は黒を返し、それ以外の場合は、発散の度合いによって色を変えています。マンデルブロ集合の描画サンプルが綺麗なのは、後者の見せ方の工夫によるところが大きく、今回は比較的地味な描画になっています。

color(cRe, cIm) {
    let re = 0
    let im = 0
    for (let n = 0; n < 50; n++) {
        let reNew = re * re - im * im + cRe
        let imNew = 2 * re * im + cIm
        re = reNew
        im = imNew
        let distance = re * re + im * im
        if (distance > 4) {
            return [n / 50 * 255, 255 / (1 + Math.exp(-Math.log(distance))), 0] // diverge
        }
    }
    return [0, 0, 0] // not diverge
}

マンデルブロ集合の定義としては、$n \rightarrow \infty$ の時に発散するか否かでしたが、$c = 0$ と $c = 1 + i$ の計算例に見るように、比較的小さな $n$ でも発散することが多いので、ここでは 50回までの計算にとどめています。そして、発散するかどうかも、複素平面上で原点からの距離の二乗が4を超えるようであれば、そこで計算を打ち切っています。4 を超える場合、どのみち、原点近くまで戻ってくることは、もうないでしょう。

最後は実際の描画処理です。Canvas の context に画像ピクセルを書き込んでいます。

drawPoint(x, y, color) {
    let image = this.ctx.createImageData(1, 1)
    image.data[0] = color[0]
    image.data[1] = color[1]
    image.data[2] = color[2]
    image.data[3] = 255
    this.ctx.putImageData(image, x, y)
}

描画に9秒かかりました。https://github.com/cslarsen/mandelbrot-j では、より高速な描画を実現しているので、パフォーマンスを上げる余地はありそうです。

Rust で描画する

次に、Rust + WebAssembly で描画してみます。ロジックは JS版と同じです。

まずは $(x, y)$ を for で回すところ。

for x in 0..pixel {
    for y in 0..pixel {
        let c_re: f64 = (x as f64)*scale/(pixel as f64) - scale/2.0 + center_x;
        let c_im = (y as f64)*scale/(pixel as f64) - scale/2.0 + center_y;
        let color = color(c_re, c_im);
        let pixel = image.get_pixel_mut(x, y);
        *pixel = image::Rgb(color);
    }
}

let mut buffer = Vec::new();
let mut writer = Cursor::new(&mut buffer);
let _ret = image.write_to(&mut writer, ImageOutputFormat::Png);
return buffer;

次に色を出し分ける、つまり収束判定をする部分。

pub fn color(c_re: f64, c_im: f64) -> [u8; 3] {
    let mut n = 0;
    let mut re: f64 = 0.0;
    let mut im: f64 = 0.0;
    while n < 50 {
        let new_re = re*re - im*im + c_re;
        let new_im = 2.0*re*im + c_im;
        re = new_re;
        im = new_im;
        let distance = re*re + im*im;
        if distance > 4.0 {
            return [
                (((n as f64)/50.0) * 255.0) as u8,
                (255.0 / (1.0+((-distance.ln()).exp()))) as u8,
                0
            ];
        }
        n += 1;
    }
    return [0, 0, 0];
}

16ms で描画できました。このように、Pure JS と Rust では、マンデルブロ集合の描画に大きなパフォーマンスの差が出ることがわかります。

コード

今回作成したコードはこちら。