数值分析 高斯—勒让德积分公式

1、高斯勒让德积分公式摘要：高斯勒让德积分公式可以用较少节点数得到高精度的计算结果,是现在现实生活中经常运用到的数值积分法。然而，当积分区间较大时，积分精度并不理想。The adva ntage of Gauss-Legendre integral formula is tend to get high-precision calculational result by using fewer Gauss-points, real life is now often applied numerical integration method. But the precision is not goo

2、d when the length of integral interval is longer.关键字： 积分计算，积分公式，高斯勒让德积分公式，MATLABKeyword:Integral Calculation , Integral formula ,Gauss-Legendre integral formula, Matlab引言：众所周知，微积分的两大部分是微分与积分。微分实际上是求一函数的导数，而积分是已知一函数的导数，求这一函数。所以，微分与积分互为逆运算。实际上，积分还可以分为两部分。第一种，是单纯的积分，也就是已知导数求原函数，称为不定积分。相对而言，另一种就是定积分了，之所

3、以称其为定积分，是因为它积分后得出的值是确定的，是一个数，而不是一个函数。计算定积分的方法很多，而高斯勒让德公式就是其中之一。高斯积分法是精度最高的插值型数值积分，具有2n+1阶精度，并且高斯积分总是稳定。而高斯求积系数，可以由Lagrange多项式插值系数进行积分得到。高斯勒让德求积公式是构造高精度差值积分的最好方法之一。他是通过让节点和积分系数待定让函数f(x)以此取i=0,1,2.n次多项式使其尽可能多的能够精确成立来求出积分节点和积分系数。高斯积分的代数精度是2n-1，而且是最高的。通常运用的是（-1，1）的积分节点和积分系数，其他积分域是通过变换x=(b-a)t/2 

4、60; +(a+b)/2 变换到-1到1之间积分。1. 现有的方法和理论1.1高斯 勒让德求积公式在高斯求积公式()中，若取权函数，区间为，则得公式 我们知道勒让德多项式是区间上的正交多项式，因此，勒让德多项式的零点就是求积公式(上式)的高斯点形如(上式)的高斯公式特别地称为高斯勒让德求积公式若取的零点做节点构造求积公式令它对准确成立，即可定出这样构造出的一点高斯勒让德求积公式是中矩形公式再取的两个零点构造求积公式令它对都准确成立，有由此解出，从而得到两点高斯勒让德求积公式三点高斯勒让德求积公式的形式是如表列出高斯勒让德求积公式的节点和系数00.00000002.000000010.57735

5、031.000000020.77459670.00000000.55555560.888888930.86113630.33998100.34785480.652145240.90617980.53846930.00000000.23692690.47862870.5688889公式()的余项由(4.5.8)得，这里是最高项系数为的勒让德多项式，由()及(3.2.7)得 当时，有它比辛普森公式余项还小，且比辛普森公式少算一个函数值当积分区间不是，而是一般的区间时，只要做变换可将化为，这时 对等式右端的积分即可使用高斯勒让德求积公式1.2复化Gauss-Legendre求积公式将被积区间m等分,

6、 记, 作变换在每个小区间上应用Gauss-Legendre公式, 累加即得复化Gauss-Legendre求积公式不妨设则有:Gauss点个数时,Gauss点个数时,总结复化Gauss-Legendre求积过程如下:1. 分割区间, 记录区间端点值；2. 通过查表或求解非线性方程组, 在所有小区间上, 将Gauss系数和Gauss点的值代入变量替换后的公式；3. 将所有区间的结果累加, 即得到整个区间上的积分近似值.针对Gauss点个数和的复化Gauss-Legendre求积公式编写的一个简单的MATLAB函数 compgauss() 如下: function = compgauss(a,

7、b, n)% Composite Gauss Integration% Equation Type: n=2, n=3% Coded by Nan.Xiao 2010-05-25% Step.1 Divide Interval% Step.2 Calculate% Step.3 Sum Resultsformat longf = (x) exp(x).*sin(x);h=(b-a)/n;xk=zeros(n+1,1);xk(1,1)=a;xk(n+1,1)=b;fk1=zeros(n,1);fk2=zeros(n,1);for i=1:n-1 xk(i+1,1)=a+h*i;endfor j=

