数值分析_第四章_数值积分与数值微分

\qquad实际问题当中常常计算积分，有些数值方法如微分方程和积分方程的求解，也都和积分计算相联系。使用牛顿-莱布尼兹($Newton-Leibniz$)公式求解往往比较困难，有时原函数不能用初等函数表示或者求解过程十分困难，因此有必要研究积分的数值计算问题。
\qquad由积分中值定理，在积分区间$[a,b]$中存在一点$\xi$使得
$${\int }_b^{a}f(x)dx=(b-a)f(\xi )$$\qquad但点$\xi$的具体位置一般难以确定，称$f(\xi )$为区间$[a,b]$上的平均高度，只要对平均高度$f(\xi )$提出一种算法，相应的获得一种数值积分方法。
\qquad如果将两端的高度$f(a)、f(b)$的算术平均值作为平均高度$f(\xi )$的近似值，则求积公式为$${\int }_a^{b}f(x)dx\approx \frac{b-a}{2}[f(a)+f(b)]$$称为梯形公式(n=1)。若使用区间中间点$c=\frac{a+b}{2}$作为$\xi$，则可得公式为
$${\int }_a^{b}f(x)dx\approx (b-a)f(\frac{a+b}{2})$$称为中矩形公式或矩形公式。
\qquad设将积分区间$[a,b]$划分为$n$等份，步长$h=\frac{b-a}{n}$，选去等距节点$x_k=a+kh$构造的插值型求积公式为$!!!$
$$I_n=(b-a)\sum _{k=0}^n{C}_k^{(n)}f(x_k)$$称为牛顿-柯特斯(Newton-Cotes)公式，式中$C_k^{(n)}$称为柯特斯系数，引进变换$x=a+th$，则有
$$C_k^{(n)}=\frac{(-1{)}^{n-k}}{nk!(n-k)!}{\int }_0^n\underset{j̸=k}{\overset{n}{\prod _{j=0}}}(t-j)dt$$\qquad辛普森(Simpson)公式(n=2)
$$I_n=(b-a)[\frac{1}{6}f(a)+\frac{4}{6}f(\frac{a+b}{2})+\frac{1}{6}f(b)]$$代数精确度为3，余项表示为
$$R[f]=K{f}^{(4)}(\eta ),\eta \in (a,b)$$\qquad更一般的，在区间$[a,b]$上适当选取某些节点$x_k$，然后用$f(x_k)$的加权平均的到平均高度的$f(\xi )$的近似值，构造的求积公式如下
$${\int }_a^{b}f(x)dx\approx \sum _{k=0}^n{A}_{k}f(x_k)$$式中$x_k$称为求积节点，$A_k$称为求积系数，亦称为伴随节点$x_k$的权。权$A_k$仅与节点$x_k$的选择有关，不依赖于被积函数$f(s)$的具体形式。
\qquad代数精度 定义1 $!!!$如果某个求积公式对于次数不超过$m$的多项式均能准确的成立，但对于$m+1$次多项式就不准确成立，则称该求积公式具有$m$次代数精度。欲使求积公式具有$m$次代数精度，只要令其对$f(x)=1,x,\cdots ,x^m$都能准确成立，即
$$\{\begin{array}{l}\sum A_k=b-a\\ \sum A_k{x}_k=\frac{1}{2}(b^2-a^2),\\ ⋮\\ \sum A_k{x}_k^m=\frac{1}{m+1}(b^{m+1}-a^{m+1})\end{array}$$\qquad求积公式的收敛性与稳定性 定义2在求积公式中，若
$$\underset{h\to 0}{\underset{n\to \infty }{\lim }}\sum _{k=0}^n{A}_{k}f(x_k)-{\int }_a^{b}f(x)dx$$其中h=max1≤i≤n{xi−xi−1}h=\mathop{max}\limits_{1\leq i\leq n}\{x_{i}-x_{i-1}\}h=1≤i≤nmax​{xi​−xi−1​}，则称求积公式收敛。在求积公式中，由于计算$f(x_k)$可能产生误差${\delta }_k$实际的得到的为$\widetilde{f_k}$，即$f(x_k)=\widetilde{f_k}+{\delta }_k$，记
$$I_n(f)=\sum _{k=0}^n{A}_{k}f(x_k)I_n(\tilde{f})=\sum _{k=0}^n{A}_k\widetilde{f_k}$$如果对任给小正数$\epsilon >0$，只要误差$|{\delta }_k|$充分小就有，
$$|I_n(f)-I_n(\tilde{f})|=|\sum _{k=0}^n{A}_k[f(x_k)-\widetilde{f_k}]|\le \epsilon$$$表明求积公式计算是稳定的$。
\qquad定义3 对任意给定的$\epsilon >0$，若$\exists \delta >0$，只要$|f(x_k)-\widetilde{f_k}|\le \delta (k=0,1,\cdots ,n)$就有上式成立，则称求积公式是稳定的。
\qquad求积公式的余项 若求积公式的代数精度为$m$，则由求积公式余项的表达式可将余项表达为
$$R[f]={\int }_a^{b}f(x)dx-\sum _{k=0}^n{A}_{k}f(x_k)=K{f}^{(m+1)}(\eta )$$\qquad$K$为不依赖$f(x)$的待定参数，$\eta \in (a,b)$。结果表明当$f(x)$是次数小于等于$m$的多项式时，由于$f^{(m+1)}(x)=0$，故此时$R[f]=0$，即求积公式精确成立。当$f(x)=x^{m+1}$时，$f^{(m+1)}(x)=(m+1)!$，又$R_n(f)̸=0$，可求得
$$\begin{gathered} K=\frac{1}{(m+1)!}[{\int }_a^b{x}^{m+1}dx-\sum _{k=0}^n{A}_k{x}_k^{m+1}] \\ =\frac{1}{(m+1)!}[\frac{1}{m+2}(b^{m+2}-a^{m+2})-\sum _{k=0}^n{A}_k{x}_k^{m+1}] \end{gathered}$$带入余项可以的到具体的余项表达式。梯形公式和中矩形公式的代数精度为1，其余项表达式为
$$R[f]=K{f}^{{}^{\prime \prime }}(\eta )$$

复合求积公式

\qquad为了提高精度通常将积分区域分成若干子区间(通常等分)，在每个子区间上用低阶求积公式，称为复合求积法。
\qquad复合梯形公式 将区间$[a,b]$分为$n$等份，分点$x_k=a+kh,h=\frac{b-a}{n},k=0,1,\cdots ,n$在每个子区间上$[x_k,x_{k+1}](k=0,1,\cdots ,n-1)$上采用梯形公式则得
$$I={\int }_a^{b}f(x)dx=\sum _{k=0}^{n-1}{\int }_{x_k}^{x_{k+1}}f(x)dx=\frac{h}{2}\sum _{k=0}^{n-1}[f(x_k+f(x_{k+1}))]+R_n(f)$$记
$$T_n=\frac{h}{2}\sum _{k=0}^{n-1}[f(x_k+f(x_{k+1}))]$$称为复合梯形公式。其余项表示为
$$R_n(f)=I-T_n=\sum _{k=0}^{n-1}[-\frac{h^3}{12}{f}^{{}^{\prime \prime }}({\eta }_k)],{\eta }_k\in (x_k,x_{k+1})$$复合辛普森求积公式与此相似，子区间上展开方式变为辛普森公式。
\qquad定理 设$f(x)\in C^{\infty }[a,b]$，则有
$$T(h)=I+{\alpha }_1{h}^2+{\alpha }_2{h}^4+\cdots +{\alpha }_l{h}^{2l}+\cdots$$其中系数${\alpha }_l(l=1,2,\cdots )$与$h$无关(由泰勒展开可得)。只要真值与近似值得误差能表示为$h$的幂级数，都能使用外推算法(将步长变为$\frac{1}{2}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)$$

自适应积分方法

\qquad在复合求积过程中，当被积函数变化剧烈时，使用较小步长，当被积函数变化平缓时，使用较大步长。