本文由CocoaChina翻譯組成員leon（社區ID）翻譯自raywenderlich
原文：Numerical Algorithms using Playgrounds

---

中學時，沒有什麼比數學圖紙更恐怖的了。

許多問題沒有現成解析方案，還有一些問題雖然可以解決，但需要大量的計算工作。這種情況下，你需要算法來輔助。

希望你沒被算法整吐。萬一你真的吐了，也要看到好的一面：又可以再吃頓午餐了:]

在這個教程中，你將學習到算法的概念，以及如何使用它們來解決難以分析的問題。通過playground，很容易使解決方案可視化，易於查看。

即使你不是個數學愛好者，對物理或計算機科學也不大感興趣，你仍會在這個教程中找到一些有意思的東西。你只需要有一些微積分和基礎物理學的知識。

你將學到如何使用數值算法解決兩個困難問題。但是為了講的更清楚點，這兩個問題也可通過解析法解決。算法在解析法無法工作時更為理想，即便如此，通過比較兩種方法，也能更加容易的理解它們的本質。

數值算法是什麼？
簡單說來，數值算法是一類解決數學問題的方法，它們不依賴於閉式解析解。

問題來了，什麼是閉式解析解？

若有一個公式，我們可以使用已知數，通過有限的操作步驟，可以獲得一個准確的結果，即稱之為閉式解析解。

再簡單點來理解一下：如果你可以使用代數的方式，找到一個表達式來解決一個未知量問題，代入某些已知數即可得到結果，就說明你已經得到了一個閉式解析方法。

何時使用數值算法？

許多問題沒有現成解析方案，還有一些問題雖然可以解決，但需要大量的計算工作。這種情況下，你需要算法來輔助。

例如，假定你需要編寫一個物理引擎，用來計算許多物體在有限時間內的行為。這時你就可以使用數值算法來更快地得到結果。

這樣做也有副作用，更快的計算速度意味著結果不會十分精確。然而，大多數情況下，這樣的結果已經夠用了。

天氣預報就得益於數值計算。天氣變化的速度、影響因素的數量都是十分驚人的。這是一個非常復雜的系統，只有數值模擬才能完成預知未來的任務。

可能正是因為缺乏這些算法，你的iPhone總是告訴你要下雨了，而外面還是艷陽高照。

開始

作為熱身活動，我們來玩個游戲，接下來，你將計算出一個給定數字的平方根。這兩個方法都將使用二分法（bisection method）。神奇的是，你可能已經知道了這個方法，但是不知道它的名字。

回想一下，在兒童時代，我們可能玩過這樣的游戲：選中100以內的一個數字，另外一個人來猜。你只會提示他猜的數字大了或者小了。

游戲開始。小明告訴小強開始猜，小強說我猜1，小明說小了，小強又猜2，小明說小了。小強接下來再猜3,4...最後選中了5，終於對了。

5步就猜中了，不錯。但是如果小明選的是78，猜中就需要花點時間。

這個游戲如果使用二分法（bisection method）來玩的話，猜中的速度會快很多。

二分法

我們知道數字介於1和100之間，因此除了一個一個的猜，或者隨便瞎猜外。我們把數字分為兩個區間：a = [1,50], b = [51, 100]。

接下來我們判斷數字是介於哪一個區間，這很簡單，只用問數字是不是50。如果數字小於50，那麼區間b就不用考慮了。接下來我們繼續把a區間分成兩半，再重復這個步驟。

例如：

數字是60.區間為：a = [1,50], b = [51, 100]。

第一步，小強說50（也即是a區間的上限），小明說小了。這時小強就知道了數字肯定在b區間上。

第二步，分解b區間為： c=[51,75]，d=[76,100]。還是選擇c區間的上限75，小強告訴他大了。這就意味著數字肯定在c區間上。

繼續分解。。。

使用這個方法，7次嘗試即可得到正確答案，一個一個試則需要60次。

1. 50 -> 小了

2. 75 -> 大了

3. 62 -> 大了

4. 56 -> 小了

5. 59 -> 小了

6. 61 -> 大了

7. 60 -> 對了！

計算x的平方根，過程也類似。數字x的平方根介於0和x之間。也可以記做：(0,x]。如果數字x大於或等於1，可以記做：[1, x]。

分解區間，得到a=(0, x/2]，b=(x/2, x]。

如果數字x是9，區間是[1, 9]，分解後的區間為a=[1, 5]，b=(5, 9]，中間值m為(1+9)/2 = 5。

