使用泰勒級數來加速計算

Question

在數學中，Taylor series對於使用小函數多項式來近似函數是很重要的。

我想知道這樣的近似如何有用，例如爲了加速計算。讓我們用著名的泰勒級數：

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

道義上，計算等級2的多項式的值應比計算log快得多。

因此，一個碼來測試此：

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 

爲什麼是多項式方法比實際慢日誌10次？我預計相反。

PS：這個問題可能是在SO和數學中間。

Answer 1

隨着Python + Numpy，它可能在這裏和那裏優化，所以它是不可能真正基準log(1+x)與x + 0.5 * x^2。 所以我轉向了C++。

結果：每操作

時間與日誌：19.57納秒每個操作
時間與日誌的順序-2泰勒展開：3。73納秒

所以大致x5因素！

---

#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); 

} 

Answer 2

使用numpy的矢量化操作幾乎總是比你自己的代碼中任何嘗試的優化都快。正如@Divakar提到的那樣，vectorize實際上只是一個編寫for循環的簡便方法，所以你的代碼將比numpy的本地代碼慢。

用標準Python代碼替換numpy的優化例程，表明您的方法速度大致相同。

import math, numpy, time 


def f(t): 
    return t + 0.5 * t ** 2 

x = numpy.random.rand(1000000) 

s = time.time() 
for num in x: 
    math.log(1 + num) 
print (time.time() - s ) 

s = time.time() 
for num in x: 
    f(num) 
print (time.time() - s)  

結果：

1.1951053142547607 
1.3485901355743408 

逼近只是稍微慢一些，但冪是非常昂貴的。與t*t更換t ** 2給出了一個很好的改善，並允許近似略強於大盤python的log

1.1818947792053223 
0.8402454853057861 

編輯：另外，由於大課在這裏進行了優化科學圖書館幾乎任何一天勝過一週的手工編寫的解決方案，這裏是泰勒系列近似與numpy的矢量化操作，這是迄今爲止最快的。請注意，唯一的重大變化是vectorize未在近似函數上調用，因此默認情況下使用numpy的矢量化操作。

import numpy, time 

def f(t): 
    return t + 0.5 * t ** 2 

x = numpy.random.rand(1000000) 
s = time.time() 
numpy.log(1 + x) 
print (time.time() - s) 

s = time.time() 
x = numpy.random.rand(100000) 
f(x) 
print (time.time() - s ) 

結果：

0.07202601432800293 
0.0019881725311279297 

有你有它，矢量化近似是在一個數量級的速度比numpy的真實矢量log。