蒙特卡洛算法的MATLAB实现

一 解决问题

蒙特卡洛算法，即随机取样法，也称为计算机随机模拟方法，源于世界著名赌城——Monte Carlo，算法由此得名。

它是基于概率论的一种计算方法：即通过随机数或者伪随机数来计算比较难以解决的问题，因为设计概率问题，那么结果就不是一个准确的值，而是一个近似值，并且所使用的样本（统计点数）越多，所得到的结果越准确。例如：不规则图形的面积问题等。

二 基本原理

我们以面积计算的问题讲解蒙特卡洛算法的实现步骤如下：

1.画出图像

求要计算图像的基本形状。此处以y1=3x；y2=8-x为例，画出图像如下：

2.确定边界

通常为矩形，要求矩形必须完整包括所求图像。

如上图，因为图像交点为（2,6），因此，确定边界如红框所示：（0,0），（8,0），（8,6），（0,6）

3.随机统计

在（矩形）边界范围内随机产生点，并统计落在所求图像中的点。所用到的函数为unifrnd函数。本例子以10^7（7次方在保证运行速度的情况下，基本可以满足准确度）为例。

4.确定面积

用比值法求面积，即落在图像上的点：整个范围内的点=所求图像面积：矩形边界面积。

即所有图像I面积=矩形边界面积*落在图像上的点/整个范围内的点。

三 相关的知识

本例中使用到的unifrnd函数的相关知识。

R = unifrnd(A,B,M,N,...) or R = unifrnd(A,B,[M,N,...])

即产生m*n阶［a，b］均匀分布U（a，b）的随机数矩阵：unifrnd (a,b,m, n)

如y=(0,8,[1,10000000])，即产生1*10000000阶位于（0,8）之间的数。

四 代码实现

%画函数图像
x=0:0.25:10;
y1=3*x;
y2=8-x;
plot(x,y1,x,y2)

axis([0 10 0 10]);
legend('y1=3x','y2=8-x');
title('蒙特卡洛算法');
text(2,6,'交点');
grid on
%产生随机点
x=unifrnd(0,8,[1,10000000]);
y=unifrnd(0,6,[1,10000000]);
%统计所在所求图形中的点
frequency=sum(y<3*x & x<=2)+sum(y<8-x & x>2);
%计算面积
area=6*8*frequency/10^7;

五 结果测验

由函数图像可知，所求区域面积为：24单位。

算法运行结果如下：

由结果可知，算法每次运行结果都不一样，这是有概率统计的特性所致。但是结果都稳定在24左右，这与理论值是一样的。

 

以上实例以规则图形为例，对于不规则图像，本算法同样适用，且更便利。

今天的分享就到这里，如有不足，欢迎大家评论指导。