数值分析







巴比伦泥板 YBC 7289
(c. 1800–1600 BCE),泥板下有根號2的六十进制近似值,精確值大約接近十進位的小數下第6位泥板,1 + 24/60 + 51/602 + 10/603 = 1.41421296...[1]


数值分析英语:numerical analysis),是指在数学分析(区别于离散数学)问题中,对使用数值近似(相对于一般化的符号运算)算法的研究。


巴比伦泥板YBC 7289是关于数值分析的最早数学作品之一,它给出了 2{displaystyle {sqrt {2}}}{sqrt {2}} 在六十进制下的一个数值逼近,2{displaystyle {sqrt {2}}}{sqrt {2}}是一個邊長為1的正方形的對角線,在西元前1800年巴比倫人也已在巴比倫泥板上計算勾股數(3,4,5){displaystyle (3,4,5)}{displaystyle (3,4,5)},即直角三角形的三邊長比。


数值分析延續了實務上數學計算的傳統。巴比倫人利用巴比伦泥板計算2{displaystyle {sqrt {2}}}{sqrt {2}}的近似值,而不是精確值。在許多實務的問題中,精確值往往無法求得,或是無法用有理數表示(如2{displaystyle {sqrt {2}}}{sqrt {2}})。数值分析的目的不在求出正確的答案,而是在其誤差在一合理範圍的條件下找到近似解。


在所有工程及科學的領域中都會用到数值分析。像天體力學研究中會用到常微分方程,最優化會用在资产组合管理中,數值線性代數是資料分析中重要的一部份,而隨機微分方程及馬可夫鏈是在醫藥或生物學中生物細胞模擬的基礎。


在電腦發明之前,数值分析主要是依靠大型的函數表及人工的內插法,但在二十世紀中被電腦的計算所取代。不過電腦的內插演算法仍然是数值分析軟體中重要的一部份。




目录






  • 1 簡介


    • 1.1 直接法和迭代法


    • 1.2 離散化




  • 2 誤差的產生及傳播


    • 2.1 捨入誤差


    • 2.2 截尾及離散化誤差


    • 2.3 數值穩定性及良置問題




  • 3 领域研究


    • 3.1 函數求值


    • 3.2 內插法、外推法、曲線擬合及回歸


    • 3.3 求解方程及方程組


    • 3.4 求解特徵值或奇異值問題


    • 3.5 最优化


    • 3.6 積分計算


    • 3.7 微分方程




  • 4 軟體


  • 5 註解


  • 6 参考文献


  • 7 外部链接


  • 8 参阅





簡介


数值分析的目的是設計及分析一些計算的方式,可針對一些問題得到近似但夠精確的結果。以下是一些會用利用数值分析處理的問題:




  • 數值天氣預報中會用到許多先進的数值分析方法。

  • 計算太空船的軌跡需要求出常微分方程的數值解。

  • 汽車公司會利用電腦模擬汽車撞擊來提昇汽車受到撞擊時的安全性。電腦的模擬會需要求出偏微分方程的數值解。


  • 对冲基金會利用各種数值分析的工具來計算股票的市值及其變異程度。

  • 航空公司會利用複雜的最佳化演算法決定票價、飛機、人員分配及用油量。此領域也稱為作業研究。

  • 保險公司會利用数值軟體進行精算分析。



直接法和迭代法



直接法和迭代法


考慮以下問題


3x3 + 4 = 28

要求解未知數x






















直接法
3x3 + 4 = 28.
減 4 3x3 = 24.
除 3
x3 = 8.
開立方
x = 2.

若是用迭代法,可用迭代法求解f(x) = 3x3 − 24=0,初值為a = 0, b = 3, f(a) = −24, f(b) = 57。




































迭代法
a b 中點
f(中點)
0 3 1.5 −13.875
1.5 3 2.25 10.17...
1.5 2.25 1.875 −4.22...
1.875 2.25 2.0625 2.32...

計算到目前為止,問題的解是界於1.875及2.0625之間,若繼續往下算,可以得到更精確的答案。



直接法利用固定次數的步驟求出問題的解。這些方式包括求解线性方程组的高斯消去法及QR演算法英语QR algorithm,求解線性規劃的单纯形法等。若利用無限精度算術的計算方式,有些問題可以得到其精確的解。不過有些問題不存在解析解(如五次方程),也就無法用直接法求解。在電腦中會使用浮點數進行運算,在假設運算方式稳定的前提下,所求得的結果可以視為是精確解的近似值。


