1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
|
//------------------------------------------------------------------------------
// GB_bitmap_assign_IxJ_template: iterate over all of C(I,J)
//------------------------------------------------------------------------------
// SuiteSparse:GraphBLAS, Timothy A. Davis, (c) 2017-2022, All Rights Reserved.
// SPDX-License-Identifier: Apache-2.0
//------------------------------------------------------------------------------
// Iterate over all positions in the IxJ Cartesian product. This is all
// entries C(i,j) where i is in the list I and j is in the list J. This
// traversal occurs whether or not C(i,j) is an entry present in C.
// The C matrix is accessed at C(I,J). The A matrix is size |I|-by-|J|.
// For bitmap assignent, C(I,J)=A is being computed. For bitmap extraction,
// C=A(I,J) so the roles of A and C are swapped (see GB_bitmap_subref.c).
{
//--------------------------------------------------------------------------
// create the tasks to iterate over IxJ
//--------------------------------------------------------------------------
int ntasks = 0, nthreads ;
GB_task_struct *TaskList = NULL ; size_t TaskList_size = 0 ;
GB_OK (GB_subassign_IxJ_slice (&TaskList, &TaskList_size, &ntasks,
&nthreads, /* I, */ nI, /* Ikind, Icolon, J, */ nJ,
/* Jkind, Jcolon, */ Context)) ;
//--------------------------------------------------------------------------
// iterate over all IxJ
//--------------------------------------------------------------------------
int taskid ;
#pragma omp parallel for num_threads(nthreads) schedule(dynamic,1) \
reduction(+:cnvals)
for (taskid = 0 ; taskid < ntasks ; taskid++)
{
//----------------------------------------------------------------------
// get the task descriptor
//----------------------------------------------------------------------
int64_t kfirst = TaskList [taskid].kfirst ;
int64_t klast = TaskList [taskid].klast ;
int64_t task_cnvals = 0 ;
bool fine_task = (klast == -1) ;
int64_t iA_start = 0, iA_end = nI ;
if (fine_task)
{
// a fine task operates on a slice of a single vector
klast = kfirst ;
iA_start = TaskList [taskid].pA ;
iA_end = TaskList [taskid].pA_end ;
}
//----------------------------------------------------------------------
// compute all vectors in this task
//----------------------------------------------------------------------
for (int64_t jA = kfirst ; jA <= klast ; jA++)
{
//------------------------------------------------------------------
// get jC, the corresponding vector of C
//------------------------------------------------------------------
int64_t jC = GB_ijlist (J, jA, Jkind, Jcolon) ;
int64_t pC0 = jC * vlen ; // first entry in C(:,jC)
int64_t pA0 = jA * nI ; // first entry in A(:,jA)
//------------------------------------------------------------------
// operate on C (I(iA_start,iA_end-1),jC)
//------------------------------------------------------------------
for (int64_t iA = iA_start ; iA < iA_end ; iA++)
{
int64_t iC = GB_ijlist (I, iA, Ikind, Icolon) ;
int64_t pC = iC + pC0 ;
int64_t pA = iA + pA0 ;
// operate on C(iC,jC) at pC (if C is bitmap or full)
// and A(iA,jA) or M(iA,jA) at pA, if A and/or M are
// bitmap or full. M(iA,jA) is accessed only for the
// subassign method when M is bitmap or full.
GB_IXJ_WORK (pC, pA) ;
}
}
cnvals += task_cnvals ;
}
//--------------------------------------------------------------------------
// free workpace
//--------------------------------------------------------------------------
GB_FREE_WORK (&TaskList, TaskList_size) ;
}
|
