Pythonでマンデルブロ集合を美しく描画する(基礎編)

# はじめに

マンデルブロ集合を描画することを主旨とした記事は多数あるので、これ自体に目新しさはないのだが、タイトルにある通り"美しく"描画することに関しては、色々と工夫の余地があるのでそれを記事にする。

"美しさ"といっても人それぞれではあるが、ここではWikiPediaのマンデルブロ集合の記事に掲載されている画像と同等のものを生成できることを目的にする。

https://ja.wikipedia.org/wiki/%E3%83%9E%E3%83%B3%E3%83%87%E3%83%AB%E3%83%96%E3%83%AD%E9%9B%86%E5%90%88

これはUltra Fractalというアプリで生成されたものなので、そちらの仕様に準拠するように作成することも意味する。

https://www.ultrafractal.com/

Pythonをそのまま使うと時間がかかりすぎるため、Numbaを使って高速化している。そのため静止画作成の描画時間はそこそこ短くなっている。CUDA環境にも対応してあるので、最後にCPU vs GPUでの描画時間等の性能評価も行なっている。

[Google Colabで動作するコード](https://colab.research.google.com/gist/T-STAR-LTD/b689665b7776bbbcaac75958e91bbffa/mandelbrot-basic.ipynb)


# 単純な実装

まずは、発散判定するまでのループ回数で視覚化するという、シンプルに最低限の実装をした場合の画像がこちら。

![Mandelbrot(Simple).png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/824193/6f8ea70f-cd1a-9c0a-a180-6ec01c39e018.png)

# 勾配をなめらかにする

上記のように普通に実装すると、ループ回数という整数値で明るさが決まるため段階的な変化がはっきり見える。これはこれで味があるかもしれないが、これを[なめらかにする手法](https://en.wikipedia.org/wiki/Plotting_algorithms_for_the_Mandelbrot_set#Continuous_(smooth)_coloring)が存在する。以下がそれを実装した場合の画像。この時、発散の条件(bailout)は理論値(4.0)より大きめにした方が勾配が綺麗になる。

![Mandelbrot(Smooth).png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/824193/13b6d4d7-f878-03a3-8d88-43a4d11787ec.png)

# 色をつける

WikiPediaの画像と大体同じになるように色をつける。色は三角関数を使って似たようなパターンになるように調整した。logを使ったり周期を0.35倍しているのもWikipediaの画像に合わせるため。

![Mandelbrot(Color).png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/824193/05e02c43-fc01-636e-1f2f-40e902e99a56.png)

# スーパーサンプリングを行う

画像のエッジ部分をなめらかにするために、[スーパーサンプリング](https://ja.wikipedia.org/wiki/%E3%82%B9%E3%83%BC%E3%83%91%E3%83%BC%E3%82%B5%E3%83%B3%E3%83%97%E3%83%AA%E3%83%B3%E3%82%B0)を行う。これはUltraFractalではAnti-aliasingをOnにした場合の処理に相当し、簡単に言えばピクセルの色を決定する際に位置を微妙にずらした計算結果を複数混ぜ合わせる。当然余計に時間がかかるようになるため、ここでは宣言時にparallel=Trueとして並列化させて、CPU数が多いほど早く実行できるようにした。

縦横倍で作成した画像を縮小して表示しても同じような結果になるが、[ShaderToy](https://www.shadertoy.com/results?query=mandelbrot)などGLSLを使ってリアルタイムでアニメーションさせる場合でも大抵この手法を使っている。
2x2の4点でサンプリングした場合の画像がこちら。

![Mandelbrot(SS).png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/824193/a4ada788-3c13-2ce5-50c5-d59d65fc5b0f.png)




[WikiPediaの記事内の拡大14](https://ja.wikipedia.org/wiki/%E3%83%9E%E3%83%B3%E3%83%87%E3%83%AB%E3%83%96%E3%83%AD%E9%9B%86%E5%90%88#/media/%E3%83%95%E3%82%A1%E3%82%A4%E3%83%AB:Mandel_zoom_14_satellite_julia_island.jpg)と同じ位置で作成した画像がこちら。色味が若干違うがほぼ同じ画像が生成できている。

![Mandelbrot.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/824193/56b2e793-04cb-0549-cd91-04b19c25cece.png)

# 性能評価

## 描画時間

上に掲載した画像(1600x1200で2x2サンプリング)の場合、Google ColabのCPUでは44秒程度かかる。GoogleColab環境のCPUは結構クロック数が低くコア数も少ないため、これは普通のPCのCPUで実行した方が速いだろう。筆者所有のM1 Mac Miniでは9秒程度だった。

CUDA対応したコードで、T4を使うと1.4秒程度にはなる。だいぶ短くなるがCUDAを使っているならもっと短くなっていてほしい気がする。Kaggle環境のP100を使うと、0.08秒という文字通り桁違いの速さになる。T4はコンシューマー向けのGPUと同様にFloat64性能が低めであるらしく、これが性能差の主因ではないかと思われる。

## 計算精度による限界

フラクタル画像は原理的には無限に拡大できるはずなのだが、Float64での計算が行われているので、計算精度によって限界がくる。以下は倍率を1E15にした際の出力画像。

![mandelbrot(1e15).png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/824193/2d614deb-0729-36ac-7315-97ce3ba2f3cc.png)

Youtube等にあるマンデルブロ集合のZoom系動画では1E1000越えも結構あって、それに比べると全く精度が足りない。倍率を上げるためにはより精度の高い演算が必要で、計算にものすごく時間がかかる。と思っていたのだが、実は[摂動論使ったアルゴリズム](https://en.wikipedia.org/wiki/Plotting_algorithms_for_the_Mandelbrot_set#Perturbation_theory_and_series_approximation)を用いることで、Float64の精度でも大幅に倍率を上げることが可能であり、おそらく最近の高精度に対応する描画アプリは全てこれを実装していると思われる。


# おわりに

まあまあ実用的な速度で綺麗なマンデルブロ集合を描画する実装ができた。
精度的な問題であまり深いズームができないが、これをベースに色付けを工夫したりすれば結構印象的な画像が生成できる。

## 関連記事

https://qiita.com/T-STAR/items/2ef76940f181acbc90f8

https://qiita.com/T-STAR/items/11229df2b6db462fb743