diff --git a/.gitignore b/.gitignore
index b25c15b81fae06e1c55946ac6270bfdb293870e8..d83979c63c40955eeb6e390f928cc25a18561dfe 100644
--- a/.gitignore
+++ b/.gitignore
@@ -1 +1,3 @@
 *~
+.DS_Store
+.vscode/
diff --git a/src/DynareMain.cc b/src/DynareMain.cc
index 4e9ff55aef86d93e355572b45ade8c5b97025e1e..23300aaa5b654b5100e9890a24a4ec66e176ce72 100644
--- a/src/DynareMain.cc
+++ b/src/DynareMain.cc
@@ -34,6 +34,7 @@
 #include "Configuration.hh"
 #include "ExtendedPreprocessorTypes.hh"
 #include "ModFile.hh"
+#include "OutputWriter.hh"
 #include "ParsingDriver.hh"
 
 /* Prototype for the function that handles the macro-expansion of the .mod file
@@ -56,7 +57,7 @@ usage()
           "[conffile=path_to_config_file] [parallel_follower_open_mode] "
           "[parallel_test] [parallel_use_psexec=true|false]"
        << " [-D<variable>[=<value>]] [-I/path] [nostrict] [stochastic] [fast] [minimal_workspace] "
-          "[compute_xrefs] [output=second|third] [language=matlab|julia]"
+          "[compute_xrefs] [output=second|third] [language=matlab|julia|python]"
        << " [params_derivs_order=0|1|2] [transform_unary_ops] "
           "[exclude_eqs=<equation_tag_list_or_file>] [include_eqs=<equation_tag_list_or_file>]"
        << " [json=parse|check|transform|compute] [jsonstdout] [onlyjson] [jsonderivsimple] "
@@ -344,6 +345,8 @@ main(int argc, char** argv)
             language = LanguageOutputType::matlab;
           else if (s == "julia")
             language = LanguageOutputType::julia;
+          else if (s == "python")
+            language = LanguageOutputType::python;
           else
             {
               cerr << "Incorrect syntax for language option" << endl;
@@ -533,15 +536,34 @@ main(int argc, char** argv)
   // Do computations
   mod_file->computingPass(no_tmp_terms, output_mode, params_derivs_order);
   if (json == JsonOutputPointType::computingpass)
-    mod_file->writeJsonOutput(basename, json, json_output_mode, onlyjson, jsonderivsimple);
-
-  // Write output files
-  if (language == LanguageOutputType::julia)
-    mod_file->writeJuliaOutput(basename);
-  else
-    mod_file->writeMOutput(basename, clear_all, clear_global, no_warn, console, nograph,
-                           nointeractive, config, check_model_changes, minimal_workspace,
-                           compute_xrefs, mexext, matlabroot, onlymodel, gui, notime);
+    {
+      unique_ptr<OutputWriter> jsonWriter = make_unique<JsonOutputWriter>(
+          *mod_file, basename, json, json_output_mode, onlyjson, jsonderivsimple);
+      jsonWriter->writeOutput();
+    }
+
+  // Write output files using appropriate OutputWriter
+  unique_ptr<OutputWriter> outputWriter;
+
+  switch (language)
+    {
+    case LanguageOutputType::matlab:
+      outputWriter = make_unique<MatlabOutputWriter>(
+          *mod_file, basename, clear_all, clear_global, no_warn, console, nograph, nointeractive,
+          config, check_model_changes, minimal_workspace, compute_xrefs, mexext, matlabroot,
+          onlymodel, gui, notime);
+      break;
+
+    case LanguageOutputType::julia:
+      outputWriter = make_unique<JuliaOutputWriter>(*mod_file, basename);
+      break;
+
+    case LanguageOutputType::python:
+      outputWriter = make_unique<PythonOutputWriter>(*mod_file, basename);
+      break;
+    }
+
+  outputWriter->writeOutput();
 
   /* Ensures that workers are not destroyed before they finish compiling.
      Also ensures that the preprocessor final message is printed after the end of
diff --git a/src/ExtendedPreprocessorTypes.hh b/src/ExtendedPreprocessorTypes.hh
index 79f9e9a6581e7fbf9bf532d18bd14d72d1ee099d..3dc87f9ee3ce51817e8408ca7251475447813db6 100644
--- a/src/ExtendedPreprocessorTypes.hh
+++ b/src/ExtendedPreprocessorTypes.hh
@@ -33,6 +33,7 @@ enum class LanguageOutputType
 {
   matlab, // outputs files for MATLAB/Octave processing
   julia,  // outputs files for Julia
+  python, // outputs files for Python
 };
 
 enum class JsonFileOutputType
diff --git a/src/ModFile.cc b/src/ModFile.cc
index 928ec611b51b1c724c52a1c7c22c13fc15d49dc1..9696baec927a003ccc38debad8c6cfbeafe434a9 100644
--- a/src/ModFile.cc
+++ b/src/ModFile.cc
@@ -282,26 +282,6 @@ ModFile::checkPass(bool nostrict, bool stochastic)
         }
     }
 
-  if (mod_file_struct.dsge_prior_weight_in_estimated_params)
-    {
-      if (!mod_file_struct.dsge_var_estimated && !mod_file_struct.dsge_var_calibrated.empty())
-        {
-          cerr << "ERROR: If dsge_prior_weight is in the estimated_params block, the prior weight "
-                  "cannot be calibrated "
-               << "via the dsge_var option in the estimation statement." << endl;
-          exit(EXIT_FAILURE);
-        }
-      else if (!mod_file_struct.dsge_var_estimated && !symbol_table.exists("dsge_prior_weight"))
-        {
-          cerr << "ERROR: If dsge_prior_weight is in the estimated_params block, it must either be "
-                  "declared as a parameter "
-               << "(deprecated) or the dsge_var option must be passed to the estimation statement "
-                  "(preferred)."
-               << endl;
-          exit(EXIT_FAILURE);
-        }
-    }
-
   /* This check must come before checking that there are as many static-only as dynamic-only
      equations, see #103 */
   dynamic_model.checkOccbinRegimes();
@@ -1005,679 +985,6 @@ ModFile::remove_directory_with_matlab_lock(const filesystem::path& dir)
   remove_all(tmp);
 }
 
