2013-06-12 10:53:40 +02:00
function [err,ss,tt,w,sdim,eigval,info] = mjdgges ( e,d,qz_criterium, fake)
2010-09-20 17:51:25 +02:00
%function [err,ss,tt,w,sdim,eigval,info] = mjdgges(e,d,qz_criterium)
2008-04-25 14:47:30 +02:00
% QZ decomposition, Sims' codes are used.
%
% INPUTS
% e [double] real square (n*n) matrix.
% d [double] real square (n*n) matrix.
% qz_criterium [double] scalar (1+epsilon).
%
% OUTPUTS
2010-09-20 17:51:25 +02:00
% err [double] scalar: 1 indicates failure, 0 indicates success
2008-04-25 14:47:30 +02:00
% ss [complex] (n*n) matrix.
% tt [complex] (n*n) matrix.
% w [complex] (n*n) matrix.
% sdim [integer] scalar.
% eigval [complex] (n*1) vector.
% info [integer] scalar.
%
% ALGORITHM
% Sims's qzdiv routine is used.
%
% SPECIAL REQUIREMENTS
% none.
2008-08-01 14:40:33 +02:00
2010-09-22 11:40:13 +02:00
% Copyright (C) 1996-2010 Dynare Team
2008-08-01 14:40:33 +02:00
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
2009-12-16 18:17:34 +01:00
2010-09-20 17:51:25 +02:00
% Check number of inputs and outputs.
2013-06-12 10:53:40 +02:00
if nargin > 4 || nargin < 2 || nargout > 7 || nargout == 0
2010-09-20 17:51:25 +02:00
error ( ' MJDGGES: takes 2 or 3 input arguments and between 1 and 7 output arguments.' )
2008-04-25 14:47:30 +02:00
end
% Check the first two inputs.
[ me , ne ] = size ( e ) ;
[ md , nd ] = size ( d ) ;
2010-09-22 11:40:13 +02:00
if ( ~ isreal ( e ) || ~ isreal ( d ) || me ~= ne || md ~= nd || me ~= nd )
2008-04-25 14:47:30 +02:00
% info should be negative in this case, see dgges.f.
2008-06-06 14:50:47 +02:00
error ( ' MJDGGES requires two square real matrices of the same dimension.' )
2008-04-25 14:47:30 +02:00
end
2009-12-08 10:36:03 +01:00
2008-04-25 14:47:30 +02:00
% Set default value of qz_criterium.
if nargin < 3
qz_criterium = 1 + 1e-6 ;
end
2009-12-08 10:36:03 +01:00
info = 0 ;
% qz() function doesn't behave the same way under Octave and MATLAB:
% - under MATLAB, complex decomposition by default, real is also available
% as an option
% - under Octave, only real decomposition available, but grouping of
% eigenvalues <= 1 is implemented as an option (criterium can't be changed)
2013-11-04 10:54:45 +01:00
if isoctave
2009-12-08 10:36:03 +01:00
[ ss , tt , w , eigval ] = qz ( e , d , ' S' ) ;
sdim = sum ( abs ( eigval ) < = 1.0 ) ;
if any ( abs ( eigval ) > 1.0 & abs ( eigval ) < = qz_criterium )
2009-12-09 10:36:54 +01:00
warning ( ' Some eigenvalues are > 1.0 but <= qz_criterium in modulus. They have nevertheless been considered as explosive, because of a limitation of Octave. To solve this, you should compile the MEX files for Octave.' )
2009-12-08 10:36:03 +01:00
end
else
% Initialization of the output arguments.
ss = zeros ( ne , ne ) ;
tt = zeros ( ne , ne ) ;
w = zeros ( ne , ne ) ;
sdim = 0 ;
eigval = zeros ( ne , 1 ) ;
% Computational part.
try
[ ss , tt , qq , w ] = qz ( e , d ) ;
[ tt , ss , qq , w ] = qzdiv ( qz_criterium , tt , ss , qq , w ) ;
warning_old_state = warning ;
warning off ;
eigval = diag ( ss ) ./ diag ( tt ) ;
warning ( warning_old_state ) ;
sdim = sum ( abs ( eigval ) < qz_criterium ) ;
catch
info = 1 ; % Not as precise as lapack's info!
end
2009-12-07 13:25:24 +01:00
end
2010-09-20 17:51:25 +02:00
err = 0 ;