# Block diagonal and sparse matrices

**URL:** <https://discourse.julialang.org/t/block-diagonal-and-sparse-matrices/105551>\
**Category:** Numerics\
**Tags:** linearalgebra\
**Created:** [October 29, 2023, 11:13pm UTC](https://discourse.julialang.org/t/block-diagonal-and-sparse-matrices/105551 "2023-10-29T23:13:12Z")\
**Posts on this page:** 1\
**Showing post:** 14

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [October 30, 2023, 5:07pm UTC](https://discourse.julialang.org/t/block-diagonal-and-sparse-matrices/105551/14 "2023-10-30T17:07:23Z")

</div>

> [@loisel](#):
>
> A while back, I was writing a preconditioner for PetSc and I was using the various matrix primitives built into PetSc. I think I multiplied some sparse matrices together. I knew that theoretically this could be done efficiently, and I could have written special-purpose code to make sure it was done efficiently, but I naively assumed PetSc’s primitive would be clever enough and unfortunately it didn’t work out. That preconditioner ended up being slow.

Side note: If you have a preconditioner ABC that is a product of sparse matrices, normally I would express this as a linear operator x \mapsto A(B(Cx)) to avoid the risk of losing sparsity in the product. (I think in Petsc you can do this as a “composite” preconditioner, or by defining a custom “shell” preconditioner with your own matvec. In Julia, you can implement such “lazy” products easily with [LinearMaps.jl](https://github.com/JuliaLinearAlgebra/LinearMaps.jl) or [LinearOperators.jl](https://github.com/JuliaSmoothOptimizers/LinearOperators.jl).)

---

_[View the full topic](https://discourse.julialang.org/t/block-diagonal-and-sparse-matrices/105551)._