python fonction intégrée pour faire de la matrice de réduction

Si vous pouvez utiliser sympy, de la Matrice.rref() peut le faire:

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])

OriginalL'auteur

Je suis d'accord avec @Mile commentaire @WinstonEwert répondre Il n'y a pas de raison un ordinateur n'a pas pu effectuer RREF avec précision donnée.

La réalisation de RREF ne devrait pas être très compliqué, et matlab réussi à avoir cette fonction, de sorte numpy devrait avoir trop.

J'ai fait un très simple et état de l'avant réalisation, de provoquer de très innefficient. Mais pour de simples matrices, ça fonctionne très 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   

Voici quelques tests simples

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.]]

OriginalL'auteur