%DRIVER 
%
%  This file is the main file to be run to replicate numerical results 
%  reported in Chari, Kehoe, and McGrattan, ``Can Sticky Price Models 
%  Generate Volatile and Persistent Real Exchange Rates?''  
%
%  Users can run their own examples by changing the parameter values
%  described below.
%
%  Parameters related to price and wage staggered contracting:
%
%      T     = fraction of year prices are fixed (T=1 ==> fixed 1 year)
%      N     = number of intermediate-good-firm cohorts
%      Nw    = number of household cohorts
%
%  Parameters of the utility function U(C,L,M)
%
%      Utype = form for utility
%            = 1 if U = [CES(C,M)*(1-L)^psi]^(1-sigma)/(1-sigma)
%            = 2 if U = [CES(C,M)]^(1-sigma)/(1-sigma)+psi*(1-L)^(1-xi)/(1-xi) 
%            = 3 if U = V^(1-sigma)/(1-sigma)
%                   V = C^(1-kappa)/(1-kappa)-psi*L^(1+xi)/(1+xi)+omega*M^nu/nu
%            = 4 if U = V^(1-sigma)/(1-sigma)
%                   V = C^(1-kappa)/(1-kappa)+psi*(1-L)^(1-xi)/(1-xi)
%                         +omega*M^nu/nu
%
%            with CES(C,M) = [omega*C^nu+(1-omega)*M^nu]^(1/nu)
%                       nu = (eta-1)/eta
%                 
%
%      eta   = real >= 0 (interest elasicity of MD)
%      kappa = real >= 0
%      omega = real in [0,1]
%      psi   = real >= 0
%      sigma = real >= 0
%      xi    = real >= 0 
%      beta  = the discount factor, real in (0,1)
%
%  Parameters of technology:
%                           Y  = (ome1*Yh^rho+(1-ome1)*Yf^rho)^(1/rho)
%                           Yh = (int Yh(i)^theta di)^(1/theta)
%                           Yh(i) + Yh*(i) = k(i)^alpha l(i)^(1-alpha)
%                           k[i,t] = (1-delta)k[i,t-1] + x[i,t]
%                                      -b/2(x[i,t]/k[i,t-1]-delta)^2*k[i,t-1]
%
%      ome1  = weight on home goods 
%      rho   = parameter related to elasticity of substitution between
%              home and foreign goods (1/(1-rho)), real <= 1
%      theta = inverse of markup 
%      alpha = capital share, real in (0,1)
%      delta = depreciation rate, real in [0,1]
%      b     = adjustment cost parameter, real >= 0
%
%  Parameters of monetary policy
%
%      Ptype = type of policy assumed
%            = 1 if Fed changes monetary growth rates 
%                (log mu[t+1] = (1-rho_m)*mu+rho_m*log mu[t] + eps[t+1])
%            = 2 if Taylor rule used 
%                (r[t] = irule' [    r[t-1]
%                                    r[t-2]
%                                    r[t-3]
%                                E_t inflation[t,t+1]
%                                    inflation[t-1,t]
%                                    inflation[t-2,t-1]
%                                    inflation[t-3,t-2]
%                                    gdp[t]
%                                    gdp[t-1]
%                                    gdp[t-2]    ]  )
%
%
%      mu    = real >= 1      
%      rhom  = real in (-1,1) 
%      sigm  = real >= 0
%      corm  = real in [0,1)
%      irule = real 10 x 1 coefficient matrix
%      rhoe  = real in (-1,1)
%      sige  = real >= 0
%      core  = real in [0,1)
%
%  Parameters of AR(1) shocks (g = government spending, a=technology)
%
%      gs    = real >=0
%      rhog  = real in (-1,1)
%      sigg  = real >=0
%      corg  = real in [0,1)
%      rhoa  = real in (-1,1)
%      siga  = real >=0
%      cora  = real in [0,1)
%
%  Parameters related to Sticky Wage case
%
%      Wtype = 1 if no sticky wages assumed
%            = 2 if sticky wages assumed
%      v     = inverse of markup of real wage over -UL/UC
%
%  Timing of shocks
%
%      timing= 1, monopolists observe only lagged shocks when choosing p[t-1]
%            = 2, monopolists observe lagged shocks and a[t]
%            = 3, monopolists observe lagged shocks and g[t],a[t]
%
%  Market structure 
%
%      Mtype = 1 if complete markets assumed
%            = 2 if incomplete markets assumed
%
%  Parameters for simulations (skip if 0)
%     
%      nti   = number of time periods for impulse responses
%      nts   = number of time periods per stochastic simulation
%      nsim  = number of stochstic simulations
%
%
%  The outputs are impulse responses (imp) for
%
%      c:           home consumption
%      l:           home employment
%      gdp:         home GDP
%      c*:          foreign consumption
%      l*:          foreign employment
%      gdp*:        foreign GDP
%      M/mu^t:      normalized money supply at home
%      P/mu^{t-1}:  normalized aggregate price at home
%      e:           nominal exchange rate
%      q:           real exchange rate
%      r:           nominal interest rate at home
%      r-pi:        real interest rate at home
%
%  and time series (sim) for 
%
%      c:           home consumption
%      x:           home investment 
%      l:           home employment
%      gdp:         home real GDP
%      ngdp:        home nominal GDP
%      c*:          foreign consumption
%      x*:          foreign investment 
%      l*:          foreign employment
%      gdp*:        foreign real GDP
%      P*/P:        price ratio
%      e:           nominal exchange rate
%      q:           real exchange rate
%      rnxy:        home real net exports/real gdp
%      adj:         home adjustment costs on capital
%       


