# Parallel sparse solve and inverse

**URL:** <https://discourse.julialang.org/t/parallel-sparse-solve-and-inverse/81538>\
**Category:** Numerics\
**Created:** [May 23, 2022, 8:32pm UTC](https://discourse.julialang.org/t/parallel-sparse-solve-and-inverse/81538 "2022-05-23T20:32:41Z")\
**Posts on this page:** 1\
**Page:** 1

<div class="post-metadata">

**Author:** ![sobhan](https://avatars.discourse-cdn.com/v4/letter/s/258eb7/32.png) [@sobhan](https://discourse.julialang.org/u/sobhan)\
**Post date:** [May 23, 2022, 8:32pm UTC](https://discourse.julialang.org/t/parallel-sparse-solve-and-inverse/81538/1 "2022-05-23T20:32:42Z")

</div>

So i need to solve a really big (and sparse) system of equation and to my knowledge there are no parallel sparse solvers for julia. I tried the pardiso for MKL but that wasn’t parallel and had performance similar to the one used by julia. What i do right now is to do an lu factorization of my lhs (lu is the factorization chosen by factorize) and then use something like

```julia
"solve Ax = b"
function par_solv!(A, b, x=similar(b))
    blocks = Threads.nthreads()
    block_size = ceil(Int, size(x, 2) / blocks)
    @threads for i in 1:blocks
            s = (i-1)*block_size+1:min(i * block_size, size(x, 2))
            ldiv!(@view(x[:, s]), A, @view(b[:, s]))
    end
    return x
end

```

This works pretty well, but can i do better? Sometimes i need the inverse of lhs and i give a dense identity matrix to this function. Can i do better?

I only with work with SparseMatrixCSC and standard arrays. i don’t care, nor plan to support any other format.
