如何优化矩阵乘法运算[重复]答案

【问题标题】：How to optimize matrix multiplication operation [duplicate]如何优化矩阵乘法运算[重复]
【发布时间】：2011-08-29 01:10:05
【问题描述】：

我需要在我的应用程序中执行大量矩阵运算。最耗时的是矩阵乘法。我是这样实现的

template<typename T>
Matrix<T> Matrix<T>::operator * (Matrix& matrix)
{


    Matrix<T> multipliedMatrix = Matrix<T>(this->rows,matrix.GetColumns(),0);

    for (int i=0;i<this->rows;i++)
    {
        for (int j=0;j<matrix.GetColumns();j++)
        {
            multipliedMatrix.datavector.at(i).at(j) = 0;
            for (int k=0;k<this->columns ;k++)
            {
                multipliedMatrix.datavector.at(i).at(j) +=  datavector.at(i).at(k) * matrix.datavector.at(k).at(j);
            }
            //cout<<(*multipliedMatrix)[i][j]<<endl;
        }
    }
    return multipliedMatrix;
}

有没有更好的写法？？到目前为止，矩阵乘法运算在我的应用程序中花费了大部分时间。也许有做这种事情的好/快图书馆吗？但是我宁愿不能使用使用图形卡进行数学运算的库，因为我在带有集成图形卡的笔记本电脑上工作。

【问题讨论】：

一般来说，如果你需要一个相当标准的高效完成的事情，你应该寻找一个做得好的库。它们是由通常比您更了解操作的人编写的，并且他们可以投入比您更多的时间来提高效率。
这些是密集矩阵吗？如果不是（在很多应用程序中你有稀疏矩阵），这个简单的算法将浪费大量时间乘以零。当然，为了有效地处理稀疏矩阵，您需要彻底改变数据的存储方式。
查看Blaze library，它可以很好地处理矩阵-矩阵乘法并提供易于使用的语法。

标签： c++ matrix matrix-multiplication

【解决方案1】：

Boost uBLAS 我认为这绝对是处理这类事情的方法。 Boost 设计精良，经过良好测试并在许多应用程序中使用。

【讨论】：

与其他可用的选项（Eigen、CUBLAS、MKL）相比，Boost uBLAS 的性能非常糟糕。特别是 Eigen 集成了 Boost uBLAS 使用的许多概念，并且可以为各种 CPU 架构发出 SIMD 指令。

【解决方案2】：

考虑GNU Scientific Library，或MV++

如果您对 C 没问题，BLAS 是一个低级库，它结合了 C 和 C 包装的 FORTRAN 指令，并用于大量高级数学库。

我对此一无所知，但另一种选择可能是Meschach 哪个seems to have decent performance。

编辑：关于您关于不想使用使用您的显卡的库的评论，我会指出，在许多情况下，使用您的显卡的库是标准（非 GPU）库的专门实现.例如，BLAS 的各种实现是listed on it's Wikipedia page，只有一些旨在利用您的 GPU。

【讨论】：

【解决方案3】：

有一本书叫Introduction to Algorithms。您可能想查看动态编程的章节。它具有使用动态编程的出色矩阵乘法算法。它值得一读。好吧，这个信息是为了以防你想编写自己的逻辑而不是使用库。

【讨论】：

【解决方案4】：

有很多高效的矩阵乘法算法。

Algorithms for efficient matrix multiplication

查看算法，找到实现。

你也可以为它做一个多线程的实现。

【讨论】：

【解决方案5】：

Eigen 是迄今为止最快的线性代数库之一，如果不是最快的话。它写得很好，而且质量很高。此外，它使用表达式模板，使编写的代码更具可读性。刚刚发布的版本 3 使用 OpenMP 进行数据并行。

#include <iostream>
#include <Eigen/Dense>

using Eigen::MatrixXd;

int main()
{
  MatrixXd m(2,2);
  m(0,0) = 3;
  m(1,0) = 2.5;
  m(0,1) = -1;
  m(1,1) = m(1,0) + m(0,1);
  std::cout << m << std::endl;
}

【讨论】：