Coverage for /pythoncovmergedfiles/medio/medio/usr/local/lib/python3.8/site-packages/scipy/sparse/linalg/

1"""Abstract linear algebra library.

3This module defines a class hierarchy that implements a kind of "lazy"

4matrix representation, called the ``LinearOperator``. It can be used to do

5linear algebra with extremely large sparse or structured matrices, without

6representing those explicitly in memory. Such matrices can be added,

7multiplied, transposed, etc.

9As a motivating example, suppose you want have a matrix where almost all of

10the elements have the value one. The standard sparse matrix representation

11skips the storage of zeros, but not ones. By contrast, a LinearOperator is

12able to represent such matrices efficiently. First, we need a compact way to

13represent an all-ones matrix::

15 >>> import numpy as np

16 >>> class Ones(LinearOperator):

17 ... def __init__(self, shape):

18 ... super().__init__(dtype=None, shape=shape)

19 ... def _matvec(self, x):

20 ... return np.repeat(x.sum(), self.shape[0])

22Instances of this class emulate ``np.ones(shape)``, but using a constant

23amount of storage, independent of ``shape``. The ``_matvec`` method specifies

24how this linear operator multiplies with (operates on) a vector. We can now

25add this operator to a sparse matrix that stores only offsets from one::

27 >>> from scipy.sparse import csr_matrix

28 >>> offsets = csr_matrix([[1, 0, 2], [0, -1, 0], [0, 0, 3]])

29 >>> A = aslinearoperator(offsets) + Ones(offsets.shape)

30 >>> A.dot([1, 2, 3])

31 array([13, 4, 15])

33The result is the same as that given by its dense, explicitly-stored

34counterpart::

36 >>> (np.ones(A.shape, A.dtype) + offsets.toarray()).dot([1, 2, 3])

37 array([13, 4, 15])

39Several algorithms in the ``scipy.sparse`` library are able to operate on

40``LinearOperator`` instances.

41"""

43import warnings

45import numpy as np

47from scipy.sparse import isspmatrix

48from scipy.sparse._sputils import isshape, isintlike, asmatrix, is_pydata_spmatrix

50__all__ = ['LinearOperator', 'aslinearoperator']

53class LinearOperator:

54 """Common interface for performing matrix vector products

56 Many iterative methods (e.g. cg, gmres) do not need to know the

57 individual entries of a matrix to solve a linear system A*x=b.

58 Such solvers only require the computation of matrix vector

59 products, A*v where v is a dense vector. This class serves as

60 an abstract interface between iterative solvers and matrix-like

61 objects.

63 To construct a concrete LinearOperator, either pass appropriate

64 callables to the constructor of this class, or subclass it.

66 A subclass must implement either one of the methods ``_matvec``

67 and ``_matmat``, and the attributes/properties ``shape`` (pair of

68 integers) and ``dtype`` (may be None). It may call the ``__init__``

69 on this class to have these attributes validated. Implementing

70 ``_matvec`` automatically implements ``_matmat`` (using a naive

71 algorithm) and vice-versa.

73 Optionally, a subclass may implement ``_rmatvec`` or ``_adjoint``

74 to implement the Hermitian adjoint (conjugate transpose). As with

75 ``_matvec`` and ``_matmat``, implementing either ``_rmatvec`` or

76 ``_adjoint`` implements the other automatically. Implementing

77 ``_adjoint`` is preferable; ``_rmatvec`` is mostly there for

78 backwards compatibility.

80 Parameters

81 ----------

82 shape : tuple

83 Matrix dimensions (M, N).

84 matvec : callable f(v)

85 Returns returns A * v.

86 rmatvec : callable f(v)

87 Returns A^H * v, where A^H is the conjugate transpose of A.

88 matmat : callable f(V)

89 Returns A * V, where V is a dense matrix with dimensions (N, K).

90 dtype : dtype

91 Data type of the matrix.

92 rmatmat : callable f(V)

93 Returns A^H * V, where V is a dense matrix with dimensions (M, K).

95 Attributes

96 ----------

97 args : tuple

98 For linear operators describing products etc. of other linear

99 operators, the operands of the binary operation.

100 ndim : int

101 Number of dimensions (this is always 2)

102

103 See Also

104 --------

105 aslinearoperator : Construct LinearOperators

106

107 Notes

108 -----

109 The user-defined matvec() function must properly handle the case

110 where v has shape (N,) as well as the (N,1) case. The shape of

111 the return type is handled internally by LinearOperator.

112

113 LinearOperator instances can also be multiplied, added with each

114 other and exponentiated, all lazily: the result of these operations

115 is always a new, composite LinearOperator, that defers linear

116 operations to the original operators and combines the results.

117

118 More details regarding how to subclass a LinearOperator and several

119 examples of concrete LinearOperator instances can be found in the

120 external project `PyLops <https://pylops.readthedocs.io>`_.

121

122

123 Examples

124 --------

125 >>> import numpy as np

126 >>> from scipy.sparse.linalg import LinearOperator

127 >>> def mv(v):

128 ... return np.array([2*v[0], 3*v[1]])

129 ...

130 >>> A = LinearOperator((2,2), matvec=mv)

131 >>> A

132 <2x2 _CustomLinearOperator with dtype=float64>

133 >>> A.matvec(np.ones(2))

134 array([ 2., 3.])

135 >>> A * np.ones(2)

136 array([ 2., 3.])

137

138 """

