saxpy.hpp

Go to the documentation of this file.

// -------------------------------------------------------//
//
// SHAMROCK code for hydrodynamics
// Copyright (c) 2021-2026 Timothée David--Cléris <tim.shamrock@proton.me>
// SPDX-License-Identifier: CeCILL Free Software License Agreement v2.1
// Shamrock is licensed under the CeCILL 2.1 License, see LICENSE for more information
//
// -------------------------------------------------------//
 
#pragma once

 
#include "shambase/assert.hpp"
#include "shambase/time.hpp"
#include "shambackends/DeviceBuffer.hpp"
#include "shambackends/DeviceScheduler.hpp"
#include "shambackends/math.hpp"
 
namespace sham::benchmarks {

    template<class T>
    inline void saxpy(u32 i, int n, T a, T *__restrict x, T *__restrict y) {
        if (i < n)
            y[i] = a * x[i] + y[i];
    }

    struct saxpy_result {
        std::string func_name;
        f64 seconds;
        f64 bandwidth;
        u64 byte_used;
    };

    template<class T>
    inline saxpy_result saxpy_bench(
        DeviceScheduler_ptr sched,
        int N,
        T init_x,
        T init_y,
        T a,
        int load_size,
        bool check_correctness) {
 
        sham::DeviceQueue &q = sched->get_queue();
 
        double seconds = shambase::get_max<double>();
 
        sham::DeviceBuffer<T> x{size_t(N), sched};
        sham::DeviceBuffer<T> y{size_t(N), sched};
 
        x.fill(init_x);
        y.fill(init_y);
 
        if (x.get_size() < N) {
            shambase::make_except_with_loc<std::runtime_error>(shambase::format(
                "x.get_size() < N\n  x.get_size() = {},\n  N = {}", x.get_size(), N));
        }
 
        if (y.get_size() < N) {
            shambase::make_except_with_loc<std::runtime_error>(shambase::format(
                "y.get_size() < N\n  y.get_size() = {},\n  N = {}", y.get_size(), N));
        }
 
        std::vector<T> y_res = {};
 
        for (int i = 0; i < 5; i++) {
 
            sham::EventList depends_list;
 
            auto x_ptr = x.get_write_access(depends_list);
            auto y_ptr = y.get_write_access(depends_list);
 
            depends_list.wait();
 
            sham::EventList empty_list{};
 
            shambase::Timer t;
            t.start();
            auto e = q.submit(empty_list, [&](sycl::handler &cgh) {
                cgh.parallel_for(sycl::range<1>{size_t(N)}, [=](sycl::item<1> item) {
                    // printf("%d\n", item.get_linear_id());
                    saxpy(item.get_linear_id(), N, a, x_ptr, y_ptr);
                });
            });
            e.wait();
            t.stop();
 
            x.complete_event_state(sycl::event{});
            y.complete_event_state(sycl::event{});
 
            seconds = sham::min(seconds, t.elapsed_sec());
 
            if (i == 0) {
                y_res = y.copy_to_stdvec();
            }
        }
 
        T expected = a * init_x + init_y;
 
        if (check_correctness) {
            T maxError = {};
            for (int i = 0; i < N; i++) {
                T delt = y_res[i] - expected;
 
                if constexpr (std::is_same_v<T, sycl::marray<float, 3>>) {
                    maxError[0] = sham::max(maxError[0], sham::abs(delt[0]));
                    maxError[1] = sham::max(maxError[1], sham::abs(delt[1]));
                    maxError[2] = sham::max(maxError[2], sham::abs(delt[2]));
                } else if constexpr (std::is_same_v<T, sycl::marray<float, 4>>) {
                    maxError[0] = sham::max(maxError[0], sham::abs(delt[0]));
                    maxError[1] = sham::max(maxError[1], sham::abs(delt[1]));
                    maxError[2] = sham::max(maxError[2], sham::abs(delt[2]));
                    maxError[3] = sham::max(maxError[3], sham::abs(delt[3]));
                } else {
                    maxError = sham::max(maxError, sham::abs(delt));
                }
            }
 
            SHAM_ASSERT(sham::equals(maxError, T{}));
        }
 
        return {
            .func_name = SourceLocation{}.loc.function_name(),
            .seconds   = seconds,
            .bandwidth = double(N) * load_size * 3 / seconds / 1e9,
            .byte_used = u64(N) * u64(load_size) * 2_u64};
    }
 
} // namespace sham::benchmarks