Python sparse.block_diag方法代码示例

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def cellGradBC(self):
        """
        The cell centered Gradient boundary condition matrix
        """
        if getattr(self, '_cellGradBC', None) is None:
            BC = self.setCellGradBC(self._cellGradBC_list)
            n = self.vnC
            if self.dim == 1:
                G = ddxCellGradBC(n[0], BC[0])
            elif self.dim == 2:
                G1 = sp.kron(speye(n[1]), ddxCellGradBC(n[0], BC[0]))
                G2 = sp.kron(ddxCellGradBC(n[1], BC[1]), speye(n[0]))
                G = sp.block_diag((G1, G2), format="csr")
            elif self.dim == 3:
                G1 = kron3(speye(n[2]), speye(n[1]), ddxCellGradBC(n[0], BC[0]))
                G2 = kron3(speye(n[2]), ddxCellGradBC(n[1], BC[1]), speye(n[0]))
                G3 = kron3(ddxCellGradBC(n[2], BC[2]), speye(n[1]), speye(n[0]))
                G = sp.block_diag((G1, G2, G3), format="csr")
            # Compute areas of cell faces & volumes
            S = self.area
            V = self.aveCC2F*self.vol  # Average volume between adjacent cells
            self._cellGradBC = sdiag(S/V)*G
        return self._cellGradBC 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def aveF2CCV(self):
        """
        averaging operator of x-faces (radial) to cell centered vectors

        Returns
        -------
        scipy.sparse.csr_matrix
            matrix that averages from faces to cell centered vectors
        """
        if getattr(self, '_aveF2CCV', None) is None:
            # n = self.vnC
            if self.isSymmetric is True:
                self._aveF2CCV = sp.block_diag(
                    (self.aveFx2CC, self.aveFz2CC), format="csr"
                )
            else:
                self._aveF2CCV = sp.block_diag(
                    (self.aveFx2CC, self.aveFy2CC, self.aveFz2CC),
                    format="csr"
                )
        return self._aveF2CCV

    ####################################################
    # Deflation Matrices
    #################################################### 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def compute_dr_wrt(self, wrt):
        if wrt is not self.a:
            return None
    
        Ainv = self.r

        if Ainv.ndim <= 2:
            return -np.kron(Ainv, Ainv.T)
        else:
            Ainv = np.reshape(Ainv,  (-1, Ainv.shape[-2], Ainv.shape[-1]))
            AinvT = np.rollaxis(Ainv, -1, -2)
            AinvT = np.reshape(AinvT, (-1, AinvT.shape[-2], AinvT.shape[-1]))
            result = np.dstack([-np.kron(Ainv[i], AinvT[i]).T for i in range(Ainv.shape[0])]).T
            result = sp.block_diag(result)

        return result 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def batch_of_relational_supports_to_support(batched_relational_supports,
                                            ordering=None):
    if ordering is not None:
        relational_supports_combined = OrderedDict([
            (k, sparse.block_diag([s[k] for s in batched_relational_supports]))
            for k in ordering])
    else:
        first_dict = batched_relational_supports[0]

        if not isinstance(first_dict, OrderedDict):
            ValueError("If you do not provide an ordering, then "
                       "`batched_relational_supports` shguld be list of "
                       "`OrderedDict`. It is a "
                       "list of {}".format(type(first_dict)))

        first_dict_order = list(first_dict.keys())

        relational_supports_combined = OrderedDict([
            (k, sparse.block_diag([s[k] for s in batched_relational_supports]))
            for k in first_dict_order])

    return relational_supports_to_support(relational_supports_combined,
                                          ordering=ordering) 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def _lazy_build_elements(self):
        self._elements = []

        for vel in self.vector_elements:
            vstart = 0
            if self.sparse:  # build block-diagonal sparse mx
                diagBlks = []
                for compbasis in self.component_bases:
                    cs = compbasis.elshape
                    comp_vel = vel[vstart:vstart + compbasis.dim]
                    diagBlks.append(comp_vel.reshape(cs))
                    vstart += compbasis.dim
                el = _sps.block_diag(diagBlks, "csr", 'complex')

            else:
                start = [0] * self.elndim
                el = _np.zeros(self.elshape, 'complex')
                for compbasis in self.component_bases:
                    cs = compbasis.elshape
                    comp_vel = vel[vstart:vstart + compbasis.dim]
                    slc = tuple([slice(start[k], start[k] + cs[k]) for k in range(self.elndim)])
                    el[slc] = comp_vel.reshape(cs)
                    vstart += compbasis.dim
                    for k in range(self.elndim): start[k] += cs[k]

            self._elements.append(el)
        if not self.sparse:  # _elements should be an array rather than a list
            self._elements = _np.array(self._elements) 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def collate_fn(batch):
        graphs, pmpds, labels = zip(*batch)
        batched_graphs = graph_batch(graphs)
        batched_pmpds = sp.block_diag(pmpds)
        batched_labels = np.concatenate(labels, axis=0)
        return batched_graphs, batched_pmpds, batched_labels 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def numpy_to_disjoint(X_list, A_list, E_list=None):
    """
    Converts a batch of graphs stored in lists (X, A, and optionally E) to the
    [disjoint mode](https://danielegrattarola.github.io/spektral/data/#disjoint-mode).

    Each entry i of the lists should be associated to the same graph, i.e.,
    `X_list[i].shape[0] == A_list[i].shape[0] == E_list[i].shape[0]`.

    :param X_list: a list of np.arrays of shape `(N, F)`;
    :param A_list: a list of np.arrays or sparse matrices of shape `(N, N)`;
    :param E_list: a list of np.arrays of shape `(N, N, S)`;
    :return:
        -  `X_out`: a rank 2 array of shape `(n_nodes, F)`;
        -  `A_out`: a rank 2 array of shape `(n_nodes, n_nodes)`;
        -  `E_out`: (only if `E_list` is given) a rank 2 array of shape
        `(n_edges, S)`;
    """
    X_out = np.vstack(X_list)
    A_list = [sp.coo_matrix(a) for a in A_list]
    if E_list is not None:
        if E_list[0].ndim == 3:
            E_list = [e[a.row, a.col] for e, a in zip(E_list, A_list)]
        E_out = np.vstack(E_list)
    A_out = sp.block_diag(A_list)
    n_nodes = np.array([x.shape[0] for x in X_list])
    I_out = np.repeat(np.arange(len(n_nodes)), n_nodes)
    if E_list is not None:
        return X_out, A_out, E_out, I_out
    else:
        return X_out, A_out, I_out 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def _test_disjoint_mode(layer, **kwargs):
    A = sp.block_diag([np.ones((N1, N1)), np.ones(
        (N2, N2)), np.ones((N3, N3))]).todense()
    X = np.random.normal(size=(N, F))
    I = np.array([0] * N1 + [1] * N2 + [2] * N3).astype(int)
    sparse = kwargs.pop('sparse', None) is not None

    A_in = Input(shape=(None,), sparse=sparse)
    X_in = Input(shape=(F,))
    I_in = Input(shape=(), dtype=tf.int32)

    layer_instance = layer(**kwargs)
    output = layer_instance([X_in, A_in, I_in])
    model = Model([X_in, A_in, I_in], output)
    output = model([X, A, I])
    X_pool, A_pool, I_pool, mask = output

    N_pool_expected = np.ceil(kwargs['ratio'] * N1) + \
                      np.ceil(kwargs['ratio'] * N2) + \
                      np.ceil(kwargs['ratio'] * N3)
    N_pool_expected = int(N_pool_expected)
    N_pool_true = A_pool.shape[0]

    _check_number_of_nodes(N_pool_expected, N_pool_true)

    assert X_pool.shape == (N_pool_expected, F)
    assert A_pool.shape == (N_pool_expected, N_pool_expected)
    assert I_pool.shape == (N_pool_expected,)

    output_shape = [o.shape for o in output]
    _check_output_and_model_output_shapes(output_shape, model.output_shape) 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def sparse_factor_solve_kkt(Q_tilde, D_tilde, A_, C_tilde, rx, rs, rz, ry, ns):
    nineq, nz, neq, nBatch = ns

    # H_ = csc_matrix((nz + nineq, nz + nineq))
    # H_[:nz, :nz] = Q_tilde
    # H_[-nineq:, -nineq:] = D_tilde
    H_ = block_diag([Q_tilde, D_tilde], format='csc')
    if neq > 0:
        g_ = torch.cat([rx, rs], 1).squeeze(0).numpy()
        h_ = torch.cat([rz, ry], 1).squeeze(0).numpy()
    else:
        g_ = torch.cat([rx, rs], 1).squeeze(0).numpy()
        h_ = rz.squeeze(0).numpy()

    H_LU = splu(H_)

    invH_A_ = csc_matrix(H_LU.solve(A_.todense().transpose()))
    invH_g_ = H_LU.solve(g_)

    S_ = A_.dot(invH_A_) + C_tilde
    S_LU = splu(S_)
    # t_ = invH_g_[np.newaxis].dot(A_.transpose()).squeeze(0) - h_
    t_ = A_.dot(invH_g_) - h_
    w_ = -S_LU.solve(t_)
    # t_ = -g_ - w_[np.newaxis].dot(A_).squeeze(0)
    t_ = -g_ - A_.transpose().dot(w_)
    v_ = H_LU.solve(t_)

    dx = v_[:nz]
    ds = v_[nz:]
    dz = w_[:nineq]
    dy = w_[nineq:] if neq > 0 else None

    dx = torch.DoubleTensor(dx).unsqueeze(0)
    ds = torch.DoubleTensor(ds).unsqueeze(0)
    dz = torch.DoubleTensor(dz).unsqueeze(0)
    dy = torch.DoubleTensor(dy).unsqueeze(0) if neq > 0 else None

    return dx, ds, dz, dy 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def compute_d1(self):
        # To stay consistent with numpy, we must upgrade 1D arrays to 2D
        mtx1r = row(self.mtx1.r) if len(self.mtx1.r.shape)<2 else self.mtx1.r
        mtx2r = col(self.mtx2.r) if len(self.mtx2.r.shape)<2 else self.mtx2.r

        if mtx1r.ndim <= 2:
            return sp.kron(sp.eye(mtx1r.shape[0], mtx1r.shape[0]),mtx2r.T)
        else:
            mtx2f = mtx2r.reshape((-1, mtx2r.shape[-2], mtx2r.shape[-1]))
            mtx2f = np.rollaxis(mtx2f, -1, -2) #transpose basically            
            result = sp.block_diag([np.kron(np.eye(mtx1r.shape[-2], mtx1r.shape[-2]),m2) for m2 in mtx2f])
            assert(result.shape[0] == self.r.size)
            return result 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def compute_d2(self):
        
        # To stay consistent with numpy, we must upgrade 1D arrays to 2D
        mtx1r = row(self.mtx1.r) if len(self.mtx1.r.shape)<2 else self.mtx1.r
        mtx2r = col(self.mtx2.r) if len(self.mtx2.r.shape)<2 else self.mtx2.r

        if mtx2r.ndim <= 1:
            return self.mtx1r
        elif mtx2r.ndim <= 2:
            return sp.kron(mtx1r, sp.eye(mtx2r.shape[1],mtx2r.shape[1]))
        else:
            mtx1f = mtx1r.reshape((-1, mtx1r.shape[-2], mtx1r.shape[-1]))            
            result = sp.block_diag([np.kron(m1, np.eye(mtx2r.shape[-1],mtx2r.shape[-1])) for m1 in mtx1f])
            assert(result.shape[0] == self.r.size)
            return result 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def duplicate_boundary_data(shape, axis: int) -> sparse.csc.csc_matrix:
    '''
    duplicate_boundary_data creates a sparse matrix that duplicated the vectorized data at
    a boundary along the axis.
    For example: if you have a space A that is 4x5 represented by vector with 20 elements. Then
        duplicate_boundary_data(x_in, y_in, axis: 0 ) will generate a 25x20 matrix that will make
        a vector representation of a 5x5 matrix where the first row is duplicated as last row.
    (used to implement periodicity)
    '''
    if axis == 0:
        block_a = sparse.eye(shape[0] - 1).tocsr()
        block_b = sparse.csr_matrix((1, shape[0] - 1))
        block_b[0, 0] = 1
        block_matrix = sparse.vstack([block_a, block_b])
        per = sparse.block_diag(tuple(shape[1] * [block_matrix]))
        block_c = sparse.csr_matrix((shape[0] - 1, 1))
        block_matrix = sparse.hstack([block_a, block_c])
        per_reverse = sparse.block_diag(tuple(shape[1] * [block_matrix]))
    elif axis == 1:
        block_1 = sparse.eye(int(shape[0] * (shape[1] - 1)))
        block_b = sparse.eye(shape[0]).tocsr()
        block_2 = sparse.hstack([
            block_b,
            sparse.csr_matrix((shape[0], shape[0] * int(shape[1] - 2)))
        ])
        per = sparse.vstack([block_1, block_2])
        block_1 = sparse.eye(int(shape[0] * (shape[1] - 1)))
        block_2 = sparse.csr_matrix((shape[0] * int(shape[1] - 1), shape[0]))
        per_reverse = sparse.hstack([block_1, block_2])
    return per, per_reverse


