任意方程求根

1.簡介

方程和函數(shù)是代數(shù)數(shù)學(xué)中最為重要的內(nèi)容之一，從初中直到大學(xué)，我們都在研究著方程與函數(shù)，甚至我們將圖形代數(shù)化，從而發(fā)展出了代數(shù)幾何、解析幾何的內(nèi)容。而在方程與函數(shù)中，我們研究其性質(zhì)最多的，往往就是方程的根（零點(diǎn)），即使是研究方程的極值點(diǎn)、鞍點(diǎn)等，我們無非也只是研究其微商的零點(diǎn)。
我們在初等數(shù)學(xué)中已經(jīng)學(xué)過許多簡單初等函數(shù)、線性方程的求解方法，在本文中，我們重點(diǎn)討論任意方程，尤其是計算困難的非線性方程的求根方法。

2.方程

2.1分類和介紹

方程就是指含有未知數(shù)的等式。是表示兩個數(shù)學(xué)式（如兩個數(shù)、函數(shù)、量、運(yùn)算）之間相等關(guān)系的一種等式，使等式成立的未知數(shù)的值稱為“解”或“根”。在這里，根據(jù)一些性質(zhì)的不同，我們將方程分成以下幾類：

單個方程
- 線性方程：本質(zhì)是等式兩邊乘以任何相同的非零數(shù)，方程的本質(zhì)都不受影響。通常認(rèn)為只含有一次項(xiàng)的方程。
- 非線性方程：是因變量與自變量之間的關(guān)系不是線性的關(guān)系的方程。
  - 多項(xiàng)式方程
  - 超越方程：指含有未知量的超越式（指數(shù)、對數(shù)、三角函數(shù)、反三角函數(shù)等）的方程。換言之，超越方程中都有無法用含有未知數(shù)的多項(xiàng)式、分式或開方表示的式子。
多個方程
- 線性方程組
- 非線性方程組

2.2方程的零點(diǎn)（根、解）

若有一個值或一些值能夠使得方程 $f(x)=0$ 成立，那么這個值就被成為方程的解，也常常被叫做零點(diǎn)和根。

若方程有且只有一個解 $x^*$ ，那么我們稱方程有單根 $x^*$ 。

若對于方程 $f(x)=0$ ，有 $f(x^*) = 0,f^{'}(x^*)=f^{''}(x^*)=\cdots=f^{(k)}(x^*)=0,f^{(k+1)}(x^*)\neq0$ ,那么稱 $x^*$ 為方程的k+1重根

PS：若方程是簡單冪函數(shù)多項(xiàng)式組成，那么方程的解的數(shù)量應(yīng)和最高此項(xiàng)的數(shù)值一致，因?yàn)榇嬖谔摳?/p>

3.求根方法

求根的方法基本上大同小異，都是通過區(qū)間去逼近方程的根的點(diǎn)。

首先我們說一個定理1：對于實(shí)函數(shù)方程 $f(x)=0$ ，當(dāng) $x\in(a,b)$ ，且 $f(x)$ 在 $x\in(a,b)$ 時單調(diào)且連續(xù)，若 $f(a)\cdot f(b)<0$ ，則方程在 $x\in(a,b)$ 有且只有一個根。

3.1二分法

3.1.1普通二分法的原理及操作

二分法和算法中的二分搜索法非常的類似，取定有根區(qū)間 $[a,b]$ 進(jìn)行對分，求得 $mid = \frac{a+b}{2}$ 進(jìn)行判斷含根區(qū)間，如果 $f(a)\cdot f(mid)<0$ ,則令 $b=mid$ ；反之若 $f(b)\cdot f(mid)<0$ ,則令 $a=mid$ 。當(dāng) $|b_n-a_n|<\epsilon$ 停止計算返回結(jié)果。

產(chǎn)生的截斷誤差為 $|e_{n-1}| = x_{n+1} - x^*\leq[b_n - a_n] = \frac{b_0 - a_0}{2^n}$ 。

可以計算出最小迭代次數(shù)為 $n = \frac{lg(c_0-a_0)-lg\epsilon}{lg2}$

代碼實(shí)現(xiàn)（更多語言的代碼見倉庫中Code文件夾）：

private static double epsilon = 0.001;
// func為函數(shù)，寫法如x=>x*x+2*x-1，a，b必須為有效的含根開區(qū)間
public static double Binary(Func<double, double> func, double a, double b)
{
    var f1 = func.Invoke(a);
    var f2 = func.Invoke(b);
    if (f1 * f2 > 0)
        throw new ArgumentException("此區(qū)間無根或根不唯一");
    double mid = (a + b) / (double)2;
    var fm = func.Invoke(mid);
    if (fm == 0)
        return fm;
    if (f1 * fm < 0)
        b = mid;
    else if (f2 * fm < 0)
        a = mid;
    if (Math.Abs(b - a) <= epsilon)
        return (a + b) / (double)2;
    return Binary(func, a, b);
}

3.1.2普通二分法準(zhǔn)確度及速度分析

假設(shè) $\left[a,b\right]$ 是二分法的初始區(qū)間，在進(jìn)行n此二分之后，得到的最終根的分布區(qū)間 $[a_n,b_n]$ 的長度為 $(b-a)/2^n$ ，我們將其中點(diǎn)作為根的最優(yōu)估計值，與真實(shí)值之間的誤差不超過區(qū)間長度的一半。

$e = \left|x_c-r\right|<\frac{b-a}{2^{n+1}}$
若 $e<\frac{1}{2}\times10^{-p}$ ,則精確到小數(shù)點(diǎn)后p位，

