2013-04-09 9 views
12

Ho creato un piccolo codice che voglio utilizzare per risolvere problemi agli autovalori che coinvolgono matrici sparse di grandi dimensioni. Funziona bene, tutto quello che voglio fare ora è impostare alcuni elementi nella matrice sparsa a zero, cioè quelli nella riga più in alto (che corrisponde all'implementazione delle condizioni al contorno). Posso semplicemente regolare i vettori delle colonne (C0, C1 e C2) sotto per ottenere quello. Tuttavia, mi sono chiesto se c'è un modo più diretto. Evidentemente, l'indicizzazione NumPy non funziona con il pacchetto sparse di SciPy.Come modificare gli elementi nella matrice sparsa in Python's SciPy?

import scipy.sparse as sp 
import scipy.sparse.linalg as la 
import numpy as np 
import matplotlib.pyplot as plt 

#discretize x-axis 
N = 11 
x = np.linspace(-5,5,N) 
print(x) 
V = x * x/2 
h = len(x)/(N) 
hi2 = 1./(h**2) 
#discretize Schroedinger Equation, i.e. build 
#banded matrix from difference equation 
C0 = np.ones(N)*30. + V 
C1 = np.ones(N) * -16. 
C2 = np.ones(N) * 1. 
diagonals = np.array([-2,-1,0,1,2]) 
H = sp.spdiags([C2, C1, C0,C1,C2],[-2,-1,0,1,2], N, N) 
H *= hi2 * (- 1./12.) * (- 1./2.) 
#solve for eigenvalues 
EV = la.eigsh(H,return_eigenvectors = False) 

#check structure of H 
plt.figure() 
plt.spy(H) 
plt.show() 

Questa è una visualizzazione della matrice creata dal codice sopra. Voglio quindi impostare gli elementi nella prima riga zero. enter image description here

+0

Ho trovato un lavoro in giro. Il formato che sto usando (dia_matrix) non va bene per il desiderio che voglio raggiungere. Userò invece csr_matrix. Dovrei chiudere questo post, allora? – seb

+0

è una domanda ben scritta e potrebbe essere utile ad altri in futuro. Che ne dici di pubblicare ciò che hai trovato come risposta? – YXD

+0

ok, lo farò – seb

risposta

13

Come suggerito nei commenti, posterò la risposta che ho trovato alla mia stessa domanda. Esistono diverse classi di matrici nel pacchetto sparse di SciPy, sono elencate nello here. Si possono convertire matrici sparse da una classe all'altra. Così, per quello che ho bisogno di fare, ho scelto di convertire il mio matrici sparse al csr_matrix di classe, semplicemente

H = sp.csr_matrix(H) 

Allora posso impostare gli elementi della prima fila a 0 utilizzando la normale notazione NumPy:

H[0,0] = 0 
H[0,1] = 0 
H[0,2] = 0 

Per completezza, inserisco di seguito il frammento di codice completo modificato.

#SciPy Sparse linear algebra takes care of sparse matrix computations 
#http://docs.scipy.org/doc/scipy/reference/sparse.linalg.html 
import scipy.sparse as sp 
import scipy.sparse.linalg as la 

import numpy as np 
import matplotlib.pyplot as plt 

#discretize x-axis 
N = 1100 
x = np.linspace(-100,100,N) 
V = x * x/2. 
h = len(x)/(N) 
hi2 = 1./(h**2) 

#discretize Schroedinger Equation, i.e. build 
#banded matrix from difference equation 
C0 = np.ones(N)*30. + V 
C1 = np.ones(N) * -16. 
C2 = np.ones(N) * 1. 

H = sp.spdiags([C2, C1, C0, C1, C2],[-2,-1,0,1,2], N, N) 
H *= hi2 * (- 1./12.) * (- 1./2.) 
H = sp.csr_matrix(H) 
H[0,0] = 0 
H[0,1] = 0 
H[0,2] = 0 

#check structure of H 
plt.figure() 
plt.spy(H) 
plt.show() 

EV = la.eigsh(H,return_eigenvectors = False) 
+1

Se hai più righe che colonne, csr è veloce, ma se hai più colonne che file, csc è più veloce. –