下一步，檢查m * m - x，是否大於期望的精度。在這裡，我們檢查m * m 大於或小於等於x。如果大於，我們使用(0, x/2]區間繼續檢查，否則，使用(x/2, x]區間。

看一下實際步驟，初始化m=5, 期望准確值為0.1:

1. 計算m * m - x: 5 * 5 - 9 = 11

2. 檢查結果是否小於等於期望值：25 - 11 <= 0.1?

3. 顯然不滿足，繼續下一個區間：m * m是否大於9？是。接下來使用區間[1, 5]，測試值為(1+5)/2=3。

4. 計算m * m - x：３* 3 - 9 = 0

5. 查查是否滿足期望值：9 - 9 <= 0.1?

6. 搞定。

備注：你可能想知道括號是什麼意思？區間有固定的格式（下限和上限）。不同的括號代表不同的區間邊界。圓括號表示邊界不在區間范圍內（即開區間），而方括號表示邊界在區間范圍內（閉區間）。如(0, a] 包含a而不包含0.在上面的例子中：區間 (0, a/2]包含a/2而不包含0；區間(a/2, a]表示所有大於a/2, 小於a的數。

在Playground上使用二分法

現在，是時候應用學到的理論了。我們來自己實現二分算法。創建一個新的playground文件，添加如下的代碼：

func bisection(x: Double) -> Double {
//1
    var lower = 1.0
    //2
    var upper = x
    //3
    var m = (lower + upper) / 2
    var epsilon = 1e-10
    //4
    while (fabs(m * m - x) > epsilon) {
    //5
        m = (lower + upper) / 2
        if m * m > x {
            upper = m
        } else {
            lower = m
        }
    }
 
    return m
}

看一下代碼各部分的含義：

1. 設置區間下限為1

2. 設置區間上限為x

3. 找到中間值，並定義期望精確度

4. 檢查操作是否能滿足精確度

5. 如果不滿足，找到新的中間值，定義新的區間，繼續查找。

添加如下代碼來測試該函數：

let bis = bisection(2.5)

在 m = (lower + upper) / 2這一行的右邊，可以看到代碼執行了35次，意味著找到結果需要35步。

馬上我們就能看到playground一個可愛功能帶來的好處：數值的完整歷史都可以查看。

二分法可以成功的算出實際結果的近似值。通過數值的歷史記錄圖，就可以看到數值算法是如何逼近正確的解。

按下option+cmd+enter打開輔助編輯器，點擊m = (lower + upper) / 2代碼行右邊的圓按鈕在輔助編輯器中添加歷史記錄。

你會看到方法一點點的跳轉到正確結果上。

古典數學也搞的定

下一個算法需要追溯到古代，它起源於公元前1750年的古巴比倫，在公元100年前亞歷山大的希羅所著《度量論》（Metrica ）一書中有所描述。這就是為何它被稱作“希羅方法”。可以通過這個 鏈接 了解更多關於希羅的知識。

這個方法使用函數，這裡你想要計算a的平方根。如果你能找到函數曲線的“0值點（zero point）”，這個點上的函數值為0，那麼求出x值就能得到a的平方根。

要完成這項工作，我們首先選擇任意x值作為起始值，計算該值對應點的tangent值（函數的切線），然後查看tangent線上的0點（即該直線和x軸的交點）。然後我們使用這個0點再次作為起始值，重復多次直到滿足精度。

因為每計算一次tangent值，都會更加逼近真正的0值，

這個過程會逐漸逼近真正的結果。下圖展示了求解a=9時，a的平方根，且起始值為1.

起始點x0 = 1，生成一條紅色的tangent線，接著生成了x1點，又生成了紫色的線；又生成了x2點，連接了藍色的線，最終找到了正確答案。

當古典數學邂逅了Playground

我們有很多古巴比倫人沒有的好東西，比如：playground。在playground上添加如下代碼，並查看：

func heron(x: Double) -> Double {
  //1
  var xOld = 0.0
  var xNew = (x + 1.0) / 2.0
  var epsilon = 1e-10
 
  //2
  while (fabs(xNew - xOld) > epsilon) {
    //3
    xOld = xNew
    xNew = (xOld + x / xOld) / 2
  }
 
  return xNew
}

這段代碼做了什麼？

1. 創建一個變量來存儲當前結果。xOld是上次計算的結果，而xNew是實際結果。

• 找到初始化xNew的方法是一個良好的起點

