2.3 數(shù)值分析

數(shù)值分析是計(jì)算數(shù)學(xué)的一個(gè)重要部分，它研究用計(jì)算機(jī)求解各種數(shù)學(xué)問(wèn)題的數(shù)值計(jì)算方法，以及其理論與軟件實(shí)現(xiàn)。用計(jì)算機(jī)求解數(shù)學(xué)問(wèn)題通常經(jīng)歷以下5個(gè)步驟：

實(shí)際問(wèn)題→數(shù)學(xué)模型→數(shù)值計(jì)算方法→程序設(shè)計(jì)→上機(jī)計(jì)算求出結(jié)果。

根據(jù)實(shí)際問(wèn)題建立數(shù)學(xué)模型往往是應(yīng)用數(shù)學(xué)的任務(wù)，計(jì)算數(shù)學(xué)更關(guān)注的是如何給出數(shù)值計(jì)算方法，并根據(jù)計(jì)算方法編制算法程序，從而求得最終的計(jì)算結(jié)果。

計(jì)算機(jī)及科學(xué)技術(shù)的快速發(fā)展使求解各種數(shù)學(xué)問(wèn)題的數(shù)值方法也越來(lái)越多，解決問(wèn)題的速度和效率也得到了很大的提升。數(shù)值分析涉計(jì)數(shù)學(xué)的各個(gè)分支，所包含的內(nèi)容十分廣泛，本節(jié)只介紹基本的數(shù)值計(jì)算方法及其理論內(nèi)容，包括：插值與數(shù)據(jù)逼近；數(shù)值微分與積分；線性方程組的數(shù)值求解；非線性方程與方程組求解等。它以數(shù)學(xué)問(wèn)題為導(dǎo)向，將理論與計(jì)算相結(jié)合，重點(diǎn)關(guān)注數(shù)學(xué)問(wèn)題的數(shù)值算法及其理論，是一門依賴于計(jì)算機(jī)技術(shù)的實(shí)用性課程。

2.3.1 插值法與數(shù)值逼近

插值法有著悠久的歷史，它來(lái)自古人的生產(chǎn)實(shí)踐活動(dòng)，可以追溯到一千多年前的隋唐時(shí)期，人們利用二次插值法制定了歷法，進(jìn)行天文計(jì)算。17世紀(jì)，微積分的產(chǎn)生給插值理論的建立提供了條件。近半個(gè)世紀(jì)以來(lái)，由于計(jì)算機(jī)的推廣，插值法逐步延伸到造船、航空航天、機(jī)械加工等實(shí)際工程問(wèn)題中，獲得了更廣泛的應(yīng)用。常見(jiàn)的插值法有拉格朗日插值法、牛頓的等距節(jié)點(diǎn)插值法，以及均差插值公式、埃爾米特插值法、分段低次插值法、三次樣條插值法等。

（1）拉格朗日插值。

若n次多項(xiàng)式lj(x)(j=0,1,…,n)在n+1個(gè)節(jié)點(diǎn)x0＜x1＜…＜xn上滿足式（2.1）：

我們稱這n+1個(gè)n次多項(xiàng)式l0(x),l1(x),…,ln(x)為節(jié)點(diǎn)x0,x1,…,xn上的n次插值基函數(shù)，可以表示為式（2.2）：

則n次插值多項(xiàng)式可以表示為式（2.3）：

由lk(x)的定義可得式（2.4）：

此時(shí)的插值多項(xiàng)式Ln(x)稱為拉格朗日插值多項(xiàng)式。

（2）均差與牛頓插值多項(xiàng)式。

利用插值基函數(shù)可以得到拉格朗日插值多項(xiàng)式，但當(dāng)插值節(jié)點(diǎn)數(shù)發(fā)生增減變化時(shí)，需要重新進(jìn)行計(jì)算。為了計(jì)算方便，可重新設(shè)計(jì)一種逐次生成插值多項(xiàng)式的方法。

已知f在插值點(diǎn)xi(i=0,1,…,n)上的值為f(xi)(i=0,1,…,n)，要求n次插值多項(xiàng)式Pn(x)滿足條件式（2.5）：

