[Scipy-svn] r3910 - in trunk/scipy/sparse: . tests
scipy-svn at scipy.org
scipy-svn at scipy.org
Sat Feb 9 16:02:39 EST 2008
Author: wnbell
Date: 2008-02-09 15:02:27 -0600 (Sat, 09 Feb 2008)
New Revision: 3910
Modified:
trunk/scipy/sparse/construct.py
trunk/scipy/sparse/tests/test_base.py
trunk/scipy/sparse/tests/test_construct.py
Log:
added sparse.kronsum and sparse.kron
deprecated sparse.spkron
Modified: trunk/scipy/sparse/construct.py
===================================================================
--- trunk/scipy/sparse/construct.py 2008-02-09 07:26:50 UTC (rev 3909)
+++ trunk/scipy/sparse/construct.py 2008-02-09 21:02:27 UTC (rev 3910)
@@ -2,7 +2,7 @@
"""
-__all__ = [ 'spdiags','speye','spidentity','spkron', 'bmat', 'lil_eye', 'lil_diags' ]
+__all__ = [ 'spdiags','speye','spidentity', 'spkron', 'kron', 'kronsum', 'bmat', 'lil_eye', 'lil_diags' ]
from itertools import izip
from warnings import warn
@@ -81,7 +81,7 @@
diags = ones((1, m), dtype = dtype)
return spdiags(diags, k, m, n).asformat(format)
-def spkron(A, B, format=None):
+def kron(A, B, format=None):
"""kronecker product of sparse matrices A and B
Parameters
@@ -100,13 +100,13 @@
>>> A = csr_matrix(array([[0,2],[5,0]]))
>>> B = csr_matrix(array([[1,2],[3,4]]))
- >>> spkron(A,B).todense()
+ >>> kron(A,B).todense()
matrix([[ 0, 0, 2, 4],
[ 0, 0, 6, 8],
[ 5, 10, 0, 0],
[15, 20, 0, 0]])
- >>> spkron(A,[[1,2],[3,4]]).todense()
+ >>> kron(A,[[1,2],[3,4]]).todense()
matrix([[ 0, 0, 2, 4],
[ 0, 0, 6, 8],
[ 5, 10, 0, 0],
@@ -159,78 +159,49 @@
return coo_matrix((data,(row,col)), shape=output_shape).asformat(format)
+def kronsum(A, B, format=None):
+ """kronecker sum of sparse matrices A and B
+ Kronecker sum of two sparse matrices is a sum of two Kronecker
+ products kron(I_n,A) + kron(B,I_m) where A has shape (m,m)
+ and B has shape (n,n) and I_m and I_n are identity matrices
+ of shape (m,m) and (n,n) respectively.
-
-def lil_eye((r,c), k=0, dtype='d'):
- """Generate a lil_matrix of dimensions (r,c) with the k-th
- diagonal set to 1.
-
Parameters
==========
- - r,c : int
- - row and column-dimensions of the output.
- - k : int
- - diagonal offset. In the output matrix,
- - out[m,m+k] == 1 for all m.
- - dtype : dtype
- - data-type of the output array.
+ A,B : squared dense or sparse matrices
+ format : format of the result (e.g. "csr")
+ - By default (format=None) an appropriate sparse matrix
+ format is returned. This choice is subject to change.
- """
- warn("lil_eye is deprecated. use speye(... , format='lil') instead", \
- DeprecationWarning)
- return speye(r,c,k,dtype=dtype,format='lil')
+ Returns
+ =======
+ kronecker sum in a sparse matrix format
+ Examples
+ ========
+
+ """
+ A = coo_matrix(A)
+ B = coo_matrix(B)
-#TODO remove this function
-def lil_diags(diags,offsets,(m,n),dtype='d'):
- """Generate a lil_matrix with the given diagonals.
+ if A.shape[0] != A.shape[1]:
+ raise ValueError('A is not square')
+
+ if B.shape[0] != B.shape[1]:
+ raise ValueError('B is not square')
- Parameters
- ==========
- - diags : list of list of values e.g. [[1,2,3],[4,5]]
- - values to be placed on each indicated diagonal.
- - offsets : list of ints
- - diagonal offsets. This indicates the diagonal on which
- the given values should be placed.
- - (r,c) : tuple of ints
- - row and column dimensions of the output.
- - dtype : dtype
- - output data-type.
+ dtype = upcast(A.dtype,B.dtype)
- Example
- =======
+ L = kron(spidentity(B.shape[0],dtype=dtype), A, format=format)
+ R = kron(B, spidentity(A.shape[0],dtype=dtype), format=format)
- >>> lil_diags([[1,2,3],[4,5],[6]],[0,1,2],(3,3)).todense()
- matrix([[ 1., 4., 6.],
- [ 0., 2., 5.],
- [ 0., 0., 3.]])
+ return (L+R).asformat(format) #since L + R is not always same format
- """
- offsets_unsorted = list(offsets)
- diags_unsorted = list(diags)
- if len(diags) != len(offsets):
- raise ValueError("Number of diagonals provided should "
- "agree with offsets.")
- sort_indices = numpy.argsort(offsets_unsorted)
- diags = [diags_unsorted[k] for k in sort_indices]
- offsets = [offsets_unsorted[k] for k in sort_indices]
- for i,k in enumerate(offsets):
- if len(diags[i]) < m-abs(k):
- raise ValueError("Not enough values specified to fill "
- "diagonal %s." % k)
- out = lil_matrix((m,n),dtype=dtype)
- for k,diag in izip(offsets,diags):
- for ix,c in enumerate(xrange(clip(k,0,n),clip(m+k,0,n))):
- out.rows[c-k].append(c)
- out.data[c-k].append(diag[ix])
- return out
-
-
def bmat( blocks, format=None, dtype=None ):
"""
Build a sparse matrix from sparse sub-blocks
@@ -332,3 +303,79 @@
return coo_matrix( (data, (row, col)), shape=shape ).asformat(format)
+
+#################################
+# Deprecated functions
+################################
+from numpy import deprecate
+
+spkron = deprecate(kron, oldname='spkron', newname='kron')
+
+def lil_eye((r,c), k=0, dtype='d'):
+ """Generate a lil_matrix of dimensions (r,c) with the k-th
+ diagonal set to 1.
+
+ Parameters
+ ==========
+ - r,c : int
+ - row and column-dimensions of the output.
+ - k : int
+ - diagonal offset. In the output matrix,
+ - out[m,m+k] == 1 for all m.
+ - dtype : dtype
+ - data-type of the output array.
+
+ """
+ warn("lil_eye is deprecated. use speye(... , format='lil') instead", \
+ DeprecationWarning)
+ return speye(r,c,k,dtype=dtype,format='lil')
+
+
+
+#TODO remove this function
+def lil_diags(diags,offsets,(m,n),dtype='d'):
+ """Generate a lil_matrix with the given diagonals.
+
+ Parameters
+ ==========
+ - diags : list of list of values e.g. [[1,2,3],[4,5]]
+ - values to be placed on each indicated diagonal.
+ - offsets : list of ints
+ - diagonal offsets. This indicates the diagonal on which
+ the given values should be placed.
+ - (r,c) : tuple of ints
+ - row and column dimensions of the output.
+ - dtype : dtype
+ - output data-type.
+
+ Example
+ =======
+
+ >>> lil_diags([[1,2,3],[4,5],[6]],[0,1,2],(3,3)).todense()
+ matrix([[ 1., 4., 6.],
+ [ 0., 2., 5.],
+ [ 0., 0., 3.]])
+
+ """
+ offsets_unsorted = list(offsets)
+ diags_unsorted = list(diags)
+ if len(diags) != len(offsets):
+ raise ValueError("Number of diagonals provided should "
+ "agree with offsets.")
+
+ sort_indices = numpy.argsort(offsets_unsorted)
+ diags = [diags_unsorted[k] for k in sort_indices]
+ offsets = [offsets_unsorted[k] for k in sort_indices]
+
+ for i,k in enumerate(offsets):
+ if len(diags[i]) < m-abs(k):
+ raise ValueError("Not enough values specified to fill "
+ "diagonal %s." % k)
+
+ out = lil_matrix((m,n),dtype=dtype)
+ for k,diag in izip(offsets,diags):
+ for ix,c in enumerate(xrange(clip(k,0,n),clip(m+k,0,n))):
+ out.rows[c-k].append(c)
+ out.data[c-k].append(diag[ix])
+ return out
+
Modified: trunk/scipy/sparse/tests/test_base.py
===================================================================
--- trunk/scipy/sparse/tests/test_base.py 2008-02-09 07:26:50 UTC (rev 3909)
+++ trunk/scipy/sparse/tests/test_base.py 2008-02-09 21:02:27 UTC (rev 3910)
@@ -17,14 +17,15 @@
import numpy
from numpy import arange, zeros, array, dot, ones, matrix, asmatrix, \
- asarray, vstack, ndarray, kron, transpose, diag
+ asarray, vstack, ndarray, transpose, diag
import random
from scipy.testing import *
+import scipy.sparse as sparse
from scipy.sparse import csc_matrix, csr_matrix, dok_matrix, \
coo_matrix, lil_matrix, dia_matrix, bsr_matrix, \
- speye, spkron, SparseEfficiencyWarning
+ speye, SparseEfficiencyWarning
from scipy.sparse.sputils import supported_dtypes
from scipy.splinalg import splu
@@ -84,12 +85,12 @@
mats.append( [[0,1],[0,2],[0,3]] )
mats.append( [[0,0,1],[0,0,2],[0,3,0]] )
- mats.append( kron(mats[0],[[1,2]]) )
- mats.append( kron(mats[0],[[1],[2]]) )
- mats.append( kron(mats[1],[[1,2],[3,4]]) )
- mats.append( kron(mats[2],[[1,2],[3,4]]) )
- mats.append( kron(mats[3],[[1,2],[3,4]]) )
- mats.append( kron(mats[3],[[1,2,3,4]]) )
+ mats.append( numpy.kron(mats[0],[[1,2]]) )
+ mats.append( numpy.kron(mats[0],[[1],[2]]) )
+ mats.append( numpy.kron(mats[1],[[1,2],[3,4]]) )
+ mats.append( numpy.kron(mats[2],[[1,2],[3,4]]) )
+ mats.append( numpy.kron(mats[3],[[1,2],[3,4]]) )
+ mats.append( numpy.kron(mats[3],[[1,2,3,4]]) )
for m in mats:
assert_equal(self.spmatrix(m).diagonal(),diag(m))
@@ -345,7 +346,7 @@
assert_equal( result, dot(a,b) )
def test_formatconversions(self):
- A = spkron([[1,0,1],[0,1,1],[1,0,0]], [[1,1],[0,1]] )
+ A = sparse.kron([[1,0,1],[0,1,1],[1,0,0]], [[1,1],[0,1]] )
D = A.todense()
A = self.spmatrix(A)
@@ -371,7 +372,7 @@
def test_tocompressedblock(self):
x = array([[1,0,2,0],[0,0,0,0],[0,0,4,5]])
y = array([[0,1,2],[3,0,5]])
- A = kron(x,y)
+ A = numpy.kron(x,y)
Asp = self.spmatrix(A)
for format in ['bsr']:
fn = getattr(Asp, 'to' + format )
@@ -1277,7 +1278,7 @@
data[3] = array([[ 0, 5, 10],
[15, 0, 25]])
- A = kron( [[1,0,2,0],[0,0,0,0],[0,0,4,5]], [[0,1,2],[3,0,5]] )
+ A = numpy.kron( [[1,0,2,0],[0,0,0,0],[0,0,4,5]], [[0,1,2],[3,0,5]] )
Asp = bsr_matrix((data,indices,indptr),shape=(6,12))
assert_equal(Asp.todense(),A)
@@ -1296,7 +1297,7 @@
assert_equal(bsr_matrix(A,blocksize=(2,2)).todense(),A)
assert_equal(bsr_matrix(A,blocksize=(2,3)).todense(),A)
- A = kron( [[1,0,2,0],[0,0,0,0],[0,0,4,5]], [[0,1,2],[3,0,5]] )
+ A = numpy.kron( [[1,0,2,0],[0,0,0,0],[0,0,4,5]], [[0,1,2],[3,0,5]] )
assert_equal(bsr_matrix(A).todense(),A)
assert_equal(bsr_matrix(A,shape=(6,12)).todense(),A)
assert_equal(bsr_matrix(A,blocksize=(1,1)).todense(),A)
@@ -1306,11 +1307,11 @@
assert_equal(bsr_matrix(A,blocksize=(3,12)).todense(),A)
assert_equal(bsr_matrix(A,blocksize=(6,12)).todense(),A)
- A = kron( [[1,0,2,0],[0,1,0,0],[0,0,0,0]], [[0,1,2],[3,0,5]] )
+ A = numpy.kron( [[1,0,2,0],[0,1,0,0],[0,0,0,0]], [[0,1,2],[3,0,5]] )
assert_equal(bsr_matrix(A,blocksize=(2,3)).todense(),A)
def test_eliminate_zeros(self):
- data = kron([1, 0, 0, 0, 2, 0, 3, 0], [[1,1],[1,1]]).T
+ data = numpy.kron([1, 0, 0, 0, 2, 0, 3, 0], [[1,1],[1,1]]).T
data = data.reshape(-1,2,2)
indices = array( [1, 2, 3, 4, 5, 6, 7, 8] )
indptr = array( [0, 3, 8] )
Modified: trunk/scipy/sparse/tests/test_construct.py
===================================================================
--- trunk/scipy/sparse/tests/test_construct.py 2008-02-09 07:26:50 UTC (rev 3909)
+++ trunk/scipy/sparse/tests/test_construct.py 2008-02-09 21:02:27 UTC (rev 3910)
@@ -1,6 +1,7 @@
"""test sparse matrix construction functions"""
-from numpy import array, matrix, kron
+import numpy
+from numpy import array, matrix
from scipy.testing import *
@@ -77,7 +78,7 @@
b = array([[1, 0, 0], [0, 1, 0], [0, 0, 1]], dtype='d')
assert_array_equal(a.toarray(), b)
- def test_spkron(self):
+ def test_kron(self):
cases = []
cases.append(array([[ 0]]))
@@ -96,11 +97,30 @@
for a in cases:
for b in cases:
- result = spkron(csr_matrix(a),csr_matrix(b)).todense()
- expected = kron(a,b)
+ result = kron(csr_matrix(a),csr_matrix(b)).todense()
+ expected = numpy.kron(a,b)
+ assert_array_equal(result,expected)
+ def test_kronsum(self):
+ cases = []
+
+ cases.append(array([[ 0]]))
+ cases.append(array([[-1]]))
+ cases.append(array([[ 4]]))
+ cases.append(array([[10]]))
+ cases.append(array([[1,2],[3,4]]))
+ cases.append(array([[0,2],[5,0]]))
+ cases.append(array([[0,2,-6],[8,0,14],[0,3,0]]))
+ cases.append(array([[1,0,0],[0,5,-1],[4,-2,8]]))
+
+ for a in cases:
+ for b in cases:
+ result = kronsum(csr_matrix(a),csr_matrix(b)).todense()
+ expected = numpy.kron(numpy.eye(len(b)), a) + \
+ numpy.kron(b, numpy.eye(len(a)))
assert_array_equal(result,expected)
+
def test_bmat(self):
A = coo_matrix([[1,2],[3,4]])
More information about the Scipy-svn
mailing list