# Symmetric Indefinite Solvers (with block diagonal pivoting)

**URL:** <https://discourse.julialang.org/t/symmetric-indefinite-solvers-with-block-diagonal-pivoting/8541>\
**Category:** Numerics\
**Tags:** fortran, linearalgebra, compilation\
**Created:** [January 23, 2018, 5:00am UTC](https://discourse.julialang.org/t/symmetric-indefinite-solvers-with-block-diagonal-pivoting/8541 "2018-01-23T05:00:29Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![jacob-roth](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jacob-roth/32/1862_2.png) [@jacob-roth](https://discourse.julialang.org/u/jacob-roth)\
**Post date:** [January 23, 2018, 5:00am UTC](https://discourse.julialang.org/t/symmetric-indefinite-solvers-with-block-diagonal-pivoting/8541/1 "2018-01-23T05:00:29Z")

</div>

Are there any modules that provide access to (sparse) symmetric indefinite solvers with Bunch-Kaufman-Parlett pivoting where I can access the `L` and `D` matrices, or would it be useful to build a module that does this? Right now, I have to do something like [this](https://stackoverflow.com/questions/48494948/julia-converting-cholmod-factor-to-sparse-matrix-and-back-again/48506723?noredirect=1#comment84032952_48506723).

As of now, is compiling IPOPT with the HSL linear solvers (for example, HSL\_MA57) the easiest way to currently use the HSL routines? I saw a previous [post](http://julia-programming-language.2336112.n4.nabble.com/Using-Custom-Compiled-IPOPT-with-JuMP-td13128.html) build a custom IPOPT with these solvers, but if I only want to call the Bunch-Kaufman-Parlett `LDLt` factorization (MA57), can I access it through IPOPT, or is there a more direct way? Btw, as another alternative, Pardiso doesn’t yet allow you to access the `L` and `D` factors.

Thanks!

---

<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:** [February 8, 2018, 5:34am UTC](https://discourse.julialang.org/t/symmetric-indefinite-solvers-with-block-diagonal-pivoting/8541/2 "2018-02-08T05:34:14Z")

</div>

Have you seen `bkfact`?

```julia
help?> bkfact
search: bkfact bkfact!

  bkfact(A, uplo::Symbol=:U, symmetric::Bool=issymmetric(A), rook::Bool=false) -> BunchKaufman

  Compute the Bunch-Kaufman [^Bunch1977] factorization of a symmetric or
  Hermitian matrix A and return a BunchKaufman object. uplo indicates which
  triangle of matrix A to reference. If symmetric is true, A is assumed to be
  symmetric. If symmetric is false, A is assumed to be Hermitian. If rook is
  true, rook pivoting is used. If rook is false, rook pivoting is not used.
  The following functions are available for BunchKaufman objects: size, \,
  inv, issymmetric, ishermitian.

  | [^Bunch1977]
  |
  | J R Bunch and L Kaufman, Some stable methods for calculating
  | inertia and solving symmetric linear systems, Mathematics of
  | Computation 31:137 (1977), 163-179. url.

```

---

<div class="post-metadata">

**Author:** ![jacob-roth](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jacob-roth/32/1862_2.png) [@jacob-roth](https://discourse.julialang.org/u/jacob-roth)\
**Post date:** [February 8, 2018, 5:50am UTC](https://discourse.julialang.org/t/symmetric-indefinite-solvers-with-block-diagonal-pivoting/8541/3 "2018-02-08T05:50:02Z")

</div>

Thanks, sorry, I wasn’t clear. I’m aware of `bkfact` and `ldlfact`, but neither allow me to access `L` and `D` separately (i.e., if I need to do a Hessian modification). I’ll clarify my question.

---

<div class="post-metadata">

**Author:** ![andreasnoack](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/andreasnoack/32/27_2.png) [@andreasnoack](https://discourse.julialang.org/u/andreasnoack)\
**Post date:** [February 8, 2018, 7:28am UTC](https://discourse.julialang.org/t/symmetric-indefinite-solvers-with-block-diagonal-pivoting/8541/4 "2018-02-08T07:28:28Z")

</div>

> but neither allow me to access L and D

This is fixed on master where you can extract the factors with

```julia
julia> F = bkfact(A);

julia> F.U
10×10 UnitUpperTriangular{Float64,Array{Float64,2}}:
 1.0 -0.989314 1.76372 0.710251 -0.499019 … -0.127294 0.546951 -0.374852
  ⋅ 1.0 0.837403 -0.227568 0.715229 0.762732 0.216349 0.743333
  ⋅ ⋅ 1.0 -0.958846 -0.0750025 -0.308083 0.115865 -0.00222697
  ⋅ ⋅ ⋅ 1.0 0.00597047 0.876363 -0.129265 0.377775
  ⋅ ⋅ ⋅ ⋅ 1.0 0.552132 0.909079 0.333753
  ⋅ ⋅ ⋅ ⋅ ⋅ … -1.18999 -0.69973 0.0101069
  ⋅ ⋅ ⋅ ⋅ ⋅ -1.14554 0.269799 0.313673
  ⋅ ⋅ ⋅ ⋅ ⋅ 1.0 -0.756874 -0.263642
  ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ 1.0 0.0
  ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ 1.0

julia> F.D
10×10 Tridiagonal{Float64,Array{Float64,1}}:
 2.03091 0.0 ⋅ ⋅ ⋅ … ⋅ ⋅ ⋅ ⋅
 0.0 -4.87095 0.0 ⋅ ⋅ ⋅ ⋅ ⋅ ⋅
  ⋅ 0.0 -0.707841 0.0 ⋅ ⋅ ⋅ ⋅ ⋅
  ⋅ ⋅ 0.0 1.13327 0.0 ⋅ ⋅ ⋅ ⋅
  ⋅ ⋅ ⋅ 0.0 10.324 ⋅ ⋅ ⋅ ⋅
  ⋅ ⋅ ⋅ ⋅ 0.0 … 0.0 ⋅ ⋅ ⋅
  ⋅ ⋅ ⋅ ⋅ ⋅ 7.7207 0.0 ⋅ ⋅
  ⋅ ⋅ ⋅ ⋅ ⋅ 0.0 -2.39402 0.0 ⋅
  ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ 0.0 -0.130385 -3.29071
  ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ -3.29071 -0.464019

```

---

<div class="post-metadata">

**Author:** ![jacob-roth](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jacob-roth/32/1862_2.png) [@jacob-roth](https://discourse.julialang.org/u/jacob-roth)\
**Post date:** [February 8, 2018, 2:55pm UTC](https://discourse.julialang.org/t/symmetric-indefinite-solvers-with-block-diagonal-pivoting/8541/5 "2018-02-08T14:55:23Z")

</div>

Thanks, sorry I missed this. Does this do fill-reducing permutations / can it handle sparse matrices?

---

<div class="post-metadata">

**Author:** ![andreasnoack](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/andreasnoack/32/27_2.png) [@andreasnoack](https://discourse.julialang.org/u/andreasnoack)\
**Post date:** [February 9, 2018, 10:20am UTC](https://discourse.julialang.org/t/symmetric-indefinite-solvers-with-block-diagonal-pivoting/8541/6 "2018-02-09T10:20:34Z")

</div>

We only have a dense Bunch-Kaufman. The `ldlfact` works for sparse matrices and you can extract the factors but it only doesn’t fill-reducing permutations, i.e. no permutations during the numerical factorization.

---

<div class="post-metadata">

**Author:** ![dpo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpo/32/3335_2.png) [@dpo](https://discourse.julialang.org/u/dpo)\
**Post date:** [June 20, 2018, 2:14am UTC](https://discourse.julialang.org/t/symmetric-indefinite-solvers-with-block-diagonal-pivoting/8541/7 "2018-06-20T02:14:16Z")

</div>

See [HSL.jl](https://github.com/JuliaSmoothOptimizers/HSL.jl). I wrote an interface to HSL\_MA97 (essentially a multi-core version of HSL\_MA57) that gives access to the factors and allows you to modify them. In the `ma57` branch, there’s a similar interface to HSL\_MA57 but a patch is required to access the factors. I’m not quite sure I’m allowed to distribute it (I should submit it to the HSL folks). In the mean time, I could send it privately. Feel free to open an issue there if anything doesn’t work right.