Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Input inconsistent with MKL Pardiso and example in README.md #73

Open
Alexander-Barth opened this issue Feb 23, 2021 · 7 comments
Open

Comments

@Alexander-Barth
Copy link

I just tried to run the example in the README.md with MKLPardisoSolver.
I use MKLSparse v1.1.0 and Pardiso v0.5.1 in Julia 1.5.1. MKL is installed automatically by MKLSparse.
Does somebody have any idea what could be the issue? The input matrix sizes seem to be correct.

Thank you for this nice package!

Example script:

using MKLSparse                                                                                                                                                       
using Pardiso                                                                                                                                                         
using SparseArrays                                                                                                                                                    
using LinearAlgebra                                                                                                                                                   
                                                                                                                                                                      
ps = MKLPardisoSolver()                                                                                                                                               
                                                                                                                                                                      
A = sparse(rand(10, 10))                                                                                                                                              
B = rand(10, 2)                                                                                                                                                       
X = zeros(10, 2)                                                                                                                                                      
solve!(ps, X, A, B)                                                                                                                                                   
                                                                                                                                                                      
#same error with                                                                                                                                                      
#B = rand(10)                                                                                                                                                         
#X = solve(ps, A, B)                                                                                                                                                  

Error message:

julia> solve!(ps, X, A, B)                                                                                                                                            
ERROR: Input inconsistent.                                                                                                                                            
Stacktrace:                                                                                                                                                           
 [1] check_error at /home/abarth/.julia/packages/Pardiso/yZsYO/src/mkl_pardiso.jl:76 [inlined]                                                                        
 [2] ccall_pardiso(::MKLPardisoSolver, ::Int64, ::Array{Float64,1}, ::Array{Int64,1}, ::Array{Int64,1}, ::Int64, ::Array{Float64,2}, ::Array{Float64,2}) at /home/aba\
rth/.julia/packages/Pardiso/yZsYO/src/mkl_pardiso.jl:71                                                                                                               
 [3] pardiso(::MKLPardisoSolver, ::Array{Float64,2}, ::SparseMatrixCSC{Float64,Int64}, ::Array{Float64,2}) at /home/abarth/.julia/packages/Pardiso/yZsYO/src/Pardiso.\
jl:337                                                                                                                                                                
 [4] solve!(::MKLPardisoSolver, ::Array{Float64,2}, ::SparseMatrixCSC{Float64,Int64}, ::Array{Float64,2}, ::Symbol) at /home/abarth/.julia/packages/Pardiso/yZsYO/src\
/Pardiso.jl:273                                                                                                                                                       
 [5] solve!(::MKLPardisoSolver, ::Array{Float64,2}, ::SparseMatrixCSC{Float64,Int64}, ::Array{Float64,2}) at /home/abarth/.julia/packages/Pardiso/yZsYO/src/Pardiso.j\
l:238                                                                                                                                                                 
 [6] top-level scope at REPL[71]:1             

Environment:

julia> versioninfo()                                                                                                                                                  
Julia Version 1.5.1                                                                                                                                                   
Commit 697e782ab8 (2020-08-25 20:08 UTC)                                                                                                                              
Platform Info:                                                                                                                                                        
  OS: Linux (x86_64-pc-linux-gnu)                                                                                                                                     
  CPU: Intel(R) Core(TM) i7-7700 CPU @ 3.60GHz                                                                                                                        
  WORD_SIZE: 64                                                                                                                                                       
  LIBM: libopenlibm                                                                                                                                                   
  LLVM: libLLVM-9.0.1 (ORCJIT, skylake)              
@KristofferC
Copy link
Member

KristofferC commented Feb 23, 2021

It seems loading MKLSparse breaks things. My guess is that https://github.com/JuliaSparse/MKLSparse.jl/blob/6060cd08a2586b275c95aac650c379ff48130ef7/src/MKLSparse.jl#L6 sets MKL in 64 bit mode but Pardiso.jl assumes that integers should be passed to MKL as 32 bits (unless

if LinearAlgebra.BLAS.vendor() === :mkl && LinearAlgebra.BlasInt == Int64
const MklInt = Int64
is true, but this is not the case here).

@Alexander-Barth
Copy link
Author

Alexander-Barth commented Feb 23, 2021

Perfect! Indeed without loading MKLSparse the code works! Thank you very much!

@KristofferC
Copy link
Member

Well, not perfect :P. We should figure out a way to make the two packages collaborate.

@Alexander-Barth
Copy link
Author

It is indeed surprising that pardiso still uses 32-bit integers.

@Alexander-Barth
Copy link
Author

For your information, on my system LinearAlgebra.BLAS.vendor() is :openblas64 and LinearAlgebra.BlasInt is Int64.
I am wondering if we can use the same logic than in MKLSparse to determine MklInt.

const MklInt = (Base.USE_BLAS64 ? Int64 : Int32)     
# or simply?
const MklInt =  LinearAlgebra.BlasInt                                                                                                             

@KristofferC
Copy link
Member

The problem is that SparseMatrices in Julia has Int64 indices in almost all cases. So the code in MKLSparse doesn't really make sense (LinearAlgebra.BlasInt) is pretty much irrelevant to sparse matrices.

It's just kinda hard how to do this properly because the MKL Interface state is global and packages use it like it is their own local setting to play with.

@lucasbanting
Copy link

Any work around for using both Pardiso and MKLSparse? Using SparseMatrices with Ti=Int32 only?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants