1 //------------------------------------------------------------------------------
2 // GB_I_inverse: invert an index list
3 //------------------------------------------------------------------------------
4
5 // SuiteSparse:GraphBLAS, Timothy A. Davis, (c) 2017-2021, All Rights Reserved.
6 // SPDX-License-Identifier: Apache-2.0
7
8 //------------------------------------------------------------------------------
9
10 // I is a large list relative to the vector length, avlen, and it is not
11 // contiguous. Scatter I into the I inverse buckets (Mark and Inext) for quick
12 // lookup.
13
14 // FUTURE:: this code is sequential. Constructing the I inverse buckets in
15 // parallel would require synchronization (a critical section for each bucket,
16 // or atomics). A more parallel approach might use qsort first, to find
17 // duplicates in I, and then construct the buckets in parallel after the qsort.
18 // But the time complexity would be higher.
19
20 #include "GB_subref.h"
21
GB_I_inverse(const GrB_Index * I,int64_t nI,int64_t avlen,int64_t * restrict * p_Mark,size_t * p_Mark_size,int64_t * restrict * p_Inext,size_t * p_Inext_size,int64_t * p_ndupl,GB_Context Context)22 GrB_Info GB_I_inverse // invert the I list for C=A(I,:)
23 (
24 const GrB_Index *I, // list of indices, duplicates OK
25 int64_t nI, // length of I
26 int64_t avlen, // length of the vectors of A
27 // outputs:
28 int64_t *restrict *p_Mark, // head pointers for buckets, size avlen
29 size_t *p_Mark_size,
30 int64_t *restrict *p_Inext, // next pointers for buckets, size nI
31 size_t *p_Inext_size,
32 int64_t *p_ndupl, // number of duplicate entries in I
33 GB_Context Context
34 )
35 {
36
37 //--------------------------------------------------------------------------
38 // get inputs
39 //--------------------------------------------------------------------------
40
41 int64_t *Mark = NULL ; size_t Mark_size = 0 ;
42 int64_t *Inext = NULL ; size_t Inext_size = 0 ;
43 int64_t ndupl = 0 ;
44
45 (*p_Mark ) = NULL ; (*p_Mark_size ) = 0 ;
46 (*p_Inext) = NULL ; (*p_Inext_size) = 0 ;
47 (*p_ndupl) = 0 ;
48
49 //--------------------------------------------------------------------------
50 // allocate workspace
51 //--------------------------------------------------------------------------
52
53 Mark = GB_CALLOC_WERK (avlen, int64_t, &Mark_size) ;
54 Inext = GB_MALLOC_WERK (nI, int64_t, &Inext_size) ;
55 if (Inext == NULL || Mark == NULL)
56 {
57 // out of memory
58 GB_FREE_WERK (&Mark, Mark_size) ;
59 GB_FREE_WERK (&Inext, Inext_size) ;
60 return (GrB_OUT_OF_MEMORY) ;
61 }
62
63 //--------------------------------------------------------------------------
64 // scatter the I indices into buckets
65 //--------------------------------------------------------------------------
66
67 // at this point, Mark is all zero, so Mark [i] < 1 for all i in
68 // the range 0 to avlen-1.
69
70 // O(nI) time; not parallel
71 for (int64_t inew = nI-1 ; inew >= 0 ; inew--)
72 {
73 int64_t i = I [inew] ;
74 ASSERT (i >= 0 && i < avlen) ;
75 int64_t ihead = (Mark [i] - 1) ;
76 if (ihead < 0)
77 {
78 // first time i has been seen in the list I
79 ihead = -1 ;
80 }
81 else
82 {
83 // i has already been seen in the list I
84 ndupl++ ;
85 }
86 Mark [i] = inew + 1 ; // (Mark [i] - 1) = inew
87 Inext [inew] = ihead ;
88 }
89
90 // indices in I are now in buckets. An index i might appear more than once
91 // in the list I. inew = (Mark [i] - 1) is the first position of i in I (i
92 // will be I [inew]), (Mark [i] - 1) is the head of a link list of all
93 // places where i appears in I. inew = Inext [inew] traverses this list,
94 // until inew is -1.
95
96 // to traverse all entries in bucket i, do:
97 // GB_for_each_index_in_bucket (inew,i)) { ... }
98
99 #define GB_for_each_index_in_bucket(inew,i) \
100 for (int64_t inew = Mark [i] - 1 ; inew >= 0 ; inew = Inext [inew])
101
102 // If Mark [i] < 1, then the ith bucket is empty and i is not in I.
103 // Otherise, the first index in bucket i is (Mark [i] - 1).
104
105 #ifdef GB_DEBUG
106 for (int64_t i = 0 ; i < avlen ; i++)
107 {
108 GB_for_each_index_in_bucket (inew, i)
109 {
110 ASSERT (inew >= 0 && inew < nI) ;
111 ASSERT (i == I [inew]) ;
112 }
113 }
114 #endif
115
116 //--------------------------------------------------------------------------
117 // return result
118 //--------------------------------------------------------------------------
119
120 (*p_Mark ) = Mark ; (*p_Mark_size ) = Mark_size ;
121 (*p_Inext) = Inext ; (*p_Inext_size) = Inext_size ;
122 (*p_ndupl) = ndupl ;
123 return (GrB_SUCCESS) ;
124 }
125
126