From bfb1efe9aa916620669d1d583b4c497b41cc2389 Mon Sep 17 00:00:00 2001 From: Maksudul Alam Date: Thu, 15 May 2025 17:46:11 +0000 Subject: [PATCH 1/4] Debug Fix + GPU Code Fix --- .../IpReSolveSolverInterface.cpp | 40 ++- src/Algorithm/LinearSolvers/test-csr.h | 257 ++++++++++++++++++ 2 files changed, 288 insertions(+), 9 deletions(-) create mode 100644 src/Algorithm/LinearSolvers/test-csr.h diff --git a/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp b/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp index 81f832a8..f0d4b92a 100644 --- a/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp +++ b/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp @@ -203,12 +203,14 @@ ESymSolverStatus ReSolveSolverInterface::InitializeStructure(Index dim, Index no if (method_ == resolve_klu) { + printf("Resolve::KLU Setup\n"); workspace_CPU_ = new ReSolve::LinAlgWorkspaceCpu(); matrix_handler_ = new ReSolve::MatrixHandler(workspace_CPU_); vector_handler_ = new ReSolve::VectorHandler(workspace_CPU_); } #if RESOLVE_WITH_GPU + printf("Resolve::GPU Setup\n"); workspace_GPU_ = new workspace_type(); workspace_GPU_->initializeHandles(); @@ -218,17 +220,20 @@ ESymSolverStatus ReSolveSolverInterface::InitializeStructure(Index dim, Index no if (method_ == resolve_rf || method_ == resolve_rf_fgmres) { - resolve_Rf_ = new rf_solver(workspace_GPU_); + printf("Resolve::RF Setup\n"); + resolve_Rf_ = new rf_solver(); } # if RESOLVE_WITH_CUDA else if (method_ == resolve_glu) { + printf("Resolve::GLU Setup\n"); resolve_GLU_ = new ReSolve::LinSolverDirectCuSolverGLU(workspace_GPU_); } # endif if (method_ == resolve_rf_fgmres) { + printf("Resolve::FGMRES Setup\n"); GS_ = new ReSolve::GramSchmidt(vector_handler_, ReSolve::GramSchmidt::CGS2); resolve_FGMRES_ = new ReSolve::LinSolverIterativeFGMRES(matrix_handler_, vector_handler_, GS_); } @@ -283,7 +288,14 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index bool full_factor_done = false; // Get Data from CPU and update the A Matrix - A_->copyDataFrom(A_->getRowData(ReSolve::memory::HOST), A_->getColData(ReSolve::memory::HOST), A_->getValues(ReSolve::memory::HOST), ReSolve::memory::HOST, ReSolve::memory::DEVICE); + + A_->setDataPointers(const_cast(ia), const_cast(ja), const_cast(A_->getValues(ReSolve::memory::HOST)), ReSolve::memory::HOST); + +#if RESOLVE_WITH_GPU + if (method_ == resolve_rf || method_ == resolve_rf_fgmres || method_ == resolve_glu) { + A_->syncData(ReSolve::memory::DEVICE); + } +#endif // FACTORIZE @@ -294,10 +306,8 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index // re_factorize_ = true; } - if (factorize_ && (new_matrix || re_factorize_)) - { - // printf("Iteration: %d: Performing KLU Factorization\n", n_iteration_); - + if( n_iteration_ == 0){ + printf("First Iteration: %d: Performing KLU Factorization\n", n_iteration_); // Symbolic Factorization if (HaveIpData()) { @@ -314,6 +324,11 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index { IpData().TimingStats().LinearSystemSymbolicFactorization().End(); } + } + + if (factorize_ && (new_matrix || re_factorize_)) + { + // printf("Iteration: %d: Performing KLU Factorization\n", n_iteration_); // perform the factorization if (HaveIpData()) @@ -322,7 +337,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index } // First Factorization is always done by KLU - std::cout << "%" << n_iteration_ << "%" << "FULL FACTORIZATIOM" << std::endl; + std::cout << "%" << n_iteration_ << "%" << "KLU FULL FACTORIZATION" << std::endl; status = resolve_KLU_->factorize(); full_factor_done = true; @@ -397,6 +412,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index // Actual Refactorize if (method_ == resolve_glu) { + std::cout << "%" << n_iteration_ << "%" << "GLU->refactorize()" << std::endl; status = resolve_GLU_->refactorize(); if (status != 0) { @@ -405,6 +421,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index } else if (method_ == resolve_rf || method_ == resolve_rf_fgmres) { + std::cout << "%" << n_iteration_ << "%" << "RF->refactorize()" << std::endl; status_refactor = resolve_Rf_->refactorize(); if (status != 0) { @@ -414,6 +431,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index # else if (method_ == resolve_rf || method_ == resolve_rf_fgmres) { + std::cout << "%" << n_iteration_ << "%" << "RF->refactorize()" << std::endl; int status = resolve_Rf_->refactorize(); if (status != 0) { @@ -424,7 +442,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index if (method_ == resolve_klu) { - std::cout << "%" << n_iteration_ << "%" << "RE-FACTORIZATIOM" << std::endl; + std::cout << "%" << n_iteration_ << "%" << "KLU->refactorize()" << std::endl; status = resolve_KLU_->refactorize(); if (status != 0) { @@ -476,6 +494,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index // Copy rhs_vals to vec_rhs vec_rhs_->copyDataFrom(rhs_vals, ReSolve::memory::HOST, ReSolve::memory::HOST); status = resolve_KLU_->solve(vec_rhs_, vec_x_); + printf("Solving using KLU!\n"); if (status != 0) { std::cout << "KLU solve status: " << status << std::endl; @@ -489,10 +508,13 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index if (method_ == resolve_rf || method_ == resolve_rf_fgmres) { printf("Iteration: %d: Setting up %s\n", n_iteration_, method_.c_str()); + ReSolve::matrix::Csc* L_csc = (ReSolve::matrix::Csc*)resolve_KLU_->getLFactor(); ReSolve::matrix::Csc* U_csc = (ReSolve::matrix::Csc*)resolve_KLU_->getUFactor(); ReSolve::matrix::Csr* L = new ReSolve::matrix::Csr(L_csc->getNumRows(), L_csc->getNumColumns(), L_csc->getNnz()); ReSolve::matrix::Csr* U = new ReSolve::matrix::Csr(U_csc->getNumRows(), U_csc->getNumColumns(), U_csc->getNnz()); + L_csc->syncData(ReSolve::memory::DEVICE); + U_csc->syncData(ReSolve::memory::DEVICE); matrix_handler_->csc2csr(L_csc, L, ReSolve::memory::DEVICE); matrix_handler_->csc2csr(U_csc, U, ReSolve::memory::DEVICE); if (L == nullptr) @@ -615,7 +637,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index if (use_rcond_) { Number rcond_val = resolve_KLU_->getMatrixConditionNumber(); - printf("RCond: %12.8e\n", rcond_val); + //printf("RCond: %12.8e\n", rcond_val); if (rcond_val < rcond_val_) { if (full_factor_done) diff --git a/src/Algorithm/LinearSolvers/test-csr.h b/src/Algorithm/LinearSolvers/test-csr.h new file mode 100644 index 00000000..0ddeacad --- /dev/null +++ b/src/Algorithm/LinearSolvers/test-csr.h @@ -0,0 +1,257 @@ +#ifndef __TEST_CSR_H__ +#define __TEST_CSR_H__ + +#include +#include +#include +#include +#include +#include +#include + +using namespace std; + +namespace debug +{ + +struct CSR +{ + vector val; + vector col_ind; + vector row_ptr; +}; + +class Debug +{ +public: + Debug() + { + } + + template static T* load_b_data(std::string prefix, int iter, int& n) + { + std::fstream file; + + std::stringstream ss; + ss << prefix << "-" << iter << ".txt"; + std::string filename = ss.str(); + + file.open(filename, std::ios::in | std::fstream::binary); + file.read((char*)&n, sizeof(n)); + T* d = new T[n]; + + T p; + for (int i = 0; i < n; i++) + { + file.read((char*)&p, sizeof(T)); + d[i] = p; + } + + file.close(); + return d; + } + + template static void save_b_data(std::string prefix, int iter, int n, const T* data) + { + std::fstream file; + std::stringstream ss; + ss << prefix << "-" << iter << ".bin"; + std::string filename = ss.str(); + + file.open(filename, std::ios::out | std::fstream::binary); + + file.write((char*)&n, sizeof(n)); + // std::cout.precision(std::numeric_limits::max_digits10 - 1); + file.write((char*)(data), n * sizeof(T)); + file.close(); + } + + template static T* load_data(std::string prefix, int iter, int& n) + { + std::fstream file; + + std::stringstream ss; + ss << prefix << "-" << iter << ".txt"; + std::string filename = ss.str(); + + file.open(filename, std::ios::in); + file >> n; + T* d = new T[n]; + + T p; + for (int i = 0; i < n; i++) + { + file >> p; + d[i] = p; + } + + file.close(); + return d; + } + + template static void save_data(std::string prefix, int iter, int n, const T* data) + { + std::fstream file; + std::stringstream ss; + ss << prefix << "-" << iter << ".txt"; + std::string filename = ss.str(); + + file.open(filename, std::ios::out); + file << n << std::endl; + + // std::cout.precision(std::numeric_limits::max_digits10 - 1); + + for (int i = 0; i < n; i++) + { + // file << std::scientific << data[i] << std::endl; + file << data[i] << std::endl; + } + // file.write((char *)(data), n * sizeof(T)); + file.close(); + } + + static void save_coo_as_mm(std::string prefix, int iter, int dim, int nnz, const int* row_ptr, const int* col_ind, double* val) + { + printf("%s\t%d\n", __FILE__, __LINE__); + int cont = 0; + + std::stringstream ss; + ss << prefix << "-matrix-" << iter << ".mtx"; + std::string filename = ss.str(); + std::cout << "Filename: " << filename << std::endl; + + std::ofstream mat_file(filename); + + mat_file << "%%MatrixMarket matrix coordinate real general" << std::endl; + mat_file << "% ID: " << iter << std::endl; + mat_file << dim << "\t" << dim << "\t" << nnz << std::endl; + + int nnz_c = 0; + + for (int i = 0; i < nnz; i++) + { + mat_file << row_ptr[i] << "\t" << col_ind[i] << "\t" << val[i] << std::endl; + } + mat_file.close(); + printf("%d, %d==%d\n", dim, nnz_c, nnz); + } + + static void save_csc_as_csv(std::string prefix, int iter, int dim, int nnz, const int* row_ptr, const int* col_ind, double* val) + { + printf("%s\t%d\n", __FILE__, __LINE__); + int cont = 0; + + std::stringstream ss; + ss << prefix << "-matrix-" << iter << ".csv"; + std::string filename = ss.str(); + std::cout << "Filename: " << filename << std::endl; + + std::ofstream mat_file(filename); + + int nnz_c = 0; + + for (int i = 1; i <= dim; i++) + { + int row_start = row_ptr[i - 1]; + int row_end = row_ptr[i]; + + nnz_c += row_end - row_start; + // printf("%d, %d, %d\n", row_start, row_end, row_end - row_start); + + for (int jj = 0; jj < dim; jj++) + { + bool found = false; + for (int j = row_start; j < row_end; j++) + { + if (jj == col_ind[j]) + { + found = true; + break; + } + } + + if (found) + { + mat_file << val[cont] << ","; + cont++; + } + else + { + mat_file << "0,"; + } + } + mat_file << std::endl; + } + mat_file.close(); + printf("%d, %d==%d\n", dim, nnz_c, nnz); + } + + static void save_csc_as_mm(std::string prefix, int iter, int dim, int nnz, const int* row_ptr, const int* col_ind, double* val) + { + printf("%s\t%d\n", __FILE__, __LINE__); + int cont = 0; + + std::stringstream ss; + ss << prefix << "-matrix-" << iter << ".mtx"; + std::string filename = ss.str(); + std::cout << "Filename: " << filename << std::endl; + + std::ofstream mat_file(filename); + + mat_file << "%%MatrixMarket matrix coordinate real general" << std::endl; + mat_file << "% ID: " << iter << std::endl; + mat_file << dim << "\t" << dim << "\t" << nnz << std::endl; + + int nnz_c = 0; + + for (int i = 1; i <= dim; i++) + { + int row_start = row_ptr[i - 1]; + int row_end = row_ptr[i]; + + nnz_c += row_end - row_start; + // printf("%d, %d, %d\n", row_start, row_end, row_end - row_start); + + for (int j = row_start; j < row_end; j++) + { + // printf("%d ", col_ind[j]); + mat_file << i << "\t" << col_ind[j] + 1 << "\t" << val[cont] << std::endl; + cont++; + } + } + mat_file.close(); + printf("%d, %d==%d\n", dim, nnz_c, nnz); + } + + static void save_vec_as_csv(std::string prefix, int iter, int dim, int nnz, const int* row_ptr, const int* col_ind, double* val) + { + } + + static void save_vec_as_mm(std::string prefix, int iter, int dim, double* val) + { + printf("%s\t%d\n", __FILE__, __LINE__); + int cont = 0; + + std::stringstream ss; + ss << prefix << "-rhs-" << iter << ".mtx"; + std::string filename = ss.str(); + std::cout << "Filename: " << filename << std::endl; + + std::ofstream vec_file(filename); + + vec_file << "%%MatrixMarket matrix array real general" << std::endl; + vec_file << "% ID: " << iter << std::endl; + vec_file << dim << "\t" << 1 << std::endl; + + int nnz_c = 0; + + for (int i = 0; i < dim; i++) + { + vec_file << val[i] << std::endl; + } + vec_file.close(); + } +}; + +} // namespace debug +#endif From 6e243317ae29ec8247d874ca97d0fa08cd4dab1d Mon Sep 17 00:00:00 2001 From: Maksudul Alam Date: Thu, 31 Jul 2025 16:14:20 +0000 Subject: [PATCH 2/4] Cleanup --- .../IpReSolveSolverInterface.cpp | 28 +- .../IpReSolveSolverInterface.hpp | 4 - src/Algorithm/LinearSolvers/NVMLHelper.hpp | 70 ----- src/Algorithm/LinearSolvers/test-csr.h | 257 ------------------ 4 files changed, 9 insertions(+), 350 deletions(-) delete mode 100644 src/Algorithm/LinearSolvers/NVMLHelper.hpp delete mode 100644 src/Algorithm/LinearSolvers/test-csr.h diff --git a/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp b/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp index f0d4b92a..0819dfd9 100644 --- a/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp +++ b/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp @@ -28,7 +28,6 @@ ReSolveSolverInterface::ReSolveSolverInterface() : val_(NULL) printf("Resolve with GPU.\n"); # if RESOLVE_WITH_CUDA printf("Resolve with CUDA.\n"); - NVMLHelper::getAvailableGPUMemory(); # else printf("Resolve with HIP.\n"); # endif @@ -39,10 +38,6 @@ ReSolveSolverInterface::ReSolveSolverInterface() : val_(NULL) ReSolveSolverInterface::~ReSolveSolverInterface() { - -// printf("Begin of Destructor\n"); -// NVMLHelper::getAvailableGPUMemory(); - DBG_START_METH("ReSolveSolverInterface::~ReSolveSolverInterface()", dbg_verbosity); delete[] val_; @@ -74,8 +69,6 @@ ReSolveSolverInterface::~ReSolveSolverInterface() delete vec_rhs_; delete vec_x_; delete A_; - -// printf("End of Destructor\n"); } void ReSolveSolverInterface::RegisterOptions(SmartPtr roptions) @@ -263,7 +256,6 @@ ESymSolverStatus ReSolveSolverInterface::InitializeStructure(Index dim, Index no vec_x_ = new ReSolve::vector::Vector(A_->getNumRows()); vec_x_->allocate(ReSolve::memory::HOST); // for KLU - // vec_x_->allocate(ReSolve::memory::DEVICE); } factorize_ = true; @@ -272,9 +264,6 @@ ESymSolverStatus ReSolveSolverInterface::InitializeStructure(Index dim, Index no initialized_ = true; pivtol_changed_ = false; -// printf("After Initialize\n"); -// NVMLHelper::getAvailableGPUMemory(); - return retval; } @@ -337,7 +326,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index } // First Factorization is always done by KLU - std::cout << "%" << n_iteration_ << "%" << "KLU FULL FACTORIZATION" << std::endl; + // std::cout << "%" << n_iteration_ << "%" << "KLU FULL FACTORIZATION" << std::endl; status = resolve_KLU_->factorize(); full_factor_done = true; @@ -412,7 +401,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index // Actual Refactorize if (method_ == resolve_glu) { - std::cout << "%" << n_iteration_ << "%" << "GLU->refactorize()" << std::endl; + //std::cout << "%" << n_iteration_ << "%" << "GLU->refactorize()" << std::endl; status = resolve_GLU_->refactorize(); if (status != 0) { @@ -421,7 +410,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index } else if (method_ == resolve_rf || method_ == resolve_rf_fgmres) { - std::cout << "%" << n_iteration_ << "%" << "RF->refactorize()" << std::endl; + //std::cout << "%" << n_iteration_ << "%" << "RF->refactorize()" << std::endl; status_refactor = resolve_Rf_->refactorize(); if (status != 0) { @@ -431,7 +420,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index # else if (method_ == resolve_rf || method_ == resolve_rf_fgmres) { - std::cout << "%" << n_iteration_ << "%" << "RF->refactorize()" << std::endl; + //std::cout << "%" << n_iteration_ << "%" << "RF->refactorize()" << std::endl; int status = resolve_Rf_->refactorize(); if (status != 0) { @@ -442,7 +431,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index if (method_ == resolve_klu) { - std::cout << "%" << n_iteration_ << "%" << "KLU->refactorize()" << std::endl; + //std::cout << "%" << n_iteration_ << "%" << "KLU->refactorize()" << std::endl; status = resolve_KLU_->refactorize(); if (status != 0) { @@ -470,7 +459,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index if (use_rcond_) { Number rcond_val = resolve_KLU_->getMatrixConditionNumber(); - printf("RCond: %12.8e\n", rcond_val); + //printf("RCond: %12.8e\n", rcond_val); if (rcond_val < rcond_val_) { if (full_factor_done) @@ -494,7 +483,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index // Copy rhs_vals to vec_rhs vec_rhs_->copyDataFrom(rhs_vals, ReSolve::memory::HOST, ReSolve::memory::HOST); status = resolve_KLU_->solve(vec_rhs_, vec_x_); - printf("Solving using KLU!\n"); + //printf("Solving using KLU!\n"); if (status != 0) { std::cout << "KLU solve status: " << status << std::endl; @@ -507,7 +496,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index # if RESOLVE_WITH_CUDA if (method_ == resolve_rf || method_ == resolve_rf_fgmres) { - printf("Iteration: %d: Setting up %s\n", n_iteration_, method_.c_str()); + printf("CUDA: Iteration: %d: Setting up %s\n", n_iteration_, method_.c_str()); ReSolve::matrix::Csc* L_csc = (ReSolve::matrix::Csc*)resolve_KLU_->getLFactor(); ReSolve::matrix::Csc* U_csc = (ReSolve::matrix::Csc*)resolve_KLU_->getUFactor(); @@ -538,6 +527,7 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index if (method_ == resolve_rf || method_ == resolve_rf_fgmres) { + printf("HIP: Iteration: %d: Setting up %s\n", n_iteration_, method_.c_str()); ReSolve::matrix::Csc* L = (ReSolve::matrix::Csc*)resolve_KLU_->getLFactor(); ReSolve::matrix::Csc* U = (ReSolve::matrix::Csc*)resolve_KLU_->getUFactor(); ReSolve::index_type* P = resolve_KLU_->getPOrdering(); diff --git a/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.hpp b/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.hpp index ebcafa08..92d10bce 100644 --- a/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.hpp +++ b/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.hpp @@ -13,10 +13,6 @@ #include "IpoptConfig.h" -#if RESOLVE_WITH_CUDA -#include "NVMLHelper.hpp" -#endif - #include #include #include diff --git a/src/Algorithm/LinearSolvers/NVMLHelper.hpp b/src/Algorithm/LinearSolvers/NVMLHelper.hpp deleted file mode 100644 index bc10329e..00000000 --- a/src/Algorithm/LinearSolvers/NVMLHelper.hpp +++ /dev/null @@ -1,70 +0,0 @@ -#pragma once - -#include -#include - -class NVMLHelper -{ -public: - static int getAvailableGPUMemory() - { - nvmlReturn_t result; - unsigned int device_count, i; - nvmlDevice_t device; - - // Initialize NVML - result = nvmlInit(); - if (result != NVML_SUCCESS) - { - printf("Failed to initialize NVML: %s\n", nvmlErrorString(result)); - return 1; - } - - // Get the number of GPUs - result = nvmlDeviceGetCount(&device_count); - if (result != NVML_SUCCESS) - { - printf("Failed to get device count: %s\n", nvmlErrorString(result)); - nvmlShutdown(); - return 1; - } - - // Loop through each GPU and get its memory information - for (i = 0; i < device_count; i++) - { - result = nvmlDeviceGetHandleByIndex(i, &device); - if (result != NVML_SUCCESS) - { - printf("Failed to get handle for device %d: %s\n", i, nvmlErrorString(result)); - nvmlShutdown(); - return 1; - } - - // Get memory information - nvmlMemory_t memory_info; - result = nvmlDeviceGetMemoryInfo(device, &memory_info); - if (result != NVML_SUCCESS) - { - printf("Failed to get memory info for device %d: %s\n", i, nvmlErrorString(result)); - nvmlShutdown(); - return 1; - } - - // Print total, free, and used memory - printf("GPU %d:\n", i); - printf(" Total memory: %ld\n", memory_info.total / (1024 * 1024)); - printf(" Free memory: %ld\n", memory_info.free / (1024 * 1024)); - printf(" Used memory: %ld\n", memory_info.used / (1024 * 1024)); - } - - // Shutdown NVML - result = nvmlShutdown(); - if (result != NVML_SUCCESS) - { - printf("Failed to shutdown NVML: %s\n", nvmlErrorString(result)); - return 1; - } - - return 0; - } -}; \ No newline at end of file diff --git a/src/Algorithm/LinearSolvers/test-csr.h b/src/Algorithm/LinearSolvers/test-csr.h deleted file mode 100644 index 0ddeacad..00000000 --- a/src/Algorithm/LinearSolvers/test-csr.h +++ /dev/null @@ -1,257 +0,0 @@ -#ifndef __TEST_CSR_H__ -#define __TEST_CSR_H__ - -#include -#include -#include -#include -#include -#include -#include - -using namespace std; - -namespace debug -{ - -struct CSR -{ - vector val; - vector col_ind; - vector row_ptr; -}; - -class Debug -{ -public: - Debug() - { - } - - template static T* load_b_data(std::string prefix, int iter, int& n) - { - std::fstream file; - - std::stringstream ss; - ss << prefix << "-" << iter << ".txt"; - std::string filename = ss.str(); - - file.open(filename, std::ios::in | std::fstream::binary); - file.read((char*)&n, sizeof(n)); - T* d = new T[n]; - - T p; - for (int i = 0; i < n; i++) - { - file.read((char*)&p, sizeof(T)); - d[i] = p; - } - - file.close(); - return d; - } - - template static void save_b_data(std::string prefix, int iter, int n, const T* data) - { - std::fstream file; - std::stringstream ss; - ss << prefix << "-" << iter << ".bin"; - std::string filename = ss.str(); - - file.open(filename, std::ios::out | std::fstream::binary); - - file.write((char*)&n, sizeof(n)); - // std::cout.precision(std::numeric_limits::max_digits10 - 1); - file.write((char*)(data), n * sizeof(T)); - file.close(); - } - - template static T* load_data(std::string prefix, int iter, int& n) - { - std::fstream file; - - std::stringstream ss; - ss << prefix << "-" << iter << ".txt"; - std::string filename = ss.str(); - - file.open(filename, std::ios::in); - file >> n; - T* d = new T[n]; - - T p; - for (int i = 0; i < n; i++) - { - file >> p; - d[i] = p; - } - - file.close(); - return d; - } - - template static void save_data(std::string prefix, int iter, int n, const T* data) - { - std::fstream file; - std::stringstream ss; - ss << prefix << "-" << iter << ".txt"; - std::string filename = ss.str(); - - file.open(filename, std::ios::out); - file << n << std::endl; - - // std::cout.precision(std::numeric_limits::max_digits10 - 1); - - for (int i = 0; i < n; i++) - { - // file << std::scientific << data[i] << std::endl; - file << data[i] << std::endl; - } - // file.write((char *)(data), n * sizeof(T)); - file.close(); - } - - static void save_coo_as_mm(std::string prefix, int iter, int dim, int nnz, const int* row_ptr, const int* col_ind, double* val) - { - printf("%s\t%d\n", __FILE__, __LINE__); - int cont = 0; - - std::stringstream ss; - ss << prefix << "-matrix-" << iter << ".mtx"; - std::string filename = ss.str(); - std::cout << "Filename: " << filename << std::endl; - - std::ofstream mat_file(filename); - - mat_file << "%%MatrixMarket matrix coordinate real general" << std::endl; - mat_file << "% ID: " << iter << std::endl; - mat_file << dim << "\t" << dim << "\t" << nnz << std::endl; - - int nnz_c = 0; - - for (int i = 0; i < nnz; i++) - { - mat_file << row_ptr[i] << "\t" << col_ind[i] << "\t" << val[i] << std::endl; - } - mat_file.close(); - printf("%d, %d==%d\n", dim, nnz_c, nnz); - } - - static void save_csc_as_csv(std::string prefix, int iter, int dim, int nnz, const int* row_ptr, const int* col_ind, double* val) - { - printf("%s\t%d\n", __FILE__, __LINE__); - int cont = 0; - - std::stringstream ss; - ss << prefix << "-matrix-" << iter << ".csv"; - std::string filename = ss.str(); - std::cout << "Filename: " << filename << std::endl; - - std::ofstream mat_file(filename); - - int nnz_c = 0; - - for (int i = 1; i <= dim; i++) - { - int row_start = row_ptr[i - 1]; - int row_end = row_ptr[i]; - - nnz_c += row_end - row_start; - // printf("%d, %d, %d\n", row_start, row_end, row_end - row_start); - - for (int jj = 0; jj < dim; jj++) - { - bool found = false; - for (int j = row_start; j < row_end; j++) - { - if (jj == col_ind[j]) - { - found = true; - break; - } - } - - if (found) - { - mat_file << val[cont] << ","; - cont++; - } - else - { - mat_file << "0,"; - } - } - mat_file << std::endl; - } - mat_file.close(); - printf("%d, %d==%d\n", dim, nnz_c, nnz); - } - - static void save_csc_as_mm(std::string prefix, int iter, int dim, int nnz, const int* row_ptr, const int* col_ind, double* val) - { - printf("%s\t%d\n", __FILE__, __LINE__); - int cont = 0; - - std::stringstream ss; - ss << prefix << "-matrix-" << iter << ".mtx"; - std::string filename = ss.str(); - std::cout << "Filename: " << filename << std::endl; - - std::ofstream mat_file(filename); - - mat_file << "%%MatrixMarket matrix coordinate real general" << std::endl; - mat_file << "% ID: " << iter << std::endl; - mat_file << dim << "\t" << dim << "\t" << nnz << std::endl; - - int nnz_c = 0; - - for (int i = 1; i <= dim; i++) - { - int row_start = row_ptr[i - 1]; - int row_end = row_ptr[i]; - - nnz_c += row_end - row_start; - // printf("%d, %d, %d\n", row_start, row_end, row_end - row_start); - - for (int j = row_start; j < row_end; j++) - { - // printf("%d ", col_ind[j]); - mat_file << i << "\t" << col_ind[j] + 1 << "\t" << val[cont] << std::endl; - cont++; - } - } - mat_file.close(); - printf("%d, %d==%d\n", dim, nnz_c, nnz); - } - - static void save_vec_as_csv(std::string prefix, int iter, int dim, int nnz, const int* row_ptr, const int* col_ind, double* val) - { - } - - static void save_vec_as_mm(std::string prefix, int iter, int dim, double* val) - { - printf("%s\t%d\n", __FILE__, __LINE__); - int cont = 0; - - std::stringstream ss; - ss << prefix << "-rhs-" << iter << ".mtx"; - std::string filename = ss.str(); - std::cout << "Filename: " << filename << std::endl; - - std::ofstream vec_file(filename); - - vec_file << "%%MatrixMarket matrix array real general" << std::endl; - vec_file << "% ID: " << iter << std::endl; - vec_file << dim << "\t" << 1 << std::endl; - - int nnz_c = 0; - - for (int i = 0; i < dim; i++) - { - vec_file << val[i] << std::endl; - } - vec_file.close(); - } -}; - -} // namespace debug -#endif From 4706b2c5f2c59d6a617336d0be4db293953e400c Mon Sep 17 00:00:00 2001 From: Maksudul Alam Date: Thu, 31 Jul 2025 20:04:13 +0000 Subject: [PATCH 3/4] Added GramSchmidt include --- src/Algorithm/LinearSolvers/IpReSolveSolverInterface.hpp | 1 + 1 file changed, 1 insertion(+) diff --git a/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.hpp b/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.hpp index 92d10bce..1d01ef0d 100644 --- a/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.hpp +++ b/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.hpp @@ -20,6 +20,7 @@ #include #include #include +#include #include #include From 85767a62ea7b1c8ec35746b50f988b32bb7d5d81 Mon Sep 17 00:00:00 2001 From: Maksudul Alam Date: Thu, 7 Aug 2025 19:17:12 +0000 Subject: [PATCH 4/4] Added parameter for rf_solver --- src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp b/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp index 0819dfd9..dd750d38 100644 --- a/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp +++ b/src/Algorithm/LinearSolvers/IpReSolveSolverInterface.cpp @@ -214,7 +214,7 @@ ESymSolverStatus ReSolveSolverInterface::InitializeStructure(Index dim, Index no if (method_ == resolve_rf || method_ == resolve_rf_fgmres) { printf("Resolve::RF Setup\n"); - resolve_Rf_ = new rf_solver(); + resolve_Rf_ = new rf_solver(workspace_GPU_); } # if RESOLVE_WITH_CUDA else if (method_ == resolve_glu) @@ -281,9 +281,15 @@ ESymSolverStatus ReSolveSolverInterface::MultiSolve(bool new_matrix, const Index A_->setDataPointers(const_cast(ia), const_cast(ja), const_cast(A_->getValues(ReSolve::memory::HOST)), ReSolve::memory::HOST); #if RESOLVE_WITH_GPU +#if RESOLVE_WITH_CUDA if (method_ == resolve_rf || method_ == resolve_rf_fgmres || method_ == resolve_glu) { A_->syncData(ReSolve::memory::DEVICE); } +#else + if (method_ == resolve_rf || method_ == resolve_rf_fgmres) { + A_->syncData(ReSolve::memory::DEVICE); + } +#endif #endif // FACTORIZE