diff --git a/preprocessor/StaticDllModel.cc b/preprocessor/StaticDllModel.cc
new file mode 100644
index 0000000000000000000000000000000000000000..da6cf82389fd7fd4a9e70eef0531ff3976337a09
--- /dev/null
+++ b/preprocessor/StaticDllModel.cc
@@ -0,0 +1,1434 @@
+/*
+ * Copyright (C) 2003-2009 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/>.
+ */
+
+#include <iostream>
+#include <cmath>
+#include <cstdlib>
+#include <cassert>
+#include <cstdio>
+#include <cerrno>
+#include "StaticDllModel.hh"
+
+// For mkdir() and chdir()
+#ifdef _WIN32
+# include <direct.h>
+#else
+# include <unistd.h>
+# include <sys/stat.h>
+# include <sys/types.h>
+#endif
+
+StaticDllModel::StaticDllModel(SymbolTable &symbol_table_arg,
+ NumericalConstants &num_constants_arg) :
+ ModelTree(symbol_table_arg, num_constants_arg),
+ max_lag(0), max_lead(0),
+ max_endo_lag(0), max_endo_lead(0),
+ max_exo_lag(0), max_exo_lead(0),
+ max_exo_det_lag(0), max_exo_det_lead(0),
+ dynJacobianColsNbr(0),
+ mode(eStandardMode),
+ cutoff(1e-15),
+ markowitz(0.7),
+ mfs(0),
+ block_triangular(symbol_table_arg, num_constants_arg)
+{
+}
+
+NodeID
+StaticDllModel::AddVariable(const string &name, int lag)
+{
+ return AddVariableInternal(name, lag);
+}
+
+void
+StaticDllModel::compileDerivative(ofstream &code_file, int eq, int symb_id, int lag, map_idx_type &map_idx) const
+ {
+ //first_derivatives_type::const_iterator it = first_derivatives.find(make_pair(eq, getDerivID(symb_id, lag)));
+ first_derivatives_type::const_iterator it = first_derivatives.find(make_pair(eq, getDerivID(symbol_table.getID(eEndogenous, symb_id), lag)));
+ if (it != first_derivatives.end())
+ (it->second)->compile(code_file, false, temporary_terms, map_idx, false);
+ else
+ code_file.write(&FLDZ, sizeof(FLDZ));
+ }
+
+
+void
+StaticDllModel::compileChainRuleDerivative(ofstream &code_file, int eqr, int varr, int lag, map_idx_type &map_idx) const
+{
+ map<pair<int, pair<int, int> >, NodeID>::const_iterator it = first_chain_rule_derivatives.find(make_pair(eqr, make_pair(varr, lag)));
+ if (it != first_chain_rule_derivatives.end())
+ (it->second)->compile(code_file, false, temporary_terms, map_idx, false);
+ else
+ code_file.write(&FLDZ, sizeof(FLDZ));
+}
+
+
+void
+StaticDllModel::BuildIncidenceMatrix()
+{
+ set<pair<int, int> > endogenous, exogenous;
+ for (int eq = 0; eq < (int) equations.size(); eq++)
+ {
+ BinaryOpNode *eq_node = equations[eq];
+ endogenous.clear();
+ NodeID Id = eq_node->get_arg1();
+ Id->collectEndogenous(endogenous);
+ Id = eq_node->get_arg2();
+ Id->collectEndogenous(endogenous);
+ for (set<pair<int, int> >::iterator it_endogenous=endogenous.begin();it_endogenous!=endogenous.end();it_endogenous++)
+ {
+ block_triangular.incidencematrix.fill_IM(eq, it_endogenous->first, 0, eEndogenous);
+ }
+ exogenous.clear();
+ Id = eq_node->get_arg1();
+ Id->collectExogenous(exogenous);
+ Id = eq_node->get_arg2();
+ Id->collectExogenous(exogenous);
+ for (set<pair<int, int> >::iterator it_exogenous=exogenous.begin();it_exogenous!=exogenous.end();it_exogenous++)
+ {
+ block_triangular.incidencematrix.fill_IM(eq, it_exogenous->first, 0, eExogenous);
+ }
+ }
+}
+
+void
+StaticDllModel::computeTemporaryTermsOrdered(Model_Block *ModelBlock)
+{
+ map<NodeID, pair<int, int> > first_occurence;
+ map<NodeID, int> reference_count;
+ int i, j, eqr, varr, lag;
+ temporary_terms_type vect;
+ ostringstream tmp_output;
+ BinaryOpNode *eq_node;
+ first_derivatives_type::const_iterator it;
+ first_chain_rule_derivatives_type::const_iterator it_chr;
+ ostringstream tmp_s;
+
+ temporary_terms.clear();
+ map_idx.clear();
+ for (j = 0;j < ModelBlock->Size;j++)
+ {
+ // Compute the temporary terms reordered
+ for (i = 0;i < ModelBlock->Block_List[j].Size;i++)
+ {
+ if (ModelBlock->Block_List[j].Equation_Type[i] == E_EVALUATE_S && i<ModelBlock->Block_List[j].Nb_Recursives && ModelBlock->Block_List[j].Equation_Normalized[i])
+ ModelBlock->Block_List[j].Equation_Normalized[i]->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, i, map_idx);
+ else
+ {
+ eq_node = equations[ModelBlock->Block_List[j].Equation[i]];
+ eq_node->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, i, map_idx);
+ }
+ }
+ for(i=0; i<(int)ModelBlock->Block_List[j].Chain_Rule_Derivatives->size();i++)
+ {
+ pair< pair<int, pair<int, int> >, pair<int, int> > it = ModelBlock->Block_List[j].Chain_Rule_Derivatives->at(i);
+ lag=it.first.first;
+ int eqr=it.second.first;
+ int varr=it.second.second;
+ it_chr=first_chain_rule_derivatives.find(make_pair(eqr, make_pair( varr, lag)));
+ it_chr->second->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, ModelBlock->Block_List[j].Size-1, map_idx);
+ }
+
+ }
+ for (j = 0;j < ModelBlock->Size;j++)
+ {
+ // Collecte the temporary terms reordered
+ for (i = 0;i < ModelBlock->Block_List[j].Size;i++)
+ {
+ if (ModelBlock->Block_List[j].Equation_Type[i] == E_EVALUATE_S && i<ModelBlock->Block_List[j].Nb_Recursives && ModelBlock->Block_List[j].Equation_Normalized[i])
+ ModelBlock->Block_List[j].Equation_Normalized[i]->collectTemporary_terms(temporary_terms, ModelBlock, j);
+ else
+ {
+ eq_node = equations[ModelBlock->Block_List[j].Equation[i]];
+ eq_node->collectTemporary_terms(temporary_terms, ModelBlock, j);
+ }
+
+ }
+ for(i=0; i<(int)ModelBlock->Block_List[j].Chain_Rule_Derivatives->size();i++)
+ {
+ pair< pair<int, pair<int, int> >, pair<int, int> > it = ModelBlock->Block_List[j].Chain_Rule_Derivatives->at(i);
+ lag=it.first.first;
+ eqr=it.second.first;
+ varr=it.second.second;
+ it_chr=first_chain_rule_derivatives.find(make_pair(eqr, make_pair( varr, lag)));
+ it_chr->second->collectTemporary_terms(temporary_terms, ModelBlock, j);
+ }
+ }
+ // Add a mapping form node ID to temporary terms order
+ j=0;
+ for (temporary_terms_type::const_iterator it = temporary_terms.begin();
+ it != temporary_terms.end(); it++)
+ map_idx[(*it)->idx]=j++;
+}
+
+void
+StaticDllModel::writeModelEquationsOrdered_M( Model_Block *ModelBlock, const string &static_basename) const
+ {
+ int i,j,k,m;
+ string tmp_s, sps;
+ ostringstream tmp_output, tmp1_output, global_output;
+ NodeID lhs=NULL, rhs=NULL;
+ BinaryOpNode *eq_node;
+ map<NodeID, int> reference_count;
+ ofstream output;
+ int nze, nze_exo, nze_other_endo;
+ vector<int> feedback_variables;
+ //For each block
+ for (j = 0;j < ModelBlock->Size;j++)
+ {
+ //recursive_variables.clear();
+ feedback_variables.clear();
+ //For a block composed of a single equation determines wether we have to evaluate or to solve the equation
+ nze = nze_exo = nze_other_endo = 0;
+ for (m=0;m<=ModelBlock->Block_List[j].Max_Lead+ModelBlock->Block_List[j].Max_Lag;m++)
+ nze+=ModelBlock->Block_List[j].IM_lead_lag[m].size;
+ tmp1_output.str("");
+ tmp1_output << static_basename << "_" << j+1 << ".m";
+ output.open(tmp1_output.str().c_str(), ios::out | ios::binary);
+ output << "%\n";
+ output << "% " << tmp1_output.str() << " : Computes static model for Dynare\n";
+ output << "%\n";
+ output << "% Warning : this file is generated automatically by Dynare\n";
+ output << "% from model file (.mod)\n\n";
+ output << "%/\n";
+ if (ModelBlock->Block_List[j].Simulation_Type==EVALUATE_BACKWARD
+ ||ModelBlock->Block_List[j].Simulation_Type==EVALUATE_FORWARD)
+ output << "function y = " << static_basename << "_" << j+1 << "(y, x, params)\n";
+ else if (ModelBlock->Block_List[j].Simulation_Type==SOLVE_FORWARD_COMPLETE
+ || ModelBlock->Block_List[j].Simulation_Type==SOLVE_BACKWARD_COMPLETE
+ || ModelBlock->Block_List[j].Simulation_Type==SOLVE_BACKWARD_SIMPLE
+ || ModelBlock->Block_List[j].Simulation_Type==SOLVE_FORWARD_SIMPLE)
+ output << "function [residual, y, g1] = " << static_basename << "_" << j+1 << "(y, x, params)\n";
+ output << " % ////////////////////////////////////////////////////////////////////////" << endl
+ << " % //" << string(" Block ").substr(int(log10(j + 1))) << j + 1 << " " << BlockTriangular::BlockType0(ModelBlock->Block_List[j].Type)
+ << " //" << endl
+ << " % // Simulation type "
+ << BlockTriangular::BlockSim(ModelBlock->Block_List[j].Simulation_Type) << " //" << endl
+ << " % ////////////////////////////////////////////////////////////////////////" << endl;
+ //The Temporary terms
+ if (ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_BACKWARD
+ && ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_FORWARD)
+ output << " g1 = spalloc(" << ModelBlock->Block_List[j].Size-ModelBlock->Block_List[j].Nb_Recursives
+ << ", " << ModelBlock->Block_List[j].Size-ModelBlock->Block_List[j].Nb_Recursives << ", " << nze << ");\n";
+
+ if (ModelBlock->Block_List[j].Temporary_InUse->size())
+ {
+ tmp_output.str("");
+ for (temporary_terms_inuse_type::const_iterator it = ModelBlock->Block_List[j].Temporary_InUse->begin();
+ it != ModelBlock->Block_List[j].Temporary_InUse->end(); it++)
+ tmp_output << " T" << *it;
+ output << " global" << tmp_output.str() << ";\n";
+ }
+ if (ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_BACKWARD && ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_FORWARD)
+ output << " residual=zeros(" << ModelBlock->Block_List[j].Size-ModelBlock->Block_List[j].Nb_Recursives << ",1);\n";
+ // The equations
+ for (i = 0;i < ModelBlock->Block_List[j].Size;i++)
+ {
+ temporary_terms_type tt2;
+ tt2.clear();
+ if (ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->size())
+ output << " " << sps << "% //Temporary variables" << endl;
+ for (temporary_terms_type::const_iterator it = ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->begin();
+ it != ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->end(); it++)
+ {
+ output << " " << sps;
+ (*it)->writeOutput(output, oMatlabStaticModelSparse, temporary_terms);
+ output << " = ";
+ (*it)->writeOutput(output, oMatlabStaticModelSparse, tt2);
+ // Insert current node into tt2
+ tt2.insert(*it);
+ output << ";" << endl;
+ }
+ string sModel = symbol_table.getName(symbol_table.getID(eEndogenous, ModelBlock->Block_List[j].Variable[i])) ;
+ eq_node = equations[ModelBlock->Block_List[j].Equation[i]];
+ lhs = eq_node->get_arg1();
+ rhs = eq_node->get_arg2();
+ tmp_output.str("");
+ /*if((ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_BACKWARD or ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_FORWARD) and (i<ModelBlock->Block_List[j].Nb_Recursives))
+ lhs->writeOutput(tmp_output, oMatlabStaticModelSparse, temporary_terms);
+ else*/
+ lhs->writeOutput(tmp_output, oMatlabStaticModelSparse, temporary_terms);
+ switch (ModelBlock->Block_List[j].Simulation_Type)
+ {
+ case EVALUATE_BACKWARD:
+ case EVALUATE_FORWARD:
+evaluation: if (ModelBlock->Block_List[j].Simulation_Type==SOLVE_TWO_BOUNDARIES_COMPLETE || ModelBlock->Block_List[j].Simulation_Type==SOLVE_TWO_BOUNDARIES_SIMPLE)
+ output << " % equation " << ModelBlock->Block_List[j].Equation[i]+1 << " variable : " << sModel
+ << " (" << ModelBlock->Block_List[j].Variable[i]+1 << ") " << block_triangular.c_Equation_Type(ModelBlock->Block_List[j].Equation_Type[i]) << endl;
+ output << " ";
+ if (ModelBlock->Block_List[j].Equation_Type[i] == E_EVALUATE)
+ {
+ output << tmp_output.str();
+ output << " = ";
+ rhs->writeOutput(output, oMatlabStaticModelSparse, temporary_terms);
+ }
+ else if (ModelBlock->Block_List[j].Equation_Type[i] == E_EVALUATE_S)
+ {
+ output << "%" << tmp_output.str();
+ output << " = ";
+ if (ModelBlock->Block_List[j].Equation_Normalized[i])
+ {
+ rhs->writeOutput(output, oMatlabStaticModelSparse, temporary_terms);
+ output << "\n ";
+ tmp_output.str("");
+ eq_node = (BinaryOpNode *)ModelBlock->Block_List[j].Equation_Normalized[i];
+ lhs = eq_node->get_arg1();
+ rhs = eq_node->get_arg2();
+ lhs->writeOutput(output, oMatlabStaticModelSparse, temporary_terms);
+ output << " = ";
+ rhs->writeOutput(output, oMatlabStaticModelSparse, temporary_terms);
+ }
+ }
+ else
+ {
+ cerr << "Type missmatch for equation " << ModelBlock->Block_List[j].Equation[i]+1 << "\n";
+ exit(EXIT_FAILURE);
+ }
+ output << ";\n";
+ break;
+ case SOLVE_BACKWARD_SIMPLE:
+ case SOLVE_FORWARD_SIMPLE:
+ case SOLVE_BACKWARD_COMPLETE:
+ case SOLVE_FORWARD_COMPLETE:
+ if (i<ModelBlock->Block_List[j].Nb_Recursives)
+ goto evaluation;
+ feedback_variables.push_back(ModelBlock->Block_List[j].Variable[i]);
+ output << " % equation " << ModelBlock->Block_List[j].Equation[i]+1 << " variable : " << sModel
+ << " (" << ModelBlock->Block_List[j].Variable[i]+1 << ") " << block_triangular.c_Equation_Type(ModelBlock->Block_List[j].Equation_Type[i]) << endl;
+ output << " " << "residual(" << i+1-ModelBlock->Block_List[j].Nb_Recursives << ") = (";
+ goto end;
+ default:
+end:
+ output << tmp_output.str();
+ output << ") - (";
+ rhs->writeOutput(output, oMatlabStaticModelSparse, temporary_terms);
+ output << ");\n";
+ }
+ }
+ // The Jacobian if we have to solve the block
+ output << " " << sps << "% Jacobian " << endl;
+ switch (ModelBlock->Block_List[j].Simulation_Type)
+ {
+ case EVALUATE_BACKWARD:
+ case EVALUATE_FORWARD:
+ break;
+ case SOLVE_BACKWARD_SIMPLE:
+ case SOLVE_FORWARD_SIMPLE:
+ case SOLVE_BACKWARD_COMPLETE:
+ case SOLVE_FORWARD_COMPLETE:
+ for(i=0; i<(int)ModelBlock->Block_List[j].Chain_Rule_Derivatives->size();i++)
+ {
+ pair< pair<int, pair<int, int> >, pair<int, int> > it = ModelBlock->Block_List[j].Chain_Rule_Derivatives->at(i);
+ k=it.first.first;
+ int eq=it.first.second.first;
+ int var=it.first.second.second;
+ int eqr=it.second.first;
+ int varr=it.second.second;
+ output << " g1(" << eq+1-ModelBlock->Block_List[j].Nb_Recursives << ", "
+ << var+1-ModelBlock->Block_List[j].Nb_Recursives << ") = ";
+ writeChainRuleDerivative(output, eqr, varr, k, oMatlabStaticModelSparse, temporary_terms);
+ output << "; % variable=" << symbol_table.getName(symbol_table.getID(eEndogenous, varr))
+ << " " << varr+1 << ", equation=" << eqr+1 << endl;
+ }
+ break;
+ default:
+ break;
+ }
+ output.close();
+ }
+ }
+
+void
+StaticDllModel::writeModelEquationsCodeOrdered(const string file_name, const Model_Block *ModelBlock, const string bin_basename, map_idx_type map_idx) const
+ {
+ struct Uff_l
+ {
+ int u, var, lag;
+ Uff_l *pNext;
+ };
+
+ struct Uff
+ {
+ Uff_l *Ufl, *Ufl_First;
+ };
+
+ int i,j,k,v;
+ string tmp_s;
+ ostringstream tmp_output;
+ ofstream code_file;
+ NodeID lhs=NULL, rhs=NULL;
+ BinaryOpNode *eq_node;
+ Uff Uf[symbol_table.endo_nbr()];
+ map<NodeID, int> reference_count;
+ vector<int> feedback_variables;
+ bool file_open=false;
+ string main_name=file_name;
+ main_name+=".cod";
+ code_file.open(main_name.c_str(), ios::out | ios::binary | ios::ate );
+ if (!code_file.is_open())
+ {
+ cout << "Error : Can't open file \"" << main_name << "\" for writing\n";
+ exit(EXIT_FAILURE);
+ }
+ //Temporary variables declaration
+ code_file.write(&FDIMST, sizeof(FDIMST));
+ k=temporary_terms.size();
+ code_file.write(reinterpret_cast<char *>(&k),sizeof(k));
+
+ for (j = 0; j < ModelBlock->Size ;j++)
+ {
+ feedback_variables.clear();
+ if (j>0)
+ code_file.write(&FENDBLOCK, sizeof(FENDBLOCK));
+ code_file.write(&FBEGINBLOCK, sizeof(FBEGINBLOCK));
+ v=ModelBlock->Block_List[j].Size - ModelBlock->Block_List[j].Nb_Recursives;
+ code_file.write(reinterpret_cast<char *>(&v),sizeof(v));
+ v=ModelBlock->Block_List[j].Simulation_Type;
+ code_file.write(reinterpret_cast<char *>(&v),sizeof(v));
+ int count_u;
+ for (i=ModelBlock->Block_List[j].Nb_Recursives; i < ModelBlock->Block_List[j].Size;i++)
+ {
+ code_file.write(reinterpret_cast<char *>(&ModelBlock->Block_List[j].Variable[i]),sizeof(ModelBlock->Block_List[j].Variable[i]));
+ code_file.write(reinterpret_cast<char *>(&ModelBlock->Block_List[j].Equation[i]),sizeof(ModelBlock->Block_List[j].Equation[i]));
+ code_file.write(reinterpret_cast<char *>(&ModelBlock->Block_List[j].Own_Derivative[i]),sizeof(ModelBlock->Block_List[j].Own_Derivative[i]));
+ }
+ if (ModelBlock->Block_List[j].Simulation_Type==SOLVE_BACKWARD_COMPLETE || ModelBlock->Block_List[j].Simulation_Type==SOLVE_FORWARD_COMPLETE)
+ {
+ int u_count_int=0;
+ Write_Inf_To_Bin_File(file_name, bin_basename, j, u_count_int,file_open);
+ code_file.write(reinterpret_cast<char *>(&ModelBlock->Block_List[j].is_linear),sizeof(ModelBlock->Block_List[j].is_linear));
+ v = u_count_int ;
+ code_file.write(reinterpret_cast<char *>(&v),sizeof(v));
+ v=symbol_table.endo_nbr();
+ code_file.write(reinterpret_cast<char *>(&v),sizeof(v));
+ v=block_triangular.ModelBlock->Block_List[j].Max_Lag;
+ code_file.write(reinterpret_cast<char *>(&v),sizeof(v));
+ v=block_triangular.ModelBlock->Block_List[j].Max_Lead;
+ code_file.write(reinterpret_cast<char *>(&v),sizeof(v));
+
+ v=u_count_int;
+ code_file.write(reinterpret_cast<char *>(&v),sizeof(v));
+ file_open=true;
+ }
+ // The equations
+ //cout << block_triangular.BlockSim(ModelBlock->Block_List[j].Simulation_Type) << " j=" << j << endl;
+ for (i = 0;i < ModelBlock->Block_List[j].Size;i++)
+ {
+ //The Temporary terms
+ //cout << "equation = " << ModelBlock->Block_List[j].Equation[i] << " variable = " << ModelBlock->Block_List[j].Variable[i] << " r[" << i << "] " << block_triangular.c_Equation_Type(ModelBlock->Block_List[j].Equation_Type[i]) << endl;
+ temporary_terms_type tt2;
+ tt2.clear();
+ for (temporary_terms_type::const_iterator it = ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->begin();
+ it != ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->end(); it++)
+ {
+ (*it)->compile(code_file, false, tt2, map_idx, false);
+ code_file.write(&FSTPST, sizeof(FSTPST));
+ map_idx_type::const_iterator ii=map_idx.find((*it)->idx);
+ v=(int)ii->second;
+ code_file.write(reinterpret_cast<char *>(&v), sizeof(v));
+ // Insert current node into tt2
+ tt2.insert(*it);
+ }
+ switch (ModelBlock->Block_List[j].Simulation_Type)
+ {
+evaluation:
+ case EVALUATE_BACKWARD:
+ case EVALUATE_FORWARD:
+ if (ModelBlock->Block_List[j].Equation_Type[i] == E_EVALUATE)
+ {
+ eq_node = equations[ModelBlock->Block_List[j].Equation[i]];
+ lhs = eq_node->get_arg1();
+ rhs = eq_node->get_arg2();
+ rhs->compile(code_file, false, temporary_terms, map_idx, false);
+ lhs->compile(code_file, true, temporary_terms, map_idx, false);
+ }
+ else if (ModelBlock->Block_List[j].Equation_Type[i] == E_EVALUATE_S)
+ {
+ eq_node = (BinaryOpNode*)ModelBlock->Block_List[j].Equation_Normalized[i];
+ lhs = eq_node->get_arg1();
+ rhs = eq_node->get_arg2();
+ rhs->compile(code_file, false, temporary_terms, map_idx, false);
+ lhs->compile(code_file, true, temporary_terms, map_idx, false);
+ }
+ break;
+ case SOLVE_BACKWARD_COMPLETE:
+ case SOLVE_FORWARD_COMPLETE:
+ if (i<ModelBlock->Block_List[j].Nb_Recursives)
+ goto evaluation;
+ feedback_variables.push_back(ModelBlock->Block_List[j].Variable[i]);
+ v=ModelBlock->Block_List[j].Equation[i];
+ Uf[v].Ufl=NULL;
+ goto end;
+ default:
+end:
+ eq_node = equations[ModelBlock->Block_List[j].Equation[i]];
+ lhs = eq_node->get_arg1();
+ rhs = eq_node->get_arg2();
+ lhs->compile(code_file, false, temporary_terms, map_idx, false);
+ rhs->compile(code_file, false, temporary_terms, map_idx, false);
+ code_file.write(&FBINARY, sizeof(FBINARY));
+ int v=oMinus;
+ code_file.write(reinterpret_cast<char *>(&v),sizeof(v));
+ code_file.write(&FSTPR, sizeof(FSTPR));
+ v = i - ModelBlock->Block_List[j].Nb_Recursives;
+ code_file.write(reinterpret_cast<char *>(&v), sizeof(v));
+ }
+ }
+ code_file.write(&FENDEQU, sizeof(FENDEQU));
+ // The Jacobian if we have to solve the block
+ if (ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_BACKWARD
+ && ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_FORWARD)
+ {
+ switch (ModelBlock->Block_List[j].Simulation_Type)
+ {
+ case SOLVE_BACKWARD_SIMPLE:
+ case SOLVE_FORWARD_SIMPLE:
+ compileDerivative(code_file, ModelBlock->Block_List[j].Equation[0], ModelBlock->Block_List[j].Variable[0], 0, map_idx);
+ code_file.write(&FSTPG, sizeof(FSTPG));
+ v=0;
+ code_file.write(reinterpret_cast<char *>(&v), sizeof(v));
+ break;
+
+ case SOLVE_BACKWARD_COMPLETE:
+ case SOLVE_FORWARD_COMPLETE:
+ count_u = feedback_variables.size();
+ for(i=0; i<(int)ModelBlock->Block_List[j].Chain_Rule_Derivatives->size();i++)
+ {
+ pair< pair<int, pair<int, int> >, pair<int, int> > it = ModelBlock->Block_List[j].Chain_Rule_Derivatives->at(i);
+ k=it.first.first;
+ int eq=it.first.second.first;
+ int var=it.first.second.second;
+ int eqr=it.second.first;
+ int varr=it.second.second;
+ int v=ModelBlock->Block_List[j].Equation[eq];
+ if(eq>=ModelBlock->Block_List[j].Nb_Recursives and var>=ModelBlock->Block_List[j].Nb_Recursives)
+ {
+ if (!Uf[v].Ufl)
+ {
+ Uf[v].Ufl=(Uff_l*)malloc(sizeof(Uff_l));
+ Uf[v].Ufl_First=Uf[v].Ufl;
+ }
+ else
+ {
+ Uf[v].Ufl->pNext=(Uff_l*)malloc(sizeof(Uff_l));
+ Uf[v].Ufl=Uf[v].Ufl->pNext;
+ }
+ Uf[v].Ufl->pNext=NULL;
+ Uf[v].Ufl->u=count_u;
+ Uf[v].Ufl->var=varr;
+ Uf[v].Ufl->lag=k;
+ compileChainRuleDerivative(code_file, eqr, varr, k, map_idx);
+ code_file.write(&FSTPSU, sizeof(FSTPSU));
+ code_file.write(reinterpret_cast<char *>(&count_u), sizeof(count_u));
+ count_u++;
+ }
+ }
+ for (i = 0;i < ModelBlock->Block_List[j].Size;i++)
+ {
+ if(i>=ModelBlock->Block_List[j].Nb_Recursives)
+ {
+ code_file.write(&FLDR, sizeof(FLDR));
+ v = i-ModelBlock->Block_List[j].Nb_Recursives;
+ code_file.write(reinterpret_cast<char *>(&v), sizeof(v));
+ code_file.write(&FLDZ, sizeof(FLDZ));
+ v=ModelBlock->Block_List[j].Equation[i];
+ for (Uf[v].Ufl=Uf[v].Ufl_First; Uf[v].Ufl; Uf[v].Ufl=Uf[v].Ufl->pNext)
+ {
+ code_file.write(&FLDSU, sizeof(FLDSU));
+ code_file.write(reinterpret_cast<char *>(&Uf[v].Ufl->u), sizeof(Uf[v].Ufl->u));
+ code_file.write(&FLDSV, sizeof(FLDSV));
+ char vc=eEndogenous;
+ code_file.write(reinterpret_cast<char *>(&vc), sizeof(vc));
+ int v1=Uf[v].Ufl->var;
+ code_file.write(reinterpret_cast<char *>(&v1), sizeof(v1));
+ code_file.write(&FBINARY, sizeof(FBINARY));
+ v1=oTimes;
+ code_file.write(reinterpret_cast<char *>(&v1), sizeof(v1));
+ code_file.write(&FCUML, sizeof(FCUML));
+ }
+ Uf[v].Ufl=Uf[v].Ufl_First;
+ while (Uf[v].Ufl)
+ {
+ Uf[v].Ufl_First=Uf[v].Ufl->pNext;
+ free(Uf[v].Ufl);
+ Uf[v].Ufl=Uf[v].Ufl_First;
+ }
+ code_file.write(&FBINARY, sizeof(FBINARY));
+ v=oMinus;
+ code_file.write(reinterpret_cast<char *>(&v), sizeof(v));
+ code_file.write(&FSTPSU, sizeof(FSTPSU));
+ v = i - ModelBlock->Block_List[j].Nb_Recursives;
+ code_file.write(reinterpret_cast<char *>(&v), sizeof(v));
+ }
+ }
+ break;
+ default:
+ break;
+ }
+ }
+ }
+ code_file.write(&FENDBLOCK, sizeof(FENDBLOCK));
+ code_file.write(&FEND, sizeof(FEND));
+ code_file.close();
+ }
+
+
+
+void
+StaticDllModel::Write_Inf_To_Bin_File(const string &static_basename, const string &bin_basename, const int &num,
+ int &u_count_int, bool &file_open) const
+ {
+ int j;
+ std::ofstream SaveCode;
+ if (file_open)
+ SaveCode.open((bin_basename + "_static.bin").c_str(), ios::out | ios::in | ios::binary | ios ::ate );
+ else
+ SaveCode.open((bin_basename + "_static.bin").c_str(), ios::out | ios::binary);
+ if (!SaveCode.is_open())
+ {
+ cout << "Error : Can't open file \"" << bin_basename << "_static.bin\" for writing\n";
+ exit(EXIT_FAILURE);
+ }
+ u_count_int=0;
+ int Size = block_triangular.ModelBlock->Block_List[num].Size - block_triangular.ModelBlock->Block_List[num].Nb_Recursives;
+ for(int i=0; i<(int)block_triangular.ModelBlock->Block_List[num].Chain_Rule_Derivatives->size();i++)
+ {
+ //Chain_Rule_Derivatives.insert(make_pair( make_pair(eq, eqr), make_pair(var, make_pair(varr, lag))));
+ pair< pair<int, pair<int, int> >, pair<int, int> > it = block_triangular.ModelBlock->Block_List[num].Chain_Rule_Derivatives->at(i);
+ int k=it.first.first;
+ int eq=it.first.second.first;
+
+ int var_init=it.first.second.second;
+ /*int eqr=it.second.first;
+ int varr=it.second.second;*/
+ if(eq>=block_triangular.ModelBlock->Block_List[num].Nb_Recursives and var_init>=block_triangular.ModelBlock->Block_List[num].Nb_Recursives)
+ {
+ int v=eq-block_triangular.ModelBlock->Block_List[num].Nb_Recursives;
+ SaveCode.write(reinterpret_cast<char *>(&v), sizeof(v));
+ int var=it.first.second.second-block_triangular.ModelBlock->Block_List[num].Nb_Recursives + k * Size;
+ SaveCode.write(reinterpret_cast<char *>(&var), sizeof(var));
+ SaveCode.write(reinterpret_cast<char *>(&k), sizeof(k));
+ int u = u_count_int + Size;
+ SaveCode.write(reinterpret_cast<char *>(&u), sizeof(u));
+ //cout << "eq=" << v << ", var=" << var << ", lag=" << k << " u=" << u << "\n";
+ u_count_int++;
+ }
+ }
+ /*cout << "u_count_int=" << u_count_int << endl;
+ cout << "block_triangular.ModelBlock->Block_List[" << num << "].Nb_Recursives=" << block_triangular.ModelBlock->Block_List[num].Nb_Recursives << " block_triangular.ModelBlock->Block_List[" << num << "].Size=" << block_triangular.ModelBlock->Block_List[num].Size << endl;*/
+ for (j=block_triangular.ModelBlock->Block_List[num].Nb_Recursives;j<block_triangular.ModelBlock->Block_List[num].Size;j++)
+ {
+ int varr=block_triangular.ModelBlock->Block_List[num].Variable[j];
+ //cout << "j=" << j << " varr=" << varr << "\n";
+ SaveCode.write(reinterpret_cast<char *>(&varr), sizeof(varr));
+ }
+ for (j=block_triangular.ModelBlock->Block_List[num].Nb_Recursives;j<block_triangular.ModelBlock->Block_List[num].Size;j++)
+ {
+ int eqr1=block_triangular.ModelBlock->Block_List[num].Equation[j];
+ SaveCode.write(reinterpret_cast<char *>(&eqr1), sizeof(eqr1));
+ }
+ SaveCode.close();
+ }
+
+void
+StaticDllModel::writeSparseStaticMFile(const string &static_basename, const string &basename, const int mode) const
+ {
+ string sp;
+ ofstream mStaticModelFile;
+ ostringstream tmp, tmp1, tmp_eq;
+ int prev_Simulation_Type;
+ //SymbolicGaussElimination SGE;
+ bool OK;
+ chdir(basename.c_str());
+ string filename = static_basename + ".m";
+ mStaticModelFile.open(filename.c_str(), ios::out | ios::binary);
+ if (!mStaticModelFile.is_open())
+ {
+ cerr << "Error: Can't open file " << filename << " for writing" << endl;
+ exit(EXIT_FAILURE);
+ }
+ mStaticModelFile << "%\n";
+ mStaticModelFile << "% " << filename << " : Computes static model for Dynare\n";
+ mStaticModelFile << "%\n";
+ mStaticModelFile << "% Warning : this file is generated automatically by Dynare\n";
+ mStaticModelFile << "% from model file (.mod)\n\n";
+ mStaticModelFile << "%/\n";
+
+ int i, k;
+
+ mStaticModelFile << "function [varargout] = " << static_basename << "(varargin)\n";
+ mStaticModelFile << " global oo_ options_ M_ ;\n";
+ //Temporary variables declaration
+ OK=true;
+ ostringstream tmp_output;
+ for (temporary_terms_type::const_iterator it = temporary_terms.begin();
+ it != temporary_terms.end(); it++)
+ {
+ if (OK)
+ OK=false;
+ else
+ tmp_output << " ";
+ (*it)->writeOutput(tmp_output, oMatlabStaticModelSparse, temporary_terms);
+ }
+ if (tmp_output.str().length()>0)
+ mStaticModelFile << " global " << tmp_output.str() << " M_ ;\n";
+
+ mStaticModelFile << " T_init=0;\n";
+ tmp_output.str("");
+ for (temporary_terms_type::const_iterator it = temporary_terms.begin();
+ it != temporary_terms.end(); it++)
+ {
+ tmp_output << " ";
+ (*it)->writeOutput(tmp_output, oMatlabStaticModel, temporary_terms);
+ tmp_output << "=T_init;\n";
+ }
+ if (tmp_output.str().length()>0)
+ mStaticModelFile << tmp_output.str();
+
+ mStaticModelFile << " y_size=M_.endo_nbr;\n";
+ mStaticModelFile << " %it is the deterministic simulation of the block decomposed static model\n";
+ mStaticModelFile << " maxit_=options_.maxit_;\n";
+ mStaticModelFile << " solve_tolf=options_.solve_tolf;\n";
+ mStaticModelFile << " y=oo_.steady_state';\n";
+ mStaticModelFile << " x=[oo_.exo_steady_state oo_.exo_det_steady_state];\n";
+
+ prev_Simulation_Type=-1;
+ mStaticModelFile << " params=M_.params;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.status = 0;\n";
+ for (i = 0;i < block_triangular.ModelBlock->Size;i++)
+ {
+ k = block_triangular.ModelBlock->Block_List[i].Simulation_Type;
+ if ((k == EVALUATE_FORWARD) && (block_triangular.ModelBlock->Block_List[i].Size))
+ {
+ mStaticModelFile << " oo_.deterministic_simulation.status = 1;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.error = 0;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.iterations = 0;\n";
+ mStaticModelFile << " if(isfield(oo_.deterministic_simulation,'block'))\n";
+ mStaticModelFile << " blck_num = length(oo_.deterministic_simulation.block)+1;\n";
+ mStaticModelFile << " else\n";
+ mStaticModelFile << " blck_num = 1;\n";
+ mStaticModelFile << " end;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.block(blck_num).status = 1;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.block(blck_num).error = 0;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.block(blck_num).iterations = 0;\n";
+ //mStaticModelFile << " for it_ = y_kmin+1:(periods+y_kmin)\n";
+ mStaticModelFile << " y=" << static_basename << "_" << i + 1 << "(y, x, params);\n";
+ mStaticModelFile << " tmp = y(:,M_.block_structure.block(" << i + 1 << ").variable);\n";
+ mStaticModelFile << " if(isnan(tmp) | isinf(tmp))\n";
+ mStaticModelFile << " disp(['Inf or Nan value during the evaluation of block " << i <<"']);\n";
+ mStaticModelFile << " return;\n";
+ mStaticModelFile << " end;\n";
+ }
+ else if ((k == EVALUATE_BACKWARD) && (block_triangular.ModelBlock->Block_List[i].Size))
+ {
+ mStaticModelFile << " oo_.deterministic_simulation.status = 1;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.error = 0;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.iterations = 0;\n";
+ mStaticModelFile << " if(isfield(oo_.deterministic_simulation,'block'))\n";
+ mStaticModelFile << " blck_num = length(oo_.deterministic_simulation.block)+1;\n";
+ mStaticModelFile << " else\n";
+ mStaticModelFile << " blck_num = 1;\n";
+ mStaticModelFile << " end;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.block(blck_num).status = 1;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.block(blck_num).error = 0;\n";
+ mStaticModelFile << " oo_.deterministic_simulation.block(blck_num).iterations = 0;\n";
+ mStaticModelFile << " " << static_basename << "_" << i + 1 << "(y, x, params);\n";
+ mStaticModelFile << " tmp = y(:,M_.block_structure.block(" << i + 1 << ").variable);\n";
+ mStaticModelFile << " if(isnan(tmp) | isinf(tmp))\n";
+ mStaticModelFile << " disp(['Inf or Nan value during the evaluation of block " << i <<"']);\n";
+ mStaticModelFile << " return;\n";
+ mStaticModelFile << " end;\n";
+ }
+ else if ((k == SOLVE_FORWARD_COMPLETE || k == SOLVE_FORWARD_SIMPLE) && (block_triangular.ModelBlock->Block_List[i].Size))
+ {
+ /*mStaticModelFile << " g1=0;\n";
+ mStaticModelFile << " r=0;\n";*/
+ tmp.str("");
+ for (int ik=block_triangular.ModelBlock->Block_List[i].Nb_Recursives ;ik<block_triangular.ModelBlock->Block_List[i].Size;ik++)
+ {
+ tmp << " " << block_triangular.ModelBlock->Block_List[i].Variable[ik]+1;
+ }
+ mStaticModelFile << " y_index = [" << tmp.str() << "];\n";
+ /*int nze, m;
+ for (nze=0,m=0;m<=block_triangular.ModelBlock->Block_List[i].Max_Lead+block_triangular.ModelBlock->Block_List[i].Max_Lag;m++)*/
+ int nze=block_triangular.ModelBlock->Block_List[i].IM_lead_lag[0].size;
+ mStaticModelFile << " if(isfield(oo_.deterministic_simulation,'block'))\n";
+ mStaticModelFile << " blck_num = length(oo_.deterministic_simulation.block)+1;\n";
+ mStaticModelFile << " else\n";
+ mStaticModelFile << " blck_num = 1;\n";
+ mStaticModelFile << " end;\n";
+ mStaticModelFile << " y = solve_one_boundary('" << static_basename << "_" << i + 1 << "'" <<
+ ", y, x, params, y_index, " << nze <<
+ ", 1, " << block_triangular.ModelBlock->Block_List[i].is_linear <<
+ ", blck_num, 0, options_.maxit_, options_.solve_tolf, options_.slowc, options_.cutoff, options_.simulation_method, 1, 0, 0);\n";
+ mStaticModelFile << " tmp = y(:,M_.block_structure.block(" << i + 1 << ").variable);\n";
+ mStaticModelFile << " if(isnan(tmp) | isinf(tmp))\n";
+ mStaticModelFile << " disp(['Inf or Nan value during the resolution of block " << i <<"']);\n";
+ mStaticModelFile << " return;\n";
+ mStaticModelFile << " end;\n";
+ }
+ else if ((k == SOLVE_BACKWARD_COMPLETE || k == SOLVE_BACKWARD_SIMPLE) && (block_triangular.ModelBlock->Block_List[i].Size))
+ {
+ /*mStaticModelFile << " g1=0;\n";
+ mStaticModelFile << " r=0;\n";*/
+ tmp.str("");
+ for (int ik=block_triangular.ModelBlock->Block_List[i].Nb_Recursives ;ik<block_triangular.ModelBlock->Block_List[i].Size;ik++)
+ {
+ tmp << " " << block_triangular.ModelBlock->Block_List[i].Variable[ik]+1;
+ }
+ mStaticModelFile << " y_index = [" << tmp.str() << "];\n";
+ /*int nze, m;
+ for (nze=0,m=0;m<=block_triangular.ModelBlock->Block_List[i].Max_Lead+block_triangular.ModelBlock->Block_List[i].Max_Lag;m++)*/
+ int nze=block_triangular.ModelBlock->Block_List[i].IM_lead_lag[0].size;
+ mStaticModelFile << " if(isfield(oo_.deterministic_simulation,'block'))\n";
+ mStaticModelFile << " blck_num = length(oo_.deterministic_simulation.block)+1;\n";
+ mStaticModelFile << " else\n";
+ mStaticModelFile << " blck_num = 1;\n";
+ mStaticModelFile << " end;\n";
+ mStaticModelFile << " y = solve_one_boundary('" << static_basename << "_" << i + 1 << "'" <<
+ ", y, x, params, y_index, " << nze <<
+ ", 0, " << block_triangular.ModelBlock->Block_List[i].is_linear <<
+ ", blck_num, y_kmin, options_.maxit_, options_.solve_tolf, options_.slowc, options_.cutoff, options_.simulation_method, 1, 0, 0);\n";
+ mStaticModelFile << " tmp = y(:,M_.block_structure.block(" << i + 1 << ").variable);\n";
+ mStaticModelFile << " if(isnan(tmp) | isinf(tmp))\n";
+ mStaticModelFile << " disp(['Inf or Nan value during the resolution of block " << i <<"']);\n";
+ mStaticModelFile << " return;\n";
+ mStaticModelFile << " end;\n";
+ }
+ }
+ mStaticModelFile << " oo_.steady_state = y';\n";
+ mStaticModelFile << "return;\n";
+
+ mStaticModelFile.close();
+
+ writeModelEquationsOrdered_M( block_triangular.ModelBlock, static_basename);
+
+ chdir("..");
+ }
+
+
+void
+StaticDllModel::writeOutput(ostream &output, const string &basename) const
+{
+ output << "options_.model_mode = " << mode << ";" << endl;
+
+ // Erase possible remnants of previous runs
+ if (mode != eStandardMode)
+ output << "delete('" << basename << "_static.m');" << endl;
+ if (mode != eDLLMode)
+ output << "erase_compiled_function('" + basename + "_static');" << endl;
+
+ // Special setup for DLL or Sparse modes
+ if (mode == eDLLMode)
+ output << "mex -O LDFLAGS='-pthread -shared -Wl,--no-undefined' " << basename << "_static.c" << endl;
+ if (mode == eSparseMode || mode == eSparseDLLMode)
+ output << "addpath " << basename << ";" << endl;
+}
+
+void
+StaticDllModel::writeOutputPostComputing(ostream &output, const string &basename) const
+{
+ if (mode == eSparseMode || mode == eSparseDLLMode)
+ output << "rmpath " << basename << ";" << endl;
+}
+
+void
+StaticDllModel::evaluateJacobian(const eval_context_type &eval_context, jacob_map *j_m, bool dynamic)
+{
+ int i=0;
+ int j=0;
+ bool *IM=NULL;
+ int a_variable_lag=-9999;
+ for (first_derivatives_type::iterator it = first_derivatives.begin();
+ it != first_derivatives.end(); it++)
+ {
+ //cout << "it->first.second=" << it->first.second << " variable_table.getSymbolID(it->first.second)=" << variable_table.getSymbolID(it->first.second) << " Type=" << variable_table.getType(it->first.second) << " eEndogenous=" << eEndogenous << " eExogenous=" << eExogenous << " variable_table.getLag(it->first.second)=" << variable_table.getLag(it->first.second) << "\n";
+ if (getTypeByDerivID(it->first.second) == eEndogenous)
+ {
+ NodeID Id = it->second;
+ double val = 0;
+ try
+ {
+ val = Id->eval(eval_context);
+ }
+ catch (ExprNode::EvalException &e)
+ {
+ cout << "evaluation of Jacobian failed for equation " << it->first.first+1 << " and variable " << symbol_table.getName(getSymbIDByDerivID(it->first.second)) << "(" << getLagByDerivID(it->first.second) << ") [" << getSymbIDByDerivID(it->first.second) << "] !" << endl;
+ Id->writeOutput(cout, oMatlabStaticModelSparse, temporary_terms);
+ cout << "\n";
+ cerr << "StaticDllModel::evaluateJacobian: evaluation of Jacobian failed for equation " << it->first.first+1 << " and variable " << symbol_table.getName(getSymbIDByDerivID(it->first.second)) << "(" << getLagByDerivID(it->first.second) << ")!" << endl;
+ }
+ int eq=it->first.first;
+ int var = symbol_table.getTypeSpecificID(getSymbIDByDerivID(it->first.second));///symbol_table.getID(eEndogenous,it->first.second);//variable_table.getSymbolID(it->first.second);
+ int k1 = getLagByDerivID(it->first.second);
+ if (a_variable_lag!=k1)
+ {
+ IM=block_triangular.incidencematrix.Get_IM(k1, eEndogenous);
+ a_variable_lag=k1;
+ }
+ if (k1==0 or !dynamic)
+ {
+ j++;
+ (*j_m)[make_pair(eq,var)]+=val;
+ }
+ if (IM[eq*symbol_table.endo_nbr()+var] && (fabs(val) < cutoff))
+ {
+ if (block_triangular.bt_verbose)
+ cout << "the coefficient related to variable " << var << " with lag " << k1 << " in equation " << eq << " is equal to " << val << " and is set to 0 in the incidence matrix (size=" << symbol_table.endo_nbr() << ")\n";
+ block_triangular.incidencematrix.unfill_IM(eq, var, k1, eEndogenous);
+ i++;
+ }
+ }
+ }
+ //Get ride of the elements of the incidence matrix equal to Zero
+ IM=block_triangular.incidencematrix.Get_IM(0, eEndogenous);
+ for (int i=0;i<symbol_table.endo_nbr();i++)
+ for (int j=0;j<symbol_table.endo_nbr();j++)
+ if (IM[i*symbol_table.endo_nbr()+j])
+ if (first_derivatives.find(make_pair(i,getDerivID(symbol_table.getID(eEndogenous, j), 0)))==first_derivatives.end())
+ block_triangular.incidencematrix.unfill_IM(i, j, 0, eEndogenous);
+ if (i>0)
+ {
+ cout << i << " elements among " << first_derivatives.size() << " in the incidence matrices are below the cutoff (" << cutoff << ") and are discarded\n";
+ cout << "the contemporaneous incidence matrix has " << j << " elements\n";
+ }
+}
+
+void
+StaticDllModel::BlockLinear(Model_Block *ModelBlock)
+{
+ int i,j,l,m,ll;
+ for (j = 0;j < ModelBlock->Size;j++)
+ {
+ if (ModelBlock->Block_List[j].Simulation_Type==SOLVE_BACKWARD_COMPLETE ||
+ ModelBlock->Block_List[j].Simulation_Type==SOLVE_FORWARD_COMPLETE)
+ {
+ ll=ModelBlock->Block_List[j].Max_Lag;
+ for (i=0;i<ModelBlock->Block_List[j].IM_lead_lag[ll].size;i++)
+ {
+ int eq=ModelBlock->Block_List[j].IM_lead_lag[ll].Equ_Index[i];
+ int var=ModelBlock->Block_List[j].IM_lead_lag[ll].Var_Index[i];
+ //first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,variable_table.getID(var,0)));
+ first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,getDerivID(symbol_table.getID(eEndogenous, var),0)));
+ if (it!= first_derivatives.end())
+ {
+ NodeID Id = it->second;
+ set<pair<int, int> > endogenous;
+ Id->collectEndogenous(endogenous);
+ if (endogenous.size() > 0)
+ {
+ for (l=0;l<ModelBlock->Block_List[j].Size;l++)
+ {
+ if (endogenous.find(make_pair(ModelBlock->Block_List[j].Variable[l], 0)) != endogenous.end())
+ {
+ ModelBlock->Block_List[j].is_linear=false;
+ goto follow;
+ }
+ }
+ }
+ }
+ }
+ }
+ else if (ModelBlock->Block_List[j].Simulation_Type==SOLVE_TWO_BOUNDARIES_COMPLETE || ModelBlock->Block_List[j].Simulation_Type==SOLVE_TWO_BOUNDARIES_SIMPLE)
+ {
+ for (m=0;m<=ModelBlock->Block_List[j].Max_Lead+ModelBlock->Block_List[j].Max_Lag;m++)
+ {
+ int k1=m-ModelBlock->Block_List[j].Max_Lag;
+ for (i=0;i<ModelBlock->Block_List[j].IM_lead_lag[m].size;i++)
+ {
+ int eq=ModelBlock->Block_List[j].IM_lead_lag[m].Equ_Index[i];
+ int var=ModelBlock->Block_List[j].IM_lead_lag[m].Var_Index[i];
+ //first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,variable_table.getID(var,k1)));
+ first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,getDerivID(symbol_table.getID(eEndogenous, var),k1)));
+ NodeID Id = it->second;
+ if (it!= first_derivatives.end())
+ {
+ set<pair<int, int> > endogenous;
+ Id->collectEndogenous(endogenous);
+ if (endogenous.size() > 0)
+ {
+ for (l=0;l<ModelBlock->Block_List[j].Size;l++)
+ {
+ if (endogenous.find(make_pair(ModelBlock->Block_List[j].Variable[l], k1)) != endogenous.end())
+ {
+ ModelBlock->Block_List[j].is_linear=false;
+ goto follow;
+ }
+ }
+ }
+ }
+ }
+ }
+ }
+follow:
+ i=0;
+ }
+}
+
+
+map<pair<int, pair<int, int > >, NodeID>
+StaticDllModel::collect_first_order_derivatives_endogenous()
+{
+ map<pair<int, pair<int, int > >, NodeID> endo_derivatives;
+ for (first_derivatives_type::iterator it2 = first_derivatives.begin();
+ it2 != first_derivatives.end(); it2++)
+ {
+ if (getTypeByDerivID(it2->first.second)==eEndogenous)
+ {
+ int eq = it2->first.first;
+ int var=symbol_table.getTypeSpecificID(getSymbIDByDerivID(it2->first.second));
+ int lag=getLagByDerivID(it2->first.second);
+ //if (lag==0)
+ endo_derivatives[make_pair(eq, make_pair(var, lag))] = it2->second;
+ }
+ }
+ return endo_derivatives;
+}
+
+
+
+void
+StaticDllModel::computingPass(const eval_context_type &eval_context, bool no_tmp_terms)
+{
+
+ // Computes static jacobian columns
+ computeStatJacobianCols();
+
+ // Compute derivatives w.r. to all endogenous, and possibly exogenous and exogenous deterministic
+ set<int> vars;
+ for (deriv_id_table_t::const_iterator it = deriv_id_table.begin();
+ it != deriv_id_table.end(); it++)
+ {
+ SymbolType type = symbol_table.getType(it->first.first);
+ if (type == eEndogenous)
+ vars.insert(it->second);
+ }
+
+ // Launch computations
+ cout << "Computing static model derivatives:" << endl
+ << " - order 1" << endl;
+ computeJacobian(vars);
+ //cout << "mode=" << mode << " eSparseDLLMode=" << eSparseDLLMode << " eSparseMode=" << eSparseMode << "\n";
+ if (mode == eSparseDLLMode || mode == eSparseMode)
+ {
+ BuildIncidenceMatrix();
+
+ jacob_map j_m;
+ evaluateJacobian(eval_context, &j_m, true);
+
+
+ if (block_triangular.bt_verbose)
+ {
+ cout << "The gross incidence matrix \n";
+ block_triangular.incidencematrix.Print_IM(eEndogenous);
+ }
+ t_etype equation_simulation_type;
+ map<pair<int, pair<int, int> >, NodeID> first_order_endo_derivatives = collect_first_order_derivatives_endogenous();
+
+ block_triangular.Normalize_and_BlockDecompose_Static_0_Model(j_m, equations, equation_simulation_type, first_order_endo_derivatives, mfs, cutoff);
+ /*for (int j = 0;j < block_triangular.ModelBlock->Size;j++)
+ {
+ for (int i = 0;i < block_triangular.ModelBlock->Block_List[j].Size;i++)
+ {
+ if (i<block_triangular.ModelBlock->Block_List[j].Nb_Recursives )
+ cout << "block=" << j << " R i=" << i << " equation=" << block_triangular.ModelBlock->Block_List[j].Equation[i]+1 << " variable=" << block_triangular.ModelBlock->Block_List[j].Variable[i]+1 << endl;
+ else
+ cout << "block=" << j << " S i=" << i << " equation=" << block_triangular.ModelBlock->Block_List[j].Equation[i]+1 << " variable=" << block_triangular.ModelBlock->Block_List[j].Variable[i]+1 << endl;
+ }
+ }*/
+
+ BlockLinear(block_triangular.ModelBlock);
+
+ computeChainRuleJacobian(block_triangular.ModelBlock);
+
+ if (!no_tmp_terms)
+ computeTemporaryTermsOrdered(block_triangular.ModelBlock);
+
+ }
+ else
+ if (!no_tmp_terms)
+ computeTemporaryTerms(mode == eStandardMode);
+}
+
+void
+StaticDllModel::writeStaticFile(const string &basename) const
+ {
+ int r;
+ switch (mode)
+ {
+ case eStandardMode:
+ break;
+ case eSparseMode:
+#ifdef _WIN32
+ r = mkdir(basename.c_str());
+#else
+ r = mkdir(basename.c_str(), 0777);
+#endif
+ if (r < 0 && errno != EEXIST)
+ {
+ perror("ERROR");
+ exit(EXIT_FAILURE);
+ }
+ writeSparseStaticMFile(basename + "_static", basename, mode);
+ block_triangular.Free_Block(block_triangular.ModelBlock);
+ block_triangular.incidencematrix.Free_IM();
+ //block_triangular.Free_IM_X(block_triangular.First_IM_X);
+ break;
+ case eDLLMode:
+ break;
+ case eSparseDLLMode:
+ // create a directory to store all the files
+#ifdef _WIN32
+ r = mkdir(basename.c_str());
+#else
+ r = mkdir(basename.c_str(), 0777);
+#endif
+ if (r < 0 && errno != EEXIST)
+ {
+ perror("ERROR");
+ exit(EXIT_FAILURE);
+ }
+ writeModelEquationsCodeOrdered(basename + "_static", block_triangular.ModelBlock, basename, map_idx);
+ block_triangular.Free_Block(block_triangular.ModelBlock);
+ block_triangular.incidencematrix.Free_IM();
+ //block_triangular.Free_IM_X(block_triangular.First_IM_X);
+ break;
+ }
+ }
+
+
+int
+StaticDllModel::computeDerivID(int symb_id, int lag)
+{
+ // Check if static variable already has a deriv_id
+ pair<int, int> key = make_pair(symb_id, lag);
+ deriv_id_table_t::const_iterator it = deriv_id_table.find(key);
+ if (it != deriv_id_table.end())
+ return it->second;
+
+ // Create a new deriv_id
+ int deriv_id = deriv_id_table.size();
+
+ deriv_id_table[key] = deriv_id;
+ inv_deriv_id_table.push_back(key);
+
+ SymbolType type = symbol_table.getType(symb_id);
+
+ if (type == eEndogenous)
+ dynJacobianColsNbr++;
+
+ return deriv_id;
+}
+
+SymbolType
+StaticDllModel::getTypeByDerivID(int deriv_id) const throw (UnknownDerivIDException)
+{
+ return symbol_table.getType(getSymbIDByDerivID(deriv_id));
+}
+
+int
+StaticDllModel::getLagByDerivID(int deriv_id) const throw (UnknownDerivIDException)
+{
+ if (deriv_id < 0 || deriv_id >= (int) inv_deriv_id_table.size())
+ throw UnknownDerivIDException();
+
+ return inv_deriv_id_table[deriv_id].second;
+}
+
+int
+StaticDllModel::getSymbIDByDerivID(int deriv_id) const throw (UnknownDerivIDException)
+{
+ if (deriv_id < 0 || deriv_id >= (int) inv_deriv_id_table.size())
+ throw UnknownDerivIDException();
+
+ return inv_deriv_id_table[deriv_id].first;
+}
+
+int
+StaticDllModel::getDerivID(int symb_id, int lag) const throw (UnknownDerivIDException)
+{
+ deriv_id_table_t::const_iterator it = deriv_id_table.find(make_pair(symb_id, lag));
+ if (it == deriv_id_table.end())
+ throw UnknownDerivIDException();
+ else
+ return it->second;
+}
+
+void
+StaticDllModel::computeStatJacobianCols()
+{
+ /* Sort the static endogenous variables by lexicographic order over (lag, type_specific_symbol_id)
+ and fill the static columns for exogenous and exogenous deterministic */
+ map<pair<int, int>, int> ordered_dyn_endo;
+
+ for (deriv_id_table_t::const_iterator it = deriv_id_table.begin();
+ it != deriv_id_table.end(); it++)
+ {
+ const int &symb_id = it->first.first;
+ const int &lag = it->first.second;
+ const int &deriv_id = it->second;
+ SymbolType type = symbol_table.getType(symb_id);
+ int tsid = symbol_table.getTypeSpecificID(symb_id);
+
+ switch (type)
+ {
+ case eEndogenous:
+ ordered_dyn_endo[make_pair(lag, tsid)] = deriv_id;
+ break;
+ case eExogenous:
+ // At this point, dynJacobianColsNbr contains the number of static endogenous
+ break;
+ case eExogenousDet:
+ // At this point, dynJacobianColsNbr contains the number of static endogenous
+ break;
+ case eParameter:
+ // We don't assign a static jacobian column to parameters
+ break;
+ default:
+ // Shut up GCC
+ cerr << "StaticDllModel::computeStatJacobianCols: impossible case" << endl;
+ exit(EXIT_FAILURE);
+ }
+ }
+
+ // Fill in static jacobian columns for endogenous
+ int sorted_id = 0;
+ for (map<pair<int, int>, int>::const_iterator it = ordered_dyn_endo.begin();
+ it != ordered_dyn_endo.end(); it++)
+ dyn_jacobian_cols_table[it->second] = sorted_id++;
+
+}
+
+int
+StaticDllModel::getDynJacobianCol(int deriv_id) const throw (UnknownDerivIDException)
+{
+ map<int, int>::const_iterator it = dyn_jacobian_cols_table.find(deriv_id);
+ if (it == dyn_jacobian_cols_table.end())
+ throw UnknownDerivIDException();
+ else
+ return it->second;
+}
+
+
+void
+StaticDllModel::computeChainRuleJacobian(Model_Block *ModelBlock)
+{
+ map<int, NodeID> recursive_variables;
+ first_chain_rule_derivatives.clear();
+ for(int blck = 0; blck<ModelBlock->Size; blck++)
+ {
+ recursive_variables.clear();
+ if (ModelBlock->Block_List[blck].Simulation_Type==SOLVE_TWO_BOUNDARIES_COMPLETE or ModelBlock->Block_List[blck].Simulation_Type==SOLVE_TWO_BOUNDARIES_SIMPLE)
+ {
+ ModelBlock->Block_List[blck].Chain_Rule_Derivatives->clear();
+ for(int i = 0; i < ModelBlock->Block_List[blck].Nb_Recursives; i++)
+ {
+ if (ModelBlock->Block_List[blck].Equation_Type[i] == E_EVALUATE_S)
+ recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = ModelBlock->Block_List[blck].Equation_Normalized[i];
+ else
+ recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = equations[ModelBlock->Block_List[blck].Equation[i]];
+ }
+ map<pair<pair<int, pair<int, int> >, pair<int, int> >, int> Derivatives = block_triangular.get_Derivatives(ModelBlock, blck);
+
+ map<pair<pair<int, pair<int, int> >, pair<int, int> >, int>::const_iterator it = Derivatives.begin();
+ //#pragma omp parallel for shared(it, blck)
+ for(int i=0; i<(int)Derivatives.size(); i++)
+ {
+ int Deriv_type = it->second;
+ pair<pair<int, pair<int, int> >, pair<int, int> > it_l(it->first);
+ it++;
+ int lag = it_l.first.first;
+ int eq = it_l.first.second.first;
+ int var = it_l.first.second.second;
+ int eqr = it_l.second.first;
+ int varr = it_l.second.second;
+ if(Deriv_type == 0)
+ {
+ first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = first_derivatives[make_pair(eqr, getDerivID(symbol_table.getID(eEndogenous, varr), lag))];
+ }
+ else if (Deriv_type == 1)
+ {
+ first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = ModelBlock->Block_List[blck].Equation_Normalized[eq]->getChainRuleDerivative(getDerivID(symbol_table.getID(eEndogenous, varr), lag), recursive_variables);
+ }
+ else if (Deriv_type == 2)
+ {
+ if(ModelBlock->Block_List[blck].Equation_Type[eq] == E_EVALUATE_S && eq<ModelBlock->Block_List[blck].Nb_Recursives)
+ first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = ModelBlock->Block_List[blck].Equation_Normalized[eq]->getChainRuleDerivative(getDerivID(symbol_table.getID(eEndogenous, varr), lag), recursive_variables);
+ else
+ first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = equations[eqr]->getChainRuleDerivative(getDerivID(symbol_table.getID(eEndogenous, varr), lag), recursive_variables);
+ }
+ ModelBlock->Block_List[blck].Chain_Rule_Derivatives->push_back(make_pair( make_pair(lag, make_pair(eq, var)), make_pair(eqr, varr)));
+ }
+ }
+ else if( ModelBlock->Block_List[blck].Simulation_Type==SOLVE_BACKWARD_SIMPLE or ModelBlock->Block_List[blck].Simulation_Type==SOLVE_FORWARD_SIMPLE
+ or ModelBlock->Block_List[blck].Simulation_Type==SOLVE_BACKWARD_COMPLETE or ModelBlock->Block_List[blck].Simulation_Type==SOLVE_FORWARD_COMPLETE)
+ {
+ ModelBlock->Block_List[blck].Chain_Rule_Derivatives->clear();
+ for(int i = 0; i < ModelBlock->Block_List[blck].Nb_Recursives; i++)
+ {
+ if (ModelBlock->Block_List[blck].Equation_Type[i] == E_EVALUATE_S)
+ recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = ModelBlock->Block_List[blck].Equation_Normalized[i];
+ else
+ recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = equations[ModelBlock->Block_List[blck].Equation[i]];
+ }
+ for(int eq = ModelBlock->Block_List[blck].Nb_Recursives; eq < ModelBlock->Block_List[blck].Size; eq++)
+ {
+ int eqr = ModelBlock->Block_List[blck].Equation[eq];
+ for(int var = ModelBlock->Block_List[blck].Nb_Recursives; var < ModelBlock->Block_List[blck].Size; var++)
+ {
+ int varr = ModelBlock->Block_List[blck].Variable[var];
+ NodeID d1 = equations[eqr]->getChainRuleDerivative(getDerivID(symbol_table.getID(eEndogenous, varr), 0), recursive_variables);
+ if (d1 == Zero)
+ continue;
+ first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, 0))] = d1;
+ ModelBlock->Block_List[blck].Chain_Rule_Derivatives->push_back(make_pair( make_pair(0, make_pair(eq, var)), make_pair(eqr, varr)));
+ }
+ }
+ }
+ }
+}
+
+
+
+void
+StaticDllModel::computeParamsDerivatives()
+{
+ for (deriv_id_table_t::const_iterator it = deriv_id_table.begin();
+ it != deriv_id_table.end(); it++)
+ {
+ if (symbol_table.getType(it->first.first) != eParameter)
+ continue;
+
+ int param = it->second;
+
+ for (first_derivatives_type::const_iterator it2 = first_derivatives.begin();
+ it2 != first_derivatives.end(); it2++)
+ {
+ int eq = it2->first.first;
+ int var = it2->first.second;
+ NodeID d1 = it2->second;
+
+ NodeID d2 = d1->getDerivative(param);
+ if (d2 == Zero)
+ continue;
+ params_derivatives[make_pair(eq, make_pair(var, param))] = d2;
+ }
+ }
+}
+
+void
+StaticDllModel::computeParamsDerivativesTemporaryTerms()
+{
+ map<NodeID, int> reference_count;
+ params_derivs_temporary_terms.clear();
+
+ for (second_derivatives_type::iterator it = params_derivatives.begin();
+ it != params_derivatives.end(); it++)
+ it->second->computeTemporaryTerms(reference_count, params_derivs_temporary_terms, true);
+}
+
+void
+StaticDllModel::writeParamsDerivativesFile(const string &basename) const
+ {
+ if (!params_derivatives.size())
+ return;
+
+ string filename = basename + "_params_derivs.m";
+
+ ofstream paramsDerivsFile;
+ paramsDerivsFile.open(filename.c_str(), ios::out | ios::binary);
+ if (!paramsDerivsFile.is_open())
+ {
+ cerr << "ERROR: Can't open file " << filename << " for writing" << endl;
+ exit(EXIT_FAILURE);
+ }
+ paramsDerivsFile << "function gp = " << basename << "_params_derivs(y, x, params, it_)" << endl
+ << "%" << endl
+ << "% Warning : this file is generated automatically by Dynare" << endl
+ << "% from model file (.mod)" << endl << endl;
+
+
+ writeTemporaryTerms(params_derivs_temporary_terms, paramsDerivsFile, oMatlabStaticModel);
+
+ paramsDerivsFile << "gp = zeros(" << equation_number() << ", " << dynJacobianColsNbr << ", "
+ << symbol_table.param_nbr() << ");" << endl;
+
+ for (second_derivatives_type::const_iterator it = params_derivatives.begin();
+ it != params_derivatives.end(); it++)
+ {
+ int eq = it->first.first;
+ int var = it->first.second.first;
+ int param = it->first.second.second;
+ NodeID d2 = it->second;
+
+ int var_col = getDynJacobianCol(var) + 1;
+ int param_col = symbol_table.getTypeSpecificID(getSymbIDByDerivID(param)) + 1;
+
+ paramsDerivsFile << "gp(" << eq+1 << ", " << var_col << ", " << param_col << ") = ";
+ d2->writeOutput(paramsDerivsFile, oMatlabStaticModel, params_derivs_temporary_terms);
+ paramsDerivsFile << ";" << endl;
+ }
+
+ paramsDerivsFile.close();
+ }
+
+
+
+void
+StaticDllModel::writeChainRuleDerivative(ostream &output, int eqr, int varr, int lag,
+ ExprNodeOutputType output_type,
+ const temporary_terms_type &temporary_terms) const
+{
+ map<pair<int, pair<int, int> >, NodeID>::const_iterator it = first_chain_rule_derivatives.find(make_pair(eqr, make_pair(varr, lag)));
+ if (it != first_chain_rule_derivatives.end())
+ (it->second)->writeOutput(output, output_type, temporary_terms);
+ else
+ output << 0;
+}
+
+
+void
+StaticDllModel::writeLatexFile(const string &basename) const
+ {
+ writeLatexModelFile(basename + "_static.tex", oLatexStaticModel);
+ }
+
+void
+StaticDllModel::jacobianHelper(ostream &output, int eq_nb, int col_nb, ExprNodeOutputType output_type) const
+{
+ output << LEFT_ARRAY_SUBSCRIPT(output_type);
+ if (IS_MATLAB(output_type))
+ output << eq_nb + 1 << ", " << col_nb + 1;
+ else
+ output << eq_nb + col_nb * equations.size();
+ output << RIGHT_ARRAY_SUBSCRIPT(output_type);
+}
+
+void
+StaticDllModel::hessianHelper(ostream &output, int row_nb, int col_nb, ExprNodeOutputType output_type) const
+{
+ output << LEFT_ARRAY_SUBSCRIPT(output_type);
+ if (IS_MATLAB(output_type))
+ output << row_nb + 1 << ", " << col_nb + 1;
+ else
+ output << row_nb + col_nb * NNZDerivatives[1];
+ output << RIGHT_ARRAY_SUBSCRIPT(output_type);
+}
+
+
diff --git a/preprocessor/StaticDllModel.hh b/preprocessor/StaticDllModel.hh
new file mode 100644
index 0000000000000000000000000000000000000000..dd067710ee3ef504affc504467e174a5547e8484
--- /dev/null
+++ b/preprocessor/StaticDllModel.hh
@@ -0,0 +1,192 @@
+/*
+ * Copyright (C) 2003-2009 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/>.
+ */
+
+#ifndef _STATICDLLMODEL_HH
+#define _STATICDLLMODEL_HH
+
+using namespace std;
+
+#include <fstream>
+
+#include "StaticModel.hh"
+#include "BlockTriangular.hh"
+
+//! Stores a static model
+class StaticDllModel : public ModelTree
+{
+public:
+ //! The modes in which StaticDllModel can work
+ enum mode_t
+ {
+ eStandardMode, //!< Standard mode (static file in Matlab)
+ eSparseMode, //!< Sparse mode (static file in Matlab with block decomposition)
+ eDLLMode, //!< DLL mode (static file in C)
+ eSparseDLLMode //!< Sparse DLL mode (static file in C with block decomposition plus a binary file)
+ };
+private:
+ typedef map<pair<int, int>, int> deriv_id_table_t;
+ //! Maps a pair (symbol_id, lag) to a deriv ID
+ deriv_id_table_t deriv_id_table;
+ //! Maps a deriv ID to a pair (symbol_id, lag)
+ vector<pair<int, int> > inv_deriv_id_table;
+
+ //! Maps a deriv_id to the column index of the static Jacobian
+ /*! Contains only endogenous, exogenous and exogenous deterministic */
+ map<int, int> dyn_jacobian_cols_table;
+
+ //! Maximum lag and lead over all types of variables (positive values)
+ /*! Set by computeDerivID() */
+ int max_lag, max_lead;
+ //! Maximum lag and lead over endogenous variables (positive values)
+ /*! Set by computeDerivID() */
+ int max_endo_lag, max_endo_lead;
+ //! Maximum lag and lead over exogenous variables (positive values)
+ /*! Set by computeDerivID() */
+ int max_exo_lag, max_exo_lead;
+ //! Maximum lag and lead over deterministic exogenous variables (positive values)
+ /*! Set by computeDerivID() */
+ int max_exo_det_lag, max_exo_det_lead;
+
+ //! Number of columns of static jacobian
+ /*! Set by computeDerivID() and computeDynJacobianCols() */
+ int dynJacobianColsNbr;
+
+ //! Derivatives of the jacobian w.r. to parameters
+ /*! First index is equation number, second is endo/exo/exo_det variable, and third is parameter.
+ Only non-null derivatives are stored in the map.
+ Variable and parameter indices are those of the getDerivID() method.
+ */
+ second_derivatives_type params_derivatives;
+
+ //! Temporary terms for the file containing parameters dervicatives
+ temporary_terms_type params_derivs_temporary_terms;
+
+ typedef map< pair< int, pair< int, int> >, NodeID> first_chain_rule_derivatives_type;
+ first_chain_rule_derivatives_type first_chain_rule_derivatives;
+
+
+ //! Writes static model file (Matlab version)
+ void writeStaticMFile(const string &static_basename) const;
+ //! Writes static model file (C version)
+ /*! \todo add third derivatives handling */
+ void writeStaticCFile(const string &static_basename) const;
+ //! Writes static model file when SparseDLL option is on
+ void writeSparseStaticMFile(const string &static_basename, const string &basename, const int mode) const;
+ //! Writes the Block reordred structure of the model in M output
+ void writeModelEquationsOrdered_M(Model_Block *ModelBlock, const string &static_basename) const;
+ //! Writes the code of the Block reordred structure of the model in virtual machine bytecode
+ void writeModelEquationsCodeOrdered(const string file_name, const Model_Block *ModelBlock, const string bin_basename, map_idx_type map_idx) const;
+ //! Computes jacobian and prepares for equation normalization
+ /*! Using values from initval/endval blocks and parameter initializations:
+ - computes the jacobian for the model w.r. to contemporaneous variables
+ - removes edges of the incidence matrix when derivative w.r. to the corresponding variable is too close to zero (below the cutoff)
+ */
+ void evaluateJacobian(const eval_context_type &eval_context, jacob_map *j_m, bool dynamic);
+ void BlockLinear(Model_Block *ModelBlock);
+ map_idx_type map_idx;
+ //! Build The incidence matrix form the modeltree
+ void BuildIncidenceMatrix();
+
+ void computeTemporaryTermsOrdered(Model_Block *ModelBlock);
+ //! Write derivative code of an equation w.r. to a variable
+ void compileDerivative(ofstream &code_file, int eq, int symb_id, int lag, map_idx_type &map_idx) const;
+ //! Write chain rule derivative code of an equation w.r. to a variable
+ void compileChainRuleDerivative(ofstream &code_file, int eq, int var, int lag, map_idx_type &map_idx) const;
+
+ virtual int computeDerivID(int symb_id, int lag);
+ //! Get the type corresponding to a derivation ID
+ SymbolType getTypeByDerivID(int deriv_id) const throw (UnknownDerivIDException);
+ //! Get the lag corresponding to a derivation ID
+ int getLagByDerivID(int deriv_id) const throw (UnknownDerivIDException);
+ //! Get the symbol ID corresponding to a derivation ID
+ int getSymbIDByDerivID(int deriv_id) const throw (UnknownDerivIDException);
+ //! Compute the column indices of the static Jacobian
+ void computeStatJacobianCols();
+ //! Computes chain rule derivatives of the Jacobian w.r. to endogenous variables
+ void computeChainRuleJacobian(Model_Block *ModelBlock);
+ //! Computes derivatives of the Jacobian w.r. to parameters
+ void computeParamsDerivatives();
+ //! Computes temporary terms for the file containing parameters derivatives
+ void computeParamsDerivativesTemporaryTerms();
+ //! Collect only the first derivatives
+ map<pair<int, pair<int, int> >, NodeID> collect_first_order_derivatives_endogenous();
+
+ //! Helper for writing the Jacobian elements in MATLAB and C
+ /*! Writes either (i+1,j+1) or [i+j*no_eq] */
+ void jacobianHelper(ostream &output, int eq_nb, int col_nb, ExprNodeOutputType output_type) const;
+
+ //! Helper for writing the sparse Hessian elements in MATLAB and C
+ /*! Writes either (i+1,j+1) or [i+j*NNZDerivatives[1]] */
+ void hessianHelper(ostream &output, int row_nb, int col_nb, ExprNodeOutputType output_type) const;
+
+ //! Write chain rule derivative of a recursive equation w.r. to a variable
+ void writeChainRuleDerivative(ostream &output, int eq, int var, int lag, ExprNodeOutputType output_type, const temporary_terms_type &temporary_terms) const;
+
+
+public:
+ StaticDllModel(SymbolTable &symbol_table_arg, NumericalConstants &num_constants);
+ //! Mode in which the ModelTree is supposed to work (Matlab, DLL or SparseDLL)
+ mode_t mode;
+ //! Adds a variable node
+ /*! This implementation allows for non-zero lag */
+ virtual NodeID AddVariable(const string &name, int lag = 0);
+ //! Absolute value under which a number is considered to be zero
+ double cutoff;
+ //! The weight of the Markowitz criteria to determine the pivot in the linear solver (simul_NG1 and simul_NG from simulate.cc)
+ double markowitz;
+ //! Compute the minimum feedback set in the static model:
+ /*! 0 : all endogenous variables are considered as feedback variables
+ 1 : the variables belonging to a non linear equation are considered as feedback variables
+ 2 : the variables belonging to a non normalizable non linear equation are considered as feedback variables
+ default value = 0 */
+ int mfs;
+ //! the file containing the model and the derivatives code
+ ofstream code_file;
+ //! Execute computations (variable sorting + derivation)
+ /*!
+ \param jacobianExo whether derivatives w.r. to exo and exo_det should be in the Jacobian (derivatives w.r. to endo are always computed)
+ \param hessian whether 2nd derivatives w.r. to exo, exo_det and endo should be computed (implies jacobianExo = true)
+ \param thirdDerivatives whether 3rd derivatives w.r. to endo/exo/exo_det should be computed (implies jacobianExo = true)
+ \param paramsDerivatives whether 2nd derivatives w.r. to a pair (endo/exo/exo_det, parameter) should be computed (implies jacobianExo = true)
+ \param eval_context evaluation context for normalization
+ \param no_tmp_terms if true, no temporary terms will be computed in the static files
+ */
+ void computingPass(const eval_context_type &eval_context, bool no_tmp_terms);
+ //! Writes model initialization and lead/lag incidence matrix to output
+ void writeOutput(ostream &output, const string &basename) const;
+ //! Write statements to be added to the main M-file, after computational tasks
+ void writeOutputPostComputing(ostream &output, const string &basename) const;
+ //! Complete set to block decompose the model
+ BlockTriangular block_triangular;
+ //! Adds informations for simulation in a binary file
+ void Write_Inf_To_Bin_File(const string &static_basename, const string &bin_basename,
+ const int &num, int &u_count_int, bool &file_open) const;
+ //! Writes static model file
+ void writeStaticFile(const string &basename) const;
+ //! Writes file containing parameters derivatives
+ void writeParamsDerivativesFile(const string &basename) const;
+
+ //! Writes LaTeX file with the equations of the static model
+ void writeLatexFile(const string &basename) const;
+
+ virtual int getDerivID(int symb_id, int lag) const throw (UnknownDerivIDException);
+ virtual int getDynJacobianCol(int deriv_id) const throw (UnknownDerivIDException);
+};
+
+#endif