# LinearAlgebra.det and symbolic variables

**URL:** <https://discourse.julialang.org/t/linearalgebra-det-and-symbolic-variables/83054>\
**Category:** Performance\
**Tags:** question, linearalgebra, symbolics\
**Created:** [June 20, 2022, 9:52am UTC](https://discourse.julialang.org/t/linearalgebra-det-and-symbolic-variables/83054 "2022-06-20T09:52:18Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![DanDoe](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dandoe/32/52717_2.png) [@DanDoe](https://discourse.julialang.org/u/DanDoe)\
**Post date:** [June 20, 2022, 9:52am UTC](https://discourse.julialang.org/t/linearalgebra-det-and-symbolic-variables/83054/1 "2022-06-20T09:52:18Z")

</div>

I am interested in the computation of the determinant of a (very sparse) matrix involving a symbolic variable.

My code is

```julia
using Symbolics, LinearAlgebra, SparseArrays

N = 11

FloatType = Float64

c = zeros(FloatType, N)
for k in 2:N
  c[k] = (k - 1)/(2.0*(N - 1))
end

A = zeros(FloatType, N, N)
A[2, 1] = c[2]
for k in 3:N
  A[k, 1] = c[k] - 1
  A[k, k-1] = 1
end

A = sparse(A)
display(A); println()

b = zeros(FloatType, N)
b[N] = 1.0
b = sparse(b)

@variables z
display(I - z * A + z * ones(FloatType, N) * b'); println()
println("Start computing determinant")
Determinant = det(I - z * A + z * ones(FloatType, N) * b')
println("Start expanding")
Determinant = expand(Determinant)

println(Determinant)

```

So I use the [`det`](https://docs.julialang.org/en/v1/stdlib/LinearAlgebra/#LinearAlgebra.det) from Julia Linear Algebra and the symbolic variable from [`Symbolics.jl`](https://github.com/JuliaSymbolics/Symbolics.jl).

For small `N` the program runs somewhat quick, but already for the given `N=11` it takes literally ages. Is there any hope to speed this up?

I can run the same problem for instance in Mathematica in basically no time, even for `N=100` it returns instantly.

---

<div class="post-metadata">

**Author:** ![greatpet](https://avatars.discourse-cdn.com/v4/letter/g/e495f1/32.png) [@greatpet](https://discourse.julialang.org/u/greatpet)\
**Post date:** [June 21, 2022, 10:01am UTC](https://discourse.julialang.org/t/linearalgebra-det-and-symbolic-variables/83054/2 "2022-06-21T10:01:35Z")

</div>

Here’s a solution using [Nemo.jl](https://nemocas.github.io/Nemo.jl/stable/), which finishes under 3 seconds for N=100. Nemo has its own matrix type, which is a little different from Julia’s own `Matrix`. I also use Nemo’s rational number type `QQ` instead of `Float64`. Unfortunately, Nemo does not support sparse matrices, so I’m computing with a dense matrix, which will cause memory usage problems if `N` is really large.

```julia
using Nemo, LinearAlgebra

N = 100

r, z = PolynomialRing(QQ, "z")

c = fill(r(0), N) # `r(0)` is the zero polynomial. `c` is a vector of polynomials
for k in 2:N
  c[k] = r((k - 1)//(2*(N - 1)))
end

A = MatrixSpace(r, N, N)() # This is using Nemo's matrix type. Intialized to zero matrix
A[2, 1] = c[2]
for k in 3:N
  A[k, 1] = c[k] - r(1)
  A[k, k-1] = r(1)
end

b = MatrixSpace(r, 1, N)()
b[1, N] = r(1)

matrix = MatrixSpace(r, N, N)(1) - z * A + z * MatrixSpace(r, N, 1)(r(1)) * b

display(Matrix(matrix)) # convert back to Julia's native Matrix type before calling `display`, as Nemo's printing of matrices is not as pretty

@time Determinant = det(matrix) # Nemo expands polynomials by default. No need to expand further

```

---

<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, 10:03am UTC](https://discourse.julialang.org/t/linearalgebra-det-and-symbolic-variables/83054/3 "2022-06-21T10:03:06Z")

</div>

This doesn’t have a specialization. It’s just using a basic fallback. It’s worth an issue to document that a specialization should be made.

---

<div class="post-metadata">

**Author:** ![DanDoe](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dandoe/32/52717_2.png) [@DanDoe](https://discourse.julialang.org/u/DanDoe)\
**Post date:** [June 22, 2022, 9:01am UTC](https://discourse.julialang.org/t/linearalgebra-det-and-symbolic-variables/83054/4 "2022-06-22T09:01:40Z")

</div>

Thanks, do you also know a function how to convert the coefficients of the determinant of type `fmpq` back to something like `Float64` ? I couldn’t find something immediately in the docs.

---

<div class="post-metadata">

**Author:** ![greatpet](https://avatars.discourse-cdn.com/v4/letter/g/e495f1/32.png) [@greatpet](https://discourse.julialang.org/u/greatpet)\
**Post date:** [June 22, 2022, 10:51am UTC](https://discourse.julialang.org/t/linearalgebra-det-and-symbolic-variables/83054/5 "2022-06-22T10:51:32Z")

</div>

Nemo provides conversion between `fmpq` and Julia’s `Rational{BigInt}`, so you can manually write a conversion function to `Float64`:

```julia
function Base.convert(::Type{Float64}, a::Nemo.fmpq)
    convert(Float64, convert(Rational{BigInt}, a))
end

```

Usage example:

```julia
julia> convert(Float64, fmpq(1//2))
0.5

```

---

<div class="post-metadata">

**Author:** ![DanDoe](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dandoe/32/52717_2.png) [@DanDoe](https://discourse.julialang.org/u/DanDoe)\
**Post date:** [June 22, 2022, 3:22pm UTC](https://discourse.julialang.org/t/linearalgebra-det-and-symbolic-variables/83054/6 "2022-06-22T15:22:50Z")

</div>

I found a small bug. `MatrixSpace(r, N, 1)(r(1))` gives not a whole vector of ones, but initializes only the first entry to one. My solution to this is

```julia
ones = MatrixSpace(r, N, 1)()
for k = 1:N
  ones[k, 1] = r(1)
end

matrix = MatrixSpace(r, N, N)(1) - z * AMat1 + z * ones * b

```

---

<div class="post-metadata">

**Author:** ![greatpet](https://avatars.discourse-cdn.com/v4/letter/g/e495f1/32.png) [@greatpet](https://discourse.julialang.org/u/greatpet)\
**Post date:** [June 22, 2022, 3:34pm UTC](https://discourse.julialang.org/t/linearalgebra-det-and-symbolic-variables/83054/7 "2022-06-22T15:34:08Z")

</div>

You’re right. The constructor only fills the diagonal elements, which is just one element in this special case. You could also first construct a Julia matrix and then convert it.

```julia
ones = MatrixSpace(r, N, 1)(fill(r(1), N, 1))

```
