ExprNode.cc 185 KB
Newer Older
1
/*
2
 * Copyright (C) 2007-2013 Dynare Team
sebastien's avatar
sebastien committed
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
 *
 * 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/>.
 */

sebastien's avatar
sebastien committed
20
21
22
23
#include <iostream>
#include <iterator>
#include <algorithm>

24
25
26
27
28
// For select1st()
#ifdef __GNUC__
# include <ext/functional>
using namespace __gnu_cxx;
#endif
sebastien's avatar
sebastien committed
29

30
#include <cassert>
31
#include <cmath>
sebastien's avatar
sebastien committed
32

sebastien's avatar
sebastien committed
33
34
#include "ExprNode.hh"
#include "DataTree.hh"
35
#include "ModFile.hh"
sebastien's avatar
sebastien committed
36

sebastien's avatar
sebastien committed
37
ExprNode::ExprNode(DataTree &datatree_arg) : datatree(datatree_arg), preparedForDerivation(false)
sebastien's avatar
sebastien committed
38
39
40
41
42
43
44
45
46
47
48
49
{
  // Add myself to datatree
  datatree.node_list.push_back(this);

  // Set my index and increment counter
  idx = datatree.node_counter++;
}

ExprNode::~ExprNode()
{
}

50
expr_t
51
ExprNode::getDerivative(int deriv_id)
sebastien's avatar
sebastien committed
52
{
sebastien's avatar
sebastien committed
53
54
  if (!preparedForDerivation)
    prepareForDerivation();
55

sebastien's avatar
sebastien committed
56
  // Return zero if derivative is necessarily null (using symbolic a priori)
57
  set<int>::const_iterator it = non_null_derivatives.find(deriv_id);
sebastien's avatar
sebastien committed
58
59
  if (it == non_null_derivatives.end())
    return datatree.Zero;
60

sebastien's avatar
sebastien committed
61
  // If derivative is stored in cache, use the cached value, otherwise compute it (and cache it)
62
  map<int, expr_t>::const_iterator it2 = derivatives.find(deriv_id);
sebastien's avatar
sebastien committed
63
64
65
66
  if (it2 != derivatives.end())
    return it2->second;
  else
    {
67
      expr_t d = computeDerivative(deriv_id);
68
      derivatives[deriv_id] = d;
sebastien's avatar
sebastien committed
69
70
71
72
      return d;
    }
}

ferhat's avatar
ferhat committed
73
int
74
ExprNode::precedence(ExprNodeOutputType output_type, const temporary_terms_t &temporary_terms) const
75
76
77
78
{
  // For a constant, a variable, or a unary op, the precedence is maximal
  return 100;
}
ferhat's avatar
ferhat committed
79

sebastien's avatar
sebastien committed
80
int
81
ExprNode::cost(const temporary_terms_t &temporary_terms, bool is_matlab) const
82
83
84
85
{
  // For a terminal node, the cost is null
  return 0;
}
sebastien's avatar
sebastien committed
86

sebastien's avatar
sebastien committed
87
88
89
90
91
void
ExprNode::collectEndogenous(set<pair<int, int> > &result) const
{
  set<pair<int, int> > symb_ids;
  collectVariables(eEndogenous, symb_ids);
92
93
  for (set<pair<int, int> >::const_iterator it = symb_ids.begin();
       it != symb_ids.end(); it++)
sebastien's avatar
sebastien committed
94
95
96
97
98
99
100
101
    result.insert(make_pair(datatree.symbol_table.getTypeSpecificID(it->first), it->second));
}

void
ExprNode::collectExogenous(set<pair<int, int> > &result) const
{
  set<pair<int, int> > symb_ids;
  collectVariables(eExogenous, symb_ids);
102
103
  for (set<pair<int, int> >::const_iterator it = symb_ids.begin();
       it != symb_ids.end(); it++)
sebastien's avatar
sebastien committed
104
105
106
107
108
109
110
111
112
    result.insert(make_pair(datatree.symbol_table.getTypeSpecificID(it->first), it->second));
}

void
ExprNode::collectModelLocalVariables(set<int> &result) const
{
  set<pair<int, int> > symb_ids;
  collectVariables(eModelLocalVariable, symb_ids);
  transform(symb_ids.begin(), symb_ids.end(), inserter(result, result.begin()),
113
            select1st<pair<int, int> >());
sebastien's avatar
sebastien committed
114
115
}

sebastien's avatar
sebastien committed
116
void
117
ExprNode::computeTemporaryTerms(map<expr_t, int> &reference_count,
118
                                temporary_terms_t &temporary_terms,
sebastien's avatar
sebastien committed
119
                                bool is_matlab) const
120
121
122
{
  // Nothing to do for a terminal node
}
sebastien's avatar
sebastien committed
123

sebastien's avatar
sebastien committed
124
void
125
ExprNode::computeTemporaryTerms(map<expr_t, int> &reference_count,
126
                                temporary_terms_t &temporary_terms,
127
                                map<expr_t, pair<int, int> > &first_occurence,
sebastien's avatar
sebastien committed
128
                                int Curr_block,
129
                                vector<vector<temporary_terms_t> > &v_temporary_terms,
130
                                int equation) const
131
132
133
{
  // Nothing to do for a terminal node
}
sebastien's avatar
sebastien committed
134

135
136
pair<int, expr_t >
ExprNode::normalizeEquation(int var_endo, vector<pair<int, pair<expr_t, expr_t> > > &List_of_Op_RHS) const
137
{
138
  /* nothing to do */
139
  return (make_pair(0, (expr_t) NULL));
140
}
141

