精华内容
下载资源
问答
  • matlab数值积分
    2020-11-06 22:16:33

    求定积分
    ∫ 1 0 x x 2 + 4 \int_{1}^{0} \frac{x}{x^{2}+4} 10x2+4x

    fun=@(x) x./(x.^2+4);
    q=quadl(fun,0,1)
    q =
        0.1116
    

    求二重定积分
    ∫ 0 1 ∫ − 1 1 y 2 sin ⁡ x d x d y \int_{0}^{1}\int_{-1}^{1} y^{2}\sin xdxdy 0111y2sinxdxdy

    fun=@(x,y) y.^2.*sin(x);
    q=dblquad(fun,-1,1,0,1);
    q =
      -1.4626e-17
    

    在数值积分情况下,此值约等于0。另外,不推荐使用 quad,改用 integral。

    更多相关内容
  • 分数阶matlab工具箱,matlab数值积分方法,matlab源码.zip,fotf,c8mnlf2.slx,optimPID,pidctrl_modelR2011a.mdl,multi_step.m,models,mod_2R2011b.mdl,mod_1R2010a.mdl,mod_1R2011b.mdl,mod_ltiR2011b.mdl,mod_3R2011a...
  • 第8章 MATLAB数值积分与微分.ppt
  • 4IntSimpson 用辛普森系列公式求积分NewtonCotes 用牛顿-科茨系列公式求积分IntGauss 用高斯公式求积分IntGaussLada 用高斯拉道公式求积分IntGaussLobato 用高斯—洛巴托公式求积分IntSample 用三次样条插值求积分...
  • 学号 班级 统计1001 姓名 指导教师 易昆南 实验题目 用多种方法计算数值积分 评 分 1设计实习目的 了解MATLAB在实际问题中的应用 通过实践加深对这门语言中M文件的了解 熟悉简单程序结构如循环结构for循环while循环...
  • Matlab非常经典的新手教程进阶提高-第8章 MATLAB数值积分与微分.ppt 对于新手,是由易而难的好教程,对于有一定基础的后面章节也有很好的提高,是个全面的ppt教程,希望能帮到大家!
  • 实验10 数值积分 实验目的 1了解数值积分的基本原理 2熟练掌握数值积分MATLAB 实现 3会用数值积分方法解决一些实际问题 实验内容 积分是数学中的一个基本概念在实际问题中也有很广泛的应用同微分一样在微 积分中它...
  • matlab数值积分以及数值求导函数源码
  • MATLAB数值积分求值实验报告范文数值积分实验报告范文.pdfMATLAB数值积分求值实验报告范文数值积分实验报告范文.pdfMATLAB数值积分求值实验报告范文数值积分实验报告范文.pdfMATLAB数值积分求值实验报告范文数值积分...
  • MATLAB数值积分求值实验报告范文数值积分实验报告范文.docxMATLAB数值积分求值实验报告范文数值积分实验报告范文.docxMATLAB数值积分求值实验报告范文数值积分实验报告范文.docxMATLAB数值积分求值实验报告范文数值...
  • M a t l a b 数 值 积 分 与 数 值 微 分 Matlab 数值积分 1. 一重数值积分的实现方法 变步长辛普森法高斯 - 克朗罗德法梯形积分法 1.1 变步长辛普森法 Matlab 提供了 quad 函数和 quadl 函数用于实现变步长 辛普森...
  • 4IntSimpson 用辛普森系列公式求积分NewtonCotes 用牛顿-科茨系列公式求积分IntGauss 用高斯公式求积分IntGaussLada 用高斯拉道公式求积分IntGaussLobato 用高斯—洛巴托公式求积分IntSample 用三次样条插值求积分...
  • MATLAB程序设计教程 MATLAB与高等数学 第08章 MATLAB数值积分与微分(共14页).ppt MATLAB程序设计教程 MATLAB与高等数学 第09章 MATLAB符号计算(共23页).ppt MATLAB程序设计教程 MATLAB与高等数学 第10章 MATLAB...
  • Matlab数值积分问题 求和命令sum调用格式. 如果x是向量,则sum(x) 给出x的各个元素的累加和;如果x是矩阵,则sum(x)是一个元素为x的每列列和的行向量. 例3.1 调用命令sum 求向量x的各个元素的累加和 解:输入 x=[1,2,3...
  • matlab数值积分.doc

    2021-09-06 19:22:23
    matlab代码,数值积分,一维积分,二维积分,三维积分
  • Matlab数值计算方法程序源代码 Matlab数值积分 共44页.pdf
  • 杨贤邦 学 号 指导教师 吴明芬 实验时间 2014.4.16 地 点 综合实验大楼 203 实验题目 数值积分方法 实 利用复化梯形辛普森公式和龙贝格数值积分公式计算定积分的 验 目 近似植 的 梯形辛普森柯特斯法及其Matlab 实现...
  • matlab数值积分的计算

    千次阅读 2018-10-03 11:00:36
    文章目录标签(空格分隔): matlab 积分 数值 计算@[toc]matlab数值积分1 Gauss-Hermite积分1.1 测试Gauss-Hermite积分函数gaussHermiteIntegral()1.2 gaussHermiteIntegral()2 Gauss-Laguerre积分2.1 测试Gauss-...

    标签(空格分隔): matlab 积分 数值 计算

    matlab数值积分

    1 Gauss-Hermite积分

    1.1 测试Gauss-Hermite积分函数gaussHermiteIntegral()

    clc,clear
    format long
    f=@(x)(exp(-x.^4))
    % % % % f=@(x)((1/sqrt(2*pi))*exp(-x.^2/2)) %标准正太分布
    % % area=[];nn=[];
    % % for n=1:10
    % %     nn=[nn,n];
    % %     area=[area,gaussHermiteIntegral(f,-inf,inf,n)];
    % % end,nn,area
    % % darea=abs(diff(area));
    % % %% 系统计算的积分精确值
    % % f_=@(x)(f(-x))
    % % area_system=quadgk(f,0,inf,'RelTol',1e-8,'AbsTol',1e-12)+...
    % %             quadgk(f_,0,inf,'RelTol',1e-8,'AbsTol',1e-12)
    % % d_my_system=abs(area-area_system);
    % % log_darea=round(log10(darea)*10)
    % % log_d_my_system=round(log10(d_my_system)*10)
    %    1.812804954110955
    n=5
    [area,gErrer]=gaussHermiteIntegral(f,-inf,inf,n,true)
    %% 从以下数据看(前向差商,系统误差):前向差商公式完全能够作为事后误差估计
    % log_darea =   -19   -28   -36   -45   -52   -60   -69   -90   -81
    % log_d_my_system = -19   -29   -37   -45   -53   -61   -69   -82   -83   -86
    %% 从以下数据看:Gauss Hermite积分的精确度e随积分节点个数n呈现关系:log10(e)=n/10
            %因此采用10n个节点的gauss积分作为程序中的误差阶数
    % f =  @(x)(exp(-x.^4))
    % n =100
    % area = 1.812804956399165
    % area_system = 1.812804954110955
    % d_my_system = 2.288210509959754e-09
    % n=41
    % d_my_system = 2.403635556458283e-05
    % n=20
    % d_my_system =   0.001299981246872
    % n=10
    % d_my_system =  0.014102477783775
    

    运行结果:

    f = 
        @(x)(exp(-x.^4))
    n =
         5
    area =
       1.812809811622279
    gErrer =
         3.459687060525241e-05
    

    1.2 gaussHermiteIntegral()

    function [area,gErrer]=gaussHermiteIntegral(func,a,b,n,withError)
    %高斯积分
    %理论依据:数值分析方法 奚梅成
    %func 被积函数;
    %[a,b] 积分区间
    %n 高斯积分的阶数
    %%
    %name:邓能财     Date: 2013/12/23
    %% 默认参数
    if nargin<5 withError=false; gErrer=0; end
    % 查错
    assert(a==-inf && b==inf,'积分区间必须为(-inf,inf)!');
    assert(ismember(n,[1:10]),'积分精确度数量级必须在1~10之间!');
    %% 函数变换
    expfunc=@(x)(exp(x^2)*func(x));
    %% 数据:10,20,30,...,maxn*10阶高斯积分的积分系数
    maxn=10; %   (±积分节点node,积分系数coef),‘±’被省略
    data= struct('node',cell(1,maxn),'coef',cell(1,maxn));
    %data()= struct('node',{[]},'coef',{[]});
    data(1)= struct('node',{[3.4361591188377e+00   
    2.5327316742328e+00   
    1.7566836492999e+00   
    1.0366108297895e+00   
    3.4290132722370e-01   ]},'coef',{[7.6404328552326e-06   
    1.3436457467812e-03   
    3.3874394455481e-02   
    2.4013861108231e-01   
    6.1086263373533e-01   ]});
    data(2)= struct('node',{[5.3874808900112e+00   
    4.6036824495507e+00   
    3.9447640401156e+00   
    3.3478545673832e+00   
    2.7888060584281e+00   
    2.2549740020893e+00   
    1.7385377121166e+00   
    1.2340762153953e+00   
    7.3747372854539e-01   
    2.4534070830090e-01   ]},'coef',{[2.2293936455341e-13   
    4.3993409922731e-10   
    1.0860693707693e-07   
    7.8025564785321e-06   
    2.2833863601635e-04   
    3.2437733422379e-03   
    2.4810520887464e-02   
    1.0901720602002e-01   
    2.8667550536283e-01   
    4.6224366960061e-01   ]});
    data(3)= struct('node',{[6.8633452935299e+00   
    6.1382792201239e+00   
    5.5331471515675e+00   
    4.9889189685899e+00   
    4.4830553570925e+00   
    4.0039086038612e+00   
    3.5444438731553e+00   
    3.0999705295864e+00   
    2.6671321245356e+00   
    2.2433914677615e+00   
    1.8267411436037e+00   
    1.4155278001982e+00   
    1.0083382710467e+00   
    6.0392105862555e-01   
    2.0112857654887e-01   ]},'coef',{[2.9082547163104e-21   
    2.8103336037565e-17   
    2.8786070806795e-14   
    8.1061862974598e-12   
    9.1785804243669e-10   
    5.1085224507758e-08   
    1.5790948873248e-06   
    2.9387252289230e-05   
    3.4831012431868e-04   
    2.7379224730677e-03   
    1.4703829704827e-02   
    5.5144176870234e-02   
    1.4673584754089e-01   
    2.8013093083921e-01   
    3.8639488954181e-01   ]});
    data(4)= struct('node',{[8.0987611392509e+00   
    7.4115825314855e+00   
    6.8402373052494e+00   
    6.3282553512201e+00   
    5.8540950560304e+00   
    5.4066542479701e+00   
    4.9792609785453e+00   
    4.5675020728444e+00   
    4.1682570668325e+00   
    3.7792067534352e+00   
    3.3985582658596e+00   
    3.0248798839013e+00   
    2.6569959984429e+00   
    2.2939171418751e+00   
    1.9347914722823e+00   
    1.5788698949316e+00   
    1.2254801090463e+00   
    8.7400661235709e-01   
    5.2387471383228e-01   
    1.7453721459758e-01   ]},'coef',{[2.5910578133142e-29   
    8.5440614848689e-25   
    2.5675935123395e-21   
    1.9891810183765e-18   
    6.0083587891597e-16   
    8.8057076445574e-14   
    7.1565280526724e-12   
    3.5256207913729e-10   
    1.1212360832276e-08   
    2.4111441636703e-07   
    3.6315761506931e-06   
    3.9369339810925e-05   
    3.1385359454133e-04   
    1.8714968295980e-03   
    8.4608880082581e-03   
    2.9312565536172e-02   
    7.8474605865404e-02   
    1.6337873271327e-01   
    2.6572825187738e-01   
    3.3864327742559e-01   ]});
    data(5)= struct('node',{[9.1824069581293e+00   
    8.5227710309178e+00   
    7.9756223682056e+00   
    7.4864094298642e+00   
    7.0343235097706e+00   
    6.6086479738554e+00   
    6.2029525192747e+00   
    5.8129946754204e+00   
    5.4357860872249e+00   
    5.0691175849172e+00   
    4.7112936661690e+00   
    4.3609731604546e+00   
    4.0170681728581e+00   
    3.6786770625153e+00   
    3.3450383139379e+00   
    3.0154977695745e+00   
    2.6894847022677e+00   
    2.3664939042987e+00   
    2.0460719686864e+00   
    1.7278065475159e+00   
    1.4113177548983e+00   
    1.0962511289577e+00   
    7.8227172955461e-01   
    4.6905905667824e-01   
    1.5630254688947e-01   ]},'coef',{[1.8283576319454e-37   
    1.6732639420678e-32   
    1.2152047767042e-28   
    2.1376568843541e-25   
    1.4170933459735e-22   
    4.4709843726175e-20   
    7.7423829723092e-18   
    8.0942618938896e-16   
    5.4659440315878e-14   
    2.5066555238998e-12   
    8.1118773649361e-11   
    1.9090405438119e-09   
    3.3467934040213e-08   
    4.4570299668179e-07   
    4.5816827079555e-06   
    3.6840190537807e-05   
    2.3426989210925e-04   
    1.1890117817496e-03   
    4.8532638261719e-03   
    1.6031941068412e-02   
    4.3079159156766e-02   
    9.4548935477086e-02   
    1.7003245567716e-01   
    2.5113085633200e-01   
    3.0508512920440e-01   ]});
    data(6)= struct('node',{[1.0159109246180e+01   
    9.5209036770133e+00   
    8.9923980014049e+00   
    8.5205692841176e+00   
    8.0851886542490e+00   
    7.6758399375049e+00   
    7.2862765943956e+00   
    6.9123815321893e+00   
    6.5512591670629e+00   
    6.2007735579934e+00   
    5.8592901963942e+00   
    5.5255210861387e+00   
    5.1984265345763e+00   
    4.8771500774732e+00   
    4.5609737579358e+00   
    4.2492864359560e+00   
    3.9415607339262e+00   
    3.6373358761707e+00   
    3.3362046535476e+00   
    3.0378033382307e+00   
    2.7418037480697e+00   
    2.4479069023077e+00   
    2.1558378712292e+00   
    1.8653415312330e+00   
    1.5761790119750e+00   
    1.2881246748689e+00   
    1.0009634995607e+00   
    7.1448878167258e-01   
    4.2850006422063e-01   
    1.4280123870344e-01   ]},'coef',{[3.7989555024173e-45   
    3.2023914909894e-40   
    3.9179142693654e-36   
    1.3335820885372e-32   
    1.7163780212922e-29   
    1.0294248187970e-26   
    3.3457568162026e-24   
    6.5125663921413e-22   
    8.1536403865806e-20   
    6.9232479198565e-18   
    4.1524441111527e-16   
    1.8166245761912e-14   
    5.9484305160636e-13   
    1.4889573490645e-11   
    2.8993590128044e-10   
    4.4568227752248e-09   
    5.4755546192770e-08   
    5.4335161342049e-07   
    4.3942869362670e-06   
    2.9187419041555e-05   
    1.6027733468185e-04   
    7.3177355696551e-04   
    2.7913248289531e-03   
    8.9321783603078e-03   
    2.4061272766109e-02   
    5.4718970932183e-02   
    1.0529876369779e-01   
    1.7177615691888e-01   
    2.3786890495866e-01   
    2.7985311752283e-01   ]});
    data(7)= struct('node',{[1.1055240743138e+01   
    1.0434459776321e+01   
    9.9210064725726e+00   
    9.4631239608462e+00   
    9.0410623133341e+00   
    8.6446518655005e+00   
    8.2677951152318e+00   
    7.9064749419688e+00   
    7.5578680649896e+00   
    7.2198940023019e+00   
    6.8909628903904e+00   
    6.5698241787567e+00   
    6.2554707188783e+00   
    5.9470752244714e+00   
    5.6439466223452e+00   
    5.3454991415733e+00   
    5.0512298489716e+00   
    4.7607019531714e+00   
    4.4735321501382e+00   
    4.1893808634926e+00   
    3.9079445988681e+00   
    3.6289498686373e+00   
    3.3521483007724e+00   
    3.0773126524498e+00   
    2.8042335229303e+00   
    2.5327166122907e+00   
    2.2625804097899e+00   
    1.9936542226231e+00   
    1.7257764756010e+00   
    1.4587932269484e+00   
    1.1925568563654e+00   
    9.2692488970662e-01   
    6.6175893081257e-01   
    3.9692367564299e-01   
    1.3228598727032e-01   ]},'coef',{[5.9928891485445e-50   
    1.3088219454484e-44   
    1.9775381832792e-41   
    5.1459796204458e-40   
    1.2902253872744e-36   
    1.3859983185653e-33   
    7.5869522154891e-31   
    2.5172152573746e-28   
    5.3429813146915e-26   
    7.6615884430059e-24   
    7.7454741231880e-22   
    5.7095124891098e-20   
    3.1533210598112e-18   
    1.3342434683070e-16   
    4.4061175397096e-15   
    1.1534952560119e-13   
    2.4259623022220e-12   
    4.1457910332392e-11   
    5.8137344484274e-10   
    6.7473312881327e-09   
    6.5292934189487e-08   
    5.3024856335052e-07   
    3.6344904047183e-06   
    2.1131233267306e-05   
    1.0466994285007e-04   
    4.4339999663333e-04   
    1.6117353953559e-03   
    5.0416367476379e-03   
    1.3605327560619e-02   
    3.1741278190051e-02   
    6.4133640985369e-02   
    1.1238816501982e-01   
    1.7101015293604e-01   
    2.2612844383734e-01   
    2.5999310620316e-01   ]});
    data(8)= struct('node',{[1.1887863560471e+01   
    1.1281694164989e+01   
    1.0780796472469e+01   
    1.0334491039538e+01   
    9.9234351145692e+00   
    9.5376679163021e+00   
    9.1712174787265e+00   
    8.8201501286611e+00   
    8.4817021873216e+00   
    8.1538381574635e+00   
    7.8350036783293e+00   
    7.5239772908069e+00   
    7.2197764697589e+00   
    6.9215953729258e+00   
    6.6287620808576e+00   
    6.3407083213300e+00   
    6.0569474734510e+00   
    5.7770582280616e+00   
    5.5006722123152e+00   
    5.2274644551007e+00   
    4.9571459285134e+00   
    4.6894576329528e+00   
    4.4241658477652e+00   
    4.1610582741258e+00   
    3.8999408693874e+00   
    3.6406352232272e+00   
    3.3829763625124e+00   
    3.1268108983731e+00   
    2.8719954485287e+00   
    2.6183952824718e+00   
    2.3658831480742e+00   
    2.1143382465070e+00   
    1.8636453287449e+00   
    1.6136938918511e+00   
    1.3643774570540e+00   
    1.1155929146033e+00   
    8.6723992270397e-01   
    6.1922034962757e-01   
    3.7143774948304e-01   
    1.2379686317313e-01   ]},'coef',{[3.8685383727931e-49   
    2.0252403350315e-45   
    5.6114941460751e-43   
    2.7312591323502e-40   
    8.0059947621974e-39   
    4.8548677202834e-38   
    1.6512006592794e-36   
    6.3922427517268e-35   
    1.8964852968049e-32   
    4.3114531467056e-30   
    6.8829012692744e-28   
    7.9853336688990e-26   
    6.9341381005532e-24   
    4.6134698831073e-22   
    2.3978402197188e-20   
    9.8965290426435e-19   
    3.2891620153310e-17   
    8.9094610374559e-16   
    1.9875208348038e-14   
    3.6848498201918e-13   
    5.7232797184518e-12   
    7.4997199146582e-11   
    8.3430559260916e-10   
    7.9229232891485e-09   
    6.4544629501197e-08   
    4.5305355658615e-07   
    2.7507246274012e-06   
    1.4496444191130e-05   
    6.6517412418152e-05   
    2.6647788669334e-04   
    9.3431768696117e-04   
    2.8731975330897e-03   
    7.7640394053997e-03   
    1.8465751154428e-02   
    3.8708599914838e-02   
    7.1600747691843e-02   
    1.1698073380430e-01   
    1.6893796927590e-01   
    2.1577494627433e-01   
    2.4383585380721e-01   ]});
    data(9)= struct('node',{[1.2668764375434e+01   
    1.2075130816911e+01   
    1.1584960252805e+01   
    1.1148508777694e+01   
    1.0746787108514e+01   
    1.0370015493271e+01   
    1.0012330403392e+01   
    9.6698698743375e+00   
    9.3399210724884e+00   
    9.0204865095226e+00   
    8.7100414635652e+00   
    8.4073884130035e+00   
    8.1115647551405e+00   
    7.8217816701739e+00   
    7.5373821269924e+00   
    7.2578111510466e+00   
    6.9825942253569e+00   
    6.7113212484547e+00   
    6.4436343875168e+00   
    6.1792187234935e+00   
    5.9177949371442e+00   
    5.6591135131207e+00   
    5.4029500908349e+00   
    5.1491016937770e+00   
    4.8973836402140e+00   
    4.6476269884214e+00   
    4.3996764055654e+00   
    4.1533883754886e+00   
    3.9086296798929e+00   
    3.6652761017604e+00   
    3.4232113106648e+00   
    3.1823258978579e+00   
    2.9425165353484e+00   
    2.7036852380959e+00   
    2.4657387122777e+00   
    2.2285877756048e+00   
    1.9921468380504e+00   
    1.7563334332534e+00   
    1.5210677923813e+00   
    1.2862724534456e+00   
    1.0518719000410e+00   
    8.1779222425490e-01   
    5.8396080911088e-01   
    3.5030602639582e-01   
    1.1675694609164e-01   ]},'coef',{[1.3926602246895e-51   
    1.9568360107978e-47   
    1.0639144183847e-44   
    7.9089094484501e-43   
    7.1933465481042e-41   
    3.4996773126588e-39   
    9.2942987391649e-38   
    5.0159398402133e-38   
    3.1216254599340e-37   
    2.8172065374851e-36   
    3.9144464842405e-34   
    5.9935373474292e-32   
    7.7735858553266e-30   
    7.7182478105170e-28   
    5.9828467035793e-26   
    3.6826786261206e-24   
    1.8269175754861e-22   
    7.3963427902401e-21   
    2.4705506323895e-19   
    6.8737797964624e-18   
    1.6064742921236e-16   
    3.1773129595293e-15   
    5.3533724860157e-14   
    7.7293031086162e-13   
    9.6137828692418e-12   
    1.0350233726258e-10   
    9.6863849248741e-10   
    7.9104166149095e-09   
    5.6567852900933e-08   
    3.5533045829498e-07   
    1.9661408494659e-06   
    9.6077760741792e-06   
    4.1557921102543e-05   
    1.5944214368481e-04   
    5.4359224952963e-04   
    1.6496055361812e-03   
    4.4622998595035e-03   
    1.0773786624642e-02   
    2.3243228173755e-02   
    4.4849988132390e-02   
    7.7467627605513e-02   
    1.1985655590239e-01   
    1.6619497260996e-01   
    2.0661430370392e-01   
    2.3035797018195e-01   ]});
    data(10)= struct('node',{[1.3406487338145e+01   
    1.2823799749488e+01   
    1.2342964222860e+01   
    1.1915061943114e+01   
    1.1521415400787e+01   
    1.1152404385585e+01   
    1.0802260753685e+01   
    1.0467185421343e+01   
    1.0144509941293e+01   
    9.8322698077780e+00   
    9.5289658233901e+00   
    9.2334208902192e+00   
    8.9446892173255e+00   
    8.6619961681345e+00   
    8.3846969404163e+00   
    8.1122473111628e+00   
    7.8441823844608e+00   
    7.5801008078575e+00   
    7.3196528223045e+00   
    7.0625310602489e+00   
    6.8084633528588e+00   
    6.5572070319215e+00   
    6.3085443611121e+00   
    6.0622788326143e+00   
    5.8182321352035e+00   
    5.5762416493299e+00   
    5.3361583601384e+00   
    5.0978451050891e+00   
    4.8611750917912e+00   
    4.6260306357872e+00   
    4.3923020786827e+00   
    4.1598868551310e+00   
    3.9286886834277e+00   
    3.6986168593185e+00   
    3.4695856364186e+00   
    3.2415136796310e+00   
    3.0143235803312e+00   
    2.7879414239820e+00   
    2.5622964023726e+00   
    2.3373204639069e+00   
    2.1129479963712e+00   
    1.8891155374270e+00   
    1.6657615087415e+00   
    1.4428259702159e+00   
    1.2202503912190e+00   
    9.9797743609810e-01   
    7.7595076154015e-01   
    5.5411482359162e-01   
    3.3241469234223e-01   
    1.1079587242244e-01   ]},'coef',{[4.1890734835153e-53   
    6.6685952840498e-49   
    1.7149047609884e-45   
    4.6329593470698e-43   
    2.3803897130070e-41   
    3.1234305446425e-39   
    2.4962135456919e-38   
    2.5988685504385e-37   
    1.3109926053898e-36   
    1.4501664855704e-36   
    5.1295692004930e-37   
    1.3477869587572e-35   
    2.8227154103471e-35   
    5.9582141086172e-34   
    7.7461504141821e-32   
    7.1102466014413e-30   
    5.0387971358429e-28   
    2.9171459549279e-26   
    1.3948232985713e-24   
    5.5610285249618e-23   
    1.8649975920043e-21   
    5.3023160356162e-20   
    1.2868329542238e-18   
    2.6824921883091e-17   
    4.8298353231531e-16   
    7.5488968761282e-15   
    1.0288749372616e-13   
    1.2278785143771e-12   
    1.2879038257424e-11   
    1.1913006349369e-10   
    9.7479212538641e-10   
    7.0758572838838e-09   
    4.5681275084851e-08   
    2.6290974837540e-07   
    1.3517971591104e-06   
    6.2215248177779e-06   
    2.5676159384549e-05   
    9.5171627785510e-05   
    3.1729197104330e-04   
    9.5269218854862e-04   
    2.5792732600591e-03   
    6.3030002856081e-03   
    1.3915665220232e-02   
    2.7779127385933e-02   
    5.0175812677429e-02   
    8.2051827391224e-02   
    1.2153798684410e-01   
    1.6313003050278e-01   
    1.9846285025419e-01   
    2.1889262958744e-01   ]});
    %% 计算
    area=0;
    noden=data(n).node;
    coefn=data(n).coef;
    k=n;
    n=10*n;  %n转换为表示Gauss积分的积分节点或积分系数的个数
    if n/2-floor(n/2)==1/2 %当阶数为奇数
        for i=1:floor(n/2)
            area=area+coefn(i)*(expfunc(noden(i))+expfunc(-noden(i)));
        end
        i=floor(n/2)+1;
        area=area+coefn(i)*expfunc(noden(i));
    elseif n/2-floor(n/2)==0 %当阶数为偶数
        for i=1:n/2
            area=area+coefn(i)*(expfunc(noden(i))+expfunc(-noden(i)));
        end
    end
    %% 误差:前向(n-1阶)差商作为事后误差估计
    if withError
        area_=gaussHermiteIntegral(func,a,b,k-1,false);
        gErrer=abs(area-area_);
    end
    end
    

    2 Gauss-Laguerre积分

    2.1 测试Gauss-Laguerre积分函数gaussLaguerreIntegral()

    clc,clear
    format long
    % f=@(x)((x+1).^(-2))
    % f=@(x)(exp(-x.^2))
    % f=@(x)(exp(-x))
    % f=@(x)(exp(-x.^4))
    f=@(x)((1/sqrt(2*pi))*exp(-x.^2/2)) %标准正太分布
    area=[];nn=[];
    for n=1:10
        nn=[nn,n];
        area=[area,gaussLaguerreIntegral(f,0,inf,n)];
    end,nn
    area_str=printArray(area)
    darea=abs(diff(area));
    %% 系统计算的积分精确值
    % area_system=quadgk(f,0,inf,'RelTol',1e-8,'AbsTol',1e-12)
    d_my_system=abs(area-.5);
    d_my_system_str=printArray(d_my_system)
    log_darea=round(log10(darea)*10);
    log_d_my_system=round(log10(d_my_system)*10);
    log_darea_str=printArray(log_darea)
    log_d_my_system_str=printArray(log_d_my_system)
    %% 测试误差项
    n=5
    [area,gErrer]=gaussLaguerreIntegral(f,0,inf,n,true)
    

    运行结果:

    f = 
        @(x)((1/sqrt(2*pi))*exp(-x.^2/2))
    n =
         5
    area =
       0.500000012302732
    gErrer =
         9.329705445981773e-07
    

    2.2 gaussLaguerreIntegral()

    function [area,gErrer]=gaussLaguerreIntegral(func,a,b,n,withError)
    %高斯积分
    %理论依据:数值分析方法 奚梅成
    %func 被积函数;
    %[a,b] 积分区间
    %n 高斯积分的阶数
    %%
    %name:邓能财     Date: 2013/12/21
    %% 默认参数
    if nargin<5 withError=false; gErrer=0; end
    %% 查错
    maxn=10;
    % assert(a==0 && b==inf,'积分区间必须为[0,inf)!');
    % assert(ismember(n,[1:maxn]),['积分精确度数量级必须在1~',num2str(maxn),'之间!']);
    %% 函数变换
    expfunc=@(x)(exp(x)*func(x));
    %% 数据:h,2h,3h,...,maxn*h阶高斯积分的积分系数
    h=5;
        %   (积分节点node,积分系数coef)
    data= struct('node',cell(1,maxn),'coef',cell(1,maxn));
    %data()= struct('node',{[]},'coef',{[]});
    data(1)= struct('node',{[2.6356031971814e-01
        1.4134030591065e+00
        3.5964257710407e+00
        7.0858100058588e+00
        1.2640800844276e+01   ]},'coef',{[5.2175561058281e-01
        3.9866681108317e-01
        7.5942449681707e-02
        3.6117586799220e-03
        2.3369972385776e-05   ]});
    data(2)= struct('node',{[1.3779347054049e-01
        7.2945454950317e-01
        1.8083429017403e+00
        3.4014336978549e+00
        5.5524961400638e+00
        8.3301527467645e+00
        1.1843785837900e+01
        1.6279257831378e+01
        2.1996585811981e+01
        2.9920697012274e+01   ]},'coef',{[3.0844111576502e-01
        4.0111992915527e-01
        2.1806828761181e-01
        6.2087456098678e-02
        9.5015169751811e-03
        7.5300838858754e-04
        2.8259233495996e-05
        4.2493139849627e-07
        1.8395648239796e-09
        9.9118272196090e-13   ]});
    data(3)= struct('node',{[9.3307812017283e-02
        4.9269174030189e-01
        1.2155954120710e+00
        2.2699495262037e+00
        3.6676227217514e+00
        5.4253366274136e+00
        7.5659162266131e+00
        1.0120228568019e+01
        1.3130282482176e+01
        1.6654407708330e+01
        2.0776478899449e+01
        2.5623894226729e+01
        3.1407519169754e+01
        3.8530683306486e+01
        4.8026085572686e+01   ]},'coef',{[2.1823488594009e-01
        3.4221017792288e-01
        2.6302757794168e-01
        1.2642581810593e-01
        4.0206864921001e-02
        8.5638778036118e-03
        1.2124361472143e-03
        1.1167439234425e-04
        6.4599267620229e-06
        2.2263169070963e-07
        4.2274303849794e-09
        3.9218972670411e-11
        1.4565152640731e-13
        1.4830270511133e-16
        1.6005949062111e-20   ]});
    data(4)= struct('node',{[7.0539889691990e-02
        3.7212681800161e-01
        9.1658210248327e-01
        1.7073065310283e+00
        2.7491992553094e+00
        4.0489253138509e+00
        5.6151749708616e+00
        7.4590174536711e+00
        9.5943928695811e+00
        1.2038802546964e+01
        1.4814293442631e+01
        1.7948895520519e+01
        2.1478788240285e+01
        2.5451702793187e+01
        2.9932554631701e+01
        3.5013434240479e+01
        4.0833057056729e+01
        4.7619994047347e+01
        5.5810795750064e+01
        6.6524416525616e+01   ]},'coef',{[1.6874680185112e-01
        2.9125436200607e-01
        2.6668610286700e-01
        1.6600245326951e-01
        7.4826064668792e-02
        2.4964417309283e-02
        6.2025508445722e-03
        1.1449623864769e-03
        1.5574177302781e-04
        1.5401440865225e-05
        1.0864863665180e-06
        5.3301209095567e-08
        1.7579811790506e-09
        3.7255024025123e-11
        4.7675292515782e-13
        3.3728442433624e-15
        1.1550143395004e-17
        1.5395221405823e-20
        5.2864427255692e-24
        1.6564566124990e-28   ]});
    data(5)= struct('node',{[5.6704775452705e-02
        2.9901089858699e-01
        7.3590955543502e-01
        1.3691831160352e+00
        2.2013260537215e+00
        3.2356758035580e+00
        4.4764966150738e+00
        5.9290837627004e+00
        7.5998993099567e+00
        9.4967492209324e+00
        1.1629014911779e+01
        1.4007957976545e+01
        1.6647125597289e+01
        1.9562898011469e+01
        2.2775241986835e+01
        2.6308772390969e+01
        3.0194291163316e+01
        3.4471097571922e+01
        3.9190608803937e+01
        4.4422349336162e+01
        5.0264574993834e+01
        5.6864967173940e+01
        6.4466670615954e+01
        7.3534234792100e+01
        8.5260155562496e+01   ]},'coef',{[1.3752601422934e-01
        2.5164527376491e-01
        2.5617600280975e-01
        1.8621549036244e-01
        1.0319984810752e-01
        4.4714161129934e-02
        1.5305232886396e-02
        4.1524146328771e-03
        8.9209907325968e-04
        1.5115601916424e-04
        2.0065531801933e-05
        2.0677743964319e-06
        1.6346520222911e-07
        9.7660150621244e-09
        4.3277207941849e-10
        1.3896009633895e-11
        3.1389227925400e-13
        4.8026148226043e-15
        4.7358853648073e-17
        2.8142053798431e-19
        9.1649543959912e-22
        1.4189400094973e-24
        8.2736519440991e-28
        1.1688817115428e-31
        1.3158315000591e-36   ]});data(6)= struct('node',{[4.7407180540804e-02   
    2.4992391675316e-01   
    6.1483345439277e-01   
    1.1431958256661e+00   
    1.8364545546226e+00   
    2.6965218745572e+00   
    3.7258145077795e+00   
    4.9272937658499e+00   
    6.3045155909651e+00   
    7.8616932933703e+00   
    9.6037759854793e+00   
    1.1536546597956e+01   
    1.3666744693064e+01   
    1.6002221188981e+01   
    1.8552134840143e+01   
    2.1327204321783e+01   
    2.4340035764533e+01   
    2.7605554796781e+01   
    3.1141586701111e+01   
    3.4969652008249e+01   
    3.9116084949068e+01   
    4.3613652908485e+01   
    4.8503986163804e+01   
    5.3841385406508e+01   
    5.9699121859235e+01   
    6.6180617794438e+01   
    7.3441238595560e+01   
    8.1736810506728e+01   
    9.1556466522537e+01   
    1.0415752443106e+02   ]},'coef',{[1.1604408602039e-01   
    2.2085112475070e-01   
    2.4139982758787e-01   
    1.9463676844642e-01   
    1.2372841596688e-01   
    6.3678780368988e-02   
    2.6860475273381e-02   
    9.3380708816042e-03   
    2.6806968913369e-03   
    6.3512912194088e-04   
    1.2390745990689e-04   
    1.9828788438953e-05   
    2.5893509291317e-06   
    2.7409428405359e-07   
    2.3328311650249e-08   
    1.5807455747823e-09   
    8.4274791230481e-11   
    3.4851612347541e-12   
    1.0990180599189e-13   
    2.5883126679993e-15   
    4.4378380184999e-17   
    5.3659179405818e-19   
    4.3939492584790e-21   
    2.3114466007031e-23   
    7.2761720424752e-26   
    1.2409628392253e-28   
    9.9460139172742e-32   
    3.0921279934282e-35   
    3.1849751877368e-39   
    1.7664292620495e-43   ]});
    data(7)= struct('node',{[4.0729209061714e-02   
    2.1468745273514e-01   
    5.2801038431934e-01   
    9.8138617345908e-01   
    1.5757259475775e+00   
    2.3122282965132e+00   
    3.1923979394931e+00   
    4.2180641250297e+00   
    5.3914027346738e+00   
    6.7149632799157e+00   
    8.1917017588376e+00   
    9.8250204970927e+00   
    1.1618816388907e+01   
    1.3577539357308e+01   
    1.5706263395763e+01   
    1.8010773287835e+01   
    2.0497671107147e+01   
    2.3174507997799e+01   
    2.6049948710370e+01   
    2.9133979209544e+01   
    3.2438171837678e+01   
    3.5976028769894e+01   
    3.9763434104800e+01   
    4.3819260118607e+01   
    4.8166197974019e+01   
    5.2831925058016e+01   
    5.7850795022134e+01   
    6.3266373675385e+01   
    6.9135413883648e+01   
    7.5534435134730e+01   
    8.2571405523583e+01   
    9.0408523864401e+01   
    9.9312973659441e+01   
    1.0979599094491e+02   
    1.2317325317538e+02   ]},'coef',{[1.0036223741924e-01   
    1.9648238336171e-01   
    2.2601456693434e-01   
    1.9627166789644e-01   
    1.3760074201313e-01   
    8.0030400464140e-02   
    3.9125771428883e-02   
    1.6186726526126e-02   
    5.6853255717743e-03   
    1.6972232991260e-03   
    4.3048506814168e-04   
    9.2634053063654e-05   
    1.6870264869314e-05   
    2.5916721474280e-06   
    3.3446109572665e-07   
    3.6077854386912e-08   
    3.2335690322103e-09   
    2.3913229815959e-10   
    1.4473323990302e-11   
    7.1013552776266e-13   
    2.7933695098016e-14   
    8.6948402151367e-16   
    2.1088459299447e-17   
    3.9129108212979e-19   
    5.4327218369250e-21   
    5.4936982290015e-23   
    3.9126007648144e-25   
    1.8797445108978e-27   
    5.7353329163766e-30   
    1.0118890306195e-32   
    6.4972453567684e-36   
    1.0249186153468e-41   
    1.0144829299011e-40   
    1.8914407768105e-43   
    2.6963083454085e-47   ]});
    data(8)= struct('node',{[3.5700394308884e-02   
    1.8816228315870e-01   
    4.6269428131458e-01   
    8.5977296397294e-01   
    1.3800108205273e+00   
    2.0242091359228e+00   
    2.7933693535068e+00   
    3.6887026779083e+00   
    4.7116411465550e+00   
    5.8638508783437e+00   
    7.1472479081023e+00   
    8.5640170175862e+00   
    1.0116634048452e+01   
    1.1807892294005e+01   
    1.3640933712537e+01   
    1.5619285893339e+01   
    1.7746905950096e+01   
    2.0028232834575e+01   
    2.2468249983498e+01   
    2.5072560772426e+01   
    2.7847480009169e+01   
    3.0800145739445e+01   
    3.3938657084914e+01   
    3.7272245880476e+01   
    4.0811492823887e+01   
    4.4568603175334e+01   
    4.8557763533060e+01   
    5.2795611187217e+01   
    5.7301863323394e+01   
    6.2100179072775e+01   
    6.7219370927127e+01   
    7.2695158847612e+01   
    7.8572802911571e+01   
    8.4911231135705e+01   
    9.1789874671236e+01   
    9.9320808717447e+01   
    1.0767244063939e+02   
    1.1712230951269e+02   
    1.2820184198826e+02   
    1.4228004446916e+02   ]},'coef',{[8.8412106190338e-02   
    1.7681473909573e-01   
    2.1136311701596e-01   
    1.9408119531860e-01   
    1.4643428242412e-01   
    9.3326798435770e-02   
    5.0932204361044e-02   
    2.3976193015685e-02   
    9.7746252467145e-03   
    3.4579399930185e-03   
    1.0622468938969e-03   
    2.8327168532433e-04   
    6.5509405003248e-05   
    1.3116069073268e-05   
    2.2684528787793e-06   
    3.3796264822014e-07   
    4.3228213222805e-08   
    4.7284937709896e-09   
    4.4031741042515e-10   
    3.4724414847684e-11   
    2.3053815448611e-12   
    1.2797725979698e-13   
    5.8941771706254e-15   
    2.2322175699210e-16   
    6.8803366293934e-18   
    1.7056039068355e-19   
    3.3537088111640e-21   
    5.1461768068696e-23   
    6.0453446308727e-25   
    5.3139346735449e-27   
    3.3657690232258e-29   
    1.2178445191183e-31   
    4.2890086271870e-35   
    2.3067506710852e-34   
    2.3302702800589e-35   
    7.9126811843002e-37   
    1.4714096228146e-38   
    6.8575263257425e-41   
    1.0375271306839e-43   
    1.2955873193971e-47   ]});
    data(9)= struct('node',{[3.1776956867006e-02   
    1.6747251657590e-01   
    4.1176935189346e-01   
    7.6501433412688e-01   
    1.2276394084198e+00   
    1.8002069428945e+00   
    2.4834175141916e+00   
    3.2781154527750e+00   
    4.1852949337855e+00   
    5.2061071301655e+00   
    6.3418686461714e+00   
    7.5940714118218e+00   
    8.9643942355385e+00   
    1.0454716247742e+01   
    1.2067132515935e+01   
    1.3803972172016e+01   
    1.5667819467697e+01   
    1.7661538267980e+01   
    1.9788300611284e+01   
    2.2051620115564e+01   
    2.4455391203194e+01   
    2.7003935367666e+01   
    2.9702056032584e+01   
    3.2555103986001e+01   
    3.5569055951601e+01   
    3.8750609641046e+01   
    4.2107299705812e+01   
    4.5647640501793e+01   
    4.9381303694811e+01   
    5.3319341780002e+01   
    5.7474473058814e+01   
    6.1861450326356e+01   
    6.6497545839445e+01   
    7.1403201447971e+01   
    7.6602919389220e+01   
    8.2126514283411e+01   
    8.8010926374694e+01   
    9.4302943609710e+01   
    1.0106347092232e+02   
    1.0837460196798e+02   
    1.1635218467520e+02   
    1.2517034745017e+02   
    1.3511619015017e+02   
    1.4673979649750e+02   
    1.6145944790909e+02   ]},'coef',{[7.9003863218578e-02   
    1.6064969814262e-01   
    1.9788792774995e-01   
    1.8978490012961e-01   
    1.5161633179829e-01   
    1.0374664065632e-01   
    6.1658521650479e-02   
    3.2073009065685e-02   
    1.4666691055305e-02   
    5.9109404105920e-03   
    2.1021480717925e-03   
    6.5997635763383e-04   
    1.8287173922646e-04   
    4.4687548597413e-05   
    9.6188526495207e-06   
    1.8207784506475e-06   
    3.0249820529759e-07   
    4.4003996132627e-08   
    5.5894745239333e-09   
    6.1800526258555e-10   
    5.9265751381681e-11   
    4.9096834353441e-12   
    3.4975836508939e-13   
    2.1317103205750e-14   
    1.1051884160140e-15   
    4.8425951880898e-17   
    1.7802221360264e-18   
    5.4452758047095e-20   
    1.3728113205342e-21   
    2.8219449500093e-23   
    4.6711191673266e-25   
    6.1366267898675e-27   
    6.3012553677705e-29   
    4.9199469546978e-31   
    2.8553923595618e-33   
    2.6513859772935e-35   
    1.7302431662816e-35   
    3.7028624824552e-36   
    1.3728993488078e-37   
    1.3069964719048e-38   
    3.7845759572530e-40   
    6.4931338105968e-44   
    1.0904747480551e-45   
    5.8188542963473e-48   
    7.4815330585462e-52   ]});
    data(10)= struct('node',{[2.8630518339375e-02   
    1.5088293567693e-01   
    3.7094878153490e-01   
    6.8909069988105e-01   
    1.1056250235399e+00   
    1.6209617511025e+00   
    2.2356103759152e+00   
    2.9501833666418e+00   
    3.7653997744058e+00   
    4.6820893875593e+00   
    5.7011975747849e+00   
    6.8237909097945e+00   
    8.0510636693908e+00   
    9.3843453082584e+00   
    1.0825109031549e+01   
    1.2374981608757e+01   
    1.4035754599830e+01   
    1.5809397197845e+01   
    1.7698070933350e+01   
    1.9704146535462e+01   
    2.1830223306578e+01   
    2.4079151444412e+01   
    2.6454057841253e+01   
    2.8958376011937e+01   
    3.1595880956623e+01   
    3.4370729963090e+01   
    3.7287510610551e+01   
    4.0351297573586e+01   
    4.3567720269995e+01   
    4.6943043991603e+01   
    5.0484267963130e+01   
    5.4199244880169e+01   
    5.8096828017249e+01   
    6.2187054175689e+01   
    6.6481373878445e+01   
    7.0992944826619e+01   
    7.5737011547727e+01   
    8.0731404802478e+01   
    8.5997211136463e+01   
    9.1559690412534e+01   
    9.7449565614851e+01   
    1.0370489123669e+02   
    1.1037385880764e+02   
    1.1751919820311e+02   
    1.2522547013347e+02   
    1.3361202792273e+02   
    1.4285832548925e+02   
    1.5326037197260e+02   
    1.6538564331668e+02   
    1.8069834370921e+02   ]},'coef',{[7.1404726135184e-02   
    1.4714860696459e-01   
    1.8567162757483e-01   
    1.8438538252735e-01   
    1.5420116860636e-01   
    1.1168536990227e-01   
    7.1052885490196e-02   
    4.0020276911508e-02   
    2.0050623080072e-02   
    8.9608512036464e-03   
    3.5781124153157e-03   
    1.2776171567891e-03   
    4.0803024498372e-04   
    1.1652883223097e-04   
    2.9741704936942e-05   
    6.7778425265418e-06   
    1.3774795031715e-06   
    2.4928861817200e-07   
    4.0103543504264e-08   
    5.7233317481479e-09   
    7.2294342491785e-10   
    8.0617101422267e-11   
    7.9133931001156e-12   
    6.8157366176493e-13   
    5.1324267151836e-14   
    3.3656247656300e-15   
    1.9134763294156e-16   
    9.3855895940761e-18   
    3.9500701020349e-19   
    1.4177501120800e-20   
    4.3099605395263e-22   
    1.1012486685628e-23   
    2.3448759094638e-25   
    4.1195425783619e-27   
    5.8517301056039e-29   
    6.7131788910435e-31   
    1.2697958123263e-32   
    6.9071714060483e-34   
    3.9707498882293e-34   
    1.2868077986345e-33   
    1.1748415289134e-33   
    1.8729029481275e-34   
    2.9809110912453e-35   
    2.1474665129166e-36   
    5.6599382245170e-38   
    1.2115127694324e-39   
    8.4123926853140e-42   
    1.1227629730641e-44   
    9.4963283290577e-48   
    4.9608486792055e-52   ]});
    %% 计算
    area=0;
    noden=data(n).node;
    coefn=data(n).coef;
    k=h*n;  %k表示Gauss积分的积分节点或积分系数的个数
    for i=1:k
        area=area+coefn(i)*expfunc(noden(i));
    end
    %% 误差:前向(n-1阶)差商作为事后误差估计
    if withError
        area_=gaussLaguerreIntegral(func,a,b,n-1,false);
        gErrer=abs(area-area_);
    end
    end
    

    输出函数:

    function str=printArray(a)
    str=sprintf('%.13e\t',a);
    

    3 Gauss-Legendre积分W(x)=1

    3.1 测试Gauss-Legendre积分W(x)=1函数gaussLegendreIntegral()

    clc,clear
    format long
    % f=@(x)(exp(x))
    % a=1,b=10
    f=@(x)(cos(x))
    % a=0,b=pi/2
    a=-pi/2,b=pi/2
    % a=-pi,b=pi
    % a=-1,b=1
    disp('各个阶数(h*)的gaussLegendreIntegral积分')
    h=5
    max=5
    n_=[2:max];area_=[];gErrer_=[];
    for n=2:max
    %% 适合区间长为数量级10左右的积分
        [area,gErrer]=gaussLegendreIntegral(f,a,b,n,true);
        area_=[area_,area];gErrer_=[gErrer_,gErrer];
    end
    n_
    area__str=printArray(area_)
    gErrer_str=printArray(gErrer_)
    area_system=quadl(f,a,b)
    

    运行结果:

    f = 
        @(x)(cos(x))
    a =
      -1.570796326794897
    b =
       1.570796326794897
    各个阶数(h*)的gaussLegendreIntegral积分
    h =
         5
    max =
         5
    n_ =
         2     3     4     5
    area__str =
    2.0000000000000e+00	2.0000000000000e+00	2.0000000000000e+00	2.0000000000000e+00	
    gErrer_str =
    1.1028448909656e-07	3.0864200084579e-14	3.0864200084579e-14	2.7089441800854e-14	
    area_system =
       1.999999977471133
    >> 
    

    3.2 gaussLegendreIntegral()

    function [area,gErrer]=gaussLegendreIntegral(func,a,b,n,withError)
    %高斯积分
    %理论依据:数值分析方法 奚梅成
    %func 被积函数;
    %[a,b] 积分区间
    %n 高斯积分的阶数
    %% 
    %name:邓能财     Date: 2013/12/21
    %% 默认参数
    if nargin<5 withError=false; gErrer=0; end
    %% 数据:1-8阶高斯积分的积分系数
    h=5;
    maxn=10; %   (±积分节点node,积分系数coef),‘±’被省略
    data= struct('node',cell(1,maxn),'coef',cell(1,maxn));
    % data(1)= struct('node',{[0           ]},'coef',{[2]});
    % data(2)= struct('node',{[0.5773502692]},'coef',{[1]});
    % data(3)= struct('node',{[0.7745966692 0           ]},'coef',{[0.5555555556 0.8888888889]});
    % data(4)= struct('node',{[0.8611363116 0.3399810436]},'coef',{[0.3478548451 0.6521451549]});
    % data(5)= struct('node',{[0.9061798459 0.5384693101 0           ]},'coef',{[0.2369268851 0.4786286705 0.5688888889]});
    % data(6)= struct('node',{[0.9324695142 0.6612093865 0.2386191861]},'coef',{[0.1713244924 0.3607615730 0.4679139346]});
    % data(7)= struct('node',{[0.9491079123 0.7415311856 0.4058451514 0           ]},'coef',{[0.1294849662 0.2797053915 0.3818300505 0.4179591837]});
    % data(8)= struct('node',{[0.9602898565 0.7966664774 0.5255324099 0.1834346425]},'coef',{[0.1012285363 0.2223810345 0.3137066459 0.3626837834]});
    data(1)= struct('node',{[9.0617984593866e-01   
    5.3846931010568e-01   
    -4.7505928254367e-17   ]},'coef',{[2.3692688505619e-01   
    4.7862867049937e-01   
    5.6888888888889e-01   ]});
    data(2)= struct('node',{[9.7390652851717e-01   
    8.6506336668898e-01   
    6.7940956829902e-01   
    4.3339539412925e-01   
    1.4887433898163e-01   ]},'coef',{[6.6671344308688e-02   
    1.4945134915058e-01   
    2.1908636251598e-01   
    2.6926671931000e-01   
    2.9552422471475e-01   ]});
    data(3)= struct('node',{[9.8799251802049e-01   
    9.3727339240071e-01   
    8.4820658341043e-01   
    7.2441773136017e-01   
    5.7097217260854e-01   
    3.9415134707756e-01   
    2.0119409399743e-01   
    1.0426645532607e-16   ]},'coef',{[3.0753241996117e-02   
    7.0366047488109e-02   
    1.0715922046717e-01   
    1.3957067792615e-01   
    1.6626920581699e-01   
    1.8616100001556e-01   
    1.9843148532711e-01   
    2.0257824192556e-01   ]});
    data(4)= struct('node',{[9.9312859918510e-01   
    9.6397192727791e-01   
    9.1223442825133e-01   
    8.3911697182222e-01   
    7.4633190646015e-01   
    6.3605368072652e-01   
    5.1086700195083e-01   
    3.7370608871542e-01   
    2.2778585114164e-01   
    7.6526521133497e-02   ]},'coef',{[1.7614007139152e-02   
    4.0601429800386e-02   
    6.2672048334110e-02   
    8.3276741576705e-02   
    1.0193011981724e-01   
    1.1819453196152e-01   
    1.3168863844918e-01   
    1.4209610931838e-01   
    1.4917298647260e-01   
    1.5275338713073e-01   ]});
    data(5)= struct('node',{[9.9555696979050e-01   
    9.7666392145952e-01   
    9.4297457122897e-01   
    8.9499199787828e-01   
    8.3344262876083e-01   
    7.5925926303736e-01   
    6.7356636847347e-01   
    5.7766293024122e-01   
    4.7300273144572e-01   
    3.6117230580939e-01   
    2.4386688372099e-01   
    1.2286469261071e-01   
    2.3477116549018e-17   ]},'coef',{[1.1393798501026e-02   
    2.6354986615032e-02   
    4.0939156701307e-02   
    5.4904695975835e-02   
    6.8038333812357e-02   
    8.0140700335000e-02   
    9.1028261982964e-02   
    1.0053594906705e-01   
    1.0851962447426e-01   
    1.1485825914571e-01   
    1.1945576353578e-01   
    1.2224244299031e-01   
    1.2317605372672e-01   ]});
    data(6)= struct('node',{[9.9689348407465e-01   
    9.8366812327975e-01   
    9.6002186496831e-01   
    9.2620004742927e-01   
    8.8256053579205e-01   
    8.2956576238277e-01   
    7.6777743210483e-01   
    6.9785049479332e-01   
    6.2052618298924e-01   
    5.3662414814202e-01   
    4.4703376953809e-01   
    3.5270472553088e-01   
    2.5463692616789e-01   
    1.5386991360858e-01   
    5.1471842555318e-02   ]},'coef',{[7.9681924961671e-03   
    1.8466468311090e-02   
    2.8784707883324e-02   
    3.8799192569627e-02   
    4.8402672830593e-02   
    5.7493156217619e-02   
    6.5974229882181e-02   
    7.3755974737705e-02   
    8.0755895229420e-02   
    8.6899787201083e-02   
    9.2122522237786e-02   
    9.6368737174644e-02   
    9.9593420586795e-02   
    1.0176238974840e-01   
    1.0285265289356e-01   ]});
    data(7)= struct('node',{[9.9770656909960e-01   
    9.8793576444385e-01   
    9.7043761603923e-01   
    9.4534514820783e-01   
    9.1285426135932e-01   
    8.7321912502522e-01   
    8.2674989909223e-01   
    7.7381025228691e-01   
    7.1481450155663e-01   
    6.5022436466589e-01   
    5.8054534474976e-01   
    5.0632277324149e-01   
    4.2813754151781e-01   
    3.4660155443081e-01   
    2.6235294120930e-01   
    1.7605106116599e-01   
    8.8371343275659e-02   
    -8.0123445265982e-18   ]},'coef',{[5.8834334204433e-03   
    1.3650828348362e-02   
    2.1322979911483e-02   
    2.8829260108895e-02   
    3.6110115863463e-02   
    4.3108422326170e-02   
    4.9769370401354e-02   
    5.6040816212370e-02   
    6.1873671966080e-02   
    6.7222285269087e-02   
    7.2044794772560e-02   
    7.6303457155442e-02   
    7.9964942242325e-02   
    8.3000593728857e-02   
    8.5386653392099e-02   
    8.7104446997184e-02   
    8.8140530430276e-02   
    8.8486794907104e-02   ]});
    data(8)= struct('node',{[9.9823770971056e-01   
    9.9072623869946e-01   
    9.7725994998377e-01   
    9.5791681921379e-01   
    9.3281280827868e-01   
    9.0209880696887e-01   
    8.6595950321226e-01   
    8.2461223083331e-01   
    7.7830565142652e-01   
    7.2731825518993e-01   
    6.7195668461418e-01   
    6.1255388966798e-01   
    5.4946712509513e-01   
    4.8307580168618e-01   
    4.1377920437160e-01   
    3.4199409082576e-01   
    2.6815218500725e-01   
    1.9269758070137e-01   
    1.1608407067526e-01   
    3.8772417506051e-02   ]},'coef',{[4.5212770985336e-03   
    1.0498284531152e-02   
    1.6421058381909e-02   
    2.2245849194167e-02   
    2.7937006980023e-02   
    3.3460195282548e-02   
    3.8782167974472e-02   
    4.3870908185673e-02   
    4.8695807635072e-02   
    5.3227846983937e-02   
    5.7439769099392e-02   
    6.1306242492929e-02   
    6.4804013456601e-02   
    6.7912045815234e-02   
    7.0611647391287e-02   
    7.2886582395805e-02   
    7.4723169057968e-02   
    7.6110361900626e-02   
    7.7039818164248e-02   
    7.7505947978425e-02   ]});
    data(9)= struct('node',{[9.9860364518194e-01   
    9.9264999844720e-01   
    9.8196871503454e-01   
    9.6660831039689e-01   
    9.4664169099563e-01   
    9.2216393671900e-01   
    8.9329167175324e-01   
    8.6016247596066e-01   
    8.2293422050209e-01   
    7.8178431259391e-01   
    7.3690884894549e-01   
    6.8852168077120e-01   
    6.3685339445322e-01   
    5.8215021256935e-01   
    5.2467282046292e-01   
    4.6469512391964e-01   
    4.0250294385854e-01   
    3.3839265425060e-01   
    2.7266976975238e-01   
    2.0564748978326e-01   
    1.3764520598325e-01   
    6.8986980163144e-02   
    1.0816665110908e-16   ]},'coef',{[3.5826631552839e-03   
    8.3231892962177e-03   
    1.3031104991583e-02   
    1.7677535257937e-02   
    2.2239847550579e-02   
    2.6696213967578e-02   
    3.1025374934516e-02   
    3.5206692201609e-02   
    3.9220236729302e-02   
    4.3046880709165e-02   
    4.6668387718374e-02   
    5.0067499237952e-02   
    5.3228016731269e-02   
    5.6134878759786e-02   
    5.8774232718842e-02   
    6.1133500831067e-02   
    6.3201440073820e-02   
    6.4968195750723e-02   
    6.6425348449842e-02   
    6.7565954163608e-02   
    6.8384577378670e-02   
    6.8877316977661e-02   
    6.9041824829232e-02   ]});
    data(10)= struct('node',{[9.9886640442007e-01   
    9.9403196943209e-01   
    9.8535408404801e-01   
    9.7286438510669e-01   
    9.5661095524281e-01   
    9.3665661894488e-01   
    9.1307855665579e-01   
    8.8596797952361e-01   
    8.5542976942995e-01   
    8.2158207085934e-01   
    7.8455583290040e-01   
    7.4449430222607e-01   
    7.0155246870682e-01   
    6.5589646568544e-01   
    6.0770292718495e-01   
    5.5715830451465e-01   
    5.0445814490746e-01   
    4.4980633497404e-01   
    3.9341431189757e-01   
    3.3550024541944e-01   
    2.7628819377953e-01   
    2.1600723687604e-01   
    1.5489058999815e-01   
    9.3174701560086e-02   
    3.1098338327189e-02   ]},'coef',{[2.9086225531552e-03   
    6.7597991957456e-03   
    1.0590548383651e-02   
    1.4380822761485e-02   
    1.8115560713489e-02   
    2.1780243170125e-02   
    2.5360673570013e-02   
    2.8842993580535e-02   
    3.2213728223579e-02   
    3.5459835615146e-02   
    3.8568756612588e-02   
    4.1528463090147e-02   
    4.4327504338803e-02   
    4.6955051303949e-02   
    4.9400938449466e-02   
    5.1655703069581e-02   
    5.3710621888997e-02   
    5.5557744806213e-02   
    5.7189925647728e-02   
    5.8600849813222e-02   
    5.9785058704266e-02   
    6.0737970841770e-02   
    6.1455899590316e-02   
    6.1936067420683e-02   
    6.2176616655347e-02   ]});
    %% 积分变量的线性变换
    y=@(x)(((b+a)+(b-a)*x)/2);
    %% 计算
    area=0;
    noden=data(n).node;
    coefn=data(n).coef;
    k=n;
    n=n*h;
    if n/2-floor(n/2)==1/2 %当阶数为奇数   
       for i=1:floor(n/2)
          area=area+coefn(i)*(func(y(noden(i)))+func(y(-1*noden(i))));
       end
       i=floor(n/2)+1;
       area=area+coefn(i)*func(y(noden(i)));
    elseif n/2-floor(n/2)==0 %当阶数为偶数  
       for i=1:n/2
          area=area+coefn(i)*(func(y(noden(i)))+func(y(-1*noden(i))));
       end
    end
    area=area*(b-a)/2;
    %% 误差:前向(n-1阶)差商作为事后误差估计
    if withError
        area_=gaussLegendreIntegral(func,a,b,k-1,false);
        gErrer=abs(area-area_);
    end
    end
    

    输出函数:

    function str=printArray(a)
    str=sprintf('%.13e\t',a);
    

    4 Gauss积分的积分系数计算

    4.1 测试Gauss积分的积分系数计算的函数

    clc,clear
    % GaussLegendre积分:(积分节点    积分系数)
    % disp('GaussLegendre积分:(积分节点+-node    积分系数coef)')
    % for n=1:16
    %     sprintf('************阶数 %d',n)
    %     [x, w] = GaussLegendre_2(n);
    %     k=floor(n/2);
    %     if k~=n/2, spac=1:k+1;else, spac=1:n/2; end
    %     node=-x(spac),coef=w(spac)
    % end
    %% GaussLaguerre积分:(积分节点    积分系数)
    % disp('GaussLaguerre积分:(积分节点    积分系数)')
    % alpha=1e-30
    % for n=1:16
    %     sprintf('************阶数 %d',n)
    %     [node,coef]=GaussLaguerre_2(n,alpha)
    % end
    %% GaussHermite 积分:(积分节点    积分系数)
    % disp('GaussHermite积分:(积分节点    积分系数)')
    % for n=1:20
    %     sprintf('************阶数 %d',n)
    %     [x, w] =GaussHermite_2(n);
    %     k=floor(n/2);
    %     if k~=n/2, spac=1:k+1;else, spac=1:n/2; end
    %     node=-x(spac),coef=w(spac)
    % end
    %%
    printGaussIntData(1,5); 
    

    printGaussIntData()函数(打印Gauss积分数据表):

    function printGaussIntData(type,h)
    %% 打印Gauss积分数据表
    str='';
    for n=h:h:h*10
        if type==2
            [x, w] =GaussLaguerre_2(n,1e-100);
        else
            if type==1
                [x, w] =GaussLegendre_2(n);
            elseif type==3
                [x, w] =GaussHermite_2(n);
            end
            k=floor(n/2);
            if k~=n/2, spac=1:k+1;else, spac=1:n/2; end
            x=-x(spac); w=w(spac);
        end
        %data()= struct('node',{[]},'coef',{[]});sprintf('%d',n)
        str=[str,'data(',num2str((n/h)),...
            ')= struct(#&#node#&#,{[',printMatrix(x,'%.13e',false),...
            ']},#&#coef#&#,{[',printMatrix(w,'%.13e',false),']});',sprintf('\n')];
    end,str
    end
    

    输出矩阵的函数:

    function s=printMatrix(m,style,with_endl)
    if nargin<2,  style='%.13e'; end
    if nargin<3,  with_endl=true; end
    [dim1,dim2]=size(m);
    s='';
    for i=1:dim1
        s=[s,sprintf([style,'   '],m(i,:)),sprintf('\n')];
    end
    if ~with_endl
        s(end)='';
    end
    end
    

    4.2 Gauss-Hermite积分系数求解函数

    function [x, w] = GaussHermite_2(n)
    %name:邓能财     Date: 2013/12/23
    i = 1:n-1; 
    a = sqrt(i/2); 
    CM = diag(a,1) + diag(a,-1);
    [V L] = eig(CM); 
    [x ind] = sort(diag(L)); 
    V = V(:,ind)'; 
    w = sqrt(pi) * V(:,1).^2;
    end
    

    4.3 Gauss-Laguerre积分系数求解函数

    function [x, w] = GaussLaguerre_2(n, alpha)
    %name:邓能财     Date: 2013/12/23
    i = 1:n; 
    a = (2*i-1) + alpha; 
    b = sqrt( i(1:n-1) .* ((1:n-1) + alpha) ); 
    CM = diag(a) + diag(b,1) + diag(b,-1);
    [V L] = eig(CM); 
    [x ind] = sort(diag(L)); 
    V = V(:,ind)'; 
    w = gamma(alpha+1) .* V(:,1).^2;
    end
    

    4.4 Gauss-Legendre积分系数求解函数

    function [x, w] = GaussLegendre_2(n)
    % http://www.mathworks.cn/matlabcentral/fileexchange/8067-gauss-laguerre
    %对应《数值分析方法 奚梅成的Gauss-Legendre积分》
    %name:邓能财     Date: 2013/12/23
    i = 1:n-1; 
    a = i./sqrt(4*i.^2-1); 
    CM = diag(a,1) + diag(a,-1);
    [V L] = eig(CM); 
    [x ind] = sort(diag(L)); 
    V = V(:,ind)'; 
    w = 2 * V(:,1).^2;
    end
    

    5 NewtonCotes积分

    5.1 测试NewtonCotes积分函数newtonCotesIntegral()

    clc,clear
    format long
    disp('各个阶数的Newton-Cote`s积分')
    % f=@(x)((x+1).^(-2))
    f =  @(x)(cos(x))
    a = -1.57079632679490
    b = 1.57079632679490
    % f=@(x)(exp(-x.^2))
    % f=@(x)(exp(-x))
    % f=@(x)(exp(-x.^4))
    % f=@(x)((1/sqrt(2*pi))*exp(-x.^2/2)) %标准正太分布
    % f=@(x)(sin(x)./x)
    % f=@(x)(x)
    % a=0,b=1
    area=[];nn=[];gErrer=[];
    for n=2:8
        nn=[nn,n];
        [arean,gErrern]=newtonCotesIntegral(f,a,b,n,true);
        area=[area,arean];
        gErrer=[gErrer,gErrern];
    end,nn,
    area_str=printArray(area)
    gErrer_str=printArray(gErrer)
    %% 系统计算的积分精确值
    area_system=quad(f,a,b,1e-13)
    d_my_system=abs(area-area_system);
    log_darea=round(log10(gErrer)*10)
    log_d_my_system=round(log10(d_my_system)*10)
    % %% 测试误差项
    % a=0,b=.01
    % n=3
    % [area,gErrer]=newtonCotesIntegral(f,a,b,n,true)
    % area_system=quad(f,a,b,1e-13)
    

    输出数组的函数:

    function str=printArray(a)
    str=sprintf('%.13e\t',a);
    

    运行结果:

    各个阶数的Newton-Cote`s积分
    f = 
        @(x)(cos(x))
    a =
      -1.570796326794900
    b =
       1.570796326794900
    nn =
         2     3     4     5     6     7     8
    area_str =
    1.5707963267949e+00	2.2672492052928e+00	2.2888179796358e+00	2.2609034070862e+00	2.2389648201115e+00	2.2178021611519e+00	2.1987231887715e+00	
    gErrer_str =
    3.1341507153454e+00	1.7494105356846e+00	6.8206453262057e-02	1.1112991494925e-01	1.0995339715508e-01	1.3352735099341e-01	1.5154966451227e-01	
    area_system =
         2
    

    5.2 newtonCotesIntegral()

    function [area,gErrer]=newtonCotesIntegral(func,a,b,n,withError)
    %Newton-Cotes积分
    %理论依据:数值分析方法 奚梅成
    %func 被积函数;
    %[a,b] 积分区间
    %n Newton-Cotes积分的阶数
    %%
    %name:邓能财     Date: 2013/12/21
    %% 默认参数
    if nargin<5 withError=false; gErrer=0; end
    %% 查错
    maxn=8;
    % assert(a~=inf && b~=inf,'积分区间不能包含inf!');
    % assert(ismember(n,[1:maxn]),['积分精确度数量级必须在1~',num2str(maxn),'之间!']);
    %% 数据:1,2,...,maxn阶Newton-Cotes积分的积分系数
        %   (积分系数coef)
    data= struct('coef',cell(1,maxn));
    %data()= struct('coef',{[]});
    data(1)= struct('coef',{[1]});
    data(2)= struct('coef',{[ 1/2          1/2  ]});
    data(3)= struct('coef',{[ 1/6          4/6        1/6  ]});
    data(4)= struct('coef',{[ 1/8          3/8        3/8        1/8 ]});
    data(5)= struct('coef',{[ 7/90        16/45       2/15       16/45      7/90 ]});
    data(6)= struct('coef',{[ 19/288      25/96      25/144      25/144     25/96      19/288 ]});
    data(7)= struct('coef',{[ 41/840      9/35       9/280       34/105     9/280       9/35      41/840 ]});
    data(8)= struct('coef',{[ 751/17280 3577/17280 1323/17280  2989/17280 2989/17280  1323/17280  3577/17280  751/17280 ]});
    %% 计算
    area=0;
    coefn=data(n).coef;
    xi=a; h=(b-a)/n;
    for i=1:n
        xi=xi+h;
        area=area+coefn(i)*func(xi);
    end
    area=area*(b-a);
    %% 误差:前向(n-1阶)差商作为事后误差估计
    if withError
        area_=newtonCotesIntegral(func,a,b,n-1,false);
        gErrer=abs(area-area_);
        gErrer=gErrer*(10^(.1*(n+1))); %根据试验值调整误差估计
    end
    % %% 误差:两分段(n-1阶)积分和与原积分之差作为事后误差估计
    % if withError
    %     area_=(newtonCotesIntegral(func,a,(a+b)/2,n,false)+...
    %            newtonCotesIntegral(func,(a+b)/2,b,n,false));
    %     gErrer=abs(area-area_);
    % end
    end
    

    6 Romberg积分

    6.1 测试Romberg积分函数rombergIntegral()

    clc,clear
    format long
    f =  @(x)(cos(x))
    a = -1.57079632679490
    b = 1.57079632679490
    '各种精度的Romberg积分'
    % f=@(x)(exp(x))
    % a=1,b=2
    area_=[];pre_=[];k_=[];
    for i=2:5
        pre=10^(-2*i)
        [area,k]=rombergIntegral(f,a,b,pre);
        pre_=[pre_,pre];area_=[area_,area];k_=[k_,k];
    end,
    k_
    area_str=printArray(area_)
    pre_str=printArray(pre_)
    area_system=quadl(f,a,b)
    

    输出数组的函数:

    function str=printArray(a)
    str=sprintf('%.13e\t',a);
    

    运行结果:

    f = 
        @(x)(cos(x))
    a =
      -1.570796326794900
    b =
       1.570796326794900
    ans =
    各种精度的Romberg积分
    pre =
         1.000000000000000e-04
    pre =
         1.000000000000000e-06
    pre =
         1.000000000000000e-08
    pre =
         1.000000000000000e-10
    k_ =
         4     5     6     6
    area_str =
    2.0000055499797e+00	1.9999999945873e+00	2.0000000000013e+00	2.0000000000013e+00	
    pre_str =
    1.0000000000000e-04	1.0000000000000e-06	1.0000000000000e-08	1.0000000000000e-10	
    area_system =
       1.999999977471133
    

    6.2 rombergIntegral()

    function [area,k]=rombergIntegral(func,a,b,prec)
    %Romberg积分
    %理论依据:数值分析方法 奚梅成
    %func 函数;
    %[a,b] 积分区间
    %pre 精度
    %程序特点:密切结合理论,将内存和运算量达到极小化
        %只记录最后一排需用的数字
    %%
    %name:邓能财     Date: 2013/12/20
    %这几乎是我编写的第一个【调试只出一个错误】的程序
    %% 计算
    %取初始分段数为1
    h=b-a;
    t=[];%记录最后一排需用的数字
    %Romberg积分的阶数:k
        %即所储存的T的长度
    k=1;
    %计算梯形积分T0(h)
    t=(func(a)+func(b))/2;
    while true
        t=[0,t];
        k=k+1;
        %h1用于T0(h)与T0(h/2)之间的递推计算
        h1=0;
        for xi=a+h/2:h:b-h/2
            h1=h1+func(xi);
        end, h1=h1*h;
        %计算梯形积分T0(h/2)
        t(1)=(t(2)+h1)/2;
        h=h/2;
        %计算所得积分最精确值
        bit=1;%记录2^(2k)
        for i=2:k
            bit=bit*4;
            t(i)=t(i-1)+(t(i-1)-t(i))/(bit-1);
        end
        %判断是否达到精度要求
        if abs(t(end)-t(end-1))<prec
            break;
        end
    end
    area=t(end);
    end
    

    7 复化积分

    7.1 测试复化积分函数complexSimpson()

    clc,clear
    % f=@(x)(exp(x))
    % a=1,b=2
    % m2=f(2)
    % m4=m2
    f =   @(x)(cos(x))
    a = -1.57079632679490
    b = 1.57079632679490
    '各种精度的复化梯形积分积分'
    % m2=1
    m4=1
    area_=[];pre_=[];n_=[];
    for i=2:5
        pre=10^(-2*i);
        [area,n]=complexSimpson(f,a,b,pre,m4);
        pre_=[pre_,pre];area_=[area_,area];n_=[n_,n];
    end,
    n_
    area_str=printArray(area_)
    pre_str=printArray(pre_)
    area_system=quadl(f,a,b)
    % pre=1e-6
    % ss=complexTrangleInt(f,a,b,pre,m2)
    % ss=complexSimpson(f,a,b,pre,m4)
    

    输出数组的函数:

    function str=printArray(a)
    str=sprintf('%.13e\t',a);
    

    运行结果:

    f = 
        @(x)(cos(x))
    a =
      -1.570796326794900
    b =
       1.570796326794900
    ans =
    各种精度的复化梯形积分积分
    m4 =
         1
    n_ =
         3     9    28    89
    area_str =
    2.0008631896735e+00	2.0000103477058e+00	2.0000001100950e+00	2.0000000010782e+00	
    pre_str =
    1.0000000000000e-04	1.0000000000000e-06	1.0000000000000e-08	1.0000000000000e-10	
    area_system =
       1.999999977471133
    

    7.2 complexSimpson()

    function [ss,n]=complexSimpson(func,a,b,pre,m4)
    %func函数
    %a,b积分区间
    %pre精度
    %m4区间上4阶导的最大值
    %name:邓能财     Date: 2013/12/23
    %%%%%%%%%%%%%%%%%%%%%%:%【表需在尾部添加分号‘;’
    n=floor((sqrt((b-a)^5)*m4/(2880*pre))^(1/4))+1%【
    h=(b-a)/(2*n);
    ss=func(a)+4*func(a+h)+func(b);
    ss1=0; %n为奇数的项
    ss2=0; %n为偶数的项
    for xi=a+2*h:2*h:a+2*(n-1)*h
        ss1=ss1+func(xi+h);
        ss2=ss2+func(xi);
    end
    ss1,ss2
    ss=ss+4*ss1+2*ss2;
    ss=ss*(h/3);
    end
    

    7.3 complexTrangleInt()

    function [ss,n]=complexTrangleInt(func,a,b,pre,m2)
    %func函数
    %a,b积分区间
    %pre精度
    %m2区间上二阶导的最大值
    %name:邓能财     Date: 2013/12/23
    %%%%%%%%%%%%%%%%%%%%%%:%【表需在尾部添加分号‘;’
    n=floor((sqrt((b-a)^3)*m2/(12*pre))^.5)+1;%【
    h=(b-a)/n;
    ss=.5*(func(a)+func(b));
    for xi=a+h:h:a+(n-1)*h
        ss=ss+func(xi);
    end
    ss=ss*h;
    end
    

    8 自适应运算积分

    8.1 测试自适应运算积分函数selfAdaptIntegral()

    clc,clear
    format long
    % f=@(x)(exp(x))
    % a=1,b=10
    % pre=1e-13
    f =  @(x)(cos(x))
    a = -1.57079632679490
    b = 1.57079632679490
    % ss=complexTrangleInt(f,a,b,pre,m2)
    % [area,deep]=selfAdaptIntegral(f,a,b,pre)
    '各种精度的自适应Simpson积分'
    area_=[];pre_=[];deep_=[];
    for i=2:5
        pre=10^(-2*i);
        [area,deep]=selfAdaptIntegral(f,a,b,pre);
        pre_=[pre_,pre];area_=[area_,area];deep_=[deep_,deep];
    end,
    area_str=printArray(area_)
    pre_str=printArray(pre_)
    deep_
    area_system=quadl(f,a,b)
    

    输出数组的函数:

    function str=printArray(a)
    str=sprintf('%.13e\t',a);
    

    运行结果:

    f = 
        @(x)(cos(x))
    a =
      -1.570796326794900
    b =
       1.570796326794900
    ans =
    各种精度的自适应Simpson积分
    area_str =
    2.0000526243412e+00	2.0000002039922e+00	2.0000000017060e+00	2.0000000000498e+00	
    pre_str =
    1.0000000000000e-04	1.0000000000000e-06	1.0000000000000e-08	1.0000000000000e-10	
    deep_ =
         1     3     5     6
    area_system =
       1.999999977471133
    

    8.2 selfAdaptIntegral()

    function [area,deep]=selfAdaptIntegral(func,a,b,prec)
    %Simpson公式的自适应积分
    %理论依据:数值分析方法 奚梅成
    %name:邓能财     Date: 2013/12/23
    %func 函数;
    %[a,b] 积分区间
    %pre 精度
    %算法:递归
    %程序特点:密切结合理论,将内存和运算量达到极小化
    %% 计算
    %%disp('##Start')
    %取分段次数为n
    %n=3时,Simpson积分系数序列为: 1424241
    			%,已经能够发挥Simpson公式的作用,且分段数最少
    n=3;
    h=(b-a)/n;
    %计算梯形积分T_n
    s1=(func(a)+func(b))/2;
    for xi=a+h:h:b-h
       s1=s1+func(xi);
    end, s1=s1*h;
    %H_n(f)
    h1=0; 
    for xi=a+h/2:h:b-h/2
       h1=h1+func(xi);
    end, h1=h1*h;
    %计算梯形积分T_2n
    s2=(s1+h1)/2; 
    %计算Simpson积分S_n
    s1=(4*s2-s1)/3;
    %H_2n(f)
    h=h/2; h2=0; 
    for xi=a+h/2:h:b-h/2
       h2=h2+func(xi);
    end, h2=h2*h;
    %计算Simpson积分S_2n
    s2=s1/2+(4*h2-h1)/6;
    %% 判断是否达到要求精度
    deep=1;
    if(abs(s1-s2)<15*prec)
       area=s2;
       %disp('##End')
       return;
    else %不满足精度,进一步细分
        [area1,deep1]=selfAdaptIntegral(func,a,(a+b)/2,prec/2);
       	[area2,deep2]=selfAdaptIntegral(func,(a+b)/2,b,prec/2);
        area=area1+area2;
        deep=deep+max(deep1,deep2);
       %disp('##End')
       return;
    end
    end
    

    9 积分系数——二

    9.1 Gauss-Legendre积分系数

    格式:序号【 积分节点 ******积分系数
    1【
     0 ********* 2
    2【
     ±0.5773502692  ********* 1
    3【
     ±0.7745966692
    0
     *********
    0.5555555556
    0.8888888889
    4【
    ±0.8611363116
    ±0.3399810436
    *********
    0.3478548451
    0.6521451549
    5【
     ±0.9061798459
    ±0.5384693101
    0
     *********
    0.2369268851
    0.4786286705
    0.5688888889
    6【
    ±0.9324695142
    ±0.6612093865
    ±0.2386191861
     *********
    0.1713244924
    0.3607615730
    0.4679139346
    7【
    ±0.9491079123
    ±0.7415311856
    ±0.4058451514
    0
     *********
    0.1294849662
    0.2797053915
    0.3818300505
    0.4179591837
    8【
    ±0.9602898565
    ±0.7966664774
    ±0.5255324099
    ±0.1834346425
     *********
    0.1012285363
    0.2223810345
    0.3137066459
    0.3626837834
    

    9.2 newton cotes积分系数

    表6—1
    【】表错误
    n  Ci(n-1)
     
    2  1/2          1/2  
    3  1/6          4/6        1/6  
    4  1/8          3/8        3/8        1/8 
    5  7/90        16/45       2/15       16/45      7/90 
    6  19/288      25/96      25/144      25/144     25/96      19/288 
    7  41/840      9/35       9/280       34/105     9/280       9/35      41/840 
    8  751/17280 3577/17280 1323/17280  2989/17280 2989/17280  1323/17280  3577/17280  751/17280 
    9  989/28350 5888/28350 -928/28350 10496/28350 -4540/28350 10496/28350 -928/28350 【5888/28350】 989/28350 
    

    联系作者 definedone@163.com

    end

    展开全文
  • 本人收集的MATLAB基础教程-第8章 MATLAB数值积分与微分.ppt 第13章 在Word环境下使用MATLAB.ppt 第12章 Simulink动态仿真集成环境.ppt 第11章 MATLAB图形...
  • 【简介】 该教程旨在帮助数学建模和需要学习MATLAB的朋友学习MATLAB的相关知识,由浅入深...数学建模竞赛培训课程 第8章 MATLAB数值积分与微分(共14页).ppt 数学建模竞赛培训课程 第9章 MATLAB符号计算(共23页).ppt
  • Matlab数值积分P44.pdf

    2022-04-15 22:25:10
    Matlab数值积分P44.pdf
  • matlab数值积分函数

    2021-04-18 04:07:43
    相关函数:%符号积分int(f,v)int(f,v,a,b)%数值积分trapz(x,y)%梯形法沿列方向求函数Y关于自变量X的积分cumtrapz(x,y)%梯形法沿列方向求函数Y关于自变量X的累计积分quad(fun,a,b,tol)%采用递推自适应Simpson法计算...

    一.相关函数:

    %符号积分

    int(f,v)

    int(f,v,a,b)

    %数值积分

    trapz(x,y)%梯形法沿列方向求函数Y关于自变量X的积分

    cumtrapz(x,y)%梯形法沿列方向求函数Y关于自变量X的累计积分

    quad(fun,a,b,tol)%采用递推自适应Simpson法计算积分

    quad1(fun,a,b,tol)%采用递推自适应Lobatto法求数值积分

    dbquad(fun,xmin,xmax,ymin,ymax,zmin,zmax,tol)%二重(闭型)数值积分指令

    triplequad(fun,xmin,xmax,ymin,ymax,zmin,zmax,tol)%三重(闭型)数值积分指令

    二.示例:

    例1:计算f(t)=exp(-t^2)在[0,1]上的定积分

    本例演示:计算定积分常用方法

    >>symsx

    int(exp(-x^2),0,1)

    ans=

    1/2*erf(1)*pi^(1/2) %erf为误差函数

    >>vpa(int(exp(-x^2),0,1))

    ans=

    .7468241328124270

    >>d=0.001;x=0:d:1;d*trapz(exp(-x.^2))

    ans=

    0.7468

    >>quad('exp(-x.^2)',0,1,1e-8)

    ans=

    0.7468

    例2:计算f(t)=1/log(t)在[0,x],0

    注意:被积函数于x=0无义,在x-->1^-处为负无穷

    本例演示:用特殊函数表示的积分结果,如何用mfun指令

    (1)

    symstx

    ft=1/log(t);

    sx=int(ft,t,0,x)

    sx=

    -Ei(1,-log(x)) %完全椭圆函数

    (2)

    x=0.5:0.1:0.9

    sx_n=-mfun('Ei',1,-log(x))

    x=

    0.5000  0.6000  0.7000  0.8000  0.9000

    sx_n=

    -0.3787 -0.5469 -0.7809 -1.1340 -1.7758

    (3)%图示被函数和积分函数

    clf

    ezplot('1/log(t)',[0.1,0.9])

    gridon

    holdon

    plot(x,sx_n,'LineWidth',3)

    Char1='1/ln(t)';

    Char2='{int_0^x}1/ln(t)dt';

    title([Char1,' and  ',Char2])

    legend(Char1,Char2,'Location','SouthWest')

    例3:计算f(t)=exp(-sin(t))在[0,4]上的定积分

    注意:本题被函数之原函数无"封闭解析表达式",符号计算无法解题!

    本例演示:符号计算有限性

    (1)符号计算解法

    symstx

    ft=exp(-sin(t))

    sx=int(ft,t,0,4)

    ft=exp(-sin(t))

    Warning:Explicitintegralcouldnotbefound.

    >Insym.intat58

    sx=

    int(exp(-sin(t)),t=0..4)

    (2)数值计算解法

    dt=0.05;          %采样间隔

    t=0:dt:4;           %数值计算适合于有限区间上,取有限个采样点

    Ft=exp(-sin(t));

    Sx=dt*cumtrapz(Ft);      %计算区间内曲线下图形面积,为小矩形面积累加得

    Sx(end)        %所求定积分值

    %图示

    plot(t,Ft,'*r','MarkerSize',4)

    holdon

    plot(t,Sx,'.k','MarkerSize',15)

    holdoff

    xlabel('x')

    legend('Ft','Sx')

    >>ans=

    3.0632

    例4:绘制积分图形,y=2/3*exp(-t/2)*cos(sqrt(3)/2*t);积分s(x)=int(y,t,0,x)于[0,4*pi]上

    symsttao

    y=2/3*exp(-t/2)*cos(sqrt(3)/2*t);

    s=subs(int(y,t,0,tao),tao,t);  %获得积分函数

    subplot(2,1,1)

    %

    ezplot(y,[0,4*pi]),ylim([-0.2,0.7]) %单变量符号函数可视化,多变量用ezsurf

    gridon

    subplot(2,1,2)

    ezplot(s,[0,4*pi])

    gridon

    title('s=inty(t)dt')

    int的积分可以是定积分,也可以是不定积分(即有没有积分上下限都可以积)可以得到解析的解,比如你对x^2积分,得到的结果是1/3*x^3,这是通过解析的方法来解的。如果int(x^2,x,1,2)得到的结果是7/3

    quad是数值积分,它只能是定积分(就是有积分上下限的积分),它是通过simpson数值积分来求得的(并不是通过解析的方法得到解析解,再将上下限代入,而是用小梯形的面积求和得到的)。如果f=inline('x.^2');quad(f,1,2)得到的结果是2.333333,这个数并不是7/3

    %% 符号变量与符号表达式

    %%%%%%%%%%%%%%%%%%%%%%%%%%%

    %1.符号变量与符号表达式

    %%%%%%%%%%%%%%%%%%%%%%%%%%%

    clear all ;

    clc;

    close all;

    % f =sym( 'sin(x)+5x')

    % f —— 符号变量名

    % sin(x)+5x—— 符号表达式

    % ' '—— 符号标识

    % 符号表达式一定要用' ' 单引号括起来matlab才能识别

    % ' ' 的内容可以是符号表达式,也可以是符号方程。

    % 例: % f1=sym('a*x^2+b*x+c') —— 二次三项式

    % f2=sym('a*x^2+b*x+c=0' )—— 方程

    % f3=sym('Dy+y^2=1') ——微分方程

    % 符号表达式或符号方程可以赋给符号变量,以后调用方便;也可以不赋给符号变量直接参与运算

    % syms 命令用来建立多个符号量,一般调用格式为:

    % syms 变量1 变量2 ... 变量n %% 符号矩阵的创建

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    %2.符号矩阵的创建

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    % 数值矩阵A=[1,2;3,4]

    % A=[a,b;c,d] —— 不识别

    % @1.用matlab函数sym创建矩阵(symbolic的缩写)

    % 命令格式:A=sym('[ ]') % ※ 符号矩阵内容同数值矩阵

    % ※ 需用sym指令定义

    % ※ 需用' '标识

    % 例如:

    A = sym('[a , 2*b ; 3*a , 0]')

    % A =

    % [ a, 2*b]

    % [3*a, 0]

    % 这就完成了一个符号矩阵的创建。

    % 注意:符号矩阵的每一行的两端都有方括号,这是与 matlab数值矩阵的一个重要区别。

    %@2.用字符串直接创建矩阵(这种方法创建的没有什么用处)

    % ※模仿matlab数值矩阵的创建方法

    % ※需保证同一列中各元素字符串有相同的长度。

    % 例:

    A =['[ a,2*b]'; '[3*a, 0]'] % A =

    % [ a, 2*b]

    % [3*a, 0]

    %@3.符号矩阵的修改

    % a.直接修改

    % 可用光标键找到所要修改的矩阵,直接修改

    % b.指令修改

    % ※用A1=sym(A,*,*,'new') 来修改。 这个经过测试,不能运行

    % ※用A1=subs(A, 'new', 'old')来修改

    % % 例如:A =[ a, 2*b]

    % [3*a, 0]

    A = sym('[a , 2*b ; 3*a , 0]')

    % A1=sym(A,2,2,'4*b') %%等效于A(2,2)='4*b';

    % A1 =[ a, 2*b]

    % [3*a, 4*b] A1=subs(A,'0','4*b')

    A2=subs(A1, 'c', 'b') % A2 =[ a, 2*c] % [3*a, 4*c] %@4.符号矩阵与数值矩阵的转换

    % ※将数值矩阵转化为符号矩阵

    % 函数调用格式:sym(A)

    A=[1/3,2.5;1/0.7,2/5]

    % A =

    % 0.3333 2.5000

    % 1.4286 0.4000

    B=sym(A)

    % ans =

    % [ 1/3, 5/2]

    % [10/7, 2/5]

    % ※将符号矩阵转化为数值矩阵

    % 函数调用格式: numeric(A)

    % B =

    % [ 1/3, 5/2]

    % [10/7, 2/5]

    %numeric(B) 这个函数不存在了

    VPA(B,4) %发现这个函数可用

    % R = VPA(S) numerically evaluates each element of the double

    matrix

    % S using variable precision floating point arithmetic with D

    decimal % digit accuracy, where D is the current setting of

    DIGITS. % The resulting R is a SYM.

    %  % VPA(S,D) uses D digits, instead of the current setting of

    DIGITS.

    % D is an integer or the SYM representation of a number.

    % ans = % [ .3333, 2.500]

    % [ 1.429, .4000]

    %% 符号运算

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    %%3. 符号运算

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    % 例1:

    f=sym( '2*x^2+3*x-5'); g=sym( 'x^2+x-7');

    h= f+g

    % h= % 3*x^2+4*x-12

    % 例2:

    f=sym('cos(x)');g=sym('sin(2*x)');

    f/g+f*g

    % ans =

    % cos(x)/sin(2*x)+cos(x)*sin(2*x)

    %% 查找符号变量

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    %%4.查找符号变量

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    % % findsym(expr) 按字母顺序列出符号表达式 expr 中的所有符号变量

    % % findsym(expr, N) 列出 expr 中离 x 最近的 N 个符号变量

    % 若表达式中有两个符号变量与 x 的距离相等,则ASCII 码大者优先。

    % ※常量 pi, i, j 不作为符号变量

    % 例:

    f=sym('2*w-3*y+z^2+5*a');

    findsym(f)

    % ans =

    % a, w, y, z

    findsym(f,3)

    % ans =

    % y,w,z

    findsym(f,1)

    % ans =

    % y

    %% 计算极限

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    %%5.计算极限

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    % limit(f,x,a): 计算f(x)当x趋向于a的极限

    % limit(f,a): 当默认变量趋向于 a 时的极限

    % limit(f): 计算 a=0 时的极限

    % limit(f,x,a,'right'): 计算右极限

    % limit(f,x,a,'left'): 计算左极限

    % 例:计算

    syms x h n; L=limit((log(x+h)-log(x))/h,h,0)

    % L = % 1/x

    M=limit((1-x/n)^n,n,inf)

    % M = % exp(-x)

    %% 计算导数

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    %%6.计算导数

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    % g=diff(f,v):求符号表达式 f 关于 v 的导数

    % g=diff(f):求符号表达式 f 关于默认变量的导数

    % g=diff(f,v,n):求 f 关于 v 的 n 阶导数

    syms x;

    f=sin(x)+3*x^2; g=diff(f,x)

    % g = % cos(x)+6*x

    %%计算积分

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    %%7.计算积分

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    % int(f,v,a,b): 计算定积分f(v)从a到b

    % int(f,a,b): 计算关于默认变量的定积分

    % int(f,v): 计算不定积分f(v)

    % int(f): 计算关于默认变量的不定积分

    f=(x^2+1)/(x^2-2*x+2)^2;

    I=int(f,x)

    % I = % 3/2*atan(x-1)+1/4*(2*x-6)/(x^2-2*x+2)

    K=int(exp(-x^2),x,0,inf)

    % K = % 1/2*pi^(1/2)

    %%函数运算

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    %%8.函数运算

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    % 1.合并、化简、展开等函数

    % collect函数:将表达式中相同幂次的项合并;

    % factor函数:将表达式因式分解;

    % simplify函数:利用代数中的函数规则对表达式进行化简;

    % numden函数:将表示式从有理数形式转变成分子与分母形式。

    % 2.反函数

    % finverse(f,v) 对指定自变量为v的函数f(v)求反函数

    % 3.复合函数

    % compose(f,g) 求f=f(x)和g=g(y)的复合函数f(g(y))

    % compose(f,g,z) 求 f=f(x)和g=g(y)的复合函数f(g(z))

    % 4.表达式替换函数(前面讲到了)

    % subs(s) 用赋值语句中给定值替换表达式中所有同名变量

    % subs (s, old, new) 用符号或数值变量new替换s中的符号变量old

    %%

    % mtaylor(f,n) —— 泰勒级数展开

    % ztrans(f) —— Z变换

    % Invztrans(f) —— 反Z变换

    % Laplace(f) —— 拉氏变换

    % Invlaplace(f) —— 反拉氏变换

    % fourier(f) —— 付氏变换

    % Invfourier(f) —— 反付氏变换

    %%

    clear

    f1 =sym('(exp(x)+x)*(x+2)');

    f2 = sym('a^3-1');

    f3 = sym('1/a^4+2/a^3+3/a^2+4/a+5');

    f4 = sym('sin(x)^2+cos(x)^2');

    collect(f1)

    % ans =

    % x^2+(exp(x)+2)*x+2*exp(x)

    expand(f1)

    % ans = % exp(x)*x+2*exp(x)+x^2+2*x

    factor(f2)

    % ans = % (a-1)*(a^2+a+1)

    [m,n]=numden(f3)

    %m为分子,n为分母

    % m =

    % 1+2*a+3*a^2+4*a^3+5*a^4

    % n =

    % a^4

    simplify(f4)

    % ans =

    % 1

    clear

    syms x y

    finverse(1/tan(x)) %求反函数,自变量为x

    % ans =

    % atan(1/x)

    f = x^2+y;

    finverse(f,y) %求反函数,自变量为y

    % ans =

    % -x^2+y clear

    syms x y z t u;

    f = 1/(1 + x^2); g = sin(y); h = x^t; p = exp(-y/u);

    compose(f,g) %求f = f(x) 和 g = g(y)的复合函数f(g(y))

    % ans =

    % 1/(1+sin(y)^2) clear

    syms a b

    subs(a+b,a,4) %用4替代a+b中的a

    % ans =

    % 4+b

    subs(cos(a)+sin(b),{a,b},{sym('alpha'),2}) %多重替换

    % ans =

    % cos(alpha)+sin(2) f=sym('x^2+3*x+2')

    % f =

    % x^2+3*x+2

    subs(f, 'x', 2) %求解f当x=2时的值

    % ans =

    % 12

    %% 方程求解

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    %%9.方程求解

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    % 1代数方程

    % 代数方程的求解由函数solve实现:

    % solve(f) 求解符号方程式f % solve(f1,…,fn) 求解由f1,…,fn组成的代数方程组 % % 2常微分方程

    % 使用函数dsolve来求解常微分方程:

    % dsolve('eq1, eq2, ...', 'cond1, cond2, ...', 'v')

    clear

    syms a b c x

    f=sym('a*x*x+b*x+c=0')

    solve(f)

    % ans =

    % [ 1/2/a*(-b+(b^2-4*c*a)^(1/2))]

    % [ 1/2/a*(-b-(b^2-4*c*a)^(1/2))]

    solve('1+x=sin(x)')

    % ans =

    % -1.9345632107520242675632614537689

    dsolve( ' Dy=x ','x') %求微分方程y'=x的通解,指定x为自变量。

    % ans =

    % 1/2*x^2+C1

    dsolve(' D2y=1+Dy ','y(0)=1','Dy(0)=0' )

    %求微分方程y''=1+y'的解,加初始条件

    % ans = % -t+exp(t)

    [x,y]=dsolve('Dx=y+x,Dy=2*x') %微分方程组的通解

    % x =

    % -1/2*C1*exp(-t)+C2*exp(2*t) % y = % C1*exp(-t)+C2*exp(2*t)

    % ezplot(y)方程解y(t)的时间曲线图

    %% funtool

    funtool %该命令将生成三个图形窗口,Figure

    No.1用于显示函数f的图形,

    % Figure No.2用于显示函数g的图形,

    % Figure No.3为一可视化的、可操作与显示一元函数的计算器界面。

    % 在该界面上由许多按钮,可以显示两个由用户输入的函数的计算结果:

    % 加、乘、微分等。funtool还有一函数存储器,允许用户将函数存入,

    % 以便后面调用。在开始时,

    % funtool显示两个函数f(x) = x与g(x) = 1在区间[-2*pi, 2*pi]上的图形。

    % Funtool同时在下面显示一控制面板,

    % 允许用户对函数f、g进行保存、更正、重新输入、联合与转换等操作。

    %% taylortool %该命令生成一图形用户界面,显示缺省函数f=x*cos(x)

    % 在区间[-2*pi,2*pi]内的图形,同时显示函数f

    % 的前N=7项的Taylor多项式级数和(在a=0附近的)图形,

    % 通过更改f(x)项可得不同的函数图形。

    % taylortool('f') %对指定的函数f,用图形用户界面显示出Taylor展开式

    %% maple内核访问函数

    % % 可以访问maple内核的matlab函数:

    % maple ——— 访问maple内核函数

    % mapleinit —— maple函数初始化

    % mpa ———— maple函数定义

    % mhelp ——— maple函数帮助命令

    % procread —— maple函数程序安装

    % 具体的操作参看相关说明

    展开全文
  • MATLAB数值积分求值实验报告.pdfMATLAB数值积分求值实验报告.pdfMATLAB数值积分求值实验报告.pdfMATLAB数值积分求值实验报告.pdfMATLAB数值积分求值实验报告.pdfMATLAB数值积分求值实验报告.pdfMATLAB数值积分求值...
  • MATLAB数值积分求值实验报告.docxMATLAB数值积分求值实验报告.docxMATLAB数值积分求值实验报告.docxMATLAB数值积分求值实验报告.docxMATLAB数值积分求值实验报告.docxMATLAB数值积分求值实验报告.docxMATLAB数值积分...
  • matlab数值积分.ppt

    2021-10-31 22:57:55
    目录如下: §5.0 引言 §5.1 插值型求积公式和 Newton-Cotes求积公式 §5.2 变步长求积算法** §5.3 Romberg算法
  • Matlab数值积分法汇总,数值积分法MATLAB,C,C++源码.zip
  • 数值分析快速复习(1)——Matlab数值积分序言数值积分基础知识为什么要进行数值积分?数值积分的几种基本方法?插值型复化求积公式Romberg公式作业:Matlab实现Romberg积分 序言 由于最近的项目需要用到数值分析的...

    序言

    由于最近的项目需要用到数值分析的相关知识,需要简单回顾以下数值积分的知识,复习内容以上学期数值分析课本内容和课后Matlab作业为主,在这里将复习的感悟分享给大家。

    数值积分基础知识

    为什么要进行数值积分?

    在工程中,自变量和因变量往往不能通过简单的函数映射关系去表达,例如最简单的光/声波的传播,也会因为噪声而变得不能用简单的三角函数去表达。因此f(x)往往是通过列表函数表达,对于列表函数,求导和积分的公式都不能使用,这就要求我们要学会数值积分这个强有力的工具。

    另外,在高等数学里我们学过,一些函数根本不能积分,找不到初等函数的原函数,例如:1/ln(x),所以对于这些函数,也必须使用数值积分进行求解。

    数值积分的几种基本方法?

    插值型

    对于列表函数或一些复杂函数,我们一般可以知道一些节点的具体坐标,根据这些坐标拟合出一个多项式,这个过程叫做插值,得到拟合的多项式再进行积分就非常简单了。

    • 两个插值节点得到的求积公式——梯形公式
      在这里插入图片描述
      在这里插入图片描述

    • 三个等距插值节点得到的求积公式——Simpson公式

    在这里插入图片描述

    在这里插入图片描述

    • 五个等距插值节点得到的求积公式——cotes公式

    复化求积公式

    插值型求积公式在较小的积分区间中非常好用,但是随着积分区间的增大,截断误差也在不断增大,因此我们很自然的想到,如果将一个很大的积分区间划分成无限小段,再对每一小段使用前面的插值型求积公式,是不是就非常好了呢?

    于是有了:

    • 复化梯形公式
      在这里插入图片描述

    • 复化Simpson公式
      在这里插入图片描述

    • 复化Cotes公式
      在这里插入图片描述

    Romberg公式

    从上面可以看到,复化Simpson公式计算简单,但是精度相比于复化cotes来说比较差,收敛速度也相对较慢。怎么样发扬优点,避免缺点呢?

    数学家从复化梯形公式的误差估计式入手,推导梯形公式、Simpson公式、Cotes公式之间的关系,这样在计算Cotes公式时就不需要上面那么复杂的公式了,而是通过求梯形公式的值——Simposon公式的值——Cotes公式的值,这就是Romberg公式。
    在这里插入图片描述
    在这里插入图片描述
    在这里插入图片描述
    在这里插入图片描述

    所以只需要计算最简单的梯形公式就可以得到很精确地Romberg积分值了!

    (想看具体推导的同学,可以翻一下书本《数值分析(东南大学出版社)》215页)

    作业:Matlab实现Romberg积分

    贴一下我们老师布置的matlab作业的代码,实现了一个小的龙伯格积分,给大家作为参考。

    %%Romberg integration method
    clear,clc;format long;
    f=@(x)1/(1+100*x^2);
    syms x;            %定义函数
    e=0.5*10^(-7);
    for c=4:100
        for i=0:c
            n=2^i;
            h=2/n;
            x=zeros(1,n+1);
            x(1)=-1;
            x(n+1)=1;
            t=0;
            if i~=0
                for j=2:n
                    x(j)=-1+2/n*(j-1);
                end
                t=f(x(2));
                if i~=1
                    for j=3:n
                        t=t+f(x(j));
                    end
                end
            end
            T(1,i+1)=h/2*(f(-1)+f(1)+2*t);
        end
        for i=1:c
            T(2,i)=4/3*T(1,i+1)-T(1,i)/3;
        end
        for i=1:c-1
            T(3,i)=16/15*T(2,i+1)-T(2,i)/15;
        end
        for i=1:c-2
            T(4,i)=64/63*T(3,i+1)-T(3,i)/63;
        end
        if abs(T(4,i)-T(4,i-1))<e
                break
        end
    end
    disp('____________________________________________________')
    fprintf('区间等分数                   %f\n',n);
    fprintf('积分值                       %f\n',T(4,i));
    disp('____________________________________________________')
    
    
    

    以上内容如有误,还请各位批评指正。

    展开全文

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 11,193
精华内容 4,477
关键字:

matlab数值积分