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"""Dictionary Of Keys based matrix""" 

2 

3__docformat__ = "restructuredtext en" 

4 

5__all__ = ['dok_matrix', 'isspmatrix_dok'] 

6 

7import itertools 

8import numpy as np 

9 

10from .base import spmatrix, isspmatrix 

11from ._index import IndexMixin 

12from .sputils import (isdense, getdtype, isshape, isintlike, isscalarlike, 

13 upcast, upcast_scalar, get_index_dtype, check_shape) 

14 

15try: 

16 from operator import isSequenceType as _is_sequence 

17except ImportError: 

18 def _is_sequence(x): 

19 return (hasattr(x, '__len__') or hasattr(x, '__next__') 

20 or hasattr(x, 'next')) 

21 

22 

23class dok_matrix(spmatrix, IndexMixin, dict): 

24 """ 

25 Dictionary Of Keys based sparse matrix. 

26 

27 This is an efficient structure for constructing sparse 

28 matrices incrementally. 

29 

30 This can be instantiated in several ways: 

31 dok_matrix(D) 

32 with a dense matrix, D 

33 

34 dok_matrix(S) 

35 with a sparse matrix, S 

36 

37 dok_matrix((M,N), [dtype]) 

38 create the matrix with initial shape (M,N) 

39 dtype is optional, defaulting to dtype='d' 

40 

41 Attributes 

42 ---------- 

43 dtype : dtype 

44 Data type of the matrix 

45 shape : 2-tuple 

46 Shape of the matrix 

47 ndim : int 

48 Number of dimensions (this is always 2) 

49 nnz 

50 Number of nonzero elements 

51 

52 Notes 

53 ----- 

54 

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

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

57 

58 Allows for efficient O(1) access of individual elements. 

59 Duplicates are not allowed. 

60 Can be efficiently converted to a coo_matrix once constructed. 

61 

62 Examples 

63 -------- 

64 >>> import numpy as np 

65 >>> from scipy.sparse import dok_matrix 

66 >>> S = dok_matrix((5, 5), dtype=np.float32) 

67 >>> for i in range(5): 

68 ... for j in range(5): 

69 ... S[i, j] = i + j # Update element 

70 

71 """ 

72 format = 'dok' 

73 

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

75 dict.__init__(self) 

76 spmatrix.__init__(self) 

77 

78 self.dtype = getdtype(dtype, default=float) 

79 if isinstance(arg1, tuple) and isshape(arg1): # (M,N) 

80 M, N = arg1 

81 self._shape = check_shape((M, N)) 

82 elif isspmatrix(arg1): # Sparse ctor 

83 if isspmatrix_dok(arg1) and copy: 

84 arg1 = arg1.copy() 

85 else: 

86 arg1 = arg1.todok() 

87 

88 if dtype is not None: 

89 arg1 = arg1.astype(dtype, copy=False) 

90 

91 dict.update(self, arg1) 

92 self._shape = check_shape(arg1.shape) 

93 self.dtype = arg1.dtype 

94 else: # Dense ctor 

95 try: 

96 arg1 = np.asarray(arg1) 

97 except Exception: 

98 raise TypeError('Invalid input format.') 

99 

100 if len(arg1.shape) != 2: 

101 raise TypeError('Expected rank <=2 dense array or matrix.') 

102 

103 from .coo import coo_matrix 

104 d = coo_matrix(arg1, dtype=dtype).todok() 

105 dict.update(self, d) 

106 self._shape = check_shape(arg1.shape) 

107 self.dtype = d.dtype 

108 

109 def update(self, val): 

110 # Prevent direct usage of update 

111 raise NotImplementedError("Direct modification to dok_matrix element " 

112 "is not allowed.") 

113 

114 def _update(self, data): 

115 """An update method for dict data defined for direct access to 

116 `dok_matrix` data. Main purpose is to be used for effcient conversion 

117 from other spmatrix classes. Has no checking if `data` is valid.""" 

118 return dict.update(self, data) 

119 

120 def set_shape(self, shape): 

121 new_matrix = self.reshape(shape, copy=False).asformat(self.format) 

122 self.__dict__ = new_matrix.__dict__ 

123 dict.clear(self) 

124 dict.update(self, new_matrix) 

125 

126 shape = property(fget=spmatrix.get_shape, fset=set_shape) 

