Coverage for /usr/lib/python3/dist-packages/scipy/sparse/_bsr.py: 16%

316 statements  

« prev     ^ index     » next       coverage.py v7.4.4, created at 2025-06-14 15:25 +0200

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

2 

3__docformat__ = "restructuredtext en" 

4 

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

6 

7from warnings import warn 

8 

9import numpy as np 

10 

11from ._matrix import spmatrix, _array_doc_to_matrix 

12from ._data import _data_matrix, _minmax_mixin 

13from ._compressed import _cs_matrix 

14from ._base import issparse, _formats, _spbase, sparray 

15from ._sputils import (isshape, getdtype, getdata, to_native, upcast, 

16 check_shape) 

17from . import _sparsetools 

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

19 bsr_matmat, bsr_transpose, bsr_sort_indices, 

20 bsr_tocsr) 

21 

22 

23class _bsr_base(_cs_matrix, _minmax_mixin): 

24 """Block Sparse Row format sparse array. 

25 

26 This can be instantiated in several ways: 

27 bsr_array(D, [blocksize=(R,C)]) 

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

29 

30 bsr_array(S, [blocksize=(R,C)]) 

31 with another sparse array S (equivalent to S.tobsr()) 

32 

33 bsr_array((M, N), [blocksize=(R,C), dtype]) 

34 to construct an empty sparse array with shape (M, N) 

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

36 

37 bsr_array((data, ij), [blocksize=(R,C), shape=(M, N)]) 

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

39 

40 bsr_array((data, indices, indptr), [shape=(M, N)]) 

41 is the standard BSR representation where the block column 

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

43 and their corresponding block values are stored in 

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

45 supplied, the array dimensions are inferred from the index arrays. 

46 

47 Attributes 

48 ---------- 

49 dtype : dtype 

50 Data type of the array 

51 shape : 2-tuple 

52 Shape of the array 

53 ndim : int 

54 Number of dimensions (this is always 2) 

55 nnz 

56 Number of stored values, including explicit zeros 

57 data 

58 Data array 

59 indices 

60 BSR format index array 

61 indptr 

62 BSR format index pointer array 

63 blocksize 

64 Block size 

65 has_sorted_indices 

66 Whether indices are sorted 

67 

68 Notes 

69 ----- 

70 Sparse arrays can be used in arithmetic operations: they support 

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

72 

73 **Summary of BSR format** 

74 

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

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

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

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

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

80 operations. 

81 

82 **Blocksize** 

83 

84 The blocksize (R,C) must evenly divide the shape of the sparse array (M,N). 

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

86 ``N % C = 0``. 

87 

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

89 an appropriate blocksize. 

90 

91 **Canonical Format** 

92 

93 In canonical format, there are no duplicate blocks and indices are sorted 

94 per row. 

95 

96 Examples 

97 -------- 

98 >>> from scipy.sparse import bsr_array 

99 >>> import numpy as np 

100 >>> bsr_array((3, 4), dtype=np.int8).toarray() 

101 array([[0, 0, 0, 0], 

102 [0, 0, 0, 0], 

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

104 

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

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

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

108 >>> bsr_array((data, (row, col)), shape=(3, 3)).toarray() 

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

110 [0, 0, 3], 

111 [4, 5, 6]]) 

112 

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

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

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

116 >>> bsr_array((data,indices,indptr), shape=(6, 6)).toarray() 

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

118 [1, 1, 0, 0, 2, 2], 

119 [0, 0, 0, 0, 3, 3], 

120 [0, 0, 0, 0, 3, 3], 

121 [4, 4, 5, 5, 6, 6], 

122 [4, 4, 5, 5, 6, 6]]) 

123 

124 """ 

125 _format = 'bsr' 

126 

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

128 _data_matrix.__init__(self) 

129 

130 if issparse(arg1): 

131 if arg1.format == self.format and copy: 

132 arg1 = arg1.copy() 

133 else: 

134 arg1 = arg1.tobsr(blocksize=blocksize) 

135 self._set_self(arg1) 

136 

137 elif isinstance(arg1,tuple): 

138 if isshape(arg1): 

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

140 self._shape = check_shape(arg1) 

141 M,N = self.shape 

142 # process blocksize 

143 if blocksize is None: 

144 blocksize = (1,1) 

145 else: 

146 if not isshape(blocksize): 

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

148 blocksize = tuple(blocksize) 

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

150 

151 R,C = blocksize 

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

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

154 

155 # Select index dtype large enough to pass array and 

156 # scalar parameters to sparsetools 

