Ali*_*eza 11 c++ random multithreading openmp
我正在使用rand()函数生成0,1之间的伪随机数用于模拟目的,但是当我决定使我的C++代码并行运行时(通过OpenMP),我注意到rand()它不是线程安全的,也不是很均匀.
所以我转而使用在其他问题的许多答案中提出的(所谓的)更均匀的生成器.看起来像这样
double rnd(const double & min, const double & max) {
static thread_local mt19937* generator = nullptr;
if (!generator) generator = new mt19937(clock() + omp_get_thread_num());
uniform_real_distribution<double> distribution(min, max);
return fabs(distribution(*generator));
}
Run Code Online (Sandbox Code Playgroud)
但是我在我模拟的原始问题中看到了许多科学错误.既反对结果rand()也反对常识的问题.
所以我编写了一个代码,用这个函数生成500k随机数,计算它们的平均值并做200次并绘制结果.
double SUM=0;
for(r=0; r<=10; r+=0.05){
#pragma omp parallel for ordered schedule(static)
for(w=1; w<=500000; w++){
double a;
a=rnd(0,1);
SUM=SUM+a;
}
SUM=SUM/w_max;
ft<<r<<'\t'<<SUM<<'\n';
SUM=0;
}
Run Code Online (Sandbox Code Playgroud)
我们知道如果不是500k,我可以无限次地做它,它应该是一个值为0.5的简单线.但是有了500k,我们的波动将在0.5左右.
使用单个线程运行代码时,结果是可以接受的:
但这是2个线程的结果:
3个主题:
4个主题:
我现在没有我的8线程CPU,但结果甚至值得.
正如你所看到的,它们都不均匀,并且在平均值附近波动很大.
这个伪随机生成器线程也不安全吗?
或者我在某个地方犯了错误?
Max*_*hof 11
关于我的测试输出有三个观察结果:
它具有比良好的随机源平均应该提供的更强的方差.您通过比较单线程结果自己观察到了这一点.
计算的平均值随线程数减少,并且永远不会达到原始的0.5(即,它不仅仅是更高的方差,而且还改变了平均值).
数据中存在时间关系,特别是在4线程情况下可见.
所有这些都可以通过代码中存在的竞争条件来解释:您SUM从多个线程分配.增加double不是原子操作(即使在x86上,你可能会在寄存器上获得原子读写).两个线程可以读取当前值(例如10),递增它(例如,两者都加0.5)然后将值写回存储器.现在你有两个线程编写10.5而不是正确的11.
尝试同时写入的线程越多SUM(没有同步),它们的更改就会丢失得越多.这解释了所有观察:
线程在每次单独运行中相互竞争的难度决定了丢失的结果数量,这可能因运行而异.
随着更多比赛(例如更多线程),平均值会降低,因为会丢失更多结果.你永远不会超过统计0.5的平均值,因为你只会丢失写入.
当线程和调度程序"安顿下来"时,方差减小.这就是为什么在进行基准测试时应该"预热"测试的原因.
不用说,这是未定义的行为.它只显示x86 CPU上的良性行为,但这不是C++标准保证的.如您所知,double的各个字节可能会被不同的线程同时写入,从而导致完全垃圾.
正确的解决方案是在本地添加双线程,然后(同步)将它们最终添加到一起.OMP为此特定目的有减少条款.
对于整数类型,您可以使用std::atomic<IntegralType>::fetch_add().std::atomic<double>存在但(在C++ 20之前)所提到的函数(和其他函数)仅适用于整数类型.
问题不在你的RNG中,而在于你的总和.只有竞争条件SUM.要解决此问题,请使用缩减,例如将pragma更改为:
#pragma omp parallel for ordered schedule(static) reduction(+:SUM)
Run Code Online (Sandbox Code Playgroud)
请注意,使用thread_localOpenMP在技术上是未定义的行为.它可能在实践中起作用,但是OpenMP和C++ 11线程概念之间的交互没有很好地定义(也参见这个问题).因此,安全的OpenMP替代方案将是:
static mt19937* generator = nullptr;
#pragma omp threadprivate(generator)
Run Code Online (Sandbox Code Playgroud)