【c/c++】Microsoft PPL/並列STL(C++17)/OpenMPを少しだけ比較【C++17】

Visual Studioc++で利用できる並列化手法としては、Microsoft PPL、並列STL(C++17)、OpenMPなどが有ります。この内どれが速いのが知りたくなったので実行速度を比較しました。
手法ごとに使う機能が異なっている(後述)ため、比較としてはかなりよろしくない部分がありますが、まあ大体の傾向は出たかな?という感じです。
CPUはi7 6700(8スレッド)で検証を行いました。

比較に用いた問題

モンテカルロ法による円周率の導出です。
問題の性質として、ループ内での処理の量は均等です。

実装で利用したもの

各並列化手法の実装では「その手法の中で完結すること」に重点を置きました。
例えば結果の集計方法では、OpenMPのreductionやPPLのcombinableなどがありますが、OpenMPではOpenMPの機能を、PPLではPPLの機能を用いて実装を行いました。
そのためどの機能が実行速度に影響しているのか上手く絞りきれていません。というかどの手法もそれぞれに詰めていくと際限が無い*1のでゆるふわメモとして残します。

各種法ごとに利用したもの

Microsoft PPL

concurrency::parallel_for_eachで並列化しました。結果の集計はconcurrency::combinableで行いました。

並列STL

algorithmのstd::for_eachを並列化しました。結果の集計はstd::atomicで行いました

OpenMP

#pragma omp parallel forで並列化しました。結果の集計は#pragma omp atomicで行いました。

コード

#include <chrono> //時間計測その他諸々
#include <cstdio>
#include <iostream>
#include <random>

#include <ppl.h> //parallel_for_each用

#include <algorithm> //for_each用
#include <execution> //algorithmの並列実行ポリシー

#include <omp.h> //OpenMP用

#define num 1000000000
#define threadnum 8;

inline double sqr(double d) { return d * d; }
inline bool check(double x, double y) { return (sqr(x) + sqr(y) > 1.0); }

double Single() { //単一スレッド
	using namespace std;

	mt19937 m = mt19937((int)time(NULL));
	constexpr int tnum = threadnum;
	mt19937_64 mts[tnum];
	for (mt19937_64 &mt : mts) { mt = mt19937_64(m()); }
	uniform_real_distribution<> uniform = uniform_real_distribution<>(0.0, 1.0);

	constexpr long L = num / tnum;
	long count = 0;
	for (int i = 0; i < tnum; i++) {
		for (long j = 0; j < L; j++) {
			if (check(uniform(mts[i]), uniform(mts[i]))) { count++; }
		}
	}

	return (double)(num - count) / (double)(num >> 2);
}

double PPL() { //concurrency::parallel_for_each
	using namespace std;

	mt19937 m = mt19937((int)time(NULL));
	constexpr int tnum = threadnum;
	mt19937_64 mts[tnum];
	for (mt19937_64 &mt : mts) { mt = mt19937_64(m()); }
	uniform_real_distribution<> uniform = uniform_real_distribution<>(0.0, 1.0);

	constexpr long L = num / threadnum;
	concurrency::combinable<long> count;
	concurrency::parallel_for_each(begin(mts), end(mts), [&count, uniform, L](mt19937_64 mt) {
	//concurrency::parallel_for(0, tnum, [&count, &mts, uniform, L](int i) { //こっちのが4~5秒ぐらい遅い
		long c = 0;
		for (long j = 0; j < L; j++) {
			if (check(uniform(mt), uniform(mt))) { c++; }
		}
		count.local() += c;
	});

	return ((double)(num - count.combine(std::plus<long>())) / (double)(num >> 2));
}

double STL() { //c++17の並列for_each
	using namespace std;

	mt19937 m = mt19937((int)time(NULL));
	constexpr int tnum = threadnum;
	mt19937_64 mts[tnum];
	for (mt19937_64 &mt : mts) { mt = mt19937_64(m()); }
	uniform_real_distribution<> uniform = uniform_real_distribution<>(0.0, 1.0);

	constexpr long L = num / threadnum;
	atomic_long count = 0;
	for_each(execution::par_unseq, begin(mts), end(mts), [&count, uniform, L](mt19937_64 mt) {
		long c = 0;
		for (long j = 0; j < L; j++) {
			if (check(uniform(mt), uniform(mt))) { c++; }
		}
		count += c;
	});

	return ((double)(num - count) / (double)(num >> 2));
}

