1 //------------------------------------------------------------------------------ 2 // GB_positional_op_ijp: C = positional_op (A), depending j 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 // TODO: rename, and use #ifdef instead of offset = 0 or 1. 11 // TODO: use this kernel for GrB_extractTuples, to create J array. 12 13 // A can be jumbled. If A is jumbled, so is C. 14 // if A and C are bitmap, not all of Cx need to be written to, but it's faster 15 // just to write to all of it. C->b is copied from A->b in the caller. 16 17 { 18 19 //-------------------------------------------------------------------------- 20 // slice the entries for each task 21 //-------------------------------------------------------------------------- 22 23 GB_WERK_DECLARE (A_ek_slicing, int64_t) ; 24 int A_ntasks, A_nthreads ; 25 GB_SLICE_MATRIX (A, 32, chunk) ; 26 27 //-------------------------------------------------------------------------- 28 // Cx = positional_op (A) 29 //-------------------------------------------------------------------------- 30 31 int tid ; 32 #pragma omp parallel for num_threads(A_nthreads) schedule(dynamic,1) 33 for (tid = 0 ; tid < A_ntasks ; tid++) 34 { 35 36 // if kfirst > klast then task tid does no work at all 37 int64_t kfirst = kfirst_Aslice [tid] ; 38 int64_t klast = klast_Aslice [tid] ; 39 40 //---------------------------------------------------------------------- 41 // C(:,kfirst:klast) = op (A(:,kfirst:klast)) 42 //---------------------------------------------------------------------- 43 44 for (int64_t k = kfirst ; k <= klast ; k++) 45 { 46 47 //------------------------------------------------------------------ 48 // find the part of A(:,k) and Cx to be operated on by this task 49 //------------------------------------------------------------------ 50 51 int64_t j = GBH (Ah, k) ; 52 int64_t pA_start, pA_end ; 53 GB_get_pA (&pA_start, &pA_end, tid, k, 54 kfirst, klast, pstart_Aslice, Ap, avlen) ; 55 56 //------------------------------------------------------------------ 57 // C(:,j) = op (A(:,j)) 58 //------------------------------------------------------------------ 59 60 GB_PRAGMA_SIMD 61 for (int64_t p = pA_start ; p < pA_end ; p++) 62 { 63 // GB_POSITION is j or j+1 64 Cx_int [p] = GB_POSITION ; 65 } 66 } 67 } 68 69 //-------------------------------------------------------------------------- 70 // free workspace 71 //-------------------------------------------------------------------------- 72 73 GB_WERK_POP (A_ek_slicing, int64_t) ; 74 } 75 76 #undef GB_POSITION 77 78