diff --git a/matlab/dynare_config.m b/matlab/dynare_config.m
index 97b18745dab06319675787eac815ac16da73ccbf..c38c254707106da52f787f63c239786f63975575 100644
--- a/matlab/dynare_config.m
+++ b/matlab/dynare_config.m
@@ -267,4 +267,17 @@ if verbose
disp(' ')
end
+% Test if getPowerDeriv DLL is present
+if exist('getPowerDeriv', 'file') == 3
+ message = '[mex] ';
+else
+ message = '[no] ';
+ addpath([dynareroot 'getPowerDeriv']);
+end
+if verbose
+ disp([ message 'getPowerDeriv routine.' ])
+ disp(' ')
+end
+
+
cd(origin);
diff --git a/matlab/getPowerDeriv.m b/matlab/getPowerDeriv/getPowerDeriv.m
similarity index 99%
rename from matlab/getPowerDeriv.m
rename to matlab/getPowerDeriv/getPowerDeriv.m
index e87eabe1e026cebbf2361d8b8b749ba92892310c..ab6e6b6d54ee765be31578d0c76ebf827d54dd63 100644
--- a/matlab/getPowerDeriv.m
+++ b/matlab/getPowerDeriv/getPowerDeriv.m
@@ -39,4 +39,4 @@ else
p = p-1;
end
end
-end
+end
\ No newline at end of file
diff --git a/mex/build/getPowerDeriv.am b/mex/build/getPowerDeriv.am
new file mode 100644
index 0000000000000000000000000000000000000000..4f7f4fd8322455713a3eaa217a364216701e63f7
--- /dev/null
+++ b/mex/build/getPowerDeriv.am
@@ -0,0 +1,5 @@
+vpath %.cc $(top_srcdir)/../../sources/get_power_deriv
+
+noinst_PROGRAMS = getPowerDeriv
+
+nodist_getPowerDeriv_SOURCES = getPowerDeriv.cc
diff --git a/mex/build/matlab/Makefile.am b/mex/build/matlab/Makefile.am
index 58dbd369f8eb5883fb3448f9a8ccfa051dc86165..8038a7ed9aa0eb2e6a66d2b4d31eaf953fae1592 100644
--- a/mex/build/matlab/Makefile.am
+++ b/mex/build/matlab/Makefile.am
@@ -2,7 +2,7 @@ ACLOCAL_AMFLAGS = -I ../../../m4
# libdynare++ must come before gensylv, k_order_perturbation, dynare_simul_
if DO_SOMETHING
-SUBDIRS = mjdgges kronecker bytecode libdynare++ gensylv k_order_perturbation dynare_simul_ estimation block_kalman_filter sobol
+SUBDIRS = mjdgges kronecker bytecode libdynare++ gensylv k_order_perturbation dynare_simul_ estimation block_kalman_filter sobol get_power_deriv
if HAVE_GSL
SUBDIRS += ms_sbvar
diff --git a/mex/build/matlab/configure.ac b/mex/build/matlab/configure.ac
index 4267a3158cbdecac349b0a29a3a22d07131aa46c..ab97d1970f92fb674738d97178e53ed44589470f 100644
--- a/mex/build/matlab/configure.ac
+++ b/mex/build/matlab/configure.ac
@@ -144,6 +144,7 @@ AC_CONFIG_FILES([Makefile
kalman_steady_state/Makefile
ms_sbvar/Makefile
block_kalman_filter/Makefile
- sobol/Makefile])
+ sobol/Makefile
+ get_power_deriv/Makefile])
AC_OUTPUT
diff --git a/mex/build/matlab/get_power_deriv/Makefile.am b/mex/build/matlab/get_power_deriv/Makefile.am
new file mode 100644
index 0000000000000000000000000000000000000000..e04608481a3fdfa5df8dc199ee6ed84dda7ef494
--- /dev/null
+++ b/mex/build/matlab/get_power_deriv/Makefile.am
@@ -0,0 +1,2 @@
+include ../mex.am
+include ../../getPowerDeriv.am
diff --git a/mex/build/octave/Makefile.am b/mex/build/octave/Makefile.am
index a6d164c9012eaf9d2ac26e2343e5a56eb5bdfb00..e11aa0423ac950b5c76aef2fd494fd4de241c74a 100644
--- a/mex/build/octave/Makefile.am
+++ b/mex/build/octave/Makefile.am
@@ -2,7 +2,7 @@ ACLOCAL_AMFLAGS = -I ../../../m4
# libdynare++ must come before gensylv, k_order_perturbation, dynare_simul_
if DO_SOMETHING
-SUBDIRS = mjdgges kronecker bytecode libdynare++ gensylv k_order_perturbation dynare_simul_ qzcomplex ordschur block_kalman_filter sobol
+SUBDIRS = mjdgges kronecker bytecode libdynare++ gensylv k_order_perturbation dynare_simul_ qzcomplex ordschur block_kalman_filter sobol get_power_deriv
if HAVE_GSL
SUBDIRS += ms_sbvar
diff --git a/mex/build/octave/configure.ac b/mex/build/octave/configure.ac
index daa13a6af00bdb92dbe6bbe4674bec58b9a12f3d..15a92f5ec0dde3c8f1e8cf32cab4b3896d319fd9 100644
--- a/mex/build/octave/configure.ac
+++ b/mex/build/octave/configure.ac
@@ -131,6 +131,7 @@ AC_CONFIG_FILES([Makefile
kalman_steady_state/Makefile
ms_sbvar/Makefile
block_kalman_filter/Makefile
- sobol/Makefile])
+ sobol/Makefile
+ get_power_deriv/Makefile])
AC_OUTPUT
diff --git a/mex/build/octave/get_power_deriv/Makefile.am b/mex/build/octave/get_power_deriv/Makefile.am
new file mode 100644
index 0000000000000000000000000000000000000000..1eaa678234a25799c50dee99f23e995ca6ba67ce
--- /dev/null
+++ b/mex/build/octave/get_power_deriv/Makefile.am
@@ -0,0 +1,3 @@
+EXEEXT = .mex
+include ../mex.am
+include ../../getPowerDeriv.am
diff --git a/mex/sources/Makefile.am b/mex/sources/Makefile.am
index 7444143f907b4d09cef614d06b33f3631dfc7bf1..963a6e0e1dd4a7d49dfb2fda2044d500822c0811 100644
--- a/mex/sources/Makefile.am
+++ b/mex/sources/Makefile.am
@@ -14,7 +14,8 @@ EXTRA_DIST = \
kalman_steady_state \
ms-sbvar \
block_kalman_filter \
- sobol
+ sobol \
+ get_power_deriv
clean-local:
rm -rf `find mex/sources -name *.o`
diff --git a/mex/sources/get_power_deriv/getPowerDeriv.cc b/mex/sources/get_power_deriv/getPowerDeriv.cc
new file mode 100644
index 0000000000000000000000000000000000000000..a1634b46a47aec356b85aa0f413aefbe30513346
--- /dev/null
+++ b/mex/sources/get_power_deriv/getPowerDeriv.cc
@@ -0,0 +1,87 @@
+/*
+** mex file for matlab getPowerDeriv.m.
+**
+** Copyright (C) 2012 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/>.
+**
+** AUTHOR(S): stephane DOT adjemian AT univ DASH lemans DOT fr
+**/
+
+
+#include <math.h>
+#include <dynmex.h>
+#include "mex.h"
+
+#ifdef _MSC_VER
+# define nearbyint(x) (fabs((x)-floor(x)) < fabs((x)-ceil(x)) ? floor(x) : ceil(x))
+#endif
+
+double getPowerDeriv(double *x, double *p, int k)
+{
+ if ( fabs(*x) < 1e-12 && *p > 0 && k >= *p && fabs(*p-nearbyint(*p)) < 1e-12 )
+ return 0.0;
+ else
+ {
+ //int i = 0;
+ double dxp = pow(*x, *p-k);
+ for (int i=0; i<k; i++)
+ dxp *= (*p)--;
+ return dxp;
+ }
+}
+
+void mexFunction( int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
+{
+ /*
+ ** INPUTS:
+ ** prhs[0] [double] scalar, x.
+ ** prhs[1] [double] scalar, p.
+ ** prhs[2] [integer] scalar, k.
+ **
+ ** OUTPUTS:
+ ** plhs[0] [double] scalar, The k-th derivative of x^p.
+ ** plhs[1] [double] scalar, info variable (equal to zero if call to getPowerDeriv is successfull).
+ */
+ if ( !( nrhs==3) )
+ {
+ DYN_MEX_FUNC_ERR_MSG_TXT("getPowerDeriv:: Three input arguments are required (x,p and k)!");
+ }
+ if (nlhs>2)
+ {
+ DYN_MEX_FUNC_ERR_MSG_TXT("getPowerDeriv:: I only return one or two output arguments!");
+ }
+ if ( !(mxGetN(prhs[0])==1 && (mxGetM(prhs[0]))==1) )
+ {
+ DYN_MEX_FUNC_ERR_MSG_TXT("getPowerDeriv:: First argument, x, must be a double scalar!");
+ }
+ if ( !(mxGetN(prhs[1])==1 && (mxGetM(prhs[1]))==1) )
+ {
+ DYN_MEX_FUNC_ERR_MSG_TXT("getPowerDeriv:: Second argument, p, must be a double scalar!");
+ }
+ if ( !(mxGetN(prhs[2])==1 && (mxGetM(prhs[2]))==1) )
+ {
+ DYN_MEX_FUNC_ERR_MSG_TXT("getPowerDeriv:: Third argument, k, must be an integer scalar!");
+ }
+ double *x = mxGetPr(prhs[0]);
+ double *p = mxGetPr(prhs[1]);
+ int k = (int) mxGetScalar(prhs[2]);
+ plhs[0] = mxCreateDoubleScalar(getPowerDeriv(x,p,k));
+ if (nlhs>1)
+ {
+ plhs[1] = mxCreateDoubleScalar(0.0);
+ }
+}