double OMP() { //OpenMP
	using namespace std;

	mt19937 m = mt19937((int)time(NULL));
	constexpr int tnum = threadnum;
	mt19937_64 mts[tnum];
	for (mt19937_64 &mt : mts) { mt = mt19937_64(m()); }
	uniform_real_distribution<> uniform = uniform_real_distribution<>(0.0, 1.0);

	constexpr long L = num / tnum;
	long count = 0;

	#pragma omp parallel for private(uniform)
	for (int i = 0; i < tnum; i++) {
		long c = 0;
		for (long j = 0; j < L; j++) {
			if (check(uniform(mts[i]), uniform(mts[i]))) {
				c++;
			}
		}
		#pragma omp atomic
		count += c;
	}

	return (double)(num - count) / (double)(num >> 2);
}

int main() {
	using namespace std;

	clock_t s, e;
	double pi;

	s = clock();
	pi = Single();
	e = clock();
	printf("StNs\t: %f sec.\tpi = %f\n", (float)(e - s) / CLOCKS_PER_SEC, pi);

	s = clock();
	pi = PPL();
	e = clock();
	printf("PPL\t: %f sec.\tpi = %f\n", (float)(e - s) / CLOCKS_PER_SEC, pi);

	s = clock();
	pi = STL();
	e = clock();
	printf("STL\t: %f sec.\tpi = %f\n", (float)(e - s) / CLOCKS_PER_SEC, pi);

	s = clock();
	pi = OMP();
	e = clock();
	printf("OMP\t: %f sec.\tpi = %f\n", (float)(e - s) / CLOCKS_PER_SEC, pi);

	cin >> pi;

	return 0;
}

実行結果

何度か試してみましたが、傾向としてはPPLが最も高速でした。

StNs    : 57.617001 sec.        pi = 3.141573
PPL     : 10.830000 sec.        pi = 3.141636
STL     : 10.840000 sec.        pi = 3.141597
OMP     : 14.084000 sec.        pi = 3.141583

考察

PPL

parallel_forも試しましたが、OpenMPより若干速い程度でした。
今回は最速となりましたが、parallel_forを用いてループを分割せずに全部回す条件ではOpenMPより5秒程度遅くなりました。
この原因は、今回の問題ではループ内の計算量が均一でしたが、恐らくダイナミックに負荷を調整しようとした結果ではないかと思います。

STL

とくに無いです。

OpenMP

最遅となりました。インデックスを使うアクセス速度が足を引っ張っていると思われ、問題との相性が悪かった結果でしょう。

その他

前回の記事から引き継いでSIMD命令も入れて試しましたが、なんと全部の手法で使わない方が速いという結果になりました。
恐らくキャッシュが足りない辺りが足を引っ張っているものだと思います。

検証できてないこと

PPLの項で触れましたが、ループ内の負荷が均一でない場合どうなるかは見ることができていません。
今回はインデックスが必須な問題ではなかったためOpenMPが不利であったりもしています。
もうちょっと別の問題でも確認してみたいですね。

まとめ

各種法の利点と欠点をまとめます。

PPL

Windows環境では条件付きで最速。
インデックス有りでも十分速い。
機能的にも様々な並列コンテナーやオブジェクトなどが揃っているが、利用可能な環境が少ない。

STL

C++17が使えるのなら、一番何も考えずに使うことができ、高速性も十分である。
インデックスが必要なループの並列化をどうするかが悩みどころ。

OpenMP

今回は最遅となったが、問題によっては最適化のしやすさから十分速くなる可能性がある。
扱える環境の幅広さも魅力。

*1:チョチョイと変えたら「SIMD使った方が遅い」まで行った辺りで泣きそうになりました。