C: Accessing lookup tables faster?

Question

I have a piece of code that traces 4 sines at a time.

My original code was making roughly 12000 sin() function calls per frame and was running at 30 fps.

I tried optimizing it by generating lookup tables. I ended up with 16 different lookup tables. I declared and load them in a separate header file at the top of my program. Each table is declared like so:

static const float d4_lookup[800] {...};

Now, with this new method I actually lost fps?! I'm running at 20 fps now instead of 30. Each frame now only has to do 8 sin / cos calls and 19200 lookup calls vs 12000 sin() calls. I compile using gcc with -O3 flag on. At the moment, the lookup tables are included at the top and are part of the global scope of the program.

I assume I'm not loading them in the right memory or something to that effect. How can I speed up the lookup time?

** EDIT 1 **

As requested, here's the function that uses the lookup calls, it is called once per frame:

void
update_sines(void)
{
    static float c1_sin, c1_cos;
    static float c2_sin, c2_cos;
    static float c3_sin, c3_cos;
    static float c4_sin, c4_cos;

    clock_gettime(CLOCK_MONOTONIC, &spec);
    s = spec.tv_sec;
    ms = spec.tv_nsec * 0.0000001;
    etime = concatenate((long)s, ms);

    c1_sin = sinf(etime * 0.00525);
    c1_cos = cosf(etime * 0.00525);
    c2_sin = sinf(etime * 0.007326);
    c2_cos = cosf(etime * 0.007326);
    c3_sin = sinf(etime * 0.0046);
    c3_cos = cosf(etime * 0.0046);
    c4_sin = sinf(etime * 0.007992);
    c4_cos = cosf(etime * 0.007992);

    int k;
    for (k = 0; k < 800; ++k)
    {       
        sine1[k] = a1_lookup[k] * ((bx1_sin_lookup[k] * c1_cos) + (c1_sin * bx1_cos_lookup[k])) + d1_lookup[k];
        sine2[k] = a2_lookup[k] * ((bx2_sin_lookup[k] * c2_cos) + (c2_sin * bx2_cos_lookup[k])) + d2_lookup[k] + 50;
        sine3[k] = a3_lookup[k] * ((bx3_sin_lookup[k] * c3_cos) + (c3_sin * bx3_cos_lookup[k])) + d3_lookup[k];
        sine4[k] = a4_lookup[k] * ((bx4_sin_lookup[k] * c4_cos) + (c4_sin * bx4_cos_lookup[k])) + d4_lookup[k] + 50;
    }

}

** UPDATE **

For anyone reading this thread, I gave up on this problem. I tried using OpenCL kernels, structs, SIMD instructions as well as all the solutions shown here. In the end the original code that computed the sinf() 12800 per frame worked faster than the lookup tables since the lookup tables didn't fit into the cache. Yet it was still only doing 30 fps. It just had too much going on to keep up with my 60fps expectations. I've decided to take a different direction. Thanks to everyone who contributed to this thread. Most of these solutions would probably work to get some half decent speed improvements but nothing like the 200% speed up I needed here to have the lookup tables work the way I wanted.

Answer 1

Sometimes it's hard to know what's slowing you down, but potentially you are going to ruin your cache hits, you could try a lookup of a struct

typedef struct 
{
  float bx1_sin;
  float bx2_sin;
  float bx3_sin;
  float bx4_sin;
  float bx1_cos;
 etc etc
 including  sine1,2,3,4 as well

} lookup_table

then

lookup_table  lookup[800]

now everything at the kth lookup will be in the same small chunk of memory.

also, if you use a macro that takes k as a parameter to do do the contents of the loop lets say SINE_CALC(k), or an inline function...

you can do

for (k = 0; k < 800; ++k)
{
  SINE_CALC(k); k++;
  SINE_CALC(k); k++;
  SINE_CALC(k); k++;
  SINE_CALC(k); k++;
  SINE_CALC(k); k++;
}

if you do a macro, make sure the k++ is outside the macro call like shown

Answer 2

Try unrolling your loops like this:

for (k = 0; k < 800; ++k)
{       
    sine1[k] = a1_lookup[k];
    sine2[k] = a2_lookup[k];
    sine3[k] = a3_lookup[k];
    sine4[k] = a4_lookup[k];
}
for (k = 0; k < 800; ++k)
{       
    sine1[k] *= ((bx1_sin_lookup[k] * c1_cos) + (c1_sin * bx1_cos_lookup[k]));
    sine2[k] *= ((bx2_sin_lookup[k] * c2_cos) + (c2_sin * bx2_cos_lookup[k]));
    sine3[k] *= ((bx3_sin_lookup[k] * c3_cos) + (c3_sin * bx3_cos_lookup[k]));
    sine4[k] *= ((bx4_sin_lookup[k] * c4_cos) + (c4_sin * bx4_cos_lookup[k]));
}
for (k = 0; k < 800; ++k)
{       
    sine1[k] += d1_lookup[k];
    sine2[k] += d2_lookup[k] + 50;
    sine3[k] += d3_lookup[k];
    sine4[k] += d4_lookup[k] + 50;
}

By accessing fewer lookup tables in each loop, you should be able to stay in the cache. The middle loop could be split up as well, but you'll need to create an intermediate table for one of the sub-expressions.

Answer 3

Intel processors can predict serial access (and perform prefetch) for up to 4 arrays both for forward and backward traverse. At least this was true in Core 2 Duo days. Split your for in:

