Python numpy.triu() Examples

The following are 30 code examples of numpy.triu(). 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 also want to check out all available functions/classes of the module numpy , or try the search function .
Example #1
Source File: test_blas.py    From Computable with MIT License 6 votes vote down vote up
def test_syrk(self):
        for f in _get_func('syrk'):
            c = f(a=self.a, alpha=1.)
            assert_array_almost_equal(np.triu(c), np.triu(self.t))

            c = f(a=self.a, alpha=1., lower=1)
            assert_array_almost_equal(np.tril(c), np.tril(self.t))

            c0 = np.ones(self.t.shape)
            c = f(a=self.a, alpha=1., beta=1., c=c0)
            assert_array_almost_equal(np.triu(c), np.triu(self.t+c0))

            c = f(a=self.a, alpha=1., trans=1)
            assert_array_almost_equal(np.triu(c), np.triu(self.tt))

    #prints '0-th dimension must be fixed to 3 but got 5', FIXME: suppress?
    # FIXME: how to catch the _fblas.error? 
Example #2
Source File: statesp_array_test.py    From python-control with BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
def test_dc_gain_integrator(self):
        """DC gain when eigenvalue at DC returns appropriately sized array of nan."""
        # the SISO case is also tested in test_dc_gain_{cont,discr}
        import itertools
        # iterate over input and output sizes, and continuous (dt=None) and discrete (dt=True) time
        for inputs, outputs, dt in itertools.product(range(1, 6), range(1, 6), [None, True]):
            states = max(inputs, outputs)

            # a matrix that is singular at DC, and has no "useless" states as in
            # _remove_useless_states
            a = np.triu(np.tile(2, (states, states)))
            # eigenvalues all +2, except for ...
            a[0, 0] = 0 if dt is None else 1
            b = np.eye(max(inputs, states))[:states, :inputs]
            c = np.eye(max(outputs, states))[:outputs, :states]
            d = np.zeros((outputs, inputs))
            sys = StateSpace(a, b, c, d, dt)
            dc = np.squeeze(np.tile(np.nan, (outputs, inputs)))
            np.testing.assert_array_equal(dc, sys.dcgain()) 
Example #3
Source File: statesp_test.py    From python-control with BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
def test_dc_gain_integrator(self):
        """DC gain when eigenvalue at DC returns appropriately sized array of nan."""
        # the SISO case is also tested in test_dc_gain_{cont,discr}
        import itertools
        # iterate over input and output sizes, and continuous (dt=None) and discrete (dt=True) time
        for inputs, outputs, dt in itertools.product(range(1, 6), range(1, 6), [None, True]):
            states = max(inputs, outputs)

            # a matrix that is singular at DC, and has no "useless" states as in
            # _remove_useless_states
            a = np.triu(np.tile(2, (states, states)))
            # eigenvalues all +2, except for ...
            a[0, 0] = 0 if dt is None else 1
            b = np.eye(max(inputs, states))[:states, :inputs]
            c = np.eye(max(outputs, states))[:outputs, :states]
            d = np.zeros((outputs, inputs))
            sys = StateSpace(a, b, c, d, dt)
            dc = np.squeeze(np.tile(np.nan, (outputs, inputs)))
            np.testing.assert_array_equal(dc, sys.dcgain()) 
Example #4
Source File: blow.py    From blow with Apache License 2.0 6 votes vote down vote up
def __init__(self,in_channel):
        super(InvConv,self).__init__()

        weight=np.random.randn(in_channel,in_channel)
        q,_=linalg.qr(weight)
        w_p,w_l,w_u=linalg.lu(q.astype(np.float32))
        w_s=np.diag(w_u)
        w_u=np.triu(w_u,1)
        u_mask=np.triu(np.ones_like(w_u),1)
        l_mask=u_mask.T

        self.register_buffer('w_p',torch.from_numpy(w_p))
        self.register_buffer('u_mask',torch.from_numpy(u_mask))
        self.register_buffer('l_mask',torch.from_numpy(l_mask))
        self.register_buffer('l_eye',torch.eye(l_mask.shape[0]))
        self.register_buffer('s_sign',torch.sign(torch.from_numpy(w_s)))
        self.w_l=torch.nn.Parameter(torch.from_numpy(w_l))
        self.w_s=torch.nn.Parameter(torch.log(1e-7+torch.abs(torch.from_numpy(w_s))))
        self.w_u=torch.nn.Parameter(torch.from_numpy(w_u))

        self.weight=None
        self.invweight=None

        return 
