IVSparse  v1.0
A sparse matrix compression library.
CSC_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  if (vals != nullptr) {
17  free(vals);
18  }
19  if (innerIdx != nullptr) {
20  free(innerIdx);
21  }
22  if (outerPtr != nullptr) {
23  free(outerPtr);
24  }
25  if (metadata != nullptr) {
26  delete[] metadata;
27  }
28 }
29 
30 // Eigen Constructor
31 template <typename T, typename indexT, bool columnMajor>
33 
34  // Make sure the matrix is compressed before reading it in
35  mat.makeCompressed();
36 
37  // set the dimensions and nnz
38  innerDim = mat.innerSize();
39  outerDim = mat.outerSize();
40  numRows = mat.rows();
41  numCols = mat.cols();
42  nnz = mat.nonZeros();
43 
44  encodeValueType();
45  index_t = sizeof(indexT);
46 
47  // allocate the memory
48  try {
49  vals = (T *)malloc(nnz * sizeof(T));
50  innerIdx = (indexT *)malloc(nnz * sizeof(indexT));
51  outerPtr = (indexT *)malloc((outerDim + 1) * sizeof(indexT));
52  } catch (std::bad_alloc &e) {
53  std::cerr << "Allocation failed: " << e.what() << '\n';
54  }
55 
56  // set the metadata
57  metadata = new uint32_t[NUM_META_DATA];
58  metadata[0] = 1;
59  metadata[1] = innerDim;
60  metadata[2] = outerDim;
61  metadata[3] = nnz;
62  metadata[4] = val_t;
63  metadata[5] = index_t;
64 
65  // copy the data
66  memcpy(vals, mat.valuePtr(), sizeof(T) * nnz);
67  memcpy(innerIdx, mat.innerIndexPtr(), sizeof(indexT) * nnz);
68  memcpy(outerPtr, mat.outerIndexPtr(), sizeof(indexT) * (outerDim + 1));
69 
70  // calculate the compressed size and run the user checks
71  calculateCompSize();
72 
73  // run the user checks
74  #ifdef IVSPARSE_DEBUG
75  userChecks();
76  #endif
77 }
78 
79 // eigen sparse matrix constructor (row major)
80 template <typename T, typename indexT, bool columnMajor>
81 SparseMatrix<T, indexT, 1, columnMajor>::SparseMatrix(Eigen::SparseMatrix<T, Eigen::RowMajor> &other) {
82 
83  other.makeCompressed();
84 
85  // set the dimensions and nnz
86  innerDim = other.innerSize();
87  outerDim = other.outerSize();
88  numRows = other.rows();
89  numCols = other.cols();
90  nnz = other.nonZeros();
91 
92  encodeValueType();
93  index_t = sizeof(indexT);
94 
95  // allocate the memory
96  try {
97  vals = (T *)malloc(nnz * sizeof(T));
98  innerIdx = (indexT *)malloc(nnz * sizeof(indexT));
99  outerPtr = (indexT *)malloc((outerDim + 1) * sizeof(indexT));
100  } catch (std::bad_alloc &e) {
101  std::cerr << "Allocation failed: " << e.what() << '\n';
102  }
103 
104  // set the metadata
105  metadata = new uint32_t[NUM_META_DATA];
106  metadata[0] = 1;
107  metadata[1] = innerDim;
108  metadata[2] = outerDim;
109  metadata[3] = nnz;
110  metadata[4] = val_t;
111  metadata[5] = index_t;
112 
113  // copy the data
114  memcpy(vals, other.valuePtr(), sizeof(T) * nnz);
115  memcpy(innerIdx, other.innerIndexPtr(), sizeof(indexT) * nnz);
116  memcpy(outerPtr, other.outerIndexPtr(), sizeof(indexT) * (outerDim + 1));
117 
118  // calculate the compressed size and run the user checks
119  calculateCompSize();
120 
121  // run the user checks
122  #ifdef IVSPARSE_DEBUG
123  userChecks();
124  #endif
125 }
126 
127 // Deep Copy Constructor
128 template <typename T, typename indexT, bool columnMajor>
130  *this = other;
131 }
132 
133 // General Conversion Constructor
134 template <typename T, typename indexT, bool columnMajor>
135 template <uint8_t compressionLevel2>
137 
138  // if already CSC then just copy
139  if constexpr (compressionLevel2 == 1) {
140  *this = other;
141  return;
142  }
143 
144  // make a temporary CSC matrix
146 
147  if constexpr (compressionLevel2 == 2) {
148  temp = other.toCSC();
149  } else if constexpr (compressionLevel2 == 3) {
150  temp = other.toCSC();
151  }
152 
153  // copy the temporary matrix
154  *this = temp;
155 
156  // run the user checks
157  #ifdef IVSPARSE_DEBUG
158  userChecks();
159  #endif
160 }
161 
162 // Raw CSC Constructor
163 template <typename T, typename indexT, bool columnMajor>
164 template <typename T2, typename indexT2>
166 T2 *vals, indexT2 *innerIndices, indexT2 *outerPtr, uint32_t num_rows, uint32_t num_cols, uint32_t nnz) {
167 
168  #ifdef IVSPARSE_DEBUG
169  assert(vals != nullptr && "Values array is null");
170  assert(innerIndices != nullptr && "Inner indices array is null");
171  assert(outerPtr != nullptr && "Outer pointer array is null");
172  assert(num_rows > 0 && "Number of rows must be greater than 0");
173  assert(num_cols > 0 && "Number of columns must be greater than 0");
174  assert(nnz > 0 && "Number of nonzeros must be greater than 0");
175  #endif
176 
177  // set the dimensions and nnz
178  if constexpr (columnMajor) {
179  innerDim = num_rows;
180  outerDim = num_cols;
181  } else {
182  innerDim = num_cols;
183  outerDim = num_rows;
184  }
185 
186  numRows = num_rows;
187  numCols = num_cols;
188  this->nnz = nnz;
189 
190  encodeValueType();
191  index_t = sizeof(indexT);
192 
193  // allocate the memory
194  try {
195  this->vals = (T *)malloc(nnz * sizeof(T));
196  innerIdx = (indexT *)malloc(nnz * sizeof(indexT));
197  this->outerPtr = (indexT *)malloc((outerDim + 1) * sizeof(indexT));
198  } catch (std::bad_alloc &e) {
199  std::cerr << "Allocation failed: " << e.what() << '\n';
200  }
201 
202  // set the metadata
203  metadata = new uint32_t[NUM_META_DATA];
204  metadata[0] = 1;
205  metadata[1] = innerDim;
206  metadata[2] = outerDim;
207  metadata[3] = nnz;
208  metadata[4] = val_t;
209  metadata[5] = index_t;
210 
211  // copy the data
212  memcpy(this->vals, vals, sizeof(T) * nnz);
213  memcpy(innerIdx, innerIndices, sizeof(indexT) * nnz);
214  memcpy(this->outerPtr, outerPtr, sizeof(indexT) * (outerDim + 1));
215 
216  // calculate the compressed size and run the user checks
217  calculateCompSize();
218 
219  // run the user checks
220  #ifdef IVSPARSE_DEBUG
221  userChecks();
222  #endif
223 }
224 
225 // Vector Constructor
226 template <typename T, typename indexT, bool columnMajor>
228 
229  // set the dimensions and nnz
230  if constexpr (columnMajor) {
231  innerDim = vec.getLength();
232  outerDim = 1;
233  numRows = vec.getLength();
234  numCols = 1;
235  } else {
236  innerDim = 1;
237  outerDim = vec.getLength();
238  numRows = 1;
239  numCols = vec.getLength();
240  }
241 
242  nnz = vec.nonZeros();
243 
244  encodeValueType();
245  index_t = sizeof(indexT);
246 
247  metadata = new uint32_t[NUM_META_DATA];
248  metadata[0] = 1;
249  metadata[1] = innerDim;
250  metadata[2] = outerDim;
251  metadata[3] = nnz;
252  metadata[4] = val_t;
253  metadata[5] = index_t;
254 
255  // if the vector is empty, return
256  if (vec.byteSize() == 0) {
257  vals = nullptr;
258  innerIdx = nullptr;
259 
260  // allocate outerPtr
261  try {
262  outerPtr = (indexT *)malloc((outerDim + 1) * sizeof(indexT));
263  } catch (std::bad_alloc &e) {
264  std::cerr << "Allocation failed: " << e.what() << '\n';
265  }
266 
267  outerPtr[0] = 0;
268  return;
269  }
270 
271  // allocate the memory
272  try {
273  vals = (T *)malloc(nnz * sizeof(T));
274  innerIdx = (indexT *)malloc(nnz * sizeof(indexT));
275  outerPtr = (indexT *)malloc((outerDim + 1) * sizeof(indexT));
276  } catch (std::bad_alloc &e) {
277  std::cerr << "Allocation failed: " << e.what() << '\n';
278  }
279 
280  // copy the data and set the outerPtr
281  memcpy(vals, vec.getValues(), sizeof(T) * nnz);
282  memcpy(innerIdx, vec.getInnerIndices(), sizeof(indexT) * nnz);
283  outerPtr[0] = 0;
284  outerPtr[1] = nnz;
285 
286  // calculate the compressed size and run the user checks
287  calculateCompSize();
288 
289  #ifdef IVSPARSE_DEBUG
290  userChecks();
291  #endif
292 }
293 
294 // Array of Vectors Constructor
295 template <typename T, typename indexT, bool columnMajor>
298 
299  // Make a temporary CSC matrix
301 
302  // append on each vector in the array
303  for (size_t i = 1; i < vecs.size(); i++) { temp.append(vecs[i]); }
304 
305  // copy the temporary matrix
306  *this = temp;
307 
308  // run the user checks and calculate the compressed size
309  calculateCompSize();
310 
311  #ifdef IVSPARSE_DEBUG
312  userChecks();
313  #endif
314 }
315 
316 // File Constructor
317 template <typename T, typename indexT, bool columnMajor>
319 
320  FILE *fp = fopen(filename, "rb");
321 
322  #ifdef IVSPARSE_DEBUG
323  // make sure the file exists
324  if (fp == NULL) {
325  std::cerr << "Error: Could not open file " << filename << std::endl;
326  exit(1);
327  }
328  #endif
329 
330  // read the metadata and set the dimensions/nnz
331  metadata = new uint32_t[NUM_META_DATA];
332  fread(metadata, sizeof(uint32_t), NUM_META_DATA, fp);
333  innerDim = metadata[1];
334  outerDim = metadata[2];
335  nnz = metadata[3];
336  val_t = metadata[4];
337  index_t = metadata[5];
338 
339  if constexpr (columnMajor) {
340  numRows = innerDim;
341  numCols = outerDim;
342  } else {
343  numRows = outerDim;
344  numCols = innerDim;
345  }
346 
347  // allocate the memory
348  try {
349  vals = (T *)malloc(nnz * sizeof(T));
350  innerIdx = (indexT *)malloc(nnz * sizeof(indexT));
351  outerPtr = (indexT *)malloc((outerDim + 1) * sizeof(indexT));
352  } catch (std::bad_alloc &e) {
353  std::cerr << "Error: Could not allocate memory for IVSparse matrix"
354  << std::endl;
355  exit(1);
356  }
357 
358  // read the data
359  fread(vals, sizeof(T), nnz, fp);
360  fread(innerIdx, sizeof(indexT), nnz, fp);
361  fread(outerPtr, sizeof(indexT), outerDim + 1, fp);
362 
363  // close the file
364  fclose(fp);
365 
366  // run the user checks
367  #ifdef IVSPARSE_DEBUG
368  userChecks();
369  #endif
370 
371  calculateCompSize();
372 }
373 
374 // COO -> CSC constructor
375 template <typename T, typename indexT, bool columnMajor>
376 template <typename T2, typename indexT2>
378 std::vector<std::tuple<indexT2, indexT2, T2>> &entries, uint32_t num_rows, uint32_t num_cols, uint32_t nnz) {
379 
380  #ifdef IVSPARSE_DEBUG
381  assert(entries.size() > 0 && "Entries array is empty");
382  assert(num_rows > 0 && "Number of rows must be greater than 0");
383  assert(num_cols > 0 && "Number of columns must be greater than 0");
384  assert(nnz > 0 && "Number of nonzeros must be greater than 0");
385  #endif
386 
387  if constexpr (columnMajor) {
388  innerDim = num_rows;
389  outerDim = num_cols;
390  } else {
391  innerDim = num_cols;
392  outerDim = num_rows;
393  }
394 
395  numRows = num_rows;
396  numCols = num_cols;
397 
398  this->nnz = nnz;
399 
400  encodeValueType();
401  index_t = sizeof(indexT);
402 
403  try {
404  vals = (T *)malloc(nnz * sizeof(T));
405  innerIdx = (indexT *)malloc(nnz * sizeof(indexT));
406  outerPtr = (indexT *)calloc(outerDim + 1, sizeof(indexT));
407  } catch (std::bad_alloc &e) {
408  std::cerr << "Allocation failed: " << e.what() << '\n';
409  }
410 
411  metadata = new uint32_t[NUM_META_DATA];
412 
413  metadata[0] = 1;
414  metadata[1] = innerDim;
415  metadata[2] = outerDim;
416  metadata[3] = nnz;
417  metadata[4] = val_t;
418  metadata[5] = index_t;
419 
420  // sort the tuples by first by column then by row
421  std::sort(entries.begin(), entries.end(),
422  [](const std::tuple<indexT2, indexT2, T2> &a,
423  const std::tuple<indexT2, indexT2, T2> &b) {
424  if (std::get<1>(a) == std::get<1>(b)) {
425  return std::get<0>(a) < std::get<0>(b);
426  } else {
427  return std::get<1>(a) < std::get<1>(b);
428  }
429  });
430 
431  int OuterIndex = -1;
432  int count = 0;
433 
434  for (size_t i = 0; i < entries.size(); i++) {
435  vals[i] = std::get<2>(entries[i]);
436  innerIdx[i] = std::get<0>(entries[i]);
437  count++;
438 
439  if (std::get<1>(entries[i]) != OuterIndex) {
440  if (OuterIndex != -1) {
441  outerPtr[OuterIndex + 1] = count - 1;
442  }
443  OuterIndex = std::get<1>(entries[i]);
444  outerPtr[OuterIndex] = count - 1;
445  }
446  }
447 
448  outerPtr[OuterIndex + 1] = count;
449 
450  // run the user checks
451  #ifdef IVSPARSE_DEBUG
452  userChecks();
453  #endif
454 
455  calculateCompSize();
456 }
457 
458 } // namespace IVSparse
Definition: CSC_SparseMatrix.hpp:24
void append(typename SparseMatrix< T, indexT, 1, columnMajor >::Vector &vec)
Definition: CSC_Methods.hpp:211
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, 1, columnMajor > toCSC()
Definition: IVCSC_Methods.hpp:126
SparseMatrix()
Definition: IVCSC_SparseMatrix.hpp:99
~SparseMatrix()
Destroy the Sparse Matrix object.
Definition: IVCSC_Constructors.hpp:15