• epsilon是期望的精確度

• 這個例子中，我們計算平方根精確度為小數點後10位。

2. 使用while循環查找是否達到期望的精確度。

3. 如果達不到精確度，設置xNew為xOld，繼續下一次迭代。

使用如下代碼驗證該函數的作用：

let her = heron(2.5)

希羅方法只需要5次迭代即可找到正確結果。

在代碼行xNew = (xOld + x/xOld)/2右邊點擊圓角button，添加值歷史，就能看到第一次迭代就能找到一個不錯的近似值。

模擬諧振子運動

在這個章節中，我們來看看如何使用數值積分算法來模擬簡諧系統運動——一種基本的動態系統。

這個系統可以描述很多現象，比如鐘擺的擺動、彈簧的振動。特別說來，它可以描述某段時間內物體移動了一定偏移量的場景。

這個例子中，假定有一個有質量的物體連接在彈簧上。為了簡化問題，我們忽略掉阻力和重力。因此唯一的作用力就是彈簧的彈力，它將物體拉回到原來的位置。

在這樣的假定下，你只需要使用兩個力學定律來處理問題：

• 牛頓第二運動定律, 它描述了物體上的作用力和加速度的關系。

• 胡克定律，它描述了彈力和物體偏移量之間的比例關系。這裡k是彈力系數而x是物體的偏移量。

簡振方程

因為彈力是物體上唯一的作用力，我們將上面的兩個方程式整理：

簡化後寫作：

也被記作，也即是共振頻率的平方。

更精確的方程式如下：

其中A是振幅，這裡即是物體的偏移量，被稱為相位差。這兩個值初始化時都被設置為定值。

如果設置時間t=0，則，且，你可以簡單的計算出振幅和相位差：

來看一個例子，我們有一個質量為2kg的物體連接在彈簧上，彈簧的彈力系數為k=196N/m。在初始時間t=0時，彈簧移動了0.1米。我們可以通過公式計算振幅、相差和共振頻率。

在Playground上實驗一下：

使用該公式計算任意時間點對應的值之前，我們需要添加一些代碼。

回到playground文件，在最後添加如下代碼：

//1
typealias Solver = (Double, Double, Double, Double, Double) -> Void
 
//2
struct HarmonicOscillator {
  var kSpring = 0.0
  var mass = 0.0
  var phase = 0.0
  var amplitude = 0.0
  var deltaT = 0.0
 
  init(kSpring: Double, mass: Double, phase: Double, amplitude: Double, deltaT: Double) {
    self.kSpring = kSpring
    self.mass = mass
    self.phase = phase
    self.amplitude = amplitude
    self.deltaT = deltaT
  }
 
  //3
  func solveUsingSolver(solver: Solver) {
    solver(kSpring, mass, phase, amplitude, deltaT)
  }
}

這些代碼塊中做了什麼？

1. 定義一個函數類型的的typealias，函數有5個Double類型的參數，返回為空。

2. 創建一個struct來定義諧振子。

3. 這個函數只是簡單的創建用來解振動問題的Clousure。（並無真實的計算代碼）

再來看一下精確方案

代碼如下：

func solveExact(amplitude: Double, phase: Double, kSpring: Double, mass: Double, t: Double) {
  var x = 0.0
  //1
  let omega = sqrt(kSpring / mass)
 
  var i = 0.0
 
  while i < 100.0 {
//2
    x =  amplitude * sin(omega * i + phase)
    i += t
  }
}

這個方法包含了所有解決運動方程需要的參數。

1. 計算共振頻率

2. 在while循環中計算物體當前的位置，i用作下次計算的自增變量。

添加如下測試代碼：

let osci = HarmonicOscillator(kSpring: 0.5, mass: 10, phase: 10, amplitude: 50, deltaT: 0.1)
osci.solveUsingSolver(solveExact)

這個方案中的方法有點奇怪：它有參數傳入，但不返回數據，也沒有顯示任何東西。

這樣做有什麼好處？

使用這個函數的目的在於，在while循環中，可以算出振動過程中具體的動態值。在Playground中，可以觀察這些動態值的歷史記錄。

在x = amplitude * sin(omega * i + phase) 處添加值歷史記錄，我們就能看到運動的軌跡。

既然第一個精確算法已經驗證通過，下面我們看一下數值計算方案。

歐拉方法

歐拉方法是用來求數值積分的一種方法。該方法1768年再Leonhard Euler所著Institutiones Calculi Integralis《積分學原理》一書中首次提出。要查看該方法的詳情，請參考：http://en.wikipedia.org/wiki/Euler_method

