2020年复化梯形公式和复化Simpson公式.pdf
(完整版)复化积分法(复化梯形求积,复化Simpson公式,变步长求积法)MATLAB编程实验报告
syms t;
y=sin(t)/t;
a=0;b=1;
R0=0.5*10^(-3);
TN(a,b,y,R0);
结论如下：
三、结论：
复化梯形求积：
复化Simpson公式
变步长求积法
h=(b-a)/n;
f=0;
for k=1:n
x(k)=a+h*(k-1);
w(k)=x(k)+h/2;
z(k)=subs(y,t,w(k));
f=f+z(k);
end
T(i)=T(i-1)/2+h/2*f
if ((T(i)-T(i-1))/3)<=R0
break;
end
end
tl=T(i)
clc,clear;
q=1;
end
p=subs(y,t,b);
T(1)=(b-a)/2*(q+p);
i=2;
n=i-1;
h=(b-a)/n;
z1=a+h/2;
z2=subs(y,t,z1);
T(2)=T(1)/2+h/2*z2;
while ((T(i)-T(i-1))/3)>R0
i=i+1
n=i-1;
n=2^(n-1)
x(k)=a+(k-1)*h
z(k)=subs(y,t,x(k));
end
for i=2:n
f=f+z(i);
end
q=subs(y,t,a);
if y=='sin(t)/t'&&a==0
q=1;
end
复合梯形公式和复合Simpson公式用法的比较
复合梯形公式和复合Simpson 公式用法的比较一、问题叙述由曲线y ()f x = ,a x b ≤≤绕x 轴旋转得到立体，其表面积计算公式为2(baS f x π=⎰分别用N=10的复合梯形公式和N=5的复合Simpson 公式求解下列每个的曲线绕x 轴旋转得到立方体的表面积。
（1）3(),01f x x x =≤≤ （2）()sin ,0/4f x x x π=≤≤ ； （3）()e ,01x f x x -=≤≤。
二、问题分析由叙述可知，三个函数在各自的定义域上函数值都大于0，且函数在各自的定义域上都是单调的。
（1）复合梯形公式取正整数n,令(b a)/n h =-,则点(k 0,1,,n)k x a kh =+=将[a,b]分为n 个小区间，利用定积分的性质可得11110011(x)(x)(x )f(x )2 =(a)f(b)2(x )2j jn n bx j j ax j j n j j h f dx f dx f h f f +--+==-=⎡⎤=≈+⎣⎦⎡⎤++⎢⎥⎣⎦∑∑⎰⎰∑所以，复合梯形公式的近似面积为：11(a)f(b)2(x )2n j j h S f f -=⎡⎤=++⎢⎥⎣⎦∑（2）复合Simpson 公式取正整数n,令h=(b-a)/2n,(j 0,1,,2n)j x a jh =+=.利用定积分的性质，以及Simpson 公式[]012(x)(x )4f(x )f(x )(f)3bahf dx f R =+++⎰在小区间222,j j x x -⎡⎤⎣⎦上的应用，有 2222221211(x)(x)(x )4f(x )f(x )3j nnbj j j j ax j j hf dx f dx f ---==⎡⎤=≈++⎣⎦∑∑⎰⎰整理可得，复合Simpson 公式的近似面积为：121211[f(a)f(b)4(x )2(x )]3n n j j j j hS f f --===+++∑∑。
复合梯形公式与复合辛普森公式对比
下面利用余项公式来估计两种算法的误差。首先需要求
f (x) sin x 1cos(xt)dt ,
x
0
所以有
的高阶导数。由于
于是
f k (x)
1 0
dk dxk
(cos
xt)dt
1tk cos(xt k )dt ,
0
2
max f k (x) 1 cos(xt k )tk dt 1tk dt 1 .
输出积分值 结束
图 4-1 复合梯形公式算法流程图 4.2 复合辛普森公式算法流程图
开始 输入区间断点 a, b
及等分数 n 求出步长各h节，各点节函点数值，相邻节求点和中su点m1, 及
各相相邻应节的点函中数心值输点出函结积数束,分值值
k=1求,2,和..n-s1um2
k=1,2,..n-1
图 4-2 复合辛普5.森M公A式TL算AB法源流程程序
1. 对多个算例进行分析，保证计算量基本相同的情况下去比较计算精度，验证复合辛普森
公式具有更高精度的结论。
2. 对多个算例进行 MATLAB 编程分析，在要求相同计算精度的前提下去比较计算量的大小，
从而分析复合梯形公式与复合辛普森公式的优劣。
参考文献
[1]穆耶赛尔·艾合买提,阿布都热西提·阿布都外力. 改进复合梯形求积公式[J]. 首都师范大学学报(自然科学 版),2016,06:1-4. [2] 刘 智 颖 , 王 向 公 , 任 威龙 , 龙 耀 萍 . 辛 普 森 公 式的 推 广 形 式 及 应 用 [J]. 山 东 理 工 大 学 学 报 ( 自 然科 学 版),2014,01:29-31. [3]陈凤忠. 辛普森公式的统一美[J]. 长春理工大学学报(综合版),2005,03:59-60+80. [4]李庆扬, 王能超, 易大义. 数值分析.第 5 版[M]. 清华大学出版社, 2008,12:106-109.
复化梯形公式和复化Simpson公式
数值计算方法上机题目3计算定积分的近似值：221x e xe dx =⎰ 要求：（1）若用复化梯形公式和复化Simpson 公式计算，要求误差限71021-⨯=ε，分别利用他们的余项估计对每种算法做出步长的事前估计；（2）分别利用复化梯形公式和复化Simpson 公式计算定积分；（3）将计算结果与精确解比较，并比较两种算法的计算量。
1.复化梯形公式程序：程序1（求f （x ）的n 阶导数：syms xf=x*exp(x) %定义函数f （x ）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n 阶导数结果1输入n=2f2 =2*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(2*exp(x) + x*exp(x))','x') %定义f(x)的二阶导数,输入程序1里求出的f2即可。
f3='-(2*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，以便求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的二阶导数的最小值点，也就是求二阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/12*((b-a)/n)^2*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hTn1=0for k=1:n-1 %求连加和xk=a+k*hTn1=Tn1+f(xk)endTn=h/2*((f(a)+2*Tn1+f(b)))z=exp(2)R=Tn-z %求已知值与计算值的差fprintf('用复化梯形算法计算的结果 Tn=')disp(Tn)fprintf('等分数 n=')disp(n) %输出等分数fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用复化梯形算法计算的结果 Tn= 7.3891等分数 n= 7019已知值与计算值的误差 R= 2.8300e-0082. Simpson公式程序：程序1：（求f（x）的n阶导数）：syms xf=x*exp(x) %定义函数f（x）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n阶导数结果1输入n=4f2 =4*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(4*exp(x) + x*exp(x))','x') %定义f(x)的四阶导数，输入程序1里求出的f2即可f3='-(4*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，一边求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的四阶导数的最小值点，也就是求四阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/180*((b-a)/(2*n))^4*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hSn1=0Sn2=0for k=0:n-1 %求两组连加和xk=a+k*hxk1=xk+h/2Sn1=Sn1+f(xk1)Sn2=Sn2+f(xk)endSn=h/6*(f(a)+4*Sn1+2*(Sn2-f(a))+f(b)) %因Sn2多加了k=0时的值，故减去f(a)z=exp(2)R=Sn-z %求已知值与计算值的差fprintf('用Simpson公式计算的结果 Sn=')disp(Sn)fprintf('等分数 n=')disp(n)fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用Simpson公式计算的结果 Sn= 7.3891等分数 n= 24已知值与计算值的误差 R= 2.7284e-008用复化梯形公式计算的结果为：7.3891，与精确解的误差为：2.8300e-008。
复化梯形公式和复化辛普生公式
1)复化梯形
#include<iostream.h>
#include<math.h>
class tixing
{
private:
int n,k;//n是等分数，
double h,f;//h步长
double a,b,integral,sum;//a上限，b下限，integral积分值
public:
《数值分析》实验报告
姓名
学号
日期
2012．11.20
实验室
设备编号
实验题目
用复化梯形公式和复化辛普生公式求I=∫01sinx/xdx
一实验目的
1．了解复化梯形公式和复化辛普生公式。
2.用复化梯形公式和复化辛普生公式求I=∫01sinx/xdx。
二实验内容
算法：复化梯形公式是Tn=∑h/2[f(xi)+ f(xi+1)]=(b-a)/2n[f(a)+2∑f(xi)+f(b)]记子段[xi,xi+1]的中点为xi+1/2,则复化Simpson公式为Sn=∑h/6[f(xi)+4f(xi+1/2)+ f(xi+1)]=b-a/6n[f(a)+4∑f(xi+1/2)+2f(xi)+f(b)]
double a,b,f,integral,interval;//a上限，b下限，integral积分值，interval每次应用法则的区间长度（2h）
double sum_odd_terms,x,sum_even_terms;//sum_even_terms偶数结点处的函数值之和，sum_odd_terms奇数结点处的函数值之和
数值分析复化Simpson积分公式和复化梯形积分公式计算积分的通用程序
数值分析第五次程序作业PB09001057 孙琪【问题】分别编写用复化Simpson 积分公式和复化梯形积分公式计算积分的通用程序；用如上程序计算积分： I (f )=∫sin (x )dx 40取节点x i ， i =0,…,N,N 为2k ,k =0,1,…,12，并分析误差；简单分析你得到的数据。
【复化Simpson 积分公式】Simpson 法则：∫f (x )dx ≈b −a 6[f (a )+4f (a +b 2)+f (b )]b a 使用偶数个子区间上的复合Simpson 法则：设n 是偶数，x i =a +ih , h =b−a n ,(0≤i ≤n) 则有∫f (x )dx =∫f (x )dx +∫f (x )dx +⋯+∫f (x )dx =∑∫f (x )dx x 2i x 2i−2n 2i=1x n x n−2x 4x 2x 2x 0b a 将Simpson 法则应用于每一个区间，得到复合Simpson 法则：∫f (x )dx ≈h 3b a [f (x 0)+2∑f (x 2i−2)n 2i=2+4∑f (x 2i−1)n 2i=1+f (x n )] 公式的误差项为：−1180(b −a )h 4f (4)(δ) 其中δ∈(a,b)【复化梯形积分公式】梯形法则：对两个节点相应的积分法则称为梯形法则：∫f (x )dx ≈b −a 2b a [f (a )+f (b )] 如果划分区间[a,b]为：a =x 0<x 1<⋯<x n =b那么在每个区间上可应用梯形法则，此时节点未必是等距的，由此得到复合梯形法则：∫f (x )dx =∑∫f (x )dx x i x i−1n i=1b a ≈12∑(x i −x i−1)[f (x i−1)+f (x i )]ni=1 对等间距h=(b-a)/n 及节点x i =a +ih ，复合梯形法则具有形式：∫f (x )dx ≈h 2[f (a )+2∑f (a +ih )n−1i=1+f (b )]b a 误差项为：−112(b −a )h 2f ′′(δ)【算法分析】复合Simpson 法则和复合梯形法则的算法上述描述中都已介绍了，在此不多做叙述。
复合梯形公式和复合辛普森公式
复合梯形公式和复合辛普森公式1.复合梯形公式步骤1：将积分区间[a,b]均匀地分成n个小区间，每个小区间的长度为h=(b-a)/n，其中n为正整数。
步骤2：定义一个函数f(x)，并在每个小区间上求出函数在小区间两个端点的函数值，记作f(x0), f(x1), f(x2), ..., f(xn)。
步骤3：根据梯形公式，每个小区间上的定积分近似值为(h/2) * (f(x0) + 2f(x1) + 2f(x2) + ... + 2f(xn-1) + f(xn))。
步骤4：将所有小区间上的定积分近似值相加，得到最终的近似值。
复合辛普森公式是通过将积分区间划分成若干个小区间，然后在每个小区间上应用辛普森公式来逼近定积分的值。
具体的步骤如下：步骤1：将积分区间[a,b]均匀地分成n个小区间，每个小区间的长度为h=(b-a)/n，其中n为正偶数。
步骤2：定义一个函数f(x)，并在每个小区间上求出函数在小区间三个点的函数值，记作f(x0),f(x1),f(x2)。
步骤3：根据辛普森公式，每个小区间上的定积分近似值为(h/3)*(f(x0)+4f(x1)+f(x2))。
步骤4：将所有小区间上的定积分近似值相加，得到最终的近似值。
复合辛普森公式的误差主要取决于小区间的数量和函数在每个小区间上的变化情况。
与复合梯形公式相比，复合辛普森公式的精度更高，但计算复杂度也更高。
综上所述，复合梯形公式和复合辛普森公式是数值积分中常用的近似计算方法。
它们通过将积分区间划分成小区间，并在每个小区间上应用相应的公式来逼近定积分的值。
这两种方法都可以通过增加小区间的数量来提高近似的精度，但相应地也会增加计算的复杂度。
根据实际情况，我们可以选择合适的方法来计算需要求解的定积分。
复化梯形公式和复化辛普森公式的精度比较(共7页)
- 让每一个人同等地提高(tí gāo)自我实验四、复化梯形(tīxíng)公式和复化Simpson 公式的精度(jīnɡ dù)比较（2 学时(xuéshí)）一、实验目的(mùdì)与要求1、熟习复化 Simpson 公式和复化梯形公式的结构原理；2、熟习并掌握两者的余项表达式；3、分别求出正确值，复化梯形的近似值，复化Simpson 的近似值，并比较后两者的精度；4、从余项表达式，即偏差曲线，来察看两者的精度，看哪个更靠近于正确值。
二、实验内容：关于函数 f (x) sin x，试利用下表计算积分 I 1sin xdx 。
x x表格以下：x 0 1/8 1/4 3/8 1/2 5/8 3/4 7/8 1 f ( x)1注：分别利用复化梯形公式和复化Simpson 公式计算，比较哪个精度更好。
此中：积分的正确值I0.9460831。
三、实验步骤1、熟习理论知识，并编写相应的程序；2、上机操作，从偏差图形上察看偏差，并与正确值对比较，看哪个精度更好；3、得出结论，并整理实验报告。
四、实验注意事项11- 让每一个人同等地提高(tí gāo)自我1、复化梯形公式，程序主体(zhǔtǐ)部分：for n=2:10T(n)=*T(n-1)for i=1:2^(n-2)；endend2、复化 Simpson 公式(gōngshì)，程序主体部分：for i=1:10n=2.^ix=0:1/n:1f=sin(x)./xf(1)=1s=0for j=1:n/2s=s+f(2*j)endt=0for j=1:(n/2-1)t=t+f(2*j-1)endS(i)=1/3/n*(f(1)+4*s+2*t+f(n+1))end五．实验(shíyàn)内容复化梯形(tīxíng)公式和复化辛普森公式的引入22- 让每一个(yī ɡè)人同等地提高自我复化梯形(tīxíng)公式：T n n 1 h [ f (x k f ( x k 1 )] ；k 0 2复化辛普森公式：S n n 1 h4 f ( x 1 ) f (x k 1 )] ；[ f (x kk 0 6k2依据题意(tí yì)和复化梯形公式、复化辛普森公式的原理编写程序求解代码以下：Matlab 代码(dài mǎ)clcs=quad( 'sin(x)./x',0,1)p1=zeros(10,1);p2=zeros(10,1);for k=6:15s1=0;s2=0;x=linspace(0,1,k);y=sin(x)./x;z=(1/(2*(k-1))):(1/(k-1)):1;sz=sin(z)./z;y(1)=1;for i=1:(k-1)s1=s1+*(x(i+1)-x(i))*(y(i)+y(i+1));endfor j=1:(k-1)s2=s2+(1/6)*(x(j+1)-x(j))*(y(j)+y(j+1)+4*sz(j));endp1(k-5)=s1-s;p2(k-5)=s2-s;endp1;33- 让每一个人同等地提高(tí gāo)自我p2;s1=s+p1(4)s2=s+p2(4)format longfor k=1:length(p1)p1(k)=abs(p1(k));p2(k)=abs(p2(k));endp1p2plot(6:1:15,p1,'-r' )hold onplot(6:1:15,10000*(p2),'-c' )hold off部分(bù fen)程序结果输出：s=s1 =s2 =结果剖析1 sin( x)；依据结果输出可知：积分 I dx 的正确值为： I=x44- 让每一个人同等(tóngděng)地提高自我经过复化梯形公式(gōngshì)和复化辛普森公式获得的积分值为：s1 =:s2 =;相对偏差为：S 1 I 100 4.15 10 4 ； 1I 100S 2 I 100 1.62 10 8 ； 2 I 100明显，从相对偏差可知经过辛普森公式获得(huòdé)的结果偏差小精度高。
复合梯形公式与复合辛普森公式对比
SHANGHAI JIAO TONG UNIVERSITY 题目名称：复合梯形公式与复合辛普森公式对比学生姓名:学生学号:班级:学院(系):目录1.概述 (3)2.问题提出 (4)3.算法推导 (5)4.算法框图 (6)复合梯形公式算法流程图 (6)复合辛普森公式算法流程图 (7)5.MATLAB源程序 (8)6.结论与展望 (9)图表目录图 4-1 复合梯形公式算法流程图 (6)图 4-2 复合辛普森公式算法流程图 (7)图 6-1 MATLAB计算结果 (9)表 2-1函数计算结果表 (4)1.概述梯形求积公式和辛普森求积公式分别是牛顿-科斯特公式中n=1和n=2时的情形。
其中梯形求积公式可表示为其公式左端是以[a,b]区间上积分，右端为b-a为高、端点函数值为上下底的梯形的面积值，故通称为梯形公式，具有1次代数精确度。
类似的，辛普森求积公式可以表示为该公式一般在立体几何中用来求拟柱体的体积，由于偶数n阶牛顿-科特斯求积公式至少具有n+1次代数精确度，所以辛普森公式实际上具有3次代数精确度。
由于牛顿-科斯特公式在n≥8时不具有稳定性，故不可能通过提高阶的方法来提高求积精度。
为了提高精度通常可把积分区间分成若干子区间（通常是等分），再在每个子区间上用低阶求积公式。
这种方法称为复合求积法。
本文主要讨论复合梯形公式和复合辛普森公式在同一数学问题中的应用。
首先给出了复合梯形公式和复合辛普森公式的推导过程以及其余项的表达形式，然后用流程图的形式介绍算法思路，再运用MATLAB编写代码计算结果，最后对结果进行对比讨论。
希望通过两个算法在同一个算例中的应用对比，更好的理解和掌握复合梯形公式和复合辛普森公式的适用范围和适用条件。
并且能够熟悉MATLAB编程求解问题的流程，掌握编程化的思想方法。
同时对两种方法的计算结果对比分析，讨论两种求积方法的计算精度。
2.问题提出对于函数f(x)=sinxx给出的函数表如下，试用复合梯形公式和复合辛普森公式计算积分I=∫sinxx dx1。
复化梯形法 复化矩形法 变步长梯形 变步长辛普森
陕西科技大学机械教改班用C++的积分其实积分的思想就是,微分—>求和—>取极限，如果是用纯手工法那就是先对一个函数微分，再求出它的面积，在取极限，因为我们的计算速度和计算量有限，现在有了计算机这个速度很快的机器，我们可以把微分后的每个小的面积加起来，为了满足精度，我们可以加大分区，即使实现不了微分出无限小的极限情况，我们也至少可以用有限次去接近他，下面我分析了四种不同的积分方法，和一个综合通用程序。
一．积分的基本思想1、思路：微分—>求和—>取极限。
2、Newton —Leibniz 公式 ⎰-=ba a Fb F dx x f )()()( 其中，)(x F 被积函数)(x f的原函数。
3、用计算机积分的思路在积分区间内“微分—>求和—>控制精度”。
因为计算机求和不可以取极限，也就是不可以无限次的加下去，所以要控制精度。
二．现有的理论1、一阶求积公式---梯形公式⎰=+-=b a T b f a f a b dx x f )]()([2)( 他只能精确计算被积函数为0、1次多项式时的积分。
2、二阶求积分公式——牛顿、科特斯公式 ⎰=+++-=ba Sb f a b f a f a b dx x f )]()2(4)([6)(他只能精确计算被积函数为0、1、2、3次多项式时的积分。
三．四种实现方法1.复化矩形法将积分区间[a,b]等分成n 个子区间:],[],[],[],[],[112322110n n n n x x x x x x x x x x ---、、、 则h=(b-a)/n,区间端点值k x =a+kh)hf(x ))f(x x (x I 11121=-=)()()x (22232x hf x f x I =-=............................)()()(111n ---=-=n n n n x hf x f x x I∑==ni i x hf T 1n )(源程序：#include <iostream.h>#include<math.h>double f(double x) //计算被积函数{double y;y=log(1+x)/(1+x*x); //被积函数return y;}double Tn(double a,double b,int n) //求Tn{double t=0.0;double xk; //区间端点值double t1,t2; //用来判断精度do{double h=(b-a)/n;for(int k=1;k<=n-1;k++) //每一小段的矩形叠加 {t1=t;xk=a+k*h;t+=h*f(xk);t2=t;}n++; //如果精度不够就对区间再次细分，直到达到精度要求 }while(fabs(t1-t2)<=1e-7); //判断计算精度return t;}void main(){double a=0.0; //积分下线double b=2.0; //积分上限int n=1024; //把区间分为1024段cout<<Tn(a,b,n)<<endl; //输出积分结果}执行结果：2.复化梯形法方法和复化矩形法类似，只是把原来的矩形小面积变成了梯形小面积，但是精确度明显提高了，也就是说达到同样的精度需要的时间少了。
复化梯形公式和复化辛普森公式
复化梯形公式和复化辛普森公式1. 引言嘿，大家好！今天我们来聊聊数学里那些看似高深莫测的公式，尤其是复化梯形公式和复化辛普森公式。
这些名字听起来就像是从某部科幻片里蹦出来的角色，但其实它们是我们在数值积分中不可或缺的好帮手。
你知道吗？它们就像是数学世界里的“超能英雄”，让我们轻松搞定积分，简直是妙不可言。
2. 复化梯形公式2.1 你知道什么是梯形吗？首先，咱们得聊聊复化梯形公式。
说白了，就是把一个复杂的积分任务，分解成几个小的梯形来求解。
想象一下，你在河边钓鱼，河水流得可欢了。
为了找一个合适的钓鱼点，你可能得把河分成几段，然后每一段的宽度就是你的小梯形。
你看，这就是复化梯形的魅力所在！2.2 如何运用复化梯形公式？用这个公式的时候，你只需把整个区间分成N个小区间，每个区间的宽度都是一样的。
然后，把每个小区间的函数值拿来加一加，再乘上宽度的一半，最后再把头尾的函数值加上。
这听起来是不是很简单？比如，你想算从0到1的某个函数的积分，只要把这个区间分成若干段，像切蛋糕一样，每一片都求个函数值，然后把结果合起来就行了。
简单得就像吃个冰淇淋，大家都喜欢。
3. 复化辛普森公式3.1 辛普森是谁？接下来，让我们来看看复化辛普森公式。
辛普森这个名字，大家可能都听过，或者说过“这是辛普森家的事儿”。
其实，他是一位牛逼的数学家，专门研究如何让积分变得更加简单。
辛普森公式就像是对梯形公式的一次升级，像换了个新款手机，功能更强大，效果更好。
3.2 如何运用复化辛普森公式？用复化辛普森公式的时候，我们也是把整个区间分成N个小区间，不过这里的N必须是偶数哦！每个小区间的宽度仍然是一样的。
然后，用函数值的加权平均法来计算。
换句话说，你把每个小区间的头尾和中间的函数值都考虑进来，像是为你的冰淇淋加上各种口味的配料。
最后，你的结果就会比单纯用梯形公式得来的要精准多了，仿佛一口下去，味蕾都在舞蹈。
4. 比较与应用4.1 谁更强？说到这儿，很多人就会问，复化梯形公式和复化辛普森公式，谁更厉害呢？其实，这就像问“苹果和橘子，哪个更好吃”。
复化梯形公式和复化Simpson公式
数值计算方法上机题目3一、计算定积分的近似值：221x e xe dx =⎰ 要求：（1）若用复化梯形公式和复化Simpson 公式计算，要求误差限71021-⨯=ε，分别利用他们的余项估计对每种算法做出步长的事前估计；（2）分别利用复化梯形公式和复化Simpson 公式计算定积分；（3）将计算结果与精确解比较，并比较两种算法的计算量。
1.复化梯形公式程序：程序1（求f （x ）的n 阶导数：syms xf=x*exp(x) %定义函数f （x ）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n 阶导数结果1输入n=2f2 =2*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(2*exp(x) + x*exp(x))','x') %定义f(x)的二阶导数,输入程序1里求出的f2即可。
f3='-(2*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，以便求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的二阶导数的最小值点，也就是求二阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/12*((b-a)/n)^2*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hTn1=0for k=1:n-1 %求连加和xk=a+k*hTn1=Tn1+f(xk)endTn=h/2*((f(a)+2*Tn1+f(b)))z=exp(2)R=Tn-z %求已知值与计算值的差fprintf('用复化梯形算法计算的结果 Tn=')disp(Tn)fprintf('等分数 n=')disp(n) %输出等分数fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用复化梯形算法计算的结果Tn= 7.3891等分数n=7019已知值与计算值的误差R= 2.8300e-0082. Simpson公式程序：程序1：（求f（x）的n阶导数）：syms xf=x*exp(x) %定义函数f（x）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n阶导数结果1输入n=4f2 =4*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(4*exp(x) + x*exp(x))','x') %定义f(x)的四阶导数，输入程序1里求出的f2即可f3='-(4*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，一边求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的四阶导数的最小值点，也就是求四阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/180*((b-a)/(2*n))^4*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hSn1=0Sn2=0for k=0:n-1 %求两组连加和xk=a+k*hxk1=xk+h/2Sn1=Sn1+f(xk1)Sn2=Sn2+f(xk)endSn=h/6*(f(a)+4*Sn1+2*(Sn2-f(a))+f(b)) %因Sn2多加了k=0时的值，故减去f(a)z=exp(2)R=Sn-z %求已知值与计算值的差fprintf('用Simpson公式计算的结果 Sn=')disp(Sn)fprintf('等分数 n=')disp(n)fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用Simpson公式计算的结果Sn= 7.3891等分数n=24已知值与计算值的误差R= 2.7284e-008用复化梯形公式计算的结果为：7.3891，与精确解的误差为：2.8300e-008。
复化梯形求积公式
复化梯形求积公式复化梯形求积法是一种用于数值积分的常见方法，它可以帮助我们求解一元定积分（即 integrals 或者是积分函数）。
在复化梯形求积法中，我们需要对函数 f(x) 在某一区间 [a,b]内的值进行插值，并将插值后的函数值相加，以计算该函数在该区间上的积分值。
复化梯形求积法的步骤如下：第一步：分解积分区间[a,b]，将其分解为 n 个子区间，其中区间长度为 h=(b-a)/n。
第二步：对积分区间[a,b] 中每个子区间使用复化梯形公式来计算它们的积分值。
根据复化梯形求积公式，积分区间 [a,b] 中每个子区间的积分值如下：Ii = h/2 * [f(xi) + f(xi+1)] (i=0,1,...,n-1)其中，xi 是第 i 个（0≤i≤n-1）子区间的左端点，xi+1 是第 i 个子区间的右端点，f(xi) 是第 i 个子区间左端点处的函数值，f(xi+1) 是第 i 个子区间右端点处的函数值。
第三步：将所有子区间积分值进行累加，得到积分区间[a,b]的总积分值。
按照复化梯形求积法，积分区间 [a,b] 的总积分值计算公式如下：I=h/2 * [f(x0)+f(xn)+∑(i=1,2,...,n-1)[f(xi)+f(xi+1)]]以上就是复化梯形求积法的基本原理及其用于求解一元定积分的步骤：首先将积分区间分解为子区间，然后使用复化梯形公式计算每个子区间的积分值，最后将子区间积分值相加得到积分区间的总积分值。
这种求积方法容易理解，计算量小，计算速度快，因此复化梯形求积法在积分计算中有着重要的应用价值。
在复化梯形求积法中，可以使用高斯公式来加快积分值的计算。
高斯公式的定义如下：I=h/2 * [f(x0)+f(xn)+c1f(x1)+c2f(x2)+...+cnf(xn-1)]其中，c1,...,cn 是某些常数，x1,...,xn-1 为积分区间上的中点。
应用高斯公式，可以将积分计算量减小一半以上，同时也大大减少计算精度上的误差。
复化梯形公式和复化Simpson公式
数值计算方法上机题目3一、计算定积分的近似值：221x e xe dx =⎰ 要求：（1）若用复化梯形公式和复化Simpson 公式计算，要求误差限71021-⨯=ε，分别利用他们的余项估计对每种算法做出步长的事前估计；（2）分别利用复化梯形公式和复化Simpson 公式计算定积分；（3）将计算结果与精确解比较，并比较两种算法的计算量。
1.复化梯形公式程序：程序1（求f （x ）的n 阶导数：syms xf=x*exp(x) %定义函数f （x ）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n 阶导数结果1输入n=2f2 =2*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(2*exp(x) + x*exp(x))','x') %定义f(x)的二阶导数,输入程序1里求出的f2即可。
f3='-(2*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，以便求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的二阶导数的最小值点，也就是求二阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/12*((b-a)/n)^2*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hTn1=0for k=1:n-1 %求连加和xk=a+k*hTn1=Tn1+f(xk)endTn=h/2*((f(a)+2*Tn1+f(b)))z=exp(2)R=Tn-z %求已知值与计算值的差fprintf('用复化梯形算法计算的结果 Tn=')disp(Tn)fprintf('等分数 n=')disp(n) %输出等分数fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用复化梯形算法计算的结果 Tn= 7.3891等分数 n=7019已知值与计算值的误差 R= 2.8300e-0082. Simpson公式程序：程序1：（求f（x）的n阶导数）：syms xf=x*exp(x) %定义函数f（x）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n阶导数结果1输入n=4f2 =4*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(4*exp(x) + x*exp(x))','x') %定义f(x)的四阶导数，输入程序1里求出的f2即可f3='-(4*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，一边求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的四阶导数的最小值点，也就是求四阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/180*((b-a)/(2*n))^4*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hSn1=0Sn2=0for k=0:n-1 %求两组连加和xk=a+k*hxk1=xk+h/2Sn1=Sn1+f(xk1)Sn2=Sn2+f(xk)endSn=h/6*(f(a)+4*Sn1+2*(Sn2-f(a))+f(b)) %因Sn2多加了k=0时的值，故减去f(a)z=exp(2)R=Sn-z %求已知值与计算值的差fprintf('用Simpson公式计算的结果 Sn=')disp(Sn)fprintf('等分数 n=')disp(n)fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用Simpson公式计算的结果 Sn= 7.3891等分数 n=24已知值与计算值的误差 R= 2.7284e-008用复化梯形公式计算的结果为：7.3891，与精确解的误差为：2.8300e-008。
复化梯形公式和复化Simpson公式
数值计算方法上机题目3一、计算定积分的近似值：要求：（1）若用复化梯形公式和复化Simpson 公式计算，要求误差限71021-⨯=ε，分别利用他们的余项估计对每种算法做出步长的事前估计；（2）分别利用复化梯形公式和复化Simpson 公式计算定积分；（3）将计算结果与精确解比较，并比较两种算法的计算量。
1.复化梯形公式程序：程序1（求f （x ）的n 阶导数：syms xf=x*exp(x) %定义函数f （x ）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n 阶导数 结果1输入n=2f2 =2*exp(x) + x*exp(x)程序2：clcclearsyms x %定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可 f2=inline('(2*exp(x) + x*exp(x))','x') %定义f(x)的二阶导数,输入程序1里求出的f2即可。
f3='-(2*exp(x) + x*exp(x))' %因fminbnd （）函数求的是表达式的最小值，且要求表达式带引号，故取负号，以便求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的二阶导数的最小值点，也就是求二阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/12*((b-a)/n)^2*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hTn1=0for k=1:n-1 %求连加和xk=a+k*hTn1=Tn1+f(xk)endTn=h/2*((f(a)+2*Tn1+f(b)))z=exp(2)R=Tn-z %求已知值与计算值的差fprintf('用复化梯形算法计算的结果 Tn=')disp(Tn)fprintf('等分数 n=')disp(n) %输出等分数fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用复化梯形算法计算的结果Tn= 7.3891等分数n=7019已知值与计算值的误差R= 2.8300e-0082. Simpson公式程序：程序1：（求f（x）的n阶导数）：syms xf=x*exp(x) %定义函数f（x）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n阶导数结果1输入n=4f2 =4*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(4*exp(x) + x*exp(x))','x') %定义f(x)的四阶导数，输入程序1里求出的f2即可f3='-(4*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，一边求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的四阶导数的最小值点，也就是求四阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/180*((b-a)/(2*n))^4*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hSn1=0Sn2=0for k=0:n-1 %求两组连加和xk=a+k*hxk1=xk+h/2Sn1=Sn1+f(xk1)Sn2=Sn2+f(xk)endSn=h/6*(f(a)+4*Sn1+2*(Sn2-f(a))+f(b)) %因Sn2多加了k=0时的值，故减去f(a)z=exp(2)R=Sn-z %求已知值与计算值的差fprintf('用Simpson公式计算的结果 Sn=')disp(Sn)fprintf('等分数 n=')disp(n)fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用Simpson公式计算的结果Sn= 7.3891等分数n=24已知值与计算值的误差R= 2.7284e-008用复化梯形公式计算的结果为：7.3891，与精确解的误差为：2.8300e-008。
利用复化梯形公式、复化simpson 公式计算积分
教师签名：
年月日
上图是利用复化simpson公式所画出的误差。其中：红线是计算误差，‘＊’号是实际误差。
注：纵轴是0.0001。
0.5208、0.1029、0.0326、0.0133、0.0064、0.0035、0.0020、0.0013、0.0008是n值分别为2到10的实际误差，-0.0083是计算误差。
成
绩
评
注：由于被积函数是x.^4，它的四阶倒数是24，所以它的积分余项是：(-(b-a)/180*((b-a)/2).^4*24)
实
验
结
果
分
析
及
心
得
体
会
上图是利用复化梯形公式所画出的误差。其中：红线是计算误差，‘＊’号是实际误差。-0.0017是计算误差。0.0417、0.0185、0.0104、0.0067 0.0046、0.0034、0.0026、0.0021、0.0017是n值分别为2到10的实际误差。
x=linspace(0,1,9);
plot(x,A,'r',x,R,'*')％将计算误差与实际误差用图像画出来
注：由于被积函数是x.^2，它的二阶倒数为2，所以积分余项为：(-(b-a)/12*h.^ 2*2)
实
验
原
理
（
算
法
流
程
图
或
者
含
注
释
的
源
代
码
）
利用复化simpson公式的程序代码如下：
同样首先建立被积函数的函数文件：
f=x.^2;％首先建立被积函数，以便于计算真实值。
a=0;％积分下线
b=1;％积分上线
复化梯形公式和复化Simpson公式
数值计算方法上机题目3一、计算定积分的近似值：221x e xe dx =⎰ 要求：（1）若用复化梯形公式和复化Simpson 公式计算，要求误差限71021-⨯=ε，分别利用他们的余项估计对每种算法做出步长的事前估计；（2）分别利用复化梯形公式和复化Simpson 公式计算定积分；（3）将计算结果与精确解比较，并比较两种算法的计算量。
1.复化梯形公式程序：程序1（求f （x ）的n 阶导数：syms xf=x*exp(x) %定义函数f （x ）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n 阶导数结果1输入n=2f2 =2*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(2*exp(x) + x*exp(x))','x') %定义f(x)的二阶导数,输入程序1里求出的f2即可。
f3='-(2*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，以便求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的二阶导数的最小值点，也就是求二阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/12*((b-a)/n)^2*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hTn1=0for k=1:n-1 %求连加和xk=a+k*hTn1=Tn1+f(xk)endTn=h/2*((f(a)+2*Tn1+f(b)))z=exp(2)R=Tn-z %求已知值与计算值的差fprintf('用复化梯形算法计算的结果 Tn=')disp(Tn)fprintf('等分数 n=')disp(n) %输出等分数fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用复化梯形算法计算的结果Tn= 7.3891等分数n=7019已知值与计算值的误差R= 2.8300e-0082. Simpson公式程序：程序1：（求f（x）的n阶导数）：syms xf=x*exp(x) %定义函数f（x）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n阶导数结果1输入n=4f2 =4*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(4*exp(x) + x*exp(x))','x') %定义f(x)的四阶导数，输入程序1里求出的f2即可f3='-(4*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，一边求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的四阶导数的最小值点，也就是求四阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/180*((b-a)/(2*n))^4*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hSn1=0Sn2=0for k=0:n-1 %求两组连加和xk=a+k*hxk1=xk+h/2Sn1=Sn1+f(xk1)Sn2=Sn2+f(xk)endSn=h/6*(f(a)+4*Sn1+2*(Sn2-f(a))+f(b)) %因Sn2多加了k=0时的值，故减去f(a)z=exp(2)R=Sn-z %求已知值与计算值的差fprintf('用Simpson公式计算的结果 Sn=')disp(Sn)fprintf('等分数 n=')disp(n)fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用Simpson公式计算的结果Sn= 7.3891等分数n=24已知值与计算值的误差R= 2.7284e-008用复化梯形公式计算的结果为：7.3891，与精确解的误差为：2.8300e-008。
数值计算考题五——复合梯形求积公式与复合辛普森求积公式求积分
数值计算考题五1. 分别用复合梯形求积公式与复合辛普森求积公式求积分I=⎰102x e sinx dx 的近似值，要求误差不超过ε=0.5⨯10-5.解：方法一： 复合梯形求积公式复合梯形求积公式是将积分区间划分为n 个很小的区间，然后将各个小区间的面积相加而得到在整个积分区间上的积分，当分成的小区间数n →∞时，求得的面积就等于积分的精确值。
由复合梯形求积公式的余项R n T 可得满足精度要求≤ε0.5⨯10-5时区间()b a ,被分成的区间数n 的最小值为700，所以在编程时循环次数应大于等于这个值，方可满足精度要求。
以下是编写的C 语言程序：#include<stdio.h>#include<math.h>void main(){int n=700,i;double x,f=0.0,t,h,T=0.0,c=2.0,a=0.0,b=1.0;h=(b-a)/n;for(i=0;i<n;i++){x=a+i*h;f=f+exp(pow(x,c))*sin(x);}t=(h/2)*(2*f+sin(1)*exp(1));printf("T=%f\n",t);}输出结果为T=0.778746.方法二：复合辛普森求积公式：复合辛普森求积法是将积分区间分割之后，在每个小区间[x i ,x i+1]上运用辛普森求积公式。
以下是编写的c 语言程序：#include<stdio.h>#include<math.h>void main(){int n=700,i;double x1,x2,f1=0.0,f2=0.0,t,h,T=0.0,c=2.0,a=0.0,b=1.0;h=(b-a)/n;for(i=0;i<n;i++){x1=a+i*h;x2=a+(i+0.5)*h;f1=f1+exp(pow(x1,c))*sin(x1);f2=f2+exp(pow(x2,c))*sin(x2); }t=(h/6)*(2*f1+sin(1)*exp(1)+4*f2); printf("T=%f\n",t);}程序输出结果为0.778745.2. 用高斯求积法求上述积分的近似值。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
数值计算方法上机题目3一、计算定积分的近似值：221x e xe dx =⎰ 要求：（1）若用复化梯形公式和复化Simpson 公式计算，要求误差限71021−⨯=ε，分别利用他们的余项估计对每种算法做出步长的事前估计；（2）分别利用复化梯形公式和复化Simpson 公式计算定积分；（3）将计算结果与精确解比较，并比较两种算法的计算量。
1.复化梯形公式程序：程序1（求f （x ）的n 阶导数：syms xf=x*exp(x) %定义函数f （x ）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n 阶导数结果1输入n=2f2 =2*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(2*exp(x) + x*exp(x))','x') %定义f(x)的二阶导数,输入程序1里求出的f2即可。
f3='-(2*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，以便求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的二阶导数的最小值点，也就是求二阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/12*((b-a)/n)^2*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hTn1=0for k=1:n-1 %求连加和xk=a+k*hTn1=Tn1+f(xk)endTn=h/2*((f(a)+2*Tn1+f(b)))z=exp(2)R=Tn-z %求已知值与计算值的差fprintf('用复化梯形算法计算的结果 Tn=')disp(Tn)fprintf('等分数 n=')disp(n) %输出等分数fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用复化梯形算法计算的结果Tn= 7.3891等分数n=7019已知值与计算值的误差R= 2.8300e-0082. Simpson公式程序：程序1：（求f（x）的n阶导数）：syms xf=x*exp(x) %定义函数f（x）n=input('输入所求导数阶数:')f2=diff(f,x,n) %求f(x)的n阶导数结果1输入n=4f2 =4*exp(x) + x*exp(x)程序2：clcclearsyms x%定义自变量xf=inline('x*exp(x)','x') %定义函数f(x)=x*exp(x)，换函数时只需换该函数表达式即可f2=inline('(4*exp(x) + x*exp(x))','x') %定义f(x)的四阶导数，输入程序1里求出的f2即可f3='-(4*exp(x) + x*exp(x))'%因fminbnd（）函数求的是表达式的最小值，且要求表达式带引号，故取负号，一边求最大值e=5*10^(-8) %精度要求值a=1 %积分下限b=2 %积分上限x1=fminbnd(f3,1,2) %求负的四阶导数的最小值点，也就是求四阶导数的最大值点对应的x值for n=2:1000000 %求等分数nRn=-(b-a)/180*((b-a)/(2*n))^4*f2(x1) %计算余项if abs(Rn)<e %用余项进行判断break% 符合要求时结束endendh=(b-a)/n %求hSn1=0Sn2=0for k=0:n-1 %求两组连加和xk=a+k*hxk1=xk+h/2Sn1=Sn1+f(xk1)Sn2=Sn2+f(xk)endSn=h/6*(f(a)+4*Sn1+2*(Sn2-f(a))+f(b)) %因Sn2多加了k=0时的值，故减去f(a)z=exp(2)R=Sn-z %求已知值与计算值的差fprintf('用Simpson公式计算的结果 Sn=')disp(Sn)fprintf('等分数 n=')disp(n)fprintf('已知值与计算值的误差 R=')disp(R)输出结果显示：用Simpson公式计算的结果Sn= 7.3891等分数n=24已知值与计算值的误差R= 2.7284e-008用复化梯形公式计算的结果为：7.3891，与精确解的误差为：2.8300e-008。
等分数n=7019用复化Simpson公式计算的结果为：7.3891，与精确解的误差为：2.7284e-008。
等分数n=243、柯斯特公式求积分：程序代码：（1）function [y,Ck,Ak]=NewtonCotes(fun,a,b,n)if nargin==1[mm,nn]=size(fun);if mm>=8error('为了保证NewtonCotes积分的稳定性，最多只能有9个等距节点！')elseif nn~=2error('fun构成应为：第一列为x，第二列为y，并且个数为小于10的等距节点！')endxk=fun(1,:);fk=fun(2,:);a=min(xk);b=max(xk);n=mm-1;elseif nargin==4xk=linspace(a,b,n+1);if isa(fun,'function_handle')fx=fun(xk);elseerror('fun积分函数的句柄，且必须能够接受矢量输入！') endelseerror('输入参数错误，请参考函数帮助！')endCk=cotescoeff(n);Ak=(b-a)*Ck;y=Ak*fx';（2）function Ck=cotescoeff(n)for i=1:n+1k=i-1;Ck(i)=(-1)^(n-k)/factorial(k)/factorial(n-k)/n*quadl(@(t)int fun(t,n,k),0,n);end（3）function f=intfun(t,n,k)f=1;for i=[0:k-1,k+1:n]f=f.*(t-i);end代码解释：function [y,Ck,Ak]=NewtonCotes(fun,a,b,n)% y=NewtonCotes(fun,a,b,n)% 牛顿-科特斯数值积分公式% 参数说明：% fun，积分表达式，这里有两种选择%(1)积分函数句柄，必须能够接受矢量输入，比如fun=@(x)sin(x).*cos(x)% (2)x,y坐标的离散点，第一列为x，第二列为y，必须等距，且节点的个数小于9，比如：fun=[1:8;sin(1:8)]'% 如果fun的表采用第二种方式，那么只需要输入第一个参数即可，否则还要输入a,b,n三个参数% a，积分下限% b，积分上限% n，牛顿-科特斯数公式的阶数，必须满足1<n<7，因为n>=8时不能保证公式的稳定性% (1)n=1，即梯形公式% (2)n=2，即辛普森公式% (3)n=4，即科特斯公式% y，数值积分结果% Ck，科特斯系数% Ak，求积系数%% Example% fun1=@(x)sin(x);%必须可以接受矢量输入% fun2=[0:0.1:0.5;sin(0:0.1:0.5)];%最多8个点，必须等距% y1=NewtonCotes(fun1,0,0.5,6)% y2==NewtonCotes(fun2)if nargin==1[mm,nn]=size(fun);if mm>=8error('为了保证NewtonCotes积分的稳定性，最多只能有9个等距节点！')elseif nn~=2error('fun构成应为：第一列为x，第二列为y，并且个数为小于10的等距节点！')endxk=fun(1,:);fk=fun(2,:);a=min(xk);b=max(xk);n=mm-1;elseif nargin==4% 计算积分节点xk和节点函数值fxxk=linspace(a,b,n+1);if isa(fun,'function_handle')fx=fun(xk);elseerror('fun积分函数的句柄，且必须能够接受矢量输入！')endelseerror('输入参数错误，请参考函数帮助！')end% 计算科特斯系数Ck=cotescoeff(n);% 计算求积系数Ak=(b-a)*Ck;% 求和算积分y=Ak*fx';function Ck=cotescoeff(n)% 由于科特斯系数最多7阶，为了方便我们可以直接使用，省得每次都计算% A1=[1,1]/2% A2=[1,4,1]/6% A3=[1,3,3,1]/8% A4=[7,32,12,32,1]/90% A5=[19,75,50,50,75,19]/288% A6=[41,216,27,272,27,216,41]/840% A7=[751,3577,1323,2989,2989,1323,3577,751]/17280% 当时为了体现公式，我们使用程序计算n阶科特斯系数for i=1:n+1k=i-1;Ck(i)=(-1)^(n-k)/factorial(k)/factorial(n-k)/n*quadl(@(t)intfun( t,n,k),0,n);endfunction f=intfun(t,n,k)% 科特斯系数中的积分表达式f=1;for i=[0:k-1,k+1:n]f=f.*(t-i);end输出结果：fun=@(x)exp(x);a=-1;b=1;n=4;NewtonCotes(fun,a,b,n)学海无涯ans =2.3505二、三点数值微分。