function[ys,check]=BFS_model_steadystate(ys,exe)
global M_
 
NumberOfParameters = M_.param_nbr;                            % Number of deep parameters.
for i = 1:NumberOfParameters                                  % Loop...
  paramname = deblank(M_.param_names(i,:));                   %    Get the name of parameter i. 
  eval([ paramname ' = M_.params(' int2str(i) ');']);         %    Get the value of parameter i.
end                                                           % End of the loop.  
check = 0;

x0 = [1.3;0.95;1.3]; 
options = optimset('Display','off');
[x,fval] = fsolve(@ss_fun,x0,options);

H_borr_bar = x(1);
tauH = x(2);
H_bar = x(3);

C_borr_bar=gama/(1-epsilon_borr)*W_bar*H_borr_bar^(-phi);
D_borr_bar=C_borr_bar/upsilon_borr;
L_borr_bar_=Pi_bar*omega_bar*D_borr_bar/RL_bar;
H_s_borr_bar=H_borr_bar;
C_s_borr_bar=C_borr_bar;
D_s_borr_bar=C_s_borr_bar/upsilon_borr_s;
L_s_borr_bar=Pi_bar*omega_bar*D_s_borr_bar/RL_s_bar;
Inv_borr_bar = deltta*D_borr_bar;
LF_by_K = 1/K_by_LF;
H_s_bar = H_bar;
C_bar = gama/(1-epsilon)*W_bar*H_bar^(-phi); 
C_s_bar = C_bar;
D_bar=C_bar/upsilon;
D_s_bar=C_s_bar/upsilon_s;
C_tot_bar = lambda*C_bar+(1-lambda)*C_borr_bar;
Inv_bar = deltta*D_bar;
C_H_bar = tauC*(lambda*C_bar+(1-lambda)*C_borr_bar);
Hc_borr_bar = tauH*H_borr_bar;
Hd_borr_bar = (1-tauH)*H_borr_bar;
Hc_s_borr_bar = tauH*H_s_borr_bar;
Hd_s_borr_bar = (1-tauH)*H_s_borr_bar;
Hc_bar = tauH*H_bar;
Hd_bar = (1-tauH)*H_bar;
Hc_s_bar = tauH*H_s_bar;
Hd_s_bar = (1-tauH)*H_s_bar;
Hc_tot_bar = lambda*Hc_bar + (1-lambda)*Hc_borr_bar;
Hd_tot_bar = lambda*Hd_bar + (1-lambda)*Hd_borr_bar;
Hc_s_tot_bar = lambda*Hc_s_bar + (1-lambda)*Hc_s_borr_bar;
Hd_s_tot_bar = lambda*Hd_s_bar + (1-lambda)*Hd_s_borr_bar;
Kc_bar = K_by_H*Hc_tot_bar;
Kd_bar = K_by_H*Hd_tot_bar;
Kc_s_bar = K_by_H_s*Hc_s_tot_bar;
Kd_s_bar = K_by_H_s*Hd_s_tot_bar;
K_bar = Kc_bar + Kd_bar;
K_s_bar = Kc_s_bar + Kd_s_bar;
Ycw_bar = K_by_H^alppha*Hc_tot_bar;
Ydw_bar = K_by_H^alppha*Hd_tot_bar;
Ycw_s_bar = K_by_H_s^alppha*Hc_s_tot_bar;
Ydw_s_bar = K_by_H_s^alppha*Hd_s_tot_bar;

Yc_bar = Ycw_bar/dispc_bar;
Yd_bar = Ydw_bar/dispd_bar;
Yc_s_bar = Ycw_s_bar/dispc_s_bar;
Yd_s_bar = Ydw_s_bar/dispd_s_bar;

