%% runsim
%%
format short

restoredefaultpath

setpathdynare4

clear
set(0,'DefaultLineLineWidth',2)
close all

paramfile_baby00


for indi=1:2
    
    
    
    STD_J      =0.0737;
    STD_K      =0.0360;
    STD_P      =0.0030;
    STD_R      =0.0013;
    STD_W      =0.0100;
    STD_Z      =0.0163;
    
    
    
    modnam_00 = 'baby00'; % base model (constraint 1 and 2 below don't bind)
    modnam_10 = 'baby10'; % first constraint is true
    modnam_01 = 'baby01'; % second constraint is true
    modnam_11 = 'baby11'; % both constraints bind
    
   % compile constraint
    constraint1 = 'lm<-lm_ss';
    constraint_relax1 = 'maxlev > 0';
    constraint2 = 'rnot < -log(r_ss)';
    constraint_relax2 = 'r > -log(r_ss)';
    
    
    curb_retrench = 0;           % if 1 slow down relaxation of constraints
    maxiter = 40;                % number of iterations a5llowed to look for a solution
    
    irfshock =char('eps_j','eps_k','eps_p','eps_r','eps_w','eps_z', 'eps_g');

   if indi==1  
     sequence1 = [ zeros(1,7)
     0.46 0 0 0 0 0 0   
     0 0 0 0 0 0 0
     0 0 0 0 0 0 0 
     0 0 0 0 0 0 0 
     0 0 0 0 0 0 0 
     0 0 0 0 0 0 0 
     0 0 0 0 0 0 0 
     0 0 0 0 0 0 0
   zeros(35,7) ];
   end    
        

  if indi==2  
    sequence1 = [ zeros(1,7)
    -0.46 0 0 0 0 0 0 
    0 0 0 0 0 0 0
    0 0 0 0 0 0 0
    0 0 0 0 0 0 0
    0 0 0 0 0 0 0
    0 0 0 0 0 0 0
    0 0 0 0 0 0 0
    0 0 0 0 0 0 0
  zeros(35,7) ];
  end
     
  nperiods = 44;  
    
    [zdatal zdatap zdatass oobase_ Mbase_] = solve_two_constraints(...
        modnam_00,modnam_10,modnam_01,modnam_11,...
        constraint1, constraint2,...
        constraint_relax1, constraint_relax2,...
        sequence1,irfshock,nperiods+50,curb_retrench,maxiter);
    
       
    
    for i=1:Mbase_.endo_nbr
        eval([deblank(Mbase_.endo_names(i,:)),'_l=zdatal(1:nperiods,i);']);
        eval([deblank(Mbase_.endo_names(i,:)),'_p=zdatap(1:nperiods,i);']);
        eval([deblank(Mbase_.endo_names(i,:)),'_ss=zdatass(i);']);
    end
    
    for i = 1:size(Mbase_.param_names,1)
        eval([Mbase_.param_names(i,:),'= Mbase_.params(i);']);
    end
    
    ctot_ss = c_ss + c1_ss;
    ctot_l = c_l + c1_l;
    ctot_p = c_p + c1_p;
    
    
    titlelist = char(...
        'Interest Rate',...
        'GDP',...
        'House Prices',...
        'Total Consumption',...
        'Consumption Saver',...
        'Consumption Borrowers',...
        'Investment',...
        'Multiplier (level)',...
        'Hours',...
        'Houses Borrowers');
    ylabels = char('Percentage Point, AR','Percent','Percent',...
        'Percent','Percent','Percent',...
        'Percent','Level','Percent', 'Percent');
    
    
    
    
    line1(:,:,indi)=100*[
        4*(r_p+r_ss-1),...
        y_p/y_ss,...
        q_p/q_ss,...
        ctot_p/ctot_ss,...
        c_p/c_ss,...
        c1_p/c1_ss,...
        ik_p/ik_ss,...
        (lm_p+lm_ss),...
        h1_p/h1_ss,...
        g_p/g_ss,...
        n_p/n_ss,...
        n1_p/n1_ss];
    line2(:,:,indi)=100*[
        4*(r_l+r_ss-1),...
        y_l/y_ss,...
        q_l/q_ss,...
        ctot_l/ctot_ss,...
        c_l/c_ss,...
        c1_l/c1_ss,...
        ik_l/ik_ss,...
        (lm_l+lm_ss),...
        h1_l/h1_ss,...
        g_l/g_ss,...
        n_l/n_ss,...
        n1_l/n1_ss];
    
    
    
    
    
    
    
    
    
    
    
    
end

%%
set(0,'defaulttextinterpreter','latex')
set(0,'DefaultAxesFontSize', 16)
set(0,'DefaultTextFontSize', 16)

nperiods = 20;

figure (1)

subplot(3,2,1)
plot(squeeze(line1(:,3,1)),'b');  hold on
plot(-squeeze(line1(:,3,2)),'r--'); 
xlim([0 nperiods])
ylim([0 15])
grid on
title({'House Price'})
ylabel('Pct. from st.st.','Rotation',90)
hold on; plot(0*r_l,'k','Linewidth',1)
legend1 = legend('Positive shock','Negative shock');
set(legend1,'Position',[0.5,0.33,0,0],'Orientation','horizontal','Interpreter','latex','Fontsize',16)


subplot(3,2,2)
plot(squeeze(line1(:,4,1)),'b');  hold on
plot(-squeeze(line1(:,4,2)),'r--'); 
xlim([0 nperiods])
ylim([0 0.8])
grid on
title({'Consumption (total)'})
ylabel('Pct. from st.st.','Rotation',90)
hold on; plot(0*r_l,'k','Linewidth',1)

subplot(3,2,3)
plot(squeeze(line1(:,8,1)),'b');  hold on
plot(squeeze(line1(:,8,2)),'r--'); 
xlim([0 nperiods])
ylim([-4 4])
grid on
title({'Borrowing Multiplier'})
ylabel('Level','Rotation',90)
hold on; plot(0*r_l,'k','Linewidth',1)

subplot(3,2,4)
plot(squeeze(line1(:,6,1)),'b');  hold on
plot(-squeeze(line1(:,6,2)),'r--'); 
xlim([0 nperiods])
ylim([0 6])
grid on
title({'Consumption (impatient)'})
ylabel('Pct. from st.st.','Rotation',90)
hold on; plot(0*r_l,'k','Linewidth',1)

printFile2(['graphs_runsim_IRFs\Fig3'],14,19)