127 

128 def getnnz(self, axis=None): 

129 if axis is not None: 

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

131 "for DOK format.") 

132 return dict.__len__(self) 

133 

134 def count_nonzero(self): 

135 return sum(x != 0 for x in self.values()) 

136 

137 getnnz.__doc__ = spmatrix.getnnz.__doc__ 

138 count_nonzero.__doc__ = spmatrix.count_nonzero.__doc__ 

139 

140 def __len__(self): 

141 return dict.__len__(self) 

142 

143 def get(self, key, default=0.): 

144 """This overrides the dict.get method, providing type checking 

145 but otherwise equivalent functionality. 

146 """ 

147 try: 

148 i, j = key 

149 assert isintlike(i) and isintlike(j) 

150 except (AssertionError, TypeError, ValueError): 

151 raise IndexError('Index must be a pair of integers.') 

152 if (i < 0 or i >= self.shape[0] or j < 0 or j >= self.shape[1]): 

153 raise IndexError('Index out of bounds.') 

154 return dict.get(self, key, default) 

155 

156 def _get_intXint(self, row, col): 

157 return dict.get(self, (row, col), self.dtype.type(0)) 

158 

159 def _get_intXslice(self, row, col): 

160 return self._get_sliceXslice(slice(row, row+1), col) 

161 

162 def _get_sliceXint(self, row, col): 

163 return self._get_sliceXslice(row, slice(col, col+1)) 

164 

165 def _get_sliceXslice(self, row, col): 

166 row_start, row_stop, row_step = row.indices(self.shape[0]) 

167 col_start, col_stop, col_step = col.indices(self.shape[1]) 

168 row_range = range(row_start, row_stop, row_step) 

169 col_range = range(col_start, col_stop, col_step) 

170 shape = (len(row_range), len(col_range)) 

171 # Switch paths only when advantageous 

172 # (count the iterations in the loops, adjust for complexity) 

173 if len(self) >= 2 * shape[0] * shape[1]: 

174 # O(nr*nc) path: loop over <row x col> 

175 return self._get_columnXarray(row_range, col_range) 

176 # O(nnz) path: loop over entries of self 

177 newdok = dok_matrix(shape, dtype=self.dtype) 

178 for key in self.keys(): 

179 i, ri = divmod(int(key[0]) - row_start, row_step) 

180 if ri != 0 or i < 0 or i >= shape[0]: 

181 continue 

182 j, rj = divmod(int(key[1]) - col_start, col_step) 

183 if rj != 0 or j < 0 or j >= shape[1]: 

184 continue 

185 x = dict.__getitem__(self, key) 

186 dict.__setitem__(newdok, (i, j), x) 

187 return newdok 

188 

189 def _get_intXarray(self, row, col): 

190 return self._get_columnXarray([row], col) 

191 

192 def _get_arrayXint(self, row, col): 

193 return self._get_columnXarray(row, [col]) 

194 

195 def _get_sliceXarray(self, row, col): 

196 row = list(range(*row.indices(self.shape[0]))) 

197 return self._get_columnXarray(row, col) 

198 

199 def _get_arrayXslice(self, row, col): 

200 col = list(range(*col.indices(self.shape[1]))) 

201 return self._get_columnXarray(row, col) 

202 

203 def _get_columnXarray(self, row, col): 

204 # outer indexing 

205 newdok = dok_matrix((len(row), len(col)), dtype=self.dtype) 

206 

207 for i, r in enumerate(row): 

208 for j, c in enumerate(col): 

209 v = dict.get(self, (r, c), 0) 

210 if v: 

211 dict.__setitem__(newdok, (i, j), v) 

212 return newdok 

213 

214 def _get_arrayXarray(self, row, col): 

215 # inner indexing 

216 i, j = map(np.atleast_2d, np.broadcast_arrays(row, col)) 

217 newdok = dok_matrix(i.shape, dtype=self.dtype) 

218 

219 for key in itertools.product(range(i.shape[0]), range(i.shape[1])): 

220 v = dict.get(self, (i[key], j[key]), 0) 

221 if v: 

222 dict.__setitem__(newdok, key, v) 

223 return newdok 

224 

225 def _set_intXint(self, row, col, x): 

226 key = (row, col) 

227 if x: 

228 dict.__setitem__(self, key, x) 

