 

% testqz
sz=100;
rep=50;
tic
for jj=1:rep;
a=randn(sz,sz);b=(randn(sz,sz));
[s,t,q,z]=qzf(a,b,1);
res(jj,1)=max(max(q*t*z'-b));
res(jj,2)=max(max(q*s*z'-a));

end
disp('max fortran residual');
disp(max(abs(res)));
disp('fortran mean time per qz')
toc/rep


tic
for jj=1:rep;
a=randn(sz,sz);b=(randn(sz,sz));
[s,t,q,z]=qz(a,b);
[s,t,q,z]=ordqz(s,t,q,z,'udi');
res(jj,1)=max(max(q'*t*z'-b));
res(jj,2)=max(max(q'*s*z'-a));

end
disp('max matlab residual');
disp(max(abs(res)));
disp('matlab mean time per qz')
toc/rep