Merge pull request #524 from nouiz/sparse

doc/library/sparse/sandbox.txt

@@ -17,5 +17,7 @@ API
 
 .. automodule:: theano.sparse.sandbox.sp
     :members:
+.. automodule:: theano.sparse.sandbox.sp2
+    :members:
 .. automodule:: theano.sparse.sandbox.truedot
     :members:

theano/sandbox/cuda/tests/test_var.py

@@ -16,12 +16,13 @@ from theano.sandbox.cuda import CudaNdarrayType, cuda_available
 if cuda_available == False:
     raise SkipTest('Optional package cuda disabled')
 
+
 def test_float32_shared_constructor():
 
-    npy_row = numpy.zeros((1,10), dtype='float32')
+    npy_row = numpy.zeros((1, 10), dtype='float32')
 
-    def eq(a,b):
-        return a==b
+    def eq(a, b):
+        return a == b
 
     # test that we can create a CudaNdarray
     assert (f32sc(npy_row).type == CudaNdarrayType((False, False)))
@@ -40,37 +41,41 @@ def test_float32_shared_constructor():
 
     # test that we can make non-matrix shared vars
     assert eq(
-            f32sc(numpy.zeros((2,3,4,5), dtype='float32')).type,
-            CudaNdarrayType((False,)*4))
+            f32sc(numpy.zeros((2, 3, 4, 5), dtype='float32')).type,
+            CudaNdarrayType((False,) * 4))
+
 
 def test_givens():
     # Test that you can use a TensorType expression to replace a
     # CudaNdarrayType in the givens dictionary.
     # This test case uses code mentionned in #757
-    data = numpy.float32([1,2,3,4])
+    data = numpy.float32([1, 2, 3, 4])
     x = f32sc(data)
-    y = x**2
-    f = theano.function([], y, givens={x:x+1})
+    y = x ** 2
+    f = theano.function([], y, givens={x: x + 1})
+    f()
+
 
 class T_updates(unittest.TestCase):
     # Test that you can use a TensorType expression to update a
     # CudaNdarrayType in the updates dictionary.
 
     def test_1(self):
-        data = numpy.float32([1,2,3,4])
+        data = numpy.float32([1, 2, 3, 4])
         x = f32sc(data)
-        y = x**2
-        f = theano.function([], y, updates={x:x+1})
+        y = x ** 2
+        f = theano.function([], y, updates={x: x + 1})
+        f()
 
     def test_2(self):
         # This test case uses code mentionned in #698
-        data = numpy.random.rand(10,10).astype('float32')
+        data = numpy.random.rand(10, 10).astype('float32')
         output_var = f32sc(name="output",
-                value=numpy.zeros((10,10), 'float32'))
+                value=numpy.zeros((10, 10), 'float32'))
 
         x = tensor.fmatrix('x')
-        output_updates = {output_var:x**2}
-        output_givens = {x:data}
+        output_updates = {output_var: x ** 2}
+        output_givens = {x: data}
         output_func = theano.function(inputs=[], outputs=[],
                 updates=output_updates, givens=output_givens)
         output_func()
@@ -78,16 +83,16 @@ class T_updates(unittest.TestCase):
     def test_3(self):
         # Test that broadcastable dimensions don't screw up
         # update expressions.
-        data = numpy.random.rand(10,10).astype('float32')
-        output_var = f32sc(name="output",
-                value=numpy.zeros((10,10), 'float32'))
+        data = numpy.random.rand(10, 10).astype('float32')
+        output_var = f32sc(name="output", value=data)
 
         # the update_var has type matrix, and the update expression
         # is a broadcasted scalar, and that should be allowed.
         output_func = theano.function(inputs=[], outputs=[],
-                updates={output_var:output_var.sum().dimshuffle('x', 'x')})
+                updates={output_var: output_var.sum().dimshuffle('x', 'x')})
         output_func()
 
+
 class T_ifelse(unittest.TestCase):
     def setUp(self):
         utt.seed_rng()
@@ -111,10 +116,10 @@ class T_ifelse(unittest.TestCase):
         f = theano.function([cond], out1)
         g = theano.function([cond], out2)
 
-        assert numpy.all(f(0) == data+1)
+        assert numpy.all(f(0) == data + 1)
         assert numpy.all(f(1) == data)
         assert numpy.all(g(0) == data)
-        assert numpy.all(g(1) == data+1)
+        assert numpy.all(g(1) == data + 1)
 
     def test_dtype_mismatch(self):
         data = self.rng.rand(5).astype('float32')
@@ -135,7 +140,7 @@ class T_ifelse(unittest.TestCase):
         self.assertRaises(TypeError, ifelse, cond, y, x)
 
     def test_broadcast_mismatch(self):
-        data = self.rng.rand(2,3).astype('float32')
+        data = self.rng.rand(2, 3).astype('float32')
         x = f32sc(data)
         print x.broadcastable
         y = tensor.frow('y')
