MATLAB数值分析: Gauss求积

Gauss求积公式为：[大谦MATLAB，dqmatlab点com]

\[\int_{-1}^{1} {f(x)dx}-\sum_{k=1}^{n} {{A}_{k}f({x}_{k})}\]

公式中的xk被称为Gauss求积点，Ak被称为求积系数。对于在一般区间[a b]中的求积，如果用Gauss求积公式，则必须进行变量替换

x=0.5*(a+b)+0.5*(b-a)t

以使[a b]→[-1 1]，并有

\[\int_{a}^{b} {f(x)dx}=\frac{b-a}{2}\int_{-1}^{1} {f[\frac{1}{2}(a+b)+\frac{1}{2}(b-a)t]dt}\]

对于上式等号右边的部分，可以应用Gauss求积公式进行数值积分。

Gauss求积公式中的Gauss求积点xk及求积系数Ak如表7-1所示。

表7-1 Gauss求积点xk及求积系数Ak

n	xk	Ak	n	xk	Ak
1	0.0000000	2.0000000	4	±0.8611363 ±0.3398810	0.3478548 0.6521452
2	±0.5773403	1.0000000	5	±0.9061793 ±0.5384693 0.00000000	0.2369269 0.4786287 0.5688889
3	±0.7745967 ±0.0000000	0.555555556 0.8888888889	6	±0.9324695 ±0.66120939 ±0.23861919	0.1713245 0.3607616 0.4679139

根据Gauss求积点与求积系数的关系，编写4点Gauss求积公式的M文件GaussLegendre.m，如下所示。

code.matlab

function s=GaussLegendre(a,b)
  % 用4点Gauss-求积公式求数值积分，其中a和b分别为积分下限和上限
  x=[0.8611363116 -0.8611363116 0.3399810436 -0.3399810436];
  u=[0.3478548451 0.3478548451 0.6521451549 0.6521451549];
  t=0.0;
  for i=1:4
      y=x(i)*(b-a)*0.5+(a+b)*0.5;
      t=t+u(i)*ff(y);
  end
  s=t*(b-a)*0.5;
  return

【例24】应用4点Gauss求积公式计算定积分

\[I=\int_{0}^{1} {{x}^{2}{e}^{x}dx}\]

在MATLAB中的实现介绍如下。

首先，编写求积函数ff.m。

code.matlab

function y=ff(x)
  y=x*x*exp(x);

然后，可以利用编写的函数进行求积：

code.matlab

>> GaussLegendre(0,1)
ans =
    0.7183