legendre多项式与gauss-legendre积分

写作目的

最近在学谱方法解偏微分方程，顺便学一下legendre多项式的一些推导，参考了北京大学蓝以中老师的《高等代数学习指南》，发博客以记之。

Legendre多项式定义

P 0 ( x ) = 1 P k ( x ) = 1 2 k k ! d k [ ( x 2 − 1 ) k ] d x k , k = { 1 , 2 , ⋯ } (1) {P_0}\left( x \right) = 1\newline{P_k}\left( x \right) = \frac{1}{{{2^k}k!}}\frac{{{d^k}\left[ {{{\left( {{x^2} - 1} \right)}^k}} \right]}}{{d{x^k}}},k = \left\{ {1,2, \cdots } \right\}\tag{1} P0​(x)=1Pk​(x)=2kk!1​dxkdk[(x2−1)k]​,k={1,2,⋯}(1)

五个引理

1. 莱布尼茨求导法则:
   $${\left(uv\right)}^{\left(k\right)}=\sum _{i=0}^k{C}_k^i{u}^{\left(i\right)}{v}^{\left(k-i\right)}\text{\hspace{1em}(2)}$$
2. 利用分部积分,有:$$\begin{gathered} \underset{a}{\overset{b}{\int }}{u}^{\left(m\right)}vdx=u^{\left(m-1\right)}v{|}_a^b-\underset{a}{\overset{b}{\int }}{u}^{\left(m-1\right)}{v}^{\left(1\right)}dx=\cdots \\ =\left[u^{\left(m-1\right)}v-u^{\left(m-2\right)}{v}^{\left(1\right)}+u^{\left(m-3\right)}{v}^{\left(2\right)}+\cdots + \\ {\left(-1\right)}^{m-1}u{v}^{\left(m-1\right)}\right]{|}_a^b \\ +{\left(-1\right)}^m\underset{a}{\overset{b}{\int }}u{v}^{\left(m\right)}dx\text{\hspace{1em}(3)} \end{gathered}$$
3. 引理三:(x^2-1)整除关系
   $$1\le l<k$$时$$\left(x^2-1\right)|\frac{d^l\left[{\left(x^2-1\right)}^k\right]}{d{x}^l}\text{\hspace{1em}(4)}$$
   这里蓝以中老师的书中用了归纳法证明，实际上不需要，很容易看出$${\left(x^2-1\right)}^k={\left(x-1\right)}^k{\left(x+1\right)}^k$$并利用引理2直接推出
   $$\frac{d^l\left[{\left(x^2-1\right)}^k\right]}{d{x}^l}=\sum _{i=0}^l{C}_l^i\frac{d^i\left[{\left(x+1\right)}^k\right]}{d{x}^i}\frac{d^{l-i}\left[{\left(x-1\right)}^k\right]}{d{x}^{l-i}}\text{\hspace{1em}(4.5)}$$
   直接看出求和项中每一项都可以被(x^2-1)整除即证。
   4.引理4：这个好像挺基础吧
   $$\underset{0}{\overset{\frac{\pi }{2}}{\int }}{\cos }^{2n+1}\left(x\right)dx=\frac{2n}{2n+1}\times \frac{2n-2}{2n-1}\times \cdots \times \frac{2}{3}\text{\hspace{1em}(5)}$$
   5.引理5：这个也显然:
   $$P_k\left(1\right)=1$$
   事实上,运用展开式(4.5)直接看出

勒让德多项式正交性证明:

$$1\le l<k$$时，可以得到:
$$\underset{-1}{\overset{1}{\int }}{P}_l\left(x\right)P_k\left(x\right)dx=\frac{1}{2^{k+l}k!l!}\underset{-1}{\overset{1}{\int }}\frac{d^k\left[{\left(x^2-1\right)}^k\right]}{d{x}^k}\frac{d^l\left[{\left(x^2-1\right)}^l\right]}{d{x}^l}dx$$
利用引理2可以得到:
$$\begin{gathered} 2^{k+l}k!l!\underset{-1}{\overset{1}{\int }}{P}_l\left(x\right)P_k\left(x\right)dx=\frac{d^{k-1}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{k-1}}\times \frac{d^l\left[{\left(x^2-1\right)}^l\right]}{d{x}^l}{|}_{-1}^1 \\ -\frac{d^{k-2}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{k-2}}\times \frac{d^{l+1}\left[{\left(x^2-1\right)}^l\right]}{d{x}^{l+1}}{|}_{-1}^1+\cdots \\ +{\left(-1\right)}^l\frac{d^{k-l-1}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{k-l-1}}\times \frac{d^{2l}\left[{\left(x^2-1\right)}^l\right]}{d{x}^{2l}}{|}_{-1}^1 \\ +{\left(-1\right)}^{l+1}\underset{-1}{\overset{1}{\int }}\frac{d^{k-l-1}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{k-l-1}}\times \frac{d^{2l+1}\left[{\left(x^2-1\right)}^l\right]}{d{x}^{2l+1}}dx \end{gathered}$$
利用引理3，除了积分号之外的结果都是0，积分号下第二项是对一个2l次的多项式求2l+1次导数，故为零，所以证明了legendre多项式的正交性。

