# Mathematica 8 can not solve Eigenvalue problem

I cannot explain this eigenvalue problem. Both matrices are positive definite and symmetric, but some solutions are complex.

But if it is matrices lower order, solutions are real. As I increase the number of order the matrices, I got more and more complex solutions (and in same time matrices are positive definite and symmetric.) What is the problem? I try to increase precision, but nothing.

From here you can copy and paste in Mathematica worksheet to see.

http://pastebin.com/raw.php?i=hM6DBVk5

Regards

-
Oh my. It would really have been useful to have pasted in stuff in InputForm[]; I know Mathematica enough to analyze snippets even without the program, but pasting in a whole notebook is sorta confusing in plaintext. – J. M. Aug 24 '11 at 2:47

The problem here is that elements of the matrix $B$ are floating point numbers. And may be somewhere along the way loss of accuracy happens. If to make $B$ a matrix with rational elements, then the answers will be real. For example,

B = Round[B 10^30]/10^30;
D0 := -x*A + B;
Solve[Det[D0] == 0, x]//N


gives

{{x -> 14847.993867199206}, {x -> 148848.79011442242}, {x ->
609570.78149666}, {x -> 1.8969502426287285*^6},
{x -> 8.338192814798254*^6}, {x -> 2.2248467905151808*^8}, {x ->
2.3271962509336177*^8}, {x -> 2.4808572853336135*^8},
{x -> 2.746203941496255*^8}, {x -> 3.9842849588286966*^8}, {x ->
1978.2943756360853}, {x -> 55725.53921427025},
{x -> 321863.65320953313}, {x -> 1.1383725875194645*^6}, {x ->
5.278252037812979*^6}, {x -> 2.1976933210891974*^8},
{x -> 2.2745356114392006*^8}, {x -> 2.4019361153335622*^8}, {x ->
2.6138647573316988*^8}, {x -> 3.547388273684647*^8}}

-
For this case is ok, but if I increase the order of matrices, I had again the complex solutions. See here pastebin.com/raw.php?i=2xxyD9vV – derdack Aug 24 '11 at 1:54
I suspect really bad conditioning; I'd try out SchurDecomposition[] (with the Pivoting option enabled if need be) on the original matrix and check for any 2-by-2 blocks... – J. M. Aug 24 '11 at 2:48
@derdack Here is again is a loss of accuracy since //N uses floating points arithmetics. This will give real answers: N[Solve[ Det[D0] == 0, x], 20]. – Andrew Aug 24 '11 at 7:28
Congratulations. Thank you very much Mr Andrew! – derdack Aug 24 '11 at 10:10
Can I decrease the time for computing because I used just one core for solving this determinant and when matrices are high order, I had a problem. – derdack Aug 24 '11 at 11:58

I am not sure if this the appropriate place to be asking this question. In the future, I would try maybe linking to a CDF document or giving a summary of the problem you are solving so that people can weigh in without being Mathematica gurus.

The problem here is that you are doing floating point arithmetic. The results are not accurate because of it. Maybe the matrices have a bad condition number. Mathematica is fortunately capable of working with numbers symbolically which is what you will want to do here.

In particular notice that the matrix b is filled with floating point numbers. These can be converted into symbolically processable numbers with the Rationalize command.

After defining the matrix b evaluate:

B = Rationalize[B, 0]

B will be filled with rational numbers instead of floating point numbers.

The results from Solve are all now real numbers.

-
You will probably want to use NSolve instead of Solve. NSolve gives out numeric answers. – Searke Aug 23 '11 at 22:13