matlab的数值积分

忘是亡心i 2022-05-29 01:20 251阅读 0赞

![Image 1][]

matlab实现标准正太分布表：

function   zhengtaifenbu                                         
    syms t x a                         
    fun=exp(-t^2/2)/sqrt(2*pi);         %定义被积函数
    Fun=int(fun,t,-inf,t);              %从-∞到t积分
    fprintf('   %c ','x');               
    for a=0.00:0.01:0.09             
     fprintf('%.2f   ',a);         
    end
     
     fprintf(' ');
    fprintf('\n');                     
    for x=0.0:0.1:3.0
        fprintf(' ');
        fprintf('%.1f ',x);            %输出除首行首列外的x的原始值
       
       for t=x+(0.00:0.01:0.09);        
       fprintf('%.4f ',eval(Fun));      %以4位小数格式输出标准正态分布函数的每个确定点处的值
       end
    fprintf('\n');   
    end

梯形公式：

function res = TiXingint(a,b,n)
        h=(b-a)/n;
        x=a:h:b;
        res=0.0;
        y=exp(-x.^2/2)/sqrt(2*pi);
    for j=1:n
        res=res+y(j)+y(j+1);
    end
        res=res*h/2;
    end

辛普森公式：

function res = Simpson(a,b,n)
        h=(b-a)/n;
        x=a:h:b;
        res=0.0;
        y=exp(-x.^2/2)/sqrt(2*pi);
        y(1) = 0;
        z=exp(-(x+h/2).^2/2)/sqrt(2*pi);
        for j=1:n
            res=res + y(j) + y(j+1) + 4*z(j); 
        end
        res=res*h/6;
    end

牛顿科特斯公式：

function res=NewtonCotes(a,b,n)
        h=(b-a)/n;
        x=a:h:b;
        res=0.0;
        y=exp(-x.^2/2)/sqrt(2*pi);
        y(1)=0;
        z0=exp(-(x+h/4).^2/2)/sqrt(2*pi);
        z1=exp(-(x+h/2).^2/2)/sqrt(2*pi);
        z2=exp(-(x+3*h/4).^2/2)/sqrt(2*pi);
        for j=1:n
            res=res+7*y(j)+7*y(j+1)+32*z0(j)+12*z1(j)+32*z2(j);
        end
        res=res*h./90;
    end

实验总结：

根据实验数据可以看出，科特斯公式的结果更为接近标准正态分布表。其中，梯形公式n=1划分，辛普森公式n=2划分，科特斯公式n=4划分，也就是说(a,b)划分为n等分，n越大，求解的实验结果更为准确。

[Image 1]: