Scalar parameter handling

core/CMakeLists.txt

@@ -9,6 +9,7 @@ target_sources(ginkgo
     base/executor.cpp
     base/mtx_io.cpp
     base/perturbation.cpp
+    base/utils.cpp
     base/version.cpp
     factorization/ic.cpp
     factorization/ilu.cpp

core/base/utils.cpp

@@ -0,0 +1,63 @@
+/*******************************<GINKGO LICENSE>******************************
+Copyright (c) 2017-2021, 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_BASE_UTILS_HPP_
+#define GKO_CORE_BASE_UTILS_HPP_
+
+
+#include "utils.hpp"
+
+
+#include <memory>
+
+
+#include <ginkgo/core/base/executor.hpp>
+#include <ginkgo/core/matrix/dense.hpp>
+
+
+namespace gko {
+
+template <typename ValueType>
+std::unique_ptr<matrix::Dense<ValueType>> scalar_to_dense(
+    const ValueType val, std::shared_ptr<const Executor> exec)
+{
+    return initialize<matrix::Dense<ValueType>>({val}, exec);
+}
+#define GKO_DECLARE_SCALAR_TO_DENSE(_type)                 \
+    std::unique_ptr<matrix::Dense<_type>> scalar_to_dense( \
+        const _type, std::shared_ptr<const Executor>)
+GKO_INSTANTIATE_FOR_EACH_VALUE_TYPE(GKO_DECLARE_SCALAR_TO_DENSE);
+
+}  // namespace gko
+
+
+#endif  // GKO_CORE_BASE_UTILS_HPP_

core/test/base/CMakeLists.txt

@@ -21,4 +21,5 @@ ginkgo_create_test(range_accessors)
 ginkgo_create_thread_test(sanitizers)
 ginkgo_create_test(types)
 ginkgo_create_test(utils)
+ginkgo_create_test(value_typed_lin_op)
 ginkgo_create_test(version)

core/test/base/value_typed_lin_op.cpp

@@ -0,0 +1,176 @@
+/*******************************<GINKGO LICENSE>******************************
+Copyright (c) 2017-2021, 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>*******************************/
+
+#include <ginkgo/core/base/lin_op.hpp>
+
+
+#include <complex>
+#include <memory>
+#include <type_traits>
+
+
+#include <gtest/gtest.h>
+
+
+#include <ginkgo/core/base/math.hpp>
+#include <ginkgo/core/matrix/dense.hpp>
+
+
+#include "core/test/utils.hpp"
+
+
+namespace {
+
+
+template <typename T>
+T get_value_from_linop(const gko::LinOp *val);
+
+
+template <typename T>
+class DummyValueTypedLinOp
+    : public gko::EnableValueTypedLinOp<DummyValueTypedLinOp<T>, T>,
+      public gko::EnableCreateMethod<DummyValueTypedLinOp<T>> {
+    using Self = DummyValueTypedLinOp<T>;
+
+public:
+    explicit DummyValueTypedLinOp(std::shared_ptr<const gko::Executor> exec,
+                                  gko::dim<2> size = gko::dim<2>{})
+        : gko::EnableValueTypedLinOp<Self, T>(exec, size), value_()
+    {}
+    DummyValueTypedLinOp(std::shared_ptr<const gko::Executor> exec,
+                         gko::dim<2> size, T value)
+        : gko::EnableValueTypedLinOp<Self, T>(exec, size), value_(value)
+    {}
+
+    T get_value() const { return value_; }
+
+protected:
+    void apply_impl(const gko::LinOp *b, gko::LinOp *x) const override {}
+
+    void apply_impl(const gko::LinOp *alpha, const gko::LinOp *b,
+                    const gko::LinOp *beta, gko::LinOp *x) const override
+    {
+        T alpha_v = get_value_from_linop<T>(alpha);
+        T beta_v = get_value_from_linop<T>(beta);
+
+        gko::as<Self>(x)->value_ =
+            alpha_v * gko::as<Self>(this)->value_ * gko::as<Self>(b)->value_ +
+            beta_v * gko::as<Self>(x)->value_;
+    }
+
+    T value_;
+};
+
+
+template <typename T>
+T get_value_from_linop(const gko::LinOp *val)
+{
+    if (auto *dense = dynamic_cast<const gko::matrix::Dense<T> *>(val)) {
+        return dense->at(0, 0);
+    } else {
+        return (dynamic_cast<const DummyValueTypedLinOp<T> *>(val))
+            ->get_value();
+    }
+}
+
+
+template <typename T>
+class EnableValueTypedLinOp : public ::testing::Test {
+protected:
+    using dummy_type = DummyValueTypedLinOp<T>;
+    using value_type = T;
+
+    gko::dim<2> dim{1, 1};
+    std::shared_ptr<const gko::ReferenceExecutor> ref{
+        gko::ReferenceExecutor::create()};
+    std::unique_ptr<dummy_type> op{dummy_type::create(ref, dim, T{1.0})};
+    T alpha_v{2.0};
+    std::unique_ptr<dummy_type> alpha{dummy_type::create(ref, dim, alpha_v)};
+    T beta_v{3.0};
+    std::unique_ptr<dummy_type> beta{dummy_type::create(ref, dim, beta_v)};
+    std::unique_ptr<dummy_type> b{dummy_type::create(ref, dim, T{4.0})};
+    std::unique_ptr<dummy_type> x{dummy_type::create(ref, dim, T{5.0})};
+};
+
+
+TYPED_TEST_SUITE(EnableValueTypedLinOp, gko::test::ValueTypes);
+
+
+TYPED_TEST(EnableValueTypedLinOp, CanCallExtendedApplyImplLinopLinop)
+{
+    using value_type = typename TestFixture::value_type;
+    this->op->apply(lend(this->alpha), lend(this->b), lend(this->beta),
+                    lend(this->x));
+
+    ASSERT_EQ(this->x->get_value(), value_type{23.0});
+}
+
+
+TYPED_TEST(EnableValueTypedLinOp, CanCallExtendedApplyImplValueValue)
+{
+    auto reference_result = gko::clone(this->x);
+
+    this->op->apply(this->alpha_v, lend(this->b), this->beta_v, lend(this->x));
+    this->op->apply(lend(this->alpha), lend(this->b), lend(this->beta),
+                    lend(reference_result));
+
+    ASSERT_EQ(reference_result->get_value(), this->x->get_value());
+}
+
+
+TYPED_TEST(EnableValueTypedLinOp, CanCallExtendedApplyImplValueLinop)
+{
+    auto reference_result = gko::clone(this->x);
+
+    this->op->apply(this->alpha_v, lend(this->b), lend(this->beta),
+                    lend(this->x));
+    this->op->apply(lend(this->alpha), lend(this->b), lend(this->beta),
+                    lend(reference_result));
+
+    ASSERT_EQ(reference_result->get_value(), this->x->get_value());
+}
+
+
+TYPED_TEST(EnableValueTypedLinOp, CanCallExtendedApplyImplLinopValue)
+{
+    auto reference_result = gko::clone(this->x);
+
+    this->op->apply(lend(this->alpha), lend(this->b), this->beta_v,
+                    lend(this->x));
+    this->op->apply(lend(this->alpha), lend(this->b), lend(this->beta),
+                    lend(reference_result));
+
+    ASSERT_EQ(reference_result->get_value(), this->x->get_value());
+}
+
+
+}  // namespace

examples/adaptiveprecision-blockjacobi/adaptiveprecision-blockjacobi.cpp

@@ -100,10 +100,8 @@ int main(int argc, char *argv[])
     b->copy_from(host_x.get());
 
     // Calculate initial residual by overwriting b
-    auto one = gko::initialize<vec>({1.0}, exec);
-    auto neg_one = gko::initialize<vec>({-1.0}, exec);
-    auto initres = gko::initialize<real_vec>({0.0}, exec);
-    A->apply(lend(one), lend(x), lend(neg_one), lend(b));
+    auto initres = gko::initialize<real_vec>({0.0}, exec->get_master());
+    A->apply(1.0, lend(x), -1.0, lend(b));
     b->compute_norm2(lend(initres));
 
     // copy b again
@@ -146,8 +144,8 @@ int main(int argc, char *argv[])
     time += std::chrono::duration_cast<std::chrono::nanoseconds>(toc - tic);
 
     // Calculate residual
-    auto res = gko::initialize<real_vec>({0.0}, exec);
-    A->apply(lend(one), lend(x), lend(neg_one), lend(b));
+    auto res = gko::initialize<real_vec>({0.0}, exec->get_master());
+    A->apply(1.0, lend(x), -1.0, lend(b));
     b->compute_norm2(lend(res));
     auto impl_res = gko::as<real_vec>(logger->get_implicit_sq_resnorm());
 

examples/cb-gmres/cb-gmres.cpp

@@ -138,7 +138,7 @@ int main(int argc, char *argv[])
         b_host->at(i, 0) =
             ValueType{1} / std::sqrt(static_cast<ValueType>(A_size[0]));
     }
