# How to find zeros of 16th degree polynomial with coefficients which contain one symbolic parameter?

I'm trying to find eigenvalues of matrix which is 16x16. Here is a part of matrix:

{{1/4 (1. + (1.16573*10^-15 - 3.81639*10^-16 I) W^(-0.0184391 -
9.99111 I) + <<12>> + (1.05471*10^-15 +
7.80626*10^-17 I) W^(-0.00601834 + 10.0408 I) + 1./
W^0.0034662), 0, 0, -0.125011, 0, 0, -0.125011, 0, 0, 0, 0, 0, 0,
0, 0, 0.125}, <<14>>, {0.125, <<14>>,  1/4 (<<19>> + <<18>> +   <<19>>/<<1>>)}}


Matrix is symmetric, with complex numbers and one symbolic-W.

My strategy was, in order to facilitate computation, to find determinant first

f[λ_] = Det[\[ScriptCapitalS] - λ IdentityMatrix[16]]


and then to solve for roots

Solve[f[λ ] == 0, λ].


But, it seems to me that computation takes too much time -- about 15 hours at the time of writing.

Is there any efficient way to take this problem in Mathematica? I can think of NSolve, FindRoot but I am not sure that will work at all (Of course, Eigenvalues is already tested, but without results). Any help is appreciated in relation to this issue.

Remark: Computation is carried out on a computer which has 4 Gb RAM and dual core processor. The software in question is Mathematica 9.0 .

Well, here is 3x3 version of my matrix, for the sake of ilustration:

{{1/4 (1. + (0.0274419 - 0.0676919 I) W^(-0.013233 -
0.0501142 I) + (0.0274419 + 0.0676919 I) W^(-0.013233 +
0.0501142 I) + 1./W^0.0131338 + 1./W^0.00896667 + 1./
W^0.0034662), 0, 0}, {0,
1/4 (0.0000907999 + 0.818533 W^(-0.0184391 - 9.99111 I) +
0.818533 W^(-0.0184391 +
9.99111 I) - (0.0298889 - 0.0435723 I) W^(-0.013233 -
0.0501142 I) - (0.0298889 + 0.0435723 I) W^(-0.013233 +
0.0501142 I) + 0.145971/W^0.0131338 +
0.810095 W^(-0.0110314 - 9.99044 I) +
0.810095 W^(-0.0110314 + 9.99044 I) - 1./
W^0.00896667 - (0.14229 - 0.112623 I) W^(-0.0075773 -
10.0401 I) - (0.14229 + 0.112623 I) W^(-0.0075773 +
10.0401 I) - (0.180361 - 0.0594458 I) W^(-0.00601834 -
10.0408 I) - (0.180361 + 0.0594458 I) W^(-0.00601834 +
10.0408 I) - 0.785802/W^0.0034662),
1/4 ((0. + 0.405084 I) W^(-0.0110314 -
9.99044 I) - (0. + 0.405084 I) W^(-0.0110314 +
9.99044 I) - (0.0297256 + 0.0901888 I) W^(-0.00601834 -
10.0408 I) - (0.0297256 - 0.0901888 I) W^(-0.00601834 +
10.0408 I))}, {0,
1/4 ((0. + 0.405084 I) W^(-0.0110314 -
9.99044 I) - (0. + 0.405084 I) W^(-0.0110314 +
9.99044 I) - (0.0297256 + 0.0901888 I) W^(-0.00601834 -
10.0408 I) - (0.0297256 - 0.0901888 I) W^(-0.00601834 +
10.0408 I)),
1/4 (0.500091 - 0.818533 W^(-0.0184391 - 9.99111 I) -
0.818533 W^(-0.0184391 +
9.99111 I) - (0.28763 + 0.102747 I) W^(-0.013233 -
0.0501142 I) - (0.28763 - 0.102747 I) W^(-0.013233 +
0.0501142 I) + 0.213796/W^0.0131338 -
0.405011 W^(-0.0110314 - 9.99044 I) -
0.405011 W^(-0.0110314 + 9.99044 I) - 1./
W^0.00896667 + (0.14229 - 0.112623 I) W^(-0.0075773 -
10.0401 I) + (0.14229 + 0.112623 I) W^(-0.0075773 +
10.0401 I) + (0.0901724 - 0.0297202 I) W^(-0.00601834 -
10.0408 I) + (0.0901724 + 0.0297202 I) W^(-0.00601834 +
10.0408 I) - 0.303692/W^0.0034662)}}