157 idx_dtype = self._get_index_dtype(maxval=max(M//R, N//C, R, C)) 

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

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

160 

161 elif len(arg1) == 2: 

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

163 self._set_self( 

164 self._coo_container(arg1, dtype=dtype, shape=shape).tobsr( 

165 blocksize=blocksize 

166 ) 

167 ) 

168 

169 elif len(arg1) == 3: 

170 # (data,indices,indptr) format 

171 (data, indices, indptr) = arg1 

172 

173 # Select index dtype large enough to pass array and 

174 # scalar parameters to sparsetools 

175 maxval = 1 

176 if shape is not None: 

177 maxval = max(shape) 

178 if blocksize is not None: 

179 maxval = max(maxval, max(blocksize)) 

180 idx_dtype = self._get_index_dtype((indices, indptr), maxval=maxval, 

181 check_contents=True) 

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

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

184 self.data = getdata(data, copy=copy, dtype=dtype) 

185 if self.data.ndim != 3: 

186 raise ValueError( 

187 'BSR data must be 3-dimensional, got shape={}'.format( 

188 self.data.shape,)) 

189 if blocksize is not None: 

190 if not isshape(blocksize): 

191 raise ValueError(f'invalid blocksize={blocksize}') 

192 if tuple(blocksize) != self.data.shape[1:]: 

193 raise ValueError('mismatching blocksize={} vs {}'.format( 

194 blocksize, self.data.shape[1:])) 

195 else: 

196 raise ValueError('unrecognized bsr_array constructor usage') 

197 else: 

198 # must be dense 

199 try: 

200 arg1 = np.asarray(arg1) 

201 except Exception as e: 

202 raise ValueError("unrecognized form for" 

203 " %s_matrix constructor" % self.format) from e 

204 arg1 = self._coo_container( 

205 arg1, dtype=dtype 

206 ).tobsr(blocksize=blocksize) 

207 self._set_self(arg1) 

208 

209 if shape is not None: 

210 self._shape = check_shape(shape) 

211 else: 

212 if self.shape is None: 

213 # shape not already set, try to infer dimensions 

214 try: 

215 M = len(self.indptr) - 1 

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

217 except Exception as e: 

218 raise ValueError('unable to infer matrix dimensions') from e 

219 else: 

220 R,C = self.blocksize 

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

222 

223 if self.shape is None: 

224 if shape is None: 

225 # TODO infer shape here 

226 raise ValueError('need to infer shape') 

227 else: 

228 self._shape = check_shape(shape) 

229 

230 if dtype is not None: 

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

232 

233 self.check_format(full_check=False) 

234 

235 def check_format(self, full_check=True): 

236 """check whether the matrix format is valid 

237 

238 *Parameters*: 

239 full_check: 

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

241 False - basic check, O(1) operations 

242 

243 """ 

244 M,N = self.shape 

245 R,C = self.blocksize 

246 

247 # index arrays should have integer data types 

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

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

250 % self.indptr.dtype.name) 

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

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

253 % self.indices.dtype.name) 

254 

255 idx_dtype = self._get_index_dtype((self.indices, self.indptr)) 

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

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

258 self.data = to_native(self.data) 

259 

260 # check array shapes 

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

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

263 if self.data.ndim != 3: 

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

265 

266 # check index pointer 

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

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

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

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

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

272 

273 # check index and data arrays 

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

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

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

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

278 "the size of index and data arrays") 

279 

280 self.prune() 

281 

282 if full_check: 

283 # check format validity (more expensive) 

284 if self.nnz > 0: 

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

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

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

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

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

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

291 "non-decreasing sequence") 

292 

293 # if not self.has_sorted_indices(): 

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

295 # self.sort_indices(check_first=False) 

296 

297 def _get_blocksize(self): 

298 return self.data.shape[1:] 

299 blocksize = property(fget=_get_blocksize) 

300 

301 def _getnnz(self, axis=None): 

302 if axis is not None: 

303 raise NotImplementedError("_getnnz over an axis is not implemented " 

304 "for BSR format") 

305 R,C = self.blocksize 

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

307 

308 _getnnz.__doc__ = _spbase._getnnz.__doc__ 

309 

310 def __repr__(self): 

311 format = _formats[self.format][1] 

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

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

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

315 (format,))) 

316 

317 def diagonal(self, k=0): 

318 rows, cols = self.shape 

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

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

321 R, C = self.blocksize 

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

323 dtype=upcast(self.dtype)) 

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

325 self.indptr, self.indices, 

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

327 return y 

328 

329 diagonal.__doc__ = _spbase.diagonal.__doc__ 

330 

331 ########################## 

332 # NotImplemented methods # 

333 ########################## 

334 

335 def __getitem__(self,key): 

336 raise NotImplementedError 

337 

338 def __setitem__(self,key,val): 

339 raise NotImplementedError 

340 

341 ###################### 

