2014-10-11 46 views
2

我目前正在爲C++編寫一個linalg庫,用於教育目的和個人使用。作爲它的一部分,我使用自定義行和列迭代器實現了一個自定義矩陣類。雖然提供了使用std :: algorithm和std :: numeric函數的非常好的功能,但我在索引和iterator/std :: inner_product方法之間執行了矩陣乘法的速度比較。結果顯著不同:C++速度比較迭代器與索引

// used later on for the custom iterator 
template<class U> 
struct EveryNth { 
    bool operator()(const U&) { return m_count++ % N == 0; } 
    EveryNth(std::size_t i) : m_count(0), N(i) {} 
    EveryNth(const EveryNth& element) : m_count(0), N(element.N) {} 

private: 
    int m_count; 
    std::size_t N; 
}; 

template<class T, 
     std::size_t rowsize, 
     std::size_t colsize> 
class Matrix 
{ 

private: 

    // Data is stored in a MVector, a modified std::vector 
    MVector<T> matrix; 

    std::size_t row_dim;     
    std::size_t column_dim; 

public: 

    // other constructors, this one is for matrix in the computation 
    explicit Matrix(MVector<T>&& s): matrix(s), 
            row_dim(rowsize), 
            column_dim(colsize){ 
    }  

    // other code... 

    typedef boost::filter_iterator<EveryNth<T>, 
            typename std::vector<T>::iterator> FilterIter; 

    // returns an iterator that skips elements in a range 
    // if "to" is to be specified, then from has to be set to a value 
    // @ param "j" - j'th column to be requested 
    // @ param "from" - starts at the from'th element 
    // @ param "to" - goes from the from'th element to the "to'th" element 
    FilterIter begin_col(std::size_t j, 
          std::size_t from = 0, 
          std::size_t to = rowsize){ 
     return boost::make_filter_iterator<EveryNth<T> >(
      EveryNth<T>(cols()), 
      matrix.Begin() + index(from, j), 
      matrix.Begin() + index(to, j) 
      ); 
    } 

    // specifies then end of the iterator 
    // so that the iterator can not "jump" past the last element into undefines behaviour 
    FilterIter end_col(std::size_t j, 
         std::size_t to = rowsize){ 
     return boost::make_filter_iterator<EveryNth<T> >(
      EveryNth<T>(cols()), 
      matrix.Begin() + index(to, j), 
      matrix.Begin() + index(to, j) 
      ); 
    } 

    FilterIter begin_row(std::size_t i, 
          std::size_t from = 0, 
          std::size_t to = colsize){ 
     return boost::make_filter_iterator<EveryNth<T> >(
      EveryNth<T>(1), 
      matrix.Begin() + index(i, from), 
      matrix.Begin() + index(i, to) 
      ); 
    } 

    FilterIter end_row(std::size_t i, 
         std::size_t to = colsize){ 
     return boost::make_filter_iterator<EveryNth<T> >(
      EveryNth<T>(1), 
      matrix.Begin() + index(i, to), 
      matrix.Begin() + index(i, to) 
      ); 
    } 

    // other code... 

    // allows to access an element of the matrix by index expressed 
    // in terms of rows and columns 
    // @ param "r" - r'th row of the matrix 
    // @ param "c" - c'th column of the matrix 
    std::size_t index(std::size_t r, std::size_t c) const { 
     return r*cols()+c; 
    } 

    // brackets operator 
    // return an elements stored in the matrix 
    // @ param "r" - r'th row in the matrix 
    // @ param "c" - c'th column in the matrix 
    T& operator()(std::size_t r, std::size_t c) { 
     assert(r < rows() && c < matrix.size()/rows()); 
     return matrix[index(r,c)]; 
    } 

    const T& operator()(std::size_t r, std::size_t c) const { 
     assert(r < rows() && c < matrix.size()/rows()); 
     return matrix[index(r,c)]; 
    } 

    // other code... 

    // end of class 
}; 

現在在運行的主要功能如下:

int main(int argc, char *argv[]){ 


    Matrix<int, 100, 100> a = Matrix<int, 100, 100>(range<int>(10000)); 


    std::clock_t begin = clock(); 
    double b = 0; 
    for(std::size_t i = 0; i < a.rows(); i++){ 
     for (std::size_t j = 0; j < a.cols(); j++) { 
       std::inner_product(a.begin_row(i), a.end_row(i), 
            a.begin_column(j),0);   
     } 
    } 

    // double b = 0; 
    // for(std::size_t i = 0; i < a.rows(); i++){ 
    //  for (std::size_t j = 0; j < a.cols(); j++) { 
    //   for (std::size_t k = 0; k < a.rows(); k++) { 
    //    b += a(i,k)*a(k,j); 
    //   } 
    //  } 
    // } 


    std::clock_t end = clock(); 
    double elapsed_secs = double(end - begin)/CLOCKS_PER_SEC; 
    std::cout << elapsed_secs << std::endl; 

    std::cout << "--- End of test ---" << std::endl; 

    std::cout << std::endl; 
    return 0; 
} 

