ComputingTasks.cc 99.1 KB
Newer Older
1
/*
2
 * Copyright (C) 2003-2014 Dynare Team
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/>.
 */

20
#include <cstdlib>
21
#include <cassert>
22
23
24
25
26
27
28
29
#include <iostream>
#include <sstream>

using namespace std;

#include "ComputingTasks.hh"
#include "Statement.hh"

30
31
32
33
34
#include <boost/algorithm/string/trim.hpp>
#include <boost/algorithm/string/split.hpp>
#include <boost/lexical_cast.hpp>
#include <boost/tokenizer.hpp>

35
36
SteadyStatement::SteadyStatement(const OptionsList &options_list_arg) :
  options_list(options_list_arg)
37
38
39
{
}

40
void
41
SteadyStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
42
{
43
  mod_file_struct.steady_present = true;
44
45
}

46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
void
SteadyStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);
  output << "steady;\n";
}

CheckStatement::CheckStatement(const OptionsList &options_list_arg) :
  options_list(options_list_arg)
{
}

void
CheckStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);
62
  output << "oo_.dr.eigval = check(M_,options_,oo_);" << endl;
63
64
65
}

void
66
CheckStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
67
68
69
70
{
  mod_file_struct.check_present = true;
}

ferhat's avatar
ferhat committed
71
ModelInfoStatement::ModelInfoStatement(const OptionsList &options_list_arg) :
72
73
74
75
  options_list(options_list_arg)
{
}

76
void
77
ModelInfoStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
78
{
ferhat's avatar
ferhat committed
79
  //mod_file_struct.model_info_present = true;
80
81
}

82
83
void
ModelInfoStatement::writeOutput(ostream &output, const string &basename) const
84
85
86
87
88
{
  options_list.writeOutput(output);
  output << "model_info();\n";
}

89
90
SimulStatement::SimulStatement(const OptionsList &options_list_arg) :
  options_list(options_list_arg)
91
92
93
94
{
}

void
95
SimulStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
96
97
{
  mod_file_struct.simul_present = true;
98
99
100

  // The following is necessary to allow shocks+endval+simul in a loop
  mod_file_struct.shocks_present_but_simul_not_yet = false;
101
102
103
104
105
106
}

void
SimulStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);
107
  output << "simul();\n";
108
109
}

sebastien's avatar
sebastien committed
110
StochSimulStatement::StochSimulStatement(const SymbolList &symbol_list_arg,
111
                                         const OptionsList &options_list_arg) :
sebastien's avatar
sebastien committed
112
  symbol_list(symbol_list_arg),
113
  options_list(options_list_arg)
114
115
116
117
{
}

void
118
StochSimulStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
119
{
sebastien's avatar
sebastien committed
120
  mod_file_struct.stoch_simul_present = true;
121
122

  // Fill in option_order of mod_file_struct
123
  OptionsList::num_options_t::const_iterator it = options_list.num_options.find("order");
124
  if (it != options_list.num_options.end())
125
    mod_file_struct.order_option = max(mod_file_struct.order_option, atoi(it->second.c_str()));
sebastien's avatar
trunk:    
sebastien committed
126

127
128
129
130
  // Fill in mod_file_struct.partial_information
  it = options_list.num_options.find("partial_information");
  if (it != options_list.num_options.end() && it->second == "1")
    mod_file_struct.partial_information = true;
sebastien's avatar
sebastien committed
131
132
133
134
135
136

  // Option k_order_solver (implicit when order >= 3)
  it = options_list.num_options.find("k_order_solver");
  if ((it != options_list.num_options.end() && it->second == "1")
      || mod_file_struct.order_option >= 3)
    mod_file_struct.k_order_solver = true;
137

138
139
140
141
142
143
}

void
StochSimulStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);
sebastien's avatar
sebastien committed
144
  symbol_list.writeOutput("var_list_", output);
145
  output << "info = stoch_simul(var_list_);" << endl;
146
147
}

148
ForecastStatement::ForecastStatement(const SymbolList &symbol_list_arg,
149
                                     const OptionsList &options_list_arg) :
150
151
152
153
154
155
156
157
158
159
  symbol_list(symbol_list_arg),
  options_list(options_list_arg)
{
}

void
ForecastStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);
  symbol_list.writeOutput("var_list_", output);
160
  output << "info = dyn_forecast(var_list_,'simul');" << endl;
161
162
}

163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
RamseyModelStatement::RamseyModelStatement(const SymbolList &symbol_list_arg,
                                             const OptionsList &options_list_arg) :
  symbol_list(symbol_list_arg),
  options_list(options_list_arg)
{
}

