巴拿赫不动点定理
(Banach Fixed-Point Theorem)

不动点理论(Fixed-Point Theory)是数学分析、拓扑学、泛函中的重要研究课题，主要研究方程 $x = g(x)$ 的解的存在性、唯一性及计算方法。满足该方程，保持到自身映射的点 $x$ 就被称为函数 $g$ 的不动点。它不仅是纯数学的基石，还被广泛应用于博弈论(如纳什均衡)、计算机科学(如程序语义)、数据科学(如图像处理)等诸多领域。

本文将介绍的巴拿赫不动点定理(Banach Fixed-Point Theorem) 描述了在完备度量空间中的收缩映射，必然存在不动点并且唯一，所以又被称为收缩映射定理。从一个初始点 $x_0$ 出发，不断迭代 $x_{n+1} = g(x_n)$ 就可以收敛到不动点上。它为数值求解提供了理论依据。上一节中介绍的牛顿法，实际上就是这种不动点的迭代过程。

1. 不动点

对于函数 $g(x)$，如果有一个点 $p$ 满足 $p = g(p)$，那么我们称 $p$ 是函数 $g(x)$ 的不动点。巴拿赫不动点理论研究的是完备度量空间 $\mathcal{X}, d$ 下的自映射 $g: \mathcal{X} \rightarrow \mathcal{X}$。所谓的度量空间，是指该空间中的任何两个点之间都有距离定义。该定义必须满足以下条件，我们用符号 $d$ 来表示。

非负数: 距离不能为负数，如果两点的距离为 0 那么它们在空间上是重合的。
对称性: 任意两点 $A$ 到 $B$ 的距离，应当等于 $B$ 到 $A$ 的距离。
三角不等式: 空间中任意三个点 $A,B,C$，两两互联形成三条边，任意两条边之和大于第三条边。

如果一个度量空间中的每个柯西序列，都能收敛到空间中的某个点上，我们就说这个空间是完备的。所谓的柯西序列是指，一个数列 $\{x_n\}$ 对于任意的正数 $\varepsilon > 0$，都存在一个正整数 $N$，使得当 $m,n \> N$ 时，都有 $d(x_m, x_n) < \varepsilon$。我们常说的有理数集 $Q$ 就是不完备的，实数集 $R$ 是完备的。开区间 $(0, 1)$ 是不完备的，因为数列 $\{\frac{1}{n}\}$ 收敛到 0，但 0 就不在开区间 $(0, 1)$ 中。但是闭区间 $[0, 1]$ 就是完备的。

如果函数 $g(x)$ 在闭区间 $[a, b]$ 内有定义，并且其值域也在 $[a, b]$ 内，即 $g([a,b]) \subset [a, b]$，那么函数 $g(x)$ 在闭区间 $[a, b]$ 内一定存在一个不动点。如果 $g(a) = a$ 或者 $g(b) = b$，那么闭区间的两个端点上就是它的不动点。否则，因为 $g([a,b]) \subset [a, b]$，一定有 $g(a) > a, g(b) < b$。我们可以构造函数 $h(x) = g(x) - x$，有$h(a) = g(a) - a < 0$，$h(b) = g(b) - b > 0$。根据介值定理，一定有 $p \in (a,b)$，使得 $h(p) = 0 = g(p) - p$，即 $g(p) = p$。该点 $p$ 就是函数 $g(x)$ 的不动点。

需要注意的是，$g([a,b]) \subset [a, b]$ 只是不动点存在的一个充分非必要条件。比如 $g(x) = x^2$ 的不动点，有方程 $x = x^2 \Longrightarrow x^2 - x = x(x - 1) = 0$。显然在闭区间 $[0.5, 1.5]$ 上存在一个不动点 $x = 1$，但是 $g(1.5) = 2.25 \notin [0.5, 1.5]$。

对于满足 $g([a,b]) \subset [a, b]$ 的函数 $g(x)$，如果在 $(a, b)$ 上存在导数 $g'(x)$，对于一个正数 $k < 1$ 有 $|g'(x)| ≤ k < 1, \forall x \in (a, b)$ 那么不动点就是唯一的。假设 $p, q$ 都是 g(x) 在区间 $[a, b]$ 内的不动点，如果 $p \neq q$，不妨设 $p < q$，那么根据中值定理，一定存在一个数 $\xi \in (p, q) \subset [a, b]$，满足 $$ \frac{g(p) - g(q)}{p - q} = g'(\xi) $$ 结合不动点的定义，有 $|p -q| = |g(p) - g(q)| = |g'(\xi)||p - q| ≤ k | p - q | < |p - q|$。矛盾。所以 $p = q$。因此区间 $[a,b]$ 内的不动点唯一。需要注意的是，条件 $|g'(x)| ≤ k < 1, \forall x \in (a, b)$ 也只是唯一性的充分非必要条件。

