16 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
21 assert(row < numRows && col < numCols &&
"Invalid row and column!");
22 assert(row >= 0 && col >= 0 &&
"Invalid row and column!");
25 return (*
this)(row, col);
29 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
35 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
39 assert(vec < outerDim && vec >= 0 &&
"Invalid vector!");
46 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
49 assert(vec < outerDim && vec >= 0 &&
"Invalid vector!");
56 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
59 assert(vec < outerDim && vec >= 0 &&
"Invalid vector!");
62 if (data[vec] == endPointers[vec]) {
65 return (
char *)endPointers[vec] - (
char *)data[vec];
71 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
75 FILE *fp = fopen(filename,
"wb");
78 fwrite(metadata, 1, NUM_META_DATA *
sizeof(uint32_t), fp);
81 for (uint32_t i = 0; i < outerDim; i++) {
82 uint64_t size = (uint8_t *)endPointers[i] - (uint8_t *)data[i];
83 fwrite(&size, 1,
sizeof(uint64_t), fp);
87 for (uint32_t i = 0; i < outerDim; i++) {
88 fwrite(data[i], 1, (
char *)endPointers[i] - (
char *)data[i], fp);
96 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
99 std::cout << std::endl;
100 std::cout <<
"IVSparse Matrix" << std::endl;
103 if (numRows < 100 && numCols < 100) {
105 for (uint32_t i = 0; i < numRows; i++) {
106 for (uint32_t j = 0; j < numCols; j++) {
107 std::cout << coeff(i, j) <<
" ";
109 std::cout << std::endl;
111 }
else if (numRows > 100 && numCols > 100) {
113 for (uint32_t i = 0; i < 100; i++) {
114 for (uint32_t j = 0; j < 100; j++) {
115 std::cout << coeff(i, j) <<
" ";
117 std::cout << std::endl;
121 std::cout << std::endl;
125 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
129 Eigen::SparseMatrix<T, columnMajor ? Eigen::ColMajor : Eigen::RowMajor> eigenMatrix(numRows, numCols);
132 for (uint32_t i = 0; i < outerDim; ++i) {
137 eigenMatrix.insert(it.row(), it.col()) = it.value();
142 eigenMatrix.makeCompressed();
152 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
155 if (compressionLevel == 2) {
162 T *values = (T *)malloc(nnz *
sizeof(T));
163 indexT *indices = (indexT *)malloc(nnz *
sizeof(indexT));
164 indexT *colPtrs = (indexT *)malloc((outerDim + 1) *
sizeof(indexT));
169 std::map<indexT, T> dict[outerDim];
172 for (uint32_t i = 0; i < outerDim; ++i) {
178 dict[i][it.getIndex()] = it.value();
182 colPtrs[i + 1] = colPtrs[i] + count;
188 for (uint32_t i = 0; i < outerDim; ++i) {
189 for (
auto &pair : dict[i]) {
190 values[count] = pair.second;
191 indices[count] = pair.first;
198 values, indices, colPtrs, numRows, numCols, nnz);
209 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
211 #ifdef IVSPARSE_DEBUG
213 assert(outerDim > 0 &&
"Cannot convert an empty matrix to an Eigen matrix!");
217 Eigen::SparseMatrix<T, columnMajor ? Eigen::ColMajor : Eigen::RowMajor> eigenMatrix(numRows, numCols);
220 for (uint32_t i = 0; i < outerDim; ++i) {
222 if (data[i] ==
nullptr) {
230 eigenMatrix.insert(it.row(), it.col()) = it.value();
235 eigenMatrix.makeCompressed();
244 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
248 #ifdef IVSPARSE_DEBUG
249 assert(vec.
getLength() == innerDim &&
"Vector must be the same size as the inner dimension!");
253 if (compSize == 0) [[unlikely]] {
257 if (vec.
begin() == vec.
end()) [[unlikely]] {
266 metadata[2] = outerDim;
270 data = (
void **)realloc(data, outerDim *
sizeof(
void *));
271 endPointers = (
void **)realloc(endPointers, outerDim *
sizeof(
void *));
272 }
catch (std::bad_alloc &e) {
273 throw std::bad_alloc();
277 data[outerDim - 1] =
nullptr;
278 endPointers[outerDim - 1] =
nullptr;
286 throw std::invalid_argument(
287 "The vector must be the same size as the outer dimension of the "
299 metadata[2] = outerDim;
304 data = (
void **)realloc(data, outerDim *
sizeof(
void *));
305 endPointers = (
void **)realloc(endPointers, outerDim *
sizeof(
void *));
306 }
catch (std::bad_alloc &e) {
307 throw std::bad_alloc();
312 data[outerDim - 1] = malloc(vec.
byteSize());
313 endPointers[outerDim - 1] = (
char *)data[outerDim - 1] + vec.
byteSize();
314 }
catch (std::bad_alloc &e) {
315 throw std::bad_alloc();
328 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
332 std::unordered_map<T, std::vector<indexT>> mapsT[innerDim];
335 for (uint32_t i = 0; i < outerDim; ++i) {
340 if constexpr (columnMajor) {
341 mapsT[it.row()][it.value()].push_back(it.col());
343 mapsT[it.col()][it.value()].push_back(it.row());
348 for (
auto &row : mapsT) {
349 for (
auto &col : row) {
351 size_t max = col.second[0];
354 for (uint32_t i = col.second.size() - 1; i > 0; --i) {
355 col.second[i] -= col.second[i - 1];
356 if ((
size_t)col.second[i] > max) max = col.second[i];
359 max = byteWidth(max);
361 col.second.push_back(max);
373 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
376 std::unordered_map<T, std::vector<indexT>> mapsT[innerDim];
379 for (uint32_t i = 0; i < outerDim; ++i) {
384 if constexpr (columnMajor) {
385 mapsT[it.row()][it.value()].push_back(it.col());
387 mapsT[it.col()][it.value()].push_back(it.row());
393 for (
auto &row : mapsT) {
394 for (
auto &col : row) {
396 size_t max = col.second[0];
399 for (uint32_t i = col.second.size() - 1; i > 0; --i) {
400 col.second[i] -= col.second[i - 1];
401 if ((
size_t)col.second[i] > max) max = col.second[i];
404 max = byteWidth(max);
406 col.second.push_back(max);
414 template <
typename T,
typename indexT, u
int8_t compressionLevel,
bool columnMajor>
415 std::vector<typename IVSparse::SparseMatrix<T, indexT, compressionLevel, columnMajor>::Vector>
418 #ifdef IVSPARSE_DEBUG
419 assert(start < outerDim && end <= outerDim && start < end &&
420 "Invalid start and end values!");
429 for (uint32_t i = start; i < end; ++i) {
435 vecs[i - start] = temp;
Definition: IVCSC_Iterator.hpp:25
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: CSC_SparseMatrix.hpp:24
Definition: VCSC_SparseMatrix.hpp:21
Definition: IVCSC_SparseMatrix.hpp:29
IVSparse::SparseMatrix< T, indexT, 2, columnMajor > toVCSC()
Definition: IVCSC_Methods.hpp:153
IVSparse::SparseMatrix< T, indexT, compressionLevel, columnMajor >::Vector getVector(uint32_t vec)
Definition: IVCSC_Methods.hpp:47
IVSparse::SparseMatrix< T, indexT, compressionLevel, columnMajor > transpose()
Definition: IVCSC_Methods.hpp:329
void append(typename SparseMatrix< T, indexT, compressionLevel, columnMajor >::Vector &vec)
Definition: IVCSC_Methods.hpp:245
IVSparse::SparseMatrix< T, indexT, 1, columnMajor > toCSC()
Definition: IVCSC_Methods.hpp:126
T coeff(uint32_t row, uint32_t col)
Definition: IVCSC_Methods.hpp:17
void write(const char *filename)
Definition: IVCSC_Methods.hpp:72
void print()
Definition: IVCSC_Methods.hpp:97
Eigen::SparseMatrix< T, columnMajor ? Eigen::ColMajor :Eigen::RowMajor > toEigen()
Definition: IVCSC_Methods.hpp:210
bool isColumnMajor() const
Definition: IVCSC_Methods.hpp:30
void * vectorPointer(uint32_t vec)
Definition: IVCSC_Methods.hpp:36
std::vector< typename IVSparse::SparseMatrix< T, indexT, compressionLevel, columnMajor >::Vector > slice(uint32_t start, uint32_t end)
Definition: IVCSC_Methods.hpp:416
void inPlaceTranspose()
Definition: IVCSC_Methods.hpp:374
size_t getVectorSize(uint32_t vec) const
Definition: IVCSC_Methods.hpp:57