Вот цикл, который я пробовал с std :: vector
и с простым старым double *
.
Для 10 миллионов элементов векторная версия стабильно выполняется примерно в 80% случаев, когда double *
версия принимает; практически любое значение N
, ve ctor заметно быстрее.
Заглянув в исходный код GCC на STL, я не заметил, чтобы std :: vector
делал что-то существенно более интересное, чем то, что делает идиома double *
(т.е. с обычным старым новым []
, оператором []
разыменовывается смещение). Этот вопрос говорит и об этом.
Есть идеи, почему векторная версия быстрее?
Compiler: GCC 4.6.1
Example compile line: g++ -Ofast -march=native -DNDEBUG \
-ftree-vectorizer-verbose=2 -o vector.bin \
vector.cpp -lrt
OS: CentOS 5
CPU: Opteron 8431
RAM: 128 GB
Результаты качественно такие же, если я использую icpc 11.1 или работаю на Xeon. Кроме того, дамп векторизатора говорит, что векторизовалась только операция заполнения в конструкторе std :: vector
.
Векторная версия:
#include
#include
#include
#include "util.h"
#include "rkck_params.h"
using namespace std;
int main( int argc, char* argv[] )
{
const size_t N = boost::lexical_cast( argv[1] );
vector y_old( N );
vector y_new( N );
vector y_err( N );
vector k0( N );
vector k1( N );
vector k2( N );
vector k3( N );
vector k4( N );
vector k5( N );
const double h = 0.5;
const timespec start = clock_tick();
for ( size_t i = 0 ; i < N ; ++i )
{
y_new[ i ] = y_old[ i ]
+ h
*(
rkck::c[0]*k0[ i ]
+ rkck::c[2]*k2[ i ]
+ rkck::c[3]*k3[ i ]
+ rkck::c[5]*k5[ i ]
);
y_err[ i ] = h
*(
rkck::cdiff[0]*k0[ i ]
+ rkck::cdiff[2]*k2[ i ]
+ rkck::cdiff[3]*k3[ i ]
+ rkck::cdiff[4]*k4[ i ]
+ rkck::cdiff[5]*k5[ i ]
);
}
const timespec stop = clock_tick();
const double total_time = seconds( start, stop );
// Output
cout << "vector\t" << N << "\t" << total_time << endl;
return 0;
}
Версия double *
:
#include
#include
#include "util.h"
#include "rkck_params.h"
using namespace std;
int main( int argc, char* argv[] )
{
const size_t N = boost::lexical_cast( argv[1] );
double* y_old = new double[ N ];
double* y_new = new double[ N ];
double* y_err = new double[ N ];
double* k0 = new double[ N ];
double* k1 = new double[ N ];
double* k2 = new double[ N ];
double* k3 = new double[ N ];
double* k4 = new double[ N ];
double* k5 = new double[ N ];
const double h = 0.5;
const timespec start = clock_tick();
for ( size_t i = 0 ; i < N ; ++i )
{
y_new[ i ]
= y_old[ i ]
+ h
*(
rkck::c[0]*k0[ i ]
+ rkck::c[2]*k2[ i ]
+ rkck::c[3]*k3[ i ]
+ rkck::c[5]*k5[ i ]
);
y_err[ i ]
= h
*(
rkck::cdiff[0]*k0[ i ]
+ rkck::cdiff[2]*k2[ i ]
+ rkck::cdiff[3]*k3[ i ]
+ rkck::cdiff[4]*k4[ i ]
+ rkck::cdiff[5]*k5[ i ]
);
}
const timespec stop = clock_tick();
const double total_time = seconds( start, stop );
delete [] y_old;
delete [] y_new;
delete [] y_err;
delete [] k0;
delete [] k1;
delete [] k2;
delete [] k3;
delete [] k4;
delete [] k5;
// Output
cout << "plain\t" << N << "\t" << total_time << endl;
return 0;
}
rkck_params.h
:
#ifndef RKCK_PARAMS_H
#define RKCK_PARAMS_H
namespace rkck
{
// C.f. $c_i$ in Ch. 16.2 of NR in C++, 2nd ed.
const double c[6]
= { 37.0/378.0,
0.0,
250.0/621.0,
125.0/594,
0.0,
512.0/1771.0 };
// C.f. $( c_i - c_i^* )$ in Ch. 16.2 of NR in C++, 2nd ed.
const double cdiff[6]
= { c[0] - 2825.0/27648.0,
c[1] - 0.0,
c[2] - 18575.0/48384.0,
c[3] - 13525.0/55296.0,
c[4] - 277.0/14336.0,
c[5] - 1.0/4.0 };
}
#endif
util.h
:
#ifndef UTIL_H
#define UTIL_H
#include
#include
inline timespec clock_tick()
{
timespec tick;
clock_gettime( CLOCK_REALTIME, &tick );
return tick;
}
// \cite{www.guyrutenberg.com/2007/09/22/profiling-code-using-clock_gettime}
inline double seconds( const timespec& earlier, const timespec& later )
{
double seconds_diff = -1.0;
double nano_diff = -1.0;
if ( later.tv_nsec < earlier.tv_nsec )
{
seconds_diff = later.tv_sec - earlier.tv_sec - 1;
nano_diff = ( 1.0e9 + later.tv_nsec - earlier.tv_nsec )*1.0e-9;
}
else
{
seconds_diff = later.tv_sec - earlier.tv_sec;
nano_diff = ( later.tv_nsec - earlier.tv_nsec )*1.0e-9;
}
return seconds_diff + nano_diff;
}
#endif