# 雅可比法

## 描述

$A\mathbf x = \mathbf b$

$A=\begin{bmatrix} a_{11} & a_{12} & \cdots & a_{1n} \\ a_{21} & a_{22} & \cdots & a_{2n} \\ \vdots & \vdots & \ddots & \vdots \\a_{n1} & a_{n2} & \cdots & a_{nn} \end{bmatrix}, \qquad \mathbf{x} = \begin{bmatrix} x_{1} \\ x_2 \\ \vdots \\ x_n \end{bmatrix} , \qquad \mathbf{b} = \begin{bmatrix} b_{1} \\ b_2 \\ \vdots \\ b_n \end{bmatrix}.$

A 可以分解成对角部分D和剩余部分R:

$A=D+R \qquad \text{其中} \qquad D = \begin{bmatrix} a_{11} & 0 & \cdots & 0 \\ 0 & a_{22} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\0 & 0 & \cdots & a_{nn} \end{bmatrix}, \qquad R = \begin{bmatrix} 0 & a_{12} & \cdots & a_{1n} \\ a_{21} & 0 & \cdots & a_{2n} \\ \vdots & \vdots & \ddots & \vdots \\a_{n1} & a_{n2} & \cdots & 0 \end{bmatrix}$

$D \mathbf{x} = \mathbf{b} - R \mathbf{x}$

$\mathbf{x}^{(k+1)} = D^{-1} (\mathbf{b} - R \mathbf{x}^{(k)}).$

$x^{(k+1)}_i = \frac{1}{a_{ii}} \left(b_i -\sum_{j\ne i}a_{ij}x^{(k)}_j\right),\quad i=1,2,\ldots,n.$

## 算法

while 未收敛
for i := 1 step until n do
$\sigma = 0$
for j := 1 step until n do
if j != i then
$\sigma = \sigma + a_{ij} x_j^{(k-1)}$
end if
end (j-loop)
$x_i^{(k)} = {{\left( {b_i - \sigma } \right)} \over {a_{ii} }}$
end (i-loop)

end (未收敛时继续循环)

## 收敛

$\rho(D^{-1}R) < 1.$[1]

$\left | a_{ii} \right | > \sum_{j \ne i} {\left | a_{ij} \right |}.$

## 举例

$A= \begin{bmatrix} 2 & 1 \\ 5 & 7 \\ \end{bmatrix}, \ b= \begin{bmatrix} 11 \\ 13 \\ \end{bmatrix} , \quad x^{(0)} = \begin{bmatrix} 1 \\ 1 \\ \end{bmatrix} .$

$D^{-1}= \begin{bmatrix} 1/2 & 0 \\ 0 & 1/7 \\ \end{bmatrix}, \ L= \begin{bmatrix} 0 & 0 \\ 5 & 0 \\ \end{bmatrix} , \quad U = \begin{bmatrix} 0 & 1 \\ 0 & 0 \\ \end{bmatrix} .$

$T=-D^{-1}(L+U)$ as

$T= \begin{bmatrix} 1/2 & 0 \\ 0 & 1/7 \\ \end{bmatrix} \left\{ \begin{bmatrix} 0 & 0 \\ -5 & 0 \\ \end{bmatrix} + \begin{bmatrix} 0 & -1 \\ 0 & 0 \\ \end{bmatrix}\right\} = \begin{bmatrix} 0 & -0.5 \\ -0.714 & 0 \\ \end{bmatrix} .$

$C = \begin{bmatrix} 1/2 & 0 \\ 0 & 1/7 \\ \end{bmatrix} \begin{bmatrix} 11 \\ 13 \\ \end{bmatrix} = \begin{bmatrix} 5.5 \\ 1.857 \\ \end{bmatrix} .$

$x^{(1)}= \begin{bmatrix} 0 & -0.5 \\ -0.714 & 0 \\ \end{bmatrix} \begin{bmatrix} 1 \\ 1 \\ \end{bmatrix} + \begin{bmatrix} 5.5 \\ 1.857 \\ \end{bmatrix} = \begin{bmatrix} 5.0 \\ 1.143 \\ \end{bmatrix} .$

$x^{(2)}= \begin{bmatrix} 0 & -0.5 \\ -0.714 & 0 \\ \end{bmatrix} \begin{bmatrix} 5.0 \\ 1.143 \\ \end{bmatrix} + \begin{bmatrix} 5.5 \\ 1.857 \\ \end{bmatrix} = \begin{bmatrix} 4.929 \\ -1.713 \\ \end{bmatrix} .$

$x=\begin{bmatrix} 7.111\\ -3.222 \end{bmatrix} .$

## 一個用Fortran解的例子

subroutine solve(A,b,x,x0,n)
implicit real*8(a-z)
integer::n,imax=200
integer::i,j,k
real*8::tol=1d-15
real*8::A(n,n),b(n),x(n),x0(n),x1(n),x2(n)

write(102,501)
501 format('Jacobi iteration',/)

x1=x0

do k=1,imax
do i=1,n
s=0
do j=1,n
if (j==i) cycle
s=s+A(i,j)*x1(j)
enddo
x2(i)=(b(i)-s)/A(i,i)
enddo
! the following step check if convergence is reached
dx2=0
do i=1,n
dx2=dx2+(x1(i)-x2(i))**2
enddo
dx2=dsqrt(dx2)
if (dx2<tol) exit
x1=x2
write(102,502)k,x1 ! record the iteration process
502 format(1x,i3,<n>(1x,1pd25.15))
enddo
x=x2
end subroutine solve

program  main
implicit real*8(a-z)
integer,parameter::n=3
real*8 ::A(n,n),b(n),x(n),x0(n)

open(unit=101,file='result.txt')
open(unit=102,file='it_result.txt')

x0=(/0d0,0d0,0d0/) ! initial guess
b=(/9d0,7d0,6d0/)
A=reshape((/10,-1,0,-1,10,-4,0,-2,10/),(/3,3/))

call solve(A,b,x,x0,n) ! solve Ax=b

write(101,501)'x(1)','x(2)','x(3)',x
501 format('jacobi result',//,<n>(1x,a25),/,<n>(1x,1pd25.15))

end program main


## 参考文献

1. ^ 1.0 1.1 [1],解线性方程组的迭代法