diff --git a/configure.ac b/configure.ac
index 712636d2dbe2479acd13831fb13e2679a5b2e0fb..1dbc85b2b69bb44947163235164eafb9a559640d 100644
--- a/configure.ac
+++ b/configure.ac
@@ -31,6 +31,7 @@ case ${host_os} in
     # On mingw32, we don't want dynamic libgcc
     # Note that static-libstdc++ is only supported since GCC 4.5 (but generates no error on older versions)
     LDFLAGS="$LDFLAGS -static-libgcc -static-libstdc++"
+    have_windows="yes"
     ;;
   *cygwin*)
     AC_MSG_WARN([You are compiling for the Cygwin target. This means that the preprocessor will])
@@ -41,6 +42,7 @@ case ${host_os} in
       # And by default, the AC_PROG_F77 will pick up g77 if it is present (even if gfortran is also here)
       F77=gfortran
     fi
+    have_windows="yes"
     ;;
 esac
 
@@ -196,8 +198,10 @@ AC_CONFIG_FILES([Makefile
 AC_ARG_ENABLE([matlab], AS_HELP_STRING([--disable-matlab], [disable compilation of MEX files for MATLAB]), [], [enable_matlab=yes])
 if test "x$enable_matlab" = "xyes"; then
   AC_CONFIG_SUBDIRS([mex/build/matlab])
+  AX_MATLAB
 fi
 AM_CONDITIONAL([ENABLE_MATLAB], [test "x$enable_matlab" = "xyes"])
+AM_CONDITIONAL([HAVE_CMD_LINE_MATLAB], [test "x$ax_enable_matlab" = "xyes" -a "x$have_windows" = "x"])
 
 AC_ARG_ENABLE([octave], AS_HELP_STRING([--disable-octave], [disable compilation of MEX files for Octave]), [], [enable_octave=yes])
 if test "x$enable_octave" = "xyes"; then
diff --git a/tests/.gitignore b/tests/.gitignore
index cfe6100aa26c0b83e229a71f0371234d12c8d0dd..f331946c03f0e3816d25c76bb9c23be3d4f73dfa 100644
--- a/tests/.gitignore
+++ b/tests/.gitignore
@@ -10,18 +10,27 @@
 *.dat
 *_simul
 
-/block_bytecode/ws
+wsOct
+/run_test_octave_output.txt
+/run_test_matlab_output.txt
+
 /block_bytecode/ls2003_tmp.mod
 /partial_information/PItest3aHc0PCLsimModPiYrVarobsAll_PCL*
 /partial_information/PItest3aHc0PCLsimModPiYrVarobsCNR_PCL*
 
+!/internals/tests.m
+!/fs2000_ssfile_aux.m
+!/printMakeCheckMatlabErrMsg.m
+!/printMakeCheckOctaveErrMsg.m
+!/ramst_initval_file_data.m
+!/run_test.m
+!/run_test_matlab.m
+!/run_test_octave.m
+!/test.m
 !/AIM/data_ca1.m
 !/AIM/fs2000_b1L1L_AIM_steadystate.m
 !/AIM/fs2000_b1L1L_steadystate.m
 !/AIM/fsdat.m
-!/block_bytecode/MARK3_endo.dat
-!/block_bytecode/MARK3_exo.dat
-!/block_bytecode/run_block_bytecode_tests.m
 !/block_bytecode/run_ls2003.m
 !/bvar_a_la_sims/bvar_sample.m
 !/external_function/extFunDeriv.m
@@ -29,6 +38,16 @@
 !/external_function/extFunWithFirstAndSecondDerivs.m
 !/fs2000/fs2000a_steadystate.m
 !/fs2000/fsdat_simul.m
+!/gsa/data_ca1.m
+!/gsa/data_ca1.mat
+!/gsa/ls2003_mean.mat
+!/gsa/ls2003_mode.mat
+!/gsa/ls2003_results.mat
+!/gsa/ls2003scr_mean.mat
+!/gsa/ls2003scr_mode.mat
+!/gsa/ls2003scr_results.mat
+!/identification/kim/kim2_steadystate.m
+!/identification/as2007/as2007_steadystate.m
 !/kalman/likelihood/compare_kalman_routines.m
 !/kalman/likelihood/simul_state_space_model.m
 !/kalman/likelihood/test1.m
@@ -37,19 +56,15 @@
 !/kalman_filter_smoother/test.mat
 !/kalman_filter_smoother/testsmoother.m
 !/k_order_perturbation/run_fs2000kplusplus.m
-!/ls2003/data_ca1.m
 !/objectives/sgu_ex1.mat
+!/ls2003/data_ca1.m
+!/measurement_errors/data_ca1.m
 !/missing/simulate_data_with_missing_observations.m
 !/parallel/data_ca1.m
 !/pi2004/idata.m
 !/pi2004/ych.dat
 !/practicing/cagan_data.mat
-!/practicing/dataHST.mat
 !/practicing/data_consRicardoypg.mat
+!/practicing/dataHST.mat
 !/practicing/datasaver.m
-!/ramst_initval_file_data.m
 !/recursive/data_ca1.m
-!/run_test.m
-!/run_test_octave.m
-!/test.m
-!/fs2000_ssfile_aux.m
diff --git a/tests/Makefile.am b/tests/Makefile.am
index c7e09ed38c2326e8173635be37d8b5a2bbb56ad3..d1e4f236fef917317698f5d84be6f198d2cf9638 100644
--- a/tests/Makefile.am
+++ b/tests/Makefile.am
@@ -1,12 +1,10 @@
-DYNARE_ROOT = $(abs_top_srcdir)/matlab
-
-MODS = \
+MODFILES = \
 	ramst.mod \
 	ramst_a.mod \
 	example1.mod \
 	example2.mod \
 	example1_use_dll.mod \
-	example1_with_tags.mod\
+	example1_with_tags.mod \
 	t_sgu_ex1.mod \
 	osr_example.mod \
 	optimal_policy/ramsey.mod \
@@ -20,7 +18,15 @@ MODS = \
 	comments.mod \
 	histval_sto.mod \
 	histval_det.mod \
-	expectation.mod \
+	auxiliary_variables/test1.mod \
+	expectations/expectation.mod \
+	expectations/expectation_ss.mod \
+	expectations/expectation_ss_old.mod \
+	steady_state/walsh1_initval.mod \
+	steady_state/walsh1_old_ss.mod \
+	steady_state/walsh1_ssm.mod \
+	steady_state/walsh1_ssm_block.mod \
+	steady_state/multi_leads.mod \
 	steady_state_operator/standard.mod \
 	steady_state_operator/use_dll.mod \
 	steady_state_operator/block.mod \
@@ -47,6 +53,7 @@ MODS = \
 	fs2000/fs2000.mod \
 	fs2000/fs2000a.mod \
 	fs2000/fs2000c.mod \
+	fs2000/fs2000d.mod \
 	homotopy/homotopy1_test.mod \
 	homotopy/homotopy2_test.mod \
 	homotopy/homotopy3_test.mod \
@@ -73,51 +80,60 @@ MODS = \
 	external_function/example1_no_deriv_functions_provided.mod \
 	external_function/example1_no_deriv_functions_provided_dll.mod \
 	seeds.mod \
-	simul/example1.mod
-
-EXTRA_DIST = $(MODS) \
+	gsa/ls2003.mod \
+	identification/kim/kim2.mod \
+	identification/as2007/as2007.mod \
+	simul/example1.mod \
+	conditional_forecasts/fs2000_est.mod
+
+EXTRA_DIST = \
+	$(MODFILES) \
+	homotopy/common.mod \
+	block_bytecode/ls2003.mod \
+	fs2000_ssfile_aux.m \
+	printMakeCheckMatlabErrMsg.m \
+	printMakeCheckOctaveErrMsg.m \
+	ramst_initval_file_data.m \
+	run_test.m \
+	run_test_matlab.m \
 	run_test_octave.m \
-	block_bytecode/run_block_bytecode_tests.m \
+	test.m \
+	AIM/data_ca1.m \
+	AIM/fs2000_b1L1L_AIM_steadystate.m \
+	AIM/fs2000_b1L1L_steadystate.m \
+	AIM/fsdat.m \
 	block_bytecode/run_ls2003.m \
-	block_bytecode/ls2003.mod \
+	bvar_a_la_sims/bvar_sample.m \
 	external_function/extFunDeriv.m \
 	external_function/extFunNoDerivs.m \
 	external_function/extFunWithFirstAndSecondDerivs.m \
-	fs2000/fsdat_simul.m \
+	expectations/expectation_ss_old_steadystate.m \
 	fs2000/fs2000a_steadystate.m \
-	identification/kim/kim2_steadystate.m \
-	identification/as2007/as2007_steadystate.m \
-	AIM/fs2000_b1L1L_steadystate.m \
-	AIM/fs2000_b1L1L_AIM_steadystate.m \
-	test.m \
-	objectives \
-	ramst_initval_file_data.m \
-	homotopy/common.mod \
-	bvar_a_la_sims/bvar_sample.m \
-	fs2000_ssfile_aux.m
+	fs2000/fsdat_simul.m \
+	k_order_perturbation/run_fs2000kplusplus.m \
+	ls2003/data_ca1.m \
+	measurement_errors/data_ca1.m \
+	missing/simulate_data_with_missing_observations.m \
+	objectives/sgu_ex1.mat \
+	conditional_forecasts/fsdat_simul.m
 
-TARGETS = check-matlab
+TARGETS =
+
+if HAVE_CMD_LINE_MATLAB
+TARGETS += check-matlab
+endif
 
 if HAVE_OCTAVE
-TARGETS += check-octave check-block-bytecode
+TARGETS += check-octave
 endif
 
 check-local: $(TARGETS)
 
-check-octave: $(MODS)
-	@set -e; \
-		for modfile in $(MODS); do \
-			$(OCTAVE) --norc --silent --no-history run_test_octave.m $$modfile $(DYNARE_ROOT) $(PACKAGE_VERSION); \
-		done
-.PHONY: check-octave
-
-check-block-bytecode:
-	cd block_bytecode && $(OCTAVE) --norc --silent --no-history run_block_bytecode_tests.m $(DYNARE_ROOT) $(PACKAGE_VERSION)
-.PHONY: check-block-bytecode
+check-matlab:
+	DYNARE_VERSION="$(PACKAGE_VERSION)" MODFILES="$(MODFILES)" $(MATLAB)/bin/matlab -nodesktop -nosplash -r run_test_matlab
 
-check-matlab: $(MODS)
-# MATLAB stuff to be added here
-.PHONY: check-matlab
+check-octave:
+	DYNARE_VERSION="$(PACKAGE_VERSION)" MODFILES="$(MODFILES)" $(OCTAVE) --norc --silent --no-history run_test_octave.m
 
 clean-local:
 	rm -f $(patsubst %.mod, %.m, $(MODS)) \
@@ -159,4 +175,9 @@ clean-local:
 
 	rm -rf partial_information/PItest3aHc0PCLsimModPiYrVarobsAll_PCL* partial_information/PItest3aHc0PCLsimModPiYrVarobsCNR_PCL*
 
-	rm -rf block_bytecode/ws block_bytecode/ls2003_tmp*
+	rm -rf block_bytecode/ls2003_tmp*
+
+	rm -f $(shell find -name wsOct) \
+		$(shell find -name wsMat.mat)
+
+	rm -f run_test_matlab_output.txt run_test_octave_output.txt
diff --git a/tests/auxiliary_variables/test1.mod b/tests/auxiliary_variables/test1.mod
new file mode 100644
index 0000000000000000000000000000000000000000..c4c01c61dfcfa72487038e583a575b2f175b619a
--- /dev/null
+++ b/tests/auxiliary_variables/test1.mod
@@ -0,0 +1,14 @@
+var x, y;
+varexo e;
+
+model;
+expectation(-2)(x(+4)*y(+1)) = 1;
+y = exp(e);
+end;
+
+initval;
+y = 1;
+x = 1;
+end;
+
+steady;
diff --git a/tests/block_bytecode/run_block_bytecode_tests.m b/tests/block_bytecode/run_block_bytecode_tests.m
deleted file mode 100644
index 86d3175894ef459c3ad697f19158d55994f882c8..0000000000000000000000000000000000000000
--- a/tests/block_bytecode/run_block_bytecode_tests.m
+++ /dev/null
@@ -1,61 +0,0 @@
-## Copyright (C) 2010 Dynare Team
-##
-## This file is part of Dynare.
-##
-## Dynare is free software: you can redistribute it and/or modify
-## it under the terms of the GNU General Public License as published by
-## the Free Software Foundation, either version 3 of the License, or
-## (at your option) any later version.
-##
-## Dynare is distributed in the hope that it will be useful,
-## but WITHOUT ANY WARRANTY; without even the implied warranty of
-## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
-## GNU General Public License for more details.
-##
-## You should have received a copy of the GNU General Public License
-## along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
-
-addpath(argv(){1})
-
-if !strcmp(dynare_version(), argv(){2})
-  error("Incorrect version of Dynare is being tested")
-endif
-
-## Ask gnuplot to create graphics in text mode
-## Note that setenv() was introduced in Octave 3.0.2, for compatibility
-## with MATLAB
-putenv("GNUTERM", "dumb")
-
-for block = 0:1
-  for bytecode = 0:1
-    ## Recall that solve_algo=7 and stack_solve_algo=2 are not supported
-    ## under Octave
-    default_solve_algo = 2;
-    default_stack_solve_algo = 0;
-    if !block && !bytecode
-      solve_algos = 0:4;
-      stack_solve_algos = 0;
-    elseif block && !bytecode
-      solve_algos = [0:4 6 8];
-      stack_solve_algos = [0 1 3 4];
-    else
-      solve_algos = [0:6 8];
-      stack_solve_algos = [0 1 3:5];
-    endif
-
-    for i = 1:length(solve_algos)
-      save ws
-      run_ls2003(block, bytecode, solve_algos(i), default_stack_solve_algo)
-      load ws
-    endfor
-    for i = 1:length(stack_solve_algos)
-      save ws
-      run_ls2003(block, bytecode, default_solve_algo, stack_solve_algos(i))
-      load ws
-    endfor
-  endfor
-endfor
-
-## Local variables:
-## mode: Octave
-## End:
diff --git a/tests/block_bytecode/run_ls2003.m b/tests/block_bytecode/run_ls2003.m
index 61d7fd838aa5a63b1f1e7d7a07dea31b83860b9d..a5f898519059594a808c5a4077c5dcbfbbeffc11 100644
--- a/tests/block_bytecode/run_ls2003.m
+++ b/tests/block_bytecode/run_ls2003.m
@@ -1,29 +1,31 @@
-## Copyright (C) 2010 Dynare Team
-##
-## This file is part of Dynare.
-##
-## Dynare is free software: you can redistribute it and/or modify
-## it under the terms of the GNU General Public License as published by
-## the Free Software Foundation, either version 3 of the License, or
-## (at your option) any later version.
-##
-## Dynare is distributed in the hope that it will be useful,
-## but WITHOUT ANY WARRANTY; without even the implied warranty of
-## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
-## GNU General Public License for more details.
-##
-## You should have received a copy of the GNU General Public License
-## along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
-
 function run_ls2003(block, bytecode, solve_algo, stack_solve_algo)
-  printf("\nTEST: ls2003 (block=%d, bytecode=%d, solve_algo=%d, stack_solve_algo=%d)...\n", block, bytecode, solve_algo, stack_solve_algo);
-  fid = fopen("ls2003_tmp.mod", "w");
+
+% Copyright (C) 2010-2011 Dynare Team
+%
+% This file is part of Dynare.
+%
+% Dynare is free software: you can redistribute it and/or modify
+% it under the terms of the GNU General Public License as published by
+% the Free Software Foundation, either version 3 of the License, or
+% (at your option) any later version.
+%
+% Dynare is distributed in the hope that it will be useful,
+% but WITHOUT ANY WARRANTY; without even the implied warranty of
+% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+% GNU General Public License for more details.
+%
+% You should have received a copy of the GNU General Public License
+% along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
+
+  disp(['TEST: ls2003 (block=' num2str(block) ', bytecode=' ...
+      num2str(bytecode) ', solve_algo=' num2str(solve_algo) ...
+      ', stack_solve_algo=' num2str(stack_solve_algo) ')...']);
+  fid = fopen('ls2003_tmp.mod', 'w');
   assert(fid > 0);
-  fprintf(fid, "@#define block = %d\n@#define bytecode = %d\n@#define solve_algo = %d\n@#define stack_solve_algo = %d\n@#include \"ls2003.mod\"\n", block, bytecode, solve_algo, stack_solve_algo)
+  fprintf(fid, ['@#define block = %d\n@#define bytecode = %d\n' ...
+      '@#define solve_algo = %d\n@#define stack_solve_algo = %d\n' ...
+      '@#include \"ls2003.mod\"\n'], block, bytecode, ...
+      solve_algo, stack_solve_algo);
   fclose(fid);
-  dynare("ls2003_tmp.mod")
-endfunction
-
-## Local variables:
-## mode: Octave
-## End:
+  dynare('ls2003_tmp.mod','console')
+end
diff --git a/tests/conditional_forecasts/fs2000_est.mod b/tests/conditional_forecasts/fs2000_est.mod
new file mode 100644
index 0000000000000000000000000000000000000000..6f07765f76c83382c9d545128ce3a9c763913a41
--- /dev/null
+++ b/tests/conditional_forecasts/fs2000_est.mod
@@ -0,0 +1,135 @@
+// See fs2000.mod in the examples/ directory for details on the model
+
+var m P c e W R k d n l gy_obs gp_obs y dA;
+varexo e_a e_m;
+
+parameters alp bet gam mst rho psi del;
+
+alp = 0.33;
+bet = 0.99;
+gam = 0.003;
+mst = 1.011;
+rho = 0.7;
+psi = 0.787;
+del = 0.02;
+
+model;
+dA = exp(gam+e_a);
+log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m;
+-P/(c(+1)*P(+1)*m)+bet*P(+1)*(alp*exp(-alp*(gam+log(e(+1))))*k^(alp-1)*n(+1)^(1-alp)+(1-del)*exp(-(gam+log(e(+1)))))/(c(+2)*P(+2)*m(+1))=0;
+W = l/n;
+-(psi/(1-psi))*(c*P/(1-n))+l/n = 0;
+R = P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(-alp)/W;
+1/(c*P)-bet*P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)/(m*l*c(+1)*P(+1)) = 0;
+c+k = exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)+(1-del)*exp(-(gam+e_a))*k(-1);
+P*c = m;
+m-1+d = l;
+e = exp(e_a);
+y = k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a));
+gy_obs = dA*y/y(-1);
+gp_obs = (P/P(-1))*m(-1)/dA;
+end;
+
+initval;
+k = 6;
+m = mst;
+P = 2.25;
+c = 0.45;
+e = 1;
+W = 4;
+R = 1.02;
+d = 0.85;
+n = 0.19;
+l = 0.86;
+y = 0.6;
+gy_obs = exp(gam);
+gp_obs = exp(-gam);
+dA = exp(gam);
+end;
+
+shocks;
+var e_a; stderr 0.014;
+var e_m; stderr 0.005;
+end;
+
+steady_state_model;
+  dA = exp(gam);
+  gst = 1/dA;
+  m = mst;
+  
+  khst = ( (1-gst*bet*(1-del)) / (alp*gst^alp*bet) )^(1/(alp-1));
+  xist = ( ((khst*gst)^alp - (1-gst*(1-del))*khst)/mst )^(-1);
+  nust = psi*mst^2/( (1-alp)*(1-psi)*bet*gst^alp*khst^alp );
+  n  = xist/(nust+xist);
+  P  = xist + nust;
+  k  = khst*n;
+
+  l  = psi*mst*n/( (1-psi)*(1-n) );
+  c  = mst/P;
+  d  = l - mst + 1;
+  y  = k^alp*n^(1-alp)*gst^alp;
+  R  = mst/bet;
+
+  W = l/n;
+  e = 1;
+
+  gp_obs = m/dA;
+  gy_obs = dA;
+end;
+
+steady;
+
+check;
+
+estimated_params;
+alp, beta_pdf, 0.356, 0.02;
+bet, beta_pdf, 0.993, 0.002;
+gam, normal_pdf, 0.0085, 0.003;
+mst, normal_pdf, 1.0002, 0.007;
+rho, beta_pdf, 0.129, 0.223;
+psi, beta_pdf, 0.65, 0.05;
+del, beta_pdf, 0.01, 0.005;
+stderr e_a, inv_gamma_pdf, 0.035449, inf;
+stderr e_m, inv_gamma_pdf, 0.008862, inf;
+end;
+
+varobs gp_obs gy_obs;
+
+options_.solve_tolf = 1e-12;
+
+// Metropolis replications are too few, this is only for testing purpose
+estimation(order=1,datafile=fsdat_simul,nobs=192,loglinear,mh_replic=10000,mh_nblocks=1,mh_jscale=0.8);
+
+conditional_forecast_paths;
+var gy_obs;
+periods  1  2  3:5;
+values   0.01 -0.02 0;
+var gp_obs;
+periods 1:5;
+values  0.05;
+end;
+
+conditional_forecast(parameter_set=prior_mode, controlled_varexo=(e_a,e_m));
+
+if ~(exist('OCTAVE_VERSION') && octave_ver_less_than('3.4.0'))
+plot_conditional_forecast(periods=10) gy_obs gp_obs;
+end
+
+conditional_forecast(parameter_set=posterior_mode, controlled_varexo=(e_a,e_m));
+
+if ~(exist('OCTAVE_VERSION') && octave_ver_less_than('3.4.0'))
+plot_conditional_forecast(periods=10) gy_obs gp_obs;
+end
+
+conditional_forecast(parameter_set=posterior_mean, controlled_varexo=(e_a,e_m));
+
+if ~(exist('OCTAVE_VERSION') && octave_ver_less_than('3.4.0'))
+plot_conditional_forecast(periods=10) gy_obs gp_obs;
+end
+
+conditional_forecast(parameter_set=posterior_median, controlled_varexo=(e_a,e_m));
+
+if ~(exist('OCTAVE_VERSION') && octave_ver_less_than('3.4.0'))
+plot_conditional_forecast(periods=10) gy_obs gp_obs;
+end
+
diff --git a/tests/expectations/expectation.mod b/tests/expectations/expectation.mod
new file mode 100644
index 0000000000000000000000000000000000000000..97252145609e1182c2c98b29e53c67dadffdd62b
--- /dev/null
+++ b/tests/expectations/expectation.mod
@@ -0,0 +1,44 @@
+// Example 1 from Collard's guide to Dynare
+var y, c, k, a, h, b;
+varexo e, u;
+
+parameters beta, rho, alpha, delta, theta, psi, tau;
+
+alpha = 0.36;
+rho   = 0.95;
+tau   = 0.025;
+beta  = 0.99;
+delta = 0.025;
+psi   = 0;
+theta = 2.95;
+
+phi   = 0.1;
+
+model;
+c*theta*h^(1+psi)=expectation(1)((1-alpha)*y)+expectation(-2)((1-alpha)*y);
+k = beta*(((exp(b)*c)/(exp(b(+1))*c(+1)))
+    *(exp(b(+1))*alpha*y(+1)+(1-delta)*k));
+y = exp(a)*(k(-1)^alpha)*(h^(1-alpha));
+k = exp(b)*(y-c)+(1-delta)*k(-1);
+a = rho*a(-1)+tau*b(-1) + e;
+b = tau*a(-1)+rho*b(-1) + u;
+end;
+
+initval;
+y = 1.08068253095672;
+c = 0.80359242014163;
+h = 0.29175631001732;
+k = 11.08360443260358;
+a = 0;
+b = 0;
+e = 0;
+u = 0;
+end;
+
+shocks;
+var e; stderr 0.009;
+var u; stderr 0.009;
+var e, u = phi*0.009*0.009;
+end;
+
+stoch_simul(order=1,irf=0);
diff --git a/tests/expectations/expectation_ss.mod b/tests/expectations/expectation_ss.mod
new file mode 100644
index 0000000000000000000000000000000000000000..592905af6eb9f9865f2069152e74dc58a1e1276e
--- /dev/null
+++ b/tests/expectations/expectation_ss.mod
@@ -0,0 +1,46 @@
+// Example 1 from Collard's guide to Dynare
+var y, c, k, a, h, b;
+varexo e, u;
+
+parameters beta, rho, alpha, delta, theta, psi, tau;
+
+alpha = 0.36;
+rho   = 0.95;
+tau   = 0.025;
+beta  = 0.99;
+
+delta = 0.025;
+psi   = 0;
+theta = 2.95;
+
+phi   = 0.1;
+
+model;
+c*theta*h^(1+psi)=(expectation(1)((1-alpha)*y)+expectation(-2)((1-alpha)*y))/2;
+k = beta*(((exp(b)*c)/(exp(b(+1))*c(+1)))
+    *(exp(b(+1))*alpha*y(+1)+(1-delta)*k));
+y = exp(a)*(k(-1)^alpha)*(h^(1-alpha));
+k = exp(b)*(y-c)+(1-delta)*k(-1);
+a = rho*a(-1)+tau*b(-1) + e;
+b = tau*a(-1)+rho*b(-1) + u;
+end;
+
+steady_state_model;
+a = 0;
+b = 0;
+h = 1/3;
+k = ((1/beta-(1-delta))/(alpha*h^(1-alpha)))^(1/(alpha-1));
+y = k^alpha*h^(1-alpha);
+c = y - delta*k;
+theta = (1-alpha)*y/(c*h^(1+psi));
+end;
+
+steady;
+
+shocks;
+var e; stderr 0.009;
+var u; stderr 0.009;
+var e, u = phi*0.009*0.009;
+end;
+
+stoch_simul(order=1,irf=0);
diff --git a/tests/expectations/expectation_ss_old.mod b/tests/expectations/expectation_ss_old.mod
new file mode 100644
index 0000000000000000000000000000000000000000..3d25ba8352c40116f77ab7c96bb31d0757ba233f
--- /dev/null
+++ b/tests/expectations/expectation_ss_old.mod
@@ -0,0 +1,36 @@
+// Example 1 from Collard's guide to Dynare
+var y, c, k, a, h, b;
+varexo e, u;
+
+parameters beta, rho, alpha, delta, theta, psi, tau;
+
+alpha = 0.36;
+rho   = 0.95;
+tau   = 0.025;
+beta  = 0.99;
+
+delta = 0.025;
+psi   = 0;
+theta = 2.95;
+
+phi   = 0.1;
+
+model;
+c*theta*h^(1+psi)=(expectation(1)((1-alpha)*y)+expectation(-2)((1-alpha)*y))/2;
+k = beta*(((exp(b)*c)/(exp(b(+1))*c(+1)))
+    *(exp(b(+1))*alpha*y(+1)+(1-delta)*k));
+y = exp(a)*(k(-1)^alpha)*(h^(1-alpha));
+k = exp(b)*(y-c)+(1-delta)*k(-1);
+a = rho*a(-1)+tau*b(-1) + e;
+b = tau*a(-1)+rho*b(-1) + u;
+end;
+
+steady;
+
+shocks;
+var e; stderr 0.009;
+var u; stderr 0.009;
+var e, u = phi*0.009*0.009;
+end;
+
+stoch_simul(order=1,irf=0);
diff --git a/tests/fs2000/fs2000d.mod b/tests/fs2000/fs2000d.mod
new file mode 100644
index 0000000000000000000000000000000000000000..04f1228053ffe5d1f28878be21e29dc6bb4f28f0
--- /dev/null
+++ b/tests/fs2000/fs2000d.mod
@@ -0,0 +1,75 @@
+// See fs2000.mod in the examples/ directory for details on the model
+
+var m P c e W R k d n l gy_obs gp_obs y dA;
+varexo e_a e_m;
+
+parameters alp bet gam mst rho psi del;
+
+alp = 0.33;
+bet = 0.99;
+gam = 0.003;
+mst = 1.011;
+rho = 0.7;
+psi = 0.787;
+del = 0.02;
+
+model;
+dA = exp(gam+e_a);
+log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m;
+-P/(c(+1)*P(+1)*m)+bet*P(+1)*(alp*exp(-alp*(gam+log(e(+1))))*k^(alp-1)*n(+1)^(1-alp)+(1-del)*exp(-(gam+log(e(+1)))))/(c(+2)*P(+2)*m(+1))=0;
+W = l/n;
+-(psi/(1-psi))*(c*P/(1-n))+l/n = 0;
+R = P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(-alp)/W;
+1/(c*P)-bet*P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)/(m*l*c(+1)*P(+1)) = 0;
+c+k = exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)+(1-del)*exp(-(gam+e_a))*k(-1);
+P*c = m;
+m-1+d = l;
+e = exp(e_a);
+y = k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a));
+gy_obs = dA*y/y(-1);
+gp_obs = (P/P(-1))*m(-1)/dA;
+end;
+
+initval;
+k = 6;
+m = mst;
+P = 2.25;
+c = 0.45;
+e = 1;
+W = 4;
+R = 1.02;
+d = 0.85;
+n = 0.19;
+l = 0.86;
+y = 0.6;
+gy_obs = exp(gam);
+gp_obs = exp(-gam);
+dA = exp(gam);
+end;
+
+shocks;
+var e_a; stderr 0.014;
+var e_m; stderr 0.005;
+end;
+
+steady;
+
+check;
+
+estimated_params;
+alp, beta_pdf, 0.356, 0.02;
+bet, beta_pdf, 0.993, 0.002;
+gam, normal_pdf, 0.0085, 0.003;
+mst, normal_pdf, 1.0002, 0.007;
+rho, beta_pdf, 0.129, 0.223;
+psi, beta_pdf, 0.65, 0.05;
+del, beta_pdf, 0.01, 0.005;
+stderr e_a, inv_gamma_pdf, 0.035449, inf;
+stderr e_m, inv_gamma_pdf, 0.008862, inf;
+end;
+
+varobs gp_obs gy_obs;
+
+options_.solve_tolf = 1e-12;
+
+estimation(order=1,datafile=fsdat_simul,nobs=192,mode_compute=5,loglinear,mh_replic=0);
diff --git a/tests/gsa/ls2003.mod b/tests/gsa/ls2003.mod
new file mode 100644
index 0000000000000000000000000000000000000000..fd46ed043f2686578e99e3c299beeefa1f525d5f
--- /dev/null
+++ b/tests/gsa/ls2003.mod
@@ -0,0 +1,224 @@
+var y y_s R pie dq pie_s de A y_obs pie_obs R_obs;
+varexo e_R e_q e_ys e_pies e_A;
+
+parameters psi1 psi2 psi3 rho_R tau alpha rr k rho_q rho_A rho_ys rho_pies;
+
+psi1 = 1.54;
+psi2 = 0.25;
+psi3 = 0.25;
+rho_R = 0.5;
+alpha = 0.3;
+rr = 2.51;
+k = 0.5;
+tau = 0.5;
+rho_q = 0.4;
+rho_A = 0.2;
+rho_ys = 0.9;
+rho_pies = 0.7;
+
+
+model(linear);
+y = y(+1) - (tau +alpha*(2-alpha)*(1-tau))*(R-pie(+1))-alpha*(tau +alpha*(2-alpha)*(1-tau))*dq(+1) + alpha*(2-alpha)*((1-tau)/tau)*(y_s-y_s(+1))-A(+1);
+pie = exp(-rr/400)*pie(+1)+alpha*exp(-rr/400)*dq(+1)-alpha*dq+(k/(tau+alpha*(2-alpha)*(1-tau)))*y+k*alpha*(2-alpha)*(1-tau)/(tau*(tau+alpha*(2-alpha)*(1-tau)))*y_s;
+pie = de+(1-alpha)*dq+pie_s;
+R = rho_R*R(-1)+(1-rho_R)*(psi1*pie+psi2*(y+alpha*(2-alpha)*((1-tau)/tau)*y_s)+psi3*de)+e_R;
+dq = rho_q*dq(-1)+e_q;
+y_s = rho_ys*y_s(-1)+e_ys;
+pie_s = rho_pies*pie_s(-1)+e_pies;
+A = rho_A*A(-1)+e_A;
+y_obs = y-y(-1)+A;
+pie_obs = 4*pie;
+R_obs = 4*R;
+end;
+
+shocks;
+	var e_R = 1.25^2;
+	var e_q = 2.5^2;
+	var e_A = 1.89;
+	var e_ys = 1.89;
+	var e_pies = 1.89;
+end;
+
+varobs y_obs R_obs pie_obs dq de;
+
+estimated_params;
+psi1 , gamma_pdf,1.5,0.5;
+psi2 , gamma_pdf,0.25,0.125;
+psi3 , gamma_pdf,0.25,0.125;
+rho_R ,beta_pdf,0.5,0.2;
+alpha ,beta_pdf,0.3,0.1;
+rr ,gamma_pdf,2.5,1;
+k , gamma_pdf,0.5,0.25;
+tau ,gamma_pdf,0.5,0.2;
+rho_q ,beta_pdf,0.4,0.2;
+rho_A ,beta_pdf,0.5,0.2;
+rho_ys ,beta_pdf,0.8,0.1;
+rho_pies,beta_pdf,0.7,0.15;
+stderr e_R,inv_gamma_pdf,1.2533,0.6551;
+stderr e_q,inv_gamma_pdf,2.5066,1.3103;
+stderr e_A,inv_gamma_pdf,1.2533,0.6551;
+stderr e_ys,inv_gamma_pdf,1.2533,0.6551;
+stderr e_pies,inv_gamma_pdf,1.88,0.9827;
+end;
+
+  
+disp(' ');
+disp('NOW I DO STABILITY MAPPING and prepare sample for Reduced form Mapping');
+disp(' ');
+disp('Press ENTER to continue'); pause(5);
+
+dynare_sensitivity(redform=1, //create sample of reduced form coefficients
+alpha2_stab=0.4,
+ksstat=0);
+// NOTE: since namendo is emppty by default, 
+// this call does not perform the mapping of reduced form coefficient: just prepares the sample
+
+/*
+disp(' ');
+disp('ANALYSIS OF REDUCED FORM COEFFICIENTS');
+disp(' ');
+disp('Press ENTER to continue'); pause(5);
+
+dynare_sensitivity(load_stab=1,  // loead previously generated sample analysed for stability
+redform=1,  // do the reduced form mapping
+logtrans_redform=1,  // estimate log-transformed reduced form coefficients (default=0)
+namendo=(pie,R),  // evaluate relationships for pie and R (namendo=(:) for all variables)
+namexo=(e_R),     // evaluate relationships with exogenous e_R (use namexo=(:) for all shocks)
+namlagendo=(R),   // evaluate relationships with lagged R (use namlagendo=(:) for all lagged endogenous)
+stab=0 // don't repeat again the stability mapping
+);
+*/
+
+disp(' ');
+disp('THE PREVIOUS TWO CALLS COULD BE DONE TOGETHER');
+disp('BY USING THE COMBINED CALL');
+disp(' ');
+disp('dynare_sensitivity(alpha2_stab=0.4, ksstat=0, redform=1,')
+disp('logtrans_redform=1, namendo=(pie,R), namexo=(e_R), namlagendo=(R));')
+disp(' ');
+disp('Press ENTER to continue'); pause(5);
+//dynare_sensitivity(
+//alpha2_stab=0.4,
+//ksstat=0,
+//redform=1, //create sample of reduced form coefficients
+//logtrans_redform=1,  // estimate log-transformed reduced form coefficients (default=0)
+//namendo=(pie,R),  // evaluate relationships for pie and R (namendo=(:) for all variables)
+//namexo=(e_R),     // evaluate relationships with exogenous e_R (use namexo=(:) for all shocks)
+//namlagendo=(R)   // evaluate relationships with lagged R (use namlagendo=(:) for all lagged endogenous)
+//);
+
+
+
+disp(' ');
+disp('MC FILTERING(rmse=1), TO MAP THE FIT FROM PRIORS');
+disp('Press ENTER to continue'); pause(5);
+
+dynare_sensitivity(datafile=data_ca1,first_obs=8,nobs=79,prefilter=1, // also presample=2,loglinear, are admissible
+load_stab=1,     // load prior sample
+istart_rmse=2,   //start computing rmse from second observation (i.e. rmse does not inlude initial big error)
+stab=0,          // don't  plot again stability analysis results
+rmse=1,          // do rmse analysis
+pfilt_rmse=0.1,  // filtering criterion, i.e. filter the best 10 percent rmse's
+alpha2_rmse=0.3, // critical value for correlations in the rmse filterting analysis: 
+                 // if ==1, means no corrleation analysis done
+alpha_rmse=1     // critical value for smirnov statistics of filtered samples
+                 // if ==1, no smornov analysis is done
+);
+
+disp(' ');
+disp('THE PREVIOUS THREE CALLS COULD BE DONE TOGETHER');
+disp('BY USING THE COMBINED CALL');
+disp(' ');
+disp('dynare_sensitivity(alpha2_stab=0.4, ksstat=0, redform=1,')
+disp('logtrans_redform=1, namendo=(pie,R), namexo=(e_R), namlagendo=(R),')   
+disp('datafile=data_ca1,first_obs=8,nobs=79,prefilter=1,')
+disp('istart_rmse=2, rmse=1, pfilt_rmse=0.1, alpha2_rmse=0.3, alpha_rmse=1);')
+disp(' ');
+disp('Press ENTER to continue'); pause(5);
+//dynare_sensitivity(
+//alpha2_stab=0.4,
+//ksstat=0,
+//redform=1, //create sample of reduced form coefficients
+//logtrans_redform=1,  // estimate log-transformed reduced form coefficients (default=0)
+//namendo=(pie,R),  // evaluate relationships for pie and R (namendo=(:) for all variables)
+//namexo=(e_R),     // evaluate relationships with exogenous e_R (use namexo=(:) for all shocks)
+//namlagendo=(R),   // evaluate relationships with lagged R (use namlagendo=(:) for all lagged endogenous)
+//datafile=data_ca1,first_obs=8,nobs=79,prefilter=1, 
+//istart_rmse=2,   //start computing rmse from second observation (i.e. rmse does not inlude initial big error)
+//rmse=1,          // do rmse analysis
+//pfilt_rmse=0.1,  // filtering criterion, i.e. filter the best 10 percent rmse's
+//alpha2_rmse=0.3, // critical value for correlations in the rmse filterting analysis: 
+//                 // if ==1, means no corrleation analysis done
+//alpha_rmse=1     // critical value for smirnov statistics of filtered samples
+//                 // if ==1, no smirnov sensitivity analysis is done
+//);
+
+
+
+disp(' ');
+disp('I ESTIMATE THE MODEL');
+disp(' ');
+disp('Press ENTER to continue'); pause(5);
+
+// run this to generate posterior mode and Metropolis files if not yet done
+//estimation(datafile=data_ca1,first_obs=8,nobs=79,mh_nblocks=2,
+//   prefilter=1,mh_jscale=0.5,mh_replic=100000, mode_compute=4, nograph, mh_drop=0.6);
+
+
+// run this to produce posterior samples of filtered, smoothed and irf variables, if not yet done
+//estimation(datafile=data_ca1,first_obs=8,nobs=79,mh_nblocks=2,prefilter=1,mh_jscale=0.3,
+//          mh_replic=0, mode_file=ls2003_mode, mode_compute=0, nograph, load_mh_file, bayesian_irf,
+//		  filtered_vars, smoother, mh_drop=0.6);
+
+
+disp(' ');
+disp('WE DO STABILITY MAPPING AGAIN, BUT FOR MULTIVARIATE SAMPLE AT THE POSTERIOR MODE (or ML) and Hessian (pprior=0 & ppost=0)');
+disp('Typical for ML estimation, also feasible for posterior mode');
+disp(' ');
+disp('Press ENTER to continue'); pause(5);
+
+dynare_sensitivity(pprior=0,Nsam=2048,alpha2_stab=0.4,
+mode_file=ls2003_mode  // specifies the mode file where the mode and Hessian are stored
+);
+
+
+disp(' ');
+disp('RMSE ANALYSIS FOR MULTIVARIATE SAMPLE AT THE POSTERIOR MODE');
+disp(' ');
+disp('Press ENTER to continue'); pause(5);
+dynare_sensitivity(mode_file=ls2003_mode,
+datafile=data_ca1,first_obs=8,nobs=79,prefilter=1,
+pprior=0,
+stab=0,
+rmse=1,
+alpha2_rmse=1, // no correlation analysis
+alpha_rmse=1  // no Smirnov sensitivity analysis
+);
+
+disp(' ');
+disp('THE LAST TWO CALLS COULD BE DONE TOGETHER');
+disp('BY USING THE COMBINED CALL');
+disp(' ');
+disp('dynare_sensitivity(pprior=0,Nsam=2048,alpha2_stab=0.4,mode_file=ls2003_mode,')
+disp('datafile=data_ca1,first_obs=8,nobs=79,prefilter=1,')
+disp('rmse=1, alpha2_rmse=1, alpha_rmse=1);')
+disp(' ');
+disp('Press ENTER to continue'); pause(5);
+//dynare_sensitivity(pprior=0,Nsam=2048,alpha2_stab=0.4,mode_file=ls2003_mode,
+//datafile=data_ca1,first_obs=8,nobs=79,prefilter=1,
+//rmse=1,
+//alpha2_rmse=1, // no correlation analysis
+//alpha_rmse=1  // no Smirnov sensitivity analysis
+//);
+
+disp(' ');
+disp('RMSE ANALYSIS FOR POSTERIOR MCMC sample (ppost=1)');
+disp('Needs a call to dynare_estimation to load all MH environment');
+disp('Press ENTER to continue'); pause(5);
+estimation(datafile=data_ca1,first_obs=8,nobs=79,mh_nblocks=2, mode_file=ls2003_mode, load_mh_file,
+  prefilter=1,mh_jscale=0.5,mh_replic=0, mode_compute=0, nograph, mh_drop=0.6);
+
+dynare_sensitivity(stab=0, // no need for stability analysis since the posterior sample is surely OK
+rmse=1,ppost=1,alpha2_rmse=1,alpha_rmse=1);
+
+
diff --git a/tests/printMakeCheckMatlabErrMsg.m b/tests/printMakeCheckMatlabErrMsg.m
new file mode 100644
index 0000000000000000000000000000000000000000..a72b55401be47b9e76475d9c69f66f18d809dfb2
--- /dev/null
+++ b/tests/printMakeCheckMatlabErrMsg.m
@@ -0,0 +1,9 @@
+function printMakeCheckMatlabErrMsg(modfilename, exception)
+    fprintf('\n********************************************\n');
+    disp('*** DYNARE-TEST-MATLAB ERROR ENCOUNTERED ***');
+    disp('********************************************');
+    disp(['  WHILE RUNNING MODFILE: ' modfilename]);
+    fprintf('\n');
+    disp(getReport(exception));
+    fprintf('*************************************\n\n\n');
+end
diff --git a/tests/printMakeCheckOctaveErrMsg.m b/tests/printMakeCheckOctaveErrMsg.m
new file mode 100644
index 0000000000000000000000000000000000000000..84e19d00edf1827559ee954ca16bc9eda6875e3e
--- /dev/null
+++ b/tests/printMakeCheckOctaveErrMsg.m
@@ -0,0 +1,14 @@
+function printMakeCheckOctaveErrMsg(modfilename, err)
+    printf("\n");
+    printf("********************************************\n");
+    printf("*** DYNARE-TEST-OCTAVE ERROR ENCOUNTERED ***\n");
+    printf("********************************************\n");
+    printf("  WHILE RUNNING MODFILE: %s\n", modfilename);
+    printf("                    MSG: %s\n", err.message);
+    if (isfield(err, 'stack'))
+        printf("                IN FILE: %s\n", err.stack(1).file);
+        printf("            IN FUNCTION: %s\n", err.stack(1).name);
+        printf("     ON LINE and COLUMN: %d and %d\n",err.stack(1).line,err.stack(1).column);
+    end
+    printf("*************************************\n\n\n");
+end
diff --git a/tests/run_test_matlab.m b/tests/run_test_matlab.m
new file mode 100644
index 0000000000000000000000000000000000000000..6a9cdedd7814bbeea716c8c2006de200e09094ec
--- /dev/null
+++ b/tests/run_test_matlab.m
@@ -0,0 +1,188 @@
+% Copyright (C) 2011 Dynare Team
+%
+% This file is part of Dynare.
+%
+% Dynare is free software: you can redistribute it and/or modify
+% it under the terms of the GNU General Public License as published by
+% the Free Software Foundation, either version 3 of the License, or
+% (at your option) any later version.
+%
+% Dynare is distributed in the hope that it will be useful,
+% but WITHOUT ANY WARRANTY; without even the implied warranty of
+% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+% GNU General Public License for more details.
+%
+% You should have received a copy of the GNU General Public License
+% along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
+
+% Implementation notes:
+%
+% Before every call to Dynare, the contents of the workspace is saved in
+% 'wsMat.mat', and reloaded after Dynare has finished (this is necessary since
+% Dynare does a 'clear -all').  Also note that we take care of clearing the
+% 'exception' variable in all 'catch' block, because otherwise the next 'load
+% wsMat' within a 'catch' block will overwrite the last exception.
+
+top_test_dir = pwd;
+addpath(top_test_dir);
+addpath([top_test_dir '/../matlab']);
+
+% Test Dynare Version
+if ~strcmp(dynare_version(), getenv('DYNARE_VERSION'))
+  error('Incorrect version of Dynare is being tested')
+end
+
+% Test MOD files listed in Makefile.am
+name=getenv('MODFILES');
+num_modfiles = 0;
+
+failedBase = {};
+
+while ~isempty(name)
+    [modfile, name] = strtok(name);
+    num_modfiles = num_modfiles + 1;
+    [directory, testfile, ext] = fileparts([top_test_dir '/' modfile]);
+    cd(directory);
+    disp('');
+    disp(['***  TESTING: ' modfile ' ***']);
+    try
+        save wsMat
+        dynare([testfile ext],'console')
+        load wsMat
+    catch exception
+        load wsMat
+        failedBase{size(failedBase,2)+1} = modfile;
+        printMakeCheckMatlabErrMsg(modfile, exception);
+        clear exception
+    end
+    close all
+    delete('wsMat.mat')
+
+    cd(top_test_dir);
+end
+
+% Test block_bytecode/ls2003.mod with various combinations of
+% block/bytecode/solve_algo/stack_solve_algo
+failedBlock = {};
+num_block_tests = 0;
+cd([top_test_dir '/block_bytecode']);
+for blockFlag = 0:1
+    for bytecodeFlag = 0:1
+        default_solve_algo = 2;
+        default_stack_solve_algo = 0;
+        if ~blockFlag && ~bytecodeFlag
+            solve_algos = 1:4;
+            stack_solve_algos = 0;
+        elseif blockFlag && ~bytecodeFlag
+            solve_algos = [1:4 6:8];
+            stack_solve_algos = 0:4;
+        else
+            solve_algos = 1:8;
+            stack_solve_algos = 0:5;
+        end
+        if license('test', 'optimization_toolbox')
+            solve_algos = [ solve_algos 0 ];
+        end
+
+        for i = 1:length(solve_algos)
+            num_block_tests = num_block_tests + 1;
+            if ~blockFlag && ~bytecodeFlag && (i == 1)
+                % This is the reference simulation path against which all
+                % other simulations will be tested
+                try
+                    save wsMat
+                    run_ls2003(blockFlag, bytecodeFlag, solve_algos(i), default_stack_solve_algo)
+                    load wsMat
+                    y_ref = oo_.endo_simul;
+                    save('test.mat','y_ref');
+                catch exception
+                    load wsMat
+                    failedBlock{size(failedBlock,2)+1} = ['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'];
+                    printMakeCheckMatlabErrMsg(['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'], exception);
+                    clear exception
+                end
+            else
+                try
+                    save wsMat
+                    run_ls2003(blockFlag, bytecodeFlag, solve_algos(i), default_stack_solve_algo)
+                    load wsMat
+                    % Test against the reference simulation path
+                    load('test.mat','y_ref');
+                    diff = oo_.endo_simul - y_ref;
+                    if(abs(diff) > options_.dynatol)
+                        failedBlock{size(failedBlock,2)+1} = ['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'];
+                        exception = MException('ERROR: simulation path differs from the reference path');
+                        printMakeCheckMatlabErrMsg(['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'], exception);
+                        clear exception
+                    end
+                catch exception
+                    load wsMat
+                    failedBlock{size(failedBlock,2)+1} = ['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'];
+                    printMakeCheckMatlabErrMsg(['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'], exception);
+                    clear exception
+                end
+            end
+        end
+        for i = 1:length(stack_solve_algos)
+            num_block_tests = num_block_tests + 1;
+            try
+                save wsMat
+                run_ls2003(blockFlag, bytecodeFlag, default_solve_algo, stack_solve_algos(i))
+                load wsMat
+                % Test against the reference simulation path
+                load('test.mat','y_ref');
+                diff = oo_.endo_simul - y_ref;
+                if(abs(diff) > options_.dynatol)
+                    failedBlock{size(failedBlock,2)+1} = ['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(default_solve_algo) ', ' num2str(stack_solve_algos(i)) ')'];
+                    exception = MException('ERROR: simulation path difers from the reference path');
+                    printMakeCheckMatlabErrMsg(['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(default_solve_algo) ', ' num2str(stack_solve_algos(i)) ')'], exception);
+                    clear exception
+                end
+            catch exception
+                load wsMat
+                failedBlock{size(failedBlock,2)+1} = ['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'];
+                printMakeCheckMatlabErrMsg(['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'], exception);
+                clear exception
+            end
+        end
+    end
+end
+delete('wsMat.mat')
+
+cd(top_test_dir);
+
+total_tests = num_modfiles+num_block_tests;
+
+% print output to screen and to file
+fid = fopen('run_test_matlab_output.txt', 'w');
+
+fprintf('\n\n\n');
+fprintf(fid,'\n\n\n');
+disp('***************************************');
+fprintf(fid,'***************************************\n');
+disp('*         DYNARE TEST RESULTS         *');
+fprintf(fid,'*         DYNARE TEST RESULTS         *\n');
+disp('*        for make check-matlab        *');
+fprintf(fid,'*        for make check-matlab        *\n');
+disp('***************************************');
+fprintf(fid,'***************************************\n');
+disp(['  ' num2str(total_tests-size(failedBase,2)-size(failedBlock,2)) ' tests PASSED out of ' num2str(total_tests) ' tests run']);
+fprintf(fid,' %d tests PASSED out of %d tests run\n', total_tests-size(failedBase,2)-size(failedBlock,2), total_tests);
+disp('***************************************');
+fprintf(fid,'***************************************\n');
+if size(failedBase,2) > 0 || size(failedBlock,2) > 0
+    disp(['List of ' num2str(size(failedBase,2)+size(failedBlock,2)) ' tests FAILED:']);
+    fprintf(fid,'List of %d tests FAILED:\n', size(failedBase,2)+size(failedBlock,2));
+    for i=1:size(failedBase,2)
+        disp(['   * ' failedBase{i}]);
+        fprintf(fid,'   * %s\n', failedBase{i});
+    end
+    for i=1:size(failedBlock,2)
+        disp(['   * ' failedBlock{i}]);
+        fprintf(fid,'   * %s\n', failedBlock{i});
+    end
+    fprintf('***************************************\n\n');
+    fprintf(fid,'***************************************\n\n');
+end
+fclose(fid);
+exit;
diff --git a/tests/run_test_octave.m b/tests/run_test_octave.m
index dc15ab08643a68c250cd6bfbd0658e327882e414..5f7bb9654283d49df96acfae55d09cf4531df4db 100644
--- a/tests/run_test_octave.m
+++ b/tests/run_test_octave.m
@@ -1,4 +1,4 @@
-## Copyright (C) 2009 Dynare Team
+## Copyright (C) 2009-2011 Dynare Team
 ##
 ## This file is part of Dynare.
 ##
@@ -15,28 +15,170 @@
 ## You should have received a copy of the GNU General Public License
 ## along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
 
-## First argument is MOD file (possibly with a path prefix)
-## Second argument is path to Dynare installation to be checked
-## Third argument is Dynare version to be checked
+## Implementation notes:
+##
+## Before every call to Dynare, the contents of the workspace is saved in
+## 'wsOct', and reloaded after Dynare has finished (this is necessary since
+## Dynare does a 'clear -all').
+
+top_test_dir = pwd;
+addpath(top_test_dir);
+addpath([top_test_dir '/../matlab']);
+
+## Test Dynare Version
+if !strcmp(dynare_version(), getenv("DYNARE_VERSION"))
+  error("Incorrect version of Dynare is being tested")
+endif
 
 ## Ask gnuplot to create graphics in text mode
 ## Note that setenv() was introduced in Octave 3.0.2, for compatibility
 ## with MATLAB
 putenv("GNUTERM", "dumb")
 
-[directory, name, ext] = fileparts(argv(){1});
+## Test MOD files listed in Makefile.am
+name = strsplit(getenv("MODFILES"), " ");
 
-printf("TEST: %s...\n", name)
+failedBase = {};
 
-addpath(argv(){2})
+for i=1:size(name,2)
+  [directory, testfile, ext] = fileparts([top_test_dir '/' name{i}]);
+  cd(directory);
+  printf("\n***  TESTING: %s ***\n", name{i});
+  try
+    save wsOct
+    dynare([testfile ext])
+    load wsOct
+  catch
+    load wsOct
+    failedBase{size(failedBase,2)+1} = name{i};
+    printMakeCheckOctaveErrMsg(name{i}, lasterror);
+  end_try_catch
+  delete('wsOct');
+  cd(top_test_dir);
+end
 
-if !strcmp(dynare_version(), argv(){3})
-  error("Incorrect version of Dynare is being tested")
-endif
+## Test block_bytecode/ls2003.mod with various combinations of
+## block/bytecode/solve_algo/stack_solve_algo
+failedBlock = {};
+num_block_tests = 0;
+cd([top_test_dir '/block_bytecode']);
+for blockFlag = 0:1
+  for bytecodeFlag = 0:1
+    ## Recall that solve_algo=7 and stack_solve_algo=2 are not supported
+    ## under Octave
+    default_solve_algo = 2;
+    default_stack_solve_algo = 0;
+    if !blockFlag && !bytecodeFlag
+      solve_algos = 0:4;
+      stack_solve_algos = 0;
+    elseif blockFlag && !bytecodeFlag
+      solve_algos = [0:4 6 8];
+      stack_solve_algos = [0 1 3 4];
+    else
+      solve_algos = [0:6 8];
+      stack_solve_algos = [0 1 3:5];
+    endif
+
+    for i = 1:length(solve_algos)
+      num_block_tests = num_block_tests + 1;
+      if !blockFlag && !bytecodeFlag && (i == 1)
+        ## This is the reference simulation path against which all
+        ## other simulations will be tested
+        try
+          save wsOct
+          run_ls2003(blockFlag, bytecodeFlag, solve_algos(i), default_stack_solve_algo)
+          load wsOct
+          y_ref = oo_.endo_simul;
+          save('test.mat','y_ref');
+        catch
+          load wsOct
+          failedBlock{size(failedBlock,2)+1} = ['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'];
+          printMakeCheckOctaveErrMsg(['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'], lasterror);
+        end_try_catch
+      else
+        try
+          save wsOct
+          run_ls2003(blockFlag, bytecodeFlag, solve_algos(i), default_stack_solve_algo)
+          load wsOct
+          ## Test against the reference simulation path
+          load('test.mat','y_ref');
+          diff = oo_.endo_simul - y_ref;
+          if(abs(diff) > options_.dynatol)
+            failedBlock{size(failedBlock,2)+1} = ['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'];
+            differr.message = ["ERROR: simulation path differs from the reference path" ];
+            printMakeCheckOctaveErrMsg(['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'], differr);
+          endif
+        catch
+          load wsOct
+          failedBlock{size(failedBlock,2)+1} = ['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'];
+          printMakeCheckOctaveErrMsg(['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'], lasterror);
+        end_try_catch
+      endif
+    endfor
+    for i = 1:length(stack_solve_algos)
+      num_block_tests = num_block_tests + 1;
+      try
+        save wsOct
+        run_ls2003(blockFlag, bytecodeFlag, default_solve_algo, stack_solve_algos(i))
+        load wsOct
+        ## Test against the reference simulation path
+        load('test.mat','y_ref');
+        diff = oo_.endo_simul - y_ref;
+        if(abs(diff) > options_.dynatol)
+          failedBlock{size(failedBlock,2)+1} = ['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(default_solve_algo) ', ' num2str(stack_solve_algos(i)) ')'];
+          differr.message = ["ERROR: simulation path differs from the reference path" ];
+          printMakeCheckOctaveErrMsg(['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(default_solve_algo) ', ' num2str(stack_solve_algos(i)) ')'], differr);
+        endif
+      catch
+        load wsOct
+        failedBlock{size(failedBlock,2)+1} = ['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'];
+        printMakeCheckOctaveErrMsg(['block_bytecode/run_ls2003.m(' num2str(blockFlag) ', ' num2str(bytecodeFlag) ', ' num2str(solve_algos(i)) ', ' num2str(default_stack_solve_algo) ')'], lasterror);
+      end_try_catch
+    endfor
+  endfor
+endfor
+
+delete('wsOct');
 
-cd(directory)
+cd(top_test_dir);
 
-dynare(name)
+total_tests = size(name,2)+num_block_tests;
+
+% print output to screen and to file
+fid = fopen("run_test_octave_output.txt", "w");
+
+printf("\n\n\n");
+fprintf(fid,'\n\n\n');
+printf("***************************************\n");
+fprintf(fid,"***************************************\n");
+printf("*         DYNARE TEST RESULTS         *\n");
+fprintf(fid,"*         DYNARE TEST RESULTS         *\n");
+printf("*        for make check-octave        *\n");
+fprintf(fid,"*        for make check-octave        *\n");
+printf("***************************************\n");
+fprintf(fid,"***************************************\n");
+printf("  %d tests PASSED out of %d tests run\n", total_tests-size(failedBase,2)-size(failedBlock,2), total_tests);
+fprintf(fid," %d tests PASSED out of %d tests run\n", total_tests-size(failedBase,2)-size(failedBlock,2), total_tests);
+printf("***************************************\n");
+fprintf(fid,"***************************************\n");
+if size(failedBase,2) > 0 || size(failedBlock,2) > 0
+  printf("List of %d tests FAILED:\n", size(failedBase,2)+size(failedBlock,2));
+  fprintf(fid,"List of %d tests FAILED:\n", size(failedBase,2)+size(failedBlock,2));
+  for i=1:size(failedBase,2)
+    printf("   * %s\n",failedBase{i});
+    fprintf(fid,"   * %s\n", failedBase{i});
+  end
+  for i=1:size(failedBlock,2)
+    printf("   * %s\n",failedBlock{i});
+    fprintf(fid,"   * %s\n", failedBlock{i});
+  end
+  printf("***************************************\n\n");
+  fprintf(fid,"***************************************\n\n");
+  fclose(fid);
+  error("make check-octave FAILED");
+else
+  fclose(fid);
+endif
 
 ## Local variables:
 ## mode: Octave
diff --git a/tests/steady_state/multi_leads.mod b/tests/steady_state/multi_leads.mod
new file mode 100644
index 0000000000000000000000000000000000000000..8d6d3c8774635298e1cbad5df89c56c040753a19
--- /dev/null
+++ b/tests/steady_state/multi_leads.mod
@@ -0,0 +1,27 @@
+    @#define N = 4
+    var A B;
+    varexo eB;
+    parameters Bss rho;
+
+    Bss=1;
+    rho=0.9;
+
+    model;
+    A = B(@{N});
+    B = (1-rho)*Bss +rho*B(-1) +eB;
+    end;
+
+    steady_state_model;
+    A = Bss;
+    B = Bss;
+    end;
+
+    resid;
+    steady;
+    check;
+
+    shocks;
+    var eB ; stderr 1 ;
+    end;
+
+    stoch_simul(order=1,irf=0) ;
diff --git a/tests/steady_state/walsh1_initval.mod b/tests/steady_state/walsh1_initval.mod
new file mode 100644
index 0000000000000000000000000000000000000000..b99016724fbe995a02525f7e2adabdc94a180f6d
--- /dev/null
+++ b/tests/steady_state/walsh1_initval.mod
@@ -0,0 +1,68 @@
+var y c k m n R pi z u;
+varexo	e sigma;	 
+// sigma stands for phi in the eq 2.37 p.69
+
+parameters alpha beta delta gamm phi1 eta a b rho  phi2 Psi thetass;  
+//phi1 stands for capital phi in eq.2.68 and 2.69
+//phi2 stands for lowercase phi in eq. 2.66
+
+change_type(var) Psi;
+change_type(parameters) n;
+
+alpha = 0.36;
+beta = 0.989; 
+gamm = 0.5;
+delta = 0.019;
+phi1 = 2;
+phi2 = 0;
+eta = 1;
+a = 0.95;
+b = 2.56;
+rho = 0.95;
+//Psi = 1.47630583;
+thetass = 1.0125;
+
+model;
+
+(a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*a*exp(c)^(-b) = (a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*(1-a)*exp(m)^(-b)+beta*(a*exp(c(+1))^(1-b)+(1-a)*exp(m(+1))^(1-b))^((b-phi1)/(1-b))*a*exp(c(+1))^(-b)/(1+pi(+1));
+
+Psi*(1-exp(n))^(-eta)/(a*exp(c)^(-b)*(a*exp(c)^(1-b) + (1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))) = (1-alpha)*exp(y)/exp(n);
+
+(a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*a*exp(c)^(-b) = beta*exp(R(+1))*(a*exp(c(+1))^(1-b)+(1-a)*exp(m(+1))^(1-b))^((b-phi1)/(1-b))*a*exp(c(+1))^(-b);
+
+exp(R) = alpha*exp(y)/exp(k(-1)) + 1-delta;
+
+exp(k) = (1-delta)*exp(k(-1))+exp(y)-exp(c);
+
+exp(y) = exp(z)*exp(k(-1))^alpha*exp(n)^(1-alpha);
+
+exp(m) = exp(m(-1))*(u+thetass)/(1+pi);
+
+z = rho*z(-1) + e;
+
+u = gamm*u(-1) + phi2*z(-1) + sigma;
+
+end;
+
+shocks;
+var e; stderr 0.007;
+var sigma;stderr 0.0089;
+end;
+
+n=log(1/3);
+
+initval;
+y =  0.2;
+c =  0.03;
+k =  2.7;
+m =  0.3;
+Psi = 1.5;
+R =  0.01;
+pi = 0.01;
+z = 0;
+u = 0;
+end;
+
+
+steady;
+
diff --git a/tests/steady_state/walsh1_old_ss.mod b/tests/steady_state/walsh1_old_ss.mod
new file mode 100644
index 0000000000000000000000000000000000000000..39eb8910d83e9406554a5891482f923ab21837e5
--- /dev/null
+++ b/tests/steady_state/walsh1_old_ss.mod
@@ -0,0 +1,51 @@
+var y c k m n R pi z u;
+varexo	e sigma;	 
+// sigma stands for phi in the eq 2.37 p.69
+
+parameters alpha beta delta gamm phi1 eta a b rho  phi2 Psi thetass;  
+//phi1 stands for capital phi in eq.2.68 and 2.69
+//phi2 stands for lowercase phi in eq. 2.66
+
+alpha = 0.36;
+beta = 0.989; 
+gamm = 0.5;
+delta = 0.019;
+phi1 = 2;
+phi2 = 0;
+eta = 1;
+a = 0.95;
+b = 2.56;
+rho = 0.95;
+Psi = 1.47630583;
+thetass = 1.0125;
+
+model;
+
+(a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*a*exp(c)^(-b) = (a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*(1-a)*exp(m)^(-b)+beta*(a*exp(c(+1))^(1-b)+(1-a)*exp(m(+1))^(1-b))^((b-phi1)/(1-b))*a*exp(c(+1))^(-b)/(1+pi(+1));
+
+Psi*(1-exp(n))^(-eta)/(a*exp(c)^(-b)*(a*exp(c)^(1-b) + (1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))) = (1-alpha)*exp(y)/exp(n);
+
+(a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*a*exp(c)^(-b) = beta*exp(R(+1))*(a*exp(c(+1))^(1-b)+(1-a)*exp(m(+1))^(1-b))^((b-phi1)/(1-b))*a*exp(c(+1))^(-b);
+
+exp(R) = alpha*exp(y)/exp(k(-1)) + 1-delta;
+
+exp(k) = (1-delta)*exp(k(-1))+exp(y)-exp(c);
+
+exp(y) = exp(z)*exp(k(-1))^alpha*exp(n)^(1-alpha);
+
+exp(m) = exp(m(-1))*(u+thetass)/(1+pi);
+
+z = rho*z(-1) + e;
+
+u = gamm*u(-1) + phi2*z(-1) + sigma;
+
+end;
+
+shocks;
+var e; stderr 0.007;
+var sigma;stderr 0.0089;
+end;
+
+
+steady;
+
diff --git a/tests/steady_state/walsh1_ssm.mod b/tests/steady_state/walsh1_ssm.mod
new file mode 100644
index 0000000000000000000000000000000000000000..8d1de749e0f4ba2449e0ca4cc7d010b430632871
--- /dev/null
+++ b/tests/steady_state/walsh1_ssm.mod
@@ -0,0 +1,72 @@
+var y c k m n R pi z u;
+varexo	e sigma;	 
+// sigma stands for phi in the eq 2.37 p.69
+
+parameters alpha beta delta gamm phi1 eta a b rho  phi2 Psi thetass;  
+//phi1 stands for capital phi in eq.2.68 and 2.69
+//phi2 stands for lowercase phi in eq. 2.66
+
+alpha = 0.36;
+beta = 0.989; 
+gamm = 0.5;
+delta = 0.019;
+phi1 = 2;
+phi2 = 0;
+eta = 1;
+a = 0.95;
+b = 2.56;
+rho = 0.95;
+Psi = 1.47630583;
+thetass = 1.0125;
+
+model;
+
+(a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*a*exp(c)^(-b) = (a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*(1-a)*exp(m)^(-b)+beta*(a*exp(c(+1))^(1-b)+(1-a)*exp(m(+1))^(1-b))^((b-phi1)/(1-b))*a*exp(c(+1))^(-b)/(1+pi(+1));
+
+Psi*(1-exp(n))^(-eta)/(a*exp(c)^(-b)*(a*exp(c)^(1-b) + (1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))) = (1-alpha)*exp(y)/exp(n);
+
+(a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*a*exp(c)^(-b) = beta*exp(R(+1))*(a*exp(c(+1))^(1-b)+(1-a)*exp(m(+1))^(1-b))^((b-phi1)/(1-b))*a*exp(c(+1))^(-b);
+
+exp(R) = alpha*exp(y)/exp(k(-1)) + 1-delta;
+
+exp(k) = (1-delta)*exp(k(-1))+exp(y)-exp(c);
+
+exp(y) = exp(z)*exp(k(-1))^alpha*exp(n)^(1-alpha);
+
+exp(m) = exp(m(-1))*(u+thetass)/(1+pi);
+
+z = rho*z(-1) + e;
+
+u = gamm*u(-1) + phi2*z(-1) + sigma;
+
+end;
+
+shocks;
+var e; stderr 0.007;
+var sigma;stderr 0.0089;
+end;
+
+steady_state_model;
+// solving in levels
+// calibrating n = 1/3 and recovering the value of Psi
+// adapting solution Walsh (2003) p. 84
+pi = thetass-1;
+en = 1/3;
+eR = 1/beta;
+y_k = (1/alpha)*(1/beta-1+delta);
+ek = en*y_k^(-1/(1-alpha));
+ec = ek*(y_k-delta);
+em = ec*(a/(1-a))^(-1/b)*((thetass-beta)/thetass)^(-1/b);
+ey = ek*y_k;
+Xss = a*ec^(1-b)*(1+(a/(1-a))^(-1/b)*((thetass-beta)/thetass)^((b-1)/b));
+Psi = (1-alpha)*(ey/en)*Xss^((b-phi1)/(1-b))*a*ec^(-b)*(1-en)^eta;
+n = log(en);
+k = log(ek);
+m = log(em);
+c = log(ec);
+y = log(ey);
+R = log(eR);
+end;
+
+steady;
+
diff --git a/tests/steady_state/walsh1_ssm_block.mod b/tests/steady_state/walsh1_ssm_block.mod
new file mode 100644
index 0000000000000000000000000000000000000000..05bd18125390f4b9829b493ef9300a63d918ae74
--- /dev/null
+++ b/tests/steady_state/walsh1_ssm_block.mod
@@ -0,0 +1,72 @@
+var y c k m n R pi z u;
+varexo	e sigma;	 
+// sigma stands for phi in the eq 2.37 p.69
+
+parameters alpha beta delta gamm phi1 eta a b rho  phi2 Psi thetass;  
+//phi1 stands for capital phi in eq.2.68 and 2.69
+//phi2 stands for lowercase phi in eq. 2.66
+
+alpha = 0.36;
+beta = 0.989; 
+gamm = 0.5;
+delta = 0.019;
+phi1 = 2;
+phi2 = 0;
+eta = 1;
+a = 0.95;
+b = 2.56;
+rho = 0.95;
+Psi = 1.47630583;
+thetass = 1.0125;
+
+model(block);
+
+(a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*a*exp(c)^(-b) = (a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*(1-a)*exp(m)^(-b)+beta*(a*exp(c(+1))^(1-b)+(1-a)*exp(m(+1))^(1-b))^((b-phi1)/(1-b))*a*exp(c(+1))^(-b)/(1+pi(+1));
+
+Psi*(1-exp(n))^(-eta)/(a*exp(c)^(-b)*(a*exp(c)^(1-b) + (1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))) = (1-alpha)*exp(y)/exp(n);
+
+(a*exp(c)^(1-b)+(1-a)*exp(m)^(1-b))^((b-phi1)/(1-b))*a*exp(c)^(-b) = beta*exp(R(+1))*(a*exp(c(+1))^(1-b)+(1-a)*exp(m(+1))^(1-b))^((b-phi1)/(1-b))*a*exp(c(+1))^(-b);
+
+exp(R) = alpha*exp(y)/exp(k(-1)) + 1-delta;
+
+exp(k) = (1-delta)*exp(k(-1))+exp(y)-exp(c);
+
+exp(y) = exp(z)*exp(k(-1))^alpha*exp(n)^(1-alpha);
+
+exp(m) = exp(m(-1))*(u+thetass)/(1+pi);
+
+z = rho*z(-1) + e;
+
+u = gamm*u(-1) + phi2*z(-1) + sigma;
+
+end;
+
+shocks;
+var e; stderr 0.007;
+var sigma;stderr 0.0089;
+end;
+
+steady_state_model;
+// solving in levels
+// calibrating n = 1/3 and recovering the value of Psi
+// adapting solution Walsh (2003) p. 84
+pi = thetass-1;
+en = 1/3;
+eR = 1/beta;
+y_k = (1/alpha)*(1/beta-1+delta);
+ek = en*y_k^(-1/(1-alpha));
+ec = ek*(y_k-delta);
+em = ec*(a/(1-a))^(-1/b)*((thetass-beta)/thetass)^(-1/b);
+ey = ek*y_k;
+Xss = a*ec^(1-b)*(1+(a/(1-a))^(-1/b)*((thetass-beta)/thetass)^((b-1)/b));
+Psi = (1-alpha)*(ey/en)*Xss^((b-phi1)/(1-b))*a*ec^(-b)*(1-en)^eta;
+n = log(en);
+k = log(ek);
+m = log(em);
+c = log(ec);
+y = log(ey);
+R = log(eR);
+end;
+
+steady;
+