diff options
author | Pauli Virtanen <pav@iki.fi> | 2013-10-07 20:51:31 +0300 |
---|---|---|
committer | Pauli Virtanen <pav@iki.fi> | 2013-10-07 23:24:20 +0300 |
commit | f2ce6351314c228a3b80cced2e2f1746596d1ddc (patch) | |
tree | 048361893da884565693a948ee54b69d7502bad9 /numpy/linalg/tests | |
parent | ec3603f43aaf4056bdd78401c23c937c5760f673 (diff) | |
download | numpy-f2ce6351314c228a3b80cced2e2f1746596d1ddc.tar.gz |
TST: linalg: add more comprehensive test cases for linalg funcs
Also remove TestCase subclassing, so that generator tests work. Also
fix bugs in the existing generator tests that were not actually run
previously.
Diffstat (limited to 'numpy/linalg/tests')
-rw-r--r-- | numpy/linalg/tests/test_linalg.py | 555 |
1 files changed, 332 insertions, 223 deletions
diff --git a/numpy/linalg/tests/test_linalg.py b/numpy/linalg/tests/test_linalg.py index 2c003a072..aeeedadf9 100644 --- a/numpy/linalg/tests/test_linalg.py +++ b/numpy/linalg/tests/test_linalg.py @@ -5,6 +5,8 @@ from __future__ import division, absolute_import, print_function import os import sys +import itertools +import traceback import numpy as np from numpy import array, single, double, csingle, cdouble, dot, identity @@ -12,8 +14,9 @@ from numpy import multiply, atleast_2d, inf, asarray, matrix from numpy import linalg from numpy.linalg import matrix_power, norm, matrix_rank from numpy.testing import ( - TestCase, assert_, assert_equal, assert_raises, assert_array_equal, - assert_almost_equal, run_module_suite, dec + assert_, assert_equal, assert_raises, assert_array_equal, + assert_almost_equal, assert_allclose, run_module_suite, + dec ) @@ -45,162 +48,287 @@ def get_complex_dtype(dtype): return {single: csingle, double: cdouble, csingle: csingle, cdouble: cdouble}[dtype] +def get_rtol(dtype): + return 0.1 * np.sqrt(np.finfo(dtype).eps) -class LinalgTestCase(object): - def test_single(self): - a = array([[1., 2.], [3., 4.]], dtype=single) - b = array([2., 1.], dtype=single) - self.do(a, b) - - def test_double(self): - a = array([[1., 2.], [3., 4.]], dtype=double) - b = array([2., 1.], dtype=double) - self.do(a, b) - - def test_double_2(self): - a = array([[1., 2.], [3., 4.]], dtype=double) - b = array([[2., 1., 4.], [3., 4., 6.]], dtype=double) - self.do(a, b) - - def test_csingle(self): - a = array([[1.+2j, 2+3j], [3+4j, 4+5j]], dtype=csingle) - b = array([2.+1j, 1.+2j], dtype=csingle) - self.do(a, b) - - def test_cdouble(self): - a = array([[1.+2j, 2+3j], [3+4j, 4+5j]], dtype=cdouble) - b = array([2.+1j, 1.+2j], dtype=cdouble) - self.do(a, b) - - def test_cdouble_2(self): - a = array([[1.+2j, 2+3j], [3+4j, 4+5j]], dtype=cdouble) - b = array([[2.+1j, 1.+2j, 1+3j], [1-2j, 1-3j, 1-6j]], dtype=cdouble) - self.do(a, b) +class LinalgCase(object): + def __init__(self, name, a, b, exception_cls=None): + assert isinstance(name, str) + self.name = name + self.a = a + self.b = b + self.exception_cls = exception_cls - def test_empty(self): - a = atleast_2d(array([], dtype = double)) - b = atleast_2d(array([], dtype = double)) + def check(self, do): + if self.exception_cls is None: + do(self.a, self.b) + else: + assert_raises(self.exception_cls, do, self.a, self.b) + + def __repr__(self): + return "<LinalgCase: %s>" % (self.name,) + + +# +# Base test cases +# + +SQUARE_CASES = [ + LinalgCase("single", + array([[1., 2.], [3., 4.]], dtype=single), + array([2., 1.], dtype=single)), + LinalgCase("double", + array([[1., 2.], [3., 4.]], dtype=double), + array([2., 1.], dtype=double)), + LinalgCase("double_2", + array([[1., 2.], [3., 4.]], dtype=double), + array([[2., 1., 4.], [3., 4., 6.]], dtype=double)), + LinalgCase("csingle", + array([[1.+2j, 2+3j], [3+4j, 4+5j]], dtype=csingle), + array([2.+1j, 1.+2j], dtype=csingle)), + LinalgCase("cdouble", + array([[1.+2j, 2+3j], [3+4j, 4+5j]], dtype=cdouble), + array([2.+1j, 1.+2j], dtype=cdouble)), + LinalgCase("cdouble_2", + array([[1.+2j, 2+3j], [3+4j, 4+5j]], dtype=cdouble), + array([[2.+1j, 1.+2j, 1+3j], [1-2j, 1-3j, 1-6j]], dtype=cdouble)), + LinalgCase("empty", + atleast_2d(array([], dtype = double)), + atleast_2d(array([], dtype = double)), + linalg.LinAlgError), + LinalgCase("nonarray", + [[1, 2], [3, 4]], + [2, 1]), + LinalgCase("matrix_b_only", + array([[1., 2.], [3., 4.]]), + matrix([2., 1.]).T), + LinalgCase("matrix_a_and_b", + matrix([[1., 2.], [3., 4.]]), + matrix([2., 1.]).T), +] + +NONSQUARE_CASES = [ + LinalgCase("single_nsq_1", + array([[1., 2., 3.], [3., 4., 6.]], dtype=single), + array([2., 1.], dtype=single)), + LinalgCase("single_nsq_2", + array([[1., 2.], [3., 4.], [5., 6.]], dtype=single), + array([2., 1., 3.], dtype=single)), + LinalgCase("double_nsq_1", + array([[1., 2., 3.], [3., 4., 6.]], dtype=double), + array([2., 1.], dtype=double)), + LinalgCase("double_nsq_2", + array([[1., 2.], [3., 4.], [5., 6.]], dtype=double), + array([2., 1., 3.], dtype=double)), + LinalgCase("csingle_nsq_1", + array([[1.+1j, 2.+2j, 3.-3j], [3.-5j, 4.+9j, 6.+2j]], dtype=csingle), + array([2.+1j, 1.+2j], dtype=csingle)), + LinalgCase("csingle_nsq_2", + array([[1.+1j, 2.+2j], [3.-3j, 4.-9j], [5.-4j, 6.+8j]], dtype=csingle), + array([2.+1j, 1.+2j, 3.-3j], dtype=csingle)), + LinalgCase("cdouble_nsq_1", + array([[1.+1j, 2.+2j, 3.-3j], [3.-5j, 4.+9j, 6.+2j]], dtype=cdouble), + array([2.+1j, 1.+2j], dtype=cdouble)), + LinalgCase("cdouble_nsq_2", + array([[1.+1j, 2.+2j], [3.-3j, 4.-9j], [5.-4j, 6.+8j]], dtype=cdouble), + array([2.+1j, 1.+2j, 3.-3j], dtype=cdouble)), + LinalgCase("cdouble_nsq_1_2", + array([[1.+1j, 2.+2j, 3.-3j], [3.-5j, 4.+9j, 6.+2j]], dtype=cdouble), + array([[2.+1j, 1.+2j], [1-1j, 2-2j]], dtype=cdouble)), + LinalgCase("cdouble_nsq_2_2", + array([[1.+1j, 2.+2j], [3.-3j, 4.-9j], [5.-4j, 6.+8j]], dtype=cdouble), + array([[2.+1j, 1.+2j], [1-1j, 2-2j], [1-1j, 2-2j]], dtype=cdouble)), +] + +HERMITIAN_CASES = [ + LinalgCase("hsingle", + array([[1., 2.], [2., 1.]], dtype=single), + None), + LinalgCase("hdouble", + array([[1., 2.], [2., 1.]], dtype=double), + None), + LinalgCase("hcsingle", + array([[1., 2+3j], [2-3j, 1]], dtype=csingle), + None), + LinalgCase("hcdouble", + array([[1., 2+3j], [2-3j, 1]], dtype=cdouble), + None), + LinalgCase("hempty", + atleast_2d(array([], dtype = double)), + None, + linalg.LinAlgError), + LinalgCase("hnonarray", + [[1, 2], [2, 1]], + None), + LinalgCase("matrix_b_only", + array([[1., 2.], [2., 1.]]), + None), + LinalgCase("hmatrix_a_and_b", + matrix([[1., 2.], [2., 1.]]), + None) +] + + +# +# Gufunc test cases +# + +GENERALIZED_SQUARE_CASES = [] +GENERALIZED_NONSQUARE_CASES = [] +GENERALIZED_HERMITIAN_CASES = [] + +for tgt, src in ((GENERALIZED_SQUARE_CASES, SQUARE_CASES), + (GENERALIZED_NONSQUARE_CASES, NONSQUARE_CASES), + (GENERALIZED_HERMITIAN_CASES, HERMITIAN_CASES)): + for case in src: + if not isinstance(case.a, np.ndarray): + continue + + a = np.array([case.a, 2*case.a, 3*case.a]) + if case.b is None: + b = None + else: + b = np.array([case.b, 7*case.b, 6*case.b]) + new_case = LinalgCase(case.name + "_tile3", a, b, + case.exception_cls) + tgt.append(new_case) + + a = np.array([case.a]*2*3).reshape((3, 2) + case.a.shape) + if case.b is None: + b = None + else: + b = np.array([case.b]*2*3).reshape((3, 2) + case.b.shape) + new_case = LinalgCase(case.name + "_tile213", a, b, + case.exception_cls) + tgt.append(new_case) + +# +# Generate stride combination variations of the above +# + +def _stride_comb_iter(x): + """ + Generate cartesian product of strides for all axes + """ + + if not isinstance(x, np.ndarray): + yield x, "nop" + return + + stride_set = [(1,)]*x.ndim + stride_set[-1] = (1, 3, -4) + if x.ndim > 1: + stride_set[-2] = (1, 3, -4) + if x.ndim > 2: + stride_set[-3] = (1, -4) + + for repeats in itertools.product(*tuple(stride_set)): + new_shape = [abs(a*b) for a, b in zip(x.shape, repeats)] + slices = tuple([slice(None, None, repeat) for repeat in repeats]) + + # new array with different strides, but same data + xi = np.empty(new_shape, dtype=x.dtype) + xi.view(np.uint32).fill(0xdeadbeef) + xi = xi[slices] + xi[...] = x + xi = xi.view(x.__class__) + assert np.all(xi == x) + yield xi, "stride_" + "_".join(["%+d" % j for j in repeats]) + +for src in (SQUARE_CASES, + NONSQUARE_CASES, + HERMITIAN_CASES, + GENERALIZED_SQUARE_CASES, + GENERALIZED_NONSQUARE_CASES, + GENERALIZED_HERMITIAN_CASES): + + new_cases = [] + for case in src: + for a, a_tag in _stride_comb_iter(case.a): + for b, b_tag in _stride_comb_iter(case.b): + new_case = LinalgCase(case.name + "_" + a_tag + "_" + b_tag, a, b, + exception_cls=case.exception_cls) + new_cases.append(new_case) + src.extend(new_cases) + + +# +# Test different routines against the above cases +# + +def _check_cases(func, cases): + for case in cases: try: - self.do(a, b) - raise AssertionError("%s should fail with empty matrices", self.__name__[5:]) - except linalg.LinAlgError as e: - pass + case.check(func) + except: + msg = "In test case: %r\n\n" % case + msg += traceback.format_exc() + raise AssertionError(msg) - def test_nonarray(self): - a = [[1, 2], [3, 4]] - b = [2, 1] - self.do(a, b) +class LinalgTestCase(object): + def test_sq_cases(self): + _check_cases(self.do, SQUARE_CASES) - def test_matrix_b_only(self): - """Check that matrix type is preserved.""" - a = array([[1., 2.], [3., 4.]]) - b = matrix([2., 1.]).T - self.do(a, b) - def test_matrix_a_and_b(self): - """Check that matrix type is preserved.""" - a = matrix([[1., 2.], [3., 4.]]) - b = matrix([2., 1.]).T - self.do(a, b) +class LinalgNonsquareTestCase(object): + def test_sq_cases(self): + _check_cases(self.do, NONSQUARE_CASES) -class LinalgNonsquareTestCase(object): - def test_single_nsq_1(self): - a = array([[1., 2., 3.], [3., 4., 6.]], dtype=single) - b = array([2., 1.], dtype=single) - self.do(a, b) - - def test_single_nsq_2(self): - a = array([[1., 2.], [3., 4.], [5., 6.]], dtype=single) - b = array([2., 1., 3.], dtype=single) - self.do(a, b) - - def test_double_nsq_1(self): - a = array([[1., 2., 3.], [3., 4., 6.]], dtype=double) - b = array([2., 1.], dtype=double) - self.do(a, b) - - def test_double_nsq_2(self): - a = array([[1., 2.], [3., 4.], [5., 6.]], dtype=double) - b = array([2., 1., 3.], dtype=double) - self.do(a, b) - - def test_csingle_nsq_1(self): - a = array([[1.+1j, 2.+2j, 3.-3j], [3.-5j, 4.+9j, 6.+2j]], dtype=csingle) - b = array([2.+1j, 1.+2j], dtype=csingle) - self.do(a, b) - - def test_csingle_nsq_2(self): - a = array([[1.+1j, 2.+2j], [3.-3j, 4.-9j], [5.-4j, 6.+8j]], dtype=csingle) - b = array([2.+1j, 1.+2j, 3.-3j], dtype=csingle) - self.do(a, b) - - def test_cdouble_nsq_1(self): - a = array([[1.+1j, 2.+2j, 3.-3j], [3.-5j, 4.+9j, 6.+2j]], dtype=cdouble) - b = array([2.+1j, 1.+2j], dtype=cdouble) - self.do(a, b) - - def test_cdouble_nsq_2(self): - a = array([[1.+1j, 2.+2j], [3.-3j, 4.-9j], [5.-4j, 6.+8j]], dtype=cdouble) - b = array([2.+1j, 1.+2j, 3.-3j], dtype=cdouble) - self.do(a, b) - - def test_cdouble_nsq_1_2(self): - a = array([[1.+1j, 2.+2j, 3.-3j], [3.-5j, 4.+9j, 6.+2j]], dtype=cdouble) - b = array([[2.+1j, 1.+2j], [1-1j, 2-2j]], dtype=cdouble) - self.do(a, b) - - def test_cdouble_nsq_2_2(self): - a = array([[1.+1j, 2.+2j], [3.-3j, 4.-9j], [5.-4j, 6.+8j]], dtype=cdouble) - b = array([[2.+1j, 1.+2j], [1-1j, 2-2j], [1-1j, 2-2j]], dtype=cdouble) - self.do(a, b) - - -def _generalized_testcase(new_cls_name, old_cls): - def get_method(old_name, new_name): - def method(self): - base = old_cls() - def do(a, b): - a = np.array([a, a, a]) - b = np.array([b, b, b]) - self.do(a, b) - base.do = do - getattr(base, old_name)() - method.__name__ = new_name - return method - - dct = dict() - for old_name in dir(old_cls): - if old_name.startswith('test_'): - new_name = old_name + '_generalized' - dct[new_name] = get_method(old_name, new_name) - - return type(new_cls_name, (object,), dct) - -LinalgGeneralizedTestCase = _generalized_testcase( - 'LinalgGeneralizedTestCase', LinalgTestCase) -LinalgGeneralizedNonsquareTestCase = _generalized_testcase( - 'LinalgGeneralizedNonsquareTestCase', LinalgNonsquareTestCase) +class LinalgGeneralizedTestCase(object): + @dec.slow + def test_generalized_sq_cases(self): + _check_cases(self.do, GENERALIZED_SQUARE_CASES) + + +class LinalgGeneralizedNonsquareTestCase(object): + @dec.slow + def test_generalized_nonsq_cases(self): + _check_cases(self.do, GENERALIZED_NONSQUARE_CASES) + + +class HermitianTestCase(object): + def test_herm_cases(self): + _check_cases(self.do, HERMITIAN_CASES) + + +class HermitianGeneralizedTestCase(object): + @dec.slow + def test_generalized_herm_cases(self): + _check_cases(self.do, GENERALIZED_HERMITIAN_CASES) def dot_generalized(a, b): a = asarray(a) - if a.ndim == 3: - return np.array([dot(ax, bx) for ax, bx in zip(a, b)]) - elif a.ndim > 3: - raise ValueError("Not implemented...") - return dot(a, b) + if a.ndim >= 3: + if a.ndim == b.ndim: + # matrix x matrix + new_shape = a.shape[:-1] + b.shape[-1:] + elif a.ndim == b.ndim + 1: + # matrix x vector + new_shape = a.shape[:-1] + else: + raise ValueError("Not implemented...") + r = np.empty(new_shape, dtype=np.common_type(a, b)) + for c in itertools.product(*map(range, a.shape[:-2])): + r[c] = dot(a[c], b[c]) + return r + else: + return dot(a, b) + def identity_like_generalized(a): a = asarray(a) - if a.ndim == 3: - return np.array([identity(a.shape[-2]) for ax in a]) - elif a.ndim > 3: - raise ValueError("Not implemented...") - return identity(a.shape[0]) + if a.ndim >= 3: + r = np.empty(a.shape, dtype=a.dtype) + for c in itertools.product(*map(range, a.shape[:-2])): + r[c] = identity(a.shape[-2]) + return r + else: + return identity(a.shape[0]) -class TestSolve(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): +class TestSolve(LinalgTestCase, LinalgGeneralizedTestCase): def do(self, a, b): x = linalg.solve(a, b) assert_almost_equal(b, dot_generalized(a, x)) @@ -265,7 +393,7 @@ class TestSolve(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): assert_(isinstance(result, ArraySubclass)) -class TestInv(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): +class TestInv(LinalgTestCase, LinalgGeneralizedTestCase): def do(self, a, b): a_inv = linalg.inv(a) assert_almost_equal(dot_generalized(a, a_inv), @@ -295,7 +423,7 @@ class TestInv(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): assert_equal(a.shape, res.shape) -class TestEigvals(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): +class TestEigvals(LinalgTestCase, LinalgGeneralizedTestCase): def do(self, a, b): ev = linalg.eigvals(a) evalues, evectors = linalg.eig(a) @@ -311,13 +439,12 @@ class TestEigvals(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): yield check, dtype -class TestEig(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): +class TestEig(LinalgTestCase, LinalgGeneralizedTestCase): def do(self, a, b): evalues, evectors = linalg.eig(a) - if evectors.ndim == 3: - assert_almost_equal(dot_generalized(a, evectors), evectors * evalues[:, None,:]) - else: - assert_almost_equal(dot(a, evectors), multiply(evectors, evalues)) + assert_allclose(dot_generalized(a, evectors), + np.asarray(evectors) * np.asarray(evalues)[...,None,:], + rtol=get_rtol(evalues.dtype)) assert_(imply(isinstance(a, matrix), isinstance(evectors, matrix))) def test_types(self): @@ -333,16 +460,15 @@ class TestEig(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): assert_equal(v.dtype, get_complex_dtype(dtype)) for dtype in [single, double, csingle, cdouble]: - yield dtype + yield check, dtype -class TestSVD(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): +class TestSVD(LinalgTestCase, LinalgGeneralizedTestCase): def do(self, a, b): u, s, vt = linalg.svd(a, 0) - if u.ndim == 3: - assert_almost_equal(a, dot_generalized(u * s[:, None,:], vt)) - else: - assert_almost_equal(a, dot(multiply(u, s), vt)) + assert_allclose(a, dot_generalized(np.asarray(u) * np.asarray(s)[...,None,:], + np.asarray(vt)), + rtol=get_rtol(u.dtype)) assert_(imply(isinstance(a, matrix), isinstance(u, matrix))) assert_(imply(isinstance(a, matrix), isinstance(vt, matrix))) @@ -360,34 +486,34 @@ class TestSVD(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): yield check, dtype -class TestCondSVD(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): +class TestCondSVD(LinalgTestCase, LinalgGeneralizedTestCase): def do(self, a, b): c = asarray(a) # a might be a matrix s = linalg.svd(c, compute_uv=False) old_assert_almost_equal(s[0]/s[-1], linalg.cond(a), decimal=5) -class TestCond2(LinalgTestCase, TestCase): +class TestCond2(LinalgTestCase): def do(self, a, b): c = asarray(a) # a might be a matrix s = linalg.svd(c, compute_uv=False) old_assert_almost_equal(s[0]/s[-1], linalg.cond(a, 2), decimal=5) -class TestCondInf(TestCase): +class TestCondInf(object): def test(self): A = array([[1., 0, 0], [0, -2., 0], [0, 0, 3.]]) assert_almost_equal(linalg.cond(A, inf), 3.) -class TestPinv(LinalgTestCase, TestCase): +class TestPinv(LinalgTestCase): def do(self, a, b): a_ginv = linalg.pinv(a) assert_almost_equal(dot(a, a_ginv), identity(asarray(a).shape[0])) assert_(imply(isinstance(a, matrix), isinstance(a_ginv, matrix))) -class TestDet(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): +class TestDet(LinalgTestCase, LinalgGeneralizedTestCase): def do(self, a, b): d = linalg.det(a) (s, ld) = linalg.slogdet(a) @@ -421,12 +547,15 @@ class TestDet(LinalgTestCase, LinalgGeneralizedTestCase, TestCase): def test_types(self): def check(dtype): x = np.array([[1, 0.5], [0.5, 1]], dtype=dtype) - assert_equal(np.linalg.det(x), get_real_dtype(dtype)) + assert_equal(np.linalg.det(x).dtype, dtype) + ph, s = np.linalg.slogdet(x) + assert_equal(s.dtype, get_real_dtype(dtype)) + assert_equal(ph.dtype, dtype) for dtype in [single, double, csingle, cdouble]: yield check, dtype -class TestLstsq(LinalgTestCase, LinalgNonsquareTestCase, TestCase): +class TestLstsq(LinalgTestCase, LinalgNonsquareTestCase): def do(self, a, b): arr = np.asarray(a) m, n = arr.shape @@ -506,70 +635,38 @@ class TestMatrixPower(object): lambda: matrix_power(self.noninv, -1)) -class TestBoolPower(TestCase): +class TestBoolPower(object): def test_square(self): A = array([[True, False], [True, True]]) assert_equal(matrix_power(A, 2), A) -class HermitianTestCase(object): - def test_single(self): - a = array([[1., 2.], [2., 1.]], dtype=single) - self.do(a, None) - - def test_double(self): - a = array([[1., 2.], [2., 1.]], dtype=double) - self.do(a, None) - - def test_csingle(self): - a = array([[1., 2+3j], [2-3j, 1]], dtype=csingle) - self.do(a, None) - - def test_cdouble(self): - a = array([[1., 2+3j], [2-3j, 1]], dtype=cdouble) - self.do(a, None) - - def test_empty(self): - a = atleast_2d(array([], dtype = double)) - assert_raises(linalg.LinAlgError, self.do, a, None) - - def test_nonarray(self): - a = [[1, 2], [2, 1]] - self.do(a, None) - - def test_matrix_b_only(self): - """Check that matrix type is preserved.""" - a = array([[1., 2.], [2., 1.]]) - self.do(a, None) - - def test_matrix_a_and_b(self): - """Check that matrix type is preserved.""" - a = matrix([[1., 2.], [2., 1.]]) - self.do(a, None) - - -HermitianGeneralizedTestCase = _generalized_testcase( - 'HermitianGeneralizedTestCase', HermitianTestCase) - -class TestEigvalsh(HermitianTestCase, HermitianGeneralizedTestCase, TestCase): +class TestEigvalsh(HermitianTestCase, HermitianGeneralizedTestCase): def do(self, a, b): # note that eigenvalue arrays must be sorted since # their order isn't guaranteed. - ev = linalg.eigvalsh(a) + ev = linalg.eigvalsh(a, 'L') evalues, evectors = linalg.eig(a) ev.sort(axis=-1) evalues.sort(axis=-1) - assert_almost_equal(ev, evalues) + assert_allclose(ev, evalues, + rtol=get_rtol(ev.dtype)) + + ev2 = linalg.eigvalsh(a, 'U') + ev2.sort(axis=-1) + assert_allclose(ev2, evalues, + rtol=get_rtol(ev.dtype)) def test_types(self): def check(dtype): x = np.array([[1, 0.5], [0.5, 1]], dtype=dtype) - assert_equal(np.linalg.eigvalsh(x), get_real_dtype(dtype)) + w = np.linalg.eigvalsh(x) + assert_equal(w.dtype, get_real_dtype(dtype)) for dtype in [single, double, csingle, cdouble]: yield check, dtype -class TestEigh(HermitianTestCase, HermitianGeneralizedTestCase, TestCase): +class TestEigh(HermitianTestCase, HermitianGeneralizedTestCase): def do(self, a, b): # note that eigenvalue arrays must be sorted since # their order isn't guaranteed. @@ -579,17 +676,29 @@ class TestEigh(HermitianTestCase, HermitianGeneralizedTestCase, TestCase): evalues.sort(axis=-1) assert_almost_equal(ev, evalues) + assert_allclose(dot_generalized(a, evc), + np.asarray(ev)[...,None,:] * np.asarray(evc), + rtol=get_rtol(ev.dtype)) + + ev2, evc2 = linalg.eigh(a, 'U') + ev2.sort(axis=-1) + assert_almost_equal(ev2, evalues) + + assert_allclose(dot_generalized(a, evc2), + np.asarray(ev2)[...,None,:] * np.asarray(evc2), + rtol=get_rtol(ev.dtype)) + def test_types(self): def check(dtype): x = np.array([[1, 0.5], [0.5, 1]], dtype=dtype) - w, v = np.linalg.eig(x) - assert_equal(w, get_real_dtype(dtype)) - assert_equal(v, dtype) + w, v = np.linalg.eigh(x) + assert_equal(w.dtype, get_real_dtype(dtype)) + assert_equal(v.dtype, dtype) for dtype in [single, double, csingle, cdouble]: yield check, dtype -class _TestNorm(TestCase): +class _TestNorm(object): dt = None dec = None @@ -640,9 +749,9 @@ class _TestNorm(TestCase): assert_almost_equal(norm(A, 2), 9.1231056256176615) assert_almost_equal(norm(A, -2), 0.87689437438234041) - self.assertRaises(ValueError, norm, A, 'nofro') - self.assertRaises(ValueError, norm, A, -3) - self.assertRaises(ValueError, norm, A, 0) + assert_raises(ValueError, norm, A, 'nofro') + assert_raises(ValueError, norm, A, -3) + assert_raises(ValueError, norm, A, 0) def test_axis(self): # Vector norms. @@ -687,20 +796,20 @@ class _TestNorm(TestCase): # Using `axis=<integer>` or passing in a 1-D array implies vector # norms are being computed, so also using `ord='fro'` raises a # ValueError. - self.assertRaises(ValueError, norm, A, 'fro', 0) - self.assertRaises(ValueError, norm, [3, 4], 'fro', None) + assert_raises(ValueError, norm, A, 'fro', 0) + assert_raises(ValueError, norm, [3, 4], 'fro', None) # Similarly, norm should raise an exception when ord is any finite # number other than 1, 2, -1 or -2 when computing matrix norms. for order in [0, 3]: - self.assertRaises(ValueError, norm, A, order, None) - self.assertRaises(ValueError, norm, A, order, (0, 1)) - self.assertRaises(ValueError, norm, B, order, (1, 2)) + assert_raises(ValueError, norm, A, order, None) + assert_raises(ValueError, norm, A, order, (0, 1)) + assert_raises(ValueError, norm, B, order, (1, 2)) # Invalid axis - self.assertRaises(ValueError, norm, B, None, 3) - self.assertRaises(ValueError, norm, B, None, (2, 3)) - self.assertRaises(ValueError, norm, B, None, (0, 1, 2)) + assert_raises(ValueError, norm, B, None, 3) + assert_raises(ValueError, norm, B, None, (2, 3)) + assert_raises(ValueError, norm, B, None, (0, 1, 2)) class TestNormDouble(_TestNorm): @@ -751,7 +860,7 @@ def test_reduced_rank(): assert_equal(matrix_rank(X), 8) -class TestQR(TestCase): +class TestQR(object): def check_qr(self, a): # This test expects the argument `a` to be an ndarray or @@ -793,7 +902,7 @@ class TestQR(TestCase): def test_qr_empty(self): a = np.zeros((0, 2)) - self.assertRaises(linalg.LinAlgError, linalg.qr, a) + assert_raises(linalg.LinAlgError, linalg.qr, a) def test_mode_raw(self): # The factorization is not unique and varies between libraries, |