迭代法是通過從一個初始估計出發尋找一系列近似解來解決問題的數學過程。和直接法不同,用迭代法求解問題時,其步驟沒有固定的次數,而且只能求得問題的近似解,所找到的一系列近似解會收敛到問題的精確解。會利用審斂法來判別所得到的近似解是否會收斂。一般而言,即使使用無限精度算術的計算方式,迭代法也無法在有限次數內得到問題的精確解。


在數值分析中用到迭代法的情形會比直接法要多。例如像牛頓法、二分法、雅可比法、廣義最小殘量方法(GMRES)及共軛梯度法等。在計算矩陣代數中,大型的問題一般會需要用迭代法來求解。



離散化


許多時候需要將連續模型的問題轉換為一個離散形式的問題,而離散形式的解可以近似原來的連續模型的解,此轉換過程稱為離散化。例如求一個函數的積分是一個連續模型的問題,也就是求一曲線以下的面積若將其離散化變成數值積分,就變成將上述面積用許多較簡單的形狀(如長方形、梯形)近似,因此只要求出這些形狀的面積再相加即可。


例如在二小時的賽車比賽中,記錄了三個不同時間點的賽車速度,如下表















時間
0:20 1:00 1:40
km/h
140 150 180

利用離散化的方式,可以假設賽車在0:00到0:40之間的速度、0:40到1:20之間的速度及1:20到2:00之間的速度分別為三個定值,因此前40分鐘的總位移可近似為(2/3h × 140 km/h) = 93.3 公里。可依此方式近似二小時內的總位移為93.3 公里 + 100 公里 + 120 公里 = 313.3 公里。位移是速度的積分,而上述的作法是用黎曼和進行數值積分的一個例子。




誤差的產生及傳播


誤差是数值分析的重要主題之一。誤差的形成可分為幾種不同的原因。



捨入誤差


當進行數值分析的設備只能用有限位數來表示一個實數時,就會出現捨入誤差(Round-off error),例如用可顯示十位數字的計算器計算1/3,所得到的結果0.333333333,和實際數值的誤差就是捨入誤差。即使進行數值分析的設備用浮點數來表示實數,仍無法完全避免捨入誤差的問題。



截尾及離散化誤差


若迭代法的数值分析算到某一程度就中止計算,或是使用一些近似的數學程序,程序所得結果和精準解不同,就會出現截尾(Truncation)誤差。將問題離散化後,由於離散化問題的解不會和原問題的解完全一様,因此會出現離散化誤差英语discretization error。例如用迭代法計算3x3+4=28{displaystyle 3x^{3}+4=28}3x^{3}+4=28的解,在計算幾次後我們認為其解為1.99,就會有0.01的截尾誤差。


一旦有了誤差,誤差就會藉著計算繼續的擴散。例如一個計算機中的加法是不準的,則a+b+c+d+e的計算也一定不準。例如剛剛計算3x3+4=28{displaystyle 3x^{3}+4=28}3x^{3}+4=28的解為1.99,若後續的運算需要用到3x3+4=28{displaystyle 3x^{3}+4=28}3x^{3}+4=28的解,用1.99代入所得的結果也會不準。


當用近似的方式處理數學式時就會出現截尾誤差。以積分為例,完全精準的積分需要求出曲線下方無限個梯形的面積和,但用在數值分析中會用有限個梯形的面積和來近似無限個梯形的面積和,此時就會出現截尾誤差。若要對一個函數進行微分,其微分量需要趨近於0,但實務上只能選擇很小的微分量。



數值穩定性及良置問題



非良置問題:考慮一函數f(x) = 1/(x − 1),f(1.1) = 10,f(1.001) = 1000。當x只改變小於0.1的數值,f(x)的變化將近1000。因此在x = 1的附近計算f(x)是一個非良置的問題。


良置問題:相反的,函數f(x)=x{displaystyle f(x)={sqrt {x}}}f(x)={sqrt  {x}}x不接近0時,其值的計算就是一個良置的問題。



數值穩定性是數值分析中一個重要的主題。若一演算法中不論什麼原因產生了誤差,此誤差不會在運算中明顯增加,此演算法為數值穩定的演算法。若問題為良置(well-conditioned)的,就會符合上述的特性,也就是問題數據微小的變化只會造成其解的微小變化。相反的,若問題數據微小的變化會造成其解的巨大變化,會稱問題為非良置或病態(ill-conditioned)。


原始問題及求解問題演算法都可以分為良置及非良置,任何的組合都是允許的。


一個求解良置問題的演算法可能是數值穩定的,也可能是數值不穩定的。數值分析的重點就是找到適定性問題的數值穩定演算法。例如,計算2的平方根(大約是1.41421)本身是一個適定性問題。許多求解的演算法都是從一個初始的近似值x1開始去求解,例如x1=1.4,再繼續計算x2x3等。巴比倫法就是一個具有此特性的演算法。另一個方法,先稱之為X方法,演算法為xk + 1 = (xk2−2)2 + xk[註 1]。以下分別用初始值 x1 = 1.4及x1 = 1.42,用二種方式進行幾次迭代。







































