gsa_speed

 [tadj, iff] = gsa_speed(A,B,mf,p),

 Written by Marco Ratto
 Joint Research Centre, The European Commission,
 (http://eemc.jrc.ec.europa.eu/),
 marco.ratto@jrc.it 

 Reference:
 M. Ratto, Global Sensitivity Analysis for Macroeconomic models, MIMEO, 2006.

0001 function [tadj, iff] = gsa_speed(A,B,mf,p),
0002 % [tadj, iff] = gsa_speed(A,B,mf,p),
0003 %
0004 % Written by Marco Ratto
0005 % Joint Research Centre, The European Commission,
0006 % (http://eemc.jrc.ec.europa.eu/),
0007 % marco.ratto@jrc.it
0008 %
0009 % Reference:
0010 % M. Ratto, Global Sensitivity Analysis for Macroeconomic models, MIMEO, 2006.
0011 
0012 % Copyright (C) 2012 Dynare Team
0013 %
0014 % This file is part of Dynare.
0015 %
0016 % Dynare is free software: you can redistribute it and/or modify
0017 % it under the terms of the GNU General Public License as published by
0018 % the Free Software Foundation, either version 3 of the License, or
0019 % (at your option) any later version.
0020 %
0021 % Dynare is distributed in the hope that it will be useful,
0022 % but WITHOUT ANY WARRANTY; without even the implied warranty of
0023 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0024 % GNU General Public License for more details.
0025 %
0026 % You should have received a copy of the GNU General Public License
0027 % along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
0028 
0029 nvar=length(mf);
0030 nstate= size(A,1);
0031 nshock = size(B,2);
0032 nrun=size(B,3);
0033 
0034 iff=zeros(nvar,nshock,nrun);
0035 tadj=iff;
0036 disp('Computing speed of adjustement ...')
0037 h = dyn_waitbar(0,'Speed of adjustement...');
0038 
0039 for i=1:nrun,
0040   irf=zeros(nvar,nshock);
0041   a=squeeze(A(:,:,i));
0042   b=squeeze(B(:,:,i));
0043   IFF=inv(eye(nstate)-a)*b;
0044   iff(:,:,i)=IFF(mf,:);
0045   IF=IFF-b;
0046   
0047   t=0;
0048   while any(any(irf<0.5))
0049     t=t+1;
0050     IFT=((eye(nstate)-a^(t+1))*inv(eye(nstate)-a))*b-b;
0051     irf=IFT(mf,:)./(IF(mf,:)+eps);
0052     irf = irf.*(abs(IF(mf,:))>1.e-7)+(abs(IF(mf,:))<=1.e-7);
0053     %irf=ft(mf,:);
0054     tt=(irf>0.5).*t;
0055     tadj(:,:,i)=((tt-tadj(:,:,i))==tt).*tt+tadj(:,:,i);
0056   end
0057   dyn_waitbar(i/nrun,h)
0058 end
0059 disp(' ')
0060 disp('.. done !')
0061 dyn_waitbar_close(h)