使用蒙特卡罗模拟计算 π 的局限性

我问过一个与此非常相似的问题，所以我会在最后提到以前的解决方案，我有一个website https://monte-carlo.sjorsvanholst.nl使用客户端的 CPU 计算 π，同时将其存储在服务器上，到目前为止我已经得到：

“701.766.448.388”指向圆内，总共“893.547.800.000”，这些数字是使用此代码计算的。 （工作示例位于：https://jsfiddle.net/d47zwvh5/2/ https://jsfiddle.net/d47zwvh5/2/)

let inside = 0;
let size = 500;

for (let i = 0; i < iterations; i++) {
  var Xpos = Math.random() * size;
  var Ypos = Math.random() * size;

  var dist = Math.hypot(Xpos - size / 2, Ypos - size / 2);

  if (dist < size / 2) {
    inside++;
  }
}

问题

(4 * 701.766.448.388) / 893.547.800.000 = 3,141483638

这是我们得到的结果，直到第四位数字为止都是正确的，4 应该是 5。

以前的问题：

1. 我搞砸了距离计算。
2. 我将圆的位置设置为 0...499，应该是 0...500
3. 我没有使用浮动，这降低了“分辨率”

免责声明

可能只是我已经达到了极限，但这个演示使用了 100 万点，得到了 3.16。考虑到我有大约 9000 亿，我认为可能更准确。

我确实明白，如果我想计算 π 这不是正确的方法，但我只是想确保一切都正确，所以我希望任何人都能发现错误，或者我只需要更多“点” '。

EDIT:有很多提到这些数字是多么不切实际，这些提到的地方是正确的，我现在已经将它们更新为正确的。

您可以轻松估计应该得到什么样的误差（误差线），这就是蒙特卡罗的美妙之处。为此，您必须计算第二动量并估计方差和标准偏差。好的一点是，收集到的值将与您收集到的平均值相同，因为您只是将 1 接 1 接 1 相加。

然后您可以获得模拟西格玛的估计以及所需值的误差线。抱歉，我不太懂 Javascript，所以这里的代码是 C# 的：

using System;

namespace Pi
{
    class Program
    {
        static void Main(string[] args)
        {
            ulong N = 1_000_000_000UL; // number of samples
            var rng = new Random(312345); // RNG

            ulong v  = 0UL; // collecting mean values here
            ulong v2 = 0UL; // collecting squares, should be the same as mean
            for (ulong k = 0; k != N; ++k) {
                double x = rng.NextDouble();
                double y = rng.NextDouble();

                var r = (x * x + y * y < 1.0) ? 1UL : 0UL;

                v  += r;
                v2 += r * r;
            }

            var mean = (double)v / (double)N;
            var varc = ((double)v2 / (double)N - mean * mean ) * ((double)N/(N-1UL)); // variance
            var stdd = Math.Sqrt(varc); // std.dev, should be sqrt(Pi/4 (1-Pi/4))
            var errr = stdd / Math.Sqrt(N);

            Console.WriteLine($"Mean = {mean}, StdDev = {stdd}, Err = {errr}");

            mean *= 4.0;
            errr *= 4.0;

            Console.WriteLine($"PI (1 sigma) = {mean - 1.0 * errr}...{mean + 1.0 * errr}");
            Console.WriteLine($"PI (2 sigma) = {mean - 2.0 * errr}...{mean + 2.0 * errr}");
            Console.WriteLine($"PI (3 sigma) = {mean - 3.0 * errr}...{mean + 3.0 * errr}");
        }
    }
}

After 10⁹ samples I've got

Mean = 0.785405665, StdDev = 0.410540627166729, Err = 1.29824345388086E-05
PI (1 sigma) = 3.14157073026184...3.14167458973816
PI (2 sigma) = 3.14151880052369...3.14172651947631
PI (3 sigma) = 3.14146687078553...3.14177844921447

看起来差不多。很容易看出，在理想情况下，方差将等于 (Pi/4)*(1-Pi/4)。确实没必要计算v2，只需将其设置为v模拟后。

坦率地说，我不知道为什么你没有得到预期的结果。求和中的精度损失可能是答案，或者我怀疑，由于播种和重叠序列，您的模拟不会产生独立的样本（因此实际的 N 远低于 900 万亿）。

但是使用这种方法您可以控制错误并检查计算的进行情况。

UPDATE

我已经插入了您的数字，以表明您显然低估了价值。代码

    N  = 893_547_800_000UL;
    v  = 701_766_448_388UL;
    v2 = v;

    var mean = (double)v / (double)N;
    var varc = ((double)v2 / (double)N - mean * mean ) * ((double)N/(N-1UL)); 
    var stdd = Math.Sqrt(varc); // should be sqrt(Pi/4 (1-Pi/4))
    var errr = stdd / Math.Sqrt(N);

    Console.WriteLine($"Mean = {mean}, StdDev = {stdd}, Err = {errr}");

    mean *= 4.0;
    errr *= 4.0;

    Console.WriteLine($"PI (1 sigma) = {mean - 1.0 * errr}...{mean + 1.0 * errr}");
    Console.WriteLine($"PI (2 sigma) = {mean - 2.0 * errr}...{mean + 2.0 * errr}");
    Console.WriteLine($"PI (3 sigma) = {mean - 3.0 * errr}...{mean + 3.0 * errr}");

并输出

Mean = 0.785370909522692, StdDev = 0.410564786603016, Err = 4.34332975349809E-07
PI (1 sigma) = 3.14148190075886...3.14148537542267
PI (2 sigma) = 3.14148016342696...3.14148711275457
PI (3 sigma) = 3.14147842609506...3.14148885008647

所以，显然你在某个地方遇到了问题（代码？表示中丢失了准确性？求和中丢失了准确性？重复/非独立采样？）