Use of Taylor series to speed up computation

Question

In math, Taylor series are important to get approximations of functions, with polynomials of small degree.

I wanted to see how such an approximation can be helpful, for example in order to speed up computations. Let's use the famous Taylor series:

log(1+x) = x + 0.5 * x^2 + (error term)

Morally, computing the value of a polynomial of degree 2 should be much faster than computing a log.

Thus a code to test this:

import numpy, time

def f(t):
    return t + 0.5 * t ** 2
f = numpy.vectorize(f)  

s = time.time()
for i in range(100):
    x = numpy.random.rand(100000) 
    numpy.log(1 + x)
print time.time() - s          # 0.556999921799 seconds

s = time.time()
for i in range(100):
    x = numpy.random.rand(100000)
    f(x)
print time.time() - s          # arghh! 4.81500005722 seconds

Why is the polynomial method 10 times slower than the actual log? I expected the contrary.

PS: This question is probably in the middle of SO and math.SE.

Answer 1

With Python+Numpy, it's probably optimized here and there, and so it's impossible to really benchmark log(1+x) vs x + 0.5 * x^2. So I moved to C++.

Result:

Time per operation with log: 19.57 ns
Time per operation with order-2 Taylor expansion of log: 3.73 ns

So roughly a x5 factor !

---

#include <iostream>
#include <math.h>
#include <time.h>
#define N (1000*1000*100)
#define NANO (1000*1000*1000)

int main()
{
  float *x = (float*) malloc(N * sizeof(float));
  float y;
  float elapsed1, elapsed2;
  clock_t begin, end;
  int i;

  for (i = 0; i < N; i++) 
    x[i] = (float) (rand() + 1) / (float)(RAND_MAX);

  begin = clock();
  for (i = 0; i < N; i++) 
    y = logf(x[i]);
  end = clock();
  elapsed1 = float(end - begin) / CLOCKS_PER_SEC / N * NANO;

  begin = clock();
  for (i = 0; i < N; i++) 
    y = x[i] + 0.5 * x[i] * x[i];  
  end = clock();
  elapsed2 = float(end - begin) / CLOCKS_PER_SEC / N * NANO;

  std::cout << "Time per operation with log: " << elapsed1 << " ns\n";  
  std::cout << "Time per operation with order-2 Taylor epansion: " << elapsed2 << " ns";

  free(x);

}