巴比倫法
巴比倫法
X方法
X方法

x1 = 1.4

x1 = 1.42

x1 = 1.4

x1 = 1.42

x2 = 1.4142857...

x2 = 1.41422535...

x2 = 1.4016

x2 = 1.42026896

x3 = 1.414213564...

x3 = 1.41421356242...

x3 = 1.4028614...

x3 = 1.42056...


...
...



x1000000 = 1.41421...

x28 = 7280.2284...

可觀察到不論初始值多少,巴比倫法都可以快速的收斂,但X方法在初始值為1.4時收斂的很慢,在初始值為1.42時X方法會發散。因此巴比倫法是數值穩定的方法,而X方法是數值不穩定的方法。



领域研究


數值分析依其待求解的問題不同,分為不同的領域。



內插法:假設一點鐘的氣溫為20度,三點鐘時為14度,可以用線性內插法推測一點半及二點鐘時的氣溫分別是18.5度及17度。


外推法:假設某國家國內生產總值平均每年成長百分之五,去年國內生產總值為一百萬元,可推測今年的國內生產總值為一百零五萬元。


A line through 20 points

回歸分析:給定幾個二維座標上的點,回歸分析就是設法找到一條最接近這些點的直線。


每杯飲料要多少錢呢?

最佳化:有一個賣飲料的小販,若每杯飲料100元,每天可以賣197杯飲料,若飲料單價增加1元,每天就會少賣1杯飲料。飲料定價為148.5元時,其每天的收入為最大值。不過由於飲料單價需為正整數,因此飲料定價可定為149元,對應每天的收入為22,052元。


圖中藍色的是風的方向,黑色的是實際軌跡,紅色的是欧拉方法所得的結果

微分方程:假設在一房間中的不同位置放置一百個風扇,然後在房間中放置一根羽毛,羽毛會依房間中氣流而移動,而房間中的氣流可能相當複雜。不過每一秒量測一次羽毛附近空氣的速度,假設羽毛下一秒是等速的直線運動,即可求得下一秒時羽毛的位置,再量測當時羽毛附近空氣的速度,……。這種方法稱為欧拉方法,常使用在常微分方程的數值分析。




函數求值


数值分析中最簡單的問題就是求出函數在某一特定數值下的值。最直覺的方法是將數值代入函數中計算,不過有時此方式的效率不佳。像針對多項式函數的求值,較有效率的方式是秦九韶算法,可以減少乘法及加法的次數。若是使用浮点数,很重要的是是估計及控制捨入誤差。



內插法、外推法、曲線擬合及回歸


內插法求解以下的問題:有一未知函數在一些特定位置下的值,求未知函數在已知數值的點之間某一點的值。


外推法類似內插法,但需要知道數值的點是在其他已知數值點的範圍以外。一般而言外推法的誤差會大於內插法。


曲線擬合是在已知一些數據的條件下,找到一條曲線完全符合現有的數據,數據可能是一些特定位置及其對應的值,也可能是其他資料,例如角度或曲率等。


回歸分析類似曲線擬合,也是根據一些特定位置及其對應的值,要找到對應曲線。但回歸分析考慮到數據可能有誤差,因此所得的的曲線不需要和數據完全符合。一般會使用最小方差法來進行回歸分析。



求解方程及方程組


另一种常見的問題是求特定方程式的解。首先會依方程式是否線性來區分,例如方程式 2x+5=3{displaystyle 2x+5=3}2x+5=3是線性方程式,而2x2+5=3{displaystyle 2x^{2}+5=3}2x^{2}+5=3是非線性方程式。


此領域許多的研究都和求解線性方程組有關。直接法是線性方程組的係數以矩陣來表示,再利用矩陣分解的方式求解,這些方法包括高斯消去法、LU分解,對於對稱矩陣(或埃爾米特矩陣)及正定矩陣可以用喬萊斯基分解英语Cholesky decomposition,非方陣的矩陣則可以用QR分解。迭代法包括有雅可比法、高斯–塞德迭代法、逐次超松驰法英语successive over-relaxation(SOR)及共轭梯度法,一般會用在大型的線性方程組中。


求根演算法是要解一非線性方程,其名稱是因為函數的根就是使其值為零的點。若函數本身可微且其導數是已知的,可以用牛頓法求解,其他的方法包括二分法、割線法等。線性化英语Linearization則是另一種求解非線性方程的方法。



求解特徵值或奇異值問題


