diff --git a/tests/.gitignore b/tests/.gitignore
index 7b5c438826aa63ae6af5d7cbf06ab5aa013605b8..3d7ccec14ddbb1e17ab97eb5ed9618ffab9cf8c6 100644
--- a/tests/.gitignore
+++ b/tests/.gitignore
@@ -96,9 +96,6 @@ wsOct
 !/kronecker/test_kron.m
 !/load_octave_packages.m
 !/ls2003/data_ca1.m
-!/matrix_notation/extFunNoDerivsMatrix.m
-!/matrix_notation/extFunFirstDerivMatrix.m
-!/matrix_notation/extFunWithFirstAndSecondDerivsMatrix.m
 !/measurement_errors/data_ca1.m
 !/measurement_errors/fs2000_corr_me_ml_mcmc/fsdat_simul.m
 !/missing/simulate_data_with_missing_observations.m
diff --git a/tests/Makefile.am b/tests/Makefile.am
index ee55632067859ee6627e93f755f5cf67b2b06b94..ef844e1bdc90b554daf1793c1585feb31669af1b 100644
--- a/tests/Makefile.am
+++ b/tests/Makefile.am
@@ -416,18 +416,7 @@ MODFILES = \
 	bgp/nk-1/nk.mod \
 	bgp/ramsey-1/ramsey.mod \
 	dynare-command-options/ramst.mod \
-	particle/local_state_space_iteration_k_test.mod \
-	matrix_notation/deterministic_matrix.mod \
-	matrix_notation/deterministic_scalar.mod \
-	matrix_notation/stochastic_matrix.mod \
-	matrix_notation/stochastic_scalar.mod \
-	matrix_notation/benchmark.mod \
-	matrix_notation/no_deriv_given_matrix_dll.mod \
-	matrix_notation/no_deriv_given_matrix.mod \
-	matrix_notation/first_deriv_given_matrix_dll.mod \
-	matrix_notation/first_deriv_given_matrix.mod \
-	matrix_notation/first_and_2nd_deriv_given_matrix_dll.mod \
-	matrix_notation/first_and_2nd_deriv_given_matrix.mod
+	particle/local_state_space_iteration_k_test.mod
 
 ECB_MODFILES = \
 	var-expectations/1/example.mod \
@@ -817,23 +806,6 @@ discretionary_policy/dennis_1_estim.o.trs: discretionary_policy/dennis_1.o.trs
 discretionary_policy/Gali_2015_chapter_3_nonlinear.m.trs: discretionary_policy/Gali_2015_chapter_3.m.trs
 discretionary_policy/Gali_2015_chapter_3_nonlinear.o.trs: discretionary_policy/Gali_2015_chapter_3.o.trs
 
-matrix_notation/deterministic_scalar.m.trs: matrix_notation/deterministic_matrix.m.trs
-matrix_notation/deterministic_scalar.o.trs: matrix_notation/deterministic_matrix.o.trs
-matrix_notation/stochastic_scalar.m.trs: matrix_notation/stochastic_matrix.m.trs
-matrix_notation/stochastic_scalar.o.trs: matrix_notation/stochastic_matrix.o.trs
-matrix_notation/first_deriv_given_matrix.m.trs: matrix_notation/benchmark.m.trs
-matrix_notation/first_deriv_given_matrix.o.trs: matrix_notation/benchmark.o.trs
-matrix_notation/first_deriv_given_matrix_dll.m.trs: matrix_notation/benchmark.m.trs
-matrix_notation/first_deriv_given_matrix_dll.o.trs: matrix_notation/benchmark.o.trs
-matrix_notation/first_and_2nd_deriv_given_matrix.m.trs: matrix_notation/benchmark.m.trs
-matrix_notation/first_and_2nd_deriv_given_matrix.o.trs: matrix_notation/benchmark.o.trs
-matrix_notation/first_and_2nd_deriv_given_matrix_dll.m.trs: matrix_notation/benchmark.m.trs
-matrix_notation/first_and_2nd_deriv_given_matrix_dll.o.trs: matrix_notation/benchmark.o.trs
-matrix_notation/no_deriv_given_matrix.m.trs: matrix_notation/benchmark.m.trs
-matrix_notation/no_deriv_given_matrix.o.trs: matrix_notation/benchmark.o.trs
-matrix_notation/no_deriv_given_matrix_dll.m.trs: matrix_notation/benchmark.m.trs
-matrix_notation/no_deriv_given_matrix_dll.o.trs: matrix_notation/benchmark.o.trs
-
 observation_trends_and_prefiltering/MCMC: m/observation_trends_and_prefiltering/MCMC o/observation_trends_and_prefiltering/MCMC
 m/observation_trends_and_prefiltering/MCMC: $(patsubst %.mod, %.m.trs, $(filter observation_trends_and_prefiltering/MCMC/%.mod, $(MODFILES)))
 o/observation_trends_and_prefiltering/MCMC: $(patsubst %.mod, %.o.trs, $(filter observation_trends_and_prefiltering/MCMC/%.mod, $(MODFILES)))
