%% Main code is chapter1v64.mod

function res=chapter2_fsteady(xx)

global p_ss rn_ss zave_ss xi_ss varphistar_ss ...
       w_ss varphiave_ss l_ss c_ss cn_ss cm_ss yn_ss kn_ss ...
       in_ss shn_ss h_ss he_ss kx_ss yx_ss ...
       kxave_ss lxave_ss yxave_ss;

global An phi0 sigma r am;
global alphan alpha nu phi Ome epsilon etay etae eta lambda rho;
global pmk_ss pmc_ss px_ss pn_ss;

[rxx, cxx] = size(xx);

for j = 1:rxx
    if xx(j) < 0
        for j = 1:rxx
            res(j) = 1000000;
        end
        return
    end
end

if xx(1)>1 | xx(2)>1 | xx(4)>1 | xx(5)>1 
    for j = 1:rxx
        res(j)=1000000;
    end
    return
end
 
delta = xx(1);
kappa = xx(2);
varsigma = xx(3);
l_ss = xx(4);
shn_ss = xx(5);

%--------------

p_ss = (((1-am)*(pn_ss^(1-rho)))+(am*(pmc_ss^(1-rho))))^(1/(1-rho));
rn_ss = pn_ss*(r+delta);
zave_ss = (phi/(epsilon-phi))*(px_ss*etay+pn_ss*eta);
xi_ss = (etae/zave_ss)*(r+lambda)*pn_ss;
varphistar_ss = phi0*(xi_ss^(1/(-epsilon)));
varphiave_ss = ((epsilon/(epsilon-phi))^(1/phi))*varphistar_ss;

w_ss = ((px_ss*varphistar_ss/(pmk_ss^nu))*...
             ((Ome/(px_ss*etay+pn_ss*eta))^(1/phi)))^(1/alpha);
w1 = pn_ss*alphan*((An*(((1-alphan)/(r+delta))^(1-alphan)))^(1/alphan)); %CORRECTION!!

kn_ss = (shn_ss*l_ss)*(1-alphan)*w_ss/(alphan*rn_ss);
in_ss = delta*kn_ss;
yn_ss = An*((shn_ss*l_ss)^alphan)*(kn_ss^(1-alphan));

h_ss = ((1-shn_ss)*l_ss)*...
          ((xi_ss*((px_ss*varphiave_ss)^phi)*...
             ((alpha/w_ss)^((1-nu)*phi))*((nu/pmk_ss)^(nu*phi)))^(-1));
he_ss = lambda*h_ss;
kx_ss = ((1-shn_ss)*l_ss)*w_ss*nu/(alpha*pmk_ss);
yx_ss = ((w_ss*((1-shn_ss)*l_ss)) + (pmk_ss*kx_ss) + (xi_ss*zave_ss*h_ss) ...
             + (px_ss*etay*xi_ss*h_ss) + (pn_ss*eta*xi_ss*h_ss))/px_ss;

kxave_ss = kx_ss/(xi_ss*h_ss);
lxave_ss = ((1-shn_ss)*l_ss)/(xi_ss*h_ss);
yxave_ss = yx_ss/(xi_ss*h_ss);
        
        
%tb_ss = ((px_ss*yx_ss)-(px_ss*etay*xi_ss*h_ss)) - ((pmk_ss*kx_ss)+(pmc_ss*cm_ss));    

c_ss = (1-l_ss)*w_ss*kappa/(p_ss*(1-kappa));
cn_ss = (1-am)*((pn_ss/p_ss)^(-rho))*c_ss;
cm_ss = am*((pmc_ss/p_ss)^(-rho))*c_ss;

cn1 = yn_ss - in_ss - (etae*he_ss) - (xi_ss*h_ss*eta);
cm1 = (((px_ss*yx_ss)-(px_ss*etay*xi_ss*h_ss)) - (pmk_ss*kx_ss))/pmc_ss;

res = zeros(rxx,1);
res(1,1)  = w_ss - w1;
res(2,1)  = (p_ss*c_ss) - (pmc_ss*cm1) - (pn_ss*cn1);
res(3,1)  = cn_ss - cn1;
res(4,1)  = cm_ss - cm1;
res(5,1) = ((px_ss*yx_ss)-(px_ss*etay*xi_ss*h_ss)) - ((pmk_ss*kx_ss)+(pmc_ss*cm_ss));