Use 64-bit Pardiso when BLAS vendor is LBT #121
Conversation
Attempt at fixing JuliaSparse#110. Note I don't really understand the internals, just trying this out and hopefully the CI tests will let us know if that was not worth it.
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## master #121 +/- ##
==========================================
- Coverage 45.74% 45.61% -0.14%
==========================================
Files 4 4
Lines 470 467 -3
==========================================
- Hits 215 213 -2
+ Misses 255 254 -1 ☔ View full report in Codecov by Sentry. 🚀 New features to boost your workflow:
|
|
Not sure this was very useful, but at least it's not using the deprecated How would I go about adding tests with |
|
Can this choice really be made at compile time now? Since with LBT you change things dynamically, there might need to be a switch at runtime to choose which ccall to make? |
|
FWIW I tried enforcing |
|
@KristofferC could you clarify how to make this switch at runtime? I'm really just looking for a MWE that uses Pardiso.jl with |
|
I mean, when we want to do a pardiso call we have to check with LBT if we are using 64 or 32 bit MKL and decide what to do at that point. |
|
You mean the isMkl() = if VERSION ≥ v"1.7.0-beta2"
config = BLAS.get_config().loaded_libs
occursin("mkl_rt", config[1].libname)
else
BLAS.vendor() === :mkl
end
isMkl64() = isMkl() && LinearAlgebra.BlasInt == Int64
MklInt() = isMkl64() ? Int64 : Int32
pardisoFunc() = isMkl64() ? :pardiso_64 : :pardisoAgain this is not my expertise so I may be misunderstanding what you mean. |
|
Yes, but all the ccalls need to be split up into two, for example: Lines 51 to 53 in 650f566 need to be something like: etc. |
|
But we are not going through LBT here though so we need to be consistent with the MKL that is loaded... |
|
After thinking about it more, I think what of what I have said here is more or less nonsense. Updating the check to use |
|
Do you have Julia code for what fails? We should look at what the actual ccall becomes and figure out what is wrong with it. |
|
Yes but I'll need a little time to cook a MWE, but essentially I commented out the checks here so my own dev'ed Pardiso.jl contains just const MklInt = Int64
const PARDISO_FUNC = :pardiso_64Then I loaded these packages (I used the LinearSolve.jl interface) using MKL
using LinearAlgebra
import Pardiso
using LinearSolveAnd finally I created a simple matrix prob = init(LinearProblem(A, b), MKLPardisoIterate())
sol = solve!(prob).u |
|
I made a branch of Pardiso.jl that forces MklInt64 there, and a little project that adds that branch and MKL.jl there. When running the simple example in the ReadMe, which is just: using Pkg
Pkg.activate(".")
Pkg.instantiate()
using MKL
using LinearAlgebra
using SparseArrays
# Load my branch of Pardiso that forces MklInt = Int64
# https://github.com/briochemc/Pardiso.jl/tree/MKL_force_Int64
using Pardiso
ps = MKLPardisoSolver()
A = sparse(rand(10, 10))
B = rand(10, 2)
X = zeros(10, 2)
solve!(ps, X, A, B)I checked that it uses MKLParido64"LinearAlgebra.BlasInt = Int64
Base.USE_BLAS64 = true
Pardiso.MklInt = Int64
Pardiso.PARDISO_FUNC = :pardiso_64But this segfaults (on my cluster) with: EDIT: So I'm likely doing something wrong of course but any help to fix this much appreciated, so that I can try with big matrices that actually requires |
2 participants
Attempt at fixing #110.
Note I don't really understand the internals, just trying this out and hopefully the CI tests will let us know if that was not worth it.