@@ -1046,10 +1018,6 @@ estimation/univariate: m/estimation/univariate o/estimation/univariate
 m/estimation/univariate: $(patsubst %.mod, %.m.trs, $(filter estimation/univariate/%.mod, $(MODFILES)))
 o/estimation/univariate: $(patsubst %.mod, %.o.trs, $(filter estimation/univariate/%.mod, $(MODFILES)))
 
-matrix_notation: m/matrix_notation o/matrix_notation
-m/matrix_notation: $(patsubst %.mod, %.m.trs, $(filter matrix_notation/%.mod, $(MODFILES)))
-o/matrix_notation: $(patsubst %.mod, %.o.trs, $(filter matrix_notation/%.mod, $(MODFILES)))
-
 # ECB files
 M_ECB_TRS_FILES = $(patsubst %.mod, %.m.trs, $(ECB_MODFILES))
 
@@ -1235,10 +1203,7 @@ EXTRA_DIST = \
 	histval_initval_file/my_assert.m \
 	histval_initval_file/ramst_data.xls \
 	histval_initval_file/ramst_data.xlsx \
-	histval_initval_file/ramst_initval_file_data.m \
-	matrix_notation/extFunFirstDerivMatrix.m \
-	matrix_notation/extFunNoDerivsMatrix.m \
-	matrix_notation/extFunWithFirstAndSecondDerivsMatrix.m
+	histval_initval_file/ramst_initval_file_data.m
 
 if ENABLE_MATLAB
 check-local: check-matlab
