diff --git a/src/CodeInterpreter.hh b/src/CodeInterpreter.hh
index 4b9911c0b909370b88ddc631bfe84828cc60e4ac..ee4ad6a21ef0bc10ea60f8b37aa5a7cd5ee7c8e2 100644
--- a/src/CodeInterpreter.hh
+++ b/src/CodeInterpreter.hh
@@ -1,5 +1,5 @@
 /*
- * Copyright © 2007-2019 Dynare Team
+ * Copyright © 2007-2020 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -121,17 +121,17 @@ enum EquationType
    E_SOLVE //!< No simple evaluation of the equation, it has to be solved
   };
 
-enum BlockSimulationType
+enum class BlockSimulationType
   {
-   UNKNOWN, //!< Unknown simulation type
-   EVALUATE_FORWARD, //!< Simple evaluation, normalized variable on left-hand side, forward
-   EVALUATE_BACKWARD, //!< Simple evaluation, normalized variable on left-hand side, backward
-   SOLVE_FORWARD_SIMPLE, //!< Block of one equation, newton solver needed, forward
-   SOLVE_BACKWARD_SIMPLE, //!< Block of one equation, newton solver needed, backward
-   SOLVE_TWO_BOUNDARIES_SIMPLE, //!< Block of one equation, newton solver needed, forward & ackward
-   SOLVE_FORWARD_COMPLETE, //!< Block of several equations, newton solver needed, forward
-   SOLVE_BACKWARD_COMPLETE, //!< Block of several equations, newton solver needed, backward
-   SOLVE_TWO_BOUNDARIES_COMPLETE //!< Block of several equations, newton solver needed, forward and backwar
+   unknown, //!< Unknown simulation type
+   evaluateForward, //!< Simple evaluation, normalized variable on left-hand side, forward
+   evaluateBackward, //!< Simple evaluation, normalized variable on left-hand side, backward
+   solveForwardSimple, //!< Block of one equation, newton solver needed, forward
+   solveBackwardSimple, //!< Block of one equation, newton solver needed, backward
+   solveTwoBoundariesSimple, //!< Block of one equation, newton solver needed, forward & ackward
+   solveForwardComplete, //!< Block of several equations, newton solver needed, forward
+   solveBackwardComplete, //!< Block of several equations, newton solver needed, backward
+   solveTwoBoundariesComplete //!< Block of several equations, newton solver needed, forward and backwar
   };
 
 //! Enumeration of possible symbol types
@@ -1428,7 +1428,7 @@ private:
   unsigned int nb_col_det_exo_jacob, nb_col_exo_jacob, nb_col_other_endo_jacob;
 public:
   inline
-  FBEGINBLOCK_() : size{0}, type{UNKNOWN},
+  FBEGINBLOCK_() : size{0}, type{static_cast<uint8_t>(BlockSimulationType::unknown)},
                    is_linear{false}, endo_nbr{0}, Max_Lag{0}, Max_Lead{0}, u_count_int{0}, nb_col_jacob{0}
   {
   }
@@ -1577,8 +1577,10 @@ public:
         CompileCode.write(reinterpret_cast<char *>(&variable[i]), sizeof(variable[0]));
         CompileCode.write(reinterpret_cast<char *>(&equation[i]), sizeof(equation[0]));
       }
-    if (type == SOLVE_TWO_BOUNDARIES_SIMPLE || type == SOLVE_TWO_BOUNDARIES_COMPLETE
-        || type == SOLVE_BACKWARD_COMPLETE || type == SOLVE_FORWARD_COMPLETE)
+    if (type == static_cast<uint8_t>(BlockSimulationType::solveTwoBoundariesSimple)
+        || type == static_cast<uint8_t>(BlockSimulationType::solveTwoBoundariesComplete)
+        || type == static_cast<uint8_t>(BlockSimulationType::solveBackwardComplete)
+        || type == static_cast<uint8_t>(BlockSimulationType::solveForwardComplete))
       {
         CompileCode.write(reinterpret_cast<char *>(&is_linear), sizeof(is_linear));
         CompileCode.write(reinterpret_cast<char *>(&endo_nbr), sizeof(endo_nbr));
@@ -1617,8 +1619,10 @@ public:
         memcpy(&bc.Equation, code, sizeof(bc.Equation)); code += sizeof(bc.Equation);
         Block_Contain_.push_back(bc);
       }
-    if (type == SOLVE_TWO_BOUNDARIES_SIMPLE || type == SOLVE_TWO_BOUNDARIES_COMPLETE
-        || type == SOLVE_BACKWARD_COMPLETE || type == SOLVE_FORWARD_COMPLETE)
+    if (type == static_cast<uint8_t>(BlockSimulationType::solveTwoBoundariesSimple)
+        || type == static_cast<uint8_t>(BlockSimulationType::solveTwoBoundariesComplete)
+        || type == static_cast<uint8_t>(BlockSimulationType::solveBackwardComplete)
+        || type == static_cast<uint8_t>(BlockSimulationType::solveForwardComplete))
       {
         memcpy(&is_linear, code, sizeof(is_linear)); code += sizeof(is_linear);
         memcpy(&endo_nbr, code, sizeof(endo_nbr)); code += sizeof(endo_nbr);
diff --git a/src/DynamicModel.cc b/src/DynamicModel.cc
index 7af8eab996dbc37a85fe3c297389b0f8470ee36a..c60320ef0175a0393fee61d85b863f91251c875d 100644
--- a/src/DynamicModel.cc
+++ b/src/DynamicModel.cc
@@ -419,27 +419,36 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
              << "% Warning : this file is generated automatically by Dynare" << endl
              << "%           from model file (.mod)" << endl << endl
              << "%/" << endl;
-      if (simulation_type == EVALUATE_BACKWARD || simulation_type == EVALUATE_FORWARD)
+      if (simulation_type == BlockSimulationType::evaluateBackward
+          || simulation_type == BlockSimulationType::evaluateForward)
         {
           output << "function [y, g1, g2, g3, varargout] = dynamic_" << block+1 << "(y, x, params, steady_state, jacobian_eval, y_kmin, periods)" << endl;
         }
-      else if (simulation_type == SOLVE_FORWARD_COMPLETE || simulation_type == SOLVE_BACKWARD_COMPLETE)
+      else if (simulation_type == BlockSimulationType::solveForwardComplete
+               || simulation_type == BlockSimulationType::solveBackwardComplete)
         output << "function [residual, y, g1, g2, g3, varargout] = dynamic_" << block+1 << "(y, x, params, steady_state, it_, jacobian_eval)" << endl;
-      else if (simulation_type == SOLVE_BACKWARD_SIMPLE || simulation_type == SOLVE_FORWARD_SIMPLE)
+      else if (simulation_type == BlockSimulationType::solveBackwardSimple
+               || simulation_type == BlockSimulationType::solveForwardSimple)
         output << "function [residual, y, g1, g2, g3, varargout] = dynamic_" << block+1 << "(y, x, params, steady_state, it_, jacobian_eval)" << endl;
       else
         output << "function [residual, y, g1, g2, g3, b, varargout] = dynamic_" << block+1 << "(y, x, params, steady_state, periods, jacobian_eval, y_kmin, y_size, Periods)" << endl;
       BlockType block_type;
-      if (simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE)
+      if (simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+          || simulation_type == BlockSimulationType::solveTwoBoundariesSimple)
         block_type = BlockType::simultan;
-      else if (simulation_type == SOLVE_FORWARD_COMPLETE || simulation_type == SOLVE_BACKWARD_COMPLETE)
+      else if (simulation_type == BlockSimulationType::solveForwardComplete
+               || simulation_type == BlockSimulationType::solveBackwardComplete)
         block_type = BlockType::simultans;
-      else if ((simulation_type == SOLVE_FORWARD_SIMPLE || simulation_type == SOLVE_BACKWARD_SIMPLE
-                || simulation_type == EVALUATE_BACKWARD || simulation_type == EVALUATE_FORWARD)
+      else if ((simulation_type == BlockSimulationType::solveForwardSimple
+                || simulation_type == BlockSimulationType::solveBackwardSimple
+                || simulation_type == BlockSimulationType::evaluateBackward
+                || simulation_type == BlockSimulationType::evaluateForward)
                && getBlockFirstEquation(block) < prologue)
         block_type = BlockType::prologue;
-      else if ((simulation_type == SOLVE_FORWARD_SIMPLE || simulation_type == SOLVE_BACKWARD_SIMPLE
-                || simulation_type == EVALUATE_BACKWARD || simulation_type == EVALUATE_FORWARD)
+      else if ((simulation_type == BlockSimulationType::solveForwardSimple
+                || simulation_type == BlockSimulationType::solveBackwardSimple
+                || simulation_type == BlockSimulationType::evaluateBackward
+                || simulation_type == BlockSimulationType::evaluateForward)
                && getBlockFirstEquation(block) >= equations.size() - epilogue)
         block_type = BlockType::epilogue;
       else
@@ -451,7 +460,8 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
              << BlockSim(simulation_type) << "  //" << endl
              << "  % ////////////////////////////////////////////////////////////////////////" << endl;
       //The Temporary terms
-      if (simulation_type == EVALUATE_BACKWARD || simulation_type == EVALUATE_FORWARD)
+      if (simulation_type == BlockSimulationType::evaluateBackward
+          || simulation_type == BlockSimulationType::evaluateForward)
         {
           output << "  if(jacobian_eval)" << endl
                  << "    g1 = spalloc(" << block_mfs  << ", " << count_col_endo << ", " << nze << ");" << endl
@@ -468,7 +478,8 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
                  << "    g1_xd=spalloc(" << block_size << ", " << count_col_exo_det  << ", " << nze_exo_det << ");" << endl
                  << "    g1_o=spalloc(" << block_size << ", " << count_col_other_endo << ", " << nze_other_endo << ");" << endl
                  << "  else" << endl;
-          if (simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE)
+          if (simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+              || simulation_type == BlockSimulationType::solveTwoBoundariesSimple)
             output << "    g1 = spalloc(" << block_mfs << "*Periods, "
                    << block_mfs << "*(Periods+" << max_leadlag_block[block].first+max_leadlag_block[block].second+1 << ")"
                    << ", " << nze << "*Periods);" << endl;
@@ -486,7 +497,8 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
             tmp_output << " T" << it;
           output << "  global" << tmp_output.str() << ";" << endl;
         }
-      if (simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE)
+      if (simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+          || simulation_type == BlockSimulationType::solveTwoBoundariesSimple)
         {
           temporary_terms_t tt2;
           for (int i = 0; i < static_cast<int>(block_size); i++)
@@ -505,16 +517,21 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
                 }
             }
         }
-      if (simulation_type == SOLVE_BACKWARD_SIMPLE || simulation_type == SOLVE_FORWARD_SIMPLE || simulation_type == SOLVE_BACKWARD_COMPLETE || simulation_type == SOLVE_FORWARD_COMPLETE)
+      if (simulation_type == BlockSimulationType::solveBackwardSimple
+          || simulation_type == BlockSimulationType::solveForwardSimple
+          || simulation_type == BlockSimulationType::solveBackwardComplete
+          || simulation_type == BlockSimulationType::solveForwardComplete)
         output << "  residual=zeros(" << block_mfs << ",1);" << endl;
-      else if (simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE)
+      else if (simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+               || simulation_type == BlockSimulationType::solveTwoBoundariesSimple)
         output << "  residual=zeros(" << block_mfs << ",y_kmin+periods);" << endl;
-      if (simulation_type == EVALUATE_BACKWARD)
+      if (simulation_type == BlockSimulationType::evaluateBackward)
         output << "  for it_ = (y_kmin+periods):y_kmin+1" << endl;
-      if (simulation_type == EVALUATE_FORWARD)
+      if (simulation_type == BlockSimulationType::evaluateForward)
         output << "  for it_ = y_kmin+1:(y_kmin+periods)" << endl;
 
-      if (simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE)
+      if (simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+          || simulation_type == BlockSimulationType::solveTwoBoundariesSimple)
         {
           output << "  b = zeros(periods*y_size,1);" << endl
                  << "  for it_ = y_kmin+1:(periods+y_kmin)" << endl
@@ -524,7 +541,8 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
           sps = "  ";
         }
       else
-        if (simulation_type == EVALUATE_BACKWARD || simulation_type == EVALUATE_FORWARD)
+        if (simulation_type == BlockSimulationType::evaluateBackward
+            || simulation_type == BlockSimulationType::evaluateForward)
           sps = "  ";
         else
           sps = "";
@@ -562,9 +580,10 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
           lhs->writeOutput(tmp_output, local_output_type, local_temporary_terms, {});
           switch (simulation_type)
             {
-            case EVALUATE_BACKWARD:
-            case EVALUATE_FORWARD:
-            evaluation:     if (simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE)
+            case BlockSimulationType::evaluateBackward:
+            case BlockSimulationType::evaluateForward:
+            evaluation:     if (simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+                                || simulation_type == BlockSimulationType::solveTwoBoundariesSimple)
                 output << "    % equation " << getBlockEquationID(block, i)+1 << " variable : " << sModel
                        << " (" << variable_ID+1 << ") " << c_Equation_Type(equ_type) << endl;
               output << "    ";
@@ -598,10 +617,10 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
                 }
               output << ";" << endl;
               break;
-            case SOLVE_BACKWARD_SIMPLE:
-            case SOLVE_FORWARD_SIMPLE:
-            case SOLVE_BACKWARD_COMPLETE:
-            case SOLVE_FORWARD_COMPLETE:
+            case BlockSimulationType::solveBackwardSimple:
+            case BlockSimulationType::solveForwardSimple:
+            case BlockSimulationType::solveBackwardComplete:
+            case BlockSimulationType::solveForwardComplete:
               if (i < block_recursive)
                 goto evaluation;
               feedback_variables.push_back(variable_ID);
@@ -609,8 +628,8 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
                      << " (" << variable_ID+1 << ") " << c_Equation_Type(equ_type) << " symb_id=" << symbol_table.getID(SymbolType::endogenous, variable_ID) << endl;
               output << "  " << "residual(" << i+1-block_recursive << ") = (";
               goto end;
-            case SOLVE_TWO_BOUNDARIES_COMPLETE:
-            case SOLVE_TWO_BOUNDARIES_SIMPLE:
+            case BlockSimulationType::solveTwoBoundariesComplete:
+            case BlockSimulationType::solveTwoBoundariesSimple:
               if (i < block_recursive)
                 goto evaluation;
               feedback_variables.push_back(variable_ID);
@@ -628,17 +647,21 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
               rhs->writeOutput(output, local_output_type, local_temporary_terms, {});
               output << ");" << endl;
 #ifdef CONDITION
-              if (simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE)
+              if (simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+                  || simulation_type == BlockSimulationType::solveTwoBoundariesSimple)
                 output << "  condition(" << i+1 << ")=0;" << endl;
 #endif
             }
         }
       // The Jacobian if we have to solve the block
-      if (simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE || simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE)
+      if (simulation_type == BlockSimulationType::solveTwoBoundariesSimple
+          || simulation_type == BlockSimulationType::solveTwoBoundariesComplete)
         output << "  " << sps << "% Jacobian  " << endl << "    if jacobian_eval" << endl;
       else
-        if (simulation_type == SOLVE_BACKWARD_SIMPLE || simulation_type == SOLVE_FORWARD_SIMPLE
-            || simulation_type == SOLVE_BACKWARD_COMPLETE || simulation_type == SOLVE_FORWARD_COMPLETE)
+        if (simulation_type == BlockSimulationType::solveBackwardSimple
+            || simulation_type == BlockSimulationType::solveForwardSimple
+            || simulation_type == BlockSimulationType::solveBackwardComplete
+            || simulation_type == BlockSimulationType::solveForwardComplete)
           output << "  % Jacobian  " << endl << "  if jacobian_eval" << endl;
         else
           output << "    % Jacobian  " << endl << "    if jacobian_eval" << endl;
@@ -744,15 +767,15 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
 
       switch (simulation_type)
         {
-        case EVALUATE_FORWARD:
-        case EVALUATE_BACKWARD:
+        case BlockSimulationType::evaluateForward:
+        case BlockSimulationType::evaluateBackward:
           output << "    end;" << endl
                  << "  end;" << endl;
           break;
-        case SOLVE_BACKWARD_SIMPLE:
-        case SOLVE_FORWARD_SIMPLE:
-        case SOLVE_BACKWARD_COMPLETE:
-        case SOLVE_FORWARD_COMPLETE:
+        case BlockSimulationType::solveBackwardSimple:
+        case BlockSimulationType::solveForwardSimple:
+        case BlockSimulationType::solveBackwardComplete:
+        case BlockSimulationType::solveForwardComplete:
           output << "  else" << endl;
           for (const auto &it : blocks_derivatives[block])
             {
@@ -775,8 +798,8 @@ DynamicModel::writeModelEquationsOrdered_M(const string &basename) const
             }
           output << "  end;" << endl;
           break;
-        case SOLVE_TWO_BOUNDARIES_SIMPLE:
-        case SOLVE_TWO_BOUNDARIES_COMPLETE:
+        case BlockSimulationType::solveTwoBoundariesSimple:
+        case BlockSimulationType::solveTwoBoundariesComplete:
           output << "    else" << endl;
           for (const auto &it : blocks_derivatives[block])
             {
@@ -892,13 +915,13 @@ DynamicModel::writeModelEquationsCode(const string &basename, const map_idx_t &m
   int u_count_int = 0;
   BlockSimulationType simulation_type;
   if ((max_endo_lag > 0) && (max_endo_lead > 0))
-    simulation_type = SOLVE_TWO_BOUNDARIES_COMPLETE;
+    simulation_type = BlockSimulationType::solveTwoBoundariesComplete;
   else if ((max_endo_lag >= 0) && (max_endo_lead == 0))
-    simulation_type = SOLVE_FORWARD_COMPLETE;
+    simulation_type = BlockSimulationType::solveForwardComplete;
   else
-    simulation_type = SOLVE_BACKWARD_COMPLETE;
+    simulation_type = BlockSimulationType::solveBackwardComplete;
 
-  Write_Inf_To_Bin_File(basename + "/model/bytecode/dynamic.bin", u_count_int, file_open, simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE, symbol_table.endo_nbr());
+  Write_Inf_To_Bin_File(basename + "/model/bytecode/dynamic.bin", u_count_int, file_open, simulation_type == BlockSimulationType::solveTwoBoundariesComplete, symbol_table.endo_nbr());
   file_open = true;
 
   //Temporary variables declaration
@@ -1178,11 +1201,13 @@ DynamicModel::writeModelEquationsCode_Block(const string &basename, const map_id
       int block_max_lag = max_leadlag_block[block].first;
       int block_max_lead = max_leadlag_block[block].second;
 
-      if (simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE || simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE
-          || simulation_type == SOLVE_BACKWARD_COMPLETE || simulation_type == SOLVE_FORWARD_COMPLETE)
+      if (simulation_type == BlockSimulationType::solveTwoBoundariesSimple
+          || simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+          || simulation_type == BlockSimulationType::solveBackwardComplete
+          || simulation_type == BlockSimulationType::solveForwardComplete)
         {
           Write_Inf_To_Bin_File_Block(basename, block, u_count_int, file_open,
-                                      simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE, linear_decomposition);
+                                      simulation_type == BlockSimulationType::solveTwoBoundariesComplete || simulation_type == BlockSimulationType::solveTwoBoundariesSimple, linear_decomposition);
           file_open = true;
         }
       map<tuple<int, int, int>, expr_t> tmp_block_endo_derivative;
@@ -1311,8 +1336,8 @@ DynamicModel::writeModelEquationsCode_Block(const string &basename, const map_id
           switch (simulation_type)
             {
             evaluation:
-            case EVALUATE_BACKWARD:
-            case EVALUATE_FORWARD:
+            case BlockSimulationType::evaluateBackward:
+            case BlockSimulationType::evaluateForward:
               equ_type = getBlockEquationType(block, i);
               {
                 FNUMEXPR_ fnumexpr(ModelEquation, getBlockEquationID(block, i));
@@ -1335,10 +1360,10 @@ DynamicModel::writeModelEquationsCode_Block(const string &basename, const map_id
                   lhs->compile(code_file, instruction_number, true, temporary_terms, map_idx, true, false);
                 }
               break;
-            case SOLVE_BACKWARD_COMPLETE:
-            case SOLVE_FORWARD_COMPLETE:
-            case SOLVE_TWO_BOUNDARIES_COMPLETE:
-            case SOLVE_TWO_BOUNDARIES_SIMPLE:
+            case BlockSimulationType::solveBackwardComplete:
+            case BlockSimulationType::solveForwardComplete:
+            case BlockSimulationType::solveTwoBoundariesComplete:
+            case BlockSimulationType::solveTwoBoundariesSimple:
               if (i < static_cast<int>(block_recursive))
                 goto evaluation;
               variable_ID = getBlockVariableID(block, i);
@@ -1371,13 +1396,13 @@ DynamicModel::writeModelEquationsCode_Block(const string &basename, const map_id
       fjmp_if_eval.write(code_file, instruction_number);
       int prev_instruction_number = instruction_number;
       // The Jacobian if we have to solve the block determinsitic block
-      if (simulation_type != EVALUATE_BACKWARD
-          && simulation_type != EVALUATE_FORWARD)
+      if (simulation_type != BlockSimulationType::evaluateBackward
+          && simulation_type != BlockSimulationType::evaluateForward)
         {
           switch (simulation_type)
             {
-            case SOLVE_BACKWARD_SIMPLE:
-            case SOLVE_FORWARD_SIMPLE:
+            case BlockSimulationType::solveBackwardSimple:
+            case BlockSimulationType::solveForwardSimple:
               {
                 FNUMEXPR_ fnumexpr(FirstEndoDerivative, getBlockEquationID(block, 0), getBlockVariableID(block, 0), 0);
                 fnumexpr.write(code_file, instruction_number);
@@ -1389,10 +1414,10 @@ DynamicModel::writeModelEquationsCode_Block(const string &basename, const map_id
               }
               break;
 
-            case SOLVE_BACKWARD_COMPLETE:
-            case SOLVE_FORWARD_COMPLETE:
-            case SOLVE_TWO_BOUNDARIES_COMPLETE:
-            case SOLVE_TWO_BOUNDARIES_SIMPLE:
+            case BlockSimulationType::solveBackwardComplete:
+            case BlockSimulationType::solveForwardComplete:
+            case BlockSimulationType::solveTwoBoundariesComplete:
+            case BlockSimulationType::solveTwoBoundariesSimple:
               count_u = feedback_variables.size();
               for (const auto &it : blocks_derivatives[block])
                 {
@@ -1403,7 +1428,9 @@ DynamicModel::writeModelEquationsCode_Block(const string &basename, const map_id
                   unsigned int varr = getBlockVariableID(block, var);
                   if (eq >= block_recursive and var >= block_recursive)
                     {
-                      if (lag != 0 && (simulation_type == SOLVE_FORWARD_COMPLETE || simulation_type == SOLVE_BACKWARD_COMPLETE))
+                      if (lag != 0
+                          && (simulation_type == BlockSimulationType::solveForwardComplete
+                              || simulation_type == BlockSimulationType::solveBackwardComplete))
                         continue;
                       if (!Uf[eqr].Ufl)
                         {
@@ -1927,7 +1954,8 @@ DynamicModel::writeSparseDynamicMFile(const string &basename) const
         block_recursive = block_size - block_mfs;
       BlockSimulationType simulation_type = getBlockSimulationType(block);
 
-      if (simulation_type == EVALUATE_FORWARD || simulation_type == EVALUATE_BACKWARD)
+      if (simulation_type == BlockSimulationType::evaluateForward
+          || simulation_type == BlockSimulationType::evaluateBackward)
         {
           for (unsigned int ik = 0; ik < block_size; ik++)
             {
@@ -1948,23 +1976,23 @@ DynamicModel::writeSparseDynamicMFile(const string &basename) const
 
       switch (simulation_type)
         {
-        case EVALUATE_FORWARD:
-        case EVALUATE_BACKWARD:
+        case BlockSimulationType::evaluateForward:
+        case BlockSimulationType::evaluateBackward:
           mDynamicModelFile << "    [y, dr(" << count_call << ").g1, dr(" << count_call << ").g2, dr(" << count_call << ").g3, dr(" << count_call << ").g1_x, dr(" << count_call << ").g1_xd, dr(" << count_call << ").g1_o]=" << basename << ".block.dynamic_" << block + 1 << "(y, x, params, steady_state, 1, it_-1, 1);" << endl
                             << "    residual(y_index_eq)=ys(y_index)-y(it_, y_index);" << endl;
           break;
-        case SOLVE_FORWARD_SIMPLE:
-        case SOLVE_BACKWARD_SIMPLE:
+        case BlockSimulationType::solveForwardSimple:
+        case BlockSimulationType::solveBackwardSimple:
           mDynamicModelFile << "    [r, y, dr(" << count_call << ").g1, dr(" << count_call << ").g2, dr(" << count_call << ").g3, dr(" << count_call << ").g1_x, dr(" << count_call << ").g1_xd, dr(" << count_call << ").g1_o]=" << basename << ".block.dynamic_" << block + 1 << "(y, x, params, steady_state, it_, 1);" << endl
                             << "    residual(y_index_eq)=r;" << endl;
           break;
-        case SOLVE_FORWARD_COMPLETE:
-        case SOLVE_BACKWARD_COMPLETE:
+        case BlockSimulationType::solveForwardComplete:
+        case BlockSimulationType::solveBackwardComplete:
           mDynamicModelFile << "    [r, y, dr(" << count_call << ").g1, dr(" << count_call << ").g2, dr(" << count_call << ").g3, dr(" << count_call << ").g1_x, dr(" << count_call << ").g1_xd, dr(" << count_call << ").g1_o]=" << basename << ".block.dynamic_" << block + 1 << "(y, x, params, steady_state, it_, 1);" << endl
                             << "    residual(y_index_eq)=r;" << endl;
           break;
-        case SOLVE_TWO_BOUNDARIES_COMPLETE:
-        case SOLVE_TWO_BOUNDARIES_SIMPLE:
+        case BlockSimulationType::solveTwoBoundariesComplete:
+        case BlockSimulationType::solveTwoBoundariesSimple:
           mDynamicModelFile << "    [r, y, dr(" << count_call << ").g1, dr(" << count_call << ").g2, dr(" << count_call << ").g3, b, dr(" << count_call << ").g1_x, dr(" << count_call << ").g1_xd, dr(" << count_call << ").g1_o]=" << basename << ".block.dynamic_" <<  block + 1 << "(y, x, params, steady_state, it_-" << max_lag << ", 1, " << max_lag << ", " << block_recursive << "," << "options_.periods" << ");" << endl
                             << "    residual(y_index_eq)=r(:,M_.maximum_lag+1);" << endl;
           break;
@@ -2017,7 +2045,7 @@ DynamicModel::writeSparseDynamicMFile(const string &basename) const
       unsigned int block_recursive = block_size - block_mfs;
       BlockSimulationType simulation_type = getBlockSimulationType(block);
 
-      if (simulation_type == EVALUATE_FORWARD && block_size)
+      if (simulation_type == BlockSimulationType::evaluateForward && block_size)
         {
           if (open_par)
             mDynamicModelFile << "  end" << endl;
@@ -2043,7 +2071,7 @@ DynamicModel::writeSparseDynamicMFile(const string &basename) const
                             << "    return;" << endl
                             << "  end;" << endl;
         }
-      else if (simulation_type == EVALUATE_BACKWARD && block_size)
+      else if (simulation_type == BlockSimulationType::evaluateBackward && block_size)
         {
           if (open_par)
             mDynamicModelFile << "  end" << endl;
@@ -2069,7 +2097,8 @@ DynamicModel::writeSparseDynamicMFile(const string &basename) const
                             << "    return;" << endl
                             << "  end;" << endl;
         }
-      else if ((simulation_type == SOLVE_FORWARD_COMPLETE || simulation_type == SOLVE_FORWARD_SIMPLE) && block_size)
+      else if ((simulation_type == BlockSimulationType::solveForwardComplete
+                || simulation_type == BlockSimulationType::solveForwardSimple) && block_size)
         {
           if (open_par)
             mDynamicModelFile << "  end" << endl;
@@ -2099,7 +2128,8 @@ DynamicModel::writeSparseDynamicMFile(const string &basename) const
                             << "    return;" << endl
                             << "  end;" << endl;
         }
-      else if ((simulation_type == SOLVE_BACKWARD_COMPLETE || simulation_type == SOLVE_BACKWARD_SIMPLE) && block_size)
+      else if ((simulation_type == BlockSimulationType::solveBackwardComplete
+                || simulation_type == BlockSimulationType::solveBackwardSimple) && block_size)
         {
           if (open_par)
             mDynamicModelFile << "  end" << endl;
@@ -2130,7 +2160,8 @@ DynamicModel::writeSparseDynamicMFile(const string &basename) const
                             << "    return;" << endl
                             << "  end;" << endl;
         }
-      else if ((simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE) && block_size)
+      else if ((simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+                || simulation_type == BlockSimulationType::solveTwoBoundariesSimple) && block_size)
         {
           if (open_par)
             mDynamicModelFile << "  end" << endl;
@@ -3173,7 +3204,7 @@ DynamicModel::writeOutput(ostream &output, const string &basename, bool block_de
           for (const auto &it : other_endo_block[block])
             for (int it1 : it.second)
               other_endogenous.insert(it1);
-          output << "block_structure.block(" << block+1 << ").Simulation_Type = " << simulation_type << ";" << endl
+          output << "block_structure.block(" << block+1 << ").Simulation_Type = " << static_cast<int>(simulation_type) << ";" << endl
                  << "block_structure.block(" << block+1 << ").maximum_lag = " << max_lag << ";" << endl
                  << "block_structure.block(" << block+1 << ").maximum_lead = " << max_lead << ";" << endl
                  << "block_structure.block(" << block+1 << ").maximum_endo_lag = " << max_lag_endo << ";" << endl
@@ -5015,7 +5046,8 @@ DynamicModel::get_Derivatives(int block)
   int max_lag, max_lead;
   map<tuple<int, int, int, int, int>, int> Derivatives;
   BlockSimulationType simulation_type = getBlockSimulationType(block);
-  if (simulation_type == EVALUATE_BACKWARD || simulation_type == EVALUATE_FORWARD)
+  if (simulation_type == BlockSimulationType::evaluateBackward
+      || simulation_type == BlockSimulationType::evaluateForward)
     {
       max_lag = 1;
       max_lead = 1;
diff --git a/src/ModelTree.cc b/src/ModelTree.cc
index e9b08abb578ac0446a4fc25530112fbece0c0d7c..43624697b129dc1bc5708ffa2074aec614b42e0b 100644
--- a/src/ModelTree.cc
+++ b/src/ModelTree.cc
@@ -1053,7 +1053,9 @@ ModelTree::printBlockDecomposition(const vector<pair<int, int>> &blocks) const
   for (unsigned int block = 0; block < Nb_TotalBlocks; block++)
     {
       BlockSimulationType simulation_type = getBlockSimulationType(block);
-      if (simulation_type == SOLVE_FORWARD_COMPLETE || simulation_type == SOLVE_BACKWARD_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE)
+      if (simulation_type == BlockSimulationType::solveForwardComplete
+          || simulation_type == BlockSimulationType::solveBackwardComplete
+          || simulation_type == BlockSimulationType::solveTwoBoundariesComplete)
         {
           Nb_SimulBlocks++;
           int size = getBlockSize(block);
@@ -1080,7 +1082,7 @@ ModelTree::reduceBlocksAndTypeDetermination(const vector<pair<int, int>> &blocks
   int Blck_Size, MFS_Size;
   int Lead, Lag;
   block_type_firstequation_size_mfs.clear();
-  BlockSimulationType Simulation_Type, prev_Type = UNKNOWN;
+  BlockSimulationType Simulation_Type, prev_Type = BlockSimulationType::unknown;
   int eq = 0;
   unsigned int l_n_static = 0, l_n_forward = 0, l_n_backward = 0, l_n_mixed = 0;
   for (i = 0; i < static_cast<int>(prologue+blocks.size()+epilogue); i++)
@@ -1140,23 +1142,23 @@ ModelTree::reduceBlocksAndTypeDetermination(const vector<pair<int, int>> &blocks
       if (Lag > 0 && Lead > 0)
         {
           if (Blck_Size == 1)
-            Simulation_Type = SOLVE_TWO_BOUNDARIES_SIMPLE;
+            Simulation_Type = BlockSimulationType::solveTwoBoundariesSimple;
           else
-            Simulation_Type = SOLVE_TWO_BOUNDARIES_COMPLETE;
+            Simulation_Type = BlockSimulationType::solveTwoBoundariesComplete;
         }
       else if (Blck_Size > 1)
         {
           if (Lead > 0)
-            Simulation_Type = SOLVE_BACKWARD_COMPLETE;
+            Simulation_Type = BlockSimulationType::solveBackwardComplete;
           else
-            Simulation_Type = SOLVE_FORWARD_COMPLETE;
+            Simulation_Type = BlockSimulationType::solveForwardComplete;
         }
       else
         {
           if (Lead > 0)
-            Simulation_Type = SOLVE_BACKWARD_SIMPLE;
+            Simulation_Type = BlockSimulationType::solveBackwardSimple;
           else
-            Simulation_Type = SOLVE_FORWARD_SIMPLE;
+            Simulation_Type = BlockSimulationType::solveForwardSimple;
         }
       l_n_static = n_static[i];
       l_n_forward = n_forward[i];
@@ -1166,18 +1168,19 @@ ModelTree::reduceBlocksAndTypeDetermination(const vector<pair<int, int>> &blocks
         {
           if (Equation_Type[equation_reordered[eq]].first == E_EVALUATE || Equation_Type[equation_reordered[eq]].first == E_EVALUATE_S)
             {
-              if (Simulation_Type == SOLVE_BACKWARD_SIMPLE)
-                Simulation_Type = EVALUATE_BACKWARD;
-              else if (Simulation_Type == SOLVE_FORWARD_SIMPLE)
-                Simulation_Type = EVALUATE_FORWARD;
+              if (Simulation_Type == BlockSimulationType::solveBackwardSimple)
+                Simulation_Type = BlockSimulationType::evaluateBackward;
+              else if (Simulation_Type == BlockSimulationType::solveForwardSimple)
+                Simulation_Type = BlockSimulationType::evaluateForward;
             }
           if (i > 0)
             {
               bool is_lead = false, is_lag = false;
               int c_Size = get<2>(block_type_firstequation_size_mfs[block_type_firstequation_size_mfs.size()-1]);
               int first_equation = get<1>(block_type_firstequation_size_mfs[block_type_firstequation_size_mfs.size()-1]);
-              if (c_Size > 0 && ((prev_Type == EVALUATE_FORWARD && Simulation_Type == EVALUATE_FORWARD && !is_lead)
-                                 || (prev_Type == EVALUATE_BACKWARD && Simulation_Type == EVALUATE_BACKWARD && !is_lag)))
+              if (c_Size > 0
+                  && ((prev_Type == BlockSimulationType::evaluateForward && Simulation_Type == BlockSimulationType::evaluateForward && !is_lead)
+                      || (prev_Type == BlockSimulationType::evaluateBackward && Simulation_Type == BlockSimulationType::evaluateBackward && !is_lag)))
                 {
                   for (int j = first_equation; j < first_equation+c_Size; j++)
                     {
@@ -1189,8 +1192,8 @@ ModelTree::reduceBlocksAndTypeDetermination(const vector<pair<int, int>> &blocks
                         is_lead = true;
                     }
                 }
-              if ((prev_Type == EVALUATE_FORWARD && Simulation_Type == EVALUATE_FORWARD && !is_lead)
-                  || (prev_Type == EVALUATE_BACKWARD && Simulation_Type == EVALUATE_BACKWARD && !is_lag))
+              if ((prev_Type == BlockSimulationType::evaluateForward && Simulation_Type == BlockSimulationType::evaluateForward && !is_lead)
+                  || (prev_Type == BlockSimulationType::evaluateBackward && Simulation_Type == BlockSimulationType::evaluateBackward && !is_lag))
                 {
                   //merge the current block with the previous one
                   BlockSimulationType c_Type = get<0>(block_type_firstequation_size_mfs[block_type_firstequation_size_mfs.size()-1]);
@@ -1259,7 +1262,8 @@ ModelTree::determineLinearBlocks()
       int block_size = getBlockSize(block);
       block_derivatives_equation_variable_laglead_nodeid_t derivatives_block = blocks_derivatives[block];
       int first_variable_position = getBlockFirstEquation(block);
-      if (simulation_type == SOLVE_BACKWARD_COMPLETE || simulation_type == SOLVE_FORWARD_COMPLETE)
+      if (simulation_type == BlockSimulationType::solveBackwardComplete
+          || simulation_type == BlockSimulationType::solveForwardComplete)
         for (const auto &[ignore, ignore2, lag, d1] : derivatives_block)
           {
             if (lag == 0)
@@ -1275,7 +1279,8 @@ ModelTree::determineLinearBlocks()
                       }
               }
           }
-      else if (simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE)
+      else if (simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+               || simulation_type == BlockSimulationType::solveTwoBoundariesSimple)
         for (const auto &[ignore, ignore2, lag, d1] : derivatives_block)
           {
             set<pair<int, int>> endogenous;
diff --git a/src/ModelTree.hh b/src/ModelTree.hh
index 9203b79d79b3b61b33d99a98b38725b72d827443..cf9a1fe8f6063e339b09d8ab034aac6f26b720fe 100644
--- a/src/ModelTree.hh
+++ b/src/ModelTree.hh
@@ -468,25 +468,25 @@ public:
   };
 
   inline static string
-  BlockSim(int type)
+  BlockSim(BlockSimulationType type)
   {
     switch (type)
       {
-      case EVALUATE_FORWARD:
+      case BlockSimulationType::evaluateForward:
         return "EVALUATE FORWARD             ";
-      case EVALUATE_BACKWARD:
+      case BlockSimulationType::evaluateBackward:
         return "EVALUATE BACKWARD            ";
-      case SOLVE_FORWARD_SIMPLE:
+      case BlockSimulationType::solveForwardSimple:
         return "SOLVE FORWARD SIMPLE         ";
-      case SOLVE_BACKWARD_SIMPLE:
+      case BlockSimulationType::solveBackwardSimple:
         return "SOLVE BACKWARD SIMPLE        ";
-      case SOLVE_TWO_BOUNDARIES_SIMPLE:
+      case BlockSimulationType::solveTwoBoundariesSimple:
         return "SOLVE TWO BOUNDARIES SIMPLE  ";
-      case SOLVE_FORWARD_COMPLETE:
+      case BlockSimulationType::solveForwardComplete:
         return "SOLVE FORWARD COMPLETE       ";
-      case SOLVE_BACKWARD_COMPLETE:
+      case BlockSimulationType::solveBackwardComplete:
         return "SOLVE BACKWARD COMPLETE      ";
-      case SOLVE_TWO_BOUNDARIES_COMPLETE:
+      case BlockSimulationType::solveTwoBoundariesComplete:
         return "SOLVE TWO BOUNDARIES COMPLETE";
       default:
         return "UNKNOWN                      ";
diff --git a/src/StaticModel.cc b/src/StaticModel.cc
index b73030cdf0307305cd25d23721a6de145c8b3453..38693849cccb0a3e9d9e805aa3823c43e071150b 100644
--- a/src/StaticModel.cc
+++ b/src/StaticModel.cc
@@ -305,20 +305,26 @@ StaticModel::writeModelEquationsOrdered_M(const string &basename) const
              << "% Warning : this file is generated automatically by Dynare" << endl
              << "%           from model file (.mod)" << endl << endl
              << "%/" << endl;
-      if (simulation_type == EVALUATE_BACKWARD || simulation_type == EVALUATE_FORWARD)
+      if (simulation_type == BlockSimulationType::evaluateBackward
+          || simulation_type == BlockSimulationType::evaluateForward)
         output << "function y = static_" << block+1 << "(y, x, params)" << endl;
       else
         output << "function [residual, y, g1] = static_" << block+1 << "(y, x, params)" << endl;
 
       BlockType block_type;
-      if (simulation_type == SOLVE_FORWARD_COMPLETE || simulation_type == SOLVE_BACKWARD_COMPLETE)
+      if (simulation_type == BlockSimulationType::solveForwardComplete
+          || simulation_type == BlockSimulationType::solveBackwardComplete)
         block_type = BlockType::simultans;
-      else if ((simulation_type == SOLVE_FORWARD_SIMPLE || simulation_type == SOLVE_BACKWARD_SIMPLE
-                || simulation_type == EVALUATE_BACKWARD || simulation_type == EVALUATE_FORWARD)
+      else if ((simulation_type == BlockSimulationType::solveForwardSimple
+                || simulation_type == BlockSimulationType::solveBackwardSimple
+                || simulation_type == BlockSimulationType::evaluateBackward
+                || simulation_type == BlockSimulationType::evaluateForward)
                && getBlockFirstEquation(block) < prologue)
         block_type = BlockType::prologue;
-      else if ((simulation_type == SOLVE_FORWARD_SIMPLE || simulation_type == SOLVE_BACKWARD_SIMPLE
-                || simulation_type == EVALUATE_BACKWARD || simulation_type == EVALUATE_FORWARD)
+      else if ((simulation_type == BlockSimulationType::solveForwardSimple
+                || simulation_type == BlockSimulationType::solveBackwardSimple
+                || simulation_type == BlockSimulationType::evaluateBackward
+                || simulation_type == BlockSimulationType::evaluateForward)
                && getBlockFirstEquation(block) >= equations.size() - epilogue)
         block_type = BlockType::epilogue;
       else
@@ -331,7 +337,8 @@ StaticModel::writeModelEquationsOrdered_M(const string &basename) const
              << "  % ////////////////////////////////////////////////////////////////////////" << endl;
       output << "  global options_;" << endl;
       //The Temporary terms
-      if (simulation_type != EVALUATE_BACKWARD && simulation_type != EVALUATE_FORWARD)
+      if (simulation_type != BlockSimulationType::evaluateBackward
+          && simulation_type != BlockSimulationType::evaluateForward)
         output << " g1 = spalloc("  << block_mfs << ", " << block_mfs << ", " << derivative_endo[block].size() << ");" << endl;
 
       if (v_temporary_terms_inuse[block].size())
@@ -342,7 +349,8 @@ StaticModel::writeModelEquationsOrdered_M(const string &basename) const
           output << "  global" << tmp_output.str() << ";" << endl;
         }
 
-      if (simulation_type != EVALUATE_BACKWARD && simulation_type != EVALUATE_FORWARD)
+      if (simulation_type != BlockSimulationType::evaluateBackward
+          && simulation_type != BlockSimulationType::evaluateForward)
         output << "  residual=zeros(" << block_mfs << ",1);" << endl;
 
       // The equations
@@ -381,8 +389,8 @@ StaticModel::writeModelEquationsOrdered_M(const string &basename) const
           lhs->writeOutput(tmp_output, local_output_type, local_temporary_terms, {});
           switch (simulation_type)
             {
-            case EVALUATE_BACKWARD:
-            case EVALUATE_FORWARD:
+            case BlockSimulationType::evaluateBackward:
+            case BlockSimulationType::evaluateForward:
             evaluation:
               output << "  % equation " << getBlockEquationID(block, i)+1 << " variable : " << sModel
                      << " (" << variable_ID+1 << ") " << c_Equation_Type(equ_type) << endl;
@@ -417,10 +425,10 @@ StaticModel::writeModelEquationsOrdered_M(const string &basename) const
                 }
               output << ";" << endl;
               break;
-            case SOLVE_BACKWARD_SIMPLE:
-            case SOLVE_FORWARD_SIMPLE:
-            case SOLVE_BACKWARD_COMPLETE:
-            case SOLVE_FORWARD_COMPLETE:
+            case BlockSimulationType::solveBackwardSimple:
+            case BlockSimulationType::solveForwardSimple:
+            case BlockSimulationType::solveBackwardComplete:
+            case BlockSimulationType::solveForwardComplete:
               if (i < block_recursive)
                 goto evaluation;
               feedback_variables.push_back(variable_ID);
@@ -437,15 +445,17 @@ StaticModel::writeModelEquationsOrdered_M(const string &basename) const
             }
         }
       // The Jacobian if we have to solve the block
-      if (simulation_type == SOLVE_BACKWARD_SIMPLE || simulation_type == SOLVE_FORWARD_SIMPLE
-          || simulation_type == SOLVE_BACKWARD_COMPLETE || simulation_type == SOLVE_FORWARD_COMPLETE)
+      if (simulation_type == BlockSimulationType::solveBackwardSimple
+          || simulation_type == BlockSimulationType::solveForwardSimple
+          || simulation_type == BlockSimulationType::solveBackwardComplete
+          || simulation_type == BlockSimulationType::solveForwardComplete)
         output << "  " << sps << "% Jacobian  " << endl;
       switch (simulation_type)
         {
-        case SOLVE_BACKWARD_SIMPLE:
-        case SOLVE_FORWARD_SIMPLE:
-        case SOLVE_BACKWARD_COMPLETE:
-        case SOLVE_FORWARD_COMPLETE:
+        case BlockSimulationType::solveBackwardSimple:
+        case BlockSimulationType::solveForwardSimple:
+        case BlockSimulationType::solveBackwardComplete:
+        case BlockSimulationType::solveForwardComplete:
           for (const auto &it : blocks_derivatives[block])
             {
               unsigned int eq, var;
@@ -497,7 +507,7 @@ StaticModel::writeModelEquationsCode(const string &basename, map_idx_t map_idx)
   FDIMST_ fdimst(temporary_terms.size());
   fdimst.write(code_file, instruction_number);
   FBEGINBLOCK_ fbeginblock(symbol_table.endo_nbr(),
-                           SOLVE_FORWARD_COMPLETE,
+                           BlockSimulationType::solveForwardComplete,
                            0,
                            symbol_table.endo_nbr(),
                            variable_reordered,
@@ -681,8 +691,10 @@ StaticModel::writeModelEquationsCode_Block(const string &basename, map_idx_t map
       unsigned int block_mfs = getBlockMfs(block);
       unsigned int block_recursive = block_size - block_mfs;
 
-      if (simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE || simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE
-          || simulation_type == SOLVE_BACKWARD_COMPLETE || simulation_type == SOLVE_FORWARD_COMPLETE)
+      if (simulation_type == BlockSimulationType::solveTwoBoundariesSimple
+          || simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+          || simulation_type == BlockSimulationType::solveBackwardComplete
+          || simulation_type == BlockSimulationType::solveForwardComplete)
         {
           Write_Inf_To_Bin_File_Block(basename, block, u_count_int, file_open);
           file_open = true;
@@ -736,8 +748,8 @@ StaticModel::writeModelEquationsCode_Block(const string &basename, map_idx_t map
           switch (simulation_type)
             {
             evaluation:
-            case EVALUATE_BACKWARD:
-            case EVALUATE_FORWARD:
+            case BlockSimulationType::evaluateBackward:
+            case BlockSimulationType::evaluateForward:
               equ_type = getBlockEquationType(block, i);
               {
                 FNUMEXPR_ fnumexpr(ModelEquation, getBlockEquationID(block, i));
@@ -760,8 +772,8 @@ StaticModel::writeModelEquationsCode_Block(const string &basename, map_idx_t map
                   lhs->compile(code_file, instruction_number, true, temporary_terms, map_idx, false, false);
                 }
               break;
-            case SOLVE_BACKWARD_COMPLETE:
-            case SOLVE_FORWARD_COMPLETE:
+            case BlockSimulationType::solveBackwardComplete:
+            case BlockSimulationType::solveForwardComplete:
               if (i < static_cast<int>(block_recursive))
                 goto evaluation;
               variable_ID = getBlockVariableID(block, i);
@@ -790,13 +802,13 @@ StaticModel::writeModelEquationsCode_Block(const string &basename, map_idx_t map
       fendequ.write(code_file, instruction_number);
 
       // The Jacobian if we have to solve the block
-      if (simulation_type != EVALUATE_BACKWARD
-          && simulation_type != EVALUATE_FORWARD)
+      if (simulation_type != BlockSimulationType::evaluateBackward
+          && simulation_type != BlockSimulationType::evaluateForward)
         {
           switch (simulation_type)
             {
-            case SOLVE_BACKWARD_SIMPLE:
-            case SOLVE_FORWARD_SIMPLE:
+            case BlockSimulationType::solveBackwardSimple:
+            case BlockSimulationType::solveForwardSimple:
               {
                 FNUMEXPR_ fnumexpr(FirstEndoDerivative, 0, 0);
                 fnumexpr.write(code_file, instruction_number);
@@ -808,8 +820,8 @@ StaticModel::writeModelEquationsCode_Block(const string &basename, map_idx_t map
               }
               break;
 
-            case SOLVE_BACKWARD_COMPLETE:
-            case SOLVE_FORWARD_COMPLETE:
+            case BlockSimulationType::solveBackwardComplete:
+            case BlockSimulationType::solveForwardComplete:
               count_u = feedback_variables.size();
               for (const auto &it : blocks_derivatives[block])
                 {
@@ -928,8 +940,8 @@ StaticModel::writeModelEquationsCode_Block(const string &basename, map_idx_t map
           switch (simulation_type)
             {
             evaluation_l:
-            case EVALUATE_BACKWARD:
-            case EVALUATE_FORWARD:
+            case BlockSimulationType::evaluateBackward:
+            case BlockSimulationType::evaluateForward:
               equ_type = getBlockEquationType(block, i);
               {
                 FNUMEXPR_ fnumexpr(ModelEquation, getBlockEquationID(block, i));
@@ -952,8 +964,8 @@ StaticModel::writeModelEquationsCode_Block(const string &basename, map_idx_t map
                   lhs->compile(code_file, instruction_number, true, tt2, map_idx2[block], false, false);
                 }
               break;
-            case SOLVE_BACKWARD_COMPLETE:
-            case SOLVE_FORWARD_COMPLETE:
+            case BlockSimulationType::solveBackwardComplete:
+            case BlockSimulationType::solveForwardComplete:
               if (i < static_cast<int>(block_recursive))
                 goto evaluation_l;
               variable_ID = getBlockVariableID(block, i);
@@ -984,8 +996,8 @@ StaticModel::writeModelEquationsCode_Block(const string &basename, map_idx_t map
       // The Jacobian if we have to solve the block determinsitic bloc
       switch (simulation_type)
         {
-        case SOLVE_BACKWARD_SIMPLE:
-        case SOLVE_FORWARD_SIMPLE:
+        case BlockSimulationType::solveBackwardSimple:
+        case BlockSimulationType::solveForwardSimple:
           {
             FNUMEXPR_ fnumexpr(FirstEndoDerivative, 0, 0);
             fnumexpr.write(code_file, instruction_number);
@@ -996,10 +1008,10 @@ StaticModel::writeModelEquationsCode_Block(const string &basename, map_idx_t map
             fstpg2.write(code_file, instruction_number);
           }
           break;
-        case EVALUATE_BACKWARD:
-        case EVALUATE_FORWARD:
-        case SOLVE_BACKWARD_COMPLETE:
-        case SOLVE_FORWARD_COMPLETE:
+        case BlockSimulationType::evaluateBackward:
+        case BlockSimulationType::evaluateForward:
+        case BlockSimulationType::solveBackwardComplete:
+        case BlockSimulationType::solveForwardComplete:
           count_u = feedback_variables.size();
           for (const auto &it : blocks_derivatives[block])
             {
@@ -1994,7 +2006,8 @@ StaticModel::writeStaticBlockMFSFile(const string &basename) const
 
       BlockSimulationType simulation_type = getBlockSimulationType(b);
 
-      if (simulation_type == EVALUATE_BACKWARD || simulation_type == EVALUATE_FORWARD)
+      if (simulation_type == BlockSimulationType::evaluateBackward
+          || simulation_type == BlockSimulationType::evaluateForward)
         {
           output << "      y_tmp = " << basename << ".block.static_" << b+1 << "(y, x, params);" << endl;
           ostringstream tmp;
@@ -2037,7 +2050,7 @@ StaticModel::writeOutput(ostream &output, bool block) const
           tmp_s << " " << getBlockVariableID(b, i)+1;
           tmp_s_eq << " " << getBlockEquationID(b, i)+1;
         }
-      output << "block_structure_stat.block(" << b+1 << ").Simulation_Type = " << simulation_type << ";" << endl
+      output << "block_structure_stat.block(" << b+1 << ").Simulation_Type = " << static_cast<int>(simulation_type) << ";" << endl
              << "block_structure_stat.block(" << b+1 << ").endo_nbr = " << block_size << ";" << endl
              << "block_structure_stat.block(" << b+1 << ").mfs = " << getBlockMfs(block) << ";" << endl
              << "block_structure_stat.block(" << b+1 << ").equation = [" << tmp_s_eq.str() << "];" << endl
@@ -2181,7 +2194,8 @@ StaticModel::computeChainRuleJacobian()
       int block_size = getBlockSize(block);
       int block_nb_mfs = getBlockMfs(block);
       int block_nb_recursives = block_size - block_nb_mfs;
-      if (simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE)
+      if (simulation_type == BlockSimulationType::solveTwoBoundariesComplete
+          || simulation_type == BlockSimulationType::solveTwoBoundariesSimple)
         {
           for (int i = 0; i < block_nb_recursives; i++)
             {