1 //------------------------------------------------------------------------------
2 // gbextract: extract entries into a GraphBLAS 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 // gbextract is an interface to GrB_Matrix_extract and GrB_Matrix_extract_[TYPE],
11 // computing the GraphBLAS expression:
12
13 // C<#M,replace> = accum (C, A (I,J)) or
14 // C<#M,replace> = accum (C, AT (I,J))
15
16 // Usage:
17
18 // C = gbextract (Cin, M, accum, A, I, J, desc)
19
20 // A is required. See GrB.m for more details.
21 // If accum or M is used, then Cin must appear.
22
23 #include "gb_matlab.h"
24 #include "GB_ij.h"
25
26 #define USAGE "usage: C = GrB.extract (Cin, M, accum, A, I, J, desc)"
27
mexFunction(int nargout,mxArray * pargout[],int nargin,const mxArray * pargin[])28 void mexFunction
29 (
30 int nargout,
31 mxArray *pargout [ ],
32 int nargin,
33 const mxArray *pargin [ ]
34 )
35 {
36
37 //--------------------------------------------------------------------------
38 // check inputs
39 //--------------------------------------------------------------------------
40
41 gb_usage (nargin >= 1 && nargin <= 7 && nargout <= 2, USAGE) ;
42
43 //--------------------------------------------------------------------------
44 // find the arguments
45 //--------------------------------------------------------------------------
46
47 mxArray *Matrix [4], *String [2], *Cell [2] ;
48 base_enum_t base ;
49 kind_enum_t kind ;
50 GxB_Format_Value fmt ;
51 int nmatrices, nstrings, ncells, sparsity ;
52 GrB_Descriptor desc ;
53 gb_get_mxargs (nargin, pargin, USAGE, Matrix, &nmatrices, String, &nstrings,
54 Cell, &ncells, &desc, &base, &kind, &fmt, &sparsity) ;
55
56 CHECK_ERROR (nmatrices < 1 || nmatrices > 3 || nstrings > 1, USAGE) ;
57
58 //--------------------------------------------------------------------------
59 // get the matrices
60 //--------------------------------------------------------------------------
61
62 GrB_Type atype, ctype = NULL ;
63 GrB_Matrix C = NULL, M = NULL, A ;
64
65 if (nmatrices == 1)
66 {
67 A = gb_get_shallow (Matrix [0]) ;
68 }
69 else if (nmatrices == 2)
70 {
71 C = gb_get_deep (Matrix [0]) ;
72 A = gb_get_shallow (Matrix [1]) ;
73 }
74 else // if (nmatrices == 3)
75 {
76 C = gb_get_deep (Matrix [0]) ;
77 M = gb_get_shallow (Matrix [1]) ;
78 A = gb_get_shallow (Matrix [2]) ;
79 }
80
81 OK (GxB_Matrix_type (&atype, A)) ;
82 if (C != NULL)
83 {
84 OK (GxB_Matrix_type (&ctype, C)) ;
85 }
86
87 //--------------------------------------------------------------------------
88 // get the operator
89 //--------------------------------------------------------------------------
90
91 GrB_BinaryOp accum = NULL ;
92
93 if (nstrings == 1)
94 {
95 // if accum appears, then Cin must also appear
96 CHECK_ERROR (C == NULL, USAGE) ;
97 accum = gb_mxstring_to_binop (String [0], ctype, ctype) ;
98 }
99
100 //--------------------------------------------------------------------------
101 // get the size of A
102 //--------------------------------------------------------------------------
103
104 GrB_Desc_Value in0 ;
105 OK (GxB_Desc_get (desc, GrB_INP0, &in0)) ;
106 GrB_Index anrows, ancols ;
107 bool A_transpose = (in0 == GrB_TRAN) ;
108 if (A_transpose)
109 {
110 // T = AT (I,J) is to be extracted where AT = A'
111 OK (GrB_Matrix_nrows (&ancols, A)) ;
112 OK (GrB_Matrix_ncols (&anrows, A)) ;
113 }
114 else
115 {
116 // T = A (I,J) is to be extracted
117 OK (GrB_Matrix_nrows (&anrows, A)) ;
118 OK (GrB_Matrix_ncols (&ancols, A)) ;
119 }
120
121 //--------------------------------------------------------------------------
122 // get I and J
123 //--------------------------------------------------------------------------
124
125 GrB_Index *I = (GrB_Index *) GrB_ALL ;
126 GrB_Index *J = (GrB_Index *) GrB_ALL ;
127 GrB_Index ni = anrows, nj = ancols ;
128 bool I_allocated = false, J_allocated = false ;
129
130 if (anrows == 1 && ncells == 1)
131 {
132 // only J is present
133 J = gb_mxcell_to_index (Cell [0], base, ancols, &J_allocated, &nj) ;
134 }
135 else if (ncells == 1)
136 {
137 // only I is present
138 I = gb_mxcell_to_index (Cell [0], base, anrows, &I_allocated, &ni) ;
139 }
140 else if (ncells == 2)
141 {
142 // both I and J are present
143 I = gb_mxcell_to_index (Cell [0], base, anrows, &I_allocated, &ni) ;
144 J = gb_mxcell_to_index (Cell [1], base, ancols, &J_allocated, &nj) ;
145 }
146
147 //--------------------------------------------------------------------------
148 // construct C if not present on input
149 //--------------------------------------------------------------------------
150
151 if (C == NULL)
152 {
153 // Cin is not present: determine its size, same type as A.
154 // T = A(I,J) or AT(I,J) will be extracted.
155 // accum must be null
156 int I_kind, J_kind ;
157 int64_t I_colon [3], J_colon [3] ;
158 GrB_Index cnrows, cncols ;
159 GB_ijlength (I, ni, anrows, (int64_t *) &cnrows, &I_kind, I_colon) ;
160 GB_ijlength (J, nj, ancols, (int64_t *) &cncols, &J_kind, J_colon) ;
161 ctype = atype ;
162
163 // create the matrix C and set its format and sparsity
164 fmt = gb_get_format (cnrows, cncols, A, NULL, fmt) ;
165 sparsity = gb_get_sparsity (A, NULL, sparsity) ;
166 C = gb_new (ctype, cnrows, cncols, fmt, sparsity) ;
167 }
168
169 //--------------------------------------------------------------------------
170 // compute C<M> += A(I,J) or AT(I,J)
171 //--------------------------------------------------------------------------
172
173 OK1 (C, GrB_Matrix_extract (C, M, accum, A, I, ni, J, nj, desc)) ;
174
175 //--------------------------------------------------------------------------
176 // free shallow copies
177 //--------------------------------------------------------------------------
178
179 OK (GrB_Matrix_free (&M)) ;
180 OK (GrB_Matrix_free (&A)) ;
181 OK (GrB_Descriptor_free (&desc)) ;
182 if (I_allocated) gb_mxfree (&I) ;
183 if (J_allocated) gb_mxfree (&J) ;
184
185 //--------------------------------------------------------------------------
186 // export the output matrix C back to MATLAB
187 //--------------------------------------------------------------------------
188
189 pargout [0] = gb_export (&C, kind) ;
190 pargout [1] = mxCreateDoubleScalar (kind) ;
191 GB_WRAPUP ;
192 }
193
194