1 //------------------------------------------------------------------------------
2 // gbbuild: build a GraphBLAS matrix or a MATLAB sparse matrix
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 // Usage:
11
12 // A = gbbuild (I, J, X)
13 // A = gbbuild (I, J, X, desc)
14 // A = gbbuild (I, J, X, m, desc)
15 // A = gbbuild (I, J, X, m, n, desc)
16 // A = gbbuild (I, J, X, m, n, dup, desc) ;
17 // A = gbbuild (I, J, X, m, n, dup, type, desc) ;
18
19 // m and n default to the largest index in I and J, respectively.
20
21 // dup is a string that defaults to 'plus.xtype' where xtype is the type of X.
22 // If dup is given by without a type, type of dup defaults to the type of X.
23
24 // type is a string that defines is the type of A, which defaults to the type
25 // of X.
26
27 // The descriptor is optional; if present, it must be the last input parameter.
28 // desc.kind is the only part used from the descriptor, and it defaults to
29 // desc.kind = 'GrB'.
30
31 #include "gb_matlab.h"
32
33 #define USAGE "usage: A = GrB.build (I, J, X, m, n, dup, type, desc)"
34
mexFunction(int nargout,mxArray * pargout[],int nargin,const mxArray * pargin[])35 void mexFunction
36 (
37 int nargout,
38 mxArray *pargout [ ],
39 int nargin,
40 const mxArray *pargin [ ]
41 )
42 {
43
44 //--------------------------------------------------------------------------
45 // check inputs
46 //--------------------------------------------------------------------------
47
48 gb_usage (nargin >= 3 && nargin <= 8 && nargout <= 2, USAGE) ;
49
50 //--------------------------------------------------------------------------
51 // get the descriptor
52 //--------------------------------------------------------------------------
53
54 base_enum_t base ;
55 kind_enum_t kind ;
56 GxB_Format_Value fmt ;
57 int sparsity ;
58 GrB_Descriptor desc = NULL ;
59 desc = gb_mxarray_to_descriptor (pargin [nargin-1], &kind, &fmt,
60 &sparsity, &base) ;
61
62 // if present, remove the descriptor from consideration
63 if (desc != NULL) nargin-- ;
64
65 OK (GrB_Descriptor_free (&desc)) ;
66
67 //--------------------------------------------------------------------------
68 // get I and J
69 //--------------------------------------------------------------------------
70
71 GrB_Index ni, nj ;
72 bool I_allocated, J_allocated ;
73 int64_t Imax = -1, Jmax = -1 ;
74
75 GrB_Index *I = (GrB_Index *) gb_mxarray_to_list (pargin [0], base,
76 &I_allocated, (int64_t *) &ni, &Imax) ;
77
78 GrB_Index *J = (GrB_Index *) gb_mxarray_to_list (pargin [1], base,
79 &J_allocated, (int64_t *) &nj, &Jmax) ;
80
81 //--------------------------------------------------------------------------
82 // get X
83 //--------------------------------------------------------------------------
84
85 GrB_Type xtype = gb_mxarray_type (pargin [2]) ;
86 GrB_Index nx = mxGetNumberOfElements (pargin [2]) ;
87
88 //--------------------------------------------------------------------------
89 // check the sizes of I, J, and X, and the type of X
90 //--------------------------------------------------------------------------
91
92 GrB_Index nvals = MAX (ni, nj) ;
93 nvals = MAX (nvals, nx) ;
94
95 if (!(ni == 1 || ni == nvals) ||
96 !(nj == 1 || nj == nvals) ||
97 !(nx == 1 || nx == nvals))
98 {
99 ERROR ("I, J, and X must have the same length") ;
100 }
101
102 CHECK_ERROR (!(mxIsNumeric (pargin [2]) || mxIsLogical (pargin [2])),
103 "X must be a numeric or logical array") ;
104 CHECK_ERROR (mxIsSparse (pargin [2]), "X cannot be sparse") ;
105
106 //--------------------------------------------------------------------------
107 // expand any scalars
108 //--------------------------------------------------------------------------
109
110 if (ni == 1 && ni < nvals)
111 {
112 GrB_Index *I2 = (GrB_Index *) mxMalloc (nvals * sizeof (GrB_Index)) ;
113 GB_matlab_helper8 ((GB_void *) I2, (GB_void *) I, nvals,
114 sizeof (GrB_Index)) ;
115 if (I_allocated) gb_mxfree (&I) ;
116 I_allocated = true ;
117 I = I2 ;
118 }
119
120 if (nj == 1 && nj < nvals)
121 {
122 GrB_Index *J2 = (GrB_Index *) mxMalloc (nvals * sizeof (GrB_Index)) ;
123 GB_matlab_helper8 ((GB_void *) J2, (GB_void *) J, nvals,
124 sizeof (GrB_Index)) ;
125 if (J_allocated) gb_mxfree (&J) ;
126 J_allocated = true ;
127 J = J2 ;
128 }
129
130 //--------------------------------------------------------------------------
131 // get m and n if present
132 //--------------------------------------------------------------------------
133
134 GrB_Index nrows = 0, ncols = 0 ;
135
136 if (nargin < 4)
137 {
138 // nrows = max entry in I + 1
139 if (Imax > -1)
140 {
141 // Imax already computed
142 nrows = Imax ;
143 }
144 else
145 {
146 // nrows = max entry in I+1
147 bool ok = GB_matlab_helper4 (I, ni, &nrows) ;
148 CHECK_ERROR (!ok, "out of memory") ;
149 }
150 }
151 else
152 {
153 // m is provided on input
154 CHECK_ERROR (!gb_mxarray_is_scalar (pargin [3]), "m must be a scalar") ;
155 nrows = (GrB_Index) mxGetScalar (pargin [3]) ;
156 }
157
158 if (nargin < 5)
159 {
160 if (Jmax > -1)
161 {
162 // Jmax already computed
163 ncols = Jmax ;
164 }
165 else
166 {
167 // ncols = max entry in J+1
168 bool ok = GB_matlab_helper4 (J, nj, &ncols) ;
169 CHECK_ERROR (!ok, "out of memory") ;
170 }
171 }
172 else
173 {
174 // n is provided on input
175 CHECK_ERROR (!gb_mxarray_is_scalar (pargin [4]), "n must be a scalar") ;
176 ncols = (GrB_Index) mxGetScalar (pargin [4]) ;
177 }
178
179 //--------------------------------------------------------------------------
180 // get the dup operator
181 //--------------------------------------------------------------------------
182
183 GrB_BinaryOp dup = NULL ;
184 if (nargin > 5)
185 {
186 dup = gb_mxstring_to_binop (pargin [5], xtype, xtype) ;
187 }
188
189 // if dup is NULL, defaults to plus.xtype, below.
190
191 //--------------------------------------------------------------------------
192 // get the output matrix type
193 //--------------------------------------------------------------------------
194
195 GrB_Type type = NULL ;
196 if (nargin > 6)
197 {
198 type = gb_mxstring_to_type (pargin [6]) ;
199 CHECK_ERROR (type == NULL, "unknown type") ;
200 }
201 else
202 {
203 type = xtype ;
204 }
205
206 //--------------------------------------------------------------------------
207 // build the matrix
208 //--------------------------------------------------------------------------
209
210 fmt = gb_get_format (nrows, ncols, NULL, NULL, fmt) ;
211 sparsity = gb_get_sparsity (NULL, NULL, sparsity) ;
212 GrB_Matrix A = gb_new (type, nrows, ncols, fmt, sparsity) ;
213
214 // expandx is true if X must be expanded from a scalar to a vector
215 void *X2 = NULL ;
216 bool expandx = (nx == 1 && nx < nvals) ;
217
218 if (xtype == GrB_BOOL)
219 {
220 bool empty = 0 ;
221 bool *X = (nvals == 0) ? &empty : mxGetData (pargin [2]) ;
222 if (dup == NULL) dup = GrB_LOR ;
223 if (expandx)
224 {
225 X2 = mxMalloc (nvals * sizeof (bool)) ;
226 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
227 sizeof (bool)) ;
228 X = (bool *) X2 ;
229 }
230 OK1 (A, GrB_Matrix_build_BOOL (A, I, J, X, nvals, dup)) ;
231 }
232 else if (xtype == GrB_INT8)
233 {
234 int8_t empty = 0 ;
235 int8_t *X = (nvals == 0) ? &empty : mxGetInt8s (pargin [2]) ;
236 if (dup == NULL) dup = GrB_PLUS_INT8 ;
237 if (expandx)
238 {
239 X2 = mxMalloc (nvals * sizeof (int8_t)) ;
240 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
241 sizeof (int8_t)) ;
242 X = (int8_t *) X2 ;
243 }
244 OK1 (A, GrB_Matrix_build_INT8 (A, I, J, X, nvals, dup)) ;
245 }
246 else if (xtype == GrB_INT16)
247 {
248 int16_t empty = 0 ;
249 int16_t *X = (nvals == 0) ? &empty : mxGetInt16s (pargin [2]) ;
250 if (dup == NULL) dup = GrB_PLUS_INT16 ;
251 if (expandx)
252 {
253 X2 = mxMalloc (nvals * sizeof (int16_t)) ;
254 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
255 sizeof (int16_t)) ;
256 X = (int16_t *) X2 ;
257 }
258 OK1 (A, GrB_Matrix_build_INT16 (A, I, J, X, nvals, dup)) ;
259 }
260 else if (xtype == GrB_INT32)
261 {
262 int32_t empty = 0 ;
263 int32_t *X = (nvals == 0) ? &empty : mxGetInt32s (pargin [2]) ;
264 if (dup == NULL) dup = GrB_PLUS_INT32 ;
265 if (expandx)
266 {
267 X2 = mxMalloc (nvals * sizeof (int32_t)) ;
268 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
269 sizeof (int32_t)) ;
270 X = (int32_t *) X2 ;
271 }
272 OK1 (A, GrB_Matrix_build_INT32 (A, I, J, X, nvals, dup)) ;
273 }
274 else if (xtype == GrB_INT64)
275 {
276 int64_t empty = 0 ;
277 int64_t *X = (nvals == 0) ? &empty : mxGetInt64s (pargin [2]) ;
278 if (dup == NULL) dup = GrB_PLUS_INT64 ;
279 if (expandx)
280 {
281 X2 = mxMalloc (nvals * sizeof (int64_t)) ;
282 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
283 sizeof (int64_t)) ;
284 X = (int64_t *) X2 ;
285 }
286 OK1 (A, GrB_Matrix_build_INT64 (A, I, J, X, nvals, dup)) ;
287 }
288 else if (xtype == GrB_UINT8)
289 {
290 uint8_t empty = 0 ;
291 uint8_t *X = (nvals == 0) ? &empty : mxGetUint8s (pargin [2]) ;
292 if (dup == NULL) dup = GrB_PLUS_UINT8 ;
293 if (expandx)
294 {
295 X2 = mxMalloc (nvals * sizeof (uint8_t)) ;
296 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
297 sizeof (uint8_t)) ;
298 X = (uint8_t *) X2 ;
299 }
300 OK1 (A, GrB_Matrix_build_UINT8 (A, I, J, X, nvals, dup)) ;
301 }
302 else if (xtype == GrB_UINT16)
303 {
304 uint16_t empty = 0 ;
305 uint16_t *X = (nvals == 0) ? &empty : mxGetUint16s (pargin [2]) ;
306 if (dup == NULL) dup = GrB_PLUS_UINT16 ;
307 if (expandx)
308 {
309 X2 = mxMalloc (nvals * sizeof (uint16_t)) ;
310 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
311 sizeof (uint16_t)) ;
312 X = (uint16_t *) X2 ;
313 }
314 OK1 (A, GrB_Matrix_build_UINT16 (A, I, J, X, nvals, dup)) ;
315 }
316 else if (xtype == GrB_UINT32)
317 {
318 uint32_t empty = 0 ;
319 uint32_t *X = (nvals == 0) ? &empty : mxGetUint32s (pargin [2]) ;
320 if (dup == NULL) dup = GrB_PLUS_UINT32 ;
321 if (expandx)
322 {
323 X2 = mxMalloc (nvals * sizeof (uint32_t)) ;
324 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
325 sizeof (uint32_t)) ;
326 X = (uint32_t *) X2 ;
327 }
328 OK1 (A, GrB_Matrix_build_UINT32 (A, I, J, X, nvals, dup)) ;
329 }
330 else if (xtype == GrB_UINT64)
331 {
332 uint64_t empty = 0 ;
333 uint64_t *X = (nvals == 0) ? &empty : mxGetUint64s (pargin [2]) ;
334 if (dup == NULL) dup = GrB_PLUS_UINT64 ;
335 if (expandx)
336 {
337 X2 = mxMalloc (nvals * sizeof (uint64_t)) ;
338 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
339 sizeof (uint64_t)) ;
340 X = (uint64_t *) X2 ;
341 }
342 OK1 (A, GrB_Matrix_build_UINT64 (A, I, J, X, nvals, dup)) ;
343 }
344 else if (xtype == GrB_FP32)
345 {
346 float empty = 0 ;
347 float *X = (nvals == 0) ? &empty : mxGetSingles (pargin [2]) ;
348 if (dup == NULL) dup = GrB_PLUS_FP32 ;
349 if (expandx)
350 {
351 X2 = mxMalloc (nvals * sizeof (float)) ;
352 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
353 sizeof (float)) ;
354 X = (float *) X2 ;
355 }
356 OK1 (A, GrB_Matrix_build_FP32 (A, I, J, X, nvals, dup)) ;
357 }
358 else if (xtype == GrB_FP64)
359 {
360 double empty = 0 ;
361 double *X = (nvals == 0) ? &empty : mxGetDoubles (pargin [2]) ;
362 if (dup == NULL) dup = GrB_PLUS_FP64 ;
363 if (expandx)
364 {
365 X2 = mxMalloc (nvals * sizeof (double)) ;
366 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
367 sizeof (double)) ;
368 X = (double *) X2 ;
369 }
370 OK1 (A, GrB_Matrix_build_FP64 (A, I, J, X, nvals, dup)) ;
371 }
372 else if (xtype == GxB_FC32)
373 {
374 GxB_FC32_t empty = GxB_CMPLXF (0,0) ;
375 GxB_FC32_t *X = &empty ;
376 if (nvals > 0) X = (GxB_FC32_t *) mxGetComplexSingles (pargin [2]) ;
377 if (dup == NULL) dup = GxB_PLUS_FC32 ;
378 if (expandx)
379 {
380 X2 = mxMalloc (nvals * sizeof (GxB_FC32_t)) ;
381 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
382 sizeof (GxB_FC32_t)) ;
383 X = (GxB_FC32_t *) X2 ;
384 }
385 OK1 (A, GxB_Matrix_build_FC32 (A, I, J, X, nvals, dup)) ;
386 }
387 else if (xtype == GxB_FC64)
388 {
389 GxB_FC64_t empty = GxB_CMPLX (0,0) ;
390 GxB_FC64_t *X = &empty ;
391 if (nvals > 0) X = (GxB_FC64_t *) mxGetComplexDoubles (pargin [2]) ;
392 if (dup == NULL) dup = GxB_PLUS_FC64 ;
393 if (expandx)
394 {
395 X2 = mxMalloc (nvals * sizeof (GxB_FC64_t)) ;
396 GB_matlab_helper8 ((GB_void *) X2, (GB_void *) X, nvals,
397 sizeof (GxB_FC64_t)) ;
398 X = (GxB_FC64_t *) X2 ;
399 }
400 OK1 (A, GxB_Matrix_build_FC64 (A, I, J, X, nvals, dup)) ;
401 }
402 else
403 {
404 ERROR ("unsupported type") ;
405 }
406
407 //--------------------------------------------------------------------------
408 // free workspace
409 //--------------------------------------------------------------------------
410
411 if (X2 != NULL ) gb_mxfree (&X2) ;
412 if (I_allocated) gb_mxfree (&I) ;
413 if (J_allocated) gb_mxfree (&J) ;
414
415 //--------------------------------------------------------------------------
416 // export the output matrix A back to MATLAB
417 //--------------------------------------------------------------------------
418
419 pargout [0] = gb_export (&A, kind) ;
420 pargout [1] = mxCreateDoubleScalar (kind) ;
421 GB_WRAPUP ;
422 }
423
424