歐拉方法的核心思想是通過使用折線逼近曲線。

具體方法是計算一個給定點的斜率，然後繪制一條同樣斜率的折線。在這條折線的末尾，繼續計算斜率，繪制另外一條線。正如你所見，該算法的精確度取決於折線的長度。

你想知道deltaT做了什麼嗎？

該數值算法中需要使用一個步長，這對算法的精確度很重要。大步長導致低精確度，但是執行速度塊。反之，精確度會提高，速度會降低。

deltaT變量就代表了步長(step size)。我們初始化這個值為0.1， 意味著我們計算每0.1秒物體移動的位置。在歐拉算法中，意味著折線在x軸上的投影長度為0.1。

在編寫代碼前，你需要再看一下這個公式：

二階微分方程可以轉化為兩個一階微分方程。可以被寫為和

我們可以用差商來完成轉換：

也會得到：

有了這些等式，我們可以直接實現歐拉方法。

在solveExact方法後添加代碼：

func solveEuler(amplitude: Double, phase: Double, kSpring: Double, mass: Double, t: Double) {
  //1
  var x = amplitude * sin(phase)
  let omega = sqrt(kSpring / mass)
  var i = 0.0
  //2
  var v = amplitude * omega * cos(phase);
  var vold = v
  var xoldEuler = x
 
  while i < 100 {
    //3
    v -= omega * omega  * x * t
    //4
    x += vold * t
    xoldEuler = x
 
    vold = v
    i+=t
    }
}

代碼的內容：

1. 設置當前的位置，和omega的值。

2. 計算當前速度。

3. 在while循環中，通過一階函數計算出新的速度。

4. 使用速度計算新的位置，在結束處，使用步長t增加i的值。

現在，添加以下代碼測試：

osci.solveUsingSolver(solveEuler)

在xoldEuler = x位置添加值歷史，並查看。你會看到這個方法顯示一個振幅增加的正弦函數。因為歐拉方法並不是一個精確算法，而這裡過大的步長0.1也導致了計算結果不精確。

以下是另外一個函數圖像，步長為0.01，這樣明顯更好。因此，使用歐拉方法時你要想到，步長越小，結果越好。但是使用折中的步長，執行起來更為簡單。

速度Verlet算法(Velocity Verlet)

最後討論的算法叫做速度Verlet。它和歐拉方法的思路一樣，但是計算新位置的方式有些許差異。

歐拉在計算新位置時，忽略實際的加速度，公式為：,而速度Verlet算法在計算時考慮到了加速度：。這再步長相等的時候，結果更優。

在solveEuler方法後添加新的代碼：

func solveVerlet(amplitude: Double, phase: Double, kSpring: Double, mass: Double, t: Double) {
  //1
  var x = amplitude * sin(phase)
  var xoldVerlet = x
  let omega = sqrt(kSpring / mass)
 
  var v = amplitude * omega * cos(phase)
  var vold = v
  var i = 0.0
 
  while i < 100 {
    //2
    x = xoldVerlet + v * t + 0.5 * omega * omega * t * t
    v -= omega * omega * x * t
    xoldVerlet = x
    i+=t
  }
}

代碼的內容：

1. 循環前的所有代碼和歐拉方法中的一樣。

2. 根據速度計算出新的位置，增加i的值。

在Playground中測試代碼:

osci.solveUsingSolver(solveVerlet)

別忘了在xoldVerlet = x代碼行上添加值歷史記錄：

接下來做什麼？

你可以通過 此鏈接 下載完整工程.

希望你在這場數值世界旅行中獲得樂趣。了解一些算法，即便只是一些古典數學的趣味歷史，都會對你帶來幫助。

繼續探索之前，你可以再改改deltaT，看看有什麼不同的結果。但是要記住playground速度並不快，選擇一個小的步長可能導致一段時間內Xcode無法操作。在我的Macbook上，設置deltaT到0.01就讓我的機器定住了好幾分鐘。

使用你學到的算法，可以解決許多問題，比如說：n-body 問題。

如果你覺得自己的算法基礎不牢固。可以在YouTube上看看MIT的 算法簡介課程。

最後，有時間可以在本站點查看更多關於這個課題的教程。

如果你對於使用數值算法有任何的疑問或建議，可以在下面評論。期望看到你的想法，並討論在app中可以用到的很酷的數學知識。