算符
设给定位移 $h$，定义位移算子 $E$：
$$
E^sp(x)=p(x+hs)
$$
可以用位移算子定义差分算子：
- 一阶向前差分算子：$\Delta_+=\Delta=E-I$
- 一阶向后差分算子：$\Delta_-=I-E^{-1}$
- 一步中心差分算子：$\Delta_0=E-E^{-1}$
- 半步中心差分算子：$\delta=E^{1/2}-E^{-1/2}$
- 二阶中心差分算子：$\delta^2=E-2I+E^{-1}$
- 平均算子：$\mu=\frac{1}{2}\qty(E^{1/2}+E^{-1/2})$
- 微分算符：$\mathop{D}=\dv{x}$
- 积分算符：$\int_{x}^{x+h} \qty(\cdot)(t)\dd{t}$
如果间隔 $h$ 需要被指明，会使用下标的形式进行表达，即 $E_h,E_\tau$
实际上，算符 $\Delta,\nabla, \delta, \mu, E, D, J$ 满足交换律、分配律、结合律。特殊的，有如下关系成立：
$$
DJ=JD=\Delta
$$
对于算子的可逆性，有以下结论：
$E$ aninhilated nothing;
$\Delta,\nabla,\delta$ annihilate any function of period $h$;
$\delta$ annihilate function of period h;
$J$ annihilates any function(with constant C);
$D$ annihilates any constant;
$\mu$ annihilates odd-harmonic function $f(x)$ for which $f(x+h)=-f(x)$.
算符的性质
我们可以使用算符 $E$ 来统一表示以上算符：
$$
\begin{aligned}
& \Delta =E-I \\
& \nabla =I-E^{-1} \\
& \delta=E^{1/2}-E^{-1/2} \\
&\mu=\frac{1}{2}\qty(E^{1/2}+E^{-1/2}) \\
& DJ=JD=\Delta=E-I
\end{aligned}
$$
若 $s$ 为非负整数，则：
$$
\begin{aligned}
\Delta^s
=&E^s\nabla^s \\
=&E^{s/2}\delta^{s} \\
=&(E-I)^{s} \\
=&\sum_{k=0}^{s}\binom{s}{k}(-1)^{s-k}E^{k}
\end{aligned}
$$
于是有：
$$
\begin{aligned}
&\Delta^s=\sum_{k=0}^{s}\binom{s}{k}(-1)^{s-k}E^{k} \\
&\nabla^s=\sum_{k=0}^{s}\binom{s}{k}(-1)^{s-k}E^{k-s} \\
&\delta^s=\sum_{k=0}^{s}\binom{s}{k}(-1)^{s-k}E^{k-s/2}
\end{aligned}
$$
进一步的，可以得到：
$$
\begin{aligned}
& E=I+\Delta \\
& E(1-\nabla)=I \\
& \qty(E^{1/2}-\frac{1}{2}\delta)^2-\frac{1}{4}\delta^2=I \\
& E^{1/2}-\frac{1}{2}\delta-\mu=0
\end{aligned}
$$
可以反解出：
$$
\begin{aligned}
& E=I+\Delta \\
& E=\frac{1}{1-\nabla} \\
& E^{1/2}=\qty(1+\frac{1}{4}\delta^2)^{1/2}+\frac{1}{2}\delta \\
& \mu=\qty(1+\frac{1}{4}\delta^2)^{1/2}
\end{aligned}
$$
其中，蕴含着互逆关系：$(I+\Delta)(I-\nabla)=I$。
$$
\begin{aligned}
&\Delta^s=(E-I)^{s} \\
=&E^s\nabla^s=E^{s/2}\delta^{s} \\
=&\sum_{k=0}^{s}\binom{s}{k}(-1)^{s-k}E^{k}
\end{aligned}
$$
所以：
$$
\begin{aligned}
&\Delta^s=\sum_{k=0}^{s}\binom{s}{k}(-1)^{s-k}E^{k} \\
&\nabla^s=\sum_{k=0}^{s}\binom{s}{k}(-1)^{s-k}E^{k-s} \\
&\delta^s=\sum_{k=0}^{s}\binom{s}{k}(-1)^{s-k}E^{k-s/2}
\end{aligned}
$$
所以：
$$
\begin{aligned}
& E=I+\Delta \\
& E(1-\nabla)=I \\
& \qty(E^{1/2}-\frac{1}{2}\delta)^2-\frac{1}{4}\delta^2=I \\
& E^{1/2}-\frac{1}{2}\delta-\mu=0
\end{aligned}
$$
所以：
$$
\begin{aligned}
& E=I+\Delta \\
& E=\frac{1}{1-\nabla} \\
& E^{1/2}=\qty(1+\frac{1}{4}\delta^2)^{1/2}+\frac{1}{2}\delta \\
& \mu=\qty(1+\frac{1}{4}\delta^2)^{1/2}
\end{aligned}
$$
所以：
$$
(I+\delta)(I-\nabla)=I
$$
考虑一个 $n$ 阶多项式 $p_n(x)$，其差分 $n+1$ 次的结果为 $0$，即：
$$
\nabla^{n+1}p_n(x)\equiv 0
$$
所以在宗量为 $n$ 阶多项式意义下，有：
$$
\begin{aligned}
&(I-\nabla)\sum_{k=0}^{n}\nabla^k=I-\nabla^{n+1}\triangleq I \\
\Rightarrow
& (I-\nabla)^{-1}=\sum_{k=0}^{n}\nabla^{k}
\end{aligned}
$$
上式成立时，要求多项式阶数不超过 $n$。
若 $n \to \infty$，则可以得到更一般的情况：
$$
E=(I-\nabla)^{-1}=\sum_{k=0}^{\infty} \nabla^{k}
$$
对于一般的光滑函数 $f(x)$，由牛顿插值可以得到 $f(x)$ 在 $x=x_0$ 处的展开：
$$
\begin{aligned}
f(x_0+sh)
=&\sum_{k=0}^{\infty}\binom{s}{k}\Delta^{k}f(x_0), & \qqtext{forward} \\
=&\sum_{k=0}^{\infty}(-1)^{k}{-s \choose k}\nabla^{k}f(x_0), & \qqtext{backward}
\end{aligned}
$$
对于微分算符，容易发现：
$$
\begin{aligned}
&Ep(x)=p(x+h)=\sum_{k = 0}^{\infty} \frac{h^k}{k!}D^kp(x)=e^{hD}p(x) \\
\Rightarrow &\boxed{E=e^{hD}} \\
\Rightarrow &
hD=\ln E
=\ln(I+\Delta)
=-\ln(I-\nabla) \\
& \quad
=2\ln(\qty(1+\frac{1}{4}\delta^2)^{1/2}+\frac{1}{2}\delta) \\
& \quad
=2\sinh^{-1}\qty(\frac{\delta}{2})
\\
& \qty(\text{note:}\sinh^{-1}(x)=\ln(x+\sqrt{x^2+1}))
\end{aligned}
$$
也可以用 $hD$ 来表示 $\delta,\mu$：
$$
\begin{aligned}
& \delta=2\sinh \frac{hD}{2} \quad
& \mu=\cosh \frac{hD}{2} \quad
& \mu\delta=\sinh hD
\end{aligned}
$$
如果采用 $\ln(1+z)$ 的有理逼近技术 $\frac{P_m(x)}{Q_n(x)}$，可得：
$$
D=\frac{1}{h}\frac{\Delta}{1+\frac{1}{2}\Delta}+O(h^2)
$$
斯特林数
对于一阶导数的生成函数形式，有：
$$
D=\frac{\ln(I+\Delta)}{h}
$$
对于 $n$ 阶导数的生成函数形式，有：
$$
D^{n}=\frac{\qty(\ln(I+\Delta))^n}{h^n}
$$
对于 $n$ 阶后向差分的生成函数形式（注意到 $I+\Delta=E=e^{hD}$），有：$I+\Delta=E=e^{hD}$，即:
$$
\Delta^{n}=\qty(e^{hD}-I)^{n}
$$
对于第一类斯特林数，有：
$$
\frac{\qty(\ln(1+x))^{n}}{n!}
=\sum_{k=n}^{\infty}\frac{1}{k!}\begin{bmatrix}k \\ n\end{bmatrix} x^k
$$
对于第二类斯特林数，有：
$$
\frac{\qty(e^{x}-1)^{n}}{n!}
=\sum_{k=n}^{\infty}\frac{1}{k!}\begin{Bmatrix}k \\ n\end{Bmatrix} x^k
$$
$(\ln(1+x))^m$
考虑到：
$$
\begin{aligned}
&\qty(\ln(1+x))^m=\qty(\sum_{n=1}^{\infty}\frac{(-1)^{n-1}}{n}x^n)^m \\
=&\sum_{n=m}^{\infty} \qty(
\sum_{\sum_{p=1}^{m}n_p=n}\prod_{p=1}^{m} \frac{(-1)^{n_p-1}}{n_p}
) x^n \\
=&\sum_{n=m}^{\infty} \qty(
\sum_{\sum_{p=1}^{m}n_p=n}\prod_{p=1}^{m} \frac{1}{n_p}
) (-1)^{n-m}x^n
\end{aligned}
$$
将 $n$ 个不同小球排列成 $m$ 个圆排列，其方案数为：
$$
\begin{bmatrix}n \\ m\end{bmatrix}
=\sum_{\sum_{p=1}^{m}n_p=n}\frac{1}{m!}\qty(n!\prod_{p=1}^{m}\frac{1}{n_p})
$$
定义有符号斯特林数：
$$
\begin{bmatrix}n \\ m\end{bmatrix}{s}
= (-1)^{n+m}\begin{bmatrix}n \\ m\end{bmatrix}
$$
则系数展开问题的解为：
$$
\frac{\qty(\ln(1+x))^m}{m!}=\sum{n=m}^{\infty}\frac{1}{n!}\begin{bmatrix}n \\ m\end{bmatrix}_{s}x^{n}
$$
$(e^x-1)^m$
考虑到：
$$
(e^{x}-1)^m
=\qty(\sum_{n=1}^{\infty} \frac{x^n}{n!})^m
=\sum_{n=m}^{\infty}\qty(\sum_{\sum_{p=1}^{m}n_p=n}\prod_{p=1}^{m}\frac{1}{n_p!})x^{n}
$$
观察 $x^{n}$ 的系数的意义，相当于从 $m$ 种无限多物品中，每种物品至少选取 $1$ 个，选取 $n_p$ 个物品对贡献系数 $\frac{1}{n_p!}$，这个系数意味着每个种类的物品本质相同。
如果把每一种物品看作盒子，物品总个数 $n=\sum n_p$ 看作小球，那么将 $n$ 个不同的小球装进 $m$ 个相同的盒子（盒子不为空）里的方案数为：
$$
\begin{Bmatrix}n \\ m\end{Bmatrix}
=\sum_{\sum_{p=1}^{m}n_p=n}\frac{1}{m!}\qty(n! \prod_{p=1}^{m}\frac{1}{n_p!})
$$
这个方案数成为第二类斯特林数，所以系数展开问题的解为：
$$
\begin{aligned}
&\qty(e^{x}-1)^{m}=\sum_{n=m}^{\infty}\frac{m!}{n!}\begin{Bmatrix}n \\ m\end{Bmatrix} x^{n} \\
\Rightarrow
&\frac{\qty(e^{x}-1)^{m}}{m!}=\sum_{n=m}^{\infty}\frac{1}{n!}\begin{Bmatrix}n \\ m\end{Bmatrix} x^{n}
\end{aligned}
$$
积分算符
考虑积分形式：
$$
J=\Delta D^{-1}
$$
则有如下分段积分：
$$
\int_{x_0}^{x_0+sh}f(x)\dd{x}=\qty(\sum_{k=0}^{s-1}E^{k})Jf(x)
=\frac{E^{s}-I}{E-I}Jf(x)
$$
考虑到 $E=I+\Delta, hD=\ln E,J=\Delta D^{-1}$，于是有：
$$
\begin{aligned}
\frac{E^{s}-I}{E-I}J
=&h\frac{(I+\Delta)^{s}-I}{\Delta}\frac{\Delta}{\ln(I+\Delta)} \\
=&h\qty(\sum_{i=0}^{\infty}{r \choose i+1}\Delta^i)
\sum_{j=0}^{\infty} c_j\Delta^{j} \\
=&h\sum_{i=0}^{\infty}\sum_{j=0}^{\infty}hc_j{r\choose i+1}\Delta^{i+j} \\
=&\sum_{k=0}^{\infty}h
\underbrace{\sum_{i=0}^{k}c_{k-i}{s \choose i+1}}_{\alpha _ {k}^{(r)}}
\Delta^{k}
\end{aligned}
$$
可以得到积分的差分形式：
$$
\int_{x_0}^{x_0+sh}f(x)\dd{x}=h\qty(\sum_{k=0}^{\infty}\alpha_{k}^{(s)}\Delta^{k})f(x)
$$
对于中心差分格式，类似的，有：
$$
\begin{aligned}
\int_{x_0-sh}^{x_0+sh}f(x)\dd{x}
=&\sum_{k=-m}^{m-1}E^{k}Jf(x) \\
=&\frac{E^{s}-E^{-s}}{E-I}Jf(x) \\
=&\frac{e^{shD}-e^{-shD}}{\Delta}\Delta D^{-1}f(x) \\
=&\frac{2\sinh(mhD)}{D}f(x)
\end{aligned}
$$
可以用波矢为 $k$ 的行波 $e^{ikx}$ 检验差分的相容阶。
$$
\begin{aligned}
&\Delta e^{ikx} \\
=&e^{ik(x+h)}-e^{ikx} \\
=&e^{ikx}\qty(e^{ikh}-1) \\
=&e^{ikx}\qty(ikh+O(k^2h^2)) \\
&=hDe^{ikx}+O(h^2) \\
&\delta^2 e^{ikx} \\
=&e^{ikx}\qty(e^{ikh}-1+e^{-ikh}-1) \\
=&e^{ikx}\qty(ikh-k^2h^2/2-ikh-k^2h^2/2+O(k^4h^4)) \\
=&e^{ikx}\qty(-(kh)^2+O((kh)^4)) \\
=&h^2D^2e^{ikx}+O(h^4)
\end{aligned}
$$
$$
\begin{aligned}
& J(\cdot)=\int_{x_0}^{x_0+h}(\cdot)\dd{x} \\
\Rightarrow
& JD=DJ=\Delta \\
\Rightarrow
& J=\Delta D^{-1} \\
\Rightarrow
& \int_{x_0}^{x_0+rh}(\cdot)\dd{x}
=\sum_{s=0}^{r-1}E^{s}J
=\frac{E^{r}-1}{E-1}J
=h\frac{(1+\Delta)^r-1}{\Delta}\frac{\Delta}{\ln(1+\Delta)}
\end{aligned}
$$
第一项，有：
$$
\frac{(1+\Delta)^r-1}{\Delta}
=\sum_{k=1}^{\infty}{r \choose k}\Delta^{k-1}
=\sum_{k=0}^{\infty}{r \choose k+1}\Delta^{k}
$$
第二项，有：
$$
\begin{aligned}
\frac{\Delta}{\ln(1+\Delta)}
&=\frac{\Delta}{\sum_{n=1}^{\infty}\frac{(-1)^{n-1}\Delta^{n}}{n}} \\
&=\frac{1}{1-
\underbrace{\sum_{n=1}^{\infty}\frac{(-1)^{n+1}\Delta^{n}}{n+1}} _ {p(\Delta)}} \\
&=1+\sum_{k=1}^{\infty}p^k(\Delta)
\end{aligned}
$$
其中，$p(\Delta)$ 为：
$$
p(\Delta)=\frac{\Delta}{2}-\frac{\Delta^2}{3}+\frac{\Delta^3}{4}+\cdots
$$
所以：
$$
\begin{aligned}
\int_{x_0}^{x_0+rh}(\cdot)\dd{x}
=&\sum_{s=0}^{r-1}E^{s}J \\
=&\frac{E^{r}-1}{E-1}J \\
=&h\frac{(1+\Delta)^r-1}{\Delta}\frac{\Delta}{\ln(1+\Delta)} \\
=&h\sum_{i=0}^{\infty}{r \choose i+1}\Delta^{i}
\sum_{j=0}^{\infty}
\qty(\sum_{k=0}^{\infty}\frac{(-1)^{k}\Delta^k}{k+1})^j
\end{aligned}
$$
另一方面：
$$
\begin{aligned}
J&=\frac{\Delta}{D}=\frac{e^{hD}-1}{D} \\
&=\sum_{n=1}^{\infty}\frac{h^n}{n!}D^{n-1} \\
&=\sum_{n=1}^{\infty}\frac{h}{n!}(\ln(1+\Delta))^{n-1}
\end{aligned}
$$
伯努利数的生成函数形式为：
$$
\frac{z}{e^{z}-1}=\sum_{n=0}^{\infty}\frac{B_n}{n!}z^n
$$
所以：
$$
\begin{aligned}
& DJ\Delta^{-1}=1 \\
\Rightarrow
& \frac{DJ}{E-1}=1 \\
\Rightarrow
& \frac{hD}{e^{hD}-1}J=h \\
\Rightarrow
& h=\sum_{n=0}^{\infty}\frac{B_n}{n!}h^n D^nJ \\
\Rightarrow
& h=J+\sum_{n=1}^{\infty}\frac{B_nh^n}{n!}D^{n-1}\Delta \\
\Rightarrow
& J=h-\sum_{n=1}^{\infty} \frac{B_nh}{n!}(\ln(1+\Delta))^{n-1}\Delta
\end{aligned}
$$
