首页/文章/ 详情

大数定理的小陷

3年前浏览4508

image.png

 过冷水诚挚邀请你加入Matlab仿真秀官方交流群进行Matlab学习、问题咨询、 Matlab相关资料下载,群号:927550334 

QQ图片20210424105303.png

        MATLAB数值计算中,积分和微分是很重要的组成部分,积分可以梯形法做计算,也可以使用Monte Carlo方法做积分计算。过冷水在具有使用Monte  Carlo方法做计算的时候取犯了一点小毛病,特此和大家分享,希望对大家精进Monte Carlo方法有帮助。

y=[0 0 0 0 0 0 0 0 1.71 3.24 4.59 5.76 6.75 7.56 8.19 8.64 8.91 9 8.91 8.64  8.19 7.56 6.75 5.76 4.59 3.24 1.71 0]';
x =[ -8 -7 -6 -5 -4.5 -4.3 -4.1 -4 -3.7 -3.4 -3.1 -2.8 -2.5 -2.2 -1.9 -1.6 -1.3 -1 -0.7 -0.4 -0.1 0.2 0.5 0.8 1.1 1.4 1.7 2 ]';
f= fit( x, y,'smoothingspline' );
figure1 = figure;
axes1 = axes('Parent',figure1,'Position',[0.13 0.11 0.8 0.8]);
hold(axes1,'on');
plot(x,y,'LineWidth',3,'Color',[1 0 0]);
ylabel('$f(x)$','FontSize',24,'Interpreter','latex');
xlabel('$x$','FontSize',24,'Interpreter','latex');
 xlim(axes1,[-10 6]);
 ylim(axes1,[-1 14]);
hold(axes1,'off');
set(axes1,'FontSize',16);
annotation(figure1,'arrow',[0.13 0.97],[0.11 0.11],'LineWidth',2);
annotation(figure1,'arrow',[0.58 0.58],[0.11 0.95],'LineWidth',2);


    该函数是一个简单函数用符号法、数值积分都可以获得精确解,实际案例函数很复杂,过冷水在此进行了模型简化,是为了说明问题,不要觉得过冷水在“杀鸡牛刀”。

    需要用大数定理求积分


    用大数定理求如图所示函数的具体表达式则为(N=1000)

1=linspace(-8,2,1000);
y1=f(x1);
I1=10*sum(y1)./1000
I1 =
   35.9138

注意我们函数图像其实在(-8~-4)区间为0所以也可以等于:

x2=linspace(-4,2,1000);
y2=f(x2);
I2=6*sum(y2)./1000
I2 =
   35.9806


        从理论上分析没有问题,实际操作 it is error!  。这是过冷水犯的第一个错误。问题是时间需求是要从f(x)不为零的地方开始积分,过冷水为了省事做整个区间的积分,就是抱着对零积分不影响结果。

        实际问题要稍微做做一下变形,要求分段积分。


        当然如果你微积分还会一点点自然而然知道:

        这个时候过冷水就沙雕了,用大数定理的或分别在(-8,-1),(-1,2),取了N(1000)个点,所有就有2N个点,所以应用大数定理就是:

x3=linspace(-8,-1,1000);
x4=linspace(-1,2,1000);
y3=f(x3);y4=f(x4);
I3=10*sum([y3;y4])/(2000)
I3 =
   42.8187

        再一次 it is error!到现在为止过冷水都没想明白为什么是错的,有对该问题感兴趣的可以留言讨论。“快就是慢,省事就是费事”简单问题还得从原理入手才可靠。

x3=linspace(-8,-1,1000);
x4=linspace(-1,2,1000);
y3=f(x3);y4=f(x4);
I4=7*sum(y3)/1000 3*sum(y4)/1000
I4 =
   35.9588

经过简单对比自然而然会发现

        在早期过冷水做计算的时候一直使用的是I3,直到出现非常明显的错误的时候,一步一步做核查才知道是最基础的公式用错了,之前怀疑是自己代码什么地方敲错了。花了一晚上去检查错误最终排查完后非常有开心なるほど!

        再回过头来讲讲为什么I1(35.9138)I2(35.9806)由于(-8~-4)函数值为零,在该区间取点分占一部分权重,(-8~-4)取点多了 (-4~2)取点自然就少一点 平均值就小一点。通常小数点后两位的误差我们是可以接受的。好巧不不巧过冷水在自己的问题计算误差放大的不合适,该错误可以用重要性抽样来解决。理论没问题,实际操作的时候就有问题了,还是不要省事严格按要求来。

        这两处错误让我对积分和大数定理有了深刻的了解,在简单的问题都要“勤实践,多思考”。总能从中学到新知识。本期的过冷水的分享就这么多,希望都能对大家有帮助,读者在实践中有遇到相似的问题,都可以留言分享,大家共同学习进步。

 

图片

        过冷水发表于 仿真秀 平台原创文章,未经授权禁止私自转载,如需转载请需要和作者沟通表明授权声明,未授权文章皆视为侵权行为,必将追责。如果您希望加入Matlab仿真秀官方交流群进行Matlab学习、问题咨询、 Matlab相关资料下载均可加群:927550334。

精品回顾

 matlab绘制农夫过河动态图

分子动力学的原子空间运动轨迹演示编程

过冷水带你用matlab制作演示动画

python批量移动文件&重命名代码分享

过冷水和你分享 matlab读取存储各种文件的方法 文末有独家金曲分享

image.png

冶金OriginMathematicaMATLAB
著作权归作者所有,欢迎分享,未经许可,不得转载
首次发布时间:2021-04-14
最近编辑:3年前
过冷水
博士 | 讲师 讨论号:927550334
获赞 361粉丝 184文章 107课程 11
点赞
收藏
未登录
还没有评论
课程
培训
服务
行家
VIP会员 学习 福利任务 兑换礼品
下载APP
联系我们
帮助与反馈