SOR iteration in numerical

1.7k Views Asked by At

In the following code I have made an SOR iteration in python. I seem to be getting the wrong output. It seems to be just outputting the b array when it should be solving the system. Is there something wrong with my function or call of my function? Here is the code below! I don't think my values of x should be all 1's and I don't think it should just take one iteration! When I run this code for smaller x, b, xo it works! but I am trying to run it for x which is a 10x10, b which is a 10x1 and xo which is a 10x1

import numpy as np 
import math 

x = np.array([[3.0, 1.0, 0., 0., 0., 0., 0., 0., 0., 0.],[1.0, 3.0, 1.0, 0., 0., 0., 0., 0., 0., 0.], [0., 1.0, 3.0, 1.0, 0., 0., 0., 0., 0., 0.], [0., 0, 1.0, 3.0, 1.0, 0., 0., 0., 0., 0.], [0., 0., 0., 1.0, 3.0, 1.0, 0., 0., 0., 0.], [0., 0., 0., 0., 1.0, 3.0, 1.0, 0., 0., 0.], [0., 0., 0., 0., 0., 1.0, 3.0, 1.0, 0., 0.], [0., 0., 0., 0., 0., 0., 1.0, 3.0, 1.0, 0.], [0., 0., 0., 0., 0., 0., 0., 1.0, 3.0, 1.0], [0., 0., 0., 0., 0., 0., 0., 0., 1.0, 3.0]])
b = np.array([1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0])
x0 = np.array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0.])
tol =  10 ** (-15)
max_iter = 20
w = 1.5

def SOR(A, b, x0, tol, max_iter, w): 
    if (w<=1 or w>2): 
        print('w should be inside [1, 2)'); 
        step = -1; 
        x = float('nan') 
        return 
    n = b.shape 
    x = x0 

    for step in range (1, max_iter): 
        for i in range(n[0]): 
            new_values_sum = np.dot(A[i, 1 : (i - 1)], x[1 : (i - 1)]) 
            for j in range(i + 1, n[0]): 
                old_values_sum = np.dot(A[i, j], x0[j]) 
            x[i] = b[i] - (old_values_sum + new_values_sum) / A[i, i] 
            x[i] = np.dot(x[i], w) + np.dot(x0[i], (1 - w))  

        if (np.linalg.norm(x - x0) < tol): 
            print(step) 
            break 

       x0 = x 


    print("X = {}".format(x)) 
    print("The number of iterations is: {}".format(step)) 
SOR(x, b, x0, tol, max_iter, w)

Which gives the following output, which is not what I am looking for

1
X = [1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
The number of iterations is: 1
1

There are 1 best solutions below

0
On BEST ANSWER

You set x = x0. This means that both variables reference the same object, changes to one are automatically changes to the other. You need to force a copy, this happens in one of the variants

x = 1.0*x0
x = 0.0+x0
x = x0.copy()

where the last is obviously the preferred version for a descriptive programming style.

Later you assign x0=x which has the same problem. Here you might not want to create a new copy and delete the previous object referenced by x0, so you can use

x0[:] = x

to force copying of the data into the existing object.