Returns the matrix square root of a positive semidefinite matrix,
    truncating negative eigenvalues.
172909
def sqrtm_psd(A, est_error=True, check_finite=True):
    
    w, v = eigh(A, check_finite=check_finite)
    mask = w <= 0
    w[mask] = 0
    np.sqrt(w, out=w)
    A_sqrt = (v * w).dot(v.conj().T)

    if est_error:
        return A_sqrt, np.linalg.norm(np.dot(A_sqrt, A_sqrt) - A, 'fro')
    else:
        return A_sqrt
953081
public int rank() {
        if (rows == 0 || columns == 0) {
            return 0;
        }

        // TODO:
        // handle small (1x1, 1xn, nx1, 2x2, 2xn, nx2, 3x3, 3xn, nx3)
        // matrices without SVD

        MatrixDecompositor decompositor = withDecompositor(LinearAlgebra.SVD);
        Matrix[] usv = decompositor.decompose();
        // TODO: Where is my pattern matching?
        Matrix s = usv[1];
        double tolerance = Math.max(rows, columns) * s.get(0, 0) * Matrices.EPS;

        int result = 0;
        for (int i = 0; i < s.rows(); i++) {
            if (s.get(i, i) > tolerance) {
                result++;
            }
        }

        return result;
    }
324694
public Matrix solve(Matrix rhs) {
        if (M != N || rhs.M != N || rhs.N != 1) {
            throw new RuntimeException("Illegal matrix dimensions.");
        }

        // create copies of the data
        Matrix A = new Matrix(this);
        Matrix b = new Matrix(rhs);

        // Gaussian elimination with partial pivoting
        for (int i = 0; i < N; i++) {

            // find pivot row and swap
            int max = i;
            for (int j = i + 1; j < N; j++) {
                if (Math.abs(A.data[j][i]) > Math.abs(A.data[max][i])) {
                    max = j;
                }
            }
            A.swap(i, max);
            b.swap(i, max);

            // singular
            if (A.data[i][i] == 0.0) {
                throw new RuntimeException("Matrix is singular.");
            }

            // pivot within b
            for (int j = i + 1; j < N; j++) {
                b.data[j][0] -= b.data[i][0] * A.data[j][i] / A.data[i][i];
            }

            // pivot within A
            for (int j = i + 1; j < N; j++) {
                double m = A.data[j][i] / A.data[i][i];
                for (int k = i + 1; k < N; k++) {
                    A.data[j][k] -= A.data[i][k] * m;
                }
                A.data[j][i] = 0.0;
            }
        }

        // back substitution
        Matrix x = new Matrix(N, 1);
        for (int j = N - 1; j >= 0; j--) {
            double t = 0.0;
            for (int k = j + 1; k < N; k++) {
                t += A.data[j][k] * x.data[k][0];
            }
            x.data[j][0] = (b.data[j][0] - t) / A.data[j][j];
        }
        return x;

    }
129669
def symbolic_Rz_matrix(symbolic_theta):
    """"""
    return sympy.Matrix([
        [sympy.cos(symbolic_theta), -sympy.sin(symbolic_theta), 0],
        [sympy.sin(symbolic_theta), sympy.cos(symbolic_theta), 0],
        [0, 0, 1]
    ])
728440
def Tmatrix(X):
    
    T = [[0., 0., 0.], [0., 0., 0.], [0., 0., 0.]]
    for row in X:
        for k in range(3):
            for l in range(3):
                T[k][l] += row[k] * row[l]
    return T
737523
def get_distance_matrix(x):
    """"""
    square = nd.sum(x ** 2.0, axis=1, keepdims=True)
    distance_square = square + square.transpose() - (2.0 * nd.dot(x, x.transpose()))
    return nd.sqrt(distance_square)
1300178
def diagonalize(self):
        ''''''
        self.eigvals, self.eigvecs = np.linalg.eig(
            (self.tensor.transpose() + self.tensor) / 2.0)
        self.eigvals = np.diag(np.dot(
            np.dot(self.eigvecs.transpose(), self.tensor), self.eigvecs))
