RBC code simulation error

[yazanjer]

Hello guys,

I am trying to simulate RBC model. I have the dynare code below:
var c k y i w Rent n r a;
varexo e;
parameters alpha beta delta sigma teta rho sigmae;
beta=0.99;
alpha=0.36;
sigma=2;
delta=0.025;
sigmae = 0.01;
rho=0.95;
teta=1;
model;
1/exp(c(-1))=beta*((1/(exp(c)))*((a*alpha*exp(k)^(alpha-1)*exp(n)^(1-alpha))+(1-delta)));
teta/(1-exp(n))=(1/exp(c))*(exp(k)^(alpha)*a*((1-alpha)*exp(n)^(-alpha)));
exp(k)=(a*exp(k(-1))^(alpha)*exp(n(-1))^(1-alpha))-exp(c)+((1-delta)*exp(k(-1)));
exp(y)=(a*exp(k(-1))^(alpha)*exp(n)^(1-alpha));
exp(i)=exp(y)-exp(c);
1/exp(c(-1))=beta*(1/exp(c))*(1+exp(r));
exp(w)=(a*exp(k)^(alpha)*((1-alpha)*exp(n)^(-alpha)));
exp(Rent)=(a*alpha*exp(k)^(alpha-1)*exp(n)^(1-alpha));
a=rho*a(-1)+e;
end;
initval;
w =log(1);
Rent =log(0.1);
r =log(0.002);
n =log(0.6);
i =log(0.6);
y =log(1.5);
c =log(1.5);
k =log(16);
end;
shocks;
var e=sigmae^2;
end;
steady;
stoch_simul(hp_filter=1600,order=1,irf=40);


I am getting the following error:

Error using lnsrch1 (line 53)
Some element of Newton direction isn't finite. Jacobian maybe
singular or there is a problem with initial values


Could you help me to find the error in my code.

Thanks

[jpfeifer]

There are a bunch of errors. Your TFP in the equations should be exp(a) as otherwise output in steady state is 0 as a has steady state 0. Moreover, your timing is wrong. For example, your Euler equation is not forward looking.

[yazanjer]

Thanks jpfeifer for your help, I considered your corrections and it worked. Do you think the following code is accepted as standard RBC model simulation.

var c k y i w Rent n r a;
varexo e;
parameters alpha beta delta sigma teta rho sigmae;
beta=0.99;
alpha=0.36;
sigma=2;
delta=0.025;
sigmae = 0.01;
rho=0.95;
teta=1;
model;
1/exp(c)=beta*((1/(exp(c(+1))))*((exp(a(+1))*alpha*exp(k(+1))^(alpha-1)*exp(n(+1))^(1-alpha))+(1-delta)));
teta/(1-exp(n))=(1/exp(c))*(exp(k)^(alpha)*exp(a)*((1-alpha)*exp(n)^(-alpha)));
exp(k)=(exp(a)*exp(k(-1))^(alpha)*exp(n(-1))^(1-alpha))-exp(c)+((1-delta)*exp(k(-1)));
exp(y)=(exp(a)*exp(k(-1))^(alpha)*exp(n)^(1-alpha));
exp(i)=exp(y)-exp(c);
1/exp(c(-1))=beta*(1/exp(c))*(1+exp(r));
exp(w)=(exp(a)*exp(k)^(alpha)*((1-alpha)*exp(n)^(-alpha)));
exp(Rent)=(exp(a)*alpha*exp(k)^(alpha-1)*exp(n)^(1-alpha));
exp(a)=rho*exp(a(-1))+e;
end;
initval;
w =1;
Rent =1;
r =1;
n =1;
i =1;
y =1;
c =1;
k =1;
a=0;
end;
shocks;
var e=sigmae^2;
end;
steady;
stoch_simul(hp_filter=1600,order=1,irf=40);

Many thanks,

[jpfeifer]

The timing for capital is wrong in some equations and correct in others.

[bpyemba]

Here is my code, I spent one two days by trying to find the solution. Someone can help?

