龙贝格算法数值实验摘要.docVIP

数值计算实验报告 题目：Romberg算法实验报告 学 院 数学与信息科学学院 专 业 13信计 学 号 2013433020 姓 名 齐学东 指导教师 高少芹 2016年 1月21日 实验三 Romberg算法 学号：2013433020 姓名：齐学东 专业：13信计 时间：2016-1-21 实验目的 通过本实验可以使学生理解如何在计算机上使用数值方法计算定积分∫?f(x)dx的近似值。 用此算法解决实际问题。主要锻炼分析问题及对数值方法产生的误差在实际执行算法中的变化，体会事前误差和事后误差的不同之处。 实验方法 龙贝格方法的基本思路： 龙贝格方法是在积分区间逐次二分的过程中，通过对梯形之值进行加速处理，从而获得高精度的积分值。 三．实验内容 (1) 求定积分 的近似值， 建立M文件sy3.m程序代码为： function [ output_args ] = Untitled1( input_args ) %UNTITLED1 Summary of this function goes here Detailed explanation goes here a=0; b=1; h=b-a; k=1; T(1,1)=(b-a)*(f(a)+f(b))/2; err=1; while(erreps) new=0; for j=1:2^(k-1) x=a+(2*j-1)*h/2; new=new+f(x); end T(k+1,1)=(T(k,1)+h*new)/2; for m=1:k T(k+1,m+1)=(4^m*T(k+1,m)-T(k,m))/(4^m-1); end err=abs(T(k+1,m+1)-T(k,m)); k=k+1; h=h/2; end k T; Vpa(T,8) %精度为10^(-8)次方 m function z=f(x) if(x~=0) z=sin(x)/x; else z=1; end 在sy3.m文件中直接运行的结果为： (2)计算积分∫2π0(1-b2cos2t)1/2dt的近似值,其中b分别取0.1,0.3,0.5,0.7,0.9。 （a）Romberg算法： 建立M文件fun.m程序代码为： function [T,result]=fun(a,b,eps,q) digits(8); h=b-a; k=1; T(1,1)=(b-a)*(f2(a,q)+f2(b,q))/2; err=1; while(erreps) new=0; for j=1:2^(k-1) x=a+(2*j-1)*h/2; new=new+f2(x,q); end T(k+1,1)=(T(k,1)+h*new)/2; for m=1:k T(k+1,m+1)=(4^m*T(k+1,m)-T(k,m))/(4^m-1); end err=abs(T(k+1,m+1)-T(k,m)); k=k+1; h=h/2; end T=vpa(T); res=vpa(T(k,m+1)); （b）复合Simpson公式（n=16） Simpson程序为： function S=S(a,b,n,q) h=(b-a)/n; X0=f2(a,q)+f2(b,q); X1=0; X2=0; for i=1:n-1 x=a+i*h; if (rem(i,2)==0) X2=X2+f2(x,q); else X1=X1+f2(x,q); end end X=h*(X0+4*X1+2*X2)/3; S=vpa(X,8); 建立M函数文件sinpson.m程序代码为： function S=simpson(a,b,n,q) h=(b-a)/n; X0=f2(a,q)+f2(b,q); X1=0; X2=0; for i=1:n-1 x=a+i*h; if (rem(i,2)==0) X2=X2+f2(x,q); else X1=X1+f2(x,q)