sebastien's avatar
sebastien committed
142
void
143
ExprNode::writeOutput(ostream &output) const
sebastien's avatar
sebastien committed
144
{
145
  writeOutput(output, oMatlabOutsideModel, temporary_terms_t());
sebastien's avatar
sebastien committed
146
147
}

148
149
150
void
ExprNode::writeOutput(ostream &output, ExprNodeOutputType output_type) const
{
151
  writeOutput(output, output_type, temporary_terms_t());
152
153
}

154
void
155
ExprNode::writeOutput(ostream &output, ExprNodeOutputType output_type, const temporary_terms_t &temporary_terms) const
156
{
157
  deriv_node_temp_terms_t tef_terms;
158
159
160
  writeOutput(output, output_type, temporary_terms, tef_terms);
}

161
162
void
ExprNode::compile(ostream &CompileCode, unsigned int &instruction_number,
163
164
                  bool lhs_rhs, const temporary_terms_t &temporary_terms,
                  const map_idx_t &map_idx, bool dynamic, bool steady_dynamic) const
165
166
167
168
169
{
  deriv_node_temp_terms_t tef_terms;
  compile(CompileCode, instruction_number, lhs_rhs, temporary_terms, map_idx, dynamic, steady_dynamic, tef_terms);
}

170
171
void
ExprNode::writeExternalFunctionOutput(ostream &output, ExprNodeOutputType output_type,
172
173
                                      const temporary_terms_t &temporary_terms,
                                      deriv_node_temp_terms_t &tef_terms) const
174
175
176
177
{
  // Nothing to do
}

178
179
void
ExprNode::compileExternalFunctionOutput(ostream &CompileCode, unsigned int &instruction_number,
180
181
182
                                        bool lhs_rhs, const temporary_terms_t &temporary_terms,
                                        const map_idx_t &map_idx, bool dynamic, bool steady_dynamic,
                                        deriv_node_temp_terms_t &tef_terms) const
183
184
185
186
{
  // Nothing to do
}

sebastien's avatar
sebastien committed
187
VariableNode *
188
ExprNode::createEndoLeadAuxiliaryVarForMyself(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
sebastien's avatar
sebastien committed
189
190
191
192
193
194
195
196
{
  int n = maxEndoLead();
  assert(n >= 2);

  subst_table_t::const_iterator it = subst_table.find(this);
  if (it != subst_table.end())
    return const_cast<VariableNode *>(it->second);

197
  expr_t substexpr = decreaseLeadsLags(n-1);
sebastien's avatar
sebastien committed
198
199
200
201
  int lag = n-2;

  // Each iteration tries to create an auxvar such that auxvar(+1)=expr(-lag)
  // At the beginning (resp. end) of each iteration, substexpr is an expression (possibly an auxvar) equivalent to expr(-lag-1) (resp. expr(-lag))
202
  while (lag >= 0)
sebastien's avatar
sebastien committed
203
    {
204
      expr_t orig_expr = decreaseLeadsLags(lag);
sebastien's avatar
sebastien committed
205
206
207
      it = subst_table.find(orig_expr);
      if (it == subst_table.end())
        {
208
          int symb_id = datatree.symbol_table.addEndoLeadAuxiliaryVar(orig_expr->idx);
sebastien's avatar
sebastien committed
209
210
211
212
213
214
215
216
217
218
219
220
221
222
          neweqs.push_back(dynamic_cast<BinaryOpNode *>(datatree.AddEqual(datatree.AddVariable(symb_id, 0), substexpr)));
          substexpr = datatree.AddVariable(symb_id, +1);
          assert(dynamic_cast<VariableNode *>(substexpr) != NULL);
          subst_table[orig_expr] = dynamic_cast<VariableNode *>(substexpr);
        }
      else
        substexpr = const_cast<VariableNode *>(it->second);

      lag--;
    }

  return dynamic_cast<VariableNode *>(substexpr);
}

sebastien's avatar
sebastien committed
223
224
225
226
227
228
229
230
231
232
VariableNode *
ExprNode::createExoLeadAuxiliaryVarForMyself(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
  int n = maxExoLead();
  assert(n >= 1);

  subst_table_t::const_iterator it = subst_table.find(this);
  if (it != subst_table.end())
    return const_cast<VariableNode *>(it->second);

233
  expr_t substexpr = decreaseLeadsLags(n);
sebastien's avatar
sebastien committed
234
235
236
237
  int lag = n-1;

  // Each iteration tries to create an auxvar such that auxvar(+1)=expr(-lag)
  // At the beginning (resp. end) of each iteration, substexpr is an expression (possibly an auxvar) equivalent to expr(-lag-1) (resp. expr(-lag))
238
  while (lag >= 0)
sebastien's avatar
sebastien committed
239
    {
240
      expr_t orig_expr = decreaseLeadsLags(lag);
sebastien's avatar
sebastien committed
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
      it = subst_table.find(orig_expr);
      if (it == subst_table.end())
        {
          int symb_id = datatree.symbol_table.addExoLeadAuxiliaryVar(orig_expr->idx);
          neweqs.push_back(dynamic_cast<BinaryOpNode *>(datatree.AddEqual(datatree.AddVariable(symb_id, 0), substexpr)));
          substexpr = datatree.AddVariable(symb_id, +1);
          assert(dynamic_cast<VariableNode *>(substexpr) != NULL);
          subst_table[orig_expr] = dynamic_cast<VariableNode *>(substexpr);
        }
      else
        substexpr = const_cast<VariableNode *>(it->second);

      lag--;
    }

  return dynamic_cast<VariableNode *>(substexpr);
}

259
260
261
262
263
264
265
266
267
268
269
270
bool
ExprNode::isNumConstNodeEqualTo(double value) const
{
  return false;
}

bool
ExprNode::isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const
{
  return false;
}

sebastien's avatar
sebastien committed
271
NumConstNode::NumConstNode(DataTree &datatree_arg, int id_arg) :
272
273
  ExprNode(datatree_arg),
  id(id_arg)
sebastien's avatar
sebastien committed
274
275
276
{
  // Add myself to the num const map
  datatree.num_const_node_map[id] = this;
sebastien's avatar
sebastien committed
277
}
sebastien's avatar
sebastien committed
278

sebastien's avatar
sebastien committed
279
280
281
282
void
NumConstNode::prepareForDerivation()
{
  preparedForDerivation = true;
sebastien's avatar
sebastien committed
283
284
285
  // All derivatives are null, so non_null_derivatives is left empty
}

286
expr_t
287
NumConstNode::computeDerivative(int deriv_id)
sebastien's avatar
sebastien committed
288
289
290
291
{
  return datatree.Zero;
}

292
void
293
NumConstNode::collectTemporary_terms(const temporary_terms_t &temporary_terms, temporary_terms_inuse_t &temporary_terms_inuse, int Curr_Block) const
294
{
295
  temporary_terms_t::const_iterator it = temporary_terms.find(const_cast<NumConstNode *>(this));
296
297
298
  if (it != temporary_terms.end())
    temporary_terms_inuse.insert(idx);
}
299

sebastien's avatar
sebastien committed
300
void
sebastien's avatar
sebastien committed
301
NumConstNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
302
303
                          const temporary_terms_t &temporary_terms,
                          deriv_node_temp_terms_t &tef_terms) const
