File: GB_extract.c

package info (click to toggle)
suitesparse 1%3A7.10.1%2Bdfsg-1
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid, trixie
  • size: 254,920 kB
  • sloc: ansic: 1,134,743; cpp: 46,133; makefile: 4,875; fortran: 2,087; java: 1,826; sh: 996; ruby: 725; python: 495; asm: 371; sed: 166; awk: 44
file content (193 lines) | stat: -rw-r--r-- 7,126 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
//------------------------------------------------------------------------------
// GB_extract: C<M> = accum(C,A(I,J))
//------------------------------------------------------------------------------

// SuiteSparse:GraphBLAS, Timothy A. Davis, (c) 2017-2025, All Rights Reserved.
// SPDX-License-Identifier: Apache-2.0

//------------------------------------------------------------------------------

// Not user-callable.  Implements the user-callable GrB_*_extract functions.
//
//      C<M> = accum (C, A (Rows,Cols)) or
//      C<M> = accum (C, AT(Rows,Cols)) where AT = A'
//
// equivalently:
//
//      C<M> = accum (C, A(Rows,Cols) )
//      C<M> = accum (C, A(Cols,Rows)')

#define GB_FREE_ALL GrB_Matrix_free (&T) ;

#include "extract/GB_extract.h"
#include "extract/GB_subref.h"
#include "mask/GB_accum_mask.h"

