Replicating Forbes, Hjortsoe and Nenova (2015): “The shocks matter: improving our estimates of exchange rate pass-through.”

This model replicates the standard NK DSGE model for a small open economy presented in Forbes, Hjortsoe and Nenova (2015): “The shocks matter: improving our estimates of exchange rate pass-through.”

Someone could help me?

Dynare answers the BK’s rank condition is not satisfied:
"There are 23 eigenvalue(s) larger than 1 in modulus
for 25 forward-looking variable(s)

The rank condition ISN’T verified!"

Can someone help me?

parameters thetah thetaf sigma n beta alphah alphaf tauh tauf ah af phih phif eta gammahpcp gammafpcp ibar pibar ahR ahPI afR afPI delta b rhoC rhoCstar rhoS;

var Yh Yf C Cstar L Lstar Ch Chstar Cf Cfstar Q bf bfstar i istar
ph pf phpcp phlcp pflcp pfstar phstar pfpcpstar phlcpstar pflcpstar
pih pihpcp pihlcp pihlcpstar pifstar pifpcpstar piflcpstar piflcp w wstar pi pistar mkp
x1pcp x2pcp x1lcp x1lcpstar x2lcp x2lcpstar x1fpcp x2fpcp x1flcp x1flcpstar x2flcp x2flcpstar
gammaC gammaCstar gammaS EPSILONS;

varexo epsilonC epsilonCstar epsilonS gammaI gammaIstar;

n = 0.05;
beta = 0.99;
delta = 0.01;
b = 117.486685;
op = 0.3;
ah = 0.715; //1-(1-n)*op
af = 0.715; //1-(1-n)*op
eta = 2;
sigma = 1.1;
alphah = 0.6;
alphaf = 0.6;
thetah = 7;
thetaf = 7;
tauh = 1/(1-thetah);
tauf = 1/(1-thetaf);
phih = 0.75;
phif = 0.75;
gammahpcp = 0.4;
gammafpcp = 0.4;
ahR = 0.7;
ahPI = 1.5;
afR = 0.7;
afPI = 1.5;
ibar = 0.010101;
pibar = 1;
rhoC = 0.9;
rhoCstar = 0.9;
rhoS = 0.9;

initval;
Yh = 0.991727466;
Yf = 0.991727466;
C = 2.166593153;
Cstar = 0.20442523;
L = 0.991727466;
Lstar = 0.991727466;
Ch = 0.828855014;
Chstar = 0.008572234;
Cf = 8.424094644;
Cfstar = 0.548354064;
Q = 0.178824091;
bf = 656.9958472;
bfstar = -34.5787288;
i = 0.010101;
istar = 0.010101;
ph = 2.302182697;
pf = 0.030675995;
phpcp = 1;
phlcp = 1;
pflcp = 1;
pfstar = 0.171542852;
phstar = 12.87400757;
pfpcpstar = 1;
phlcpstar = 1;
pflcpstar = 1;
pi = 1;
pistar = 1;
pih = 1;
pihpcp = 1;
pihlcp = 1;
pihlcpstar = 1;
pifstar = 1;
pifpcpstar = 1;
piflcpstar = 1;
piflcp = 1;
w = 2.302182697;
wstar = 0.171542852;
x1pcp = 2.8028366;
x2pcp = 2.8028366;
x1lcp = 2.342523776;
x2lcp = 2.342523776;
x1lcpstar = 0.024226991;
x2lcpstar = 0.024226991;
x1fpcp = 2.8028366;
x2fpcp = 2.8028366;
x1flcp = 23.80831589;
x2flcp = 23.80831589;
x1flcpstar = 1.549767343;
x2flcpstar = 1.549767343;
mkp = 1;
gammaC = 1;
gammaCstar = 1;
gammaS = 1;
end;

model;

// domestic LCP and PCP firms

