octave: libinterp/corefcn/qz.cc annotate

annotate libinterp/corefcn/qz.cc @ 27919:1891570abac8

update Octave Project Developers copyright for the new year In files that have the "Octave Project Developers" copyright notice, update for 2020.

author	John W. Eaton <jwe@octave.org>
date	Mon, 06 Jan 2020 22:29:51 -0500
parents	b442ec6dda5c
children	bd51beb6205e

rev	line source
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1 /*
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	2
27919 1891570abac8 update Octave Project Developers copyright for the new year John W. Eaton <jwe@octave.org> parents: 27918 diff changeset	3 Copyright (C) 1998-2020 The Octave Project Developers
27918 b442ec6dda5c use centralized file for copyright info for individual contributors John W. Eaton <jwe@octave.org> parents: 26376 diff changeset	4
b442ec6dda5c use centralized file for copyright info for individual contributors John W. Eaton <jwe@octave.org> parents: 26376 diff changeset	5 See the file COPYRIGHT.md in the top-level directory of this distribution
b442ec6dda5c use centralized file for copyright info for individual contributors John W. Eaton <jwe@octave.org> parents: 26376 diff changeset	6 or <https://octave.org/COPYRIGHT.html/>.
b442ec6dda5c use centralized file for copyright info for individual contributors John W. Eaton <jwe@octave.org> parents: 26376 diff changeset	7
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	8
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	9 This file is part of Octave.
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	10
24534 194eb4bd202b maint: Update punctuation for GPL v3 license text. Rik <rik@octave.org> parents: 24441 diff changeset	11 Octave is free software: you can redistribute it and/or modify it
22755 3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22407 diff changeset	12 under the terms of the GNU General Public License as published by
24534 194eb4bd202b maint: Update punctuation for GPL v3 license text. Rik <rik@octave.org> parents: 24441 diff changeset	13 the Free Software Foundation, either version 3 of the License, or
22755 3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22407 diff changeset	14 (at your option) any later version.
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	15
22755 3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22407 diff changeset	16 Octave is distributed in the hope that it will be useful, but
3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22407 diff changeset	17 WITHOUT ANY WARRANTY; without even the implied warranty of
3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22407 diff changeset	18 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22407 diff changeset	19 GNU General Public License for more details.
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	20
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	21 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: 5823 diff changeset	22 along with Octave; see the file COPYING. If not, see
24534 194eb4bd202b maint: Update punctuation for GPL v3 license text. Rik <rik@octave.org> parents: 24441 diff changeset	23 <https://www.gnu.org/licenses/>.
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	24
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	25 */
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	26
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	27 // Generalized eigenvalue balancing via LAPACK
3911 8389e78e67d4 [project @ 2002-04-28 02:15:38 by jwe] jwe parents: 3887 diff changeset	28
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	29 // Originally written by A. S. Hodel <scotte@eng.auburn.edu>, but is
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	30 // substantially different with the change to use LAPACK.
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	31
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	32 #undef DEBUG
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	33 #undef DEBUG_SORT
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	34 #undef DEBUG_EIG
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	35
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	36 #if defined (HAVE_CONFIG_H)
21301 40de9f8f23a6 Use '#include "config.h"' rather than <config.h>. Rik <rik@octave.org> parents: 21279 diff changeset	37 # include "config.h"
9786 2c279308f6ab fix includes in some src/DLD-FUNCTIONS files John W. Eaton <jwe@octave.org> parents: 9758 diff changeset	38 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	39
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	40 #include <cctype>
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	41 #include <cmath>
4051 b79da8779a0e [project @ 2002-08-17 19:38:32 by jwe] jwe parents: 3911 diff changeset	42
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	43 #if defined (DEBUG_EIG)
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	44 # include <iomanip>
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	45 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	46
4153 6b96ce9f5743 [project @ 2002-11-06 20:38:49 by jwe] jwe parents: 4051 diff changeset	47 #include "f77-fcn.h"
22322 93b3cdd36854 move most f77 function decls to separate header files John W. Eaton <jwe@octave.org> parents: 22317 diff changeset	48 #include "lo-lapack-proto.h"
21279 eb1524b07fe3 better use of templates for qr classes John W. Eaton <jwe@octave.org> parents: 21200 diff changeset	49 #include "qr.h"
4153 6b96ce9f5743 [project @ 2002-11-06 20:38:49 by jwe] jwe parents: 4051 diff changeset	50 #include "quit.h"
6b96ce9f5743 [project @ 2002-11-06 20:38:49 by jwe] jwe parents: 4051 diff changeset	51
15039 e753177cde93 maint: Move non-dynamically linked functions from DLD-FUNCTIONS/ to corefcn/ directory Rik <rik@octave.org> parents: 14876 diff changeset	52 #include "defun.h"
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	53 #include "error.h"
21100 e39e05d90788 Switch gripe_XXX to either err_XXX or warn_XXX naming scheme. Rik <rik@octave.org> parents: 21024 diff changeset	54 #include "errwarn.h"
20940 48b2ad5ee801 maint: Rename oct-obj.[cc\|h] to ovl.[cc\|h] for clarity. Rik <rik@octave.org> parents: 20892 diff changeset	55 #include "ovl.h"
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	56 #if defined (DEBUG) \|\| defined (DEBUG_SORT)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	57 # include "pager.h"
21200 fcac5dbbf9ed maint: Indent #ifdef blocks in libinterp. Rik <rik@octave.org> parents: 21110 diff changeset	58 # include "pr-output.h"
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	59 #endif
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	60
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	61 // FIXME: Matlab does not produce lambda as the first output argument.
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	62 // Compatibility problem?
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	63
15039 e753177cde93 maint: Move non-dynamically linked functions from DLD-FUNCTIONS/ to corefcn/ directory Rik <rik@octave.org> parents: 14876 diff changeset	64 DEFUN (qz, args, nargout,
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	65 doc: /* -- texinfo --
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	66 @deftypefn {} {@var{lambda} =} qz (@var{A}, @var{B})
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	67 @deftypefnx {} {[@var{AA}, @var{BB}, @var{Q}, @var{Z}, @var{V}, @var{W}, @var{lambda}] =} qz (@var{A}, @var{B})
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	68 @deftypefnx {} {[@var{AA}, @var{BB}, @var{Z}] =} qz (@var{A}, @var{B}, @var{opt})
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	69 @deftypefnx {} {[@var{AA}, @var{BB}, @var{Z}, @var{lambda}] =} qz (@var{A}, @var{B}, @var{opt})
25106 d7ad543255c5 doc: Shorten very long first sentences of docstrings (bug #53388). Rik <rik@octave.org> parents: 25103 diff changeset	70 Compute the QZ@tie{}decomposition of a generalized eigenvalue problem.
d7ad543255c5 doc: Shorten very long first sentences of docstrings (bug #53388). Rik <rik@octave.org> parents: 25103 diff changeset	71
d7ad543255c5 doc: Shorten very long first sentences of docstrings (bug #53388). Rik <rik@octave.org> parents: 25103 diff changeset	72 The generalized eigenvalue problem is defined as
d7ad543255c5 doc: Shorten very long first sentences of docstrings (bug #53388). Rik <rik@octave.org> parents: 25103 diff changeset	73
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	74 @tex
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	75 $$A x = \lambda B x$$
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	76 @end tex
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	77 @ifnottex
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	78
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	79 @math{A x = @var{lambda} B x}
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	80
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	81 @end ifnottex
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	82
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	83 There are three calling forms of the function:
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	84
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	85 @enumerate
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	86 @item @code{@var{lambda} = qz (@var{A}, @var{B})}
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	87
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	88 Compute the generalized eigenvalues
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	89 @tex
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	90 $\lambda.$
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	91 @end tex
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	92 @ifnottex
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	93 @var{lambda}.
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	94 @end ifnottex
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	95
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	96 @item @code{[@var{AA}, @var{BB}, @var{Q}, @var{Z}, @var{V}, @var{W}, @var{lambda}] = qz (@var{A}, @var{B})}
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	97
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	98 Compute QZ@tie{}decomposition, generalized eigenvectors, and generalized
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	99 eigenvalues.
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	100 @tex
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	101 $$ AV = BV{ \rm diag }(\lambda) $$
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	102 $$ W^T A = { \rm diag }(\lambda)W^T B $$
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	103 $$ AA = Q^T AZ, BB = Q^T BZ $$
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	104 @end tex
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	105 @ifnottex
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	106
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	107 @example
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	108 @group
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	109
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	110 @var{A} * @var{V} = @var{B} * @var{V} * diag (@var{lambda})
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	111 @var{W}' * @var{A} = diag (@var{lambda}) * @var{W}' * @var{B}
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	112 @var{AA} = @var{Q} * @var{A} * @var{Z}, @var{BB} = @var{Q} * @var{B} * @var{Z}
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	113
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	114 @end group
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	115 @end example
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	116
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	117 @end ifnottex
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	118 with @var{Q} and @var{Z} orthogonal (unitary for complex case).
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	119
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	120 @item @code{[@var{AA}, @var{BB}, @var{Z} @{, @var{lambda}@}] = qz (@var{A}, @var{B}, @var{opt})}
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	121
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	122 As in form 2 above, but allows ordering of generalized eigenpairs for, e.g.,
25003 2365c2661b3c doc: Spellcheck documentation ahead of 4.4 release. Rik <rik@octave.org> parents: 24534 diff changeset	123 solution of discrete time algebraic @nospell{Riccati} equations. Form 3 is not
2365c2661b3c doc: Spellcheck documentation ahead of 4.4 release. Rik <rik@octave.org> parents: 24534 diff changeset	124 available for complex matrices, and does not compute the generalized
2365c2661b3c doc: Spellcheck documentation ahead of 4.4 release. Rik <rik@octave.org> parents: 24534 diff changeset	125 eigenvectors @var{V}, @var{W}, nor the orthogonal matrix @var{Q}.
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	126
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	127 @table @var
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	128 @item opt
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	129 for ordering eigenvalues of the @nospell{GEP} pencil. The leading block of
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	130 the revised pencil contains all eigenvalues that satisfy:
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	131
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	132 @table @asis
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	133 @item @qcode{"N"}
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	134 unordered (default)
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	135
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	136 @item @qcode{"S"}
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	137 small: leading block has all
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	138 @tex
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	139 $\|\lambda\| < 1$
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	140 @end tex
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	141 @ifnottex
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	142 \|@var{lambda}\| < 1
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	143 @end ifnottex
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	144
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	145 @item @qcode{"B"}
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	146 big: leading block has all
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	147 @tex
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	148 $\|\lambda\| \geq 1$
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	149 @end tex
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	150 @ifnottex
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	151 \|@var{lambda}\| @geq{} 1
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	152 @end ifnottex
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	153
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	154 @item @qcode{"-"}
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	155 negative real part: leading block has all eigenvalues in the open left
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	156 half-plane
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	157
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	158 @item @qcode{"+"}
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	159 non-negative real part: leading block has all eigenvalues in the closed right
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	160 half-plane
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	161 @end table
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	162 @end table
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	163 @end enumerate
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	164
112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	165 Note: @code{qz} performs permutation balancing, but not scaling
24440 f8bbaacefc33 doc: Fix various inconsistencies in manual (bug #52712). Rik <rik@octave.org> parents: 23428 diff changeset	166 (@pxref{XREFbalance,,balance}), which may be lead to less accurate results than
23428 1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	167 @code{eig}. The order of output arguments was selected for compatibility with
1bc0e610e293 doc: Redo docstring for qz (bug #50846). Rik <rik@octave.org> parents: 23219 diff changeset	168 @sc{matlab}.
25734 c7095a755185 New function "ordeig" (patch #9670) Sébastien Villemot <sebastien@debian.org> parents: 25570 diff changeset	169 @seealso{eig, ordeig, balance, lu, chol, hess, qr, qzhess, schur, svd}
21966 112b20240c87 move docstrings in C++ files out of C strings and into comments John W. Eaton <jwe@octave.org> parents: 21751 diff changeset	170 @end deftypefn */)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	171 {
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	172 int nargin = args.length ();
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	173
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	174 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	175 octave_stdout << "qz: nargin = " << nargin
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	176 << ", nargout = " << nargout << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	177 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	178
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	179 if (nargin < 2 \|\| nargin > 3 \|\| nargout > 7)
20802 8bb38ba1bad6 eliminate return statements after calls to print_usage John W. Eaton <jwe@octave.org> parents: 20700 diff changeset	180 print_usage ();
8bb38ba1bad6 eliminate return statements after calls to print_usage John W. Eaton <jwe@octave.org> parents: 20700 diff changeset	181
8bb38ba1bad6 eliminate return statements after calls to print_usage John W. Eaton <jwe@octave.org> parents: 20700 diff changeset	182 if (nargin == 3 && (nargout < 3 \|\| nargout > 4))
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	183 error ("qz: invalid number of output arguments for form 3 call");
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	184
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	185 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	186 octave_stdout << "qz: determine ordering option" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	187 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	188
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	189 // Determine ordering option.
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	190
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	191 char ord_job = 'N';
24116 e1020353fff5 qz: initialize optional ordering flag safely Mike Miller <mtmiller@octave.org> parents: 23826 diff changeset	192 double safmin = 0.0;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	193
24116 e1020353fff5 qz: initialize optional ordering flag safely Mike Miller <mtmiller@octave.org> parents: 23826 diff changeset	194 if (nargin == 3)
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	195 {
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	196 std::string opt = args(2).xstring_value ("qz: OPT must be a string");
8927 d75f4ee0538d qz.cc (Fqz): avoid maybe clobbred by vfork warning from GCC John W. Eaton <jwe@octave.org> parents: 8920 diff changeset	197
24116 e1020353fff5 qz: initialize optional ordering flag safely Mike Miller <mtmiller@octave.org> parents: 23826 diff changeset	198 if (opt.empty ())
e1020353fff5 qz: initialize optional ordering flag safely Mike Miller <mtmiller@octave.org> parents: 23826 diff changeset	199 error ("qz: OPT must be a non-empty string");
e1020353fff5 qz: initialize optional ordering flag safely Mike Miller <mtmiller@octave.org> parents: 23826 diff changeset	200
e1020353fff5 qz: initialize optional ordering flag safely Mike Miller <mtmiller@octave.org> parents: 23826 diff changeset	201 ord_job = std::toupper (opt[0]);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	202
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	203 std::string valid_opts = "NSB+-";
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	204
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	205 if (valid_opts.find_first_of (ord_job) == std::string::npos)
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	206 error ("qz: invalid order option '%c'", ord_job);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	207
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	208 // overflow constant required by dlag2
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	209 F77_FUNC (xdlamch, XDLAMCH) (F77_CONST_CHAR_ARG2 ("S", 1),
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	210 safmin
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	211 F77_CHAR_ARG_LEN (1));
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	212
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	213 #if defined (DEBUG_EIG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	214 octave_stdout << "qz: initial value of safmin="
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	215 << setiosflags (std::ios::scientific)
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	216 << safmin << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	217 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	218
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	219 // Some machines (e.g., DEC alpha) get safmin = 0;
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	220 // for these, use eps instead to avoid problems in dlag2.
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	221 if (safmin == 0)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	222 {
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	223 #if defined (DEBUG_EIG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	224 octave_stdout << "qz: DANGER WILL ROBINSON: safmin is 0!"
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	225 << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	226 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	227
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	228 F77_FUNC (xdlamch, XDLAMCH) (F77_CONST_CHAR_ARG2 ("E", 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	229 safmin
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	230 F77_CHAR_ARG_LEN (1));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	231
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	232 #if defined (DEBUG_EIG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	233 octave_stdout << "qz: safmin set to "
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	234 << setiosflags (std::ios::scientific)
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	235 << safmin << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	236 #endif
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	237 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	238 }
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	239
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	240 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	241 octave_stdout << "qz: check matrix A" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	242 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	243
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	244 // Matrix A: check dimensions.
22988 cd33c785e80e put to_f77_int inside octave namespace John W. Eaton <jwe@octave.org> parents: 22964 diff changeset	245 F77_INT nn = octave::to_f77_int (args(0).rows ());
cd33c785e80e put to_f77_int inside octave namespace John W. Eaton <jwe@octave.org> parents: 22964 diff changeset	246 F77_INT nc = octave::to_f77_int (args(0).columns ());
20892 c07bee629973 2015 Code Sprint: use ovl (). Rik <rik@octave.org> parents: 20853 diff changeset	247
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	248 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	249 octave_stdout << "Matrix A dimensions: (" << nn << ',' << nc << ')'
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	250 << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	251 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	252
23577 80c42f4cca13 maint: Deprecate is_empty and replace with isempty. Rik <rik@octave.org> parents: 23477 diff changeset	253 if (args(0).isempty ())
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	254 {
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	255 warn_empty_arg ("qz: A");
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	256 return octave_value_list (2, Matrix ());
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	257 }
22964 0c12642be005 use F77_INT instead of octave_idx_type for libinterp QZ function John W. Eaton <jwe@octave.org> parents: 22860 diff changeset	258 else if (nc != nn)
21110 3d0d84305600 Use err_square_matrix_required more widely. Rik <rik@octave.org> parents: 21100 diff changeset	259 err_square_matrix_required ("qz", "A");
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	260
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	261 // Matrix A: get value.
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	262 Matrix aa;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	263 ComplexMatrix caa;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	264
23581 c3075ae020e1 maint: Deprecate is_complex_type and replace with iscomplex. Rik <rik@octave.org> parents: 23577 diff changeset	265 if (args(0).iscomplex ())
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	266 caa = args(0).complex_matrix_value ();
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	267 else
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	268 aa = args(0).matrix_value ();
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	269
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	270 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	271 octave_stdout << "qz: check matrix B" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	272 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	273
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	274 // Extract argument 2 (bb, or cbb if complex).
22988 cd33c785e80e put to_f77_int inside octave namespace John W. Eaton <jwe@octave.org> parents: 22964 diff changeset	275 F77_INT b_nr = octave::to_f77_int (args(1).rows ());
cd33c785e80e put to_f77_int inside octave namespace John W. Eaton <jwe@octave.org> parents: 22964 diff changeset	276 F77_INT b_nc = octave::to_f77_int (args(1).columns ());
22964 0c12642be005 use F77_INT instead of octave_idx_type for libinterp QZ function John W. Eaton <jwe@octave.org> parents: 22860 diff changeset	277
0c12642be005 use F77_INT instead of octave_idx_type for libinterp QZ function John W. Eaton <jwe@octave.org> parents: 22860 diff changeset	278 if (nn != b_nc \|\| nn != b_nr)
21100 e39e05d90788 Switch gripe_XXX to either err_XXX or warn_XXX naming scheme. Rik <rik@octave.org> parents: 21024 diff changeset	279 err_nonconformant ();
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	280
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	281 Matrix bb;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	282 ComplexMatrix cbb;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	283
23581 c3075ae020e1 maint: Deprecate is_complex_type and replace with iscomplex. Rik <rik@octave.org> parents: 23577 diff changeset	284 if (args(1).iscomplex ())
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	285 cbb = args(1).complex_matrix_value ();
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	286 else
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	287 bb = args(1).matrix_value ();
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	288
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	289 // Both matrices loaded, now check whether to calculate complex or real val.
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	290
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	291 bool complex_case = (args(0).iscomplex () \|\| args(1).iscomplex ());
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	292
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	293 if (nargin == 3 && complex_case)
20831 35241c4b696c eliminate return statements after calls to error John W. Eaton <jwe@octave.org> parents: 20802 diff changeset	294 error ("qz: cannot re-order complex qz decomposition");
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	295
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	296 // First, declare variables used in both the real and complex cases.
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	297 // FIXME: There are a lot of excess variables declared.
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	298 // Probably a better way to handle this.
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	299 Matrix QQ (nn,nn), ZZ (nn,nn), VR (nn,nn), VL (nn,nn);
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	300 RowVector alphar (nn), alphai (nn), betar (nn);
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	301 ComplexRowVector xalpha (nn), xbeta (nn);
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	302 ComplexMatrix CQ (nn,nn), CZ (nn,nn), CVR (nn,nn), CVL (nn,nn);
22964 0c12642be005 use F77_INT instead of octave_idx_type for libinterp QZ function John W. Eaton <jwe@octave.org> parents: 22860 diff changeset	303 F77_INT ilo, ihi, info;
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	304 char comp_q = (nargout >= 3 ? 'V' : 'N');
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	305 char comp_z = ((nargout >= 4 \|\| nargin == 3)? 'V' : 'N');
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	306
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	307 // Initialize Q, Z to identity matrix if either is needed
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	308 if (comp_q == 'V' \|\| comp_z == 'V')
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	309 {
23477 3530b956d707 maint: Strip trailing whitespace from source files. John W. Eaton <jwe@octave.org> parents: 23457 diff changeset	310 double *QQptr = QQ.fortran_vec ();
3530b956d707 maint: Strip trailing whitespace from source files. John W. Eaton <jwe@octave.org> parents: 23457 diff changeset	311 double *ZZptr = ZZ.fortran_vec ();
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	312 std::fill_n (QQptr, QQ.numel (), 0.0);
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	313 std::fill_n (ZZptr, ZZ.numel (), 0.0);
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	314 for (F77_INT i = 0; i < nn; i++)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	315 {
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	316 QQ(i,i) = 1.0;
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	317 ZZ(i,i) = 1.0;
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	318 }
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	319 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	320
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	321 // Always perform permutation balancing.
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	322 const char bal_job = 'P';
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	323 RowVector lscale (nn), rscale (nn), work (6 * nn), rwork (nn);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	324
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	325 if (complex_case)
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	326 {
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	327 #if defined (DEBUG)
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	328 if (comp_q == 'V')
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	329 octave_stdout << "qz: performing balancing; CQ =\n" << CQ << std::endl;
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	330 #endif
23582 0cc2011d800e maint: Deprecate is_real_type and replace with isreal. Rik <rik@octave.org> parents: 23581 diff changeset	331 if (args(0).isreal ())
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	332 caa = ComplexMatrix (aa);
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	333
23582 0cc2011d800e maint: Deprecate is_real_type and replace with isreal. Rik <rik@octave.org> parents: 23581 diff changeset	334 if (args(1).isreal ())
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	335 cbb = ComplexMatrix (bb);
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	336
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	337 if (comp_q == 'V')
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	338 CQ = ComplexMatrix (QQ);
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	339
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	340 if (comp_z == 'V')
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	341 CZ = ComplexMatrix (ZZ);
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	342
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	343 F77_XFCN (zggbal, ZGGBAL,
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	344 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
22407 34ce5be04942 maint: Style check C++ code in libinterp/. Rik <rik@octave.org> parents: 22372 diff changeset	345 nn, F77_DBLE_CMPLX_ARG (caa.fortran_vec ()), nn,
34ce5be04942 maint: Style check C++ code in libinterp/. Rik <rik@octave.org> parents: 22372 diff changeset	346 F77_DBLE_CMPLX_ARG (cbb.fortran_vec ()),
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	347 nn, ilo, ihi, lscale.fortran_vec (),
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	348 rscale.fortran_vec (), work.fortran_vec (), info
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	349 F77_CHAR_ARG_LEN (1)));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	350 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	351 else
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	352 {
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	353 #if defined (DEBUG)
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	354 if (comp_q == 'V')
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	355 octave_stdout << "qz: performing balancing; QQ =\n" << QQ << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	356 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	357
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	358 F77_XFCN (dggbal, DGGBAL,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	359 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	360 nn, aa.fortran_vec (), nn, bb.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	361 nn, ilo, ihi, lscale.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	362 rscale.fortran_vec (), work.fortran_vec (), info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	363 F77_CHAR_ARG_LEN (1)));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	364 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	365
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	366 // Only permutation balance above is done. Skip scaling balance.
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	367
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	368 #if 0
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	369 // Since we just want the balancing matrices, we can use dggbal
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	370 // for both the real and complex cases; left first
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	371
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	372 if (comp_q == 'V')
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	373 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	374 F77_XFCN (dggbak, DGGBAK,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	375 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	376 F77_CONST_CHAR_ARG2 ("L", 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	377 nn, ilo, ihi, lscale.data (), rscale.data (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	378 nn, QQ.fortran_vec (), nn, info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	379 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	380 F77_CHAR_ARG_LEN (1)));
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	381
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	382 #if defined (DEBUG)
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	383 if (comp_q == 'V')
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	384 octave_stdout << "qz: balancing done; QQ =\n" << QQ << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	385 #endif
25103 078b795c5219 maint: style check C++ ahead of 4.4 release. Rik <rik@octave.org> parents: 25054 diff changeset	386 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	387
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	388 // then right
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	389 if (comp_z == 'V')
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	390 {
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	391 F77_XFCN (dggbak, DGGBAK,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	392 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	393 F77_CONST_CHAR_ARG2 ("R", 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	394 nn, ilo, ihi, lscale.data (), rscale.data (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	395 nn, ZZ.fortran_vec (), nn, info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	396 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	397 F77_CHAR_ARG_LEN (1)));
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	398
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	399 #if defined (DEBUG)
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	400 if (comp_z == 'V')
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	401 octave_stdout << "qz: balancing done; ZZ=\n" << ZZ << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	402 #endif
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	403 }
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	404 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	405
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	406 char qz_job = (nargout < 2 ? 'E' : 'S');
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	407
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	408 if (complex_case)
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	409 {
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	410 // Complex case.
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	411
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	412 // The QR decomposition of cbb.
22317 6ca3acf5fad8 move some new numeric classes to namespace octave::math John W. Eaton <jwe@octave.org> parents: 22137 diff changeset	413 octave::math::qr<ComplexMatrix> cbqr (cbb);
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	414 // The R matrix of QR decomposition for cbb.
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	415 cbb = cbqr.R ();
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	416 // (Q*)caa for following work.
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	417 caa = (cbqr.Q ().hermitian ()) * caa;
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	418 CQ = CQ * cbqr.Q ();
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	419
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	420 F77_XFCN (zgghrd, ZGGHRD,
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	421 (F77_CONST_CHAR_ARG2 (&comp_q, 1),
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	422 F77_CONST_CHAR_ARG2 (&comp_z, 1),
22135 407c66ae1e20 reduce warnings from GCC's link-time optimization feature (bug #48531) John W. Eaton <jwe@octave.org> parents: 22133 diff changeset	423 nn, ilo, ihi, F77_DBLE_CMPLX_ARG (caa.fortran_vec ()),
22407 34ce5be04942 maint: Style check C++ code in libinterp/. Rik <rik@octave.org> parents: 22372 diff changeset	424 nn, F77_DBLE_CMPLX_ARG (cbb.fortran_vec ()), nn,
34ce5be04942 maint: Style check C++ code in libinterp/. Rik <rik@octave.org> parents: 22372 diff changeset	425 F77_DBLE_CMPLX_ARG (CQ.fortran_vec ()), nn,
22135 407c66ae1e20 reduce warnings from GCC's link-time optimization feature (bug #48531) John W. Eaton <jwe@octave.org> parents: 22133 diff changeset	426 F77_DBLE_CMPLX_ARG (CZ.fortran_vec ()), nn, info
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	427 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	428 F77_CHAR_ARG_LEN (1)));
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	429
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	430 ComplexRowVector cwork (nn);
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	431
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	432 F77_XFCN (zhgeqz, ZHGEQZ,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	433 (F77_CONST_CHAR_ARG2 (&qz_job, 1),
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	434 F77_CONST_CHAR_ARG2 (&comp_q, 1),
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	435 F77_CONST_CHAR_ARG2 (&comp_z, 1),
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11570 diff changeset	436 nn, ilo, ihi,
22135 407c66ae1e20 reduce warnings from GCC's link-time optimization feature (bug #48531) John W. Eaton <jwe@octave.org> parents: 22133 diff changeset	437 F77_DBLE_CMPLX_ARG (caa.fortran_vec ()), nn,
407c66ae1e20 reduce warnings from GCC's link-time optimization feature (bug #48531) John W. Eaton <jwe@octave.org> parents: 22133 diff changeset	438 F77_DBLE_CMPLX_ARG (cbb.fortran_vec ()), nn,
407c66ae1e20 reduce warnings from GCC's link-time optimization feature (bug #48531) John W. Eaton <jwe@octave.org> parents: 22133 diff changeset	439 F77_DBLE_CMPLX_ARG (xalpha.fortran_vec ()),
407c66ae1e20 reduce warnings from GCC's link-time optimization feature (bug #48531) John W. Eaton <jwe@octave.org> parents: 22133 diff changeset	440 F77_DBLE_CMPLX_ARG (xbeta.fortran_vec ()),
407c66ae1e20 reduce warnings from GCC's link-time optimization feature (bug #48531) John W. Eaton <jwe@octave.org> parents: 22133 diff changeset	441 F77_DBLE_CMPLX_ARG (CQ.fortran_vec ()), nn,
407c66ae1e20 reduce warnings from GCC's link-time optimization feature (bug #48531) John W. Eaton <jwe@octave.org> parents: 22133 diff changeset	442 F77_DBLE_CMPLX_ARG (CZ.fortran_vec ()), nn,
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	443 F77_DBLE_CMPLX_ARG (cwork.fortran_vec ()), nn,
25103 078b795c5219 maint: style check C++ ahead of 4.4 release. Rik <rik@octave.org> parents: 25054 diff changeset	444 rwork.fortran_vec (), info
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	445 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	446 F77_CHAR_ARG_LEN (1)
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11570 diff changeset	447 F77_CHAR_ARG_LEN (1)));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	448
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	449 if (comp_q == 'V')
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	450 {
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	451 // Left eigenvector.
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	452 F77_XFCN (zggbak, ZGGBAK,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	453 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	454 F77_CONST_CHAR_ARG2 ("L", 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	455 nn, ilo, ihi, lscale.data (), rscale.data (),
22135 407c66ae1e20 reduce warnings from GCC's link-time optimization feature (bug #48531) John W. Eaton <jwe@octave.org> parents: 22133 diff changeset	456 nn, F77_DBLE_CMPLX_ARG (CQ.fortran_vec ()), nn, info
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	457 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	458 F77_CHAR_ARG_LEN (1)));
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11570 diff changeset	459 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	460
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	461 if (comp_z == 'V')
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	462 {
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	463 // Right eigenvector.
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	464 F77_XFCN (zggbak, ZGGBAK,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	465 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	466 F77_CONST_CHAR_ARG2 ("R", 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	467 nn, ilo, ihi, lscale.data (), rscale.data (),
22135 407c66ae1e20 reduce warnings from GCC's link-time optimization feature (bug #48531) John W. Eaton <jwe@octave.org> parents: 22133 diff changeset	468 nn, F77_DBLE_CMPLX_ARG (CZ.fortran_vec ()), nn, info
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	469 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	470 F77_CHAR_ARG_LEN (1)));
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11570 diff changeset	471 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	472
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	473 }
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	474 else
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	475 {
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	476 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	477 octave_stdout << "qz: performing qr decomposition of bb" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	478 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	479
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	480 // Compute the QR factorization of bb.
22317 6ca3acf5fad8 move some new numeric classes to namespace octave::math John W. Eaton <jwe@octave.org> parents: 22137 diff changeset	481 octave::math::qr<Matrix> bqr (bb);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	482
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	483 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	484 octave_stdout << "qz: qr (bb) done; now performing qz decomposition"
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	485 << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	486 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	487
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	488 bb = bqr.R ();
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	489
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	490 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	491 octave_stdout << "qz: extracted bb" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	492 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	493
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	494 aa = (bqr.Q ()).transpose () * aa;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	495
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	496 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	497 octave_stdout << "qz: updated aa " << std::endl;
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	498 octave_stdout << "bqr.Q () =\n" << bqr.Q () << std::endl;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	499
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	500 if (comp_q == 'V')
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	501 octave_stdout << "QQ =" << QQ << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	502 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	503
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	504 if (comp_q == 'V')
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	505 QQ = QQ * bqr.Q ();
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	506
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	507 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	508 octave_stdout << "qz: precursors done..." << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	509 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	510
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	511 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	512 octave_stdout << "qz: comp_q = " << comp_q << ", comp_z = " << comp_z
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	513 << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	514 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	515
21562 6c2fd62db1f7 maint: Eliminate accidental double spaces in code. Rik <rik@octave.org> parents: 21301 diff changeset	516 // Reduce to generalized Hessenberg form.
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	517 F77_XFCN (dgghrd, DGGHRD,
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	518 (F77_CONST_CHAR_ARG2 (&comp_q, 1),
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	519 F77_CONST_CHAR_ARG2 (&comp_z, 1),
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	520 nn, ilo, ihi, aa.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	521 nn, bb.fortran_vec (), nn, QQ.fortran_vec (), nn,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	522 ZZ.fortran_vec (), nn, info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	523 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	524 F77_CHAR_ARG_LEN (1)));
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	525
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	526 // Check if just computing generalized eigenvalues,
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	527 // or if we're actually computing the decomposition.
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	528
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	529 // Reduce to generalized Schur form.
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	530 F77_XFCN (dhgeqz, DHGEQZ,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	531 (F77_CONST_CHAR_ARG2 (&qz_job, 1),
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	532 F77_CONST_CHAR_ARG2 (&comp_q, 1),
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	533 F77_CONST_CHAR_ARG2 (&comp_z, 1),
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	534 nn, ilo, ihi, aa.fortran_vec (), nn, bb.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	535 nn, alphar.fortran_vec (), alphai.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	536 betar.fortran_vec (), QQ.fortran_vec (), nn,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	537 ZZ.fortran_vec (), nn, work.fortran_vec (), nn, info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	538 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	539 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	540 F77_CHAR_ARG_LEN (1)));
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	541
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	542 if (comp_q == 'V')
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	543 {
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	544 F77_XFCN (dggbak, DGGBAK,
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	545 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	546 F77_CONST_CHAR_ARG2 ("L", 1),
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	547 nn, ilo, ihi, lscale.data (), rscale.data (),
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	548 nn, QQ.fortran_vec (), nn, info
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	549 F77_CHAR_ARG_LEN (1)
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	550 F77_CHAR_ARG_LEN (1)));
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	551
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	552 #if defined (DEBUG)
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	553 if (comp_q == 'V')
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	554 octave_stdout << "qz: balancing done; QQ=\n" << QQ << std::endl;
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	555 #endif
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	556 }
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	557
17787 175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	558 // then right
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	559 if (comp_z == 'V')
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	560 {
17787 175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	561 F77_XFCN (dggbak, DGGBAK,
175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	562 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	563 F77_CONST_CHAR_ARG2 ("R", 1),
175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	564 nn, ilo, ihi, lscale.data (), rscale.data (),
175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	565 nn, ZZ.fortran_vec (), nn, info
175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	566 F77_CHAR_ARG_LEN (1)
175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	567 F77_CHAR_ARG_LEN (1)));
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	568
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	569 #if defined (DEBUG)
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	570 if (comp_z == 'V')
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	571 octave_stdout << "qz: balancing done; ZZ=\n" << ZZ << std::endl;
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	572 #endif
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	573 }
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	574
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	575 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	576
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	577 // Order the QZ decomposition?
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	578 if (ord_job != 'N')
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	579 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	580 if (complex_case)
20831 35241c4b696c eliminate return statements after calls to error John W. Eaton <jwe@octave.org> parents: 20802 diff changeset	581 // Probably not needed, but better be safe.
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	582 error ("qz: cannot re-order complex QZ decomposition");
20831 35241c4b696c eliminate return statements after calls to error John W. Eaton <jwe@octave.org> parents: 20802 diff changeset	583
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	584 #if defined (DEBUG_SORT)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	585 octave_stdout << "qz: ordering eigenvalues: ord_job = "
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	586 << ord_job << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	587 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	588
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	589 Array<F77_LOGICAL> select (dim_vector (nn, 1));
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	590
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	591 for (int j = 0; j < nn; j++)
20831 35241c4b696c eliminate return statements after calls to error John W. Eaton <jwe@octave.org> parents: 20802 diff changeset	592 {
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	593 switch (ord_job)
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	594 {
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	595 case 'S':
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	596 select(j) = alphar(j)alphar(j) + alphai(j)alphai(j) < betar(j)*betar(j);
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	597 break;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	598
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	599 case 'B':
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	600 select(j) = alphar(j)alphar(j) + alphai(j)alphai(j) >= betar(j)*betar(j);
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	601 break;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	602
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	603 case '+':
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	604 select(j) = alphar(j) * betar(j) >= 0;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	605 break;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	606
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	607 case '-':
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	608 select(j) = alphar(j) * betar(j) < 0;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	609 break;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	610
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	611 default:
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	612 // Invalid order option
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	613 // (should never happen since options were checked at the top).
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	614 panic_impossible ();
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	615 break;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	616 }
20831 35241c4b696c eliminate return statements after calls to error John W. Eaton <jwe@octave.org> parents: 20802 diff changeset	617 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	618
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	619 F77_LOGICAL wantq = 0, wantz = (comp_z == 'V');
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	620 F77_INT ijob = 0, mm, lrwork3 = 4*nn+16, liwork = nn;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	621 F77_DBLE pl, pr;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	622 RowVector rwork3(lrwork3);
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	623 Array<F77_INT> iwork (dim_vector (liwork, 1));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	624
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	625 F77_XFCN (dtgsen, DTGSEN,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	626 (ijob, wantq, wantz,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	627 select.fortran_vec (), nn,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	628 aa.fortran_vec (), nn,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	629 bb.fortran_vec (), nn,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	630 alphar.fortran_vec (),
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	631 alphai.fortran_vec (),
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	632 betar.fortran_vec (),
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	633 nullptr, nn,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	634 ZZ.fortran_vec (), nn,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	635 mm,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	636 pl, pr,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	637 nullptr,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	638 rwork3.fortran_vec (), lrwork3,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	639 iwork.fortran_vec (), liwork,
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	640 info));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	641
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	642 #if defined (DEBUG_SORT)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	643 octave_stdout << "qz: back from dtgsen: aa =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	644 octave_print_internal (octave_stdout, aa);
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	645 octave_stdout << "\nbb =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	646 octave_print_internal (octave_stdout, bb);
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	647 if (comp_z == 'V')
20831 35241c4b696c eliminate return statements after calls to error John W. Eaton <jwe@octave.org> parents: 20802 diff changeset	648 {
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	649 octave_stdout << "\nZZ =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	650 octave_print_internal (octave_stdout, ZZ);
20831 35241c4b696c eliminate return statements after calls to error John W. Eaton <jwe@octave.org> parents: 20802 diff changeset	651 }
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	652 octave_stdout << "\nqz: info=" << info;
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	653 octave_stdout << "\nalphar =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	654 octave_print_internal (octave_stdout, Matrix (alphar));
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	655 octave_stdout << "\nalphai =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	656 octave_print_internal (octave_stdout, Matrix (alphai));
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	657 octave_stdout << "\nbeta =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	658 octave_print_internal (octave_stdout, Matrix (betar));
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	659 octave_stdout << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	660 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	661 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	662
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	663 // Compute the generalized eigenvalues as well?
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	664 ComplexColumnVector gev;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	665
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	666 if (nargout < 2 \|\| nargout == 7 \|\| (nargin == 3 && nargout == 4))
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	667 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	668 if (complex_case)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	669 {
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	670 ComplexColumnVector tmp (nn);
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	671
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	672 for (F77_INT i = 0; i < nn; i++)
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	673 tmp(i) = xalpha(i) / xbeta(i);
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	674
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	675 gev = tmp;
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	676 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	677 else
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	678 {
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	679 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	680 octave_stdout << "qz: computing generalized eigenvalues" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	681 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	682
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	683 // Return finite generalized eigenvalues.
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	684 ComplexColumnVector tmp (nn);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	685
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	686 F77_INT cnt = 0;
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	687 for (F77_INT i = 0; i < nn; i++)
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	688 if (betar(i) != 0)
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	689 tmp(cnt++) = Complex (alphar(i), alphai(i)) / betar(i);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	690
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	691 tmp.resize (cnt); // Trim vector to number of return values
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	692
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	693 gev = tmp;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	694 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	695 }
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	696
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	697 // Right, left eigenvector matrices.
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	698 if (nargout >= 5)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	699 {
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	700 // Which side to compute?
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	701 char side = (nargout == 5 ? 'R' : 'B');
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	702 // Compute all of them and backtransform
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	703 char howmany = 'B';
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	704 // Dummy pointer; select is not used.
23457 21baad6b35c4 maint: Use C++11 nullptr rather than 0 or NULL when possible. Rik <rik@octave.org> parents: 23449 diff changeset	705 F77_INT *select = nullptr;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	706
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	707 if (complex_case)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	708 {
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	709 CVL = CQ;
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	710 CVR = CZ;
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	711 ComplexRowVector cwork2 (2 * nn);
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	712 RowVector rwork2 (8 * nn);
22964 0c12642be005 use F77_INT instead of octave_idx_type for libinterp QZ function John W. Eaton <jwe@octave.org> parents: 22860 diff changeset	713 F77_INT m;
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	714
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	715 F77_XFCN (ztgevc, ZTGEVC,
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	716 (F77_CONST_CHAR_ARG2 (&side, 1),
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	717 F77_CONST_CHAR_ARG2 (&howmany, 1),
22407 34ce5be04942 maint: Style check C++ code in libinterp/. Rik <rik@octave.org> parents: 22372 diff changeset	718 select, nn, F77_DBLE_CMPLX_ARG (caa.fortran_vec ()), nn,
34ce5be04942 maint: Style check C++ code in libinterp/. Rik <rik@octave.org> parents: 22372 diff changeset	719 F77_DBLE_CMPLX_ARG (cbb.fortran_vec ()),
34ce5be04942 maint: Style check C++ code in libinterp/. Rik <rik@octave.org> parents: 22372 diff changeset	720 nn, F77_DBLE_CMPLX_ARG (CVL.fortran_vec ()), nn,
34ce5be04942 maint: Style check C++ code in libinterp/. Rik <rik@octave.org> parents: 22372 diff changeset	721 F77_DBLE_CMPLX_ARG (CVR.fortran_vec ()), nn, nn,
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	722 m, F77_DBLE_CMPLX_ARG (cwork2.fortran_vec ()),
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	723 rwork2.fortran_vec (), info
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	724 F77_CHAR_ARG_LEN (1)
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11570 diff changeset	725 F77_CHAR_ARG_LEN (1)));
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	726 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	727 else
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	728 {
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	729 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	730 octave_stdout << "qz: computing generalized eigenvectors" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	731 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	732
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	733 VL = QQ;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	734 VR = ZZ;
22964 0c12642be005 use F77_INT instead of octave_idx_type for libinterp QZ function John W. Eaton <jwe@octave.org> parents: 22860 diff changeset	735 F77_INT m;
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	736
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	737 F77_XFCN (dtgevc, DTGEVC,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	738 (F77_CONST_CHAR_ARG2 (&side, 1),
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	739 F77_CONST_CHAR_ARG2 (&howmany, 1),
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	740 select, nn, aa.fortran_vec (), nn, bb.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	741 nn, VL.fortran_vec (), nn, VR.fortran_vec (), nn, nn,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	742 m, work.fortran_vec (), info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	743 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	744 F77_CHAR_ARG_LEN (1)));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	745
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	746 // Now construct the complex form of VV, WW.
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	747 F77_INT j = 0;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	748
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	749 while (j < nn)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	750 {
22860 0b1e25cc4457 eliminate use of OCTAVE_QUIT macro in C++ sources John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	751 octave_quit ();
4153 6b96ce9f5743 [project @ 2002-11-06 20:38:49 by jwe] jwe parents: 4051 diff changeset	752
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	753 // See if real or complex eigenvalue.
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	754
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	755 // Column increment; assume complex eigenvalue.
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	756 int cinc = 2;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	757
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	758 if (j == (nn-1))
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	759 // Single column.
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	760 cinc = 1;
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	761 else if (aa(j+1,j) == 0)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	762 cinc = 1;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	763
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	764 // Now copy the eigenvector (s) to CVR, CVL.
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	765 if (cinc == 1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	766 {
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	767 for (F77_INT i = 0; i < nn; i++)
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	768 CVR(i,j) = VR(i,j);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	769
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	770 if (side == 'B')
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	771 for (F77_INT i = 0; i < nn; i++)
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	772 CVL(i,j) = VL(i,j);
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	773 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	774 else
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	775 {
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	776 // Double column; complex vector.
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	777
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	778 for (F77_INT i = 0; i < nn; i++)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	779 {
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	780 CVR(i,j) = Complex (VR(i,j), VR(i,j+1));
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	781 CVR(i,j+1) = Complex (VR(i,j), -VR(i,j+1));
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	782 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	783
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	784 if (side == 'B')
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	785 for (F77_INT i = 0; i < nn; i++)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	786 {
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	787 CVL(i,j) = Complex (VL(i,j), VL(i,j+1));
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	788 CVL(i,j+1) = Complex (VL(i,j), -VL(i,j+1));
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	789 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	790 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	791
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	792 // Advance to next eigenvectors (if any).
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	793 j += cinc;
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	794 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	795 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	796 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	797
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	798 octave_value_list retval (nargout);
32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	799
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	800 switch (nargout)
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	801 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	802 case 7:
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	803 retval(6) = gev;
23826 d69021d58a61 avoid fallthrough warnings John W. Eaton <jwe@octave.org> parents: 23807 diff changeset	804 OCTAVE_FALLTHROUGH;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	805
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	806 case 6:
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	807 // Return left eigenvectors.
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	808 retval(5) = CVL;
23826 d69021d58a61 avoid fallthrough warnings John W. Eaton <jwe@octave.org> parents: 23807 diff changeset	809 OCTAVE_FALLTHROUGH;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	810
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	811 case 5:
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	812 // Return right eigenvectors.
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	813 retval(4) = CVR;
23826 d69021d58a61 avoid fallthrough warnings John W. Eaton <jwe@octave.org> parents: 23807 diff changeset	814 OCTAVE_FALLTHROUGH;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	815
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	816 case 4:
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	817 if (nargin == 3)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	818 {
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	819 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	820 octave_stdout << "qz: sort: retval(3) = gev =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	821 octave_print_internal (octave_stdout, ComplexMatrix (gev));
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	822 octave_stdout << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	823 #endif
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	824 retval(3) = gev;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	825 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	826 else
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	827 {
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	828 if (complex_case)
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	829 retval(3) = CZ;
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	830 else
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	831 retval(3) = ZZ;
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	832 }
23826 d69021d58a61 avoid fallthrough warnings John W. Eaton <jwe@octave.org> parents: 23807 diff changeset	833 OCTAVE_FALLTHROUGH;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	834
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	835 case 3:
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	836 if (nargin == 3)
14876 54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	837 {
54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	838 if (complex_case)
54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	839 retval(2) = CZ;
54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	840 else
54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	841 retval(2) = ZZ;
54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	842 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	843 else
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	844 {
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	845 if (complex_case)
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	846 retval(2) = CQ.hermitian ();
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	847 else
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	848 retval(2) = QQ.transpose ();
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	849 }
23826 d69021d58a61 avoid fallthrough warnings John W. Eaton <jwe@octave.org> parents: 23807 diff changeset	850 OCTAVE_FALLTHROUGH;
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	851
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	852 case 2:
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	853 {
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	854 if (complex_case)
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	855 {
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	856 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	857 octave_stdout << "qz: retval(1) = cbb =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	858 octave_print_internal (octave_stdout, cbb);
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	859 octave_stdout << "\nqz: retval(0) = caa =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	860 octave_print_internal (octave_stdout, caa);
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	861 octave_stdout << std::endl;
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11570 diff changeset	862 #endif
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	863 retval(1) = cbb;
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	864 retval(0) = caa;
a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	865 }
17787 175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	866 else
175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	867 {
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	868 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	869 octave_stdout << "qz: retval(1) = bb =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	870 octave_print_internal (octave_stdout, bb);
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	871 octave_stdout << "\nqz: retval(0) = aa =\n";
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	872 octave_print_internal (octave_stdout, aa);
acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	873 octave_stdout << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	874 #endif
17787 175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	875 retval(1) = bb;
175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	876 retval(0) = aa;
175b392e91fe Use GNU style coding conventions for code in libinterp/ Rik <rik@octave.org> parents: 17744 diff changeset	877 }
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11570 diff changeset	878 }
10311 a217e1d74353 untabify qz.cc John W. Eaton <jwe@octave.org> parents: 10307 diff changeset	879 break;
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	880
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	881 case 1:
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	882 case 0:
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	883 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	884 octave_stdout << "qz: retval(0) = gev = " << gev << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	885 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	886 retval(0) = gev;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	887 break;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	888
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	889 default:
11553 01f703952eff Improve docstrings for functions in DLD-FUNCTIONS directory. Rik <octave@nomad.inbox5.com> parents: 11523 diff changeset	890 error ("qz: too many return arguments");
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	891 break;
25103 078b795c5219 maint: style check C++ ahead of 4.4 release. Rik <rik@octave.org> parents: 25054 diff changeset	892 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	893
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21562 diff changeset	894 #if defined (DEBUG)
26001 acb4689aa5f2 * qz.cc: Use octave_stdout instead of std::cout for debugging messages. John W. Eaton <jwe@octave.org> parents: 25734 diff changeset	895 octave_stdout << "qz: exiting (at long last)" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	896 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	897
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	898 return retval;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	899 }
13074 0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	900
0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	901 /*
0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	902 %!shared a, b, c
14501 60e5cf354d80 Update %!tests in DLD-FUNCTIONS/ directory with Octave coding conventions. Rik <octave@nomad.inbox5.com> parents: 14360 diff changeset	903 %! a = [1 2; 0 3];
60e5cf354d80 Update %!tests in DLD-FUNCTIONS/ directory with Octave coding conventions. Rik <octave@nomad.inbox5.com> parents: 14360 diff changeset	904 %! b = [1 0; 0 0];
60e5cf354d80 Update %!tests in DLD-FUNCTIONS/ directory with Octave coding conventions. Rik <octave@nomad.inbox5.com> parents: 14360 diff changeset	905 %! c = [0 1; 0 0];
60e5cf354d80 Update %!tests in DLD-FUNCTIONS/ directory with Octave coding conventions. Rik <octave@nomad.inbox5.com> parents: 14360 diff changeset	906 %!assert (qz (a,b), 1)
60e5cf354d80 Update %!tests in DLD-FUNCTIONS/ directory with Octave coding conventions. Rik <octave@nomad.inbox5.com> parents: 14360 diff changeset	907 %!assert (isempty (qz (a,c)))
13074 0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	908
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	909 ## Example 7.7.3 in Golub & Van Loan
13074 0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	910 %!test
0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	911 %! a = [ 10 1 2;
0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	912 %! 1 2 -1;
0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	913 %! 1 1 2];
14501 60e5cf354d80 Update %!tests in DLD-FUNCTIONS/ directory with Octave coding conventions. Rik <octave@nomad.inbox5.com> parents: 14360 diff changeset	914 %! b = reshape (1:9,3,3);
13074 0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	915 %! [aa, bb, q, z, v, w, lambda] = qz (a, b);
14501 60e5cf354d80 Update %!tests in DLD-FUNCTIONS/ directory with Octave coding conventions. Rik <octave@nomad.inbox5.com> parents: 14360 diff changeset	916 %! sz = length (lambda);
60e5cf354d80 Update %!tests in DLD-FUNCTIONS/ directory with Octave coding conventions. Rik <octave@nomad.inbox5.com> parents: 14360 diff changeset	917 %! observed = (b * v * diag ([lambda;0])) (:, 1:sz);
18678 6113e0c6920b maint: Clean up extra spaces before/after parentheses. Rik <rik@octave.org> parents: 17787 diff changeset	918 %! assert ((av)(:, 1:sz), observed, norm (observed) 1e-14);
13088 4ffea87ad71b codesprint: Fix tolerance for qz.cc tests. Ben Abbott <bpabbott@mac.com> parents: 13074 diff changeset	919 %! observed = (diag ([lambda;0]) * w' * b) (1:sz, :);
18678 6113e0c6920b maint: Clean up extra spaces before/after parentheses. Rik <rik@octave.org> parents: 17787 diff changeset	920 %! assert ((w'a)(1:sz, :) , observed, norm (observed) 1e-13);
13088 4ffea87ad71b codesprint: Fix tolerance for qz.cc tests. Ben Abbott <bpabbott@mac.com> parents: 13074 diff changeset	921 %! assert (q * a * z, aa, norm (aa) * 1e-14);
4ffea87ad71b codesprint: Fix tolerance for qz.cc tests. Ben Abbott <bpabbott@mac.com> parents: 13074 diff changeset	922 %! assert (q * b * z, bb, norm (bb) * 1e-14);
13074 0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	923
14876 54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	924 %!test
54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	925 %! A = [0, 0, -1, 0; 1, 0, 0, 0; -1, 0, -2, -1; 0, -1, 1, 0];
54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	926 %! B = [0, 0, 0, 0; 0, 1, 0, 0; 0, 0, 1, 0; 0, 0, 0, 1];
54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	927 %! [AA, BB, Q, Z1] = qz (A, B);
23430 32d532b8e7d0 qz.cc: Overhaul qz function. Rik <rik@octave.org> parents: 23429 diff changeset	928 %! [AA, BB, Z2] = qz (A, B, "-");
14876 54a386f2ac4e correcly compute Z for 3-output call to qz (bug #36728) John W. Eaton <jwe@octave.org> parents: 14844 diff changeset	929 %! assert (Z1, Z2);
25570 da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	930
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	931 %!test
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	932 %! A = [ -1.03428 0.24929 0.43205 -0.12860;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	933 %! 1.16228 0.27870 2.12954 0.69250;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	934 %! -0.51524 -0.34939 -0.77820 2.13721;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	935 %! -1.32941 2.11870 0.72005 1.00835 ];
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	936 %! B = [ 1.407302 -0.632956 -0.360628 0.068534;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	937 %! 0.149898 0.298248 0.991777 0.023652;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	938 %! 0.169281 -0.405205 -1.775834 1.511730;
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	939 %! 0.717770 1.291390 -1.766607 -0.531352 ];
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	940 %! [AA, BB, Z, lambda] = qz (A, B, "+");
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	941 %! assert (all (real (lambda(1:3)) >= 0))
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	942 %! assert (real (lambda(4) < 0))
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	943 %! [AA, BB, Z, lambda] = qz (A, B, "-");
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	944 %! assert (real (lambda(1) < 0))
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	945 %! assert (all (real (lambda(2:4)) >= 0))
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	946 %! [AA, BB, Z, lambda] = qz (A, B, "B");
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	947 %! assert (all (abs (lambda(1:3)) >= 1))
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	948 %! assert (abs (lambda(4) < 1))
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	949 %! [AA, BB, Z, lambda] = qz (A, B, "S");
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	950 %! assert (abs (lambda(1) < 1))
da2b38b7a077 Use LAPACK DTGSEN subrouting for ordered-qz computation (bug #53761) Sébastien Villemot <sebastien@debian.org> parents: 25106 diff changeset	951 %! assert (all (abs (lambda(2:4)) >= 1))
13074 0b789a03bde1 codesprint: Add 3 tests for qz.cc Jordi Gutiérrez Hermoso <jordigh@gmail.com> parents: 12642 diff changeset	952 */

Mercurial > octave

annotate libinterp/corefcn/qz.cc @ 27919:1891570abac8