139

140 ndim = 2

141

142 def __new__(cls, *args, **kwargs):

143 if cls is LinearOperator:

144 # Operate as _CustomLinearOperator factory.

145 return super(LinearOperator, cls).__new__(_CustomLinearOperator)

146 else:

147 obj = super(LinearOperator, cls).__new__(cls)

148

149 if (type(obj)._matvec == LinearOperator._matvec

150 and type(obj)._matmat == LinearOperator._matmat):

151 warnings.warn("LinearOperator subclass should implement"

152 " at least one of _matvec and _matmat.",

153 category=RuntimeWarning, stacklevel=2)

154

155 return obj

156

157 def __init__(self, dtype, shape):

158 """Initialize this LinearOperator.

159

160 To be called by subclasses. ``dtype`` may be None; ``shape`` should

161 be convertible to a length-2 tuple.

162 """

163 if dtype is not None:

164 dtype = np.dtype(dtype)

165

166 shape = tuple(shape)

167 if not isshape(shape):

168 raise ValueError("invalid shape %r (must be 2-d)" % (shape,))

169

170 self.dtype = dtype

171 self.shape = shape

172

173 def _init_dtype(self):

174 """Called from subclasses at the end of the __init__ routine.

175 """

176 if self.dtype is None:

177 v = np.zeros(self.shape[-1])

178 self.dtype = np.asarray(self.matvec(v)).dtype

179

180 def _matmat(self, X):

181 """Default matrix-matrix multiplication handler.

182

183 Falls back on the user-defined _matvec method, so defining that will

184 define matrix multiplication (though in a very suboptimal way).

185 """

186

187 return np.hstack([self.matvec(col.reshape(-1,1)) for col in X.T])

188

189 def _matvec(self, x):

190 """Default matrix-vector multiplication handler.

191

192 If self is a linear operator of shape (M, N), then this method will

193 be called on a shape (N,) or (N, 1) ndarray, and should return a

194 shape (M,) or (M, 1) ndarray.

195

196 This default implementation falls back on _matmat, so defining that

197 will define matrix-vector multiplication as well.

198 """

199 return self.matmat(x.reshape(-1, 1))

200

201 def matvec(self, x):