void
RamseyModelStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
{
  mod_file_struct.ramsey_model_present = true;

  /* Fill in option_order of mod_file_struct
     Since ramsey model needs one further order of derivation (for example, for 1st order
     approximation, it needs 2nd derivatives), we add 1 to the order declared by user */
  OptionsList::num_options_t::const_iterator it = options_list.num_options.find("order");
  if (it != options_list.num_options.end())
    {
      int order = atoi(it->second.c_str());
      if (order > 2)
        {
          cerr << "ERROR: ramsey_model: order > 2 is not  implemented" << endl;
          exit(EXIT_FAILURE);
        }
      mod_file_struct.order_option = max(mod_file_struct.order_option, order + 1);
    }

  // Fill in mod_file_struct.partial_information
  it = options_list.num_options.find("partial_information");
  if (it != options_list.num_options.end() && it->second == "1")
    mod_file_struct.partial_information = true;

  // Option k_order_solver (implicit when order >= 3)
  it = options_list.num_options.find("k_order_solver");
  if ((it != options_list.num_options.end() && it->second == "1")
      || mod_file_struct.order_option >= 3)
    mod_file_struct.k_order_solver = true;
}

void
RamseyModelStatement::writeOutput(ostream &output, const string &basename) const
{
205
206
207
208
  // options_.ramsey_policy indicates that a Ramsey model is present in the *.mod file
  // this affects the computation of the steady state that uses a special algorithm
  // It should probably rather be a M_ field, but we leave it in options_ for historical reason
  output << "options_.ramsey_policy = 1;\n";
209
210
}

sebastien's avatar
sebastien committed
211
RamseyPolicyStatement::RamseyPolicyStatement(const SymbolList &symbol_list_arg,
212
                                             const OptionsList &options_list_arg) :
sebastien's avatar
sebastien committed
213
  symbol_list(symbol_list_arg),
214
215
216
217
218
  options_list(options_list_arg)
{
}

void
219
RamseyPolicyStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
220
{
221
222
223
224
  // ramsey_model_present indicates that the model is augmented with the FOC of the planner problem
  mod_file_struct.ramsey_model_present = true;
  // ramsey_policy_present indicates that ramsey_policy instruction for computation of first order approximation 
  // of  a stochastic Ramsey problem if present in the *.mod file
sebastien's avatar
sebastien committed
225
  mod_file_struct.ramsey_policy_present = true;
226
227
228
229

  /* Fill in option_order of mod_file_struct
     Since ramsey policy needs one further order of derivation (for example, for 1st order
     approximation, it needs 2nd derivatives), we add 1 to the order declared by user */
230
  OptionsList::num_options_t::const_iterator it = options_list.num_options.find("order");
231
  if (it != options_list.num_options.end())
232
233
    {
      int order = atoi(it->second.c_str());
234
      if (order > 2)
235
        {
236
          cerr << "ERROR: ramsey_policy: order > 2 is not  implemented" << endl;
237
238
239
240
          exit(EXIT_FAILURE);
        }
      mod_file_struct.order_option = max(mod_file_struct.order_option, order + 1);
    }
sebastien's avatar
trunk:    
sebastien committed
241

242
243
244
245
  // Fill in mod_file_struct.partial_information
  it = options_list.num_options.find("partial_information");
  if (it != options_list.num_options.end() && it->second == "1")
    mod_file_struct.partial_information = true;
sebastien's avatar
sebastien committed
246
247
248
249
250
251

  // Option k_order_solver (implicit when order >= 3)
  it = options_list.num_options.find("k_order_solver");
  if ((it != options_list.num_options.end() && it->second == "1")
      || mod_file_struct.order_option >= 3)
    mod_file_struct.k_order_solver = true;
252
253
254
255
256
257
}

void
RamseyPolicyStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);
sebastien's avatar
sebastien committed
258
  symbol_list.writeOutput("var_list_", output);
259
260
261
  output << "ramsey_policy(var_list_);\n";
}

262
263
264
265
266
267
268
269
DiscretionaryPolicyStatement::DiscretionaryPolicyStatement(const SymbolList &symbol_list_arg,
							   const OptionsList &options_list_arg) :
  symbol_list(symbol_list_arg),
  options_list(options_list_arg)
{
}

void
270
DiscretionaryPolicyStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
271
272
273
{
  mod_file_struct.discretionary_policy_present = true;

274
275
276
277
278
279
  if (options_list.symbol_list_options.find("instruments") == options_list.symbol_list_options.end())
    {
      cerr << "ERROR: discretionary_policy: the instruments option is required." << endl;
      exit(EXIT_FAILURE);
    }

280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
  /* Fill in option_order of mod_file_struct
     Since discretionary policy needs one further order of derivation (for example, for 1st order
     approximation, it needs 2nd derivatives), we add 1 to the order declared by user */
  OptionsList::num_options_t::const_iterator it = options_list.num_options.find("order");
  if (it != options_list.num_options.end())
    {
      int order = atoi(it->second.c_str());
      if (order > 1)
        {
          cerr << "ERROR: discretionary_policy: order > 1 is not yet implemented" << endl;
          exit(EXIT_FAILURE);
        }
      mod_file_struct.order_option = max(mod_file_struct.order_option, order + 1);
    }

  // Fill in mod_file_struct.partial_information
  it = options_list.num_options.find("partial_information");
  if (it != options_list.num_options.end() && it->second == "1")
    mod_file_struct.partial_information = true;

  // Option k_order_solver (implicit when order >= 3)
  it = options_list.num_options.find("k_order_solver");
  if ((it != options_list.num_options.end() && it->second == "1")
      || mod_file_struct.order_option >= 3)
    mod_file_struct.k_order_solver = true;
}

