Hi Warp,
In K3DSurf, I'm using mainly three steps to draw Isosurfaces:
1) Grid initialisation:
Preparing the Grid: the grid is a 3D matrix that holds values obtained from
parsing the Isosurface math formula. This process can be simplified to that
code :
++++++++++++++++++++++++++++++++++++++++++++++++++++
void InitGrig(){
     FunctionParser Isosurface;
     double Var[3];
  double X[100], y[100], Z[100];
  // init the parser:
     Isosurface.Parse("cos(x)+cos(y)+cos(z)+log(abs(cos(x))+1)", "x,y,z");
     // Init the axes to a Cube [-1,1][-1,1][-1,1]
  for(int i=0; i<100; i++) X[i] = y[i] = Z[i] = -1+i*(2/100);
  // Compute the Grid:
     for(i=0; i<limitX; i++) {
         Var[0] = X[i];  // Axe x
      for(j=0; j<limitY; j++) {
         Var[2] = Y[j];  // Axe y
       for(k=0; k<limitZ; k++) {
         Var[3] = Z[k];  // Axe z
         Grid[i][j][k] = Isosurface.Eval(Var);
       }}}
}
+++++++++++++++++++++++++++++++++++++++++++++++++++++
That's it. This simple example show how to compute the grid for the
isosurface : F(x, y, z) = cos(x)+cos(y)+cos(z)+log(abs(cos(x))+1)
2) Generate Triangles from the grid:
Use the marching cube algorithm to generate triangles from the matrix Grid.
3) Draw the image in the window:
Use OpenGl to draw thoses triangles.

The steps 2 and 3 are now extremely fast in K3DSurf: This part of code is
able to treat
(as Thorsten said) millions of triangles in less than a second.
Infortunatly, the first step is very slow compared the the others.
PovRay is not using the marching cube algorithm but it's also using a
mathematical parser in the process of generating the 3D images: even if
PovRay holds the fastest implementation of the raytracing algo, it's going
to be very slow because of it's math parser and this is especially true
when computing big grids and complicated Isosurfaces formulas.
I know that claiming to make a new parser "technologie" is quite crazy but
I'm here to show you that it's really what it looks like. The only
requirement to understand easily what I'm doing is just to have some
understanding of the 3D drawing technics because this solution is inspired
mainly from some basic geometry thoughts.
Now I'll try to explain what's the relationship between the parser and
Isosurfaces:
The original parser you created was calculating a single point at a time.
This is abvious when you look at the code shown in the top and specialy the
line:
>> Grid[i][j][k] = Isosurface.Eval(Var); // Var holds p(x,y,z) of a single point
Imagine now that we want to define a new Eval() fct that can take as
parameter sets of N points defined like this: P(X[N], Y[N], Z[N]) where X, Y
and Z and vectors.
For N =64,the code of the newly defined fct NewEval() will looks like this:

The original code:
-------------------------------------------
    for(IP=0; IP<ByteCodeSize; ++IP){
        switch(ByteCode[IP])
        {
         case   cCos: Stack[SP] = cos(Stack[SP]);
                      break;
........

// Variables:
          default:
              Stack[++SP] = Vars[ByteCode[IP]-VarBegin];
        }
    }

Where Stack = new double[StackSize];
and   Vars  = new double[3];
---------------------------------------------

Changed to :
---------------------------------------------

    for(IP=0; IP<ByteCodeSize; ++IP)
    {
        switch(ByteCode[IP])
        {
         case   cCos:
NewStack = &StackArray[SP*64];
for(i=0; i<64; i++) NewStack[i] = cos(NewStack[i]);
break;
........

// Variables:
          default:
    ++SP;
    NewStack = &StackArray[SP*64];
     if((ByteCode[IP]-VarBegin) == 0) memcpy(NewStack, &Vars[0]  ,
64*sizeof(double));
else if((ByteCode[IP]-VarBegin) == 1) memcpy(NewStack, &Vars[64] ,
64*sizeof(double));
else if((ByteCode[IP]-VarBegin) == 2) memcpy(NewStack, &Vars[124],
64*sizeof(double));

Where StackArray = new double[64*StackSize];
and   Vars       = new double[3*64];
-------------------------------------------------
This code holds almost all the idea behind this new kind of parsers.
Now the best part is that I'm going to use my knowlege of the isosurfaces
geometry to make the newly created function NewEval() working for real:
Imgine we have an isosurface to evaluate for a grid of 100x100x100, your
original code was used like I showd in code on the top.
To use the second newly created fct, the code on the top should be changed
to something like that :

++++++++++++++++++++++++++++++++++++++++++++++++++++
void InitGrig(){
     FunctionParser Isosurface;
     double Var[3*64];
  double X[100], y[100], Z[100];
  // init the parser:
     Isosurface.Parse("cos(x)+cos(y)+cos(z)+log(abs(cos(x))+1)", "x,y,z");
     // Init the axes to a Cube [-1,1][-1,1][-1,1]
  for(int i=0; i<100; i++) X[i] = y[i] = Z[i] = -1+i*(2/100);
  // Compute the Grid:

for(StartI=0; StartI<100; StartI+=4){
    l=0;
    for(i=0; i<4; i++)
       for(j=0; j<4; j++)
          for(k=0; k<4; k++)
                Vars[l++]      = X[i+StartI];
for(StartJ=0; StartJ<100; StartJ+=4){
    l=0;
    for(i=0; i<4; i++)
       for(j=0; j<4; j++)
          for(k=0; k<4; k++)
                Vars[64+(l++)]  = Y[j+StartJ];
for(StartK=0; StartK<100; StartK+=4){
    l=0;
    for(i=0; i<4; i++)
       for(j=0; j<4; j++)
          for(k=0; k<4; k++)
                Vars[128+(l++)]  = Z[k+StartK];


   Grid[i][j][k] = Isosurface.NewEval(Var);
}}}
}
+++++++++++++++++++++++++++++++++++++++++++++++++++++

That's it that's all! This code can be resumed to this : I splited my
isosurfaces calculations to chunks of 64 operations and that amount of
operations is equivalent to a small grid of 4x4x4.
To calculate the hole grid of 100x100x100, I have to use NewEval() : n =
(100x100x100)/(4x4x4) times.
So, now we can use the parser to evaluate more than one operations in the
loop.
I didn't finish yet...I just discovered (and implemented) a very efficient
way of sampling the calculations by using cubic symetrys.
I'll talk about it later.
Cheers,
Taha

Post a reply to this message