Mathematica works this example fine and gives eigenvalues as a function of W.

• The incomplete matrix is not really helpful with regard to actually evaluating your code. Commented Jun 19, 2015 at 14:11
• What's your ultimate goal here? Do you want a function f[W] that returns the eigenvalues for a given W? Is a closed-form expression for W needed, or are numerical answers for specific W values OK? Do you want a graph of the real and/or imaginary roots as a function of W? A parametric plot of the root loci in the complex plane? Any guidance as to what your final goal is will be helpful in figuring out how to help you. Commented Jun 19, 2015 at 14:12
• @YvesKlett, matrix is too big to be presented here, unfortunately. The point is that entries of the matrix contain that parameter W, and consequently that W goes to the polynomial coefficients, which in my opinion complicates calculation. But, W is unavoidable.
– Moki
Commented Jun 19, 2015 at 14:43
• @MichaelSeifert, ultimate goal is to find zeros (roots) of given polynomial as functions of W. The next step would be finding eigenvectors. Actually, behind W is time and this is connected with dynamics. I hope now it is clearer where the problem lies.
– Moki
Commented Jun 19, 2015 at 14:52
• Is the matrix "symmetric", or is it Hermitian? That will make a difference as to whether you'll have to deal with complex eigenvalues & eigenvectors or not. Commented Jun 19, 2015 at 15:19

There is, of course, no closed-form solution for the eigenvalues and eigenvectors of a 16 x 16 matrix (or, more generally, for the roots of an $n$th-order polynomial with $n > 5$.) So you're going to be stuck with numerical methods no matter what you do. The easiest way to create a graph of the eigenvalues as a function of W will be to define a function using SetDelayed (i.e., :=) and then graph that function, like so:

mat1[W_] = RandomComplex[{-1 - I, 1 + I}, {16, 16}] + W RandomComplex[{-1 - I, 1 + I}, {16, 16}];
\[ScriptCapitalS][W_] = mat1[W] + Transpose[mat1[W]];
evals[W_] := Sort[Eigenvalues[\[ScriptCapitalS][W]]]
Plot[Re[evals[W]], {W, 0, 1}]


Creating this graph takes about 5 seconds on my machine; it might take a bit longer on yours, since your machine is older than mine, and also your definition for $\mathcal{S}$ involves complex powers of W (whereas mine is just linear in W.) Still, I doubt it will take 15 hours.

Note the use of the := in the definition of evals[W]. This tells Mathematica that you don't want the result evaluated immediately; you only want it evaluated when you call the function later. In this case, this means that Mathematica doesn't waste time trying to store the full expressions for the eigenvalues in terms of W (which is what I suspect your machine is getting stuck on.) instead, it just applies the numbers fed to it by Plot when it's called later on.

I hope this helps, but I admit that I'm not 100% from your question what your exact problem is. If this doesn't help, feel free to tell me so in the comments and I'll see what I can do to address it.

• "no closed-form solution" - there are closed forms, but they are somewhat unwieldy, which is why Mathematica uses Root[] instead. Commented Jun 19, 2015 at 18:04
• @Guesswhoitis, even solutions with Root are welcome. I do not expect that polynomial of this degree has "simple" roots.
– Moki
Commented Jun 20, 2015 at 11:43
• I have found your answer very interesting, but it seems that is more than I looking for, in some respects.
– Moki
Commented Jun 20, 2015 at 11:46
• Please find 3x3 matrix in my edited post. Does it case of 16x16 matrix seems solvable? Or other strategies must be employed?
– Moki
Commented Jun 20, 2015 at 12:03