模长证明

$$\underset{-1}{\overset{1}{\int }}{P}_k^2\left(x\right)dx={\left(\frac{1}{k!2^k}\right)}^2\underset{-1}{\overset{1}{\int }}\frac{d^k\left[{\left(x^2-1\right)}^k\right]}{d{x}^k}\times \frac{d^k\left[{\left(x^2-1\right)}^k\right]}{d{x}^k}dx$$
同样运用引理2得到
$$\begin{gathered} {\left(k!2^k\right)}^2\underset{-1}{\overset{1}{\int }}{P}_k^2\left(x\right)dx=\frac{d^{k-1}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{k-1}}\times \frac{d^k\left[{\left(x^2-1\right)}^k\right]}{d{x}^k}{|}_{-1}^1 \\ -\frac{d^{k-2}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{k-2}}\times \frac{d^{k+1}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{k+1}}{|}_{-1}^1+\cdots + \\ {\left(-1\right)}^k\frac{d\left[{\left(x^2-1\right)}^k\right]}{dx}\times \frac{d^{2k-2}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{2k-2}}{|}_{-1}^1 \\ +{\left(-1\right)}^{k+1}\left(x^2-1\right)\times \frac{d^{2k-1}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{2k-1}}{|}_{-1}^1 \\ +{\left(-1\right)}^{k+2}\underset{-1}{\overset{1}{\int }}{\left(x^2-1\right)}^2\frac{d^{2k}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{2k}}dx\text{\hspace{1em}(6)} \end{gathered}$$
同样(6)等式右端只有积分号下不为0,得到:
$$\begin{gathered} {\left(k!2^k\right)}^2\underset{-1}{\overset{1}{\int }}{P}_k^2\left(x\right)dx={\left(-1\right)}^{k+2}\underset{-1}{\overset{1}{\int }}{\left(x^2-1\right)}^2\frac{d^{2k}\left[{\left(x^2-1\right)}^k\right]}{d{x}^{2k}}dx= \\ {\left(-1\right)}^k\left(2k\right)!\underset{-1}{\overset{1}{\int }}{\left(x^2-1\right)}^{2}dx\text{\hspace{1em}(7)} \end{gathered}$$
(7)中最后式子利用三角积分换元与对称性得到:
$$\begin{gathered} {\left(-1\right)}^k\left(2k\right)!\underset{-1}{\overset{1}{\int }}{\left(x^2-1\right)}^{2}dx= \\ 2\times {\left(-1\right)}^{2k}\times \left(2k\right)!\underset{0}{\overset{\frac{\pi }{2}}{\int }}{\cos }^{2k+1}\left(x\right)dx\text{\hspace{1em}(8)} \end{gathered}$$
(8)利用引理4得到:
$$\begin{gathered} \underset{-1}{\overset{1}{\int }}{P}_k^2\left(x\right)dx=\frac{1}{{\left(k!2^k\right)}^2}\times \left(2\times \left(2k\right)!\right)\times \\ \frac{2k}{2k+1}\times \frac{2k-2}{2k-1}\times \cdots \times \frac{2}{3}=\frac{2}{2k+1}\text{\hspace{1em}(9)} \end{gathered}$$

递推公式证明

