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