void
DiscretionaryPolicyStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);
  symbol_list.writeOutput("var_list_", output);
  output << "discretionary_policy(var_list_);\n";
}

sebastien's avatar
sebastien committed
315
EstimationStatement::EstimationStatement(const SymbolList &symbol_list_arg,
316
                                         const OptionsList &options_list_arg) :
sebastien's avatar
sebastien committed
317
  symbol_list(symbol_list_arg),
318
  options_list(options_list_arg)
319
320
321
322
{
}

void
323
EstimationStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
324
{
sebastien's avatar
sebastien committed
325
  mod_file_struct.estimation_present = true;
326
327

  // Fill in option_order of mod_file_struct
328
  OptionsList::num_options_t::const_iterator it = options_list.num_options.find("order");
329
  if (it != options_list.num_options.end())
330
331
332
333
334
335
336
337
338
339
340
341
    {
      int order = atoi(it->second.c_str());
          
      if (order > 2)
        {
          cerr << "ERROR: order > 2 is not supported in estimation" << endl;
          exit(EXIT_FAILURE);
        }
      
      mod_file_struct.order_option = max(mod_file_struct.order_option, order);
    }
  
342
343
344
345
  // Fill in mod_file_struct.partial_information
  it = options_list.num_options.find("partial_information");
  if (it != options_list.num_options.end() && it->second == "1")
    mod_file_struct.partial_information = true;
346

347
348
349
350
351
  // Fill in mod_file_struct.estimation_analytic_derivation
  it = options_list.num_options.find("analytic_derivation");
  if (it != options_list.num_options.end() && it->second == "1")
    mod_file_struct.estimation_analytic_derivation = true;

352
353
  it = options_list.num_options.find("dsge_var");
  if (it != options_list.num_options.end())
354
355
356
357
358
359
360
361
362
    // Ensure that irf_shocks & dsge_var have not both been passed
    if (options_list.symbol_list_options.find("irf_shocks") != options_list.symbol_list_options.end())
      {
        cerr << "The irf_shocks and dsge_var options may not both be passed to estimation." << endl;
        exit(EXIT_FAILURE);
      }
    else
      // Fill in mod_file_struct.dsge_var_calibrated
      mod_file_struct.dsge_var_calibrated = it->second;
363
364

  // Fill in mod_file_struct.dsge_var_estimated
365
  OptionsList::string_options_t::const_iterator it_str = options_list.string_options.find("dsge_var");
366
367
368
369
370
371
372
373
374
375
  if (it_str != options_list.string_options.end())
    mod_file_struct.dsge_var_estimated = true;

  // Fill in mod_file_struct.bayesian_irf_present
  it = options_list.num_options.find("bayesian_irf");
  if (it != options_list.num_options.end() && it->second == "1")
    mod_file_struct.bayesian_irf_present = true;

  it = options_list.num_options.find("dsge_varlag");
  if (it != options_list.num_options.end())
376
377
    if (mod_file_struct.dsge_var_calibrated.empty()
        && !mod_file_struct.dsge_var_estimated)
378
379
380
381
382
383
      {
        cerr << "ERROR: The estimation statement requires a dsge_var option to be passed "
             << "if the dsge_varlag option is passed." << endl;
        exit(EXIT_FAILURE);
      }

384
385
  if (!mod_file_struct.dsge_var_calibrated.empty()
      && mod_file_struct.dsge_var_estimated)
386
387
388
389
    {
      cerr << "ERROR: An estimation statement cannot take more than one dsge_var option." << endl;
      exit(EXIT_FAILURE);
    }
390
391
392
393

  if (options_list.string_options.find("datafile") == options_list.string_options.end() &&
      !mod_file_struct.estimation_data_statement_present)
    {
394
      cerr << "ERROR: The estimation statement requires a data file to be supplied via the datafile option." << endl;
395
396
      exit(EXIT_FAILURE);
    }
397
398
399
400
401
402
403

  if (options_list.string_options.find("mode_file") != options_list.string_options.end() &&
      mod_file_struct.estim_params_use_calib)
    {
      cerr << "ERROR: The mode_file option of the estimation statement is incompatible with the use_calibration option of the estimated_params_init block." << endl;
      exit(EXIT_FAILURE);
    }
404
405
406
407
408
409
}

void
EstimationStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);
410
411
412
413
414
415
416
417

  // Special treatment for order option and particle filter
  OptionsList::num_options_t::const_iterator it = options_list.num_options.find("order");
  if (it == options_list.num_options.end())
    output << "options_.order = 1;" << endl;
  else if (atoi(it->second.c_str()) == 2)
    output << "options_.particle.status = 1;" << endl;

