高斯-勒让德求积公式及Matlab实现

初学Markdown编辑器，编排及内容错误请指正，谢谢。

1. 引言

       众所周知，微积分的两大部分是微分与积分。微分实际上是求一函数的导数，而积分是已知一函数的导数，求这一函数。所以，微分与积分互为逆运算。实际上，积分还可以分为两部分。第一种，是单纯的积分，也就是已知导数求原函数，称为不定积分。相对而言，另一种就是定积分了，之所以称其为定积分，是因为它积分后得出的值是确定的，是一个数，而不是一个函数。计算定积分的方法很多，而高斯—勒让德公式就是其中之一。
       高斯积分法是精度最高的插值型数值积分，具有$2n+1$阶精度，并且高斯积分总是稳定。而高斯求积系数，可以由Lagrange多项式插值系数进行积分得到。
       高斯-勒让德求积公式是构造高精度差值积分的最好方法之一。他是通过让节点和积分系数待定让函数$f(x)$以此取$i=0,1,2....n$次多项式使其尽可能多的能够精确成立来求出积分节点和积分系数。高斯积分的代数精度是$2n-1$，而且是最高的。通常运用的是$(-1,1)$的积分节点和积分系数，其他积分域是通过一定的变换变换到-1到1之间积分。

2. 高斯-勒让德求积公式

       在区间$[-1,1]$上，高斯-勒让德求积公式为

$$\int^{1}_{-1}{f(x)}dx\approx\sum^{x}_{k=0}A_{k}f(x_k)$$ 我们知道勒让德多项式 $P_{n+1}(x)$ 的零点就是求积公式的高斯点，形如上式的高斯公式特别的称为高斯-勒让德公式。
       若取 $P_1(x)=x$ 的零点 $x_0=0$ 做节点构造求积公式
$$\int^{1}_{-1}{f(x)}dx\approx A_0f(0)$$ 令它对 $f(x)=1$准确成立，即可定出 $A_0=2$ 。这样构造出的一点高斯-勒让德求积公式是中矩形公式，再取 $P_2(x)=\frac{1}{2}(3x^2-1)$的两个零点 $\pm\frac{1}{\sqrt{3}}$构造求积公式
$$\int^{1}_{-1}{f(x)}dx\approx A_0f(-\frac{1}{\sqrt{3}})+A_1f(\frac{1}{\sqrt{3}})$$ 令它对 $f(x)=1,x$ 都准确成立，有
$$\left\{ \begin{aligned} &A_0+A_1=2, \\ &A_0(-\frac{1}{\sqrt{3}})+A_1(\frac{1}{\sqrt{3}})=0 \end{aligned} \right.$$
由此解出 $A_0=A_1=1$，从而得到两点高斯-勒让德求积公式
$$\int^{1}_{-1}{f(x)}dx\approx f(-\frac{1}{\sqrt{3}})+f(\frac{1}{\sqrt{3}})$$
三点高斯-勒让德公式的形式是
$$\int^{1}_{-1}{f(x)}dx\approx \frac{5}{9}f(-\frac{\sqrt{15}}{5})+\frac{8}{9}f(0)+ \frac{5}{9}f(\frac{\sqrt{15}}{5})$$ 下表给出常用的高斯-勒让德求积公式的节点和系数

节点数	$x_k$	$A_k$
1	0.0000000	2.0000000
2	$\pm$0.5773503	1.0000000
3	$\pm$0.7745967 0.0000000	0.5555556 0.8888889
4	$\pm$0.8611363 $\pm$0.3399810	0.3478548 0.6521452

当积分区间不是$[-1,1]$，而是一般的区间$[a,b]$，只要做变换

$$x=\frac{b-a}{2}t+\frac{a+b}{2}$$ 可将 $[a,b]$化为 $[-1,1]$，这时
$$\int^{b}_{a}{f(x)}dx=\frac{b-a}{2}\int^{1}_{-1}{f(\frac{b-a}{2}t+\frac{a+b}{2})}dt$$

3. Matlab数值求解高斯-勒让德积分

例1：用两点高斯-勒让德积分公式求 $\int_{0}^1{x^3}dx$ 。
       我们可以很容易知道结果为 $0.25$ 接下来用matlab数值求解：

GaussP=[-0.5773503 0.5773503];                              %高斯点
GaussA=[1 1];                                               %高斯系数
h = 0.1;                                                    %剖分步长
x = 0:h:1;                                                  %区间[0,1]
for i=1:length(x)-1
    points = h/2*GaussP + (x(i+1)+x(i))/2;                  %区间变换
    f(i) = 0;
    for k=1:2
        f(i) = f(i) + h/2*points(k)^2*GaussA(k);
    end
end
result = sum(f)

可以得到结果为：0.250000000133413。
例2：用三点高斯-勒让德积分公式求 $\int_{0}^1{x^2(1-x)^2}dx$。
       我们可以很容易知道结果为 $\frac{1}{30}\approx0.033333333333333$ 接下来用matlab数值求解：

GaussP=[-0.7745967 0 0.7745967];                            %高斯点
GaussA=[0.5555556 0.8888889 0.5555556];                     %高斯系数
h = 0.1;                                                    %剖分步长
x = 0:h:1;                                                  %区间[0,1]
for i=1:length(x)-1
    points = h/2*GaussP + (x(i+1)+x(i))/2;                  %区间变换
    f(i) = 0;
    for k=1:3
        f(i) = f(i) + h/2*(points(k)^2*(1-points(k))^2)*GaussA(k);
    end
end
result = sum(f)

结果为：0.033333334999780