diff --git a/CHANGELOG.md b/CHANGELOG.md index 3c760b73..fdeb1ade 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -4,6 +4,10 @@ - Added classes and tests for permutation, Ruiz scaling, Cholesky factorization, Schur complement conjugate gradient and matrix multiplication and addition. +- Changed random number generation int tests to be C++ style and fixed-seed, to avoid random failures. + +- Added Schur Complement Conjugate Gradient class. + ## Changes to Re::Solve since release 0.99.2 - Added cmake-format. diff --git a/resolve/Common.hpp b/resolve/Common.hpp index af30adad..5d5ae88a 100644 --- a/resolve/Common.hpp +++ b/resolve/Common.hpp @@ -14,11 +14,12 @@ namespace ReSolve namespace constants { - constexpr real_type ZERO = 0.0; - constexpr real_type ONE = 1.0; - constexpr real_type TWO = 2.0; - constexpr real_type HALF = 0.5; - constexpr real_type MINUS_ONE = -1.0; + constexpr real_type ZERO = 0.0; + constexpr real_type ONE = 1.0; + constexpr real_type TWO = 2.0; + constexpr real_type HALF = 0.5; + constexpr real_type MINUS_ONE = -1.0; + constexpr index_type SEED = 12345; constexpr real_type MACHINE_EPSILON = std::numeric_limits::epsilon(); } // namespace constants diff --git a/resolve/hykkt/cholesky/CholeskySolverHip.cpp b/resolve/hykkt/cholesky/CholeskySolverHip.cpp index cc8ae4fa..358ebea9 100644 --- a/resolve/hykkt/cholesky/CholeskySolverHip.cpp +++ b/resolve/hykkt/cholesky/CholeskySolverHip.cpp @@ -60,6 +60,7 @@ namespace ReSolve } A_chol_ = convertToCholmod(A); A_ = A; + mem_.deviceSynchronize(); } void CholeskySolverHip::symbolicAnalysis() @@ -70,6 +71,7 @@ namespace ReSolve { out::error() << "Cholesky symbolic analysis failed with status: " << Common_.status << "\n"; } + mem_.deviceSynchronize(); } /** @@ -156,6 +158,7 @@ namespace ReSolve out::error() << "Refactorization step failed with status: " << status << "\n"; } L_->setUpdated(memory::DEVICE); + mem_.deviceSynchronize(); } } diff --git a/tests/unit/hykkt/HykktCholeskyTests.hpp b/tests/unit/hykkt/HykktCholeskyTests.hpp index ec214381..cf1a9b8a 100644 --- a/tests/unit/hykkt/HykktCholeskyTests.hpp +++ b/tests/unit/hykkt/HykktCholeskyTests.hpp @@ -4,6 +4,7 @@ #include #include #include +#include #include #include #include @@ -26,8 +27,8 @@ namespace ReSolve class HykktCholeskyTests : public TestBase { public: - HykktCholeskyTests(memory::MemorySpace memspace, MatrixHandler& matrixHandler) - : memspace_(memspace), matrixHandler_(matrixHandler) + HykktCholeskyTests(memory::MemorySpace memspace, MatrixHandler& matrixHandler, std::mt19937& generator) + : memspace_(memspace), matrixHandler_(matrixHandler), generator_(generator) { cholmod_start(&Common); } @@ -127,7 +128,6 @@ namespace ReSolve solver.addMatrixInfo(A); solver.symbolicAnalysis(); solver.numericalFactorization(); - // Generate a random vector x_expected and compute b = A * x_expected vector::Vector* x_expected = randomVector(n); @@ -242,7 +242,8 @@ namespace ReSolve // reset values for (size_t j = 0; j < L->nzmax; j++) { - static_cast(L->x)[j] = 2.0 * static_cast(rand()) / RAND_MAX - 1.0; + std::uniform_real_distribution distribution(-1.0, 1.0); + static_cast(L->x)[j] = 2.0 * distribution(generator_) - 1.0; } cholmod_free_sparse(&L_tr, &Common); cholmod_free_sparse(&A_chol, &Common); @@ -266,6 +267,7 @@ namespace ReSolve private: ReSolve::memory::MemorySpace memspace_; MatrixHandler& matrixHandler_; + std::mt19937& generator_; cholmod_common Common; @@ -281,14 +283,16 @@ namespace ReSolve L_p[i + 1] = L_p[i]; for (size_t j = i; j < n; ++j) { - // with probability 'density', add a non-zero entry - if (i == j || static_cast(rand()) / RAND_MAX < density) + std::uniform_real_distribution rand_dist(0.0, 1.0); + if (i == j || rand_dist(generator_) < density) { L_i.push_back((int) j); // force diagonal entry to be non-zero double value = 2.0; if (i != j) { + std::uniform_real_distribution value_dist(-1.0, 1.0); + value = 2.0 * value_dist(generator_); value = 2.0 * static_cast(rand()) / RAND_MAX - 1.0; } L_x.push_back(value); @@ -311,9 +315,10 @@ namespace ReSolve { vector::Vector* v = new vector::Vector(n); v->allocate(memory::HOST); + std::uniform_real_distribution distribution(0.0, 1.0); for (index_type i = 0; i < n; ++i) { - v->getData(memory::HOST)[i] = static_cast(rand()) / RAND_MAX; + v->getData(memory::HOST)[i] = distribution(generator_); } v->setDataUpdated(memory::HOST); if (memspace_ == memory::DEVICE) diff --git a/tests/unit/hykkt/HykktSCCGTests.hpp b/tests/unit/hykkt/HykktSCCGTests.hpp index 2d853474..035c82cd 100644 --- a/tests/unit/hykkt/HykktSCCGTests.hpp +++ b/tests/unit/hykkt/HykktSCCGTests.hpp @@ -6,6 +6,7 @@ #pragma once #include +#include #include #include @@ -35,9 +36,10 @@ namespace ReSolve * @param[in] memspace Memory space for the test (HOST or DEVICE). * @param[in] matrix_handler Reference to a matrix handler for the selected backend. * @param[in] vector_handler Reference to a vector handler for the selected backend. + * @param[in] generator Reference to a C++ random number generator. */ - HykktSchurComplementConjugateGradientTests(memory::MemorySpace memspace, MatrixHandler& matrix_handler, VectorHandler& vector_handler) - : memspace_(memspace), matrix_handler_(matrix_handler), vector_handler_(vector_handler) + HykktSchurComplementConjugateGradientTests(memory::MemorySpace memspace, MatrixHandler& matrix_handler, VectorHandler& vector_handler, std::mt19937& generator) + : memspace_(memspace), matrix_handler_(matrix_handler), vector_handler_(vector_handler), generator_(generator) { } @@ -130,16 +132,18 @@ namespace ReSolve memory::MemorySpace memspace_; ///< Memory space used by the test. MatrixHandler& matrix_handler_; ///< Backend-specific matrix handler. VectorHandler& vector_handler_; ///< Backend-specific vector handler. + std::mt19937& generator_; ///< C++ random number generator. /** - * @brief Generate a random vector of doubles between 0 and RAND_MAX. Copied from HykktCholeskyTests.hpp. + * @brief Generate a random vector of doubles between 0 and 1. Copied from HykktCholeskyTests.hpp. * @param[in] vec Target vector to write to. */ void randomVector(vector::Vector* vec) { + std::uniform_real_distribution distribution(0.0, 1.0); for (index_type i = 0; i < vec->getSize(); ++i) { - vec->getData(memory::HOST)[i] = static_cast(rand()) / RAND_MAX; + vec->getData(memory::HOST)[i] = distribution(generator_); } vec->setDataUpdated(memory::HOST); if (memspace_ == memory::DEVICE) diff --git a/tests/unit/hykkt/runHykktCholeskyTests.cpp b/tests/unit/hykkt/runHykktCholeskyTests.cpp index d05fb560..55ebeb7d 100644 --- a/tests/unit/hykkt/runHykktCholeskyTests.cpp +++ b/tests/unit/hykkt/runHykktCholeskyTests.cpp @@ -7,8 +7,10 @@ */ #include #include +#include #include +#include "resolve/Common.hpp" #include "tests/unit/hykkt/HykktCholeskyTests.hpp" /** @@ -24,9 +26,9 @@ void runTests(const std::string& backend, ReSolve::memory::MemorySpace memspace, WorkspaceType workspace; workspace.initializeHandles(); - ReSolve::MatrixHandler handler(&workspace); - - ReSolve::tests::HykktCholeskyTests test(memspace, handler); + ReSolve::MatrixHandler handler(&workspace); + std::mt19937 generator(ReSolve::constants::SEED); // set random seed for reproducibility + ReSolve::tests::HykktCholeskyTests test(memspace, handler, generator); result += test.minimalCorrectness(); handler.setValuesChanged(true, memspace); diff --git a/tests/unit/hykkt/runHykktSCCGTests.cpp b/tests/unit/hykkt/runHykktSCCGTests.cpp index b3462676..33f0e385 100644 --- a/tests/unit/hykkt/runHykktSCCGTests.cpp +++ b/tests/unit/hykkt/runHykktSCCGTests.cpp @@ -5,6 +5,7 @@ */ #include #include +#include #include #include @@ -35,7 +36,8 @@ void runTests(const std::string& backend, ReSolve::memory::MemorySpace memspace, workspace.initializeHandles(); ReSolve::MatrixHandler matrix_handler(&workspace); ReSolve::VectorHandler vector_handler(&workspace); - ReSolve::tests::HykktSchurComplementConjugateGradientTests test(memspace, matrix_handler, vector_handler); + std::mt19937 generator(ReSolve::constants::SEED); + ReSolve::tests::HykktSchurComplementConjugateGradientTests test(memspace, matrix_handler, vector_handler, generator); result += test.SCCGTest();