From f520fe5f0384c793c40cd2a2fd9e257eb7bac718 Mon Sep 17 00:00:00 2001 From: Emily Boudreaux Date: Mon, 9 Nov 2026 13:25:52 -0500 Subject: [PATCH] feat(solver/linear): Implemented a few linear solver frontends As with much of the code, the actual solvers are implemented (through MFEM) in the backend. This provides an interface so we can more easily swap out the linear solvers on the backend in the future --- src/include/serif/solver/linear/failure.hpp | 19 +++++++ src/include/serif/solver/linear/linear.hpp | 55 +++++++++++++++++++++ src/lib/serif/solver/linear/linear.cpp | 32 ++++++++++++ 3 files changed, 106 insertions(+) create mode 100644 src/include/serif/solver/linear/failure.hpp create mode 100644 src/include/serif/solver/linear/linear.hpp create mode 100644 src/lib/serif/solver/linear/linear.cpp diff --git a/src/include/serif/solver/linear/failure.hpp b/src/include/serif/solver/linear/failure.hpp new file mode 100644 index 00000000..c2fd99c3 --- /dev/null +++ b/src/include/serif/solver/linear/failure.hpp @@ -0,0 +1,19 @@ +#pragma once + +#include +#include "serif/utils/error/nonfinite.hpp" + +namespace serif::solver::linear { + /** + * @brief Bridge a failed checked inner solve through MFEM's void operator-application interface. + * Only this numerical failure is intercepted by the checked backend. Dimension, ownership, and + * other callback errors retain their original exception types and reach the library's caller. + */ + class ApplicationFailure : public std::runtime_error { + public: + using std::runtime_error::runtime_error; + }; + + /** @brief Stop Krylov before a nonfinite application reaches MFEM internal arithmetic assertions. */ + using NonFiniteApplication = utils::errors::NonFiniteValue; +} diff --git a/src/include/serif/solver/linear/linear.hpp b/src/include/serif/solver/linear/linear.hpp new file mode 100644 index 00000000..12d2a884 --- /dev/null +++ b/src/include/serif/solver/linear/linear.hpp @@ -0,0 +1,55 @@ +#pragma once + +#include + +namespace serif::solver::linear { + /** @brief Numerical termination of a checked Krylov solve. Invalid configuration remains a programming error. */ + enum class Status { converged, maximum_iterations, breakdown, non_finite, backend_failure }; + + /** @brief Accuracy is measured against the current right-hand side, independently of the initial guess. */ + struct Control { + double relative_tolerance{1.0e-8}; + double absolute_tolerance{0.0}; + int maximum_iterations{100}; + + /** @brief Reject invalid tolerances or iteration budgets before numerical work begins. */ + void validate() const; + + /** @brief Return max(absolute tolerance, relative tolerance times the right-hand-side norm). */ + [[nodiscard]] double convergence_threshold(double right_hand_side_norm) const; + }; + + /** @brief Flexible GMRES accepts both fixed linear inverses and variable inner preconditioners. */ + struct FGMRES { + int restart_length{160}; + int print_level{-1}; + + /** @brief Check the Krylov restart length and MFEM diagnostic verbosity. */ + void validate() const; + }; + + /** + * @brief Bounded diagnostics from one solve, including an independently evaluated true residual. + * @note Operator counts include warm-start safeguarding and final verification. Iterations count + * completed preconditioner applications; restarts count additional Krylov cycles entered. + */ + struct Report { + Status status{Status::backend_failure}; + Control control{}; + int iterations{0}; + int restarts{0}; + double right_hand_side_norm{0.0}; + double initial_residual_norm{0.0}; + double reported_residual_norm{0.0}; + double true_residual_norm{0.0}; + double relative_true_residual_norm{0.0}; + std::uint64_t operator_applications{0}; + std::uint64_t preconditioner_applications{0}; + double solve_seconds{0.0}; + double operator_seconds{0.0}; + double preconditioner_seconds{0.0}; + + /** @brief Whether the final correction satisfied the independently checked residual tolerance. */ + [[nodiscard]] bool converged() const noexcept { return status == Status::converged; } + }; +} diff --git a/src/lib/serif/solver/linear/linear.cpp b/src/lib/serif/solver/linear/linear.cpp new file mode 100644 index 00000000..5930c13f --- /dev/null +++ b/src/lib/serif/solver/linear/linear.cpp @@ -0,0 +1,32 @@ +#include "serif/solver/linear/linear.hpp" + +#include +#include +#include +#include "serif/utils/error/errors.hpp" + +namespace serif::solver::linear { + void Control::validate() const { + utils::errors::serif_ensure( + std::isfinite(relative_tolerance) && relative_tolerance >= 0.0 && + std::isfinite(absolute_tolerance) && absolute_tolerance >= 0.0 && maximum_iterations > 0, + "A linear solve needs finite nonnegative tolerances and a positive iteration budget; received {}, {}, and {}.", + relative_tolerance, absolute_tolerance, maximum_iterations); + } + + double Control::convergence_threshold(double right_hand_side_norm) const { + validate(); + const double relative_threshold = relative_tolerance * right_hand_side_norm; + utils::errors::serif_ensure( + std::isfinite(right_hand_side_norm) && right_hand_side_norm >= 0.0 && std::isfinite(relative_threshold), + "A linear convergence threshold needs a finite nonnegative right-hand-side norm and finite product with the relative tolerance; received {}.", + right_hand_side_norm); + return std::max(absolute_tolerance, relative_threshold); + } + + void FGMRES::validate() const { + utils::errors::serif_ensure( + restart_length > 0 && print_level >= -1 && print_level <= 3, + "FGMRES needs a positive restart length and print level between -1 and 3; received {} and {}.", restart_length, print_level); + } +}