数学实验“几种常见的求积分近似解的方法”实验报告(内含matlab程序).doc

数学实验“几种常见的求积分近似解的方法”实验报告(内含matlab程序)

———————————————————————————————— 作者：
———————————————————————————————— 日期：

2

个人收集整理 勿做商业用途
西京学院数学软件实验任务书
课程名称
数学软件实验
班级
数0901
学号
0912020107
姓名
李亚强
实验课题
Romberg积分法,Gauss型积分法,高斯—勒让德积分法，高斯—切比雪夫积分法，高斯-拉盖尔积分法,高斯—埃尔米特积分法
实验目的
熟悉 Romberg积分法，Gauss型积分法，高斯—勒让德积分法，高斯-切比雪夫积分法，高斯—拉盖尔积分法，高斯-埃尔米特积分法
实验要求
运用Matlab/C/C++/Java/Maple/Mathematica等其中一种语言完成
实验内容
Romberg积分法，Gauss型积分法,高斯—勒让德积分法，高斯—切比雪夫积分法，高斯—拉盖尔积分法，高斯—埃尔米特积分法
成绩

教师

实验二十一实验报告
一、 实验名称：Romberg积分法，Gauss型积分法，高斯-勒让德积分法，高斯-切比雪夫积分法，高斯-拉盖尔积分法，高斯-埃尔米特积分法.
二、 实验目的:进一步熟悉Romberg积分法，Gauss型积分法,高斯-勒让德积分法,高斯—切比雪夫积分法，高斯-拉盖尔积分法，高斯-埃尔米特积分法。
三、 实验要求：运用Matlab/C/C++/Java/Maple/Mathematica等其中一种语言完成程序设计。
四、 实验原理：
1． Romberg积分法：
龙贝格积分法是用里查森外推算法来加快复合梯形求积公式的收敛速度，它的算法如下，其中是通过一系列逼近原定积分的龙贝格分值。
计算
对，计算下列各步：

对和,计算
随着计算的步骤的增加，越来越逼近积分。
2． Gauss型积分法：
高斯积分公式的思想是用个不等距的节点对被积函数进行插值，然后对插值后的函数进行积分，其积分公式为:

如果积分区间不是，则需转换到此区间：

其中系数、节点与的关系如下表所示:
3． 高斯-切比雪夫积分法:
第一类切比雪夫积分形式为：

其中，
4． 高斯—拉盖尔积分法：
高斯—拉盖尔公式有两种形式：

下面编制的程序是针对第一种形式的高斯—拉盖尔公式,即

因此程序的第一个输入参数—-被积函数，是上式中的。
5． 高斯—埃尔米特积分法：
高斯-埃尔米特公式有以下两种形式：

下面编制的程序是针对第一种形式的高斯—埃尔米特公式，即

