annotate src/qr.cc @ 704:250fc1c93fe2

[project @ 1994-09-15 02:44:42 by jwe]
author jwe
date Thu, 15 Sep 1994 02:47:09 +0000
parents 0a81458ef677
children 36ba0576bd1b
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
7
d68036bcad4c [project @ 1993-08-08 03:53:43 by jwe]
jwe
parents: 1
diff changeset
1 // f-qr.cc -*- C++ -*-
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
2 /*
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
3
453
393e95f46b51 [project @ 1994-06-06 00:05:20 by jwe]
jwe
parents: 240
diff changeset
4 Copyright (C) 1993, 1994 John W. Eaton
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
5
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
6 This file is part of Octave.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
7
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
8 Octave is free software; you can redistribute it and/or modify it
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
9 under the terms of the GNU General Public License as published by the
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
10 Free Software Foundation; either version 2, or (at your option) any
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
11 later version.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
12
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
13 Octave is distributed in the hope that it will be useful, but WITHOUT
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
16 for more details.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
17
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
18 You should have received a copy of the GNU General Public License
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
19 along with Octave; see the file COPYING. If not, write to the Free
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
20 Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
21
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
22 */
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
23
240
a99f28f5e351 [project @ 1993-11-30 20:24:36 by jwe]
jwe
parents: 162
diff changeset
24 #ifdef HAVE_CONFIG_H
a99f28f5e351 [project @ 1993-11-30 20:24:36 by jwe]
jwe
parents: 162
diff changeset
25 #include "config.h"
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
26 #endif
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
27
453
393e95f46b51 [project @ 1994-06-06 00:05:20 by jwe]
jwe
parents: 240
diff changeset
28 #include "dbleQR.h"
393e95f46b51 [project @ 1994-06-06 00:05:20 by jwe]
jwe
parents: 240
diff changeset
29 #include "CmplxQR.h"
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
30
539
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
31 #include "dbleQRP.h"
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
32 #include "CmplxQRP.h"
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
33
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
34 #include "tree-const.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
35 #include "user-prefs.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
36 #include "gripes.h"
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
37 #include "utils.h"
544
20fbad23ae51 [project @ 1994-07-22 05:04:44 by jwe]
jwe
parents: 539
diff changeset
38 #include "help.h"
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
39 #include "defun-dld.h"
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
40
701
0a81458ef677 [project @ 1994-09-15 02:23:24 by jwe]
jwe
parents: 636
diff changeset
41 DEFUN_DLD_BUILTIN ("qr", Fqr, Sqr, 2, 2,
539
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
42 "[Q, R] = qr (X): form Q unitary and R upper triangular such\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
43 that Q * R = X\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
44 \n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
45 [Q, R] = qr (X, 0): form the economy decomposition such that if X is\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
46 if X is m by n then only the first n columns of Q\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
47 are computed.\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
48 \n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
49 [Q, R, P] = qr (X): form QRP factorization of X where\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
50 P is a permutation matrix such that\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
51 A * P = Q * R\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
52 \n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
53 [Q, R, P] = qr (X, 0): form the economy decomposition with \n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
54 permutation vector P such that Q * R = X (:, P)\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
55 \n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
56 qr (X) alone returns the output of the LAPACK routine dgeqrf, such\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
57 that R = triu (qr (X))")
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
58 {
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
59 Octave_object retval;
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
60
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
61 int nargin = args.length ();
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
62
539
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
63 if (nargin != 2 && nargin != 3 || nargout > 3)
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
64 {
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
65 print_usage ("qr");
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
66 return retval;
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
67 }
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
68
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
69 tree_constant arg = args(1);
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
70
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
71 if (empty_arg ("qr", arg.rows (), arg.columns ()) < 0)
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
72 return retval;
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
73
539
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
74 QR::type type = nargout == 1 ? QR::raw
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
75 : (nargin == 3 ? QR::economy : QR::std);
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
76
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
77 if (arg.is_real_type ())
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
78 {
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
79 Matrix m = arg.matrix_value ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
80
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
81 if (! error_state)
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
82 {
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
83 if (nargout < 3)
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
84 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
85 QR fact (m, type);
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
86 retval(1) = fact.R ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
87 retval(0) = fact.Q ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
88 }
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
89 else
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
90 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
91 QRP fact (m, type);
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
92 retval(2) = fact.P ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
93 retval(1) = fact.R ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
94 retval(0) = fact.Q ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
95 }
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
96 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
97 }
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
98 else if (arg.is_complex_type ())
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
99 {
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
100 ComplexMatrix m = arg.complex_matrix_value ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
101
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
102 if (! error_state)
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
103 {
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
104 if (nargout < 3)
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
105 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
106 ComplexQR fact (m, type);
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
107 retval(1) = fact.R ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
108 retval(0) = fact.Q ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
109 }
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
110 else
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
111 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
112 ComplexQRP fact (m, type);
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
113 retval(2) = fact.P ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
114 retval(1) = fact.R ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
115 retval(0) = fact.Q ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
116 }
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
117 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
118 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
119 else
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
120 {
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
121 gripe_wrong_type_arg ("qr", arg);
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
122 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
123
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
124 return retval;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
125 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
126
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
127 /*
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
128 ;;; Local Variables: ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
129 ;;; mode: C++ ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
130 ;;; page-delimiter: "^/\\*" ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
131 ;;; End: ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
132 */