Example #5
Source File: test_twodim_base.py    From recruit with Apache License 2.0 6 votes vote down vote up
def test_tril_triu_ndim3():
    for dtype in np.typecodes['AllFloat'] + np.typecodes['AllInteger']:
        a = np.array([
            [[1, 1], [1, 1]],
            [[1, 1], [1, 0]],
            [[1, 1], [0, 0]],
            ], dtype=dtype)
        a_tril_desired = np.array([
            [[1, 0], [1, 1]],
            [[1, 0], [1, 0]],
            [[1, 0], [0, 0]],
            ], dtype=dtype)
        a_triu_desired = np.array([
            [[1, 1], [0, 1]],
            [[1, 1], [0, 0]],
            [[1, 1], [0, 0]],
            ], dtype=dtype)
        a_triu_observed = np.triu(a)
        a_tril_observed = np.tril(a)
        assert_array_equal(a_triu_observed, a_triu_desired)
        assert_array_equal(a_tril_observed, a_tril_desired)
        assert_equal(a_triu_observed.dtype, a.dtype)
        assert_equal(a_tril_observed.dtype, a.dtype) 
Example #6
Source File: test_twodim_base.py    From recruit with Apache License 2.0 6 votes vote down vote up
def test_tril_triu_dtype():
    # Issue 4916
    # tril and triu should return the same dtype as input
    for c in np.typecodes['All']:
        if c == 'V':
            continue
        arr = np.zeros((3, 3), dtype=c)
        assert_equal(np.triu(arr).dtype, arr.dtype)
        assert_equal(np.tril(arr).dtype, arr.dtype)

    # check special cases
    arr = np.array([['2001-01-01T12:00', '2002-02-03T13:56'],
                    ['2004-01-01T12:00', '2003-01-03T13:45']],
                   dtype='datetime64')
    assert_equal(np.triu(arr).dtype, arr.dtype)
    assert_equal(np.tril(arr).dtype, arr.dtype)

    arr = np.zeros((3,3), dtype='f4,f4')
    assert_equal(np.triu(arr).dtype, arr.dtype)
    assert_equal(np.tril(arr).dtype, arr.dtype) 
Example #7
Source File: test_tools.py    From tenpy with GNU General Public License v3.0 6 votes vote down vote up
def test_qr_li():
    cutoff = 1.e-10
    for shape in [(5, 4), (4, 5)]:
        print('shape =', shape)
        A = np.arange(20).reshape(shape)  # linearly dependent: only two rows/columns independent
        A[3, :] = np.random.random() * (cutoff / 100)  # nearly linear dependent
        q, r = tools.math.qr_li(A)
        assert np.linalg.norm(r - np.triu(r)) == 0.
        qdq = q.T.conj().dot(q)
        assert np.linalg.norm(qdq - np.eye(len(qdq))) < 1.e-13
        assert np.linalg.norm(q.dot(r) - A) < cutoff * 20
        r, q = tools.math.rq_li(A)
        assert np.linalg.norm(r - np.triu(r, r.shape[1] - r.shape[0])) == 0.
        qqd = q.dot(q.T.conj())
        assert np.linalg.norm(qqd - np.eye(len(qqd))) < 1.e-13
        assert np.linalg.norm(r.dot(q) - A) < cutoff * 20 