Y_bar = Yc_bar + Yd_bar;
Y_s_bar = Yc_s_bar + Yd_s_bar;
Y_emu_bar = n*Y_bar +(1-n)*Y_s_bar;
Gov_bar = g*Yc_bar;
Gov_s_bar = g_s*Yc_s_bar;
LF_bar = LF_by_K*K_bar;
N_bar = K_bar-LF_bar;
LF_s_bar = LF_bar;
C_exit_bar = Theta*(1-gammaB_bar)*N_bar/(gammaB_bar);
varrhoF1_bar = varrhoF1_by_LF*LF_bar;
varrhoF2_bar = varrhoF2_by_LF*LF_bar;
tax_exiting_entrepreneur_bar = (1-Theta)*(1-gammaB_bar)*N_bar/(gammaB_bar);
T_borr_bar = (Gov_bar -  tax_exiting_entrepreneurs_on*tax_exiting_entrepreneur_bar)/(1-lambda); 
Bb_bar_ = LF_bar + (1-lambda)*L_borr_bar_;
Bb_s_bar = Bb_bar_;
Kb_bar = Kb_by_Bb*Bb_bar_;
bank_profit = Rb_bar*Bb_bar_-RR_bar*(Bb_bar_-Kb_bar)-.5*kappaK*(Kb_by_Bb-etab_bar)^2*Kb_bar;
Ib_bar = varphib*bank_profit;
I_bar = deltaK*K_bar;
I_s_bar = I_bar;
H_total_bar=lambda*H_bar+(1-lambda)*H_borr_bar;
H_s_total_bar=H_total_bar;
UC_bar = gama/(C_bar*(1-epsilon));
UC_s_bar = gamma_s/(C_s_bar*(1-epsilon)); 
spread_R_Rb_bar = Rb_bar-RR_bar;
rl = log(RL_bar);
rl_s = log(RL_s_bar);

betta = betta_bar;
beta_borr = beta_borr_bar;
q = log(1);
c = log(C_bar);
d = log(D_bar);
c_borr = log(C_borr_bar);
d_borr = log(D_borr_bar);
dpd = log(Pi_bar);
dpc = log(Pi_bar);
dph = log(Pi_bar);
hc = log(Hc_bar);
hd = log(Hd_bar);
hc_borr = log(Hc_borr_bar);
hd_borr = log(Hd_borr_bar);
varrho = log(gama/(betta*(1-epsilon)*C_bar));
inv = log(deltta*D_bar);
varrho_borr = log(gama/(beta_borr*(1-epsilon_borr)*C_borr_bar));
inv_borr = log(deltta*D_borr_bar);
lead_Varrho_borr=exp(varrho_borr);
lead_Varrho=exp(varrho);
l_borr = log(L_borr_bar_);

q_s = log(1);
c_s = log(C_s_bar);
d_s = log(D_s_bar);
c_s_borr = log(C_s_borr_bar);
d_s_borr = log(D_s_borr_bar);
dpd_s = log(Pi_bar);
dpc_s = log(Pi_bar);
dpf = log(Pi_bar);
hc_s = log(Hc_s_bar);
hd_s = log(Hd_s_bar);
hc_s_borr = log(Hc_s_borr_bar);
hd_s_borr = log(Hd_s_borr_bar);
varrho_s = log(gamma_s/(betta*(1-epsilon)*C_s_bar));
inv_s = log(deltta*D_s_bar);
varrho_s_borr = log(gamma_s/(beta_borr*(1-epsilon_borr)*C_s_borr_bar));
inv_s_borr = log(deltta*D_s_borr_bar);
lead_Varrho_borr_s=exp(varrho_s_borr);
lead_Varrho_s=exp(varrho_s);
l_s_borr = log(L_s_borr_bar);

lag_Inv_ = deltta*D_bar;
lag_Inv_s_ = deltta*D_s_bar;
lag_Inv_borr_ = deltta*D_borr_bar;
lag_Inv_borr_s_ = deltta*D_s_borr_bar;

r = log(RR_bar);
omegaa = log(1-chi);
omegaa_s = log(1-chi);
omegap = log(1-chi);
omegap_s = log(1-chi);
terms = log(1);
b = log(1);

k = log(K_bar);
lag_K_ = K_bar;
kd = log(Kd_bar);
k_s = log(K_s_bar);
lag_K_s_ = K_s_bar;
kd_s = log(Kd_s_bar);
tobinq = log(tobinQ_bar);
tobinq_s = log(tobinQ_s_bar);
rlf_h = log(RLF_H_bar);
rlf_h_s = log(RLF_H_s_bar);
rlf_f = log(RLF_F_bar);
rlf_f_s = log(RLF_F_s_bar);
lf = log(LF_bar);
lf_s = log(LF_s_bar);
varrhoF2 = varrhoF2_bar;
varrhoF2_s = varrhoF2_bar;
I = I_bar;
I_s = I_s_bar;
omegaFa = log(OmegaF_bar);
omegaFp = log(OmegaF_bar);
omegaFa_s = log(OmegaF_bar);
omegaFp_s = log(OmegaF_bar);

