MATLAB already comes with Intel MKL for its BLAS implementation. There's no reason to replace that.
As for using GPUs, if you make your array a gpuArray
(to do that, just do gpuArray(A)
), then you can use MATLAB's matrix multiplication and it will use optimized kernals from MAGMA to perform the computation. You can Google around to reason some people saying this outperforms CUBLAS by like 10%, but the comments are usually old (2013) and blablabla: it's fast enough that it's likely the best option if you're in MATLAB (though if you really want performance, you should look at Julia with CUBLAS, which will have a lower interop overhead and faster user-compiled kernals).
The reason why you won't find a BLAS implementation that "has both" is because the implementations have to be completely different to fully leverage the GPU, and so at that point they might as well be different libraries since the reason to bundle code is usually for some form of code reuse. You can try swapping out backend BLAS implementations in MATLAB to learn, but it likely won't cause a performance change it may give difficulties because it's very undocumented. That's just one problem among many with closed-source software. If you want to be able to modify every detail, swap BLAS libraries and libm implementations, etc., you may want to look into using an open-source software instead, which as above I'd recommend Julia (or you may be able to do this with Octave, though I don't know and will refer to whatever documentation they have).