Help us understand the problem. What is going on with this article?

C/C++で書いたループを高速化する(その2) openACC

目次

C/C++で書いたループを高速化する(その1) ベクトル化
C/C++で書いたループを高速化する(その2) openACC <- ここ

はじめに

NvidiaコパイラでopenACCを使ってみました。コンパイラの入手方法はその1を見てください。そこそこの計算量がないと差が出にくいので、XYZの3次元座標1024点の相互距離の平均値を計算するプログラムを作ってみました。average_distance1とaverage_distance2の二つの関数を下に書いてあります。データの並びが違うだけです。その1に触れていますが、後者はCPUで実行するとSIMDを使える並びです。

CPUで実行する場合には
nvc++ -fast -Minfo=all distance1.cpp average_distance1.cpp
でコンパイルしています。関数の中にpragmaを記入していますがここでは無視されます。SIMDの使用の有無で速度の違いがでています。

次にGPUで実行する場合です。
nvc++ -acc -Minfo=all distance1.cpp average_distance1.cpp
両関数とものGPUが使われます。
CPUもGPUもだいぶ古いものですが、速くなっていることがわかります。データの並び方で差はほとんどないようです。
このようにループの直前に#pragma~の1行加えただけでGPUが使えるようになります。

結果

CPU SIMD CPU 実行時間(ms) GPU(ms)
average_distance1 X 2.909 0.907
average_distance2 O 1.458 0.963

(実行時間はループ一回分に換算しています)

CPU Intel Core i7-2600 CPU 3.40GHz
GPU GTX 750
CentOS 8.2

使用したプログラム

average_distance1.cpp
#include <cmath>

inline double distance(double x1, double y1, double z1, double x2, double y2, double z2){
    return(sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)+(z1-z2)*(z1-z2)));
}

double average_distance1(int n, const double xyz[][3]){
    double sum=0.0;
    #pragma acc kernels
    for(int i=0; i<n; i++)
        for(int j=i+1; j<n; j++)
            sum+=distance(xyz[i][0], xyz[i][1], xyz[i][2], xyz[j][0], xyz[j][1], xyz[j][2]);
    return(sum/(n*(n-1.0)/2.0));
}

double average_distance2(int n, const double xyz2[][1024]){
    double sum=0.0;
    #pragma acc kernels
    for(int i=0; i<n; i++)
        for(int j=i+1; j<n; j++)
            sum+=distance(xyz2[0][i], xyz2[1][i], xyz2[2][i], xyz2[0][j], xyz2[1][j], xyz2[2][j]);
    return(sum/(n*(n-1.0)/2.0));
}
distance1.cpp
#include <iostream>
#include <chrono>
#include <random>

using namespace std;

double average_distance1(int n, const double xyz[][3]);
double average_distance2(int n, const double xyz2[][1024]);

int main(){
    std::mt19937 mt(123);
    std::uniform_real_distribution<double> u_rand(-10.0, 10.0);

    double data1[1024][3];
    for(int i=0; i<1024; i++){
         data1[i][0]=u_rand(mt);
         data1[i][1]=u_rand(mt);
         data1[i][2]=u_rand(mt);
    }

    double data2[3][1024];
    for(int i=0; i<1024; i++){
         data2[0][i]=data1[i][0];
         data2[1][i]=data1[i][1];
         data2[2][i]=data1[i][2];
    }

    chrono::system_clock::time_point start, end;
    double time;
    volatile double result;

    start=chrono::system_clock::now();
    for(int i=0; i<1000; i++){
        result=average_distance1(1024, data1);
    }
    end=chrono::system_clock::now();
    time=static_cast<double>(chrono::duration_cast<chrono::microseconds>(end - start).count()/1000.0);
    cout << result << endl;
    cout << time << endl;

    start=chrono::system_clock::now();
    for(int i=0; i<1000; i++){
        result=average_distance2(1024, data2);
    }
    end=chrono::system_clock::now();
    time=static_cast<double>(chrono::duration_cast<chrono::microseconds>(end - start).count()/1000.0);
    cout << result << endl;
    cout << time << endl;
}
ki073
業務の手抜き?をできるようにRubyやRを中心に使っています。今は高速化をねらってHaskellにも手を出しています。最適化問題で困っている人が多いので最近ではGLPKのプログラムを多く公開しています。質問や要望などがありましたらコメント欄に書き込んでください。
Why not register and get more from Qiita?
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away