2013-04-09 11 views
12

Zbudowałem mały kod, którego chcę użyć do rozwiązywania problemów wartości własnej z udziałem dużych rzadkich macierzy. Działa to dobrze, wszystko, co chcę teraz zrobić, to ustawić niektóre elementy w rzadkiej macierzy na zero, tj. Te w najwyższym rzędzie (co odpowiada implementacji warunków brzegowych). Mogę tylko dostosować wektory kolumn (C0, C1 i C2) poniżej, aby to osiągnąć. Zastanawiałem się jednak, czy istnieje bardziej bezpośredni sposób. Oczywiście indeksowanie NumPy nie działa z rzadkim pakietem SciPy.Jak zmienić elementy w rzadkiej macierzy w SciPy Pythona?

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

Jest to wizualizacja macierzy zbudowanej na podstawie powyższego kodu. Chcę tak ustawić elementy w pierwszym wierszu zero. enter image description here

+0

Znalazłem pracę. Format, którego używam (dia_matrix) nie jest dobry, jeśli chcę, chcę go osiągnąć. Zamiast tego użyję csr_matrix. Czy powinienem wtedy zamknąć ten post? – seb

+0

to dobrze napisane pytanie i może być przydatne dla innych osób w przyszłości. Co powiesz na umieszczenie tego, co znalazłeś jako odpowiedź? – YXD

+0

ok, zrobię to – seb

Odpowiedz

13

Jak zasugerowałem w komentarzach, opublikuję odpowiedź, którą znalazłem na własne pytanie. Istnieje kilka klas macierzy w rzadkim pakiecie SciPy, są one wymienione here. Można konwertować rzadkie macierze z jednej klasy na drugą. Więc co trzeba zrobić, to wybrać, aby przekształcić moje rzadki matrycy do csr_matrix klasy, po prostu

H = sp.csr_matrix(H) 

Następnie można ustawić elementy w pierwszym rzędzie do 0 za pomocą zwykłej notacji NumPy:

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

Dla kompletności, publikuję pełny zmodyfikowany fragment kodu poniżej.

#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

Jeśli masz więcej wierszy niż kolumn, csr jest szybkie, ale jeśli masz więcej kolumn niż wierszy, csc jest szybsze. –

Powiązane problemy