8、1:n fk1(j)=f(xk(j)+xk(j+1)/2+(h/2)*(-1/sqrt(3)+. f(xk(j)+xk(j+1)/2+(h/2)*(1/sqrt(3);endfor r=1:n fk2(r)=(5/9)*f(xk(r)+xk(r+1)/2+(h/2)*(-sqrt(15)/5)+. (8/9)*f(xk(r)+xk(r+1)/2+(h/2)*(0)+. (5/9)*f(xk(r)+xk(r+1)/2+(h/2)*(sqrt(15)/5);endmysum1=h*sum(fk1)/2;mysum2=h*sum(fk2)/2;disp('Result of 2 Nodes:

9、')disp(mysum1);disp('Result of 3 Nodes:')disp(mysum2);end1.3龙贝格，三点，五点以及变步长高斯勒让德求积法以下是关于龙贝格，三点，五点以及变步长高斯勒让德之间精度的相互比较#include <iostream.h> #include <math.h> #include <iomanip.h> #define Precision1 0.000000000001 # define e  

10、60;      2.71828183 #define  MAXRepeat 10   double function (double x)   double s; s=1/x; return s; double Romberg(double a,double b,double f(double x)   int m,n,k;  

11、;double yMAXRepeat,h,ep,p,xk,s,q; h=b-a;  y0=h*(f(a)+f(b)/2.0;/计算T1(h)=1/2(b-a)(f(a)+f(b);  m=1; n=1; ep=Precision1+1;  while(ep>=Precision1)&&(m<MAXRepeat)  p=0.0; for(k=0;k<n;k+)  xk=a+(k+0.5)*h;  p=p+f(xk);       &#

12、160;            p=(y0+h*p)/2.0;  /Tm(h/2),变步长梯形求积公式  s=1.0; for(k=1;k<=m;k+) s=4.0*s;/ pow(4,m) q=(s*p-yk-1)/(s-1.0); yk-1=p; p=q;  ep=fabs(q-ym-1);  m=m+1;         &#

13、160;   ym-1=q; n=n+n;   /  2 4 8 16    h=h/2.0;/二倍分割区间 return q; double ThreePointGaussLegendre(double a,double b,double f(double x)  double x,w;  static double X3=-sqrt(15)/5.0,0,

14、sqrt(15)/5.0;   static double L3=5/9.0,8/9.0,5/9.0;     w=0.0;   for(int i=0;i<3;i+)                         x=(b-a)*Xi+(b+a)/2

15、.0;                           w=w+f(x)*Li;                return w; double FivePointGaussL

16、egendre(double a,double b,double f(double x)   double x,w;     static double X5=-0.9061798459,-0.5384693101,0,0.5384693101,0.9061798459;     static double L5=0.2369268851,0.4786286705,0.5688888889,0.478628

17、6705,0.2369268851;   w=0.0;     for(int i=0;i<5;i+)                         x=(b-a)*Xi+(b+a)/2.0;        &#

18、160;                     w=w+f(x)*Li;/每一次小区间利用勒让德公式计算的结果             return w; double FivePointPrecisionGaussLegendre(double a,dou

19、ble b,double f(double x)   int m,i,j;     double s,p,ep,h,aa,bb,w,x,g;     static double X5=-0.9061798459,-0.5384693101,0,0.5384693101,0.9061798459;   m=1;   h=b-a;     s=fabs

20、(0.001*h); p=1.0e+35;    ep=Precision1+1;    while(ep>=Precision1)&&(fabs(h)>s)           g=0.0;         for(i=0;i<m;i+)        

21、             aa=a+i*h;            bb=aa+h;             w=0.0;          for(

22、j=0;j<=4;j+)                           x=(bb-aa)*Xj+(bb+aa)/2.0;                w=w+f(x)*Lj;

23、                       g=g+w;/各个区间计算结果之和相加                     g=g*h/2.0;     &#

24、160;       ep=fabs(g-p)/(1.0+fabs(g);/计算精度             p=g;             m=m+1;            

25、h=(b-a)/m;/分割区间             return g; main()     double a,b,s;     cout<<"请输入积分下限:"     cin>>a;     cout<<"请输入积分上限:"

26、    cin>>b;     cout<<"的真值为:"<<endl;     cout<<"1.098612289"<<endl;     /*龙贝格求积*/     s=Romberg( a, b, function);    

27、60;cout<<"龙贝格求积公式:"<<endl;     cout<<setiosflags(ios:fixed)<<setprecision(14)<<s<<endl;     /*三点求积公式*/     s=ThreePointGaussLegendre( a, b, function);     cout

28、<<"三点求积公式:"<<endl;     cout<<setiosflags(ios:fixed)<<setprecision(14)<<s<<endl;     /*五点求积公式*/     s=FivePointGaussLegendre( a, b, function);     cout<<&

29、quot;五点求积公式"<<endl;     cout<<setiosflags(ios:fixed)<<setprecision(14)<<s<<endl;      s=FivePointPrecisionGaussLegendre(a, b,function);     cout<<"控制精度五点求积公式"<<endl;  &

30、#160;  cout<<setiosflags(ios:fixed)<<setprecision(14)<<s<<endl;     return 0; 2. 高斯勒让德求积的程序2.1三点高斯勒让德公式的代码function gl=f(str,a,b)x=zeros(3,1);y=zeros(3,1);x(1)=-sqrt(15)/5;x(2)=0;x(3)=sqrt(15)/5;for i=1:3 t=(b-a)/2*x(i)+(a+b)/2; y(i)=eval(str);%exp(t)*sin(t);%此处为求积的函数,t为自变量endgl=5/9*y(1)+8/9*y(2)+5/9*y(3);上面的代码保存为f.m文件，调用的时候如下f('t*2',-1,1)f('exp(t)*sin(t)',1,3)其中第一个参数为求积