202 """Matrix-vector multiplication.

203

204 Performs the operation y=A*x where A is an MxN linear

205 operator and x is a column vector or 1-d array.

206

207 Parameters

208 ----------

209 x : {matrix, ndarray}

210 An array with shape (N,) or (N,1).

211

212 Returns

213 -------

214 y : {matrix, ndarray}

215 A matrix or ndarray with shape (M,) or (M,1) depending

216 on the type and shape of the x argument.

217

218 Notes

219 -----

220 This matvec wraps the user-specified matvec routine or overridden

221 _matvec method to ensure that y has the correct shape and type.

222

223 """

224

225 x = np.asanyarray(x)

226

227 M,N = self.shape

228

229 if x.shape != (N,) and x.shape != (N,1):

230 raise ValueError('dimension mismatch')

231

232 y = self._matvec(x)

233

234 if isinstance(x, np.matrix):

235 y = asmatrix(y)

236 else:

237 y = np.asarray(y)

238

239 if x.ndim == 1:

240 y = y.reshape(M)

241 elif x.ndim == 2:

242 y = y.reshape(M,1)

243 else:

244 raise ValueError('invalid shape returned by user-defined matvec()')

245

246 return y

247

248 def rmatvec(self, x):

249 """Adjoint matrix-vector multiplication.

250

251 Performs the operation y = A^H * x where A is an MxN linear

252 operator and x is a column vector or 1-d array.

253

254 Parameters

255 ----------

256 x : {matrix, ndarray}

257 An array with shape (M,) or (M,1).

258

259 Returns

260 -------

261 y : {matrix, ndarray}

262 A matrix or ndarray with shape (N,) or (N,1) depending

263 on the type and shape of the x argument.

264

265 Notes

266 -----

267 This rmatvec wraps the user-specified rmatvec routine or overridden

268 _rmatvec method to ensure that y has the correct shape and type.

269

270 """

271

272 x = np.asanyarray(x)

273

274 M,N = self.shape

275

276 if x.shape != (M,) and x.shape != (M,1):

277 raise ValueError('dimension mismatch')

278

279 y = self._rmatvec(x)

280

281 if isinstance(x, np.matrix):

282 y = asmatrix(y)

283 else:

284 y = np.asarray(y)

285

286 if x.ndim == 1:

287 y = y.reshape(N)

288 elif x.ndim == 2:

289 y = y.reshape(N,1)

290 else:

291 raise ValueError('invalid shape returned by user-defined rmatvec()')

292

293 return y

294

295 def _rmatvec(self, x):

296 """Default implementation of _rmatvec; defers to adjoint."""

297 if type(self)._adjoint == LinearOperator._adjoint:

298 # _adjoint not overridden, prevent infinite recursion

299 raise NotImplementedError

300 else:

301 return self.H.matvec(x)

302

303 def matmat(self, X):

304 """Matrix-matrix multiplication.

305

306 Performs the operation y=A*X where A is an MxN linear

307 operator and X dense N*K matrix or ndarray.

308

309 Parameters

310 ----------

311 X : {matrix, ndarray}

312 An array with shape (N,K).

313

314 Returns

315 -------

316 Y : {matrix, ndarray}

317 A matrix or ndarray with shape (M,K) depending on

318 the type of the X argument.

319

320 Notes

321 -----

322 This matmat wraps any user-specified matmat routine or overridden

323 _matmat method to ensure that y has the correct type.

324

325 """

326

327 X = np.asanyarray(X)

328

329 if X.ndim != 2:

330 raise ValueError('expected 2-d ndarray or matrix, not %d-d'

331 % X.ndim)

332

333 if X.shape[0] != self.shape[1]:

334 raise ValueError('dimension mismatch: %r, %r'

335 % (self.shape, X.shape))

336

337 Y = self._matmat(X)

338

339 if isinstance(Y, np.matrix):

340 Y = asmatrix(Y)

341

342 return Y

343

344 def rmatmat(self, X):

