蒙特卡罗模拟计算π\piπ

Q11.蒙特卡罗模拟是一种以概率和统计理论为基础的计算方法，它能通过随机数来解决很多计算问题。蒙特卡罗方法将所求解的问题同一定的概率模型相联系，用计算机实现统计模拟或抽样，以获得问题的近似解。本题将利用蒙特卡罗模拟来计算圆周率。

有同学可能会好奇为什么这个关系式会成立，实际上这里用到了概率论与数理统计中几何概率和大数定律的思想，我们这里不深入探究，大家有个印象即可。

现在请大家使用MATLAB生成N个随机点（N取10000），每个随机点的横纵坐标都分别在区间[0,1]上均匀分布，接下来统计出位于1/4圆内的随机点的数量n，然后套用上面的公式即可计算出一个近似的圆周率。将你的N分别变成之前的10倍、100倍和1000倍，你计算出来的圆周率的精度有何变化？

解:随着N取得越大，则越来越逼近π\piπ。但取得太大，数量太大导致计算机运算时间更长，需要更大的内存。

clear,clc
x=rand(1,10000);
y=rand(1,10000);
N=length(x);
scatter(x,y,'.')
hold on
t=0:pi/200:pi/2;
x1=cos(t);
y1=sin(t);
plot(x1,y1,'r','Linewidth',2)      %画1/4圆
for i=1:length(x)if sqrt(x(i)^2+y(i)^2)<=1X(i)=x(i);Y(i)=y(i);        %在1/4圆内的点elseX1(i)=x(i);Y1(i)=y(i);       %不在1/4圆内的点end
end
scatter(X,Y,'b.')         %在1/4圆内的点的散点图
scatter(X1,Y1,'k.')
ind=find(X~=0);
n=length(ind);
format long
pai=4*n/N
e=abs(pai-pi)             %这是误差pai =3.143200000000000e =0.001607346410207

蒙特卡罗模拟求定积分

解：学过微积分的都知道，f(x)f(x)f(x)在[a,b][a,b][a,b]上的积分就是曲线y=f(x)y=f(x)y=f(x)在积分区间是与xxx轴围成的曲边梯形的面积，即图中的S1S1S1部分。用蒙特卡洛的思想就是：曲线下方与xxx轴围成的曲边梯形的面积(S1S1S1)与以a，b边为矩形的面积(a，b边为矩形的面积(a，b边为矩形的面积(S1+S2S1+S2S1+S2)之比，即)之比，即)之比，即P=S1S1+S2=落在区域S1上的点数落在矩形上全部点数P=\frac{S1}{S1+S2}=\frac{落在区域S1上的点数}{落在矩形上全部点数}P=S1+S2S1=落在矩形上全部点数落在区域S1上的点数
设落在区域S1S1S1上的点数为nnn,落在矩形上全部点数为NNN,则
S1S1+S2=落在区域S1上的点数落在矩形上全部点数=nN\frac{S1}{S1+S2}=\frac{落在区域S1上的点数}{落在矩形上全部点数}=\frac{n}{N}S1+S2S1=落在矩形上全部点数落在区域S1上的点数=Nn即∫abf(x)dx=S1=n(S1+S2)N\int_a^b{f(x)}\,{\rm d}x=S1=\frac{n(S1+S2)}{N}∫abf(x)dx=S1=Nn(S1+S2)

clear,clc
x=2*rand(1,10000);
y=2*rand(1,10000);
S=2*2;
N=length(x);
y1=x./(exp(x)-1);
for i=1:length(x)if y(i)<=y1(i)Y(i)=y(i);end
end
ind=find(Y~=0);
n=length(ind);
disp('积分的估计值为：')
jf=S*n/N
积分的估计值为：jf =1.2164

蒙特卡罗模拟计算概率

