- expressions can now be read from streams; the input expression can contain

[ginac.git] / ginac / matrix.cpp

diff --git a/ginac/matrix.cpp b/ginac/matrix.cpp
index 37716290f77317e4916b1705073ba3834332b1ce..f219f27de0c3dcb1432bcc7d944b4ef8bcd63dcc 100644 (file)
--- a/ginac/matrix.cpp
+++ b/ginac/matrix.cpp
@@ -26,9 +26,10 @@
 
 #include "matrix.h"
 #include "archive.h"
+#include "numeric.h"
+#include "lst.h"
 #include "utils.h"
 #include "debugmsg.h"
-#include "numeric.h"
 
 #ifndef NO_NAMESPACE_GINAC
 namespace GiNaC {
@@ -469,6 +470,16 @@ ex matrix::determinant(void) const
     if (numeric_flag)
         return determinant_numeric();
     
+    if (5*sparse_count<row*col) {     // MAGIC, maybe 10 some bright day?
+        matrix M(*this);
+        // int sign = M.division_free_elimination();
+        int sign = M.fraction_free_elimination();
+        if (normal_flag)
+            return sign*M(row-1,col-1).normal();
+        else
+            return sign*M(row-1,col-1).expand();
+    }
+    
     // Now come the minor expansion schemes.  We always develop such that the
     // smallest minors (i.e, the trivial 1x1 ones) are on the rightmost column.
     // For this to be efficient it turns out that the emptiest columns (i.e.
@@ -498,11 +509,6 @@ ex matrix::determinant(void) const
             result[r*col+c] = m[r*col+(*i)];
     }
     
-    if (0*sparse_count>row*col) {     // MAGIC, maybe 10 some day?
-        if (normal_flag)
-            return sign*matrix(row,col,result).determinant_minor_sparse().normal();
-        return sign*matrix(row,col,result).determinant_minor_sparse();
-    }
     if (normal_flag)
         return sign*matrix(row,col,result).determinant_minor_dense().normal();
     return sign*matrix(row,col,result).determinant_minor_dense();
@@ -1213,6 +1219,29 @@ int matrix::pivot(unsigned ro, bool symbolic)
     return 0;
 }
 
+/** Convert list of lists to matrix. */
+ex lst_to_matrix(const ex &l)
+{
+       if (!is_ex_of_type(l, lst))
+               throw(std::invalid_argument("argument to lst_to_matrix() must be a lst"));
+
+       // Find number of rows and columns
+       unsigned rows = l.nops(), cols = 0, i, j;
+       for (i=0; i<rows; i++)
+               if (l.op(i).nops() > cols)
+                       cols = l.op(i).nops();
+
+       // Allocate and fill matrix
+       matrix &m = *new matrix(rows, cols);
+       for (i=0; i<rows; i++)
+               for (j=0; j<cols; j++)
+                       if (l.op(i).nops() > j)
+                               m.set(i, j, l.op(i).op(j));
+                       else
+                               m.set(i, j, ex(0));
+       return m;
+}
+
 //////////
 // global constants
 //////////