Mercurial > octave
annotate liboctave/numeric/chol.cc @ 21270:230e186e292d
make building without qrupdate work again
* dMatrix.h (Matrix::hermitian): New function.
* fMatrix.h (FloatMatrix::hermitian): New function.
* liboctave/numeric/chol.cc: Fix function declarations, definition of
zero, and names of imag and conj functions.
author | John W. Eaton <jwe@octave.org> |
---|---|
date | Tue, 16 Feb 2016 11:13:55 -0500 |
parents | 3c8a3d35661a |
children | d3b265a83adc |
rev | line source |
---|---|
457 | 1 /* |
2 | |
19697
4197fc428c7d
maint: Update copyright notices for 2015.
John W. Eaton <jwe@octave.org>
parents:
17769
diff
changeset
|
3 Copyright (C) 1994-2015 John W. Eaton |
11523 | 4 Copyright (C) 2008-2009 Jaroslav Hajek |
457 | 5 |
6 This file is part of Octave. | |
7 | |
8 Octave is free software; you can redistribute it and/or modify it | |
9 under the terms of the GNU General Public License as published by the | |
7016 | 10 Free Software Foundation; either version 3 of the License, or (at your |
11 option) any later version. | |
457 | 12 |
13 Octave is distributed in the hope that it will be useful, but WITHOUT | |
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |
16 for more details. | |
17 | |
18 You should have received a copy of the GNU General Public License | |
7016 | 19 along with Octave; see the file COPYING. If not, see |
20 <http://www.gnu.org/licenses/>. | |
457 | 21 |
22 */ | |
23 | |
24 #ifdef HAVE_CONFIG_H | |
21202
f7121e111991
maint: indent #ifdef blocks in liboctave and src directories.
Rik <rik@octave.org>
parents:
21136
diff
changeset
|
25 # include <config.h> |
457 | 26 #endif |
27 | |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
28 #include <vector> |
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
29 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
30 #include "CColVector.h" |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
31 #include "CMatrix.h" |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
32 #include "CRowVector.h" |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
33 #include "chol.h" |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
34 #include "dColVector.h" |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
35 #include "dMatrix.h" |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
36 #include "dRowVector.h" |
1847 | 37 #include "f77-fcn.h" |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
38 #include "fCColVector.h" |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
39 #include "fCMatrix.h" |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
40 #include "fCRowVector.h" |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
41 #include "fColVector.h" |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
42 #include "fMatrix.h" |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
43 #include "fRowVector.h" |
457 | 44 #include "lo-error.h" |
8377
25bc2d31e1bf
improve OCTAVE_LOCAL_BUFFER
Jaroslav Hajek <highegg@gmail.com>
parents:
7725
diff
changeset
|
45 #include "oct-locbuf.h" |
9862
c0aeedd8fb86
improve chol Matlab compatibility
Jaroslav Hajek <highegg@gmail.com>
parents:
8562
diff
changeset
|
46 #include "oct-norm.h" |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
47 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
48 #if ! defined (HAVE_QRUPDATE) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
49 # include "CmplxQR.h" |
21202
f7121e111991
maint: indent #ifdef blocks in liboctave and src directories.
Rik <rik@octave.org>
parents:
21136
diff
changeset
|
50 # include "dbleQR.h" |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
51 # include "fCmplxQR.h" |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
52 # include "floatQR.h" |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
53 #endif |
457 | 54 |
55 extern "C" | |
56 { | |
4552 | 57 F77_RET_T |
11518 | 58 F77_FUNC (dpotrf, DPOTRF) (F77_CONST_CHAR_ARG_DECL, |
59 const octave_idx_type&, double*, | |
60 const octave_idx_type&, octave_idx_type& | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
61 F77_CHAR_ARG_LEN_DECL); |
5340 | 62 |
63 F77_RET_T | |
11518 | 64 F77_FUNC (dpotri, DPOTRI) (F77_CONST_CHAR_ARG_DECL, |
65 const octave_idx_type&, double*, | |
66 const octave_idx_type&, octave_idx_type& | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
67 F77_CHAR_ARG_LEN_DECL); |
6486 | 68 |
69 F77_RET_T | |
11518 | 70 F77_FUNC (dpocon, DPOCON) (F77_CONST_CHAR_ARG_DECL, |
71 const octave_idx_type&, double*, | |
72 const octave_idx_type&, const double&, | |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11570
diff
changeset
|
73 double&, double*, octave_idx_type*, |
11518 | 74 octave_idx_type& |
75 F77_CHAR_ARG_LEN_DECL); | |
8547
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
76 #ifdef HAVE_QRUPDATE |
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
77 |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
78 F77_RET_T |
11495 | 79 F77_FUNC (dch1up, DCH1UP) (const octave_idx_type&, double*, |
80 const octave_idx_type&, double*, double*); | |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
81 |
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
82 F77_RET_T |
11495 | 83 F77_FUNC (dch1dn, DCH1DN) (const octave_idx_type&, double*, |
84 const octave_idx_type&, double*, double*, | |
85 octave_idx_type&); | |
86 | |
87 F77_RET_T | |
88 F77_FUNC (dchinx, DCHINX) (const octave_idx_type&, double*, | |
89 const octave_idx_type&, const octave_idx_type&, | |
8547
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
90 double*, double*, octave_idx_type&); |
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
91 |
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
92 F77_RET_T |
11495 | 93 F77_FUNC (dchdex, DCHDEX) (const octave_idx_type&, double*, |
94 const octave_idx_type&, const octave_idx_type&, | |
95 double*); | |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
96 |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
97 F77_RET_T |
11495 | 98 F77_FUNC (dchshx, DCHSHX) (const octave_idx_type&, double*, |
99 const octave_idx_type&, const octave_idx_type&, | |
8547
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
100 const octave_idx_type&, double*); |
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
101 #endif |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
102 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
103 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
104 F77_FUNC (spotrf, SPOTRF) (F77_CONST_CHAR_ARG_DECL, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
105 const octave_idx_type&, float*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
106 const octave_idx_type&, octave_idx_type& |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
107 F77_CHAR_ARG_LEN_DECL); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
108 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
109 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
110 F77_FUNC (spotri, SPOTRI) (F77_CONST_CHAR_ARG_DECL, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
111 const octave_idx_type&, float*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
112 const octave_idx_type&, octave_idx_type& |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
113 F77_CHAR_ARG_LEN_DECL); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
114 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
115 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
116 F77_FUNC (spocon, SPOCON) (F77_CONST_CHAR_ARG_DECL, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
117 const octave_idx_type&, float*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
118 const octave_idx_type&, const float&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
119 float&, float*, octave_idx_type*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
120 octave_idx_type& |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
121 F77_CHAR_ARG_LEN_DECL); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
122 #ifdef HAVE_QRUPDATE |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
123 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
124 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
125 F77_FUNC (sch1up, SCH1UP) (const octave_idx_type&, float*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
126 const octave_idx_type&, float*, float*); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
127 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
128 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
129 F77_FUNC (sch1dn, SCH1DN) (const octave_idx_type&, float*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
130 const octave_idx_type&, float*, float*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
131 octave_idx_type&); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
132 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
133 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
134 F77_FUNC (schinx, SCHINX) (const octave_idx_type&, float*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
135 const octave_idx_type&, const octave_idx_type&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
136 float*, float*, octave_idx_type&); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
137 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
138 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
139 F77_FUNC (schdex, SCHDEX) (const octave_idx_type&, float*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
140 const octave_idx_type&, const octave_idx_type&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
141 float*); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
142 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
143 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
144 F77_FUNC (schshx, SCHSHX) (const octave_idx_type&, float*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
145 const octave_idx_type&, const octave_idx_type&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
146 const octave_idx_type&, float*); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
147 #endif |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
148 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
149 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
150 F77_FUNC (zpotrf, ZPOTRF) (F77_CONST_CHAR_ARG_DECL, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
151 const octave_idx_type&, Complex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
152 const octave_idx_type&, octave_idx_type& |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
153 F77_CHAR_ARG_LEN_DECL); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
154 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
155 F77_FUNC (zpotri, ZPOTRI) (F77_CONST_CHAR_ARG_DECL, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
156 const octave_idx_type&, Complex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
157 const octave_idx_type&, octave_idx_type& |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
158 F77_CHAR_ARG_LEN_DECL); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
159 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
160 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
161 F77_FUNC (zpocon, ZPOCON) (F77_CONST_CHAR_ARG_DECL, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
162 const octave_idx_type&, Complex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
163 const octave_idx_type&, const double&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
164 double&, Complex*, double*, octave_idx_type& |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
165 F77_CHAR_ARG_LEN_DECL); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
166 #ifdef HAVE_QRUPDATE |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
167 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
168 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
169 F77_FUNC (zch1up, ZCH1UP) (const octave_idx_type&, Complex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
170 const octave_idx_type&, Complex*, double*); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
171 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
172 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
173 F77_FUNC (zch1dn, ZCH1DN) (const octave_idx_type&, Complex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
174 const octave_idx_type&, Complex*, double*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
175 octave_idx_type&); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
176 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
177 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
178 F77_FUNC (zchinx, ZCHINX) (const octave_idx_type&, Complex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
179 const octave_idx_type&, const octave_idx_type&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
180 Complex*, double*, octave_idx_type&); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
181 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
182 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
183 F77_FUNC (zchdex, ZCHDEX) (const octave_idx_type&, Complex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
184 const octave_idx_type&, const octave_idx_type&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
185 double*); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
186 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
187 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
188 F77_FUNC (zchshx, ZCHSHX) (const octave_idx_type&, Complex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
189 const octave_idx_type&, const octave_idx_type&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
190 const octave_idx_type&, Complex*, double*); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
191 #endif |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
192 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
193 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
194 F77_FUNC (cpotrf, CPOTRF) (F77_CONST_CHAR_ARG_DECL, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
195 const octave_idx_type&, FloatComplex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
196 const octave_idx_type&, octave_idx_type& |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
197 F77_CHAR_ARG_LEN_DECL); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
198 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
199 F77_FUNC (cpotri, CPOTRI) (F77_CONST_CHAR_ARG_DECL, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
200 const octave_idx_type&, FloatComplex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
201 const octave_idx_type&, octave_idx_type& |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
202 F77_CHAR_ARG_LEN_DECL); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
203 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
204 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
205 F77_FUNC (cpocon, CPOCON) (F77_CONST_CHAR_ARG_DECL, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
206 const octave_idx_type&, FloatComplex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
207 const octave_idx_type&, const float&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
208 float&, FloatComplex*, float*, octave_idx_type& |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
209 F77_CHAR_ARG_LEN_DECL); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
210 #ifdef HAVE_QRUPDATE |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
211 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
212 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
213 F77_FUNC (cch1up, CCH1UP) (const octave_idx_type&, FloatComplex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
214 const octave_idx_type&, FloatComplex*, float*); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
215 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
216 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
217 F77_FUNC (cch1dn, CCH1DN) (const octave_idx_type&, FloatComplex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
218 const octave_idx_type&, FloatComplex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
219 float*, octave_idx_type&); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
220 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
221 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
222 F77_FUNC (cchinx, CCHINX) (const octave_idx_type&, FloatComplex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
223 const octave_idx_type&, const octave_idx_type&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
224 FloatComplex*, float*, octave_idx_type&); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
225 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
226 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
227 F77_FUNC (cchdex, CCHDEX) (const octave_idx_type&, FloatComplex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
228 const octave_idx_type&, const octave_idx_type&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
229 float*); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
230 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
231 F77_RET_T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
232 F77_FUNC (cchshx, CCHSHX) (const octave_idx_type&, FloatComplex*, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
233 const octave_idx_type&, const octave_idx_type&, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
234 const octave_idx_type&, FloatComplex*, float*); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
235 #endif |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
236 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
237 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
238 static Matrix |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
239 chol2inv_internal (const Matrix& r, bool is_upper = true) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
240 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
241 Matrix retval; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
242 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
243 octave_idx_type r_nr = r.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
244 octave_idx_type r_nc = r.cols (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
245 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
246 if (r_nr != r_nc) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
247 (*current_liboctave_error_handler) ("chol2inv requires square matrix"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
248 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
249 octave_idx_type n = r_nc; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
250 octave_idx_type info = 0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
251 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
252 Matrix tmp = r; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
253 double *v = tmp.fortran_vec (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
254 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
255 if (info == 0) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
256 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
257 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
258 F77_XFCN (dpotri, DPOTRI, (F77_CONST_CHAR_ARG2 ("U", 1), n, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
259 v, n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
260 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
261 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
262 F77_XFCN (dpotri, DPOTRI, (F77_CONST_CHAR_ARG2 ("L", 1), n, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
263 v, n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
264 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
265 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
266 // If someone thinks of a more graceful way of doing this (or |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
267 // faster for that matter :-)), please let me know! |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
268 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
269 if (n > 1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
270 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
271 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
272 for (octave_idx_type j = 0; j < r_nc; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
273 for (octave_idx_type i = j+1; i < r_nr; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
274 tmp.xelem (i, j) = tmp.xelem (j, i); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
275 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
276 for (octave_idx_type j = 0; j < r_nc; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
277 for (octave_idx_type i = j+1; i < r_nr; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
278 tmp.xelem (j, i) = tmp.xelem (i, j); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
279 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
280 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
281 retval = tmp; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
282 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
283 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
284 return retval; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
285 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
286 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
287 static FloatMatrix |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
288 chol2inv_internal (const FloatMatrix& r, bool is_upper = true) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
289 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
290 FloatMatrix retval; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
291 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
292 octave_idx_type r_nr = r.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
293 octave_idx_type r_nc = r.cols (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
294 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
295 if (r_nr != r_nc) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
296 (*current_liboctave_error_handler) ("chol2inv requires square matrix"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
297 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
298 octave_idx_type n = r_nc; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
299 octave_idx_type info = 0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
300 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
301 FloatMatrix tmp = r; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
302 float *v = tmp.fortran_vec (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
303 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
304 if (info == 0) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
305 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
306 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
307 F77_XFCN (spotri, SPOTRI, (F77_CONST_CHAR_ARG2 ("U", 1), n, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
308 v, n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
309 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
310 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
311 F77_XFCN (spotri, SPOTRI, (F77_CONST_CHAR_ARG2 ("L", 1), n, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
312 v, n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
313 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
314 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
315 // If someone thinks of a more graceful way of doing this (or |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
316 // faster for that matter :-)), please let me know! |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
317 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
318 if (n > 1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
319 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
320 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
321 for (octave_idx_type j = 0; j < r_nc; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
322 for (octave_idx_type i = j+1; i < r_nr; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
323 tmp.xelem (i, j) = tmp.xelem (j, i); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
324 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
325 for (octave_idx_type j = 0; j < r_nc; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
326 for (octave_idx_type i = j+1; i < r_nr; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
327 tmp.xelem (j, i) = tmp.xelem (i, j); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
328 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
329 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
330 retval = tmp; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
331 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
332 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
333 return retval; |
457 | 334 } |
335 | |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
336 static ComplexMatrix |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
337 chol2inv_internal (const ComplexMatrix& r, bool is_upper = true) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
338 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
339 ComplexMatrix retval; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
340 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
341 octave_idx_type r_nr = r.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
342 octave_idx_type r_nc = r.cols (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
343 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
344 if (r_nr != r_nc) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
345 (*current_liboctave_error_handler) ("chol2inv requires square matrix"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
346 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
347 octave_idx_type n = r_nc; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
348 octave_idx_type info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
349 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
350 ComplexMatrix tmp = r; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
351 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
352 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
353 F77_XFCN (zpotri, ZPOTRI, (F77_CONST_CHAR_ARG2 ("U", 1), n, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
354 tmp.fortran_vec (), n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
355 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
356 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
357 F77_XFCN (zpotri, ZPOTRI, (F77_CONST_CHAR_ARG2 ("L", 1), n, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
358 tmp.fortran_vec (), n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
359 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
360 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
361 // If someone thinks of a more graceful way of doing this (or |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
362 // faster for that matter :-)), please let me know! |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
363 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
364 if (n > 1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
365 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
366 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
367 for (octave_idx_type j = 0; j < r_nc; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
368 for (octave_idx_type i = j+1; i < r_nr; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
369 tmp.xelem (i, j) = tmp.xelem (j, i); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
370 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
371 for (octave_idx_type j = 0; j < r_nc; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
372 for (octave_idx_type i = j+1; i < r_nr; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
373 tmp.xelem (j, i) = tmp.xelem (i, j); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
374 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
375 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
376 retval = tmp; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
377 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
378 return retval; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
379 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
380 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
381 static FloatComplexMatrix |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
382 chol2inv_internal (const FloatComplexMatrix& r, bool is_upper = true) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
383 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
384 FloatComplexMatrix retval; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
385 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
386 octave_idx_type r_nr = r.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
387 octave_idx_type r_nc = r.cols (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
388 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
389 if (r_nr != r_nc) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
390 (*current_liboctave_error_handler) ("chol2inv requires square matrix"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
391 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
392 octave_idx_type n = r_nc; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
393 octave_idx_type info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
394 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
395 FloatComplexMatrix tmp = r; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
396 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
397 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
398 F77_XFCN (cpotri, CPOTRI, (F77_CONST_CHAR_ARG2 ("U", 1), n, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
399 tmp.fortran_vec (), n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
400 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
401 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
402 F77_XFCN (cpotri, CPOTRI, (F77_CONST_CHAR_ARG2 ("L", 1), n, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
403 tmp.fortran_vec (), n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
404 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
405 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
406 // If someone thinks of a more graceful way of doing this (or |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
407 // faster for that matter :-)), please let me know! |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
408 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
409 if (n > 1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
410 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
411 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
412 for (octave_idx_type j = 0; j < r_nc; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
413 for (octave_idx_type i = j+1; i < r_nr; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
414 tmp.xelem (i, j) = tmp.xelem (j, i); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
415 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
416 for (octave_idx_type j = 0; j < r_nc; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
417 for (octave_idx_type i = j+1; i < r_nr; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
418 tmp.xelem (j, i) = tmp.xelem (i, j); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
419 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
420 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
421 retval = tmp; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
422 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
423 return retval; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
424 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
425 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
426 template <typename T> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
427 T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
428 chol2inv (const T& r) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
429 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
430 return chol2inv_internal (r); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
431 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
432 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
433 // Compute the inverse of a matrix using the Cholesky factorization. |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
434 template <typename T> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
435 T |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
436 chol<T>::inverse (void) const |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
437 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
438 return chol2inv_internal (chol_mat, is_upper); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
439 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
440 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
441 template <typename T> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
442 void |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
443 chol<T>::set (const T& R) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
444 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
445 if (! R.is_square ()) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
446 (*current_liboctave_error_handler) ("chol: requires square matrix"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
447 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
448 chol_mat = R; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
449 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
450 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
451 #if ! defined (HAVE_QRUPDATE) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
452 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
453 template <typename T> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
454 void |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
455 chol<T>::update (const VT& u) |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
456 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
457 warn_qrupdate_once (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
458 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
459 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
460 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
461 if (u.numel () != n) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
462 (*current_liboctave_error_handler) ("cholupdate: dimension mismatch"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
463 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
464 init (chol_mat.hermitian () * chol_mat + T (u) * T (u).hermitian (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
465 true, false); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
466 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
467 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
468 template <typename T> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
469 static bool |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
470 singular (const T& a) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
471 { |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
472 static typename T::element_type zero (0); |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
473 for (octave_idx_type i = 0; i < a.rows (); i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
474 if (a(i,i) == zero) return true; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
475 return false; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
476 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
477 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
478 template <typename T> |
5275 | 479 octave_idx_type |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
480 chol<T>::downdate (const VT& u) |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
481 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
482 warn_qrupdate_once (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
483 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
484 octave_idx_type info = -1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
485 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
486 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
487 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
488 if (u.numel () != n) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
489 (*current_liboctave_error_handler) ("cholupdate: dimension mismatch"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
490 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
491 if (singular (chol_mat)) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
492 info = 2; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
493 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
494 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
495 info = init (chol_mat.hermitian () * chol_mat |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
496 - T (u) * T (u).hermitian (), true, false); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
497 if (info) info = 1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
498 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
499 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
500 return info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
501 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
502 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
503 template <typename T> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
504 octave_idx_type |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
505 chol<T>::insert_sym (const VT& u, octave_idx_type j) |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
506 { |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
507 static typename T::element_type zero (0); |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
508 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
509 warn_qrupdate_once (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
510 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
511 octave_idx_type info = -1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
512 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
513 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
514 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
515 if (u.numel () != n + 1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
516 (*current_liboctave_error_handler) ("cholinsert: dimension mismatch"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
517 if (j < 0 || j > n) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
518 (*current_liboctave_error_handler) ("cholinsert: index out of range"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
519 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
520 if (singular (chol_mat)) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
521 info = 2; |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
522 else if (imag (u(j)) != zero) |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
523 info = 3; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
524 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
525 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
526 T a = chol_mat.hermitian () * chol_mat; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
527 T a1 (n+1, n+1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
528 for (octave_idx_type k = 0; k < n+1; k++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
529 for (octave_idx_type l = 0; l < n+1; l++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
530 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
531 if (l == j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
532 a1(k, l) = u(k); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
533 else if (k == j) |
21270
230e186e292d
make building without qrupdate work again
John W. Eaton <jwe@octave.org>
parents:
21269
diff
changeset
|
534 a1(k, l) = conj (u(l)); |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
535 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
536 a1(k, l) = a(k < j ? k : k-1, l < j ? l : l-1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
537 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
538 info = init (a1, true, false); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
539 if (info) info = 1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
540 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
541 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
542 return info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
543 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
544 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
545 template <typename T> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
546 void |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
547 chol<T>::delete_sym (octave_idx_type j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
548 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
549 warn_qrupdate_once (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
550 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
551 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
552 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
553 if (j < 0 || j > n-1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
554 (*current_liboctave_error_handler) ("choldelete: index out of range"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
555 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
556 T a = chol_mat.hermitian () * chol_mat; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
557 a.delete_elements (1, idx_vector (j)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
558 a.delete_elements (0, idx_vector (j)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
559 init (a, true, false); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
560 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
561 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
562 template <typename T> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
563 void |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
564 chol<T>::shift_sym (octave_idx_type i, octave_idx_type j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
565 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
566 warn_qrupdate_once (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
567 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
568 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
569 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
570 if (i < 0 || i > n-1 || j < 0 || j > n-1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
571 (*current_liboctave_error_handler) ("cholshift: index out of range"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
572 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
573 T a = chol_mat.hermitian () * chol_mat; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
574 Array<octave_idx_type> p (dim_vector (n, 1)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
575 for (octave_idx_type k = 0; k < n; k++) p(k) = k; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
576 if (i < j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
577 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
578 for (octave_idx_type k = i; k < j; k++) p(k) = k+1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
579 p(j) = i; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
580 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
581 else if (j < i) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
582 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
583 p(j) = i; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
584 for (octave_idx_type k = j+1; k < i+1; k++) p(k) = k-1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
585 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
586 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
587 init (a.index (idx_vector (p), idx_vector (p)), true, false); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
588 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
589 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
590 #endif |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
591 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
592 // Specializations. |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
593 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
594 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
595 octave_idx_type |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
596 chol<Matrix>::init (const Matrix& a, bool upper, bool calc_cond) |
457 | 597 { |
5275 | 598 octave_idx_type a_nr = a.rows (); |
599 octave_idx_type a_nc = a.cols (); | |
1944 | 600 |
457 | 601 if (a_nr != a_nc) |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
602 (*current_liboctave_error_handler) ("chol: requires square matrix"); |
457 | 603 |
5275 | 604 octave_idx_type n = a_nc; |
605 octave_idx_type info; | |
457 | 606 |
20462
5ce959c55cc0
Propagate 'lower' in chol(a, 'lower') to underlying library function.
PrasannaKumar Muralidharan <prasannatsmkumar@gmail.com>
parents:
20232
diff
changeset
|
607 is_upper = upper; |
5ce959c55cc0
Propagate 'lower' in chol(a, 'lower') to underlying library function.
PrasannaKumar Muralidharan <prasannatsmkumar@gmail.com>
parents:
20232
diff
changeset
|
608 |
9862
c0aeedd8fb86
improve chol Matlab compatibility
Jaroslav Hajek <highegg@gmail.com>
parents:
8562
diff
changeset
|
609 chol_mat.clear (n, n); |
20462
5ce959c55cc0
Propagate 'lower' in chol(a, 'lower') to underlying library function.
PrasannaKumar Muralidharan <prasannatsmkumar@gmail.com>
parents:
20232
diff
changeset
|
610 if (is_upper) |
20628
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
611 for (octave_idx_type j = 0; j < n; j++) |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
612 { |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
613 for (octave_idx_type i = 0; i <= j; i++) |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
614 chol_mat.xelem (i, j) = a(i, j); |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
615 for (octave_idx_type i = j+1; i < n; i++) |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
616 chol_mat.xelem (i, j) = 0.0; |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
617 } |
20462
5ce959c55cc0
Propagate 'lower' in chol(a, 'lower') to underlying library function.
PrasannaKumar Muralidharan <prasannatsmkumar@gmail.com>
parents:
20232
diff
changeset
|
618 else |
20628
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
619 for (octave_idx_type j = 0; j < n; j++) |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
620 { |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
621 for (octave_idx_type i = 0; i < j; i++) |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
622 chol_mat.xelem (i, j) = 0.0; |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
623 for (octave_idx_type i = j; i < n; i++) |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
624 chol_mat.xelem (i, j) = a(i, j); |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
625 } |
1944 | 626 double *h = chol_mat.fortran_vec (); |
457 | 627 |
6486 | 628 // Calculate the norm of the matrix, for later use. |
629 double anorm = 0; | |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11570
diff
changeset
|
630 if (calc_cond) |
9862
c0aeedd8fb86
improve chol Matlab compatibility
Jaroslav Hajek <highegg@gmail.com>
parents:
8562
diff
changeset
|
631 anorm = xnorm (a, 1); |
6486 | 632 |
20462
5ce959c55cc0
Propagate 'lower' in chol(a, 'lower') to underlying library function.
PrasannaKumar Muralidharan <prasannatsmkumar@gmail.com>
parents:
20232
diff
changeset
|
633 if (is_upper) |
20628
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
634 F77_XFCN (dpotrf, DPOTRF, (F77_CONST_CHAR_ARG2 ("U", 1), n, h, n, info |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
635 F77_CHAR_ARG_LEN (1))); |
20462
5ce959c55cc0
Propagate 'lower' in chol(a, 'lower') to underlying library function.
PrasannaKumar Muralidharan <prasannatsmkumar@gmail.com>
parents:
20232
diff
changeset
|
636 else |
20628
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
637 F77_XFCN (dpotrf, DPOTRF, (F77_CONST_CHAR_ARG2 ("L", 1), n, h, n, info |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
638 F77_CHAR_ARG_LEN (1))); |
457 | 639 |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7017
diff
changeset
|
640 xrcond = 0.0; |
9862
c0aeedd8fb86
improve chol Matlab compatibility
Jaroslav Hajek <highegg@gmail.com>
parents:
8562
diff
changeset
|
641 if (info > 0) |
c0aeedd8fb86
improve chol Matlab compatibility
Jaroslav Hajek <highegg@gmail.com>
parents:
8562
diff
changeset
|
642 chol_mat.resize (info - 1, info - 1); |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11570
diff
changeset
|
643 else if (calc_cond) |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7017
diff
changeset
|
644 { |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7017
diff
changeset
|
645 octave_idx_type dpocon_info = 0; |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7017
diff
changeset
|
646 |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7017
diff
changeset
|
647 // Now calculate the condition number for non-singular matrix. |
11570
57632dea2446
attempt better backward compatibility for Array constructors
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
648 Array<double> z (dim_vector (3*n, 1)); |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7017
diff
changeset
|
649 double *pz = z.fortran_vec (); |
11570
57632dea2446
attempt better backward compatibility for Array constructors
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
650 Array<octave_idx_type> iz (dim_vector (n, 1)); |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7017
diff
changeset
|
651 octave_idx_type *piz = iz.fortran_vec (); |
20462
5ce959c55cc0
Propagate 'lower' in chol(a, 'lower') to underlying library function.
PrasannaKumar Muralidharan <prasannatsmkumar@gmail.com>
parents:
20232
diff
changeset
|
652 if (is_upper) |
20628
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
653 F77_XFCN (dpocon, DPOCON, (F77_CONST_CHAR_ARG2 ("U", 1), n, h, |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
654 n, anorm, xrcond, pz, piz, dpocon_info |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
655 F77_CHAR_ARG_LEN (1))); |
20462
5ce959c55cc0
Propagate 'lower' in chol(a, 'lower') to underlying library function.
PrasannaKumar Muralidharan <prasannatsmkumar@gmail.com>
parents:
20232
diff
changeset
|
656 else |
20628
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
657 F77_XFCN (dpocon, DPOCON, (F77_CONST_CHAR_ARG2 ("L", 1), n, h, |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
658 n, anorm, xrcond, pz, piz, dpocon_info |
48fedd8fbff7
maint: Apply Octave coding style to Cholesky classes
Mike Miller <mtmiller@octave.org>
parents:
20531
diff
changeset
|
659 F77_CHAR_ARG_LEN (1))); |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7017
diff
changeset
|
660 |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11570
diff
changeset
|
661 if (dpocon_info != 0) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
662 info = -1; |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7017
diff
changeset
|
663 } |
457 | 664 |
665 return info; | |
666 } | |
667 | |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
668 #if defined (HAVE_QRUPDATE) |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
669 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
670 template <> |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
671 void |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
672 chol<Matrix>::update (const ColumnVector& u) |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
673 { |
7559
07522d7dcdf8
fixes to QR and Cholesky updating code
Jaroslav Hajek <highegg@gmail.com>
parents:
7554
diff
changeset
|
674 octave_idx_type n = chol_mat.rows (); |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
675 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
676 if (u.numel () != n) |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
677 (*current_liboctave_error_handler) ("cholupdate: dimension mismatch"); |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
678 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
679 ColumnVector utmp = u; |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
680 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
681 OCTAVE_LOCAL_BUFFER (double, w, n); |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
682 |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
683 F77_XFCN (dch1up, DCH1UP, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
684 utmp.fortran_vec (), w)); |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
685 } |
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
686 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
687 template <> |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
688 octave_idx_type |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
689 chol<Matrix>::downdate (const ColumnVector& u) |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
690 { |
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
691 octave_idx_type info = -1; |
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
692 |
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
693 octave_idx_type n = chol_mat.rows (); |
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
694 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
695 if (u.numel () != n) |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
696 (*current_liboctave_error_handler) ("cholupdate: dimension mismatch"); |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
697 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
698 ColumnVector utmp = u; |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
699 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
700 OCTAVE_LOCAL_BUFFER (double, w, n); |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
701 |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
702 F77_XFCN (dch1dn, DCH1DN, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
703 utmp.fortran_vec (), w, info)); |
7554
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
704 |
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
705 return info; |
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
706 } |
40574114c514
implement Cholesky factorization updating
Jaroslav Hajek <highegg@gmail.com>
parents:
7482
diff
changeset
|
707 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
708 template <> |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
709 octave_idx_type |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
710 chol<Matrix>::insert_sym (const ColumnVector& u, octave_idx_type j) |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
711 { |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
712 octave_idx_type info = -1; |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
713 |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
714 octave_idx_type n = chol_mat.rows (); |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11570
diff
changeset
|
715 |
20232
a9574e3c6e9e
Deprecate Array::length() and Sparse::length() in favour of ::numel().
Carnë Draug <carandraug@octave.org>
parents:
19697
diff
changeset
|
716 if (u.numel () != n + 1) |
8547
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
717 (*current_liboctave_error_handler) ("cholinsert: dimension mismatch"); |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
718 if (j < 0 || j > n) |
8547
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
719 (*current_liboctave_error_handler) ("cholinsert: index out of range"); |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
720 |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
721 ColumnVector utmp = u; |
8547
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
722 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
723 OCTAVE_LOCAL_BUFFER (double, w, n); |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
724 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
725 chol_mat.resize (n+1, n+1); |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
726 |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
727 F77_XFCN (dchinx, DCHINX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
728 j + 1, utmp.fortran_vec (), w, info)); |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
729 |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
730 return info; |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
731 } |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
732 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
733 template <> |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
734 void |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
735 chol<Matrix>::delete_sym (octave_idx_type j) |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
736 { |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
737 octave_idx_type n = chol_mat.rows (); |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11570
diff
changeset
|
738 |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
739 if (j < 0 || j > n-1) |
8547
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
740 (*current_liboctave_error_handler) ("choldelete: index out of range"); |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
741 |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
742 OCTAVE_LOCAL_BUFFER (double, w, n); |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
743 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
744 F77_XFCN (dchdex, DCHDEX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
745 j + 1, w)); |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
746 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
747 chol_mat.resize (n-1, n-1); |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
748 } |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
749 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
750 template <> |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
751 void |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
752 chol<Matrix>::shift_sym (octave_idx_type i, octave_idx_type j) |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
753 { |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
754 octave_idx_type n = chol_mat.rows (); |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11570
diff
changeset
|
755 |
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11570
diff
changeset
|
756 if (i < 0 || i > n-1 || j < 0 || j > n-1) |
8547
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
757 (*current_liboctave_error_handler) ("cholshift: index out of range"); |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
758 |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
759 OCTAVE_LOCAL_BUFFER (double, w, 2*n); |
8547
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
760 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
761 F77_XFCN (dchshx, DCHSHX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
762 i + 1, j + 1, w)); |
7700
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
763 } |
efccca5f2ad7
more QR & Cholesky updating functions
Jaroslav Hajek <highegg@gmail.com>
parents:
7559
diff
changeset
|
764 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
765 #endif |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
766 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
767 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
768 octave_idx_type |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
769 chol<FloatMatrix>::init (const FloatMatrix& a, bool upper, bool calc_cond) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
770 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
771 octave_idx_type a_nr = a.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
772 octave_idx_type a_nc = a.cols (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
773 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
774 if (a_nr != a_nc) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
775 (*current_liboctave_error_handler) ("chol: requires square matrix"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
776 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
777 octave_idx_type n = a_nc; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
778 octave_idx_type info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
779 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
780 is_upper = upper; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
781 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
782 chol_mat.clear (n, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
783 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
784 for (octave_idx_type j = 0; j < n; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
785 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
786 for (octave_idx_type i = 0; i <= j; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
787 chol_mat.xelem (i, j) = a(i, j); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
788 for (octave_idx_type i = j+1; i < n; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
789 chol_mat.xelem (i, j) = 0.0f; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
790 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
791 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
792 for (octave_idx_type j = 0; j < n; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
793 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
794 for (octave_idx_type i = 0; i < j; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
795 chol_mat.xelem (i, j) = 0.0f; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
796 for (octave_idx_type i = j; i < n; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
797 chol_mat.xelem (i, j) = a(i, j); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
798 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
799 float *h = chol_mat.fortran_vec (); |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
800 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
801 // Calculate the norm of the matrix, for later use. |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
802 float anorm = 0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
803 if (calc_cond) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
804 anorm = xnorm (a, 1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
805 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
806 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
807 F77_XFCN (spotrf, SPOTRF, (F77_CONST_CHAR_ARG2 ("U", 1), n, h, n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
808 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
809 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
810 F77_XFCN (spotrf, SPOTRF, (F77_CONST_CHAR_ARG2 ("L", 1), n, h, n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
811 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
812 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
813 xrcond = 0.0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
814 if (info > 0) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
815 chol_mat.resize (info - 1, info - 1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
816 else if (calc_cond) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
817 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
818 octave_idx_type spocon_info = 0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
819 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
820 // Now calculate the condition number for non-singular matrix. |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
821 Array<float> z (dim_vector (3*n, 1)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
822 float *pz = z.fortran_vec (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
823 Array<octave_idx_type> iz (dim_vector (n, 1)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
824 octave_idx_type *piz = iz.fortran_vec (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
825 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
826 F77_XFCN (spocon, SPOCON, (F77_CONST_CHAR_ARG2 ("U", 1), n, h, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
827 n, anorm, xrcond, pz, piz, spocon_info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
828 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
829 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
830 F77_XFCN (spocon, SPOCON, (F77_CONST_CHAR_ARG2 ("L", 1), n, h, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
831 n, anorm, xrcond, pz, piz, spocon_info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
832 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
833 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
834 if (spocon_info != 0) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
835 info = -1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
836 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
837 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
838 return info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
839 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
840 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
841 #ifdef HAVE_QRUPDATE |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
842 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
843 template <> |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
844 void |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
845 chol<FloatMatrix>::update (const FloatColumnVector& u) |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
846 { |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
847 octave_idx_type n = chol_mat.rows (); |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
848 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
849 if (u.numel () != n) |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
850 (*current_liboctave_error_handler) ("cholupdate: dimension mismatch"); |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
851 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
852 FloatColumnVector utmp = u; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
853 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
854 OCTAVE_LOCAL_BUFFER (float, w, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
855 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
856 F77_XFCN (sch1up, SCH1UP, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
857 utmp.fortran_vec (), w)); |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
858 } |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
859 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
860 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
861 octave_idx_type |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
862 chol<FloatMatrix>::downdate (const FloatColumnVector& u) |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
863 { |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
864 octave_idx_type info = -1; |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
865 |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
866 octave_idx_type n = chol_mat.rows (); |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
867 |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
868 if (u.numel () != n) |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
869 (*current_liboctave_error_handler) ("cholupdate: dimension mismatch"); |
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
870 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
871 FloatColumnVector utmp = u; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
872 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
873 OCTAVE_LOCAL_BUFFER (float, w, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
874 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
875 F77_XFCN (sch1dn, SCH1DN, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
876 utmp.fortran_vec (), w, info)); |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
877 |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
878 return info; |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
879 } |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
880 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
881 template <> |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
882 octave_idx_type |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
883 chol<FloatMatrix>::insert_sym (const FloatColumnVector& u, octave_idx_type j) |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
884 { |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
885 octave_idx_type info = -1; |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
886 |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
887 octave_idx_type n = chol_mat.rows (); |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
888 |
20232
a9574e3c6e9e
Deprecate Array::length() and Sparse::length() in favour of ::numel().
Carnë Draug <carandraug@octave.org>
parents:
19697
diff
changeset
|
889 if (u.numel () != n + 1) |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
890 (*current_liboctave_error_handler) ("cholinsert: dimension mismatch"); |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
891 if (j < 0 || j > n) |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
892 (*current_liboctave_error_handler) ("cholinsert: index out of range"); |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
893 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
894 FloatColumnVector utmp = u; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
895 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
896 OCTAVE_LOCAL_BUFFER (float, w, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
897 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
898 chol_mat.resize (n+1, n+1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
899 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
900 F77_XFCN (schinx, SCHINX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
901 j + 1, utmp.fortran_vec (), w, info)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
902 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
903 return info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
904 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
905 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
906 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
907 void |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
908 chol<FloatMatrix>::delete_sym (octave_idx_type j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
909 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
910 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
911 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
912 if (j < 0 || j > n-1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
913 (*current_liboctave_error_handler) ("choldelete: index out of range"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
914 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
915 OCTAVE_LOCAL_BUFFER (float, w, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
916 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
917 F77_XFCN (schdex, SCHDEX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
918 j + 1, w)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
919 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
920 chol_mat.resize (n-1, n-1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
921 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
922 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
923 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
924 void |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
925 chol<FloatMatrix>::shift_sym (octave_idx_type i, octave_idx_type j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
926 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
927 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
928 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
929 if (i < 0 || i > n-1 || j < 0 || j > n-1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
930 (*current_liboctave_error_handler) ("cholshift: index out of range"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
931 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
932 OCTAVE_LOCAL_BUFFER (float, w, 2*n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
933 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
934 F77_XFCN (schshx, SCHSHX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
935 i + 1, j + 1, w)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
936 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
937 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
938 #endif |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
939 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
940 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
941 octave_idx_type |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
942 chol<ComplexMatrix>::init (const ComplexMatrix& a, bool upper, bool calc_cond) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
943 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
944 octave_idx_type a_nr = a.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
945 octave_idx_type a_nc = a.cols (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
946 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
947 if (a_nr != a_nc) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
948 (*current_liboctave_error_handler) ("chol: requires square matrix"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
949 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
950 octave_idx_type n = a_nc; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
951 octave_idx_type info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
952 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
953 is_upper = upper; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
954 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
955 chol_mat.clear (n, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
956 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
957 for (octave_idx_type j = 0; j < n; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
958 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
959 for (octave_idx_type i = 0; i <= j; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
960 chol_mat.xelem (i, j) = a(i, j); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
961 for (octave_idx_type i = j+1; i < n; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
962 chol_mat.xelem (i, j) = 0.0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
963 } |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
964 else |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
965 for (octave_idx_type j = 0; j < n; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
966 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
967 for (octave_idx_type i = 0; i < j; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
968 chol_mat.xelem (i, j) = 0.0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
969 for (octave_idx_type i = j; i < n; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
970 chol_mat.xelem (i, j) = a(i, j); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
971 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
972 Complex *h = chol_mat.fortran_vec (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
973 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
974 // Calculate the norm of the matrix, for later use. |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
975 double anorm = 0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
976 if (calc_cond) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
977 anorm = xnorm (a, 1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
978 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
979 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
980 F77_XFCN (zpotrf, ZPOTRF, (F77_CONST_CHAR_ARG2 ("U", 1), n, h, n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
981 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
982 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
983 F77_XFCN (zpotrf, ZPOTRF, (F77_CONST_CHAR_ARG2 ("L", 1), n, h, n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
984 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
985 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
986 xrcond = 0.0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
987 if (info > 0) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
988 chol_mat.resize (info - 1, info - 1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
989 else if (calc_cond) |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
990 { |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
991 octave_idx_type zpocon_info = 0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
992 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
993 // Now calculate the condition number for non-singular matrix. |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
994 Array<Complex> z (dim_vector (2*n, 1)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
995 Complex *pz = z.fortran_vec (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
996 Array<double> rz (dim_vector (n, 1)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
997 double *prz = rz.fortran_vec (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
998 F77_XFCN (zpocon, ZPOCON, (F77_CONST_CHAR_ARG2 ("U", 1), n, h, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
999 n, anorm, xrcond, pz, prz, zpocon_info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1000 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1001 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1002 if (zpocon_info != 0) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1003 info = -1; |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1004 } |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1005 |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1006 return info; |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1007 } |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1008 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1009 #ifdef HAVE_QRUPDATE |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1010 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1011 template <> |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1012 void |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1013 chol<ComplexMatrix>::update (const ComplexColumnVector& u) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1014 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1015 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1016 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1017 if (u.numel () != n) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1018 (*current_liboctave_error_handler) ("cholupdate: dimension mismatch"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1019 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1020 ComplexColumnVector utmp = u; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1021 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1022 OCTAVE_LOCAL_BUFFER (double, rw, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1023 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1024 F77_XFCN (zch1up, ZCH1UP, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1025 utmp.fortran_vec (), rw)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1026 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1027 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1028 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1029 octave_idx_type |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1030 chol<ComplexMatrix>::downdate (const ComplexColumnVector& u) |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1031 { |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1032 octave_idx_type info = -1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1033 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1034 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1035 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1036 if (u.numel () != n) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1037 (*current_liboctave_error_handler) ("cholupdate: dimension mismatch"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1038 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1039 ComplexColumnVector utmp = u; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1040 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1041 OCTAVE_LOCAL_BUFFER (double, rw, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1042 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1043 F77_XFCN (zch1dn, ZCH1DN, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1044 utmp.fortran_vec (), rw, info)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1045 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1046 return info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1047 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1048 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1049 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1050 octave_idx_type |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1051 chol<ComplexMatrix>::insert_sym (const ComplexColumnVector& u, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1052 octave_idx_type j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1053 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1054 octave_idx_type info = -1; |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1055 |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1056 octave_idx_type n = chol_mat.rows (); |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1057 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1058 if (u.numel () != n + 1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1059 (*current_liboctave_error_handler) ("cholinsert: dimension mismatch"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1060 if (j < 0 || j > n) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1061 (*current_liboctave_error_handler) ("cholinsert: index out of range"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1062 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1063 ComplexColumnVector utmp = u; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1064 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1065 OCTAVE_LOCAL_BUFFER (double, rw, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1066 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1067 chol_mat.resize (n+1, n+1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1068 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1069 F77_XFCN (zchinx, ZCHINX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1070 j + 1, utmp.fortran_vec (), rw, info)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1071 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1072 return info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1073 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1074 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1075 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1076 void |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1077 chol<ComplexMatrix>::delete_sym (octave_idx_type j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1078 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1079 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1080 |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1081 if (j < 0 || j > n-1) |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1082 (*current_liboctave_error_handler) ("choldelete: index out of range"); |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
1083 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1084 OCTAVE_LOCAL_BUFFER (double, rw, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1085 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1086 F77_XFCN (zchdex, ZCHDEX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1087 j + 1, rw)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1088 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1089 chol_mat.resize (n-1, n-1); |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1090 } |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1091 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1092 template <> |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1093 void |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1094 chol<ComplexMatrix>::shift_sym (octave_idx_type i, octave_idx_type j) |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1095 { |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1096 octave_idx_type n = chol_mat.rows (); |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1097 |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11570
diff
changeset
|
1098 if (i < 0 || i > n-1 || j < 0 || j > n-1) |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1099 (*current_liboctave_error_handler) ("cholshift: index out of range"); |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20629
diff
changeset
|
1100 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1101 OCTAVE_LOCAL_BUFFER (Complex, w, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1102 OCTAVE_LOCAL_BUFFER (double, rw, n); |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1103 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1104 F77_XFCN (zchshx, ZCHSHX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1105 i + 1, j + 1, w, rw)); |
8562
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1106 } |
a6edd5c23cb5
use replacement methods if qrupdate is not available
Jaroslav Hajek <highegg@gmail.com>
parents:
8547
diff
changeset
|
1107 |
8547
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
1108 #endif |
d66c9b6e506a
imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents:
8377
diff
changeset
|
1109 |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1110 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1111 octave_idx_type |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1112 chol<FloatComplexMatrix>::init (const FloatComplexMatrix& a, bool upper, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1113 bool calc_cond) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1114 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1115 octave_idx_type a_nr = a.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1116 octave_idx_type a_nc = a.cols (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1117 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1118 if (a_nr != a_nc) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1119 (*current_liboctave_error_handler) ("chol: requires square matrix"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1120 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1121 octave_idx_type n = a_nc; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1122 octave_idx_type info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1123 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1124 is_upper = upper; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1125 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1126 chol_mat.clear (n, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1127 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1128 for (octave_idx_type j = 0; j < n; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1129 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1130 for (octave_idx_type i = 0; i <= j; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1131 chol_mat.xelem (i, j) = a(i, j); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1132 for (octave_idx_type i = j+1; i < n; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1133 chol_mat.xelem (i, j) = 0.0f; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1134 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1135 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1136 for (octave_idx_type j = 0; j < n; j++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1137 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1138 for (octave_idx_type i = 0; i < j; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1139 chol_mat.xelem (i, j) = 0.0f; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1140 for (octave_idx_type i = j; i < n; i++) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1141 chol_mat.xelem (i, j) = a(i, j); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1142 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1143 FloatComplex *h = chol_mat.fortran_vec (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1144 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1145 // Calculate the norm of the matrix, for later use. |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1146 float anorm = 0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1147 if (calc_cond) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1148 anorm = xnorm (a, 1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1149 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1150 if (is_upper) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1151 F77_XFCN (cpotrf, CPOTRF, (F77_CONST_CHAR_ARG2 ("U", 1), n, h, n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1152 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1153 else |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1154 F77_XFCN (cpotrf, CPOTRF, (F77_CONST_CHAR_ARG2 ("L", 1), n, h, n, info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1155 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1156 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1157 xrcond = 0.0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1158 if (info > 0) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1159 chol_mat.resize (info - 1, info - 1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1160 else if (calc_cond) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1161 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1162 octave_idx_type cpocon_info = 0; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1163 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1164 // Now calculate the condition number for non-singular matrix. |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1165 Array<FloatComplex> z (dim_vector (2*n, 1)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1166 FloatComplex *pz = z.fortran_vec (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1167 Array<float> rz (dim_vector (n, 1)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1168 float *prz = rz.fortran_vec (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1169 F77_XFCN (cpocon, CPOCON, (F77_CONST_CHAR_ARG2 ("U", 1), n, h, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1170 n, anorm, xrcond, pz, prz, cpocon_info |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1171 F77_CHAR_ARG_LEN (1))); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1172 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1173 if (cpocon_info != 0) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1174 info = -1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1175 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1176 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1177 return info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1178 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1179 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1180 #ifdef HAVE_QRUPDATE |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1181 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1182 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1183 void |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1184 chol<FloatComplexMatrix>::update (const FloatComplexColumnVector& u) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1185 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1186 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1187 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1188 if (u.numel () != n) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1189 (*current_liboctave_error_handler) ("cholupdate: dimension mismatch"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1190 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1191 FloatComplexColumnVector utmp = u; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1192 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1193 OCTAVE_LOCAL_BUFFER (float, rw, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1194 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1195 F77_XFCN (cch1up, CCH1UP, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1196 utmp.fortran_vec (), rw)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1197 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1198 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1199 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1200 octave_idx_type |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1201 chol<FloatComplexMatrix>::downdate (const FloatComplexColumnVector& u) |
5340 | 1202 { |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1203 octave_idx_type info = -1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1204 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1205 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1206 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1207 if (u.numel () != n) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1208 (*current_liboctave_error_handler) ("cholupdate: dimension mismatch"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1209 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1210 FloatComplexColumnVector utmp = u; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1211 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1212 OCTAVE_LOCAL_BUFFER (float, rw, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1213 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1214 F77_XFCN (cch1dn, CCH1DN, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1215 utmp.fortran_vec (), rw, info)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1216 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1217 return info; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1218 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1219 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1220 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1221 octave_idx_type |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1222 chol<FloatComplexMatrix>::insert_sym (const FloatComplexColumnVector& u, |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1223 octave_idx_type j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1224 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1225 octave_idx_type info = -1; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1226 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1227 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1228 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1229 if (u.numel () != n + 1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1230 (*current_liboctave_error_handler) ("cholinsert: dimension mismatch"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1231 if (j < 0 || j > n) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1232 (*current_liboctave_error_handler) ("cholinsert: index out of range"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1233 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1234 FloatComplexColumnVector utmp = u; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1235 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1236 OCTAVE_LOCAL_BUFFER (float, rw, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1237 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1238 chol_mat.resize (n+1, n+1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1239 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1240 F77_XFCN (cchinx, CCHINX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1241 j + 1, utmp.fortran_vec (), rw, info)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1242 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1243 return info; |
5340 | 1244 } |
21269
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1245 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1246 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1247 void |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1248 chol<FloatComplexMatrix>::delete_sym (octave_idx_type j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1249 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1250 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1251 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1252 if (j < 0 || j > n-1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1253 (*current_liboctave_error_handler) ("choldelete: index out of range"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1254 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1255 OCTAVE_LOCAL_BUFFER (float, rw, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1256 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1257 F77_XFCN (cchdex, CCHDEX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1258 j + 1, rw)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1259 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1260 chol_mat.resize (n-1, n-1); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1261 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1262 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1263 template <> |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1264 void |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1265 chol<FloatComplexMatrix>::shift_sym (octave_idx_type i, octave_idx_type j) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1266 { |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1267 octave_idx_type n = chol_mat.rows (); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1268 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1269 if (i < 0 || i > n-1 || j < 0 || j > n-1) |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1270 (*current_liboctave_error_handler) ("cholshift: index out of range"); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1271 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1272 OCTAVE_LOCAL_BUFFER (FloatComplex, w, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1273 OCTAVE_LOCAL_BUFFER (float, rw, n); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1274 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1275 F77_XFCN (cchshx, CCHSHX, (n, chol_mat.fortran_vec (), chol_mat.rows (), |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1276 i + 1, j + 1, w, rw)); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1277 } |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1278 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1279 #endif |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1280 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1281 // Instantiations we need. |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1282 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1283 template class chol<Matrix>; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1284 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1285 template class chol<FloatMatrix>; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1286 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1287 template class chol<ComplexMatrix>; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1288 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1289 template class chol<FloatComplexMatrix>; |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1290 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1291 template Matrix |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1292 chol2inv<Matrix> (const Matrix& r); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1293 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1294 template ComplexMatrix |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1295 chol2inv<ComplexMatrix> (const ComplexMatrix& r); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1296 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1297 template FloatMatrix |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1298 chol2inv<FloatMatrix> (const FloatMatrix& r); |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1299 |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1300 template FloatComplexMatrix |
3c8a3d35661a
better use of templates for Cholesky factorization
John W. Eaton <jwe@octave.org>
parents:
21202
diff
changeset
|
1301 chol2inv<FloatComplexMatrix> (const FloatComplexMatrix& r); |