const unsigned int size=256;
#define RANGE 1.2
const double xlow=-RANGE;
const double xhigh=RANGE;
const double ylow=-RANGE;
const double yhigh=RANGE;
const double zlow=-RANGE;
const double zhigh=RANGE;

const unsigned int maxiterations=80;
const double mandpow=8.0;
double r;

double valInRange(double low, double high, unsigned int size, unsigned int off)
{
  return low+((high-low)/(double)size)*(double)off;
}

unsigned int doPoint(double cx, double cy, double cz)
{
  // program from http://www.treblig.org/3dbrot/3dbrot.c
  double x,y,z;
  double newx,newy,newz;
  double theta,phi,rpow;
  //double r;
  unsigned int i;
  x=0.0;
  y=0.0;
  z=0.0;

  for(i=0;(i<maxiterations) && ((x*x+y*y+z*z) < 2.0);i++)
  {
   /* These maths from http://www.skytopia.com/project/fractal/mandelbulb.html */
    r = sqrt(x*x + y*y + z*z );
    theta = atan2(sqrt(x*x + y*y) , z);
    phi = atan2(y,x);
    rpow = pow(r,mandpow);
   
    newx = rpow * sin(theta*mandpow) * cos(phi*mandpow);
    newy = rpow * sin(theta*mandpow) * sin(phi*mandpow);
    newz = rpow * cos(theta*mandpow);
   
    x=newx+cx;
    y=newy+cy;
    z=newz+cz;
  }
  return i;
}


void createMandelbulb(PolyVox::SimpleVolume<PolyVox::Material8>& volData)
{
   PolyVox::Material8 uValue=1;
   //This three-level for loop iterates over every voxel in the volume
   for (int z = 0; z < volData.getWidth(); z++)
   {
      double fz=valInRange(zlow, zhigh, size, z);
      for (int y = 0; y < volData.getHeight(); y++)
      {
         double fy=valInRange(ylow, yhigh, size, y);
         for (int x = 0; x < volData.getDepth(); x++)
         {
            double fx=valInRange(xlow, xhigh, size, x);
            unsigned int val=doPoint(fx,fy,fz);
            //printf("%i ",val);
            if (val>=maxiterations-1)
            {
               uValue=r*12; //for example use r for material
               volData.setVoxelAt(x, y, z, uValue);
            }
         }
      }
   }
}
...
createMandelbulb(volData);
...