%  Ellen McGrattan, 9-15-99
%  Revised, ERM, 10-1-04


%  CASES reported in the paper:
%      1-16:  Vary N
%     17-24:  Sensitivity Analysis,
%             e.g., High Exports, Taylor rule, Incomplete Markets,
%                   Sticky Wages, Uncorrelated Money, Persistent Money
%                   Growth, Real Shocks added, Nonseparable Preferences
%     25-33:  Vary Sigma
%     34-42:  Vary Sigma, Fix adjustment costs
%     43-47:  Vary Correlation of money growth

T         = [  .25,  0.5, 0.75,  1.0, 1.25,  1.5, 1.75,  2.0, ...
              2.25,  2.5, 2.75,  3.0, 3.25,  3.5, 3.75,  4.0, ...
               1.0,  1.0,  1.0,  1.0,  1.0,  1.0,  1.0,  1.0, ...
               1.0,  1.0,  1.0,  1.0,  1.0,  1.0,  1.0,  1.0,  1.0, ...
               1.0,  1.0,  1.0,  1.0,  1.0,  1.0,  1.0,  1.0,  1.0, ...
               1.0,  1.0,  1.0,  1.0,  1.0];
ncase     = length(T);
N         = [    1,    2,    3,    4,    5,    6,    7,    8, ...
                 9,   10,   11,   12,   13,   14,   15,   16, ...
                 4,    4,    4,    4,    4,    4,    4,    4, ...
                 4,    4,    4,    4,    4,    4,    4,    4,   4,...
                 4,    4,    4,    4,    4,    4,    4,    4,   4,...
                 4,    4,    4,    4,    4];
Nw        =  4  * ones(1,ncase);
Utype     =  2  * ones(1,ncase);
eta       = .39 * ones(1,ncase);
kappa     = .99 * ones(1,ncase);
omega     = .94 * ones(1,ncase);
Utype(18) = 4;
kappa(18) = 6;
omega(18) = 1-.97^(T(8)/N(8));
Utype(24) = 1;
psi       = [   10,   10,   10,   10,   10,   10,   10,   10,  ...
                10,   10,   10,   10,   10,   10,   10,   10,  ...
                80,  100,   10,    9,   10,   10,   40, 2.25,  ...
               2.3,  3.3,    5,    7,   15,   20,   30,   45,   80,  ...
               2.3,  3.3,    5,    7,   15,   20,   30,   45,   80,  ...
                10,   10,   10,   10,   10];
sigma     = [    5,    5,    5,    5,    5,    5,    5,    5,  ...
                 5,    5,    5,    5,    5,    5,    5,    5,  ...
                 5,.0001,    5,    5,    5,    5,    5,    5,  ...
              1.01,    2,    3,    4,    6,    7,    8,    9,   10,  ...
              1.01,    2,    3,    4,    6,    7,    8,    9,   10,  ...
                 5,    5,    5,    5,    5];
