IVSparse  v1.0
A sparse matrix compression library.
VCSC_Constructors.hpp
Go to the documentation of this file.
1 
9 #pragma once
10 
11 namespace IVSparse {
12 
13 // Destructor
14 template <typename T, typename indexT, bool columnMajor>
16  // delete the meta data
17  if (metadata != nullptr) {
18  delete[] metadata;
19  }
20 
21  // delete the values
22  if (values != nullptr) {
23  for (size_t i = 0; i < outerDim; i++) {
24  if (values[i] != nullptr) {
25  free(values[i]);
26  }
27  }
28  free(values);
29  }
30  if (counts != nullptr) {
31  for (size_t i = 0; i < outerDim; i++) {
32  if (counts[i] != nullptr) {
33  free(counts[i]);
34  }
35  }
36  free(counts);
37  }
38  if (indices != nullptr) {
39  for (size_t i = 0; i < outerDim; i++) {
40  if (indices[i] != nullptr) {
41  free(indices[i]);
42  }
43  }
44  free(indices);
45  }
46 
47  if (valueSizes != nullptr) {
48  free(valueSizes);
49  }
50  if (indexSizes != nullptr) {
51  free(indexSizes);
52  }
53 }
54 
55 // Eigen Constructor
56 template <typename T, typename indexT, bool columnMajor>
58  // make sure the matrix is compressed
59  mat.makeCompressed();
60 
61  // get the number of rows and columns
62  numRows = mat.rows();
63  numCols = mat.cols();
64 
65  outerDim = columnMajor ? numCols : numRows;
66  innerDim = columnMajor ? numRows : numCols;
67 
68  // get the number of non-zero elements
69  nnz = mat.nonZeros();
70 
71  // call the compression function
72  compressCSC(mat.valuePtr(), mat.innerIndexPtr(), mat.outerIndexPtr());
73 }
74 
75 // Eigen Row Major Constructor
76 template <typename T, typename indexT, bool columnMajor>
77 SparseMatrix<T, indexT, 2, columnMajor>::SparseMatrix(Eigen::SparseMatrix<T, Eigen::RowMajor> &mat) {
78  // make sure the matrix is compressed
79  mat.makeCompressed();
80 
81  // get the number of rows and columns
82  numRows = mat.rows();
83  numCols = mat.cols();
84 
85  outerDim = numRows;
86  innerDim = numCols;
87 
88  // get the number of non-zero elements
89  nnz = mat.nonZeros();
90 
91  // call the compression function
92  compressCSC(mat.valuePtr(), mat.innerIndexPtr(), mat.outerIndexPtr());
93 }
94 
95 // Deep Copy Constructor
96 template <typename T, typename indexT, bool columnMajor>
98  *this = other;
99 }
100 
101 // Conversion Constructor
102 template <typename T, typename indexT, bool columnMajor>
103 template <uint8_t otherCompressionLevel>
105  // if already the right compression level
106  if constexpr (otherCompressionLevel == 2) {
107  *this = other;
108  return;
109  }
110 
111  // make a temporary matrix of the correct compression level
113 
114  // convert other to the right compression level
115  if constexpr (otherCompressionLevel == 1) {
116  temp = other.toVCSC();
117  } else if constexpr (otherCompressionLevel == 3) {
118  temp = other.toVCSC();
119  }
120 
121  // other should be the same compression level as this now
122  *this = temp;
123 
124  // run the user checks and calculate the compression size
125  calculateCompSize();
126 
127  #ifdef IVSPARSE_DEBUG
128  userChecks();
129  #endif
130 }
131 
132 // Raw CSC Constructor
133 template <typename T, typename indexT, bool columnMajor>
134 template <typename T2, typename indexT2>
136  T2 *vals, indexT2 *innerIndices, indexT2 *outerPtr, uint32_t num_rows, uint32_t num_cols, uint32_t nnz) {
137 
138  #ifdef IVSPARSE_DEBUG
139  assert(num_rows > 0 && num_cols > 0 && nnz > 0 &&
140  "Error: Matrix dimensions must be greater than 0");
141  assert(innerIndices != nullptr && outerPtr != nullptr && vals != nullptr &&
142  "Error: Pointers cannot be null");
143  #endif
144 
145 
146  // set the dimensions
147  if (columnMajor) {
148  innerDim = num_rows;
149  outerDim = num_cols;
150  } else {
151  innerDim = num_cols;
152  outerDim = num_rows;
153  }
154  numRows = num_rows;
155  numCols = num_cols;
156  this->nnz = nnz;
157 
158  // call the compression function
159  compressCSC(vals, innerIndices, outerPtr);
160 }
161 
162 // COO Constructor
163 template <typename T, typename indexT, bool columnMajor>
164 template <typename T2, typename indexT2>
166  std::vector<std::tuple<indexT2, indexT2, T2>> &entries, uint32_t num_rows, uint32_t num_cols, uint32_t nnz) {
167 
168  #ifdef IVSPARSE_DEBUG
169  assert(num_rows > 0 && num_cols > 0 && nnz > 0 &&
170  "Error: Matrix dimensions must be greater than 0");
171  #endif
172 
173  // see if the matrix is empty
174  if (nnz == 0) {
176  }
177 
178  // set the dimensions
179  if (columnMajor) {
180  innerDim = num_rows;
181  outerDim = num_cols;
182  } else {
183  innerDim = num_cols;
184  outerDim = num_rows;
185  }
186 
187  numRows = num_rows;
188  numCols = num_cols;
189  this->nnz = nnz;
190  encodeValueType();
191  index_t = sizeof(indexT);
192 
193  metadata = new uint32_t[NUM_META_DATA];
194  metadata[0] = 2;
195  metadata[1] = innerDim;
196  metadata[2] = outerDim;
197  metadata[3] = nnz;
198  metadata[4] = val_t;
199  metadata[5] = index_t;
200 
201  // allocate the vectors
202  try {
203  values = (T **)malloc(sizeof(T *) * outerDim);
204  counts = (indexT **)malloc(sizeof(indexT *) * outerDim);
205  indices = (indexT **)malloc(sizeof(indexT *) * outerDim);
206  valueSizes = (indexT *)malloc(sizeof(indexT) * outerDim);
207  indexSizes = (indexT *)malloc(sizeof(indexT) * outerDim);
208  } catch (std::bad_alloc &ba) {
209  std::cerr << "bad_alloc caught: " << ba.what() << '\n';
210  throw std::runtime_error("Error: Could not allocate memory");
211  }
212 
213  // sort the tuples by first by column then by row
214  std::sort(entries.begin(), entries.end(),
215  [](const std::tuple<indexT2, indexT2, T2> &a,
216  const std::tuple<indexT2, indexT2, T2> &b) {
217  if (std::get<1>(a) == std::get<1>(b)) {
218  return std::get<0>(a) < std::get<0>(b);
219  } else {
220  return std::get<1>(a) < std::get<1>(b);
221  }
222  });
223 
224  std::map<T2, std::vector<indexT2>> maps[outerDim];
225 
226  // loop through the tuples
227  for (size_t i = 0; i < nnz; i++) {
228  // get the column
229  indexT2 row = std::get<0>(entries[i]);
230  indexT2 col = std::get<1>(entries[i]);
231  T2 val = std::get<2>(entries[i]);
232 
233  // check if the value is already in the map
234  if (maps[col].find(val) != maps[col].end()) {
235  // value found positive delta encode it
236  maps[col][val].push_back(row);
237  } else {
238  // value not found
239  maps[col][val] = std::vector<indexT2>{row};
240  }
241 
242  } // end of loop through tuples
243 
244  // loop through the array
245  #ifdef IVSPARSE_HAS_OPENMP
246  #pragma omp parallel for
247  #endif
248  for (size_t i = 0; i < outerDim; i++) {
249  // check if the column is empty
250  if (maps[i].empty()) {
251  values[i] = nullptr;
252  counts[i] = nullptr;
253  indices[i] = nullptr;
254  valueSizes[i] = 0;
255  indexSizes[i] = 0;
256  continue;
257  }
258 
259  size_t performanceVecSize = 0;
260  size_t numInidces = 0;
261 
262  // loop through the vectors of the map
263  for (auto &val : maps[i]) {
264  performanceVecSize++;
265  numInidces += val.second.size();
266  }
267 
268  try {
269  values[i] = (T *)malloc(sizeof(T) * maps[i].size());
270  counts[i] = (indexT *)malloc(sizeof(indexT) * maps[i].size());
271  indices[i] = (indexT *)malloc(sizeof(indexT) * numInidces);
272  } catch (std::bad_alloc &ba) {
273  std::cerr << "bad_alloc caught: " << ba.what() << '\n';
274  throw std::runtime_error("Error: Could not allocate memory");
275  }
276 
277  valueSizes[i] = maps[i].size();
278  indexSizes[i] = numInidces;
279  nnz += numInidces;
280 
281  size_t index = 0;
282  size_t valIndex = 0;
283 
284  for (auto &val : maps[i]) {
285  values[i][valIndex] = val.first;
286  counts[i][valIndex] = val.second.size();
287 
288  for (auto &indexVal : val.second) {
289  indices[i][index] = indexVal;
290  index++;
291  }
292 
293  valIndex++;
294  }
295 
296  } // end of loop through the array
297 
298  // run the user checks and calculate the compression size
299  calculateCompSize();
300 
301  #ifdef IVSPARSE_DEBUG
302  userChecks();
303  #endif
304 }
305 
306 // IVSparse Vector Constructor
307 template <typename T, typename indexT, bool columnMajor>
310 
311  // Get the dimensions and metadata
312  if (columnMajor) {
313  numRows = vec.getLength();
314  numCols = 1;
315  innerDim = numRows;
316  outerDim = numCols;
317  } else {
318  numRows = 1;
319  numCols = vec.getLength();
320  innerDim = numCols;
321  outerDim = numRows;
322  }
323  nnz = vec.nonZeros();
324  encodeValueType();
325  index_t = sizeof(indexT);
326 
327  metadata = new uint32_t[NUM_META_DATA];
328  metadata[0] = 2;
329  metadata[1] = innerDim;
330  metadata[2] = outerDim;
331  metadata[3] = nnz;
332  metadata[4] = val_t;
333  metadata[5] = index_t;
334 
335  // allocate the vectors
336  try {
337  values = (T **)malloc(sizeof(T *));
338  counts = (indexT **)malloc(sizeof(indexT *));
339  indices = (indexT **)malloc(sizeof(indexT *));
340  valueSizes = (indexT *)malloc(sizeof(indexT));
341  indexSizes = (indexT *)malloc(sizeof(indexT));
342  } catch (std::bad_alloc &ba) {
343  std::cerr << "bad_alloc caught: " << ba.what() << '\n';
344  throw std::runtime_error("Error: Could not allocate memory");
345  }
346 
347  // check if the vector is empty
348  if (vec.byteSize() == 0) [[unlikely]] {
349  values[0] = nullptr;
350  counts[0] = nullptr;
351  indices[0] = nullptr;
352  valueSizes[0] = 0;
353  indexSizes[0] = 0;
354  return;
355  }
356 
357  // set the sizes
358  valueSizes[0] = vec.uniqueVals();
359  indexSizes[0] = vec.nonZeros();
360 
361  // allocate the memory for the one vector
362  try {
363  values[0] = (T *)malloc(sizeof(T) * valueSizes[0]);
364  counts[0] = (indexT *)malloc(sizeof(indexT) * valueSizes[0]);
365  indices[0] = (indexT *)malloc(sizeof(indexT) * indexSizes[0]);
366  } catch (std::bad_alloc &ba) {
367  std::cerr << "bad_alloc caught: " << ba.what() << '\n';
368  throw std::runtime_error("Error: Could not allocate memory");
369  }
370 
371  // copy the values
372  memcpy(values[0], vec.getValues(), sizeof(T) * valueSizes[0]);
373  memcpy(counts[0], vec.getCounts(), sizeof(indexT) * valueSizes[0]);
374  memcpy(indices[0], vec.getIndices(), sizeof(indexT) * indexSizes[0]);
375 
376  // run the user checks and calculate the compression size
377  calculateCompSize();
378  #ifdef IVSPARSE_DEBUG
379  userChecks();
380  #endif
381 } // end of IVSparse Vector Constructor
382 
383 // Array of Vectors Constructor
384 template <typename T, typename indexT, bool columnMajor>
386  std::vector<typename IVSparse::SparseMatrix<T, indexT, 2,columnMajor>::Vector> &vecs) {
387 
388  // Construct a one vector matrix to append to
390 
391  // append the rest of the vectors
392  for (size_t i = 1; i < vecs.size(); i++) {
393  temp.append(vecs[i]);
394  }
395 
396  // copy the temp matrix to this
397  *this = temp;
398 
399  // run the user checks and calculate the compression size
400  calculateCompSize();
401 
402  #ifdef IVSPARSE_DEBUG
403  userChecks();
404  #endif
405 }
406 
407 // File Constructor
408 template <typename T, typename indexT, bool columnMajor>
410  // open the file
411  FILE *fp = fopen(filename, "rb");
412 
413  #ifdef IVSPARSE_DEBUG
414  // check if the file was opened
415  if (fp == nullptr) {
416  throw std::runtime_error("Error: Could not open file");
417  }
418  #endif
419 
420  // read the metadata
421  metadata = new uint32_t[NUM_META_DATA];
422  fread(metadata, sizeof(uint32_t), NUM_META_DATA, fp);
423 
424  // set the matrix info
425  innerDim = metadata[1];
426  outerDim = metadata[2];
427  nnz = metadata[3];
428  val_t = metadata[4];
429  index_t = metadata[5];
430  numRows = columnMajor ? innerDim : outerDim;
431  numCols = columnMajor ? outerDim : innerDim;
432 
433  #ifdef IVSPARSE_DEBUG
434  // if the compression level of the file is different than the compression
435  // level of the class
436  if (metadata[0] != 2) {
437  // throw an error
438  throw std::runtime_error(
439  "Error: Compression level of file does not match compression level of "
440  "class");
441  }
442  #endif
443 
444  // allocate the vectors
445  try {
446  values = (T **)malloc(sizeof(T *) * outerDim);
447  counts = (indexT **)malloc(sizeof(indexT *) * outerDim);
448  indices = (indexT **)malloc(sizeof(indexT *) * outerDim);
449  valueSizes = (indexT *)malloc(sizeof(indexT) * outerDim);
450  indexSizes = (indexT *)malloc(sizeof(indexT) * outerDim);
451  } catch (std::bad_alloc &ba) {
452  std::cerr << "bad_alloc caught: " << ba.what() << '\n';
453  throw std::runtime_error("Error: Could not allocate memory");
454  }
455 
456  // read in the value sizes
457  for (size_t i = 0; i < outerDim; i++) {
458  fread(&valueSizes[i], sizeof(indexT), 1, fp);
459  }
460 
461  // read in the index sizes
462  for (size_t i = 0; i < outerDim; i++) {
463  fread(&indexSizes[i], sizeof(indexT), 1, fp);
464  }
465 
466  // read in the values
467  for (size_t i = 0; i < outerDim; i++) {
468  try {
469  values[i] = (T *)malloc(sizeof(T) * valueSizes[i]);
470  } catch (std::bad_alloc &ba) {
471  std::cerr << "bad_alloc caught: " << ba.what() << '\n';
472  throw std::runtime_error("Error: Could not allocate memory");
473  }
474  fread(values[i], sizeof(T), valueSizes[i], fp);
475  }
476 
477  // read in the counts
478  for (size_t i = 0; i < outerDim; i++) {
479  try {
480  counts[i] = (indexT *)malloc(sizeof(indexT) * valueSizes[i]);
481  } catch (std::bad_alloc &ba) {
482  std::cerr << "bad_alloc caught: " << ba.what() << '\n';
483  throw std::runtime_error("Error: Could not allocate memory");
484  }
485  fread(counts[i], sizeof(indexT), valueSizes[i], fp);
486  }
487 
488  // read in the indices
489  for (size_t i = 0; i < outerDim; i++) {
490  try {
491  indices[i] = (indexT *)malloc(sizeof(indexT) * indexSizes[i]);
492  } catch (std::bad_alloc &ba) {
493  std::cerr << "bad_alloc caught: " << ba.what() << '\n';
494  throw std::runtime_error("Error: Could not allocate memory");
495  }
496  fread(indices[i], sizeof(indexT), indexSizes[i], fp);
497  }
498 
499  // close the file
500  fclose(fp);
501 
502  // calculate the compresssion size
503  calculateCompSize();
504 
505  // run the user checks
506  #ifdef IVSPARSE_DEBUG
507  userChecks();
508  #endif
509 } // end of File Constructor
510 
511 //* Private Constructors *//
512 
513 // Private Tranpose Constructor
514 template <typename T, typename indexT, bool columnMajor>
516  std::unordered_map<T, std::vector<indexT>> maps[], uint32_t num_rows, uint32_t num_cols) {
517 
518  // set class variables
519  if constexpr (columnMajor) {
520  innerDim = num_cols;
521  outerDim = num_rows;
522  } else {
523  innerDim = num_rows;
524  outerDim = num_cols;
525  }
526 
527  numRows = num_cols;
528  numCols = num_rows;
529  encodeValueType();
530  index_t = sizeof(indexT);
531 
532  // allocate the vectors
533  try {
534  values = (T **)malloc(sizeof(T *) * outerDim);
535  counts = (indexT **)malloc(sizeof(indexT *) * outerDim);
536  indices = (indexT **)malloc(sizeof(indexT *) * outerDim);
537  valueSizes = (indexT *)malloc(sizeof(indexT) * outerDim);
538  indexSizes = (indexT *)malloc(sizeof(indexT) * outerDim);
539  } catch (std::bad_alloc &ba) {
540  std::cerr << "bad_alloc caught: " << ba.what() << '\n';
541  throw std::runtime_error("Error: Could not allocate memory");
542  }
543 
544  // loop through the array
545  #ifdef IVSPARSE_HAS_OPENMP
546  #pragma omp parallel for
547  #endif
548  for (size_t i = 0; i < outerDim; i++) {
549  // check if the column is empty
550  if (maps[i].empty()) [[unlikely]] {
551  values[i] = nullptr;
552  counts[i] = nullptr;
553  indices[i] = nullptr;
554  valueSizes[i] = 0;
555  indexSizes[i] = 0;
556  continue;
557  }
558 
559  size_t byteSize = 0;
560  size_t numInidces = 0;
561 
562  // loop through the vectors of the map
563  for (auto &val : maps[i]) {
564  // add the size of the vector to the byteSize
565  byteSize += (sizeof(indexT) * val.second.size());
566 
567  // add the size of the vector to the numIndices
568  numInidces += val.second.size();
569  }
570 
571  try {
572  values[i] = (T *)malloc(sizeof(T) * maps[i].size());
573  counts[i] = (indexT *)malloc(sizeof(indexT) * maps[i].size());
574  indices[i] = (indexT *)malloc(sizeof(indexT) * numInidces);
575  } catch (std::bad_alloc &ba) {
576  std::cerr << "bad_alloc caught: " << ba.what() << '\n';
577  throw std::runtime_error("Error: Could not allocate memory");
578  }
579 
580  valueSizes[i] = maps[i].size();
581  indexSizes[i] = numInidces;
582  nnz += numInidces;
583 
584  size_t index = 0;
585  size_t valIndex = 0;
586 
587  for (auto &val : maps[i]) {
588  values[i][valIndex] = val.first;
589  counts[i][valIndex] = val.second.size();
590 
591  for (auto &indexVal : val.second) {
592  indices[i][index] = indexVal;
593  index++;
594  }
595 
596  valIndex++;
597  }
598 
599  } // end of loop through the array
600 
601  // set the metadata
602  metadata = new uint32_t[NUM_META_DATA];
603  metadata[0] = 2;
604  metadata[1] = innerDim;
605  metadata[2] = outerDim;
606  metadata[3] = nnz;
607  metadata[4] = val_t;
608  metadata[5] = index_t;
609 
610  // run the user checks and calculate the compression size
611  calculateCompSize();
612 
613  #ifdef IVSPARSE_DEBUG
614  userChecks();
615  #endif
616 } // end of Private Tranpose Constructor
617 
618 } // namespace IVSparse
Definition: VCSC_SparseMatrix.hpp:21
void append(typename SparseMatrix< T, indexT, 2, columnMajor >::Vector &vec)
Definition: VCSC_Methods.hpp:241
uint32_t nonZeros() const
Definition: IVSparse_Base_Methods.hpp:39
size_t byteSize() const
Definition: IVSparse_Base_Methods.hpp:42
Definition: IVCSC_SparseMatrix.hpp:29
IVSparse::SparseMatrix< T, indexT, 2, columnMajor > toVCSC()
Definition: IVCSC_Methods.hpp:153
SparseMatrix()
Definition: IVCSC_SparseMatrix.hpp:99
~SparseMatrix()
Destroy the Sparse Matrix object.
Definition: IVCSC_Constructors.hpp:15