Q13.清风订了一份报纸，送报人可能在早上6：30至7：30之间把报纸送到他家门口，而清风出门的时间在早上7：00到8：00之间（假设送报人到达的时间和清风出发的时间在对应的时间区间内都是均匀分布的）。请使用蒙特卡罗模拟来估计清风在离开家之前能拿到报纸的概率。注意：这是概率论中几何概型的一道经典题目，准确的概率是7/8.

解：如下图所示，在红线框内即是清风在离开家之前能拿到报纸的概率。通过分析，我们知道它是S1S1S1,则有P=S1S1+S2=落在区域S1上的点数落在小正方形上全部点P=\frac{S1}{S1+S2}=\frac{落在区域S1上的点数}{落在小正方形上全部点}P=S1+S2S1=落在小正方形上全部点落在区域S1上的点数这里为了方便计算，我们记落在小三角形S1S1S1上的nnn,全部点为NNN,则有P=S1S1+S2=N−nNP=\frac{S1}{S1+S2}=\frac{N-n}{N}P=S1+S2S1=NN−n

我们把小正方形按同一比例大小方缩到0~1之间，则有

clear,clc
x=rand(1,10000);
y=rand(1,10000);
t=0.5:0.01:1;
y1=t-0.5;
plot(t,y1,'r','Linewidth',2)
hold on
scatter(x,y,'k.')
for i=1:length(x)if x(i)>=0.5&y(i)<=x(i)-0.5     %找在S1和S2的点X(i)=x(i);Y(i)=y(i);elseX1(i)=x(i);Y1(i)=y(i);end
end
ind=find(X~=0);
ind1=find(X1~=0);
scatter(X,Y,'b.')
scatter(X1,Y1,'k.')
n=length(ind);
N=length(x);
format long
p=(N-n)/N
e=abs(p-7/8)   %误差p =0.876100000000000e =0.001100000000000

枚举法与网格搜索法

Q18.枚举法与网格搜索法：

解：枚举法和网格搜索法不知道的同学可以自行去查一下。
(1)这个很简单，只要给的xxx向量步长足够的小，那么得到的结果就会更精确。但不能太小，前面说了，会耗时间，占内存，只是要理解就可以。

clear,clcx=0:0.0001:1;
y1=2*x;
y2=exp(x);
Y=y1-y2;
plot(x,y1,'b',x,y2,'r')
legend('\ity=2*x','\ity=e(x)')
disp('最大值为：')
M=max(Y)
[ind,location]=find(Y==max(Y))
x(location)
最小值为：M =-0.613705641106080ind =1location =6932ans =0.693100000000000

大家可以看到曲线y=2xy=2xy=2x始终位于曲线y=exy=e^xy=ex的下方，所以最大值固然是个负数，所以通过图形来看，图像离得越近，上面他们的差值就越大。

(2)这里可以去了解一下到数值分析中的共轭梯度法里面的最速下降法，取yyy的相反数f=−y=−(2x−ex)f=-y=-(2x-e^x)f=−y=−(2x−ex)即f=ex−2xf=e^x-2xf=ex−2x。求yyy的最大值即求fff的最小值，感兴趣的同学可以去了解一下。最速下降法之后我会单独讲。
(3)这里说一下曼哈顿距离¹。
d(x,y)=∑i=1n∣xi−yi∣d(x,y)=\sum_{i=1}^n{|x_{i}-y_{i}|}d(x,y)=i=1∑n∣xi−yi∣

这个图中，红、蓝、黄这些轨迹都是曼哈顿距离，那条直线是欧氏距离。举个例子，如果点A(x1,y1),B(x2,y2)A(x1,y1),B(x2,y2)A(x1,y1),B(x2,y2),则它们的曼哈顿距离为d(A,B)=∣x1−x2∣+∣y1−y2∣d(A,B)=|x1-x2|+|y1-y2|d(A,B)=∣x1−x2∣+∣y1−y2∣

clear,clc
x=[2 5 9 12 15];
y=[8 -4 -9 3 -1];
x1=0:0.01:20;
for i=1:length(x1)
s=0;
for j=1:5d=abs(x1(i)-x(j))+abs(0-y(j));s=s+d;
ends1(i)=s;
end
s1;
[m,n]=min(s1)
x0=x1(n)m =45n =901x0 =9

