[Octave-bug-tracker] [bug #46103] inversion of a singular diagonal matrix A produces inv(A)==A

[jan]

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



With all due respect, inverse is mathematically a different concept to
pseudo-inverse. And no function should knowlingly give such a misleading
result that  can give very odd behaviour. 

And since the diagonal matrix is now produced by e.g. diag(), already the
bw-compatibility requires this. 

It seems that simply the following patch would do, but as this is my first
peak into octave source, I might be wrong. 


address@hidden:/opt/octave/latest/octave$ hg diff
diff -r 6256f6e366ac liboctave/array/fCDiagMatrix.cc
--- a/liboctave/array/fCDiagMatrix.cc   Fri Oct 02 05:50:43 2015 +0200
+++ b/liboctave/array/fCDiagMatrix.cc   Fri Oct 02 23:35:33 2015 +0300
@@ -375,10 +375,7 @@
   for (octave_idx_type i = 0; i < length (); i++)
     {
       if (elem (i, i) == 0.0f)
-        {
-          info = -1;
-          return *this;
-        }
+        retval.elem (i, i) = std::numeric_limits<double>::infinity ();
       else
         retval.elem (i, i) = 1.0f / elem (i, i);
     }
diff -r 6256f6e366ac liboctave/array/fDiagMatrix.cc
--- a/liboctave/array/fDiagMatrix.cc    Fri Oct 02 05:50:43 2015 +0200
+++ b/liboctave/array/fDiagMatrix.cc    Fri Oct 02 23:35:33 2015 +0300
@@ -280,10 +280,7 @@
   for (octave_idx_type i = 0; i < len; i++)
     {
       if (elem (i, i) == 0.0)
-        {
-          info = -1;
-          return *this;
-        }
+        retval.elem (i, i) = std::numeric_limits<double>::infinity ();
       else
         retval.elem (i, i) = 1.0 / elem (i, i);
     }

 


    _______________________________________________________

Reply to this item at:

  <http://savannah.gnu.org/bugs/?46103>

_______________________________________________
  Meddelandet skickades via/av Savannah
  http://savannah.gnu.org/