文档库

最新最全的文档下载
当前位置:文档库 > matlab实现复化NewtonCotes公式求积分的程序应用和代码

matlab实现复化NewtonCotes公式求积分的程序应用和代码

执行函数为1、使用方法:

Step1:在MATLAB命令窗口输入被积函数

2

1

2

t

t e dt 。

输入应为:。

Step2:执行函数。输入形式为mymulNewtonCotes(ft,a,b,m,n);

其中ft—被积函数,此体重,已经在第一步赋值;

a—积分下限,本题中为0;

b—积分上限,本题中为1;

m—将区间[a,b]等分的子区间数量,本题可选为10;

n—采用的Newton-Cotes公式的阶数,必须满足n<8,否则积分没法保证稳定性。

当n=1时,即为复化梯形公式;n=2时,即为复化复化辛普森公式。

所以,分别输入mymulNewtonCotes(ft,0,1,10,1)和mymulNewtonCotes(ft,0,1,10,2)就可以得到两种方法的积分计算结果。

2、计算结果

而根据积分运算,可得:

说明复化梯形和复化辛普森公式计算出的结果基本一致,与实际结果相符。

3、程序代码

functionyy=mymulNewtonCotes(ft,a,b,m,n)

%复化Newton-Cotes数值积分公式,即在每个子区间上使用Newton-Cotes公式,然后求和,

%参考的输入形式为mymulNewtonCotes(ft,0,1,10,2)

%参数说明:

%ft——被积函数,此题中ft=@(t)t.*exp(t^2/2)

%a——积分下限

%b——积分上限

%m——将区间[a,b]等分的子区间数量

%n——采用的Newton-Cotes公式的阶数,必须满足n<8,否则积分没法保证稳定性

%(1)n=1时为复化梯形公式

%(2)n=2时为复化辛普森公式

xx=linspace(a,b,m+1);

forl=1:m

s(l)=myNewtonCotes(ft,xx(l),xx(l+1),n);

end

yy=sum(s);

function[y,Ck,Ak]=myNewtonCotes(ft,a,b,n)

%牛顿-科特斯数值积分公式

%Ck——科特斯系数

%Ak——求积系数

%y——牛顿-科特斯数值积分结果

xk=linspace(a,b,n+1);

forj=1:n+1

ff(j)=ft(xk(j));

end

%计算科特斯系数

fori=1:n+1

k=i-1;

Ck(i)=(-1)^(n-k)/factorial(k)/factorial(n-k)/n*quadl(@(t)intfun(t,n,k),0,n); end

%计算求积系数

Ak=(b-a)*Ck;

%求和算积分

y=Ak*ff';

functionf=intfun(t,n,k)

%科特斯系数中的积分表达式

f=1;

fori=[0:k-1,k+1:n]

f=f.*(t-i);

end