File: GB_AxB_dot4_cij.c

package info (click to toggle)
suitesparse-graphblas 7.4.0%2Bdfsg-1
  • links: PTS, VCS
  • area: main
  • in suites: bookworm
  • size: 67,112 kB
  • sloc: ansic: 1,072,243; cpp: 8,081; sh: 512; makefile: 503; asm: 369; python: 125; awk: 10
file content (140 lines) | stat: -rw-r--r-- 4,487 bytes parent folder | download | duplicates (3)
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
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
//------------------------------------------------------------------------------
// GB_AxB_dot4_cij.c: C(i,j) = A(:,i)'*B(:,j) for dot4 method
//------------------------------------------------------------------------------

// SuiteSparse:GraphBLAS, Timothy A. Davis, (c) 2017-2022, All Rights Reserved.
// SPDX-License-Identifier: Apache-2.0

//------------------------------------------------------------------------------

// A is sparse or hypersparse, B is full or bitmap, and C is full

{

    //--------------------------------------------------------------------------
    // get C(i,j)
    //--------------------------------------------------------------------------

    const int64_t pC = i + pC_start ;   // C(i,j) is at Cx [pC]
    GB_CTYPE GB_GET4C (cij, pC) ;       // cij = Cx [pC]

    //--------------------------------------------------------------------------
    // C(i,j) += A (:,i)*B(:,j): a single dot product
    //--------------------------------------------------------------------------

    #if ( GB_B_IS_FULL )
    {

        //----------------------------------------------------------------------
        // A is sparse/hyper and B is full
        //----------------------------------------------------------------------

        #if GB_IS_PAIR_MULTIPLIER
        { 
            #if GB_IS_EQ_MONOID
            // EQ_PAIR semiring
            cij = (cij == 1) ;
            #elif (GB_CTYPE_BITS > 0)
            // PLUS, XOR monoids: A(:,i)'*B(:,j) is nnz(A(:,i)),
            // for bool, 8-bit, 16-bit, or 32-bit integer
            uint64_t t = ((uint64_t) cij) + ainz ;
            cij = (GB_CTYPE) (t & GB_CTYPE_BITS) ;
            #elif GB_IS_PLUS_FC32_MONOID
            // PLUS monoid for float complex
            cij = GxB_CMPLXF (crealf (cij) + (float) ainz, 0) ;
            #elif GB_IS_PLUS_FC64_MONOID
            // PLUS monoid for double complex
            cij = GxB_CMPLX (creal (cij) + (double) ainz, 0) ;
            #else
            // PLUS monoid for float, double, or 64-bit integers 
            cij += (GB_CTYPE) ainz ;
            #endif
        }
        #elif GB_IS_MIN_FIRSTJ_SEMIRING
        {
            // MIN_FIRSTJ semiring: take the 1st entry in A(:,i)
            if (ainz > 0)
            { 
                int64_t k = Ai [pA] + GB_OFFSET ;
                cij = GB_IMIN (cij, k) ;
            }
        }
        #elif GB_IS_MAX_FIRSTJ_SEMIRING
        {
            // MAX_FIRSTJ semiring: take last entry in A(:,i)
            if (ainz > 0)
            { 
                int64_t k = Ai [pA_end-1] + GB_OFFSET ;
                cij = GB_IMAX (cij, k) ;
            }
        }
        #else
        {
            GB_PRAGMA_SIMD_DOT (cij)
            for (int64_t p = pA ; p < pA_end ; p++)
            { 
                int64_t k = Ai [p] ;
                GB_DOT (k, p, pB+k) ;   // cij += A(k,i)*B(k,j)
            }
        }
        #endif

    }
    #else
    {

        //----------------------------------------------------------------------
        // A is sparse/hyper and B is bitmap
        //----------------------------------------------------------------------

        #if GB_IS_MIN_FIRSTJ_SEMIRING
        {
            // MIN_FIRSTJ semiring: take the first entry
            for (int64_t p = pA ; p < pA_end ; p++)
            {
                int64_t k = Ai [p] ;
                if (Bb [pB+k])
                { 
                    cij = GB_IMIN (cij, k + GB_OFFSET) ;
                    break ;
                }
            }
        }
        #elif GB_IS_MAX_FIRSTJ_SEMIRING
        {
            // MAX_FIRSTJ semiring: take the last entry
            for (int64_t p = pA_end-1 ; p >= pA ; p--)
            {
                int64_t k = Ai [p] ;
                if (Bb [pB+k])
                { 
                    cij = GB_IMAX (cij, k + GB_OFFSET) ;
                    break ;
                }
            }
        }
        #else
        {
            GB_PRAGMA_SIMD_DOT (cij)
            for (int64_t p = pA ; p < pA_end ; p++)
            {
                int64_t k = Ai [p] ;
                if (Bb [pB+k])
                { 
                    GB_DOT (k, p, pB+k) ; // cij+=A(k,i)*B(k,j)
                }
            }

        }
        #endif

    }
    #endif

    //--------------------------------------------------------------------------
    // save C(i,j)
    //--------------------------------------------------------------------------

    Cx [pC] = cij ;
}