diff --git a/mex/sources/korderpert/src/dynamic_dll.cpp b/mex/sources/korderpert/src/dynamic_dll.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..9eec4b53da76cd11e7f2dc12475bc30ca58a31f9
--- /dev/null
+++ b/mex/sources/korderpert/src/dynamic_dll.cpp
@@ -0,0 +1,244 @@
+/*
+ * Copyright (C) 2008-2009 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 <http://www.gnu.org/licenses/>.
+ */
+
+
+#include "k_ord_dynare.h"
+#include "dynamic_dll.h"
+#include "math.h"
+#include <cstring>
+
+//////////////////////////////////////////////////////
+// Convert Matlab Dynare endo and exo names array to C type array of string pointers
+// Poblem is that Matlab mx function returns a long string concatenated by columns rather than rows
+// hence a rather low level approach is needed
+///////////////////////////////////////////////////////
+const char **
+DynareMxArrayToString(const mxArray *mxFldp, const int len, const int width)
+{
+  char *cNamesCharStr = mxArrayToString(mxFldp);
+  const char **ret = DynareMxArrayToString(cNamesCharStr, len, width);
+  return ret;
+}
+
+const char **
+DynareMxArrayToString(const char *cNamesCharStr, const int len, const int width)
+{
+  char cNamesMX[len][width+1]; //
+#ifdef DEBUG
+  mexPrintf("loop DynareMxArrayToString cNamesCharStr = %s \n", cNamesCharStr);
+#endif
+  for (int i = 0; i < width; i++)
+    {
+      for (int j = 0; j < len; j++)
+        {
+          // Allow alphanumeric and underscores "_" only:
+          if (isalnum(cNamesCharStr[j+i*len]) || ('_' == cNamesCharStr[j+i*len]))
+            {
+              cNamesMX[j][i] = cNamesCharStr[j+i*len];
+            }
+          else cNamesMX[j][i] = '\0';
+        }
+    }
+  const char **ret = (const char **) mxCalloc(len, sizeof(char *));
+  for (int j = 0; j < len; j++)
+    {
+      cNamesMX[j][width] = '\0';
+#ifdef DEBUG
+      //				mexPrintf("String [%d]= %s \n", j, cNamesMX[j]);
+#endif
+      char *token = (char *) mxCalloc(strlen(cNamesMX[j])+1, sizeof(char));
+      strcpy(token, cNamesMX[j]);
+      ret[j] = token;
+#ifdef DEBUG
+      mexPrintf("ret [%d]= %s \n", j, ret[j]);
+#endif
+    }
+  return ret;
+}
+
+/***********************************
+ * Members of DynamicModelDLL for handling loading and calling
+ * <model>_dynamic () function
+ **************************************/
+DynamicModelDLL::DynamicModelDLL(const char *modName, const int y_length, const int j_cols,
+                                 const int n_max_lag, const int n_exog, const char *sExt)
+  : length(y_length), jcols(j_cols), nMax_lag(n_max_lag), nExog(n_exog)
+{
+  char fName[MAX_MODEL_NAME];
+  strcpy(fName, modName);
+  using namespace std;
+  strcat(fName, "_dynamic");
+#ifdef DEBUG
+  mexPrintf("MexPrintf: Call Load run DLL %s .\n", fName);
+#endif
+
+  try
+    {
+#ifdef WINDOWS
+      if (sExt == NULL) sExt = (".dll");
+      HINSTANCE dynamicHinstance;
+      //		dynamicHinstance=::LoadLibraryEx(strcat(fNname,"_.dll"),NULL,DONT_RESOLVE_DLL_REFERENCES);//sExt); //"_.dll");
+      dynamicHinstance = ::LoadLibrary(strcat(fName, sExt)); //.dll); //"_.dll");
+      if (dynamicHinstance == NULL)
+        throw 1;                  //alt: return;
+      //		(DynamicFn*)	typedef void * (__stdcall *DynamicFn)();
+# ifdef DEBUG
+      mexPrintf("MexPrintf: Call GetProcAddress  %s .\n", fName);
+# endif
+      Dynamic = (DynamicFn *) ::GetProcAddress(dynamicHinstance, "Dynamic");
+
+#else // __linux__
+      if (sExt == NULL) sExt = (".so");
+      dynamicHinstance = dlopen(strcat(fName, sExt), RTLD_NOW);
+      if ((dynamicHinstance == NULL) || dlerror())
+        {
+          cerr << dlerror() << endl;
+          mexPrintf("MexPrintf:Error loading DLL: %s", dlerror);
+          throw 1;
+        }
+      Dynamic = (DynamicFn) dlsym(dynamicHinstance, "Dynamic");
+      if ((Dynamic  == NULL) || dlerror())
+        {
+          cerr << dlerror() << endl;
+          mexPrintf("MexPrintf:Error finding DLL function: %s", dlerror);
+          throw 2;
+        }
+#endif
+
+    }
+  catch (int i)
+    {
+      mexPrintf("MexPrintf: error in Load and run DLL %s , %d.\n", fName, i);
+      mexErrMsgTxt("Err: An error in Load and run DLL  .\n");
+      return;
+
+    }
+  catch (...)
+    {
+      mexPrintf("MexPrintf: Unknown error in Call MATLAB function %s.\n", fName);
+      mexErrMsgTxt("Err: Unknown error in Load and run DLL  .\n");
+      return;
+    }
+}
+
+// close DLL: If the referenced object was successfully closed,
+// close() returns 0, non 0 otherwise
+int
+DynamicModelDLL::close()
+{
+#ifdef WINDOWS
+  // MS FreeLibrary returns non 0 if OK, 0 if fails.
+  bool rb = FreeLibrary(dynamicHinstance);
+  if (rb)
+    return 0;
+  else
+    return 1;
+#else // linux
+  //If OK, dlclose() returns 0, non 0 otherwise
+  return dlclose(dynamicHinstance);
+#endif
+};
+
+void
+DynamicModelDLL::eval(const Vector &y, const TwoDMatrix &x, const  Vector *modParams,
+                      int it_, Vector &residual, TwoDMatrix *g1, TwoDMatrix *g2)
+{
+
+  double  *dresidual, *dg1 = NULL, *dg2 = NULL;
+  //int length=y.length(); // not!
+  if ((jcols-nExog) != y.length())
+    {
+      // throw DLL Error
+      mexPrintf(" DLL Error: (jcols-nExog)!=ys.length() \n");
+      return;
+    }
+  if (residual.length() < length)  // dummy or insufficient
+    {
+      Vector *tempv = new Vector(length);
+      residual = *tempv;
+      delete tempv;
+      residual.zeros();
+    }
+  if (g1 != NULL)
+    {
+      if (g1->nrows() != length)  // dummy
+        {
+          delete g1;
+          g1 =     new TwoDMatrix(length, jcols); // and get a new one
+          g1->zeros();
+        }
+      dg1 = const_cast<double *>(g1->base());
+    }
+  if (g2 != NULL)
+    {
+      dg2 = const_cast<double *>(g2->base());
+    }
+  dresidual = const_cast<double *>(residual.base());
+  double *dy = const_cast<double *>(y.base());
+  double *dx = const_cast<double *>(x.base());
+  double *dbParams = const_cast<double *>(modParams->base());
+#ifdef DEBUG
+  mexPrintf(" try eval Dynamic with ne g1: cols=%d , rows=%d\n",
+            g1->ncols(), g1->nrows());
+  for (int i = 0; i < modParams->length(); i++)
+    {
+      mexPrintf("k_ord_perturbation: Params[%d]= %g.\n", i, (*modParams)[i]);
+    }
+  for (int i = 0; i < jcols-nExog; i++)
+    {
+      mexPrintf("k_ord_perturbation: Ys[%d]= %g.\n", i, dy[i]);
+    }
+  mexPrintf("k_order_perturbation: call <model> Dynamic dParams= %g ,  , dy = %g dx = %f .\n",
+            dbParams[0], dy[0], dx[0]);
+
+#endif
+  try
+    {
+      Dynamic(dy, dx, nExog, dbParams, it_, dresidual, dg1, dg2);
+    }
+  catch (...)
+    {
+      mexPrintf("MexPrintf: error in run Dynamic DLL \n");
+    }
+};
+
+void
+DynamicModelDLL::eval(const Vector &y, const TwoDMatrix &x, const Vector *modParams,
+                      Vector &residual, TwoDMatrix *g1, TwoDMatrix *g2)
+{
+
+  eval(y, x, modParams, nMax_lag, residual, g1, g2);
+};
+
+void
+DynamicModelDLL::eval(const Vector &y, const Vector &x, const Vector *modParams,
+                      Vector &residual, TwoDMatrix *g1, TwoDMatrix *g2)
+{
+
+  /** ignore given exogens and create new 2D x matrix since
+   * when calling <model>_dynamic(z,x,params,it_) x must be equal to
+   * zeros(M_.maximum_lag+1,M_.exo_nbr)
+   **/
+  TwoDMatrix &mx = *(new TwoDMatrix(nMax_lag+1, nExog));
+  mx.zeros(); // initialise shocks to 0s
+
+  eval(y, mx, modParams, nMax_lag, residual, g1, g2);
+};
+
+
diff --git a/mex/sources/korderpert/src/k_order_perturbation.h b/mex/sources/korderpert/src/dynamic_dll.h
similarity index 100%
rename from mex/sources/korderpert/src/k_order_perturbation.h
rename to mex/sources/korderpert/src/dynamic_dll.h
diff --git a/mex/sources/korderpert/src/k_ord_dynare.cpp b/mex/sources/korderpert/src/k_ord_dynare.cpp
index b0a7a5bf16d910d6072c4a63109212b90d6f8820..1206305273d0162e745f4606d530d39a3493a5df 100644
--- a/mex/sources/korderpert/src/k_ord_dynare.cpp
+++ b/mex/sources/korderpert/src/k_ord_dynare.cpp
@@ -22,6 +22,7 @@
 #include <vector>
 #include "first_order.h"
 #include "k_ord_dynare.h"