xi        = sigma;
beta      = 0.97 .^[T./N];
rho       = 1/3 * ones(1,ncase);
ome1      = (.984/.016).^(1-rho)./(1+(.984/.016).^(1-rho));
ome1(17)  = (.850/.150)^(2/3)    /(1+(.850/.150)^(2/3));
theta     = 0.9 * ones(1,ncase);
alpha     = 1/3 * ones(1,ncase);
delta     = 1-.92  .^[T./N];
b         = [  124,  129,  135,  138,  140,  143,  149,  152, ...
               156,  161,  165,  171,  177,  181,  183,  185, ...
               139,  155,  135,  133,  135,  135,   93,   .1, ...
                27,   52,   80,  107,  166,  204,  242,  288,  335, ...
               138,  138,  138,  138,  138,  138,  138,  138,  138, ...
               138,  138,  138,  138,  138];
Ptype     = ones(1,ncase);
Ptype(18) = 2;
mu        = 1.0   .^[T./N];
rhom      = 0.683 .^(4*T./N);
rhom(22)  = 0.99  .^(4*T(22)/N(22));
corm      = [   .5,   .5,   .5,   .5,   .5,   .5,   .5,   .5, ...
                .5,   .5,   .5,   .5,   .5,   .5,   .5,   .5, ...
                .5,   .5,   .5,   .5,   .0,   .5,   .5,   .5, ...
                .5,   .5,   .5,   .5,   .5,   .5,   .5,   .5,   .5,  ...
                .5,   .5,   .5,   .5,   .5,   .5,   .5,   .5,   .5,  ...
                .2,   .4,   .6,   .8,  .95];
sigm      = [.0226,.0198,.0181,.01685,.0160,.0154,.0151,.0148, ...
             .0146,.0144,.0143,.0142,.0142,.0141,.0140,.01385, ...
             .0170,.0171,.0167,.0151,.0192,.00147,.0168,.000192, ...
            .00365,.0068,.0101,.0134,.0202,.0240,.0277,.0317,.0356, ...
            .0050,.00865,.0118,.0145,.0189,.0208,.0224,.0239,.0252, ...
             .0190,.0178,.0156,.01175,.00613];
rhor      = 0.66    * ones(1,ncase);
ar        = 1.8     * ones(1,ncase);
br        = 0.12/4  * ones(1,ncase);
rhoe      = zeros(1,ncase);
sige      = 0.0273   * ones(1,ncase);
sige      = 0.0281   * ones(1,ncase);  % use for R/W
sige      = 0.0228   * ones(1,ncase);  % use for CGG
core      = zeros(1,ncase);
core(18)  = .5;
gs        = 10^(-7) * ones(1,ncase);
rhog      = zeros(1,ncase) ;
sigg      = 10^(-7) * ones(1,ncase);
corg      = zeros(1,ncase);
rhoa      = zeros(1,ncase);
siga      = 10^(-7) * ones(1,ncase);
cora      = zeros(1,ncase);
timing    = ones(1,ncase);
gs(23)    = .13;
rhog(23)  = .97;
sigg(23)  = .01;
rhoa(23)  = .95;
siga(23)  = .007;
cora(23)  = .25;
timing(23)= 3;
Mtype     = ones(1,ncase);
Mtype(19) = 2;
nti       = 30  * ones(1,ncase);
nts       = 88  * ones(1,ncase);
nsim      = 100 * ones(1,ncase);
Wtype     = ones(1,ncase);
Wtype(20) = 2;
v         = 0.87 * ones(1,ncase);
irule     = [.66;0;0;0;.34*1.8;0;0;0;.34*.03;0]* ones(1,ncase);    %CGG 1998
irule     = [.507;.023;.158;0;-.062;.427;.3; 
             .624/4;-.059/4;-.418/4]* ones(1,ncase);               %RW  1997
irule     = [.79;0;0;.21*2.15;0;0;0;.21*.93/4;0;0]* ones(1,ncase); %CGG 2000
irule     = irule+1e-10;  % needed to aid computation of generalized eigs

