extract the split residual update and update test

Co-authored-by: Marcel Koch <[email protected]>

core/solver/chebyshev.cpp

@@ -40,6 +40,7 @@ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 
 #include "core/distributed/helpers.hpp"
 #include "core/solver/ir_kernels.hpp"
+#include "core/solver/residual_update.hpp"
 #include "core/solver/solver_base.hpp"
 #include "core/solver/solver_boilerplate.hpp"
 
@@ -266,52 +267,20 @@ void Chebyshev<ValueType>::apply_dense_impl(const VectorType* dense_b,
     int iter = -1;
     while (true) {
         ++iter;
-        if (iter == 0) {
-            // In iter 0, the iteration and residual are updated.
-            bool all_stopped = stop_criterion->update()
-                                   .num_iterations(iter)
-                                   .residual(residual_ptr)
-                                   .solution(dense_x)
-                                   .check(relative_stopping_id, true,
-                                          &stop_status, &one_changed);
-            this->template log<log::Logger::iteration_complete>(
-                this, dense_b, dense_x, iter, residual_ptr, nullptr, nullptr,
-                &stop_status, all_stopped);
-            if (all_stopped) {
-                break;
-            }
-        } else {
-            // In the other iterations, the residual can be updated separately.
-            bool all_stopped = stop_criterion->update()
-                                   .num_iterations(iter)
-                                   .solution(dense_x)
-                                   // we have the residual check later
-                                   .ignore_residual_check(true)
-                                   .check(relative_stopping_id, false,
-                                          &stop_status, &one_changed);
-            if (all_stopped) {
-                this->template log<log::Logger::iteration_complete>(
-                    this, dense_b, dense_x, iter, nullptr, nullptr, nullptr,
-                    &stop_status, all_stopped);
-                break;
-            }
-            residual_ptr = residual;
-            // residual = b - A * x
-            residual->copy_from(dense_b);
-            this->get_system_matrix()->apply(neg_one_op, dense_x, one_op,
-                                             residual);
-            all_stopped = stop_criterion->update()
-                              .num_iterations(iter)
-                              .residual(residual_ptr)
-                              .solution(dense_x)
-                              .check(relative_stopping_id, true, &stop_status,
-                                     &one_changed);
+        auto log_func = [this](auto solver, auto dense_b, auto dense_x,
+                               auto iter, auto residual_ptr,
+                               array<stopping_status>& stop_status,
+                               bool all_stopped) {
             this->template log<log::Logger::iteration_complete>(
-                this, dense_b, dense_x, iter, residual_ptr, nullptr, nullptr,
+                solver, dense_b, dense_x, iter, residual_ptr, nullptr, nullptr,
                 &stop_status, all_stopped);
-            if (all_stopped) {
-                break;
-            }
+        };
+        bool all_stopped = residual_update(
+            this, iter, one_op, neg_one_op, dense_b, dense_x, residual,
+            residual_ptr, stop_criterion, relative_stopping_id, stop_status,
+            one_changed, log_func);
+        if (all_stopped) {
+            break;
         }
 
         if (solver_->apply_uses_initial_guess()) {

core/solver/ir.cpp

@@ -40,6 +40,7 @@ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 
 #include "core/distributed/helpers.hpp"
 #include "core/solver/ir_kernels.hpp"
+#include "core/solver/residual_update.hpp"
 #include "core/solver/solver_base.hpp"
 #include "core/solver/solver_boilerplate.hpp"
 
@@ -223,54 +224,23 @@ void Ir<ValueType>::apply_dense_impl(const VectorType* dense_b,
     while (true) {
         ++iter;
 
-        if (iter == 0) {
-            // In iter 0, the iteration and residual are updated.
-            bool all_stopped = stop_criterion->update()
-                                   .num_iterations(iter)
-                                   .residual(residual_ptr)
-                                   .solution(dense_x)
-                                   .check(relative_stopping_id, true,
-                                          &stop_status, &one_changed);
+        auto log_func = [this](auto solver, auto dense_b, auto dense_x,
+                               auto iter, auto residual_ptr,
+                               array<stopping_status>& stop_status,
+                               bool all_stopped) {
             this->template log<log::Logger::iteration_complete>(
-                this, dense_b, dense_x, iter, residual_ptr, nullptr, nullptr,
+                solver, dense_b, dense_x, iter, residual_ptr, nullptr, nullptr,
                 &stop_status, all_stopped);
-            if (all_stopped) {
-                break;
-            }
-        } else {
-            // In the other iterations, the residual can be updated separately.
-            bool all_stopped = stop_criterion->update()
-                                   .num_iterations(iter)
-                                   .solution(dense_x)
-                                   // we have the residual check later
-                                   .ignore_residual_check(true)
-                                   .check(relative_stopping_id, false,
-                                          &stop_status, &one_changed);
-            if (all_stopped) {
-                this->template log<log::Logger::iteration_complete>(
-                    this, dense_b, dense_x, iter, nullptr, nullptr, nullptr,
-                    &stop_status, all_stopped);
-                break;
-            }
-            residual_ptr = residual;
-            // residual = b - A * x
-            residual->copy_from(dense_b);
-            this->get_system_matrix()->apply(neg_one_op, dense_x, one_op,
-                                             residual);
-            all_stopped = stop_criterion->update()
-                              .num_iterations(iter)
-                              .residual(residual_ptr)
-                              .solution(dense_x)
-                              .check(relative_stopping_id, true, &stop_status,
-                                     &one_changed);
-            this->template log<log::Logger::iteration_complete>(
-                this, dense_b, dense_x, iter, residual_ptr, nullptr, nullptr,
-                &stop_status, all_stopped);
-            if (all_stopped) {
-                break;
-            }
+        };
+        bool all_stopped = residual_update(
+            this, iter, one_op, neg_one_op, dense_b, dense_x, residual,
+            residual_ptr, stop_criterion, relative_stopping_id, stop_status,
+            one_changed, log_func);
+        if (all_stopped) {
+            break;
         }
 
+
         if (solver_->apply_uses_initial_guess()) {
             // Use the inner solver to solve
             // A * inner_solution = residual

core/solver/residual_update.hpp

@@ -0,0 +1,103 @@
+/*******************************<GINKGO LICENSE>******************************
+Copyright (c) 2017-2023, the Ginkgo authors
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions
+are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above copyright
+notice, this list of conditions and the following disclaimer in the
+documentation and/or other materials provided with the distribution.
+
+3. Neither the name of the copyright holder nor the names of its
+contributors may be used to endorse or promote products derived from
+this software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
+IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
+HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
+DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
+THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
+(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
+OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+******************************<GINKGO LICENSE>*******************************/
+
+#ifndef GKO_CORE_SOLVER_RESIDUAL_UPDATE_HPP_
+#define GKO_CORE_SOLVER_RESIDUAL_UPDATE_HPP_
+
+
+#include <ginkgo/core/base/array.hpp>
+#include <ginkgo/core/matrix/dense.hpp>
+#include <ginkgo/core/stop/criterion.hpp>
+
+
+namespace gko {
+namespace solver {
+
+
+template <typename SolverType, typename VectorType, typename ScalarType,
+          typename LogFunc>
+bool residual_update(SolverType* solver, int iter, const ScalarType* one_op,
+                     const ScalarType* neg_one_op, const VectorType* dense_b,
+                     VectorType* dense_x, VectorType* residual,
+                     const VectorType*& residual_ptr,
+                     std::unique_ptr<gko::stop::Criterion>& stop_criterion,
+                     uint8 relative_stopping_id,
+                     array<stopping_status>& stop_status, bool& one_changed,
+                     LogFunc log)
+{
+    if (iter == 0) {
+        // In iter 0, the iteration and residual are updated.
+        bool all_stopped =
+            stop_criterion->update()
+                .num_iterations(iter)
+                .residual(residual_ptr)
+                .solution(dense_x)
+                .check(relative_stopping_id, true, &stop_status, &one_changed);
+        log(solver, dense_b, dense_x, iter, residual_ptr, stop_status,
+            all_stopped);
+        return all_stopped;
+    } else {
+        // In the other iterations, the residual can be updated separately.
+        bool all_stopped =
+            stop_criterion->update()
+                .num_iterations(iter)
+                .solution(dense_x)
+                // we have the residual check later
+                .ignore_residual_check(true)
+                .check(relative_stopping_id, false, &stop_status, &one_changed);
+        if (all_stopped) {
+            log(solver, dense_b, dense_x, iter, nullptr, stop_status,
+                all_stopped);
+            return all_stopped;
+        }
+        residual_ptr = residual;
+        // residual = b - A * x
+        residual->copy_from(dense_b);
+        solver->get_system_matrix()->apply(neg_one_op, dense_x, one_op,
+                                           residual);
+        all_stopped =
+            stop_criterion->update()
+                .num_iterations(iter)
+                .residual(residual_ptr)
+                .solution(dense_x)
+                .check(relative_stopping_id, true, &stop_status, &one_changed);
+        log(solver, dense_b, dense_x, iter, residual_ptr, stop_status,
+            all_stopped);
+        return all_stopped;
+    }
+}
+
+
+}  // namespace solver
+}  // namespace gko
+
+#endif  // GKO_CORE_SOLVER_RESIDUAL_UPDATE_HPP_

core/test/solver/chebyshev.cpp

@@ -195,10 +195,7 @@ TYPED_TEST(Chebyshev, CanSetInnerSolverInFactory)
     auto chebyshev_factory =
         Solver::build()
             .with_criteria(
-                gko::stop::Iteration::build().with_max_iters(3u).on(this->exec),
-                gko::stop::ResidualNorm<value_type>::build()
-                    .with_reduction_factor(r<value_type>::value)
-                    .on(this->exec))
+                gko::stop::Iteration::build().with_max_iters(3u).on(this->exec)
             .with_solver(
                 Solver::build()
                     .with_criteria(

include/ginkgo/core/solver/chebyshev.hpp

@@ -158,8 +158,8 @@ class Chebyshev : public EnableLinOp<Chebyshev<ValueType>>,
             GKO_FACTORY_PARAMETER_VECTOR(criteria, nullptr);
 
         /**
-         * Inner solver factory. If not provided this will result in a
-         * non-preconditioned Chebyshev iteration.
+         * Inner solver (preconditioner) factory. If not provided this will
+         * result in a non-preconditioned Chebyshev iteration.
          */
         std::shared_ptr<const LinOpFactory> GKO_FACTORY_PARAMETER_SCALAR(
             solver, nullptr);
@@ -172,8 +172,9 @@ class Chebyshev : public EnableLinOp<Chebyshev<ValueType>>,
             generated_solver, nullptr);
 
         /**
-         * The pair of foci of ellipse. It is usually be {lower bound of eigval,
-         * upper bound of eigval} for real matrices.
+         * The pair of foci of ellipse, which covers the eigenvalues of
+         * preconditioned system. It is usually be {lower bound of eigval, upper
+         * bound of eigval} of preconditioned real matrices.
          */
         std::pair<value_type, value_type> GKO_FACTORY_PARAMETER_VECTOR(
             foci, value_type{0}, value_type{1});