test_arithmetic1d.py 12 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341
  1. """Test of 1D arithmetic operations"""
  2. import pytest
  3. import numpy as np
  4. from numpy.testing import assert_equal, assert_allclose
  5. from scipy.sparse import coo_array, csr_array
  6. from scipy.sparse._sputils import isscalarlike
  7. spcreators = [coo_array, csr_array]
  8. math_dtypes = [np.int64, np.float64, np.complex128]
  9. def toarray(a):
  10. if isinstance(a, np.ndarray) or isscalarlike(a):
  11. return a
  12. return a.toarray()
  13. @pytest.fixture
  14. def dat1d():
  15. return np.array([3, 0, 1, 0], 'd')
  16. @pytest.fixture
  17. def datsp_math_dtypes(dat1d):
  18. dat_dtypes = {dtype: dat1d.astype(dtype) for dtype in math_dtypes}
  19. return {
  20. sp: [(dtype, dat, sp(dat)) for dtype, dat in dat_dtypes.items()]
  21. for sp in spcreators
  22. }
  23. @pytest.mark.parametrize("spcreator", spcreators)
  24. class TestArithmetic1D:
  25. def test_empty_arithmetic(self, spcreator):
  26. shape = (5,)
  27. for mytype in [
  28. np.dtype('int32'),
  29. np.dtype('float32'),
  30. np.dtype('float64'),
  31. np.dtype('complex64'),
  32. np.dtype('complex128'),
  33. ]:
  34. a = spcreator(shape, dtype=mytype)
  35. b = a + a
  36. c = 2 * a
  37. assert isinstance(a @ a.tocsr(), np.ndarray)
  38. assert isinstance(a @ a.tocoo(), np.ndarray)
  39. for m in [a, b, c]:
  40. assert m @ m == a.toarray() @ a.toarray()
  41. assert m.dtype == mytype
  42. assert toarray(m).dtype == mytype
  43. def test_abs(self, spcreator):
  44. A = np.array([-1, 0, 17, 0, -5, 0, 1, -4, 0, 0, 0, 0], 'd')
  45. assert_equal(abs(A), abs(spcreator(A)).toarray())
  46. def test_round(self, spcreator):
  47. A = np.array([-1.35, 0.56, 17.25, -5.98], 'd')
  48. Asp = spcreator(A)
  49. assert_equal(np.around(A, decimals=1), round(Asp, ndigits=1).toarray())
  50. def test_elementwise_power(self, spcreator):
  51. A = np.array([-4, -3, -2, -1, 0, 1, 2, 3, 4], 'd')
  52. Asp = spcreator(A)
  53. assert_equal(np.power(A, 2), Asp.power(2).toarray())
  54. # element-wise power function needs a scalar power
  55. with pytest.raises(NotImplementedError, match='input is not scalar'):
  56. spcreator(A).power(A)
  57. def test_real(self, spcreator):
  58. D = np.array([1 + 3j, 2 - 4j])
  59. A = spcreator(D)
  60. assert_equal(A.real.toarray(), D.real)
  61. def test_imag(self, spcreator):
  62. D = np.array([1 + 3j, 2 - 4j])
  63. A = spcreator(D)
  64. assert_equal(A.imag.toarray(), D.imag)
  65. def test_mul_scalar(self, spcreator, datsp_math_dtypes):
  66. for dtype, dat, datsp in datsp_math_dtypes[spcreator]:
  67. assert_equal(dat * 2, (datsp * 2).toarray())
  68. assert_equal(dat * 17.3, (datsp * 17.3).toarray())
  69. def test_rmul_scalar(self, spcreator, datsp_math_dtypes):
  70. for dtype, dat, datsp in datsp_math_dtypes[spcreator]:
  71. assert_equal(2 * dat, (2 * datsp).toarray())
  72. assert_equal(17.3 * dat, (17.3 * datsp).toarray())
  73. def test_sub(self, spcreator, datsp_math_dtypes):
  74. for dtype, dat, datsp in datsp_math_dtypes[spcreator]:
  75. if dtype == np.dtype('bool'):
  76. # boolean array subtraction deprecated in 1.9.0
  77. continue
  78. assert_equal((datsp - datsp).toarray(), np.zeros(4))
  79. assert_equal((datsp - 0).toarray(), dat)
  80. A = spcreator([1, -4, 0, 2], dtype='d')
  81. assert_equal((datsp - A).toarray(), dat - A.toarray())
  82. assert_equal((A - datsp).toarray(), A.toarray() - dat)
  83. # test broadcasting
  84. assert_equal(datsp.toarray() - dat[0], dat - dat[0])
  85. def test_add0(self, spcreator, datsp_math_dtypes):
  86. for dtype, dat, datsp in datsp_math_dtypes[spcreator]:
  87. # Adding 0 to a sparse matrix
  88. assert_equal((datsp + 0).toarray(), dat)
  89. # use sum (which takes 0 as a starting value)
  90. sumS = sum([k * datsp for k in range(1, 3)])
  91. sumD = sum([k * dat for k in range(1, 3)])
  92. assert_allclose(sumS.toarray(), sumD)
  93. def test_elementwise_multiply(self, spcreator):
  94. # real/real
  95. A = np.array([4, 0, 9])
  96. B = np.array([0, 7, -1])
  97. Asp = spcreator(A)
  98. Bsp = spcreator(B)
  99. assert_allclose(Asp.multiply(Bsp).toarray(), A * B) # sparse/sparse
  100. assert_allclose(Asp.multiply(B).toarray(), A * B) # sparse/dense
  101. # complex/complex
  102. C = np.array([1 - 2j, 0 + 5j, -1 + 0j])
  103. D = np.array([5 + 2j, 7 - 3j, -2 + 1j])
  104. Csp = spcreator(C)
  105. Dsp = spcreator(D)
  106. assert_allclose(Csp.multiply(Dsp).toarray(), C * D) # sparse/sparse
  107. assert_allclose(Csp.multiply(D).toarray(), C * D) # sparse/dense
  108. # real/complex
  109. assert_allclose(Asp.multiply(Dsp).toarray(), A * D) # sparse/sparse
  110. assert_allclose(Asp.multiply(D).toarray(), A * D) # sparse/dense
  111. def test_elementwise_multiply_broadcast(self, spcreator):
  112. A = np.array([4])
  113. B = np.array([[-9]])
  114. C = np.array([1, -1, 0])
  115. D = np.array([[7, 9, -9]])
  116. E = np.array([[3], [2], [1]])
  117. F = np.array([[8, 6, 3], [-4, 3, 2], [6, 6, 6]])
  118. G = [1, 2, 3]
  119. H = np.ones((3, 4))
  120. J = H.T
  121. K = np.array([[0]])
  122. L = np.array([[[1, 2], [0, 1]]])
  123. # Some arrays can't be cast as spmatrices (A, C, L) so leave
  124. # them out.
  125. Asp = spcreator(A)
  126. Csp = spcreator(C)
  127. Gsp = spcreator(G)
  128. # 2d arrays
  129. Bsp = spcreator(B)
  130. Dsp = spcreator(D)
  131. Esp = spcreator(E)
  132. Fsp = spcreator(F)
  133. Hsp = spcreator(H)
  134. Hspp = spcreator(H[0, None])
  135. Jsp = spcreator(J)
  136. Jspp = spcreator(J[:, 0, None])
  137. Ksp = spcreator(K)
  138. matrices = [A, B, C, D, E, F, G, H, J, K, L]
  139. spmatrices = [Asp, Bsp, Csp, Dsp, Esp, Fsp, Gsp, Hsp, Hspp, Jsp, Jspp, Ksp]
  140. sp1dmatrices = [Asp, Csp, Gsp]
  141. # sparse/sparse
  142. for i in sp1dmatrices:
  143. for j in spmatrices:
  144. try:
  145. dense_mult = i.toarray() * j.toarray()
  146. except ValueError:
  147. with pytest.raises(ValueError, match='inconsistent shapes'):
  148. i.multiply(j)
  149. continue
  150. sp_mult = i.multiply(j)
  151. assert_allclose(sp_mult.toarray(), dense_mult)
  152. # sparse/dense
  153. for i in sp1dmatrices:
  154. for j in matrices:
  155. try:
  156. dense_mult = i.toarray() * j
  157. except TypeError:
  158. continue
  159. except ValueError:
  160. matchme = 'broadcast together|inconsistent shapes'
  161. with pytest.raises(ValueError, match=matchme):
  162. i.multiply(j)
  163. continue
  164. try:
  165. sp_mult = i.multiply(j)
  166. except ValueError:
  167. continue
  168. assert_allclose(toarray(sp_mult), dense_mult)
  169. def test_elementwise_divide(self, spcreator, dat1d):
  170. datsp = spcreator(dat1d)
  171. expected = np.array([1, np.nan, 1, np.nan])
  172. actual = datsp / datsp
  173. # need assert_array_equal to handle nan values
  174. np.testing.assert_array_equal(actual, expected)
  175. denom = spcreator([1, 0, 0, 4], dtype='d')
  176. expected = [3, np.nan, np.inf, 0]
  177. np.testing.assert_array_equal(datsp / denom, expected)
  178. # complex
  179. A = np.array([1 - 2j, 0 + 5j, -1 + 0j])
  180. B = np.array([5 + 2j, 7 - 3j, -2 + 1j])
  181. Asp = spcreator(A)
  182. Bsp = spcreator(B)
  183. assert_allclose(Asp / Bsp, A / B)
  184. # integer
  185. A = np.array([1, 2, 3])
  186. B = np.array([0, 1, 2])
  187. Asp = spcreator(A)
  188. Bsp = spcreator(B)
  189. with np.errstate(divide='ignore'):
  190. assert_equal(Asp / Bsp, A / B)
  191. # mismatching sparsity patterns
  192. A = np.array([0, 1])
  193. B = np.array([1, 0])
  194. Asp = spcreator(A)
  195. Bsp = spcreator(B)
  196. with np.errstate(divide='ignore', invalid='ignore'):
  197. assert_equal(Asp / Bsp, A / B)
  198. def test_pow(self, spcreator):
  199. A = np.array([1, 0, 2, 0])
  200. B = spcreator(A)
  201. # unusual exponents
  202. with pytest.raises(ValueError, match='negative integer powers'):
  203. B**-1
  204. with pytest.raises(NotImplementedError, match='zero power'):
  205. B**0
  206. for exponent in [1, 2, 3, 2.2]:
  207. ret_sp = B**exponent
  208. ret_np = A**exponent
  209. assert_equal(ret_sp.toarray(), ret_np)
  210. assert_equal(ret_sp.dtype, ret_np.dtype)
  211. def test_dot_scalar(self, spcreator, dat1d):
  212. A = spcreator(dat1d)
  213. scalar = 10
  214. actual = A.dot(scalar)
  215. expected = A * scalar
  216. assert_allclose(actual.toarray(), expected.toarray())
  217. def test_matmul(self, spcreator):
  218. Msp = spcreator([2, 0, 3.0])
  219. B = spcreator(np.array([[0, 1], [1, 0], [0, 2]], 'd'))
  220. col = np.array([[1, 2, 3]]).T
  221. # check sparse @ dense 2d column
  222. assert_allclose(Msp @ col, Msp.toarray() @ col)
  223. # check sparse1d @ sparse2d, sparse1d @ dense2d, dense1d @ sparse2d
  224. assert_allclose((Msp @ B).toarray(), (Msp @ B).toarray())
  225. assert_allclose(Msp.toarray() @ B, (Msp @ B).toarray())
  226. assert_allclose(Msp @ B.toarray(), (Msp @ B).toarray())
  227. # check sparse1d @ dense1d, sparse1d @ sparse1d
  228. V = np.array([0, 0, 1])
  229. assert_allclose(Msp @ V, Msp.toarray() @ V)
  230. Vsp = spcreator(V)
  231. Msp_Vsp = Msp @ Vsp
  232. assert isinstance(Msp_Vsp, np.ndarray)
  233. assert Msp_Vsp.shape == ()
  234. # output is 0-dim ndarray
  235. assert_allclose(np.array(3), Msp_Vsp)
  236. assert_allclose(np.array(3), Msp.toarray() @ Vsp)
  237. assert_allclose(np.array(3), Msp @ Vsp.toarray())
  238. assert_allclose(np.array(3), Msp.toarray() @ Vsp.toarray())
  239. # check error on matrix-scalar
  240. with pytest.raises(ValueError, match='Scalar operands are not allowed'):
  241. Msp @ 1
  242. with pytest.raises(ValueError, match='Scalar operands are not allowed'):
  243. 1 @ Msp
  244. def test_sub_dense(self, spcreator, datsp_math_dtypes):
  245. # subtracting a dense matrix to/from a sparse matrix
  246. for dtype, dat, datsp in datsp_math_dtypes[spcreator]:
  247. if dtype == np.dtype('bool'):
  248. # boolean array subtraction deprecated in 1.9.0
  249. continue
  250. # Manually add to avoid upcasting from scalar
  251. # multiplication.
  252. sum1 = (dat + dat + dat) - datsp
  253. assert_equal(sum1, dat + dat)
  254. sum2 = (datsp + datsp + datsp) - dat
  255. assert_equal(sum2, dat + dat)
  256. def test_size_zero_matrix_arithmetic(self, spcreator):
  257. # Test basic matrix arithmetic with shapes like 0, (1, 0), (0, 3), etc.
  258. mat = np.array([])
  259. a = mat.reshape(0)
  260. d = mat.reshape((1, 0))
  261. f = np.ones([5, 5])
  262. asp = spcreator(a)
  263. dsp = spcreator(d)
  264. # bad shape for addition
  265. with pytest.raises(ValueError, match='inconsistent shapes'):
  266. asp.__add__(dsp)
  267. # matrix product.
  268. assert_equal(asp.dot(asp), np.dot(a, a))
  269. # bad matrix products
  270. with pytest.raises(ValueError, match='dimension mismatch|shapes.*not aligned'):
  271. asp.dot(f)
  272. # elemente-wise multiplication
  273. assert_equal(asp.multiply(asp).toarray(), np.multiply(a, a))
  274. assert_equal(asp.multiply(a).toarray(), np.multiply(a, a))
  275. assert_equal(asp.multiply(6).toarray(), np.multiply(a, 6))
  276. # bad element-wise multiplication
  277. with pytest.raises(ValueError, match='inconsistent shapes'):
  278. asp.multiply(f)
  279. # Addition
  280. assert_equal(asp.__add__(asp).toarray(), a.__add__(a))