function res=steady_state(X)
global MULT Y K K_L C C_Y K_Y Lh L PHI Kh Ks Pstar_P Wh_P Wstar_P W_P R Rk_P gammaw_1 gammaw_2 MC_P sp sw pi gammap_1 gammap_2 r premium G B_P C1 Lh1 C1_Y C2 Lh2 C2_Y;

MULT=X(1);
Y=X(2);
K=X(3);
K_L=X(4);
C=X(5);
C_Y=X(6);
K_Y=X(7);
Lh=X(8);
L=X(9);
PHI=X(10); 
Kh=X(11);
Ks=X(12);
Pstar_P=X(13);
Wh_P=X(14);
Wstar_P=X(15);
W_P=X(16);
R=X(17);
Rk_P=X(18);
premium=X(19);
gammaw_1=X(20);
gammaw_2=X(21); 
MC_P=X(22);
r=X(23);
gammap_1=X(24);
gammap_2=X(25);
G = X(26);
B_P = X(27);
C1 = X(28);
Lh1 = X(29);
C1_Y = X(30);
C2 = X(31);
Lh2 = X(32);
C2_Y = X(33);

sp=1;
sw=1;
pi=1;
delta=0.025;
lambda_w=0.5;
lambda_p=0.62;
sigma_c=1.39;
lambda=0.71;
sigma_l=1.92;
theta_p=0.65;
theta_w=0.73;
PHI_w=1.62;
alpha=0.19;
beta=0.998;
tau_w = 0.3;
% rho_R=0.8;
% phi_pi=1.5;
% phi_y=0.5;
% rho_a=0.95;
% rho_b=0.18;
% rho_r=0.12;
% sigma_a=0.45;
% sigma_b=0.23;
% sigma_r=0.24;
G_ss = 0.1678;
nu = 0.1;

%% Steady State

res(1,1) = PHI-(Y*(PHI_w -1));
res(2,1) = r-(log(R)-log(pi));   % if prices stay the same, then there is not going to be inflation!                 
res(3,1) = premium-(Rk_P-r);   % differential in the expected real net returns
res(4,1) = Lh-L;
res(5,1) = R-beta^(-1);
res(6,1) = MC_P-(1+lambda_p)^(-1);
res(7,1) = Rk_P-(1/beta-1+delta);
res(8,1) = W_P-(MC_P*(1-alpha)*((alpha/(1-alpha))^(alpha))*Rk_P^(-alpha))^(1/(1-alpha));
res(9,1) = K_L-(alpha/(1-alpha))*W_P*Rk_P^(-1);
res(10,1) = K_Y-K_L^(1-alpha)*PHI_w;
res(11,1) = C_Y-(1-delta*K_Y) + G/Y;
res(12,1) = Wh_P-(1+lambda_w)^(-1)*W_P;
res(13,1) = Lh1-((1-lambda)*(C1_Y*K_Y^(-1)*K_L*L/Lh1)/((1-tau_w)*Wh_P))^(-1/(1+sigma_l));
res(14,1) = Y-(K_L^alpha)*L-PHI;
res(15,1) = C-C_Y*Y;
res(16,1) = K-K_L*L;
res(17,1) = gammaw_1-W_P^((1+lambda_w)/lambda_w)*L/(1-beta*theta_w);
res(18,1) = gammaw_2-W_P^((1+lambda_w)/lambda_w)*L*Wh_P/(1-beta*theta_w);
res(19,1) = gammap_1-Y/(1-beta*theta_p);
res(20,1) = Pstar_P-1;
res(21,1) = Wstar_P-W_P;
res(22,1) = MULT-1/C1;
res(23,1) = gammap_2-Y/(1-beta*theta_p)*MC_P;
res(24,1) = Kh-K;
res(25,1) = Ks-K;
res(26,1) = G + B_P*(1/pi - 1/R) - tau_w*Wh_P*Lh;
res(27,1) = G - G_ss;
res(28,1) = C1-C1_Y*Y;
res(29,1) = C- (1-nu)*C1 - nu*C2;
res(30,1) = Lh - (1-nu)*Lh1 - nu*Lh2;
res(31,1) = Lh2-((1-lambda)*(C2_Y*K_Y^(-1)*K_L*L/Lh2)/((1-tau_w)*Wh_P))^(-1/(1+sigma_l));
res(32,1) = C2 - (1-tau_w)*Wh_P*Lh2;
res(33,1) = C2-C2_Y*Y;

res=res';



% X = [MULT, Y, K, K_L, C, C_Y, K_Y, Lh, L, PHI, Kh, Ks, Pstar_P, Wh_P, Wstar_P, W_P, R, Rk_P, pi, gammaw_1, gammaw_2,...
%     sw MC_P, gammap_1, gammap_2, sp, r, premium]';
end