爲什麼Matlab比C++快11倍

Question

我正在比較Matlab和C++之間Vanilla調用選項的Monte Carlo定價算法的速度。這與Why is MATLAB so fast in matrix multiplication?並不相同，因爲加速不是由於矩陣乘法（只有一個點積快速完成），但似乎是由於其高效的高斯隨機數生成器。

在Matlab中的代碼已經被矢量化和代碼如下

function [ value ] = OptionMCValue(yearsToExpiry, spot, strike, riskFreeRate, dividendYield, volatility, numPaths ) 

    sd = volatility*sqrt(yearsToExpiry); 
    sAdjusted = spot * exp((riskFreeRate - dividendYield - 0.5*volatility*volatility) * yearsToExpiry); 

    g = randn(1,numPaths); 
    sT = sAdjusted * exp(g * sd); 
    values = max(sT-strike,0);` 
    value = mean(values); 
    value = value * exp(-riskFreeRate * yearsToExpiry); 

end 

如果我運行這個擁有1000萬點的路徑如下

strike = 100.0; 
yearsToExpiry = 2.16563; 
spot = 100.0; 
volatility = 0.20; 
dividendYield = 0.03; 
riskFreeRate = 0.05; 
oneMillion = 1000000; 
numPaths = 10*oneMillion; 

tic 
value = OptionMCValue(yearsToExpiry, spot, strike, riskFreeRate, dividendYield, volatility, numPaths ); 
toc 

我得到

Elapsed time is 0.359304 seconds. 
    12.8311 

現在我在VS2013中用C++做同樣的事情

我的代碼是在OptionMC類，並且是如下

double OptionMC::value(double yearsToExpiry, 
        double spot, 
        double riskFreeRate, 
        double dividendYield, 
        double volatility, 
        unsigned long numPaths) 
{ 
    double sd = volatility*sqrt(yearsToExpiry); 
    double sAdjusted = spot * exp((riskFreeRate - dividendYield - 0.5*volatility*volatility) * yearsToExpiry); 
    double value = 0.0; 
    double g, sT; 

    for (unsigned long i = 0; i < numPaths; i++) 
    { 
     g = GaussianRVByBoxMuller(); 
     sT = sAdjusted * exp(g * sd); 
     value += Max(sT - m_strike, 0.0); 
    } 

    value = value * exp(-riskFreeRate * yearsToExpiry); 
    value /= (double) numPaths; 
    return value; 
} 

的BM代碼如下

double GaussianRVByBoxMuller() 
{ 
double result; 
double x; double y;; 
double w; 

do 
{ 
    x = 2.0*rand()/static_cast<double>(RAND_MAX)-1; 
    y = 2.0*rand()/static_cast<double>(RAND_MAX)-1; 
    w = x*x + y*y; 
} while (w >= 1.0); 

w = sqrt(-2.0 * log(w)/w); 
result = x*w; 

return result; 
} 

我已設置的優化選項以優化在Visual Studio中的速度。

對於10米路徑需要4.124秒。

這比Matlab慢11倍。

任何人都可以解釋區別？

編輯：在進一步測試減速確實似乎是調用GaussianRVByBoxMuller。 Matlab似乎有一個非常有效的實現 - Ziggurat方法。請注意，BM在這裏是次優的，因爲它會產生2個RV，而我只使用1個。解決這個問題會導致2倍加速。

Answer 1

就在現在，您正在生成單線程代碼。猜測，Matlab正在使用多線程代碼。這允許它以大約N倍的速度運行，其中N = CPU中的內核數量。

還有一點比這個故事還多。出現的另一個問題是您正在使用rand()，它使用隱藏的全局狀態。因此，如果您對代碼進行簡單重寫以使其成爲多線程，那麼很有可能衝突超過rand()的內部狀態將阻止您獲得更多速度改進（並且可能會運行得更慢 - 也許相當比較慢）。

爲了解決這個問題，您可以考慮（例如）使用在C++ 11中添加的新隨機數生成（可能是分佈式）類。有了這些，您可以爲每個線程創建一個單獨的隨機數生成器實例，從而避免與其內部狀態發生衝突。

我改寫了你的代碼一點點地使用這些，並調用該函數，得到這個：

double m_strike = 100.0; 

class generator { 
    std::normal_distribution<double> dis; 
    std::mt19937_64 gen; 
public: 
    generator(double lower = 0.0, double upper = 1.0) 
     : gen(std::random_device()()), dis(lower, upper) {} 

    double operator()() { 
     return dis(gen); 
    } 
}; 

double value(double yearsToExpiry, 
    double spot, 
    double riskFreeRate, 
    double dividendYield, 
    double volatility, 
    unsigned long numPaths) 
{ 
    double sd = volatility*sqrt(yearsToExpiry); 
    double sAdjusted = spot * exp((riskFreeRate - dividendYield - 0.5*volatility*volatility) * yearsToExpiry); 
    double value = 0.0; 
    double g, sT; 

    generator gen; 

// run iterations in parallel, with a private random number generator for each thread: 
#pragma omp parallel for reduction(+:value) private(gen) 
    for (long i = 0; i < numPaths; i++) 
    { 
     g = gen(); // GaussianRVByBoxMuller(); 
     sT = sAdjusted * exp(g * sd); 
     value += std::max(sT - m_strike, 0.0); 
    } 

    value = value * exp(-riskFreeRate * yearsToExpiry); 
    value /= (double)numPaths; 
    return value; 
} 

int main() { 
    std::cout << "value: " << value(2.16563, 100.0, 0.05, 0.03, 0.2, 10'000'000) << "\n"; 
} 

我用VC++ 2015年編制本，使用以下命令行：

cl -openmp -Qpar -arch:AVX -O2b2 -GL test.cpp 

在AMD A8-7600上運行時間約爲0.31秒。
在Intel i7處理器上，運行時間爲~16秒。

當然，如果你的CPU擁有更多的核心，那麼你仍然有一個很好的機會運行它。現在看來，我的代碼需要VC++ 2015而不是2013年，但我懷疑它會影響性能。這主要是方便，就像使用10'000'000而不是10000000（但我不打算在這臺機器上安裝2013的副本，只是爲了弄清楚我需要改變以適應它）。

另請注意，與最近的英特爾處理器相比，您可能（或可能不會）通過將arch:AVX更改爲arch:AVX2來獲得一些改進。

單線程代碼的快速檢查表明您的Box-Muller分發代碼可能比標準庫的正常分發代碼要快一點，因此切換到線程友好的版本可能會獲得更多的速度（並且優化版本應該也是的約兩倍）。