# Geometry matrix functions
############################ 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def symmetry_matrix2D_sign_correction(x_vector: np.array,
                                      y_vector: np.array,
                                      axis: int,
                                      sign: int = 1,
                                      periodicity: bool = False):
    '''
    SymmetryMatrix2D_sign_correction generates a matrixon3 that corrects the sign of every
        parameter component (f, fx, fy, fxy) according to the parametrization being symmetric
        or anti-symmetric.

        The correct symmetry_matrix for [f, fx, fy, fxy] will be of the form:
            sign_matrix@block_diag(4*symmetry matrix)

    Args:
     x_vector: x vector of the fine grid
     y_vector: y vector of the fine grid
     axis: Symmetry axis.
     sign: 1 or -1 depanding whether or not you have symmetry or anti-symmetry.
     periodicity: Indicating whether or not there is periodicity.

    Return:
     sign_correction matrix


    '''
    n_x = len(x_vector)
    n_y = len(y_vector)
    shape = (n_x, n_y)
    return symmetry_matrix_sign_correction(shape, axis, sign, periodicity) 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def aveF2CCV(self):
        "Construct the averaging operator on cell faces to cell centers."
        if getattr(self, '_aveF2CCV', None) is None:
            if self.dim == 1:
                self._aveF2CCV = self.aveFx2CC
            elif self.dim == 2:
                self._aveF2CCV = sp.block_diag((
                    self.aveFx2CC, self.aveFy2CC
                ), format="csr")
            elif self.dim == 3:
                self._aveF2CCV = sp.block_diag((
                    self.aveFx2CC, self.aveFy2CC, self.aveFz2CC
                ), format="csr")
        return self._aveF2CCV 

