(A, Cnodes, AggOp, T)
return {'P_I': P_I, 'I_F': I_F, 'I_C': I_C, 'Cpts': Cpts, 'Fpts': Fpts}
Return C and F pts.

    Helper function that returns a dictionary of sparse matrices and arrays
    which allow us to easily operate on Cpts and Fpts separately.

    Parameters
    ----------
    A : {csr_matrix, bsr_matrix}
        Operator
    Cnodes : {array}
        Array of all root node indices.  This is an array of nodal indices,
        not degree-of-freedom indices.  If the blocksize of T is 1, then
        nodal indices and degree-of-freedom indices coincide.
    AggOp : {csr_matrix}
        Aggregation operator corresponding to A
    T : {bsr_matrix}
        Tentative prolongator based on AggOp

    Returns
    -------
    Dictionary containing these parameters:

    P_I : {bsr_matrix}
        Interpolation operator that carries out only simple injection from the
        coarse grid to fine grid Cpts nodes
    I_F : {bsr_matrix}
        Identity operator on Fpts, i.e., the action of this matrix zeros
        out entries in a vector at all Cpts, leaving Fpts untouched
    I_C : {bsr_matrix}
        Identity operator on Cpts nodes, i.e., the action of this matrix zeros
        out entries in a vector at all Fpts, leaving Cpts untouched
    Cpts : {array}
        An array of all root node dofs, corresponding to the F/C splitting
    Fpts : {array}
        An array of all non root node dofs, corresponding to the F/C splitting

    Examples
    --------
    >>> from numpy import array
    >>> from pyamg.util.utils import get_Cpt_params
    >>> from pyamg.gallery import poisson
    >>> from scipy.sparse import csr_matrix, bsr_matrix
    >>> A = poisson((10,), format='csr')
    >>> Cpts = array([3, 7])
    >>> AggOp = ([[ 1., 0.], [ 1., 0.],
    ...           [ 1., 0.], [ 1., 0.],
    ...           [ 1., 0.], [ 0., 1.],
    ...           [ 0., 1.], [ 0., 1.],
    ...           [ 0., 1.], [ 0., 1.]])
    >>> AggOp = csr_matrix(AggOp)
    >>> T = AggOp.copy().tobsr()
    >>> params = get_Cpt_params(A, Cpts, AggOp, T)
    >>> params['P_I'].toarray()
    array([[0., 0.],
           [0., 0.],
           [0., 0.],
           [1., 0.],
           [0., 0.],
           [0., 0.],
           [0., 0.],
           [0., 1.],
           [0., 0.],
           [0., 0.]])

    Notes
    -----
    The principal calling routine is
    aggregation.smooth.energy_prolongation_smoother,
    which uses the Cpt_param dictionary for root-node style
    prolongation smoothing
