【C++】【vectorクラス】高速なベクトル足し算

Question

###C++のstd::vectorクラスにおいての足し算実装

初歩的なものですが、C++を用いてベクトル同士の足し算を実装中です。
今現在、以下のような単純なもの浮かばないのですが
より高速なアルゴリズムをご存知の方はいませんか？

// ベクトルA,Bのサイズは同じ(サイズをNとする)
vector<int> Addition(vector<int> A, vector<int> B, int N){
    vector<int> temp(N);
    for(int i = 0 ; i < N ; i++){
        temp[i] = A[i] + B[i];
    }
    return temp
}

/* 
例：
A = [1,2,3,4], B = [5,6,7,8]
print Addition(A,B)
> [6,8,10,12]
*/

Answer 1

高速化のための tips は皆さんの助言どおり。
ですがより速い"アルゴリズム"はなさそうです。
ハードウェアの助けを借りて力技でゴリ押すならマルチスレッド(SIMDもアリ)、
もっと欲張るならOpenCLとかCUDAとか。

[追記] ...ってことで、マルチスレッドとCUDAでやってみたわけよ。

#include <thread>
#include <vector>
#include <iostream>
#include <chrono>
#include <cassert>
#include <thread>
#include <cuda_runtime.h>
#include <device_launch_parameters.h>

template<typename Function>
float time_call(Function&& f) {
  using namespace std::chrono;
  auto start = high_resolution_clock::now();
  f();
  auto stop =  high_resolution_clock::now();
  return duration_cast<microseconds>(stop - start).count() / 1000.0f;
}

void Addition_single(std::vector<int>& C, const std::vector<int>& A, const std::vector<int>& B, int N) {
  for ( int i = 0; i < N; ++i ) {
    C[i] = A[i] + B[i];
  }
}

void Addition_multi(std::vector<int>& C, const std::vector<int>& A, const std::vector<int>& B, int N) {
  const int nthr = 4;
  std::thread threads[nthr];
  for ( int t = 0; t < nthr; ++t ) {
    threads[t] = std::thread(
      [&](int begin, int end, int stride) {
        for ( int i = begin; i < end; i += stride ) {
          C[i] = A[i] + B[i];
        }
      },
      t, N, nthr);
  }
  for ( auto& thr : threads ) thr.join();
}

__global__ void kernel(int* pC, const int* pA, const int* pB, int size) {
  int i = blockDim.x * blockIdx.x + threadIdx.x;
  if ( i < size ) { 
    pC[i] = pA[i] + pB[i];
  }
}

void Addition_cuda(std::vector<int>& C, const std::vector<int>& A, const std::vector<int>& B, int N) {
  int* pA;
  int* pB;
  int* pC;
  cudaMalloc(&pA, N*sizeof(int));
  cudaMalloc(&pB, N*sizeof(int));
  cudaMalloc(&pC, N*sizeof(int));
  cudaMemcpy(pA, A.data(), N*sizeof(int), cudaMemcpyHostToDevice);
  cudaMemcpy(pB, B.data(), N*sizeof(int), cudaMemcpyHostToDevice);
  kernel<<<(N+255)/256,256>>>(pC, pA, pB, N);
  cudaMemcpy(C.data(), pC, N*sizeof(int), cudaMemcpyDeviceToHost);
}

int main() {
  int N = 2000000;
  std::vector<int> A(N);
  std::vector<int> B(N);
  for ( int i = 0; i < N; ++i ) {
    A[i] = i; B[i] = -i;
  }
  std::vector<int> C(N);
  
  float single = time_call([&]() { Addition_single(C, A, B, N); });
  for ( int i = 0; i < N; ++i ) assert( T[i] == 0 );
  std::cout << "single = " << single << "[ms]\n";

  float multi  = time_call([&]() { Addition_multi(C, A, B, N); });
  for ( int i = 0; i < N; ++i ) assert( T[i] == 0 );
  std::cout << "multi  = " << multi  << "[ms]\n" ;

  float cuda   = time_call([&]() { Addition_cuda(C, A, B, N); });
  for ( int i = 0; i < N; ++i ) assert( T[i] == 0 );
  std::cout << "cuda   = " << cuda  << "[ms]\n" ;
}

結果: ぜーんぜん速くならんです。
thread/CUDAのオーバヘッドが速くなった分を食ってしまう。
時間計算量 O(N) なのを頑張って速くしようと考えん方がむしろ幸せなんじゃないでしょか。

Answer 2

vectorを使ってるので要素数を渡す必要は無いですね。size()で得ることが出来ます。
それと、引数は参照にすべきです。そうしないとvectorのコピーが発生するのでオーバーヘッドがかかります。
また、std::transformを使えばもう少し簡単に(forが不要)かけます。
高速かどうかは計ってないのでわかりません。

