# Example of direct linear solve, Sparse Matrix in cuda

**URL:** <https://discourse.julialang.org/t/example-of-direct-linear-solve-sparse-matrix-in-cuda/22000>\
**Category:** General Usage\
**Tags:** gpu, linearalgebra\
**Created:** [March 18, 2019, 1:50pm UTC](https://discourse.julialang.org/t/example-of-direct-linear-solve-sparse-matrix-in-cuda/22000 "2019-03-18T13:50:20Z")\
**Posts on this page:** 3\
**Page:** 1

<div class="post-metadata">

**Author:** ![rveltz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rveltz/32/2707_2.png) [@rveltz](https://discourse.julialang.org/u/rveltz)\
**Post date:** [March 18, 2019, 1:50pm UTC](https://discourse.julialang.org/t/example-of-direct-linear-solve-sparse-matrix-in-cuda/22000/1 "2019-03-18T13:50:20Z")

</div>

Hi,

In hope to write a **simple** example for my package [PseudoArcLengthContinuation.jl](https://github.com/rveltz/PseudoArcLengthContinuation.jl), I tried to look for examples of direct sparse linear solves **on the GPU** in the context of PDE.

So far I have not found an example where a direct method for solving `A * x = rhs` (with `A` being the discretisation of a differential operator) outperforms in term of computation time its CPU equivalent. I used the `CUSOLVER` part of `CuArrays`. I noticed that `cuda 10.0` has much better performance with regard to this problem. I tried, among other things, the matrix `L1` in this [example](https://rveltz.github.io/PseudoArcLengthContinuation.jl/dev/#Example-2:-Snaking-with-2d-Swift-Hohenberg-equation-1).

Hence, my question is “Does anybody have an example of fast **direct** sparse linear solve on the GPU in the context of PDE?”

I guess I would also be happy with an iterative method simple enough to avoid the use of preconditionning.

Thank you a lot for your help,

Best

---

<div class="post-metadata">

**Author:** ![mohamed82008](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mohamed82008/32/18171_2.png) [@mohamed82008](https://discourse.julialang.org/u/mohamed82008)\
**Post date:** [March 18, 2019, 2:38pm UTC](https://discourse.julialang.org/t/example-of-direct-linear-solve-sparse-matrix-in-cuda/22000/2 "2019-03-18T14:38:32Z")

</div>

If your matrix is SPD, `IterativeSolvers.cg` supports arbitrary array types. You can use `CuVector` for the rhs and `CuArrays.CUSPARSE.CuSparseMatrixCSC` for the matrix or the csr equivalent.

---

<div class="post-metadata">

**Author:** ![rveltz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rveltz/32/2707_2.png) [@rveltz](https://discourse.julialang.org/u/rveltz)\
**Post date:** [March 18, 2019, 3:17pm UTC](https://discourse.julialang.org/t/example-of-direct-linear-solve-sparse-matrix-in-cuda/22000/3 "2019-03-18T15:17:47Z")

</div>

Thank you! I looked at the test example for `cg` and for `A = laplace_matrix(Float64, 1000, 2)` there is no CV convergence in `1000` steps. Also, the preconditionner is not improving this…