@@ -146,7 +151,7 @@ class T_ifelse(unittest.TestCase):
         self.assertRaises(TypeError, ifelse, cond, y, x)
 
     def test_sparse_tensor_error(self):
-        data = self.rng.rand(2,3).astype('float32')
+        data = self.rng.rand(2, 3).astype('float32')
         x = f32sc(data)
         y = sparse.matrix('csc', dtype='float32', name='y')
         z = sparse.matrix('csr', dtype='float32', name='z')
@@ -160,5 +165,3 @@ class T_ifelse(unittest.TestCase):
         self.assertRaises((TypeError, ValueError), ifelse, cond, z, x)
         self.assertRaises((TypeError, ValueError), ifelse, cond, y, z)
         self.assertRaises((TypeError, ValueError), ifelse, cond, z, y)
-
-

theano/sparse/sandbox/sp2.py

-from theano.sparse.basic import * # To facilitate later merge into sparse module
+import numpy
+
+from theano import gof, tensor, scalar
+from theano.tensor import blas
+
 from theano.sparse.basic import (
-    _is_sparse, _is_sparse_variable, _is_dense_variable,
-    _is_sparse, _is_dense, _kmap_eq, _kmap_hash)
+    as_sparse_variable, SparseType, add_s_s, neg,
+    mul_s_s, mul_s_d,
+    CSMProperties, CSM, register_specialize,
+    _is_sparse_variable, CSC, CSR,
+    csm_data, csm_indices, csm_indptr, csm_shape,
+    _is_sparse)
 
 
 class Cast(gof.op.Op):
@@ -33,19 +41,21 @@ def local_add_s_s(node):
     """
     If two matrices are known to have the same sparsity pattern,
     optimize  the addition by only adding their data vector.
-    
+
     Very special case optimization. Activate when for add(x, y),
     y is an expression like sp_ones_like(x) * another_matrix.
-    
+
     This is useful for sparse weight updates.
 
-    Work also for add(x, neg(y)) in the same case. 
-    As of this writting sub is only implemented as x + neg(y) for sparse matrix.
+    Work also for add(x, neg(y)) in the same case.
+
+    As of this writting sub is only implemented as x + neg(y) for
+    sparse matrix.
     """
     if node.op == add_s_s:
         x, y = node.inputs
-        
-        # In case addition was transformed to subtraction        
+
+        # In case addition was transformed to subtraction
         if hasattr(y.owner, 'op') and y.owner.op == neg:
             y_ = y.owner.inputs[0]
         else:
@@ -54,38 +64,46 @@ def local_add_s_s(node):
             return False
         if hasattr(y_.owner, 'op') and y_.owner.op not in [mul_s_s, mul_s_d]:
             return False
-        
+
         def same_pattern(node):
             """Check node has same sparsity as x."""
-            # In case the sparse matrix is multiplied by a scalar (ex: learning rate)
+            # In case the sparse matrix is multiplied by a scalar (ex:
+            # learning rate)
             if hasattr(node.owner, 'op') and node.owner.op == mul_scalar:
                 node = node.owner.inputs[1]
-            
+
             # Check node creates a matrix
-            if not hasattr(node.owner, 'op') or not isinstance(node.owner.op, CSM):
-              return False
-            
+            if not hasattr(node.owner, 'op') or not isinstance(node.owner.op,
+                                                               CSM):
+                return False
+
             # Check matrix is creates from CSMProperties
-            if filter(lambda i: not hasattr(i.owner, 'op') or not isinstance(i.owner.op, CSMProperties), node.owner.inputs[1:]):
-              return False
-            
+            if filter(lambda i: not hasattr(i.owner, 'op') or
+                      not isinstance(i.owner.op, CSMProperties),
+                      node.owner.inputs[1:]):
+                return False
+
             # Verify indices, indptr and shape are the same as x
             if filter(lambda i: i.owner.inputs[0] != x, node.owner.inputs[1:]):
-              return False
-                                
+                return False
+
             return True
-        
+
         if filter(same_pattern, y_.owner.inputs):
             return [add_s_s_data(x, y)]
     return False
 register_specialize(local_add_s_s)
 
+
 class AddSSData(gof.op.Op):
-    '''Add two sparse matrices assuming they have the same sparsity pattern. '''
+    '''Add two sparse matrices assuming they have the same sparsity
+    pattern. '''
     def __eq__(self, other):
         return (type(self) == type(other))
+
     def __hash__(self):
         return hash(type(self))
+
     def make_node(self, x, y):
         x, y = map(as_sparse_variable, [x, y])
         if x.type.dtype != y.type.dtype:
@@ -94,46 +112,50 @@ class AddSSData(gof.op.Op):
             raise NotImplementedError()
         return gof.Apply(self,
                          [x, y],
-                         [SparseType(dtype = x.type.dtype,
-                                 format = x.type.format).make_variable()])
-    def perform(self, node, (x, y), (out, )): 
+                         [SparseType(dtype=x.type.dtype,
+                                 format=x.type.format).make_variable()])
+
+    def perform(self, node, (x, y), (out, )):
         assert _is_sparse(x) and _is_sparse(y)
         assert x.shape == y.shape
         out[0] = x.copy()
         out[0].data += y.data
 add_s_s_data = AddSSData()
 