Example #8
Source File: test_twodim_base.py    From lambda-packs with MIT License 6 votes vote down vote up
def test_tril_triu_ndim3():
    for dtype in np.typecodes['AllFloat'] + np.typecodes['AllInteger']:
        a = np.array([
            [[1, 1], [1, 1]],
            [[1, 1], [1, 0]],
            [[1, 1], [0, 0]],
            ], dtype=dtype)
        a_tril_desired = np.array([
            [[1, 0], [1, 1]],
            [[1, 0], [1, 0]],
            [[1, 0], [0, 0]],
            ], dtype=dtype)
        a_triu_desired = np.array([
            [[1, 1], [0, 1]],
            [[1, 1], [0, 0]],
            [[1, 1], [0, 0]],
            ], dtype=dtype)
        a_triu_observed = np.triu(a)
        a_tril_observed = np.tril(a)
        yield assert_array_equal, a_triu_observed, a_triu_desired
        yield assert_array_equal, a_tril_observed, a_tril_desired
        yield assert_equal, a_triu_observed.dtype, a.dtype
        yield assert_equal, a_tril_observed.dtype, a.dtype 
Example #9
Source File: test_twodim_base.py    From lambda-packs with MIT License 6 votes vote down vote up
def test_tril_triu_dtype():
    # Issue 4916
    # tril and triu should return the same dtype as input
    for c in np.typecodes['All']:
        if c == 'V':
            continue
        arr = np.zeros((3, 3), dtype=c)
        assert_equal(np.triu(arr).dtype, arr.dtype)
        assert_equal(np.tril(arr).dtype, arr.dtype)

    # check special cases
    arr = np.array([['2001-01-01T12:00', '2002-02-03T13:56'],
                    ['2004-01-01T12:00', '2003-01-03T13:45']],
                   dtype='datetime64')
    assert_equal(np.triu(arr).dtype, arr.dtype)
    assert_equal(np.tril(arr).dtype, arr.dtype)

    arr = np.zeros((3,3), dtype='f4,f4')
    assert_equal(np.triu(arr).dtype, arr.dtype)
    assert_equal(np.tril(arr).dtype, arr.dtype) 
Example #10
Source File: transformer.py    From ITDD with MIT License 6 votes vote down vote up
def _get_attn_subsequent_mask(self, size):
        """
        Get an attention mask to avoid using the subsequent info.

        Args:
            size: int

        Returns:
            (`LongTensor`):

            * subsequent_mask `[1 x size x size]`
        """
        attn_shape = (1, size, size)
        subsequent_mask = np.triu(np.ones(attn_shape), k=1).astype('uint8')
        subsequent_mask = torch.from_numpy(subsequent_mask)
        return subsequent_mask 
Example #11
Source File: ktransformer.py    From ITDD with MIT License 6 votes vote down vote up
def _get_attn_subsequent_mask(self, size):
        """
        Get an attention mask to avoid using the subsequent info.

        Args:
            size: int

        Returns:
            (`LongTensor`):

            * subsequent_mask `[1 x size x size]`
        """
        attn_shape = (1, size, size)
        subsequent_mask = np.triu(np.ones(attn_shape), k=1).astype('uint8')
        subsequent_mask = torch.from_numpy(subsequent_mask)
        return subsequent_mask 
Example #12
Source File: mtransformer.py    From ITDD with MIT License 6 votes vote down vote up
def _get_attn_subsequent_mask(self, size):
        """
        Get an attention mask to avoid using the subsequent info.

        Args:
            size: int

        Returns:
            (`LongTensor`):

            * subsequent_mask `[1 x size x size]`
        """
        attn_shape = (1, size, size)
        subsequent_mask = np.triu(np.ones(attn_shape), k=1).astype('uint8')
        subsequent_mask = torch.from_numpy(subsequent_mask)
        return subsequent_mask 
Example #13
Source File: test_twodim_base.py    From vnpy_crypto with MIT License 6 votes vote down vote up
def test_tril_triu_ndim3():
    for dtype in np.typecodes['AllFloat'] + np.typecodes['AllInteger']:
        a = np.array([
            [[1, 1], [1, 1]],
            [[1, 1], [1, 0]],
            [[1, 1], [0, 0]],
            ], dtype=dtype)
        a_tril_desired = np.array([
            [[1, 0], [1, 1]],
            [[1, 0], [1, 0]],
            [[1, 0], [0, 0]],
            ], dtype=dtype)
        a_triu_desired = np.array([
            [[1, 1], [0, 1]],
            [[1, 1], [0, 0]],
            [[1, 1], [0, 0]],
            ], dtype=dtype)
        a_triu_observed = np.triu(a)
        a_tril_observed = np.tril(a)
        yield assert_array_equal, a_triu_observed, a_triu_desired
        yield assert_array_equal, a_tril_observed, a_tril_desired
        yield assert_equal, a_triu_observed.dtype, a.dtype
        yield assert_equal, a_tril_observed.dtype, a.dtype 
