POV Ray 3.7 code
//=========================================================
// Quaternion 3D Julia-Julia-Julia fractal set
//  generated randomly
// (POV-Ray point by point generation)
//=========================================================

#include "colors.inc"    // Standard Color definitions
#include "glass.inc"    // Glass pigment effect
#include "metals.inc"   // Metal pigment effect

  global_settings {assumed_gamma 1.0}
  background { color rgb <00,0.0,0.0> }

 camera {							// set view point (camera) location
    location <-20, 20, -300>
    look_at  <-5, 0,  0>
  }

   light_source {						// set point light sources location
	< -50, 20, -10>
	rgb <1.000000, 1.000000, 1.000000> * 2.0
 }

 light_source {
	< 20, 20, -10>
	rgb <1.000000, 1.000000, 1.000000> * 2.0
 }

#declare R =  .1;    // set radius value
#declare L = 2;     // set square area side
#declare n = 200;     // set number of pixels per area side
#declare st = 1;
#declare Nr = 30;  // Nr = 40
#declare Th = 0;
#declare Ph = 30;

#declare Cx  = -0.7454294;   // JuliaX parameter
#declare Cy = 0;		 // JuliaY parameter
#declare Cz  = 0;  		// JuliaZ parameter

#declare K = array[2*n+1][2*n+1][2*n+1];
  	#for (p, -n, n, st)
     #declare IncX = p*L/n;   // JuliaX increment

    #for (q, -n, n, st)
    #declare IncY = q*L/n;    // JuliaY increment

    #for (r, -n, n, st)
    #declare IncZ = r*L/n;  // JuliaZ increment

      #declare X = IncX;   // start JuliaX
      #declare Y = IncY;   // start JuliaY
      #declare Z = IncZ;   // start  JuliaZ
     #declare K[p+n][q+n][r+n] = 0;

      #for (k,0,Nr)
      #declare XX = X*X - Y*Y - Z*Z + Cx;   // cycle JuliaX
     #declare YY = 2*X*Y + Cy;       		      // cycle JuliaY
     #declare ZZ = 2*X*Z + Cz;       		      // cycle JuliaZ
       #declare X = XX;
       #declare Y = YY;
       #declare Z = ZZ;

       #declare W = X*X +Y*Y + Z*Z;

   #if ( W > R)    		  // escape if
     	#if (W < R + 0.02)     // escape if
   	#declare K[p+n][q+n][r+n] = k;
#end  // end if
#end  // end if
#end  // end for k
#end  // end for q
#end  // end for p
#end  // end for r

#declare Rnd_1 = seed (1153);
#declare Rnd_2 = seed (553);
#declare Rnd_3 = seed (876);

 union{
//  replace n*n*n by  n*n*n*clock for animation
 #for(j,0,n*n*n)  // alternative   n*n*n/16 n*n*n/64  n*n*n/256
   #declare p = int(2*n*rand(Rnd_1));
   #declare q = int(2*n*rand(Rnd_2));
   #declare r = int(2*n*rand(Rnd_3));
    #if(K[p][q][r] > 0)
 	 	sphere {
  		  < -n+p, -n+q, -n+r >, 1
     			 texture {
    			  pigment { color Col_Glass_Yellow }
      		 }
      		 finish { ambient rgb <0.3,0.1,0.1>
     			  diffuse .3
     			 reflection .3
     			 specular 1 }   // plot sphere

translate < 0, -10, 0 >
rotate < 0, Th, Ph > }
#end // end if
#end}  // end for j