From cf77c3f24519976ce111ada29eba4b1a7b627f58 Mon Sep 17 00:00:00 2001
From: sebastien <sebastien@ac1d8469-bf42-47a9-8791-bf33cf982152>
Date: Mon, 14 Dec 2009 11:06:29 +0000
Subject: [PATCH] * Fix USE_DLL option at order 2 (closes #64) * Added a
 corresponding test

git-svn-id: https://www.dynare.org/svn/dynare/trunk@3236 ac1d8469-bf42-47a9-8791-bf33cf982152
---
 matlab/dr1.m                   |  5 +++
 matlab/global_initialization.m |  3 ++
 preprocessor/ModFile.cc        |  3 +-
 tests/Makefile.am              |  1 +
 tests/example1_use_dll.mod     | 73 ++++++++++++++++++++++++++++++++++
 5 files changed, 84 insertions(+), 1 deletion(-)
 create mode 100644 tests/example1_use_dll.mod

diff --git a/matlab/dr1.m b/matlab/dr1.m
index 0be5b97c6..1d2ec1b92 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 365658374..89021f92f 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 0695a5370..42c829b84 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 e8d809d6d..dc215e0c1 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 000000000..4e630c21e
--- /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;