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

Fix : 

"Bytecode does not enforce positivity constraint on irreversible investment model #439 "
when block model option is used.
When cutoff is set to 0 the block decomposition is performed on the incidence matrix instead of the static Jacobian matrix. This correction avoids to use derivatives on logical operators that are always equal to 0.
issue#70
Ferhat Mihoubi 2016-03-25 15:38:49 +01:00
parent d52c378460
commit 4419b9e9e6
1 changed files with 33 additions and 7 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

40
ModelTree.cc
View File

@ -429,9 +429,21 @@ ModelTree::computePrologueAndEpilogue(const jacob_map_t &static_jacobian_arg, ve
eq2endo[*it] = i;
equation_reordered[i] = i;
variable_reordered[*it] = i;
}
for (jacob_map_t::const_iterator it = static_jacobian_arg.begin(); it != static_jacobian_arg.end(); it++)
IM[it->first.first * n + endo2eq[it->first.second]] = true;
}
if (cutoff == 0)
{
set<pair<int, int> > endo;
for (int i = 0; i < n; i++)
{
endo.clear();
equations[i]->collectEndogenous(endo);
for (set<pair<int, int> >::const_iterator it = endo.begin(); it != endo.end(); it++)
IM[i * n + endo2eq[it->first]] = true;
}
}
else
for (jacob_map_t::const_iterator it = static_jacobian_arg.begin(); it != static_jacobian_arg.end(); it++)
IM[it->first.first * n + endo2eq[it->first.second]] = true;
bool something_has_been_done = true;
prologue = 0;
int k = 0;
@ -631,9 +643,23 @@ ModelTree::computeBlockDecompositionAndFeedbackVariablesForEachBlock(const jacob
{
reverse_equation_reordered[equation_reordered[i]] = i;
reverse_variable_reordered[variable_reordered[i]] = i;
}
}
jacob_map_t tmp_normalized_contemporaneous_jacobian;
if (cutoff == 0)
{
set<pair<int, int> > endo;
for (int i = 0; i < nb_var; i++)
{
endo.clear();
equations[i]->collectEndogenous(endo);
for (set<pair<int, int> >::const_iterator it = endo.begin(); it != endo.end(); it++)
tmp_normalized_contemporaneous_jacobian[make_pair(i, it->first)] = 1;
for (jacob_map_t::const_iterator it = static_jacobian.begin(); it != static_jacobian.end(); it++)
}
}
else
tmp_normalized_contemporaneous_jacobian = static_jacobian;
for (jacob_map_t::const_iterator it = tmp_normalized_contemporaneous_jacobian.begin(); it != tmp_normalized_contemporaneous_jacobian.end(); it++)
if (reverse_equation_reordered[it->first.first] >= (int) prologue && reverse_equation_reordered[it->first.first] < (int) (nb_var - epilogue)
&& reverse_variable_reordered[it->first.second] >= (int) prologue && reverse_variable_reordered[it->first.second] < (int) (nb_var - epilogue)
&& it->first.first != endo2eq[it->first.second])
@ -1116,7 +1142,7 @@ ModelTree::computeJacobian(const set<int> &vars)
continue;
first_derivatives[make_pair(eq, *it)] = d1;
++NNZDerivatives[0];
}
}
}
}
@ -1635,7 +1661,7 @@ ModelTree::computeParamsDerivatives()
{
set<int> deriv_id_set;
addAllParamDerivId(deriv_id_set);
for (set<int>::const_iterator it = deriv_id_set.begin();
it != deriv_id_set.end(); it++)
{