通过分析可以看到，该桥修在x(9,0)x(9,0)x(9,0)处，使得它到5个村庄的总的曼哈尔顿距离最短。最短为45。

内容来自: von Neumann ↩︎

数学建模清风微信公众号的习题答案(挑战篇-蒙特卡罗思想、枚举法和网格搜索法)相关推荐

数学建模清风微信公众号的习题答案(挑战篇-完结)
以下题目是来自微信公众号数学建模清风老师的题目以下是本人结合在微信公众号上学到的知识去做的,如有不正确或不足,欢迎指正! Q15.在本章3.3.5小节介绍sort函数时,我们留下了一个问题:如果存在 ...
数学建模清风微信公众号的习题答案(挑战篇-数据异常值处理)
以下题目是来自微信公众数学建模清风老师的题目 1.0数据异常处理 1.1正态分布的3σ\sigmaσ原则 1.2箱线图识别异常值 2.0 数据异常处理实例 3.0 总结分析以下是个人结合在微信公众号 ...
数学建模清风微信公众号的习题答案(挑战篇2)
以下题目是来自微信公众号数学建模清风老师的题目以下是个人结合在微信公众号上学到的知识去做的,如有不正确或不足之处,欢迎指正! Q6.最近短视频上有一个有趣的街头抽奖游戏,规则如下:摆摊的店家准备了2 ...
数学建模清风微信公众号的习题答案(挑战篇1)
以下题目是来自微信公众号数学建模清风老师的题目以下是个人结合在微信公众号上学到的知识去做的,若有不正确或不足之处,欢迎指正! Q1.A是一个矩阵,A(:)可以按照列方向的顺序(线性索引)将A中每个元 ...
数学建模清风微信公众号的习题答案(提高篇2）
以下题目是来自微信公众号数学建模清风老师的题目以下是个人结合在微信公众号上学到的知识去做的,若有不正确或不足,欢迎指正! Q11.怎样将一个方阵A的主对角线元素重新赋值为0.例如: 解:这个主要考察 ...
数学建模清风微信公众号的习题答案(基础篇-填空题）
以下题目是来自微信公众号数学建模清风老师的题目如有不正确或不足的地方,欢迎指正! Q1.填空题 (1) MATLAB中矩阵的元素应包括在方括号括号中,矩阵的同行元素之间由(,)分隔,行与行之间用(; ...
微信公众号开发入门教程第一篇
微信公众号开发入门教程第一篇关键字:微信公众平台开发作者:方倍工作室在这篇微信公众平台开发教程中,我们假定你已经有了PHP语言程序.MySQL数据库.计算机网络通讯.及HTTP/XML/CSS/ ...
php早起打卡,微信公众号早起打卡挑战应用制作教程
原标题:微信公众号早起打卡挑战应用制作教程最近有很多公众号都推送了"早起打卡挑战"应用,粉丝可支付一定的金额参与早起打卡,若次日打卡则可获得支付金额+奖励金额,若忘记打卡则会被扣 ...
微信公众号排版多少钱一篇？
公众号运营可以拆解为文章撰写和图文排版,公众号推文可以整体来做,也可以拆开来做,当你写好公众号文章以后,再找专业的人士进行排版美化. 微信公众号排版多少钱一篇?今天伯乐网络传媒就来给大家聊聊这个话题. ...

数学建模清风微信公众号的习题答案(挑战篇-蒙特卡罗思想、枚举法和网格搜索法)

以下题目来自微信公众数学建模清风老师的题目

蒙特卡罗模拟计算π\piπ

蒙特卡罗模拟求定积分

蒙特卡罗模拟计算概率

枚举法与网格搜索法

数学建模清风微信公众号的习题答案(挑战篇-蒙特卡罗思想、枚举法和网格搜索法)相关推荐

最新文章

热门文章