diff --git a/matlab/dr1.m b/matlab/dr1.m
index 0be5b97c697db7d2a100617b473b0caf329185d6..1d2ec1b9252bc28720c3cc1c3e079eb397e17cb8 100644
--- a/matlab/dr1.m
+++ b/matlab/dr1.m
@@ -240,6 +240,11 @@ function [dr,info,M_,options_,oo_] = dr1(dr,task,M_,options_,oo_)
             [junk,jacobia_,hessian1] = feval([M_.fname '_dynamic'],z,...
                                             [oo_.exo_simul ...
                                 oo_.exo_det_simul], M_.params, it_);
+            if options_.use_dll
+                % In USE_DLL mode, the hessian is in the 3-column sparse representation
+                hessian1 = sparse(hessian1(:,1), hessian1(:,2), hessian1(:,3), ...
+                                  size(jacobia_, 1), size(jacobia_, 2)*size(jacobia_, 2));
+            end
         end
     end
     
diff --git a/matlab/global_initialization.m b/matlab/global_initialization.m
index 3656583749a4798653f68fbc381bd7dead0d5e53..89021f92fcf26672b5692824c2d05c6467379c64 100644
--- a/matlab/global_initialization.m
+++ b/matlab/global_initialization.m
@@ -226,6 +226,9 @@ function global_initialization()
   % did model undergo block decomposition + minimum feedback set computation ?
   options_.block = 0;
 
+  % model evaluated using a compiled MEX
+  options_.use_dll = 0;
+  
   % model evaluated using bytecode.dll
   options_.bytecode = 0;
   
diff --git a/preprocessor/ModFile.cc b/preprocessor/ModFile.cc
index 0695a53703ed11ccc7bac9a8fc9a3a2eb40d0ca6..42c829b84ccf4fbed6d6bafd5be313177d52e119 100644
--- a/preprocessor/ModFile.cc
+++ b/preprocessor/ModFile.cc
@@ -287,7 +287,8 @@ ModFile::writeOutputFiles(const string &basename, bool clear_all
     mOutputFile << "options_.linear = 1;" << endl;
 
   mOutputFile << "options_.block=" << block << ";" << endl
-              << "options_.bytecode=" << byte_code << ";" << endl;
+              << "options_.bytecode=" << byte_code << ";" << endl
+              << "options_.use_dll=" << use_dll << ";" << endl;
 
   if (byte_code)
     mOutputFile << "if exist('bytecode') ~= 3" << endl
diff --git a/tests/Makefile.am b/tests/Makefile.am
index e8d809d6d044fc8addd13c3c1af87ef4a02502eb..dc215e0c1b87e09ad4ce232426898c86eaf0533f 100644
--- a/tests/Makefile.am
+++ b/tests/Makefile.am
@@ -5,6 +5,7 @@ OCTAVE_MODS = \
 	ramst.mod \
 	ramst_a.mod \
 	example1.mod \
+	example1_use_dll.mod \
 	t_sgu_ex1.mod \
 	ramsey.mod \
 	ramst_initval_file.mod \
diff --git a/tests/example1_use_dll.mod b/tests/example1_use_dll.mod
new file mode 100644
index 0000000000000000000000000000000000000000..4e630c21ece20035e56432f22a768e0f31ddd5c0
--- /dev/null
+++ b/tests/example1_use_dll.mod
@@ -0,0 +1,73 @@
+// Test USE_DLL option at order 2
+
+var y, c, k, a, h, b;
+varexo e,u;
+
+parameters beta, rho, alpha, delta, theta, psi, tau, phi;
+
+alpha = 0.36;
+rho   = 0.95;
+tau   = 0.025;
+beta  = 0.99;
+delta = 0.025;
+psi   = 0;
+theta = 2.95;
+
+phi   = 0.1;
+
+model(use_dll);
+c*theta*h^(1+psi)=(1-alpha)*y;
+k = beta*(((exp(b)*c)/(exp(b(+1))*c(+1)))
+    *(exp(b(+1))*alpha*y(+1)+(1-delta)*k));
+y = exp(a)*(k(-1)^alpha)*(h^(1-alpha));
+k = exp(b)*(y-c)+(1-delta)*k(-1);
+a = rho*a(-1)+tau*b(-1) + e;
+b = tau*a(-1)+rho*b(-1) + u;
+end;
+
+initval;
+y = 1.08068253095672;
+c = 0.80359242014163;
+h = 0.29175631001732;
+k = 5;
+a = 0;
+b = 0;
+e = 0;
+u = 0;
+end;
+
+shocks;
+var e; stderr 0.009;
+var u; stderr 0.009;
+var e, u = phi*0.009*0.009;
+end;
+
+stoch_simul(nograph);
+
+if ~exist('example1_results.mat','file');
+   error('example1 must be run first');
+end;
+
+oo1 = load('example1_results','oo_');
+
+dr0 = oo1.oo_.dr;
+dr = oo_.dr;
+
+if max(max(abs(dr0.ghx - dr.ghx))) > 1e-12;
+   error('error in ghx');
+end;
+if max(max(abs(dr0.ghu - dr.ghu))) > 1e-12;
+   error('error in ghu');
+end;
+if max(max(abs(dr0.ghxx - dr.ghxx))) > 1e-12;
+   error('error in ghxx');
+end;
+if max(max(abs(dr0.ghuu - dr.ghuu))) > 1e-12;
+   error('error in ghuu');
+end;
+if max(max(abs(dr0.ghxu - dr.ghxu))) > 1e-12;
+   error('error in ghxu');
+end;
+if max(max(abs(dr0.ghs2 - dr.ghs2))) > 1e-12;
+   error('error in ghs2');
+end;