Storage-level mul! gives wrong results or aborts for SubArray operands #700

The storage-level mul! methods in lib/mkl/linalg.jl pass any strided SubArray of a oneArray to oneMKL, and the wrappers in lib/mkl/wrappers_blas.jl use lda = max(1, stride(A, 2)) and incx = stride(x, 1) without checking the layout. Depending on the view, that gives silently wrong results or aborts the process:

using oneAPI, LinearAlgebra
A = oneArray(rand(Float32, 8, 4)); B = oneArray(rand(Float32, 4, 5)); x = oneArray(rand(Float32, 4))

# row-stepping view: wrong result, for gemm and gemv alike
C = mul!(oneArray(zeros(Float32, 4, 5)), 'N', 'N', view(A, 1:2:8, :), B, true, false)
Array(C) ≈ Array(A)[1:2:8, :] * Array(B)   # false
y = mul!(oneArray(zeros(Float32, 4)), 'N', view(A, 1:2:8, :), x, true, false)
Array(y) ≈ Array(A)[1:2:8, :] * Array(x)   # false

# reversed vector: wrong result (BLAS reads a negative increment from the end of the vector)
y = mul!(oneArray(zeros(Float32, 8)), 'N', A, view(x, 4:-1:1), true, false)
Array(y) ≈ Array(A) * reverse(Array(x))    # false

# reversed columns: the process aborts
mul!(oneArray(zeros(Float32, 8, 5)), 'N', 'N', view(A, :, 4:-1:1), B, true, false)
# terminate called after throwing an instance of 'oneapi::mkl::invalid_argument'
#   what():  oneapi::mkl::oneapi::mkl::blas::sgemm: invalid argument: Illegal value supplied for parameter lda
# signal 6 (-6): Aborted

The abort also means onemklSgemm in deps/src/onemkl.cpp lets oneMKL exceptions escape through the C interface instead of turning them into an error code.

The public mul!(C, A, B) isn't affected, because GPUArrays intercepts products with views before they reach these methods. JuliaGPU/GPUArrays.jl#800 starts forwarding matrix-matrix products with views here again, but only N/T/C products where every matrix has a unit row stride and stride(X, 2) >= size(X, 1). Calling the storage-level mul! directly is public API since Julia 1.13 (JuliaLang/LinearAlgebra.jl#1671), though. CUDA.jl added the equivalent guards in JuliaGPU/CUDA.jl#3316:

• only call oneMKL when every matrix has stride(X, 1) == 1 and stride(X, 2) >= max(1, size(X, 1)), and every vector has a positive stride;
• otherwise fall back to GPUArrays.generic_matmatmul!.

That applies to gemm, gemv, symm, hemm, symv and hemv.

Reproduced with oneAPI.jl 2.9.2 on an Iris Xe, Julia 1.12.7. The mul! methods and the gemm/gemv/symm wrappers are unchanged on main.