418
419
420
421
422
  // Do not check for the steady state in diffuse filter mode (#400)
  it = options_list.num_options.find("diffuse_filter");
  if (it != options_list.num_options.end() && it->second == "1")
    output << "options_.steadystate.nocheck = 1;" << endl;

sebastien's avatar
sebastien committed
423
  symbol_list.writeOutput("var_list_", output);
424
425
426
427
428
429
430
431
  output << "dynare_estimation(var_list_);\n";
}

DynareSensitivityStatement::DynareSensitivityStatement(const OptionsList &options_list_arg) :
  options_list(options_list_arg)
{
}

sebastien's avatar
sebastien committed
432
void
433
DynareSensitivityStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
sebastien's avatar
sebastien committed
434
{
435
  OptionsList::num_options_t::const_iterator it = options_list.num_options.find("identification");
sebastien's avatar
sebastien committed
436
437
438
439
440
  if (it != options_list.num_options.end()
      && it->second == "1")
    mod_file_struct.identification_present = true;
}

441
442
443
void
DynareSensitivityStatement::writeOutput(ostream &output, const string &basename) const
{
444
  options_list.writeOutput(output, "options_gsa");
445
446
447
448
449
450
451

  /* Ensure that nograph, nodisplay and graph_format are also set in top-level
     options_.
     \todo factorize this code between identification and dynare_sensitivity,
     and provide a generic mechanism for this situation (maybe using regexps) */
  OptionsList::num_options_t::const_iterator it = options_list.num_options.find("nodisplay");
  if (it != options_list.num_options.end())
452
    output << "options_.nodisplay = " << it->second << ";" << endl;
453
454
  it = options_list.num_options.find("nograph");
  if (it != options_list.num_options.end())
455
    output << "options_.nograph = " << it->second << ";" << endl;
456
457
  OptionsList::string_options_t::const_iterator it2 = options_list.string_options.find("graph_format");
  if (it2 != options_list.string_options.end())
458
    output << "options_.graph_format = '" << it2->second << "';" << endl;
459
  
460
  output << "dynare_sensitivity(options_gsa);" << endl;
461
462
}

sebastien's avatar
sebastien committed
463
RplotStatement::RplotStatement(const SymbolList &symbol_list_arg,
464
                               const OptionsList &options_list_arg) :
sebastien's avatar
sebastien committed
465
  symbol_list(symbol_list_arg),
466
467
468
469
470
471
472
473
  options_list(options_list_arg)
{
}

void
RplotStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);
sebastien's avatar
sebastien committed
474
  symbol_list.writeOutput("var_list_", output);
475
476
477
  output << "rplot(var_list_);\n";
}

478
479
480
481
482
483
484
485
486
487
488
UnitRootVarsStatement::UnitRootVarsStatement(void)
{
}

void
UnitRootVarsStatement::writeOutput(ostream &output, const string &basename) const
{
  output << "options_.diffuse_filter = 1;" << endl
	 << "options_.steadystate.nocheck = 1;" << endl;
}

489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
PeriodsStatement::PeriodsStatement(int periods_arg) : periods(periods_arg)
{
}

void
PeriodsStatement::writeOutput(ostream &output, const string &basename) const
{
  output << "options_.periods = " << periods << ";" << endl;
}

DsampleStatement::DsampleStatement(int val1_arg) : val1(val1_arg), val2(-1)
{
}

DsampleStatement::DsampleStatement(int val1_arg, int val2_arg) : val1(val1_arg), val2(val2_arg)
{
}

void
DsampleStatement::writeOutput(ostream &output, const string &basename) const
{
  if (val2 < 0)
511
    output << "dsample(" << val1 << ");" << endl;
512
  else
513
    output << "dsample(" << val1 << ", " << val2 << ");" << endl;
514
515
516
517
518
519
520
521
522
}

EstimatedParamsStatement::EstimatedParamsStatement(const vector<EstimationParams> &estim_params_list_arg,
                                                   const SymbolTable &symbol_table_arg) :
  estim_params_list(estim_params_list_arg),
  symbol_table(symbol_table_arg)
{
}