則Pn(x)可以表示為式（2.6）：

式中，a0,a1,…,an為待定系數(shù)。與拉格朗日插值不同的是，這里的Pn(x)是由基函數(shù){1,x-x0,…,(x-x0)…(x-xn-1)}逐次遞推得到的。為了確定a0,a1,…,an的表達(dá)式，需要引入均差的定義。

一般地，稱式

為f(x)的k階均差。

借助均差的定義，可以將x看作[a,b]上的一點(diǎn)，可得式（2.8）。

式中，

其系數(shù)ak滿足式（2.9）：

我們稱Pn(x)為牛頓均差插值多項(xiàng)式，它比拉格朗日插值計(jì)算量小，更便于程序設(shè)計(jì)。

在數(shù)值計(jì)算中經(jīng)常要計(jì)算函數(shù)值，當(dāng)函數(shù)只在有限點(diǎn)集上給定函數(shù)值，要在包含該點(diǎn)集的區(qū)間上用公式給出函數(shù)的簡(jiǎn)單表達(dá)式，這就涉及在區(qū)間[a,b]上用簡(jiǎn)單函數(shù)逼近已知復(fù)雜函數(shù)的問(wèn)題，這就是函數(shù)逼近問(wèn)題。

2.3.2 數(shù)值積分與數(shù)值微分

積分和微分是用極限來(lái)定義的兩種分析運(yùn)算的，處理數(shù)值分析和數(shù)值微分的基本方法是逼近法：設(shè)構(gòu)造某個(gè)簡(jiǎn)單函數(shù)P(x)近似f(x)，然后對(duì)P(x)求積（求導(dǎo)）得到f(x)的積分（導(dǎo)數(shù)）近似值。插值求積公式分為牛頓-柯特斯公式和高斯公式兩類，實(shí)際計(jì)算宜采用復(fù)合求積的方法。高斯公式精度高、計(jì)算穩(wěn)定，但節(jié)點(diǎn)選取較困難。帶權(quán)高斯求積方法能把復(fù)雜積分簡(jiǎn)單化，還可以直接計(jì)算奇異積分。基于理查森外推的龍貝格求積方法由于計(jì)算程序簡(jiǎn)單、精度較高，是一種在計(jì)算機(jī)上求積的有效算法。在數(shù)值微分中也有相似的算法，但由于計(jì)算的不穩(wěn)定性，在數(shù)值微分中，步長(zhǎng)的選取很重要。

2.3.3 解線性方程的直接方法與迭代法

關(guān)于線性方程組的數(shù)值解法一般有直接法和迭代法兩類。直接法就是經(jīng)過(guò)有限步算術(shù)運(yùn)算，求得方程組的精確解，但由于實(shí)際計(jì)算中存在誤差，這種方法有一定的局限性，只能求得線性方程組的近似解。迭代法是利用某種極限過(guò)程去逐步逼近線性方程組的精確解。迭代法要求計(jì)算機(jī)的存儲(chǔ)單元較少、程序設(shè)計(jì)簡(jiǎn)單、原始系數(shù)矩陣在計(jì)算過(guò)程中始終不變，但存在收斂性及收斂速度的問(wèn)題。迭代法是解決大型稀疏矩陣方程組的重要方法。

1.直接方法

（1）高斯消去法。

設(shè)有線性方程組：

或?qū)憺榫仃囆问剑?/p>

簡(jiǎn)記為Ax=b。

高斯消去法解方程組的基本思想，是用逐次消去未知數(shù)的方法把原方程組Ax=b化為與其等價(jià)的三角形方程組，而求解三角形方程組可用回代的方法求解。

下面我們討論如何用算法來(lái)實(shí)現(xiàn)高斯消去的過(guò)程：

設(shè)A∈Rm×n(m＞1),s=m:n(m-1,n)，如果，利用高斯方法將A化為上梯形，且A(k)覆蓋A，乘數(shù)mik覆蓋aik。

對(duì)于k=1,2,…,s：

① 如果akk=0，則停止計(jì)算。

② 對(duì)于i=k+1,…,m

i.aik←mik=aik/akk

ii.對(duì)于j=k+1,…,n