#include <algorithm>

using namespace std;

vector<int> Addition(vector<int>& v1, vector<int>& v2)
{
  vector<int> v3 = vector<int>(v1.size());
  transform(v1.begin(), v1.end(), v2.begin(), v3.begin(), plus<int>());
  return v3;
}

上記は2つのvectorのサイズがv1のほうがサイズが大きい場合は例外が出ると思います。サイズが違ってたら空っぽのvectorを返すとか、小さいサイズを変換の元にするなどが必要かと思います。その辺は自分なりに工夫してください。
（boost.rangeを使えばtransformはもう少し簡潔に記述できます。）

戻り値を新たに返すのではなく、例えば引数のv1に返すという手もあります。そうすれば戻り値を返す時のコピーを防げます。

#include <algorithm>

using namespace std;

void Addition(vector<int>& v1, const vector<int>& v2)
{
  transform(v1.begin(), v1.end(), v2.begin(), v1.begin(), plus<int>());
}

ただ、これなら1行で書けるので関数にする必要もなくなるという考えもあります。

更に追記です。

vectorはデバッグビルドだとかなり速度が遅くなります（VC++ではそうです）。なので、リリースビルドすれば高速に動作する可能性があります。

Answer 3

高速化の要点をまとめます。

• vectorオブジェクトは参照渡しで

PineMatsuさんが指摘するとおり、値渡しだと各要素のコピーが発生し、パフォーマンス低下の要因となります。もしvectorのサイズが100万だった場合100万件分のコピーが発生します。

• vectorオブジェクトはreturnで返さない

Chironianさんが書かれているように、場合によってはRVOの仕組みによりコピーが省かれることもありますが、常にRVOが働くとは限らないので、参照渡しの引数で結果を返した方がパフォーマンス的には有利です。

• 「X = A + B」ではなく、「A += B」とする

オブジェクトの作成と破棄には相応のコストがかかります。前者ではオブジェクトを1個余分に作る必要があるため、パフォーマンスに影響が出ます。場合によってはメモリへの負荷も気になります。

• vectorの各要素のアクセスは、インデックスではなくイテレーターを使う

組み込み配列のインデックス指定だと、今時のコンパイラーならかなりの最適化が期待できますが、vectorのインデックス指定だと最適化があまり効かず、毎回アドレス計算を行ってしまい、高速化の足かせになってしまいます。
その点、イテレーターは中身はポインタなのでアドレス計算のコストが省けます。

• いっそのこと生ポインタを使う

とにかく最大限に高速化したいという場合は、最終的にはこれになります。
vectorの要素はメモリ内で連続して配置されていることが保証されているので、イテレーターの代わりに生ポインタによる連続アクセスが可能です。生ポインタによるアクセスは最適化が効きやすいので高速化が期待できます。
ただし、バッファオーバーランなどの問題を引き起こしやすいので、慎重にコードを書かないといけません。

Answer 4

こんにちは。

高速化ポイントとしては２つ思いつきますが、その両方とも最適化により自動的に高速化されてしまうだろうと思います。

1. for文

最適化されない場合、indexではなくてiteratorを使った方が、微妙に速くなるCPUもあると思います。

1. return文

文法通りに解釈するとこの返し方の場合、std::vector<int>のコピーが発生するので遅いのですが、最近のPC用コンパイラでは、RVO(Return Value Optimization)と言う最適化が働いてコピーされなくなると思います。

Comments

[Chironian]

PineMatsuさんが書かれている引数側の問題を見落としてました。
確かに要素数は態々渡す必要はない筈です。またstd::vector<>の引数は参照かポインタにしないとこちらでコピーが発生します。それは痛いです。

[PineMatsu]

C++11だったら範囲forを使う手もありますが速度的にイテレータとくらべてどうなんでしょうね？処理系に依存するのかな?

[Chironian]

範囲ベースforも中身はイテレータです。ただ、毎回end()関数を呼ぶことになっているので、そこも最適化したい場合はイテレータの方が良いだろうと思います。

[PineMatsu]

なるほど。終了条件をforの外に出してしまうということですね。ありがとうございます。

Answer 5

vectorの要素が100万くらいで、要素の型がintなら実際にわずかな時間しかかからないし、普通のプログラムであれば、全体に対してもほんの僅かなのでは？オイラの環境では、10^7要素で0.1秒です。究極の方法使って最大0.1秒節約。全体が10秒であれば、1%の節約。

profilerを使ってどこで時間がかかっているか調べましょう。要素の型がもっと複雑であれば、その型に応じた最適化があるかも です。

g++ -O3 使って実行時間が変わるか、openMPをloop部分に使って、

<p>#pragma omp parallel for for (int i = 0; i < size; ++i) </p> g++ -fopenmp で試す。 というあたりからがいいのでは？