Natürlich kann ich auch "linalg.svc(matrix)" benutzen und die Anzahl der nicht-zu-kleinen Diagonalelemente zählen. Aber das wäre ein Riesenaufwand bei dieser Anzahl Matrizen, da
(1) die Berechnung singulärer Werte überhaupt nicht benötigt wird, und
(2) numpy mit float arbeitet und wegen der Rundungsfehler massenhaft Werte unterhalb der Toleranzgrenze findet.
Hat jemand eine Idee, wie man den Rang wirklich *schnell* berechnen kann?
NACHTRAG (5.10.2011):
Ein Algorithmus zur Rankberechnung ist im wesentlichen dasselbe wie ein Algorithmus zur Determinantenberechnung, der die Matrix wie bei einer LU-Zerlegung zeilenweise (oder spaltenweise) reduziert und dann den Schritt ausgibt, bei welchem dieser Vorgang abbricht (dieser Schritt ist der Rang). Leider ist diese einfache Erweiterung bei numpy.linalg.det() nicht vorgesehen.
Ich habe inzwischen die folgende ganzzahlige Rangsuche selber programmiert, aber sie ist für 0-1-Matrizen der Ordnung 27x12 ungefähr 5mal-10mal so langsam wie meine (weiter unten ebenfalls aufgelistete) svd-basierte Alternative. Ausgabe ist (rank,quasidet), wobei quasidet nur bei quadratischen Matrizen der Determinant dieser Matrix ist, bei nichtquadratischen Matrizen ist quasidet beliebig.
Code: Alles auswählen
import numpy as sc
def irank_direkt(matrix):
matnum = sc.matrix(matrix)
assert matnum.dtype==int
#
quasidet = 1
(zlz,spz) = matnum.shape
for schritt in range(min(zlz,spz)):
#Elemente der Matrix immer *zeilenweise* absuchen:
for (pos,wert) in enumerate(matnum.flat):
if not wert==0: break
else: return (schritt,0)
(zl,sp) = divmod(pos,spz) #(zl,sp)==(pos//spz,pos%spz)
#
#(hier sollten nichtverschwindende Elemente bereits gefunden werden,
#aber dann verzichtet man auf die numpy-Matrixoperationen)
xpiv = matnum[zl,sp]
matnum = (xpiv*matnum - matnum[:,sp]*matnum[zl,:]) // quasidet
quasidet = xpiv
#
#(folgende Operation ist anscheinend sehr langsam)
matnum = sc.hstack( (matnum[zl+1:,:sp], matnum[zl+1:,sp+1:]) )
(zlz,spz) = matnum.shape
else: return (schritt+1,quasidet)
Code: Alles auswählen
import numpy as sc
def irank_svd(matrix, tol=1.e-9):
matrix = sc.matrix(matrix)
sv_lst = sc.linalg.svd(matrix, compute_uv=False)
rank = len([sv for sv in sv_lst if sv>tol])
return rank