304
{
305
  temporary_terms_t::const_iterator it = temporary_terms.find(const_cast<NumConstNode *>(this));
306
307
308
  if (it != temporary_terms.end())
    if (output_type == oMatlabDynamicModelSparse)
      output << "T" << idx << "(it_)";
ferhat's avatar
ferhat committed
309
    else
310
311
312
313
      output << "T" << idx;
  else
    output << datatree.num_constants.get(id);
}
sebastien's avatar
sebastien committed
314

sebastien's avatar
sebastien committed
315
double
316
NumConstNode::eval(const eval_context_t &eval_context) const throw (EvalException, EvalExternalFunctionException)
sebastien's avatar
sebastien committed
317
{
318
  return (datatree.num_constants.getDouble(id));
sebastien's avatar
sebastien committed
319
320
}

321
void
322
323
324
325
NumConstNode::compile(ostream &CompileCode, unsigned int &instruction_number,
                      bool lhs_rhs, const temporary_terms_t &temporary_terms,
                      const map_idx_t &map_idx, bool dynamic, bool steady_dynamic,
                      deriv_node_temp_terms_t &tef_terms) const
326
327
{
  FLDC_ fldc(datatree.num_constants.getDouble(id));
328
  fldc.write(CompileCode, instruction_number);
329
}
330

sebastien's avatar
sebastien committed
331
void
sebastien's avatar
sebastien committed
332
333
334
NumConstNode::collectVariables(SymbolType type_arg, set<pair<int, int> > &result) const
{
}
ferhat's avatar
ferhat committed
335

336
337
pair<int, expr_t >
NumConstNode::normalizeEquation(int var_endo, vector<pair<int, pair<expr_t, expr_t> > > &List_of_Op_RHS) const
338
{
339
  /* return the numercial constant */
340
  return (make_pair(0, datatree.AddNonNegativeConstant(datatree.num_constants.get(id))));
341
}
ferhat's avatar
ferhat committed
342

343
344
expr_t
NumConstNode::getChainRuleDerivative(int deriv_id, const map<int, expr_t> &recursive_variables)
ferhat's avatar
ferhat committed
345
{
ferhat's avatar
ferhat committed
346
  return datatree.Zero;
ferhat's avatar
ferhat committed
347
348
}

349
expr_t
sebastien's avatar
sebastien committed
350
NumConstNode::toStatic(DataTree &static_datatree) const
351
{
352
  return static_datatree.AddNonNegativeConstant(datatree.num_constants.get(id));
353
}
sebastien's avatar
sebastien committed
354

355
356
357
expr_t
NumConstNode::cloneDynamic(DataTree &dynamic_datatree) const
{
358
  return dynamic_datatree.AddNonNegativeConstant(datatree.num_constants.get(id));
359
360
}

sebastien's avatar
sebastien committed
361
362
363
364
365
int
NumConstNode::maxEndoLead() const
{
  return 0;
}
ferhat's avatar
ferhat committed
366

sebastien's avatar
sebastien committed
367
368
369
370
371
372
int
NumConstNode::maxExoLead() const
{
  return 0;
}

373
374
375
376
377
378
379
380
381
382
383
384
int
NumConstNode::maxEndoLag() const
{
  return 0;
}

int
NumConstNode::maxExoLag() const
{
  return 0;
}

385
386
387
388
389
390
int
NumConstNode::maxLead() const
{
  return 0;
}

391
expr_t
sebastien's avatar
sebastien committed
392
393
394
395
396
NumConstNode::decreaseLeadsLags(int n) const
{
  return const_cast<NumConstNode *>(this);
}