for (k = 0; k < 800; ++k)
    sine1[k] = a1_lookup[k] * ((bx1_sin_lookup[k] * c1_cos) + (c1_sin * bx1_cos_lookup[k])) + d1_lookup[k];
for (k = 0; k < 800; ++k)
    sine2[k] = a2_lookup[k] * ((bx2_sin_lookup[k] * c2_cos) + (c2_sin * bx2_cos_lookup[k])) + d2_lookup[k] + 50;
for (k = 0; k < 800; ++k)
    sine3[k] = a3_lookup[k] * ((bx3_sin_lookup[k] * c3_cos) + (c3_sin * bx3_cos_lookup[k])) + d3_lookup[k];
for (k = 0; k < 800; ++k)
    sine4[k] = a4_lookup[k] * ((bx4_sin_lookup[k] * c4_cos) + (c4_sin * bx4_cos_lookup[k])) + d4_lookup[k] + 50;

I guess you have more cache load than benchmarks in other answers so this does matters. I recommend you not to unroll loops, compilers do it well.

Answer 4

Using a simple sin lookup table will yields >20% speed increase on my linux machine (vm, gcc, 64bit). Interestingly, the size of lookup table (within reasonable < L1 cache size values) does not influence the speed of execution.

Using a fastsin simple implementation from here I got >45% improvement.

Code:

#include <math.h>
#include <stdio.h>
#include <stdint.h>
#include <sys/time.h>
#include <time.h>

#define LOOKUP_SIZE 628


uint64_t currentTimestampUs( void )
{
  struct timeval tv;
  time_t localTimeRet;
  uint64_t timestamp = 0;
  //time_t tzDiff = 0;
  struct tm when;
  int64_t localeOffset = 0;

  {
    localTimeRet = time(NULL);
    localtime_r ( &localTimeRet, &when );
    localeOffset = when.tm_gmtoff * 1000000ll;
  }

  gettimeofday ( &tv, NULL );
  timestamp = ((uint64_t)((tv.tv_sec) * 1000000ll) ) + ( (uint64_t)(tv.tv_usec) );
  timestamp+=localeOffset;

  return timestamp;
}


const double PI   = 3.141592653589793238462;
const double PI2  = 3.141592653589793238462 * 2;

static float sinarr[LOOKUP_SIZE];

void initSinArr() {
  int a =0;
  for (a=0; a<LOOKUP_SIZE; a++) {
    double arg = (1.0*a/LOOKUP_SIZE)*((double)PI * 0.5);
    float sinval_f = sin(arg); // double computation earlier to avoid losing precision on value
    sinarr[a] = sinval_f;
  }
}

float sinlookup(float val) {
  float normval = val;
  while (normval < 0) {
    normval += PI2;
  }
  while (normval > PI2) {
    normval -= PI2;
  }
  int index = LOOKUP_SIZE*(2*normval/PI);

  if (index > 3*LOOKUP_SIZE) {
    index = -index + 4*LOOKUP_SIZE;//LOOKUP_SIZE - (index-3*LOOKUP_SIZE);
    return -sinarr[index];
  } else if (index > 2*LOOKUP_SIZE) {
    index = index - 2*LOOKUP_SIZE;
    return -sinarr[index];
  } else if (index > LOOKUP_SIZE) {
    index = 2*LOOKUP_SIZE - index;
    return sinarr[index];
  } else {
    return sinarr[index];
  }
}


float sin_fast(float x) {
  while (x < -PI)
      x += PI2;

  while (x >  PI)
      x -= PI2;

  //compute sine
  if (x < 0)
      return 1.27323954 * x + .405284735 * x * x;
  else
      return 1.27323954 * x - 0.405284735 * x * x;

}

int main(void) {
  initSinArr();
  int a = 0;
  float val = 0;
  const int num_tries = 100000;

  uint64_t startLookup = currentTimestampUs();

  for (a=0; a<num_tries; a++) {
    for (val=0; val<PI2; val+=0.01) {
      float compval = sinlookup(val);
      (void)compval;
    }
  }

  uint64_t startSin = currentTimestampUs();

  for (a=0; a<num_tries; a++) {
    for (val=0; val<PI2; val+=0.01) {
      float compval = sin(val);
      (void)compval;
    }
  }

  uint64_t startFastSin = currentTimestampUs();

  for (a=0; a<num_tries; a++) {
    for (val=0; val<PI2; val+=0.01) {
      float compval = sin_fast(val);
      (void)compval;
    }
  }
  uint64_t end = currentTimestampUs();

  int64_t lookupMs = (startSin - startLookup)/1000;
  int64_t sinMs = (startFastSin - startSin)/1000;
  int64_t fastSinMs = (end - startFastSin)/1000;
  printf(" lookup: %lld ms\n", lookupMs );
  printf("    sin: %lld ms\n", sinMs );
  printf("   diff: %lld ms\n", sinMs-lookupMs);
  printf("  diff%: %lld %\n", 100*(sinMs-lookupMs)/sinMs);

  printf("fastsin: %lld ms\n", fastSinMs );
  printf("    sin: %lld ms\n", sinMs );
  printf("   diff: %lld ms\n", sinMs-fastSinMs);
  printf("  diff%: %lld %\n", 100*(sinMs-fastSinMs)/sinMs);
}

Sample result:

 lookup: 2276 ms
    sin: 3004 ms
   diff: 728 ms
  diff%: 24 %
fastsin: 1500 ms
    sin: 3004 ms
   diff: 1504 ms
  diff%: 50 %