I am looking at the Jacobi Method which is the simplest (and slowest) iterative method for solving a large system of linear equations. Such systems are too large for traditional methods to work. I am solving:

A x = b

The Jacobi method uses the following matrix form for its iterations:

D x = b – (L + U) x

Where D, L and U are the diagonal, strictly left lower triangular, and strictly right upper triangular matrices, respectively. The algorithm is guaranteed to work, with any initial guess for x, if matrix A has diagonal dominant elements.
In case matrix A is not diagonal dominant, I am trying to find a trick to add a preconditioning matrix P to make A diagonal dominant. P is a diagonal matrix with zeros in non-diagonal elements. Each diagonal value of matrix P is equal to the sum of absolute values in the corresponding row of matrix A. Thus I have:

(A + P) x = b + P x

The iteration formula becomes:

(D + P) x = b + P x – (L + U) x

If I use a guess for x on the right hand side, the iterations do not converge. However, if I solve for x using a traditional method to get a good approximation for x and use that in evaluating the term P x, then the above equation converges to a more accurate solution. In fact, using this approach the expression b + P x can be regarded as equal to b' which is a new constants vector. So basically I am changing the original equation to a more favorable form by updating matrix A and vector b.

The question to is this. Am I cheating in the above method that I described? Can you offer any better solutions?

Namir

PS: whatever approach works with the Jacobi method can also be applied to the Gauss-Seidel method which is a better method. There are algorithms (variants of the conjugate gradient method) that solve large linear systems of equations where matrix A is not diagonal dominant or symmetrical. I am trying to go back to basics and see if the Jacobi and Gauss-Seidel methods can get a face lift.

Edited: 6 Nov 2012, 10:41 a.m.