許多重要的問題可以用奇異值分解或特徵分解來表示。例如有些图像压缩演算法[2]就是以奇異值分解為基礎。統計學中對應的工具稱為主成分分析。



最优化



最优化問題的目的是要找到使特定目標函數有最大值(或最小值)的點,一般而言這個點需符合一些約束。


依目標函數及約束條件的不同,最佳化又可以再細分:例如線性規劃處理目標函數及約束條件均為線性的情形,常用單純形法來求解。若目標函數及約束條件其中有一項為非線性,就是非線性規劃的範圍。


有約束條件的問題可以利用拉格朗日乘数轉換為沒有約束條件的問題。



積分計算



數值積分的目的是在求一定積分的值。一般常用牛頓-寇次公式,包括辛普森積分法、高斯求積等。上述方式是利用分治法來處理積分問題,也就是將大範圍的積分切割成許多小範圍的積分,再進行計算。不過在高維度時,上述作法可能會因為要作許多的計算而變得不實用(也就是維數之咒所描述的情形),此時可以採用蒙地卡羅方法或半蒙地卡羅方法。(可參照蒙地卡羅積分英语Monte Carlo integration,或是適用於高維度的稀疏网格英语sparse grid法。)



微分方程


数值分析也會用近似的方式計算微分方程的解,包括常微分方程及偏微分方程。


常微分方程往往會使用迭代法,已知曲線的一點,設法算出其斜率,找到下一點,再推出下一點的資料。歐拉方法是其中最簡單的方式,較常使用的是龍格-庫塔法。


偏微分方程的数值分析解法一般都會先將問題離散化,轉換成有限元素的次空間。可以透過有限元素法、有限差分法及有限體積法英语finite volume method,這些方法可將偏微分方程轉換為代數方程,但其理論論證往往和泛函分析的定理有關。另一種偏微分方程的数值分析解法則是利用離散傅立葉變換或快速傅立葉變換。



軟體


20世紀末,大部份数值分析的演算法都已用許多不同的程式語言實現。Netlib英语Netlib软件库包含了許多数值分析演算法的程式,大部份是Fortran及C語言的程式。商業產品也實現了許多不同的数值分析演算法,包括國際數學及統計程序庫數字型檔及英商纳格资讯英语Numerical Algorithms Group软件库,GNU科学数值库則是自由軟體的数值分析演算法软件库。


数值分析的商用應用程式包括MATLAB、S-PLUS英语S-PLUS、LabVIEW及互動式數據語言(IDL)等,自由軟體或開源軟體的数值分析應用程式則包括FreeMat、Scilab、GNU Octave (類似Matlab)、IT++(C++函式庫連 library)、R語言 (類似S-PLUS)及一些Python的衍生版本。各應用程式的性能有很大的差異:一般而言向量及矩陣的運算都很快,而各應用程式純量運算的速度差異則可能會超過10倍以上[3][4]


許多計算機代數系統的軟體(像Mathematica及Maple)由於使用無限精度算術的計算方式,可以得到比一般軟體更準確的結果。


電子試算表的軟體也可以處理一些簡單的數值分析問題。



註解




  1. ^ 這是一個針對方程式x=(x2−2)2+x=f(x){displaystyle x=(x^{2}-2)^{2}+x=f(x)}x=(x^{2}-2)^{2}+x=f(x)定点迭代法英语fixed point iteration,其解包括2{displaystyle {sqrt {2}}}{sqrt {2}}。由於f(x)≥x{displaystyle f(x)geq x}f(x)geq x,每次迭代會使數值增加,因此x1=1.4<2{displaystyle x_{1}=1.4<{sqrt {2}}}x_{1}=1.4<{sqrt  {2}}會收斂,而x1=1.42>2{displaystyle x_{1}=1.42>{sqrt {2}}}x_{1}=1.42>{sqrt  {2}}會發散。



参考文献





  1. ^ Photograph, illustration, and description of the root(2) tablet from the Yale Babylonian Collection


  2. ^ The Singular Value Decomposition and Its Applications in Image Compression 互联网档案馆的存檔,存档日期2006-10-04.


  3. ^ Speed comparison of various number crunching packages 互联网档案馆的存檔,存档日期2006-10-05.


  4. ^ Comparison of mathematical programs for data analysis Stefan Steinhaus, ScientificWeb.com




外部链接



参阅




  • 算法

  • 計算科學

  • 数值分析主題列表英语List of numerical analysis topics

  • 格拉姆-施密特正交化

  • 數值微分

  • 符号数值计算英语Symbolic-numeric computation

  • 算法分析

  • Numerical Recipes英语Numerical Recipes









Popular posts from this blog

數位音樂下載

When can things happen in Etherscan, such as the picture below?

格利澤436b