riskFa = log(sigmaF_bar);
riskFa_s = log(sigmaF_bar);
riskFp = log(sigmaF_bar);
riskFp_s = log(sigmaF_bar);
riskp = log(sigma_bar);
riskp_s = log(sigma_bar);
riska = log(sigma_bar);
riska_s = log(sigma_bar);

phtilde = log(Phtilde_bar);
pftilde = log(Pftilde_bar);
pdtilde = log(Pdtilde_bar);
pdtilde_s = log(Pdtilde_s_bar);
vc = log(dispc_bar);
vc_s = log(dispc_s_bar);
vd = log(dispd_bar);
vd_s = log(dispd_s_bar);

xc1 = log(UC_bar*MCc_bar*Ycw_bar/(1-theta_c*betta*Pi_bar^(1+sigmma))); 
xc2 = log(UC_bar*Ycw_bar/(1-theta_c*betta*Pi_bar^(sigmma)));
xc1_s = log(UC_s_bar*MCc_s_bar*Ycw_s_bar/(1-theta_c_s*betta*Pi_bar^(1+sigmma))); 
xc2_s = log(UC_s_bar*Ycw_s_bar/(1-theta_c_s*betta*Pi_bar^(sigmma))); 
xd1 = log(UC_bar*MCd_bar*Ydw_bar/(1-theta_d*betta*Pi_bar^(1+sigmma))); 
xd2 = log(UC_bar*Ydw_bar/(1-theta_d*betta*Pi_bar^(sigmma)));
xd1_s = log(UC_s_bar*MCd_s_bar*Ydw_s_bar/(1-theta_d_s*betta*Pi_bar^(1+sigmma)));
xd2_s = log(UC_s_bar*Ydw_s_bar/(1-theta_d_s*betta*Pi_bar^(sigmma)));

prefd = log(1);
prefd_s = log(1);
premium = log(1);
prefc = log(1);
prefc_s = log(1);
zc = log(1);
zd = log(1);
zc_s = log(1);
zd_s = log(1);

kb = log(Kb_bar);
kb_s = log(Kb_bar);
lag_Kb_ = Kb_bar;
lag_Kb_s_ = Kb_bar;

lead_LAMBDA = betta;
lead_LAMBDA_s = betta;

A = 1;
Gov = g*Yc_bar;
Gov_s = g_s*Yc_s_bar;
gammaB = gammaB_bar;
etab = etab_bar;
etab_s = etab_s_bar;


Y_emu_a = Y_emu_bar;
R_a = RR_bar;
R_ib_a = RR_bar;
B_a = log(1);
Pi_emu_a = Pi_bar;
C_tot_a = C_tot_bar;
I_a = I_bar;
H_tot_a = H_total_bar;
K_by_B_a = Kb_by_Bb;
Bb_a = Bb_bar_;
Loans_to_firms = LF_bar;
Prop_Loans_to_Firms = LF_bar/(LF_bar+(1-lambda)*L_borr_bar_);
Prop_business_loans_crossborder = 1-tauL;
TPS = 1;
wages = W_bar;
return_capital = RK_bar;
capital_rental_rate = rK_bar;
bank_capital_investment =Ib_bar;
res_inv = Yd_bar;
net_worth = N_bar;
tobQ = tobinQ_bar;
marg_cost_C = MC_bar;
marg_cost_D = MC_bar;
spread_R_Rb = spread_R_Rb_bar;
mortgages = L_borr_bar_;

 Xc1 = exp(xc1);
 Xc2 = exp(xc2);
 Xc1_s = exp(xc1_s);
 Xc2_s = exp(xc2_s);
 Xd2 = exp(xd2);
 Xd1_s = exp(xd1_s);
 Xd2_s = exp(xd2_s);
fprintf('\n Residual of Equation 65: %d',sigmma/(sigmma-1) * Xc1 - Xc2 * exp(phtilde))
    
NumberOfEndogenousVariables = M_.endo_nbr;                    % Number of endogenous variables.
ys = zeros(NumberOfEndogenousVariables,1);                    % Initialization of ys (steady state).
for i = 1:NumberOfEndogenousVariables                         % Loop...
  
  varname = deblank(M_.endo_names(i,:));                      %    Get the name of endogenous variable i.                     
  eval(['ys(' int2str(i) ') = ' varname ';']);                %    Get the steady state vZNue of this variable.
  
end