params    = [T;N;Nw;Utype;eta;kappa;omega;psi;sigma;xi;beta;ome1;rho;theta;
             alpha;delta;b;Ptype;mu;rhom;sigm;corm;irule;rhoe;sige;core;
             gs;rhog;sigg;corg;rhoa;siga;cora;timing;Mtype;Wtype;v;
             nti;nts;nsim];
%            1 2 3   4     5   6     7    8     9   10  11   12  13   14
%             15    16  17  18   19  20   21   22   23:32 33   34  35
%            36  37   38   39   40   41   42   43     44    45   46
%            47  48   49  

mti       = max(nti);
mts       = max(nts);
msim      = max(nsim);


randn('seed',91262)
mshocks   = randn(mts,msim);
mashocks  = randn(mts,msim);
gshocks   = randn(mts,msim);
gashocks  = randn(mts,msim);
ashocks   = randn(mts,msim);
aashocks  = randn(mts,msim);
Imp       = zeros(mti*ncase,12);
T1m       = zeros(10,ncase);
T1s       = zeros(10,ncase);
T2m       = zeros(15,ncase);
T2s       = zeros(15,ncase);
T3        = zeros(60,ncase);
T4        = zeros(3,ncase);
Adj       = zeros(2,ncase);
ifig      = 0;
iplot     = 0;

ibase     = 4;   % baseline case

diary temp.out
disp(sprintf('Computing and simulating %g economies...',ncase))
disp(       ['Results appear on the screen and in temp.out ', ...
             'when the program is complete.'])