事實(shí)上我們也可以簡單的計算出函數(shù)執(zhí)行的次數(shù)為n+2次。

3.2淺談迭代

3.2.1 迭代是什么

迭代是重復(fù)反饋過程的活動，其目的通常是為了逼近所需目標(biāo)或結(jié)果。每一次對過程的重復(fù)稱為一次“迭代”，而每一次迭代得到的結(jié)果會作為下一次迭代的初始值。
重復(fù)執(zhí)行一系列運(yùn)算步驟，從前面的量依次求出后面的量的過程。此過程的每一次結(jié)果，都是由對前一次所得結(jié)果施行相同的運(yùn)算步驟得到的。這篇文章中對于迭代有一個非常不錯的概述究竟什么是迭代？

3.2.2 不動點(diǎn)的定義

各位可以嘗試以下操作，

隨意輸入一個數(shù)字 $\lambda$
然后對其進(jìn)行 $cos\lambda$ 運(yùn)算，
將運(yùn)算結(jié)果作為新的值傳回 $cosx$ 函數(shù)之中

當(dāng)你一直重復(fù)以上三個操作，你會發(fā)現(xiàn)數(shù)字最后會定格在0.73908513左右不動了。這是一個非常有趣的現(xiàn)象，事實(shí)上我們的這個操作就是 $x_n = cosx_{n-1}$ ，而最后不變化的數(shù)值實(shí)際上就是 $x = cosx$ 這個方程的解。

這就是我們不動點(diǎn)的一種實(shí)際情況，不動點(diǎn)原理是數(shù)學(xué)上一個重要的原理，也叫壓縮映像原理或巴拿赫（Banach）不動點(diǎn)定理，完整的表達(dá)：完備的度量空間上，到自身的一個壓縮映射存在唯一的不動點(diǎn)。用初等數(shù)學(xué)可以這么理解：連續(xù)映射 $f$ 的定義域包含值域，則存在一個 $x$ 使得 $f(x)=x$ 。

若某函數(shù)滿足 $f(\lambda)=\lambda,\lambda \in R$ ，我們就稱 $\lambda$ 為函數(shù)的一階不動點(diǎn)。同樣的，我們推廣一下，若 $f(f(\lambda))=\lambda,\lambda \in R$ ，則稱為二階不動點(diǎn)，一階不動點(diǎn)必定是二階不動點(diǎn)。

不動點(diǎn)的存在定理：若某函數(shù) $y=f(x)$ 與 $y=x$ 存在至少一個交點(diǎn)，那么函數(shù)必然存在不動點(diǎn)。

3.2.3 函數(shù)的相似性

我們常常說圖形之間的相似，事實(shí)上函數(shù)也有其相似的定義。若函數(shù) $f(x)$ 對其換元，令 $t = \varphi(x),x=\varphi^{-1}(t),f(x)=>f(\varphi^{-1}(t))$ 。我們的不動點(diǎn)是一個基于迭代的過程，迭代就是讓 $f(x_0) = x_1 = t,f(x_1)=f(f(x_0))=f(t)=x_2...$ 的過程，如果我們換個角度去思考，迭代不就是和我們剛剛提到的換元是一個道理嗎？那么我們思考一下，如果對于二次迭代 $f(f(x))$ ，如果說我們把函數(shù)完全看成是一個復(fù)合函數(shù)，那么我們現(xiàn)在需要做的就是試圖將函數(shù)寫成 $f(\varphi^{-1}(t))$ 的形式。

現(xiàn)在做出如下變換，將外層的 $f(x)$ 換元為 $f(\varphi^{-1}(t))$ ，那么二次迭代式就變成了 $g(t) = \varphi(f(\varphi^{-1}(t)))$ 。隨后對我們的 $g(x)$ 再進(jìn)行迭代就會得到 $g(g(t)) = \varphi(f(f(\varphi^{-1}(t)))$ ，利用數(shù)學(xué)歸納法計算n次迭代式，那么就會得到 $g_n(t) = \varphi(f_n(\varphi^{-1}(t)))$ 。

我們在這給出總結(jié)和定義：若函數(shù) $g(x),f(x),\varphi(x)$ ，若有 $g(x)=\varphi(f(\varphi^{-1}(x)))$ ，那么稱函數(shù) $g(x)$ 和 $f(x)$ 通過 $\varphi(x)$ 相似，記作 $f\sim g$ ， $\varphi(x)$ 稱為相似的橋函數(shù)，且之前我們的過程得出了一個重要的結(jié)論就是 $g\sim f=>g^n\sim f^n$ ，同時，相似函數(shù)的不動點(diǎn)是完全一致的，證明如下：

$若函數(shù)g(x)的不動點(diǎn)為x_0，則有$
$g(x_0)=x_0=\varphi(f(\varphi^{-1}(x_0)))$

$根據(jù)反函數(shù)性質(zhì)有\(zhòng)varphi(\varphi^{-1}(x))=x$

$故f(\varphi^{-1}(x_0)) = \varphi^{-1}(x_0)=x_0,\varphi^{-1}(x_0)是f(x)的一個不動點(diǎn)$

通過上述得出的不動點(diǎn)，很容易通過剛才講的內(nèi)容得出 $\varphi^{-1}(x_0)$ 實(shí)際上也是 $g(x)$ 換元后得到的不動點(diǎn)的位置，因此相似函數(shù)的不動點(diǎn)是一一對應(yīng)的。

3.2.4 迭代收斂速度

迭代是一種逼近、利用數(shù)列、函數(shù)收斂的性質(zhì)進(jìn)行求解的方法，那么對于收斂的速度，我們很難直接從數(shù)值看出速度，因?yàn)槭諗窟^程中，誤差的變化是越來越小的，因此我們再次進(jìn)行一個比階。

