python - transpuesta - metodo de gauss jordan

Python función incorporada para hacer la reducción de matriz

Question

¿Python tiene una función incorporada que convierte una matriz en una forma escalonada (también conocida como triangular superior)?

Answer 1

Estoy de acuerdo con un comentario de @Mile para @WinstonEwert No hay ninguna razón para que una computadora no pueda ejecutar RREF con una precisión dada .

La realización de RREF no debería ser muy complicada, y matlab de alguna manera logró tener esta función, por lo que el número también debería tenerla.

Hice una realización hacia adelante muy simple y estatal, de causa muy ineficiente. Pero para matrices simples funciona bastante bien:

def rref(mat,precision=0,GJ=False):
    m,n = mat.shape
    p,t = precision, 1e-1**precision
    A = around(mat.astype(float).copy(),decimals=p )
    if GJ:
        A = hstack((A,identity(n)))
    pcol = -1 #pivot colum
    for i in xrange(m):
        pcol += 1
        if pcol >= n : break
        #pivot index
        pid = argmax( abs(A[i:,pcol]) )
        #Row exchange
        A[i,:],A[pid+i,:] = A[pid+i,:].copy(),A[i,:].copy()
        #pivot with given precision
        while pcol < n and abs(A[i,pcol]) < t:
            #pivot index
            pid = argmax( abs(A[i:,pcol]) )
            #Row exchange
            A[i,:],A[pid+i,:] = A[pid+i,:].copy(),A[i,:].copy()
            pcol += 1
        if pcol >= n : break
        pivot = float(A[i,pcol])
        for j in xrange(m):
            if j == i: continue
            mul = float(A[j,pcol])/pivot
            A[j,:] = around(A[j,:] - A[i,:]*mul,decimals=p)
        A[i,:] /= pivot
        A[i,:] = around(A[i,:],decimals=p)

    if GJ:
        return A[:,:n].copy(),A[:,n:].copy()
    else:
        return A   

Aquí algunas pruebas sencillas.

print "/*--------------------------------------/"
print "/             Simple TEST               /"
print "/--------------------------------------*/"
A = array([[1,2,3],[4,5,6],[-7,8,9]])
print "A:/n",R
R = rref(A,precision=6)
print "R:/n",R
print
print "With GJ "
R,E =   rref(A,precision=6,GJ=True)
print "R:/n",R
print "E:/n",E
print "AdotE:/n",around( dot(A,E),decimals=0)

/*--------------------------------------/
/             Simple TEST               /
/--------------------------------------*/
A:
[[ 1.  0.  1.]
 [-0.  1.  1.]
 [ 0.  0.  0.]
 [ 0.  0.  0.]]
R:
[[ 1.  0.  0.]
 [ 0.  1.  0.]
 [ 0.  0.  1.]]

With GJ 
R:
[[ 1.  0.  0.]
 [ 0.  1.  0.]
 [ 0.  0.  1.]]
E:
[[-0.071428  0.142857 -0.071429]
 [-1.857142  0.714285  0.142857]
 [ 1.595237 -0.523809 -0.071428]]
AdotE:
[[ 1.  0.  0.]
 [ 0.  1.  0.]
 [-0.  0.  1.]]

print "/*--------------------------------------/"
print "/        Not Invertable TEST            /"
print "/--------------------------------------*/"
A = array([
    [2,2,4, 4],
    [3,1,6, 2],
    [5,3,10,6]])
print "A:/n",A
R = rref(A,precision=2)
print "R:/n",R
print
print "A^{T}:/n",A.T
R = rref(A.T,precision=10)
print "R:/n",R
/*--------------------------------------/
/        Not Invertable TEST            /
/--------------------------------------*/
A:
[[ 2  2  4  4]
 [ 3  1  6  2]
 [ 5  3 10  6]]
R:
[[ 1.  0.  2.  0.]
 [-0.  1. -0.  2.]
 [ 0.  0.  0.  0.]]

A^{T}:
[[ 2  3  5]
 [ 2  1  3]
 [ 4  6 10]
 [ 4  2  6]]
R:
[[ 1.  0.  1.]
 [-0.  1.  1.]
 [ 0.  0.  0.]
 [ 0.  0.  0.]]

Answer 2

Sí. En scipy.linalg , lu hace la descomposición de LU que esencialmente le dará forma de escalón de fila.

Hay otras factorizaciones, como qr , rq , svd y más, si está interesado.

Documentation.

Answer 3

ver http://mail.scipy.org/pipermail/numpy-discussion/2008-November/038705.html

Básicamente: no lo hagas.

El algoritmo rref produce demasiada inexactitud cuando se implementa en una computadora. Entonces, o bien quiere resolver el problema de otra manera, o usar símbolos como @aix sugirió.

Answer 4

Si puedes usar sympy , Matrix.rref() puede hacerlo:

In [8]: sympy.Matrix(np.random.random((4,4))).rref()
Out[8]: 
([1, 1.42711055402454e-17, 0, -1.38777878078145e-17]
[0,                  1.0, 0,  2.22044604925031e-16]
[0, -2.3388341405089e-16, 1, -2.22044604925031e-16]
[0, 3.65674099486992e-17, 0,                   1.0],
 [0, 1, 2, 3])