Return C and F pts.
Return C and F pts .
def get_Cpt_params(A, Cnodes, AggOp, T):
    """Return C and F pts.

    Helper function that returns a dictionary of sparse matrices and arrays
    which allow us to easily operate on Cpts and Fpts separately.

    Parameters
    ----------
    A : {csr_matrix, bsr_matrix}
        Operator
    Cnodes : {array}
        Array of all root node indices.  This is an array of nodal indices,
        not degree-of-freedom indices.  If the blocksize of T is 1, then
        nodal indices and degree-of-freedom indices coincide.
    AggOp : {csr_matrix}
        Aggregation operator corresponding to A
    T : {bsr_matrix}
        Tentative prolongator based on AggOp

    Returns
    -------
    Dictionary containing these parameters:

    P_I : {bsr_matrix}
        Interpolation operator that carries out only simple injection from the
        coarse grid to fine grid Cpts nodes
    I_F : {bsr_matrix}
        Identity operator on Fpts, i.e., the action of this matrix zeros
        out entries in a vector at all Cpts, leaving Fpts untouched
    I_C : {bsr_matrix}
        Identity operator on Cpts nodes, i.e., the action of this matrix zeros
        out entries in a vector at all Fpts, leaving Cpts untouched
    Cpts : {array}
        An array of all root node dofs, corresponding to the F/C splitting
    Fpts : {array}
        An array of all non root node dofs, corresponding to the F/C splitting

    Examples
    --------
    >>> from numpy import array
    >>> from pyamg.util.utils import get_Cpt_params
    >>> from pyamg.gallery import poisson
    >>> from scipy.sparse import csr_matrix, bsr_matrix
    >>> A = poisson((10,), format='csr')
    >>> Cpts = array([3, 7])
    >>> AggOp = ([[ 1., 0.], [ 1., 0.],
    ...           [ 1., 0.], [ 1., 0.],
    ...           [ 1., 0.], [ 0., 1.],
    ...           [ 0., 1.], [ 0., 1.],
    ...           [ 0., 1.], [ 0., 1.]])
    >>> AggOp = csr_matrix(AggOp)
    >>> T = AggOp.copy().tobsr()
    >>> params = get_Cpt_params(A, Cpts, AggOp, T)
    >>> params['P_I'].toarray()
    array([[0., 0.],
           [0., 0.],
           [0., 0.],
           [1., 0.],
           [0., 0.],
           [0., 0.],
           [0., 0.],
           [0., 1.],
           [0., 0.],
           [0., 0.]])

    Notes
    -----
    The principal calling routine is
    aggregation.smooth.energy_prolongation_smoother,
    which uses the Cpt_param dictionary for root-node style
    prolongation smoothing

    """
    if not isspmatrix_bsr(A) and not isspmatrix_csr(A):
        raise TypeError('Expected BSR or CSR matrix A')
    if not isspmatrix_csr(AggOp):
        raise TypeError('Expected CSR matrix AggOp')
    if not isspmatrix_bsr(T):
        raise TypeError('Expected BSR matrix T')
    if T.blocksize[0] != T.blocksize[1]:
        raise TypeError('Expected square blocksize for BSR matrix T')
    if A.shape[0] != A.shape[1]:
        raise TypeError('Expected square matrix A')
    if T.shape[0] != A.shape[0]:
        raise TypeError('Expected compatible dimensions for T and A,\
                         T.shape[0] = A.shape[0]')
    if Cnodes.shape[0] != AggOp.shape[1]:
        if AggOp.shape[1] > 1:
            raise TypeError('Number of columns in AggOp must equal number\
                             of Cnodes')

    if isspmatrix_bsr(A) and A.blocksize[0] > 1:
        # Expand the list of Cpt nodes to a list of Cpt dofs
        blocksize = A.blocksize[0]
        Cpts = np.repeat(blocksize*Cnodes, blocksize)
        for k in range(1, blocksize):
            Cpts[list(range(k, Cpts.shape[0], blocksize))] += k
    else:
        blocksize = 1
        Cpts = Cnodes
    Cpts = np.array(Cpts, dtype=int)

    # More input checking
    if Cpts.shape[0] != T.shape[1]:
        if T.shape[1] > blocksize:
            raise ValueError('Expected number of Cpts to match T.shape[1]')
    if blocksize != T.blocksize[0]:
        raise ValueError('Expected identical blocksize in A and T')
    if AggOp.shape[0] != int(T.shape[0]/blocksize):
        raise ValueError('Number of rows in AggOp must equal number of\
                          fine-grid nodes')

    # Create two maps, one for F points and one for C points
    ncoarse = T.shape[1]
    I_C = eye(A.shape[0], A.shape[1], format='csr')
    I_F = I_C.copy()
    I_F.data[Cpts] = 0.0
    I_F.eliminate_zeros()
    I_C = I_C - I_F
    I_C.eliminate_zeros()

    # Find Fpts, the complement of Cpts
    Fpts = I_F.indices.copy()

    # P_I only injects from Cpts on the coarse grid to the fine grid, but
    # because of it's later uses, it must have the CSC indices ordered as
    # in Cpts
    if I_C.nnz > 0:
        indices = Cpts.copy()
        indptr = np.arange(indices.shape[0]+1)
    else:
        indices = np.zeros((0,), dtype=T.indices.dtype)
        indptr = np.zeros((ncoarse+1,), dtype=T.indptr.dtype)

    P_I = csc_matrix((I_C.data.copy(), indices, indptr),
                     shape=(I_C.shape[0], ncoarse))
    P_I = P_I.tobsr(T.blocksize)

    # Use same blocksize as A
    if isspmatrix_bsr(A):
        I_C = I_C.tobsr(A.blocksize)
        I_F = I_F.tobsr(A.blocksize)
    else:
        I_C = I_C.tobsr(blocksize=(1, 1))
        I_F = I_F.tobsr(blocksize=(1, 1))

    return {'P_I': P_I, 'I_F': I_F, 'I_C': I_C, 'Cpts': Cpts, 'Fpts': Fpts}