345 """Adjoint matrix-matrix multiplication.

346

347 Performs the operation y = A^H * x where A is an MxN linear

348 operator and x is a column vector or 1-d array, or 2-d array.

349 The default implementation defers to the adjoint.

350

351 Parameters

352 ----------

353 X : {matrix, ndarray}

354 A matrix or 2D array.

355

356 Returns

357 -------

358 Y : {matrix, ndarray}

359 A matrix or 2D array depending on the type of the input.

360

361 Notes

362 -----

363 This rmatmat wraps the user-specified rmatmat routine.

364

365 """

366

367 X = np.asanyarray(X)

368

369 if X.ndim != 2:

370 raise ValueError('expected 2-d ndarray or matrix, not %d-d'

371 % X.ndim)

372

373 if X.shape[0] != self.shape[0]:

374 raise ValueError('dimension mismatch: %r, %r'

375 % (self.shape, X.shape))

376

377 Y = self._rmatmat(X)

378 if isinstance(Y, np.matrix):

379 Y = asmatrix(Y)

380 return Y

381

382 def _rmatmat(self, X):

383 """Default implementation of _rmatmat defers to rmatvec or adjoint."""

384 if type(self)._adjoint == LinearOperator._adjoint:

385 return np.hstack([self.rmatvec(col.reshape(-1, 1)) for col in X.T])

386 else:

387 return self.H.matmat(X)

388

389 def __call__(self, x):

390 return self*x

391

392 def __mul__(self, x):

393 return self.dot(x)

394

395 def dot(self, x):

396 """Matrix-matrix or matrix-vector multiplication.

397

398 Parameters

399 ----------

400 x : array_like

401 1-d or 2-d array, representing a vector or matrix.

402

403 Returns

404 -------

405 Ax : array

406 1-d or 2-d array (depending on the shape of x) that represents

407 the result of applying this linear operator on x.

408

409 """

410 if isinstance(x, LinearOperator):

411 return _ProductLinearOperator(self, x)

412 elif np.isscalar(x):

413 return _ScaledLinearOperator(self, x)

414 else:

415 x = np.asarray(x)

416

417 if x.ndim == 1 or x.ndim == 2 and x.shape[1] == 1:

418 return self.matvec(x)

419 elif x.ndim == 2:

420 return self.matmat(x)

421 else:

422 raise ValueError('expected 1-d or 2-d array or matrix, got %r'

423 % x)

424

425 def __matmul__(self, other):

426 if np.isscalar(other):

427 raise ValueError("Scalar operands are not allowed, "

428 "use '*' instead")

429 return self.__mul__(other)

430

431 def __rmatmul__(self, other):

432 if np.isscalar(other):

433 raise ValueError("Scalar operands are not allowed, "

434 "use '*' instead")

435 return self.__rmul__(other)

436

437 def __rmul__(self, x):

438 if np.isscalar(x):

439 return _ScaledLinearOperator(self, x)

440 else:

441 return NotImplemented

442

443 def __pow__(self, p):

444 if np.isscalar(p):

445 return _PowerLinearOperator(self, p)

446 else:

447 return NotImplemented

448

449 def __add__(self, x):

450 if isinstance(x, LinearOperator):

451 return _SumLinearOperator(self, x)

452 else:

453 return NotImplemented

454

455 def __neg__(self):

456 return _ScaledLinearOperator(self, -1)

457

458 def __sub__(self, x):

459 return self.__add__(-x)

460

461 def __repr__(self):

462 M,N = self.shape

463 if self.dtype is None:

464 dt = 'unspecified dtype'

465 else:

466 dt = 'dtype=' + str(self.dtype)

467

468 return '<%dx%d %s with %s>' % (M, N, self.__class__.__name__, dt)

469

470 def adjoint(self):

471 """Hermitian adjoint.

472

473 Returns the Hermitian adjoint of self, aka the Hermitian

474 conjugate or Hermitian transpose. For a complex matrix, the

475 Hermitian adjoint is equal to the conjugate transpose.

476

477 Can be abbreviated self.H instead of self.adjoint().

478

479 Returns

480 -------

481 A_H : LinearOperator

482 Hermitian adjoint of self.

483 """

