diff --git a/matlab/run_simulations.m b/matlab/run_simulations.m
new file mode 100644
index 0000000000000000000000000000000000000000..5bb027a168ae6a5ae2756ab13958910c70d42fec
--- /dev/null
+++ b/matlab/run_simulations.m
@@ -0,0 +1,272 @@
+function run_simulations(model, stack_solve_algo, static_solve_algo, niter, ccrit)
+
+% INPUTS
+% - model               [char]      1×n array, name of the mod file (without extension).
+% - stack_solve_algo    [integer]   1×p array, identifiers for the solver of the stack model.
+% - static_solve_algo   [integer]   1×q array, identifiers for the solver of the static model (only used with stack_solve_algo==7)
+% - niter               [integer]   scalar, number of simulation of each configuration.
+% - ccrit               [double]    scalar, precision threshold.
+%
+% OUTPUTS
+% None
+
+% Copyright © 2024 Dynare Team
+
+id = 1:length(stack_solve_algo);
+jd = 1:length(static_solve_algo);
+
+%
+% Initialise the arrys storing the timings
+%
+Timings = matfile(sprintf('simulations-%s.mat', model));
+if exist(sprintf('simulations-%s.mat', model), 'file')
+    Timings.Properties.Writable = true;
+end
+if ~ismember('MATLAB_computing', {whos(Timings).name})
+    Timings.MATLAB_preprocessing = NaN(length(stack_solve_algo), length(static_solve_algo));
+    Timings.MATLAB_computing = NaN(length(stack_solve_algo), length(static_solve_algo));
+end
+if ~ismember('USE_DLL_computing', {whos(Timings).name})
+    Timings.USE_DLL_preprocessing = NaN(length(stack_solve_algo), length(static_solve_algo));
+    Timings.USE_DLL_computing = NaN(length(stack_solve_algo), length(static_solve_algo));
+end
+if ~ismember('BYTECODE_computing', {whos(Timings).name})
+    Timings.BYTECODE_preprocessing = NaN(length(stack_solve_algo), 1);
+    Timings.BYTECODE_computing = NaN(length(stack_solve_algo), 1);
+end
+if ~ismember('MATLAB_WITH_BLOCKS_computing', {whos(Timings).name})
+    Timings.MATLAB_WITH_BLOCKS_preprocessing = NaN(length(stack_solve_algo), 4);
+    Timings.MATLAB_WITH_BLOCKS_computing = NaN(length(stack_solve_algo), 4);
+end
+if ~ismember('USE_DLL_WITH_BLOCKS_computing', {whos(Timings).name})
+    Timings.USE_DLL_WITH_BLOCKS_preprocessing = NaN(length(stack_solve_algo), 4);
+    Timings.USE_DLL_WITH_BLOCKS_computing = NaN(length(stack_solve_algo), 4);
+end
+if ~ismember('BYTECODE_WITH_BLOCKS_computing', {whos(Timings).name})
+    Timings.BYTECODE_WITH_BLOCKS_preprocessing = NaN(length(stack_solve_algo), 4);
+    Timings.BYTECODE_WITH_BLOCKS_computing = NaN(length(stack_solve_algo), 4);
+end
+%
+% Matlab
+%
+for algo = [0, 1]
+    if algo<7
+        if isnan(Timings.MATLAB_computing(id(stack_solve_algo==algo), 1))
+            % Try to simulate the model if this configuration has not been tried before (with or without success).
+            try
+                disp_title(sprintf('| MATLAB (solve_algo=%u, no blocks, model=%s) |', algo, model));
+                dprintf('dynare %s onlyclearglobals -DITERATIONS=%u -DSTACK_SOLVE_ALGO_VALUE=%u', model, niter, algo);
+                skipline()
+                Timings.MATLAB_computing(id(stack_solve_algo==algo), 1) = -3; % Provision for MATLAB crash
+                info = dynare(sprintf('%s', model), 'onlyclearglobals', sprintf('-DITERATIONS=%u', niter), sprintf('-DSTACK_SOLVE_ALGO_VALUE=%u', algo));
+                if oo_.deterministic_simulation.error<ccrit
+                    Timings.MATLAB_preprocessing(id(stack_solve_algo==algo), 1) = info.time.preprocessor;
+                    Timings.MATLAB_computing(id(stack_solve_algo==algo), 1) = info.time.compute/niter;
+                else
+                    Timings.MATLAB_preprocessing(id(stack_solve_algo==algo), 1) = info.time.preprocessor;
+                    Timings.MATLAB_computing(id(stack_solve_algo==algo), 1) = -2;
+                end
+            catch E
+                fprintf(1, 'Simulation failed:\n\n%s', E.message);
+                Timings.MATLAB_preprocessing(id(stack_solve_algo==algo), 1) = info.time.preprocessor;
+                Timings.MATLAB_computing(id(stack_solve_algo==algo), 1) = -1;
+            end
+        end
+    else
+        for solver = static_solve_algo
+            if isnan(Timings.MATLAB_computing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)))
+                % Try to simulate the model if this configuration has not been tried before (with or without success).
+                try
+                    disp_title(sprintf('| MATLAB (solve_algo=%u, no blocks, model=%s) |', algo, model));
+                    dprintf('dynare %s onlyclearglobals -DITERATIONS=%u -DSTACK_SOLVE_ALGO_VALUE=%u -DSTEADY_SOLVE_ALGO_VALUE=%u', model, niter, algo, solver);
+                    skipline()
+                    Timings.MATLAB_computing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = -3; % Provision for MATLAB crash
+                    info = dynare(sprintf('%s', model), 'onlyclearglobals', sprintf('-DITERATIONS=%u', niter), sprintf('-DSTACK_SOLVE_ALGO_VALUE=%u', algo), sprintf('-DSTEADY_SOLVE_ALGO_VALUE=%u', solver));
+                    if oo_.deterministic_simulation.error<ccrit
+                        Timings.MATLAB_preprocessing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = info.time.preprocessor;
+                        Timings.MATLAB_computing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = info.time.compute/niter;
+                    else
+                        Timings.MATLAB_preprocessing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = info.time.preprocessor;
+                        Timings.MATLAB_computing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = -2;
+                    end
+                catch E
+                    fprintf(1, 'Simulation failed:\n\n%s', E.message);
+                    Timings.MATLAB_preprocessing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = info.time.preprocessor;
+                    Timings.MATLAB_computing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = -1;
+                end
+            end
+        end
+    end
+end
+for algo = [2, 3, 4, 5]
+    % These algorithms are not available
+    Timings.MATLAB_preprocessing(id(stack_solve_algo==algo), 1) = -10;
+    Timings.MATLAB_computing(id(stack_solve_algo==algo), 1) = -10;
+end
+%
+% use_dll
+%
+for algo = [0, 1]
+    if algo<7
+        try
+            if isnan(Timings.USE_DLL_computing(id(stack_solve_algo==algo), 1))
+                disp_title(sprintf('| USE_DLL (solve_algo=%u, no blocks, model=%s) |', algo, model));
+                dprintf('dynare %s onlyclearglobals fast -DUSE_DLL=true -DITERATIONS=%u -DSTACK_SOLVE_ALGO_VALUE=%u', model, niter, algo);
+                skipline()
+                Timings.USE_DLL_computing(id(stack_solve_algo==algo), 1) = -3; % Provision for MATLAB crash
+                info = dynare(sprintf('%s', model), 'onlyclearglobals', 'fast', '-DUSE_DLL=true', sprintf('-DITERATIONS=%u', niter), sprintf('-DSTACK_SOLVE_ALGO_VALUE=%u', algo));
+                if oo_.deterministic_simulation.error<ccrit
+                    Timings.USE_DLL_preprocessing(id(stack_solve_algo==algo), 1) = info.time.preprocessor;
+                    Timings.USE_DLL_computing(id(stack_solve_algo==algo), 1) = info.time.compute/niter;
+                else
+                    Timings.USE_DLL_preprocessing(id(stack_solve_algo==algo), 1) = info.time.preprocessor;
+                    Timings.USE_DLL_computing(id(stack_solve_algo==algo), 1) = -2;
+                end
+            end
+        catch E
+            fprintf(1, 'Simulation failed:\n\n%s', E.message);
+            Timings.USE_DLL_preprocessing(id(stack_solve_algo==algo), 1) = info.time.preprocessor;
+            Timings.USE_DLL_computing(id(stack_solve_algo==algo), 1) = -1;
+        end
+    else
+        for solver = static_solve_algo
+            try
+                if isnan(Timings.USE_DLL_computing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)))
+                    disp_title(sprintf('| USE_DLL (solve_algo=%u, no blocks, model=%s) |', algo, model));
+                    dprintf('dynare %s onlyclearglobals fast -DUSE_DLL=true -DITERATIONS=%u -DSTACK_SOLVE_ALGO_VALUE=%u -DSTEADY_SOLVE_ALGO_VALUE=%u', model, niter, algo, solver);
+                    skipline()
+                    Timings.USE_DLL_computing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = -3; % Provision for MATLAB crash
+                    info = dynare(sprintf('%s', model), 'onlyclearglobals', 'fast', '-DUSE_DLL=true', sprintf('-DITERATIONS=%u', niter), sprintf('-DSTACK_SOLVE_ALGO_VALUE=%u', algo), sprintf('-DSTEADY_SOLVE_ALGO_VALUE=%u', solver));
+                    if oo_.deterministic_simulation.error<ccrit
+                        Timings.USE_DLL_preprocessing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = info.time.preprocessor;
+                        Timings.USE_DLL_computing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = info.time.compute/niter;
+                    else
+                        Timings.USE_DLL_preprocessing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = info.time.preprocessor;
+                        Timings.USE_DLL_computing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = -2;
+                    end
+                end
+            catch E
+                fprintf(1, 'Simulation failed:\n\n%s', E.message);
+                Timings.USE_DLL_preprocessing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = info.time.preprocessor;
+                Timings.USE_DLL_computing(id(stack_solve_algo==algo), jd(solver==static_solve_algo)) = -1;
+            end
+        end
+    end
+end
+for algo = [2, 3, 4, 5]
+    % These algorithms are not available
+    Timings.USE_DLL_preprocessing(id(stack_solve_algo==algo), 1) = -10;
+    Timings.USE_DLL_computing(id(stack_solve_algo==algo), 1) = -10;
+end
+%
+% Bytecode
+%
+for algo = [0, 1, 2, 3, 4, 5]
+    try
+        if isnan(Timings.BYTECODE_computing(id(stack_solve_algo==algo), 1))
+            disp_title(sprintf('| BYTECODE (solve_algo=%u, no blocks, model=%s) |', algo, model));
+            dprintf('dynare %s onlyclearglobals -DBYTECODE=true -DITERATIONS=%u -DSTACK_SOLVE_ALGO_VALUE=%u', model, niter, algo);
+            skipline()
+            Timings.BYTECODE_computing(id(stack_solve_algo==algo), 1) = -3; % Provision for MATLAB crash
+            info = dynare(sprintf('%s', model), 'onlyclearglobals', sprintf('-DITERATIONS=%u', niter), '-DBYTECODE=true', sprintf('-DSTACK_SOLVE_ALGO_VALUE=%u', algo));
+            if oo_.deterministic_simulation.error<ccrit
+                Timings.BYTECODE_preprocessing(id(stack_solve_algo==algo), 1) = info.time.preprocessor;
+                Timings.BYTECODE_computing(id(stack_solve_algo==algo), 1) = info.time.compute/niter;
+            else
+                Timings.BYTECODE_preprocessing(id(stack_solve_algo==algo), 1) = info.time.preprocessor;
+                Timings.BYTECODE_computing(id(stack_solve_algo==algo), 1) = -2;
+            end
+        end
+    catch E
+        fprintf(1, 'Simulation failed:\n\n%s', E.message);
+        Timings.BYTECODE_preprocessing(id(stack_solve_algo==algo), 1) = info.time.preprocessor;
+        Timings.BYTECODE_computing(id(stack_solve_algo==algo), 1) = -1;
+    end
+end
+%
+% Matlab with blocks
+%
+for mfs = 0:3
+    for algo = [0, 1, 2, 3, 4]
+        if isnan(Timings.MATLAB_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1))
+            try
+                disp_title(sprintf('| MATLAB (solve_algo=%u, with blocks, mfs=%u, model=%s) |', algo, mfs, model));
+                dprintf('dynare %s onlyclearglobals -DITERATIONS=%u -DSTACK_SOLVE_ALGO_VALUE=%u -DMFS_VALUE=%u -DBLOCKS=true', model, niter, algo, mfs);
+                skipline()
+                Timings.MATLAB_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = -3; % Provision for MATLAB crash
+                info = dynare(sprintf('%s', model), 'onlyclearglobals', sprintf('-DITERATIONS=%u', niter), sprintf('-DMFS_VALUE=%u', mfs), '-DBLOCKS=true', sprintf('-DSTACK_SOLVE_ALGO_VALUE=%u', algo));
+                if oo_.deterministic_simulation.error<ccrit
+                    Timings.MATLAB_WITH_BLOCKS_preprocessing(id(stack_solve_algo==algo), mfs+1) = info.time.preprocessor;
+                    Timings.MATLAB_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = info.time.compute/niter;
+                else
+                    Timings.MATLAB_WITH_BLOCKS_preprocessing(id(stack_solve_algo==algo), mfs+1) = info.time.preprocessor;
+                    Timings.MATLAB_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = -2;
+                end
+            catch E
+                fprintf(1, 'Simulation failed:\n\n%s', E.message);
+                Timings.MATLAB_WITH_BLOCKS_preprocessing(id(stack_solve_algo==algo), mfs+1) = info.time.preprocessor;
+                Timings.MATLAB_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = -1;
+            end
+        end
+    end
+end
+%
+% USE_DLL with blocks
+%
+for mfs = 0:3
+    for algo = [0, 1, 2, 3, 4]
+        if isnan(Timings.USE_DLL_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1))
+            try
+                disp_title(sprintf('| USE_DLL (solve_algo=%u, with blocks, mfs=%u, model=%s) |', algo, mfs, model));
+                dprintf('dynare %s onlyclearglobals fast -DUSE_DLL=true -DITERATIONS=%u -DSTACK_SOLVE_ALGO_VALUE=%u -DMFS_VALUE=%u -DBLOCKS=true', model, niter, algo, mfs);
+                skipline()
+                Timings.USE_DLL_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = -3; % Provision for MATLAB crash
+                info = dynare(sprintf('%s', model), 'onlyclearglobals', 'fast', sprintf('-DITERATIONS=%u', niter), '-DBLOCKS=true', sprintf('-DMFS_VALUE=%u', mfs), '-DUSE_DLL=true', sprintf('-DSTACK_SOLVE_ALGO_VALUE=%u', algo));
+                if oo_.deterministic_simulation.error<ccrit
+                    Timings.USE_DLL_WITH_BLOCKS_preprocessing(id(stack_solve_algo==algo), mfs+1) = info.time.preprocessor;
+                    Timings.USE_DLL_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = info.time.compute/niter;
+                else
+                    Timings.USE_DLL_WITH_BLOCKS_preprocessing(id(stack_solve_algo==algo), mfs+1) = info.time.preprocessor;
+                    Timings.USE_DLL_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = -2;
+                end
+            catch E
+                fprintf(1, 'Simulation failed:\n\n%s', E.message);
+                Timings.USE_DLL_WITH_BLOCKS_preprocessing(id(stack_solve_algo==algo), mfs+1) = info.time.preprocessor;
+                Timings.USE_DLL_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = -1;
+            end
+        end
+    end
+    system(sprintf('rm -rf %s', model));
+    system(sprintf('rm -rf +%s', model));
+end
+%
+% Bytecode with blocks
+%
+for mfs = 0:3
+    for algo = [0, 1, 2, 3, 4, 5]
+        if isnan(Timings.BYTECODE_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1))
+            try
+                disp_title(sprintf('| BYTECODE (solve_algo=%u, with blocks, mfs=%u, model=%s) |', algo, mfs, model));
+                dprintf('dynare %s onlyclearglobals -DBYTECODE=true -DITERATIONS=%u -DSTACK_SOLVE_ALGO_VALUE=%u -DMFS_VALUE=%u -DBLOCKS=true', model, niter, algo, mfs);
+                skipline()
+                Timings.BYTECODE_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = -3; % Provision for MATLAB crash
+                info = dynare(sprintf('%s', model), 'onlyclearglobals', sprintf('-DITERATIONS=%u', niter), '-DBLOCKS=true', sprintf('-DMFS_VALUE=%u', mfs), '-DBYTECODE=true', sprintf('-DSTACK_SOLVE_ALGO_VALUE=%u', algo));
+                if oo_.deterministic_simulation.error<ccrit
+                    Timings.BYTECODE_WITH_BLOCKS_preprocessing(id(stack_solve_algo==algo), mfs+1) = info.time.preprocessor;
+                    Timings.BYTECODE_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = info.time.compute/niter;
+                else
+                    Timings.BYTECODE_WITH_BLOCKS_preprocessing(id(stack_solve_algo==algo), mfs+1) = info.time.preprocessor;
+                    Timings.BYTECODE_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = -2;
+                end
+            catch E
+                fprintf(1, 'Simulation failed:\n\n%s', E.message);
+                Timings.BYTECODE_WITH_BLOCKS_preprocessing(id(stack_solve_algo==algo), mfs+1) = info.time.preprocessor;
+                Timings.BYTECODE_WITH_BLOCKS_computing(id(stack_solve_algo==algo), mfs+1) = -1;
+            end
+        end
+    end
+end
+%
+% Touch a file if everything is done.
+%
+system(sprintf('touch %s-done.info', model))