changeset 8954:97c84c4c2247

Make the column permutation vector in sparse LU cols()-long.
author Jason Riedy <jason@acm.org>
date Tue, 10 Mar 2009 21:54:34 -0400
parents a6945f92b868
children 6d3fcbf89267
files liboctave/ChangeLog liboctave/sparse-base-lu.cc test/ChangeLog test/build_sparse_tests.sh
diffstat 4 files changed, 56 insertions(+), 5 deletions(-) [+]
line wrap: on
line diff
--- a/liboctave/ChangeLog
+++ b/liboctave/ChangeLog
@@ -1,3 +1,8 @@
+2009-03-10  Jason Riedy  <jason@acm.org>
+
+	* sparse-base-lu.cc (Pc_vec): The column permutation should be
+	Ufact.cols ()-long, not Lfact.rows ()-long.
+
 2009-03-10  Jason Riedy  <jason@acm.org>
 
 	* dSparse.cc (SparseMatrix::SparseMatrix (const PermMatrix&)):
--- a/liboctave/sparse-base-lu.cc
+++ b/liboctave/sparse-base-lu.cc
@@ -121,11 +121,11 @@
 sparse_base_lu <lu_type, lu_elt_type, p_type, p_elt_type> :: Pc_vec (void) const
 {
 
-  octave_idx_type nr = Lfact.rows ();
+  octave_idx_type nc = Ufact.cols ();
 
-  ColumnVector Pout (nr);
+  ColumnVector Pout (nc);
 
-  for (octave_idx_type i = 0; i < nr; i++)
+  for (octave_idx_type i = 0; i < nc; i++)
     Pout.xelem (i) = static_cast<double> (Q(i) + 1);
 
   return Pout;
--- a/test/ChangeLog
+++ b/test/ChangeLog
@@ -1,3 +1,7 @@
+2009-03-10  Jason Riedy  <jason@acm.org>
+
+	* build_sparse_tests.sh: Add LU tests to the rectangular tests.
+
 2009-03-10  Jason Riedy  <jason@acm.org>
 
 	* test_diag_perm.m: Add a test for conversion to sparse form.
--- a/test/build_sparse_tests.sh
+++ b/test/build_sparse_tests.sh
@@ -665,8 +665,8 @@
 %! assert(j-i>=0);
 
 %!testif HAVE_UMFPACK ;# LU with vector permutations
-%! [L,U,P] = lu(bs,'vector');
-%! assert(L(P,:)*U,bs,1e-10);
+%! [L,U,P,Q] = lu(bs,'vector');
+%! assert(L(P,:)*U(:,Q),bs,1e-10);
 %! # triangularity
 %! [i,j,v]=find(L);
 %! assert(i-j>=0);
@@ -745,6 +745,48 @@
     # gen_divop_tests # Disable rectangular \ and / for now
     gen_matrixdiag_tests
     gen_matrixreshape_tests
+    cat >>$TESTS <<EOF
+%!testif HAVE_UMFPACK ;# permuted LU
+%! [L,U] = lu(bs);
+%! assert(L*U,bs,1e-10);
+
+%!testif HAVE_UMFPACK ;# simple LU + row permutations
+%! [L,U,P] = lu(bs);
+%! assert(P'*L*U,bs,1e-10);
+%! # triangularity
+%! [i,j,v]=find(L);
+%! assert(i-j>=0);
+%! [i,j,v]=find(U);
+%! assert(j-i>=0);
+
+%!testif HAVE_UMFPACK ;# simple LU + row/col permutations
+%! [L,U,P,Q] = lu(bs);
+%! assert(P'*L*U*Q',bs,1e-10);
+%! # triangularity
+%! [i,j,v]=find(L);
+%! assert(i-j>=0);
+%! [i,j,v]=find(U);
+%! assert(j-i>=0);
+
+%!testif HAVE_UMFPACK ;# LU with vector permutations
+%! [L,U,P,Q] = lu(bs,'vector');
+%! assert(L(P,:)*U(:,Q),bs,1e-10);
+%! # triangularity
+%! [i,j,v]=find(L);
+%! assert(i-j>=0);
+%! [i,j,v]=find(U);
+%! assert(j-i>=0);
+
+%!testif HAVE_UMFPACK ;# LU with scaling
+%! [L,U,P,Q,R] = lu(bs);
+%! assert(R*P'*L*U*Q',bs,1e-10);
+%! # triangularity
+%! [i,j,v]=find(L);
+%! assert(i-j>=0);
+%! [i,j,v]=find(U);
+%! assert(j-i>=0);
+
+EOF
 }