上述命题的证明过程中，我们看到不动点 $x = g(x)$ 与我们想要求解的方程根 $f(x) = 0$ 之间是很大关联的。我们可以通过恒等变形从形式上将一种问题转换成另一种问题，有时甚至可以通过数值上的近似来完成转换。从形式上看，牛顿法就是这样的一个恒等变形: $$ f(x) = 0 \Longleftrightarrow x = g(x) = x - \frac{f(x)}{f'(x)} $$

对于一个求方程根 $f(x) = 0$ 的问题，我们可以有无数种方案，将其转换成一个不动点 $x = g(x)$ 的问题。那么究竟那种转换方式更好呢？

2. 单点迭代的收敛速度

先来看一下方程 $(x - 1)(x + 1)(x - 2) = 0$，显然它有三个根 $x = ±1, x = 2$。将其展开有 $x^3 - 2x^2 - x + 2 = 0$。通过一些数学形式上的恒等变换，我们可以写出如下四种不动点方程。

$$ g_1(x) = x = x^3 - 2x^2 + 2 $$ $$ g_2(x) = x = \sqrt{\frac{1}{2}(x^3 - x + 2)} $$ $$ g_3(x) = x = 2 + \frac{1}{x} - \frac{2}{x^2} $$ $$ g_4(x) = x = x - \frac{x^3 - 2x^2 - x + 2}{3x^2 - 4x - 1} $$

	$g_1(x)$	$g_2(x)$	$g_3(x)$	$g_4(x)$
1	0.990101	1.005062	1.029507	0.999948
2	1.009996	1.002547	1.084341	1.000000
3	0.990105	1.001278	1.221242	1.000000
4	1.009992	1.000640	1.477846	1.000000
5	0.990109	1.000320	1.760921	1.000000
6	1.009988	1.000160	1.922899	1.000000
7	0.990113	1.000080	1.979148	1.000000
8	1.009984	1.000040	1.994677	1.000000
9	0.990117	1.000020	1.998662	1.000000

根 $x = 1$ 显然在闭区间 $[0.9, 1.1]$ 内，我们用相同的初值 $x_0 = 1.01$ 对这四个不动点方程分别进行了 9 次迭代，如上面右侧的表格所示。从表中，我们可以看到 $g_1(x)$ 的序列在 0.99+ 与 1.00+ 来回振荡，振幅有在不断变小，但速度很慢。$g_2(x)$ 是收敛的，但相比于 $g_4(x)$ 来说它的收敛速度比较慢。我们是按照牛顿法的迭代公式构造的 $g_4(x)$ 所以它收敛很快。 $g_3(x)$ 干脆就发散了距离解 $x = 1$ 越来越远。

从上面的例子可以看出，选择合理的不动点方程，对于算法效率有很大影响。下面我们来研究一下不动点方程的迭代收敛速度。

假设 $\{p_n\}$ 是一个收敛到 $p$ 的数列，并且 $p_n \neq p, \forall n \in N^*$。如果存在正数 $\lambda$ 和 $\alpha$ 使得如下极限成立

$$ \lim_{n \to \infty} \frac{| p_{n+1} - p |}{| p_n - p |^{\alpha}} = \lambda $$

那么，我们说数列 $\{p_n\}$ 以 $\alpha$ 的阶数收敛到 $p$，渐进误差常数为 $\lambda$。对于不动点迭代序列 $p_{n+1} = g(p_n)$ 而言，如果 $\{p_n\}$ 以 $\alpha$ 阶收敛于 $p$，那么该迭代的收敛过程就是 $\alpha$ 阶的。阶数越高，收敛速度就越快。

假设函数 $g(x)$ 在闭区间 $[a, b]$ 上满足不动点存在且唯一的两个充分条件。此外，我们还假设 $g'$ 在开区间 $a,b$ 上连续并且有一个正数 $k < 1$ 使得 $|g'(x)| ≤ k, \forall x \in (a, b)$。如果 $g'(p) \neq 0$，那么对于任意 $p_0 \neq p, p_0 \in [a, b]$，序列 $p_{n+1} = g(p_n)$ 只能线性收敛到唯一不动点 $p$。即收敛阶数为 1。根据中值定理，有

$$ p_{n+1} - p = g(p_n) - g(p) = g'(\xi_n)(p_n - p) $$

其中 $\xi_n$ 是 $p_n$ 与 $p$ 之间的一个数。因为 $\{p_n\}$ 收敛于 $p$，所以 $\{\xi_n\}$ 也收敛于 $p$，又 $g'$ 在 $(a,b)$ 上连续，有 $\lim_{n \to \infty} g'(\xi_n) = g'(p)$。因此

$$ \lim_{n \to \infty} \frac{p_{n+1} - p}{p_n - p} = \lim_{n \to \infty}g'(\xi_n) = g'(p) \Longrightarrow \lim_{n \to \infty} \frac{|p_{n+1} - p|}{|p_n - p|} = |g'(p)| $$

因此，只要 $g'(p) \neq 0$，那么序列 $\{p_n\}$ 只能线性收敛到 $p$。

如果 $g'(p) = 0$，$g''$ 在包含 $p$ 的一个区间 $(p - \delta, p + \delta)$ 上连续且 $|g''(x)| < M$。我们在区间 $(p - \delta, p + \delta)$ 上对 $g(x)$ 泰勒展开有:

$$ g(x) = g(p) + g'(p)(x - p) + \frac{g''(\xi)}{2}(x - p)^2 $$

代入 $g(p) = p, g'(p) = 0$ 有:

$$ g(x) = p + \frac{g''(\xi)}{2}(x - p)^2 $$

特别的 $g(p_n) = p_{n+1}$，有

$$ p_{n+1} = p + \frac{g''(\xi_n)}{2}(p_n - p)^2 \Longrightarrow \lim_{n \to \infty} \frac{|p_{n+1} - p|}{|p_n - p|^2} = \frac{|g''(p)|}{2} < \frac{M}{2} \Longrightarrow |p_{n+1} - p| < \frac{M}{2}|p_n - p|^2 $$

所以，如果 $g'(p) = 0$ 那么 $\{p_n\}$ 至少二阶收敛。

3. 加速收敛

根据前文的推导过程，$g'(p) \neq 0$ 时有极限 $\lim_{n \to \infty}\frac{p_{n+1} - p}{p_n - p} = g'(p)$。所以存在一个正整数 $N > 0$，使得 $n > N$ 时有近似关系:

$$ \frac{p_{n+1} - p}{p_n - p} \approx \frac{p_{n+2} - p}{p_{n+1} - p} \Longrightarrow (p_{n+1} - p)^2 \approx (p_{n+2} - p)(p_n - p) $$

上式展开有

$$ \begin{aligned} p_{n+1}^2 - 2 p_{n+1}p + p^2 & \approx p_{n+2}p_n - (p_n + p_{n+2})p + p^2 \\ (p_{n+2} + p_n - 2p_{n+1})p & \approx p_{n+2}p_n - p_{n+1}^2 \\ p & \approx \frac{p_{n+2}p_n - p_{n+1}^2}{p_{n+2} + p_n - 2p_{n+1}} \\ p & \approx p_n - \frac{(p_{n+1} - p_n)^2}{p_{n+2} - 2p_{n+1} + p_n } \end{aligned} $$

我们将上式记为 $\hat{p_n}$。它认为上式是对 $p$ 的更好的估计，它认为数列 $\{\hat{p_n}\}$ 比 $\{p_n\}$ 收敛的更快

$$ \hat{p}_n \approx p_n - \frac{(p_{n+1} - p_n)^2}{p_{n+2} - 2p_{n+1} + p_n } $$

上述方法被称为Aitken $\mathbb{\Delta^2}$ 方法，给定一个收敛速度缓慢的序列，可以用它生成一个速度更快的新序列。后来 Steffensen 将它与不动点迭代本身结合起来。每迭代两步就用 Aitken 插入一步预测值，并把预测值作为下一次迭代的起点，从而衍生出了具有二次收敛速度的 Steffensen 迭代法。其实现过程很简单，如下代码所示，为了防止计算 Aitken $\mathbb{\Delta^2}$ 出现除以 0 的情况，我们增加了 7-8 两行，检查连续两次不动点迭代的值是否足够接近。

        template <typename DataType>
        DataType AitkenSteffensen(std::function<DataType(DataType)> g, DataType p0, int max_iter = 100, DataType tol = SMALL_VALUE)
        {
            for (int i = 0; i < max_iter; ++i) {
                DataType p1 = g(p0);
                DataType p2 = g(p1);
                if (std::abs(p2 - p1) < tol)
                    return p2;
                DataType p = p0 - (p1 - p0)*(p1 - p0)/(p2 - 2*p1 + p0);
                if (std::abs(p - p0) < tol)
                    return p;
                p0 = p;
            }
            return p0;
        }

4. 完

对于一个求方程根 $f(x) = 0$ 的问题，我们可以有无数种方案，将其转换成一个不动点 $x = g(x)$ 的问题。如果可以，我们希望得到的不动点函数 $g(x)$ 在 $p$ 处的导数为0，即 $g'(p) = 0$。这样序列 $p_{n+1} = g(p_n)$ 就会以至少二阶的速度收敛到 $p$。牛顿法就是这样一种不动点迭代方法。对于只有一阶线性收敛速度的 $g(x)$ 我们也可以通过 Aitken $\Delta^2$ 方法和 Steffenson 迭代得到二次收敛速度。