About those eight seconds, part II: mea culpa

Monday, August 05, 2013

The images presented in the previous two posts have been mostly reconstructions. The only piece of software I have from that project twenty-some years ago is the piece that drew the chaotic attractors, and that was specific to the SGI and had to be rewritten (extensively) for a modern system.

The other program—the one that generated the maps? I have no clue as to where that code went. Or rather, why I neglected to keep a copy of it. So I had to basically recreate that program from scratch. The original reconstruction started out as:

int mainloop(const double A,const double B,const double C,const double D)
{
  double xn,xn1;
  double yn,yn1;
  int    ix;
  int    iy;
  int    count;
  bool   pix[DIM][DIM];
  
  memset(pix,0,sizeof(pix));
  
  for (count = 0 , xn = yn = 0.5 ; count < MAXLOOP ; count++)
  {
    xn1 = ((A * yn) + B) * xn * (1.0 - xn);
    yn1 = ((C * xn) + D) * yn * (1.0 - yn);    
    xn  = xn1;
    yn  = yn1;

But since I was setting memory and not plotting points, I needed to be very careful about the values of xn and yn. I remembered that the values needed to be between 0 and 1, so attempt one I had:

assert(xn >= 0.0);
assert(xn <  1.0);
assert(yn >= 0.0);
assert(yn <  1.0);

But those assumptions were proven false rather quickly. Attempt two I decided to print the values of xn and yn and as it happened, very quickly (like on the first pass) the values hit +inf. I interpreted that to mean I could bail early, and thus finished the reconstruction with attempt three:

    if (xn <  0.0) return MAXLOOP;
    if (xn >= 1.0) return MAXLOOP;
    if (yn <  0.0) return MAXLOOP;
    if (yn >= 1.0) return MAXLOOP;
    
    ix  = (int)(xn * DIM);
    iy  = (int)(yn * DIM);
    
    if (pix[ix][iy])
      break;
    pix[ix][iy] = true;
  }
  
  return(count);
}

But the sharp edges in the maps were bothering me. So much so that I rewrote the code to continue the loop if the values of xn and yn were out of range, thinking that maybe they might come back into range.

And yes, they did:

[As odd as it sounds, that looks much better]

The downside—this revised code now takes a bit longer to run. The LuaJIT version is still faster than the C version, but the difference (on average 1 minute for the LuaJIT version; 1 minute, 30 seconds for the C version) isn't as jaw-droppingly great.

Sigh—it's not quite as fast as I made it out to be.

Update on Tuesday, August 6^th, 2013

Yes, it is as fast as I think it was …

The Boston Diaries

Monday, August 05, 2013