Welcome to Fractal Forums

vec4   is a 4 dimensional vector

pixel_resolution_3d =  pixel resolution of 3d array;
   //if pixel_resolution_3d=.05      (-.71,.18,.23) =  (-.70,.20,.25)

arrayElements = elements of x,y, and z in the array;  //set odd because of origin

arrayBoundingCubeEdgeLength =  pixel_resolution_3d * (arrayElements-1); 

// outside of the array bounding cube, distance is calculated the old way

x_int, y_int, z_int = integer value of pixels with range of  [1 to arrayElements]

  if  (arrayElements==513)  
     x_int, y_int, and z_int    range from 1 to 513  

     (0.,0.,0.) = (x_int=257, y_int=257, z_int=257) 

  halfBCEL= arrayBoundingCubeEdgeLength/2;   //half the bounding cube edge length

     (-halfBCEL, -halfBCEL, -halfBCEL) = (x_int=1, y_int=1, z_int=1)    =  (1,1,1)  
     ( halfBCEL,  halfBCEL,  halfBCEL) = (x_int=513, y_int=513, z_int=513)   =  (513,513,513)


  At every (x_int, y_int, z_int)  store orbitTrap and DE information for each iteration 
(or application) of formulas so that we can use the data if we want to start recalculating
 from a slightly earlier iteration:

vec4   pointData [x_int][y_int][z_int][iteration] =  
    vec4 (x_orbit [iteration], y_orbit [iteration], z_orbit[iteration], DE [iteration])

   //note we can recover our running derivative, dr (or whatever you call it), from our DE using 
the x,y,and z orbit values for the last iteration, so we have ALL of the information available to 
us to us in the above 4-dimensional array of a 4 dimensional vector

pointData [1][1][1][8]  would be the point data for the 8th iteration of formulas 
applied at the point (-halfBCEL, -halfBCEL, -halfBCEL) 

  1)  Streamline (reduce compute cost of) fractal calculation (ray paths cross the same areas over and over, although not one of them follows the EXACT same path).  

    condition a:  we take our initial DE, DEi, of point (x,y,z) and the fractal is much greater distance away than pixel_resolution_3d, so we assign DE to many points in the area:

         -assign DE (DEi) or an average DE to all points within  (x +/- [DEi * sqrt (3/8)] , y +/- [DEi * sqrt(3/8)] , z +/- [DEi *sqrt(3/8)]) 

                  - average DE, ADEi, of point DEi, using the center of 6 faces of a cube, and its center:  

                 ADEi=  [DE  (x+DEi*sqrt(3/8), y, z)  +  DE  (x - DEi *sqrt(3/8), y, z) + DE (x, y+ DEi *sqrt(3/8) +... DE (x, y, z - DEi * sqrt(3/8) +  DEi] / 7;
         
  We don't have to recalculate DE for every point in every ray- we check to see what point we are closest to and see if it has DE assigned or not- if not, we calculate DE for it, and assign DE to points in the vicinity that do not already have DE assigned.  Because our space is quantized, more than one ray will go through the same point too:

  At resolution_3d= .05, a ray going through  (-.71,.18,.23) will use the DE for (-.70,.20,.25), as will a ray going through (-.69,.21,.26)


  2)  If we apply additional transforms or skews to the fractal, we don't have to recalculate everything, we just apply the transforms to the existing data.  

 // Outside of trace(...) function:

#group Quantization
uniform float CubeCircumsphereRadius;slider[0.01,2.,4.]
uniform int CubePixelsPerRadiusPower2;slider[5,6,10]
float CubePixelsPerRadius=pow(2.,float(CubePixelsPerRadiusPower2));
float resolution3d=CubeCircumsphereRadius/(CubePixelsPerRadius);
float inverseResolution3d=1./resolution3d;

vec3 quantize (vec3 point) {
	point*=inverseResolution3d;
		  //round is not in version 120!!! is version 130!!!!
         //floor is good enough, although round would be better...  
		point = vec3 (floor(point.x),floor(point.y),floor(point.z));
		point *= resolution3d;

	return point;
}

//Inside of trace(..) function's DE check loop:

vec3 p = from + totalDist * direction;
		p=quantize(p);
		dist = DE(p);