for i=1:ncase;
  p      = params(:,i);
  shocks = [mshocks(1:nts(i),1:nsim(i)),mashocks(1:nts(i),1:nsim(i)), ...
            gshocks(1:nts(i),1:nsim(i)),gashocks(1:nts(i),1:nsim(i)), ...
            ashocks(1:nts(i),1:nsim(i)),aashocks(1:nts(i),1:nsim(i))];
  
  if p(44) == 1;             % complete markets
    if p(18) == 1;           % exogenous money growth
      if p(45) == 1;         % sticky prices only
        if p(2) == 1;        % N=1

           [imp,sim] = iscmem1_ms(p,shocks);     

        else                 

           [imp,sim] = iscmem_ms(p,shocks);

        end


      else

        if p(2) == 1;

           %[imp,sim] = iscmemsw1_ms(p,shocks);     
           disp(' This option is not yet available')

        else

           [imp,sim] = iscmemsw_ms(p,shocks);

        end

      end
    else
      if p(2) == 1;

         [imp,sim] = iscmtr1_ms(p,shocks);

      else

         [imp,sim] = iscmtr_ms(p,shocks);

      end;
    end
  elseif p(44) == 2;

    if p(2) == 1;

       [imp,sim] = isimem1_ms(p,shocks);

    else

       [imp,sim] = isimem_ms(p,shocks);

    end

  end


  % GENERATE FIGURES
  %
  %    Figure 1.  Impulse response of home C,L,GDP
  %    Figure 2.  Impulse response of foreign C,L,GDP
  %    Figure 3.  Impulse response of home M,P
  %    Figure 4.  Impulse response of exchange rates
  %    Figure 5.  Impulse response of interest rates
  %

  if p(47) > 0;

    Imp((i-1)*mti+1:(i-1)*mti+p(47),:) = imp;

    if iplot;
      ifig  = ifig + 1;
      istr  = int2str(i);
      figure(ifig)
      orient tall
      subplot(511), 
      plot(0:p(47)-1,imp(:,1:3)),   
      title(['Case ',istr,': Responses of C,L,GDP'])

      subplot(512), 
      plot(0:p(47)-1,imp(:,4:6)),   
      title(['Case ',istr,': Responses of C*,L*,GDP*'])

      subplot(513), 
      plot(0:p(47)-1,imp(:,7:8)),   
      title(['Case ',istr,': Responses of M/mu^t,P/mu^t'])

      subplot(514), 
      plot(0:p(47)-1,imp(:,9:10)),  
      title(['Case ',istr,': Responses of q,e'])

      subplot(515), 
      plot(0:p(47)-1,imp(:,11:12)), 
      title(['Case ',istr,': Responses of Real R, Nominal R'])
    end;

  end;

  % GENERATE TABLES
  %  
  %    * Exchange Rates and Prices
  %    * Business Cycle Statistics
  %

  n        = nts(i);
  ns       = nsim(i);
  o        = zeros(n,ns);
  c        = o; c(:)     = sim(:,1);
  x        = o; x(:)     = max(sim(:,2),.0001);
  l        = o; l(:)     = sim(:,3);
  gdp      = o; gdp(:)   = sim(:,4);
  ngdp     = o; ngdp(:)  = sim(:,5);
  ca       = o; ca(:)    = sim(:,6);
  xa       = o; xa(:)    = max(sim(:,7),.0001);
  la       = o; la(:)    = sim(:,8);
  gdpa     = o; gdpa(:)  = sim(:,9);
  prat     = o; prat(:)  = sim(:,10);
  e        = o; e(:)     = sim(:,11);
  q        = o; q(:)     = sim(:,12);
  rnxy     = o; rnxy(:)  = sim(:,13);
  adj      = o; adj(:)   = sim(:,14);

  dc       = log(c)    - hptrend(log(c),1600);
  dx       = log(x)    - hptrend(log(x),1600);
  dl       = log(l)    - hptrend(log(l),1600);
  dgdp     = log(gdp)  - hptrend(log(gdp),1600);
  dngdp    = log(ngdp) - hptrend(log(ngdp),1600);
  dca      = log(ca)   - hptrend(log(ca),1600);
  dxa      = log(xa)   - hptrend(log(xa),1600);
  dla      = log(la)   - hptrend(log(la),1600);
  dgdpa    = log(gdpa) - hptrend(log(gdpa),1600);
  dprat    = log(prat) - hptrend(log(prat),1600);
  de       = log(e)    - hptrend(log(e),1600);
  dq       = log(q)    - hptrend(log(q),1600);
  drnxy    = rnxy      - hptrend(rnxy,1600);
  ddc      = dc-dca;
  ddl      = dl-dla;

  A        = [dq, dq,   dq,    dq,  dgdp,  dc,  dx,  dl,  drnxy];
  B        = [de, dgdp, drnxy, ddc, dgdpa, dca, dxa, dla, dgdp];
  C        = (sum(A.*B)-mean(A).*mean(B)*n)/(n-1);
  D        = (sum(A.*A)-mean(A).*mean(A)*n)/(n-1);
  E        = (sum(B.*B)-mean(B).*mean(B)*n)/(n-1);

  cc       = zeros(ns,9);
  cc(:)    = C'./sqrt(D'.*E');      
  if ns>1;
    cc       = [mean(cc)',std(cc)'];     % = cross-correlations
  else
    cc       = [cc',0*cc'];
  end;

  A        = [dprat, de, dq, dgdp, dc, dx, dl, drnxy];
  sd       = zeros(ns,8);
  sd(:)    = std(A)' *100;
  sd       = [sd(:,1)./sd(:,4),sd(:,2)./sd(:,4),sd(:,3)./sd(:,4),sd(:,4), ...
              sd(:,8),sd(:,5)./sd(:,4),sd(:,6)./sd(:,4),sd(:,7)./sd(:,4)];
  if ns>1;
    sd       = [mean(sd)',std(sd)'];     % = standard deviations
  else
    sd       = [sd',0*sd'];
  end;

  B        = A(2:n,:);
  A        = A(1:n-1,:);
  C        = (sum(A.*B)-mean(A).*mean(B)*(n-1))/(n-2);
  D        = (sum(A.*A)-mean(A).*mean(A)*(n-1))/(n-2);
  E        = (sum(B.*B)-mean(B).*mean(B)*(n-1))/(n-2);
  
  ac       = zeros(ns,8);
  ac(:)    = C'./sqrt(D'.*E');      
  if ns>1;
    ac       = [mean(ac)',std(ac)'];     % = autocorrelations
  else
    ac       = [ac',0*ac'];
  end;

  T1m(:,i) = [sd(1:3,1);ac(1:3,1);cc(1:4,1)];
  T1s(:,i) = [sd(1:3,2);ac(1:3,2);cc(1:4,2)];

  T2m(:,i) = [sd(4:8,1);ac(4:8,1);cc(5:9,1)];
  T2s(:,i) = [sd(4:8,2);ac(4:8,2);cc(5:9,2)];

  Adj(:,i) = [mean(mean(adj')),std(mean(adj'))]';

  A        = [ddc, ddl, dgdp];
  sd       = zeros(ns,3);
  sd(:)    = std(A)' *100;
  sd       = [sd(:,1)./sd(:,3),sd(:,2)./sd(:,3)];
  if ns>1;
    sd       = mean(sd.^2)';            % = autocorrelations
  else
    sd       = (sd.^2)';
  end;

  T4(1:2,i)= sd;

  A        = ddc;
  B        = ddl;
  C        = (sum(A.*B)-mean(A).*mean(B)*n)/(n-1);
  T4(3,i)  = mean(C./std(dgdp).^2);   % = var(x)/var(dgdp), x=dc-dc*,dl-dl*

  nex1     = o; nex1(:)  = sim(:,15);
  nex2     = o; nex2(:)  = sim(:,16);
  nim1     = o; nim1(:)  = sim(:,17);
  nim2     = o; nim2(:)  = sim(:,18);
  rex1     = o; rex1(:)  = sim(:,19);
  rim1     = o; rim1(:)  = sim(:,20);
  tot1     = o; tot1(:)  = sim(:,21);
  tot2     = o; tot2(:)  = sim(:,22);
  tot3     = o; tot3(:)  = sim(:,23);

  dnxy1    = (nex1-nim1)./ngdp - hptrend((nex1-nim1)./ngdp,1600);
  dnxy2    = (nex2-nim2)./ngdp - hptrend((nex2-nim2)./ngdp,1600);
  dnxy3    = (nex2-nim1)./ngdp - hptrend((nex2-nim1)./ngdp,1600);
  dtot1    = log(tot1) - hptrend(log(tot1),1600);
  dtot2    = log(tot2) - hptrend(log(tot2),1600);
  dtot3    = log(tot3) - hptrend(log(tot3),1600);

  A        = [dnxy1, dnxy2, dnxy3, dnxy1, dnxy2, dnxy3, ...
              dnxy1, dnxy2, dnxy3, dtot1, dtot2, dtot3, ...
              dtot1, dtot2, dtot3, dtot1, dtot2, dtot3];
  B        = [ dgdp,  dgdp,  dgdp,    de,    de,    de, ...
                 dq,    dq,    dq,  dgdp,  dgdp,  dgdp, ...
                 de,    de,    de,    dq,    dq,    dq];
  C        = (sum(A.*B)-mean(A).*mean(B)*n)/(n-1);
  D        = (sum(A.*A)-mean(A).*mean(A)*n)/(n-1);
  E        = (sum(B.*B)-mean(B).*mean(B)*n)/(n-1);

  cc       = zeros(ns,18);
  cc(:)    = C'./sqrt(D'.*E');      
  if ns>1;
    cc       = [mean(cc)',std(cc)'];     % = cross-correlations
  else
    cc       = [cc',0*cc'];
  end;

  A        = [dnxy1,dnxy2,dnxy3,dtot1,dtot2,dtot3];
  sd       = zeros(ns,6);
  sd(:)    = std(A)' *100;
  if ns>1;
    sd       = [mean(sd)',std(sd)'];     % = standard deviations
  else
    sd       = [sd',0*sd'];
  end;

  B        = A(2:n,:);
  A        = A(1:n-1,:);
  C        = (sum(A.*B)-mean(A).*mean(B)*(n-1))/(n-2);
  D        = (sum(A.*A)-mean(A).*mean(A)*(n-1))/(n-2);
  E        = (sum(B.*B)-mean(B).*mean(B)*(n-1))/(n-2);
  
  ac       = zeros(ns,6);
  ac(:)    = C'./sqrt(D'.*E');      
  if ns>1;
    ac       = [mean(ac)',std(ac)'];     % = autocorrelations
  else
    ac       = [ac',0*ac'];
  end;

  T3(:,i)  = [sd(1:3,1);sd(1:3,2);sd(4:6,1);sd(4:6,2);
              ac(1:3,1);ac(1:3,2);ac(4:6,1);ac(4:6,2);
              cc(1:3,1);cc(1:3,2);cc(4:6,1);cc(4:6,2);
              cc(7:9,1);cc(7:9,2);cc(10:12,1);cc(10:12,2);
              cc(13:15,1);cc(13:15,2);cc(16:18,1);cc(16:18,2)];

  disp(sprintf(' .... Case %g complete',i))
