Mandelbulb (3D Fractal) on Codea

Hey all,

So I’ve been looking into raytracing lately and decided to try the famous Mandelbulb 3D fractal.
Needless to say, this is a pretty heavy task to render, but I eventually managed to make it run rather smoothly at 1/4 resolution on my 3rd gen iPad.

Here is a video showing a couple different ones running with Codea. Note that the recording lowered the frame rate (very irregularly, but especially during zooming), and it is actually quite better on the device.
http://www.youtube.com/watch?v=9vGyDp5DQ-w

Again, I will release the source as soon as 1.5 lands.

Cheers

I’d really love to see this, but when I try to play the YouTube video it tells me:

This video contains content from UMG, who has blocked it in your country on copyright grounds.

I’m guessing the soundtrack you chose may not be allowed on YouTube in Australia. Any chance of uploading an alternate version?

Not allowed in France too. Please, show us…!

Blocked in US also.

oops, reuploading video (had some daft punk audio from the Tron movie…)

@Jmv38 moi aussi je suis Français ^^

Anyway, I think fractals are really cool ! Being able to zoom in and not reveal pixels or polygons heh… Rendering it split in half (as in, starting the raytracing from inside the fractal) looks is really cool , like a brain scan ^^

Impressive! Is this a volumic 3d computation? That generate transparent/opaque voxels? Then you just plot and shade all the voxels that are not transparent?

@jmv38 - The renderer is distance based using this function I found here


// Distance Estimated mandelbulb
float DE(vec3 pos) {
	vec3 z = pos;
	float dr = 1.0;
	float r = 0.0;
	for (int i = 0; i < Iterations ; i++) {
		r = length(z);
		if (r>Bailout) break;
		
		// convert to polar coordinates
		float theta = acos(z.z/r);
		float phi = atan(z.y,z.x);
		dr =  pow( r, Power-1.0)*Power*dr + 1.0;
		
		// scale and rotate the point
		float zr = pow( r,Power);
		theta = theta*Power;
		phi = phi*Power;
		
		// convert back to cartesian coordinates
		z = zr*vec3(sin(theta)*cos(phi), sin(phi)*sin(theta), cos(theta)) + pos;
	}
	return 0.5*log(r)*r/dr;
}

Thanks.

Sup, added distance based LoD (more detail as you zoom in) and progressive quality (image gets sharper when “standing” still).

Also got rid of all the trig functions, so it renders around twice as fast now !

Since I like zooming into them, I made a final video (might want to hit the mute button).
I apologize for the low quality, messed up with the picture in picture :frowning:

http://www.youtube.com/watch?v=7Eiz3XMhqEU

I’m about done with this now, going back to my other projects, but this was really interesting !

Cheers

Once again, @Xavier, I’m just blown away both by the quality and the performance that you’ve managed with these shaders.

I’ve been beating my head against a 2D game in progress and haven’t had much time to really work with the shaders, but you’ve made me very jealous of your progress.

Man you are a pro!

One word Amazing

@Xavier Would you be able to share the code now even if it can’t be used until 1.5 is out. I would like to see the math involved and what it takes to calculate these shapes. I downloaded the app Mandelbulb HD just to play around with something similar to your example.

@Mark - I’m quite surprised at the power of our tablets ^^ I’m curious to see how an iPad 4 would run it (or my old iPad 1 for that matter:P)

@dave1707 - Here is a code that will work on previous Codea versions (slow as hell, but really means to be readable)


displayMode(FULLSCREEN)
function setup()
    img = image(WIDTH/8, HEIGHT/8)
    mandelbulb = MandelBulb(img)
    mandelbulb:draw()
end

function draw()
    background(0)
    sprite(img, WIDTH/2, HEIGHT/2, WIDTH)
end

--------------------------------------------------------------------------
MandelBulb = class()

function MandelBulb:init(img)
    local x = 0*math.pi/180
    local y = -60*math.pi/180

    self.camPos = vec3(math.cos(x)*math.cos(y), math.sin(y), math.cos(y)*math.sin(x))*2.0

    self.POWER = 8
    self.MAX_MARCH = 30
    self.DETAIL = 0.0025
    self.MAX_ITER = 4
    self.ITER_BAILOUT = 1.25331

    self.MIN_DIST = 0.9
    self.MAX_DIST = 2.3
    self.COLOR_STEP = 255/self.MAX_MARCH

    self.width, self.height = spriteSize(img)
    self.surface = img
    
    self.camTarget = vec3(0, 0, 0)
    self.camUp = vec3(0, 1, 0)
    self.camDir = (self.camTarget - self.camPos):normalize()
    self.camSide = self.camDir:cross(self.camUp):normalize()
    self.camCDS = self.camDir:cross(self.camSide)    
    
    self.ray = vec3(0, 0, 0)