484 return self._adjoint()

485

486 H = property(adjoint)

487

488 def transpose(self):

489 """Transpose this linear operator.

490

491 Returns a LinearOperator that represents the transpose of this one.

492 Can be abbreviated self.T instead of self.transpose().

493 """

494 return self._transpose()

495

496 T = property(transpose)

497

498 def _adjoint(self):

499 """Default implementation of _adjoint; defers to rmatvec."""

500 return _AdjointLinearOperator(self)

501

502 def _transpose(self):

503 """ Default implementation of _transpose; defers to rmatvec + conj"""

504 return _TransposedLinearOperator(self)

505

506

507class _CustomLinearOperator(LinearOperator):

508 """Linear operator defined in terms of user-specified operations."""

509

510 def __init__(self, shape, matvec, rmatvec=None, matmat=None,

511 dtype=None, rmatmat=None):

512 super().__init__(dtype, shape)

513

514 self.args = ()

515

516 self.__matvec_impl = matvec

517 self.__rmatvec_impl = rmatvec

518 self.__rmatmat_impl = rmatmat

519 self.__matmat_impl = matmat

520

521 self._init_dtype()

522

523 def _matmat(self, X):

524 if self.__matmat_impl is not None:

525 return self.__matmat_impl(X)

526 else:

527 return super()._matmat(X)

528

529 def _matvec(self, x):

530 return self.__matvec_impl(x)

531

532 def _rmatvec(self, x):

533 func = self.__rmatvec_impl

534 if func is None:

535 raise NotImplementedError("rmatvec is not defined")

536 return self.__rmatvec_impl(x)

537

538 def _rmatmat(self, X):

539 if self.__rmatmat_impl is not None:

540 return self.__rmatmat_impl(X)

541 else:

542 return super()._rmatmat(X)

543

544 def _adjoint(self):

545 return _CustomLinearOperator(shape=(self.shape[1], self.shape[0]),

546 matvec=self.__rmatvec_impl,

547 rmatvec=self.__matvec_impl,

548 matmat=self.__rmatmat_impl,

549 rmatmat=self.__matmat_impl,

550 dtype=self.dtype)

551

552

553class _AdjointLinearOperator(LinearOperator):

554 """Adjoint of arbitrary Linear Operator"""

555

556 def __init__(self, A):

557 shape = (A.shape[1], A.shape[0])

558 super().__init__(dtype=A.dtype, shape=shape)

559 self.A = A

560 self.args = (A,)

561

562 def _matvec(self, x):

563 return self.A._rmatvec(x)

564

565 def _rmatvec(self, x):

566 return self.A._matvec(x)

567

568 def _matmat(self, x):

569 return self.A._rmatmat(x)

570

571 def _rmatmat(self, x):

572 return self.A._matmat(x)

573

574class _TransposedLinearOperator(LinearOperator):

575 """Transposition of arbitrary Linear Operator"""

576

577 def __init__(self, A):

578 shape = (A.shape[1], A.shape[0])

579 super().__init__(dtype=A.dtype, shape=shape)

580 self.A = A

581 self.args = (A,)

582

583 def _matvec(self, x):

584 # NB. np.conj works also on sparse matrices

585 return np.conj(self.A._rmatvec(np.conj(x)))

586

587 def _rmatvec(self, x):

588 return np.conj(self.A._matvec(np.conj(x)))

589

590 def _matmat(self, x):

591 # NB. np.conj works also on sparse matrices

592 return np.conj(self.A._rmatmat(np.conj(x)))

593

594 def _rmatmat(self, x):

595 return np.conj(self.A._matmat(np.conj(x)))

596

597def _get_dtype(operators, dtypes=None):

598 if dtypes is None:

599 dtypes = []

600 for obj in operators:

