Implementing Gauss Seidel Method in Java

Gauss-Seidel (also known as successive displacement method) is a mathematical computational method mainly used to find the solution of a System of Linear Algebra. Ever heard of the Jacobi Method? Well, the Gauss-Seidel is nothing but a better version of the Jacobi.  

Example :

Input:
Enter the number of variables in the equation:
2
Enter the augmented matrix:
1 2 3
6 5 4
 
Output: 
6.0 5.0 4.0 
1.0 2.0 3.0 
X0 = {0.6666666666666666 1.1666666666666667 }
X1 = {-0.30555555555555564 1.652777777777778 }
X2 = {-0.7106481481481481 1.855324074074074 }
X3 = {-0.8794367283950617 1.9397183641975309 }
X4 = {-0.9497653034979425 1.9748826517489713 }
X5 = {-0.9790688764574759 1.9895344382287379 }
X6 = {-0.9912786985239483 1.9956393492619742 }
X7 = {-0.9963661243849785 1.9981830621924892 }
X8 = {-0.9984858851604077 1.9992429425802039 }
X9 = {-0.9993691188168363 1.999684559408418 }
X10 = {-0.9997371328403484 1.999868566420174 }
X11 = {-0.9998904720168117 1.9999452360084058 }
X12 = {-0.999954363340338 1.999977181670169 }
X13 = {-0.9999809847251406 1.9999904923625702 }
X14 = {-0.9999920769688085 1.9999960384844042 }
X15 = {-0.9999966987370034 1.9999983493685016 }
X16 = {-0.9999986244737512 1.9999993122368755 }
X17 = {-0.9999994268640631 1.9999997134320315 }
X18 = {-0.9999997611933598 1.9999998805966799 }
X19 = {-0.9999999004972331 1.9999999502486165 }
X20 = {-0.9999999585405137 1.9999999792702567 }
X21 = {-0.999999982725214 1.999999991362607 }
X22 = {-0.9999999928021724 1.9999999964010862 }
X23 = {-0.999999997000905 1.9999999985004524 }
X24 = {-0.999999998750377 1.9999999993751885 }
X25 = {-0.9999999994793237 1.9999999997396618 }
X26 = {-0.9999999997830514 1.9999999998915257 }
X27 = {-0.9999999999096048 1.9999999999548024 }
X28 = {-0.9999999999623352 1.9999999999811675 }
X29 = {-0.9999999999843061 1.999999999992153 }
X30 = {-0.9999999999934606 1.9999999999967302 }
X31 = {-0.9999999999972751 1.9999999999986375 }
X32 = {-0.9999999999988646 1.9999999999994322 }
X33 = {-0.9999999999995268 1.9999999999997633 }
X34 = {-0.9999999999998028 1.9999999999999014 }
X35 = {-0.9999999999999176 1.9999999999999587 }
X36 = {-0.9999999999999656 1.9999999999999827 }
X37 = {-0.9999999999999855 1.9999999999999927 }
X38 = {-0.9999999999999938 1.999999999999997 }
X39 = {-0.9999999999999973 1.9999999999999987 }
X40 = {-0.9999999999999988 1.9999999999999993 }
X41 = {-0.9999999999999993 1.9999999999999996 }

A basic mathematical algorithm for the Gauss Seidel:

Given a set of n equations and n unknowns:

a11x1+a12x2+a13x3+......+a1nxn=c1

a21x1+a22x2+a23x3+......+a2nxn=c2
.                        .
.                        .
.                        .
a1nx1+a2nx2+a3nx3+......+annxn=cn

If the diagonal elements are non-zero, each equation is rewritten for the corresponding unknown, that is, the first equation is rewritten with x₁ on the left hand side, the second equation is rewritten with x₂ on the left hand side and so on as follows:

x1=(c1-a12x2-a13x3-....-a1nxn)/a11
x2=(c2-a21x1-a23x3-....-a2nxn)/a22
.
.
.
xn=(cn-an1x1-an2x2-....-an.n-1xn-1)/ann

