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