342 # Arithmetic methods # 

343 ###################### 

344 

345 def _add_dense(self, other): 

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

347 

348 def _mul_vector(self, other): 

349 M,N = self.shape 

350 R,C = self.blocksize 

351 

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

353 

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

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

356 other, result) 

357 

358 return result 

359 

360 def _mul_multivector(self,other): 

361 R,C = self.blocksize 

362 M,N = self.shape 

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

364 

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

366 

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

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

369 other.ravel(), result.ravel()) 

370 

371 return result 

372 

373 def _mul_sparse_matrix(self, other): 

374 M, K1 = self.shape 

375 K2, N = other.shape 

376 

377 R,n = self.blocksize 

378 

379 # convert to this format 

380 if other.format == "bsr": 

381 C = other.blocksize[1] 

382 else: 

383 C = 1 

384 

385 if other.format == "csr" and n == 1: 

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

387 else: 

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

389 

390 idx_dtype = self._get_index_dtype((self.indptr, self.indices, 

391 other.indptr, other.indices)) 

392 

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

394 self.indptr.astype(idx_dtype), 

395 self.indices.astype(idx_dtype), 

396 other.indptr.astype(idx_dtype), 

397 other.indices.astype(idx_dtype)) 

398 

399 idx_dtype = self._get_index_dtype((self.indptr, self.indices, 

400 other.indptr, other.indices), 

401 maxval=bnnz) 

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

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

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

405 

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

407 self.indptr.astype(idx_dtype), 

408 self.indices.astype(idx_dtype), 

409 np.ravel(self.data), 

410 other.indptr.astype(idx_dtype), 

411 other.indices.astype(idx_dtype), 

412 np.ravel(other.data), 

413 indptr, 

414 indices, 

415 data) 

416 

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

418 

419 # TODO eliminate zeros 

420 

421 return self._bsr_container( 

422 (data, indices, indptr), shape=(M, N), blocksize=(R, C) 

423 ) 

424 

425 ###################### 

426 # Conversion methods # 

427 ###################### 

428 

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

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

431 

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

433 matrix and the resultant bsr_array. 

434 

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

436 block size of the bsr_array. 

437 """ 

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

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

440 if copy: 

441 return self.copy() 

442 else: 

443 return self 

444 

445 def tocsr(self, copy=False): 

446 M, N = self.shape 

447 R, C = self.blocksize 

448 nnz = self.nnz 

449 idx_dtype = self._get_index_dtype((self.indptr, self.indices), 

450 maxval=max(nnz, N)) 

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

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

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

454 

455 bsr_tocsr(M // R, # n_brow 

456 N // C, # n_bcol 

457 R, C, 

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

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

460 self.data, 

461 indptr, 

462 indices, 

463 data) 

464 return self._csr_container((data, indices, indptr), shape=self.shape) 

465 

466 tocsr.__doc__ = _spbase.tocsr.__doc__ 

467 

468 def tocsc(self, copy=False): 

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

470 

471 tocsc.__doc__ = _spbase.tocsc.__doc__ 

472 

473 def tocoo(self, copy=True): 

474 """Convert this matrix to COOrdinate format. 

475 

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

477 this matrix and the resultant coo_matrix. 

478 """ 

479 

480 M,N = self.shape 

481 R,C = self.blocksize 

482 

483 indptr_diff = np.diff(self.indptr) 

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

485 # Check for potential overflow 

486 indptr_diff_limited = indptr_diff.astype(np.intp) 

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

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

489 indptr_diff = indptr_diff_limited 

490 

491 idx_dtype = self._get_index_dtype(maxval=max(M, N)) 

492 row = (R * np.arange(M//R, dtype=idx_dtype)).repeat(indptr_diff) 

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

494 row += np.tile(np.arange(R, dtype=idx_dtype).reshape(-1,1), (1,C)) 

495 row = row.reshape(-1) 

496 

497 col = (C * self.indices).astype(idx_dtype, copy=False).repeat(R*C).reshape(-1,R,C) 

498 col += np.tile(np.arange(C, dtype=idx_dtype), (R,1)) 

499 col = col.reshape(-1) 

500 

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

502 

503 if copy: 

504 data = data.copy() 

505 

506 return self._coo_container( 

507 (data, (row, col)), shape=self.shape 

508 ) 

509 

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

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

512 

513 toarray.__doc__ = _spbase.toarray.__doc__ 

514 

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

516 if axes is not None: 

517 raise ValueError("Sparse matrices do not support " 

518 "an 'axes' parameter because swapping " 

519 "dimensions is the only logical permutation.") 

520 

521 R, C = self.blocksize 

522 M, N = self.shape 

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

524 

525 if self.nnz == 0: 

526 return self._bsr_container((N, M), blocksize=(C, R), 

527 dtype=self.dtype, copy=copy) 

528 

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

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

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

532 

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

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

535 indptr, indices, data.ravel()) 

536 

537 return self._bsr_container((data, indices, indptr), 

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

539 

540 transpose.__doc__ = _spbase.transpose.__doc__ 

541 

542 ############################################################## 

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

544 ############################################################## 

545 

546 def eliminate_zeros(self): 

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

548 

549 if not self.nnz: 

550 return # nothing to do 

551 

552 R,C = self.blocksize 

553 M,N = self.shape 

554 

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

556 

557 nonzero_blocks = mask.nonzero()[0] 

558 

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

560 

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

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

563 self.indices, mask) 

564 self.prune() 

565 

566 def sum_duplicates(self): 

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

568 

569 The is an *in place* operation 

570 """ 

