Bytecode: the block and non-block versions now coexist in parallel

As a consequence, a new “block_decomposed” option of the bytecode MEX has been introduced to explicitly select the block-decomposed version. Note that we do not always use the “block_decomposed” option even when the “block” option has been passed to the user, in situations where the block decomposition brings nothing (e.g. when evaluating the residuals of the whole model).
2023-01-09 15:19:36 +01:00 · 2023-01-09 15:19:36 +01:00 · f84753025d
parent f769420888
commit f84753025d
14 changed files with 78 additions and 68 deletions
--- a/matlab/block_bytecode_mfs_steadystate.m
+++ b/matlab/block_bytecode_mfs_steadystate.m
@ -2,7 +2,7 @@ function [r, g1] = block_bytecode_mfs_steadystate(y, b, y_all, exo, params, T, M
 % Wrapper around the *_static.m file, for use with dynare_solve,
 % when block_mfs option is given to steady.

-% Copyright © 2009-2022 Dynare Team
+% Copyright © 2009-2023 Dynare Team
 %
 % This file is part of Dynare.
 %
@ -21,4 +21,4 @@ function [r, g1] = block_bytecode_mfs_steadystate(y, b, y_all, exo, params, T, M

 indx = M.block_structure_stat.block(b).variable;
 y_all(indx) = y;
-[r, g1] = bytecode(y_all, exo, params, y_all, 1, y_all, T, 'evaluate', 'static', ['block = ' int2str(b) ]);
+[r, g1] = bytecode(y_all, exo, params, y_all, 1, y_all, T, 'evaluate', 'static', 'block_decomposed', ['block = ' int2str(b) ]);
--- a/matlab/dr_block.m
+++ b/matlab/dr_block.m
@ -36,7 +36,7 @@ function [dr,info,M_,oo_] = dr_block(dr,task,M_,options_,oo_,varargin)
 %   none.
 %

-% Copyright © 2010-2022 Dynare Team
+% Copyright © 2010-2023 Dynare Team
 %
 % This file is part of Dynare.
 %
@ -79,7 +79,7 @@ end
 data = M_.block_structure.block;

 if options_.bytecode
-    [~, data]= bytecode('dynamic','evaluate', z, zx, M_.params, dr.ys, 1, data);
+    [~, data]= bytecode('dynamic', 'evaluate', 'block_decomposed', z, zx, M_.params, dr.ys, 1, data);
 else
    T=NaN(M_.block_structure.dyn_tmp_nbr, 1);
    it_=M_.maximum_lag+1;
--- a/matlab/dynare_solve_block_or_bytecode.m
+++ b/matlab/dynare_solve_block_or_bytecode.m
@ -1,6 +1,6 @@
 function [x,info] = dynare_solve_block_or_bytecode(y, exo, params, options, M)

-% Copyright © 2010-2022 Dynare Team
+% Copyright © 2010-2023 Dynare Team
 %
 % This file is part of Dynare.
 %
@ -54,7 +54,11 @@ if options.block && ~options.bytecode
 elseif options.bytecode
    if options.solve_algo >= 5 && options.solve_algo <= 8
        try
-            x = bytecode('static', x, exo, params);
+            if options.block
+                x = bytecode('static', 'block_decomposed', x, exo, params);
+            else
+                x = bytecode('static', x, exo, params);
+            end
        catch ME
            disp(ME.message);
            info = 1;
@ -79,7 +83,7 @@ elseif options.bytecode
            % Also update the temporary terms vector (needed for the dynare_solve case)
            try
                [~, ~, x, T] = bytecode(x, exo, params, x, 1, x, T, 'evaluate', 'static', ...
-                                        ['block = ' int2str(b)]);
+                                        'block_decomposed', ['block = ' int2str(b)]);
            catch ME
                disp(ME.message);
                info = 1;
--- a/matlab/model_diagnostics.m
+++ b/matlab/model_diagnostics.m
@ -16,7 +16,7 @@ function model_diagnostics(M,options,oo)
 %   none.
 %

-% Copyright © 1996-2020 Dynare Team
+% Copyright © 1996-2023 Dynare Team
 %
 % This file is part of Dynare.
 %
@ -153,7 +153,7 @@ for b=1:nb
                                    'evaluate', 'static');
        else
            [res, jacob] = bytecode(dr.ys, exo, M.params, dr.ys, 1, exo, ...
-                                    'evaluate', 'static',['block=' ...
+                                    'evaluate', 'static', 'block_decomposed', ['block=' ...
                                int2str(b)]);
        end
    else
--- a/matlab/perfect-foresight-models/det_cond_forecast.m
+++ b/matlab/perfect-foresight-models/det_cond_forecast.m
@ -12,7 +12,7 @@ function data_set = det_cond_forecast(varargin)
 %  dataset                [dseries]     Returns a dseries containing the forecasted endgenous variables and shocks
 %

-% Copyright © 2013-2020 Dynare Team
+% Copyright © 2013-2023 Dynare Team
 %
 % This file is part of Dynare.
 %
@ -161,7 +161,11 @@ else
            if options_.bytecode
                save_options_dynatol_f = options_.dynatol.f;
                options_.dynatol.f = 1e-7;
-                [endo, exo] = bytecode('extended_path', plan, oo_.endo_simul, oo_.exo_simul, M_.params, oo_.steady_state, options_.periods);
+                if options_.block
+                    [endo, exo] = bytecode('extended_path', 'block_decomposed', plan, oo_.endo_simul, oo_.exo_simul, M_.params, oo_.steady_state, options_.periods);
+                else
+                    [endo, exo] = bytecode('extended_path', plan, oo_.endo_simul, oo_.exo_simul, M_.params, oo_.steady_state, options_.periods);
+                end
                options_.dynatol.f = save_options_dynatol_f;

                oo_.endo_simul = endo;
@ -466,16 +470,13 @@ if pf && ~surprise
        indx_x = [];
        for k = 1 : per
            if k == 1
-                if (isfield(M_,'block_structure'))
-                    data1 = M_.block_structure.block;
-                    Size = length(M_.block_structure.block);
-                else
-                    data1 = M_;
-                    Size = 1;
-                end
                data1 = M_;
                if (options_.bytecode)
-                    [zz, data1]= bytecode('dynamic','evaluate', z, zx, M_.params, oo_.steady_state, k, data1);
+                    if options_.block
+                        [zz, data1]= bytecode('dynamic','block_decomposed','evaluate', z, zx, M_.params, oo_.steady_state, k, data1);
+                    else
+                        [zz, data1]= bytecode('dynamic','evaluate', z, zx, M_.params, oo_.steady_state, k, data1);
+                    end
                else
                    [zz, g1b] = feval([M_.fname '.dynamic'], z', zx, M_.params, oo_.steady_state, k);
                    data1.g1_x = g1b(:,end - M_.exo_nbr + 1:end);
@ -736,16 +737,13 @@ else
            indx_x = [];
            for k = 1 : per
                if k == 1
-                    if (isfield(M_,'block_structure'))
-                        data1 = M_.block_structure.block;
-                        Size = length(M_.block_structure.block);
-                    else
-                        data1 = M_;
-                        Size = 1;
-                    end
                    data1 = M_;
                    if (options_.bytecode)
-                        [zz, data1]= bytecode('dynamic','evaluate', z, zx, M_.params, oo_.steady_state, k, data1);
+                        if options_.block
+                            [zz, data1]= bytecode('dynamic','block_decomposed','evaluate', z, zx, M_.params, oo_.steady_state, k, data1);
+                        else
+                            [zz, data1]= bytecode('dynamic','evaluate', z, zx, M_.params, oo_.steady_state, k, data1);
+                        end
                    else
                        [zz, g1b] = feval([M_.fname '.dynamic'], z', zx, M_.params, oo_.steady_state, k);
                        data1.g1_x = g1b(:,end - M_.exo_nbr + 1:end);
--- a/matlab/perfect-foresight-models/perfect_foresight_solver_core.m
+++ b/matlab/perfect-foresight-models/perfect_foresight_solver_core.m
@ -11,7 +11,7 @@ function [oo_, maxerror] = perfect_foresight_solver_core(M_, options_, oo_)
 % - oo_                 [struct] contains results
 % - maxerror            [double] contains the maximum absolute error

-% Copyright © 2015-2022 Dynare Team
+% Copyright © 2015-2023 Dynare Team
 %
 % This file is part of Dynare.
 %
@ -55,7 +55,7 @@ if options_.block
    end
    if options_.bytecode
        try
-            oo_.endo_simul = bytecode('dynamic', oo_.endo_simul, oo_.exo_simul, M_.params, repmat(oo_.steady_state,1, periods+2), periods);
+            oo_.endo_simul = bytecode('dynamic', 'block_decomposed', oo_.endo_simul, oo_.exo_simul, M_.params, repmat(oo_.steady_state,1, periods+2), periods);
            oo_.deterministic_simulation.status = true;
        catch ME
            disp(ME.message)
--- a/matlab/print_bytecode_dynamic_model.m
+++ b/matlab/print_bytecode_dynamic_model.m
@ -11,7 +11,7 @@ function print_bytecode_dynamic_model()
 % SPECIAL REQUIREMENTS
 %   none

-% Copyright © 2001-2017 Dynare Team
+% Copyright © 2001-2023 Dynare Team
 %
 % This file is part of Dynare.
 %
@ -29,7 +29,11 @@ function print_bytecode_dynamic_model()
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 global options_
 if options_.bytecode
-    bytecode('print','dynamic');
+    if options_.block
+        bytecode('print','dynamic','block_decomposed');
+    else
+        bytecode('print','dynamic');
+    end
 else
    disp('You have to use bytecode option in model command to use print_bytecode_dynamic_model');
 end
--- a/matlab/print_bytecode_static_model.m
+++ b/matlab/print_bytecode_static_model.m
@ -11,7 +11,7 @@ function print_bytecode_static_model()
 % SPECIAL REQUIREMENTS
 %   none

-% Copyright © 2001-2017 Dynare Team
+% Copyright © 2001-2023 Dynare Team
 %
 % This file is part of Dynare.
 %
@ -29,7 +29,11 @@ function print_bytecode_static_model()
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 global options_
 if options_.bytecode
-    bytecode('print','static');
+    if options_.block
+        bytecode('print','static','block_decomposed');
+    else
+        bytecode('print','static');
+    end
 else
    disp('You have to use bytecode option in model command to use print_bytecode_static_model');
 end
--- a/mex/sources/bytecode/Interpreter.cc
+++ b/mex/sources/bytecode/Interpreter.cc
@ -1,5 +1,5 @@
 /*
- * Copyright © 2007-2022 Dynare Team
+ * Copyright © 2007-2023 Dynare Team
 *
 * This file is part of Dynare.
 *
@ -20,6 +20,7 @@
 #include <sstream>
 #include <algorithm>
 #include <cstring>
+#include <filesystem>

 #include "Interpreter.hh"

@ -31,8 +32,8 @@ Interpreter::Interpreter(double *params_arg, double *y_arg, double *ya_arg, doub
                         int maxit_arg_, double solve_tolf_arg, size_t size_of_direction_arg, int y_decal_arg, double markowitz_c_arg,
                         string &filename_arg, int minimal_solving_periods_arg, int stack_solve_algo_arg, int solve_algo_arg,
                         bool global_temporary_terms_arg, bool print_arg, bool print_error_arg, mxArray *GlobalTemporaryTerms_arg,
-                         bool steady_state_arg, bool print_it_arg, int col_x_arg, int col_y_arg, BasicSymbolTable &symbol_table_arg)
-: dynSparseMatrix {y_size_arg, y_kmin_arg, y_kmax_arg, print_it_arg, steady_state_arg, periods_arg, minimal_solving_periods_arg, symbol_table_arg}
+                         bool steady_state_arg, bool block_decomposed_arg, bool print_it_arg, int col_x_arg, int col_y_arg, BasicSymbolTable &symbol_table_arg)
+: dynSparseMatrix {y_size_arg, y_kmin_arg, y_kmax_arg, print_it_arg, steady_state_arg, block_decomposed_arg, periods_arg, minimal_solving_periods_arg, symbol_table_arg}
 {
  params = params_arg;
  y = y_arg;
@ -551,18 +552,16 @@ Interpreter::print_a_block()
 }

 void
-Interpreter::ReadCodeFile(string file_name, CodeLoad &code)
+Interpreter::ReadCodeFile(const string &file_name, CodeLoad &code)
 {
-  if (steady_state)
-    file_name += "/model/bytecode/static";
-  else
-    file_name += "/model/bytecode/dynamic";
+  filesystem::path codfile {file_name + "/model/bytecode/" + (block_decomposed ? "block/" : "")
+    + (steady_state ? "static" : "dynamic") + ".cod"};

  //First read and store in memory the code
-  code_liste = code.get_op_code(file_name);
+  code_liste = code.get_op_code(codfile);
  EQN_block_number = code.get_block_number();
  if (!code_liste.size())
-    throw FatalException{"In compute_blocks, " + file_name + ".cod cannot be opened"};
+    throw FatalException{"In compute_blocks, " + codfile.string() + " cannot be opened"};
  if (block >= code.get_block_number())
    throw FatalException{"In compute_blocks, input argument block = " + to_string(block+1)
                         + " is greater than the number of blocks in the model ("
--- a/mex/sources/bytecode/Interpreter.hh
+++ b/mex/sources/bytecode/Interpreter.hh
@ -1,5 +1,5 @@
 /*
- * Copyright © 2007-2021 Dynare Team
+ * Copyright © 2007-2023 Dynare Team
 *
 * This file is part of Dynare.
 *
@ -51,12 +51,12 @@ public:
              int maxit_arg_, double solve_tolf_arg, size_t size_of_direction_arg, int y_decal_arg, double markowitz_c_arg,
              string &filename_arg, int minimal_solving_periods_arg, int stack_solve_algo_arg, int solve_algo_arg,
              bool global_temporary_terms_arg, bool print_arg, bool print_error_arg, mxArray *GlobalTemporaryTerms_arg,
-              bool steady_state_arg, bool print_it_arg, int col_x_arg, int col_y_arg, BasicSymbolTable &symbol_table_arg);
+              bool steady_state_arg, bool block_decomposed_arg, bool print_it_arg, int col_x_arg, int col_y_arg, BasicSymbolTable &symbol_table_arg);
  bool extended_path(const string &file_name, const string &bin_basename, bool evaluate, int block, int &nb_blocks, int nb_periods, const vector<s_plan> &sextended_path, const vector<s_plan> &sconstrained_extended_path, const vector<string> &dates, const table_conditional_global_type &table_conditional_global);
  bool compute_blocks(const string &file_name, const string &bin_basename, bool evaluate, int block, int &nb_blocks);
  void check_for_controlled_exo_validity(FBEGINBLOCK_ *fb, const vector<s_plan> &sconstrained_extended_path);
  bool MainLoop(const string &bin_basename, const CodeLoad &code, bool evaluate, int block, bool last_call, bool constrained, const vector<s_plan> &sconstrained_extended_path, const vector_table_conditional_local_type &vector_table_conditional_local);
-  void ReadCodeFile(string file_name, CodeLoad &code);
+  void ReadCodeFile(const string &file_name, CodeLoad &code);

  inline mxArray *
  get_jacob(int block_num) const
--- a/mex/sources/bytecode/SparseMatrix.cc
+++ b/mex/sources/bytecode/SparseMatrix.cc
@ -1,5 +1,5 @@
 /*
- * Copyright © 2007-2022 Dynare Team
+ * Copyright © 2007-2023 Dynare Team
 *
 * This file is part of Dynare.
 *
@ -19,12 +19,14 @@

 #include <sstream>
 #include <algorithm>
+#include <filesystem>

 #include "SparseMatrix.hh"

-dynSparseMatrix::dynSparseMatrix(int y_size_arg, int y_kmin_arg, int y_kmax_arg, bool print_it_arg, bool steady_state_arg, int periods_arg,
+dynSparseMatrix::dynSparseMatrix(int y_size_arg, int y_kmin_arg, int y_kmax_arg, bool print_it_arg, bool steady_state_arg, bool block_decomposed_arg, int periods_arg,
                                 int minimal_solving_periods_arg, BasicSymbolTable &symbol_table_arg) :
-  Evaluate {y_size_arg, y_kmin_arg, y_kmax_arg, print_it_arg, steady_state_arg, periods_arg, minimal_solving_periods_arg, symbol_table_arg}
+  Evaluate {y_size_arg, y_kmin_arg, y_kmax_arg, print_it_arg, steady_state_arg, periods_arg, minimal_solving_periods_arg, symbol_table_arg},
+  block_decomposed {block_decomposed_arg}
 {
  pivotva = nullptr;
  g_save_op = nullptr;
@ -288,19 +290,11 @@ dynSparseMatrix::Read_SparseMatrix(const string &file_name, int Size, int period
  mem_mngr.fixe_file_name(file_name);
  if (!SaveCode.is_open())
    {
-      if (steady_state)
-        SaveCode.open(file_name + "/model/bytecode/static.bin", ios::in | ios::binary);
-      else
-        SaveCode.open(file_name + "/model/bytecode/dynamic.bin", ios::in | ios::binary);
+      filesystem::path binfile {file_name + "/model/bytecode/" + (block_decomposed ? "block/" : "")
+        + (steady_state ? "static" : "dynamic") + ".bin"};
+      SaveCode.open(binfile, ios::in | ios::binary);
      if (!SaveCode.is_open())
-        {
-          ostringstream tmp;
-          if (steady_state)
-            tmp << "In Read_SparseMatrix, " << file_name << "/model/bytecode/static.bin cannot be opened";
-          else
-            tmp << "In Read_SparseMatrix, " << file_name << "/model/bytecode/dynamic.bin cannot be opened";
-          throw FatalException{tmp.str()};
-        }
+        throw FatalException{"In Read_SparseMatrix, " + binfile.string() + " cannot be opened"};
    }
  IM_i.clear();
  if (two_boundaries)
--- a/mex/sources/bytecode/SparseMatrix.hh
+++ b/mex/sources/bytecode/SparseMatrix.hh
@ -1,5 +1,5 @@
 /*
- * Copyright © 2007-2022 Dynare Team
+ * Copyright © 2007-2023 Dynare Team
 *
 * This file is part of Dynare.
 *
@ -68,7 +68,7 @@ constexpr double mem_increasing_factor = 1.1;
 class dynSparseMatrix : public Evaluate
 {
 public:
-  dynSparseMatrix(int y_size_arg, int y_kmin_arg, int y_kmax_arg, bool print_it_arg, bool steady_state_arg, int periods_arg, int minimal_solving_periods_arg, BasicSymbolTable &symbol_table_arg);
+  dynSparseMatrix(int y_size_arg, int y_kmin_arg, int y_kmax_arg, bool print_it_arg, bool steady_state_arg, bool block_decomposed_arg, int periods_arg, int minimal_solving_periods_arg, BasicSymbolTable &symbol_table_arg);
  void Simulate_Newton_Two_Boundaries(int blck, int y_size, int y_kmin, int y_kmax, int Size, int periods, bool cvg, int minimal_solving_periods, int stack_solve_algo, const vector_table_conditional_local_type &vector_table_conditional_local);
  void Simulate_Newton_One_Boundary(bool forward);
  void fixe_u(double **u, int u_count_int, int max_lag_plus_max_lead_plus_1);
@ -142,6 +142,9 @@ private:
  // Computes A−B where A and B are dense. The result is dense.
  static mxArray *subtract_A_B(const mxArray *A_m, const mxArray *B_m);
 protected:
+  // Whether to use the block-decomposed version of the bytecode file
+  bool block_decomposed;
+
  stack<double> Stack;
  int nb_prologue_table_u, nb_first_table_u, nb_middle_table_u, nb_last_table_u;
  int nb_prologue_table_y, nb_first_table_y, nb_middle_table_y, nb_last_table_y;
--- a/mex/sources/bytecode/bytecode.cc
+++ b/mex/sources/bytecode/bytecode.cc
@ -1,5 +1,5 @@
 /*
- * Copyright © 2007-2022 Dynare Team
+ * Copyright © 2007-2023 Dynare Team
 *
 * This file is part of Dynare.
 *
@ -58,7 +58,8 @@ Get_Arguments_and_global_variables(int nrhs,
                                   double *steady_yd[], size_t &steady_row_y, size_t &steady_col_y,
                                   unsigned int &periods,
                                   mxArray *block_structur[],
-                                   bool &steady_state, bool &evaluate, int &block,
+                                   bool &steady_state, bool &block_decomposed,
+                                   bool &evaluate, int &block,
                                   mxArray *M_[], mxArray *oo_[], mxArray *options_[], bool &global_temporary_terms,
                                   bool &print,
                                   bool &print_error,
@ -116,6 +117,8 @@ Get_Arguments_and_global_variables(int nrhs,
          steady_state = true;
        else if (Get_Argument(prhs[i]) == "dynamic")
          steady_state = false;
+        else if (Get_Argument(prhs[i]) == "block_decomposed")
+          block_decomposed = true;
        else if (Get_Argument(prhs[i]) == "evaluate")
          evaluate = true;
        else if (Get_Argument(prhs[i]) == "global_temporary_terms")
@ -204,6 +207,7 @@ mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
  unsigned int periods = 1;
  double *direction;
  bool steady_state = false;
+  bool block_decomposed {false};
  bool evaluate = false;
  int block = -1;
  double *params = nullptr;
@ -237,7 +241,7 @@ mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
                                         &steady_yd, steady_row_y, steady_col_y,
                                         periods,
                                         &block_structur,
-                                         steady_state, evaluate, block,
+                                         steady_state, block_decomposed, evaluate, block,
                                         &M_, &oo_, &options_, global_temporary_terms,
                                         print, print_error, &GlobalTemporaryTerms,
                                         &plan, &pfplan, &extended_path, &extended_path_struct);
@ -715,7 +719,7 @@ mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
                          periods, y_kmin, y_kmax, maxit_, solve_tolf, size_of_direction, y_decal,
                          markowitz_c, file_name, minimal_solving_periods, stack_solve_algo,
                          solve_algo, global_temporary_terms, print, print_error, GlobalTemporaryTerms,
-                          steady_state, print_it, col_x, col_y, symbol_table};
+                          steady_state, block_decomposed, print_it, col_x, col_y, symbol_table};
  string f(fname);
  mxFree(fname);
  int nb_blocks = 0;
--- a/2
+++ b/2
@ -1 +1 @@
-Subproject commit 4348f4d57f83aa8f96c3ec23899c33971e847428
+Subproject commit a4f299c0b781676afe3d9d073a7e2cd0551af997