diff --git a/src/DataTree.cc b/src/DataTree.cc
index cb1931b8326ace2d622e2bfdfb86030283a2d2d5..37d429f52de6ba78baf4f7f24c95cfed67fc67b3 100644
--- a/src/DataTree.cc
+++ b/src/DataTree.cc
@@ -904,13 +904,6 @@ DataTree::minLagForSymbol(int symb_id) const
   return r;
 }
 
-void
-DataTree::writePowerDerivCHeader(ostream &output) const
-{
-  if (isBinaryOpUsed(BinaryOpcode::powerDeriv))
-    output << "double getPowerDeriv(double, double, int);" << endl;
-}
-
 void
 DataTree::writePowerDeriv(ostream &output) const
 {
diff --git a/src/DataTree.hh b/src/DataTree.hh
index d0dde84662a29435edc068f169fc523a2c73cd55..23ce9c11ff8eeb062ae36c28aa7174edf89e0c7b 100644
--- a/src/DataTree.hh
+++ b/src/DataTree.hh
@@ -278,8 +278,6 @@ public:
   //! Returns the minimum lag (as a negative number) of the given symbol in the whole data tree (and not only in the equations !!)
   /*! Returns 0 if the symbol is not used */
   int minLagForSymbol(int symb_id) const;
-  //! Write the C Header for getPowerDeriv when use_dll is used
-  void writePowerDerivCHeader(ostream &output) const;
   //! Write getPowerDeriv in C
   void writePowerDeriv(ostream &output) const;
   //! Thrown when trying to access an unknown variable by deriv_id
diff --git a/src/DynamicModel.cc b/src/DynamicModel.cc
index aa3637011a652153ea89da95be185ecb9ab223f8..567142d54459fe444659fb0edfb90854ecebf739 100644
--- a/src/DynamicModel.cc
+++ b/src/DynamicModel.cc
@@ -1160,140 +1160,88 @@ DynamicModel::writeDynamicCFile(const string &basename) const
 {
   filesystem::create_directories(basename + "/model/src");
   string filename = basename + "/model/src/dynamic.c";
-  string filename_mex = basename + "/model/src/dynamic_mex.c";
-  ofstream mDynamicModelFile, mDynamicMexFile;
 
   int ntt = temporary_terms_mlv.size() + temporary_terms_derivatives[0].size() + temporary_terms_derivatives[1].size() + temporary_terms_derivatives[2].size() + temporary_terms_derivatives[3].size();
 
-  mDynamicModelFile.open(filename, ios::out | ios::binary);
-  if (!mDynamicModelFile.is_open())
+  ofstream output;
+  output.open(filename, ios::out | ios::binary);
+  if (!output.is_open())
     {
       cerr << "Error: Can't open file " << filename << " for writing" << endl;
       exit(EXIT_FAILURE);
     }
-  mDynamicModelFile << "/*" << endl
-                    << " * " << filename << " : Computes dynamic model for Dynare" << endl
-                    << " *" << endl
-                    << " * Warning : this file is generated automatically by Dynare" << endl
-                    << " *           from model file (.mod)" << endl
-                    << " */" << endl
-                    << "#include <math.h>" << endl;
-
-  mDynamicModelFile << "#include <stdlib.h>" << endl;
-
-  if (external_functions_table.get_total_number_of_unique_model_block_external_functions())
-    // External Matlab function, implies Dynamic function will call mex
-    mDynamicModelFile << R"(#include "mex.h")" << endl;
-
-  mDynamicModelFile << "#define max(a, b) (((a) > (b)) ? (a) : (b))" << endl
-                    << "#define min(a, b) (((a) > (b)) ? (b) : (a))" << endl;
+  output << "/*" << endl
+         << " * " << filename << " : Computes dynamic model for Dynare" << endl
+         << " *" << endl
+         << " * Warning : this file is generated automatically by Dynare" << endl
+         << " *           from model file (.mod)" << endl
+         << " */" << endl
+         << endl
+         << "#include <math.h>" << endl
+         << "#include <stdlib.h>" << endl
+         << R"(#include "mex.h")" << endl
+         << endl
+         << "#define max(a, b) (((a) > (b)) ? (a) : (b))" << endl
+         << "#define min(a, b) (((a) > (b)) ? (b) : (a))" << endl
+         << endl;
 
   // Write function definition if BinaryOpcode::powerDeriv is used
-  writePowerDerivCHeader(mDynamicModelFile);
-
-  mDynamicModelFile << endl;
+  writePowerDeriv(output);
 
-  // Writing the function body
-  writeDynamicModel(mDynamicModelFile, true, false);
+  output << endl;
 
-  mDynamicModelFile << endl;
+  writeDynamicModel(output, true, false);
 
-  writePowerDeriv(mDynamicModelFile);
-  mDynamicModelFile.close();
-
-  mDynamicMexFile.open(filename_mex, ios::out | ios::binary);
-  if (!mDynamicMexFile.is_open())
-    {
-      cerr << "Error: Can't open file " << filename_mex << " for writing" << endl;
-      exit(EXIT_FAILURE);
-    }
+  output << "void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])" << endl
+         << "{" << endl
+         << "  if (nlhs > " << computed_derivs_order + 1 << ")" << endl
+         << R"(    mexErrMsgTxt("Derivatives of higher order than computed have been requested");)" << endl
+         << endl
+         << "  double *y = mxGetPr(prhs[0]);" << endl
+         << "  double *x = mxGetPr(prhs[1]);" << endl
+         << "  double *params = mxGetPr(prhs[2]);" << endl
+         << "  double *steady_state = mxGetPr(prhs[3]);" << endl
+         << "  int it_ = (int) mxGetScalar(prhs[4]) - 1;" << endl
+         << "  int nb_row_x = mxGetM(prhs[1]);" << endl
+         << endl
+         << "  double *T = (double *) malloc(sizeof(double)*" << ntt << ");" << endl
+         << endl
+         << "  if (nlhs >= 1)" << endl
+         << "    {" << endl
+         << "       plhs[0] = mxCreateDoubleMatrix(" << equations.size() << ",1, mxREAL);" << endl
+         << "       double *residual = mxGetPr(plhs[0]);" << endl
+         << "       dynamic_resid_tt(y, x, nb_row_x, params, steady_state, it_, T);" << endl
+         << "       dynamic_resid(y, x, nb_row_x, params, steady_state, it_, T, residual);" << endl
+         << "    }" << endl
+         << endl
+         << "  if (nlhs >= 2)" << endl
+         << "    {" << endl
+         << "       plhs[1] = mxCreateDoubleMatrix(" << equations.size() << ", " << dynJacobianColsNbr << ", mxREAL);" << endl
+         << "       double *g1 = mxGetPr(plhs[1]);" << endl
+         << "       dynamic_g1_tt(y, x, nb_row_x, params, steady_state, it_, T);" << endl
+         << "       dynamic_g1(y, x, nb_row_x, params, steady_state, it_, T, g1);" << endl
+         << "    }" << endl
+         << endl
+         << "  if (nlhs >= 3)" << endl
+         << "    {" << endl
+         << "      plhs[2] = mxCreateDoubleMatrix(" << NNZDerivatives[2] << ", " << 3 << ", mxREAL);" << endl
+         << "      double *v2 = mxGetPr(plhs[2]);" << endl
+         << "      dynamic_g2_tt(y, x, nb_row_x, params, steady_state, it_, T);" << endl
+         << "      dynamic_g2(y, x, nb_row_x, params, steady_state, it_, T, v2);" << endl
+         << "    }" << endl
+         << endl
+         << "  if (nlhs >= 4)" << endl
+         << "    {" << endl
+         << "      plhs[3] = mxCreateDoubleMatrix(" << NNZDerivatives[3] << ", " << 3 << ", mxREAL);" << endl
+         << "      double *v3 = mxGetPr(plhs[3]);" << endl
+         << "      dynamic_g3_tt(y, x, nb_row_x, params, steady_state, it_, T);" << endl
+         << "      dynamic_g3(y, x, nb_row_x, params, steady_state, it_, T, v3);" << endl
+         << "    }" << endl
+         << endl
+         << "  free(T);"
+         << "}" << endl;
 
-  // Writing the gateway routine
-  mDynamicMexFile << "/*" << endl
-                  << " * " << filename_mex << " : The gateway routine used to call the Dynamic function "
-                  << "located in " << filename << endl
-                  << " *" << endl
-                  << " * Warning : this file is generated automatically by Dynare" << endl
-                  << " *           from model file (.mod)" << endl
-                  << endl
-                  << " */" << endl
-                  << endl
-                  << "#include <stdlib.h>" << endl
-                  << R"(#include "mex.h")" << endl
-                  << endl
-                  << "void dynamic_resid_tt(const double *y, const double *x, int nb_row_x, const double *params, const double *steady_state, int it_, double *T);" << endl
-                  << "void dynamic_resid(const double *y, const double *x, int nb_row_x, const double *params, const double *steady_state, int it_, const double *T, double *residual);" << endl
-                  << "void dynamic_g1_tt(const double *y, const double *x, int nb_row_x, const double *params, const double *steady_state, int it_, double *T);" << endl
-                  << "void dynamic_g1(const double *y, const double *x, int nb_row_x, const double *params, const double *steady_state, int it_, const double *T, double *g1);" << endl
-                  << "void dynamic_g2_tt(const double *y, const double *x, int nb_row_x, const double *params, const double *steady_state, int it_, double *T);" << endl
-                  << "void dynamic_g2(const double *y, const double *x, int nb_row_x, const double *params, const double *steady_state, int it_, const double *T, double *v2);" << endl
-                  << "void dynamic_g3_tt(const double *y, const double *x, int nb_row_x, const double *params, const double *steady_state, int it_, double *T);" << endl
-                  << "void dynamic_g3(const double *y, const double *x, int nb_row_x, const double *params, const double *steady_state, int it_, const double *T, double *v3);" << endl
-                  << "void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])" << endl
-                  << "{" << endl
-                  << "  /* Check that no derivatives of higher order than computed are being requested */" << endl
-                  << "  if (nlhs > " << computed_derivs_order + 1 << ")" << endl
-                  << R"(    mexErrMsgTxt("Derivatives of higher order than computed have been requested");)" << endl
-                  << "  /* Create a pointer to the input matrix y. */" << endl
-                  << "  double *y = mxGetPr(prhs[0]);" << endl
-                  << endl
-                  << "  /* Create a pointer to the input matrix x. */" << endl
-                  << "  double *x = mxGetPr(prhs[1]);" << endl
-                  << endl
-                  << "  /* Create a pointer to the input matrix params. */" << endl
-                  << "  double *params = mxGetPr(prhs[2]);" << endl
-                  << endl
-                  << "  /* Create a pointer to the input matrix steady_state. */" << endl
-                  << "  double *steady_state = mxGetPr(prhs[3]);" << endl
-                  << endl
-                  << "  /* Fetch time index */" << endl
-                  << "  int it_ = (int) mxGetScalar(prhs[4]) - 1;" << endl
-                  << endl
-                  << "  /* Gets number of rows of matrix x. */" << endl
-                  << "  int nb_row_x = mxGetM(prhs[1]);" << endl
-                  << endl
-                  << "  double *T = (double *) malloc(sizeof(double)*" << ntt << ");"
-                  << endl
-                  << "  if (nlhs >= 1)" << endl
-                  << "  {" << endl
-                  << "     /* Set the output pointer to the output matrix residual. */" << endl
-                  << "     plhs[0] = mxCreateDoubleMatrix(" << equations.size() << ",1, mxREAL);" << endl
-                  << "     double *residual = mxGetPr(plhs[0]);" << endl
-                  << "     dynamic_resid_tt(y, x, nb_row_x, params, steady_state, it_, T);" << endl
-                  << "     dynamic_resid(y, x, nb_row_x, params, steady_state, it_, T, residual);" << endl
-                  << "  }" << endl
-                  << endl
-                  << "  if (nlhs >= 2)" << endl
-                  << "  {" << endl
-                  << "     /* Set the output pointer to the output matrix g1. */" << endl
-                  << "     plhs[1] = mxCreateDoubleMatrix(" << equations.size() << ", " << dynJacobianColsNbr << ", mxREAL);" << endl
-                  << "     double *g1 = mxGetPr(plhs[1]);" << endl
-                  << "     dynamic_g1_tt(y, x, nb_row_x, params, steady_state, it_, T);" << endl
-                  << "     dynamic_g1(y, x, nb_row_x, params, steady_state, it_, T, g1);" << endl
-                  << "  }" << endl
-                  << endl
-                  << " if (nlhs >= 3)" << endl
-                  << "  {" << endl
-                  << "     /* Set the output pointer to the output matrix v2. */" << endl
-                  << "     plhs[2] = mxCreateDoubleMatrix(" << NNZDerivatives[2] << ", " << 3
-                  << ", mxREAL);" << endl
-                  << "     double *v2 = mxGetPr(plhs[2]);" << endl
-                  << "     dynamic_g2_tt(y, x, nb_row_x, params, steady_state, it_, T);" << endl
-                  << "     dynamic_g2(y, x, nb_row_x, params, steady_state, it_, T, v2);" << endl
-                  << "  }" << endl
-                  << endl
-                  << " if (nlhs >= 4)" << endl
-                  << "  {" << endl
-                  << "     /* Set the output pointer to the output matrix v3. */" << endl
-                  << "     plhs[3] = mxCreateDoubleMatrix(" << NNZDerivatives[3] << ", " << 3 << ", mxREAL);" << endl
-                  << "     double *v3 = mxGetPr(plhs[3]);" << endl
-                  << "     dynamic_g3_tt(y, x, nb_row_x, params, steady_state, it_, T);" << endl
-                  << "     dynamic_g3(y, x, nb_row_x, params, steady_state, it_, T, v3);" << endl
-                  << "  }" << endl
-                  << endl
-                  << " free(T);"
-                  << "}" << endl;
-  mDynamicMexFile.close();
+  output.close();
 }
 
 string
@@ -4541,7 +4489,7 @@ DynamicModel::writeDynamicFile(const string &basename, bool block, bool linear_d
   else if (use_dll)
     {
       writeDynamicCFile(basename);
-      compileDll(basename, "dynamic", mexext, matlabroot, dynareroot);
+      compileMEX(basename, "dynamic", mexext, matlabroot, dynareroot);
     }
   else if (julia)
     writeDynamicJuliaFile(basename);
diff --git a/src/ModelTree.cc b/src/ModelTree.cc
index c64e0330079c08edf5a6d3207687b741b6a1ea09..1dec13509ad36eab636e21813b07507f7bff206e 100644
--- a/src/ModelTree.cc
+++ b/src/ModelTree.cc
@@ -1978,7 +1978,7 @@ ModelTree::matlab_arch(const string &mexext)
 }
 
 void
-ModelTree::compileDll(const string &basename, const string &static_or_dynamic, const string &mexext, const filesystem::path &matlabroot, const filesystem::path &dynareroot) const
+ModelTree::compileMEX(const string &basename, const string &funcname, const string &mexext, const filesystem::path &matlabroot, const filesystem::path &dynareroot) const
 {
   const string opt_flags = "-O3 -g0 --param ira-max-conflict-table-size=1 -fno-forward-propagate -fno-gcse -fno-dce -fno-dse -fno-tree-fre -fno-tree-pre -fno-tree-cselim -fno-tree-dse -fno-tree-dce -fno-tree-pta -fno-gcse-after-reload";
 
@@ -2058,11 +2058,9 @@ ModelTree::compileDll(const string &basename, const string &static_or_dynamic, c
     }
 
   auto model_dir = filesystem::path{basename} / "model" / "src";
-  filesystem::path main_src{model_dir / (static_or_dynamic + ".c")},
-    mex_src{model_dir / (static_or_dynamic + "_mex.c")};
-
+  filesystem::path src{model_dir / (funcname + ".c")};
   filesystem::path mex_dir{"+" + basename};
-  filesystem::path binary{mex_dir / (static_or_dynamic + "." + mexext)};
+  filesystem::path binary{mex_dir / (funcname + "." + mexext)};
 
   ostringstream cmd;
 
@@ -2092,7 +2090,7 @@ ModelTree::compileDll(const string &basename, const string &static_or_dynamic, c
   if (!user_set_add_flags.empty())
     cmd << user_set_add_flags << " ";
 
-  cmd << main_src << " " << mex_src << " -o " << binary << " ";
+  cmd << src << " -o " << binary << " ";
 
   if (user_set_subst_libs.empty())
     cmd << libs;
@@ -2106,7 +2104,7 @@ ModelTree::compileDll(const string &basename, const string &static_or_dynamic, c
   cmd << '"';
 #endif
 
-  cout << "Compiling " << static_or_dynamic << " MEX..." << endl << cmd.str() << endl;
+  cout << "Compiling " << funcname << " MEX..." << endl << cmd.str() << endl;
 
   if (system(cmd.str().c_str()))
     {
diff --git a/src/ModelTree.hh b/src/ModelTree.hh
index adde143464e24cae7f71f405b72bad3961f112f0..77741fa930b8c82cc5cf49605859550b4c8f754d 100644
--- a/src/ModelTree.hh
+++ b/src/ModelTree.hh
@@ -405,8 +405,8 @@ private:
   void copyHelper(const ModelTree &m);
   //! Returns the name of the MATLAB architecture given the extension used for MEX files
   static string matlab_arch(const string &mexext);
-  //! Compiles the MEX file
-  void compileDll(const string &basename, const string &static_or_dynamic, const string &mexext, const filesystem::path &matlabroot, const filesystem::path &dynareroot) const;
+  //! Compiles a MEX file
+  void compileMEX(const string &basename, const string &funcname, const string &mexext, const filesystem::path &matlabroot, const filesystem::path &dynareroot) const;
 
 public:
   ModelTree(SymbolTable &symbol_table_arg,
diff --git a/src/StaticModel.cc b/src/StaticModel.cc
index aa6be689412660844915f279eda58fd108b81f6e..ebfb747d6458fc84deceff3a8938b384532a9016 100644
--- a/src/StaticModel.cc
+++ b/src/StaticModel.cc
@@ -1562,7 +1562,6 @@ StaticModel::writeStaticCFile(const string &basename) const
   // Writing comments and function definition command
   filesystem::create_directories(basename + "/model/src");
   string filename = basename + "/model/src/static.c";
-  string filename_mex = basename + "/model/src/static_mex.c";
 
   int ntt = temporary_terms_mlv.size() + temporary_terms_derivatives[0].size() + temporary_terms_derivatives[1].size() + temporary_terms_derivatives[2].size() + temporary_terms_derivatives[3].size();
 
@@ -1580,74 +1579,33 @@ StaticModel::writeStaticCFile(const string &basename) const
          << " * Warning : this file is generated automatically by Dynare" << endl
          << " *           from model file (.mod)" << endl << endl
          << " */" << endl
-         << "#include <math.h>" << endl;
-
-  output << "#include <stdlib.h>" << endl;
-
-  if (external_functions_table.get_total_number_of_unique_model_block_external_functions())
-    // External Matlab function, implies Static function will call mex
-    output << R"(#include "mex.h")" << endl;
-
-  output << "#define max(a, b) (((a) > (b)) ? (a) : (b))" << endl
-         << "#define min(a, b) (((a) > (b)) ? (b) : (a))" << endl;
+         << endl
+         << "#include <math.h>" << endl
+         << "#include <stdlib.h>" << endl
+         << R"(#include "mex.h")" << endl
+         << endl
+         << "#define max(a, b) (((a) > (b)) ? (a) : (b))" << endl
+         << "#define min(a, b) (((a) > (b)) ? (b) : (a))" << endl
+         << endl;
 
   // Write function definition if BinaryOpcode::powerDeriv is used
-  writePowerDerivCHeader(output);
+  writePowerDeriv(output);
 
   output << endl;
 
-  // Writing the function body
   writeStaticModel(output, true, false);
 
-  output << endl;
-
-  writePowerDeriv(output);
-  output.close();
-
-  output.open(filename_mex, ios::out | ios::binary);
-  if (!output.is_open())
-    {
-      cerr << "ERROR: Can't open file " << filename_mex << " for writing" << endl;
-      exit(EXIT_FAILURE);
-    }
-
-  // Writing the gateway routine
-  output << "/*" << endl
-         << " * " << filename_mex << " : The gateway routine used to call the Static function "
-         << "located in " << filename << endl
-         << " *" << endl
-         << " * Warning : this file is generated automatically by Dynare" << endl
-         << " *           from model file (.mod)" << endl << endl
-         << " */" << endl
-         << endl
-         << "#include <stdlib.h>" << endl
-         << R"(#include "mex.h")" << endl
-         << endl
-         << "void static_resid_tt(const double *y, const double *x, int nb_row_x, const double *params, double *T);" << endl
-         << "void static_resid(const double *y, const double *x, int nb_row_x, const double *params, const double *T, double *residual);" << endl
-         << "void static_g1_tt(const double *y, const double *x, int nb_row_x, const double *params, double *T);" << endl
-         << "void static_g1(const double *y, const double *x, int nb_row_x, const double *params, const double *T, double *g1);" << endl
-         << "void static_g2_tt(const double *y, const double *x, int nb_row_x, const double *params, double *T);" << endl
-         << "void static_g2(const double *y, const double *x, int nb_row_x, const double *params, const double *T, double *v2);" << endl
-         << "void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])" << endl
+  output << "void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])" << endl
          << "{" << endl
-         << "  /* Create a pointer to the input matrix y. */" << endl
          << "  double *y = mxGetPr(prhs[0]);" << endl
-         << endl
-         << "  /* Create a pointer to the input matrix x. */" << endl
          << "  double *x = mxGetPr(prhs[1]);" << endl
-         << endl
-         << "  /* Create a pointer to the input matrix params. */" << endl
          << "  double *params = mxGetPr(prhs[2]);" << endl
-         << endl
-         << "  /* Gets number of rows of matrix x. */" << endl
          << "  int nb_row_x = mxGetM(prhs[1]);" << endl
          << endl
-         << "  double *T = (double *) malloc(sizeof(double)*" << ntt << ");"
+         << "  double *T = (double *) malloc(sizeof(double)*" << ntt << ");" << endl
          << endl
          << "  if (nlhs >= 1)" << endl
          << "    {" << endl
-         << "      /* Set the output pointer to the output matrix residual. */" << endl
          << "      plhs[0] = mxCreateDoubleMatrix(" << equations.size() << ",1, mxREAL);" << endl
          << "      double *residual = mxGetPr(plhs[0]);" << endl
          << "      static_resid_tt(y, x, nb_row_x, params, T);" << endl
@@ -1656,7 +1614,6 @@ StaticModel::writeStaticCFile(const string &basename) const
          << endl
          << "  if (nlhs >= 2)" << endl
          << "    {" << endl
-         << "      /* Set the output pointer to the output matrix g1. */" << endl
          << "      plhs[1] = mxCreateDoubleMatrix(" << equations.size() << ", " << symbol_table.endo_nbr() << ", mxREAL);" << endl
          << "      double *g1 = mxGetPr(plhs[1]);" << endl
          << "      static_g1_tt(y, x, nb_row_x, params, T);" << endl
@@ -1665,9 +1622,7 @@ StaticModel::writeStaticCFile(const string &basename) const
          << endl
          << "  if (nlhs >= 3)" << endl
          << "    {" << endl
-         << "      /* Set the output pointer to the output matrix v2. */" << endl
-         << "      plhs[2] = mxCreateDoubleMatrix(" << NNZDerivatives[2] << ", " << 3
-         << ", mxREAL);" << endl
+         << "      plhs[2] = mxCreateDoubleMatrix(" << NNZDerivatives[2] << ", " << 3 << ", mxREAL);" << endl
          << "      double *v2 = mxGetPr(plhs[2]);" << endl
          << "      static_g2_tt(y, x, nb_row_x, params, T);" << endl
          << "      static_g2(y, x, nb_row_x, params, T, v2);" << endl
@@ -1675,6 +1630,7 @@ StaticModel::writeStaticCFile(const string &basename) const
          << endl
          << "  free(T);" << endl
          << "}" << endl;
+
   output.close();
 }
 
@@ -1699,7 +1655,7 @@ StaticModel::writeStaticFile(const string &basename, bool block, bool bytecode,
   else if (use_dll)
     {
       writeStaticCFile(basename);
-      compileDll(basename, "static", mexext, matlabroot, dynareroot);
+      compileMEX(basename, "static", mexext, matlabroot, dynareroot);
     }
   else if (julia)
     writeStaticJuliaFile(basename);