octave: src/DLD-FUNCTIONS/qz.cc annotate

annotate src/DLD-FUNCTIONS/qz.cc @ 10307:4e4270ab70d6

support complex case in qz

author	Jyh-miin Lin <jyhmiinlin@gmail.com>
date	Thu, 11 Feb 2010 08:13:49 +0100
parents	d0ce5e973937
children	a217e1d74353

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
8920 eb63fbe60fab update copyright notices John W. Eaton <jwe@octave.org> parents: 8517 diff changeset	3 Copyright (C) 1998, 1999, 2000, 2002, 2003, 2004, 2005, 2006, 2007,
eb63fbe60fab update copyright notices John W. Eaton <jwe@octave.org> parents: 8517 diff changeset	4 2008, 2009 A. S. Hodel
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	5
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	6 This file is part of Octave.
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	7
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	8 Octave is free software; you can redistribute it and/or modify it
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	9 under the terms of the GNU General Public License as published by the
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5823 diff changeset	10 Free Software Foundation; either version 3 of the License, or (at your
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5823 diff changeset	11 option) any later version.
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	12
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	13 Octave is distributed in the hope that it will be useful, but WITHOUT
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	16 for more details.
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	17
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	18 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	19 along with Octave; see the file COPYING. If not, see
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5823 diff changeset	20 <http://www.gnu.org/licenses/>.
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	21
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	22 */
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	23
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	24 // Generalized eigenvalue balancing via LAPACK
3911 8389e78e67d4 [project @ 2002-04-28 02:15:38 by jwe] jwe parents: 3887 diff changeset	25
8389e78e67d4 [project @ 2002-04-28 02:15:38 by jwe] jwe parents: 3887 diff changeset	26 // Author: A. S. Hodel <scotte@eng.auburn.edu>
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	27
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	28 #undef DEBUG
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	29 #undef DEBUG_SORT
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	30 #undef DEBUG_EIG
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	31
9786 2c279308f6ab fix includes in some src/DLD-FUNCTIONS files John W. Eaton <jwe@octave.org> parents: 9758 diff changeset	32 #ifdef HAVE_CONFIG_H
2c279308f6ab fix includes in some src/DLD-FUNCTIONS files John W. Eaton <jwe@octave.org> parents: 9758 diff changeset	33 #include <config.h>
2c279308f6ab fix includes in some src/DLD-FUNCTIONS files John W. Eaton <jwe@octave.org> parents: 9758 diff changeset	34 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	35
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	36 #include <cfloat>
4051 b79da8779a0e [project @ 2002-08-17 19:38:32 by jwe] jwe parents: 3911 diff changeset	37
3523 b80bbb43a1a9 [project @ 2000-02-02 10:25:52 by jwe] jwe parents: 3500 diff changeset	38 #include <iostream>
4051 b79da8779a0e [project @ 2002-08-17 19:38:32 by jwe] jwe parents: 3911 diff changeset	39 #include <iomanip>
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	40
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	41 #include "CmplxQRP.h"
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	42 #include "CmplxQR.h"
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	43 #include "dbleQR.h"
4153 6b96ce9f5743 [project @ 2002-11-06 20:38:49 by jwe] jwe parents: 4051 diff changeset	44 #include "f77-fcn.h"
7231 2eb392d058bb [project @ 2007-11-30 18:53:29 by jwe] jwe parents: 7017 diff changeset	45 #include "lo-math.h"
4153 6b96ce9f5743 [project @ 2002-11-06 20:38:49 by jwe] jwe parents: 4051 diff changeset	46 #include "quit.h"
6b96ce9f5743 [project @ 2002-11-06 20:38:49 by jwe] jwe parents: 4051 diff changeset	47
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	48 #include "defun-dld.h"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	49 #include "error.h"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	50 #include "gripes.h"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	51 #include "oct-obj.h"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	52 #include "oct-map.h"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	53 #include "ov.h"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	54 #include "pager.h"
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	55 #if defined (DEBUG) \|\| defined (DEBUG_SORT)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	56 #include "pr-output.h"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	57 #endif
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	58 #include "symtab.h"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	59 #include "utils.h"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	60 #include "variables.h"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	61
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	62 typedef octave_idx_type (*sort_function) (const octave_idx_type& LSIZE, const double& ALPHA,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	63 const double& BETA, const double& S,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	64 const double& P);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	65
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	66 extern "C"
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	67 {
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	68 F77_RET_T
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	69 F77_FUNC (dggbal, DGGBAL) (F77_CONST_CHAR_ARG_DECL,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	70 const octave_idx_type& N, double* A, const octave_idx_type& LDA,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	71 double* B, const octave_idx_type& LDB, octave_idx_type& ILO,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	72 octave_idx_type& IHI, double* LSCALE, double* RSCALE,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	73 double* WORK, octave_idx_type& INFO
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	74 F77_CHAR_ARG_LEN_DECL);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	75
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	76 F77_RET_T
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	77 F77_FUNC (zggbal, ZGGBAL) (F77_CONST_CHAR_ARG_DECL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	78 const octave_idx_type& N, Complex* A, const octave_idx_type& LDA,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	79 Complex* B, const octave_idx_type& LDB, octave_idx_type& ILO,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	80 octave_idx_type& IHI, double* LSCALE, double* RSCALE,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	81 double* WORK, octave_idx_type& INFO
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	82 F77_CHAR_ARG_LEN_DECL);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	83
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	84 F77_RET_T
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	85 F77_FUNC (dggbak, DGGBAK) (F77_CONST_CHAR_ARG_DECL,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	86 F77_CONST_CHAR_ARG_DECL,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	87 const octave_idx_type& N, const octave_idx_type& ILO,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	88 const octave_idx_type& IHI, const double* LSCALE,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	89 const double* RSCALE, octave_idx_type& M, double* V,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	90 const octave_idx_type& LDV, octave_idx_type& INFO
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	91 F77_CHAR_ARG_LEN_DECL
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	92 F77_CHAR_ARG_LEN_DECL);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	93
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	94 F77_RET_T
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	95 F77_FUNC (zggbak, ZGGBAK) (F77_CONST_CHAR_ARG_DECL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	96 F77_CONST_CHAR_ARG_DECL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	97 const octave_idx_type& N, const octave_idx_type& ILO,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	98 const octave_idx_type& IHI, const double* LSCALE,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	99 const double* RSCALE, octave_idx_type& M, Complex* V,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	100 const octave_idx_type& LDV, octave_idx_type& INFO
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	101 F77_CHAR_ARG_LEN_DECL
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	102 F77_CHAR_ARG_LEN_DECL);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	103
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	104 F77_RET_T
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	105 F77_FUNC (dgghrd, DGGHRD) (F77_CONST_CHAR_ARG_DECL,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	106 F77_CONST_CHAR_ARG_DECL,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	107 const octave_idx_type& N, const octave_idx_type& ILO,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	108 const octave_idx_type& IHI, double* A,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	109 const octave_idx_type& LDA, double* B,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	110 const octave_idx_type& LDB, double* Q,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	111 const octave_idx_type& LDQ, double* Z,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	112 const octave_idx_type& LDZ, octave_idx_type& INFO
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	113 F77_CHAR_ARG_LEN_DECL
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	114 F77_CHAR_ARG_LEN_DECL);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	115
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	116 F77_RET_T
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	117 F77_FUNC (zgghrd, ZGGHRD) (F77_CONST_CHAR_ARG_DECL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	118 F77_CONST_CHAR_ARG_DECL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	119 const octave_idx_type& N, const octave_idx_type& ILO,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	120 const octave_idx_type& IHI, Complex* A,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	121 const octave_idx_type& LDA, Complex* B,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	122 const octave_idx_type& LDB, Complex* Q,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	123 const octave_idx_type& LDQ, Complex* Z,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	124 const octave_idx_type& LDZ, octave_idx_type& INFO
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	125 F77_CHAR_ARG_LEN_DECL
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	126 F77_CHAR_ARG_LEN_DECL);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	127
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	128 F77_RET_T
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	129 F77_FUNC (dhgeqz, DHGEQZ) (F77_CONST_CHAR_ARG_DECL,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	130 F77_CONST_CHAR_ARG_DECL,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	131 F77_CONST_CHAR_ARG_DECL,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	132 const octave_idx_type& N, const octave_idx_type& ILO, const octave_idx_type& IHI,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	133 double* A, const octave_idx_type& LDA, double* B,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	134 const octave_idx_type& LDB, double* ALPHAR,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	135 double* ALPHAI, double* BETA, double* Q,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	136 const octave_idx_type& LDQ, double* Z,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	137 const octave_idx_type& LDZ, double* WORK,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	138 const octave_idx_type& LWORK, octave_idx_type& INFO
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	139 F77_CHAR_ARG_LEN_DECL
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	140 F77_CHAR_ARG_LEN_DECL
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	141 F77_CHAR_ARG_LEN_DECL);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	142
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	143 F77_RET_T
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	144 F77_FUNC (zhgeqz, ZHGEQZ) (F77_CONST_CHAR_ARG_DECL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	145 F77_CONST_CHAR_ARG_DECL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	146 F77_CONST_CHAR_ARG_DECL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	147 const octave_idx_type& N, const octave_idx_type& ILO, const octave_idx_type& IHI,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	148 Complex* A, const octave_idx_type& LDA, Complex* B,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	149 const octave_idx_type& LDB, Complex* ALPHA, Complex* BETA, Complex* CQ, const octave_idx_type& LDQ,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	150 Complex* CZ, const octave_idx_type& LDZ,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	151 Complex* WORK,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	152 const octave_idx_type& LWORK, double* RWORK,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	153 octave_idx_type& INFO
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	154 F77_CHAR_ARG_LEN_DECL
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	155 F77_CHAR_ARG_LEN_DECL
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	156 F77_CHAR_ARG_LEN_DECL);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	157
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	158 F77_RET_T
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	159 F77_FUNC (dlag2, DLAG2) (const double* A, const octave_idx_type& LDA, const double* B,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	160 const octave_idx_type& LDB, const double& SAFMIN,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	161 double& SCALE1, double& SCALE2,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	162 double& WR1, double& WR2, double& WI);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	163
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	164 // Van Dooren's code (netlib.org: toms/590) for reordering
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	165 // GEP. Only processes Z, not Q.
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	166 F77_RET_T
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	167 F77_FUNC (dsubsp, DSUBSP) (const octave_idx_type& NMAX, const octave_idx_type& N, double* A,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	168 double* B, double* Z, sort_function,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	169 const double& EPS, octave_idx_type& NDIM, octave_idx_type& FAIL,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	170 octave_idx_type* IND);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	171
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	172 // documentation for DTGEVC incorrectly states that VR, VL are
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	173 // complex*16; they are declared in DTGEVC as double precision
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	174 // (probably a cut and paste problem fro ZTGEVC)
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	175 F77_RET_T
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	176 F77_FUNC (dtgevc, DTGEVC) (F77_CONST_CHAR_ARG_DECL,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	177 F77_CONST_CHAR_ARG_DECL,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	178 octave_idx_type* SELECT, const octave_idx_type& N, double* A,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	179 const octave_idx_type& LDA, double* B,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	180 const octave_idx_type& LDB, double* VL,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	181 const octave_idx_type& LDVL, double* VR,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	182 const octave_idx_type& LDVR, const octave_idx_type& MM,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	183 octave_idx_type& M, double* WORK, octave_idx_type& INFO
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	184 F77_CHAR_ARG_LEN_DECL
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	185 F77_CHAR_ARG_LEN_DECL);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	186
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	187 F77_RET_T
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	188 F77_FUNC (ztgevc, ZTGEVC) (F77_CONST_CHAR_ARG_DECL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	189 F77_CONST_CHAR_ARG_DECL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	190 octave_idx_type* SELECT, const octave_idx_type& N,const Complex* A,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	191 const octave_idx_type& LDA,const Complex* B,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	192 const octave_idx_type& LDB, Complex* xVL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	193 const octave_idx_type& LDVL, Complex* xVR,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	194 const octave_idx_type& LDVR, const octave_idx_type& MM,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	195 octave_idx_type& M, Complex* CWORK, double* RWORK, octave_idx_type& INFO
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	196 F77_CHAR_ARG_LEN_DECL
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	197 F77_CHAR_ARG_LEN_DECL);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	198
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	199 F77_RET_T
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	200 F77_FUNC (xdlamch, XDLAMCH) (F77_CONST_CHAR_ARG_DECL,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	201 double& retval
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	202 F77_CHAR_ARG_LEN_DECL);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	203
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	204 F77_RET_T
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	205 F77_FUNC (xdlange, XDLANGE) (F77_CONST_CHAR_ARG_DECL,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	206 const octave_idx_type&, const octave_idx_type&, const double*,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	207 const octave_idx_type&, double*, double&
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	208 F77_CHAR_ARG_LEN_DECL);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	209 }
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	210
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	211 // fcrhp, fin, fout, folhp:
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	212 // routines for ordering of generalized eigenvalues
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	213 // return 1 if test is passed, 0 otherwise
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	214 // fin: \|lambda\| < 1
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	215 // fout: \|lambda\| >= 1
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	216 // fcrhp: real(lambda) >= 0
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	217 // folhp: real(lambda) < 0
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	218
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	219 static octave_idx_type
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	220 fcrhp (const octave_idx_type& lsize, const double& alpha,
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	221 const double& beta, const double& s, const double&)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	222 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	223 if (lsize == 1)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	224 return (alpha*beta >= 0 ? 1 : -1);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	225 else
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	226 return (s >= 0 ? 1 : -1);
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	227 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	228
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	229 static octave_idx_type
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	230 fin (const octave_idx_type& lsize, const double& alpha,
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	231 const double& beta, const double&, const double& p)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	232 {
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	233 octave_idx_type retval;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	234
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	235 if (lsize == 1)
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	236 retval = (fabs (alpha) < fabs (beta) ? 1 : -1);
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	237 else
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	238 retval = (fabs (p) < 1 ? 1 : -1);
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	239
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	240 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	241 std::cout << "qz: fin: retval=" << retval << 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
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	244 return retval;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	245 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	246
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	247 static octave_idx_type
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	248 folhp (const octave_idx_type& lsize, const double& alpha,
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	249 const double& beta, const double& s, const double&)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	250 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	251 if (lsize == 1)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	252 return (alpha*beta < 0 ? 1 : -1);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	253 else
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	254 return (s < 0 ? 1 : -1);
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	255 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	256
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	257 static octave_idx_type
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	258 fout (const octave_idx_type& lsize, const double& alpha,
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	259 const double& beta, const double&, const double& p)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	260 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	261 if (lsize == 1)
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	262 return (fabs (alpha) >= fabs (beta) ? 1 : -1);
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	263 else
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	264 return (fabs (p) >= 1 ? 1 : -1);
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	265 }
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	266
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	267 DEFUN_DLD (qz, args, nargout,
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	268 "-- texinfo --\n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	269 @deftypefn {Loadable Function} {@var{lambda} =} qz (@var{a}, @var{b})\n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	270 Generalized eigenvalue problem @math{A x = s B x},\n\
9064 7c02ec148a3c Check grammar on all .cc files Rik <rdrider0-list@yahoo.com> parents: 8927 diff changeset	271 @var{QZ} decomposition. There are three ways to call this function:\n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	272 @enumerate\n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	273 @item @code{lambda = qz(A,B)}\n\
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	274 \n\
5016 bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	275 Computes the generalized eigenvalues\n\
bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	276 @tex\n\
bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	277 $\\lambda$\n\
bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	278 @end tex\n\
8517 81d6ab3ac93c Allow documentation tobe built for other formats than tex and info sh@sh-laptop parents: 8286 diff changeset	279 @ifnottex\n\
5016 bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	280 @var{lambda}\n\
8517 81d6ab3ac93c Allow documentation tobe built for other formats than tex and info sh@sh-laptop parents: 8286 diff changeset	281 @end ifnottex\n\
5016 bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	282 of @math{(A - s B)}.\n\
3500 7923abdeb4e5 [project @ 2000-01-31 06:35:00 by jwe] jwe parents: 3427 diff changeset	283 @item @code{[AA, BB, Q, Z, V, W, lambda] = qz (A, B)}\n\
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	284 \n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	285 Computes qz decomposition, generalized eigenvectors, and \n\
5481 65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	286 generalized eigenvalues of @math{(A - sB)}\n\
5016 bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	287 @tex\n\
bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	288 $$ AV = BV{ \\rm diag }(\\lambda) $$\n\
bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	289 $$ W^T A = { \\rm diag }(\\lambda)W^T B $$\n\
bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	290 $$ AA = Q^T AZ, BB = Q^T BZ $$\n\
bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	291 @end tex\n\
8517 81d6ab3ac93c Allow documentation tobe built for other formats than tex and info sh@sh-laptop parents: 8286 diff changeset	292 @ifnottex\n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	293 @example\n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	294 @group\n\
5481 65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	295 \n\
65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	296 AV = BV*diag(lambda)\n\
65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	297 W'A = diag(lambda)W'*B\n\
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	298 AA = QAZ, BB = QBZ\n\
5481 65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	299 \n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	300 @end group\n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	301 @end example\n\
8517 81d6ab3ac93c Allow documentation tobe built for other formats than tex and info sh@sh-laptop parents: 8286 diff changeset	302 @end ifnottex\n\
5016 bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	303 with @var{Q} and @var{Z} orthogonal (unitary)= @var{I}\n\
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	304 \n\
5016 bdbee5282954 [project @ 2004-09-22 02:50:35 by jwe] jwe parents: 4629 diff changeset	305 @item @code{[AA,BB,Z@{, lambda@}] = qz(A,B,opt)}\n\
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	306 \n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	307 As in form [2], but allows ordering of generalized eigenpairs\n\
5481 65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	308 for (e.g.) solution of discrete time algebraic Riccati equations.\n\
65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	309 Form 3 is not available for complex matrices, and does not compute\n\
65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	310 the generalized eigenvectors @var{V}, @var{W}, nor the orthogonal matrix @var{Q}.\n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	311 @table @var\n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	312 @item opt\n\
9064 7c02ec148a3c Check grammar on all .cc files Rik <rdrider0-list@yahoo.com> parents: 8927 diff changeset	313 for ordering eigenvalues of the GEP pencil. The leading block\n\
5481 65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	314 of the revised pencil contains all eigenvalues that satisfy:\n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	315 @table @code\n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	316 @item \"N\"\n\
5481 65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	317 = unordered (default) \n\
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	318 \n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	319 @item \"S\"\n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	320 = small: leading block has all \|lambda\| <=1 \n\
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	321 \n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	322 @item \"B\"\n\
5481 65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	323 = big: leading block has all \|lambda\| >= 1 \n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	324 \n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	325 @item \"-\"\n\
5481 65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	326 = negative real part: leading block has all eigenvalues\n\
65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	327 in the open left half-plane\n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	328 \n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	329 @item \"+\"\n\
9209 923c7cb7f13f Simplify TeXinfo files by eliminating redundant @iftex followed by @tex construction. Rik <rdrider0-list@yahoo.com> parents: 9064 diff changeset	330 = non-negative real part: leading block has all eigenvalues\n\
5481 65e02518f4f5 [project @ 2005-10-04 17:12:08 by jwe] jwe parents: 5333 diff changeset	331 in the closed right half-plane\n\
9064 7c02ec148a3c Check grammar on all .cc files Rik <rdrider0-list@yahoo.com> parents: 8927 diff changeset	332 @end table\n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	333 @end table\n\
f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	334 @end enumerate\n\
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	335 \n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	336 Note: qz performs permutation balancing, but not scaling (see balance).\n\
9064 7c02ec148a3c Check grammar on all .cc files Rik <rdrider0-list@yahoo.com> parents: 8927 diff changeset	337 Order of output arguments was selected for compatibility with @sc{matlab}\n\
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	338 \n\
8286 6f2d95255911 fix @seealso references to point to existing anchors Thorsten Meyer <thorsten.meyier@gmx.de> parents: 7520 diff changeset	339 @seealso{balance, eig, schur}\n\
3372 f16c2ce14886 [project @ 1999-11-23 19:07:09 by jwe] jwe parents: 3187 diff changeset	340 @end deftypefn")
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	341 {
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	342 octave_value_list retval;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	343 int nargin = args.length ();
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	344
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	345 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	346 std::cout << "qz: nargin = " << nargin << ", nargout = " << nargout << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	347 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	348
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	349 if (nargin < 2 \|\| nargin > 3 \|\| nargout > 7)
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	350 {
5823 080c08b192d8 [project @ 2006-05-19 05:32:17 by jwe] jwe parents: 5775 diff changeset	351 print_usage ();
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	352 return retval;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	353 }
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	354 else if (nargin == 3 && (nargout < 3 \|\| nargout > 4))
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	355 {
3427 e098ebb77023 [project @ 2000-01-13 09:25:53 by jwe] jwe parents: 3372 diff changeset	356 error ("qz: invalid number of output arguments for form [3] call");
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	357 return retval;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	358 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	359
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	360 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	361 std::cout << "qz: determine ordering option" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	362 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	363
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	364 // Determine ordering option
8927 d75f4ee0538d qz.cc (Fqz): avoid maybe clobbred by vfork warning from GCC John W. Eaton <jwe@octave.org> parents: 8920 diff changeset	365 volatile char ord_job = 0;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	366 static double safmin;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	367
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	368 if (nargin == 2)
8927 d75f4ee0538d qz.cc (Fqz): avoid maybe clobbred by vfork warning from GCC John W. Eaton <jwe@octave.org> parents: 8920 diff changeset	369 ord_job = 'N';
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	370 else if (!args(2).is_string ())
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	371 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	372 error ("qz: argument 3 must be a string");
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	373 return retval;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	374 }
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	375 else
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	376 {
8927 d75f4ee0538d qz.cc (Fqz): avoid maybe clobbred by vfork warning from GCC John W. Eaton <jwe@octave.org> parents: 8920 diff changeset	377 std::string tmp = args(2).string_value ();
d75f4ee0538d qz.cc (Fqz): avoid maybe clobbred by vfork warning from GCC John W. Eaton <jwe@octave.org> parents: 8920 diff changeset	378
d75f4ee0538d qz.cc (Fqz): avoid maybe clobbred by vfork warning from GCC John W. Eaton <jwe@octave.org> parents: 8920 diff changeset	379 if (! tmp.empty ())
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	380 ord_job = tmp[0];
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	381
8927 d75f4ee0538d qz.cc (Fqz): avoid maybe clobbred by vfork warning from GCC John W. Eaton <jwe@octave.org> parents: 8920 diff changeset	382 if (! (ord_job == 'N' \|\| ord_job == 'n'
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	383 \|\| ord_job == 'S' \|\| ord_job == 's'
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	384 \|\| ord_job == 'B' \|\| ord_job == 'b'
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	385 \|\| ord_job == '+' \|\| ord_job == '-'))
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	386 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	387 error ("qz: invalid order option");
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	388 return retval;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	389 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	390
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	391 // overflow constant required by dlag2
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	392 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	393 safmin
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	394 F77_CHAR_ARG_LEN (1));
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	395
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	396 #ifdef DEBUG_EIG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	397 std::cout << "qz: initial value of safmin=" << setiosflags (std::ios::scientific)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	398 << safmin << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	399 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	400
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	401 // some machines (e.g., DEC alpha) get safmin = 0;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	402 // for these, use eps instead to avoid problems in dlag2
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	403 if (safmin == 0)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	404 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	405 #ifdef DEBUG_EIG
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	406 std::cout << "qz: DANGER WILL ROBINSON: safmin is 0!" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	407 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	408
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	409 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	410 safmin
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	411 F77_CHAR_ARG_LEN (1));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	412
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	413 #ifdef DEBUG_EIG
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	414 std::cout << "qz: safmin set to " << setiosflags (std::ios::scientific)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	415 << safmin << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	416 #endif
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	417 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	418 }
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	419
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	420 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	421 std::cout << "qz: check argument 1" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	422 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	423
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	424 // Argument 1: check if it's o.k. dimensioned
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	425 octave_idx_type nn = args(0).rows ();
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	426
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	427 #ifdef DEBUG
3531 97cf542676e1 [project @ 2000-02-02 11:30:40 by jwe] jwe parents: 3523 diff changeset	428 std::cout << "argument 1 dimensions: (" << nn << "," << args(0).columns () << ")"
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	429 << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	430 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	431
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	432 int arg_is_empty = empty_arg ("qz", nn, args(0).columns ());
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	433
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	434 if (arg_is_empty < 0)
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	435 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	436 gripe_empty_arg ("qz: parameter 1", 0);
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	437 return retval;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	438 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	439 else if (arg_is_empty > 0)
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	440 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	441 gripe_empty_arg ("qz: parameter 1; continuing", 0);
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	442 return octave_value_list (2, Matrix ());
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	443 }
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	444 else if (args(0).columns () != nn)
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	445 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	446 gripe_square_matrix_required ("qz");
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	447 return retval;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	448 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	449
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	450 // Argument 1: dimensions look good; get the value
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	451 Matrix aa;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	452 ComplexMatrix caa;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	453
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	454 if (args(0).is_complex_type ())
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	455 caa = args(0).complex_matrix_value ();
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	456 else
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	457 aa = args(0).matrix_value ();
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	458
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	459 if (error_state)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	460 return retval;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	461
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	462 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	463 std::cout << "qz: check argument 2" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	464 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	465
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	466 // Extract argument 2 (bb, or cbb if complex)
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	467 if ((nn != args(1).columns ()) \|\| (nn != args(1).rows ()))
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	468 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	469 gripe_nonconformant ();
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	470 return retval;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	471 }
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	472
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	473 Matrix bb;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	474 ComplexMatrix cbb;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	475
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	476 if (args(1).is_complex_type ())
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	477 cbb = args(1).complex_matrix_value ();
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	478 else
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	479 bb = args(1).matrix_value ();
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	480
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	481 if (error_state)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	482 return retval;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	483
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	484 // Both matrices loaded, now let's check what kind of arithmetic:
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	485 //declared static to avoid compiler warnings about long jumps, vforks.
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	486
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	487 int complex_case
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	488 = (args(0).is_complex_type () \|\| args(1).is_complex_type ());
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	489 // static complex_case causing random segfault, so it is removed
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	490 if (nargin == 3 && complex_case)
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	491 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	492 error ("qz: cannot re-order complex qz decomposition.");
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	493 return retval;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	494 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	495
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	496 // first, declare variables used in both the real and complex case
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	497 Matrix QQ(nn,nn), ZZ(nn,nn), VR(nn,nn), VL(nn,nn);
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	498 RowVector alphar(nn), alphai(nn), betar(nn);
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	499 ComplexRowVector xalpha(nn), xbeta(nn);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	500 ComplexMatrix CQ(nn,nn), CZ(nn,nn), CVR(nn,nn), CVL(nn,nn);
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	501 octave_idx_type ilo, ihi, info;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	502 char compq = (nargout >= 3 ? 'V' : 'N');
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	503 char compz = (nargout >= 4 ? 'V' : 'N');
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	504
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	505 // initialize Q, Z to identity if we need either of them
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	506 if (compq == 'V' \|\| compz == 'V')
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	507 for (octave_idx_type ii = 0; ii < nn; ii++)
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5016 diff changeset	508 for (octave_idx_type jj = 0; jj < nn; jj++)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	509 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	510 OCTAVE_QUIT;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	511 QQ(ii,jj) = ZZ(ii,jj) = (ii == jj ? 1.0 : 0.0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	512 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	513
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	514 // always perform permutation balancing
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	515 const char bal_job = 'P';
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	516 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	517
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	518 if (complex_case)
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	519 {
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	520 #ifdef DEBUG
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	521 if (compq == 'V')
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	522 std::cout << "qz: performing balancing; CQ=" << std::endl << CQ << std::endl;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	523 #endif
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	524 if (args(0).is_real_type ())
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	525 caa = ComplexMatrix (aa);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	526
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	527 if (args(1).is_real_type ())
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	528 cbb = ComplexMatrix (bb);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	529
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	530 if (compq == 'V')
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	531 CQ = ComplexMatrix (QQ);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	532
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	533 if (compz == 'V')
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	534 CZ = ComplexMatrix (ZZ);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	535
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	536 F77_XFCN (zggbal, ZGGBAL,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	537 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	538 nn, caa.fortran_vec (), nn, cbb.fortran_vec (),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	539 nn, ilo, ihi, lscale.fortran_vec (),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	540 rscale.fortran_vec (), work.fortran_vec (), info
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	541 F77_CHAR_ARG_LEN (1)));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	542 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	543 else
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	544 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	545 #ifdef DEBUG
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	546 if (compq == 'V')
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	547 std::cout << "qz: performing balancing; QQ=" << std::endl << QQ << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	548 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	549
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	550 F77_XFCN (dggbal, DGGBAL,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	551 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	552 nn, aa.fortran_vec (), nn, bb.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	553 nn, ilo, ihi, lscale.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	554 rscale.fortran_vec (), work.fortran_vec (), info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	555 F77_CHAR_ARG_LEN (1)));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	556 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	557
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	558 // Since we just want the balancing matrices, we can use dggbal
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	559 // for both the real and complex cases;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	560 // left first
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	561
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	562 /* if (compq == 'V')
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	563 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	564 F77_XFCN (dggbak, DGGBAK,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	565 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	566 F77_CONST_CHAR_ARG2 ("L", 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	567 nn, ilo, ihi, lscale.data (), rscale.data (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	568 nn, QQ.fortran_vec (), nn, info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	569 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	570 F77_CHAR_ARG_LEN (1)));
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	571
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	572 #ifdef DEBUG
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	573 if (compq == 'V')
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	574 std::cout << "qz: balancing done; QQ=" << std::endl << QQ << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	575 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	576 }
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	577
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	578 // then right
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	579 if (compz == 'V')
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	580 {
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4153 diff changeset	581 F77_XFCN (dggbak, DGGBAK,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	582 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	583 F77_CONST_CHAR_ARG2 ("R", 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	584 nn, ilo, ihi, lscale.data (), rscale.data (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	585 nn, ZZ.fortran_vec (), nn, info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	586 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	587 F77_CHAR_ARG_LEN (1)));
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	588
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	589 #ifdef DEBUG
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	590 if (compz == 'V')
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	591 std::cout << "qz: balancing done; ZZ=" << std::endl << ZZ << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	592 #endif
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	593 } */
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	594
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	595 static char qz_job;
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	596 qz_job = (nargout < 2 ? 'E' : 'S');
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	597
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	598 if (complex_case)
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	599 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	600 // complex case
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	601 ComplexQR cbqr (cbb); // declare cbqr as the QR decomposition of cbb
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	602 cbb = cbqr.R (); // the R matrix of QR decomposition for cbb
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	603 caa = (cbqr.Q ().hermitian ())caa; // (Q)caa for following work
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	604 //if (compq == 'V')
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	605 CQ = CQ*cbqr.Q ();
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	606 F77_XFCN (zgghrd, ZGGHRD,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	607 (F77_CONST_CHAR_ARG2 (&compq, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	608 F77_CONST_CHAR_ARG2 (&compz, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	609 nn, ilo, ihi, caa.fortran_vec (),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	610 nn, cbb.fortran_vec (), nn, CQ.fortran_vec (), nn,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	611 CZ.fortran_vec (), nn, info
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	612 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	613 F77_CHAR_ARG_LEN (1)));
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	614 ComplexRowVector cwork(1*nn);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	615 F77_XFCN (zhgeqz, ZHGEQZ,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	616 (F77_CONST_CHAR_ARG2 (&qz_job, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	617 F77_CONST_CHAR_ARG2 (&compq, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	618 F77_CONST_CHAR_ARG2 (&compz, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	619 nn, ilo, ihi,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	620 caa.fortran_vec (), nn,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	621 cbb.fortran_vec (),nn,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	622 xalpha.fortran_vec (), xbeta.fortran_vec (),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	623 CQ.fortran_vec (), nn,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	624 CZ.fortran_vec (), nn,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	625 cwork.fortran_vec (), nn, rwork.fortran_vec (), info
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	626 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	627 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	628 F77_CHAR_ARG_LEN (1)));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	629
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	630 if (compq == 'V')
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	631 {
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	632 // Left eigenvector
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	633 F77_XFCN (zggbak, ZGGBAK,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	634 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	635 F77_CONST_CHAR_ARG2 ("L", 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	636 nn, ilo, ihi, lscale.data (), rscale.data (),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	637 nn, CQ.fortran_vec (), nn, info
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	638 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	639 F77_CHAR_ARG_LEN (1)));
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	640 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	641
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	642 // then right
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	643 if (compz == 'V')
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	644 {
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	645 F77_XFCN (zggbak, ZGGBAK,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	646 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	647 F77_CONST_CHAR_ARG2 ("R", 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	648 nn, ilo, ihi, lscale.data (), rscale.data (),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	649 nn, CZ.fortran_vec (), nn, info
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	650 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	651 F77_CHAR_ARG_LEN (1)));
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	652 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	653
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	654 }
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	655 else // real matrices case
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	656 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	657 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	658 std::cout << "qz: peforming qr decomposition of bb" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	659 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	660
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	661 // compute the QR factorization of bb
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	662 QR bqr (bb);
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	663
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	664 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	665 std::cout << "qz: qr (bb) done; now peforming qz decomposition" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	666 #endif
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 bb = bqr.R ();
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	669
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	670 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	671 std::cout << "qz: extracted bb" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	672 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	673
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	674 aa = (bqr.Q ()).transpose ()*aa;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	675
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	676 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	677 std::cout << "qz: updated aa " << std::endl;
0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	678 std::cout << "bqr.Q () = " << std::endl << bqr.Q () << std::endl;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	679
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	680 if (compq == 'V')
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	681 std::cout << "QQ =" << QQ << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	682 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	683
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	684 if (compq == 'V')
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	685 QQ = QQ*bqr.Q ();
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	686
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	687 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	688 std::cout << "qz: precursors done..." << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	689 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	690
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	691 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	692 std::cout << "qz: compq = " << compq << ", compz = " << compz << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	693 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	694
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	695 // reduce to generalized hessenberg form
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	696 F77_XFCN (dgghrd, DGGHRD,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	697 (F77_CONST_CHAR_ARG2 (&compq, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	698 F77_CONST_CHAR_ARG2 (&compz, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	699 nn, ilo, ihi, aa.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	700 nn, bb.fortran_vec (), nn, QQ.fortran_vec (), nn,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	701 ZZ.fortran_vec (), nn, info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	702 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	703 F77_CHAR_ARG_LEN (1)));
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	704
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	705 // check if just computing generalized eigenvalues or if we're
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	706 // actually computing the decomposition
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	707
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	708 // reduce to generalized Schur form
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	709 F77_XFCN (dhgeqz, DHGEQZ,
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	710 (F77_CONST_CHAR_ARG2 (&qz_job, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	711 F77_CONST_CHAR_ARG2 (&compq, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	712 F77_CONST_CHAR_ARG2 (&compz, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	713 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	714 nn, alphar.fortran_vec (), alphai.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	715 betar.fortran_vec (), QQ.fortran_vec (), nn,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	716 ZZ.fortran_vec (), nn, work.fortran_vec (), nn, info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	717 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	718 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	719 F77_CHAR_ARG_LEN (1)));
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	720 if (compq == 'V')
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	721 {
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	722 F77_XFCN (dggbak, DGGBAK,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	723 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	724 F77_CONST_CHAR_ARG2 ("L", 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	725 nn, ilo, ihi, lscale.data (), rscale.data (),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	726 nn, QQ.fortran_vec (), nn, info
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	727 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	728 F77_CHAR_ARG_LEN (1)));
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	729
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	730 #ifdef DEBUG
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	731 if (compq == 'V')
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	732 std::cout << "qz: balancing done; QQ=" << std::endl << QQ << std::endl;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	733 #endif
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	734 }
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	735
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	736 // then right
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	737 if (compz == 'V')
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	738 {
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	739 F77_XFCN (dggbak, DGGBAK,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	740 (F77_CONST_CHAR_ARG2 (&bal_job, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	741 F77_CONST_CHAR_ARG2 ("R", 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	742 nn, ilo, ihi, lscale.data (), rscale.data (),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	743 nn, ZZ.fortran_vec (), nn, info
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	744 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	745 F77_CHAR_ARG_LEN (1)));
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	746
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	747 #ifdef DEBUG
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	748 if (compz == 'V')
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	749 std::cout << "qz: balancing done; ZZ=" << std::endl << ZZ << std::endl;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	750 #endif
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	751 }
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	752
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	753 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	754
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	755 // order the QZ decomposition?
8927 d75f4ee0538d qz.cc (Fqz): avoid maybe clobbred by vfork warning from GCC John W. Eaton <jwe@octave.org> parents: 8920 diff changeset	756 if (! (ord_job == 'N' \|\| ord_job == 'n'))
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	757 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	758 if (complex_case)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	759 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	760 // probably not needed, but better be safe
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	761 error ("qz: cannot re-order complex qz decomposition.");
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	762 return retval;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	763 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	764 else
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	765 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	766 #ifdef DEBUG_SORT
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	767 std::cout << "qz: ordering eigenvalues: ord_job = "
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	768 << ord_job << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	769 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	770
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	771 // declared static to avoid vfork/long jump compiler complaints
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	772 static sort_function sort_test;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	773 sort_test = 0;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	774
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	775 switch (ord_job)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	776 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	777 case 'S':
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	778 case 's':
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	779 sort_test = &fin;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	780 break;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	781
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	782 case 'B':
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	783 case 'b':
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	784 sort_test = &fout;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	785 break;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	786
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	787 case '+':
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	788 sort_test = &fcrhp;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	789 break;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	790
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	791 case '-':
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	792 sort_test = &folhp;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	793 break;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	794
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	795 default:
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	796 // invalid order option (should never happen, since we
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	797 // checked the options at the top).
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	798 panic_impossible ();
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	799 break;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	800 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	801
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	802 octave_idx_type ndim, fail;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	803 double inf_norm;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	804
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	805 F77_XFCN (xdlange, XDLANGE,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	806 (F77_CONST_CHAR_ARG2 ("I", 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	807 nn, nn, aa.data (), nn, work.fortran_vec (), inf_norm
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	808 F77_CHAR_ARG_LEN (1)));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	809
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	810 double eps = DBL_EPSILONinf_normnn;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	811
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	812 #ifdef DEBUG_SORT
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	813 std::cout << "qz: calling dsubsp: aa=" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	814 octave_print_internal (std::cout, aa, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	815 std::cout << std::endl << "bb=" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	816 octave_print_internal (std::cout, bb, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	817 if (compz == 'V')
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	818 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	819 std::cout << std::endl << "ZZ=" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	820 octave_print_internal (std::cout, ZZ, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	821 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	822 std::cout << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	823 std::cout << "alphar = " << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	824 octave_print_internal (std::cout, (Matrix) alphar, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	825 std::cout << std::endl << "alphai = " << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	826 octave_print_internal (std::cout, (Matrix) alphai, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	827 std::cout << std::endl << "beta = " << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	828 octave_print_internal (std::cout, (Matrix) betar, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	829 std::cout << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	830 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	831
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	832 Array<octave_idx_type> ind (nn);
3550 bc492f4a94cb [project @ 2000-02-03 01:36:30 by jwe] jwe parents: 3538 diff changeset	833
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	834 F77_XFCN (dsubsp, DSUBSP,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	835 (nn, nn, aa.fortran_vec (), bb.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	836 ZZ.fortran_vec (), sort_test, eps, ndim, fail,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	837 ind.fortran_vec ()));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	838
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	839 #ifdef DEBUG
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	840 std::cout << "qz: back from dsubsp: aa=" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	841 octave_print_internal (std::cout, aa, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	842 std::cout << std::endl << "bb=" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	843 octave_print_internal (std::cout, bb, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	844 if (compz == 'V')
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	845 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	846 std::cout << std::endl << "ZZ=" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	847 octave_print_internal (std::cout, ZZ, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	848 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	849 std::cout << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	850 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	851
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	852 // manually update alphar, alphai, betar
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	853 static int jj;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	854
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	855 jj=0;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	856 while (jj < nn)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	857 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	858 #ifdef DEBUG_EIG
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	859 std::cout << "computing gen eig #" << jj << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	860 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	861
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	862 static int zcnt; // number of zeros in this block
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	863
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	864 if (jj == (nn-1))
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	865 zcnt = 1;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	866 else if (aa(jj+1,jj) == 0)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	867 zcnt = 1;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	868 else zcnt = 2;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	869
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	870 if (zcnt == 1) // real zero
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	871 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	872 #ifdef DEBUG_EIG
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	873 std::cout << " single gen eig:" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	874 std::cout << " alphar(" << jj << ") = " << aa(jj,jj) << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	875 std::cout << " betar( " << jj << ") = " << bb(jj,jj) << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	876 std::cout << " alphai(" << jj << ") = 0" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	877 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	878
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	879 alphar(jj) = aa(jj,jj);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	880 alphai(jj) = 0;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	881 betar(jj) = bb(jj,jj);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	882 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	883 else
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	884 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	885 // complex conjugate pair
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	886 #ifdef DEBUG_EIG
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	887 std::cout << "qz: calling dlag2:" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	888 std::cout << "safmin="
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	889 << setiosflags (std::ios::scientific) << safmin << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	890
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	891 for (int idr = jj; idr <= jj+1; idr++)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	892 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	893 for (int idc = jj; idc <= jj+1; idc++)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	894 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	895 std::cout << "aa(" << idr << "," << idc << ")="
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	896 << aa(idr,idc) << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	897 std::cout << "bb(" << idr << "," << idc << ")="
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	898 << bb(idr,idc) << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	899 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	900 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	901 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	902
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	903 // FIXME -- probably should be using
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	904 // fortran_vec instead of &aa(jj,jj) here.
4566 30ba814d6700 [project @ 2003-10-29 15:06:29 by jwe] jwe parents: 4552 diff changeset	905
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	906 double scale1, scale2, wr1, wr2, wi;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	907 const double aa_ptr = aa.data () + jjnn+jj;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	908 const double bb_ptr = bb.data () + jjnn+jj;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	909 F77_XFCN (dlag2, DLAG2,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	910 (aa_ptr, nn, bb_ptr, nn, safmin,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	911 scale1, scale2, wr1, wr2, wi));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	912
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	913 #ifdef DEBUG_EIG
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	914 std::cout << "dlag2 returns: scale1=" << scale1
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	915 << "\tscale2=" << scale2 << std::endl
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	916 << "\twr1=" << wr1 << "\twr2=" << wr2
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	917 << "\twi=" << wi << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	918 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	919
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	920 // just to be safe, check if it's a real pair
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	921 if (wi == 0)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	922 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	923 alphar(jj) = wr1;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	924 alphai(jj) = 0;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	925 betar(jj) = scale1;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	926 alphar(jj+1) = wr2;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	927 alphai(jj+1) = 0;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	928 betar(jj+1) = scale2;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	929 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	930 else
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	931 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	932 alphar(jj) = alphar(jj+1)=wr1;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	933 alphai(jj) = -(alphai(jj+1) = wi);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	934 betar(jj) = betar(jj+1) = scale1;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	935 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	936 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	937
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	938 // advance past this block
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	939 jj += zcnt;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	940 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	941
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	942 #ifdef DEBUG_SORT
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	943 std::cout << "qz: back from dsubsp: aa=" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	944 octave_print_internal (std::cout, aa, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	945 std::cout << std::endl << "bb=" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	946 octave_print_internal (std::cout, bb, 0);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	947
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	948 if (compz == 'V')
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	949 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	950 std::cout << std::endl << "ZZ=" << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	951 octave_print_internal (std::cout, ZZ, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	952 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	953 std::cout << std::endl << "qz: ndim=" << ndim << std::endl
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	954 << "fail=" << fail << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	955 std::cout << "alphar = " << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	956 octave_print_internal (std::cout, (Matrix) alphar, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	957 std::cout << std::endl << "alphai = " << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	958 octave_print_internal (std::cout, (Matrix) alphai, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	959 std::cout << std::endl << "beta = " << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	960 octave_print_internal (std::cout, (Matrix) betar, 0);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	961 std::cout << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	962 #endif
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	963 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	964 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	965
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	966 // compute generalized eigenvalues?
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	967 ComplexColumnVector gev;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	968
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	969 if (nargout < 2 \|\| nargout == 7 \|\| (nargin == 3 && nargout == 4))
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	970 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	971 if (complex_case)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	972 {
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	973 int cnt = 0;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	974
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	975 for (int ii = 0; ii < nn; ii++)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	976 // if (cbetar(ii) != 0)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	977 cnt++;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	978
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	979 ComplexColumnVector tmp(cnt);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	980
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	981 cnt = 0;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	982 for (int ii = 0; ii < nn; ii++)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	983 // if (cbetar(ii) != 0)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	984 tmp(cnt++) = xalpha(ii)/xbeta(ii);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	985 gev = tmp;
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	986 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	987 else
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	988 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	989 #ifdef DEBUG
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	990 std::cout << "qz: computing generalized eigenvalues" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	991 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	992
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	993 // return finite generalized eigenvalues
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	994 int cnt = 0;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	995
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	996 for (int ii = 0; ii < nn; ii++)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	997 if (betar(ii) != 0)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	998 cnt++;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	999
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1000 ComplexColumnVector tmp(cnt);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1001
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1002 cnt = 0;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1003 for (int ii = 0; ii < nn; ii++)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1004 if (betar(ii) != 0)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1005 tmp(cnt++) = Complex(alphar(ii), alphai(ii))/betar(ii);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1006 gev = tmp;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1007 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1008 }
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1009
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1010 // right, left eigenvector matrices
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1011 if (nargout >= 5)
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1012 {
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1013 char side = (nargout == 5 ? 'R' : 'B'); // which side to compute?
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1014 char howmny = 'B'; // compute all of them and backtransform
7520 b166043585a8 use 0 instead of NULL John W. Eaton <jwe@octave.org> parents: 7482 diff changeset	1015 octave_idx_type *select = 0; // dummy pointer; select is not used.
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1016
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1017 if (complex_case)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1018 {
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1019 octave_idx_type m;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1020 CVL=CQ;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1021 CVR=CZ;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1022 ComplexRowVector cwork2(2*nn);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1023 RowVector rwork2(8*nn);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1024
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1025 //octave_idx_type n=nn;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1026 F77_XFCN (ztgevc, ZTGEVC,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1027 (F77_CONST_CHAR_ARG2 (&side, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1028 F77_CONST_CHAR_ARG2 (&howmny, 1),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1029 select, nn, caa.fortran_vec (), nn, cbb.fortran_vec (),
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1030 nn, CVL.fortran_vec (), nn, CVR.fortran_vec (), nn, nn,
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1031 m, cwork2.fortran_vec (), rwork2.fortran_vec (), info
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1032 F77_CHAR_ARG_LEN (1)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1033 F77_CHAR_ARG_LEN (1)));
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1034 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1035 else
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1036 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1037 #ifdef DEBUG
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1038 std::cout << "qz: computing generalized eigenvectors" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1039 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1040
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1041 VL = QQ;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1042 VR = ZZ;
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1043 octave_idx_type m;
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1044 F77_XFCN (dtgevc, DTGEVC,
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1045 (F77_CONST_CHAR_ARG2 (&side, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1046 F77_CONST_CHAR_ARG2 (&howmny, 1),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1047 select, nn, aa.fortran_vec (), nn, bb.fortran_vec (),
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1048 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	1049 m, work.fortran_vec (), info
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1050 F77_CHAR_ARG_LEN (1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1051 F77_CHAR_ARG_LEN (1)));
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1052
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1053 // now construct the complex form of VV, WW
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1054 int jj = 0;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1055
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1056 while (jj < nn)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1057 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1058 OCTAVE_QUIT;
4153 6b96ce9f5743 [project @ 2002-11-06 20:38:49 by jwe] jwe parents: 4051 diff changeset	1059
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1060 // see if real or complex eigenvalue
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1061 int cinc = 2; // column increment; assume complex eigenvalue
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1062
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1063 if (jj == (nn-1))
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1064 cinc = 1; // single column
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1065 else if (aa(jj+1,jj) == 0)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1066 cinc = 1;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1067
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1068 // now copy the eigenvector (s) to CVR, CVL
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1069 if (cinc == 1)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1070 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1071 for (int ii = 0; ii < nn; ii++)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1072 CVR(ii,jj) = VR(ii,jj);
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1073
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1074 if (side == 'B')
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1075 for (int ii = 0; ii < nn; ii++)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1076 CVL(ii,jj) = VL(ii,jj);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1077 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1078 else
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1079 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1080 // double column; complex vector
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1081
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1082 for (int ii = 0; ii < nn; ii++)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1083 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1084 CVR(ii,jj) = Complex (VR(ii,jj), VR(ii,jj+1));
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1085 CVR(ii,jj+1) = Complex (VR(ii,jj), -VR(ii,jj+1));
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1086 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1087
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1088 if (side == 'B')
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1089 for (int ii = 0; ii < nn; ii++)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1090 {
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1091 CVL(ii,jj) = Complex (VL(ii,jj), VL(ii,jj+1));
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1092 CVL(ii,jj+1) = Complex (VL(ii,jj), -VL(ii,jj+1));
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1093 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1094 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1095
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1096 // advance to next eigenvectors (if any)
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1097 jj += cinc;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1098 }
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1099 }
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1100 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1101
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1102 switch (nargout)
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1103 {
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1104 case 7:
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1105 retval(6) = gev;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1106
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1107 case 6: // return eigenvectors
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1108 retval(5) = CVL;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1109
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1110 case 5: // return eigenvectors
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1111 retval(4) = CVR;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1112
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1113 case 4:
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1114 if (nargin == 3)
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1115 {
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1116 #ifdef DEBUG
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1117 std::cout << "qz: sort: retval(3) = gev = " << std::endl;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1118 octave_print_internal (std::cout, gev);
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1119 std::cout << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1120 #endif
10154 40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1121 retval(3) = gev;
40dfc0c99116 DLD-FUNCTIONS/.cc: untabify John W. Eaton <jwe@octave.org>* parents: 9786 diff changeset	1122 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1123 else
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1124 {if (complex_case)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1125 retval(3) = CZ;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1126 else
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1127 retval(3) = ZZ;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1128 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1129
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1130 case 3:
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1131 if (nargin == 3)
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1132 retval(2) = CZ;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1133 else
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1134 {if (complex_case)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1135 retval(2) = CQ.hermitian(); // compabible with MATLAB output
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1136 else
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1137 retval(2) = QQ.transpose();
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1138 }
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1139 case 2:
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1140 {if (complex_case)
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1141 {
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1142 #ifdef DEBUG
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1143 std::cout << "qz: retval (1) = cbb = " << std::endl;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1144 octave_print_internal (std::cout, cbb, 0);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1145 std::cout << std::endl << "qz: retval(0) = caa = " <<std::endl;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1146 octave_print_internal (std::cout, caa, 0);
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1147 std::cout << std::endl;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1148 #endif
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1149 retval(1) = cbb;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1150 retval(0) = caa;
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1151 }
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1152 else
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1153 { // real case
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1154 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	1155 std::cout << "qz: retval (1) = bb = " << std::endl;
3531 97cf542676e1 [project @ 2000-02-02 11:30:40 by jwe] jwe parents: 3523 diff changeset	1156 octave_print_internal (std::cout, bb, 0);
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	1157 std::cout << std::endl << "qz: retval(0) = aa = " <<std::endl;
3531 97cf542676e1 [project @ 2000-02-02 11:30:40 by jwe] jwe parents: 3523 diff changeset	1158 octave_print_internal (std::cout, aa, 0);
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	1159 std::cout << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1160 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1161 retval(1) = bb;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1162 retval(0) = aa;
10307 4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1163 }
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1164 break;}
4e4270ab70d6 support complex case in qz Jyh-miin Lin <jyhmiinlin@gmail.com> parents: 10155 diff changeset	1165
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1166
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1167 case 1:
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1168 case 0:
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1169 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	1170 std::cout << "qz: retval(0) = gev = " << gev << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1171 #endif
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1172 retval(0) = gev;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1173 break;
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1174
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1175 default:
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1176 error ("qz: too many return arguments.");
9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1177 break;
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1178 }
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1179
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1180 #ifdef DEBUG
3538 0ff7323dab8b [project @ 2000-02-02 12:36:25 by jwe] jwe parents: 3531 diff changeset	1181 std::cout << "qz: exiting (at long last)" << std::endl;
3185 9580887dd160 [project @ 1998-09-26 02:45:55 by jwe] jwe parents: 3183 diff changeset	1182 #endif
3183 5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1183
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1184 return retval;
5edc539c2f80 [project @ 1998-09-25 03:37:22 by jwe] jwe parents: diff changeset	1185 }

Mercurial > octave

annotate src/DLD-FUNCTIONS/qz.cc @ 10307:4e4270ab70d6