1 //------------------------------------------------------------------------------
2 // gbmxm: sparse matrix-matrix multiplication
3 //------------------------------------------------------------------------------
4 
5 // SuiteSparse:GraphBLAS, Timothy A. Davis, (c) 2017-2021, All Rights Reserved.
6 // SPDX-License-Identifier: GPL-3.0-or-later
7 
8 //------------------------------------------------------------------------------
9 
10 // gbmxm is an interface to GrB_mxm.
11 
12 // Usage:
13 
14 // C = gbmxm (semiring, A, B)
15 // C = gbmxm (semiring, A, B, desc)
16 // C = gbmxm (Cin, accum, semiring, A, B, desc)
17 // C = gbmxm (Cin, M, semiring, A, B, desc)
18 // C = gbmxm (Cin, M, accum, semiring, A, B, desc)
19 
20 // If Cin is not present then it is implicitly a matrix with no entries, of the
21 // right size (which depends on A, B, and the descriptor).
22 
23 #include "gb_matlab.h"
24 
25 #define USAGE "usage: C = GrB.mxm (Cin, M, accum, semiring, A, B, desc)"
26 
mexFunction(int nargout,mxArray * pargout[],int nargin,const mxArray * pargin[])27 void mexFunction
28 (
29     int nargout,
30     mxArray *pargout [ ],
31     int nargin,
32     const mxArray *pargin [ ]
33 )
34 {
35 
36     //--------------------------------------------------------------------------
37     // check inputs
38     //--------------------------------------------------------------------------
39 
40     gb_usage (nargin >= 3 && nargin <= 7 && nargout <= 2, USAGE) ;
41 
42     //--------------------------------------------------------------------------
43     // find the arguments
44     //--------------------------------------------------------------------------
45 
46     mxArray *Matrix [4], *String [2], *Cell [2] ;
47     base_enum_t base ;
48     kind_enum_t kind ;
49     GxB_Format_Value fmt ;
50     int nmatrices, nstrings, ncells, sparsity ;
51     GrB_Descriptor desc ;
52     gb_get_mxargs (nargin, pargin, USAGE, Matrix, &nmatrices, String, &nstrings,
53         Cell, &ncells, &desc, &base, &kind, &fmt, &sparsity) ;
54 
55     CHECK_ERROR (nmatrices < 2 || nstrings < 1 || ncells > 0, USAGE) ;
56 
57     // ensure the descriptor is present, and set GxB_SORT to true
58     if (desc == NULL)
59     {
60         OK (GrB_Descriptor_new (&desc)) ;
61     }
62     OK (GxB_Desc_set (desc, GxB_SORT, true)) ;
63 
64     //--------------------------------------------------------------------------
65     // get the matrices
66     //--------------------------------------------------------------------------
67 
68     GrB_Type atype, btype, ctype = NULL ;
69     GrB_Matrix C = NULL, M = NULL, A, B ;
70 
71     if (nmatrices == 2)
72     {
73         A = gb_get_shallow (Matrix [0]) ;
74         B = gb_get_shallow (Matrix [1]) ;
75     }
76     else if (nmatrices == 3)
77     {
78         C = gb_get_deep    (Matrix [0]) ;
79         A = gb_get_shallow (Matrix [1]) ;
80         B = gb_get_shallow (Matrix [2]) ;
81     }
82     else // if (nmatrices == 4)
83     {
84         C = gb_get_deep    (Matrix [0]) ;
85         M = gb_get_shallow (Matrix [1]) ;
86         A = gb_get_shallow (Matrix [2]) ;
87         B = gb_get_shallow (Matrix [3]) ;
88     }
89 
90     OK (GxB_Matrix_type (&atype, A)) ;
91     OK (GxB_Matrix_type (&btype, B)) ;
92     if (C != NULL)
93     {
94         OK (GxB_Matrix_type (&ctype, C)) ;
95     }
96 
97     //--------------------------------------------------------------------------
98     // get the operators
99     //--------------------------------------------------------------------------
100 
101     GrB_BinaryOp accum = NULL ;
102     GrB_Semiring semiring ;
103 
104     if (nstrings == 1)
105     {
106         semiring = gb_mxstring_to_semiring (String [0], atype, btype) ;
107     }
108     else
109     {
110         // if accum appears, then Cin must also appear
111         CHECK_ERROR (C == NULL, USAGE) ;
112         accum    = gb_mxstring_to_binop    (String [0], ctype, ctype) ;
113         semiring = gb_mxstring_to_semiring (String [1], atype, btype) ;
114     }
115 
116     //--------------------------------------------------------------------------
117     // construct C if not present on input
118     //--------------------------------------------------------------------------
119 
120     // If C is NULL, then it is not present on input.
121     // Construct C of the right size and type.
122 
123     if (C == NULL)
124     {
125         // get the descriptor contents to determine if A and B are transposed
126         GrB_Desc_Value in0, in1 ;
127         OK (GxB_Desc_get (desc, GrB_INP0, &in0)) ;
128         OK (GxB_Desc_get (desc, GrB_INP1, &in1)) ;
129         bool A_transpose = (in0 == GrB_TRAN) ;
130         bool B_transpose = (in1 == GrB_TRAN) ;
131 
132         // get the size of A and B
133         GrB_Index anrows, ancols, bnrows, bncols ;
134         OK (GrB_Matrix_nrows (&anrows, A)) ;
135         OK (GrB_Matrix_ncols (&ancols, A)) ;
136         OK (GrB_Matrix_nrows (&bnrows, B)) ;
137         OK (GrB_Matrix_ncols (&bncols, B)) ;
138 
139         // determine the size of C
140         GrB_Index cnrows = (A_transpose) ? ancols : anrows ;
141         GrB_Index cncols = (B_transpose) ? bnrows : bncols ;
142 
143         // use the semiring's additive monoid as the type of C
144         GrB_Monoid add_monoid ;
145         GrB_BinaryOp add ;
146         OK (GxB_Semiring_add (&add_monoid, semiring)) ;
147         OK (GxB_Monoid_operator (&add, add_monoid)) ;
148         OK (GxB_BinaryOp_ztype (&ctype, add)) ;
149 
150         // create the matrix C and set its format and sparsity
151         fmt = gb_get_format (cnrows, cncols, A, B, fmt) ;
152         sparsity = gb_get_sparsity (A, B, sparsity) ;
153         C = gb_new (ctype, cnrows, cncols, fmt, sparsity) ;
154     }
155 
156     //--------------------------------------------------------------------------
157     // compute C<M> += A*B
158     //--------------------------------------------------------------------------
159 
160     OK1 (C, GrB_mxm (C, M, accum, semiring, A, B, desc)) ;
161 
162     //--------------------------------------------------------------------------
163     // free shallow copies
164     //--------------------------------------------------------------------------
165 
166     OK (GrB_Matrix_free (&M)) ;
167     OK (GrB_Matrix_free (&A)) ;
168     OK (GrB_Matrix_free (&B)) ;
169     OK (GrB_Descriptor_free (&desc)) ;
170 
171     //--------------------------------------------------------------------------
172     // export the output matrix C back to MATLAB
173     //--------------------------------------------------------------------------
174 
175     pargout [0] = gb_export (&C, kind) ;
176     pargout [1] = mxCreateDoubleScalar (kind) ;
177     GB_WRAPUP ;
178 }
179 
180