用于蒙特卡洛集成的线程安全随机数生成 [英] Thread-safe random number generation for Monte-Carlo integration

问题描述

我正在尝试写一些可以很快计算出随机数并且可以应用于多个线程的东西.我当前的代码是:

Im trying to write something which very quickly calculates random numbers and can be applied on multiple threads. My current code is:

/* Approximating PI using a Monte-Carlo method. */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include <omp.h>
#define N 1000000000  /* As lareg as possible for increased accuracy */

double random_function(void);

int main(void)
{
   int i = 0;
    double X, Y;
   double count_inside_temp = 0.0, count_inside = 0.0;
   unsigned int th_id = omp_get_thread_num();
   #pragma omp parallel private(i, X, Y) firstprivate(count_inside_temp)
   {
      srand(th_id);
      #pragma omp for schedule(static)
      for (i = 0; i <= N; i++) {
         X = 2.0 * random_function() - 1.0;
         Y = 2.0 * random_function() - 1.0;
         if ((X * X) + (Y * Y) < 1.0) {
        count_inside_temp += 1.0;
     }
  }
  #pragma omp atomic
      count_inside += count_inside_temp;
   }
   printf("Approximation to PI is = %.10lf\n", (count_inside * 4.0)/ N);
   return 0;
}

double random_function(void)
{
   return ((double) rand() / (double) RAND_MAX);
}

这有效，但是通过观察资源管理器，我知道它没有使用所有线程. rand()是否适用于多线程代码?如果没有，还有一个很好的选择吗?非常感谢.杰克

This works but from observing a resource manager I know its not using all the threads. Does rand() work for multithreaded code? And if not is there a good alternative? Many Thanks. Jack

推荐答案

rand()线程安全吗?也许(也许不是):

Is rand() thread safe? Maybe, maybe not:

rand()函数不需要重新输入.不需要重入的函数不需要是线程安全的."

一种测试和良好的学习习惯是，用一个固定的整数替换对rand()的调用，然后看看会发生什么.

One test and good learning exercise would be to replace the call to rand() with, say, a fixed integer and see what happens.

我认为伪随机数生成器的方式是一个黑匣子，它以整数作为输入，并以整数作为输出.对于任何给定的输入，输出始终是相同的，但是数字序列中没有模式，并且该序列在可能的输出范围内均匀分布. (此模型并不完全准确，但是可以.)使用此黑盒的方式是选择一个凝视数字(种子)，并使用应用程序中的输出值并将其作为下一次调用的输入随机数生成器.设计API的常见方法有两种:

The way I think of pseudo-random number generators is as a black box which take an integer as input and return an integer as output. For any given input the output is always the same, but there is no pattern in the sequence of numbers and the sequence is uniformly distributed over the range of possible outputs. (This model isn't entirely accurate, but it'll do.) The way you use this black box is to choose a staring number (the seed) use the output value in your application and as the input for the next call to the random number generator. There are two common approaches to designing an API:

1. 两个函数，一个用于设置初始种子(例如srand(seed))，另一个用于从序列中检索下一个值(例如rand()). PRNG的状态以某种全局变量的形式存储在内部.生成新的随机数要么不是线程安全的(很难说，但是输出流将不可重现)，要么在多线程代码中将很慢(最终会在状态值周围进行一些序列化).
2. PRNG状态暴露给应用程序程序员的接口.在这里，您通常具有三个函数:init_prng(seed)，它返回PRNG状态的某些不透明表示； get_prng(state)，它返回随机数并更改状态变量；以及destroy_peng(state)，它仅清除分配的内存，依此类推. .带有这种类型的API的PRNG应该都是线程安全的，并且应该在没有锁定的情况下并行运行(因为您要管理(现在是线程本地的)状态变量.

1. Two functions, one to set the initial seed (e.g. srand(seed)) and one to retrieve the next value from the sequence (e.g. rand()). The state of the PRNG is stored internally in sort of global variable. Generating a new random number either will not be thread safe (hard to tell, but the output stream won't be reproducible) or will be slow in multithreded code (you end up with some serialization around the state value).
2. A interface where the PRNG state is exposed to the application programmer. Here you typically have three functions: init_prng(seed), which returns some opaque representation of the PRNG state, get_prng(state), which returns a random number and changes the state variable, and destroy_peng(state), which just cleans up allocated memory and so on. PRNGs with this type of API should all be thread safe and run in parallel with no locking (because you are in charge of managing the (now thread local) state variable.

我通常使用Fortran编写并使用 Ladd的实现Mersenne Twister PRNG(该链接值得一读).在C中有很多合适的PRNG可以将状态暴露给您的控件. PRNG 看起来不错，并使用它(在并行区域和私有状态内进行初始化和销毁​​调用)变量)应该可以让您获得不错的加速效果.

I generally write in Fortran and use Ladd's implementation of the Mersenne Twister PRNG (that link is worth reading). There are lots of suitable PRNG's in C which expose the state to your control. PRNG looks good and using this (with initialization and destroy calls inside the parallel region and private state variables) should give you a decent speedup.

最后，通常情况下，如果您一口气要求完整的随机数序列，可以使PRNG表现更好(例如，编译器可以对PRNG内部进行矢量化处理).因此，这些库通常具有类似get_prng_array(state)函数的功能，这些函数可以将充满随机数的数组返回给您，就像您将get_prng放入填充数组元素的循环中一样-它们可以更快地完成操作.这将是第二次优化(并且需要在并行的for循环内添加一个for循环.显然，您不想这样做会耗尽每个线程的堆栈空间！

Finally, it's often the case that PRNGs can be made to perform better if you ask for a whole sequence of random numbers in one go (e.g. the compiler can vectorize the PRNG internals). Because of this libraries often have something like get_prng_array(state) functions which give you back an array full of random numbers as if you put get_prng in a loop filling the array elements - they just do it more quickly. This would be a second optimization (and would need an added for loop inside the parallel for loop. Obviously, you don't want to run out of per-thread stack space doing this!

这篇关于用于蒙特卡洛集成的线程安全随机数生成的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！