@RJDennis the deterministic simulations look great now!

I’m having some issues w/ simulations for stochastic models.

```
#Stochastic NCGM
states:
k, z
end
jumps:
c
end
shocks:
ε
end
parameters:
β = 0.99
σ = 1.1
δ = 0.025
α = 0.30
ρ = 0.8
σ_ε = 0.01
end
equations:
k(+1) = (1.0 - δ)*k + exp(z)*k^α - c
c^(-σ) = β*c(+1)^(-σ)*(1.0 - δ + α*exp(z(+1))*k(+1)^(α - 1.0))
z(+1) = ρ*z + σ_ε*ε
end
```

Julia program:

```
filename = "StochasticNCGM.txt";
path = joinpath(@__DIR__,filename); m = get_model(path);
xg = [0.05, 3.05, 0.7]; tol = 1e-8; maxiters = 1000;
ss = compute_steady_state(m,xg,tol,maxiters); #Steady state #zss, kss, css
```

Third order perturbation gives an error message when simulating

```
P3 = PerturbationScheme(ss,1.0,"third");
s_P3 = solve_model(m,P3)
sim_data = simulate(s_P3, x0b, T)
UndefVarError: hx not defined
```

Next

```
P2 = PerturbationScheme(ss,1.0,"second");
s_P2 = solve_model(m,P2);
C1 = ChebyshevSchemeStoch(ss,chebyshev_nodes,[21,21],9,3,[0.07 23.0; -0.07 19.5],tol,1e-6,maxiters);
s_C1 = solve_model(m,C1);
L1 = SmolyakSchemeStoch(ss,chebyshev_gauss_lobatto,9,4,[0.07 23.0; -0.07 19.5],tol,1e-6,maxiters);
s_L1 = solve_model(m,L1);
T=100; x0a= ss[1:2]*2;
#
sim_data = simulate(s_P2, x0a, T)
plot(1:T, sim_data', label = hcat(m.variables...))
plot!(ss, seriestype = :hline, label=false)
#
sim_data = simulate(s_C1, x0a, T)
plot(1:T, sim_data', label = hcat(m.variables...))
plot!(ss, seriestype = :hline, label=false)
#
sim_data = simulate(s_L1, x0a, T)
plot(1:T, sim_data', label = hcat(m.variables...))
plot!(ss, seriestype = :hline, label=false)
#
T_irf=50; ix_shock=1; n_MC= 10000;
#
pos_imps, neg_imps = impulses(s_P2,T_irf,ix_shock,n_MC)
plot(1:50, pos_imps', title = "Positive IRF", label = hcat(m.variables...))
plot(1:50, neg_imps', title = "Negative IRF", label = hcat(m.variables...))
#
pos_imps, neg_imps = impulses(s_C1,T_irf,ix_shock,n_MC)
plot(1:50, pos_imps', title = "Positive IRF", label = hcat(m.variables...))
plot(1:50, neg_imps', title = "Negative IRF", label = hcat(m.variables...))
```

The IRFs look great but the simulations don’t seem to start above steady state.