octave-nkf: liboctave/SparseQR.cc annotate

annotate liboctave/SparseQR.cc @ 12312:b10ea6efdc58 release-3-4-x ss-3-3-91

version is now 3.3.91

author	John W. Eaton <jwe@octave.org>
date	Mon, 31 Jan 2011 08:36:58 -0500
parents	12df7854fa7c
children	72c96de7a403

rev	line source
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	1 /*
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	2
11523 fd0a3ac60b0e update copyright notices John W. Eaton <jwe@octave.org> parents: 11510 diff changeset	3 Copyright (C) 2005-2011 David Bateman
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	4
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 6924 diff changeset	5 This file is part of Octave.
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 6924 diff changeset	6
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	7 Octave is free software; you can redistribute it and/or modify it
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	8 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: 6924 diff changeset	9 Free Software Foundation; either version 3 of the License, or (at your
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 6924 diff changeset	10 option) any later version.
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	11
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	12 Octave is distributed in the hope that it will be useful, but WITHOUT
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	15 for more details.
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	16
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	17 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: 6924 diff changeset	18 along with Octave; see the file COPYING. If not, see
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 6924 diff changeset	19 <http://www.gnu.org/licenses/>.
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	20
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	21 */
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	22
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	23 #ifdef HAVE_CONFIG_H
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	24 #include <config.h>
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	25 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	26 #include <vector>
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	27
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	28 #include "lo-error.h"
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	29 #include "SparseQR.h"
8377 25bc2d31e1bf improve OCTAVE_LOCAL_BUFFER Jaroslav Hajek <highegg@gmail.com> parents: 7789 diff changeset	30 #include "oct-locbuf.h"
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	31
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	32 SparseQR::SparseQR_rep::SparseQR_rep (const SparseMatrix& a, int order)
11528 feb90c7cfc16 Build successfully without HAVE_CXSPARSE David Grundberg <individ@acc.umu.se> parents: 11523 diff changeset	33 : count (1), nrows (0)
feb90c7cfc16 Build successfully without HAVE_CXSPARSE David Grundberg <individ@acc.umu.se> parents: 11523 diff changeset	34 #ifdef HAVE_CXSPARSE
feb90c7cfc16 Build successfully without HAVE_CXSPARSE David Grundberg <individ@acc.umu.se> parents: 11523 diff changeset	35 , S (0), N (0)
feb90c7cfc16 Build successfully without HAVE_CXSPARSE David Grundberg <individ@acc.umu.se> parents: 11523 diff changeset	36 #endif
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	37 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	38 #ifdef HAVE_CXSPARSE
5648 69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	39 CXSPARSE_DNAME () A;
10506 bdf5d85cfc5e replace nzmax by nnz where appropriate in liboctave Jaroslav Hajek <highegg@gmail.com> parents: 10350 diff changeset	40 A.nzmax = a.nnz ();
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	41 A.m = a.rows ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	42 A.n = a.cols ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	43 nrows = A.m;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	44 // Cast away const on A, with full knowledge that CSparse won't touch it
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	45 // Prevents the methods below making a copy of the data.
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	46 A.p = const_cast<octave_idx_type *>(a.cidx ());
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	47 A.i = const_cast<octave_idx_type *>(a.ridx ());
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	48 A.x = const_cast<double *>(a.data ());
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	49 A.nz = -1;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	50 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	51 #if defined(CS_VER) && (CS_VER >= 2)
eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	52 S = CXSPARSE_DNAME (_sqr) (order, &A, 1);
eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	53 #else
eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	54 S = CXSPARSE_DNAME (_sqr) (&A, order - 1, 1);
eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	55 #endif
eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	56
5648 69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	57 N = CXSPARSE_DNAME (_qr) (&A, S);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	58 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	59 if (!N)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	60 (*current_liboctave_error_handler)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	61 ("SparseQR: sparse matrix QR factorization filled");
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	62 count = 1;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	63 #else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	64 (*current_liboctave_error_handler)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	65 ("SparseQR: sparse matrix QR factorization not implemented");
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	66 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	67 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	68
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	69 SparseQR::SparseQR_rep::~SparseQR_rep (void)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	70 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	71 #ifdef HAVE_CXSPARSE
5648 69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	72 CXSPARSE_DNAME (_sfree) (S);
69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	73 CXSPARSE_DNAME (_nfree) (N);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	74 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	75 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	76
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	77 SparseMatrix
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	78 SparseQR::SparseQR_rep::V (void) const
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	79 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	80 #ifdef HAVE_CXSPARSE
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	81 // Drop zeros from V and sort
5775 ace8d8d26933 [project @ 2006-04-24 19:13:06 by jwe] jwe parents: 5681 diff changeset	82 // FIXME Is the double transpose to sort necessary?
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	83 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5648 69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	84 CXSPARSE_DNAME (_dropzeros) (N->L);
69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	85 CXSPARSE_DNAME () *D = CXSPARSE_DNAME (_transpose) (N->L, 1);
69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	86 CXSPARSE_DNAME (_spfree) (N->L);
69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	87 N->L = CXSPARSE_DNAME (_transpose) (D, 1);
69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	88 CXSPARSE_DNAME (_spfree) (D);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	89 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	90
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	91 octave_idx_type nc = N->L->n;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	92 octave_idx_type nz = N->L->nzmax;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	93 SparseMatrix ret (N->L->m, nc, nz);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	94 for (octave_idx_type j = 0; j < nc+1; j++)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	95 ret.xcidx (j) = N->L->p[j];
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	96 for (octave_idx_type j = 0; j < nz; j++)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	97 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	98 ret.xridx (j) = N->L->i[j];
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	99 ret.xdata (j) = N->L->x[j];
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	100 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	101 return ret;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	102 #else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	103 return SparseMatrix ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	104 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	105 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	106
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	107 ColumnVector
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	108 SparseQR::SparseQR_rep::Pinv (void) const
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	109 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	110 #ifdef HAVE_CXSPARSE
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	111 ColumnVector ret(N->L->m);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	112 for (octave_idx_type i = 0; i < N->L->m; i++)
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	113 #if defined(CS_VER) && (CS_VER >= 2)
eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	114 ret.xelem(i) = S->pinv[i];
eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	115 #else
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	116 ret.xelem(i) = S->Pinv[i];
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	117 #endif
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	118 return ret;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	119 #else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	120 return ColumnVector ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	121 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	122 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	123
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	124 ColumnVector
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	125 SparseQR::SparseQR_rep::P (void) const
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	126 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	127 #ifdef HAVE_CXSPARSE
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	128 ColumnVector ret(N->L->m);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	129 for (octave_idx_type i = 0; i < N->L->m; i++)
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	130 #if defined(CS_VER) && (CS_VER >= 2)
eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	131 ret.xelem(S->pinv[i]) = i;
eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	132 #else
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	133 ret.xelem(S->Pinv[i]) = i;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	134 #endif
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	135 return ret;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	136 #else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	137 return ColumnVector ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	138 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	139 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	140
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	141 SparseMatrix
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	142 SparseQR::SparseQR_rep::R (const bool econ) const
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	143 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	144 #ifdef HAVE_CXSPARSE
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	145 // Drop zeros from R and sort
5775 ace8d8d26933 [project @ 2006-04-24 19:13:06 by jwe] jwe parents: 5681 diff changeset	146 // FIXME Is the double transpose to sort necessary?
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	147 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5648 69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	148 CXSPARSE_DNAME (_dropzeros) (N->U);
69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	149 CXSPARSE_DNAME () *D = CXSPARSE_DNAME (_transpose) (N->U, 1);
69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	150 CXSPARSE_DNAME (_spfree) (N->U);
69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	151 N->U = CXSPARSE_DNAME (_transpose) (D, 1);
69a4f320d95a [project @ 2006-03-08 20:17:37 by dbateman] dbateman parents: 5610 diff changeset	152 CXSPARSE_DNAME (_spfree) (D);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	153 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	154
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	155 octave_idx_type nc = N->U->n;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	156 octave_idx_type nz = N->U->nzmax;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	157
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	158 SparseMatrix ret ((econ ? (nc > nrows ? nrows : nc) : nrows), nc, nz);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	159
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	160 for (octave_idx_type j = 0; j < nc+1; j++)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	161 ret.xcidx (j) = N->U->p[j];
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	162 for (octave_idx_type j = 0; j < nz; j++)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	163 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	164 ret.xridx (j) = N->U->i[j];
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	165 ret.xdata (j) = N->U->x[j];
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	166 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	167 return ret;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	168 #else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	169 return SparseMatrix ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	170 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	171 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	172
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	173 Matrix
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	174 SparseQR::SparseQR_rep::C (const Matrix &b) const
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	175 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	176 #ifdef HAVE_CXSPARSE
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	177 octave_idx_type b_nr = b.rows();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	178 octave_idx_type b_nc = b.cols();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	179 octave_idx_type nc = N->L->n;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	180 octave_idx_type nr = nrows;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	181 const double *bvec = b.fortran_vec();
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	182 Matrix ret (b_nr, b_nc);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	183 double *vec = ret.fortran_vec();
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	184 if (nr < 0 \|\| nc < 0 \|\| nr != b_nr)
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	185 (*current_liboctave_error_handler) ("matrix dimension mismatch");
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	186 else if (nr == 0 \|\| nc == 0 \|\| b_nc == 0)
be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	187 ret = Matrix (nc, b_nc, 0.0);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	188 else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	189 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	190 OCTAVE_LOCAL_BUFFER (double, buf, S->m2);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	191 for (volatile octave_idx_type j = 0, idx = 0; j < b_nc; j++, idx+=b_nr)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	192 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	193 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	194 for (octave_idx_type i = nr; i < S->m2; i++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	195 buf[i] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	196 volatile octave_idx_type nm = (nr < nc ? nr : nc);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	197 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	198 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	199 CXSPARSE_DNAME (_ipvec) (S->pinv, bvec + idx, buf, b_nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	200 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	201 CXSPARSE_DNAME (_ipvec) (b_nr, S->Pinv, bvec + idx, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	202 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	203 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	204
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	205 for (volatile octave_idx_type i = 0; i < nm; i++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	206 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	207 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	208 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	209 CXSPARSE_DNAME (_happly) (N->L, i, N->B[i], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	210 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	211 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	212 for (octave_idx_type i = 0; i < b_nr; i++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	213 vec[i+idx] = buf[i];
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	214 }
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	215 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	216 return ret;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	217 #else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	218 return Matrix ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	219 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	220 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	221
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	222 Matrix
7505 f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	223 SparseQR::SparseQR_rep::Q (void) const
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	224 {
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	225 #ifdef HAVE_CXSPARSE
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	226 octave_idx_type nc = N->L->n;
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	227 octave_idx_type nr = nrows;
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	228 Matrix ret (nr, nr);
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	229 double *vec = ret.fortran_vec();
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	230 if (nr < 0 \|\| nc < 0)
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	231 (*current_liboctave_error_handler) ("matrix dimension mismatch");
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	232 else if (nr == 0 \|\| nc == 0)
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	233 ret = Matrix (nc, nr, 0.0);
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	234 else
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	235 {
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	236 OCTAVE_LOCAL_BUFFER (double, bvec, nr + 1);
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	237 for (octave_idx_type i = 0; i < nr; i++)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	238 bvec[i] = 0.;
7505 f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	239 OCTAVE_LOCAL_BUFFER (double, buf, S->m2);
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	240 for (volatile octave_idx_type j = 0, idx = 0; j < nr; j++, idx+=nr)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	241 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	242 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	243 bvec[j] = 1.0;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	244 for (octave_idx_type i = nr; i < S->m2; i++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	245 buf[i] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	246 volatile octave_idx_type nm = (nr < nc ? nr : nc);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	247 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
7505 f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	248 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	249 CXSPARSE_DNAME (_ipvec) (S->pinv, bvec, buf, nr);
7505 f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	250 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	251 CXSPARSE_DNAME (_ipvec) (nr, S->Pinv, bvec, buf);
7505 f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	252 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	253 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
7505 f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	254
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	255 for (volatile octave_idx_type i = 0; i < nm; i++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	256 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	257 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	258 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	259 CXSPARSE_DNAME (_happly) (N->L, i, N->B[i], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	260 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	261 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	262 for (octave_idx_type i = 0; i < nr; i++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	263 vec[i+idx] = buf[i];
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	264 bvec[j] = 0.0;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	265 }
7505 f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	266 }
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	267 return ret.transpose ();
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	268 #else
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	269 return Matrix ();
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	270 #endif
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	271 }
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	272
f5005d9510f4 Remove dispatched sparse functions and treat in the generic versions of the functions David Bateman <dbateman@free.fr> parents: 7017 diff changeset	273 Matrix
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	274 qrsolve(const SparseMatrix&a, const Matrix &b, octave_idx_type& info)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	275 {
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	276 info = -1;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	277 #ifdef HAVE_CXSPARSE
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	278 octave_idx_type nr = a.rows();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	279 octave_idx_type nc = a.cols();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	280 octave_idx_type b_nc = b.cols();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	281 octave_idx_type b_nr = b.rows();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	282 const double *bvec = b.fortran_vec();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	283 Matrix x;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	284
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	285 if (nr < 0 \|\| nc < 0 \|\| nr != b_nr)
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	286 (*current_liboctave_error_handler)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	287 ("matrix dimension mismatch in solution of minimum norm problem");
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	288 else if (nr == 0 \|\| nc == 0 \|\| b_nc == 0)
be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	289 x = Matrix (nc, b_nc, 0.0);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	290 else if (nr >= nc)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	291 {
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	292 SparseQR q (a, 3);
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	293 if (! q.ok ())
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	294 return Matrix();
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	295 x.resize(nc, b_nc);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	296 double *vec = x.fortran_vec();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	297 OCTAVE_LOCAL_BUFFER (double, buf, q.S()->m2);
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	298 for (volatile octave_idx_type i = 0, idx = 0, bidx = 0; i < b_nc;
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	299 i++, idx+=nc, bidx+=b_nr)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	300 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	301 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	302 for (octave_idx_type j = nr; j < q.S()->m2; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	303 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	304 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	305 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	306 CXSPARSE_DNAME (_ipvec) (q.S()->pinv, bvec + bidx, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	307 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	308 CXSPARSE_DNAME (_ipvec) (nr, q.S()->Pinv, bvec + bidx, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	309 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	310 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	311 for (volatile octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	312 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	313 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	314 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	315 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	316 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	317 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	318 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	319 CXSPARSE_DNAME (_usolve) (q.N()->U, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	320 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	321 CXSPARSE_DNAME (_ipvec) (q.S()->q, buf, vec + idx, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	322 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	323 CXSPARSE_DNAME (_ipvec) (nc, q.S()->Q, buf, vec + idx);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	324 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	325 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	326 }
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	327 info = 0;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	328 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	329 else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	330 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	331 SparseMatrix at = a.hermitian();
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	332 SparseQR q (at, 3);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	333 if (! q.ok ())
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	334 return Matrix();
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	335 x.resize(nc, b_nc);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	336 double *vec = x.fortran_vec();
5681 233d98d95659 [project @ 2006-03-16 17:48:55 by dbateman] dbateman parents: 5648 diff changeset	337 volatile octave_idx_type nbuf = (nc > q.S()->m2 ? nc : q.S()->m2);
233d98d95659 [project @ 2006-03-16 17:48:55 by dbateman] dbateman parents: 5648 diff changeset	338 OCTAVE_LOCAL_BUFFER (double, buf, nbuf);
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	339 for (volatile octave_idx_type i = 0, idx = 0, bidx = 0; i < b_nc;
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	340 i++, idx+=nc, bidx+=b_nr)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	341 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	342 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	343 for (octave_idx_type j = nr; j < nbuf; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	344 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	345 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	346 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	347 CXSPARSE_DNAME (_pvec) (q.S()->q, bvec + bidx, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	348 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	349 CXSPARSE_DNAME (_pvec) (nr, q.S()->Q, bvec + bidx, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	350 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	351 CXSPARSE_DNAME (_utsolve) (q.N()->U, buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	352 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	353 for (volatile octave_idx_type j = nr-1; j >= 0; j--)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	354 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	355 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	356 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	357 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	358 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	359 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	360 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	361 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	362 CXSPARSE_DNAME (_pvec) (q.S()->pinv, buf, vec + idx, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	363 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	364 CXSPARSE_DNAME (_pvec) (nc, q.S()->Pinv, buf, vec + idx);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	365 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	366 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	367 }
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	368 info = 0;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	369 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	370
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	371 return x;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	372 #else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	373 return Matrix ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	374 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	375 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	376
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	377 SparseMatrix
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	378 qrsolve(const SparseMatrix&a, const SparseMatrix &b, octave_idx_type &info)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	379 {
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	380 info = -1;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	381 #ifdef HAVE_CXSPARSE
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	382 octave_idx_type nr = a.rows();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	383 octave_idx_type nc = a.cols();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	384 octave_idx_type b_nr = b.rows();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	385 octave_idx_type b_nc = b.cols();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	386 SparseMatrix x;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	387 volatile octave_idx_type ii, x_nz;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	388
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	389 if (nr < 0 \|\| nc < 0 \|\| nr != b_nr)
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	390 (*current_liboctave_error_handler)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	391 ("matrix dimension mismatch in solution of minimum norm problem");
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	392 else if (nr == 0 \|\| nc == 0 \|\| b_nc == 0)
be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	393 x = SparseMatrix (nc, b_nc);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	394 else if (nr >= nc)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	395 {
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	396 SparseQR q (a, 3);
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	397 if (! q.ok ())
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	398 return SparseMatrix();
10506 bdf5d85cfc5e replace nzmax by nnz where appropriate in liboctave Jaroslav Hajek <highegg@gmail.com> parents: 10350 diff changeset	399 x = SparseMatrix (nc, b_nc, b.nnz());
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	400 x.xcidx(0) = 0;
10506 bdf5d85cfc5e replace nzmax by nnz where appropriate in liboctave Jaroslav Hajek <highegg@gmail.com> parents: 10350 diff changeset	401 x_nz = b.nnz();
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	402 ii = 0;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	403 OCTAVE_LOCAL_BUFFER (double, Xx, (b_nr > nc ? b_nr : nc));
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	404 OCTAVE_LOCAL_BUFFER (double, buf, q.S()->m2);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	405 for (volatile octave_idx_type i = 0, idx = 0; i < b_nc; i++, idx+=nc)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	406 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	407 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	408 for (octave_idx_type j = 0; j < b_nr; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	409 Xx[j] = b.xelem(j,i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	410 for (octave_idx_type j = nr; j < q.S()->m2; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	411 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	412 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	413 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	414 CXSPARSE_DNAME (_ipvec) (q.S()->pinv, Xx, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	415 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	416 CXSPARSE_DNAME (_ipvec) (nr, q.S()->Pinv, Xx, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	417 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	418 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	419 for (volatile octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	420 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	421 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	422 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	423 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	424 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	425 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	426 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	427 CXSPARSE_DNAME (_usolve) (q.N()->U, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	428 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	429 CXSPARSE_DNAME (_ipvec) (q.S()->q, buf, Xx, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	430 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	431 CXSPARSE_DNAME (_ipvec) (nc, q.S()->Q, buf, Xx);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	432 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	433 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	434
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	435 for (octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	436 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	437 double tmp = Xx[j];
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	438 if (tmp != 0.0)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	439 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	440 if (ii == x_nz)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	441 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	442 // Resize the sparse matrix
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	443 octave_idx_type sz = x_nz * (b_nc - i) / b_nc;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	444 sz = (sz > 10 ? sz : 10) + x_nz;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	445 x.change_capacity (sz);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	446 x_nz = sz;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	447 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	448 x.xdata(ii) = tmp;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	449 x.xridx(ii++) = j;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	450 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	451 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	452 x.xcidx(i+1) = ii;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	453 }
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	454 info = 0;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	455 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	456 else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	457 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	458 SparseMatrix at = a.hermitian();
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	459 SparseQR q (at, 3);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	460 if (! q.ok ())
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	461 return SparseMatrix();
10506 bdf5d85cfc5e replace nzmax by nnz where appropriate in liboctave Jaroslav Hajek <highegg@gmail.com> parents: 10350 diff changeset	462 x = SparseMatrix (nc, b_nc, b.nnz());
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	463 x.xcidx(0) = 0;
10506 bdf5d85cfc5e replace nzmax by nnz where appropriate in liboctave Jaroslav Hajek <highegg@gmail.com> parents: 10350 diff changeset	464 x_nz = b.nnz();
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	465 ii = 0;
5681 233d98d95659 [project @ 2006-03-16 17:48:55 by dbateman] dbateman parents: 5648 diff changeset	466 volatile octave_idx_type nbuf = (nc > q.S()->m2 ? nc : q.S()->m2);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	467 OCTAVE_LOCAL_BUFFER (double, Xx, (b_nr > nc ? b_nr : nc));
5681 233d98d95659 [project @ 2006-03-16 17:48:55 by dbateman] dbateman parents: 5648 diff changeset	468 OCTAVE_LOCAL_BUFFER (double, buf, nbuf);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	469 for (volatile octave_idx_type i = 0, idx = 0; i < b_nc; i++, idx+=nc)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	470 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	471 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	472 for (octave_idx_type j = 0; j < b_nr; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	473 Xx[j] = b.xelem(j,i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	474 for (octave_idx_type j = nr; j < nbuf; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	475 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	476 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	477 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	478 CXSPARSE_DNAME (_pvec) (q.S()->q, Xx, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	479 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	480 CXSPARSE_DNAME (_pvec) (nr, q.S()->Q, Xx, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	481 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	482 CXSPARSE_DNAME (_utsolve) (q.N()->U, buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	483 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	484 for (volatile octave_idx_type j = nr-1; j >= 0; j--)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	485 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	486 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	487 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	488 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	489 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	490 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	491 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	492 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	493 CXSPARSE_DNAME (_pvec) (q.S()->pinv, buf, Xx, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	494 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	495 CXSPARSE_DNAME (_pvec) (nc, q.S()->Pinv, buf, Xx);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	496 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	497 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	498
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	499 for (octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	500 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	501 double tmp = Xx[j];
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	502 if (tmp != 0.0)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	503 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	504 if (ii == x_nz)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	505 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	506 // Resize the sparse matrix
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	507 octave_idx_type sz = x_nz * (b_nc - i) / b_nc;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	508 sz = (sz > 10 ? sz : 10) + x_nz;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	509 x.change_capacity (sz);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	510 x_nz = sz;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	511 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	512 x.xdata(ii) = tmp;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	513 x.xridx(ii++) = j;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	514 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	515 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	516 x.xcidx(i+1) = ii;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	517 }
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	518 info = 0;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	519 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	520
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	521 x.maybe_compress ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	522 return x;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	523 #else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	524 return SparseMatrix ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	525 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	526 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	527
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	528 ComplexMatrix
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	529 qrsolve(const SparseMatrix&a, const ComplexMatrix &b, octave_idx_type &info)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	530 {
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	531 info = -1;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	532 #ifdef HAVE_CXSPARSE
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	533 octave_idx_type nr = a.rows();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	534 octave_idx_type nc = a.cols();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	535 octave_idx_type b_nc = b.cols();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	536 octave_idx_type b_nr = b.rows();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	537 ComplexMatrix x;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	538
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	539 if (nr < 0 \|\| nc < 0 \|\| nr != b_nr)
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	540 (*current_liboctave_error_handler)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	541 ("matrix dimension mismatch in solution of minimum norm problem");
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	542 else if (nr == 0 \|\| nc == 0 \|\| b_nc == 0)
be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	543 x = ComplexMatrix (nc, b_nc, Complex (0.0, 0.0));
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	544 else if (nr >= nc)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	545 {
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	546 SparseQR q (a, 3);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	547 if (! q.ok ())
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	548 return ComplexMatrix();
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	549 x.resize(nc, b_nc);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	550 Complex *vec = x.fortran_vec();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	551 OCTAVE_LOCAL_BUFFER (double, Xx, (b_nr > nc ? b_nr : nc));
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	552 OCTAVE_LOCAL_BUFFER (double, Xz, (b_nr > nc ? b_nr : nc));
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	553 OCTAVE_LOCAL_BUFFER (double, buf, q.S()->m2);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	554 for (volatile octave_idx_type i = 0, idx = 0; i < b_nc; i++, idx+=nc)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	555 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	556 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	557 for (octave_idx_type j = 0; j < b_nr; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	558 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	559 Complex c = b.xelem (j,i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	560 Xx[j] = std::real (c);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	561 Xz[j] = std::imag (c);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	562 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	563 for (octave_idx_type j = nr; j < q.S()->m2; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	564 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	565 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	566 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	567 CXSPARSE_DNAME (_ipvec) (q.S()->pinv, Xx, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	568 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	569 CXSPARSE_DNAME (_ipvec) (nr, q.S()->Pinv, Xx, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	570 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	571 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	572 for (volatile octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	573 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	574 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	575 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	576 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	577 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	578 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	579 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	580 CXSPARSE_DNAME (_usolve) (q.N()->U, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	581 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	582 CXSPARSE_DNAME (_ipvec) (q.S()->q, buf, Xx, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	583 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	584 CXSPARSE_DNAME (_ipvec) (nc, q.S()->Q, buf, Xx);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	585 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	586 for (octave_idx_type j = nr; j < q.S()->m2; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	587 buf[j] = 0.;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	588 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	589 CXSPARSE_DNAME (_ipvec) (q.S()->pinv, Xz, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	590 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	591 CXSPARSE_DNAME (_ipvec) (nr, q.S()->Pinv, Xz, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	592 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	593 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	594 for (volatile octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	595 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	596 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	597 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	598 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	599 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	600 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	601 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	602 CXSPARSE_DNAME (_usolve) (q.N()->U, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	603 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	604 CXSPARSE_DNAME (_ipvec) (q.S()->q, buf, Xz, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	605 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	606 CXSPARSE_DNAME (_ipvec) (nc, q.S()->Q, buf, Xz);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	607 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	608 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	609 for (octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	610 vec[j+idx] = Complex (Xx[j], Xz[j]);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	611 }
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	612 info = 0;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	613 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	614 else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	615 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	616 SparseMatrix at = a.hermitian();
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	617 SparseQR q (at, 3);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	618 if (! q.ok ())
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	619 return ComplexMatrix();
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	620 x.resize(nc, b_nc);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	621 Complex *vec = x.fortran_vec();
5681 233d98d95659 [project @ 2006-03-16 17:48:55 by dbateman] dbateman parents: 5648 diff changeset	622 volatile octave_idx_type nbuf = (nc > q.S()->m2 ? nc : q.S()->m2);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	623 OCTAVE_LOCAL_BUFFER (double, Xx, (b_nr > nc ? b_nr : nc));
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	624 OCTAVE_LOCAL_BUFFER (double, Xz, (b_nr > nc ? b_nr : nc));
5681 233d98d95659 [project @ 2006-03-16 17:48:55 by dbateman] dbateman parents: 5648 diff changeset	625 OCTAVE_LOCAL_BUFFER (double, buf, nbuf);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	626 for (volatile octave_idx_type i = 0, idx = 0; i < b_nc; i++, idx+=nc)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	627 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	628 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	629 for (octave_idx_type j = 0; j < b_nr; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	630 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	631 Complex c = b.xelem (j,i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	632 Xx[j] = std::real (c);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	633 Xz[j] = std::imag (c);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	634 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	635 for (octave_idx_type j = nr; j < nbuf; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	636 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	637 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	638 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	639 CXSPARSE_DNAME (_pvec) (q.S()->q, Xx, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	640 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	641 CXSPARSE_DNAME (_pvec) (nr, q.S()->Q, Xx, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	642 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	643 CXSPARSE_DNAME (_utsolve) (q.N()->U, buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	644 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	645 for (volatile octave_idx_type j = nr-1; j >= 0; j--)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	646 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	647 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	648 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	649 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	650 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	651 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	652 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	653 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	654 CXSPARSE_DNAME (_pvec) (q.S()->pinv, buf, Xx, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	655 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	656 CXSPARSE_DNAME (_pvec) (nc, q.S()->Pinv, buf, Xx);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	657 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	658 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	659 for (octave_idx_type j = nr; j < nbuf; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	660 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	661 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	662 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	663 CXSPARSE_DNAME (_pvec) (q.S()->q, Xz, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	664 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	665 CXSPARSE_DNAME (_pvec) (nr, q.S()->Q, Xz, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	666 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	667 CXSPARSE_DNAME (_utsolve) (q.N()->U, buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	668 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	669 for (volatile octave_idx_type j = nr-1; j >= 0; j--)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	670 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	671 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	672 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	673 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	674 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	675 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	676 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	677 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	678 CXSPARSE_DNAME (_pvec) (q.S()->pinv, buf, Xz, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	679 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	680 CXSPARSE_DNAME (_pvec) (nc, q.S()->Pinv, buf, Xz);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	681 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	682 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	683 for (octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	684 vec[j+idx] = Complex (Xx[j], Xz[j]);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	685 }
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	686 info = 0;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	687 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	688
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	689 return x;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	690 #else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	691 return ComplexMatrix ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	692 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	693 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	694
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	695 SparseComplexMatrix
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	696 qrsolve(const SparseMatrix&a, const SparseComplexMatrix &b, octave_idx_type &info)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	697 {
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	698 info = -1;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	699 #ifdef HAVE_CXSPARSE
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	700 octave_idx_type nr = a.rows();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	701 octave_idx_type nc = a.cols();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	702 octave_idx_type b_nr = b.rows();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	703 octave_idx_type b_nc = b.cols();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	704 SparseComplexMatrix x;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	705 volatile octave_idx_type ii, x_nz;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	706
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	707 if (nr < 0 \|\| nc < 0 \|\| nr != b_nr)
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	708 (*current_liboctave_error_handler)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	709 ("matrix dimension mismatch in solution of minimum norm problem");
6924 be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	710 else if (nr == 0 \|\| nc == 0 \|\| b_nc == 0)
be176b7e110a [project @ 2007-09-25 19:08:19 by dbateman] dbateman parents: 5797 diff changeset	711 x = SparseComplexMatrix (nc, b_nc);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	712 else if (nr >= nc)
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	713 {
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	714 SparseQR q (a, 3);
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	715 if (! q.ok ())
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	716 return SparseComplexMatrix();
10506 bdf5d85cfc5e replace nzmax by nnz where appropriate in liboctave Jaroslav Hajek <highegg@gmail.com> parents: 10350 diff changeset	717 x = SparseComplexMatrix (nc, b_nc, b.nnz());
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	718 x.xcidx(0) = 0;
10506 bdf5d85cfc5e replace nzmax by nnz where appropriate in liboctave Jaroslav Hajek <highegg@gmail.com> parents: 10350 diff changeset	719 x_nz = b.nnz();
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	720 ii = 0;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	721 OCTAVE_LOCAL_BUFFER (double, Xx, (b_nr > nc ? b_nr : nc));
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	722 OCTAVE_LOCAL_BUFFER (double, Xz, (b_nr > nc ? b_nr : nc));
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	723 OCTAVE_LOCAL_BUFFER (double, buf, q.S()->m2);
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	724 for (volatile octave_idx_type i = 0, idx = 0; i < b_nc; i++, idx+=nc)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	725 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	726 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	727 for (octave_idx_type j = 0; j < b_nr; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	728 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	729 Complex c = b.xelem (j,i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	730 Xx[j] = std::real (c);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	731 Xz[j] = std::imag (c);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	732 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	733 for (octave_idx_type j = nr; j < q.S()->m2; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	734 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	735 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	736 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	737 CXSPARSE_DNAME (_ipvec) (q.S()->pinv, Xx, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	738 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	739 CXSPARSE_DNAME (_ipvec) (nr, q.S()->Pinv, Xx, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	740 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	741 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	742 for (volatile octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	743 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	744 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	745 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	746 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	747 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	748 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	749 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	750 CXSPARSE_DNAME (_usolve) (q.N()->U, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	751 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	752 CXSPARSE_DNAME (_ipvec) (q.S()->q, buf, Xx, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	753 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	754 CXSPARSE_DNAME (_ipvec) (nc, q.S()->Q, buf, Xx);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	755 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	756 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	757 for (octave_idx_type j = nr; j < q.S()->m2; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	758 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	759 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	760 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	761 CXSPARSE_DNAME (_ipvec) (q.S()->pinv, Xz, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	762 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	763 CXSPARSE_DNAME (_ipvec) (nr, q.S()->Pinv, Xz, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	764 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	765 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	766 for (volatile octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	767 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	768 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	769 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	770 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	771 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	772 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	773 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	774 CXSPARSE_DNAME (_usolve) (q.N()->U, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	775 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	776 CXSPARSE_DNAME (_ipvec) (q.S()->q, buf, Xz, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	777 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	778 CXSPARSE_DNAME (_ipvec) (nc, q.S()->Q, buf, Xz);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	779 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	780 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	781
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	782 for (octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	783 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	784 Complex tmp = Complex (Xx[j], Xz[j]);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	785 if (tmp != 0.0)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	786 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	787 if (ii == x_nz)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	788 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	789 // Resize the sparse matrix
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	790 octave_idx_type sz = x_nz * (b_nc - i) / b_nc;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	791 sz = (sz > 10 ? sz : 10) + x_nz;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	792 x.change_capacity (sz);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	793 x_nz = sz;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	794 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	795 x.xdata(ii) = tmp;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	796 x.xridx(ii++) = j;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	797 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	798 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	799 x.xcidx(i+1) = ii;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	800 }
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	801 info = 0;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	802 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	803 else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	804 {
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	805 SparseMatrix at = a.hermitian();
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	806 SparseQR q (at, 3);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	807 if (! q.ok ())
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	808 return SparseComplexMatrix();
10506 bdf5d85cfc5e replace nzmax by nnz where appropriate in liboctave Jaroslav Hajek <highegg@gmail.com> parents: 10350 diff changeset	809 x = SparseComplexMatrix (nc, b_nc, b.nnz());
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	810 x.xcidx(0) = 0;
10506 bdf5d85cfc5e replace nzmax by nnz where appropriate in liboctave Jaroslav Hajek <highegg@gmail.com> parents: 10350 diff changeset	811 x_nz = b.nnz();
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	812 ii = 0;
5681 233d98d95659 [project @ 2006-03-16 17:48:55 by dbateman] dbateman parents: 5648 diff changeset	813 volatile octave_idx_type nbuf = (nc > q.S()->m2 ? nc : q.S()->m2);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	814 OCTAVE_LOCAL_BUFFER (double, Xx, (b_nr > nc ? b_nr : nc));
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	815 OCTAVE_LOCAL_BUFFER (double, Xz, (b_nr > nc ? b_nr : nc));
5681 233d98d95659 [project @ 2006-03-16 17:48:55 by dbateman] dbateman parents: 5648 diff changeset	816 OCTAVE_LOCAL_BUFFER (double, buf, nbuf);
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	817 for (volatile octave_idx_type i = 0, idx = 0; i < b_nc; i++, idx+=nc)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	818 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	819 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	820 for (octave_idx_type j = 0; j < b_nr; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	821 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	822 Complex c = b.xelem (j,i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	823 Xx[j] = std::real (c);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	824 Xz[j] = std::imag (c);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	825 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	826 for (octave_idx_type j = nr; j < nbuf; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	827 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	828 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	829 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	830 CXSPARSE_DNAME (_pvec) (q.S()->q, Xx, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	831 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	832 CXSPARSE_DNAME (_pvec) (nr, q.S()->Q, Xx, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	833 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	834 CXSPARSE_DNAME (_utsolve) (q.N()->U, buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	835 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	836 for (volatile octave_idx_type j = nr-1; j >= 0; j--)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	837 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	838 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	839 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	840 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	841 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	842 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	843 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	844 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	845 CXSPARSE_DNAME (_pvec) (q.S()->pinv, buf, Xx, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	846 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	847 CXSPARSE_DNAME (_pvec) (nc, q.S()->Pinv, buf, Xx);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	848 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	849 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	850 for (octave_idx_type j = nr; j < nbuf; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	851 buf[j] = 0.;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	852 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	853 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	854 CXSPARSE_DNAME (_pvec) (q.S()->q, Xz, buf, nr);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	855 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	856 CXSPARSE_DNAME (_pvec) (nr, q.S()->Q, Xz, buf);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	857 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	858 CXSPARSE_DNAME (_utsolve) (q.N()->U, buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	859 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	860 for (volatile octave_idx_type j = nr-1; j >= 0; j--)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	861 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	862 octave_quit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	863 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	864 CXSPARSE_DNAME (_happly) (q.N()->L, j, q.N()->B[j], buf);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	865 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	866 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	867 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	868 #if defined(CS_VER) && (CS_VER >= 2)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	869 CXSPARSE_DNAME (_pvec) (q.S()->pinv, buf, Xz, nc);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	870 #else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	871 CXSPARSE_DNAME (_pvec) (nc, q.S()->Pinv, buf, Xz);
5792 eb90c83b4f91 [project @ 2006-05-04 20:14:49 by dbateman] dbateman parents: 5775 diff changeset	872 #endif
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	873 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	874
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	875 for (octave_idx_type j = 0; j < nc; j++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	876 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	877 Complex tmp = Complex (Xx[j], Xz[j]);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	878 if (tmp != 0.0)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	879 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	880 if (ii == x_nz)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	881 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	882 // Resize the sparse matrix
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	883 octave_idx_type sz = x_nz * (b_nc - i) / b_nc;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	884 sz = (sz > 10 ? sz : 10) + x_nz;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	885 x.change_capacity (sz);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	886 x_nz = sz;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	887 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	888 x.xdata(ii) = tmp;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	889 x.xridx(ii++) = j;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	890 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	891 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	892 x.xcidx(i+1) = ii;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	893 }
5797 11fcab4c461d [project @ 2006-05-09 06:15:17 by dbateman] dbateman parents: 5792 diff changeset	894 info = 0;
5610 9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	895 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	896
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	897 x.maybe_compress ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	898 return x;
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	899 #else
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	900 return SparseComplexMatrix ();
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	901 #endif
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	902 }
9761b7d24e9e [project @ 2006-02-09 09:12:02 by dbateman] dbateman parents: diff changeset	903
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	904 Matrix
12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	905 qrsolve(const SparseMatrix &a, const MArray<double> &b,
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	906 octave_idx_type &info)
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	907 {
12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	908 return qrsolve (a, Matrix (b), info);
7789 82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: 7505 diff changeset	909 }
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: 7505 diff changeset	910
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	911 ComplexMatrix
12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	912 qrsolve(const SparseMatrix &a, const MArray<Complex> &b,
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	913 octave_idx_type &info)
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11528 diff changeset	914 {
7789 82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: 7505 diff changeset	915 return qrsolve (a, ComplexMatrix (b), info);
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: 7505 diff changeset	916 }

Mercurial > octave-nkf

annotate liboctave/SparseQR.cc @ 12312:b10ea6efdc58 release-3-4-x ss-3-3-91