# Algebraic equation in matrix form

**URL:** <https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621>\
**Category:** New to Julia\
**Tags:** question\
**Created:** [November 22, 2022, 3:11am UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621 "2022-11-22T03:11:41Z")\
**Posts on this page:** 14\
**Page:** 1

<div class="post-metadata">

**Author:** ![TrumeAAA](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/trumeaaa/32/44541_2.png) [@TrumeAAA](https://discourse.julialang.org/u/TrumeAAA)\
**Post date:** [November 22, 2022, 3:11am UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/1 "2022-11-22T03:11:42Z")

</div>

Hi, I want to solve a Algebraic equation in matrix form:

 ![Inked屏幕截图 2022-11-22 110800_LI](https://global.discourse-cdn.com/julialang/original/3X/0/0/00def1d5b9f5c89189925f54f2e8bd0a9dd95d87.jpeg)  
Actually, I need to solve AX=0, that is to say I need to solve det(A)=0 with respect to k’ when k in (-1,1).  
This is easy in matlab by running:

```julia
...
syms k;
kk = double(solve(det(A)));
...

```

But it takes a few hours to solve the whole problem(need to solve 10^3 same eqs, and in matlab for in for, very slow). So I tried to do it in Julia. However, I have no idea how to do it in Julia.  
So anyone knows how to do it?  
Any suggestion will be helpful!

---

<div class="post-metadata">

**Author:** ![sylvaticus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sylvaticus/32/203883_2.png) [@sylvaticus](https://discourse.julialang.org/u/sylvaticus)\
**Post date:** [November 22, 2022, 8:22am UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/2 "2022-11-22T08:22:54Z")

</div>

crosspost on SO: [julia - Algebraic equation in matrix form - Stack Overflow](https://stackoverflow.com/questions/74527017/algebraic-equation-in-matrix-form)

(I believe the policy is that crossposting if fine, but should be declared)

Edit: replaced the wrong link

---

<div class="post-metadata">

**Author:** ![baggepinnen](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/baggepinnen/32/693_2.png) [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)\
**Post date:** [November 22, 2022, 8:49am UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/3 "2022-11-22T08:49:30Z")

</div>

Do you have to solve it symbolically or would a numerical solution be okay? A numerical solution could likely be obtained in milliseconds.

---

<div class="post-metadata">

**Author:** ![TrumeAAA](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/trumeaaa/32/44541_2.png) [@TrumeAAA](https://discourse.julialang.org/u/TrumeAAA)\
**Post date:** [November 22, 2022, 1:03pm UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/4 "2022-11-22T13:03:07Z")

</div>

Hi, a numerical solution is ok. It is true that a numerical solution could be fast. However, I need to solve the eq for every k in (-1,1) and at every moment t in (0,1000) and this will cost a lot of time in matlab. But I am not familiar with Julia, so I dont know how to do it in Julia.

---

<div class="post-metadata">

**Author:** ![TrumeAAA](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/trumeaaa/32/44541_2.png) [@TrumeAAA](https://discourse.julialang.org/u/TrumeAAA)\
**Post date:** [November 22, 2022, 1:03pm UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/5 "2022-11-22T13:03:40Z")

</div>

Sorry, I dont understand…

---

<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:** [November 22, 2022, 1:14pm UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/6 "2022-11-22T13:14:37Z")

</div>

> [@TrumeAAA](#):
>
> Actually, I need to solve AX=0, that is to say I need to solve det(A)=0 with respect to k’ when k in (-1,1)

This is a [nonlinear eigenproblem](https://en.wikipedia.org/wiki/Nonlinear_eigenproblem); there are much more efficient ways to solve it than black-box root solvers for \det A. See e.g. [NonlinearEigenproblems.jl](https://github.com/nep-pack/NonlinearEigenproblems.jl).

Furthermore, if you are trying to find the roots as a function of the continuous parameter k, that is a [numerical continuation](https://en.wikipedia.org/wiki/Numerical_continuation) problem, and there are better algorithms than re-starting the root-finding from scratch for each k. See e.g. [BifurcationKit.jl](https://github.com/bifurcationkit/BifurcationKit.jl).

---

<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:** [November 22, 2022, 1:20pm UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/7 "2022-11-22T13:20:42Z")

</div>

> [@stevengj](#):
>
> This is a [nonlinear eigenproblem](https://en.wikipedia.org/wiki/Nonlinear_eigenproblem); there are much more efficient ways to solve it than black-box root solvers for \det A. See e.g. [NonlinearEigenproblems.jl](https://github.com/nep-pack/NonlinearEigenproblems.jl).

In fact, it’s even better than that — looking at your matrix, it’s actually a _linear_ eigenproblem in -k'^2 (assuming all of your \eta\_\ell factors are constants): if you just move the k'^2 terms to the other side, the equation you wrote takes the form Bc = -k'^2 c, so you can just use `sqrt.(.- eigvals(B))` to find the roots k'.

If your matrix is 7 \times 7 as shown in your post, you should be able to solve 1000 of these eigenproblems in a fraction of a second.

---

<div class="post-metadata">

**Author:** ![antoine-levitt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/antoine-levitt/32/4008_2.png) [@antoine-levitt](https://discourse.julialang.org/u/antoine-levitt)\
**Post date:** [November 22, 2022, 2:09pm UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/8 "2022-11-22T14:09:24Z")

</div>

Not to mention it’s diagonal plus convolution, so you can use iterative eigensolvers+FFTs+diagonal preconditioning to solve it at a cost linear in the number of unknowns. I imagine this comes from a band structure computation in quantum mechanics? You might be able to adapt [GitHub - JuliaMolSim/DFTK.jl: Density-functional toolkit](https://github.com/JuliaMolSim/DFTK.jl/) to your purposes (see eg [Custom potential · DFTK.jl](https://docs.dftk.org/dev/examples/custom_potential/) for how to incorporate a custom potential), although that might be overkill for this relatively simple problem.

---

<div class="post-metadata">

**Author:** ![baptiste](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/baptiste/32/9813_2.png) [@baptiste](https://discourse.julialang.org/u/baptiste)\
**Post date:** [November 22, 2022, 11:54pm UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/9 "2022-11-22T23:54:05Z")

</div>

What do you mean? This is my question, and I did not post on SO.

---

<div class="post-metadata">

**Author:** ![sylvaticus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sylvaticus/32/203883_2.png) [@sylvaticus](https://discourse.julialang.org/u/sylvaticus)\
**Post date:** [November 23, 2022, 4:10am UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/10 "2022-11-23T04:10:26Z")

</div>

Sorry, copypasted the wrong link, edited.

---

<div class="post-metadata">

**Author:** ![TrumeAAA](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/trumeaaa/32/44541_2.png) [@TrumeAAA](https://discourse.julialang.org/u/TrumeAAA)\
**Post date:** [November 24, 2022, 2:36am UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/11 "2022-11-24T02:36:19Z")

</div>

Thank you for your kind reply, It really helps a lot! Have a nice day!

---

<div class="post-metadata">

**Author:** ![TrumeAAA](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/trumeaaa/32/44541_2.png) [@TrumeAAA](https://discourse.julialang.org/u/TrumeAAA)\
**Post date:** [November 24, 2022, 2:39am UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/12 "2022-11-24T02:39:59Z")

</div>

> [@antoine-levitt](#):
>
> although that might be overkill for this relatively simple problem.

Hi, I do tried to solve a band structure. Actually it is a band structure of a photonic crystal of 1d in normal incident case. Thank you for your kind reply, It helps a lot!

---

<div class="post-metadata">

**Author:** ![TrumeAAA](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/trumeaaa/32/44541_2.png) [@TrumeAAA](https://discourse.julialang.org/u/TrumeAAA)\
**Post date:** [November 24, 2022, 2:40am UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/13 "2022-11-24T02:40:17Z")

</div>

Have a nice day!

---

<div class="post-metadata">

**Author:** ![TrumeAAA](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/trumeaaa/32/44541_2.png) [@TrumeAAA](https://discourse.julialang.org/u/TrumeAAA)\
**Post date:** [November 26, 2022, 7:27am UTC](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/14 "2022-11-26T07:27:05Z")

</div>

I tried to solve the problem directly as an algebraic equation and it took a lot of time in matlab. However time cost can be reduced greatly as I solved it in matlab as a eigenproblem as stenengj said.Total time was reduced to 20s-40s in matlab and 8s or so in Julia. Still a lot to learn…

> [@Algebraic equation in matrix form](https://discourse.julialang.org/t/algebraic-equation-in-matrix-form/90621/7):
>
> In fact, it’s even better than that — looking at your matrix, it’s actually a linear eigenproblem in -k'^2 (assuming all of your \eta\_\ell factors are constants): if you just move the k'^2 terms to the other side, the equation you wrote takes the form Bc = -k'^2 c, so you can just use sqrt.(.- eigvals(B)) to find the roots k'. If your matrix is 7 \times 7 as shown in your post, you should be able to solve 1000 of these eigenproblems in a fraction of a second.
