IVSparse  v1.0
A sparse matrix compression library.
IVCSC_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, uint8_t compressionLevel, bool columnMajor>
16  // delete the meta data
17  if (metadata != nullptr) {
18  delete[] metadata;
19  }
20 
21  if (data != nullptr) {
22  // free the data
23  for (size_t i = 0; i < outerDim; i++) {
24  if (data[i] != nullptr) {
25  free(data[i]);
26  }
27  }
28  free(data);
29  }
30 
31  if (endPointers != nullptr) {
32  free(endPointers);
33  }
34 }
35 
36 // Row and Column Constructor
37 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
39 
40  #ifdef IVSPARSE_DEBUG
41  // check that the number of rows and columns is greater than 0
42  assert(num_rows > 0 && num_cols > 0 &&
43  "The number of rows and columns must be greater than 0");
44  #endif
45 
46  // set class variables
47  if constexpr (columnMajor) {
48  innerDim = num_cols;
49  outerDim = num_rows;
50  } else {
51  innerDim = num_rows;
52  outerDim = num_cols;
53  }
54  numRows = num_rows;
55  numCols = num_cols;
56  encodeValueType();
57  index_t = sizeof(indexT);
58  nnz = 0;
59 
60  // allocate memory for the data
61  try {
62  data = (void **)malloc(outerDim * sizeof(void *));
63  endPointers = (void **)malloc(outerDim * sizeof(void *));
64  } catch (const std::exception &e) {
65  std::cerr << e.what() << '\n';
66  }
67 
68  // set all data and endpointers to the nullptr
69  for (size_t i = 0; i < outerDim; i++) {
70  data[i] = nullptr;
71  endPointers[i] = nullptr;
72  }
73 
74  // set the metadata
75  metadata = new uint32_t[NUM_META_DATA];
76  metadata[0] = compressionLevel;
77  metadata[1] = innerDim;
78  metadata[2] = outerDim;
79  metadata[3] = nnz;
80  metadata[4] = val_t;
81  metadata[5] = index_t;
82 
83  // calculate the compression size
84  calculateCompSize();
85 
86  // run the user checks
87  #ifdef IVSPARSE_DEBUG
88  userChecks();
89  #endif
90 }
91 
92 // Eigen Constructor
93 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
95 
96  // check if the matrix is empty
97  if (mat.nonZeros() == 0) {
98  val_t = 0;
99  index_t = 0;
100 
101  data = nullptr;
102  endPointers = nullptr;
103  metadata = nullptr;
104  return;
105  }
106 
107  // make sure the matrix is compressed
108  mat.makeCompressed();
109 
110  // get the number of rows and columns
111  numRows = mat.rows();
112  numCols = mat.cols();
113 
114  outerDim = columnMajor ? numCols : numRows;
115  innerDim = columnMajor ? numRows : numCols;
116 
117  // get the number of non-zero elements
118  nnz = mat.nonZeros();
119 
120  // call the compression function
121  compressCSC(mat.valuePtr(), mat.innerIndexPtr(), mat.outerIndexPtr());
122 }
123 
124 // Eigen Row Major Constructor
125 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
126 SparseMatrix<T, indexT, compressionLevel, columnMajor>::SparseMatrix(Eigen::SparseMatrix<T, Eigen::RowMajor> &mat) {
127 
128  // check if the matrix is empty
129  if (mat.nonZeros() == 0) {
130  val_t = 0;
131  index_t = 0;
132 
133  data = nullptr;
134  endPointers = nullptr;
135  metadata = nullptr;
136  return;
137  }
138 
139  // make sure the matrix is compressed
140  mat.makeCompressed();
141 
142  // get the number of rows and columns
143  numRows = mat.rows();
144  numCols = mat.cols();
145 
146  outerDim = numRows;
147  innerDim = numCols;
148 
149  // get the number of non-zero elements
150  nnz = mat.nonZeros();
151 
152  // call the compression function
153  compressCSC(mat.valuePtr(), mat.innerIndexPtr(), mat.outerIndexPtr());
154 }
155 
156 // Deep Copy Constructor
157 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
160 
161  *this = other;
162 }
163 
164 // Conversion Constructor
165 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
166 template <uint8_t otherCompressionLevel>
169 
170  // if already the right compression level
171  if constexpr (otherCompressionLevel == compressionLevel) {
172  *this = other;
173  return;
174  }
175 
176  // make a temporary matrix of the correct compression level
178 
179  // convert other to the right compression level
180  if constexpr (otherCompressionLevel == 1) {
181  temp = other.toIVCSC();
182  } else if constexpr (otherCompressionLevel == 2) {
183  temp = other.toIVCSC();
184  }
185 
186  // other should be the same compression level as this now
187  *this = temp;
188 
189  // run the user checks
190  #ifdef IVSPARSE_DEBUG
191  userChecks();
192  #endif
193 }
194 
195 // Raw CSC Constructor
196 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
197 template <typename T2, typename indexT2>
199  T2 *vals, indexT2 *innerIndices, indexT2 *outerPtr, uint32_t num_rows, uint32_t num_cols, uint32_t nnz) {
200 
201  #ifdef IVSPARSE_DEBUG
202  assert(num_rows > 0 && num_cols > 0 &&
203  "The number of rows and columns must be greater than 0");
204  assert(nnz > 0 && "The number of non-zero elements must be greater than 0");
205  #endif
206 
207  // see if the matrix is empty
208  if (nnz == 0) [[unlikely]] {
209  val_t = 0;
210  index_t = 0;
211  data = nullptr;
212  endPointers = nullptr;
213  metadata = nullptr;
214  return;
215  }
216 
217  // set class variables
218  if (columnMajor) {
219  innerDim = num_rows;
220  outerDim = num_cols;
221  } else {
222  innerDim = num_cols;
223  outerDim = num_rows;
224  }
225  numRows = num_rows;
226  numCols = num_cols;
227  this->nnz = nnz;
228 
229  // call the compression function
230  compressCSC(vals, innerIndices, outerPtr);
231 }
232 
233 // COO Constructor
234 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
235 template <typename T2, typename indexT2>
237  std::vector<std::tuple<indexT2, indexT2, T2>> &entries, uint32_t num_rows, uint32_t num_cols, uint32_t nnz) {
238 
239  #ifdef IVSPARSE_DEBUG
240  assert(num_rows > 0 && num_cols > 0 &&
241  "The number of rows and columns must be greater than 0");
242  assert(nnz > 0 && "The number of non-zero elements must be greater than 0");
243  #endif
244 
245  // see if the matrix is empty
246  if (nnz == 0) [[unlikely]] {
247  *this = SparseMatrix<T, indexT, compressionLevel, columnMajor>(num_rows, num_cols);
248  return;
249  }
250 
251  // set class variables
252  if (columnMajor) {
253  innerDim = num_rows;
254  outerDim = num_cols;
255  } else {
256  innerDim = num_cols;
257  outerDim = num_rows;
258  }
259 
260  numRows = num_rows;
261  numCols = num_cols;
262  this->nnz = nnz;
263  encodeValueType();
264  index_t = sizeof(indexT);
265 
266  metadata = new uint32_t[NUM_META_DATA];
267  metadata[0] = compressionLevel;
268  metadata[1] = innerDim;
269  metadata[2] = outerDim;
270  metadata[3] = nnz;
271  metadata[4] = val_t;
272  metadata[5] = index_t;
273 
274  // allocate memory for the data
275  try {
276  data = (void **)malloc(outerDim * sizeof(void *));
277  endPointers = (void **)malloc(outerDim * sizeof(void *));
278  } catch (std::bad_alloc &e) {
279  std::cerr << "Error: Could not allocate memory for IVSparse matrix"
280  << std::endl;
281  exit(1);
282  }
283 
284  // set all data and endpointers to the nullptr
285  for (size_t i = 0; i < outerDim; i++) {
286  data[i] = nullptr;
287  endPointers[i] = nullptr;
288  }
289 
290  // sort the tuples by first by column then by row
291  std::sort(entries.begin(), entries.end(),
292  [](const std::tuple<indexT2, indexT2, T2> &a,
293  const std::tuple<indexT2, indexT2, T2> &b) {
294  if (std::get<1>(a) == std::get<1>(b)) {
295  return std::get<0>(a) < std::get<0>(b);
296  } else {
297  return std::get<1>(a) < std::get<1>(b);
298  }
299  });
300 
301  std::map<T2, std::vector<indexT2>> maps[outerDim];
302 
303  // loop through the tuples
304  for (size_t i = 0; i < nnz; i++) {
305  // get the column
306  indexT2 row = std::get<0>(entries[i]);
307  indexT2 col = std::get<1>(entries[i]);
308  T2 val = std::get<2>(entries[i]);
309 
310  // check if the value is already in the map
311  if (maps[col].find(val) != maps[col].end()) {
312  // value found positive delta encode it
313  maps[col][val].push_back(row - maps[col][val][1]);
314 
315  // update the last index
316  maps[col][val][1] = row;
317 
318  // update the maximum delta
319  if (maps[col][val][maps[col][val].size() - 1] > maps[col][val][0])
320  maps[col][val][0] = maps[col][val][maps[col][val].size() - 1];
321  } else {
322  // value not found
323  maps[col][val] = std::vector<indexT2>{row};
324 
325  // add maximum delta and last index placeholders
326  maps[col][val].push_back(row);
327  maps[col][val].push_back(row);
328  }
329  } // end of loop through tuples
330 
331  // loop through the maps
332  #ifdef IVSPARSE_HAS_OPENMP
333  #pragma omp parallel for
334  #endif
335  for (uint32_t i = 0; i < outerDim; i++) {
336  size_t outerByteSize = 0;
337 
338  for (auto &pair : maps[i]) {
339  // change first value to be byte width of the maximum delta
340  pair.second[0] = byteWidth(pair.second[0]);
341 
342  // add the size of the run to the size of the column
343  //* value + index width + indices * index width + delimiter (index width)
344  outerByteSize += sizeof(T) + 1 +
345  (pair.second[0] * (pair.second.size() - 2)) +
346  pair.second[0];
347  }
348 
349  // if column is empty set the data and endpointer to nullptr
350  if (outerByteSize == 0) {
351  data[i] = nullptr;
352  endPointers[i] = nullptr;
353  continue;
354  }
355 
356  // allocate space for the column
357  try {
358  data[i] = malloc(outerByteSize);
359  } catch (std::bad_alloc &e) {
360  std::cout << "Error: " << e.what() << std::endl;
361  exit(1);
362  }
363 
364  // get a help pointer for moving through raw memory
365  void *helpPtr = data[i];
366 
367  // loop through the dictionary and write to memory
368  for (auto &pair : maps[i]) {
369  // Write the value to memory
370  *(T *)helpPtr = (T)pair.first;
371  helpPtr = (T *)helpPtr + 1;
372 
373  // also write the index width
374  *(uint8_t *)helpPtr = (uint8_t)pair.second[0];
375  helpPtr = (uint8_t *)helpPtr + 1;
376 
377  // loop through the indices and write them to memory
378  for (size_t k = 0; k < pair.second.size(); k++) {
379  // if compression level 3 skip the first two indices and cast the index
380  if (k == 0 || k == 1) {
381  continue;
382  }
383 
384  // create a type of the correct width
385  switch (pair.second[0]) {
386  case 1:
387  *(uint8_t *)helpPtr = (uint8_t)pair.second[k];
388  helpPtr = (uint8_t *)helpPtr + 1;
389  break;
390  case 2:
391  *(uint16_t *)helpPtr = (uint16_t)pair.second[k];
392  helpPtr = (uint16_t *)helpPtr + 1;
393  break;
394  case 4:
395  *(uint32_t *)helpPtr = (uint32_t)pair.second[k];
396  helpPtr = (uint32_t *)helpPtr + 1;
397  break;
398  case 8:
399  *(uint64_t *)helpPtr = (uint64_t)pair.second[k];
400  helpPtr = (uint64_t *)helpPtr + 1;
401  break;
402  }
403 
404  } // End of index loop
405 
406  // write a delimiter of the correct width
407  switch (pair.second[0]) {
408  case 1:
409  *(uint8_t *)helpPtr = (uint8_t)DELIM;
410  helpPtr = (uint8_t *)helpPtr + 1;
411  break;
412  case 2:
413  *(uint16_t *)helpPtr = (uint16_t)DELIM;
414  helpPtr = (uint16_t *)helpPtr + 1;
415  break;
416  case 4:
417  *(uint32_t *)helpPtr = (uint32_t)DELIM;
418  helpPtr = (uint32_t *)helpPtr + 1;
419  break;
420  case 8:
421  *(uint64_t *)helpPtr = (uint64_t)DELIM;
422  helpPtr = (uint64_t *)helpPtr + 1;
423  break;
424  }
425  // Set a pointer to the end of the data
426  endPointers[i] = helpPtr;
427 
428  } // End of dictionary loop
429  }
430 
431  // calculate the compression size
432  calculateCompSize();
433 }
434 
435 // IVSparse Vector Constructor
436 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
439 
440 
441  // if it is the matrix is empty and we need to construct it
442  if (columnMajor) {
443  numRows = vec.getLength();
444  numCols = 1;
445  innerDim = numRows;
446  outerDim = numCols;
447  } else {
448  numRows = 1;
449  numCols = vec.getLength();
450  innerDim = numCols;
451  outerDim = numRows;
452  }
453 
454  nnz = vec.nonZeros();
455 
456  encodeValueType();
457  index_t = sizeof(indexT);
458 
459  // malloc the data
460  try {
461  data = (void **)malloc(sizeof(void *));
462  endPointers = (void **)malloc(sizeof(void *));
463  } catch (std::bad_alloc &e) {
464  throw std::bad_alloc();
465  }
466 
467  // malloc the vector
468  try {
469  // if vector is empty set the data to null
470  if (vec.begin() == vec.end()) {
471  data[0] = nullptr;
472  endPointers[0] = nullptr;
473  } else {
474  data[0] = malloc(vec.byteSize());
475  endPointers[0] = (char *)data[0] + vec.byteSize();
476  }
477  } catch (std::bad_alloc &e) {
478  throw std::bad_alloc();
479  }
480 
481  // copy the vector into the matrix if the vector is not empty
482  if (vec.begin() != vec.end()) memcpy(data[0], vec.begin(), vec.byteSize());
483 
484  // set the metadata
485  metadata = new uint32_t[NUM_META_DATA];
486  metadata[0] = compressionLevel;
487  metadata[1] = innerDim;
488  metadata[2] = outerDim;
489  metadata[3] = nnz;
490  metadata[4] = val_t;
491  metadata[5] = index_t;
492 
493  // run the user checks and calculate the compression size
494  calculateCompSize();
495  #ifdef IVSPARSE_DEBUG
496  userChecks();
497  #endif
498 }
499 
500 // Array of Vectors Constructor
501 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
504 
506 
507  for (size_t i = 1; i < vecs.size(); i++) {
508  temp.append(vecs[i]);
509  }
510 
511  *this = temp;
512 
513  calculateCompSize();
514 
515  // run the user checks
516  #ifdef IVSPARSE_DEBUG
517  userChecks();
518  #endif
519 }
520 
521 // File Constructor
522 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
524 
525  FILE *fp = fopen(filename, "rb");
526 
527  #ifdef IVSPARSE_DEBUG
528  if (fp == nullptr) {
529  throw std::runtime_error("Error: Could not open file");
530  }
531  #endif
532 
533  // read the metadata
534  metadata = new uint32_t[NUM_META_DATA];
535  fread(metadata, sizeof(uint32_t), NUM_META_DATA, fp);
536 
537  // set the matrix info
538  innerDim = metadata[1];
539  outerDim = metadata[2];
540  nnz = metadata[3];
541  val_t = metadata[4];
542  index_t = metadata[5];
543 
544  numRows = columnMajor ? innerDim : outerDim;
545  numCols = columnMajor ? outerDim : innerDim;
546 
547  #ifdef IVSPARSE_DEBUG
548  // if the compression level of the file is different than the compression
549  // level of the class
550  if (metadata[0] != compressionLevel) {
551  // throw an error
552  throw std::runtime_error(
553  "Error: Compression level of file does not match compression level of "
554  "class");
555  }
556  #endif
557 
558  // allocate the memory
559  try {
560  data = (void **)malloc(outerDim * sizeof(void *));
561  endPointers = (void **)malloc(outerDim * sizeof(void *));
562  } catch (std::bad_alloc &e) {
563  std::cerr << "Error: Could not allocate memory for IVSparse matrix"
564  << std::endl;
565  exit(1);
566  }
567 
568  // get the vector sizes and allocate memory
569  for (size_t i = 0; i < outerDim; i++) {
570  // get the size of the column
571  uint64_t size;
572  fread(&size, sizeof(uint64_t), 1, fp);
573 
574  // if the size is 0, set the data and endpointer to nullptr
575  if (size == 0) {
576  data[i] = nullptr;
577  endPointers[i] = nullptr;
578  continue;
579  }
580 
581  // malloc the column
582  try {
583  data[i] = malloc(size);
584  endPointers[i] = (char *)data[i] + size;
585  } catch (std::bad_alloc &e) {
586  throw std::bad_alloc();
587  }
588  }
589 
590  // read the data
591  for (size_t i = 0; i < outerDim; i++) {
592  fread(data[i], 1, (uint8_t *)endPointers[i] - (uint8_t *)data[i], fp);
593  }
594 
595  // close the file
596  fclose(fp);
597 
598  // calculate the compresssion size
599  calculateCompSize();
600 
601  // run the user checks
602  #ifdef IVSPARSE_DEBUG
603  userChecks();
604  #endif
605 
606 } // end of file constructor
607 
608 //* Private Constructors *//
609 
610 // Private Tranpose Constructor
611 template <typename T, typename indexT, uint8_t compressionLevel, bool columnMajor>
613  std::unordered_map<T, std::vector<indexT>> maps[], uint32_t num_rows, uint32_t num_cols) {
614 
615  // set class variables
616  if constexpr (columnMajor) {
617  innerDim = num_cols;
618  outerDim = num_rows;
619  } else {
620  innerDim = num_rows;
621  outerDim = num_cols;
622  }
623 
624  numRows = num_cols;
625  numCols = num_rows;
626  encodeValueType();
627  index_t = sizeof(indexT);
628 
629  // allocate memory for the data
630  try {
631  data = (void **)malloc(outerDim * sizeof(void *));
632  endPointers = (void **)malloc(outerDim * sizeof(void *));
633  } catch (const std::exception &e) {
634  std::cerr << e.what() << '\n';
635  }
636 
637  // set all data and endpointers to the nullptr
638  for (size_t i = 0; i < outerDim; i++) {
639  data[i] = nullptr;
640  endPointers[i] = nullptr;
641  }
642 
643  //* logic here
644 
645  // loop through the array
646  #ifdef IVSPARSE_HAS_OPENMP
647  #pragma omp parallel for
648  #endif
649  for (size_t i = 0; i < outerDim; i++) {
650  // check if the column is empty
651  if (maps[i].empty()) [[unlikely]] {
652  data[i] = nullptr;
653  endPointers[i] = nullptr;
654  continue;
655  }
656 
657  size_t byteSize = 0;
658 
659  // loop through the vectors of the map
660  for (auto &val : maps[i]) {
661  // add the size of the vector to the byteSize
662  byteSize += sizeof(T) + 1 +
663  (val.second[val.second.size() - 1] * (val.second.size() - 1) +
664  val.second[val.second.size() - 1]);
665  }
666 
667  // allocate memory for the vector
668  try {
669  data[i] = malloc(byteSize);
670  } catch (const std::exception &e) {
671  std::cerr << e.what() << '\n';
672  }
673 
674  // set the end pointer
675  endPointers[i] = (char *)data[i] + byteSize;
676 
677  // compressCSC the column
678  void *helpPtr = data[i];
679 
680  for (auto &val : maps[i]) {
681  nnz += val.second.size() - 1;
682 
683  // set the value
684  *(T *)helpPtr = val.first;
685  helpPtr = (char *)helpPtr + sizeof(T);
686  *(uint8_t *)helpPtr = (uint8_t)val.second[val.second.size() - 1];
687  helpPtr = (uint8_t *)helpPtr + 1;
688 
689  // write the indices
690  for (size_t k = 0; k < val.second.size(); k++) {
691  if (k == val.second.size() - 1) break;
692 
693  switch (val.second[val.second.size() - 1]) {
694  case 1:
695  *(uint8_t *)helpPtr = (uint8_t)val.second[k];
696  helpPtr = (uint8_t *)helpPtr + 1;
697  break;
698  case 2:
699  *(uint16_t *)helpPtr = (uint16_t)val.second[k];
700  helpPtr = (uint16_t *)helpPtr + 1;
701  break;
702  case 4:
703  *(uint32_t *)helpPtr = (uint32_t)val.second[k];
704  helpPtr = (uint32_t *)helpPtr + 1;
705  break;
706  case 8:
707  *(uint64_t *)helpPtr = (uint64_t)val.second[k];
708  helpPtr = (uint64_t *)helpPtr + 1;
709  break;
710  }
711  }
712 
713  // write a delimiter of the correct width
714  switch (val.second[val.second.size() - 1]) {
715  case 1:
716  *(uint8_t *)helpPtr = (uint8_t)DELIM;
717  helpPtr = (uint8_t *)helpPtr + 1;
718  break;
719  case 2:
720  *(uint16_t *)helpPtr = (uint16_t)DELIM;
721  helpPtr = (uint16_t *)helpPtr + 1;
722  break;
723  case 4:
724  *(uint32_t *)helpPtr = (uint32_t)DELIM;
725  helpPtr = (uint32_t *)helpPtr + 1;
726  break;
727  case 8:
728  *(uint64_t *)helpPtr = (uint64_t)DELIM;
729  helpPtr = (uint64_t *)helpPtr + 1;
730  break;
731  }
732  }
733  }
734 
735  //* end logic
736 
737  metadata = new uint32_t[NUM_META_DATA];
738 
739  // Set the meta data
740  metadata[0] = compressionLevel;
741  metadata[1] = innerDim;
742  metadata[2] = outerDim;
743  metadata[3] = nnz;
744  metadata[4] = val_t;
745  metadata[5] = index_t;
746 
747  // calculate the compression size
748  calculateCompSize();
749 
750  // run the user checks
751  #ifdef IVSPARSE_DEBUG
752  userChecks();
753  #endif
754 
755 } // end of private transpose constructor
756 
757 } // namespace IVSparse
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: IVCSC_SparseMatrix.hpp:29
void append(typename SparseMatrix< T, indexT, compressionLevel, columnMajor >::Vector &vec)
Definition: IVCSC_Methods.hpp:245
SparseMatrix()
Definition: IVCSC_SparseMatrix.hpp:99
~SparseMatrix()
Destroy the Sparse Matrix object.
Definition: IVCSC_Constructors.hpp:15