229 elif dict.__contains__(self, key): 

230 del self[key] 

231 

232 def _set_arrayXarray(self, row, col, x): 

233 row = list(map(int, row.ravel())) 

234 col = list(map(int, col.ravel())) 

235 x = x.ravel() 

236 dict.update(self, zip(zip(row, col), x)) 

237 

238 for i in np.nonzero(x == 0)[0]: 

239 key = (row[i], col[i]) 

240 if dict.__getitem__(self, key) == 0: 

241 # may have been superseded by later update 

242 del self[key] 

243 

244 def __add__(self, other): 

245 if isscalarlike(other): 

246 res_dtype = upcast_scalar(self.dtype, other) 

247 new = dok_matrix(self.shape, dtype=res_dtype) 

248 # Add this scalar to every element. 

249 M, N = self.shape 

250 for key in itertools.product(range(M), range(N)): 

251 aij = dict.get(self, (key), 0) + other 

252 if aij: 

253 new[key] = aij 

254 # new.dtype.char = self.dtype.char 

255 elif isspmatrix_dok(other): 

256 if other.shape != self.shape: 

257 raise ValueError("Matrix dimensions are not equal.") 

258 # We could alternatively set the dimensions to the largest of 

259 # the two matrices to be summed. Would this be a good idea? 

260 res_dtype = upcast(self.dtype, other.dtype) 

261 new = dok_matrix(self.shape, dtype=res_dtype) 

262 dict.update(new, self) 

263 with np.errstate(over='ignore'): 

264 dict.update(new, 

265 ((k, new[k] + other[k]) for k in other.keys())) 

266 elif isspmatrix(other): 

267 csc = self.tocsc() 

268 new = csc + other 

269 elif isdense(other): 

270 new = self.todense() + other 

271 else: 

272 return NotImplemented 

273 return new 

274 

275 def __radd__(self, other): 

276 if isscalarlike(other): 

277 new = dok_matrix(self.shape, dtype=self.dtype) 

278 M, N = self.shape 

279 for key in itertools.product(range(M), range(N)): 

280 aij = dict.get(self, (key), 0) + other 

281 if aij: 

282 new[key] = aij 

283 elif isspmatrix_dok(other): 

284 if other.shape != self.shape: 

285 raise ValueError("Matrix dimensions are not equal.") 

286 new = dok_matrix(self.shape, dtype=self.dtype) 

287 dict.update(new, self) 

288 dict.update(new, 

289 ((k, self[k] + other[k]) for k in other.keys())) 

290 elif isspmatrix(other): 

291 csc = self.tocsc() 

292 new = csc + other 

293 elif isdense(other): 

294 new = other + self.todense() 

295 else: 

296 return NotImplemented 

297 return new 

298 

299 def __neg__(self): 

300 if self.dtype.kind == 'b': 

301 raise NotImplementedError('Negating a sparse boolean matrix is not' 

302 ' supported.') 

303 new = dok_matrix(self.shape, dtype=self.dtype) 

304 dict.update(new, ((k, -self[k]) for k in self.keys())) 

305 return new 

306 

307 def _mul_scalar(self, other): 

308 res_dtype = upcast_scalar(self.dtype, other) 

309 # Multiply this scalar by every element. 

310 new = dok_matrix(self.shape, dtype=res_dtype) 

311 dict.update(new, ((k, v * other) for k, v in self.items())) 

312 return new 

313 

314 def _mul_vector(self, other): 

315 # matrix * vector 

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

317 for (i, j), v in self.items(): 

318 result[i] += v * other[j] 

319 return result 

320 

321 def _mul_multivector(self, other): 

322 # matrix * multivector 

323 result_shape = (self.shape[0], other.shape[1]) 

324 result_dtype = upcast(self.dtype, other.dtype) 

325 result = np.zeros(result_shape, dtype=result_dtype) 

326 for (i, j), v in self.items(): 

327 result[i,:] += v * other[j,:] 

328 return result 

329 

330 def __imul__(self, other): 

331 if isscalarlike(other): 

332 dict.update(self, ((k, v * other) for k, v in self.items())) 

333 return self 

334 return NotImplemented 

335 

336 def __truediv__(self, other): 

337 if isscalarlike(other): 

338 res_dtype = upcast_scalar(self.dtype, other) 