601 if obj is not None and hasattr(obj, 'dtype'):

602 dtypes.append(obj.dtype)

603 return np.result_type(*dtypes)

604

605

606class _SumLinearOperator(LinearOperator):

607 def __init__(self, A, B):

608 if not isinstance(A, LinearOperator) or \

609 not isinstance(B, LinearOperator):

610 raise ValueError('both operands have to be a LinearOperator')

611 if A.shape != B.shape:

612 raise ValueError('cannot add %r and %r: shape mismatch'

613 % (A, B))

614 self.args = (A, B)

615 super().__init__(_get_dtype([A, B]), A.shape)

616

617 def _matvec(self, x):

618 return self.args[0].matvec(x) + self.args[1].matvec(x)

619

620 def _rmatvec(self, x):

621 return self.args[0].rmatvec(x) + self.args[1].rmatvec(x)

622

623 def _rmatmat(self, x):

624 return self.args[0].rmatmat(x) + self.args[1].rmatmat(x)

625

626 def _matmat(self, x):

627 return self.args[0].matmat(x) + self.args[1].matmat(x)

628

629 def _adjoint(self):

630 A, B = self.args

631 return A.H + B.H

632

633

634class _ProductLinearOperator(LinearOperator):

635 def __init__(self, A, B):

636 if not isinstance(A, LinearOperator) or \

637 not isinstance(B, LinearOperator):

638 raise ValueError('both operands have to be a LinearOperator')

639 if A.shape[1] != B.shape[0]:

640 raise ValueError('cannot multiply %r and %r: shape mismatch'

641 % (A, B))

642 super().__init__(_get_dtype([A, B]),

643 (A.shape[0], B.shape[1]))

644 self.args = (A, B)

645

646 def _matvec(self, x):

647 return self.args[0].matvec(self.args[1].matvec(x))

648

649 def _rmatvec(self, x):

650 return self.args[1].rmatvec(self.args[0].rmatvec(x))

651

652 def _rmatmat(self, x):

653 return self.args[1].rmatmat(self.args[0].rmatmat(x))

654

655 def _matmat(self, x):

656 return self.args[0].matmat(self.args[1].matmat(x))

657

658 def _adjoint(self):

659 A, B = self.args

660 return B.H * A.H

661

662

663class _ScaledLinearOperator(LinearOperator):

664 def __init__(self, A, alpha):

665 if not isinstance(A, LinearOperator):

666 raise ValueError('LinearOperator expected as A')

667 if not np.isscalar(alpha):

668 raise ValueError('scalar expected as alpha')

669 dtype = _get_dtype([A], [type(alpha)])

670 super().__init__(dtype, A.shape)

671 self.args = (A, alpha)

672

673 def _matvec(self, x):

674 return self.args[1] * self.args[0].matvec(x)

675

676 def _rmatvec(self, x):

677 return np.conj(self.args[1]) * self.args[0].rmatvec(x)

678

679 def _rmatmat(self, x):

680 return np.conj(self.args[1]) * self.args[0].rmatmat(x)

681

682 def _matmat(self, x):

683 return self.args[1] * self.args[0].matmat(x)

684

685 def _adjoint(self):

686 A, alpha = self.args

687 return A.H * np.conj(alpha)

688

689

690class _PowerLinearOperator(LinearOperator):

691 def __init__(self, A, p):

692 if not isinstance(A, LinearOperator):

693 raise ValueError('LinearOperator expected as A')

694 if A.shape[0] != A.shape[1]:

695 raise ValueError('square LinearOperator expected, got %r' % A)

696 if not isintlike(p) or p < 0:

697 raise ValueError('non-negative integer expected as p')

698

699 super().__init__(_get_dtype([A]), A.shape)

700 self.args = (A, p)

701

702 def _power(self, fun, x):

703 res = np.array(x, copy=True)

704 for i in range(self.args[1]):

705 res = fun(res)

706 return res

707