523
void
524
EstimatedParamsStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
525
526
527
{
  for (vector<EstimationParams>::const_iterator it = estim_params_list.begin();
       it != estim_params_list.end(); it++)
528
529
530
531
    {
      if (it->name == "dsge_prior_weight")
        mod_file_struct.dsge_prior_weight_in_estimated_params = true;

532
      // Handle case of degenerate beta prior
533
      if (it->prior == eBeta)
534
        try
535
          {
536
537
            if (it->mean->eval(eval_context_t()) == 0.5
                && it->std->eval(eval_context_t()) == 0.5)
538
539
540
541
542
543
544
545
              {
                cerr << "ERROR: The prior density is not defined for the beta distribution when the mean = standard deviation = 0.5." << endl;
                exit(EXIT_FAILURE);
              }
          }
        catch (ExprNode::EvalException &e)
          {
            // We don't have enough information to compute the numerical value, skip the test
546
547
          }
    }
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578

  // Check that no parameter/endogenous is declared twice in the block
  set<string> already_declared;
  set<pair<string, string> > already_declared_corr;
  for (vector<EstimationParams>::const_iterator it = estim_params_list.begin();
       it != estim_params_list.end(); it++)
    {
      if (it->type == 3) // Correlation
        {
          // Use lexical ordering for the pair of symbols
          pair<string, string> x = it->name < it->name2 ? make_pair(it->name, it->name2) : make_pair(it->name2, it->name);

          if (already_declared_corr.find(x) == already_declared_corr.end())
            already_declared_corr.insert(x);
          else
            {
              cerr << "ERROR: in `estimated_params' block, the correlation between " << it->name << " and " << it->name2 << " is declared twice." << endl;
              exit(EXIT_FAILURE);
            }
        }
      else
        {
          if (already_declared.find(it->name) == already_declared.end())
            already_declared.insert(it->name);
          else
            {
              cerr << "ERROR: in `estimated_params' block, the symbol " << it->name << " is declared twice." << endl;
              exit(EXIT_FAILURE);
            }
        }
    }
579
580
581
582

  // Fill in mod_file_struct.estimated_parameters (related to #469)
  for (vector<EstimationParams>::const_iterator it = estim_params_list.begin();
       it != estim_params_list.end(); it++)
583
    if (it->type == 2 && it->name != "dsge_prior_weight")
584
      mod_file_struct.estimated_parameters.insert(symbol_table.getID(it->name));
585
586
}

587
588
589
void
EstimatedParamsStatement::writeOutput(ostream &output, const string &basename) const
{
590
591
592
593
594
  output << "global estim_params_" << endl
         << "estim_params_.var_exo = [];" << endl
         << "estim_params_.var_endo = [];" << endl
         << "estim_params_.corrx = [];" << endl
         << "estim_params_.corrn = [];" << endl
Sébastien Villemot's avatar
Sébastien Villemot committed
595
         << "estim_params_.param_vals = [];" << endl;
596
597
598

  vector<EstimationParams>::const_iterator it;

599
  for (it = estim_params_list.begin(); it != estim_params_list.end(); it++)
600
    {
sebastien's avatar
sebastien committed
601
      int symb_id = symbol_table.getTypeSpecificID(it->name) + 1;
602
603
      SymbolType symb_type = symbol_table.getType(it->name);

604
      switch (it->type)
605
606
        {
        case 1:
607
          if (symb_type == eExogenous)
608
            output << "estim_params_.var_exo = [estim_params_.var_exo; ";
609
          else if (symb_type == eEndogenous)
610
            output << "estim_params_.var_endo = [estim_params_.var_endo; ";
611
          output << symb_id;
612
613
          break;
        case 2:
614
615
          output << "estim_params_.param_vals = [estim_params_.param_vals; "
                 << symb_id;
616
617
          break;
        case 3:
618
          if (symb_type == eExogenous)
619
            output << "estim_params_.corrx = [estim_params_.corrx; ";
620
          else if (symb_type == eEndogenous)
621
            output << "estim_params_.corrn = [estim_params_.corrn; ";
sebastien's avatar
sebastien committed
622
          output << symb_id << " " << symbol_table.getTypeSpecificID(it->name2)+1;
623
624
          break;
        }
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
      output << ", ";
      it->init_val->writeOutput(output);
      output << ", ";
      it->low_bound->writeOutput(output);
      output << ", ";
      it->up_bound->writeOutput(output);
      output << ", "
             << it->prior << ", ";
      it->mean->writeOutput(output);
      output << ", ";
      it->std->writeOutput(output);
      output << ", ";
      it->p3->writeOutput(output);
      output << ", ";
      it->p4->writeOutput(output);
      output << ", ";
      it->jscale->writeOutput(output);
      output << " ];" << endl;
643
644
645
646
    }
}

EstimatedParamsInitStatement::EstimatedParamsInitStatement(const vector<EstimationParams> &estim_params_list_arg,
647
648
                                                           const SymbolTable &symbol_table_arg,
                                                           const bool use_calibration_arg) :
649
  estim_params_list(estim_params_list_arg),
650
651
  symbol_table(symbol_table_arg),
  use_calibration(use_calibration_arg)
652
653
654
{
}

655
656
657
658
659
660
661
void
EstimatedParamsInitStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
{
  if (use_calibration)
    mod_file_struct.estim_params_use_calib = true;
}

