高斯积分

2020年6月11日2022年3月3日 numanal 0 个评论 Gauss-Legendre, Gauss-Lobatto, 积分法则, 高斯积分

积分法则

考虑实数区间 $[-1,1]$ 上的积分 $$ \int_{-1}^1 f(\xi)\, d\xi. \tag{1} $$ 一个数值积分法则 $ (\xi_i,\ w_i),i=1,\dots,n$, 即积分点和积分权重, 具有如下形式 $$ \int_{-1}^1 f(\xi) \,d\xi \approx \sum_{i=1}^{n} \omega_i f \left( \xi_i\right). \tag{2} $$ 我们希望选取的积分法则具有尽可能高的数值精度, 即数值积分 (2) 对尽可能高阶的多项式准确成立. 一般来说, $2n$ 个参数可以唯一确定一个 $2n-1$ 次多项式. 因此我们希望可以找到积分法则使得下面的等式成立 $$ \int_{-1}^{1} p(\xi) \,d\xi = \sum_{i=1}^{n} \omega_i p \left( \xi_i\right), \quad \forall\, p\in \mathbb{P}_{2n-1}. \tag{3} $$

一个简单的例子

问题: 寻找积分法则 $(\xi_i,\ w_i),i=1,2$ 使如下等式成立 \[ \int_{-1}^{1} p(\xi) \,d\xi = \omega_1 p(\xi_1)+\omega_2 p(\xi_2), \quad \forall\, p\in \mathbb{P}_{3}([-1,1]). \tag{4} \] 求解: 在等式 (4) 中分别另 $p(x)=1,\ \xi,\ \xi^2,\ \xi^3$, 得到 4 个方程: $$ \omega_1+\omega_2=\int_{-1}^{1} 1 \,d\xi =2,\quad \omega_1 \xi_1+\omega_2\xi_2=\int_{-1}^{1} \xi \,d\xi =0,\\ \omega_1 \xi_1^2+\omega_2\xi_2^2=\int_{-1}^{1} \xi^2 \,d\xi =\frac{3}{2}, \quad \omega_1 \xi_1^3+\omega_2\xi_2^3=\int_{-1}^{1} \xi^3 \,d\xi =0. $$ 求解上面的非线性方程组 (可以利用 Python 符号求解) 可以得到积分法则: \[ (\xi_1,\ \omega_1)=(-\frac{\sqrt{3}}{3},\ 1)=(-0.577350269189626,\ 1)\\ (\xi_2,\ \omega_2)=(\frac{\sqrt{3}}{3},\ 1)=(0.577350269189626,\ 1). \] 可以验证上面得到的积分法则对所有的 3-阶多项式精确成立.

# 在线运行如下代码求解非线性方程组 https://live.sympy.org/
from sympy import Symbol, nonlinsolve

w1, w2 = Symbol('w1'), Symbol('w2')
x1, x2 = Symbol('x1'), Symbol('x2')
eq1 = w1 + w2 - 2
eq2 = w1 * x1 + w2 * x2
eq3 = w1 * x1**2 + w2 * x2**2 - 2/3
eq4 = w1 * x1**3 + w2 * x2**3

nonlinsolve([eq1, eq2, eq3, eq4], [w1, w2, x1, x2])
# print(nonlinsolve([eq1, eq2, eq3, eq4], [w1, w2, x1, x2]))

Gauss 积分

更高阶的积分法则可以通过上面的方法计算得到, 但是最常用的 Gauss 积分可以很容易地通过正交多项式 (比如 Legendre 多项式) 来计算.

Legendre 多项式

在 $[-1,1]$ 上满足如下性质的多项式被称为 Legendre 多项式:

对任一 $n\in \mathbb{N}$, $P_n(x)$ 是 $n$-阶首一多项式 (最高阶项系数为 1);
对所有阶数小于 $n$ 的多项式 $P(x)$, $\int_{-1}^{1}P(x)P_n(x)\ dx=0$.

可以证明每一个 Legendre 多项式的根都是单根且对称分布在 $(-1, 1)$ 内. 重要地是, 如果我们以 $n$-阶 Legendre 多项式 $P_n(x)$ 的根作为积分节点, 再加上合适的权重, 就可以得到 Gauss-Legendre 积分法则 (具有 $2n-1$ 阶精度).

一个重要定理

