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