From fe6e715c4c4eba2af06674aebf25eb81654b8a46 Mon Sep 17 00:00:00 2001 From: Tobias Ribizel Date: Wed, 15 Nov 2023 13:19:43 -0700 Subject: [PATCH] add example --- examples/CMakeLists.txt | 1 + .../CMakeLists.txt | 14 +++ .../reordered-preconditioned-solver/build.sh | 17 +++ .../data/A.mtx | 114 ++++++++++++++++++ .../data/b.mtx | 21 ++++ .../data/x0.mtx | 21 ++++ .../doc/builds-on | 1 + .../doc/intro.dox | 5 + .../reordered-preconditioned-solver/doc/kind | 1 + .../doc/results.dox | 35 ++++++ .../doc/short-intro | 1 + .../doc/tooltip | 1 + .../reordered-preconditioned-solver.cpp | 111 +++++++++++++++++ 13 files changed, 343 insertions(+) create mode 100644 examples/reordered-preconditioned-solver/CMakeLists.txt create mode 100755 examples/reordered-preconditioned-solver/build.sh create mode 100644 examples/reordered-preconditioned-solver/data/A.mtx create mode 100644 examples/reordered-preconditioned-solver/data/b.mtx create mode 100644 examples/reordered-preconditioned-solver/data/x0.mtx create mode 100644 examples/reordered-preconditioned-solver/doc/builds-on create mode 100644 examples/reordered-preconditioned-solver/doc/intro.dox create mode 100644 examples/reordered-preconditioned-solver/doc/kind create mode 100644 examples/reordered-preconditioned-solver/doc/results.dox create mode 100644 examples/reordered-preconditioned-solver/doc/short-intro create mode 100644 examples/reordered-preconditioned-solver/doc/tooltip create mode 100644 examples/reordered-preconditioned-solver/reordered-preconditioned-solver.cpp diff --git a/examples/CMakeLists.txt b/examples/CMakeLists.txt index 33e3bab735a..5c98f3b8332 100644 --- a/examples/CMakeLists.txt +++ b/examples/CMakeLists.txt @@ -30,6 +30,7 @@ set(EXAMPLES_LIST par-ilu-convergence performance-debugging preconditioner-export + reordered-preconditioned-solver simple-solver-logging) if(GINKGO_BUILD_CUDA AND GINKGO_BUILD_OMP) diff --git a/examples/reordered-preconditioned-solver/CMakeLists.txt b/examples/reordered-preconditioned-solver/CMakeLists.txt new file mode 100644 index 00000000000..17c17f2b71c --- /dev/null +++ b/examples/reordered-preconditioned-solver/CMakeLists.txt @@ -0,0 +1,14 @@ +cmake_minimum_required(VERSION 3.16) +project(reordered-preconditioned-solver) + +# We only need to find Ginkgo if we build this example stand-alone +if (NOT GINKGO_BUILD_EXAMPLES) + find_package(Ginkgo 1.8.0 REQUIRED) +endif() +add_executable(reordered-preconditioned-solver reordered-preconditioned-solver.cpp) +target_link_libraries(reordered-preconditioned-solver Ginkgo::ginkgo) + +# Copy the data files to the execution directory +configure_file(data/A.mtx data/A.mtx COPYONLY) +configure_file(data/b.mtx data/b.mtx COPYONLY) +configure_file(data/x0.mtx data/x0.mtx COPYONLY) diff --git a/examples/reordered-preconditioned-solver/build.sh b/examples/reordered-preconditioned-solver/build.sh new file mode 100755 index 00000000000..cd35af13dc9 --- /dev/null +++ b/examples/reordered-preconditioned-solver/build.sh @@ -0,0 +1,17 @@ +#!/bin/bash + +# set up script +if [ $# -ne 1 ]; then + echo -e "Usage: $0 GINKGO_BUILD_DIRECTORY" + exit 1 +fi +BUILD_DIR=$1 +THIS_DIR=$( cd "$( dirname "${BASH_SOURCE[0]}" )" &>/dev/null && pwd ) + +source ${THIS_DIR}/../build-setup.sh + +# build +${CXX} -std=c++14 -o ${THIS_DIR}/preconditioned-solver \ + ${THIS_DIR}/reordered-preconditioned-solver.cpp \ + -I${THIS_DIR}/../../include -I${BUILD_DIR}/include \ + -L${THIS_DIR} ${LINK_FLAGS} diff --git a/examples/reordered-preconditioned-solver/data/A.mtx b/examples/reordered-preconditioned-solver/data/A.mtx new file mode 100644 index 00000000000..c67437da567 --- /dev/null +++ b/examples/reordered-preconditioned-solver/data/A.mtx @@ -0,0 +1,114 @@ +%%MatrixMarket matrix coordinate integer symmetric +%------------------------------------------------------------------------------- +% UF Sparse Matrix Collection, Tim Davis +% http://www.cise.ufl.edu/research/sparse/matrices/JGD_Trefethen/Trefethen_20b +% name: JGD_Trefethen/Trefethen_20b +% [Diagonal matrices with primes, Nick Trefethen, Oxford Univ.] +% id: 2203 +% date: 2008 +% author: N. Trefethen +% ed: J.-G. Dumas +% fields: name title A id date author ed kind notes +% kind: combinatorial problem +%------------------------------------------------------------------------------- +% notes: +% Diagonal matrices with primes, Nick Trefethen, Oxford Univ. +% From Jean-Guillaume Dumas' Sparse Integer Matrix Collection, +% http://ljk.imag.fr/membres/Jean-Guillaume.Dumas/simc.html +% +% Problem 7 of the Hundred-dollar, Hundred-digit Challenge Problems, +% SIAM News, vol 35, no. 1. +% +% 7. Let A be the 20,000 x 20,000 matrix whose entries are zero +% everywhere except for the primes 2, 3, 5, 7, . . . , 224737 along the +% main diagonal and the number 1 in all the positions A(i,j) with +% |i-j| = 1,2,4,8, . . . ,16384. What is the (1,1) entry of inv(A)? +% +% http://www.siam.org/news/news.php?id=388 +% +% Filename in JGD collection: Trefethen/trefethen_20__19_minor.sms +%------------------------------------------------------------------------------- +19 19 83 +1 1 3 +2 1 1 +3 1 1 +5 1 1 +9 1 1 +17 1 1 +2 2 5 +3 2 1 +4 2 1 +6 2 1 +10 2 1 +18 2 1 +3 3 7 +4 3 1 +5 3 1 +7 3 1 +11 3 1 +19 3 1 +4 4 11 +5 4 1 +6 4 1 +8 4 1 +12 4 1 +5 5 13 +6 5 1 +7 5 1 +9 5 1 +13 5 1 +6 6 17 +7 6 1 +8 6 1 +10 6 1 +14 6 1 +7 7 19 +8 7 1 +9 7 1 +11 7 1 +15 7 1 +8 8 23 +9 8 1 +10 8 1 +12 8 1 +16 8 1 +9 9 29 +10 9 1 +11 9 1 +13 9 1 +17 9 1 +10 10 31 +11 10 1 +12 10 1 +14 10 1 +18 10 1 +11 11 37 +12 11 1 +13 11 1 +15 11 1 +19 11 1 +12 12 41 +13 12 1 +14 12 1 +16 12 1 +13 13 43 +14 13 1 +15 13 1 +17 13 1 +14 14 47 +15 14 1 +16 14 1 +18 14 1 +15 15 53 +16 15 1 +17 15 1 +19 15 1 +16 16 59 +17 16 1 +18 16 1 +17 17 61 +18 17 1 +19 17 1 +18 18 67 +19 18 1 +19 19 71 diff --git a/examples/reordered-preconditioned-solver/data/b.mtx b/examples/reordered-preconditioned-solver/data/b.mtx new file mode 100644 index 00000000000..05d92ecc6f7 --- /dev/null +++ b/examples/reordered-preconditioned-solver/data/b.mtx @@ -0,0 +1,21 @@ +%%MatrixMarket matrix array real general +19 1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 diff --git a/examples/reordered-preconditioned-solver/data/x0.mtx b/examples/reordered-preconditioned-solver/data/x0.mtx new file mode 100644 index 00000000000..91d470cdbcd --- /dev/null +++ b/examples/reordered-preconditioned-solver/data/x0.mtx @@ -0,0 +1,21 @@ +%%MatrixMarket matrix array real general +19 1 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 diff --git a/examples/reordered-preconditioned-solver/doc/builds-on b/examples/reordered-preconditioned-solver/doc/builds-on new file mode 100644 index 00000000000..9b64c9bfd28 --- /dev/null +++ b/examples/reordered-preconditioned-solver/doc/builds-on @@ -0,0 +1 @@ +preconditioned-solver diff --git a/examples/reordered-preconditioned-solver/doc/intro.dox b/examples/reordered-preconditioned-solver/doc/intro.dox new file mode 100644 index 00000000000..aba0c05d0ed --- /dev/null +++ b/examples/reordered-preconditioned-solver/doc/intro.dox @@ -0,0 +1,5 @@ + +

