新邊界條件下的三次樣條插值函數(shù)

于洋++袁健華++錢江++王美玲

摘要：本文討論了三次樣條插值函數(shù)（cubic spline）邊界條件的更一般情形。將邊界條件的“端點”導數(shù)條件換成“中間插值節(jié)點”的導數(shù)，從而將求樣條函數(shù)的“三對角矩陣”進行了降階并向“上（下）三角矩陣”的轉(zhuǎn)化。在理論上證明了新邊界條件下樣條函數(shù)構造的唯一性，并通過數(shù)值實驗驗證了新邊界條件下的樣條函數(shù)與元函數(shù)有較好地擬合度。

關鍵詞：計算數(shù)學；三次樣條插值；邊界條件；三對角矩陣；降階

中圖分類號：0241.3

文獻標識碼：A

DOI：10.3969/j.issn.1003-6970.2016.02.006

引言

在科學研究與工程設計的諸多問題中，經(jīng)常會要求計算函數(shù)值，或者給出表征某一現(xiàn)象的函數(shù)解析式。當遇到復雜問題時，常常找不到簡單函數(shù)顯式地表出物理量間的關系，這時根據(jù)觀測的若干離散數(shù)據(jù)點構造簡單的函數(shù)代替未知函數(shù)或復雜函數(shù)是高效的解決辦法，即插值法。插值法是函數(shù)逼近的一個重要方法，是計算數(shù)學中常用的技巧，在實際應用問題中也有著廣泛的應用。本文討論的核心內(nèi)容是在處理實際問題中常用的樣條插值函數(shù)。樣條函數(shù)是計算數(shù)學以及計算機輔助幾何設計的重要工具。近年來，隨著計算機技術的不斷進步，樣條函數(shù)的應用研究得到了迅速的發(fā)展和廣泛的應用。

樣條（spline）是繪圖員用來描繪光滑曲線的一種簡單工具。為了把一些指定點（樣點）連接成一條光滑曲線，往往用一條富有彈性的細長材料，如木條（樣條）固定在樣點上，其他地方讓它自由彎曲。這樣連接樣點的曲線稱為樣條曲線。樣條插值函數(shù)就是對樣條曲線進行數(shù)學模擬得到的。下文我們主要針對三次樣條插值函數(shù)進行分析和論證。

設f（x）是區(qū)間[a，b]上的一個二次連續(xù)可微函數(shù)。在區(qū)間[a，b]上給定一組基點：

是二次連續(xù)可微的，Si（x），（i=1，2，…，n）均是不高于三次的多項式或零多項式，且滿足條件

則稱S（x）為函數(shù)f（x）的三次樣條插值函數(shù)，簡稱三次樣條。

記mi=s"（xi），f（xi）=fi，根據(jù)三次樣條的定義知，S（x）的二階導數(shù)S"（x）在每個子區(qū)間[xi，xi+1]上都是線性函數(shù)。所以在[xi，Xi+l]上S（x）=Si（x）的二階導數(shù)

對（4）式連續(xù)積分兩次，并根據(jù)條件（2）和（3）可得

hi-1mi-1+2（hi-1+hi）mi+mi+1=6{f[xi，Xi+1]-f[xi-1，xi]}（5）

其中i=2，…，n，這是一個含有n+l個未知量m1，m2，…，mn+1，而只有n-l個方程的線性方程組。通常要補加兩個約束條件，稱為邊界條件或端點條件，使（5）式成為由n+l個方程構成的方程組。

一般，為使上述方程組有唯一解，在區(qū)間[a，b]的端點a，b（即x1，xn+l）處有以下三種邊界條件：

（a.）給定f（x）在兩端點d，6的二階導數(shù)；

（b.）給定f（x）在兩端點a，b的一階導數(shù)；

（c.）未給任何導數(shù)時，求S'（x）在兩端點的近似值。以X1，x2，x3，X4為基點做一個三次Newton插值多項式Na（X）；再以點xn+1，xn，xn-1，xn-2可得Nb（x），令S'（a）=Na（a），s（b）=Nb（b），剩下的工作同（b.），這種情況也稱作Lagrange三次樣條插值函數(shù)。

1 新邊值條件及理論證明

根據(jù)上述論述知：（5）式是一個含有n+l個未知量m1，m2，…，mn+1，而只有n-l個方程的線性方程組。由矩陣知識易得，該方程組的解包含兩個自由未知量。不同于（a.）（b.）（c.）三種“端點處的邊界條件”，我們給出以下四組“關于中間插值節(jié)點的新邊界條件”：

對比新舊邊界條件易得：（A）（B.）分別是對（a.）（b.）的一般化推廣，（D.）是對（C.）的一般化推廣。對于（A.）（B.）（C.）（D.）中后k，j的關系，有如下12中情形：對（A.）邊界情況有三類：2

為了避免過多贅述，下面僅給出以上12種情況中比較有代表性的情形的理論證明：即（D.）類中3

由（4）式Si（x）求-次積分并化簡可得：