end;
params = [params;1:ncase];
disp('PARAMETERS')
disp('----------')
disp(        '                                                ')
disp(sprintf('    Baseline case:  %3g                         ',ibase))
disp(        '                                                ')
disp(        '    Parameters related to staggered contracting:')
disp(        '                                                ')
disp(        '        Case   N    Nw    T                     ')
disp(sprintf('         %3g   %2g   %2g    %4.2f\n',params([50,2,3,1],:)))
disp(        '                                                ')
disp(        '    Parameters of the utility function U(C,L,M) ')
disp(        '                                                ')
disp(        '        Case  Utype eta  kappa omega psi   sigma xi   beta')
disp(sprintf(['         %3g   %2g   %4.2f %5.3f %4.2f  %5.1f %4.1f', ...
              '  %4.1f %6.4f\n'],params([50,4:11],:)))
disp(        '                                                ')
disp(        '    Parameters of technology Y, Yh, Yh(i), k(i,t)')
disp(        '                                                ')
disp(        '        Case  ome1   rho    theta  alpha  delta  b')
disp(sprintf('         %3g  %6.4f %6.4f %4.2f   %4.2f   %6.4f %6.2f\n', ...
                                params([50,12:17],:)))
disp(        '                                                ')  
disp(        '    Parameters of exogenous money policy        ')
disp(        '                                                ')  
disp(        '        Case  Ptype mu    rhom sigm     corm    ')
disp(sprintf('         %3g   %2g   %5.3f %4.2f %8.2e %4.2f\n', ...
                                params([50,18:22],:)))