Introduction

+ +

About the example

+ diff --git a/examples/reordered-preconditioned-solver/doc/kind b/examples/reordered-preconditioned-solver/doc/kind new file mode 100644 index 00000000000..53a96d5771f --- /dev/null +++ b/examples/reordered-preconditioned-solver/doc/kind @@ -0,0 +1 @@ +preconditioners diff --git a/examples/reordered-preconditioned-solver/doc/results.dox b/examples/reordered-preconditioned-solver/doc/results.dox new file mode 100644 index 00000000000..66dcde2c83b --- /dev/null +++ b/examples/reordered-preconditioned-solver/doc/results.dox @@ -0,0 +1,35 @@ +

Results

+This is the expected output: + +@code{.cpp} + +Solution (x): +%%MatrixMarket matrix array real general +19 1 +0.252218 +0.108645 +0.0662811 +0.0630433 +0.0384088 +0.0396536 +0.0402648 +0.0338935 +0.0193098 +0.0234653 +0.0211499 +0.0196413 +0.0199151 +0.0181674 +0.0162722 +0.0150714 +0.0107016 +0.0121141 +0.0123025 +Residual norm sqrt(r^T r): +%%MatrixMarket matrix array real general +1 1 +4.82005e-08 + +@endcode + +

Comments about programming and debugging

