使用直接方法求解线性系统
解线性系统方法
矩阵分解技术
迭代和最小二乘法解线性系统
迭代方法
雅可比法
给定方程的线性系统:
a 11 x 1 + a 12 x 2 + ⋯ + a 1 n x n = b 1 a 21 x 1 + a 22 x 2 + ⋯ + a 2 n x n = b 2 ⋮ a i 1 x 1 + a i 2 x 2 + ⋯ + a i n x n = b i ⋮ a n 1 x 1 + a n 2 x 2 + ⋯ + a n n x n = b b \boldsymbol{a}{11}\boldsymbol{x}1+\boldsymbol{a}{12}\boldsymbol{x}2+\cdots +\boldsymbol{a}{1\boldsymbol{n}}\boldsymbol{x}{\boldsymbol{n}}=\boldsymbol{b}1\\\boldsymbol{a}{21}\boldsymbol{x}1+\boldsymbol{a}{22}\boldsymbol{x}2+\cdots +\boldsymbol{a}{2\boldsymbol{n}}\boldsymbol{x}{\boldsymbol{n}}=\boldsymbol{b}2\\\vdots \\\boldsymbol{a}{\boldsymbol{i}1}\boldsymbol{x}1+\boldsymbol{a}{\boldsymbol{i}2}\boldsymbol{x}2+\cdots +\boldsymbol{a}{\boldsymbol{in}}\boldsymbol{x}{\boldsymbol{n}}=\boldsymbol{b}{\boldsymbol{i}}\\\vdots \\\boldsymbol{a}{\boldsymbol{n}1}\boldsymbol{x}1+\boldsymbol{a}{\boldsymbol{n}2}\boldsymbol{x}2+\cdots +\boldsymbol{a}{\boldsymbol{nn}}\boldsymbol{x}{\boldsymbol{n}}=\boldsymbol{b}{\boldsymbol{b}} a 11 x 1 + a 12 x 2 + ⋯ + a 1 n x n = b 1 a 21 x 1 + a 22 x 2 + ⋯ + a 2 n x n = b 2 ⋮ a i 1 x 1 + a i 2 x 2 + ⋯ + a in x n = b i ⋮ a n 1 x 1 + a n 2 x 2 + ⋯ + a nn x n = b b
从上面的等式可以得出:
x 1 = 1 a 11 ( b 1 − a 12 x 2 − ⋯ − a n 1 x n ) x 2 = 1 a 22 ( b 2 − a 21 x 1 − a 23 x 3 − ⋯ − a n 1 x n ) ⋮ ⋮ x i = 1 a i i ( b i − a i 1 x 1 − ⋯ − a i i − 1 x i − 1 − a i i + 1 x i + 1 − ⋯ − a i n x n ) ⋮ ⋮ x n = 1 a n n ( b n − a n 1 x 1 − ⋯ − a n n − 1 x n − 1 ) \boldsymbol{x}1=\frac{1}{\boldsymbol{a}{11}}\left( \boldsymbol{b}1-\boldsymbol{a}{12}\boldsymbol{x}2-\cdots -\boldsymbol{a}{\boldsymbol{n}1}\boldsymbol{x}{\boldsymbol{n}} \right) \\\boldsymbol{x}2=\frac{1}{\boldsymbol{a}{22}}\left( \boldsymbol{b}2-\boldsymbol{a}{21}\boldsymbol{x}1-\boldsymbol{a}{23}\boldsymbol{x}3-\cdots -\boldsymbol{a}{\boldsymbol{n}1}\boldsymbol{x}{\boldsymbol{n}} \right) \\\vdots \,\, \vdots \\\boldsymbol{x}{\boldsymbol{i}}=\frac{1}{\boldsymbol{a}{\boldsymbol{ii}}}\left( \boldsymbol{b}{\boldsymbol{i}}-\boldsymbol{a}{\boldsymbol{i}1}\boldsymbol{x}1-\cdots -\boldsymbol{a}{\boldsymbol{ii}-1}\boldsymbol{x}{\boldsymbol{i}-1}-\boldsymbol{a}{\boldsymbol{ii}+1}\boldsymbol{x}{\boldsymbol{i}+1}-\cdots -\boldsymbol{a}{\boldsymbol{in}}\boldsymbol{x}{\boldsymbol{n}} \right) \\\vdots \,\, \vdots \\\boldsymbol{x}{\boldsymbol{n}}=\frac{1}{\boldsymbol{a}{\boldsymbol{nn}}}\left( \boldsymbol{b}{\boldsymbol{n}}-\boldsymbol{a}{\boldsymbol{n}1}\boldsymbol{x}1-\cdots -\boldsymbol{a}{\boldsymbol{nn}-1}\boldsymbol{x}{\boldsymbol{n}-1} \right) x 1 = a 11 1 ( b 1 − a 12 x 2 − ⋯ − a n 1 x n ) x 2 = a 22 1 ( b 2 − a 21 x 1 − a 23 x 3 − ⋯ − a n 1 x n ) ⋮ ⋮ x i = a ii 1 ( b i − a i 1 x 1 − ⋯ − a ii − 1 x i − 1 − a ii + 1 x i + 1 − ⋯ − a in x n ) ⋮ ⋮ x n = a nn 1 ( b n − a n 1 x 1 − ⋯ − a nn − 1 x n − 1 )
雅可比法是一种迭代方法,它从对解[ x 1 ( 0 ) , x 2 ( 0 ) , ⋯ , x n ( 0 ) ] T \left[ \boldsymbol{x}{1}^{\left( 0 \right)},\boldsymbol{x}{2}^{\left( 0 \right)},\cdots ,\boldsymbol{x}_{\boldsymbol{n}}^{\left( 0 \right)} \right] ^{\boldsymbol{T}} [ x 1 ( 0 ) , x 2 ( 0 ) , ⋯ , x n ( 0 ) ] T 的初始猜测开始。 然后,使用迭代k k k 中的解来找到迭代k + 1 k+1 k + 1 中系统解的近似值。 完成如下:
x 1 ( k + 1 ) = 1 a 11 ( b 1 − a 12 x 2 ( k ) − ⋯ − a n 1 x n ( k ) ) x 2 ( k + 1 ) = 1 a 22 ( b 2 − a 21 x 1 ( k ) − a 23 x 3 ( k ) − ⋯ − a n 1 x n ( k ) ) ⋮ ⋮ x i ( k + 1 ) = 1 a i i ( b i − a i 1 x 1 ( k ) − ⋯ − a i i − 1 x i − 1 ( k ) − a i i + 1 x i + 1 ( k ) − ⋯ − a i n x n ( k ) ) ⋮ ⋮ x n ( k + 1 ) = 1 a n n ( b n − a n 1 x 1 ( k ) − ⋯ − a n n − 1 x n − 1 ( k ) ) \boldsymbol{x}{1}^{\left( \boldsymbol{k}+1 \right)}=\frac{1}{\boldsymbol{a}{11}}\left( \boldsymbol{b}1-\boldsymbol{a}{12}\boldsymbol{x}{2}^{\left( \boldsymbol{k} \right)}-\cdots -\boldsymbol{a}{\boldsymbol{n}1}\boldsymbol{x}{\boldsymbol{n}}^{\left( \boldsymbol{k} \right)} \right) \\\boldsymbol{x}{2}^{\left( \boldsymbol{k}+1 \right)}=\frac{1}{\boldsymbol{a}{22}}\left( \boldsymbol{b}2-\boldsymbol{a}{21}\boldsymbol{x}{1}^{\left( \boldsymbol{k} \right)}-\boldsymbol{a}{23}\boldsymbol{x}{3}^{\left( \boldsymbol{k} \right)}-\cdots -\boldsymbol{a}{\boldsymbol{n}1}\boldsymbol{x}{\boldsymbol{n}}^{\left( \boldsymbol{k} \right)} \right) \\\vdots \,\, \vdots \\\boldsymbol{x}{\boldsymbol{i}}^{\left( \boldsymbol{k}+1 \right)}=\frac{1}{\boldsymbol{a}{\boldsymbol{ii}}}\left( \boldsymbol{b}{\boldsymbol{i}}-\boldsymbol{a}{\boldsymbol{i}1}\boldsymbol{x}{1}^{\left( \boldsymbol{k} \right)}-\cdots -\boldsymbol{a}{\boldsymbol{ii}-1}\boldsymbol{x}{\boldsymbol{i}-1}^{\left( \boldsymbol{k} \right)}-\boldsymbol{a}{\boldsymbol{ii}+1}\boldsymbol{x}{\boldsymbol{i}+1}^{\left( \boldsymbol{k} \right)}-\cdots -\boldsymbol{a}{\boldsymbol{in}}\boldsymbol{x}{\boldsymbol{n}}^{\left( \boldsymbol{k} \right)} \right) \\\vdots \,\, \vdots \\\boldsymbol{x}{\boldsymbol{n}}^{\left( \boldsymbol{k}+1 \right)}=\frac{1}{\boldsymbol{a}{\boldsymbol{nn}}}\left( \boldsymbol{b}{\boldsymbol{n}}-\boldsymbol{a}{\boldsymbol{n}1}\boldsymbol{x}{1}^{\left( \boldsymbol{k} \right)}-\cdots -\boldsymbol{a}{\boldsymbol{nn}-1}\boldsymbol{x}{\boldsymbol{n}-1}^{\left( \boldsymbol{k} \right)} \right) x 1 ( k + 1 ) = a 11 1 ( b 1 − a 12 x 2 ( k ) − ⋯ − a n 1 x n ( k ) ) x 2 ( k + 1 ) = a 22 1 ( b 2 − a 21 x 1 ( k ) − a 23 x 3 ( k ) − ⋯ − a n 1 x n ( k ) ) ⋮ ⋮ x i ( k + 1 ) = a ii 1 ( b i − a i 1 x 1 ( k ) − ⋯ − a ii − 1 x i − 1 ( k ) − a ii + 1 x i + 1 ( k ) − ⋯ − a in x n ( k ) ) ⋮ ⋮ x n ( k + 1 ) = a nn 1 ( b n − a n 1 x 1 ( k ) − ⋯ − a nn − 1 x n − 1 ( k ) )
通常,迭代x i ( k + 1 ) \boldsymbol{x}_{\boldsymbol{i}}^{\left( \boldsymbol{k}+1 \right)} x i ( k + 1 ) 中的解可以表示为:
x i ( k + 1 ) = 1 a i i ( b i − ∑ j = 1 , j ≠ i n a i j x j ( k ) ) , i = 1 , ⋯ , n \boldsymbol{x}{\boldsymbol{i}}^{\left( \boldsymbol{k}+1 \right)}=\frac{1}{\boldsymbol{a}{\boldsymbol{ii}}}\left( \boldsymbol{b}{\boldsymbol{i}}-\sum{\boldsymbol{j}=1,\boldsymbol{j}\ne \boldsymbol{i}}^{\boldsymbol{n}}{\boldsymbol{a}{\boldsymbol{ij}}\boldsymbol{x}{\boldsymbol{j}}^{\left( \boldsymbol{k} \right)}} \right) ,\boldsymbol{i}=1,\cdots ,\boldsymbol{n} x i ( k + 1 ) = a ii 1 ( b i − ∑ j = 1 , j = i n a ij x j ( k ) ) , i = 1 , ⋯ , n
对于任意ε > 0 \boldsymbol{\varepsilon }>0 ε > 0 ,当∥ x i ( k + 1 ) − x i ( k ) ∥ < ε \left\| \boldsymbol{x}{\boldsymbol{i}}^{\left( \boldsymbol{k}+1 \right)}-\boldsymbol{x}{\boldsymbol{i}}^{\left( \boldsymbol{k} \right)} \right\| <\boldsymbol{\varepsilon } x i ( k + 1 ) − x i ( k ) < ε 时,雅可比迭代停止。
**示例1,**为线性系统编写雅可比方法的前三个迭代:
( 2 − 1 1 − 2 5 − 1 1 − 2 4 ) ( x 1 x 2 x 3 ) = ( − 1 1 3 ) \left( \begin{matrix} 2& -1& 1\\ -2& 5& -1\\ 1& -2& 4\\\end{matrix} \right) \left( \begin{array}{c} \boldsymbol{x}_1\\ \boldsymbol{x}_2\\ \boldsymbol{x}_3\\\end{array} \right) =\left( \begin{array}{c} -1\\ 1\\ 3\\\end{array} \right) 2 − 2 1 − 1 5 − 2 1 − 1 4 x 1 x 2 x 3 = − 1 1 3
从零向量开始
x ( 0 ) = ( 0 0 0 ) \boldsymbol{x}^{\left( 0 \right)}=\left( \begin{array}{c} 0\\ 0\\ 0\\\end{array} \right) x ( 0 ) = 0 0 0
解决
我们写成:
x 1 ( k + 1 ) = 1 2 ( − 1 + x 2 ( k ) − x 3 ( k ) ) x 2 ( k + 1 ) = 1 5 ( 1 + 2 x 1 ( k ) + x 3 ( k ) ) x 3 ( k + 1 ) = 1 4 ( 3 − x 1 ( k ) + 2 x 2 ( k ) ) \boldsymbol{x}{1}^{\left( \boldsymbol{k}+1 \right)}=\frac{1}{2}\left( -1+\boldsymbol{x}{2}^{\left( \boldsymbol{k} \right)}-\boldsymbol{x}{3}^{\left( \boldsymbol{k} \right)} \right) \\\boldsymbol{x}{2}^{\left( \boldsymbol{k}+1 \right)}=\frac{1}{5}\left( 1+2\boldsymbol{x}{1}^{\left( \boldsymbol{k} \right)}+\boldsymbol{x}{3}^{\left( \boldsymbol{k} \right)} \right) \\\boldsymbol{x}{3}^{\left( \boldsymbol{k}+1 \right)}=\frac{1}{4}\left( 3-\boldsymbol{x}{1}^{\left( \boldsymbol{k} \right)}+2\boldsymbol{x}_{2}^{\left( \boldsymbol{k} \right)} \right) x 1 ( k + 1 ) = 2 1 ( − 1 + x 2 ( k ) − x 3 ( k ) ) x 2 ( k + 1 ) = 5 1 ( 1 + 2 x 1 ( k ) + x 3 ( k ) ) x 3 ( k + 1 ) = 4 1 ( 3 − x 1 ( k ) + 2 x 2 ( k ) )
第一次迭代k = 0 k=0 k = 0
x 1 ( 1 ) = 1 2 ( − 1 + x 2 ( 0 ) − x 3 ( 0 ) ) = 1 2 ( − 1 + 0 − 0 ) = − 1 2 x 2 ( 1 ) = 1 5 ( 1 + 2 x 1 ( 0 ) + x 3 ( 0 ) ) = 1 5 ( 1 + 2 ( 0 ) + 0 ) = 1 5 x 3 ( 1 ) = 1 4 ( 3 − x 1 ( 0 ) + 2 x 2 ( 0 ) ) = 1 4 ( 3 − 0 + 2 ( 0 ) ) = 3 4 \boldsymbol{x}{1}^{\left( 1 \right)}=\frac{1}{2}\left( -1+\boldsymbol{x}{2}^{\left( 0 \right)}-\boldsymbol{x}{3}^{\left( 0 \right)} \right) =\frac{1}{2}\left( -1+0-0 \right) =-\frac{1}{2}\\\boldsymbol{x}{2}^{\left( 1 \right)}=\frac{1}{5}\left( 1+2\boldsymbol{x}{1}^{\left( 0 \right)}+\boldsymbol{x}{3}^{\left( 0 \right)} \right) =\frac{1}{5}\left( 1+2\left( 0 \right) +0 \right) =\frac{1}{5}\\\boldsymbol{x}{3}^{\left( 1 \right)}=\frac{1}{4}\left( 3-\boldsymbol{x}{1}^{\left( 0 \right)}+2\boldsymbol{x}_{2}^{\left( 0 \right)} \right) =\frac{1}{4}\left( 3-0+2\left( 0 \right) \right) =\frac{3}{4} x 1 ( 1 ) = 2 1 ( − 1 + x 2 ( 0 ) − x 3 ( 0 ) ) = 2 1 ( − 1 + 0 − 0 ) = − 2 1 x 2 ( 1 ) = 5 1 ( 1 + 2 x 1 ( 0 ) + x 3 ( 0 ) ) = 5 1 ( 1 + 2 ( 0 ) + 0 ) = 5 1 x 3 ( 1 ) = 4 1 ( 3 − x 1 ( 0 ) + 2 x 2 ( 0 ) ) = 4 1 ( 3 − 0 + 2 ( 0 ) ) = 4 3
第二次迭代k = 1 k=1 k = 1
x 1 ( 2 ) = 1 2 ( − 1 + x 2 ( 1 ) − x 3 ( 1 ) ) = 1 2 ( − 1 + 1 5 − 3 4 ) = − 31 40 x 2 ( 2 ) = 1 5 ( 1 + 2 x 1 ( 1 ) + x 3 ( 1 ) ) = 1 5 ( 1 + 2 ⋅ − 1 2 + 3 4 ) = 3 20 x 3 ( 2 ) = 1 4 ( 3 − x 1 ( 1 ) + 2 x 2 ( 1 ) ) = 1 4 ( 3 − − 1 2 + 2 1 5 ) = 39 40 \boldsymbol{x}{1}^{\left( 2 \right)}=\frac{1}{2}\left( -1+\boldsymbol{x}{2}^{\left( 1 \right)}-\boldsymbol{x}{3}^{\left( 1 \right)} \right) =\frac{1}{2}\left( -1+\frac{1}{5}-\frac{3}{4} \right) =-\frac{31}{40}\\\boldsymbol{x}{2}^{\left( 2 \right)}=\frac{1}{5}\left( 1+2\boldsymbol{x}{1}^{\left( 1 \right)}+\boldsymbol{x}{3}^{\left( 1 \right)} \right) =\frac{1}{5}\left( 1+2\cdot \frac{-1}{2}+\frac{3}{4} \right) =\frac{3}{20}\\\boldsymbol{x}{3}^{\left( 2 \right)}=\frac{1}{4}\left( 3-\boldsymbol{x}{1}^{\left( 1 \right)}+2\boldsymbol{x}_{2}^{\left( 1 \right)} \right) =\frac{1}{4}\left( 3-\frac{-1}{2}+2\frac{1}{5} \right) =\frac{39}{40} x 1 ( 2 ) = 2 1 ( − 1 + x 2 ( 1 ) − x 3 ( 1 ) ) = 2 1 ( − 1 + 5 1 − 4 3 ) = − 40 31 x 2 ( 2 ) = 5 1 ( 1 + 2 x 1 ( 1 ) + x 3 ( 1 ) ) = 5 1 ( 1 + 2 ⋅ 2 − 1 + 4 3 ) = 20 3 x 3 ( 2 ) = 4 1 ( 3 − x 1 ( 1 ) + 2 x 2 ( 1 ) ) = 4 1 ( 3 − 2 − 1 + 2 5 1 ) = 40 39
第三次迭代$k=2$
x 1 ( 3 ) = 1 2 ( − 1 + x 2 ( 2 ) − x 3 ( 2 ) ) = 1 2 ( − 1 + 3 20 − 39 40 ) = − 73 80 x 2 ( 3 ) = 1 5 ( 1 + 2 x 1 ( 2 ) + x 3 ( 2 ) ) = 1 5 ( 1 + 2 ⋅ − 31 40 − − 73 80 ) = 17 200 x 3 ( 3 ) = 1 4 ( 3 − x 1 ( 2 ) + 2 x 2 ( 2 ) ) = 1 4 ( 3 − 31 40 + 2 3 20 ) = 163 160 \boldsymbol{x}{1}^{\left( 3 \right)}=\frac{1}{2}\left( -1+\boldsymbol{x}{2}^{\left( 2 \right)}-\boldsymbol{x}{3}^{\left( 2 \right)} \right) =\frac{1}{2}\left( -1+\frac{3}{20}-\frac{39}{40} \right) =-\frac{73}{80}\\\boldsymbol{x}{2}^{\left( 3 \right)}=\frac{1}{5}\left( 1+2\boldsymbol{x}{1}^{\left( 2 \right)}+\boldsymbol{x}{3}^{\left( 2 \right)} \right) =\frac{1}{5}\left( 1+2\cdot \frac{-31}{40}-\frac{-73}{80} \right) =\frac{17}{200}\\\boldsymbol{x}{3}^{\left( 3 \right)}=\frac{1}{4}\left( 3-\boldsymbol{x}{1}^{\left( 2 \right)}+2\boldsymbol{x}_{2}^{\left( 2 \right)} \right) =\frac{1}{4}\left( 3-\frac{31}{40}+2\frac{3}{20} \right) =\frac{163}{160} x 1 ( 3 ) = 2 1 ( − 1 + x 2 ( 2 ) − x 3 ( 2 ) ) = 2 1 ( − 1 + 20 3 − 40 39 ) = − 80 73 x 2 ( 3 ) = 5 1 ( 1 + 2 x 1 ( 2 ) + x 3 ( 2 ) ) = 5 1 ( 1 + 2 ⋅ 40 − 31 − 80 − 73 ) = 200 17 x 3 ( 3 ) = 4 1 ( 3 − x 1 ( 2 ) + 2 x 2 ( 2 ) ) = 4 1 ( 3 − 40 31 + 2 20 3 ) = 160 163
**示例2,**雅可比方法将用于求解线性系统:
( − 5 1 − 2 1 6 3 2 − 1 − 4 ) ( x 1 x 2 x 3 ) = ( 13 1 − 1 ) \left( \begin{matrix} -5& 1& -2\\ 1& 6& 3\\ 2& -1& -4\\\end{matrix} \right) \left( \begin{array}{c} \boldsymbol{x}_1\\ \boldsymbol{x}_2\\ \boldsymbol{x}_3\\\end{array} \right) =\left( \begin{array}{c} 13\\ 1\\ -1\\\end{array} \right) − 5 1 2 1 6 − 1 − 2 3 − 4 x 1 x 2 x 3 = 13 1 − 1
MATLAB代码
Copy function x = JacobiSolve(A, b, Eps)
n = length(b) ;
x0 = zeros(3, 1) ;
x = ones(size(x0)) ;
while norm(x-x0, inf) ≥ Eps
x0 = x ;
for i = 1 : n
x(i) = b(i) ;
for j = 1 : n
if j , i
x(i) = x(i) - A(i, j)*x0(j) ;
end
end
x(i) = x(i) / A(i, i) ;
end
Copy >> A = [-5 1 -2; 1 6 3; 2 -1 -4] ;
>> b = [13; 1; -1] ;
>> JacobiSolveLinSystem(A, b)
-2.0000
1.0000
-1.0000
使用Python,功能JacobiSolve具有以下代码:
Copy def JacobiSolve(A, b, Eps):
import numpy as np
n = len(b)
x0, x = np.zeros((n, 1), 'float'), np.ones((n, 1), 'float')
while np.linalg.norm(x-x0, np.inf) ≥ Eps:
x0 = x.copy()
for i in range(n):
x[i] = b[i]
for j in range(n):
if j != i:
x[i] -= A[i][j]*x0[j]
x[i] /= A[i][i]
return x
通过调用函数JacobiSolve求解给定的线性系统,我们获得:
Copy In [3]: A = np.array([[-5, 1, -2], [1, 6, 3], [2, -1, -4]])
In [4]: b = np.array([[13], [1], [-1]])
In [5]: x = JacobiSolve(A, b, Eps)
In [6]: print(’x = \\n’, x)
Out[6]:
x =
[[-2. ],
[ 1.00000002],
[-1.00000002]]
最小二乘解
线性系统解中的病态和正则化技术
解非线性方程组
数据插补
微分和积分
解非线性常微分方程
常微分方程的非标准有限差分法
线性规划和二次规划
非线性规划
最优控制问题