Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1"""Compressed Block Sparse Row matrix format""" 

2 

3__docformat__ = "restructuredtext en" 

4 

5__all__ = ['bsr_matrix', 'isspmatrix_bsr'] 

6 

7from warnings import warn 

8 

9import numpy as np 

10 

11from .data import _data_matrix, _minmax_mixin 

12from .compressed import _cs_matrix 

13from .base import isspmatrix, _formats, spmatrix 

14from .sputils import (isshape, getdtype, to_native, upcast, get_index_dtype, 

15 check_shape) 

16from . import _sparsetools 

17from ._sparsetools import (bsr_matvec, bsr_matvecs, csr_matmat_maxnnz, 

18 bsr_matmat, bsr_transpose, bsr_sort_indices, 

19 bsr_tocsr) 

20 

21 

22class bsr_matrix(_cs_matrix, _minmax_mixin): 

23 """Block Sparse Row matrix 

24 

25 This can be instantiated in several ways: 

26 bsr_matrix(D, [blocksize=(R,C)]) 

27 where D is a dense matrix or 2-D ndarray. 

28 

29 bsr_matrix(S, [blocksize=(R,C)]) 

30 with another sparse matrix S (equivalent to S.tobsr()) 

31 

32 bsr_matrix((M, N), [blocksize=(R,C), dtype]) 

33 to construct an empty matrix with shape (M, N) 

34 dtype is optional, defaulting to dtype='d'. 

35 

36 bsr_matrix((data, ij), [blocksize=(R,C), shape=(M, N)]) 

37 where ``data`` and ``ij`` satisfy ``a[ij[0, k], ij[1, k]] = data[k]`` 

38 

39 bsr_matrix((data, indices, indptr), [shape=(M, N)]) 

40 is the standard BSR representation where the block column 

41 indices for row i are stored in ``indices[indptr[i]:indptr[i+1]]`` 

42 and their corresponding block values are stored in 

43 ``data[ indptr[i]: indptr[i+1] ]``. If the shape parameter is not 

44 supplied, the matrix dimensions are inferred from the index arrays. 

45 

46 Attributes 

47 ---------- 

48 dtype : dtype 

49 Data type of the matrix 

50 shape : 2-tuple 

51 Shape of the matrix 

52 ndim : int 

53 Number of dimensions (this is always 2) 

54 nnz 

55 Number of stored values, including explicit zeros 

56 data 

57 Data array of the matrix 

58 indices 

59 BSR format index array 

60 indptr 

61 BSR format index pointer array 

62 blocksize 

63 Block size of the matrix 

64 has_sorted_indices 

65 Whether indices are sorted 

66 

67 Notes 

68 ----- 

69 Sparse matrices can be used in arithmetic operations: they support 

70 addition, subtraction, multiplication, division, and matrix power. 

71 

72 **Summary of BSR format** 

73 

74 The Block Compressed Row (BSR) format is very similar to the Compressed 

75 Sparse Row (CSR) format. BSR is appropriate for sparse matrices with dense 

76 sub matrices like the last example below. Block matrices often arise in 

77 vector-valued finite element discretizations. In such cases, BSR is 

78 considerably more efficient than CSR and CSC for many sparse arithmetic 

79 operations. 

80 

81 **Blocksize** 

82 

83 The blocksize (R,C) must evenly divide the shape of the matrix (M,N). 

84 That is, R and C must satisfy the relationship ``M % R = 0`` and 

85 ``N % C = 0``. 

86 

87 If no blocksize is specified, a simple heuristic is applied to determine 

88 an appropriate blocksize. 

89 

90 Examples 

91 -------- 

92 >>> from scipy.sparse import bsr_matrix 

93 >>> bsr_matrix((3, 4), dtype=np.int8).toarray() 

94 array([[0, 0, 0, 0], 

95 [0, 0, 0, 0], 

96 [0, 0, 0, 0]], dtype=int8) 

97 

98 >>> row = np.array([0, 0, 1, 2, 2, 2]) 

99 >>> col = np.array([0, 2, 2, 0, 1, 2]) 

100 >>> data = np.array([1, 2, 3 ,4, 5, 6]) 

101 >>> bsr_matrix((data, (row, col)), shape=(3, 3)).toarray() 

102 array([[1, 0, 2], 

103 [0, 0, 3], 

104 [4, 5, 6]]) 

105 

106 >>> indptr = np.array([0, 2, 3, 6]) 

107 >>> indices = np.array([0, 2, 2, 0, 1, 2]) 

108 >>> data = np.array([1, 2, 3, 4, 5, 6]).repeat(4).reshape(6, 2, 2) 

109 >>> bsr_matrix((data,indices,indptr), shape=(6, 6)).toarray() 

110 array([[1, 1, 0, 0, 2, 2], 

111 [1, 1, 0, 0, 2, 2], 

112 [0, 0, 0, 0, 3, 3], 

113 [0, 0, 0, 0, 3, 3], 

114 [4, 4, 5, 5, 6, 6], 

115 [4, 4, 5, 5, 6, 6]]) 

116 

117 """ 