# 需要导入模块: from scipy import sparse [as 别名]
# 或者: from scipy.sparse import block_diag [as 别名]
def aveCCV2F(self):
        """
        Construct the averaging operator on cell centers to
        faces as a vector.
        """
        if getattr(self, '_aveCCV2F', None) is None:
            if self.dim == 1:
                self._aveCCV2F = self.aveCC2F
            elif self.dim == 2:
                aveCCV2Fx = sp.kron(speye(self.nCy), av_extrap(self.nCx))
                aveCC2VFy = sp.kron(av_extrap(self.nCy), speye(self.nCx))
                self._aveCCV2F = sp.block_diag((
                    aveCCV2Fx, aveCC2VFy
                ), format="csr")
            elif self.dim == 3:
                aveCCV2Fx = kron3(
                    speye(self.nCz), speye(self.nCy), av_extrap(self.nCx)
                )
                aveCC2VFy = kron3(
                    speye(self.nCz), av_extrap(self.nCy), speye(self.nCx)
                )
                aveCC2BFz = kron3(
                    av_extrap(self.nCz), speye(self.nCy), speye(self.nCx)
                )
                self._aveCCV2F = sp.block_diag((
                        aveCCV2Fx, aveCC2VFy, aveCC2BFz
                ), format="csr")
        return self._aveCCV2F