Mercurial > octave

changeset 22241:71d86e88589f

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
chol2inv: fix support for sparse matrices by fixing MatrixType (bug #36437) * liboctave/array/MatrixType.cc: a Banded matrix can also be Upper or Lower. In such case, set type to Lower/Upper and not to Banded. The check for Banded matrix does not exist when input is non-sparse, which is why this only affected sparse matrices. * libinterp/dldfcn/chol.cc: add tests for sparse matrices.
author Barbara Locsi <locsi.barbara@gmail.com>
date Sat, 06 Aug 2016 19:38:32 +0200
parents 669fc8cf1fdd
children 5417dad26a25
files libinterp/dldfcn/chol.cc liboctave/array/MatrixType.cc
diffstat 2 files changed, 35 insertions(+), 2 deletions(-) [+]
line wrap: on
line diff
--- a/libinterp/dldfcn/chol.cc	Tue Aug 09 15:07:30 2016 -0700
+++ b/libinterp/dldfcn/chol.cc	Sat Aug 06 19:38:32 2016 +0200
@@ -570,6 +570,35 @@
   return retval;
 }
 
+/*
+
+## Test for bug #36437
+%!function sparse_chol2inv (T, tol)
+%!  iT = inv (T);
+%!  ciT = chol2inv (chol (T));
+%!  assert (ciT, iT, tol);
+%!  assert (chol2inv (chol ( full (T))), ciT, tol*2);
+%!endfunction
+
+%!test
+%! A = gallery ("poisson", 3);
+%! sparse_chol2inv (A, eps);
+
+%!test
+%! n = 10;
+%! B = spdiags (ones (n, 1) * [1 2 1], [-1 0 1], n, n);
+%! sparse_chol2inv (B, eps*100);
+
+%!test
+%! C = gallery("tridiag", 5);
+%! sparse_chol2inv (C, eps*10);
+
+%!test
+%! D = gallery("wathen", 1, 1);
+%! sparse_chol2inv (D, eps*10^4);
+
+*/
+
 DEFUN_DLD (cholupdate, args, nargout,
            doc: /* -*- texinfo -*-
 @deftypefn {} {[@var{R1}, @var{info}] =} cholupdate (@var{R}, @var{u}, @var{op})
--- a/liboctave/array/MatrixType.cc	Tue Aug 09 15:07:30 2016 -0700
+++ b/liboctave/array/MatrixType.cc	Sat Aug 06 19:38:32 2016 +0200
@@ -351,7 +351,9 @@
               else
                 dense = false;
             }
-          else if (upper_band == 0)
+
+          // If a matrix is Banded but also Upper/Lower, set to the latter.
+          if (upper_band == 0)
             typ = MatrixType::Lower;
           else if (lower_band == 0)
             typ = MatrixType::Upper;
@@ -668,7 +670,9 @@
               else
                 dense = false;
             }
-          else if (upper_band == 0)
+
+          // If a matrix is Banded but also Upper/Lower, set to the latter.
+          if (upper_band == 0)
             typ = MatrixType::Lower;
           else if (lower_band == 0)
             typ = MatrixType::Upper;