求助，刚刚学习Dynare，想要解一个新凯恩斯菲利普斯曲线，在estimation那一步之前做了steady；check；stoch_simul(order = 1); 都没有报错。跑到estimation那一行就总会报错 BK不满足，没有稳定的均衡。
下面是我的code
var pi_p, pi_w, x_tilde, pi_p_bar, r, r_bar, mu_p, g, s, U, pi_w_bar;

varexo e_pi_p e_pi_w e_mu e_s e_g e_U e_x e_r e_pi_p_bar e_r_bar;


parameters theta1 theta2 theta3 theta4 phix_x phix_p phix_w phix_r psix_x1
           psix_x2 psix_x3 psix_p1 psix_p2 psix_p3 psix_w1 psix_w2 psix_w3
           psix_r1 psix_r2 psix_r3 phir_x phir_p phir_r psir_x psir_x1
           psir_x2 psir_p psir_p1 psir_p2 psir_r1 psir_r2 psir_r3;

parameters gamma_p kappa_p gamma_w alpha_w; % parameters to be estimated


% parameter values
theta1 = 0.42737;
theta2 = -0.13702;
theta3 = -0.34640;
theta4 = 0.08971;

phix_x = -0.027322;
phix_p = -0.051434;
phix_w = 0.005238;
phix_r = -0.056990;

psix_x1 = 0.367390;
psix_x2 = 0.218561;
psix_x3 = 0.009072;

psix_p1 = 0.006256;
psix_p2 = -0.101948;
psix_p3 = -0.251117;

psix_w1 = -0.101734;
psix_w2 = 0.083261;
psix_w3 = 0.042870;

psix_r1 = -0.036025;
psix_r2 = -0.195659;
psix_r3 = -0.051661;

phir_x = 0.009420;
phir_p = 0.102174;
phir_r = -0.024188;
psir_x = 0.203054;
psir_x1 = 0.449247;
psir_x2 = -0.038923;
psir_p = 0.114563;
psir_p1 = -0.029460;
psir_p2 = -0.005519;
psir_r1 = 0.452309;
psir_r2 = -0.006560;
psir_r3 = 0.017668;


gamma_p = 0.737;
kappa_p = 0.0035;
gamma_w = 1.0;
alpha_w = -0.034;


model;
%local variables
#kappa_w = 0.5*kappa_p;


% price equation
[name = 'price equation']
pi_p = (0.25*gamma_p*(0.02*(pi_p(-1) + pi_p(-2) + pi_p(-3)) + pi_p(-4)) + 0.98*EXPECTATION(-1)(pi_p(+1)) + 0.02*(1-gamma_p)*pi_p_bar(-1)+kappa_p*mu_p(-1))/(1+0.245*gamma_p)+e_pi_p;

% wage equation
[name = 'wage equation']
pi_w = (0.25*gamma_w*(0.02*(pi_w(-1) + pi_w(-2) + pi_w(-3)) + pi_w(-4)) + 0.98*EXPECTATION(-1)(pi_w(+1))+ 0.02*(1-gamma_w)*(pi_p_bar(-1)+g(-1))+alpha_w*U(-1)-kappa_w*s(-1))/(1+0.245*gamma_w)+e_pi_w;

% auxiliary equations
[name = 'mu']
mu_p=mu_p(-1)+ 0.0025 *(pi_w - pi_p - g)+e_mu;

[name = 's']
s = 1.25* mu_p + e_s;

[name = 'g']
g = g(-1) + e_g;

[name = 'U']
U = theta1*U(-1)+theta2*(U(-1)-U(-2))+theta3*x_tilde+theta4*(x_tilde-x_tilde(-1))+e_U;


[name = 'gdp gap']
x_tilde-x_tilde(-1) = phix_x*x_tilde(-1)+phix_p*(pi_p(-1)-pi_p_bar(-1))+phix_w * (pi_w(-1) - pi_w_bar(-1)) + phix_r*(r(-1) - r_bar(-1))+ (psix_x1*(x_tilde(-1)-x_tilde(-2)) + psix_x2*(x_tilde(-2)-x_tilde(-3))+psix_x3*(x_tilde(-3)-x_tilde(-4))+psix_p1*(pi_p(-1)-pi_p(-2))+psix_p2*(pi_p(-2)-pi_p(-3))+psix_p3*(pi_p(-3)-pi_p(-4))+ psix_w1*(pi_w(-1)-pi_w(-2))+psix_w2*(pi_w(-2)-pi_w(-3))+psix_w3*(pi_w(-3)-pi_w(-4))+psix_r1*(r(-1)-r(-2))+psix_r2*(r(-2)-r(-3))+psix_r3*(r(-3)-r(-4))) + e_x;


[name = 'r']
r - r(-1) = phir_x * x_tilde(-1) + phir_p*(pi_p(-1)-pi_p_bar(-1)) + phir_r*(r(-1)-r_bar(-1)) + psir_x * (x_tilde - x_tilde(-1)) + psir_x1*(x_tilde(-1)-x_tilde(-2)) + psir_x2*(x_tilde(-2)-x_tilde(-3)) + psir_p*(pi_p-pi_p(-1)) + psir_p1*(pi_p(-1)-pi_p(-2)) + psir_p2*(pi_p(-2)-pi_p(-3)) + psir_r1*(r(-1) - r(-2)) + psir_r2*(r(-2)-r(-3)) + psir_r3*(r(-3)-r(-4)) + e_r;


[name = 'pi_p_bar']
pi_p_bar = pi_p_bar(-1) + e_pi_p_bar;


[name = 'r_bar']
r_bar = r_bar(-1) + e_r_bar;


[name = 'pi_w']
pi_w_bar = pi_p_bar+g;
end;


estimated_params;
    gamma_p, 0.2, 0, 1;
    kappa_p, 0.3, 0, Inf;
    gamma_w, 0.4, 0, 1;
    alpha_w, -0.03, -Inf, 0;
    stderr e_pi_p, 0.1;
    stderr e_pi_w, 0.1;
    stderr e_mu, 0.1;
    stderr e_s, 0.1;
    stderr e_g, 0.1;
    stderr e_U, 0.044;
    stderr e_x, 0.5055;
    stderr e_r, 0.8404;
    stderr e_pi_p_bar, 0.1;
    stderr e_r_bar, 0.1;
    corr e_pi_p, e_pi_w, 0.8;
end;

steady;
check;
stoch_simul(order = 1, nograph);

varobs pi_p pi_w x_tilde pi_p_bar r r_bar mu_p g s U;

estimated_params_init(use_calibration);
end;
estimation(datafile=load_nkpc_data)pi_p pi_w x_tilde pi_p_bar r r_bar mu_p g s U;