Example #14
Source File: test_twodim_base.py    From vnpy_crypto with MIT License 6 votes vote down vote up
def test_tril_triu_dtype():
    # Issue 4916
    # tril and triu should return the same dtype as input
    for c in np.typecodes['All']:
        if c == 'V':
            continue
        arr = np.zeros((3, 3), dtype=c)
        assert_equal(np.triu(arr).dtype, arr.dtype)
        assert_equal(np.tril(arr).dtype, arr.dtype)

    # check special cases
    arr = np.array([['2001-01-01T12:00', '2002-02-03T13:56'],
                    ['2004-01-01T12:00', '2003-01-03T13:45']],
                   dtype='datetime64')
    assert_equal(np.triu(arr).dtype, arr.dtype)
    assert_equal(np.tril(arr).dtype, arr.dtype)

    arr = np.zeros((3,3), dtype='f4,f4')
    assert_equal(np.triu(arr).dtype, arr.dtype)
    assert_equal(np.tril(arr).dtype, arr.dtype) 
Example #15
Source File: _bksf_test.py    From OpenFermion with Apache License 2.0 6 votes vote down vote up
def test_bravyi_kitaev_fast_edgeoperator_Bi(self):
        # checking the edge operators
        edge_matrix = numpy.triu(numpy.ones((4, 4)))
        edge_matrix_indices = numpy.array(
            numpy.nonzero(numpy.triu(edge_matrix) -
                          numpy.diag(numpy.diag(edge_matrix))))

        correct_operators_b0 = ((0, 'Z'), (1, 'Z'), (2, 'Z'))
        correct_operators_b1 = ((0, 'Z'), (3, 'Z'), (4, 'Z'))
        correct_operators_b2 = ((1, 'Z'), (3, 'Z'), (5, 'Z'))
        correct_operators_b3 = ((2, 'Z'), (4, 'Z'), (5, 'Z'))

        qterm_b0 = QubitOperator(correct_operators_b0, 1)
        qterm_b1 = QubitOperator(correct_operators_b1, 1)
        qterm_b2 = QubitOperator(correct_operators_b2, 1)
        qterm_b3 = QubitOperator(correct_operators_b3, 1)
        self.assertTrue(qterm_b0 ==
                        _bksf.edge_operator_b(edge_matrix_indices, 0))
        self.assertTrue(qterm_b1 ==
                        _bksf.edge_operator_b(edge_matrix_indices, 1))
        self.assertTrue(qterm_b2 ==
                        _bksf.edge_operator_b(edge_matrix_indices, 2))
        self.assertTrue(qterm_b3 ==
                        _bksf.edge_operator_b(edge_matrix_indices, 3)) 
Example #16
Source File: models.py    From fragile with MIT License 6 votes vote down vote up
def _cov_matrix_diagonalization(self):
        # Decomposition of cov_matrix into coords_matrix*diag(scaling_diag.^2)*coords_matrix'
        # (diagonalization)
        if (
            self._count_eval - self.n_eigen_eval
            > self.pop_size / (self.lr_covrank1_const + self.lr_mu_const) / self.n_dims / 10
        ):
            self.n_eigen_eval = self._count_eval
            self.cov_matrix = numpy.triu(self.cov_matrix) + numpy.triu(self.cov_matrix, 1).T
            eigvals, eigvects = numpy.linalg.eig(self.cov_matrix)
            self.scaling_diag = numpy.diag(eigvals)  # [::-1])
            self.coords_matrix = eigvects  # [:, ::-1]
            assert numpy.abs(numpy.imag(self.coords_matrix).sum()) == 0, self.coords_matrix
            assert numpy.abs(numpy.imag(self.scaling_diag).sum()) == 0, self.scaling_diag
            self.scaling_diag = numpy.sqrt(numpy.diag(self.scaling_diag)).reshape(-1, 1)
            self.invsqrtC = numpy.matmul(
                numpy.matmul(self.coords_matrix, numpy.diag(self.scaling_diag.flatten() ** -1)),
                self.coords_matrix.T,
            ) 