662
663
664
void
EstimatedParamsInitStatement::writeOutput(ostream &output, const string &basename) const
{
665
666
667
  if (use_calibration)
    output << "options_.use_calibration_initialization = 1;" << endl;

668
669
  vector<EstimationParams>::const_iterator it;

670
  for (it = estim_params_list.begin(); it != estim_params_list.end(); it++)
671
    {
sebastien's avatar
sebastien committed
672
      int symb_id = symbol_table.getTypeSpecificID(it->name) + 1;
673
674
      SymbolType symb_type = symbol_table.getType(it->name);

675
676
      if (it->type < 3)
        {
677
          if (symb_type == eExogenous)
678
            {
679
680
681
682
              output << "tmp1 = find(estim_params_.var_exo(:,1)==" << symb_id << ");" << endl;
              output << "estim_params_.var_exo(tmp1,2) = ";
              it->init_val->writeOutput(output);
              output << ";" << endl;
683
            }
684
          else if (symb_type == eEndogenous)
685
            {
686
687
688
689
              output << "tmp1 = find(estim_params_.var_endo(:,1)==" << symb_id << ");" << endl;
              output << "estim_params_.var_endo(tmp1,2) = ";
              it->init_val->writeOutput(output);
              output << ";" << endl;
690
            }
691
          else if (symb_type == eParameter)
692
            {
693
694
695
696
              output << "tmp1 = find(estim_params_.param_vals(:,1)==" << symb_id << ");" << endl;
              output << "estim_params_.param_vals(tmp1,2) = ";
              it->init_val->writeOutput(output);
              output << ";" << endl;
697
698
699
700
            }
        }
      else
        {
701
          if (symb_type == eExogenous)
702
            {
703
704
              output << "tmp1 = find((estim_params_.corrx(:,1)==" << symb_id << " & estim_params_.corrx(:,2)==" << symbol_table.getTypeSpecificID(it->name2)+1 << ") | "
                     <<             "(estim_params_.corrx(:,2)==" << symb_id << " & estim_params_.corrx(:,1)==" << symbol_table.getTypeSpecificID(it->name2)+1 << "));" << endl;
705
706
707
              output << "estim_params_.corrx(tmp1,3) = ";
              it->init_val->writeOutput(output);
              output << ";" << endl;
708
            }
709
          else if (symb_type == eEndogenous)
710
            {
711
712
              output << "tmp1 = find((estim_params_.corrn(:,1)==" << symb_id << " & estim_params_.corrn(:,2)==" << symbol_table.getTypeSpecificID(it->name2)+1 << ") | "
                     <<             "(estim_params_.corrn(:,2)==" << symb_id << " & estim_params_.corrn(:,1)==" << symbol_table.getTypeSpecificID(it->name2)+1 << "));" << endl;
713
714
715
              output << "estim_params_.corrn(tmp1,3) = ";
              it->init_val->writeOutput(output);
              output << ";" << endl;
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
            }
        }
    }
}

EstimatedParamsBoundsStatement::EstimatedParamsBoundsStatement(const vector<EstimationParams> &estim_params_list_arg,
                                                               const SymbolTable &symbol_table_arg) :
  estim_params_list(estim_params_list_arg),
  symbol_table(symbol_table_arg)
{
}

void
EstimatedParamsBoundsStatement::writeOutput(ostream &output, const string &basename) const
{
  vector<EstimationParams>::const_iterator it;

733
  for (it = estim_params_list.begin(); it != estim_params_list.end(); it++)
734
    {
sebastien's avatar
sebastien committed
735
      int symb_id = symbol_table.getTypeSpecificID(it->name) + 1;
736
737
      SymbolType symb_type = symbol_table.getType(it->name);

738
739
      if (it->type < 3)
        {
740
          if (symb_type == eExogenous)
741
            {
742
743
744
745
746
747
748
749
750
              output << "tmp1 = find(estim_params_.var_exo(:,1)==" << symb_id << ");" << endl;

              output << "estim_params_.var_exo(tmp1,3) = ";
              it->low_bound->writeOutput(output);
              output << ";" << endl;

              output << "estim_params_.var_exo(tmp1,4) = ";
              it->up_bound->writeOutput(output);
              output << ";" << endl;
751
            }
752
          else if (symb_type == eEndogenous)
753
            {
754
755
756
757
758
759
760
761
762
              output << "tmp1 = find(estim_params_.var_endo(:,1)==" << symb_id << ");" << endl;

              output << "estim_params_.var_endo(tmp1,3) = ";
              it->low_bound->writeOutput(output);
              output << ";" << endl;

              output << "estim_params_.var_endo(tmp1,4) = ";
              it->up_bound->writeOutput(output);
              output << ";" << endl;
763
            }
764
          else if (symb_type == eParameter)
765
            {
766
767
768
769
770
771
772
773
774
              output << "tmp1 = find(estim_params_.param_vals(:,1)==" << symb_id << ");" << endl;

              output << "estim_params_.param_vals(tmp1,3) = ";
              it->low_bound->writeOutput(output);
              output << ";" << endl;

              output << "estim_params_.param_vals(tmp1,4) = ";
              it->up_bound->writeOutput(output);
              output << ";" << endl;
775
776
777
778
            }
        }
      else
        {
779
          if (symb_type == eExogenous)
780
            {
781
782
              output << "tmp1 = find((estim_params_.corrx(:,1)==" << symb_id << " & estim_params_.corrx(:,2)==" << symbol_table.getTypeSpecificID(it->name2)+1 << ") | "
                     <<             "(estim_params_.corrx(:,2)==" << symb_id << " & estim_params_.corrx(:,1)==" << symbol_table.getTypeSpecificID(it->name2)+1 << "));" << endl;
783
784
785
786
787
788
789
790

              output << "estim_params_.corrx(tmp1,4) = ";
              it->low_bound->writeOutput(output);
              output << ";" << endl;

              output << "estim_params_.corrx(tmp1,5) = ";
              it->up_bound->writeOutput(output);
              output << ";" << endl;
791
            }
792
          else if (symb_type == eEndogenous)
793
            {
794
795
              output << "tmp1 = find((estim_params_.corrn(:,1)==" << symb_id << " & estim_params_.corrn(:,2)==" << symbol_table.getTypeSpecificID(it->name2)+1 << ") | "
                     <<             "(estim_params_.corrn(:,2)==" << symb_id << " & estim_params_.corrn(:,1)==" << symbol_table.getTypeSpecificID(it->name2)+1 << "));" << endl;
796
797
798
799
800
801
802
803

              output << "estim_params_.corrn(tmp1,4) = ";
              it->low_bound->writeOutput(output);
              output << ";" << endl;

              output << "estim_params_.corrn(tmp1,5) = ";
              it->up_bound->writeOutput(output);
              output << ";" << endl;
804
805
806
807
808
            }
        }
    }
}

