数值积分与微分
代数精度
基本定义,只要满足
∫
a
b
X
m
+
1
d
x
≠
∑
k
=
0
n
A
k
X
k
m
+
1
\int_a^bX^{m+1}dx\neq\sum_{k=0}^nA_kX_k^{m+1}
∫abXm+1dx=k=0∑nAkXkm+1
在m次以内都是相等的那么就称求积公式具有m次的代数精度
注:代数精度与特殊函数没有关系,只和选取的点以及这些点的权值有关
插值型求积公式
插值型求积函数其权系数 A i = ∫ a b l i d x A_i=\int_{a}^{b}l_idx Ai=∫ablidx
插值型求积公式就是利用拉个朗日或者牛顿插值多项式构造插值函数,而且我们知道拉格朗日和牛顿插值多项的余项可以表示为
R
n
(
x
)
=
f
n
+
1
(
ξ
)
(
n
+
1
)
!
ω
n
+
1
(
x
)
ω
n
+
1
(
x
)
=
(
x
−
x
0
)
(
x
−
x
1
)
.
.
.
(
x
−
x
n
)
R_n(x)=\frac{f^{n+1}(\xi)}{(n+1)!}\omega_{n+1}(x) \\[10pt]\omega_{n+1}(x)=(x-x_0)(x-x_1)...(x-x_n)
Rn(x)=(n+1)!fn+1(ξ)ωn+1(x)ωn+1(x)=(x−x0)(x−x1)...(x−xn)
有这个表达式可以推出,插值型求积公式的余项:
R
[
f
]
=
∫
a
b
f
n
+
1
(
ξ
)
(
n
+
1
)
!
ω
n
+
1
(
x
)
d
x
R[f]=\int_a^b\frac{f^{n+1}(\xi)}{(n+1)!}\omega_{n+1}(x)dx
R[f]=∫ab(n+1)!fn+1(ξ)ωn+1(x)dx
对于
f
(
x
)
f(x)
f(x)我们很容易得出当
f
(
x
)
f(x)
f(x)的次数等于n的时候我们经过n+1次求导积分内的表达式为0,故插值型求积公式至少具有n次代数精度
=>定理: I n = ∑ k = 0 n A k f ( x k ) I_n=\sum_{k=0}^{n}A_{k}f(x_k) In=∑k=0nAkf(xk)至少具有n次代数精度的充要条件是他是插值型的
求积公式的稳定性与收敛性
收敛性
基本定义随着步数增多,误差缩小,机械求积的值逐渐趋向于原本函数的值,那么这个求积公式就可以称之为收敛的
lim
h
→
0
∑
k
=
0
n
A
k
f
(
x
k
)
=
∫
a
b
f
(
x
)
d
x
\lim_{h\to0}\sum_{k=0}^nA_kf(x_k)=\int_a^bf(x)dx
h→0limk=0∑nAkf(xk)=∫abf(x)dx
其中
h
=
max
1
≤
i
≤
n
(
x
i
−
x
i
−
1
)
h=\max_{1\leq{i}\leq{n}}(x_i-x_{i-1})
h=1≤i≤nmax(xi−xi−1)
注:可以对比插值多项式,对于多个点的插值,当点过多的时候可能会产生龙格现象,误差随着插值点增加逐渐增加。但使用分段插值的时候可以避免误差的累积和龙格现象的产生。
稳定性
定义
∀
ϵ
,
∃
δ
对
∀
x
K
,
使
得
f
(
x
k
)
−
f
k
∣
<
δ
时
有
∣
∑
k
=
0
n
A
k
f
(
x
k
)
−
∑
k
=
0
n
A
k
f
k
∣
<
ϵ
\forall\epsilon,\exists\delta对\forall{x_K},使得f(x_k)-f_k|<\delta时有|\sum_{k=0}^{n}A_kf(x_k)-\sum_{k=0}^{n}A_kf_k|\lt\epsilon
∀ϵ,∃δ对∀xK,使得f(xk)−fk∣<δ时有∣k=0∑nAkf(xk)−k=0∑nAkfk∣<ϵ
则称这个机械求积是稳定的
=>当 A k > 0 A_k>0 Ak>0对任意k都成立的时候求积一定是稳定的,简单推理即可
牛顿-柯斯特公式
把积分区间等分成n个小区间,求的每个点的
f
(
x
)
f(x)
f(x)的值,利用拉格朗日插值公式求的插值函数,再对插值函数求积分
L
n
(
x
)
=
∑
k
=
0
n
f
(
x
k
)
l
k
(
x
)
余
项
:
R
n
(
x
)
=
f
(
n
+
1
)
(
x
)
(
n
+
1
)
!
ω
n
+
1
(
x
)
故
可
以
推
出
,
牛
顿
−
柯
斯
特
的
基
本
形
式
:
I
n
=
∑
k
=
0
n
f
(
x
)
∫
a
b
l
k
(
x
)
d
x
误
差
项
:
R
(
I
n
)
=
∫
a
b
R
n
(
x
)
d
x
L_n(x)=\sum_{k=0}^{n}f(x_k)l_k(x) \\[10pt] 余项:R_n(x)=\frac{f^{(n+1)}(x)}{(n+1)!}\omega_{n+1}(x) \\[10pt] 故可以推出,牛顿-柯斯特的基本形式:\\[10pt]I_n=\sum_{k=0}^{n}f(x)\int_{a}^{b}l_k(x)dx \\[10pt]误差项:R(I_n)=\int_{a}^{b}R_n(x)dx
Ln(x)=k=0∑nf(xk)lk(x)余项:Rn(x)=(n+1)!f(n+1)(x)ωn+1(x)故可以推出,牛顿−柯斯特的基本形式:In=k=0∑nf(x)∫ablk(x)dx误差项:R(In)=∫abRn(x)dx
由于牛顿-柯斯特公式是等距节点, x k , x h x_k,x_h xk,xh可以进行替换
最终形式:
I
n
(
f
)
=
∑
k
=
0
n
A
k
f
(
x
k
)
=
(
b
−
a
)
∑
k
=
0
n
C
k
(
n
)
f
(
x
k
)
C
k
n
=
(
−
1
)
n
−
k
n
∗
k
!
∗
(
n
−
k
)
!
∫
0
n
∏
0
≤
j
≤
n
,
j
≠
k
(
t
−
j
)
d
t
I_n(f)=\sum_{k=0}{n}A_kf(x_k)=(b-a)\sum_{k=0}^nC_k^{(n)}f(x_k) \\[10pt]C_k^{n}=\frac{(-1)^{n-k}}{n*k!*(n-k)!}\int_0^n\prod_{0\leq{j}\leq{n},j\neq{k}}(t-j)dt
In(f)=k=0∑nAkf(xk)=(b−a)k=0∑nCk(n)f(xk)Ckn=n∗k!∗(n−k)!(−1)n−k∫0n0≤j≤n,j=k∏(t−j)dt
=> 牛顿-柯斯特公式至少具有n+1次代数精度
低阶牛顿-柯斯特公式,及其余项
低阶牛顿-柯斯特公式就是把n取低阶的时候得到的公式,牛顿-柯斯特公式在原有的插值求积的基础上增加了,等分区间这一条件
再复习一遍,牛顿-柯斯特公式的求积公式的系数
C
k
(
n
)
=
(
−
1
)
(
n
−
k
)
n
∗
k
!
∗
(
n
−
k
)
!
∫
0
n
∏
0
≤
j
≤
n
,
j
≠
k
(
t
−
j
)
d
t
C_k^{(n)}=\frac{(-1)^{(n-k)}}{n*k!*(n-k)!}\int_0^n\prod_{0\leq{j}\leq{n},j\neq{k}}(t-j)dt
Ck(n)=n∗k!∗(n−k)!(−1)(n−k)∫0n0≤j≤n,j=k∏(t−j)dt
梯形公式 n=1
n = 1 , C 0 = 1 2 , C 1 = 1 2 I n = ( b − a ) ( 1 2 f ( a ) + 1 2 f ( x b ) ) R ( I n ) = ∫ a b f ′ ′ ( x ) 2 ( x − a ) ( x − b ) d x = − 1 12 ( b − a ) 3 f ′ ′ ( x ) n=1 , C_0=\frac{1}{2}, C_1=\frac{1}{2} \\[10pt] I_n=(b-a)(\frac{1}{2}f(a)+\frac{1}{2}f(x_b)) \\[10pt]R(I_n)=\int_a^b\frac{f^{''}(x)}{2}(x-a)(x-b)dx=-\frac{1}{12}(b-a)^3f^{''}(x) n=1,C0=21,C1=21In=(b−a)(21f(a)+21f(xb))R(In)=∫ab2f′′(x)(x−a)(x−b)dx=−121(b−a)3f′′(x)
可以看出梯形公式具有一次代数精度
Simpson公式(三点公式) n=2
n = 3 , C 0 ( 2 ) = 1 6 , C 1 ( 2 ) = 1 4 , C 2 ( 2 ) = 1 6 I n = ( b − a ) ( 1 6 f ( a ) + 4 6 f ( a + b 2 ) + 1 6 f ( b ) ) R ( I n ) = − b − a 180 ( b − a 2 ) 4 f ( 4 ) ( η ) n=3,C_0^{(2)}=\frac{1}{6},C_1^{(2)}=\frac{1}{4},C_2^{(2)}=\frac{1}{6} \\[10pt]I_n=(b-a)(\frac{1}{6}f(a)+\frac{4}{6}f(\frac{a+b}{2})+\frac{1}{6}{f(b)}) \\[10pt]R(I_n)=-\frac{b-a}{180}(\frac{b-a}{2})^4f^{(4)}(\eta) n=3,C0(2)=61,C1(2)=41,C2(2)=61In=(b−a)(61f(a)+64f(2a+b)+61f(b))R(In)=−180b−a(2b−a)4f(4)(η)
Simpson三点公式具有三次代数精度
柯斯特求积公式及其余项 n=4
n = 4 , C 0 ( 4 ) = 7 90 , C 1 ( 4 ) = 32 90 , C 2 ( 4 ) = 12 90 , C 3 ( 4 ) = 32 90 , C 4 ( 4 ) = 7 90 I 4 = b − a 90 [ 7 f ( x 0 ) + 12 f ( x 1 ) + 32 f ( x 2 ) + 12 f ( x 3 ) + 7 f ( x 4 ) ] R ( I n ) = − 2 ( b − a ) 945 ( b − a 4 ) 6 f ( 6 ) ( η ) n=4,C_0^{(4)}=\frac{7}{90},C_1^{(4)}=\frac{32}{90},C_2^{(4)}=\frac{12}{90},C_3^{(4)}=\frac{32}{90},C_4^{(4)}=\frac{7}{90} \\[10pt]I_4=\frac{b-a}{90}[7f(x_0)+12f(x_1)+32f(x_2)+12f(x_3)+7f(x_4)] \\[10pt]R(I_n)=-\frac{2(b-a)}{945}(\frac{b-a}{4})^6f^{(6)}(\eta) n=4,C0(4)=907,C1(4)=9032,C2(4)=9012,C3(4)=9032,C4(4)=907I4=90b−a[7f(x0)+12f(x1)+32f(x2)+12f(x3)+7f(x4)]R(In)=−9452(b−a)(4b−a)6f(6)(η)
柯斯特公式具有五次代数精度
牛顿-柯斯特公式的稳定性
由插值公式稳定性推理可得到只要满足权函数 ∀ i ∈ { 0 , 1 , . . . , n } , A i > 0 \forall{i\in\lbrace0,1,...,n\rbrace},A_i>0 ∀i∈{0,1,...,n},Ai>0插值求积公式都是稳定的
=> n < 8 n\lt8 n<8的时候牛顿克斯特公式均是稳定的
但牛顿-柯斯特公式中 C k C_k Ck不保证大于0的时候稳定性无法得到保证,不使用不稳定的高阶牛顿-柯斯特公式
复化求积公式
复化求积公式,把积分区间等分成n个小区间,每个小区间使用牛顿-柯斯特公式进行积分,最后将所有小段的积分累加,近似总区间段的积分
针对于每个小区间等分成l份,那么步长为
h
l
\frac{h}{l}
lh,对应每个节点的值分别为
x
k
,
x
k
+
h
l
,
.
.
.
,
x
k
+
l
h
l
=
x
k
+
1
x_k,x_k+\frac{h}{l},...,x_k+\frac{lh}{l}=x_{k+1}
xk,xk+lh,...,xk+llh=xk+1
[
x
k
,
x
k
+
1
]
[x_k,x_{k+1}]
[xk,xk+1]区间上的的积分:
∫
x
k
x
k
+
1
f
(
x
)
d
x
=
h
∑
i
=
0
l
C
i
(
n
)
f
(
x
k
+
i
h
l
)
\int_{x_k}^{x_{k+1}}f(x)dx=h\sum_{i=0}^{l}C_i^{(n)}f(x_k+\frac{ih}{l})
∫xkxk+1f(x)dx=hi=0∑lCi(n)f(xk+lih)
I
n
I_n
In等于n个小区间积分相加的总和
I
n
=
h
∑
k
=
0
n
−
1
∑
i
=
0
l
C
i
(
l
)
f
(
x
k
+
i
h
l
)
I_n=h\sum_{k=0}^{n-1}\sum_{i=0}^{l}C_i^{(l)}f(x_k+\frac{ih}{l})
In=hk=0∑n−1i=0∑lCi(l)f(xk+lih)
注意
C
k
(
n
)
=
(
−
1
)
(
n
−
k
)
n
∗
k
!
∗
(
n
−
k
)
!
∫
0
n
∏
0
≤
j
≤
n
,
j
≠
k
(
t
−
j
)
d
t
C_k^{(n)}=\frac{(-1)^{(n-k)}}{n*k!*(n-k)!}\int_0^n\prod_{0\leq{j}\leq{n},j\neq{k}}(t-j)dt
Ck(n)=n∗k!∗(n−k)!(−1)(n−k)∫0n0≤j≤n,j=k∏(t−j)dt
C
k
(
n
)
C_k^{(n)}
Ck(n)与区间的起始点
x
k
x_k
xk和
f
(
x
)
f(x)
f(x)的值无关,每段上面的
C
k
(
n
)
C_k^{(n)}
Ck(n)都满足同意序列,故可以先求
C
k
(
n
)
C_k^{(n)}
Ck(n)再带入公式求解
复化梯形公式
在每个小区间上使用梯形公式求积,l=1
T
n
=
b
−
a
2
n
[
f
(
a
)
+
∑
k
=
1
n
−
1
f
(
x
k
)
+
f
(
b
)
]
T_n=\frac{b-a}{2n}[f(a)+\sum_{k=1}^{n-1}f(x_k)+f(b)]
Tn=2nb−a[f(a)+k=1∑n−1f(xk)+f(b)]
复化辛普森公式
在每个小区间上使用辛普森公式求积,l=2
S
n
=
b
−
a
6
n
[
f
(
a
)
+
4
∑
k
=
1
n
−
1
f
(
x
k
+
1
2
)
+
2
∑
k
=
1
n
−
1
f
(
x
k
)
+
f
(
b
)
]
S_n=\frac{b-a}{6n}[f(a)+4\sum_{k=1}^{n-1}f(x_{k+\frac{1}{2}})+2\sum_{k=1}^{n-1}f(x_k)+f(b)]
Sn=6nb−a[f(a)+4k=1∑n−1f(xk+21)+2k=1∑n−1f(xk)+f(b)]
复化柯斯特公式
在每个小区间上使用柯斯特公式,l=4
S
n
=
b
−
a
90
n
[
7
f
(
a
)
+
∑
k
=
0
n
−
1
[
32
f
(
x
k
+
1
4
)
+
12
f
(
x
k
+
2
4
)
+
32
f
(
x
k
+
3
4
)
]
+
14
∑
k
=
1
n
−
1
f
(
x
k
)
+
7
f
(
b
)
]
S_n=\frac{b-a}{90n}[7f(a)+\sum_{k=0}^{n-1}[32f(x_{k+\frac{1}{4}})+12f(x_{k+\frac{2}{4}})+32f(x_{k+\frac{3}{4}})]+14\sum_{k=1}^{n-1}f(x_k)+7f(b)]
Sn=90nb−a[7f(a)+k=0∑n−1[32f(xk+41)+12f(xk+42)+32f(xk+43)]+14k=1∑n−1f(xk)+7f(b)]
注:柯斯特复化求积公式收敛速度最快且精度最高,复化辛普森公式次之,复化梯形公式最坏
复化求积公式的余项分别为
复化梯形求积公式
I − T n = − b − a 12 f ′ ′ ( η ) ≈ − 1 12 h 2 [ f ′ ( b ) − f ′ ( a ) ] = O ( h 2 ) I-T_n=-\frac{b-a}{12}f''(\eta)\approx-\frac{1}{12}h^2[f'(b)-f'(a)]=O(h^2) I−Tn=−12b−af′′(η)≈−121h2[f′(b)−f′(a)]=O(h2)
复化辛普森公式
I − S n = − b − a 180 ( h 2 ) 4 f ( 4 ) ( η ) ≈ − 1 180 ( h 2 ) 4 [ f ′ ′ ′ ( b ) − f ′ ′ ′ ( a ) ] = O ( h 4 ) I-S_n=-\frac{b-a}{180}(\frac{h}{2})^4f^{(4)}(\eta)\approx-\frac{1}{180}(\frac{h}{2})^4[f'''(b)-f'''(a)]=O(h^4) I−Sn=−180b−a(2h)4f(4)(η)≈−1801(2h)4[f′′′(b)−f′′′(a)]=O(h4)
复化柯斯特公式
I − C n = − 2 ( b − a ) 945 ( h 4 ) 6 f ( 6 ) ( η ) ≈ − 2 945 ( h 4 ) 6 [ f ( 5 ) ( b ) − f ( 5 ) ( a ) ] = O ( h 6 ) I-C_n=-\frac{2(b-a)}{945}(\frac{h}{4})^6f^{(6)}(\eta)\approx-\frac{2}{945}(\frac{h}{4})^6[f^{(5)}(b)-f^{(5)}(a)]=O(h^6) I−Cn=−9452(b−a)(4h)6f(6)(η)≈−9452(4h)6[f(5)(b)−f(5)(a)]=O(h6)
收敛阶
如果满足
lim
h
→
0
I
−
I
n
h
p
=
c
\lim_{h\to0}\frac{I-I_n}{h^p}=c
h→0limhpI−In=c
则称复化求积公式是
p
p
p阶收敛的,
p
p
p称为
I
n
I_n
In的收敛阶
对于刚才的复化求积余项我们可以得到复化梯形公式,复化辛普森公式,复化柯斯特公式的收敛阶分别为2,4,6阶
( p s ps ps:这里可以观察到,对于牛顿-柯斯特求积公式的基本形式,梯形公式,辛普森公式,柯斯特公式的的代数精度分别是1,3,5)
龙贝格求积公式
梯形求积公式,实在复化求积公式的基础上采取变步长来确定原积分区间应该化为多少等分
设
原
有
区
间
分
为
n
个
等
份
,
具
有
n
+
1
个
节
点
对
于
指
定
区
间
[
x
k
,
x
k
+
1
]
增
加
中
分
点
x
k
+
1
2
,
对
分
割
后
的
区
间
使
用
复
化
梯
形
求
积
公
式
S
[
x
k
,
x
k
+
1
]
=
h
4
[
f
(
x
k
)
+
2
f
(
x
k
+
1
2
)
+
f
(
x
k
+
1
)
]
则
总
区
间
积
分
经
过
化
简
以
后
有
T
2
n
=
1
2
T
n
+
h
2
∑
k
=
0
n
−
1
f
[
a
+
(
2
k
+
1
)
h
2
]
设原有区间分为n个等份,具有n+1个节点 \\[10pt] 对于指定区间[x_k,x_{k+1}]增加中分点x_{k+\frac{1}{2}},对分割后的区间使用复化梯形求积公式 \\[10pt] S_{[x_k,x_{k+1}]}=\frac{h}{4}[f(x_k)+2f(x_{k+\frac{1}{2}})+f(x_{k+1})] \\[10pt] 则总区间积分经过化简以后有 \\[10pt]T_{2n}=\frac{1}{2}T_n+\frac{h}{2}\sum_{k=0}^{n-1}f[a+(2k+1)\frac{h}{2}]
设原有区间分为n个等份,具有n+1个节点对于指定区间[xk,xk+1]增加中分点xk+21,对分割后的区间使用复化梯形求积公式S[xk,xk+1]=4h[f(xk)+2f(xk+21)+f(xk+1)]则总区间积分经过化简以后有T2n=21Tn+2hk=0∑n−1f[a+(2k+1)2h]
其中a表示的是积分区间的其实后面的表达式,表示所有新添加的点的函数值累加
∗
h
2
*\frac{h}{2}
∗2h加上原梯形公式的二分之一,h会随着公式的迭代不断改变
改进的梯形求积公式
已知梯形求积公式的余项
R
n
(
f
)
=
−
b
−
a
12
h
2
f
′
′
(
η
)
R_n(f)=-\frac{b-a}{12}h^2f''(\eta)
Rn(f)=−12b−ah2f′′(η)
在把每个小区间等分后,区间个数翻倍,对应的余项
R
2
n
(
f
)
=
−
b
−
a
12
(
h
2
)
2
f
′
′
(
η
‾
)
R_{2n}(f)=-\frac{b-a}{12}(\frac{h}{2})^2f''(\overline\eta)
R2n(f)=−12b−a(2h)2f′′(η)
我们令
f
′
′
(
η
)
=
f
′
′
(
η
‾
)
f''(\eta)=f''(\overline\eta)
f′′(η)=f′′(η)就有
R
n
=
4
R
2
n
R
2
n
=
T
n
−
T
2
n
3
R_n=4R_{2n} \\[10pt]R_{2n}=\frac{T_n-T_{2n}}{3}
Rn=4R2nR2n=3Tn−T2n
而
R
2
n
R_{2n}
R2n表示2n区间的积分误差
再将所求的的误差加到
T
2
n
T_{2n}
T2n中那么就得到改进的梯形求积公式
T
‾
2
n
=
T
2
n
+
1
3
(
T
n
−
T
2
n
)
=
4
3
T
2
n
−
T
n
=
h
6
[
f
(
a
)
+
4
∑
k
=
1
n
−
1
f
(
x
k
+
1
2
)
+
2
∑
k
=
1
n
−
1
f
(
x
k
)
+
f
(
b
)
]
\overline{T}_{2n}=T_{2n}+\frac{1}{3}(T_n-T_{2n})=\frac{4}{3}T_{2n}-T_n\\=\frac{h}{6}[f(a)+4\sum_{k=1}^{n-1}f(x_{k+\frac{1}{2}})+2\sum_{k=1}^{n-1}f(x_k)+f(b)]
T2n=T2n+31(Tn−T2n)=34T2n−Tn=6h[f(a)+4k=1∑n−1f(xk+21)+2k=1∑n−1f(xk)+f(b)]
在误差加上后,改进的梯形求积公式实际上得到的就是复化辛普森公式
同理,利用复化辛普森公式前后两个值的关系我们可以得到复化柯斯特公式
C
n
=
16
15
S
2
n
−
1
15
S
n
C_n=\frac{16}{15}S_{2n}-\frac{1}{15}S_n
Cn=1516S2n−151Sn
龙贝格求积公式
龙贝格求积公式就是通过复化柯斯特公式,二分前后两个只得到的
R
n
=
64
63
S
2
n
−
1
63
S
n
R_n=\frac{64}{63}S_{2n}-\frac{1}{63}S_n
Rn=6364S2n−631Sn
综上所述我们拥有了一种从梯形公式逐步求到龙贝格公示的方法
k | 区间等分数 n = 2 k n=2^k n=2k | 梯形公式 | 辛普森公式 | 柯斯特公式 | 龙贝格公式 |
---|---|---|---|---|---|
0 | 2 0 = 1 2^0=1 20=1 | T 1 T_1 T1 | |||
1 | 2 1 = 2 2^1=2 21=2 | T 2 T_2 T2 | S 1 S_1 S1 | ||
2 | 2 2 = 4 2^2=4 22=4 | T 4 T_4 T4 | S 2 S_2 S2 | C 1 C_1 C1 | |
3 | 2 3 = 8 2^3=8 23=8 | T 8 T_8 T8 | S 4 S_4 S4 | C 2 C_2 C2 | R 1 R_1 R1 |
2 4 = 16 2^4=16 24=16 | T 16 T_{16} T16 | S 8 S_8 S8 | C 4 C_4 C4 | R 2 R_2 R2 | |
5 | 2 5 = 32 2^5=32 25=32 | T 32 T_{32} T32 | S 16 S_{16} S16 | C 8 C_8 C8 | R 4 R_4 R4 |
… |
理查森外推加速法
理查森外推加速法的基本公式
T
m
(
h
)
=
4
m
4
m
−
1
T
m
−
1
(
h
2
)
−
1
4
m
−
1
T
m
−
1
(
h
)
T_m(h)^=\frac{4^m}{4^m-1}T_{m-1}(\frac{h}{2})-\frac{1}{4^m-1}T_{m-1}(h)
Tm(h)=4m−14mTm−1(2h)−4m−11Tm−1(h)
这里m表示第m次加速,如
T
m
T_m
Tm表示经过m次理查森外推法加速梯形的值
同时我们也知道 T m ( h 2 ) T_{m}(\frac{h}{2}) Tm(2h)意味着在将原等分长度h在折半,故可以得到更加常见的表达式
T
m
(
k
)
=
4
m
4
m
−
1
T
m
−
1
(
k
+
1
)
−
1
4
m
−
1
T
m
−
1
(
k
)
T_m^{(k)}=\frac{4^m}{4^m-1}T_{m-1}^{(k+1)}-\frac{1}{4^m-1}T_{m-1}^{(k)}
Tm(k)=4m−14mTm−1(k+1)−4m−11Tm−1(k)
上标k表示经过第k次二分
下标m表示经过第m次加速
高斯求积公式
高斯公式的引言
高斯求积公式一般性理论
对于 n + 1 n+1 n+1个节点的求积公式,最高可得 2 n + 1 2n+1 2n+1代数精度的求积公式,这类具有最高求积精度的求积公式称为高斯求积公式
注:高斯求积公式,用到的高斯点具有特殊性
高斯点以及求积系数
=> 高斯点定理:对于节点
a
≤
x
0
≤
x
1
≤
.
.
.
≤
x
n
≤
b
a\leq{x_0}\leq{x_{1}}\leq{...}\leq{x_n}\leq{b}
a≤x0≤x1≤...≤xn≤b是高斯点的充要条件是以这些节点为零点的零点多项式
ω
n
+
1
(
x
)
\omega_{n+1}(x)
ωn+1(x)与次数不超过n的多项式都正交
ω
n
+
1
(
x
)
=
(
x
−
x
0
)
(
x
−
x
1
)
⋯
(
x
−
x
n
)
\omega_{n+1}(x)=(x-x_0)(x-x_1)\cdots(x-x_n)
ωn+1(x)=(x−x0)(x−x1)⋯(x−xn)
又已知高斯求积公式具有2n+1次代数精度对于
f
(
x
)
=
x
,
⋯
,
x
2
n
+
1
f(x)=x,\cdots,x^{2n+1}
f(x)=x,⋯,x2n+1都能准确成立故可以列出2n+1个方程组,又因为对于求积系数A有n+1个,只要其中的n+1个等式既可以列出方程组
矩阵形式
x
k
,
k
∈
{
0
,
1
,
⋯
,
n
}
x_k,k\in\{0,1,\cdots,n\}
xk,k∈{0,1,⋯,n}已知,
A
k
,
k
∈
{
0
,
1
⋯
,
n
}
A_k,k\in\{0,1\cdots,n\}
Ak,k∈{0,1⋯,n}是待求量
(
1
1
⋯
1
x
1
x
2
⋯
x
n
(
x
1
)
2
(
x
2
)
2
⋯
(
x
n
)
2
⋮
⋮
⋯
⋮
(
x
1
)
n
(
x
2
)
n
⋯
(
x
n
)
2
)
∗
(
A
0
A
1
A
2
⋮
A
n
)
=
(
∫
a
b
1
d
x
∫
a
b
x
d
x
∫
a
b
x
2
d
x
⋮
∫
a
b
x
n
d
x
)
\begin{pmatrix} 1 & 1 & \cdots & 1 \\ x_1 & x_2 & \cdots & x_n \\ (x_1)^2 & (x_2)^2 & \cdots &(x_n)^2\\ \vdots & \vdots & \cdots &\vdots \\ (x_1)^n & (x_2)^n & \cdots & (x_n)^2 \\ \end{pmatrix}* \begin{pmatrix} A_0\\ A_1 \\ A_2 \\ \vdots \\ A_n \end{pmatrix}= \begin{pmatrix} \int_a^b1dx\\ \int_a^bxdx\\ \int_a^bx^2dx\\ \vdots\\ \int_a^bx^ndx \end{pmatrix}
⎝⎜⎜⎜⎜⎜⎛1x1(x1)2⋮(x1)n1x2(x2)2⋮(x2)n⋯⋯⋯⋯⋯1xn(xn)2⋮(xn)2⎠⎟⎟⎟⎟⎟⎞∗⎝⎜⎜⎜⎜⎜⎛A0A1A2⋮An⎠⎟⎟⎟⎟⎟⎞=⎝⎜⎜⎜⎜⎜⎜⎛∫ab1dx∫abxdx∫abx2dx⋮∫abxndx⎠⎟⎟⎟⎟⎟⎟⎞
高斯求积公式的余项
R n [ f ] = f 2 n + 2 ( η ) ( 2 n + 2 ) ! ∫ a b ω n + 1 2 ( x ) ρ ( x ) d x R_n[f]=\frac{f^{2n+2}(\eta)}{(2n+2)!}\int_a^b\omega^2_{n+1}(x)\rho(x)dx Rn[f]=(2n+2)!f2n+2(η)∫abωn+12(x)ρ(x)dx
高斯求积公式的稳定性与收敛性
=> 高斯求积公式的求积系数都是正的,高斯求积公式是稳定的
=> 高斯求积公式是收敛的
常用的高斯求积公式
高斯-勒让德求积公式
勒让德正交多项式的基本表达式
P 0 ( x ) = 1 P n ( x ) = 1 2 n n ! d n d x n ( x 2 − 1 ) n 或 P n ( x ) = n ! ( 2 n ) ! d n d x n [ ( x 2 − 1 ) n ] P_0(x)=1\\ P_n(x)=\frac{1}{2^nn!}\frac{d^n}{dx^n}(x^2-1)^n\\ 或P_n(x)=\frac{n!}{(2n)!}\frac{d^n}{dx^n}[(x^2-1)^n] P0(x)=1Pn(x)=2nn!1dxndn(x2−1)n或Pn(x)=(2n)!n!dxndn[(x2−1)n]
勒让德正交多项式的特性
=>1. 正交性 , 勒让德正交多项式满足
KaTeX parse error: No such environment: equation at position 8: \begin{̲e̲q̲u̲a̲t̲i̲o̲n̲}̲\int_{-1}^{1}P_…
=>2. 奇偶性
P
n
(
x
)
=
(
−
1
)
n
P
n
(
x
)
P_n(x)=(-1)^nP_n(x)
Pn(x)=(−1)nPn(x)
=>3. 递推关系 ( n + 1 ) P n + 1 ( x ) = ( 2 n + 1 ) x P n ( x ) − n P n − 1 ( x ) (n+1)P_{n+1}(x)=(2n+1)xP_n(x)-nP_{n-1}(x) (n+1)Pn+1(x)=(2n+1)xPn(x)−nPn−1(x)
=>4. P n ( x ) P_n(x) Pn(x)在区间[-1,1]内有n个不同的实0点
常用的勒让德正交多项式
P 0 ( x ) = 1 P 1 ( x ) = x P 2 ( x ) = ( 3 x 2 − 1 ) 2 P 3 ( x ) = ( 5 x 3 − 3 x ) 2 P 4 ( x ) = ( 35 x 4 − 30 x 2 + 3 ) 8 P 5 ( x ) = 63 x 5 − 70 x 3 + 15 x 8 P 6 ( x ) = ( 231 x 6 − 351 x 4 + 105 x 2 − 5 ) 16 P_0(x)=1\\ P_1(x)=x\\ P_2(x)=\frac{(3x^2-1)}{2}\\ P_3(x)=\frac{(5x^3-3x)}{2}\\ P_4(x)=\frac{(35x^4-30x^2+3)}{8}\\ P_5(x)=\frac{63x^5-70x^3+15x}{8}\\ P_6(x)=\frac{(231x^6-351x^4+105x^2-5)}{16}\\ P0(x)=1P1(x)=xP2(x)=2(3x2−1)P3(x)=2(5x3−3x)P4(x)=8(35x4−30x2+3)P5(x)=863x5−70x3+15xP6(x)=16(231x6−351x4+105x2−5)
因为 P o ( x ) = 1 P_o(x)=1 Po(x)=1对于区 P 1 ( x ) = x P_1(x)=x P1(x)=x与任意次数不超过1的多项式都正交
对于 P 2 ( x ) 与 P 1 ( x ) , P 0 ( x ) P_2(x)与P_1(x),P_0(x) P2(x)与P1(x),P0(x)都正交,任意次数不超过2的多项式都可由 P 1 , P 0 P_1,P_0 P1,P0线性表出均与 P 2 P_2 P2正交,并且由性质4可以知道 P 2 ( x ) P_2(x) P2(x)具有两个不同的实零点, P 2 ( x ) = ( x − x 0 ) ( x — x 1 ) P_2(x)=(x-x_0)(x—x_1) P2(x)=(x−x0)(x—x1)成立, x 0 , x 1 x_0,x_1 x0,x1为高斯点
同理 P 3 P_3 P3可以得到3个高斯点
常用的高斯勒让德公式
1 点 的 高 斯 − 勒 让 德 公 式 ∫ − 1 1 f ( x ) d x = 2 f ( 0 ) 2 点 的 高 斯 − 勒 让 德 公 式 ∫ − 1 1 f ( x ) d x = f ( − 1 3 ) + f ( 1 3 ) 3 点 的 高 斯 − 勒 让 德 公 式 ∫ − 1 1 f ( x ) d x = 5 9 f ( − 15 5 ) + 8 9 f ( 0 ) + 5 9 f ( 15 5 ) \\1点的高斯-勒让德公式 \quad \int_{-1}^{1}f(x)dx=2f(0) \\2点的高斯-勒让德公式 \quad \int_{-1}^1f(x)dx=f(-\frac{1}{\sqrt{3}})+f(\frac{1}{\sqrt{3}}) \\3点的高斯-勒让德公式 \quad \int_{-1}^1f(x)dx=\frac{5}{9}f(-\frac{\sqrt{15}}{5})+\frac{8}{9}f(0)+\frac{5}{9}f(\frac{\sqrt{15}}{5}) 1点的高斯−勒让德公式∫−11f(x)dx=2f(0)2点的高斯−勒让德公式∫−11f(x)dx=f(−31)+f(31)3点的高斯−勒让德公式∫−11f(x)dx=95f(−515)+98f(0)+95f(515)
ps:无论是高斯插值点或者高斯求积系数都与具体的函数f(x)无关
n=高斯点数-1
高斯勒让德公式的余项
R n [ f ] = f ( 2 n + 2 ) ( η ) ( 2 n + 2 ) ! ∫ − 1 1 P n + 1 2 ( x ) d x η ∈ [ − 1 , 1 ] R_n[f]=\frac{f^{(2n+2)}(\eta)}{(2n+2)!}\int_{-1}^1P_{n+1}^2(x)dx \quad \eta\in[-1,1] Rn[f]=(2n+2)!f(2n+2)(η)∫−11Pn+12(x)dxη∈[−1,1]
高斯-切比雪夫积分公式
高斯-切比雪夫积分公式的基本表达式
T n ( x ) = cos ( n arccos x ) T_n(x)=\cos(n\arccos{x}) Tn(x)=cos(narccosx)
高斯-切比雪夫多项式的特性
=>1. 递推关系 T n + 1 ( x ) = 2 x T n ( x ) − T n − 1 ( x ) , 其 中 T 0 ( x ) = 1 , T 1 ( x ) = x T_{n+1}(x)=2xT_n(x)-T_{n-1}(x),其中T_0(x)=1,T_1(x)=x Tn+1(x)=2xTn(x)−Tn−1(x),其中T0(x)=1,T1(x)=x
=>2. 切比雪夫正交多项式的值
KaTeX parse error: No such environment: equation at position 8: \begin{̲e̲q̲u̲a̲t̲i̲o̲n̲}̲\int_{-1}^1\fra…
=>3.
T
2
k
(
x
)
T_{2k}(x)
T2k(x)只包含
x
x
x的偶次幂,
T
2
k
+
1
(
x
)
T_{2k+1}(x)
T2k+1(x)只包含
x
x
x的奇次幂
=>4. T n ( x ) T_n(x) Tn(x)在[-1,1]区间上具有n个零点
常用的切比雪夫多项式
T 0 ( x ) = 1 T 1 ( x ) = x T 2 ( x ) = 2 x 2 − 1 T 3 ( x ) = 4 x 3 − 3 x T 4 ( x ) = 8 x 4 − 8 x 2 + 1 T 5 ( x ) = 16 x 5 − 20 x 3 + 5 x T 6 ( x ) = 32 x 6 − 48 x 4 + 18 x 2 − 1 T_0(x)=1\\ T_1(x)=x\\ T_2(x)=2x^2-1\\ T_3(x)=4x^3-3x\\ T_4(x)=8x^4-8x^2+1\\ T_5(x)=16x^5-20x^3+5x\\ T_6(x)=32x^6-48x^4+18x^2-1\\ T0(x)=1T1(x)=xT2(x)=2x2−1T3(x)=4x3−3xT4(x)=8x4−8x2+1T5(x)=16x5−20x3+5xT6(x)=32x6−48x4+18x2−1
高斯-切比雪夫多项式的零点
x k = cos ( 2 k + 1 2 n + 2 π ) k = 0 , 1 , ⋯ , n 求 积 公 式 ∫ − 1 1 f ( x ) 1 − x 2 d x ≈ ∑ k = 0 n A k f ( x k ) 求 积 系 数 A k = π n + 1 x_k=\cos(\frac{2k+1}{2n+2}\pi) \quad k=0,1,\cdots,n \\求积公式\quad \int_{-1}^{1}\frac{f(x)}{\sqrt{1-x^2}}dx\approx\sum_{k=0}^{n}A_kf(x_k) \\求积系数 \quad A_k=\frac{\pi}{n+1} xk=cos(2n+22k+1π)k=0,1,⋯,n求积公式∫−111−x2f(x)dx≈k=0∑nAkf(xk)求积系数Ak=n+1π
高斯切比雪夫多项式余项
R [ f ] = 2 π 2 2 n ( 2 n ) ! f ( 2 n ) ( η ) η ∈ ( − 1 , 1 ) R[f]=\frac{2\pi}{2^{2n}(2n)!}f^{(2n)}(\eta) \quad \eta\in{(-1,1)} R[f]=22n(2n)!2πf(2n)(η)η∈(−1,1)
数值微分
中点方法与误差分析
f + ′ ( a ) ≈ f ( a + h ) − f ( a ) h f − ′ ( a ) ≈ f ( a − h ) − f ( a ) h f ′ ( a ) ≈ f ( a + h ) − f ( a − h ) 2 h 中 点 公 式 f'_+(a)\approx\frac{f(a+h)-f(a)}{h}\\ f'_-(a)\approx\frac{f(a-h)-f(a)}{h}\\ f'(a)\approx\frac{f(a+h)-f(a-h)}{2h} \quad 中点公式 f+′(a)≈hf(a+h)−f(a)f−′(a)≈hf(a−h)−f(a)f′(a)≈2hf(a+h)−f(a−h)中点公式
其中h表示步长,步长越小计算越精确但步长过小,会导致两个相近的数相减,有效数字损失
最优步长 h o p t h_{opt} hopt
h o p t = 3 3 ϵ M ϵ = max { ∣ ϵ 1 ∣ , ∣ ϵ 2 ∣ } ϵ 1 , ϵ 2 分 别 表 示 f ( a + h ) , f ( a − h ) 的 舍 入 误 差 M = max ∣ x − a ∣ ≤ h ∣ f ′ ′ ′ ( x ) ∣ h_{opt}=^3\sqrt{\frac{3\epsilon}{M}} \\ \epsilon=\max\{|\epsilon_1|,|\epsilon_2|\} \quad \epsilon_1,\epsilon_2分别表示f(a+h),f(a-h)的舍入误差 \\M=\max_{|x-a|\leq{h}}|f'''(x)| hopt=3M3ϵϵ=max{∣ϵ1∣,∣ϵ2∣}ϵ1,ϵ2分别表示f(a+h),f(a−h)的舍入误差M=∣x−a∣≤hmax∣f′′′(x)∣
插值型求导公式
插值型求导公式根据已知点构造插值多项式拉格朗日插值法,牛顿插值法
在对插值函数求导,求得导函数,计算已知点的导数即可得到微分
插值求导公式的微分公式的余项
R ( f n ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ω n + 1 ′ ( x k ) R(f_n)=\frac{f^{(n+1)}(\xi)}{(n+1)!}\omega'_{n+1}(x_k) R(fn)=(n+1)!f(n+1)(ξ)ωn+1′(xk)
两点公式
x
0
,
x
1
x_0,x_1
x0,x1进行插值,对应求导
P
1
′
(
x
)
=
1
h
[
−
f
(
x
0
)
+
f
(
x
1
)
]
P
1
′
(
x
0
)
=
1
h
[
f
(
x
1
)
−
f
(
x
0
)
]
P'_{1}(x)=\frac{1}{h}[-f(x_0)+f(x_1)] \\ P'_1(x_0)=\frac{1}{h}[f(x_1)-f(x_0)]
P1′(x)=h1[−f(x0)+f(x1)]P1′(x0)=h1[f(x1)−f(x0)]
x
0
,
x
1
x_0,x_1
x0,x1的导数值
P
1
′
(
x
0
)
=
1
h
[
f
(
x
1
)
−
f
(
x
0
)
]
P
1
′
(
x
1
)
=
1
h
[
f
(
x
1
)
−
f
(
x
0
)
]
P'_1(x_0)=\frac{1}{h}[f(x_1)-f(x_0)] \\ P'_1(x_1)=\frac{1}{h}[f(x_1)-f(x_0)]
P1′(x0)=h1[f(x1)−f(x0)]P1′(x1)=h1[f(x1)−f(x0)]
x
0
,
x
1
x_0,x_1
x0,x1的余项分别是
f
′
(
x
0
)
=
1
h
[
f
(
x
1
)
−
f
(
x
0
)
]
−
h
2
f
′
′
(
ξ
)
f
′
(
x
1
)
=
1
h
[
f
(
x
1
)
−
f
(
x
0
)
]
+
h
2
f
′
′
(
ξ
)
f'(x_0)=\frac{1}{h}[f(x_1)-f(x_0)]-\frac{h}{2}f''(\xi) \\ f'(x_1)=\frac{1}{h}[f(x_1)-f(x_0)]+\frac{h}{2}f''(\xi) \\
f′(x0)=h1[f(x1)−f(x0)]−2hf′′(ξ)f′(x1)=h1[f(x1)−f(x0)]+2hf′′(ξ)
三点公式
P 2 ′ ( x 0 + t h ) = 1 2 h [ ( 2 t − 3 ) f ( x 0 ) − ( 4 t − 4 ) f ( x 1 ) + ( 2 t − 1 ) f ( x 2 ) ] P'_2(x_0+th)=\frac{1}{2h}[(2t-3)f(x_0)-(4t-4)f(x_1)+(2t-1)f(x_2)] P2′(x0+th)=2h1[(2t−3)f(x0)−(4t−4)f(x1)+(2t−1)f(x2)]
x
0
,
x
1
,
x
2
x_0,x_1,x_2
x0,x1,x2的值
P
2
′
(
x
0
)
=
1
2
h
[
−
3
f
(
x
0
)
+
4
f
(
x
1
)
−
f
(
x
2
)
]
P
2
′
(
x
1
)
=
1
2
h
[
−
f
(
x
0
)
+
f
(
x
2
)
]
P
2
′
(
x
2
)
=
1
2
h
[
f
(
x
0
)
−
4
f
(
x
1
)
+
3
f
(
x
2
)
]
P'_2(x_0)=\frac{1}{2h} [-3f(x_0)+4f(x_1)-f(x_2)] \\ P'_2(x_1)=\frac{1}{2h} [-f(x_0)+f(x_2)]\\ P'_2(x_2)=\frac{1}{2h} [f(x_0)-4f(x_1)+3f(x_2)]\\
P2′(x0)=2h1[−3f(x0)+4f(x1)−f(x2)]P2′(x1)=2h1[−f(x0)+f(x2)]P2′(x2)=2h1[f(x0)−4f(x1)+3f(x2)]
对应三点的余项分别为
f
′
(
x
0
)
=
1
2
h
[
−
3
f
(
x
0
)
+
4
f
(
x
1
)
−
f
(
x
2
)
]
+
h
2
3
f
′
′
′
(
ξ
)
f
′
(
x
)
=
1
2
[
−
f
(
x
0
)
+
f
(
x
2
)
]
−
h
2
6
f
′
′
′
(
ξ
)
f
′
(
x
)
=
1
2
[
f
(
x
0
)
−
4
f
(
x
1
)
+
3
f
(
x
2
)
]
+
h
2
3
f
′
′
′
(
ξ
)
f'(x_0)=\frac{1}{2h}[-3f(x_0)+4f(x_1)-f(x_2)]+\frac{h^2}{3}f'''(\xi) \\ f'(x)=\frac{1}{2}[-f(x_0)+f(x_2)]-\frac{h^2}{6}f'''(\xi) \\ f'(x)=\frac{1}{2}[f(x_0)-4f(x_1)+3f(x_2)]+\frac{h^2}{3}f'''(\xi) \\
f′(x0)=2h1[−3f(x0)+4f(x1)−f(x2)]+3h2f′′′(ξ)f′(x)=21[−f(x0)+f(x2)]−6h2f′′′(ξ)f′(x)=21[f(x0)−4f(x1)+3f(x2)]+3h2f′′′(ξ)
可在插值积分多项式的基础上建立高阶数值微分公式
f
(
k
)
≈
P
m
(
k
)
(
x
)
f^{(k)}\approx{}P_m^{(k)}(x)
f(k)≈Pm(k)(x)
已知
P
2
P_2
P2关于变量t的表达式可对t求导得二阶数值微分公式
注:
- 使用两点公式的时候应该选择步长较小的函数值
- 一般情况下两点公式的精确性没有三点公式精确