function [ys,check] = Thesis_nlnbl_steadystate(ys,exo)
% function [ys,check] = Thesis_nlnbl_steadystate(ys,exo)
% computes the steady state for the NK_baseline.mod and uses a numerical
% solver to do so
% Inputs: 
%   - ys        [vector] vector of initial values for the steady state of
%                   the endogenous variables
%   - exo       [vector] vector of values for the exogenous variables
%
% Output: 
%   - ys        [vector] vector of steady state values fpr the the endogenous variables
%   - check     [scalar] set to 0 if steady state computation worked and to
%                    1 of not (allows to impos restriction on parameters)

global M_ 

% read out parameters to access them with their name
NumberOfParameters = M_.param_nbr;
for ii = 1:NumberOfParameters
  paramname = deblank(M_.param_names(ii,:));
  eval([ paramname ' = M_.params(' int2str(ii) ');']);
end
% initialize indicator
check = 0;


%% Enter model equations here

options=optimset(); % set options for numerical solver

% 1. Steady state policy interest rate
r_ib=(1+pi_bar)/betta-1;

% Set parameter r_ib_ss
r_ib_bar = r_ib;

% 10-13. Interest rates
r_f    = rstar+pbar;
R_b    = r_ib;
r_b    = R_b+muy_b;

% exchange rate depreciation
ge = (1+r_ib)/(1+r_f)-1;

% 4. FOC investment
q      = 1;

% exogenous processes
% 33. Technology Shock
a      = 1;
xia    = 0;
% 34. Foreign output Shock
ystar  = 1;
xiy    = 0;

%Set parameter ystar_ss
ystar_ss = ystar;

% 35. Monetary shock
u      = 1;
xiu    = 0;

% 22. CPI inflation
pi_cpi = pi_bar;

% 23-24. Domestic goods & imported goods inflation
pi_h   = pi_bar;
pi_f   = pi_bar;
pih   = pi_bar;
pif   = pi_bar;

%8-9. Domestic & Imported good optimal price
pi_hr   = (((1+pi_cpi)^(1-epsl)-thetta_h)/(1-thetta_h))^(1/(1-epsl));
pi_fr   = (((1+pi_cpi)^(1-epsl)-thetta_f)/(1-thetta_f))^(1/(1-epsl));

%6-7. Domestic & import price dispersion
vh = ((1-thetta_h)*((1+pi_h)/(1+pi_hr))^epsl)/(1-thetta_h*(1+pi_h)^(epsl));
vf = ((1-thetta_f)*((1+pi_f)/(1+pi_fr))^epsl)/(1-thetta_f*(1+pi_f)^(epsl));

% Real marginal cost 
mc = (1-betta*thetta_h*(1+pi_h)^(epsl))*(1+pi_hr)*(epsl-1)/((1-betta*thetta_h*(1+pi_h)^(epsl-1))*(1+pi_h)*epsl);

% lopg
lopg = (1-betta*thetta_f*(1+pi_f)^(epsl))*(1+pi_fr)*(epsl-1)/((1-betta*thetta_f*(1+pi_f)^(epsl-1))*(1+pi_f)*epsl);

% Nominal wage
w      = (alppha/(r_b+delt_k))^(alppha/(1-alppha));

% Foreign demand for domestic goods
c_hstar=gam*lopg^etta*ystar;

% {solve numerically for n}
[n,fval,exitflag]=fzero(@(n)n^(varrphi+alppha)/w^alppha*gam-n*((w/(r_b+delt_k))*(alppha/(1-alppha)))^alppha+c_hstar,0.6,options);
if exitflag <1
    %indicate the SS computation was not sucessful; this would also be detected by Dynare
    %setting the indicator here shows how to use this functionality to
    %filter out parameter draws
    check=1; %set failure indicator
    return; %return without updating steady states
end


% capital
k      = n*(w/(r_b+delt_k))*(alppha/(1-alppha));

% 26. Production function
y      = a*k^alppha*n^(1-alppha);

% 16. Local demand for domestic goods
c_h     = y-c_hstar;

% Total consumption
c = c_h/gam;

% Local demand for imported goods
c_f     = (1-gam)*c;

% 27. Law of motion for capital of domestic producers
invest  = delt_k*k;

% 28. New investment financed by bank loans (credit market clearing) CHECK if necessary
b      = k;

% Leverage
k_b    = b*nuy;

% 29. Law of motion for capital of banks
PI_b   = k_b*delt_b;

% 30. Bank balance sheet
m      = b*(1-nuy);

% Bank profit
j_b = delt_b*k_b;

% Debt
d = r_ib*(1-nuy)*k+(1+pi_bar)*c_hstar-mc*(1+pi_bar)*y*vh-lopg*(1+pi_bar)*c_f*vf;

% Set parameter d_ss
d_bar   = d;

x1=c^(-siggma)*y*mc/(1-betta*thetta_h*(1+pi_bar)^epsl);

x2=c^(-siggma)*y/(1-betta*thetta_h*(1+pi_bar)^(epsl-1));

g1=c^(-siggma)*c_f*lopg/(1-betta*thetta_f*(1+pi_bar)^epsl);

g2=c^(-siggma)*c_f/(1-betta*thetta_f*(1+pi_bar)^(epsl-1));


check=0;

ys = [log(y);
      log(k);
      log(invest);
      log(c);
      log(c_h);
      log(c_f);
      log(c_hstar);
      log(vh);
      log(vf);
      d;
      log(n);
      log(w);
      log(mc);
      pi_h;
      pi_hr;
      log(x1);
      log(x2);
      log(lopg);
      pi_f;
      pi_fr;
      log(g1);
      log(g2);
      pi_cpi;
      log(pih);
      log(pif);
      ge;
      log(q);
      r_b;
      R_b;
      r_ib;
      r_f;
      log(b);
      log(k_b);
      log(m);
      log(j_b);
      log(a);
      log(ystar);
      log(u)];