Hi Dave,
The plugin is not is the database yet, and I made some small changes to 3DFractalRaytrace which I am testing along with the new plugin. They should be in the database sometime tomorrow. The biggest challenge was that my normal distance estimate method can't be used as I can't seem to get the right derivative. I have an alternate method which doesn't require the derivative, but it is a little tricky to use. My brute force method also doesn't seem to behave well.
Hi Ron, here's my code for the relevant formulas:
<snip>
elseif @fractaltype==1 ; Truly 3D Mandelbrot formulas
zri = @startri
zj = @startj
cri = x1 + flip(y1)
cj = z1
magn = |@startri| + sqr(@startj)
if @useDE && (@fractal3D<4 || @fractal3D>7)
dzri = (1,0)
dzj = 0.0
endif
elseif @fractaltype==2 ; Truly 3D Julia Formulas
zri = x1 + flip(y1)
zj = z1
cri = @constri
cj = @constj
magn = |zri| + sqr(zj)
if @useDE && (@fractal3D<4 || @fractal3D>7)
dzri = (1,0)
dzj = 0.0
endif
<snip>
elseif @fractaltype<3 ; Truly 3D
if (@fractal3D<4 || @fractal3D>5) && magn<1e-100
zri = cri
zj = cj
elseif @fractal3D==0 ; White/Nylander
if @useDE
magn = sqrt(magn)
r = sqrt(|dzri| + sqr(dzj))
th = atan2(zri)
ph = asin(zj/magn)
dr = @mpwr*r*magn^(@mpwr-1.0)
dph = (@mpwr-1.0)*ph + asin(dzj/r)
dth = (@mpwr-1.0)*th + atan2(dzri)
dzj = -dr*sin(dph)
dzri = dr*cos(dph)*(cos(dth) + flip(sin(dth)))
if @fractaltype==1
dzri = dzri + 1.0
endif
r = magn^@mpwr
dph = @mpwr*ph
dth = @mpwr*th
zj = -r*sin(dph) + cj
zri = r*cos(dph)*(cos(dth) + flip(sin(dth))) + cri
else
r = (magn=sqrt(magn))^@mpwr
ph = @mpwr*asin(zj/magn)
th = @mpwr*atan2(zri)
zj = -r*sin(ph) + cj
zri = r*cos(ph)*(cos(th) + flip(sin(th))) + cri
endif
elseif @fractal3D==1 ; Thornton 1
if @useDE
magn = sqrt(magn)
r = sqrt(|dzri| + sqr(dzj))
th = atan2(zri)
ph = acos(zj/magn)
dr = @mpwr*r*magn^(@mpwr-1.0)
dph = (@mpwr-1.0)*ph + acos(dzj/r)
dth = (@mpwr-1.0)*th + atan2(dzri)
dzj = dr*cos(dph)
dzri = dr*sin(dph)*(cos(dth) + flip(sin(dth)))
if @fractaltype==1
dzri = dzri + 1.0
endif
r = magn^@mpwr
dph = @mpwr*ph
dth = @mpwr*th
zj = r*cos(dph) + cj
zri = r*sin(dph)*(cos(dth) + flip(sin(dth))) + cri
else
r = (magn=sqrt(magn))^@mpwr
ph = @mpwr*acos(zj/magn)
th = @mpwr*atan2(zri)
zj = r*cos(ph) + cj
zri = r*sin(ph)*(cos(th) + flip(sin(th))) + cri
endif
elseif @fractal3D==2 ; Thornton 2
if @useDE
magn = sqrt(magn)
r = sqrt(|dzri| + sqr(dzj))
th = atan2(zri)
if (ph = acos(zj/magn))>0.5*#pi
ph = #pi - ph
endif
if (dph = acos(dzj/r))>0.5*#pi
dph = #pi - dph
endif
dr = @mpwr*r*magn^(@mpwr-1.0)
dph = (@mpwr-1.0)*ph + dph
dth = (@mpwr-1.0)*th + atan2(dzri)
dzj = dr*cos(dph)
dzri = dr*sin(dph)*(cos(dth) + flip(sin(dth)))
if @fractaltype==1
dzri = dzri + 1.0
endif
r = magn^@mpwr
dth = @mpwr*th
dph = @mpwr*ph
zri = r*sin(dph)*(cos(dth) + flip(sin(dth))) + cri
zj = r*cos(dph) + cj
else
r = (magn=sqrt(magn))^@mpwr
if (ph = acos(zj/magn))>0.5*#pi
ph = #pi - ph
endif
th = @mpwr*atan2(zri)
ph = @mpwr*ph
zj = r*cos(ph) + cj
zri = r*sin(ph)*(cos(th) + flip(sin(th))) + cri
endif
elseif @fractal3D==3 ; Rucker
if @useDE
magn = sqrt(magn)
r = sqrt(|dzri| + sqr(dzj))
th = atan2(zri)
ph = atan2(real(zri)+flip(zj))
dr = @mpwr*r*magn^(@mpwr-1.0)
dph = (@mpwr-1.0)*ph + atan2(real(dzri)+flip(dzj))
dth = (@mpwr-1.0)*th + atan2(dzri)
dzj = dr*sin(dph)
dzri = dr*cos(dph)*(cos(dth) + flip(sin(dth)))
if @fractaltype==1
dzri = dzri + 1.0
endif
r = magn^@mpwr
dph = @mpwr*ph
dth = @mpwr*th
zj = r*sin(dph) + cj
zri = r*cos(dph)*(cos(dth) + flip(sin(dth))) + cri
else
r = (magn=sqrt(magn))^@mpwr
ph = @mpwr*atan2(real(zri)+flip(zj))
th = @mpwr*atan2(zri)
zj = r*sin(ph) + cj
zri = r*cos(ph)*(cos(th) + flip(sin(th))) + cri
endif
elseif @fractal3D==4 ; Makin 3D
r = zj
zj = 2.0*zj*(real(zri) - imag(zri)) + cj
zri = zri*zri - r*r + cri
; r = zj
; zj = zj*(2.0*(real(zri) + imag(zri)) - zj) + cj
; zri = real(zri)*(2.0*(imag(zri) + r) - real(zri)) \
; + flip(imag(zri)*(2.0*(real(zri) + r) - imag(zri))) + cri
; * | r i j
; ----------------
; r | r i j
; i | i -r -j
; j | j -j -r
elseif @fractal3D==5 ; Makin 3D 2
r = zj
zj = -zj*zj + 2.0*real(zri)*imag(zri) + cj
zri = real(zri)*real(zri) - flip(imag(zri)*imag(zri)) \
+ 2.0*r*(imag(zri) - flip(real(zri))) + cri
; * | r i j
; ----------------
; r | r j -i
; i | j -i r
; j | -i r -j
elseif @fractal3D==6 ; Makin 3D 3
ztemp = zri
zri = sqr(zri) - zj*zj
zri = real(zri) + flip(imag(zri)*(1.0 - sqr(zj)/|ztemp|)) + cri
zj = 2.0*real(ztemp)*zj \
*(1.0 - sqr(imag(ztemp))/|real(ztemp)+flip(zj)|) + cj
elseif @fractal3D==7 ; Makin 3D 4
ztemp = zri
zri = sqr(zri) - zj*zj
zri = real(zri) + flip(imag(zri)*(1.0 + sqr(zj)/|ztemp|)) + cri
zj = 2.0*real(ztemp)*zj \
*(1.0 + sqr(imag(ztemp))/|real(ztemp)+flip(zj)|) + cj
elseif @fractal3D==8 ; White
;r = sqrt(x^2 + y^2 + z^2 )
;theta = atan2(sqrt(x^2 + y^2) , z)
;phi = atan2(y,x) ;
;
;newx = r^n * sin(theta*n) * cos(phi*n)
;newy = r^n * sin(theta*n) * sin(phi*n)
;newz = r^n * cos(theta*n)
if @useDE
magn = sqrt(magn)
r = sqrt(|dzri| + sqr(dzj))
ph = atan2(zri)
th = asin(zj/magn)
dr = @mpwr*r*magn^(@mpwr-1.0)
dth = (@mpwr-1.0)*th + asin(dzj/r)
dph = (@mpwr-1.0)*ph + atan2(dzri)
dzj = dr*cos(dth)
dzri = dr*sin(dth)*(cos(dph) + flip(sin(dph)))
if @fractaltype==1
dzri = dzri + 1.0
endif
r = magn^@mpwr
dth = @mpwr*th
dph = @mpwr*ph
zj = r*cos(dth) + cj
zri = r*sin(dth)*(cos(dph) + flip(sin(dph))) + cri
else
r = (magn=sqrt(magn))^@mpwr
th = @mpwr*asin(zj/magn)
ph = @mpwr*atan2(zri)
zj = r*cos(th) + cj
zri = r*sin(th)*(cos(ph) + flip(sin(ph))) + cri
endif
elseif @fractal3D==9 ; White/Nylander quick z^2+c
if @useDE
th = atan2(zri)
ph = asin(zj/(magn=sqrt(magn)))
dr = 2.0*(r=sqrt(|dzri|+sqr(dzj)))*magn
dth = th + atan2(dzri)
dph = ph + asin(dzj/r)
dzri = dr*(cos(dph)*cos(dth) + flip(cos(dph)*sin(dth)))
dzj = -dr*sin(dph)
if @fractaltype==1
dzri = dzri + 1.0
endif
endif
zri = sqr(zri)*(1.0 - sqr(zj)/(r=|zri|)) + cri
zj = -2.0*zj*sqrt(r) + cj
elseif @fractal3D==10 ; White/Nylander quick z^3+c
if @useDE
th = atan2(zri)
ph = asin(zj/sqrt(magn))
dr = 3.0*(r=sqrt(|dzri|+sqr(dzj)))*magn
dth = 2.0*th + atan2(dzri)
dph = 2.0*ph + asin(dzj/r)
dzri = dr*(cos(dph)*cos(dth) + flip(cos(dph)*sin(dth)))
dzj = -dr*sin(dph)
if @fractaltype==1
dzri = dzri + 1.0
endif
endif
zri = (zri^3)*(1.0 - 3.0*sqr(zj)/(r=|zri|)) + cri
zj = zj^3 - 3.0*zj*r + cj
elseif @fractal3D==11 ; White/Nylander quick z^4+c
if @useDE
th = atan2(zri)
ph = asin(zj/(magn=sqrt(magn)))
dr = 4.0*(r=sqrt(|dzri|+sqr(dzj)))*magn^3
dth = 3.0*th + atan2(dzri)
dph = 3.0*ph + asin(dzj/r)
dzri = dr*(cos(dph)*cos(dth) + flip(cos(dph)*sin(dth)))
dzj = -dr*sin(dph)
if @fractaltype==1
dzri = dzri + 1.0
endif
endif
r = |zri|
zri = (zri^4)*(1.0 - (6.0*r*sqr(zj) - zj^4)/sqr(r)) + cri
zj = -4.0*zj*sqrt(r)*(r - sqr(zj)) + cj
elseif @fractal3D==12 ; White/Nylander quick z^5+c
if @useDE
th = atan2(zri)
ph = asin(zj/(magn=sqrt(magn)))
dr = 5.0*(r=sqrt(|dzri|+sqr(dzj)))*magn^4
dth = 4.0*th + atan2(dzri)
dph = 4.0*ph + asin(dzj/r)
dzri = dr*(cos(dph)*cos(dth) + flip(cos(dph)*sin(dth)))
dzj = -dr*sin(dph)
if @fractaltype==1
dzri = dzri + 1.0
endif
endif
r = |zri|
zri = (zri^5)*(1.0 + 5.0*(zj^4 - 2.0*r*sqr(zj))/sqr(r)) + cri
zj = -zj*(zj^4 - 10.0*r*sqr(zj) + 5.0*sqr(r)) + cj
elseif @fractal3D==13 ; White/Nylander quick z^6+c
if @useDE
th = atan2(zri)
ph = asin(zj/(magn=sqrt(magn)))
dr = 6.0*(r=sqrt(|dzri|+sqr(dzj)))*magn^5
dth = 5.0*th + atan2(dzri)
dph = 5.0*ph + asin(dzj/r)
dzri = dr*(cos(dph)*cos(dth) + flip(cos(dph)*sin(dth)))
dzj = -dr*sin(dph)
if @fractaltype==1
dzri = dzri + 1.0
endif
endif
r = |zri|
zri = (zri^6)*(1.0 - (zj^6 - 15.0*r*zj^4 + 15.0*sqr(r)*sqr(zj)) \
/r^3) + cri
zj = -2.0*zj*sqrt(r)*(3.0*zj^4 - 10.0*r*sqr(zj) + 3.0*sqr(r)) + cj
elseif @fractal3D==14 ; White/Nylander quick z^p+c
if @useDE
th = atan2(zri)
ph = asin(zj/(magn=sqrt(magn)))
dr = @mpwr*(r=sqrt(|dzri|+sqr(dzj)))*magn^(@mpwr-1.0)
dth = (@mpwr-1.0)*th + atan2(dzri)
dph = (@mpwr-1.0)*ph + asin(dzj/r)
dzri = dr*(cos(dph)*cos(dth) + flip(cos(dph)*sin(dth)))
dzj = -dr*sin(dph)
if @fractaltype==1
dzri = dzri + 1.0
endif
endif
ztemp = ((r=cabs(zri)) + flip(zj))^@mpwr
zri = real(ztemp)*(zri/r)^@mpwr + cri
zj = -imag(ztemp) + cj
endif
magn = |zri| + sqr(zj)
Note the code for handling the case when (x,y,z) is near zero, it's a while since I put that in and I can't remember whether I checked to see if it's really required or not - in any case it doesn't handle the derivative correctly except maybe the initial iteration.
Also note the fact that magn = |x+iy+jz| is precalculated before entry into the iteration loop.
As you can see several of the formulas (mine) do not yet have analytical DE versions and in the others I've only used the trig version for the derivative but it is possible to get the derivative without the trig.
Also they aren't all fully optimised yet
Also I'm not sure the analytical DE is correct for the "Rucker" formula - it may be that Jos Leys' method for the derivative is not applicable in that case.
April 01, 2012, 09:17:50 PM
3D Fractals with Ultrafractal
UltraFractal
« 1 2 »
Dilber MC Theme by HarzeM