用一根线模拟浦丰投针实验—

用一根线模拟浦丰投针实验——Java

浦丰投针实验（Buffon’s needle problem）
本文先给出正确的推理过程，错误的推导过程放在文章末尾，方便大家探讨。
不多废话，直奔主题。
浦丰投针实验的内容：
博物学家浦丰1777年邀请朋友在自己家做实验，要求朋友们向桌子上随意抛出准备好的针，在桌子上有等间距的平行线，针的长度小于平行线的间距。朋友们一共抛出了n=2212根针，其中与任意一根平行线相交的针m=704，于是普丰得出π=nm≈3.142\pi=\frac{n}{m}\approx3.142π=mn≈3.142。
抛出的针的图形抽象如下

原理是什么呢？

首先针长小于平行线距离，所以针只有可能与一根线相交。如上图，设针与线相交的角度为α，那么从针的中点做平行线的垂线，并且在离相交平行线最近的针的一端做平行线，与垂线向交，形成直角三角形，形成的夹角也为α，设平行线的间距为a，针的长度为b，因此可以得到结论，要使针与线相交：
x≤12bsin⁡αx \le \frac{1}{2}b\sin\alphax≤21bsinα

建立横坐标为α取值，纵坐标为x的坐标系统，如图：

由定义可知α的取值范围是[0,π\piπ]的闭区间，x的取值区间为[0,12a\frac{1}{2}a21a]，而针以任意角度落在桌子上的概率是相同的，所以所有针投射在桌面的总的概率（表现为矩形的面积）为：

S0=12aπS_0=\frac{1}{2}a\piS0=21aπ

针与线相交的概率是定积分:

SA=∫0π12bsin⁡αdα=bS_A=\int_0^\pi \frac{1}{2}b\sin\alpha d\alpha\,=bSA=∫0π21bsinαdα=b

我们知道，针与线相交的概率：

ρ=mn=SaS0=2baπ\rho=\frac{m}{n}=\frac{S_a}{S_0}=\frac{2b}{a\pi}ρ=nm=S0Sa=aπ2b

此时，我们令针长b=12ab=\frac{1}{2}ab=21a，则可得到结论π=nm\pi=\frac{n}{m}π=mn
以上便是浦丰投针的证明过程，下面给出计算机建模的过程。
假设有一根与坐标轴x重合无限长的线，我们投的针长为b=0.5b=0.5b=0.5，每次投针先用random函数确定针的A端纵坐标y1∈[−0.5,0.5]y_1\in\left [-0.5,0.5 \right ]y1∈[−0.5,0.5]，然后再取随机角度得到B端纵坐标y2y_2y2，只要y1y_1y1和y2y_2y2异号就说明该针与线相交了。
至此，一根针就抛投完毕了，剩下的针抛投方式以此类推。最后只要用总针数除以相交的针数就得到了圆周率的近似值。
那么这一根线的模型是什么样子呢？

注意图中“针1”和“针2”的位置，两根针所处的相对线位置完全一样，即
“针1”相对纵坐标0的线，从中点到线0的长度为x=12bsin⁡α=12∗0.5∗sin⁡90°=0.25x = \frac{1}{2}b\sin\alpha=\frac{1}{2}*0.5*\sin90\degree=0.25x=21bsinα=21∗0.5∗sin90°=0.25

“针2”相对纵坐标1的线，从中点到线1的长度为x=12bsin⁡α=12∗0.5∗sin⁡90°=0.25x = \frac{1}{2}b\sin\alpha=\frac{1}{2}*0.5*\sin90\degree=0.25x=21bsinα=21∗0.5∗sin90°=0.25

因此我们得出结论，在上图中，所有与线相交针的物理状态都可以从坐标为0的线与针相交的集合中找到对应的关系。通俗来说就是在[−0.5,0.5]\left[-0.5,0.5 \right ][−0.5,0.5]区间内与线相交的针可以完全模拟多根线的情况。
最终执行代码得到圆周率近似值3.14。
代码片段如下：

static Random random = new SecureRandom();static double m;static double n = 9999999;// 默认除法运算精度private static final Integer DEF_DIV_SCALE = 50;public static void play() {for (int i = 0; i < n; i++) {throwNeedle();}double pi = divide(n, m, DEF_DIV_SCALE);System.out.println("圆周率PI：" + pi);}public static void throwNeedle() {Needle needle = generateNeedle();if (needle.crossLine) {m++;}}private static Needle generateNeedle() {double y1 = random.nextDouble() - 0.5;double angle = 180 * random.nextDouble();double y2 = y1 + Math.sin(angle) * 0.5;boolean crossLine = false;//判断针与横线是否相交if ((y1 < 0 && y2 > 0) || (y1 > 0 && y2 < 0)) {crossLine = true;}return new Needle(y1, y2, crossLine);}private static Double divide(Double dividend, Double divisor, Integer scale) {if (scale < 0) {throw new IllegalArgumentException("The scale must be a positive integer or zero");}BigDecimal b1 = new BigDecimal(Double.toString(dividend));BigDecimal b2 = new BigDecimal(Double.toString(divisor));return b1.divide(b2, scale, RoundingMode.HALF_UP).doubleValue();}public static void main(String[] args) {play();}

项目代码地址：https://github.com/ChenjiQiute/Buffon-needle-problem.git
下面是实验的模型建立过程，有兴趣的朋友可以了解一下。
无意间接触到了计算圆周率的浦丰投针实验，被数学家深深的折服，因此总想自己去动手去做这个实验，但是总不能学前辈买几千根针去随机抛针吧。正好作为Java开发，用计算机模拟投针，这不是妥妥的么，于是有了以上这篇文章。
原本以为比较简单，但是当把数学模型落实到计算机代码，还是需要一些思考过程。
最初我建立了多根线的模型，确定纵坐标y1后，然后再根据0.5的针长，随机确定y2的坐标，想当然的认为只要y1-y2绝对值不大于0.5，模型就完整了。
这样实际是错误的，我盲目的用现实世界的观察结果去建立计算机模型，导致改变了概率模型，得到了错误的结果，让我一度怀疑程序伪随机造成结果错误。
后来参考了一下网上的代码，发现如果画多根线，会造成循环多次计算，增加算法复杂度，然后由多根线变为两根线，最后精简到一根线也是为了减少计算次数，才算顺利结束。
另外网络上的代码大多只有晦涩的代码没有建立模型的理论说明，看的人一头雾水，不由的想要详细的分享出来。

引申思考：

1、作为一个程序员，大多数时间都是在实现业务，对于数学模型敏感度太低，思考的太少，对自身发展无益。而且如果你想将别人的算法模型转化为生产实践，用计算机建模也是必不可少，所幸我能接触到算法的同事，让我了解了不足。
2、关于编程语言的精度问题，在计算无理数时的劣势就体现出来了，让我又对matlab这个商业计算软件有了新的认识。不管在什么领域，如果你能深耕，总有一片蓝海。
3、对于一根线的证明，缺乏数学的证明思路，不够严谨，这个是在我写文章的时候发现的，有时你认为你明白的道理，让你证明时，才发现其实自己也是一知半解，误打误撞。感叹数学家证明1+1=2的功底。数学思维要加强！
先想这么多，还有很多不足，望各位留言赐教！