diff --git a/examples/reordered-preconditioned-solver/doc/short-intro b/examples/reordered-preconditioned-solver/doc/short-intro new file mode 100644 index 00000000000..6edaeec3931 --- /dev/null +++ b/examples/reordered-preconditioned-solver/doc/short-intro @@ -0,0 +1 @@ +The reordered preconditioned solver example. diff --git a/examples/reordered-preconditioned-solver/doc/tooltip b/examples/reordered-preconditioned-solver/doc/tooltip new file mode 100644 index 00000000000..70701fc58a8 --- /dev/null +++ b/examples/reordered-preconditioned-solver/doc/tooltip @@ -0,0 +1 @@ +Use a reordering in Ginkgo. Solve a linear system with a preconditioner. diff --git a/examples/reordered-preconditioned-solver/reordered-preconditioned-solver.cpp b/examples/reordered-preconditioned-solver/reordered-preconditioned-solver.cpp new file mode 100644 index 00000000000..104ebcb9ea6 --- /dev/null +++ b/examples/reordered-preconditioned-solver/reordered-preconditioned-solver.cpp @@ -0,0 +1,111 @@ +// SPDX-FileCopyrightText: 2017-2023 The Ginkgo authors +// +// SPDX-License-Identifier: BSD-3-Clause + + +#include + + +#include +#include +#include +#include +#include "ginkgo/core/matrix/permutation.hpp" + + +int main(int argc, char* argv[]) +{ + // Some shortcuts + using ValueType = double; + using RealValueType = gko::remove_complex; + using IndexType = int; + + using vec = gko::matrix::Dense; + using real_vec = gko::matrix::Dense; + using mtx = gko::matrix::Csr; + using cg = gko::solver::Cg; + using bj = gko::preconditioner::Jacobi; + + // Print version information + std::cout << gko::version_info::get() << std::endl; + + // Figure out where to run the code + if (argc == 2 && (std::string(argv[1]) == "--help")) { + std::cerr << "Usage: " << argv[0] << " [executor]" << std::endl; + std::exit(-1); + } + + // Figure out where to run the code + const auto executor_string = argc >= 2 ? argv[1] : "reference"; + std::map()>> + exec_map{ + {"omp", [] { return gko::OmpExecutor::create(); }}, + {"cuda", + [] { + return gko::CudaExecutor::create(0, + gko::OmpExecutor::create()); + }}, + {"hip", + [] { + return gko::HipExecutor::create(0, gko::OmpExecutor::create()); + }}, + {"dpcpp", + [] { + return gko::DpcppExecutor::create(0, + gko::OmpExecutor::create()); + }}, + {"reference", [] { return gko::ReferenceExecutor::create(); }}}; + + // executor where Ginkgo will perform the computation + const auto exec = exec_map.at(executor_string)(); // throws if not valid + + // Read data + auto A = share(gko::read(std::ifstream("data/A.mtx"), exec)); + auto b = gko::read(std::ifstream("data/b.mtx"), exec); + auto x = gko::read(std::ifstream("data/x0.mtx"), exec); + + auto reordering = + gko::experimental::reorder::Rcm::build().on(exec)->generate( + A); + // Permute matrix and vectors + auto A_reordered = share(A->permute(reordering)); + auto b_reordered = b->permute(reordering, gko::matrix::permute_mode::rows); + // this reordering is not necessary, but it maps the initial guess to the + // unreordered case + auto x_reordered = x->permute(reordering, gko::matrix::permute_mode::rows); + + const RealValueType reduction_factor{1e-7}; + // Create solver factory + auto solver_gen = + cg::build() + .with_criteria(gko::stop::Iteration::build().with_max_iters(20u), + gko::stop::ResidualNorm::build() + .with_reduction_factor(reduction_factor)) + // Add preconditioner, these 2 lines are the only + // difference from the simple solver example + .with_preconditioner(bj::build().with_max_block_size(8u)) + .on(exec); + // Create solver + auto solver = solver_gen->generate(A_reordered); + + // Solve system + solver->apply(b_reordered, x_reordered); + + // Revert permutation to get the unpermuted solution + x_reordered->permute(reordering, x, + gko::matrix::permute_mode::inverse_rows); + + // Print solution + std::cout << "Solution (x):\n"; + write(std::cout, x); + + // Calculate residual + auto one = gko::initialize({1.0}, exec); + auto neg_one = gko::initialize({-1.0}, exec); + auto res = gko::initialize({0.0}, exec); + A->apply(one, x, neg_one, b); + b->compute_norm2(res); + + std::cout << "Residual norm sqrt(r^T r):\n"; + write(std::cout, res); +}