【数值分析】高斯型求积公式,任意区间三点gauss求积公式,matlab实现
你哥同学 2024-07-02 11:35:03 阅读 53
Gauss型求积公式
Gauss型求积公式定义
∫
a
b
ρ
(
x
)
f
(
x
)
d
x
≈
∑
i
=
1
n
A
i
f
(
x
i
)
\int_{ a }^{b} \rho(x)f(x) \mathrm dx \approx \sum_{i=1}^{ n}A_if(x_i)
∫abρ(x)f(x)dx≈i=1∑nAif(xi)
如果求积公式具有
2
n
−
1
{2n-1}
2n−1 次代数精度,则称对应的节点
x
1
,
x
2
,
⋯
,
x
n
{x_1,x_2,\cdots ,x_n}
x1,x2,⋯,xn 为Gauss点,此时求积公式称为Gauss型求积公式。
为了讨论方便,本节取
n
{n}
n 个节点,并记节点为
x
1
,
x
2
,
⋯
x
n
{x_1,x_2,\cdots x_n}
x1,x2,⋯xn ,从
1
{1}
1 开始取!同时,所讨论的积分均为带有权函数
ρ
(
x
)
{\rho(x)}
ρ(x) 的积分。
插值型求积公式
∫
a
b
ρ
(
x
)
f
(
x
)
d
x
≈
∑
i
=
1
n
A
i
f
(
x
i
)
\int_{ a }^{b} \rho(x)f(x) \mathrm dx \approx \sum_{i=1}^{ n}A_if(x_i)
∫abρ(x)f(x)dx≈i=1∑nAif(xi)
的代数精度最高不超过
2
n
−
1
{2n-1}
2n−1 ,且达到
2
n
−
1
{2n-1}
2n−1 时,所有求积系数为正。所以高斯公式是给定节点数下代数精度最高的求积公式。
构造Gauss型求积公式的步骤
对给定区间
[
a
,
b
]
{[a,b]}
[a,b] 及权函数
ρ
(
x
)
{\rho(x)}
ρ(x) ,由Schmidt正交化过程构造正交多项式
P
0
(
x
)
,
P
1
(
x
)
,
⋯
,
P
n
(
x
)
{P_0(x),P_1(x), \cdots ,P_n(x)}
P0(x),P1(x),⋯,Pn(x)求出
P
n
(
x
)
{P_n(x)}
Pn(x) 的
n
{n}
n 个零点
x
1
,
x
2
,
⋯
,
x
n
{ x_1,x_2, \cdots, x_n}
x1,x2,⋯,xn 即为Gauss点计算求积系数
A
i
=
∫
a
b
ρ
(
x
)
l
i
(
x
)
d
x
{A_i= \int_{ a }^{b} \rho(x)l_i(x) \mathrm dx }
Ai=∫abρ(x)li(x)dx,
i
=
1
,
2
,
⋯
,
n
{ i=1,2,\cdots,n }
i=1,2,⋯,n,
l
i
{l_i}
li 为拉格朗日基函数
[!example]-
求计算积分
∫
−
1
1
x
2
f
(
x
)
d
x
{ \int_{ -1 }^{1} x^2f(x) \mathrm dx}
∫−11x2f(x)dx 的两点Gauss公式。
解:
ρ
(
x
)
=
x
2
,
n
=
2
\rho(x)=x^2, n=2
ρ(x)=x2,n=2
首先按Schmidt正交化求出正交多项式:
P
0
(
x
)
=
1
P
1
(
x
)
=
x
−
(
x
,
P
0
(
x
)
)
(
P
0
(
x
)
,
P
0
(
x
)
)
P
0
(
x
)
=
x
−
∫
−
1
1
x
3
d
x
∫
−
1
1
x
2
d
x
=
x
P
2
(
x
)
=
x
2
−
(
x
2
,
P
0
(
x
)
)
(
P
0
(
x
)
,
P
0
(
x
)
)
P
0
(
x
)
−
(
x
2
,
P
1
(
x
)
)
(
P
1
(
x
)
,
P
1
(
x
)
)
P
1
(
x
)
=
x
2
−
∫
−
1
1
x
4
d
x
∫
−
1
1
x
2
d
x
−
∫
−
1
1
x
5
d
x
∫
−
1
1
x
4
d
x
x
=
x
2
−
3
5
\begin{align*} P_0(x)=&1 \\ \\ P_1(x)=&x- \frac{(x,P_0(x))}{(P_0(x),P_0(x))}P_0(x)=x- \frac{\int_{ -1 }^{1} x^3 \mathrm dx}{\int_{ -1 }^{1} x^2 \mathrm dx}=x \\ \\ P_2(x)=&x^2- \frac{(x^2,P_0(x))}{(P_0(x),P_0(x))}P_0(x)-\frac{(x^2,P_1(x))}{(P_1(x),P_1(x))}P_1(x) \\ \\ =&x^2- \frac{\int_{ -1 }^{1} x^4 \mathrm dx}{\int_{ -1 }^{1} x^2 \mathrm dx}- \frac{\int_{ -1 }^{1} x^5 \mathrm dx}{\int_{ -1 }^{1} x^4 \mathrm dx}x=x^2- \frac{3}{5} \end{align*}
P0(x)=P1(x)=P2(x)==1x−(P0(x),P0(x))(x,P0(x))P0(x)=x−∫−11x2dx∫−11x3dx=xx2−(P0(x),P0(x))(x2,P0(x))P0(x)−(P1(x),P1(x))(x2,P1(x))P1(x)x2−∫−11x2dx∫−11x4dx−∫−11x4dx∫−11x5dxx=x2−53
再令
P
2
(
x
)
=
0
{P_2(x)=0}
P2(x)=0 求出Gauss点:
x
1
=
−
3
5
,
x
2
=
3
5
x_1=-\sqrt{\frac{3}{5}} \,\,,\,\, x_2=\sqrt{\frac{3}{5}}
x1=−53
,x2=53
最后计算求积系数:
A
1
=
∫
−
1
1
x
2
l
1
(
x
)
d
x
=
∫
−
1
1
x
2
x
−
x
2
x
1
−
x
2
d
x
=
1
3
A
2
=
∫
−
1
1
x
2
l
2
(
X
)
d
x
=
∫
−
1
1
x
2
x
−
x
1
x
1
−
x
2
d
x
=
1
3
\begin{align*} A_1=& \int_{ -1 }^{1} x^2l_1(x) \mathrm dx= \int_{ -1 }^{1} x^2 \frac{x-x_2}{x_1-x_2} \mathrm dx= \frac{1}{3}\\ \\ A_2=& \int_{ -1 }^{1} x^2l_2(X) \mathrm dx= \int_{ -1 }^{1} x^2 \frac{x-x_1}{x_1-x_2} \mathrm dx= \frac{1}{3} \end{align*}
A1=A2=∫−11x2l1(x)dx=∫−11x2x1−x2x−x2dx=31∫−11x2l2(X)dx=∫−11x2x1−x2x−x1dx=31
所以两点Gauss公式为
∫
−
1
1
x
2
f
(
x
)
d
x
≈
1
3
[
f
(
−
3
5
)
+
f
(
3
5
)
]
\int_{ -1 }^{1} x^2f(x) \mathrm dx \approx \frac{1}{3}[f(- \sqrt{\frac{3}{5}})+f(\sqrt{\frac{3}{5}})]
∫−11x2f(x)dx≈31[f(−53
)+f(53
)]
4.1 Gauss-legendre勒让德求积公式
区间
[
−
1
,
1
]
{[-1,1]}
[−1,1] ,权函数
ρ
(
x
)
=
1
{\rho(x)=1}
ρ(x)=1 。
Gauss求积公式中在区间
[
−
1
,
1
]
{[-1,1]}
[−1,1] 上权函数为
1
{1}
1 的Gauss点其实是确定的,所以可以通过把待求区间伸缩变换到
[
−
1
,
1
]
{[-1,1]}
[−1,1] 间,再通过Gauss求积公式求解。
一点Gauss求积公式
∫
−
1
1
f
(
x
)
d
x
≈
2
f
(
0
)
\int_{ -1 }^{1} f(x) \mathrm dx \approx 2f(0)
∫−11f(x)dx≈2f(0)
两点Gauss求积公式
∫
−
1
1
f
(
x
)
d
x
≈
f
(
−
1
3
)
+
f
(
1
3
)
\int_{ -1 }^{1} f(x) \mathrm dx \approx f(- \frac{1}{\sqrt{3}})+f(\frac{1}{\sqrt{3}})
∫−11f(x)dx≈f(−3
1)+f(3
1)
三点Gauss求积公式
∫
−
1
1
f
(
x
)
d
x
≈
5
9
f
(
−
3
5
)
+
8
9
f
(
0
)
+
5
9
f
(
3
5
)
\int_{ -1 }^{1} f(x) \mathrm dx \approx \frac{5}{9} f(- \sqrt{\frac{3}{5}})+ \frac{8}{9}f(0) +\frac{5}{9}f(\sqrt{\frac{3}{5}})
∫−11f(x)dx≈95f(−53
)+98f(0)+95f(53
)
任意区间的三点高斯求积公式:
∫
a
b
f
(
x
)
d
x
≈
(
b
−
a
)
2
(
5
9
f
(
(
a
+
b
)
2
−
3
5
(
b
−
a
)
2
)
+
8
9
f
(
(
a
+
b
)
2
)
+
5
9
f
(
(
a
+
b
)
2
+
3
5
(
b
−
a
)
2
)
)
\int_{ a }^{b} f(x) \mathrm dx \approx \frac{(b-a)}{2} \bigg( \frac{5}{9} f(\frac{(a+b)}{2} - \sqrt{\frac{3}{5}}\frac{(b-a)}{2})+ \frac{8}{9}f(\frac{(a+b)}{2}) +\frac{5}{9}f(\frac{(a+b)}{2} + \sqrt{\frac{3}{5}}\frac{(b-a)}{2}) \bigg)
∫abf(x)dx≈2(b−a)(95f(2(a+b)−53
2(b−a))+98f(2(a+b))+95f(2(a+b)+53
2(b−a)))
matlab实现
%% 三点高斯勒让德求积公式
% 输入函数,积分上界,积分下界
function I = gaussL3P(f,a,b)
x = (a+b)/2+(b-a)/2*[-sqrt(0.6) 0 sqrt(0.6)];
I = (b-a)/2*f(x)*[5 8 5]'/9;
end
四点Gauss求积公式
高斯点
:
±
0.8611363
,
高斯系数
:
0.3478548
高斯点:\pm0.8611363 \,\,,\,\, 高斯系数:0.3478548
高斯点:±0.8611363,高斯系数:0.3478548
高斯点
:
±
0.3399810
,
高斯系数
:
0.6521452
高斯点:\pm0.3399810 \,\,,\,\, 高斯系数:0.6521452
高斯点:±0.3399810,高斯系数:0.6521452
[!example]-
已知三点Gauss公式
∫
−
1
1
f
(
x
)
d
x
≈
5
9
f
(
0.6
)
+
8
9
f
(
0
)
+
5
9
f
(
−
0.6
)
\int_{ -1 }^{1} f(x) \mathrm dx \approx \frac{5}{9}f(\sqrt{0.6})+ \frac{8}{9}f(0)+ \frac{5}{9}f(- \sqrt{0.6})
∫−11f(x)dx≈95f(0.6
)+98f(0)+95f(−0.6
)
试用如上公式计算
∫
0.5
1
x
d
x
{ \int_{ 0.5 }^{1} \sqrt{x} \mathrm dx}
∫0.51x
dx 的值。
解:换限,设
t
=
a
x
+
b
t=ax+b
t=ax+b
{
0.5
a
+
b
=
−
1
a
+
b
=
1
⇒
{
a
=
4
b
=
−
3
\begin{cases} 0.5a+b=-1 \\ \\ a+b=1 \end{cases} \Rightarrow \begin{cases} a=4 \\ \\ b=-3 \end{cases}
⎩
⎨
⎧0.5a+b=−1a+b=1⇒⎩
⎨
⎧a=4b=−3
∫
0.5
1
f
(
x
)
d
x
=
1
4
∫
−
1
1
f
(
t
+
3
4
)
d
t
≈
1
4
[
5
9
f
(
0.6
+
3
4
)
+
8
9
f
(
0
+
3
4
)
+
5
9
f
(
−
0.6
+
3
4
)
]
=
1
4
[
5
9
0.6
+
3
4
+
8
9
0
+
3
4
+
5
9
−
0.6
+
3
4
]
=
0.4310
\begin{align*} \int_{ 0.5 }^{1} f(x) \mathrm dx=& \frac{1}{4} \int_{ -1 }^{1} f(\frac{t+3}{4}) \mathrm dt \\ \\ \approx& \frac{1}{4}[\frac{5}{9} f(\frac{\sqrt{0.6}+3}{4})+ \frac{8}{9}f(\frac{0+3}{4}) + \frac{5}{9}f(\frac{- \sqrt{0.6}+3}{4}) ] \\ \\ =&\frac{1}{4}[\frac{5}{9} \sqrt{\frac{\sqrt{0.6}+3}{4}}+ \frac{8}{9}\sqrt{\frac{0+3}{4}} + \frac{5}{9}\sqrt{\frac{- \sqrt{0.6}+3}{4}} ] \\ \\ =&0.4310 \end{align*}
∫0.51f(x)dx=≈==41∫−11f(4t+3)dt41[95f(40.6
+3)+98f(40+3)+95f(4−0.6
+3)]41[9540.6
+3
+9840+3
+954−0.6
+3
]0.4310
4.2 变形的Gauss型求积公式
Gauss-Laguerre拉盖尔求积公式
区间
[
0
,
+
∞
)
{[0,+\infty)}
[0,+∞) ,权函数
ρ
(
x
)
=
e
−
x
{\rho(x)=e^{-x}}
ρ(x)=e−x
Gauss-Hermite求积公式
两点Gauss-Hermite求积公式
∫
−
∞
+
∞
e
−
x
2
f
(
x
)
d
x
≈
π
2
f
(
−
2
2
)
+
π
2
f
(
2
2
)
\int_{ -\infty }^{+\infty} e^{-x^2}f(x) \mathrm dx \approx \frac{\sqrt{\pi}}{2} f(-\frac{\sqrt{2}}{2})+\frac{\sqrt{\pi}}{2} f(\frac{\sqrt{2}}{2})
∫−∞+∞e−x2f(x)dx≈2π
f(−22
)+2π
f(22
)
三点Gauss-Hermite求积公式
∫
−
∞
+
∞
e
−
x
2
f
(
x
)
d
x
≈
π
6
f
(
−
6
2
)
+
2
π
3
f
(
0
)
+
π
6
f
(
6
2
)
\int_{ -\infty }^{+\infty} e^{-x^2}f(x) \mathrm dx \approx \frac{\sqrt{\pi}}{6} f(-\frac{\sqrt{6}}{2})+ \frac{2 \sqrt{\pi}}{3} f(0)+\frac{\sqrt{\pi}}{6} f(\frac{\sqrt{6}}{2})
∫−∞+∞e−x2f(x)dx≈6π
f(−26
)+32π
f(0)+6π
f(26
)
区间
(
−
∞
,
+
∞
)
{(-\infty,+\infty)}
(−∞,+∞) ,权函数
ρ
(
x
)
=
e
−
x
2
{\rho(x)=e^{-x^2}}
ρ(x)=e−x2
上一篇: 【C++】万字详解IO流(输入输出流+文件流+字符串流)
本文标签
声明
本文内容仅代表作者观点,或转载于其他网站,本站不以此文作为商业用途
如有涉及侵权,请联系本站进行删除
转载本站原创文章,请注明来源及作者。