IVSparse  v1.0
A sparse matrix compression library.
VCSC_Methods.hpp
1 
9 #pragma once
10 
11 namespace IVSparse
12 {
13 
14  //* Getters *//
15 
16  // Gets the element stored at the given row and column
17  template <typename T, typename indexT, bool columnMajor>
18  T SparseMatrix<T, indexT, 2, columnMajor>::coeff(uint32_t row, uint32_t col) { return (*this)(row, col); }
19 
20  // Check for Column Major
21  template <typename T, typename indexT, bool columnMajor>
22  bool SparseMatrix<T, indexT, 2, columnMajor>::isColumnMajor() const { return columnMajor; }
23 
24  // get the values vector
25  template <typename T, typename indexT, bool columnMajor>
26  T *SparseMatrix<T, indexT, 2, columnMajor>::getValues(uint32_t vec) const { return values[vec]; }
27 
28  // get the counts vector
29  template <typename T, typename indexT, bool columnMajor>
30  indexT *SparseMatrix<T, indexT, 2, columnMajor>::getCounts(uint32_t vec) const { return counts[vec]; }
31 
32  // get the indices vector
33  template <typename T, typename indexT, bool columnMajor>
34  indexT *SparseMatrix<T, indexT, 2, columnMajor>::getIndices(uint32_t vec) const { return indices[vec]; }
35 
36  // get the number of unique values in a vector
37  template <typename T, typename indexT, bool columnMajor>
39  {
40  if (valueSizes == nullptr) { return 0; }
41  return valueSizes[vec];
42  }
43 
44  // get the number of indices in a vector
45  template <typename T, typename indexT, bool columnMajor>
47  {
48  if (indexSizes == nullptr) { return 0; }
49  return indexSizes[vec];
50  }
51 
52  // get the vector at the given index
53  template <typename T, typename indexT, bool columnMajor>
55 
56  //* Utility Methods *//
57 
58  // Writes the matrix to file
59  template <typename T, typename indexT, bool columnMajor>
61  {
62  // Open the file
63  FILE *fp = fopen(filename, "wb");
64 
65  // Write the metadata
66  fwrite(metadata, 1, NUM_META_DATA * sizeof(uint32_t), fp);
67 
68  // write the lengths of the vectors
69  for (uint32_t i = 0; i < outerDim; ++i) { fwrite(&valueSizes[i], 1, sizeof(indexT), fp); }
70  for (uint32_t i = 0; i < outerDim; ++i) { fwrite(&indexSizes[i], 1, sizeof(indexT), fp); }
71 
72  // write the values
73  for (uint32_t i = 0; i < outerDim; ++i) { fwrite(values[i], 1, valueSizes[i] * sizeof(T), fp); }
74 
75  // write the counts
76  for (uint32_t i = 0; i < outerDim; ++i) { fwrite(counts[i], 1, valueSizes[i] * sizeof(indexT), fp); }
77 
78  // write the indices
79  for (uint32_t i = 0; i < outerDim; ++i) { fwrite(indices[i], 1, indexSizes[i] * sizeof(indexT), fp); }
80 
81  // close the file
82  fclose(fp);
83  }
84 
85  // Prints the matrix dense to console
86  template <typename T, typename indexT, bool columnMajor>
88  {
89  std::cout << std::endl;
90  std::cout << "IVSparse Matrix" << std::endl;
91 
92  // if the matrix is less than 100 rows and columns print the whole thing
93  if (numRows < 100 && numCols < 100)
94  {
95  // print the matrix
96  for (uint32_t i = 0; i < numRows; i++)
97  {
98  for (uint32_t j = 0; j < numCols; j++)
99  {
100  std::cout << coeff(i, j) << " ";
101  }
102  std::cout << std::endl;
103  }
104  }
105  else if (numRows > 100 && numCols > 100)
106  {
107  // print the first 100 rows and columns
108  for (uint32_t i = 0; i < 100; i++)
109  {
110  for (uint32_t j = 0; j < 100; j++)
111  {
112  std::cout << coeff(i, j) << " ";
113  }
114  std::cout << std::endl;
115  }
116  }
117 
118  std::cout << std::endl;
119  }
120 
121  // Convert a IVCSC matrix to CSC
122  template <typename T, typename indexT, bool columnMajor>
124  {
125  // create a new sparse matrix
126  Eigen::SparseMatrix<T, columnMajor ? Eigen::ColMajor : Eigen::RowMajor> eigenMatrix(numRows, numCols);
127 
128  // iterate over the matrix
129  for (uint32_t i = 0; i < outerDim; ++i)
130  {
131  for (typename SparseMatrix<T, indexT, 2>::InnerIterator it(*this, i); it; ++it)
132  {
133  // add the value to the matrix
134  eigenMatrix.insert(it.row(), it.col()) = it.value();
135  }
136  }
137 
138  // finalize the matrix
139  eigenMatrix.makeCompressed();
140 
141  // make a CSC matrix
143 
144  // return the matrix
145  return CSCMatrix;
146  }
147 
148  // Convert a IVCSC matrix to a VCSC matrix
149  template <typename T, typename indexT, bool columnMajor>
151  {
152 
153  // make a pointer for the CSC pointers
154  T *values = (T *)malloc(nnz * sizeof(T));
155  indexT *indices = (indexT *)malloc(nnz * sizeof(indexT));
156  indexT *colPtrs = (indexT *)malloc((outerDim + 1) * sizeof(indexT));
157 
158  colPtrs[0] = 0;
159 
160  // make an array of ordered maps to hold the data
161  std::map<indexT, T> dict[outerDim];
162 
163  // iterate through the data using the iterator
164  for (uint32_t i = 0; i < outerDim; ++i)
165  {
166  size_t count = 0;
167 
168  for (typename SparseMatrix<T, indexT, 2>::InnerIterator it(*this, i); it; ++it)
169  {
170  dict[i][it.getIndex()] = it.value();
171  count++;
172  }
173  colPtrs[i + 1] = colPtrs[i] + count;
174  }
175  size_t count = 0;
176 
177  // loop through the dictionary and populate values and indices
178  for (uint32_t i = 0; i < outerDim; ++i)
179  {
180  for (auto &pair : dict[i])
181  {
182  values[count] = pair.second;
183  indices[count] = pair.first;
184  count++;
185  }
186  }
187 
188  // return a IVCSC matrix from the CSC vectors
189  IVSparse::SparseMatrix<T, indexT, 3, columnMajor> mat(values, indices, colPtrs, numRows, numCols, nnz);
190 
191  // free the CSC vectors
192  free(values);
193  free(indices);
194  free(colPtrs);
195 
196  return mat;
197  }
198 
199  // converts the csf matrix to an eigen one and returns it
200  template <typename T, typename indexT, bool columnMajor>
201  Eigen::SparseMatrix<T, columnMajor ? Eigen::ColMajor : Eigen::RowMajor> SparseMatrix<T, indexT, 2, columnMajor>::toEigen()
202  {
203 
204  #ifdef CSF_DEBUG
205  // assert that the matrix is not empty
206  assert(outerDim > 0 && "Cannot convert an empty matrix to an Eigen matrix!");
207  #endif
208 
209  // create a new sparse matrix
210  Eigen::SparseMatrix<T, columnMajor ? Eigen::ColMajor : Eigen::RowMajor> eigenMatrix(numRows, numCols);
211 
212  // iterate over the matrix
213  for (uint32_t i = 0; i < outerDim; ++i)
214  {
215  for (typename SparseMatrix<T, indexT, 2>::InnerIterator it(*this, i); it; ++it)
216  {
217  // add the value to the matrix
218  eigenMatrix.insert(it.row(), it.col()) = it.value();
219  }
220  }
221 
222  // finalize the matrix
223  eigenMatrix.makeCompressed();
224 
225  // return the matrix
226  return eigenMatrix;
227  }
228 
229  //* Conversion/Transformation Methods *//
230 
231  // appends a vector to the back of the storage order of the matrix
232  template <typename T, typename indexT, bool columnMajor>
234  {
235  // check if the matrix is empty
236  if (numRows < 1 && numCols < 1) [[unlikely]]
237  {
239  }
240  else
241  {
242  // check if the vector is empty, if so change the implementation details
243  if (vec.nonZeros() == 0)
244  {
245 
246  if (columnMajor)
247  {
248  numCols++;
249  }
250  else
251  {
252  numRows++;
253  }
254  outerDim++;
255 
256  // update metadata
257  metadata[2] = outerDim;
258 
259  // realloc the vectors
260  try
261  {
262  values = (T **)realloc(values, outerDim * sizeof(T *));
263  counts = (indexT **)realloc(counts, outerDim * sizeof(indexT *));
264  indices = (indexT **)realloc(indices, outerDim * sizeof(indexT *));
265  valueSizes = (indexT *)realloc(valueSizes, outerDim * sizeof(indexT));
266  indexSizes = (indexT *)realloc(indexSizes, outerDim * sizeof(indexT));
267  }
268  catch (std::bad_alloc &e)
269  {
270  std::cerr << "Error: " << e.what() << std::endl;
271  exit(1);
272  }
273 
274  // set the last vector to be empty
275  values[outerDim - 1] = nullptr;
276  counts[outerDim - 1] = nullptr;
277  indices[outerDim - 1] = nullptr;
278  valueSizes[outerDim - 1] = 0;
279  indexSizes[outerDim - 1] = 0;
280 
281  calculateCompSize();
282  return;
283  }
284  else
285  {
286 
287  // check that the vector is the correct size
288  if ((vec.getLength() != innerDim))
289  throw std::invalid_argument("The vector must be the same size as the outer dimension of the matrix!");
290 
291  outerDim++;
292  nnz += vec.nonZeros();
293 
294  if (columnMajor)
295  {
296  numCols++;
297  }
298  else
299  {
300  numRows++;
301  }
302 
303  // update metadata
304  metadata[2] = outerDim;
305  metadata[3] = nnz;
306 
307  // realloc the vectors
308  try
309  {
310  values = (T **)realloc(values, outerDim * sizeof(T *));
311  counts = (indexT **)realloc(counts, outerDim * sizeof(indexT *));
312  indices = (indexT **)realloc(indices, outerDim * sizeof(indexT *));
313  valueSizes = (indexT *)realloc(valueSizes, outerDim * sizeof(indexT));
314  indexSizes = (indexT *)realloc(indexSizes, outerDim * sizeof(indexT));
315  }
316  catch (std::bad_alloc &e)
317  {
318  std::cerr << "Error: " << e.what() << std::endl;
319  exit(1);
320  }
321 
322  // set the sizes of the new vector
323  valueSizes[outerDim - 1] = vec.uniqueVals();
324  indexSizes[outerDim - 1] = vec.nonZeros();
325 
326  // allocate the new vectors
327  try
328  {
329  values[outerDim - 1] = (T *)malloc(valueSizes[outerDim - 1] * sizeof(T));
330  counts[outerDim - 1] = (indexT *)malloc(sizeof(indexT) * valueSizes[outerDim - 1]);
331  indices[outerDim - 1] = (indexT *)malloc(indexSizes[outerDim - 1] * sizeof(indexT));
332  }
333  catch (std::bad_alloc &e)
334  {
335  std::cerr << "Error: " << e.what() << std::endl;
336  exit(1);
337  }
338 
339  // copy the data from the vector to the new vectors
340  memcpy(values[outerDim - 1], vec.getValues(), valueSizes[outerDim - 1] * sizeof(T));
341  memcpy(counts[outerDim - 1], vec.getCounts(), valueSizes[outerDim - 1] * sizeof(indexT));
342  memcpy(indices[outerDim - 1], vec.getIndices(), indexSizes[outerDim - 1] * sizeof(indexT));
343 
344  // update the compressed size
345  calculateCompSize();
346  }
347  }
348 
349  } // end append
350 
351  // tranposes the csf matrix
352  template <typename T, typename indexT, bool columnMajor>
354  {
355  // make a data structure to store the tranpose
356  std::unordered_map<T, std::vector<indexT>> mapsT[innerDim];
357 
358  // populate the transpose data structure
359  for (uint32_t i = 0; i < outerDim; ++i)
360  {
361  for (typename SparseMatrix<T, indexT, 2>::InnerIterator it(*this, i); it; ++it)
362  {
363  // add the value to the map
364  if constexpr (columnMajor)
365  {
366  mapsT[it.row()][it.value()].push_back(it.col());
367  }
368  else
369  {
370  mapsT[it.col()][it.value()].push_back(it.row());
371  }
372  }
373  }
374 
375  // create a new matrix passing in transposedMap
376  IVSparse::SparseMatrix<T, indexT, 2, columnMajor> temp(mapsT, numRows, numCols);
377 
378  // return the new matrix
379  return temp;
380  }
381 
382  // Transpose In Place Method
383  template <typename T, typename indexT, bool columnMajor>
385  {
386 
387  // make a data structure to store the tranpose
388  std::unordered_map<T, std::vector<indexT>> mapsT[innerDim];
389 
390  // populate the transpose data structure
391  for (uint32_t i = 0; i < outerDim; ++i)
392  {
393  for (typename SparseMatrix<T, indexT, 2>::InnerIterator it(*this, i); it; ++it)
394  {
395  // add the value to the map
396  if constexpr (columnMajor)
397  {
398  mapsT[it.row()][it.value()].push_back(it.col());
399  }
400  else
401  {
402  mapsT[it.col()][it.value()].push_back(it.row());
403  }
404  }
405  }
406 
407  // set this to the transposed matrix
408  *this = IVSparse::SparseMatrix<T, indexT, 2, columnMajor>(mapsT, numRows, numCols);
409  }
410 
411  // slice method that returns a vector of IVSparse vectors
412  template <typename T, typename indexT, bool columnMajor>
413  std::vector<typename IVSparse::SparseMatrix<T, indexT, 2, columnMajor>::Vector> SparseMatrix<T, indexT, 2, columnMajor>::slice(uint32_t start, uint32_t end)
414  {
415  #ifdef CSF_DEBUG
416  assert(start < outerDim && end <= outerDim && start < end && "Invalid start and end values!");
417  #endif
418 
419  // make a vector of IVSparse vectors
420  std::vector<typename IVSparse::SparseMatrix<T, indexT, 2, columnMajor>::Vector> vecs(end - start);
421 
422  // grab the vectors and add them to vecs
423  for (uint32_t i = start; i < end; ++i)
424  {
425  // make a temp vector
427 
428  // add the vector to vecs
429  vecs[i - start] = temp;
430  }
431 
432  // return the vector
433  return vecs;
434  }
435 
436 } // end namespace IVSparse
Definition: IVCSC_Iterator.hpp:25
Definition: CSC_SparseMatrix.hpp:24
Definition: VCSC_SparseMatrix.hpp:22
uint32_t nonZeros() const
Definition: IVSparse_Base_Methods.hpp:34
Definition: IVCSC_SparseMatrix.hpp:27
IVSparse::SparseMatrix< T, indexT, compressionLevel, columnMajor >::Vector getVector(uint32_t vec)
Definition: IVCSC_Methods.hpp:29
IVSparse::SparseMatrix< T, indexT, compressionLevel, columnMajor > transpose()
Definition: IVCSC_Methods.hpp:269
void append(typename SparseMatrix< T, indexT, compressionLevel, columnMajor >::Vector &vec)
Definition: IVCSC_Methods.hpp:202
IVSparse::SparseMatrix< T, indexT, 1, columnMajor > toCSC()
Definition: IVCSC_Methods.hpp:92
T coeff(uint32_t row, uint32_t col)
Definition: IVCSC_Methods.hpp:17
void write(const char *filename)
Definition: IVCSC_Methods.hpp:42
void print()
Definition: IVCSC_Methods.hpp:64
Eigen::SparseMatrix< T, columnMajor ? Eigen::ColMajor :Eigen::RowMajor > toEigen()
Definition: IVCSC_Methods.hpp:169
bool isColumnMajor() const
Definition: IVCSC_Methods.hpp:21
std::vector< typename IVSparse::SparseMatrix< T, indexT, compressionLevel, columnMajor >::Vector > slice(uint32_t start, uint32_t end)
Definition: IVCSC_Methods.hpp:350
void inPlaceTranspose()
Definition: IVCSC_Methods.hpp:311