概要
先日、VRChatで以下のようなワールドを作成しPublic化しました。シェーダのソースコードを解説します。
「Fractal Room」Public化しました！ こんな感じで、拡大・縮小を続ける巨大なマンデルブロ集合を眺められるワールドです。全部シェーダで書いてるので、ワールドの容量は正真正銘の0MB！遊びに来てね～ #VRC #VRChat_world紹介 pic.twitter.com/DGL4FSsmIM
— Hibit (@hibit_at) 2018年9月22日
一部引用の明記
マンデルブロ集合の計算については、この記事を参考にさせていただきました（@notargsさん、ありがとうございます）。本コードではこれの一部を改変して、
- 時間変化によって任意の座標を中心に拡大・縮小する
- 発散速度によってカラフルなグラデーションをつける
という要素を付け加えています。
ソースコード
とりあえずソースを見せろという方へ。
Shader "Custom/Mandelbrot"
{
Properties{
_MaxIteration("MaxIteration", Range(1, 2048)) = 512
_Threshold("Threshold", Range(1, 100)) = 2
_PosX("PosX",Range(-2,2)) = -0.685
_PosY("PosY",Range(-2,2)) = -0.3
}
SubShader
{
Tags { "RenderType" = "Opaque" }
LOD 100
Pass
{
CGPROGRAM
#pragma vertex vert
#pragma fragment frag
float _Threshold;
int _MaxIteration;
float _PosX;
float _PosY;
float2 cpow(float2 v)
{
float a = 3.14;
if (v.x != 0){
a = atan2(v.y, v.x) * 2;
}
return float2(cos(a), sin(a)) * pow(length(v), 2);
}
half mandelbrot(half2 c)
{
half2 z = half2(0, 0);
for (int i = 0; i < _MaxIteration; i++)
{
z = cpow(z);
z += c;
if (length(z) > _Threshold){
return (float)i / _MaxIteration;
}
}
return 1.0;
}
struct appdata
{
float4 vertex : POSITION;
float2 uv : TEXCOORD0;
};
struct v2f
{
float2 uv : TEXCOORD0;
float4 vertex : SV_POSITION;
};
v2f vert (appdata v)
{
v2f o;
o.vertex = UnityObjectToClipPos(v.vertex);
o.uv = v.uv;
return o;
}
fixed4 frag (v2f i) : SV_Target
{
float scale = 2.1+pow(sin(_Time.x),5)*2;
scale = pow(10,scale)-1;
i.uv.x = i.uv.x + _PosX * scale;
i.uv.y = i.uv.y + _PosY * scale;
fixed col = mandelbrot((i.uv - 0.5)/scale);
fixed col_r = 0.5+sin(col*6.28)/2;
fixed col_g = 0.5+sin(col*6.28+2.09)/2;
fixed col_b = 0.5+sin(col*6.28+4.19)/2;
if (col == 1){
col_r = 0;
col_g = 0;
col_b = 0;
}
return fixed4(col_r,col_g,col_b,0);
}
ENDCG
}
}
Fallback "Diffuse"
}
複素数の2乗を定義
マンデルブロ集合の計算には、複素数の2乗を計算する必要がありますが、その計算関数を自前で用意します。ある複素数$z = x + yi$が複素数平面の極座標で$z(r,\theta)$として表せるとき、
z(r,\theta)^2 = z(r^2,2\theta)
が成り立ちます。$r$は$\sqrt{x^2+y^2}^{※1}$、$\theta$は$tan^{-1}(\frac{y}{x})^{※2}$で求められるので、コードとしては以下のようになります。
float2 cpow(float2 v)
{
float a = 3.14;
if (v.x != 0){
a = atan2(v.y, v.x) * 2;
}
return float2(cos(a), sin(a)) * pow(length(v), 2);
}
※1 … この計算はlength関数としてデフォルトで用意されています。
※2 … x = 0の場合は0除算でエラーになってしまうのでif文で個別処理します。
マンデルブロ集合を計算
ある複素数$c$について、
z_0 = 0\\
z_{n+1} = z_n^2 + c
という計算をn回おこなった時、$z_n$の長さが一定以上になる時を「発散」とします。規定の繰り返しまでに「発散」すればその回数の逆数（≒ 1をMAXとした発散しづらさ）を、しなければあきらめて1を返すような関数Mandelbrotを定義します。これによってマンデルブロ集合が描かれます。何回計算して（MaxIteration）何を閾値とするか（Threshold）はプロパティで調整しますが、それぞれ512、2ぐらいがちょうど良いかなと思います。
half mandelbrot(half2 c)
{
half2 z = half2(0, 0);
for (int i = 0; i < _MaxIteration; i++)
{
z = cpow(z);
z += c;
if (length(z) > _Threshold){
return (float)i / _MaxIteration;
}
}
return 1.0;
}
時間によって拡大・縮小 + 色付け
fixed4 frag (v2f i) : SV_Target
{
float scale = 2.1+pow(sin(_Time.x),5)*2;
scale = pow(10,scale)-1;
i.uv.x = i.uv.x + _PosX * scale;
i.uv.y = i.uv.y + _PosY * scale;
fixed col = mandelbrot((i.uv - 0.5)/scale);
fixed col_r = 0.5+sin(col*6.28)/2;
fixed col_g = 0.5+sin(col*6.28+2.09)/2;
fixed col_b = 0.5+sin(col*6.28+4.19)/2;
if (col == 1){
col_r = 0;
col_g = 0;
col_b = 0;
}
return fixed4(col_r,col_g,col_b,0);
}
この部分で
- 時間変化による、任意の座標を中心とした拡大・縮小
- 発散速度を元にしたグラデーション色付け
を行っています。
拡大・縮小
float scale = 2.1+pow(sin(_Time.y),5)*2;
scale = pow(10,scale)-1;
i.uv.x = i.uv.x + _PosX * scale;
i.uv.y = i.uv.y + _PosY * scale;
fixed col = mandelbrot((i.uv - 0.5)/scale);
倍率（scale）は、1倍から10,000倍ぐらいまでが周期変化するようにしたいのですが、線形的に変化させると前半が速くなりすぎる&後半が遅くなりすぎるので、指数的に変化するようにします。$10^n$の右肩にある指数部を、三角関数を用いて0.1～4.1で周期変化するようにしています。なぜ+0.1しているのかについては、0にすると$10^0-1 = 0$になってエラーが起こるからです。
拡大・縮小の中心となる座標については、プロパティ（PosXとPosY)で調整します。縮小した時はともかく、拡大した時は0.01動くだけですごく動いてしまうので、ミリ単位での調整が必要になってきます。私はPosX = -0.685,PosY = -0.3がいい感じなのを見つけたのでそれにしてますが、各自エモい組み合わせを見つけてください。
グラデーション
fixed col_r = 0.5+sin(col*6.28)/2;
fixed col_g = 0.5+sin(col*6.28+2.09)/2;
fixed col_b = 0.5+sin(col*6.28+4.19)/2;
if (col == 1){
col_r = 0;
col_g = 0;
col_b = 0;
}
return fixed4(col_r,col_g,col_b,0);
先程、colという変数に各座標におけるMandelbrot関数を適用した後の値をいれました。これはその座標の「発散しづらさ」を示しています。0だとすぐに発散して、何回計算しても発散しなかったら1になります。
この0～1の間の数値にうまく色をあてはめてやればグラデーションができます。このコードではRGBの3色について三角関数を(2/3)πずつずらして足し合わせることによりそれを表現しています。また、発散しなかった場合はその領域をハッキリと明示したいので（巷でよくみるマンデルブロ集合の絵もよくそうなっています）そこだけはif文の個別処理で黒色にしています。
ここの味付けも各自調節して、最も美しいと感じる組み合わせを見つけるのもいいと思います。本記事ではソースコードまで貼りませんが、グラデーション用のテクスチャ画像して、それを読み込むようにすればより自由にグラデーションを定義できます。