Mercurial > octave-nkf

comparison liboctave/SparsedbleCHOL.cc @ 5506:b4cfbb0ec8c4

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 2005-10-23 19:09:32 by dbateman]
author dbateman
date Sun, 23 Oct 2005 19:09:33 +0000
parents
children 6b9cec830d72
comparison
equal deleted inserted replaced
5505:17682e3fba2a 5506:b4cfbb0ec8c4
1 /*
2
3 Copyright (C) 2005 David Bateman
4 Copyright (C) 1998-2005 Andy Adler
5
6 Octave is free software; you can redistribute it and/or modify it
7 under the terms of the GNU General Public License as published by the
8 Free Software Foundation; either version 2, or (at your option) any
9 later version.
10
11 Octave is distributed in the hope that it will be useful, but WITHOUT
12 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14 for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with this program; see the file COPYING. If not, write to the
18 Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
19 Boston, MA 02110-1301, USA.
20
21 */
22
23 #ifdef HAVE_CONFIG_H
24 #include <config.h>
25 #endif
26
27 #include "SparsedbleCHOL.h"
28
29 // Instantiate the base CHOL class for the type we need
30 #define OCTAVE_CHOLMOD_TYPE CHOLMOD_REAL
31 #include "sparse-base-chol.h"
32 #include "sparse-base-chol.cc"
33 template class sparse_base_chol <SparseMatrix, double, SparseMatrix>;
34
35 // Compute the inverse of a matrix using the Cholesky factorization.
36 SparseMatrix
37 chol2inv (const SparseMatrix& r)
38 {
39 octave_idx_type r_nr = r.rows ();
40 octave_idx_type r_nc = r.cols ();
41 SparseMatrix retval;
42
43 if (r_nr == r_nc)
44 {
45 SparseType mattype (r);
46 int typ = mattype.type (false);
47 double rcond;
48 octave_idx_type info;
49 SparseMatrix rinv;
50
51 if (typ == SparseType::Upper)
52 {
53 rinv = r.inverse(mattype, info, rcond, true, false);
54 retval = rinv.transpose() * rinv;
55 }
56 else if (typ == SparseType::Lower)
57 {
58 rinv = r.transpose().inverse(mattype, info, rcond, true, false);
59 retval = rinv.transpose() * rinv;
60 }
61 else
62 (*current_liboctave_error_handler)
63 ("spchol2inv requires triangular matrix");
64 }
65 else
66 (*current_liboctave_error_handler) ("spchol2inv requires square matrix");
67
68 return retval;
69 }
70
71 /*
72 ;;; Local Variables: ***
73 ;;; mode: C++ ***
74 ;;; End: ***
75 */