定理: 令 $\xi_i,i=1,2,\dots,n$ 为 $P_n(\xi)$ 的 $n$ 个根, 选取权重 $$ \omega_i = \int_{-1}^{1} L_i(\xi)\ d\xi, $$ 这里 $L_i(\xi) = \sum_{j=1,j\neq i}^n \frac{\xi-\xi_j} {\xi_i-\xi_j}$ 为 Lagrange 插值基函数. 那么, 对所有小于 $2n-1$ 阶的多项式 $P(\xi)$, 我们有 $$ \int_{-1}^{1} P(\xi) \ d\xi = \sum_{i=1}^{n} \omega_i P(\xi_i). $$ 证明: 由多项式分解理论可知, 存在 $Q(\xi),R(\xi)\in \mathbb{P}{n-1}$ 使得 $P(\xi)$ 具有如下分解 $$ P(\xi) = Q(\xi)\, P_n(\xi)+R(\xi). $$ 因为 $P_n(\xi_i)=0$, 所以 $P(\xi_i)=R(\xi_i)$. 由 Legendre 多项式的正交性质, 有 $$ \int_{-1}^{1} P(\xi)\ d\xi = \int_{-1}^{1} Q(\xi)\, P_n(\xi)+R(\xi) \ dx = \int_{-1}^{1} R(\xi) \ d\xi. \tag{5} $$ 因为 $\mathbb{P}_{n-1}\ni R(\xi)=\sum_{i=1}^n R(\xi_i)L_i(\xi) =\sum_{i=1}^n P(\xi_i)L_i(\xi)$, 可知 $$ \int_{-1}^{1} R(\xi) \ d\xi = \int_{-1}^{1} \left( \sum_{i=1}^n P(\xi_i)L_i(\xi) \right) \ d\xi \\ =\sum_{i=1}^n \left(\int_{-1}^{1} L_i(\xi) \ d\xi \right) P(\xi_i) = \sum_{i=1}^n \omega_i P(\xi_i). $$ 代入 (5) 完成定理证明.

一般区间上的数值积分

通常的做法是将一般区间 $[a,b], a,b\in \mathbb{R}$ 上的积分线性变换到参考单元 $[-1,1]$ 上. 具体来说, 积分 (1) 有如下数值积分公式 $$ \int_{a}^b f(x) dx=\frac{b-a}{2} \int_{-1}^{1}f \left( \frac{b-a}{2}\xi+\frac{a+b}{2} \right) \, d\xi \\ \approx \frac{b-a}{2} \sum_{i}\omega_i f \left( \frac{b-a}{2}\xi_i+\frac{a+b}{2} \right) $$ 因此, 区间 $[a,b]$ 上的积分法则 $(x_i, w_i)$ 为 $$ x_i = \frac{b-a}{2}\xi_i+\frac{a+b}{2},\quad w_i = \frac{b-a}{2} \omega_i, $$ 即 $$ \int_{a}^b f(x) dx \approx \sum_{i}w_i f (x_i) . $$

附录 A: Gauss-Legendre

n	积分点 $\pm\xi_i$	权重 $\omega_i$
2	1.0	1.000000000000000
3	0	1.333333333333333
	1.0	0.333333333333333
4	0.447213595499958	0.833333333333333
	1.0	0.166666666666667
5	0.0	0.56888888888888889
	0.5384693101056830	0.47862867049936647
	0.90617984593866399	0.23692688505618909
6	0.2386191860831969	0.46791393457269105
	0.6612093864662645	0.36076157304813860
	0.9324695142031520	0.17132449237917035
7	0.0	0.41795918367346939
	0.4058451513773972	0.38183005050511894
	0.7415311855993944	0.27970539148927667
	0.9491079123427585	0.12948496616886969
8	0.1834346424956498	0.36268378337836198
	0.5255324099163290	0.31370664587788729
	0.7966664774136267	0.22238103445337447
	0.9602898564975362	0.10122853629037626

附录 B: Gauss-Lobatto

n	积分点 $\pm\xi_i$	权重 $\omega_i$
2	1.0	1.0
3	0.0	1.333333333333333
	1.0	0.333333333333333
4	0.447213595499958	0.833333333333333
	1.0	0.166666666666667
5	0.0	0.711111111111111
	0.654653670707977	0.544444444444444
	1.0	0.100000000000000
6	0.285231516480645	0.554858377035486
	0.765055323929465	0.378474956297847
	1.0	0.066666666666667
7	0.0	0.487619047619048
	0.468848793470714	0.431745381209863
	0.830223896278567	0.276826047361566
	1.0	0.047619047619048
8	0.209299217902479	0.412458794658704
	0.591700181433142	0.341122692483504
	0.871740148509607	0.210704227143506
	1.0	0.035714285714286

参考文献

分享到: