oneMKL: transposed sparse gemv reads y when beta is zero
Maintainer thường phản hồi trong vòng 1 ngày
Chưa có ai nhận issue này.
Đánh giá
- Độ khó
- 4/5
- Thời gian dự kiến
- 1-2 ngày
- Mức phù hợp với người mới
- 22/100
Hướng nghiên cứu
The reproducer is the standalone repro.cpp in the issue, built with icpx -fsycl -qmkl=sequential and run on the GPU device; the calls under suspicion are sparse::gemv with trans and conjtrans and beta set to zero. The wrapper calls live in deps/src/onemkl.cpp, which is where the beta-zero workaround from #702 would sit. Done means the upstream oneMKL report is filed on the Intel forum and linked from this issue, since the fix belongs in oneMKL and not in oneAPI.jl.
Do mô hình lập chỉ mục viết ra từ nội dung của issue.
Mô tả
On the GPU, oneapi::mkl::sparse::gemv with a CSR handle and transpose::trans or transpose::conjtrans reads y when beta is zero, so NaNs in an uninitialized output buffer end up in the result. The same call with transpose::nontrans ignores y, as do CSR gemm with any op and the CPU device. Prefilling y with 42 gives the correct result, so y is scaled by beta rather than skipped. BLAS convention is that y need not be set on input when beta is zero, and Julia's mul! relies on that: * passes uninitialized memory from similar.
In oneAPI.jl this made oneSparseMatrixCSC * x return NaNs depending on what the allocator handed out, because CSC matrices are stored as the CSR handle of their transpose and multiplied with a transposed op. That is the intermittent oneMKL.oneSparseMatrixCSC test failure in GPUArrays.jl's downstream CI. #702 works around it by zeroing the output before oneMKL calls that take a beta when beta is zero.
Reproducer
The calls mirror what our support library does (deps/src/onemkl.cpp): init_matrix_handle, set_csr_data with one-based indices, then gemv/gemm with USM pointers and no optimize_* call.
repro.cpp
// oneMKL sparse::gemv on a CSR handle with op = trans/conjtrans reads y when beta == 0
#include <sycl/sycl.hpp>
#include <oneapi/mkl.hpp>
#include <mkl.h>
#include <algorithm>
#include <cmath>
#include <cstdint>
#include <cstdio>
#include <utility>
using oneapi::mkl::transpose;
namespace sparse = oneapi::mkl::sparse;
int main() {
sycl::queue q{sycl::default_selector_v, sycl::property::queue::in_order()};
char ver[256];
mkl_get_version_string(ver, sizeof(ver));
auto dev = q.get_device();
printf("%s\n%s (driver %s)\n\n", ver, dev.get_info<sycl::info::device::name>().c_str(),
dev.get_info<sycl::info::device::driver_version>().c_str());
// 4x3 matrix in one-based CSR; row 2 and column 2 are empty
// [1 0 4]
// [0 0 0]
// [0 0 2]
// [3 0 0]
const int64_t m = 4, n = 3, nnz = 4;
int32_t *row_ptr = sycl::malloc_shared<int32_t>(m + 1, q);
int32_t *col_ind = sycl::malloc_shared<int32_t>(nnz, q);
float *val = sycl::malloc_shared<float>(nnz, q);
const int32_t rp[] = {1, 3, 3, 4, 5}, ci[] = {1, 3, 3, 1};
const float v[] = {1, 4, 2, 3};
std::copy(rp, rp + m + 1, row_ptr);
std::copy(ci, ci + nnz, col_ind);
std::copy(v, v + nnz, val);
sparse::matrix_handle_t A;
sparse::init_matrix_handle(&A);
sparse::set_csr_data(q, A, m, n, nnz, oneapi::mkl::index_base::one, row_ptr, col_ind, val).wait();
// y = 1 * op(A) * x + 0 * y with x = ones and y prefilled; gemm uses 2 columns
const float expect_N[] = {5, 0, 2, 3}, expect_T[] = {4, 0, 6};
float *x = sycl::malloc_shared<float>(2 * m, q), *y = sycl::malloc_shared<float>(2 * m, q);
const std::pair<transpose, const char *> ops[] = {
{transpose::nontrans, "nontrans"}, {transpose::trans, "trans"}, {transpose::conjtrans, "conjtrans"}};
for (int64_t ncols : {1, 2})
for (auto [op, opname] : ops)
for (float fill : {NAN, 42.0f}) {
const bool N = op == transpose::nontrans;
const int64_t nx = N ? n : m, ny = N ? m : n;
std::fill(x, x + nx * ncols, 1.0f);
std::fill(y, y + ny * ncols, fill);
if (ncols == 1)
sparse::gemv(q, op, 1.0f, A, x, 0.0f, y).wait();
else
sparse::gemm(q, oneapi::mkl::layout::col_major, op, transpose::nontrans, 1.0f, A,
x, ncols, nx, 0.0f, y, ny).wait();
bool ok = true;
printf("%s %-9s y_in=%-3g y =", ncols == 1 ? "gemv" : "gemm", opname, fill);
for (int64_t i = 0; i < ny * ncols; i++) {
printf(" %g", y[i]);
ok &= y[i] == (N ? expect_N : expect_T)[i % ny];
}
printf(" %s\n", ok ? "OK" : "WRONG");
}
sparse::release_matrix_handle(q, &A).wait();
for (void *p : {(void *)row_ptr, (void *)col_ind, (void *)val, (void *)x, (void *)y})
sycl::free(p, q);
}
icpx -fsycl -qmkl=sequential -O2 repro.cpp -o repro
ONEAPI_DEVICE_SELECTOR=level_zero:gpu ./repro
Intel(R) oneAPI Math Kernel Library Version 2026.1-Product Build 20260612 for Intel(R) 64 architecture applications
Intel(R) Iris(R) Xe Graphics (driver 1.15.38308)
gemv nontrans y_in=nan y = 5 0 2 3 OK
gemv nontrans y_in=42 y = 5 0 2 3 OK
gemv trans y_in=nan y = nan nan nan WRONG
gemv trans y_in=42 y = 4 0 6 OK
gemv conjtrans y_in=nan y = nan nan nan WRONG
gemv conjtrans y_in=42 y = 4 0 6 OK
gemm nontrans y_in=nan y = 5 0 2 3 5 0 2 3 OK
gemm nontrans y_in=42 y = 5 0 2 3 5 0 2 3 OK
gemm trans y_in=nan y = 4 0 6 4 0 6 OK
gemm trans y_in=42 y = 4 0 6 4 0 6 OK
gemm conjtrans y_in=nan y = 4 0 6 4 0 6 OK
gemm conjtrans y_in=42 y = 4 0 6 4 0 6 OK
With ONEAPI_DEVICE_SELECTOR=opencl:cpu every case is OK.
Variations:
- Every output entry becomes NaN, not just the one for the empty column, and a matrix without empty rows or columns fails the same way.
- Zero-based indices fail the same way.
- Calling
sparse::optimize_gemv(q, op, A)before each gemv makes all cases correct. Calling it only fornontransand then running atransgemv still fails.
On the older oneMKL of the Aurora LTS stack, sparse gemm with beta = 0 also propagated NaN, for every op (first CI run of #702). I haven't been able to reproduce that on 2026.1.
Versions
- oneMKL 2026.1 (Product Build 20260612), the release
oneAPI_Support_jll0.11.0 is built from - Intel oneAPI DPC++/C++ Compiler 2026.1.0 (2026.1.0.20260617)
- Intel Iris Xe Graphics (TigerLake-LP GT2), Level Zero driver 1.15.38308 (NEO 26.18.38308 from
NEO_jll), Level Zero loader 1.29.0 - Debian 13
Float64 is untested: this GPU has no FP64 support.
To do: report this to Intel on the oneMKL community forum, and link the report here.
- Ngôn ngữ chính
- Julia
- Star
- 217
- Fork
- 37
- Merge trung bình
- 1 ngày 5 giờ
- Pull request đã merge (30 ngày)
- 40
Chuẩn bị môi trường
Dự án này không cung cấp dev container, Dockerfile hay hướng dẫn đóng góp, nên bạn cần tự thiết lập môi trường: hãy bắt đầu từ README và xem hướng dẫn đóng góp lần đầu của chúng tôi để biết các bước chung.
Bắt đầu từ đâu
- Đọc hết issue, rồi đọc hướng dẫn đóng góp của dự án.
- Bình luận trên issue rằng bạn sẽ nhận — tránh hai người làm cùng một việc.
- Fork repository và làm thay đổi trên một nhánh.
- Mở pull request có tham chiếu số hiệu của issue.
Issue khác của JuliaGPU/oneAPI.jl
-
Độ khó 3/5 1-2 ngày Mức phù hợp với người mới 62/100
Maintainer thường phản hồi trong vòng 1 ngày
-
Make blocking Level Zero calls GC-safeCó thể đã có người làm @maleadt đã nhận 1 ngày trước. Đang mở
Độ khó 3/5 1-2 ngày Mức phù hợp với người mới 55/100
Maintainer thường phản hồi trong vòng 1 ngày
-
Array construction queries device properties on every allocationCó thể đã có người làm @maleadt đã nhận 1 ngày trước. Đang mở
Độ khó 3/5 Nửa ngày Mức phù hợp với người mới 65/100
Maintainer thường phản hồi trong vòng 1 ngày
-
Độ khó 5/5 Hơn một tuần Mức phù hợp với người mới 35/100
JuliaGPU/oneAPI.jl#663 · 1 reaction ·
Maintainer thường phản hồi trong vòng 1 ngày
-
Độ khó 4/5 3-5 ngày Mức phù hợp với người mới 35/100
JuliaGPU/oneAPI.jl#607 · 2 bình luận ·
Maintainer thường phản hồi trong vòng 1 ngày
Tất cả issue của JuliaGPU/oneAPI.jl
Issue tương tự
-
Add more snippets to R and JuliaĐang mởenhancement good first issue help wanted
Độ khó 2/5 Nửa ngày Mức phù hợp với người mới 82/100
-
authority-vechtstromen ready-for-agent
Độ khó 1/5 Dưới một giờ Mức phù hợp với người mới 75/100
Deltares/Ribasim-NL#890 ·
Maintainer thường phản hồi trong vòng 1 ngày
-
Độ khó 1/5 Dưới một giờ Mức phù hợp với người mới 88/100
Maintainer thường phản hồi trong vòng 1 ngày
-
Precompilation fails on Julia 1.14 nightly: HessenbergQ lmul! methods overwrite LinearAlgebra'sĐang mở
Độ khó 2/5 1-3 giờ Mức phù hợp với người mới 78/100
JuliaLinearAlgebra/GenericLinearAlgebra.jl#183 · 1 bình luận ·
-
Độ khó 2/5 1-3 giờ Mức phù hợp với người mới 72/100
JuliaGPU/KernelAbstractions.jl#897 ·
Maintainer thường phản hồi trong vòng 1 ngày