disp(        '                                                ')  
disp(        '    Parameters of Taylor rule                   ')
disp(        '                                                ')
disp(        '        Case  rhoe sige     core   irule        ')
disp(sprintf(['         %3g  %4.2f %8.2e %4.2f   %5.2f %5.2f %5.2f %5.2f', ...
              ' %5.2f %5.2f %5.2f %5.2f %5.2f %5.2f\n'], ...
                                params([50,33:35,23:32],:)))
disp(        '                                                ')  
disp(        '    Parameters of exogenous real shocks         ')
disp(        '                                                ')  
disp(        ['        Case   gs   rhog sigg     corg', ...
              ' rhoa siga     cora timing'])
disp(sprintf(['         %3g   %4.2f %4.2f %8.2e %4.2f ', ...
              '%4.2f %8.2e %4.2f %1g\n'],params([50,36:43],:)))
disp(        '                                                ')  
disp(        '    Parameter indicating Market structure       ')
disp(        '                                                ')  
disp(        '        Case  Mtype                             ')
disp(sprintf('         %3g   %2g         \n',params([50,44],:)))
disp(        '                                                ')  
disp(        '    Parameters related to Sticky Wage cases     ')
disp(        '                                                ')  
disp(        '        Case  Wtype  v                          ')
disp(sprintf('         %3g   %2g    %4.2f\n',params([50,45,46],:)))
disp(        '                                                ')  
disp(        '    Parameters for simulations                  ')
disp(        '                                                ')  
disp(        '        Case  nti  nts  nsim                    ')
disp(sprintf('         %3g  %3g  %3g  %3g\n',params([50,47:49],:)))
disp(        '                                                ')  
disp(        '                                                ')  


disp('RESULTS REPORTED IN THE TEXT')
disp('----------------------------')

cases=char(' N = 1 ',' N = 2 ',' N = 3 ',' N = 4 ', ...
           ' N = 5 ',' N = 6 ',' N = 7 ',' N = 8 ');

