Skip to content

Avoid reallocation for self residual norm calculation #1898

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

Open
yhmtsai wants to merge 1 commit into
base: develop
Choose a base branch
from
Open

Avoid reallocation for self residual norm calculation #1898

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
18 changes: 14 additions & 4 deletions core/stop/residual_norm.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,10 +1,14 @@
// SPDX-FileCopyrightText: 2017 - 2024 The Ginkgo authors
// SPDX-FileCopyrightText: 2017 - 2025 The Ginkgo authors
//
// SPDX-License-Identifier: BSD-3-Clause

#include "ginkgo/core/stop/residual_norm.hpp"

#include <type_traits>

#include <ginkgo/core/base/precision_dispatch.hpp>
#include <ginkgo/core/distributed/vector.hpp>
#include <ginkgo/core/matrix/dense.hpp>

#include "core/base/dispatch_helper.hpp"
#include "core/components/fill_array_kernels.hpp"
Expand Down Expand Up @@ -183,9 +187,15 @@ bool ResidualNormBase<ValueType>::check_impl(
auto exec = this->get_executor();
norm_dispatch<ValueType>(
[&](auto dense_b, auto dense_x) {
auto dense_r = dense_b->clone();
system_matrix_->apply(neg_one_, dense_x, one_, dense_r);
dense_r->compute_norm2(u_dense_tau_, reduction_tmp_);
using dense_b_type = std::remove_cv_t<
std::remove_reference_t<decltype(*dense_b)>>;
Comment on lines +190 to +191
Copy link
Member

@MarcelKoch MarcelKoch Jul 21, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Isn't this the same as

Suggested change
using dense_b_type = std::remove_cv_t<
std::remove_reference_t<decltype(*dense_b)>>;
using dense_b_type = std::decay_t<decltype(*dense_b)>>;

if (!this->rhs_) {
this->rhs_ = dense_b_type::create_with_config_of(dense_b);
}
auto rhs = as<dense_b_type>(this->rhs_);
rhs->copy_from(dense_b);
system_matrix_->apply(neg_one_, dense_x, one_, rhs);
rhs->compute_norm2(u_dense_tau_, reduction_tmp_);
},
b_.get(), updater.solution_);
dense_tau = u_dense_tau_.get();
Expand Down
4 changes: 3 additions & 1 deletion include/ginkgo/core/stop/residual_norm.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
// SPDX-FileCopyrightText: 2017 - 2024 The Ginkgo authors
// SPDX-FileCopyrightText: 2017 - 2025 The Ginkgo authors
//
// SPDX-License-Identifier: BSD-3-Clause

Expand Down Expand Up @@ -86,6 +86,8 @@ class ResidualNormBase
std::shared_ptr<const Vector> neg_one_{};
// workspace for reduction
mutable gko::array<char> reduction_tmp_;
// temporary rhs for residual computation
mutable std::shared_ptr<LinOp> rhs_{};
Comment on lines 88 to +90
Copy link
Member

@upsj upsj Jul 21, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think we should avoid adding more mutable functionality, and instead either find a way to include the stopping criteria in the workspace setup, or hide it behind an opaque pointer, so at least there is no mutability in the public interface. We could do this by adding a workspace parameter to the check function and having the old version create a temporary workspace.

};


Expand Down
126 changes: 120 additions & 6 deletions reference/test/stop/residual_norm_kernels.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
// SPDX-FileCopyrightText: 2017 - 2024 The Ginkgo authors
// SPDX-FileCopyrightText: 2017 - 2025 The Ginkgo authors
//
// SPDX-License-Identifier: BSD-3-Clause

Expand All @@ -13,7 +13,19 @@
#include "core/test/utils.hpp"


namespace {
class AllocationLogger : public gko::log::Logger {
public:
mutable int count = 0;

protected:
void on_allocation_completed(const gko::Executor* exec,
const gko::size_type& num_bytes,
const gko::uintptr& location) const override
{
std::cout << num_bytes << std::endl;
++count;
}
};


template <typename T>
Expand Down Expand Up @@ -42,7 +54,7 @@ class ResidualNorm : public ::testing::Test {
std::unique_ptr<typename gko::stop::ResidualNorm<T>::Factory> rhs_factory_;
std::unique_ptr<typename gko::stop::ResidualNorm<T>::Factory> rel_factory_;
std::unique_ptr<typename gko::stop::ResidualNorm<T>::Factory> abs_factory_;
std::shared_ptr<const gko::Executor> exec_;
std::shared_ptr<gko::Executor> exec_;
};

TYPED_TEST_SUITE(ResidualNorm, gko::test::ValueTypes, TypenameNameGenerator);
Expand Down Expand Up @@ -417,6 +429,111 @@ TYPED_TEST(ResidualNorm, SelfCalculatesAndWaitsTillResidualGoal)
}


