1 //------------------------------------------------------------------------------
2 // GB_mex_reduce_terminal: [c,flag] = sum(A), reduce to scalar
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 // Compute c=max(A,x) where all entries in A are known a priori to be <= x.
11 // x becomes the terminal value of a user-defined max monoid.
12 
13 #include "GB_mex.h"
14 
15 #define USAGE "c = GB_mex_reduce_terminal (A, terminal)"
16 
17 #define FREE_ALL                        \
18 {                                       \
19     GrB_Matrix_free_(&A) ;              \
20     GrB_BinaryOp_free_(&Max) ;          \
21     GrB_Monoid_free_(&Max_Terminal) ;   \
22     GB_mx_put_global (true) ;           \
23 }
24 
25 void maxdouble (double *z, const double *x, const double *y) ;
26 
maxdouble(double * z,const double * x,const double * y)27 void maxdouble (double *z, const double *x, const double *y)
28 {
29     // this is not safe with NaNs
30     (*z) = ((*x) > (*y)) ? (*x) : (*y) ;
31 }
32 
mexFunction(int nargout,mxArray * pargout[],int nargin,const mxArray * pargin[])33 void mexFunction
34 (
35     int nargout,
36     mxArray *pargout [ ],
37     int nargin,
38     const mxArray *pargin [ ]
39 )
40 {
41 
42     bool malloc_debug = GB_mx_get_global (true) ;
43     GrB_Matrix A = NULL ;
44     GrB_BinaryOp Max = NULL;
45     GrB_Monoid Max_Terminal = NULL ;
46     GrB_Info info ;
47 
48     // check inputs
49     if (nargout > 1 || nargin < 1 || nargin > 2)
50     {
51         mexErrMsgTxt ("Usage: " USAGE) ;
52     }
53 
54     #define GET_DEEP_COPY ;
55     #define FREE_DEEP_COPY ;
56 
57     // get A (shallow copy)
58     A = GB_mx_mxArray_to_Matrix (pargin [0], "A input", false, true) ;
59     if (A == NULL)
60     {
61         FREE_ALL ;
62         mexErrMsgTxt ("A failed") ;
63     }
64 
65     if (A->type != GrB_FP64)
66     {
67         FREE_ALL ;
68         mexErrMsgTxt ("A must be double precision") ;
69     }
70 
71     // get the terminal value, if present.  Default is 1.
72     double GET_SCALAR (1, double, terminal, 1) ;
73 
74     // create the Max operator
75     info = GrB_BinaryOp_new (&Max, maxdouble, GrB_FP64, GrB_FP64, GrB_FP64);
76     if (info != GrB_SUCCESS)
77     {
78         mexErrMsgTxt ("Max failed") ;
79     }
80 
81     // create the Max monoid
82     info = GxB_Monoid_terminal_new_FP64_(&Max_Terminal, Max, (double) 0,
83         terminal) ;
84     if (info != GrB_SUCCESS)
85     {
86         mexErrMsgTxt ("Max_Terminal failed") ;
87     }
88 
89     // reduce to a scalar
90     double c ;
91     info = GrB_Matrix_reduce_FP64_(&c, NULL, Max_Terminal, A, NULL) ;
92     if (info != GrB_SUCCESS)
93     {
94         printf ("error: %d\n", info) ;
95         mexErrMsgTxt ("reduce failed") ;
96     }
97 
98     // printf ("result %g\n", c) ;
99 
100     // return C to MATLAB as a scalar
101     pargout [0] = mxCreateDoubleScalar (c) ;
102 
103     FREE_ALL ;
104 }
105 
106