Mercurial > octave-dspies

annotate liboctave/floatQR.cc @ 8547:d66c9b6e506a

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
imported patch qrupdate.diff
author Jaroslav Hajek <highegg@gmail.com>
date Tue, 20 Jan 2009 21:16:42 +0100
parents 4976f66d469b
children a6edd5c23cb5
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
1 /*
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
2
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
3 Copyright (C) 1994, 1995, 1996, 1997, 2002, 2003, 2004, 2005, 2007
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
4 John W. Eaton
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
5 Copyright (C) 2008, 2009 Jaroslav Hajek
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
6
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
7 This file is part of Octave.
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
8
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
9 Octave is free software; you can redistribute it and/or modify it
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
10 under the terms of the GNU General Public License as published by the
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
11 Free Software Foundation; either version 3 of the License, or (at your
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
12 option) any later version.
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
13
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
14 Octave is distributed in the hope that it will be useful, but WITHOUT
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
15 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
16 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
17 for more details.
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
18
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
19 You should have received a copy of the GNU General Public License
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
20 along with Octave; see the file COPYING. If not, see
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
21 <http://www.gnu.org/licenses/>.
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
22
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
23 */
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
24
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
25 #ifdef HAVE_CONFIG_H
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
26 #include <config.h>
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
27 #endif
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
28
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
29 #include "floatQR.h"
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
30 #include "f77-fcn.h"
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
31 #include "lo-error.h"
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
32 #include "Range.h"
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
33 #include "idx-vector.h"
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
34 #include "oct-locbuf.h"
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
35
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
36 extern "C"
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
37 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
38 F77_RET_T
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
39 F77_FUNC (sgeqrf, SGEQRF) (const octave_idx_type&, const octave_idx_type&, float*, const octave_idx_type&,
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
40 float*, float*, const octave_idx_type&, octave_idx_type&);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
41
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
42 F77_RET_T
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
43 F77_FUNC (sorgqr, SORGQR) (const octave_idx_type&, const octave_idx_type&, const octave_idx_type&, float*,
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
44 const octave_idx_type&, float*, float*, const octave_idx_type&, octave_idx_type&);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
45
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
46 #ifdef HAVE_QRUPDATE
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
47
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
48 F77_RET_T
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
49 F77_FUNC (sqr1up, SQR1UP) (const octave_idx_type&, const octave_idx_type&, const octave_idx_type&,
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
50 float*, const octave_idx_type&, float*, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
51 float*, float*, float*);
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
52
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
53 F77_RET_T
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
54 F77_FUNC (sqrinc, SQRINC) (const octave_idx_type&, const octave_idx_type&, const octave_idx_type&,
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
55 float*, const octave_idx_type&, float*, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
56 const octave_idx_type&, const float*, float*);
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
57
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
58 F77_RET_T
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
59 F77_FUNC (sqrdec, SQRDEC) (const octave_idx_type&, const octave_idx_type&, const octave_idx_type&,
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
60 float*, const octave_idx_type&, float*, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
61 const octave_idx_type&, float*);
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
62
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
63 F77_RET_T
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
64 F77_FUNC (sqrinr, SQRINR) (const octave_idx_type&, const octave_idx_type&,
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
65 float*, const octave_idx_type&, float*, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
66 const octave_idx_type&, const float*, float*);
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
67
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
68 F77_RET_T
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
69 F77_FUNC (sqrder, SQRDER) (const octave_idx_type&, const octave_idx_type&,
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
70 float*, const octave_idx_type&, float*, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
71 const octave_idx_type&, float*);
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
72
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
73 F77_RET_T
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
74 F77_FUNC (sqrshc, SQRSHC) (const octave_idx_type&, const octave_idx_type&, const octave_idx_type&,
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
75 float*, const octave_idx_type&, float*, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
76 const octave_idx_type&, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
77 float*);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
78
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
79 #endif
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
80 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
81
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
82 FloatQR::FloatQR (const FloatMatrix& a, QR::type qr_type)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
83 : q (), r ()
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
84 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
85 init (a, qr_type);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
86 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
87
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
88 void
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
89 FloatQR::init (const FloatMatrix& a, QR::type qr_type)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
90 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
91 octave_idx_type m = a.rows ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
92 octave_idx_type n = a.cols ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
93
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
94 if (m == 0 || n == 0)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
95 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
96 (*current_liboctave_error_handler) ("QR must have non-empty matrix");
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
97 return;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
98 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
99
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
100 octave_idx_type min_mn = m < n ? m : n;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
101 Array<float> tau (min_mn);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
102 float *ptau = tau.fortran_vec ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
103
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
104 octave_idx_type lwork = 32*n;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
105 Array<float> work (lwork);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
106 float *pwork = work.fortran_vec ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
107
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
108 octave_idx_type info = 0;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
109
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
110 FloatMatrix A_fact = a;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
111 if (m > n && qr_type != QR::economy)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
112 A_fact.resize (m, m, 0.0);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
113
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
114 float *tmp_data = A_fact.fortran_vec ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
115
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
116 F77_XFCN (sgeqrf, SGEQRF, (m, n, tmp_data, m, ptau, pwork, lwork, info));
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
117
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
118 if (qr_type == QR::raw)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
119 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
120 for (octave_idx_type j = 0; j < min_mn; j++)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
121 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
122 octave_idx_type limit = j < min_mn - 1 ? j : min_mn - 1;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
123 for (octave_idx_type i = limit + 1; i < m; i++)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
124 A_fact.elem (i, j) *= tau.elem (j);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
125 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
126
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
127 r = A_fact;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
128
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
129 if (m > n)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
130 r.resize (m, n);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
131 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
132 else
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
133 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
134 octave_idx_type n2 = (qr_type == QR::economy) ? min_mn : m;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
135
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
136 if (qr_type == QR::economy && m > n)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
137 r.resize (n, n, 0.0);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
138 else
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
139 r.resize (m, n, 0.0);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
140
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
141 for (octave_idx_type j = 0; j < n; j++)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
142 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
143 octave_idx_type limit = j < min_mn-1 ? j : min_mn-1;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
144 for (octave_idx_type i = 0; i <= limit; i++)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
145 r.elem (i, j) = tmp_data[m*j+i];
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
146 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
147
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
148 lwork = 32 * n2;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
149 work.resize (lwork);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
150 float *pwork2 = work.fortran_vec ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
151
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
152 F77_XFCN (sorgqr, SORGQR, (m, n2, min_mn, tmp_data, m, ptau,
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
153 pwork2, lwork, info));
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
154
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
155 q = A_fact;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
156 q.resize (m, n2);
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
157 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
158 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
159
7924
4976f66d469b miscellaneous cleanup
John W. Eaton <jwe@octave.org>
parents: 7792
diff changeset
160 FloatQR::FloatQR (const FloatMatrix& q_arg, const FloatMatrix& r_arg)
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
161 {
7924
4976f66d469b miscellaneous cleanup
John W. Eaton <jwe@octave.org>
parents: 7792
diff changeset
162 if (q_arg.columns () != r_arg.rows ())
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
163 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
164 (*current_liboctave_error_handler) ("QR dimensions mismatch");
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
165 return;
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
166 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
167
7924
4976f66d469b miscellaneous cleanup
John W. Eaton <jwe@octave.org>
parents: 7792
diff changeset
168 this->q = q_arg;
4976f66d469b miscellaneous cleanup
John W. Eaton <jwe@octave.org>
parents: 7792
diff changeset
169 this->r = r_arg;
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
170 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
171
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
172 #ifdef HAVE_QRUPDATE
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
173
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
174 void
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
175 FloatQR::update (const FloatColumnVector& u, const FloatColumnVector& v)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
176 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
177 octave_idx_type m = q.rows ();
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
178 octave_idx_type n = r.columns ();
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
179 octave_idx_type k = q.columns ();
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
180
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
181 if (u.length () == m && v.length () == n)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
182 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
183 FloatColumnVector utmp = u, vtmp = v;
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
184 OCTAVE_LOCAL_BUFFER (float, w, 2*k);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
185 F77_XFCN (sqr1up, SQR1UP, (m, n, k, q.fortran_vec (), m, r.fortran_vec (), k,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
186 utmp.fortran_vec (), vtmp.fortran_vec (), w));
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
187 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
188 else
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
189 (*current_liboctave_error_handler) ("QR update dimensions mismatch");
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
190 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
191
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
192 void
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
193 FloatQR::update (const FloatMatrix& u, const FloatMatrix& v)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
194 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
195 octave_idx_type m = q.rows ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
196 octave_idx_type n = r.columns ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
197 octave_idx_type k = q.columns ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
198
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
199 if (u.rows () == m && v.rows () == n && u.cols () == v.cols ())
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
200 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
201 OCTAVE_LOCAL_BUFFER (float, w, 2*k);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
202 for (octave_idx_type i = 0; i < u.cols (); i++)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
203 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
204 FloatColumnVector utmp = u.column (i), vtmp = v.column (i);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
205 F77_XFCN (sqr1up, SQR1UP, (m, n, k, q.fortran_vec (), m, r.fortran_vec (), k,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
206 utmp.fortran_vec (), vtmp.fortran_vec (), w));
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
207 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
208 }
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
209 else
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
210 (*current_liboctave_error_handler) ("qrupdate: dimensions mismatch");
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
211 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
212
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
213 void
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
214 FloatQR::insert_col (const FloatColumnVector& u, octave_idx_type j)
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
215 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
216 octave_idx_type m = q.rows ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
217 octave_idx_type n = r.columns ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
218 octave_idx_type k = q.columns ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
219
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
220 if (u.length () != m)
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
221 (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
222 else if (j < 0 || j > n)
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
223 (*current_liboctave_error_handler) ("qrinsert: index out of range");
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
224 else
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
225 {
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
226 if (k < m)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
227 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
228 q.resize (m, k+1);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
229 r.resize (k+1, n+1);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
230 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
231 else
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
232 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
233 r.resize (k, n+1);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
234 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
235
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
236 FloatColumnVector utmp = u;
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
237 OCTAVE_LOCAL_BUFFER (float, w, k);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
238 F77_XFCN (sqrinc, SQRINC, (m, n, k, q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
239 r.fortran_vec (), r.rows (), j + 1,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
240 utmp.data (), w));
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
241 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
242 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
243
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
244 void
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
245 FloatQR::insert_col (const FloatMatrix& u, const Array<octave_idx_type>& j)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
246 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
247 octave_idx_type m = q.rows ();
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
248 octave_idx_type n = r.columns ();
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
249 octave_idx_type k = q.columns ();
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
250
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
251 Array<octave_idx_type> jsi;
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
252 Array<octave_idx_type> js = j.sort (jsi, ASCENDING);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
253 octave_idx_type nj = js.length ();
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
254 bool dups = false;
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
255 for (octave_idx_type i = 0; i < nj - 1; i++)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
256 dups = dups && js(i) == js(i+1);
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
257
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
258 if (dups)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
259 (*current_liboctave_error_handler) ("qrinsert: duplicate index detected");
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
260 else if (u.length () != m || u.columns () != nj)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
261 (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
262 else if (nj > 0 && (js(0) < 0 || js(nj-1) > n))
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
263 (*current_liboctave_error_handler) ("qrinsert: index out of range");
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
264 else if (nj > 0)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
265 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
266 octave_idx_type kmax = std::min (k + nj, m);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
267 if (k < m)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
268 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
269 q.resize (m, kmax);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
270 r.resize (kmax, n + nj);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
271 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
272 else
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
273 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
274 r.resize (k, n + nj);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
275 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
276
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
277 OCTAVE_LOCAL_BUFFER (float, w, kmax);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
278 for (octave_idx_type i = 0; i < js.length (); i++)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
279 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
280 FloatColumnVector utmp = u.column (jsi(i));
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
281 F77_XFCN (sqrinc, SQRINC, (m, n + i, std::min (kmax, k + i),
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
282 q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
283 r.fortran_vec (), r.rows (), js(i) + 1,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
284 utmp.data (), w));
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
285 }
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
286 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
287 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
288
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
289 void
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
290 FloatQR::delete_col (octave_idx_type j)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
291 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
292 octave_idx_type m = q.rows ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
293 octave_idx_type k = r.rows ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
294 octave_idx_type n = r.columns ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
295
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
296 if (j < 0 || j > n-1)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
297 (*current_liboctave_error_handler) ("qrdelete: index out of range");
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
298 else
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
299 {
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
300 OCTAVE_LOCAL_BUFFER (float, w, k);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
301 F77_XFCN (sqrdec, SQRDEC, (m, n, k, q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
302 r.fortran_vec (), r.rows (), j + 1, w));
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
303
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
304 if (k < m)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
305 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
306 q.resize (m, k-1);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
307 r.resize (k-1, n-1);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
308 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
309 else
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
310 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
311 r.resize (k, n-1);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
312 }
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
313 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
314 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
315
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
316 void
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
317 FloatQR::delete_col (const Array<octave_idx_type>& j)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
318 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
319 octave_idx_type m = q.rows ();
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
320 octave_idx_type n = r.columns ();
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
321 octave_idx_type k = q.columns ();
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
322
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
323 Array<octave_idx_type> jsi;
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
324 Array<octave_idx_type> js = j.sort (jsi, DESCENDING);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
325 octave_idx_type nj = js.length ();
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
326 bool dups = false;
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
327 for (octave_idx_type i = 0; i < nj - 1; i++)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
328 dups = dups && js(i) == js(i+1);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
329
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
330 if (dups)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
331 (*current_liboctave_error_handler) ("qrinsert: duplicate index detected");
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
332 else if (nj > 0 && (js(0) > n-1 || js(nj-1) < 0))
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
333 (*current_liboctave_error_handler) ("qrinsert: index out of range");
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
334 else if (nj > 0)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
335 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
336 OCTAVE_LOCAL_BUFFER (float, w, k);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
337 for (octave_idx_type i = 0; i < js.length (); i++)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
338 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
339 F77_XFCN (sqrdec, SQRDEC, (m, n - i, k == m ? k : k - i,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
340 q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
341 r.fortran_vec (), r.rows (), js(i) + 1, w));
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
342 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
343 if (k < m)
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
344 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
345 q.resize (m, k - nj);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
346 r.resize (k - nj, n - nj);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
347 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
348 else
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
349 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
350 r.resize (k, n - nj);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
351 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
352
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
353 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
354 }
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
355
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
356 void
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
357 FloatQR::insert_row (const FloatRowVector& u, octave_idx_type j)
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
358 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
359 octave_idx_type m = r.rows ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
360 octave_idx_type n = r.columns ();
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
361 octave_idx_type k = std::min (m, n);
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
362
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
363 if (! q.is_square () || u.length () != n)
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
364 (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
365 else if (j < 0 || j > m)
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
366 (*current_liboctave_error_handler) ("qrinsert: index out of range");
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
367 else
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
368 {
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
369 q.resize (m + 1, m + 1);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
370 r.resize (m + 1, n);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
371 FloatRowVector utmp = u;
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
372 OCTAVE_LOCAL_BUFFER (float, w, k);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
373 F77_XFCN (sqrinr, SQRINR, (m, n, q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
374 r.fortran_vec (), r.rows (),
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
375 j + 1, utmp.fortran_vec (), w));
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
376
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
377 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
378 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
379
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
380 void
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
381 FloatQR::delete_row (octave_idx_type j)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
382 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
383 octave_idx_type m = r.rows ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
384 octave_idx_type n = r.columns ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
385
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
386 if (! q.is_square ())
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
387 (*current_liboctave_error_handler) ("qrdelete: dimensions mismatch");
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
388 else if (j < 0 || j > m-1)
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
389 (*current_liboctave_error_handler) ("qrdelete: index out of range");
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
390 else
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
391 {
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
392 OCTAVE_LOCAL_BUFFER (float, w, 2*m);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
393 F77_XFCN (sqrder, SQRDER, (m, n, q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
394 r.fortran_vec (), r.rows (), j + 1,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
395 w));
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
396
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
397 q.resize (m - 1, m - 1);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
398 r.resize (m - 1, n);
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
399 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
400 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
401
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
402 void
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
403 FloatQR::shift_cols (octave_idx_type i, octave_idx_type j)
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
404 {
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
405 octave_idx_type m = q.rows ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
406 octave_idx_type k = r.rows ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
407 octave_idx_type n = r.columns ();
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
408
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
409 if (i < 0 || i > n-1 || j < 0 || j > n-1)
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
410 (*current_liboctave_error_handler) ("qrshift: index out of range");
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
411 else
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
412 {
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
413 OCTAVE_LOCAL_BUFFER (float, w, 2*k);
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
414 F77_XFCN (sqrshc, SQRSHC, (m, n, k,
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
415 q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
416 r.fortran_vec (), r.rows (),
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
417 i + 1, j + 1, w));
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
418 }
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
419 }
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
420
8547
d66c9b6e506a imported patch qrupdate.diff
Jaroslav Hajek <highegg@gmail.com>
parents: 7924
diff changeset
421 #endif
7792
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
422
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
423
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
424 /*
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
425 ;;; Local Variables: ***
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
426 ;;; mode: C++ ***
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
427 ;;; End: ***
39c1026191e9 add missing files from single-precision merge
John W. Eaton <jwe@octave.org>
parents:
diff changeset
428 */