[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[CHANGESET] Add the amd function
From: |
David Bateman |
Subject: |
[CHANGESET] Add the amd function |
Date: |
Fri, 21 Mar 2008 13:33:04 +0100 |
User-agent: |
Thunderbird 2.0.0.12 (X11/20080306) |
As the amd code is already used with umfpack, which Octave uses for the
sparse LU function, and Tim Davis supplies a mex file for the equivalent
of the amd function in matlab under the GPL, it is trivial to write an
Octave version of the same code.. Attached is a matlab compatible
version of the amd function for Octave.
D.
--
David Bateman address@hidden
Motorola Labs - Paris +33 1 69 35 48 04 (Ph)
Parc Les Algorithmes, Commune de St Aubin +33 6 72 01 06 33 (Mob)
91193 Gif-Sur-Yvette FRANCE +33 1 69 35 77 01 (Fax)
The information contained in this communication has been classified as:
[x] General Business Information
[ ] Motorola Internal Use Only
[ ] Motorola Confidential Proprietary
# HG changeset patch
# User David Bateman <address@hidden>
# Date 1206102654 -3600
# Node ID 5429d6a96b4e39a1db5fbaf3692e059bed23be28
# Parent 6e930aebeebffd2fcd91fe94d7c47269df38d76a
Add the amd function
diff --git a/ChangeLog b/ChangeLog
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,7 @@ 2008-03-18 David Bateman <address@hidden
+2008-03-21 David Bateman <address@hidden>
+
+ * configure.in (HAVE_AMD): Complete test for presence of amd.
+
2008-03-18 David Bateman <address@hidden>
* configure.in (AC_CHECK_FUNCS): Also check lgamma_r.
diff --git a/configure.in b/configure.in
--- a/configure.in
+++ b/configure.in
@@ -734,7 +734,26 @@ AC_SUBST(LAPACK_DIR)
# Check for AMD library
AMD_LIBS=
AC_SUBST(AMD_LIBS)
-AC_CHECK_LIB(amd, amd_postorder, [AMD_LIBS="-lamd";
with_amd=yes],[with_amd=no])
+
+AC_ARG_WITH(amd,
+ [AS_HELP_STRING([--without-amd],
+ [don't use AMD, disable some sparse functionality])],
+ with_amd=$withval, with_amd=yes)
+
+warn_amd="AMD not found. This will result in some lack of functionality for
sparse matrices."
+if test "$with_amd" = yes; then
+ with_amd=no
+ AC_CHECK_HEADERS([suitesparse/amd.h ufsparse/amd.h amd/amd.h amd.h], [
+ AC_CHECK_LIB(amd, amd_postorder, [AMD_LIBS="-lamd"; with_amd=yes])
+ if test "$with_amd" = yes; then
+ AC_DEFINE(HAVE_AMD, 1, [Define if the AMD library is used.])
+ warn_amd=
+ fi
+ break])
+fi
+if test -n "$warn_amd"; then
+ AC_MSG_WARN($warn_amd)
+fi
# Check for CAMD library
CAMD_LIBS=
@@ -1948,6 +1967,11 @@ if test -n "$warn_umfpack"; then
warn_msg_printed=true
fi
+if test -n "$warn_amd"; then
+ AC_MSG_WARN($warn_amd)
+ warn_msg_printed=true
+fi
+
if test -n "$warn_colamd"; then
AC_MSG_WARN($warn_colamd)
warn_msg_printed=true
diff --git a/doc/ChangeLog b/doc/ChangeLog
--- a/doc/ChangeLog
+++ b/doc/ChangeLog
@@ -1,3 +1,7 @@ 2008-03-19 Michael D. Godfrey <godfrey
+2008-03-21 David Bateman <address@hidden>
+
+ * interpreter/sparse.txi: Document amd function.
+
2008-03-19 Michael D. Godfrey <address@hidden>
* interpreter/plot.txi: Reorder symbol character table.
diff --git a/doc/interpreter/sparse.txi b/doc/interpreter/sparse.txi
--- a/doc/interpreter/sparse.txi
+++ b/doc/interpreter/sparse.txi
@@ -469,7 +469,7 @@ used.
@c @dfn{treelayout}
@item Sparse matrix reordering:
- @dfn{ccolamd}, @dfn{colamd}, @dfn{colperm}, @dfn{csymamd},
+ @dfn{amd}, @dfn{ccolamd}, @dfn{colamd}, @dfn{colperm}, @dfn{csymamd},
@dfn{dmperm}, @dfn{symamd}, @dfn{randperm}, @dfn{symrcm}
@item Linear algebra:
@@ -621,7 +621,7 @@ Several functions are available to reord
Several functions are available to reorder depending on the type of the
matrix to be factorized. If the matrix is symmetric positive-definite,
then @dfn{symamd} or @dfn{csymamd} should be used. Otherwise
address@hidden or @dfn{ccolamd} should be used. For completeness
address@hidden, @dfn{colamd} or @dfn{ccolamd} should be used. For completeness
the reordering functions @dfn{colperm} and @dfn{randperm} are
also available.
@@ -706,6 +706,8 @@ Finally, Octave implicitly reorders the
Finally, Octave implicitly reorders the matrix when using the div (/)
and ldiv (\) operators, and so no the user does not need to explicitly
reorder the matrix to maximize performance.
+
address@hidden(amd)
@DOCSTRING(ccolamd)
diff --git a/liboctave/ChangeLog b/liboctave/ChangeLog
--- a/liboctave/ChangeLog
+++ b/liboctave/ChangeLog
@@ -1,3 +1,7 @@ 2008-03-20 David Bateman <address@hidden
+2008-03-21 David Bateman <address@hidden>
+
+ * oct-sparse.h: Add headers for amd.h.
+
2008-03-20 David Bateman <address@hidden>
* Array.cc (Array<T> Array<T>::diag (octave_idx_type) const): New
diff --git a/liboctave/oct-sparse.h b/liboctave/oct-sparse.h
--- a/liboctave/oct-sparse.h
+++ b/liboctave/oct-sparse.h
@@ -25,6 +25,16 @@ along with Octave; see the file COPYING.
#ifdef HAVE_CONFIG_H
#include <config.h>
+#endif
+
+#if defined (HAVE_SUITESPARSE_AMD_H)
+#include <suitesparse/amd.h>
+#elif defined (HAVE_UFSPARSE_AMD_H)
+#include <ufsparse/amd.h>
+#elif defined (HAVE_AMD_AMD_H)
+#include <amd/amd.h>
+#elif defined (HAVE_AMD_H)
+#include <amd.h>
#endif
#if defined (HAVE_SUITESPARSE_UMFPACK_H)
diff --git a/src/ChangeLog b/src/ChangeLog
--- a/src/ChangeLog
+++ b/src/ChangeLog
@@ -1,3 +1,8 @@ 2008-03-20 David Bateman <address@hidden
+2008-03-21 David Bateman <address@hidden>
+
+ * DLD-FUNCTIONS/amd.cc: New file.
+ * Makefile.in (DLD_XSRC): Add amd.cc.
+
2008-03-20 David Bateman <address@hidden>
* Cell.cc (Cell Cell::diag (void) const): delete.
diff --git a/src/DLD-FUNCTIONS/amd.cc b/src/DLD-FUNCTIONS/amd.cc
new file mode 100644
--- /dev/null
+++ b/src/DLD-FUNCTIONS/amd.cc
@@ -0,0 +1,212 @@
+/*
+
+Copyright (C) 2008 David Bateman
+
+This file is part of Octave.
+
+Octave is free software; you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation; either version 3 of the License, or (at your
+option) any later version.
+
+Octave is distributed in the hope that it will be useful, but WITHOUT
+ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+for more details.
+
+You should have received a copy of the GNU General Public License
+along with Octave; see the file COPYING. If not, see
+<http://www.gnu.org/licenses/>.
+
+*/
+
+// This is the octave interface to amd, which bore the copyright given
+// in the help of the functions.
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include <cstdlib>
+
+#include <string>
+#include <vector>
+
+#include "ov.h"
+#include "defun-dld.h"
+#include "pager.h"
+#include "ov-re-mat.h"
+
+#include "ov-re-sparse.h"
+#include "ov-cx-sparse.h"
+#include "oct-map.h"
+
+#include "oct-sparse.h"
+
+#ifdef IDX_TYPE_LONG
+#define AMD_NAME(name) amd_l ## name
+#else
+#define AMD_NAME(name) amd ## name
+#endif
+
+DEFUN_DLD (amd, args, nargout,
+ "-*- texinfo -*-\n\
address@hidden {Loadable Function} address@hidden =} amd (@var{s})\n\
address@hidden {Loadable Function} address@hidden =} amd (@var{s},
@var{opts})\n\
+\n\
+Returns the approximate minimum degree permutation of a matrix. This\n\
+permutation such that the Cholesky factorization of @address@hidden
(@var{p},\n\
address@hidden)} tends to be sparser than the Cholesky factorization of
@var{s}\n\
+itself. @code{amd} is typically faster than @code{symamd} but serves a\n\
+similar purpose.\n\
+\n\
+The optional parameter @var{opts} is a structure that controls the\n\
+behavior of @code{amd}. The fields of these structure are\n\
+\n\
address@hidden @asis\n\
address@hidden opts.dense\n\
+Determines what @code{amd} considers to be a dense row or column of the\n\
+input matrix. Rows or columns with more that @code{max(16, (dense *\n\
+sqrt (@var{n})} entries, where @var{n} is the order of the matrix @var{s},\n\
+are igorned by @code{amd} during the calculation of the permutation\n\
+The value of dense must be a positive scalar and its default value is 10.0\n\
+\n\
address@hidden opts.aggressive\n\
+If this value is a non zero scalar, then @code{amd} performs agressive\n\
+absorption. The default is not to perform agressive absorption.\n\
address@hidden table\n\
+\n\
+The author of the code itself is Timothy A. Davis (davis@@cise.ufl.edu),\n\
+University of Florida (see
@url{http://www.cise.ufl.edu/research/sparse/amd}).\n\
address@hidden, colamd}\n\
address@hidden deftypefn")
+{
+ octave_value_list retval;
+
+#ifdef HAVE_AMD
+ int nargin = args.length ();
+
+ if (nargin < 1 || nargin > 2)
+ print_usage ();
+ else
+ {
+ octave_idx_type n_row, n_col, nnz;
+ const octave_idx_type *ridx, *cidx;
+ SparseMatrix sm;
+ SparseComplexMatrix scm;
+
+ if (args(0).is_sparse_type ())
+ {
+ if (args(0).is_complex_type ())
+ {
+ scm = args(0).sparse_complex_matrix_value ();
+ n_row = scm.rows ();
+ n_col = scm.cols ();
+ nnz = scm.nzmax ();
+ ridx = scm.xridx ();
+ cidx = scm.xcidx ();
+ }
+ else
+ {
+ sm = args(0).sparse_matrix_value ();
+ n_row = sm.rows ();
+ n_col = sm.cols ();
+ nnz = sm.nzmax ();
+ ridx = sm.xridx ();
+ cidx = sm.xcidx ();
+ }
+ }
+ else
+ {
+ if (args(0).is_complex_type ())
+ sm = SparseMatrix (real (args(0).complex_matrix_value ()));
+ else
+ sm = SparseMatrix (args(0).matrix_value ());
+
+ n_row = sm.rows ();
+ n_col = sm.cols ();
+ nnz = sm.nzmax ();
+ ridx = sm.xridx ();
+ cidx = sm.xcidx ();
+ }
+
+ if (!error_state && n_row != n_col)
+ error ("amd: input matrix must be square");
+
+ if (!error_state)
+ {
+ OCTAVE_LOCAL_BUFFER (double, Control, AMD_CONTROL);
+ AMD_NAME (_defaults) (Control) ;
+ if (nargin > 1)
+ {
+ Octave_map arg1 = args(1).map_value ();
+
+ if (!error_state)
+ {
+ if (arg1.contains ("dense"))
+ {
+ Cell c = arg1.contents ("dense");
+ if (c.length() == 1)
+ Control[AMD_DENSE] = c.elem(0).double_value ();
+ else
+ error ("amd: invalid options structure");
+ }
+ if (arg1.contains ("aggressive"))
+ {
+ Cell c = arg1.contents ("aggressive");
+ if (c.length() == 1)
+ Control[AMD_AGGRESSIVE] = c.elem(0).double_value ();
+ else
+ error ("amd: invalid options structure");
+ }
+ }
+ }
+
+ if (!error_state)
+ {
+ OCTAVE_LOCAL_BUFFER (octave_idx_type, P, n_col);
+ Matrix xinfo (AMD_INFO, 1);
+ double *Info = xinfo.fortran_vec ();
+
+ // FIXME: How can we manage the memory allocation of amd in
+ // a cleaner manner?
+ amd_malloc = malloc;
+ amd_free = free;
+ amd_calloc = calloc;
+ amd_realloc = realloc;
+ amd_printf = printf;
+
+ octave_idx_type result = AMD_NAME (_order) (n_col, cidx, ridx, P,
+ Control, Info);
+
+ switch (result)
+ {
+ case AMD_OUT_OF_MEMORY:
+ error ("amd: out of memory");
+ break;
+ case AMD_INVALID:
+ error ("amd: input matrix is corrupted");
+ break;
+ default:
+ {
+ if (nargout > 1)
+ retval(1) = xinfo;
+
+ Matrix Pout (1, n_col);
+ for (octave_idx_type i = 0; i < n_col; i++)
+ Pout.xelem (i) = P[i] + 1;
+
+ retval (0) = Pout;
+ }
+ }
+ }
+ }
+ }
+#else
+
+ error ("amd: not available in this version of Octave");
+
+#endif
+
+ return retval;
+}
diff --git a/src/Makefile.in b/src/Makefile.in
--- a/src/Makefile.in
+++ b/src/Makefile.in
@@ -62,8 +62,8 @@ OPT_IN := $(addprefix ../liboctave/, $(a
OPT_IN := $(addprefix ../liboctave/, $(addsuffix .in, $(OPT_BASE)))
OPT_INC := $(addprefix ../liboctave/, $(addsuffix .h, $(OPT_BASE)))
-DLD_XSRC := balance.cc besselj.cc betainc.cc bsxfun.cc cellfun.cc chol.cc \
- ccolamd.cc colamd.cc colloc.cc conv2.cc convhulln.cc daspk.cc \
+DLD_XSRC := amd.cc balance.cc besselj.cc betainc.cc bsxfun.cc cellfun.cc \
+ chol.cc ccolamd.cc colamd.cc colloc.cc conv2.cc convhulln.cc daspk.cc \
dasrt.cc dassl.cc det.cc dispatch.cc dlmread.cc dmperm.cc eig.cc \
expm.cc fft.cc fft2.cc fftn.cc fftw.cc filter.cc find.cc fsolve.cc \
gammainc.cc gcd.cc getgrent.cc getpwent.cc getrusage.cc \
- [CHANGESET] Add the amd function,
David Bateman <=