[Octave-bug-tracker] [bug #53897] poly of conjugate complex pairs should return a real polynomial

[Rik]

Follow-up Comment #3, bug #53897 (project octave):

Yes, but is it on the right track?  Using your example from comment #2 I find
that the largest imaginary value is 10 eps which is quite small.  At least
there are no horrendous errors like 256i as in the first example.

It seems that sort followed by a test to replace small values with 0 would
work.  The algorithm in poly.m is


  for j = 1:n;
    y(2:(j+1)) = y(2:(j+1)) - v(j) .* y(1:j);
  endfor


which has a multiplication and a subtraction per round of the algorithm.  To
me that says worst case of 2 eps per cycle or 2 x n x eps for the algorithm. 
It's not perfect, but for eps one could use the largest value in v.


epsmax = 2* n * eps (max (abs (v(:))));
idx = abs (imag (y)) < epsmax;
y(idx) = real (y(idx));


Of course, this is starting to look like a fair amount of work so if detecting
the conjugate pairs is easier then we could do that as you suggested in
comment #0.


    _______________________________________________________

Reply to this item at:

  <https://savannah.gnu.org/bugs/?53897>

_______________________________________________
  Message sent via Savannah
  https://savannah.gnu.org/