-    auto b_norm = gko::initialize<real_vec>({0.0}, exec);
+    auto b_norm = gko::initialize<real_vec>({0.0}, exec->get_master());
     b_host->compute_norm2(lend(b_norm));
     auto b = clone(exec, lend(b_host));
 
@@ -194,25 +194,20 @@ int main(int argc, char *argv[])
 
     // To measure if your solution has actually converged, the error of the
     // solution is measured.
-    // one, neg_one are objects that represent the numbers which allow for a
-    // uniform interface when computing on any device. To compute the residual,
-    // the (advanced) apply method is used.
-    auto one = gko::initialize<vec>({1.0}, exec);
-    auto neg_one = gko::initialize<vec>({-1.0}, exec);
-
-    auto res_norm_keep = gko::initialize<real_vec>({0.0}, exec);
-    auto res_norm_reduce = gko::initialize<real_vec>({0.0}, exec);
+    // To compute the residual, the (advanced) apply method is used.
+    auto res_norm_keep = gko::initialize<real_vec>({0.0}, exec->get_master());
+    auto res_norm_reduce = gko::initialize<real_vec>({0.0}, exec->get_master());
     auto tmp = gko::clone(gko::lend(b));
 
     // tmp = Ax - tmp
-    A->apply(lend(one), lend(x_keep), lend(neg_one), lend(tmp));
+    A->apply(1.0, lend(x_keep), -1.0, lend(tmp));
     tmp->compute_norm2(lend(res_norm_keep));
 
     std::cout << "\nResidual norm without compression:\n";
     write(std::cout, lend(res_norm_keep));
 
     tmp->copy_from(lend(b));
