octave-lyh: liboctave/dbleQR.cc annotate

annotate liboctave/dbleQR.cc @ 8966:1bba53c0a38d

Implement diag + sparse, diag - sparse, sparse + diag, sparse - diag. Date: Mon, 9 Mar 2009 17:45:22 -0400 This does not use the typical sparse-mx-ops generator. I suspect the sematics of elementwise multiplication and division to be rather controversial, so they are not included. If comparison operations are added, the implementation should be shifted over to use the typical generator. The template in Sparse-diag-op-defs.h likely could use const bools rather than functional argument operations. I haven't measured which is optimized more effectively. Also, the Octave binding layer in op-dm-scm.cc likely could use all sorts of macro or template trickery, but it's far easier to let Emacs handle it for now. That would be worth revisiting if further elementwise sparse and diagonal operations are added.

author	Jason Riedy <jason@acm.org>
date	Mon, 09 Mar 2009 17:49:14 -0400
parents	20dfb885f877
children	b03953732530

rev	line source
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	1 /*
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	2
7017 a1dbe9d80eee [project @ 2007-10-12 21:27:11 by jwe] jwe parents: 7016 diff changeset	3 Copyright (C) 1994, 1995, 1996, 1997, 2002, 2003, 2004, 2005, 2007
a1dbe9d80eee [project @ 2007-10-12 21:27:11 by jwe] jwe parents: 7016 diff changeset	4 John W. Eaton
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	5 Copyright (C) 2008, 2009 Jaroslav Hajek
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	6 Copyright (C) 2009 VZLU Prague
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	7
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	8 This file is part of Octave.
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	9
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	10 Octave is free software; you can redistribute it and/or modify it
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	11 under the terms of the GNU General Public License as published by the
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5307 diff changeset	12 Free Software Foundation; either version 3 of the License, or (at your
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5307 diff changeset	13 option) any later version.
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	14
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	15 Octave is distributed in the hope that it will be useful, but WITHOUT
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	17 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	18 for more details.
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	19
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	20 You should have received a copy of the GNU General Public License
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5307 diff changeset	21 along with Octave; see the file COPYING. If not, see
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5307 diff changeset	22 <http://www.gnu.org/licenses/>.
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	23
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	24 */
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	25
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	26 #ifdef HAVE_CONFIG_H
1192 b6360f2d4fa6 [project @ 1995-03-30 21:38:35 by jwe] jwe parents: 1011 diff changeset	27 #include <config.h>
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	28 #endif
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	29
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	30 #include "dbleQR.h"
1847 2ffe49eb95a5 [project @ 1996-02-03 12:47:55 by jwe] jwe parents: 1368 diff changeset	31 #include "f77-fcn.h"
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	32 #include "lo-error.h"
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	33 #include "Range.h"
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	34 #include "idx-vector.h"
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	35 #include "oct-locbuf.h"
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	36
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	37 extern "C"
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	38 {
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4192 diff changeset	39 F77_RET_T
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	40 F77_FUNC (dgeqrf, DGEQRF) (const octave_idx_type&, const octave_idx_type&, double*, const octave_idx_type&,
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	41 double, double, const octave_idx_type&, octave_idx_type&);
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	42
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4192 diff changeset	43 F77_RET_T
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	44 F77_FUNC (dorgqr, DORGQR) (const octave_idx_type&, const octave_idx_type&, const octave_idx_type&, double*,
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	45 const octave_idx_type&, double, double, const octave_idx_type&, octave_idx_type&);
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	46
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	47 #ifdef HAVE_QRUPDATE
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	48
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	49 F77_RET_T
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	50 F77_FUNC (dqr1up, DQR1UP) (const octave_idx_type&, const octave_idx_type&, const octave_idx_type&,
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	51 double, const octave_idx_type&, double, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	52 double, double, double*);
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	53
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	54 F77_RET_T
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	55 F77_FUNC (dqrinc, DQRINC) (const octave_idx_type&, const octave_idx_type&, const octave_idx_type&,
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	56 double, const octave_idx_type&, double, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	57 const octave_idx_type&, const double, double);
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	58
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	59 F77_RET_T
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	60 F77_FUNC (dqrdec, DQRDEC) (const octave_idx_type&, const octave_idx_type&, const octave_idx_type&,
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	61 double, const octave_idx_type&, double, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	62 const octave_idx_type&, double*);
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	63
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	64 F77_RET_T
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	65 F77_FUNC (dqrinr, DQRINR) (const octave_idx_type&, const octave_idx_type&,
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	66 double, const octave_idx_type&, double, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	67 const octave_idx_type&, const double, double);
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	68
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	69 F77_RET_T
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	70 F77_FUNC (dqrder, DQRDER) (const octave_idx_type&, const octave_idx_type&,
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	71 double, const octave_idx_type&, double, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	72 const octave_idx_type&, double*);
7700 efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	73
efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	74 F77_RET_T
efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	75 F77_FUNC (dqrshc, DQRSHC) (const octave_idx_type&, const octave_idx_type&, const octave_idx_type&,
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	76 double, const octave_idx_type&, double, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	77 const octave_idx_type&, const octave_idx_type&,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	78 double*);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	79
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	80 #endif
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	81 }
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	82
539 5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe] jwe parents: 457 diff changeset	83 QR::QR (const Matrix& a, QR::type qr_type)
2763 d9d00d7e271e [project @ 1997-03-01 02:14:33 by jwe] jwe parents: 2244 diff changeset	84 : q (), r ()
d9d00d7e271e [project @ 1997-03-01 02:14:33 by jwe] jwe parents: 2244 diff changeset	85 {
d9d00d7e271e [project @ 1997-03-01 02:14:33 by jwe] jwe parents: 2244 diff changeset	86 init (a, qr_type);
d9d00d7e271e [project @ 1997-03-01 02:14:33 by jwe] jwe parents: 2244 diff changeset	87 }
d9d00d7e271e [project @ 1997-03-01 02:14:33 by jwe] jwe parents: 2244 diff changeset	88
d9d00d7e271e [project @ 1997-03-01 02:14:33 by jwe] jwe parents: 2244 diff changeset	89 void
d9d00d7e271e [project @ 1997-03-01 02:14:33 by jwe] jwe parents: 2244 diff changeset	90 QR::init (const Matrix& a, QR::type qr_type)
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	91 {
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	92 octave_idx_type m = a.rows ();
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	93 octave_idx_type n = a.cols ();
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	94
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	95 octave_idx_type min_mn = m < n ? m : n;
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	96 OCTAVE_LOCAL_BUFFER (double, tau, min_mn);
1921 ce0db0ca0729 [project @ 1996-02-11 02:31:58 by jwe] jwe parents: 1882 diff changeset	97
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	98 octave_idx_type info = 0;
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	99
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	100 Matrix afact = a;
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	101 if (m > n && qr_type == QR::std)
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	102 afact.resize (m, m);
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	103
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	104 if (m > 0)
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	105 {
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	106 // workspace query.
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	107 double rlwork;
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	108 F77_XFCN (dgeqrf, DGEQRF, (m, n, afact.fortran_vec (), m, tau, &rlwork, -1, info));
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	109
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	110 // allocate buffer and do the job.
8811 20dfb885f877 int -> octave_idx fixes John W. Eaton <jwe@octave.org> parents: 8806 diff changeset	111 octave_idx_type lwork = rlwork;
20dfb885f877 int -> octave_idx fixes John W. Eaton <jwe@octave.org> parents: 8806 diff changeset	112 lwork = std::max (lwork, static_cast<octave_idx_type> (1));
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	113 OCTAVE_LOCAL_BUFFER (double, work, lwork);
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	114 F77_XFCN (dgeqrf, DGEQRF, (m, n, afact.fortran_vec (), m, tau, work, lwork, info));
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	115 }
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	116
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	117 form (n, afact, tau, qr_type);
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	118 }
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	119
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	120 void QR::form (octave_idx_type n, Matrix& afact,
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	121 double *tau, QR::type qr_type)
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	122 {
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	123 octave_idx_type m = afact.rows (), min_mn = std::min (m, n);
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	124 octave_idx_type info;
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	125
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7017 diff changeset	126 if (qr_type == QR::raw)
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7017 diff changeset	127 {
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7017 diff changeset	128 for (octave_idx_type j = 0; j < min_mn; j++)
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7017 diff changeset	129 {
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7017 diff changeset	130 octave_idx_type limit = j < min_mn - 1 ? j : min_mn - 1;
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7017 diff changeset	131 for (octave_idx_type i = limit + 1; i < m; i++)
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	132 afact.elem (i, j) *= tau[j];
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7017 diff changeset	133 }
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7017 diff changeset	134
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	135 r = afact;
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7017 diff changeset	136 }
539 5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe] jwe parents: 457 diff changeset	137 else
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	138 {
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	139 // Attempt to minimize copying.
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	140 if (m >= n)
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	141 {
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	142 // afact will become q.
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	143 q = afact;
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	144 octave_idx_type k = qr_type == QR::economy ? n : m;
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	145 r = Matrix (k, n);
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	146 for (octave_idx_type j = 0; j < n; j++)
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	147 {
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	148 octave_idx_type i = 0;
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	149 for (; i <= j; i++)
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	150 r.xelem (i, j) = afact.xelem (i, j);
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	151 for (;i < k; i++)
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	152 r.xelem (i, j) = 0;
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	153 }
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	154 afact = Matrix (); // optimize memory
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	155 }
539 5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe] jwe parents: 457 diff changeset	156 else
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	157 {
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	158 // afact will become r.
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	159 q = Matrix (m, m);
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	160 for (octave_idx_type j = 0; j < m; j++)
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	161 for (octave_idx_type i = j + 1; i < m; i++)
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	162 {
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	163 q.xelem (i, j) = afact.xelem (i, j);
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	164 afact.xelem (i, j) = 0;
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	165 }
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	166 r = afact;
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	167 }
1921 ce0db0ca0729 [project @ 1996-02-11 02:31:58 by jwe] jwe parents: 1882 diff changeset	168
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	169
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	170 if (m > 0)
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	171 {
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	172 octave_idx_type k = q.columns ();
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	173 // workspace query.
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	174 double rlwork;
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	175 F77_XFCN (dorgqr, DORGQR, (m, k, min_mn, q.fortran_vec (), m, tau,
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	176 &rlwork, -1, info));
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	177
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	178 // allocate buffer and do the job.
8811 20dfb885f877 int -> octave_idx fixes John W. Eaton <jwe@octave.org> parents: 8806 diff changeset	179 octave_idx_type lwork = rlwork;
20dfb885f877 int -> octave_idx fixes John W. Eaton <jwe@octave.org> parents: 8806 diff changeset	180 lwork = std::max (lwork, static_cast<octave_idx_type> (1));
8597 c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	181 OCTAVE_LOCAL_BUFFER (double, work, lwork);
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	182 F77_XFCN (dorgqr, DORGQR, (m, k, min_mn, q.fortran_vec (), m, tau,
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	183 work, lwork, info));
c86718093c1b improve & fix QR classes Jaroslav Hajek <highegg@gmail.com> parents: 8562 diff changeset	184 }
1921 ce0db0ca0729 [project @ 1996-02-11 02:31:58 by jwe] jwe parents: 1882 diff changeset	185 }
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	186 }
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	187
7924 4976f66d469b miscellaneous cleanup John W. Eaton <jwe@octave.org> parents: 7725 diff changeset	188 QR::QR (const Matrix& q_arg, const Matrix& r_arg)
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	189 {
8562 a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	190 octave_idx_type qr = q_arg.rows (), qc = q_arg.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	191 octave_idx_type rr = r_arg.rows (), rc = r_arg.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	192 if (qc == rr && (qr == qc \|\| (qr > qc && rr == rc)))
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	193 {
8562 a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	194 q = q_arg;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	195 r = r_arg;
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	196 }
8562 a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	197 else
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	198 (*current_liboctave_error_handler) ("QR dimensions mismatch");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	199 }
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	200
8562 a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	201 QR::type
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	202 QR::get_type (void) const
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	203 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	204 QR::type retval;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	205 if (!q.is_empty () && q.is_square ())
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	206 retval = QR::std;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	207 else if (q.rows () > q.columns () && r.is_square ())
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	208 retval = QR::economy;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	209 else
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	210 retval = QR::raw;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	211 return retval;
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	212 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	213
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	214 #ifdef HAVE_QRUPDATE
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	215
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	216 void
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	217 QR::update (const ColumnVector& u, const ColumnVector& v)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	218 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	219 octave_idx_type m = q.rows ();
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	220 octave_idx_type n = r.columns ();
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	221 octave_idx_type k = q.columns ();
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	222
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	223 if (u.length () == m && v.length () == n)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	224 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	225 ColumnVector utmp = u, vtmp = v;
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	226 OCTAVE_LOCAL_BUFFER (double, w, 2*k);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	227 F77_XFCN (dqr1up, DQR1UP, (m, n, k, q.fortran_vec (), m, r.fortran_vec (), k,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	228 utmp.fortran_vec (), vtmp.fortran_vec (), w));
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	229 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	230 else
8562 a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	231 (*current_liboctave_error_handler) ("qrupdate: dimensions mismatch");
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	232 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	233
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	234 void
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	235 QR::update (const Matrix& u, const Matrix& v)
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	236 {
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	237 octave_idx_type m = q.rows ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	238 octave_idx_type n = r.columns ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	239 octave_idx_type k = q.columns ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	240
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	241 if (u.rows () == m && v.rows () == n && u.cols () == v.cols ())
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	242 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	243 OCTAVE_LOCAL_BUFFER (double, w, 2*k);
8806 c7864bb74914 avoid some GCC warnings John W. Eaton <jwe@octave.org> parents: 8721 diff changeset	244 for (volatile octave_idx_type i = 0; i < u.cols (); i++)
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	245 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	246 ColumnVector utmp = u.column (i), vtmp = v.column (i);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	247 F77_XFCN (dqr1up, DQR1UP, (m, n, k, q.fortran_vec (), m, r.fortran_vec (), k,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	248 utmp.fortran_vec (), vtmp.fortran_vec (), w));
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	249 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	250 }
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	251 else
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	252 (*current_liboctave_error_handler) ("qrupdate: dimensions mismatch");
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	253 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	254
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	255 void
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	256 QR::insert_col (const ColumnVector& u, octave_idx_type j)
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	257 {
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	258 octave_idx_type m = q.rows ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	259 octave_idx_type n = r.columns ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	260 octave_idx_type k = q.columns ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	261
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	262 if (u.length () != m)
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	263 (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
7560 0ef0f9802a37 modify QR updating methods to use 0-based indexing Jaroslav Hajek <highegg@gmail.com> parents: 7559 diff changeset	264 else if (j < 0 \|\| j > n)
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	265 (*current_liboctave_error_handler) ("qrinsert: index out of range");
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	266 else
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	267 {
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	268 if (k < m)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	269 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	270 q.resize (m, k+1);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	271 r.resize (k+1, n+1);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	272 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	273 else
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	274 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	275 r.resize (k, n+1);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	276 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	277
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	278 ColumnVector utmp = u;
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	279 OCTAVE_LOCAL_BUFFER (double, w, k);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	280 F77_XFCN (dqrinc, DQRINC, (m, n, k, q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	281 r.fortran_vec (), r.rows (), j + 1,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	282 utmp.data (), w));
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	283 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	284 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	285
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	286 void
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	287 QR::insert_col (const Matrix& u, const Array<octave_idx_type>& j)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	288 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	289 octave_idx_type m = q.rows ();
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	290 octave_idx_type n = r.columns ();
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	291 octave_idx_type k = q.columns ();
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	292
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	293 Array<octave_idx_type> jsi;
8721 e9cb742df9eb imported patch sort3.diff Jaroslav Hajek <highegg@gmail.com> parents: 8597 diff changeset	294 Array<octave_idx_type> js = j.sort (jsi, 0, ASCENDING);
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	295 octave_idx_type nj = js.length ();
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	296 bool dups = false;
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	297 for (octave_idx_type i = 0; i < nj - 1; i++)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	298 dups = dups && js(i) == js(i+1);
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	299
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	300 if (dups)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	301 (*current_liboctave_error_handler) ("qrinsert: duplicate index detected");
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	302 else if (u.length () != m \|\| u.columns () != nj)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	303 (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	304 else if (nj > 0 && (js(0) < 0 \|\| js(nj-1) > n))
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	305 (*current_liboctave_error_handler) ("qrinsert: index out of range");
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	306 else if (nj > 0)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	307 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	308 octave_idx_type kmax = std::min (k + nj, m);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	309 if (k < m)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	310 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	311 q.resize (m, kmax);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	312 r.resize (kmax, n + nj);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	313 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	314 else
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	315 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	316 r.resize (k, n + nj);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	317 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	318
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	319 OCTAVE_LOCAL_BUFFER (double, w, kmax);
8806 c7864bb74914 avoid some GCC warnings John W. Eaton <jwe@octave.org> parents: 8721 diff changeset	320 for (volatile octave_idx_type i = 0; i < js.length (); i++)
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	321 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	322 ColumnVector utmp = u.column (jsi(i));
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	323 F77_XFCN (dqrinc, DQRINC, (m, n + i, std::min (kmax, k + i),
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	324 q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	325 r.fortran_vec (), r.rows (), js(i) + 1,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	326 utmp.data (), w));
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	327 }
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	328 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	329 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	330
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	331 void
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	332 QR::delete_col (octave_idx_type j)
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	333 {
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	334 octave_idx_type m = q.rows ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	335 octave_idx_type k = r.rows ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	336 octave_idx_type n = r.columns ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	337
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	338 if (j < 0 \|\| j > n-1)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	339 (*current_liboctave_error_handler) ("qrdelete: index out of range");
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	340 else
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	341 {
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	342 OCTAVE_LOCAL_BUFFER (double, w, k);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	343 F77_XFCN (dqrdec, DQRDEC, (m, n, k, q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	344 r.fortran_vec (), r.rows (), j + 1, w));
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	345
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	346 if (k < m)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	347 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	348 q.resize (m, k-1);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	349 r.resize (k-1, n-1);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	350 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	351 else
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	352 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	353 r.resize (k, n-1);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	354 }
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	355 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	356 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	357
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	358 void
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	359 QR::delete_col (const Array<octave_idx_type>& j)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	360 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	361 octave_idx_type m = q.rows ();
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	362 octave_idx_type n = r.columns ();
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	363 octave_idx_type k = q.columns ();
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	364
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	365 Array<octave_idx_type> jsi;
8721 e9cb742df9eb imported patch sort3.diff Jaroslav Hajek <highegg@gmail.com> parents: 8597 diff changeset	366 Array<octave_idx_type> js = j.sort (jsi, 0, DESCENDING);
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	367 octave_idx_type nj = js.length ();
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	368 bool dups = false;
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	369 for (octave_idx_type i = 0; i < nj - 1; i++)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	370 dups = dups && js(i) == js(i+1);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	371
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	372 if (dups)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	373 (*current_liboctave_error_handler) ("qrinsert: duplicate index detected");
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	374 else if (nj > 0 && (js(0) > n-1 \|\| js(nj-1) < 0))
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	375 (*current_liboctave_error_handler) ("qrinsert: index out of range");
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	376 else if (nj > 0)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	377 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	378 OCTAVE_LOCAL_BUFFER (double, w, k);
8806 c7864bb74914 avoid some GCC warnings John W. Eaton <jwe@octave.org> parents: 8721 diff changeset	379 for (volatile octave_idx_type i = 0; i < js.length (); i++)
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	380 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	381 F77_XFCN (dqrdec, DQRDEC, (m, n - i, k == m ? k : k - i,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	382 q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	383 r.fortran_vec (), r.rows (), js(i) + 1, w));
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	384 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	385 if (k < m)
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	386 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	387 q.resize (m, k - nj);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	388 r.resize (k - nj, n - nj);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	389 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	390 else
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	391 {
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	392 r.resize (k, n - nj);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	393 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	394
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	395 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	396 }
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	397
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	398 void
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	399 QR::insert_row (const RowVector& u, octave_idx_type j)
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	400 {
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	401 octave_idx_type m = r.rows ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	402 octave_idx_type n = r.columns ();
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	403 octave_idx_type k = std::min (m, n);
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	404
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	405 if (! q.is_square () \|\| u.length () != n)
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	406 (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
7560 0ef0f9802a37 modify QR updating methods to use 0-based indexing Jaroslav Hajek <highegg@gmail.com> parents: 7559 diff changeset	407 else if (j < 0 \|\| j > m)
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	408 (*current_liboctave_error_handler) ("qrinsert: index out of range");
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	409 else
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	410 {
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	411 q.resize (m + 1, m + 1);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	412 r.resize (m + 1, n);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	413 RowVector utmp = u;
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	414 OCTAVE_LOCAL_BUFFER (double, w, k);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	415 F77_XFCN (dqrinr, DQRINR, (m, n, q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	416 r.fortran_vec (), r.rows (),
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	417 j + 1, utmp.fortran_vec (), w));
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	418
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	419 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	420 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	421
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	422 void
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	423 QR::delete_row (octave_idx_type j)
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	424 {
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	425 octave_idx_type m = r.rows ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	426 octave_idx_type n = r.columns ();
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	427
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	428 if (! q.is_square ())
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	429 (*current_liboctave_error_handler) ("qrdelete: dimensions mismatch");
7560 0ef0f9802a37 modify QR updating methods to use 0-based indexing Jaroslav Hajek <highegg@gmail.com> parents: 7559 diff changeset	430 else if (j < 0 \|\| j > m-1)
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	431 (*current_liboctave_error_handler) ("qrdelete: index out of range");
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	432 else
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	433 {
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	434 OCTAVE_LOCAL_BUFFER (double, w, 2*m);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	435 F77_XFCN (dqrder, DQRDER, (m, n, q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	436 r.fortran_vec (), r.rows (), j + 1,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	437 w));
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	438
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	439 q.resize (m - 1, m - 1);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	440 r.resize (m - 1, n);
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	441 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	442 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	443
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	444 void
7700 efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	445 QR::shift_cols (octave_idx_type i, octave_idx_type j)
efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	446 {
efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	447 octave_idx_type m = q.rows ();
efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	448 octave_idx_type k = r.rows ();
efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	449 octave_idx_type n = r.columns ();
efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	450
efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	451 if (i < 0 \|\| i > n-1 \|\| j < 0 \|\| j > n-1)
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	452 (*current_liboctave_error_handler) ("qrshift: index out of range");
7700 efccca5f2ad7 more QR & Cholesky updating functions Jaroslav Hajek <highegg@gmail.com> parents: 7560 diff changeset	453 else
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	454 {
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	455 OCTAVE_LOCAL_BUFFER (double, w, 2*k);
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	456 F77_XFCN (dqrshc, DQRSHC, (m, n, k,
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	457 q.fortran_vec (), q.rows (),
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	458 r.fortran_vec (), r.rows (),
d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	459 i + 1, j + 1, w));
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	460 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	461 }
56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	462
8562 a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	463 #else
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	464
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	465 // Replacement update methods.
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	466
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	467 void
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	468 QR::update (const ColumnVector& u, const ColumnVector& v)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	469 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	470 warn_qrupdate_once ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	471
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	472 octave_idx_type m = q.rows ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	473 octave_idx_type n = r.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	474
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	475 if (u.length () == m && v.length () == n)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	476 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	477 init(qr + Matrix (u) Matrix (v).transpose (), get_type ());
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	478 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	479 else
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	480 (*current_liboctave_error_handler) ("qrupdate: dimensions mismatch");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	481 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	482
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	483 void
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	484 QR::update (const Matrix& u, const Matrix& v)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	485 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	486 warn_qrupdate_once ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	487
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	488 octave_idx_type m = q.rows ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	489 octave_idx_type n = r.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	490
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	491 if (u.rows () == m && v.rows () == n && u.cols () == v.cols ())
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	492 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	493 init(qr + u v.transpose (), get_type ());
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	494 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	495 else
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	496 (*current_liboctave_error_handler) ("qrupdate: dimensions mismatch");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	497 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	498
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	499 static
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	500 Matrix insert_col (const Matrix& a, octave_idx_type i,
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	501 const ColumnVector& x)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	502 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	503 Matrix retval (a.rows (), a.columns () + 1);
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	504 retval.assign (idx_vector::colon, idx_vector (0, i),
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	505 a.index (idx_vector::colon, idx_vector (0, i)));
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	506 retval.assign (idx_vector::colon, idx_vector (i), x);
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	507 retval.assign (idx_vector::colon, idx_vector (i+1, retval.columns ()),
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	508 a.index (idx_vector::colon, idx_vector (i, a.columns ())));
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	509 return retval;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	510 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	511
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	512 static
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	513 Matrix insert_row (const Matrix& a, octave_idx_type i,
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	514 const RowVector& x)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	515 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	516 Matrix retval (a.rows () + 1, a.columns ());
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	517 retval.assign (idx_vector (0, i), idx_vector::colon,
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	518 a.index (idx_vector (0, i), idx_vector::colon));
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	519 retval.assign (idx_vector (i), idx_vector::colon, x);
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	520 retval.assign (idx_vector (i+1, retval.rows ()), idx_vector::colon,
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	521 a.index (idx_vector (i, a.rows ()), idx_vector::colon));
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	522 return retval;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	523 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	524
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	525 static
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	526 Matrix delete_col (const Matrix& a, octave_idx_type i)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	527 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	528 Matrix retval = a;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	529 retval.delete_elements (1, idx_vector (i));
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	530 return retval;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	531 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	532
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	533 static
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	534 Matrix delete_row (const Matrix& a, octave_idx_type i)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	535 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	536 Matrix retval = a;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	537 retval.delete_elements (0, idx_vector (i));
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	538 return retval;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	539 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	540
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	541 static
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	542 Matrix shift_cols (const Matrix& a,
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	543 octave_idx_type i, octave_idx_type j)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	544 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	545 octave_idx_type n = a.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	546 Array<octave_idx_type> p (n);
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	547 for (octave_idx_type k = 0; k < n; k++) p(k) = k;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	548 if (i < j)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	549 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	550 for (octave_idx_type k = i; k < j; k++) p(k) = k+1;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	551 p(j) = i;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	552 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	553 else if (j < i)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	554 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	555 p(j) = i;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	556 for (octave_idx_type k = j+1; k < i+1; k++) p(k) = k-1;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	557 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	558
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	559 return a.index (idx_vector::colon, idx_vector (p));
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	560 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	561
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	562 void
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	563 QR::insert_col (const ColumnVector& u, octave_idx_type j)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	564 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	565 warn_qrupdate_once ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	566
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	567 octave_idx_type m = q.rows ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	568 octave_idx_type n = r.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	569
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	570 if (u.length () != m)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	571 (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	572 else if (j < 0 \|\| j > n)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	573 (*current_liboctave_error_handler) ("qrinsert: index out of range");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	574 else
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	575 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	576 init (::insert_col (q*r, j, u), get_type ());
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	577 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	578 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	579
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	580 void
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	581 QR::insert_col (const Matrix& u, const Array<octave_idx_type>& j)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	582 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	583 warn_qrupdate_once ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	584
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	585 octave_idx_type m = q.rows ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	586 octave_idx_type n = r.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	587
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	588 Array<octave_idx_type> jsi;
8721 e9cb742df9eb imported patch sort3.diff Jaroslav Hajek <highegg@gmail.com> parents: 8597 diff changeset	589 Array<octave_idx_type> js = j.sort (jsi, 0, ASCENDING);
8562 a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	590 octave_idx_type nj = js.length ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	591 bool dups = false;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	592 for (octave_idx_type i = 0; i < nj - 1; i++)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	593 dups = dups && js(i) == js(i+1);
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	594
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	595 if (dups)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	596 (*current_liboctave_error_handler) ("qrinsert: duplicate index detected");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	597 else if (u.length () != m \|\| u.columns () != nj)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	598 (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	599 else if (nj > 0 && (js(0) < 0 \|\| js(nj-1) > n))
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	600 (*current_liboctave_error_handler) ("qrinsert: index out of range");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	601 else if (nj > 0)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	602 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	603 Matrix a = q*r;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	604 for (octave_idx_type i = 0; i < js.length (); i++)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	605 a = ::insert_col (a, js(i), u.column (i));
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	606 init (a, get_type ());
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	607 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	608 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	609
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	610 void
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	611 QR::delete_col (octave_idx_type j)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	612 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	613 warn_qrupdate_once ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	614
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	615 octave_idx_type m = q.rows ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	616 octave_idx_type n = r.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	617
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	618 if (j < 0 \|\| j > n-1)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	619 (*current_liboctave_error_handler) ("qrdelete: index out of range");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	620 else
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	621 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	622 init (::delete_col (q*r, j), get_type ());
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	623 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	624 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	625
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	626 void
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	627 QR::delete_col (const Array<octave_idx_type>& j)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	628 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	629 warn_qrupdate_once ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	630
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	631 octave_idx_type m = q.rows ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	632 octave_idx_type n = r.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	633
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	634 Array<octave_idx_type> jsi;
8721 e9cb742df9eb imported patch sort3.diff Jaroslav Hajek <highegg@gmail.com> parents: 8597 diff changeset	635 Array<octave_idx_type> js = j.sort (jsi, 0, DESCENDING);
8562 a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	636 octave_idx_type nj = js.length ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	637 bool dups = false;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	638 for (octave_idx_type i = 0; i < nj - 1; i++)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	639 dups = dups && js(i) == js(i+1);
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	640
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	641 if (dups)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	642 (*current_liboctave_error_handler) ("qrinsert: duplicate index detected");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	643 else if (nj > 0 && (js(0) > n-1 \|\| js(nj-1) < 0))
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	644 (*current_liboctave_error_handler) ("qrinsert: index out of range");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	645 else if (nj > 0)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	646 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	647 Matrix a = q*r;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	648 for (octave_idx_type i = 0; i < js.length (); i++)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	649 a = ::delete_col (a, js(i));
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	650 init (a, get_type ());
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	651 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	652 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	653
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	654 void
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	655 QR::insert_row (const RowVector& u, octave_idx_type j)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	656 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	657 warn_qrupdate_once ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	658
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	659 octave_idx_type m = r.rows ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	660 octave_idx_type n = r.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	661
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	662 if (! q.is_square () \|\| u.length () != n)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	663 (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	664 else if (j < 0 \|\| j > m)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	665 (*current_liboctave_error_handler) ("qrinsert: index out of range");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	666 else
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	667 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	668 init (::insert_row (q*r, j, u), get_type ());
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	669 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	670 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	671
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	672 void
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	673 QR::delete_row (octave_idx_type j)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	674 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	675 octave_idx_type m = r.rows ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	676 octave_idx_type n = r.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	677
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	678 if (! q.is_square ())
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	679 (*current_liboctave_error_handler) ("qrdelete: dimensions mismatch");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	680 else if (j < 0 \|\| j > m-1)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	681 (*current_liboctave_error_handler) ("qrdelete: index out of range");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	682 else
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	683 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	684 init (::delete_row (q*r, j), get_type ());
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	685 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	686 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	687
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	688 void
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	689 QR::shift_cols (octave_idx_type i, octave_idx_type j)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	690 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	691 warn_qrupdate_once ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	692
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	693 octave_idx_type m = q.rows ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	694 octave_idx_type n = r.columns ();
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	695
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	696 if (i < 0 \|\| i > n-1 \|\| j < 0 \|\| j > n-1)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	697 (*current_liboctave_error_handler) ("qrshift: index out of range");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	698 else
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	699 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	700 init (::shift_cols (q*r, i, j), get_type ());
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	701 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	702 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	703
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	704 void warn_qrupdate_once (void)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	705 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	706 static bool warned = false;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	707 if (! warned)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	708 {
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	709 (*current_liboctave_warning_handler)
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	710 ("In this version of Octave, QR & Cholesky updating routines\n"
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	711 "simply update the matrix and recalculate factorizations.\n"
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	712 "To use fast algorithms, link Octave with the qrupdate library.\n"
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	713 "See <http://sourceforge.net/projects/qrupdate>.\n");
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	714 warned = true;
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	715 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	716 }
a6edd5c23cb5 use replacement methods if qrupdate is not available Jaroslav Hajek <highegg@gmail.com> parents: 8547 diff changeset	717
8547 d66c9b6e506a imported patch qrupdate.diff Jaroslav Hajek <highegg@gmail.com> parents: 8367 diff changeset	718 #endif
7553 56be6f31dd4e implementation of QR factorization updating Jaroslav Hajek <highegg@gmail.com> parents: 7482 diff changeset	719
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	720 /*
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	721 ;;; Local Variables: ***
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	722 ;;; mode: C++ ***
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	723 ;;; End: ***
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	724 */

Mercurial > hg > octave-lyh

annotate liboctave/dbleQR.cc @ 8966:1bba53c0a38d