+#include "dynamic_dll.h"
 
 #include "mex.h"
 
@@ -40,12 +41,12 @@ class KordpJacobian;
 KordpDynare::KordpDynare(const char **endo,  int num_endo,
                          const char **exo, int nexog, int npar, //const char** par,
                          Vector *ysteady, TwoDMatrix *vcov, Vector *inParams, int nstat,
-                         int npred, int nforw, int nboth, const int jcols, const Vector *nnzd, 
+                         int npred, int nforw, int nboth, const int jcols, const Vector *nnzd, 
                          const int nsteps, int norder, //const char* modName,
                          Journal &jr, DynamicModelDLL &dynamicDLL, double sstol, 
                          const vector<int> *var_order, const TwoDMatrix *llincidence, double criterium)
   : nStat(nstat), nBoth(nboth), nPred(npred), nForw(nforw), nExog(nexog), nPar(npar),
-  nYs(npred + nboth), nYss(nboth + nforw), nY(num_endo), nJcols(jcols), NNZD(nnzd), nSteps(nsteps), 
+  nYs(npred + nboth), nYss(nboth + nforw), nY(num_endo), nJcols(jcols), NNZD(nnzd), nSteps(nsteps), 
   nOrder(norder), journal(jr),  dynamicDLL(dynamicDLL), ySteady(ysteady), vCov(vcov), params(inParams),
   md(1), dnl(NULL), denl(NULL), dsnl(NULL), ss_tol(sstol), varOrder(var_order),
   ll_Incidence(llincidence), qz_criterium(criterium)