708 def _matvec(self, x):

709 return self._power(self.args[0].matvec, x)

710

711 def _rmatvec(self, x):

712 return self._power(self.args[0].rmatvec, x)

713

714 def _rmatmat(self, x):

715 return self._power(self.args[0].rmatmat, x)

716

717 def _matmat(self, x):

718 return self._power(self.args[0].matmat, x)

719

720 def _adjoint(self):

721 A, p = self.args

722 return A.H ** p

723

724

725class MatrixLinearOperator(LinearOperator):

726 def __init__(self, A):

727 super().__init__(A.dtype, A.shape)

728 self.A = A

729 self.__adj = None

730 self.args = (A,)

731

732 def _matmat(self, X):

733 return self.A.dot(X)

734

735 def _adjoint(self):

736 if self.__adj is None:

737 self.__adj = _AdjointMatrixOperator(self)

738 return self.__adj

739

740class _AdjointMatrixOperator(MatrixLinearOperator):

741 def __init__(self, adjoint):

742 self.A = adjoint.A.T.conj()

743 self.__adjoint = adjoint

744 self.args = (adjoint,)

745 self.shape = adjoint.shape[1], adjoint.shape[0]

746

747 @property

748 def dtype(self):

749 return self.__adjoint.dtype

750

751 def _adjoint(self):

752 return self.__adjoint

753

754

755class IdentityOperator(LinearOperator):

756 def __init__(self, shape, dtype=None):

757 super().__init__(dtype, shape)

758

759 def _matvec(self, x):

760 return x

761

762 def _rmatvec(self, x):

763 return x

764

765 def _rmatmat(self, x):

766 return x

767

768 def _matmat(self, x):

769 return x

770

771 def _adjoint(self):

772 return self

773

774

775def aslinearoperator(A):

776 """Return A as a LinearOperator.

777

778 'A' may be any of the following types:

779 - ndarray

780 - matrix

781 - sparse matrix (e.g. csr_matrix, lil_matrix, etc.)

782 - LinearOperator

783 - An object with .shape and .matvec attributes

784

785 See the LinearOperator documentation for additional information.

786

787 Notes

788 -----

789 If 'A' has no .dtype attribute, the data type is determined by calling

790 :func:`LinearOperator.matvec()` - set the .dtype attribute to prevent this

791 call upon the linear operator creation.

792

793 Examples

794 --------

795 >>> import numpy as np

796 >>> from scipy.sparse.linalg import aslinearoperator

797 >>> M = np.array([[1,2,3],[4,5,6]], dtype=np.int32)

798 >>> aslinearoperator(M)

799 <2x3 MatrixLinearOperator with dtype=int32>

800 """

801 if isinstance(A, LinearOperator):

802 return A

803

804 elif isinstance(A, np.ndarray) or isinstance(A, np.matrix):

805 if A.ndim > 2:

806 raise ValueError('array must have ndim <= 2')

807 A = np.atleast_2d(np.asarray(A))

808 return MatrixLinearOperator(A)

809

810 elif isspmatrix(A) or is_pydata_spmatrix(A):

811 return MatrixLinearOperator(A)

812

813 else:

814 if hasattr(A, 'shape') and hasattr(A, 'matvec'):

815 rmatvec = None

816 rmatmat = None

817 dtype = None

818

819 if hasattr(A, 'rmatvec'):

820 rmatvec = A.rmatvec

821 if hasattr(A, 'rmatmat'):

822 rmatmat = A.rmatmat

823 if hasattr(A, 'dtype'):

824 dtype = A.dtype

825 return LinearOperator(A.shape, A.matvec, rmatvec=rmatvec,

826 rmatmat=rmatmat, dtype=dtype)

827

828 else:

829 raise TypeError('type not understood')

Coverage for /pythoncovmergedfiles/medio/medio/usr/local/lib/python3.8/site-packages/scipy/sparse/linalg/_interface.py: 29%

346 statements