disp(' ')
disp(' ')
disp('EXCHANGE RATES AND PRICES -- VARYING N')
disppr([T1m(:,1:8),T1s(:,1:8)]',cases)
disp(' ')
disp(' ')
disp('BUSINESS CYCLE STATISTICS -- VARYING N')
dispbc([T2m(:,1:8),T2s(:,1:8)]',cases)

cases=char(' N = 9 ', ' N = 10 ',' N = 11 ',' N = 12 ', ...
           ' N = 13 ',' N = 14 ',' N = 15 ',' N = 16 ');

disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp('EXCHANGE RATES AND PRICES -- VARYING N')
disppr([T1m(:,9:16),T1s(:,9:16)]',cases)
disp(' ')
disp(' ')
disp('BUSINESS CYCLE STATISTICS -- VARYING N')
dispbc([T2m(:,9:16),T2s(:,9:16)]',cases)


cases=char(' Baseline ',' High Exports ',' Taylor Rule ', ...
           ' Inc. Markets ',' Sticky Wages ',' Uncorr. Money ', ...
           ' Pers. Money ',' All Shocks ',' Nonsep. Pref.');
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp('EXCHANGE RATES AND PRICES - SENSITIVITY ANALYSIS')
disppr([T1m(:,[ibase,17:24]),T1s(:,[ibase,17:24])]',cases)
disp(' ')
disp(' ')
disp('BUSINESS CYCLE STATISTICS - SENSITIVITY ANALYSIS')
dispbc([T2m(:,[ibase,17:24]),T2s(:,[ibase,17:24])]',cases)

cases=char(' Baseline ',[' Sigma=',num2str(sigma(25))], ...
                        [' Sigma=',num2str(sigma(26))], ...
                        [' Sigma=',num2str(sigma(27))], ...
                        [' Sigma=',num2str(sigma(28))], ...
                        [' Sigma=',num2str(sigma(29))], ...
                        [' Sigma=',num2str(sigma(30))], ...
                        [' Sigma=',num2str(sigma(31))], ...
                        [' Sigma=',num2str(sigma(32))], ...
                        [' Sigma=',num2str(sigma(33))]);
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp('EXCHANGE RATES AND PRICES - VARYING SIGMA')
disppr([T1m(:,[ibase,25:33]),T1s(:,[ibase,25:33])]',cases)
disp(' ')
disp(' ')
disp('BUSINESS CYCLE STATISTICS - VARYING SIGMA')
dispbc([T2m(:,[ibase,25:33]),T2s(:,[ibase,25:33])]',cases)

disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp('EXCHANGE RATES AND PRICES - VARY SIGMA, FIXING AC')
disppr([T1m(:,[ibase,34:42]),T1s(:,[ibase,34:42])]',cases)
disp(' ')
disp(' ')
disp('BUSINESS CYCLE STATISTICS - VARY SIGMA, FIXING AC')
dispbc([T2m(:,[ibase,34:42]),T2s(:,[ibase,34:42])]',cases)


cases=char(' Baseline ',[' Cor M=',num2str(corm(43))], ...
                        [' Cor M=',num2str(corm(44))], ...
                        [' Cor M=',num2str(corm(45))], ...
                        [' Cor M=',num2str(corm(46))], ...
                        [' Cor M=',num2str(corm(47))]);
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp('EXCHANGE RATES AND PRICES - VARYING CORR(M,M*)')
disppr([T1m(:,[ibase,43:47]),T1s(:,[ibase,43:47])]',cases)
disp(' ')
disp(' ')
disp('BUSINESS CYCLE STATISTICS - VARYING CORR(M,M*)')
dispbc([T2m(:,[ibase,43:47]),T2s(:,[ibase,43:47])]',cases)

cases=char(' ePh*,Pf ',' Ph,ePf* ', ' Ph,Pf ');
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp(' ')
disp('ALTERNATIVE DECENTRALIZATIONS')
dispde(T3(:,ibase),cases)


cases=char(' Baseline ',[' Sigma=',num2str(sigma(25))], ...
                        [' Sigma=',num2str(sigma(26))], ...
                        [' Sigma=',num2str(sigma(27))], ...
                        [' Sigma=',num2str(sigma(28))], ...
                        [' Sigma=',num2str(sigma(29))], ...
                        [' Sigma=',num2str(sigma(30))], ...
                        [' Sigma=',num2str(sigma(31))], ...
                        [' Sigma=',num2str(sigma(32))], ...
                        [' Sigma=',num2str(sigma(33))]);
disp(' ')
disp(' ')
disp('COMPONENTS OF VAR(Q)/VAR(Y) CALCULATION')
dispvq(T4(:,[ibase,25:33])',cases)

diary
