From c785a81d03dbd4ca6fdaa683304042eae7ad7698 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?S=C3=A9bastien=20Villemot?= <sebastien@dynare.org>
Date: Tue, 18 Feb 2020 12:18:39 +0100
Subject: [PATCH] k-order MEX: fix crash on Octave when higher-order
 derivatives have numerical zero

In Octave, when some values given to the sparse() function are numerically
zero, then the nzmax of the generated sparse matrix is shrinked accordingly;
while under MATLAB, the nzmax is the length of the vector of values, zeros
included.

The check at the top of
DynamicModelMFile::unpackSparseMatrixAndCopyIntoTwoDMatData() would then fail
under Octave if some higher-derivatives had an element which is symbolically
non-zero but numerically zero.

We therefore relax the check, and accordingly adapt the code that handles
numerical zeros.

This bug was uncovered by tests/pruning/AnSchorfheide_pruned_state_space.mod,
which was failing under Octave.

(cherry picked from commit 144bdf34bab11a4d0d51ba0115ba98b3cc29d979)
---
 mex/sources/k_order_perturbation/dynamic_m.cc | 12 +++++++-----
 1 file changed, 7 insertions(+), 5 deletions(-)

diff --git a/mex/sources/k_order_perturbation/dynamic_m.cc b/mex/sources/k_order_perturbation/dynamic_m.cc
index 04473cfe24..c0ff5c2ca5 100644
--- a/mex/sources/k_order_perturbation/dynamic_m.cc
+++ b/mex/sources/k_order_perturbation/dynamic_m.cc
@@ -1,5 +1,5 @@
 /*
- * Copyright © 2010-2019 Dynare Team
+ * Copyright © 2010-2020 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -35,11 +35,13 @@ DynamicModelMFile::unpackSparseMatrixAndCopyIntoTwoDMatData(mxArray *sparseMat,
 {
   int totalCols = mxGetN(sparseMat);
   mwIndex *rowIdxVector = mxGetIr(sparseMat);
-  mwSize sizeRowIdxVector = mxGetNzmax(sparseMat);
   mwIndex *colIdxVector = mxGetJc(sparseMat);
 
   assert(tdm.ncols() == 3);
-  assert(tdm.nrows() == sizeRowIdxVector);
+  /* Under MATLAB, the following check always holds at equality; under Octave,
+     there may be an inequality, because Octave diminishes nzmax if one gives
+     zeros in the values vector when calling sparse(). */
+  assert(tdm.nrows() >= mxGetNzmax(sparseMat));
 
   double *ptr = mxGetPr(sparseMat);
 
@@ -55,11 +57,11 @@ DynamicModelMFile::unpackSparseMatrixAndCopyIntoTwoDMatData(mxArray *sparseMat,
         output_row++;
       }
 
-  /* If there are less elements than Nzmax (that might happen if some
+  /* If there are less elements than expected (that might happen if some
      derivative is symbolically not zero but numerically zero at the evaluation
      point), then fill in the matrix with empty entries, that will be
      recognized as such by KordpDynare::populateDerivativesContainer() */
-  while (output_row < static_cast<int>(sizeRowIdxVector))
+  while (output_row < tdm.nrows())
     {
       tdm.get(output_row, 0) = 0;
       tdm.get(output_row, 1) = 0;
-- 
GitLab