809
ObservationTrendsStatement::ObservationTrendsStatement(const trend_elements_t &trend_elements_arg,
810
811
812
813
814
815
816
817
818
819
820
                                                       const SymbolTable &symbol_table_arg) :
  trend_elements(trend_elements_arg),
  symbol_table(symbol_table_arg)
{
}

void
ObservationTrendsStatement::writeOutput(ostream &output, const string &basename) const
{
  output << "options_.trend_coeff_ = {};" << endl;

821
  trend_elements_t::const_iterator it;
822

823
  for (it = trend_elements.begin(); it != trend_elements.end(); it++)
824
    {
825
      SymbolType type = symbol_table.getType(it->first);
826
827
828
829
830
831
832
833
834
835
836
837
      if (type == eEndogenous)
        {
          output << "tmp1 = strmatch('" << it->first << "',options_.varobs,'exact');\n";
          output << "options_.trend_coeffs{tmp1} = '";
          it->second->writeOutput(output);
          output << "';" << endl;
        }
      else
        cout << "Error : Non-variable symbol used in TREND_COEFF: " << it->first << endl;
    }
}

sebastien's avatar
sebastien committed
838
839
OsrParamsStatement::OsrParamsStatement(const SymbolList &symbol_list_arg) :
  symbol_list(symbol_list_arg)
840
841
842
{
}

843
void
844
OsrParamsStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
845
846
847
848
{
  mod_file_struct.osr_params_present = true;
}

849
850
851
void
OsrParamsStatement::writeOutput(ostream &output, const string &basename) const
{
sebastien's avatar
sebastien committed
852
  symbol_list.writeOutput("osr_params_", output);
853
854
}

sebastien's avatar
sebastien committed
855
OsrStatement::OsrStatement(const SymbolList &symbol_list_arg,
856
                           const OptionsList &options_list_arg) :
sebastien's avatar
sebastien committed
857
  symbol_list(symbol_list_arg),
858
859
860
861
862
  options_list(options_list_arg)
{
}

void
863
OsrStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
864
{
sebastien's avatar
sebastien committed
865
  mod_file_struct.osr_present = true;
866
867

  // Fill in option_order of mod_file_struct
868
  OptionsList::num_options_t::const_iterator it = options_list.num_options.find("order");
869
  if (it != options_list.num_options.end())
870
    mod_file_struct.order_option = max(mod_file_struct.order_option, atoi(it->second.c_str()));
sebastien's avatar
trunk:    
sebastien committed
871

872
873
874
875
  // Fill in mod_file_struct.partial_information
  it = options_list.num_options.find("partial_information");
  if (it != options_list.num_options.end() && it->second == "1")
    mod_file_struct.partial_information = true;
sebastien's avatar
sebastien committed
876
877
878
879
880
881

  // Option k_order_solver (implicit when order >= 3)
  it = options_list.num_options.find("k_order_solver");
  if ((it != options_list.num_options.end() && it->second == "1")
      || mod_file_struct.order_option >= 3)
    mod_file_struct.k_order_solver = true;
882
883
884
885
886
887
}

void
OsrStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);
sebastien's avatar
sebastien committed
888
  symbol_list.writeOutput("var_list_", output);
889
  output << "oo_.osr = osr(var_list_,osr_params_,obj_var_,optim_weights_);\n";
890
891
}

892
893
OptimWeightsStatement::OptimWeightsStatement(const var_weights_t &var_weights_arg,
                                             const covar_weights_t &covar_weights_arg,
894
895
896
897
898
899
900
                                             const SymbolTable &symbol_table_arg) :
  var_weights(var_weights_arg),
  covar_weights(covar_weights_arg),
  symbol_table(symbol_table_arg)
{
}

