Write a program to solve a linear system of the form Ax b b

Write a program to solve a linear system of the form Ax = b by Gauss-Seidel iteration. You should pass the matrix A. solution vector b, the initial guess vector x, the size of the system N, and the relative error tolerance in the solution. Overwrite the x vector with the solution. Make sure you have a bound on the number of iterations and return an informative message if the maximum number of iterations has been exceeded. Return the solution and the number of iterations it took to converge. Test your program on the system in exercise #7 above to verify the accuracy of your code. Run your program on the system of equations 9x_1 - 3x_2 = 6 -2X_1 + 8X_2 = -4 with starting vector x = (0 0)^T. Change the coefficient of x_1 in the first equation from 9 to 1 see what happens to your program.

Solution

CODE:

import numpy as np

N = 1000
#Initialize the given coefficients of equation
A = np.array([[9,-3],
[-2,8]])
#Initialize the RHS coefficients
b = np.array([6,-45])

for i in range(A.shape[0]):

row = [\"{}*x{}\".format(A[i, j], j + 1) for j in range(A.shape[1])]
print(\" + \".join(row), \"=\", b[i])
x = np.zeros_like(b)
for it_count in range(N):
print(\"Current solution:\", x)
x_new = np.zeros_like(x)

for i in range(A.shape[0]):
s1 = np.dot(A[i, :i], x_new[:i])
s2 = np.dot(A[i, i + 1:], x[i + 1:])
x_new[i] = (b[i] - s1 - s2) / A[i, i]

if np.allclose(x, x_new, rtol=1e-8):
break

x = x_new

print(\"Solution:\")
print(x)

OUTPUT:

(\'9*x1 + -3*x2\', \'=\', 6)
(\'-2*x1 + 8*x2\', \'=\', -45)
(\'Current solution:\', array([0, 0]))
(\'Current solution:\', array([ 0, -6]))
(\'Current solution:\', array([-2, -7]))
Solution:
[-2 -7]

EXPLANATION:

First for loop is formating the given array into this fromat and displaying :

(\'9*x1 + -3*x2\', \'=\', 6)

Then, we are taking Current solution as 0,0 and iterating till the final solution is achieved.

We are obtiaining the dot product by using method dot by passing both A array and in each iteration x_new array.

And then putting that in the final equation of the Guass-Seidel Method

x_new[i] = (b[i] - s1 - s2) / A[i, i]

And printing the Current Result each time.

We are checking with allclose() method that two arrays are equal within a tolerance , if so it returns true

Then we are breaking out of the nested loop and printing the Final result.

A.shape[0] : shape attribute returns the dimensions of the array.If A has n rows and m columns, then Y.shape is (n,m). So A.shape[0] is n.

x = np.zeros_like(b): Return an array of zeros with the same shape and type as a given array(b).