end


function MandelBulb:draw()
    for x=1, self.width do
        for y=1, self.height do
    
            -- make our ray
            local px = (x*2 - self.width)/self.height
            local py = (y*2 - self.height)/self.height
            self.ray = (self.camSide*px + self.camCDS*py + self.camDir):normalize()
    
            -- check if ray goes near our target
            if self:rayInCircle() then
                local col = 255
                local dist = 0
                local total_dist = self.MIN_DIST
                
                for i=1, self.MAX_MARCH do
                
                    total_dist = total_dist + dist
                    dist = self:DE(self.camPos + self.ray*total_dist)
                    
                    -- the farther away, the darker
                    col = col - self.COLOR_STEP

                    if dist < self.DETAIL or total_dist > self.MAX_DIST then
                        break
                    end
                end
		        
                if total_dist > self.MAX_DIST then col = 0 end
                
                self.surface:set(x, y, col, col, col)
            end
        end
    end
end


-- Distance Estimated mandelbulb
function MandelBulb:DE(pos)
    local sin = math.sin
    local cos = math.cos
    local atan2 = math.atan2
    local acos = math.acos
    local pow = math.pow

    local w = vec3(pos.x, pos.y, pos.z)
    local dr = 1
    local r = 0
    
    for i=1, self.MAX_ITER do
        r = w:len()
        if r > self.ITER_BAILOUT then
            break
        end

        dr = pow( r, self.POWER - 1) * self.POWER * dr + 1
        
        -- convert to polar coordinates
        local theta = acos(w.y/r)
        local phi = atan2(w.x,w.z)
        
        -- scale and rotate the point
        r = pow(r, self.POWER)
        theta = theta * self.POWER
        phi = phi * self.POWER
        
        -- convert back to cartesian coordinates
        w = vec3(sin(theta)*sin(phi), cos(theta), sin(theta)*cos(phi))
        w = w*r + pos

    end
    return 0.5*math.log(r)*r/dr
end

-- Check to see if ray near mandelbulb
function MandelBulb:rayInCircle()
    local rdt = self.camPos:dot(self.ray)
    local rdr = self.camPos:dot(self.camPos) - 1.25331 -- sqrt(PI/2)
    return (rdt*rdt)-rdr > 0
end

Cheers

.@Xavier Could you post the code with shaders as well? I have an iPad 4 and am in the beta so could run it.

@Andrew_Stacey - here is the shader version http://pastebin.com/PMq3jBHQ
Warning : Not clean, not documented

Note to others, this will not run on 1.4.6

@Xavier Thanks for the code. Speed doesn’t matter because I was just interested in seeing the math. I down loaded a Mandelbulb program on my PC that allows all kinds of Mandelbulb views based on different parameters, so that’s what I use for speed.

Wow, the shader version is really nice :slight_smile:

Thanks for sharing. Besides the great mandelbuld program, your example of shader-in-project is great!

@Jmv38 thanks, I really like the option of being able to build my shader code from a pool of various modules. It’s used a lot in webgl (well at least I use it a lot)

This is how my fragment shader is built on my raytracer

    local fragmentSource =
    header..
    intersectFunctions..
    inShadowFunctions..
    getColor..
    main
    
    return fragmentSource

```


This is how the getColor string is built
~~~
--------------------------
---  get color of ray  ---
--------------------------
local getColor =
[[
vec3 getColor(vec3 vo, vec3 vd)
{
    float coef = 1.0;
    vec3 col = vec3(0., 0., 0.);
    for (int i=0; i<5; i++)
    {
    float t = INFINITY;
]]..
self:getIntersectCode()..
self:getMinimumIntersectCode()..
[[
    if (t == INFINITY) return col;
    
    vec3 intersection = vo + vd*t;
    vec3 normal;
]]..
self:getNormalCode()..
[[
    vec3 lo = intersection;
    vec3 ld;
    vec3 lightColor;
]]..
self:loopLights()..
self:getShadowCode()..
[[
    float n = dot(normal,ld)*coef;
    vec3 blinn = normalize(ld-vd);
    float b = max(dot(blinn, normal), 0.0);
    b = pow(b, 60.)*coef;
]]..
self:getObjectColor()..
[[
    }
    
    // bounce ray
    coef *= reflection;
    if (coef == 0.0) return col;
    float reflect = 2.* dot(normal, vd);
    vo = intersection;
    vd = vd - reflect*normal;
    }
    
    return col;
}
]]
~~~

There are a few downsides obviously, like the debugging which is harder since you get no error messages, especially if you're not familiar with glsl (or if you mistype flaot somewhere in your 150 lines of code :P)

Cheers