%%% Welfare-Maximizing Monetary and Fiscal Policy Rules in Closed Economy%%%
var c gc ig y b w l r rk theta gk k mc tauk tauc taul t pi pstar n j ip def;
varexo e_gc e_ig e_theta;

parameters beta, alpha, omega, gamma, psi, delta, r_ss, tauc_ss, tauk_ss, l_ss, taul_ss, b_ss, y_ss, theta_ss, gc_ss, ig_ss, pi_ss, eta, nu, rhogc, rhoig, rhotheta, rhopi, rhoy, rhodef, phikb, phikgc, phikig, phiktheta,
phicb, phicgc, phicig, phictheta, philb, philgc, philig, philtheta;

beta = .99;
alpha = .24;
omega = .38;
gamma = .39;
delta = 0.025;
r_ss = 1.1;
tauc_ss = 0.05;
tauk_ss = 0.15;
l_ss = .33;
taul_ss = 0.26;
b_ss = 0.5;
theta_ss = 1;
y_ss = 1.12923626;
gc_ss = 0.1179;
ig_ss = 0.0996;
pi_ss = 1;
eta = 0.75;
nu = 6;
rhogc = 0.95;
rhoig = 0.95;
rhotheta = 0.95;
rhopi = 1.53;
rhoy = 0;
rhodef = 0;
phikb = 8.56;
phikgc = 0;
phikig = 0;
phiktheta = 0;
phicb = 8.56;
phicgc = 0;
phicig = 0;
phictheta = 0;
philb = 8.56;
philgc = 0;
philig = 0;
philtheta = 0;
b_ss = 0;
gk_ss = (ig_ss/delta);
mup = (nu/(nu-1));
mc_ss = (1/mup);
r_ss = (1/beta);
pstar_ss = ((1-eta*pi_ss^(nu-1))/(1-eta))^(1/(1-nu));
k_ss = (alpha*mc_ss*theta_ss*(gk_ss^omega)*(l_ss^(1-alpha)))/((r_ss - (delta*tauk_ss) - (1-delta))/(1-tauk_ss));
rk_ss = alpha*mc_ss*theta_ss*(gk_ss^omega)*(l_ss^(1-alpha))*(k_ss^(alpha-1));
ip_ss = delta*k_ss;
w_ss = (1-alpha)*mc_ss*theta_ss*(gk_ss^omega)*(k_ss^alpha)*(l_ss^(-alpha));
n_ss = (1/(1-eta*beta))*theta_ss*(gk_ss^omega)*(k_ss^alpha)*(l_ss^(1-alpha))*mc_ss;
j_ss = (1/(1-eta*beta))*theta_ss*(gk_ss^omega)*(k_ss^alpha)*(l_ss^(1-alpha))*pstar_ss;
j_ss = (nu/(nu-1))*n_ss;
y_ss = theta_ss*(gk_ss^omega)*(k_ss^alpha)*(l_ss^(1-alpha));
c_ss = y_ss - ip_ss - ig_ss - gc_ss;
psi = (1/(c_ss + gamma*gc_ss))*((1-taul_ss)/(1+tauc_ss))*w_ss;
t_ss = gc_ss + ig_ss;
t_ss = tauc_ss*c_ss + tauk_ss*rk_ss*k_ss - delta*tauk_ss*k_ss + taul_ss*w_ss*l_ss;
def_ss = b_ss;

model;

