From 9ecb9ff200731d7edf42ccdbcf368017441ecddb Mon Sep 17 00:00:00 2001 From: shakedregev Date: Fri, 29 Aug 2025 13:27:47 +0000 Subject: [PATCH 1/8] fixed examples --- examples/gpuRefactor.cpp | 6 +++--- examples/sysRefactor.cpp | 2 +- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/examples/gpuRefactor.cpp b/examples/gpuRefactor.cpp index ca11e7134..597238eff 100644 --- a/examples/gpuRefactor.cpp +++ b/examples/gpuRefactor.cpp @@ -267,8 +267,8 @@ int gpuRefactor(int argc, char* argv[]) if (i == 1) { // Extract factors and configure refactorization solver - matrix::Csc* L = (matrix::Csc*) KLU.getLFactor(); - matrix::Csc* U = (matrix::Csc*) KLU.getUFactor(); + matrix::Csr* L = (matrix::Csr*) KLU.getLFactorCsr(); + matrix::Csr* U = (matrix::Csr*) KLU.getUFactorCsr(); if (L == nullptr || U == nullptr) { std::cout << "Factor extraction from KLU failed!\n"; @@ -276,7 +276,7 @@ int gpuRefactor(int argc, char* argv[]) index_type* P = KLU.getPOrdering(); index_type* Q = KLU.getQOrdering(); - Rf.setup(A, L, U, P, Q, vec_rhs); + Rf.setupCsr(A, L, U, P, Q, vec_rhs); // Setup iterative refinement solver if (is_iterative_refinement) diff --git a/examples/sysRefactor.cpp b/examples/sysRefactor.cpp index 1cf3d48c5..66afbdc0c 100644 --- a/examples/sysRefactor.cpp +++ b/examples/sysRefactor.cpp @@ -81,7 +81,7 @@ int main(int argc, char* argv[]) else { std::cout << "Re::Solve is not built with support for " << opt->second; - std::cout << "backend.\n"; + std::cout << " backend.\n"; return 1; } From 52d23afcd34321f162c3be7e631ea56c921edeec Mon Sep 17 00:00:00 2001 From: shakedregev Date: Tue, 2 Sep 2025 14:47:57 +0000 Subject: [PATCH 2/8] kalmarek example --- examples/CMakeLists.txt | 5 + examples/gpuRefactorSame.cpp | 248 +++++++++++++++++++++++++++++++++++ 2 files changed, 253 insertions(+) create mode 100644 examples/gpuRefactorSame.cpp diff --git a/examples/CMakeLists.txt b/examples/CMakeLists.txt index e97ec0531..d5eee7a68 100644 --- a/examples/CMakeLists.txt +++ b/examples/CMakeLists.txt @@ -30,6 +30,11 @@ if(RESOLVE_USE_KLU) target_link_libraries(gpuRefactor.exe PRIVATE ReSolve) endif(RESOLVE_USE_GPU) + if(RESOLVE_USE_HIP) #HIP specific example + add_executable(gpuRefactorSame.exe gpuRefactorSame.cpp) + target_link_libraries(gpuRefactorSame.exe PRIVATE ReSolve) + endif(RESOLVE_USE_HIP) + # Create KLU+CUDA examples if(RESOLVE_USE_CUDA) diff --git a/examples/gpuRefactorSame.cpp b/examples/gpuRefactorSame.cpp new file mode 100644 index 000000000..c8c2ef271 --- /dev/null +++ b/examples/gpuRefactorSame.cpp @@ -0,0 +1,248 @@ +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include + +#include "ExampleHelper.hpp" + +int main() +{ + std::cout << "Simple ReSolve Example (following gpuRefactor pattern)\n"; + std::cout << "=====================================================\n"; + + // Matrix data in CSC format (converted from 1-based to 0-based indexing) + // Julia CSC data: colptr, rowval, nzval + std::vector csc_col_ptr = {0, 1, 2, 3, 6, 8, 10}; // Julia colptr - 1 + std::vector csc_row_ind = {0, 1, 2, 0, 1, 3, 0, 4, 1, 5}; // Julia rowval - 1 + std::vector csc_values = { + 0.00141521219668486, + 0.6984313028075056, + -10.0, + -0.8307014463420237, + 0.1846074461471301, + -10.0, + -0.8408959803799315, + -10.0, + -0.835721417925707, + -10.0 + }; + using namespace ReSolve::examples; + // Convert CSC to CSR manually for this example + int n = 6; // 6x6 matrix + int nnz = csc_values.size(); // 10 non-zeros + + std::cout << "Matrix size: " << n << "x" << n << "\n"; + std::cout << "Number of non-zeros: " << nnz << "\n"; + + // Convert CSC to CSR + std::vector csr_row_ptr(n + 1, 0); + std::vector csr_col_ind(nnz); + std::vector csr_values(nnz); + + // Count entries per row + for (int i = 0; i < nnz; ++i) { + csr_row_ptr[csc_row_ind[i] + 1]++; + } + + // Convert counts to pointers + for (int i = 1; i <= n; ++i) { + csr_row_ptr[i] += csr_row_ptr[i - 1]; + } + + // Fill CSR arrays + std::vector temp_row_ptr = csr_row_ptr; + for (int col = 0; col < n; ++col) { + for (int idx = csc_col_ptr[col]; idx < csc_col_ptr[col + 1]; ++idx) { + int row = csc_row_ind[idx]; + int pos = temp_row_ptr[row]++; + csr_col_ind[pos] = col; + csr_values[pos] = csc_values[idx]; + } + } + + std::cout << "Matrix data converted to CSR format\n"; + + // Create workspace + ReSolve::LinAlgWorkspaceHIP workspace; + workspace.initializeHandles(); + ExampleHelper helper(workspace); + + + // Create matrix and vector handlers + ReSolve::MatrixHandler matrix_handler(&workspace); + ReSolve::VectorHandler vector_handler(&workspace); + + // Direct solvers instantiation + ReSolve::LinSolverDirectKLU KLU; + ReSolve::LinSolverDirectRocSolverRf Rf(&workspace); + + // Create CSR matrix + ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(n, n, nnz, false, false); + A->allocateMatrixData(ReSolve::memory::HOST); + A->copyDataFrom(csr_row_ptr.data(), csr_col_ind.data(), csr_values.data(), + ReSolve::memory::HOST, ReSolve::memory::HOST); + A->allocateMatrixData(ReSolve::memory::DEVICE); + A->syncData(ReSolve::memory::DEVICE); + + std::cout << "CSR matrix created and synced to device\n"; + + // A->print(std::cout, 0); + + // Create vectors + ReSolve::vector::Vector* vec_rhs = new ReSolve::vector::Vector(n); + ReSolve::vector::Vector* vec_x = new ReSolve::vector::Vector(n); + + vec_rhs->allocate(ReSolve::memory::HOST); + vec_rhs->allocate(ReSolve::memory::DEVICE); + vec_x->allocate(ReSolve::memory::HOST); + vec_x->allocate(ReSolve::memory::DEVICE); + + // Set right-hand side vector + std::vector rhs_data = { + 0.04204480190421101, + 0.1868729984154292, + -0.1581138884334949, + 0.04939382906591484, + -1.0, + 0.1118034015563139 + }; + vec_rhs->copyDataFrom(rhs_data.data(), ReSolve::memory::HOST, ReSolve::memory::HOST); + vec_rhs->syncData(ReSolve::memory::DEVICE); + + std::cout << "Right-hand side vector set\n"; + + // Setup KLU solver + KLU.setup(A); + matrix_handler.setValuesChanged(true, ReSolve::memory::DEVICE); + + // Analysis (symbolic factorization) + int status = KLU.analyze(); + std::cout << "KLU analysis status: " << status << std::endl; + if (status != 0) { + std::cerr << "KLU analysis failed!" << std::endl; + return 1; + } + + // Numeric factorization + status = KLU.factorize(); + std::cout << "KLU factorization status: " << status << std::endl; + if (status != 0) { + std::cerr << "KLU factorization failed!" << std::endl; + return 1; + } + + // Triangular solve + status = KLU.solve(vec_rhs, vec_x); + std::cout << "KLU solve status: " << status << std::endl; + if (status != 0) { + std::cerr << "KLU solve failed!" << std::endl; + return 1; + } + + std::cout << "System solved successfully with KLU\n"; + + // Get solution + double* solution = vec_x->getData(ReSolve::memory::HOST); + + std::cout << "Solution vector:\n"; + std::cout << "["; + for (int i = 0; i < n; ++i) { + std::cout << solution[i] << (i < n - 1 ? ", " : ""); + } + std::cout << "]" << std::endl; + + helper.resetSystem(A, vec_rhs, vec_x); + helper.printShortSummary(); + // Extract factors and configure refactorization solver + // ReSolve::matrix::Csc* L = (ReSolve::matrix::Csc*) KLU.getLFactor(); + // ReSolve::matrix::Csc* U = (ReSolve::matrix::Csc*) KLU.getUFactor(); + // ReSolve::matrix::Csc* F = (ReSolve::matrix::Csc*) KLU.getFFactor(); + ReSolve::matrix::Csr* L = (ReSolve::matrix::Csr*) KLU.getLFactorCsr(); + ReSolve::matrix::Csr* U = (ReSolve::matrix::Csr*) KLU.getUFactorCsr(); + if (L == nullptr || U == nullptr) { + std::cout << "Factor extraction from KLU failed!\n"; + } else { + std::cout << "L and U factors extracted successfully\n"; + // Print L and U factors + std::cout << "L factor:\n"; + L->print(std::cout, 0); + std::cout << "U factor:\n"; + U->print(std::cout, 0); + if (F != nullptr) { + std::cout << "F factor:\n"; + F->print(std::cout, 0); + } else { + std::cout << "F factor is not available.\n"; + } + + ReSolve::index_type* P = KLU.getPOrdering(); + ReSolve::index_type* Q = KLU.getQOrdering(); + double* rhs_before = vec_rhs->getData(ReSolve::memory::HOST); + std::cout << "["; + for (int i = 0; i < n; ++i) { + std::cout << rhs_before[i] << (i < n - 1 ? ", " : ""); + } + std::cout << "]" << std::endl; + //Rf.setupCsr(A, L, U, P, Q, vec_rhs); + Rf.setup(A, L, U, P, Q, vec_rhs); + std::cout << "RocSolverRf setup completed\n"; + + // Test refactorization with the same matrix (in practice, you'd change matrix values) + std::cout << "\nTesting refactorization with same matrix...\n"; + + // Refactorize + status = Rf.refactorize(); + std::cout << "RocSolverRf refactorization status: " << status << std::endl; + if (status != 0) { + std::cerr << "Refactorization failed!" << std::endl; + } else { + // Solve with refactorization + status = Rf.solve(vec_rhs, vec_x); + std::cout << "RocSolverRf solve status: " << status << std::endl; + helper.resetSystem(A, vec_rhs, vec_x); + + if (status == 0) { + std::cout << "System solved successfully with refactorization\n"; + + // Get solution + vec_x->syncData(ReSolve::memory::HOST); + solution = vec_x->getData(ReSolve::memory::HOST); + + std::cout << "Refactorization solution vector:\n"; + std::cout << "["; + for (int i = 0; i < n; ++i) { + std::cout << solution[i] << (i < n - 1 ? ", " : ""); + } + std::cout << "]" << std::endl; + + // Get the rhs + double* rhs_solution = vec_rhs->getData(ReSolve::memory::HOST); + std::cout << "Right-hand side vector after solve:\n"; + std::cout << "["; + for (int i = 0; i < n; ++i) { + std::cout << rhs_solution[i] << (i < n - 1 ? ", " : ""); + } + std::cout << "]" << std::endl; + } + } + } + + // Cleanup + delete A; + delete vec_x; + delete vec_rhs; + + std::cout << "Example completed successfully!\n"; + return 0; +} \ No newline at end of file From ca8d179e51725c7f6e09ab938a02c7912c841c11 Mon Sep 17 00:00:00 2001 From: shakedregev Date: Tue, 2 Sep 2025 11:07:34 -0400 Subject: [PATCH 3/8] moved to experimental, kalmarek example works --- examples/CMakeLists.txt | 5 ----- examples/experimental/CMakeLists.txt | 3 +++ examples/{ => experimental}/gpuRefactorSame.cpp | 11 ++--------- 3 files changed, 5 insertions(+), 14 deletions(-) rename examples/{ => experimental}/gpuRefactorSame.cpp (96%) diff --git a/examples/CMakeLists.txt b/examples/CMakeLists.txt index d5eee7a68..e97ec0531 100644 --- a/examples/CMakeLists.txt +++ b/examples/CMakeLists.txt @@ -30,11 +30,6 @@ if(RESOLVE_USE_KLU) target_link_libraries(gpuRefactor.exe PRIVATE ReSolve) endif(RESOLVE_USE_GPU) - if(RESOLVE_USE_HIP) #HIP specific example - add_executable(gpuRefactorSame.exe gpuRefactorSame.cpp) - target_link_libraries(gpuRefactorSame.exe PRIVATE ReSolve) - endif(RESOLVE_USE_HIP) - # Create KLU+CUDA examples if(RESOLVE_USE_CUDA) diff --git a/examples/experimental/CMakeLists.txt b/examples/experimental/CMakeLists.txt index 5d22d2fce..2e6c4cb28 100644 --- a/examples/experimental/CMakeLists.txt +++ b/examples/experimental/CMakeLists.txt @@ -32,6 +32,9 @@ if(RESOLVE_USE_KLU) add_executable(klu_rocsolverrf_check_redo.exe r_KLU_rocsolverrf_redo_factorization.cpp) target_link_libraries(klu_rocsolverrf_check_redo.exe PRIVATE ReSolve) + add_executable(gpuRefactorSame.exe gpuRefactorSame.cpp) # Example with refactorization on GPU with the same system (from kalmarek) + target_link_libraries(gpuRefactorSame.exe PRIVATE ReSolve) + endif(RESOLVE_USE_HIP) endif(RESOLVE_USE_KLU) diff --git a/examples/gpuRefactorSame.cpp b/examples/experimental/gpuRefactorSame.cpp similarity index 96% rename from examples/gpuRefactorSame.cpp rename to examples/experimental/gpuRefactorSame.cpp index c8c2ef271..18c2e4893 100644 --- a/examples/gpuRefactorSame.cpp +++ b/examples/experimental/gpuRefactorSame.cpp @@ -14,7 +14,7 @@ #include #include -#include "ExampleHelper.hpp" +#include "examples/ExampleHelper.hpp" int main() { @@ -179,12 +179,6 @@ int main() L->print(std::cout, 0); std::cout << "U factor:\n"; U->print(std::cout, 0); - if (F != nullptr) { - std::cout << "F factor:\n"; - F->print(std::cout, 0); - } else { - std::cout << "F factor is not available.\n"; - } ReSolve::index_type* P = KLU.getPOrdering(); ReSolve::index_type* Q = KLU.getQOrdering(); @@ -194,8 +188,7 @@ int main() std::cout << rhs_before[i] << (i < n - 1 ? ", " : ""); } std::cout << "]" << std::endl; - //Rf.setupCsr(A, L, U, P, Q, vec_rhs); - Rf.setup(A, L, U, P, Q, vec_rhs); + Rf.setupCsr(A, L, U, P, Q, vec_rhs); std::cout << "RocSolverRf setup completed\n"; // Test refactorization with the same matrix (in practice, you'd change matrix values) From afdecb2750967a7a966f9ec70f2a6bef516257fb Mon Sep 17 00:00:00 2001 From: shakedregev Date: Tue, 2 Sep 2025 11:10:39 -0400 Subject: [PATCH 4/8] removed commented code --- examples/experimental/gpuRefactorSame.cpp | 4 ---- 1 file changed, 4 deletions(-) diff --git a/examples/experimental/gpuRefactorSame.cpp b/examples/experimental/gpuRefactorSame.cpp index 18c2e4893..d84a52cba 100644 --- a/examples/experimental/gpuRefactorSame.cpp +++ b/examples/experimental/gpuRefactorSame.cpp @@ -164,10 +164,6 @@ int main() helper.resetSystem(A, vec_rhs, vec_x); helper.printShortSummary(); - // Extract factors and configure refactorization solver - // ReSolve::matrix::Csc* L = (ReSolve::matrix::Csc*) KLU.getLFactor(); - // ReSolve::matrix::Csc* U = (ReSolve::matrix::Csc*) KLU.getUFactor(); - // ReSolve::matrix::Csc* F = (ReSolve::matrix::Csc*) KLU.getFFactor(); ReSolve::matrix::Csr* L = (ReSolve::matrix::Csr*) KLU.getLFactorCsr(); ReSolve::matrix::Csr* U = (ReSolve::matrix::Csr*) KLU.getUFactorCsr(); if (L == nullptr || U == nullptr) { From f437d4c64606b2397724cb61fe76ffecf9892e5c Mon Sep 17 00:00:00 2001 From: shakedregev Date: Tue, 2 Sep 2025 15:27:04 +0000 Subject: [PATCH 5/8] Apply pre-commmit fixes --- examples/experimental/gpuRefactorSame.cpp | 437 +++++++++++----------- 1 file changed, 225 insertions(+), 212 deletions(-) diff --git a/examples/experimental/gpuRefactorSame.cpp b/examples/experimental/gpuRefactorSame.cpp index d84a52cba..86c146df6 100644 --- a/examples/experimental/gpuRefactorSame.cpp +++ b/examples/experimental/gpuRefactorSame.cpp @@ -1,11 +1,11 @@ -#include -#include #include -#include #include +#include #include +#include +#include -#include +#include "examples/ExampleHelper.hpp" #include #include #include @@ -13,225 +13,238 @@ #include #include #include - -#include "examples/ExampleHelper.hpp" +#include int main() { - std::cout << "Simple ReSolve Example (following gpuRefactor pattern)\n"; - std::cout << "=====================================================\n"; - - // Matrix data in CSC format (converted from 1-based to 0-based indexing) - // Julia CSC data: colptr, rowval, nzval - std::vector csc_col_ptr = {0, 1, 2, 3, 6, 8, 10}; // Julia colptr - 1 - std::vector csc_row_ind = {0, 1, 2, 0, 1, 3, 0, 4, 1, 5}; // Julia rowval - 1 - std::vector csc_values = { - 0.00141521219668486, - 0.6984313028075056, - -10.0, - -0.8307014463420237, - 0.1846074461471301, - -10.0, - -0.8408959803799315, - -10.0, - -0.835721417925707, - -10.0 - }; - using namespace ReSolve::examples; - // Convert CSC to CSR manually for this example - int n = 6; // 6x6 matrix - int nnz = csc_values.size(); // 10 non-zeros - - std::cout << "Matrix size: " << n << "x" << n << "\n"; - std::cout << "Number of non-zeros: " << nnz << "\n"; - - // Convert CSC to CSR - std::vector csr_row_ptr(n + 1, 0); - std::vector csr_col_ind(nnz); - std::vector csr_values(nnz); - - // Count entries per row - for (int i = 0; i < nnz; ++i) { - csr_row_ptr[csc_row_ind[i] + 1]++; - } - - // Convert counts to pointers - for (int i = 1; i <= n; ++i) { - csr_row_ptr[i] += csr_row_ptr[i - 1]; - } - - // Fill CSR arrays - std::vector temp_row_ptr = csr_row_ptr; - for (int col = 0; col < n; ++col) { - for (int idx = csc_col_ptr[col]; idx < csc_col_ptr[col + 1]; ++idx) { - int row = csc_row_ind[idx]; - int pos = temp_row_ptr[row]++; - csr_col_ind[pos] = col; - csr_values[pos] = csc_values[idx]; - } - } - - std::cout << "Matrix data converted to CSR format\n"; - - // Create workspace - ReSolve::LinAlgWorkspaceHIP workspace; - workspace.initializeHandles(); - ExampleHelper helper(workspace); - - - // Create matrix and vector handlers - ReSolve::MatrixHandler matrix_handler(&workspace); - ReSolve::VectorHandler vector_handler(&workspace); - - // Direct solvers instantiation - ReSolve::LinSolverDirectKLU KLU; - ReSolve::LinSolverDirectRocSolverRf Rf(&workspace); - - // Create CSR matrix - ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(n, n, nnz, false, false); - A->allocateMatrixData(ReSolve::memory::HOST); - A->copyDataFrom(csr_row_ptr.data(), csr_col_ind.data(), csr_values.data(), - ReSolve::memory::HOST, ReSolve::memory::HOST); - A->allocateMatrixData(ReSolve::memory::DEVICE); - A->syncData(ReSolve::memory::DEVICE); - - std::cout << "CSR matrix created and synced to device\n"; - - // A->print(std::cout, 0); - - // Create vectors - ReSolve::vector::Vector* vec_rhs = new ReSolve::vector::Vector(n); - ReSolve::vector::Vector* vec_x = new ReSolve::vector::Vector(n); - - vec_rhs->allocate(ReSolve::memory::HOST); - vec_rhs->allocate(ReSolve::memory::DEVICE); - vec_x->allocate(ReSolve::memory::HOST); - vec_x->allocate(ReSolve::memory::DEVICE); - - // Set right-hand side vector - std::vector rhs_data = { - 0.04204480190421101, - 0.1868729984154292, - -0.1581138884334949, - 0.04939382906591484, - -1.0, - 0.1118034015563139 - }; - vec_rhs->copyDataFrom(rhs_data.data(), ReSolve::memory::HOST, ReSolve::memory::HOST); - vec_rhs->syncData(ReSolve::memory::DEVICE); - - std::cout << "Right-hand side vector set\n"; - - // Setup KLU solver - KLU.setup(A); - matrix_handler.setValuesChanged(true, ReSolve::memory::DEVICE); - - // Analysis (symbolic factorization) - int status = KLU.analyze(); - std::cout << "KLU analysis status: " << status << std::endl; - if (status != 0) { - std::cerr << "KLU analysis failed!" << std::endl; - return 1; - } - - // Numeric factorization - status = KLU.factorize(); - std::cout << "KLU factorization status: " << status << std::endl; - if (status != 0) { - std::cerr << "KLU factorization failed!" << std::endl; - return 1; - } - - // Triangular solve - status = KLU.solve(vec_rhs, vec_x); - std::cout << "KLU solve status: " << status << std::endl; - if (status != 0) { - std::cerr << "KLU solve failed!" << std::endl; - return 1; + std::cout << "Simple ReSolve Example (following gpuRefactor pattern)\n"; + std::cout << "=====================================================\n"; + + // Matrix data in CSC format (converted from 1-based to 0-based indexing) + // Julia CSC data: colptr, rowval, nzval + std::vector csc_col_ptr = {0, 1, 2, 3, 6, 8, 10}; // Julia colptr - 1 + std::vector csc_row_ind = {0, 1, 2, 0, 1, 3, 0, 4, 1, 5}; // Julia rowval - 1 + std::vector csc_values = { + 0.00141521219668486, + 0.6984313028075056, + -10.0, + -0.8307014463420237, + 0.1846074461471301, + -10.0, + -0.8408959803799315, + -10.0, + -0.835721417925707, + -10.0}; + using namespace ReSolve::examples; + // Convert CSC to CSR manually for this example + int n = 6; // 6x6 matrix + int nnz = csc_values.size(); // 10 non-zeros + + std::cout << "Matrix size: " << n << "x" << n << "\n"; + std::cout << "Number of non-zeros: " << nnz << "\n"; + + // Convert CSC to CSR + std::vector csr_row_ptr(n + 1, 0); + std::vector csr_col_ind(nnz); + std::vector csr_values(nnz); + + // Count entries per row + for (int i = 0; i < nnz; ++i) + { + csr_row_ptr[csc_row_ind[i] + 1]++; + } + + // Convert counts to pointers + for (int i = 1; i <= n; ++i) + { + csr_row_ptr[i] += csr_row_ptr[i - 1]; + } + + // Fill CSR arrays + std::vector temp_row_ptr = csr_row_ptr; + for (int col = 0; col < n; ++col) + { + for (int idx = csc_col_ptr[col]; idx < csc_col_ptr[col + 1]; ++idx) + { + int row = csc_row_ind[idx]; + int pos = temp_row_ptr[row]++; + csr_col_ind[pos] = col; + csr_values[pos] = csc_values[idx]; } - - std::cout << "System solved successfully with KLU\n"; - - // Get solution - double* solution = vec_x->getData(ReSolve::memory::HOST); - - std::cout << "Solution vector:\n"; + } + + std::cout << "Matrix data converted to CSR format\n"; + + // Create workspace + ReSolve::LinAlgWorkspaceHIP workspace; + workspace.initializeHandles(); + ExampleHelper helper(workspace); + + // Create matrix and vector handlers + ReSolve::MatrixHandler matrix_handler(&workspace); + ReSolve::VectorHandler vector_handler(&workspace); + + // Direct solvers instantiation + ReSolve::LinSolverDirectKLU KLU; + ReSolve::LinSolverDirectRocSolverRf Rf(&workspace); + + // Create CSR matrix + ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(n, n, nnz, false, false); + A->allocateMatrixData(ReSolve::memory::HOST); + A->copyDataFrom(csr_row_ptr.data(), csr_col_ind.data(), csr_values.data(), ReSolve::memory::HOST, ReSolve::memory::HOST); + A->allocateMatrixData(ReSolve::memory::DEVICE); + A->syncData(ReSolve::memory::DEVICE); + + std::cout << "CSR matrix created and synced to device\n"; + + // A->print(std::cout, 0); + + // Create vectors + ReSolve::vector::Vector* vec_rhs = new ReSolve::vector::Vector(n); + ReSolve::vector::Vector* vec_x = new ReSolve::vector::Vector(n); + + vec_rhs->allocate(ReSolve::memory::HOST); + vec_rhs->allocate(ReSolve::memory::DEVICE); + vec_x->allocate(ReSolve::memory::HOST); + vec_x->allocate(ReSolve::memory::DEVICE); + + // Set right-hand side vector + std::vector rhs_data = { + 0.04204480190421101, + 0.1868729984154292, + -0.1581138884334949, + 0.04939382906591484, + -1.0, + 0.1118034015563139}; + vec_rhs->copyDataFrom(rhs_data.data(), ReSolve::memory::HOST, ReSolve::memory::HOST); + vec_rhs->syncData(ReSolve::memory::DEVICE); + + std::cout << "Right-hand side vector set\n"; + + // Setup KLU solver + KLU.setup(A); + matrix_handler.setValuesChanged(true, ReSolve::memory::DEVICE); + + // Analysis (symbolic factorization) + int status = KLU.analyze(); + std::cout << "KLU analysis status: " << status << std::endl; + if (status != 0) + { + std::cerr << "KLU analysis failed!" << std::endl; + return 1; + } + + // Numeric factorization + status = KLU.factorize(); + std::cout << "KLU factorization status: " << status << std::endl; + if (status != 0) + { + std::cerr << "KLU factorization failed!" << std::endl; + return 1; + } + + // Triangular solve + status = KLU.solve(vec_rhs, vec_x); + std::cout << "KLU solve status: " << status << std::endl; + if (status != 0) + { + std::cerr << "KLU solve failed!" << std::endl; + return 1; + } + + std::cout << "System solved successfully with KLU\n"; + + // Get solution + double* solution = vec_x->getData(ReSolve::memory::HOST); + + std::cout << "Solution vector:\n"; + std::cout << "["; + for (int i = 0; i < n; ++i) + { + std::cout << solution[i] << (i < n - 1 ? ", " : ""); + } + std::cout << "]" << std::endl; + + helper.resetSystem(A, vec_rhs, vec_x); + helper.printShortSummary(); + ReSolve::matrix::Csr* L = (ReSolve::matrix::Csr*) KLU.getLFactorCsr(); + ReSolve::matrix::Csr* U = (ReSolve::matrix::Csr*) KLU.getUFactorCsr(); + if (L == nullptr || U == nullptr) + { + std::cout << "Factor extraction from KLU failed!\n"; + } + else + { + std::cout << "L and U factors extracted successfully\n"; + // Print L and U factors + std::cout << "L factor:\n"; + L->print(std::cout, 0); + std::cout << "U factor:\n"; + U->print(std::cout, 0); + + ReSolve::index_type* P = KLU.getPOrdering(); + ReSolve::index_type* Q = KLU.getQOrdering(); + double* rhs_before = vec_rhs->getData(ReSolve::memory::HOST); std::cout << "["; - for (int i = 0; i < n; ++i) { - std::cout << solution[i] << (i < n - 1 ? ", " : ""); + for (int i = 0; i < n; ++i) + { + std::cout << rhs_before[i] << (i < n - 1 ? ", " : ""); } std::cout << "]" << std::endl; - - helper.resetSystem(A, vec_rhs, vec_x); - helper.printShortSummary(); - ReSolve::matrix::Csr* L = (ReSolve::matrix::Csr*) KLU.getLFactorCsr(); - ReSolve::matrix::Csr* U = (ReSolve::matrix::Csr*) KLU.getUFactorCsr(); - if (L == nullptr || U == nullptr) { - std::cout << "Factor extraction from KLU failed!\n"; - } else { - std::cout << "L and U factors extracted successfully\n"; - // Print L and U factors - std::cout << "L factor:\n"; - L->print(std::cout, 0); - std::cout << "U factor:\n"; - U->print(std::cout, 0); - - ReSolve::index_type* P = KLU.getPOrdering(); - ReSolve::index_type* Q = KLU.getQOrdering(); - double* rhs_before = vec_rhs->getData(ReSolve::memory::HOST); + Rf.setupCsr(A, L, U, P, Q, vec_rhs); + std::cout << "RocSolverRf setup completed\n"; + + // Test refactorization with the same matrix (in practice, you'd change matrix values) + std::cout << "\nTesting refactorization with same matrix...\n"; + + // Refactorize + status = Rf.refactorize(); + std::cout << "RocSolverRf refactorization status: " << status << std::endl; + if (status != 0) + { + std::cerr << "Refactorization failed!" << std::endl; + } + else + { + // Solve with refactorization + status = Rf.solve(vec_rhs, vec_x); + std::cout << "RocSolverRf solve status: " << status << std::endl; + helper.resetSystem(A, vec_rhs, vec_x); + + if (status == 0) + { + std::cout << "System solved successfully with refactorization\n"; + + // Get solution + vec_x->syncData(ReSolve::memory::HOST); + solution = vec_x->getData(ReSolve::memory::HOST); + + std::cout << "Refactorization solution vector:\n"; std::cout << "["; - for (int i = 0; i < n; ++i) { - std::cout << rhs_before[i] << (i < n - 1 ? ", " : ""); + for (int i = 0; i < n; ++i) + { + std::cout << solution[i] << (i < n - 1 ? ", " : ""); } std::cout << "]" << std::endl; - Rf.setupCsr(A, L, U, P, Q, vec_rhs); - std::cout << "RocSolverRf setup completed\n"; - - // Test refactorization with the same matrix (in practice, you'd change matrix values) - std::cout << "\nTesting refactorization with same matrix...\n"; - - // Refactorize - status = Rf.refactorize(); - std::cout << "RocSolverRf refactorization status: " << status << std::endl; - if (status != 0) { - std::cerr << "Refactorization failed!" << std::endl; - } else { - // Solve with refactorization - status = Rf.solve(vec_rhs, vec_x); - std::cout << "RocSolverRf solve status: " << status << std::endl; - helper.resetSystem(A, vec_rhs, vec_x); - - if (status == 0) { - std::cout << "System solved successfully with refactorization\n"; - - // Get solution - vec_x->syncData(ReSolve::memory::HOST); - solution = vec_x->getData(ReSolve::memory::HOST); - - std::cout << "Refactorization solution vector:\n"; - std::cout << "["; - for (int i = 0; i < n; ++i) { - std::cout << solution[i] << (i < n - 1 ? ", " : ""); - } - std::cout << "]" << std::endl; - - // Get the rhs - double* rhs_solution = vec_rhs->getData(ReSolve::memory::HOST); - std::cout << "Right-hand side vector after solve:\n"; - std::cout << "["; - for (int i = 0; i < n; ++i) { - std::cout << rhs_solution[i] << (i < n - 1 ? ", " : ""); - } - std::cout << "]" << std::endl; - } + + // Get the rhs + double* rhs_solution = vec_rhs->getData(ReSolve::memory::HOST); + std::cout << "Right-hand side vector after solve:\n"; + std::cout << "["; + for (int i = 0; i < n; ++i) + { + std::cout << rhs_solution[i] << (i < n - 1 ? ", " : ""); } + std::cout << "]" << std::endl; + } } + } - // Cleanup - delete A; - delete vec_x; - delete vec_rhs; + // Cleanup + delete A; + delete vec_x; + delete vec_rhs; - std::cout << "Example completed successfully!\n"; - return 0; -} \ No newline at end of file + std::cout << "Example completed successfully!\n"; + return 0; +} From bfc1e27ed402464fb6dcd5690b54b9af416e89f7 Mon Sep 17 00:00:00 2001 From: Shaked Regev <35384901+shakedregev@users.noreply.github.com> Date: Tue, 2 Sep 2025 11:28:40 -0400 Subject: [PATCH 6/8] Update examples/experimental/CMakeLists.txt --- examples/experimental/CMakeLists.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/experimental/CMakeLists.txt b/examples/experimental/CMakeLists.txt index 2e6c4cb28..055b86252 100644 --- a/examples/experimental/CMakeLists.txt +++ b/examples/experimental/CMakeLists.txt @@ -32,7 +32,7 @@ if(RESOLVE_USE_KLU) add_executable(klu_rocsolverrf_check_redo.exe r_KLU_rocsolverrf_redo_factorization.cpp) target_link_libraries(klu_rocsolverrf_check_redo.exe PRIVATE ReSolve) - add_executable(gpuRefactorSame.exe gpuRefactorSame.cpp) # Example with refactorization on GPU with the same system (from kalmarek) + add_executable(gpuRefactorSame.exe gpuRefactorSame.cpp) # Example with refactorization on GPU with the same system (from kalmarek). target_link_libraries(gpuRefactorSame.exe PRIVATE ReSolve) endif(RESOLVE_USE_HIP) From 7cb5d9609104972f8647287cb60d0c998cbdde0b Mon Sep 17 00:00:00 2001 From: shakedregev Date: Thu, 4 Sep 2025 12:15:41 -0400 Subject: [PATCH 7/8] addressed change requests --- .../{gpuRefactorSame.cpp => r_KLU_rocsolverrf_asym6x6.cpp} | 5 +++++ 1 file changed, 5 insertions(+) rename examples/experimental/{gpuRefactorSame.cpp => r_KLU_rocsolverrf_asym6x6.cpp} (97%) diff --git a/examples/experimental/gpuRefactorSame.cpp b/examples/experimental/r_KLU_rocsolverrf_asym6x6.cpp similarity index 97% rename from examples/experimental/gpuRefactorSame.cpp rename to examples/experimental/r_KLU_rocsolverrf_asym6x6.cpp index 86c146df6..1360f52b7 100644 --- a/examples/experimental/gpuRefactorSame.cpp +++ b/examples/experimental/r_KLU_rocsolverrf_asym6x6.cpp @@ -1,3 +1,8 @@ +/* +This is a simple example demonstrating the use of ReSolve's direct solvers KLU with RocSolverRf on the same system. +Previously, this example was flagged as failing in issue #332. +*/ + #include #include #include From a35406ef3c0ce8b0c6491e72a8e3da85fdb4c5f5 Mon Sep 17 00:00:00 2001 From: shakedregev Date: Thu, 4 Sep 2025 13:13:47 -0400 Subject: [PATCH 8/8] changed cmake --- examples/experimental/CMakeLists.txt | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/examples/experimental/CMakeLists.txt b/examples/experimental/CMakeLists.txt index 055b86252..8b8518521 100644 --- a/examples/experimental/CMakeLists.txt +++ b/examples/experimental/CMakeLists.txt @@ -32,8 +32,8 @@ if(RESOLVE_USE_KLU) add_executable(klu_rocsolverrf_check_redo.exe r_KLU_rocsolverrf_redo_factorization.cpp) target_link_libraries(klu_rocsolverrf_check_redo.exe PRIVATE ReSolve) - add_executable(gpuRefactorSame.exe gpuRefactorSame.cpp) # Example with refactorization on GPU with the same system (from kalmarek). - target_link_libraries(gpuRefactorSame.exe PRIVATE ReSolve) + add_executable(r_KLU_rocsolverrf_asym6x6.exe r_KLU_rocsolverrf_asym6x6.cpp) # Example with refactorization on GPU with the same system (from kalmarek). + target_link_libraries(r_KLU_rocsolverrf_asym6x6.exe PRIVATE ReSolve) endif(RESOLVE_USE_HIP)