線性代數在數值分析的應用 (一)：兩點邊值問題

Posted on 02/24/2010 by ccjou

本文的閱讀等級：初級

數值分析 (numerical analysis) 是一門研究連續數學問題演算法的科目，內容包羅萬象，如解非線性方程、解線性或非線性聯立方程組、插值 (interpolation) 與曲線配適 (curve fitting)、數值微分與積分、常微分方程與偏微分方程的數值解、邊值 (boundary-value) 與特徵值 (characteristic-value) 問題等。線性代數和數值分析有極密切的關係，一方面數值分析涵蓋一些數值線性代數問題，另一方面數值分析的部分演算法建立於線性代數之上。本文介紹線性代數在數值分析中的一個簡單應用──將兩點邊值 (two-point boundary-value) 問題轉換為線性聯立方程組的求解問題。

 
考慮下面的二階微分方程：

{x}''-x=t

兩邊值為 x(0)=2，x(1)=-1。我們的目標是計算 x 於區間 [0,1] 中 n 個等距間隔點 t_i 的近似值，t_i=i/(n+1)，i=1,2,\ldots,n。

 
求兩點邊值問題數值解的第一步驟是將連續函數的微分以差分取代。設函數所定義的區間為 [a,b]，將此區間分割為 n+1 個等份，令 h=\frac{b-a}{n+1}，切割的等間隔點為 t_i=a+ih，i=0,1,\ldots,n+1。考慮在 t_i+h 和 t_i-h 的泰勒展開式 (Taylor series expansion)：

\begin{aligned}\displaystyle x(t_i+h)&=x(t_i)+x'(t_i)h+ \frac{h^2}{2!}x'' (t_i)+\frac{h^3}{3!}x''' (t_i)+\cdots\\ x(t_i-h)&=x(t_i)-x'(t_i)h+ \frac{h^2}{2!}x'' (t_i)-\frac{h^3}{3!}x''' (t_i)+\cdots\end{aligned}

將上面二式相減與相加，分別得到

\begin{aligned}\displaystyle x'(t_i)&=\frac{x(t_i+h)-x(t_i-h)}{2h}+O(h^2)\\ x''(t_i)&=\frac{x(t_i+h)-2x(t_i)+x(t_i-h)}{h^2}+O(h^2)\end{aligned}

其中符號 O(h^k) 代表包含 k 次冪和更高次冪的項。為簡化書寫，令 x_{i}=x(t_i)。如果捨去 O(h^2)，結果稱為中央差分近似 (central difference approximation)，如下：

\begin{aligned}\displaystyle x'(t_i)&\approx\frac{x_{i+1}-x_{i-1}}{2h}\\ x''(t_i)&\approx\frac{x_{i+1}-2x_i+x_{i-1}}{h^2}\end{aligned}

已知 a=0，b=1，故 t_i=ih，h=1/(n+1)。當 t=t_i 時，將 x''(t_i) 的近似式代入原方程式：

\displaystyle \frac{x_{i+1}-2x_i+x_{i-1}}{h^2}-x_i=t_i

整理為

x_{i-1}-(2+h^2)x_i+x_{i+1}=h^2t_i=ih^3,~~i=1,2,\ldots,n

另外還需考慮已知邊值 x_0=x(t_0)=x(0)=2 和 x_{n+1}=x(t_{n+1})=x(1)=-1。

 
設未知向量 \mathbf{x}=(x_1,x_2,\ldots,x_n)^T，上述 n 個聯立方程式可表示為矩陣形式 A\mathbf{x}=\mathbf{b}，A 為 n\times n 階，如下：

\begin{aligned} &\begin{bmatrix} -2-h^2&1&0&\cdots&0&0&0\\ 1&-2-h^2&1&\cdots&0&0&0\\ 0&1&-2-h^2&\cdots&0&0&0\\ \vdots&\vdots&\vdots&\ddots&\vdots&\vdots\\ 0&0&0&\cdots&-2-h^2&1&0\\ 0&0&0&\cdots&1&-2-h^2&1\\ 0&0&0&\cdots&0&1&-2-h^2 \end{bmatrix}\begin{bmatrix} x_1\\ x_2\\ x_3\\ \vdots\\ x_{n-2}\\ x_{n-1}\\ x_n \end{bmatrix}\\ &=\begin{bmatrix} h^3-2\\ 2h^3\\ 3h^3\\ \vdots\\ (n-2)h^3\\ (n-1)h^3\\ nh^3+1 \end{bmatrix}\end{aligned}

解出此線性方程即可得到 x(t) 的近似值，如果增加區間切割數 n+1，h 值將隨之減小，可得到更緊密且精確的結果。當然，代價是必須解出一個尺寸較大的 n 階線性方程。觀察上面的係數矩陣形式發現它是一個三對角矩陣，以高斯消去法化簡三對角矩陣僅需消耗 O(n) 計算量 (見“特殊矩陣 (11)：三對角矩陣”)。

Advertisement
This entry was posted in 線性代數專欄, 應用之道 and tagged , , , , . Bookmark the permalink.

Leave a Reply

Fill in your details below or click an icon to log in:

Gravatar
WordPress.com Logo

You are commenting using your WordPress.com account. ( Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. ( Log Out /  Change )

Cancel

Connecting to %s