# HG changeset patch
# User Jaroslav Hajek <highegg@gmail.com>
# Date 1254754993 -7200
# Node ID 9fba7e1da785baa3ec1eea0b14668d5184308736
# Parent  50db3c5175b567e783d580374962a91be6d4fdc8
correct algorithm for perm matrix det (sign)

diff --git a/liboctave/ChangeLog b/liboctave/ChangeLog
--- a/liboctave/ChangeLog
+++ b/liboctave/ChangeLog
@@ -1,3 +1,8 @@
+2009-10-05  Jaroslav Hajek  <highegg@gmail.com>
+
+	* PermMatrix.cc (PermMatrix::determinant): Implement a (hopefully)
+	working algorithm.
+
 2009-10-05  Jaroslav Hajek  <highegg@gmail.com>
 
 	* dim-vector.h (operator ==): Include fast case.
diff --git a/liboctave/PermMatrix.cc b/liboctave/PermMatrix.cc
--- a/liboctave/PermMatrix.cc
+++ b/liboctave/PermMatrix.cc
@@ -28,6 +28,7 @@
 #include "idx-vector.h"
 #include "error.h"
 #include "Array-util.h"
+#include "oct-locbuf.h"
 
 static void
 gripe_invalid_permutation (void)
@@ -100,16 +101,30 @@
 octave_idx_type 
 PermMatrix::determinant (void) const
 {
-  Array<octave_idx_type> pa = *this;
-  octave_idx_type len = pa.length (), *p = pa.fortran_vec ();
-  bool neg = false;
+  // Determine the sign of a permutation in linear time.
+  // Is this widely known?
+
+  octave_idx_type len = perm_length ();
+  const octave_idx_type *pa = data ();
+
+  OCTAVE_LOCAL_BUFFER (octave_idx_type, p, len);
+  OCTAVE_LOCAL_BUFFER (octave_idx_type, q, len);
+
   for (octave_idx_type i = 0; i < len; i++)
     {
-      octave_idx_type j = p[i];
+      p[i] = pa[i];
+      q[p[i]] = i;
+    }
+
+  bool neg = false;
+
+  for (octave_idx_type i = 0; i < len; i++)
+    {
+      octave_idx_type j = p[i], k = q[i];
       if (j != i)
         {
-          p[i] = p[j];
-          p[j] = j;
+          p[k] = p[i];
+          q[j] = q[i];
           neg = ! neg;
         }
     }