這里，且方程的系數(shù)矩陣是主對角線占優(yōu)的線性方程組，即存在唯一解。則由（II）式可解mj，mj+1…，mk-2，mk-1。

接下來，問題轉(zhuǎn)化成：一方面，已知S"（xj=mj，S"（xj+1）=mj+1求關于未知量m1，m2，…，mj-1的三線型上三角方程組；另一方面，已知s"（Xk-1）=mk-1，s"（Xk）=mk求關于未知量mk+1，mk+2…，mn+1的三線型下三角方程組。

其中（12）（13）式子中Ψi，δi記法同（II）式。

綜上，即得到了所有mi（i=1，2，…，n+l）的值，從而得到唯一確定的樣條函數(shù)。

2 數(shù)值算例

下面給出新邊界條件下樣條函數(shù)的構造算法和計算實例。應用matlab編寫函數(shù)程序，得到樣條函數(shù)對原函數(shù)的擬合對比曲線。

2.1 算法

輸入：進行插值的全部節(jié)點的橫坐標xi（i=1，…，n+l），對應節(jié)點的縱坐標fi（i=l，…，n+1）和邊界條件。

輸出：三次樣條插值函數(shù)si（x）（i=1，…，n）的系數(shù)。

Step l：計算步長hi=Xi+1-xi，其中i=l，…，n；

Step 2：根據(jù)Step l的步長計算差商并形成三對角方程組（或三線型上三角方程組，三線型下三角方程組）；

Step 3：解Step 2中的方程組，找出所有mi（i=1，…，n+l）；

Step 4：將mi（i=l，…，n+l）帶回到樣條函數(shù)Si（x）（i=1，…，n）中；

Step 5：將預測點x帶入f（x）的樣條函數(shù)S（x）得近似值，matlab作圖得插值曲線。

2.2 算例

計算f（x）=X?在以下集合中的點{-5，-4，-3，-2，-1，0，1，2，3，4，5}處的三次樣條插值函數(shù)，并且滿足邊界條件f（-2）=-4，f”（2）=2。

分析 ：同上我們有記法δi=hi+hi+1，Ψi= f[x，+1，xi+2]-f[Xi，xi+1]。題目已知M4=-4，m8=2，步長hi=h=l，首先需要解如下三對角方程組：

利用追趕法懈（14）式的三對角方程組，可得m4，m5，m6，m77的值。

下面分別解m1，m2，m3對應的三線型上三角方程組：

以及解m8，m9，m10，m11對應的三線型下三角方程組：解（15）和（16）式即得剩下的所有mi的值。最后，將Si（x）整理成關于（x-xi）的三次多項式：

由（17）式中（x-Xi）?，（x-Xi）?，（X-xi）的系數(shù)及常數(shù)項關于mi，hi，fi的表達式，可計算出S（x）的系數(shù)矩陣，帶回（17）式即可得插值函數(shù)的表達式。

最后，通過matlab程序得圖1：新邊界條件下樣條函數(shù)與原函數(shù)的擬合對比圖，從圖1可見新邊界條件下的插值虛線與原函數(shù)圖形的走向一致，曲線在坐標系下的位置也相近。

與圖1相對應地，通過matlab軟件得表1：新邊界條件下樣條函數(shù)與原函數(shù)的擬合對比數(shù)據(jù)，由表1中數(shù)據(jù)可知：當步長變小時，相對誤差變小。由此可見新邊界條件下構造樣條函數(shù)的可行性和相對可靠性。

3 結論

本文給出了三次樣條插值函數(shù)更一般化的構造條件，并且由上述新邊界條件的論證過程及數(shù)值實驗可知，任給樣條插值節(jié)點的一階或者二階導數(shù)的值中的兩個，便可以將整個樣條函數(shù)確定。特別地，由（B.）情況可知，可以用“中間節(jié)點的導數(shù)值”取代“邊界節(jié)點的單側導數(shù)值”。進一步看，（c.）情況中選擇Newton插值節(jié)點時可將“端點”換成“中間節(jié)點”，即插值區(qū)間的中部，從而降低了插值誤差，提高了精度。精度的提高，會直接影響算法對生產(chǎn)實際問題處理的有效性。工程中模擬機翼下輪廓曲線等的插值問題可能會要求中間某點的曲率是固定值，或二階導數(shù)為固定值的情況，則本文提出的新邊界條件能有效處理這一情形。此外，邊界條件取中間節(jié)點相對于取端點的情形，對線性方程組的系數(shù)矩陣“三對角矩陣”進行了維數(shù)上的降階和形式上向“三線型上三角”（或“三線型下三角”）的轉(zhuǎn)化。同階的“三對角矩陣”的計算量是O（n?），而“三線型上三角”（或“三線型下三角”）的計算量是O（n），這對大型矩陣計算是有益的。

正如很多實際問題需要從一維推廣到二維、三維情況，我們期望進一步地工作可以在高維樣條插值函數(shù)的構造上展開，將本文提出的新邊界條件構造三次樣條函數(shù)的方法應用于三維情況。綜上，本文對三次樣條函數(shù)的構造提供了新思路，為工業(yè)樣條的計算提供了理論依據(jù)。