Pardiso 6 significantly slower than MKL Pardiso

[ma-sadeghi]

Hi,

I just tried solving the 2D Laplacian on a 2,000 by 2,000 grid using Pardiso 6 and MKL Pardiso and MKL Pardiso is almost 2X faster. Based on the Pardiso website, it should be the opposite by a significant margin (they have reported ~10X faster factorization using Pardiso 6).

Have you also experienced the same issue and is that normal?

Thanks, Amin

[KristofferC]

Can you post the code that generates the matrix and the code you use to factorize it. Would be good to have for eventual future investigation.

[ma-sadeghi]

There you go:

using Pardiso
using LinearAlgebra # for norm
using SparseArrays
using Random

function spdiagm_nonsquare(m, n, args...)
    I, J, V = SparseArrays.spdiagm_internal(args...)
    return sparse(I, J, V, m, n)
end

function ∇²(n₁,n₂)
    o₁ = ones(n₁)
    ∂₁ = spdiagm_nonsquare(n₁+1,n₁,-1=>-o₁,0=>o₁)
    o₂ = ones(n₂)
    ∂₂ = spdiagm_nonsquare(n₂+1,n₂,-1=>-o₂,0=>o₂)
    return kron(sparse(I,n₂,n₂), ∂₁'*∂₁) + kron(∂₂'*∂₂, sparse(I,n₁,n₁))
end

@info("Assembling A and b...")
A = ∇²(1000, 1000)
b = rand(A.n)

ps = MKLPardisoSolver()

@info("Solving...")
@time X = solve(ps, A, b)

[RoyiAvital]

On Pardiso's site they show graphs which suggests they are much much faster. Though they show results for Pardiso 7.

Is there a way to use the package with updated version of Pardiso?

[ma-sadeghi]

@RoyiAvital The site has recently been updated. They showed similar performance gains for Pardiso v6 too, but I tried using both Julia's wrapper (Pardiso.jl) and Python's wrapper (pypardiso), and the performance gains didn't realize.

Just to be clear, I'm not suggesting they "lie" on their website, but for sure the performance gains are for their specific problem. One possibility is that their problem involves partial refactorization that might not be available on MKL Pardiso (?), or some other explanation.

[j-fu]

For Pardiso 7 see #75 - in the moment (just checked) they still offer Pardiso 6 when it comes to download.