-void
-ModFile::writeMOutput(const string& basename, bool clear_all, bool clear_global, bool no_warn,
-                      bool console, bool nograph, bool nointeractive, const Configuration& config,
-                      bool check_model_changes, bool minimal_workspace, bool compute_xrefs,
-                      const string& mexext, const filesystem::path& matlabroot, bool onlymodel,
-                      bool gui, bool notime) const
-{
-  if (basename.empty())
-    {
-      cerr << "ERROR: Missing file name" << endl;
-      exit(EXIT_FAILURE);
-    }
-
-  auto plusfolder {DataTree::packageDir(basename)};
-
-  if (check_model_changes && !heterogeneity_table.empty())
-    {
-      cerr << "ERROR: the 'fast' option is not supported for heterogeneous models" << endl;
-      exit(EXIT_FAILURE);
-    }
-  bool hasModelChanged = !dynamic_model.isChecksumMatching(basename) || !check_model_changes;
-  if (hasModelChanged)
-    {
-      // Erase possible remnants of previous runs
-
-      /* Under MATLAB+Windows (but not under Octave nor under GNU/Linux or
-         macOS), if we directly remove the "+" subdirectory, then the
-         preprocessor is not able to recreate it afterwards (presumably because
-         MATLAB maintains some sort of lock on it). So we use a hack. */
-      remove_directory_with_matlab_lock(plusfolder);
-
-      filesystem::remove_all(basename + "/model/src");
-      filesystem::remove_all(basename + "/model/bytecode");
-      // Do not remove basename/model/julia/, otherwise it would break calls to
-      // writeToFileIfModified()
-    }
-
-  create_directory(plusfolder);
-  filesystem::path fname {plusfolder / "driver.m"};
-  ofstream mOutputFile {fname, ios::out | ios::binary};
-  if (!mOutputFile.is_open())
-    {
-      cerr << "ERROR: Can't open file " << fname.string() << " for writing" << endl;
-      exit(EXIT_FAILURE);
-    }
-
-  mOutputFile << "%" << endl
-              << "% Status : main Dynare file" << endl
-              << "%" << endl
-              << "% Warning : this file is generated automatically by Dynare" << endl
-              << "%           from model file (.mod)" << endl
-              << endl;
-
-  if (no_warn)
-    mOutputFile << "warning off" << endl; // This will be executed *after* function warning_config()
-
-  if (clear_all)
-    mOutputFile << "clearvars -global" << endl
-                << "clear_persistent_variables(fileparts(which('dynare')), false)" << endl;
-  else if (clear_global)
-    mOutputFile << "clearvars -global M_ options_ oo_ estim_params_ bayestopt_ dataset_ "
-                   "dataset_info estimation_info;"
-                << endl;
-
-  if (!notime)
-    mOutputFile << "tic0 = tic;" << endl;
-
-  mOutputFile
-      << "% Define global variables." << endl
-      << "global M_ options_ oo_ estim_params_ bayestopt_ dataset_ dataset_info estimation_info"
-      << endl
-      << "options_ = [];" << endl
-      << "M_.fname = '" << basename << "';" << endl
-      << "M_.dynare_version = '" << PACKAGE_VERSION << "';" << endl
-      << "oo_.dynare_version = '" << PACKAGE_VERSION << "';" << endl
-      << "options_.dynare_version = '" << PACKAGE_VERSION << "';" << endl
-      << "%" << endl
-      << "% Some global variables initialization" << endl
-      << "%" << endl;
-  if (!onlymodel)
-    config.writeHooks(mOutputFile);
-  mOutputFile << "global_initialization;" << endl;
-
-  if (console)
-    mOutputFile << "options_.console_mode = true;" << endl << "options_.nodisplay = true;" << endl;
-  if (nograph)
-    mOutputFile << "options_.nograph = true;" << endl;
-
-  if (nointeractive)
-    mOutputFile << "options_.nointeractive = true;" << endl;
-
-  if (param_used_with_lead_lag)
-    mOutputFile << "M_.parameter_used_with_lead_lag = true;" << endl;
-
-  symbol_table.writeOutput(mOutputFile);
-  heterogeneity_table.writeOutput(mOutputFile);
-
-  var_model_table.writeOutput(basename, mOutputFile);
-  trend_component_model_table.writeOutput(basename, mOutputFile);
-  var_expectation_model_table.writeOutput(mOutputFile);
-  pac_model_table.writeOutput(mOutputFile);
-
-  // Initialize M_.Sigma_e, M_.Correlation_matrix, M_.H, and M_.Correlation_matrix_ME
-  mOutputFile << "M_.Sigma_e = zeros(" << symbol_table.exo_nbr() << ", " << symbol_table.exo_nbr()
-              << ");" << endl
-              << "M_.Correlation_matrix = eye(" << symbol_table.exo_nbr() << ", "
-              << symbol_table.exo_nbr() << ");" << endl;
-  for (int hd {0}; hd < heterogeneity_table.size(); hd++)
-    mOutputFile << "M_.heterogeneity(" << hd + 1 << ").Sigma_e = zeros("
-                << symbol_table.het_exo_nbr(hd) << ", " << symbol_table.het_exo_nbr(hd) << ");"
-                << endl;
-
-  if (mod_file_struct.calibrated_measurement_errors)
-    mOutputFile << "M_.H = zeros(" << symbol_table.observedVariablesNbr() << ", "
-                << symbol_table.observedVariablesNbr() << ");" << endl
-                << "M_.Correlation_matrix_ME = eye(" << symbol_table.observedVariablesNbr() << ", "
-                << symbol_table.observedVariablesNbr() << ");" << endl;
-  else
-    mOutputFile << "M_.H = 0;" << endl << "M_.Correlation_matrix_ME = 1;" << endl;
-
-  // May be later modified by a shocks block
-  mOutputFile << "M_.sigma_e_is_diagonal = true;" << endl;
-
-  /* Initialize the structures created for several blocks, as part of the implementation of the
-     “overwrite” option */
-  mOutputFile << "M_.det_shocks = [];" << endl
-              << "M_.surprise_shocks = [];" << endl
-              << "M_.learnt_shocks = [];" << endl
-              << "M_.learnt_endval = [];" << endl
-              << "M_.heteroskedastic_shocks.Qvalue_orig = [];" << endl
-              << "M_.heteroskedastic_shocks.Qscale_orig = [];" << endl
-              << "M_.matched_irfs = {};" << endl
-              << "M_.matched_irfs_weights = {};" << endl
-              << "M_.perfect_foresight_controlled_paths = [];" << endl;
-
-  // NB: options_.{ramsey,discretionary}_policy should rather be fields of M_
-  mOutputFile << boolalpha << "options_.linear = " << linear << ";" << endl
-              << "options_.block = " << block << ";" << endl
-              << "options_.bytecode = " << bytecode << ";" << endl
-              << "options_.use_dll = " << use_dll << ";" << endl
-              << "options_.ramsey_policy = " << mod_file_struct.ramsey_model_present << ";" << endl
-              << "options_.discretionary_policy = " << mod_file_struct.discretionary_policy_present
-              << ";" << endl;
-
-  if (mod_file_struct.discretionary_policy_present)
-    mOutputFile << "M_.discretionary_orig_eq_nbr = " << original_model.equation_number() << ";"
-                << endl;
-
-  if (parallel_local_files.size() > 0)
-    {
-      mOutputFile << "options_.parallel_info.local_files = {" << endl;
-      for (const auto& parallel_local_file : parallel_local_files)
-        {
-          size_t j = parallel_local_file.find_last_of(R"(/\)");
-          if (j == string::npos)
-            mOutputFile << "'', '" << parallel_local_file << "';" << endl;
-          else
-            mOutputFile << "'" << parallel_local_file.substr(0, j + 1) << "', '"
-                        << parallel_local_file.substr(j + 1) << "';" << endl;
-        }
-      mOutputFile << "};" << endl;
-    }
-
-  if (dynamic_model.isHessianComputed())
-    {
-      mOutputFile << "M_.nonzero_hessian_eqs = ";
-      dynamic_model.printNonZeroHessianEquations(mOutputFile);
-      mOutputFile << ";" << endl << "M_.hessian_eq_zero = isempty(M_.nonzero_hessian_eqs);" << endl;
-    }
-
-  if (!onlymodel)
-    config.writeCluster(mOutputFile);
-
-  if (bytecode)
-    mOutputFile << "if exist('bytecode') ~= 3" << endl
-                << "  error('DYNARE: Can''t find bytecode DLL. Please compile it or remove the "
-                   "''bytecode'' option.')"
-                << endl
-                << "end" << endl;
-
-  mOutputFile << "M_.eq_nbr = " << dynamic_model.equation_number() << ";" << endl
-              << "M_.ramsey_orig_eq_nbr = " << mod_file_struct.ramsey_orig_eq_nbr << ";" << endl
-              << "M_.ramsey_orig_endo_nbr = " << mod_file_struct.ramsey_orig_endo_nbr << ";" << endl
-              << "M_.set_auxiliary_variables = exist(['./+' M_.fname "
-                 "'/set_auxiliary_variables.m'], 'file') == 2;"
-              << endl;
-
-  epilogue.writeOutput(mOutputFile);
-
-  if (dynamic_model.equation_number() > 0)
-    {
-      dynamic_model.writeDriverOutput(mOutputFile, compute_xrefs);
-      if (!no_static)
-        {
-          static_model.writeDriverOutput(mOutputFile);
-          if (mod_file_struct.ramsey_model_present)
-            static_model.writeDriverRamseyMultipliersDerivativesSparseIndices(mOutputFile);
-        }
-    }
-
-  for (const auto& hm : heterogeneous_models)
-    hm.writeDriverOutput(mOutputFile);
-
-  if (onlymodel || gui)
-    for (const auto& statement : statements)
-      {
-        /* Special treatment for initval block: insert initial values for the
-           auxiliary variables and initialize exo det */
-        if (auto ivs = dynamic_cast<InitValStatement*>(statement.get()); ivs)
-          {
-            ivs->writeOutput(mOutputFile, basename, minimal_workspace);
-            static_model.writeAuxVarInitval(mOutputFile, ExprNodeOutputType::matlabOutsideModel);
-            ivs->writeOutputPostInit(mOutputFile);
-          }
-
-        // Special treatment for endval block: insert initial values for the auxiliary variables
-        if (auto evs = dynamic_cast<EndValStatement*>(statement.get()); evs)
-          {
-            evs->writeOutput(mOutputFile, basename, minimal_workspace);
-            static_model.writeAuxVarInitval(mOutputFile, ExprNodeOutputType::matlabOutsideModel);
-          }
-
-        if (auto ips = dynamic_cast<InitParamStatement*>(statement.get()); ips)
-          ips->writeOutput(mOutputFile, basename, minimal_workspace);
-
-        if (auto ss = dynamic_cast<ShocksStatement*>(statement.get()); ss)
-          ss->writeOutput(mOutputFile, basename, minimal_workspace);
-
-        if (auto eps = dynamic_cast<EstimatedParamsStatement*>(statement.get()); eps)
-          eps->writeOutput(mOutputFile, basename, minimal_workspace);
-
-        if (auto sgs = dynamic_cast<ShockGroupsStatement*>(statement.get()); sgs)
-          sgs->writeOutput(mOutputFile, basename, minimal_workspace);
-
-        if (gui)
-          if (auto it = dynamic_cast<NativeStatement*>(statement.get()); it)
-            it->writeOutput(mOutputFile, basename, minimal_workspace);
-      }
-  else
-    {
-      for (const auto& statement : statements)
-        {
-          statement->writeOutput(mOutputFile, basename, minimal_workspace);
-
-          /* Special treatment for initval block: insert initial values for the
-             auxiliary variables and initialize exo det */
-          if (auto ivs = dynamic_cast<InitValStatement*>(statement.get()); ivs)
-            {
-              static_model.writeAuxVarInitval(mOutputFile, ExprNodeOutputType::matlabOutsideModel);
-              ivs->writeOutputPostInit(mOutputFile);
-            }
-
-          // Special treatment for endval block: insert initial values for the auxiliary variables
-          if (auto evs = dynamic_cast<EndValStatement*>(statement.get()); evs)
-            static_model.writeAuxVarInitval(mOutputFile, ExprNodeOutputType::matlabOutsideModel);
-
-          // Special treatment for load params and steady state statement: insert initial values for
-          // the auxiliary variables
-          if (auto lpass = dynamic_cast<LoadParamsAndSteadyStateStatement*>(statement.get());
-              lpass && !no_static)
-            static_model.writeAuxVarInitval(mOutputFile, ExprNodeOutputType::matlabOutsideModel);
-        }
-
-      if (!notime)
-        mOutputFile << endl
-                    << endl
-                    << "oo_.time = toc(tic0);" << endl
-                    << "disp(['Total computing time : ' dynsec2hms(oo_.time) ]);" << endl;
-
-      mOutputFile << "if ~exist([M_.dname filesep 'Output'],'dir')" << endl
-                  << "    mkdir(M_.dname,'Output');" << endl
-                  << "end" << endl
-                  << "save([M_.dname filesep 'Output' filesep '" << basename
-                  << "_results.mat'], 'oo_', 'M_', 'options_');" << endl
-                  << "if exist('estim_params_', 'var') == 1" << endl
-                  << "  save([M_.dname filesep 'Output' filesep '" << basename
-                  << "_results.mat'], 'estim_params_', '-append');" << endl
-                  << "end" << endl
-                  << "if exist('bayestopt_', 'var') == 1" << endl
-                  << "  save([M_.dname filesep 'Output' filesep '" << basename
-                  << "_results.mat'], 'bayestopt_', '-append');" << endl
-                  << "end" << endl
-                  << "if exist('dataset_', 'var') == 1" << endl
-                  << "  save([M_.dname filesep 'Output' filesep '" << basename
-                  << "_results.mat'], 'dataset_', '-append');" << endl
-                  << "end" << endl
-                  << "if exist('estimation_info', 'var') == 1" << endl
-                  << "  save([M_.dname filesep 'Output' filesep '" << basename
-                  << "_results.mat'], 'estimation_info', '-append');" << endl
-                  << "end" << endl
-                  << "if exist('dataset_info', 'var') == 1" << endl
-                  << "  save([M_.dname filesep 'Output' filesep '" << basename
-                  << "_results.mat'], 'dataset_info', '-append');" << endl
-                  << "end" << endl
-                  << "if exist('oo_recursive_', 'var') == 1" << endl
-                  << "  save([M_.dname filesep 'Output' filesep '" << basename
-                  << "_results.mat'], 'oo_recursive_', '-append');" << endl
-                  << "end" << endl
-                  << "if exist('options_mom_', 'var') == 1" << endl
-                  << "  save([M_.dname filesep 'Output' filesep '" << basename
-                  << "_results.mat'], 'options_mom_', '-append');" << endl
-                  << "end" << endl;
-
-      config.writeEndParallel(mOutputFile);
-
-      if (!no_warn)
-        {
-          if (int num_warnings {warnings.numWarnings()}; num_warnings > 0)
-            mOutputFile << "disp('Note: " << num_warnings
-                        << " warning(s) encountered in the preprocessor')" << endl;
-
-          mOutputFile << "if ~isempty(lastwarn)" << endl
-                      << "  disp('Note: warning(s) encountered in MATLAB/Octave code')" << endl
-                      << "end" << endl;
-        }
-    }
-
-  mOutputFile.close();
-
-  if (hasModelChanged)
-    {
-      // Create static and dynamic files
-      if (dynamic_model.equation_number() > 0)
-        {
-          if (!no_static)
-            {
-              static_model.writeStaticFile(basename, use_dll, mexext, matlabroot, false);
-              static_model.writeParamsDerivativesFile<false>(basename);
-              if (mod_file_struct.ramsey_model_present)
-                {
-                  if (use_dll)
-                    static_model.writeRamseyMultipliersDerivativesCFile(
-                        basename, mexext, matlabroot, mod_file_struct.ramsey_orig_endo_nbr);
-                  else
-                    static_model.writeRamseyMultipliersDerivativesMFile(
-                        basename, mod_file_struct.ramsey_orig_endo_nbr);
-                }
-            }
-
-          dynamic_model.writeDynamicFile(basename, use_dll, mexext, matlabroot, false);
-
-          dynamic_model.writeParamsDerivativesFile<false>(basename);
-        }
-
-      // Create steady state file
-      steady_state_model.writeSteadyStateFile(basename, false);
-
-      // Create epilogue file
-      epilogue.writeEpilogueFile(basename);
-
-      pac_model_table.writeTargetCoefficientsFile(basename);
-
-      for (const auto& hm : heterogeneous_models)
-        hm.writeModelFiles(basename, false);
-    }
-}
-
-void
-ModFile::writeJuliaOutput(const string& basename) const
-{
-  if (dynamic_model.equation_number() > 0)
-    {
-      if (!no_static)
-        {
-          static_model.writeStaticFile(basename, false, "", {}, true);
-          static_model.writeParamsDerivativesFile<true>(basename);
-        }
-      dynamic_model.writeDynamicFile(basename, use_dll, "", {}, true);
-      dynamic_model.writeParamsDerivativesFile<true>(basename);
-    }
-  steady_state_model.writeSteadyStateFile(basename, true);
-}
-
-void
-ModFile::writeJsonOutput(const string& basename, JsonOutputPointType json,
-                         JsonFileOutputType json_output_mode, bool onlyjson, bool jsonderivsimple)
-{
-  if (json == JsonOutputPointType::nojson)
-    return;
-
-  if (json == JsonOutputPointType::parsing || json == JsonOutputPointType::checkpass)
-    symbol_table.freeze();
-
-  if (json_output_mode == JsonFileOutputType::standardout)
-    cout << "//-- BEGIN JSON --// " << endl << "{" << endl;
-
-  writeJsonOutputParsingCheck(basename, json_output_mode,
-                              json == JsonOutputPointType::transformpass,
-                              json == JsonOutputPointType::computingpass);
-
-  if (json == JsonOutputPointType::parsing || json == JsonOutputPointType::checkpass)
-    symbol_table.unfreeze();
-
-  if (json == JsonOutputPointType::computingpass)
-    writeJsonComputingPassOutput(basename, json_output_mode, jsonderivsimple);
-
-  if (json_output_mode == JsonFileOutputType::standardout)
-    cout << "}" << endl << "//-- END JSON --// " << endl;
-
-  switch (json)
-    {
-    case JsonOutputPointType::parsing:
-      cout << "JSON written after Parsing step." << endl;
-      break;
-    case JsonOutputPointType::checkpass:
-      cout << "JSON written after Check step." << endl;
-      break;
-    case JsonOutputPointType::transformpass:
-      cout << "JSON written after Transform step." << endl;
-      break;
-    case JsonOutputPointType::computingpass:
-      cout << "JSON written after Computing step." << endl;
-      break;
-    case JsonOutputPointType::nojson:
-      cerr << "ModFile::writeJsonOutput: should not arrive here." << endl;
-      exit(EXIT_FAILURE);
-    }
-
-  if (onlyjson)
-    exit(EXIT_SUCCESS);
-}
-
-void
-ModFile::writeJsonOutputParsingCheck(const string& basename, JsonFileOutputType json_output_mode,
-                                     bool transformpass, bool computingpass) const
-{
-  ostringstream output;
-  output << "{" << endl;
-
-  symbol_table.writeJsonOutput(output);
-  output << ", ";
-  if (!heterogeneity_table.empty())
-    {
-      heterogeneity_table.writeJsonOutput(output);
-      output << ", ";
-    }
-  dynamic_model.writeJsonOutput(output);
-  output << ", ";
-  static_model.writeJsonOutput(output);
-
-  if (!statements.empty() || !var_model_table.empty() || !trend_component_model_table.empty())
-    {
-      output << R"(, "statements": [)";
-      if (!var_model_table.empty())
-        {
-          var_model_table.writeJsonOutput(output);
-          output << ", ";
-        }
-
-      if (!trend_component_model_table.empty())
-        {
-          trend_component_model_table.writeJsonOutput(output);
-          output << ", ";
-        }
-
-      if (!var_expectation_model_table.empty())
-        {
-          var_expectation_model_table.writeJsonOutput(output);
-          output << ", ";
-        }
-
-      if (!pac_model_table.empty())
-        {
-          pac_model_table.writeJsonOutput(output);
-          output << ", ";
-        }
-
-      for (bool printed_something {false}; auto& it : statements)
-        {
-          if (exchange(printed_something, true))
-            output << ", " << endl;
-          it->writeJsonOutput(output);
-        }
-      output << "]" << endl;
-    }
-
-  if (computingpass)
-    {
-      output << ",";
-      dynamic_model.writeJsonDynamicModelInfo(output);
-    }
-
-  output << R"(, "steady_state_model": )" << mod_file_struct.steady_state_model_present << endl;
-
-  output << "}" << endl;
-
-  ostringstream original_model_output;
-  original_model_output << "";
-  if (transformpass || computingpass)
-    {
-      original_model_output << "{";
-      original_model.writeJsonOriginalModelOutput(original_model_output);
-      if (!statements.empty() || !var_model_table.empty() || !trend_component_model_table.empty())
-        {
-          original_model_output << endl << R"(, "statements": [)";
-          if (!var_model_table.empty())
-            {
-              var_model_table.writeJsonOutput(original_model_output);
-              original_model_output << ", ";
-            }
-          if (!trend_component_model_table.empty())
-            {
-              trend_component_model_table.writeJsonOutput(original_model_output);
-              original_model_output << ", ";
-            }
-          if (!pac_model_table.empty())
-            {
-              pac_model_table.writeJsonOutput(original_model_output);
-              original_model_output << ", ";
-            }
-          for (bool printed_something {false}; const auto& it : statements)
-            {
-              original_model_output << (exchange(printed_something, true) ? "," : "") << endl;
-              it->writeJsonOutput(original_model_output);
-            }
-          original_model_output << "]" << endl;
-        }
-      original_model_output << "}" << endl;
-    }
-
-  ostringstream steady_state_model_output;
-  steady_state_model_output << "";
-  if (dynamic_model.equation_number() > 0)
-    steady_state_model.writeJsonSteadyStateFile(steady_state_model_output,
-                                                transformpass || computingpass);
-
-  if (json_output_mode == JsonFileOutputType::standardout)
-    {
-      if (transformpass || computingpass)
-        cout << R"("transformed_modfile": )";
-      else
-        cout << R"("modfile": )";
-      cout << output.str();
-      if (!original_model_output.str().empty())
-        cout << R"(, "original_model": )" << original_model_output.str();
-      if (!steady_state_model_output.str().empty())
-        cout << R"(, "steady_state_model": )" << steady_state_model_output.str();
-    }
-  else
-    {
-      if (!basename.size())
-        {
-          cerr << "ERROR: Missing file name" << endl;
-          exit(EXIT_FAILURE);
-        }
-
-      filesystem::create_directories(basename + "/model/json");
-      const filesystem::path fname {basename + "/model/json/modfile.json"};
-      ofstream jsonOutputFile {fname, ios::out | ios::binary};
-      if (!jsonOutputFile.is_open())
-        {
-          cerr << "ERROR: Can't open file " << fname.string() << " for writing" << endl;
-          exit(EXIT_FAILURE);
-        }
-
-      jsonOutputFile << output.str();
-      jsonOutputFile.close();
-
-      if (!original_model_output.str().empty())
-        {
-          if (basename.size())
-            {
-              const filesystem::path fname {basename + "/model/json/modfile-original.json"};
-              jsonOutputFile.open(fname, ios::out | ios::binary);
-              if (!jsonOutputFile.is_open())
-                {
-                  cerr << "ERROR: Can't open file " << fname.string() << " for writing" << endl;
-                  exit(EXIT_FAILURE);
-                }
-            }
-          else
-            {
-              cerr << "ERROR: Missing file name" << endl;
-              exit(EXIT_FAILURE);
-            }
-
-          jsonOutputFile << original_model_output.str();
-          jsonOutputFile.close();
-        }
-      if (!steady_state_model_output.str().empty())
-        {
-          if (basename.size())
-            {
-              const filesystem::path fname {basename + "/model/json/steady_state_model.json"};
-              jsonOutputFile.open(fname, ios::out | ios::binary);
-              if (!jsonOutputFile.is_open())
-                {
-                  cerr << "ERROR: Can't open file " << fname.string() << " for writing" << endl;
-                  exit(EXIT_FAILURE);
-                }
-            }
-          else
-            {
-              cerr << "ERROR: Missing file name" << endl;
-              exit(EXIT_FAILURE);
-            }
-
-          jsonOutputFile << steady_state_model_output.str();
-          jsonOutputFile.close();
-        }
-    }
-}
-
-void
-ModFile::writeJsonComputingPassOutput(const string& basename, JsonFileOutputType json_output_mode,
-                                      bool jsonderivsimple) const
-{
-  if (basename.empty() && json_output_mode != JsonFileOutputType::standardout)
-    {
-      cerr << "ERROR: Missing file name" << endl;
-      exit(EXIT_FAILURE);
-    }
-
-  ostringstream tmp_out, static_output, dynamic_output, static_paramsd_output,
-      dynamic_paramsd_output;
-
-  static_output << "{";
-  static_model.writeJsonComputingPassOutput(static_output, !jsonderivsimple);
-  static_output << "}";
-
-  dynamic_output << "{";
-  dynamic_model.writeJsonComputingPassOutput(dynamic_output, !jsonderivsimple);
-  dynamic_output << "}";
-
-  static_model.writeJsonParamsDerivatives(tmp_out, !jsonderivsimple);
-  if (!tmp_out.str().empty())
-    static_paramsd_output << "{" << tmp_out.str() << "}" << endl;
-
-  tmp_out.str("");
-  dynamic_model.writeJsonParamsDerivatives(tmp_out, !jsonderivsimple);
-  if (!tmp_out.str().empty())
-    dynamic_paramsd_output << "{" << tmp_out.str() << "}" << endl;
-
-  if (json_output_mode == JsonFileOutputType::standardout)
-    {
-      cout << R"(, "static_model": )" << static_output.str() << endl
-           << R"(, "dynamic_model": )" << dynamic_output.str() << endl;
-
-      if (!static_paramsd_output.str().empty())
-        cout << R"(, "static_params_deriv": )" << static_paramsd_output.str() << endl;
-
-      if (!dynamic_paramsd_output.str().empty())
-        cout << R"(, "dynamic_params_deriv": )" << dynamic_paramsd_output.str() << endl;
-    }
-  else
-    {
-      filesystem::create_directories(basename + "/model/json");
-
-      writeJsonFileHelper(basename + "/model/json/static.json", static_output);
-      writeJsonFileHelper(basename + "/model/json/dynamic.json", dynamic_output);
-
-      if (!static_paramsd_output.str().empty())
-        writeJsonFileHelper(basename + "/model/json/static_params_derivs.json",
-                            static_paramsd_output);
-
-      if (!dynamic_paramsd_output.str().empty())
-        writeJsonFileHelper(basename + "/model/json/params_derivs.json", dynamic_paramsd_output);
-    }
-}
-
-void
-ModFile::writeJsonFileHelper(const filesystem::path& fname, ostringstream& output) const
-{
-  ofstream jsonOutput {fname, ios::out | ios::binary};
-  if (!jsonOutput.is_open())
-    {
-      cerr << "ERROR: Can't open file " << fname.string() << " for writing" << endl;
-      exit(EXIT_FAILURE);
-    }
-  jsonOutput << output.str();
-  jsonOutput.close();
-}
-
 filesystem::path
 ModFile::unique_path()
 {
diff --git a/src/ModFile.hh b/src/ModFile.hh
index 97423f164423c84b31fe488f06c58bb9101e6d21..84523276bf4c85b3fbc2ee3d40329b75c113cdf1 100644
--- a/src/ModFile.hh
+++ b/src/ModFile.hh
@@ -35,6 +35,7 @@
 #include "ModelEquationBlock.hh"
 #include "NumericalConstants.hh"
 #include "NumericalInitialization.hh"
+#include "OutputWriter.hh"
 #include "Statement.hh"
 #include "StaticModel.hh"
 #include "SubModel.hh"
@@ -46,6 +47,12 @@ using namespace std;
 //! The abstract representation of a "mod" file
 class ModFile
 {
+  friend class OutputWriter;
+  friend class MatlabOutputWriter;
+  friend class JuliaOutputWriter;
+  friend class PythonOutputWriter;
+  friend class JsonOutputWriter;
+
 public:
   explicit ModFile(WarningConsolidation& warnings_arg);
   // For heterogeneity dimensions
@@ -135,23 +142,6 @@ private:
   ModFileStructure mod_file_struct;
   //! Warnings Encountered
   WarningConsolidation& warnings;
-  //! Functions used in writing of JSON outut. See writeJsonOutput
-  void writeJsonOutputParsingCheck(const string& basename, JsonFileOutputType json_output_mode,
-                                   bool transformpass, bool computingpass) const;
-  void writeJsonComputingPassOutput(const string& basename, JsonFileOutputType json_output_mode,
-                                    bool jsonderivsimple) const;
-  void writeJsonFileHelper(const filesystem::path& fname, ostringstream& output) const;
-  /* Generate a random temporary path, in the current directory. Equivalent to
-     boost::filesystem::unique_path(). Both are insecure, but currently there
-     is no better portable solution. Maybe in a later C++ standard? */
-  static filesystem::path unique_path();
-
-  /* Hack for removing a directory that is locked by MATLAB/Windows. The
-     directory is renamed before being deleted. The renaming must occur in the
-     same directory, otherwise it may fail if the destination is not on the
-     same filesystem. This technique is applied recursively to subdirectories.
-     Works even if the argument does not exist or is not a directory. */
-  static void remove_directory_with_matlab_lock(const filesystem::path& dir);
 
 public:
   //! Add a statement
@@ -174,35 +164,24 @@ public:
    * files */
   /*! \param params_derivs_order compute this order of derivs wrt parameters */
   void computingPass(bool no_tmp_terms, OutputType output, int params_derivs_order);
-  //! Writes Matlab/Octave output files
-  /*!
-    \param basename The base name used for writing output files. Should be the name of the mod file
-    without its extension \param clear_all Should a "clear all" instruction be written to output ?
-    \param console Are we in console mode ?
-    \param nograph Should we build the figures?
-    \param nointeractive Should Dynare request user input?
-    \param cygwin Should the MEX command of use_dll be adapted for Cygwin?
-    \param msvc Should the MEX command of use_dll be adapted for MSVC?
-    \param mingw Should the MEX command of use_dll be adapted for MinGW?
-    \param compute_xrefs if true, equation cross references will be computed
-  */
-  void writeMOutput(const string& basename, bool clear_all, bool clear_global, bool no_warn,
-                    bool console, bool nograph, bool nointeractive, const Configuration& config,
-                    bool check_model_changes, bool minimal_workspace, bool compute_xrefs,
-                    const string& mexext, const filesystem::path& matlabroot, bool onlymodel,
-                    bool gui, bool notime) const;
-
-  void writeJuliaOutput(const string& basename) const;
-
   void computeChecksum();
-  //! Write JSON representation of ModFile object
-  //! Initially created to enable Julia to work with .mod files
-  //! Potentially outputs ModFile after the various parts of processing (parsing, checkPass,
-  //! transformPass, computingPass) Allows user of other host language platforms (python, fortran,
-  //! etc) to provide support for dynare .mod files
-  void writeJsonOutput(const string& basename, JsonOutputPointType json,
-                       JsonFileOutputType json_output_mode, bool onlyjson,
-                       bool jsonderivsimple = false);
+  //! JSON helper methods that can be accessed by JsonOutputWriter
+  void writeJsonOutputParsingCheck(const string& basename, JsonFileOutputType json_output_mode,
+                                   bool transformpass, bool computingpass) const;
+  void writeJsonComputingPassOutput(const string& basename, JsonFileOutputType json_output_mode,
+                                    bool jsonderivsimple) const;
+  void writeJsonFileHelper(const filesystem::path& fname, ostringstream& output) const;
+  /* Generate a random temporary path, in the current directory. Equivalent to
+     boost::filesystem::unique_path(). Both are insecure, but currently there
+     is no better portable solution. Maybe in a later C++ standard? */
+  static filesystem::path unique_path();
+
+  /* Hack for removing a directory that is locked by MATLAB/Windows. The
+     directory is renamed before being deleted. The renaming must occur in the
+     same directory, otherwise it may fail if the destination is not on the
+     same filesystem. This technique is applied recursively to subdirectories.
+     Works even if the argument does not exist or is not a directory. */
+  static void remove_directory_with_matlab_lock(const filesystem::path& dir);
 };
 
 #endif
diff --git a/src/OutputWriter.cc b/src/OutputWriter.cc
new file mode 100644
index 0000000000000000000000000000000000000000..02389246a8f82fb8e38d2a24ba49471fc7c45d1e
--- /dev/null
+++ b/src/OutputWriter.cc
@@ -0,0 +1,706 @@
+/*
+ * Copyright © 2006-2025 Dynare Team
+ *
+ * This file is part of Dynare.
+ *
+ * Dynare is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Dynare is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
+ */
+
+#include "OutputWriter.hh"
+#include "ModFile.hh"
+#include "StaticModel.hh"
+#include <fstream>
+#include <iostream>
+
+// Base OutputWriter implementation
+OutputWriter::OutputWriter(const ModFile& modFile, const std::string& basename) :
+    modFile(modFile), basename(basename)
+{
+}
+
+// MATLAB OutputWriter implementation
+MatlabOutputWriter::MatlabOutputWriter(
+    const ModFile& modFile, const std::string& basename, bool clear_all, bool clear_global,
+    bool no_warn, bool console, bool nograph, bool nointeractive, const Configuration& config,
+    bool check_model_changes, bool minimal_workspace, bool compute_xrefs, const std::string& mexext,
+    const std::filesystem::path& matlabroot, bool onlymodel, bool gui, bool notime) :
+    OutputWriter(modFile, basename),
+    clear_all(clear_all),
+    clear_global(clear_global),
+    no_warn(no_warn),
+    console(console),
+    nograph(nograph),
+    nointeractive(nointeractive),
+    config(config),
+    check_model_changes(check_model_changes),
+    minimal_workspace(minimal_workspace),
+    compute_xrefs(compute_xrefs),
+    mexext(mexext),
+    matlabroot(matlabroot),
+    onlymodel(onlymodel),
+    gui(gui),
+    notime(notime)
+{
+}
+
+void
+MatlabOutputWriter::writeOutput() const
+{
+  if (basename.empty())
+    {
+      std::cerr << "ERROR: Missing file name" << std::endl;
+      exit(EXIT_FAILURE);
+    }
+
+  auto plusfolder = DataTree::packageDir(basename);
+
+  if (check_model_changes && !modFile.heterogeneity_table.empty())
+    {
+      std::cerr << "ERROR: the 'fast' option is not supported for heterogeneous models"
+                << std::endl;
+      exit(EXIT_FAILURE);
+    }
+  bool hasModelChanged
+      = !modFile.dynamic_model.isChecksumMatching(basename) || !check_model_changes;
+  if (hasModelChanged)
+    {
+      // Erase possible remnants of previous runs
+
+      /* Under MATLAB+Windows (but not under Octave nor under GNU/Linux or
+         macOS), if we directly remove the "+" subdirectory, then the
+         preprocessor is not able to recreate it afterwards (presumably because
+         MATLAB maintains some sort of lock on it). So we use a hack. */
+      ModFile::remove_directory_with_matlab_lock(plusfolder);
+
+      std::filesystem::remove_all(basename + "/model/src");
+      std::filesystem::remove_all(basename + "/model/bytecode");
+      // Do not remove basename/model/julia/, otherwise it would break calls to
+      // writeToFileIfModified()
+    }
+
+  create_directory(plusfolder);
+  std::filesystem::path fname = plusfolder / "driver.m";
+  std::ofstream mOutputFile(fname, std::ios::out | std::ios::binary);
+  if (!mOutputFile.is_open())
+    {
+      std::cerr << "ERROR: Can't open file " << fname.string() << " for writing" << std::endl;
+      exit(EXIT_FAILURE);
+    }
+
+  mOutputFile << "%" << std::endl
+              << "% Status : main Dynare file" << std::endl
+              << "%" << std::endl
+              << "% Warning : this file is generated automatically by Dynare" << std::endl
+              << "%           from model file (.mod)" << std::endl
+              << std::endl;
+
+  if (no_warn)
+    mOutputFile << "warning off"
+                << std::endl; // This will be executed *after* function warning_config()
+
+  if (clear_all)
+    mOutputFile << "clearvars -global" << std::endl
+                << "clear_persistent_variables(fileparts(which('dynare')), false)" << std::endl;
+  else if (clear_global)
+    mOutputFile << "clearvars -global M_ options_ oo_ estim_params_ bayestopt_ dataset_ "
+                   "dataset_info estimation_info;"
+                << std::endl;
+
+  if (!notime)
+    mOutputFile << "tic0 = tic;" << std::endl;
+
+  mOutputFile
+      << "% Define global variables." << std::endl
+      << "global M_ options_ oo_ estim_params_ bayestopt_ dataset_ dataset_info estimation_info"
+      << std::endl
+      << "options_ = [];" << std::endl
+      << "M_.fname = '" << basename << "';" << std::endl
+      << "M_.dynare_version = '" << PACKAGE_VERSION << "';" << std::endl
+      << "oo_.dynare_version = '" << PACKAGE_VERSION << "';" << std::endl
+      << "options_.dynare_version = '" << PACKAGE_VERSION << "';" << std::endl
+      << "%" << std::endl
+      << "% Some global variables initialization" << std::endl
+      << "%" << std::endl;
+  if (!onlymodel)
+    config.writeHooks(mOutputFile);
+  mOutputFile << "global_initialization;" << std::endl;
+
+  if (console)
+    mOutputFile << "options_.console_mode = true;" << std::endl
+                << "options_.nodisplay = true;" << std::endl;
+  if (nograph)
+    mOutputFile << "options_.nograph = true;" << std::endl;
+
+  if (nointeractive)
+    mOutputFile << "options_.nointeractive = true;" << std::endl;
+
+  if (modFile.param_used_with_lead_lag)
+    mOutputFile << "M_.parameter_used_with_lead_lag = true;" << std::endl;
+
+  modFile.symbol_table.writeOutput(mOutputFile);
+  modFile.heterogeneity_table.writeOutput(mOutputFile);
+
+  modFile.var_model_table.writeOutput(basename, mOutputFile);
+  modFile.trend_component_model_table.writeOutput(basename, mOutputFile);
+  modFile.var_expectation_model_table.writeOutput(mOutputFile);
+  modFile.pac_model_table.writeOutput(mOutputFile);
+
+  // Initialize M_.Sigma_e, M_.Correlation_matrix, M_.H, and M_.Correlation_matrix_ME
+  mOutputFile << "M_.Sigma_e = zeros(" << modFile.symbol_table.exo_nbr() << ", "
+              << modFile.symbol_table.exo_nbr() << ");" << std::endl
+              << "M_.Correlation_matrix = eye(" << modFile.symbol_table.exo_nbr() << ", "
+              << modFile.symbol_table.exo_nbr() << ");" << std::endl;
+  for (int hd = 0; hd < modFile.heterogeneity_table.size(); hd++)
+    mOutputFile << "M_.heterogeneity(" << hd + 1 << ").Sigma_e = zeros("
+                << modFile.symbol_table.het_exo_nbr(hd) << ", "
+                << modFile.symbol_table.het_exo_nbr(hd) << ");" << std::endl;
+
+  if (modFile.mod_file_struct.calibrated_measurement_errors)
+    mOutputFile << "M_.H = zeros(" << modFile.symbol_table.observedVariablesNbr() << ", "
+                << modFile.symbol_table.observedVariablesNbr() << ");" << std::endl
+                << "M_.Correlation_matrix_ME = eye(" << modFile.symbol_table.observedVariablesNbr()
+                << ", " << modFile.symbol_table.observedVariablesNbr() << ");" << std::endl;
+  else
+    mOutputFile << "M_.H = 0;" << std::endl << "M_.Correlation_matrix_ME = 1;" << std::endl;
+
+  // May be later modified by a shocks block
+  mOutputFile << "M_.sigma_e_is_diagonal = true;" << std::endl;
+
+  /* Initialize the structures created for several blocks, as part of the implementation of the
+     "overwrite" option */
+  mOutputFile << "M_.det_shocks = [];" << std::endl
+              << "M_.surprise_shocks = [];" << std::endl
+              << "M_.learnt_shocks = [];" << std::endl
+              << "M_.learnt_endval = [];" << std::endl
+              << "M_.heteroskedastic_shocks.Qvalue_orig = [];" << std::endl
+              << "M_.heteroskedastic_shocks.Qscale_orig = [];" << std::endl
+              << "M_.matched_irfs = {};" << std::endl
+              << "M_.matched_irfs_weights = {};" << std::endl
+              << "M_.perfect_foresight_controlled_paths = [];" << std::endl;
+
+  // NB: options_.{ramsey,discretionary}_policy should rather be fields of M_
+  mOutputFile << std::boolalpha << "options_.linear = " << modFile.linear << ";" << std::endl
+              << "options_.block = " << modFile.block << ";" << std::endl
+              << "options_.bytecode = " << modFile.bytecode << ";" << std::endl
+              << "options_.use_dll = " << modFile.use_dll << ";" << std::endl
+              << "options_.ramsey_policy = " << modFile.mod_file_struct.ramsey_model_present << ";"
+              << std::endl
+              << "options_.discretionary_policy = "
+              << modFile.mod_file_struct.discretionary_policy_present << ";" << std::endl;
+
+  if (modFile.mod_file_struct.discretionary_policy_present)
+    mOutputFile << "M_.discretionary_orig_eq_nbr = " << modFile.original_model.equation_number()
+                << ";" << std::endl;
+
+  if (modFile.parallel_local_files.size() > 0)
+    {
+      mOutputFile << "options_.parallel_info.local_files = {" << std::endl;
+      for (const auto& parallel_local_file : modFile.parallel_local_files)
+        {
+          size_t j = parallel_local_file.find_last_of(R"(/\)");
+          if (j == std::string::npos)
+            mOutputFile << "'', '" << parallel_local_file << "';" << std::endl;
+          else
+            mOutputFile << "'" << parallel_local_file.substr(0, j + 1) << "', '"
+                        << parallel_local_file.substr(j + 1) << "';" << std::endl;
+        }
+      mOutputFile << "};" << std::endl;
+    }
+
+  if (modFile.dynamic_model.isHessianComputed())
+    {
+      mOutputFile << "M_.nonzero_hessian_eqs = ";
+      modFile.dynamic_model.printNonZeroHessianEquations(mOutputFile);
+      mOutputFile << ";" << std::endl
+                  << "M_.hessian_eq_zero = isempty(M_.nonzero_hessian_eqs);" << std::endl;
+    }
+
+  if (!onlymodel)
+    config.writeCluster(mOutputFile);
+
+  if (modFile.bytecode)
+    mOutputFile << "if exist('bytecode') ~= 3" << std::endl
+                << "  error('DYNARE: Can''t find bytecode DLL. Please compile it or remove the "
+                   "'bytecode' option.')"
+                << std::endl
+                << "end" << std::endl;
+
+  mOutputFile << "M_.eq_nbr = " << modFile.dynamic_model.equation_number() << ";" << std::endl
+              << "M_.ramsey_orig_eq_nbr = " << modFile.mod_file_struct.ramsey_orig_eq_nbr << ";"
+              << std::endl
+              << "M_.ramsey_orig_endo_nbr = " << modFile.mod_file_struct.ramsey_orig_endo_nbr << ";"
+              << std::endl
+              << "M_.set_auxiliary_variables = exist(['./+' M_.fname "
+                 "'/set_auxiliary_variables.m'], 'file') == 2;"
+              << std::endl;
+
+  modFile.epilogue.writeOutput(mOutputFile);
+
+  if (modFile.dynamic_model.equation_number() > 0)
+    {
+      modFile.dynamic_model.writeDriverOutput(mOutputFile, compute_xrefs);
+      if (!modFile.no_static)
+        {
+          modFile.static_model.writeDriverOutput(mOutputFile);
+          if (modFile.mod_file_struct.ramsey_model_present)
+            modFile.static_model.writeDriverRamseyMultipliersDerivativesSparseIndices(mOutputFile);
+        }
+    }
+
+  for (const auto& hm : modFile.heterogeneous_models)
+    hm.writeDriverOutput(mOutputFile);
+
+  if (onlymodel || gui)
+    for (const auto& statement : modFile.statements)
+      {
+        /* Special treatment for initval block: insert initial values for the
+           auxiliary variables and initialize exo det */
+        if (auto ivs = dynamic_cast<InitValStatement*>(statement.get()); ivs)
+          {
+            ivs->writeOutput(mOutputFile, basename, minimal_workspace);
+            modFile.static_model.writeAuxVarInitval(mOutputFile,
+                                                    ExprNodeOutputType::matlabOutsideModel);
+            ivs->writeOutputPostInit(mOutputFile);
+          }
+
+        // Special treatment for endval block: insert initial values for the auxiliary variables
+        if (auto evs = dynamic_cast<EndValStatement*>(statement.get()); evs)
+          {
+            evs->writeOutput(mOutputFile, basename, minimal_workspace);
+            modFile.static_model.writeAuxVarInitval(mOutputFile,
+                                                    ExprNodeOutputType::matlabOutsideModel);
+          }
+
+        if (auto ips = dynamic_cast<InitParamStatement*>(statement.get()); ips)
+          ips->writeOutput(mOutputFile, basename, minimal_workspace);
+
+        if (auto ss = dynamic_cast<ShocksStatement*>(statement.get()); ss)
+          ss->writeOutput(mOutputFile, basename, minimal_workspace);
+
+        if (auto eps = dynamic_cast<EstimatedParamsStatement*>(statement.get()); eps)
+          eps->writeOutput(mOutputFile, basename, minimal_workspace);
+
+        if (auto sgs = dynamic_cast<ShockGroupsStatement*>(statement.get()); sgs)
+          sgs->writeOutput(mOutputFile, basename, minimal_workspace);
+
+        if (gui)
+          if (auto it = dynamic_cast<NativeStatement*>(statement.get()); it)
+            it->writeOutput(mOutputFile, basename, minimal_workspace);
+      }
+  else
+    {
+      for (const auto& statement : modFile.statements)
+        {
+          statement->writeOutput(mOutputFile, basename, minimal_workspace);
+
+          /* Special treatment for initval block: insert initial values for the
+             auxiliary variables and initialize exo det */
+          if (auto ivs = dynamic_cast<InitValStatement*>(statement.get()); ivs)
+            {
+              modFile.static_model.writeAuxVarInitval(mOutputFile,
+                                                      ExprNodeOutputType::matlabOutsideModel);
+              ivs->writeOutputPostInit(mOutputFile);
+            }
+
+          // Special treatment for endval block: insert initial values for the auxiliary variables
+          if (auto evs = dynamic_cast<EndValStatement*>(statement.get()); evs)
+            modFile.static_model.writeAuxVarInitval(mOutputFile,
+                                                    ExprNodeOutputType::matlabOutsideModel);
+
+          // Special treatment for load params and steady state statement: insert initial values for
+          // the auxiliary variables
+          if (auto lpass = dynamic_cast<LoadParamsAndSteadyStateStatement*>(statement.get());
+              lpass && !modFile.no_static)
+            modFile.static_model.writeAuxVarInitval(mOutputFile,
+                                                    ExprNodeOutputType::matlabOutsideModel);
+        }
+
+      if (!notime)
+        mOutputFile << std::endl
+                    << std::endl
+                    << "oo_.time = toc(tic0);" << std::endl
+                    << "disp(['Total computing time : ' dynsec2hms(oo_.time) ]);" << std::endl;
+
+      mOutputFile << "if ~exist([M_.dname filesep 'Output'],'dir')" << std::endl
+                  << "    mkdir(M_.dname,'Output');" << std::endl
+                  << "end" << std::endl
+                  << "save([M_.dname filesep 'Output' filesep '" << basename
+                  << "_results.mat'], 'oo_', 'M_', 'options_');" << std::endl
+                  << "if exist('estim_params_', 'var') == 1" << std::endl
+                  << "  save([M_.dname filesep 'Output' filesep '" << basename
+                  << "_results.mat'], 'estim_params_', '-append');" << std::endl
+                  << "end" << std::endl
+                  << "if exist('bayestopt_', 'var') == 1" << std::endl
+                  << "  save([M_.dname filesep 'Output' filesep '" << basename
+                  << "_results.mat'], 'bayestopt_', '-append');" << std::endl
+                  << "end" << std::endl
+                  << "if exist('dataset_', 'var') == 1" << std::endl
+                  << "  save([M_.dname filesep 'Output' filesep '" << basename
+                  << "_results.mat'], 'dataset_', '-append');" << std::endl
+                  << "end" << std::endl
+                  << "if exist('estimation_info', 'var') == 1" << std::endl
+                  << "  save([M_.dname filesep 'Output' filesep '" << basename
+                  << "_results.mat'], 'estimation_info', '-append');" << std::endl
+                  << "end" << std::endl
+                  << "if exist('dataset_info', 'var') == 1" << std::endl
+                  << "  save([M_.dname filesep 'Output' filesep '" << basename
+                  << "_results.mat'], 'dataset_info', '-append');" << std::endl
+                  << "end" << std::endl
+                  << "if exist('oo_recursive_', 'var') == 1" << std::endl
+                  << "  save([M_.dname filesep 'Output' filesep '" << basename
+                  << "_results.mat'], 'oo_recursive_', '-append');" << std::endl
+                  << "end" << std::endl
+                  << "if exist('options_mom_', 'var') == 1" << std::endl
+                  << "  save([M_.dname filesep 'Output' filesep '" << basename
+                  << "_results.mat'], 'options_mom_', '-append');" << std::endl
+                  << "end" << std::endl;
+
+      config.writeEndParallel(mOutputFile);
+
+      if (!no_warn)
+        {
+          if (int num_warnings = modFile.warnings.numWarnings(); num_warnings > 0)
+            mOutputFile << "disp('Note: " << num_warnings
+                        << " warning(s) encountered in the preprocessor')" << std::endl;
+
+          mOutputFile << "if ~isempty(lastwarn)" << std::endl
+                      << "  disp('Note: warning(s) encountered in MATLAB/Octave code')" << std::endl
+                      << "end" << std::endl;
+        }
+    }
+
+  mOutputFile.close();
+
+  if (hasModelChanged)
+    {
+      // Create static and dynamic files
+      if (modFile.dynamic_model.equation_number() > 0)
+        {
+          if (!modFile.no_static)
+            {
+              this->writeStaticFile(modFile.static_model, modFile.use_dll, mexext, matlabroot);
+              modFile.static_model.writeParamsDerivativesFile<false>(basename);
+              if (modFile.mod_file_struct.ramsey_model_present)
+                {
+                  if (modFile.use_dll)
+                    modFile.static_model.writeRamseyMultipliersDerivativesCFile(
+                        basename, mexext, matlabroot, modFile.mod_file_struct.ramsey_orig_endo_nbr);
+                  else
+                    modFile.static_model.writeRamseyMultipliersDerivativesMFile(
+                        basename, modFile.mod_file_struct.ramsey_orig_endo_nbr);
+                }
+            }
+
+          modFile.dynamic_model.writeDynamicFile(basename, modFile.use_dll, mexext, matlabroot,
+                                                 false);
+
+          modFile.dynamic_model.writeParamsDerivativesFile<false>(basename);
+        }
+
+      // Create steady state file
+      modFile.steady_state_model.writeSteadyStateFile(basename, false);
+
+      // Create epilogue file
+      modFile.epilogue.writeEpilogueFile(basename);
+
+      modFile.pac_model_table.writeTargetCoefficientsFile(basename);
+
+      for (const auto& hm : modFile.heterogeneous_models)
+        hm.writeModelFiles(basename, false);
+    }
+}
+
+void
+MatlabOutputWriter::writeStaticFile(const StaticModel& staticModel, bool use_dll,
+                                    const std::string& mexext,
+                                    const std::filesystem::path& matlabroot) const
+{
+  // Implementation for MATLAB-specific static file writing
+  std::filesystem::path model_dir {basename};
+  model_dir /= "model";
+  if (use_dll)
+    {
+      create_directories(model_dir / "src" / "sparse");
+      if (staticModel.block_decomposed)
+        create_directories(model_dir / "src" / "sparse" / "block");
+    }
+
+  // Create the needed MATLAB directories
+  auto plusfolder {StaticModel::packageDir(basename)};
+  /* The following is not a duplicate of the same call from
+     ModFile::writeMOutput(), because of planner_objective which needs its
+     +objective subdirectory */
+  create_directories(plusfolder);
+
+  auto sparsefolder {plusfolder / "+sparse"};
+  create_directories(sparsefolder);
+  if (staticModel.block_decomposed)
+    create_directories(sparsefolder / "+block");
+
+  create_directories(plusfolder / "+debug");
+  create_directories(model_dir / "bytecode" / "block");
+
+  /* PlannerObjective subclass or discretionary optimal policy models don't
+     have as many variables as equations; bytecode does not support that
+     case */
+  if (static_cast<int>(staticModel.equations.size()) == staticModel.symbol_table.endo_nbr())
+    staticModel.writeStaticBytecode(basename);
+
+  if (staticModel.block_decomposed)
+    staticModel.writeStaticBlockBytecode(basename);
+
+  // Sparse representation
+  if (use_dll)
+    staticModel.writeSparseModelCFiles<false>(basename, mexext, matlabroot);
+  else // MATLAB/Octave
+    staticModel.writeSparseModelMFiles<false>(basename);
+
+  staticModel.writeSetAuxiliaryVariablesFile<false>(basename, false);
+  staticModel.writeComplementarityConditionsFile<false>(basename);
+
+  // Support for model debugging
+  staticModel.writeDebugModelMFiles<false>(basename);
+}
+
+// Julia OutputWriter implementation
+JuliaOutputWriter::JuliaOutputWriter(const ModFile& modFile, const std::string& basename) :
+    OutputWriter(modFile, basename)
+{
+}
+
+void
+JuliaOutputWriter::writeOutput() const
+{
+  if (modFile.dynamic_model.equation_number() > 0)
+    {
+      if (!modFile.no_static)
+        {
+          modFile.static_model.writeStaticFile(basename, false, "", {}, true);
+          modFile.static_model.writeParamsDerivativesFile<true>(basename);
+        }
+      modFile.dynamic_model.writeDynamicFile(basename, modFile.use_dll, "", {}, true);
+      modFile.dynamic_model.writeParamsDerivativesFile<true>(basename);
+    }
+  modFile.steady_state_model.writeSteadyStateFile(basename, true);
+}
+
+void
+JuliaOutputWriter::writeStaticFile(const StaticModel& staticModel, bool use_dll,
+                                   const std::string& mexext,
+                                   const std::filesystem::path& matlabroot) const
+{
+  // Implementation for Julia-specific static file writing
+  std::filesystem::path model_dir {basename};
+  model_dir /= "model";
+  if (use_dll)
+    {
+      create_directories(model_dir / "src" / "sparse");
+      if (staticModel.block_decomposed)
+        create_directories(model_dir / "src" / "sparse" / "block");
+    }
+
+  // Create the needed Julia directories
+  create_directories(model_dir / "julia");
+  create_directories(model_dir / "bytecode" / "block");
+
+  /* PlannerObjective subclass or discretionary optimal policy models don't
+     have as many variables as equations; bytecode does not support that
+     case */
+  if (static_cast<int>(staticModel.equations.size()) == staticModel.symbol_table.endo_nbr())
+    staticModel.writeStaticBytecode(basename);
+
+  if (staticModel.block_decomposed)
+    staticModel.writeStaticBlockBytecode(basename);
+
+  // Sparse representation
+  if (use_dll)
+    staticModel.writeSparseModelCFiles<false>(basename, mexext, matlabroot);
+  else // Julia files
+    staticModel.writeSparseModelJuliaFiles<false>(basename);
+
+  staticModel.writeSetAuxiliaryVariablesFile<false>(basename, true);
+}
+
+// Python OutputWriter implementation
+PythonOutputWriter::PythonOutputWriter(const ModFile& modFile, const std::string& basename) :
+    OutputWriter(modFile, basename)
+{
+}
+
+void
+PythonOutputWriter::writeOutput() const
+{
+  if (modFile.dynamic_model.equation_number() > 0)
+    {
+      if (!modFile.no_static)
+        {
+          modFile.static_model.writeStaticFile(basename, false, "", {}, true);
+          modFile.static_model.writeParamsDerivativesFile<true>(basename);
+        }
+      modFile.dynamic_model.writeDynamicFile(basename, modFile.use_dll, "", {}, true);
+      modFile.dynamic_model.writeParamsDerivativesFile<true>(basename);
+    }
+  modFile.steady_state_model.writeSteadyStateFile(basename, true);
+
+  // Create Python-specific interface files
+  std::filesystem::path model_dir(basename);
+  model_dir /= "model";
+  std::filesystem::path python_dir = model_dir / "python";
+  create_directories(python_dir);
+
+  // Write Python initialization file
+  std::ofstream init_file(python_dir / "__init__.py");
+  if (init_file.is_open())
+    {
+      init_file << "# Generated by Dynare preprocessor\n\n";
+      init_file << "from .model import *\n";
+      init_file.close();
+    }
+
+  // Write Python model interface file
+  std::ofstream model_file(python_dir / "model.py");
+  if (model_file.is_open())
+    {
+      model_file << "# Generated by Dynare preprocessor\n\n";
+      model_file << "import numpy as np\n";
+      model_file << "from scipy.optimize import root\n\n";
+
+      // Add wrapper functions to call Julia-generated model code
+      model_file << "def load_model():\n";
+      model_file << "    \"\"\"Load and initialize the model\"\"\"\n";
+      model_file << "    # TODO: Add Python-specific model loading code\n";
+      model_file << "    return {'M_': {}, 'options_': {}, 'oo_': {}}\n\n";
+
+      // Add other necessary functions similar to Julia interface
+      model_file << "def steady_state(model):\n";
+      model_file << "    \"\"\"Compute steady state of the model\"\"\"\n";
+      model_file << "    # TODO: Implement steady state computation\n";
+      model_file << "    pass\n\n";
+
+      model_file.close();
+    }
+}
+
+void
+PythonOutputWriter::writeStaticFile(const StaticModel& staticModel, bool use_dll,
+                                    const std::string& mexext,
+                                    const std::filesystem::path& matlabroot) const
+{
+  // Implementation for Python-specific static file writing
+  std::filesystem::path model_dir {basename};
+  model_dir /= "model";
+  if (use_dll)
+    {
+      create_directories(model_dir / "src" / "sparse");
+      if (staticModel.block_decomposed)
+        create_directories(model_dir / "src" / "sparse" / "block");
+    }
+
+  // Create needed Python directories
+  create_directories(model_dir / "python");
+  create_directories(model_dir / "bytecode" / "block");
+
+  /* PlannerObjective subclass or discretionary optimal policy models don't
+     have as many variables as equations; bytecode does not support that
+     case */
+  if (static_cast<int>(staticModel.equations.size()) == staticModel.symbol_table.endo_nbr())
+    staticModel.writeStaticBytecode(basename);
+
+  if (staticModel.block_decomposed)
+    staticModel.writeStaticBlockBytecode(basename);
+
+  // For Python, we can leverage the Julia writer and then add Python-specific files
+  if (use_dll)
+    staticModel.writeSparseModelCFiles<false>(basename, mexext, matlabroot);
+  else
+    staticModel.writeSparseModelJuliaFiles<false>(basename);
+
+  staticModel.writeSetAuxiliaryVariablesFile<false>(basename, true);
+
+  // Additional Python-specific code would go here
+}
+
+// JSON OutputWriter implementation
+JsonOutputWriter::JsonOutputWriter(const ModFile& modFile, const std::string& basename,
+                                   JsonOutputPointType json, JsonFileOutputType json_output_mode,
+                                   bool onlyjson, bool jsonderivsimple) :
+    OutputWriter(modFile, basename),
+    json(json),
+    json_output_mode(json_output_mode),
+    onlyjson(onlyjson),
+    jsonderivsimple(jsonderivsimple)
+{
+}
+
+void
+JsonOutputWriter::writeOutput() const
+{
+  if (json == JsonOutputPointType::nojson)
+    return;
+
+  if (json_output_mode == JsonFileOutputType::standardout)
+    std::cout << "//-- BEGIN JSON --// " << std::endl << "{" << std::endl;
+
+  modFile.writeJsonOutputParsingCheck(basename, json_output_mode,
+                                      json == JsonOutputPointType::transformpass,
+                                      json == JsonOutputPointType::computingpass);
+
+  if (json == JsonOutputPointType::computingpass)
+    modFile.writeJsonComputingPassOutput(basename, json_output_mode, jsonderivsimple);
+
+  if (json_output_mode == JsonFileOutputType::standardout)
+    std::cout << "}" << std::endl << "//-- END JSON --// " << std::endl;
+
+  switch (json)
+    {
+    case JsonOutputPointType::parsing:
+      std::cout << "JSON written after Parsing step." << std::endl;
+      break;
+    case JsonOutputPointType::checkpass:
+      std::cout << "JSON written after Check step." << std::endl;
+      break;
+    case JsonOutputPointType::transformpass:
+      std::cout << "JSON written after Transform step." << std::endl;
+      break;
+    case JsonOutputPointType::computingpass:
+      std::cout << "JSON written after Computing step." << std::endl;
+      break;
+    case JsonOutputPointType::nojson:
+      std::cerr << "ModFile::writeJsonOutput: should not arrive here." << std::endl;
+      exit(EXIT_FAILURE);
+    }
+
+  if (onlyjson)
+    exit(EXIT_SUCCESS);
+}
+
+void
+JsonOutputWriter::writeStaticFile(const StaticModel& staticModel, bool a, const std::string& mexext,
+                                  const std::filesystem::path& matlabroot) const
+{
+  // Implementation for JSON-specific static file writing
+  std::filesystem::path model_dir {basename};
+  model_dir /= "model";
+
+  // Create required directories
+  create_directories(model_dir / "bytecode" / "block");
+
+  if (static_cast<int>(staticModel.equations.size()) == staticModel.symbol_table.endo_nbr())
+    staticModel.writeStaticBytecode(basename);
+
+  if (staticModel.block_decomposed)
+    staticModel.writeStaticBlockBytecode(basename);
+
+  // JSON doesn't need all the output files, just the JSON representation
+}
diff --git a/src/OutputWriter.hh b/src/OutputWriter.hh
new file mode 100644
index 0000000000000000000000000000000000000000..c57307fdae0faf81b7a5245ae7cd8f79ca706581
--- /dev/null
+++ b/src/OutputWriter.hh
@@ -0,0 +1,140 @@
+/*
+ * Copyright © 2006-2025 Dynare Team
+ *
+ * This file is part of Dynare.
+ *
+ * Dynare is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Dynare is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
+ */
+
+#ifndef _OUTPUTWRITER_HH
+#define _OUTPUTWRITER_HH
+
+#include "Configuration.hh"
+#include "ExtendedPreprocessorTypes.hh"
+#include <filesystem>
+#include <string>
+
+class ModFile;
+class StaticModel;
+
+/**
+ * Abstract OutputWriter class following the Builder pattern
+ * Defines the interface for various output writers (MATLAB, Julia, Python, JSON)
+ */
+class OutputWriter
+{
+protected:
+  const ModFile& modFile;
+  const std::string& basename;
+
+public:
+  OutputWriter(const ModFile& modFile, const std::string& basename);
+  virtual ~OutputWriter() = default;
+
+  /**
+   * Main method to generate output files
+   */
+  virtual void writeOutput() const = 0;
+
+  /**
+   * Abstract method to write static file - implemented by subclasses
+   */
+  virtual void writeStaticFile(const StaticModel& staticModel, bool use_dll,
+                               const std::string& mexext,
+                               const std::filesystem::path& matlabroot) const
+      = 0;
+};
+
+/**
+ * Concrete builder for MATLAB output
+ */
+class MatlabOutputWriter : public OutputWriter
+{
+private:
+  bool clear_all;
+  bool clear_global;
+  bool no_warn;
+  bool console;
+  bool nograph;
+  bool nointeractive;
+  const Configuration& config;
+  bool check_model_changes;
+  bool minimal_workspace;
+  bool compute_xrefs;
+  const std::string& mexext;
+  const std::filesystem::path& matlabroot;
+  bool onlymodel;
+  bool gui;
+  bool notime;
+
+public:
+  MatlabOutputWriter(const ModFile& modFile, const std::string& basename, bool clear_all,
+                     bool clear_global, bool no_warn, bool console, bool nograph,
+                     bool nointeractive, const Configuration& config, bool check_model_changes,
+                     bool minimal_workspace, bool compute_xrefs, const std::string& mexext,
+                     const std::filesystem::path& matlabroot, bool onlymodel, bool gui,
+                     bool notime);
+
+  void writeOutput() const override;
+  void writeStaticFile(const StaticModel& staticModel, bool use_dll, const std::string& mexext,
+                       const std::filesystem::path& matlabroot) const override;
+};
+
+/**
+ * Concrete builder for Julia output
+ */
+class JuliaOutputWriter : public OutputWriter
+{
+public:
+  JuliaOutputWriter(const ModFile& modFile, const std::string& basename);
+
+  void writeOutput() const override;
+  void writeStaticFile(const StaticModel& staticModel, bool use_dll, const std::string& mexext,
+                       const std::filesystem::path& matlabroot) const override;
+};
+
+/**
+ * Concrete builder for Python output
+ */
+class PythonOutputWriter : public OutputWriter
+{
+public:
+  PythonOutputWriter(const ModFile& modFile, const std::string& basename);
+
+  void writeOutput() const override;
+  void writeStaticFile(const StaticModel& staticModel, bool use_dll, const std::string& mexext,
+                       const std::filesystem::path& matlabroot) const override;
+};
+
+/**
+ * Concrete builder for JSON output
+ */
+class JsonOutputWriter : public OutputWriter
+{
+private:
+  JsonOutputPointType json;
+  JsonFileOutputType json_output_mode;
+  bool onlyjson;
+  bool jsonderivsimple;
+
+public:
+  JsonOutputWriter(const ModFile& modFile, const std::string& basename, JsonOutputPointType json,
+                   JsonFileOutputType json_output_mode, bool onlyjson, bool jsonderivsimple);
+
+  void writeOutput() const override;
+  void writeStaticFile(const StaticModel& staticModel, bool use_dll, const std::string& mexext,
+                       const std::filesystem::path& matlabroot) const override;
+};
+
+#endif // _OUTPUTWRITER_HH
diff --git a/src/StaticModel.cc b/src/StaticModel.cc
index 93944b89d5825476a71db128b47b32fd8303e5e9..1c16914a0f3b879ce0139920c488e3aee06be1a4 100644
--- a/src/StaticModel.cc
+++ b/src/StaticModel.cc
@@ -256,63 +256,6 @@ StaticModel::computingPass(int derivsOrder, int paramsDerivsOrder,
   computeMCPEquationsReordering();
 }
 
-void
-StaticModel::writeStaticFile(const string& basename, bool use_dll, const string& mexext,
-                             const filesystem::path& matlabroot, bool julia) const
-{
-  filesystem::path model_dir {basename};
-  model_dir /= "model";
-  if (use_dll)
-    {
-      create_directories(model_dir / "src" / "sparse");
-      if (block_decomposed)
-        create_directories(model_dir / "src" / "sparse" / "block");
-    }
-  if (julia)
-    create_directories(model_dir / "julia");
-  else
-    {
-      auto plusfolder {packageDir(basename)};
-      /* The following is not a duplicate of the same call from
-         ModFile::writeMOutput(), because of planner_objective which needs its
-         +objective subdirectory */
-      create_directories(plusfolder);
-
-      auto sparsefolder {plusfolder / "+sparse"};
-      create_directories(sparsefolder);
-      if (block_decomposed)
-        create_directories(sparsefolder / "+block");
-
-      create_directories(plusfolder / "+debug");
-    }
-  create_directories(model_dir / "bytecode" / "block");
-
-  /* PlannerObjective subclass or discretionary optimal policy models don’t
-     have as many variables as equations; bytecode does not support that
-     case */
-  if (static_cast<int>(equations.size()) == symbol_table.endo_nbr())
-    writeStaticBytecode(basename);
-  if (block_decomposed)
-    writeStaticBlockBytecode(basename);
-
-  // Sparse representation
-  if (use_dll)
-    writeSparseModelCFiles<false>(basename, mexext, matlabroot);
-  else if (julia)
-    writeSparseModelJuliaFiles<false>(basename);
-  else // MATLAB/Octave
-    writeSparseModelMFiles<false>(basename);
-
-  writeSetAuxiliaryVariablesFile<false>(basename, julia);
-
-  if (!julia)
-    {
-      writeComplementarityConditionsFile<false>(basename);
-      // Support for model debugging
-      writeDebugModelMFiles<false>(basename);
-    }
-}
-
 bool
 StaticModel::exoPresentInEqs() const
 {
diff --git a/src/StaticModel.hh b/src/StaticModel.hh
index eb6ebe8bf6b0fa890a99a644d2876358cf0602f0..12429383103e1f8e5cb5fa4fccb7b9d51aabfe5e 100644
--- a/src/StaticModel.hh
+++ b/src/StaticModel.hh
@@ -25,10 +25,12 @@
 
 #include "Bytecode.hh"
 #include "ModelTree.hh"
+#include "OutputWriter.hh"
 
 using namespace std;
 
 class DynamicModel;
+class OutputWriter;
 
 //! Stores a static model, as derived from the "model" block when leads and lags have been removed
 class StaticModel : public ModelTree
@@ -144,10 +146,6 @@ public:
   void computingPass(int derivsOrder, int paramsDerivsOrder, const eval_context_t& eval_context,
                      bool no_tmp_terms, bool block, bool use_dll);
 
-  //! Writes static model file (+ bytecode)
-  void writeStaticFile(const string& basename, bool use_dll, const string& mexext,
-                       const filesystem::path& matlabroot, bool julia) const;
-
   //! Write JSON Output (used by PlannerObjectiveStatement)
   void writeJsonOutput(ostream& output) const;