$\exists p\in N \And C>0,\displaystyle \lim_{k \to \infty}\frac{e_{k+1}}{e_k^p} = C$

其中：

$e_k = |x_k-x^*|$

$p= \begin{cases} p=1,線性收斂\\ p=2,平方收斂\\ p>2,高階收斂 \end{cases}$

并且有一個定理：若 $\varphi(x)$ 在所求的根 $x^*$ 鄰域有連續(xù)的二階導(dǎo)數(shù)，并且有 $|\varphi^{'}(x)|<1$ ，有以下特點(diǎn)：

當(dāng) $\varphi^{'}(x)\neq0$ ，那么迭代過程線性收斂
當(dāng) $\varphi^{'}(x)=0,\varphi^{''}(x)\neq0$ ，那么迭代是平方收斂的。

證明過程留給讀者，只需要利用泰勒展開是便可以證明該定理。

3.3不動點(diǎn)迭代法

3.3.1不動點(diǎn)迭代法的原理及操作

迭代法是將方程求根問題轉(zhuǎn)換成了函數(shù)求交點(diǎn)問題再轉(zhuǎn)換成數(shù)列求極限的問題，這個過程中，對方程進(jìn)行一個巧妙的變換之后，方程就可以在若干次迭代之后解出一個近似解。

操作方法如下：首先將方程 $f(x)=0$ 改寫成 $x = g(x)$ 的形式，這樣就可以將方程的解看成是函數(shù) $y=x$ 和 $y=g(x)$ 的交點(diǎn)。給定初始值 $x_0$ 后，則 $x_1 = g(x_0)$ ，不斷重復(fù)這個過程，若 $\displaystyle \lim_{k \to \infty}x_k$ 存在，則迭代便可以達(dá)到使得 $x_k$ 趨于交點(diǎn)，而這個交點(diǎn)就是我們剛才講了那么久的不動點(diǎn)。

我們即使構(gòu)造出了 $x=g(x)$ 的迭代式，往往由于我們的變換只是一些簡單的移項(xiàng)、通分等四則運(yùn)算獲得，使得迭代式也并不是一個很容易被計算的函數(shù)，這個時候，我們之前講到的那么多函數(shù)相似就派上了用場，我們可以利用函數(shù)相似性去尋找一個更為簡單的函數(shù) $\varphi(x)\sim g(x)$ ，由于相似的性質(zhì)，不動點(diǎn)并不會發(fā)生變化。不過，我們需要保證最終的迭代函數(shù)在指定的含根區(qū)間內(nèi)存在一階導(dǎo)數(shù)且導(dǎo)數(shù)值 $|\varphi^{'}(x)|<1$ ，否則迭代函數(shù)將會不收斂。

迭代過程收斂的情況如下圖所示：

avatar

如果說我們的導(dǎo)數(shù)不滿足條件，那么我們的迭代將會發(fā)散：

avatar

代碼樣例

// func是迭代函數(shù)而不是實(shí)際函數(shù)
public static double Iterative(Func<double, double> func, double x)
{
    double temp = func.Invoke(x);
    if (Math.Abs(temp - x) <= epsilon)
    {
        return temp;
    }
    x = temp;
    return Iterative(func, x);
}

3.3.2迭代法的收斂性證明

在這里，我們將證明迭代法求根的合理性和可行性。

前提條件：設(shè)函數(shù) $\varphi(x), x\in[a,b]$ 有連續(xù)的一階導(dǎo)數(shù)，并且滿足以下條件：

$\forall x\in [a,b],\varphi(x)\in[a,b]$
$\exists L \in (0,1),\forall x\in[a,b],|\varphi^{'}(x)|\leq L$

要證明和解決以下命題和問題：

在 $x\in[a,b],\exists x^*,\varphi(x^*) = 0$
$\forall x_0\in[a,b]$ ，迭代過程 $x_{k+1} = \varphi(x_k)$ 均收斂與 $x^*$
求解誤差分析式

現(xiàn)在開始證明

1：證明在區(qū)間內(nèi)有且只有一個根存在：

證：在 $x \in [a,b]$ 時， $\varphi^{'}(x)$ 存在，所以有 $\varphi(x)$ 連續(xù)，于是可以作 $g(x) = x - \varphi(x)$ ，易知 $g(x)$ 連續(xù)。

又因?yàn)?img class="math-inline" src="https://math.jianshu.com/math?formula=%5Cvarphi(x)%20%5Cin%20%5Ba%2Cb%5D" alt="\varphi(x) \in [a,b]" mathimg="1">，且 $g(a)*g(b)<0$ ,故存在實(shí)根，使得 $x=\varphi(x)$ 。

利用反證法：若在[a,b]上還有一實(shí)根 $\bar{x}$ ，那么通過拉格朗日中值定理必定有：

$x^*-\bar{x} = \varphi(x^*)-\varphi(\bar{x}) = \varphi^{'}(\xi)(x^*-\bar{x})\Longrightarrow \varphi^{'}(\xi) = 1$

顯然與假設(shè)不符合。

2：證明這個根收斂于 $x^*$

根據(jù)拉格朗日中值定理，容易知道

$x^*-x_{k+1} = \varphi(x^*)-\varphi(x_k) = \varphi^{'}(\xi)(x^*-x_k)$

又因 $|\varphi^{'}(x)|\leq L$ ，故易得：

$|x^*-x_{k+1}|\leq|L(x^*-x_k)|=|L^2(x^*-x_{k-1})=\cdots=|L^k(x^*-x_0)|$

因?yàn)?img class="math-inline" src="https://math.jianshu.com/math?formula=%5Cdisplaystyle%20%5Clim_%7Bk%20%5Cto%20%5Cinfty%7DL%5Ek%3D0" alt="\displaystyle \lim_{k \to \infty}L^k=0" mathimg="1">，故 $\displaystyle \lim_{k \to \infty}|x^*-x_{k+1}|=0$ (絕對值必定非負(fù))，得 $x^*=x_{k+1}(k\to\infty)$

3:迭代法的誤差式：

設(shè)某一次迭代后誤差為 $\epsilon$ ，則可以推出：

$|x_{k+1}-x_{k}| = |(x^*-x_k)-(x^*-x_{k+1})\geq|x^*-x_k|-|x^*-x_{k+1}|\geq|x^*-x_k|-L|x^*-x_k|$

其中從 $|x^*-x_{k+1}| \leq L|x^*-x_k|$ 則是因?yàn)閷τ谝粋€可以收斂的函數(shù)而言，每一次迭代后，誤差總是減少的。

故可以輕松的計算出誤差估計式為：

$|x^*-x_k|\leq\frac{1}{1-L}|x_{k+1}-x_k|$

通過中值定理，又可以推出：

$|x_{k+1}-x_k| = |\varphi'(\xi)(x_k-x_{k-1})|$

因?yàn)?img class="math-inline" src="https://math.jianshu.com/math?formula=%7C%5Cvarphi'(%5Cxi)%7C%5Cleq%20L" alt="|\varphi'(\xi)|\leq L" mathimg="1">，將上式遞推后放縮成第二種誤差估計式：

$|x^*-x_k|\leq\frac{L^k}{1-L}|x_1-x_0|$

其中L可以近似的看作是一個足夠接近方程解的函數(shù)導(dǎo)數(shù)的值。

3.3.3不動點(diǎn)迭代法的收斂速度

首先我們先直接給出一個結(jié)論：在不動點(diǎn)迭代過程中，第i步迭代和第i+1步迭代的表達(dá)式分別為 $e_i,e_{i+1}$ ，且有 $\displaystyle \lim_{ i \to \infty}\frac{e_{i+1}}{e_i}=S<1\neq0$ ，因此不動點(diǎn)迭代法是一個線性收斂的算法。

如何能證明我們的算法是一個線性收斂的函數(shù)呢？實(shí)際上非常簡單。證明過程如下：

$若\varphi(x)連續(xù)可微,不動點(diǎn)的一個足夠近似的值為x_0，x_{i+1}和x_i是兩次迭代的估計值$

$根據(jù)中值定理，\exists r\in (x_i,r),\frac{\varphi(x_i)-x_0}{x_i-x_0}=\varphi'(r)$

$\varphi(x_i) = x_{i+1},故x_{i+1}-r = \varphi'(r)(x_i-r)$

$e_{i+1}=\left|\varphi'(r)\right|e_i，即S=\left|\varphi'(r)\right|為常數(shù)$

我們可以在這里舉幾個例子，例如我們希望求函數(shù) $f(x) = x^3+x-1=0$ 的根，倘若我們?nèi)?img class="math-inline" src="https://math.jianshu.com/math?formula=%5Cvarphi(x)%20%3D%201-x%5E3" alt="\varphi(x) = 1-x^3" mathimg="1">，根據(jù)我們誤差分析的式子，包括我們收斂分析的公式來說，后一步的誤差將會是 $e_{i+1}=-3r^2e_i$ ，我們預(yù)估的根大概是0.6823左右，可以發(fā)現(xiàn)我們的收斂速度的絕對值已經(jīng)大于1，意味著每一步的誤差都在擴(kuò)大，因此函數(shù)不會收斂。如果我們?nèi)?img class="math-inline" src="https://math.jianshu.com/math?formula=%5Cvarphi(x)%20%3D%20%5Csqrt%7B1-x%7D" alt="\varphi(x) = \sqrt{1-x}" mathimg="1">，我們繼續(xù)計算出的收斂速度大約是0.7左右，很明顯誤差在以線性縮小。如果再次修改我們的遞推式 $\varphi(x) = \frac{1+3x^3}{1+3x^2}$ ，我們的S幾乎是等于0的，因此也解釋了這種方法為何是最快收斂到我們的根。

//todo 圖片

再舉一個例子，函數(shù) $f(x) = cosx-sinx=0$ 的迭代求根過程，我先給定我們的近似值大約是0.7853982。這里似乎很難構(gòu)造我們的迭代遞推式，事實(shí)上只需要方程兩邊同時加x就可以湊出 $x=cosx-sinx+x=\varphi(x)$ 。我們利用計算機(jī)編寫好程序并計算，最后得到了這樣一組迭代過程中的一些數(shù)據(jù)：

誤差圖

不難看出我們的 $\frac{e_i}{e_{i-1}}$ 在一個恒定值0.414保持了相當(dāng)?shù)拇螖?shù)，這也是因?yàn)?img class="math-inline" src="https://math.jianshu.com/math?formula=S%3D%5Cleft%7C%5Cvarphi'(r)%5Cright%7C%3D%7C1-sinr-cosr%7C%5Capprox0.414" alt="S=\left|\varphi'(r)\right|=|1-sinr-cosr|\approx0.414" mathimg="1">，印證了我們的推論是完全正確的，當(dāng)我們不動點(diǎn)迭代在一個區(qū)間內(nèi)，往往會以線性模式進(jìn)行誤差的縮小。

