I have the impression that this has to do with the precision and size of your coefficients,
and more precisely of the magnitude difference between them.
Your constant coefficient is around 4e18.
If I scale the solutions by a factor 10, by dividing coefficient 3 by 100 and coefficient 5 by 10000, I obtain the correct result (to be interpreted 10x larger):
VNL roots are
4554.6 +i 4.62223e-32
-4457.61 +i -1.44445e-34
-4554.6 +i -4.62223e-32
4457.61 +i 1.44445e-34
So try reducing too extreme values first (e.g. by dividing coefficient n by 10^n) before calling the solver.
I presume that the original algorithm (by Kriegman and Ponce) has this same instability.
Not sure whether this scaling solution should be introduced into the VNL interface layer, or whether this should better be left to the user who calls the method. Suggestions welcome, of course.
--- Den ons 2010-09-01 skrev somi <email@example.com>: