error with matrix in python (method of gauss-jordan) (numpy)

2

Good day estimated, implement the gauss-jordan method in python in the following way and it is working:

#La idea de este metodo es que el usuario ingrese una matriz MxM y un vector de tamano M
import numpy
def gaussJordan(m):
    #creamos una matriz, un vector lleno de ceros, y el vector solucion con la misma cantidad de zeros
    matrix = numpy.zeros((m,m))
    vector = numpy.zeros((m))
    x = numpy.zeros((m))
    #se llena la matriz y el vector
    for r in range(0, m):
        for c in range(0, m):
            matrix[(r), (c)] =(input("Elemento a[" + str(r+1)+"," +str(c+1)+"] "))
        vector[(r)]=(input('b[]'+str(r+1)+']: '))
    #asi funciona el metodo
    for k in range(0, m):
        for r in range(k+1, m):
            factor=(matrix[r,k]/matrix[k,k])
            vector[r]=vector[r]-(factor*vector[k])
            for c in range(0,m):
                matrix[r,c]=matrix[r,c]-(factor*matrix[k,c])

    x[m-1]=vector[m-1]/matrix[m-1, m-1]

    for r in range(m-2, -1, -1):
        suma = 0
        for c in range(0,m):
            suma=suma+matrix[r,c]*x[c]
        x[r]=(vector[r]-suma)/matrix[r, r]  
    return x

m = int(input('Valor de m:'))
print(gaussJordan(m))

The problem is that a colleague asked me if it would be possible that instead of entering m in the method and then building the matrix, I could enter the matrix and the vector already made, and it seemed like a simple change but I came across a couple of errors that I do not understand the truth, first when entering the matrix [[2,6,1], [1,2, -1], [5,7, -4]] I returned an error in which apparently did not recognize it as a matrix so use numpy.matrix and the same with the vector just in case numpy.array. After this I had no problems in executing the method but for some reason I returned a different result than it should be ... when I execute it by filling the matrix within the method, if it returns the result it should be, it should return [10, - 3, 5] but it returns me [21.5, -8, 12] if someone knows how to tell me why I really appreciate it, I do not think it's an algorithm problem, because it would throw an incorrect result in which if it asks to fill the matrix internally, I think what is the problem of how the matrix is reading the program, THE DATA I USE FOR THE TEST ARE MATRIX = [[2,6,1], [1,2, -1], [5,7, -4]] VECTOR = [7, -1,9] here is the code:

#La idea de este metodo es que el usuario ingrese una matriz MxM y un vector de tamano M
import numpy
def gaussJordan(matriz, vector):

    matrix = numpy.matrix(matriz)
    vector = numpy.array(vector)
    m = len(vector)
    x = numpy.zeros((m))
    #se llena la matriz y el vector
##    for r in range(0, m):
##        for c in range(0, m):
##            matrix[(r), (c)] =(input("Elemento a[" + str(r+1)+"," +str(c+1)+"] "))
##        vector[(r)]=(input('b[]'+str(r+1)+']: '))

    #asi funciona el metodo
    for k in range(0, m):
        for r in range(k+1, m):
            factor=(matrix[r,k]/matrix[k,k])
            vector[r]=vector[r]-(factor*vector[k])
            for c in range(0,m):
                matrix[r,c]=matrix[r,c]-(factor*matrix[k,c])

    x[m-1]=vector[m-1]/matrix[m-1, m-1]
    #print (x[m-1])

    for r in range(m-2, -1, -1):
        suma = 0
        for c in range(0,m):
            suma=suma+matrix[r,c]*x[c]
        x[r]=(vector[r]-suma)/matrix[r, r]  
    return x

print(gaussJordan([[2,6,1],[1,2,-1],[5,7,-4]],[7,-1,9]))
    
asked by felpex 30.05.2017 в 12:33
source

1 answer

2

The problem is that NumPy has static and homogeneous typing. I mean, if you declare an array of integers this can only store integers and if it's float32 you can only store 32 bit floats. This is important because a Python list can store any variety of objects even if they are of different types within the same list and change the type of a given element at any time (it really only saves objects without caring what they are).

When you use input your arrays vector and matrix are created by numpy.zeros that returns a array of float 64 bits . During your algorithm, both vector and matrix may end up having to store floats, in this case there is no problem and the method works.

When you use the list to build both arrays an array of a given type is created according to the type of data contained in the input lists , that is, two integer arrays are created in this case. The problem is that these arrays can not store decimals so if you do something like:

vector[0] = 4/3 

vector[0] does not store 1.33333 as we can expect, but only the whole part (1). It is limited to storing the whole part but does not throw an exception for the difference in types. This causes important cumulative errors in the algorithm and you end up with important errors in the result.

The solution is simple, force the type of the arrays to be floats even if the input lists are integers using the argument dtype :

import numpy
def gaussJordan(matriz, vector):

    matrix = numpy.array(matriz, dtype=numpy.float64)
    vector = numpy.array(vector, dtype=numpy.float64)

    m = len(vector)
    x = numpy.zeros(m)

    for k in range(0, m):
        for r in range(k+1, m):
            factor=(matrix[r,k]/matrix[k,k])
            vector[r]=vector[r]-(factor*vector[k])
            for c in range(0,m):
                matrix[r,c]=matrix[r,c]-(factor*matrix[k,c])

    x[m-1]=vector[m-1]/matrix[m-1, m-1]

    for r in range(m-2, -1, -1):
        suma = 0
        for c in range(0,m):
            suma=suma+matrix[r,c]*x[c]
        x[r]=(vector[r]-suma)/matrix[r, r]  
    return x

print(gaussJordan([[2,6,1],[1,2,-1],[5,7,-4]],[7,-1,9]))

Exit:

  

[10. -3. 5.]

    
answered by 30.05.2017 / 13:48
source