Skip to content

Add sysGmres example #401

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
pelesh merged 7 commits into from
Dec 18, 2025
Merged

Add sysGmres example #401

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
5377927
Add sysGmres to CMakeLists.txt
KakeruUeda Dec 14, 2025
7c08771
Add getter for memspace
KakeruUeda Dec 15, 2025
8cd7ba0
Add sysGmres example
KakeruUeda Dec 15, 2025
f04f381
Fix comments
KakeruUeda Dec 16, 2025
4caae5a
Apply pre-commmit fixes
KakeruUeda Dec 17, 2025
500b3e5
Update CHANGELOG.md
KakeruUeda Dec 18, 2025
d204ada
Apply pre-commmit fixes
KakeruUeda Dec 18, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions CHANGELOG.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -61,3 +61,5 @@ It is seamless from the user perspective and fixed many bugs.
13. Added LUSOL direct solver, which can factorize matrices and extract factors independently of KLU.

14. Various Spack updates.

15. Added examples/sysGmres.cpp to demonstrate how to use SystemSolver with GMRES.
6 changes: 6 additions & 0 deletions examples/CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,10 @@
add_executable(rand_gmres.exe randGmres.cpp)
target_link_libraries(rand_gmres.exe PRIVATE ReSolve)

# Build an example with a system GMRES solver
add_executable(sysGmres.exe sysGmres.cpp)
target_link_libraries(sysGmres.exe PRIVATE ReSolve)

if(RESOLVE_USE_KLU)

# Build example with KLU factorization on CPU
Expand Down Expand Up @@ -46,6 +50,8 @@ set(installable_executables "")
# Install all examples in bin directory
list(APPEND installable_executables rand_gmres.exe)

list(APPEND installable_executables sysGmres.exe)

if(RESOLVE_USE_KLU)
list(APPEND installable_executables kluFactor.exe kluRefactor.exe
sysRefactor.exe
Expand Down
7 changes: 7 additions & 0 deletions examples/ExampleHelper.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -99,11 +99,18 @@ namespace ReSolve
}
}

/// Returns the configured hardware backend
std::string getHardwareBackend() const
{
return hardware_backend_;
}

/// Returns the configured memory space
ReSolve::memory::MemorySpace getMemspace() const
{
return memspace_;
}

/**
* @brief Set the new linear system together with its computed solution
* and compute solution error and residual norms.
Expand Down
319 changes: 319 additions & 0 deletions examples/sysGmres.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,319 @@
/**
* @file sysGmres.cpp
* @author Kakeru Ueda (ueda.k.2290@m.isct.ac.jp)
* @brief Example of solving a linear system with GMRES using SystemSolver.
*/

#include <iomanip>
#include <iostream>
#include <string>

#include "ExampleHelper.hpp"
#include <resolve/SystemSolver.hpp>
#include <resolve/matrix/Csr.hpp>
#include <resolve/matrix/io.hpp>
#include <resolve/utilities/params/CliOptions.hpp>
#include <resolve/vector/Vector.hpp>
#include <resolve/workspace/LinAlgWorkspace.hpp>

// Uses ReSolve data types
using namespace ReSolve;
using namespace ReSolve::examples;
using namespace ReSolve::memory;
using vector_type = ReSolve::vector::Vector;

/// Prints help message describing system usage
void printHelpInfo()
{
std::cout << "\nsysGmres.exe loads a linear system from files and solves it using GMRES.\n\n";
std::cout << "Usage:\n\t./";
std::cout << "sysGmres.exe -m <matrix file> -r <rhs file>\n\n";
std::cout << "Optional features:\n";
std::cout << "\t-b <cpu|cuda|hip> \tSelects hardware backend.\n";
std::cout << "\t-h \tPrints this message.\n";
std::cout << "\t-i <iter method> \tIterative method: randgmres or fgmres (default 'randgmres').\n";
std::cout << "\t-g <gs method> \tGram-Schmidt method: cgs1, cgs2, or mgs (default 'cgs2').\n";
std::cout << "\t-s <sketching method> \tSketching method: count or fwht (default 'count')\n";
std::cout << "\t-x <flexible> \tEnable flexible: yes or no (default 'yes')\n\n";
}

//
// Forward declarations of functions
//

/**
* @brief Example of solving a linear system with GMRES using SystemSolver.
*
* @tparam workspace_type - Type of the workspace to use
* @param[in] argc - Number of command line arguments
* @param[in] argv - Command line arguments
* @return 0 if the example ran successfully, 1 otherwise
*/
template <class workspace_type>
static int sysGmres(int argc, char* argv[]);

/// Checks if inputs for GMRES are valid, otherwise sets defaults
static void processInputs(std::string& method,
std::string& gs,
std::string& sketch,
std::string& flexible);

/// Main function selects example to be run
int main(int argc, char* argv[])
{
CliOptions options(argc, argv);

bool is_help = options.hasKey("-h");
if (is_help)
{
printHelpInfo();
return 0;
}

auto opt = options.getParamFromKey("-b");
if (!opt)
{
std::cout << "No backend option provided. Defaulting to CPU.\n";
return sysGmres<ReSolve::LinAlgWorkspaceCpu>(argc, argv);
}
#ifdef RESOLVE_USE_CUDA
else if (opt->second == "cuda")
{
return sysGmres<ReSolve::LinAlgWorkspaceCUDA>(argc, argv);
}
#endif
#ifdef RESOLVE_USE_HIP
else if (opt->second == "hip")
{
return sysGmres<ReSolve::LinAlgWorkspaceHIP>(argc, argv);
}
#endif
else if (opt->second == "cpu")
{
return sysGmres<ReSolve::LinAlgWorkspaceCpu>(argc, argv);
}
else
{
std::cout << "Re::Solve is not build with support for " << opt->second;
std::cout << " backend.\n";
return 1;
}

return 0;
}

