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