+
 # register a specialization to replace MulSD -> MulSDCSX
 @gof.local_optimizer([mul_s_d])
 def local_mul_s_d(node):
     if node.op == mul_s_d:
         x, y = node.inputs
-        
+
         x_is_sparse_variable = _is_sparse_variable(x)
-        y_is_sparse_variable = _is_sparse_variable(y)
-        
+        # y_is_sparse_variable = _is_sparse_variable(y)
+
         if x_is_sparse_variable:
-          svar = x
-          dvar = y
+            svar = x
+            dvar = y
         else:
-          svar = y
-          dvar = x
-        
+            svar = y
+            dvar = x
+
         if dvar.type.ndim != 2:
             return False
         if svar.type.format == 'csc':
-          CSx = CSC
-          mul_s_d_csx = mul_s_d_csc
+            CSx = CSC
+            mul_s_d_csx = mul_s_d_csc
         elif svar.type.format == 'csr':
-          CSx = CSR
-          mul_s_d_csx = mul_s_d_csr
+            CSx = CSR
+            mul_s_d_csx = mul_s_d_csr
         else:
             raise NotImplemented()
-        
-        c_data = mul_s_d_csx(csm_data(svar), csm_indices(svar), csm_indptr(svar), dvar)
-        
-        return [CSx(c_data, csm_indices(svar), csm_indptr(svar), csm_shape(svar))]
-    
+
+        c_data = mul_s_d_csx(csm_data(svar), csm_indices(svar),
+                             csm_indptr(svar), dvar)
+
+        return [CSx(c_data, csm_indices(svar), csm_indptr(svar),
+                    csm_shape(svar))]
+
     return False
 register_specialize(local_mul_s_d)
 
@@ -141,15 +163,19 @@ register_specialize(local_mul_s_d)
 class MulSDCSC(gof.Op):
     def __eq__(self, other):
         return (type(self) == type(other))
+
     def __hash__(self):
         return hash(type(self))
+
     def make_node(self, a_data, a_indices, a_indptr, b):
         assert b.type.ndim == 2
         return gof.Apply(self, [a_data, a_indices, a_indptr, b],
                                [tensor.tensor(b.dtype, (False,))])
+
     #def perform(self, node, (a_data, a_indices, a_indptr, b), (out,)):
     #    return NotImplementedError()
-    def c_code(self, node, name, (_data, _indices, _indptr, _b,), (_zout, ), sub):
+    def c_code(self, node, name, (_data, _indices, _indptr, _b,),
+               (_zout, ), sub):
 
         if node.inputs[0].type.dtype in ('complex64', 'complex128'):
             raise NotImplementedError('Complex types are not supported for a')
@@ -209,22 +235,26 @@ class MulSDCSC(gof.Op):
             }
         }
 
-        """% dict(locals(), **sub)
+        """ % dict(locals(), **sub)
 mul_s_d_csc = MulSDCSC()
 
 
 class MulSDCSR(gof.Op):
     def __eq__(self, other):
         return (type(self) == type(other))
+
     def __hash__(self):
         return hash(type(self))
+
     def make_node(self, a_data, a_indices, a_indptr, b):
         assert b.type.ndim == 2
         return gof.Apply(self, [a_data, a_indices, a_indptr, b],
                                [tensor.tensor(b.dtype, (False,))])
+
     #def perform(self, node, (a_data, a_indices, a_indptr, b), (out,)):
     #    return NotImplemented()
-    def c_code(self, node, name, (_data, _indices, _indptr, _b,), (_zout, ), sub):
+    def c_code(self, node, name, (_data, _indices, _indptr, _b,),
+               (_zout, ), sub):
 
         if node.inputs[0].type.dtype in ('complex64', 'complex128'):
             raise NotImplementedError('Complex types are not supported for a')
@@ -284,9 +314,10 @@ class MulSDCSR(gof.Op):
             }
         }
 
-        """% dict(locals(), **sub)
+        """ % dict(locals(), **sub)
 mul_s_d_csr = MulSDCSR()
 
+
 class Poisson(gof.op.Op):
     def __eq__(self, other):
         return (type(self) == type(other))
@@ -306,6 +337,7 @@ class Poisson(gof.op.Op):
         out[0].eliminate_zeros()
 poisson = Poisson()
 
+
 class Multinomial(gof.op.Op):
     def __eq__(self, other):
         return (type(self) == type(other))
@@ -713,7 +745,6 @@ class SamplingDotCsr(gof.Op):
         return blas.blas_header_text()
 
     def c_libraries(self):
-        import pdb; pdb.set_trace()
         return blas.ldflags()
 
     def c_compile_args(self):