exp(c) = exp(y) - exp(ip) - exp(ig) - exp(gc);
exp(k) = exp(ip) + exp(k(-1));
exp(gk) = exp(ig) + exp(gk(-1));
beta*(((exp(c) + gamma*exp(gc))/(exp(c(+1)) + gamma*exp(gc(+1))))*((1 + exp(tauc))/(1 + exp(tauc(+1))))*(1/exp(pi(+1)))*exp(r)) = 1;
beta*(((exp(c) + gamma*exp(gc))/(exp(c(+1)) + gamma*exp(gc(+1))))*((1 + exp(tauc))/(1 + exp(tauc(+1))))*(((1-exp(tauk(+1)))*exp(rk(+1)) + (delta*exp(tauk(+1))) + (1-delta)))) = 1;
psi = (1/(exp(c) + gamma*exp(gc))) *((1 - exp(taul))/(1 + exp(tauc)))*exp(w);
exp(y) = exp(theta)*(exp(gk)^omega)*(exp(k)^alpha)*(exp(l)^(1-alpha));
exp(rk) = alpha*exp(mc)*exp(y)/exp(k);
exp(w) = (1-alpha)*exp(mc)*(exp(y)/exp(l));
exp(j) = (nu/(nu-1))*exp(n);
exp(pstar) = ((1 - eta*(exp(pstar)^(nu - 1)))/(1 - eta))^(1/(1-nu));
exp(n) = exp(y)*exp(mc) + (eta*beta*((exp(c) + gamma*exp(gc))/(exp(c(+1)) + gamma*exp(gc(+1))))*(exp(pi(+1)^nu)*exp(n(+1))));
exp(j) = exp(y)*exp(pstar) + (eta*beta*((exp(c) + gamma*exp(gc))/(exp(c(+1)) + gamma*exp(gc(+1))))*(exp(pstar)/exp(pstar(+1)))*(exp(pi(+1)^(nu-1))*exp(j(+1))));
exp(gc) + exp(ig) + exp(r(-1))*b(-1) = b + exp(t);
exp(t) = exp(tauc)*exp(c) + exp(tauk)*exp(rk)*exp(k(-1)) - delta*exp(tauk)*exp(k(-1)) + exp(w)*exp(l)*exp(taul);
exp(gc - STEADY_STATE(gc)) = exp(rhogc*(gc(-1) - STEADY_STATE(gc)) + e_gc);
exp(ig - STEADY_STATE(ig)) = exp(rhogc*(ig(-1) - STEADY_STATE(ig)) + e_ig);
exp(theta - STEADY_STATE(theta)) = exp(rhotheta*(theta(-1) - STEADY_STATE(theta)) + e_theta);
def = b - b(-1);
exp(tauc - STEADY_STATE(tauc)) = phicb*(b - STEADY_STATE(b)) + exp((phicgc*(gc - STEADY_STATE(gc))) + phicig*(ig - STEADY_STATE(ig)) + (phictheta*(theta - STEADY_STATE(theta))));
exp(tauk - STEADY_STATE(tauk)) = phikb*(b - STEADY_STATE(b)) + exp((phikgc*(gc - STEADY_STATE(gc))) + phikig*(ig - STEADY_STATE(ig)) + (phiktheta*(theta - STEADY_STATE(theta))));
exp(taul - STEADY_STATE(taul)) = philb*(b - STEADY_STATE(b)) + exp((philgc*(gc - STEADY_STATE(gc))) + philig*(ig - STEADY_STATE(ig)) + (philtheta*(theta - STEADY_STATE(theta))));
exp(r - STEADY_STATE(r)) = exp((rhopi*(pi - STEADY_STATE(pi))) + (rhoy*(y - STEADY_STATE(y)))) + (rhodef*(def - STEADY_STATE(def)));


end;

initval;

c = log(c_ss);
gc = log(gc_ss);
ig = log(ig_ss);
y = log(y_ss);
b = b_ss;
w = log(w_ss);
l = log(l_ss);
r = log(r_ss);
rk = log(rk_ss);
theta = log(theta_ss);
gk = log(gk_ss);
k = log(k_ss);
mc = log(mc_ss);
tauk = log(tauk_ss);
tauc = log(tauc_ss);
taul = log(taul_ss);
pi = log(pi_ss);
pstar = log(pstar_ss);
n = log(n_ss);
j = log(j_ss);
ip = log(ip_ss);
def = def_ss;


end;

steady;
check;
shocks;
var e_gc = 0.01;
var e_ig = 0.01;
var e_theta = 0.01;

end;

options_.printing=1;
stoch_simul(order=2, irf=25) c y l ig gc pi gk k; %, hp_filter = 1600);