PySpark, 利用蒙特卡洛测算圆柱体体积

发表于 | 分类于 |

蒙特卡洛算法

刚开始我以为蒙特卡洛是一个人的名字，后来才发现是一个城市的名字，而且是一个赌城的名字。有趣的是还有个
算法叫做拉斯维加斯算法。另外，不知道有没有澳门算法。

先简单说说蒙特算法的原理吧。就比如说你要在一堆苹果中找出一个最大的苹果。一个一个地拿过来，第一个拿出来
然后它就是当前最大的，然后拿第二个，和刚拿的第一个比较，看看谁大，如果比第一个大的话，就讲第二个苹果作为
最大的苹果，否则第一个依然视为最大的。

就这样一个一个苹果不断地拿出来比较、迭代。拿出的苹果越多就越接近所有苹果中最大的那个苹果。

蒙特卡洛计算Pi

先拿Spark官方的一个例子，用pyspark来计算Pi的值

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
import sys
from random import random
from operator import add
from pyspark.sql import SparkSession
if __name__ == "__main__":
    """
        Usage: pi [partitions]
    """
    spark = SparkSession\
        .builder\
        .appName("PythonPi")\
        .getOrCreate()
    partitions = int(sys.argv[1]) if len(sys.argv) > 1 else 2
    n = 100000 * partitions
    def f(_):
        x = random() * 2 - 1
        y = random() * 2 - 1
        return 1 if x ** 2 + y ** 2 <= 1 else 0
    count = spark.sparkContext.parallelize(range(1, n + 1), partitions).map(f).reduce(add)
    print("Pi is roughly %f" % (4.0 * count / n))
    spark.stop()

解析下这段代码。关键算法是在f(_)。就是在一定的数值范围内，随机生成x,y 两个数。如果x,y 符合x ** 2 + y ** 2 <= 1
这公式的话就返回1，否则返回0

通过pyspark迭代100000次，执行100000次f(_)这个方法，这样就会生成一个以1,0组合成，length为100000
的数组。然后sum这个大数组，把数组里面的1全部累加起来。然后，看累加起来的数与100000相除，也就是算出1在
数组中出现的比例，再用这个比例乘以4就可以近似算出Pi的值。执行次数越多，越接近Pi的准确值。

蒙特卡洛算体积

算体积的原理也是类似。比如我们算球的体积，可以假设球被一个立方体包住，然后随机生成很多个点（x,y,z) 遍布在
立方体的方位内，如果满足球的体积公式则返回1，否则返回0。然后看看遍历后占所有随机点的比例有多大。再将这个比例乘以
立方体的体积，就能得出球的体积

代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
radius = 4
def f1(_):
    x = uniform(0 - radius, 0 + radius)
    y = uniform(1 - radius, 1 + radius)
    z = uniform(-2 - radius, -2 + radius)
    return 1 if x ** 2 + (y - 1) ** 2 + (z + 2) ** 2 <= 16 else 0
N = 10000000
partitions = 1
cube_volume = (radius * 2) ** 3
count = spark.sparkContext.parallelize(range(1, N + 1), partitions).map(f1).reduce(add)
print(count)
print("Volume is roughly %f" % (count / N * cube_volume))

总结

只要掌握蒙特卡洛原理，实现起来并不难。有问题欢迎留言，或联系我。