//
// Definitions of functions
//

template <class workspace_type>
int sysGmres(int argc, char* argv[])
{
// return_code is used as a failure flag.
int return_code = 0;
int status = 0;

// Collect all CLI
CliOptions options(argc, argv);

bool is_help = options.hasKey("-h");
if (is_help)
{
printHelpInfo();
return 0;
}

// Read matrix file
auto opt = options.getParamFromKey("-m");
std::string matrix_pathname("");
if (opt)
{
matrix_pathname = opt->second;
}
else
{
std::cout << "Incorrect input!\n";
printHelpInfo();
return 1;
}

// Read right-hand-side vector file
std::string rhs_pathname("");
opt = options.getParamFromKey("-r");
if (opt)
{
rhs_pathname = opt->second;
}
else
{
std::cout << "Incorrect input!\n";
printHelpInfo();
return 1;
}

// Read GMRES-related options
opt = options.getParamFromKey("-i");
std::string method = opt ? (*opt).second : "randgmres";

opt = options.getParamFromKey("-g");
std::string gs = opt ? (*opt).second : "cgs2";

opt = options.getParamFromKey("-s");
std::string sketch = opt ? (*opt).second : "count";

opt = options.getParamFromKey("-x");
std::string flexible = opt ? (*opt).second : "yes";

processInputs(method, gs, sketch, flexible);

std::cout << "Matrix file: " << matrix_pathname << "\n"
<< "RHS file: " << rhs_pathname << "\n";

// Create workspace
workspace_type workspace;
workspace.initializeHandles();

// Create a helper object (computing errors, printing summaries, etc.)
ExampleHelper<workspace_type> helper(workspace);
std::string hw_backend = helper.getHardwareBackend();
std::cout << "sysGmres with " << hw_backend << " backend\n";

// Set memory space
MemorySpace memspace = helper.getMemspace();

// Set solver
SystemSolver solver(&workspace,
"none",
"none",
method,
"ilu0",
"none");

solver.setGramSchmidtMethod(gs);

// Read and open matrix and right-hand-side vector
std::ifstream mat_file(matrix_pathname);
if (!mat_file.is_open())
{
std::cout << "Failed to open matrix file: " << matrix_pathname << "\n";
return 1;
}
std::ifstream rhs_file(rhs_pathname);
if (!rhs_file.is_open())
{
std::cout << "Failed to open RHS file: " << rhs_pathname << "\n";
return 1;
}

bool is_expand_symmetric = true;

// Load system matrix and RHS vector from input files
matrix::Csr* A = io::createCsrFromFile(mat_file, is_expand_symmetric);
vector_type* vec_rhs = io::createVectorFromFile(rhs_file);

// Create solution vector
vector_type* vec_x = new vector_type(A->getNumRows());
vec_x->allocate(memspace);

if (memspace == memory::DEVICE)
{
// Copy data to the device
A->syncData(memspace);
vec_rhs->syncData(memspace);
}

mat_file.close();
rhs_file.close();

// Set iterative solver options
solver.getIterativeSolver().setCliParam("maxit", "2500");
solver.getIterativeSolver().setCliParam("tol", "1e-12");

status = solver.setMatrix(A);
std::cout << "solver.setMatrix returned status: " << status << "\n";
if (status != 0)
{
return_code = 1;
}

// Set GMRES solver options
if (method == "randgmres")
{
solver.setSketchingMethod(sketch);
}
solver.getIterativeSolver().setCliParam("flexible", flexible);
solver.getIterativeSolver().setCliParam("restart", "200");

// Set up the preconditioner
if (return_code == 0)
{
status = solver.preconditionerSetup();
std::cout << "solver.preconditionerSetup returned status: " << status << "\n";
if (status != 0)
{
return_code = 1;
}
}

// Solve the system
if (return_code == 0)
{
status = solver.solve(vec_rhs, vec_x);
std::cout << "solver.solve returned status: " << status << "\n";
if (status != 0)
{
return_code = 1;
}
}

if (return_code == 0)
{
helper.resetSystem(A, vec_rhs, vec_x);

// Get reference to iterative solver and print results
LinSolverIterative& iter_solver = solver.getIterativeSolver();
helper.printIterativeSolverSummary(&iter_solver);
}

// Free matrix/vector objects.
delete A;
delete vec_rhs;
delete vec_x;

return return_code;
}

void processInputs(std::string& method, std::string& gs, std::string& sketch, std::string& flexible)
{
if (method == "randgmres")
{
if ((sketch != "count") && (sketch != "fwht"))
{
std::cout << "Sketching method " << sketch << " not recognized.\n";
std::cout << "Setting sketch to the default (count).\n\n";
sketch = "count";
}
}

if ((method != "randgmres") && (method != "fgmres"))
{
std::cout << "Iterative method " << method << " not recognized.\n";
std::cout << "Setting iterative solver method to the default (RANDGMRES).\n\n";
method = "randgmres";
}

if (gs != "cgs1" && gs != "cgs2" && gs != "mgs" && gs != "mgs_two_sync"
&& gs != "mgs_pm")
{
std::cout << "Orthogonalization method " << gs << " not recognized.\n";
std::cout << "Setting orthogonalization to the default (CGS2).\n\n";
gs = "cgs2";
}

if ((flexible != "yes") && (flexible != "no"))
{
std::cout << "Flexible option " << flexible << " not recognized.\n";
std::cout << "Setting flexible to the default (yes).\n\n";
flexible = "yes";
}
}