% The MATLAB subroutine employ.m recalled by the function 
% fsolve that computes the steady-state of each variable is

function F = employ(n_ss, alpha, omega, beta, gamma, phi, rho, delta, a_ss, B, chi);
 k_ss   = ((1+delta*beta-beta)/(alpha*beta*a_ss*(n_ss^(1-alpha))))^(1/(alpha-1));
 Ub_ss  = 1 - (1-rho)*n_ss;
 y_ss   = a_ss*(k_ss^alpha)*(n_ss^(1-alpha));
 h_ss   = rho*n_ss;
 x_ss   = h_ss/Ub_ss;
 g_ss   = (a_ss^gamma)*B*(x_ss^omega);
 u_ss   = (1 - n_ss);
 in_ss  = delta*k_ss;
 c_ss   = y_ss-(g_ss*h_ss)-in_ss;
 F      = (1-alpha)*a_ss*k_ss^alpha*n_ss^(-alpha)-a_ss^gamma*B*(omega+1)*x_ss^omega-chi*c_ss*n_ss^phi+beta*a_ss^gamma*B*(1-rho)*((omega+1)*x_ss^omega-omega*x_ss^(omega+1));