function [ys,params,check] = Mendoza(ys,exo,M_,options_)
% function [ys,params,check] = NK_baseline_steadystate(ys,exo,M_,options_)
% 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
%   - M_        [structure] Dynare model structure
%   - options   [structure] Dynare options structure
%
% Output: 
%   - ys        [vector] vector of steady state values for the the endogenous variables
%   - params    [vector] vector of parameter values
%   - check     [scalar] set to 0 if steady state computation worked and to
%                    1 of not (allows to impose restrictions on parameters)

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


%% Enter model equations here

Z=1;
Z2=1;
Z3=1;
mu=0;
lambda(1-Z*R)=mu;


GDP=Z*(K^(1-aplhaa))*(N^aplhaa); % needs N
C_N = GDP-T_N;
ppi=GDP-W;
C_T=((1+tauu)*(-(N)^(deltaa-1))/omegaa*W*(omegaa*(C_T)^(-etaa)+(1-omegaa)*(C_N)^(-etaa))^(-(1/etaa)-1))^(1/(-etaa-1));   
                                                                                                                                                                                                                      
P_N = (omegaa/(1-omegaa))*(C_N/C_T)^(-etaa-1); 
P=P_N/((omegaa*(C_T)^(-etaa))+(1-omegaa)*(C_N)^(-etaa))^(-(1/etaa)-1);
C=((omegaa*(C_T)^(-etaa))+(1-omegaa)*(C_N)^(-etaa))^(-(1/etaa)-1); 

W=P_N*Z3*aplhaa*(K^(1-aplhaa))*(N^aplhaa);          
D=-(1+tauu)*C_T+P_N*(C_N)+Z2*Y_T+D(-1)*Z*R-T;
N=(ppi-D-((1-tauu)*(C_T + P_N*C_N)) + (Z*R*D)-T-(P_N*T_N))/W;  
 
%% end own model equations

params=NaN(NumberOfParameters,1);
for iter = 1:length(M_.params) %update parameters set in the file
  eval([ 'params(' num2str(iter) ') = ' M_.param_names{iter} ';' ])
end
 
NumberOfEndogenousVariables = M_.orig_endo_nbr; %auxiliary variables are set automatically
for ii = 1:NumberOfEndogenousVariables
  varname = M_.endo_names{ii};
  eval(['ys(' int2str(ii) ') = ' varname ';']);
end
