我们从Python开源项目中，提取了以下50个代码示例，用于说明如何使用scipy.sparse.csc_matrix()。
def load(group): gene_ids = list(getattr(group, cr_constants.H5_GENE_IDS_ATTR).read()) if hasattr(group, cr_constants.H5_GENE_NAMES_ATTR): gene_names = list(getattr(group, cr_constants.H5_GENE_NAMES_ATTR).read()) else: gene_names = gene_ids assert len(gene_ids) == len(gene_names) genes = [cr_constants.Gene(id, name, None, None, None) for id, name in itertools.izip(gene_ids, gene_names)] bcs = list(getattr(group, cr_constants.H5_BCS_ATTR).read()) matrix = GeneBCMatrix(genes, bcs) shape = getattr(group, cr_constants.H5_MATRIX_SHAPE_ATTR).read() data = getattr(group, cr_constants.H5_MATRIX_DATA_ATTR).read() indices = getattr(group, cr_constants.H5_MATRIX_INDICES_ATTR).read() indptr = getattr(group, cr_constants.H5_MATRIX_INDPTR_ATTR).read() # quick check to make sure indptr increases monotonically (to catch overflow bugs) assert np.all(np.diff(indptr)>=0) matrix.m = sp_sparse.csc_matrix((data, indices, indptr), shape=shape) return matrix
def _local_add_sparse(ltenss): """Computes the local tensors of a sum of MPArrays (except for the boundary tensors). Works only for products right now :param ltenss: Raveled local tensors :returns: Correct local tensor representation """ dim = len(ltenss[0]) nr_summands = len(ltenss) indptr = np.arange(nr_summands * dim + 1) indices = np.concatenate((np.arange(nr_summands),) * dim) data = np.concatenate([lt[None, :] for lt in ltenss]) data = np.rollaxis(data, 1).ravel() return ssp.csc_matrix((data, indices, indptr), shape=(nr_summands, dim * nr_summands))
def write_data(self, result_dict): for key, result in six.iteritems(result_dict): if ss.isspmatrix(result): if np.isnan(result.data).any(): raise ValueError("data {} have nan".format(key)) elif np.isnan(result).any(): raise ValueError("data {} have nan".format(key)) with SimpleTimer("Writing generated data {} to hdf5 file" .format(key), end_in_new_line=False): if key in self.h5f: # self.h5f[key][...] = result raise NotImplementedError("Overwriting not supported.") else: if (isinstance(result, ss.csc_matrix) or isinstance(result, ss.csr_matrix)): # sparse matrix h5sparse.Group(self.h5f).create_dataset(key, data=result) else: self.h5f.create_dataset(key, data=result) self.h5f.flush()
def solve_linear2(model:Model.fem_model): K_bar,F_bar,index=model.K_,model.F_,model.index Dvec=model.D Logger.info('Solving linear model with %d DOFs...'%model.DOF) n_nodes=model.node_count #sparse matrix solution delta_bar = sl.spsolve(sp.csc_matrix(K_bar),F_bar) #delta_bar=linalg.solve(K_bar,F_bar,sym_pos=True) delta = delta_bar #fill original displacement vector prev = 0 for idx in index: gap=idx-prev if gap>0: delta=np.insert(delta,prev,[0]*gap) prev = idx + 1 if idx==index[-1] and idx!=n_nodes-1: delta = np.insert(delta,prev, [0]*(n_nodes*6-prev)) delta += Dvec model.is_solved=True return delta
def theta_matrix(coord, adj, preload=True, train=True): print "creating adjacent theta matrix ..." if preload is True: if train is True: theta_matrix = np.load('../data/theta_matrix_train_100.npy') else: theta_matrix = np.load('../data/theta_matrix_test_100.npy') else: theta_matrix = [] for i in tqdm(range(coord.shape[0])): for j in range(coord.shape[1]): theta_row = angle(coord[i,adj[i][j].nonzero()[1],:] - coord[i,j,:]) col_indice = adj[i][j].nonzero()[1] row_indice = (np.zeros(col_indice.shape[0])).astype(int32) if j == 0: theta_matrix_tmp = sp.csc_matrix((theta_row, (row_indice, col_indice)), shape=(1,coord.shape[1])) else: theta_matrix_tmp = sp.vstack((theta_matrix_tmp, sp.csc_matrix((theta_row, (row_indice, col_indice)), shape=(1,coord.shape[1])))) theta_matrix.append(theta_matrix_tmp) theta_matrix = np.array(theta_matrix) return theta_matrix
def compute_dr_wrt(self, wrt): if wrt is not self.v: return None v = self.v.r.reshape(-1, 3) blocks = -np.einsum('ij,ik->ijk', v, v) * (self.ss**(-3./2.)).reshape((-1, 1, 1)) for i in range(3): blocks[:, i, i] += self.s_inv if True: # pylint: disable=using-constant-test data = blocks.ravel() indptr = np.arange(0, (self.v.r.size+1)*3, 3) indices = col(np.arange(0, self.v.r.size)) indices = np.hstack([indices, indices, indices]) indices = indices.reshape((-1, 3, 3)) indices = indices.transpose((0, 2, 1)).ravel() result = sp.csc_matrix((data, indices, indptr), shape=(self.v.r.size, self.v.r.size)) return result else: matvec = lambda x: np.einsum('ijk,ik->ij', blocks, x.reshape((blocks.shape[0], 3))).ravel() return sp.linalg.LinearOperator((self.v.r.size, self.v.r.size), matvec=matvec)
def compute_dr_wrt(self, wrt): if wrt is not self.v: return None cplus = self.cplus cminus = self.cminus vplus = self.f[:, cplus] vminus = self.f[:, cminus] vplus3 = row(np.hstack([col(vplus*3), col(vplus*3+1), col(vplus*3+2)])) vminus3 = row(np.hstack([col(vminus*3), col(vminus*3+1), col(vminus*3+2)])) IS = row(np.arange(0, vplus3.size)) ones = np.ones(vplus3.size) shape = (self.f.size, self.v.r.size) dr_vplus, dr_vminus = [ sp.csc_matrix((ones, np.vstack([IS, item])), shape=shape) for item in vplus3, vminus3 # FIXME change item to a DAMP ] return dr_vplus - dr_vminus
def lchol(A): import copy import numpy as np import scipy.sparse as sp from blmath.numerics.linalg import cholmod # pylint: disable=no-name-in-module # the cholmod port needs 64b data, while scipy sparse is 32b A64 = copy.copy(A) A64.indices = A.indices.astype(np.int64) A64.indptr = A.indptr.astype(np.int64) [L_ind, L_iptr, L_data, L_nonpsd, q] = cholmod.lchol_c(A64) L = sp.csc_matrix((L_data, L_ind, L_iptr), shape=A64.shape) S = sp.csc_matrix((np.ones(len(q)), (q, range(len(q)))), shape=(max(q)+1, len(q))) return L, L_nonpsd, S
def setIWMult(self, mult): m = mult / float(self.iwMult) if self.sparse: v = spsparse.lil_matrix((self.hw.shape[0], self.hw.shape[0])) v.setdiag(np.ones(v.shape[0], dtype=self.dtype)) v.setdiag([m,]*(self.nIn+1)) self.hw = v*self.hw # good for debugging above #w = self.hw.todense() #w[:self.nIn+1,:] *= m #self.hw = spsparse.csc_matrix(w) else: self.hw[:self.nIn+1,:] *= m self.iwMult = mult
def setRWScale(self, x, scale): # why does this method not work? XXX - idfah m = scale / float(self.rwScale) if self.sparse: v = spsparse.lil_matrix((self.hw.shape[0], self.hw.shape[0])) d = np.ones(v.shape[0], dtype=self.dtype) d[self.nIn+1:] = m v.setdiag(d) self.hw = v*self.hw # good for debugging above #w = self.hw.todense() #w[self.nIn+1:,:] *= m #self.hw = spsparse.csc_matrix(w) else: self.hw[self.nIn+1:,:] *= m self.scaleIW(x)
def __get_prolongation_matrix(ndofs_coarse, ndofs_fine): """Helper routine for the prolongation operator Args: ndofs_fine (int): number of DOFs on the fine grid ndofs_coarse (int): number of DOFs on the coarse grid Returns: scipy.sparse.csc_matrix: sparse prolongation matrix of size `ndofs_fine` x `ndofs_coarse` """ # This is a workaround, since I am not aware of a suitable way to do # this directly with sparse matrices. P = np.zeros((ndofs_fine, ndofs_coarse)) np.fill_diagonal(P[1::2, :], 1) np.fill_diagonal(P[0::2, :], 1.0/2.0) np.fill_diagonal(P[2::2, :], 1.0/2.0) return sp.csc_matrix(P)
def __init__(self, ndofs, A, rhs, *args, **kwargs): """Initialization routine for a problem Args: ndofs (int): number of degrees of freedom A (scipy.sparse.csc_matrix): a sparse system matrix rhs (numpy.ndarray): right-hand side vector *args: Variable length argument list **kwargs: Arbitrary keyword arguments """ # we use private attributes with getters and setters here to avoid overwriting by accident self._ndofs = None self._A = None self._rhs = None # now set the attributes using potential validation through the setters defined below self.ndofs = ndofs self.A = A self.rhs = rhs
def load(): '''Load ML-100k data Returns the review matrix as a numpy array''' import numpy as np from scipy import sparse from os import path if not path.exists('data/ml-100k/u.data'): raise IOError("Data has not been downloaded.\nTry the following:\n\n\tcd data\n\t./download.sh") # The input is in the form of a CSC sparse matrix, so it's a natural fit to # load the data, but we then convert to a more traditional array before # returning data = np.loadtxt('data/ml-100k/u.data') ij = data[:, :2] ij -= 1 # original data is in 1-based system values = data[:, 2] reviews = sparse.csc_matrix((values, ij.T)).astype(float) return reviews.toarray()
def generateCosineNeighborGraph(hin,kNeighbors=10,tf_param={'word':True, 'entity':False, 'we_weight':1}): X, newIds, entIds = GraphGenerator.getTFVectorX(hin,param=tf_param) cosX = cosine_similarity(X) #return sparse.csc_matrix(X.dot(X.transpose())),newIds n = cosX.shape[0] graph = np.zeros((n,n)) tic = time.time() for i in range(n): for j in np.argpartition(-cosX[i],kNeighbors)[:kNeighbors]: if j == i: continue #graph[i, j] += cosX[i, j] #graph[j, i] += cosX[i, j] graph[i, j] += 1 graph[j, i] += 1 toc = time.time() - tic return sparse.csc_matrix(graph), newIds
def generateCosineNeighborGraphfromX(X, kNeighbors=10): cosX = cosine_similarity(X) # return sparse.csc_matrix(X.dot(X.transpose())),newIds #print cosX.shape n = cosX.shape[0] graph = np.zeros((n, n)) tic = time.time() for i in range(n): for j in np.argpartition(-cosX[i], kNeighbors)[:kNeighbors]: if j == i: continue # graph[i, j] += cosX[i, j] # graph[j, i] += cosX[i, j] graph[i, j] += 1 graph[j, i] += 1 toc = time.time() - tic #print 'graph generation done in %f seconds.' % toc return sparse.csc_matrix(graph)
def semihin_experiment(scope, scope_name, count, X, newIds, label_num=5): experiment_path = 'data/local/split/' + scope_name + '/' with open('data/local/' + scope_name + '.dmp') as f: hin = pk.load(f) n = X.shape[0] e = X.shape[1] if not type(X) is np.ndarray: X = X.toarray() graph = np.zeros((n + e, n + e)) graph[0:n, n:n + e] = X graph[n:n + e, 0:n] = X.transpose() graph = sparse.csc_matrix(graph) newLabel = GraphGenerator.getNewLabels(hin) lp_param = {'alpha': 0.98, 'normalization_factor': 5, 'method': 'variant'} ssl = SSLClassifier(graph, newLabel, scope, lp_param, repeatTimes=50, trainNumbers=label_num, classCount=count) ssl.repeatedFixedExperimentwithNewIds(pathPrefix=experiment_path + 'lb' + str(label_num).zfill(3) + '_', newIds=newIds) return ssl.get_mean()
def saveSparseHDF5(matrix, prefix, fname): """ matrix: sparse matrix prefix: prefix of dataset fname : name of h5py file where matrix will be saved """ assert matrix.__class__==csr_matrix or matrix.__class__==csc_matrix,'Expecting csc/csr' with h5py.File(fname,mode='a') as f: for info in ['data','indices','indptr','shape']: key = '%s_%s'%(prefix,info) try: data = getattr(matrix, info) except: assert False,'Expecting attribute '+info+' in matrix' """ For empty arrays, data, indicies and indptr will be [] To deal w/ this use np.nan in its place """ if len(data)==0: f.create_dataset(key,data=np.array([np.nan])) else: f.create_dataset(key,data= data) key = prefix+'_type' val = matrix.__class__.__name__ f.attrs[key] = np.string_(val)
def _testSparse(): m1 = np.random.randn(4,12) m1[:,3:8] = 0 m1[0:2,:] = 0 csr_1 = csr_matrix(m1) csc_1 = csc_matrix(m1) m2 = np.random.randn(4,12) csr_2 = csr_matrix(m2) csc_2 = csc_matrix(m2) fname = 'tmp.h5' saveSparseHDF5(csc_1, 'c1', fname) saveSparseHDF5(csr_1, 'r1', fname) saveSparseHDF5(csc_2, 'c2', fname) saveSparseHDF5(csr_2, 'r2', fname) l_csc_1 = loadSparseHDF5('c1', fname) l_csr_1 = loadSparseHDF5('r1', fname) l_csc_2 = loadSparseHDF5('c2', fname) l_csr_2 = loadSparseHDF5('r2', fname) result = 0. for init,final in zip([csc_1,csr_1,csc_2,csr_2],[l_csc_1,l_csr_1,l_csc_2,l_csr_2]): result += (init-final).sum() +(init.toarray()-final.toarray()).sum() print 'Diff b/w saved vs loaded matrices',result os.unlink(fname)
def _predict_output_derivatives(self, x): n = x.shape[0] nt = self.nt ny = self.training_points[None][0][1].shape[1] num = self.num dy_dstates = np.empty(n * num['dof']) self.rbfc.compute_jac(n, x.flatten(), dy_dstates) dy_dstates = dy_dstates.reshape((n, num['dof'])) dstates_dytl = np.linalg.inv(self.mtx) ones = np.ones(self.nt) arange = np.arange(self.nt) dytl_dyt = csc_matrix((ones, (arange, arange)), shape=(num['dof'], self.nt)) dy_dyt = (dytl_dyt.T.dot(dstates_dytl.T).dot(dy_dstates.T)).T dy_dyt = np.einsum('ij,k->ijk', dy_dyt, np.ones(ny)) return {None: dy_dyt}
def get_ground_state(sparse_operator): """Compute lowest eigenvalue and eigenstate. Returns: eigenvalue: The lowest eigenvalue, a float. eigenstate: The lowest eigenstate in scipy.sparse csc format. """ if not is_hermitian(sparse_operator): raise ValueError('sparse_operator must be Hermitian.') values, vectors = scipy.sparse.linalg.eigsh( sparse_operator, 2, which='SA', maxiter=1e7) eigenstate = scipy.sparse.csc_matrix(vectors[:, 0]) eigenvalue = values[0] return eigenvalue, eigenstate.getH()
def epclosure(A, source, target, B=None, retpath=False, kind='ultrametric'): """ Source target "epistemic" closure. Python implementation. See `epclosuress` for parameters. """ # ensure A is CSR A = sp.csr_matrix(A) if B is None: B = A.tocsc() else: # ensure B is CSC B = sp.csc_matrix(B) _caps, _paths = cclosuress(A, source, retpaths=retpath, kind=kind) row = A.getrow(source) s_neighbors = row.indices[row.data > 0] s_reachables = set(np.where(_caps)[0]) return _epclosureone(A, source, target, B, retpath, _caps, _paths, s_neighbors, s_reachables)
def constr(self): def fun(x): x_coord, y_coord, z_coord = self._get_cordinates(x) return x_coord**2 + y_coord**2 + z_coord**2 - 1 def jac(x): x_coord, y_coord, z_coord = self._get_cordinates(x) Jx = 2 * np.diag(x_coord) Jy = 2 * np.diag(y_coord) Jz = 2 * np.diag(z_coord) return csc_matrix(np.hstack((Jx, Jy, Jz))) def hess(x, v): D = 2 * np.diag(v) return block_diag(D, D, D) return NonlinearConstraint(fun, ("less",), jac, hess)
def test_trust_region_barely_feasible(self): H = csc_matrix([[6, 2, 1, 3], [2, 5, 2, 4], [1, 2, 4, 5], [3, 4, 5, 7]]) A = csc_matrix([[1, 0, 1, 0], [0, 1, 1, 1]]) c = np.array([-2, -3, -3, 1]) b = -np.array([3, 0]) trust_radius = 2.32379000772445021283 Z, _, Y = projections(A) x, info = projected_cg(H, c, Z, Y, b, tol=0, trust_radius=trust_radius) assert_equal(info["stop_cond"], 2) assert_equal(info["hits_boundary"], True) assert_array_almost_equal(np.linalg.norm(x), trust_radius) assert_array_almost_equal(x, -Y.dot(b))
def test_hits_boundary(self): H = csc_matrix([[6, 2, 1, 3], [2, 5, 2, 4], [1, 2, 4, 5], [3, 4, 5, 7]]) A = csc_matrix([[1, 0, 1, 0], [0, 1, 1, 1]]) c = np.array([-2, -3, -3, 1]) b = -np.array([3, 0]) trust_radius = 3 Z, _, Y = projections(A) x, info = projected_cg(H, c, Z, Y, b, tol=0, trust_radius=trust_radius) assert_equal(info["stop_cond"], 2) assert_equal(info["hits_boundary"], True) assert_array_almost_equal(np.linalg.norm(x), trust_radius)
def test_negative_curvature(self): H = csc_matrix([[1, 2, 1, 3], [2, 0, 2, 4], [1, 2, 0, 2], [3, 4, 2, 0]]) A = csc_matrix([[1, 0, 1, 0], [0, 1, 0, 1]]) c = np.array([-2, -3, -3, 1]) b = -np.array([3, 0]) Z, _, Y = projections(A) trust_radius = 1000 x, info = projected_cg(H, c, Z, Y, b, tol=0, trust_radius=trust_radius) assert_equal(info["stop_cond"], 3) assert_equal(info["hits_boundary"], True) assert_array_almost_equal(np.linalg.norm(x), trust_radius) # The box contraints are inactive at the solution but # are active during the iterations.
def test_inactive_box_constraints(self): H = csc_matrix([[6, 2, 1, 3], [2, 5, 2, 4], [1, 2, 4, 5], [3, 4, 5, 7]]) A = csc_matrix([[1, 0, 1, 0], [0, 1, 1, 1]]) c = np.array([-2, -3, -3, 1]) b = -np.array([3, 0]) Z, _, Y = projections(A) x, info = projected_cg(H, c, Z, Y, b, tol=0, lb=[0.5, -np.inf, -np.inf, -np.inf], return_all=True) x_kkt, _ = eqp_kktfact(H, c, A, b) assert_equal(info["stop_cond"], 1) assert_equal(info["hits_boundary"], False) assert_array_almost_equal(x, x_kkt) # The box contraints active and the termination is # by maximum iterations (infeasible iteraction).
def test_active_box_constraints_hits_boundaries(self): H = csc_matrix([[6, 2, 1, 3], [2, 5, 2, 4], [1, 2, 4, 5], [3, 4, 5, 7]]) A = csc_matrix([[1, 0, 1, 0], [0, 1, 1, 1]]) c = np.array([-2, -3, -3, 1]) b = -np.array([3, 0]) trust_radius = 3 Z, _, Y = projections(A) x, info = projected_cg(H, c, Z, Y, b, tol=0, ub=[np.inf, np.inf, 1.6, np.inf], trust_radius=trust_radius, return_all=True) assert_equal(info["stop_cond"], 2) assert_equal(info["hits_boundary"], True) assert_array_almost_equal(x[2], 1.6) # The box contraints are active and the termination is # because it hits boundary (infeasible iteraction).
def test_active_box_constraints_hits_boundaries_infeasible_iter(self): H = csc_matrix([[6, 2, 1, 3], [2, 5, 2, 4], [1, 2, 4, 5], [3, 4, 5, 7]]) A = csc_matrix([[1, 0, 1, 0], [0, 1, 1, 1]]) c = np.array([-2, -3, -3, 1]) b = -np.array([3, 0]) trust_radius = 4 Z, _, Y = projections(A) x, info = projected_cg(H, c, Z, Y, b, tol=0, ub=[np.inf, 0.1, np.inf, np.inf], trust_radius=trust_radius, return_all=True) assert_equal(info["stop_cond"], 2) assert_equal(info["hits_boundary"], True) assert_array_almost_equal(x[1], 0.1) # The box contraints are active and the termination is # because it hits boundary (no infeasible iteraction).
def test_active_box_constraints_negative_curvature(self): H = csc_matrix([[1, 2, 1, 3], [2, 0, 2, 4], [1, 2, 0, 2], [3, 4, 2, 0]]) A = csc_matrix([[1, 0, 1, 0], [0, 1, 0, 1]]) c = np.array([-2, -3, -3, 1]) b = -np.array([3, 0]) Z, _, Y = projections(A) trust_radius = 1000 x, info = projected_cg(H, c, Z, Y, b, tol=0, ub=[np.inf, np.inf, 100, np.inf], trust_radius=trust_radius) assert_equal(info["stop_cond"], 3) assert_equal(info["hits_boundary"], True) assert_array_almost_equal(x[2], 100)
def test_nullspace_and_least_squares_sparse(self): A_dense = np.array([[1, 2, 3, 4, 0, 5, 0, 7], [0, 8, 7, 0, 1, 5, 9, 0], [1, 0, 0, 0, 0, 1, 2, 3]]) At_dense = A_dense.T A = csc_matrix(A_dense) test_points = ([1, 2, 3, 4, 5, 6, 7, 8], [1, 10, 3, 0, 1, 6, 7, 8], [1.12, 10, 0, 0, 100000, 6, 0.7, 8]) for method in available_sparse_methods: Z, LS, _ = projections(A, method) for z in test_points: # Test if x is in the null_space x = Z.matvec(z) assert_array_almost_equal(A.dot(x), 0) # Test orthogonality assert_array_almost_equal(orthogonality(A, x), 0) # Test if x is the least square solution x = LS.matvec(z) x2 = scipy.linalg.lstsq(At_dense, z)[0] assert_array_almost_equal(x, x2)
def test_iterative_refinements_sparse(self): A_dense = np.array([[1, 2, 3, 4, 0, 5, 0, 7], [0, 8, 7, 0, 1, 5, 9, 0], [1, 0, 0, 0, 0, 1, 2, 3]]) A = csc_matrix(A_dense) test_points = ([1, 2, 3, 4, 5, 6, 7, 8], [1, 10, 3, 0, 1, 6, 7, 8], [1.12, 10, 0, 0, 100000, 6, 0.7, 8], [1, 0, 0, 0, 0, 1, 2, 3+1e-10]) for method in available_sparse_methods: Z, LS, _ = projections(A, method, orth_tol=1e-18, max_refin=100) for z in test_points: # Test if x is in the null_space x = Z.matvec(z) assert_array_almost_equal(A.dot(x), 0, decimal=14) # Test orthogonality assert_array_almost_equal(orthogonality(A, x), 0, decimal=16)
def test_compare_dense_and_sparse2(self): D1 = np.diag([-1.7, 1, 0.5]) D2 = np.diag([1, -0.6, -0.3]) D3 = np.diag([-0.3, -1.5, 2]) A = np.hstack([D1, D2, D3]) A_sparse = csc_matrix(A) np.random.seed(0) Z, LS, Y = projections(A) Z_sparse, LS_sparse, Y_sparse = projections(A_sparse) for k in range(1): z = np.random.normal(size=(9,)) assert_array_almost_equal(Z.dot(z), Z_sparse.dot(z)) assert_array_almost_equal(LS.dot(z), LS_sparse.dot(z)) x = np.random.normal(size=(3,)) assert_array_almost_equal(Y.dot(x), Y_sparse.dot(x))
def test_sparse_matrix(self): A = np.array([[1, 2, 3, 4, 0, 5, 0, 7], [0, 8, 7, 0, 1, 5, 9, 0], [1, 0, 0, 0, 0, 1, 2, 3]]) A = csc_matrix(A) test_vectors = ([-1.98931144, -1.56363389, -0.84115584, 2.2864762, 5.599141, 0.09286976, 1.37040802, -0.28145812], [697.92794044, -4091.65114008, -3327.42316335, 836.86906951, 99434.98929065, -1285.37653682, -4109.21503806, 2935.29289083]) test_expected_orth = (0, 0) for i in range(len(test_vectors)): x = test_vectors[i] orth = test_expected_orth[i] assert_array_almost_equal(orthogonality(A, x), orth)
def test_local_csm_properties_csm(): data = tensor.vector() indices, indptr, shape = (tensor.ivector(), tensor.ivector(), tensor.ivector()) mode = theano.compile.mode.get_default_mode() mode = mode.including("specialize", "local_csm_properties_csm") for CS, cast in [(sparse.CSC, sp.csc_matrix), (sparse.CSR, sp.csr_matrix)]: f = theano.function([data, indices, indptr, shape], sparse.csm_properties( CS(data, indices, indptr, shape)), mode=mode) assert not any( isinstance(node.op, (sparse.CSM, sparse.CSMProperties)) for node in f.maker.fgraph.toposort()) v = cast(random_lil((10, 40), config.floatX, 3)) f(v.data, v.indices, v.indptr, v.shape)
def test_local_csm_grad_c(): raise SkipTest("Opt disabled as it don't support unsorted indices") if not theano.config.cxx: raise SkipTest("G++ not available, so we need to skip this test.") data = tensor.vector() indices, indptr, shape = (tensor.ivector(), tensor.ivector(), tensor.ivector()) mode = theano.compile.mode.get_default_mode() if theano.config.mode == 'FAST_COMPILE': mode = theano.compile.Mode(linker='c|py', optimizer='fast_compile') mode = mode.including("specialize", "local_csm_grad_c") for CS, cast in [(sparse.CSC, sp.csc_matrix), (sparse.CSR, sp.csr_matrix)]: cost = tensor.sum(sparse.DenseFromSparse()(CS(data, indices, indptr, shape))) f = theano.function( [data, indices, indptr, shape], tensor.grad(cost, data), mode=mode) assert not any(isinstance(node.op, sparse.CSMGrad) for node in f.maker.fgraph.toposort()) v = cast(random_lil((10, 40), config.floatX, 3)) f(v.data, v.indices, v.indptr, v.shape)
def test_equality_case(self): """ Test assuring normal behaviour when values in the matrices are equal """ scipy_ver = [int(n) for n in scipy.__version__.split('.')[:2]] if (bool(scipy_ver < [0, 13])): raise SkipTest("comparison operators need newer release of scipy") x = sparse.csc_matrix() y = theano.tensor.matrix() m1 = sp.csc_matrix((2, 2), dtype=theano.config.floatX) m2 = numpy.asarray([[0, 0], [0, 0]], dtype=theano.config.floatX) for func in self.testsDic: op = func(y, x) f = theano.function([y, x], op) self.assertTrue(numpy.array_equal(f(m2, m1), self.testsDic[func](m2, m1)))
def test_csm_properties_grad(self): sp_types = {'csc': sp.csc_matrix, 'csr': sp.csr_matrix} for format in ['csc', 'csr']: for dtype in ['float32', 'float64']: spmat = sp_types[format](random_lil((4, 3), dtype, 3)) verify_grad_sparse(lambda *x: CSMProperties()(*x)[0], [spmat], structured=True) verify_grad_sparse(lambda *x: CSMProperties()(*x)[1], [spmat], structured=True) verify_grad_sparse(lambda *x: CSMProperties()(*x)[2], [spmat], structured=True) verify_grad_sparse(lambda *x: CSMProperties()(*x)[2], [spmat], structured=True)
def test_csm_unsorted(self): """ Test support for gradients of unsorted inputs. """ sp_types = {'csc': sp.csc_matrix, 'csr': sp.csr_matrix} for format in ['csr', 'csc', ]: for dtype in ['float32', 'float64']: x = tensor.tensor(dtype=dtype, broadcastable=(False,)) y = tensor.ivector() z = tensor.ivector() s = tensor.ivector() # Sparse advanced indexing produces unsorted sparse matrices a = sparse_random_inputs(format, (4, 3), out_dtype=dtype, unsorted_indices=True)[1][0] # Make sure it's unsorted assert not a.has_sorted_indices def my_op(x): y = tensor.constant(a.indices) z = tensor.constant(a.indptr) s = tensor.constant(a.shape) return tensor.sum( dense_from_sparse(CSM(format)(x, y, z, s) * a)) verify_grad_sparse(my_op, [a.data])
def test_csm(self): sp_types = {'csc': sp.csc_matrix, 'csr': sp.csr_matrix} for format in ['csc', 'csr']: for dtype in ['float32', 'float64']: x = tensor.tensor(dtype=dtype, broadcastable=(False,)) y = tensor.ivector() z = tensor.ivector() s = tensor.ivector() f = theano.function([x, y, z, s], CSM(format)(x, y, z, s)) spmat = sp_types[format](random_lil((4, 3), dtype, 3)) res = f(spmat.data, spmat.indices, spmat.indptr, numpy.asarray(spmat.shape, 'int32')) assert numpy.all(res.data == spmat.data) assert numpy.all(res.indices == spmat.indices) assert numpy.all(res.indptr == spmat.indptr) assert numpy.all(res.shape == spmat.shape)
def test_dot_sparse_sparse(self): # test dot for 2 input sparse matrix sparse_dtype = 'float64' sp_mat = {'csc': sp.csc_matrix, 'csr': sp.csr_matrix, 'bsr': sp.csr_matrix} for sparse_format_a in ['csc', 'csr', 'bsr']: for sparse_format_b in ['csc', 'csr', 'bsr']: a = SparseType(sparse_format_a, dtype=sparse_dtype)() b = SparseType(sparse_format_b, dtype=sparse_dtype)() d = theano.dot(a, b) f = theano.function([a, b], theano.Out(d, borrow=True)) topo = f.maker.fgraph.toposort() for M, N, K, nnz in [(4, 3, 2, 3), (40, 30, 20, 3), (40, 30, 20, 30), (400, 3000, 200, 6000), ]: a_val = sp_mat[sparse_format_a]( random_lil((M, N), sparse_dtype, nnz)) b_val = sp_mat[sparse_format_b]( random_lil((N, K), sparse_dtype, nnz)) f(a_val, b_val)
def setUp(self): super(DotTests, self).setUp() x_size = (10, 100) y_size = (100, 1000) utt.seed_rng() self.x_csr = scipy.sparse.csr_matrix( numpy.random.binomial(1, 0.5, x_size), dtype=theano.config.floatX) self.x_csc = scipy.sparse.csc_matrix( numpy.random.binomial(1, 0.5, x_size), dtype=theano.config.floatX) self.y = numpy.asarray(numpy.random.uniform(-1, 1, y_size), dtype=theano.config.floatX) self.y_csr = scipy.sparse.csr_matrix( numpy.random.binomial(1, 0.5, y_size), dtype=theano.config.floatX) self.y_csc = scipy.sparse.csc_matrix( numpy.random.binomial(1, 0.5, y_size), dtype=theano.config.floatX) self.v_10 = numpy.asarray(numpy.random.uniform(-1, 1, 10), dtype=theano.config.floatX) self.v_100 = numpy.asarray(numpy.random.uniform(-1, 1, 100), dtype=theano.config.floatX)
def test_csc_dense(self): x = theano.sparse.csc_matrix('x') y = theano.tensor.matrix('y') v = theano.tensor.vector('v') for (x, y, x_v, y_v) in [(x, y, self.x_csc, self.y), (x, v, self.x_csc, self.v_100), (v, x, self.v_10, self.x_csc)]: f_a = theano.function([x, y], theano.sparse.dot(x, y)) f_b = lambda x, y: x * y utt.assert_allclose(f_a(x_v, y_v), f_b(x_v, y_v)) # Test infer_shape self._compile_and_check([x, y], [theano.sparse.dot(x, y)], [x_v, y_v], (Dot, Usmm, UsmmCscDense))
def test_size(): """ Ensure the `size` attribute of sparse matrices behaves as in numpy. """ for sparse_type in ('csc_matrix', 'csr_matrix'): x = getattr(theano.sparse, sparse_type)() y = getattr(scipy.sparse, sparse_type)((5, 7)).astype(config.floatX) get_size = theano.function([x], x.size) def check(): assert y.size == get_size(y) # We verify that the size is correctly updated as we store more data # into the sparse matrix (including zeros). check() y[0, 0] = 1 check() y[0, 1] = 0 check()
def test_infer_shape(self): mat = (numpy.arange(12) + 1).reshape((4, 3)) mat[0, 1] = mat[1, 0] = mat[2, 2] = 0 x_csc = theano.sparse.csc_matrix(dtype=theano.config.floatX) mat_csc = sp.csc_matrix(mat, dtype=theano.config.floatX) self._compile_and_check([x_csc], [Remove0()(x_csc)], [mat_csc], self.op_class) x_csr = theano.sparse.csr_matrix(dtype=theano.config.floatX) mat_csr = sp.csr_matrix(mat, dtype=theano.config.floatX) self._compile_and_check([x_csr], [Remove0()(x_csr)], [mat_csr], self.op_class)
def test_GetItemList(self): a, A = sparse_random_inputs('csr', (4, 5)) b, B = sparse_random_inputs('csc', (4, 5)) y = a[0][[0, 1, 2, 3, 1]] z = b[0][[0, 1, 2, 3, 1]] fa = theano.function([a[0]], y) fb = theano.function([b[0]], z) t_geta = fa(A[0]).todense() t_getb = fb(B[0]).todense() s_geta = scipy.sparse.csr_matrix(A[0])[[0, 1, 2, 3, 1]].todense() s_getb = scipy.sparse.csc_matrix(B[0])[[0, 1, 2, 3, 1]].todense() utt.assert_allclose(t_geta, s_geta) utt.assert_allclose(t_getb, s_getb)
def test_mul_s_v(self): sp_types = {'csc': sp.csc_matrix, 'csr': sp.csr_matrix} for format in ['csr', 'csc']: for dtype in ['float32', 'float64']: x = theano.sparse.SparseType(format, dtype=dtype)() y = tensor.vector(dtype=dtype) f = theano.function([x, y], mul_s_v(x, y)) spmat = sp_types[format](random_lil((4, 3), dtype, 3)) mat = numpy.asarray(numpy.random.rand(3), dtype=dtype) out = f(spmat, mat) utt.assert_allclose(spmat.toarray() * mat, out.toarray())
def test_structured_add_s_v(self): sp_types = {'csc': sp.csc_matrix, 'csr': sp.csr_matrix} for format in ['csr', 'csc']: for dtype in ['float32', 'float64']: x = theano.sparse.SparseType(format, dtype=dtype)() y = tensor.vector(dtype=dtype) f = theano.function([x, y], structured_add_s_v(x, y)) spmat = sp_types[format](random_lil((4, 3), dtype, 3)) spones = spmat.copy() spones.data = numpy.ones_like(spones.data) mat = numpy.asarray(numpy.random.rand(3), dtype=dtype) out = f(spmat, mat) utt.assert_allclose(as_ndarray(spones.multiply(spmat + mat)), out.toarray())
def perform(self, node, inputs, outputs): # for efficiency, if remap does nothing, then do not apply it (data, indices, indptr, shape) = inputs (out,) = outputs if len(shape) != 2: raise ValueError('Shape should be an array of length 2') if data.shape != indices.shape: errmsg = ('Data (shape ' + repr(data.shape) + ' must have the same number of elements ' + 'as indices (shape' + repr(indices.shape) + ')') raise ValueError(errmsg) if self.format == 'csc': out[0] = scipy.sparse.csc_matrix((data, indices.copy(), indptr.copy()), numpy.asarray(shape), copy=False) else: assert self.format == 'csr' out[0] = scipy.sparse.csr_matrix((data, indices.copy(), indptr.copy()), shape.copy(), copy=False)
def perform(self, node, inp, outputs): (out,) = outputs x = inp[0] ind1 = inp[1] ind2 = inp[2] gz = inp[3] if x.format in ["csr"]: y = scipy.sparse.csr_matrix((x.shape[0], x.shape[1])) else: y = scipy.sparse.csc_matrix((x.shape[0], x.shape[1])) z = 0 for z in range(0, len(ind1)): y[(ind1[z], ind2[z])] = gz[z] out[0] = y
