1function test14
2%TEST14 test GrB_reduce
3
4% SuiteSparse:GraphBLAS, Timothy A. Davis, (c) 2017-2021, All Rights Reserved.
5% SPDX-License-Identifier: Apache-2.0
6
7fprintf ('\ntest14: reduce to column and scalar\n') ;
8
9[~, ~, add_ops, types, ~, ~] = GB_spec_opsall ;
10types = types.all ;
11
12rng ('default') ;
13
14m = 8 ;
15n = 4 ;
16dt = struct ('inp0', 'tran') ;
17
18for k1 = 1:length(types)
19    atype = types {k1} ;
20    fprintf ('.') ;
21    A = GB_spec_random (m, n, 0.3, 100, atype) ;
22    B = GB_spec_random (n, m, 0.3, 100, atype) ;
23    w = GB_spec_random (m, 1, 0.3, 100, atype) ;
24    cin = GB_mex_cast (0, atype) ;
25    mask = GB_random_mask (m, 1, 0.5, true, false) ;
26
27    if (isequal (atype, 'logical'))
28        ops = {'or', 'and', 'xor', 'eq', 'any'} ;
29    else
30        ops = {'min', 'max', 'plus', 'times', 'any'} ;
31    end
32
33    if (isequal (atype, 'double'))
34        hrange = [0 1] ;
35        crange = [0 1] ;
36    else
37        hrange = 0 ;
38        crange = 1 ;
39    end
40
41    is_float = contains (atype, 'single') || contains (atype, 'double') ;
42
43    for A_is_hyper = 0:1
44    for A_is_csc   = 0:1
45
46    A.is_csc = A_is_csc ; A.is_hyper = A_is_hyper ;
47    B.is_csc = A_is_csc ; B.is_hyper = A_is_hyper ;
48
49    for k2 = 1:length(add_ops)
50        op = add_ops {k2} ;
51
52        if (isequal (op, 'any'))
53            tol = [ ] ;
54        elseif (contains (atype, 'single'))
55            tol = 1e-5 ;
56        elseif (contains (atype, 'double'))
57            tol = 1e-12 ;
58        else
59            tol = 0 ;
60        end
61
62        try
63            GB_spec_operator (op, atype) ;
64            identity = GB_spec_identity (op, atype) ;
65        catch
66            continue
67        end
68
69        % no mask
70        w1 = GB_spec_reduce_to_vector (w, [], [], op, A, []) ;
71        w2 = GB_mex_reduce_to_vector  (w, [], [], op, A, []) ;
72        GB_spec_compare (w1, w2, identity, tol) ;
73
74        % no mask, with accum
75        w1 = GB_spec_reduce_to_vector (w, [], 'plus', op, A, []) ;
76        w2 = GB_mex_reduce_to_vector  (w, [], 'plus', op, A, []) ;
77        GB_spec_compare (w1, w2, identity, tol) ;
78
79        % with mask
80        w1 = GB_spec_reduce_to_vector (w, mask, [], op, A, []) ;
81        w2 = GB_mex_reduce_to_vector  (w, mask, [], op, A, []) ;
82        GB_spec_compare (w1, w2, identity, tol) ;
83
84        % with mask and accum
85        w1 = GB_spec_reduce_to_vector (w, mask, 'plus', op, A, []) ;
86        w2 = GB_mex_reduce_to_vector  (w, mask, 'plus', op, A, []) ;
87        GB_spec_compare (w1, w2, identity, tol) ;
88
89        % no mask, transpose
90        w1 = GB_spec_reduce_to_vector (w, [], [], op, B, dt) ;
91        w2 = GB_mex_reduce_to_vector  (w, [], [], op, B, dt) ;
92        GB_spec_compare (w1, w2, identity, tol) ;
93
94        % no mask, with accum, transpose
95        w1 = GB_spec_reduce_to_vector (w, [], 'plus', op, B, dt) ;
96        w2 = GB_mex_reduce_to_vector  (w, [], 'plus', op, B, dt) ;
97        GB_spec_compare (w1, w2, identity, tol) ;
98
99        % with mask, transpose
100        w1 = GB_spec_reduce_to_vector (w, mask, [], op, B, dt) ;
101        w2 = GB_mex_reduce_to_vector  (w, mask, [], op, B, dt) ;
102        GB_spec_compare (w1, w2, identity, tol) ;
103
104        % with mask and accum, transpose
105        w1 = GB_spec_reduce_to_vector (w, mask, 'plus', op, B, dt) ;
106        w2 = GB_mex_reduce_to_vector  (w, mask, 'plus', op, B, dt) ;
107        GB_spec_compare (w1, w2, identity, tol) ;
108
109        % GB_spec_reduce_to_scalar always operates column-wise, but GrB_reduce
110        % operates in whatever order it is given: by column if CSC or by row if
111        % CSR.  The result can vary slightly because of different round off
112        % errors.  A_flip causes GB_spec_reduce_to_scalar to operate in the
113        % same order as GrB_reduce.
114
115        A_flip = A ;
116        if (~A.is_csc && is_float)
117            A_flip.matrix = A_flip.matrix.' ;
118            A_flip.pattern = A_flip.pattern' ;
119            A_flip.is_csc = true ;
120        end
121
122        % Parallel reduction leads to different roundoff.  So even with A_flip,
123        % c1 and c2 can only be compared to within round-off error.
124
125        % to scalar
126        c2 = GB_mex_reduce_to_scalar  (cin, [ ], op, A) ;
127        if (isequal (op, 'any'))
128            X = GB_mex_cast (full (A.matrix (A.pattern)), A.class) ;
129            assert (any (X == c2)) ;
130        else
131            c1 = GB_spec_reduce_to_scalar (cin, [ ], op, A_flip) ;
132            if (is_float)
133                assert (abs (c1-c2) < tol *  (abs(c1) + 1))
134            else
135                assert (isequal (c1, c2)) ;
136            end
137        end
138
139        % to scalar, with accum
140        c2 = GB_mex_reduce_to_scalar (cin, 'plus', op, A) ;
141        if (~isequal (op, 'any'))
142            c1 = GB_spec_reduce_to_scalar (cin, 'plus', op, A_flip) ;
143            if (is_float)
144                assert (abs (c1-c2) < tol *  (abs(c1) + 1))
145            else
146                assert (isequal (c1, c2)) ;
147            end
148        end
149    end
150    end
151    end
152end
153
154fprintf ('\ntest14: all tests passed\n') ;
155
156