Line data Source code
1 : //----------------------------------------------------------------------------
2 : // LAGraph/experimental/test/LG_check_argminmax.c: simple arg min/max method
3 : // ----------------------------------------------------------------------------
4 :
5 : // LAGraph, (c) 2019-2025 by The LAGraph Contributors, All Rights Reserved.
6 : // SPDX-License-Identifier: BSD-2-Clause
7 : //
8 : // For additional details (including references to third party source code and
9 : // other files) see the LICENSE file or contact permission@sei.cmu.edu. See
10 : // Contributors.txt for a full list of contributors. Created, in part, with
11 : // funding and support from the U.S. Government (see Acknowledgments.txt file).
12 : // DM22-0790
13 :
14 : // Contributed by Tim Davis, Texas A&M University
15 :
16 : // ----------------------------------------------------------------------------
17 :
18 : #include "LG_Xtest.h"
19 : #include "LG_internal.h"
20 :
21 : #undef LG_FREE_WORK
22 : #define LG_FREE_WORK \
23 : { \
24 : LAGraph_Free ((void **) &I, msg) ; \
25 : LAGraph_Free ((void **) &J, msg) ; \
26 : LAGraph_Free ((void **) &X, msg) ; \
27 : LAGraph_Free ((void **) &v, msg) ; \
28 : LAGraph_Free ((void **) &p, msg) ; \
29 : }
30 :
31 : #undef LG_FREE_ALL
32 : #define LG_FREE_ALL \
33 : { \
34 : LG_FREE_WORK ; \
35 : GrB_free (x_result) ; \
36 : GrB_free (p_result) ; \
37 : }
38 :
39 : //------------------------------------------------------------------------------
40 : // LG_check_argminmax: compute argmin/max of each row/column of A
41 : //------------------------------------------------------------------------------
42 :
43 : // This method is single threaded and thus slow (like many LG_check_* methods).
44 : // For simplicity, this method does all its computations in FP64, and then
45 : // typecasts its results when done. This can cause loss of precision if the
46 : // input matrix is of type int64 or uint64 and has entries larger than about
47 : // 2^52. It can compute a different result for the position p if dim==0;
48 : // otherwise, its results should match LAGraph_argminmax.
49 :
50 265 : int LG_check_argminmax
51 : (
52 : // output
53 : GrB_Vector *x_result, // min/max value in each row/col of A
54 : GrB_Vector *p_result, // index of min/max value in each row/col of A
55 : // input
56 : GrB_Matrix A,
57 : int dim, // dim=1: cols of A, dim=2: rows of A
58 : // dim=0: return a scalar A(i,j) in x(0) and
59 : // its row and col in p.
60 : bool is_min,
61 : char *msg
62 : )
63 : {
64 :
65 : //--------------------------------------------------------------------------
66 : // check inputs
67 : //--------------------------------------------------------------------------
68 :
69 265 : LG_CLEAR_MSG ;
70 :
71 265 : uint64_t *I = NULL, *J = NULL, *p = NULL ;
72 265 : double *X = NULL, *v = NULL ;
73 :
74 265 : LG_ASSERT (A != NULL, GrB_NULL_POINTER) ;
75 265 : LG_ASSERT (x_result != NULL, GrB_NULL_POINTER) ;
76 265 : LG_ASSERT (p_result != NULL, GrB_NULL_POINTER) ;
77 :
78 265 : (*x_result) = NULL ;
79 265 : (*p_result) = NULL ;
80 :
81 : //--------------------------------------------------------------------------
82 : // extract the entries from the matrix
83 : //--------------------------------------------------------------------------
84 :
85 : GrB_Index nrows, ncols, nvals, n, np ;
86 265 : GRB_TRY (GrB_Matrix_nrows (&nrows, A)) ;
87 265 : GRB_TRY (GrB_Matrix_ncols (&ncols, A)) ;
88 265 : GRB_TRY (GrB_Matrix_nvals (&nvals, A)) ;
89 :
90 265 : switch (dim)
91 : {
92 88 : case 0: n = 1 ; np = 2 ; break ;
93 88 : case 1: n = ncols ; np = n ; break ;
94 88 : case 2: n = nrows ; np = n ; break ;
95 1 : default: LG_ASSERT (false, GrB_INVALID_VALUE) ;
96 : }
97 :
98 : //--------------------------------------------------------------------------
99 : // extract the entries from the matrix
100 : //--------------------------------------------------------------------------
101 :
102 264 : LG_TRY (LAGraph_Malloc ((void **) &I, nvals, sizeof (uint64_t), msg)) ;
103 264 : LG_TRY (LAGraph_Malloc ((void **) &J, nvals, sizeof (uint64_t), msg)) ;
104 264 : LG_TRY (LAGraph_Malloc ((void **) &X, nvals, sizeof (double), msg)) ;
105 :
106 264 : GRB_TRY (GrB_Matrix_extractTuples_FP64 (I, J, X, &nvals, A)) ;
107 :
108 264 : GrB_Type type = NULL ;
109 : #if LAGRAPH_SUITSPARSE
110 : GRB_TRY (GxB_Matrix_type (&type, A)) ;
111 : #else
112 : char typename [LAGRAPH_MAX_NAME_LEN+1] ;
113 264 : LG_TRY (LAGraph_Matrix_TypeName (typename, A, msg)) ;
114 264 : LG_TRY (LAGraph_TypeFromName (&type, typename, msg)) ;
115 : #endif
116 :
117 : //--------------------------------------------------------------------------
118 : // allocate temporary arrays for the results
119 : //--------------------------------------------------------------------------
120 :
121 264 : LG_TRY (LAGraph_Calloc ((void **) &v, n, sizeof (double), msg)) ;
122 264 : LG_TRY (LAGraph_Malloc ((void **) &p, np, sizeof (uint64_t), msg)) ;
123 93324 : for (int64_t k = 0 ; k < np ; k++)
124 : {
125 : // kth entry is not yet present
126 93060 : p [k] = -1 ;
127 : }
128 :
129 : //--------------------------------------------------------------------------
130 : // compute the argmin/argmax
131 : //--------------------------------------------------------------------------
132 :
133 264 : if (dim == 0)
134 : {
135 :
136 : //----------------------------------------------------------------------
137 : // global argmin/argmax
138 : //----------------------------------------------------------------------
139 :
140 88 : if (nvals > 0)
141 : {
142 : // (v,p) = first entry in the matrix
143 88 : v [0] = X [0] ;
144 88 : p [0] = I [0] ;
145 88 : p [1] = J [0] ;
146 :
147 1855590 : for (int64_t k = 1 ; k < nvals ; k++)
148 : {
149 1855502 : if (is_min)
150 : {
151 : // argmin
152 927751 : if (X [k] < v [0])
153 : {
154 98 : v [0] = X [k] ;
155 98 : p [0] = I [k] ;
156 98 : p [1] = J [k] ;
157 : }
158 : }
159 : else
160 : {
161 : // argmax
162 927751 : if (X [k] > v [0])
163 : {
164 101 : v [0] = X [k] ;
165 101 : p [0] = I [k] ;
166 101 : p [1] = J [k] ;
167 : }
168 : }
169 : }
170 :
171 : }
172 :
173 : }
174 176 : else if (dim == 1)
175 : {
176 :
177 : //----------------------------------------------------------------------
178 : // argmin/argmax of each column of A
179 : //----------------------------------------------------------------------
180 :
181 1855678 : for (int64_t k = 0 ; k < nvals ; k++)
182 : {
183 1855590 : int64_t i = I [k] ;
184 1855590 : int64_t j = J [k] ;
185 1855590 : double x = X [k] ;
186 1855590 : if (p [j] == -1)
187 : {
188 : // first entry seen in column j
189 46442 : v [j] = x ;
190 46442 : p [j] = i ;
191 : }
192 1809148 : else if (is_min)
193 : {
194 904574 : if (x < v [j])
195 : {
196 30132 : v [j] = x ;
197 30132 : p [j] = i ;
198 : }
199 : }
200 : else
201 : {
202 904574 : if (x > v [j])
203 : {
204 44540 : v [j] = x ;
205 44540 : p [j] = i ;
206 : }
207 : }
208 : }
209 :
210 : }
211 : else // (dim == 2)
212 : {
213 :
214 : //----------------------------------------------------------------------
215 : // argmin/argmax of each row of A
216 : //----------------------------------------------------------------------
217 :
218 1855678 : for (int64_t k = 0 ; k < nvals ; k++)
219 : {
220 1855590 : int64_t i = I [k] ;
221 1855590 : int64_t j = J [k] ;
222 1855590 : double x = X [k] ;
223 1855590 : if (p [i] == -1)
224 : {
225 : // first entry seen in row i
226 46442 : v [i] = x ;
227 46442 : p [i] = j ;
228 : }
229 1809148 : else if (is_min)
230 : {
231 904574 : if (x < v [i])
232 : {
233 30064 : v [i] = x ;
234 30064 : p [i] = j ;
235 : }
236 : }
237 : else
238 : {
239 904574 : if (x > v [i])
240 : {
241 44176 : v [i] = x ;
242 44176 : p [i] = j ;
243 : }
244 : }
245 : }
246 : }
247 :
248 : //--------------------------------------------------------------------------
249 : // create the output vectors and fill with the results
250 : //--------------------------------------------------------------------------
251 :
252 264 : LG_TRY (GrB_Vector_new (x_result, type, n)) ;
253 264 : LG_TRY (GrB_Vector_new (p_result, GrB_INT64, np)) ;
254 :
255 264 : if (dim == 0)
256 : {
257 88 : if (nvals > 0)
258 : {
259 88 : GRB_TRY (GrB_Vector_setElement_FP64 (*x_result, v [0], 0)) ;
260 88 : GRB_TRY (GrB_Vector_setElement_INT64 (*p_result, p [0], 0)) ;
261 88 : GRB_TRY (GrB_Vector_setElement_INT64 (*p_result, p [1], 1)) ;
262 : }
263 : }
264 : else
265 : {
266 93060 : for (int64_t k = 0 ; k < n ; k++)
267 : {
268 92884 : if (p [k] != -1)
269 : {
270 92884 : GRB_TRY (GrB_Vector_setElement_FP64 (*x_result, v [k], k)) ;
271 92884 : GRB_TRY (GrB_Vector_setElement_INT64 (*p_result, p [k], k)) ;
272 : }
273 : }
274 : }
275 :
276 264 : GRB_TRY (GrB_wait (*x_result, GrB_MATERIALIZE)) ;
277 264 : GRB_TRY (GrB_wait (*p_result, GrB_MATERIALIZE)) ;
278 :
279 : //--------------------------------------------------------------------------
280 : // free workspace and return result
281 : //--------------------------------------------------------------------------
282 :
283 264 : LG_FREE_WORK ;
284 264 : return (GrB_SUCCESS) ;
285 : }
286 :
|