GrB_Info GB_extract                 // C<M> = accum (C, A(I,J))
(
    GrB_Matrix C,                   // input/output matrix for results
    const bool C_replace,           // C matrix descriptor
    const GrB_Matrix M,             // optional mask for C, unused if NULL
    const bool Mask_comp,           // mask descriptor
    const bool Mask_struct,         // if true, use the only structure of M
    const GrB_BinaryOp accum,       // optional accum for Z=accum(C,T)
    const GrB_Matrix A,             // input matrix
    const bool A_transpose,         // A matrix descriptor
    const void *Rows,               // row indices
    const bool Rows_is_32,          // if true, Rows is 32-bit; else 64-bit
    const uint64_t nRows_in,        // number of row indices
    const void *Cols,               // column indices
    const bool Cols_is_32,          // if true, Rows is 32-bit; else 64-bit
    const uint64_t nCols_in,        // number of column indices
    GB_Werk Werk
)
{

    //--------------------------------------------------------------------------
    // check inputs
    //--------------------------------------------------------------------------

    // C may be aliased with M and/or A

    GrB_Info info ;
    struct GB_Matrix_opaque T_header ;
    GrB_Matrix T = NULL ;
    GB_RETURN_IF_NULL (Rows) ;
    GB_RETURN_IF_NULL (Cols) ;
    GB_RETURN_IF_FAULTY_OR_POSITIONAL (accum) ;

    ASSERT_MATRIX_OK (C, "C input for GB_Matrix_extract", GB0) ;
    ASSERT_MATRIX_OK_OR_NULL (M, "M for GB_Matrix_extract", GB0) ;
    ASSERT_BINARYOP_OK_OR_NULL (accum, "accum for GB_Matrix_extract", GB0) ;
    ASSERT_MATRIX_OK (A, "A input for GB_Matrix_extract", GB0) ;

    // check domains and dimensions for C<M> = accum (C,T)
    GB_OK (GB_compatible (C->type, C, M, Mask_struct, accum, A->type,
        Werk)) ;

    // check the dimensions of C
    int64_t cnrows = GB_NROWS (C) ;
    int64_t cncols = GB_NCOLS (C) ;

    int64_t nRows, nCols, RowColon [3], ColColon [3] ;
    int rkind, ckind ;

    if (!A_transpose)
    { 
        // T = A(Rows,Cols)
        GB_ijlength (Rows, Rows_is_32, nRows_in, GB_NROWS (A), &nRows, &rkind,
            RowColon) ;
        GB_ijlength (Cols, Cols_is_32, nCols_in, GB_NCOLS (A), &nCols, &ckind,
            ColColon) ;
    }
    else
    { 
        // T = A(Cols,Rows)
        GB_ijlength (Rows, Rows_is_32, nRows_in, GB_NCOLS (A), &nRows, &rkind,
            RowColon) ;
        GB_ijlength (Cols, Cols_is_32, nCols_in, GB_NROWS (A), &nCols, &ckind,
            ColColon) ;
    }

    if (cnrows != nRows || cncols != nCols)
    { 
        GB_ERROR (GrB_DIMENSION_MISMATCH,
            "Dimensions not compatible:\n"
            "required size of output is " GBd "-by-" GBd "\n"
            "but actual size output is  " GBd "-by-" GBd "\n",
            nRows, nCols, cnrows, cncols) ;
    }

    // quick return if an empty mask is complemented
    GB_RETURN_IF_QUICK_MASK (C, C_replace, M, Mask_comp, Mask_struct) ;

    // delete any lingering zombies and assemble any pending tuples
    GB_MATRIX_WAIT (M) ;        // cannot be jumbled
    GB_MATRIX_WAIT (A) ;        // cannot be jumbled

    GB_BURBLE_DENSE (C, "(C %s) ") ;
    GB_BURBLE_DENSE (M, "(M %s) ") ;
    GB_BURBLE_DENSE (A, "(A %s) ") ;

    //--------------------------------------------------------------------------
    // handle the CSR/CSC format and transpose; T = A (I,J) or T = A (J,I)
    //--------------------------------------------------------------------------

    const void *I, *J ;
    int64_t ni, nj ;
    bool T_is_csc ;
    bool I_is_32, J_is_32 ;

    if (A->is_csc)
    {
        if (!A_transpose)
        { 
            // T = A(Rows,Cols) where both A and T are in CSC format
            I = Rows ; ni = nRows_in ;  // indices into the vectors
            J = Cols ; nj = nCols_in ;  // vectors
            T_is_csc = true ;           // return T in CSC format
            I_is_32 = Rows_is_32 ;
            J_is_32 = Cols_is_32 ;
        }
        else
        { 
            // T = A(Cols,Rows) where A is CSC and T is returned as CSR format
            I = Cols ; ni = nCols_in ;  // indices into the vectors
            J = Rows ; nj = nRows_in ;  // vectors
            T_is_csc = false ;          // return T in CSR format
            I_is_32 = Cols_is_32 ;
            J_is_32 = Rows_is_32 ;
        }
    }
    else
    {
        if (!A_transpose)
        { 
            // T = A(Rows,Cols) where both A and T are in CSR format
            I = Cols ; ni = nCols_in ;  // indices into the vectors
            J = Rows ; nj = nRows_in ;  // vectors
            T_is_csc = false ;          // return T in CSR format
            I_is_32 = Cols_is_32 ;
            J_is_32 = Rows_is_32 ;
        }
        else
        { 
            // T = A(Cols,Rows) where A is CSR but T is returned as CSC format
            I = Rows ; ni = nRows_in ;  // indices into the vectors
            J = Cols ; nj = nCols_in ;  // vectors
            T_is_csc = true ;           // return T in CSC format
            I_is_32 = Rows_is_32 ;
            J_is_32 = Cols_is_32 ;
        }
    }

    // T has T->vdim = |J|, each vector of length T->vlen = |J|, regardless of
    // its CSR/CSC format.

    // J is a list of length |J| of vectors in the range 0:A->vdim-1
    // I is a list of length |I| of indices in the range 0:A->vlen-1

    // |I| and |J| are either nRows or nCols, depending on the 4 cases above.

    // T has the same hypersparsity as A.

    //--------------------------------------------------------------------------
    // T = A (I,J)
    //--------------------------------------------------------------------------

    // TODO::: iso:  if accum is PAIR, extract T as iso

    GB_CLEAR_MATRIX_HEADER (T, &T_header) ;
    GB_OK (GB_subref (T, false, T_is_csc, A,
        I, I_is_32, ni, J, J_is_32, nj, false, Werk)) ;
    ASSERT_MATRIX_OK (T, "T extracted", GB0) ;
    ASSERT (GB_JUMBLED_OK (T)) ;

    //--------------------------------------------------------------------------
    // C<M> = accum (C,T): accumulate the results into C via the mask M
    //--------------------------------------------------------------------------

    return (GB_accum_mask (C, M, NULL, accum, &T, C_replace, Mask_comp,
        Mask_struct, Werk)) ;
}