我正在为具有不同数据结构和技术(向量,数组和OpenMP)的矩阵实现C乘法,我发现一个奇怪的情况…我的动态数组版本工作更好:
时间:
openmp mult_1: time: 5.882000 s
array mult_2: time: 1.478000 s
我的汇编标志是:
/usr/bin/g++ -fopenmp -pthread -std=c++1y -O3
C矢量版本
- typedef std::vector<std::vector<float>> matrix_f;
- void mult_1 (const matrix_f & matrixOne,const matrix_f & matrixTwo,matrix_f & result) {
- const int matrixSize = (int)result.size();
- #pragma omp parallel for simd
- for (int rowResult = 0; rowResult < matrixSize; ++rowResult) {
- for (int colResult = 0; colResult < matrixSize; ++colResult) {
- for (int k = 0; k < matrixSize; ++k) {
- result[rowResult][colResult] += matrixOne[rowResult][k] * matrixTwo[k][colResult];
- }
- }
- }
- }
动态阵列版本
- void mult_2 ( float * matrixOne,float * matrixTwo,float * result,int size) {
- for (int row = 0; row < size; ++row) {
- for (int col = 0; col < size; ++col) {
- for (int k = 0; k < size; ++k) {
- (*(result+(size*row)+col)) += (*(matrixOne+(size*row)+k)) * (*(matrixTwo+(size*k)+col));
- }
- }
- }
- }
测试:
C矢量版本
- utils::ChronoTimer timer;
- /* set Up simple matrix */
- utils::matrix::matrix_f matr1 = std::vector<std::vector<float>>(size,std::vector<float>(size));
- fillRandomMatrix(matr1);
- utils::matrix::matrix_f matr2 = std::vector<std::vector<float>>(size,std::vector<float>(size));
- fillRandomMatrix(matr2);
- utils::matrix::matrix_f result = std::vector<std::vector<float>>(size,std::vector<float>(size));
- timer.init();
- utils::matrix::mult_1(matr1,matr2,result);
- std::printf("openmp mult_1: time: %f ms\n",timer.now() / 1000);
动态阵列版本
- utils::ChronoTimer timer;
- float *p_matr1 = new float[size*size];
- float *p_matr2 = new float[size*size];
- float *p_result = new float[size*size];
- fillRandomMatrixArray(p_matr1,size);
- fillRandomMatrixArray(p_matr2,size);
- timer.init();
- utils::matrix::mult_2(p_matr1,p_matr2,p_result,size);
- std::printf("array mult_2: time: %f ms\n",timer.now() / 1000);
- delete [] p_matr1;
- delete [] p_matr2;
- delete [] p_result;
我正在检查一些以前的帖子,但我找不到任何与我的问题相关的link,link2,link3:
更新:
我重新测试了答案,而矢量作品稍微好一些:
vector mult: time: 1.194000 s
array mult_2: time: 1.202000 s
C矢量版本
- void mult (const std::vector<float> & matrixOne,const std::vector<float> & matrixTwo,std::vector<float> & result,int size) {
- for (int row = 0; row < size; ++row) {
- for (int col = 0; col < size; ++col) {
- for (int k = 0; k <size; ++k) {
- result[(size*row)+col] += matrixOne[(size*row)+k] * matrixTwo[(size*k)+col];
- }
- }
- }
- }
动态阵列版本
- void mult_2 ( float * matrixOne,int size) {
- for (int row = 0; row < size; ++row) {
- for (int col = 0; col < size; ++col) {
- for (int k = 0; k < size; ++k) {
- (*(result+(size*row)+col)) += (*(matrixOne+(size*row)+k)) * (*(matrixTwo+(size*k)+col));
- }
- }
- }
- }
此外,我的矢量化版本工作更好(0.803秒);
解决方法
向量向量类似于锯齿状数组 – 每个条目都是一个指针的数组,每个指针指向另一个浮点数组.
相比之下,原始阵列版本是一个内存块,您可以在数学中找到元素.
使用单个向量,而不是向量向量,并手动进行数学运算.或者,使用fix-sized std ::数组的向量.或者写一个需要浮点数(一维)向量的帮助器类型,并给出一个2维视图.
连续缓冲区中的数据比一系列断开连接的缓冲区中的数据更有缓存和优化.
- template<std::size_t Dim,class T>
- struct multi_dim_array_view_helper {
- std::size_t const* dims;
- T* t;
- std::size_t stride() const {
- return
- multi_dim_array_view_helper<Dim-1,T>{dims+1,nullptr}.stride()
- * *dims;
- }
- multi_dim_array_view_helper<Dim-1,T> operator[](std::size_t i)const{
- return {dims+1,t+i*stride()};
- }
- };
- template<class T>
- struct multi_dim_array_view_helper<1,T> {
- std::size_t stride() const{ return 1; }
- T* t;
- T& operator[](std::size_t i)const{
- return t[i];
- }
- multi_dim_array_view_helper( std::size_t const*,T* p ):t(p) {}
- };
- template<std::size_t Dims>
- using dims_t = std::array<std::size_t,Dims-1>;
- template<std::size_t Dims,class T>
- struct multi_dim_array_view_storage
- {
- dims_t<Dims> storage;
- };
- template<std::size_t Dims,class T>
- struct multi_dim_array_view:
- multi_dim_array_view_storage<Dims,T>,multi_dim_array_view_helper<Dims,T>
- {
- multi_dim_array_view( dims_t<Dims> d,T* t ):
- multi_dim_array_view_storage<Dims,T>{std::move(d)},T>{
- this->storage.data(),t
- }
- {}
- };
现在你可以这样做:
- std::vector<float> blah = {
- 01.f,02.f,03.f,11.f,12.f,13.f,21.f,22.f,23.f,};
- multi_dim_array_view<2,float> view = { {3},blah.data() };
- for (std::size_t i = 0; i < 3; ++i )
- {
- std::cout << "[";
- for (std::size_t j = 0; j < 3; ++j )
- std::cout << view[i][j] << ",";
- std::cout << "]\n";
- }
在视图类中没有复制数据.它只是提供一个多维数组的平面数组的视图.