1 //------------------------------------------------------------------------------
2 // GB_subassign_13: C(I,J)<!M> = scalar ; using S
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 // Method 13: C(I,J)<!M> = scalar ; using S
11
12 // M: present
13 // Mask_comp: true
14 // C_replace: false
15 // accum: NULL
16 // A: scalar
17 // S: constructed
18
19 // C: not bitmap, but can be full since no zombies are inserted in that case
20 // M: not bitmap
21
22 #include "GB_subassign_methods.h"
23
GB_subassign_13(GrB_Matrix C,const GrB_Index * I,const int64_t ni,const int64_t nI,const int Ikind,const int64_t Icolon[3],const GrB_Index * J,const int64_t nj,const int64_t nJ,const int Jkind,const int64_t Jcolon[3],const GrB_Matrix M,const bool Mask_struct,const void * scalar,const GrB_Type atype,GB_Context Context)24 GrB_Info GB_subassign_13
25 (
26 GrB_Matrix C,
27 // input:
28 const GrB_Index *I,
29 const int64_t ni,
30 const int64_t nI,
31 const int Ikind,
32 const int64_t Icolon [3],
33 const GrB_Index *J,
34 const int64_t nj,
35 const int64_t nJ,
36 const int Jkind,
37 const int64_t Jcolon [3],
38 const GrB_Matrix M,
39 const bool Mask_struct,
40 const void *scalar,
41 const GrB_Type atype,
42 GB_Context Context
43 )
44 {
45
46 //--------------------------------------------------------------------------
47 // check inputs
48 //--------------------------------------------------------------------------
49
50 ASSERT (!GB_IS_BITMAP (C)) ;
51 ASSERT (!GB_aliased (C, M)) ; // NO ALIAS of C==M
52
53 //--------------------------------------------------------------------------
54 // S = C(I,J)
55 //--------------------------------------------------------------------------
56
57 GB_EMPTY_TASKLIST ;
58 GB_OK (GB_subassign_symbolic (S, C, I, ni, J, nj, true, Context)) ;
59
60 //--------------------------------------------------------------------------
61 // get inputs
62 //--------------------------------------------------------------------------
63
64 GB_MATRIX_WAIT_IF_JUMBLED (M) ;
65
66 GB_GET_C ; // C must not be bitmap
67 const int64_t Cnvec = C->nvec ;
68 const int64_t *restrict Ch = C->h ;
69 const int64_t *restrict Cp = C->p ;
70 const bool C_is_hyper = (Ch != NULL) ;
71 GB_GET_MASK ;
72 GB_GET_SCALAR ;
73 GB_GET_S ;
74 GrB_BinaryOp accum = NULL ;
75
76 //--------------------------------------------------------------------------
77 // Method 13: C(I,J)<!M> = scalar ; using S
78 //--------------------------------------------------------------------------
79
80 // Time: Close to optimal; must visit all IxJ, so Omega(|I|*|J|) is
81 // required. The sparsity of !M cannot be exploited.
82
83 // Methods 13, 15, 17, and 19 are very similar.
84
85 //--------------------------------------------------------------------------
86 // Parallel: all IxJ (Methods 01, 03, 13, 15, 17, 19)
87 //--------------------------------------------------------------------------
88
89 GB_SUBASSIGN_IXJ_SLICE ;
90
91 //--------------------------------------------------------------------------
92 // phase 1: create zombies, update entries, and count pending tuples
93 //--------------------------------------------------------------------------
94
95 #pragma omp parallel for num_threads(nthreads) schedule(dynamic,1) \
96 reduction(+:nzombies)
97 for (taskid = 0 ; taskid < ntasks ; taskid++)
98 {
99
100 //----------------------------------------------------------------------
101 // get the task descriptor
102 //----------------------------------------------------------------------
103
104 GB_GET_IXJ_TASK_DESCRIPTOR_PHASE1 (iA_start, iA_end) ;
105
106 //----------------------------------------------------------------------
107 // compute all vectors in this task
108 //----------------------------------------------------------------------
109
110 for (int64_t j = kfirst ; j <= klast ; j++)
111 {
112
113 //------------------------------------------------------------------
114 // get jC, the corresponding vector of C
115 //------------------------------------------------------------------
116
117 int64_t jC = GB_ijlist (J, j, Jkind, Jcolon) ;
118
119 //------------------------------------------------------------------
120 // get S(iA_start:end,j) and M(iA_start:end,j)
121 //------------------------------------------------------------------
122
123 GB_GET_VECTOR_FOR_IXJ (S, iA_start) ;
124 GB_GET_VECTOR_FOR_IXJ (M, iA_start) ;
125
126 //------------------------------------------------------------------
127 // C(I(iA_start,iA_end-1),jC)<!M,repl> = scalar
128 //------------------------------------------------------------------
129
130 for (int64_t iA = iA_start ; iA < iA_end ; iA++)
131 {
132
133 //--------------------------------------------------------------
134 // Get the indices at the top of each list.
135 //--------------------------------------------------------------
136
137 int64_t iS = (pS < pS_end) ? GBI (Si, pS, Svlen) : INT64_MAX ;
138 int64_t iM = (pM < pM_end) ? GBI (Mi, pM, Mvlen) : INT64_MAX ;
139
140 //--------------------------------------------------------------
141 // find the smallest index of [iS iA iM] (always iA)
142 //--------------------------------------------------------------
143
144 int64_t i = iA ;
145
146 //--------------------------------------------------------------
147 // get M(i,j)
148 //--------------------------------------------------------------
149
150 bool mij ;
151 if (i == iM)
152 {
153 // mij = (bool) M [pM]
154 mij = GBB (Mb, pM) && GB_mcast (Mx, pM, msize) ;
155 GB_NEXT (M) ;
156 }
157 else
158 {
159 // mij not present, implicitly false
160 ASSERT (i < iM) ;
161 mij = false ;
162 }
163
164 // complement the mask entry mij since Mask_comp is true
165 mij = !mij ;
166
167 //--------------------------------------------------------------
168 // assign the entry
169 //--------------------------------------------------------------
170
171 if (i == iS)
172 {
173 ASSERT (i == iA) ;
174 {
175 // both S (i,j) and A (i,j) present
176 if (mij)
177 {
178 // ----[C A 1] or [X A 1]---------------------------
179 // [C A 1]: action: ( =A ): copy A, no accum
180 // [X A 1]: action: ( undelete ): zombie lives
181 GB_C_S_LOOKUP ;
182 GB_noaccum_C_A_1_scalar ;
183 }
184 GB_NEXT (S) ;
185 }
186 }
187 else
188 {
189 ASSERT (i == iA) ;
190 {
191 // S (i,j) is not present, A (i,j) is present
192 if (mij)
193 {
194 // ----[. A 1]--------------------------------------
195 // [. A 1]: action: ( insert )
196 task_pending++ ;
197 }
198 }
199 }
200 }
201 }
202
203 GB_PHASE1_TASK_WRAPUP ;
204 }
205
206 //--------------------------------------------------------------------------
207 // phase 2: insert pending tuples
208 //--------------------------------------------------------------------------
209
210 GB_PENDING_CUMSUM ;
211
212 #pragma omp parallel for num_threads(nthreads) schedule(dynamic,1) \
213 reduction(&&:pending_sorted)
214 for (taskid = 0 ; taskid < ntasks ; taskid++)
215 {
216
217 //----------------------------------------------------------------------
218 // get the task descriptor
219 //----------------------------------------------------------------------
220
221 GB_GET_IXJ_TASK_DESCRIPTOR_PHASE2 (iA_start, iA_end) ;
222
223 //----------------------------------------------------------------------
224 // compute all vectors in this task
225 //----------------------------------------------------------------------
226
227 for (int64_t j = kfirst ; j <= klast ; j++)
228 {
229
230 //------------------------------------------------------------------
231 // get jC, the corresponding vector of C
232 //------------------------------------------------------------------
233
234 int64_t jC = GB_ijlist (J, j, Jkind, Jcolon) ;
235
236 //------------------------------------------------------------------
237 // get S(iA_start:end,j) and M(iA_start:end,j)
238 //------------------------------------------------------------------
239
240 GB_GET_VECTOR_FOR_IXJ (S, iA_start) ;
241 GB_GET_VECTOR_FOR_IXJ (M, iA_start) ;
242
243 //------------------------------------------------------------------
244 // C(I(iA_start,iA_end-1),jC)<!M,repl> = scalar
245 //------------------------------------------------------------------
246
247 for (int64_t iA = iA_start ; iA < iA_end ; iA++)
248 {
249
250 //--------------------------------------------------------------
251 // Get the indices at the top of each list.
252 //--------------------------------------------------------------
253
254 int64_t iS = (pS < pS_end) ? GBI (Si, pS, Svlen) : INT64_MAX ;
255 int64_t iM = (pM < pM_end) ? GBI (Mi, pM, Mvlen) : INT64_MAX ;
256
257 //--------------------------------------------------------------
258 // find the smallest index of [iS iA iM] (always iA)
259 //--------------------------------------------------------------
260
261 int64_t i = iA ;
262
263 //--------------------------------------------------------------
264 // get M(i,j)
265 //--------------------------------------------------------------
266
267 bool mij ;
268 if (i == iM)
269 {
270 // mij = (bool) M [pM]
271 mij = GBB (Mb, pM) && GB_mcast (Mx, pM, msize) ;
272 GB_NEXT (M) ;
273 }
274 else
275 {
276 // mij not present, implicitly false
277 ASSERT (i < iM) ;
278 mij = false ;
279 }
280
281 // complement the mask entry mij since Mask_comp is true
282 mij = !mij ;
283
284 //--------------------------------------------------------------
285 // assign the entry
286 //--------------------------------------------------------------
287
288 if (i == iS)
289 {
290 ASSERT (i == iA) ;
291 {
292 GB_NEXT (S) ;
293 }
294 }
295 else
296 {
297 ASSERT (i == iA) ;
298 {
299 // S (i,j) is not present, A (i,j) is present
300 if (mij)
301 {
302 // ----[. A 1]--------------------------------------
303 // [. A 1]: action: ( insert )
304 int64_t iC = GB_ijlist (I, iA, Ikind, Icolon) ;
305 GB_PENDING_INSERT (scalar) ;
306 }
307 }
308 }
309 }
310 }
311
312 GB_PHASE2_TASK_WRAPUP ;
313 }
314
315 //--------------------------------------------------------------------------
316 // finalize the matrix and return result
317 //--------------------------------------------------------------------------
318
319 GB_SUBASSIGN_WRAPUP ;
320 }
321
322