1 //------------------------------------------------------------------------------
2 // GB_concat_hyper: concatenate an array of matrices into a hypersparse matrix
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 #define GB_FREE_ALL \
11 GB_FREE (&Wi, Wi_size) ; \
12 GB_FREE_WERK (&Wj, Wj_size) ; \
13 GB_FREE_WERK (&Wx, Wx_size) ; \
14 GB_phbix_free (C) ;
15
16 #include "GB_concat.h"
17
GB_concat_hyper(GrB_Matrix C,const int64_t cnz,const GrB_Matrix * Tiles,const GrB_Index m,const GrB_Index n,const int64_t * restrict Tile_rows,const int64_t * restrict Tile_cols,GB_Context Context)18 GrB_Info GB_concat_hyper // concatenate into a hypersparse matrix
19 (
20 GrB_Matrix C, // input/output matrix for results
21 const int64_t cnz, // # of entries in C
22 const GrB_Matrix *Tiles, // 2D row-major array of size m-by-n,
23 const GrB_Index m,
24 const GrB_Index n,
25 const int64_t *restrict Tile_rows, // size m+1
26 const int64_t *restrict Tile_cols, // size n+1
27 GB_Context Context
28 )
29 {
30
31 //--------------------------------------------------------------------------
32 // allocate triplet workspace to construct C as hypersparse
33 //--------------------------------------------------------------------------
34
35 GrB_Info info ;
36 GrB_Matrix A = NULL ;
37
38 int64_t *restrict Wi = NULL ; size_t Wi_size = 0 ;
39 int64_t *restrict Wj = NULL ; size_t Wj_size = 0 ;
40 GB_void *restrict Wx = NULL ; size_t Wx_size = 0 ;
41
42 GrB_Type ctype = C->type ;
43 int64_t cvlen = C->vlen ;
44 int64_t cvdim = C->vdim ;
45 bool csc = C->is_csc ;
46 size_t csize = ctype->size ;
47 GB_Type_code ccode = ctype->code ;
48
49 float hyper_switch = C->hyper_switch ;
50 float bitmap_switch = C->bitmap_switch ;
51 int sparsity_control = C->sparsity ;
52
53 bool static_header = C->static_header ;
54 GB_phbix_free (C) ;
55
56 Wi = GB_MALLOC (cnz, int64_t, &Wi_size) ; // becomes C->i
57 Wj = GB_MALLOC_WERK (cnz, int64_t, &Wj_size) ; // freed below
58 Wx = GB_MALLOC_WERK (cnz * csize, GB_void, &Wx_size) ; // freed below
59 if (Wi == NULL || Wj == NULL || Wx == NULL)
60 {
61 // out of memory
62 GB_FREE_ALL ;
63 return (GrB_OUT_OF_MEMORY) ;
64 }
65
66 GB_GET_NTHREADS_MAX (nthreads_max, chunk, Context) ;
67
68 int64_t nouter = csc ? n : m ;
69 int64_t ninner = csc ? m : n ;
70
71 //--------------------------------------------------------------------------
72 // concatenate all matrices into the list of triplets
73 //--------------------------------------------------------------------------
74
75 int64_t pC = 0 ;
76 for (int64_t outer = 0 ; outer < nouter ; outer++)
77 {
78 for (int64_t inner = 0 ; inner < ninner ; inner++)
79 {
80
81 //------------------------------------------------------------------
82 // get the tile A
83 //------------------------------------------------------------------
84
85 A = csc ? GB_TILE (Tiles, inner, outer)
86 : GB_TILE (Tiles, outer, inner) ;
87 ASSERT (!GB_ANY_PENDING_WORK (A)) ;
88
89 //------------------------------------------------------------------
90 // determine where to place the tile in C
91 //------------------------------------------------------------------
92
93 // The tile A appears in vectors cvstart:cvend-1 of C, and indices
94 // cistart:ciend-1.
95
96 int64_t cvstart, cistart ;
97 if (csc)
98 {
99 // C is held by column
100 // Tiles is row-major and accessed in column order
101 cvstart = Tile_cols [outer] ;
102 cistart = Tile_rows [inner] ;
103 }
104 else
105 {
106 // C is held by row
107 // Tiles is row-major and accessed in row order
108 cvstart = Tile_rows [outer] ;
109 cistart = Tile_cols [inner] ;
110 }
111
112 //------------------------------------------------------------------
113 // extract the tuples from tile A
114 //------------------------------------------------------------------
115
116 int64_t anz = GB_NNZ (A) ;
117 GB_OK (GB_extractTuples (
118 (GrB_Index *) ((csc ? Wi : Wj) + pC),
119 (GrB_Index *) ((csc ? Wj : Wi) + pC),
120 Wx + pC * csize, (GrB_Index *) (&anz), ccode, A, Context)) ;
121
122 //------------------------------------------------------------------
123 // adjust the indices to reflect their new place in C
124 //------------------------------------------------------------------
125
126 int nth = GB_nthreads (anz, chunk, nthreads_max) ;
127 if (cistart > 0 && cvstart > 0)
128 {
129 int64_t pA ;
130 #pragma omp parallel for num_threads(nth) schedule(static)
131 for (pA = 0 ; pA < anz ; pA++)
132 {
133 Wi [pC + pA] += cistart ;
134 Wj [pC + pA] += cvstart ;
135 }
136 }
137 else if (cistart > 0)
138 {
139 int64_t pA ;
140 #pragma omp parallel for num_threads(nth) schedule(static)
141 for (pA = 0 ; pA < anz ; pA++)
142 {
143 Wi [pC + pA] += cistart ;
144 }
145 }
146 else if (cvstart > 0)
147 {
148 int64_t pA ;
149 #pragma omp parallel for num_threads(nth) schedule(static)
150 for (pA = 0 ; pA < anz ; pA++)
151 {
152 Wj [pC + pA] += cvstart ;
153 }
154 }
155
156 //------------------------------------------------------------------
157 // advance the tuple counter
158 //------------------------------------------------------------------
159
160 pC += anz ;
161 }
162 }
163
164 //--------------------------------------------------------------------------
165 // build C from the triplets
166 //--------------------------------------------------------------------------
167
168 GB_OK (GB_builder
169 (
170 C, // create C using a static or dynamic header
171 ctype, // C->type
172 cvlen, // C->vlen
173 cvdim, // C->vdim
174 csc, // C->is_csc
175 (int64_t **) &Wi, // Wi becomes C->i on output, or freed on error
176 &Wi_size,
177 (int64_t **) &Wj, // Wj, free on output
178 &Wj_size,
179 (GB_void **) &Wx, // Wx, free on output
180 &Wx_size,
181 false, // tuples need to be sorted
182 true, // no duplicates
183 cnz, // size of Wi and Wj in # of tuples
184 true, // is_matrix: unused
185 NULL, NULL, NULL, // original I,J,S tuples, not used here
186 cnz, // # of tuples
187 NULL, // op for assembling duplicates (there are none)
188 ccode, // type of Wx
189 Context
190 )) ;
191
192 C->hyper_switch = hyper_switch ;
193 C->bitmap_switch = bitmap_switch ;
194 C->sparsity = sparsity_control ;
195 ASSERT (C->static_header == static_header) ;
196 ASSERT (GB_IS_HYPERSPARSE (C)) ;
197
198 // workspace has been freed by GB_builder, or transplanted into C
199 ASSERT (Wi == NULL) ;
200 ASSERT (Wj == NULL) ;
201 ASSERT (Wx == NULL) ;
202
203 return (GrB_SUCCESS) ;
204 }
205
206