# Solving Sparse Linear Systems fast

**URL:** <https://discourse.julialang.org/t/solving-sparse-linear-systems-fast/83071>\
**Category:** Performance\
**Tags:** sparse, linearsolve\
**Created:** [June 20, 2022, 2:11pm UTC](https://discourse.julialang.org/t/solving-sparse-linear-systems-fast/83071 "2022-06-20T14:11:41Z")\
**Posts on this page:** 1\
**Showing post:** 9

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [June 21, 2022, 9:01am UTC](https://discourse.julialang.org/t/solving-sparse-linear-systems-fast/83071/9 "2022-06-21T09:01:22Z")

</div>

> [@pushkar\_khandare](#):
>
> Unfortunately when I applied for an academic license to use Pardiso I got an email stating that academic licences are not being given out starting January 2022.  
> Haven’t been able to use Pardiso.jl yet.

You don’t need it for MKLPardiso.

MWE with LinearSolve.jl:

```julia
using LinearSolve, LinearSolvePardiso, SparseArrays, Random

A = sparse([1.0 0 -2 3
    0 5 1 2
    -2 1 4 -7
    3 2 -7 5])
b = rand(4)
prob = LinearProblem(A, b)
u = solve(prob, MKLPardisoFactorize(); cache_kwargs...).u

```

should work and automatically install any binaries needed.

> [@pushkar\_khandare](#):
>
> This is for a smaller matrix. The first few and last few rows have diagonal entries only, and the rest come from the fact that I am solving a FEM problem on a square mesh.

This has a regular sparsity pattern, so it should be a case where UMFPACK-type algorithms do well because they should be able to exploit the parallelism of the BLAS3 kernel. I’d try:

```julia
using LinearSolve, LinearSolvePardiso, SparseArrays
A = # your matrix
b = # test b

prob = LinearProblem(A, b)

for alg in (
    MKLPardisoFactorize(),
    MKLPardisoIterate(),
    UMFPACKFactorization(),
    KLUFactorization())

    u = solve(prob, alg).u
end

using AlgebraicMultigrid
ml = ruge_stuben(A) # Construct a Ruge-Stuben solver
pl = aspreconditioner(ml)
solve(prob1, KrylovJL_GMRES(), Pl = pl).u

using IncompleteLU
pl = ilu(A, τ = 0.1) # τ needs to be tuned per problem
solve(prob1, KrylovJL_GMRES(), Pl = pl).u

```

I should make that into a tutorial. You can see this in action in [https://diffeq.sciml.ai/stable/tutorials/advanced\_ode\_example/](https://diffeq.sciml.ai/stable/tutorials/advanced_ode_example/)

---

_[View the full topic](https://discourse.julialang.org/t/solving-sparse-linear-systems-fast/83071)._