[Scipy-svn] r3687 - in trunk/scipy/sparse: . tests
scipy-svn at scipy.org
scipy-svn at scipy.org
Tue Dec 18 18:20:29 EST 2007
Author: wnbell
Date: 2007-12-18 17:20:23 -0600 (Tue, 18 Dec 2007)
New Revision: 3687
Modified:
trunk/scipy/sparse/base.py
trunk/scipy/sparse/compressed.py
trunk/scipy/sparse/csc.py
trunk/scipy/sparse/csr.py
trunk/scipy/sparse/data.py
trunk/scipy/sparse/tests/test_base.py
Log:
made __pow__ work like dense matrices
consolidated other functionality
Modified: trunk/scipy/sparse/base.py
===================================================================
--- trunk/scipy/sparse/base.py 2007-12-18 19:56:14 UTC (rev 3686)
+++ trunk/scipy/sparse/base.py 2007-12-18 23:20:23 UTC (rev 3687)
@@ -6,7 +6,7 @@
from numpy import asarray, asmatrix, asanyarray, ones
-from sputils import isdense, isscalarlike
+from sputils import isdense, isscalarlike, isintlike
# The formats that we might potentially understand.
@@ -229,14 +229,9 @@
def __div__(self, other):
# Always do true division
return self.__truediv__(other)
-
- def __pow__(self, other):
- csc = self.tocsc()
- return csc ** other
-
+
def __neg__(self):
- csc = self.tocsc()
- return -csc
+ return -self.tocsr()
def __iadd__(self, other):
raise NotImplementedError
@@ -253,6 +248,31 @@
def __itruediv__(self, other):
raise NotImplementedError
+ def __pow__(self, other):
+ if self.shape[0] != self.shape[1]:
+ raise TypeError,'matrix is not square'
+
+ if isintlike(other):
+ other = int(other)
+ if other < 0:
+ raise ValueError,'exponent must be >= 0'
+
+ if other == 0:
+ from construct import spidentity
+ return spidentity( self.shape[0], dtype=self.dtype )
+ elif other == 1:
+ return self.copy()
+ else:
+ result = self
+ for i in range(1,other):
+ result = result*self
+ return result
+ elif isscalarlike(other):
+ raise ValueError,'exponent must be an integer'
+ else:
+ raise NotImplementedError
+
+
def __getattr__(self, attr):
if attr == 'A':
return self.toarray()
Modified: trunk/scipy/sparse/compressed.py
===================================================================
--- trunk/scipy/sparse/compressed.py 2007-12-18 19:56:14 UTC (rev 3686)
+++ trunk/scipy/sparse/compressed.py 2007-12-18 23:20:23 UTC (rev 3687)
@@ -45,14 +45,16 @@
if rank(arg1) == 2:
from coo import coo_matrix
self.shape = arg1.shape
- self._set_self( self._tothis(coo_matrix(arg1)) )
+ self._set_self( self.__class__(coo_matrix(arg1)) )
else:
raise ValueError, "dense array must have rank 1 or 2"
elif isspmatrix(arg1):
- if copy:
+ if arg1.format == self.format and copy:
arg1 = arg1.copy()
- self._set_self( self._tothis(arg1) )
+ else:
+ arg1 = getattr(arg1,'to' + self.format)()
+ self._set_self( arg1 )
elif isinstance(arg1, tuple):
if isshape(arg1):
@@ -67,7 +69,7 @@
if len(arg1) == 2:
# (data, ij) format
from coo import coo_matrix
- other = self._tothis(coo_matrix(arg1, shape=shape))
+ other = self.__class__( coo_matrix(arg1, shape=shape) )
self._set_self( other )
elif len(arg1) == 3:
# (data, indices, indptr) format
@@ -181,13 +183,6 @@
"non-decreasing sequence'
-
- def __imul__(self, other): #self *= other
- if isscalarlike(other):
- self.data *= other
- return self
- else:
- raise NotImplementedError
def __add__(self,other):
# First check if argument is a scalar
@@ -271,29 +266,6 @@
else:
raise NotImplementedError
- def __itruediv__(self, other): #self *= other
- if isscalarlike(other):
- recip = 1.0 / other
- self.data *= recip
- return self
- else:
- raise NotImplementedError
-
- def __pow__(self, other):
- """ Element-by-element power (unless other is a scalar, in which
- case return the matrix power.)
- """
- if isscalarlike(other):
- return self._with_data(self.data**other)
- elif isspmatrix(other):
- if (other.shape != self.shape):
- raise ValueError, "inconsistent shapes"
-
- warn("use .multiply(other) for elementwise multiplication", \
- DeprecationWarning)
- return self._binopt(other,'_elmul_')
- else:
- raise NotImplementedError
def multiply(self, other):
"""Point-wise multiplication by another matrix
@@ -319,7 +291,7 @@
major_axis = self._swap((M,N))[0]
indptr = empty( major_axis + 1, dtype=intc )
- other = self._tothis(other)
+ other = self.__class__(other) #convert to this format
fn = getattr(sparsetools, self.format + '_matmat_pass1')
fn( M, N, self.indptr, self.indices, \
other.indptr, other.indices, \
@@ -420,7 +392,9 @@
def _get_slice(self, i, start, stop, stride, shape):
- """Returns a view of the elements [i, myslice.start:myslice.stop].
+ """Returns a copy of the elements
+ [i, start:stop:string] for row-oriented matrices
+ [start:stop:string, i] for column-oriented matrices
"""
if stride != 1:
raise ValueError, "slicing with step != 1 not supported"
@@ -433,17 +407,14 @@
if self.indices[ind] >= start and self.indices[ind] < stop:
indices.append(ind)
- index = self.indices[indices] - start
+ index = self.indices[indices] - start
data = self.data[indices]
indptr = numpy.array([0, len(indices)])
return self.__class__((data, index, indptr), shape=shape, \
dtype=self.dtype)
- def _transpose(self, cls, copy=False):
- M, N = self.shape
- return cls((self.data,self.indices,self.indptr),(N,M),copy=copy)
-
+ # conversion methods
def todia(self):
return self.tocoo(copy=False).todia()
@@ -472,8 +443,15 @@
from coo import coo_matrix
return coo_matrix((data,(row,col)), self.shape)
+
+ def toarray(self):
+ A = self.tocoo(copy=False)
+ M = zeros(self.shape, dtype=self.dtype)
+ M[A.row, A.col] = A.data
+ return M
-
+
+ # methods that modify the internal data structure
def sorted_indices(self):
"""Return a copy of this matrix with sorted indices
"""
@@ -483,7 +461,7 @@
# an alternative that has linear complexity is the following
# typically the previous option is faster
- #return self._toother()._toother()
+ #return self.toother().toother()
def sort_indices(self):
"""Sort the indices of this matrix *in place*
@@ -520,6 +498,38 @@
self.data = self.data[:self.nnz]
self.indices = self.indices[:self.nnz]
+
+ # needed by _data_matrix
+ def _with_data(self,data,copy=True):
+ """Returns a matrix with the same sparsity structure as self,
+ but with different data. By default the structure arrays
+ (i.e. .indptr and .indices) are copied.
+ """
+ if copy:
+ return self.__class__((data,self.indices.copy(),self.indptr.copy()), \
+ shape=self.shape,dtype=data.dtype)
+ else:
+ return self.__class__((data,self.indices,self.indptr), \
+ shape=self.shape,dtype=data.dtype)
+
+ # utility functions
+ def _binopt(self, other, op, in_shape=None, out_shape=None):
+ """apply the binary operation fn to two sparse matrices"""
+ other = self.__class__(other)
+
+ if in_shape is None:
+ in_shape = self.shape
+ if out_shape is None:
+ out_shape = self.shape
+
+ # e.g. csr_plus_csr, cscmucsc, etc.
+ fn = getattr(sparsetools, self.format + op + self.format)
+
+ indptr, ind, data = fn(in_shape[0], in_shape[1], \
+ self.indptr, self.indices, self.data,
+ other.indptr, other.indices, other.data)
+ return self.__class__((data, ind, indptr), shape=out_shape)
+
def _get_submatrix( self, shape0, shape1, slice0, slice1 ):
"""Return a submatrix of this matrix (new matrix is created)."""
def _process_slice( sl, num ):
@@ -564,35 +574,3 @@
data, indices, indptr = aux[2], aux[1], aux[0]
return data, indices, indptr, i1 - i0, j1 - j0
- # needed by _data_matrix
- def _with_data(self,data,copy=True):
- """Returns a matrix with the same sparsity structure as self,
- but with different data. By default the structure arrays
- (i.e. .indptr and .indices) are copied.
- """
- if copy:
- return self.__class__((data,self.indices.copy(),self.indptr.copy()), \
- shape=self.shape,dtype=data.dtype)
- else:
- return self.__class__((data,self.indices,self.indptr), \
- shape=self.shape,dtype=data.dtype)
-
- # utility functions
- def _binopt(self, other, op, in_shape=None, out_shape=None):
- """apply the binary operation fn to two sparse matrices"""
- other = self._tothis(other)
-
- if in_shape is None:
- in_shape = self.shape
- if out_shape is None:
- out_shape = self.shape
-
- # e.g. csr_plus_csr, cscmucsc, etc.
- fn = getattr(sparsetools, self.format + op + self.format)
-
- indptr, ind, data = fn(in_shape[0], in_shape[1], \
- self.indptr, self.indices, self.data,
- other.indptr, other.indices, other.data)
- return self.__class__((data, ind, indptr), shape=out_shape)
-
-
Modified: trunk/scipy/sparse/csc.py
===================================================================
--- trunk/scipy/sparse/csc.py 2007-12-18 19:56:14 UTC (rev 3686)
+++ trunk/scipy/sparse/csc.py 2007-12-18 23:20:23 UTC (rev 3687)
@@ -81,15 +81,16 @@
else:
return _cs_matrix.__getattr__(self, attr)
+ def transpose(self, copy=False):
+ from csr import csr_matrix
+ M,N = self.shape
+ return csr_matrix((self.data,self.indices,self.indptr),(N,M),copy=copy)
+
def __iter__(self):
csr = self.tocsr()
for r in xrange(self.shape[0]):
yield csr[r,:]
- def transpose(self, copy=False):
- from csr import csr_matrix
- return _cs_matrix._transpose(self, csr_matrix, copy)
-
def __getitem__(self, key):
#TODO unify these in _cs_matrix
if isinstance(key, tuple):
@@ -217,12 +218,6 @@
from csr import csr_matrix
return csr_matrix((data, indices, indptr), self.shape)
- def toarray(self):
- A = self.tocoo(copy=False)
- M = zeros(self.shape, dtype=self.dtype)
- M[A.row, A.col] = A.data
- return M
-
def get_submatrix( self, slice0, slice1 ):
"""Return a submatrix of this matrix (new matrix is created).
Contigous range of rows and columns can be selected using:
@@ -241,15 +236,7 @@
"""
return (x[1],x[0])
- def _otherformat(self):
- return "csr"
- def _toother(self):
- return self.tocsr()
-
- def _tothis(self, other):
- return other.tocsc()
-
from sputils import _isinstance
def isspmatrix_csc(x):
Modified: trunk/scipy/sparse/csr.py
===================================================================
--- trunk/scipy/sparse/csr.py 2007-12-18 19:56:14 UTC (rev 3686)
+++ trunk/scipy/sparse/csr.py 2007-12-18 23:20:23 UTC (rev 3687)
@@ -73,10 +73,6 @@
"""
-
-
-
-
def __getattr__(self, attr):
if attr == 'colind':
warn("colind attribute no longer in use. Use .indices instead",
@@ -87,7 +83,8 @@
def transpose(self, copy=False):
from csc import csc_matrix
- return _cs_matrix._transpose(self, csc_matrix, copy)
+ M,N = self.shape
+ return csc_matrix((self.data,self.indices,self.indptr),(N,M),copy=copy)
def __getitem__(self, key):
#TODO unify these in _cs_matrix
@@ -227,12 +224,6 @@
from csc import csc_matrix
return csc_matrix((data, indices, indptr), self.shape)
- def toarray(self):
- A = self.tocoo(copy=False)
- M = zeros(self.shape, dtype=self.dtype)
- M[A.row, A.col] = A.data
- return M
-
def get_submatrix( self, slice0, slice1 ):
"""Return a submatrix of this matrix (new matrix is created).
Contigous range of rows and columns can be selected using:
@@ -251,16 +242,7 @@
"""
return (x[0],x[1])
- def _otherformat(self):
- return "csc"
- def _toother(self):
- return self.tocsc()
-
- def _tothis(self, other):
- return other.tocsr()
-
-
from sputils import _isinstance
def isspmatrix_csr(x):
Modified: trunk/scipy/sparse/data.py
===================================================================
--- trunk/scipy/sparse/data.py 2007-12-18 19:56:14 UTC (rev 3686)
+++ trunk/scipy/sparse/data.py 2007-12-18 23:20:23 UTC (rev 3687)
@@ -9,6 +9,7 @@
__all__ = []
from base import spmatrix
+from sputils import isscalarlike
class _data_matrix(spmatrix):
def __init__(self):
@@ -32,6 +33,21 @@
def __neg__(self):
return self._with_data(-self.data)
+ def __imul__(self, other): #self *= other
+ if isscalarlike(other):
+ self.data *= other
+ return self
+ else:
+ raise NotImplementedError
+
+ def __itruediv__(self, other): #self /= other
+ if isscalarlike(other):
+ recip = 1.0 / other
+ self.data *= recip
+ return self
+ else:
+ raise NotImplementedError
+
def astype(self, t):
return self._with_data(self.data.astype(t))
Modified: trunk/scipy/sparse/tests/test_base.py
===================================================================
--- trunk/scipy/sparse/tests/test_base.py 2007-12-18 19:56:14 UTC (rev 3686)
+++ trunk/scipy/sparse/tests/test_base.py 2007-12-18 23:20:23 UTC (rev 3687)
@@ -159,7 +159,23 @@
denom = self.spmatrix(matrix([[1,0,0,4],[-1,0,0,0],[0,8,0,-5]],'d'))
res = matrix([[1,0,0,0.5],[-3,0,numpy.inf,0],[0,0.25,0,0]],'d')
assert_array_equal((self.datsp / denom).todense(),res)
+
+ def check_pow(self):
+ A = matrix([[1,0,2,0],[0,3,4,0],[0,5,0,0],[0,6,7,8]])
+ B = self.spmatrix( A )
+ for exponent in [0,1,2,3]:
+ assert_array_equal((B**exponent).todense(),A**exponent)
+
+ #invalid exponents
+ for exponent in [-1, 2.2, 1 + 3j]:
+ self.assertRaises( Exception, B.__pow__, exponent )
+
+ #nonsquare matrix
+ B = self.spmatrix(A[:3,:])
+ self.assertRaises( Exception, B.__pow__, 1 )
+
+
def check_rmatvec(self):
M = self.spmatrix(matrix([[3,0,0],[0,1,0],[2,0,3.0],[2,3,0]]))
assert_array_almost_equal([1,2,3,4]*M, dot([1,2,3,4], M.toarray()))
More information about the Scipy-svn
mailing list