571 if self.has_canonical_format: 

572 return 

573 self.sort_indices() 

574 R, C = self.blocksize 

575 M, N = self.shape 

576 

577 # port of _sparsetools.csr_sum_duplicates 

578 n_row = M // R 

579 nnz = 0 

580 row_end = 0 

581 for i in range(n_row): 

582 jj = row_end 

583 row_end = self.indptr[i+1] 

584 while jj < row_end: 

585 j = self.indices[jj] 

586 x = self.data[jj] 

587 jj += 1 

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

589 x += self.data[jj] 

590 jj += 1 

591 self.indices[nnz] = j 

592 self.data[nnz] = x 

593 nnz += 1 

594 self.indptr[i+1] = nnz 

595 

596 self.prune() # nnz may have changed 

597 self.has_canonical_format = True 

598 

599 def sort_indices(self): 

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

601 """ 

602 if self.has_sorted_indices: 

603 return 

604 

605 R,C = self.blocksize 

606 M,N = self.shape 

607 

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

609 

610 self.has_sorted_indices = True 

611 

612 def prune(self): 

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

614 """ 

615 

616 R,C = self.blocksize 

617 M,N = self.shape 

618 

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

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

621 

622 bnnz = self.indptr[-1] 

623 

624 if len(self.indices) < bnnz: 

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

626 if len(self.data) < bnnz: 

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

628 

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

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

631 

632 # utility functions 

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

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

635 

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

637 # and explode self and other to match. 

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

639 

640 # e.g. bsr_plus_bsr, etc. 

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

642 

643 R,C = self.blocksize 

644 

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

646 idx_dtype = self._get_index_dtype((self.indptr, self.indices, 

647 other.indptr, other.indices), 

648 maxval=max_bnnz) 

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

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

651 

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

653 if op in bool_ops: 

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

655 else: 

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

657 

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

659 self.indptr.astype(idx_dtype), 

660 self.indices.astype(idx_dtype), 

661 self.data, 

662 other.indptr.astype(idx_dtype), 

663 other.indices.astype(idx_dtype), 

664 np.ravel(other.data), 

665 indptr, 

666 indices, 

667 data) 

668 

669 actual_bnnz = indptr[-1] 

670 indices = indices[:actual_bnnz] 

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

672 

673 if actual_bnnz < max_bnnz/2: 

674 indices = indices.copy() 

675 data = data.copy() 

676 

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

678 

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

680 

681 # needed by _data_matrix 

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

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

684 but with different data. By default the structure arrays 

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

686 """ 

687 if copy: 

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

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

690 else: 

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

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

693 

694# # these functions are used by the parent class 

695# # to remove redudancy between bsc_matrix and bsr_matrix 

696# def _swap(self,x): 

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

698# """ 

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

700 

701 

702def isspmatrix_bsr(x): 

703 """Is `x` of a bsr_matrix type? 

704 

705 Parameters 

706 ---------- 

707 x 

708 object to check for being a bsr matrix 

709 

710 Returns 

711 ------- 

712 bool 

713 True if `x` is a bsr matrix, False otherwise 

714 

715 Examples 

716 -------- 

717 >>> from scipy.sparse import bsr_array, bsr_matrix, csr_matrix, isspmatrix_bsr 

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

719 True 

720 >>> isspmatrix_bsr(bsr_array([[5]])) 

721 False 

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

723 False 

724 """ 

725 return isinstance(x, bsr_matrix) 

726 

727 

728# This namespace class separates array from matrix with isinstance 

729class bsr_array(_bsr_base, sparray): 

730 pass 

731 

732bsr_array.__doc__ = _bsr_base.__doc__ 

733 

734class bsr_matrix(spmatrix, _bsr_base): 

735 pass 

736 

737bsr_matrix.__doc__ = _array_doc_to_matrix(_bsr_base.__doc__)