stepan
/
dynare
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

adding possibility to pass derivative matrices to k_order_perturbation 

MEX function
time-shift
Michel Juillard 2012-07-07 21:21:28 +02:00
parent 50ef5377e4
commit 84164dc57e
3 changed files with 88 additions and 28 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

74
mex/sources/k_order_perturbation/k_ord_dynare.cc
View File

@ -46,7 +46,30 @@ KordpDynare::KordpDynare(const vector<string> &endo, int num_endo,
nYs(npred + nboth), nYss(nboth + nforw), nY(num_endo), nJcols(jcols), NNZD(nnzd), nSteps(nsteps),
nOrder(norder), journal(jr), ySteady(ysteady), params(inParams), vCov(vcov),
md(1), dnl(*this, endo), denl(*this, exo), dsnl(*this, dnl, denl), ss_tol(sstol), varOrder(var_order),
ll_Incidence(llincidence), qz_criterium(criterium), dynamicModelFile(dynamicModelFile_arg)
ll_Incidence(llincidence), qz_criterium(criterium), dynamicModelFile(dynamicModelFile_arg), g1p(NULL),
g2p(NULL), g3p(NULL)
{
ReorderDynareJacobianIndices();
// Initialise ModelDerivativeContainer(*this, this->md, nOrder);
for (int iord = 1; iord <= nOrder; iord++)
md.insert(new FSSparseTensor(iord, nY+nYs+nYss+nExog, nY));
}
KordpDynare::KordpDynare(const vector<string> &endo, int num_endo,
const vector<string> &exo, int nexog, int npar,
Vector &ysteady, TwoDMatrix &vcov, Vector &inParams, int nstat,
int npred, int nforw, int nboth, const int jcols, const Vector &nnzd,
const int nsteps, int norder,
Journal &jr, DynamicModelAC *dynamicModelFile_arg, double sstol,
const vector<int> &var_order, const TwoDMatrix &llincidence, double criterium,
TwoDMatrix *g1_arg, TwoDMatrix *g2_arg, TwoDMatrix *g3_arg) throw (TLException) :
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),
nOrder(norder), journal(jr), ySteady(ysteady), params(inParams), vCov(vcov),
md(1), dnl(*this, endo), denl(*this, exo), dsnl(*this, dnl, denl), ss_tol(sstol), varOrder(var_order),
ll_Incidence(llincidence), qz_criterium(criterium), dynamicModelFile(dynamicModelFile_arg),
g1p(g1_arg), g2p(g2_arg), g3p(g3_arg)
{
ReorderDynareJacobianIndices();
@ -89,32 +112,35 @@ KordpDynare::evaluateSystem(Vector &out, const Vector &yym, const Vector &yy,
void
KordpDynare::calcDerivativesAtSteady() throw (DynareException)
{
g1p = new TwoDMatrix(nY, nJcols);
g1p->zeros();
if (nOrder > 1)
if (g1p == NULL)
{
// allocate space for sparse Hessian
g2p = new TwoDMatrix((int) NNZD[1], 3);
g2p->zeros();
g1p = new TwoDMatrix(nY, nJcols);
g1p->zeros();
if (nOrder > 1)
{
// allocate space for sparse Hessian
g2p = new TwoDMatrix((int) NNZD[1], 3);
g2p->zeros();
}
if (nOrder > 2)
{
g3p = new TwoDMatrix((int) NNZD[2], 3);
g3p->zeros();
}
Vector xx(nexog());
xx.zeros();
Vector out(nY);
out.zeros();
Vector llxSteady(nJcols-nExog);
LLxSteady(ySteady, llxSteady);
dynamicModelFile->eval(llxSteady, xx, params, ySteady, out, g1p, g2p, g3p);
}
if (nOrder > 2)
{
g3p = new TwoDMatrix((int) NNZD[2], 3);
g3p->zeros();
}
Vector xx(nexog());
xx.zeros();
Vector out(nY);
out.zeros();
Vector llxSteady(nJcols-nExog);
LLxSteady(ySteady, llxSteady);
dynamicModelFile->eval(llxSteady, xx, params, ySteady, out, g1p, g2p, g3p);
populateDerivativesContainer(*g1p, 1, JacobianIndices);
delete g1p;

8
mex/sources/k_order_perturbation/k_ord_dynare.hh
View File

@ -126,6 +126,14 @@ public:
Journal &jr, DynamicModelAC *dynamicModelFile_arg, double sstol,
const vector<int> &varOrder, const TwoDMatrix &ll_Incidence,
double qz_criterium) throw (TLException);
KordpDynare(const vector<string> &endo, int num_endo,
const vector<string> &exo, int num_exo, int num_par,
Vector &ySteady, TwoDMatrix &vCov, Vector &params, int nstat, int nPred,
int nforw, int nboth, const int nJcols, const Vector &NNZD,
const int nSteps, const int ord,
Journal &jr, DynamicModelAC *dynamicModelFile_arg, double sstol,
const vector<int> &varOrder, const TwoDMatrix &ll_Incidence,
double qz_criterium, TwoDMatrix *g1_arg, TwoDMatrix *g2_arg, TwoDMatrix *g3_arg) throw (TLException);
virtual ~KordpDynare();
int

34
mex/sources/k_order_perturbation/k_order_perturbation.cc
View File

@ -68,8 +68,8 @@ extern "C" {
mexFunction(int nlhs, mxArray *plhs[],
int nrhs, const mxArray *prhs[])
{
if (nrhs != 3 || nlhs < 2)
DYN_MEX_FUNC_ERR_MSG_TXT("Must have exactly 3 input parameters and take at least 2 output parameters.");
if (nrhs < 3 || nlhs < 2)
DYN_MEX_FUNC_ERR_MSG_TXT("Must have at least 3 input parameters and takes at least 2 output parameters.");
const mxArray *dr = prhs[0];
const mxArray *M_ = prhs[1];
@ -189,7 +189,32 @@ extern "C" {
if ((nEndo != nendo) || (nExog != nexo))
DYN_MEX_FUNC_ERR_MSG_TXT("Incorrect number of input parameters.");
/* Fetch time index */
TwoDMatrix *g1m=NULL;
TwoDMatrix *g2m=NULL;
TwoDMatrix *g3m=NULL;
// derivatives passed as arguments */
if (nrhs > 3)
{
const mxArray *g1 = prhs[3];
int m = (int) mxGetM(g1);
int n = (int) mxGetN(g1);
g1m = new TwoDMatrix(m, n, mxGetPr(g1));
if (nrhs > 4)
{
const mxArray *g2 = prhs[4];
int m = (int) mxGetM(g2);
int n = (int) mxGetN(g2);
g2m = new TwoDMatrix(m, n, mxGetPr(g2));
if (nrhs > 5)
{
const mxArray *g3 = prhs[5];
int m = (int) mxGetM(g3);
int n = (int) mxGetN(g3);
g3m = new TwoDMatrix(m, n, mxGetPr(g3));
}
}
}
const int nSteps = 0; // Dynare++ solving steps, for time being default to 0 = deterministic steady state
const double sstol = 1.e-13; //NL solver tolerance from
@ -216,7 +241,8 @@ extern "C" {
KordpDynare dynare(endoNames, nEndo, exoNames, nExog, nPar,
ySteady, vCov, modParams, nStat, nPred, nForw, nBoth,
jcols, NNZD, nSteps, kOrder, journal, dynamicModelFile,
sstol, var_order_vp, llincidence, qz_criterium);
sstol, var_order_vp, llincidence, qz_criterium,
g1m, g2m, g3m);
// construct main K-order approximation class