annotate liboctave/numeric/schur.cc @ 21663:79414f8e6f24

schur.cc: Declare selector variable as volatile to suppress gcc warning. schur.cc (init): Declare selector variable as volatile to suppress gcc warning.
author Rik <rik@octave.org>
date Mon, 02 May 2016 14:27:52 -0700
parents 40de9f8f23a6
children aba2e6293dd8
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
1 /*
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
2
19697
4197fc428c7d maint: Update copyright notices for 2015.
John W. Eaton <jwe@octave.org>
parents: 17769
diff changeset
3 Copyright (C) 1994-2015 John W. Eaton
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
4
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
5 This file is part of Octave.
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
6
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
7 Octave is free software; you can redistribute it and/or modify it
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
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: 5307
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: 5307
diff changeset
10 option) any later version.
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
11
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
12 Octave is distributed in the hope that it will be useful, but WITHOUT
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
15 for more details.
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
16
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
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: 5307
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: 5307
diff changeset
19 <http://www.gnu.org/licenses/>.
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
20
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
21 */
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
22
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
23 #ifdef HAVE_CONFIG_H
21301
40de9f8f23a6 Use '#include "config.h"' rather than <config.h>.
Rik <rik@octave.org>
parents: 21274
diff changeset
24 # include "config.h"
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
25 #endif
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
26
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
27 #include "CMatrix.h"
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
28 #include "dMatrix.h"
1847
2ffe49eb95a5 [project @ 1996-02-03 12:47:55 by jwe]
jwe
parents: 1756
diff changeset
29 #include "f77-fcn.h"
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
30 #include "fCMatrix.h"
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
31 #include "fMatrix.h"
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
32 #include "lo-error.h"
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
33 #include "schur.h"
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
34
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
35 typedef octave_idx_type (*double_selector) (const double&, const double&);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
36 typedef octave_idx_type (*float_selector) (const float&, const float&);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
37 typedef octave_idx_type (*complex_selector) (const Complex&);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
38 typedef octave_idx_type (*float_complex_selector) (const FloatComplex&);
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
39
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
40 extern "C"
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
41 {
4552
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe]
jwe
parents: 4192
diff changeset
42 F77_RET_T
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe]
jwe
parents: 4192
diff changeset
43 F77_FUNC (dgeesx, DGEESX) (F77_CONST_CHAR_ARG_DECL,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
44 F77_CONST_CHAR_ARG_DECL,
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
45 double_selector,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
46 F77_CONST_CHAR_ARG_DECL,
11495
8a5e980da6aa style fixes
John W. Eaton <jwe@octave.org>
parents: 10822
diff changeset
47 const octave_idx_type&, double*,
8a5e980da6aa style fixes
John W. Eaton <jwe@octave.org>
parents: 10822
diff changeset
48 const octave_idx_type&, octave_idx_type&,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
49 double*, double*, double*, const octave_idx_type&,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
50 double&, double&, double*, const octave_idx_type&,
11495
8a5e980da6aa style fixes
John W. Eaton <jwe@octave.org>
parents: 10822
diff changeset
51 octave_idx_type*, const octave_idx_type&,
8a5e980da6aa style fixes
John W. Eaton <jwe@octave.org>
parents: 10822
diff changeset
52 octave_idx_type*, octave_idx_type&
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
53 F77_CHAR_ARG_LEN_DECL
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
54 F77_CHAR_ARG_LEN_DECL
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
55 F77_CHAR_ARG_LEN_DECL);
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
56
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
57 F77_RET_T
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
58 F77_FUNC (sgeesx, SGEESX) (F77_CONST_CHAR_ARG_DECL,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
59 F77_CONST_CHAR_ARG_DECL,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
60 float_selector,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
61 F77_CONST_CHAR_ARG_DECL,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
62 const octave_idx_type&, float*,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
63 const octave_idx_type&, octave_idx_type&,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
64 float*, float*, float*, const octave_idx_type&,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
65 float&, float&, float*, const octave_idx_type&,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
66 octave_idx_type*, const octave_idx_type&,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
67 octave_idx_type*, octave_idx_type&
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
68 F77_CHAR_ARG_LEN_DECL
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
69 F77_CHAR_ARG_LEN_DECL
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
70 F77_CHAR_ARG_LEN_DECL);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
71
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
72 F77_RET_T
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
73 F77_FUNC (zgeesx, ZGEESX) (F77_CONST_CHAR_ARG_DECL,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
74 F77_CONST_CHAR_ARG_DECL,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
75 complex_selector,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
76 F77_CONST_CHAR_ARG_DECL,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
77 const octave_idx_type&, Complex*,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
78 const octave_idx_type&, octave_idx_type&,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
79 Complex*, Complex*, const octave_idx_type&,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
80 double&, double&, Complex*,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
81 const octave_idx_type&, double*,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
82 octave_idx_type*, octave_idx_type&
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
83 F77_CHAR_ARG_LEN_DECL
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
84 F77_CHAR_ARG_LEN_DECL
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
85 F77_CHAR_ARG_LEN_DECL);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
86
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
87 F77_RET_T
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
88 F77_FUNC (zrsf2csf, ZRSF2CSF) (const octave_idx_type&, Complex *,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
89 Complex *, double *, double *);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
90 F77_RET_T
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
91 F77_FUNC (cgeesx, CGEESX) (F77_CONST_CHAR_ARG_DECL,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
92 F77_CONST_CHAR_ARG_DECL,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
93 float_complex_selector,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
94 F77_CONST_CHAR_ARG_DECL,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
95 const octave_idx_type&, FloatComplex*,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
96 const octave_idx_type&, octave_idx_type&,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
97 FloatComplex*, FloatComplex*,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
98 const octave_idx_type&, float&, float&,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
99 FloatComplex*, const octave_idx_type&,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
100 float*, octave_idx_type*, octave_idx_type&
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
101 F77_CHAR_ARG_LEN_DECL
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
102 F77_CHAR_ARG_LEN_DECL
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
103 F77_CHAR_ARG_LEN_DECL);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
104
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
105 F77_RET_T
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
106 F77_FUNC (crsf2csf, CRSF2CSF) (const octave_idx_type&, FloatComplex *,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
107 FloatComplex *, float *, float *);
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
108 }
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
109
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
110 // For real types.
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
111
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
112 template <typename T>
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
113 static octave_idx_type
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
114 select_ana (const T& a, const T&)
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
115 {
17769
49a5a4be04a1 maint: Use GNU style coding conventions for code in liboctave/
Rik <rik@octave.org>
parents: 17744
diff changeset
116 return (a < 0.0);
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
117 }
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
118
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
119 template <typename T>
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
120 static octave_idx_type
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
121 select_dig (const T& a, const T& b)
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
122 {
1251
97eac19837dc [project @ 1995-04-11 15:58:32 by jwe]
jwe
parents: 1192
diff changeset
123 return (hypot (a, b) < 1.0);
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
124 }
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
125
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
126 // For complex types.
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
127
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
128 template <typename T>
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
129 static octave_idx_type
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
130 select_ana (const T& a)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
131 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
132 return a.real () < 0.0;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
133 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
134
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
135 template <typename T>
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
136 static octave_idx_type
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
137 select_dig (const T& a)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
138 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
139 return (abs (a) < 1.0);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
140 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
141
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
142 template <>
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
143 octave_idx_type
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
144 schur<Matrix>::init (const Matrix& a, const std::string& ord, bool calc_unitary)
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
145 {
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
146 octave_idx_type a_nr = a.rows ();
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
147 octave_idx_type a_nc = a.cols ();
1929
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
148
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
149 if (a_nr != a_nc)
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
150 (*current_liboctave_error_handler) ("schur: requires square matrix");
21136
7cac4e7458f2 maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents: 19697
diff changeset
151
7cac4e7458f2 maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents: 19697
diff changeset
152 if (a_nr == 0)
10607
f7501986e42d make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents: 10350
diff changeset
153 {
f7501986e42d make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents: 10350
diff changeset
154 schur_mat.clear ();
f7501986e42d make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents: 10350
diff changeset
155 unitary_mat.clear ();
f7501986e42d make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents: 10350
diff changeset
156 return 0;
f7501986e42d make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents: 10350
diff changeset
157 }
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
158
3334
5187390bfde6 [project @ 1999-11-03 20:41:16 by jwe]
jwe
parents: 2847
diff changeset
159 // Workspace requirements may need to be fixed if any of the
5187390bfde6 [project @ 1999-11-03 20:41:16 by jwe]
jwe
parents: 2847
diff changeset
160 // following change.
5187390bfde6 [project @ 1999-11-03 20:41:16 by jwe]
jwe
parents: 2847
diff changeset
161
5008
c2bb27ada496 [project @ 2004-09-17 14:45:39 by jwe]
jwe
parents: 4552
diff changeset
162 char jobvs;
1930
d20ab06301e8 [project @ 1996-02-11 22:30:18 by jwe]
jwe
parents: 1929
diff changeset
163 char sense = 'N';
d20ab06301e8 [project @ 1996-02-11 22:30:18 by jwe]
jwe
parents: 1929
diff changeset
164 char sort = 'N';
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
165
5008
c2bb27ada496 [project @ 2004-09-17 14:45:39 by jwe]
jwe
parents: 4552
diff changeset
166 if (calc_unitary)
c2bb27ada496 [project @ 2004-09-17 14:45:39 by jwe]
jwe
parents: 4552
diff changeset
167 jobvs = 'V';
c2bb27ada496 [project @ 2004-09-17 14:45:39 by jwe]
jwe
parents: 4552
diff changeset
168 else
c2bb27ada496 [project @ 2004-09-17 14:45:39 by jwe]
jwe
parents: 4552
diff changeset
169 jobvs = 'N';
c2bb27ada496 [project @ 2004-09-17 14:45:39 by jwe]
jwe
parents: 4552
diff changeset
170
1756
1af643fa00e3 [project @ 1996-01-22 04:55:17 by jwe]
jwe
parents: 1631
diff changeset
171 char ord_char = ord.empty () ? 'U' : ord[0];
1af643fa00e3 [project @ 1996-01-22 04:55:17 by jwe]
jwe
parents: 1631
diff changeset
172
1af643fa00e3 [project @ 1996-01-22 04:55:17 by jwe]
jwe
parents: 1631
diff changeset
173 if (ord_char == 'A' || ord_char == 'D' || ord_char == 'a' || ord_char == 'd')
1930
d20ab06301e8 [project @ 1996-02-11 22:30:18 by jwe]
jwe
parents: 1929
diff changeset
174 sort = 'S';
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
175
21663
79414f8e6f24 schur.cc: Declare selector variable as volatile to suppress gcc warning.
Rik <rik@octave.org>
parents: 21301
diff changeset
176 volatile double_selector selector = 0;
1929
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
177 if (ord_char == 'A' || ord_char == 'a')
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
178 selector = select_ana<double>;
1929
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
179 else if (ord_char == 'D' || ord_char == 'd')
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
180 selector = select_dig<double>;
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
181
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
182 octave_idx_type n = a_nc;
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
183 octave_idx_type lwork = 8 * n;
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
184 octave_idx_type liwork = 1;
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
185 octave_idx_type info;
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
186 octave_idx_type sdim;
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
187 double rconde;
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
188 double rcondv;
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
189
1929
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
190 schur_mat = a;
5008
c2bb27ada496 [project @ 2004-09-17 14:45:39 by jwe]
jwe
parents: 4552
diff changeset
191
c2bb27ada496 [project @ 2004-09-17 14:45:39 by jwe]
jwe
parents: 4552
diff changeset
192 if (calc_unitary)
10607
f7501986e42d make schur more Matlab compatible
Jaroslav Hajek <highegg@gmail.com>
parents: 10350
diff changeset
193 unitary_mat.clear (n, n);
1929
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
194
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
195 double *s = schur_mat.fortran_vec ();
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
196 double *q = unitary_mat.fortran_vec ();
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
197
11570
57632dea2446 attempt better backward compatibility for Array constructors
John W. Eaton <jwe@octave.org>
parents: 11523
diff changeset
198 Array<double> wr (dim_vector (n, 1));
1929
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
199 double *pwr = wr.fortran_vec ();
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
200
11570
57632dea2446 attempt better backward compatibility for Array constructors
John W. Eaton <jwe@octave.org>
parents: 11523
diff changeset
201 Array<double> wi (dim_vector (n, 1));
1929
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
202 double *pwi = wi.fortran_vec ();
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
203
11570
57632dea2446 attempt better backward compatibility for Array constructors
John W. Eaton <jwe@octave.org>
parents: 11523
diff changeset
204 Array<double> work (dim_vector (lwork, 1));
1929
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
205 double *pwork = work.fortran_vec ();
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
206
3334
5187390bfde6 [project @ 1999-11-03 20:41:16 by jwe]
jwe
parents: 2847
diff changeset
207 // BWORK is not referenced for the non-ordered Schur routine.
11570
57632dea2446 attempt better backward compatibility for Array constructors
John W. Eaton <jwe@octave.org>
parents: 11523
diff changeset
208 octave_idx_type ntmp = (ord_char == 'N' || ord_char == 'n') ? 0 : n;
57632dea2446 attempt better backward compatibility for Array constructors
John W. Eaton <jwe@octave.org>
parents: 11523
diff changeset
209 Array<octave_idx_type> bwork (dim_vector (ntmp, 1));
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
210 octave_idx_type *pbwork = bwork.fortran_vec ();
1929
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
211
11570
57632dea2446 attempt better backward compatibility for Array constructors
John W. Eaton <jwe@octave.org>
parents: 11523
diff changeset
212 Array<octave_idx_type> iwork (dim_vector (liwork, 1));
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5105
diff changeset
213 octave_idx_type *piwork = iwork.fortran_vec ();
1929
908f5b6676d7 [project @ 1996-02-11 22:05:08 by jwe]
jwe
parents: 1882
diff changeset
214
4552
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe]
jwe
parents: 4192
diff changeset
215 F77_XFCN (dgeesx, DGEESX, (F77_CONST_CHAR_ARG2 (&jobvs, 1),
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
216 F77_CONST_CHAR_ARG2 (&sort, 1),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
217 selector,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
218 F77_CONST_CHAR_ARG2 (&sense, 1),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
219 n, s, n, sdim, pwr, pwi, q, n, rconde, rcondv,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
220 pwork, lwork, piwork, liwork, pbwork, info
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
221 F77_CHAR_ARG_LEN (1)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
222 F77_CHAR_ARG_LEN (1)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
223 F77_CHAR_ARG_LEN (1)));
457
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
224
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
225 return info;
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
226 }
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe]
jwe
parents:
diff changeset
227
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
228 template <>
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
229 octave_idx_type
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
230 schur<FloatMatrix>::init (const FloatMatrix& a, const std::string& ord,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
231 bool calc_unitary)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
232 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
233 octave_idx_type a_nr = a.rows ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
234 octave_idx_type a_nc = a.cols ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
235
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
236 if (a_nr != a_nc)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
237 (*current_liboctave_error_handler) ("SCHUR requires square matrix");
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
238
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
239 if (a_nr == 0)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
240 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
241 schur_mat.clear ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
242 unitary_mat.clear ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
243 return 0;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
244 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
245
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
246 // Workspace requirements may need to be fixed if any of the
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
247 // following change.
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
248
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
249 char jobvs;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
250 char sense = 'N';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
251 char sort = 'N';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
252
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
253 if (calc_unitary)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
254 jobvs = 'V';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
255 else
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
256 jobvs = 'N';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
257
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
258 char ord_char = ord.empty () ? 'U' : ord[0];
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
259
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
260 if (ord_char == 'A' || ord_char == 'D' || ord_char == 'a' || ord_char == 'd')
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
261 sort = 'S';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
262
21663
79414f8e6f24 schur.cc: Declare selector variable as volatile to suppress gcc warning.
Rik <rik@octave.org>
parents: 21301
diff changeset
263 volatile float_selector selector = 0;
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
264 if (ord_char == 'A' || ord_char == 'a')
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
265 selector = select_ana<float>;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
266 else if (ord_char == 'D' || ord_char == 'd')
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
267 selector = select_dig<float>;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
268
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
269 octave_idx_type n = a_nc;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
270 octave_idx_type lwork = 8 * n;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
271 octave_idx_type liwork = 1;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
272 octave_idx_type info;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
273 octave_idx_type sdim;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
274 float rconde;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
275 float rcondv;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
276
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
277 schur_mat = a;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
278
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
279 if (calc_unitary)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
280 unitary_mat.clear (n, n);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
281
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
282 float *s = schur_mat.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
283 float *q = unitary_mat.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
284
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
285 Array<float> wr (dim_vector (n, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
286 float *pwr = wr.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
287
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
288 Array<float> wi (dim_vector (n, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
289 float *pwi = wi.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
290
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
291 Array<float> work (dim_vector (lwork, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
292 float *pwork = work.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
293
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
294 // BWORK is not referenced for the non-ordered Schur routine.
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
295 octave_idx_type ntmp = (ord_char == 'N' || ord_char == 'n') ? 0 : n;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
296 Array<octave_idx_type> bwork (dim_vector (ntmp, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
297 octave_idx_type *pbwork = bwork.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
298
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
299 Array<octave_idx_type> iwork (dim_vector (liwork, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
300 octave_idx_type *piwork = iwork.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
301
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
302 F77_XFCN (sgeesx, SGEESX, (F77_CONST_CHAR_ARG2 (&jobvs, 1),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
303 F77_CONST_CHAR_ARG2 (&sort, 1),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
304 selector,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
305 F77_CONST_CHAR_ARG2 (&sense, 1),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
306 n, s, n, sdim, pwr, pwi, q, n, rconde, rcondv,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
307 pwork, lwork, piwork, liwork, pbwork, info
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
308 F77_CHAR_ARG_LEN (1)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
309 F77_CHAR_ARG_LEN (1)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
310 F77_CHAR_ARG_LEN (1)));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
311
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
312 return info;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
313 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
314
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
315 template <>
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
316 octave_idx_type
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
317 schur<ComplexMatrix>::init (const ComplexMatrix& a, const std::string& ord,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
318 bool calc_unitary)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
319 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
320 octave_idx_type a_nr = a.rows ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
321 octave_idx_type a_nc = a.cols ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
322
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
323 if (a_nr != a_nc)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
324 (*current_liboctave_error_handler) ("SCHUR requires square matrix");
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
325
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
326 if (a_nr == 0)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
327 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
328 schur_mat.clear ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
329 unitary_mat.clear ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
330 return 0;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
331 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
332
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
333 // Workspace requirements may need to be fixed if any of the
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
334 // following change.
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
335
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
336 char jobvs;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
337 char sense = 'N';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
338 char sort = 'N';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
339
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
340 if (calc_unitary)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
341 jobvs = 'V';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
342 else
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
343 jobvs = 'N';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
344
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
345 char ord_char = ord.empty () ? 'U' : ord[0];
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
346
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
347 if (ord_char == 'A' || ord_char == 'D' || ord_char == 'a' || ord_char == 'd')
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
348 sort = 'S';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
349
21663
79414f8e6f24 schur.cc: Declare selector variable as volatile to suppress gcc warning.
Rik <rik@octave.org>
parents: 21301
diff changeset
350 volatile complex_selector selector = 0;
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
351 if (ord_char == 'A' || ord_char == 'a')
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
352 selector = select_ana<Complex>;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
353 else if (ord_char == 'D' || ord_char == 'd')
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
354 selector = select_dig<Complex>;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
355
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
356 octave_idx_type n = a_nc;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
357 octave_idx_type lwork = 8 * n;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
358 octave_idx_type info;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
359 octave_idx_type sdim;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
360 double rconde;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
361 double rcondv;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
362
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
363 schur_mat = a;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
364 if (calc_unitary)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
365 unitary_mat.clear (n, n);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
366
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
367 Complex *s = schur_mat.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
368 Complex *q = unitary_mat.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
369
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
370 Array<double> rwork (dim_vector (n, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
371 double *prwork = rwork.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
372
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
373 Array<Complex> w (dim_vector (n, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
374 Complex *pw = w.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
375
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
376 Array<Complex> work (dim_vector (lwork, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
377 Complex *pwork = work.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
378
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
379 // BWORK is not referenced for non-ordered Schur.
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
380 octave_idx_type ntmp = (ord_char == 'N' || ord_char == 'n') ? 0 : n;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
381 Array<octave_idx_type> bwork (dim_vector (ntmp, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
382 octave_idx_type *pbwork = bwork.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
383
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
384 F77_XFCN (zgeesx, ZGEESX, (F77_CONST_CHAR_ARG2 (&jobvs, 1),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
385 F77_CONST_CHAR_ARG2 (&sort, 1),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
386 selector,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
387 F77_CONST_CHAR_ARG2 (&sense, 1),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
388 n, s, n, sdim, pw, q, n, rconde, rcondv,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
389 pwork, lwork, prwork, pbwork, info
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
390 F77_CHAR_ARG_LEN (1)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
391 F77_CHAR_ARG_LEN (1)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
392 F77_CHAR_ARG_LEN (1)));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
393
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
394 return info;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
395 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
396
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
397 template <>
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
398 schur<ComplexMatrix>
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
399 rsf2csf<ComplexMatrix, Matrix> (const Matrix& s_arg, const Matrix& u_arg)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
400 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
401 ComplexMatrix s (s_arg);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
402 ComplexMatrix u (u_arg);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
403
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
404 octave_idx_type n = s.rows ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
405
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
406 if (s.columns () != n || u.rows () != n || u.columns () != n)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
407 (*current_liboctave_error_handler)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
408 ("rsf2csf: inconsistent matrix dimensions");
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
409
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
410 if (n > 0)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
411 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
412 OCTAVE_LOCAL_BUFFER (double, c, n-1);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
413 OCTAVE_LOCAL_BUFFER (double, sx, n-1);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
414
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
415 F77_XFCN (zrsf2csf, ZRSF2CSF, (n, s.fortran_vec (),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
416 u.fortran_vec (), c, sx));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
417 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
418
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
419 return schur<ComplexMatrix> (s, u);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
420 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
421
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
422 template <>
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
423 octave_idx_type
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
424 schur<FloatComplexMatrix>::init (const FloatComplexMatrix& a,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
425 const std::string& ord, bool calc_unitary)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
426 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
427 octave_idx_type a_nr = a.rows ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
428 octave_idx_type a_nc = a.cols ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
429
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
430 if (a_nr != a_nc)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
431 (*current_liboctave_error_handler) ("SCHUR requires square matrix");
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
432
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
433 if (a_nr == 0)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
434 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
435 schur_mat.clear ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
436 unitary_mat.clear ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
437 return 0;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
438 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
439
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
440 // Workspace requirements may need to be fixed if any of the
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
441 // following change.
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
442
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
443 char jobvs;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
444 char sense = 'N';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
445 char sort = 'N';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
446
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
447 if (calc_unitary)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
448 jobvs = 'V';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
449 else
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
450 jobvs = 'N';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
451
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
452 char ord_char = ord.empty () ? 'U' : ord[0];
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
453
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
454 if (ord_char == 'A' || ord_char == 'D' || ord_char == 'a' || ord_char == 'd')
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
455 sort = 'S';
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
456
21663
79414f8e6f24 schur.cc: Declare selector variable as volatile to suppress gcc warning.
Rik <rik@octave.org>
parents: 21301
diff changeset
457 volatile float_complex_selector selector = 0;
21266
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
458 if (ord_char == 'A' || ord_char == 'a')
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
459 selector = select_ana<FloatComplex>;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
460 else if (ord_char == 'D' || ord_char == 'd')
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
461 selector = select_dig<FloatComplex>;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
462
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
463 octave_idx_type n = a_nc;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
464 octave_idx_type lwork = 8 * n;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
465 octave_idx_type info;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
466 octave_idx_type sdim;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
467 float rconde;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
468 float rcondv;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
469
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
470 schur_mat = a;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
471 if (calc_unitary)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
472 unitary_mat.clear (n, n);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
473
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
474 FloatComplex *s = schur_mat.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
475 FloatComplex *q = unitary_mat.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
476
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
477 Array<float> rwork (dim_vector (n, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
478 float *prwork = rwork.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
479
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
480 Array<FloatComplex> w (dim_vector (n, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
481 FloatComplex *pw = w.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
482
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
483 Array<FloatComplex> work (dim_vector (lwork, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
484 FloatComplex *pwork = work.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
485
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
486 // BWORK is not referenced for non-ordered Schur.
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
487 octave_idx_type ntmp = (ord_char == 'N' || ord_char == 'n') ? 0 : n;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
488 Array<octave_idx_type> bwork (dim_vector (ntmp, 1));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
489 octave_idx_type *pbwork = bwork.fortran_vec ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
490
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
491 F77_XFCN (cgeesx, CGEESX, (F77_CONST_CHAR_ARG2 (&jobvs, 1),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
492 F77_CONST_CHAR_ARG2 (&sort, 1),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
493 selector,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
494 F77_CONST_CHAR_ARG2 (&sense, 1),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
495 n, s, n, sdim, pw, q, n, rconde, rcondv,
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
496 pwork, lwork, prwork, pbwork, info
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
497 F77_CHAR_ARG_LEN (1)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
498 F77_CHAR_ARG_LEN (1)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
499 F77_CHAR_ARG_LEN (1)));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
500
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
501 return info;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
502 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
503
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
504 template <>
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
505 schur<FloatComplexMatrix>
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
506 rsf2csf<FloatComplexMatrix, FloatMatrix> (const FloatMatrix& s_arg, const FloatMatrix& u_arg)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
507 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
508 FloatComplexMatrix s (s_arg);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
509 FloatComplexMatrix u (u_arg);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
510
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
511 octave_idx_type n = s.rows ();
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
512
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
513 if (s.columns () != n || u.rows () != n || u.columns () != n)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
514 (*current_liboctave_error_handler)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
515 ("rsf2csf: inconsistent matrix dimensions");
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
516
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
517 if (n > 0)
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
518 {
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
519 OCTAVE_LOCAL_BUFFER (float, c, n-1);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
520 OCTAVE_LOCAL_BUFFER (float, sx, n-1);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
521
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
522 F77_XFCN (crsf2csf, CRSF2CSF, (n, s.fortran_vec (),
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
523 u.fortran_vec (), c, sx));
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
524 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
525
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
526 return schur<FloatComplexMatrix> (s, u);
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
527 }
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
528
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
529 // Instantiations we need.
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
530
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
531 template class schur<ComplexMatrix>;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
532
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
533 template class schur<FloatComplexMatrix>;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
534
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
535 template class schur<FloatMatrix>;
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
536
e69eaee28737 make better use of templates for Schur decomposition
John W. Eaton <jwe@octave.org>
parents: 21202
diff changeset
537 template class schur<Matrix>;