diff --git a/tests/matrix_notation/benchmark.mod b/tests/matrix_notation/benchmark.mod
deleted file mode 100644
index c33ae1358ffa1d51f63fbbee12077a2fab7decc8..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/benchmark.mod
+++ /dev/null
@@ -1,44 +0,0 @@
-/* This file is used as a benchmark against with the other tests with external
-   functions are compared.
-   It is almost the same as example1.mod, except that the shocks process is
-   nonlinear (in order to have a non-zero Hessian of the
-   external function) */
-
-var y, c, k, h, a, b;
-varexo e, u;
-parameters beta, alpha, delta, theta, psi, rho, 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)=(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)^2+5*tau*b(-1)^2 + e;
-b = 3*tau/2*a(-1)^2+7*rho/3*b(-1)^2 + u;
-end;
-
-initval;
-y = 1.08068253095672;
-c = 0.80359242014163;
-h = 0.29175631001732;
-k = 11.08360443260358;
-end;
-
-shocks;
-var e; stderr 0.009;
-var u; stderr 0.009;
-var e, u = phi*0.009*0.009;
-end;
-
-stoch_simul;
diff --git a/tests/matrix_notation/deterministic_matrix.mod b/tests/matrix_notation/deterministic_matrix.mod
deleted file mode 100644
index ff72d56c030f3ca8edff99e9bd1f88ac129d9658..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/deterministic_matrix.mod
+++ /dev/null
@@ -1,38 +0,0 @@
-// Do not forget to update deterministic_scalar.mod when this file is modified
-
-var(rows=3) X, Y, Z;
-var t;
-varexo(rows=3) E, U, V;
-parameters(rows=3, cols=3) P;
-
-P = [.1, .2, -.3; -.1, .5, .1; .2, -.3, .9];
-
-model;
-  X = P*X(-1)+E;
-  transpose(Y) = transpose(Y(-1))*transpose(P)+transpose(U);
-  Z = V;
-  t = Y[1](-1) + Y[2](+1);
-end;
-
-shocks;
-  var E;
-  periods 1 2;
-  values [0.5; 0.6; 0.7] [0.2; 0.3; 0.4];
-
-  var U[1];
-  periods 3;
-  values 0.8;
-
-  var U[2];
-  periods 3;
-  values 0.7;
-
-  var V[:];
-  periods 1;
-  values 0.1;
-end;
-
-perfect_foresight_setup(periods=50);
-perfect_foresight_solver;
-
-write_latex_dynamic_model;
diff --git a/tests/matrix_notation/deterministic_scalar.mod b/tests/matrix_notation/deterministic_scalar.mod
deleted file mode 100644
index 3a893acc49c066e99fe029d63705d4481cf0b1f3..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/deterministic_scalar.mod
+++ /dev/null
@@ -1,70 +0,0 @@
-// Scalar translation of deterministic_matrix.mod
-// The declaration orders are the same as those of the expanded matrix version
-var t X_1 X_2 X_3 Y_1 Y_2 Y_3 Z_1 Z_2 Z_3;
-varexo E_1 E_2 E_3 U_1 U_2 U_3 V_1 V_2 V_3;
-parameters P_1_1 P_1_2 P_1_3 P_2_1 P_2_2 P_2_3 P_3_1 P_3_2 P_3_3;
-
-P_1_1 = .1;
-P_1_2 = .2;
-P_1_3 = -.3;
-P_2_1 = -.1;
-P_2_2 = .5;
-P_2_3 = .1;
-P_3_1 = .2;
-P_3_2 = -.3;
-P_3_3 = .9;
-
-model;
-  X_1 = P_1_1*X_1(-1) + P_1_2*X_2(-1) + P_1_3*X_3(-1) + E_1;
-  X_2 = P_2_1*X_1(-1) + P_2_2*X_2(-1) + P_2_3*X_3(-1) + E_2;
-  X_3 = P_3_1*X_1(-1) + P_3_2*X_2(-1) + P_3_3*X_3(-1) + E_3;
-  Y_1 = P_1_1*Y_1(-1) + P_1_2*Y_2(-1) + P_1_3*Y_3(-1) + U_1;
-  Y_2 = P_2_1*Y_1(-1) + P_2_2*Y_2(-1) + P_2_3*Y_3(-1) + U_2;
-  Y_3 = P_3_1*Y_1(-1) + P_3_2*Y_2(-1) + P_3_3*Y_3(-1) + U_3;
-  Z_1 = V_1;
-  Z_2 = V_2;
-  Z_3 = V_3;
-  t = Y_1(-1) + Y_2(+1);
-end;
-
-shocks;
-  var E_1;
-  periods 1 2;
-  values 0.5 0.2;
-
-  var E_2;
-  periods 1 2;
-  values 0.6 0.3;
-
-  var E_3;
-  periods 1 2;
-  values 0.7 0.4;
-
-  var U_1;
-  periods 3;
-  values 0.8;
-
-  var U_2;
-  periods 3;
-  values 0.7;
-
-  var V_1;
-  periods 1;
-  values 0.1;
-
-  var V_2;
-  periods 1;
-  values 0.1;
-
-  var V_3;
-  periods 1;
-  values 0.1;
-end;
-
-perfect_foresight_setup(periods=50);
-perfect_foresight_solver;
-
-L = load('deterministic_matrix_results.mat');
-if max(max(abs(L.oo_.endo_simul - oo_.endo_simul))) > 1e-12
-  error('Failure in matrix expansion')
-end
diff --git a/tests/matrix_notation/extFunFirstDerivMatrix.m b/tests/matrix_notation/extFunFirstDerivMatrix.m
deleted file mode 100644
index 68c27d3a883cc179615b364d8bcee17405ba7b04..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/extFunFirstDerivMatrix.m
+++ /dev/null
@@ -1,4 +0,0 @@
-function dy=extFunFirstDerivMatrix(a,b,c)
-dy = [ c(1)^2 0        5*c(2)^2 0        2*a(1)*c(1) 10*b(1)*c(2);
-       0      3*c(1)^2 0        7*c(2)^2 6*a(2)*c(1) 14*b(2)*c(2) ];
-end
diff --git a/tests/matrix_notation/extFunNoDerivsMatrix.m b/tests/matrix_notation/extFunNoDerivsMatrix.m
deleted file mode 100644
index e7541cee70badbb3424555eee5c1a45608049e1b..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/extFunNoDerivsMatrix.m
+++ /dev/null
@@ -1,3 +0,0 @@
-function y=extFunNoDerivsMatrix(a,b,c)
-y=[a(1), 5*b(1); 3*a(2) 7*b(2)]*c.^2;
-end
diff --git a/tests/matrix_notation/extFunWithFirstAndSecondDerivsMatrix.m b/tests/matrix_notation/extFunWithFirstAndSecondDerivsMatrix.m
deleted file mode 100644
index e7b227909c2ef3dfa7461db60fff6e0e3295c6f8..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/extFunWithFirstAndSecondDerivsMatrix.m
+++ /dev/null
@@ -1,20 +0,0 @@
-function [y dy d2y]=extFunWithFirstAndSecondDerivsMatrix(a,b,c)
-y=[a(1), 5*b(1); 3*a(2), 7*b(2)]*c.^2;
-
-dy = [ c(1)^2 0        5*c(2)^2 0        2*a(1)*c(1) 10*b(1)*c(2);
-       0      3*c(1)^2 0        7*c(2)^2 6*a(2)*c(1) 14*b(2)*c(2) ];
-
-d2y=zeros(2,6,6);
-d2y(1,1,5) = 2*c(1);
-d2y(1,3,6) = 10*c(2);
-d2y(1,5,1) = 2*c(1);
-d2y(1,5,5) = 2*a(1);
-d2y(1,6,3) = 10*c(2);
-d2y(1,6,6) = 10*b(1);
-d2y(2,2,5) = 6*c(1);
-d2y(2,4,6) = 14*c(2);
-d2y(2,5,2) = 6*c(1);
-d2y(2,5,5) = 6*a(2);
-d2y(2,6,4) = 14*c(2);
-d2y(2,6,6) = 14*b(2);
-end
diff --git a/tests/matrix_notation/first_and_2nd_deriv_given_matrix.mod b/tests/matrix_notation/first_and_2nd_deriv_given_matrix.mod
deleted file mode 100644
index 120ef1043a2c0ec03f482f72ab022816acb092bb..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/first_and_2nd_deriv_given_matrix.mod
+++ /dev/null
@@ -1,71 +0,0 @@
-var y, c, k, h;
-var(rows = 2) A; // A == [a; b]
-varexo(rows = 2) E; // E == [e; u]
-parameters beta, alpha, delta, theta, psi;
-parameters(rows=2,cols=2) P;
-
-alpha = 0.36;
-rho   = 0.95;
-tau   = 0.025;
-beta  = 0.99;
-delta = 0.025;
-psi   = 0;
-theta = 2.95;
-
-phi   = 0.1;
-
-P = [rho, tau; tau/2, rho/3];
-
-external_function(nargs=3, name=extFunWithFirstAndSecondDerivsMatrix, in_arg_dim=[2,2,2],
-                  out_arg_dim=2, first_deriv_provided, second_deriv_provided);
-
-model;
-[endogenous='c',name='law of motion of capital']
-c*theta*h^(1+psi)=(1-alpha)*y;
-k = beta*(((exp(A[2])*c)/(exp(A[2](+1))*c(+1)))
-          *(exp(A[2](+1))*alpha*y(+1)+(1-delta)*k));
-y = exp(A[1])*(k(-1)^alpha)*(h^(1-alpha));
-k = exp(A[2])*(y-c)+(1-delta)*k(-1) ;
-A = extFunWithFirstAndSecondDerivsMatrix(P[:,1], transpose(P[:,2]), A(-1)) + E;
-end;
-
-initval;
-y = 1.08068253095672;
-c = 0.80359242014163;
-h = 0.29175631001732;
-k = 11.08360443260358;
-
-A[1] = 0;
-A[2] = 0;
-E[:] = 0;
-end;
-
-shocks;
-var E[:]; stderr 0.009;
-var E[1], E[2] = phi*0.009*0.009;
-end;
-
-stoch_simul;
-
-L = load('benchmark_results.mat');
-if max(max(abs(L.oo_.dr.ghu - oo_.dr.ghu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghx - oo_.dr.ghx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxu - oo_.dr.ghxu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxx - oo_.dr.ghxx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghuu - oo_.dr.ghuu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghs2 - oo_.dr.ghs2))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.var - oo_.var))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
diff --git a/tests/matrix_notation/first_and_2nd_deriv_given_matrix_dll.mod b/tests/matrix_notation/first_and_2nd_deriv_given_matrix_dll.mod
deleted file mode 100644
index 570cf021dedddde43c01f11f0d8e5e3761a79f98..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/first_and_2nd_deriv_given_matrix_dll.mod
+++ /dev/null
@@ -1,71 +0,0 @@
-var y, c, k, h;
-var(rows = 2) A; // A == [a; b]
-varexo(rows = 2) E; // E == [e; u]
-parameters beta, alpha, delta, theta, psi;
-parameters(rows=2,cols=2) P;
-
-alpha = 0.36;
-rho   = 0.95;
-tau   = 0.025;
-beta  = 0.99;
-delta = 0.025;
-psi   = 0;
-theta = 2.95;
-
-phi   = 0.1;
-
-P = [rho, tau; tau/2, rho/3];
-
-external_function(nargs=3, name=extFunWithFirstAndSecondDerivsMatrix, in_arg_dim=[2,2,2],
-                  out_arg_dim=2, first_deriv_provided, second_deriv_provided);
-
-model(use_dll);
-[endogenous='c',name='law of motion of capital']
-c*theta*h^(1+psi)=(1-alpha)*y;
-k = beta*(((exp(A[2])*c)/(exp(A[2](+1))*c(+1)))
-          *(exp(A[2](+1))*alpha*y(+1)+(1-delta)*k));
-y = exp(A[1])*(k(-1)^alpha)*(h^(1-alpha));
-k = exp(A[2])*(y-c)+(1-delta)*k(-1) ;
-A = extFunWithFirstAndSecondDerivsMatrix(P[:,1], transpose(P[:,2]), A(-1)) + E;
-end;
-
-initval;
-y = 1.08068253095672;
-c = 0.80359242014163;
-h = 0.29175631001732;
-k = 11.08360443260358;
-
-A[1] = 0;
-A[2] = 0;
-E[:] = 0;
-end;
-
-shocks;
-var E[:]; stderr 0.009;
-var E[1], E[2] = phi*0.009*0.009;
-end;
-
-stoch_simul;
-
-L = load('benchmark_results.mat');
-if max(max(abs(L.oo_.dr.ghu - oo_.dr.ghu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghx - oo_.dr.ghx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxu - oo_.dr.ghxu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxx - oo_.dr.ghxx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghuu - oo_.dr.ghuu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghs2 - oo_.dr.ghs2))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.var - oo_.var))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
diff --git a/tests/matrix_notation/first_deriv_given_matrix.mod b/tests/matrix_notation/first_deriv_given_matrix.mod
deleted file mode 100644
index 431df3d64395ef2c97e1a4a6209a8c401eb7c136..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/first_deriv_given_matrix.mod
+++ /dev/null
@@ -1,72 +0,0 @@
-// --+ options: json=compute +--
-
-var y, c, k, h;
-var(rows = 2) A; // A == [a; b]
-varexo(rows = 2) E; // E == [e; u]
-parameters beta, alpha, delta, theta, psi;
-parameters(rows=2,cols=2) P;
-
-alpha = 0.36;
-rho   = 0.95;
-tau   = 0.025;
-beta  = 0.99;
-delta = 0.025;
-psi   = 0;
-theta = 2.95;
-
-phi   = 0.1;
-
-P = [rho, tau; tau/2, rho/3];
-
-external_function(nargs=3, name=extFunNoDerivsMatrix, in_arg_dim=[2,2,2],
-                  out_arg_dim=2, first_deriv_provided=extFunFirstDerivMatrix);
-
-model;
-c*theta*h^(1+psi)=(1-alpha)*y;
-k = beta*(((exp(A[2])*c)/(exp(A[2](+1))*c(+1)))
-          *(exp(A[2](+1))*alpha*y(+1)+(1-delta)*k));
-y = exp(A[1])*(k(-1)^alpha)*(h^(1-alpha));
-k = exp(A[2])*(y-c)+(1-delta)*k(-1) ;
-A = extFunNoDerivsMatrix(P[:,1], transpose(P[:,2]), A(-1)) + E;
-end;
-
-initval;
-y = 1.08068253095672;
-c = 0.80359242014163;
-h = 0.29175631001732;
-k = 11.08360443260358;
-
-A[1] = 0;
-A[2] = 0;
-E[:] = 0;
-end;
-
-shocks;
-var E[:]; stderr 0.009;
-var E[1], E[2] = phi*0.009*0.009;
-end;
-
-stoch_simul;
-
-L = load('benchmark_results.mat');
-if max(max(abs(L.oo_.dr.ghu - oo_.dr.ghu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghx - oo_.dr.ghx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxu - oo_.dr.ghxu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxx - oo_.dr.ghxx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghuu - oo_.dr.ghuu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghs2 - oo_.dr.ghs2))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.var - oo_.var))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
diff --git a/tests/matrix_notation/first_deriv_given_matrix_dll.mod b/tests/matrix_notation/first_deriv_given_matrix_dll.mod
deleted file mode 100644
index 070195c210c4d1c4db8a3b3a6de874f4189050e9..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/first_deriv_given_matrix_dll.mod
+++ /dev/null
@@ -1,73 +0,0 @@
-// --+ options: json=compute +--
-
-var y, c, k, h;
-var(rows = 2) A; // A == [a; b]
-varexo(rows = 2) E; // E == [e; u]
-parameters beta, alpha, delta, theta, psi;
-parameters(rows=2,cols=2) P;
-
-alpha = 0.36;
-rho   = 0.95;
-tau   = 0.025;
-beta  = 0.99;
-delta = 0.025;
-psi   = 0;
-theta = 2.95;
-
-phi   = 0.1;
-
-P = [rho, tau; tau/2, rho/3];
-
-external_function(nargs=3, name=extFunNoDerivsMatrix, in_arg_dim=[2,2,2],
-                  out_arg_dim=2, first_deriv_provided=extFunFirstDerivMatrix);
-
-model(use_dll);
-c*theta*h^(1+psi)=(1-alpha)*y;
-k = beta*(((exp(A[2])*c)/(exp(A[2](+1))*c(+1)))
-          *(exp(A[2](+1))*alpha*y(+1)+(1-delta)*k));
-y = exp(A[1])*(k(-1)^alpha)*(h^(1-alpha));
-k = exp(A[2])*(y-c)+(1-delta)*k(-1) ;
-A = extFunNoDerivsMatrix(P[:,1], transpose(P[:,2]), A(-1)) + E;
-end;
-
-initval;
-y = 1.08068253095672;
-c = 0.80359242014163;
-h = 0.29175631001732;
-k = 11.08360443260358;
-
-A[1] = 0;
-A[2] = 0;
-E[:] = 0;
-
-end;
-
-shocks;
-var E[:]; stderr 0.009;
-var E[1], E[2] = phi*0.009*0.009;
-end;
-
-stoch_simul;
-
-L = load('benchmark_results.mat');
-if max(max(abs(L.oo_.dr.ghu - oo_.dr.ghu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghx - oo_.dr.ghx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxu - oo_.dr.ghxu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxx - oo_.dr.ghxx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghuu - oo_.dr.ghuu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghs2 - oo_.dr.ghs2))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.var - oo_.var))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
diff --git a/tests/matrix_notation/no_deriv_given_matrix.mod b/tests/matrix_notation/no_deriv_given_matrix.mod
deleted file mode 100644
index 08c22198ffa5c313da9c1f24db60884ad4ee2ed0..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/no_deriv_given_matrix.mod
+++ /dev/null
@@ -1,71 +0,0 @@
-// --+ options: json=compute +--
-
-var y, c, k, h;
-var(rows = 2) A; // A == [a; b]
-varexo(rows = 2) E; // E == [e; u]
-parameters beta, alpha, delta, theta, psi;
-parameters(rows=2,cols=2) P;
-
-alpha = 0.36;
-rho   = 0.95;
-tau   = 0.025;
-beta  = 0.99;
-delta = 0.025;
-psi   = 0;
-theta = 2.95;
-
-phi   = 0.1;
-
-P = [rho, tau; tau/2, rho/3];
-
-external_function(nargs=3, name=extFunNoDerivsMatrix, in_arg_dim=[2,2,2], out_arg_dim=2);
-
-model;
-c*theta*h^(1+psi)=(1-alpha)*y;
-k = beta*(((exp(A[2])*c)/(exp(A[2](+1))*c(+1)))
-          *(exp(A[2](+1))*alpha*y(+1)+(1-delta)*k));
-y = exp(A[1])*(k(-1)^alpha)*(h^(1-alpha));
-k = exp(A[2])*(y-c)+(1-delta)*k(-1) ;
-A = extFunNoDerivsMatrix(P[:,1], transpose(P[:,2]), A(-1)) + E;
-end;
-
-initval;
-y = 1.08068253095672;
-c = 0.80359242014163;
-h = 0.29175631001732;
-k = 11.08360443260358;
-
-A[1] = 0;
-A[2] = 0;
-E[:] = 0;
-end;
-
-shocks;
-var E[:]; stderr 0.009;
-var E[1], E[2] = phi*0.009*0.009;
-end;
-
-stoch_simul;
-
-L = load('benchmark_results.mat');
-if max(max(abs(L.oo_.dr.ghu - oo_.dr.ghu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghx - oo_.dr.ghx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxu - oo_.dr.ghxu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxx - oo_.dr.ghxx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghuu - oo_.dr.ghuu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghs2 - oo_.dr.ghs2))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.var - oo_.var))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
diff --git a/tests/matrix_notation/no_deriv_given_matrix_dll.mod b/tests/matrix_notation/no_deriv_given_matrix_dll.mod
deleted file mode 100644
index ccbedff1c778bf66748e7e9c5ec1074ea190fab3..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/no_deriv_given_matrix_dll.mod
+++ /dev/null
@@ -1,71 +0,0 @@
-// --+ options: json=compute +--
-
-var y, c, k, h;
-var(rows = 2) A; // A == [a; b]
-varexo(rows = 2) E; // E == [e; u]
-parameters beta, alpha, delta, theta, psi;
-parameters(rows=2,cols=2) P;
-
-alpha = 0.36;
-rho   = 0.95;
-tau   = 0.025;
-beta  = 0.99;
-delta = 0.025;
-psi   = 0;
-theta = 2.95;
-
-phi   = 0.1;
-
-P = [rho, tau; tau/2, rho/3];
-
-external_function(nargs=3, name=extFunNoDerivsMatrix, in_arg_dim=[2,2,2], out_arg_dim=2);
-
-model(use_dll);
-c*theta*h^(1+psi)=(1-alpha)*y;
-k = beta*(((exp(A[2])*c)/(exp(A[2](+1))*c(+1)))
-          *(exp(A[2](+1))*alpha*y(+1)+(1-delta)*k));
-y = exp(A[1])*(k(-1)^alpha)*(h^(1-alpha));
-k = exp(A[2])*(y-c)+(1-delta)*k(-1) ;
-A = extFunNoDerivsMatrix(P[:,1], transpose(P[:,2]), A(-1)) + E;
-end;
-
-initval;
-y = 1.08068253095672;
-c = 0.80359242014163;
-h = 0.29175631001732;
-k = 11.08360443260358;
-
-A[1] = 0;
-A[2] = 0;
-E[:] = 0;
-end;
-
-shocks;
-var E[:]; stderr 0.009;
-var E[1], E[2] = phi*0.009*0.009;
-end;
-
-stoch_simul;
-
-L = load('benchmark_results.mat');
-if max(max(abs(L.oo_.dr.ghu - oo_.dr.ghu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghx - oo_.dr.ghx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxu - oo_.dr.ghxu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghxx - oo_.dr.ghxx))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghuu - oo_.dr.ghuu))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.dr.ghs2 - oo_.dr.ghs2))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
-if max(max(abs(L.oo_.var - oo_.var))) > 1e-12
-  error('Failure in matrix notation with external function')
-end
diff --git a/tests/matrix_notation/stochastic_matrix.mod b/tests/matrix_notation/stochastic_matrix.mod
deleted file mode 100644
index d055f5e60c4c5ee028b968792e63c76889589c61..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/stochastic_matrix.mod
+++ /dev/null
@@ -1,46 +0,0 @@
-// Do not forget to update stochastic_scalar.mod when this file is modified
-
-var(rows=3) X, Y, Z, T;
-var r;
-varexo(rows=3) E, U, V, W;
-varexo f;
-parameters(rows=3, cols=3) P;
-varexo_det(rows=3) D;
-
-P = [.1, .2, -.3; -.1, .5, .1; .2, -.3, .9];
-
-model;
-  X = P*X(-1)+E;
-  transpose(Y) = transpose(Y(-1))*transpose(P)+transpose(U);
-  Z[1] = V[1];
-  Z[2] = V[2];
-  Z[3] = V[3];
-  T = D + W;
-  r = 0.9*r(+1)+f;
-end;
-
-shocks;
-  var E = [ 0.1, 0.2, 0.3 ];
-  var E[1], E[2] = 0.08;
-
-  var U[1] = 0.4;
-  var U[2]; stderr 0.5;
-  var U[3]; stderr 0.6;
-  corr U[1], U[2] = 0.09;
-
-  var V; stderr [ 1.1, 1.2, 1.3 ];
-
-  vcov W = [
-      1, 0.2, 0.3;
-    0.2, 5,   0.6;
-    0.3, 0.6,   9
-  ];
-
-  var f = 3;
-  var f, E[2] = 0.5;
-  corr U[1], f = 0.7;
-end;
-
-stoch_simul(order=1, nodecomposition, irf=0);
-
-write_latex_dynamic_model;
diff --git a/tests/matrix_notation/stochastic_scalar.mod b/tests/matrix_notation/stochastic_scalar.mod
deleted file mode 100644
index 3c5be1857c21938e5398bffb8f341eba5241d8bf..0000000000000000000000000000000000000000
--- a/tests/matrix_notation/stochastic_scalar.mod
+++ /dev/null
@@ -1,72 +0,0 @@
-// Scalar translation of stochastic_matrix.mod
-// The declaration orders are the same as those of the expanded matrix version
-var r X_1 X_2 X_3 Y_1 Y_2 Y_3 Z_1 Z_2 Z_3 T_1 T_2 T_3;
-varexo f E_1 E_2 E_3 U_1 U_2 U_3 V_1 V_2 V_3 W_1 W_2 W_3;
-parameters P_1_1 P_1_2 P_1_3 P_2_1 P_2_2 P_2_3 P_3_1 P_3_2 P_3_3;
-varexo_det D_1 D_2 D_3;
-
-P_1_1 = .1;
-P_1_2 = .2;
-P_1_3 = -.3;
-P_2_1 = -.1;
-P_2_2 = .5;
-P_2_3 = .1;
-P_3_1 = .2;
-P_3_2 = -.3;
-P_3_3 = .9;
-
-model;
-  X_1 = P_1_1*X_1(-1) + P_1_2*X_2(-1) + P_1_3*X_3(-1) + E_1;
-  X_2 = P_2_1*X_1(-1) + P_2_2*X_2(-1) + P_2_3*X_3(-1) + E_2;
-  X_3 = P_3_1*X_1(-1) + P_3_2*X_2(-1) + P_3_3*X_3(-1) + E_3;
-  Y_1 = P_1_1*Y_1(-1) + P_1_2*Y_2(-1) + P_1_3*Y_3(-1) + U_1;
-  Y_2 = P_2_1*Y_1(-1) + P_2_2*Y_2(-1) + P_2_3*Y_3(-1) + U_2;
-  Y_3 = P_3_1*Y_1(-1) + P_3_2*Y_2(-1) + P_3_3*Y_3(-1) + U_3;
-  Z_1 = V_1;
-  Z_2 = V_2;
-  Z_3 = V_3;
-  T_1 = D_1 + W_1;
-  T_2 = D_2 + W_2;
-  T_3 = D_3 + W_3;
-  r = 0.9*r(+1)+f;
-end;
-
-shocks;
-  var E_1 = 0.1;
-  var E_2 = 0.2;
-  var E_3 = 0.3;
-  var E_1, E_2 = 0.08;
-
-  var U_1 = 0.4;
-  var U_2; stderr 0.5;
-  var U_3; stderr 0.6;
-  corr U_1, U_2 = 0.09;
-
-  var V_1; stderr 1.1;
-  var V_2; stderr 1.2;
-  var V_3; stderr 1.3;
-
-  var W_1 = 1;
-  var W_2 = 5;
-  var W_3 = 9;
-  var W_1, W_2 = 0.2;
-  var W_1, W_3 = 0.3;
-  var W_2, W_3 = 0.6;
-
-  var f = 3;
-  var f, E_2 = 0.5;
-  corr U_1, f = 0.7;
-end;
-
-stoch_simul(order=1, nodecomposition, irf=0);
-
-L = load('stochastic_matrix_results.mat');
-if max(max(abs(L.oo_.dr.ghu - oo_.dr.ghu))) > 1e-12
-  error('Failure in matrix expansion')
-end
-if max(max(abs(L.oo_.dr.ghx - oo_.dr.ghx))) > 1e-12
-  error('Failure in matrix expansion')
-end
-if max(max(abs(L.oo_.var - oo_.var))) > 1e-12
-  error('Failure in matrix expansion')
-end