-    A->apply(lend(one), lend(x_reduce), lend(neg_one), lend(tmp));
+    A->apply(1.0, lend(x_reduce), -1.0, lend(tmp));
     tmp->compute_norm2(lend(res_norm_reduce));
 
     std::cout << "\nResidual norm with compression:\n";

examples/custom-logger/custom-logger.cpp

@@ -343,18 +343,14 @@ int main(int argc, char *argv[])
     logger->write();
 
     // To measure if your solution has actually converged, you can measure
-    // the error of the solution. one, neg_one are objects that represent
-    // the numbers which allow for a uniform interface when computing on any
-    // device. To compute the residual, all you need to do is call the apply
-    // method, which in this case is an spmv and equivalent to the LAPACK
-    // z_spmv routine. Finally, you compute the euclidean 2-norm with the
-    // compute_norm2 function.
-    auto one = gko::initialize<vec>({1.0}, exec);
-    auto neg_one = gko::initialize<vec>({-1.0}, exec);
-    auto res = gko::initialize<real_vec>({0.0}, exec);
-    A->apply(gko::lend(one), gko::lend(x), gko::lend(neg_one), gko::lend(b));
-    b->compute_norm2(gko::lend(res));
+    // the error of the solution. To compute the residual, all you need to do is
+    // call the apply method, which in this case is an spmv and equivalent to
+    // the LAPACK z_spmv routine. Finally, you compute the euclidean 2-norm with
+    // the compute_norm2 function.
+    auto res = gko::initialize<real_vec>({0.0}, exec->get_master());
+    A->apply(1.0, gko::lend(x), -1.0, gko::lend(b));
+    b->compute_norm2(lend(res));
 
     std::cout << "Residual norm sqrt(r^T r):\n";
-    write(std::cout, gko::lend(res));
+    write(std::cout, lend(res));
 }

examples/custom-stopping-criterion/custom-stopping-criterion.cpp

@@ -127,10 +127,8 @@ void run_solver(volatile bool *stop_iteration_process,
     write(std::cout, lend(x));
 
     // Calculate residual
-    auto one = gko::initialize<vec>({1.0}, exec);
-    auto neg_one = gko::initialize<vec>({-1.0}, exec);
-    auto res = gko::initialize<real_vec>({0.0}, exec);
-    A->apply(lend(one), lend(x), lend(neg_one), lend(b));
+    auto res = gko::initialize<real_vec>({0.0}, exec->get_master());
+    A->apply(1.0, lend(x), -1.0, lend(b));
     b->compute_norm2(lend(res));
 
     std::cout << "Residual norm sqrt(r^T r): \n";

examples/heat-equation/heat-equation.cpp

@@ -158,7 +158,7 @@ int main(int argc, char *argv[])
     // create output vector with initial guess for
     auto out_vector = in_vector->clone();
     // create scalar for source update
-    auto tau_source_scalar = gko::initialize<vec>({source_scale * tau}, exec);
+    auto tau_source_scalar = source_scale * tau;
     // create stencil matrix as shared_ptr for solver
     auto stencil_matrix = gko::share(mtx::create(exec));
     // assemble matrix
@@ -212,7 +212,7 @@ int main(int argc, char *argv[])
                     ->get_const_values());
         }
         // add heat source contribution
-        in_vector->add_scaled(gko::lend(tau_source_scalar), gko::lend(source));
+        in_vector->add_scaled(tau_source_scalar, gko::lend(source));
         // execute Euler step
         solver->apply(gko::lend(in_vector), gko::lend(out_vector));
         // swap input and output

examples/ilu-preconditioned-solver/ilu-preconditioned-solver.cpp

@@ -134,12 +134,10 @@ int main(int argc, char *argv[])
     write(std::cout, gko::lend(x));
 
     // Calculate residual
-    auto one = gko::initialize<vec>({1.0}, exec);
-    auto neg_one = gko::initialize<vec>({-1.0}, exec);
-    auto res = gko::initialize<real_vec>({0.0}, exec);
-    A->apply(gko::lend(one), gko::lend(x), gko::lend(neg_one), gko::lend(b));
-    b->compute_norm2(gko::lend(res));
+    auto res = gko::initialize<real_vec>({0.0}, exec->get_master());
+    A->apply(1.0, gko::lend(x), -1.0, gko::lend(b));
+    b->compute_norm2(lend(res));
 
     std::cout << "Residual norm sqrt(r^T r):\n";
-    write(std::cout, gko::lend(res));
+    write(std::cout, lend(res));
 }