Monte Carlo模拟误差分析课程设计.doc

Monte Carlo模拟误差分析课程设计1. 实验目的1.1 学习并掌握MATLAB软件的基本功能和使用。1.2 学习并掌握基于Monte Carlo Method(MCM)分析的不确定度计算方法。1.3 研究Guide to the expression of Uncertainty in Measurement(GUM)法与MCM法的区别与联系和影响因素，自适应MCM方法，基于最短包含区间的MCM法。2. MATLAB软件介绍实验内容2.1 介绍MATLAB软件的基本知识MATLAB名字由MATrix和LABoratory 两词的前三个字母组合而成。20世纪七十年代，时任美国新墨西哥大学计算机科学系主任的Cleve Moller出于减轻学生编程负担的动机，为学生设计了一组调用LINPACK和EISPACK矩阵软件工具包库程序的“通俗易用”的接口，此即用FORTRAN编写的萌芽状态的MATLABMATLAB语言的主要特点(1). 具有丰富的数学功能(2). 具有很好的图视系统(3). 可以直接处理声言和图形文件。(4). 具有若干功能强大的应用工具箱。(5). 使用方便，具有很好的扩张功能。(6). 具有很好的帮助功能演示内容：(1). MATLAB的数值计算功能在“命令行”Command提示窗口中键入：“A=eye(5,5);A=zeros(5,5);A=ones(5,5)”等命令生成各类矩阵；在“命令行”Command提示窗口中键入：“v,d=eig (A)”生成特征矩阵和特征向量；在“命令行”Command提示窗口中键入：“expm(A)”对矩阵A求幂；在“命令行”Command提示窗口中键入：x=1 3 5;y=2 4 6;z=conv(x,y)；显示结果：z = 2 10 28 38 30(2). MATLAB的符号计算功能在“命令行”Command提示窗口中键入：syms a x；f=sin(a*x); df=diff(f,x); dfa=diff(f,a)；Command提示窗口显示结果： df =cos(a*x)*a； dfa =cos(a*x)*x；2.2 MATLAB软件画图特性(1). MATLAB二维绘图命令函数：plot 参数：线型、颜色、多重线、网格和标记、画面窗口分割、其他方式、隐函数的描绘)(2). MATLAB三维画图曲面与网格图命令函数：mesh三维带阴影曲面图：surf三维曲线命令：plot3演示内容：(1). MATLAB的二维绘图功能在命令行Command提示窗口中键入：close all; x=linspace(0, 2*pi, 100); % 100个点的x座标 y=sin(x); % 对应的y座标 plot(x,y); 得到如下的结果：图1在命令行Command提示窗口中键入：“plot(x, sin(x), x, cos(x);”得到如下的结果：图2在命令行Command提示窗口中键入：plot(x, sin(x), co, x, cos(x), g*); 得到如下的结果：图3在命令行Command提示窗口中键入：xlabel(Input Value); % x轴注解 ylabel(Function Value); % y轴注解 title(Two Trigonometric Functions); % 图形标题 legend(y = sin(x),y = cos(x); % 图形注解 grid on; % 显示格线 得到如下的结果：图4(2). MATLAB的多维绘图功能在命令行Command提示窗口中键入：X,Y = meshgrid(-3:0.125:3); % 生成二维网格点Z = peaks(X,Y); % 生成某种内置函数mesh(X,Y,Z);得到如下的结果图5其他的演示功能详见“MATLAB画图文档”3. Monte Carlo模拟误差分析的实验原理在误差分析的过程中，常用的方法是通过测量方程推导出误差传递方程，再通过不确定度的合成公式获得间接测量量的标准不确定度和扩展不确定度(GUM)。在有些场合下，测量方程较难获得，在这种情况下研究误差的特性就需要借助于模拟统计的方式进行计算。Monte Carlo(MCM)法就是较为常用的数学工具，具体原理相见相关资料。此次课程设计中按照实验要求产生的随机数可以模拟测量误差，通过对这些随机数的概率密度分布函数的面积、包络线和概率特征点的求取，可以获得随机误差的标准不确定度(MCM)，并与理论上估计标准不确定度的Bessel公式、极差法作(GUM)比较，完成实验内容。并以此作为基础，分析GUM法与MCM法的区别与联系，影响MCM法的参数，自适应MCM法和基于最短包含区间的MCM法。已知两项误差分量服从正态分布，标准不确定度分别为mV， mV，用统计模拟分析法给出两项误差和的分布(误差分布的统计直方图，合成的标准差，合成的置信概率 P为99.73%的扩展不确定度)。4. Monte Carlo模拟误差分析的实验内容4.1 MCM法与GUM法进行模拟误差分析和不确定度计算(1). 利用MATLAB软件生成0，1区间的均匀分布的随机数；(2). 给出误差分量的随机值：利用MATLAB，由均匀分布随机数生成标准正态分布随机数，误差分量随机数可表示为mV；同理得 mV(3). 求和的随机数：误差和的随机数；(4). 重复以上步骤，得误差和的随机数系列： ；(5). 作误差和的统计直方图：以误差数值为横坐标，以频率为纵坐标作图。作图区间应包含所有数据，按数值将区间等分为组（尽可能大），每组间隔为，记数各区间的随机数的数目，以为底，以为高作第()区间的矩形，最终的组矩形构成误差和的分布直方图，该图包络线线即为实验的误差分布曲线。(6). 以频率为界划定区间，该区间半宽即为测量总误差的置信概率为99.73%的扩展不确定度。(7). 合成的标准不确定度：A. 总误差随机数平均值 B. 各误差随机数的残差 C. 按照Bessel公式估计标准不确定度 实验流程图：图6实验示例程序：tic;clear;clc;close all;bdclose all;%设定参数值%随机信号点数N，均值Mu，标准差Sigma%N=105;Mu=1;Sigma=2;M=N/10;x=0:1:M;x_=1:M;u1=0.005;u2=0.007;%产生两个在(0,1)上服从均匀分布的,种子为0,每一次都相同的随机数X1和X2%rand(state,0);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);grid on;figure(2),axis(0,N,-max(5*Y2),max(5*Y2)plot(Y2);grid on;% hold on% plot(x,0,k);grid on;% plot(x,1,r-);grid on;% plot(x,-1,r-);grid on;% hold on%变换为任意均值和方差的正态分布%Z1=Sigma*Y1+Mu;%作图%高斯随机信号% subplot(2,2,2)% axis(0,N,-6,6)% plot(Z1);grid on;% hold on% plot(x,Mu,k);% plot(x,Mu+Sigma,r-);grid on;% plot(x,Mu-Sigma,r-);grid on;% hold on%正态分布误差1幅度直方图%figure(3)axis(-1,1,0,N)hist(delta1,M);grid on;%正态分布误差2幅度直方图%figure(4)axis(-1,1,0,N)hist(delta2,M);grid on;%合成误差幅度直方图%figure(5)axis(-1,1,0,N)H=hist(delta,M);hist(delta,M);grid on;%画包络线%figure(6)HH=envelope(x_,H);plot(delta_n,HH,b:);grid on;hold on;%计算直方图的面积%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);for ii=1:M/2-1area(ii+1)=area(ii)+s_(ii);end% 计算99.73%的直方图面积for iii=1:M/2; area(iii);if (area(iii)-0.9973*S)=0; breakendendplot(delta_n(M/2-iii+1),delta_n(M/2+iii),H(M/2-iii),H(M/2+iii),ro);grid on; 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;实验需要讨论的问题：(1). N(采样点数)，M(划分的区间数)与直方图的关系(平滑，Y轴的高度)。(2). Bessel公式计算的标准不确定度(GUM)与99.73%直方图面积的扩展不确定度(MCM)两者之间会存在误差，这个误差与哪些因素有关(N,M,N与M的关系)。4.2 自适应MCM法在执行自适应蒙特卡洛方法的基本过程中，蒙特卡洛试验次数不断增加，直至所需要的各种结果达到统计意义上的稳定。如果某结果的两倍标准偏差小于标准不确定度的数值容差时，则认定该数值结果稳定。(1). 基于前一个实验，构建衡量Monte Carlo法和GUM法计算得到标准不确定度差值的函数。(2). 将随机数个数N，分割区间数M分别作为该函数的自变量，定义自变量的取值范围，从而获得相应的函数值。(3). 分别进行三维网格作图和三维曲线作图，通过观察曲线获得最佳的N，M组合。实验示例程序：tic;warning off;a,b=meshgrid(logspace(1,6);for j=1:max(size(a); for jj=1:max(size(b); Result1(j,jj)=shiyan(a(j),b(jj); endendfigure(1),surfc(a,b,Result1);c=logspace(1,6);d=logspace(1,6);for jjj=1:max(size(c); Result2(jjj)=shiyan(c(jjj),d(jjj);endfigure(2),plot3(c,d,Result2);grid on;toc;实验需要讨论的问题：如何根据三维网格曲线和三维曲线获得最佳