迭代的終止

3.4牛頓法

牛頓在整個微積分、數(shù)值分析中都有著舉足輕重的地位，這里闡述的牛頓法，全名牛頓-拉夫遜法，就是Taylor展開式的一部分。牛頓迭代法的核心思想就是：設(shè)法將一個非線性方程 $f(x)=0$ 轉(zhuǎn)化為某種線性方程去求解。在數(shù)學(xué)意義上，我們知道泰勒公式可以將任意函數(shù)以簡單冪函數(shù)的形式表示出來，而在幾何意義上，我們是利用切線與X軸的交點(diǎn)去進(jìn)行迭代，在根處，切線必過零點(diǎn)。若設(shè) $f(x)=0$ 的近似解為 $x_k$ ，則方程可以用一階Taylor展開進(jìn)行近似表示，牛頓迭代法的核心公式如下：

$p_1(x) = f(x_k)+f^{'}(x_k)(x-x_k)$

3.4.1普通牛頓法

從上述公式中我們知道，其中 $p_1(x)$ 就是泰勒多項(xiàng)式的表達(dá)，將其看成是方程 $f(x)=0$ ，通過迭代的思想，從而得到了一個線性方程：

$x_{k+1} = x_k - \frac{f(x_k)}{f^{'}(x_k)}$

這就是普通牛頓法的迭代公式。在幾何意義上，他就是一個切線與x軸的交點(diǎn)去逼近零點(diǎn)，如圖：

newtown

我們不斷的通過迭代這個過程，就能讓我們的值越來越精確。

3.4.2牛頓下山法

所有的迭代法都有一個無法逃過的缺點(diǎn)，如果選取的初始值離根太遠(yuǎn)，則會導(dǎo)致迭代過程次數(shù)過多，并且有可能導(dǎo)致迭代發(fā)散，因?yàn)榈ㄖ痪哂芯植渴諗啃浴?/p>

為了避免迭代失敗或時間過長，我們加上這樣一個條件用于保證數(shù)值穩(wěn)定下降

$|f(x_{k+1})|<|f(x_k)|$

將這個條件和牛頓法結(jié)合再一起，即再保證數(shù)值穩(wěn)定下降的過程中，利用牛頓法加快迭代，這就是牛頓下山法。具體的操作如下：

將牛頓法的結(jié)果 $\bar{x}_{k+1}$ 與前一步的迭代值 $x_k$ 進(jìn)行加權(quán)平均作為新的迭代值 $x_{k+1}$ ，

則有：

$x_{k+1} = \lambda\bar{x}_{k+1} + (1-\lambda)x_k$
或
$x_{k+1} = x_k - \lambda\frac{f(x_k)}{f^{'}(x_k)}$

其中 $\lambda(0\leq\lambda<1)$ 稱為下山因子，它的值是一個逐步試探的過程，可以從1開始取值，一旦滿足 $|f(x_{k+1})|<|f(x_k)|$ 則稱為下山成功，否則需要另選初始值 $x_0$ 進(jìn)行試算。

3.4.3簡單牛頓法

牛頓法可以說是一個非常有效的計算方法，準(zhǔn)確度和迭代次數(shù)上都要比普通的迭代法要好得多，但是牛頓法最大的問題是我們需要求方程的導(dǎo)數(shù)，對于某些極其復(fù)雜的函數(shù)而言，導(dǎo)數(shù)是無法通過人工的方式計算，假如我們使用微積分的方式去求解導(dǎo)數(shù)，這對整個程序的性能會有比較大的影響，因此我們可以利用一個常數(shù)值 $\lambda$ 來代替方程中的導(dǎo)數(shù)項(xiàng)。此時迭代公式為：

$x_{k+1} =x_k - \frac{f(x_k)}{\lambda}$

不過對于常數(shù) $\lambda$ 的取值是有限制的，因?yàn)槲覀冃枰ＷC迭代函數(shù)的收斂性，如果函數(shù)不收斂于 $x^*$ ，那么一切都沒有意義。于是有：

$\varphi(x) = x - \frac{f(x)}{\lambda}\Longrightarrow\varphi^{'}(x) = 1-\frac{f^{'}(x)}{\lambda}$

牛頓迭代法的收斂性遵循前文中普通迭代法的收斂性，于是可以得到：

$|\varphi^{'}(x)| = |1-\frac{f^{'}(x)}{\lambda}|\Longrightarrow0<\frac{f^{'}(x)}{c}<2$

這樣我們就可以很輕松的確定下c的值了。

簡單牛頓法的幾何意義就簡單許多了，和我們之前討論的普通迭代法一致，只不過普通迭代法是將函數(shù)值和 $y=x$ 進(jìn)行處理變換，而簡單牛頓法則是和 $y= \lambda(x-x_k)+f(x_k)$ 進(jìn)行變換，本質(zhì)是一致的。

///這里只寫普通牛頓法，另外的函數(shù)由讀者自己思考
// 其中f1x為導(dǎo)數(shù)
public static double Newtown(Func<double, double> fx, Func<double, double> f1x, double x)
{
    var temp = f1x.Invoke(x);
    if (temp == 0)
    {
        throw new ArgumentException();
    }
    x = x - fx.Invoke(x) / temp;
    if (Math.Abs(fx.Invoke(x)) <= epsilon)
    {
        return x;
    }
    return Newtown(fx, f1x, x);
}

3.4.4 改進(jìn)牛頓法——弦截法

弦截法是牛頓法的一種改進(jìn)，保留了牛頓法中收斂速度快的優(yōu)點(diǎn)，還克服了再牛頓法中需要計算函數(shù)導(dǎo)數(shù)值 $f^{'}$ 的缺點(diǎn)。

弦截法中，我們用差商

$\frac{f(x_k)-f(x_{k-1})}{x_k-x_{k-1}}$

去代替牛頓法中的導(dǎo)數(shù)值，于是可以得到以下離散化的迭代

$x_{k+1} = x_k - \frac{f(x_k)}{f(x_k)-f(x_{k-1})}(x_k-x_{k-1})$

這種方法叫做雙點(diǎn)弦截法，如圖所示：

弦截法

從圖中可以知道，弦截法一直利用兩點(diǎn)之間的連線作為迭代的內(nèi)容，那么，他的合理性在哪里呢？

整個迭代法都離不開中值定理，這里也是這樣，事實(shí)上，這個差商之所以可以對導(dǎo)數(shù)值進(jìn)行替代，是因?yàn)橹兄刀ɡ碇姓f過，連續(xù)函數(shù)中兩函數(shù)上的點(diǎn)的連線的斜率必為兩點(diǎn)之間某一點(diǎn)的導(dǎo)數(shù)值，并且迭代過程中，這兩點(diǎn)的中值會越來越逼近函數(shù)零點(diǎn)，事實(shí)上這已經(jīng)說明了弦截法是牛頓法的改進(jìn)方法了。

不過，如果將函數(shù)中 $x_{k-1}$ 用一個定點(diǎn) $x_0$ 代替，這種方法叫做單點(diǎn)弦截法，幾何意義如圖所示：

單點(diǎn)弦截法

//這里就對雙點(diǎn)弦截法進(jìn)行編程
public static double StringCut(Func<double, double> func, double x1, double x2)
{
    var temp = x1 - (func.Invoke(x1) / (func.Invoke(x1) - func.Invoke(x2))) * (x1 - x2);
    x2 = x1;
    x1 = temp;
    if (Math.Abs(func.Invoke(x1)) <= epsilon)
    {
        return x1;
    }
    return StringCut(func, x1, x2);
}

3.4.5 牛頓法的收斂性

牛頓法的收斂性有可能是一階收斂也有可能是二階收斂，具體要看我們的迭代函數(shù)。我們首先講講牛頓法二次收斂的情形。

假定我們牛頓法的迭代函數(shù)為 $g(x) = x-\frac{f(x)}{f'(x)}$ ，那么我們迭代函數(shù)的導(dǎo)數(shù)為 $g'(x) = \frac{f(x)f''(x)}{f'(x)^2}$ ，設(shè)根為r，那么易知 $g'(r)=0$ ，迭代函數(shù)是收斂的。

我們再將 $f(r)$ 進(jìn)行泰勒展開，取 $x_i$ 是一個足夠接近根的一次迭代結(jié)果，并且代入 $f(r)=0$ 則有

$f(r) = f(x_i)+(r-x_i)f'(x_i)+\frac{(r-x_i)^2}{2}f''(\xi),\xi \in (x_i,r)$

$f(r) = 0 =>x_i-\frac{f(x_i)}{f'(x_i)} -r = \frac{(r-x_i)^2}{2}\frac{f''(\xi)}{f'(x_i)}$

若 $f'(x_i)\neq0$ 對此函數(shù)迭代式在進(jìn)行一次迭代可以得到以下式子

$x_{i+1}-r=e_i^2\frac{f''(\xi)}{2f'(x_i)}$

$e_{i+1} = e_i^2\left|\frac{f''(\xi)}{2f'(x_i)}\right|$

這樣我們就得到了牛頓法的二階收斂的定義，對于二階收斂，我們要求是需要在該處的一階導(dǎo)數(shù)值不為0。

不過我們的牛頓法并不是總為二次收斂的，在面臨多重根的時候，牛頓法往往會倒退到線性收斂的速度。，例如函數(shù) $f(x)=x^m$ 利用牛頓法求根的時候，牛頓公式如下：

$x_{i+1} = x_i-\frac{x_i^m}{mx_i^{m-1}}=\frac{m-1}{m}x_i$

唯一的根是0，但是0是一個m-1重根，因此得到

$e_{i+1} = Se_i,S=\frac{m-1}{m}$

我們這里可以給出一個定義，若m+1階連續(xù)可微函數(shù)有一個m重根，那么利用牛頓法會線性局部收斂，誤差式如上。

收斂迭代的加速與精度

利用迭代法進(jìn)行收斂的計算時，往往面臨一個很麻煩的問題，就是當(dāng)運(yùn)行足夠多次迭代的時候，我們的精度往往會不降反升，這是因?yàn)橐环N類似搖擺的問題出現(xiàn)，例如我規(guī)定此刻你站在坐標(biāo)1上，我希望你走到坐標(biāo)0上，你每次可以向左或者向右走三步，那么第一次你向左走三步之后，得到的迭代結(jié)果時-2，反而離精確值更遠(yuǎn)了。這是一個令人糟心的問題。這里我們對這種現(xiàn)象進(jìn)行一次解釋。

收斂的精度

y還是x？

事實(shí)上我們的誤差分為了前向誤差和后向誤差，什么是前向誤差呢？實(shí)際上就是我們要求解的根與實(shí)際的差值，假設(shè) $f(r)=0$ 這個方程求根，前向誤差定義為 $\left|r-x_i\right|$ ，而后向誤差指代的就是我們在函數(shù)值上的接近程度，也就是y軸的接近程度，我們定義為 $\left|f(x_i)|$ 。

對于計算機(jī)而言，存儲的精度是有極限的，假定我們現(xiàn)在有一臺古老的計算機(jī)，他的存儲精度最多就到小數(shù)點(diǎn)后4位，那么假設(shè)我們的迭代函數(shù)是 $f(x)=x^2$ ，我們知道在精度范圍內(nèi)，0.001是他能取到最小具有有效數(shù)字的精度解，如果我們現(xiàn)在通過迭代取到的數(shù)字是0.001的話，我們會發(fā)現(xiàn)一個驚人的問題，再進(jìn)行一次迭代之后，得到 $f(0.001)=0.000001$ ，由于我們的計算機(jī)無法存儲足夠的精度，也就是說這個值會被計算機(jī)看成是0，但我們知道這并不是我們的最優(yōu)解。

這其實(shí)就是我們后向誤差足夠小，但是前向誤差并不夠小導(dǎo)致的。我們在實(shí)際的代碼編寫的時候，一定要注意的一個地方就是我們的誤差在什么時候應(yīng)當(dāng)利用前向誤差進(jìn)行終止，什么時候又應(yīng)該利用后向誤差進(jìn)行求解。

敏感性

讀者可以嘗試將迭代初始值設(shè)置成我們方程實(shí)際的根進(jìn)行迭代，按著我們正常人的理解，他應(yīng)該立即停止迭代直接返回我們的初始值即可。但是算法往往是不盡人意的，在絕大多數(shù)工具或語言的內(nèi)置求根函數(shù)中，即使我們傳入的初始值是根，往往得到的結(jié)果反而是一個不精確的解。例如著名的威爾金森多項(xiàng)式：

$W(x) = (x-1)(x-2)...(x-20)$

如果我們利用這個已經(jīng)因式分解好了的式子進(jìn)行求根運(yùn)算，這似乎沒有必要，我們可以得到一個非常工整的方程根，他就是1-20之間所有的自然數(shù)組成，一旦我們將其展開后投入運(yùn)算，由于存在許多大數(shù)字，很多大數(shù)字會因?yàn)榇鎯?、計算等過程中產(chǎn)生損失和誤差，這些誤差實(shí)際上是很小的數(shù)字，但是卻對我們的計算結(jié)果產(chǎn)生了很大的誤差，最終我們難以得到一個精確解。這就是對于求根過程中對數(shù)字?jǐn)_動的敏感性問題。

我們進(jìn)一步去探索一下究竟是什么導(dǎo)致了誤差放大，并且能粗略的計算出這些誤差究竟有多少，我們假設(shè)是尋找一個函數(shù) $f(x)=0$ 的根r，同時對我們的函數(shù)做出一個很小的擾動 $\epsilon g(x)$ ，其中我們的系數(shù) $\epsilon$ 很小，用 $\Delta r$ 表示根的變化。那么會有以下式子。

$f(r+\Delta r)+\epsilon g(r+\Delta r) = 0$

將上式統(tǒng)統(tǒng)泰勒展開

$f(r)+(\Delta r)f'(r)+\epsilon g(r)+\epsilon(\Delta r)g'(r)+O((\Delta r)^2) = 0,其中\(zhòng)displaystyle O((\Delta r)^2) \to 0,\Delta r \to 0$

忽略后面的高階無窮小，同時r又是根，于是可以得到一個關(guān)于r的變化值的近似

$\Delta r (f'(r)+\epsilon g'(r)) \approx -f(r)-\epsilon g(r) =-\epsilon g(r)$

又因?yàn)?img class="math-inline" src="https://math.jianshu.com/math?formula=%5Cepsilon" alt="\epsilon" mathimg="1">也是一個很小的值（需要遠(yuǎn)小于 $f'(r)$ ），我們直接認(rèn)為他就是0，得到了 $\Delta r$ 的最終估計:

$\Delta r \approx -\epsilon\frac{g(r)}{f'(r)}$

這就是我們根的敏感性公式，利用這個公式，我們可以得出，如果出現(xiàn)了某些數(shù)字的擾動，根會如何變化。

這里舉一個黃皮書上的例子，如果函數(shù) $P(x) = (x-1)(x-2)...(x-6)-10^{-6}x^7$ ，求出它的最大根。我們?nèi)庋劭梢娮畲蟾鶓?yīng)該是在6附近，如果沒有最后那一項(xiàng)，我們的最大根就是6，那么問題就在于如果我去除最后一項(xiàng)，對根會有多大的影響呢？

那么我們直接設(shè)近似函數(shù)為 $f(x)=(x-1)(x-2)...(x-6)$ ，令我們的 $\epsilon g(x)=-10^{-6}x^7$ 。通過我們的根敏感性公式計算一番之后，得到了

$\Delta r \approx -\epsilon\frac{g(r)}{f'(r)} = 2332.8\times 10^{-6}$

因此我們估計的最大根應(yīng)當(dāng)是6.0023328左右，讀者可以自行使用迭代法去求解最大根，最后得出的結(jié)果與這個估計值是極度接近的。這里我們定義一個新東西誤差放大因子，也就是我們的 $\left|\frac{g(r)}{rf'(r)}\right|$ 。

普通迭代加速

對于迭代過程，利用中值定理有

$x^*-x_{k+1} = \varphi(x^*)-\varphi(x_k) = \varphi^{'}(\xi)(x^*-x_k)$

當(dāng) $\Delta x \rightarrow0$ ，我們將 $\varphi^{'}(\xi)$ 看成定值 $\lambda(\lambda<1)$

于是有

$\lambda(x^*-x_k) = x^*-x_{k+1}\longrightarrow x^* = \frac{1}{1-\lambda}\bar{x}_{k+1}-\frac{\lambda}{1-\lambda}x_k$

故可以推出最終的迭代公式為

$\begin{cases} \bar{x}_{k+1} = \varphi(x_k)\\ x_{k+1} = \bar{x}_{k+1} + \frac{\lambda}{1-\lambda}(\bar{x}_{k+1}-x_k) \end{cases}$

利用這種方式的好處就是再求得 $x_k$ 時，利用加權(quán)的方式，使得迭代法變得有點(diǎn)類似像牛頓迭代法一樣變成切線性質(zhì)的線而不是x=c或y=c的線，你經(jīng)過畫圖和簡單的代數(shù)運(yùn)算后，你會發(fā)現(xiàn) $\bar{x}_{k+1}<x_{k+1}$ ，也就是說達(dá)到了我們的加速目的。

Atiken加速法

Atiken加速法其實(shí)就是再普通加速迭代公式上在進(jìn)行一次迭代，這里直接寫出公式：

$\begin{cases} \bar{x}_{k+1} = \varphi(x_k)\\ \bar{\bar{x}}_{k+1} = \varphi(\bar{x}_{k+1})\\ x_{k+1} = \bar{\bar{x}}_{k+1} - \frac{(\bar{\bar{x}}_{k+1}-\bar{x}_{k+1})}{\bar{\bar{x}}_{k+1} - 2\bar{x}_{k+1}+x_k} \end{cases}$

練習(xí)題

試著不借助計算器計算 $\sqrt{2}$ 的值
試著證明 $f(x)=ax+b$ 在牛頓法中會一步收斂
編程題：由一個高度為10m的圓柱構(gòu)成的發(fā)射井頂部是一個半球，體積 $400m^3$ ，確定發(fā)射井底部班級
設(shè)函數(shù)為 $f(x)=x^n-ax^{n-1},g(x)=x^n$ ，利用敏感性公式手動估計 $f_\epsilon(x) = x^n-ax^{n-1}+\epsilon x^n$ 的非零根

Reference

《數(shù)值分析》（原書第二版） —— Timothy Sauer

《數(shù)值計算方法》（清華大學(xué)出版社） —— 呂同富等

色偷偷精品伊人,欧洲久久精品,欧美综合婷婷骚逼,国产AV主播,国产最新探花在线,九色在线视频一区,伊人大交九欧美,1769亚洲,黄色成人av

【機(jī)器學(xué)習(xí)】數(shù)值分析（1）—— 任意方程求根

【機(jī)器學(xué)習(xí)】數(shù)值分析（1）—— 任意方程求根

任意方程求根

1.簡介

2.方程

2.1分類和介紹

2.2方程的零點(diǎn)（根、解）

3.求根方法

3.1二分法

3.1.1普通二分法的原理及操作

3.1.2普通二分法準(zhǔn)確度及速度分析

3.2淺談迭代

3.2.1 迭代是什么

3.2.2 不動點(diǎn)的定義

3.2.3 函數(shù)的相似性

3.2.4 迭代收斂速度

3.3不動點(diǎn)迭代法

3.3.1不動點(diǎn)迭代法的原理及操作

3.3.2迭代法的收斂性證明

3.3.3不動點(diǎn)迭代法的收斂速度

迭代的終止

3.4牛頓法

3.4.1普通牛頓法

3.4.2牛頓下山法

3.4.3簡單牛頓法

3.4.4 改進(jìn)牛頓法——弦截法

3.4.5 牛頓法的收斂性

收斂迭代的加速與精度

收斂的精度

y還是x？

敏感性

普通迭代加速

Atiken加速法

練習(xí)題

Reference

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容

色偷偷精品伊人,欧洲久久精品,欧美综合婷婷骚逼,国产AV主播,国产最新探花在线,九色在线视频一区,伊人大交九 欧美,1769亚洲,黄色成人av

【機(jī)器學(xué)習(xí)】數(shù)值分析（1）—— 任意方程求根

任意方程求根

1.簡介

2.方程

2.1分類和介紹

2.2方程的零點(diǎn)（根、解）

3.求根方法

3.1二分法

3.1.1普通二分法的原理及操作

3.1.2普通二分法準(zhǔn)確度及速度分析

3.2淺談迭代

3.2.1 迭代是什么

3.2.2 不動點(diǎn)的定義

3.2.3 函數(shù)的相似性

3.2.4 迭代收斂速度

3.3不動點(diǎn)迭代法

3.3.1不動點(diǎn)迭代法的原理及操作

3.3.2迭代法的收斂性證明

3.3.3不動點(diǎn)迭代法的收斂速度

迭代的終止

3.4牛頓法

3.4.1普通牛頓法

3.4.2牛頓下山法

3.4.3簡單牛頓法

3.4.4 改進(jìn)牛頓法——弦截法

3.4.5 牛頓法的收斂性

收斂迭代的加速與精度

收斂的精度

y還是x？

敏感性

普通迭代加速

Atiken加速法

練習(xí)題

Reference

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容

色偷偷精品伊人,欧洲久久精品,欧美综合婷婷骚逼,国产AV主播,国产最新探花在线,九色在线视频一区,伊人大交九欧美,1769亚洲,黄色成人av

2.2方程的零點(diǎn)（根、解）

y還是x？