IVSparse  v1.0
A sparse matrix compression library.
VCSC_Private_Methods.hpp
Go to the documentation of this file.
1 
9 #pragma once
10 
11 namespace IVSparse {
12 
13 // Encodes the value type of the matrix in a uint32_t
14 template <typename T, typename indexT, bool columnMajor>
15 void SparseMatrix<T, indexT, 2, columnMajor>::encodeValueType() {
16  uint8_t byte0 = sizeof(T);
17  uint8_t byte1 = std::is_floating_point<T>::value ? 1 : 0;
18  uint8_t byte2 = std::is_signed_v<T> ? 1 : 0;
19  uint8_t byte3 = columnMajor ? 1 : 0;
20 
21  val_t = (byte3 << 24) | (byte2 << 16) | (byte1 << 8) | byte0;
22 }
23 
24 // Checks if the value type is correct for the matrix
25 template <typename T, typename indexT, bool columnMajor>
26 void SparseMatrix<T, indexT, 2, columnMajor>::checkValueType() {
27  uint8_t byte0 = val_t & 0xFF;
28  uint8_t byte1 = (val_t >> 8) & 0xFF;
29  uint8_t byte2 = (val_t >> 16) & 0xFF;
30  uint8_t byte3 = (val_t >> 24) & 0xFF;
31  assert(byte0 == sizeof(T) && "Value type size does not match");
32  assert(byte1 == std::is_floating_point_v<T> &&
33  "Value type is not floating point");
34  assert(byte2 == std::is_signed_v<T> && "Value type is not signed");
35  assert(byte3 == columnMajor && "Major direction does not match");
36 }
37 
38 // performs some simple user checks on the matrices metadata
39 template <typename T, typename indexT, bool columnMajor>
40 void SparseMatrix<T, indexT, 2, columnMajor>::userChecks() {
41  assert((innerDim > 1 || outerDim > 1 || nnz > 1) &&
42  "The matrix must have at least one row, column, and nonzero value");
43  assert(std::is_floating_point<indexT>::value == false &&
44  "The index type must be a non-floating point type");
45  assert((std::is_arithmetic<T>::value && std::is_arithmetic<indexT>::value) &&
46  "The value and index types must be numeric types");
47  assert((std::is_same<indexT, bool>::value == false) &&
48  "The index type must not be bool");
49  assert((innerDim < std::numeric_limits<indexT>::max() &&
50  outerDim < std::numeric_limits<indexT>::max()) &&
51  "The number of rows and columns must be less than the maximum value "
52  "of the index type");
53  checkValueType();
54 }
55 
56 // Calculates the current byte size of the matrix in memory
57 template <typename T, typename indexT, bool columnMajor>
58 void SparseMatrix<T, indexT, 2, columnMajor>::calculateCompSize() {
59  // set compSize to zero
60  compSize = 0;
61 
62  // add the size of the metadata
63  compSize += META_DATA_SIZE;
64 
65  // add the performance vectors
66  compSize += (sizeof(T *) * outerDim); // values
67  compSize += (sizeof(indexT *) * outerDim); // counts
68  compSize += (sizeof(indexT *) * outerDim); // indices
69 
70  compSize += (sizeof(indexT) * outerDim); // valueSizes
71  compSize += (sizeof(indexT) * outerDim); // indexSizes
72  for (uint32_t i = 0; i < outerDim; i++) {
73  compSize += (sizeof(T) * valueSizes[i]); // values
74  compSize += (sizeof(indexT) * valueSizes[i]); // counts
75  compSize += (sizeof(indexT) * indexSizes[i]); // indices
76  }
77 }
78 
79 // Compression Algorithm for going from CSC to IVCSC
80 template <typename T, typename indexT, bool columnMajor>
81 template <typename T2, typename indexT2>
82 void SparseMatrix<T, indexT, 2, columnMajor>::compressCSC(T2 *vals, indexT2 *innerIndices, indexT2 *outerPointers) {
83 
84  // ---- Stage 1: Setup the Matrix ---- //
85 
86  // set the value and index types of the matrix
87  encodeValueType();
88  index_t = sizeof(indexT);
89 
90  // allocate space for metadata
91  metadata = new uint32_t[NUM_META_DATA];
92  metadata[0] = 2;
93  metadata[1] = innerDim;
94  metadata[2] = outerDim;
95  metadata[3] = nnz;
96  metadata[4] = val_t;
97  metadata[5] = index_t;
98 
99  // run the user checks on the metadata
100  #ifdef IVSPARSE_DEBUG
101  userChecks();
102  #endif
103 
104  // allocate space for the 2D Run lenngth encoded CSC matrix
105  try {
106  values = (T **)malloc(sizeof(T *) * outerDim);
107  counts = (indexT **)malloc(sizeof(indexT *) * outerDim);
108  indices = (indexT **)malloc(sizeof(indexT *) * outerDim);
109 
110  valueSizes = (indexT *)malloc(sizeof(indexT) * outerDim);
111  indexSizes = (indexT *)malloc(sizeof(indexT) * outerDim);
112  } catch (std::bad_alloc &e) {
113  std::cerr << "Error: Could not allocate memory for the matrix" << std::endl;
114  exit(1);
115  }
116 
117  // ---- Stage 2: Construct the Dictionary For Each Column ---- //
118 
119  // Loop through each column and construct a middle data structre for the matrix
120  #ifdef IVSPARSE_HAS_OPENMP
121  #pragma omp parallel for
122  #endif
123  for (uint32_t i = 0; i < outerDim; i++) {
124  // create the data structure to temporarily hold the data
125  std::map<T2, std::vector<indexT2>>
126  dict; // Key = value, Value = vector of indices
127 
128  // check if the current column is empty
129  if (outerPointers[i] == outerPointers[i + 1]) {
130  valueSizes[i] = 0;
131  indexSizes[i] = 0;
132 
133  values[i] = nullptr;
134  counts[i] = nullptr;
135  indices[i] = nullptr;
136  continue;
137  }
138 
139  // create a variable to hold the size of the column
140  size_t numIndices = 0;
141 
142  // loop through each value in the column and add it to dict
143  for (indexT2 j = outerPointers[i]; j < outerPointers[i + 1]; j++) {
144  // check if the value is already in the dictionary or not
145  if (dict.find(vals[j]) != dict.end()) {
146  // add the index
147  dict[vals[j]].push_back(innerIndices[j]);
148 
149  numIndices++;
150  } else {
151  // if value not already in the dictionary add it
152 
153  // create a new vector for the indices
154  dict[vals[j]] = std::vector<indexT2>{innerIndices[j]};
155 
156  numIndices++;
157  }
158 
159  } // end value loop
160 
161  // ---- Stage 3: Allocate Size of Column Data ---- //
162 
163  try {
164  values[i] = (T *)malloc(sizeof(T) * dict.size());
165  counts[i] = (indexT *)malloc(sizeof(indexT) * dict.size());
166  indices[i] = (indexT *)malloc(sizeof(indexT) * numIndices);
167  } catch (std::bad_alloc &e) {
168  std::cerr << "Error: Could not allocate memory for the matrix"
169  << std::endl;
170  exit(1);
171  }
172 
173  // set the size of the column
174  valueSizes[i] = dict.size();
175  indexSizes[i] = numIndices;
176  size_t performanceVecSize = 0;
177  size_t indexSize = 0;
178 
179  // ---- Stage 4: Populate the Column Data ---- //
180 
181  for (auto &pair : dict) {
182  values[i][performanceVecSize] = pair.first;
183  counts[i][performanceVecSize] = pair.second.size();
184 
185  for (uint32_t j = 0; j < pair.second.size(); j++) {
186  indices[i][indexSize] = pair.second[j];
187  indexSize++;
188  }
189  performanceVecSize++;
190  }
191 
192  } // end column loop
193 
194  calculateCompSize();
195 
196 } // end compressCSC
197 
198 } // end namespace IVSparse