@@ -94,8 +95,8 @@ KordpDynare::KordpDynare(const char **endo,  int num_endo,
 KordpDynare::KordpDynare(const KordpDynare &dynare)
   : nStat(dynare.nStat), nBoth(dynare.nBoth),       nPred(dynare.nPred),
   nForw(dynare.nForw), nExog(dynare.nExog),  nPar(dynare.nPar),
-  nYs(dynare.nYs), nYss(dynare.nYss), nY(dynare.nY), nJcols(dynare.nJcols),
-  NNZD(dynare.NNZD), nSteps(dynare.nSteps), nOrder(dynare.nOrder), 
+  nYs(dynare.nYs), nYss(dynare.nYss), nY(dynare.nY), nJcols(dynare.nJcols),
+  NNZD(dynare.NNZD), nSteps(dynare.nSteps), nOrder(dynare.nOrder), 
   journal(dynare.journal), dynamicDLL(dynare.dynamicDLL), //modName(dynare.modName),
   ySteady(NULL), params(NULL), vCov(NULL), md(dynare.md),
   dnl(NULL), denl(NULL), dsnl(NULL), ss_tol(dynare.ss_tol),
@@ -203,13 +204,13 @@ KordpDynare::calcDerivatives(const Vector &yy, const Vector &xx)
   // Hessian TwoDMatrix *g2;
   if (nOrder > 1)
     {
-      //g2 = new TwoDMatrix(nY, nJcols*nJcols); // generate g2 for Hessian
-      // generate g2 space for sparse Hessian 3x NNZH = 3x NNZD[1]
+      //g2 = new TwoDMatrix(nY, nJcols*nJcols); // generate g2 for Hessian
+      // generate g2 space for sparse Hessian 3x NNZH = 3x NNZD[1]
       g2 = new TwoDMatrix((int) (*NNZD)[1],3);
-      g2->zeros();
-#ifdef DEBUG
-      mexPrintf(" g2 cols %d rows %d \n", g2->numCols(), g2->numRows());
-      //g2->print();
+      g2->zeros();
+#ifdef DEBUG
+      mexPrintf(" g2 cols %d rows %d \n", g2->numCols(), g2->numRows());
+      //g2->print();
 #endif
     }
   Vector &out = *(new Vector(nY));
@@ -249,12 +250,12 @@ KordpDynare::calcDerivatives(const Vector &yy, const Vector &xx)
   populateDerivativesContainer(g1, 1, JacobianIndices);
   if (nOrder > 1)
     {
-      //		ReorderBlocks(g2,JacobianIndices);
-#ifdef DEBUG
-      mexPrintf(" post dll g2 cols %d rows %d \n", g2->numCols(), g2->numRows());
-      for (int ii=0;ii<g2->numRows(); ii++)
-        mexPrintf(" g2[%d]:  %d %d %f \n", ii, (int)g2->get(ii,0),(int)g2->get(ii,1),g2->get(ii,2));
-      //g2->print();
+      //		ReorderBlocks(g2,JacobianIndices);
+#ifdef DEBUG
+      mexPrintf(" post dll g2 cols %d rows %d \n", g2->numCols(), g2->numRows());
+      for (int ii=0;ii<g2->numRows(); ii++)
+        mexPrintf(" g2[%d]:  %d %d %f \n", ii, (int)g2->get(ii,0),(int)g2->get(ii,1),g2->get(ii,2));
+      //g2->print();
 #endif
       populateDerivativesContainer(g2, 2, JacobianIndices);
     }
@@ -336,7 +337,7 @@ KordpDynare::populateDerivativesContainer(TwoDMatrix *g, int ord, const vector<i
     }
   else if (ord == 2) 
     {
-    int nJcols1 = nJcols-nExog;
+    int nJcols1 = nJcols-nExog;
     vector<int> revOrder(nJcols1);
     for (int i = 0; i < nJcols1; i++)
       revOrder[(*vOrder)[i]] = i;
@@ -344,7 +345,7 @@ KordpDynare::populateDerivativesContainer(TwoDMatrix *g, int ord, const vector<i
     	{
 	    int j = (int)g->get(i,0)-1; // hessian indices start with 1
 	    int i1 = (int)g->get(i,1) -1;
-	    int s0 = (int)floor(i1/nJcols);
+	    int s0 = (int)floor(i1/nJcols);
       int s1 = i1- (nJcols*s0);
 	    if (s0 < nJcols1)
 	      s[0] = revOrder[s0];
@@ -355,20 +356,20 @@ KordpDynare::populateDerivativesContainer(TwoDMatrix *g, int ord, const vector<i
 	    else
 	      s[1] = s1;
 	    if (s[1] >= s[0])
-	      {
+	      {
 	      double x = g->get(i,2);
     		mdTi->insert(s, j, x);
-#ifdef DEBUG
-        mexPrintf(" s[0]=%d  s[1]=%d j=%d  x=%f \n", s[0],s[1],j,x);
-	      s.print();
-	      std::cout << s[0] << " " << s[1] << "\n";
-#endif
+#ifdef DEBUG
+        mexPrintf(" s[0]=%d  s[1]=%d j=%d  x=%f \n", s[0],s[1],j,x);
+	      s.print();
+	      std::cout << s[0] << " " << s[1] << "\n";
+#endif
 	      }
     	}
-    }
+    }
 #ifdef DEBUG
-  mexPrintf("k_ord_dynare.cpp: END populate FSSparseTensor in calcDerivatives: ord=%d \n",ord);
-  mdTi->print();
+  mexPrintf("k_ord_dynare.cpp: END populate FSSparseTensor in calcDerivatives: ord=%d \n",ord);
+  mdTi->print();
 #endif
   // md container
   //md.clear();// this is to be used only for 1st order!!
diff --git a/mex/sources/korderpert/src/k_ord_dynare.h b/mex/sources/korderpert/src/k_ord_dynare.h
index 2bf505c38cf608e7b566fff953111357dc902af5..f57d0668fdbdef2e46849da6a4081baaa30380e1 100644
--- a/mex/sources/korderpert/src/k_ord_dynare.h
+++ b/mex/sources/korderpert/src/k_ord_dynare.h
@@ -37,7 +37,6 @@
 #include "kord_exception.h"
 #include "nlsolve.h"
 #include "approximation.h"
-#include "k_order_perturbation.h"
 
 class KordpDynare;
 
@@ -123,8 +122,8 @@ class KordpDynare : public DynamicModel
   const int nYs; // ={npred + nboth ; }
   const int nYss; // nyss ={ nboth + nforw ; }
   const int nY;  // = num_endo={ nstat + npred + nboth + nforw ; }
-  const int nJcols; // no of jacobian columns= nExog+nEndo+nsPred+nsForw
-  const Vector *NNZD;  //the total number of non-zero derivative elements
+  const int nJcols; // no of jacobian columns= nExog+nEndo+nsPred+nsForw
+  const Vector *NNZD;  //the total number of non-zero derivative elements
                       //  where hessian is 2nd : NZZD(order=2)
   const int nSteps;
   const int nOrder;
@@ -147,10 +146,10 @@ public:
   KordpDynare(const char **endo, int num_endo,
               const char **exo, int num_exo, int num_par, //const char** par,
               Vector *ySteady, TwoDMatrix *vCov, Vector *params, int nstat, int nPred,
-              int nforw, int nboth, const int nJcols, const Vector *NNZD, 
+              int nforw, int nboth, const int nJcols, const Vector *NNZD, 
               const int nSteps, const int ord,    //const char* modName,
               Journal &jr, DynamicModelDLL &dynamicDLL, double sstol,
-              const vector<int> *varOrder, const TwoDMatrix *ll_Incidence, 
+              const vector<int> *varOrder, const TwoDMatrix *ll_Incidence, 
               double qz_criterium);
 
   /** Makes a deep copy of the object. */
diff --git a/mex/sources/korderpert/src/k_order_perturbation.cpp b/mex/sources/korderpert/src/k_order_perturbation.cpp
index 12794bbac2e82bbe8d63cca1c08cfa32866ce345..d9f67d72c59b2ff3b689c5b8c75a1048c62f9d1e 100644
--- a/mex/sources/korderpert/src/k_order_perturbation.cpp
+++ b/mex/sources/korderpert/src/k_order_perturbation.cpp
@@ -39,6 +39,7 @@
  **********************************************************/
 
 #include "k_ord_dynare.h"
+#include "dynamic_dll.h"
 #include "math.h"
 #include <cstring>
 
@@ -82,6 +83,7 @@ CK_order_perturbation::CK_order_perturbation()
 
 #endif // _MSC_VER && WINDOWS
 
+#ifdef MATLAB_MEX_FILE  // exclude mexFunction for other applications
 extern "C" {
 
   // mexFunction: Matlab Inerface point and the main application driver
@@ -176,8 +178,8 @@ extern "C" {
     const int nPar = (int) mxGetScalar(mxFldp);
     // it_ should be set to M_.maximum_lag
     mxFldp = mxGetField(M_, 0, "maximum_lag");
-    const int nMax_lag = (int) mxGetScalar(mxFldp);
-
+    const int nMax_lag = (int) mxGetScalar(mxFldp);
+
     nPred -= nBoth; // correct nPred for nBoth.
 
     mxFldp      = mxGetField(dr, 0, "order_var");
@@ -221,15 +223,15 @@ extern "C" {
           }
       }
 #endif
-
-    //get NNZH =NNZD(2) = the total number of non-zero Hessian elements 
-    mxFldp = mxGetField(M_, 0, "NNZDerivatives");
-    dparams = (double *) mxGetData(mxFldp);
-    Vector *NNZD =  new Vector (dparams, (int) mxGetM(mxFldp));
-#ifdef DEBUG
-    mexPrintf("NNZH=%d, \n", (int) (*NNZD)[1]);
-#endif
-
+
+    //get NNZH =NNZD(2) = the total number of non-zero Hessian elements 
+    mxFldp = mxGetField(M_, 0, "NNZDerivatives");
+    dparams = (double *) mxGetData(mxFldp);
+    Vector *NNZD =  new Vector (dparams, (int) mxGetM(mxFldp));
+#ifdef DEBUG
+    mexPrintf("NNZH=%d, \n", (int) (*NNZD)[1]);
+#endif
+
     const int jcols = nExog+nEndo+nsPred+nsForw; // Num of Jacobian columns
     mexPrintf("k_order_perturbation: jcols= %d .\n", jcols);
 
@@ -308,7 +310,7 @@ extern "C" {
         // make KordpDynare object
         KordpDynare dynare(endoNamesMX,  nEndo, exoNamesMX,  nExog, nPar, // paramNames,
                            ySteady, vCov, modParams, nStat, nPred, nForw, nBoth,
-                           jcols, NNZD, nSteps, kOrder, journal, dynamicDLL, 
+                           jcols, NNZD, nSteps, kOrder, journal, dynamicDLL, 
                            sstol, var_order_vp, llincidence, qz_criterium);
 
         // construct main K-order approximation class
@@ -469,222 +471,10 @@ extern "C" {
       }  //catch
   }; // end of mexFunction()
 }; // end of extern C
+#endif // ifdef MATLAB_MEX_FILE  to exclude mexFunction for other applications
 
-//////////////////////////////////////////////////////
-// Convert Matlab Dynare endo and exo names array to C type array of string pointers
-// Poblem is that Matlab mx function returns a long string concatenated by columns rather than rows
-// hence a rather low level approach is needed
-///////////////////////////////////////////////////////
-const char **
-DynareMxArrayToString(const mxArray *mxFldp, const int len, const int width)
-{
-  char *cNamesCharStr = mxArrayToString(mxFldp);
-  const char **ret = DynareMxArrayToString(cNamesCharStr, len, width);
-  return ret;
-}
-
-const char **
-DynareMxArrayToString(const char *cNamesCharStr, const int len, const int width)
-{
-  char cNamesMX[len][width+1]; //
-#ifdef DEBUG
-  mexPrintf("loop DynareMxArrayToString cNamesCharStr = %s \n", cNamesCharStr);
-#endif
-  for (int i = 0; i < width; i++)
-    {
-      for (int j = 0; j < len; j++)
-        {
-          // Allow alphanumeric and underscores "_" only:
-          if (isalnum(cNamesCharStr[j+i*len]) || ('_' == cNamesCharStr[j+i*len]))
-            {
-              cNamesMX[j][i] = cNamesCharStr[j+i*len];
-            }
-          else cNamesMX[j][i] = '\0';
-        }
-    }
-  const char **ret = (const char **) mxCalloc(len, sizeof(char *));
-  for (int j = 0; j < len; j++)
-    {
-      cNamesMX[j][width] = '\0';
-#ifdef DEBUG
-      //				mexPrintf("String [%d]= %s \n", j, cNamesMX[j]);
-#endif
-      char *token = (char *) mxCalloc(strlen(cNamesMX[j])+1, sizeof(char));
-      strcpy(token, cNamesMX[j]);
-      ret[j] = token;
-#ifdef DEBUG
-      mexPrintf("ret [%d]= %s \n", j, ret[j]);
-#endif
-    }
-  return ret;
-}
-
-/***********************************
- * Members of DynamicModelDLL for handling loading and calling
- * <model>_dynamic () function
- **************************************/
-DynamicModelDLL::DynamicModelDLL(const char *modName, const int y_length, const int j_cols,
-                                 const int n_max_lag, const int n_exog, const char *sExt)
-  : length(y_length), jcols(j_cols), nMax_lag(n_max_lag), nExog(n_exog)
-{
-  char fName[MAX_MODEL_NAME];
-  strcpy(fName, modName);
-  using namespace std;
-  strcat(fName, "_dynamic");
-#ifdef DEBUG
-  mexPrintf("MexPrintf: Call Load run DLL %s .\n", fName);
-#endif
-
-  try
-    {
-#ifdef WINDOWS
-      if (sExt == NULL) sExt = (".dll");
-      HINSTANCE dynamicHinstance;
-      //		dynamicHinstance=::LoadLibraryEx(strcat(fNname,"_.dll"),NULL,DONT_RESOLVE_DLL_REFERENCES);//sExt); //"_.dll");
-      dynamicHinstance = ::LoadLibrary(strcat(fName, sExt)); //.dll); //"_.dll");
-      if (dynamicHinstance == NULL)
-        throw 1;                  //alt: return;
-      //		(DynamicFn*)	typedef void * (__stdcall *DynamicFn)();
-# ifdef DEBUG
-      mexPrintf("MexPrintf: Call GetProcAddress  %s .\n", fName);
-# endif
-      Dynamic = (DynamicFn *) ::GetProcAddress(dynamicHinstance, "Dynamic");
-
-#else // __linux__
-      if (sExt == NULL) sExt = (".so");
-      dynamicHinstance = dlopen(strcat(fName, sExt), RTLD_NOW);
-      if ((dynamicHinstance == NULL) || dlerror())
-        {
-          cerr << dlerror() << endl;
-          mexPrintf("MexPrintf:Error loading DLL: %s", dlerror);
-          throw 1;
-        }
-      Dynamic = (DynamicFn) dlsym(dynamicHinstance, "Dynamic");
-      if ((Dynamic  == NULL) || dlerror())
-        {
-          cerr << dlerror() << endl;
-          mexPrintf("MexPrintf:Error finding DLL function: %s", dlerror);
-          throw 2;
-        }
-#endif
-
-    }
-  catch (int i)
-    {
-      mexPrintf("MexPrintf: error in Load and run DLL %s , %d.\n", fName, i);
-      mexErrMsgTxt("Err: An error in Load and run DLL  .\n");
-      return;
-
-    }
-  catch (...)
-    {
-      mexPrintf("MexPrintf: Unknown error in Call MATLAB function %s.\n", fName);
-      mexErrMsgTxt("Err: Unknown error in Load and run DLL  .\n");
-      return;
-    }
-}
-
-// close DLL: If the referenced object was successfully closed,
-// close() returns 0, non 0 otherwise
-int
-DynamicModelDLL::close()
-{
-#ifdef WINDOWS
-  // MS FreeLibrary returns non 0 if OK, 0 if fails.
-  bool rb = FreeLibrary(dynamicHinstance);
-  if (rb)
-    return 0;
-  else
-    return 1;
-#else // linux
-  //If OK, dlclose() returns 0, non 0 otherwise
-  return dlclose(dynamicHinstance);
-#endif
-};
-
-void
-DynamicModelDLL::eval(const Vector &y, const TwoDMatrix &x, const  Vector *modParams,
-                      int it_, Vector &residual, TwoDMatrix *g1, TwoDMatrix *g2)
-{
-
-  double  *dresidual, *dg1 = NULL, *dg2 = NULL;
-  //int length=y.length(); // not!
-  if ((jcols-nExog) != y.length())
-    {
-      // throw DLL Error
-      mexPrintf(" DLL Error: (jcols-nExog)!=ys.length() \n");
-      return;
-    }
-  if (residual.length() < length)  // dummy or insufficient
-    {
-      Vector *tempv = new Vector(length);
-      residual = *tempv;
-      delete tempv;
-      residual.zeros();
-    }
-  if (g1 != NULL)
-    {
-      if (g1->nrows() != length)  // dummy
-        {
-          delete g1;
-          g1 =     new TwoDMatrix(length, jcols); // and get a new one
-          g1->zeros();
-        }
-      dg1 = const_cast<double *>(g1->base());
-    }
-  if (g2 != NULL)
-    {
-      dg2 = const_cast<double *>(g2->base());
-    }
-  dresidual = const_cast<double *>(residual.base());
-  double *dy = const_cast<double *>(y.base());
-  double *dx = const_cast<double *>(x.base());
-  double *dbParams = const_cast<double *>(modParams->base());
-#ifdef DEBUG
-  mexPrintf(" try eval Dynamic with ne g1: cols=%d , rows=%d\n",
-            g1->ncols(), g1->nrows());
-  for (int i = 0; i < modParams->length(); i++)
-    {
-      mexPrintf("k_ord_perturbation: Params[%d]= %g.\n", i, (*modParams)[i]);
-    }
-  for (int i = 0; i < jcols-nExog; i++)
-    {
-      mexPrintf("k_ord_perturbation: Ys[%d]= %g.\n", i, dy[i]);
-    }
-  mexPrintf("k_order_perturbation: call <model> Dynamic dParams= %g ,  , dy = %g dx = %f .\n",
-            dbParams[0], dy[0], dx[0]);
-
-#endif
-  try
-    {
-      Dynamic(dy, dx, nExog, dbParams, it_, dresidual, dg1, dg2);
-    }
-  catch (...)
-    {
-      mexPrintf("MexPrintf: error in run Dynamic DLL \n");
-    }
-};
-
-void
-DynamicModelDLL::eval(const Vector &y, const TwoDMatrix &x, const Vector *modParams,
-                      Vector &residual, TwoDMatrix *g1, TwoDMatrix *g2)
-{
 
-  eval(y, x, modParams, nMax_lag, residual, g1, g2);
-};
 
-void
-DynamicModelDLL::eval(const Vector &y, const Vector &x, const Vector *modParams,
-                      Vector &residual, TwoDMatrix *g1, TwoDMatrix *g2)
-{
 
-  /** ignore given exogens and create new 2D x matrix since
-   * when calling <model>_dynamic(z,x,params,it_) x must be equal to
-   * zeros(M_.maximum_lag+1,M_.exo_nbr)
-   **/
-  TwoDMatrix &mx = *(new TwoDMatrix(nMax_lag+1, nExog));
-  mx.zeros(); // initialise shocks to 0s
 
-  eval(y, mx, modParams, nMax_lag, residual, g1, g2);
-};
 
diff --git a/mex/sources/korderpert/src/k_order_test_main.cpp b/mex/sources/korderpert/src/k_order_test_main.cpp
index 57d4fce5dcd60b3494cedfc1368d35ed79e27f02..de438b2ac77bc5db77becead1d8a751e1b36ddff 100644
--- a/mex/sources/korderpert/src/k_order_test_main.cpp
+++ b/mex/sources/korderpert/src/k_order_test_main.cpp
@@ -26,6 +26,7 @@
 
 //#include "stdafx.h"
 #include "k_ord_dynare.h"
+#include "dynamic_dll.h"
 
 int
 main(int argc, char *argv[])
@@ -79,9 +80,9 @@ main(int argc, char *argv[])
   };
   const int nSteady = 16; //27 //31;//29, 16 (int)mxGetM(mxFldp);
   Vector *ySteady =  new Vector(dYSparams, nSteady);
-
-  double nnzd[3]={ 77,217,0};
-  const Vector *NNZD =  new Vector(nnzd, 3);
+
+  double nnzd[3]={ 77,217,0};
+  const Vector *NNZD =  new Vector(nnzd, 3);
 
   //mxFldp = mxGetField(dr, 0,"nstatic" );
   const int nStat = 7; //(int)mxGetScalar(mxFldp);
diff --git a/mex/sources/korderpert/src/linkDLL.bat b/mex/sources/korderpert/src/linkDLL.bat
index 7b528611871b608474c8a163909fc06bcc7413ea..b3ded88ed0e0b179e85574ddd4e15d216fe9eb87 100644
--- a/mex/sources/korderpert/src/linkDLL.bat
+++ b/mex/sources/korderpert/src/linkDLL.bat
@@ -1 +1 @@
-gcc -DMATLAB_MEX_FILE -mthreads -shared -DWINDOWS -DPOSIX_THREADS -I"c:/Program Files"/MATLAB_SV71/extern/include -I Dynare_pp/src -I Dynare_pp/kord -I Dynare_pp/tl/cc -I Dynare_pp/utils/cc/ -I Dynare_pp/sylv/cc/ -I"f:/Pthreads/Pre-built.2/include" -I"f:/mingw/include" -o k_order_perturbation.dll k_ord_dynare.o k_order_perturbation.o nlsolve.o "dynare_pp/extern/matlab"/dynarelib.a -Wl,-L"c:/Program Files"/MATLAB_SV71/extern/lib/win32/microsoft/ -Wl,-llibmex -Wl,-llibmx -Wl,-llibmwlapack -Wl,-llibdflapack -lg2c -lmingw32  -lstdc++ -L"f:/Pthreads/Pre-built.2/lib" -lpthreadGC2
+gcc -DMATLAB_MEX_FILE -mthreads -shared -DWINDOWS -DPOSIX_THREADS -I"c:/Program Files"/MATLAB_SV71/extern/include -I Dynare_pp/src -I Dynare_pp/kord -I Dynare_pp/tl/cc -I Dynare_pp/utils/cc/ -I Dynare_pp/sylv/cc/ -I"f:/Pthreads/Pre-built.2/include" -I"f:/mingw/include" -o k_order_perturbation.dll k_ord_dynare.o k_order_perturbation.o dynamic_dll.o nlsolve.o "dynare_pp/extern/matlab"/dynarelib.a -Wl,-L"c:/Program Files"/MATLAB_SV71/extern/lib/win32/microsoft/ -Wl,-llibmex -Wl,-llibmx -Wl,-llibmwlapack -Wl,-llibdflapack -lg2c -lmingw32  -lstdc++ -L"f:/Pthreads/Pre-built.2/lib" -lpthreadGC2
diff --git a/mex/sources/korderpert/src/linkdbgexe.bat b/mex/sources/korderpert/src/linkdbgexe.bat
index ef625380f1fcb991b3b4c5ce1a02da2495189c44..e6ea506fd11bf1885e9183b0c6ae994ad58465b8 100644
--- a/mex/sources/korderpert/src/linkdbgexe.bat
+++ b/mex/sources/korderpert/src/linkdbgexe.bat
@@ -1 +1 @@
-gcc -DMATLAB_MEX_FILE -mthreads -g -DWINDOWS -DPOSIX_THREADS -I"c:/Program Files"/MATLAB_SV71/extern/include -I Dynare_pp/src -I Dynare_pp/kord -I Dynare_pp/tl/cc -I Dynare_pp/utils/cc/ -I Dynare_pp/sylv/cc/ -I"f:/Pthreads/Pre-built.2/include" -I"f:/mingw/include" -o k_orddbgtest.exe k_order_test_main.o k_ord_dynare.o k_order_perturbation.o nlsolve.o "dynare_pp/extern/matlab"/dynarelib.a -Wl,-L"c:/Program Files"/MATLAB_SV71/extern/lib/win32/microsoft/ -Wl,-llibmex -Wl,-llibmx -Wl,-llibmwlapack -Wl,-llibdflapack -lg2c -lmingw32  -lstdc++ -L"f:/Pthreads/Pre-built.2/lib" -lpthreadGC2
+gcc -DMATLAB_MEX_FILE -mthreads -g -DWINDOWS -DPOSIX_THREADS -I"c:/Program Files"/MATLAB_SV71/extern/include -I Dynare_pp/src -I Dynare_pp/kord -I Dynare_pp/tl/cc -I Dynare_pp/utils/cc/ -I Dynare_pp/sylv/cc/ -I"f:/Pthreads/Pre-built.2/include" -I"f:/mingw/include" -o k_orddbgtest.exe k_order_test_main.o k_ord_dynare.o k_order_perturbation.o nlsolve.o  dynamic_dll.o "dynare_pp/extern/matlab"/dynarelib.a -Wl,-L"c:/Program Files"/MATLAB_SV71/extern/lib/win32/microsoft/ -Wl,-llibmex -Wl,-llibmx -Wl,-llibmwlapack -Wl,-llibdflapack -lg2c -lmingw32  -lstdc++ -L"f:/Pthreads/Pre-built.2/lib" -lpthreadGC2