考虑欧氏空间$$R{\left[x\right]}_{n+2}$$
存在一组正交基$$P_0\left(x\right),P_1\left(x\right),\cdots ,P_{n+1}\left(x\right)$$
则n+1次多项式可被正交基唯一地线性表示$$x{P}_n\left(x\right)=c_{n+1}{P}_0\left(x\right)+c_n{P}_1\left(x\right)+\cdots +c_0{P}_{n+1}\left(x\right)\text{\hspace{1em}(10)}$$
则先证明:
$$c_3=c_4=\cdots =c_{n+1}=0$$
事实上，等式(10)两边乘以P_0(x)并积分得到
$$\underset{-1}{\overset{1}{\int }}x{P}_n\left(x\right)P_0\left(x\right)dx=c_{n+1}\underset{-1}{\overset{1}{\int }}{P}_0^2\left(x\right)dx$$
且:
$$\begin{gathered} \partial \left(x{P}_0\left(x\right)\right)=1<n\Rightarrow \underset{-1}{\overset{1}{\int }}x{P}_n\left(x\right)P_0\left(x\right)dx=0\Rightarrow c_{n+1}=0 \\ \partial \left(x{P}_1\left(x\right)\right)=2<n\Rightarrow \underset{-1}{\overset{1}{\int }}x{P}_n\left(x\right)P_1\left(x\right)dx=0\Rightarrow c_n=0 \\ \cdots \\ \partial \left(x{P}_{n-2}\left(x\right)\right)=n-1<n\Rightarrow \underset{-1}{\overset{1}{\int }}x{P}_n\left(x\right)P_{n-2}\left(x\right)dx=0\Rightarrow c_3=0 \end{gathered}$$
因此
$$x{P}_n\left(x\right)=c_0{P}_{n+1}\left(x\right)+c_1{P}_n\left(x\right)+c_2{P}_{n-1}\left(x\right)\text{\hspace{1em}(11)}$$
注意到根据勒让德多项式的定义(1):
$$x^{2n-1}\notin {\left(x^2-1\right)}^n\Rightarrow x^{n-1}\notin \frac{d^n\left[{\left(x^2-1\right)}^n\right]}{d{x}^n}$$
又由于$$x^n\notin P_{n+1}\left(x\right),P_{n-1}\left(x\right)\Rightarrow c_1=0$$
因此$$x{P}_n\left(x\right)=c_0{P}_{n+1}\left(x\right)+c_2{P}_{n-1}\left(x\right)\text{\hspace{1em}(12)}$$
比较等式(12)两边x^{n+1}次系数可得:
$$\begin{gathered} \frac{2n\left(2n-1\right)\cdots \left(n+1\right)}{2^{n}n!}=c_0\frac{\left(2n+2\right)\left(2n+1\right)\cdots \left(n+2\right)}{2^{n+1}\left(n+1\right)!}\Rightarrow \\ {c}_0=\frac{n+1}{2n+1} \end{gathered}$$
最后算c_2时直接利用引理5
$$1=\frac{n+1}{2n+1}+c_2\times 1\Rightarrow c_2=\frac{n}{2n+1}$$
故得到
$$x{P}_n\left(x\right)=\frac{n+1}{2n+1}{P}_{n+1}\left(x\right)+\frac{n}{2n+1}{P}_{n-1}\left(x\right)\text{\hspace{1em}(13)}$$
即为递推式的等价形式。

零点分布

n次勒让德多项式有n个实根且分布在(-1,1)开区间内:
证明使用罗德里格斯公式(1):首先注意到${\left(x^2-1\right)}^k,k\ge 1$一定有-1,1为根,使用rolle中值定理可以得到$\frac{d{\left(x^2-1\right)}^k}{dx}$存在$-1<{\xi }_{1,1}<1$为根,同时若$k>1$则可导出$\frac{d^2{\left(x^2-1\right)}^k}{d{x}^2}$有$-1<{\xi }_{2,1}<{\xi }_{2,2}<1$为根,以此类推:$\frac{d^{k-1}{\left(x^2-1\right)}^k}{d{x}^{k-1}},k\ge 1$有$-1<{\xi }_{k-1,1}<\cdots <{\xi }_{k-1,k-1}<1$为根,最后使用的rolle定理使得在这k+1个点形成的k个段内各有一个实根,即证.

等价的三对角矩阵的谱:

注意到递推公式(13):
等价表达式:
$$P_{n+1}\left(x\right)=\frac{2n+1}{n+1}x{P}_n\left(x\right)-\frac{n}{n+1}{P}_{n-1}\left(x\right)\text{\hspace{1em}(14)}$$
因此等价的行列式为,这里$D_n$为三对角矩阵:
$$D_{n+1}=\left[\begin{array}{cc}{D}_n& \sqrt{\frac{n}{n+1}}\\ \sqrt{\frac{n}{n+1}}& \frac{2n+1}{n+1}x\end{array}\right]$$
因此推得legendre多项式$P_{n+1}$的根是下列实对称三对角阵的谱:
T = d i a g ( k 4 k 2 − 1 , 1 ) + d i a g ( k 4 k 2 − 1 , − 1 ) , k ∈ { 1 , 2 , ⋯ n } T = diag\left( {\frac{k}{{\sqrt {4{k^2} - 1} }},1} \right) + diag\left( {\frac{k}{{\sqrt {4{k^2} - 1} }}, - 1} \right),k \in \left\{ {1,2, \cdots n} \right\} T=diag(4k2−1 ​k​,1)+diag(4k2−1 ​k​,−1),k∈{1,2,⋯n}
求积系数是对应特征向量模长平方的两倍