146962
def left_eliminate!
      #    inv = Algebra.SquareMatrix(ground, rsize).unity
      k = ground.unity
      pi = 0
      each_j do |j|
        next unless i = (pi...rsize).find { |i1| !self[i1, j].zero? }
        if i != pi
          swap_r!(pi, i) # ; inv.swap_r!(pi, i)
          k = -k
        end
        c = ground.unity / self[pi, j] # this lets the entries be in ground
        multiply_r!(pi, c) # ; inv.multiply_r!(pi, c)
        k *= c
        each_i do |i0|
          next if i0 == pi
          d = self[i0, j] # / self[pi, j]
          mix_r!(i0, pi, -d) # ; inv.mix_r!(i0, pi, -d)
        end
        pi += 1
      end
      [left, k]
    end
517864
public Matrix getD () {
		Matrix X = new Matrix(n,n);
		double[][] D = X.getArray();
		for (int i = 0; i < n; i++) {
			for (int j = 0; j < n; j++) {
				D[i][j] = 0.0;
			}
			D[i][i] = d[i];
			if (e[i] > 0) {
				D[i][i+1] = e[i];
			} else if (e[i] < 0) {
				D[i][i-1] = e[i];
			}
		}
		return X;
	}
600535
def Ry_matrix(theta):
    """"""
    return np.array([
        [np.cos(theta), 0, np.sin(theta)],
        [0, 1, 0],
        [-np.sin(theta), 0, np.cos(theta)]
    ])
1283103
def get_pij_matrix(t, diag, A, A_inv):
    
    return A.dot(np.diag(np.exp(diag * t))).dot(A_inv)
516223
public static FloatMatrix SVDValues(ComplexFloatMatrix A) {
        int m = A.rows;
        int n = A.columns;
        FloatMatrix S = new FloatMatrix(min(m, n));
        float[] rwork = new float[5*min(m,n)];

        int info = NativeBlas.cgesvd('N', 'N', m, n, A.dup().data, 0, m, S.data, 0, null, 0, 1, null, 0, min(m,n), rwork, 0);

        if (info > 0) {
          throw new LapackConvergenceException("GESVD", info + " superdiagonals of an intermediate bidiagonal form failed to converge.");
        }

        return S;
    }
765172
def normalize_matrix(matrix):
  """"""
  abs_matrix = np.abs(matrix.copy())
  return abs_matrix / abs_matrix.max()
600155
public static DMatrixRMaj nullspaceQRP( DMatrixRMaj A , int totalSingular ) {
        SolveNullSpaceQRP_DDRM solver = new SolveNullSpaceQRP_DDRM();

        DMatrixRMaj nullspace = new DMatrixRMaj(1,1);

        if( !solver.process(A,totalSingular,nullspace))
            throw new RuntimeException("Solver failed. try SVD based method instead?");

        return nullspace;
    }
295917
public static double[][] pseudoInverse(double[][] matrix){
		if(isSolverUseApacheCommonsMath) {
			// Use LU from common math
			SingularValueDecomposition svd = new SingularValueDecomposition(new Array2DRowRealMatrix(matrix));
			double[][] matrixInverse = svd.getSolver().getInverse().getData();

			return matrixInverse;
		}
		else {
			return org.jblas.Solve.pinv(new org.jblas.DoubleMatrix(matrix)).toArray2();
		}
	}
699882
def funTransQuadF(k, s):
    
    sqrtk = np.sqrt(complex(k))
    a = np.cos(sqrtk * s)
    b = np.sin(sqrtk * s) / sqrtk
    c = -sqrtk * np.sin(sqrtk * s)
    d = np.cos(sqrtk * s)
    return np.matrix([[a.real, b.real], [c.real, d.real]], dtype=np.double)
928914
def _svd(cls, matrix, num_concepts=5):
        
        u, s, v = svds(matrix, k=num_concepts)
        return u, s, v