.These equations can be re-written in the form as:

x1=(c1-nj=1,j!=1∑a1jxj)/a11
x2=(c2-nj=1,j!=2∑a2jxj)/a22
.
.
.
x2=(cn-nj=1,j!=n∑anjxj)/ann

Hence, for any row i ;

xi=(ci- nj=i,j!=i∑aijxj)/aii , i=1,2,..,n

Note: Gauss-Seidel is applicable to strictly diagonally dominant or symmetric positive definite.

In short, for a linear system of equations, say: 

Ax=b

Deduced formula:

xi(k)=(bi-∑j<i a i j x j(k)-∑j>i a i j x j(k-1))/aii

In terms of matrices, the definition of Gauss-Seidel:

x(k)=(D-L)-1(Ux(k-1)+b)

Here, D represents the diagonal part of the matrix A, -L strictly lower triangle part of A, U strictly upper triangle part of A.

Programming algorithm:

1. Take the number of variables in the equation and the values of the augmented matrix 
    (formed by appending the columns of the two given matrices) as the input.
2. Checking if the given augmented matrix is diagonally dominant.
3. Attempting to transform the given augmented matrix into diagonally dominant. If proved otherwise, then an appropriate message returned to the user.
4. Computing and printing the final result as well as the iterations.
   

Implementation in Java: 

Java

// Implementing Gauss Seidel Method in Java
import java.io.BufferedReader;
import java.io.IOException;
import java.io.InputStreamReader;
import java.io.PrintWriter;
import java.util.Arrays;
import java.util.StringTokenizer;

class GFG {
    // we set a max number of iterations to
    // prevent an infinite loop
    public static final int MAX_ITERATIONS = 100;
    private double[][] M;
    public GFG(double[][] matrix) { M = matrix; }

    public void print() // printing
    {
        int n = M.length;
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n + 1; j++)
                System.out.print(M[i][j] + " ");
            System.out.println();
        }
    }
    // attempting to change a matrix to dominant
    // if proved that it is not
    public boolean transformToDominant(int r, boolean[] V,
                                       int[] R)
    {
        int n = M.length;
        if (r == M.length) {
            double[][] T = new double[n][n + 1];
            for (int i = 0; i < R.length; i++) {
                for (int j = 0; j < n + 1; j++)
                    T[i][j] = M[R[i]][j];
            }
            M = T;
            return true;
        }
        for (int i = 0; i < n; i++) {
            if (V[i])
                continue;
            double sum = 0;
            for (int j = 0; j < n; j++)
                sum += Math.abs(M[i][j]);
            if (2 * Math.abs(M[i][r]) > sum) {
                // diagonally dominant?
                V[i] = true;
                R[r] = i;
                if (transformToDominant(r + 1, V, R))
                    return true;
                V[i] = false;
            }
        }
        return false;
    }
    // method to check whether matrix is
    // diagonally dominant or not
    public boolean makeDominant()
    {
        boolean[] visited = new boolean[M.length];
        int[] rows = new int[M.length];
        Arrays.fill(visited, false);
        return transformToDominant(0, visited, rows);
    }

    // method to find the solution of the matrix
    // after all conditions are satisfied
    public void solve()
    {
        int iterations = 0;
        int n = M.length;
        double epsilon = 1e-15;
        double[] X = new double[n]; // Approximations
        double[] P = new double[n]; // Prev
        Arrays.fill(X, 0);
        while (true) {
            for (int i = 0; i < n; i++) {
                double sum = M[i][n]; // b_n
                for (int j = 0; j < n; j++)
                    if (j != i)
                        sum -= M[i][j] * X[j];
                // Update xi to use in the next
                // row calculation
                X[i] = 1 / M[i][i] * sum;
            }
            System.out.print("X" + iterations + " = {");
            for (int i = 0; i < n; i++)
                System.out.print(X[i] + " ");
            System.out.println("}");
            iterations++;
            if (iterations == 1)
                continue;
            boolean stop = true;
            for (int i = 0; i < n && stop; i++)
                if (Math.abs(X[i] - P[i]) > epsilon)
                    stop = false;
            if (stop || iterations == MAX_ITERATIONS)
                break;
            P = (double[])X.clone();
        }
    }

    public static void main(String[] args)
        throws IOException
    {
        PrintWriter writer
            = new PrintWriter(System.out, true);
        int n = 2, k = 1;
        double[][] M = new double[n][n + 1];
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n + 1; j++)
                M[i][j] = k++;
        }
        GFG gausSeidel = new GFG(M);
        if (!gausSeidel.makeDominant()) {

            // if it is found that a matrix cannot be
            // changed into diagonally dominant then we
            // return the message to the user

            writer.println(
                "The system isn't diagonally dominant: "
                + "The method cannot guarantee convergence.");
        }
        writer.println();
        gausSeidel.print();
        gausSeidel.solve();
    }
}

