1 /****************************************************************/
2 /* Parallel Combinatorial BLAS Library (for Graph Computations) */
3 /* version 1.6 -------------------------------------------------*/
4 /* date: 6/15/2017 ---------------------------------------------*/
5 /* authors: Ariful Azad, Aydin Buluc --------------------------*/
6 /****************************************************************/
7 /*
8 Copyright (c) 2010-2017, The Regents of the University of California
9
10 Permission is hereby granted, free of charge, to any person obtaining a copy
11 of this software and associated documentation files (the "Software"), to deal
12 in the Software without restriction, including without limitation the rights
13 to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
14 copies of the Software, and to permit persons to whom the Software is
15 furnished to do so, subject to the following conditions:
16
17 The above copyright notice and this permission notice shall be included in
18 all copies or substantial portions of the Software.
19
20 THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
21 IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
22 FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
23 AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
24 LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
25 OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
26 THE SOFTWARE.
27 */
28
29
30 #include "SpTuples.h"
31 #include "SpParHelper.h"
32 #include <iomanip>
33
34 namespace combblas {
35
36 template <class IT,class NT>
SpTuples(int64_t size,IT nRow,IT nCol)37 SpTuples<IT,NT>::SpTuples(int64_t size, IT nRow, IT nCol)
38 :m(nRow), n(nCol), nnz(size)
39 {
40 if(nnz > 0)
41 {
42 tuples = new std::tuple<IT, IT, NT>[nnz];
43 }
44 else
45 {
46 tuples = NULL;
47 }
48 isOperatorNew = false;
49 }
50
51 template <class IT,class NT>
SpTuples(int64_t size,IT nRow,IT nCol,std::tuple<IT,IT,NT> * mytuples,bool sorted,bool isOpNew)52 SpTuples<IT,NT>::SpTuples (int64_t size, IT nRow, IT nCol, std::tuple<IT, IT, NT> * mytuples, bool sorted, bool isOpNew)
53 :tuples(mytuples), m(nRow), n(nCol), nnz(size), isOperatorNew(isOpNew)
54 {
55 if(!sorted)
56 {
57 SortColBased();
58 }
59
60 }
61
62 /**
63 * Generate a SpTuples object from an edge list
64 * @param[in,out] edges: edge list that might contain duplicate edges. freed upon return
65 * Semantics differ depending on the object created:
66 * NT=bool: duplicates are ignored
67 * NT='countable' (such as short,int): duplicated as summed to keep count
68 **/
69 template <class IT, class NT>
SpTuples(int64_t maxnnz,IT nRow,IT nCol,std::vector<IT> & edges,bool removeloops)70 SpTuples<IT,NT>::SpTuples (int64_t maxnnz, IT nRow, IT nCol, std::vector<IT> & edges, bool removeloops):m(nRow), n(nCol)
71 {
72 if(maxnnz > 0)
73 {
74 tuples = new std::tuple<IT, IT, NT>[maxnnz];
75 }
76 for(int64_t i=0; i<maxnnz; ++i)
77 {
78 rowindex(i) = edges[2*i+0];
79 colindex(i) = edges[2*i+1];
80 numvalue(i) = (NT) 1;
81 }
82 std::vector<IT>().swap(edges); // free memory for edges
83
84 nnz = maxnnz; // for now (to sort)
85 SortColBased();
86
87 int64_t cnz = 0;
88 int64_t dup = 0; int64_t self = 0;
89 nnz = 0;
90 while(cnz < maxnnz)
91 {
92 int64_t j=cnz+1;
93 while(j < maxnnz && rowindex(cnz) == rowindex(j) && colindex(cnz) == colindex(j))
94 {
95 numvalue(cnz) += numvalue(j);
96 numvalue(j++) = 0; // mark for deletion
97 ++dup;
98 }
99 if(removeloops && rowindex(cnz) == colindex(cnz))
100 {
101 numvalue(cnz) = 0;
102 --nnz;
103 ++self;
104 }
105 ++nnz;
106 cnz = j;
107 }
108
109 std::tuple<IT, IT, NT> * ntuples = new std::tuple<IT,IT,NT>[nnz];
110 int64_t j = 0;
111 for(int64_t i=0; i<maxnnz; ++i)
112 {
113 if(numvalue(i) != 0)
114 {
115 ntuples[j++] = tuples[i];
116 }
117 }
118 assert(j == nnz);
119
120 delete [] tuples;
121 tuples = ntuples;
122 isOperatorNew = false;
123 }
124
125
126 /**
127 * Generate a SpTuples object from StackEntry array, then delete that array
128 * @param[in] multstack {value-key pairs where keys are pair<col_ind, row_ind> sorted lexicographically}
129 * \remark Since input is column sorted, the tuples are automatically generated in that way too
130 **/
131 template <class IT, class NT>
SpTuples(int64_t size,IT nRow,IT nCol,StackEntry<NT,std::pair<IT,IT>> * & multstack)132 SpTuples<IT,NT>::SpTuples (int64_t size, IT nRow, IT nCol, StackEntry<NT, std::pair<IT,IT> > * & multstack)
133 :m(nRow), n(nCol), nnz(size)
134 {
135 isOperatorNew = false;
136 if(nnz > 0)
137 {
138 tuples = new std::tuple<IT, IT, NT>[nnz];
139 }
140 for(int64_t i=0; i<nnz; ++i)
141 {
142 colindex(i) = multstack[i].key.first;
143 rowindex(i) = multstack[i].key.second;
144 numvalue(i) = multstack[i].value;
145 }
146 delete [] multstack;
147 }
148
149
150 template <class IT,class NT>
~SpTuples()151 SpTuples<IT,NT>::~SpTuples()
152 {
153 if(nnz > 0)
154 {
155 if(isOperatorNew)
156 ::operator delete(tuples);
157 else
158 delete [] tuples;
159 }
160 }
161
162 /**
163 * Hint1: copy constructor (constructs a new object. i.e. this is NEVER called on an existing object)
164 * Hint2: Base's default constructor is called under the covers
165 * Normally Base's copy constructor should be invoked but it doesn't matter here as Base has no data members
166 */
167 template <class IT,class NT>
SpTuples(const SpTuples<IT,NT> & rhs)168 SpTuples<IT,NT>::SpTuples(const SpTuples<IT,NT> & rhs): m(rhs.m), n(rhs.n), nnz(rhs.nnz)
169 {
170 tuples = new std::tuple<IT, IT, NT>[nnz];
171 isOperatorNew = false;
172 for(IT i=0; i< nnz; ++i)
173 {
174 tuples[i] = rhs.tuples[i];
175 }
176 }
177
178 //! Constructor for converting SpDCCols matrix -> SpTuples
179 template <class IT,class NT>
SpTuples(const SpDCCols<IT,NT> & rhs)180 SpTuples<IT,NT>::SpTuples (const SpDCCols<IT,NT> & rhs): m(rhs.m), n(rhs.n), nnz(rhs.nnz)
181 {
182 if(nnz > 0)
183 {
184 FillTuples(rhs.dcsc);
185 }
186 isOperatorNew = false;
187 }
188
189 template <class IT,class NT>
FillTuples(Dcsc<IT,NT> * mydcsc)190 inline void SpTuples<IT,NT>::FillTuples (Dcsc<IT,NT> * mydcsc)
191 {
192 tuples = new std::tuple<IT, IT, NT>[nnz];
193 IT k = 0;
194 for(IT i = 0; i< mydcsc->nzc; ++i)
195 {
196 for(IT j = mydcsc->cp[i]; j< mydcsc->cp[i+1]; ++j)
197 {
198 colindex(k) = mydcsc->jc[i];
199 rowindex(k) = mydcsc->ir[j];
200 numvalue(k++) = mydcsc->numx[j];
201 }
202 }
203 }
204
205
206 // Hint1: The assignment operator (operates on an existing object)
207 // Hint2: The assignment operator is the only operator that is not inherited.
208 // Make sure that base class data are also updated during assignment
209 template <class IT,class NT>
operator =(const SpTuples<IT,NT> & rhs)210 SpTuples<IT,NT> & SpTuples<IT,NT>::operator=(const SpTuples<IT,NT> & rhs)
211 {
212 if(this != &rhs) // "this" pointer stores the address of the class instance
213 {
214 if(nnz > 0)
215 {
216 // make empty
217 if(isOperatorNew)
218 ::operator delete(tuples);
219 else
220 delete [] tuples;
221 }
222 m = rhs.m;
223 n = rhs.n;
224 nnz = rhs.nnz;
225 isOperatorNew = false;
226
227 if(nnz> 0)
228 {
229 tuples = new std::tuple<IT, IT, NT>[nnz];
230 for(IT i=0; i< nnz; ++i)
231 {
232 tuples[i] = rhs.tuples[i];
233 }
234 }
235 }
236 return *this;
237 }
238
239 /**
240 * \pre {The object is either column-sorted or row-sorted, either way the identical entries will be consecutive}
241 **/
242 template <class IT,class NT>
243 template <typename BINFUNC>
RemoveDuplicates(BINFUNC BinOp)244 void SpTuples<IT,NT>::RemoveDuplicates(BINFUNC BinOp)
245 {
246 if(nnz > 0)
247 {
248 std::vector< std::tuple<IT, IT, NT> > summed;
249 summed.push_back(tuples[0]);
250
251 for(IT i=1; i< nnz; ++i)
252 {
253 if((joker::get<0>(summed.back()) == joker::get<0>(tuples[i])) && (joker::get<1>(summed.back()) == joker::get<1>(tuples[i])))
254 {
255 joker::get<2>(summed.back()) = BinOp(joker::get<2>(summed.back()), joker::get<2>(tuples[i]));
256 }
257 else
258 {
259 summed.push_back(tuples[i]);
260
261 }
262 }
263 if(isOperatorNew)
264 ::operator delete(tuples);
265 else
266 delete [] tuples;
267 tuples = new std::tuple<IT, IT, NT>[summed.size()];
268 isOperatorNew = false;
269 std::copy(summed.begin(), summed.end(), tuples);
270 nnz = summed.size();
271 }
272 }
273
274
275 //! Loads a triplet matrix from infile
276 //! \remarks Assumes matlab type indexing for the input (i.e. indices start from 1)
277 template <class IT,class NT>
getstream(std::ifstream & infile)278 std::ifstream& SpTuples<IT,NT>::getstream (std::ifstream& infile)
279 {
280 std::cout << "Getting... SpTuples" << std::endl;
281 IT cnz = 0;
282 if (infile.is_open())
283 {
284 while ( (!infile.eof()) && cnz < nnz)
285 {
286 infile >> rowindex(cnz) >> colindex(cnz) >> numvalue(cnz); // row-col-value
287
288 rowindex(cnz) --;
289 colindex(cnz) --;
290
291 if((rowindex(cnz) > m) || (colindex(cnz) > n))
292 {
293 std::cerr << "supplied matrix indices are beyond specified boundaries, aborting..." << std::endl;
294 }
295 ++cnz;
296 }
297 assert(nnz == cnz);
298 }
299 else
300 {
301 std::cerr << "input file is not open!" << std::endl;
302 }
303 return infile;
304 }
305
306 //! Output to a triplets file
307 //! \remarks Uses matlab type indexing for the output (i.e. indices start from 1)
308 template <class IT,class NT>
putstream(std::ofstream & outfile) const309 std::ofstream& SpTuples<IT,NT>::putstream(std::ofstream& outfile) const
310 {
311 outfile << m <<"\t"<< n <<"\t"<< nnz<<std::endl;
312 for (IT i = 0; i < nnz; ++i)
313 {
314 outfile << rowindex(i)+1 <<"\t"<< colindex(i)+1 <<"\t"
315 << numvalue(i) << std::endl;
316 }
317 return outfile;
318 }
319
320 template <class IT,class NT>
PrintInfo()321 void SpTuples<IT,NT>::PrintInfo()
322 {
323 std::cout << "This is a SpTuples class" << std::endl;
324
325 std::cout << "m: " << m ;
326 std::cout << ", n: " << n ;
327 std::cout << ", nnz: "<< nnz << std::endl;
328
329 for(IT i=0; i< nnz; ++i)
330 {
331 if(rowindex(i) < 0 || colindex(i) < 0)
332 {
333 std::cout << "Negative index at " << i << std::endl;
334 return;
335 }
336 else if(rowindex(i) >= m || colindex(i) >= n)
337 {
338 std::cout << "Index " << i << " too big with values (" << rowindex(i) << ","<< colindex(i) << ")" << std::endl;
339 }
340 }
341
342 if(m < 8 && n < 8) // small enough to print
343 {
344 NT ** A = SpHelper::allocate2D<NT>(m,n);
345 for(IT i=0; i< m; ++i)
346 for(IT j=0; j<n; ++j)
347 A[i][j] = 0.0;
348
349 for(IT i=0; i< nnz; ++i)
350 {
351 A[rowindex(i)][colindex(i)] = numvalue(i);
352 }
353 for(IT i=0; i< m; ++i)
354 {
355 for(IT j=0; j<n; ++j)
356 {
357 std::cout << std::setiosflags(std::ios::fixed) << std::setprecision(2) << A[i][j];
358 std::cout << " ";
359 }
360 std::cout << std::endl;
361 }
362 SpHelper::deallocate2D(A,m);
363 }
364 }
365
366 }
367