dynare/dynare++/src/main.cc

// Copyright © 2004-2011, Ondra Kamenik

#include "dynare3.hh"
#include "dynare_exception.hh"
#include "dynare_params.hh"

#include "utils/cc/exception.hh"
#include "parser/cc/parser_exception.hh"
#include "../sylv/cc/SylvException.hh"
#include "../kord/seed_generator.hh"
#include "../kord/global_check.hh"
#include "../kord/approximation.hh"

int
main(int argc, char **argv)
{
  DynareParams params(argc, argv);
  if (params.help)
    {
      params.printHelp();
      return 0;
    }
  if (params.version)
    {
      printf(u8"Dynare++ v. %s. Copyright © 2004-2011, Ondra Kamenik\n",
             DYNVERSION);
      printf("Dynare++ comes with ABSOLUTELY NO WARRANTY and is distributed under\n");
      printf("GPL: modules integ, tl, kord, sylv, src, extern and documentation\n");
      printf("LGPL: modules parser, utils\n");
      printf(" for GPL  see http://www.gnu.org/licenses/gpl.html\n");
      printf(" for LGPL see http://www.gnu.org/licenses/lgpl.html\n");
      return 0;
    }
  sthread::detach_thread_group::max_parallel_threads = params.num_threads;

  try
    {
      // make journal name and journal
      std::string jname(params.basename);
      jname += ".jnl";
      Journal journal(jname.c_str());

      // make dynare object
      Dynare dynare(params.modname, params.order, params.ss_tol, journal);
      // make list of shocks for which we will do IRFs
      std::vector<int> irf_list_ind;
      if (params.do_irfs_all)
        for (int i = 0; i < dynare.nexog(); i++)
          irf_list_ind.push_back(i);
      else
        irf_list_ind = ((const DynareNameList &) dynare.getExogNames()).selectIndices(params.irf_list);

      // write matlab files
      FILE *mfd;
      std::string mfile1(params.basename);
      mfile1 += "_f.m";
      if (nullptr == (mfd = fopen(mfile1.c_str(), "w")))
        {
          fprintf(stderr, "Couldn't open %s for writing.\n", mfile1.c_str());
          exit(1);
        }
      ogdyn::MatlabSSWriter writer0(dynare.getModel(), params.basename.c_str());
      writer0.write_der0(mfd);
      fclose(mfd);

      std::string mfile2(params.basename);
      mfile2 += "_ff.m";
      if (nullptr == (mfd = fopen(mfile2.c_str(), "w")))
        {
          fprintf(stderr, "Couldn't open %s for writing.\n", mfile2.c_str());
          exit(1);
        }
      ogdyn::MatlabSSWriter writer1(dynare.getModel(), params.basename.c_str());
      writer1.write_der1(mfd);
      fclose(mfd);

      // open mat file
      std::string matfile(params.basename);
      matfile += ".mat";
      mat_t *matfd = Mat_Create(matfile.c_str(), nullptr);
      if (matfd == nullptr)
        {
          fprintf(stderr, "Couldn't open %s for writing.\n", matfile.c_str());
          exit(1);
        }

      // write info about the model (dimensions and variables)
      dynare.writeMat(matfd, params.prefix);
      // write the dump file corresponding to the input
      dynare.writeDump(params.basename);

      seed_generator::set_meta_seed(static_cast<std::mt19937::result_type>(params.seed));

      TLStatic::init(dynare.order(),
                     dynare.nstat()+2*dynare.npred()+3*dynare.nboth()
                     +2*dynare.nforw()+dynare.nexog());

      Approximation app(dynare, journal, params.num_steps, params.do_centralize, params.qz_criterium);
      try
        {
          app.walkStochSteady();
        }
      catch (const KordException &e)
        {
          // tell about the exception and continue
          printf("Caught (not yet fatal) Kord exception: ");
          e.print();
          JournalRecord rec(journal);
          rec << "Solution routine not finished (" << e.get_message()
              << "), see what happens" << endrec;
        }

      std::string ss_matrix_name(params.prefix);
      ss_matrix_name += "_steady_states";
      ConstTwoDMatrix(app.getSS()).writeMat(matfd, ss_matrix_name.c_str());

      // check the approximation
      if (params.check_along_path || params.check_along_shocks
          || params.check_on_ellipse)
        {
          GlobalChecker gcheck(app, sthread::detach_thread_group::max_parallel_threads, journal);
          if (params.check_along_shocks)
            gcheck.checkAlongShocksAndSave(matfd, params.prefix,
                                           params.getCheckShockPoints(),
                                           params.getCheckShockScale(),
                                           params.check_evals);
          if (params.check_on_ellipse)
            gcheck.checkOnEllipseAndSave(matfd, params.prefix,
                                         params.getCheckEllipsePoints(),
                                         params.getCheckEllipseScale(),
                                         params.check_evals);
          if (params.check_along_path)
            gcheck.checkAlongSimulationAndSave(matfd, params.prefix,
                                               params.getCheckPathPoints(),
                                               params.check_evals);
        }

      // write the folded decision rule to the Mat-4 file
      app.getFoldDecisionRule().writeMat(matfd, params.prefix);

      // simulate conditional
      if (params.num_condper > 0 && params.num_condsim > 0)
        {
          SimResultsDynamicStats rescond(dynare.numeq(), params.num_condper, 0);
          Vector det_ss{app.getSS().getCol(0)};
          rescond.simulate(params.num_condsim, app.getFoldDecisionRule(), det_ss, dynare.getVcov(), journal);
          rescond.writeMat(matfd, params.prefix);
        }

      // simulate unconditional
      //const DecisionRule& dr = app.getUnfoldDecisionRule();
      const DecisionRule &dr = app.getFoldDecisionRule();
      if (params.num_per > 0 && params.num_sim > 0)
        {
          SimResultsStats res(dynare.numeq(), params.num_per, params.num_burn);
          res.simulate(params.num_sim, dr, dynare.getSteady(), dynare.getVcov(), journal);
          res.writeMat(matfd, params.prefix);

          // impulse response functions
          if (!irf_list_ind.empty())
            {
              IRFResults irf(dynare, dr, res, irf_list_ind, journal);
              irf.writeMat(matfd, params.prefix);
            }
        }

      // simulate with real-time statistics
      if (params.num_rtper > 0 && params.num_rtsim > 0)
        {
          RTSimResultsStats rtres(dynare.numeq(), params.num_rtper, params.num_burn);
          rtres.simulate(params.num_rtsim, dr, dynare.getSteady(), dynare.getVcov(), journal);
          rtres.writeMat(matfd, params.prefix);
        }

      Mat_Close(matfd);

    }
  catch (const KordException &e)
    {
      printf("Caught Kord exception: ");
      e.print();
      return e.code();
    }
  catch (const TLException &e)
    {
      printf("Caught TL exception: ");
      e.print();
      return 255;
    }
  catch (SylvException &e)
    {
      printf("Caught Sylv exception: ");
      e.printMessage();
      return 255;
    }
  catch (const DynareException &e)
    {
      printf("Caught Dynare exception: %s\n", e.message().c_str());
      return 255;
    }
  catch (const ogu::Exception &e)
    {
      printf("Caught ogu::Exception: ");
      e.print();
      return 255;
    }
  catch (const ogp::ParserException &e)
    {
      printf("Caught parser exception: %s\n", e.message());
      return 255;
    }

  return 0;
}