function [A,B,Q,Z] = qzswitch(i,A,B,Q,Z) Takes U.T. matrices A, B, orthonormal matrices Q,Z, interchanges diagonal elements i and i+1 of both A and B, while maintaining Q'AZ' and Q'BZ' unchanged. If diagonal elements of A and B are zero at matching positions, the returned A will have zeros at both positions on the diagonal. This is natural behavior if this routine is used to drive all zeros on the diagonal of A to the lower right, but in this case the qz transformation is not unique and it is not possible simply to switch the positions of the diagonal elements of both A and B.
0001 function [A,B,Q,Z] = qzswitch(i,A,B,Q,Z) 0002 %function [A,B,Q,Z] = qzswitch(i,A,B,Q,Z) 0003 % 0004 % Takes U.T. matrices A, B, orthonormal matrices Q,Z, interchanges 0005 % diagonal elements i and i+1 of both A and B, while maintaining 0006 % Q'AZ' and Q'BZ' unchanged. If diagonal elements of A and B 0007 % are zero at matching positions, the returned A will have zeros at both 0008 % positions on the diagonal. This is natural behavior if this routine is used 0009 % to drive all zeros on the diagonal of A to the lower right, but in this case 0010 % the qz transformation is not unique and it is not possible simply to switch 0011 % the positions of the diagonal elements of both A and B. 0012 0013 % Original file downloaded from: 0014 % http://sims.princeton.edu/yftp/gensys/mfiles/qzswitch.m 0015 0016 % Copyright (C) 1993-2007 Christopher Sims 0017 % Copyright (C) 2008-2011 Dynare Team 0018 % 0019 % This file is part of Dynare. 0020 % 0021 % Dynare is free software: you can redistribute it and/or modify 0022 % it under the terms of the GNU General Public License as published by 0023 % the Free Software Foundation, either version 3 of the License, or 0024 % (at your option) any later version. 0025 % 0026 % Dynare is distributed in the hope that it will be useful, 0027 % but WITHOUT ANY WARRANTY; without even the implied warranty of 0028 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0029 % GNU General Public License for more details. 0030 % 0031 % You should have received a copy of the GNU General Public License 0032 % along with Dynare. If not, see <http://www.gnu.org/licenses/>. 0033 0034 realsmall=sqrt(eps)*10; 0035 %realsmall=1e-3; 0036 a = A(i,i); d = B(i,i); b = A(i,i+1); e = B(i,i+1); 0037 c = A(i+1,i+1); f = B(i+1,i+1); 0038 % A(i:i+1,i:i+1)=[a b; 0 c]; 0039 % B(i:i+1,i:i+1)=[d e; 0 f]; 0040 if (abs(c)<realsmall && abs(f)<realsmall) 0041 if abs(a)<realsmall 0042 % l.r. coincident 0's with u.l. of A=0; do nothing 0043 return 0044 else 0045 % l.r. coincident zeros; put 0 in u.l. of a 0046 wz=[b; -a]; 0047 wz=wz/sqrt(wz'*wz); 0048 wz=[wz [wz(2)';-wz(1)'] ]; 0049 xy=eye(2); 0050 end 0051 elseif (abs(a)<realsmall && abs(d)<realsmall) 0052 if abs(c)<realsmall 0053 % u.l. coincident zeros with l.r. of A=0; do nothing 0054 return 0055 else 0056 % u.l. coincident zeros; put 0 in l.r. of A 0057 wz=eye(2); 0058 xy=[c -b]; 0059 xy=xy/sqrt(xy*xy'); 0060 xy=[[xy(2)' -xy(1)'];xy]; 0061 end 0062 else 0063 % usual case 0064 wz = [c*e-f*b, (c*d-f*a)']; 0065 xy = [(b*d-e*a)', (c*d-f*a)']; 0066 n = sqrt(wz*wz'); 0067 m = sqrt(xy*xy'); 0068 if m<eps*100 0069 % all elements of A and B proportional 0070 return 0071 end 0072 wz = n\wz; 0073 xy = m\xy; 0074 wz = [wz; -wz(2)', wz(1)']; 0075 xy = [xy;-xy(2)', xy(1)']; 0076 end 0077 A(i:i+1,:) = xy*A(i:i+1,:); 0078 B(i:i+1,:) = xy*B(i:i+1,:); 0079 A(:,i:i+1) = A(:,i:i+1)*wz; 0080 B(:,i:i+1) = B(:,i:i+1)*wz; 0081 Z(:,i:i+1) = Z(:,i:i+1)*wz; 0082 Q(i:i+1,:) = xy*Q(i:i+1,:);