def get_Cpt_params ( A , Cnodes , AggOp , T ) : if not isspmatrix_bsr ( A ) and not isspmatrix_csr ( A ) : raise TypeError ( 'Expected BSR or CSR matrix A' ) if not isspmatrix_csr ( AggOp ) : raise TypeError ( 'Expected CSR matrix AggOp' ) if not isspmatrix_bsr ( T ) : raise TypeError ( 'Expected BSR matrix T' ) if T . blocksize [ 0 ] != T . blocksize [ 1 ] : raise TypeError ( 'Expected square blocksize for BSR matrix T' ) if A . shape [ 0 ] != A . shape [ 1 ] : raise TypeError ( 'Expected square matrix A' ) if T . shape [ 0 ] != A . shape [ 0 ] : raise TypeError ( 'Expected compatible dimensions for T and A,\
                         T.shape[0] = A.shape[0]' ) if Cnodes . shape [ 0 ] != AggOp . shape [ 1 ] : if AggOp . shape [ 1 ] > 1 : raise TypeError ( 'Number of columns in AggOp must equal number\
                             of Cnodes' ) if isspmatrix_bsr ( A ) and A . blocksize [ 0 ] > 1 : # Expand the list of Cpt nodes to a list of Cpt dofs blocksize = A . blocksize [ 0 ] Cpts = np . repeat ( blocksize * Cnodes , blocksize ) for k in range ( 1 , blocksize ) : Cpts [ list ( range ( k , Cpts . shape [ 0 ] , blocksize ) ) ] += k else : blocksize = 1 Cpts = Cnodes Cpts = np . array ( Cpts , dtype = int ) # More input checking if Cpts . shape [ 0 ] != T . shape [ 1 ] : if T . shape [ 1 ] > blocksize : raise ValueError ( 'Expected number of Cpts to match T.shape[1]' ) if blocksize != T . blocksize [ 0 ] : raise ValueError ( 'Expected identical blocksize in A and T' ) if AggOp . shape [ 0 ] != int ( T . shape [ 0 ] / blocksize ) : raise ValueError ( 'Number of rows in AggOp must equal number of\
                          fine-grid nodes' ) # Create two maps, one for F points and one for C points ncoarse = T . shape [ 1 ] I_C = eye ( A . shape [ 0 ] , A . shape [ 1 ] , format = 'csr' ) I_F = I_C . copy ( ) I_F . data [ Cpts ] = 0.0 I_F . eliminate_zeros ( ) I_C = I_C - I_F I_C . eliminate_zeros ( ) # Find Fpts, the complement of Cpts Fpts = I_F . indices . copy ( ) # P_I only injects from Cpts on the coarse grid to the fine grid, but # because of it's later uses, it must have the CSC indices ordered as # in Cpts if I_C . nnz > 0 : indices = Cpts . copy ( ) indptr = np . arange ( indices . shape [ 0 ] + 1 ) else : indices = np . zeros ( ( 0 , ) , dtype = T . indices . dtype ) indptr = np . zeros ( ( ncoarse + 1 , ) , dtype = T . indptr . dtype ) P_I = csc_matrix ( ( I_C . data . copy ( ) , indices , indptr ) , shape = ( I_C . shape [ 0 ] , ncoarse ) ) P_I = P_I . tobsr ( T . blocksize ) # Use same blocksize as A if isspmatrix_bsr ( A ) : I_C = I_C . tobsr ( A . blocksize ) I_F = I_F . tobsr ( A . blocksize ) else : I_C = I_C . tobsr ( blocksize = ( 1 , 1 ) ) I_F = I_F . tobsr ( blocksize = ( 1 , 1 ) ) return { 'P_I' : P_I , 'I_F' : I_F , 'I_C' : I_C , 'Cpts' : Cpts , 'Fpts' : Fpts }