397
expr_t
sebastien's avatar
sebastien committed
398
NumConstNode::decreaseLeadsLagsPredeterminedVariables() const
399
{
sebastien's avatar
sebastien committed
400
  return const_cast<NumConstNode *>(this);
401
402
}

403
expr_t
404
NumConstNode::substituteEndoLeadGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool deterministic_model) const
sebastien's avatar
sebastien committed
405
406
407
408
{
  return const_cast<NumConstNode *>(this);
}

409
expr_t
410
411
412
413
414
NumConstNode::substituteEndoLagGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
  return const_cast<NumConstNode *>(this);
}

415
expr_t
416
NumConstNode::substituteExoLead(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool deterministic_model) const
sebastien's avatar
sebastien committed
417
418
419
420
{
  return const_cast<NumConstNode *>(this);
}

421
expr_t
sebastien's avatar
sebastien committed
422
NumConstNode::substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
sebastien's avatar
sebastien committed
423
424
425
426
{
  return const_cast<NumConstNode *>(this);
}

427
expr_t
428
429
430
431
432
NumConstNode::substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const
{
  return const_cast<NumConstNode *>(this);
}

433
434
435
436
437
438
expr_t
NumConstNode::substituteLogPow(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs1, vector<BinaryOpNode *> &neweqs2) const
{
  return const_cast<NumConstNode *>(this);
}

439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
bool
NumConstNode::isNumConstNodeEqualTo(double value) const
{
  if (datatree.num_constants.getDouble(id) == value)
    return true;
  else
    return false;
}

bool
NumConstNode::isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const
{
  return false;
}

454
455
456
457
458
459
bool
NumConstNode::containsEndogenous(void) const
{
  return false;
}

460
461
462
463
464
465
466
expr_t
NumConstNode::replaceTrendVar() const
{
  return const_cast<NumConstNode *>(this);
}

expr_t
467
NumConstNode::detrend(int symb_id, bool log_trend, expr_t trend) const
468
469
470
471
472
473
474
475
476
477
{
  return const_cast<NumConstNode *>(this);
}

expr_t
NumConstNode::removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const
{
  return const_cast<NumConstNode *>(this);
}

Houtan Bastani's avatar
Houtan Bastani committed
478
479
480
481
482
483
484
485
486
487
488
489
490
491
VariableNode::VariableNode(DataTree &datatree_arg, int symb_id_arg, int lag_arg) :
  ExprNode(datatree_arg),
  symb_id(symb_id_arg),
  type(datatree.symbol_table.getType(symb_id_arg)),
  lag(lag_arg)
{
  // Add myself to the variable map
  datatree.variable_node_map[make_pair(symb_id, lag)] = this;

  // It makes sense to allow a lead/lag on parameters: during steady state calibration, endogenous and parameters can be swapped
  assert(type != eExternalFunction
         && (lag == 0 || (type != eModelLocalVariable && type != eModFileLocalVariable)));
}

sebastien's avatar
sebastien committed
492
493
494
495
496
void
VariableNode::prepareForDerivation()
{
  if (preparedForDerivation)
    return;
497

sebastien's avatar
sebastien committed
498
  preparedForDerivation = true;
sebastien's avatar
sebastien committed
499

sebastien's avatar
sebastien committed
500
  // Fill in non_null_derivatives
ferhat's avatar
ferhat committed
501
  switch (type)
sebastien's avatar
sebastien committed
502
503
504
505
506
    {
    case eEndogenous:
    case eExogenous:
    case eExogenousDet:
    case eParameter:
507
    case eTrend:
508
    case eLogTrend:
sebastien's avatar
sebastien committed
509
      // For a variable or a parameter, the only non-null derivative is with respect to itself
510
      non_null_derivatives.insert(datatree.getDerivID(symb_id, lag));
sebastien's avatar
sebastien committed
511
      break;
sebastien's avatar
sebastien committed
512
    case eModelLocalVariable:
sebastien's avatar
sebastien committed
513
      datatree.local_variables_table[symb_id]->prepareForDerivation();
sebastien's avatar
sebastien committed
514
      // Non null derivatives are those of the value of the local parameter
sebastien's avatar
sebastien committed
515
516
517
      non_null_derivatives = datatree.local_variables_table[symb_id]->non_null_derivatives;
      break;
    case eModFileLocalVariable:
518
    case eStatementDeclaredVariable:
sebastien's avatar
sebastien committed
519
      // Such a variable is never derived
sebastien's avatar
sebastien committed
520
      break;
521
    case eExternalFunction:
sebastien's avatar
sebastien committed
522
      cerr << "VariableNode::prepareForDerivation: impossible case" << endl;
523
      exit(EXIT_FAILURE);
sebastien's avatar
sebastien committed
524
525
526
    }
}