901
void
902
OptimWeightsStatement::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings)
903
904
905
906
{
  mod_file_struct.optim_weights_present = true;
}

907
908
909
void
OptimWeightsStatement::writeOutput(ostream &output, const string &basename) const
{
sebastien's avatar
sebastien committed
910
911
912
913
914
  output << "%" << endl
         << "% OPTIM_WEIGHTS" << endl
         << "%" << endl
         << "optim_weights_ = sparse(M_.endo_nbr,M_.endo_nbr);" << endl
         << "obj_var_ = [];" << endl << endl;
915

916
  for (var_weights_t::const_iterator it = var_weights.begin();
917
       it != var_weights.end(); it++)
918
919
    {
      const string &name = it->first;
920
      const expr_t value = it->second;
sebastien's avatar
sebastien committed
921
      int id = symbol_table.getTypeSpecificID(name) + 1;
Houtan Bastani's avatar
Houtan Bastani committed
922
      output << "optim_weights_(" << id << "," << id << ") = ";
923
924
925
926
927
      value->writeOutput(output);
      output << ";" << endl;
      output << "obj_var_ = [obj_var_; " << id << "];\n";
    }

928
  for (covar_weights_t::const_iterator it = covar_weights.begin();
929
       it != covar_weights.end(); it++)
930
931
932
    {
      const string &name1 = it->first.first;
      const string &name2 = it->first.second;
933
      const expr_t value = it->second;
sebastien's avatar
sebastien committed
934
935
      int id1 = symbol_table.getTypeSpecificID(name1) + 1;
      int id2 = symbol_table.getTypeSpecificID(name2) + 1;
Houtan Bastani's avatar
Houtan Bastani committed
936
      output << "optim_weights_(" << id1 << "," << id2 << ") = ";
937
938
      value->writeOutput(output);
      output << ";" << endl;
939
      output << "obj_var_ = [obj_var_; " << id1 << "; " << id2 << "];\n";
940
941
942
    }
}

sebastien's avatar
sebastien committed
943
DynaSaveStatement::DynaSaveStatement(const SymbolList &symbol_list_arg,
sebastien's avatar
trunk:    
sebastien committed
944
                                     const string &filename_arg) :
sebastien's avatar
sebastien committed
945
  symbol_list(symbol_list_arg),
sebastien's avatar
trunk:    
sebastien committed
946
  filename(filename_arg)
947
948
949
950
951
952
{
}

void
DynaSaveStatement::writeOutput(ostream &output, const string &basename) const
{
sebastien's avatar
sebastien committed
953
  symbol_list.writeOutput("var_list_", output);
sebastien's avatar
trunk:    
sebastien committed
954
955
  output << "dynasave('" << filename
         << "',var_list_);" << endl;
956
957
}

sebastien's avatar
sebastien committed
958
DynaTypeStatement::DynaTypeStatement(const SymbolList &symbol_list_arg,
sebastien's avatar
trunk:    
sebastien committed
959
                                     const string &filename_arg) :
sebastien's avatar
sebastien committed
960
  symbol_list(symbol_list_arg),
sebastien's avatar
trunk:    
sebastien committed
961
  filename(filename_arg)
962
963
964
965
966
967
{
}

void
DynaTypeStatement::writeOutput(ostream &output, const string &basename) const
{
sebastien's avatar
sebastien committed
968
  symbol_list.writeOutput("var_list_", output);
sebastien's avatar
trunk:    
sebastien committed
969
970
  output << "dynatype('" << filename
         << "',var_list_);" << endl;
971
972
}

973
ModelComparisonStatement::ModelComparisonStatement(const filename_list_t &filename_list_arg,
974
975
976
977
978
979
980
981
982
983
984
                                                   const OptionsList &options_list_arg) :
  filename_list(filename_list_arg),
  options_list(options_list_arg)
{
}

void
ModelComparisonStatement::writeOutput(ostream &output, const string &basename) const
{
  options_list.writeOutput(output);

sebastien's avatar
trunk:    
sebastien committed
985
986
  output << "ModelNames_ = {};" << endl;
  output << "ModelPriors_ = [];" << endl;
987

988
  for (filename_list_t::const_iterator it = filename_list.begin();
989
       it != filename_list.end(); it++)
990
    {
sebastien's avatar
trunk:    
sebastien committed
991
992
      output << "ModelNames_ = { ModelNames_{:} '" << (*it).first << "'};" << endl;
      output << "ModelPriors_ = [ ModelPriors_ ; " << (*it).second << "];" << endl;
993
    }
sebastien's avatar
trunk:    
sebastien committed
994
  output << "model_comparison(ModelNames_,ModelPriors_,oo_,options_,M_.fname);" << endl;
995
996
}

sebastien's avatar
sebastien committed
997
PlannerObjectiveStatement::PlannerObjectiveStatement(StaticModel *model_tree_arg) :
998
999
1000
  model_tree(model_tree_arg)
{
}