- 1、本文档共5页,可阅读全部内容。
- 2、原创力文档(book118)网站文档一经付费(服务费),不意味着购买了该文档的版权,仅供个人/单位学习、研究之用,不得用于商业用途,未经授权,严禁复制、发行、汇编、翻译或者网络传播等,侵权必究。
- 3、本站所有内容均由合作方或网友上传,本站不对文档的完整性、权威性及其观点立场正确性做任何保证或承诺!文档内容仅供研究参考,付费前请自行鉴别。如您付费,意味着您自己接受本站规则且自行承担风险,本站不退款、不进行额外附加服务;查看《如何避免下载的几个坑》。如果您已付费下载过本站文档,您可以点击 这里二次下载。
- 4、如文档侵犯商业秘密、侵犯著作权、侵犯人身权等,请点击“版权申诉”(推荐),也可以打举报电话:400-050-0827(电话支持时间:9:00-18:30)。
查看更多
matlab实现复化梯形公式,复化simpson公式以及romberg积分
实验目的
熟悉并掌握数值积分的方法,重要训练复化梯形公式,复化simpson公式以及romberg积分。
问题描述
问题三数值积分椭圆周长的计算。考虑椭圆,为计算其周长,只要计算其第一象限的长度即可.
用参数方程可以表示为,
计算公式为
为计算方便,我们可以令,即计算下面的积分
(可以归结为上面的形式)
采用复化梯形公式,复化Simpson公式以及Romberg积分的方法计算积分
给出通用程序,该通用程序可以计算任何一个函数在任意一个区间在给定的精度下的数值积分。程序输出为计算出的数值积分值以及计算函数值的次数。
算法介绍
首先利用给出的各迭代公式,设计程序。在matlab对话框中输入要计算的函数,给出区间和精度。
复化梯形的迭代公式为:
abf(x)dx=h/2fa+2j=1n-1fxj+f(b);
复化simpson迭代公式为:
abf(x)dx=h/3fa+2j=1(n2)-1fx2j+4j=1(n2)fx2j-1+f(b);
Romberg迭代公式为:
Rk,j=Rk,j-1+Rk,j-1—Rk-1,j-14j-1-1。
程序
对于复化梯形公式和复化simpson公式,我们放在jifenn.m中。
(%标记后的程序可用来把b看为变量时的算法实现)
%复化梯形公式
function y=jifenn(f,n,a,b) (说明:f表示任一函数,n精度,a,b为区间)
fi=f(a)+f(b);
h=(b-a)/n;
d=1;
%function f=jifen(n,a,b,c)
%syms t
%y=sqrt(1+(c^2-1)*cos(t)^2);
%ya=subs(y,t,a);
%yb=subs(y,t,b);
%fi=ya+yb;
for i=1:n-1
x=a+i*h;
fi=fi+2*f(x);
d=d+1;
%yx=subs(y,t,x);
%fi=fi+2*yx;
end
f4=h/2*fi,d
%复化simposon公式
f1=0;
f2=0;
dd=1;
for i=1:n-1
dd=dd+1;
if rem(i,2)~=0;
x1=a+i*h;
f1=f1+f(x1);
else rem(i,2)==0;
x2=a+i*h;
f2=f2+f(x2) ;
end
end
f3=(h/3)*(f(a)+4*f1+2*f2+f(b)),dd
对于romberg积分,建立romberg.m文件。
function y=romberg(f,n,a,b) (说明:f表示任一函数,n精度,a,b为区间)
z=zeros(n,n);
h=b-a;
z(1,1)=(h/2)*(f(a)+f(b));
f1=0;
for i=2:n
for k=1:2^(i-2)
f1=f1+f(a+(k-0.5)*h);
end
z(i,1)=0.5*z(i-1,1)+0.5*h*f1;
h=h/2;
f1=0;
for j=2:i
z(i,j)=z(i,j-1)+(z(i,j-1)-z(i-1,j-1))/(4^(j-1)-1);
end
end
z,n
运行结果
对于复化梯形公式和复化simpson公式,我们运行下列语句并得到结果:
fun=inline(sqrt(1+(0.5^2-1).*cos(t).^2));
jifenn(fun,8,0,pi/2)
f4 =
1.2111
d =
8
f3 =
1.2111
dd =
8
1.2111*4
ans =
4.8444
1.2111*4
ans =
4.8444
(说明:在本题中将椭圆中的未知量a取为1,b取为0.5。f4为复化梯形公式得到的椭圆周长,f3为复化simpson公式得到的椭圆周长)。
对于romberg,运行下列语句并最终得到结果为:
fun=inline(sqrt(1+(0.5^2-1).*cos(t).^2));
romberg(fun,8,0,pi/0.5)
z =
3.1416 0 0 0 0 0 0 0
3.1416 3.1416 0 0 0 0 0 0
4.7124 5.2360 5.3756
文档评论(0)