527
expr_t
sebastien's avatar
sebastien committed
528
VariableNode::computeDerivative(int deriv_id)
sebastien's avatar
sebastien committed
529
{
ferhat's avatar
ferhat committed
530
  switch (type)
sebastien's avatar
sebastien committed
531
532
533
534
    {
    case eEndogenous:
    case eExogenous:
    case eExogenousDet:
sebastien's avatar
sebastien committed
535
    case eParameter:
536
    case eTrend:
537
    case eLogTrend:
sebastien's avatar
sebastien committed
538
      if (deriv_id == datatree.getDerivID(symb_id, lag))
sebastien's avatar
sebastien committed
539
540
541
        return datatree.One;
      else
        return datatree.Zero;
sebastien's avatar
sebastien committed
542
    case eModelLocalVariable:
sebastien's avatar
sebastien committed
543
      return datatree.local_variables_table[symb_id]->getDerivative(deriv_id);
sebastien's avatar
sebastien committed
544
545
    case eModFileLocalVariable:
      cerr << "ModFileLocalVariable is not derivable" << endl;
546
      exit(EXIT_FAILURE);
547
548
549
    case eStatementDeclaredVariable:
      cerr << "eStatementDeclaredVariable is not derivable" << endl;
      exit(EXIT_FAILURE);
550
    case eExternalFunction:
551
      cerr << "Impossible case!" << endl;
552
      exit(EXIT_FAILURE);
sebastien's avatar
sebastien committed
553
    }
sebastien's avatar
sebastien committed
554
  // Suppress GCC warning
555
  exit(EXIT_FAILURE);
sebastien's avatar
sebastien committed
556
557
}

558
void
559
VariableNode::collectTemporary_terms(const temporary_terms_t &temporary_terms, temporary_terms_inuse_t &temporary_terms_inuse, int Curr_Block) const
560
{
561
  temporary_terms_t::const_iterator it = temporary_terms.find(const_cast<VariableNode *>(this));
562
563
564
565
566
  if (it != temporary_terms.end())
    temporary_terms_inuse.insert(idx);
  if (type == eModelLocalVariable)
    datatree.local_variables_table[symb_id]->collectTemporary_terms(temporary_terms, temporary_terms_inuse, Curr_Block);
}
567

sebastien's avatar
sebastien committed
568
void
sebastien's avatar
sebastien committed
569
VariableNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
570
571
                          const temporary_terms_t &temporary_terms,
                          deriv_node_temp_terms_t &tef_terms) const
