Shocks.cc 7.94 KB
Newer Older
sebastien's avatar
sebastien committed
1
/*
sebastien's avatar
sebastien committed
2
 * Copyright (C) 2003-2009 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/>.
 */

20
using namespace std;
sebastien's avatar
sebastien committed
21

22
23
#include <cassert>
#include <cstdlib>
sebastien's avatar
sebastien committed
24
25
#include <iostream>

26
27
#include "Shocks.hh"

sebastien's avatar
sebastien committed
28
AbstractShocksStatement::AbstractShocksStatement(bool mshocks_arg,
29
30
                                                 const det_shocks_type &det_shocks_arg,
                                                 const SymbolTable &symbol_table_arg) :
sebastien's avatar
sebastien committed
31
32
33
  mshocks(mshocks_arg),
  det_shocks(det_shocks_arg),
  symbol_table(symbol_table_arg)
34
35
36
{
}

sebastien's avatar
sebastien committed
37
38
void
AbstractShocksStatement::writeDetShocks(ostream &output) const
39
{
sebastien's avatar
sebastien committed
40
  int exo_det_length = 0;
41

sebastien's avatar
sebastien committed
42
43
44
  for(det_shocks_type::const_iterator it = det_shocks.begin();
      it != det_shocks.end(); it++)
    {
sebastien's avatar
sebastien committed
45
      int id = symbol_table.getTypeSpecificID(it->first) + 1;
sebastien's avatar
sebastien committed
46
47
      bool exo_det = (symbol_table.getType(it->first) == eExogenousDet);
      int set_shocks_index = ((int) mshocks) + 2 * ((int) exo_det);
48

sebastien's avatar
sebastien committed
49
50
51
52
      for (unsigned int i = 0; i < it->second.size(); i++)
        {
          const int &period1 = it->second[i].period1;
          const int &period2 = it->second[i].period2;
sebastien's avatar
sebastien committed
53
          const NodeID value = it->second[i].value;
sebastien's avatar
sebastien committed
54
55

          if (period1 == period2)
sebastien's avatar
sebastien committed
56
57
58
59
60
61
            {
              output << "set_shocks(" << set_shocks_index << "," << period1
                     << ", " << id << ", ";
              value->writeOutput(output);
              output << ");" << endl;
            }
sebastien's avatar
sebastien committed
62
          else
sebastien's avatar
sebastien committed
63
64
65
66
67
68
            {
              output << "set_shocks(" << set_shocks_index << "," << period1
                     << ":" << period2 << ", " << id << ", ";
              value->writeOutput(output);
              output << ");" << endl;
            }
sebastien's avatar
sebastien committed
69
70
71
72
73
74

          if (exo_det && (period2 > exo_det_length))
            exo_det_length = period2;
        }
    }
  output << "M_.exo_det_length = " << exo_det_length << ";\n";
75
76
}

sebastien's avatar
sebastien committed
77
void
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
AbstractShocksStatement::checkPass(ModFileStructure &mod_file_struct)
{
  mod_file_struct.shocks_present = true;
}


ShocksStatement::ShocksStatement(const det_shocks_type &det_shocks_arg,
                                 const var_and_std_shocks_type &var_shocks_arg,
                                 const var_and_std_shocks_type &std_shocks_arg,
                                 const covar_and_corr_shocks_type &covar_shocks_arg,
                                 const covar_and_corr_shocks_type &corr_shocks_arg,
                                 const SymbolTable &symbol_table_arg) :
  AbstractShocksStatement(false, det_shocks_arg, symbol_table_arg),
  var_shocks(var_shocks_arg),
  std_shocks(std_shocks_arg),
  covar_shocks(covar_shocks_arg),
  corr_shocks(corr_shocks_arg)
{
}

void
ShocksStatement::writeOutput(ostream &output, const string &basename) const
{
  output << "%" << endl
         << "% SHOCKS instructions" << endl
         << "%" << endl;

  // Write instruction that initializes a shock
  output << "make_ex_;" << endl;

  writeDetShocks(output);
  writeVarAndStdShocks(output);
  writeCovarAndCorrShocks(output);
  if (covar_shocks.size()+corr_shocks.size() > 0)
    output << "M_.sigma_e_is_diagonal = 0;" << endl;
  else
    output << "M_.sigma_e_is_diagonal = 1;" << endl;
}

void
ShocksStatement::writeVarAndStdShocks(ostream &output) const
119
{
sebastien's avatar
sebastien committed
120
  var_and_std_shocks_type::const_iterator it;
121

sebastien's avatar
sebastien committed
122
  for(it = var_shocks.begin(); it != var_shocks.end(); it++)
123
    {
sebastien's avatar
sebastien committed
124
      int id = symbol_table.getTypeSpecificID(it->first) + 1;
sebastien's avatar
sebastien committed
125
126
127
128
      const NodeID value = it->second;
      output << "M_.Sigma_e(" << id << ", " << id << ") = ";
      value->writeOutput(output);
      output << ";" << endl;
129
    }
sebastien's avatar
sebastien committed
130
131

  for(it = std_shocks.begin(); it != std_shocks.end(); it++)
132
    {
sebastien's avatar
sebastien committed
133
      int id = symbol_table.getTypeSpecificID(it->first) + 1;
sebastien's avatar
sebastien committed
134
135
136
137
      const NodeID value = it->second;
      output << "M_.Sigma_e(" << id << ", " << id << ") = (";
      value->writeOutput(output);
      output << ")^2;" << endl;
138
139
140
    }
}