TYPED_TEST(ResidualNorm, SelfCalculatesWithoutReallocation)
{
using Mtx = typename TestFixture::Mtx;
using NormVector = typename TestFixture::NormVector;
using T = TypeParam;
using T_nc = gko::remove_complex<TypeParam>;
auto initial_res = gko::initialize<Mtx>({100.0}, this->exec_);
auto system_mtx = share(gko::initialize<Mtx>({1.0}, this->exec_));

T rhs_val = 10.0;
std::shared_ptr<gko::LinOp> rhs =
gko::initialize<Mtx>({rhs_val}, this->exec_);
auto rhs_criterion = this->rhs_factory_->generate(system_mtx, rhs, nullptr,
initial_res.get());
auto rel_criterion = this->rel_factory_->generate(system_mtx, rhs, nullptr,
initial_res.get());
auto abs_criterion = this->abs_factory_->generate(system_mtx, rhs, nullptr,
initial_res.get());
{
Copy link
Member

@MarcelKoch MarcelKoch Jul 21, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Either these blocks need a SCOPED_TRACE or the above is extracted into a new test fixture and the blocks are fully separated.

auto solution = gko::initialize<Mtx>({rhs_val - T{10.0}}, this->exec_);
auto rhs_norm = gko::initialize<NormVector>({100.0}, this->exec_);
gko::as<Mtx>(rhs)->compute_norm2(rhs_norm);
constexpr gko::uint8 RelativeStoppingId{1};
bool one_changed{};
gko::array<gko::stopping_status> stop_status(this->exec_, 1);
stop_status.get_data()[0].reset();
auto logger = std::make_shared<AllocationLogger>();

ASSERT_FALSE(rhs_criterion->update().solution(solution).check(
RelativeStoppingId, true, &stop_status, &one_changed));

solution->at(0) = rhs_val - r<T>::value * T{1.1} * rhs_norm->at(0);
ASSERT_FALSE(rhs_criterion->update().solution(solution).check(
RelativeStoppingId, true, &stop_status, &one_changed));
ASSERT_EQ(stop_status.get_data()[0].has_converged(), false);
ASSERT_EQ(one_changed, false);

solution->at(0) = rhs_val - r<T>::value * T{0.5} * rhs_norm->at(0);
this->exec_->add_logger(logger);
ASSERT_TRUE(rhs_criterion->update().solution(solution).check(
RelativeStoppingId, true, &stop_status, &one_changed));
ASSERT_EQ(stop_status.get_data()[0].has_converged(), true);
ASSERT_EQ(one_changed, true);
ASSERT_EQ(logger->count, 0);
this->exec_->remove_logger(logger);
}
{
T initial_norm = 100.0;
auto solution =
gko::initialize<Mtx>({rhs_val - initial_norm}, this->exec_);
constexpr gko::uint8 RelativeStoppingId{1};
bool one_changed{};
gko::array<gko::stopping_status> stop_status(this->exec_, 1);
stop_status.get_data()[0].reset();
auto logger = std::make_shared<AllocationLogger>();

ASSERT_FALSE(rel_criterion->update().solution(solution).check(
RelativeStoppingId, true, &stop_status, &one_changed));

solution->at(0) = rhs_val - r<T>::value * T{1.1} * initial_norm;
ASSERT_FALSE(rel_criterion->update().solution(solution).check(
RelativeStoppingId, true, &stop_status, &one_changed));
ASSERT_EQ(stop_status.get_data()[0].has_converged(), false);
ASSERT_EQ(one_changed, false);

solution->at(0) = rhs_val - r<T>::value * T{0.5} * initial_norm;
this->exec_->add_logger(logger);
ASSERT_TRUE(rel_criterion->update().solution(solution).check(
RelativeStoppingId, true, &stop_status, &one_changed));
ASSERT_EQ(stop_status.get_data()[0].has_converged(), true);
ASSERT_EQ(one_changed, true);
ASSERT_EQ(logger->count, 0);
this->exec_->remove_logger(logger);
}
{
auto solution = gko::initialize<Mtx>({rhs_val - T{100.0}}, this->exec_);
constexpr gko::uint8 RelativeStoppingId{1};
bool one_changed{};
gko::array<gko::stopping_status> stop_status(this->exec_, 1);
stop_status.get_data()[0].reset();
auto logger = std::make_shared<AllocationLogger>();

ASSERT_FALSE(abs_criterion->update().solution(solution).check(
RelativeStoppingId, true, &stop_status, &one_changed));

// TODO FIXME: NVHPC calculates different result of rhs - r*1.2 from
// rhs - tmp = rhs - (r * 1.2). https://godbolt.org/z/GrGE9PE67
solution->at(0) = rhs_val - r<T>::value * T{1.4};
ASSERT_FALSE(abs_criterion->update().solution(solution).check(
RelativeStoppingId, true, &stop_status, &one_changed));
ASSERT_EQ(stop_status.get_data()[0].has_converged(), false);
ASSERT_EQ(one_changed, false);

solution->at(0) = rhs_val - r<T>::value * T{0.5};
this->exec_->add_logger(logger);
ASSERT_TRUE(abs_criterion->update().solution(solution).check(
RelativeStoppingId, true, &stop_status, &one_changed));
ASSERT_EQ(stop_status.get_data()[0].has_converged(), true);
ASSERT_EQ(one_changed, true);
ASSERT_EQ(logger->count, 0);
this->exec_->remove_logger(logger);
}
}


TYPED_TEST(ResidualNorm, WaitsTillResidualGoalMultipleRHS)
{
using Mtx = typename TestFixture::Mtx;
Expand Down Expand Up @@ -1079,6 +1196,3 @@ TYPED_TEST(ResidualNormWithAbsolute, WaitsTillResidualGoalMultipleRHS)
ASSERT_EQ(stop_status.get_data()[1].has_converged(), true);
ASSERT_EQ(one_changed, true);
}


} // namespace
Loading