339 new = dok_matrix(self.shape, dtype=res_dtype) 

340 dict.update(new, ((k, v / other) for k, v in self.items())) 

341 return new 

342 return self.tocsr() / other 

343 

344 def __itruediv__(self, other): 

345 if isscalarlike(other): 

346 dict.update(self, ((k, v / other) for k, v in self.items())) 

347 return self 

348 return NotImplemented 

349 

350 def __reduce__(self): 

351 # this approach is necessary because __setstate__ is called after 

352 # __setitem__ upon unpickling and since __init__ is not called there 

353 # is no shape attribute hence it is not possible to unpickle it. 

354 return dict.__reduce__(self) 

355 

356 # What should len(sparse) return? For consistency with dense matrices, 

357 # perhaps it should be the number of rows? For now it returns the number 

358 # of non-zeros. 

359 

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

361 if axes is not None: 

362 raise ValueError("Sparse matrices do not support " 

363 "an 'axes' parameter because swapping " 

364 "dimensions is the only logical permutation.") 

365 

366 M, N = self.shape 

367 new = dok_matrix((N, M), dtype=self.dtype, copy=copy) 

368 dict.update(new, (((right, left), val) 

369 for (left, right), val in self.items())) 

370 return new 

371 

372 transpose.__doc__ = spmatrix.transpose.__doc__ 

373 

374 def conjtransp(self): 

375 """Return the conjugate transpose.""" 

376 M, N = self.shape 

377 new = dok_matrix((N, M), dtype=self.dtype) 

378 dict.update(new, (((right, left), np.conj(val)) 

379 for (left, right), val in self.items())) 

380 return new 

381 

382 def copy(self): 

383 new = dok_matrix(self.shape, dtype=self.dtype) 

384 dict.update(new, self) 

385 return new 

386 

387 copy.__doc__ = spmatrix.copy.__doc__ 

388 

389 def tocoo(self, copy=False): 

390 from .coo import coo_matrix 

391 if self.nnz == 0: 

392 return coo_matrix(self.shape, dtype=self.dtype) 

393 

394 idx_dtype = get_index_dtype(maxval=max(self.shape)) 

395 data = np.fromiter(self.values(), dtype=self.dtype, count=self.nnz) 

396 row = np.fromiter((i for i, _ in self.keys()), dtype=idx_dtype, count=self.nnz) 

397 col = np.fromiter((j for _, j in self.keys()), dtype=idx_dtype, count=self.nnz) 

398 A = coo_matrix((data, (row, col)), shape=self.shape, dtype=self.dtype) 

399 A.has_canonical_format = True 

400 return A 

401 

402 tocoo.__doc__ = spmatrix.tocoo.__doc__ 

403 

404 def todok(self, copy=False): 

405 if copy: 

406 return self.copy() 

407 return self 

408 

409 todok.__doc__ = spmatrix.todok.__doc__ 

410 

411 def tocsc(self, copy=False): 

412 return self.tocoo(copy=False).tocsc(copy=copy) 

413 

414 tocsc.__doc__ = spmatrix.tocsc.__doc__ 

415 

416 def resize(self, *shape): 

417 shape = check_shape(shape) 

418 newM, newN = shape 

419 M, N = self.shape 

420 if newM < M or newN < N: 

421 # Remove all elements outside new dimensions 

422 for (i, j) in list(self.keys()): 

423 if i >= newM or j >= newN: 

424 del self[i, j] 

425 self._shape = shape 

426 

427 resize.__doc__ = spmatrix.resize.__doc__ 

428 

429 

430def isspmatrix_dok(x): 

431 """Is x of dok_matrix type? 

432 

433 Parameters 

434 ---------- 

435 x 

436 object to check for being a dok matrix 

437 

438 Returns 

439 ------- 

440 bool 

441 True if x is a dok matrix, False otherwise 

442 

443 Examples 

444 -------- 

445 >>> from scipy.sparse import dok_matrix, isspmatrix_dok 

446 >>> isspmatrix_dok(dok_matrix([[5]])) 

447 True 

448 

449 >>> from scipy.sparse import dok_matrix, csr_matrix, isspmatrix_dok 

450 >>> isspmatrix_dok(csr_matrix([[5]])) 

451 False 

452 """ 

453 return isinstance(x, dok_matrix)