aij←aij-mik*akj

由此可見(jiàn)，算法的第k步需要做m-k次除法和(m-k)(n-k)次乘法運(yùn)算，因此，本算法大約需要s3/3-(m+n)s2/2+mns次乘法運(yùn)算。當(dāng)m=n時(shí)，總計(jì)大約需要n3/3次乘法運(yùn)算。

（2）高斯主元素消去法

由高斯主元素消去法可知，在消元過(guò)程中可能出現(xiàn)的情況，即無(wú)法進(jìn)行消去法，即使主元素，但很小時(shí)，用其作除數(shù)會(huì)導(dǎo)致其他元素?cái)?shù)量級(jí)的嚴(yán)重增長(zhǎng)和舍入誤差的擴(kuò)散，最后，計(jì)算結(jié)果不準(zhǔn)確。

目前，主要使用的是列主元素消去法，并假定式（2.11）的A∈Rn×n為非奇異的。

設(shè)方程組（2.11）的增廣矩陣為

首先，在A的第1列中選取絕對(duì)值最大的元素作為主元素，如；然后，交換B的第1行與第i1行，經(jīng)第1次消元計(jì)算得。

重復(fù)上述過(guò)程，設(shè)完成第k-1步的選主元素，交換兩行及消元計(jì)算，(A b)化為

式中，A(k)的元素仍記為aij，b(k)的元素仍記為bi。

第k步選主元素（在A(k)右下角方陣的第1列內(nèi)選），即確定ik，

交換(A(k)b(k))第k行與ik(k=1,2,…,n-1)行的元素，再進(jìn)行消元計(jì)算，最后將原線性方程組化為：

回代求解

下面，我們討論如何用算法來(lái)實(shí)現(xiàn)列主元素消去法的過(guò)程。

設(shè)Ax=b，本算法用A具有行交換的列主元素消去法。消元結(jié)果沖掉A，乘數(shù)mij沖掉aij，計(jì)算解x沖掉常數(shù)項(xiàng)b，行列式存放在det中。

① det←1。

② 對(duì)于k=1,2,…,n-1，

i.按列選取主元素。

ii.如果aik,k=0，則停止計(jì)算(det(A)=0)。

iii.如果ik=k，則轉(zhuǎn)iv，

akj?aik,j(j=k,k+1,…,n)。

換行：bk?bik，

det←det。

iv.消元計(jì)算。

對(duì)于i=k+1,…,n，

a.aik←mik=aik/akk。

b.對(duì)于j=k+1,…,n，

aij←aij-mik·akj。

c.bi←bi-mik·bk。

v.det←akk·det。

③ 如果ann=0，則停止計(jì)算[det(A)=0]。

④ 回代求解。

i.bn←bn/ann。

ii.對(duì)于i=n-1,…,2,1，

⑤ det←ann·det。

（2）迭代法。對(duì)于由工程技術(shù)中產(chǎn)生的大型稀疏矩陣方程組（A的階數(shù)n很大，但零元素較多），利用迭代法求解線性方程組Ax=b更合適。下面主要介紹迭代法中的雅克比迭代法和高斯-塞德?tīng)柕ā?/p>

1）雅克比迭代法。

解Ax=b的雅克比迭代法的計(jì)算公式為

2）高斯-塞德?tīng)柕ā?/p>

解Ax=b的雅克比迭代法的計(jì)算公式為

下面我們討論如何用算法來(lái)實(shí)現(xiàn)高斯-塞德?tīng)柕ǖ倪^(guò)程。

設(shè)Ax=b，其中A∈Rn×n為非奇異矩陣，且aii≠0(i=1,2,…,n)，本算法用高斯-塞德?tīng)柕ń?span id="sbvrg7t" class="italic">Ax=b，數(shù)組x(n)開(kāi)始存放x(0)，后存放x(k)，N0為最大迭代次數(shù)。

① xi←0.0(i=1,2,…,n)。

② 對(duì)于k=1,2,…,N0，i=1,2,…,n，

迭代一次，這個(gè)算法需要的運(yùn)算次數(shù)最多與矩陣A的非零元素的個(gè)數(shù)一樣。