Python Code for Guass-Seidel Method

Gauss-Seidel method is one of the common iterative methods to find solution of linear systems. Its working principle fairly simple, just same as normal algebraic way for finding values of unknown variables. However, in this method, when we express coefficients in a matrix form, diagonal dominant matrix is necessary to make sure our numerical solution converges. Hence, below you will see function to check diagonal dominance of matrix and function for Gauss-Seidel method. Lets go!

#Check if the matrix is diagonal dominant
def check_diaDom(A):
    r = len(A)
    # Always good idea to see if the code is working the way it should be
    # so, I have included print option to check throughout the process
    diagnl = 0
    non_diagnl = 0

    for i in range(r):
        for j in range(r):
            # print(f'Loop i = {i}, j = {j}\n')
            # print(A[i, j])
            if i ==j:
             diagnl+=abs(A[i] [j])
            elif i !=j: #here you can simply use else rather than else if
             non_diagnl+=abs(A[i] [j])
    
    # print('Sum of diagonal elements\n', diagnl)
    # print('Sum of non-diagonal elements\n', non_diagnl)

    if diagnl >=non_diagnl:
       verdict = True
       #print('Diagonal dominant Matrix')
    else:
       #print('Non-diagonal dominant Matrix')
       verdict = False
    return verdict

def GS_method(A, Y, X):
    n = len(A) # to find no of rows or col of sq matrix A
    for j in range(0, n):
        #old_val = X[i]
        summ_val = Y[j]
        for i in range(0, n):
            if (j!=i):
                summ_val-=A[j] [i] * X[i]
                #print(summ_val)
        X[j] = summ_val/A[j][j]
        #print(X[i]) 
    return X

Now, lets test our function! For example I am using set of equations given below.

Example equation sets

A = [[8, 3, -3],[-2, -8, 5],[3, 5, 10]] #coefficient of unknown variables
Y = [14, 5, -8] #terms that appears after '=' sign in the equation
X = [0, 0, 0] #initial guess for our method

if check_diaDom(A)==1:
 for i in range(1000):
    X = GS_method(A, Y, X)
    print(X)
else: 
   print('Non-Diagonal matrix')

On successful running the code above, you can expect output as given below.

The method we have just seen is more traditional way of doing but numpy package already has easier function to solve this like linear equation system.

Get kafleZ’s stories in your inbox

Join Medium for free to get updates from this writer.

Remember me for faster sign in

We just need to use numpy’s linear algebra solve function. For example:

#Solving linear algebraic equations using numpy

import numpy as np

A = np.array([[8, 3, -3], [-2, -8, 5], [3, 5, 10]])
y = np.array([14, 5, -8])

x = np.linalg.solve(A, y)
print(x)

Sample output

I hope this post is helpful! Thank you!