572
573
{
  // If node is a temporary term
574
  temporary_terms_t::const_iterator it = temporary_terms.find(const_cast<VariableNode *>(this));
575
576
577
578
579
580
581
582
  if (it != temporary_terms.end())
    {
      if (output_type == oMatlabDynamicModelSparse)
        output << "T" << idx << "(it_)";
      else
        output << "T" << idx;
      return;
    }
sebastien's avatar
sebastien committed
583

584
585
586
587
588
589
  if (IS_LATEX(output_type))
    {
      if (output_type == oLatexDynamicSteadyStateOperator)
        output << "\\bar{";
      output << datatree.symbol_table.getTeXName(symb_id);
      if (output_type == oLatexDynamicModel
590
          && (type == eEndogenous || type == eExogenous || type == eExogenousDet || type == eModelLocalVariable || type == eTrend || type == eLogTrend))
591
592
593
594
595
596
597
598
        {
          output << "_{t";
          if (lag != 0)
            {
              if (lag > 0)
                output << "+";
              output << lag;
            }
599
          output << "}";
600
601
602
603
604
        }
      else if (output_type == oLatexDynamicSteadyStateOperator)
        output << "}";
      return;
    }
sebastien's avatar
sebastien committed
605

606
607
608
609
610
  int i;
  int tsid = datatree.symbol_table.getTypeSpecificID(symb_id);
  switch (type)
    {
    case eParameter:
611
      if (output_type == oMatlabOutsideModel)
612
613
614
615
        output << "M_.params" << "(" << tsid + 1 << ")";
      else
        output << "params" << LEFT_ARRAY_SUBSCRIPT(output_type) << tsid + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type);
      break;
ferhat's avatar
ferhat committed
616

617
    case eModelLocalVariable:
Sébastien Villemot's avatar
Sébastien Villemot committed
618
      if (output_type == oMatlabDynamicModelSparse || output_type == oMatlabStaticModelSparse)
619
620
        {
          output << "(";
621
          datatree.local_variables_table[symb_id]->writeOutput(output, output_type, temporary_terms, tef_terms);
622
623
624
          output << ")";
        }
      else
625
626
627
        /* We append underscores to avoid name clashes with "g1" or "oo_" (see
           also ModelTree::writeModelLocalVariables) */
        output << datatree.symbol_table.getName(symb_id) << "__";
628
      break;
ferhat's avatar
ferhat committed
629

630
631
632
633
    case eModFileLocalVariable:
      output << datatree.symbol_table.getName(symb_id);
      break;

634
635
636
637
638
639
640
641
    case eEndogenous:
      switch (output_type)
        {
        case oMatlabDynamicModel:
        case oCDynamicModel:
          i = datatree.getDynJacobianCol(datatree.getDerivID(symb_id, lag)) + ARRAY_SUBSCRIPT_OFFSET(output_type);
          output <<  "y" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << RIGHT_ARRAY_SUBSCRIPT(output_type);
          break;
642
        case oCStaticModel:
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
        case oMatlabStaticModel:
        case oMatlabStaticModelSparse:
          i = tsid + ARRAY_SUBSCRIPT_OFFSET(output_type);
          output <<  "y" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << RIGHT_ARRAY_SUBSCRIPT(output_type);
          break;
        case oMatlabDynamicModelSparse:
          i = tsid + ARRAY_SUBSCRIPT_OFFSET(output_type);
          if (lag > 0)
            output << "y" << LEFT_ARRAY_SUBSCRIPT(output_type) << "it_+" << lag << ", " << i << RIGHT_ARRAY_SUBSCRIPT(output_type);
          else if (lag < 0)
            output << "y" << LEFT_ARRAY_SUBSCRIPT(output_type) << "it_" << lag << ", " << i << RIGHT_ARRAY_SUBSCRIPT(output_type);
          else
            output << "y" << LEFT_ARRAY_SUBSCRIPT(output_type) << "it_, " << i << RIGHT_ARRAY_SUBSCRIPT(output_type);
          break;
        case oMatlabOutsideModel:
          output << "oo_.steady_state(" << tsid + 1 << ")";
          break;
        case oMatlabDynamicSteadyStateOperator:
661
        case oMatlabDynamicSparseSteadyStateOperator:
662
          output << "steady_state(" << tsid + 1 << ")";
663
          break;
664
665
666
        case oCDynamicSteadyStateOperator:
          output << "steady_state[" << tsid << "]";
          break;
667
668
669
        case oSteadyStateFile:
          output << "ys_(" << tsid + 1 << ")";
          break;
670
        default:
671
672
          cerr << "VariableNode::writeOutput: should not reach this point" << endl;
          exit(EXIT_FAILURE);
673
674
        }
      break;
ferhat's avatar
ferhat committed
675

676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
    case eExogenous:
      i = tsid + ARRAY_SUBSCRIPT_OFFSET(output_type);
      switch (output_type)
        {
        case oMatlabDynamicModel:
        case oMatlabDynamicModelSparse:
          if (lag > 0)
            output <<  "x(it_+" << lag << ", " << i << ")";
          else if (lag < 0)
            output <<  "x(it_" << lag << ", " << i << ")";
          else
            output <<  "x(it_, " << i << ")";
          break;
        case oCDynamicModel:
          if (lag == 0)
            output <<  "x[it_+" << i << "*nb_row_x]";
          else if (lag > 0)
            output <<  "x[it_+" << lag << "+" << i << "*nb_row_x]";
          else
            output <<  "x[it_" << lag << "+" << i << "*nb_row_x]";
          break;
697
        case oCStaticModel:
698
699
700
701
702
703
704
705
706
707
708
        case oMatlabStaticModel:
        case oMatlabStaticModelSparse:
          output << "x" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << RIGHT_ARRAY_SUBSCRIPT(output_type);
          break;
        case oMatlabOutsideModel:
          assert(lag == 0);
          output <<  "oo_.exo_steady_state(" << i << ")";
          break;
        case oMatlabDynamicSteadyStateOperator:
          output <<  "oo_.exo_steady_state(" << i << ")";
          break;
709
710
711
        case oSteadyStateFile:
          output << "exo_(" << i << ")";
          break;
712
        default:
713
714
          cerr << "VariableNode::writeOutput: should not reach this point" << endl;
          exit(EXIT_FAILURE);
715
716
        }
      break;
ferhat's avatar
ferhat committed
717

718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
    case eExogenousDet:
      i = tsid + datatree.symbol_table.exo_nbr() + ARRAY_SUBSCRIPT_OFFSET(output_type);
      switch (output_type)
        {
        case oMatlabDynamicModel:
        case oMatlabDynamicModelSparse:
          if (lag > 0)
            output <<  "x(it_+" << lag << ", " << i << ")";
          else if (lag < 0)
            output <<  "x(it_" << lag << ", " << i << ")";
          else
            output <<  "x(it_, " << i << ")";
          break;
        case oCDynamicModel:
          if (lag == 0)
733
            output <<  "x[it_+" << i << "*nb_row_x]";
734
          else if (lag > 0)
735
            output <<  "x[it_+" << lag << "+" << i << "*nb_row_x]";
736
          else
737
            output <<  "x[it_" << lag << "+" << i << "*nb_row_x]";
738
          break;
739
        case oCStaticModel:
740
741
742
743
744
745
746
747
748
749
750
        case oMatlabStaticModel:
        case oMatlabStaticModelSparse:
          output << "x" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << RIGHT_ARRAY_SUBSCRIPT(output_type);
          break;
        case oMatlabOutsideModel:
          assert(lag == 0);
          output <<  "oo_.exo_det_steady_state(" << tsid + 1 << ")";
          break;
        case oMatlabDynamicSteadyStateOperator:
          output <<  "oo_.exo_det_steady_state(" << tsid + 1 << ")";
          break;
751
752
753
        case oSteadyStateFile:
          output << "exo_(" << i << ")";
          break;
754
        default:
755
756
          cerr << "VariableNode::writeOutput: should not reach this point" << endl;
          exit(EXIT_FAILURE);
757
758
        }
      break;
ferhat's avatar
ferhat committed
759

760
    case eExternalFunction:
761
    case eTrend:
762
    case eLogTrend:
763
    case eStatementDeclaredVariable:
764
765
766
767
      cerr << "Impossible case" << endl;
      exit(EXIT_FAILURE);
    }
}
sebastien's avatar
sebastien committed
768

sebastien's avatar
sebastien committed
769
double
770
VariableNode::eval(const eval_context_t &eval_context) const throw (EvalException, EvalExternalFunctionException)
sebastien's avatar
sebastien committed
771
{
772
  eval_context_t::const_iterator it = eval_context.find(symb_id);
sebastien's avatar
sebastien committed
773
  if (it == eval_context.end())
sebastien's avatar
sebastien committed
774
    throw EvalException();
sebastien's avatar
sebastien committed
775

sebastien's avatar
sebastien committed
776
  return it->second;
sebastien's avatar
sebastien committed
777
778
}