Example #17
Source File: test_twodim_base.py    From auto-alt-text-lambda-api with MIT License 6 votes vote down vote up
def test_tril_triu_ndim3():
    for dtype in np.typecodes['AllFloat'] + np.typecodes['AllInteger']:
        a = np.array([
            [[1, 1], [1, 1]],
            [[1, 1], [1, 0]],
            [[1, 1], [0, 0]],
            ], dtype=dtype)
        a_tril_desired = np.array([
            [[1, 0], [1, 1]],
            [[1, 0], [1, 0]],
            [[1, 0], [0, 0]],
            ], dtype=dtype)
        a_triu_desired = np.array([
            [[1, 1], [0, 1]],
            [[1, 1], [0, 0]],
            [[1, 1], [0, 0]],
            ], dtype=dtype)
        a_triu_observed = np.triu(a)
        a_tril_observed = np.tril(a)
        yield assert_array_equal, a_triu_observed, a_triu_desired
        yield assert_array_equal, a_tril_observed, a_tril_desired
        yield assert_equal, a_triu_observed.dtype, a.dtype
        yield assert_equal, a_tril_observed.dtype, a.dtype 
Example #18
Source File: test_twodim_base.py    From auto-alt-text-lambda-api with MIT License 6 votes vote down vote up
def test_tril_triu_dtype():
    # Issue 4916
    # tril and triu should return the same dtype as input
    for c in np.typecodes['All']:
        if c == 'V':
            continue
        arr = np.zeros((3, 3), dtype=c)
        assert_equal(np.triu(arr).dtype, arr.dtype)
        assert_equal(np.tril(arr).dtype, arr.dtype)

    # check special cases
    arr = np.array([['2001-01-01T12:00', '2002-02-03T13:56'],
                    ['2004-01-01T12:00', '2003-01-03T13:45']],
                   dtype='datetime64')
    assert_equal(np.triu(arr).dtype, arr.dtype)
    assert_equal(np.tril(arr).dtype, arr.dtype)

    arr = np.zeros((3,3), dtype='f4,f4')
    assert_equal(np.triu(arr).dtype, arr.dtype)
    assert_equal(np.tril(arr).dtype, arr.dtype) 
Example #19
Source File: slinalg.py    From D-VAE with MIT License 5 votes vote down vote up
def perform(self, node, inputs, outputs):
        """
        Implements the "reverse-mode" gradient [1]_ for the
        Cholesky factorization of a positive-definite matrix.

        References
        ----------
        .. [1] S. P. Smith. "Differentiation of the Cholesky Algorithm".
           Journal of Computational and Graphical Statistics,
           Vol. 4, No. 2 (Jun.,1995), pp. 134-147
           http://www.jstor.org/stable/1390762

        """
        x = inputs[0]
        L = inputs[1]
        dz = inputs[2]
        dx = outputs[0]
        N = x.shape[0]
        if self.lower:
            F = numpy.tril(dz)
            for k in xrange(N - 1, -1, -1):
                for j in xrange(k + 1, N):
                    for i in xrange(j, N):
                        F[i, k] -= F[i, j] * L[j, k]
                        F[j, k] -= F[i, j] * L[i, k]
                for j in xrange(k + 1, N):
                    F[j, k] /= L[k, k]
                    F[k, k] -= L[j, k] * F[j, k]
                F[k, k] /= (2 * L[k, k])
        else:
            F = numpy.triu(dz)
            for k in xrange(N - 1, -1, -1):
                for j in xrange(k + 1, N):
                    for i in xrange(j, N):
                        F[k, i] -= F[j, i] * L[k, j]
                        F[k, j] -= F[j, i] * L[k, i]
                for j in xrange(k + 1, N):
                    F[k, j] /= L[k, k]
                    F[k, k] -= L[k, j] * F[k, j]
                F[k, k] /= (2 * L[k, k])
        dx[0] = F 