Output

The system isn't diagonally dominant: The method cannot guarantee convergence.

1.0 2.0 3.0 
4.0 5.0 6.0 
X0 = {3.0 -1.2000000000000002 }
X1 = {5.4 -3.1200000000000006 }
X2 = {9.240000000000002 -6.192000000000002 }
X3 = {15.384000000000004 -11.107200000000004 }
X4 = {25.21440000000001 -18.97152000000001 }
X5 = {40.94304000000002 -31.554432000000016 }
X6 = {66.10886400000004 -51.68709120000003 }
X7 = {106.37418240000007 -83.89934592000006 }
X8 = {170.79869184000012 -135.4389534720001 }
X9 = {273.8779069440002 -217.90232555520015 }
X10 = {438.8046511104003 -349.84372088832026 }
X11 = {702.6874417766405 -560.9499534213124 }
X12 = {1124.899906842625 -898.7199254740999 }
X13 = {1800.4398509481998 -1439.1518807585599 }
X14 = {2881.3037615171197 -2303.843009213696 }
X15 = {4610.686018427392 -3687.348814741914 }
X16 = {7377.697629483828 -5900.958103587062 }
X17 = {11804.916207174125 -9442.7329657393 }
X18 = {18888.4659314786 -15109.57274518288 }
X19 = {30222.14549036576 -24176.51639229261 }
X20 = {48356.03278458522 -38683.62622766818 }
X21 = {77370.25245533636 -61895.00196426909 }
X22 = {123793.00392853818 -99033.20314283055 }
X23 = {198069.4062856611 -158454.3250285289 }
X24 = {316911.6500570578 -253528.12004564624 }
X25 = {507059.2400912925 -405646.192073034 }
X26 = {811295.384146068 -649035.1073168544 }
X27 = {1298073.2146337088 -1038457.3717069671 }
X28 = {2076917.7434139343 -1661532.9947311475 }
X29 = {3323068.989462295 -2658453.991569836 }
X30 = {5316910.983139672 -4253527.586511738 }
X31 = {8507058.173023475 -6805645.338418781 }
X32 = {1.3611293676837562E7 -1.088903374147005E7 }
X33 = {2.17780704829401E7 -1.742245518635208E7 }
X34 = {3.484491337270416E7 -2.787592949816333E7 }
X35 = {5.575186199632666E7 -4.460148839706133E7 }
X36 = {8.920297979412267E7 -7.136238263529813E7 }
X37 = {1.4272476827059627E8 -1.1417981341647702E8 }
X38 = {2.2835962983295405E8 -1.8268770266636324E8 }
X39 = {3.653754083327265E8 -2.923003254661812E8 }
X40 = {5.846006539323624E8 -4.6768052194588995E8 }
X41 = {9.353610468917799E8 -7.48288836313424E8 }
X42 = {1.496577675626848E9 -1.1972621393014784E9 }
X43 = {2.394524281602957E9 -1.9156194240823655E9 }
X44 = {3.831238851164731E9 -3.064991079731785E9 }
X45 = {6.12998216246357E9 -4.903985728770856E9 }
X46 = {9.807971460541712E9 -7.84637716723337E9 }
X47 = {1.569275433746674E10 -1.2554203468773392E10 }
X48 = {2.5108406940546783E10 -2.0086725551237427E10 }
X49 = {4.017345110547485E10 -3.2138760883179886E10 }
X50 = {6.427752176935977E10 -5.142201741428782E10 }
X51 = {1.0284403483157564E11 -8.227522786406052E10 }
X52 = {1.6455045573112103E11 -1.3164036458369684E11 }
X53 = {2.6328072917039368E11 -2.1062458333511496E11 }
X54 = {4.212491666732299E11 -3.36999333337384E11 }
X55 = {6.73998666677768E11 -5.391989333410144E11 }
X56 = {1.0783978666850288E12 -8.627182933468231E11 }
X57 = {1.7254365866966462E12 -1.3803492693561172E12 }
X58 = {2.7606985387152344E12 -2.208558830970988E12 }
X59 = {4.417117661944976E12 -3.533694129554781E12 }
X60 = {7.067388259112562E12 -5.65391060728885E12 }
X61 = {1.13078212145807E13 -9.04625697166336E12 }
X62 = {1.809251394332972E13 -1.4474011154662576E13 }
X63 = {2.8948022309328152E13 -2.3158417847461324E13 }
X64 = {4.631683569492565E13 -3.705346855593932E13 }
X65 = {7.410693711188164E13 -5.928554968950412E13 }
X66 = {1.1857109937901123E14 -9.48568795032078E13 }
X67 = {1.897137590064186E14 -1.517710072051337E14 }
X68 = {3.035420144102704E14 -2.4283361152821512E14 }
X69 = {4.8566722305643325E14 -3.8853377844514544E14 }
X70 = {7.770675568902939E14 -6.216540455122339E14 }
X71 = {1.2433080910244708E15 -9.946464728195755E14 }
X72 = {1.989292945639154E15 -1.591434356511322E15 }
X73 = {3.182868713022647E15 -2.5462949704181165E15 }
X74 = {5.092589940836236E15 -4.0740719526689875E15 }
X75 = {8.148143905337978E15 -6.518515124270381E15 }
X76 = {1.3037030248540764E16 -1.042962419883261E16 }
X77 = {2.0859248397665224E16 -1.668739871813218E16 }
X78 = {3.3374797436264364E16 -2.6699837949011492E16 }
X79 = {5.3399675898022984E16 -4.2719740718418392E16 }
X80 = {8.5439481436836784E16 -6.8351585149469432E16 }
X81 = {1.36703170298938864E17 -1.09362536239151104E17 }
X82 = {2.18725072478302208E17 -1.74980057982641792E17 }
X83 = {3.4996011596528358E17 -2.7996809277222688E17 }
X84 = {5.5993618554445376E17 -4.4794894843556301E17 }
X85 = {8.9589789687112602E17 -7.1671831749690086E17 }
X86 = {1.43343663499380173E18 -1.14674930799504141E18 }
X87 = {2.29349861599008282E18 -1.8347988927920663E18 }
X88 = {3.6695977855841326E18 -2.9356782284673065E18 }
X89 = {5.871356456934613E18 -4.697085165547691E18 }
X90 = {9.394170331095382E18 -7.5153362648763064E18 }
X91 = {1.5030672529752613E19 -1.2024538023802092E19 }
X92 = {2.4049076047604183E19 -1.9239260838083346E19 }
X93 = {3.847852167616669E19 -3.0782817340933358E19 }
X94 = {6.1565634681866715E19 -4.925250774549338E19 }
X95 = {9.850501549098675E19 -7.880401239278941E19 }
X96 = {1.5760802478557882E20 -1.2608641982846306E20 }
X97 = {2.5217283965692612E20 -2.017382717255409E20 }
X98 = {4.034765434510818E20 -3.227812347608655E20 }
X99 = {6.45562469521731E20 -5.164499756173848E20 }

Time Complexity: O(N²)