對於性病:: inner_product /迭代器,它需要的方法:

bash-3.2$ ./main 

3.78358 
--- End of test --- 

和索引(// out)方法:

bash-3.2$ ./main 

0.106173 
--- End of test --- 

這比迭代器方法快了近40倍。你看到代碼中的任何東西都會減慢迭代器計算的速度嗎?我應該提及的是,我嘗試了兩種方法,併產生了正確的結果。

謝謝你的想法。

+2

它可能屬於代碼審查......你也應該補充你的編譯器,以及所有編譯器選項。 – quantdev 2014-10-11 12:17:21

+0

看起來你並沒有使用你的迭代器迭代,而是每次都創建新的迭代器。或者我讀錯了嗎? – Galik 2014-10-11 12:21:25

+0

@Galik:更新了帖子。複製時發生錯誤。 – Vincent 2014-10-11 12:31:55

回答

2

你必須明白的是,矩陣操作是很好理解的,編譯器非常擅長對矩陣操作中涉及的事物進行優化。考慮C = AB,其中C是M×N,A是M×Q,B是Q×N。

double a[M][Q], b[Q][N], c[M][N]; 
for(unsigned i = 0; i < M; i++){ 
    for (unsigned j = 0; j < N; j++) { 
    double temp = 0.0; 
    for (unsigned k = 0; k < Q; k++) { 
     temp += a[i][k]*b[k][j]; 
    } 
    c[i][j] = temp; 
    } 
} 

(你不會相信如何誘惑我寫上面FORTRAN IV)

編譯器着眼於這一點,注意到什麼是真正發生的事情是,他通過和步行c步長爲1,b步長爲Q.他消除了下標計算中的乘法,並進行直接索引。

在這一點上,內環的形式爲:

temp += a[r1] * b[r2]; 
r1 += 1; 
r2 += Q; 

而且你周圍的循環(重新)初始化r1和r2每次穿過。

這是您可以做一個簡單的矩陣乘法的絕對最小計算。你做不到這一點,因爲你必須做這些乘法和增加以及索引調整。

你所能做的就是增加開銷。

這就是iterator和std :: inner_product()方法的作用:它增加了公制的開銷。

+0

嗯,我預期該指數的方法要快一些,只是沒有那麼多。所以我開始擔心。 – Vincent 2014-10-11 21:11:44

+0

記住,在C/C++ for循環的語義需要編譯器發射的代碼,重新評估循環的終止條件,其全部,在每次穿過環。編譯器不知道循環體沒有重構矩陣,所以必須檢查。如果您知道矩陣的大小是固定的乘法的時間(這該死的好應該的!),你可以而且應該在初始化緩存,比方說,a.rows()和a.cols()的值部分的形式。 – 2014-10-14 21:09:47

+0

你能提供一個簡短的例子嗎? – Vincent 2014-10-15 06:46:31

2

這只是對低級代碼優化的一些額外信息和一般建議。


要最終找出時間是在低級別的代碼花(緊密循環和熱點),

  1. 您必須能夠實現代碼的多個版本計算相同的結果,採用不同的實施策略。
    • 您需要廣泛的數學和計算知識才能做到這一點。
  2. 您必須檢查拆裝(機器碼)。
  3. 您還必須在指令級採樣分析器下運行您的代碼,以查看機器代碼的哪一部分執行得最爲嚴重(即熱點)。
    • 爲了收集足夠數量的樣本分析器的,你需要在一個緊湊的循環運行的代碼,在數百萬或數十億倍。
  4. 您必須比較不同版本的代碼(來自不同的實施策略)之間的熱點反彙編。
  5. 根據以上信息,您可以得出結論:某些實施策略效率較低(更浪費或多餘)。
    • 如果你到達這一步,你現在可以發佈並與他人分享你的發現。

一些可能性:

  1. 使用boost::filter_iterator爲執行該跳過每N元素是一種浪費的迭代器。內部實現必須一次遞增1。如果N很大,則通過boost::filter_iterator訪問下一個元素將變爲O(N)操作,而不是簡單的迭代器算術,該操作將是O(1)操作。
  2. 您的boost::filter_iterator實現使用模運算符。儘管現代CPU上的整數除法和模運算速度很快,但它仍然不如簡單整數算法那麼快。

簡單地說,

  • 遞增,遞減,加法和減法是最快的,對於整數和浮點運算。
  • 乘法和位移稍慢。
  • 分部和模數操作會下注較慢。
  • 最後,浮點三角函數和超越函數,尤其是那些需要調用標準數學庫函數的函數,將是最慢的。
+0

非常感謝您的評論。我仍然處於開發的早期階段,所以實際上我很滿意迭代器是按照原樣工作的。不過,我相信這個測試顯示還有很多工作還在進行中。到時候我會發布結果和代碼審查。 – Vincent 2014-10-12 11:01:56