Example #20
Source File: data.py    From controllable-text-attribute-transfer with Apache License 2.0 5 votes vote down vote up
def subsequent_mask(size):
    "Mask out subsequent positions."
    attn_shape = (1, size, size)
    subsequent_mask = np.triu(np.ones(attn_shape), k=1).astype('uint8')
    return torch.from_numpy(subsequent_mask) == 0 
Example #21
Source File: graphTools.py    From graph-neural-networks with GNU General Public License v3.0 5 votes vote down vote up
def edgeFailSampling(W, p):
    """
    edgeFailSampling: randomly delete the edges of a given graph
    
    Input:
        W (np.array): adjacency matrix
        p (float): probability of deleting an edge
    
    Output:
        W (np.array): adjacency matrix with some edges randomly deleted
        
    Obs.: The resulting graph need not be connected (even if the input graph is)
    """
    
    assert 0 <= p <= 1
    N = W.shape[0]
    assert W.shape[1] == N
    undirected = np.allclose(W, W.T, atol = zeroTolerance)
    
    maskEdges = np.random.rand(N, N)
    maskEdges = (maskEdges > p).astype(W.dtype) # Put a 1 with probability 1-p
    
    W = maskEdges * W
    if undirected:
        W = np.triu(W)
        W = W + W.T
        
    return W 
Example #22
Source File: slinalg.py    From D-VAE with MIT License 5 votes vote down vote up
def __init__(self, lower=True):
        assert lower in [True, False]
        self.lower = lower
        if lower:
            self.tri0 = numpy.tril
            self.tri1 = lambda a: numpy.triu(a, 1)
        else:
            self.tri0 = numpy.triu
            self.tri1 = lambda a: numpy.tril(a, -1) 
Example #23
Source File: model.py    From controllable-text-attribute-transfer with Apache License 2.0 5 votes vote down vote up
def subsequent_mask(size):
    "Mask out subsequent positions."
    attn_shape = (1, size, size)
    subsequent_mask = np.triu(np.ones(attn_shape), k=1).astype('uint8')
    return torch.from_numpy(subsequent_mask) == 0 
Example #24
Source File: nlinalg.py    From D-VAE with MIT License 5 votes vote down vote up
def perform(self, node, inputs, outputs):
        """
        Implements the "reverse-mode" gradient for the eigensystem of
        a square matrix.

        """
        x, w, v, W, V = inputs
        N = x.shape[0]
        outer = numpy.outer

        def G(n):
            return sum(v[:, m] * V.T[n].dot(v[:, m]) / (w[n] - w[m])
                       for m in xrange(N) if m != n)

        g = sum(outer(v[:, n], v[:, n] * W[n] + G(n))
                for n in xrange(N))

        # Numpy's eigh(a, 'L') (eigh(a, 'U')) is a function of tril(a)
        # (triu(a)) only.  This means that partial derivative of
        # eigh(a, 'L') (eigh(a, 'U')) with respect to a[i,j] is zero
        # for i < j (i > j).  At the same time, non-zero components of
        # the gradient must account for the fact that variation of the
        # opposite triangle contributes to variation of two elements
        # of Hermitian (symmetric) matrix. The following line
        # implements the necessary logic.
        out = self.tri0(g) + self.tri1(g).T

        # The call to self.tri0 in perform upcast from float32 to
        # float64 or from int* to int64 in numpy 1.6.1 but not in
        # 1.6.2. We do not want version dependent dtype in Theano.
        # We think it should be the same as the output.
        outputs[0][0] = numpy.asarray(out, dtype=node.outputs[0].dtype) 
Example #25
Source File: model.py    From controllable-text-attribute-transfer with Apache License 2.0 5 votes vote down vote up
def subsequent_mask(size):
    "Mask out subsequent positions."
    attn_shape = (1, size, size)
    subsequent_mask = np.triu(np.ones(attn_shape), k=1).astype('uint8')
    return torch.from_numpy(subsequent_mask) == 0 