因此程序的第一个输入参数——被积函数，是上式中的。
五、 实验内容：
%Romberg积分法
function [q，step］=Roberg(f,a，b,eps）
if(nargin==3）
eps=1.0e—4;
end
M=1;
tol=10；
k=0;
T=zeros(1,1);
h=b-a；
T(1,1)=(h/2）*(subs（sym(f),findsym（sym（f)),a）+subs(sym（f），findsym（sym(f))，b)）;
while tol〉eps
k=k+1;
h=h/2;
Q=0；
for i=1:M
x=a+h*(2＊i-1）；
Q=Q+subs(sym（f），findsym(sym(f)），x）;
end
T（k+1,1）=T（k,1）/2+h＊Q;
M=2＊M;
for j=1:k
T（k+1,j+1）=T（k+1,j)+(T（k+1j）—T（k,j）)/(4^j-1);
end
tol=abs（T(k+1,j+1）—T（k，j）);
end
q=T(k+1,j+1）；
step=k;
％Gauss型积分法
function q=IntGauss（f,a，b，n,AK,XK）
if(n〈5＆＆nargin==4)
AK=0；
XK=0;
else
XK1=(（b-a)/2）＊XK+（(a+b）/2）;
q=（(b-a)/2）＊sum（AK。*subs(sym（f），findsym（f),XK1))；
end
ta=（b—a)/2；
tb=(ab）/2；
switch n
case 1
q=2*ta*subs（sym（f)，findsym（sym(f)）,tb);
case 2
q=ta＊（subs(sym(f),findsym(sym（f）),ta＊0.5773503+tb）+subs（sym(f)，findsym（sym(f））,—ta＊0.5773503+tb)）;
case 3
q=ta＊(0.55555556*subs(sym（f），findsym(sym（f)）,ta＊0。7745967+tb）+0。55555556*subs（symf),findsym（sym(f))，-ta*0。7745967+tb)+0.88888889*subs（sym(f),findsym(sym(f)）,tb)）；
case 4
q=ta*（0。3478548*subs(sym(f）,findsym（sym(f)），ta＊0。8611363+tb)+0。3478548＊subs（sym(f)，findsym（sym(f)),—ta＊0。8611363+tb)+0.6521452*subs（sym（f)，findsym(sym(f)),ta*0。3398810+tb)+0.6521452＊subs(sym(f)，findsym(sym（f）),—ta＊0.3398810+tb)）；
case 5
q=ta＊(0.2369269*subs（sym(f），findsym(sym（f）），ta＊0。9061793+tb）+0.2369269＊subs(sym（f）,findsym（sym(f）），—ta*0.9061793+tb）+0.4786287*subs（sym(f),findsym（sym(f）),ta*0。5384693+tb)+0。4786287*subs(sym（f），findsym(sym(f)）,-ta＊0.5384693+tb)+0.5688889*subs（sym(f)，findsym(sym(f))，tb)）;
end
％高斯-勒让德积分法
function [ql，Ak,xk]=guasslegendre（fun,a，b,n，tol)
if nargin==1
a=—1;b=1；n=7;tol=1e-8;
elseif nargin==3
n=7；tol=1e—8；
elseif nargin==4
tol=1e-8;
elseif nargin==2｜nargin〉5
error(’The Number of Input Arguments Is Wrong!'）；
end
syms x
p=sym2poly（diff（（x^2-1）^（n+1），n+1））/(2^n*factorial(n））;
tk=roots(p)；
Ak=zeros(n+1,1);
for i=1：n+1
xkt=tk;
xkt(i）=［］;
pn=poly（xkt);
fp=@（x）polyval（pn,x)/polyval(pn,tk(i））；
Ak(i）=quadl(fp,—1，1，tol);
end
xk=（b-a）/2*tk+（b+a）/2；
fun=fcnchk（fun，'vectorize’);
fx=fun（xk)b—a）/2；
ql=sum（Ak。*fx）；
％高斯-切比雪夫积分法
function q = IntQBXF1（func,n）
format long；
pi = 3.1415926535；
q = 0；
A = zeros(n,1)；
x = zeros（n,1）；
for i=1:n
A(i） = pi/n;
x(i) = cos（pi＊(2＊i—1)/2/n）；
y(i) = subs(sym（func）, findsym（sym（func）)， x(i）);
q = q + A（i)＊y(i）;
end
%高斯—拉盖尔积分法
function q=IntGaussLager(f,n，AK，XK)
if（n<6＆＆nargin==2）
AK=0;
XK=0；
else
q=sum（AK.＊subs(sym(f）,findsym(f）,XK）);
end
switch n
case 2
q=0.853553＊subs(sym(f）,findsym（sym（f）)，—0。585786)+0.146447＊subs（sym(f)，findsym（sym（f)）,3.414214)；
case 3
q=0.711093*subs（sym(f）,findsym（sym（f）)，0.415575)+0。278518*subs(sym(f),findsym(sym(f)),2.294280)+0.0103893*subs（sym（f），findsym（sym(f））,6。289945）;
case 4
q=0.603154＊subs(sym(f)，findsym（sym(f）)，0。322548）+0.357419＊subs（sym（f)，findsym(sym(f)），1。745761)+0.0388879＊subs(sym（f），findsym（sym（f）),4。536620)+0000539295＊subs（sym（f)，findsym（sym（f）），9。395071)；
case 5
q=0。521756*subs（sym（f）,findsym（sym(f)），0。263560）+0。398667*subs(sym（f),findsym(sym(f）)，1。413403)+0。0759424*subs（sym（f），findsym(sym（f）），3。596426）+0。00361176*subs（sym（f),findsym(sym（f））,7.085810)+0.0000233700*subs（sym(f)，findsym(symf)）,12。640801);
end
%高斯—埃尔米特积分法
function q=IntGaussHermite(f，n，AK,XK）
if（n〈6&&nargin==2)
AK=0;
XK=0;
else
q=sum（AK.*subs(sym（f)，findsym（f),XK));
end
switch n
case 2
q=0。886227＊（subs（sym（f）,findsym（sym(f）），—0.707107)+subs(sym（f），findsym(sym（f）)，0。707107))；
case 3
q=1.181636＊subs（sym(f）,findsym(sym（f)）)+0。295409*subs（sym（f),findsym(sym(f）)，1.224745)+subs(sym（f）,findsym(sym（f）),—1.224745);
case 4
q=0。544444*(subs（sym(f)，findsym(sym(f))，0.524648)+subs(sym(f）,findsym（sym(f））,—0.524648））+0.100000＊（subs(sym(f),findsym(sym(f）),1.650680)+subs(sym(f),findsym（sym（f）),—1。650680）)；
case 5
q=0。945309＊subs(sym（f），findsym(sym(f))，0）+0。393619*(subs(sym(f）,findsym（sym(f)),0.958572）+subs（sym（f)，findsym（sym（f））,-0。958572））+0。199532*（subs(sym（f），findsym(sym（f）),2。020183）+subs（sym(f），findsym（sym(f)),-2。020183）);
end
- 8 -