Line data Source code
1 :
2 : //------------------------------------------------------------------------------
3 : // LAGraph_Fast_Build: Uses saxpy methods for faster builds, especially powerful
4 : // when output is bitmap.
5 : //------------------------------------------------------------------------------
6 :
7 : // LAGraph, (c) 2019-2024 by The LAGraph Contributors, All Rights Reserved.
8 : // SPDX-License-Identifier: BSD-2-Clause
9 : //
10 : // For additional details (including references to third party source code and
11 : // other files) see the LICENSE file or contact permission@sei.cmu.edu. See
12 : // Contributors.txt for a full list of contributors. Created, in part, with
13 : // funding and support from the U.S. Government (see Acknowledgments.txt file).
14 : // DM22-0790
15 :
16 : // Contributed by Gabriel Gomez, Texas A&M University
17 :
18 : //------------------------------------------------------------------------------
19 :
20 : // This method is allows for assigns or build into a vector. It is created to
21 : // target some specific mxv GraphBLAS methods that occasionally tend to be
22 : // faster than both build and assign.
23 :
24 : // This method is used for cases where you want to build a vector equivilent to
25 : // the following for loop:
26 : // for (j = 0 ; j < n ; j++)
27 : // {
28 : // uint64_t i = I_vec [j] ;
29 : // c [i] += X_vec [j] ;
30 : // }
31 :
32 : // It is fastest when the accum biop is equivalent to the dup moniod and c is a
33 : // full vector or when there is no accumulator and I_vec.nvals is
34 : // sufficiently large when compared to c.nrows.
35 :
36 : // It builds a matix P which is coustructed in O(1) time if a
37 : // filled ramp vector is passed in.
38 :
39 : // P is built by column and will contain one entry per column placed at row
40 : // I_vec [j]. So, P*X_vec will find for every row the intersections with X_vec,
41 : // apply the dup monoid and put the result with the accumulator back into c at
42 : // that row.
43 :
44 : // X_vec can be sparse if dup is a semiring (dup_second), in which case the jth
45 : // column correspoding to an empty X_vec[j] have no effect on the outcome.
46 :
47 : // desc can also be modified to affect the mask or to pass in GRB_TRAN for
48 : // GrB_INP0, this has some interesting usecases, specifically when X_vec is a
49 : // "map". That is:
50 : // LAGraph_FastAssign(c, NULL, NULL, I_vec, map, NULL, any_second, NULL, msg)
51 : // is equivalent to:
52 : // for (j = 0 ; j < n ; j++)
53 : // {
54 : // c [j] = map[I_vec[j]] ;
55 : // }
56 : // Currently, this is not implemented for LAGraph_FastAssign_Monoid so a
57 : // semiring must be passed in.
58 :
59 :
60 : #include "LG_internal.h"
61 : #include "LAGraphX.h"
62 : #if LG_SUITESPARSE_GRAPHBLAS_V10
63 :
64 : // Uncomment if you would like to use the monoid version of FastAssign.
65 : // Passing in a semiring is faster but this may be more convienient.
66 : #if 0
67 : #undef LG_FREE_ALL
68 : #define LG_FREE_ALL \
69 : { \
70 : GrB_free(&sem); \
71 : }
72 :
73 : int LAGraph_FastAssign_Monoid
74 : (
75 : // output
76 : // Vector to be built (or assigned): initialized with correct dimensions.
77 : GrB_Vector c,
78 : // inputs
79 : const GrB_Vector mask,
80 : const GrB_BinaryOp accum,
81 : const GrB_Vector I_vec, // Indecies (duplicates allowed)
82 : const GrB_Vector X_vec, // Values
83 : // Optional (Give me a ramp with size > X_vec.size for faster calculations)
84 : const GrB_Vector ramp,
85 : const GrB_Monoid dup, // Applied to duplicates
86 : const GrB_Descriptor desc,
87 : char *msg
88 : )
89 : {
90 : GrB_BinaryOp op = NULL;
91 : GrB_Semiring sem = NULL;
92 : int code = 0;
93 : GrB_get(X_vec, &code, GrB_EL_TYPE_CODE);
94 : switch (code)
95 : {
96 : case GrB_BOOL_CODE : op = GrB_SECOND_BOOL ; break ;
97 : case GrB_INT8_CODE : op = GrB_SECOND_INT8 ; break ;
98 : case GrB_INT16_CODE : op = GrB_SECOND_INT16 ; break ;
99 : case GrB_INT32_CODE : op = GrB_SECOND_INT32 ; break ;
100 : case GrB_INT64_CODE : op = GrB_SECOND_INT64 ; break ;
101 : case GrB_UINT8_CODE : op = GrB_SECOND_UINT8 ; break ;
102 : case GrB_UINT16_CODE : op = GrB_SECOND_UINT16 ; break ;
103 : case GrB_UINT32_CODE : op = GrB_SECOND_UINT32 ; break ;
104 : case GrB_UINT64_CODE : op = GrB_SECOND_UINT64 ; break ;
105 : case GrB_FP32_CODE : op = GrB_SECOND_FP32 ; break ;
106 : case GrB_FP64_CODE : op = GrB_SECOND_FP64 ; break ;
107 : case GxB_FC32_CODE : op = GxB_SECOND_FC32 ; break ;
108 : case GxB_FC64_CODE : op = GxB_SECOND_FC64 ; break ;
109 : default :
110 : LG_ERROR_MSG("LAGraph failed (file %s, line %d):" \
111 : " LAGraph_FastAssign_Monoid not implemented for UDTs",
112 : __FILE__, __LINE__);
113 : break ;
114 : }
115 : GRB_TRY (GrB_Semiring_new(&sem, dup, op)) ;
116 : // GxB_print(sem, stdout) ;
117 : LG_TRY (LAGraph_FastAssign_Semiring
118 : (c, mask, accum, I_vec, X_vec, ramp, sem, desc, msg)) ;
119 : LG_FREE_ALL ;
120 : return (GrB_SUCCESS);
121 : }
122 : #endif
123 :
124 : #undef LG_FREE_ALL
125 : #define LG_FREE_ALL \
126 : { \
127 : GrB_free(&P); \
128 : GrB_free(&con); \
129 : LAGraph_Free(&ramp_a, msg); \
130 : LAGraph_Free(&i_a, msg); \
131 : }
132 :
133 : // This method can be faster if given a builtin semiring.
134 105584 : int LAGraph_FastAssign_Semiring
135 : (
136 : // output
137 : // Vector to be built (or assigned): initialized with correct dimensions.
138 : GrB_Vector c,
139 : // inputs
140 : const GrB_Vector mask,
141 : const GrB_BinaryOp accum,
142 : const GrB_Vector I_vec, // Indecies (duplicates allowed)
143 : const GrB_Vector X_vec, // Values
144 : // Optional (Give me a ramp with size > X_vec.size for faster calculations)
145 : const GrB_Vector ramp,
146 : // monoid is applied to duplicates. Binary op should be SECOND.
147 : const GrB_Semiring semiring,
148 : const GrB_Descriptor desc,
149 : char *msg
150 : )
151 : {
152 : // TODO: Let indicies be specified by value of by index in I_vec via
153 : // descriptor (although build or assign would be better if I_vec is by
154 : // index since it is sorted and has no dups). By value could be useful if
155 : // I_vec is not full.
156 : // TODO: Ditto for X_vec
157 :
158 105584 : GrB_Matrix P = NULL;
159 : GrB_Index n, nrows;
160 105584 : GxB_Container con = NULL;
161 105584 : void *ramp_a = NULL, *i_a =NULL;
162 105584 : int ramp_h = GrB_DEFAULT, trsp = GrB_DEFAULT, i_h = GrB_DEFAULT;
163 105584 : uint64_t ramp_n = 0, ramp_size = 0, i_n = 0, i_size= 0;
164 105584 : GrB_Type i_type = NULL, ramp_type = NULL;
165 : //----------------------------------------------------------------------
166 : // Check inputs
167 : //----------------------------------------------------------------------
168 : //TODO: assert inputs are full or desc says to use by value etc.
169 105584 : LG_ASSERT (c != NULL, GrB_NULL_POINTER) ;
170 105584 : LG_ASSERT (I_vec != NULL, GrB_NULL_POINTER) ;
171 105584 : LG_ASSERT (X_vec != NULL, GrB_NULL_POINTER) ;
172 105584 : LG_ASSERT_MSG (c != X_vec && c != I_vec && c != ramp && c != mask,
173 : GrB_NOT_IMPLEMENTED, "c cannot be aliased with any input.") ;
174 :
175 : //----------------------------------------------------------------------
176 : // Find dimensions and type
177 : //----------------------------------------------------------------------
178 105584 : GRB_TRY (GrB_Vector_size(&n, I_vec)) ;
179 :
180 105584 : if(desc != NULL)
181 : {
182 20089 : GRB_TRY (GrB_get(desc, &trsp, GrB_INP0)) ;
183 : }
184 105584 : if(trsp == GrB_TRAN)
185 : {
186 20089 : GRB_TRY (GrB_Vector_size(&nrows, X_vec)) ;
187 : }
188 : else
189 : {
190 85495 : GRB_TRY (GrB_Vector_size(&nrows, c)) ;
191 : }
192 :
193 : //----------------------------------------------------------------------
194 : // Load up containers
195 : //----------------------------------------------------------------------
196 105584 : GRB_TRY (GrB_Matrix_new(&P, GrB_BOOL, nrows, n));
197 105134 : GRB_TRY (GxB_Container_new(&con));
198 :
199 104234 : if(ramp == NULL)
200 : {
201 : //FUTURE: maybe let user input a size 0 ramp and build it for them?
202 310 : GRB_TRY (GrB_free(&(con->p))) ;
203 310 : ramp_type = (n + 1 <= INT32_MAX)? GrB_UINT32: GrB_UINT64;
204 620 : GrB_IndexUnaryOp idxnum = (n + 1 <= INT32_MAX)?
205 310 : GrB_ROWINDEX_INT32: GrB_ROWINDEX_INT64;
206 310 : GRB_TRY (GrB_Vector_new(&(con->p), ramp_type, n + 1));
207 294 : GRB_TRY (GrB_assign (con->p, NULL, NULL, 0, GrB_ALL, 0, NULL)) ;
208 286 : GRB_TRY (GrB_apply (con->p, NULL, NULL, idxnum, con->p, 0, NULL)) ;
209 : }
210 : else
211 : {
212 103924 : GRB_TRY (GxB_Vector_unload(
213 : ramp, &ramp_a, &ramp_type, &ramp_n, &ramp_size, &ramp_h, NULL)) ;
214 103924 : LG_ASSERT_MSG (ramp_n > n, GrB_DIMENSION_MISMATCH, "Ramp too small!");
215 103924 : GRB_TRY (GxB_Vector_load(
216 : con->p, &ramp_a, ramp_type, n + 1, ramp_size,
217 : GxB_IS_READONLY, NULL)) ;
218 : // Since con->p won't free this array I should be safe to load it back
219 : // into ramp.
220 103924 : GRB_TRY (GxB_Vector_load(
221 : ramp, &ramp_a, ramp_type, ramp_n, ramp_size, ramp_h, NULL)) ;
222 103924 : ramp_a = NULL;
223 : }
224 : // con->i = I_vec;
225 104194 : GRB_TRY (GxB_Vector_unload(
226 : I_vec, &i_a, &i_type, &i_n, &i_size, &i_h, NULL)) ;
227 104194 : GRB_TRY (GxB_Vector_load(
228 : con->i, &i_a, i_type, i_n, i_size, GxB_IS_READONLY, NULL)) ;
229 : // Since con->i won't free this array I should be safe to load it back
230 : // into I_vec.
231 104194 : GRB_TRY (GxB_Vector_load(
232 : I_vec, &i_a, i_type, i_n, i_size, i_h, NULL)) ;
233 104194 : i_a = NULL;
234 :
235 : // con->x [0] = false, of length 1
236 104194 : GRB_TRY (GrB_free(&(con->x))) ;
237 104194 : GRB_TRY (GrB_Vector_new (&(con->x), GrB_BOOL, 1)) ;
238 103894 : GRB_TRY (GrB_assign (con->x, NULL, NULL, 0, GrB_ALL, 1, NULL)) ;
239 102844 : con->format = GxB_SPARSE;
240 102844 : con->orientation = GrB_COLMAJOR;
241 102844 : con->nrows = nrows;
242 102844 : con->ncols = n ;
243 102844 : con->nvals = n ;
244 102844 : con->nrows_nonempty = -1 ;
245 102844 : con->ncols_nonempty = n ;
246 102844 : con->iso = true ;
247 102844 : con->jumbled = false ;
248 102844 : con->format = GxB_SPARSE ;
249 102844 : con->orientation = GrB_COLMAJOR ;
250 102844 : con->Y = NULL ;
251 : //----------------------------------------------------------------------
252 : // Load P and do the mxv
253 : //----------------------------------------------------------------------
254 102844 : GRB_TRY (GxB_load_Matrix_from_Container(P, con, NULL));
255 : // GRB_TRY (GxB_fprint(P, GxB_SHORT, stdout));
256 102844 : GRB_TRY (GrB_mxv(c, mask, accum, semiring, P, X_vec, desc));
257 : //----------------------------------------------------------------------
258 : // Free work.
259 : // Note: this does not free inputs since they are marked GxB_IS_READONLY
260 : //----------------------------------------------------------------------
261 102244 : GrB_free(&P) ;
262 102244 : GrB_free(&con) ;
263 102244 : return (GrB_SUCCESS) ;
264 : }
265 : #endif
|