sebastien's avatar
sebastien committed
141
void
142
ShocksStatement::writeCovarAndCorrShocks(ostream &output) const
143
{
sebastien's avatar
sebastien committed
144
145
146
  covar_and_corr_shocks_type::const_iterator it;

  for(it = covar_shocks.begin(); it != covar_shocks.end(); it++)
147
    {
sebastien's avatar
sebastien committed
148
149
      int id1 = symbol_table.getTypeSpecificID(it->first.first) + 1;
      int id2 = symbol_table.getTypeSpecificID(it->first.second) + 1;
sebastien's avatar
sebastien committed
150
151
152
153
      const NodeID value = it->second;
      output << "M_.Sigma_e(" << id1 << ", " << id2 << ") = ";
      value->writeOutput(output);
      output << "; M_.Sigma_e(" << id2 << ", " << id1 << ") = M_.Sigma_e("
sebastien's avatar
sebastien committed
154
             << id1 << ", " << id2 << ");\n";
155
    }
sebastien's avatar
sebastien committed
156
157

  for(it = corr_shocks.begin(); it != corr_shocks.end(); it++)
158
    {
sebastien's avatar
sebastien committed
159
160
      int id1 = symbol_table.getTypeSpecificID(it->first.first) + 1;
      int id2 = symbol_table.getTypeSpecificID(it->first.second) + 1;
sebastien's avatar
sebastien committed
161
162
163
164
      const NodeID value = it->second;
      output << "M_.Sigma_e(" << id1 << ", " << id2 << ") = ";
      value->writeOutput(output);
      output << "*sqrt(M_.Sigma_e(" << id1 << ", " << id1 << ")*M_.Sigma_e("
165
             << id2 << ", " << id2 << ")); M_.Sigma_e(" << id2 << ", "
sebastien's avatar
sebastien committed
166
             << id1 << ") = M_.Sigma_e(" << id1 << ", " << id2 << ");\n";
167
168
169
    }
}

sebastien's avatar
sebastien committed
170
171
MShocksStatement::MShocksStatement(const det_shocks_type &det_shocks_arg,
                                   const SymbolTable &symbol_table_arg) :
172
  AbstractShocksStatement(true, det_shocks_arg, symbol_table_arg)
173
174
175
{
}

sebastien's avatar
sebastien committed
176
void
sebastien's avatar
sebastien committed
177
MShocksStatement::writeOutput(ostream &output, const string &basename) const
178
{
sebastien's avatar
sebastien committed
179
  output << "%" << endl
180
         << "% MSHOCKS instructions" << endl
sebastien's avatar
sebastien committed
181
         << "%" << endl;
182

sebastien's avatar
sebastien committed
183
  // Write instruction that initializes a shock
sebastien's avatar
sebastien committed
184
  output << "make_ex_;" << endl;
185

sebastien's avatar
sebastien committed
186
  writeDetShocks(output);
187
}
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240

ConditionalForecastPathsStatement::ConditionalForecastPathsStatement(const AbstractShocksStatement::det_shocks_type &paths_arg, const SymbolTable &symbol_table_arg) :
  paths(paths_arg),
  symbol_table(symbol_table_arg),
  path_length(-1)
{
}

void
ConditionalForecastPathsStatement::checkPass(ModFileStructure &mod_file_struct)
{
  for(AbstractShocksStatement::det_shocks_type::const_iterator it = paths.begin();
      it != paths.end(); it++)
    {
      int this_path_length = 0;
      const vector<AbstractShocksStatement::DetShockElement> &elems = it->second;
      for(int i = 0; i < (int) elems.size(); i++)
        // Period1 < Period2, as enforced in ParsingDriver::add_period()
        this_path_length = max(this_path_length, elems[i].period2);
      if (path_length == -1)
        path_length = this_path_length;
      else if (path_length != this_path_length)
        {
          cerr << "conditional_forecast_paths: all constrained paths must have the same length!" << endl;
          exit(EXIT_FAILURE);
        }
    }
}

void
ConditionalForecastPathsStatement::writeOutput(ostream &output, const string &basename) const
{
  assert(path_length > 0);
  output << "constrained_vars_ = [];" << endl
         << "constrained_paths_ = zeros(" << paths.size() << ", " << path_length << ");" << endl;

  int k = 1;

  for(AbstractShocksStatement::det_shocks_type::const_iterator it = paths.begin();
      it != paths.end(); it++)
    {
      output << "constrained_vars_ = strvcat(constrained_vars_, '" << it->first << "');" << endl;
      const vector<AbstractShocksStatement::DetShockElement> &elems = it->second;
      for(int i = 0; i < (int) elems.size(); i++)
        for(int j = elems[i].period1; j <= elems[i].period2; j++)
          {
            output << "constrained_paths_(" << k << "," << j << ")=";
            elems[i].value->writeOutput(output);
            output << ";" << endl;
          }
      k++;
    }
}