118 format = 'bsr' 

119 

120 def __init__(self, arg1, shape=None, dtype=None, copy=False, blocksize=None): 

121 _data_matrix.__init__(self) 

122 

123 if isspmatrix(arg1): 

124 if isspmatrix_bsr(arg1) and copy: 

125 arg1 = arg1.copy() 

126 else: 

127 arg1 = arg1.tobsr(blocksize=blocksize) 

128 self._set_self(arg1) 

129 

130 elif isinstance(arg1,tuple): 

131 if isshape(arg1): 

132 # it's a tuple of matrix dimensions (M,N) 

133 self._shape = check_shape(arg1) 

134 M,N = self.shape 

135 # process blocksize 

136 if blocksize is None: 

137 blocksize = (1,1) 

138 else: 

139 if not isshape(blocksize): 

140 raise ValueError('invalid blocksize=%s' % blocksize) 

141 blocksize = tuple(blocksize) 

142 self.data = np.zeros((0,) + blocksize, getdtype(dtype, default=float)) 

143 

144 R,C = blocksize 

145 if (M % R) != 0 or (N % C) != 0: 

146 raise ValueError('shape must be multiple of blocksize') 

147 

148 # Select index dtype large enough to pass array and 

149 # scalar parameters to sparsetools 

150 idx_dtype = get_index_dtype(maxval=max(M//R, N//C, R, C)) 

151 self.indices = np.zeros(0, dtype=idx_dtype) 

152 self.indptr = np.zeros(M//R + 1, dtype=idx_dtype) 

153 

154 elif len(arg1) == 2: 

155 # (data,(row,col)) format 

156 from .coo import coo_matrix 

157 self._set_self(coo_matrix(arg1, dtype=dtype).tobsr(blocksize=blocksize)) 

158 

159 elif len(arg1) == 3: 

160 # (data,indices,indptr) format 

161 (data, indices, indptr) = arg1 

162 

163 # Select index dtype large enough to pass array and 

164 # scalar parameters to sparsetools 

165 maxval = 1 

166 if shape is not None: 

167 maxval = max(shape) 

168 if blocksize is not None: 

169 maxval = max(maxval, max(blocksize)) 

170 idx_dtype = get_index_dtype((indices, indptr), maxval=maxval, check_contents=True) 

171 

172 self.indices = np.array(indices, copy=copy, dtype=idx_dtype) 

173 self.indptr = np.array(indptr, copy=copy, dtype=idx_dtype) 

174 self.data = np.array(data, copy=copy, dtype=getdtype(dtype, data)) 

175 else: 

176 raise ValueError('unrecognized bsr_matrix constructor usage') 

177 else: 

178 # must be dense 

179 try: 

180 arg1 = np.asarray(arg1) 

181 except Exception: 

182 raise ValueError("unrecognized form for" 

183 " %s_matrix constructor" % self.format) 

184 from .coo import coo_matrix 

185 arg1 = coo_matrix(arg1, dtype=dtype).tobsr(blocksize=blocksize) 

186 self._set_self(arg1) 

187 

188 if shape is not None: 

189 self._shape = check_shape(shape) 

190 else: 

191 if self.shape is None: 

192 # shape not already set, try to infer dimensions 

193 try: 

194 M = len(self.indptr) - 1 

195 N = self.indices.max() + 1 

196 except Exception: 

197 raise ValueError('unable to infer matrix dimensions') 

198 else: 

199 R,C = self.blocksize 

200 self._shape = check_shape((M*R,N*C)) 

201 

202 if self.shape is None: 

203 if shape is None: 

204 # TODO infer shape here 

205 raise ValueError('need to infer shape') 

206 else: 

207 self._shape = check_shape(shape) 

208 

209 if dtype is not None: 

210 self.data = self.data.astype(dtype, copy=False) 

211 

212 self.check_format(full_check=False) 

213 

214 def check_format(self, full_check=True): 

215 """check whether the matrix format is valid 

216 

217 *Parameters*: 

218 full_check: 

219 True - rigorous check, O(N) operations : default 

220 False - basic check, O(1) operations 

221 

222 """ 

223 M,N = self.shape 

224 R,C = self.blocksize 

225 

226 # index arrays should have integer data types 

227 if self.indptr.dtype.kind != 'i': 

228 warn("indptr array has non-integer dtype (%s)" 

229 % self.indptr.dtype.name) 

230 if self.indices.dtype.kind != 'i': 

231 warn("indices array has non-integer dtype (%s)" 

232 % self.indices.dtype.name) 

233 

234 idx_dtype = get_index_dtype((self.indices, self.indptr)) 

235 self.indptr = np.asarray(self.indptr, dtype=idx_dtype) 

236 self.indices = np.asarray(self.indices, dtype=idx_dtype) 

237 self.data = to_native(self.data) 

238 

239 # check array shapes 

240 if self.indices.ndim != 1 or self.indptr.ndim != 1: 

241 raise ValueError("indices, and indptr should be 1-D") 

242 if self.data.ndim != 3: 

243 raise ValueError("data should be 3-D") 

244 

245 # check index pointer 

246 if (len(self.indptr) != M//R + 1): 

247 raise ValueError("index pointer size (%d) should be (%d)" % 

248 (len(self.indptr), M//R + 1)) 

249 if (self.indptr[0] != 0): 

250 raise ValueError("index pointer should start with 0") 

251 

252 # check index and data arrays 

253 if (len(self.indices) != len(self.data)): 

254 raise ValueError("indices and data should have the same size") 

255 if (self.indptr[-1] > len(self.indices)): 

256 raise ValueError("Last value of index pointer should be less than " 

257 "the size of index and data arrays") 

258 

259 self.prune() 

260 

261 if full_check: 

262 # check format validity (more expensive) 

263 if self.nnz > 0: 

264 if self.indices.max() >= N//C: 

265 raise ValueError("column index values must be < %d (now max %d)" % (N//C, self.indices.max())) 

266 if self.indices.min() < 0: 

267 raise ValueError("column index values must be >= 0") 

268 if np.diff(self.indptr).min() < 0: 

269 raise ValueError("index pointer values must form a " 

270 "non-decreasing sequence") 

271 

272 # if not self.has_sorted_indices(): 

273 # warn('Indices were not in sorted order. Sorting indices.') 

274 # self.sort_indices(check_first=False) 

275 

276 def _get_blocksize(self): 

277 return self.data.shape[1:] 

278 blocksize = property(fget=_get_blocksize) 

279 

280 def getnnz(self, axis=None): 

281 if axis is not None: 

282 raise NotImplementedError("getnnz over an axis is not implemented " 

283 "for BSR format") 

284 R,C = self.blocksize 

285 return int(self.indptr[-1] * R * C) 

286 

287 getnnz.__doc__ = spmatrix.getnnz.__doc__ 

288 

289 def __repr__(self): 

290 format = _formats[self.getformat()][1] 

291 return ("<%dx%d sparse matrix of type '%s'\n" 

292 "\twith %d stored elements (blocksize = %dx%d) in %s format>" % 

293 (self.shape + (self.dtype.type, self.nnz) + self.blocksize + 

294 (format,))) 

295 

296 def diagonal(self, k=0): 

297 rows, cols = self.shape 

298 if k <= -rows or k >= cols: 

299 return np.empty(0, dtype=self.data.dtype) 

300 R, C = self.blocksize 

301 y = np.zeros(min(rows + min(k, 0), cols - max(k, 0)), 

302 dtype=upcast(self.dtype)) 

303 _sparsetools.bsr_diagonal(k, rows // R, cols // C, R, C, 

304 self.indptr, self.indices, 

305 np.ravel(self.data), y) 

306 return y 

307 

308 diagonal.__doc__ = spmatrix.diagonal.__doc__ 

309 

310 ########################## 

311 # NotImplemented methods # 

312 ########################## 

313 

314 def __getitem__(self,key): 

315 raise NotImplementedError 

316 

317 def __setitem__(self,key,val): 

318 raise NotImplementedError 

319 

320 ###################### 

321 # Arithmetic methods # 

322 ###################### 

323 

324 @np.deprecate(message="BSR matvec is deprecated in SciPy 0.19.0. " 

325 "Use * operator instead.") 

326 def matvec(self, other): 

327 """Multiply matrix by vector.""" 

328 return self * other 

329 

330 @np.deprecate(message="BSR matmat is deprecated in SciPy 0.19.0. " 

331 "Use * operator instead.") 

332 def matmat(self, other): 

333 """Multiply this sparse matrix by other matrix.""" 

334 return self * other 

335 

336 def _add_dense(self, other): 

337 return self.tocoo(copy=False)._add_dense(other) 

338 

339 def _mul_vector(self, other): 

340 M,N = self.shape 

341 R,C = self.blocksize 

342 

343 result = np.zeros(self.shape[0], dtype=upcast(self.dtype, other.dtype)) 

344 

345 bsr_matvec(M//R, N//C, R, C, 

346 self.indptr, self.indices, self.data.ravel(), 

347 other, result) 

348 

349 return result 

350 

351 def _mul_multivector(self,other): 

352 R,C = self.blocksize 

353 M,N = self.shape 

354 n_vecs = other.shape[1] # number of column vectors 

355 

356 result = np.zeros((M,n_vecs), dtype=upcast(self.dtype,other.dtype)) 

357 

358 bsr_matvecs(M//R, N//C, n_vecs, R, C, 

359 self.indptr, self.indices, self.data.ravel(), 

360 other.ravel(), result.ravel()) 

361 

362 return result 

363 

364 def _mul_sparse_matrix(self, other): 

365 M, K1 = self.shape 

366 K2, N = other.shape 

367 

368 R,n = self.blocksize 

369 

370 # convert to this format 

371 if isspmatrix_bsr(other): 

372 C = other.blocksize[1] 

373 else: 

374 C = 1 

375 

376 from .csr import isspmatrix_csr 

377 

378 if isspmatrix_csr(other) and n == 1: 

379 other = other.tobsr(blocksize=(n,C), copy=False) # lightweight conversion 

380 else: 

381 other = other.tobsr(blocksize=(n,C)) 

382 

383 idx_dtype = get_index_dtype((self.indptr, self.indices, 

384 other.indptr, other.indices)) 

385 

386 bnnz = csr_matmat_maxnnz(M//R, N//C, 

387 self.indptr.astype(idx_dtype), 

388 self.indices.astype(idx_dtype), 

389 other.indptr.astype(idx_dtype), 

390 other.indices.astype(idx_dtype)) 

391 

392 idx_dtype = get_index_dtype((self.indptr, self.indices, 

393 other.indptr, other.indices), 

394 maxval=bnnz) 

395 indptr = np.empty(self.indptr.shape, dtype=idx_dtype) 

396 indices = np.empty(bnnz, dtype=idx_dtype) 

397 data = np.empty(R*C*bnnz, dtype=upcast(self.dtype,other.dtype)) 

398 

399 bsr_matmat(bnnz, M//R, N//C, R, C, n, 

400 self.indptr.astype(idx_dtype), 

401 self.indices.astype(idx_dtype), 

402 np.ravel(self.data), 

403 other.indptr.astype(idx_dtype), 

404 other.indices.astype(idx_dtype), 

405 np.ravel(other.data), 

406 indptr, 

407 indices, 

408 data) 

409 

410 data = data.reshape(-1,R,C) 

411 

412 # TODO eliminate zeros 

413 

414 return bsr_matrix((data,indices,indptr),shape=(M,N),blocksize=(R,C)) 

415 

416 ###################### 

417 # Conversion methods # 

418 ###################### 

419 

420 def tobsr(self, blocksize=None, copy=False): 

421 """Convert this matrix into Block Sparse Row Format. 

422 

423 With copy=False, the data/indices may be shared between this 

424 matrix and the resultant bsr_matrix. 

425 

426 If blocksize=(R, C) is provided, it will be used for determining 

427 block size of the bsr_matrix. 

428 """ 

429 if blocksize not in [None, self.blocksize]: 

430 return self.tocsr().tobsr(blocksize=blocksize) 

431 if copy: 

432 return self.copy() 

433 else: 

434 return self 

435 

436 def tocsr(self, copy=False): 

437 M, N = self.shape 

438 R, C = self.blocksize 

439 nnz = self.nnz 

440 idx_dtype = get_index_dtype((self.indptr, self.indices), 

441 maxval=max(nnz, N)) 

442 indptr = np.empty(M + 1, dtype=idx_dtype) 

443 indices = np.empty(nnz, dtype=idx_dtype) 

444 data = np.empty(nnz, dtype=upcast(self.dtype)) 

445 

446 bsr_tocsr(M // R, # n_brow 

447 N // C, # n_bcol 

448 R, C, 

449 self.indptr.astype(idx_dtype, copy=False), 

450 self.indices.astype(idx_dtype, copy=False), 

451 self.data, 

452 indptr, 

453 indices, 

454 data) 

455 from .csr import csr_matrix 

456 return csr_matrix((data, indices, indptr), shape=self.shape) 

457 

458 tocsr.__doc__ = spmatrix.tocsr.__doc__ 

459 

460 def tocsc(self, copy=False): 

461 return self.tocsr(copy=False).tocsc(copy=copy) 

462 

463 tocsc.__doc__ = spmatrix.tocsc.__doc__ 

464 

465 def tocoo(self, copy=True): 

466 """Convert this matrix to COOrdinate format. 

467 

468 When copy=False the data array will be shared between 

469 this matrix and the resultant coo_matrix. 

470 """ 

471 

472 M,N = self.shape 

473 R,C = self.blocksize 

474 

475 indptr_diff = np.diff(self.indptr) 

476 if indptr_diff.dtype.itemsize > np.dtype(np.intp).itemsize: 

477 # Check for potential overflow 

478 indptr_diff_limited = indptr_diff.astype(np.intp) 

479 if np.any(indptr_diff_limited != indptr_diff): 

480 raise ValueError("Matrix too big to convert") 

481 indptr_diff = indptr_diff_limited 

482 

483 row = (R * np.arange(M//R)).repeat(indptr_diff) 

484 row = row.repeat(R*C).reshape(-1,R,C) 

485 row += np.tile(np.arange(R).reshape(-1,1), (1,C)) 

486 row = row.reshape(-1) 

487 

488 col = (C * self.indices).repeat(R*C).reshape(-1,R,C) 

489 col += np.tile(np.arange(C), (R,1)) 

490 col = col.reshape(-1) 

491 

492 data = self.data.reshape(-1) 

493 

494 if copy: 

495 data = data.copy() 

496 

497 from .coo import coo_matrix 

498 return coo_matrix((data,(row,col)), shape=self.shape) 

499 

500 def toarray(self, order=None, out=None): 

501 return self.tocoo(copy=False).toarray(order=order, out=out) 

502 

503 toarray.__doc__ = spmatrix.toarray.__doc__ 

504 

505 def transpose(self, axes=None, copy=False): 

506 if axes is not None: 

507 raise ValueError(("Sparse matrices do not support " 

508 "an 'axes' parameter because swapping " 

509 "dimensions is the only logical permutation.")) 

510 

511 R, C = self.blocksize 

512 M, N = self.shape 

513 NBLK = self.nnz//(R*C) 

514 

515 if self.nnz == 0: 

516 return bsr_matrix((N, M), blocksize=(C, R), 

517 dtype=self.dtype, copy=copy) 

518 

519 indptr = np.empty(N//C + 1, dtype=self.indptr.dtype) 

520 indices = np.empty(NBLK, dtype=self.indices.dtype) 

521 data = np.empty((NBLK, C, R), dtype=self.data.dtype) 

522 

523 bsr_transpose(M//R, N//C, R, C, 

524 self.indptr, self.indices, self.data.ravel(), 

525 indptr, indices, data.ravel()) 

526 

527 return bsr_matrix((data, indices, indptr), 

528 shape=(N, M), copy=copy) 

529 

530 transpose.__doc__ = spmatrix.transpose.__doc__ 

531 

532 ############################################################## 

533 # methods that examine or modify the internal data structure # 

534 ############################################################## 

535 

536 def eliminate_zeros(self): 

537 """Remove zero elements in-place.""" 

538 

539 if not self.nnz: 

540 return # nothing to do 

541 

542 R,C = self.blocksize 

543 M,N = self.shape 

544 

545 mask = (self.data != 0).reshape(-1,R*C).sum(axis=1) # nonzero blocks 

546 

547 nonzero_blocks = mask.nonzero()[0] 

548 

549 self.data[:len(nonzero_blocks)] = self.data[nonzero_blocks] 

550 

551 # modifies self.indptr and self.indices *in place* 

552 _sparsetools.csr_eliminate_zeros(M//R, N//C, self.indptr, 

553 self.indices, mask) 

554 self.prune() 

555 

556 def sum_duplicates(self): 

557 """Eliminate duplicate matrix entries by adding them together 

558 

559 The is an *in place* operation 

560 """ 

561 if self.has_canonical_format: 

562 return 

563 self.sort_indices() 

564 R, C = self.blocksize 

565 M, N = self.shape 

566 

567 # port of _sparsetools.csr_sum_duplicates 

568 n_row = M // R 

569 nnz = 0 

570 row_end = 0 

571 for i in range(n_row): 

572 jj = row_end 

573 row_end = self.indptr[i+1] 

574 while jj < row_end: 

575 j = self.indices[jj] 

576 x = self.data[jj] 

577 jj += 1 

578 while jj < row_end and self.indices[jj] == j: 

579 x += self.data[jj] 

580 jj += 1 

581 self.indices[nnz] = j 

582 self.data[nnz] = x 

583 nnz += 1 

584 self.indptr[i+1] = nnz 

585 

586 self.prune() # nnz may have changed 

587 self.has_canonical_format = True 

588 

589 def sort_indices(self): 

590 """Sort the indices of this matrix *in place* 

591 """ 

592 if self.has_sorted_indices: 

593 return 

594 

595 R,C = self.blocksize 

596 M,N = self.shape 

597 

598 bsr_sort_indices(M//R, N//C, R, C, self.indptr, self.indices, self.data.ravel()) 

599 

600 self.has_sorted_indices = True 

601 

602 def prune(self): 

603 """ Remove empty space after all non-zero elements. 

604 """ 

605 

606 R,C = self.blocksize 

607 M,N = self.shape 

608 

609 if len(self.indptr) != M//R + 1: 

610 raise ValueError("index pointer has invalid length") 

611 

612 bnnz = self.indptr[-1] 

613 

614 if len(self.indices) < bnnz: 

615 raise ValueError("indices array has too few elements") 

616 if len(self.data) < bnnz: 

617 raise ValueError("data array has too few elements") 

618 

619 self.data = self.data[:bnnz] 

620 self.indices = self.indices[:bnnz] 

621 

622 # utility functions 

623 def _binopt(self, other, op, in_shape=None, out_shape=None): 

624 """Apply the binary operation fn to two sparse matrices.""" 

625 

626 # Ideally we'd take the GCDs of the blocksize dimensions 

627 # and explode self and other to match. 

628 other = self.__class__(other, blocksize=self.blocksize) 

629 

630 # e.g. bsr_plus_bsr, etc. 

631 fn = getattr(_sparsetools, self.format + op + self.format) 

632 

633 R,C = self.blocksize 

634 

635 max_bnnz = len(self.data) + len(other.data) 

636 idx_dtype = get_index_dtype((self.indptr, self.indices, 

637 other.indptr, other.indices), 

638 maxval=max_bnnz) 

639 indptr = np.empty(self.indptr.shape, dtype=idx_dtype) 

640 indices = np.empty(max_bnnz, dtype=idx_dtype) 

641 

642 bool_ops = ['_ne_', '_lt_', '_gt_', '_le_', '_ge_'] 

643 if op in bool_ops: 

644 data = np.empty(R*C*max_bnnz, dtype=np.bool_) 

645 else: 

646 data = np.empty(R*C*max_bnnz, dtype=upcast(self.dtype,other.dtype)) 

647 

648 fn(self.shape[0]//R, self.shape[1]//C, R, C, 

649 self.indptr.astype(idx_dtype), 

650 self.indices.astype(idx_dtype), 

651 self.data, 

652 other.indptr.astype(idx_dtype), 

653 other.indices.astype(idx_dtype), 

654 np.ravel(other.data), 

655 indptr, 

656 indices, 

657 data) 

658 

659 actual_bnnz = indptr[-1] 

660 indices = indices[:actual_bnnz] 

661 data = data[:R*C*actual_bnnz] 

662 

663 if actual_bnnz < max_bnnz/2: 

664 indices = indices.copy() 

665 data = data.copy() 

666 

667 data = data.reshape(-1,R,C) 

668 

669 return self.__class__((data, indices, indptr), shape=self.shape) 

670 

671 # needed by _data_matrix 

672 def _with_data(self,data,copy=True): 

673 """Returns a matrix with the same sparsity structure as self, 

674 but with different data. By default the structure arrays 

675 (i.e. .indptr and .indices) are copied. 

676 """ 

677 if copy: 

678 return self.__class__((data,self.indices.copy(),self.indptr.copy()), 

679 shape=self.shape,dtype=data.dtype) 

680 else: 

681 return self.__class__((data,self.indices,self.indptr), 

682 shape=self.shape,dtype=data.dtype) 

683 

684# # these functions are used by the parent class 

685# # to remove redudancy between bsc_matrix and bsr_matrix 

686# def _swap(self,x): 

687# """swap the members of x if this is a column-oriented matrix 

688# """ 

689# return (x[0],x[1]) 

690 

691 

692def isspmatrix_bsr(x): 

693 """Is x of a bsr_matrix type? 

694 

695 Parameters 

696 ---------- 

697 x 

698 object to check for being a bsr matrix 

699 

700 Returns 

701 ------- 

702 bool 

703 True if x is a bsr matrix, False otherwise 

704 

705 Examples 

706 -------- 

707 >>> from scipy.sparse import bsr_matrix, isspmatrix_bsr 

708 >>> isspmatrix_bsr(bsr_matrix([[5]])) 

709 True 

710 

711 >>> from scipy.sparse import bsr_matrix, csr_matrix, isspmatrix_bsr 

712 >>> isspmatrix_bsr(csr_matrix([[5]])) 

713 False 

714 """ 

715 return isinstance(x, bsr_matrix)