哈工大误差分析课程设计MonteCarlo.docx
《哈工大误差分析课程设计MonteCarlo.docx》由会员分享，可在线阅读，更多相关《哈工大误差分析课程设计MonteCarlo.docx（38页珍藏版）》请在冰豆网上搜索。
哈工大误差分析课程设计MonteCarlo
MonteCarlo模拟误差分析课程设计
1实验目的
1.1了解MATLAB软件的基本功能和使用
1.2学习不确定度的统计模拟分析方法
1.3研究误差概率密度函数和Bessel公式获得扩展不确定度的方法和影响因素
2实验原理
在误差分析的过程中，常用的方法是通过测量方程推导出误差传递方程，再通过不确定度的合成公式获得间接测量量的标准不确定度和扩展不确定度（GUM）。
在有些场合下，测量方程较难获得，在这种情况下研究误差的特性就需要借助于模拟统计的方式进行计算。
MonteCarlo（MCM）法就是较为常用的数学工具，具体原理相见相关资料。
此次课程设计中按照实验要求产生的随机数可以模拟测量误差，通过对这些随机数的概率密度分布函数的面积、包络线和概率特征点的求取，可以获得随机误差的标准不确定度——（MCM），并与理论上估计标准不确定度的Bessel公式、极差法作——（GUM）比较，完成实验内容。
并以此作为基础，分析GUM法与MCM法的区别与联系，影响MCM法的参数，自适应MCM法和基于最短包含区间的MCM法。
已知两项误差分量服从正态分布，标准不确定度分别为
mV，
mV，用统计模拟分析法给出两项误差和的分布（误差分布的统计直方图，合成的标准差，合成的置信概率P为99.73%的扩展不确定度）。
3实验内容
（1）.利用MATLAB软件生成[0，1]区间的均匀分布的随机数
；
（2）.给出误差分量的随机值：
利用MATLAB，由均匀分布随机数
生成标准正态分布随机数
，误差分量随机数可表示为
mV；
mV
（3）.求和的随机数：
误差和的随机数
；
（4）.重复以上步骤，得误差和的随机数系列：
；
（5）.作误差和的统计直方图：
以误差数值为横坐标，以频率为纵坐标作图。
作图区间应包含所有数据，按数值将区间等分为
组（
尽可能大），每组间隔为
，记数各区间的随机数的数目
，以
为底，以
为高作第
（
）区间的矩形，最终构成误差和的分布直方图，该图包络线线即为实验的误差分布曲线。
（6）.以频率
为界划定区间，该区间半宽即为测量总误差的置信概率为95%的扩展不确定度。
（7）.合成的标准不确定度：
4.实验流程图：
一．实验1
本实验中随机数种子为014。
并使分别取N为100000点和10000点两种情况下，得到M值分别为5*N,2*N,N,N/2,N/5,N/10五种情况下的模拟图像。
1．实验1程序
tic;
clear;clc;closeall;
%%设定参数值%%
%%随机信号点数N，均值为1，标准差u1,u2%%
N=10^5;
M=N/10;
x=0:
1:
M;
x_=[1:
M];
u1=0.005;
u2=0.007;
%%产生两个在（0,1）上服从均匀分布的,种子为0,每一次都相同的随机数X1和X2%%
rand（'state',014）;
X1=rand（1,N）;
X2=rand（1,N）;
%%按照Box-Mueller变换方法产生标准正态分布Y1和Y2%%
Y1=sqrt（-2*log（X1））.*cos（2*pi*X2）;
Y2=sqrt（-2*log（X1））.*sin（2*pi*X2）;
%%为做直方图先定义好X轴的坐标数据%%
delta1=u1*Y1;
delta2=u2*Y2;
delta=delta1+delta2;
d_delta=（max（delta）-min（delta））/（M-1）;%%d_delta为误差分布的间距
delta_n=[min（delta）:
d_delta:
max（delta）];%%delta_n为误差分布序列
%%作图%%
%%高斯随机信号%%
figure
（1）,
axis（[0,N,-max（5*Y1）,max（5*Y1）]）
plot（Y1）;gridon;
figure
（2）,
axis（[0,N,-max（5*Y2）,max（5*Y2）]）
plot（Y2）;gridon;
%holdon
%plot（x,0,'k'）;gridon;
%plot（x,1,'r--'）;gridon;
%plot（x,-1,'r--'）;gridon;
%holdon
%%变换为任意均值和方差的正态分布%%
%Z1=Sigma*Y1+Mu;
%%作图%%
%%高斯随机信号%%
%subplot（2,2,2）
%axis（[0,N,-6,6]）
%plot（Z1）;gridon;
%holdon
%plot（x,Mu,'k'）;
%plot（x,Mu+Sigma,'r--'）;gridon;
%plot（x,Mu-Sigma,'r--'）;gridon;
%holdon
%%正态分布误差1幅度直方图%%
figure（3）
axis（[-1,1,0,N]）
hist（delta1,M）;gridon;
%%正态分布误差2幅度直方图%%
figure（4）
axis（[-1,1,0,N]）
hist（delta2,M）;gridon;
%%合成误差幅度直方图%%
figure（5）
axis（[-1,1,0,N]）
H=hist（delta,M）;
hist（delta,M）;gridon;
%%画包络线%%
figure（6）
HH=envelope（x_,H）;
plot（delta_n,HH,'b:
'）;gridon;
holdon;
%%计算直方图的面积%%
S=sum（d_delta*abs（H））;
%%计算直方图的面积%%
%%s_1表示正向直方图的每一个单元的面积
%%s_2表示反向直方图的每一个单元的面积
%%s_表示正反向两两对称每一对单元的面积
%%area表示以中心为对称轴的累加面积
i=1:
1:
M/2;
s_1（i）=d_delta*abs（floor（H（floor（M/2+i））））;
s_2（i）=d_delta*abs（floor（H（floor（M/2-i+1））））;
s_（i）=s_1（i）+s_2（i）;
area
（1）=s_
（1）;
forii=1:
M/2-1
area（ii+1）=area（ii）+s_（ii）;
end
%%计算99.73%的直方图面积
foriii=1:
M/2;
area（iii）;
if（area（iii）-0.9973*S）>=0;
break
end
end
plot（[delta_n（M/2-iii+1）,delta_n（M/2+iii）],[H（M/2-iii+1）,H（M/2+iii）],'ro'）;gridon;
delta_n_u=（delta_n（floor（M/2+iii））-delta_n（floor（M/2-iii+1）））/6;
%%理论计算标准不确定度%%
delta_mean=mean（delta）;
delta_cancha=delta-delta_mean;
s=sqrt（（sum（delta_cancha.^2））/（N-1））;
%%%%%%%%%%%%%%%
toc;
2.实验1程序运行结果图
（1）当M=N/10时
Figure1
Figure2
Figure3
Figure4
Figure5
Figure6
（2）当更改N与M不同的倍数关系时，可得到不同的计算结果，如以下个图所示：
图1.1N=10^5,M=N*5，s=0.0086，detla_n_u=0.0087
图1.2N=10^5,M=N*2，s=0.0086，detla_n_u=0.0087
图1.3N=10^5,M=N，s=0.0086，detla_n_u=0.0087
图1.4N=10^5,M=N/2，s=0.0086，detla_n_u=0.0087
图1.5N=10^5,M=N/5，s=0.0086，detla_n_u=0.0086
图1.6N=10^5,M=N/10，s=0.0086，detla_n_u=0.0085
图1.7N=10^4,M=N*5，s=0.0085，detla_n_u=0.0087
图1.8N=10^4,M=N*2，s=0.0085，detla_n_u=0.0087
图1.9N=10^4,M=N，s=0.0085，detla_n_u=0.0084
图1.10N=10^4,M=N/2，s=0.0085，detla_n_u=0.0082
图1.11N=10^4,M=N/5，s=0.0085，detla_n_u=0.0078
图1.12N=10^4,M=N/10，s=0.0085，detla_n_u=0.0074
表2N=10^5时，N与M成不同倍数k时，直方图计算结果与理论计算结果的差异
k=N/M
1/5
1/2
1
2
5
10
s
0.0086
0.0086
0.0086
0.0086
0.0086
0.0086
delta_n_u
0.0087
0.0087
0.0087
0.0087
0.0086
0.0085
|delta_n_u﹣s|
0.0001
0.0001
0.0001
0.0001
0
0.0001
表2N=10^4时，N与M成不同倍数k时，直方图计算结果与理论计算结果的差异
k=N/M
1/5
1/2
1
2
5
10
s
0.0085
0.0085
0.0085
0.0085
0.0085
0.0085
delta_n_u
0.0087
0.0087
0.0084
0.0082
0.0078
0.0074
|delta_n_u﹣s|
0.0002
0.0002
0.0001
0.0003
0.0007
0.0011
3实验需要讨论的问题
（1）.N（采样点数），M（划分的区间数）与直方图的关系（平滑，Y轴的高度）。
有图1.1~1.12可知：
当N固定的情况下，随着M值得增大直方图的平滑性变差，Y轴高度下降。
其中，M当N改变时，即当N增大时可使直方图更为精细，且此时不改变直方图包络的基本形状。
（2）.Bessel公式计算的标准不确定度与99.73%直方图面积的扩展不确定度两者之间会存在误差，这个误差与哪些因素有关（N,M,N>=M）
此误差的大小和M、N的相对大小值有关。
当N>=M时，由于对离散的误差值统计运算存在舍入误差导致误差，此误差随着M的增大可消除此项舍入误差。
当M>N时，增大M值，误差值稳定，且不能改善误差值。
二．实验2——自适应MCM法
在执行自适应蒙特卡洛方法的基本过程中，蒙特卡洛试验次数不断增加，直至所需要的各种结果达到统计意义上的稳定。
如果某结果的两倍标准偏差小于标准不确定度的数值容差时，则认定该数值结果稳定。
（1）.基于前一个实验，构建衡量MonteCarlo法和GUM法计算得到标准不确定度差值的函数。
（2）.将随机数个数N，分割区间数M分别作为该函数的自变量，定义自变量的取值范围，从而获得相应的函数值。
（3）.分别进行三维网格作图和三维曲线作图，通过观察曲线获得最佳的N，
升级会员 