Mercurial > octave
annotate src/DLD-FUNCTIONS/schur.cc @ 10607:f7501986e42d
make schur more Matlab compatible
author | Jaroslav Hajek <highegg@gmail.com> |
---|---|
date | Thu, 06 May 2010 09:49:36 +0200 |
parents | d0ce5e973937 |
children | 9c9e07f5eb1c |
rev | line source |
---|---|
2928 | 1 /* |
2 | |
8920 | 3 Copyright (C) 1996, 1997, 1999, 2000, 2004, 2005, 2006, 2007, 2008, 2009 |
7017 | 4 John W. Eaton |
2928 | 5 |
6 This file is part of Octave. | |
7 | |
8 Octave is free software; you can redistribute it and/or modify it | |
9 under the terms of the GNU General Public License as published by the | |
7016 | 10 Free Software Foundation; either version 3 of the License, or (at your |
11 option) any later version. | |
2928 | 12 |
13 Octave is distributed in the hope that it will be useful, but WITHOUT | |
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |
16 for more details. | |
17 | |
18 You should have received a copy of the GNU General Public License | |
7016 | 19 along with Octave; see the file COPYING. If not, see |
20 <http://www.gnu.org/licenses/>. | |
2928 | 21 |
22 */ | |
23 | |
24 #ifdef HAVE_CONFIG_H | |
25 #include <config.h> | |
26 #endif | |
27 | |
28 #include <string> | |
29 | |
30 #include "CmplxSCHUR.h" | |
31 #include "dbleSCHUR.h" | |
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
32 #include "fCmplxSCHUR.h" |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
33 #include "floatSCHUR.h" |
2928 | 34 |
35 #include "defun-dld.h" | |
36 #include "error.h" | |
37 #include "gripes.h" | |
38 #include "oct-obj.h" | |
39 #include "utils.h" | |
40 | |
41 DEFUN_DLD (schur, args, nargout, | |
3548 | 42 "-*- texinfo -*-\n\ |
3372 | 43 @deftypefn {Loadable Function} {@var{s} =} schur (@var{a})\n\ |
10607
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
44 @deftypefnx {Loadable Function} {@var{s} =} schur (@var{a}, \"complex\")\n\ |
3372 | 45 @deftypefnx {Loadable Function} {[@var{u}, @var{s}] =} schur (@var{a}, @var{opt})\n\ |
46 @cindex Schur decomposition\n\ | |
47 The Schur decomposition is used to compute eigenvalues of a\n\ | |
48 square matrix, and has applications in the solution of algebraic\n\ | |
49 Riccati equations in control (see @code{are} and @code{dare}).\n\ | |
50 @code{schur} always returns\n\ | |
51 @tex\n\ | |
52 $S = U^T A U$\n\ | |
53 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
54 @ifnottex\n\ |
3372 | 55 @code{s = u' * a * u}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
56 @end ifnottex\n\ |
3372 | 57 where\n\ |
58 @tex\n\ | |
59 $U$\n\ | |
60 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
61 @ifnottex\n\ |
3372 | 62 @code{u}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
63 @end ifnottex\n\ |
3372 | 64 is a unitary matrix\n\ |
65 @tex\n\ | |
66 ($U^T U$ is identity)\n\ | |
67 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
68 @ifnottex\n\ |
3372 | 69 (@code{u'* u} is identity)\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
70 @end ifnottex\n\ |
3372 | 71 and\n\ |
72 @tex\n\ | |
73 $S$\n\ | |
74 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
75 @ifnottex\n\ |
3372 | 76 @code{s}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
77 @end ifnottex\n\ |
3372 | 78 is upper triangular. The eigenvalues of\n\ |
79 @tex\n\ | |
80 $A$ (and $S$)\n\ | |
81 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
82 @ifnottex\n\ |
3372 | 83 @code{a} (and @code{s})\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
84 @end ifnottex\n\ |
3372 | 85 are the diagonal elements of\n\ |
86 @tex\n\ | |
5555 | 87 $S$.\n\ |
3372 | 88 @end tex\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
89 @ifnottex\n\ |
5555 | 90 @code{s}.\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
91 @end ifnottex\n\ |
3372 | 92 If the matrix\n\ |
93 @tex\n\ | |
94 $A$\n\ | |
95 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
96 @ifnottex\n\ |
3372 | 97 @code{a}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
98 @end ifnottex\n\ |
3372 | 99 is real, then the real Schur decomposition is computed, in which the\n\ |
100 matrix\n\ | |
101 @tex\n\ | |
102 $U$\n\ | |
103 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
104 @ifnottex\n\ |
3372 | 105 @code{u}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
106 @end ifnottex\n\ |
3372 | 107 is orthogonal and\n\ |
108 @tex\n\ | |
109 $S$\n\ | |
110 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
111 @ifnottex\n\ |
3372 | 112 @code{s}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
113 @end ifnottex\n\ |
3372 | 114 is block upper triangular\n\ |
115 with blocks of size at most\n\ | |
116 @tex\n\ | |
117 $2\\times 2$\n\ | |
118 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
119 @ifnottex\n\ |
3372 | 120 @code{2 x 2}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
121 @end ifnottex\n\ |
3600 | 122 along the diagonal. The diagonal elements of\n\ |
3372 | 123 @tex\n\ |
124 $S$\n\ | |
125 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
126 @ifnottex\n\ |
3372 | 127 @code{s}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
128 @end ifnottex\n\ |
3372 | 129 (or the eigenvalues of the\n\ |
130 @tex\n\ | |
131 $2\\times 2$\n\ | |
132 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
133 @ifnottex\n\ |
3372 | 134 @code{2 x 2}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
135 @end ifnottex\n\ |
3372 | 136 blocks, when\n\ |
137 appropriate) are the eigenvalues of\n\ | |
138 @tex\n\ | |
139 $A$\n\ | |
140 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
141 @ifnottex\n\ |
3372 | 142 @code{a}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
143 @end ifnottex\n\ |
3372 | 144 and\n\ |
145 @tex\n\ | |
146 $S$.\n\ | |
147 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
148 @ifnottex\n\ |
3372 | 149 @code{s}.\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
150 @end ifnottex\n\ |
2928 | 151 \n\ |
10607
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
152 A complex decomposition may be forced by passing \"complex\" as @var{opt}.\n\ |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
153 \n\ |
3372 | 154 The eigenvalues are optionally ordered along the diagonal according to\n\ |
155 the value of @code{opt}. @code{opt = \"a\"} indicates that all\n\ | |
156 eigenvalues with negative real parts should be moved to the leading\n\ | |
157 block of\n\ | |
158 @tex\n\ | |
159 $S$\n\ | |
160 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
161 @ifnottex\n\ |
3372 | 162 @code{s}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
163 @end ifnottex\n\ |
3372 | 164 (used in @code{are}), @code{opt = \"d\"} indicates that all eigenvalues\n\ |
165 with magnitude less than one should be moved to the leading block of\n\ | |
166 @tex\n\ | |
167 $S$\n\ | |
168 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
169 @ifnottex\n\ |
3372 | 170 @code{s}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
171 @end ifnottex\n\ |
3372 | 172 (used in @code{dare}), and @code{opt = \"u\"}, the default, indicates that\n\ |
173 no ordering of eigenvalues should occur. The leading\n\ | |
174 @tex\n\ | |
175 $k$\n\ | |
176 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
177 @ifnottex\n\ |
3372 | 178 @code{k}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
179 @end ifnottex\n\ |
3372 | 180 columns of\n\ |
181 @tex\n\ | |
182 $U$\n\ | |
183 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
184 @ifnottex\n\ |
3372 | 185 @code{u}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
186 @end ifnottex\n\ |
3372 | 187 always span the\n\ |
188 @tex\n\ | |
189 $A$-invariant\n\ | |
190 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
191 @ifnottex\n\ |
3372 | 192 @code{a}-invariant\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
193 @end ifnottex\n\ |
3372 | 194 subspace corresponding to the\n\ |
195 @tex\n\ | |
196 $k$\n\ | |
197 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
198 @ifnottex\n\ |
3372 | 199 @code{k}\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
200 @end ifnottex\n\ |
3372 | 201 leading eigenvalues of\n\ |
202 @tex\n\ | |
203 $S$.\n\ | |
204 @end tex\n\ | |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
205 @ifnottex\n\ |
3372 | 206 @code{s}.\n\ |
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7814
diff
changeset
|
207 @end ifnottex\n\ |
3372 | 208 @end deftypefn") |
2928 | 209 { |
210 octave_value_list retval; | |
211 | |
212 int nargin = args.length (); | |
213 | |
214 if (nargin < 1 || nargin > 2 || nargout > 2) | |
215 { | |
5823 | 216 print_usage (); |
2928 | 217 return retval; |
218 } | |
219 | |
220 octave_value arg = args(0); | |
221 | |
3523 | 222 std::string ord; |
2928 | 223 |
224 if (nargin == 2) | |
225 { | |
226 ord = args(1).string_value (); | |
227 | |
228 if (error_state) | |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
229 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
230 error ("schur: expecting string as second argument"); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
231 return retval; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
232 } |
2928 | 233 } |
234 | |
10607
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
235 bool force_complex = false; |
2928 | 236 |
10607
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
237 if (ord == "complex") |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
238 { |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
239 force_complex = true; |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
240 ord = std::string (); |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
241 } |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
242 else |
2928 | 243 { |
10607
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
244 char ord_char = ord.empty () ? 'U' : ord[0]; |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
245 |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
246 if (ord_char != 'U' && ord_char != 'A' && ord_char != 'D' |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
247 && ord_char != 'u' && ord_char != 'a' && ord_char != 'd') |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
248 { |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
249 warning ("schur: incorrect ordered schur argument `%c'", |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
250 ord.c_str ()); |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
251 return retval; |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
252 } |
2928 | 253 } |
254 | |
5275 | 255 octave_idx_type nr = arg.rows (); |
256 octave_idx_type nc = arg.columns (); | |
2928 | 257 |
258 if (nr != nc) | |
259 { | |
260 gripe_square_matrix_required ("schur"); | |
261 return retval; | |
262 } | |
263 | |
10607
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
264 if (! arg.is_numeric_type ()) |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
265 gripe_wrong_type_arg ("schur", arg); |
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
266 else if (arg.is_single_type ()) |
2928 | 267 { |
10607
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
268 if (! force_complex && arg.is_real_type ()) |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
269 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
270 FloatMatrix tmp = arg.float_matrix_value (); |
2928 | 271 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
272 if (! error_state) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
273 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
274 if (nargout == 0 || nargout == 1) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
275 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
276 FloatSCHUR result (tmp, ord, false); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
277 retval(0) = result.schur_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
278 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
279 else |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
280 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
281 FloatSCHUR result (tmp, ord, true); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
282 retval(1) = result.schur_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
283 retval(0) = result.unitary_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
284 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
285 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
286 } |
10607
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
287 else |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
288 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
289 FloatComplexMatrix ctmp = arg.float_complex_matrix_value (); |
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
290 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
291 if (! error_state) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
292 { |
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
293 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
294 if (nargout == 0 || nargout == 1) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
295 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
296 FloatComplexSCHUR result (ctmp, ord, false); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
297 retval(0) = result.schur_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
298 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
299 else |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
300 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
301 FloatComplexSCHUR result (ctmp, ord, true); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
302 retval(1) = result.schur_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
303 retval(0) = result.unitary_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
304 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
305 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
306 } |
2928 | 307 } |
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
308 else |
2928 | 309 { |
10607
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
310 if (! force_complex && arg.is_real_type ()) |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
311 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
312 Matrix tmp = arg.matrix_value (); |
2928 | 313 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
314 if (! error_state) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
315 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
316 if (nargout == 0 || nargout == 1) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
317 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
318 SCHUR result (tmp, ord, false); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
319 retval(0) = result.schur_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
320 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
321 else |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
322 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
323 SCHUR result (tmp, ord, true); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
324 retval(1) = result.schur_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
325 retval(0) = result.unitary_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
326 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
327 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
328 } |
10607
f7501986e42d
make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents:
10155
diff
changeset
|
329 else |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
330 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
331 ComplexMatrix ctmp = arg.complex_matrix_value (); |
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
332 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
333 if (! error_state) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
334 { |
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
335 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
336 if (nargout == 0 || nargout == 1) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
337 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
338 ComplexSCHUR result (ctmp, ord, false); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
339 retval(0) = result.schur_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
340 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
341 else |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
342 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
343 ComplexSCHUR result (ctmp, ord, true); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
344 retval(1) = result.schur_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
345 retval(0) = result.unitary_matrix (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
346 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
347 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9758
diff
changeset
|
348 } |
2928 | 349 } |
350 | |
351 return retval; | |
352 } | |
353 | |
354 /* | |
7814
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
355 |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
356 %!test |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
357 %! a = [1, 2, 3; 4, 5, 9; 7, 8, 6]; |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
358 %! [u, s] = schur (a); |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
359 %! assert(u' * a * u, s, sqrt (eps)); |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
360 |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
361 %!test |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
362 %! a = single([1, 2, 3; 4, 5, 9; 7, 8, 6]); |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
363 %! [u, s] = schur (a); |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
364 %! assert(u' * a * u, s, sqrt (eps('single'))); |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
365 |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
366 %!test |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
367 %! fail("schur ([1, 2; 3, 4], 2)","warning"); |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
368 |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
369 %!error <Invalid call to schur.*> schur (); |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
370 %!error schur ([1, 2, 3; 4, 5, 6]); |
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
371 |
2928 | 372 */ |