Example #26
Source File: test_twodim_base.py    From vnpy_crypto with MIT License 5 votes vote down vote up
def test_mask_indices():
    # simple test without offset
    iu = mask_indices(3, np.triu)
    a = np.arange(9).reshape(3, 3)
    assert_array_equal(a[iu], array([0, 1, 2, 4, 5, 8]))
    # Now with an offset
    iu1 = mask_indices(3, np.triu, 1)
    assert_array_equal(a[iu1], array([1, 2, 5])) 
Example #27
Source File: test_linalg.py    From vnpy_crypto with MIT License 5 votes vote down vote up
def test_dynamic_programming_logic(self):
        # Test for the dynamic programming part
        # This test is directly taken from Cormen page 376.
        arrays = [np.random.random((30, 35)),
                  np.random.random((35, 15)),
                  np.random.random((15, 5)),
                  np.random.random((5, 10)),
                  np.random.random((10, 20)),
                  np.random.random((20, 25))]
        m_expected = np.array([[0., 15750., 7875., 9375., 11875., 15125.],
                               [0.,     0., 2625., 4375.,  7125., 10500.],
                               [0.,     0.,    0.,  750.,  2500.,  5375.],
                               [0.,     0.,    0.,    0.,  1000.,  3500.],
                               [0.,     0.,    0.,    0.,     0.,  5000.],
                               [0.,     0.,    0.,    0.,     0.,     0.]])
        s_expected = np.array([[0,  1,  1,  3,  3,  3],
                               [0,  0,  2,  3,  3,  3],
                               [0,  0,  0,  3,  3,  3],
                               [0,  0,  0,  0,  4,  5],
                               [0,  0,  0,  0,  0,  5],
                               [0,  0,  0,  0,  0,  0]], dtype=int)
        s_expected -= 1  # Cormen uses 1-based index, python does not.

        s, m = _multi_dot_matrix_chain_order(arrays, return_costs=True)

        # Only the upper triangular part (without the diagonal) is interesting.
        assert_almost_equal(np.triu(s[:-1, 1:]),
                            np.triu(s_expected[:-1, 1:]))
        assert_almost_equal(np.triu(m), np.triu(m_expected)) 
Example #28
Source File: bierman.py    From vnpy_crypto with MIT License 5 votes vote down vote up
def udu(M):
    """Construct the UDU' decomposition of a positive, semidefinite matrix M

    Parameters
    ----------
    M : [n, n] array
        Matrix to factorize

    Returns
    -------
    UDU : UDU_decomposition of size n
        UDU' representation of M
    """
    assert np.allclose(M, M.T), 'M must be symmetric, positive semidefinite'
    n = M.shape[0]

    # perform Bierman's COV2UD subroutine (fucking inclusive indices)
    M = np.triu(M)
    U = np.eye(n)
    d = np.zeros(n)
    for j in reversed(range(2, n + 1)):
        d[j - 1] = M[j - 1, j - 1]
        if d[j - 1] > 0:
            alpha = 1.0 / d[j - 1]
        else:
            alpha = 0.0
        for k in range(1, j):
            beta = M[k - 1, j - 1]
            U[k - 1, j - 1] = alpha * beta
            M[0:k, k - 1] = M[0:k, k - 1] - beta * U[0:k, j - 1]

    d[0] = M[0, 0]

    return UDU_decomposition(U, d) 
Example #29
Source File: test_sample.py    From strawberryfields with Apache License 2.0 5 votes vote down vote up
def test_invalid_adjacency(self, dim):
        """Test if function raises a ``ValueError`` for a matrix that is not symmetric"""
        with pytest.raises(ValueError, match="Input must be a NumPy array"):
            adj_asym = np.triu(np.ones((dim, dim)))
            sample.sample(A=adj_asym, n_mean=1.0) 
Example #30
Source File: test_twodim_base.py    From vnpy_crypto with MIT License 5 votes vote down vote up
def test_tril_triu_with_inf():
    # Issue 4859
    arr = np.array([[1, 1, np.inf],
                    [1, 1, 1],
                    [np.inf, 1, 1]])
    out_tril = np.array([[1, 0, 0],
                         [1, 1, 0],
                         [np.inf, 1, 1]])
    out_triu = out_tril.T
    assert_array_equal(np.triu(arr), out_triu)
    assert_array_equal(np.tril(arr), out_tril)