matlab求根:

n = 10;%11次勒让德多项式
A = diag((1:n)./sqrt(4*(1:n).^2-1),1)+diag((1:n)./sqrt(4*(1:n).^2-1),-1);
    [v,e] = eig(A);
    diag(e)%11个根

画出图像(python实现)

import numpy as np
import matplotlib.pyplot as plt
n = 10;
P = np.zeros((n+1,n+1));
P[n,0] = 1;
P[n-1,1] = 1;
x = np.linspace(-1, 1, 100)
J = np.zeros((n+1,n+1));
for i in range(n):
    J[i,i+1] = 1
#print(J)
for i in range(1,n):
    P[:,i+1] = ((2*i+1)/(i+1))*np.dot(J,P[:,i]) - i/(i+1)*P[:,i-1];
#print(P)
for i in range(n+1):
    y = np.polyval(P[:,i], x)
    plt.plot(x,y)
plt.show()

积分实例

尝试计算:
$$\underset{-1}{\overset{1}{\int }}\frac{1}{x^2+1}dx=\frac{\pi }{2}$$

import numpy as np  
import matplotlib.pyplot as plt  
def legendre_polynomial_integral_coefficient(n):  
    DOF = np.zeros((n,n))  
    for iter in range(n-1):  
        DOF[iter,iter+1] = (iter+1)/np.sqrt(4*(iter+1)**2-1)  
        DOF[iter+1,iter] = DOF[iter,iter+1]  
        [xr,v] = np.linalg.eig(DOF)  
        A = 2*v[0,:]**2  
#     print(np.linalg.norm(v@np.diag(xr)@v.T - DOF))  
    return xr,A  
def gauss_legendre_integral(f,xr,A,a=-1,b=1):#高斯勒让德求积  
       xrr = ((a+b) + (b-a)*xr)/2  
       return np.sum(f(xrr)*A*(b-a)/2)  
def grid_on_gauss(f,xr,A,intervals):   #复化高斯求积  
       s = 0.0;  
       for i in range(len(intervals)-1):  
         s = s + gauss_legendre_integral(f,xr,A,intervals[i],intervals[i+1])  
       return s  
def interpolant(x):  
    return 1/(1+x**2)  
if __name__=="__main__":  
    xr,A = legendre_polynomial_integral_coefficient(2) #五次多项式具有9次代数精度  
    print('求积节点和求积系数:',end = '')  
    print(xr,A)  
    err = np.zeros((1,8))  
    for i in range(2,10):  
        intervals = np.linspace(-1,1,i)  
        print('gauss-legendre求积复化节点:',end = '|')  
        print(intervals.round(4))  
        val = grid_on_gauss(interpolant,xr,A,intervals)  
        print('val = {:.16f}'.format(val))  
        err[0,i-2] = val - np.pi/2  
#     print(err[0,1],err[0,3],err[0,7])  
    print('2段和4段复化得到误差阶次{:.4f}'.format(np.log(np.abs(err[0,1]/err[0,3]))/np.log(2)))  
    print('4段和8段复化得到误差阶次{:.4f}'.format(np.log(np.abs(err[0,3]/err[0,7]))/np.log(2))) 

求积节点和求积系数:[ 0.57735027 -0.57735027] [1. 1.]
gauss-legendre求积复化节点:|[-1. 1.]
val = 1.4999999999999998
gauss-legendre求积复化节点:|[-1. 0. 1.]
val = 1.5737704918032787
gauss-legendre求积复化节点:|[-1. -0.3333 0.3333 1. ]
val = 1.5706921944035348
gauss-legendre求积复化节点:|[-1. -0.5 0. 0.5 1. ]
val = 1.5708049465155127
gauss-legendre求积复化节点:|[-1. -0.6 -0.2 0.2 0.6 1. ]
val = 1.5707969907681516
gauss-legendre求积复化节点:|[-1. -0.6667 -0.3333 0. 0.3333 0.6667 1. ]
val = 1.5707966371919944
gauss-legendre求积复化节点:|[-1. -0.7143 -0.4286 -0.1429 0.1429 0.4286 0.7143 1. ]
val = 1.5707964458929606
gauss-legendre求积复化节点:|[-1. -0.75 -0.5 -0.25 0. 0.25 0.5 0.75 1. ]
val = 1.5707963805291514
2段和4段复化得到误差阶次8.4306
4段和8段复化得到误差阶次7.3257