copy from src/DLD-FUNCTIONS/spchol.cc
copy to src/DLD-FUNCTIONS/symbfact.cc
--- a/src/DLD-FUNCTIONS/spchol.cc
+++ b/src/DLD-FUNCTIONS/symbfact.cc
@@ -25,313 +25,19 @@
 #include <config.h>
 #endif
 
+#include "SparseCmplxCHOL.h"
+#include "SparsedbleCHOL.h"
+#include "oct-spparms.h"
+#include "sparse-util.h"
+
+#include "ov-re-sparse.h"
+#include "ov-cx-sparse.h"
 #include "defun-dld.h"
 #include "error.h"
 #include "gripes.h"
 #include "oct-obj.h"
 #include "utils.h"
 
-#include "SparseCmplxCHOL.h"
-#include "SparsedbleCHOL.h"
-#include "ov-re-sparse.h"
-#include "ov-cx-sparse.h"
-#include "oct-spparms.h"
-#include "sparse-util.h"
-
-static octave_value_list
-sparse_chol (const octave_value_list& args, const int nargout, 
-	     const std::string& name, const bool LLt)
-{
-  octave_value_list retval;
-  int nargin = args.length ();
-
-  if (nargin != 1 || nargout > 3)
-    {
-      print_usage ();
-      return retval;
-    }
-
-  octave_value arg = args(0);
-    
-  octave_idx_type nr = arg.rows ();
-  octave_idx_type nc = arg.columns ();
-  bool natural = (nargout != 3);
-
-  int arg_is_empty = empty_arg (name.c_str(), nr, nc);
-
-  if (arg_is_empty < 0)
-    return retval;
-  if (arg_is_empty > 0)
-    return octave_value (Matrix ());
-
-  if (arg.is_real_type ())
-    {
-      SparseMatrix m = arg.sparse_matrix_value ();
-
-      if (! error_state)
-	{
-	  octave_idx_type info;
-	  SparseCHOL fact (m, info, natural);
-	  if (nargout == 3)
-	    retval(2) = fact.Q();
-
-	  if (nargout > 1 || info == 0)
-	    {
-	      retval(1) = fact.P();
-	      if (LLt)
-		retval(0) = fact.L();
-	      else
-		retval(0) = fact.R();
-	    }
-	  else
-	    error ("%s: matrix not positive definite", name.c_str());
-	}
-    }
-  else if (arg.is_complex_type ())
-    {
-      SparseComplexMatrix m = arg.sparse_complex_matrix_value ();
-
-      if (! error_state)
-	{
-	  octave_idx_type info;
-	  SparseComplexCHOL fact (m, info, natural);
-
-	  if (nargout == 3)
-	    retval(2) = fact.Q();
-	  
-	  if (nargout > 1 || info == 0)
-	    {
-	      retval(1) = fact.P();
-	      if (LLt)
-		retval(0) = fact.L();
-	      else
-		retval(0) = fact.R();
-	    }
-	  else
-	    error ("%s: matrix not positive definite", name.c_str());
-	}
-    }
-  else
-    gripe_wrong_type_arg (name.c_str(), arg);
-
-  return retval;
-}
-
-// PKG_ADD: dispatch ("chol", "spchol", "sparse matrix");
-// PKG_ADD: dispatch ("chol", "spchol", "sparse complex matrix");
-// PKG_ADD: dispatch ("chol", "spchol", "sparse bool matrix");
-DEFUN_DLD (spchol, args, nargout,
-  "-*- texinfo -*-\n\
-@deftypefn {Loadable Function} {@var{r} =} spchol (@var{a})\n\
-@deftypefnx {Loadable Function} {[@var{r}, @var{p}] =} spchol (@var{a})\n\
-@deftypefnx {Loadable Function} {[@var{r}, @var{p}, @var{q}] =} spchol (@var{a})\n\
-@cindex Cholesky factorization\n\
-Compute the Cholesky factor, @var{r}, of the symmetric positive definite\n\
-sparse matrix @var{a}, where\n\
-@iftex\n\
-@tex\n\
-$ R^T R = A $.\n\
-@end tex\n\
-@end iftex\n\
-@ifinfo\n\
-\n\
-@example\n\
-r' * r = a.\n\
-@end example\n\
-@end ifinfo\n\
-\n\
-If called with 2 or more outputs @var{p} is the 0 when @var{r} is positive\n\
-definite and @var{p} is a positive integer otherwise.\n\
-\n\
-If called with 3 outputs then a sparsity preserving row/column permutation\n\
-is applied to @var{a} prior to the factorization. That is @var{r}\n\
-is the factorization of @code{@var{a}(@var{q},@var{q})} such that\n\
-@iftex\n\
-@tex\n\
-$ R^T R = Q A Q^T$.\n\
-@end tex\n\
-@end iftex\n\
-@ifinfo\n\
-\n\
-@example\n\
-r' * r = q * a * q'.\n\
-@end example\n\
-@end ifinfo\n\
-\n\
-Note that @code{splchol} factorization is faster and uses less memory.\n\
-@seealso{spcholinv, spchol2inv, splchol}\n\
-@end deftypefn")
-{
-  return sparse_chol (args, nargout, "spchol", false);
-}
-
-// PKG_ADD: dispatch ("lchol", "splchol", "sparse matrix");
-// PKG_ADD: dispatch ("lchol", "splchol", "sparse complex matrix");
-// PKG_ADD: dispatch ("lchol", "splchol", "sparse bool matrix");
-DEFUN_DLD (splchol, args, nargout,
-  "-*- texinfo -*-\n\
-@deftypefn {Loadable Function} {@var{l} =} splchol (@var{a})\n\
-@deftypefnx {Loadable Function} {[@var{l}, @var{p}] =} splchol (@var{a})\n\
-@deftypefnx {Loadable Function} {[@var{l}, @var{p}, @var{q}] =} splchol (@var{a})\n\
-@cindex Cholesky factorization\n\
-Compute the Cholesky factor, @var{l}, of the symmetric positive definite\n\
-sparse matrix @var{a}, where\n\
-@iftex\n\
-@tex\n\
-$ L L^T = A $.\n\
-@end tex\n\
-@end iftex\n\
-@ifinfo\n\
-\n\
-@example\n\
-l * l' = a.\n\
-@end example\n\
-@end ifinfo\n\
-\n\
-If called with 2 or more outputs @var{p} is the 0 when @var{l} is positive\n\
-definite and @var{l} is a positive integer otherwise.\n\
-\n\
-If called with 3 outputs that a sparsity preserving row/column permutation\n\
-is applied to @var{a} prior to the factorization. That is @var{l}\n\
-is the factorization of @code{@var{a}(@var{q},@var{q})} such that\n\
-@iftex\n\
-@tex\n\
-$ L R^T = A (Q, Q)$.\n\
-@end tex\n\
-@end iftex\n\
-@ifinfo\n\
-\n\
-@example\n\
-r * r' = a (q, q).\n\
-@end example\n\
-@end ifinfo\n\
-\n\
-Note that @code{splchol} factorization is faster and uses less memory\n\
-than @code{spchol}. @code{splchol(@var{a})} is equivalent to\n\
-@code{spchol(@var{a})'}.\n\
-@seealso{spcholinv, spchol2inv, splchol}\n\
-@end deftypefn")
-{
-  return sparse_chol (args, nargout, "splchol", true);
-}
-
-// PKG_ADD: dispatch ("cholinv", "spcholinv", "sparse matrix");
-// PKG_ADD: dispatch ("cholinv", "spcholinv", "sparse complex matrix");
-// PKG_ADD: dispatch ("cholinv", "spcholinv", "sparse bool matrix");
-DEFUN_DLD (spcholinv, args, ,
-  "-*- texinfo -*-\n\
-@deftypefn {Loadable Function} {} spcholinv (@var{a})\n\
-Use the Cholesky factorization to compute the inverse of the\n\
-sparse symmetric positive definite matrix @var{a}.\n\
-@seealso{spchol, spchol2inv}\n\
-@end deftypefn")
-{
-  octave_value retval;
-
-  int nargin = args.length ();
-
-  if (nargin == 1)
-    {
-      octave_value arg = args(0);
-    
-      octave_idx_type nr = arg.rows ();
-      octave_idx_type nc = arg.columns ();
-
-      if (nr == 0 || nc == 0)
-	retval = Matrix ();
-      else
-	{
-	  if (arg.is_real_type ())
-	    {
-	      SparseMatrix m = arg.sparse_matrix_value ();
-
-	      if (! error_state)
-		{
-		  octave_idx_type info;
-		  SparseCHOL chol (m, info);
-		  if (info == 0)
-		    retval = chol.inverse ();
-		  else
-		    error ("spcholinv: matrix not positive definite");
-		}
-	    }
-	  else if (arg.is_complex_type ())
-	    {
-	      SparseComplexMatrix m = arg.sparse_complex_matrix_value ();
-
-	      if (! error_state)
-		{
-		  octave_idx_type info;
-		  SparseComplexCHOL chol (m, info);
-		  if (info == 0)
-		    retval = chol.inverse ();
-		  else
-		    error ("spcholinv: matrix not positive definite");
-		}
-	    }
-	  else
-	    gripe_wrong_type_arg ("spcholinv", arg);
-	}
-    }
-  else
-    print_usage ();
-
-  return retval;
-}
-
-// PKG_ADD: dispatch ("chol2inv", "spchol2inv", "sparse matrix");
-// PKG_ADD: dispatch ("chol2inv", "spchol2inv", "sparse complex matrix");
-// PKG_ADD: dispatch ("chol2inv", "spchol2inv", "sparse bool matrix");
-DEFUN_DLD (spchol2inv, args, ,
-  "-*- texinfo -*-\n\
-@deftypefn {Loadable Function} {} spchol2inv (@var{u})\n\
-Invert a sparse symmetric, positive definite square matrix from its\n\
-Cholesky decomposition, @var{u}.  Note that @var{u} should be an\n\
-upper-triangular matrix with positive diagonal elements.\n\
-@code{chol2inv (@var{u})} provides @code{inv (@var{u}'*@var{u})} but\n\
-it is much faster than using @code{inv}.\n\
-@seealso{spchol, spcholinv}\n\
-@end deftypefn")
-{
-  octave_value retval;
-
-  int nargin = args.length ();
-
-  if (nargin == 1)
-    {
-      octave_value arg = args(0);
-    
-      octave_idx_type nr = arg.rows ();
-      octave_idx_type nc = arg.columns ();
-
-      if (nr == 0 || nc == 0)
-	retval = Matrix ();
-      else
-	{
-	  if (arg.is_real_type ())
-	    {
-	      SparseMatrix r = arg.sparse_matrix_value ();
-
-	      if (! error_state)
-		retval = chol2inv (r);
-	    }
-	  else if (arg.is_complex_type ())
-	    {
-	      SparseComplexMatrix r = arg.sparse_complex_matrix_value ();
-
-	      if (! error_state)
-		retval = chol2inv (r);
-	    }
-	  else
-	    gripe_wrong_type_arg ("spchol2inv", arg);
-	}
-    }
-  else
-    print_usage ();
-
-  return retval;
-}
-
 DEFUN_DLD (symbfact, args, nargout,
     "-*- texinfo -*-\n\
 @deftypefn {Loadable Function} {[@var{count}, @var{h}, @var{parent}, @var{post}, @var{r}] =} symbfact (@var{s}, @var{typ}, @var{mode})\n\