マンデルブロ集合を3次元化する試み-その1-

マンデルブロ集合を3次元化したい！と思って今Juliaでソースを書いてみてます。

今のところ、、惜しい！って感じです。

今できてる部分はこんな感じ。

複素数の拡張
マンデルブロ集合の計算
アフィン変換
描画

複素数の拡張

まずは、3次元化するために複素数を1次元拡張して、実部、虚部１，虚部２っていう風にします。

クォータニオン（四元数）をお手本にやってみました。

struct Wcomprex{T<:Real}<:Number
    real::T
    imag::T
    jmag::T
end

まずはWcomprexっていう新しい型を定義して

そして、2乗と足し算と絶対値の計算の関数を作る

function ww(w::Wcomprex)
    Wcomprex(w.real^2-w.imag^2-w.jmag^2,
    2*w.real*w.imag,
    2*w.real*w.jmag)
end
function wplus(w::Wcomprex,q::Wcomprex)
    Wcomprex(w.real+q.real,
    w.imag+q.imag,
    w.jmag+q.jmag)
end
function wabs(w::Wcomprex)
    sqrt(w.real^2+w.imag^2+w.jmag^2)
end

マンデルブロ集合の計算

マンデルブロ集合は$$z_{n+1}=z^2_n+c$$っていう漸化式を計算して、発散しないようなものを描画していくんだけどこれを拡張

$$w_{n+1}=w^2_n+d$$

$$w_n=x_n+iy_n+jz_n$$

$$d=a+ib+jc$$

これを計算します。

function checkmandel3d(d::Wcomprex)
    w = Wcomprex(0.0,0.0,0.0)
    for i=1:40
        w = wplus(ww(w),d)
        if wabs(w)>4
            return false
        end
    end
    return true
end

アフィン変換

3次元のものを描くときは回転とかも考えないといけないので

回転のためのアフィン変換を実装

function affineroll(x,y,z,xrad,yrad,zrad)
    inp = [x;y;z;1;;]
    xrol = [1 0 0 0;
            0 cos(xrad) -sin(xrad) 0;
            0 sin(xrad) cos(xrad) 0;
            0 0 0 1]
    yrol = [cos(yrad) 0 sin(yrad) 0;
            0 1 0 0;
            -sin(yrad) 0 cos(yrad) 0;
            0 0 0 1]
    zrol = [cos(zrad) -sin(zrad) 0 0;
            sin(zrad) cos(zrad) 0 0;
            0 0 1 0;
            0 0 0 1]

    xrol*yrol*zrol*inp

end

描画

最後に描画

@png begin
    
    background("black")
    dif = 0.005
    dotsize = 2.5
    for re=-2:dif:2
        for i=-2:dif:2
            for j=-2:dif:2
                if checkmandel3d(Wcomprex(re,i,j))
                    zahyo = affineroll(re,i,j,pi/3,pi/3,pi/6)
                    xpos = rescale(zahyo[1],-1.5,1.5,-1150,1150)
                    ypos = rescale(zahyo[2],-1.5,1.5,-1150,1150)

                    r = rescale(zahyo[1],-1,1,0.5,0.8)
                    g = rescale(zahyo[2],-1,1,0.7,0.9)
                    b = rescale(zahyo[3],-1,1,0.6,0.9)
                    a = rescale(zahyo[3],-1,1,0.5,0.8)

                    setcolor(r,g,b,a)
                    move(Point(xpos,ypos))
                    line(Point(xpos + dotsize, ypos))
                    line(Point(
                            xpos + dotsize,
                            ypos + dotsize,
                    ))
                     line(Point(xpos, ypos + dotsize))
                    fillpath()
                end
            end
        end
    end

end 2400 2400