Two country, two goods model

Steady state
1

Can anyone help me debug this? I have been staring at this for days and I still couldn’t figure out why would my Steady State be wrong.

var c cstar n nstar k kstar z zstar r rstar w wstar y ystar B Bstar q sigma sigmastar p_a p_b pstar_a pstar_b a astar b bstar;
predetermined_variables k kstar B Bstar;
varexo e_z estar_z e_sigma estar_sigma;
parameters beta theta delta rho rhostar A Bbar pac omega omegastar alpha domshare domsharestar;

alpha = 0.67;
beta = 0.9895;
delta = 0.0255;
theta = 0.36;
A = 2;
rho = 0.98;
rhostar = 0.98;
Bbar = 0;
pac = 0.00074;
omega = 0.98;
omegastar = 0.98;
domshare = 0.5;
domsharestar = 0.5;

model;
/*
Home Economy
*/
y = exp(z) * k^theta * n^(1 - theta);
y = a + astar;

/*
Constraints for home Economy
*/
c + k(+1) + q * B(+1) + 0.5 * pac * (B(+1) - Bbar)^2 = w * n + r * k + (1 - delta) * k + B;

/*
FOC - intra temporal, capital accumulation, and bond accumulation
*/
A/(1-n) = (1/c) * w;
(1/c) = beta * ((1/c(+1)) * (r(+1) + 1 - delta));
beta * (1/c(+1)) = (q + pac * (B(+1) - Bbar)) * (1/c);

/*
Prices
*/
r = theta * exp(z) * k^(theta - 1) * n^(1-theta);
w = (1 - theta) * exp(z) * k^(theta) * n^(-theta);
p_a = (domshare * a^(1 - alpha) + (1 - domshare) * b^(1 - alpha))^(alpha / (1 - alpha)) * (domshare * a^(-alpha));
p_b = (domshare * a^(1 - alpha) + (1 - domshare) * b^(1 - alpha))^(alpha / (1 - alpha)) * ((1 - domshare) * b^(-alpha));

/*
TFP
*/
z = rho * z(-1) + exp(sigma) * e_z;
sigma = omega * sigma(-1) + e_sigma;

/*
Foreign Economy
*/
ystar = exp(zstar) * kstar^theta * nstar^(1 - theta);
ystar = bstar + b;

/*
Constraints for foreign Economy
*/
cstar + kstar(+1) + q * Bstar(+1) + 0.5 * pac * (Bstar(+1) - Bbar)^2 = wstar * nstar + rstar * kstar + (1 - delta) * kstar + Bstar;

/*
FOC - intra temporal, capital accumulation, and bond accumulation
*/
A/(1-nstar) = (1/cstar) * wstar;
(1/cstar) = beta * ((1/cstar(+1)) * (rstar(+1) + 1 - delta));
beta * (1/cstar(+1)) = (q + pac * (Bstar(+1) - Bbar)) * (1/cstar);

/*
Prices
*/
rstar = theta * exp(zstar) * kstar^(theta - 1) * nstar^(1-theta);
wstar = (1 - theta) * exp(zstar) * kstar^(theta) * nstar^(-theta);
pstar_b = (domsharestar * bstar^(1 - alpha) + (1 - domsharestar) * astar^(1 - alpha))^(alpha / (1 - alpha)) * (domsharestar * bstar^(-alpha));
pstar_a = (domsharestar * bstar^(1 - alpha) + (1 - domsharestar) * astar^(1 - alpha))^(alpha / (1 - alpha)) * ((1 - domsharestar) * astar^(-alpha));

/*
TFP
*/
zstar = rhostar * zstar(-1) + exp(sigmastar) * estar_z;
sigmastar = omegastar * sigmastar(-1) + estar_sigma;

/*
Market clearing
*/
B(+1) + Bstar(+1) = 0;
(domshare * a^(1 - alpha) + (1 - domshare) * b^(1 - alpha))^(1 / (1 - alpha)) = c + k(+1) - (1 - delta) * k;
(domsharestar * bstar^(1 - alpha) + (1 - domsharestar) * astar^(1 - alpha))^(1 / (1 - alpha)) = cstar + kstar(+1) - (1 - delta) * kstar;
end;

steady_state_model;
z = 0;
zstar = 0;
sigma = 0;
sigmastar = 0;

r = (1/beta) + delta - 1;
n = 1 / (1 + (A/(1-theta)) * (1 - delta * theta * (beta/(1 + beta * delta - beta))));
k = (r / (theta * n^(1 - theta)))^(1/(theta-1));
y = k^theta * n^(1 - theta);
w = (1 - theta) * (y / n);
c = w * n + k * (r - delta);

B = 0;
Bstar = 0;

q = beta - pac * (B - Bbar);

rstar = (1/beta) + delta - 1;
nstar = 1 / (1 + (A/(1-theta)) * (1 - delta * theta * (beta/(1 + beta * delta - beta))));
kstar = (rstar / (theta * nstar^(1 - theta)))^(1/(theta-1));
ystar = kstar^theta * nstar^(1 - theta);
wstar = (1 - theta) * (ystar / nstar);
cstar = wstar * nstar + kstar * (rstar - delta);

a = (y * (domshare / (1 - domshare))^(1 / alpha)) / (1 + ( (domshare / (1 - domshare))^(1 / alpha) ));
astar = y - a;

bstar = (ystar * (domsharestar / (1 - domsharestar))^(1 / alpha)) / (1 + ( (domsharestar / (1 - domsharestar))^(1 / alpha) ));
b = ystar - bstar;



p_a = (domshare * a^(1 - alpha) + (1 - domshare) * b^(1 - alpha))^(alpha / (1 - alpha)) * (domshare * a^(-alpha));
p_b = (domshare * a^(1 - alpha) + (1 - domshare) * b^(1 - alpha))^(alpha / (1 - alpha)) * ((1 - domshare) * b^(-alpha));
pstar_a = (domsharestar * bstar^(1 - alpha) + (1 - domsharestar) * astar^(1 - alpha))^(alpha / (1 - alpha)) * ((1 - domsharestar) * astar^(-alpha));
pstar_b = (domsharestar * bstar^(1 - alpha) + (1 - domsharestar) * astar^(1 - alpha))^(alpha / (1 - alpha)) * (domsharestar * bstar^(-alpha));




end;

model_diagnostics;
resid;

steady;


/*
shocks;
var e_z = 0.01^2;
var estar_z = 0.01^2;

var e_sigma = 0.048^2;
var estar_sigma = 0.048^2;
end;

stoch_simul(order=3, pruning, irf=300, periods=300);

//stoch_simul(order=3, pruning, irf=300, periods=300, replic=20000, simul_replic=20000);
simulated_series_raw=get_simul_replications(M_,options_);
*/
2

That often means that the underlying approach to computing the steady state is wrong. Often, one can approach it from the production and the demand side. Both approaches should deliver the same result. Sometimes they don’t, pointing to the error.