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.
Because at some point throwing exceptions from MEX files (with mexErrMsgTxt())
was not working under Windows 64-bit, we had designed a workaround to avoid
using exceptions.
Most MEX files were returning an error code as their first (or sometimes last)
argument, and that code would have to be checked from the MATLAB code.
Since this workaround is no longer needed, this commit removes it. As a
consequence, the interface of many MEX files is modified.
For some background, see https://www.dynare.org/pipermail/dev/2010-September/000895.html
It applies the approximated policy function to a set of particles, using
Dynare++ routines.
There is support for parallelization, using Dynare++ multithreading
model (itself based on C++11 threads; we don’t use OpenMP because it is
incompatible with MKL). For the time being, default to a single thread. This
should be later refined through empirical testing.
This MEX solves nonlinear systems of equations using a trust region algorithm.
The problem is subdivided in smaller problems by doing a block
triangularisation of the Jacobian at the guess value, using the
Dulmage-Mendelsohn algorithm.
The interface of the MEX is simply:
[x, info] = block_trust_region(f, guess_value);
Where f is either a function handle or a string designating a function.
f must take one argument (the evaluation point), and return either one or two
arguments (the residuals and, optionally, the Jacobian).
On success, info=0; on failure, info=1.
The logic of the dynSparseMatrix::Sparse_substract_SA_SB() routine was
incorrect.
In some cases, it would read past the last nonzero elements of the A matrix,
and consequently write past the number of allocated nonzero elements of the C
matrix.
This would lead to crashes and, probably, to wrong results under certain
circumstances.
Closes: #1652
It constructs the stacked residuals and jacobian of the perfect foresight
problem.
It is an almost perfect replacement for the perfect_foresight_problem.m
routine, while being much more efficient.
Note however that the DLL never return complex numbers (it instead puts NaNs at
the place where there would have been complex). This may create problems for
some MOD files; the algorithms will need to be adapted to use a more
line-search method.
Various modernizations and simplifications.
Also remove a workaround for a LAPACK bug in DGGES (the VSL argument was
apparently referenced even though JOBVSL="N"). Hopefully the bug has been fixed
everywhere now.
Previously there were GeneralMatrix::numRows() and TwoDMatrix::nrows() for doing
the same thing (and same for columns and Const versions).
Merge these two into GeneralMatrix::nrows().
Sébastien Villemot
Houtan Bastani