779
void
780
781
782
783
VariableNode::compile(ostream &CompileCode, unsigned int &instruction_number,
                      bool lhs_rhs, const temporary_terms_t &temporary_terms,
                      const map_idx_t &map_idx, bool dynamic, bool steady_dynamic,
                      deriv_node_temp_terms_t &tef_terms) const
784
785
{
  if (type == eModelLocalVariable || type == eModFileLocalVariable)
786
    datatree.local_variables_table[symb_id]->compile(CompileCode, instruction_number, lhs_rhs, temporary_terms, map_idx, dynamic, steady_dynamic, tef_terms);
787
788
789
790
791
792
793
794
795
796
797
798
  else
    {
      int tsid = datatree.symbol_table.getTypeSpecificID(symb_id);
      if (type == eExogenousDet)
        tsid += datatree.symbol_table.exo_nbr();
      if (!lhs_rhs)
        {
          if (dynamic)
            {
              if (steady_dynamic)  // steady state values in a dynamic model
                {
                  FLDVS_ fldvs(type, tsid);
799
                  fldvs.write(CompileCode, instruction_number);
800
801
802
803
804
805
                }
              else
                {
                  if (type == eParameter)
                    {
                      FLDV_ fldv(type, tsid);
806
                      fldv.write(CompileCode, instruction_number);
807
808
809
810
                    }
                  else
                    {
                      FLDV_ fldv(type, tsid, lag);
811
                      fldv.write(CompileCode, instruction_number);
812
813
814
815
816
817
                    }
                }
            }
          else
            {
              FLDSV_ fldsv(type, tsid);
818
              fldsv.write(CompileCode, instruction_number);
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
            }
        }
      else
        {
          if (dynamic)
            {
              if (steady_dynamic)  // steady state values in a dynamic model
                {
                  cerr << "Impossible case: steady_state in rhs of equation" << endl;
                  exit(EXIT_FAILURE);
                }
              else
                {
                  if (type == eParameter)
                    {
                      FSTPV_ fstpv(type, tsid);
835
                      fstpv.write(CompileCode, instruction_number);
836
837
838
839
                    }
                  else
                    {
                      FSTPV_ fstpv(type, tsid, lag);
840
                      fstpv.write(CompileCode, instruction_number);
841
842
843
844
845
846
                    }
                }
            }
          else
            {
              FSTPSV_ fstpsv(type, tsid);
847
              fstpsv.write(CompileCode, instruction_number);
848
849
850
851
            }
        }
    }
}
852

853
void
854
VariableNode::computeTemporaryTerms(map<expr_t, int> &reference_count,
855
                                    temporary_terms_t &temporary_terms,
856
                                    map<expr_t, pair<int, int> > &first_occurence,
857
                                    int Curr_block,
858
                                    vector<vector<temporary_terms_t> > &v_temporary_terms,
859
                                    int equation) const
860
861
862
863
{
  if (type == eModelLocalVariable)
    datatree.local_variables_table[symb_id]->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, v_temporary_terms, equation);
}
864

sebastien's avatar
sebastien committed
865
void
sebastien's avatar
sebastien committed
866
867
868
869
870
871
872
VariableNode::collectVariables(SymbolType type_arg, set<pair<int, int> > &result) const
{
  if (type == type_arg)
    result.insert(make_pair(symb_id, lag));
  if (type == eModelLocalVariable)
    datatree.local_variables_table[symb_id]->collectVariables(type_arg, result);
}
ferhat's avatar
ferhat committed
873

874
875
pair<int, expr_t>
VariableNode::normalizeEquation(int var_endo, vector<pair<int, pair<expr_t, expr_t> > > &List_of_Op_RHS) const
876
{
877
878
879
880
881
882
883
884
885
886
887
888
  /* The equation has to be normalized with respect to the current endogenous variable ascribed to it. 
     The two input arguments are : 
        - The ID of the endogenous variable associated to the equation.
        - The list of operators and operands needed to normalize the equation*
        
     The pair returned by NormalizeEquation is composed of 
      - a flag indicating if the expression returned contains (flag = 1) or not (flag = 0) 
        the endogenous variable related to the equation.
        If the expression contains more than one occurence of the associated endogenous variable, 
        the flag is equal to 2. 
      - an expression equal to the RHS if flag = 0 and equal to NULL elsewhere
  */
889
890
891
  if (type == eEndogenous)
    {
      if (datatree.symbol_table.getTypeSpecificID(symb_id) == var_endo && lag == 0)
892
        /* the endogenous variable */
893
        return (make_pair(1, (expr_t) NULL));
894
895
896
897
898
899
900
901
902
903
904
      else
        return (make_pair(0, datatree.AddVariableInternal(symb_id, lag)));
    }
  else
    {
      if (type == eParameter)
        return (make_pair(0, datatree.AddVariableInternal(symb_id, 0)));
      else
        return (make_pair(0, datatree.AddVariableInternal(symb_id, lag)));
    }
}
ferhat's avatar
ferhat committed
905

