CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Wiki > Jacobi method

Jacobi method

From CFD-Wiki

Revision as of 07:54, 8 April 2006 by Jasond (Talk | contribs)
Jump to: navigation, search

Contents

Introduction

We seek the solution to set of linear equations:

A \phi = b

In matrix terms, the definition of the Jacobi method can be expressed as :

\phi^{(k+1)} = D^{ - 1} \left[\left( {L + U} \right)\phi^{(k)} + b\right]

where D, L, and U represent the diagonal, lower triangular, and upper triangular parts of the coefficient matrix A and k is the iteration count. This matrix expression is mainly of academic interest, and is not used to program the method. Rather, an element-based approach is used:

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

Note that the computation of \phi^{(k+1)}_i requires each element in \phi^{(k)} except itself. Then, unlike in the Gauss-Seidel method, we can't overwrite \phi^{(k)}_i with \phi^{(k+1)}_i, as that value will be needed by the rest of the computation. This is the most meaningful difference between the Jacobi and Gauss-Seidel methods. The minimum ammount of storage is two vectors of size n, and explicit copying will need to take place.

Algorithm

Chose an initial guess \phi^{0} to the solution

for k := 1 step 1 untill convergence do
for i := 1 step until n do
\sigma = 0
for j := 1 step until n do
if j != i then
\sigma = \sigma + a_{ij} \phi_j^{(k-1)}
end if
end (j-loop)
\phi_i^{(k)} = {{\left( {b_i - \sigma } \right)} \over {a_{ii} }}
end (i-loop)
check if convergence is reached
end (k-loop)

Convergence

It is proven that if the absolute value of the diagonal term is always greater than the sum of the absolute values of other term in the row:

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

then the method always converge.

Usually, but not always, the method converges even if this condition is not satisfied, but the diagonal terms in the matrix are greater by the absolute values than the other terms.

Example Calculation

As with Gauss-Seidel, Jacobi iteration lends itself to situations in which we need not explicitly represent the matrix. Consider the simple heat equation problem

\nabla^2 T(x) = 0,\ x\in [0,1]

subject to the boundary conditions T(0)=0 and T(1)=1. The exact solution to this problem is T(x)=x. The standard second-order finite difference discretization is

T_{i-1}-2T_i+T_{i+1} = 0,

where T_i is the (discrete) solution available at uniformly spaced nodes (see the Gauss-Seidel example for the matrix form). For any given T_i for 1 < i < n, we can write

T_i = \frac{1}{2}(T_{i-1}+T_{i+1}).

Then, stepping through the solution vector from i=2 to i=n-1, we can compute the next iterate from the two surrounding values. For a proper Jacobi iteration, we'll need to use values from the previous iteration on the right-hand side:

T_i^{k+1} = \frac{1}{2}(T_{i-1}^{k}+T_{i+1}^k).

The following table gives the results of 10 iterations with 5 nodes (3 interior and 2 boundary) as well as L_2 norm error.

Jacobi Solution
Iteration T_1 T_2 T_3 T_4 T_5 L_2 error
0 0.0000E+00 0.0000E+00 0.0000E+00 0.0000E+00 1.0000E+00 1.0000E+00
1 0.0000E+00 0.0000E+00 0.0000E+00 5.0000E-01 1.0000E+00 6.1237E-01
2 0.0000E+00 0.0000E+00 2.5000E-01 5.0000E-01 1.0000E+00 4.3301E-01
3 0.0000E+00 1.2500E-01 2.5000E-01 6.2500E-01 1.0000E+00 3.0619E-01
4 0.0000E+00 1.2500E-01 3.7500E-01 6.2500E-01 1.0000E+00 2.1651E-01
5 0.0000E+00 1.8750E-01 3.7500E-01 6.8750E-01 1.0000E+00 1.5309E-01
6 0.0000E+00 1.8750E-01 4.3750E-01 6.8750E-01 1.0000E+00 1.0825E-01
7 0.0000E+00 2.1875E-01 4.3750E-01 7.1875E-01 1.0000E+00 7.6547E-02
8 0.0000E+00 2.1875E-01 4.6875E-01 7.1875E-01 1.0000E+00 5.4127E-02
9 0.0000E+00 2.3438E-01 4.6875E-01 7.3438E-01 1.0000E+00 3.8273E-02
10 0.0000E+00 2.3438E-01 4.8438E-01 7.3438E-01 1.0000E+00 2.7063E-02


External link

Retrieved from "https://openfoam.cfd-online.com/Wiki/Jacobi_method"
My wiki