Coverage for /usr/lib/python3/dist-packages/scipy/sparse/_coo.py: 16%
299 statements
« prev ^ index » next coverage.py v7.4.4, created at 2025-06-14 15:25 +0200
« prev ^ index » next coverage.py v7.4.4, created at 2025-06-14 15:25 +0200
1""" A sparse matrix in COOrdinate or 'triplet' format"""
3__docformat__ = "restructuredtext en"
5__all__ = ['coo_array', 'coo_matrix', 'isspmatrix_coo']
7from warnings import warn
9import numpy as np
11from ._matrix import spmatrix, _array_doc_to_matrix
12from ._sparsetools import coo_tocsr, coo_todense, coo_matvec
13from ._base import issparse, SparseEfficiencyWarning, _spbase, sparray
14from ._data import _data_matrix, _minmax_mixin
15from ._sputils import (upcast, upcast_char, to_native, isshape, getdtype,
16 getdata, downcast_intp_index,
17 check_shape, check_reshape_kwargs)
19import operator
22class _coo_base(_data_matrix, _minmax_mixin):
23 """
24 A sparse matrix in COOrdinate format.
26 Also known as the 'ijv' or 'triplet' format.
28 This can be instantiated in several ways:
29 coo_array(D)
30 with a dense matrix D
32 coo_array(S)
33 with another sparse matrix S (equivalent to S.tocoo())
35 coo_array((M, N), [dtype])
36 to construct an empty matrix with shape (M, N)
37 dtype is optional, defaulting to dtype='d'.
39 coo_array((data, (i, j)), [shape=(M, N)])
40 to construct from three arrays:
41 1. data[:] the entries of the matrix, in any order
42 2. i[:] the row indices of the matrix entries
43 3. j[:] the column indices of the matrix entries
45 Where ``A[i[k], j[k]] = data[k]``. When shape is not
46 specified, it is inferred from the index arrays
48 Attributes
49 ----------
50 dtype : dtype
51 Data type of the matrix
52 shape : 2-tuple
53 Shape of the matrix
54 ndim : int
55 Number of dimensions (this is always 2)
56 nnz
57 Number of stored values, including explicit zeros
58 data
59 COO format data array of the matrix
60 row
61 COO format row index array of the matrix
62 col
63 COO format column index array of the matrix
65 Notes
66 -----
68 Sparse matrices can be used in arithmetic operations: they support
69 addition, subtraction, multiplication, division, and matrix power.
71 Advantages of the COO format
72 - facilitates fast conversion among sparse formats
73 - permits duplicate entries (see example)
74 - very fast conversion to and from CSR/CSC formats
76 Disadvantages of the COO format
77 - does not directly support:
78 + arithmetic operations
79 + slicing
81 Intended Usage
82 - COO is a fast format for constructing sparse matrices
83 - Once a matrix has been constructed, convert to CSR or
84 CSC format for fast arithmetic and matrix vector operations
85 - By default when converting to CSR or CSC format, duplicate (i,j)
86 entries will be summed together. This facilitates efficient
87 construction of finite element matrices and the like. (see example)
89 Canonical format
90 - Entries and indices sorted by row, then column.
91 - There are no duplicate entries (i.e. duplicate (i,j) locations)
92 - Arrays MAY have explicit zeros.
94 Examples
95 --------
97 >>> # Constructing an empty matrix
98 >>> import numpy as np
99 >>> from scipy.sparse import coo_array
100 >>> coo_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)
105 >>> # Constructing a matrix using ijv format
106 >>> row = np.array([0, 3, 1, 0])
107 >>> col = np.array([0, 3, 1, 2])
108 >>> data = np.array([4, 5, 7, 9])
109 >>> coo_array((data, (row, col)), shape=(4, 4)).toarray()
110 array([[4, 0, 9, 0],
111 [0, 7, 0, 0],
112 [0, 0, 0, 0],
113 [0, 0, 0, 5]])
115 >>> # Constructing a matrix with duplicate indices
116 >>> row = np.array([0, 0, 1, 3, 1, 0, 0])
117 >>> col = np.array([0, 2, 1, 3, 1, 0, 0])
118 >>> data = np.array([1, 1, 1, 1, 1, 1, 1])
119 >>> coo = coo_array((data, (row, col)), shape=(4, 4))
120 >>> # Duplicate indices are maintained until implicitly or explicitly summed
121 >>> np.max(coo.data)
122 1
123 >>> coo.toarray()
124 array([[3, 0, 1, 0],
125 [0, 2, 0, 0],
126 [0, 0, 0, 0],
127 [0, 0, 0, 1]])
129 """
130 _format = 'coo'
132 def __init__(self, arg1, shape=None, dtype=None, copy=False):
133 _data_matrix.__init__(self)
135 if isinstance(arg1, tuple):
136 if isshape(arg1):
137 M, N = arg1
138 self._shape = check_shape((M, N))
139 idx_dtype = self._get_index_dtype(maxval=max(M, N))
140 data_dtype = getdtype(dtype, default=float)
141 self.row = np.array([], dtype=idx_dtype)
142 self.col = np.array([], dtype=idx_dtype)
143 self.data = np.array([], dtype=data_dtype)
144 self.has_canonical_format = True
145 else:
146 try:
147 obj, (row, col) = arg1
148 except (TypeError, ValueError) as e:
149 raise TypeError('invalid input format') from e
151 if shape is None:
152 if len(row) == 0 or len(col) == 0:
153 raise ValueError('cannot infer dimensions from zero '
154 'sized index arrays')
155 M = operator.index(np.max(row)) + 1
156 N = operator.index(np.max(col)) + 1
157 self._shape = check_shape((M, N))
158 else:
159 # Use 2 steps to ensure shape has length 2.
160 M, N = shape
161 self._shape = check_shape((M, N))
163 idx_dtype = self._get_index_dtype((row, col), maxval=max(self.shape), check_contents=True)
164 self.row = np.array(row, copy=copy, dtype=idx_dtype)
165 self.col = np.array(col, copy=copy, dtype=idx_dtype)
166 self.data = getdata(obj, copy=copy, dtype=dtype)
167 self.has_canonical_format = False
168 else:
169 if issparse(arg1):
170 if arg1.format == self.format and copy:
171 self.row = arg1.row.copy()
172 self.col = arg1.col.copy()
173 self.data = arg1.data.copy()
174 self._shape = check_shape(arg1.shape)
175 else:
176 coo = arg1.tocoo()
177 self.row = coo.row
178 self.col = coo.col
179 self.data = coo.data
180 self._shape = check_shape(coo.shape)
181 self.has_canonical_format = False
182 else:
183 #dense argument
184 M = np.atleast_2d(np.asarray(arg1))
186 if M.ndim != 2:
187 raise TypeError('expected dimension <= 2 array or matrix')
189 self._shape = check_shape(M.shape)
190 if shape is not None:
191 if check_shape(shape) != self._shape:
192 raise ValueError('inconsistent shapes: %s != %s' %
193 (shape, self._shape))
194 index_dtype = self._get_index_dtype(maxval=max(self._shape))
195 row, col = M.nonzero()
196 self.row = row.astype(index_dtype, copy=False)
197 self.col = col.astype(index_dtype, copy=False)
198 self.data = M[self.row, self.col]
199 self.has_canonical_format = True
201 if dtype is not None:
202 self.data = self.data.astype(dtype, copy=False)
204 self._check()
206 def reshape(self, *args, **kwargs):
207 shape = check_shape(args, self.shape)
208 order, copy = check_reshape_kwargs(kwargs)
210 # Return early if reshape is not required
211 if shape == self.shape:
212 if copy:
213 return self.copy()
214 else:
215 return self
217 nrows, ncols = self.shape
219 if order == 'C':
220 # Upcast to avoid overflows: the coo_array constructor
221 # below will downcast the results to a smaller dtype, if
222 # possible.
223 dtype = self._get_index_dtype(maxval=(ncols * max(0, nrows - 1) + max(0, ncols - 1)))
225 flat_indices = np.multiply(ncols, self.row, dtype=dtype) + self.col
226 new_row, new_col = divmod(flat_indices, shape[1])
227 elif order == 'F':
228 dtype = self._get_index_dtype(maxval=(nrows * max(0, ncols - 1) + max(0, nrows - 1)))
230 flat_indices = np.multiply(nrows, self.col, dtype=dtype) + self.row
231 new_col, new_row = divmod(flat_indices, shape[0])
232 else:
233 raise ValueError("'order' must be 'C' or 'F'")
235 # Handle copy here rather than passing on to the constructor so that no
236 # copy will be made of new_row and new_col regardless
237 if copy:
238 new_data = self.data.copy()
239 else:
240 new_data = self.data
242 return self.__class__((new_data, (new_row, new_col)),
243 shape=shape, copy=False)
245 reshape.__doc__ = _spbase.reshape.__doc__
247 def _getnnz(self, axis=None):
248 if axis is None:
249 nnz = len(self.data)
250 if nnz != len(self.row) or nnz != len(self.col):
251 raise ValueError('row, column, and data array must all be the '
252 'same length')
254 if self.data.ndim != 1 or self.row.ndim != 1 or \
255 self.col.ndim != 1:
256 raise ValueError('row, column, and data arrays must be 1-D')
258 return int(nnz)
260 if axis < 0:
261 axis += 2
262 if axis == 0:
263 return np.bincount(downcast_intp_index(self.col),
264 minlength=self.shape[1])
265 elif axis == 1:
266 return np.bincount(downcast_intp_index(self.row),
267 minlength=self.shape[0])
268 else:
269 raise ValueError('axis out of bounds')
271 _getnnz.__doc__ = _spbase._getnnz.__doc__
273 def _check(self):
274 """ Checks data structure for consistency """
276 # index arrays should have integer data types
277 if self.row.dtype.kind != 'i':
278 warn("row index array has non-integer dtype (%s) "
279 % self.row.dtype.name)
280 if self.col.dtype.kind != 'i':
281 warn("col index array has non-integer dtype (%s) "
282 % self.col.dtype.name)
284 idx_dtype = self._get_index_dtype((self.row, self.col), maxval=max(self.shape))
285 self.row = np.asarray(self.row, dtype=idx_dtype)
286 self.col = np.asarray(self.col, dtype=idx_dtype)
287 self.data = to_native(self.data)
289 if self.nnz > 0:
290 if self.row.max() >= self.shape[0]:
291 raise ValueError('row index exceeds matrix dimensions')
292 if self.col.max() >= self.shape[1]:
293 raise ValueError('column index exceeds matrix dimensions')
294 if self.row.min() < 0:
295 raise ValueError('negative row index found')
296 if self.col.min() < 0:
297 raise ValueError('negative column index found')
299 def transpose(self, axes=None, copy=False):
300 if axes is not None:
301 raise ValueError("Sparse matrices do not support "
302 "an 'axes' parameter because swapping "
303 "dimensions is the only logical permutation.")
305 M, N = self.shape
306 return self.__class__((self.data, (self.col, self.row)),
307 shape=(N, M), copy=copy)
309 transpose.__doc__ = _spbase.transpose.__doc__
311 def resize(self, *shape):
312 shape = check_shape(shape)
313 new_M, new_N = shape
314 M, N = self.shape
316 if new_M < M or new_N < N:
317 mask = np.logical_and(self.row < new_M, self.col < new_N)
318 if not mask.all():
319 self.row = self.row[mask]
320 self.col = self.col[mask]
321 self.data = self.data[mask]
323 self._shape = shape
325 resize.__doc__ = _spbase.resize.__doc__
327 def toarray(self, order=None, out=None):
328 """See the docstring for `_spbase.toarray`."""
329 B = self._process_toarray_args(order, out)
330 fortran = int(B.flags.f_contiguous)
331 if not fortran and not B.flags.c_contiguous:
332 raise ValueError("Output array must be C or F contiguous")
333 M,N = self.shape
334 coo_todense(M, N, self.nnz, self.row, self.col, self.data,
335 B.ravel('A'), fortran)
336 return B
338 def tocsc(self, copy=False):
339 """Convert this matrix to Compressed Sparse Column format
341 Duplicate entries will be summed together.
343 Examples
344 --------
345 >>> from numpy import array
346 >>> from scipy.sparse import coo_array
347 >>> row = array([0, 0, 1, 3, 1, 0, 0])
348 >>> col = array([0, 2, 1, 3, 1, 0, 0])
349 >>> data = array([1, 1, 1, 1, 1, 1, 1])
350 >>> A = coo_array((data, (row, col)), shape=(4, 4)).tocsc()
351 >>> A.toarray()
352 array([[3, 0, 1, 0],
353 [0, 2, 0, 0],
354 [0, 0, 0, 0],
355 [0, 0, 0, 1]])
357 """
358 if self.nnz == 0:
359 return self._csc_container(self.shape, dtype=self.dtype)
360 else:
361 M,N = self.shape
362 idx_dtype = self._get_index_dtype(
363 (self.col, self.row), maxval=max(self.nnz, M)
364 )
365 row = self.row.astype(idx_dtype, copy=False)
366 col = self.col.astype(idx_dtype, copy=False)
368 indptr = np.empty(N + 1, dtype=idx_dtype)
369 indices = np.empty_like(row, dtype=idx_dtype)
370 data = np.empty_like(self.data, dtype=upcast(self.dtype))
372 coo_tocsr(N, M, self.nnz, col, row, self.data,
373 indptr, indices, data)
375 x = self._csc_container((data, indices, indptr), shape=self.shape)
376 if not self.has_canonical_format:
377 x.sum_duplicates()
378 return x
380 def tocsr(self, copy=False):
381 """Convert this matrix to Compressed Sparse Row format
383 Duplicate entries will be summed together.
385 Examples
386 --------
387 >>> from numpy import array
388 >>> from scipy.sparse import coo_array
389 >>> row = array([0, 0, 1, 3, 1, 0, 0])
390 >>> col = array([0, 2, 1, 3, 1, 0, 0])
391 >>> data = array([1, 1, 1, 1, 1, 1, 1])
392 >>> A = coo_array((data, (row, col)), shape=(4, 4)).tocsr()
393 >>> A.toarray()
394 array([[3, 0, 1, 0],
395 [0, 2, 0, 0],
396 [0, 0, 0, 0],
397 [0, 0, 0, 1]])
399 """
400 if self.nnz == 0:
401 return self._csr_container(self.shape, dtype=self.dtype)
402 else:
403 M,N = self.shape
404 idx_dtype = self._get_index_dtype(
405 (self.row, self.col), maxval=max(self.nnz, N)
406 )
407 row = self.row.astype(idx_dtype, copy=False)
408 col = self.col.astype(idx_dtype, copy=False)
410 indptr = np.empty(M + 1, dtype=idx_dtype)
411 indices = np.empty_like(col, dtype=idx_dtype)
412 data = np.empty_like(self.data, dtype=upcast(self.dtype))
414 coo_tocsr(M, N, self.nnz, row, col, self.data,
415 indptr, indices, data)
417 x = self._csr_container((data, indices, indptr), shape=self.shape)
418 if not self.has_canonical_format:
419 x.sum_duplicates()
420 return x
422 def tocoo(self, copy=False):
423 if copy:
424 return self.copy()
425 else:
426 return self
428 tocoo.__doc__ = _spbase.tocoo.__doc__
430 def todia(self, copy=False):
431 self.sum_duplicates()
432 ks = self.col - self.row # the diagonal for each nonzero
433 diags, diag_idx = np.unique(ks, return_inverse=True)
435 if len(diags) > 100:
436 # probably undesired, should todia() have a maxdiags parameter?
437 warn("Constructing a DIA matrix with %d diagonals "
438 "is inefficient" % len(diags), SparseEfficiencyWarning)
440 #initialize and fill in data array
441 if self.data.size == 0:
442 data = np.zeros((0, 0), dtype=self.dtype)
443 else:
444 data = np.zeros((len(diags), self.col.max()+1), dtype=self.dtype)
445 data[diag_idx, self.col] = self.data
447 return self._dia_container((data, diags), shape=self.shape)
449 todia.__doc__ = _spbase.todia.__doc__
451 def todok(self, copy=False):
452 self.sum_duplicates()
453 dok = self._dok_container((self.shape), dtype=self.dtype)
454 dok._update(zip(zip(self.row,self.col),self.data))
456 return dok
458 todok.__doc__ = _spbase.todok.__doc__
460 def diagonal(self, k=0):
461 rows, cols = self.shape
462 if k <= -rows or k >= cols:
463 return np.empty(0, dtype=self.data.dtype)
464 diag = np.zeros(min(rows + min(k, 0), cols - max(k, 0)),
465 dtype=self.dtype)
466 diag_mask = (self.row + k) == self.col
468 if self.has_canonical_format:
469 row = self.row[diag_mask]
470 data = self.data[diag_mask]
471 else:
472 row, _, data = self._sum_duplicates(self.row[diag_mask],
473 self.col[diag_mask],
474 self.data[diag_mask])
475 diag[row + min(k, 0)] = data
477 return diag
479 diagonal.__doc__ = _data_matrix.diagonal.__doc__
481 def _setdiag(self, values, k):
482 M, N = self.shape
483 if values.ndim and not len(values):
484 return
485 idx_dtype = self.row.dtype
487 # Determine which triples to keep and where to put the new ones.
488 full_keep = self.col - self.row != k
489 if k < 0:
490 max_index = min(M+k, N)
491 if values.ndim:
492 max_index = min(max_index, len(values))
493 keep = np.logical_or(full_keep, self.col >= max_index)
494 new_row = np.arange(-k, -k + max_index, dtype=idx_dtype)
495 new_col = np.arange(max_index, dtype=idx_dtype)
496 else:
497 max_index = min(M, N-k)
498 if values.ndim:
499 max_index = min(max_index, len(values))
500 keep = np.logical_or(full_keep, self.row >= max_index)
501 new_row = np.arange(max_index, dtype=idx_dtype)
502 new_col = np.arange(k, k + max_index, dtype=idx_dtype)
504 # Define the array of data consisting of the entries to be added.
505 if values.ndim:
506 new_data = values[:max_index]
507 else:
508 new_data = np.empty(max_index, dtype=self.dtype)
509 new_data[:] = values
511 # Update the internal structure.
512 self.row = np.concatenate((self.row[keep], new_row))
513 self.col = np.concatenate((self.col[keep], new_col))
514 self.data = np.concatenate((self.data[keep], new_data))
515 self.has_canonical_format = False
517 # needed by _data_matrix
518 def _with_data(self,data,copy=True):
519 """Returns a matrix with the same sparsity structure as self,
520 but with different data. By default the index arrays
521 (i.e. .row and .col) are copied.
522 """
523 if copy:
524 return self.__class__((data, (self.row.copy(), self.col.copy())),
525 shape=self.shape, dtype=data.dtype)
526 else:
527 return self.__class__((data, (self.row, self.col)),
528 shape=self.shape, dtype=data.dtype)
530 def sum_duplicates(self):
531 """Eliminate duplicate matrix entries by adding them together
533 This is an *in place* operation
534 """
535 if self.has_canonical_format:
536 return
537 summed = self._sum_duplicates(self.row, self.col, self.data)
538 self.row, self.col, self.data = summed
539 self.has_canonical_format = True
541 def _sum_duplicates(self, row, col, data):
542 # Assumes (data, row, col) not in canonical format.
543 if len(data) == 0:
544 return row, col, data
545 # Sort indices w.r.t. rows, then cols. This corresponds to C-order,
546 # which we rely on for argmin/argmax to return the first index in the
547 # same way that numpy does (in the case of ties).
548 order = np.lexsort((col, row))
549 row = row[order]
550 col = col[order]
551 data = data[order]
552 unique_mask = ((row[1:] != row[:-1]) |
553 (col[1:] != col[:-1]))
554 unique_mask = np.append(True, unique_mask)
555 row = row[unique_mask]
556 col = col[unique_mask]
557 unique_inds, = np.nonzero(unique_mask)
558 data = np.add.reduceat(data, unique_inds, dtype=self.dtype)
559 return row, col, data
561 def eliminate_zeros(self):
562 """Remove zero entries from the matrix
564 This is an *in place* operation
565 """
566 mask = self.data != 0
567 self.data = self.data[mask]
568 self.row = self.row[mask]
569 self.col = self.col[mask]
571 #######################
572 # Arithmetic handlers #
573 #######################
575 def _add_dense(self, other):
576 if other.shape != self.shape:
577 raise ValueError('Incompatible shapes ({} and {})'
578 .format(self.shape, other.shape))
579 dtype = upcast_char(self.dtype.char, other.dtype.char)
580 result = np.array(other, dtype=dtype, copy=True)
581 fortran = int(result.flags.f_contiguous)
582 M, N = self.shape
583 coo_todense(M, N, self.nnz, self.row, self.col, self.data,
584 result.ravel('A'), fortran)
585 return self._container(result, copy=False)
587 def _mul_vector(self, other):
588 #output array
589 result = np.zeros(self.shape[0], dtype=upcast_char(self.dtype.char,
590 other.dtype.char))
591 coo_matvec(self.nnz, self.row, self.col, self.data, other, result)
592 return result
594 def _mul_multivector(self, other):
595 result = np.zeros((other.shape[1], self.shape[0]),
596 dtype=upcast_char(self.dtype.char, other.dtype.char))
597 for i, col in enumerate(other.T):
598 coo_matvec(self.nnz, self.row, self.col, self.data, col, result[i])
599 return result.T.view(type=type(other))
602def isspmatrix_coo(x):
603 """Is `x` of coo_matrix type?
605 Parameters
606 ----------
607 x
608 object to check for being a coo matrix
610 Returns
611 -------
612 bool
613 True if `x` is a coo matrix, False otherwise
615 Examples
616 --------
617 >>> from scipy.sparse import coo_array, coo_matrix, csr_matrix, isspmatrix_coo
618 >>> isspmatrix_coo(coo_matrix([[5]]))
619 True
620 >>> isspmatrix_coo(coo_array([[5]]))
621 False
622 >>> isspmatrix_coo(csr_matrix([[5]]))
623 False
624 """
625 return isinstance(x, coo_matrix)
628# This namespace class separates array from matrix with isinstance
629class coo_array(_coo_base, sparray):
630 pass
632coo_array.__doc__ = _coo_base.__doc__
634class coo_matrix(spmatrix, _coo_base):
635 pass
637coo_matrix.__doc__ = _array_doc_to_matrix(_coo_base.__doc__)