906
907
expr_t
VariableNode::getChainRuleDerivative(int deriv_id, const map<int, expr_t> &recursive_variables)
sebastien's avatar
sebastien committed
908
{
ferhat's avatar
ferhat committed
909
910
911
912
913
914
  switch (type)
    {
    case eEndogenous:
    case eExogenous:
    case eExogenousDet:
    case eParameter:
915
    case eTrend:
916
    case eLogTrend:
sebastien's avatar
sebastien committed
917
      if (deriv_id == datatree.getDerivID(symb_id, lag))
ferhat's avatar
ferhat committed
918
919
920
921
        return datatree.One;
      else
        {
          //if there is in the equation a recursive variable we could use a chaine rule derivation
922
          map<int, expr_t>::const_iterator it = recursive_variables.find(datatree.getDerivID(symb_id, lag));
923
          if (it != recursive_variables.end())
ferhat's avatar
ferhat committed
924
            {
925
              map<int, expr_t>::const_iterator it2 = derivatives.find(deriv_id);
926
927
928
929
              if (it2 != derivatives.end())
                return it2->second;
              else
                {
930
                  map<int, expr_t> recursive_vars2(recursive_variables);
931
                  recursive_vars2.erase(it->first);
932
                  //expr_t c = datatree.AddNonNegativeConstant("1");
933
                  expr_t d = datatree.AddUMinus(it->second->getChainRuleDerivative(deriv_id, recursive_vars2));
934
                  //d = datatree.AddTimes(c, d);
sebastien's avatar
sebastien committed
935
                  derivatives[deriv_id] = d;
936
937
                  return d;
                }
ferhat's avatar
ferhat committed
938
939
940
941
942
            }
          else
            return datatree.Zero;
        }
    case eModelLocalVariable:
sebastien's avatar
sebastien committed
943
      return datatree.local_variables_table[symb_id]->getChainRuleDerivative(deriv_id, recursive_variables);
ferhat's avatar
ferhat committed
944
945
946
    case eModFileLocalVariable:
      cerr << "ModFileLocalVariable is not derivable" << endl;
      exit(EXIT_FAILURE);
947
948
949
    case eStatementDeclaredVariable:
      cerr << "eStatementDeclaredVariable is not derivable" << endl;
      exit(EXIT_FAILURE);
950
    case eExternalFunction:
ferhat's avatar
ferhat committed
951
952
953
954
955
      cerr << "Impossible case!" << endl;
      exit(EXIT_FAILURE);
    }
  // Suppress GCC warning
  exit(EXIT_FAILURE);
sebastien's avatar
sebastien committed
956
957
}

958
expr_t
ferhat's avatar
ferhat committed
959
VariableNode::toStatic(DataTree &static_datatree) const
960
{
961
  return static_datatree.AddVariable(symb_id);
962
}
ferhat's avatar
ferhat committed
963

964
965
966
967
968
969
expr_t
VariableNode::cloneDynamic(DataTree &dynamic_datatree) const
{
  return dynamic_datatree.AddVariable(symb_id, lag);
}

sebastien's avatar
sebastien committed
970
971
972
int
VariableNode::maxEndoLead() const
{
973
  switch (type)
sebastien's avatar
sebastien committed
974
975
976
977
978
979
980
981
982
983
    {
    case eEndogenous:
      return max(lag, 0);
    case eModelLocalVariable:
      return datatree.local_variables_table[symb_id]->maxEndoLead();
    default:
      return 0;
    }
}

sebastien's avatar
sebastien committed
984
985
986
int
VariableNode::maxExoLead() const
{
987
  switch (type)
sebastien's avatar
sebastien committed
988
989
990
991
992
993
994
995
996
997
    {
    case eExogenous:
      return max(lag, 0);
    case eModelLocalVariable:
      return datatree.local_variables_table[symb_id]->maxExoLead();
    default:
      return 0;
    }
}

998
999
1000
1001
1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025
int
VariableNode::maxEndoLag() const
{
  switch (type)
    {
    case eEndogenous:
      return max(-lag, 0);
    case eModelLocalVariable:
      return datatree.local_variables_table[symb_id]->maxEndoLag();
    default:
      return 0;
    }
}

int
VariableNode::maxExoLag() const
{
  switch (type)
    {
    case eExogenous:
      return max(-lag, 0);
    case eModelLocalVariable:
      return datatree.local_variables_table[symb_id]->maxExoLag();
    default:
      return 0;
    }
}

1026
1027
1028
1029
1030
1031
1032
1033
1034
1035
1036
1037
1038
1039
1040
1041
int
VariableNode::maxLead() const
{
  switch (type)
    {
    case eEndogenous:
      return lag;
    case eExogenous:
      return lag;
    case eModelLocalVariable:
      return datatree.local_variables_table[symb_id]->maxLead();
    default:
      return 0;
    }
}

1042
expr_t
sebastien's avatar
sebastien committed
1043
1044
VariableNode::decreaseLeadsLags(int n) const
{
1045
  switch (type)
sebastien's avatar
sebastien committed
1046
1047
1048
1049
    {
    case eEndogenous:
    case eExogenous:
    case eExogenousDet:
1050
    case eTrend:
1051
    case eLogTrend:
sebastien's avatar
sebastien committed
1052
1053
1054
1055
1056
1057
1058
1059
      return datatree.AddVariable(symb_id, lag-n);
    case eModelLocalVariable:
      return datatree.local_variables_table[symb_id]->decreaseLeadsLags(n);
    default:
      return const_cast<VariableNode *>(this);
    }
}

1060
expr_t
sebastien's avatar
sebastien committed
1061
VariableNode::decreaseLeadsLagsPredeterminedVariables() const
1062
{
sebastien's avatar
sebastien committed
1063
  if (datatree.symbol_table.isPredetermined(symb_id))
1064
1065
1066
1067
1068
    return decreaseLeadsLags(1);
  else
    return const_cast<VariableNode *>(this);
}

1069
expr_t
1070
VariableNode::substituteEndoLeadGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool deterministic_model) const