x1pcp = (thetah/(thetah-1)) (C^(-sigma)) w (phpcp^(-thetah)) (Ch+((1-n)/n) (ph^(-thetah)) (Q^thetah) (phstar^thetah) Chstar)+beta alphah (pihpcp(+1)^thetah) x1pcp(+1); //1
x2pcp = (1-tauh)
(C^(-sigma)) ph (phpcp^(1-thetah)) (Ch+((1-n)/n) (ph^(-thetah)) (Q^thetah) (phstar^thetah) Chstar)+beta alphah*(pihpcp(+1)^(thetah-1))*x2pcp(+1); //2

((1-alphah*pihpcp^(thetah-1))/(1-alphah))^(1/(1-thetah))=x1pcp/x2pcp; //3

x1lcp = (thetah/(thetah-1)) (C^(-sigma)) w (phlcp^(-thetah)) Ch+beta alphah (pihlcp(+1)^thetah) x1lcp(+1); //4
x2lcp = (1-tauh)
(C^(-sigma)) ph (phlcp^(1-thetah)) Ch+beta alphah*(pihlcp(+1)^(thetah-1))*x2lcp(+1); //5

((1-alphah*pihlcp^(thetah-1))/(1-alphah))^(1/(1-thetah))=x1lcp/x2lcp; //6

x1lcpstar = (thetah/(thetah-1)) (C^(-sigma)) w (phlcpstar^(-thetah)) Chstar+beta alphah (pihlcpstar(+1)^thetah) x1lcpstar(+1); //7
x2lcpstar = (1-tauh)
(C^(-sigma)) phstar Q*(phlcpstar^(1-thetah)) Chstar+beta alphah*(pihlcpstar(+1)^(thetah-1))*x2lcpstar(+1); //8

((1-alphah*pihlcpstar^(thetah-1))/(1-alphah))^(1/(1-thetah))=x1lcpstar/x2lcpstar; //9

// foreign LCP and PCP firms

x1fpcp = (thetaf/(thetaf-1)) (Cstar^(-sigma)) wstar (pfpcpstar^(-thetaf)) (Cfstar+(n/(1-n)) ((Q pfstar)^(-thetaf)) (pf^thetaf) Cf)+beta alphaf (pifpcpstar(+1)^thetaf) x1fpcp(+1); //10
x2fpcp = (1-tauf) pfstar (Cstar^(-sigma))
(pfpcpstar^(1-thetaf)) (Cfstar+(n/(1-n)) ((Q pfstar)^(-thetaf)) (pf^thetaf) Cf)+beta alphaf*(pifpcpstar(+1)^(thetaf-1))*x2fpcp(+1); //11

((1-alphaf*pifpcpstar^(thetaf-1))/(1-alphaf))^(1/(1-thetaf))=x1fpcp/x2fpcp; //12

x1flcpstar = (thetaf/(thetaf-1)) (Cstar^(-sigma)) wstar (pflcpstar^(-thetaf)) Cfstar+beta alphaf (piflcpstar(+1)^thetaf) x1flcpstar(+1); //13
x2flcpstar = (1-tauf)
(Cstar^(-sigma)) pfstar (pflcpstar^(1-thetah)) Cfstar+beta alphaf*(piflcpstar(+1)^(thetaf-1))*x2flcpstar(+1); //14

((1-alphaf*piflcpstar^(thetaf-1))/(1-alphaf))^(1/(1-thetaf))=x1flcpstar/x2flcpstar; //15

x1flcp = (thetaf/(thetaf-1)) (Cstar^(-sigma)) wstar (pflcp^(-thetaf)) Cf+beta alphaf (piflcp(+1)^thetaf) x1flcp(+1); //16
x2flcp = (1-tauf)
(Cstar^(-sigma)) pf (pflcp^(1-thetaf)) Cf Q^(-1)+beta alphaf (piflcp(+1)^(thetaf-1))*x2flcp(+1); //17

((1-alphaf*piflcp^(thetaf-1))/(1-alphaf))^(1/(1-thetaf))=x1flcp/x2flcp; //18

// consumption demand

