Python numpy.linalg() Examples
The following are 30 code examples for showing how to use numpy.linalg(). These examples are extracted from open source projects. You can vote up the ones you like or vote down the ones you don't like, and go to the original project or source file by following the links above each example.
You may check out the related API usage on the sidebar.
You may also want to check out all available functions/classes of the module
numpy
, or try the search function
.
Example 1
Project: pyscf Author: pyscf File: 001-remove_linear_dep.py License: Apache License 2.0 | 6 votes |
def eig(h, s): from pyscf import symm nirrep = len(mol.symm_orb) h = symm.symmetrize_matrix(h, mol.symm_orb) s = symm.symmetrize_matrix(s, mol.symm_orb) cs = [] es = [] # # Linear dependency are removed by looping over different symmetry irreps. # for ir in range(nirrep): d, t = numpy.linalg.eigh(s[ir]) x = t[:,d>1e-8] / numpy.sqrt(d[d>1e-8]) xhx = reduce(numpy.dot, (x.T, h[ir], x)) e, c = numpy.linalg.eigh(xhx) cs.append(reduce(numpy.dot, (mol.symm_orb[ir], x, c))) es.append(e) e = numpy.hstack(es) c = numpy.hstack(cs) return e, c
Example 2
Project: pyscf Author: pyscf File: 42-remove_linear_dep.py License: Apache License 2.0 | 6 votes |
def eig(h, s): from pyscf import symm nirrep = len(mol.symm_orb) h = symm.symmetrize_matrix(h, mol.symm_orb) s = symm.symmetrize_matrix(s, mol.symm_orb) cs = [] es = [] # # Linear dependency are removed by looping over different symmetry irreps. # for ir in range(nirrep): d, t = numpy.linalg.eigh(s[ir]) x = t[:,d>1e-8] / numpy.sqrt(d[d>1e-8]) xhx = reduce(numpy.dot, (x.T, h[ir], x)) e, c = numpy.linalg.eigh(xhx) cs.append(reduce(numpy.dot, (mol.symm_orb[ir], x, c))) es.append(e) e = numpy.hstack(es) c = numpy.hstack(cs) return e, c
Example 3
Project: SparseSC Author: microsoft File: fit_fast.py License: MIT License | 6 votes |
def _weights(V , X_treated, X_control, w_pen): V = np.diag(V) #make square #weights = np.zeros((X_control.shape[0], X_treated.shape[0])) w_pen_mat = 2 * w_pen * np.diag(np.ones(X_control.shape[0])) A = X_control.dot(2 * V).dot(X_control.T) + w_pen_mat # 5 B = ( X_treated.dot(2 * V).dot(X_control.T).T + 2 * w_pen / X_control.shape[0] ) # 6 try: b = scipy.linalg.solve(A, B) except scipy.linalg.LinAlgError as exc: print("Unique weights not possible.") if w_pen == 0: print("Try specifying a very small w_pen rather than 0.") raise exc return b
Example 4
Project: D-VAE Author: muhanzhang File: test_slinalg.py License: MIT License | 6 votes |
def test_eigvalsh(): if not imported_scipy: raise SkipTest("Scipy needed for the geigvalsh op.") import scipy.linalg A = theano.tensor.dmatrix('a') B = theano.tensor.dmatrix('b') f = function([A, B], eigvalsh(A, B)) rng = numpy.random.RandomState(utt.fetch_seed()) a = rng.randn(5, 5) a = a + a.T for b in [10 * numpy.eye(5, 5) + rng.randn(5, 5)]: w = f(a, b) refw = scipy.linalg.eigvalsh(a, b) numpy.testing.assert_array_almost_equal(w, refw) # We need to test None separatly, as otherwise DebugMode will # complain, as this isn't a valid ndarray. b = None B = theano.tensor.NoneConst f = function([A], eigvalsh(A, B)) w = f(a) refw = scipy.linalg.eigvalsh(a, b) numpy.testing.assert_array_almost_equal(w, refw)
Example 5
Project: D-VAE Author: muhanzhang File: test_slinalg.py License: MIT License | 6 votes |
def test_perform(self): if not imported_scipy: raise SkipTest('kron tests need the scipy package to be installed') for shp0 in [(2,), (2, 3), (2, 3, 4), (2, 3, 4, 5)]: x = tensor.tensor(dtype='floatX', broadcastable=(False,) * len(shp0)) a = numpy.asarray(self.rng.rand(*shp0)).astype(config.floatX) for shp1 in [(6,), (6, 7), (6, 7, 8), (6, 7, 8, 9)]: if len(shp0) + len(shp1) == 2: continue y = tensor.tensor(dtype='floatX', broadcastable=(False,) * len(shp1)) f = function([x, y], kron(x, y)) b = self.rng.rand(*shp1).astype(config.floatX) out = f(a, b) # Newer versions of scipy want 4 dimensions at least, # so we have to add a dimension to a and flatten the result. if len(shp0) + len(shp1) == 3: scipy_val = scipy.linalg.kron( a[numpy.newaxis, :], b).flatten() else: scipy_val = scipy.linalg.kron(a, b) utt.assert_allclose(out, scipy_val)
Example 6
Project: D-VAE Author: muhanzhang File: test_nlinalg.py License: MIT License | 6 votes |
def test_numpy_compare(self): rng = numpy.random.RandomState(utt.fetch_seed()) M = tensor.matrix("A", dtype=theano.config.floatX) V = tensor.vector("V", dtype=theano.config.floatX) a = rng.rand(4, 4).astype(theano.config.floatX) b = rng.rand(4).astype(theano.config.floatX) A = ( [None, 'fro', 'inf', '-inf', 1, -1, None, 'inf', '-inf', 0, 1, -1, 2, -2], [M, M, M, M, M, M, V, V, V, V, V, V, V, V], [a, a, a, a, a, a, b, b, b, b, b, b, b, b], [None, 'fro', inf, -inf, 1, -1, None, inf, -inf, 0, 1, -1, 2, -2]) for i in range(0, 14): f = function([A[1][i]], norm(A[1][i], A[0][i])) t_n = f(A[2][i]) n_n = numpy.linalg.norm(A[2][i], A[3][i]) assert _allclose(n_n, t_n)
Example 7
Project: vnpy_crypto Author: birforce File: linalg.py License: MIT License | 6 votes |
def logdet_symm(m, check_symm=False): """ Return log(det(m)) asserting positive definiteness of m. Parameters ---------- m : array-like 2d array that is positive-definite (and symmetric) Returns ------- logdet : float The log-determinant of m. """ from scipy import linalg if check_symm: if not np.all(m == m.T): # would be nice to short-circuit check raise ValueError("m is not symmetric.") c, _ = linalg.cho_factor(m, lower=True) return 2*np.sum(np.log(c.diagonal()))
Example 8
Project: OpenFermion Author: quantumlib File: _davidson_test.py License: Apache License 2.0 | 6 votes |
def setUp(self): """Sets up all variables needed for Davidson class.""" dimension = 10 matrix = generate_matrix(dimension) def mat_vec(vec): """Trivial matvec with a numpy matrix.""" return numpy.dot(matrix, vec) self.linear_operator = scipy.sparse.linalg.LinearOperator( (dimension, dimension), matvec=mat_vec) self.diagonal = numpy.diag(matrix) self.davidson = Davidson(linear_operator=self.linear_operator, linear_operator_diagonal=self.diagonal) self.matrix = matrix self.initial_guess = numpy.eye(self.matrix.shape[0], 10) self.eigen_values = numpy.array([ 1.15675714, 1.59132505, 2.62268014, 4.44533793, 5.3722743, 5.54393114, 7.73652405, 8.50089897, 9.4229309, 15.54405993, ])
Example 9
Project: OpenFermion Author: quantumlib File: _davidson_test.py License: Apache License 2.0 | 6 votes |
def setUp(self): """Sets up all variables needed for SparseDavidson class.""" logging.basicConfig(level=logging.INFO) self.dimension = 1000 self.sparse_matrix = generate_sparse_matrix(self.dimension) self.davidson_options = DavidsonOptions(max_subspace=100, max_iterations=50, real_only=True) # Checks for built-in eigh() function. self.eigen_values, self.eigen_vectors = numpy.linalg.eigh( self.sparse_matrix) self.assertAlmostEqual(get_difference( self.sparse_matrix, self.eigen_values, self.eigen_vectors), 0) # Makes sure eigenvalues are sorted. self.eigen_values = sorted(self.eigen_values)
Example 10
Project: OpenFermion Author: quantumlib File: _sparse_tools.py License: Apache License 2.0 | 6 votes |
def variance(operator, state): """Compute variance of operator with a state. Args: operator(scipy.sparse.spmatrix or scipy.sparse.linalg.LinearOperator): The operator whose expectation value is desired. state(numpy.ndarray or scipy.sparse.spmatrix): A numpy array representing a pure state or a sparse matrix representing a density matrix. Returns: A complex number giving the variance. Raises: ValueError: Input state has invalid format. """ return (expectation(operator ** 2, state) - expectation(operator, state) ** 2)
Example 11
Project: OpenFermion Author: quantumlib File: _sparse_tools.py License: Apache License 2.0 | 6 votes |
def get_gap(sparse_operator, initial_guess=None): """Compute gap between lowest eigenvalue and first excited state. Args: sparse_operator (LinearOperator): Operator to find the ground state of. initial_guess (ndarray): Initial guess for eigenspace. A good guess dramatically reduces the cost required to converge. Returns: A real float giving eigenvalue gap. """ if not is_hermitian(sparse_operator): raise ValueError('sparse_operator must be Hermitian.') values, _ = scipy.sparse.linalg.eigsh( sparse_operator, k=2, v0=initial_guess, which='SA', maxiter=1e7) gap = abs(values[1] - values[0]) return gap
Example 12
Project: OpenFermion Author: quantumlib File: _linear_qubit_operator.py License: Apache License 2.0 | 6 votes |
def generate_linear_qubit_operator(qubit_operator, n_qubits=None, options=None): """ Generates a LinearOperator from a QubitOperator. Args: qubit_operator(QubitOperator): A qubit operator to be applied on vectors. n_qubits(int): The total number of qubits options(LinearQubitOperatorOptions): Options for the ParallelLinearQubitOperator. Returns: linear_operator(scipy.sparse.linalg.LinearOperator): A linear operator. """ if options is None: linear_operator = LinearQubitOperator(qubit_operator, n_qubits) else: linear_operator = ParallelLinearQubitOperator( qubit_operator, n_qubits, options) return linear_operator
Example 13
Project: OpenFermion Author: quantumlib File: _davidson.py License: Apache License 2.0 | 6 votes |
def generate_random_vectors(row, col, real_only=False): """Generates orthonormal random vectors with col columns. Args: row(int): Number of rows for the vectors. col(int): Number of columns for the vectors. real_only(bool): Real vectors or complex ones. Returns: random_vectors(numpy.ndarray(complex)): Orthonormal random vectors. """ random_vectors = numpy.random.rand(row, col) if not real_only: random_vectors = random_vectors + numpy.random.rand(row, col) * 1.0j random_vectors = scipy.linalg.orth(random_vectors) return random_vectors
Example 14
Project: attention-lvcsr Author: rizar File: test_slinalg.py License: MIT License | 6 votes |
def test_eigvalsh(): if not imported_scipy: raise SkipTest("Scipy needed for the geigvalsh op.") import scipy.linalg A = theano.tensor.dmatrix('a') B = theano.tensor.dmatrix('b') f = function([A, B], eigvalsh(A, B)) rng = numpy.random.RandomState(utt.fetch_seed()) a = rng.randn(5, 5) a = a + a.T for b in [10 * numpy.eye(5, 5) + rng.randn(5, 5)]: w = f(a, b) refw = scipy.linalg.eigvalsh(a, b) numpy.testing.assert_array_almost_equal(w, refw) # We need to test None separatly, as otherwise DebugMode will # complain, as this isn't a valid ndarray. b = None B = theano.tensor.NoneConst f = function([A], eigvalsh(A, B)) w = f(a) refw = scipy.linalg.eigvalsh(a, b) numpy.testing.assert_array_almost_equal(w, refw)
Example 15
Project: attention-lvcsr Author: rizar File: test_slinalg.py License: MIT License | 6 votes |
def test_perform(self): if not imported_scipy: raise SkipTest('kron tests need the scipy package to be installed') for shp0 in [(2,), (2, 3), (2, 3, 4), (2, 3, 4, 5)]: x = tensor.tensor(dtype='floatX', broadcastable=(False,) * len(shp0)) a = numpy.asarray(self.rng.rand(*shp0)).astype(config.floatX) for shp1 in [(6,), (6, 7), (6, 7, 8), (6, 7, 8, 9)]: if len(shp0) + len(shp1) == 2: continue y = tensor.tensor(dtype='floatX', broadcastable=(False,) * len(shp1)) f = function([x, y], kron(x, y)) b = self.rng.rand(*shp1).astype(config.floatX) out = f(a, b) # Newer versions of scipy want 4 dimensions at least, # so we have to add a dimension to a and flatten the result. if len(shp0) + len(shp1) == 3: scipy_val = scipy.linalg.kron( a[numpy.newaxis, :], b).flatten() else: scipy_val = scipy.linalg.kron(a, b) utt.assert_allclose(out, scipy_val)
Example 16
Project: attention-lvcsr Author: rizar File: test_nlinalg.py License: MIT License | 6 votes |
def test_numpy_compare(self): rng = numpy.random.RandomState(utt.fetch_seed()) M = tensor.matrix("A", dtype=theano.config.floatX) V = tensor.vector("V", dtype=theano.config.floatX) a = rng.rand(4, 4).astype(theano.config.floatX) b = rng.rand(4).astype(theano.config.floatX) A = ( [None, 'fro', 'inf', '-inf', 1, -1, None, 'inf', '-inf', 0, 1, -1, 2, -2], [M, M, M, M, M, M, V, V, V, V, V, V, V, V], [a, a, a, a, a, a, b, b, b, b, b, b, b, b], [None, 'fro', inf, -inf, 1, -1, None, inf, -inf, 0, 1, -1, 2, -2]) for i in range(0, 14): f = function([A[1][i]], norm(A[1][i], A[0][i])) t_n = f(A[2][i]) n_n = numpy.linalg.norm(A[2][i], A[3][i]) assert _allclose(n_n, t_n)
Example 17
Project: neural-belief-tracker Author: nmrksic File: nbt.py License: Apache License 2.0 | 6 votes |
def xavier_vector(word, D=300): """ Returns a D-dimensional vector for the word. We hash the word to always get the same vector for the given word. """ seed_value = hash_string(word) numpy.random.seed(seed_value) neg_value = - math.sqrt(6)/math.sqrt(D) pos_value = math.sqrt(6)/math.sqrt(D) rsample = numpy.random.uniform(low=neg_value, high=pos_value, size=(D,)) norm = numpy.linalg.norm(rsample) rsample_normed = rsample/norm return rsample_normed
Example 18
Project: Uranium Author: Ultimaker File: MeshData.py License: GNU Lesser General Public License v3.0 | 6 votes |
def calculateNormalsFromIndexedVertices(vertices: numpy.ndarray, indices: numpy.ndarray, face_count: int) -> numpy.ndarray: """Calculate the normals of this mesh of triagles using indexes. :param vertices: :type{narray} list of vertices as a 1D list of float triples :param indices: :type{narray} list of indices as a 1D list of integers :param face_count: :type{integer} the number of triangles defined by the indices array :return: :type{narray} list normals as a 1D array of floats, each group of 3 floats is a vector """ start_time = time() # Numpy magic! # First, reset the normals normals = numpy.zeros((face_count*3, 3), dtype=numpy.float32) for face in indices[0:face_count]: normals[face[0]] = numpy.cross(vertices[face[0]] - vertices[face[1]], vertices[face[0]] - vertices[face[2]]) length = numpy.linalg.norm(normals[face[0]]) normals[face[0]] /= length normals[face[1]] = normals[face[0]] normals[face[2]] = normals[face[0]] end_time = time() Logger.log("d", "Calculating normals took %s seconds", end_time - start_time) return normals
Example 19
Project: FreeCAD_assembly2 Author: hamish2014 File: __init__.py License: GNU Lesser General Public License v2.1 | 6 votes |
def rotation_matrix_axis_and_angle_2(R, debug=False, errorThreshold=10**-7, msg=None): w, v = numpy.linalg.eig(R) #this method is not used at the primary method as numpy.linalg.eig does not return answers in high enough precision angle, axis = None, None eigErrors = abs(w -1) #errors from 1+0j i = (eigErrors == min(eigErrors)).tolist().index(True) axis = numpy.real(v[:,i]) if i != 1: angle = arccos2( numpy.real( w[1] ) ) else: angle = arccos2( numpy.real( w[0] ) ) error = norm(axis_rotation_matrix(angle, *axis) - R) if debug: print('rotation_matrix_axis_and_angle error %1.1e' % error) if error > errorThreshold: angle = -angle error = norm(axis_rotation_matrix(angle, *axis) - R) if error > errorThreshold: R_pickle_str = pickle.dumps(R) #R_abs_minus_identity = abs(R) - numpy.eye(3) print(R*R.transpose()) raise ValueError( 'rotation_matrix_axis_and_angle_2: no solution found! locals %s' % str(locals())) return axis, angle
Example 20
Project: pyscf Author: pyscf File: 001-remove_linear_dep.py License: Apache License 2.0 | 5 votes |
def eig(h, s): d, t = numpy.linalg.eigh(s) # Removing the eigenvectors assoicated to the smallest eigenvalue, the new # basis defined by x matrix has 139 vectors. x = t[:,d>1e-8] / numpy.sqrt(d[d>1e-8]) xhx = reduce(numpy.dot, (x.T, h, x)) e, c = numpy.linalg.eigh(xhx) c = numpy.dot(x, c) # Return 139 eigenvalues and 139 eigenvectors. return e, c # # Replacing the default eig function with the above one, the HF solver # generate only 139 canonical orbitals #
Example 21
Project: pyscf Author: pyscf File: 42-remove_linear_dep.py License: Apache License 2.0 | 5 votes |
def eig(h, s): d, t = numpy.linalg.eigh(s) # Removing the eigenvectors assoicated to the smallest eigenvalue, the new # basis defined by x matrix has 139 vectors. x = t[:,d>1e-8] / numpy.sqrt(d[d>1e-8]) xhx = reduce(numpy.dot, (x.T, h, x)) e, c = numpy.linalg.eigh(xhx) c = numpy.dot(x, c) # Return 139 eigenvalues and 139 eigenvectors. return e, c # # Replacing the default eig function with the above one, the HF solver # generate only 139 canonical orbitals #
Example 22
Project: tangent Author: google File: grads.py License: Apache License 2.0 | 5 votes |
def adet(z, x): """d|A|/dA = adj(A).T See Jacobi's formula: https://en.wikipedia.org/wiki/Jacobi%27s_formula """ adjugate = numpy.linalg.det(x) * numpy.linalg.pinv(x) d[x] = d[z] * numpy.transpose(adjugate) # # Tangent adjoints #
Example 23
Project: recruit Author: Frank-qlu File: test_public_api.py License: Apache License 2.0 | 5 votes |
def test_numpy_linalg(): bad_results = check_dir(np.linalg) assert bad_results == {}
Example 24
Project: hoggorm Author: olivertomic File: statTools.py License: BSD 2-Clause "Simplified" License | 5 votes |
def ortho(arr1, arr2): """ This function orthogonalises arr1 with respect to arr2. The function then returns orthogonalised array arr1_orth. PARAMETERS ---------- arr1 : numpy array A numpy array containing some data arr2 : numpy array A numpy array containing some data RETURNS ------- numpy array A numpy array holding orthogonalised numpy array ``arr1``. Examples -------- some examples """ # Find number of rows, such that identity matrix I can be created numberRows = numpy.shape(arr1)[0] I = numpy.identity(numberRows, float) # Compute transpose of arr1 arr2_T = numpy.transpose(arr2) term1 = numpy.linalg.inv(numpy.dot(arr2_T, arr2)) term2 = numpy.dot(arr2, term1) term3 = numpy.dot(term2, arr2_T) arr1_orth = numpy.dot((I - term3), arr1) return arr1_orth
Example 25
Project: hoggorm Author: olivertomic File: statTools.py License: BSD 2-Clause "Simplified" License | 5 votes |
def matrixRank(arr, tol=1e-8): """ Computes the rank of an array/matrix, i.e. number of linearly independent variables. This is not the same as numpy.rank() which only returns the number of ways (2-way, 3-way, etc) an array/matrix has. PARAMETERS ---------- arrX : numpy array A numpy array containing the data RETURNS ------- scalar Rank of matrix. Examples -------- >>> import hoggorm as ho >>> >>> # Get the rank of the data >>> ho.matrixRank(myData) >>> 8 """ if len(arr.shape) != 2: raise ValueError('Input must be a 2-d array or Matrix object') s = numpy.linalg.svd(arr, compute_uv=0) return numpy.sum(numpy.where(s > tol, 1, 0))
Example 26
Project: D-VAE Author: muhanzhang File: test_slinalg.py License: MIT License | 5 votes |
def test_eigvalsh_grad(): if not imported_scipy: raise SkipTest("Scipy needed for the geigvalsh op.") import scipy.linalg rng = numpy.random.RandomState(utt.fetch_seed()) a = rng.randn(5, 5) a = a + a.T b = 10 * numpy.eye(5, 5) + rng.randn(5, 5) tensor.verify_grad(lambda a, b: eigvalsh(a, b).dot([1, 2, 3, 4, 5]), [a, b], rng=numpy.random)
Example 27
Project: D-VAE Author: muhanzhang File: test_slinalg.py License: MIT License | 5 votes |
def test_solve_correctness(self): if not imported_scipy: raise SkipTest("Scipy needed for the Cholesky op.") rng = numpy.random.RandomState(utt.fetch_seed()) A = theano.tensor.matrix() b = theano.tensor.matrix() y = self.op(A, b) gen_solve_func = theano.function([A, b], y) cholesky_lower = Cholesky(lower=True) L = cholesky_lower(A) y_lower = self.op(L, b) lower_solve_func = theano.function([L, b], y_lower) cholesky_upper = Cholesky(lower=False) U = cholesky_upper(A) y_upper = self.op(U, b) upper_solve_func = theano.function([U, b], y_upper) b_val = numpy.asarray(rng.rand(5, 1), dtype=config.floatX) # 1-test general case A_val = numpy.asarray(rng.rand(5, 5), dtype=config.floatX) # positive definite matrix: A_val = numpy.dot(A_val.transpose(), A_val) assert numpy.allclose(scipy.linalg.solve(A_val, b_val), gen_solve_func(A_val, b_val)) # 2-test lower traingular case L_val = scipy.linalg.cholesky(A_val, lower=True) assert numpy.allclose(scipy.linalg.solve_triangular(L_val, b_val, lower=True), lower_solve_func(L_val, b_val)) # 3-test upper traingular case U_val = scipy.linalg.cholesky(A_val, lower=False) assert numpy.allclose(scipy.linalg.solve_triangular(U_val, b_val, lower=False), upper_solve_func(U_val, b_val))
Example 28
Project: D-VAE Author: muhanzhang File: test_slinalg.py License: MIT License | 5 votes |
def test_expm(): if not imported_scipy: raise SkipTest("Scipy needed for the expm op.") rng = numpy.random.RandomState(utt.fetch_seed()) A = rng.randn(5, 5).astype(config.floatX) ref = scipy.linalg.expm(A) x = tensor.matrix() m = expm(x) expm_f = function([x], m) val = expm_f(A) numpy.testing.assert_array_almost_equal(val, ref)
Example 29
Project: D-VAE Author: muhanzhang File: test_nlinalg.py License: MIT License | 5 votes |
def test_qr_modes(): rng = numpy.random.RandomState(utt.fetch_seed()) A = tensor.matrix("A", dtype=theano.config.floatX) a = rng.rand(4, 4).astype(theano.config.floatX) f = function([A], qr(A)) t_qr = f(a) n_qr = numpy.linalg.qr(a) assert _allclose(n_qr, t_qr) for mode in ["reduced", "r", "raw", "full", "economic"]: f = function([A], qr(A, mode)) t_qr = f(a) n_qr = numpy.linalg.qr(a, mode) if isinstance(n_qr, (list, tuple)): assert _allclose(n_qr[0], t_qr[0]) assert _allclose(n_qr[1], t_qr[1]) else: assert _allclose(n_qr, t_qr) try: n_qr = numpy.linalg.qr(a, "complete") f = function([A], qr(A, "complete")) t_qr = f(a) assert _allclose(n_qr, t_qr) except TypeError as e: assert "name 'complete' is not defined" in str(e)
Example 30
Project: D-VAE Author: muhanzhang File: test_nlinalg.py License: MIT License | 5 votes |
def test_svd(): rng = numpy.random.RandomState(utt.fetch_seed()) A = tensor.matrix("A", dtype=theano.config.floatX) U, V, T = svd(A) fn = function([A], [U, V, T]) a = rng.rand(4, 4).astype(theano.config.floatX) n_u, n_v, n_t = numpy.linalg.svd(a) t_u, t_v, t_t = fn(a) assert _allclose(n_u, t_u) assert _allclose(n_v, t_v) assert _allclose(n_t, t_t)