annotate src/svd.cc @ 1192:b6360f2d4fa6

[project @ 1995-03-30 21:38:35 by jwe]
author jwe
date Thu, 30 Mar 1995 21:38:35 +0000
parents dfe01093f657
children f62e82e0a2ce
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-svd.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
1009
dfe01093f657 [project @ 1995-01-04 04:05:12 by jwe]
jwe
parents: 718
diff changeset
4 Copyright (C) 1993, 1994, 1995 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
1192
b6360f2d4fa6 [project @ 1995-03-30 21:38:35 by jwe]
jwe
parents: 1009
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 "dbleSVD.h"
393e95f46b51 [project @ 1994-06-06 00:05:20 by jwe]
jwe
parents: 240
diff changeset
29 #include "CmplxSVD.h"
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
30
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
31 #include "tree-const.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
32 #include "user-prefs.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
33 #include "gripes.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
34 #include "error.h"
628
aecbe369233b [project @ 1994-08-19 13:27:51 by jwe]
jwe
parents: 622
diff changeset
35 #include "utils.h"
544
20fbad23ae51 [project @ 1994-07-22 05:04:44 by jwe]
jwe
parents: 539
diff changeset
36 #include "help.h"
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
37 #include "defun-dld.h"
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
38
701
0a81458ef677 [project @ 1994-09-15 02:23:24 by jwe]
jwe
parents: 636
diff changeset
39 DEFUN_DLD_BUILTIN ("svd", Fsvd, Ssvd, 2, 3,
539
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
40 "S = svd (X) or [U, S, V] = svd (X [, 0])\n\
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
41 \n\
539
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
42 Compute the singular value decomposition of X. Given a second input\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
43 argument, an `economy' sized factorization is computed that omits\n\
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
44 unnecessary rows and columns of U and V")
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
45 {
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
46 Octave_object retval;
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
47
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
48 int nargin = args.length ();
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
49
712
36ba0576bd1b [project @ 1994-09-19 14:18:15 by jwe]
jwe
parents: 701
diff changeset
50 if (nargin < 1 || nargin > 2 || nargout == 2 || nargout > 3)
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
51 {
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
52 print_usage ("svd");
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
53 return retval;
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
54 }
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
55
712
36ba0576bd1b [project @ 1994-09-19 14:18:15 by jwe]
jwe
parents: 701
diff changeset
56 tree_constant arg = args(0);
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
57
718
e81d3a66725e [project @ 1994-09-21 14:58:18 by jwe]
jwe
parents: 712
diff changeset
58 int arg_is_empty = empty_arg ("svd", arg.rows (), arg.columns ());
e81d3a66725e [project @ 1994-09-21 14:58:18 by jwe]
jwe
parents: 712
diff changeset
59
e81d3a66725e [project @ 1994-09-21 14:58:18 by jwe]
jwe
parents: 712
diff changeset
60 if (arg_is_empty < 0)
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
61 return retval;
718
e81d3a66725e [project @ 1994-09-21 14:58:18 by jwe]
jwe
parents: 712
diff changeset
62 else if (arg_is_empty > 0)
e81d3a66725e [project @ 1994-09-21 14:58:18 by jwe]
jwe
parents: 712
diff changeset
63 return Octave_object (3, Matrix ());
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
64
712
36ba0576bd1b [project @ 1994-09-19 14:18:15 by jwe]
jwe
parents: 701
diff changeset
65 SVD::type type = (nargin == 2) ? SVD::economy : SVD::std;
539
5ec10a984241 [project @ 1994-07-21 22:40:04 by jwe]
jwe
parents: 519
diff changeset
66
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
67 if (arg.is_real_type ())
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
68 {
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
69 Matrix tmp = arg.matrix_value ();
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: 628
diff changeset
71 if (! error_state)
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
72 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
73 SVD result (tmp, type);
622
6636e6198f01 [project @ 1994-08-16 15:08:03 by jwe]
jwe
parents: 620
diff changeset
74
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
75 DiagMatrix sigma = result.singular_values ();
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
76
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
77 if (nargout == 0 || nargout == 1)
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
78 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
79 retval(0) = tree_constant (sigma.diag (), 1);
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
80 }
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
81 else
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
82 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
83 retval(2) = result.right_singular_matrix ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
84 retval(1) = sigma;
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
85 retval(0) = result.left_singular_matrix ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
86 }
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
87 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
88 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
89 else if (arg.is_complex_type ())
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
90 {
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
91 ComplexMatrix ctmp = arg.complex_matrix_value ();
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
92
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
93 if (! error_state)
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
94 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
95 ComplexSVD result (ctmp, type);
622
6636e6198f01 [project @ 1994-08-16 15:08:03 by jwe]
jwe
parents: 620
diff changeset
96
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
97 DiagMatrix sigma = result.singular_values ();
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
98
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
99 if (nargout == 0 || nargout == 1)
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
100 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
101 retval(0) = tree_constant (sigma.diag (), 1);
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
102 }
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
103 else
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
104 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
105 retval(2) = result.right_singular_matrix ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
106 retval(1) = sigma;
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
107 retval(0) = result.left_singular_matrix ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 628
diff changeset
108 }
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
109 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
110 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
111 else
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
112 {
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
113 gripe_wrong_type_arg ("svd", arg);
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
114 return retval;
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
115 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
116
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
117 return retval;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
118 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
119
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
120 /*
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
121 ;;; Local Variables: ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
122 ;;; mode: C++ ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
123 ;;; page-delimiter: "^/\\*" ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
124 ;;; End: ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
125 */