Ch = ah C ph^(-phih); //19
Cf = (1-ah) C pf^(-phih); //20
Chstar = (1-af) Cstar phstar^(-phif); //21
Cfstar = af Cstar pfstar^(-phif); //22

// labor supply

(L^eta)*(C^sigma) gammaC^(-1) = w; //23
(Lstar^eta)
(Cstar^sigma)*gammaCstar^(-1) = wstar; //24

// composition Euler equations

beta*((gammaC(+1) C(+1)^(-sigma))/(gammaC C^(-sigma))) ((1+i)/pi(+1))=1; //25
beta
((gammaC(+1) C(+1)^(-sigma))/(gammaC C^(-sigma))) ((1+istar)/pistar(+1)) (Q(+1)/Q) = gammaS/(1-delta*(Q bf - b)); //26
beta
((gammaCstar(+1) Cstar(+1)^(-sigma))/(gammaCstar Cstar^(-sigma)))*((1+istar)/pistar(+1))=1; //27

pi = ph(-1)*pih/ph; //28
pistar = pfstar(-1) pifstar/pfstar; //29
mkp = pf/(Q
wstar); //30

// price indices

1 = ah ph^(1-phih)+(1-ah) pf^(1-phih); //31
1 = af
pfstar^(1-phif)+(1-af) phstar^(1-phif); //32
1 = gammahpcp
phpcp^(1-thetah)+(1-gammahpcp) phlcp^(1-thetah); //33
1 = gammahpcp
(phpcp
ph*(Q^(-1)) phstar^(-1))^(1-thetah)+(1-gammahpcp) phlcpstar^(1-thetah); //34
1 = gammafpcp
pfpcpstar^(1-thetaf)+(1-gammafpcp) pflcpstar^(1-thetaf); //35
1 = gammafpcp
(pfpcpstar
(pfstar*Q/pf))^(1-thetaf)+(1-gammafpcp)*pflcp^(1-thetaf); //36

// inflation rates

pihpcp = phpcp pih/phpcp(-1); //37
pifpcpstar = pfpcpstar
pifstar/pfpcpstar(-1); //38
pihlcp = phlcp pih/phlcp(-1); //39
pihlcpstar = (phlcpstar
phstar*pfstar(-1) pifstar)/(phlcpstar(-1) phstar(-1) pfstar); //40
piflcpstar = pflcpstar
pifstar/pflcpstar(-1); //41
piflcp = (pflcp
pf
ph(-1)*pih)/(pflcp(-1)*pf(-1)*ph); //42

// production functions

Yh = L; //43
Yf = Lstar; //44

// resource constraints

C+Q bf/((1+istar) (1-delta*(Q bf - b))) = Yh+(Q bf(-1)/pistar); // //45
n*bf+(1-n)*bfstar = 0; //46

// monetary policy rules

log(i/ibar) = ahR log(i(-1)/ibar)+ahPI log(pi/pibar)+gammaI;
log(istar/ibar) = afR log(istar(-1)/ibar)+afPI log(pistar/pibar)+gammaIstar;

EPSILONS = -epsilonS;

// good markets equilibrium

Yh = ah*(ph^(-phih)) C+((1-n)/n) (1-af) (ph^(-phih)) (Q^phih) Cstar; //49
Yf = (n/(1-n))
(1-ah) (pfstar^(-phif)) (Q^(-phif)) C+af (pfstar^(-phif))*Cstar; //50

// other shocks

gammaC = 1-rhoC+rhoC gammaC(-1)+epsilonC; //51
gammaCstar = 1-rhoCstar+rhoCstar
gammaCstar(-1)+epsilonCstar; //52
gammaS = 1-rhoS+rhoS*gammaS(-1)+EPSILONS; //53

end;

model_diagnostics;
resid;

steady;
check;
shocks;
var epsilonC = 1;
var epsilonCstar = 1;
var gammaI = 1;
var gammaIstar = 1;
var epsilonS = 1;
end;

stoch_simul(order=1,irf = 100,nograph,periods=50000);

Please do not double-post. See BK rank condition not satisfied