Eigen: Langsamer Zugriff auf Spalten von Matrix 4

Question

Ich benutze Eigen für Operationen ähnlich zu Cholesky update, was eine Menge von AXPY (Summe plus Multiplikation mit einem Skalar) auf die Spalten einer Matrix mit fester Größe, typischerweise ein Matrix4d, impliziert. Kurz gesagt, ist es 3-mal teurer zu den Spalten einer Matrix zuzugreifen 4 als zu einem Vektor 4.

Typischerweise wird der Code unter:

for(int i=0;i<4;++i) L.col(0) += x*y[i]; 

ist 3-mal weniger effizient als der Code unten :

for(int i=0;i<4;++i) l4 += x*y[i]; 

wobei L 4.

Außerdem typischerweise eine Matrix der Größe 4, x, y und l4 sind Vektoren der Größe ist, die verbrachte Zeit in der ersten Zeile des Codes nicht auf dem m wird in Abhängigkeit Atrix-Speicherorganisation (entweder RowMajor von ColMajor).

Auf einem Intel i7 (2,5 GHz), dauert es etwa 0,007us für Vektor-Operationen, und 0,02us für Matrix-Operationen (Zeiten werden durch Wiederholen 100000 Mal die gleiche Operation). Meine Anwendung würde Tausende solcher Operationen in Zeiten benötigen, die hoffentlich weit unter der Millisekunde liegen.

Frage: Ich mache etwas falsch beim Zugriff auf Spalten meiner 4x4-Matrix? Gibt es etwas zu tun, um die erste Codezeile effizienter zu machen?

Voll Code für Timings verwendet wird, ist unter:

#include <iostream> 
#include <Eigen/Core> 
#include <vector> 
#include <sys/time.h> 

typedef Eigen::Matrix<double,4,1,Eigen::ColMajor> Vector4; 
//typedef Eigen::Matrix<double,4,4,Eigen::RowMajor,4,4> Matrix4; 
typedef Eigen::Matrix<double,4,4,Eigen::ColMajor,4,4> Matrix4; 

inline double operator- ( const struct timeval & t1,const struct timeval & t0) 
{ 
    /* TODO: double check the double conversion from long (on 64x). */ 
    return double(t1.tv_sec - t0.tv_sec)+1e-6*double(t1.tv_usec - t0.tv_usec); 
} 

void sumCols(Matrix4 & L, 
       Vector4 & x4, 
       Vector4 & y) 
{ 
    for(int i=0;i<4;++i) 
    { 
     L.col(0) += x4*y[i]; 
    } 
} 

void sumVec(Vector4 & L, 
      Vector4 & x4, 
      Vector4 & y) 
{ 
    for(int i=0;i<4;++i) 
    { 
     //L.tail(4-i) += x4.tail(4-i)*y[i]; 
     L   += x4   *y[i]; 
    } 
} 

int main() 
{ 
    using namespace Eigen; 

    const int NBT = 1000000; 

    struct timeval t0,t1; 

    std::vector<  Vector4> x4s(NBT); 
    std::vector<  Vector4> y4s(NBT); 
    std::vector<  Vector4> z4s(NBT); 
    std::vector<  Matrix4> L4s(NBT); 

    for(int i=0;i<NBT;++i) 
    { 
    x4s[i] = Vector4::Random(); 
    y4s[i] = Vector4::Random(); 
    L4s[i] = Matrix4::Random(); 
    } 

    int sample = int(z4s[55][2]/10*NBT); 
    std::cout << "*** SAMPLE = " << sample << std::endl; 

    gettimeofday(&t0,NULL); 
    for(int i=0;i<NBT;++i) 
    { 
     sumCols(L4s[i], x4s[i], y4s[i]); 
    } 
    gettimeofday(&t1,NULL); 
    std::cout << (t1-t0) << std::endl; 
    std::cout << "\t\t\t\t\t\t\tForce check" << L4s[sample](1,0) << std::endl; 

    gettimeofday(&t0,NULL); 
    for(int i=0;i<NBT;++i) 
    { 
     sumVec(z4s[i], x4s[i], y4s[i]); 
    } 
    gettimeofday(&t1,NULL); 
    std::cout << (t1-t0) << std::endl; 
    std::cout << "\t\t\t\t\t\t\tForce check" << z4s[sample][2] << std::endl; 

    return -1; 
} 

Answer 1

Wie ich in einem Kommentar sagte, die generierte Assembly genau die für beide Funktionen gleich ist.

Das Problem ist, dass Ihr Benchmark in dem Sinne verzerrt ist, dass L4s 4 mal größer als z4s ist, und Sie somit mehr Cache-Misses im Matrixfall als im Vektorfall bekommen.