[Octave-bug-tracker] [bug #60786] mpower: negative scalar to power of matrix containing non-integers returns nan

[John W. Eaton]

Follow-up Comment #8, bug #60786 (project octave):

Rik: This code has origins even earlier than I thought.  The first commit in
the hg archive dated August 1993 has this code for the scalar x matrix case:


tree_constant
xpow (double a, Matrix& b)
{
  tree_constant retval;

  int nr = b.rows ();
  int nc = b.columns ();

  if (nr == 0 || nc == 0 || nr != nc)
    error ("for x^A, A must be square");
  else
    {
      EIG b_eig (b);
      ComplexColumnVector lambda (b_eig.eigenvalues ());
      ComplexMatrix Q (b_eig.eigenvectors ());

      for (int i = 0; i < nr; i++)
        {
          Complex elt = lambda.elem (i);
          if (imag (elt) == 0.0)
            lambda.elem (i) = pow (a, real (elt));
          else
            lambda.elem (i) = pow (a, elt);
        }
      ComplexDiagMatrix D (lambda);

      ComplexMatrix result = Q * D * Q.inverse ();
      retval = tree_constant (result);
    }

  return retval;
}


so the check for real eigenvalues and the conversion to real has been in the
code from the beginning.  If this isn't the part of the code you were asking
about yesterday, then please let me know which function it was because I'm
curious now to know if there was a later change that introduced the complex to
real conversion that you were asking about.

And, if there is a better way to do this now using templates, I'm all in favor
of fixing it, but only if someone has time and is interested in working on it,
of course.

BTW, the hg archive apparently doesn't have ALL changes since it begins with a
date of 1993.  I guess that's the earliest change that was in the CVS archive
that we converted to hg.  I'm fairly sure I used some form of version control
from the very beginning though (probably RCS), so it is a little surprising
that there are no earlier changes recorded since I started working on Octave
in February of 1992.

    _______________________________________________________

Reply to this item at:

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

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